2020-08-09 / Blog 5

Fitting soil water retention curve data in R

The relationship between water content (θ) of the soil and its matric suction (Ψ), or water potential, is defined using the soil water retention curve (SWRC). One of the most commonly used equations to define the SWRC of a soil is proposed by van Genuchten (1980). In terms of usual geotechnical parameters, SWRC equation can be expressed as follows:

where θs and θr are the saturated and residual water content, α and n are the fitting parameters of the equation. α is related - but not equal - to the air-entry value of the soil, whereas n controls the slope of the mid-section of the curve.

You can determine SWRC in many ways. But it can be summarised in two main categories as suction-controlled and suction-monitored. Suction-controlled techniques are mainly used in the laboratory, whereas suction-monitoring techniques can be applied both under laboratory and field conditions. Whatever technique you use, what you end up as data is combinations of volumetric, or gravimetric, water content and matric suction. What you need after you obtain the data is an algorithm to fit the equation of van Genuchten (1980) to the measured data.

We start with defining a function to calculate the volumetric water content at a discrete value of matric suction.

vg <- function(suction,thetas,thetar,a,n)

{

((thetas-thetar)/((((a*suction)^n)+1)^(1-(1/n))))+thetar)

}

We will use the nls function in R to estimate the fitting parameters of the van Genuchten (1980) equation using a non-linear least squares algorithm.

fit <- nls(y~vg(x,thetas,thetar,a,n),start=list(thetas=42,thetar=5,a=0.33,n=5))

Here, y is the measured volumetric water content values, while x is the measured matric suction values. You need to enter a rough estimate of starting values for the fitting parameters. R will try to estimate your fitting parameters using a convergence criteria. If successful, you will end up with such a screen:

Formula: y ~ vg(x, thetas, thetar, a, n)

Parameters:

Estimate Std. Error t value Pr(>|t|)

thetas 4.230e+01 8.654e-02 488.8 <2e-16 ***

thetar 5.942e+00 5.255e-02 113.1 <2e-16 ***

a 2.703e-01 2.953e-04 915.4 <2e-16 ***

n 7.422e+00 3.931e-02 188.8 <2e-16 ***

---

Signif. codes: 0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1


Residual standard error: 0.7644 on 1535 degrees of freedom


Number of iterations to convergence: 9

Achieved convergence tolerance: 6.835e-06

We see that the estimation has been performed successfully in 9 iterations, and all our fitting parameters are statistically significant. We obtain the fitting parameters as θs = 42.3%, θr=5.94%, a=0.27 and n=7.42.