Original Article

Predicting land surface movement based on experimentally determined soil shrinkage behaviour

Bente R. Lexmond¹, Gilles Erkens^1,2,3, Jasper Griffioen^4,5, Sanneke Van Asselen², Laura Pedretti⁶, Claudia Meisina⁶ and Esther Stouthamer¹

¹Department of Physical Geography, Utrecht University, Utrecht, The Netherlands; ²Deltares research Institute, Delft, The Netherlands; ³Department of Biology, Utrecht University, Utrecht, The Netherlands; ⁴TNO Geological Survey of the Netherlands, Utrecht, The Netherlands; ⁵Copernicus Institute of Sustainable Development, Utrecht University, Utrecht, The Netherlands; ⁶Department of Earth and Environmental Sciences, University of Pavia, Pavia, Italy

Abstract

Understanding and predicting shallow subsurface deformation and resulting vertical land surface movement due to shrink-swell behaviour of clayey soils is crucial for infrastructure resilience and evaluating long-term land surface elevation trends. Specifically, because shrink-swell movements seem to increase with increasing drought intensity, duration and occurrence. This study aimed to identify a universally applicable and straightforward method to predict shallow subsurface deformation based on changes in volumetric water content (VWC) or soil water suction and laboratory-measured shrinkage behaviour. An existing method based on changes in VWC and a newly proposed method based on changes in soil water suction were intercompared. The newly proposed method relates the vertical- or volumetric shrinkage to the matric suction, in combination with a depth-based overburden term. The input for both prediction models consisted of the shrinkage behaviour measured in a laboratory setup. The shrinkage behaviour was described as the change in sample- volume, height, the water content and soil water suction, during drying. Samples were collected from three sites in the Netherlands, where expansive clay minerals were identified and the models could be validated using extensometer data. Both prediction models proved to be able to capture seasonal shrinkage behaviour well. However, short-term (days to weeks) deformations were more challenging to predict due to the difficulty in calibrating VWC-sensors for the full VWC range. Moreover, seasonal swell was not predicted accurately for a site with a relatively low expansive clay mineral content. The Briaud model was the most accurate prediction model for shrinkage behaviour tested in this study.

Keywords: expansive soil; extensometer measurements; seasonal soil deformation; shrinkage and swelling

 

 

Introduction

Damage to infrastructure due to shrinkage and swelling of the shallow subsurface is a global issue and occurs due to the presence of expansive soils and alternating wetting and drying conditions (Bronswijk & Evers-Vermeer, 1990; Doris et al., 2008; Lu & Dong, 2017). While often associated with arid climates, these problems also affect temperate regions like the Netherlands, where expansive clay minerals are present in the subsurface (Griffioen et al., 2016). Annual shrinkage and swelling-induced damage costs reach tens of billions of US dollars in the USA and France, and hundreds of millions in countries including the UK, China, Saudi Arabia and Australia (Adem & Vanapalli, 2015; Azam, 2003; Fityus et al., 2004; Jones & Jefferson, 2012; Puppala & Cerato, 2009; Shi et al., 2002). In the Netherlands, increased damage to foundations and infrastructure due to weather events has been observed since the dry summer of 2018, particularly in areas with clay-rich soils (Van Essen, 2025). Climate change is expected to intensify extreme weather events (Seneviratne et al., 2021). More severe and frequent droughts further increase the potential for shrink–swell behaviour.

Expansive soils are classified as soils that change in volume upon changes in water content, primarily due to the presence of expansive clay minerals, such as smectite (Jones & Jefferson, 2012). Shrinkage and swelling behaviour are governed by the soil moisture variability and soil characteristics such as bulk density, grain size distribution (GSD), clay mineralogy, soil organic matter content (SOM), carbonate content, cation exchange capacity (CEC) and Atterberg limits (Boivin et al., 2004; Davidson & Page, 1956; Greene-Kelley, 1974; Mitchell, 1993; Peng & Horn, 2005, 2013). Also, land use (tractor traffic, tillage and cattle trampling), management practices (type of ploughing, application of compost and manure affecting soil porosity) and drying history influence the shrink and swell characteristics (Alaoui et al., 2011; Basma et al., 1996; Dörner et al., 2009; Vignozzi & Pellegrini, 2004) such as hysteresis (Estabragh et al., 2018). This list of characteristics highlights the complexity of predicting shrinkage and swelling behaviour.

Numerous studies have investigated drying and wetting cycles of expansive soils (e.g. Basma et al., 1996; Day, 1994, 1995; Estabragh et al., 2018; Osipov et al., 1987). These studies consistently show that drying intensity affects swell potential. When soils are dried below the shrinkage limit, that is, the water content at which a soil stops shrinking, the swelling potential increases due to entrapped air (Estabragh et al., 2018). However, when the water content of a soil reaches unprecedented low levels, the suction stresses can yield restructuring of the soil matrix due to increased interparticle stress and ongoing suction. Restructuring of the soil matrix can yield irreversible shrinkage or physical ripening, leading to land subsidence (Basma et al., 1996). Irreversible shrinkage can either increase or decrease the swelling potential, depending on the original microstructure and the characteristics of the bonding between particles (Basma et al., 1996). Additionally, the clay mineral orientation plays a role: edge-to-face structures reduce swelling potential, while face-to-face structures increase it (Basma et al., 1996; Osipov et al., 1987).

Observed changes in the amplitude of shrinkage and swelling due to increasing drying intensity stabilise after 3–8 shrinkage and swelling cycles (Basma et al., 1996; Day, 1994, 1995; Estabragh et al., 2018; Osipov et al., 1987). Osipov et al. (1987) found that the number of wetting and drying cycles necessary to obtain maximum swelling increases as the structural bonds strengthen. Estabragh et al. (2018) show that shrinkage and swelling hysteresis disappears when a soil sample has reached the shrinkage and swelling equilibrium. The equilibrium was defined based on a relation between the sample’s suction and void ratio. Equilibrium shrinkage and swelling can be expected for soils for which no extreme weather conditions have occurred in recent years or the environment remained unchanged in another way (Jones et al., 2009). Shrinkage tests have been successfully applied to predict field shrinkage behaviour (Briaud et al., 2003). In contrast, Giziensky and Lee reported that laboratory swelling tests may overpredict field swelling behaviour. This discrepancy is attributed to differences in the constraints imposed on soil volume change under laboratory and field conditions (Dhowian & Al-Sadaan, 2010). In laboratory tests, the lateral and vertical swelling constraints do not simulate those of the field accurately, with regard to loading and lateral mobility of the testing apparatus. Additionally, the particle orientation changes during shrinkage (Basma et al., 1996), whereby randomly oriented clay particles preferably rotate to the direction perpendicular to the direction of the major consolidation stress (Kirkpatrick & Rennie, 1972; Liu et al., 2021). The follow-up swelling pressure is greater in the direction perpendicular to the particle orientation, compared to parallel. In the field, the drying and wetting stresses are partly directed by surrounding soil layers, structural loads and natural soil heterogeneity, which also influences the particle orientation and subsequent swelling behaviour. If swell is more constrained in field conditions than shrinkage, irreversible volume loss will occur, resulting in land subsidence.

Several generalised models have been developed to predict shrinkage and swelling behaviour, using volumetric water content (VWC) (Briaud et al., 2003; Overton et al., 2006), suction (Adem & Vanapalli, 2015; Wray et al., 2005) or consolidation theory (Abed, 2008; Vu & Fredlund, 2004; Zhang, 2004), often combined with geotechnical parameters. These methods were reviewed by Adem and Vanapalli (2015). Among them, the Briaud model accurately predicted soil movement and stands out for its simplicity in terms of required parameters, input data and software. A key advantage is its use of VWC rather than suction or stress to predict strain, as VWC is a parameter that is generally more easily monitored than suction or stress. Measurement equipment for VWC is more robust during dry periods compared to tensiometers used to measure suction. Additionally, several models exist to estimate VWC profiles in the unsaturated zone above the groundwater table based on environmental and site-specific conditions.

The Briaud model uses in-situ measured or modelled VWC in combination with laboratory-measured shrinkage behaviour. Shrinkage behaviour is defined as the relationship between VWC and both volumetric and vertical strain during drying. The model has yielded predictions that agree reasonably well with measured land surface movements (Briaud et al., 2003). However, it tends to underestimate shrinkage and overestimate swell (Briaud et al., 2003). When VWC is monitored using Time Domain Reflectometry (TDR) or Time Domain Transmission (TDT) sensors, calibration across the full moisture range is essential in clay-rich soils. Water bound to clay minerals alters the dielectric constant, interfering with standard linear calibration. Proper calibration enables these sensors to capture full moisture dynamics, unlike tensiometers which capture suction up to a maximum of 140 kPa (TEROS 32, MeterGROUP), which is about 10% of the suction possible in natural Dutch soils (Bronswijk & Evers-Vermeer, 1990).

The fundamental basis of the model proposed by Briaud et al. (2003) is that the water content is directly linked to suction through the soil water retention curve. As moisture flows through the unsaturated subsurface due to gradients in mechanical energy, the moisture transport is typically described using soil water suction (Wray et al., 2005). Additionally, whereas VWC is influenced by soil properties like clay content, suction shows a continuous depth profile making it easier to simulate (Wray et al., 2005). Therefore, when soil moisture content variations are limited to the measurement range of tensiometers, a suction-based model can provide an easier tool to simulate shrinkage and swelling, based on the continuous suction profile.

As the shrinkage and swelling behaviour and amplitude depend on many variables, both in subsurface characteristics and local environmental factors, there is a need for a universal, simple method to predict shrinkage and swelling. Thus, the aim of this study is to test the applicability of the Briaud model and a newly proposed method to predict land surface movement due to shrinkage and swelling for three locations in the Netherlands. These predictions are based on changes in VWC (Briaud) or matric suction (newly developed), combined with laboratory-measured shrinkage behaviour. The objectives are to (1) accurately characterise shrinkage behaviour of undisturbed clay samples, and (2) simulate shallow subsurface deformation using measured VWC and suction data. For VWC-based simulations, we apply the Briaud et al. (2003) model. For suction-based predictions, we propose a simplified approach and test whether overburden pressure should be included. Simulations are validated using in-situ measured layer deformation. We hypothesise that (1) laboratory-measured shrinkage behaviour of undisturbed clay samples can be scaled to predict field-measured shallow subsurface deformation and associated land surface movement, (2) a substantial part of the observed land surface movement at the study sites can be explained by shrinkage–swelling-induced deformation of the shallow subsurface, and (3) matric suction proves to be a more suitable parameter to simulate shrinkage and swelling.

Methodology

