Brooks–Corey Modeling by One-Dimensional Vertical Inﬁltration Method

: The laboratory methods used for the soil water retention curve (SWRC) construction and parameter estimation is time-consuming. A vertical inﬁltration method was proposed to estimate parameters α and n and to further construct the SWRC. In the present study, the relationships describing the cumulative inﬁltration and inﬁltration rate with the depth of the wetting front were established, and simpliﬁed expressions for estimating α and n parameters were proposed. The one-dimensional vertical inﬁltration experiments of four soils were conducted to verify if the proposed method would accurately estimate α and n . The ﬁtted values of α and n , obtained from the RETC software, were consistent with the calculated values obtained from the inﬁltration method. The comparison between the measured SWRCs obtained from the centrifuge method and the calculated SWRCs that were based on the inﬁltration method displayed small values of root mean square error ( RMSE ), mean absolute percentage error ( MAPE ), and mean absolute error. SWMS_2D-based simulations of cumulative inﬁltration, based on the calculated α and n , remained consistent with the measured values due to small RMSE and MAPE values. The experiments veriﬁed the proposed one-dimensional vertical inﬁltration method, which has applications in ﬁeld hydraulic parameter estimation.


Introduction
The soil water retention curve (SWRC) represents the soil's ability to store and release water as it is subjected to soil suction, and is described as the relation between soil moisture and soil suction [1,2]. The SWRC is used to evaluate the size and the distribution of soil pores and its water-holding capacity [3,4]. Soil hydraulic parameters can be obtained from the SWRC and used as key inputs for the water movement prediction that simulates the evaporation, the infiltration, and the runoff dynamics of field soils. Further, these parameters can be used to estimate the functions of various hydraulic properties of unsaturated and saturated soils to model solute transport [5,6]. The soil water retention behavior is important for many applications.
The soil hydraulic parameters primarily include θ s , θ r , α, and n, which can be estimated from the SWRC [7]. Numerous techniques are employed for the SWRC determination, such as the pressure plate method, the filter paper method, the tensiometer method, the osmotic suction control method, and the centrifuge method. Xing et al. [8] indicated that soil samples exhibit compressibility under a centrifugal force when using the centrifuge method-which can lead to clear changes in the bulk density-and similar problems are found with the tensiometer technique and the pressure plate technique [9][10][11]. The collected soil moisture data for the SWRC requires correction. Previous research [12] called for these corrections because accurate data collection during SWRC measurement enhances the accuracy of the obtained hydraulic parameters. These methods are time-consuming, requiring a measurement of the sedimentation of soil samples, whereas the change in bulk density is easily ignored [6,13,14]. Due to the shortcomings of SWRC data collection, direct SWRC measurements were avoided in this study.
The key parameters α and n were estimated through one-dimensional vertical infiltration experiments, followed by SWRC construction. The adopted method (i.e., the field-infiltration method) matched the actual field situation, where water infiltration is common. In order to simplify field operation, we used simplified visual observation and calculation for parameter estimation. The relationships describing the cumulative infiltration and the infiltration rate with the depth of wetting front are presented.
The main purposes of this study were to verify if the selected method was effective for parameter estimation during SWRC construction and to verify if the obtained SWRC was effective for the cumulative infiltration prediction. Specifically, the aims were (1) to provide a simplified and practical method to estimate soil hydraulic parameters and SWRCs and (2) to investigate the effect that the estimated SWRC obtained through the one-dimensional vertical infiltration method had on cumulative infiltration simulations when using SWMS_2D (Version 1.21, U.S. Department of Agriculture, Riverside, CA, USA), a Fortran-based open-source software package, which is simpler than HYDRUS for modeling soil water dynamics [8,15].

