# A free energy based model for water transfer in amphiphilic soils

Florian Cajot, Claude Doussan, Philippe Beltrame

► To cite this version:

Florian Cajot, Claude Doussan, Philippe Beltrame. A free energy based model for water transfer in amphiphilic soils. Advances in Water Resources, 2025, 198, pp.104915. 10.1016/j.advwatres.2025.104915. hal-04965328

HAL Id: hal-04965328

https://hal.science/hal-04965328v1

Submitted on 25 Feb 2025

HAL is a multi-disciplinary open access
archive for the deposit and dissemination of sci-
entific research documents, whether they are pub-
lished or not. The documents may come from
teaching and research institutions in France or
abroad, or from public or private research centers.

L'archive ouverte pluridisciplinaire HAL, est
destinée au dépôt et à la diffusion de documents
scientifiques de niveau recherche, publiés ou non,
émanant des établissements d'enseignement et de
recherche français ou étrangers, des laboratoires
publics ou privés.

<figure><img src="image_1.png" /></figure>

Distributed under a Creative Commons Attribution 4.0 International License
<figure><img src="image_2.png" /></figure>

<figure><img src="image_3.png" /></figure>

# A free energy based model for water transfer in amphiphilic soils

Florian Cajot<sup>ID</sup> *, Claude Doussan<sup>ID</sup>, Philippe Beltrame<sup>ID</sup>

UMR1114 EMMAH INRAE-AU, 228, Route de L'Aérodrome, Avignon, 84194, France

## ARTICLE INFO

**Keywords:**

Water transfer model

Amphiphilic matter

Wettability

Free energy

Exopolysaccharides

Rhizosphere

Soil

## ABSTRACT

A 3D macroscopic gradient-dynamics model is developed and applied to sandy soil in presence of exopolysaccharides (EPS), to mimic a soil influenced by root exuded mucilages (rhizospheric soil). Depending on water content, amphiphilic soil has a hydrophilic or hydrophobic behavior which impacts water transfer and retention. To model this saturation-dependent wettability, we propose a nonequilibrium thermodynamic approach based on the definition of the free energy of the system. The free energy functional contains gravity energy, the free surface energy of the water and the effective interaction (attractive and repellent) between water and the porous matrix with the amphiphilic matter. The latter defines the wettability of the porous medium. Water flow dynamics is derived from Onsager's variational principle leading to a non-linear fourth order PDE on the saturation generalizing the Richards equation. The new formulation reproduces a range of water flow regimes encountered in soil with EPS: (i) the stoppage of imbibition front in a homogeneous soil leading to equilibrium where moist regions coexist with a dry region, (ii) a decrease in capillary height in comparison with a sand without amphiphilic matter and (iii) the existence of a threshold of amphiphilic concentration for which the capillary rise is stopped at the dry layer containing the amphiphilic matter. After calibrating the parameters of our model, numerical simulation is in qualitative and quantitative agreement with experiments from the literature.

## 1. Introduction

Rising temperatures and climate deregulation over long periods impact soils and vegetation (Intergovernmental Panel on Climate Change (IPCC), 2014; Sinha et al., 2021). These changes influence the dynamics of the soil water status. In particular, higher evapotranspiration, longer and more frequent drought periods put an emphasis on increased importance of soil water and its use by vegetation and crops (and possible associated water stress). As a result, this will increase the pressure on available water resources because of a greater need for irrigation. This calls for a reliable management of water in soil, with a sound understanding of water transfer and retention in soil as well as a better knowledge of the interactions between climate, soil and vegetation (Rodriguez-Iturbe, 2000).

Soil and plants are intimately related in soil by the soil zone close to the roots, called the rhizosphere (Hinsinger et al., 2009; Shi et al., 2016). The activity of roots strongly influences the biogeochemical activity and processes of the soil in this zone. In addition to processes such as water and nutrient uptake or growth, roots deliver a range of organic matter to soil, among which mucilage. The mucilage is a gel substance mainly made up of Exopolysaccharides (EPS) (Bérard et al., 2020; Nazari, 2021). In recent literature, it has been demonstrated that mucilage can have a hydrophobic effect (Ahmed et al., 2015; Benard

et al., 2019; Knott et al., 2022) and that EPS and other molecules in mucilage have part of amphiphilic behavior due to re-orientation of polar/apolar part (Knott et al., 2022; Flemming et al., 2023). The change in wettability induced by EPS is assumed to be the cause of the regulating effect of water transfer in the rhizosphere.

More generally, water transfer and retention in an unsaturated porous medium depend not only on the structure of the medium (porous matrix geometry), but also on the physico-chemical interactions (hydrophilic or hydrophobic) between the solid matrix, the gas phase and the water. The polar affinity between soil minerals and water leads to hydrophilic interactions. These interactions are responsible for the movement of water toward drier zones, as shown for example in the capillary rise process (Hird and Bolton, 2017; Leverett, 1941). In presence of amphiphilic matter, this physico-chemical interaction is amplified in moist soil and then increases the water storage in the rhizosphere (Le Gall et al., 2021; Doussan et al., 2024). In contrast, in dry soil, amphiphilic matter provides a hydrophobic property leading to fingering flows or slowing down the water flow as a result of the delay of the rhizosphere rehumectation called the buffering effect (Carminati et al., 2010). This latter phenomenon is also revealed in the capillary rise experiment in a stratified sand column with an upper sand layer

<footnote>* Corresponding author.</footnote>

<footnote>E-mail addresses: florian.cajot@gmail.com (F. Cajot), claude.doussan@inrae.fr (C. Doussan), philippe.beltrame@univ-avignon.fr (P. Beltrame).</footnote>
impregnated with EPS (Kroener et al., 2014; Kroener, 2015). The higher the concentration, the slower the rise of the water saturation front in the sand-EPS layer. Beyond a threshold concentration, the water no longer infiltrates for at least 5 minutes. In the literature, the observed effects are primarily attributed to a significant reduction in hydraulic conductivity at high EPS concentrations. In contrast, the effects of the change in wettability are modeled by introducing a relaxation time associated with conformational changes in amphiphilic matter (Carminati et al., 2010; Ahmed et al., 2014; Kroener et al., 2014), which play only a minor role in slowing the capillary rise. Therefore, the impacts of variable wettability linked to amphiphilic properties remain largely unexplored and poorly modeled (Benard et al., 2019).

Wettability models in the literature are based on the definition of a capillary pressure representing water suction. There are several relationships relating capillary pressure and saturation (i.e. the water retention curve) (Brooks and Corey, 1966; van Genuchten, 1980) relating pressure and saturation. All of them are decreasing functions with respect to the water saturation. The capillary pressure term in addition to the gravity term are integrated into Darcy law for unsaturated porous media leading to the Richards equations. This water transfer model in unsaturated porous media shows behaviors characteristic of hydrophilic interaction, such as capillary rise or wetting bulb dynamic to the driest regions.

The wettability of the medium (more or less hydrophilic), with this approach, is taken into account by a contact angle through the capillary height in the model of capillary pressure. The relationship between the measured capillary height and the value used in the Richards model is not obvious, since the simulated equilibrium saturation profile does not show a clear demarcation between the wet and dry regions, unlike the capillary rise experiment. Another approach to estimate the soil wettability is based on infiltration time experiments such as Water Drop Penetration Time (WDPT) or Molarity of Ethanol Droplet (MED). The medium is qualified from very hydrophilic to very hydrophobic depending on the penetration time (Bisdom et al., 1993) then it is associated with a contact angle (Mandal and Jayaprakash, 2009).

The link between contact angle and capillary height is based on the strong assumption that the porous medium is represented by cylindrical tubes. However, it has been shown that the complex geometry of the matrix makes this relationship questionable (Bamorovat Abadi and Bahrami, 2022). Moreover, the experimental measurement of a contact angle for a soil is complex to assess. Even if the measure of this angle and the effect of hydrophobicity at pore scale have been reported (Müehl et al., 2012), the complexity of a porous medium and organic matter distribution make upscaling of microscopic wettability to the macroscale a very challenging task (Rabbani et al., 2018). On the other hand, some experiments are based on measuring the contact angle on a thin layer of porous medium on an impermeable substrate (Alghunaim et al., 2016; Heib and Schmitt, 2016). There is no proof, however, that this measure represents the contact angle in the porous medium.

Finally, the most fundamental point is that it has been shown that the Richards equation cannot reproduce phenomena related to hydrophobicity (DiCarlo, 2013) such as spontaneous fingered flows in a hydrophobic layer even if the layer is homogeneous (Nieber et al., 2005).

One of the hypotheses of the link between hydrophobicity and fingering is the existence of the hysteresis of the contact angle which is revealed through the hysteresis of the retention curves (Beliaev and Hassanizadeh, 2001; Eliassi and Glass, 2002). In the last decade, several models introduced time-dependent terms directly into a Richards-type equation in order to reproduce an unstable wetting front or the buffering effect (Rätz and Schweizer, 2012; Wang and Wallach, 2022). However, these ad-hoc terms are still related to the notion of contact angle in porous media.

In the context of water transfer in the rhizosphere, Kroener et al. (2014) introduced a relaxation term into the Richards equation to represent the delay in water content balancing with the matrix containing

mucilage. In this way, they reproduce the effect of hysteresis during rewetting/drying cycles.

Therefore, there is a consensus in the literature that the dynamics of water transfer in a porous medium containing amphiphilic matter is a non-equilibrium process that is not reproduced by the Richards model.