To predict in-situ surface and shallow subsurface vertical movement, the shrinkage behaviour of samples extracted at three different sites was studied. The shrinkage behaviour was defined as the change in sample dimensions and the VWC measured in the laboratory setup. The shrinkage behaviour was measured as the decrease in the height and radius of cylindrical samples of 250 mL. The acquired data were used to predict in-situ shrinkage and swelling behaviour, based on the VWC measured in the field over periods ranging up to several years. The field setup was designed to measure the bulk VWC, which is the ratio between the volume of water and the bulk volume of the soil. This is also the definition of VWC that will be used within this paper. The predictions were based on the approach of Briaud et al. (2003) and using the direct application of the shrinkage curve defined as the normalised- sample height and volume in relation to the measured matric suction. The sampling method, lab setup and prediction methodologies are described first (‘Sample collection and preparation’ section to ‘Testing the performance of the prediction methods’ section) and thereafter the site descriptions and sample characteristics (‘Site 1: Hazerswoude’ section to ‘Site 3: Nieuwolda’ section).

Sample collection and preparation

The sample extraction sites were located in the vicinity of extensometer monitoring stations, that include field VWC measurements in the unsaturated zone at depths with sediments containing around 10% or more clay-sized particles by volume. The samples were extracted at depths ranging between 30 and 130 cm below the land surface, using a ring soil sampling kit (Sample ring kit 0784SC, Ø 84 mm, Eijkelkamp Soil & Water B.V.) to minimise disturbance. The rings were either carefully pushed or hammered into the subsurface. From the sampling holes, groundwater was collected just below the groundwater table. The water in the hole was pumped out and after waiting for it to have refilled, the water was extracted using a peristaltic pump. The water in the hole was pumped out to minimise the amount of sediment in the water. The water was fully filled in plastic containers and stored in a dark and cool (4°C) place before saturation of samples could be started. The samples were extracted at all sites during March (Site 3), April (Site 1) 2023 and February 2024 (Site 2).

The rings containing the samples were sealed airtight, transported to the lab and saturated with the collected groundwater. The samples had a diameter of 80 mm and a height of 50 mm. The saturation was carried out by submerging the samples up to 80% of their height in glass desiccators for 3 weeks, starting maximum 3 days after sampling. The samples were not submerged fully to prevent air trapping from the top of the samples. The glass desiccators were used to prevent drying of the samples during the rewetting in the relatively dry lab rooms.

After saturation, the samples were transferred from the sampling rings into rubber membranes (Lexmond et al., 2024; Schindler et al., 2015). After transferring the samples, the dimensions of the samples were measured using a digital calliper to confirm their volume. The last step in sample preparation was to drill a small hole in the bottom of the samples using a small hand auger with a diameter of 5.3 mm and a depth of 12.5 mm to fit a tensiometer (HYPROP2, MeterGROUP). The hole was filled with water before inserting the tensiometer. After this final sample preparation step, the HYPROP2 was mounted on the samples to minimise handling of the samples. The HYPROP2 setup normally incorporates two tensiometers, but when using two, the tensiometers bend together due to the shrinkage of the samples. Therefore, the 12.5 mm tensiometer was used solely, while the other tensiometer opening was covered with tape. The top of the samples was uncovered just before the automated measurements were started.

Experimental setup

To characterise the shrinkage behaviour of the soil samples, the sample dimensions, weight and soil water suction were monitored over a prolonged period (up to 5 weeks) following the method proposed by Lexmond et al. (2024). The experimental setup was located in a climate-controlled room that was set to a temperature of 18°C and humidity of 60%. These conditions were maintained to minimise the evaporation rate and prevent desiccation cracking of the samples (Peron et al., 2009). Weight and soil water suction were monitored using the HYPROP2 setup (Meter Group, Germany), which also recorded the temperature. Sample dimensions were obtained with optical distance sensors with a beam diameter of 8 mm (Baumer FADK 14U4470/S14/IO), measuring the samples’ height (number of measurements [n] = 15) and radius (n = 2). The radius was measured 30 mm from the bottom of the sample (Lexmond et al., 2024). The sample volume was approximated assuming a truncated cone shape, calculated based on average height, top radius and bottom radius. The average height was determined based on the 15 height measurements. The volume was thus calculated, every 10 minutes, as follows:

where h is the average height (mm), r_high is the average radius (mm) at the top of the sample, and r_low is the average radius (mm) at the bottom of the sample.

Measurements were fully automated, taken every 10 minutes. Data collection was discontinued when no changes in average height or diameter were observed for 24 hours. The sample dimensions were measured using a digital calliper directly after the termination of shrinkage experiments. After confirming the final weight as well, the samples were resaturated for 2 months. After resaturation, the sample dimensions and weight were measured to determine the irreversible shrinkage that can occur. The samples (originally about 250 mL) were split into two parts: two-thirds for oven-drying to determine the water content and the other one-third for ‘wet’ analyses; the particle size distribution and Atterberg limits. Gravimetric water content was determined by oven-drying the large subsamples at 60°C for 24 hours (Dexter & Richard, 2009) after completing the shrinkage measurements. The VWC of the sample was derived from the measured weight changes combined with the density of water (temperature based) and the oven-dried weight of the sample (Eq. 2).

In which the sample’ bulk weight and the temperature of the sample were measured every 10 minutes. The ρ_water refers to the density of water and was defined as depending on the temperature (Eq. 3) and constant due to the climate-controlled conditions. The density of water was estimated using a polynomial approximation based on the formulation of Kell (1975):

Where T refers to the temperature in degrees Celsius. Thereafter, the dried subsamples were used to determine mineralogical composition and other analyses that are discussed below. A detailed overview of the measurement setup is described by Lexmond et al. (2024).

Sample characterisation

After the experiments were terminated, the GSD of the samples was determined by taking subsamples and removing cementing agents and organic matter from the samples, using respectively diluted hydrochloric acid (10%) and diluted hydrogen peroxide (30%). After cleaning the samples, the GSD was determined by laser diffraction using a Malvern Mastersizer 2000 particle sizer. The instrument settings, the detailed GSDs and a note on the uncertainty associated with optical particle sizing of sand-rich samples are provided in SI.2. These distributions are used here as supporting context only. The SOM and CaCO₃ contents were determined using thermogravimetric analysis (heating at a rate of 1 degree/minute up to 1000°C using the LECO TGA701), where SOM was set equal to weight loss between 105 and 550°C and CO₂ from CaCO₃ loss between 600 and 850°C. The liquid limit was determined using the fall cone technique and the plastic limit by the rolling technique at an external geotechnical laboratory (Geolab Wiertsema, The Netherlands). The clay mineralogical and bulk mineralogical compositions were determined using X-ray diffraction (XRD) (QMineral, Belgium). The CEC was determined based on the Co-Hexamine trichloride method (QMineral, Belgium). The XRD and CEC analyses were carried out for one sample per location instead of per sample.

Discretisation of soil layers

To be able to predict the soil layer thickness, the soil layers had to be discretised in order to incorporate the measured VWC and different sample characteristics. The discretisation of the unsaturated subsurface in simulation layers is defined based on the subsurface composition build-up, the placement of the VWC sensors and the depth of the upper two extensometer anchors. Two extensometer anchors define the boundaries of the combined simulation layers. The combined layers simulate the soil between depths of 0.06–0.07 m and 0.40–0.50 m for the three extensometer locations. For the simulation, the anchor-bound sublayers are split in two to five individual simulation layers, based on the homogeneity of the soil texture and the depth of the VWC sensors. If the subsurface build-up was relatively homogeneous, with relation to the texture, the simulation layers were solely based on the placement of the VWC sensors. However, differences in the characteristics of the soil layers yield potential jumps in the VWC gradient with depth. To understand if and where these jumps occur in the soil profile, HYDRUS 1D (Šimůnek et al., 2013) was used to simulate the VWC with depth for wet and dry periods, in order to define an accurate soil layer characterisation in relation to the measured VWC values. The simulations were based on the subsurface characteristics and average groundwater levels, under dry conditions. The simulation in- and outputs are shown in SI.3. The resulting discretisation is noted in the materials section.

The two prediction methods of the periodically unsaturated soil layer thickness change

Volumetric water content

The Shrink Test-Water Content model proposed by Briaud et al. (2003), referred to as the Briaud model in this work, is based on the relation between the volumetric- and vertical strain and the relation between the volumetric strain and change in VWC (Equations 4–6), as measured in the shrinkage experiments:

 

 

where ∆h_i is the difference in thickness compared to the initial thickness of soil layer i, which is defined at the maximum water content measured at the field site, f is the shrinkage ratio and E is the shrink-swell modulus. The vertical and volumetric strain are determined as the relative change in comparison to the dimensions at maximum field-VWC. The difference in the thickness of the soil layer is calculated for every available VWC measurement, which is every 15 minutes in this study. For this study, the method of Briaud et al. (2003) was applied incorporating a larger sample size (diameter 80 mm, height 50 mm).

The maximum and minimum VWC values were measured in the field using TDR and TDT VWC sensors (TEROS 12 and GroPoint Profile sensors). TEROS 12 (METERGroup, Germany) is a TDR sensor that measures the soil moisture based on a 1010 mL volume of soil, with a general accuracy of 0.03 m³/m³. The GroPoint Profile sensor (GroPoint, Canada) is a TDT sensor and measures the VWC at 15 cm intervals with an accuracy of 2–4%, depending on the calibration. The accuracy of both sensors depends on correct installation, meaning that the previously mentioned accuracies are the optimal accuracies of the sensors. Additionally, both sensor types function better when calibrated for the site and depth specific sediments. The TEROS 12 sensors were calibrated, but the GroPoint Profile sensor was not. The GroPoint sensor is difficult to calibrate due to its dimensions. The manufacturer provided specific relations between dielectric permittivity and VWC for this sensor. We have simulated the thickness change due to shrinkage and swelling based on both calibrated (calibrated TEROS 12 data for Site 1 and 2; the clay-rich specific mineral soil relation for GroPoint Profile Site 3) and non-calibrated (factory settings for both types of sensors and all locations) data within this work. The maximum VWC was used to determine the change in VWC per time step and the minimum was used to constrain the VWC range and define E and f within the proper range. The strain relations described by E and f (Eq. 5 and 6) were based on the shrinkage behaviour determined in the lab. Equations 4–6 were applied to the simulation layers defined on the measured VWC and textural classification. The used simulation layers and their thicknesses are noted down in ‘Site 1: Hazerswoude’ section to ‘Site 3: Nieuwolda’ section, the descriptions of the study sites. The sample properties are noted down in the ‘Site 1: Hazerswoude’ section to ‘Site 3: Nieuwolda’ section.

Soil water suction

The second prediction method relates the measured matric suction to the relative sample dimensions (the relative height and volume). The relation was defined as a polynomial function, by defining a polynomial trendline for suction and height and suction and volume (Eq. 7) for each sample. In a second step, a relation between the suction and sample weight, that were both measured with the HYPROP2 device (Eq. 8), was defined. That bulk density relation was used to determine the overburden pressure at the depth used for analysis, in relation to the depth and gravitational acceleration (Eq. 9). The relation between matric suction and relative sample dimensions was then applied to the resulting total pressure (summation of suction and overburden, Eq. 10). The matric suction was only measured for Sites 1 and 2, so this method was only carried out for these locations.

 

 

 

Testing the performance of the prediction methods

