Skip to content

GitLab

D
dorie
Commit 964606c4 authored by Lukas Riedel's avatar Lukas Riedel Committed by Santiago Ospina De Los Ríos
Browse files
Options

Reformulate TimeController and add unit test 

* Increase encapsulation of TimeController and reduce capabilities
* TimeController now only acts upon being called
* Remove validation routine
* Add unit test loading the default config time values

This commit breaks all models.
parent 2bf8045c
Hide whitespace changes
Inline Side-by-side
Showing with 679 additions and 345 deletions
cmake/modules/DorieTesting.cmake
View file @ 964606c4
......@@ -113,6 +113,13 @@ endfunction()
# Registers the created tests as unit tests, including coverage flags.
# If not specified, the tests are registered as system tests.
#
# .. cmake_param:: CUSTOM_MAIN
# :option:
#
# Write a custom `main()` function for the unit test executables instead
# of generating a default one automatically. Only applies if UNIT_TEST
# is enabled.
#
# .. cmake_param:: TARGET
# :single:
#
......@@ -158,7 +165,7 @@ endfunction()
# executable will be linked to the libraries DORiE depends on.
#
function(dorie_add_metaini_test)
set(OPTIONS UNIT_TEST)
set(OPTIONS UNIT_TEST CUSTOM_MAIN)
set(SINGLE TARGET METAINI SCRIPT BASENAME CREATED_TARGETS)
cmake_parse_arguments(SYSTEM_TEST "${OPTIONS}" "${SINGLE}" "" ${ARGN})
......@@ -166,6 +173,11 @@ function(dorie_add_metaini_test)
message(SEND_ERROR "No meta ini file given!")
endif()
if(SYSTEM_TEST_CUSTOM_MAIN AND NOT SYSTEM_TEST_UNIT_TEST)
message(WARNING "Ignoring option CUSTOM_MAIN because option UNIT_TEST "
"was not enabled")
endif()
# configure meta ini file or just copy.
get_filename_component(metaini-name ${SYSTEM_TEST_METAINI} NAME_WE)
get_filename_component(metaini-extension ${SYSTEM_TEST_METAINI} EXT)
......@@ -217,8 +229,14 @@ function(dorie_add_metaini_test)
# add dependencies and flags
if(SYSTEM_TEST_UNIT_TEST)
add_coverage_links(${created_targets})
# add_coverage_links(${created_targets})
add_dependencies(build_unit_tests ${created_targets})
if (SYSTEM_TEST_CUSTOM_MAIN)
target_link_libraries(${created_targets} gtest)
else ()
target_link_libraries(${created_targets} gtest_main)
endif()
else()
add_dependencies(build_system_tests ${created_targets})
endif()
......
dune/dorie/common/time_controller.hh deleted 100755 → 0
View file @ 2bf8045c
#ifndef DUNE_DORIE_CALCULATION_CONTROLLER_HH
#define DUNE_DORIE_CALCULATION_CONTROLLER_HH
#include <iostream>
#include <cmath>
#include <dune/common/exceptions.hh>
#include <dune/common/parametertree.hh>
#include <dune/dorie/common/logging.hh>
namespace Dune{
namespace Dorie{
/// Class to control time and maximum number of iterations for the Newton solver
template<typename R, typename TSC>
class CalculationController
{
private:
//! Typedef of logger for brevity
using Logger = std::shared_ptr<spdlog::logger>;
R time; //!< current time
R dt; //!< current time step
R suggested_dt; //!< suggestion for the next dt
bool use_dt_suggestion; //!< flag for suggested dt
int it; //!< current number of allowed iterations for Newton solver
//! logger of this instance
const Logger _log;
const R dtmax, dtmin; //!< time step limits
const R eps; //!< error margin
const R dtinc, dtdec; //!< time step mutliplicators
const R tBegin; //!< model time limits
R tEnd; //!< model time limits
const int itmax, itmin; //!< newton iteration limits
const TSC& tsc; //!< Class with Function getNextTimeStamp()
public:
/// Read relevant time values, check for consistency and adapt time step to BC change, if necessary
/** \param tsc_ Class with public function getNextTimeStamp(time)
* \param log The logger to use for this instance,
* defaults to retrieving the base logger.
*/
CalculationController(const ParameterTree& config,
const TSC& tsc_,
const Logger log=get_logger(log_base))
: time(config.get<R>("time.start"))
, dt(config.get<R>("time.startTimestep"))
, suggested_dt(dt)
, use_dt_suggestion(false)
, _log(log)
, dtmax(config.get<R>("time.maxTimestep"))
, dtmin(config.get<R>("time.minTimestep"))
, eps(dtmin/10.0)
, dtinc(config.get<R>("time.timestepIncreaseFactor"))
, dtdec(config.get<R>("time.timestepDecreaseFactor"))
, tBegin(config.get<R>("time.start"))
, tEnd(config.get<R>("time.end"))
, itmax(config.get<int>("time.maxIterations"))
, itmin(config.get<int>("time.minIterations"))
, tsc(tsc_)
{
if (not input_valid()) {
DUNE_THROW(IOError, "Invalid input for time controller");
}
// Check whether first time step needs adjustment
const R changeTime = tsc.getNextTimeStamp(time);
if(changeTime>0.0 && time+dt >= changeTime){
dt = changeTime-time-eps;
_log->trace("Adapting time step to change in boundary condition at "
"{:.2e}", changeTime);
}
it=calcIterations(dt);
}
/// suggest timestep
inline void suggest_timestep (R suggestion) {
if (FloatCmp::lt(suggestion,dt)) {
_log->trace("Adapting time step due to a suggestion: "
"{:.2e} -> {:.2e}", dt, suggestion);
dt = suggestion;
}
}
/// Return current time
inline R getTime () const { return time; }
/// Return timestep for next calculation
inline R getDT () const { return dt; }
/// Return begin time
inline R getBeginTime () const { return tBegin; }
/// Return end time
inline R getEndTime () const { return tEnd; }
/// Return maximum number of Newton iterations allowed for next calculation
inline int getIterations () const { return it; }
void set_time (const R _time)
{
time = _time;
}
/// Set a new finish time for the computation
void set_t_end (const R _t_end)
{
tEnd = _t_end;
}
/// Validation of solution and time step calculation
/** If the Newton solver throws an exception (because the solution did not converge
* in the allowed steps or the timestep is too large), exc=true is passed to this
* function and the timestep will be reduced. Returning false causes the Problem
* solver to recalculate the timestep with new parameters provided by the controller.
* The next time step is also adjusted to the next time the BC changes
* \param exc Boolean if Newton Class threw an exception
* \return Boolean if solution is valid
* \throw Dune::Exception Newton solver does not converge at dt>=dtmin and it<=itmax
*/
bool validate (const bool exc)
{
R changeTime = 0.0;
if (!exc) // solution was successfully computed
{
// time has advanced after calculation
time += dt;
if(time >= tEnd-eps) // Model has ended. Do nothing
return true;
dt = std::max(dt,dtmin);
dt = std::min(dt*dtinc,dtmax);
changeTime = tsc.getNextTimeStamp(time);
if(time >= changeTime-eps && time <= changeTime){
// small step to ensure that next time step starts with new BC
time = changeTime;
changeTime = tsc.getNextTimeStamp(time);
}
if(changeTime < 0.0 && time+dt > tEnd){
// time step adjustment to tEnd
dt = tEnd-time;
_log->trace("Adjusting time step to model end at: {:.2e}",
tEnd);
}
else if(changeTime > 0.0 && time+dt >= changeTime){
// time step adjustment to BC Change
dt = changeTime-time-eps;
_log->trace("Adapting time step to change in boundary condition at "
"{:.2e}", changeTime);
}
}
else // Newton solver threw error. Time is not advanced
{
if (dt <= dtmin) {
_log->error("Solver threw an exception at minimal time step "
"{:.2e}. Time step cannot be reduced", dt);
DUNE_THROW(Exception,"Solver fails at minimal time step");
}
auto failed_dt = dt;
dt = std::max(dt*dtdec,dtmin);
_log->trace("Adapting time step to improve solver: "
"{:.2e} -> {:.2e}", failed_dt, dt);
}
// Calculate new iterations based on adapted time step
it=calcIterations(dt);
return !exc;
}
private:
/// Calculate the allowed number of Newton iterations by logarithmic interpolation between max and min timestep
/** \param dt Timestep size
* \return Number of allowed Newton steps for inserted dt
*/
inline int calcIterations (const R dt)
{
return std::round( (itmin-itmax) * std::log10(dt-dtmin+1.0) / std::log10(dtmax-dtmin+1.0) + itmax );
}
/// Some checks if user input makes sense
/**
* \note This method does not throw exceptions by itself!
* \return True if the input variables are consistent
*/
bool input_valid () const
{
if (dtinc<=1.0) {
_log->error("timestepIncreaseFactor must be >1.0. Value is: {}",
dtinc);
return false;
}
if (dtdec>=1.0 || dtdec<=0.0) {
_log->error("timestepDecreaseFactor exceeds (0.0, 1.0). "
"Value is: {}", dtdec);
return false;
}
if (itmax<itmin) {
_log->error("maxIterations must be >=minIterations. "
"Values are: {}, {}", itmin, itmax);
return false;
}
if (dtmax<dtmin) {
_log->error("maxTimestep must be >=minTimestep. "
"Values are: {}, {}", dtmin, dtmax);
return false;
}
if (dt<dtmin || dt>dtmax) {
_log->error("startTimestep exceeds [minTimestep, maxTimestep]. "
"Value is: {}", dt);
return false;
}
return true;
}
};
} // namespace Dorie
} // namespace Dune
#endif
dune/dorie/common/time_step_controller.hh 0 → 100755
View file @ 964606c4
#ifndef DUNE_DORIE_TIME_CONTROLLER_HH
#define DUNE_DORIE_TIME_CONTROLLER_HH
#include <iostream>
#include <cmath>
#include <algorithm>
#include <dune/common/exceptions.hh>
#include <dune/common/parametertree.hh>
#include <dune/common/float_cmp.hh>
#include <dune/dorie/common/logging.hh>
namespace Dune{
namespace Dorie{
/// Provides access to time step information, as read from a config file
/** This class facilitates the management of time and time step related user
* input via config file settings. It does not control the actual simulation
* time or time step. Many settings can also be directly manipulated,
* if required.
*
* It is tailored for an implict time step scheme without specific criterion
* for a time step size. Time step suggestions are incorporated
* conservatively: Only suggestions below the current time step are accepted.
*
* The usual way of using this controller would be as follows:
* 1. Retrieve a time step suggestion from somewhere else.
* 2. Call suggest_time_step() with the suggestion. Depending on the
* current time step, the suggestion may or may not be accepted.
* 3. Adjust the time step to the simulation end by calling
* get_time_step() with the current simulation time as argument.
* 4. Perform a solver iteration with the given time step.
* 5. Increase the time step if the solver succeeded, and decrease it
* otherwise.
* 6. Use set_time_step() to choose any time step within the specified
* limits. These limits can optionally be ignored.
*
* \note set_time_step() and suggest_time_step() always return the new
* new time step. Programmers are free to ignore this value and
* retrieve it via get_time_step(), granting the option of adjusting
* the time step to the simulation end time.
*
* \tparam TF Time data type
* \ingroup Common
* \author Lukas Riedel
* \date 2019
*/
template<typename TF>
class TimeStepController
{
private:
using Logger = std::shared_ptr<spdlog::logger>;
/// The logger for this instance
const Logger _log;
const TF _time_begin; //!< Simulation start
TF _time_end; //!< Simulation end
const TF _time_step_min, //!< Minimum time step
_time_step_max, //!< Maximum time step
_time_step_inc, //!< Time step increment factor
_time_step_dec; //!< Time step decrement factor
/// Current time step
TF _time_step;
public:
/// Read settings from config file and check input parameters
/**
* \param config The configuration tree to read from
* \param log The logger for this instance
* \throw IOError if certain settings contradict each other
*/
explicit TimeStepController (const Dune::ParameterTree& config,
const Logger log=Dorie::get_logger(log_base)):
_log(log),
_time_begin(config.get<TF>("time.start")),
_time_end(config.get<TF>("time.end")),
_time_step_min(config.get<TF>("time.minTimestep")),
_time_step_max(config.get<TF>("time.maxTimestep")),
_time_step_inc(config.get<TF>("time.timestepIncreaseFactor")),
_time_step_dec(config.get<TF>("time.timestepDecreaseFactor")),
_time_step(config.get<TF>("time.startTimestep"))
{
if (not input_valid()) {
DUNE_THROW(IOError, "Invalid input for time controller");
}
}
/// Calculate a new time step from the current time and the suggestion
/** If the new time step exceeds the simulation end time (#_time_end),
* it is adjusted to hit it exactly.
*
* \param time Current time (required for adjustment to end time)
* \return The new time step
*/
TF get_time_step (const TF time)
{
using namespace Dune::FloatCmp;
if (ge(time, _time_end)) {
_log->warn("Current time exceeds the prescribed simulation "
"end time");
}
// adjust to the simulation end time
else if (ge(time + _time_step, _time_end)) {
_log->trace("Adjusting time step to simulation end time");
set_time_step(_time_end - time);
}
return get_time_step();
}
/// Return the currently stored time step
TF get_time_step () const { return _time_step; }
/// Return the minimum time step
TF get_time_step_min () const { return _time_step_min; }
/// Return the maximum time step
TF get_time_step_max () const { return _time_step_max; }
/// Return the simulation start time as stored in this class
TF get_time_begin () const { return _time_begin; }
/// Return the simulation end time as stored in this class
TF get_time_end () const { return _time_end; }
/// Set the time step to any value within the prescribed value range
/** If the value exceeds the range, a warning is issued and it is
* clamped.
*
* \param time_step New time step
* \param force Ignore time step limits when setting the value
* \return The new time step
*
* \warning Forcing the time step removes all sanity checks. In
* particular, negative time steps would be accepted.
*/
TF set_time_step (const TF time_step, const bool force=false)
{
if (not force && not time_step_valid(time_step)) {
_log->warn("Time step '{}' incompatible to chosen limits. "
"Adjusting automatically",
time_step);
_time_step = std::clamp(time_step, _time_step_min, _time_step_max);
}
else {
_time_step = time_step;
}
return get_time_step();
}
/// Suggest a new time step
/** Uses the minimum of the stored time step and the suggested time step.
*
* In the special case where the suggested time step and the current time
* step differ by less than the minimum time step, the resulting time
* step will be appropriately increased.
*
* \param time_step_suggestion External suggestion for a time step.
* \return The new time step
*
* \todo Add different policies on how to deal with suggested time steps
*/
TF suggest_time_step (const TF time_step_suggestion)
{
using namespace Dune::FloatCmp;
const auto time_step_new = std::min(_time_step,
time_step_suggestion);
// take suggestion to avoid time step smaller than dt_min
if (ge(time_step_new + _time_step_min, time_step_suggestion))
{
set_time_step(time_step_suggestion);
}
else {
set_time_step(time_step_new);
}
return get_time_step();
}
/// Set a new finish time of the simulation.
/** \param time_end The new end time
*/
void set_time_end (const TF time_end)
{
_time_end = time_end;
}
/// Decrease the time step by a certain decrement factor
/** \param time_step_decrement Decrement factor
*/
void decrease_time_step (const TF time_step_decrement)
{
_time_step = std::max(_time_step * time_step_decrement,
_time_step_min);
}
/// Decrease the time step by the default decrement factor
void decrease_time_step ()
{
decrease_time_step(_time_step_dec);
}
/// Increase the time step by a certain increment factor
/** \param time_step_increment Increment factor
*/
void increase_time_step (const TF time_step_increment)
{
_time_step = std::min(_time_step * time_step_increment,
_time_step_max);
}
/// Increase the time step by the default increment factor
void increase_time_step ()
{
increase_time_step(_time_step_inc);
}
private:
/// Some checks if user input makes sense
/**
* \note This method does not throw exceptions by itself!
* \return True if the input variables are consistent
*/
bool input_valid () const
{
if (_time_step_inc <= 1.0) {
_log->error("timestepIncreaseFactor must be >1.0. Value is: {}",
_time_step_inc);
return false;
}
if (_time_step_dec>=1.0 || _time_step_dec<=0.0) {
_log->error("timestepDecreaseFactor exceeds (0.0, 1.0). "
"Value is: {}", _time_step_dec);
return false;
}
if (_time_step_max < _time_step_min) {
_log->error("maxTimestep must be >=minTimestep. "
"Values are: {}, {}", _time_step_max, _time_step_min);
return false;
}
if (not time_step_valid(_time_step)) {
_log->error("Current time step exceeds [minTimestep, "
"maxTimestep]. Value is: {}",
_time_step);
return false;
}
return true;
}
/// Assert that a time step is within valid limits
/** \param time_step The time step to check
* \return True if the time step is valid
*/
bool time_step_valid (const TF time_step) const
{
if (time_step<_time_step_min || time_step>_time_step_max) {
return false;
}
return true;
}
};
} // namespace Dorie
} // namespace Dune
#endif // DUNE_DORIE_TIME_CONTROLLER_HH
dune/dorie/model/base.hh
View file @ 964606c4
......@@ -44,6 +44,7 @@ public:
)
: _output_policy(output_policy)
, _adapt_policy(adapt_policy)
, _time(0.0)
, _log(create_logger(log_name, helper, spdlog::level::from_str(log_level)))
{ }
......@@ -161,7 +162,7 @@ public:
*
* @return Current time of the model.
*/
virtual double current_time() const = 0;
virtual double current_time() const { return _time; }
/*-----------------------------------------------------------------------*//**
* @brief Suggest a time step to the model.
......@@ -207,6 +208,8 @@ private:
AdaptivityPolicy _adapt_policy;
protected:
/// The current time of the model
double _time;
//! The logger of this model
std::shared_ptr<spdlog::logger> _log;
};
......
dune/dorie/model/coupling/richards_coupling.cc
View file @ 964606c4
......@@ -97,7 +97,7 @@ void ModelRichardsTransportCoupling<Traits>::step()
if ( Dune::FloatCmp::ge(_time_transport, _transport->end_time()) )
{
_current_time = _time_richards;
this->_time = _time_richards;
return;
}
......@@ -109,7 +109,6 @@ void ModelRichardsTransportCoupling<Traits>::step()
// solve transport until synchronization with richards
while (Dune::FloatCmp::lt(_time_transport,_time_richards))
{
_time_transport = _transport->current_time();
_transport->suggest_timestep(_time_richards - _time_transport);
_transport->step();
// intermediate adapt if policy targets solute
......@@ -120,8 +119,10 @@ void ModelRichardsTransportCoupling<Traits>::step()
adapt_grid();
post_adapt_grid();
}
_time_transport = _transport->current_time();
}
_current_time = _time_richards;
this->_time = _time_richards;
// drop the old grid functions in transport model.
igf_water_content->pop();
......
dune/dorie/model/coupling/richards_coupling.hh
View file @ 964606c4
......@@ -89,7 +89,6 @@ private:
std::shared_ptr<typename Traits::BaseTraits::Grid> _grid;
Dune::ParameterTree _inifile_richards, _inifile_transport;
TimeField _current_time;
std::unique_ptr<ModelRichards> _richards;
std::unique_ptr<ModelTransport> _transport;
......@@ -169,7 +168,7 @@ public:
*
* @return Begin time of the model.
*/
TimeField begin_time() const override
double begin_time() const override
{
return _richards->begin_time();
}
......@@ -179,21 +178,11 @@ public:
*
* @return End time of the model.
*/
TimeField end_time() const override
double end_time() const override
{
return _richards->end_time();
}
/*------------------------------------------------------------------------*//**
* @brief Method that provides the current time of the model.
*
* @return Current time of the model.
*/