Recently, a soil water transfer model has been developed to simulate various nonequilibrium phenomena linked to the presence of hydrophobic matter, such as water trapping, intermittent flow and fingering instability (Beltrame and Cajot, 2022). Instead of a contact angle, the soil wettability is associated with the energy of attractive and repellent interactions in a similar way to substrate wettability associated with the disjoining pressure (de Gennes, 1985; Derjaguin et al., 1987). Because of the hydrophobic interaction, the effective disjoining pressure in the soil is not necessarily monotonic, unlike the capillary pressure. The governing equation is a gradient dynamics model based on the functional free energy in the soil using Onsager's variational principle (Doi, 2011). It is similar to the Richards model (Richards, 1931) but (i) the capillary pressure is replaced by a disjoining pressure deriving from the wettability energy and (ii) an effective surface tension energy is added. The latter term was first introduced in soils to explain the spontaneous appearance of preferential flows interpreted as a Saffman-Taylor instability by Cueto-Felgueroso and Juanes (2008) and Beljadid et al. (2020). In Beltrame and Cajot (2022), the front instability results mainly from hydrophobic interaction in a stratified porous medium whereas the effective surface tension scales the fingers width. The latter term is required to have a well-posed mathematical formulation in the case where the disjoining pressure is no longer monotonic (Beltrame and Cajot, 2022).

This approach appears particularly promising for the case of amphiphilic matter, as the wettability energy integrates both hydrophilic and hydrophobic interactions.

In this paper, we extend the thermodynamic approach and dynamic model of Beltrame and Cajot (2022) to the presence of amphiphilic matter in the soil or porous medium. Thus, the model derives from a functional free energy including wettability and free surface tension. The advantage of this theoretical framework is that the out-of-equilibrium processes that occur remain consistent with thermodynamics.

The paper is organized as follows. Section 2 is devoted to the construction of the functional free energy in a soil containing amphiphilic matter, resulting in a thermodynamic model of water transfer. A one-dimensional numerical simulation of horizontal water flow and capillary rise is performed in Section 3 to highlight the qualitative effects of the amphiphilic matter in the horizontal and vertical water dynamic of a homogeneous sandy medium. In order to calibrate the parameters of our model, we simulate capillary rise experiments (Kroener et al., 2014; Kroener, 2015) of water in stratified media incorporating a layer with amphiphilic matter. Finally, we provide some concluding remarks in Section 4.

## 2. Modeling

### 2.1. *Outline of the problem*

We consider water transfer in a sandy soil containing a concentration $c$ of EPS, an amphiphilic matter that is widely present in the rhizosphere. This specific choice of sand with EPS is motivated by the opportunity to use controlled experiments from the literature. The concentration $c$ is defined as the mass of dry EPS in 100 g of dry substrate. To mimic experiment of the literature, we address two water flow configurations: (a) water flow in a horizontal homogeneous domain of length $L$ (Fig. 1(a)) and (b) capillary rise in a vertical stratified sandy media of height $L$ with or without EPS (Fig. 1(b)).

Regardless the water flow configuration, we model the amphiphilic porous medium starting from the definition of the components of the free energy of the system. We then derive the equation of the dynamics of the water saturation in amphiphilic porous medium with the Onsager principle.
<figure><img src="image_4.png" /></figure>

(a) Horizontal water flow in a homogeneous sandy soil containing EPS.

<figure><img src="image_5.png" /></figure>

(b) Capillary rise experiment in a stratified sandy soil containing a layer of sand with EPS on a layer of sand without EPS.

<caption>Fig. 1. Sketch of the geometry studied. The hatched area corresponds to the medium containing EPS.</caption>

## 2.2. Free energy in amphiphilic porous medium

Let describe the fundamental property and form of the water free energy $F$ in amphiphilic porous medium. The density free energy is defined for a homogeneous and isotropic Representative Elementary Volume (REV). Spontaneous water movements are due to gradient of the density of free energy $f$ in the medium. In presence of amphiphilic matter, it depends on the saturation $S$ and the concentration of amphiphilic matter $c$,

$$F[S, c] = \int_{0}^{L} f(S, c)dx. \qquad (1)$$

The density of free energy is the sum of the energies which are in the medium: $f_w$ the wettability energy between water and porous matrix (Beltrame and Cajot, 2022), and $f_v$ a free surface energy of water in the porous medium (Cueto-Felgueroso and Juanes, 2008),

$$f(S, c) = f_w(S, c) + f_v(S). \qquad (2)$$

Each of these terms will be described in the following subsections.

The functional derivative of the free energy $\delta F/\delta S$ expresses the hydraulic pressure, also called water potential in soil literature (Richards, 1931; Emmerich, 2008; Cueto-Felgueroso and Juanes, 2009a).

### 2.2.1. Soil wettability and wetting free energy

We define soil wettability as any interaction between porous matrix and water regardless of hydrophilic or hydrophobic nature of the interaction. In this study, we neglect the effect of the solutes interactions with the porous matrix which can affect soil wettability and we assume that amphiphilic matter is adsorbed on the solid matrix. In literature, porous matrix and water interaction is usually associated with capillary pressure (Richards, 1931; Leverett, 1941; Cueto-Felgueroso and Juanes, 2009a). This interaction is then attractive. However, as we have shown in Beltrame and Cajot (2022), pressure formulations that consider only capillary pressure cannot reproduce the stoppage of imbibition front that typically occurs in a hydrophobic porous medium. To overcome this problem, a repellent interaction acting at low saturation has been added in the free energy (Beltrame and Cajot, 2022). This coexistence of attractive and repulsive interactions is analogous of partial wetting of a liquid on a substrate (de Gennes, 1985). The formulation of the pressure proposed in Beltrame and Cajot (2022) is relevant in the present study because it can represent hydrophilic character of the medium at higher water content and hydrophobicity with lower water saturation, as amphiphilic matter does.

In a porous medium close to the saturation, it is common to observe that air is trapped in the soil structure. In particular, the air can be surrounded by water in micro-pores. The contact line may be pinned due to contact angle hysteresis leading to water trapping. In this framework, variation in saturation results in a variation in air pressure. The negative of this compression work can be classically interpreted as free energy noted $f_{comp}$ (Cueto-Felgueroso and Juanes, 2009b). Since this trapping phenomenon essentially occurs only for saturation close to one, this term depends on the saturation and is non-zero in a range of values close to saturation. Although this is not really an interaction between a porous matrix and water, it is included in wettability energy.

From these premises, the wettability density energy in the porous medium $f_w$ is the sum of three terms: capillary energy $f_c$ (attractive interaction), the repellent energy $f_r$ (repellent interaction) and the air compression energy $f_{comp}$,

$$f_w(S, c) = f_c(S) + f_r(c, S) + f_{comp}(S). \qquad (3)$$

Let us specify each of these terms.

**Capillary energy $f_c(S)$**. We adopt the Brooks and Corey analytical expression for the capillary density energy relevant for sandy porous medium (Brooks and Corey, 1966; Leverett, 1941; Cueto-Felgueroso and Juanes, 2009b),

$$f_c(S) = -\rho g H_{cap} \frac{\lambda}{\lambda - 1} S^{1-\frac{1}{\lambda}}, \qquad (4)$$

with $H_{cap}$ the capillary height and $\lambda \neq 1$ the Brooks and Corey parameter (Brooks and Corey, 1966) related to the pore size distribution.

**Repellent energy $f_r(S, c)$**. We assume that repellent interaction is only due to the presence of EPS such that the hydrophobic effect increases with respect to the EPS concentration. We call $H_{rep}(c)$ a height depending on $c$ increasing the repellent energy with respect to $c$ such as:

$$f_r(S, c) = \rho g H_{rep}(c) \tilde{f}_r(S), \qquad (5)$$

where $\tilde{f}_r(S)$ is the normalized repellent energy depending only on saturation $S$. The latter is modeled in accordance with two properties: (a) apparent hydrophobicity decreases in the presence of water, and (b) it vanishes when saturation reaches about $S = 0.5$ (Kroener et al., 2018). To this end, we adopt a Gaussian distribution centered at $S_r = 0.05$ with a dispersion coefficient $\sigma = 0.2$ for the repellent pressure term $-\partial f_r/\partial S$. Note that the Gaussian is not centered at zero to avoid any singularity in the case of a dry porous medium modeling. As a
<figure><img src="image_6.png" /><figcaption>Fig. 2. Disjoining pressure in fine sandy soil for different parameter $\chi$. Parameters are $\rho = 10^3$ kg m$^{-3}$, $g = 9.81$ m s$^{-2}$, $H_{cap} = 0.26$ m, $S_r = 0.05$, $\sigma = 0.2$, $\lambda = 0.972$ and $\kappa = 200$.</figcaption></figure>

physical consequence, this assumption implies that a few percent of water always remains (Beltrame and Cajot, 2022). By integrating the repellent pressure term, we obtain the repellent density energy as a erf function of S:

$$\tilde{f}_r(S) = 2 \operatorname{erf} \left( \frac{S - S_r}{\sqrt{2\sigma}} \right). \qquad (6)$$

The parameters $S_r$ and $\sigma$ are constants fixed throughout the article.

We associate to the height $H_{rep}$ the dimensionless factor $\chi(c)$ such as $H_{rep}(c) = H_{cap} \chi(c)$:

$$\chi(c) = \frac{H_{rep}(c)}{H_{cap}}. \qquad (7)$$

The c-dependence of $H_{rep}$ and $\chi$ will be determined by an empirical relationship in Section 3.3.2.

**Air compression energy.** The air compression density energy $f_{comp}$ writes

$$f_{comp}(S) = \operatorname{sgn}(\lambda - 1) \exp[-\kappa(1-S)] f_c(S), \qquad (8)$$