The performance of the prediction methods was tested to find the best choice to predict land surface movement induced by shrinkage and swelling. As the prediction models are based on shrinkage, we expect a better performance of the models during drying, simulating shrinkage and a decrease in soil thickness. Therefore, we tested the models on land surface movement change and on drying/shrinkage predictions separately. The predictions were tested by defining the root mean square error (RMSE) of the simulated change in soil layer thickness, totalling the change in all simulation layers per location in relation to the extensometer-measured change of thickness of the upper half a meter. RMSE values that are not exceeding 10% of the monitored changes are deemed satisfactory.

For Site 1 (Figure 1), there are several years of data available to base the predictions on and compare the predictions to. The simulation methods do not incorporate irreversible shrinkage and hysteresis between shrinkage and swelling. Therefore, the predictions for Site 1 have been corrected at the start of the yearly seasonal shrinkage to be able to assess the annual prediction accuracy, by eliminating the difference between simulated and measured thickness change at the start of the seasonal shrinkage.

Figure 1. Left: Map indicating the measurement- and sample locations in the Netherlands. Right: Lithological build-up of the Holocene sequence at the measurement- and sample locations. The average highest (upper 12%) and lowest (lower 12%) groundwater (GW) levels are indicated with the dashed lines.

We have studied the predictions using the Briaud model more in depth by assessing the coefficient of determination (R²) and applying a limited-memory bound-constrained Broyden-Fletcher-Goldfarb-Shanno (L-BFGS-B) multi-dimensional optimisation (R language) based on the first year and separate seasonal shrinkage and swell periods, by using the measured changes in VWC, the topsoil deformation and the Briaud formula (Equations 4–6). The f and E parameters were constrained between the values 0.09 and 50 to enable convergence, based on trial and error. The initial values for the f and E parameters were the measured values.

Study sites and sample characteristics

Three study sites were selected based on the subsurface build-up, specifically the presence of smectite minerals and the presence of extensometers that have been operating for several years (Van Asselen et al., 2024a, 2024b). This extended period of operation allows for the evaluation of long-term model performance, offering insights beyond isolated shrink-swell cycles and leading to a more comprehensive understanding of overall behaviour, including the identification of irreversible shrinkage. The sites are located in the north and west of the Netherlands, in areas rich in fluvial and marine clay deposits of Mid- to Late-Holocene age (Figure 1). The subsurface build-up of the sites exists of Holocene clay-rich sediments alternating with peat layers on top of a semi-stable Pleistocene sand layer (Van der Meulen et al., 2007) (Figure 1). The sites differ in depositional origin: Sites 1 and 3 were influenced by tidal marine processes, while the deposits in Site 2 were shaped by fluvial activity. Site 1 contains more silt and sand in the upper 2 metres compared to the other two sites. Sites 2 and 3 are rich in peat layers between 2 and 7–8 metres depth. The upper 50 cm of the subsurface of the sites consists of 10–50% clay-sized particles (Table 1), with a relative similar mineralogical composition dominated by smectite and vermiculite (Table 2).

Table 1. The sample characteristics determined after the shrinkage experiments.

Sample number and depth	Liquid limit [WC%]	Plastic limit [WC%]	Plasticity [%]	SOM [w%]	CaCO3 [w%]	CEC [meq/100 g]	GSD Clay/silt/sand [volume % of grains]	Saturated density of the sample [kg/L]
Site 1
Sample 1.1: 40 cm	33.5	22.9	10.6	2.42	0.083	9.7	8/43/49	1.90
Sample 1.2: 40 cm	31.2	22.0	9.22	2.35	0.11		11/47/42	1.88
Sample 1.3: 40 cm	32.6	20.5	12.1	2.35	0.11		10/44/46	1.95
Site 2
Sample 2.1: 30 cm	75.1	32.2	42.9	7.56	0.21	23.6	38/59/3	1.81
Site 3
Sample 3.1: 30 cm	80.5	34.5	46.0	8.44	0.25	30.9	46/51/3	1.70
The saturated density of the sample refers to the density at the onset of the drying experiments. WC[%] refers to the water content (%); w% refers to the percentage of the weight; meq is the milliequivalent. The clay size fraction has an upper limit of 8 μm, silt of 63 μm and sand of 1000 μm.

 

Table 2. The mineralogical composition of one sample per site, as determined by XRD-analysis for the bulk and clay size fraction.

Mineral – bulk	w% Site 1 Sample 1.3	w% Site 2 Sample 2.2*	w% Site 3 Sample 3.1
Quartz	61.3	32.2	17.1
K-feldspar	5.5	3.3	2.3
Plagioclase	9.3	6.2	2
Carbonates	0.6	13.4	0.4
2:1-layer silicates*	18.0	33.4	61.0
Kaolinite group	1.0	3.0	7.5
Chlorite group	0.7	2.1	1.3
*2:1-layer silicates in detail
Dioctahedral Muscovite	1.4	4.0	8.0
Dioctahedral disordered illite	1.7	4.5	5.4
Dioctahedral Illite/Smectite/Vermiculite	9.7	19.2	35.4
Dioctahedral Smectite/Vermiculite	5.2	5.8	12.2
Mineral – clay fraction
Kaolinite	2.4	2.1	2.6
Kaolinite- Smectite interstratified	4.2	6.9	6.9
Illite	9.3	13.6	9.2
Illite/Smectite/Vermiculite interstratified	54.5	58.5	60.2
Smectite/Vermiculite interstratified	29.4	17.6	20.8
Chlorite	0.3	1.3	0.4
* The drying and shrinkage measurement of Sample 2.2 is not included in this work, as the sample was extracted below the zone of interest at 105 cm depth at Site 2. w% refers to the mass percentage of the sample after removing organic matter. XRD: X-ray diffraction.

The extensometers monitor the change in subsurface levels relative to a reference point. The reference anchor is an iron massive cone founded in the semi-stable Pleistocene sand underlying the Holocene deposits (Van Asselen et al., 2020) for all three sites. The sand layer is considered stable for the period of interest: several years (Van Asselen et al., 2020). Borros anchors were installed at several depths between the land surface and reference points, from which the level relative to the reference point was monitored at a frequency of 15 minutes. The shallowest anchor is a square perforated stainless-steel plate (0.5 × 0.5 m). The anchor at 0.4 m depth is a horizontal steel strip. The displacement sensors have an accuracy of 0.1 mm, and the anchor can be installed with a depth precision of 3 mm and the risk of poor anchor borehole-wall coupling is low (Kooi & Waardenberg, 2024). In addition to the changes in the anchor levels, the VWC, suction, temperature and electric conductivity of the porewater were also measured in the upper 50 cm of the subsurface at the sites. Although factors such as precipitation, solar radiation, phreatic and confined groundwater heads, and nearby ditch water levels were also measured at the sites, they are not considered in this analysis as the VWC was directly observed.

The setup of the extensometers was similar for the sites of interest: the surface anchors were situated a few centimetres below the actual surface, where the second anchor was placed at depths of 40–50 cm (Figure 1). The VWC, matric suction and soil temperature were monitored at 20, 35 and 50 cm below the land surface at Sites 1 and 2 (TEROS 12; Meter group), whereas Site 3 has a GroPoint Profile sensor installed, that measured the average VWC at 20–35 and 35–50 cm below the land surface. The matric suction (TEROS 32, Meter group) was also measured at depths of 20, 35 and 50 cm below the surface at Sites 1 and 2 (Van Asselen et al., 2024a).

Site 1: Hazerswoude

This site is located in the west of the Netherlands, in the Hazerswoude-Dorp polder, where the subsurface is characterised by clay-rich shallow marine (tidal) deposits overlying a Pleistocene sand layer. The surfacing layer (0–22 cm) consists of silty clay with 23.2% clay-sized particles, underlain by a loam layer containing 16.5% clay-sized particles, with a thickness of 113 cm (Van Asselen et al., 2024a). The soil is classified as entisol and the groundwater level is situated 94 cm below the surface on average. A comprehensive subsurface description is provided by Van Asselen et al. (2024a). The site, previously utilised as arable land, transitioned to grassland after the installation of the extensometer in 2021.

The two upper anchors were installed at 6 and 50 cm depth. To study the shrinkage and swelling behaviour at this site, three samples were extracted at a depth of 35–40 cm below the land surface. The sample was extracted from this depth as the upper 40 cm of the subsurface was considered to have a similar GSD. The sample was extracted from the deepest segment of this layer to avoid the dried-out conditions at the top of the layer that would not allow for undistributed sampling. The plasticity of the samples was low (Table 1), and the samples contained 2.4% organic material and 0.1% calcium carbonate. The mineralogy of the site was analysed by X-Ray diffraction. The bulk of the sample consists of quartz, K-feldspar and plagioclase (80.3%), with 18.0% of the bulk of the sample exists of 2:1-layer silicates (Table 2). The clay-sized fraction of the sample is predominantly composed of interstratified minerals containing smectite and vermiculite minerals (88.1%) (Table 2). The interstratified smectite-vermiculite has the highest expansive potential and makes up 29.4% of the clay-sized mineral content. Only 2.6% of the clay-sized mineral fraction exists of non-expansive minerals (kaolinite and chlorite). Illite has a lower expansive potential than smectite and vermiculite and made up 9.3% of the total clay-sized fraction.

The surface layer is discretised in five layers to predict the thickness change based on VWC measurements. The first layer is the clay-rich top layer (6–22 cm). This layer has been defined based on the limited gradient in the VWC that was modelled using Hydrus 1D (Šimůnek et al., 2013). Between 22 and 50 cm depth, a linear VWC gradient was simulated. Therefore, this layer was discretised in four layers of 7 cm. The VWC values measured at 35 and 50 cm are directly applied to the 29–36 cm and 43–50 cm layer. The average of the measured VWC at 35 and 50 cm was applied to 36–43 cm. The trend in VWC (declining or increasing with depth) between 35 and 50 cm depth was used to determine the VWC value for 22–29 cm, by assuming a continuation of this trend between 22 and 29 cm.

Site 2: Cabauw

This site is situated in the west of the Netherlands near the village of Cabauw in an area historically shaped by fluvial processes (Van Asselen et al., 2024a). The shallow subsurface consists of fluvial clay containing varying amounts of silt. The uppermost 15 cm contains 52% clay-sized particles, and clay content increases to 60.7% at a depth of 40–50 cm. The soil is classified as entisol. The Pleistocene sand layer is situated about 8.5 m below the surface. The groundwater level is situated at 62 cm below land surface on average. The current and long-term land use of the site is pasture. However, the area has been fenced off since the installation of the measurement site in 2021 to prevent grazing. The VWC and matric suction were measured since the end of 2023. The upper two anchors were installed at 6 and 40 cm depth. The VWC sensors are installed at 20, 35 and 50 cm depth at this site too.

To predict the land surface movement at this location, a sample was extracted from a depth 30–35 cm below the surface. This depth was selected as the SOM was comparable to the SOM deeper in the subsurface. The sample’s plasticity was high, and the sample contained 7.6% of organic material and 0.2% of calcium carbonate. In contrast to the mineralogy of Site 1, most of the bulk did exist of 2:1-layer silicates (33.4%), next to 32.2% quartz (Table 2). The expansive mineral content of the clay-sized fraction is very similar to Site 1 with 96.6%. However, the interstratified smectite-vermiculite content is significantly lower at 17.6%. The clay-sized fraction exists of 13.6% illite and 83% expansive interstratified clay minerals (Table 2). Only one sample was used to determine the shrinkage behaviour for practical reasons.

