Fitting functions of soil hydraulic properties

Easiest way to start learning how to use unsatfit is to run sample codes as instructed in this page. You can optimize parameters of WRF (water retention function) and HCF (hydraulic conductivity function) of various available models to measured data set. You can test with sample data provided in this page. For optimizing only WRF parameters, please use SWRC Fit.

**List of sample codes**

- Models
- Multiple curves
- Contour plot

where

\[\begin{eqnarray} S(h) &=& \biggl[\dfrac{1}{1+(ah)^n}\biggr]^m\\ m &=& 1-1/n \end{eqnarray}\]- Install Python 3 and unsatfit. Install pandas by
`python -m pip install pandas`

. - Prepare datafile as csv (comma-separated values) format for (h, θ) data in the filename “swrc.csv” with a header of “h, theta”. See sample data of clay (Unsoda 2362).
- Prepare (h, K) data in the filename “hcc.csv” with a header of “h, K”. See sample (h, K) data of clay (Unsoda 2362).
- Download a sample code for VGM model.
- Read the sample code and edit initial parameters, bound of parameters and figure setting as written in the comment of the sample code if necessary. For more figure options, see source code of unsatfit.
- Run the sample code at the same directory with data files (swrc.csv and hcc.csv). For running the code on Mac or unix-like system, edit the first line (shebang) and mark the file executable by
`chmod +x VGM.py`

. For running on Windows, please refer to Python on Windows FAQ. - It first optimizes WRF parameters (θ
_{s}, θ_{r}, a, n) of VG model, and then optimizes HCF parameters (K_{s}, p) of VG Mualem model, or modified VG model (h_{b}=2) Vogel et al. (2000) when n > 1.1. h_{b}value can be changed as HB value in the sample code. - Fitted parameters are shown at the standard output, where qs and qr means θ
_{s}and θ_{r}respectively, and R2 q means R^{2}for θ of water retention curve and R2 logK means R^{2}for log(K) of hydraulic conductivity curve. - Note that the program is unit independent, meaning that the unit of the parameters depends on the unit of the input data. Unit of pressure head is assumed as cm for a (cm
^{-1}) and h_{b}. - Figure file is produced as VG.png. For use in papers, pdf file can be produced as instructed in the sample code.

Result with sample data of clay (Unsoda 2362) is shown below.

Result with sample data of Gilat loam is shown below. In this case, bimodal model is appropriate, as shown below.

Bimodal model with general HCF (Seki et al., 2022) can represent water retention and hydraulic conductivity of various types of soil in a wide range of pressure head, as verified in Seki et al., 2023. You can conduct the same fitting as written in the paper by using the following sample codes.

where

\[\begin{eqnarray} S(h) &=& \begin{cases}w S_1(h) + (1-w)\left(h/H\right)^{-\lambda_2} & (h>H)\\ w S_1(h) + 1-w & (h \le H)\end{cases}\\ S_1(h) &=& Q \biggl[\dfrac{\ln(h/H)}{\sigma_1}\biggr]\\ Q(x) &=& \frac{1-\mathrm{erf}(x/\sqrt{2})}{2}\\ b_1 &=& w \exp\biggl(\frac{q^2 \sigma_1^2}{2}\biggr)\\ b_2 &=& (1-w)\biggl(\frac{q}{\lambda_2} + 1\biggr)^{-1}\\ \gamma &=& b_1+b_2 \end{eqnarray}\]- Use sample code for KBC model. See instruction in the VGM model above.
- It first optimizes WRF parameters (θ
_{s}, w, H, σ_{1}, λ_{2}) of KBC model (θ_{r}=0), and then optimizes general HCF parameters (K_{s}, p, q) of KBC model (r=1), or modified KBC model (h_{b}=2, r=1) when σ_{1}> 2. - Result with sample data of Gilat loam is shown below.

where

\[\begin{eqnarray} S(h) &=& w\bigl[1+(ah)^{n_1}\bigr]^{-m_1} + (1-w)\bigl[1+(ah)^{n_2}\bigr]^{-m_2}\\ m_i&=&1-1/{n_i}\\ \Gamma_i(h) &=& 1-\biggl[1-\big[1+(ah)^{n_i}\bigr]^{-1}\biggr]^{m_i} \end{eqnarray}\]- Use sample code for DVC model. See instruction in the VGM model above.
- It optimizes WRF parameters (θ
_{s}, w, a, n_{1}, n_{2}) , and then optimizes general HCF parameters (K_{s}, p, r) of DVC model or modified DVC model (h_{b}=2) when n_{1}or n_{2}is smaller than 1.1. - Result with sample data of Gilat loam is shown below.

where

\[\begin{eqnarray} S(h) &=& \begin{cases}w \left(h / H\right)^{-\lambda_1} + (1-w)\left(h / H\right)^{-\lambda_2} & (h>H)\\ 1 & (h \le H)\end{cases}\\ B_i &=& \biggl(\frac{q}{\lambda_i} + 1\biggr)^{-1}\\ \gamma &=& wB_1+(1-w)B_2\\ \Gamma_i(h) &=& \left(h / H\right)^{-\lambda_i-q} \end{eqnarray}\]- Use sample code for DBC model. See instruction in the VGM model above.
- It optimizes WRF parameters (θ
_{s}, w, λ_{1}, λ_{2}), and then optimizes general HCF parameters (K_{s}, p, q) of DBC model. - Result with sample data of Gilat loam is shown below.

where

\[\begin{eqnarray} S(h) &=& w S_1(h) + (1-w)S_2(h)\\ S_1(h) &=& Q \biggl[\dfrac{\ln(h/H)}{\sigma}\biggr]\\ Q(x) &=& \frac{1-\mathrm{erf}(x/\sqrt{2})}{2}\\ S_2(h) &=& \begin{cases}\dfrac{L(h_0)-L(h)}{L(h_0)-L(H)} & (h>H)\\ 1 & (h \le H)\end{cases}\\ L(h) &=& \ln(1+h/H)\\ K_1(h) &=& {S_1(h)}^p \Biggl[ Q \biggl[\dfrac{\ln(h/H)}{\sigma} + \sigma \biggr] \Biggr]^2\\ K_2(h) &=& \begin{cases}(h/H)^a & (h>H)\\ 1 & (h \le H)\end{cases}\\ \end{eqnarray}\]- Use sample code for Peters model. See instruction in the VGM model above.
- h
_{0}=6.3×10^{6}is constant and can be edited as H0 value in the sample code. - It optimizes WRF parameters (θ
_{s}, w, H, σ), and then optimizes general HCF parameters (K_{s}, p, a, ω) of Peters model or modified Peters model (h_{b}=2) when σ > 2. - Result with sample data of Gilat loam is shown below.

- Use sample code for multiple curves. See instruction in the VGM model above.
- Draw DBC, DVC, KBC and Peters model in the same figure.
- Result with sample data of Gilat loam is shown below.

- Use sample code for contour plot. See instruction in the VGM model above.
- It draws a contour plot of (p, q) for RMSE of estimated log
_{10}(K) for KBC (r=1) optimization. - Result with sample data of Gilat loam is shown below.