### 5.4.4 Conclusion

The infiltration model developed in section 5.1 reproduces different types of infiltration corresponding to the literature on powdery media [82, 73, 228, 96, 177]. It would be interesting to compare this model with spreading dynamics on brush polymers [166], or with other empirical laws such as Tanner's law [320].

Nevertheless, this model offers a new opportunity to better understand the influence of amphiphilic matter on infiltration phenomena in a porous medium. In addition, it illustrates that this transfer model in the presence of amphiphilic matter is relevant to the study of water transfer behaviors in the rhizosphere in the presence of a root.

## 5.5 Water transfer around a single root

The aim of this section is to model water transfer between bulk soil and a root, passing through the rhizosphere. To achieve this, we propose to consider a horizontal soil fraction with a root where we distinguish: bulk and rhizosphere soil (white dashed line in Figure 5.16). The rhizosphere soil contain an EPS concentration from plant due to root activity unlike the bulk soil. However, the bulk soil can contain an Extracellullar Polymeric Substances concentration due to microorganism in stress [254].

<figure><img src="image_105.png" /><figcaption>Figure 5.16: Sketch of the considered geometry for soil around a single root.</figcaption></figure>

### 5.5.1 Model of soil around a single root

In this study, we assume that Extracellullar Polymeric Substances are similar to the EPS effect. Thus, the EPS concentration in the bulk soil is lower than a critical concentration $C_{crit}$. This concentration represents a boundary between bulk and rhizosphere. In this study, we define the rhizosphere from root radius $r_{root}$ to the rhizosphere radius $r_{rhizo}$ and a bulk soil from the rhizosphere radius $r_{rhizo}$ to the bulk radius $r_{bulk}$. The rhizosphere contains a concentration of EPS that varies with the distance between the soil and the root. Then we consider an EPS concentration $c$ depending on the space $r$.

As the EPS concentration is decreasing from the root exuded concentration $c_{max}$ to the microorganism in stress $C_{crit}$ in the rhizosphere, EPS concentration decreases with increasing distance from the root (grey area in Figure 5.16). Therefore, the concentration distribution $c(r)$ in the rhizosphere
needs to satisfy:

$$c(r) = \begin{cases} c_{max}, & \text{if } r = r_{root}, \\ \ge C_{crit}, & \text{if } r_{root} < r \le r_{rhizo}, \\ \le C_{crit}, & \text{else.} \end{cases} \qquad (5.38)$$

Landl et al. [204] propose to represent this distribution in the soil depending on the age of the root. Their expression should be interesting to be considered when we study all the root. Here we only focus on a local root area, i.e we assume that the age of the root is fixed. Moreover, in this study we consider a simple situation and neglect diffusion of the concentration. This strong assumption allows us to study the behavior of the water with out time dependency of the concentration. Kroener et al. [200] propose in this situation to represent the distribution of the EPS in the soil with a decreasing exponential function

$$c(r) = A \exp(-ar), \qquad (5.39)$$

where A is an amplitude coefficient and a is a range action. These parameters are calibrated using Eq. (5.38), such that

$$A \exp(-ar_{root}) = c_{max}, \qquad (5.40a)$$

$$A \exp(-ar_{rhizo}) = C_{crit}. \qquad (5.40b)$$

Eq. (5.40a) give an expression of a according to A,

$$a = \ln \left[ \left( \frac{A}{C_{crit}} \right)^{\frac{1}{r_{rhizo}}} \right]. \qquad (5.41)$$

Inserting this equation in Eq. (5.40b), the amplitude coefficient A is given by

$$A \exp(-ar_{rhizo}) = C_{crit} \frac{r_{root}}{r_{rhizo}} A^{1-\frac{r_{root}}{r_{rhizo}}}. \qquad (5.42)$$

By introducing the variable $b = \frac{r_{root}}{r_{rhizo}}$ and Eq. (5.40a), the parameter A can be written

$$A = \left( \frac{c_{max}}{C_{crit}} \right)^{\frac{1}{b}}. \qquad (5.43)$$

Thus, we deduce the range action parameter a using Eqs. (5.43)–(5.40b),

$$\begin{aligned} a &= \ln \left[ \left( \frac{A}{C_{crit}} \right)^{\frac{1}{r_{rhizo}}} \right] \\ &= \frac{1}{r_{rhizo} - r_{root}} \ln \left[ \frac{c_{max}}{C_{crit}} \right] \end{aligned} \qquad (5.44)$$

Therefore, the distribution of EPS concentration in the soil is given by

$$c(r) = c_{max} \frac{r_{rhizo} - r}{r_{rhizo} - r_{root}} C_{crit} \frac{r - r_{root}}{r_{rhizo} - r_{root}}. \qquad (5.45)$$
Zickenrott et al. [376] estimate the concentration in the rhizosphere between $5 \times 10^{-3}$ and 5 gram of dry mucilage per 100 gram of dry soil. Thus, we consider $c_{max} = 5$ and $C_{crit} = 5 \times 10^{-3}$. Note that these values are difficult to defined experimentally and should change with respect to plants and soil. Note, these values are difficult to define experimentally and depend on the plant and soil.

To be consistent with our experimental data defined for drop infiltration, we propose to modify the EPS distribution by taking the minimum between $c(r)$ and $c = 1$, i.e. the new EPS distribution is given by $\tilde{c}$ such that,

$$ \tilde{c}(r) = \max\{c(r), 1\}. \qquad (5.46) $$

For the rest of the study, we will use the notation without the tilde to simplify the expression.

Carminati [65] then later Kroener et al. [200] propose to model rhizosphere problem with the following parameters :

$$ \begin{cases} r_{\text{root}} = 0.05\text{cm} = 5 \times 10^{-4}\text{m}, \\ r_{\text{rhizo}} = 0.2\text{cm} = 2 \times 10^{-3}\text{m}, \\ r_{\text{bulk}} = 1\text{cm} = 10^{-2}\text{m}. \end{cases} \qquad (5.47) $$

Note that the radius of the root and the rhizosphere can change from plant effect (growth, age, hydric stress) but we neglect these effect in this study. Therefore, we define the studied soil using the EPS distribution in the soil around the root as displayed in Figure 5.17.

<figure><img src="image_106.png" /><figcaption>Figure 5.17: EPS distribution in the soil around a root</figcaption></figure>

### 5.5.2 Model of water transfer around a single root

The water transfer around a root is generally considered as angular invariance around [61, 301], therefore we work on an axisymmetric system. Using the soil transfer model developed to represent the infiltration of a drop, Eq. (5.13), coupled with a root sampling $S$, we obtain the following system

$$ \begin{cases} \phi \partial_t s = \frac{1}{r} \partial_r \left[ r K(s) \partial_r \frac{\delta F}{\delta s} \right] + S, \\ \frac{\delta F}{\delta s} = - \rho g H_{\text{cap}} \Pi_w(s, c) - \Gamma \Delta s \end{cases} \qquad (5.48) $$
We use the same parameters and functions to represent energies and hydraulic conductivity as in section 5.1.3 but we add in the free energy an osmotic energy $f_{osm}$, i.e

$$f(s, c) = f_w(s, c) + f_{\nabla}(s) + f_{osm}(s, c). \qquad (5.49)$$

The osmotic energy represents the swelling property of the EPS, i.e it represents the ability to store more water in the soil in presence of EPS. This property represents the hydrophilic effects of the amphiphilic matter. We propose to model this effect as acting only at high saturation, since at low saturation we find the hydrophobic property of EPS. In addition, hydrophilic effect increases with respect to the EPS concentration. Therefore, the osmotic energy is a function such that

* the osmotic pressure $h_{osm}$ is an increasing function with respect to s and c;

* the osmotic effect are neglected for low saturation;

* there is no osmotic effect for a soil without EPS ($c = 0$).

In this study, we assume a theoretical expression of the osmotic pressure $h_{osm}$,

$$h_{osm}(s, c) = \frac{\delta f_{osm}}{\delta s}(s, c) = c(\exp(s) - 1). \qquad (5.50)$$

Note that this function require some experimental calibration model a specific situation. Moreover, since we neglect the degradation of the amphiphilic matter, the distribution of the EPS concentration remains unchanged over time and the repellent and osmotic property of EPS depend only the initial distribution of the EPS.

In the case of soil humectation near a root, the soil is model using specific boundary condition far the root/rhizosphere interface. We propose to impose a saturation and pressure on the bulk soil ($r = r_{bulk}$) corresponding to an equilibrium far from the root. Moreover, when the soil is very dry, we can consider that the root imposes no suction, i.e. we impose non-flow conditions on the root ($r = r_{root}$).

To solve PDE system Eq.(5.48), we employ a finite element approach on the weak formulation (see appendix I) using the C++ FEM library OOMPH-LIB [169]. For time-stepping, a backward differen- tiation formula (BDF) scheme of second order is used, while both the time-steps and the spatial mesh are adaptive to the dynamics.

## 5.6 Humectation of the soil

In this context of rhizosphere humectation, we propose to present two simulations: a first one without taking into account the amphiphilic matter [290, 183, 111] and a second one considering the spatial distribution of the amphiphilic matter.

In the first simulation, we start with a dry profile close to the root and with a high saturation of 80% at $r = r_{bulk}$ (Figure 5.18a). The imbibition front approaches the root over time: it arrives in the rhizosphere after 3.6 s (Figure 5.18b), then passes through an intermediate point (black dashed line in
Figure 5.18) after 6.1 s (Figure 5.18c) to finally reach the root after 7 s (Figure 5.18d). It can be seen
that the advance of the front is the same in bulk soil and in the rhizosphere.

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

(a) $t = 0$ s

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

(b) $t = 3.6$ s

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

(c) $t = 6.1$ s

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

(d) $t = 7$ s

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

(e) Equilibrium

<caption>Figure 5.18: Water transfer from the bulk soil to the root without EPS concentration.</caption>

In the second simulation, we start from the same initial state as in first simulation (Figure 5.19a).
The front advances globally in the same way as the first simulation. However, imbibition is slower: the
imbibition front reaches the root after 13.3 s (Figure 5.19d). This slower dynamic is already observed
when the rhizosphere is reached after 4.6 s, which is 1 s more than previously (Figure 5.19b). Moreover,
in the rhizosphere, the imbibition front took more time to cover the distance between the rhizosphere
at $r = r_{rhizo}$ and the intermediate point than to cross the entire bulk soil. Thus, we observe a delayed
dynamics of rhizosphere wetting in the presence of amphiphilic matter as mentioned in Carminati [65],
Kroener et al. [200] or Zarebanadkouki et al. [366]. However, this delayed dynamics in this simulation
is only of few seconds compare to realistic situation. This effect is due to the fact that we consider no modification of the hydraulic conductivity in the presence of EPS. In the case where the hydraulic conductivity decrease in presence of EPS as proposed by Kroener et al. [199] and Landl et al. [206] then our model is able to obtain a delayed dynamics of few hours as realistic situation. It should be noted that the delayed dynamics could also be a consequence of the dynamics of degradation, production and diffusion of EPS, which we neglect in this situation.

Moreover, the equilibrium profiles in the first simulation (Figure 5.19e) is uniform saturation equal

<figure><img src="image_112.png" /><figcaption>Figure 5.19: Water transfer from the bulk soil to the root by taking into account the EPS concentration.</figcaption></figure>

to soil saturation ($s = 0.8$) far from the root whereas in the second simulation (Figure 5.18e), we observe that the rhizosphere stores more water than the soil. The latter result is in agreement with the observation obtained by Kroener et al. [200] who showed that the rhizosphere retains more water
than the bulk soil at equilibrium. This phenomenon is due to the swelling of the EPS.

In order to improve the model, it should be relevant to improve the representation of the water
storage and the delayed dynamic between rhizosphere and bulk soil using a better calibration of the
osmotic and the hydraulic conductivity function.

5.7 Conclusion

Two water transfer models were proposed in this chapter: first the coupling between a water droplet and the soil, taking into account surface interactions, and then water transfer around a single root. These models are based on the same construction of water transfer in the rhizosphere in the presence of EPS. Using the same free energy construct, we were able to reproduce different behaviors found in the literature, such as the dynamics of infiltration and spreading of a drop on the soil surface, slow diffusion in the presence of EPS, capillary rise in stratified media and delayed wetting of the rhizosphere.

The free energy concept has allowed us not only to build a thermodynamically coherent model, but also to facilitate the formulation of process coupling. In particular, the coupled influence of wettability between a thin film of liquid on the surface and a porous medium could not have been taken into account by working only on the pressure terms. In this way, the notion of free energy is a key factor in the construction of transfer dynamics in soils in the presence of amphiphilic matter.

These models offer new ways of representing water transfer and coupling processes. Therefore it is
an encouraging view to better understanding water processes in presence of EPS.
Conclusion

The final objective of this work was to propose a soil transfer model taking into account the presence
of amphiphilic matter. This work was carried out in four stages:

1) The development of a water transfer model based on the free energy of water in the soil. This model provides new parameters to describe the macroscopic processes at work between water, porous matrix and organic matter.

2) The implementation of numerical tools to simulate the temporal evolution of the developed models. We also provide analysis tools to study water transfer behaviors and the detection of thresholds.

3) Characterization of hydraulic properties in the presence of amphiphilic matter in a rhizosphere analogue medium, to understand the role of this substance. We also parameterize the hydraulic functions related to this kind of medium.

4) The development of models able to represent water transfers in a rhizosphere analogue medium in the presence of amphiphilic matter, such as infiltration and humectation near the root.

We clarified first the ambiguity between capillary pressure and thermodynamic potential in soils: capillary pressure is the functional derivative of free energy density. Then we noted that capillary pressure is often associated with soil wettability by taking into account a contact angle Θ and the related capillary height. However, contact angle or capillary height do not allow to reproduce flows specific to the presence of organic matter such as preferential flows or water retention. In contrast, the free energy provides a general framework to derive a model of gradient dynamics as it is done in chapter 2. Moreover, to take into account all components acting on the water, the total free energy is the sum of the different components taken separately. But this reasoning cannot be applied to pressure. We define two new terms in the free energy: i) the wetting energy and ii) the free surface energy. i) The wetting energy represents the interaction energy between water and the soil matrix. Both attractive and repellent interactions are included. We derive a wetting pressure generalizing the capillary pressure. The integration of repellent term in free energy is not documented in the literature and it is a key term of the improvement of the model. ii) The surface energy was already introduced in the literature. In our context, it plays a regularizing role in modeling hydrophobic media.

The constructed wettability introduces a non-monotonic curve, reproducing different dynamics
during wetting/drying cycles. This allows us to reproduce hysteresis phenomenon during water transfer
in soil. Thus, we have provided an alternative to the relaxation time models for modeling hysteresis
and wetting phenomenon in the presence of organic matter [201].
Using this concept, we have developed a new model of water transfer in the presence of hydrophobic matter. This model is a highly non-linear fourth-order PDE that requires numerical methods to study the dynamics of water transfer. We began by developing our own tools in Matlab, using the proposed numerical schemes, such as taking into account explicitly the destabilizing term that occurs in the presence of hydrophobic terms. However, to reduce computation time and optimize the algorithm, we switched to a finite element library Oomph-lib. This library offers adaptive time-stepping and meshes for irregular dynamics as intermittent flow in 1D and 2D. In addition, a first version of continuation code has been developed to track water behavior for different situations and scenarios. These numerical developments associated with this free energy model highlighted several transfer dynamics: fingering, water trapping, intermittency.

Then, we carried out hydraulic experiments on rhizosphere analogue media to characterize it with respect to the concentration of amphiphilic matter. This characterization enabled us to extract a first range of parameters to model a rhizospheric soil. Further experiments are required to extract all the parameters for modeling, in particular the influence of amphiphilic matter concentration on the hydrophobicity of the medium.

Once we had calibrated our model on the basis of our experiments and literature data, we were able to obtain results applied to the rhizosphere, such as drop-soil exchanges, infiltration dynamics and humectation delay in the rhizosphere. These results are consistent with the literature, but require further experiments to be completely calibrated. We correlated the retention curve for different EPS concentration *c* with a reference medium (*c* = 0%). We adjusted the free energy of the hydrophobic part with our water retention experiments and experiment capillary rise in a stratified medium with EPS literature experiments [201]. However, we have not adjusted the hydrophilic part in relation to EPS swelling. Taking into account the calibration of this energy should better simulate the water storage capacity in the rhizosphere as defined in the literature [288, 280, 130, 302].

A first outlook in the continuity of this work is the calibration of the osmotic effect of EPS in contact with water by modifying the wetting energy density. An increasing of the hydrophilic term in the saturation *s* range near 1 could reproduce the osmotic effect near saturation. Moreover experiments could be realize to calibrate more accurately the law using to represent the hydrophobic effect of EPS concentration on the water transfer or define the EPS degradation (project MucilageMatter submit). It is also possible to define experiments to calibrate surface tension, such as measuring the typical size between two fingering. Moreover, the analogy between multiphase media and porous media can be discussed to define this tension. Indeed, a more complex expressions could be considered as proposes by Beljadid *et al* [17].

Another perspective of this work would be to provide more continuation tools to study behavior in several contexts than the current implementation (periodic boundary condition in 1D). An improvement to our continuation tools with the PDE2PATH library would be to develop a version with Dirichlet or Neumann boundary conditions in 1D and 2D. In this way, a new tool, not currently available in our laboratory, would be developed to study stationary saturation states in the bulk soil and rhizosphere around a root under different wetting conditions. Moreover, the library *trullekrul* developed under MATLAB would be used to adapt the mesh to our solution obtained with a finite element method.

Finally, this free energy approach appears to be a potential candidate for modeling coupling processes in the soil. This type of construction would make it possible to associate an energy with each
process to define a global energy for the system. For example, it could be possible to modeling the
dynamics of EPS degradation with soil transfer as biofilm model proposed by Trinschek *et al.* [336, 337]
and coupled it with the soil transfer model developed during this thesis. Finally, to couple all these
dynamics with the root uptake. Integrated into a broader framework, at the scale of a root sys-
tem [110, 299], this modelling will contribute to the analysis of rhizospheric traits likely to improve
water use by plants in a context of climate change and practices.

In another context, this model could be used to study the transport of substances to remediate soil pollution. For this, it would be interesting to introduce energies associated with diffusiophoresis, as proposed in the TRACE-it project: Controlling particle flow driven by local concentration gradients in geological porous media.

Therefore, this work opens up a new way of coupling the processes associated with soil water
transfer in the presence of amphiphilic matter.
# Appendices
<figure><img src="image_113.png" /></figure>
# Appendix A

## Numerical time scheme

The aim of this chapter is to present several time scheme and prediction-correction methods used in this manuscript. Let's define $\Omega$ a domain in $\mathbb{R}^N$ with $N$ the dimension of the domain and a time interval $[0, T]$. We consider a physical phenomenon represented by the function $y(\mathbf{x}, t)$ in the domain during the time interval, i.e the spatial coordinates are given by $\mathbf{x} \in \Omega$ and the time duration by $t \in [0, T]$. This phenomenon is express by the general formulation of differential equation,

$$ \dot{y}(\mathbf{x}, t) = \partial_t y(\mathbf{x}, t) = f(y(\mathbf{x}, t), t), \qquad (A.1) $$

where $y(\mathbf{x}, 0) = y_0(\mathbf{x})$ is the initial condition, $f$ is a smooth (linear or non-linear) function defined on $\Omega \times [0, T]$ and satisfying the Lipschitz condition.

To simplify the notations, the spatial part is neglected and we concentrate on solving the EDO from Eq. (A.1):

$$ \dot{y}(t) = f(y(t), t). \qquad (A.2) $$

### A.1 Time scheme

A numerical time scheme is used to numerically approximate the solution at each time $t$. These schemes are used to solve problems where the analytical solution is not known. The finite difference method is used to generate our numerical schemes. This method is based on Taylor series applied in a discrete domain.

**Definition A.1.1.** Let's a function $y$ is infinitely differentiable at a real or complex number $t_1$. The Taylor series of this function is a power series:

$$ y(t) = \sum_{k=0}^{\infty} \frac{y^{(k)}(t_1)}{k!} (t - t_1)^k, \qquad (A.3) $$

where $k!$ denotes the factorial of $k$, $y^{(k)}(t_1)$ denotes the $k$th derivative of $y$ evaluated at the point $t_1$.

A subdivision of time is assumed $0 = t_0 < t_1 < t_2 < \dots < t_m = T$ such that $[0, T] = \cup_{i=0}^{m-1} [t_i, t_{i+1}]$. This subdivision forms a temporal mesh in which each point constitutes a discretization node. The
discretization step between the mesh nodes on the interval $[t_i, t_{i+1}]$ is defined by $dt_n$ such that $t_{n+1} = dt_n + t_n$.

A numerical scheme allows the calculation of an approximation of the solution at time $t_{n+1}$, denoted $y_{n+1}$, and an approximation error is made between the exact formulation of $y(t_{n+1})$ and $y_{n+1}$. This error is called the local truncation error (LTE), i.e.

$$LTE = y_n - y(t_n). \qquad (A.4)$$

This error is considered of order $q$, if $y_n$ is of the form $\sum_{n, l \in [0, m], n \neq l}$.

$$y_n = \sum_{k=0}^{q-1} \frac{y_l^{(k)}}{k!} (t_n - t_l)^k. \qquad (A.5)$$

The LTE is therefore of the form

$$LTE = \sum_{k=q}^{\infty} \frac{y_l^{(k)}}{k!} (t_n - t_l)^k \approx \mathcal{O}((t_n - t_l)^{q+1}). \qquad (A.6)$$

In the rest of this chapter, we consider that the approximation is given for all $n \in [0, m]$,

$$y(t_n) \approx y_n. \qquad (A.7)$$

### A.1.1 Time scheme at one step

A one-step scheme is a time scheme that depends only on the previous instant. By applying the Taylor series given by Eq.(A.5), the exact behaviour of the physical phenomenon at time $t_{n+1}$ is naturally given as a function of time $t_n$ by the relation:

$$y(t_{n+1}) = \sum_{k=0}^{\infty} \frac{y^{(k)}(t_n)}{k!} dt_n^k, \qquad (A.8)$$

or

$$y(t_{n+1}) = y(t_n) - \sum_{k=1}^{\infty} \frac{y^{(k)}(t_{n+1})}{k!} (-dt_n)^k. \qquad (A.9)$$

The one-step schemes are linear combinations of its two exact expressions. Because we are trying to solve Eq. A.1, the schemes that we need to use are limited to the function $y$ and its first derivative. Let $\theta \in [0, 1]$, we define the family of $\theta$-schemes:

$$y_{n+1} = y_n + dt_n [(\theta - 1) \dot{y}_n + \theta \dot{y}_{n+1}]. \qquad (A.10)$$

There are two main categories of scheme:

* explicit schemes that depend only on known instants, in this case $t_n$;

* implicit schemes depending on the future moment $t_{n+1}$.
In the family of $\theta$ schemes, this differentiation is made according to the value of $\theta$. If $\theta = 0$, then the scheme is explicit, otherwise the scheme is implicit. Among these schemes, we find the classically used Euler schemes:

* Forward Euler scheme (FE), which is an explicite scheme

$$y_{n+1} = y_n + dt_n \dot{y}_n. \qquad (A.11)$$

* Backward Euler scheme (BE), which is an implicite scheme

$$y_{n+1} = y_n + dt_n \dot{y}_{n+1}. \qquad (A.12)$$

The particularity of these schemes is their simplicity to develop, but their order of accuracy is low. They are first-order methods, with a LTE,

$$LTE^{FE} = -\frac{dt_n^2}{2} \ddot{y}_n + O(dt_n^3) \qquad (A.13)$$

and

$$LTE^{BE} = \frac{dt_n^2}{2} \ddot{y}_n + O(dt_n^3). \qquad (A.14)$$

Generally, $\theta$ schemes are simple to use but have a low order of accuracy, except for the Crank-Nicolson scheme, also known as the Trapezoid rule (TR) [180]. This scheme is of order 2, when the temporal mesh is uniform, so $\theta = \frac{1}{2}$,

$$y_{n+1} = y_n + \frac{dt_n}{2} [\dot{y}_n + \dot{y}_{n+1}], \qquad (A.15)$$

with a LTE,

$$LTE^{TR} = \frac{dt_n^3}{12} \dddot{y}_n + O(dt_n^4) \qquad (A.16)$$

In order to have wider range of schemes with better accuracy, it is a smart strategy to turn to multi-step schemes.

### A.1.2 Multi-step time scheme

Multi-step schemes are schemes that depend on the previous time instants $t_{n-p}$, with $n \ge p \ge 0$. Let's introduce $K_1 \ge 1$, $K_2 \ge 0$ and $K = \max(K_1, K_2)$. We consider two real sequences $\{\alpha_j\}_{j \in I_1}$ and $\{\beta_j\}_{j \in I_2}$ where $I_1 = \{1, \dots, K_1\}$ and $I_2 = \{0, \dots, K_2\}$. The linear multi-step method is a K-step method given for all $n \ge K$ [45, 51, 305, 154]:

$$y_n = \sum_{j=1}^{K_1} \alpha_j y_{n-j} + dt \sum_{j=0}^{K_2} \beta_j \dot{y}_{n-j}. \qquad (A.17)$$

There exists two main classes of methods: Adams methods and Backward differential formula (BDF).
## Adams methods

schemes based on the Adams method are divided into two categories: Adams-Bashforth schemes and Adams-Moulton schemes.

Adams-Bashforth schemes are explicit schemes of order q written as

$$y_n = y_{n-j} + dt \sum_{j=1}^{q} \beta_j \dot{y}_{n-j}. \qquad (A.18)$$

The form of this class of scheme is in the form of a multi-step linear scheme with $K_1 = 1$ (and $\alpha_1 = 1$), $K_2 = q$ and $\beta_0 = 0$. The FE scheme is a first-order scheme from this class, denoted AB1. schemes of higher order than 1 are present, such as AB2 [305], which is a scheme of order 2, defined by:

$$\frac{y_{n+1} - y_n}{dt_n} = \left(1 + \frac{dt_n}{2dt_{n-1}}\right) \dot{y}_n - \frac{1}{2dt_{n-1}} \ddot{y}_{n-1}, \qquad (A.19)$$

with a LTE,

$$LTE^{AB2} = \frac{dt_n^2}{12} (2dt_n + 3dt_{n-1}) \dddot{y}_n + \mathcal{O}(dt_n^4). \qquad (A.20)$$

Adams-Moulton schemes are implicit schemes of order q written as

$$y_n = y_{n-j} + dt \sum_{j=0}^{q-1} \beta_j \dot{y}_{n-j}. \qquad (A.21)$$

This form corresponds to a multi-step linear scheme with $K_1 = 1$ (and $\alpha_1 = 1$), $K_2 = q-1$ and $\beta_0 > 0$. This class of scheme includes the BE and TR schemes, called AM1 and AM2 respectively.

## BDF methods

The BDF schemes [146, 172, 173] are implicit schemes of order q written

$$y_n = \sum_{j=1}^{q} \alpha_j y_{n-j}. \qquad (A.22)$$

This form corresponds to a multi-step linear scheme with $K_1 = q$, $K_2 = 0$ and $\beta_0 > 0$. The BE scheme is a first-order scheme in this class, denoted BDF1. Note that the BE scheme is a member not only of this class but also of the Adams-Moulton class. Higher-order schemes are present in this class, such as the BDF2 scheme [70], which is a second-order scheme,

$$\frac{y_{n+1} - y_n}{dt_n} = \frac{dt_n}{2dt_n + dt_{n-1}} \frac{y_n - y_{n-1}}{dt_{n-1}} + \frac{dt_n + dt_{n-1}}{2dt_n + dt_{n-1}} \dot{y}_{n+1}. \qquad (A.23)$$

## Other scheme

Linear multi-step schemes are not limited to these two main classes of schemes. In these two classes, only Adams-Bashforth schemes are explicit schemes. We construct a scheme with $K_1 > 1$, we have schemes that are not in this class. For example, with $K_1 = 2$ and $K_2 = 1$ with $\beta_0 = 0$, we obtain the
generalised Leapfrog scheme,