The composition of the top layer is relatively homogeneous, consisting of silty clay. The SOM content decreases with depth from the surface to a depth of 26 cm. Between 26 and 40 cm, the clay content varies, increasing with depth. Due to the relatively minor differences in the characteristics of this subsurface layer with depth, the surface layer is defined to consist of two model layers: 6–23 cm and 23–40 cm.

Site 3: Nieuwolda

This site is situated in the north of the Netherlands. The shallow subsurface is composed of a clay-rich tidal deposit. This 2-metre-thick deposit is situated on top of peat layers (alternating with clay and about 4 metres thick [Figure 1]) and the Pleistocene sand layer. The clay content of the upper 50 cm of the subsurface is about 42%. The groundwater level is situated about 150 cm below the land surface year-round, aside from short-time raises. The VWC was measured using GroPoint Profile TDT sensor (GroPoint, Canada) that has not been calibrated at depths 20–35 cm and 35–50 cm for the surface layer.

The two upper anchors were installed at 6 and 50 cm depth. Site 3 was sampled at a depth of 30–35 cm. The sample was extracted from this depth, as it was the only interval without an abundance of reed fragments that would result in cracking of the sample during the drying experiments. Cracks in the sample cannot be measured with the measurement set-up. The sample had a very high plasticity and contained 8.4% organic matter and 0.3% calcium carbonate (Table 1). The bulk sample contained mostly 2:1-layer silicates (61.0%) and only 17.1% quartz (Table 2). The expansive mineral content of the clay-sized fraction is very similar to the other two sites with 97.1%. The interstratified smectite-vermiculite content is in between that of the other two sites with 20.8% (Table 2). The clay-sized fraction exists of 9.2% illite and 87.9% expansive interstratified clay minerals. Only one sample was collected due to the aforementioned high amount of reed in the subsurface.

The texture of the shallow subsurface of Site 3 was homogeneous with depth and was silty clay. Only the colour of the material was found to vary slightly within the upper 50 cm. Therefore, the discretisation of the measured upper layer was solely based on the placement of the VWC sensor. Two layers were defined: 6–28 cm and 28–50 cm. The measured VWC was directly used as simulation value.

Results

This section is organised per data type: ‘Field data’ section focusses on the extensometer data to show the changes in the thickness of the Holocene deposits for the three locations. ‘Shrinkage behaviour’ section shows the measured shrinkage behaviour in soil shrinkage curves, including a summary of the parameters that are used in the Briaud et al. (2003) model, and the relations between the matric suction and relative dimension change of the samples used in the alternative method. ‘Predictions of shrinkage and swelling in-situ’ section shows the predictions based on the modelling methods, and ‘Performance of predictions’ section summarises the performance of the models.

Field data

Two trends are visible in all sites: there is a strong seasonal surface movement, and most of this vertical deformation occurs in the upper 40–50 cm for all sites (Figure 2). Permanent deformation (long-term subsidence) is foremost visible in Site 1 and to a lesser extent in Site 2. At Site 1, the maximum vertical deformation (maximum thickness – minimum thickness) of the subsurface layers diminishes with depth, when the deformation is normalised over the thickness of the subsurface. The top layer caused 36% of the total deformation of 14 mm, from February 2022 to September 2022, 48% of the maximum total deformation during 2023 (January to July) and 42% of the total maximum deformation during 2024 (January to August). Contrastingly, the thickness change of the bottom layer (650–900 cm) is in the range of 10^-3 mm. Every summer, the maximum thickness of the top layer (6–50 cm depth) decreases with 1–2 mm. Whereas the maximum thickness of the total Holocene subsurface decreased with 3 mm between 2023 and 2024 (Figure 2 Site 1).

Figure 2. The vertical deformation, defined as the thickness change in relation to the initial (time of installing) thickness, of the monitored subsurface layers for the three sites. The grey zones in the plots indicate the availability of measured volumetric water content data and thus the periods that are simulated in this work. Note the differences in the vertical axes for the different sites. From May 2022 to May 2023, the third anchor installed at -260 cm was malfunctioning at Site 1; this data is not included in this plot.

The total deformation measured at Site 2 is larger than at Site 1. The maximum deformation of the measured sequence (6–1050 cm) is 34 mm during 2022 (February to September), 27 mm during 2023 (January to July) and 17 mm during 2024 (February to August). The top layer (6–40 cm) shows a maximum deformation equalling 53–56% of the maximum total deformation during a year (Figure 2 Site 2). For Site 2, the vertical deformation of the subsurface layers diminishes with depth, again highlighting the fact that most land surface movement stems from the upper 40–50 cm of the subsurface column. There is an irreversible 1–2 mm decrease of the top layer (6–40 cm) per summer.

The vertical deformation of the shallow subsurface at Site 3 is in the same range as Site 1, with a maximum deformation of 8 mm during 2023 (January to June) and 15 mm during 2024 (July–September). The maximum deformation of the top layer (7–50 cm) accounts for 79–103% of the total deformation. During 2024, the top layer shrinks more than the maximum total strain, because the thickness of the peat-rich layers (195–780 cm) increased. The peat-rich layer increases in thickness during the start of every summer and starts to decrease during September/October, resulting in smaller total deformation over time. The increase in the thickness of the peat-rich layers is a result of the water management at the site: the groundwater levels are increased in the summer time via the surrounding ditches. The extensometer was installed later in the year than at the other two locations. Combined with the water management, this could explain the increase in the subsurface thickness during the period of May 2022–September 2022.

The deformation of the shallowest subsurface layer at Site 1 agrees well with the variation in the measured VWC in general (Figure 3). However, short time changes in the VWC (less than a week) do not seem to be related to thickness changes of the shallow subsurface layer. The VWC decreases with depth in the shallow subsurface layer, as a result of the low hydraulic conductivity of the silty clay overlying loam with higher conductivity. At Site 2, the deformation of the shallow subsurface layer exceeds that of Site 1, but a lower variation in VWC was measured. During the summer of 2024, the thickness of the shallow subsurface layer at Site 2 decreases at the same time as the VWC decreases.

Figure 3. The measured volumetric water content (VWC) at Site 1 (depth VWC 1 = 20 cm; VWC 2 = 35 cm; VWC 3 = 50 cm), Site 2 (depth VWC 1 = 20 cm; VWC 2 = 35 cm) and Site 3 (depth VWC 1 = 20–35 cm; VWC 2 = 35–50 cm) (left vertical axes), in combination with the measured thickness change of the upper subsurface layers (right vertical axes).

The deformation at Site 3 shows a short time increase in thickness during the winter of 2022–2023, due to a cold day combined with precipitation (1.4 mm), when the temperature was 0.1°C 10 cm above the surface at a local measurement station (KNMI, Nieuw Beerta, approximately 10 km from Site 3) (SI. Figure SI.3). During summer, the VWC decreases, as does the thickness of the shallow subsurface layer, but the latter delayed. The general variation in VWC is faster and has a higher amplitude than at the other two sites. Also, the limit of the VWC sensors at this site is reached during the spring of 2023. The relation between the VWC variation and change in shallow subsurface thickness is least pronounced at Site 3.

The matric suction was also measured at sites 1 and 2 at similar depths as the VWC. Air entered the shallowest tensiometer at Site 1 and the tensiometer at 35 cm depth at Site 2, yielding some gaps in the data. For Site 2, an interpolation was used based on the other two tensiometers to get an insight into the matric suction at a depth of 35 cm (Figure 4 & SI Figure S1). The layer deformation at these two sites fits the changes in matric suction over time. Especially, the change in matric suction measured at a depth of 20 cm coincides with the layer deformation.

Figure 4. The measured matric suction at Site 1 (depth matric suction 1 = 20 cm; suction = 35 cm; suction 3 = 50 cm), Site 2 (depth suction 1 = 20 cm; suction 2 = 35 cm) (left vertical axes), in combination with the measured thickness change of the upper subsurface layers (right vertical axes).

Shrinkage behaviour

This section describes the shrinkage behaviour of the samples extracted at the measurement sites. The shrinkage behaviour is defined as relative dimensions in relation to the VWC, instead of soil shrinkage characteristics curves, as the relative changes are applied in the prediction models, based on the VWC. The shrinkage of the samples extracted at Site 1 was limited, with a maximum volume loss of 5% and a height loss of maximum 2.5% (Figure 5). The initial water content of the samples (volume water in sample/sample volume) is close to 0.45 (mL/mL) on average and decreased to 0.10. All three samples from Site 1 show a relatively fast initial shrinkage that decreases over time.

Figure 5. Top row: the relative volume (sample volume/initial sample volume) in relation to the volumetric water content (volume water in sample/sample volume) of the samples. Bottom row: the relative sample height (sample height/initial sample height) at decreasing volumetric water content. Note the difference in the y-axis for top row, Site 2.

The sample extracted at Site 2 had a slightly higher VWC than the samples from Site 1 (0.51). The sample was air-dried to a VWC of 0.12. The sample from Site 2 shows the largest decrease in volume and height of all samples during drying: 36% (volume) and 14% (height) of the initial dimensions were lost during drying.

The clay-rich sample from Site 3 shows a volumetric shrinkage of 5% and a vertical shrinkage of 2%. Most of the shrinkage from this sample was thus radial shrinkage. The initial VWC was slightly higher than for the other sites at 0.58, but it did decrease to 0.10. The maximum VWC measured in these experiments is lower than the maximum VWC that was measured at the field sites.

The shrinkage behaviour was also defined with the E and f factor, according to the methodology proposed by Briaud et al. (2003). E is derived from the relation between volumetric strain and change in VWC and f is the ratio between vertical and volumetric strain (Equations 3 and 4). The defined E and f factors are shown in Table 3. The ratio of f to E is determining the amount of strain that is expected during drying or wetting of the subsurface. The ratio between f and E for Site 2 is a factor 10 larger than for the other sites and samples, indicating that the expected strain is about 10 times higher with change in water content for this location, compared to the others.

Table 3. The E (volumetric strain in relation to volumetric water content change) and f (ratio between vertical strain and volumetric strain) factors per sample.

Briaud parameters	Site 1: Sample 1.1	Site 1: Sample 1.2	Site 1: Sample 1.3	Site 2: Sample 2.1	Site 3: Sample 3.1
E [-]	6.32	6.21	9.22	0.558	4.51
f [-]	0.309	0.315	0.566	0.317	0.254
Ratio f/E [-]	0.0489	0.0507	0.0614	0.567	0.0564

The relations between the matric suction and relative sample dimensions were defined between matric suction 0 to -100 kPa (Figure 6). The relations that are indicated in Figure 6 are described with polynomial trendlines (SI: Table SI.1) to predict the layer deformation later on. The trendlines had 20 decimal numbers to ensure the proper relation was simulated. For some samples, two trendlines were defined to be able to capture the full shrinkage behaviour in relation to matric suction. Relatively more strain occurred in the radial direction, compared to the vertical for all sites (Figure 6).

