Problem

Summary

The ProblemStat and ProblemInstat are classes collecting everything necessary to assemble stationary or instationary PDEs, respectively.

A ProblemStat contains the data structures for the linear system, i.e. system matrix, solution and right-hand side vector, it contains the grid and global basis, a linear solver, marker and error estimator, and much more, see below.

A ProblemInstat is responsible for storing data needed for time discretization, like an old-solution vector, references to time, timestep width and the stationary problem to solve in each timestep.

These two classes build the basis for most AMDiS projects.

Problem classes

Class Descriptions
ProblemStat Standard container for stationary problems (class template)
ProblemInstat Standard container for instationary problems (class template)

Interfaces

Interface Descriptions
StandardProblemIteration A master problem for a single non coupled problem. (class)
ProblemStatBase Base class for stationary problems. (class)
ProblemInstatBase Base class for instationary problems. (class)
ProblemTimeInterface Interface for time dependent problems. (abstract class)
ProblemIterationInterface Interface for iterations of stationary problems (abstract class)

class ProblemStat

Defined in header <amdis/ProblemStat.hpp>

template <class Traits>
class ProblemStat
  : public ProblemStatBase
  , public StandardProblemIteration

ProblemStat is a structure storing data structures for stationary problems.

The template parameter Traits defines the global basis and with this the GridView and Grid. It is required to have at least the member types GlobalBasis and CoefficientType and may provide additionally a static create(GridView) function to contruct the global basis from a GridView object. See DefaultProblemTraits and DefaultBasisCreate for a detailed reference.

Member Types

Member Type Definition
GlobalBasis typename Traits::GlobalBasis
GridView typename GlobalBasis::GridView
Grid typename GridView::Grid
WorldMatrix FieldMatrix<K, dow, dow>
WorldVector typename Element::Geometry::GlobalCoordinate
SystemMatrix DOFMatrix<GlobalBasis, GlobalBasis, C>
SystemVector DOFVector<GlobalBasis, C>

Here, K is the coordinate type, i.e. K = typename Grid::ctype, the constant dow is the world dimension, i.e. dow = Grid::dimensionworld, and the type C is the coefficients type, given in the Traits as C = typename Traits::CoefficientType.

Member functions

Function Descriptions
(constructor) Construct the problem
initialize Initializes all data-structures in the problem
addMatrixOperator Adds an operator to the system matrix
addVectorOperator Adds an operator to the rhs vector
addDirichletBC Adds a Dirichlet boundary condition to the system
addPeriodicBC Adds a periodic boundary condition to the system
writeFiles Writes output files.

Getters and Setters

Function Descriptions
name The name of the problem
grid, setGrid Get and set the stored grid
gridView Get the GridView stored in the global basis
boundaryManager Return the boundary manager to identify boundary segments
globalBasis Return the stored global basis
solver, setSolver Get and set the linear solver
systemMatrix, solutionVector, rhsVector Return the matrix and vector of the linear system
solution Return a grid function of the solution

Implementation of the ProblemStatBase interface

Function Descriptions
buildAfterAdapt, assemble Assembles the linear system after grid operations.
markElements Marks mesh elements for refinement and coarsening.
adaptGrid Refinement/coarsening of the grid.
globalCoarsen Uniform global grid coarsening.
globalRefine Uniform global refinement.
solve Solves the assembled system.
estimate A-posteriori error estimation.

Implementation of the ProblemIterationInterface interface

Function Descriptions
oneIteration A single build-solve-adapt step

function ProblemStat::ProblemStat

explicit ProblemStat(std::string const& name)                               // (1)

ProblemStat(std::string const& name, Grid& grid)                            // (2)

ProblemStat(std::string const& name, Grid& grid, GlobalBasis& globalBasis)  // (3)

(1) Construct the problem from the grid type and basis type given in the Traits template parameter. Thereby it is assumed, that the Traits type provides a static create(GridView) function to create the basis. The grid is created by inspecting the parameter file and using the MeshCreator class.

(2) Construct the problem with given grid. Store the grid reference into a non-destroying shared_ptr and create the basis from the Traits class, see (1).

(3) Construct the problem with given grid and basis, by storing both references into non-destroying shared_ptr.

Arguments

std::string name
The name of the problem that is used to identify parameters in the initfile.
Grid grid
A grid implementing the Dune::Grid interface.
GlobalBasis globalBasis
A basis implementing the Dune::Functions::GlobalBasis interface.

Example

using Grid = Dune::YaspGrid<2>;
Grid grid({1.0, 1.0}, {2, 2});

using namespace Dune::Functions::BasisFactory;
auto basis = makeBasis(grid.leafGridView(), lagrange<2>());

// use predefined traits type `LagrangeBasis`
ProblemStat<LagrangeBasis<Grid, 2>> prob1("prob");

// use predefined traits type `LagrangeBasis` but provide a grid directly
// The grid type must match the one defined in the Traits class.
ProblemStat<LagrangeBasis<Grid, 2>> prob2("prob", grid);

