Seepage Analysis — MADEP CPT Interpreter

Scope and positioning

1. Analysis class and intended engineering use

The public seepage route is a two-dimensional steady-state saturated-flow analysis on the Stage 6 section. It is intended for cross-sectional engineering interpretation: phreatic-line refinement, pore-pressure distribution, exit-gradient checking, discharge visualization, and support of slope-stability assessment on the same section geometry.

Model class. The present public route is a Darcy-flow FEM screen on a 2D section. It is not a transient groundwater package, not an unsaturated-flow solver, and not a coupled hydro-mechanical consolidation analysis.

Primary quantities

h(x, y): Total hydraulic head [m].
u(x, y): Pore-water pressure [kPa].
q = (q_x, q_y): Darcy velocity or specific discharge vector [m/s].
∇h: Hydraulic gradient vector [-].
k = diag(k_x, k_y): Permeability tensor aligned with the section axes in the current implementation [m/s].

The default soil source is the interpreted CPT-derived layer model unless custom polygons are enabled.
The same shared measurement line used elsewhere in Stage 6 is used here for arbitrary line probing.
The result may be sampled by Bishop when FEM pore pressure is explicitly enabled.

Governing equations

2. Darcy flow, head field, and derived pore pressure

2.1 Darcy’s law

The seepage solver follows the classical Darcy formulation for saturated flow:

q = −k ∇h

q_x = −k_x ∂h/∂x

q_y = −k_y ∂h/∂y

Notation

q, q_x, q_y: Darcy discharge vector and its Cartesian components [m/s]
k, k_x, k_y: hydraulic conductivity tensor and its principal components [m/s]
h: total hydraulic head referenced to the section datum [m]
x, y: section coordinates [m]
∂h/∂x, ∂h/∂y: head gradient components [-]

Here q is the Darcy velocity or specific discharge, not the actual pore-water velocity. In the current public workflow this distinction is important because the plotted discharge field is a Darcy field; no porosity-based seepage-velocity correction is currently exposed.

2.2 Continuity and the steady-state PDE

Under steady-state conditions, continuity combined with Darcy’s law gives:

∂/∂x (k_x ∂h/∂x) + ∂/∂y (k_y ∂h/∂y) = 0

Notation

k_x, k_y: conductivity components in the section axes [m/s]
h: total hydraulic head [m]
x, y: section coordinates [m]

For homogeneous isotropic soil, this reduces to Laplace’s equation. In the app, the general anisotropic region-wise form is kept because the section can contain multiple permeability zones.

2.3 Head and pore pressure

Once the total head field is known, pore pressure follows directly:

u = γ_w(h − y)

h = y + u / γ_w

Notation

u: pore-water pressure [kPa]
γ_w: unit weight of water [kN/m³]
h: total hydraulic head [m]
y: elevation relative to the section datum [m]

The public seepage route is therefore most naturally interpreted as a head solver with pore pressure and discharge derived afterward. The phreatic line is the locus where u = 0, equivalently h = y.

Input model

3. Shared geometry, materials, and hydraulic inputs

3.1 Shared Stage 6 geometry

The seepage workspace does not introduce a second drawing model. It reuses the shared Stage 6 section geometry: terrain polyline, model base, lateral boundaries, soil polygons, retaining-wall geometry, and the interpreted layer state when custom polygons are not active.

Custom polygons can locally replace the default laterally extended CPT layering.
Retaining walls are represented hydraulically through permeability contrast, not as structural seepage elements. Per-wall material now carries anisotropic conductivity in the wall-local frame (k_⊥ across, k_∥ along); see §4.5.
The seepage mesh is geometry-consistent with the shared Stage 6 section workflow.

3.2 Material data

Each soil zone contributes hydraulic conductivity data:

k = diag(k_x, k_y)

Notation

k: axis-aligned conductivity tensor [m/s]
k_x, k_y: horizontal and vertical conductivity components in the section axes [m/s]