where $\kappa \ge 0$ is a parameter of air compression (Cueto-Felgueroso and Juanes, 2009a). Therefore, the porous medium is never completely saturated and we call $S_{comp}$ the saturated profile of water which is slightly inferior to one. This profile is a stable equilibrium saturation defined numerically by Maxwell construction and it allows to avoids having numerical values of S greater than one. This numerical phenomenon is enhanced by overshoot appearing during the advancement of a front (DiCarlo, 2004; Nieber et al., 2005; Egorov et al., 2003; DiCarlo, 2007; Cueto-Felgueroso and Juanes, 2009b), and it is also particularly present in hydrophobic media (Beltrame and Cajot, 2022).

**Disjoining pressure.** From the wetting energy density $f_w$, we derive the dimensionless pressure $\Pi_w$, called disjoining pressure in Beltrame and Cajot (2022), which is given by

$$\Pi_w(S, c) = -\frac{1}{\rho g H_{cap}} \partial_S f_w(S, c), \qquad (9a)$$

$$= -\frac{1}{\rho g H_{cap}} \partial_S [f_c(S) + f_{comp}(S)] - \chi(c) \partial_S \tilde{f}_r(S). \qquad (9b)$$

Until Section 3.3.2, we study the $\Pi_w(S, c)$ curves as a function of S and the amplitude $\chi$ instead of c as we do not know the relationship between $\chi$ and c. For simplicity, we omit the c-dependency of the $\Pi_w$ function.

The disjoining pressure expressed in Eq. (9b) results from Eqs. (3), (5) and (7). Thus, without hydrophobic interactions, i.e. $\chi = 0$ then $f_r = 0$, $\Pi_w$ corresponds to the capillary pressure plus an air compression term. The curve $\Pi_w(S)$ is decreasing, as shown in Fig. 2, for $\chi = 0$. The only minimum is reached for a value near one ($S_{comp} = 0.956$). In the absence of other interactions, water diffuses uniformly through a homogeneous porous medium. This is not necessarily true, if the repellent term is present ($\chi \neq 0$).

For $\chi > \chi_i \approx 1.73$, an inflection point appears and the curve is no longer monotonous, i.e $\Pi_w(S)$ increases on the interval $[S_m, S_M]$ where the saturations $S_m$ and $S_M$ are extrema of $\Pi_w(S)$. For $S > S_M$ hydrophobicity becomes negligible and $\Pi_w$ only slightly differs from the capillary pressure of Brooks and Corey model. For S in the neighborhood of the $S_m$ value, the action of hydrophobicity is the strongest. In the saturation range $[S_m, S_M]$, uniform water distribution is no longer a stable equilibrium. Hydrophobicity causes dewetting, giving rise to a new equilibrium in which certain dry regions spatially coexist with wet ones (Beltrame and Cajot, 2022). These equilibria are studied in Section 3.1 and a thermodynamics interpretation is given.

### 2.2.2. Water free surface energy in porous medium

In unsaturated porous media, the presence of an interface between water and air implies the existence of a free surface. This free surface can be associated to an energy due to surface tension. To estimate this energy in a porous medium, Cueto-Felgueroso and Juanes (2008) propose a phase field approach where the saturation is the phase parameter. The authors applied the Cahn-Hilliard approach (Cahn and Hilliard, 1958) to porous media where saturation takes the role of phase field in a two-phase medium consisting of water and air. The energy associated to the water/air interface is expressed as in Cueto-Felgueroso and Juanes works (Cueto-Felgueroso and Juanes, 2008, 2009a,b,c):

$$f_v = \frac{\Gamma}{2} |\nabla S|^2, \qquad (10)$$

where $\Gamma$ is called the macroscopic surface tension (Cueto-Felgueroso and Juanes, 2009b). The free energy term $f_v$ is also called non-local term since it depends on the gradient of the saturation and not only the local value of saturation. Physically, it represents the existence of a continuous free surface. Recently, Beljadid et al. (2020) proposed a more complex expression where macroscopic surface tension is considered as a function depending on the saturation in the medium, called capillary pinning function. In this present study, we consider the expression given by Eq. (10).

In presence of repellent interaction, it is crucial to consider this surface energy. Without this effective surface tension $\Gamma$, as shown in Beltrame and Cajot (2022), hydrophobic interactions lead to infinitely small watered zones and a mathematical singularity, which are not physical phenomena. Therefore, this free surface energy in the porous medium plays a regularizing role in the wetting process due to the repellent part of the wettability energy (Beltrame and Cajot, 2022).

### 2.3. Governing equations of water saturation evolution

To establish the equation for the evolution of water saturation, we use a generic construction. The equation must have the form of a continuity equation:

$$\phi \partial_t S = - \nabla \cdot \mathcal{J} \qquad (11)$$

where the hydraulic flow $\mathcal{J}(S)$ is related to the gradient of free energy. To specify its expression, we assume that the energy dissipation is a quadratic in the flux $\mathcal{J}(S)$ which constitutes the simplest dissipation functional. Based on the Onsager variational principle (Doi, 2011) then the hydraulic flux has to be related to the free energy as follows (Cajot, 2024):

$$\mathcal{J} = K(S) \nabla \left( \frac{\delta \mathcal{F}}{\delta S} - \rho g z \right), \qquad (12)$$

where $K(S)$ is a positive function and $\delta \mathcal{F}/\delta S$ is the variational derivative. Finally, the time evolution equation reads

$$\phi \partial_t S = \nabla \cdot \left[ K(S) \nabla \left( \frac{\delta \mathcal{F}}{\delta S} - \rho g z \right) \right], \qquad (13a)$$

$$\frac{\delta \mathcal{F}}{\delta S} = \partial_s f_w(S, c) - \Gamma \Delta S. \qquad (13b)$$

We retrieve the classical description of the water flux in soil literature with the Richards equation by identifying $K(S)$ as the hydraulic conductivity and replacing $-\delta \mathcal{F}/\delta S$ without the macroscopic surface
<caption>**Table 1**</caption>

Brooks–Corey parameters for sandy soil extracted from Kroener et al. (2014) and Kroener et al. (2018).

<table><thead><tr><th>Soil</th><th>H<sub>cap</sub> [mm]</th><th>λ</th></tr></thead><tbody><tr><td>Fine sand</td><td>260</td><td>0.972</td></tr><tr><td>Medium sand</td><td>130</td><td>0.77</td></tr><tr><td>Coarse sand</td><td>10.2</td><td>0.561</td></tr></tbody></table>

tension by the capillary pressure. Therefore, Eq. (13a) is a generalization of the Richards equation classically used to simulate water transfer in unsaturated porous media (Richards, 1931).

Let us introduce the relative conductivity $k_r$ such as

$$K(S) = \frac{K_{sat}}{\rho g} k_r(S), \qquad (14)$$

where $K_{sat}$ is the saturated hydraulic conductivity. The function $k_r(S)$ is an increasing and generally convex function (Richards, 1931; Bear, 1988). The hydraulic conductivity $K$ is assumed to be a function invariant with respect to $c$ in order to highlight the effect of wettability on the system. In the following, we adopt the van Genuchten–Mualem formulation to represent the relative permeability (Mualem, 1976; van Genuchten, 1980),

$$k_r(S) = \sqrt{S} \left[ 1 - \left( 1 - S^{\frac{1}{m}} \right)^m \right]^2, \qquad (15)$$