Theory
The one-dimensional infiltration process can be described by Equation (1). The initial and boundary conditions are shown in Equations (2) and (3), respectively.
where D represents unsaturated soil water diffusivity (cm 2 /min), K represents unsaturated hydraulic conductivity (cm/min), z represents the depth from the soil surface (cm, taken positive downward), H z represents the depth of simulation domain (cm), and t represents time (min). Parlange [16] stated that the value of The infiltration rate i at any time equals the water flux at the soil surface (Equation (5)): Equation (5) is described as Equation (6) in the form of soil suction, with combinations of K(θ) = C(θ)D(θ) and C(θ) = −dθ/dh.
From Equation (8), the K corresponding to any depth z could be obtained Equation (9). When z equaled the depth of the wetting front (i.e., z = z f ), θ denoted the initial soil moisture. K could be ignored under a low soil moisture condition. Under this circumstance, Equation (8) was simplified into Equation (10).
The cumulative infiltration corresponding to z f is expressed as: where I represents the cumulative infiltration (cm) and θ 0 represents the initial soil moisture (cm 3 /cm 3 ). The Taylor series expansion (first order) of Equation (14) was approximately obtained, displayed as Equation (16). Equation (16) was substituted into the Equation (15) in order to obtain Equation (17).
Linear relationships between i and 1/z f and I and z f were taken from Equations (18) and (19).
The parameters α and n were estimated through the following Equation:

Soil Sample Collection and SWRC Measurement
Loam and sandy loam ( Table 1) were selected from the cultivated field soil layers in Liaoning, Shandong, Shaanxi, and Xinjiang Provinces in China. Liaoning Province is located on the northeast plain, where cultivated farmland soil is primarily nutrient-rich black soil. Shandong Province has extensive river systems, which feed into the Bohai Sea and the Yellow Sea, resulting in moist and brown soils. Shaanxi Province is located on the Loess Plateau, which is full of loose soils. Xinjiang region has strong atmospheric evaporation and little rainfall, creating a saline-alkali soil. Soil samples were selected to represent their respective regions.
A moderate amount of scattered soil from each sample was compacted into cutting rings (100 cm 3 ) at the designated bulk density (Table 1), and then saturated with distilled water for 48 h. The soil samples were placed into a centrifuge (CR21G II, HITACHI, Tokyo, Japan) to obtain SWRC data for the four soils. In order to avoid compressibility during centrifugation, some low suctions were set, including 100 cm, 200 cm, 300 cm, 400 cm, 500 cm, 600 cm, 700 cm, 800 cm, and 1000 cm. After the equilibrium time for a particular suction was obtained, the soil samples were then removed from the centrifuge and weighed on an electronic balance with a precision of 0.01 g. The soil samples were oven dried at 105 • C until a constant weight was obtained. The SWRCs for the four experimental soils were constructed in accordance with the θ-h parameter planes. Note: a Soils collected from Liaoning Province; b soils collected from Xinjiang Province; c soils collected from Shandong Province; d soils collected from Shaanxi Province; e θ s represents the saturated soil moisture; f θ r represents the residual soil moisture; g γ d represents the initial dry bulk density; h SOC represents the soil organic carbon.

Laboratory One-Dimensional Infiltration
The soils were compacted into the Plexiglas columns at the designated dry bulk density in order to simulate the in situ bulk density of a local cultivated field, with an initial soil moisture of approximately 0.25 cm 3 /cm 3 . The depth of soil compaction was 20 cm. The soils were compacted into soil columns and left for 24 h until the infiltration experiment was conducted. Each treatment was performed thrice. A graduated Marriott bottle filled with distilled water was connected to the top of each soil column to maintain a steady infiltration head at 5 cm, as well as a continuous water supply. The water level in the Marriott bottle and the depth of wetting front were recorded. The time interval was 1, 2, 3, and 5 min during the processes of durations 1~10, 10~30, 30~60, and >60 min, respectively. The infiltration finished when the wetting front reached a depth of 20 cm. The cumulative infiltration and the depth of the wetting front were used to estimate the α and n parameters and for further constructing the SWRC.

Error Analysis
Root mean square error (RMSE), the mean absolute percentage error (MAPE), and the mean absolute error (MAE)-Equations (21)-(23), respectively-were selected to conduct a differential analysis between the measured and the calculated soil moisture. The first two indicators evaluated the differences between the measured cumulative infiltration and the calculated measurement using the parameters obtained via the infiltration method. A small value of RMSE, MAPE, or MAE indicated a high predictive accuracy [19].
where p is the sample size, V i is the ith measured soil moisture or cumulative infiltration, andV i is the ith calculated soil moisture or the cumulative infiltration.

