OGS 6.3.0-179-g962fdcd4e.dirty.20200403132553
ODE Solver Library

Detailed Description

This ODE solver library has been designed with an implicit first-order quasilinear ODE in mind. However, it is in principle not restricted to such a kind of equation, but can be extended to also solve other equation types. In particular it is possible to introduce equation types that are no ODEs, but are nonlinear equations without time derivative terms.

The aim of the library's design is being able to formulate FEM processes without having to care which time discretization scheme and which nonlinear iteration method will be used to solve it.

The library offers different time discretization schemes, cf. the conceputal remarks on them, namely the forward and backward Euler and Crank-Nicolson methods and Backward Differentiation Formulas. The design follows Gear's method, which also underlies the DASSL algorithm of Petzold et al., cf. Differential-algebraic equations article, section Numerical methods/Direct discretization [3].

A rough overview over the interplay between the various parts of this library is given in the image below. Therein red symbols indicate which classes own which matrices or vectors (for the meaning of the different symbols refer to the documentation of the respective classes). The word own should not be taken too strict in the sense of C++ member ownership, although currently it is implemented as such; but in the future this implementation detail might change. Rather own means that the class is in charge of the respective matrix or vector, i.e., it can read from and write to it or pass it on to functions; or in other words: The class knows about the meaning of that matrix/vector. Note that only matrices and vectors that describe some proper state of the respective classes are shown; those storing only intermediate computations have been omitted from the image.

Interplay of the different parts of the ODE solver library at the example of a first-order implicit quasilinear ODE.
In the ODE solver library the instances of some classes can work together with several instances of other classes throughout there lifetime, whereas they have a strict one-to-one correspondence to objects of some different type. Those relations are given in the following table.

Class 1 Class 2 Relation Remarks
TimeDiscretizedODESystem ODESystem 1:1 the ODESystem represents part of the state of the TimeDiscretizedODESystem
TimeDiscretizedODESystem TimeDiscretization 1:1 analogous for the TimeDiscretization
TimeDiscretizedODESystem MatrixTranslator 1:1 analogous for the MatrixTranslator
NonlinearSolver NonlinearSystem 1:n a nonlinear solver can solve various equations, one after the other
LinearSolver NonlinearSolver 1:n various NonlinearSolver's can share the same LinearSolver


class  NumLib::EquationSystem
class  NumLib::MatrixTranslator< ODETag >
class  NumLib::MatrixTranslator< ODESystemTag::FirstOrderImplicitQuasilinear >
class  NumLib::MatrixTranslatorGeneral< ODETag >
class  NumLib::MatrixTranslatorGeneral< ODESystemTag::FirstOrderImplicitQuasilinear >
class  NumLib::MatrixTranslatorForwardEuler< ODETag >
class  NumLib::MatrixTranslatorForwardEuler< ODESystemTag::FirstOrderImplicitQuasilinear >
class  NumLib::MatrixTranslatorCrankNicolson< ODETag >
class  NumLib::MatrixTranslatorCrankNicolson< ODESystemTag::FirstOrderImplicitQuasilinear >
class  NumLib::NonlinearSolver< NLTag >
class  NumLib::NonlinearSolver< NonlinearSolverTag::Newton >
class  NumLib::NonlinearSolver< NonlinearSolverTag::Picard >
class  NumLib::NonlinearSystem< NLTag >
class  NumLib::NonlinearSystem< NonlinearSolverTag::Newton >
class  NumLib::NonlinearSystem< NonlinearSolverTag::Picard >
class  NumLib::ODESystem< ODETag, NLTag >
class  NumLib::ODESystem< ODESystemTag::FirstOrderImplicitQuasilinear, NonlinearSolverTag::Picard >
class  NumLib::ODESystem< ODESystemTag::FirstOrderImplicitQuasilinear, NonlinearSolverTag::Newton >
class  NumLib::InternalMatrixStorage
class  NumLib::TimeDiscretization
class  NumLib::BackwardEuler
 Backward Euler scheme. More...
class  NumLib::ForwardEuler
 Forward Euler scheme. More...
class  NumLib::CrankNicolson
 Generalized Crank-Nicolson scheme. More...
class  NumLib::BackwardDifferentiationFormula
 Backward differentiation formula. More...
class  NumLib::TimeDiscretizedODESystemBase< NLTag >
class  NumLib::TimeDiscretizedODESystem< ODETag, NLTag >
class  NumLib::TimeDiscretizedODESystem< ODESystemTag::FirstOrderImplicitQuasilinear, NonlinearSolverTag::Newton >
class  NumLib::TimeDiscretizedODESystem< ODESystemTag::FirstOrderImplicitQuasilinear, NonlinearSolverTag::Picard >


enum  NumLib::IterationResult : char { NumLib::IterationResult::SUCCESS, NumLib::IterationResult::FAILURE, NumLib::IterationResult::REPEAT_ITERATION }
 Status flags telling the NonlinearSolver if an iteration succeeded. More...
