Drawdown of urban drain trenches triggering 2-D transient seepage in soil massifs subject to managed aquifer discharge: sandbox experiments, analytical and HYDRUS2D modeling

ABSTRACT

Managed Aquifer Discharge (MAD) through trenches lowers the water table of ‘shallow perched aquifers’ in urban areas (Muscat-Oman). Drawdown of the water table in periodic trenches or a single trench is experimentally studied in a sandbox, and modelled numerically (by HYDRUS2D) and analytically. Evaluation of transient flow fields of pore pressure and hydraulic gradients is necessary for assessing the stability of trench slopes to seepage-induced erosion. A controlled emptying of trench’s water level is examined. An analytical solution is obtained for an unconfined seepage at an early stage of drainage. Numerical simulations are conducted for early and late stages of seepage towards a partially filled trench in soils with capillarity. For a steady-state Darcian flow limit, comparisons of HYDRUS results with the Charny analytical formula for capillarity-free soils are done. The flow rates, isobars, isotachs and streamlines are determined, with assessment of zones of high hydraulic gradients susceptible to seepage erosion.

KEYWORDS:

• Urban drainage by trenches
• sandbox drawdown experiments
• complex potential-conformal mappings
• boundary-value problem and the keldysh-sedov formula
• HYDRUS2D simulations
• isobars-isohumes-isotachs-streamlines

Acknowledgements

This work was supported by SQU, grants IG/VC/WRC/21/01, DR/RG/17. The development program of the Scientific and Educational Mathematical Center of the Volga Federal District, Russia, agreement No. No. 075-02-2021-1393 is acknowledged. Helpful comments by an anonymous referee are highly appreciated.

Disclosure statement

No potential conflict of interest was reported by the author(s).

List on Main Notations

a₀, b =vertical and horizontal sizes of the drained massif (m)

$h^{\ast }=$ piezometric head (m)

$h_{0s}=$ water level in the downstream trench after sudden emptying (m)

k_s = saturated hydraulic conductivity (m/s)

$p^{\ast }$= pressure head (m)

Q = seepage flow rate (m²/s)

$t^{\ast }$= time (s)

$V^{\ast }\left[u^{\ast },v^{\ast }\right]$ = Darcian velocity vector (m/s)

$w^{\ast }=$ complex potential

$\left(x^{\ast },y^{\ast }\right)=$ Cartesian coordinates (m)

$z^{\ast }=$ complex physical coordinate

(${\theta }_r$, ${\theta }_s$, α, n) = the Van Genuchten soil parameters (1,1,1/m,1)

${\varphi }^{\ast }$= velocity potential (m²/s)

${\psi }^{\ast }=$ stream function (m²/s)

$\zeta =\zeta +i\eta$ complex reference variable

Notes

1. The early stage of Experiments R1 and R2 can be modeled in the same manner.

2. Stoker (1957) wrote: “One class of problems is solved by assuming a solution in the form of power series in the

time, which implies that initial motions and motions for a short time only can be determined in general. Nevertheless, some interesting cases can be dealt with, even rather easily, by using the so-called Lagrange representation, rather than the Euler representation which is used otherwise throughout the book. The problem of the breaking of a dam, and, more generally, problems of the collapse of columns of a liquid resting on a rigid horizontal plane can be treated in this way.” In our problem, at early time the pressure head, the piezometric head, and the stream function are harmonic in the rectangle G_z, see Figure 4. Therefore, our complex potential is a holomorphic function in G_z. In Stoker’s dam break problem (see his book, 513–522), his 2D pressure p₀ is a harmonic function of his space variables a and b (Stoker’s Lagrangian coordinates). Stoker’s characteristic function Z₂ (a complexified acceleration) is holomorphic in his half-strip (which is a special case of our rectangle G_z). Stoker (1957) proceeded in the same way as we are doing in this MS. Specifically, he solved a boundary value problem (the Dirichlet one) for Z₂ in his half-strip. For this purpose, he mapped the half-strip onto a reference half-plane by an elementary function. He used the Schwarz-Christoffel formula for this mapping. In our problem, the conformal mapping of G_z onto a reference half-plane (Figure 4b) is determined by an elliptic function, Equation (2)(2) $z=\frac{1}{2\mathrm{K}(\lambda )}\underset{0}{\overset{\zeta }{\int }}\frac{\mathrm{d}\tau }{\sqrt{\left(1-{\tau }^2\right)\left(1-{\lambda }^2{\tau }^2\right)}}-1/2,$(2) . Our boundary-value problem is not as simple as one of Stoker (1957), viz. we have to solve a mixed boundary-value problem, Equation (5)(5) $\begin{array}{rl}& \varphi (\xi )=0,\xi \in (-\mathrm{\infty },-1/\lambda )\cup (1,\mathrm{\infty }),\\ & \varphi (\xi )=a_0-y(\xi ),\xi \in (-1/\lambda ,-\mu ),\\ & \varphi (\xi )=h_d,\xi \in (-\mu ,-1),\\ & \psi (\xi )=0,\xi \in (-1,1).\end{array}$(5) . That results in an integral solution (6) with the kernel of the integral gien by Equations (7)(7) $g(s)=a_0+\frac{1}{2K(\lambda )}\underset{-1}{\overset{s}{\int }}\frac{\mathrm{d}\tau }{\sqrt{({\tau }^2-1)(1-{\lambda }^2{\tau }^2)}}=a_0-\frac{1}{2K(\lambda )}F\left(arcsin\sqrt{\frac{1-1/s^2}{1-{\lambda }^2}},\sqrt{1-{\lambda }^2}\right)$(7) -(8(8) $g\left(s\right)\equiv h_d,s\in \left(-\mu ,-1\right).$(8) ). Unfortunately, our solution, Equation (6)(6) $w\left(\zeta \right)=\frac{\sqrt{{\zeta }^2-1}}{\pi i}\underset{-1/\lambda }{\overset{-1}{\int }}\frac{g\left(\tau \right)}{\sqrt{\left({\tau }^2-1\right)}}\frac{\mathrm{d}\tau }{\tau -\zeta },$(6) , can not be expressed in elementary functions, as Stoker’s solution, i.e. we have to tackle singular integrals in Equation (6)(6) $w\left(\zeta \right)=\frac{\sqrt{{\zeta }^2-1}}{\pi i}\underset{-1/\lambda }{\overset{-1}{\int }}\frac{g\left(\tau \right)}{\sqrt{\left({\tau }^2-1\right)}}\frac{\mathrm{d}\tau }{\tau -\zeta },$(6) .

3. The vertical axis in HYDRUS is z and default notation of the pressure head is different from adopted in the paper.

4. All our HYDRUS2D project files are available upon request.

5. Best fitting of the dyad (α, n) can be done as in Al-Mayahi et al. (2020), viz. by experimental determination of the water retention curve of the sandbox filling and sensitivity analysis.

6. We used the option ‘Boundary fluxes’ in HYDRUS ‘Results’ at t = 24 hours, which gives practically identical values to ones for t = 18 hours that ensures a steady state seepage regime. We evaluated the exfiltration flux through the seepage face Q_nCsCv = 8 cm²/hour. Charny’s formula can not evaluate separately the fluxes through the seepage face C_sC_v and constant head segment AC_s in Figure 8.

Additional information

Funding

This work was supported by the Sultan Qaboos University [DR/RG/17,IG/VC/WRC/21/01].