Comparison of the Fitted and Calculated Values of α and n
The measured data from the thrice-replicated samples of all four soils demonstrated a similar tendency between i and 1/z f and I and z f , showing that i and I increased as 1/z f and z f increased. The mean values were used to construct Figure 1 to highlight the relationships. Equations (18) and (19) revealed that the non-zero-axial and zero-axial linear functions were used to fit i-1/z f and I-z f curves. The values of parameters A 1 , A 2 , and A 3 (listed in Table 2) displayed a good linear relationship between i and 1/z f and I and z f due to a high R 2 (>0.95). The results justified Equations (18) and (19).
The parameter α and n values were calculated using Equation (20). The fitted values were obtained from the RETC software (Agricultural Research Service, Riverside, CA, USA) using the measured data from the centrifuge method. Figure 2 shows that the α and n datasets were around the 1:1 line. The mean values of the relative error between the fitted and the calculated α were 8.27%, 9.73%, 3.33%, and 7.00% for the soil from Liaoning, Xinjiang, Shandong, and Shaanxi. The values for samples between the fitted and the calculated n were 3.49%, 1.91%, 7.21%, and 2.80%.

Comparison of the Measured and Calculated SWRCs
The SWRC generated via the measured θ-h data planes, obtained from centrifuge method, was defined as the measured SWRC. The generated SWRC, combined with the α and n obtained from the infiltration method, was defined as the calculated SWRC. The measured and calculated SWRCs for the four soils (plotted in Figure 3) demonstrated a good identity of the two types of the SWRC. The measured soil moisture for each tested suction was similar to the fitted soil moisture. Further error analysis showed that the correlation coefficient r between the measured and the calculated soil moisture for all four soils were larger than 0.99. The RMSE and MAE values were slightly larger than 0, demonstrating a near perfect estimation result. This result was also due to small MAPE values ( Table 3). The observed SWRCs of the experimental soils and the error values indicated that using the one-dimensional vertical infiltration method for parameter estimation (i.e., α and n) and for the SWRC estimation was reasonable.
Water 2018, 10, x FOR PEER REVIEW 8 of 11

Comparison of the Measured and Calculated SWRCs
The SWRC generated via the measured θ-h data planes, obtained from centrifuge method, was defined as the measured SWRC. The generated SWRC, combined with the α and n obtained from the infiltration method, was defined as the calculated SWRC. The measured and calculated SWRCs for the four soils (plotted in Figure 3) demonstrated a good identity of the two types of the SWRC. The measured soil moisture for each tested suction was similar to the fitted soil moisture. Further error analysis showed that the correlation coefficient r between the measured and the calculated soil moisture for all four soils were larger than 0.99. The RMSE and MAE values were slightly larger than 0, demonstrating a near perfect estimation result. This result was also due to small MAPE values ( Table 3). The observed SWRCs of the experimental soils and the error values indicated that using the one-dimensional vertical infiltration method for parameter estimation (i.e., α and n) and for the SWRC estimation was reasonable.  Table 3. Error analysis for the measured and calculated soil moisture contents.

Performance of Cumulative Infiltration Simulation
The cumulative infiltration was simulated via SWMS_2D software, using the α and n parameters obtained by the infiltration method. The changes of simulated cumulative infiltration over time remained consistent with the observed values, although some gaps were exposed (Figure 4). The errors demonstrated that the RMSE and MAPE values between the observed and simulated values were 0.0698 cm and 2.0234% for the sandy loam from Liaoning, 0.1453 cm and 4.3688% for the sandy loam from Xinjiang, 0.1139 cm and 3.5108% for the loam from Shandong, and 0.1529 cm and 2.7700% for the loam from Shaanxi. The small RMSE and MAPE values indicated that the infiltration method had potential for parameter estimation (i.e., α and n), which could be considered as the inputs for the SWMS_2D.
Water 2018, 10, x FOR PEER REVIEW 9 of 11 The cumulative infiltration was simulated via SWMS_2D software, using the α and n parameters obtained by the infiltration method. The changes of simulated cumulative infiltration over time remained consistent with the observed values, although some gaps were exposed ( Figure  4). The errors demonstrated that the RMSE and MAPE values between the observed and simulated values were 0.0698 cm and 2.0234% for the sandy loam from Liaoning, 0.1453 cm and 4.3688% for the sandy loam from Xinjiang, 0.1139 cm and 3.5108% for the loam from Shandong, and 0.1529 cm and 2.7700% for the loam from Shaanxi. The small RMSE and MAPE values indicated that the infiltration method had potential for parameter estimation (i.e., α and n), which could be considered as the inputs for the SWMS_2D.

