Modeling of soil water retention curves based on two programs

We compared two methods of modeling soil water retention curves used by scientists in Slovakia. The first modeling was done using the GENRET program and the second using the RETC program. Samples of pure sandy soil and sandy soil with applied biochar in three different particle sizes were used for the simulation. Sandy soil has a very low retention capacity therefore the modeling of soil water retention curves is not easy. Our results showed that the GENRET program can model the curve even at high pressures, which the RETC program modeled only in one variant, but the RETC program had better agreement with measured data. The most significant differences between the programs were at the lowest and highest pressures.


Introduction
The soil water retention curve (SWRC) is the relationship between water content and matric potential. It is one of the most significant hydraulic functions for modeling flow transport in porous media (Zha et al., 2015). SWRC is one of the most important and complex hydro-physical characteristics of the soil, because we can determine: a) energetic characteristics of soil water, b) qualitative and quantitative characteristics of soil pores, c) hydrolimits, d) changes in physical and hydro-physical characteristics due to the action of natural factors or human activity, e) design parameters for irrigation and soil drainage and f) input data for the calculation of other soil hydro-physical characteristics (Antal and Fídler, 1989). Therefore, SWRC is considered as one of the most fundamental soil hydraulic properties, and the accurate acquisition of SWRC and its parameterization have great significance in understanding soil moisture dynamic and soil hydrology. There are two main approaches to obtaining SWRC. The first is experimental determination and the second is derivation from basic soil properties by using pedotransfer functions (PTF). Although, experimental approach is time consuming and costly, it is undoubtedly more precise and reliable for SWRC of specific soils (Pan et al., 2019). SWRC is affected by many soil environmental factors, and thus it is very difficult to accurately determinates relationship with these factors. The course of SWRC depends on the granular and mineralogical composition of the soil, on the content of humus, exchangeable cations, structure and bulk density.
For that reason, it is established for each soil, or for each horizon separately (Fulajtár, 2006). The SWRC can be measured by a variety of techniques. Among the available techniques, the most commonly used are the Richards' pressure apparatus (Dane and Hopmans, 2002;Richards, 1948;Richards, 1965), also called pressure plates (and hereafter called pressure plates apparatus). We have used this technique for our experiment, as well. The SWRCs are usually measured during the drainage process, starting with soil saturated with water (this is the preferred case), or during soil wetting, then starting with air-dry soil. Both SWRC main branches are typical of a particular soil; the shape and position of the SWRC reflect porous space dimensions, and porous space reflect soil texture (Novák and Hlaváčiková, 2019). Several mathematical formulations for describing the SWRC are available (Brooks and Corey, 1966;van Genuchten, 1980;Vogel and Cislerova, 1988). Early empirical models primarily describe the wet end of the SWRC. For example, the Brooks and Corey (1966) and van Genuchten (1980) models (abbreviated the BC and VG models, respectively) are the most popular for predicting SWRC under wet conditions. However, the two models do not predict the SWRC at oven dryness, because they assume the matric suction to be infinite when the soil moisture approaches the residual water content (Silvia and Grifoll, 2007). According to Bittelli and Flury (2009) directly observed soil suction is limited from 0 up to 1500 kPa (15 bars) and the measurements are time consuming and laborious at high suction values. In our experiment we used a sandy soil, as not very typical agricultural soil, but there are several countries and its locations where this type of soil has to be used for agricultural utilization. The biochar is one of tools to improve hydro-physical parameters and fertilization of sandy soils (Dokoohaki et al., 2017;Toková et al., 2022), that's why we present three variants of SWRCs with biochar and one variant without biochar. However, in this study we present only SWRC modeling differences and program selection.

Preparation of soil and soil-biochar mixtures
The sandy soil used in this experiment was taken from Plavecký Štvrtok area at Záhorská lowland. Particle size distribution was measured by the hydrometer method (Velebný, 1981). It consists of 91% sand, 7.5% silt and 1.5% clay, so it is classified as sand, based on the USDA classification. The soil was sieved to a fraction with particle diameter size ≤2 mm. Used biochar was produced from willow trees by pyrolysis temperature at 300°C. The size of biochar was 0-10 mm, so it was homogenized using a hammer mill and sieved to a fraction with a particle size 125 µm and 2 mm. Samples of soil-biochar mixtures were prepared in laboratory conditions with a biochar application rate of 20 t ha -1 in Kopecky cylinders with a volume of 100 cm 3 . Overall four variants were prepared: pure sandy soil (control), mixture of soil and biochar with particle size <125µm, mixture of soil and biochar with particle size 125µm-2mm and mixture of soil and biochar with particle size >2mm. For each variant were prepared three repetitions.

Overpressure equipment (pressure plate method)
To determine the SWRC of the soil, we used 9 measurement points at pressure potentials of 0; 0.002; 0.06; 0.1; 0.3; 0.56; 1; 3 and 4.8 bars. For the analytical construction of the SWRC, it was necessary to calculate the residual soil moisture -θr, because the relationship between soil moisture and soil water potential is determined only in the range of θr and θs (saturated water content) (Skalová et al., 2015). We determined the residual moisture of the soil using the measured points with the trend line and saturated soil moisture was measured. In this study, we present the average soil moisture values for each variant.