Figure 6. Top row: the relative volume (sample volume/initial sample volume) in relation to the matric suction of the samples. Bottom row: the relative sample height (sample height/initial sample height) in relation to the matric suction of the samples. The absolute value of the matric suction is noted in these graphs.

Predictions of shrinkage and swelling in-situ

The relations between strain and water content presented in the previous sections have been applied to the measured water content at the three sites to predict shrinkage and swelling. The fit of the Briaud-based predictions and measured deformation differs per site (Figure 7). Site 1 shows a lot of similarities between the predicted and measured deformation, but the short-term deformation is slightly overestimated by the predictions. There is also little difference between the predictions based on calibrated and non-calibrated VWC data. The best fit for the seasonal behaviour is based on Sample 1.3. However, swell is overestimated for this location, creating a larger cumulative offset with every wet period.

Figure 7. The change in the thickness of the top layer measured using an extensometer at the three different locations in combination with the Briaud-based thickness predictions where nc indicates the predictions based on non-calibrated volumetric water content (VWC) data and c the calibrated version. At Site 1, the thickness of the top layer, i.e. the thickness of the layer between the first and the second anchor below the surface at the moment of installation of the extensometers, is 44 cm, at Site 2 34 cm and at Site 3, 43 cm. The nc data series refer to the simulated changes based on uncalibrated VWC measurements, whereas n data series refer to simulations based on calibrated VWC measurements. The numbers refer to the sample numbers.

The Briaud-based predictions for Site 2 fit the measured data less than for Site 1 (Figure 7). The swell that occurs during the 2023–2024 winter is not simulated. The simulated seasonal negative strain, starting in the summer of 2024, has the same dimensions as what is measured. The positive deformation that follows during the wet autumn is also in line with the measured deformation. However, a highly dynamic subsurface deformation is simulated, but not measured. Also, there is a large difference in outcome based on either the calibrated and non-calibrated VWC input.

The Briaud models for Site 3 do not fit the measured data well (Figure 7), both for the short- and long-term behaviour. Most striking is the positive deformation that is measured in December 2022. The swell occurred after precipitation and some days with temperatures below 0°C. Thereafter, the measured deformation is slightly negative, whereas the Briaud models simulate a positive deformation during this period. The simulated and measured negative deformation do not align well temporally. From June 2023 on, the positive deformation aligns better temporally and dimensionally.

The newly proposed applications of the shrinkage curve yield different results with regard to the fit between predicted and measured layer thickness changes (Figure 8 compared to Figure 7). The prediction based on the relative height and volume in relation to the matric suction do show negative and positive deformation with similar timing as the Briaud-based predictions, but an increase in the amplitude of the short-term deformation for Site 1 (Figure 8), compared to the Briaud-based predictions. For Site 2, there is little similarity between simulated and measured deformation dynamics. For both sites, the best fit for the seasonal shrinkage is the simulation based on volumetric shrinkage, including a simple overburden pressure term. Whereas, the short-term behaviour is simulated best with the simulation based on only vertical deformation in combination with overburden pressure. Swell is overestimated by every simulation based on the relation between matric suction and dimensional changes.

Figure 8. The change in the thickness of the top layer measured using an extensometer at the three different locations in combination with the matric suction-relative dimension thickness predictions. In this, V indicates the relations between matric suction and volumetric deformation and H the relations between matric suction and vertical deformation. At Site 1, the initial thickness of the top layer was 44 cm, and at Site 2 34 cm.

Performance of predictions

The performance of the predictions is assessed based on the difference between predicted and measured thickness for each site, the RMSE, and the R² of the predicted thickness in relation to the measured seasonal shrinkage. The RMSE of the Briaud models differs an order of magnitude between the three sites. For Site 1 the average RMSE of the Briaud-based prediction is 0.472 mm, based on a period of 3 years, where the average maximum shrinkage was 4.1 mm (Table 4). The maximum shrinkage is defined as the maximum minus minimum thickness of the top layer during the dry season. When the predictions for this site based on the Briaud model were corrected for the overestimation of seasonal swell, the RMSE decreased to 0.102 mm. Thus, the RMSE of the Briaud-based predictions of Site 1 is between 10 and 2.5% of the average seasonal measured shrinkage. The yearly-fitted data underestimates shrinkage with an average of 0.88 mm, which is about 20% of the average seasonal shrinkage. The R² of the Briaud predictions was 0.69 when assessing the full simulation period. The maximum shrinkage simulated with the Briaud model was very close to the measured maximum shrinkage for the studied Site 1.

Table 4. Summary of model performance determined based on the root mean square error.

Model	Site 1	Site 2	Site 3
Measured maximum shrinkage (average of 3 years for Site 1 only)	4.1 mm	9.0 mm	5.9 mm
Briaud-simulated maximum shrinkage	4 mm	10 mm	6 mm
Vertical-shrinkage + overburden matric suction model maximum shrinkage	1.8 mm	2 mm	-
Volumetric-shrinkage + overburden matric suction model maximum shrinkage	7 mm	6 mm	-
RMSE Briaud model	0.472 mm	1.12 mm	0.692 mm
RMSE Vertical-shrinkage + overburden matric suction model	10.9 mm (full simulation range) 1.84 mm (year 1)	4.86 mm	-
RMSE Volumetric-shrinkage + overburden matric suction model	11.7 mm (full simulation range) 2.19 mm (year 1)	4.64 mm	-
Root Mean Square Error.

The RMSE for Site 2 is about two to three times larger than for Site 1 (1.12 mm), corresponding to a maximum shrinkage of 9.0 mm (June–September 2024). The RMSE is thus 12.4% of the seasonal change, which is above the 10% limit. Yet, the predicted and measured thickness align very well at the end of the wet period (difference of 0.01 mm) and the R² was only slightly lower than for Site 1 at 0.64.

The RMSE of the Briaud-based predictions of Site 3 was 0.692 mm, where the measured maximum shrinkage was 5.9 mm. This RMSE was thus 12% of the measured seasonal shrinkage. The predicted shrinkage seemed to be delayed compared to the measured shrinkage and overestimated by 0.1 mm (1.7% of the measured shrinkage). The R² value of 0.04 indicated a lack of correlation between simulated and measured deformation for this site.

The similarity between the Briaud-based modelled and measured deformation dynamics are indicative that most of the deformation results from changes in water content. However, the newly proposed method relating suction and overburden pressure to shrinkage and swelling, yielded higher RMSE values and thus indicating less of a fit to the measured behaviour. This high RMSE is partially explained by the overestimation of swell in the predictions for Site 1. As, the offset between the simulation and the measured deformation increases over time. The results do highlight a better fit for the vertical-shrinkage prediction for Site 1 and the volumetric-shrinkage prediction for Site 2. Also, the maximum predicted shrinkage is closer to the measured maximum shrinkage, based on the vertical-shrinkage for Site 1 and volumetric-shrinkage for Site 2.

The RMSE is a valuable tool for assessing model performance, but temporal trends, for instance due to an overestimation of swell, do interfere with the assessment of the short-term quality of the predictions. Therefore, we have also tested the quality of the fit when the Briaud model was optimised based on the in-situ measured data. In other words, we have tested what the best f and E parameters would be to simulate the top layer deformation, by optimising the f and E parameters based on the measured data (Figure 9). The optimised versions of the Briaud model show a high-quality fit to the measured top layer deformation for the year 2022. The swell is overestimated from September 2022 on, but overall, the simulation fits the measurements well. The overestimation of swell behaviour was also found in the non-optimised versions of the Briaud models and is, therefore, not related to the E and f parameters. The optimised f/E ratio is also very close to the f/E ratio determined in the shrinkage experiments (Table 5). The f/E ratios that are noted down in Table 5 are averages based on the different simulation layers that make up the total simulation surface layer. The f/E ratio for the upper simulation layer (6–22 cm) was 10 (2022) to 100 (2023 and 2024) times higher than the other and lower four simulation layers, indicating that most deformation occurs in the upper 20 cm of the shallow subsurface layer.

Figure 9. Left: the simulated layer deformation at Site 1, for year 1, 2022, based on optimised f and E parameters. Right: the optimised simulation for Site 2 for the shrinkage period. The values of the parameters are summarised in Table 5.

Table 5. The measured and optimised E and f parameters per optimisation for all three sites.

Scenario	E [-]	f [-]	f/E [-]
Measured S1.1	6.32	0.309	0.0489
Optimised S1.1; year 1	31.3	0.457	0.0501
Optimised S1.1; shrink 1	6.31	0.369	0.0588
Optimised S1.1; shrink 2	25.6	0.295	0.0391
Optimised S1.1; shrink 3	31.5	0.344	0.0462
Measured S1.2	6.21	0.315	0.0507
Optimised S1.2; year 1	27.0	0.362	0.0503
Optimised S1.2; shrink 1	6.20	0.361	0.0585
Optimised S1.2; shrink 2	29.5	0.289	0.0388
Optimised S1.2; shrink 3	37.9	0.227	0.0263
Measured S1.3	9.22	0.566	0.0614
Optimised S1.3; year 1	29.0	0.513	0.0502
Optimised S1.3; shrink 1	9.21	0.522	0.0570
Optimised S1.3; shrink 2	30.8	0.402	0.0388
Optimised S1.3; shrink 3	32.8	0.476	0.0462
Measured S2.1	0.558	0.317	0.567
Optimised S2.1; full period	Boundaries of optimisation parameters 50.0 and 0.09 for E and f
Optimised S2.1; shrinkage period	E and f not separately determined (no convergence)		0.33
Measured S3.1	4.51	0.254	0.0564
Optimised S3.1; full period	Did not converge
Optimised S3.1; shrinkage period

For Site 2, the optimalisation of the f and E factors was carried out similarly for the full period, yielding f and E parameters that had the values of the optimalisation boundaries (Table 5). The optimised f/E ratio was well below the experimentally defined one (Table 5). The simulated layer deformation stemming from the optimisation shows behaviour very different from what is measured (Figure 9) and is off by about 1 mm at the starting point of the simulation/measurement period. Optimising the Briaud model for Site 3 was not possible, as the model was not able to converge.

For Site 1, the simulation was optimised per seasonal shrinkage and swelling period as well (Figure 10). The f/E ratios that are used are noted down in the SI Table 2. When applying this seasonal optimisation, the periods of short deformation are still overestimated, but it is possible to render the seasonal changes properly. To do so, the f/E ratio for the swell periods has to diminish every year. The average f/E ratio for the first swell period (2022) is 0.033 and decreases to 0.21 and 0.17 for respectively the second (2023) and third (2024) swell periods. The f/E ratio for the seasonal shrinkage phases does not show a trend.

Figure 10. Comparison between the Briaud model results optimised seasonally and measurements for Site 1, Hazerswoude.

