2. Infiltration in Homogeneous Medium

One of the most simple DORiE simulations is the case of constant infiltration on a 2D homogeneous medium. In particular, we will create a similar simulation to the one conducted by Kurt Roth in his Soil Physics Lecture Notes, Chapter 6.3.1.

2.1. Study Case

Consider a uniform and isotropic soil with constant water table (\(h_m = 0\,\text{m}\)) at height \(y=0\,\text{m}\) and vertical flux in the vadose zone. We choose the \(y\)-axis to point upwards. For times \(t < 0\), the water phase is in static equilibrium, i.e., \(j_w = 0\,\text{ms}^{-1}\) in the entire domain. The soil surface is located at \(y=4\,\text{m}\). For \(t \leq 0\), the water flux through the surface boundary is set to \(j_w = 5.56 \cdot 10^{−6}\,\text{ms}^{-1} = 20\,\text{mm}\,\text{h}^{-1}\), equivalent to heavy rainfall.

2.2. Configure DORiE Input Files

In what follows, we will set up input files and run the simulation step-by-step reviewing the most important parts of the DORiE work-flow.

2.2.1. Simulation Mode

The first decision to make is to choose the mode of your simulation. In this case, we are only interested in the water flow movement. Hence, the richards mode in the configuration file is well suited for our purpose.

mode = richards

2.2.2. Grid Creation

For any simulation, the grid settings have to be setup exactly once since all models share the same grid.

In this case, we will create a rectangular grid by specifying the number of cells and the extensions of the domain. Grids of this type are directly created within DORiE, thus, the keyword gridFile is ignored.


The original simulation is one dimensional, but DORiE only supports two and three dimensions. Hence, we use a 2D simulation with 1 cell for the \(x\)-direction, as the simulation is symmetrical along this axis.

We set 1 4 as the extensions of the domain. This means that the rectangular grid will be generated with an extension of \(1\,\text{m}\) in the \(x\)-direction and an extension of \(4\,\text{m}\) in the \(y\)-direction. Notice that the \(x\)-direction points to the right, and the \(y\)-direction upwards. The point of origin in DORiE’s reference frame is located at the lower left point.

Now we have to fill a domain of this size with rectangular grid cells by specifying the number of cells into each direction. For the \(x\)-direction, we will set this to 1. For the \(y\)-direction, we choose a reasonable resolution of \(2\,\text{cm}\) per cell, meaning that we need 200 cells in total. That is, we set the pair 1 200 of cells in the config file.

gridType = rectangular
dimensions = 2
extensions = 1 4
cells = 1 200

2.2.3. Soil Parameterization

First, in the configuration file we set the parameterization file that we want to use. In this case, we select the file richards_param.yml provided by the dorie create command.

file = richards_param.yml

Now, for homogeneous materials, the key grid.mapping.volume in the config file refers to the keyword index to use in the parameterization file for the whole domain. That said, if the parameterization file looks like this:

    index: 0
      type: MvG
        alpha: -2.3
        n: 4.17
        k0: 2.2E-5
        theta_r: 0.03
        theta_s: 0.31
        tau: -1.1
    index: 1
      type: MvG
        alpha: -0.7
        n: 1.3
        k0: 1.0E-5
        theta_r: 0.01
        theta_s: 0.41
        tau: 0.0

  type: None

then, a volume set to 0 will have a parameterization for sand while a volume set to 1 will have a parameterization for silt. For now, let’s say we want to simulate an homogeneous sand.

volume = 0

2.2.4. Initial Condition

The initial condition can be fed with HDF5 data files or with analytic functions. In this case, we set an initial condition with the water table at \(y = 0\,\text{m}\) with a fluid phase in hydrostatic equilibrium. This can be represented by an analytic function where the matricHead is simply set to -h. See the documentation of analytic initial conditions for details.

type = analytic
quantity = matricHead
equation = -h

2.2.5. Boundary Condition

The boundary conditions file has to be specified by the keyword richards.boundary.file in the configuration file. For now, we will use the boundary condition file richards_bc.yml provided by the command dorie create. By default, this file sets a constant infiltration of \(j_w = -5.55\cdot 10^{-6}\,\text{ms}^{-1}\) at the top, a constant matric head of \(h_m = 0\,\text{m}\) at the bottom, and a no-flux condition on the sides of the domain.

file = richards_bc.yml

2.2.6. Output

Finally, we give a name and a path to the output files of the simulation:

fileName = infiltr_homogeneous_sand_2D
outputPath = ./

2.3. Run DORiE

Once everything is set up, we run DORiE by calling the command dorie run followed by the configuration file confing.ini:

dorie run config.ini

By doing this, the simulation should start and provide information about the status of the simulation depending on the on the logLevel keyword on the configuration file. A typical DORiE simulation has the following output:

[18:27:32.762 I]  Starting DORiE
[18:27:32.762 I]  Reading configuration file: config.ini
[18:27:32.768 I]  Creating output directory: ./
[18:27:32.768 I]  Creating a rectangular grid in 2D
[18:27:32.776 I] [RIC]  Loading parameter data file: richards_param.yml
[18:27:32.779 I] [RIC]  Reading boundary condition data file: richards_bc.yml
[18:27:32.784 I] [RIC]  Setup complete
[18:27:33.296 I] [RIC]  Time Step 0: 0.00e+00 + 1.00e+01 -> 1.00e+01
[18:27:33.581 I] [RIC]  Time Step 1: 1.00e+01 + 1.50e+01 -> 2.50e+01
[18:27:33.899 I] [RIC]  Time Step 2: 2.50e+01 + 2.25e+01 -> 4.75e+01
[18:27:34.177 I] [RIC]  Time Step 3: 4.75e+01 + 3.38e+01 -> 8.12e+01
[18:27:46.863 I] [RIC]  Time Step 51: 7.67e+05 + 1.00e+05 -> 8.67e+05
[18:27:46.894 I] [RIC]  Time Step 52: 8.67e+05 + 1.00e+05 -> 9.67e+05
[18:27:46.923 I] [RIC]  Time Step 53: 9.67e+05 + 3.34e+04 -> 1.00e+06
[18:27:46.938 I]  DORiE finished after 1.47e+01s :)

2.4. Results

The results should have been written in several .vtu files, one for each time step, and gathered by a .pvd file. By opening the latter in Paraview (or VisIt) it is possible to visualize the dynamics of the matric head, water content, and water flux as shown below.


Input files