Estimating Soil Water Retention Curve by Inverse Modelling from Combination of In Situ Dynamic Soil Water Content and Soil Potential Data

Soil water retention curves (SWRCs) are crucial for characterizing soil moisture dynamics, and are particularly relevant in the context of irrigation management. Inverse modelling is one of the methods used to parameterize models representing these curves, which are closest to the field reality. The objective of this study is to estimate the soil hydraulic properties through inverse modelling using the HYDRUS-1D code based on soil moisture and potential data acquired in the field. The in situ SWRCs acquired every 30 min are based on simultaneous soil water content and soil water potential measurements with 10HS and MPS-2 sensors, respectively, in five experimental fields. The fields were planted with drip-irrigated lettuces from February to March 2016 in the Chrey Bak catchment located in the Tonlé Sap Lake region, Cambodia. After calibration of the van Genuchten soil water retention model parameters, we used them to evaluate the performance of HYDRUS-1D to predict soil moisture dynamics in the studied fields. Water flow was reasonably well reproduced in all sites covering a range of soil types (loamy sand and loamy soil) with root mean square errors ranging from 0.02 to 0.03 cm3 cm−3.


Introduction
In agriculture, soil water in unsaturated porous soil media is crucial for crop development [1].Adequate characterization of soil water movement can improve agricultural water management fundamentally and economically [2][3][4].The soil water retention curve (SWRC) is considered to be a paramount and a priori property of the hydraulic behavior of soils [2,[5][6][7][8][9].SWRCs are essentially required to estimate the soil water availability for plant use, and to simulate water and solute flow in vadose zones [10].A SWRC is defined as the relationship between soil water matric potential (h) and volumetric soil water content (θ) [8].
SWRC is considered to be a difficult-to-measure soil hydraulic property [11].There are many direct and indirect methods to characterize the SWRC [12].Direct measurement can be conducted in the laboratory with a small soil sample or in a small-scale in situ experiment using some devices [5].Generally, SWRC is obtained directly by laboratory experiments using porous media-based methods such as the sand box, hanging water column, pressure cells, pressure plate extractors, and centrifuge [13].Tension disc infiltrometer experiments are commonly used for SWRC field experiments [12,[14][15][16][17][18][19].Soil moisture and soil water potential probes can also be used to obtain in situ SWRCs [20].Both laboratory and filed experiments have their own advantages and limitations.
The field experiments can minimize the disturbance caused by soil sampling and unchanged soil boundary conditions, whereas laboratory experiments can measure SWRCs in a larger range, especially in the wettest and the driest events [21].However, these direct measurement methods are time consuming and expensive [9].
Because of the disadvantages of direct measurement, alternative indirect methods using computer modelling have been developed for fast and accurate prediction [22].A rapid approach is the use of pedotransfer functions (PTFs) based on easy-to-measure soil data such as soil texture, bulk density, and soil organicity [23,24].However, the reliability of applying these relationships in other contexts is uncertain and requires cautious validation for different regions [25].
In recent decades, an inverse solution has appeared as an attractive procedure to obtain SWRCs [26,27].Inversion estimation can be an easy and reliable procedure [28].Several existing simulation model software tools include inverse estimation as embedded functions within the program [29].This approach involves estimating a set of limited unknown model parameters by using easily measurable variables (model output) such as water flow to compare to observation through the process of objective function optimization [9,26].The soil hydraulic functions for inversion are described principally by analytical functions [30] such as Gardner [31], Durner [32], Campbell [33], and van Genuchten [34].van Genuchten has gained wide recognition in generating a reasonable SWRC from various laboratory and field experiments [35].However, this model may not properly describe the hydrodynamic functioning of clay soils [36].Obtaining accurately the important parameters of the function is still a great challenge [37].The limitations of inverse modelling are mostly related to the solution uniqueness, and insufficient data due to limited measurement range of instruments used [9,22].Minimization of the these problems can be obtained with the following precautions: (i) it should have inverse input data for objective function with a wide range of water content, soil pressure head and additional retention curve data from simultaneously measuring pressure head and water content data in the soil profile [38][39][40][41]; (ii) initial soil parameter values should be reasonably close to their true values [38].Some studies proved that including water content or infiltration rate alone will not provide uniqueness of the optimized parameters in the inverse solution [40,41].Schelle et al. [42] recommended that instrumentation of a lysimeter equipped with pressure head sensors can significantly improve parameter identifiability of the inverse solution.
The inverse solution has been successfully applied to laboratory experiments with quick and precise results [5,[41][42][43].A recent laboratory method for SWRC determination using transient water release and imbibitions (TWRI) method was developed by Wayllace et al. [44].This transient method involves using physical tests using electrical balance, and numerical tests using inverse modelling with HYDRUS and the van Genuchten model [45].The procedure presents a fast, accurate, and simple testing tool for obtaining SWRC and hydraulic conductivity functions of various types of soils under drying and wetting conditions.To date, this approach has been applied in geotechnical engineering with few experimental results establishing its validity and generality [45][46][47][48][49].
In the field conditions, numerous studies of the inverse estimation technique have been applied [23,41,50].Notably, the HYDRUS model has increasingly and successfully been used for inverse estimation particularly in in situ conditions [5,18,27,51,52].Lai and Ren [4] determined soil hydraulic parameters at field scale by inverse modelling using combined HYDRUS-1D and PEST models.The authors confirmed that there are no unique effective average properties for a heterogeneous field to simulate field water content.Their results showed that the inverse modelling approach simulated soil water dynamics well, but still needs to be improved for layered soils, especially with fine-textured soil layer.Inverse modelling with the HYDRUS model was successfully performed by Filipović et al. [52] who used infiltration experiments with tension disc infiltration data to estimate soil hydraulic properties.Similarly, Rashid et al. [12] conclude that the HYDRUS inverse solution approach applied to infiltration data measured with a tension disc infiltrometer is a useful method to characterize soil hydraulic properties.Le Bourgeois et al. [5] investigated an inverse modelling method using HYDRUS-1D and the NSGA-II algorithm to estimate soil hydraulic properties from in situ water content measurement.Their calibrated Mualem-van Genuchten parameters had very low uncertainty and could still result in a good simulation of water flow.
Previous studies on earlier inverse modelling studies mostly using HYDRUS focused on identifying soil hydraulic properties using infiltration experiments based on using tension disc infiltrometer data as discussed above.Thus, it is more interesting to inverse the soil hydraulic properties using field retention curve data with a fine timestep of dynamic soil moisture content and pressure head measured simultaneously.
The objectives of this study are twofold: (i) to estimate soil properties for retention curve using inverse modelling with HYDRUS-1D with input field data of dynamic water content and SWRCs; and (ii) to evaluate its ability to simulate the water flow process in a case study of Cambodian soils.

Experimental Sites
The five experimental sites (T1, T2, T3, T4 and T5) are in the Chrey Bak catchment, known as an important agricultural development area situated in Kampong Chhnang Province in north-western Cambodia (Figure 1) [53].The study area is influenced by the tropical monsoon climate, having two seasons, the dry season from December to May and the rainy season from June to November.Climate data was collected from a weather station downstream of the catchment.Data recorded from 2012-2014 showed a mean annual temperature of 27 • C, and a mean annual minimum and maximum humidity ranging from 51 to 92 per cent, respectively.Average solar radiation was from 16.8 to 18 mega joules per square meter per day.Mean annual precipitation ranged from 1600 mm to 1788 mm.situ water content measurement.Their calibrated Mualem-van Genuchten parameters had very low uncertainty and could still result in a good simulation of water flow.
Previous studies on earlier inverse modelling studies mostly using HYDRUS focused on identifying soil hydraulic properties using infiltration experiments based on using tension disc infiltrometer data as discussed above.Thus, it is more interesting to inverse the soil hydraulic properties using field retention curve data with a fine timestep of dynamic soil moisture content and pressure head measured simultaneously.
The objectives of this study are twofold: (i) to estimate soil properties for retention curve using inverse modelling with HYDRUS-1D with input field data of dynamic water content and SWRCs; and (ii) to evaluate its ability to simulate the water flow process in a case study of Cambodian soils.

Experimental Sites
The five experimental sites (T1, T2, T3, T4 and T5) are in the Chrey Bak catchment, known as an important agricultural development area situated in Kampong Chhnang Province in north-western Cambodia (Figure 1) [53].The study area is influenced by the tropical monsoon climate, having two seasons, the dry season from December to May and the rainy season from June to November.Climate data was collected from a weather station downstream of the catchment.Data recorded from 2012-2014 showed a mean annual temperature of 27 °C, and a mean annual minimum and maximum humidity ranging from 51 to 92 per cent, respectively.Average solar radiation was from 16.8 to 18 mega joules per square meter per day.Mean annual precipitation ranged from 1600 mm to 1788 mm.

Experimental Set-Up
The experimental sites 400 m 2 in area had lettuce planted with a density of 16 plants m −2 .Transplantation of the lettuce seedling and harvesting was on 30 January and 4 March 2016, respectively.Irrigation was applied using the drip irrigation system.The drip line has a 30 cm emitter space with maximum manufactory capacity of 3 L h −1 .The soil bed was raised to 0.3 m in height with the top width of 1 m.There was no rain during the experiment.

Soil Measurement
The soil characteristics of the sites are described in Table 1 using the United Nations' Food and Agriculture Organization (FAO) soil description guide [54].The soils of each field were sampled to determine the soil texture and bulk density with three replications at the laboratory.Table 2

Experimental Set-Up
The experimental sites 400 m 2 in area had lettuce planted with a density of 16 plants m −2 .Transplantation of the lettuce seedling and harvesting was on 30 January and 4 March 2016, respectively.Irrigation was applied using the drip irrigation system.The drip line has a 30 cm emitter space with maximum manufactory capacity of 3 L h −1 .The soil bed was raised to 0.3 m in height with the top width of 1 m.There was no rain during the experiment.

Soil Measurement
The soil characteristics of the sites are described in Table 1 using the United Nations' Food and Agriculture Organization (FAO) soil description guide [54].The soils of each field were sampled to determine the soil texture and bulk density with three replications at the laboratory.Table 2 presents the basic soil physical properties.Soil texture and bulk density measurement methods have been described in a previous study of Ket et al. [55].Porosity (∅) and void ratio (e) were calculated from bulk density (Bd) assuming an average mineral density of 2.65 g cm −3 [56].pH of the water was measured by potentiometry [57].
In situ SWRC data collection was conducted in the five experimental fields.The SWRCs were measured using coupled inexpensive Decagon Devices sensors, e.g., soil moisture sensor 10HS (Decagon Devices, Pullman, WA, USA) and a soil water matric potential sensor MPS-2 (Decagon Devices, Pullman, WA, USA) connected to Em50 data logger (Decagon Devices, Pullman, WA, USA) [59].The sensors were installed horizontally near each other at depths of 10 cm and 20 cm below the soil surface between irrigated lettuces at sites T1, T2, and T3 and below the bare soil at sites T4 and T5.The sensor record was from 1 to 21 March 2016 at site T1, from 22 February to 21 March at site T2, from 14 February to 21 March at site T3, from 22 February to 21 March at sites T4 and from 7 to 21 March at site T5.Data from the sensors were recorded at 30 min intervals.The sensors of 10HS and MPS-2 for SWRCs measurement were available only at 20 cm depth at T1 and 10 cm depth at T5. ), ∅: porosity (-), e: void ratio (-), θ ini : initial water content (cm 3 cm −3 ), Ks: saturated hydraulic conductivity (mm s −1 ), θ s : saturated soil moisture (cm 3 cm −3 ).

Instrumentation
MPS-2 is the dielectric water matric potential sensor [60].MPS-2 measures the water content of porous ceramic discs using a capacitive reading and converts the measured water content to water potential using the moisture characteristic curve of the ceramic.The range of the measurement is from −10 to −500 kPa (pF 2.01 to pF 3.71).Sensor accuracy is ±25% of the reading between −9 and −100 kPa.MPS-2 has been confirmed to have good reliability in its respective range [20].Accuracy can decrease up to approximately ±35% and ±50% at −300 kPa and −500 kPa respectively [61].The MPS-2 pre-calibrated by the manufacturer is not affected by soil type and requires correct installation with adequate hydraulic contact.10HS is a capacitance sensor related to dielectric permittivity [62].10HS has a measurement range of 0-0.57cm 3 cm −3 and an accuracy of ±0.03 cm 3 cm −3 .It operates at a frequency 70 MHz and it can be used at temperatures between 0 and +50 • C with permittivity measurement volume of 1 dm 3 [63].Mittelbach et al. [64] mentioned that the 10HS sensor fails to measure soil moisture (θ) above 0.4 cm 3 cm −3 and presents a decreasing sensitivity in measuring θ with increasing θ, resulting in a poor ability to represent the variability in θ for moist conditions.The 10HS sensors in this study were calibrated following Decagon step-by-step instructions [65] with the calibrated equation in Table 3.

5-Min Timestep Reference Evapotranspiration Computation
The reference evapotranspiration (ETo) was calculated based on standardized American Society of Civil Engineers (ASCE) Penman-Monteith equation (ASCE-PM) [66,67] using meteorological data collected in 5-min timesteps (e.g., precipitation, air temperature, solar radiation, relative humidity, and wind speed) from the weather station (Figure 1).The standardized ASCE-PM equation is: where ET o is the standardized grass-reference evapotranspiration( ET) (mm 5-min −1 ), ∆ is the slope of saturation vapour pressure versus air temperature curve (kPa • C −1 ), Rn is the calculated net radiation at the crop surface (Mega Joule (MJ) m −2 5-min −1 ), G is the heat flux density at the soil surface (MJ m −2 5-min −1 ), T is the air temperature at 1.5 to 2.5 m height ( • C), U 2 is the wind speed at 2-m height, e s is the saturation vapour pressure (kPa), e a is the actual vapour pressure (kPa), γ is the psychrometric constant (kPa • C −1 ).C n is the numerator constant that changes with reference surface and calculation timestep (900 • C mm s 3 Mg −1 day −1 for 24 h timesteps, and 37 • C mm s 3 Mg −1 h −1 for hourly timesteps for the grass-reference surface, so 3.083 • C mm s 3 Mg −1 h −1 for 5-min timestep).Cd is the denominator constant that changes with reference surface and calculation timestep (0.34 s m −1 for 24 h timesteps, 0.24 s m −1 for hourly timesteps during daytime, and 0.96 s m −1 for hourly night time for hourly or shorter time steps for the grass-reference surface) [68].The values for Cn and Cd are associated with bulk surface resistance (r s ) and aerodynamic roughness of the surface (r a ) [66,69].The values r s-daytime = 50 s m −1 and r s-nighttime = 200 s m −1 were adapted in this study.The 5-min step G was a function of R n (G 5min-day time = 0.1 R n and G 5min-night time = 0.5 R n ).

Partitioning Crop Evapotranspiration
When the soil has full canopy cover, crop evapotranspiration (ET c ) is often assumed to be similar to transpiration [70].The soil not shadowed by crops is exposed to radiation, which is considered to be soil evaporation [71].Therefore, crop transpiration is closely related to canopy cover.Partitioning of crop evapotranspiration (ET c , mm min −1 ) into soil evaporation (E, mm min −1 ) and crop transpiration (T, mm min −1 ) was based on the method proposed and described by Gallardo et al. [72].Crop transpiration (T) was estimated according to reference evapotranspiration (ET o ) and ground canopy cover (G) using the following equation.
where T is the crop transpiration (mm min −1 ), ET o is the reference evapotranspiration (mm min −1 ).K c max is the maximum crop coefficient (K c ) value.When ground cover is less than 10%, a K c of lettuce is about 1.05 if it is well irrigated [73].When a canopy cover reaching about 75%, lettuce has a K c > 1.05 and K c max = 1.10 was adapted for lettuce during the mid-season.G is the ground cover percentage, ]. G x is maximum ground canopy cover (%).The coefficients a = 6.58 and b = −10.02were adapted in this study.N is normalized accumulative reference evapotranspiration.Soil evaporation (E) was partitioned with below equation.
At site T4 and T5, sensors were placed in the soil profile where plants did not grow.Therefore, there was no partitioning ET c at these sites.The crop input data for partitioning are shown in Table 4.The maximum crop canopy covers were measured using a smartphone camera at 1 m height at harvest time and converted into canopy cover percentage using the Canopeo smartphone app.
where θ is the volumetric water content (L 3 L −3 ), t is the time (T), z is the positive downward vertical space coordinate (L), K(h) is the unsaturated hydraulic conductivity function (LT −1 ), h is the pressure head (L).
where S e is the effective water content (-), θ r and θ s are the residual and saturated water content (L 3 L −3 ), α (>0, in L −1 ) is related to the inverse of the air-entry suction, m and n is the curve shape parameters, m = 1 − 1/n, n (>1) is a measure of the pore-size distribution, l is the pore-connectivity/tortuosity parameter (-), K s is the saturated hydraulic conductivity (L T −1 ).Of these, θ r , θ s and K s have a physical meaning, and α, n and l are empirical curve shape parameters [34,78].l was set to be 0.5 recommended by Mualem [77].The parameter n affects the steepness of the curve.Large values of n result in a steeper curve [79].The period of simulation was 21, 28, 36, 29 and 12 days at sites T1, T2, T3, T4 and T5 respectively during February and March 2016 according to the available data from the soil sensors, 10HS and MPS-2.

Air-Entry Value
The air-entry value (AEV) is a critical and commonly used variable obtained from the SWRC for estimating other unsaturated soil properties, i.e., permeability [80,81].AEV is defined as the suction at which drainage of the soil pores begins [79,82].AEV occurs usually between -h = 1 to 10 kPa [79].The AEV reflects the maximum soil pore size and the soil texture [83].A decrease in soil grain size leads to an increase in AEV and flattening of the SWRC slope [82].A coarse-grained soil has a lower air-entry value, and lower residual suction than a fine-grained soil [84].Besides soil texture, other effects listed by Wang et al. [85] are soil structure, initial water content, contact angle, organic matter, clay content, and bulk density.Gallage and Uchimura [84] found that soils with low density have lower AEV and residual suction than soils with a high density.Numerically, AEV is related to the α and n parameters of the vG model [83].Lower values of α indicate that the air-entry region is broad [79].Soltani [80] has proposed simplified methods to determine AEV as below.

Time-Variable Boundary and Initial Conditions
HYDRUS-1D input data components in this study include input data for the variable boundary condition (i.e., evapotranspiration, transpiration, and irrigation), soil hydraulic properties and inverse input data (soil moisture in time and retention curve data).The variable boundary condition of this study is illustrated in Figure 2.
The upper boundary condition was set to an atmospheric condition influenced by irrigation supply, and crop evapotranspiration with a surface layer as indicated in the following equation: where q o (t) is the difference between irrigation, transpiration, and evaporation rate.
The free/zero-gradient drainage boundary condition is suitable for water flow simulation of unsaturated soil in which the soil domain of interest is not affected by groundwater [86,87].The process allows water to leave a flow domain by gravity, assuming a unit vertical hydraulic gradient without external forced drainage conditions [87].Based on these conditions and the assumptions of no influence of the groundwater to the studied soil profile, the free drainage was considered at the bottom of the soil domain.The free drainage condition is represented as follows: Where L sp is the soil profile assumed at 40 cm depth.HYDRUS-1D can give water flow at any specific soil depth [88].The observation nodes were set at 10 and 20 cm following the location of the 10HS and MPS-2 sensors.The soil moisture at the beginning of the simulation in the soil profile domain was set to the initial observation from these 10HS sensors.

Time-Variable Boundary and Initial Conditions
HYDRUS-1D input data components in this study include input data for the variable boundary condition (i.e., evapotranspiration, transpiration, and irrigation), soil hydraulic properties and inverse input data (soil moisture in time and retention curve data).The variable boundary condition of this study is illustrated in Figure 2.
The upper boundary condition was set to an atmospheric condition influenced by irrigation supply, and crop evapotranspiration with a surface layer as indicated in the following equation: where q o (t) is the difference between irrigation, transpiration, and evaporation rate.
The free/zero-gradient drainage boundary condition is suitable for water flow simulation of unsaturated soil in which the soil domain of interest is not affected by groundwater [86,87].The process allows water to leave a flow domain by gravity, assuming a unit vertical hydraulic gradient without external forced drainage conditions [87].Based on these conditions and the assumptions of no influence of the groundwater to the studied soil profile, the free drainage was considered at the bottom of the soil domain.The free drainage condition is represented as follows: Where L  is the soil profile assumed at 40 cm depth.HYDRUS-1D can give water flow at any specific soil depth [88].The observation nodes were set at 10 and 20 cm following the location of the 10HS and MPS-2 sensors.The soil moisture at the beginning of the simulation in the soil profile domain was set to the initial observation from these 10HS sensors.The initial soil hydraulic parameters (e.g., θ r , θ s , α, n and l) were determined by the Rosetta method of Schaap et al. [89] based on soil texture, and initial K s parameters were obtained from tension infiltrometer measurement (Table 5).The initial soil hydraulic parameters (e.g., θ r , θ s , α, n and l) were determined by the Rosetta method of Schaap et al. [89] based on soil texture, and initial K s parameters were obtained from tension infiltrometer measurement (Table 5).

Inverse Solution
In HYDRUS-1D, the Marquart-Levenberge method is used to optimize the soil hydraulic parameters [23,[90][91][92].The optimization process to generate vG parameters by the model in this study is to minimize the difference between simulated and observed values of water content, θ(t), and soil water retention data, h(θ) through an objective function, Φ(θ(t), h(θ)) as described below [93].
where m is the two types of dataset, i.e., soil water content θ(t), and retention curve h(θ), n j is the number of measurements of the jth dataset, p * ij and p ij (b) are the observations and predictions for the jth measurement set, b (e.g., θ r , θ s , α, n, and K s ) is the vector of optimsed parameters.The first set of measurements is water content in time, θ(t) using 10HS measurement, and the second is the retention curve h(θ) from simultaneous measurement of 10HS soil water content and MPS-2 soil water potential.The hydraulic parameters for two soil depths at 10 and 20 cm were optimized simultaneously in the inversion process.
To evaluate the model performance, we used three fitted statistical indicators, e.g., the root mean square error (RMSE), the Nash and Sutcliffe model efficiency (NSE), and the coefficient of determination (R 2 ) with the following expressions.
where p * ij and p ij (b) are the observation and simulation as described above, n j is total amount of the data.R 2 ranges between 0 and 1. Value 1 indicates that the simulation dispersion is equal to the observation, whereas R 2 equals to 0, there is no any correlation between them [94].The lower RMSE value is, the better model performs [95].NSE values are between −∞ and 1.0 (1 perfect fit) [95].If 0 < NSE < 1.0, it is considered to be an acceptable performance, otherwise a negative NSE indicates unacceptable performance [95][96][97].NSE < 1.0, it is considered to be an acceptable performance, otherwise a negative NSE indicates unacceptable performance [95][96][97].

Partition of Crop Evapotranspiration
Results of estimation of transpiration and evaporation from partitioning crop evapotranspiration (ET  ) based on the method proposed by Gallardo et al. [72] at sites T1, T2 and T3 are presented in Figure 4.The transpiration estimation at site T1 and T2 was done during mid-season with full canopy cover of 35% and 57% respectively and during development stage at T3 reaching the full canopy cover of 44%.The results show that soil evaporation was the main part of ET  during the growing season.The soil evaporation ranged from 47% to 67% at the maximum canopy cover stage in the three sites.

Partition of Crop Evapotranspiration
Results of estimation of transpiration and evaporation from partitioning crop evapotranspiration (ET c ) based on the method proposed by Gallardo et al. [72] at sites T1, T2 and T3 are presented in Figure 4.The transpiration estimation at site T1 and T2 was done during mid-season with full canopy cover of 35% and 57% respectively and during development stage at T3 reaching the full canopy cover of 44%.The results show that soil evaporation was the main part of ET c during the growing season.The soil evaporation ranged from 47% to 67% at the maximum canopy cover stage in the three sites.

Estimated Parameters and Model Evaluation
The inverse estimation technique using HYDRUS-1D has been successfully used to determine the soil hydraulic function of van Genuchten (vG) parameters in the two soil depths of 10 cm and 20 cm in the five experimental sites.Table 6 shows the results of optimized vG parameter sets in all experimental sites.The estimated θ s ranged from 0.35 to 0.43 cm 3 cm −3 for loamy soils (T1 and T5), 0.30 to 0.35 cm 3 cm −3 for sand soils (T2 and T3) and 0.32 to 0.37 cm 3 cm −3 for loam soil (T4).It was noted that these estimations are similar to the measured values.At site T4, the estimated θ s (ranging from 0.32-0.37cm 3 cm −3 ) is much lower than the measured value (0.43 cm 3 cm −3 ).Commonly, the laboratory-measured θ s can be smaller than in situ values because of incomplete saturation and air entrapment of the soil [18,99], while it is possibly over/underestimated in the laboratory-measured θ s [100].
For the computation of saturated hydraulic conductivity parameters, it was noted that inputting low initial measured saturated conductivity in sand soil during the optimization process with K s = 0.43 mm min −1 at site T2 and K s = 1.45 mm min −1 at site T3 resulted in a divergence of the model.In that case, initial K s were estimated using Rosetta method of Schaap et al. [89] based on soil texture.

Estimated Parameters and Model Evaluation
The inverse estimation technique using HYDRUS-1D has been successfully used to determine the soil hydraulic function of van Genuchten (vG) parameters in the two soil depths of 10 cm and 20 cm in the five experimental sites.Table 6 shows the results of optimized vG parameter sets in all experimental sites.The estimated θ s ranged from 0.35 to 0.43 cm 3 cm −3 for loamy soils (T1 and T5), 0.30 to 0.35 cm 3 cm −3 for sand soils (T2 and T3) and 0.32 to 0.37 cm 3 cm −3 for loam soil (T4).It was noted that these estimations are similar to the measured values.At site T4, the estimated θ s (ranging from 0.32-0.37cm 3 cm −3 ) is much lower than the measured value (0.43 cm 3 cm −3 ).Commonly, the laboratory-measured θ s can be smaller than in situ values because of incomplete saturation and air entrapment of the soil [18,99], while it is possibly over/underestimated in the laboratory-measured θ s [100].
For the computation of saturated hydraulic conductivity parameters, it was noted that inputting low initial measured saturated conductivity in sand soil during the optimization process with K s = 0.43 mm min −1 at site T2 and K s = 1.45 mm min −1 at site T3 resulted in a divergence of the model.In that case, initial K s were estimated using Rosetta method of Schaap et al. [89] based on soil texture.Therefore, K s = 4.95 mm min −1 from Rosetta method at sites T2 and T3 was selected as the initial value and resulted in convergence.It highlights the need for reasonable input values before an inversion process.
The fitted K s in Table 6 were quite higher than the measured values from the tension infiltrometer (TI) in sites T1, T2, T3.The difference between in situ measured and simulated values of K s is most probably due to the different time and space at the measurement from the scale at which the processes are modelled [18].It is commonly known that agricultural soils frequently exhibit extensive spatial and temporal changes in pore characteristics that cause variation change in soil hydraulic conductivity, K(h).On the other hand, it was noted that K s were measured by TI at the soil surface.Thus, it caused the limitation both on K(h) extrapolation towards saturation and on its application to extended depth (10 and 20 cm) [101].Another potential reason can be due to the limitations of TI method in measurement in the coarse soil texture.The result study of March [102] suggested that TI method underestimated K s under high permeability conditions.Similarly, Yoon et al. [101] observed that there was greater fluctuation in measuring hydraulic conductivity near the saturated condition especially for sandy soil affected by the presence of irregular distribution of macropore possibly with entrapped air.The other limitations were described by Roulier et al. [103] mainly associated with the simplifying assumptions of the analysis methods and instrumental restriction [101].Indeed, K s obtained from TI was estimated indirectly and assuming only matrix flow, that might lead to inaccurate K s [104].
Therefore, the inverse model estimation suggested a reasonable value of K s using the observed soil water dynamic.There were low estimated K s of coarse soils at this soil depth at site T3 at 20 cm and T5 at 10 cm.The soil compaction at 20 cm soil depth could be the reason of the low estimated permeability at site T3.It was noted that at site T5 the second soil layer below 15 cm depth was highly compacted and that could lead to the low K s estimation for this loamy sand at 10 cm depth.In contrast, the measured K s in T4 and T5 are similar to the optimized values.This suggests that a tension infiltrometer might be useful in determining K s for inverse estimation for fine soil type.
These results also suggest that for further work, a sensitivity analysis of the vG parameters could be interesting, to investigate the constrains of its impact on the simulation.Tables 7 and 8 show the evaluation model performances of the inverse estimation to generate SWRCs and water flow, respectively.After inversion, the results showed the improvement of model performance in SWRC simulation with lower RMSE (ranging from 0.01 to 0.04 cm 3 cm −3 ) and higher coefficients of NSE (ranged from 0.53 to 0.99) and determination (R 2 ) (ranging from 0.78 to 0.99) (Table 7).Similarly, water flow simulation resulted in satisfied performance with RMSE ranged from 0.02 to 0.03 cm 3 cm −3 , NSE from 0.64 to 0.83 and R 2 from 0.85 to 0.99 (Table 8).Overall, these results confirm that the inverse modelling with HYDRUS-1D can predict soil water retention curve and dynamic water content in the soil profile with a reasonable degree of accuracy.Figure 5 shows the soil water retention curves before and after inversion.Simulated SWRCs are close to the observation in all sites and depths.The excellent fits were at site T1 at 20 cm depth (loamy sand soil) (RMSE = 0.01 cm 3 cm −3 ) and at site T2 at 20 cm depth (RMSE = 0.002 cm 3 cm −3 ).Otherwise, the other SWRCs in the other sites and depths fitted fairly well with RMSE, ranging from 0.016 to 0.048 cm 3 cm −3 .
The data of field SWRC at near-saturated and driest points in all experimental sites were not presented because of data inaccuracy and limitation from the sensors (10HS and MPS-2) at these wet and dry ranges.MPS-2 is unable to provide accurate pressures at the wet end of the SWRC because of the air-entry limit [20].The acceptable measured ranges of the sensors were selected for inverse input data as illustrated in Figure 5.Despite its limitations, MPS-2 can be useful in a drier range and to pilot irrigation [20].Similarly, the errors from the 10HS sensor can be due to its accuracy limitations in measuring the water content in moist conditions [64].On the other hand, the effect of temperature on the 10HS sensor response could have caused significant data errors in the fine soil and for high water content, especially at upper soil layers where there is higher temperature fluctuation [105].The estimated AEV of the tested soils are shown in Table 6.AEV ranged from 0.86 to 12.64 kPa.AEV of loamy sand were from 2.28 to 9.60 kPa, sand were from 0.86 to 12.64 kPa, and loam were 3.37 to 5.16 kPa.Reference [106] showed that the soil with higher sand content has smaller AEV using laboratory pressure plate extractor test.They also found the effect of initial water content on the AEV.The sand soil with higher initial water content has larger AEV (i.e., the soil with sand content of 20% to 60% resulted in AEV of 6.04 to 1.94 kPa respectively).Similarly, Konyai et al. [107], using pressure plates for defining SWRCs, found the AEV of loamy sand ranged from 1.30 to 2.00 kPa and loam soil of 0.90 kPa.However, those results were from the laboratory.It is seen that the generated SWRC at site T3 with high sand content proposed a higher AEV of 10.09 and 12.64 kPa.The physical environmental effects to the SWRC of sand soil are complex in their observations at near real AEV.The limitation of sensors to catch the accuracy of the low range of suction can be the main reason for the high estimated AEV of the sand soil at T3.In Figure 5, the generated initial SWRC from the Rosetta method are often obtained from laboratory SWRC.The inverse estimated SWRC reflected the effect of the dynamic physical environmental.Therefore, it is interesting to investigate the field SWRC with the influence of the in situ environment.

Dynamic Soil Water Content
Observed and simulated soil water content distributions in the soil profile at 10 and 20 cm at the 5 experimental sites are shown in Figures 6 and 7, respectively.The fluctuation of soil water content was coming from only irrigation applied at sites T1, T2, T3 and T4.The water flow simulation in all sites gave good and excellent correlations between observed and simulated water contents in the soil profile (i.e., R 2 values ranged from 0.85 to 0.99).It is noted that the water flow simulation was similar to the soil moisture observation in the wet range but widely over/underestimated in the end dry events after the inversion (i.e., RMSE ranged from 0.02 to 0.03 cm 3 cm −3 and NSE ranged from 0.67 to 0.84).The important factor influencing the simulation model performance is the assumption of homogeneous soil properties through the entire simulation.This assumption is not flexible to real phenomena e.g., the dynamic change of soil hydraulic properties during the growing season, the macropore flow, and lateral flux of soil moisture.Schwen et al. [2] revealed that the accuracy of the soil water flow simulation, especially of that near the surface, can be improved by using time-variable hydraulic parameters.Despite several simplifying assumptions of HYDRUS-1D model, the results of the inversion yet provided reasonable estimates of water dynamics and a better improvement from the initial simulation.

Conclusions
In this study, inversion approach using HYDRUS-1D based on only in situ measurement of soil hydraulic properties were used to estimate the van Genuchten (vG) parameters of the soil properties in experimental fields with different soil types including loamy sand, sand, and loam soil, in Chrey Bak Catchment, Cambodia.The in situ soil hydraulic measured data that were used as inverse input data included saturated hydraulic conductivity using tension infiltrometer, observed dynamic soil water content and SWRC from simultaneous measures of Decagon 10HS soil moisture and MPS-2 soil water potential sensors.To apply this approach, first we computed 5-min timestep reference

Conclusions
In this study, inversion approach using HYDRUS-1D based on only in situ measurement of soil hydraulic properties were used to estimate the van Genuchten (vG) parameters of the soil properties in experimental fields with different soil types including loamy sand, sand, and loam soil, in Chrey Bak Catchment, Cambodia.The in situ soil hydraulic measured data that were used as inverse input data included saturated hydraulic conductivity using tension infiltrometer, observed dynamic soil water content and SWRC from simultaneous measures of Decagon 10HS soil moisture and MPS-2 soil water potential sensors.To apply this approach, first we computed 5-min timestep reference evapotranspiration based on the ASCE Penman-Monteith equation and obtained reasonable results.Then we partitioned the crop evapotranspiration based on the Gallardo et al. [72] method in the field after lettuce was planted.
Analyzing the SWRC shows that the simulated SWRC using optimum vG parameters of the inversion closely and fairly matched the field SWRC data, with RMSE ranging from 0.002 to 0.048 cm 3 cm −3 .The fair match was due to the limitation in measurement of the sensors in the wet and the driest range of moisture condition in the field.In testing water flow simulation, generally the model resulted in a good simulation performance for predicting the soil water content in the experimental period in all the studied soil types, with the RMSE ranging from 0.02 to 0.03 cm 3 cm −3 .The results proposed that soil hydraulic properties can be estimated effectively and rapidly with inverse modelling using only the information from the field soil hydraulic measurement without the required expensive laboratory data.We found that a tension infiltrometer might be useful in determining saturated hydraulic conductivity (K s ) for inverse estimation for fine (loam) soil type, as the generated K s was slightly changed from the initial measured values after inversion.We confirmed that the inverse process in the model is sensitive to initial hydraulic soil parameter input, and can lead to model divergence.For further work, validation on another data set should be investigated.
It is recommended that more additional experimental data, particularly accurate field dynamic SWRC in the wet and dry encompassing new technics and methods, should be included to improve the inversion approach in further research.Additionally, the dynamic change of soil properties should be considered in the inversion model to improve the accuracy of water flow prediction.

Figure 1 .
Figure 1.Map of experimental sites.

Figure 1 .
Figure 1.Map of experimental sites.

Figure 2 .
Figure 2. Scheme of boundary condition of HYDRUS in present study.

Figure 2 .
Figure 2. Scheme of boundary condition of HYDRUS in present study.

3. 1 . 1 . 5 -Figure 3
Figure 3 illustrates the result of reference evapotranspiration (ET o ) computed based on the Penman-Monteith equation (Equation (1)) by using meteorological data collected at 5 min resolution over the growing season from February to March 2016.5-min ET o between February and March ranged from 0.000143 to 0.075 mm 5-min −1 , with a mean value of 0.017 mm 5-min −1 (Figure 3).The accumulative 5-min ET o in daily estimates ranged from 4 to 6 mm day −1 , which are reasonable values for this study area climate.Nobuhiro et al. [98] found that the average daily evapotranspiration levels during the late rainy season and the middle of the dry season in central Cambodia were 4.3 and 4.6 mm day −1 , respectively, and the maximum daily ET o levels were 5.2 and 5.7 mm day −1 respectively.

3. 1 . 1 . 5 -Figure 3
Figure 3 illustrates the result of reference evapotranspiration (ET  ) computed based on the Penman-Monteith equation (Equation (1)) by using meteorological data collected at 5 min resolution over the growing season from February to March 2016.5-min ET  between February and March ranged from 0.000143 to 0.075 mm 5-min −1 , with a mean value of 0.017 mm 5-min −1 (Figure3).The accumulative 5-min ET  in daily estimates ranged from 4 to 6 mm day −1 , which are reasonable values for this study area climate.Nobuhiro et al.[98] found that the average daily evapotranspiration levels during the late rainy season and the middle of the dry season in central Cambodia were 4.3 and 4.6 mm day −1 , respectively, and the maximum daily ET  levels were 5.2 and 5.7 mm day −1 respectively.

Figure 5 .
Figure 5. SWRC in before and after inversion solution.(a) SWRC at 20 cm depth at site T1; (b) SWRC at 10 cm depth at site T2; (c) SWRC at 20 cm depth at site T2; (d) SWRC at 10 cm depth at site T3, (e) SWRC at 20 cm depth at site T3; (f) SWRC at 10 cm depth at site T4; (g) SWRC at 20 cm depth at site T4; (h) SWRC at 10 cm depth at site T5.The measured soil water potential (h) ranged from −20 to −95 kPa at site T1, from −13 to −99 kPa and −15 to −48 kPa at soil depth of 10 and 20 cm respectively at site T2, from −14 to −476 kPa and −48 to −608 kPa at 10 and 20 cm respectively at site T3, from −11 to −69 and −10 to −111 kPa at 10 and 20 cm respectively at T4, from −48 to −212 kPa at 10 cm at site T5.

Figure 5 .
Figure 5. SWRC in before and after inversion solution.(a) SWRC at 20 cm depth at site T1; (b) SWRC at 10 cm depth at site T2; (c) SWRC at 20 cm depth at site T2; (d) SWRC at 10 cm depth at site T3, (e) SWRC at 20 cm depth at site T3; (f) SWRC at 10 cm depth at site T4; (g) SWRC at 20 cm depth at site T4; (h) SWRC at 10 cm depth at site T5.The measured soil water potential (h) ranged from −20 to −95 kPa at site T1, from −13 to −99 kPa and −15 to −48 kPa at soil depth of 10 and 20 cm respectively at site T2, from −14 to −476 kPa and −48 to −608 kPa at 10 and 20 cm respectively at site T3, from −11 to −69 and −10 to −111 kPa at 10 and 20 cm respectively at T4, from −48 to −212 kPa at 10 cm at site T5.

23 Figure 6 .
Figure 6.Observed and simulated dynamic soil moisture in 30-min timestep using inversion approach during the growing season in 2016 at different sites at 10 cm depth.(a) at site T1; (b) at site T2; (c) at site T3; (d) at site T4; (e) at site T5.

Figure 6 .
Figure 6.Observed and simulated dynamic soil moisture in 30-min timestep using inversion approach during the growing season in 2016 at different sites at 10 cm depth.(a) at site T1; (b) at site T2; (c) at site T3; (d) at site T4; (e) at site T5.

Figure 7 .
Figure 7. Observed and simulated dynamic soil moisture in 30-min timestep using inversion approach during the growing season in 2016 at different sites at 20 cm depth.(a) at site T1; (b) at site T2; (c) at site T3; (d) at site T4.

Figure 7 .
Figure 7. Observed and simulated dynamic soil moisture in 30-min timestep using inversion approach during the growing season in 2016 at different sites at 20 cm depth.(a) at site T1; (b) at site T2; (c) at site T3; (d) at site T4.

Table 1 .
Soil description of the study sites.

Table 3 .
Calibrated equation of 10HS in different sites.

Table 4 .
Crop data for partitioning crop evapotranspiration.

Table 6 .
Soil hydraulic vG parameter sets derived from inversion estimation.

Table 7 .
Evaluation of model performance in simulation of soil water retention curve (SWRC) before and after inversion.

Table 8 .
Evaluation of model performance in simulation of water flow before and after inversion.