Using the above optimisation, the predicted swell compared better to the measured swell. However, the shrinkage simulated in 2023 and 2024 does not fit the measured shrinkage better, compared to the predictions based on the lab measurements. The simulation of the seasonal shrinkage works adequately without optimalisation. The same lab-based simulation can be used for several years when the overestimation of swell is corrected for at the start of the seasonal shrinkage period (Figure 11). Roughly double as much swell is predicted compared to measured swell for Site 1. For the other two sites this comparison is not possible due to lack of data. For Site 2, there is no indication that swell is overestimated.

Figure 11. The predictions of layer deformation for Site 1 using the Briaud model and corrected at the onset of seasonal shrinkage.

Discussion

The goal of this study was to predict vertical land surface movement due to soil shrinkage and swelling based on the VWC (Briaud et al., 2003) and matric suction for three different locations in the Netherlands. The Briaud model relates the change in VWC to the relative vertical and volumetric strain. In addition to this method, we have related the relative vertical and volumetric strain to matric suction and predicted the shallow deformation for two measurement sites. The performance of the models was quantified using the RMSE value and the Briaud predictions by running optimised versions of the simulation.

In general, the predictions of land surface movement defined as layer deformation, based on in-situ measured VWC (Briaud et al., 2003), performed adequately with regard to predicting soil deformation due to shrinkage, during a period of prolonged lower water content in the topsoil layer (up to a depth of 50 cm). The seasonal swell was also predicted accurately for Site 2, but overestimated for Site 1 and Site 3. The method performs relatively well on short-term behaviour for Site 1, but not for the other sites (based on the RMSE values). The direct application of relative sample dimensions based on matric suction highlights the complexity of the shrinkage behaviour: the ratio between vertical and volumetric shrinkage seems to depend on the duration or intensity of a dry period (high matric suction). The predictions improve when the overburden pressure is added to the matric suction, before applying the relation determined in the lab to predict the shallow deformation. The results are discussed per site first and integrated thereafter.

Site 1: Hazerswoude

The shrinkage measured in the laboratory was completely reversible for the samples extracted at Site 1. Contrastingly, the measured in-situ swell of the top layer (6–50 cm) is less than predicted based on these swell measurements. This top layer lost between 1.25 and 2 mm of thickness per year for the period 2022–2025. The decrease in swelling potential and the progressive loss of surface layer thickness could be the result of a relatively dry summer, such as in 2018 just preceding the measurement period. If the surface layer was dried below the shrinkage limit at that time, the swelling potential has been increased directly afterwards, due to entrapment of air and restructuring of the matrix (Basma, 1996; Estabragh et al., 2018). The shrinkage and swelling will equilibrate in the years following the dry event. However, this explanation is not likely, as the samples were extracted during 2023 well after 2018, and especially as the loss of the layer thickness is constant per year. Therefore, the overestimation of the swell behaviour in the simulations must be the result of the different conditions between the measurement setup in the laboratory and in the field, such as lateral restraining or overburden pressure, or disturbance of the samples during extraction. With regard to the latter, the deformation of the pore space is known to affect shrinkage less than swelling behaviour (Makki et al., 2024). Also, swell is a slower process than shrinkage, due to the smaller pore sizes in a dry soil that is being wetted and the hydration forces (Bronswijk, 1991; Kodikara et al., 1999; Makki et al., 2024). This could explain the lack of fully reversible shrinkage, but the wet period during the winter of 2023–2024 does not show swell, whereas the VWC values remain high for an extended period. These factors seem to point in the direction that either the swell pressure is too low (due to the relative low amount of expansive minerals) to increase the layer thickness in relation to the overburden pressure, or the constraining lateral force, that results from the surrounding subsurface material, limits the swell (Kodikara et al., 1999).

The shrinkage behaviour was well simulated based on the Briaud model and the VWC values, independent of the overestimation of the swell. A simple correction consisting of subtracting the overestimated swell, results in a near-perfect fit for the seasonal shrinkage simulation (Figure 11). The shrinkage behaviour was not affected by the irreversible volume loss.

The measured matric suction dynamics at Site 1 shows a clear seasonal trend that fits the surface layer deformation (Figure 4). The predictions based on the relation between the vertical deformation and matric suction, show the same build-up in offset as the VWC-based simulations and overestimate the seasonal swell (Figure 8). The addition of an overburden term in the simulations prevents most of the swell overestimations. The short-term shrinkage and swelling are predicted best with the vertical-deformation model in this study, whereas the seasonal shrinkage is predicted better with the volumetric-deformation model for 2022 and 2023. For 2024, the vertical-deformation model is also better suited for the seasonal shrinkage prediction. During the wet year 2024, the matric suction decreased mainly in the upper 20 cm and decreased with only 15 kPa and 5 kPa at 35 and 50 cm depth, respectively. This is significantly less than in the previous years. When soils dry slowly, and the soil water suction gradient is low, no desiccation cracks are formed (Peron et al., 2009). The results presented in this work indicate that the shrinkage is mostly directed vertically in these circumstances, especially when the drying does not only occur at the land surface, but also at 35 cm depth.

Site 2: Cabauw

The predicted surface layer deformation (6–40 cm) with both methods does not agree well with the measured deformation of Site 2 (Figures 7 and 8). The Briaud models show a maximum thickness change that agrees with the measured maximum thickness change. Yet, the short-term surface layer deformation at this site is not described well based on the VWC data. The VWC values are relatively high (between 40 and 50% most of the time, Figure 2). These values are outside the calibrated domain of the VWC sensor. Due to the high content of expansive clay minerals at Site 2, the VWC values and thus changes in the VWC values might be overestimated in this wet and not-calibrated domain of the VWC sensors, as the dielectric constant that is used as measurement proxy is influenced by the amount of expansive clays’ interaction with the water. This could partly explain the mismatch between the measured and Briaud-based simulated deformation dynamics.

The matric suction at Site 2 was monitored from November 2023 until the beginning of 2025, which was a relatively wet period (Figure SI.1), with relatively low matric suction. The matric suction was around 0 kPa until July 2024, indicating a saturated soil. For the same timeframe, the VWC values indicate a saturated soil too. The groundwater level was also monitored at the measurement site (SI: Figure SI.4) and shows a similar result: the subsurface layer of interest (initially 6–40 cm depth) is partly to almost fully saturated during the simulation period (Figure 12). When the groundwater level rises, swelling occurs, due to a decrease in effective stress and poro-elastic response within the shallow subsurface layer. Every time the groundwater level is above 20 cm depth, there is a slight increase in the thickness of the shallow subsurface layer (Figure 12). This highlights the threshold between the normal shrinkage-swelling range and swell due to increased pore pressure. Additionally, it is possible that the water was migrating from macro- and mesopores to micropores and the interlayer domain at these saturated conditions, resulting in swell/positive layer deformation. The simulation is not taking the redistribution of water in the soil into account and therefore the positive layer deformation was not simulated. From July 2024 on, the matric suction increased significantly at all measured depths and the magnitude and timing of the shrinkage during this period was simulated correctly (Figure 8).

Figure 12. The groundwater level measured at Site 2: Cabauw, in combination with the measured layer deformation of the initial 6–40 cm depth layer.

Site 3: Nieuwolda

The simulation of surface layer (7–50 cm) deformation based on changes in VWC for Site 3 does not agree to the measured deformation (Figure 7). This may be explained by the fact that the in-situ VWC measurement equipment was not calibrated. Because the matric suction was not monitored at this site it is very hard to assess the quality of the measured VWC data, due to the non-linearity in the relation between the dielectric constant and water content resulting from the interaction between the water and expansive clay minerals. Site 3 is very rich in expansive clay minerals (>40 w%), which influence the dielectric constant due to their interaction with water. The ‘plateau’ at a VWC of 0.66 in March 2023 (Figure 3) also indicates that the VWC values are at the sensor operational limits. The sensors are best fit to measure at lower VWC values and their accuracy diminishes when the VWC is high (GroPoint, 2025). The measured surface layer deformation shows a striking pattern that seems to be related to freezing events: during December 2022 and January 2024 there is an increase in the surface layer thickness of 2–4 mm (Figure 2) for a few days and thereafter the thickness of the surface layer slowly diminishes again. These short bursts of increased surface layer thickness occur during days when the temperature becomes 1°C or lower at a depth of 20–35 cm (SI Figure SI.3) and freezing occurred at the land surface. A temperature of 1°C or lower was not measured below the surface at the other two sites. Freezing and thawing cycles create microcracks in the soil structure, resulting in a decrease of shrinkage (Zhao et al., 2021), which also could partly explain the lower thickness changes that have been measured (3 mm swell), compared to the predicted (1 mm shrinkage) changes (Figure 7).

Integration: performance of the tested models

This section aims to summarise and integrate the results discussed previously. The importance of proper VWC sensor calibration, the underestimation of some shrinkage events and swell constraints in relation to matrix composition and interference with other subsurface processes are discussed.

The performance of the tested models is not uniform among the three sites of interest. The lack of uniform performance is due the different soil characteristics and groundwater dynamics at the sites and probably due to the calibration of the in-situ VWC measurement equipment at the sites. The measurement setup at Sites 1 and 2 is similar, but the calibration of the VWC sensors had to be carried out for different VWC ranges. For Site 2, the calibration was not carried out for the full in-situ VWC range. The clay-sized fraction and expansive clay (interlayered smectite minerals) contents at Site 2 are 30 and 10% higher, respectively, than at Site 1. Clay minerals affect the measurement of the VWC, when using time-domain reflectometry to measure the VWC (TDR; METERgroup, 2024). The presence of clay minerals, and expansive clay minerals more specifically, influences the relation between dielectric permittivity and VWC such that VWC is overestimated at high water contents and underestimated at low water contents (Bobrov et al., 2019; Kaviratna et al., 1996; Schwartz et al., 2013). For Site 2 we expect the effect of the calibration to be small, as the sensors have been calibrated for most of the measurement range and the amount of expansive clay minerals is about 25%. However, the effect cannot be ignored for Site 3, having an expansive clay content of 47%, while the TDT sensors were not calibrated at all. The measured VWCs are also relatively high, exceeding 0.60. This even outranges the maximum VWC measured in the drying experiments (0.50), which occurred after saturating the sample for 3 weeks. The lack of calibration across the full VWC range may have caused the VWC measurements, and consequently the estimated changes in VWC used to predict shrinkage and swelling effects, to be overestimated. This effect, in combination with the different restraining conditions of mostly swelling in field conditions compared to laboratory conditions, results in less predicted land surface movement than expected from the VWC measurements.

Simulating the layer deformation based on the matric suction, as proposed in this work, shows an interesting shift in the shrinkage behaviour with the drying intensity. For short term (days to week) changes in the matric suction, the best fit is found for the simulations based on the vertical deformation in combination with overburden pressure. However, when the matric suction drops significantly at a depth of 50 cm for Site 1 and 35 cm for Site 2, the simulation based on the volumetric shrinkage in combination with overburden pressure fits the measured deformation best. Both the newly proposed and the Briaud model incorporate a specific dynamic relation between the vertical strain and VWC. Briaud et al. (2003) found that overburden pressure does not influence the relation between VWC change and vertical strain. The underestimation of shrinkage is thus most likely due to low water content in the shallowest part of the layer of interest. The VWC was not measured at a depth of 6–10 cm, which creates the need to simulate the VWC at this depth. The simulation probably overestimates the VWC around a depth of 10 cm, where thus more shrinkage occurs in-situ than is predicted using the models.