where the van Genuchten–Mualem $m \in [0; 1[$ and $n > 1$ parameters are linked to the pore size distribution and the relationship between $m$ and $n$ is $m = 1 - 1/n$.

To solve the PDE system (Eq. (13)), we use a Finite Element Method implemented by the C++ FEM library OOMPH-LIB (Heil and Hazel, 2006). For time-stepping, a Backward Differentiation Formula (BDF) scheme of second order is used, while both the time-steps and the spatial mesh are adaptive to the dynamics (Cajot, 2024).

In the next section, the PDE system (Eq. (13)) is projected onto horizontal $x$ axis or the vertical $z$ axis depending on the problem studied.

### 3. Numerical simulation

In this section, the dynamics of water transfer for the two configurations of Fig. 1 are numerically simulated. Different EPS concentrations and sand grain sizes are studied. We consider three grain sizes: fine, medium and coarse. The corresponding retention parameters are summarized in Table 1.

#### 3.1. Horizontal homogeneous fine sandy soil

We study the water flow in a horizontal homogeneous fine sandy layer of length $L = 400$ mm (Fig. 1(a)). We represent the water saturation as a 1D model, i.e. spatial saturation profile depends only on $x$. At the boundaries $x = 0$ and $x = L$, we consider no water flux and no pressure gradient. In this context, the total water saturation is conserved. These physical boundary conditions are mathematically expressed as Neumann boundary conditions for $S$ and $\delta F/\delta S$. We project on the $x$-axis the PDE system (Eq. (13)) with Neumann boundary conditions:

$$\phi \partial_t S = \partial_x \left[ K(S) \partial_x \left( \frac{\delta F}{\delta S} \right) \right], \qquad (16a)$$

$$\frac{\delta F}{\delta S} = \partial_S f_w(S, c) - \Gamma \partial_{xx}^2 S, \qquad (16b)$$

$$\frac{\partial S}{\partial n} = 0 \qquad \text{at } x = 0, \text{ and } x = L, \qquad (16c)$$

$$\frac{\partial}{\partial n} \left( \frac{\delta F}{\delta S} \right) = 0, \qquad \text{at } x = 0 \text{ and } x = L, \qquad (16d)$$

where $n$ denotes the exterior normal of the domain boundary ($x = 0$ and $x = L$). Note, that in our model, four boundary conditions are required, in contrast to the Darcy–Richards model, for which only two boundary

<figure><img src="image_7.png" /><figcaption>**Fig. 3.** Horizontal saturation $S(x)$ profile in homogeneous fine sandy soil for different values of repellent factor $\chi$. Black dashed line is the initial profile of saturation $S_0$ defined in Eq. (17). Parameters are given for fine sandy soil in Table 1 and the other are $\rho = 10^3$ kg m$^{-3}$, $g = 9.81$ m s$^{-2}$, $S_r = 0.05$, $\sigma = 0.2$, $\phi = 0.43$, $K_{sat} = 10^{-4}$ m s$^{-1}$, $\Gamma = 1.7$ mN and $\kappa = 200$.</figcaption></figure>

conditions are required.

All numerical simulations are performed with the same initial water distribution $S_0(x)$. Water is initially concentrated on $1/4$ of the length $L$ at the middle of the domain with a saturation $S = 0.92$ and the rest of the domain is practically dry ($S = 0.02$). We model this initial distribution $s_0(x)$ with a continuous profile (black dashed line in Fig. 3) given by

$$S_0(x) = 0.9 \frac{1 + \tanh \left[ 10 \cos \left( \pi \tanh \left\{ 0.011 \left( x - \frac{L}{2} \right) \right\} \right) \right]}{2} + 0.02 \qquad (17)$$

Fig. 3 displays the equilibrium profile of horizontal saturation $S(x)$ for different values of the repellent factor $\chi$ whose disjoining pressure is given by Eq. (3) for fine sand in Fig. 2. Two types of equilibrium saturation profile are observed: (i) the profile uniform $P_u$, when $\chi \le 2$ and (ii) the profile nonuniform $P_n$ where two different values of saturation coexist, when $\chi \ge 2.5$.

These states have been described in Beltrame and Cajot (2022) using continuation of steady-states with the mean saturation as bifurcation parameter. Uniform equilibrium $P_u$ becomes unstable when $\Pi_u(S)$ is no longer monotonous, i.e. if $\chi > \chi_i$. If the mean saturation $S_m$ is in a range where $\Pi_u(S)$ is increasing, then uniform state becomes linearly unstable. As numerically shown in Beltrame and Cajot (2022), the nonuniform equilibrium $P_n$ not only exists for the instability condition of $P_u$, but may coexist with the uniform equilibrium for a mean saturation where $P_u$ is stable. Then, we have a bi-stability state with a hysteresis phenomenon.

Without numerical investigation, the precedent result can be interpreted by analogy with the Van der Waals (VdW) isotherms for the liquid-vapor change. The state curve gives the pressure as a function of volume for a temperature. For a critical temperature, the VdW isotherm is no longer monotonous. The single gas phase is then no longer stable and liquid and gas coexist. The fraction of gas and liquid is given by the Maxwell construction (de Gennes et al., 2004; Doi, 2013).

In the soil, the curve $\Pi_u(s)$, for a value of $\chi$, plays the role of state curve of the VdW isotherm at temperature $T$. The order parameter is then $\chi$ instead $T$. The gas or liquid states of the VdW equilibria are associated with wet and dry regions in the soil. A necessary condition for two different saturations to coexist is pressure equilibrium.

Therefore, for a decreasing curve, only one state can exist, leading to the uniform profile state $P_u$; it is the case (i).

If the curve is no longer monotonic, then for a given pressure interval, there are three possible saturation values. For instance, for $\chi = 6$ in Fig. 4, the horizontal line associated to the dimensionless expression of the pressure $\Pi_u$ equal to $-0.094$ intersects the curve in three points $S_1 \approx 0.092 < S_2 \approx 0.396 < S_3 \approx 0.956$. Note that $S_3$ corresponds to the
<figure><img src="image_8.png" /><figcaption>Fig. 4. Disjoining pressure in fine sandy soil for $\chi = 6$. The horizontal black dashed line is the dimensionless equilibrium pressure $\Pi_w$ equal to -0.094. The red stars are the cross point between disjoining pressure and the equilibrium pressure at saturation $S_1 \approx 0.092$, $S_2 \approx 0.396$ and $S_3 \approx 0.956$. The green area $A_1$ between $S_1$ and $S_2$ is equal in absolute value to the area $A_2$ between $S_2$ and $S_3$ ($|A_1| = |A_2| = 0.69$).</figcaption></figure>

compressed saturation $S_{comp}$. The equality of pressure is not sufficient to ensure equilibrium, the three saturation values must be such that the areas $A_1$ and $A_2$, formed between the curve and the horizontal line, are equal (see Fig. 4). This later condition of the energy constitutes the Maxwell construction. Saturation $S_2$ corresponds to a uniform but unstable equilibrium. The system will undergo “phase separation”, with wet region of saturation $S_3$ and dry region of saturation $S_1$ coexisting.

The critical $\chi_c$ value for which we go from case (i) to (ii) is the one for which the curve admits an inflection point, i.e. $\chi_c = \chi_l \approx 1.73$. However, this critical value is valid only for an infinite domain. For finite domain, the value $\chi_c$ is slightly larger than $\chi_l$ because the uniform state with $S_2$ saturation can be still stable (Beltrame and Cajot, 2022). That is why we observe a uniform profile for $\chi = 2 > \chi_l$. Moreover, the existence of surface tension $\Gamma \neq 0$ prevents a discontinuity, so the saturation profile is continuous and the saturation plateaus are no longer exactly flat. For $\chi = 6$, the two saturation extrema $S_1 \approx 0.092$ and $S_3 \approx 0.956$ in Fig. 3 are well approximated by Maxwell construction in Fig. 4.

The model reproduces a stoppage of imbibition front for large enough concentration of amphiphilic matter in horizontal homogeneous fine sandy soil. We determine the onset of this phenomenon related to the amplitude $\chi$. The values of saturation of nonuniform equilibria for the case of a stopped imbibition front are approximated by the Maxwell construction. In the following, we focus on a vertical water flow.

## 3.2. Vertical homogeneous fine sandy soil

We consider the 1D model of a vertical water transfer in a sandy soil where we impose a saturation $S = 0.96$ and a suction of 250 Pa at the bottom $z = 0$. At the top, $z = L = 550$ mm, a no-flux condition for the water and the energy is imposed. Then PDE system with the boundary conditions reads:

$$ \phi \partial_t S = \partial_z \left[ K(S) \partial_z \left( \frac{\delta F}{\delta S} - \rho g z \right) \right], \quad (18a) $$

$$ \frac{\delta F}{\delta S} = \partial_S f_w(S, c) - \Gamma \partial_{zz}^2 S, \quad (18b) $$

$$ S = 0.96 \text{ and } \frac{\delta F}{\delta S} = -250 \text{ Pa at } z = 0, \quad (18c) $$

$$ \frac{\partial S}{\partial n} = 0, \quad \text{at } z = L, \quad (18d) $$

$$ \frac{\partial}{\partial n} \left( \frac{\delta F}{\delta S} \right) = 0, \quad \text{at } z = L. \quad (18e) $$

We study the dynamic of capillary rise in a sandy soil by varying the repellent factor $\chi$ from 0 to 6. Fig. 5 displays the saturation profiles

at $t = 10, 60, 120$ s and at equilibrium for each factor $\chi$. Black dashed line is the capillary height ($H_{cap} = 260$ mm) in the fine sandy soil for $\chi = 0$ according to the Brooks and Corey model (Brooks and Corey, 1966). In the Richards model, there is a continuum from high to low saturation, as shows the equilibrium profile for $\chi = 0$ in Fig. 5(d). Then we arbitrarily define capillary height, the height at which saturation reaches $s_{ch} = 0.91$. This definition is consistent with the capillary height for the case $\chi = 0$ of the Richards model.

Capillary rise dynamics is slower when the repellent factor $\chi$ is greater. At the beginning of the capillary rise, namely $t = 10$ s, the saturation profiles are almost identical, with a delay as $\chi$ increases. The case $\chi = 6$ is an exception: saturation does not exceed 0.1 as soon as $z > 0$ and is less than 0.01 at $z = 30$ mm. In other words, there is no capillary rise because of hydrophobic behavior of the dry amphiphilic matter. At $t = 60$ s, capillary rise is complete, i.e. the height at which saturation reaches 0.91 no longer evolves. Nevertheless, water continues to rise, as shows the profile evolution. For $t \ge 120$ s, there is a clear difference between the profile with $\chi = 0$ and the other with $\chi \ge 2$. While for $\chi = 0$ the curve is convex for $s < 0.8$, all the saturation profiles for $\chi \ge 2$ show a plateau. The height of this plateau decreases when $\chi$ increases. The plateau indicates an abrupt change in water saturation, from moist to almost dry within a small height range. The equilibrium profile ($S, z$) shows a “rebound” for low saturation but the water saturation is less than 0.01, i.e. the soil is dry. For instance, for $\chi = 3$, the plateau occurs in the range [284, 297] mm and the water is less than 0.01 for $z = 361$ mm. This means that there is a sharp separation between the wet and dry zones, unlike in the case of $\chi = 0$. It is analogous of the coexistence of two states: dry and saturated for the horizontal substrate.

According to the definition of the capillary height, the height at which saturation reaches $S_{ch} = 0.91$, we represent in Fig. 6 the capillary height as a function of $\chi$. In general, capillary height is a decreasing function of $\chi$.

For $\chi \le 3$ the capillary height variation is negligible. In the range $3.5 \le \chi \le 5.5$ capillary height decreases almost linearly. Capillary rise vanishes for a critical value between 5.5 and 6.0. This fact can be explained by the disjoining pressure curve as displayed in Fig. 4. According to the Maxwell construction, saturation $S = 0.956$ coexists with saturation in very dry substrate for $\chi = 0.6$. Since we have set $S = 0.96$ at bottom boundary condition, it is very close to the equilibrium saturation of the Maxwell construction and then it can coexist with a dry region and capillary rise does not occur. On the contrary, for lower value of $\chi$, the horizontal line from $\Pi_w(S = 0.96)$ intersects the $\Pi_w$ curve at three points such that $A_2 > A_1$ implying that the equilibrium condition is not met and that there is a potential energy for water to rise.

Note that the curve of Fig. 6 depends on the choice of the saturation $S_{ch}$ criterium to determine the capillary height. For instance, if we define capillary height, the height at which saturation reaches $S_{ch} = 0.6$ instead of 0.91, then the capillary height decreases significantly for $\chi = 0, 2$ and 3. Regardless of the choice of $S_{ch}$, however, capillary height decreases with the hydrophobic parameter $\chi$.

The wettability model based on the competition of attractive and repellent energies allows to reproduce capillary rise and a stoppage of imbibition front in soil without reducing hydraulic conductivity.

## 3.3. Vertical stratified fine sandy soil

The goal of this section is to compare our model with an existing experiment in literature of capillary rise in a stratified medium (Kroener, 2015; Kroener et al., 2014). This will enable the model to be calibrated.
<figure><img src="image_9.png" /></figure>

<figure><img src="image_10.png" /></figure>

<figure><img src="image_11.png" /></figure>

<figure><img src="image_12.png" /><figcaption>Fig. 5. Capillary rise dynamic in homogeneous fine sandy soil with different repellent factor $\chi$. The horizontal black dashed line is the capillary height in this soil with $\chi = 0$. Parameters are given for fine sandy soil in Table 1 and the other are $\rho = 10^3$ kg m$^{-3}$, $g = 9.81$ m s$^{-2}$, $S_r = 0.05$, $\sigma = 0.2$, $\phi = 0.43$, $K_{sat} = 10^{-4}$ m s$^{-1}$, $\Gamma = 1.7$ mN and $\kappa = 200$.</figcaption></figure>

<figure><img src="image_13.png" /><figcaption>Fig. 6. Capillary height in homogeneous fine sandy soil with respect to the repellent factor $\chi$. The horizontal black dashed line is the capillary height in this soil with $\chi = 0$.</figcaption></figure>

<figure><img src="image_14.png" /><figcaption>Fig. 7. Equilibrium capillary rise profile depending to the repellent factor in a stratified fine sandy soil. The first layer ($0 \le z \le 38$ mm) does not contain repellent interaction ($\chi = 0$) whereas the second layer ($z > 38$ mm) contain a $\chi \ne 0$. The horizontal black dashed line delimits the two layers and the black line corresponds to the initial saturation. Parameters are given for fine sandy soil in Table 1 and the other are $\rho = 10^3$ kg m$^{-3}$, $g = 9.81$ m s$^{-2}$, $S_r = 0.05$, $\sigma = 0.2$, $\phi = 0.43$, $K_{sat} = 10^{-4}$ m s$^{-1}$, $\Gamma = 1.7$ mN and $\kappa = 200$.</figcaption></figure>

### 3.3.1. Equilibria in stratified layers

We consider the 1D model of a vertical stratification of two kinds of fine sandy soil such that the heterogeneity depends only on $z$. The vertical domain length is fixed to $L = 55$ mm and we assume the same boundary conditions as those given in Section 3.2.

The density free energy of wettability given by Eq. (2) depends now on $z$. As in Beltrame and Cajot (2022), only the repellent interaction varies and the density free energy reads:

$$f_w(S, z) = f_c(S) + \xi_{het}(z) \rho g H_{cap} \chi f_r(S) + f_{comp}(S) \quad (19)$$

The heterogeneous function $\xi_{het}(z)$ varies from 0 to 1 approximating the step function:

$$\xi_{het}(z) = \frac{1 + \tanh(100(z - 38.5 \text{ mm}))}{2} \quad (20)$$

The function is almost zero in the vertical range $z \in [0, 38.5 \text{ mm}]$ corresponding to negligible repellent interaction in the fine sand layer. The function is almost 1 in the remaining top layer corresponding to a repellent free energy of amplitude $\chi$ for the fine sand layer containing amphiphilic matter. Black dashed line in Fig. 7 represents this heterogeneous function where $z$ is the ordinate axis. Equilibrium profiles are displayed in Fig. 7. For $\chi \le 5.3$, water rises until the domain

boundary and all the medium is saturated. For $\chi = 5.4$ water infiltrates only partially the layer with amphiphilic matter. Then, a sharp limit appears and the remaining layer is almost dry with a saturation about $S = 0.10$. For $\chi \le 5.5$ capillary rises stops at the amphiphilic sand layer. The transition from saturated to dry is sharper when $\chi$ increases. Therefore, there is a threshold value of $\chi$ about 5.5 for which the water does not infiltrate the layer. This threshold is similar to the previous case of homogeneous layer for which for $\chi$ about 6 the capillary rise does not occur. Here, the critical value of $\chi$ is slightly lower because a part of potential energy is lost by the suction at boundary condition and by the gravity during the capillary rise in the sand without amphiphilic matter.

The threshold of imbibition in the amphiphilic layer is sharply defined and then we chose this experiment of capillary rise to accurately calibrate the wettability properties of the model.
<figure><img src="image_15.png" /><figcaption>Fig. 8. Repellent factor $H_{rep}(c) = H_{cap}\chi(c)$ with respect to the concentration c of amphiphilic matter. The fitted curve, based on Eq. (21), is calibrated for a sandy soil containing EPS using fine sand, medium sand and coarse sand. The R² error of the fitting parameters is 0.998 for $\theta = 1418$ mm and $\nu = 1.402$. Other parameters are given in Table 1.</figcaption></figure>

<figure><img src="image_16.png" /><figcaption>Fig. 9. Disjoining pressure in fine sandy soil for different amphiphilic concentrations. Parameters are $\rho = 10^3$ kg m⁻³, $g = 9.81$ m s⁻², $H_{cap} = 0.26$ m, $S_r = 0.05$, $\sigma = 0.2$, $\lambda = 0.972$, $\kappa = 200$, $\theta = 1418$ mm and $\nu = 1.402$.</figcaption></figure>

### 3.3.2. Comparison with capillary rise experiment and model calibration

In this section, we consider the experimental set-up presented in Kroener (2015) and Kroener et al. (2014). The experiment consists of placing a 5 mm layer of sand mixed with EPS (analogous to a rhizosphere) between two layers of dry clean sand (analogous to the bulk soil, EPS free in rhizosphere and root). The samples were placed on suction tables with a suction of 2.5 cm (hydraulic head) in order to impose a constant pressure at their base. Experiments were carried out on samples of sand with different grain size: fine (0.125–0.2 mm grain diameter), medium (0.2–0.355 mm) and coarse (0.63–1.0 mm).

This experiment highlights a concentration threshold $c_{crit}$ where water is not able to pass the sand and EPS stratum during a certain period when the concentration is higher than $c_{crit}$, after 5 min of capillary rise, then passes by the border of the sample. The phenomenon of capillary rise stalling is reminiscent of the dynamics presented in Section 3.3.2. Thus, the repellent factor $\chi$ can be calibrated from this experiment such that the critical concentration of the model matches the experiment.

The critical concentration depends on the grain size of the sand ($c_{crit} = 1$ for fine sand, 0.75 for medium sand and 0.05 for coarse sand) and we fit the height $H_{rep}$ such as the threshold appears. In this way, we extract a couple ($H_{rep}(c), c$) for each medium. These couples are displayed in red dot in Fig. 8. In addition, we assume that $H_{rep}$ follows a power law:

$$H_{rep}(c) = \theta c^{\nu}. \qquad (21)$$

The fitted curve, blue line in Fig. 8, leads to the values of the prefactor $\theta = 1418$ mm and the exponent $\nu = 1.402$ where the R² error of this fit is 0.998. Finally, we deduce the function $\chi(c)$ from $H_{rep}(c)$ using Eq. (7).

Inserting this amplitude factor in Eq. (9a), the disjoining pressure $\Pi_w$ is defined with respect to the EPS concentration c. Fig. 9 displays the disjoining pressure for fine sand with $c = \{0, 0.1, 0.6, 0.71, 0.8, 1\}$. This figure is similar to Fig. 2, only the parametrization with c instead of $\chi$ changes. The curve $\Pi_w(s, c)$ is no longer monotonic for $c > 0.44$. For higher concentration, the curve $\Pi_w(s, c)$ increases on the interval $[S_m, S_M]$. The saturation $S_m$ and $S_M$ vary with respect to the value of c. The repellent effect increases at low saturation as c increases, namely the range $[S_m, S_M]$ and the amplitude of extrema increase. For example, for $c = 0.8$ the dimensionless disjoining pressure at $S_m(0.8) = 0.214$ is $\Pi_w(S_m(0.8), 0.8) = -0.815$ whereas for $c = 0$ at the same saturation the pressure is much higher $\Pi_w(S_m(0.8), 0) = 5.119$.

To model capillary rise and buffering effect due to the amphiphilic sand layer, only the first two layers are considered in the model. Namely, the first layer of pure sand, 38 mm high, overlaid by a layer of sand-EPS mixture, 17 mm high, with a total domain length $L = 55$ mm. We apply Dirichlet boundary condition at $z = 0$ and homogeneous Neumann boundary condition at $z = L$.

The heterogeneous function, defined in Section 3.3.1 (Eq. (20)), displayed in dashed line in Fig. 10, is almost equal to 1 in the region

containing EPS otherwise it is almost zero.

The numerical simulation of the capillary rise dynamics, for different sandy soil with the fitted $\chi$ function for a range of different EPS concentration, is shown in Fig. 10. All simulations are performed using the same initial profile indicated with a black line in Fig. 10.

In fine sand, water rises in both layers for an EPS concentration strictly below 1, and there is no noticeable difference in water saturation between the two layers at equilibrium. For $c = 1$, the water rises a few millimeters in the EPS layer and then stops. At equilibrium, there is only about 10% water in the EPS layer (Fig. 10(a)).

For a medium sand (Fig. 10(b)) a similar behavior occurs but for a lower value $c_{crit} = 0.75$. Note that the capillary rise at the interface is even slightly higher, whereas for fine sand, the capillary rise was a few millimeters above the interface.

For coarse sand, the saturation profiles differ because the capillary height ($H_{cap} \approx 10$ mm) is much lower than the height of the interface (38.5 mm). Fig. 10(c) shows that the soil is almost saturated over the capillary height ($\approx 10$ mm) and then decreases gradually. As long as the water is far from the layers interface, all water contents are similar for a given z. Thus, at $z = 30$ mm, the saturation reaches $S = 0.5$ regardless of concentration. In the EPS layer, an infiltration threshold is still observed. For concentrations above 0.05, saturation does not exceed 0.05 at the interface. For $c = 0.05$ and $c = 0$, saturation is around 0.4 at the top of the EPS layer. Saturation slightly decreases between $c = 0$ and $c = 0.05$. This non-negligible effect of EPS can be explained by the fact that the disjoining pressure curves differ slightly at this saturation range.

For any grain size, a critical concentration determines an infiltration threshold in the EPS layer. The critical concentration of EPS increases for smaller grain sizes. This result is consistent with the fact that sand grains are associated with polar, i.e. hydrophilic, interactions, so decreasing grain size increases hydrophilic interaction energy per unit volume. Whereas hydrophobic interactions, essentially due to amphiphilic matter, do not vary with grain size. So, for finer sand grain sizes, hydrophilic interactions are relatively more important.

To conclude, the model of amphiphilic soil, defined by Eq. (13) and based on wettability energy depending on c, highlights a phenomenon of capillary rise threshold and it is able to reproduce both scenarios of capillary rise and the same dependence of $c_{crit}$ on grain size than the literature experiment.

## 4. Conclusion

In this paper, we have proposed a model of out-of-equilibrium water transfer in an amphiphilic soil. The modeling is consistent with a thermodynamic description of the wettability of a porous medium in presence of hydrophobic matter (Beltrame and Cajot, 2022). For
<figure><img src="image_17.png" /><figcaption>**Fig. 10.** Equilibrium profiles of water saturation after capillary rise in sandy stratified media made up of three different grain sizes: (a) fine grain, (b) medium grain and (c) coarse grain. The first layer of each sand medium ($0 \le z \le 38$ mm) does not contain EPS while the second layer ($z > 38$ mm) is a mixture of sand containing EPS at different concentrations $c$ indicated in legend. The horizontal black dashed line delimits the two layers and the black line corresponds to the initial saturation. The capillary height and exponent (pore size distribution) factors of the Brooks and Corey law are indicated in each panel. Other parameters are $\rho = 10^3$ kg m$^{-3}$, $g = 9.81$ m s$^{-2}$, $S_r = 0.05$, $\sigma = 0.2$, $\phi = 0.43$, $K_{sat} = 10^{-4}$ m s$^{-1}$, $\Gamma = 1.7$ mN and $\kappa = 200$.</figcaption></figure>

a better account the amphiphilic matter effect on water-soil interactions, we propose a paradigm shift in the modeling of soil wettability. Instead of associating only a contact angle with wettability, we define a functional energy representing the effective interactions between the soil matrix (Beltrame and Cajot, 2022) and the free water surface (Cueto-Felgueroso and Juanes, 2009b; Beljadid et al., 2020). Defining a hydraulic conductivity and a dissipation energy term, water flow dynamics is derived from the generic Onsager variational principle leading to a non-linear fourth order PDE on the saturation. This equation generalizes the Darcy-Richards equation (Richards, 1931) and its improvement (Beljadid et al., 2020) where notably a disjoining pressure term derived from the wettability free energy replaces the capillary pressure.

This new formulation enables to simulate phenomena specific to amphiphilic wettability properties. A first phenomenon is the simulation of the stoppage of imbibition front in a homogeneous soil appearing for a large enough concentration of amphiphilic matter: at equilibrium, a moist soil coexists with a dry soil. This type of solution is made possible with this model by the fact that the disjoining pressure is no longer a monotonic function of saturation due to the hydrophobic property of amphiphilic matter at low saturation. This coexistence of wet/dry states is analogous to the liquid-vapor equilibrium for a Van der Waals gas. Secondly, for low amphiphilic matter concentration, numerical simulation demonstrates a slower capillary rise dynamics and a decrease in capillary height in comparison with a sand without amphiphilic matter. It is remarkable that the equilibrium profile displays a clear distinction between the region close to saturation and the dry zone in contrast to the equilibrium profile modeled by the Richards equation with the Brooks and Corey model. As a third result, we have captured a threshold of amphiphilic concentration for which the capillary rise does not occur in a dry homogeneous amphiphilic medium or is stopped at the dry layer containing the amphiphilic matter without reducing hydraulic conductivity. Comparison between simulations and experiments of capillary rise in fine, medium and coarse sand enabled the calibration of the relationship between amphiphilic matter concentration and the repellent pressure term in the disjoining pressure. Once the model has been calibrated, the numerical simulation of capillary rise is in agreement with experimental observations (Kroener, 2015).

Delayed capillary rise or delayed rehumectation have been interpreted differently in the literature devoted to the rhizosphere (Ahmed et al., 2014; Kroener et al., 2014, 2016, 2018). In presence of amphiphilic matter, hydraulic conductivity decreases and a relaxation term representing the out-of-equilibrium rehumectation dynamics is added. The slowdown effect is interpreted as being mainly due to the low conductivity (Kroener et al., 2016). In contrast, in our model, only wettability is responsible for the delay. Both interpretations are

relevant.

Nevertheless, the experiment presented in Kroener (2015) shows a threshold concentration for which the capillary rise is stopped after 5 minutes then continues to rise slowly in the EPS layer. A low conductivity cannot be the reason for stopping capillary rise or the buffering effect. This demonstrates that wettability plays a crucial role in infiltration dynamics.

Furthermore, saturated hydraulic conductivity measurements do not show a systematic decrease with the increase in amphiphilic matter concentration (Kroener et al., 2018). In order to clearly distinguish the respective roles of wettability and hydraulic conductivity, further experiments need to be carried out.

In this article, we have neglected the osmotic pressure, which occurs in nearly saturated porous media with amphiphilic matter. To take into account this phenomenon, a hydrophilic component was added in the free energy applied to a model of rhizosphere with EPS (Cajot, 2024). Numerical simulation of rehumectation in the neighborhood of the root demonstrated an impact on both the rehumectation duration and the final water content.

Our construction of the free energy in amphiphilic soil based on wettability energy provides a new way to understand water transfer in presence of amphiphilic matter. It is remarkable that our model requires the surface tension term introduced by Cueto-Felgueroso and Juanes (2009a), Beljadid et al. (2020) and Moure et al. (2023) for front instability in homogeneous soil without organic matter. This suggests that this term should be taken into account in soil infiltration models. Similarly, the hydrophobicity term is not only useful in the context of amphiphilic matter but is also a determining factor in the stability of fronts in stratified media (see Beltrame and Cajot, 2022). Our approach with a new definition of wettability is therefore promising for modeling environments containing amphiphilic matter and, more generally, for understanding the impact of organic matter on water transfer such as drop infiltration into a soil (Cajot et al., 2025).

### CRediT authorship contribution statement

**Florian Cajot:** Writing – review & editing, Writing – original draft, Validation, Software, Methodology, Formal analysis, Conceptualization. **Claude Doussan:** Writing – review & editing, Supervision, Funding acquisition. **Philippe Beltrame:** Writing – review & editing, Writing – original draft, Supervision, Funding acquisition, Conceptualization.

### Declaration of competing interest

The authors declare that they have no known competing financial interests or personal relationships that could have appeared to influence the work reported in this paper.
## Acknowledgments

The authors acknowledge support by INRAE and Avignon University which provided the PhD grant of FC. FC and PB acknowledge fruitful discussions with Uwe Thiele and greatly appreciate his valuable advices, aided by a Perdiguier program awarded to FC by Avignon Université. The authors acknowledge fruitful discussion with Laurent Legendre and Laurent Limat and greatly appreciate their valuable advices during the scientist expert meeting of FC. The authors thank Annette Bérard, Clémentine Lapie, Léo Corridor, Emy Garcia and Arnaud Chapelet for technical assistance and help in the design of experiments of drop infiltration in an amphiphilic sandy medium.

## Data availability

No data was used for the research described in the article.

## References

Ahmed, M.A., Holz, M., Woche, S.K., Bachmann, J., Carminati, A., 2015. Effect of soil drying on mucilage exudation and its water repellency: a new method to collect mucilage. J. Plant Nutr. Soil Sci. 178 (6), 821-824. http://dx.doi.org/10.1002/jpln.201500177, URL: https://onlinelibrary.wiley.com/doi/abs/10.1002/jpln.201500177, _eprint: https://onlinelibrary.wiley.com/doi/pdf/10.1002/jpln.201500177.

Ahmed, M.A., Kroener, E., Holz, M., Zarebanadkouki, M., Carminati, A., 2014. Mucilage exudation facilitates root water uptake in dry soils. Funct. Plant Biol. 41 (11), 1129. http://dx.doi.org/10.1071/FP13330, URL: http://www.publish.csiro.au/?paper=FP13330.

Alghunaim, A., Kirdponpattara, S., Newby, B.-m.Z., 2016. Techniques for determining contact angle and wettability of powders. Powder Technol. 287, 201-215. http://dx.doi.org/10.1016/j.powtec.2015.10.002, URL: https://www.sciencedirect.com/science/article/pii/S0032591015300917.

Bamorovat Abadi, G., Bahrami, M., 2022. The effect of surface roughness on capillary rise in micro-grooves. Sci. Rep. 12 (1), 14867. http://dx.doi.org/10.1038/s41598-022-19111-w.

Bear, J., 1988. Dynamics of Fluids in Porous Media. In: Dover Books on Physics and Chemistry, Dover, New York.

Belaiev, A.Y., Hassanizadeh, S.M., 2001. A theoretical model of hysteresis and dynamic effects in the capillary relation for two-phase flow in porous media. Transp. Porous Media 43 (3), 487-510. http://dx.doi.org/10.1023/A:1010736108256.

Beljadid, A., Cueto-Felgueroso, L., Juanes, R., 2020. A continuum model of unstable infiltration in porous media endowed with an entropy function. Adv. Water Resour. 144, 103684. http://dx.doi.org/10.1016/j.adwatres.2020.103684, URL: https://linkinghub.elsevier.com/retrieve/pii/S0309170819308954.

Beltrame, P., Cajot, F., 2022. Model of hydrophobic porous media applied to stratified media: Water trapping, intermittent flow and fingering instability. Europhys. Lett. 138 (5), 53004. http://dx.doi.org/10.1209/0295-5075/ac71c1, Publisher: EDP Sciences, IOP Publishing and Società Italiana di Fisica.

Benard, P., Zarebanadkouki, M., Carminati, A., 2019. Physics and hydraulics of the rhizosphere network. J. Plant Nutr. Soil Sci. 182 (1), 5-8. http://dx.doi.org/10.1002/jpln.201800042, URL: https://onlinelibrary.wiley.com/doi/abs/10.1002/jpln.201800042, _eprint: https://onlinelibrary.wiley.com/doi/pdf/10.1002/jpln.201800042.

Bérard, A., Clavel, T., Le Bourvellec, C., Davoine, A., Le Gall, S., Doussan, C., Bureau, S., 2020. Exopolysaccharides in the rhizosphere: A comparative study of extraction methods. Application to their quantification in Mediterranean soils. Soil Biol. Biochem. 149, 107961. http://dx.doi.org/10.1016/j.soilbio.2020.107961, URL: https://www.sciencedirect.com/science/article/pii/S0038071720302571.

Bisdom, E.B.A., Dekker, L.W., Schoute, J.F.T., 1993. Water repellency of sieve fractions from sandy soils and relationships with organic material and soil structure. In: Brussaard, L., Kooistra, M.J. (Eds.), Soil Structure/Soil Biota Inter-relationships. Elsevier, Amsterdam, pp. 105-118. http://dx.doi.org/10.1016/B978-0-444-81490-6.50013-3, URL: https://www.sciencedirect.com/science/article/pii/B9780444814906500133.

Brooks, R.H., Corey, A.T., 1966. Properties of porous media affecting fluid flow. J. Irrig. Drain. Div. 92 (2), 61-88. http://dx.doi.org/10.1061/JRCEA4.0000425, URL: https://ascelibrary.org/doi/10.1061/JRCEA4.0000425. Publisher: American Society of Civil Engineers.

Cahn, J.W., Hilliard, J.E., 1958. Free energy of a nonuniform system. i. interfacial free energy. J. Chem. Phys. 28 (2), 258-267. http://dx.doi.org/10.1063/1.1744102.

Cajot, F., 2024. Modeling Water Transfer in Soil in the Presence of Amphiphilic Matter: Application to the Rhizosphere (These de doctorat). Avignon, Avignon.

Cajot, F., Doussan, C., Hartmann, S., Beltrame, P., 2025. Model of drop infiltration into a thin amphiphilic porous medium. J. Colloid Interface Sci. 684, 35-46. http://dx.doi.org/10.1016/j.jcis.2024.12.216, URL: https://www.sciencedirect.com/science/article/pii/S0021979724030704.

Carminati, A., Moradi, A.B., Vetterlein, D., Vontobel, P., Lehmann, E., Weller, U., Vogel, H.-J., Oswald, S.E., 2010. Dynamics of soil water content in the rhizosphere. Plant Soil 332 (1-2), 163-176. http://dx.doi.org/10.1007/s11104-010-0283-8, URL: http://link.springer.com/10.1007/s11104-010-0283-8.

Cueto-Felgueroso, L., Juanes, R., 2008. Nonlocal interface dynamics and pattern formation in gravity-driven unsaturated flow through porous media. Phys. Rev. Lett. 101 (24), 244504. http://dx.doi.org/10.1103/PhysRevLett.101.244504, URL: https://link.aps.org/doi/10.1103/PhysRevLett.101.244504. Publisher: American Physical Society.

Cueto-Felgueroso, L., Juanes, R., 2009a. Adaptive rational spectral methods for the linear stability analysis of nonlinear fourth-order problems. J. Comput. Phys. 228 (17), 6536-6552. http://dx.doi.org/10.1016/j.jcp.2009.05.045, URL: https://linkinghub.elsevier.com/retrieve/pii/S0021999109003179.

Cueto-Felgueroso, L., Juanes, R., 2009b. A phase field model of unsaturated flow. Water Resour. Res. 45 (10), http://dx.doi.org/10.1029/2009WR007945, URL: https://onlinelibrary.wiley.com/doi/abs/10.1029/2009WR007945, _eprint: https://onlinelibrary.wiley.com/doi/pdf/10.1029/2009WR007945.

Cueto-Felgueroso, L., Juanes, R., 2009c. Stability analysis of a phase-field model of gravity-driven unsaturated flow through porous media. Phys. Rev. E 79 (3), 036301. http://dx.doi.org/10.1103/PhysRevE.79.036301, URL: https://link.aps.org/doi/10.1103/PhysRevE.79.036301.

de Gennes, P.G., 1985. Wetting: statics and dynamics. Rev. Modern Phys. 57 (3), 827-863. http://dx.doi.org/10.1103/RevModPhys.57.827, URL: https://link.aps.org/doi/10.1103/RevModPhys.57.827. Publisher: American Physical Society.

de Gennes, P.-G., Brochard-Wyart, F., Quéré, D., 2004. Wetting and long-range forces. In: de Gennes, P.-G., Brochard-Wyart, F., Quéré, D. (Eds.), Capillarity and Wetting Phenomena: Drops, Bubbles, Pearls, Waves. Springer, New York, NY, pp. 87-105. http://dx.doi.org/10.1007/978-0-387-21656-0-4.

Derjaguin, B.V., Churaev, N.V., Muller, V.M., 1987. Surface Forces. Springer US, Boston, MA, http://dx.doi.org/10.1007/978-1-4757-6639-4, URL: http://link.springer.com/10.1007/978-1-4757-6639-4.

DiCarlo, D.A., 2004. Experimental measurements of saturation overshoot on infiltration. Water Resour. Res. 40 (4), http://dx.doi.org/10.1029/2003WR002670, URL: https://onlinelibrary.wiley.com/doi/abs/10.1029/2003WR002670, _eprint: https://onlinelibrary.wiley.com/doi/pdf/10.1029/2003WR002670.

DiCarlo, D.A., 2007. Capillary pressure overshoot as a function of imbibition flux and initial water content. Water Resour. Res. 43 (8), http://dx.doi.org/10.1029/2006WR005550, URL: https://onlinelibrary.wiley.com/doi/abs/10.1029/2006WR005550, _eprint: https://onlinelibrary.wiley.com/doi/pdf/10.1029/2006WR005550.

DiCarlo, D.A., 2013. Stability of gravity-driven multiphase flow in porous media: 40 Years of advancements. Water Resour. Res. 49 (8), 4531-4544. http://dx.doi.org/10.1002/wrcr.20359, URL: https://onlinelibrary.wiley.com/doi/abs/10.1002/wrcr.20359, _eprint: https://onlinelibrary.wiley.com/doi/pdf/10.1002/wrcr.20359.

Doi, M., 2011. Onsager's variational principle in soft matter. J. Phys.: Condens. Matter. 23 (28), 284118. http://dx.doi.org/10.1088/0953-8984/23/28/284118.

Doi, M., 2013. Soft Matter Physics, first ed. Oxford University Press, Oxford ; New York, OCLC: ocn857650756.

Doussan, C., Gall, S.L., Ruy, S., Bérard, A., 2024. Utiliser les racines pour moduler les impacts des déficits en eau et améliorer la gestion de l'eau dans les agrosystèmes. Sci. Eau Territ. (45), http://dx.doi.org/10.20870/Revue-SET.2024.45.8144, 8144-8144. URL: https://revue-set.fr/article/view/8144. Number: 45.

Egorov, A.G., Dautov, R.Z., Nieber, J.L., Sheshukov, A.Y., 2003. Stability analysis of gravity-driven infiltrating flow. Water Resour. Res. 39 (9), http://dx.doi.org/10.1029/2002WR001886, URL: https://onlinelibrary.wiley.com/doi/abs/10.1029/2002WR001886, _eprint: https://onlinelibrary.wiley.com/doi/pdf/10.1029/2002WR001886.

Eliassi, M., Glass, R.J., 2002. On the porous-continuum modeling of gravity-driven fingers in unsaturated materials: Extension of standard theory with a hold-back-pile-up effect. Water Resour. Res. 38 (11), 16-1-16-11. http://dx.doi.org/10.1029/2001WR001131.

Emmerich, H., 2008. Advances of and by phase-field modelling in condensed-matter physics. Adv. Phys. 57 (1), 1-87. http://dx.doi.org/10.1080/00018730701822522, Publisher: Taylor & Francis _eprint: https://doi.org/10.1080/00018730701822522.

Flemming, H.-C., van Hullebusch, E.D., Neu, T.R., Nielsen, P.H., Seviour, T., Stoodley, P., Wingender, J., Wuerzt, S., 2023. The biofilm matrix: multitasking in a shared space. Nat. Rev. Microbiol. 21 (2), 70-86. http://dx.doi.org/10.1038/s41579-022-00791-0, URL: https://www.nature.com/articles/s41579-022-00791-0. Publisher: Nature Publishing Group.

Heib, F., Schmitt, M., 2016. Statistical contact angle analyses with the high-precision drop shape analysis (HPDSA) approach: Basic principles and applications. Coatings 6 (4), 57. http://dx.doi.org/10.3390/coatings6040057, URL: https://mdpi.com/2079-6412/6/4/57. Number: 4 Publisher: Multidisciplinary Digital Publishing Institute.
Heil, M., Hazel, A.L., 2006. Oomph-lib – an object-oriented multi-physics finite-element library. In: Bungartz, H.-J., Schäfer, M. (Eds.), Fluid-Structure Interaction. In: Lecture Notes in Computational Science and Engineering, Springer, Berlin, Heidelberg, pp. 19–49. http://dx.doi.org/10.1007/3-540-34596-5_2.

Hinsinger, P., Bengough, A.G., Vetterlein, D., Young, I.M., 2009. Rhizosphere: bio-physics, biogeochemistry and ecological relevance. Plant Soil 321 (1), 117–152. http://dx.doi.org/10.1007/s11104-008-9885-9.

Hird, R., Bolton, M.D., 2017. Clarification of capillary rise in dry sand. Eng. Geol. 230, 77–83. http://dx.doi.org/10.1016/j.enggeo.2017.09.023, URL: https://www.sciencedirect.com/science/article/pii/S0013795217302302.

Intergovernmental Panel on Climate Change (IPCC), 2014. Climate Change 2013 – The Physical Science Basis: Working Group I Contribution to the Fifth Assessment Report of the Intergovernmental Panel on Climate Change. Cambridge University Press, Cambridge, http://dx.doi.org/10.1017/CBO9781107415324, URL: https://www.cambridge.org/core/books/climate-change-2013-the-physical-science-basis/BE9453E500DEF3640B383BADDC332C3E.

Knott, M., Ani, M., Kroener, E., Diehl, D., 2022. Effect of changing chemical environment on physical properties of maize root mucilage. Plant Soil 478 (1), 85–101. http://dx.doi.org/10.1007/s11104-022-05577-0.

Kroener, E., 2015. Roots at the percolation threshold. Phys. Rev. E 91 (4), http://dx.doi.org/10.1103/PhysRevE.91.042706.

Kroener, E., Holz, M., Zarebanadkouki, M., Ahmed, M., Carminati, A., 2018. Effects of mucilage on rhizosphere hydraulic functions depend on soil particle size. Vadose Zone J. 17 (1), 170056. http://dx.doi.org/10.2136/vzj2017.03.0056, URL: https://onlinelibrary.wiley.com/doi/abs/10.2136/vzj2017.03.0056, _eprint: https://onlinelibrary.wiley.com/doi/pdf/10.2136/vzj2017.03.0056.

Kroener, E., Zarebanadkouki, M., Bittelli, M., Carminati, A., 2016. Simulation of root water uptake under consideration of nonequilibrium dynamics in the rhizosphere. Water Resour. Res. 52 (8), 5755–5770. http://dx.doi.org/10.1002/2015WR018579, URL: http://doi.wiley.com/10.1002/2015WR018579.

Kroener, E., Zarebanadkouki, M., Kaestner, A., Carminati, A., 2014. Nonequilibrium water dynamics in the rhizosphere: How mucilage affects water flow in soils. Water Resour. Res. 50 (8), 6479–6495. http://dx.doi.org/10.1002/2013WR014756, URL: https://onlinelibrary.wiley.com/doi/abs/10.1002/2013WR014756, _eprint: https://onlinelibrary.wiley.com/doi/pdf/10.1002/2013WR014756.

Le Gall, S., Bérard, A., Page, D., Lanoe, L., Bertin, N., Doussan, C., 2021. Increased exopolysaccharide production and microbial activity affect soil water retention and field performance of tomato under water deficit. Rhizosphere 19, 100408. http://dx.doi.org/10.1016/j.rhisp.2021.100408, URL: https://www.sciencedirect.com/science/article/pii/S245221982100104X.

Leverett, M., 1941. Capillary behavior in porous solids. Trans. AIME 142 (01), 152–169. http://dx.doi.org/10.2118/941152-G.

Mandal, D., Jayaprakash, J., 2009. Water repellency of soils in the lower Himalayan regions of India: impact of land use. Current Sci. 96 (1), 148–152, URL: https://www.jstor.org/stable/24104743. Publisher: Current Science Association.

Moure, A., Jones, N., Pawlak, J., Meyer, C., Fu, X., 2023. A thermodynamic nonequilibrium model for preferential infiltration and refreezing of melt in snow. Water Resour. Res. 59 (5), e02022WR034035. http://dx.doi.org/10.1029/2022WR034035, URL: https://onlinelibrary.wiley.com/doi/abs/10.1029/2022WR034035, _eprint: https://onlinelibrary.wiley.com/doi/pdf/10.1029/2022WR034035.

Mualem, Y., 1976. A new model for predicting the hydraulic conductivity of unsaturated porous media. Water Resour. Res. 12 (3), 513–522. http://dx.doi.org/10.1029/WR0121i003p00513, URL: https://onlinelibrary.wiley.com/doi/abs/10.1029/WR0121i003p00513, _eprint: https://onlinelibrary.wiley.com/doi/pdf/10.1029/WR0121i003p00513.

Müehl, G.J.H., Rüehlmann, J., Goebel, M.-O., Bachmann, J., 2012. Application of confocal laser scanning microscopy (CLSM) to visualize the effect of porous media wettability on unsaturated pore water configuration. J. Soils Sediments 12 (1), 75–85. http://dx.doi.org/10.1007/s11368-011-0395-7.

Nazari, M., 2021. Plant mucilage components and their functions in the rhizosphere. Rhizosphere 18, 100344. http://dx.doi.org/10.1016/j.rhisp.2021.100344, URL: https://www.sciencedirect.com/science/article/pii/S2452219821000409.

Nieber, J.L., Dautov, R.Z., Egorov, A.G., Sheshukov, A.Y., 2005. Dynamic capillary pressure mechanism for instability in gravity-driven flows; review and extension to very dry conditions. Transp. Porous Media 58 (1), 147–172. http://dx.doi.org/10.1007/s11242-004-5473-5.

Rabbani, H.S., Zhao, B., Juanes, R., Shokri, N., 2018. Pore geometry control of apparent wetting in porous media. Sci. Rep. 8 (1), 15729. http://dx.doi.org/10.1038/s41598-018-34146-8.

Rätz, A., Schweizer, B., 2012. Hysteresis models and gravity fingering in porous media. ZAMM - J. Appl. Math. Mech. http://dx.doi.org/10.1002/zamm.201200052, Publisher: WILEY-VCH Verlag.

Richards, L.A., 1931. Capillary conduction of liquids through porous mediums. Physics 1 (5), 318–333. http://dx.doi.org/10.1063/1.1745010, URL: http://aip.scitation.org/doi/10.1063/1.1745010.

Rodríguez-Iturbe, I., 2000. Ecohydrology: A hydrologic perspective of climate-soil-vegetation dynamics. Water Resour. Res. 36 (1), 3–9. http://dx.doi.org/10.1029/1999WR900210, URL: https://onlinelibrary.wiley.com/doi/abs/10.1029/1999WR900210, _eprint: https://onlinelibrary.wiley.com/doi/pdf/10.1029/1999WR900210.

Shi, S., Nuccio, E.E., Shi, Z.J., He, Z., Zhou, J., Firestone, M.K., 2016. The interconnected rhizosphere: High network complexity dominates rhizosphere assemblages. Ecol. Lett. 19 (8), 926–936. http://dx.doi.org/10.1111/ele.12630, URL: https://onlinelibrary.wiley.com/doi/abs/10.1111/ele.12630, _eprint: https://onlinelibrary.wiley.com/doi/pdf/10.1111/ele.12630.

Sinha, R., Fritschi, F.B., Zandalinas, S.I., Mittler, R., 2021. The impact of stress combination on reproductive processes in crops. Plant Sci. 311, 111007. http://dx.doi.org/10.1016/j.plantsci.2021.111007, URL: https://www.sciencedirect.com/science/article/pii/S016894522100203X.

van Genuchten, M.T., 1980. A closed-form equation for predicting the hydraulic conductivity of unsaturated soils. Soil Sci. Am. J. 44 (5), 892–898. http://dx.doi.org/10.2136/sssaj1980.03615995004400050002x, URL: https://onlinelibrary.wiley.com/doi/abs/10.2136/sssaj1980.03615995004400050002x, _eprint: https://onlinelibrary.wiley.com/doi/pdf/10.2136/sssaj1980.03615995004400050002x.

Wang, Z., Wallach, R., 2022. Modeling gravity-driven unstable flow in subcritical water-repellent soils with a time-dependent contact angle. Water Resour. Res. 58 (6), e2021WR031859. http://dx.doi.org/10.1029/2021WR031859, URL: https://onlinelibrary.wiley.com/doi/abs/10.1029/2021WR031859, _eprint: https://onlinelibrary.wiley.com/doi/pdf/10.1029/2021WR031859.