GENRET
The GENRET program is a part of the GLOBAL model, which was developed at the Institute of Hydrology SAS (Majerčák and Novák, 1995). Non-linear least-squares analysis of the SWRC is approximated by van Genuchten method (1980) and SWRC hysteresis is considered based on Kool-Parker theory (1987). Input data are saturated hydraulic conductivity, saturated water content and measured points of SWRCs. Residual water content could be computed or fitted; parameters alpha and n can be computed based on Mualem's or Burdine's theory or fitted. In our study were parameters alpha and n determinated based on the Mualem's theory. The unsaturated hydraulic conductivity is computed based on Mualem's theory (1976).

RETC
The RETC (RETention Curve) is a program for analyzing the hydraulic conductivity properties of unsaturated soils (van Genuchten et al., 1991). The parametric models of Brooks-Corey and van Genuchten are used to represent the SWRC. The theoretical pore-size distribution models of Mualem and Burdine predict the unsaturated hydraulic conductivity function. The simulation can be generated from observed soil water retention data assuming that one observed conductivity value (not necessarily at saturation) is available (van Ganuchten et al., 1991). The program also permits users to fit analytical functions simultaneously to observed water retention and hydraulic conductivity data. The RETC program assign fixed parameters alpha and n for sandy soil. In our study we used van Genuchten model for SWRC simulation. It is widely used expression for describing the soil water retention function (1): is soil water potential, r and sdenote residual and saturated water contents, respectively, n and mare dimensionless empirical shape factors, 1 is a pore-size distribution index. For practical purposes it is also denoted as n in the program RETC, α is also empirical shape parameter.

Statistical analysis
Comparison of two programs for SWRC simulations was evaluated using the Microsoft Excel software. The statistical analysis was done by using standard box plot to show distribution of a set of data.

Results and discussion
Our results show that there are deviations in SWRC modeling when using two different programs ( Fig. 1-4). The biggest differences were detected at the highest modeled pressuresa significant deviation between the curves was detected at pressures of pF 2.5 and above. It follows from the modeled SWRCs that at pF 3 soil moisture was 50% lower in the GENRET program compared to the RETC program, and at pF 3.5 it was up to 76% (Fig. 1). Higher differences in soil moisture were also detected for the variant with a particle size of biochar >2mm (Fig. 4), when at pF 3 it was 40% for the GENRET program and at pF 3.5 to 66% compared to the RETC  Average SWRC for sandy soil and biochar with particle size <125µm modeled by GENRET and RETC programs. Fig. 3. Average SWRC for sandy soil and biochar with particle size 125µm-2mm modeled by GENRET and RETC programs.

GENRET RETC
program. At saturated soil values and low pressures, the difference was negligiblemaximum to 6% (Fig. 3  and 4). Almost identical SWRCs were simulated between GENRET and RETC programs for the variant of soil and biochar with particle size <125µm. The most significant difference was 13% at pF 3.5 (Fig. 2). From Fig. 1-4 it is clear that the GENRET program is able to model soil moisture even at higher pressures, despite the fact that the input data for sandy soil modeling deflect from the standard values. The statistical analysis (Fig. 5) shows that while the average soil moisture was modeled at 10% vol. with the GENRET program, it was around 18% vol. with the RETC program. For better visualization are differences in measured and modeled SWRC points shown in Table 1. Better agreement was between measured data and modeled data by RETC program. In higher pressures were differences between programs and measured data higher, sometimes underestimated (more often GENRET program) and sometimes overestimated (more often RETC program). Soil moisture for pF 3.7 was not modeled by RETC program for variant soil + biochar >2mm. The shape differences of the SWRCs were due to the different alpha and n parameters ( Table 2). The fitted alpha parameter in RETC program was lower than the alpha values calculated by GENRET program; the parameter n was significantly higher in RETC compared to the calculated n value for all variants in GENRET.

Fig. 4.
Average SWRC for sandy soil and biochar with particle size >2mm modeled by GENRET and RETC programs.

Conclusion
Sandy soil has a very low retention capacity; therefore, it is not completely suitable for agricultural use. By using biochar, it is possible to improve its basic hydro-physical properties. SWRC is one of the most important parameters that can be used to assess soil retention capacity. Measuring SWRC points for sandy soil is problematic because the infiltrated water drains away very quickly (even at the lowest pressures). In our experiment, soil moisture was significantly reduced already at a pressure of 0.06 bar (6 kPa), and from a pressure of 0.56 bar (56 kPa) it changed indistinctly. Therefore, SWRC modeling is also more complicated, and not every model can handle non-standard input data. We compared the SWRC simulations by using the GENRET and RETC programs. They showed dispersion mainly at higher pressures and also when applying biochar with the largest particles. The RETC program failed to complete the modeling at the highest pressures and only for the soil and biochar <125µm variant modeled SWRC above pF 5 (10 000 kPa). It can be explained by the fact, that the smallest biochar particles (<125µm) can fill the pores between the sandy soil particles and thereby increase its retention capacity, which was also reflected in the values of the input data to the programs. Other way round, the GENRET program had a problem only with the variant with pure soil, when the modeling ended at a pF 4.8 (6 309 kPa). For the other variants, it reached a value of pF 5.5 (31 622 kPa). Measured SWRC points had better agreement with RETC program. We can conclude that both programs have their limitations for simulations the SWRC of sandy soil. Our next research is continuing by using different method of θr determination for the GENRET program and calculating parameters alpha and n for RETC program.