// use predefined traits type `LagrangeBasis` but provide a grid and basis directly
// The grid and basis type must match those defined in the Traits class.
ProblemStat<LagrangeBasis<Grid, 2>> prob3("prob", grid, basis);

// using c++17 class template argument deduction
ProblemStat prob4("prob", grid, basis);

function ProblemStat::initialize

void initialize(Flag initFlag, ProblemStat* adoptProblem = nullptr, Flag adoptFlag = INIT_NOTHING)

Initialization of the problem and its data members.

Constructs the grid and basis (if not yet provided otherwise) and creates a linear system data-structures, i.e. system matrix, solution vector and right-hand side vector, based on the size of the basis.

Arguments

Flag initFlat
A flag indicating what to initialize, use INIT_ALL to initialize everything. See below for possible values.
ProblemStat* adoptProblem
A problem from which to adopt data not initialized in this problem.
Flag adoptFlag
A flag indicating what to adopt from the adoptProblem.

Init Flags

Flag Descriptions
INIT_FE_SPACE Initialize the global basis
INIT_MESH Initialize the grid
CREATE_MESH Create the grid
INIT_SYSTEM Initialize system matrix, solution vector and right-hand side vector
INIT_SOLVER Initialize a linear solver
INIT_ESTIMATOR Initialize an error estimator
INIT_MARKER Initialize a grid marker
INIT_FILEWRITER Initialize a file writer
INIT_NOTHING Do not initialize anything
INIT_ALL Initialize everything

function ProblemStat::addMatrixOperator

// (1)
template <class Operator, class RowTreePath = RootTreePath, class ColTreePath = RootTreePath>
void addMatrixOperator(Operator const& op, RowTreePath row = {}, ColTreePath col = {})

// (2)
template <class Operator, class RowTreePath = RootTreePath, class ColTreePath = RootTreePath>
void addMatrixOperator(BoundaryType b, Operator const& op, RowTreePath row = {}, ColTreePath col = {})

(1) Add a local operator, evaluated on the elements, to the system matrix.

(2) Add a local operator, evaluated on boundary intersections, to the system matrix.

Arguments

Operator op
A (pre-) local operator to assemble on the local context, see LocalOperator and GridFunctionOperator.
RowTreePath row
TreePath identifying the sub-basis in the global basis tree corresponding to the row basis, see makeTreePath()
ColTreePath col
TreePath identifying the sub-basis in the global basis tree corresponding to the column basis, see makeTreePath()
BoundaryType b
Boundary indentifier/number to set on which part of the boundary to assemble this operator. Can be constructed from an integer, see BoundaryType.

Note

If no row or col tree-path is given, the root tree-path is assumed, identifying the root node in the basis tree.

The tree-paths can be constructed from integers or integral-constants, or by using the makeTreePath() function. It is then transformed into a Dune::TypeTree::HybridTreePath<...> type.

Examples

auto op1 = makeOperator(tag::gradtest_gradtrial{}, alpha);
prob.addMatrixOperator(op1, 0, 0);

auto op2 = makeOperator(tag::test_trial{}, beta);
prob.addMatrixOperator(BoundaryType{1}, op2, 0, 0);

See Also

function ProblemStat::addVectorOperator

// (1)
template <class Operator, class TreePath = RootTreePath>
void addVectorOperator(Operator const& op, TreePath path = {})

// (2)
template <class Operator, class TreePath = RootTreePath>
void addVectorOperator(BoundaryType b, Operator const& op, TreePath path = {})

(1) Add a local operator, evaluated on the elements, to the right-hand side vector of the linear system.

(2) Add a local operator, evaluated on boundary intersections, to the right-hand side vector of the linear system

Arguments

Operator op
A (pre-) local operator to assemble on the local context, see LocalOperator and GridFunctionOperator.
TreePath path
TreePath identifying the sub-basis in the global basis tree corresponding to the row basis, see makeTreePath()
BoundaryType b
Boundary indentifier/number to set on which part of the boundary to assemble this operator. Can be constructed from an integer, see BoundaryType.

Note

If no path tree-path is given, the root tree-path is assumed, identifying the root node in the basis tree.

The tree-path can be constructed from integers or integral-constants, or by using the makeTreePath() function. It is then transformed into a Dune::TypeTree::HybridTreePath<...> type.

Examples

auto op = makeOperator(tag::test{}, [g](auto const& x) { return g(x); });
prob.addVectorOperator(BoundaryType{1}, op, _0);

See Also

function ProblemStat::addDirichletBC

// (1)
template <class Predicate, class RowTreePath, class ColTreePath, class Values>
void addDirichletBC(Predicate const& predicate,
                    RowTreePath row, ColTreePath col,
                    Values const& values)

// (2)
template <class RowTreePath, class ColTreePath, class Values>
void addDirichletBC(BoundaryType id,
                    RowTreePath row, ColTreePath col,
                    Values const& values)

Add a Dirichlet boundary condition to the linear system.