$$y_{n+1} = \left(1 - \left(\frac{dt_n}{dt_{n-1}}\right)^2\right) y_n + \left(\frac{dt_n}{dt_{n-1}}\right)^2 y_{n-1} + dt_n \left(1 + \frac{dt_n}{dt_{n-1}}\right) \dot{y}_n, \quad (A.24)$$

with a LTE,

$$LTE^{Leapfrog} = \frac{dt_n^3}{6} \left(1 - \frac{dt_{n-1}}{dt_n}\right) \ddot{y}_n + \mathcal{O}(dt_n^4). \quad (A.25)$$

This scheme is an explicit scheme of order 2.

## A.2 Prediction-correction methods

Prediction-correction methods are iterative methods used to define the future time step once $y_{n+1}$ is obtained. In order to better follow the physical dynamics represented by Eq. (A.1), these methods are based on the prediction and correction of the time step $dt_{n+1}$.

### A.2.1 General principle

The principle of these methods consists in the prediction of a time step using the estimation of the expression of an implicit scheme's LTE depending on the expression of an explicite scheme of the same order. The prediction is considered correct if the matrix system associated with Eq. (A.1) converges at the next iteration, otherwise the time step must be corrected by reducing it. To prevent multiple corrections, it is important to predict the future time step. We explain below how the prediction principle works.

An explicit scheme is used to predict the solution of Eq. (A.1), this solution is denoted $y_{n+1}^P$. This solution is used to solve the matrix system associated with the problem based on an implicit scheme and define the solution $y_{n+1}$. However, an approximation error is generated using this implicit scheme. In order to reduce the error at the next iteration, the LTE of the implicit scheme is expressed as a function of the predicted solution $y_{n+1}^P$ and the solution $y_{n+1}$, at time $t_{n+1}$. A constraint is defined in order to establish our estimate and remain within a finite and bounded error formulation,

$$\|d_n\| \le \epsilon y_{max}, \quad (A.26)$$

where $|| \cdot ||$ is a norm, $\epsilon$ is a tolerance and $y_{max}$ is the maximum of $y_{n+1}$, $\forall n \ge 0$.

In addition to this constraint, the ratio of the norms of the LTEs over two iterations is expressed as a function of the previous discretisation steps by a function $h$,

$$\frac{\|d_n\|}{\|d_{n-1}\|} = h(dt_{n+1}, dt_n, \dots, dt_1). \quad (A.27)$$

Using these two expressions, we determine the relationship

$$\|d_n\| = h(dt_{n+1}, dt_n, \dots, dt_1) \|d_{n-1}\| \le \epsilon y_{max}.$$

From this last relationship, we estimate the future discretisation step $dt_{n+1}$.
## A.2.2 Order 1

One of the first-order prediction-correction methods is based on the FE (Eq. (A.11)) and BE (Eq. (A.12)) schemes. The solutions are defined by

$$ 
\begin{aligned}
y_{n+1}^{P} &= y_n + dt_n \dot{y}_n, \\
y_{n+1} &= y_n + dt_n \dot{y}_{n+1}.
\end{aligned}
$$

From Eq. (A.14), the LTE associated with the BE scheme is

$$ d_n = \frac{dt_n^2}{2} \ddot{y}_n + \mathcal{O}(dt_n^3). $$

Using the expression for the predicted solution $y_{n+1}^P$ and its LTE (Eq. (A.13)), we determine an expression for $d_n$ as a function of the two solutions

$$ d_n = \frac{y_{n+1} - y_{n+1}^P}{2}. $$

A recurrence relation is defined between the LTEs from previous iterations,

$$ \frac{d_n}{d_{n-1}} = \left( \frac{dt_n}{dt_{n-1}} \right)^2 + \mathcal{O}(dt). \qquad (A.28) $$

Combining the constraint given by Eq. (A.26) and the recurrence relation, Eq. (A.28), the LTE is bounded according to the previous time steps,

$$ ||d_n|| = \left( \frac{dt_n}{dt_{n-1}} \right)^2 ||d_{n-1}|| \le \epsilon y_{max}. $$

In other words, the actual time step is bounded as a function of the previous steps,

$$ dt_n \le \left( \frac{\epsilon y_{max}}{||d_n||} \right)^{\frac{1}{2}} dt_{n-1}, $$

where the norm of $d_n$ is given by

$$ ||d_n|| = \frac{||y_{n+1} - y_{n+1}^P||}{2}. $$

Therefore the future time step is estimated as:

$$ dt_{n+1} = \left( \frac{\epsilon y_{max}}{||d_n||} \right)^{\frac{1}{2}} dt_n. \qquad (A.29) $$

## A.2.3 Order 2

Similar approach than the previous section A.2.2 allow to define the future time step.

With $dt$ is the actual time step and $dt_p$ the previous time step, we select the Backward differential
formulation (BDF2) scheme as the time scheme

$$ \dot{y}_{n+1} = \left( \frac{1}{dt} + \frac{1}{dt + dt_p} \right) y_{n+1} - \frac{dt + dt_p}{dt dt_p} y_n + \frac{dt}{dt_p (dt + dt_p)} y_{n-1}, $$

and the Generalized Leapfrog scheme to predict the solution $y_{n+1}^P$,

$$ y_{n+1}^P = \left(1 - \frac{dt^2}{dt_p^2}\right) y_n + \frac{dt^2}{dt_p^2} y_{n-1} + \left(1 + \frac{dt}{dt_p}\right) dt \dot{y}_n. $$

The local truncation error (LTE) of the BDF2 scheme is

$$ d_n = \frac{\left(1 + \frac{dt_p}{dt}\right)^2}{1 + 3\frac{dt_p}{dt} + 4\left(\frac{dt_p}{dt}\right)^2 + 2\left(\frac{dt_p}{dt}\right)^3} (y_{n+1} - y_{n+1}^P) + O(dt^4). $$

Soderlind and Wang [317] mention that the potential new time step is given by

$$ dt_n = dt \left( \frac{\epsilon}{\|d_n\|} \right)^{\frac{1}{3}}. $$
# Appendix B

## Matrix time system

The aim of this chapter is to formulate the matrix forms $A \in \mathcal{M}_2(\mathbb{R})$, $U \in \mathbb{R}^2$ and $b \in \mathbb{R}^2$ associated with the following problem,

$$ egin{cases} 
\phi \partial_t s = \nabla \cdot \left[ K(s) \nabla \frac{\delta \mathcal{F}}{\delta s} \right], \\ 
\frac{\delta \mathcal{F}}{\delta s} = -\Pi_w(s) - \Gamma \Delta - \rho g z. 
\end{cases} 
\quad (B.1) $$

Eq.(B.1) is expressed in weak form as follows

$$ egin{cases} 
\int_{\Omega} \partial_t sv = - \int_{\Omega} K(s) \nabla \frac{\delta \mathcal{F}}{\delta s} \nabla v \\ 
\int_{\Omega} \frac{\delta \mathcal{F}}{\delta s} w = - \int_{\Omega} \Pi(s) w + \int_{\Omega} T \nabla s \nabla w - \int_{\Omega} \rho g z w 
\end{cases} 
\quad (B.2) $$

This weak formulation, Eq.(B.2), can be expressed in matrix form $AU = b$. To facilitate matrix construction, let introduce $F_1(s, \frac{\delta \mathcal{F}}{\delta s}) = K(s) \nabla \frac{\delta \mathcal{F}}{\delta s}$ and $F_2(s, \frac{\delta \mathcal{F}}{\delta s}) = \Pi(s)$. Since $K(s)$ and $\Pi(s)$ are non-linear functions, we need to linearize them during the construction of the time scheme in order to obtain the matrix form associated with Eq.(B.2).

## B.1 Linear operator

In this section, we assume that the indices $i$ and $j$ represent the indices of discrete variables, and $\phi$ are the test functions resulting from the weak formulation. We define the usual linear operators to construct the matrix system.

The matrix $M$ is called the mass matrix and its elements are

$$ M_{ij} = \int_{\Omega} \phi_i \phi_j, \quad (B.3) $$

The matrix $K$ is called the stiffness matrix and its elements are

$$ K_{ij} = \int_{\Omega} \nabla \phi_i \nabla \phi_j \quad (B.4) $$
The stiffness matrix is generalized in the presence of a function $f$ by the matrix operator $L^f$ whose matrix elements are

$$L_{ij}^f = \int_{\Omega} f \nabla \phi_i \nabla \phi_j \qquad (B.5)$$

This last operator is useful for formulating linearized expressions.

## B.2 Backward Euler scheme

Using the Backward Euler scheme defined in Eq. (A.12), the system Eq. B.2 is written

