*2.4. Stage 2*

Once the wetting front reaches the interface between layers 1 and 2, the low permeability of layer 2 means that moisture begins to accumulate about the interface. The additional moisture creates a moisture gradient both below and above the interface in different directions. These gradients allow more or less moisture than would occur by gravity alone at that moisture content and allows the flux and soil suction to be continuous at the interface. A mass balance argumen<sup>t</sup> means that the difference in flux at the upper and lower ends of the segmen<sup>t</sup> is equal to the rate of accumulation in that segment. The time for saturation to occur is given by:

$$dt\_2 = \int (\partial - \partial\_b) dz \langle (A\_n - A\_o), \tag{7}$$

where θ*b* is the background water content. For layer 1, this background value is the new moisture content and for layer 2, this is the old moisture content. We will return to the calculation of integral (7) later.

## *2.5. Stage 3*

The dimensionless equations describing the mass balance of the perched layer are:

$$
\beta \frac{\partial \mathfrak{h}}{\partial t} = B(\frac{\partial}{\partial \mathfrak{x}}(h \frac{\partial \mathfrak{h}}{\partial \mathfrak{x}})) + A - q \ 0 < \mathfrak{x} < 1 \tag{8}
$$

$$
\beta \frac{\partial \mathbf{h}}{\partial t} = B(\frac{\partial}{\partial \mathbf{x}}(h \frac{\partial \mathbf{h}}{\partial \mathbf{x}})) - q \text{ 1} < \mathbf{x} < \mathbf{x}\_1 \tag{9}
$$

where

$$B = \mathbf{K}\_{\mathbb{S}1}{}^h l\_2 \, ^2/(\mathbf{K}\_{\mathbb{S}2}{}^v \mathbf{x}\_0 \, ^2) \tag{10}$$

$$A = I A\_{\rm u} \%\_{s2} \text{ v} \tag{11}$$

$$\beta = (\theta\_s^{-1} - \theta\_n^{-1}) \langle (\theta\_s^{-2} - \theta\_o^{-2}) \rangle \tag{12}$$

*h* is the head of the perched water table, *x* = 1 represents the edge of the irrigation field; *x* = *x*1 is the edge of the wetted zone outside of the irrigation field; x0 is the half-width of the irrigated field; and *q* is the vertical flux into layer 2. Continuity in *h* and the flux of water (and therefore gradient in *h*) is assumed to occur at *x* = 1. At *x* = *x*1*, h* is zero.

For the above equations, *A* is a dimensionless parameter that reflects the degree of perching, with perching not occurring for smaller *A* and interception of the perched head with the upper boundary condition for larger *A*. The dimensionless parameter, *B,* reflects the degree to which lateral movement occurs. As *B* approaches zero, there is no lateral movement and the system behaves as a 1D system. For very large *B*, the perched layer spreads thinly across the impeding layer. We also shall assume for this stage, that the head of the perched layer is lower than the upper boundary, i.e.,

$$h \le l\_1 / l\_2 \tag{13}$$

Darcy's Law across the saturated zone implies that:

$$q = 1 + h/z\_{\text{sat}} \tag{14}$$

where *zsat* is the depth of the saturation front. This assumes that the main hydraulic impedance is in the second layer. Under the wetting front model,

$$dz\_{w\theta}dt = q\tag{15}$$

*Water* **2020**, *12*, 944

This assumes that the new flux is much greater than the old flux. During Stage 3, the wetting front has not reached the base of layer 2.

In addition to the equations above, some further assumptions are added in order to estimate recharge:

1. We shall ignore the effect of the ponded head outside the irrigated field on the infiltration into the impeding layer, i.e.,

$$
\eta = dz\_{w\eta} dt = 1,\ 1 < \mathbf{x} < \mathbf{x}\_1 \tag{16}
$$

2. We shall assume quasi-steady-state Depuit–Forchheimer equations for this area, which leads to the following equations:

$$h = (1 + \text{sqrt}(B)h\_1 - \mathbf{x}) / \text{sqrt}(B), \quad 1 < \mathbf{x} < \mathbf{x}\_1 \tag{17}$$

$$\mathbf{x}\_1 = h\_1 \text{sqrt}(\mathbf{B}) + \mathbf{1}\_\text{.} \text{ and} \tag{18}$$

$$Q = -Bl\frac{\partial h}{\partial \mathbf{x}} = h\_1 \text{sqrt}(B) \tag{19}$$

where *Q* is the non-dimensional lateral flux at *x* = 1 and *h*1 is the head of the perched layer at *x* = 1.

3. We shall assume that the head is constant across the irrigated field. Combining Equations (8), (9), (18) and (19) gives:

$$
\beta dh\_1 / dt = A - q - h\_1 \text{sptt(B)}.\ 0 < x < x\_1\tag{20}
$$

4. We shall assume in early times of ponding that the lateral movement is small, and processes are vertical. We shall also assume that the separation between saturation fronts and wetting fronts is constant. By defining the dimensionless parameter:

$$
\alpha = h\sharp z\_{\text{sat}\_\*} \tag{21}
$$

we find that *zsat* and *h* increase linearly:

$$z\_{\text{sat}} = (1+\alpha)t \tag{22}$$

$$h = \alpha (1 + \alpha) t \tag{23}$$

$$\alpha = (-(1+\beta) + sqrt((1+\beta)^2 + 4(A-1)\beta)/(2\beta), \text{and} \tag{24}$$

$$q = 1 + \alpha \tag{25}$$

Equation (25) indicates a flux greater than the free drainage flux through a saturated clay layer (*q* = 1).

Stage 3 finishes when the wetting front reaches the bottom of layer 2. To estimate when this occurs, it is necessary to estimate the thickness of the zone, Φ*,* between the wetting front and the saturated zone (Figure 3a). If we assume that this zone is in a quasi-steady state, this can be estimated from Darcy equation to give:

$$\Phi = q\wp(q-1) = \int d\psi \ell(q\%, \zeta(\psi) - 1) l\_2\rangle,\tag{26}$$

$$d\Phi = h\_{b2}\ell(q-1)l\_2) + \int d\psi \psi(q\ell(K\_r(\psi)-1)l\_2) \tag{27}$$

where ϕ is a soil hydraulic property and the integral in Equation (26) goes from 0 (saturated) to ψ3, and in Equation (27) from *hb*2 to ψ3, where ψ3 is either (a) the matric potential relating to the pre-development drainage, where the transition zone is entirely within the clay layer or (b) the soil matric suction at the interface with layer 3. We shall assume that ϕ is constant with respect to *q.*

Once the wetting front reaches the base of layer 2, the flux at top of layer 3 increases to *q* = 1 + α and the wetting front begins to move through this layer. This occurs at non-dimensional time after perching begins:

$$t\_3 = (1 - q\rho\alpha)/(1 + \alpha),\tag{28}$$

and the ponded head will be:

$$h\_0 = \alpha - \varphi \tag{29}$$

By knowing the flux during this third stage, it is possible to use the steady-state Darcy's Law to calculate *t*2 in Equation (7).
