Analytically Supported Numerical Modeling of Horizontal and Radial Collector Wells

Today, drinking water supply is a strategic task all over the world. One of the most important drinking water sources in Central Europe, especially in Hungary, is riverbank filtration (RBF), which is based on abstraction wells located close to rivers (Grischek et al. 2012). RBF offers various advantages, such as protection against chemicals and pathogens, buffering capacity during contaminant transport and extreme hydrologic events, and a pretreatment step during water treatment (Sandhu and Grischek 2012). As a result of the bank filtration (BF) process particles, biodegradable organic compounds, pathogens, and nitrate are completely or partly removed from the surface water (Dragon et al. 2018). At the end of the BF process the quality of raw water often meets the drinking water quality standards of different sites (Hiscock and Grischek 2002). Through RBF, many communities along riverbanks pump drinking water from alluvial aquifers using different types of wells (Ray et al. 2002). In addition to the quality issues of water production, optimization of the quantity of water production should be taken into account. Properly selected and designed water production facilities can contribute to achieving higher yield with one construction. Such facilities may include horizontal wells and radial collector (also referred to as Ranney or star) wells with horizontal arms, which can provide high flow rate. RBF-based drinking water production and horizontal collector wells are used by many cities all over the world such as Warsaw, London, Budapest, and Belgrade (Babac and Babac 2009; Hunt 2003). It is important to know the hydraulics of these complex systems in order to optimize safe and sustainable drinking water production. The aim of our work is to better understand the hydraulic conditions of radial collector wells and to calculate the water levels and flow rates using different modeling methods.

One of the main objectives of hydrogeological modeling is to approximate the hydrogeological observations with the help of different modeling methods (Szucs et al. 2006). These methods can help maintain sustainable water extraction (Buday et al. 2015), especially in drinking water supply.

In this study, the steady-state and unsteady flow to wells with horizontal screens is analyzed. Nonleaky (confined or unconfined) homogeneous and isotropic aquifers are considered. Linear systems with no wellbore storage are assumed in comparative analyses. The nonlinear processes in the screen pipes (axial turbulent friction loss, kinetic energy variation) and the turbulent screen orifice loss are briefly overviewed in the next section. The saturated thickness of the aquifer is assumed as a time-invariant parameter, that is, the thickness reduction effect by depletion of the free surface at the top is ignored. The thickness, hydraulic conductivity, specific storage, and specific yield constitute the hydrogeologic parameters. The number, depth, direction, drilling radius, length, and distance of screened intervals characterize the simulated wells. Zero, one, and four fixed-head boundaries are considered in the presented simulation examples.

In the three presented analyses, infinite, semi-infinite, and square-shaped bounded aquifers are considered. Three solutions (Renard and Dupuy 1991; Morozov 2018; Milojevic 1963) are used to verify the CW software (Szekely 2021), which is utilized in all three examples. Zhan and Zlotnik (2002) and Bakker et al. (2005) provide conceptual and quantitative support for comparative analyses for infinite and bounded aquifers. More methods are used in (1) comparative simulations; and (2) the Appendix demonstrating the technique for unsteady-state drawdown called late time partitioning (LTP).

Szekely (2021) presented an analytically based method for collector wells in homogeneous, vertically anisotropic aquifers and the related CW software developed for steady-state flow conditions. This program is the base tool providing analytical support.

Fig. 1 shows the vertical section of a radial collector well with arms at two different levels. The two-level option is used in the example of the section “Unsteady Flow to the Two-Level Radial Collector Well in an Unconfined Aquifer.” The five-arm example in the section “Steady-State Flow to the Five-Arm Radial Collector Well in a Square-Shaped Aquifer” demonstrates the case where all the arms are drilled at the same level. Finally, the horizontal well of the section “Unsteady Flow to Horizontal Well in an Infinite Unconfined Aquifer” may be approximated as a two-arm radial collector well with no caisson in the center.

For the screened sections of radially located horizontal arms, the following geometric parameters are considered: number, direction, depth $d$, length $l$, drilling radius $r_d$, and distance of the midpoint of the screen to the center of the caisson $c$. The caisson is considered as part of the aquifer of thickness $b$; therefore, its diameter and impermeable wall have no hydraulic impact on the flow patterns around the well. The CW software is designed and tested for horizontal and radial collector wells under three types of constant-head boundaries generating recharge through (1) the vertical cylindrical surface (circular recharge boundary), (2) the vertical plane (RBF), or (3) the top of the aquifer (well under the river bed). This software is the main tool in accuracy analysis and analytical support of (confirming or improving) the numerical models.

The closed-form solution by Renard and Dupuy [1991, Eq. (5)] combines two analytical methods. It was used to confirm the CW simulated well drawdown for square-shaped or circular discharge area with the effect of vertical anisotropy (Szekely 2020). The comparative analysis showed a difference in the third digit only. As shown in the Appendix, this solution in combination with the LTP method provides a reasonable approximation for the late time unsteady drawdown in infinite aquifers.