The current public implementation assumes the conductivity tensor is aligned with the global section axes. Rotated permeability tensors are outside the current public scope.

Transparency note. The app preserves manual permeability overrides, uses CPT-derived representative conductivity when present, and otherwise falls back to a documented soil-type heuristic. This is a practical engineering workflow, not a substitute for laboratory or pumping-test calibration.

3.3 Hydraulic boundary inputs

The public route exposes head-type and no-flow-type semantics on the outer boundary of the section. This design keeps the model auditable: the engineer can see exactly where hydraulic conditions are imposed.

Boundary conditions

4. Boundary-condition semantics in the public workflow

4.1 Prescribed head

h = h₀ on Γ_h

Notation

h: solved total head on the boundary [m]
h₀: prescribed boundary head [m]
Γ_h: boundary segment carrying the Dirichlet head condition

This is the standard Dirichlet condition. It is used for reservoir levels, tailwater levels, water-filled excavation boundaries, and other edges where the total head is known.

4.2 No flow

q_n = 0 on Γ_q

Notation

q_n: normal Darcy flux across the boundary [m/s]
Γ_q: boundary segment carrying the Neumann no-flow condition

No-flow is the natural zero-flux boundary. In the current workflow it is the default condition for any outer edge that has not been explicitly assigned another condition.

4.3 Seepage face

A seepage face is a boundary segment where water may exit freely, so the pore pressure is atmospheric and the hydraulic head equals the elevation:

h = y on the active seepage face

Notation

h: total head on the active seepage-face segment [m]
y: boundary elevation [m]

The difficulty is that the active extent is not known a priori. The public route therefore uses an iterative active-set style update in which the candidate seepage face is adjusted until the head field and active boundary extent are mutually compatible.

4.4 Interior drains

A drain is a one-dimensional submanifold Γ_d ⊂ Ω on which the head field is controlled by a user-specified target. The polyline is enforced as a constrained edge in the Triangle PSLG (marker markerType: 'drain') so the mesh conforms to the drain trace; nodes on Γ_d are then candidate Dirichlet sites for the seepage solver.

Notation

Γ_d: drain trace, a polyline embedded in the analysis domain
H_d: prescribed drain head [m]; constant per drain in the current release
ψ_i = h_i − y_i: pressure head at drain node i [m]
R_i: Galerkin reaction at constrained drain node i; sign convention R_i < 0 = water leaving the domain (drain collecting)

Three gating modes are exposed:

always — Γ_d enforces h = H_d at every accepted Newton iterate, independent of pressure state. Appropriate when the drain is permanently below the water table and extraction direction is known a priori.
when-saturated — Dirichlet on Γ_d is gated by a hysteretic pressure-and-flux active set, in direct analogy with the outer-boundary seepage-face logic. A node activates when ψ_i ≥ ψ_activate and the predicted reaction does not indicate injection; it deactivates when ψ_i falls below −ψ_keep or the reaction crosses the injection threshold. Mirrors a relief drain above the water table.
head-cap — one-way head capping. Operationally a linear complementarity problem on Γ_d:

h_i − H_d ≤ 0, R_i ≤ 0, (h_i − H_d)·R_i = 0 ∀ i ∈ Γ_d

These are the Karush–Kuhn–Tucker conditions of an obstacle problem: the head is capped at H_d, the drain may only extract, and complementarity forces each node into either slack (h_i < H_d, R_i = 0) or active-at-cap (h_i = H_d, R_i < 0). The system is solved by a primal–dual active-set semismooth Newton loop run jointly with the outer-boundary seepage-face active set; the joint iteration terminates when the combined mask signature is stable or has entered a two-cycle, with an iteration cap as a final safeguard.

Discharge. Per-drain inflow is reported as a two-element-sided Darcy flux integral across each drain mesh edge, summed over Γ_d. The reaction sum Q_R = −Σ R_i over active nodes is computed in parallel as a self-test and must agree with the per-edge integral to FEM mass-balance accuracy. Reported values are m²/s per metre of out-of-plane length, consistent with the plane-strain convention.