(1) The boundary part is identified using a boolean predicate function that is evaluated on the barycenters of the boundary intersections

(2) The boundary part is identified using a boundary identifier/number that was set in the BoundaryManager before.

Arguments

Predicate predicate
A functor bool(WorldVector) returning true, if the boundary condition should be assembled for the DOFs on the boundary intersection.
RowTreePath row
TreePath identifying the sub-basis in the global basis tree corresponding to the row basis, see makeTreePath()
ColTreePath col
TreePath identifying the sub-basis in the global basis tree corresponding to the column basis, see makeTreePath()
Values values
A GridFunction evaluating to the Dirichlet values that should be set for the identified DOFs.
BoundaryType id
Boundary indentifier/number to set on which part of the boundary to assemble this boundary condition. Can be constructed from an integer, see BoundaryType.

Examples

// add a Dirichlet BC at the left boundary of the domain [0,1]^2
prob.addDirichletBC([](auto const& x) { return x[0] < 1.e-8; },
                    0, 0,
                    [](auto const& x) { return 0.0; });

// same as above but using the boundary id
prob.boundaryManager()->setBoxBoundary({1,2,2,2});
prob.addDirichletBC(BoundaryType{1}, 0, 0, 0.0);

class ProblemInstat

Defined in header <amdis/ProblemInstat.hpp>

template <class Traits>
class ProblemInstat
  : public ProblemInstatBase

Standard implementation of ProblemTimeInterface for a time dependent problems.

Member Functions

Function Descriptions
(constructor) Constructs a ProblemInstat from a stationary problem.
initialize Initialization of the instationary problem.
initTimestep Implementation of ProblemTimeInterface::initTimestep().
closeTimestep Implementation of ProblemTimeInterface::closeTimestep().
problemStat Returns the stored stationary problem.
oldSolutionVector Returns the data vector for the old-solution
oldSolution Returns the old-solution as discrete function
transferInitialSolution Implementation of ProblemTimeInterface::transferInitialSolution().

class StandardProblemIteration

Defined in header <amdis/StandardProblemIteration.hpp>

class StandardProblemIteration
  : public virtual ProblemIterationInterface

A master problem for a single non-coupled problem.

Member Functions

Function Descriptions
(constructor) Constructs a StandardProblemIteration
beginIteration Implementation of ProblemIterationIterface::beginIteration()
oneIteration Implementation of ProblemIterationInterface::oneIteration()
endIteration Implementation of ProblemIterationInterface::endIteration()
name Returns the name of the iteration interface.
numProblems Returns number of managed problems.
problem Returns a managed ProblemStat.

class ProblemStatBase

Defined in header <amdis/ProblemStatBase.hpp>

class ProblemStatBase

Interface for stationary problems. Concrete problems must override all pure virtual methods. Base class for ProblemStat.

Member Functions

Function Descriptions
markElements Marks mesh elements for refinement and coarsening.
buildAfterAdapt Assembling of system matrices and vectors after adaption.
adaptGrid Refinement/coarsening of the grid.
globalCoarsen Uniform global grid coarsening.
globalRefine Uniform global grid refinement.
solve Solves the assembled system.
estimate A-posteriori error estimation.
name Returns the name of the problem.

class ProblemInstatBase

Defined in header <amdis/ProblemInstatBase.hpp>

class ProblemInstatBase
  : public virtual ProblemTimeInterface

Base class for ProblemInstat.

Member Functions

Function Descriptions
(constructor) Constructs a new ProblemStatBase
solveInitialProblem Implementation of ProblemTimeInterface::solveInitialProblem().
name Returns the name of the instationary problem.
time Returns reference to current simulation time set in setTime() from AdaptInfo::time().
setTime Implementation of ProblemTimeInterface::setTime().
tau Returns reference to current simulation timestep set in setTime() from AdaptInfo::timestep().
invTau Returns reference to current simulation 1.0/timestep.

class ProblemIterationInterface

Defined in header <amdis/ProblemIterationInterface.hpp>

class ProblemIterationInterface

Interface for master problems needed by the adaption loop. A master problem can handle one single or multiple coupled problems. In the latter case, the master problem can determine the execution order of the build, solve, estimate, and adapt steps of the single problems in oneIteration().

Member Functions

Function Descriptions
beginIteration Called before each adaption loop iteration.
oneIteration Determines the execution order of the single adaption steps.
endIteration Called after each adaption loop iteration.
numProblems Returns number of managed problems.
problem Returns the problem with the given number.
name Returns the name of the iteration interface.

class ProblemTimeInterface

Defined in header <amdis/ProblemTimeInterface.hpp>

class ProblemTimeInterface

Interface for time dependent problems.

Member Functions

Function Descriptions
initTimeInterface Called at the beginning of the adaption loop before any other call.
setTime Executes all needed operations when the simulation time changes.
initTimestep Called at the beginning of each timestep.
closeTimestep Called at the end of each timestep.
solveInitialProblem Solves the initial problem.
transferInitialSolution Transfer the initial problem.