Commit 5ce20ffc authored by Lukas Riedel's avatar Lukas Riedel

Merge branch 'cherry-pick-dac2a4e1' into '46-rework-system-tests'

Merge branch '43-add-options-for-altering-dg-formulation-in-local-operator' into 'master'

See merge request !34

cherry picked from commit dac2a4e1
parents 380e9be6 118d9bce
......@@ -11,6 +11,12 @@ if(CMAKE_BUILD_TYPE_UPPER MATCHES DEBUG)
set(CMAKE_VERBOSE_MAKEFILE ON)
endif()
# option to change DG scheme via config file
option(EXPERIMENTAL_DG_FEATURES
"Enable experimental DG settings through the config file"
OFF
)
#
if(NOT (dune-common_DIR OR dune-common_ROOT OR
"${CMAKE_PREFIX_PATH}" MATCHES ".*dune-common.*"))
......
......@@ -143,6 +143,15 @@ If you installed [Anaconda](https://conda.io/docs/user-guide/install/download.ht
**Installing DUNE and DORiE via `make install` into your usual work environment is strongly discouraged!**
### Experimental Features
The local operator implementing Richards equation's discretization supports
multiple scheme settings. Setting these via the config file is disabled by
default. You can enable this feature by reconfiguring DORiE with the CMake flag
`-DEXPERIMENTAL_DG_FEATURES=On`, and rebuilding it.
The configuration settings in the section `[dg.experimental]` will then override
the default settings.
## Recommended Third-Party Software
......
......@@ -412,4 +412,41 @@ adding an empty line, make text **bold** or ``monospaced``.
</parameter>
</category>
<category name="dg.experimental">
<parameter name="method">
<definition> DG discretization method for skeleton terms.
**SIPG:** Symmetric Interior Penalty
**NIPG:** Non-Symmetric Interior Penalty
**OOB:** Oden, Babuska, Baumann: no penalty term
**IIP:** Incomplete Interior Penalty: no symmetry term
</definition>
<values> SIPG, NIPG, OOB, IIP </values>
<suggestion> SIPG </suggestion>
<comment> Experimental settings are enabled by the appropriate CMake flag.
</comment>
</parameter>
<parameter name="upwinding">
<definition> Upwinding method for skeleton terms.
**semiUpwind:** Apply upwinding to conductivity factor (only).
**fullUpwind:** Apply upwinding on numeric flux and conductivity.
</definition>
<values> none, semiUpwind, fullUpwind </values>
<suggestion> none </suggestion>
</parameter>
<parameter name="weights">
<definition> Apply harmonic weighting to skeleton term contributions.
</definition>
<values> true, false </values>
<suggestion> true </suggestion>
</parameter>
</category>
</dorie>
......@@ -3,4 +3,9 @@ add_subdirectory(solver)
add_subdirectory(impl)
add_executable("dorie" dorie.cc)
dune_target_link_libraries(dorie dorie-impl ${DUNE_LIBS})
\ No newline at end of file
dune_target_link_libraries(dorie dorie-impl ${DUNE_LIBS})
# enable setting operator scheme from config file
if(EXPERIMENTAL_DG_FEATURES)
target_compile_definitions("dorie" PUBLIC -DEXPERIMENTAL_DG_FEATURES)
endif()
\ No newline at end of file
......@@ -42,7 +42,7 @@ private:
/// Adaptivity GFS
using AGFS = typename AGFSHelper::Type;
/// Error estimator local operator
using ESTLOP = Dune::Dorie::FluxErrorEstimator<Traits,Param,Boundary>;
using ESTLOP = Dune::Dorie::Operator::FluxErrorEstimator<Traits,Param,Boundary>;
/// Empty constraints container
using NoTrafo = Dune::PDELab::EmptyTransformation;
/// Solution vector type
......
......@@ -31,7 +31,20 @@ Simulation<Traits>::Simulation (Dune::MPIHelper& _helper, std::shared_ptr<Grid>
finitial = std::make_unique<FlowInitial>(inifile,*param,gv);
// --- Local Operators ---
#ifdef EXPERIMENTAL_DG_FEATURES
// read experimental settings from inifile
namespace OP = Dune::Dorie::Operator;
const auto settings = OP::read_experimental_operator_settings(inifile);
const auto method = std::get<OP::RichardsDGMethod::Type>(settings);
const auto upwinding = std::get<OP::RichardsDGUpwinding::Type>(settings);
const auto weights = std::get<OP::RichardsDGWeights::Type>(settings);
slop = std::make_unique<SLOP>(inifile, *param, gv, *fboundary, *fsource,
method, upwinding, weights);
#else
slop = std::make_unique<SLOP>(inifile,*param,gv,*fboundary,*fsource);
#endif // EXPERIMENTAL_DG_FEATURES
tlop = std::make_unique<TLOP>(inifile,*param,gv);
controller = std::make_unique<CalculationController>(inifile,*fboundary,helper);
......@@ -128,6 +141,7 @@ bool Simulation<Traits>::compute_time_step ()
const RF t = controller->getTime();
const RF dt = controller->getDT();
bool exception = false;
try
{
// solve in parallel
......@@ -140,7 +154,7 @@ bool Simulation<Traits>::compute_time_step ()
// solve sequentially
else {
pdesolver_seq->setMaxIterations(controller->getIterations());
osm_seq->apply(t,dt,*uold,*unew);
osm_seq->apply(t, dt, *uold, *unew);
}
*uold = *unew;
......
This diff is collapsed.
......@@ -12,6 +12,7 @@
namespace Dune {
namespace Dorie {
namespace Operator {
/// Local operator for residual-based error estimation.
/*
......@@ -48,18 +49,23 @@ namespace Dune {
const Parameter& param;
const Boundary& boundary;
const RF penalty_factor;
const int intorderadd;
const int quadrature_factor;
RF penalty_factor;
RF time;
FluxErrorEstimator (const Dune::ParameterTree& config,
const Parameter& param_, const Boundary& boundary_,
RichardsDGMethod::Type method_ = RichardsDGMethod::SIPG,
int intorderadd_ = 2, int quadrature_factor_ = 2)
: param(param_), boundary(boundary_),
intorderadd(intorderadd_), quadrature_factor(quadrature_factor_),
penalty_factor(config.get<RF>("dg.penaltyFactor")),
intorderadd(intorderadd_), quadrature_factor(quadrature_factor_)
{ }
time(0.0)
{
if (method_ == RichardsDGMethod::OOB)
penalty_factor = 0.0;
}
/// skeleton integral depending on test and ansatz functions
/** each face is only visited ONCE!
......@@ -330,7 +336,8 @@ namespace Dune {
};
}
}
} // namespace Operator
} // namespace Dorie
} // namespace Dune
#endif // DUNE_DORIE_ERROR_INDICATOR_HH
\ No newline at end of file
......@@ -22,8 +22,8 @@ namespace Dune {
explicit MualemVanGenuchten(const ParameterTree& fieldConfig_, const Dune::MPIHelper& helper_,
Interpolator& interpolator_, H5File& h5file_, const int verbose_)
: PB(fieldConfig_, "vanGenuchten", parameterArray, helper_, interpolator_, h5file_, verbose_)
, thetaR("theta_r", 0.0, 1.0), thetaS("theta_s", 0.0, 1.0), alpha("alpha", -10.0, 0.0)
, n("n", 1.0, 100.0), tau("tau", -100.0, 0.0), k0("k0", 0.0, 1.0)
, thetaR("theta_r", 0.0, 1.0), thetaS("theta_s", 0.0, 1.0), alpha("alpha", -20.0, 0.0)
, n("n", 1.0, 100.0), tau("tau", -100.0, 1.0), k0("k0", 0.0, 1.0)
{
parameterArray.push_back(& thetaR);
parameterArray.push_back(& thetaS);
......
Markdown is supported
0% or
You are about to add 0 people to the discussion. Proceed with caution.
Finish editing this message first!
Please register or to comment