9.2. A well with a periodic discharge

The discharge of a well may vary through time periodically. For example, in an aquifer storage and recovery (ASR) system, water is injected into the aquifer during times of abundant water availability and pumped out of the aquifer during times of water scarcity. For an ASR well, the discharge varies around zero, for example sinusoidally as

(9.23)\[ Q(t) = A \sin(2 \pi t / \tau) \]

where \(A\) and \(\tau\) are the amplitude and period of the sinusoidal discharge. The boundary condition at the well is approximated as

(9.24)\[ \lim_{r \to 0} -2\pi r Q_r = A\sin(2\pi t / \tau) \]

while the boundary condition at infinity is ((9.8)). There is no initial condition: The well has been pumping with the sinusoidal discharge (Eq. (9.23)) for a long time. Flow in the aquifer is governed by the diffusion equation ((9.5)).
      The solution for the stated problem is obtained using separation of variables: The solution is written as a function of the radial distance \(r\) multiplied by a function of the time \(t\). The analysis is easiest conducted in terms of complex variables. The head solution corresponding to boundary condition ((9.23)) is the imaginary part of the complex solution (Bruggeman, 1999):

(9.25)\[ h = -\frac{A}{2 \pi T} \Im\left[\text{K}_0\left(\frac{r}{\mu}\sqrt{i}\right)\text{e}^{2\pi i t/\tau}\right] \]

where the characteristic length \(\mu\) is defined as

(9.26)\[ \mu = \sqrt{\frac{T}{S}\frac{\tau}{2\pi}} \]

Note that \(\mu\) has the dimension of length. The corresponding radial component of the discharge vector is

(9.27)\[ Q_r = -\frac{A}{2 \pi} \Im\left[ \frac{\sqrt{i}}{\mu} \text{K}_1\left(\frac{r}{\mu}\sqrt{i}\right)\text{e}^{2\pi i t/\tau}\right] \]

The sinusoidal head variation is largest at the well and decreases with distance from the well. Discharge variations with a short period damp out at a shorter distance from the well than variations with a longer period. It is noted that taking the real part of the complex function in Eq. ((9.25)) results in the solution for the case that the sine function in Eq. ((9.23)) is replaced by a cosine function.
      In the example below, the solution is plotted for a sinusoidal well discharge with a period of 1 year in an unconfined aquifer where the transmissivity is again approximated as constant. The function kv of the scipy.special package is used to compute the modified Bessel functions of complex argument, because the functions k0 and k1 only accept real arguments. In the left graph of Figure 9.4, the head variation is plotted at the well and at a distance \(r=2\mu\). In the right graph of Figure 9.4, the total flow in the aquifer \(Q_\text{tot}=2\pi rQ_r\) is plotted at the well and at \(r = 6\mu\). Both the amplitude of the head and of the total discharge decrease with the distance to the well, but the amplitude of the head decreases much faster. The amplitude of the head decreases to less than 5% of the amplitude of the well at \(r = 2 \mu\), while the amplitude of the total flow decreases to less than \(0.05A\) at \(r = 6 \mu\).

# parameters
T = 100 # transmissivity, m^2/d
S = 1e-1 # storage coefficient, -
tau = 365 # period, d
A = 200 # amplitude of discharge function, m^3/d
rw = 0.3 # radius of well, m
mu = np.sqrt(T * tau / (S * 2 * np.pi))
print(f'characteristic length mu: {mu:.2f} m')

characteristic length mu: 241.02 m

# solution
from scipy.special import kv
def hsinwell(r, t):
    hcomplex = -A / (2 * np.pi * T) * kv(0, r * np.sqrt(1j) / mu) * \
               np.exp(2 * np.pi * 1j * t / tau)
    return hcomplex.imag

def Qrsinwell(r, t):
    Qrcomplex = -A / (2 * np.pi) * np.sqrt(1j) / mu * \
                kv(1, r * np.sqrt(1j) / mu) * np.exp(2 * np.pi * 1j * t / tau)
    return Qrcomplex.imag

# basic plot
t = np.linspace(0, 2 * tau, 100)
plt.subplot(121)
plt.plot(t, hsinwell(rw, t), label='$r=r_w$')
plt.plot(t, hsinwell(2 * mu, t), label='$r=2\\mu$')
plt.legend()
plt.subplot(122)
plt.plot(t, Qrsinwell(rw, t) * 2 * np.pi * rw, label='$r=r_w$')
plt.plot(t, Qrsinwell(6 * mu, t) * 2 * np.pi * 6 * mu, 'C2', label='$r=6\\mu$')
plt.legend();

Figure 9.4 Sinusoidal variations vs. time at different distances for a periodic well for the head (left) and \(Q_\text{tot}\) (right).

      As an illustration, an animation of the head surface of a periodic well is shown in Figure 9.5

Exercise 9.7

Use numerical derivatives to verify that \(Q_r\) as programmed is indeed \(Q_r=-T\partial h/\partial r\) for several values of \(r\) and \(t\).

Solution

They are the same for several decimals

Once Loop Reflect

Figure 9.5 The head surface of a periodic well.