Neuman and Witherspoon Solution for Confined Two-Aquifer Systems with Leakage

A mathematical solution by Neuman and Witherspoon (1969) is useful for determining the hydraulic properties of leaky confined aquifer systems (transmissivity and storage coefficient of pumped aquifer; vertical hydraulic conductivity and storage coefficient of aquitard; and transmissivity and storage coefficient of unpumped aquifer). Analysis involves matching solution type curves to drawdown data collected during a pumping test. Observation wells may be screened in the pumped aquifer, unpumped aquifer or aquitard.

Unlike the Hantush and Jacob (1955) and Hantush (1960) solutions for leaky confined aquifers, the Neuman and Witherspoon (1969) solution accounts for drawdown in the unpumped aquifer. Neuman and Witherspoon showed that the assumption of no drawdown in the unpumped aquifer (i.e., a constant head boundary condition) can lead to significant errors when estimating the hydraulic properties of the pumped aquifer.

You are not restricted to constant-rate tests with the Neuman and Witherspoon solution. AQTESOLV incorporates the principle of superposition in time to simulate variable-rate and recovery tests with this method.

Assumptions

• aquifer has infinite areal extent
• aquifer is homogeneous, isotropic and of uniform thickness
• aquifer potentiometric surface is initially horizontal
• control well is fully penetrating
• flow to control well is horizontal
• aquifer is leaky
• flow is unsteady
• water is released instantaneously from storage with decline of hydraulic head
• diameter of control well is very small so that storage in the well can be neglected
• aquitard has infinite areal extent, uniform vertical hydraulic conductivity and storage coefficient, and uniform thickness
• flow in the aquitard is vertical

Equations

Neuman and Witherspoon (1969) derived the following equations for drawdown in a confined two-aquifer system:

$$s_1=\frac{Q}{4\pi T}{\int }_0^{\infty }\left(1-e^{-y^2{\overline{t}}_{D1}}\right)·$$ $$\left\{\left[1+G\left(y\right)\right]J_0\left({\omega }_1\left(y\right)\right)+\left[1-G\left(y\right)\right]J_0\left({\omega }_2\left(y\right)\right)\right\}\frac{dy}{y}\text{(1)}$$ $$s_1^{\prime }=\frac{Q}{4\pi T}\frac{2}{\pi }\sum _0^{\infty }\frac{1}{n}\sin \frac{n\pi z}{b_1^{\prime }}{\int }_0^{\infty }\left[1-e^{·n^2{\pi }^2{\overline{t}}_{D1}}+\frac{e^{-n^2{\pi }^2{\overline{t}}_{D1}}-e^{-y^2{\overline{t}}_{D1}}}{1-y^2/\left(n^2{\pi }^2\right)}\right]·$$ $$\{\left[\frac{2{\left(r/B_{21}\right)}^2{\left(-1\right)}^{n}y}{F\left(y\right)\sin y}-G\left(y\right)-1\right]J_0\left({\omega }_1\left(y\right)\right)-$$ $$\left[\frac{2{\left(r/B_{21}\right)}^2{\left(-1\right)}^{n}y}{F\left(y\right)\sin y}-G\left(y\right)+1\right]J_0\left({\omega }_2\left(y\right)\right)\}\frac{dy}{y}\text{(2)}$$ $$s_2=\frac{Q}{4\pi T}{\int }_0^{\infty }\left(1-e^{-y^2{\overline{t}}_{D1}}\right)\frac{2{\left(r/B_{21}\right)}^2}{F\left(y\right)}·$$ $$\left[J_0\left({\omega }_1\left(y\right)\right)-J_0\left({\omega }_2\left(y\right)\right)\right]\frac{dy}{\sin y}\text{(3)}$$ $$t_{Di}=\frac{T_{i}t}{S_i{r}^2}\text{(4)}$$ $${\overline{t}}_{D1}=t_{D1}{\left(r/B_{11}\right)}^4/{\left(4{\beta }_{11}\right)}^2\text{(5)}$$ $$r/B_{ij}=r\sqrt{\frac{K_j^{\prime }}{T_i{b}_j^{\prime }}}\text{(6)}$$ $${\beta }_{ij}=\frac{r}{4}\sqrt{\frac{K_j^{\prime }{S}_j^{\prime }}{b_j^{\prime }{T}_i{S}_i}}\text{(7)}$$ $$G\left(y\right)=M\left(y\right)/F\left(y\right)\text{(8)}$$ $${\omega }_1^2\left(y\right)=\frac{1}{2}\left[N\left(y\right)+F\left(y\right)\right]\text{(9)}$$ $${\omega }_2^2\left(y\right)=\frac{1}{2}\left[N\left(y\right)-F\left(y\right)\right]\text{(10)}$$ $$F^2\left(y\right)=M^2\left(y\right)+{\left[\frac{2\left(r/B_{11}\right)\left(r/B_{21}\right)y}{\sin y}\right]}^2\text{(11)}$$ $$M\left(y\right)=\left[\frac{{\left(r/B_{11}\right)}^4}{{\left(4{\beta }_{11}\right)}^2}-\frac{{\left(r/B_{21}\right)}^4}{{\left(4{\beta }_{21}\right)}^2}\right]y^2-$$ $$\left[{\left(\frac{r}{B_{11}}\right)}^2-{\left(\frac{r}{B_{21}}\right)}^2\right]y\cot y\text{(12)}$$ $$N\left(y\right)=\left[\frac{{\left(r/B_{11}\right)}^4}{{\left(4{\beta }_{11}\right)}^2}+\frac{{\left(r/B_{21}\right)}^4}{{\left(4{\beta }_{21}\right)}^2}\right]y^2-\left[{\left(\frac{r}{B_{11}}\right)}^2+{\left(\frac{r}{B_{21}}\right)}^2\right]y\cot y\text{(13)}$$

