Commit a8efe9a8 authored by Timo Koch's avatar Timo Koch
Browse files

[doc] Improve documentation for assembly folder

parent 72135d81
......@@ -393,7 +393,7 @@
/* ***************** Assembly ******************/
/*!
* \defgroup Assembly Assembly
* TODO: Doc me in modules.txt!
* \brief Assembly of linear systems (Jacobian and residual)
*/
/* ***************** Common ******************/
......
......@@ -18,8 +18,9 @@
*****************************************************************************/
/*!
* \file
* \brief An assembler for the global linear system for fully implicit models
* and cell-centered discretization schemes using Newton's method.
* \ingroup Assembly
* \ingroup BoxDiscretization
* \brief An assembler for Jacobian and residual contribution per element (box method)
*/
#ifndef DUMUX_BOX_LOCAL_ASSEMBLER_HH
#define DUMUX_BOX_LOCAL_ASSEMBLER_HH
......@@ -33,16 +34,23 @@
namespace Dumux {
/*!
* \ingroup ImplicitModel
* \brief An assembler for the local contributions (per element) to the global
* linear system for fully implicit models and cell-centered discretization schemes.
* \ingroup Assembly
* \ingroup BoxDiscretization
* \brief An assembler for Jacobian and residual contribution per element (box method)
* \tparam TypeTag the TypeTag
* \tparam DM the differentiation method to residual compute derivatives
* \tparam implicit if to use an implicit or explicit time discretization
*/
template<class TypeTag,
DiffMethod DM = DiffMethod::numeric,
bool implicit = true>
class BoxLocalAssembler;
/*!
* \ingroup Assembly
* \ingroup BoxDiscretization
* \brief Box local assembler using numeric differentiation and implicit time discretization
*/
template<class TypeTag>
class BoxLocalAssembler<TypeTag,
DiffMethod::numeric,
......@@ -442,6 +450,12 @@ private:
}; // implicit BoxAssembler with numeric Jacobian
/*!
* \ingroup Assembly
* \ingroup BoxDiscretization
* \brief Box local assembler using numeric differentiation and explicit time discretization
*/
template<class TypeTag>
class BoxLocalAssembler<TypeTag,
DiffMethod::numeric,
......@@ -809,6 +823,12 @@ private:
}; // explicit BoxAssembler with numeric Jacobian
/*!
* \ingroup Assembly
* \ingroup BoxDiscretization
* \brief Box local assembler using analytic (hard coded) derivatives and implicit time discretization
*/
template<class TypeTag>
class BoxLocalAssembler<TypeTag,
DiffMethod::analytic,
......@@ -1121,6 +1141,12 @@ private:
}; // implicit BoxAssembler with analytic Jacobian
/*!
* \ingroup Assembly
* \ingroup BoxDiscretization
* \brief Box local assembler using analytic (hard coded) derivatives and explicit time discretization
*/
template<class TypeTag>
class BoxLocalAssembler<TypeTag,
DiffMethod::analytic,
......
......@@ -18,7 +18,9 @@
*****************************************************************************/
/*!
* \file
* \brief Calculates the residual of models based on the box scheme element-wise.
* \ingroup Assembly
* \ingroup BoxDiscretization
* \brief Calculates the element-wise residual for the box scheme
*/
#ifndef DUMUX_BOX_LOCAL_RESIDUAL_HH
#define DUMUX_BOX_LOCAL_RESIDUAL_HH
......@@ -29,14 +31,13 @@
#include <dumux/common/properties.hh>
#include <dumux/assembly/fvlocalresidual.hh>
namespace Dumux
{
namespace Dumux {
/*!
* \ingroup BoxModel
* \brief Element-wise calculation of the residual for models
* based on the fully implicit box scheme.
*
* \todo Please doc me more!
* \ingroup Assembly
* \ingroup BoxDiscretization
* \brief The element-wise residual for the box scheme
* \tparam TypeTag the TypeTag
*/
template<class TypeTag>
class BoxLocalResidual : public FVLocalResidual<TypeTag>
......@@ -67,6 +68,7 @@ class BoxLocalResidual : public FVLocalResidual<TypeTag>
public:
using ParentType::ParentType;
//! evaluate flux residuals for one sub control volume face and add to residual
void evalFlux(ElementResidualVector& residual,
const Problem& problem,
const Element& element,
......@@ -91,7 +93,7 @@ public:
}
}
//! evaluate flux residuals for one sub control volume face
ResidualVector evalFlux(const Problem& problem,
const Element& element,
const FVElementGeometry& fvGeometry,
......
......@@ -18,8 +18,12 @@
*****************************************************************************/
/*!
* \file
* \brief An assembler for the global linear system for fully implicit models
* and cell-centered discretization schemes using Newton's method.
* \ingroup Assembly
* \ingroup CCDiscretization
* \brief An assembler for Jacobian and residual contribution per element (cell-centered methods)
* \tparam TypeTag the TypeTag
* \tparam DM the differentiation method to residual compute derivatives
* \tparam implicit if to use an implicit or explicit time discretization
*/
#ifndef DUMUX_CC_LOCAL_ASSEMBLER_HH
#define DUMUX_CC_LOCAL_ASSEMBLER_HH
......@@ -34,16 +38,23 @@
namespace Dumux {
/*!
* \ingroup ImplicitModel
* \brief An assembler for the local contributions (per element) to the global
* linear system for fully implicit models and cell-centered discretization schemes.
* \ingroup Assembly
* \ingroup CCDiscretization
* \brief An assembler for Jacobian and residual contribution per element (cell-centered methods)
* \tparam TypeTag the TypeTag
* \tparam DM the differentiation method to residual compute derivatives
* \tparam implicit if to use an implicit or explicit time discretization
*/
template<class TypeTag,
DiffMethod DM = DiffMethod::numeric,
bool implicit = true>
class CCLocalAssembler;
/*!
* \ingroup Assembly
* \ingroup CCDiscretization
* \brief Cell-centered scheme local assembler using numeric differentiation and implicit time discretization
*/
template<class TypeTag>
class CCLocalAssembler<TypeTag,
DiffMethod::numeric,
......@@ -586,7 +597,12 @@ private:
{ return gridVolVars.volVars(scv); }
};
//! Explicit assembler with numeric differentiation
/*!
* \ingroup Assembly
* \ingroup CCDiscretization
* \brief Cell-centered scheme local assembler using numeric differentiation and explicits time discretization
*/
template<class TypeTag>
class CCLocalAssembler<TypeTag,
DiffMethod::numeric,
......@@ -926,7 +942,11 @@ private:
{ return gridVolVars.volVars(scv); }
};
//! implicit assembler using analytic differentiation
/*!
* \ingroup Assembly
* \ingroup CCDiscretization
* \brief Cell-centered scheme local assembler using analytic (hand-coded) differentiation and implicit time discretization
*/
template<class TypeTag>
class CCLocalAssembler<TypeTag,
DiffMethod::analytic,
......@@ -1187,7 +1207,11 @@ private:
}
};
//! explicit assembler using analytic differentiation
/*!
* \ingroup Assembly
* \ingroup CCDiscretization
* \brief Cell-centered scheme local assembler using analytic (hand-coded) differentiation and explicit time discretization
*/
template<class TypeTag>
class CCLocalAssembler<TypeTag,
DiffMethod::analytic,
......
......@@ -18,7 +18,9 @@
*****************************************************************************/
/*!
* \file
* \brief Calculates the residual of models based on the box scheme element-wise.
* \ingroup Assembly
* \ingroup CCDiscretization
* \brief Calculates the element-wise residual for cell-centered discretization schemes
*/
#ifndef DUMUX_CC_LOCAL_RESIDUAL_HH
#define DUMUX_CC_LOCAL_RESIDUAL_HH
......@@ -28,14 +30,12 @@
#include <dumux/common/properties.hh>
#include <dumux/assembly/fvlocalresidual.hh>
namespace Dumux
{
namespace Dumux {
/*!
* \ingroup CCModel
* \brief Element-wise calculation of the residual for models
* based on the fully implicit cell-centered scheme.
*
* \todo Please doc me more!
* \ingroup Assembly
* \ingroup CCDiscretization
* \brief Calculates the element-wise residual for the cell-centered discretization schemes
*/
template<class TypeTag>
class CCLocalResidual : public FVLocalResidual<TypeTag>
......@@ -54,6 +54,7 @@ class CCLocalResidual : public FVLocalResidual<TypeTag>
public:
using ParentType::ParentType;
//! evaluate the flux residual for a sub control volume face and add to residual
void evalFlux(ElementResidualVector& residual,
const Problem& problem,
const Element& element,
......@@ -68,6 +69,7 @@ public:
residual[localScvIdx] += evalFlux(problem, element, fvGeometry, elemVolVars, elemFluxVarsCache, scvf);
}
//! evaluate the flux residual for a sub control volume face
ResidualVector evalFlux(const Problem& problem,
const Element& element,
const FVElementGeometry& fvGeometry,
......
......@@ -18,15 +18,21 @@
*****************************************************************************/
/*!
* \file
* \brief An enum class to define various methods available in order to compute
* \ingroup Assembly
* \brief An enum class to define various differentiation methods available in order to compute
the derivatives of the residual i.e. the entries in the jacobian matrix.
*/
#ifndef DUMUX_JACOBIAN_DIFFERENTIATION_METHODS_HH
#define DUMUX_JACOBIAN_DIFFERENTIATION_METHODS_HH
namespace Dumux {
/*!
* \ingroup Assembly
* \brief Differentiation methods in order to compute the derivatives
* of the residual i.e. the entries in the jacobian matrix.
* \todo automatic differentation is not yet implemented
*/
enum class DiffMethod
{
numeric, analytic, automatic
......
......@@ -18,8 +18,8 @@
*****************************************************************************/
/*!
* \file
* \brief An assembler for the global linear system
* for fully implicit models and vertex-centered discretization schemes.
* \ingroup Assembly
* \brief A linear system assembler (residual and Jacobian) for finite volume schemes
*/
#ifndef DUMUX_FV_ASSEMBLER_HH
#define DUMUX_FV_ASSEMBLER_HH
......@@ -40,9 +40,11 @@
namespace Dumux {
/*!
* \ingroup ImplicitModel
* \brief An assembler for the global linear system
* for fully implicit models and cell-centered discretization schemes.
* \ingroup Assembly
* \brief A linear system assembler (residual and Jacobian) for finite volume schemes
* \tparam TypeTag the TypeTag
* \tparam diffMethod the differentiation method to residual compute derivatives
* \tparam isImplicit if to use an implicit or explicit time discretization
*/
template<class TypeTag, DiffMethod diffMethod, bool isImplicit = true>
class FVAssembler
......
......@@ -18,7 +18,8 @@
*****************************************************************************/
/*!
* \file
* \brief Calculates the element-wise residual of finite-volume models.
* \ingroup Assembly
* \brief The element-wise residual for finite volume schemes
*/
#ifndef DUMUX_FV_LOCAL_RESIDUAL_HH
#define DUMUX_FV_LOCAL_RESIDUAL_HH
......@@ -29,15 +30,13 @@
#include <dumux/common/timeloop.hh>
#include <dumux/discretization/methods.hh>
namespace Dumux
{
namespace Dumux {
/*!
* \brief Element-wise calculation of the residual matrix for models
* using a fully implicit discretization.
*
* \ingroup Assembly
* \brief The element-wise residual for finite volume schemes
* \note This class defines the interface used by the assembler using
* static polymorphism.
* static polymorphism. Implementations are specialized for a certain discretization scheme
*/
template<class TypeTag>
class FVLocalResidual
......@@ -85,8 +84,12 @@ public:
* This can be used to figure out how much of each conservation
* quantity is inside the element.
*
* \param problem The problem to solve
* \param element The DUNE Codim<0> entity for which the storage
* term ought to be calculated
* \param fvGridGeometry The finite-volume grid geometry
* \param gridVariables The grid variables (volume and flux variables)
* \param sol The solution vector
*/
ElementResidualVector evalStorage(const Problem& problem,
const Element &element,
......@@ -127,19 +130,17 @@ public:
/*!
* \brief Compute the local residual, i.e. the deviation of the
* equations from zero for instationary problems.
*
* \param problem The problem to solve
* \param element The DUNE Codim<0> entity for which the residual
* ought to be calculated
* \param fvGeometry The finite-volume geometry of the element
* \param prevVolVars The volume averaged variables for all
* sub-control volumes of the element at the previous
* time level
* \param curVolVars The volume averaged variables for all
* sub-control volumes of the element at the current
* time level
* \param bcTypes The types of the boundary conditions for all
* vertices of the element
* \param prevElemVolVars The volume averaged variables for all
* sub-control volumes of the element at the previous time level
* \param curElemVolVars The volume averaged variables for all
* sub-control volumes of the element at the current time level
* \param bcTypes The types of the boundary conditions for all boundary entities of an element
* \param elemFluxVarsCache The flux variable caches for the element stencil
*/
ElementResidualVector eval(const Problem& problem,
const Element& element,
......@@ -176,19 +177,17 @@ public:
/*!
* \brief Compute the storage local residual, i.e. the deviation of the
* storage term from zero for instationary problems.
*
* \param problem The problem to solve
* \param element The DUNE Codim<0> entity for which the residual
* ought to be calculated
* \param fvGeometry The finite-volume geometry of the element
* \param prevVolVars The volume averaged variables for all
* sub-control volumes of the element at the previous
* time level
* \param curVolVars The volume averaged variables for all
* sub-control volumes of the element at the current
* time level
* \param bcTypes The types of the boundary conditions for all
* vertices of the element
* \param prevElemVolVars The volume averaged variables for all
* sub-control volumes of the element at the previous time level
* \param curElemVolVars The volume averaged variables for all
* sub-control volumes of the element at the current time level
* \param bcTypes The types of the boundary conditions for all boundary entities of an element
* \param elemFluxVarsCache The flux variable caches for the element stencil
*/
ElementResidualVector evalStorage(const Problem& problem,
const Element& element,
......@@ -218,16 +217,15 @@ public:
/*!
* \brief Compute the local residual, i.e. the deviation of the
* equations from zero for stationary problem.
*
* \param problem The problem to solve
* \param element The DUNE Codim<0> entity for which the residual
* ought to be calculated
* \param fvGeometry The finite-volume geometry of the element
* \param curVolVars The volume averaged variables for all
* sub-control volumes of the element at the current
* time level
* \param bcTypes The types of the boundary conditions for all
* vertices of the element
* \param curElemVolVars The volume averaged variables for all
* sub-control volumes of the element at the current time level
* \param bcTypes The types of the boundary conditions for all boundary entities of an element
* \param elemFluxVarsCache The flux variable caches for the element stencil
*/
ElementResidualVector eval(const Problem& problem,
const Element& element,
......@@ -261,14 +259,16 @@ public:
/*!
* \name Model specific interface
* \note The following method are the model specific implementations of the actual residual
* \note The following method are the model specific implementations of the local residual
*/
// \{
/*!
* \brief Calculate the source term of the equation
*
* \param scv The sub-control volume over which we integrate the source term
* \param problem The problem to solve
* \param scv The sub-control volume over which we integrate the storage term
* \param volVars The volume variables associated with the scv
* \note has to be implemented by the model specific residual class
*
*/
......@@ -282,6 +282,11 @@ public:
/*!
* \brief Calculate the source term of the equation
*
* \param problem The problem to solve
* \param element The DUNE Codim<0> entity for which the residual
* ought to be calculated
* \param fvGeometry The finite-volume geometry of the element
* \param elemVolVars The volume variables associated with the element stencil
* \param scv The sub-control volume over which we integrate the source term
* \note This is the default implementation for all models as sources are computed
* in the user interface of the problem
......@@ -307,7 +312,14 @@ public:
/*!
* \brief Calculate the source term of the equation
*
* \param scv The sub-control volume over which we integrate the source term
* \param problem The problem to solve
* \param element The DUNE Codim<0> entity for which the residual
* ought to be calculated
* \param fvGeometry The finite-volume geometry of the element
* \param elemVolVars The volume variables associated with the element stencil
* \param scvf The sub-control volume over which we integrate the flux
* \param elemFluxVarsCache the flux variable caches for the element's flux stencils
*
* \note has to be implemented by the model specific residual class
*
*/
......@@ -329,6 +341,21 @@ public:
*/
// \{
/*!
* \brief Compute the storage local residual, i.e. the deviation of the
* storage term from zero for instationary problems.
*
* \param residual The residual vector to fill
* \param problem The problem to solve
* \param element The DUNE Codim<0> entity for which the residual
* ought to be calculated
* \param fvGeometry The finite-volume geometry of the element
* \param prevElemVolVars The volume averaged variables for all
* sub-control volumes of the element at the previous time level
* \param curElemVolVars The volume averaged variables for all
* sub-control volumes of the element at the current time level
* \param scv The sub control volume the storage term is integrated over
*/
void evalStorage(ElementResidualVector& residual,
const Problem& problem,
const Element& element,
......@@ -344,7 +371,7 @@ public:
// \f$\frac{m}{\partial t}\f$ term if using implicit or explicit
// euler as time discretization.
//
// We might need a more explicit way for
// TODO: We might need a more explicit way for
// doing the time discretization...
//! Compute storage with the model specific storage residual
......@@ -361,6 +388,19 @@ public:
residual[scv.indexInElement()] += storage;
}
/*!
* \brief Compute the source local residual, i.e. the deviation of the
* source term from zero.
*
* \param residual The residual vector to fill
* \param problem The problem to solve
* \param element The DUNE Codim<0> entity for which the residual
* ought to be calculated
* \param fvGeometry The finite-volume geometry of the element
* \param curElemVolVars The volume averaged variables for all
* sub-control volumes of the element at the current time level
* \param scv The sub control volume the source term is integrated over
*/
void evalSource(ElementResidualVector& residual,
const Problem& problem,
const Element& element,
......@@ -377,6 +417,22 @@ public:
residual[scv.indexInElement()] -= source;
}
/*!
* \brief Compute the source local residual, i.e. the deviation of the
* source term from zero and add to the residual
* \note This is implemented for different discretization schemes
*
* \param residual The residual vector to fill
* \param problem The problem to solve
* \param element The DUNE Codim<0> entity for which the residual
* ought to be calculated
* \param fvGeometry The finite-volume geometry of the element
* \param elemVolVars The volume averaged variables for all
* sub-control volumes of the element at the current time level
* \param elemBcTypes the boundary types for the boundary entities of an elements
* \param elemFluxVarsCache The flux variable caches for the element stencil
* \param scvf The sub control volume face the flux term is integrated over
*/
void evalFlux(ElementResidualVector& residual,
const Problem& problem,
const Element& element,
......@@ -386,6 +442,19 @@ public:
const ElementFluxVariablesCache& elemFluxVarsCache,
const SubControlVolumeFace& scvf) const {}
/*!
* \brief Compute the flux local residual, i.e. the deviation of the
* flux term from zero.
*
* \param problem The problem to solve
* \param element The DUNE Codim<0> entity for which the residual
* ought to be calculated
* \param fvGeometry The finite-volume geometry of the element
* \param elemVolVars The volume averaged variables for all
* sub-control volumes of the element at the current time level
* \param elemFluxVarsCache The flux variable caches for the element stencil
* \param scvf The sub control volume face the flux term is integrated over
*/
ResidualVector evalFlux(const Problem& problem,
const Element& element,
const FVElementGeometry& fvGeometry,
......@@ -396,6 +465,14 @@ public:
return asImp().evalFlux(problem, element, fvGeometry, elemVolVars, elemFluxVarsCache, scvf);
}
//\}
/*!
* \name Interfaces for analytic Jacobian computation
*/
// \{
//! Compute the derivative of the storage residual
template<class PartialDerivativeMatrix>
void addStorageDerivatives(PartialDerivativeMatrix& partialDerivatives,
const Problem& problem,
......@@ -407,6 +484,7 @@ public:
DUNE_THROW(Dune::NotImplemented, "analytic storage derivative");
}
//! Compute the derivative of the source residual
template<class PartialDerivativeMatrix>
void addSourceDerivatives(PartialDerivativeMatrix& partialDerivatives,
const Problem& problem,
......@@ -418,6 +496,7 @@ public:
DUNE_THROW(Dune::NotImplemented, "analytic source derivative");
}
//! Compute the derivative of the flux residual
template<class PartialDerivativeMatrices, class T = TypeTag>
std::enable_if_t<GET_PROP_VALUE(T, DiscretizationMethod) != DiscretizationMethods::Box, void>
addFluxDerivatives(PartialDerivativeMatrices& derivativeMatrices,
......@@ -431,6 +510,7 @@ public:
DUNE_THROW(Dune::NotImplemented, "analytic flux derivative for cell-centered models");
}
//! Compute the derivative of the flux residual for the box method
template<class JacobianMatrix, class T = TypeTag>
std::enable_if_t<GET_PROP_VALUE(T, DiscretizationMethod) == DiscretizationMethods::Box, void>
addFluxDerivatives(JacobianMatrix& A,
......@@ -444,6 +524,7 @@ public:
DUNE_THROW(Dune::NotImplemented, "analytic flux derivative for box models");
}
//! Compute the derivative of the Dirichlet flux residual for cell-centered schemes
template<class PartialDerivativeMatrices>
void addCCDirichletFluxDerivatives(PartialDerivativeMatrices& derivativeMatrices,
const Problem& problem,
......@@ -456,6 +537,7 @@ public:
DUNE_THROW(Dune::NotImplemented, "analytic Dirichlet flux derivative");
}
//! Compute the derivative of Robin type boundary conditions ("solution dependent Neumann")
template<class PartialDerivativeMatrices>
void addRobinFluxDerivatives(PartialDerivativeMatrices& derivativeMatrices,
const Problem& problem,
......@@ -468,15 +550,17 @@ public:
DUNE_THROW(Dune::NotImplemented, "analytic Robin flux derivative");
}
//\}
/*!
* \brief Sets the solution from which to start the time integration. Has to be
* called prior to assembly for time-dependent problems.
* \brief Sets the solution from which to start the time integration.
* \note Has to be called prior to assembly for time-dependent problems.
*/
void setPreviousSolution(const SolutionVector& u)
{ prevSol_ = &u; }
/*!
* \brief Return the solution that has been set as the previous one.
* \brief Return the solution of the previous time level that has been set
*/