Commit 84595911 authored by Santiago Ospina De Los Ríos's avatar Santiago Ospina De Los Ríos

Merge branch '125-new-data-structure-for-storing-boundary-conditions' into 'master'

Resolve "New data structure for storing boundary conditions"

Closes #151, #136, #29, #126, and #125

See merge request !121
parents 59cf9b32 cef239ed
......@@ -36,12 +36,12 @@
* Split grid adaptivity process into marking and actual adaptation !91
* Rework VTK adapters and VTK output writer !64
* `RichardsSimulation` now implements the abstract simulation base class !89
* Swich to the stable release branch `releases/2.6` of `dune-testtools` !97
* Switch to the stable release branch `releases/2.6` of `dune-testtools` !97
* Merge Python packages into single new package `dorie` !100
* Move `dorie` Command Line Interface script into Python package !102
* Parameterization data input via YAML and H5 files !82
* Simplify H5 reader to only read datasets !109
* Extend run config file to contain data on multipe models !103
* Extend run config file to contain data on multiple models !103
* DORiE now writes vertex data by default. !128
* Switch license from MIT to GPLv3 !135
* Specifying scaling field `extensions` and `offset` is now optional !133
......@@ -53,6 +53,7 @@
* Build independent library and executable for each compile-time setting !144
* `SimulationBase` unit test now uses Google Test !159
* Deploy online documentation for each branch to private server !163
* Use YAML (instead of muPhi `.bcdat`) files for specifying BCs !121
### Fixed
* Allow meta-ini files for unit tests !101
......@@ -60,6 +61,8 @@
* Shape of input datasets was flipped when loading scaling factors !124
* `dune_add_system_test` requires target with location after bugfix !165
* `make all` would not build the solver application target `dorie` !167
* Allow Neumman BC to be applied on a different direction other than gravity !121
* Transport model option `dirichletMode` was not working properly !121
### Deprecated
......
......@@ -95,6 +95,11 @@ function(dorie_add_unit_test)
muparser::muparser hdf5 yaml-cpp spdlog)
# add_coverage_links(${UNIT_TEST_TARGET})
target_compile_definitions(${UNIT_TEST_TARGET}
PUBLIC
GTEST
)
if (UNIT_TEST_CUSTOM_MAIN)
target_link_libraries(${UNIT_TEST_TARGET} gtest)
else ()
......@@ -232,6 +237,11 @@ function(dorie_add_metaini_test)
# add_coverage_links(${created_targets})
add_dependencies(build_unit_tests ${created_targets})
target_compile_definitions(${created_targets}
PUBLIC
GTEST
)
if (SYSTEM_TEST_CUSTOM_MAIN)
target_link_libraries(${created_targets} gtest)
else ()
......
......@@ -7,8 +7,8 @@ configure_file(${CMAKE_CURRENT_BINARY_DIR}/tutorial-1.ini
COPYONLY)
# Copy files needed for the test to run
configure_file(${CMAKE_BINARY_DIR}/doc/default_files/2d_infiltr.bcdat
2d_infiltr.bcdat
configure_file(${CMAKE_BINARY_DIR}/doc/default_files/infiltration.yml
infiltration.yml
COPYONLY)
configure_file(${CMAKE_BINARY_DIR}/doc/default_files/richards_param.yml
richards_param.yml
......
......@@ -27,7 +27,7 @@ quantity = matricHead
equation = -h
[richards.boundary]
file = 2d_infiltr.bcdat
file = infiltration.yml
[richards.output]
fileName = infiltr_homogeneous_sand_2D
......
......@@ -121,7 +121,7 @@ Boundary Condition
The :doc:`boundary conditions file </manual/bcfile>` has to be specified by the
keyword ``richards.boundary.file`` in the configuration file. For now, we will
use the infiltration file for 2D, ``2d_infiltr.bcdat``, provided by the command
use the infiltration file, ``infiltration.yml``, provided by the command
``dorie create``. By default, this file sets a constant infiltration of
:math:`j_w = -5.55\cdot 10^{-6}\,\text{ms}^{-1}` at the top, a constant matric
head of :math:`h_m = 0\,\text{m}` at the bottom, and a no-flux condition on the
......@@ -194,6 +194,6 @@ shown below.
============= ======================================================================
Configuration :download:`config.ini <config.ini>`
Boundary :download:`2d_infiltr.bcdat </default_files/2d_infiltr.bcdat>`
Boundary :download:`infiltration.yml </default_files/infiltration.yml>`
Parameters :download:`richards_param.yml </default_files/richards_param.yml>`
============= ======================================================================
spatial_resolution_north 0
spatial_resolution_south 0
spatial_resolution_west -1
spatial_resolution_east -1
number_BC_change_times 1
0 neumann -5.55e-6 dirichlet 0
\ No newline at end of file
spatial_resolution_north 2 0.25 0.75
spatial_resolution_south -1
spatial_resolution_west -1
spatial_resolution_east -1
number_BC_change_times 2
0 neumann 0 dirichlet 1 neumann 0
1E3 neumann 0 neumann 0 neumann 0
spatial_resolution_north_we 0
spatial_resolution_north_fb 0
spatial_resolution_south_we 0
spatial_resolution_south_fb 0
spatial_resolution_west_sn -1
spatial_resolution_west_fb -1
spatial_resolution_east_sn -1
spatial_resolution_east_fb -1
spatial_resolution_front_sn -1
spatial_resolution_front_we -1
spatial_resolution_back_sn -1
spatial_resolution_back_we -1
number_BC_change_times 1
0 neumann -5.55e-6 dirichlet 0
\ No newline at end of file
spatial_resolution_north_we 2 0.25 0.75
spatial_resolution_north_fb 2 0.25 0.75
spatial_resolution_south_we -1
spatial_resolution_south_fb -1
spatial_resolution_west_sn -1
spatial_resolution_west_fb -1
spatial_resolution_east_sn -1
spatial_resolution_east_fb -1
spatial_resolution_front_sn -1
spatial_resolution_front_we -1
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
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,9 +117,10 @@ function(create_default_config INPUT OUTPUT SOURCE_DIR CSS)
endfunction()
# copy BC & parameter files
file(COPY 2d_infiltr.bcdat 3d_infiltr.bcdat
2d_solute.bcdat 3d_solute.bcdat
richards_param.yml transport_param.yml
file(COPY infiltration.yml
solute.yml
richards_param.yml
transport_param.yml
DESTINATION .)
# Random field generator
......
upper:
index: 1
conditions:
heavy_rainfall:
type: Neumann
time: 0
value: -5.55e-6
horizontal_projection: true
lower:
index: 0
conditions:
water_table:
type: Dirichlet
time: 0
value: 0
default:
type: Neumann
value: 0
......@@ -90,25 +90,11 @@ adding an empty line, make text **bold** or ``monospaced``.
<category name="boundary">
<parameter name="file">
<definition> Path to the boundary condition file. </definition>
<suggestion> 2d_infiltr.bcdat </suggestion>
<definition> Path to the YAML boundary condition data file.
<suggestion> infiltration.yml </suggestion>
</definition>
<values> path </values>
</parameter>
<parameter name="fileType">
<definition> Type of spatial segmentation of the boundaries specified in the BC file </definition>
<values> rectangularGrid </values>
<suggestion> rectangularGrid </suggestion>
<comment> Choose type of boundary segmentation: rectangularGrid </comment>
</parameter>
<parameter name="interpolateBCvalues">
<definition> Whether to interpolate between the boundary conditions
at different times linearly (``true``) or not at all (``false``). May require
different boundary condition files. </definition>
<values> true, false </values>
<suggestion> false </suggestion>
</parameter>
</category>
<category name="initial">
......
upper:
index: 1
conditions:
tracer:
type: Dirichlet
time: 0
value: 1
concentration_type: water_phase
default:
type: Outflow
value: 0
......@@ -101,30 +101,10 @@ adding an empty line, make text **bold** or ``monospaced``.
<category name="boundary">
<parameter name="file">
<definition> Path to the boundary condition file. </definition>
<values> path </values>
</parameter>
<parameter name="fileType">
<definition> Type of spatial segmentation of the boundaries specified in the BC file </definition>
<values> rectangularGrid </values>
<suggestion> rectangularGrid </suggestion>
<comment> Choose type of boundary segmentation: rectangularGrid </comment>
</parameter>
<parameter name="interpolateBCvalues">
<definition> Whether to interpolate between the boundary conditions
at different times linearly (``true``) or not at all (``false``). May require
different boundary condition files. </definition>
<values> true, false </values>
<suggestion> false </suggestion>
</parameter>
<parameter name="dirichletMode">
<definition> Type of the input value for the dirichlet condition.
<definition> Path to the YAML boundary condition data file.
<suggestion> solute.yml </suggestion>
</definition>
<values> soluteConcentration, totalSolute </values>
<suggestion> soluteConcentration </suggestion>
<values> path </values>
</parameter>
</category>
......
spatial_resolution_north 2 0.4 0.8
spatial_resolution_south 0
spatial_resolution_west -1
spatial_resolution_east -1
number_BC_change_times 1
0 neumann 0.0 neumann -5.55e-6 neumann 0.0 dirichlet 0
\ No newline at end of file
spatial_resolution_north_we 2 0.4 0.6
spatial_resolution_north_fb 2 0.4 0.6
spatial_resolution_south_we 0
spatial_resolution_south_fb 0
spatial_resolution_west_sn -1
spatial_resolution_west_fb -1
spatial_resolution_east_sn -1
spatial_resolution_east_fb -1
spatial_resolution_front_sn -1
spatial_resolution_front_we -1
spatial_resolution_back_sn -1
spatial_resolution_back_we -1
number_BC_change_times 1
0 neumann 0.0 neumann 0.0 neumann 0.0 neumann 0.0 neumann -5.55e-6 neumann 0.0 neumann 0.0 neumann 0.0 neumann 0.0 dirichlet 0
\ No newline at end of file
......@@ -27,10 +27,8 @@ Templates for the required files can conveniently be generated with the
| :ref:`Parameterization Data File | :file:`.yml` | Specifies the soil and solute | :download:`richards_param.yml </default_files/richards_param.yml>`, |
| <man-parameter_file>` | | parameterization type(s) and parameters. | :download:`transport_param.yml </default_files/transport_param.yml>` |
+----------------------------------------+----------------+------------------------------------------+----------------------------------------------------------------------+
| :doc:`Boundary Condition Data File | :file:`.bcdat` | Specifies the boundary segmentation and | :download:`2d_infiltr.bcdat </default_files/2d_infiltr.bcdat>`, |
| </manual/bcfile>` | | the boundary conditions. | :download:`3d_infiltr.bcdat </default_files/3d_infiltr.bcdat>`, |
| | | | :download:`2d_solute.bcdat </default_files/2d_solute.bcdat>`, |
| | | | :download:`3d_solute.bcdat </default_files/3d_solute.bcdat>` |
| :doc:`Boundary Condition Data File | :file:`.yml` | Specifies the boundary segmentation and | :download:`infiltration.yml </default_files/infiltration.yml>`, |
| </manual/bcfile>` | | the boundary conditions. | :download:`solute.yml </default_files/solute.yml>` |
+----------------------------------------+----------------+------------------------------------------+----------------------------------------------------------------------+
| :ref:`GMSH Mesh File | :file:`.gmsh` | Mesh with possibly complicated geometry | --- |
| <man-grid_gmsh>` (optional) | | to build an unstructured grid from. | |
......
This diff is collapsed.
......@@ -68,6 +68,7 @@ successful. You can deactivate the ``venv`` at any time by executing
You are now ready to use the DORiE Command Line Interface!
.. _man-cli-commands:
CLI Commands
============
......
#ifndef DUNE_DORIE_COMMON_BOUNDARY_CONDITION_BASE_HH
#define DUNE_DORIE_COMMON_BOUNDARY_CONDITION_BASE_HH
#include <string>
#include <string_view>
#include <optional>
#include <memory>
#include <dune/common/exceptions.hh>
#include <dune/dorie/common/logging.hh>
#include <dune/dorie/common/typedefs.hh>
#include <dune/common/float_cmp.hh>
namespace Dune {
namespace Dorie {
/// A generic boundary condition storage
/** Stores the boundary condition time interval and values,
* and interpolates them linearly
*
* \tparam RF Type of the boundary codition value
* \ingroup Boundary
*/
template<typename RF>
class BoundaryCondition
{
public:
/// The type used for identifying time
using Time = double;
protected:
/// The time interval in which the BC is applied
const TimeInterval<Time> _time_interval;
/// The starting (or overall) value
const RF _value_begin;
/// The optional final value
const std::optional<RF> _value_end{};
/// The name given by the user in the data file
const std::string _name;
/// Check if the given time is inside the time interval
/** \param time The point in time to check
* \return True if the time is inside the interval of this BC
*/
bool in_interval (const Time time) const
{
const auto& t_begin = _time_interval.begin;
const auto& t_end = _time_interval.end;
using namespace Dune::FloatCmp;
if (lt(time, t_begin) or gt(time, t_end))
{
auto log = Dorie::get_logger(Dorie::log_base);
log->error("Querying boundary condition '{}' outside its time "
"interval [{}, {}] with time: {}",
_name, t_begin, t_end, time);
return false;
}
return true;
}
public:
/// Construct a static boundary condition
BoundaryCondition (const std::string_view name,
const TimeInterval<Time> interval,
const RF value)
:
_time_interval(interval),
_value_begin(value),
_name(name)
{ }
/// Construct an interpolated boundary condition
BoundaryCondition (const std::string_view name,
const TimeInterval<Time> interval,
const RF value_begin,
const std::optional<RF> value_end)
:
_time_interval(interval),
_value_begin(value_begin),
_value_end(value_end),
_name(name)
{ }
/// Destroy this object
virtual ~BoundaryCondition () = default;
/// Return the name of this boundary condition
std::string_view name () const { return _name; }
/// Return the time interval of this boundary condition
const TimeInterval<Time>& time_interval () const
{
return _time_interval;
}
/// Return the type of this boundary condition
virtual Operator::BCType type () const
{
return Operator::BCType::none;
}
/// Return the value for a given point in time
/** \param time The time at which to evaluate the BC
* \return The (possibly interpolated) value
* \throw InvalidStateException if the time is not inside the interval
* \todo Maybe use `assert` on in_interval
*/
virtual RF evaluate (const Time time) const
{
if (not in_interval(time)) {
DUNE_THROW(InvalidStateException, "Evaluating boundary condition "
"with invalid time");
}
// without interpolation, return the first value
if (not _value_end)
return _value_begin;
// else: interpolate values
const auto& t_begin = _time_interval.begin;
const auto& t_end = _time_interval.end;
const auto& v_end = _value_end.value();
return _value_begin
+ (time - t_begin)
* (v_end - _value_begin)
/ (t_end - t_begin);
}
};
} // namespace Dorie
} // namespace Dune
#endif // DUNE_DORIE_COMMON_BOUNDARY_CONDITION_BASE_HH
#ifndef DUNE_DORIE_COMMON_BOUNDARY_CONDITION_DIRICHLET_HH
#define DUNE_DORIE_COMMON_BOUNDARY_CONDITION_DIRICHLET_HH
#include <dune/dorie/common/boundary_condition/boundary_condition.hh>
namespace Dune {
namespace Dorie {
/// Toggle for concentration represented by Dirichlet BC
/** Used for \ref TransportModel model only.
*
* The relation of water_phase and total is given by
* \f[
* C_w = C_t \theta_w \, ,
* \f]
* where \f$ \theta_w \f$ is the volumetric water content.
*/
enum class SoluteConcentration {
total, //!< Total concentration in the volume \f$ C_t \f$.
water_phase, //!< Concentration in water phase \f$ C_w \f$.
none //!< No type, used for \ref RichardsModel model.
};
/// Class representing a Dirichlet boundary condition (fixed solution)
/** In the \ref TransportModel model, the #_solute_concentration switch is used
* to determine which type of solute concentration the boundary condition
* refers to.
*
* \ingroup Boundary
*/
template<typename RF>
class DirichletBoundaryCondition : public BoundaryCondition<RF>
{
protected:
using Base = BoundaryCondition<RF>;
using Time = typename Base::Time;
/// Setting for the type of concentration encoded in this BC value
const SoluteConcentration _solute_concentration;
public:
/// Construct a static boundary condition
DirichletBoundaryCondition (
const std::string_view name,
const TimeInterval<Time> interval,
const RF value,
const SoluteConcentration concentration=SoluteConcentration::none)
:
Base(name, interval, value),
_solute_concentration(concentration)
{ }
/// Construct an interpolated boundary condition
DirichletBoundaryCondition (
const std::string_view name,
const TimeInterval<Time> interval,
const RF value_begin,
const std::optional<RF> value_end,
const SoluteConcentration concentration=SoluteConcentration::none)
:
Base(name, interval, value_begin, value_end),
_solute_concentration(concentration)
{ }
/// Destroy this object
virtual ~DirichletBoundaryCondition () override = default;
/// Return the type of this boundary condition
virtual Operator::BCType type () const override
{
return Operator::BCType::Dirichlet;
}
/// The type string of this boundary condition class
static constexpr std::string_view type_str = "dirichlet";
/// Return the quantity encoded in this Dirichlet value
SoluteConcentration concentration_type () const
{
return _solute_concentration;
}
};
} // namespace Dorie
} // namespace Dune
#endif // DUNE_DORIE_COMMON_BOUNDARY_CONDITION_DIRICHLET_HH
This diff is collapsed.
This diff is collapsed.
// These groups should gather information on strongly related portions of
// the code. In them we aim to explain to other developers:
// - The rationale of the design
// - How careful someone has to be when modifying it
// - Specific naming or other conventions upheld in the code
/**
@defgroup Boundary Boundary Conditions
@{
@ingroup Common
@brief Representing boundary conditions for all models
## Overview
There are several types of boundary conditions. They all store the time
interval during which they can be evaluated. Typically, one static value,
or two values for interpolation, are stored. Since boundary conditions are
implemented as derived classed of Dune::Dorie::BoundaryCondition, they
may also store additional members and define additional methods (see
Dune::Dorie::NeumannBoundaryCondition for an example).
## Implementation
Boundary conditions are handled by the
Dune::Dorie::BoundaryConditionManager. It loads the user input data and
assembles the necessary data structures. Its behavior is the same for all
models, with the difference that models may not support all types of
boundary conditions. Therefore, the manager has a template parameter for
setting its model mode, which is propagated to the
Dune::Dorie::BoundaryConditionFactory. New boundary conditions have to be
"registered" there, making sure they can be used for a certain (or all)
models.
@}
**/
#ifndef DUNE_DORIE_COMMON_BOUNDARY_CONDITION_NEUMANN_HH
#define DUNE_DORIE_COMMON_BOUNDARY_CONDITION_NEUMANN_HH
#include <dune/dorie/common/boundary_condition/boundary_condition.hh>
namespace Dune {
namespace Dorie {
/// Class representing a Neumann boundary condition (fixed flux)
/** In the \ref RichardsModel model, this class typically models irrigation and
* evaporation from soil surfaces and therefore stores an additional switch
* which can be used to consider the surface tilt when evaluating the boundary
* condition.
*
* If #_horizontal_projection is set to true, the local operator scales the
* boundary condition \f$ j \f$ with
* \f[
* j = j_N \left| \mathbf{\hat{n}}_F \cdot \mathbf{\hat{g}}_F \right|
* \f]
* where the vectors are the unit normal vector of the boundary face, and
* the unit vector along the gravitational force, respectively.
*
* #_horizontal_projection is ignored in the \ref TransportModel model.
*
* \ingroup Boundary
*/
template<typename RF>
class NeumannBoundaryCondition : public BoundaryCondition<RF>
{
protected:
using Base = BoundaryCondition<RF>;
using Time = typename Base::Time;
/// Value of the surface adjustment switch
const bool _horizontal_projection;
public:
/// Construct a static boundary condition
NeumannBoundaryCondition (const std::string_view name,
const TimeInterval<Time> interval,
const RF value,
const bool horizontal_projection)
:
Base(name, interval, value),
_horizontal_projection(horizontal_projection)
{ }
/// Construct an interpolated boundary condition
NeumannBoundaryCondition (const std::string_view name,
const TimeInterval<Time> interval,
const RF value_begin,
const std::optional<RF> value_end,
const bool horizontal_projection)
:
Base(name, interval, value_begin, value_end),
_horizontal_projection(horizontal_projection)
{ }
/// Destroy this object
virtual ~NeumannBoundaryCondition () override = default;
/// Return the type of this boundary condition
virtual Operator::BCType type () const override
{
return Operator::BCType::Neumann;
}
/// Return true the BC is to be projected horizontally
bool horizontal_projection () const { return _horizontal_projection; }
/// The type string of this boundary condition class
static constexpr std::string_view type_str = "neumann";
};
} // namespace Dorie
} // namespace Dune
#endif // DUNE_DORIE_COMMON_BOUNDARY_CONDITION_NEUMANN_HH
#ifndef DUNE_DORIE_COMMON_BOUNDARY_CONDITION_OUTFLOW_HH
#define DUNE_DORIE_COMMON_BOUNDARY_CONDITION_OUTFLOW_HH
#include <dune/dorie/common/boundary_condition/boundary_condition.hh>
namespace Dune {
namespace Dorie {
/// Class representing an Outflow BC (matter may leave through the boundary)
/** The exact implementation depends on the model.
*
* \ingroup Boundary
*/
template<typename RF>
class OutflowBoundaryCondition : public BoundaryCondition<RF>
{
protected:
using Base = BoundaryCondition<RF>;
using Time = typename Base::Time;
public:
/// Re-use base constructor
using Base::Base;
/// Destroy this object
virtual ~OutflowBoundaryCondition () override = default;
/// Return the type of this boundary condition
virtual Operator::BCType type () const override
{
return Operator::BCType::Outflow;
}
/// The type string of this boundary condition class
static constexpr std::string_view type_str = "outflow";
};
} // namespace Dorie
} // namespace Dune
#endif // DUNE_DORIE_COMMON_BOUNDARY_CONDITION_OUTFLOW_HH