where

• $b_i$ is thickness of aquifer i [L]
• $b_1^{\prime }$ is thickness of aquitard [L]
• $J_0$ is Bessel function of first kind, zero order
• $K_1^{\prime }$ is vertical hydraulic conductivity of aquitard [L/T]
• $Q$ is pumping rate [L³/T]
• $r$ is radial distance from pumping well to observation well [L]
• $s_1$ is drawdown in pumped aquifer [L]
• $s_1^{\prime }$ is drawdown in aquitard [L]
• $s_2$ is drawdown in unpumped aquifer [L]
• $S$, $S_1$ is storativity in pumped aquifer [dimensionless]
• $S_2$ is storativity in unpumped aquifer [dimensionless]
• $S_1^{\prime }$ is storativity of aquitard [dimensionless]
• $t$ is elapsed time since start of pumping [T]
• $T$, $T_1$ is transmissivity in pumped aquifer [L²/T]
• $T_2$ is transmissivity in unpumped aquifer [L²/T]

Expressions for drawdown in the aquitard for fully or partially penetrating wells can be found by averaging the drawdown over the length of the well screen.

Data Requirements

• pumping and observation well locations
• pumping rate(s)
• observation well measurements (time and displacement)
• partial penetration depths (optional, aquitard only)
• aquitard thickness, b' (for partially penetrating wells)

Solution Options

• variable pumping rates
• multiple pumping wells
• multiple observation wells
• boundaries

Estimated Parameters

• $T$ (transmissivity in pumped aquifer = T₁)
• $S$ (storativity in pumped aquifer = S₁)
• $r/B$ (leakage parameter = r/B₁₁)
• $\beta$ (leakage parameter, = β₁₁)
• $T2$ (transmissivity in unpumped aquifer = T₂)
• $S2$ (storativity in unpumped aquifer = S₂)

The Report view shows aquitard properties, $K\prime /b\prime$ and $K\prime$, computed from the leakage parameter, $r/B$.

Curve Matching Tips

• Use the Cooper and Jacob (1946) solution to obtain preliminary estimates of aquifer properties.
• Choose Match>Visual to perform visual curve matching using the procedure for type curve solutions.
• Use active type curves for more effective visual matching with variable-rate pumping tests.
• Select values of $r/B$ and $\beta$ from the Family and Curve drop-down lists on the toolbar.
• Use parameter tweaking to perform visual curve matching and sensitivity analysis.
• Perform visual curve matching prior to automatic estimation to obtain reasonable starting values for the aquifer properties.

References

Neuman, S.P. and P.A. Witherspoon, 1969. Theory of flow in a confined two aquifer system, Water Resources Research, vol. 5, no. 4, pp. 803-816.