enum  NumLib::NonlinearSolverTag : bool { NumLib::NonlinearSolverTag::Picard, NumLib::NonlinearSolverTag::Newton }
 Tag used to specify which nonlinear solver will be used. More...
enum  NumLib::ODESystemTag : char { NumLib::ODESystemTag::FirstOrderImplicitQuasilinear, NumLib::ODESystemTag::NeumannBC }
 Tag used to specify the type of ODE. More...


template<ODESystemTag ODETag>
std::unique_ptr< MatrixTranslator< ODETag > > NumLib::createMatrixTranslator (TimeDiscretization const &timeDisc)
std::pair< std::unique_ptr< NonlinearSolverBase >, NonlinearSolverTagNumLib::createNonlinearSolver (GlobalLinearSolver &linear_solver, BaseLib::ConfigTree const &config)

Enumeration Type Documentation

◆ IterationResult

enum NumLib::IterationResult : char

Status flags telling the NonlinearSolver if an iteration succeeded.


Definition at line 21 of file EquationSystem.h.

◆ NonlinearSolverTag

enum NumLib::NonlinearSolverTag : bool

Tag used to specify which nonlinear solver will be used.


Picard fixpoint iteration scheme


Newton-Raphson iteration scheme

Definition at line 20 of file Types.h.

◆ ODESystemTag

enum NumLib::ODESystemTag : char

Tag used to specify the type of ODE.


First order implicit quasi-linear ODE

This is an ODE of the form \( M(x,t)\cdot \dot x + K(x,t) \cdot x - b(x,t) =: r(\dot x, x, t) \stackrel{!}{=} 0 \)


Definition at line 26 of file Types.h.

26  : char
27 {
35  NeumannBC // Sure, that's misuse of this enum, so sue me!
36 };

Function Documentation

◆ createMatrixTranslator()

template<ODESystemTag ODETag>
std::unique_ptr<MatrixTranslator<ODETag> > NumLib::createMatrixTranslator ( TimeDiscretization const &  timeDisc)

Creates a GlobalMatrix translator suitable to work together with the given time discretization scheme.

Definition at line 291 of file MatrixTranslator.h.

293 {
294  if (auto* fwd_euler = dynamic_cast<ForwardEuler const*>(&timeDisc))
295  {
296  return std::unique_ptr<MatrixTranslator<ODETag>>(
297  new MatrixTranslatorForwardEuler<ODETag>(*fwd_euler));
298  }
299  if (auto* crank = dynamic_cast<CrankNicolson const*>(&timeDisc))
300  {
301  return std::unique_ptr<MatrixTranslator<ODETag>>(
302  new MatrixTranslatorCrankNicolson<ODETag>(*crank));
303  }
305  return std::unique_ptr<MatrixTranslator<ODETag>>(
306  new MatrixTranslatorGeneral<ODETag>(timeDisc));
307 }

◆ createNonlinearSolver()

std::pair< std::unique_ptr< NonlinearSolverBase >, NonlinearSolverTag > NumLib::createNonlinearSolver ( GlobalLinearSolver &  linear_solver,
BaseLib::ConfigTree const &  config 

Creates a new nonlinear solver from the given configuration.

linear_solverthe linear solver that will be used by the nonlinear solver
configconfiguration settings
a pair (nl_slv, tag) where nl_slv is the generated nonlinear solver instance and the tag indicates if it uses the Picard or Newton-Raphson method
Input File Parameter:
Input File Parameter:
Input File Parameter:

Definition at line 417 of file NonlinearSolver.cpp.

References BaseLib::ConfigTree::getConfigParameter(), NumLib::Newton, OGS_FATAL, and NumLib::Picard.

Referenced by ProjectData::parseNonlinearSolvers().

419 {
421  auto const type = config.getConfigParameter<std::string>("type");
423  auto const max_iter = config.getConfigParameter<int>("max_iter");
425  if (type == "Picard") {
426  auto const tag = NonlinearSolverTag::Picard;
427  using ConcreteNLS = NonlinearSolver<tag>;
428  return std::make_pair(
429  std::make_unique<ConcreteNLS>(linear_solver, max_iter), tag);
430  }
431  if (type == "Newton")
432  {
434  auto const damping = config.getConfigParameter<double>("damping", 1.0);
435  if (damping <= 0)
436  {
438  "The damping factor for the Newon method must be positive, got "
439  "%g.",
440  damping);
441  }
442  auto const tag = NonlinearSolverTag::Newton;
443  using ConcreteNLS = NonlinearSolver<tag>;
444  return std::make_pair(
445  std::make_unique<ConcreteNLS>(linear_solver, max_iter, damping),
446  tag);
447  }
448  OGS_FATAL("Unsupported nonlinear solver type");
449 }
#define OGS_FATAL(fmt,...)
Definition: Error.h:64