Drain ↔ wall geometry handles weep-hole configurations explicitly: drain segments may share boundary edges with a wall polygon, in which case the drain Dirichlet wins at the shared nodes and the wall material continues to govern element assembly off Γ_d. Drains buried inside or bisecting a wall polygon are rejected at validation.

4.5 Per-wall hydraulic conductivity

Each retaining wall carries its own material with anisotropic conductivity in the wall-local (across, along) frame:

K_wall^local = diag(k_⊥, k_∥) ≡ diag(k_across, k_along)

Notation

k_⊥: across-wall conductivity [m/s]; controls leakage normal to the wall plane
k_∥: along-wall conductivity [m/s]; controls leakage along the wall axis (sheet-pile joints, slurry-wall vertical pathways)

For the axis-aligned vertical walls supported in the current public route, the global-frame tensor is K_wall^global = diag(k_⊥, k_∥) and the element conductivity matrix is unchanged in form. General-orientation walls require the full 2×2 tensor in the element matrix and are deferred.

Five engineering presets are exposed:

Sheet pile — k_⊥ ≈ 10⁻¹⁰ m/s, k_∥ ≈ 10⁻⁸ m/s. Joint leakage along the sheet interlocks is a small but non-negligible vertical pathway.
Slurry / diaphragm wall — isotropic 10⁻⁹ m/s.
Soil-mix wall — isotropic 10⁻⁷ m/s.
Relief wall — 10⁻⁶ to 10⁻⁵ m/s. Deliberate leakage by design.
Legacy impermeable — k_⊥ = k_∥ = 10⁻¹⁰ m/s. Back-fill default for saved projects predating per-wall material; result identical to the legacy hard-coded cut-off within FEM tolerance.

Numerical floor. The element-matrix conductivity floor is 10⁻²⁰ m/s, set to keep the assembled system non-singular. In practice a wall conductivity below 10⁻¹⁰ × k_soil,min is asymptotically indistinguishable from a hard cut-off and degrades CG convergence without changing the head field; the UI flags this as the practical floor.

Finite-element formulation

5. T3 formulation, assembly, and solution strategy

5.1 Element family

The public seepage route uses three-node linear triangles. The head field is interpolated linearly over each element, and the hydraulic gradient is therefore constant within an individual triangle.

Interpretation consequence. The plotted head field is nodal and continuous, but |∇h|, |q|, and the directional discharge components are piecewise constant at the element level before contour post-processing.

5.2 Element conductivity matrix

For each T3 element, the conductivity matrix is assembled in the standard form:

K_e = ∫_Ωe B^T D B dΩ

Notation

K_e: element conductivity matrix [m/s] for unit out-of-plane thickness
Ω_e: element area [m²]
B: head-gradient matrix of the T3 shape functions [1/m]
D: element conductivity tensor [m/s]

where D is the element conductivity tensor and B is the gradient matrix associated with the linear head shape functions.

5.3 Global assembly and BC application

The global head system is assembled from the element matrices and then modified for the imposed Dirichlet conditions. The resulting linear system is symmetric positive definite under the normal well-posed seepage configurations used in the public app.

5.4 Mesh generation and region tagging

The mesh is built from the constrained section geometry. Material regions are preserved as meshing constraints, so hydraulic interfaces remain explicit in the FE model. Region tagging determines which conductivity tensor is assigned to each element.

Free surface

6. Free-surface iteration and seepage-face activation

Unconfined flow introduces nonlinearity because the extent of the saturated domain is part of the solution. The current public route avoids remeshing on every iteration by using a practical active/inactive treatment of wet and dry regions together with the seepage-face update.

The entire section is meshed once; the wet or dry interpretation is updated iteratively.
The seepage face is activated only where the solution indicates outward drainage under atmospheric pressure.
The public route is therefore a practical engineering free-surface approximation, not a fully general transient saturation algorithm.