$$\left\{
\begin{aligned}
\int_{\Omega} s_{n+1} v &= \int_{\Omega} s_n v - \int_{\Omega} dt_n F_1 \left(s_{n+1}, \frac{\delta \mathcal{F}_{n+1}}{\delta s}\right) \nabla v \\
\int_{\Omega} \frac{\delta \mathcal{F}_{n+1}}{\delta s} w &= - \int_{\Omega} F_2 \left(s_{n+1}, \frac{\delta \mathcal{F}_{n+1}}{\delta s}\right) w + \int_{\Omega} T \nabla s_{n+1} \nabla w - \int_{\Omega} \rho gzw
\end{aligned}
\right.
\qquad (B.6)$$

In linearized form, this expression becomes

$$\left\{
\begin{aligned}
\int_{\Omega} s_{n+1} v + \int_{\Omega} dt_n K(s_n) \nabla \frac{\delta \mathcal{F}_{n+1}}{\delta s} \nabla v &= \int_{\Omega} s_n v \\
\int_{\Omega} \frac{\delta \mathcal{F}_{n+1}}{\delta s} w - \int_{\Omega} T \nabla s_{n+1} \nabla w &= - \int_{\Omega} \Pi(s_n) w - \int_{\Omega} \rho gzw
\end{aligned}
\right.
\qquad (B.7)$$

This system is of the form $AU = b$ with the following matrix and vector elements:

$$A = \begin{bmatrix} M & dt_n L^{K(s_n)} \\ -L^T & M \end{bmatrix}, \qquad (B.8)$$

$$U = \begin{bmatrix} s_{n+1} \\ \frac{\delta \mathcal{F}_{n+1}}{\delta s} \end{bmatrix} \qquad (B.9)$$

and

$$b = \begin{bmatrix} M s_n \\ -M (\Pi(s_n) + \rho g z) \end{bmatrix}. \qquad (B.10)$$

## B.3 Trapezoid rule scheme

Using the Trapezoid rule scheme defined in Eq. (A.15), the system Eq. B.2 is written

$$\left\{
\begin{aligned}
\int_{\Omega} s_{n+1} v &= \int_{\Omega} s_n v - \int_{\Omega} \frac{dt_n}{2} F_1 \left(s_n, \frac{\delta \mathcal{F}_n}{\delta s}\right) \nabla v - \int_{\Omega} \frac{dt_n}{2} F_1 \left(s_{n+1}, \frac{\delta \mathcal{F}_{n+1}}{\delta s}\right) \nabla v \\
\int_{\Omega} \frac{\delta \mathcal{F}_{n+1}}{\delta s} w &= - \int_{\Omega} F_2 \left(s_{n+1}, \frac{\delta \mathcal{F}_{n+1}}{\delta s}\right) w + \int_{\Omega} T \nabla s_{n+1} \nabla w - \int_{\Omega} \rho gzw
\end{aligned}
\right.
\qquad (B.11)$$

In linearized form, this expression becomes
$$ \left\{ \begin{aligned} & \int_{\Omega} s_{n+1} v = \int_{\Omega} s_{n} v - \int_{\Omega} dt_{n} K(s_{n}) \nabla \frac{\delta \mathcal{F}_{n+1}}{\delta s} \nabla v \\ & \int_{\Omega} \frac{\delta \mathcal{F}_{n+1}}{\delta s} w - \int_{\Omega} T \nabla s_{n+1} \nabla w = - \int_{\Omega} \Pi(s_{n}) w - \int_{\Omega} \rho gzw \end{aligned} \right. \qquad (B.12) $$

This system is of the form $AU = b$ with the following matrix and vector elements:

$$ A = \begin{bmatrix} M & \frac{dt_n}{2} L^{K(s_n)} \\ -L^T & M \end{bmatrix}, \qquad (B.13) $$

$$ U = \begin{bmatrix} s_{n+1} \\ \frac{\delta \mathcal{F}_{n+1}}{\delta s} \end{bmatrix} \qquad (B.14) $$

and

$$ b = \begin{bmatrix} M s_n - \frac{dt_n}{2} L^{K(s_n)} s_n \\ -M \left( \Pi(s_n) + \rho g z \right) \end{bmatrix}. \qquad (B.15) $$

## B.4 BDF2 scheme

Using the BDF2 scheme defined in Eq. (A.23), the system Eq. B.2 is written

$$ \left\{ \begin{aligned} & \int_{\Omega} s_{n+1} v = \frac{(dt_n + dt_{n-1})^2}{2dt_n + dt_{n-1}} \int_{\Omega} s_n v - \frac{dt_n^2}{dt_{n-1}(2dt_n + dt_{n-1})} \int_{\Omega} s_{n-1} v \\ & \qquad - \frac{dt_n (dt_n + dt_{n-1})}{2dt_n + dt_{n-1}} \int_{\Omega} F_1 \left( s_{n+1}, \frac{\delta \mathcal{F}_{n+1}}{\delta s} \right) \nabla v \\ & \int_{\Omega} \frac{\delta \mathcal{F}_{n+1}}{\delta s} w = - \int_{\Omega} F_2 \left( s_{n+1}, \frac{\delta \mathcal{F}_{n+1}}{\delta s} \right) w + \int_{\Omega} T \nabla s_{n+1} \nabla w - \int_{\Omega} \rho gzw \end{aligned} \right. \qquad (B.16) $$

In linearized form, this expression becomes

$$ \left\{ \begin{aligned} & \int_{\Omega} s_{n+1} v + \frac{dt_n (dt_n + dt_{n-1})}{2dt_n + dt_{n-1}} \int_{\Omega} K(s_n) \nabla \frac{\delta \mathcal{F}_{n+1}}{\delta s} \nabla v = \frac{(dt_n + dt_{n-1})^2}{2dt_n + dt_{n-1}} \int_{\Omega} s_n v \\ & \qquad - \frac{dt_n^2}{dt_{n-1}(2dt_n + dt_{n-1})} \int_{\Omega} s_{n-1} v \qquad (B.17) \\ & \int_{\Omega} \frac{\delta \mathcal{F}_{n+1}}{\delta s} w - \int_{\Omega} T \nabla s_{n+1} \nabla w = - \int_{\Omega} \Pi(s_n) w - \int_{\Omega} \rho gzw \end{aligned} \right. $$

This system is of the form $AU = b$ with the following matrix and vector elements:

$$ A = \begin{bmatrix} M & \frac{dt_n(dt_n + dt_{n-1})}{2dt_n + dt_{n-1}} L^{K(s_n)} \\ -L^T & M \end{bmatrix}, \qquad (B.18) $$

$$ U = \begin{bmatrix} s_{n+1} \\ \frac{\delta \mathcal{F}_{n+1}}{\delta s} \end{bmatrix} \qquad (B.19) $$
and

$$b = \left[ \begin{array}{l} \frac{(dt_n + dt_{n-1})^2}{2dt_n + dt_{n-1}} Ms_n - \frac{dt_n^2}{dt_{n-1}(2dt_n + dt_{n-1})} Ms_{n-1} \\ -M(\Pi(s_n) + \rho gz) \end{array} \right]. \qquad (B.20)$$
# Appendix C

## Variational formulation of the problem

The aim of this chapter is to define the weak formulations associated to each boundary condition that can be studied.

**Remark C.0.1.** To simplify notation, we will write

* $\int_{\Omega} v = \int_{\Omega} v(x)dx$ and $\int_{\partial\Omega} v = \int_{\partial\Omega} v(x)d\lambda(x);$

* $\delta_s \mathcal{F}$ for $\frac{\delta \mathcal{F}}{\delta s};$

* $K(s)\nabla\delta_s \mathcal{F}$ for $K(s)\nabla\frac{\delta \mathcal{F}}{\delta s}.$

Therefore, the system (B.1) is equivalent to:

$$ \partial_t s = \nabla \cdot (K(s)\nabla\delta_s \mathcal{F}) + \chi(s), \qquad (C.1a) $$

$$ \delta_s \mathcal{F} = -\Pi(s) - T\Delta s. \qquad (C.1b) $$

We multiply Eq.(C.1a) by a test function $v \in H^1(\Omega)$, then integrate over $\Omega$:

$$ \int_{\Omega} \partial_t s \cdot v = \int_{\Omega} (\nabla \cdot (K(s)\nabla\delta_s \mathcal{F})) \cdot v + \chi(s) \cdot v. \qquad (C.2) $$

Using Green formula, we obtain the following expression:

$$ \int_{\Omega} \partial_t s \cdot v = \int_{\Omega} \chi(s) \cdot v - \int_{\Omega} K(s) \nabla \delta_s \mathcal{F} \cdot \nabla v + \int_{\partial\Omega} (K(s) \partial_n(\delta_s \mathcal{F}) \cdot v). \qquad (C.3) $$

Similarly, we obtain for $w \in H^1(\Omega)$:

$$ \int_{\Omega} \delta_s \mathcal{F} \omega = \int_{\Omega} T \nabla s \cdot \nabla \omega - \int_{\Omega} \Pi(s) \omega - \int_{\partial\Omega} T \partial_n(s) \omega. \qquad (C.4) $$

The equations (C.3) and (C.4) are the weak formulation (also called variational formulation) of the
(C.1) problem:

$$ \left\{ \begin{aligned} \int_{\Omega} \partial_t s \cdot v &= \int_{\Omega} \chi(s) \cdot v - \int_{\Omega} K(s) \nabla \delta_s \mathcal{F} \cdot \nabla v + \int_{\partial \Omega} K(s) \partial_n (\delta_s \mathcal{F}) v, \\ \int_{\Omega} \delta_s \mathcal{F} \omega &= - \int_{\Omega} \Pi(s) \omega + \int_{\Omega} T \nabla s \cdot \nabla \omega - \int_{\partial \Omega} T \partial_n(s) \omega. \end{aligned} \right. \qquad (C.5) $$

To complete this formulation, we need to add boundary conditions, which will be discussed in the next section.

**Remark C.0.2.** The residuals $r_0$ and $r_1$ will be defined as follows:

$$ \left\{ \begin{aligned} r_0 &= - \int_{\Omega} \partial_t s \cdot v + \int_{\Omega} \chi(s) \cdot v - \int_{\Omega} K(s) \nabla \delta_s \mathcal{F} \cdot \nabla v + \int_{\partial \Omega} K(s) \partial_n (\delta_s \mathcal{F}) v, \\ r_1 &= - \int_{\Omega} \delta_s \mathcal{F} \omega - \int_{\Omega} \Pi(s) \omega + \int_{\Omega} T \nabla s \cdot \nabla \omega - \int_{\partial \Omega} T \partial_n(s) \omega. \end{aligned} \right. \qquad (C.6) $$

## C.1 Boundary conditions

We use this notation to facilitate reading: D for Dirichlet condition, N for Neumann condition, NH for homogeneous Neumann condition, P for periodic condition and F for flux condition.

### C.1.1 General case

$$ \left\{ \begin{aligned} \partial_t s &= \nabla \cdot (K(s) \nabla \delta_s \mathcal{F}) + \chi(s), \\ \delta_s \mathcal{F} &= -T \Delta s - \Pi(s). \end{aligned} \right. $$

Using Green formula, we define the weak formulation on the domain $\Omega$:

$$ \left\{ \begin{aligned} \int_{\Omega} \partial_t s v &= \int_{\partial \Omega} K(s) \partial_n \delta_s \mathcal{F} v d\sigma - \int_{\Omega} K(s) \nabla \delta_s \mathcal{F} \cdot \nabla v + \int_{\partial \Omega} \chi(s) v, \\ \int_{\Omega} \delta_s \mathcal{F} w &= - \int_{\partial \Omega} T \partial_n s w d\sigma + \int_{\Omega} T \nabla s \cdot \nabla w - \int_{\partial \Omega} \Pi(s) w. \end{aligned} \right. $$

Let $\Gamma_D$, $\Gamma_N$, $\Gamma_F$ and $\Gamma_P$, the spaces on the edge of the domain with Dirichlet, Neumann, Flux and Periodic conditions on the edge, respectively. Thus, we define the boundaries of the domain $\partial \Omega$ such that

$$ \partial \Omega = \Gamma_D \cup \Gamma_N \cup \Gamma_F \cup \Gamma_P. $$
$$
\left\{
\begin{aligned}
\int_{\Omega} \partial_t sv &= - \int_{\Omega} K(s) \nabla \delta_s \mathcal{F} \nabla v + \int_{\Omega} \chi(s) v + \int_{\Gamma_D} K(s) \partial_n \delta_s \mathcal{F} v d\sigma + \int_{\Gamma_N} K(s) \partial_n \delta_s \mathcal{F} v d\sigma \\
&\quad + \int_{\Gamma_F} K(s) \partial_n \delta_s \mathcal{F} v d\sigma + \int_{\Gamma_P} K(s) \partial_n \delta_s \mathcal{F} v d\sigma \\
\int_{\Omega} \delta_s \mathcal{F} w &= \int_{\Omega} T \nabla s \cdot \nabla w - \int_{\Omega} \Pi(s) w, - \int_{\Gamma_D} T \partial_n sw d\sigma - \int_{\Gamma_N} T \partial_n sw d\sigma \\
&\quad - \int_{\Gamma_F} T \partial_n sw d\sigma - \int_{\Gamma_P} T \partial_n sw d\sigma.
\end{aligned}
\right.
$$

C.1.2 In dimension 1

Let $\partial\Omega = \Gamma_{in} \cup \Gamma_{out}$ with $\Gamma_{in}$ the input edge and $\Gamma_{out}$ the output edge. Let $a, b$ reals such as $a \in \Gamma_{in}$ and $b \in \Gamma_{out}$. Here are the different boundary conditions that can be encountered in one dimension:

<table><thead><tr><th colspan="2">Saturation s</th><th colspan="2">δ<sub>s</sub>F</th><th colspan="2">K(s)∇δ<sub>s</sub>F</th></tr><tr><th>Input</th><th>Output</th><th>Input</th><th>Output</th><th>Input</th><th>Output</th></tr></thead><tbody><tr><td>D</td><td>D</td><td>NH</td><td>NH</td><td>×</td><td>×</td></tr><tr><td>NH</td><td>NH</td><td>NH</td><td>NH</td><td>×</td><td>×</td></tr><tr><td>D</td><td>NH</td><td>NH</td><td>NH</td><td>×</td><td>×</td></tr><tr><td>NH</td><td>D</td><td>NH</td><td>NH</td><td>×</td><td>×</td></tr><tr><td>P</td><td>P</td><td>P</td><td>P</td><td>×</td><td>×</td></tr><tr><td>D</td><td>D</td><td>×</td><td>×</td><td>Flux</td><td>Flux</td></tr><tr><td>NH</td><td>NH</td><td>×</td><td>×</td><td>Flux</td><td>Flux</td></tr><tr><td>D</td><td>NH</td><td>×</td><td>×</td><td>Flux</td><td>Flux</td></tr><tr><td>NH</td><td>D</td><td>×</td><td>×</td><td>Flux</td><td>Flux</td></tr><tr><td>D</td><td>D</td><td>NH</td><td>×</td><td>×</td><td>Flux</td></tr><tr><td>NH</td><td>NH</td><td>NH</td><td>×</td><td>×</td><td>Flux</td></tr><tr><td>D</td><td>NH</td><td>NH</td><td>×</td><td>×</td><td>Flux</td></tr><tr><td>NH</td><td>D</td><td>NH</td><td>×</td><td>×</td><td>Flux</td></tr><tr><td>D</td><td>D</td><td>×</td><td>NH</td><td>Flux</td><td>×</td></tr><tr><td>NH</td><td>NH</td><td>×</td><td>NH</td><td>Flux</td><td>×</td></tr><tr><td>D</td><td>NH</td><td>×</td><td>NH</td><td>Flux</td><td>×</td></tr><tr><td>NH</td><td>D</td><td>×</td><td>NH</td><td>Flux</td><td>×</td></tr></tbody></table>

Therefore, we will try to find the solution space of the problem for the different cases detailed in the table.

Dirichlet input/output conditions for $s$, with NH input/output for free energy

An input/output Dirichlet condition is imposed for $s$ and $\partial_n \delta_s \mathcal{F} = 0$ on $\Gamma_{in} \cup \Gamma_{out}$:
$$ \left\{ \begin{aligned} & \int_{\Omega} \partial_t s.v = \int_{\Omega} K(s) \nabla \delta_s \mathcal{F} \cdot \nabla v + \int_{\Omega} \chi(s).v, \\ & \int_{\Omega} \delta_s \mathcal{F} \omega = - \int_{\Omega} \Pi(s) \omega - \int_{\partial \Omega} T \cdot \partial_n(s) \omega + \int_{\Omega} T \nabla s \cdot \nabla \omega. \end{aligned} \right. $$

This is equivalent to solve the following problem, with $V_1, V_2$ defined as follows:

$$ \begin{aligned} V_1 &= \{s \in H^1(\Omega) \text{ such as } s(a) = s_{in} \text{ and } s(b) = s_{out}\} \text{ and} \\ V_2 &= \{\delta_s \mathcal{F} \in H^1(\Omega) \text{ such as } \partial_n(\delta_s \mathcal{F})(a) = \partial_n(\delta_s \mathcal{F})(b) = 0\} \end{aligned} $$

$$ \left\{ \begin{aligned} & \text{Find } (s, \delta_s \mathcal{F}) \text{ in } V_1 \times V_2 \text{ such as:} \\ & \int_{\Omega} \partial_t s.v = \int_{\Omega} K(s) \nabla \delta_s \mathcal{F} \cdot \nabla v + \int_{\Omega} \chi(s).v, \quad \text{for all } v \in H_0^1(\Omega) \\ & \int_{\Omega} \delta_s \mathcal{F} \omega = - \int_{\Omega} \Pi(s) \omega - \int_{\partial \Omega} T \cdot \partial_n(s) \omega + \int_{\Omega} T \nabla s \cdot \nabla \omega, \quad \text{for all } \omega \in H^1(\Omega) \end{aligned} \right. $$

## NH input/output conditions for $s$, with NH input/output for free energy

We impose a homogeneous Neumann input/output condition for $s$ and $\partial_n \delta_s \mathcal{F} = 0$ on $\Gamma_{in} \cup \Gamma_{out}$:

$$ \left\{ \begin{aligned} & \int_{\Omega} \partial_t s.v = \int_{\Omega} K(s) \nabla \delta_s \mathcal{F} \cdot \nabla v + \int_{\Omega} \chi(s).v, \\ & \int_{\Omega} \delta_s \mathcal{F} \omega = - \int_{\Omega} \Pi(s) \omega + \int_{\Omega} T \nabla s \cdot \nabla \omega. \end{aligned} \right. $$

This is equivalent to solve the following problem, with $V_1, V_2$ defined as follows:

$$ \begin{aligned} V_1 &= \{s \in H^1(\Omega) \text{ such as } \partial_n s(a) = \partial_n s(b) = 0\} \text{ and} \\ V_2 &= \{\delta_s \mathcal{F} \in H^1(\Omega) \text{ such as } \partial_n(\delta_s \mathcal{F})(a) = \partial_n(\delta_s \mathcal{F})(b) = 0\}. \end{aligned} $$

$$ \left\{ \begin{aligned} & \text{Find } (s, \delta_s \mathcal{F}) \text{ in } V_1 \times V_2 \text{ such as:} \\ & \int_{\Omega} \partial_t s.v = \int_{\Omega} K(s) \nabla \delta_s \mathcal{F} \cdot \nabla v + \int_{\Omega} \chi(s).v, \quad \text{for all } v \in H^1(\Omega), \\ & \int_{\Omega} \delta_s \mathcal{F} \omega = - \int_{\Omega} \Pi(s) \omega + \int_{\Omega} T \nabla s \cdot \nabla \omega, \quad \text{for all } \omega \in H^1(\Omega). \end{aligned} \right. $$

## Dirichlet/Neumann conditions Homogeneous for $s$, with NH input/output for free energy

Consider first the case where we impose a Dirichlet condition on the input for $s$, a Homogeneous Neumann condition ($\partial_n s = 0$) on the output for $s$, and $\partial_n \delta_s \mathcal{F} = 0$ on $\Gamma_{in} \cup \Gamma_{out}$:

$$ \left\{ \begin{aligned} & \int_{\Omega} \partial_t s.v = \int_{\Omega} K(s) \nabla \delta_s \mathcal{F} \cdot \nabla v + \int_{\Omega} \chi(s).v, \\ & \int_{\Omega} \delta_s \mathcal{F} \omega = - \int_{\Omega} \Pi(s) \omega - \int_{\Gamma_{in}} T \cdot \partial_n(s) \omega + \int_{\Omega} T \nabla s \cdot \nabla \omega. \end{aligned} \right. $$

This is equivalent to solve the following problem, with $V_1, V_2$ defined as follows:

$$ \begin{aligned} V_1 &= \{s \in H^1(\Omega) \text{ such as } s(a) = s_{in} \text{ and } \partial_n s(b) = 0\} \text{ and} \\ V_2 &= \{\delta_s \mathcal{F} \in H^1(\Omega) \text{ such as } \partial_n(\delta_s \mathcal{F})(a) = \partial_n(\delta_s \mathcal{F})(b) = 0\}. \end{aligned} $$
$$
\left\{
\begin{array}{ll}
\text{Find } (s, \delta_s \mathcal{F}) \text{ in } V_1 \times V_2 \text{ such as:} & \\
\displaystyle \int_{\Omega} \partial_t s.v = \int_{\Omega} K(s) \nabla \delta_s \mathcal{F} \cdot \nabla v + \int_{\Omega} \chi(s).v, & \text{for all } v \in H^1(\Omega) \text{ such as } v = 0 \text{ on } \Gamma_{in}, \\
\displaystyle \int_{\Omega} \delta_s \mathcal{F} \omega = - \int_{\Omega} \Pi(s) \omega - \int_{\Gamma_{in}} T \cdot \partial_n(s) \omega + \int_{\Omega} T \nabla s \cdot \nabla \omega, & \text{for all } \omega \in H^1(\Omega).
\end{array}
\right.
$$

Similarly, for a NH input and Dirichlet output condition on saturation s, we obtain the following problem, with: $V_1 = \{s \in H^1(\Omega) \text{ such as } s(b) = s_{out} \text{ and } \partial_n s(a) = 0\}$ and $V_2 = \{\delta_s \mathcal{F} \in H^1(\Omega) \text{ such as } \partial_n(\delta_s \mathcal{F})(a) = \partial_n(\delta_s \mathcal{F})(b) = 0\}$

$$
\left\{
\begin{array}{ll}
\text{Find } (s, \delta_s \mathcal{F}) \text{ in } V_1 \times V_2 \text{ such as:} & \\
\displaystyle \int_{\Omega} \partial_t s.v = \int_{\Omega} K(s) \nabla \delta_s \mathcal{F} \cdot \nabla v + \int_{\Omega} \chi(s).v, & \text{for all } v \in H^1(\Omega) \text{ such as } v = 0 \text{ on } \Gamma_{\text{out}}, \\
\displaystyle \int_{\Omega} \delta_s \mathcal{F} \omega = - \int_{\Omega} \Pi(s) \omega - \int_{\Gamma_{\text{out}}} T \cdot \partial_n(s) \omega + \int_{\Omega} T \nabla s \cdot \nabla \omega, & \text{for all } \omega \in H^1(\Omega).
\end{array}
\right.
$$

Periodic conditions for s and for free energy $\delta_s \mathcal{F}$

We impose a periodic condition for $s$ and $\delta_s \mathcal{F}$ such that $s(a) = s(b)$ and $\delta_s \mathcal{F}(a) = \delta_s \mathcal{F}(b)$. We obtain the problem:

$$
\left\{
\begin{array}{l}
\displaystyle \int_{\Omega} \partial_t s.v = \int_{\Omega} K(s) \nabla \delta_s \mathcal{F} \cdot \nabla v + \int_{\Omega} \chi(s).v, \\
\\
\displaystyle \int_{\Omega} \delta_s \mathcal{F} \omega = - \int_{\Omega} \Pi(s) \omega + \int_{\Omega} T \nabla s \cdot \nabla \omega.
\end{array}
\right.
$$

This is equivalent to solve the following problem, with $V_1, V_2$ defined as follows:
$V_1 = \{s \in H^1(\Omega) \text{ such as } s(a) = s(b)\}$ and $V_2 = \{\delta_s \mathcal{F} \in H^1(\Omega) \text{ such as } \delta_s \mathcal{F}(a) = \delta_s \mathcal{F}(b)\}$

$$
\left\{
\begin{array}{ll}
\text{Find } (s, \delta_s \mathcal{F}) \text{ in } V_1 \times V_2 \text{ such as:} & \\
\displaystyle \int_{\Omega} \partial_t s.v = \int_{\Omega} K(s) \nabla \delta_s \mathcal{F} \cdot \nabla v + \int_{\Omega} \chi(s).v, & \text{for all } v \in H^1(\Omega), \\
\displaystyle \int_{\Omega} \delta_s \mathcal{F} \omega = - \int_{\Omega} \Pi(s) \omega + \int_{\Omega} T \nabla s \cdot \nabla \omega, & \text{for all } \omega \in H^1(\Omega).
\end{array}
\right.
$$

C.1.3 NH condition on free energy input, and a flux condition on output

Taking now the case where we impose an outgoing flow and a homogeneous Neumann condition on $\delta_s \mathcal{F}$, the system (C.5) becomes:

$$
\left\{
\begin{array}{l}
\displaystyle \int_{\Omega} \partial_t s.v = \int_{\Omega} \chi(s).v - \int_{\Omega} K(s) \nabla \delta_s \mathcal{F} \cdot \nabla v + \int_{\Gamma_{\text{out}}} (K(s) \partial_n \delta_s \mathcal{F}) v, \\
\\
\displaystyle \int_{\Omega} \delta_s \mathcal{F} \omega = - \int_{\Omega} \Pi(s) \omega + \int_{\Omega} T \nabla s \cdot \nabla \omega - \int_{\Gamma_{\text{in}} \cup \Gamma_{\text{out}}} T \cdot \partial_n(s) \omega.
\end{array}
\right.
\qquad (C.7)
$$
Now we will do the case disjunction on s:

* For Dirichlet conditions on s input and output, the system (C.7) remains unchanged. So the problem will be with: $V_1 = \{s \in H^1(\Omega) \text{ such as } s(a) = s_{in} \text{ and } s(b) = s_{out}\}$ and $V_2 = \{\delta_s \mathcal{F} \in H^1(\Omega) \text{ such as } \partial_n \delta_s \mathcal{F}(a) = 0 \text{ and } K(s)\nabla(\delta_s \mathcal{F})(b) = 0\}$

$$ \begin{cases} \text{Find } (s, \delta_s \mathcal{F}) \text{ in } V_1 \times V_2 \text{ such as:} \\ \displaystyle \int_{\Omega} \partial_t s \cdot v = \int_{\Omega} \chi(s) \cdot v - \int_{\Omega} K(s) \nabla \delta_s \mathcal{F} \cdot \nabla v + \int_{\Gamma_{out}} (K(s) \partial_n \delta_s \mathcal{F}) v, & \text{for all } v \in H_0^1(\Omega), \\ \displaystyle \int_{\Omega} \delta_s \mathcal{F} \omega = - \int_{\Omega} \Pi(s) \omega + \int_{\Omega} T \nabla s \cdot \nabla \omega - \int_{\Gamma_{in} \cup \Gamma_{out}} T \cdot \partial_n(s) \omega, & \text{for all } \omega \in H^1(\Omega). \end{cases} $$

* For homogeneous Neumann conditions on s input and output, we have:
  $V_1 = \{s \in H^1(\Omega) \text{ such as } \partial_n s(a) = 0 \text{ and } \partial_n s(b) = 0\}$ and
  $V_2 = \{\delta_s \mathcal{F} \in H^1(\Omega) \text{ such as } \partial_n \delta_s \mathcal{F}(a) = 0 \text{ and } K(s)\nabla(\delta_s \mathcal{F})(b) = 0\}$

$$ \begin{cases} \text{Find } (s, \delta_s \mathcal{F}) \text{ in } V_1 \times V_2 \text{ such as:} \\ \displaystyle \int_{\Omega} \partial_t s \cdot v = \int_{\Omega} \chi(s) \cdot v - \int_{\Omega} K(s) \nabla \delta_s \mathcal{F} \cdot \nabla v + \int_{\Gamma_{out}} (K(s) \partial_n \delta_s \mathcal{F}) v, & \text{for all } v \in H^1(\Omega), \\ \displaystyle \int_{\Omega} \delta_s \mathcal{F} \omega = - \int_{\Omega} \Pi(s) \omega + \int_{\Omega} T \nabla s \cdot \nabla \omega, & \text{for all } \omega \in H^1(\Omega). \end{cases} $$

* For a Dirichlet condition at the input and a homogeneous Neumann condition at the output:
  $V_1 = \{s \in H^1(\Omega) \text{ such as } s(a) = s_{in} \text{ and } \partial_n s(b) = 0\}$ and
  $V_2 = \{\delta_s \mathcal{F} \in H^1(\Omega) \text{ such as } \partial_n \delta_s \mathcal{F}(a) = 0 \text{ and } K(s)\nabla(\delta_s \mathcal{F})(b) = 0\}$

$$ \begin{cases} \text{Find } (s, \delta_s \mathcal{F}) \text{ in } V_1 \times V_2 \text{ such as:} \\ \displaystyle \int_{\Omega} \partial_t s \cdot v = \int_{\Omega} \chi(s) \cdot v - \int_{\Omega} K(s) \nabla \delta_s \mathcal{F} \cdot \nabla v + \int_{\Gamma_{out}} (K(s) \partial_n \delta_s \mathcal{F}) v, \\ \text{for all } v \in H^1(\Omega) \text{ such as } v = 0 \text{ on } \Gamma_{in}, \\ \displaystyle \int_{\Omega} \delta_s \mathcal{F} \omega = - \int_{\Omega} \Pi(s) \omega + \int_{\Omega} T \nabla s \cdot \nabla \omega - \int_{\Gamma_{in}} T \cdot \partial_n(s) \omega, & \text{for all } \omega \in H^1(\Omega). \end{cases} $$

**Remark C.1.1.** For a homogeneous Neumann condition at the input and a Dirichlet condition at the output, we find the same form as this one, with only the term $\int_{\Gamma_{in}} T \cdot \partial_n(s) \omega$ replaced by $\int_{\Gamma_{out}} T \cdot \partial_n(s) \omega$ and $\partial_n s(a) = 0, s(b) = s_{out}$.

**Remark C.1.2.**
* Considering a NH condition on the free energy output and a flux condition on the input, we simply repeat the results of this subsection and replace $\int_{\Gamma_{out}} (K(s) \partial_n \delta_s \mathcal{F}) v$ by $\int_{\Gamma_{in}} (K(s) \partial_n \delta_s \mathcal{F}) v$.

* Considering an input-output flow condition, we simply repeat the results of this subsection and replace $\int_{\Gamma_{out}} (K(s) \partial_n \delta_s \mathcal{F}) v$ by $\int_{\partial \Omega} (K(s) \partial_n \delta_s \mathcal{F}) v$.
# Appendix D

## Numerical resolution in time

Solving the matrix problem, Eq. (3.9), is performed using numerical solution algorithms.

There are two main categories for solving the system numerically: direct methods and iterative methods.

The direct method to solve a linear system is an algorithm which would give the solution in a finite number of operations. There are also iterative methods which involve constructing a sequence of vectors $u_n^k$ converging towards the solution $u_{n+1}$, i.e. when the number of iterations $k \to +\infty$ then $u_n^k \to u_{n+1}$.

In our case, the problem is non-linear. Iterative methods are used to define the solution. Figure 3.8 represents the diagram associated with solving the matrix system using iterative methods.

An iterative algorithm is developed to define the solution at time $t_{n+1}$, this means using a resolution loop of $k$ iterations. The initialisation of this algorithm is based on the solution at time $t_n$. The solution obtained at the first iteration $u_n^0$ is defined using an explicit scheme to approximate the solution at time $t_{n+1}$. Therefore, it is denoted $u_{n+1}^0 = u_n^P$.

The resolution loop continues to run until the profile has converged, i.e. until the residual $r = b(u) - A(u)\delta u$ obtained from the matrix resolution has a norm less than a sufficiently small $\epsilon$ tolerance.

However, limiting the algorithm to this single criterion means that it may converge after a high number of iterations, i.e imposing a very long resolution time, or diverge because the choice of time step is incorrect (section 3.3.2). In order to prevent these two situations, a finite number of iterations $k_{itermax}$ is imposed. Thus, if the norm of the residual is less than $\epsilon$ and $k < k_{itermax}$ then the iterative solution algorithm has converged and $u_{n+1} = u_n^k$. If $k > k_{itermax}$ and the norm of the residual is greater than $\epsilon$ then the system has not converged.

Now that we have presented the general principle of the iterative algorithm, it is necessary to choose a resolution algorithm to establish the value of the residual. In the remainder of this section, we propose to focus on the Picard and Newton methods.
## D.1 Picard method

Picard method is an iterative method for solving non-linear systems of the form [208]:

$$F(u) = \frac{\partial u}{\partial t} - \nabla \cdot (\alpha(u)\nabla u) - f(u) = 0, \quad \text{in } \Omega \times (0, T). \qquad (D.1)$$

This system can be expressed in a matrix form as mentioned in section 3.2,

$$A(u)u = b(u) \qquad (D.2)$$

Eq. (D.2) is solved in linearized form starting from an initial value of $u$. Let $u^{-}$ be this initial value and $u$ the solution of the system Eq. (D.1) at time $t > 0$. In other words, this method solves a problem of the form:

$$\hat{F}(u, u^{-}) = 0, \qquad (D.3)$$

where $\hat{F}$ is the linear function of $F$. It is defined by

$$\hat{F}(u, u^{-}) = \frac{\partial u}{\partial t} - \nabla \cdot (\alpha(u^{-})\nabla u) - f(u^{-}). \qquad (D.4)$$

This formulation is governed by a matrix form such that Eq. (D.3) is expressed in the form of Eq. (D.2),

$$A(u^{-})u = b(u^{-}). \qquad (D.5)$$

The matrix $A(u^{-})$ and the vector $b(u^{-})$ are the matrix forms defined in section 3.2 and described according to the time scheme chosen in Appendix B. Thus, the resolution required to use Picard method is to solve the given matrix form without modifying the system.

During the iterations of the solution, the computed solution can be very different from the solution predicted by the explicit scheme. To prevent this, we use a relaxation method. A relaxation method adapts the solution $u^{*}$ obtained at each iteration by making it depend on the solution obtained at the previous iteration. In this way, the evolution of the solution is not too violent. For this reason, the solution defined after each iteration is given by the relationship:

$$u = u^{*}\omega + (1 - \omega)u^{-}. \qquad (D.6)$$

In our study, we consider $\omega = 1$.

This means that we solve $\hat{F}(u^{*}, u^{-}) = 0$, i.e. $A(u^{-})u^{*} = b(u^{-})$ which defines the solution of the system, governed by Eq. (D.2), given by Eq. (D.6).

A general algorithm associated with Picard method is given below. This algorithm is one of the versions developed and tested to solve our system governed by Eq. (3.5).

## D.2 Newton method

Newton method, like Picard method, can be used to solve a system in the form of Eq. (D.1). However, this method requires the use of the Jacobian. The solution is approximated by estimating the tangent
Algorithm for relaxed Picard iterations

**Initialization:**

$u^{-} = u_{n+1}^{P}$: initial guess

$r = A(u^{-})u^{-} - b(u^{-})$

$\omega$ relaxed value

**while** $||r|| > \epsilon$ **do**

solve $A(u^{-})u^{\star} = b(u^{-})$

$u = u^{\star}\omega + (1 - \omega)u^{-}$

$u^{-} \leftarrow u$

$r = A(u^{-})u^{-} - b(u^{-})$

**end while**

**Output:** $u_{n+1} = u^{-}$ solution of equation (D.2)

of the function $F$ at the solution $u$. In other words, we are trying to solve the following system:

$$ \mathbf{D}_u F(u) = 0. \qquad (D.7) $$

This weakly formulated system can be expressed in matrix form,

$$ J(u)u = r(u), \qquad (D.8) $$

where $r(u) = A(u)u - b(u)$ is the residual of the system governed by Eq. (D.2) and $J$ is the Jacobian matrix associated with the differentiable form of $F$ at $u$, $\mathbf{D}_u F$.

Eq. (D.8) is solved in linearized form from a value of $u$. Using the same notations as those defined in Picard method, this amounts to solve

$$ \mathbf{D}_u \hat{F}(u, u^{-}) = 0. \qquad (D.9) $$

This formulation is governed by a matrix form such that Eq. (D.3) is expressed in the form of Eq. (D.8),

$$ J(u^{-})u = r(u^{-}). \qquad (D.10) $$

The matrix $J(u^{-})$ depends on the time scheme chosen.

The relaxation method, given by Eq. (D.6), is also used for Newton method. Like Picard method, this method allows us to remain in the neighbourhood of the solution predicted by an explicit scheme and not to diverge or converge towards a solution that does not correspond to our dynamics.

This means that we solve $\mathbf{D}_u \hat{F}(u^{\star}, u^{-}) = 0$, i.e. $J(u^{-})u^{\star} = r(u^{-})$ in order to define the solution of the system, governed by Eq. (D.8), given by Eq. (D.6).

A general algorithm associated with Newton method is given below. It is one of the versions developed and tested to solve our system governed by Eq. (3.5).
Algorithm for relaxed Newton iterations

Initialization:

$u^{-} = u_{n+1}^{P}$: initial guess

$r = A(u^{-})u^{-} - b(u^{-})$

$\omega$ relaxed value

while $||r|| > \epsilon$ do

solve $J(u^{-})u^{*} = r$

$u = u^{*}\omega + (1 - \omega)u^{-}$

$u^{-} \leftarrow u$

$r = A(u^{-})u^{-} - b(u^{-})$

end while

Output: $u_{n+1} = u^{-}$

## D.3 Synthesis

Our system governed by Eq. (3.5) is iteratively solved. The Picard and Newton algorithms have been presented to solve it. Both methods use a number of iterations to converge to the solution depending on the choice of the first element of the sequence at initialization time. Naturally, the solution must prealably approximated using an explicit scheme. In addition, both methods have advantages and disadvantages.

Picard method is fast to implement and requires relatively small storage space. The matrix form defined by Eq. (3.9) is directly reused in the Picard method. However, its resolution time is relatively long because of the large number of iterations required for the system to converge.

Newton method requires the Jacobian matrix of the system to be defined. The construction of this matrix requires more time and storage than Picard method. However, its convergence time is relatively faster due to its low number of iterations.

Depending on the user's needs for solving the system, both methods can be used. For fast tests that do not require a lot of construction time and low memory storage, Picard method is recommended, whereas for faster long-term simulation, Newton method is better. An algorithmic version giving the choice of the method has been implemented using the Picard/Newton algorithms. It is described below.
Algorithm for relaxed Picard/Newton iterations

Initialization:

$u^{-} = u_{n+1}^{P}$: initial guess

$r = A(u^{-})u^{-} - b(u^{-})$

$\gamma$

$\omega$ relaxed value

while $||r|| > \epsilon$ do

if $\gamma = 1$ then

$\tilde{b} = r$

else if $\gamma = 0$ then

$\tilde{b} = b$

end if

$\tilde{A} = A(u^{-}) + \gamma (J(u^{-}) - A(u^{-}))$

solve $\tilde{A}u^{\star} = \tilde{b}$

$u = u^{\star}\omega + (1 - \omega)u^{-}$

$u^{-} \leftarrow u$

$r = A(u^{-})u^{-} - b(u^{-})$

end while

Output: $u_{n+1} = u^{-}$
# Appendix E

Chemical diversity of crop plants root mucilages. Implications for their maximal water content and mineralization in soils

This appendix refers to the shortnote of 15 pages, accepted with minor correction, in the journal Rhizosphere. Note that the aims of this journal is to advance the frontier of our understanding of plant-soil interactions.
1 Chemical diversity of crop plants root mucilages. Implications for their

2 maximal water content and mineralization in soils.

3

4 Abstract

5 The biophysical functioning of the rhizosphere can act on plant water use efficiency through root
6 mucilage release and their stimulation of microorganisms' activity. Mucilages characteristics and roles
7 are still poorly studied for the diversity of agronomic plant species. We compared mucilages collected
8 from roots of germinated seeds (aerohydroponics technique) of eight plant species (four
9 Eudicotyledons, four Monocotyledon) in terms of total sugars content, MIR spectra, mucilage maximal
10 water content and induced respiration in contrasting soils. Our results suggest that the "chemical
11 fingerprint" of these mucilages could be discriminated according to the phylogenetic proximity of the
12 plant species. In addition, the maximal water content that mucilages retain seems to be linked to their
13 chemical composition and seems more related to the presence of high molecular weight sugars than
14 their total sugar amounts. If mucilage-induced respiration by soil microbiota appeared to be
15 independent from phylogeny, some mucilages induced more respiration than others regardless of the
16 soil studied. Microbial communities and soil physico-chemical properties interact in mineralization
17 with variations in mucilages' chemical composition.

18 **Keywords:** root, mucilage, agronomic plant species, MIR spectra, EPS, water content,
19 mineralization

20

21 Improvements in the efficiency of water use and resistance to water deficits are required for the
22 agroecological transition (Ahmed et al., 2018) and need new production systems as well as new
23 cultivars. To this aim, we need to understand the functional traits of the plant that could improve its
24 resistance to water deficit.

25 In relation with the root system (Vetterlein and Doussan, 2016), the biophysical functioning of the
26 rhizosphere (e.g variation in soil water retention, wettability, porosity, microorganisms activity, due to
27 root activity) is a characteristic that could modify plant water use efficiency. Peculiar biophysical
28 characteristics of rhizospheric soil (water retention and aggregation) are thought to be linked in part
29 to exopolymers (EPS) produced by plant roots via their mucilages and also by rhizospheric
30 microorganisms, themselves stimulated by rhizodeposits (including mucilages, Brax et al., 2019). Root
31 mucilages are mainly composed of exopolysaccharides (EPSac), with minor amount of proteins,
32 minerals and lipids (Nazari, 2021). Their fate in the soil in relation to their mineralization and according
to their plant origin has been little studied to date (Naveed et al., 2017) and no comparison between different soils (with different microbial communities) has been made to our knowledge. In relation with the diversity of agronomic plant species, their variability in characteristics and roles are however poorly studied (Galloway et al., 2020).

Our aim was to characterize the root mucilages of various crop species, through their chemistry (sugar content, MIR spectra), physical parameter (maximal water content) and mineralization potential in contrasting soils (MicroResp™ bioassays).

We used an aero-hydroponic method adapted from Zickenrott et al. (2016), Holz et al. (2018) and Brax et al. (2019) to harvest the mucilage of eight plant species of agronomic interest (four Eudicotyledons, four Monocotyledons, Tab. 1). After sterilization the seeds (Isopropanol 70%, NaOCl 1% then tween 20 0.1%, followed by a wash with ultrapure water), the seeds were germinated on a tray (3x3 mm stainless steel mesh, adding a very thin layer of glass wool for very small seeds such as tomato, rapeseed and sorghum), incubated in the dark at 25°C +/- 1°C, in an hermetic box above a solution of CaSO4 (2 mM) aerated continuously with two pumps. After 2-3 days of incubation, mucilage is collected from the tip of each root daily using a baby fly connected to a suction pump over a further 1-2 days, until the roots reach the CaSO4 solution. For each species, five growth box systems were set up and run in parallel. From four growth boxes, the mucilage was collected and pooled to obtain enough material for analysis and experimentation. Globally, the mucilage pools were obtained from 1100 to 29800 germinated seeds, depending on the species tested (Tab.1). The fifth growth box was used to get the amount of mucilage (in total sugar equivalent) produced per dry mass of roots. In this box, after the first mucilage collection, the roots were cut along the grid surface, dried and weighed.

The collected mucilage was pooled into tubes, weighed and first stored at -20°C. An aliquot was freeze-dried to measure mucilage maximal water content (MWC, calculated as mucilage wet weight minus dry weight divided by dry weight, Nazari et al., 2020) and for analyses with Medium Infrared (MIR) spectroscopy. Total sugars were analyzed using the Dubois method (Dubois et al., 1956). MIR measurements were performed with a Tensor 27 FTIR spectrometer (Bruker Optics, Wissembourg, France). The samples were scanned between 4000 and 600 cm⁻¹ and each spectrum was obtained by averaging 16 successive scans. For each sample, four replications on different aliquots were made (Bérard et al., 2020). MIR data were processed with ChemFlow chemometrics software (GalaxyProject). Principal component analysis (PCA) of the MIR spectra of the eight species mucilage showed that the first two axes explained 82 % of variance.

We also performed a MicroResp™ assay (Bérard et al., 2011; De Vries et al., 2019) of respiration induced by the different mucilages in three contrasting soils selected for their distinct microbial biomass, texture, organic C, and water holding capacity. These include a sandy clay loam from an irrigated
67 permanent grassland ("Crau"), a silty clay loam soil from a conventional agricultural field ("Lysi") and a
68 sandy soil from a pine forest ("Tavel") (Bérard et al., 2020, SuppMat1). We introduced 25 µl of fresh
69 mucilage per well of the microplate system (containing 0.48 to 0.57 g in dry mass equivalent per well
70 depending on the soil type), with 4 replications per mucilage and per soil tested (the moisture content
71 of the test soils was 40% of their water holding capacity). Respiration measurements were recorded
72 after six hours' incubation at 23°C and respiration induced by each mucilage was normalized by the total
73 sugars amount introduced per soil sample. Non-parametric tests (Kruskal Wallis, Dunn Bonferroni test,
74 p<0.05) were used for data analysis.

75 From the results obtained, the variations of the eight species root mucilage in its chemical, physical and
76 degradation characteristics enabled to address the three points below:

77 - The chemical composition of plant root mucilages seems related to their phylogenetic proximity

78 The PCA of MIR measurements shows that the species are relatively well separated from each other
79 (Fig.1). The dendrogram of relative distances between points effectively separates Monocotyledons
80 (negative part of PC1) from Eudicotyledons (positive part of PC1), Asterids (Tomato) from Rosids,
81 Malvids (Rapeseed) from Fabids, and Roinoids (Lupins) from Genisteae (Alfalfa) (SuppMat2). In the
82 Monocotyledons, Panicoideae having C4 photosynthesis strategy (Maize and Sorghum) are at the
83 bottom left of the PCA, and Pooideae with C3 photosynthesis (Wheat and Barley) at the top left. This
84 suggests that these plant mucilage chemical compositions are related to phylogenetic proximity, at
85 least up to the family scale, corroborating earlier studies by (Vančura and Hovadík, 1965) who
86 observed similar molecular compositions of root exudates whose plant species were closely related.
87 More recently, similar observations have been made on the phylogeny-chemical composition link of
88 root exudates from different plant families using GC-MS (Herz et al., 2018; Williams et al., 2022) or
89 using ESI FT-ICR MS (Miao et al., 2020).

90 According to MIR spectra acquired, dicots are located on the positive part of PC1 (Fig.1a),
91 characterized by a high absorption in a broad band of wave numbers between 1217 and 1681 cm⁻¹
92 (Fig.1b). This zone is largely represented by nitrogenous molecules of the amine, amide (I, II, III types),
93 proteins, as well as phosphates and fats (phospholipids, triglycerides, lipids, fatty acids). Benzene
94 (1328; 1489 cm⁻¹), phenyl (1500 -1510 cm⁻¹) and quinoid (1588 cm⁻¹) rings are also in this band. Wave
95 numbers between 1419 and 1444 cm⁻¹ indicate pectin- and cellulose-type polysaccharides (Talari et
96 al., 2017).

97 The negative part of PC1 is represented by Monocots-Poaceae with MIR absorption bands at 932-1064
98 cm-1 and 1140 and 1198 cm-1. These zones are mainly represented by simple carbohydrates (sucrose,
99 mannose, galactose, fructose), polysaccharides (pectin, cellulose, glycogen) and, to a lesser extent,
100 phosphate, phosphodiester, amino acid (thyrosine) and aromatic molecules (Talari et al., 2017).
Setting aside the redundancies between spectra regions, these results suggest that PC1 is structured primarily by a gradient of fat abundance (lipids, fatty acids, triglyceride; wave numbers above 1300 cm⁻¹). Nitrogenous molecules also drive PC1, with an abundance of tertiary amides to its positive part (wave numbers above 1400 cm⁻¹). This joins the observations of (Miao et al., 2020): the degree of unsaturation and the proportion of nitrogen in the molecules studied were the two main factors discriminating species mucilage chemical composition according to their phylogenetic families. So although mucilage from different phylogenetic families possess a "largely common composition" (Herz et al., 2018), particularly for their carbohydrate base (Vančura and Hovadík, 1965), the abundance of fat and nitrogenous molecules (amide III), albeit minor in mucilage composition (Nazari, 2021; Miao et al., 2020) seem the primary factors in their variability.

## 111 - The chemical composition of mucilages seems linked to their MWC

A positive correlation arises ($R^2=0.568, p=0.051$) between amount of sugar (per g of mucilage dry matter) and maximal water content of mucilage, only if maize is excluded. Nazari et al. (2020) also found no correlation for maize cultivars. Furthermore, the PC2 axis structured by a wavelength band between 1100 and 1250 cm⁻¹, corresponding to polysaccharide (pectin and cellulose) and sugar type molecules (Talari et al., 2017), appears to rank mucilages according to the quantity of high molecular weight sugar (at the top is the tomato whose exudates are rich in pectin, Jeong et al., 2018) and therefore also to their ability to retain water (Kazanskii and Dubrovskii, 1992). This is underlined by the positive correlation between the coordinates of the PC2 axis and MWC values ($R^2=0.623, p<0.0001$) (Fig.2).

## 121 - Mucilage-induced respiration depends on the type of soil and probably also on its plant origin.

Irrespective of the mucilage tested, we observed a ranking of respirations according to the 3 soils, with mucilage-induced respirations systematically higher in Crau and Lysi soils than in Tavel soil (p<0.0001). This shows an effect of microbial communities and their edaphic environment (with Crau>Lysi>Tavel for MO content, Microbial Biomass, CEC, WHC and Basicity, Bérard et al., 2020; Riffaldi et al., 1996) greater than the effect of the mucilage's chemical composition on its induced respiration. However, some mucilages are classified (in terms of induced respiration intensity) in the same way whatever the soil (wheat, lupin, maize mucilages being the least degraded, while tomato, rapeseed mucilages the most - Fig. 2). In particular, respiration induced by tomato mucilage was always much higher than that induced by maize mucilage (p<0.002). The water content of tomato mucilage is higher than that of maize mucilage (Tab.1), probably making the first more available to microorganisms. In addition, these two mucilages are very far apart on the PCA of the MIR spectra (Fig.1), suggesting different compositions. To a lesser extent than Naveed et al. (2017), we observe the same trend when comparing maize and barley mucilages. This suggests that the chemistry and/or physics of the
mucilage have a dominant role on its availability/degradation by soil microorganisms. It would be interesting to develop chemical and physical analyses to characterize, for example, the size and branching of the polymers, and their water retention and viscosity (parameters likely to influence the degradability of these mucilages) in order to confirm these hypotheses. On the other hand, plant phylogeny does not seem to be related to mucilage-induced respiration levels. Finally, in some cases, mucilages-induced respiration differently according to soil type. For example, rapeseed mucilage-induced respiration is higher than wheat or lupin mucilage in Crau soil (p<0.001), whereas this difference fades in Lysi and Tavel soils (p>0.05). The soil, through its microbial community and its physicochemical properties, can therefore interact with the degradation of the different mucilages.

The sugar amount of mucilage produced per seed is extremely variable (Tab.1) and not associated with phylogeny but strongly correlated with seed weight (R²= 0.92, p=0.0001). This suggests that species involve more or less the same fraction of seed's C reserve in sugars of mucilage production. As a corollary, this also results in a large range of amounts of mucilage-derived sugars delivered to soil among species (according to seed weights), which may impact considerably in soil the shaping and interactions in this "early stage rhizosphere" at a critical stage in plant establishment. Analyzes of C stored in the seeds, in relation with mucilage-derived sugars could strengthen these findings. We also observed that the sugars amount produced per root weight was highly variable between different species (6-fold, Tab.1), which in real plant conditions in soils could impact rhizospheric soil water retention. This exploratory study, based on mucilage harvested in aero-hydroponics from various plant species, needs to be extended to later phenological stages in real field conditions, with a variety of soils and water availability conditions, as done for maize (Nazari et al., 2023) and tomato (Le Gall et al., 2021).

**Acknowledgements:** This work was made possible thanks to the support of our institute, INRAE, for a post-doctoral grant and the funding of apprenticeship training for engineers, and to the funding of a doctoral grant by INRAE and the University of Avignon. Special thanks to the SQPOV-INRAE unit (S. Bureau and C. Le Bourvellec) for making their FTIR spectrometer available to us and for their advice on the MIR spectral measurements and data.

## References

Ahmed, M.L.A., Passioura, J., Carminati, A., 2018. Hydraulic processes in roots and the rhizosphere pertinent to increasing yield of water-limited grain crops: a critical review. J. Exp. Bot. 69, 3255–3265. https://doi.org/10.1093/jxb/ery183.
169 Galloway, A.F., Knox P., Krause K., 2020. Sticky mucilages and exudates of plants: putative
170 microenvironmental design elements with biotechnological value. New Phytologist. 225, 1461-
171 1469. https://doi.org/10.1111/nph.16144.

172 Holz, M., Leue, M., Ahmed, M.A., Benard, P., Gerke, H.H., Carminati, A., 2018. Spatial Distribution
173 of Mucilage in the Rhizosphere Measured With Infrared Spectroscopy. Front. Environ. Science.
174 6, 1-7. https://doi.org/10.3389/fenvs.2018.00087.

175 Bérard, A., Clavel, T., Le Bourvellec, C., Davoine, A., Le Gall, S., Doussan, C., Bureau, S., 2020.
176 Exopolysaccharides in the rhizosphere: A comparative study of extraction methods. Application
177 to their quantification in Mediterranean soils. Soil Biology and Biochemistry. 149, 1-10.
178 https://doi.org/10.1016/j.soilbio.2020.107961.

179 Bérard, A., Bouchet, T., Sevenier, G., Pablo, A.L., Gros, R., 2011. Resilience of soil microbial
180 communities impacted by severe drought and high temperature in the context of
181 Mediterranean heat-waves. Eur. J. Soil Biol. 47, 333-342.
182 https://doi.org/10.1016/j.ejsobi.2011.08.004.

183 Brax, M., Buchmann, C., Kenngott, K., Schaumann, G.H., Diehl, D., 2019. Influence of the physico-
184 chemical properties of root mucilage and model substances on the microstructural stability of
185 sand. Biogeochem. 147, 35-52. https://doi.org/10.1007/s10533-019-00626-w.

186 Dubois, M., Gilles, K.A., Hamilton, J.K., Rebers, P.A., Smith, F., 1956. Colorimetric Method for
187 Determination of Sugars and Related Substances. Analytical Chem. 28, 3, 350-356.
188 https://doi.org/10.1021/ac60111a017.

189 Herz, K., Dietz, S., Gorzolka, K., Haider, S., Jandt, U., Scheel, D., Bruelheide, H., 2018. Linking Root
190 Exudates to Functional Plant Traits. PLOS ONE. 14, 3.
191 https://doi.org/10.1371/journal.pone.0204128.

192 Jeong, H.Y., Nguyen, H.P., Eom, S.H., Lee, C., 2018. Integrative Analysis of Pectin Methylesterase
193 (PME) and PME Inhibitors in Tomato (Solanum Lycopersicum): Identification, Tissue-Specific
194 Expression, and Biochemical Characterization. Plant Physiol. Biochem. 132, 557-65.
195 https://doi.org/10.1016/j.plaphy.2018.10.006.

196 Kazanskii, K. S., and S. A. Dubrovskii. 1992. Chemistry and Physics of 'Agricultural' Hydrogels.
197 Polyelectrolytes Hydrogels Chromatographic Materials. Advances in Polymer Science. 97-133.
198 https://doi.org/10.1007/3-540-55109-3_3.

199 Le Gall, S., Bérard, A., Page, D., Lanoë, L., Bertin, N., Doussan, C., 2021. Increased exopolysaccharide
200 production and microbial activity affect soil water retention and field performance of tomato
201 under water deficit. Rhizosphere. 19, e100408. https://doi.org/10.1016/j.rhisph.2021.100408.
202 Miao, Y., Lv, J., Huang, H., Cao, D., Zhang, S., 2020. Molecular Characterization of Root Exudates Using Fourier Transform Ion Cyclotron Resonance Mass Spectrometry. J. Environm. Sci. 98, 22–30. https://doi.org/10.1016/j.jes.2020.05.011.

205 Naveed, M., Brown, L. K., Raffan, A. C., George, T. S., Bengough, A. G., Roosevelt, T., Sinclair, I., Koebernick, N., Cooper, L., Hackett, C. A., Hallett, P.D., 2017. Plant exudates may stabilize or weaken soil depending on species, origin and time. Eur. J. Soil Sci. 68, 806–816. https://doi.org/10.1111/ejss.12487.

209 Nazari, M., Bilyera, N., Banfield, C.C., Mason-Jones, K., Zarebanadkouki, M., Munene, R., Dippold, M.A., 2023. Soil, climate, and variety impact on quantity and quality of maize root mucilage exudation. Plant Soil. 482, 25–38. https://doi.org/10.1007/s11104-022-05669-x.

212 Nazari M., 2021. Plant mucilage components and their functions in the rhizosphere. Rhizosphere 18, 100344. https://doi.org/10.1016/j.rhisph.2021.100344

214 Nazari, M., Riebeling, S., Banfield, C.C., Akale, A., Crosta, M., Mason-Jones, K., Dippold, M. A., Ahmed, M.A., 2020. Mucilage Polysaccharide Composition and Exudation in Maize From Contrasting Climatic Regions. Front. Plant Sci. 11. https://doi.org/10.3389/fpls.2020.587610.

217 Riffaldi, R., Saviozzi, A., Levi-Minzi, R., 1996. Carbon mineralization kinetics as influenced by soil properties, Biol. Fertil. Soils. 22, 293-298. https://doi.org/10.1007/BF00334572.

219 Talari, A.C.S, Garcia Martinez, M.A., Movasaghi, Z., Rehman, S., Ur Rehman, I., 2017. Advances in Fourier transform infrared (FTIR) spectroscopy of biological tissues, App. Spectroscopy Reviews. 52, 5, 456-506. https://doi.org/10.1080/05704928.2016.1230863

222 Vančura, V., Hovadík, A., 1965. Root Exudates of Plants. Plant Soil. 22, 1, 21–32. https://doi.org/10.1007/BF01377686.

224 Vetterlein, D., Doussan, C., 2016. Root age distribution: how does it matter in plant processes? A focus on water uptake. Plant Soil. 407, 145–160. https://doi.org/10.1007/s11104-016-2849-6.

226 De Vries, F.T., Williams, A., Stringer, F., Willcocks, R., McEwing, R., Langridge, H., Straathof, A.L., 2019. Changes in root-exudate-induced respiration reveal a novel mechanism through which drought affects ecosystem carbon cycling. New Phytologist. 224, 1, 132-145 https://doi.org/10.1111/nph.16001.

230 Williams, A., Langridge, H., Straathof, A.L., Muhamadali, H., Hollywood, K.A., Goodacre, R., De Vries, F.T., 2022. Root Functional Traits Explain Root Exudation Rate and Composition across a Range of Grassland Species. J. Ecol. 110, 1, 21–33. https://doi.org/10.1111/1365-2745.13630.

233 Zickenrott, I.M., Woche, S.K., Bachmann, J., Ahmed, M.A., Vetterlein, D., 2016. An efficient method for the collection of root mucilage from different plant species—A case study on the effect of
235 mucilage on soil water repellency. J. Plant. Nutr. Soil Sci. 179, 294-302.

236 https://doi.org/10.1002/jpln.201500511.

237

238 **Figures and Table captions**

239 **Figure1:** Results of the Principal Component Analysis (PCA) performed on MIR spectra (four spectra per species) of mucilages harvested from the roots of aero-hydroponically germinated seeds of eight plant species. **a-** PCA representing the eight species. **b-** PCA eigenvector of the PC1 (black curve) and PC2 (red curve) axes of Fig.1a were built from the wavenumbers of the MIR spectra.

243 **Figure.2:** Correlation between the chemical characteristics of the different mucilages (harvested from the roots of aero-hydroponically germinated seeds of eight plant species) (x-axis: coordinates of the PC2 axis of the PCA performed with the MIR spectra, see Fig.1) and a physical characteristic (y-axis: mucilage maximal water content, see Tab.1) of these mucilages.

247 **Figure.3:** Respiration measurements induced by different mucilages (harvested from the roots of aero-hydroponically germinated seeds of eight plant species) in three contrasting soils (MicroResp™ technique). **a-** Respiration in Tavel soil. **b-** Respiration in Lysi soil. **c-** Respiration in Crau soil. Different letters indicate significant differences (four replications, Kruskal-Wallis, Dunn Bonferroni test). Tom: Tomato, Sor: Sorghum, Rap: Rape, Alf: Alfalfa, Bar: Barley, Lup: Lupin, Wht: Wheat, Mze: Maize.

252 **Table.1:** Phylogeny of the eight vegetal species studied, experiment characteristics, mucilages characteristics harvested from the roots of aero-hydroponically germinated seeds of these eight plant species. MWC= maximal water content (calculated as mucilage wet weight minus dry weight divided by dry weight).

256

257 **SuppMat.1:** Main physical, chemical and microbial properties of the three soils used for MicroResp™ bioassay with the eight different mucilages tested. SOC: Soil Organic Carbon, TN : Total Nitrogen, CEC: Cation Exchange Capacity, WHC: Water Holding Capacity, MB: Microbial Biomass, WSA : Water Stable Aggregates

261 **SupMat2.a:** Results of the Dendogram (Agglomerative Nesting Hierarchical Clustering) performed on MIR spectra (four spectra per species) of mucilages harvested from the roots of aero-hydroponically germinated seeds of eight vegetal species. **SupMat2.b:** Years of phylogenetic phyla from compilations
264 of articles from the Angiosperm Phylogeny Website hosted by the Missouri Botanical Garden:

265 http://www.mobot.org/MOBOT/Research/APweb/welcome.html.

266 **SuppMat.3** Discussion of the limitations of the mucilage collection method
<table><caption>**Table.1:** Phylogeny of the eight vegetal species studied, experiment characteristics, mucilages characteristics of these vegetal species. MWC= maximal water content (calculated as mucilage wet weight minus dry weight divided by dry weight).</caption><thead><tr><th rowspan="2">Species</th><th rowspan="2">Variety</th><th rowspan="2">Tribe</th><th rowspan="2">Sub-family</th><th rowspan="2">Family</th><th rowspan="2">Clade</th><th rowspan="2">Clade</th><th rowspan="2">Thousand seed weight (g)</th><th rowspan="2">Number of Seeds per tray</th><th rowspan="2">Germination time (day)</th><th colspan="4">Mucilage characteristics</th></tr><tr><th>Total sugar per dry root [eq glc] mg/g</th><th>Total sugar per dry matter [eq glc] mg/g</th><th>Total sugar per 1000 seed [eq glc] mg</th><th>MWC water [g] per dry weight [g]</th></tr></thead><tbody><tr><td>Maize</td><td>MAS 26.R</td><td>Triticeae</td><td>Pooideae</td><td>Poaceae</td><td>Commelinids</td><td>Monocot</td><td>286</td><td>277</td><td>3</td><td>0.89</td><td>8.16</td><td>1.97</td><td>353</td></tr><tr><td>Sorghum</td><td>Solarius</td><td>Triticeae</td><td>Pooideae</td><td>Poaceae</td><td>Commelinids</td><td>Monocot</td><td>11</td><td>2846</td><td>3</td><td>1.28</td><td>3.48</td><td>0.07</td><td>452</td></tr><tr><td>Wheat (Spring wheat)</td><td>Togado</td><td>Andropogoneae</td><td>Panicoidae</td><td>Poaceae</td><td>Commelinids</td><td>Monocot</td><td>37</td><td>865</td><td>3</td><td>0.28</td><td>6.31</td><td>0.38</td><td>566</td></tr><tr><td>Barley (Spring barley)</td><td>RGT Planet</td><td>Andropogoneae</td><td>Panicoidae</td><td>Poaceae</td><td>Commelinids</td><td>Monocot</td><td>43</td><td>794</td><td>3</td><td>0.95</td><td>5.41</td><td>0.65</td><td>761</td></tr><tr><td>Alfalfa</td><td>Europe</td><td>Robinoids</td><td>Faboideae</td><td>Fabaceae</td><td>Rosids</td><td>Eudicot</td><td>2</td><td>7456</td><td>2</td><td>0.70</td><td>1.37</td><td>0.05</td><td>140</td></tr><tr><td>Lupin</td><td>Rumba</td><td>Genistae</td><td>Faboideae</td><td>Fabaceae</td><td>Rosids</td><td>Eudicot</td><td>154</td><td>520</td><td>3</td><td>0.20</td><td>4.05</td><td>0.68</td><td>252</td></tr><tr><td>Rape</td><td>Cleopatra</td><td>Brassiceae</td><td>Brassicoideae</td><td>Brassicaceae</td><td>Rosids</td><td>Eudicot</td><td>5</td><td>3804</td><td>3</td><td>0.67</td><td>2.92</td><td>0.04</td><td>468</td></tr><tr><td>Tomato</td><td>Saint Pierre</td><td>Solanae</td><td>Solenoidae</td><td>Solanaceae</td><td>Asterids</td><td>Eudicot</td><td>2</td><td>6138</td><td>3</td><td>0.19</td><td>5.47</td><td>0.10</td><td>944</td></tr><tr><td colspan="12" style="text-align: center;">Mean for all species</td></tr><tr><td colspan="12" style="text-align: center;">standard deviation for all species</td></tr><tr><td colspan="12" style="text-align: center;">0.64</td></tr><tr><td colspan="12" style="text-align: center;">0.40</td></tr><tr><td colspan="12" style="text-align: center;">4.65</td></tr><tr><td colspan="12" style="text-align: center;">0.49</td></tr><tr><td colspan="12" style="text-align: center;">0.65</td></tr><tr><td colspan="12" style="text-align: center;">263</td></tr></tbody></table>
Figures 1-2-3

<figure><img src="image_114.png" /><figcaption>Figure 1</figcaption></figure>

Click here to access/download;Figure;Figures 1-2-3-revised.pptx
<figure><img src="image_115.png" /><figcaption>Figure 2</figcaption></figure>
<figure><img src="image_116.png" /><figcaption>Figure 3</figcaption></figure>
# Appendix F

## Material and methods : Designing a rhizosphere soil

This chapter presents in poster form the steps involved in designing a rhizosphere model medium (Figure F.1-F.2). First, we present the collection of maize mucilage, followed by the collection of chia seed mucilage. Then we present the method of washing sand with acid and finally the creation of the model medium with the sampling of a modified medium.

<figure><img src="image_117.png" /><figcaption>Figure F.1: Preparation cycle for not mixed soil</figcaption></figure>

<figure><img src="image_118.png" /><figcaption>Figure F.2: Preparation cycle for mixed soil</figcaption></figure>
Brace maize mucilage

<figure><img src="image_119.png" /></figure>
Protocole – Récolte de mucilage de graine de chia

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

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

En contact de
l'eau et par
agitation, les
graines de chia
subissent un
stress et produise
du mucilage. Ce
mucilage se dilue
dans l'eau dans
laquelle baignent
les graines. Il
forme un « gel »
lié aux graines

Après aspiration + 2 filtrations

<figure><img src="image_122.png" /><figcaption>Mucilage dilué</figcaption></figure>

La « matière fraîche » de mucilage obtenue peut ensuite être convertie en « matière sèche » en prélevant des petites quantités de matières fraîche et en les plaçant à l'étuve. La différence des masses permet d'obtenir un ratio

Matière sèche
Matière fraîche

<caption>Graines de chia</caption>
INRAO

Protocole – Lavage de sable

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

<figure><img src="image_124.png" /><figcaption>Sable de FontaineBleau</figcaption></figure>

<figure><img src="image_125.png" /><figcaption>200 microns</figcaption></figure>

<caption>400 microns</caption>

Dans un premier temps, le sable est tamisé entre 200 et 400 microns

Le sable doit subir un lavage à l'acide (HCL 12%) afin d'être entièrement purifié. Pour cela, il est placé dans un bac en plastique puis il est complètement immergé dans l'acide pendant une nuit. Pendant les 5 premières heures, le sable doit être mélangé régulièrement (au moins toutes les heures)

<figure><img src="image_126.png" /><figcaption>Sable + HCL
(12%)</figcaption></figure>

<figure><img src="image_127.png" /><figcaption>PH-mètre</figcaption></figure>

Sonde ph

PH-mètre

Après avoir retiré l'excédent d'acide, le sable peut être placé dans le tamis à 200 microns afin d'être rincé abondamment jusqu'à que le ph devienne neutre. Le ph est contrôlé à l'aide d'un ph mètre. Un dernier rinçage est ensuite effectué à l'eau distillé.

Une fois rincé, le sable est disposé dans un bac de séchage puis il est placé à l'étuve à 70° C pendant 2 jours

<figure><img src="image_128.png" /><figcaption>Etuve</figcaption></figure>
INRAO

Protocole - Mélange sable/mucilage →
création du sol artificiel

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

Au préalable, un calcul est effectué à partir du ratio matière sèche/matière fraiche afin de connaitre les quantités de sable et de matière fraiche à mélanger pour pouvoir obtenir la concentration souhaitée pour chaque sol

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

1 nuit au réfrigérateur

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

Un mélange doit être de nouveau effectué avant la mise en bac

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

Mise en échantillon (voir protocole)
INRAO

Protocole – Mise en échantillon

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

<figure><img src="image_134.png" /></figure>
# Appendix G

## Experiments on the degradation of samples of chia mucilage mixed with sand

The aim of this study was to assess the extent to which the mucilage samples mixed with sand had undergone mucilage degradation over the time scale of the water retention measurements (approximately 6 weeks).

### G.1 Materials and methods:

We tested 3 concentration levels of mucilage (dry matter) mixed with Fontainebleau sand (200-250 µm grain size sieved and acid washed): 0% (control), 0.1% and 1%. We incubated samples of mucilage mixed with sand in metal rings placed on a suction table (100% saturation, 3 replications) for 44 days in the dark at laboratory temperature (25°C ±2). Measurements were taken at 4 key dates on the water equilibrium time scale: 0 days, 12 days (saturation on suction table), 29 days (middle of retention measurement experiment) and 44 days (end of retention measurement experiment). A total of 36 samples were studied. Two types of measurement were carried out.

#### G.1.1 Microbial degradation

To estimate microbial degradation of the mucilage we carried out micro-respiration measurements based on the MicroRespTM tool [56]. To summarize: aliquots of sand + mucilage (of the order of 0.4 g precisely weighed) were distributed in the deep wells of a plate (3 to 4 analytical replications). An initial reading of the optical density at 570 nm of the detection plate (containing a gel with a coloured indicator whose colour change after incubation correlates with the release of CO₂ from the sample) is taken. The MicroRespTM system was set up (the detection plate facing the deep-well plate containing the samples, the 2 separated by a silicone seal individualising the wells, the whole being placed in a press to hermetically seal the system, Figure G.1) and incubated in the dark at 23°C for 6 hours. The optical density of the detection plate was then measured a second time and the respiration rate of the sample was calculated to obtain µgC – CO₂g⁻¹ sand-mucilage h⁻¹ [52]. This respiration is indicative of the degradation activity of the organic matter (in this case the mucilage introduced into the sand)
by the microorganisms present in the sample, and is also indicative of the active microbial biomass of the sample [9].

<figure><img src="image_135.png" /><figcaption>Figure G.1: Principle of the MicroRespTM colorimetric measurement to assess soil microbial respiration. In our case, we studied the respiration of sand samples mixed with chia mucilage.</figcaption></figure>

### G.1.2 Presence of ExoPolySaccharides (EPS) adhering to the sand

To estimate the amount of EPS adhering to the sand and therefore likely to influence the water retention of the sample, we carried out total sugar measurements of the same samples after water extraction. According to Yang et al (2012), water extraction of the sample combined with centrifugation separates degraded EPS (supernatant after centrifugation) from non-degraded EPS (pellet after centrifugation). To do this, we introduced 0.03 g of sample and 20 µL of MilliQ water into an Eppendorf tube. After vortexing for 3 minutes, the sample was centrifuged (6869 g) for 3 minutes. The supernatant was removed and we analysed the total sugars in the pellet using the method proposed by Marasco et al [226]. The pellet was mixed with phenol and sulphuric acid and the coloration was measured by spectrophotometry (2 analytical replications, [113, 209]). Results are obtained in µg (glucose equivalent) L⁻¹.

Statistical analyses were performed using XLStat software. The data did not meet the ANOVA criteria, non-parametric tests (Kruskal Wallis, Dunn Bonferroni test) were then used for data analysis.

## G.2 Results and discussion

Microbial respiration was significantly higher in the sand-mucilage samples, with a systematic gradient according to the concentration of mucilage introduced into the sand (p < 0.001, Figure G.2), which suggests a stimulation of the microbial communities present in the samples by the introduced chia mucilage.

We observed a peak in respiration from the first sampling date at 12 days of incubation, followed by a drop in respiration at 29 days of incubation, which remained stable at 44 days of incubation (Figure G.2).
<figure><img src="image_136.png" /><figcaption>Figure G.2: Microbial respiration (CO₂ released) of sand + chia seed mucilage samples incubated on a suction table for 44 days. 0% = control sand without mucilage, 0.1% = sand with 0.1% mucilage, 1% = sand with 1% mucilage. Letters indicate significant differences between measurements (Kruskal Wallis, Dunn Bonferroni test, *p* < 0.05).</figcaption></figure>

As the suction table is not sterile, nor is the chia seed mucilage used, despite washing the sand with acid before handling, we observed significant microbial activity in the sand (including without mucilage, *p* = 0.0004) from 12 days of incubation on the suction table. These results are consistent with the study by Yang *et al.* [362] (carried out, however, with sand collected *in situ* and not washed, but with a lower moisture content), which showed microbial development from the first day of incubation, as well as those of Ahmed *et al.* [4, 1], who observed rapid degradation (a few days) of maize mucilage in a soil and microbial development rates linked to the initial quantities of mucilage added, as well as the study by van Veelen *et al.* [344], who observed the degradation of a chia mucilage in an artificial soil within 2 weeks.

With regard to the quantity of EPS adhering to the soil (total sugars in glucose equivalent), we observed the same pattern as for microbial activity measurements: a difference in EPS quantities depending on the amount of mucilage introduced into the sand and a peak at 12 days (Figure G.3).

Microbial degradation of the mucilage in the sand did not, therefore, reduce these measured quantities of EPS adhering to the sand, on the contrary. According to the work of Ahmed *et al.* [4, 1], the quantity of mucilage integrated into the microbial biomass after 15 days of incubation in soil at optimum moisture (which is the case on the suction table, > 66% saturation) is of the order of 14% and can reach 40% for certain microbial groups. Our hypothesis is that the microbial biomass, which was probably low at the start of incubation (washed sand), was strongly stimulated by the addition of mucilage (Figure G.1), in an environment that was initially very low in organic matter, and that the microorganisms not only degraded this mucilage, but also incorporated it into their biomass and
<figure><img src="image_137.png" /></figure>

Figure G.3: Quantity of EPS adhering to the sand in samples of sand + chia seed mucilage incubated on a suction table for 44 days. 0% = control sand without mucilage, 0.1% = sand with 0.1% mucilage, 1% = sand with 1% mucilage. The EPS were measured using colorimetric analysis of total sugars, directly on the samples (after pre-extraction of the dissolved EPS with water). Letters indicate significant differences between measurements (Kruskal Wallis, Dunn Bonferroni test, *p* < 0.05).

their own EPS production. Our EPS measurements were carried out directly on the sand-mucilage samples (after extraction with water), but without extraction of the EPS adhering to the soil and not incorporated into the microbial biomass (additional step with an ion exchange resin, [53]). The EPS peak observed at 12 days probably corresponds to EPS produced and partially integrated into the biomass and microbial EPS.

This hypothesis would need to be confirmed with microbial biomass measurements combined with isotopic labelling methods, in order to identify the origin of the EPS incorporated into microbial biomasses [4, 1]. In addition, this reality of the rhizospheric environment evolving over time (degradation of mucilages brought by roots into the soil, microbial growth and production of other EPS molecules by microorganisms stimulated by these mucilages) should be taken into account in complementary studies to address the evolution of the physical and hydrophysical properties of the rhizosphere, likely to evolve as a result [244, 344].
# Appendix H

## Coupled system: Drop infiltration in porous medium with amphiphilic matter

### H.1 Dissipative free energy

We consider the system: water drop + porous media.

This system is isolated, if we neglect evaporation and impose zero mass and energy flux through the boundaries of the system.

Therefore it is possible to define a free energy $\mathcal{F}$ that decreases with time.

The system has the following form:

$$ \begin{cases} \frac{\partial \zeta}{\partial t} = \nabla \cdot \left[ Q(\zeta) \nabla \frac{\delta \mathcal{F}}{\delta \zeta} \right] - v \\ \frac{\partial \xi}{\partial t} = \nabla \cdot \left[ K(\xi) \nabla \frac{\delta \mathcal{F}}{\delta \xi} \right] + v \end{cases} \qquad (H.1) $$

with

$$ v = \beta \left( \frac{\delta \mathcal{F}}{\delta \zeta} - \frac{\delta \mathcal{F}}{\delta \xi} \right). \qquad (H.2) $$

We check the sign of $\frac{d\mathcal{F}}{dt}$.
$$
\begin{align*}
\frac{d\mathcal{F}}{dt}(\zeta, \xi) &= \left\langle \frac{\delta\mathcal{F}}{\delta\zeta}, \frac{d\zeta}{dt} \right\rangle + \left\langle \frac{\delta\mathcal{F}}{\delta\xi}, \frac{d\xi}{dt} \right\rangle \\
&= \int_D \frac{\delta\mathcal{F}}{\delta\zeta} \frac{d\zeta}{dt} d\Omega + \int_D \frac{\delta\mathcal{F}}{\delta\xi} \frac{d\xi}{dt} d\Omega \\
&= \int_D \frac{\delta\mathcal{F}}{\delta\zeta} \left( \nabla \cdot \left[ Q(\zeta) \nabla \frac{\delta\mathcal{F}}{\delta\zeta} \right] - v \right) d\Omega + \int_D \frac{\delta\mathcal{F}}{\delta\xi} \left( \nabla \cdot \left[ K(\xi) \nabla \frac{\delta\mathcal{F}}{\delta\xi} \right] + v \right) d\Omega \\
&= - \int_D \left[ Q(\zeta) \left|\nabla \frac{\delta\mathcal{F}}{\delta\zeta}\right|^2 + K(\xi) \left|\nabla \frac{\delta\mathcal{F}}{\delta\xi}\right|^2 \right] d\Omega - \int_D v \left( \frac{\delta\mathcal{F}}{\delta\zeta} - \frac{\delta\mathcal{F}}{\delta\xi} \right) d\Omega \\
&\quad + \int_{\partial D} \frac{\delta\mathcal{F}}{\delta\zeta} Q(\zeta) \underbrace{\partial_n \frac{\delta\mathcal{F}}{\delta\zeta}}_{=0} d\sigma + \int_{\partial D} \frac{\delta\mathcal{F}}{\delta\xi} K(\xi) \underbrace{\partial_n \frac{\delta\mathcal{F}}{\delta\xi}}_{=0} d\sigma \\
&= - \underbrace{\int_D Q(\zeta) \left|\nabla \frac{\delta\mathcal{F}}{\delta\zeta}\right|^2 d\Omega}_{\ge 0} - \underbrace{\int_D K(\xi) \left|\nabla \frac{\delta\mathcal{F}}{\delta\xi}\right|^2 d\Omega}_{\ge 0} - \underbrace{\int_D \beta \left( \frac{\delta\mathcal{F}}{\delta\zeta} - \frac{\delta\mathcal{F}}{\delta\xi} \right)^2 d\Omega}_{\ge 0} \\
&\le 0
\end{align*}
\tag{H.3}
$$

Therefore, the system minimizes its free energy with respect to the time.

Note: if $v \ge 0$ and $\frac{\delta\mathcal{F}}{\delta\zeta} - \frac{\delta\mathcal{F}}{\delta\xi} \ge 0$ then the water goes from the drop to the porous medium.

## H.2 Model parameters

<table><caption>Table H.1: Model parameters used for the simulations of the film and porous substrate model.</caption><thead><tr><th>Parameter name</th><th>Symbol</th><th>Value</th></tr></thead><tbody><tr><td>initial drop volume</td><td>V</td><td>100 &mu;L</td></tr><tr><td>liquid-gas interface energy</td><td>&gamma;</td><td>70 mN.m<sup>-1</sup></td></tr><tr><td>water density</td><td>&rho;</td><td>10<sup>3</sup> kg.m<sup>-3</sup></td></tr><tr><td>gravitational acceleration</td><td>g</td><td>9.81 m.s<sup>-2</sup></td></tr><tr><td>viscosity (T=20&deg;C)</td><td>&mu;</td><td>10<sup>-3</sup> Pa.s</td></tr><tr><td>precursor layer height</td><td>&zeta;<sub>p</sub></td><td>20 &mu;m</td></tr><tr><td>soil thickness</td><td>D</td><td>5 mm</td></tr><tr><td>soil porosity</td><td>&phi;</td><td>0.43</td></tr><tr><td>soil compression</td><td>&kappa;</td><td>200</td></tr><tr><td>macroscopic surface tension</td><td>&Gamma;</td><td>1.7 mN</td></tr><tr><td>capillary height</td><td>H<sub>cap</sub></td><td>26 cm</td></tr><tr><td>saturated hydraulic conductivity</td><td>K<sub>sat</sub></td><td>10<sup>-4</sup> m.s<sup>-1</sup></td></tr><tr><td>pore size distribution coefficient (Van Genuchten model)</td><td>n</td><td>3.31</td></tr><tr><td>pore size distribution coefficient (Brooks and Corey model)</td><td>&lambda;</td><td>0.972</td></tr><tr><td>imbibition rate coefficient</td><td>&beta;</td><td>4.4 &times; 10<sup>-8</sup> m.Pa<sup>-1</sup>.s<sup>-1</sup></td></tr><tr><td>simulation domain length</td><td>L</td><td>40 mm</td></tr><tr><td>repellent saturation</td><td>s<sub>r</sub></td><td>0.05</td></tr><tr><td>repellent range coefficient</td><td>&sigma;</td><td>0.2</td></tr><tr><td>adaption substrate wetting</td><td>c<sub>0</sub></td><td>0.4</td></tr><tr><td>amplitude factor amplitude</td><td>&zeta;</td><td>1.418 m</td></tr><tr><td>amplitude factor exponent</td><td>&upsilon;</td><td>1.402</td></tr></tbody></table>
## H.3 Stability analysis

The dimensional system (5.30) can formally be written as

$$ \partial_t U = F(U) = \begin{bmatrix} F_{\zeta}(U) \\ F_{\xi}(U) \end{bmatrix}, \qquad (H.4) $$

where

$$ U = \begin{bmatrix} \zeta \\ \xi \end{bmatrix} $$

The linear stability of the uniform steady film $\xi_0$ and the uniform height saturation $\xi_0$ is determined by the spectral analysis of the Jacobian matrix $\mathbf{J}$ of $F$

$$ \mathbf{J} = \frac{DF}{D(\zeta, \xi)}(U_0) = \begin{bmatrix} \frac{DF_{\zeta}}{D\zeta} & \frac{DF_{\zeta}}{D\xi} \\ \frac{DF_{\xi}}{D\xi} & \frac{DF_{\xi}}{D\xi} \end{bmatrix} \qquad (H.5) $$

where $U_0 = (\zeta_0, \xi_0)$ and the index 0 represents the function at the equilibrium $(\zeta_0, \xi_0)$.

Since the problem is invariant w.r.t translation, the linear system can be solved for each Fourier mode $k$, using the ansatz

$$ U_k = U_0 + \epsilon u_j \exp(\sigma t + ikx). \qquad (H.6) $$

We define $\mathbf{J}_k$ a matrix 2 × 2 in the Fourier mode $k$, as in [7],

$$ \mathbf{J}_k = k^2 \mathbf{Q} \cdot \mathbf{E} - \beta \mathbf{V}, \qquad (H.7) $$

where $\mathbf{Q}$ is the mobility matrix

$$ \mathbf{Q} = \begin{bmatrix} Q_0 & 0 \\ 0 & DK_0 \end{bmatrix}, \qquad (H.8) $$

$\mathbf{V}$ is the density flux exchange matrix

$$ \mathbf{V} = \begin{bmatrix} \partial_{\zeta\zeta} f_w - \partial_{\zeta\xi} f_w + \gamma k^2 & \partial_{\zeta\xi} f_w - \partial_{\xi\xi} f_w - Tk^2 \\ -\partial_{\zeta\zeta} f_w + \partial_{\zeta\xi} f_w - \gamma k^2 & -\partial_{\zeta\xi} f_w + \partial_{\xi\xi} f_w + Tk^2 \end{bmatrix}, \qquad (H.9) $$

$\mathbf{E}$ is the energy matrix in Fourier space

$$ \mathbf{E} = -\mathbf{W} - k^2 \mathbf{T}, \qquad (H.10) $$

with the surface tension matrix

$$ \mathbf{T} = \begin{bmatrix} \gamma & 0 \\ 0 & T \end{bmatrix}, \qquad (H.11) $$

and the wetting energy matrix

$$ \mathbf{W} = \begin{bmatrix} \partial_{\zeta\zeta} f_w & \partial_{\zeta\xi} f_w \\ \partial_{\zeta\xi} f_w & \partial_{\xi\xi} f_w \end{bmatrix}. \qquad (H.12) $$

Eq.(H.7) is used to determine the eigenvalues (growth rates $\sigma$) and the eigenvectors (amplitude
ratio $u$).

### H.3.1 Uniform variation

To study the uniform variation of film and height of saturation, we focus on mode $k = 0$. Eq. (H.7) becomes

$$ \mathbf{J}_0 = -\beta \begin{bmatrix} \partial_{\zeta\zeta} f_w - \partial_{\zeta\xi} f_w & \partial_{\zeta\xi} f_w - \partial_{\xi\xi} f_w \\ -\partial_{\zeta\zeta} f_w + \partial_{\zeta\xi} f_w & -\partial_{\zeta\xi} f_w + \partial_{\xi\xi} f_w \end{bmatrix} \quad (H.13) $$

The computation of its eigenvalues $\sigma$ and the corresponding eigenvectors $u$ gives

$$ \sigma_1 = 0 \quad \text{and} \quad u_1 = \left(-\frac{\partial_{\zeta\xi} f_w - \partial_{\xi\xi} f_w}{\partial_{\zeta\zeta} f_w - \partial_{\zeta\xi} f_w}, 1\right), \quad (H.14) $$

$$ \sigma_2 = -\beta (\partial_{\zeta\zeta} f_w - 2\partial_{\zeta\xi} f_w + \partial_{\xi\xi} f_w) \quad \text{and} \quad u_2 = (1, -1). $$

The stability condition of the varicose mode is given by

$$ \partial_{\zeta\zeta} f_w - 2\partial_{\zeta\xi} f_w + \partial_{\xi\xi} f_w > 0 \quad (H.15) $$

### H.3.2 Equilibrium stability

We consider $U_0$ the equilibrium satisfying the system of Eq. (5.35). Under the criteria of water height conservation Eq. (5.33), the couple $(\zeta_0, \xi_0)$ is an extremum or an inflection point of the free energy for flat film and uniform height of saturation.

#### H.3.2.1 Nature of critical point

To define the stability of equilibrium couple $U_0$, we study the nature of critical point using hessian matrix $\mathbf{H}$,

$$ \mathbf{H} = \begin{bmatrix} \partial_{\zeta\zeta} f_w & \partial_{\zeta\xi} f_w \\ \partial_{\zeta\xi} f_w & \partial_{\xi\xi} f_w \end{bmatrix}. \quad (H.16) $$

#### H.3.2.1.1 Minimum
The couple $U_0$ is a minimum, i.e the determinant of the hessian matrix and its trace are positive,

$$ \begin{gathered} \partial_{\zeta\zeta} f_w \partial_{\xi\xi} f_w > (\partial_{\zeta\xi} f_w)^2, \\ \partial_{\zeta\zeta} f_w + \partial_{\xi\xi} f_w > 0. \end{gathered} \quad (H.17) $$

As $\partial_{\xi} a(\xi, c) < 0$ and $\partial_{\zeta} f_w^{\zeta}(\zeta) > 0$ then $\partial_{\zeta\xi} f_w < 0$.

The equilibrium couple satisfy the stability condition Eq. (H.15), i.e $U_0$ is a stable equilibrium couple.

#### H.3.2.1.2 Maximum
The couple $U_0$ is a maximum, i.e the determinant of the hessian matrix is positive and its trace is negative,

$$ \begin{gathered} \partial_{\zeta\zeta} f_w \partial_{\xi\xi} f_w > (\partial_{\zeta\xi} f_w)^2, \\ \partial_{\zeta\zeta} f_w + \partial_{\xi\xi} f_w < 0. \end{gathered} \quad (H.18) $$
The stability of a critical couple corresponding to a maximum is not as simple as in the case of a minimum.

First of all, using Eq. (H.18) and the properties of the Hessian matrix **H**, we deduce that the derivatives $\partial_{\zeta\zeta}f_w$ and $\partial_{\xi\xi}f_w$ are negative at the equilibrium point $U_0$,

As $\partial_{\zeta\xi}f_w < 0$ and using the positive determinant of **H**, we have

$$0 > \partial_{\zeta\xi}f_w > -\sqrt{\partial_{\zeta\zeta}f_w \partial_{\xi\xi}f_w} > -\frac{|\partial_{\zeta\zeta}f_w| + |\partial_{\xi\xi}f_w|}{2} \quad (H.19)$$

**Remark H.3.1.** We can express the derivative of $f_w$, $\partial_{\zeta\zeta}f_w$ and $\partial_{\xi\xi}f_w$, as

$$\partial_{\zeta\zeta}f_w = -|\partial_{\zeta\zeta}f_w|$$

$$\partial_{\xi\xi}f_w = -|\partial_{\xi\xi}f_w|$$

Then the equilibrium couple don't satisfy the stability condition of Eq. (H.15), i.e $U_0$ is an unstable equilibrium couple.

**H.3.2.1.3 Inflection point** The couple $U_0$ is an inflection point, i.e the determinant of the hessian matrix is negative,

$$\partial_{\zeta\zeta}f_w \partial_{\xi\xi}f_w < (\partial_{\zeta\xi}f_w)^2, \quad (H.20)$$

In this situation, we can conclude that the couple is

* stable if the hessian trace is positive

* unstable if the hessian trace is negative and $\partial_{\zeta\zeta}f_w \partial_{\xi\xi}f_w \le 0$

* Else we can not conclude the stability without studying in more detail functions.

**H.3.2.2 Free energy density variation**

To define the stability of equilibrium couple $U_0$, we study the variation of the free energy density under the conservation criteria.

Let us consider a total water in the system **H** and introduce the function $\zeta(\xi) = H - \xi$ that represent the height of water in the flat film $\zeta$ depending on the height of water in the porous media $\xi$. We define a function $\varphi(\xi)$ such that it represent the free energy density under the conservation criteria

$$\varphi(\xi) = f(H - \xi, \xi).$$

The second derivative of $\varphi$ is exactly the function required to study the stability of the varicose mode,

$$\begin{aligned} \frac{d^2\varphi}{d\xi^2}(\xi) &= \frac{d^2}{d\xi^2} \left[ f_w^\zeta(\zeta(\xi)) + f_w^\xi(\xi) + f_w^{\text{int}}(\zeta(\xi), \xi) \right] \\ &= \frac{d^2 f_w^h}{d\zeta^2}(\zeta(\xi)) + \frac{d^2 f_w^\xi}{d\xi^2}(\xi) - 2 \frac{\partial^2 f_w^{\text{int}}}{\partial\zeta\partial\xi}(\zeta(\xi), \xi) + \frac{\partial^2 f_w^{\text{int}}}{\partial\zeta^2}(\zeta(\xi), \xi) + \frac{\partial^2 f_w^{\text{int}}}{\partial\xi^2}(\zeta(\xi), \xi) \\ &= \frac{\partial^2 f_w}{\partial\zeta^2}(\zeta(\xi), \xi) + \frac{\partial^2 f_w}{\partial\xi^2}(\zeta(\xi), \xi) - 2 \frac{\partial^2 f_w}{\partial\zeta\partial\xi}(\zeta(\xi), \xi). \end{aligned}$$
The stability condition of Eq.(H.15) can be rewrite as

$$ \frac{d^2\varphi}{d\xi^2}(\xi) > 0. \qquad (H.21) $$

Note : if the couple corresponds to an extremum of the free energy function, than this criteria correspond to an extremum in a convex area. In other word, the stability condition is verified only for minima of the free energy density function.

## H.4 Extrema of free energy density

For uniform saturation and flat film, we express the free energy density

$$ f(\zeta, \xi, c) = f_w^{\zeta}(\zeta) + f_w^{\xi}(\xi, c) + f_w^{\text{int}}(\zeta, \xi, c) + \Lambda(H - \xi - \zeta), \qquad (H.22) $$

where $\Lambda$ is the Lagrange multiplier.

Thus, the derivative of $f$ with respect to $\zeta$ and $\xi$ are

$$ \frac{\delta F_1}{\delta \zeta} = \partial_{\zeta} f_w^{\zeta}(\zeta) + \partial_{\zeta} f_w^{\text{int}}(\zeta, \xi, c) - \Lambda = \frac{\delta F}{\delta \zeta} - \Lambda \qquad (H.23) $$

and

$$ \frac{\delta F_1}{\delta \xi} = \partial_{\xi} f_w^{\xi}(\xi, c) + \partial_{\xi} f_w^{\text{int}}(\zeta, \xi, c) - \Lambda = \frac{\delta F}{\delta \xi} - \Lambda \qquad (H.24) $$

where $\frac{\delta F}{\delta \zeta}$ and $\frac{\delta F}{\delta \xi}$ are the functional free energy density without respect the mass conservative criteria.

We consider a steady-states $U_0 = (\zeta_0, \xi_0)$ of the system (5.30) is define by Eq.(5.35). In other words, $U_0$ is an extrema of the free energy density Eq.(5.34),

$$ \left\{ \begin{array}{l} \displaystyle \frac{\delta F_1}{\delta \zeta}(U_0) = 0, \\ \\ \displaystyle \frac{\delta F_1}{\delta \xi}(U_0) = 0, \\ \\ \displaystyle \frac{\delta F_1}{\delta \Lambda}(U_0) = 0, \end{array} \right. \qquad (H.25) $$

We consider only film height depending on $\xi$ such that $\zeta(\xi) = H - \xi$, i.e the free energy density can express without the conservative condition as a function $\varphi$ depending only of the height of water in the soil,

$$ \varphi(\xi) = f(\zeta(\xi), \xi, c) = f_w^{\zeta}(\zeta(\xi)) + f_w^{\xi}(\xi, c) + f_w^{\text{int}}(\zeta(\xi), \xi, c). \qquad (H.26) $$

Then equilibria of the system correspond to extrema of $\varphi$.
# Appendix I

## Water transfer with single root

### I.1 Finite Element Method

In this study, we assume a symmetry with respect to z axis. In other words, the angular part is invariant. Moreover, when we fix a depth z the studying system can be express in dimensional or dimensionless formulation as in the following section.

#### I.1.1 Dimensional formulation

We consider the water transfer in the soil close to a single root, i.e the water transfer is modelled using

$$ 
\left\{ 
\begin{aligned}
\phi \partial_t s &= \frac{1}{r} \partial_r \left[ r K(s) \partial_r \frac{\delta \mathcal{F}}{\delta s} \right] + S, \\
\frac{\delta \mathcal{F}}{\delta s} &= - \rho g H_{cap} \Pi_w(s, c) - \Gamma \Delta s,
\end{aligned}
\right. 
\qquad (I.1) $$

where S is a sink term representing the water uptake of the root. This term is classically a boundary condition at the interface root/soil.

Thus, this model is similar to the system Eqs (2.40). Using FEM Eqs (2.40) is express

$$ 
\left\{ 
\begin{aligned}
\phi \int_{\Omega} \partial_t svd\Omega &= \int_{\Omega} \nabla \cdot \left[ K(s) \nabla \frac{\delta \mathcal{F}}{\delta s} \right] v d\Omega, \\
\int_{\Omega} \frac{\delta \mathcal{F}}{\delta s} w d\Omega &= - \int_{\Omega} \rho g H_{cap} \Pi_w(s, c) w d\Omega - \int_{\Omega} \Gamma \Delta s w d\Omega.
\end{aligned}
\right. 
\qquad (I.2) $$

This is simplified using Green formula,

$$ 
\left\{ 
\begin{aligned}
\phi \int_{\Omega} \partial_t svd\Omega &= \int_{\partial\Omega} K(s) \partial_n \frac{\delta \mathcal{F}}{\delta s} v d\sigma - \int_{\Omega} K(s) \nabla \frac{\delta \mathcal{F}}{\delta s} \cdot \nabla v d\Omega, \\
\int_{\Omega} \frac{\delta \mathcal{F}}{\delta s} w d\Omega &= - \int_{\Omega} \rho g H_{cap} \Pi_w(s, c) w d\Omega - \int_{\partial\Omega} \Gamma \partial_n s w d\sigma + \int_{\Omega} \Gamma \nabla s \nabla w d\Omega.
\end{aligned}
\right. 
\qquad (I.3) $$
In cylindrical coordinate, this system becomes

$$ egin{cases} 
\phi \int_{\Omega_r} \partial_t svr dr = S - \int_{\Omega_r} K(s) \partial_r \frac{\delta \mathcal{F}}{\delta s} \partial_r vr dr, \\ 
\int_{\Omega_r} \frac{\delta \mathcal{F}}{\delta s} wr dr = - \int_{\Omega_r} \rho g H_{cap} \Pi_w(s, c) wr dr - \int_{\partial \Omega_r} \Gamma \partial_n swrd \sigma_r + \int_{\Omega_r} \Gamma \partial_r s \partial_r wr dr, 
\end{cases} 
\quad (1.4) $$

where the sink term $S = \int_{\partial \Omega_r} K(s) \partial_n \frac{\delta \mathcal{F}}{\delta s} vrd\sigma_r$.

### I.1.2 Dimensionless formulation

The dimensionless expression of the system Eqs (I.1) is given by

$$ egin{cases} 
\frac{\phi}{T^*} \partial_{\tilde{t}} s = \frac{1}{\tilde{r}} \partial_{\tilde{r}} \left[ \tilde{r} k_r(s) \partial_{\tilde{r}} \frac{\delta \tilde{\mathcal{F}}}{\delta s} \right] + S, \\ 
\frac{\delta \tilde{\mathcal{F}}}{\delta s} = - \Pi_w(s, c) - \frac{T}{\tilde{r}} \partial_{\tilde{r}} [\tilde{r} \partial_{\tilde{r}} s], 
\end{cases} 
\quad (1.5) $$

where $T = \frac{\Gamma}{\rho g H_{cap} L^{*2}}$ and $T^* = \frac{\phi \mu L^{*2}}{\rho g H_{cap} k} = \frac{\phi L^{*2}}{H_{cap} K_{sat}}$.

Thus, this model is similar to the dimensionless system Eqs (2.40). Using the FEM of this dimensionless system, we have

$$ egin{cases} 
\int_{\Omega} \partial_t svd\Omega = \int_{\partial\Omega} k_r(s) \partial_n \frac{\delta \mathcal{F}}{\delta s} v d\sigma - \int_{\Omega} k_r(s) \nabla \frac{\delta \mathcal{F}}{\delta s} \cdot \nabla v d\Omega, \\ 
\int_{\Omega} \frac{\delta \mathcal{F}}{\delta s} w d\Omega = - \int_{\partial\Omega} \Pi_w(s, c) w d\Omega - \int_{\partial\Omega} T \partial_n sw d\sigma + \int_{\Omega} T \nabla s \nabla w d\Omega. 
\end{cases} 
\quad (1.6) $$

In cylindrical coordinate, this system becomes

$$ egin{cases} 
\int_{\Omega_r} \partial_t svr dr = S - \int_{\Omega_r} k_r(s) \partial_r \frac{\delta \mathcal{F}}{\delta s} \partial_r vr dr, \\ 
\int_{\Omega_r} \frac{\delta \mathcal{F}}{\delta s} wr dr = - \int_{\Omega_r} \Pi_w(s, c) wr dr - \int_{\partial \Omega_r} T \partial_n swrd \sigma_r + \int_{\Omega_r} T \partial_r s \partial_r wr dr, 
\end{cases} 
\quad (1.7) $$

where the sink term $S = \int_{\partial \Omega_r} k_r(s) \partial_n \frac{\delta \mathcal{F}}{\delta s} vrd\sigma_r$.

## I.2 Simulated parameters
<table><caption>Table I.1: Model parameters used for the simulations of transfer applied to the rhizosphere.</caption><thead><tr><th>Parameter name</th><th>Symbol</th><th>Value</th></tr></thead><tbody><tr><td>maximum EPS concentration at root surface</td><td>c<sub>max</sub></td><td>5%</td></tr><tr><td>minimum EPS concentration</td><td>c<sub>crit</sub></td><td>10<sup>-3</sup>%</td></tr><tr><td>dimensional radius of the root</td><td>r<sub>root</sub></td><td>0.05 cm</td></tr><tr><td>dimensional radius of the rhizospheric soil</td><td>r<sub>rhizo</sub></td><td>0.2 cm</td></tr><tr><td>dimensional radius of the bulk soil</td><td>r<sub>bulk</sub></td><td>1 cm</td></tr><tr><td>water density</td><td>ρ</td><td>10<sup>3</sup> kg.m<sup>-3</sup></td></tr><tr><td>gravitational acceleration</td><td>g</td><td>9.81 m.s<sup>-2</sup></td></tr><tr><td>viscosity (T=20°C)</td><td>μ</td><td>10<sup>-3</sup> Pa.s</td></tr><tr><td>intrinsic permeability</td><td>k</td><td>10<sup>-11</sup> m<sup>2</sup></td></tr><tr><td>soil porosity</td><td>φ</td><td>0.43</td></tr><tr><td>soil compression</td><td>κ</td><td>200</td></tr><tr><td>sort of macroscopic surface tension</td><td>Γ</td><td>1.7 mN</td></tr><tr><td>capillary height</td><td>H<sub>cap</sub></td><td>26 cm</td></tr><tr><td>saturated hydraulic conductivity</td><td>K<sub>sat</sub></td><td>10<sup>-4</sup> m.s<sup>-1</sup></td></tr><tr><td>pore size distribution coefficient (Van Genuchten model)</td><td>n</td><td>3.31</td></tr><tr><td>pore size distribution coefficient (Brooks and Corey model)</td><td>λ</td><td>0.972</td></tr><tr><td>simulation domain length</td><td>L</td><td>0.95 cm</td></tr><tr><td>repellent saturation</td><td>s<sub>r</sub></td><td>0.05</td></tr><tr><td>repellent range coefficient</td><td>σ</td><td>0.2</td></tr><tr><td>amplitude factor amplitude</td><td>ζ</td><td>1.418 m</td></tr><tr><td>amplitude factor exponent</td><td>ν</td><td>1.402</td></tr><tr><td>dimensionless length</td><td>L<sup>*</sup></td><td>0.2 cm</td></tr><tr><td>dimensionless time</td><td>T<sup>*</sup></td><td>66.2 ms</td></tr><tr><td>dimensionless radius of the root</td><td>r̃<sub>root</sub></td><td>0.25</td></tr><tr><td>dimensionless radius of the rhizospheric soil</td><td>r̃<sub>rhizo</sub></td><td>1</td></tr><tr><td>dimensionless radius of the bulk soil</td><td>r̃<sub>bulk</sub></td><td>5</td></tr></tbody></table>
# Appendix J

## Publications and Seminars

### J.1 Publications in peer-reviewed journals

* Philippe BELTRAME, Florian CAJOT. Model of hydrophobic porous media applied to stratified media: Water trapping, intermittent flow and fingering instability. 06/2022, Europhysics Letters, Volume 138, Number 5, doi:10.1209/0295-5075/ac71c1.

* Le Gall S., Lapie C., Cajot F., Doussan C., Corridor L., Bérard A.Chemical diversity of crop plants root mucilages. Implications for their maximal water content and mineralization in soils. Accepted 19/01/2024, Rhizosphere.

### J.2 Seminars and conferences

* Florian CAJOT, Philippe BELTRAME, Annette BERARD, Claude DOUSSAN. Model of water drop infiltration in porous media with amphiphilic matter, 10/2023, Systèmes multiphases - multicomposants : Présentation oral, 16e journées d'études des milieux poreux : JEMP 2023

* Florian CAJOT, Philippe BELTRAME, Annette BERARD, Claude DOUSSAN. Model of water drop infiltration in porous media with amphiphilic matter, 05/2023, MS06-A-Physics of multi-phase flow in diverse porous media : Présentation oral, 15th Annual International Conference on Porous Media : InterPore2023

* Florian CAJOT. The role of Exopolysaccharides in water transfer within the rhizosphere, 07/2022, Présentation oral, Symposium Tersys-Implanteus, 1st Summer School Avignon

* Florian CAJOT, Philippe BELTRAME. Intermittent water transfer: A wettability model for porous media., 08/2021, FS05 - Porous Media : Poster, 25th International Congress of Theoretical and Applied Mechanics Milan (virtuel)
Bibliography

[1] M. A. Ahmed, C. C. Banfield, M. Sanaullah, A. Gunina, and M. A. Dippold. Utilisation of mucilage C by microbial communities under drought. *Biology and Fertility of Soils*, 54(1):83-94, Jan. 2018.

[2] M. A. Ahmed, M. Holz, S. K. Woche, J. Bachmann, and A. Carminati. Effect of soil drying on mucilage exudation and its water repellency: a new method to collect mucilage. *Journal of Plant Nutrition and Soil Science*, 178(6):821-824, 2015. _eprint: https://onlinelibrary.wiley.com/doi/pdf/10.1002/jpln.201500177.

[3] M. A. Ahmed, E. Kroener, M. Holz, M. Zarebanadkouki, and A. Carminati. Mucilage exudation facilitates root water uptake in dry soils. *Functional Plant Biology*, 41(11):1129, 2014.

[4] M. A. Ahmed, M. Sanaullah, E. Blagodatskaya, K. Mason-Jones, H. Jawad, Y. Kuzyakov, and M. A. Dippold. Soil microorganisms exhibit enzymatic and priming response to root mucilage under drought. *Soil Biology and Biochemistry*, 116:410-418, Jan. 2018.

[5] A. A. Albalasmeh and T. A. Ghezzehei. Interplay between soil drying and root exudation in rhizosheath development. *Plant and Soil*, 374(1):739-751, Jan. 2014.

[6] V. Allègre, L. Jouniaux, F. Lehmann, and P. Sailhac. Streaming potential dependence on water-content in Fontainebleau sand. *Geophysical Journal International*, 182(3):1248-1266, Sept. 2010.

[7] G. Amarandei, P. Beltrame, I. Clancy, C. O'Dwyer, A. Arshak, U. Steiner, D. Corcoran, and U. Thiele. Pattern formation induced by an electric field in a polymer-air-polymer thin film system. *Soft Matter*, 8(23):6333, 2012.

[8] S. Amstutz, D. Gourion, and M. Zabiba. VARIATIONAL APPROXIMATION OF INTERFACE ENERGIES AND APPLICATIONS. *Interfaces and Free Boundaries*, (23):59-102, 2021. Publisher: European Mathematical Society.

[9] J. P. E. Anderson and K. H. Domsch. A physiological method for the quantitative measurement of microbial biomass in soils. *Soil Biology and Biochemistry*, 10(3):215-221, Jan. 1978.

[10] J. Ballester, X. Rodó, and F. Giorgi. Future changes in Central Europe heat waves expected to mostly follow summer mean warming. *Climate Dynamics*, 35(7):1191-1205, Dec. 2010.

[11] J. W. Barrett, J. F. Blowey, and H. Garcke. Finite element approximation of a fourth order nonlinear degenerate parabolic equation. *Numerische Mathematik*, 80(4):525-556, Oct. 1998.
[12] J. W. Barrett, J. F. Blowey, and H. Garcke. Finite Element Approximation of the Cahn-Hilliard Equation with Degenerate Mobility. *SIAM Journal on Numerical Analysis*, 37(1):286-318, Jan. 1999. Publisher: Society for Industrial and Applied Mathematics.

[13] P. Barré and P. D. Hallett. Rheological stabilization of wet soils by model root and fungal exudates depends on clay mineralogy. *European Journal of Soil Science*, 60(4):525-538, 2009. _eprint: https://onlinelibrary.wiley.com/doi/pdf/10.1111/j.1365-2389.2009.01151.x.

[14] S. Batany. *Influence d'un macropore sur l'écoulement et le transport de solutés en milieu poreux : expérimentations sur sol modèle macroporé et simulations numériques*. phdthesis, Université Paris-Est, Nov. 2016.

[15] T. W. J. Bauters, T. S. Steenhuis, J.-Y. Parlange, and D. A. DiCarlo. Preferential Flow in Water-Repellent Sands. *Soil Science Society of America Journal*, 62(5):1185-1190, 1998.

[16] A. Beliaev and S. Hassanizadeh. A Theoretical Model of Hysteresis and Dynamic Effects in the Capillary Relation for Two-phase Flow in Porous Media. *Transport in Porous Media*, 43(3):487-510, 2001.

[17] A. Beljadid, L. Cueto-Felgueroso, and R. Juanes. A continuum model of unstable infiltration in porous media endowed with an entropy function. *Advances in Water Resources*, 144:103684, Oct. 2020.

[18] P. Beltrame. Partial and complete wetting in a micro-tube. *Europhysics Letters*, 121(6):64002, May 2018. Publisher: EDP Sciences, IOP Publishing and Società Italiana di Fisica.

[19] P. Beltrame. Soil Wettability. In *Soils as a Key Component of the Critical Zone 3*, pages 53-84. John Wiley & Sons, Ltd, 2018.

[20] P. Beltrame and F. Cajot. Model of hydrophobic porous media applied to stratified media: Water trapping, intermittent flow and fingering instability. *Europhysics Letters*, 138(5):53004, June 2022. Publisher: EDP Sciences, IOP Publishing and Società Italiana di Fisica.

[21] P. Beltrame, E. Knobloch, P. Hänggi, and U. Thiele. Rayleigh and depinning instabilities of forced liquid ridges on heterogeneous substrates. *Physical Review E*, 83(1):016305, Jan. 2011. Publisher: American Physical Society.

[22] P. Beltrame and U. Thiele. Time Integration and Steady-State Continuation for 2d Lubrication Equations. *SIAM Journal on Applied Dynamical Systems*, 9(2):484-518, Jan. 2010. Publisher: Society for Industrial and Applied Mathematics.

[23] P. Benard, M. Zarebanadkouki, M. Brax, R. Kaltenbach, I. Jerjen, F. Marone, E. Couradeau, V. J. Felde, A. Kaestner, and A. Carminati. Microhydrological Niches in Soils: How Mucilage and EPS Alter the Biophysical Properties of the Rhizosphere and Other Biological Hotspots. *Vadose Zone Journal*, 18(1):180211, 2019. _eprint: https://onlinelibrary.wiley.com/doi/pdf/10.2136/vzj2018.12.0211.
[24] P. Benard, M. Zarebanadkouki, C. Hedwig, M. Holz, M. Ahmed, and A. Carminati. Pore-Scale Distribution of Mucilage Affecting Water Repellency in the Rhizosphere. *Vadose Zone Journal*, 17(1):170013, 2018. _eprint: https://onlinelibrary.wiley.com/doi/pdf/10.2136/vzj2017.01.0013.

[25] E. Benizri, C. Nguyen, S. Piutti, S. Slezack-Deschaumes, and L. Philippot. Additions of maize root mucilage to soil changed the structure of the bacterial community. *Soil Biology and Biochemistry*, 39(5):1230-1233, May 2007.

[26] E. Beretta, M. Bertsch, and R. Passo. Nonnegative solutions of a fourth-order nonlinear degenerate parabolic equation. *Archive for Rational Mechanics and Analysis*, 129:175-200, June 1995.

[27] F. Bernis. Finite speed of propagation and continuity of the interface for thin viscous flows. *Advances in Differential Equations*, 1(3):337-368, Jan. 1996. Publisher: Khayyam Publishing, Inc.

[28] F. Bernis and A. Friedman. Higher order nonlinear degenerate parabolic equations. *Journal of Differential Equations*, 83(1):179-206, Jan. 1990.

[29] A. L. Bertozzi and M. Pugh. The lubrication approximation for thin viscous films: Regularity and long-time behavior of weak solutions. *Communications on Pure and Applied Mathematics*, 49(2):85-123, 1996.
_eprint: https://onlinelibrary.wiley.com/doi/pdf/10.1002/%28SICI%291097-0312%28199602%2949%3A2%3C85%3A%3AAID-CPA1%3E3.0.CO%3B2-2.

[30] M. Bertsch, R. D. Passo, H. Garcke, and G. Grün. The thin viscous flow equation in higher space dimensions. *Advances in Differential Equations*, 3(3):417-440, Jan. 1998. Publisher: Khayyam Publishing, Inc.

[31] K. Beven and P. Germann. Macropores and water flow in soils. *Water Resources Research*, 18(5):1311-1325, 1982. _eprint: https://onlinelibrary.wiley.com/doi/pdf/10.1029/WR018i005p01311.

[32] E. B. A. Bisdom, L. W. Dekker, and J. F. T. Schouté. Water repellency of sieve fractions from sandy soils and relationships with organic material and soil structure. In L. Brussaard and M. J. Kooistra, editors, *Soil Structure/Soil Biota Interrelationships*, pages 105-118. Elsevier, Amsterdam, Jan. 1993.

[33] D. Bonn, J. Eggers, J. Indekeu, J. Meunier, and E. Rolley. Wetting and spreading. *Reviews of Modern Physics*, 81(2):739-805, May 2009. Publisher: American Physical Society.

[34] C. Bordes, L. Jouniaux, M. Dietrich, J.-P. Pozzi, and S. Garambois. First laboratory measurements of seismo-magnetic conversions in fluid-filled Fontainebleau sand. *Geophysical Research Letters*, 33(1), 2006. _eprint: https://onlinelibrary.wiley.com/doi/pdf/10.1029/2005GL024582.

[35] E. Bormashenko and V. Starov. Impact of surface forces on wetting of hierarchical surfaces and contact angle hysteresis. *Colloid and Polymer Science*, 291(2):343-346, Feb. 2013.
[36] F. Boyer. *Écoulements diphasiques de type Cahn-Hilliard*. These de doctorat, Bordeaux 1, Jan. 2001.

[37] M. Brax, C. Buchmann, K. Kenngott, G. E. Schaumann, and D. Diehl. Influence of the physico-chemical properties of root mucilage and model substances on the microstructural stability of sand. *Biogeochemistry*, 147(1):35-52, Jan. 2020.

[38] M. Brax, C. Buchmann, and G. E. Schaumann. Biohydrogel induced soil-water interactions: how to untangle the gel effect? A review. *Journal of Plant Nutrition and Soil Science*, 180(2):121-141, 2017. _eprint: https://onlinelibrary.wiley.com/doi/pdf/10.1002/jpln.201600453.

[39] M. Brax, G. E. Schaumann, and D. Diehl. Gel formation mechanism and gel properties controlled by Ca2+ in chia seed mucilage and model substances. *Journal of Plant Nutrition and Soil Science*, 182(1):92-103, 2019. _eprint: https://onlinelibrary.wiley.com/doi/pdf/10.1002/jpln.201800430.

[40] C. J. Brinker and G. W. Scherer. *Sol-gel science: the physics and chemistry of sol-gel processing*. Academic Press, Boston, 1990.

[41] F. Brochard-Wyart, P. Nassoy, and P.-H. Puech. *Physique de la matière molle*. Sciences sup. Dunod, Malakoff, 2018.

[42] R. H. Brooks and A. T. Corey. Properties of Porous Media Affecting Fluid Flow. *Journal of the Irrigation and Drainage Division*, 92(2):61-88, June 1966. Publisher: American Society of Civil Engineers.

[43] Brooks, R.H. and Corey, A.T. HYDRAULIC PROPERTIES OF POROUS MEDIA. *Hydrology Paper*, (3), 1964.

[44] Brooks, R.H.; Corey, A.T. *Hydraulic properties of porous media*. 1964.

[45] P. N. Brown, G. D. Byrne, and A. C. Hindmarsh. VODE: A Variable-Coefficient ODE Solver. *SIAM J. Sci. and Stat. Comput.*, 10(5):1038-1051, Sept. 1989.

[46] H. Brézis. Analyse fonctionnelle: théorie et applications; [cours master, agrégation, licence 4/5]. Mathématiques appliquées pour le master. Dunod, Paris, nouvelle présentation, [nachdr.] edition, 2010.

[47] H. Brézis and H. Brézis. *Functional analysis, Sobolev spaces and partial differential equations*. Universitext. Springer, New York ; London, 2011.

[48] E. Buckingham and United States. Bureau of Soils. *Studies on the movement of soil moisture*. Washington, Govt. Print. Off., 1907.

[49] W. J. Bullied, P. R. Bullock, and R. C. Van Acker. Modelling Soil Water Retention for Weed Seed Germination Sensitivity to Water Potential. *Applied and Environmental Soil Science*, 2012:e812561, May 2012. Publisher: Hindawi.

[50] N. Burdine. Relative Permeability Calculations From Pore Size Distribution Data. *Journal of Petroleum Technology*, 5(03):71-78, Mar. 1953.
[51] K. Burrage, J. C. Butcher, and F. H. Chipman. An implementation of singly-implicit Runge-Kutta methods. *BIT*, 20(3):326-340, Sept. 1980.

[52] A. Bérard, T. Bouchet, G. Sévenier, A. L. Pablo, and R. Gros. Resilience of soil microbial communities impacted by severe drought and high temperature in the context of Mediterranean heat waves. *European Journal of Soil Biology*, 47(6):333-342, Nov. 2011.

[53] A. Bérard, T. Clavel, C. Le Bourvellec, A. Davoine, S. Le Gall, C. Doussan, and S. Bureau. Exopolysaccharides in the rhizosphere: A comparative study of extraction methods. Application to their quantification in Mediterranean soils. *Soil Biology and Biochemistry*, 149:107961, Oct. 2020.

[54] J. W. Cahn and J. E. Hilliard. Free energy of a nonuniform System. 1. Interfacial free energy. *J. Chem. Phys.*, 28:258-267, 1958.

[55] R. Calvet, C. Chenu, and S. Houot. *Les matières organiques des sols: rôles agronomiques et environnementaux*. Agri production. Éd. France agricole, Paris, 2011.

[56] C. D. Campbell, S. J. Chapman, C. M. Cameron, M. S. Davidson, and J. M. Potts. A Rapid Microtiter Plate Method To Measure Carbon Dioxide Evolved from Carbon Substrate Amendments so as To Determine the Physiological Profiles of Soil Microbial Communities by Using Whole Soil. *Applied and Environmental Microbiology*, 69(6):3593-3599, June 2003. Publisher: American Society for Microbiology.

[57] G. S. Campbell. A SIMPLE METHOD FOR DETERMINING UNSATURATED CONDUCTIVITY FROM MOISTURE RETENTION DATA: *Soil Science*, 117(6):311-314, June 1974.

[58] A. Canarini, C. Kaiser, A. Merchant, A. Richter, and W. Wanek. Root Exudation of Primary Metabolites: Mechanisms and Their Roles in Plant Responses to Environmental Stimuli. *Frontiers in Plant Science*, 10, 2019.

[59] M. I. Capitani, L. J. Corzo-Rios, L. A. Chel-Guerrero, D. A. Betancur-Ancona, S. M. Nolasco, and M. C. Tomás. Rheological properties of aqueous dispersions of chia (Salvia hispanica L.) mucilage. *Journal of Food Engineering*, 149:70-77, Mar. 2015.

[60] Y. Capowiez, N. Bottinelli, S. Sammartino, E. Michel, and P. Jouquet. Morphological and functional characterisation of the burrow systems of six earthworm species (Lumbricidae). *Biology and Fertility of Soils*, 51(7):869-877, Oct. 2015.

[61] A. Carminati. A Model of Root Water Uptake Coupled with Rhizosphere Dynamics. *Vadose Zone Journal*, 11(3):vzj2011.0106, Aug. 2012.

[62] A. Carminati. Rhizosphere wettability decreases with root age: a problem or a strategy to increase water uptake of young roots? *Frontiers in Plant Science*, 4, 2013.

[63] A. Carminati, P. Benard, M. A. Ahmed, and M. Zarebanadkouki. Liquid bridges at the root-soil interface. *Plant and Soil*, 417(1):1-15, Aug. 2017.
[64] A. Carminati, E. Kroener, M. A. Ahmed, M. Zarebanadkouki, M. Holz, and T. Ghezzehei. Water for Carbon, Carbon for Water. *Vadose Zone Journal*, 15(2):vzj2015.04.0060, Feb. 2016.

[65] A. Carminati, A. B. Moradi, D. Vetterlein, P. Vontobel, E. Lehmann, U. Weller, H.-J. Vogel, and S. E. Oswald. Dynamics of soil water content in the rhizosphere. *Plant and Soil*, 332(1-2):163-176, July 2010.

[66] A. Carminati, C. L. Schneider, A. B. Moradi, M. Zarebanadkouki, D. Vetterlein, H.-J. Vogel, A. Hildebrandt, U. Weller, L. Schüler, and S. E. Oswald. How the Rhizosphere May Favor Water Availability to RootsAll rights reserved. No part of this periodical may be reproduced or transmitted in any form or by any means, electronic or mechanical, including photocopying, recording, or any information storage and retrieval system, without permission in writing from the publisher. *Vadose Zone Journal*, 10(3):988-998, Aug. 2011.

[67] A. Carminati and D. Vetterlein. Plasticity of rhizosphere hydraulic properties as a key for efficient utilization of scarce resources. *Annals of Botany*, 112(2):277-290, July 2013.

[68] A. Carminati, M. Zarebanadkouki, E. Kroener, M. A. Ahmed, and M. Holz. Biophysical rhizosphere processes affecting root water uptake. *Annals of Botany*, 118(4):561-571, Oct. 2016.

[69] A. B. D. Cassie and S. Baxter. Wettability of porous surfaces. *Transactions of the Faraday Society*, 40:546, 1944.

[70] E. A. Celaya, J. J. A. Aguirrezabala, and P. Chatzipantelidis. Implementation of an Adaptive BDF2 Formula and Comparison with the MATLAB Ode15s. *Procedia Computer Science*, 29:1014-1026, Jan. 2014.

[71] A. Chaboud. Isolation, purification and chemical composition of maize root cap slime. *Plant and Soil*, 73(3):395-402, Oct. 1983.

[72] J. M. Chaparro, D. V. Badri, and J. M. Vivanco. Rhizosphere microbiome assemblage is affected by plant development. *The ISME Journal*, 8(4):790-803, Apr. 2014. Number: 4 Publisher: Nature Publishing Group.

[73] H. R. Charles-Williams, R. Wengeler, K. Flore, H. Feise, M. J. Hounslow, and A. D. Salman. Granule nucleation and growth: Competing drop spreading and infiltration processes. *Powder Technology*, 206(1):63-71, Jan. 2011.

[74] C. Chenu. Influence of a fungal polysaccharide, scleroglucan, on clay microstructures. *Soil Biology and Biochemistry*, 21(2):299-305, Jan. 1989.

[75] C. Chenu. Clay—or sand—polysaccharide associations as models for the interface between micro-organisms and soil: water related properties and microstructure. In L. Brussaard and M. J. Kooistra, editors, *Soil Structure/Soil Biota Interrelationships*, pages 143-156. Elsevier, Amsterdam, Jan. 1993.

[76] C. Chenu and D. Cosentino. Microbial regulation of soil structural dynamics. *The architecture and biology of soils: life in inner space*, pages 37-70, Jan. 2011.
[77] C. Chenu, Y. Le Bissonnais, and D. Arrouays. Organic Matter Influence on Clay Wettability and Soil Aggregate Stability. *Soil Science Society of America Journal*, 64(4):1479-1486, 2000. _eprint: https://onlinelibrary.wiley.com/doi/pdf/10.2136/sssaj2000.6441479x.

[78] C. Chenu and A. Plante. Clay-sized organo-mineral complexes in a cultivation chronosequence: revisiting the concept of the "primary organo-mineral complex". *European Journal of Soil Science*, 56(4):596, 2006.

[79] C. Chenu and E. B. Roberson. Diffusion of glucose in microbial extracellular polysaccharide as affected by water potential. *Soil Biology and Biochemistry*, 28(7):877-884, July 1996.

[80] M. O. Ciantia, M. Arroyo, C. O'Sullivan, A. Gens, and T. Liu. Grading evolution and critical state in a discrete numerical model of Fontainebleau sand. *Géotechnique*, 69(1):1-15, Jan. 2019. Publisher: ICE Publishing.

[81] P. G. Ciarlet. *The finite element method for elliptic problems*. Number 40 in Classics in applied mathematics. Society for Industrial and Applied Mathematics, Philadelphia, Pa, unabridged republ. of the work first publ. by north-holland, 1978 edition, 2002.

[82] A. Clarke, T. D. Blake, K. Carruthers, and A. Woodward. Spreading and Imbibition of Liquid Droplets on Porous Surfaces. *Langmuir*, 18(8):2980-2984, Apr. 2002. Publisher: American Chemical Society.

[83] H. Crockford, S. Topalidis, and D. P. Richardson. Water repellency in a dry sclerophyll eucalypt forest — measurements and processes. *Hydrological Processes*, 5(4):405-420, 1991. _eprint: https://onlinelibrary.wiley.com/doi/pdf/10.1002/hyp.3360050408.

[84] L. Cueto-Felgueroso and R. Juanes. Nonlocal Interface Dynamics and Pattern Formation in Gravity-Driven Unsaturated Flow through Porous Media. *Physical Review Letters*, 101(24):244504, Dec. 2008. Publisher: American Physical Society.

[85] L. Cueto-Felgueroso and R. Juanes. Adaptive rational spectral methods for the linear stability analysis of nonlinear fourth-order problems. *Journal of Computational Physics*, 228(17):6536-6552, Sept. 2009.

[86] L. Cueto-Felgueroso and R. Juanes. A phase field model of unsaturated flow. *Water Resources Research*, 45(10), 2009. _eprint: https://onlinelibrary.wiley.com/doi/pdf/10.1029/2009WR007945.

[87] L. Cueto-Felgueroso and R. Juanes. Stability analysis of a phase-field model of gravity-driven unsaturated flow through porous media. *Physical Review E*, 79(3):036301, Mar. 2009.

[88] S. Czarnes, P. D. Hallett, A. G. Bengough, and I. M. Young. Root- and microbial-derived mucilages affect soil structure and water transport. *European Journal of Soil Science*, 51(3):435-443, 2000. _eprint: https://onlinelibrary.wiley.com/doi/pdf/10.1046/j.1365-2389.2000.00327.x.

[89] S. Czarnes, S. Hiller, A. R. Dexter, P. D. Hallett, and F. Bartoli. Root:soil adhesion in the maize rhizosphere: the rheological approach. *Plant and Soil*, 211(1):69-86, Apr. 1999.
[90] H. Darcy. *Les fontaines publiques de la ville de Dijon: Exposition et application des principes à suivre et des formules à employer dans les questions de distribution d'eau: Ouvrage terminé par un appendice relatif aux fournitures d'eau de plusieurs villes, au filtrage des eaux et à la fabrication des tuyaux de fonte, de plomb, de tôle et de bitume.* Victor Dalmont, éditeur, 1856.
Google-Books-ID: 42EUAAAAQAAJ.

[91] J.-F. Daïan. Equilibre et transferts en milieux poreux. Oct. 2013.

[92] P. G. de Gennes. Wetting: statics and dynamics. *Reviews of Modern Physics*, 57(3):827-863, July 1985. Publisher: American Physical Society.

[93] P.-G. d. de Gennes, F. Brochard-Wyart, and D. Quéré. *Gouttes, bulles, perles et ondes*. Belin, Paris, France, 2005. OCLC: 1191071573.

[94] B. P. Degens. Macro-aggregation of soils by biological bonding and binding mechanisms and the factors affecting these: a review. *Soil Research*, 35(3):431-460, 1997. Publisher: CSIRO PUBLISHING.

[95] E. Delfosse-Ribay, I. Djeran-Maigre, R. Cabrillac, and D. Gouvenot. Shear modulus and damping ratio of grouted sand. *Soil Dynamics and Earthquake Engineering*, 24(6):461-471, Aug. 2004.

[96] M. Denesuk, G. L. Smith, B. J. J. Zelinski, N. J. Kreidl, and D. R. Uhlmann. Capillary Penetration of Liquid Droplets into Porous Materials. *Journal of Colloid and Interface Science*, 158(1):114-120, June 1993.

[97] B. V. Derjaguin and N. V. Churaev. Polymolecular adsorption and capillary condensation in narrow slit pores. *Journal of Colloid and Interface Science*, 54(2):157-175, Feb. 1976.

[98] B. V. Derjaguin, N. V. Churaev, and V. M. Muller. *Surface Forces*. Springer US, Boston, MA, 1987.

[99] B. Despres. Différences Finies et Volumes Finis Master Mathématiques et Applications.

[100] D. A. DiCarlo. Stability of gravity-driven multiphase flow in porous media: 40 Years of advancements. *Water Resources Research*, 49(8):4531-4544, 2013.

[101] D. Diehl, M. Knott, and G. E. Schaumann. Purification effects show seed and root mucilage's ability to respond to changing rhizosphere conditions. *Biopolymers*, 114(8):e23561, Aug. 2023.

[102] H. A. Dijkstra, F. W. Wubs, A. K. Cliffe, Eusebius Doedel, I. F. Dragomirescu, B. Eckhardt, A. Y. Gelfgat, A. L. Hazel, V. Lucarini, A. G. Salinger, E. T. Phipps, J. Sanchez-Umbria, H. Schuttelaars, L. S. Tuckerman, and U. Thiele. Numerical Bifurcation Methods and their Application to Fluid Dynamics: Analysis beyond Simulation. *Commun. Comput. Phys.*, 15:1-45, 2014.

[103] S. H. Doerr. On standardizing the ‘Water Drop Penetration Time’ and the ‘Molarity of an Ethanol Droplet’ techniques to classify soil hydrophobicity: A case study
using medium textured soils. *Earth Surface Processes and Landforms*, 23(7):663-668, 1998. _eprint: https://onlinelibrary.wiley.com/doi/pdf/10.1002/%28SICI%291096-9837%28199807%2923%3A7%3C663%3A%3AAID-ESP909%3E3.0.CO%3B2-6.

[104] S. H. Doerr, R. A. Shakesby, and R. P. D. Walsh. Soil hydrophobicity variations with depth and particle size fraction in burned and unburned Eucalyptus globulus and Pinus pinaster forest terrain in the Águeda Basin, Portugal. *CATENA*, 27(1):25-47, June 1996.

[105] S. H. Doerr, R. A. Shakesby, and R. P. D. Walsh. Soil water repellency: its causes, characteristics and hydro-geomorphological significance. *Earth-Science Reviews*, 51(1):33-65, Aug. 2000.

[106] M. Doi. Onsager's variational principle in soft matter. *Journal of Physics: Condensed Matter*, 23(28):284118, June 2011.

[107] M. Doi. *Soft Matter Physics*. OUP Oxford, July 2013. Google-Books-ID: ccUaBj73PZsC.

[108] C. Doussan and E. Garrigues. Measuring and Imaging the Soil-root-water System with a Light Transmission 2D Technique. *BIO-PROTOCOL*, 9(6), 2019.

[109] C. Doussan, L. Pagès, and A. Pierret. Soil Exploration and Resource Acquisition by Plant Roots: An Architectural and Modelling Point of View. In E. Lichtfouse, M. Navarrete, P. De-baeke, S. Véronique, and C. Alberola, editors, *Sustainable Agriculture*, pages 583-600. Springer Netherlands, Dordrecht, 2009.

[110] C. Doussan, L. Pagès, and G. Vercambre. Modelling of the Hydraulic Architecture of Root Systems: An Integrated Approach to Water Absorption—Model Description. *Annals of Botany*, 81(2):213-223, Feb. 1998.

[111] C. Doussan, A. Pierret, E. Garrigues, and L. Pagès. Water Uptake by Plant Roots: II – Modelling of Water Transfer in the Soil Root-system with Explicit Account of Flow within the Root System – Comparison with Experiments. *Plant and Soil*, 283(1):99-117, May 2006.

[112] C. Doussan and S. Ruy. Prediction of unsaturated soil hydraulic conductivity with electrical conductivity. *Water Resources Research*, 45(10), 2009. _eprint: https://onlinelibrary.wiley.com/doi/pdf/10.1029/2008WR007309.

[113] M. DuBois, K. A. Gilles, J. K. Hamilton, P. A. Rebers, and F. Smith. Colorimetric Method for Determination of Sugars and Related Substances. *Analytical Chemistry*, 28(3):350-356, Mar. 1956.

[114] V. Dunbabin, Z. Rengel, and A. J. Diggle. Simulating form and function of root systems: efficiency of nitrate uptake is dependent on root system architecture and the spatial and temporal variability of nitrate supply. *Functional Ecology*, 18(2):204-211, 2004. _eprint: https://onlinelibrary.wiley.com/doi/pdf/10.1111/j.0269-8463.2004.00827.x.

[115] L. Durlofsky and J. F. Brady. Analysis of the Brinkman equation as a model for flow in porous media. *The Physics of Fluids*, 30(11):3329-3341, Nov. 1987.
[116] A. P. Edwards and J. M. Bremner. Microaggregates in Soils1. *Journal of Soil Science*, 18(1):64-73, 1967. _eprint: https://onlinelibrary.wiley.com/doi/pdf/10.1111/j.1365-2389.1967.tb01488.x.

[117] A. Einstein. Folgerungen aus den Capillaritätserscheinungen [AdP 4, 513 (1901)]. *Annalen der Physik*, 14(S1):87-98, 2005. Publisher: WILEY-VCH Verlag.

[118] A. Eissa and M. Zaki. The impact of global climatic changes on the aquatic environment. *Energy Policy - ENERG POLICY*, 4:251-259, Dec. 2011.

[119] M. Eliassi and R. J. Glass. On the porous-continuum modeling of gravity-driven fingers in unsaturated materials: Extension of standard theory with a hold-back-pile-up effect. *Water Resources Research*, 38(11):16-31-16-11, 2002.

[120] C. M. Elliott and H. Garcke. On the Cahn-Hilliard Equation with Degenerate Mobility. *SIAM Journal on Mathematical Analysis*, 27(2):404-423, Mar. 1996. Publisher: Society for Industrial and Applied Mathematics.

[121] C. M. Elliott and A. M. Stuart. The Global Dynamics of Discrete Semilinear Parabolic Equations. *SIAM Journal on Numerical Analysis*, 30(6):1622-1663, Dec. 1993. Publisher: Society for Industrial and Applied Mathematics.

[122] S. Engelnkemper, S. V. Gurevich, H. Uecker, D. Wetzel, and U. Thiele. Continuation for Thin Film Hydrodynamics and Related Scalar Problems. In A. Gelfgat, editor, *Computational Modelling of Bifurcations and Instabilities in Fluid Dynamics*, Computational Methods in Applied Sciences, pages 459-501. Springer International Publishing, Cham, 2019.

[123] A. Ern. *Aide-mémoire des éléments finis*. Dunod, Aug. 2005. Google-Books-ID: 5048WMM-Nfc8C.

[124] O. Esmaelipoor Jahromi, M. Knott, R. K. Mysore Janakiram, R. Rahim, and E. Kroener. Pore-scale simulation of mucilage drainage. *Vadose Zone Journal*, 21(5):e20218, 2022. _eprint: https://onlinelibrary.wiley.com/doi/pdf/10.1002/vzj2.20218.

[125] O. Esmaelipoor Jahromi, R. A. Patel, J. A. Huisman, and J. Vanderborght. Pore-scale simulation of mucilage drainage using phase field method. pages EGU-6120, May 2023. Conference Name: EGU General Assembly Conference Abstracts ADS Bibcode: 2023EGUGA..25.6120E.

[126] D. A. Farrell and W. E. Larson. Modeling the pore structure of porous media. *Water Resources Research*, 8(3):699-706, 1972. _eprint: https://onlinelibrary.wiley.com/doi/pdf/10.1029/WR008i003p00699.

[127] M. W. Farthing and F. L. Ogden. Numerical Solution of Richards' Equation: A Review of Advances and Challenges. *Soil Science Society of America Journal*, 81(6):1257-1269, 2017. _eprint: https://onlinelibrary.wiley.com/doi/pdf/10.2136/sssaj2017.02.0058.

[128] C. Fauchard, V. Guilbert, F. Sagnard, and M. Froumentin. Mesures de teneurs en eau volumique et massique sur du sable. *Bulletin des Laboratoires des Ponts et Chaussées*, (274):pp 109-125, Jan. 2009.
[129] D. S. Feeney, J. W. Crawford, T. Daniell, P. D. Hallett, N. Nunan, K. Ritz, M. Rivers, and I. M. Young. Three-dimensional Microorganization of the Soil-Root-Microbe System. *Microbial Ecology*, 52(1):151-158, July 2006.

[130] H.-C. Flemming and J. Wingender. Relevance of microbial extracellular polymeric substances (EPSs) - Part I: Structural and ecological aspects. *Water Science and Technology*, 43(6):1-8, Mar. 2001.

[131] H.-C. Flemming and J. Wingender. The biofilm matrix. *Nature Reviews Microbiology*, 8(9):623-633, Sept. 2010. Number: 9 Publisher: Nature Publishing Group.

[132] P. Forster, V. Ramaswamy, P. Artaxo, T. Berntsen, R. Betts, D. W. Fahey, J. Haywood, J. Lean, D. C. Lowe, G. Myhre, J. Nganga, R. Prinn, G. Raga, M. Schulz, R. van Dorland, G. Bodeker, O. Boucher, W. D. Collins, T. J. Conway, E. Dlugokencky, J. W. Elkins, D. Etheridge, P. Foukal, P. Fraser, M. Geller, F. Joos, C. D. Keeling, S. Kinne, K. Lassey, U. Lohmann, A. C. Manning, S. Montzka, D. Oram, K. O'Shaughnessy, S. Piper, G.-K. Plattner, M. Ponater, N. Ramankutty, G. Reid, D. Rind, K. Rosenlof, R. Sausen, D. Schwarzkopf, S. K. Solanki, G. Stenchikov, N. Stuber, T. Takemura, C. Textor, R. Wang, R. Weiss, and T. Whorf. Changes in Atmospheric Constituents and in Radiative Forcing. In S. Solomon, D. Qin, M. Manning, Z. Chen, M. Marquis, K. B. Averyt, M. Tignor, and H. L. Miller, editors, *Climate Change 2007: The Physical Science Basis. Contribution of Working Group I to the 4th Assessment Report of the Intergovernmental Panel on Climate Change*. Cambridge University Press, Cambridge, United Kingdom and New York, USA, Oct. 2007.

[133] M. Fortier, V. Lemaitre, A. Gaudry, B. Pawlak, A. Driouich, M.-L. Follet-Gueye, and M. Vicré. A fine-tuned defense at the pea root caps: Involvement of border cells and arabinogalactan proteins against soilborne diseases. *Frontiers in Plant Science*, 14, 2023.

[134] C. M. M. Franco, M. E. Tate, and J. M. Oades. Studies on non-wetting sands. 1. The role of intrinsic particulate organic-matter in the development of water-repellency in non-wetting sands. *Soil Research*, 33(2):253-263, 1995. Publisher: CSIRO PUBLISHING.

[135] E. I. Franses, O. A. Basaran, and C.-H. Chang. Techniques to measure dynamic surface tension. *Current Opinion in Colloid & Interface Science*, 1(2):296-303, Apr. 1996.

[136] D. Fredlund and A. Xing. Equations for the soil-water characteristic curve. *Canadian Geotechnical Journal*, 31(4):521-532, Aug. 1994. Publisher: NRC Research Press.

[137] R. A. Friesner, L. S. Tuckerman, B. C. Dornblaser, and T. V. Russo. A method for exponential propagation of large systems of stiff nonlinear differential equations. *Journal of Scientific Computing*, 4(4):327-354, Dec. 1989.

[138] T. Frohoff-Hülsmann and U. Thiele. Localized states in coupled Cahn-Hilliard equations. *IMA Journal of Applied Mathematics*, 86(5):924-943, Oct. 2021. arXiv:2010.14183 [cond-mat, physics:nlin].
[139] T. Frohoff-Hülsmann and U. Thiele. Nonreciprocal Cahn-Hilliard Model Emerges as a Universal Amplitude Equation. *Physical Review Letters*, 131(10):107201, Sept. 2023. Publisher: American Physical Society.

[140] T. Fürst, R. Vodák, M. Šír, and M. Bíl. On the incompatibility of Richards' equation and finger-like infiltration in unsaturated homogeneous porous media: RICHARDS' EQUATION AND FINGER-LIKE SOLUTIONS. *Water Resources Research*, 45(3), Mar. 2009.

[141] A. F. Galloway, P. Knox, and K. Krause. Sticky mucilages and exudates of plants: putative microenvironmental design elements with biotechnological value. *New Phytologist*, 225(4):1461-1469, 2020. _eprint: https://onlinelibrary.wiley.com/doi/pdf/10.1111/nph.16144.

[142] T. Gambaryan-Roisman. Liquids on porous layers: wetting, imbibition and transport processes. *Current Opinion in Colloid & Interface Science*, 19(4):320-335, Aug. 2014.

[143] Y. Gao, Y. Yang, W. Ling, H. Kong, and X. Zhu. Gradient Distribution of Root Exudates and Polycyclic Aromatic Hydrocarbons in Rhizosphere Soil. *Soil Science Society of America Journal*, 75(5):1694-1703, 2011. _eprint: https://onlinelibrary.wiley.com/doi/pdf/10.2136/sssaj2010.0244.

[144] W. Gardner. THE CAPILLARY POTENTIAL AND ITS RELATION TO SOILMOISTURE CONSTANTS. *Soil Science*, 10(5):357, Nov. 1920.

[145] W. R. Gardner. Dynamic Aspects of Soil-Water Availability to Plants. *Annual Review of Plant Physiology*, 16(1):323-342, 1965. _eprint: https://doi.org/10.1146/annurev.pp.16.060165.001543.

[146] C. W. Gear. *Numerical initial value problems in ordinary differential equations*. Prentice-Hall series in automatic computation. Prentice-Hall, Englewood Cliffs, N.J., 1971. OCLC: 159668.

[147] H. H. Gerke and M. T. van Genuchten. A dual-porosity model for simulating the preferential movement of water and solutes in structured porous media. *Water Resources Research*, 29(2):305-319, 1993. _eprint: https://onlinelibrary.wiley.com/doi/pdf/10.1029/92WR02339.

[148] C. Gessa and S. Deiana. Fibrillar structure of Ca polygalacturonate as a model for a soil-root interface. *Plant and Soil*, 129(2):211-217, Dec. 1990.

[149] M. Giaquinta and S. Hildebrandt. *Calculus of variations*. Number 310-311 in Grundlehren der mathematischen wissenschaften. Springer, Berlin, 1996.

[150] J.-M. Gobat, M. Aragno, and W. Matthey. *Le sol vivant: bases de pédologie, biologie des sols*. Science & ingénierie de l'environnement. Presses polytechniques et universitaires romandes, Lausanne, 3e éd. revue et augmentée edition, 2010.

[151] A. Golchin, J. M. Oades, J. O. Skjemstad, and P. Clarke. Soil structure and carbon cycling. *Soil Research*, 32(5):1043-1068, 1994. Publisher: CSIRO PUBLISHING.
[152] A. A. Golovin, A. A. Nepomnyashchy, S. H. Davis, and M. A. Zaks. Convective Cahn-Hilliard Models: From Coarsening to Roughening. *Physical Review Letters*, 86(8):1550–1553, Feb. 2001. Publisher: American Physical Society.

[153] P. J. Gregory, A. G. Bengough, T. S. George, and P. D. Hallett. Rhizosphere Engineering by Plants: Quantifying Soil-Root Interactions. In *Enhancing Understanding and Quantification of Soil-Root Growth Interactions*, pages 1–30. John Wiley & Sons, Ltd, 2013. Section: 1 _eprint: https://onlinelibrary.wiley.com/doi/pdf/10.2134/advagricsystmodel4.c1.

[154] P. M. Gresho, R. L. Sani, and M. S. Engelman. *Incompressible flow and the finite element method: advection-diffusion and isothermal laminar flow*. Wiley, Chichester Weinheim, reprinted with corr edition, 1999.

[155] D. Greve, S. Hartmann, and U. Thiele. Stick-slip dynamics in the forced wetting of polymer brushes. *Soft Matter*, 19(22):4041–4061, June 2023. Publisher: The Royal Society of Chemistry.

[156] P. H. Groenevelt and C. D. Grant. A new model for the soil-water retention curve that solves the problem of residual water contents. *European Journal of Soil Science*, 55(3):479–485, 2004. _eprint: https://onlinelibrary.wiley.com/doi/pdf/10.1111/j.1365-2389.2004.00617.x.

[157] M. L. Gruwel. In Situ Magnetic Resonance Imaging of Plant Roots. *Vadose Zone Journal*, 13(3):vzj2013.08.0158, Mar. 2014.

[158] G. Grün. Degenerate Parabolic Differential Equations of Fourth Order and a Plasticity Model with Non-Local Hardening. *Zeitschrift für Analysis und ihre Anwendungen*, 14(3):541–574, Sept. 1995.

[159] G. Grün. On the convergence of entropy consistent schemes for lubrication type equations in multiple space dimensions. *Mathematics of Computation*, 72(243):1251–1279, July 2003.

[160] G. Grün and M. Rumpf. Nonnegativity preserving convergent schemes for the thin film equation. *Numerische Mathematik*, 87(1):113–152, Nov. 2000.

[161] G. Grün and M. Rumpf. Simulation of singularities and instabilities arising in thin film flow. *European Journal of Applied Mathematics*, 12(3):293–320, June 2001. Publisher: Cambridge University Press.

[162] K. P. Hapgood, J. D. Litster, S. R. Biggs, and T. Howes. Drop Penetration into Porous Powder Beds. *Journal of Colloid and Interface Science*, 253(2):353–366, Sept. 2002.

[163] W. Hardy. III. The spreading of fluids on glass. *The London, Edinburgh, and Dublin Philosophical Magazine and Journal of Science*, 38(223):49–55, July 1919. Publisher: Taylor & Francis _eprint: https://doi.org/10.1080/14786440708635928.

[164] R. J. Harper and R. J. Gilkes. Soil attributes related to water repellency and the utility of soil survey for predicting its occurrence. *Soil Research*, 32(5):1109–1124, 1994. Publisher: CSIRO PUBLISHING.
[165] S. Hartmann, C. Diddens, M. Jalaal, and U. Thiele. Sessile drop evaporation in a gap – crossover between diffusion-limited and phase transition-limited regime. *Journal of Fluid Mechanics*, 960:A32, Apr. 2023. Publisher: Cambridge University Press.

[166] S. Hartmann, J. Diekmann, D. Greve, and U. Thiele. Drops on polymer brushes – advances in thin-film modelling of adaptive substrates, Dec. 2023. arXiv:2311.07307 [cond-mat, physics:physics].

[167] M. Hassanizadeh and W. G. Gray. General conservation equations for multi-phase systems: 3. Constitutive theory for porous media flow. *Advances in Water Resources*, 3(1):25–40, Mar. 1980.

[168] S. M. Hassanizadeh and W. G. Gray. Thermodynamic basis of capillary pressure in porous media. *Water Resources Research*, 29(10):3389–3405, 1993. _eprint: https://onlinelibrary.wiley.com/doi/pdf/10.1029/93WR01495.

[169] M. Heil and A. L. Hazel. oomph-lib – An Object-Oriented Multi-Physics Finite-Element Library. In H.-J. Bungartz and M. Schäfer, editors, *Fluid-Structure Interaction, Lecture Notes in Computational Science and Engineering*, pages 19–49, Berlin, Heidelberg, 2006. Springer.

[170] D. Hillel. *L'Eau et le sol: principes et processus physiques*. Academia, Louvain-la-Neuve, 1988. OCLC: 19693229.

[171] L. Hiltner. Über nevere Erfahrungen und Probleme auf dem Gebiet der Boden Bakteriologie und unter besonderer Beurchsichtigung der Grundungung und Broche. *Arbeit. Deut. Landw. Ges. Berlin*, 98:59–78, 1904.

[172] A. C. Hindmarsh and L. R. Petzold. Algorithms and software for ordinary differential equations and differential-algebraic equations, Part I: Euler methods and error estimation. *Computer in Physics*, 9(1):34–41, Jan. 1995.

[173] A. C. Hindmarsh and L. R. Petzold. Algorithms and software for ordinary differential equations and differential-algebraic equations, Part II: Higher-order methods and software packages. *Computer in Physics*, 9(2):148–155, Mar. 1995.

[174] P. Hinsinger, A. G. Bengough, D. Vetterlein, and I. M. Young. Rhizosphere: biophysics, biogeochemistry and ecological relevance. *Plant and Soil*, 321(1):117–152, Aug. 2009.

[175] M. Hochbruck, C. Lubich, and H. Selhofer. Exponential Integrators for Large Systems of Differential Equations. *SIAM Journal on Scientific Computing*, 19(5):1552–1574, Jan. 1998. Publisher: Society for Industrial and Applied Mathematics.

[176] M. P. Holl, A. J. Archer, and U. Thiele. Efficient calculation of phase coexistence and phase diagrams: application to a binary phase-field-crystal model. *Journal of Physics: Condensed Matter*, 33(11):115401, Dec. 2020. Publisher: IOP Publishing.

[177] R. K. Holman, M. J. Cima, S. A. Uhland, and E. Sachs. Spreading and Infiltration of Inkjet-Printed Polymer Solution Droplets on a Porous Substrate. *Journal of Colloid and Interface Science*, 249(2):432–440, May 2002.
[178] M. Holz, M. Leue, M. A. Ahmed, P. Benard, H. H. Gerke, and A. Carminati. Spatial Distribution of Mucilage in the Rhizosphere Measured With Infrared Spectroscopy. *Frontiers in Environmental Science*, 6, 2018.

[179] M. Holz, M. Zarebanadkouki, A. Kaestner, Y. Kuzyakov, and A. Carminati. Rhizodeposition under drought is controlled by root growth rate and rhizosphere water content. *Plant and Soil*, 423(1):429-442, Feb. 2018.

[180] T. J. R. Hughes. *The finite element method: linear static and dynamic finite element analysis*. Dover Publications, Mineola, NY, 2000.

[181] Intergovernmental Panel on Climate Change (IPCC). *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, 2014.

[182] J. N. Israelachvili. *Intermolecular and surface forces*. Academic press, Burlington (Mass.), 3rd ed edition, 2011.

[183] M. Javaux, T. Schröder, J. Vanderborght, and H. Vereecken. Use of a Three-Dimensional Detailed Modeling Approach for Predicting Root Water UptakeAll rights reserved. No part of this periodical may be reproduced or transmitted in any form or by any means, electronic or mechanical, including photocopying, recording, or any information storage and retrieval system, without permission in writing from the publisher. *Vadose Zone Journal*, 7(3):1079–1088, Aug. 2008.

[184] C. G. Jones, J. H. Lawton, and M. Shachak. Organisms as Ecosystem Engineers. *Oikos*, 69(3):373-386, 1994. Publisher: [Nordic Society Oikos, Wiley].

[185] D. L. Jones, C. Nguyen, and R. D. Finlay. Carbon flow in the rhizosphere: carbon trading at the soil-root interface. *Plant and Soil*, 321(1):5-33, Aug. 2009.

[186] H. Jorda, K. Huber, A. Kunkel, J. Vanderborght, M. Javaux, C. Oberdörster, K. Hammel, and A. Schnepf. Mechanistic modeling of pesticide uptake with a 3D plant architecture model. *Environmental Science and Pollution Research*, 28(39):55678-55689, Oct. 2021.

[187] H. Jorda, A. Perelman, N. Lazarovitch, and J. Vanderborght. Exploring Osmotic Stress and Differences between Soil-Root Interface and Bulk Salinities. *Vadose Zone Journal*, 17(1):170029, 2018. _eprint: https://onlinelibrary.wiley.com/doi/pdf/10.2136/vzj2017.01.0029.

[188] C. JOUANY, C. CHENU, and P. CHASSIN. Détermination de la mouillabilité des constituants du sol à partir de mesures d'angles de contact : revue bibliographique. *Détermination de la mouillabilité des constituants du sol à partir de mesures d'angles de contact : revue bibliographique*, 30(1):33-47, 1992. Place: Plaisir Publisher: AFES.

[189] J. Jurin. II. An account of some experiments shown before the Royal Society; with an enquiry into the cause of the ascent and suspension of water in capillary tubes. *Philosophical Transactions of the Royal Society of London*, 30(355):739-747, Dec. 1719.
[190] S. Kalliadasis and U. Thiele, editors. *Thin films of soft matter*. Number 490 in Courses and lectures / International Centre for Mechanical Sciences. Springer, Wien, 2007.

[191] O. Kap, S. Hartmann, H. Hoek, S. de Beer, I. Siretanu, U. Thiele, and F. Mugele. Nonequilibrium configurations of swelling polymer brush layers induced by spreading drops of weakly volatile oil. *The Journal of Chemical Physics*, 158(17):174903, May 2023.

[192] D. Kashyap and T. Agarwal. Climate Change, Water Resources, and Agriculture: Impacts and Adaptation Measures. In V. Venkatramanan, S. Shah, and R. Prasad, editors, *Global Climate Change and Environmental Policy: Agriculture Perspectives*, pages 227–256. Springer, Singapore, 2020.

[193] I. G. Kevrekidis, B. Nicolaenko, and J. C. Scovel. Back in the Saddle Again: A Computer Assisted Study of the Kuramoto–Sivashinsky Equation. *SIAM Journal on Applied Mathematics*, 50(3):760–790, June 1990. Publisher: Society for Industrial and Applied Mathematics.

[194] P. M. King. Comparison of methods for measuring severity of water repellence of sandy soils and assessment of some factors that affect its measurement. *Soil Research*, 19(3):275–285, 1981. Publisher: CSIRO PUBLISHING.

[195] E. M. Knee, F.-C. Gong, M. Gao, M. Teplitski, A. R. Jones, A. Foxworthy, A. J. Mort, and W. D. Bauer. Root Mucilage from Pea and Its Utilization by Rhizosphere Bacteria as a Sole Carbon Source. *Molecular Plant-Microbe Interactions®*, 14(6):775–784, June 2001. Publisher: Scientific Societies.

[196] K. Kosugi. Lognormal Distribution Model for Unsaturated Soil Hydraulic Properties. *Water Resources Research*, 32(9):2697–2703, 1996. _eprint:_ https://onlinelibrary.wiley.com/doi/pdf/10.1029/96WR01776.

[197] B. Krauskopf, H. M. Osinga, J. Galán-Vioque, and E. Doedel, editors. *Numerical continuation methods for dynamical systems: path following and boundary value problems: dedicated to Eu-sebius J. Doedel for his 60th birthday*. Understanding Complex Systems. Springer, Dordrecht, 2007.

[198] E. Kroener, M. A. Ahmed, and A. Carminati. Roots at the percolation threshold. *Physical Review E*, 91(4):042706, Apr. 2015. Publisher: American Physical Society.

[199] E. Kroener, M. Holz, M. Zarebanadkouki, M. Ahmed, and A. Carminati. Effects of Mucilage on Rhizosphere Hydraulic Functions Depend on Soil Particle Size. *Vadose Zone Journal*, 17(1):170056, 2018. _eprint:_ https://onlinelibrary.wiley.com/doi/pdf/10.2136/vzj2017.03.0056.

[200] E. Kroener, M. Zarebanadkouki, M. Bittelli, and A. Carminati. Simulation of root water uptake under consideration of nonequilibrium dynamics in the rhizosphere. *Water Resources Research*, 52(8):5755–5770, 2016. _eprint:_ https://onlinelibrary.wiley.com/doi/pdf/10.1002/2015WR018579.
[201] E. Kroener, M. Zarebanadkouki, A. Kaestner, and A. Carminati. Nonequilibrium water dynamics in the rhizosphere: How mucilage affects water flow in soils. *Water Resources Research*, 50(8):6479-6495, 2014. _eprint: https://onlinelibrary.wiley.com/doi/pdf/10.1002/2013WR014756.

[202] F. Lafolie, I. Cousin, P.-A. Marron, A. Mollier, V. Pot, N. Moitrier, N. Moitrier, and C. Nouguier. The "VSOIL" modeling platform. *Revue forestière française*, 66:np, 2014. Publisher: AgroParis-Tech.

[203] J. Lallement. *Modélisation et simulation numérique d'écoulements de films minces avec effet de mouillage partiel*. These de doctorat, Toulouse, ISAE, Feb. 2019.

[204] M. Landl, A. Haupenthal, D. Leitner, E. Kroener, D. Vetterlein, R. Bol, H. Vereecken, J. Vanderborght, and A. Schnepf. Simulating rhizodeposition patterns around growing and exuding root systems. *in silico Plants*, 3(2):diab028, July 2021.

[205] M. Landl, K. Huber, A. Schnepf, J. Vanderborght, M. Javaux, A. Glyn Bengough, and H. Vereecken. A new model for root growth in soil with macropores. *Plant and Soil*, 415(1):99-116, June 2017.

[206] M. Landl, M. Phalempin, S. Schlüter, D. Vetterlein, J. Vanderborght, E. Kroener, and A. Schnepf. Modeling the Impact of Rhizosphere Bulk Density and Mucilage Gradients on Root Water Uptake. *Frontiers in Agronomy*, 3, 2021.

[207] I. Langmuir. Oil Lenses on Water and the Nature of Monomolecular Expanded Films. *The Journal of Chemical Physics*, 1(11):756-776, Nov. 2004.

[208] H. P. Langtangen. Solving nonlinear ODE and PDE problems, 2016.

[209] S. Le Gall, A. Bérard, D. Page, L. Lanoe, N. Bertin, and C. Doussan. Increased exopolysaccharide production and microbial activity affect soil water retention and field performance of tomato under water deficit. *Rhizosphere*, 19:100408, Sept. 2021.

[210] D. Leitner, A. Schnepf, S. Klepsch, and T. Roosevelt. Comparison of nutrient uptake between three-dimensional simulation and an averaged root system model. *Plant Biosystems - An International Journal Dealing with all Aspects of Plant Biology*, 144(2):443-447, June 2010. Publisher: Taylor & Francis _eprint: https://doi.org/10.1080/11263501003726334.

[211] E. C. Leong and H. Rahardjo. Review of Soil-Water Characteristic Curve Equations. *Journal of Geotechnical and Geoenvironmental Engineering*, 123(12):1106-1117, Dec. 1997. Publisher: American Society of Civil Engineers.

[212] M. Leverett. Capillary Behavior in Porous Solids. *Transactions of the AIME*, 142(01):152-169, Dec. 1941.

[213] R. M. Lewis G. N. *Thermodynamics and the free energy of chemical substances*. MGH, 1923.
[214] M. Li, T. Shinano, and T. Tadano. Distribution of exudates of lupin roots in the rhizosphere under phosphorus deficient conditions. *Soil Science and Plant Nutrition*, 43(1):237-245, Mar. 1997.

[215] Z. Li, S. Escoffier, and P. Kotronis. Using centrifuge tests data to identify the dynamic soil properties: Application to Fontainebleau sand. *Soil Dynamics and Earthquake Engineering*, 52:77-87, Sept. 2013.

[216] K.-Y. Lin, J. R. Daniel, and R. L. Whistler. Structure of chia seed polysaccharide exudate. *Carbohydrate Polymers*, 23(1):13-18, Jan. 1994.

[217] D. R. Linden, P. F. Hendrix, D. C. Coleman, and P. C. van Vliet. Faunal Indicators of Soil Quality. In *Defining Soil Quality for a Sustainable Environment*, pages 91-106. John Wiley & Sons, Ltd, 1994. Section: 6 _eprint: https://onlinelibrary.wiley.com/doi/pdf/10.2136/ssaspecpub35.c6.

[218] A.-S. Lissy. *Imagerie quantitative non invasive et expérimentation pour l'intégration des phénomènes d'écoulement non uniforme dans les modèles de transfert de masse en milieux poreux hétérogènes : Application aux sols structurés*. These de doctorat, Avignon, Feb. 2019.

[219] H. Liu, Z. Ju, J. Bachmann, R. Horton, and T. Ren. Moisture-Dependent Wettability of Artificial Hydrophobic Soils and Its Relevance for Soil Water Desorption Curves. *Soil Science Society of America Journal*, 76(2):342-349, 2012. _eprint: https://onlinelibrary.wiley.com/doi/pdf/10.2136/ssaj2011.0081.

[220] M. Lucas, S. Schlüter, H.-J. Vogel, and D. Vetterlein. Roots compact the surrounding soil depending on the structures they encounter. *Scientific Reports*, 9(1):16236, Nov. 2019. Number: 1 Publisher: Nature Publishing Group.

[221] S. Luo, N. Lu, C. Zhang, and W. Likos. Soil water potential: A historical perspective and recent breakthroughs. *Vadose Zone Journal*, 21(4):e20203, 2022. _eprint: https://onlinelibrary.wiley.com/doi/pdf/10.1002/vzj2.20203.

[222] R. J. Luxmoore. Micro-, meso-, and macroporosity of soil. *Soil Sci. Soc. Am. J.*; (United States), 45:3, May 1981. Institution: Oak Ridge National Lab., TN.

[223] F. H. Macdougall. *Thermodynamics and Chemistry*. John Wiley & Sons, 1939.

[224] D. Mandal and J. Jayaprakash. Water repellency of soils in the lower Himalayan regions of India: impact of land use. *Current Science*, pages 148-152, 2009. Publisher: JSTOR.

[225] V. Manet. Méthode des Éléments Finis: vulgarisation des aspects mathématiques et illustration de la méthode. page 429, Dec. 2018.

[226] R. Marasco, M. Fusi, M. Mosqueira, J. M. Booth, F. Rossi, M. Cardinale, G. Michoud, E. Rolli, G. Mugnai, L. Vergani, S. Borin, R. De Philippis, A. Cherif, and D. Daffonchio. Rhizosheath-root system changes exopolysaccharide content but stabilizes bacterial community across contrasting seasons in a desert environment. *Environmental Microbiome*, 17(1):14, Apr. 2022.
[227] A. Marchand. *Mouillage statique et dynamique : Influences géométriques aux échelles moléculaires*. phdthesis, Université Paris-Diderot - Paris VII, Nov. 2011.

[228] A. Marmur. Drop penetration into a thin porous medium. *Journal of Colloid and Interface Science*, 123(1):161-169, May 1988.

[229] B. Mary, C. Fresneau, J. L. Morel, and A. Mariotti. C and N cycling during decomposition of root mucilage, roots and glucose in soil. *Soil Biology and Biochemistry*, 25(8):1005-1014, Aug. 1993.

[230] M. Ma'shum and V. C. Farmer. Origin and assessment of water repellency of a sandy South Australian soil. *Soil Research*, 23(4):623-626, 1985. Publisher: CSIRO PUBLISHING.

[231] M. E. McCully and J. S. Boyer. The expansion of maize root-cap mucilage during hydration. 3. Changes in water potential and water content. *Physiologia Plantarum*, 99(1):169-177, 1997. _eprint: https://onlinelibrary.wiley.com/doi/pdf/10.1111/j.1399-3054.1997.tb03445.x.

[232] C. R. McKee and A. C. Bumb. Flow-Testing Coalbed Methane Production Wells in the Presence of Water and Gas. *SPE Formation Evaluation*, 2(04):599-608, Dec. 1987.

[233] McKee, C.R and Bumb, A.C. The Importance of Unsaturated Flow Parameters in Designing a Monitoring System for Hazardous Wastes and Environmental Emergencies. pages 50-58, Houston, 1984.

[234] L. M. Milne-Thomson. *The calculus of finite differences*. AMS Chelsea Pub, Providence, R.I, 2000.

[235] T. Mimmo, C. Marzadori, and C. E. Gessa. Does the degree of pectin esterification influence aluminium sorption by the root apoplast? *Plant and Soil*, 314(1):159-168, Jan. 2009.

[236] V. S. Mitlin. Dewetting of Solid Surface: Analogy with Spinodal Decomposition. *Journal of Colloid and Interface Science*, 156(2):491-497, Mar. 1993.

[237] L. Modica. The gradient theory of phase transitions and the minimal interface criterion. *Archive for Rational Mechanics and Analysis*, 98(2):123-142, June 1987.

[238] A. B. Moradi, A. Carminati, D. Vetterlein, P. Vontobel, E. Lehmann, U. Weller, J. W. Hopmans, H.-J. Vogel, and S. E. Oswald. Three-dimensional visualization and quantification of water content in the rhizosphere. *New Phytologist*, 192(3):653-663, 2011. _eprint: https://onlinelibrary.wiley.com/doi/pdf/10.1111/j.1469-8137.2011.03826.x.

[239] E. R. Morris, A. N. Cutler, S. B. Ross-Murphy, D. A. Rees, and J. Price. Concentration and shear rate dependence of viscosity in random coil polysaccharide solutions. *Carbohydrate Polymers*, 1(1):5-21, Sept. 1981.

[240] Y. Mualem. A new model for predicting the hydraulic conductivity of unsaturated porous media. *Water Resources Research*, 12(3):513-522, 1976. _eprint: https://onlinelibrary.wiley.com/doi/pdf/10.1029/WR012i003p00513.
[241] G. J. H. Müehl, J. Rüehlmann, M.-O. Goebel, and J. Bachmann. Application of confocal laser scanning microscopy (CLSM) to visualize the effect of porous media wettability on unsaturated pore water configuration. *Journal of Soils and Sediments*, 12(1):75-85, Jan. 2012.

[242] M. Müller and L. G. MacDowell. Wetting of polymer liquids: Monte Carlo simulations and self-consistent field calculations. *Journal of Physics: Condensed Matter*, 15(19):R609, May 2003.

[243] H. Naseem, M. Ahsan, M. A. Shahid, and N. Khan. Exopolysaccharides producing rhizobacteria and their role in plant growth and drought tolerance. *Journal of Basic Microbiology*, 58(12):1009-1022, Dec. 2018.

[244] M. Naveed, M. A. Ahmed, P. Benard, L. K. Brown, T. S. George, A. G. Bengough, T. Roosevelt, N. Koebernick, and P. D. Hallett. Surface tension, rheology and hydrophobicity of rhizodeposits and seed mucilage influence soil water retention and hysteresis. *Plant and Soil*, 437(1):65-81, Apr. 2019.

[245] M. Naveed, L. K. Brown, A. C. Raffan, T. S. George, A. G. Bengough, T. Roosevelt, I. Sinclair, N. Koebernick, L. Cooper, C. A. Hackett, and P. D. Hallett. Plant exudates may stabilize or weaken soil depending on species, origin and time. *European Journal of Soil Science*, 68(6):806-816, 2017. _eprint: https://onlinelibrary.wiley.com/doi/pdf/10.1111/ejss.12487.

[246] M. Nazari. Plant mucilage components and their functions in the rhizosphere. *Rhizosphere*, 18:100344, June 2021.

[247] M. Nazari, S. Riebeling, C. C. Banfield, A. Akale, M. Crosta, K. Mason-Jones, M. A. Dippold, and M. A. Ahmed. Mucilage Polysaccharide Composition and Exudation in Maize From Contrasting Climatic Regions. *Frontiers in Plant Science*, 11, 2020.

[248] G. Neumann, S. Bott, M. Ohler, H.-P. Mock, R. Lippmann, R. Grosch, and K. Smalla. Root exudation and root development of lettuce (Lactuca sativa L. cv. Tizian) as affected by different soils. *Frontiers in Microbiology*, 5, 2014.

[249] C. Nguyen, F. Froux, S. Recous, T. Morvan, and C. Robin. Net N immobilisation during the biodegradation of mucilage in soil as affected by repeated mineral and organic fertilisation. *Nutrient Cycling in Agroecosystems*, 80(1):39-47, Jan. 2008.

[250] J. L. Nieber, T. W. J. Bauters, T. S. Steenhuis, and J. Y. Parlange. Numerical simulation of experimental gravity-driven unstable flow in water repellent sand. *Journal of Hydrology*, 231-232:295-307, 2000.

[251] F. Néron. *Petit précis d'agriculture: de la politique à la technique le contexte de l'exploitation agricole, les productions végétales et animales, les différents types d'agriculture*. Agriproduction. Éd. France agricole, Paris, 2e éd edition, 2014.

[252] E. Oleghe, M. Naveed, E. M. Baggs, and P. D. Hallett. Plant exudates improve the mechanical conditions for root penetration through compacted soils. *Plant and Soil*, 421(1):19-30, Dec. 2017.
[253] C. T. Omuto. Biexponential model for water retention characteristics. *Geoderma*, 149(3):235-242, Mar. 2009.

[254] D. Or, S. Phutane, and A. Dechesne. Extracellular Polymeric Substances Affecting Pore-Scale Hydrologic Conditions for Bacterial Activity in Unsaturated SoilsAll rights reserved. No part of this periodical may be reproduced or transmitted in any form or by any means, electronic or mechanical, including photocopying, recording, or any information storage and retrieval system, without permission in writing from the publisher. *Vadose Zone Journal*, 6(2):298-305, May 2007.

[255] D. Or, B. F. Smets, J. M. Wraith, A. Dechesne, and S. P. Friedman. Physical constraints affecting bacterial habitats and activity in unsaturated porous media – a review. *Advances in Water Resources*, 30(6):1505–1527, June 2007.

[256] D. Or, M. Tuller, and J. Wraith. Water potential. In *Encyclopedia of Soils in the Environment*, pages 270–277. Jan. 2004. Journal Abbreviation: Encyclopedia of Soils in the Environment.

[257] A. Oron, S. H. Davis, and S. G. Bankoff. Long-scale evolution of thin liquid films. *Reviews of Modern Physics*, 69(3):931–980, July 1997. Publisher: American Physical Society.

[258] H. M. I. Osborn, F. Lochey, L. Mosley, and D. Read. Analysis of polysaccharides and monosaccharides in the root mucilage of maize (Zea mays L.) by gas chromatography. *Journal of Chromatography A*, 831(2):267–276, Jan. 1999.

[259] J. Pang, M. H. Ryan, K. H. M. Siddique, and R. J. Simpson. Unwrapping the rhizosheath. *Plant and Soil*, 418(1):129–139, Sept. 2017.

[260] R. I. Papendick and G. S. Campbell. Theory and Measurement of Water Potential. In *Water Potential Relations in Soil Microbiology*, pages 1–22. John Wiley & Sons, Ltd, 1981. Section: 1 _eprint: https://onlinelibrary.wiley.com/doi/pdf/10.2136/sssaspecpub9.c1.

[261] J. C. Parker and A. J. Valocchi. Constraints on the Validity of Equilibrium and First-Order Kinetic Transport Models in Structured Soils. *Water Resources Research*, 22(3):399–407, 1986. _eprint: https://onlinelibrary.wiley.com/doi/pdf/10.1029/WR022i003p00399.

[262] J.-Y. Parlange and D. E. Hill. Theoretical Analysis of wetting front instability in soils. *Soil Science*, 122(4), 1976.

[263] R. D. Passo, H. Garcke, and G. Grün. On a Fourth-Order Degenerate Parabolic Equation: Global Entropy Estimates, Existence, and Qualitative Behavior of Solutions. *SIAM Journal on Mathematical Analysis*, 29(2):321–342, Jan. 1998. Publisher: Society for Industrial and Applied Mathematics.

[264] R. D. Passo, L. Giacomelli, and G. Grün. A waiting time phenomenon for thin film equations.

[265] E. Paterson, T. Gebbing, C. Abel, A. Sim, and G. Telfer. Rhizodeposition shapes rhizosphere microbial community structure in organic soil. *New Phytologist*, 173(3):600–610, 2007. _eprint: https://onlinelibrary.wiley.com/doi/pdf/10.1111/j.1469-8137.2006.01931.x.
[266] X. Peng, P. D. Hallett, B. Zhang, and R. Horn. Physical response of rigid and non-rigid soils to analogues of biological exudates. *European Journal of Soil Science*, 62(5):676-684, 2011. _eprint: https://onlinelibrary.wiley.com/doi/pdf/10.1111/j.1365-2389.2011.01383.x.

[267] E. Perrier, C. Mullon, M. Rieu, and G. de Marsily. Computer construction of fractal soil structures: Simulation of their hydraulic and shrinkage properties. *Water Resources Research*, 31(12):2927-2943, 1995. _eprint: https://onlinelibrary.wiley.com/doi/pdf/10.1029/95WR02214.

[268] M. Phalempin, E. Lippold, D. Vetterlein, and S. Schlüter. Soil texture and structure heterogeneity predominantly governs bulk density gradients around roots. *Vadose Zone Journal*, 20(5):e20147, 2021. _eprint: https://onlinelibrary.wiley.com/doi/pdf/10.1002/vzj2.20147.

[269] J. Philip. General Method of Exact Solution of the Concentration?Dependent Diffusion Equation. *Australian Journal of Physics*, 13(1):1, 1960.

[270] J. R. Philip. ABSOLUTE THERMODYNAMIC FUNCTIONS IN SOIL-WATER STUDIES. *Soil Science*, 89(2):111, Feb. 1960.

[271] L. M. Pismen and Y. Pomeau. Disjoining potential and spreading of thin liquid layers in the diffuse-interface model coupled to hydrodynamics. *Physical Review E*, 62(2):2480-2492, Aug. 2000. Publisher: American Physical Society.

[272] V. Pot, S. Peth, O. Monga, L. E. Vogel, A. Genty, P. Garnier, L. Vieublé-Gonod, M. Ogurreck, F. Beckmann, and P. C. Baveye. Three-dimensional distribution of water and air in soil pores: Comparison of two-phase two-relaxation-times lattice-Boltzmann and morphological model outputs with synchrotron X-ray computed tomography data. *Advances in Water Resources*, 84:87-102, Oct. 2015.

[273] V. Pot, X. Portell, W. Otten, P. Garnier, O. Monga, and P. C. Baveye. Accounting for soil architecture and microbial dynamics in microscale models: Current practices in soil science and the path ahead. *European Journal of Soil Science*, 73(1):e13142, 2022. _eprint: https://onlinelibrary.wiley.com/doi/pdf/10.1111/ejss.13142.

[274] J. M. Raaijmakers, I. De Bruijn, O. Nybroe, and M. Ongena. Natural functions of lipopeptides from Bacillus and Pseudomonas: more than surfactants and antibiotics. *FEMS Microbiology Reviews*, 34(6):1037-1062, Nov. 2010.

[275] M. A. Rab, R. E. Haling, S. R. Aarons, M. Hannah, I. M. Young, and D. Gibson. Evaluation of X-ray computed tomography for quantifying macroporosity of loamy pasture soils. *Geoderma*, 213:460-470, Jan. 2014.

[276] A. Raetz and B. Schweizer. Hysteresis models and gravity fingering in porous media. ZAMM - *Journal of Applied Mathematics and Mechanics*, 2012. Publisher: WILEY-VCH Verlag.

[277] R. Rahim, O. E. Jahromi, W. Amelung, and E. Kroener. Rhizosheath formation depends on mucilage concentration and water content. *Plant and Soil*, Nov. 2023.
[278] D. P. Rasse, C. Rumpel, and M.-F. Dignac. Is soil carbon mostly root carbon? Mechanisms for a specific stabilisation. *Plant and Soil*, 269(1):341–356, Feb. 2005.

[279] N. Ray, A. Rupp, R. Schulz, and P. Knabner. Old and New Approaches Predicting the Diffusion in Porous Media. *Transport in Porous Media*, 124(3):803–824, Sept. 2018.

[280] D. Read, P. Gregory, and A. Bell. Physical properties of axenic maize root mucilage. *Plant and Soil*, 211(1):87–91, Apr. 1999.

[281] D. B. Read, A. G. Bengough, P. J. Gregory, J. W. Crawford, D. Robinson, C. M. Scrimgeour, I. M. Young, K. Zhang, and X. Zhang. Plant roots release phospholipid surfactants that modify the physical and chemical properties of soil. *New Phytologist*, 157(2):315–326, 2003. _eprint: https://onlinelibrary.wiley.com/doi/pdf/10.1046/j.1469-8137.2003.00665.x.

[282] D. B. Read and P. J. Gregory. Surface tension and viscosity of axenic maize and lupin root mucilages. *The New Phytologist*, 137(4):623–628, Dec. 1997. Publisher: Cambridge University Press.

[283] D. A. Rees and E. J. Welsh. Secondary and Tertiary Structure of Polysaccharides in Solutions and Gels. *Angewandte Chemie International Edition in English*, 16(4):214–224, 1977. _eprint: https://onlinelibrary.wiley.com/doi/pdf/10.1002/anie.197702141.

[284] L. A. Richards. CAPILLARY CONDUCTION OF LIQUIDS THROUGH POROUS MEDIUMS. *Physics*, 1(5):318–333, Nov. 1931.

[285] V. Richefeu, F. Radjai, and J.-Y. Delenne. Lattice Boltzmann modelling of liquid distribution in unsaturated granular media. *Computers and Geotechnics*, 80:353–359, Dec. 2016.

[286] M. C. Rillig and D. L. Mummey. Mycorrhizas and soil structure. *New Phytologist*, 171(1):41–53, 2006. _eprint: https://onlinelibrary.wiley.com/doi/pdf/10.1111/j.1469-8137.2006.01750.x.

[287] E. B. Roberson, C. Chenu, and M. K. Firestone. Microstructural changes in bacterial exopolysaccharides during desiccation. *Soil Biology and Biochemistry*, 25(9):1299–1301, Sept. 1993.

[288] E. B. Roberson and M. K. Firestone. Relationship between Desiccation and Exopolysaccharide Production in a Soil Pseudomonas sp. *Applied and Environmental Microbiology*, 58(4):1284–1291, Apr. 1992. Publisher: American Society for Microbiology.

[289] S. Roman, C. Soulaine, M. A. AlSaud, A. Kovscek, and H. Tchelepi. Particle velocimetry analysis of immiscible two-phase flow in micromodels. *Advances in Water Resources*, 95:199–211, Sept. 2016.

[290] T. Roosevelt and A. C. Fowler. A mathematical model for water and nutrient uptake by plant root systems. *Journal of Theoretical Biology*, 228(2):173–184, May 2004.

[291] T. Roosevelt, S. D. Keyes, K. R. Daly, A. Carminati, W. Otten, D. Vetterlein, and S. Peth. Challenges in imaging and predictive modeling of rhizosphere processes. *Plant and Soil*, 407(1):9–38, Oct. 2016.
[292] R. Rosenzweig, U. Shavit, and A. Furman. Water Retention Curves of Biofilm-Affected Soils using Xanthan as an Analogue. *Soil Science Society of America Journal*, 76(1):61-69, 2012. _eprint: https://onlinelibrary.wiley.com/doi/pdf/10.2136/sssaj2011.0155.

[293] P. J. Ross and K. R. J. Smettem. A Simple Treatment of Physical Nonequilibrium Water Flow in Soils. *Soil Science Society of America Journal*, 64(6):1926-1930, 2000. _eprint: https://onlinelibrary.wiley.com/doi/pdf/10.2136/sssaj2000.6461926x.

[294] S. Ruiz, D. McKay Fletcher, K. Williams, and T. Roosevelt. Plant-Soil Modelling. volume 4, pages 127-198. Mar. 2021.

[295] D. Russo. Determining soil hydraulic properties by parameter estimation: On the selection of a model for the hydraulic properties. *Water Resources Research*, 24(3):453-459, 1988. _eprint: https://onlinelibrary.wiley.com/doi/pdf/10.1029/WR024i003p00453.

[296] M. Rötzer, A. Prechtel, and N. Ray. Pore scale modeling of the mutual influence of roots and soil aggregation in the rhizosphere. *Frontiers in Soil Science*, 3, 2023. Publisher: Frontiers Media S.A.

[297] S. Sammartino, A.-S. Lissy, C. Bogner, R. Van Den Bogaert, Y. Capowiez, S. Ruy, and S. Cornu. Identifying the Functional Macropore Network Related to Preferential Flow in Structured Soils. *Vadose Zone Journal*, 14(10):vzj2015.05.0070, Oct. 2015.

[298] A. Schnepf, A. Carminati, M. A. Ahmed, M. Ani, P. Benard, J. Bentz, M. Bonkowski, M. Knott, D. Diehl, P. Duddek, E. Kröner, M. Javaux, M. Landl, E. Lehndorff, E. Lippold, A. Lieu, C. W. Mueller, E. Oburger, W. Otten, X. Portell, M. Phalempin, A. Prechtel, R. Schulz, J. Vanderborght, and D. Vetterlein. Linking rhizosphere processes across scales: Opinion. *Plant and Soil*, 478(1):5-42, Sept. 2022.

[299] A. Schnepf, D. Leitner, M. Landl, G. Lobet, T. H. Mai, S. Morandage, C. Sheng, M. Zörner, J. Vanderborght, and H. Vereecken. CRootBox: a structural-functional modelling framework for root systems. *Annals of Botany*, 121(5):1033-1053, Apr. 2018.

[300] N. Schröder, N. Lazarovitch, J. Vanderborght, H. Vereecken, and M. Javaux. Linking transpiration reduction to rhizosphere salinity using a 3D coupled soil-plant model. *Plant and Soil*, 377(1):277-293, Apr. 2014.

[301] N. Schwartz, A. Carminati, and M. Javaux. The impact of mucilage on root water uptake—A numerical study. *Water Resources Research*, 52(1):264-277, 2016. _eprint: https://onlinelibrary.wiley.com/doi/pdf/10.1002/2015WR018150.

[302] M. R. Segura-Campos, N. Ciau-Solís, G. Rosado-Rubio, L. Chel-Guerrero, and D. Betancur-Ancona. Chemical and Functional Properties of Chia Seed (*Salvia hispanica* L.) Gum. *International Journal of Food Science*, 2014:e241053, Mar. 2014. Publisher: Hindawi.

[303] J. S. Selker and M. H. Schroth. Evaluation of hydrodynamic scaling in porous media using finger dimensions. *Water Resources Research*, 34(8):1935-1940, 1998. _eprint: https://onlinelibrary.wiley.com/doi/pdf/10.1029/98WR00625.
[304] P. Seppecher. Second-gradient theory: Application to Cahn-Hilliard fluids. In G. A. Maugin, R. Drouot, and F. Sidoroff, editors, *Continuum Thermomechanics: The Art and Science of Modelling Material Behaviour*, Solid Mechanics and Its Applications, pages 379–388. Springer Netherlands, Dordrecht, 2002.

[305] L. F. Shampine and M. K. Gordon. *Computer solution of ordinary differential equations: the initial value problem*. W.H. Freeman, San Francisco, 1975. OCLC: 1093545.

[306] S. Shi, E. E. Nuccio, Z. J. Shi, Z. He, J. Zhou, and M. K. Firestone. The interconnected rhizosphere: High network complexity dominates rhizosphere assemblages. *Ecology Letters*, 19(8):926–936, 2016. _eprint: https://onlinelibrary.wiley.com/doi/pdf/10.1111/ele.12630.

[307] O. T. Sililo and J. H. Tellam. Fingering in Unsaturated Zone Flow: A Qualitative Review with Laboratory Experiments on Heterogeneous Systems. *Ground Water*, 38(6):864–871, Nov. 2000.

[308] R. Sinha, F. B. Fritschi, S. I. Zandalinas, and R. Mittler. The impact of stress combination on reproductive processes in crops. *Plant Science*, 311:111007, Oct. 2021.

[309] S. Solomon, D. Qin, M. Manning, Z. Chen, M. Marquis, K. Avery, M. Tignor, and H. Miller. Climate Change 2007: The Physical Science Basis. Working Group I Contribution to the Fourth Assessment Report of the IPCC. volume 1. Jan. 2007.

[310] C. Soulaine, J. Maes, and S. Roman. Computational Microfluidics for Geosciences. *Frontiers in Water*, 3, 2021.

[311] C. Soulaine, S. Roman, A. Kovscek, and H. A. Tchelepi. Pore-scale modelling of multiphase reactive flow: application to mineral dissolution with production of. *Journal of Fluid Mechanics*, 855:616–645, Nov. 2018. Publisher: Cambridge University Press.

[312] G. P. Sparling. Soil microbial biomass, activity and nutrient cycling as indicators of soil health. *Biological indicators of soil health.*, pages 97–119, 1997. Publisher: CAB INTERNATIONAL.

[313] L. Stingaciu, H. Schulz, A. Pohlmeier, S. Behnke, H. Zilken, M. Javaux, and H. Vereecken. In Situ Root System Architecture Extraction from Magnetic Resonance Imaging for Water Uptake Modeling. *Vadose Zone Journal*, 12(1):vzj2012.0019, Feb. 2013.

[314] P. Stoodley, R. Cargo, C. J. Rupp, S. Wilson, and I. Klapper. Biofilm material properties as related to shear-induced deformation and detachment phenomena. *Journal of Industrial Microbiology and Biotechnology*, 29(6):361–367, Dec. 2002.

[315] M. C. Sukop and D. Or. Lattice Boltzmann method for modeling liquid-vapor interface configurations in porous media. *Water Resources Research*, 40(1), 2004. _eprint: https://onlinelibrary.wiley.com/doi/pdf/10.1029/2003WR002333.

[316] I. W. Sutherland. The biofilm matrix – an immobilized but dynamic microbial environment. *Trends in Microbiology*, 9(5):222–227, May 2001. Publisher: Elsevier.

[317] G. Söderlind and L. Wang. Adaptive time-stepping and computational stability. *Journal of Computational and Applied Mathematics*, 185(2):225–243, Jan. 2006.