Simulation of Groundwater Flow: Computational Notes & Analysis

Detailed Guide: Simulation of Groundwater Flow

The practice of Simulation of Groundwater Flow combines principles of fluid mechanics, hydrogeology, and computational mathematics. Subsurface flow systems are intrinsically complex due to geological variations, boundary constraints, and transient recharge mechanisms. Numerical simulation allows civil and water resource engineers to evaluate aquifer behavior, predict water table drawdowns, prevent slope instabilities, and design seepage control systems for dams and dikes.

1. Fundamental Equations and Flow-Net Analysis

The foundation for any Simulation of Groundwater Flow rests upon Darcy’s Law, which states that the volumetric flow rate ($Q$) through a porous medium is directly proportional to the hydraulic gradient ($i = -\frac{dh}{dx}$) and the cross-sectional area ($A$):

Combining Darcy’s Law with the steady-state continuity equation for an incompressible fluid yields Laplace’s Equation in two dimensions for isotropic, homogeneous porous media:

Graphical solutions to Laplace’s equation yield a Flow-Net, which is an orthogonal network consisting of two sets of curves:

• Streamlines (Flow Lines): Paths along which individual particles of water move through the soil. No flow occurs across streamlines.
• Equipotential Lines: Contours along which the total hydraulic head $h$ remains constant. Piezometers installed anywhere along an equipotential line will exhibit identical water levels.
• Flow Potential ($\Phi$): Defined as $\Phi = -K h$, representing the potential energy driving the subsurface fluid motion.
• Rules for Flow-Net Construction: In an isotropic medium, equipotential lines and streamlines intersect strictly at $90^\circ$ right angles, forming elementary “square” grids where $B/L \approx 1$.

2. Finite Difference Scheme for 2D Groundwater Simulation

When physical geometry or boundary conditions become irregular, graphical flow-nets become impractical. In modern computational water resources engineering, numerical approaches such as the Finite Difference Method (FDM) are utilized for Simulation of Groundwater Flow.

By discretizing a 2D horizontal or vertical domain into a grid with spacing $\Delta x$ and $\Delta y$, the partial differential equation for steady-state 2D confined groundwater flow is approximated at grid node $(i, j)$ using central finite differences:

Where $T_x$ and $T_y$ represent the transmissivity in orthogonal directions, and $Q_{i,j}$ accounts for external pumping wells or groundwater recharge sources. Solving this system yields nodal hydraulic head distributions across the entire aquifer domain.

3. Simulation of Seepage Under a Dam (2D Model)

A fundamental practical application in IOE numerical exams is the 2D steady-state seepage simulation beneath a concrete dam with sheet piles. The boundary conditions typically include:

• Upstream Bed Boundary: Known high hydraulic head $h = H_1$ (Dirichlet boundary condition).
• Downstream Bed Boundary: Known low hydraulic head $h = H_2$ (Dirichlet boundary condition).
• Impermeable Dam Base & Impervious Basal Bed: Zero normal flux $\frac{\partial h}{\partial n} = 0$ (Neumann boundary condition).

Using 2D finite difference schemes, engineers solve for the pore water pressure field, total seepage flow rate ($Q = K \cdot H \cdot \frac{N_f}{N_d}$), and uplift pressures acting on the base of hydraulic structures to guarantee geotechnical stability against piping and sliding.

4. 1D Implicit Model for River-Stage and Aquifer Interaction

In riparian corridors, surface water bodies (rivers) dynamically exchange water with adjacent unconfined aquifers. The 1D Boussinesq equation governs this unconfined transient interaction:

Because explicit finite difference schemes suffer from strict numerical stability limits (Courant/von Neumann conditions), an Implicit Finite Difference Scheme (such as the Crank-Nicolson or Fully Implicit method) is employed. This formulation produces a tridiagonal system of linear algebraic equations solved efficiently using the Thomas Algorithm (TDMA), allowing stable numerical time-stepping during flood wave passage or prolonged drought conditions.