Richards’ Equation, transient problems

Nearly all production ATS simulations start with subsurface flow as their basis. The basic form of this is the Richards equation.

This tests and demonstrates a variety of transient Richards equation solutions, demonstrating boundary conditions and seepage face conditions. All are on a single column of cells, oriented vertically, and are initialized as the hydrostatic solution. Infilitration is turned on, and/or seepage occurs to raise or lower the water table. These boundary conditions are useful in 2D transect simulations and watershed simulations as well; this demo shows how to use them.

Author: Ethan Coon

%matplotlib inline

import sys,os
sys.path.append(os.path.join(os.environ['ATS_SRC_DIR'],'tools', 'utils'))
sys.path.append(os.path.join(os.environ['ATS_SRC_DIR'],'docs', 'documentation', 'source', 'ats_demos'))

import ats_xdmf
import plot_column_data
import numpy as np
from matplotlib import pyplot as plt

import run_demos

def plot(testname, title, force=False):
    """A generic runner and plotting function that plots saturation as a function of depth."""
    run_demos.run_demo_local(testname, force=force)
    dirname = testname+'.demo'

    fig = plt.figure()

    vis = ats_xdmf.VisFile(dirname, time_unit='d')
    vis.loadMesh(columnar=True)

    # plot the solution
    ax = fig.subplots(1,1)
    plot_column_data.plot_subsurface(vis, ['saturation_liquid',], ax, None, cmap='jet')

    ax.set_xlabel('saturation [-]')
    ax.set_ylabel('z [m]')

The first example simply looks at infiltrating into a column with no flux boundary conditions, slowly filling up the subsurface. Note that the color scheme goes from the initial condition in blue to later times in red.

run_demos.run_demo_local('infiltration')
plot("./infiltration", "infiltration example")

  Test log file : LOGS/ats-tests-2024-04-22_12-02-41.testlog
Running ATS demo problems :
02_richards          | .

----------------------------------------------------------------------
Regression test summary:
    Total run time: 2.01131 [s]
    Total tests : 1
    Tests run : 1
    All tests passed.

The next set of problems look at seepage faces, which are boundary conditions that allow water to seep out if the soil is saturated, but don’t allow water in. Starting from a water table below the seepage face at 8m, we infiltrate until the water table rises to this seepage face. The water table does not rise above the seepage face, but instead flows out of the domain reaching steady-state.

plot("seepage_infiltration", "infiltration to seepage face")

  Test log file : LOGS/ats-tests-2024-04-22_12-02-43.testlog
Running ATS demo problems :
02_richards          | .

----------------------------------------------------------------------
Regression test summary:
    Total run time: 3.01102 [s]
    Total tests : 1
    Tests run : 1
    All tests passed.

In the second example, we start from a water table above the seepage face, and let the water seep out until the water table draws down to the seepage face. The water table does not go below the seepage face at 6m.

plot("seepage_drawdown", "drawdown example")

  Test log file : LOGS/ats-tests-2024-04-22_12-02-46.testlog
Running ATS demo problems :
02_richards          | .

----------------------------------------------------------------------
Regression test summary:
    Total run time: 1.0108 [s]
    Total tests : 1
    Tests run : 1
    All tests passed.

The next example shows how a seepage face works in a more dynamic system that ebbs and flows. We set a seepage face above the water table, and infiltrate until the water table is at the seepage face. Then we draw down the water table by taking water out of the bottom of the domain. The water table drops below the seepage face.

This shows that the seepage face can turn both on and off, i.e. swapping back from Neumann to Dirichlet to Neumann again.

plot("seepage_exfiltration", "on-and-off example, FV")

  Test log file : LOGS/ats-tests-2024-04-22_12-02-47.testlog
Running ATS demo problems :
02_richards          | .

----------------------------------------------------------------------
Regression test summary:
    Total run time: 4.01255 [s]
    Total tests : 1
    Tests run : 1
    All tests passed.

Demonstration and test of a more general seepage condition in which an infiltration flux is specified on the seepage face. The condition has constraints on both sides – infiltration is less than or equal to the specified flux, and the pressure is less than atmospheric. This is useful for, for instance, spinning up coupled surface-subsurface runs. The subsurface only is run, with this type of BC everywhere and the infiltration as the mean annual rainfall rate. Everywhere there are rivers end up with atmospheric pressure, while everywhere there is not ponded water sees the infiltration condition (at steady-state).

In this problem, we set the seepage face on the top surface. Rain is continuous on this surface for all time. The water table is initially below the surface. The seepage face infiltrates until the entire column is saturated, at which point the seepage face sets the max pressure as p_atm. Then, after 1/2 year, the bottom boundary condition is set to drain the column at twice the rain rate. This draws down the water table below the surface, and the rain infiltration is again active. At 3/4 years, the drainage condition is set to equal the rainfall rate, ensuring that steady state is reached and the water table stays below the surface.

plot("infiltration_then_seepage", "infiltrate-then-seep example")

  Test log file : LOGS/ats-tests-2024-04-22_12-02-51.testlog
Running ATS demo problems :
02_richards          | .

----------------------------------------------------------------------
Regression test summary:
    Total run time: 2.01109 [s]
    Total tests : 1
    Tests run : 1
    All tests passed.