Commit 298a8fd2 authored by Lukas Riedel's avatar Lukas Riedel

Merge branch 'master' into 133-implement-a-class-for-initial-conditions-from-hdf5-data

Merge and update user docs.
parents d3f42c97 a1afde48
......@@ -5,8 +5,14 @@ variables:
CPUS_DIND: 2
DOCKER_LOGIN:
docker login -u $DOCKER_HUB_USER -p $DOCKER_HUB_PW
DUNE_ENV_IMAGE:
dorie/dune-env:2.6
# DUNE environment image
BASE_IMAGE: dorie/dune-env
# Use semantic versioning (not the version of DUNE) and bump according to
# to whether changes are backwards-compatible or not.
IMAGE_VERSION: "1.0"
DUNE_ENV_IMAGE: ${BASE_IMAGE}:img-v${IMAGE_VERSION}
CMAKE_FLAGS:
-DDUNE_PYTHON_VIRTUALENV_SETUP=True
-DDUNE_PYTHON_ALLOW_GET_PIP=True
......
......@@ -6,6 +6,12 @@
### Added
* DORiE now depends on [`yaml-cpp`](https://github.com/jbeder/yaml-cpp/), a
library for reading YAML files. The version required is >=5.2.0.
* New classes representing parameterizations. Every parameterization must now
derive from `RichardsParameterization`, and return function objects
representing parameterization functions.
* Added a abstract base class [`SimulationBase`](dune/dorie/interface/base_simulation.hh)
that models the concept of simulation steps so that they
can be later coupled with other simulations.
* Added an abstract base class
[`SimulationBase`](dune/dorie/common/simulation_base.hh)
that models the concept of simulation steps so that they
......@@ -28,6 +34,10 @@
* Add logging framework 'spdlog' as submodule !106
* Support input of boundary segmentations !119
* Reconstruction of continuous fluxes using RT finite elements !105
* Custom DG finite element map for simplices based on Pk polynomials !125
* Initial conditions expressed as analytic functions using
[muparser](http://beltoforion.de/article.php?a=muparser) !131
* Coupling between transient models for water flow and solute transport !96
### Changed
* `Simulation` is renamed `RichardsSimulation` and moved to
......@@ -35,9 +45,6 @@
* `RichardsSimulation` now has its own `RichardsSimulationTraits` derived from
`BaseTraits`, which defines all its member types. `BaseTraits` now have
reduced content and are intended to be shared between models/simulations.
* Grid adaptation now is done in two steps: (1) mark entities of the grid to
refine/coarse and (2) adapt the grid and project the degrees of freedom on
the new grid.
* Data structures for storing and accessing parameters. The new class
`FlowParameters` maps parameter sets and scaling factors to every grid cell
on the coarsest grid level, and only stores one set of parameters for
......@@ -45,6 +52,15 @@
to call `FlowParameters::bind` with a grid entity to cache the appropriate
data. All classes and functions querying parameter data have been updated.
The old parameterization classes located in
[`param_base.hh`](dune/dorie/solver/param_base.hh),
[`param_factory.hh`](dune/dorie/solver/param_factory.hh), and
[`param_van_genuchten.hh`](dune/dorie/solver/param_van_genuchten.hh)
are still used for reading in data for the new storage scheme.
The respective objects are freed once a simulation commences.
* Grid adaptation now is done in two steps: (1) mark entities of the grid to
refine/coarse and (2) adapt the grid and project the degrees of freedom on
the new grid.
The old parameterization classes have been removed.
* `RichardsSimulation` now uses one vector of coefficients instead of two, which
now is a `shared_ptr` instead a `unique_ptr`.
......
......@@ -77,6 +77,7 @@ by CI tests.
| SuperLU | 5.2 |
| [yaml-cpp](https://github.com/jbeder/yaml-cpp) | >= 5.2.0 |
| [spdlog](https://github.com/gabime/spdlog) | 1.1.0 | Included as Git Submodule
| [muparser](http://beltoforion.de/article.php?a=muparser) | master |
| [dune-common](https://gitlab.dune-project.org/core/dune-common) | releases/2.6
| [dune-geometry](https://gitlab.dune-project.org/core/dune-geometry) | releases/2.6
| [dune-grid](https://gitlab.dune-project.org/core/dune-grid) | releases/2.6
......@@ -116,14 +117,14 @@ If you installed [Anaconda](https://conda.io/docs/user-guide/install/download.ht
apt update
apt install cmake doxygen gcc g++ gfortran \
git libatlas-base-dev libfftw3-dev libfftw3-mpi-dev \
libfreetype6-dev libhdf5-mpi-dev libopenmpi-dev libpng-dev \
libsuperlu-dev libyaml-cpp-dev libxft-dev \
python3-dev python3-pip
libfreetype6-dev libhdf5-mpi-dev libmuparser-dev \
libopenmpi-dev libpng-dev libsuperlu-dev libyaml-cpp-dev \
libxft-dev python3-dev python3-pip
**macOS:**
brew update
brew install cmake doxygen gcc libpng open-mpi \
brew install cmake doxygen gcc libpng open-mpi muparser \
pkg-config python3 superlu yaml-cpp
brew install hdf5 --with-mpi
brew install fftw --with-mpi
......
......@@ -9,6 +9,7 @@ FIND_PACKAGE (FFTW REQUIRED)
FIND_PACKAGE (SuperLU REQUIRED)
FIND_PACKAGE (MPI REQUIRED)
find_package (yaml-cpp 0.5.2 REQUIRED)
FIND_PACKAGE (muparser REQUIRED)
FIND_PACKAGE (METIS)
FIND_PACKAGE (ParMETIS)
......@@ -20,7 +21,8 @@ include_directories(${FFTW_INCLUDES}
list (APPEND DUNE_LIBS
${FFTW_LIBRARIES}
${HDF5_LIBRARIES}
${YAML_CPP_LIBRARIES})
${YAML_CPP_LIBRARIES}
muparser::muparser)
MESSAGE(STATUS "DUNE Libraries: ${DUNE_LIBS}")
# Add DORiE testing functions
......
# - Find the muParser library
#
# Usage:
# find_package(muParser [REQUIRED] [QUIET] )
#
# Hints may be given by the (environment) variables
# MUPARSER_ROOT ... Path to the installation library
#
# It sets the following variables:
# MUPARSER_FOUND ... true if muParser is found on the system
# MUPARSER_LIBRARIES ... full path to muParser library
# MUPARSER_INCLUDES ... muParser include directory
#
# It defines the following targets:
# muparser::muparser ... muparser library to link against
#
find_library(MUPARSER_LIBRARY
NAMES muparser muparserd
HINTS ${MUPARSER_ROOT} ENV MUPARSER_ROOT
)
find_path(MUPARSER_INCLUDE_DIR
muParserDef.h
HINTS ${MUPARSER_ROOT} ENV MUPARSER_ROOT
)
mark_as_advanced(MUPARSER_INCLUDE_DIR MUPARSER_LIBRARY)
find_package_handle_standard_args(MUPARSER
DEFAULT_MSG
MUPARSER_LIBRARY MUPARSER_INCLUDE_DIR)
if (MUPARSER_FOUND)
set(MUPARSER_LIBRARIES ${MUPARSER_LIBRARY} )
set(MUPARSER_INCLUDES ${MUPARSER_INCLUDE_DIR} )
# add the target
add_library(muparser::muparser SHARED IMPORTED)
set_target_properties(muparser::muparser
PROPERTIES IMPORTED_LOCATION ${MUPARSER_LIBRARIES}
INTERACE_INCLUDE_DIRECTORIES ${MUPARSER_INCLUDES}
)
endif()
......@@ -4,4 +4,4 @@ spatial_resolution_west -1
spatial_resolution_east -1
number_BC_change_times 2
0 neumann 0 dirichlet 1 neumann 0
1E7 neumann 0 neumann 0 neumann 0
\ No newline at end of file
1E3 neumann 0 neumann 0 neumann 0
......@@ -12,4 +12,4 @@ spatial_resolution_back_sn -1
spatial_resolution_back_we -1
number_BC_change_times 2
0 neumann 0 neumann 0 neumann 0 neumann 0 dirichlet 1 neumann 0 neumann 0 neumann 0 neumann 0
1E7 neumann 0 neumann 0 neumann 0 neumann 0 neumann 0 neumann 0 neumann 0 neumann 0 neumann 0
\ No newline at end of file
1E3 neumann 0 neumann 0 neumann 0 neumann 0 neumann 0 neumann 0 neumann 0 neumann 0 neumann 0
\ No newline at end of file
......@@ -117,7 +117,10 @@ function(create_default_config INPUT OUTPUT SOURCE_DIR CSS)
endfunction()
# copy BC & parameter files
file(COPY 2d_infiltr.bcdat 3d_infiltr.bcdat param.yml DESTINATION .)
file(COPY 2d_infiltr.bcdat 3d_infiltr.bcdat
2d_solute.bcdat 3d_solute.bcdat
param.yml
DESTINATION .)
# Random field generator
scrape_parameters(
......@@ -150,14 +153,6 @@ scrape_parameters(
MODEL "richards"
)
# Final Dorie config
create_default_config(
${CMAKE_CURRENT_BINARY_DIR}/config.yml
"config.ini;parameters.html;parameters.rst"
${PROJECT_SOURCE_DIR}/dune/dorie
${CMAKE_CURRENT_SOURCE_DIR}/parameters.css
)
# Model: Transport
scrape_parameters(
SOURCE_DIR ${PROJECT_SOURCE_DIR}/dune/dorie/model/transport
......@@ -166,11 +161,10 @@ scrape_parameters(
MODEL "transport"
)
# Transport config (temporally hidden to the user)
# TODO: Include in API once transport is ready
# Final Dorie config
create_default_config(
${CMAKE_CURRENT_BINARY_DIR}/config.yml
"transport-config.ini;transport-parameters.html;transport-parameters.rst"
"config.ini;parameters.html;parameters.rst"
${PROJECT_SOURCE_DIR}/dune/dorie
${CMAKE_CURRENT_SOURCE_DIR}/parameters.css
)
)
\ No newline at end of file
......@@ -29,8 +29,16 @@ All attributes are optional.
The parser supports rudimentary markdown / styling. You can add a paragraph by
adding an empty line, make text **bold** or ``monospaced``.
-->
<dorie>
<category name="simulation">
<parameter name="mode">
<definition> Mode of simulation in DORiE. </definition>
<suggestion> richards </suggestion>
<values> richards, richards+transport </values>
<comment> richards, richards+transport </comment>
</parameter>
</category>
<category name="grid">
<parameter name="gridType">
<definition> Grid geometry. The grid can either be structured rectangular / cubic
......
......@@ -61,6 +61,12 @@ adding an empty line, make text **bold** or ``monospaced``.
<suggestion> true </suggestion>
</parameter>
<parameter name="policy">
<definition> Policy to write the data. </definition>
<suggestion> endOfRichardsStep </suggestion>
<values> endOfRichardsStep, none </values>
</parameter>
<parameter name="subsamplingLevel">
<definition> Plot VTK files with virtually refined grids. VTK only
supports bilinear triangulations and displays higher-order solutions
......@@ -107,11 +113,10 @@ adding an empty line, make text **bold** or ``monospaced``.
<category name="initial">
<parameter name="type">
<definition> The type of data to create the initial condition from.
</definition>
<values> data, linear </values>
<suggestion> linear </suggestion>
<comment> Choose source type: data, linear </comment>
<definition> There are currently two possible initial conditions to choose from: A HDF datafile (``data``), or analytic equations (``analytic``).</definition>
<values> data, analytic </values>
<suggestion> analytic </suggestion>
<comment> Choose initial condition type: data, analytic </comment>
</parameter>
<parameter name="quantity">
......@@ -123,20 +128,10 @@ adding an empty line, make text **bold** or ``monospaced``.
<comment> Choose quantity represented: matricHead </comment>
</parameter>
<parameter name="headLower">
<definition> Initial matric head at the lower boundary of the domain.
(``linear`` type only)
</definition>
<values> float </values>
<suggestion> 0 </suggestion>
</parameter>
<parameter name="headGradient">
<definition> Gradient of the matric head towards the upper boundary
(positive y coordinate). (``linear`` type only)
</definition>
<values> float </values>
<suggestion> -1 </suggestion>
<parameter name="equation">
<definition> Equation for the initial condition </definition>
<values> equation [x,y,z,h] </values>
<suggestion> -h </suggestion>
</parameter>
<parameter name="file">
......
......@@ -61,6 +61,12 @@ adding an empty line, make text **bold** or ``monospaced``.
<suggestion> true </suggestion>
</parameter>
<parameter name="policy">
<definition> Policy to write the data. </definition>
<suggestion> endOfRichardsStep </suggestion>
<values> endOfTransportStep, endOfRichardsStep, none </values>
</parameter>
<parameter name="subsamplingLevel">
<definition> Plot VTK files with virtually refined grids. VTK only
supports bilinear triangulations and displays higher-order solutions
......@@ -106,11 +112,40 @@ adding an empty line, make text **bold** or ``monospaced``.
</category>
<category name="initial">
<parameter name="value">
<definition> Initial constant value over the whole domain. </definition>
<values> float </values>
<parameter name="type">
<definition> There are currently two possible initial conditions to choose from: A HDF datafile (``data``), or analytic equations (``analytic``).</definition>
<values> data, analytic </values>
<suggestion> analytic </suggestion>
<comment> Choose initial condition type: data, analytic </comment>
</parameter>
<parameter name="quantity">
<definition> The physical quantity represented by the initial condition
data.
</definition>
<values> soluteConcentration </values>
<suggestion> soluteConcentration </suggestion>
<comment> Choose quantity represented: soluteConcentration </comment>
</parameter>
<parameter name="equation">
<definition> Equation for the initial condition </definition>
<values> equation [x,y,z,h] </values>
<suggestion> 0 </suggestion>
</parameter>
<parameter name="datafile">
<definition> Path to the initial condition data file.
(``data`` type only)
</definition>
<values> path </values>
</parameter>
<parameter name="dataset">
<definition> Dataset to use as initial condition. (``data`` type only)
</definition>
<values> string </values>
</parameter>
</category>
<category name="time">
......@@ -178,8 +213,8 @@ adding an empty line, make text **bold** or ``monospaced``.
</category>
<category name="parameters">
<parameter name="molecularDiffusion">
<definition> Molecular diffusion of the fluid phase </definition>
<parameter name="effHydrDips">
<definition> Effective hydrodynamic dispersion </definition>
<values> float &gt; 0 </values>
<suggestion> 2E-9 </suggestion>
</parameter>
......@@ -187,20 +222,32 @@ adding an empty line, make text **bold** or ``monospaced``.
<category name="numerics">
<parameter name="timestepMethod">
<definition> Numerical scheme to perform time steps in the simulation
<definition> Numerical scheme to perform time steps in the simulation.
``alex2`` and ``implicit_euler`` are implicit methods.
``explicit_euler`` is a explicit method.
</definition>
<values> string </values>
<suggestion> explicit_euler, implicit_euler, heun, shu3, rk4, alex2, fractional, alex3
</suggestion>
<values> explicit_euler, implicit_euler, alex2 </values>
<suggestion> alex2 </suggestion>
</parameter>
<parameter name="courant">
<definition> Courant number for explicit methods. Explicit methods have a
maximum time step according to the CFL condition.
<definition> Courant number for explicit methods. It is a scale for the
maximum stable time step according to the CFL condition. The CFL
condition penalizes high velocities, low dispersion coefficients, and
small grid elements. Hence, Courant numbers near to 1 tend to be more
instable while numbers near to 0 tend to considerably limit the
maximum time step.
</definition>
<values> float &gt; 0 </values>
<values> 0 &lt; float &lt; 1 </values>
<suggestion> 0.8 </suggestion>
</parameter>
<parameter name="FEorder">
<definition> Order of the finite element method used. Values &gt; 1 are not
thoroughly tested. </definition>
<suggestion> 0 </suggestion>
<values> 0 </values>
</parameter>
</category>
<category name="fluxReconstruction">
......@@ -230,4 +277,22 @@ adding an empty line, make text **bold** or ``monospaced``.
<suggestion> 1E-10 </suggestion>
</parameter>
</category>
<category name="solverParameters">
<parameter name="reduction">
<definition> Required relative defect reduction.
Reduce this value to increase precision.
</definition>
<suggestion> 1E-6 </suggestion>
<values> float </values>
</parameter>
<parameter name="minDefect">
<definition> Minimum absolute defect at which linear solver stops.
Reduce this value to increase precision.
</definition>
<suggestion> 1E-20 </suggestion>
<values> float </values>
</parameter>
</category>
</dorie>
......@@ -47,13 +47,29 @@ Output Files
saved into the directory from which DORiE was executed.
* VTK Files (``.vtu``): The output data in unstructured grid format.
It contains multiple datasets:
- ``head``: Matric head
- ``K0``: Saturated conductivity
- ``flux``: Water flux (not reconstructed!)
- ``theta_w``: Volumetric water content
- ``Theta``: Water saturation
It contains multiple multiple datasets depending on the type of model:
.. object:: Richards
- ``head``: Matric head :math:`h_m \, [\mathrm{m}]`
- ``K0``: Saturated conductivity :math:`K_0 \, [\mathrm{ms}^{-1}]`
- ``flux``: Water flux :math:`\mathbf{j}_w \, [\mathrm{ms}^{-1}]`
(not reconstructed!)
- ``flux_RT<k-1>``: Reconstructed water flux
:math:`\mathbf{j}_{w, \mathrm{rc}} \, [\mathrm{ms}^{-1}]`
(if enabled!)
- ``theta_w``: Volumetric water content :math:`\theta \, [\mathrm{-}]`
- ``Theta``: Water saturation :math:`\Theta \, [\mathrm{-}]`
.. object:: Transport
- ``solute``: Solute concentration in water phase
:math:`c_w \, [\mathrm{kg}/\mathrm{m}^3]`
- ``solute_total``: Total solute concentration
:math:`c_t \, [\mathrm{kg}/\mathrm{m}^3]`
- ``flux_RT<k-1>``: Reconstructed solute flux
:math:`\mathbf{j}_{s, \mathrm{rc}} \, [\mathrm{kg}/\mathrm{m}^{-2}\mathrm{s}^{-1}]`
(if enabled!)
* VTK Parallel Collection Files (``.pvtu``): Merging multiple VTK files in case
DORiE is run in parallel.
......
......@@ -36,5 +36,40 @@ soil. Depending on the soil architecture, the hydraulic parameters may be
discontinuous which would lead to water content and conductivity exhibiting
singularities.
Use the value ``richards`` in the ``simulation.mode`` keyword to run the
standalone Richards solver.
Solute Transport Solver
=======================
The solute transport equation for unsaturated media describes movement of
solute by the *solute concentration* :math:`c_w \, [\mathrm{kg}/\mathrm{m}^3]`,
.. math::
\frac{\partial}{\partial t} c_w \theta_w
- \nabla \left[ \vec{j}_w c_w + D \theta_w \nabla c_w
\right]
= 0,
where :math:`D \, [\mathrm{m}^2/\mathrm{s}]` is the effective hydrodynamic
dispersion tensor, and :math:`\vec{j}_w` is the water flux. Currently,
:math:`D` can only be used as a constant parameter by the keyword
``effHydrDips`` in the configuration file. The solute concentration per total
volume :math:`c_t \, [\mathrm{kg}/\mathrm{m}^3]` is calculated by multiplying
the concentration in the water phase with the volumetric water content,
.. math::
c_t = c_w \theta_w
.. TODO: Update effective hydrodynamic dispersion description once
parameterization is merged in master.
The solute transport requires a water flux of the unsaturated media, hence,
the Richards solver must be utilized to solve them. Use the value
``richards+transport`` in the ``simulation.mode`` keyword to run the solute
transport solver.
.. _DUNE: https://www.dune-project.org/
.. _DUNE-PDELab: https://www.dune-project.org/modules/dune-pdelab/
......@@ -21,15 +21,36 @@ Dirichlet
.. object:: dirichlet <head>
Set a fixed matric head [m] at the boundary.
Set a fixed matric head [:math:`m`] at the boundary.
.. object:: dirichlet <solute>
Set a fixed solute concentration at the boundary.
Neumann
+++++++
.. object:: neumann <flux>
.. object:: neumann <water_flux>
Set a fixed volumetric water flux [:math:`m/s`] at the boundary. Positive values
imply fluxes out of the domain, negative values imply fluxes into the domain.
.. object:: neumann <solute_flux>
Set a fixed volumetric solute flux [:math:`kg/(s\,m^3)`] at the boundary. Positive values
imply fluxes out of the domain, negative values imply fluxes into the domain.
Outflow
+++++++
Set a fixed volumetric flux [m/s] at the boundary. Positive values imply fluxes out of the domain, negative values imply fluxes into the domain.
.. object:: outflow <solute_outflux>
Set a spatial rate of change of volumetric solute flux
[:math:`kg/(s\,m\,m^3)`] at the boundary. Positive values imply higher rate
of change out of the domain, negative values imply lower rate of change into
the domain. Influx situations are not taken into account. For most of the
cases an outflow of 0. [:math:`kg/(s\,m\,m^3)`] is a reasonable situation.
Datafile Structure
==================
......@@ -83,7 +104,7 @@ These lines follow a simple grammar:
bc_line: time { group }*
time: `float`
group: ( bc_type `float` )
bc_type: "neumann" | "dirichlet"
bc_type: "neumann" | "dirichlet" | "outflow"
The boundary conditions defined here are parsed in the same order as the boundary segments have been specified. In 3D, the rectangular boundary segments are parsed in a tabular fashion, where columns run faster than rows. Columns are defined along the first direction specified in the `Boundary Segmentation`_, and rows are defined along the second direction.
......
......@@ -6,9 +6,9 @@ Depending on the application for what you want to use DORiE, there might be the
Understanding the water flux output
-----------------------------------
Firstly, we have to recall that DORiE solves a Discontinuous Galerking finite element problem with *matric head* as unknown. It means that the solution of the *matric head* (and therefore the *water flux*) is continuous only element-wise, or in other words, it is discontinuous on the intersections between elements. On the other hand, the dG method solves *numerical fluxes* on the intersections between elements composed together with a penalty term that increases with respect to the discontinuity of the *matric head*. This ensures that the local solution is conservative while keeps the discontinuity as low as possible.
Firstly, we have to recall that DORiE solves a Discontinuous Galerking finite element problem with *matric head*/*solute concentration* as unknown. It means that the solution of the *matric head*/*solute concentration* (and therefore the *water flux*/*solute flux*) is continuous only element-wise, or in other words, it is discontinuous on the intersections between elements. On the other hand, the dG method solves *numerical fluxes* on the intersections between elements composed together with a penalty term that increases with respect to the discontinuity of the *matric head*/*solute concentration*. This ensures that the local solution is conservative while keeps the discontinuity as low as possible.
From the description above one can infer that one has to distinguish between *water fluxes* at the interior of each element and at the intersections of all elements (we call these intersections skeleton of the grid). Unfortunately, there is no a standard form to write the skeleton fluxes on formats like VTK and that's the main reason why DORiE only provides the interior fluxes. However, assuming one can write both fluxes into some output format, they are still discontinuous (notice that direct use of discontinuous fluxes are useless for conservative computations since the transported quantities are very likely to get stagnated or over-transported in the nearby of intersections between elements). It means that it is needed some sort of post-processing that ensures that the *water mass* is still locally and globally conserved.
From the description above one can infer that one has to distinguish between * fluxes* at the interior of each element and at the intersections of all elements (we call these intersections skeleton of the grid). Unfortunately, there is no a standard form to write the skeleton fluxes on formats like VTK and that's the main reason why DORiE only provides the interior fluxes. However, assuming one can write both fluxes into some output format, they are still discontinuous (notice that direct use of discontinuous fluxes are useless for conservative computations since the transported quantities are very likely to get stagnated or over-transported in the nearby of intersections between elements). It means that it is needed some sort of post-processing that ensures that the *mass* is still locally and globally conserved.
.. DORiE allow finite volume computations under certain specific conditions. In such case, if generated, the raw flux output generated by DORiE has no meaning. The reason is that finite volumes are computed with finite elements of order 0 where gradients are 0.
......
......@@ -16,33 +16,70 @@ solution quantity.
Initial condition input is controlled entirely via the
:doc:`Configuration File <man-config-file>`.
.. note::
The initial condition is projected onto the actual solution function space.
Depending on grid shape and resolution, function space (order), and
interpolator (if applicable), the resulting solution may vary greatly from
the actual input data.
Input Types
-----------
This is an overview of all input types for initial conditions.
They are controlled by the ``initial.type`` key and available for every model.
.. object:: Linear Function
.. object:: Analytic
* ``type = analytic``
An analytic function :math:`f(\vec{p})` which depends on the physical
position :math:`\vec{p}`. The function must be defined via the key
``initial.equation``. For parsing the input expression, we use muparser_
which supports a set of common mathematical functions. Additionally, the
following variables can be used:
Available variables:
* ``x``: X-coordinate :math:`p_1 \, [\mathrm{m}]`.
* ``y``: Y-coordinate :math:`p_2 \, [\mathrm{m}]`.
* ``z``: Z-coordinate :math:`p_3 \, [\mathrm{m}]` (only in 3D).
* ``h``: Height above origin. Synonymous to ``y`` in 2D and ``z`` in 3D.
* ``pi``: Mathematical constant :math:`\pi`.
* ``dim``: Number of physical dimensions.
.. tip::
Assuming the target quantity is the matric head (see
:ref:`initial-transformation`), typical initial conditions for a
Richards simulation are
* ``type = linear``
* Hydrostatic equilibrium: A vertical gradient of :math:`-1` and a
fixed value ``<v>`` at height :math:`h = 0 \, \mathrm{m}`::
This type of initial condition defines an analytic linear function
:math:`f(h)` which only depends on height :math:`h`. The user specifies
the value at height zero, :math:`f_0 = f(h = 0 \, \mathrm{m})`, and a
gradient :math:`\Delta f`. The function defining the solution is then
given by
initial.equation = -h + <v>
.. math::
* Gravity flow: Constant value.
f (h) = f_0 + h \cdot \Delta f
.. tip::
The expression for a gaussian pulse of solute concentration centered at
:math:`\vec{p} = [0.5, 0.5]^T \, \mathrm{m}` is::
.. object:: HDF5_ Dataset
initial.equation = exp(-sqrt((x-0.5)^2+(y-0.5)^2)/(4.*0.002))/(4*pi*0.002)^(2/dim).
.. object:: Dataset
* ``type = data``
This type of initial condition defines function :math:`f(x,y,z)` depending
on the data included in the datafile ``initial.datafile`` and the dataset
``initial.dataset`` within the datafile. Values between the centers of the
data set are interpolated depending on the keyword ``initial.interpolator``
Load the initial condition from a data file ``initial.file`` by opening the
dataset ``initial.dataset`` inside this file. The data is interpreted as
function :math:`f(\mathbf{p})` of the physical position :math:`\mathbf{p}`
using one of the available :ref:`interpolators <man-interpolators>`, which
can be chosen using the setting ``initial.interpolation``. The input data
is automatically streched to match the grid extensions.
Supported file extensions:
* ``.h5``: H5_ data file. ``initial.dataset`` may be a file-internal path
to the target dataset.
.. _initial-transformation:
Transformation Types
--------------------
......@@ -58,15 +95,18 @@ Initial condition tranformations for the Richards solver.
* ``quantity = matricHead``
The input data is directly interpreted as matric head
:math:`h_m [\mathrm{m}]`.
:math:`f = h_m [\mathrm{m}]`.
Interpolators
-------------
Transport
^^^^^^^^^
Initial condition tranformations for the Transport solver.
.. object:: Nearest neighbor interpolator
.. object:: No Transformation
Interprets dataset values as cell values.
* ``quantity = soluteConcentration``
* ``interpolation: nearest``
The input data is directly interpreted as solute concentration,
:math:`f = c_w [\mathrm{kg}/\mathrm{m}^3]`.
.. _HDF5: https://www.h5py.org/
\ No newline at end of file
.. _H5: https://www.h5py.org/