The semianalytical solution by Morozov [2018, Figs. 2(a and b)] was used to verify the CW software for both the steady-state well drawdown and the flux distribution along the well screen. Circular discharge areas with four and eight radially positioned horizontal wells were assumed and the ratio of the flux variation to the uniform drawdown was calculated. For each well layout, three options with different aquifer thicknesses are evaluated. According to Szekely (2020), the 114 data by Morozov (personal communication) and the corresponding CW simulated values exhibit a small mean absolute relative deviation of 0.0159 and a high 0.9995 correlation coefficient.

Milojevic [1963, Eq. (9)] presented an empirical formula to define the relationship between the pumping rate and the steady-state drawdown in the arms of a single radial collector well. The well operates in the semi-infinite aquifer with constant-head boundary that simulates the effect of RBF. The referenced RBF solution for an isotropic aquifer is based on laboratory experiments and considers uniform head in the arms. The presented comparative analysis by Szekely (2020) is performed with the following parameters: aquifer thickness $b=10\mathrm{m}$, drilling depth $d=5\mathrm{m}$, drilling radius $r_d=0.25\mathrm{m}$, length of fully screened arms $l=25\mathrm{m}$, number of $\mathrm{arms}=4$, distance to the $\mathrm{river}=200\mathrm{m}$. The six listed parameters are in the validity range of the empirical equation published by Milojevic (1963, Table 3). Because no caisson is simulated, the distance required by CW between the origin and the middle of the screened section of the lateral is set to the half-screen’s length of $c=l/2=12.5\mathrm{m}$. The ratio of $Q/k=20$ is applied in the analysis where $Q$ and $k$ are the pumping rate (${\mathrm{m}}^3/\mathrm{day}$) and the hydraulic conductivity ($\mathrm{m}/\mathrm{day}$), respectively. The two methods show a low relative deviation of 0.01125. The comparative analysis was extended to wells with 8 and 12 arms. The higher but still low relative deviations of 0.01987 (eight arms) and 0.05219 (12 arms) are acceptable when comparing results obtained with experimental and theoretical analyses.

The CW software has also been successfully verified for nonlinear hydraulic processes. Following Cyr (2007), the steady-state flow around and in the horizontal well was analyzed as outlined next. A presumably infinite and homogeneous aquifer beneath a water reservoir is assumed. Rough, turbulent flow through the screen orifice and in the screen pipe is considered. The impact of the dynamic pressure by the momentum effect is also taken into account. The improved modeling option called Modified HGS-2 by Cyr (2007) and the relevant CW simulation show appropriate agreement of both the calculated heads and fluxes (Szekely 2020).

The finite-difference (FD) numerical groundwater flow modeling is used in all examples. FLOW software is used for simulating the unsteady and steady-state flow in large aquifers of infinite or semi-infinite extension. The MODFLOW software involving the Revised Multi-Node Well (MNW2) package (Konikow et al. 2009) is used to simulate the steady-state flow to the five-arm radial collector well in a square-shaped aquifer. FLOW uses the point-centered FD scheme, whereas MODFLOW employs the block-centered FD scheme.

The analytical methods and WT and WF software for vertical wells allow for unsteady drawdown analysis of single wells or well fields, respectively. These methods are applied to the replacing vertical wells and provide useful tools to confirm the results obtained for horizontal wells. The short, screened replacing vertical wells are positioned along the horizontal screens of the simulated collector well. They are located in the centers of segments or horizontal line sinks constituting the horizontal screens of the CW model. The well line operates at a discharge rate of the segments and generates independent unsteady drawdown data for analytical support.

The FW method combines the FLOW and WF concepts and simulates the drawdown in the wells of the replacing well line. The FW procedure uses the FLOW numerical modeling environment by calling the routine for well drawdown in the nodes.

The MODFLOW MNW2 and FW software use the option of calculating separate drawdown in the wellbore and for the associated nodes of the FD simulation mesh incorporating the relevant screen section. They allow for taking into account the parameters of screen loss or skin effect caused by change of hydraulic conductivity outside the screen.

In the examples presented in the following subsections, various layouts of the numerical models and pumping wells are applied. Table 1 summarizes the main features of the simulation models.

In our work, we dealt with horizontal and radial collector wells, which are one of the most important well types in RBF systems. We aimed to better understand the hydraulic conditions of this kind of important water production facility. In this study, we used various methods for modeling the hydraulic conditions of horizontal and radial collector wells. The following six points summarize our main conclusions:

In future field studies, the two presented numerical simulation methods can be used for (1) heterogeneous aquifers with irregular, curve-linear boundaries; (2) water table aquifers with drawdown-dependent saturated thickness; and (3) wells with screens operated under the riverbed common in riverbank filtration schemes.

Hantush (1964, Fig. 25) analyzed and visualized the unsteady drawdown around radial collector wells of arms’ length $l$ in infinite water table aquifers. He concluded that at distance more than $5l$ “the equation of drawdown for a collector well of at least two laterals on a line can be given by the Theis formula” (Hantush 1964). This means that the drawdown contour lines become practically circular at and beyond the distance $R_c=5·l$ from the center of the well. For fully penetrating vertical fracture or horizontal trench, the same is demonstrated graphically by Van Tonder et al. (2002).

For single vertical wells operating in infinite confined or unconfined aquifers, the $s–\log t$ plot at late time converges to the straight line and the Theis solution may be approximated by $Q/(4·\pi ·T)·\ln [(2.2458·T·t)/(r^2·S)]$, referred to as Jacob’s formula (Cooper and Jacob 1946). Here symbols $s$, $Q$, $T$, $S$, $r$, and $t$ denote the drawdown, discharge rate, transmissivity, storativity, distance, and time, respectively. An arbitrary recharge boundary $R$ is introduced by adding and subtracting the term $\ln (R)$. In this way the previously referenced single term unsteady drawdown formula can be partitioned or decomposed into the steady-state analytical core $s_1$ and the unsteady part as follows:where

Eq. (2) is valid at $r<R$. The term $s_1$ is the Dupuit equation at assumed constant-head circular recharge boundary of radius $R$, while the second term on the right-hand side of Eq. (1) defines the time-variant drawdown at $R$. Similar LTP drawdown approximation may also be introduced for vertical fractures, horizontal drains, or collector wells. It requires replacement of the Dupuit equation with the relevant steady-state solutions $s(x,y,z,R)$ for the selected abstraction scheme considering circular recharge boundary with an influence radius $R$. The following comparative analysis demonstrates and confirms this procedure.

Gringarten et al. (1974) developed analytical methods for two-dimensional (2D) unsteady drawdown in a confined homogeneous aquifer tapped by a fully penetrating vertical fracture of length $2l$. The options of uniform flux and infinite fracture conductivity or uniform head are considered. Table 1 in Gringarten et al. (1974) provides the dimensionless drawdownin the fracture of infinite conductivity. The dimensionless time is defined as

As the time increases, (1) the gradient stabilizes at $r<R$, (2) the $s–\log t$ graph becomes linear, and (3) the flow patterns converge to the radial flow beyond the distance $R_c$. At this stage, the $s_D$ function for the fracture can be approximated with the LTP method.

Citrini (1951) analyzed the 2D steady-state flow to the radial system of symmetrically located $N$ drains fully penetrating the confined aquifer. The drains of zero width and uniform length $l$ extend from the origin and the drawdown are given as (Citrini 1951; Babac and Babac 2009)

The equivalent well radius $r_e$ is to be calculated as

In the case of single drain $N=2$, the flow patterns are identical to the vertical fracture of infinite conductivity. Another formula was presented by Dvoretskii and Yarmakhov (1998, p. 131). Both solutions yield identical results. Dvoretskii and Yarmakhov (1998) also published a closed-form analytical solution for the steady-state-specific (per unit length) flux along the axes of the drains (Dvoretskii and Yarmakhov 1998).

Renard and Dupuy (1991) presented the closed-form formula for the drawdown in a fully penetrating horizontal well positioned in the center of the symmetric 2D flow domain in the form of

Following Giger (1987), $R^{*}$ denotes either the radius of constant-head circular recharge boundary or the ${(A/\pi )}^{1/2}$ parameter, where $A$ is the area of the presumably square-shaped flow domain with zero drawdown at the boundary.

Table 4 compares the analytical (Gringarten et al. 1974, Table 1), the LTP modified (Citrini 1951; Renard and Dupuy 1991) solutions at the following parameters: $T=0.001{\mathrm{m}}^2/\mathrm{s}$, $S=0.00001$, $Q=0.001{\mathrm{m}}^3/\mathrm{s}$, and $l=10\mathrm{m}$. A value of $R=100\mathrm{m}$ is selected for calculations; this assumed radius of influence yielded the following steady-state $s_1$ data: 4.7677 m (Citrini model) and 4.7637 m (Renard and Dupuy model). The small deviation from the analytical data verifies the applicability of the LTP transformation for this flow problem.

Some or all data, models, or code that support the findings of this study are available from the corresponding author upon reasonable request.

The research was carried out within the 2018-1.2.1-NKP-2018-00011 “Drinking Water: Multidisciplinary Assessment of Secure Supply from the Source to the Consumers” project. We would like to thank the reviewers for their useful comments to improve our manuscript.