Summary and Conclusions
In the present study, three sets of α and n values were calculated through a one-dimensional infiltration experiment (Figure 2). The standard deviations of α were 1.34 × 10 −4 , 8.32 × 10 −4 , 3.26 × 10 −4 , and 8.08 × 10 −4 , and n reached 0.0139, 0.0170, 0.0013, and 0.0116 for the soil from Liaoning, Xinjiang, Shandong, and Shaanxi, respectively. The small standard deviations indicated the calculated α and n values were accurate via the infiltration method, although additional soils and duplicate tests remain necessary. Ma et al. [20] and Wang et al. [21] proposed a horizontal infiltration method to estimate the parameters of the Brooks-Corey model. The present study adopted the vertical infiltration method for parameter estimation, which would be more practical since water infiltration is a more common in the cropland, particularly after irrigation events. In these studies, undisturbed 20 cm depth soil samples were obtained for the hydraulic parameters that were estimated through the vertical infiltration method, which would save more time than the direct SWRC measurements. Su et al. [22] used second-, third-, and fourth-order Taylor series expansions to obtain an approximate analytical solution for a one-dimensional, constant water head, vertical

Summary and Conclusions
In the present study, three sets of α and n values were calculated through a one-dimensional infiltration experiment (Figure 2). The standard deviations of α were 1.34 × 10 −4 , 8.32 × 10 −4 , 3.26 × 10 −4 , and 8.08 × 10 −4 , and n reached 0.0139, 0.0170, 0.0013, and 0.0116 for the soil from Liaoning, Xinjiang, Shandong, and Shaanxi, respectively. The small standard deviations indicated the calculated α and n values were accurate via the infiltration method, although additional soils and duplicate tests remain necessary. Ma et al. [20] and Wang et al. [21] proposed a horizontal infiltration method to estimate the parameters of the Brooks-Corey model. The present study adopted the vertical infiltration method for parameter estimation, which would be more practical since water infiltration is a more common in the cropland, particularly after irrigation events. In these studies, undisturbed 20 cm depth soil samples were obtained for the hydraulic parameters that were estimated through the vertical infiltration method, which would save more time than the direct SWRC measurements. Su et al. [22] used second-, third-, and fourth-order Taylor series expansions to obtain an approximate analytical solution for a one-dimensional, constant water head, vertical infiltration. This study, however, used first-order Taylor series expansion for parameter estimation with satisfying accuracy, which simplified the calculation. The SWRCs were ultimately determined based on the calculated α and n, which were considered the inputs for the accumulation infiltration simulation. This justified the obtained SWRCs via the one-dimensional vertical infiltration method.
The theoretical relationships between I and i and z f were established, and a simplified parameter for estimating expressions of α and n was proposed. One-dimensional vertical infiltration experiments of the four loam and sandy loam soils showed linear relationships between i and 1/z f and I and z f . The calculated values of α and n obtained through the infiltration method were consistent with the fitted values obtained from the RETC software (Agricultural Research Service, Riverside, CA, USA), with a relative error of 3.33-9.73% and 1.91-7.21%. The values of soil moisture for the calculated SWRCs were similar to the measured SWRCs, with a RMSE of 0.0018-0.0101 cm 3 /cm −3 , a MAE of 0.0012-0.0100 cm 3 /cm −3 , and a MAPE of 0.3806-5.1445%. The SWMS_2D-based cumulative infiltration simulations for the four soils revealed that the simulated values (on the basis of calculated α and n) were consistent with the observed values, demonstrated by the small RMSE (0.0698-0.1529 cm) and MAPE (2.0234-4.3688%).