Grid Creation and Mapping

To guarantee numeric accuracy, boundary conditions and parameterizations must exactly map to grid boundary faces and grid cells, respectively. DORiE therefore not only requires a specification of the actual grid, but also of these mappings.

DORiE supports two types of grid input, which are controlled via the configuration file parameter grid.gridType:

  1. Building a rectangular grid on the fly (rectangular). Depending on whether grid adaptivity is enabled or not, DORiE will use a structured or an unstructured grid manager. This does not change the input scheme. The user has to specify the config parameters for [grid] dimensions, spatial extensions, and grid cells in each direction. For the mapping of boundary conditions and parameters, mapping datasets are required.

  2. Building a grid according to a GMSH grid file (gmsh). This mode always uses an unstructured grid manager. The user only has to specify the [grid] dimensions and the GMSH gridFile. The mapping is read from this file and must be considered when building the mesh.


The point of origin in DORiE’s reference frame is located at the lower left (front) point of the domain. In 2D, the \(x\)-direction points to the right, and the \(y\)-direction upwards. In 3D, the \(x\)-direction points to the right, the \(y\)-direction points backwards, and the \(z\)-direction points upwards. Notice that DORiE does not support negative coordinates!

Mapping Datasets

Mappings are realized as HDF5 datasets, where every value is mapped to a single grid entity. Changing the grid configuration requires adapting these datasets! All config keys of this part refer to the config section [grid.mapping]. There is one dataset for the mapping of the soil architecture (volume) and one dataset for each segment of the boundary (boundaryXYZ). The datasets contain integer indices, which are then identified with certain boundary conditions or parameterizations given by the respective YAML input files. All mapping datasets must be stored inside the same H5 file (file). Use the Python module h5py for easily creating these data files.

The dataset layout follows the C memory layout, with the primary dimension running fastest. The coordinate system is (right-hand) cartesian, with the x-axis running to the right in both 2D and 3D. All datasets follow the general layout array[z, y, x], where unused dimensions are simply omitted. The first dimension is always the vertical axis. The volume dataset must have the same dimension as the grid itself (grid.dimensions), whereas the boundaryXYZ datasets are interpreted as projections of the boundary in their respective normal direction, reducing the dimension by one. The boundaryLower dataset would have the layout array[x] for a 2D grid, and array[y, x] for a 3D grid. Therefore, boundaryLower and boundaryUpper are “seen from above” (3D-layout: array[y, x]), and boundaryLeft and boundaryRight are “seen from the right” (3D-layout: array[z, y]).

The following Python code creates a volume mapping dataset layered in a file mapping.h5, where the lower half of the domain maps to medium 0, and the upper half to medium 1. The 3D domain contains 10 cells in each direction.

import numpy as np
import h5py

# Create dataset. Mind the data type!
size = 10
layered = np.zeros((size, size, size), dtype=np.int_)
layered[5:, ...] = 1

# Write dataset to file in 'append' mode
with h5py.File("mapping.h5", 'a') as file:
    file.create_dataset("layered", data=layered)

Homogeneous Mappings

If the entire volume, or a complete boundary, should be mapped to a single index, no dataset is required. Instead, you can set the value for the respective config file parameter to the desired index. If no dataset should be read at all, set file to none.

Even with file set to a valid H5 file, DORiE will always try to interpret the input values for grid.mapping as integers. If this succeeds, the value is interpreted as “global” mapping index for the respective part of the grid. Therefore, do not use dataset names starting with digits!

Mapping Soil Layers with GMSH

GMSH supports mapping lines, surfaces, and volumes to physical entities. These entities may combine multiple of the aforementioned geometric entities. Physical entities are assigned a non-negative index upon creation. These indices can be set by the user in any of the GMSH frontends. Indices may not occur multiple times, even if they are assigned to different types of physical entities.

The following .geo GMSH input file creates a rectangular domain that is split in half. The lower part is mapped to index 1, the upper part to index 2. Additionally, a different index for every boundary orientation is set. Notice that the left and right boundaries consist of two geometric lines each. This code can directly be transferred to a Python script using the pygmsh module. It writes a .geo file while checking for a correct syntax within your script. It is readily installed into the virtual environment.

// define geometric entities
Point(1) = {0, 0, 0, 0.1};
Point(2) = {2, 0, 0, 0.1};
Point(3) = {2, 2, 0, 0.1};
Point(4) = {0, 2, 0, 0.1};
Point(5) = {0, 1, 0, 0.1};
Point(6) = {2, 1, 0, 0.1};

Line(1) = {1, 2};
Line(2) = {2, 6};
Line(3) = {6, 3};
Line(4) = {3, 4};
Line(5) = {4, 5};
Line(6) = {5, 1};
Line(7) = {5, 6};

Line Loop(1) = {1, 2, -7, 6};
Plane Surface(1) = {1};
Line Loop(2) = {7, 3, 4, 5};
Plane Surface(2) = {2};

// define physical entities, index in round brackets
Physical Surface(1) = {1}; // lower
Physical Surface(2) = {2}; // upper

// entire set of physical entities must always be defined!
Physical Line(3) = {1}; // bottom
Physical Line(4) = {2, 3}; // right
Physical Line(5) = {4}; // top
Physical Line(6) = {5, 6}; // left

A .geo file is the basis for creating the actual mesh in GMSH. You can load it into the GMSH GUI, or perform the meshing directly using the GMSH command line interface:

gmsh <geo-file> -<dim>

Replace <geo-file> with the appropriate file, and <dim? with the spatial dimension of the intended mesh.