3.6. Variable Flux Boundary

The variable flux boundary condition simulates an infinite extent of the regional groundwater system around the model boundary. This boundary condition simulates boundary fluxes that would result if the modeled flow domain extended outward a great distance from the actual boundary. This works by attaching the analytical solution for a semi-infinite linear aquifer to the boundary of the modeled flow domain.

To ensure proper implementation of variable-flux boundary conditions, place the actual model domain boundary far enough from simulated hydraulic stresses such as mining and pumping that their effects at the variable-flux boundary become negligible. The water exchange at the boundary of the modeled domain from variable-flux boundary conditions takes the general form shown below.

(3.14)\[q_{v} = \sum_{i=1}^{n} \left[ C_{Vi}\,\bigl(H_{Vi} - H_i\bigr) + Q'_{Vi} \right] + Q_{V0}\]

Where:

  • \(q_V\) = functional representation of the groundwater discharge across the variable flux boundary [L³T⁻¹],

  • \(C_{Vi}\) = coefficient derived from the analytical solution for a semi-infinite linear aquifer (the equation) [L²T⁻¹],

  • \(H_{Vi}\) = steady-state hydraulic head [L],

  • \(H_i\) = computed hydraulic head for the current time step [L],

  • \(Q'_{Vi}\) = discharge across the boundary in the current time step due to head changes at the variable flux boundary in past time steps [L³T⁻¹],

  • \(Q_{V0}\) = initial discharge across the boundary [L³T⁻¹], and

  • n = number of nodes [L°].

The initial discharge across the variable flux boundary (\(Q_{V0}\)) must be included when calculating discharge across variable-flux boundaries. This value typically comes from the steady-state model, where a constant-head boundary condition replaced the variable-flux boundary condition. Therefore, \(Q_{V0}\) represents the discharge required to achieve the model’s specified initial head distribution.

The coefficient \(C_{Vi}\), initial discharge \(Q_{v0}\), and discharge \(Q′_{Vi}`\) have non-zero values at each node on the variable-flux boundary and zero values elsewhere. The coefficient \(C_{Vi}\), calculated using the equation below, derives from the analytical solution for discharge from a semi-infinite linear aquifer.

(3.15)\[C_{Vi} = \frac{2}{\Delta t}\,\left[ \frac{K\,B\,W}{\left(\pi K / S_S\right)^{1/2}} \right]\,\Delta t^{1/2}\]

Where:

  • K = hydraulic conductivity of the element associated with node i [LT⁻¹],

  • B = thickness of the face normal to the direction of variable flux associated with node i [L],

  • W = width of the face normal to the direction of variable flux associated with node i [L],

  • \(S_s\) = specific storage of the element associated with node i [L⁻¹], and

  • t = time [T].

The discharge at x = 0 for a semi-infinite linear aquifer extending from x = 0 to x = ∞ appears below (Carslaw and Jaeger, 1959 [3]).

(3.16)\[q_{0} = -\left[ \frac{K\,B\,W}{\sqrt{\pi K t / S_S}} \right]\,\Delta h\]

Where:

  • \(q_o\) = instantaneous discharge rate at x = 0 [L³T⁻¹], and

  • \(\Delta h\) = step change in the boundary head at x = 0 in the semi-infinite aquifer [L].

A positive qo value represents discharge from the semi-infinite aquifer into the modeled domain. The average discharge over a specified time period comes from integrating over the period and dividing by the period duration. The integral appears below and evaluates in the following equation.

(3.17)\[Q_{0} = -\frac{1}{\Delta t} \int_{t_{1}}^{t_{2}} \left[ \frac{K\,B\,W}{\sqrt{\pi K t / S_S}} \right] \,\Delta h \, dt\]
(3.18)\[Q_{0} = -\frac{2}{\Delta t}\,\left[ \frac{K\,B\,W}{\left(\pi K / S_S\right)^{1/2}} \right]\,\Delta h\,\left(t_{2}^{1/2} - t_{1}^{1/2}\right)\]

Where:

  • \(Q_0\) = average discharge during period \(t_1\) to \(t_2\) [L³T⁻¹],

  • \(t_1\) = time to start of a step head change (Δh) [T], and

  • \(t_2\) = time to end of a step head change (Δh) [T].

The equation rewrites as shown below to demonstrate its consistency with the first term in the equation.

(3.19)\[Q_{0} = -\,C_{Vi}\,\Delta h_{i}\]

Where:

\[-\Delta h_{i} = H_{Vi} - H_i\]

The second term, \(Q′_{Vi}\), accounts for discharge from head changes in all past time steps where the variable flux boundary condition applied. The equation for \(Q′_{Vi}\) below derives from the equation to account for discharge from head differences in previous time steps.

(3.20)\[Q'_{v} = \sum_{j=1}^{m-1} \frac{2}{\Delta t_{m}} \left[ \frac{K\,B\,W}{\left(\pi K / S_S\right)^{1/2}} \right] \,\Delta h_{j}\, \left[ \left(t_{m}-t_{j-1}\right)^{1/2} - \left(t_{m-1}-t_{j-1}\right)^{1/2} \right]\]

Where:

  • \(j\) = the index for a time step [L°], and

  • \(m\) = the index for the current time step [L°].