Outputs and interpretation

7. Head, pore pressure, discharge, and arbitrary line probes

7.1 Canvas fields

Field	Meaning	Interpretation note
h	Total head	Primary solved nodal field.
u	Pore pressure	Derived from `u = γw(h - y)`.
\|∇h\|	Hydraulic-gradient magnitude	Element-wise constant before contour post-processing.
\|q\|, q_x, q_y	Darcy discharge field	Specific discharge, not pore-water velocity.

7.2 Shared line probe

The shared Stage 6 measurement line doubles as a seepage probe. It can sample head, pore pressure, gradient magnitude, discharge magnitude, directional discharge components, and the line-normal discharge:

q_n = q · n

Notation

q_n: probe-normal Darcy discharge component [m/s]
q: Darcy discharge vector at the sampled point [m/s]
n: unit normal to the user-defined probe line [-]

Because q_n depends on probe orientation, it is intentionally a line quantity rather than a global contour variable.

Coupling

8. Coupling of seepage results to slope stability

The seepage result can replace the simple hydrostatic pore-pressure interpretation in the Bishop workflow. When enabled, slice-base pore pressure is sampled from the FEM field rather than from the drawn phreatic line alone.

Important. This coupling is opt-in. If no seepage result exists, or if a requested point lies outside the seepage solution domain, the stability workflow falls back to its hydrostatic interpretation rather than silently failing.

Assumptions and limitations

9. Transparent modeling boundaries

The public seepage route is steady-state only.
Unsaturated flow, suction, and Richards-equation behavior are not exposed.
No hydro-mechanical coupling is exposed in the public workflow.
Permeability anisotropy in soil regions is axis-aligned in the current implementation. Wall material now exposes an anisotropic (k_⊥, k_∥) tensor in the wall-local frame; general-orientation walls requiring a full 2×2 global tensor are deferred.
Head-prescribed interior drains with three gating modes (always, when-saturated, head-cap) are supported through a Dirichlet active-set on the conforming PSLG. Flow-prescribed line sinks (q = −Q δ(x − x_w)), finite-resistance Cauchy drains (q = β(h − H_d)), and drain-water-level coupling to an outlet elevation are not yet exposed.
The route is intended for engineering screening and coupled interpretation, not as a replacement for a dedicated regional hydrogeological model.

References

10. Reference basis

Darcy, H. (1856) — foundational law for saturated flow through porous media.
EN 1997-1 — Eurocode 7 hydraulic-heave and piping context for geotechnical interpretation.
Bear, J. Hydraulics of Groundwater. McGraw-Hill, 1979. ISBN 0-07-004170-9. Reference for anisotropic Darcy flow and the wall-local conductivity tensor framing in §4.5.
Bathe, K.-J., and Khoshgoftaar, M. R. "Finite Element Free Surface Seepage Analysis Without Mesh Iteration." International Journal for Numerical and Analytical Methods in Geomechanics, 3(1), 13–22, 1979. DOI: 10.1002/nag.1610030103. Reference for the fixed-mesh wet/dry residual-flow approach used in iterate mode and shared by the drain active-set update.
Hintermüller, M., Ito, K., and Kunisch, K. "The Primal-Dual Active Set Strategy as a Semismooth Newton Method." SIAM Journal on Optimization, 13(3), 865–888, 2002. DOI: 10.1137/S1052623401383558. Mathematical basis for the head-cap LCP solver and its joint coupling with the outer seepage-face active set in §4.4.
Standard finite-element seepage formulations — the present route follows the classical T3 Darcy/Laplace formulation documented in the project seepage specification.
GeoStudio / SEEP/W style free-surface practice — relevant as background for the practical active-set and reduced-dry-zone style interpretation used in engineering seepage workflows.

For the exhaustive derivation, pseudocode, and implementation-specific notes, continue to the full seepage reference section.