For all sites, the predictions overestimate the short-term swell activity with respect to thickness changes, compared to the actual thickness changes. The predictions based on the soil matric suction yield an insight from the lack of agreement between measured and modelled short-term swell: the addition of an overburden term minimises the overestimation of short-term swell. The results show this term is more dominant with regard to swell, compared to shrinkage. Additionally, the seasonal swell was also overestimated by the predictions for Sites 1 and 3. At Site 1, the swell potential seems to be too low to regain the soil structure before drying during the monitored period, due to the overburden pressure (Peng & Horn, 2007). Swell potential is related to the specific surface area of the soil particles, related to the clay mineralogy (Prakash & Sridharan, 2004; Ross, 1978). The swell potential of Site 1 is thus less due to textural composition and expansive mineral content, compared to Sites 2 and 3, explaining the overestimation of the swell based on shrinkage. Also, the accessibility of the micropores can be reduced during shrinkage, whereas the VWC is not affected (Basma et al., 1996; Zhao et al., 2021). Most swell occurs in the meso- and micropore size range and can thus be limited when these close during shrinkage. This explains the yearly loss of surface layer thickness for Sites 1 and 2. Periods of slow positive deformation can be distinguished for Sites 1 and 2 (Figure 1). This slow swell is not affected by the VWC dynamics during this period. The increase in porewater pressure due to increase in the groundwater level yields a faster flow from the larger to smaller pores, which explains the swell in relation to the peaks in groundwater level (Figure 12 and Figure SI.4).

The Briaud et al. (2003) and matric suction-deformation methods work well to predict the seasonal shrinkage during a dry period. The application of the Briaud method can be improved by measuring the VWC at the top of the layer of interest too. This provides a time- and cost-efficient prediction method for shrinkage for all locations we have tested. The Briaud model can also be used based on simulated VWC values instead of measured, making the method widely applicable, especially as it is simple to apply with limited testing involved. The presented methods perform less when it comes to predicting swell. To simulate the land movement fully, the pore pressure changes due to high groundwater levels should be included, especially in the shallow (0–20 cm depth) zone where overburden pressure is zero to minimal. The response of the deformation is instantaneous when the groundwater level exceeds a depth of 20 cm. Including the measurement of deformation at higher water levels in the lab would provide a proper addition to be able to simulate the system above this threshold.

Land surface movement due to shrinkage and swelling is expected to increase with the increase of wet and dry spells as extreme weather conditions under climate change. The ability to predict shrinkage and swelling will thus become more important in temperate regions as the Netherlands. More intense and prolonged dry periods can drive soils beyond previously experienced moisture deficits, enhancing shrinkage through increased suction stresses and potentially causing irreversible structural rearrangement of the soil matrix (Basma et al., 1996; Estabragh et al., 2018). These structural changes may alter subsequent swelling and shrinkage behaviour by modifying pore geometry, particle bonding, and clay mineral orientation (Basma et al., 1996; Osipov et al., 1987). As increasingly severe dry periods penetrate deeper into the subsurface, deeper soil layers are also likely to experience moisture deficits comparable to those already occurring in near-surface layers. In cases where the main distinction between upper and lower horizons is their historically experienced minimum water content, the shrinkage behaviour at depth may progressively converge toward that of the upper soil horizons. As a result, prediction methods that relate shrinkage behaviour to VWC may remain relatively robust under changing climatic conditions, provided that the governing soil properties remain unchanged, provided that the governing soil properties remain unchanged. Swelling behaviour, however, is expected to remain more difficult to predict accurately, because it is more strongly influenced by soil structural changes, confinement conditions, hysteresis, and the wetting–drying history of the soil (Basma et al., 1996; Estabragh et al., 2018). Consequently, further investigation into field-scale swelling behaviour and its controlling mechanisms is required to improve future predictions of land surface movement and subsidence risk under climate change.

Conclusion

This study aimed to test the existing Briaud et al. (2003) model and a newly proposed simple model for predicting shallow subsurface deformation of clayey soils based on changes in VWC or matric suction in combination with laboratory-measured shrinkage behaviour of soil samples. The shrinkage behaviour determined in the lab varied within and between locations. However, the differences within locations were limited. The seasonal shrinkage was predicted well with the model devised by Briaud et al. (2003), where the offset between measured and simulated seasonal shrinkage was less than 1 mm. The offset in simulated and measured shrinkage was between 2 and 11% of the measured seasonal shrinkage for the tested locations. The performance of the simulations varied per location, partly due to the calibration of the VWC sensors (TDR and TDT), as well as the high groundwater levels that influence the porewater pressure. Due to the high groundwater levels, the simulations were outside the range where the shrinkage behaviour was measured in the experimental set-up. The Briaud model consistently outperformed other approaches in capturing seasonal shrinkage at all locations, even under conditions of potentially overestimated VWC values (Sites 2 and 3). The better performance of the Briaud model as applied in this work is likely related to the use of VWC as an input parameter, which represents a bulk, volume-averaged measurement of the sample, whereas the suction-based approach relies on matric suction measured at a single point within the sample. The method would be more robust when the suction would be monitored at multiple levels or in the middle of the sample. This approach would only be possible when using a more durable and robust tensiometer in the laboratory setup. Both the application of the Briaud model and the suction-based simulations highlight that the shrinkage behaviour changes with the moisture content. The application of the Briaud method presented here can be improved by including a simulation of the shallow VWC (<10 cm depth). Overall, the Briaud model is an adequate tool to predict land surface movement due to shrinkage and swelling, when accurate VWC values are available. The method could be further improved by including swell tests to capture hysteresis too. Yet, the results indicate that the shrinkage behaviour is not affected by irreversible shrinkage on this time-scale (several years), even though irreversible shrinkage was observed. In the context of projected increases in dry periods frequency and severity in the Netherlands and other temperate regions, the applicability of the presented methods becomes particularly relevant. The results indicate that the Briaud-based approach is robust in capturing seasonal shrinkage during prolonged dry periods, suggesting that it can be effectively applied to assess land surface deformation under future drier climate conditions. However, the need to include swell and hysteresis measurements becomes more urgent.

References

Abed, A.A., 2008. Numerical modeling of expansive soil behaviour. PhD thesis, Stuttgart University, Stuttgart.

Adem, H.H. & Vanapalli, S.K., 2015. Review of methods for predicting in situ volume change movement of expansive soil over time. Journal of Rock Mechanics and Geotechnical Engineering 7(1): 73–86. DOI: https://doi.org/10.1016/j.jrmge.2014.11.002

Alaoui, E.A., Lipiec, J. & Gerke, H.H., 2011. A review of the changes in the soil pore system due to soil deformation: a hydrodynamic perspective. Soil and Tillage Research 115–116: 1–15. DOI: https://doi.org/10.1016/j.still.2011.06.002

Azam, S., 2003. Influence of mineralogy on swelling and consolidation of soils in eastern Saudi Arabia. Canadian Geotechnical Journal 40(5): 964–975. DOI: https://doi.org/10.1139/t03-047

Basma, A.A., Al-Homoud, A.S., Malkawi, A.I.H. & Al-Bashabsheh, M.A., 1996. Swelling-shrinkage behavior of natural expansive clays. Applied Clay Science 11(2–4): 211–227. DOI: https://doi.org/10.1016/S0169-1317(96)00009-9

Bobrov, P.P., Belyaeva, T.A., Kroshka, E.S. & Rodionova, O.V., 2019. Soil moisture measurement by the dielectric method. Eurasian Soil Science 52: 822–833. DOI: https://doi.org/10.1134/S106422931905003X

Boivin, P., Garnier, P. & Tessier, D., 2004. Relationship between clay content, clay type, and shrinkage properties of soil samples. Soil Science Society of America Journal 68(4): 1145–1153. DOI: https://doi.org/10.2136/sssaj2004.1145

Briaud, J.L., Zhang, X. & Moon, S., 2003. Shrink test–water content method for shrink and swell predictions. Journal of Geotechnical and Geoenvironmental Engineering 129(7): 590–600. DOI: https://doi.org/10.1061/(ASCE)1090-0241(2003)129:7(590)

Bronswijk, J.J.B., 1991. Magnitude, modeling and significance of swelling and shrinkage processes in clay soils. Wageningen: Wageningen University.

Bronswijk, J.J.B. & Evers-Vermeer, J.J., 1990. Shrinkage of Dutch clay soil aggregates. Netherlands Journal of Agricultural Science 38(2): 175–194. DOI: https://doi.org/10.18174/njas.v38i2.16603

Davidson, S.E. & Page, J.B., 1956. Factors influencing swelling and shrinking in soils. Soil Science Society of America Journal 20(3): 320–324. DOI: https://doi.org/10.2136/sssaj1956.03615995002000030007x

Day, R.W., 1994. Swell-shrink behavior of compacted clay. Journal of Geotechnical Engineering 120(3): 618–623. DOI: https://doi.org/10.1061/(ASCE)0733-9410(1994)120:3(618)

Day, R.W., 1995. Ultimate density of a compacted clay subjected to cycles of wetting and drying. Environmental & Engineering Geoscience 1(2): 229–232. DOI: https://doi.org/10.2113/gseegeosci.I.2.229

Dexter, A.R. & Richard, G., 2009. Water potentials produced by oven drying of soil samples. Soil Science Society of America Journal 73: 1646–1651. DOI: https://doi.org/10.2136/sssaj2008.0294N

Dhowian, A.W. & Al-Saadan, T.A., 2010. Swell behavior of expansive soil with free lateral movements. Journal of King Saud University-Engineering Sciences 22(2): 51–63. DOI: https://doi.org/10.1016/S1018-3639(18)30493-8

Doris, J.J., Rizzo, D.M. & Dewoolkar, M.M., 2008. Forecasting vertical ground surface movement from shrinking/swelling soils with artificial neural networks. International Journal for Numerical and Analytical Methods in Geomechanics 32(10): 1229–1245. DOI: https://doi.org/10.1002/nag.666

Dörner, J., Dec, D., Peng, X. & Horn, R., 2009. Change of shrink age behavior of an Andisol in southern Chile: effects of land use and wetting/drying cycles. Soil and Tillage Research 106(1): 45–53. DOI: https://doi.org/10.1016/j.still.2009.09.013

Estabragh, A.R., Soltani, A. & Javadi, A.A., 2018. Effect of pore water chemistry on the behaviour of a kaolin–bentonite mixture during drying and wetting cycles. European Journal of Environmental and Civil Engineering 24(7): 895–914. DOI: https://doi.org/10.1080/19648189.2018.1428691

Fityus, S.G., Smith, D.W. & Allman, M.A., 2004. Expansive soil test site near Newcastle. Journal of Geotechnical and Geoenvironmental Engineering 130(7): 686–695. DOI: https://doi.org/10.1061/(ASCE)1090-0241(2004)130:7(686)

Greene-Kelly, R., 1974. Shrinkage of clay soils: a statistical correlation with other soil properties. Geoderma 11(4): 243–257. DOI: https://doi.org/10.1016/0016-7061(74)90052-4

Griffioen, J., Klaver, G. & Westerhoff, W.E., 2016. The mineralogy of suspended matter, fresh and Cenozoic sediments in the fluvio-deltaic Rhine–Meuse–Scheldt–Ems area, the Netherlands: an overview and review. Netherlands Journal of Geosciences 95(1): 23–107. DOI: https://doi.org/10.1017/njg.2015.32

GroPoint, 2025. GroPoint™ profile multi segment soil moisture & temperature profiling probe; models 2625-S-T-x, 2625-S-x operation manual. Victoria, BC: GroPoint.

Jones, G.M., Cassidy, N.J., Thomas, P.A., Plante, S. & Pringle, J.K., 2009. Imaging and monitoring tree-induced subsidence using electrical resistivity imaging. Near Surface Geophysics 7(3): 191–206. DOI: https://doi.org/10.3997/1873-0604.2009017

Jones, L.D. & Jefferson, I., 2012. Expansive soils.Expansive soils. In: Burland, J., Chapman, T., Skinner, H. & Brown, M. (eds), ICE manual of geotechnical engineering, Volume 1. London: ICE Publishing, pp. 413–441.

Kaviratna, P.D., Pinnavaia, T.J. & Schroeder, P.A., 1996. Dielectric properties of smectite clays. Journal of Physics and Chemistry of Solids 57(12): 1897–1906. DOI: https://doi.org/10.1016/S0022-3697(96)00076-5

Kell, G.S., 1975. Density, thermal expansivity, and compressibility of liquid water from 0. deg. to 150. deg.. Correlations and tables for atmospheric pressure and saturation reviewed and expressed on 1968 temperature scale. Journal of Chemical and Engineering Data 20(1): 97–105. DOI: https://doi.org/10.1021/je60064a005

Kirkpatrick, W.M. & Rennie, I.A., 1972. Directional properties of consolidated kaolin. Geotechnique 22(1): 166–169. DOI: https://doi.org/10.1680/geot.1972.22.1.166

Kodikara, J., Barbour, S.L. & Fredlund, D.G., 1999, February. Changes in clay structure and behaviour due to wetting and drying. In: Proceedings 8th Australia New Zealand Conference on Geomechanics: Consolidating Knowledge. Hobart: Australian Geomechanics Society, pp. 179–185.

Kooi, H. & Waardenburg, E.H., 2024. Evaluation of alternative extensometer designs and methods: desk study report for Regio Deal Bodemdaling Groene Hart project 44. Delft: Deltares.

Lexmond, B.R., van Dam, B., Hockin, C.V., Erkens, G., Griffioen, J. & Stouthamer, E., 2024. Measuring shrinkage of expansive soils using a novel automated high-frequency setup. Soil Science Society of America Journal 88(6): 2343–2352. DOI: https://doi.org/10.1002/saj2.20755

Liu, Z., Zhang, R., Liu, Z. & Zhang, Y., 2021. Experimental study on swelling behavior and its anisotropic evaluation of unsaturated expansive soil. Advances in Materials Science and Engineering 2021(1): 6937240. DOI: https://doi.org/10.1155/2021/6937240

Lu, N. & Dong, Y., 2017. Correlation between soil-shrinkage curve and water-retention characteristics. Journal of Geotechnical and Geoenvironmental Engineering 143(9): 04017054. DOI: https://doi.org/10.1061/(ASCE)GT.1943-5606.0001741

Makki, L., Duc, M., Coppée, T. & Szymkiewicz, F., 2024. The swelling–shrinkage properties of intact and disturbed clayey and marly soils: the density effect. Geotechnics 4(2): 512–529. DOI: https://doi.org/10.3390/geotechnics4020028

METER Group, 2024. TEROS 11/12 manual. Accessed 20 May 2025. 3.3.1 Volumetric Water Content. 20587_TEROS11-12_Manual_Web.pdf

Mitchell, J.K., 1993. Fundamentals of soil behavior. New York: Wiley.

Osipov, V.I., Bik, N.N. & Rumjantseva, N.A., 1987. Cyclic swelling of clays. Applied Clay Science 2(4): 363–374. DOI: https://doi.org/10.1016/0169-1317(87)90042-1

Overton, D.D., Chao, K.C. & Nelson, J.D., 2006. Time rate of heave prediction for expansive soils. In: Proceedings of the GeoCongress, Atlanta, GA. Reston, VA: American Society of Civil Engineers, pp. 1–6.

Peng, X. & Horn, R., 2005. Modeling soil shrinkage curve across a wide range of soil types. Soil Science Society of America Journal 69(3): 584–592. DOI: https://doi.org/10.2136/sssaj2004.0146

Peng, X. & Horn, R., 2007. Anisotropic shrinkage and swelling of some organic and inorganic soils. European Journal of Soil Science 58(1): 98–107. DOI: https://doi.org/10.1111/j.1365-2389.2006.00808.x

Peng, X. & Horn, R., 2013. Identifying six types of soil shrinkage curves from a large set of experimental data. Soil Science Society of America Journal 77(2): 372–381. DOI: https://doi.org/10.2136/sssaj2011.0422

Peron, H., Hueckel, T., Laloui, L. & Hu, L., 2009. Fundamentals of desiccation cracking of fine-grained soils: experimental characterisation and mechanisms identification. Canadian Geotechnical Journal 46(10): 1177–1201. DOI: https://doi.org/10.1139/T09-054

Prakash, K. & Sridharan, A., 2004. Free swell ratio and clay mineralogy of fine-grained soils. Geotechnical Testing Journal 27(2): 220–225. DOI: https://doi.org/10.1520/GTJ10860

Puppala, A.J. & Cerato, A., 2009. Heave distress problems in chemically-treated sulfate-laden materials. Geo-Strata – Geo Institute of ASCE 10(2): 28.

Ross, G.J., 1978. Relationships of specific surface area and clay content to shrink-swell potential of soils having different clay mineralogical compositions. Canadian Journal of Soil Science 58(2): 159–166. DOI: https://doi.org/10.4141/cjss78-020

Schindler, U., Doerner, J. & Mueller, L., 2015. Simplified method for quantifying the hydraulic properties of shrinking soils. Journal of Plant Nutrition and Soil Science 178(1): 136–145. DOI: https://doi.org/10.1002/jpln.201300556

Schwartz, R.C., Casanova, J.J., Pelletier, M.G., Evett, S.R. & Baumhardt, R.L., 2013. Soil permittivity response to bulk electrical conductivity for selected soil water sensors. Vadose Zone Journal 12(2): vzj2012-0133. DOI: https://doi.org/10.2136/vzj2012.0133

Seneviratne, S.I., Zhang, X., Adnan, M., Badi, W., Dereczynski, C., Di Luca, A., Ghosh, S., Iskandar, I., Kossin, J., Lewis, S., Otto, F., Pinto, I., Satoh, M., Vicente-Serrano, S.M., Wehner, S. & Zhou, B., 2021. Weather and climate extreme events in a changing climate. In: V. Masson-Delmotte, P. Zhai, A. Pirani, S.L. Connors, C. Péan, S. Berger, N. Caud, Y. Chen, L. Goldfarb, M.I. Gomis, M. Huang, K. Leitzell, E. Lonnoy, J.B.R. Matthews, T.K. Maycock, T. Waterfield, O. Yelekçi, R. Yu & B. Zhou (eds.). Climate change 2021: the physical science basis. Contribution of Working Group I to the Sixth Assessment Report of the Intergovernmental Panel on Climate Change. Cambridge and New York, NY: Cambridge University Press, pp. 1513–1766. DOI: https://doi.org/10.1017/9781009157896.013

Shi, B., Jiang, H., Liu, Z. & Fang, H.Y., 2002. Engineering geological characteristics of expansive soils in China. Engineering Geology 67(1–2): 63–71. DOI: https://doi.org/10.1016/S0013-7952(02)00145-X

Šimůnek, J., Šejna, M., Saito, H., Sakai, M. & van Genuchten, M.T., 2013. The hydrus-1D software package for simulating the movement of water, heat, and multiple solutes in variably saturated media, version 4.17, HYDRUS Software Series 3. Riverside, CA: Department of Environmental Sciences, University of California Riverside, p. 342.

van Asselen, S., Erkens, G. & de Graaf, F., 2020. Monitoring shallow subsidence in cultivated peatlands. Proceedings of the International Association of Hydrological Sciences 382: 189–194. DOI: https://doi.org/10.5194/piahs-382-189-2020

Van Asselen, S., Kooi, H., Bakema, G., Van der Marel, H., Hanssen, R., Koster, K. & Hartig, R., 2024a. Een bodembewegingsmeetnetwerk in het Groene Hart. Deelrapport van Regio Deal Bodemdaling Groene Hart project 44: Bodemdaling in kaart en kijken in de bodem; meten, monitoring en ontrafelen van bodemdalingsprocessen. Delft: Deltares.

Van Asselen, S., Kooi, H., Erkens, G., Bakema, G., Van Diepen, S., Conroy, P. & Hanssen, R., 2024b. Eerste resultaten van bodemdalingsmeetlocaties en een dynamisch maaiveldhoogtemodel (D-DEM) in het Groene Hart Deelrapport van Regiodeal Bodemdaling Groene Hart project 44. Delft: Deltares.

Van der Meulen, M.J., Maljers, D., Van Gessel, S.F. & Gruijters, S.H.L.L., 2007. Clay resources in the Netherlands. Netherlands Journal of Geosciences 86(2): 117–130. DOI: https://doi.org/10.1017/S001677460002312X

Van Essen, H., 2025. Laboratory research on swelling clays interim report. Deltares with support from Holland Innovation team. Document ID: 11211478-004-BGS-0001. Delft: Deltares.

Vignozzi, N. & Pellegrini, S., 2004. Soil structure and the effect of management practices. Soil and Tillage Research 79: 131–143. DOI: https://doi.org/10.1016/j.still.2004.07.002

Vu, H.Q. & Fredlund, D.G., 2004. The prediction of one-, two-, and three-dimensional heave in expansive soils. Canadian Geotechnical Journal 41(4): 713–737. DOI: https://doi.org/10.1139/t04-023

Wray, W.K., El-Garhy, B.M. & Youssef, A.A., 2005. Three-dimensional model for moisture and volume changes prediction in expansive soils. Journal of Geotechnical and Geoenvironmental Engineering 131(3): 311–324. DOI: https://doi.org/10.1061/(ASCE)1090-0241(2005)131:3(311)

Zhang, X., 2004. Consolidation theories for saturated-unsaturated soils and numerical simulations of residential buildings on expansive soils. PhD thesis, Texas A&M University, College Station, TX.

Zhao, G.T., Zou, W.L., Han, Z., Wang, D.X. & Wang, X.Q., 2021. Evolution of soil-water and shrinkage characteristics of an expansive clay during freeze-thaw and drying-wetting cycles. Cold Regions Science and Technology 186: 103275. DOI: https://doi.org/10.1016/j.coldregions.2021.103275