version 3.9-dev
Multi-point flux approximation (Mpfa)

A cell-centered finite volume scheme with multi-point flux approximation. More...

Description

We continue from the introduction to Cell-centered Finite Volume Methods using as an example the following scalar elliptic problem for the unknown \(u\),

\begin{equation} \begin{aligned} \nabla \cdot \left( - \boldsymbol{\Lambda} \nabla u \right) &= q &&\mathrm{in} \, \Omega \\ \left( - \boldsymbol{\Lambda} \nabla u \right) \cdot \boldsymbol{n} &= v_N &&\mathrm{on} \, \Gamma_N \\ u &= u_D &&\mathrm{on} \, \Gamma_D, \label{eq:elliptic} \end{aligned} \end{equation}

and recall that after integration of the governing equation over a control volume \(K\) and application of the divergence theorem, we replace exact fluxes by numerical approximations \(F_{K, \sigma} \approx \int_{\sigma} \left( - \boldsymbol{\Lambda}_K \nabla u \right) \cdot \boldsymbol{n} \mathrm{d} \Gamma\).

For the Mpfa method, expressions for the face fluxes \(F_{K, \sigma}\) are obtained by introducing intermediate face unknowns \(u_\sigma\) in addition to the cell unknowns \(u_K\) and enforcing the physically motivated continuity of fluxes and continuity of the solution across the faces. For a face \(\sigma\) between the two polygons \(K\) and \(L\) these conditions read:

\begin{equation} \begin{aligned} &F_{K, \sigma} + F_{L, \sigma} = 0 \\ &{u}_{K,\sigma} = {u}_{L,\sigma} = {u}_{\sigma}. \label{eq:sigmaConditions} \end{aligned} \end{equation}

Using these conditions, the intermediate face unknowns \({u}_\sigma\) can be eliminated and the fluxes are expressed as a function of the cell unknowns \(u_N\) and associated transmissibilities \(t^N_{K,\sigma}\):

\begin{equation} F_{K,\sigma} = \sum_{N \in \mathcal{S}_{K,\sigma}} t^N_{K,\sigma} u_{N}. \label{eq:FVFluxExpression} \end{equation}

Interaction region for the Mpfa-O method. Three elements belong to the interaction region around vertex v. The fluxes over the sub-control volume faces therefore depend on the three cell unknowns uK, uL, uM.

The main difference between different finite volume schemes is the choice of the flux expression, i.e. the computation of the \(t^N_{K,\sigma}\) and the size of \(\mathcal{S}_{K,\sigma}\). For the Tpfa method this is presented in Two-point flux approximation (Tpfa). There, the stencil and transmissibilities are given as

\begin{equation*} \mathcal{S}_{K,\sigma} = \lbrace K,L \rbrace, \quad t^K_{K,\sigma} = \vert{\sigma}\vert \frac{t_{K,\sigma} t_{L,\sigma}}{t_{K,\sigma} + t_{L,\sigma}},\; t^L_{K,\sigma} = -\vert{\sigma}\vert \frac{t_{K,\sigma} t_{L,\sigma}}{t_{K,\sigma} + t_{L,\sigma}}, \end{equation*}

with the transmissibilities \(t_{K,\sigma}, t_{L,\sigma}\) (see Two-point flux approximation (Tpfa)).

In the following, we present a multi-point flux approximation method called Mpfa-O method introduced in [1]. The main difference to the Tpfa scheme is the fact that a consistent discrete gradient is constructed, i.e. the term \(\nabla u \cdot \mathbf{d}^{\bot}_{K,\sigma}\) is not neglected.

For this scheme, a dual grid is created by connecting the barycenters of the cells with the barycenters of the faces ( \(d=2\)) or the barycenters of the faces and edges ( \(d=3\)). This divides each primary grid face into \(n\) sub-control volume faces \(\sigma^v\), where \(n\) is the number of corners of the primary grid face and the superscript \(v\) refers to the vertex the sub-face can be associated with (see figure above). Also, it sub-divides the control volume \(K\) into regions \(K_v\) that can be associated with the vertex \(v\). These regions are used to construct a local interpolation scheme. We allow for piecewise constant \(\mathbf{\Lambda}\) (denoted as \(\mathbf{\Lambda}_K\) for each cell \(K\)) and construct discrete gradients \(\nabla_\mathcal{D}^{K_v} u\) in each region \(K_v\). In the following, we restrict our discussion to the two-dimensional setup that is shown in the figure above. Here, the discrete gradients are constructed to be consistent such that the following conditions hold:

\begin{equation} \nabla_\mathcal{D}^{K_v} u \cdot (\mathbf{x}_{\sigma^v_1}- \mathbf{x}_{K}) = u_{\sigma^v_1} - u_K, \quad \nabla_\mathcal{D}^{K_v} u \cdot (\mathbf{x}_{\sigma^v_3}- \mathbf{x}_{K}) = u_{\sigma^v_3} - u_K. \end{equation}

Thus, a discrete gradient in \(K_v\) that fulfills these conditions is given as

\begin{equation} \nabla_\mathcal{D}^{K_v} u = \mathbb{D}^{-T}_{K_v} \begin{bmatrix} u_{\sigma^v_1} - u_K \\ u_{\sigma^v_3} - u_K \end{bmatrix}, \qquad \text{ with }\; \mathbb{D}_{K_v} := \begin{bmatrix} \mathbf{x}_{\sigma^v_1}- \mathbf{x}_K & \mathbf{x}_{\sigma^v_3} - \mathbf{x}_K \end{bmatrix}. \end{equation}

This enables us to write the discrete flux across \(\sigma^v_1\) from cell \(K\) as follows:

\begin{equation} F_{K, \sigma^v_1} := - |\sigma^v_1| \mathbf{n}_{\sigma^v_1}^T \mathbf{\Lambda}_K \nabla_\mathcal{D}^{K_v} u. \label{eq:discreteFlux} \end{equation}

Inserting the discrete gradient yields

\begin{equation} F_{K, \sigma^v_1} = \omega_{K,\sigma^v_1\sigma^v_1}(u_K - u_{\sigma^v_1}) + \omega_{K,\sigma^v_1 \sigma^v_3}(u_K - u_{\sigma^v_3}), \label{eq:discreteFluxRef} \end{equation}

with \((\omega_{K,\sigma^v_1\sigma^v_1},\omega_{K,\sigma^v_1 \sigma^v_3})^T = |\sigma^v_1| \mathbb{D}^{-1}_{K_v}\mathbf{\Lambda}_K \mathbf{n}_{\sigma^v_1}\). These values are calculated in DuMux by using the function Dumux::computeMpfaTransmissibility.

To deduce a cell-centered scheme, the introduced face unknowns \(u_{\sigma^v_i}\) have to be eliminated. This is done by enforcing flux continuity for each sub-control volume face, i.e.

\begin{align} F_{K, \sigma^v_1} + F_{L, \sigma^v_1} &= 0, \\ F_{K, \sigma^v_3} + F_{M, \sigma^v_3} &= 0, \\ F_{L, \sigma^v_2} + F_{M, \sigma^v_2} &= 0. \end{align}

This results in a system of equations for the face unknowns \(\mathbf{u}_{\sigma}\)

\begin{equation} \mathbb{A}^{3\times 3} \mathbf{u}_{\sigma} = \mathbb{B}^{3\times 3} \mathbf{u}, \end{equation}

where \(\mathbf{u}\) contains the three cell unknowns \(u_K,u_L,u_M\) and \(\mathbf{u}_{\sigma}\) the three face unknowns \(u_{\sigma^v_1}, u_{\sigma^v_2}, u_{\sigma^v_3}\). Inserting these face unknowns into the flux expression yields

\begin{equation} F_{K,\sigma^v_i} = \sum_{N \in \lbrace K,L,M \rbrace } t^N_{K,\sigma^v_i} u_{N} = \mathbf{t}_{K,\sigma^v_i} \cdot \mathbf{u}, \end{equation}

for each cell \(K\) and sub-control volume face \(\sigma^v_i\).

Files

file  mpfa/computetransmissibility.hh
 This file contains free functions to evaluate the transmissibilities associated with flux evaluations across sub-control volume faces in the context of cell-centered Mpfa schemes.
 
file  cellcentered/mpfa/connectivitymap.hh
 Stores the face indices corresponding to the neighbors of an element that contribute to the derivative calculation. Depending on if an mpfa scheme leads to a symmetric/unsymmetric sparsity pattern, the adequate implementation of the connectivity map is chosen.
 
file  dualgridindexset.hh
 Class for the index sets of the dual grid in mpfa schemes.
 
file  discretization/cellcentered/mpfa/elementfluxvariablescache.hh
 The element-local object of flux variables caches.
 
file  cellcentered/mpfa/elementvolumevariables.hh
 The local (stencil) volume variables class for cell centered mpfa models.
 
file  discretization/cellcentered/mpfa/fvelementgeometry.hh
 Stencil-local finite volume geometry (scvs and scvfs) for cell-centered mpfa models This builds up the sub control volumes and sub control volume faces for each element in the local scope we are restricting to, e.g. stencil or element.
 
file  discretization/cellcentered/mpfa/fvgridgeometry.hh
 The finite volume geometry (scvs and scvfs) for cell-centered mpfa models on a grid view This builds up the sub control volumes and sub control volume faces for each element of the grid partition.
 
file  cellcentered/mpfa/fvgridgeometrytraits.hh
 Traits class to be used in conjunction with the CCMpfaFVGridGeometry.
 
file  discretization/cellcentered/mpfa/gridfluxvariablescache.hh
 Flux variable caches on a gridview.
 
file  gridinteractionvolumeindexsets.hh
 Class for the grid interaction volume index sets of mpfa schemes.
 
file  cellcentered/mpfa/gridvolumevariables.hh
 The grid volume variables class for cell centered mpfa models.
 
file  helper.hh
 Helper class to get data required for mpfa scheme.
 
file  interactionvolumebase.hh
 Base class for interaction volumes of mpfa methods.
 
file  interactionvolumedatahandle.hh
 Data handle class for interaction volumes of mpfa methods. This class is passed to interaction volumes to store the necessary data in it.
 
file  localassemblerbase.hh
 Defines the general interface of classes used for the assembly of the local systems of equations involved in the transmissibility computation in mpfa schemes.
 
file  localassemblerhelper.hh
 A class that contains helper functions as well as functionality which is common to different mpfa schemes and which solely operate on the interaction volume.
 
file  localfacedata.hh
 Data structure holding interaction volume-local information for a grid subb-control volume face embedded in it.
 
file  methods.hh
 The available mpfa schemes in Dumux.
 
file  discretization/cellcentered/mpfa/omethod/interactionvolume.hh
 Class for the interaction volume of the mpfa-o scheme.
 
file  interactionvolumeindexset.hh
 Class for the index set within an interaction volume of the mpfa-o scheme.
 
file  discretization/cellcentered/mpfa/omethod/localassembler.hh
 Class for the assembly of the local systems of equations involved in the transmissibility computation in an mpfa-o type manner.
 
file  discretization/cellcentered/mpfa/omethod/localsubcontrolentities.hh
 Classes for sub control entities of the mpfa-o method.
 
file  scvgeometryhelper.hh
 Helper class providing functionality to compute the geometry of the interaction-volume local sub-control volumes involved in the mpfa-o scheme.
 
file  staticinteractionvolume.hh
 Class for the interaction volume of the mpfa-o scheme to be used when the sizes are known at compile time, e.g. for non-boundary interaction volumes on structured grids.
 
file  scvgradients.hh
 Class providing functionality for the reconstruction of the gradients in the sub-control volumes involved in mpfa schemes.
 
file  discretization/cellcentered/mpfa/subcontrolvolumeface.hh
 The sub control volume face.
 

Classes

struct  Dumux::NodalIndexSetDefaultTraits< GV >
 Default traits to be used in conjunction with the dual grid nodal index set. More...
 
class  Dumux::CCMpfaDualGridNodalIndexSet< T >
 Nodal index set for mpfa schemes, constructed around grid vertices. More...
 
class  Dumux::CCMpfaDualGridIndexSet< NI >
 Class for the index sets of the dual grid in mpfa schemes. More...
 
struct  Dumux::InteractionVolumeDataStorage< PrimaryIV, PrimaryIVDataHandle, SecondaryIV, SecondaryIVDataHandle >
 Structure to store interaction volumes and data handles. More...
 
class  Dumux::CCMpfaElementFluxVariablesCache< GFVC, cachingEnabled >
 The flux variables caches for an element. More...
 
class  Dumux::CCMpfaElementFluxVariablesCache< GFVC, true >
 The flux variables caches for an element with caching enabled. More...
 
class  Dumux::CCMpfaElementFluxVariablesCache< GFVC, false >
 The flux variables caches for an element with caching disabled. More...
 
class  Dumux::CCMpfaElementVolumeVariables< GVV, cachingEnabled >
 The local (stencil) volume variables class for cell centered mpfa models. More...
 
class  Dumux::CCMpfaElementVolumeVariables< GVV, true >
 The local (stencil) volume variables class for cell centered mpfa models with caching. More...
 
class  Dumux::CCMpfaElementVolumeVariables< GVV, false >
 The local (stencil) volume variables class for cell centered tpfa models with caching. More...
 
class  Dumux::CCMpfaFVElementGeometry< GG, enableGridGeometryCache >
 Stencil-local finite volume geometry (scvs and scvfs) for cell-centered mpfa models This builds up the sub control volumes and sub control volume faces for each element in the local scope we are restricting to, e.g. stencil or element. More...
 
class  Dumux::CCMpfaFVElementGeometry< GG, true >
 Stencil-local finite volume geometry (scvs and scvfs) for cell-centered mpfa models Specialization for grid caching enabled. More...
 
class  Dumux::CCMpfaFVElementGeometry< GG, false >
 Stencil-local finite volume geometry (scvs and scvfs) for cell-centered TPFA models Specialization for grid caching disabled. More...
 
class  Dumux::CCMpfaFVGridGeometry< GridView, Traits, enableCache >
 The finite volume geometry (scvs and scvfs) for cell-centered mpfa models on a grid view This builds up the sub control volumes and sub control volume faces. More...
 
class  Dumux::CCMpfaFVGridGeometry< GV, Traits, true >
 The finite volume geometry (scvs and scvfs) for cell-centered mpfa models on a grid view This builds up the sub control volumes and sub control volume faces. More...
 
class  Dumux::CCMpfaFVGridGeometry< GV, Traits, false >
 The finite volume geometry (scvs and scvfs) for cell-centered mpfa models on a grid view This builds up the sub control volumes and sub control volume faces. More...
 
struct  Dumux::CCMpfaFVGridGeometryTraits< GV, NI, PIV, SIV >
 Traits class to be used for the CCMpfaFVGridGeometry. More...
 
struct  Dumux::IvDataHandlePhysicsTraits< ModelTraits >
 Data handle physics traits. More...
 
struct  Dumux::CCMpfaDefaultGridFluxVariablesCacheTraits< P, FVC, FVCF, PIV, SIV, PDH, SDH >
 Data handle physics traits. More...
 
class  Dumux::CCMpfaGridFluxVariablesCache< Traits, cachingEnabled >
 Flux variable caches on a gridview. More...
 
class  Dumux::CCMpfaGridFluxVariablesCache< TheTraits, true >
 Flux variable caches on a gridview with grid caching enabled. More...
 
class  Dumux::CCMpfaGridFluxVariablesCache< TheTraits, false >
 Flux variable caches on a gridview with grid caching disabled. More...
 
class  Dumux::CCMpfaGridInteractionVolumeIndexSets< FVG, NI, PI, SI >
 Class that holds all interaction volume index sets on a grid view. More...
 
struct  Dumux::CCMpfaDefaultGridVolumeVariablesTraits< P, VV >
 Traits for the base class for the grid volume variables. More...
 
class  Dumux::CCMpfaGridVolumeVariables< Problem, VolumeVariables, cachingEnabled, Traits >
 Base class for the grid volume variables. More...
 
class  Dumux::MpfaDimensionHelper< GridGeometry, dim, dimWorld >
 Dimension-specific helper class to get data required for mpfa scheme. More...
 
class  Dumux::MpfaDimensionHelper< GridGeometry, 2, 2 >
 Dimension-specific mpfa helper class for dim == 2 & dimWorld == 2. More...
 
class  Dumux::MpfaDimensionHelper< GridGeometry, 2, 3 >
 Dimension-specific mpfa helper class for dim == 2 & dimWorld == 3. Reuses some functionality of the specialization for dim = dimWorld = 2. More...
 
class  Dumux::MpfaDimensionHelper< GridGeometry, 3, 3 >
 Dimension-specific mpfa helper class for dim == 3 & dimWorld == 3. More...
 
class  Dumux::CCMpfaHelper< GridGeometry >
 Helper class to get the required information on an interaction volume. More...
 
class  Dumux::CCMpfaInteractionVolumeBase< T >
 Base class for the interaction volumes of mpfa methods. It defines the interface and actual implementations should derive from this class. More...
 
class  Dumux::CCMpfaDataHandleBases::SystemMatricesHandle< MVT, size1, size2 >
 Common base class to all handles. Stores arrays of the matrices involved in the interaction volume-local systems of equations. Apart from the transmissibility matrix we store those matrices that are needed e.g. for later face pressure reconstruction. The fluxes as well as the local systems of equations can be expressed as functions of the intermediate unknown face face values \(\bar{\mathbf{u}}\) and the known cell/Dirichlet values \(\mathbf{u}\) using the matrices \(\mathbf{A}\), \(\mathbf{B}\), \(\mathbf{C}\), \(\mathbf{D}\) and the vector of Neumann fluxes \(\mathbf{N}\) as follows: More...
 
class  Dumux::CCMpfaDataHandleBases::SystemVectorsHandle< MVT, size1, size2 >
 Common base class to all handles. Stores arrays of the vectors involved in the interaction volume-local systems of equations. More...
 
class  Dumux::AdvectionDataHandle< MatVecTraits, PhysicsTraits, EnableAdvection >
 Data handle for quantities related to advection. More...
 
class  Dumux::DiffusionDataHandle< MatVecTraits, PhysicsTraits, EnableDiffusion >
 Data handle for quantities related to diffusion. More...
 
class  Dumux::HeatConductionDataHandle< MatVecTraits, PhysicsTraits, enableHeatConduction >
 Data handle for quantities related to heat conduction. More...
 
class  Dumux::InteractionVolumeDataHandle< MVT, PT >
 Class for the interaction volume data handle. More...
 
class  Dumux::InteractionVolumeAssemblerBase< P, EG, EV >
 Defines the general interface of the local assembler classes for the assembly of the interaction volume-local transmissibility matrix. Specializations have to be provided for the available interaction volume implementations. these should derive from this base class. More...
 
class  Dumux::InteractionVolumeAssemblerHelper
 A class that contains helper functions as well as functionality which is common to different mpfa schemes and which solely operate on the interaction volume. More...
 
class  Dumux::InteractionVolumeLocalFaceData< GridIndexType, LocalIndexType >
 General implementation of a data structure holding interaction volume-local information for a grid sub-control volume face embedded in it. More...
 
struct  Dumux::CCMpfaODefaultInteractionVolumeTraits< NodalIndexSet, Scalar >
 The default interaction volume traits class for the mpfa-o method. This uses dynamic types types for matrices/vectors in order to work on general grids. For interaction volumes known at compile time use the static interaction volume implementation. More...
 
class  Dumux::CCMpfaOInteractionVolume< Traits >
 Forward declaration of the o-method's interaction volume. More...
 
class  Dumux::CCMpfaOInteractionVolumeIndexSet< DualGridNodalIndexSet >
 The interaction volume index set class for the mpfa-o scheme. More...
 
class  Dumux::MpfaOInteractionVolumeAssembler< P, EG, EV >
 Specialization of the interaction volume-local assembler class for the schemes using an mpfa-o type assembly. More...
 
class  Dumux::CCMpfaOInteractionVolumeLocalScv< IvIndexSet, Scalar, dim, dimWorld >
 Class for the interaction volume-local sub-control volume used in the mpfa-o scheme. More...
 
struct  Dumux::CCMpfaOInteractionVolumeLocalScvf< IvIndexSet >
 Class for the interaction volume-local sub-control volume face used in the mpfa-o scheme. More...
 
class  Dumux::CCMpfaOScvGeometryHelper< LocalScvType >
 Helper class providing functionality to compute the geometry of the interaction-volume local sub-control volumes of mpfa-o type. More...
 
struct  Dumux::CCMpfaODefaultStaticInteractionVolumeTraits< NI, S, C, F >
 The default interaction volume traits class for the mpfa-o method with known size of the interaction volumes at compile time. It uses statically sized containers for the iv-local data structures and static matrices and vectors. More...
 
class  Dumux::CCMpfaOStaticInteractionVolume< Traits >
 Forward declaration of the o-method's static interaction volume. More...
 
class  Dumux::CCMpfaScvGradients
 Class providing functionality for the reconstruction of the gradients in the sub-control volumes involved in mpfa schemes. More...
 
struct  Dumux::CCMpfaDefaultScvfGeometryTraits< GridView >
 Default traits class to be used for the sub-control volume faces for the cell-centered finite volume scheme using MPFA. More...
 
class  Dumux::CCMpfaSubControlVolumeFace< GV, T >
 Class for a sub control volume face in mpfa methods, i.e a part of the boundary of a control volume we compute fluxes on. More...
 

Enumerations

enum class  Dumux::MpfaMethods : unsigned int { Dumux::MpfaMethods::oMethod }
 The available mpfa schemes in Dumux. More...
 

Functions

template<class EG , class IVSubControlVolume , class Tensor >
Dune::FieldVector< typename Tensor::field_type, IVSubControlVolume::myDimension > Dumux::computeMpfaTransmissibility (const EG &fvGeometry, const IVSubControlVolume &scv, const typename EG::SubControlVolumeFace &scvf, const Tensor &t, typename IVSubControlVolume::ctype extrusionFactor)
 Free function to evaluate the Mpfa transmissibility associated with the flux (in the form of flux = t*gradU) across a sub-control volume face stemming from a given sub-control volume with corresponding tensor t. More...
 
template<class EG , class IVSubControlVolume , class Tensor , std::enable_if_t< Dune::IsNumber< Tensor >::value, int > = 1>
Dune::FieldVector< Tensor, IVSubControlVolume::myDimension > Dumux::computeMpfaTransmissibility (const EG &fvGeometry, const IVSubControlVolume &scv, const typename EG::SubControlVolumeFace &scvf, const Tensor &t, typename IVSubControlVolume::ctype extrusionFactor)
 Free function to evaluate the Mpfa transmissibility associated with the flux (in the form of flux = t*gradU) across a sub-control volume face stemming from a given sub-control volume with corresponding tensor t, where t is a scalar. More...
 
template<class FVElementGeometry >
std::size_t Dumux::CCMpfa::maxNumBoundaryVolVars (const FVElementGeometry &fvGeometry)
 Computes how many boundary vol vars come into play for flux calculations on an element (for a given element finite volume geometry). This number here is probably always higher than the actually needed number of volume variables. However, we want to make sure it is high enough so that enough memory is reserved in the element volume variables below. More...
 
template<class VolumeVariables , class IndexType , class Problem , class FVElemGeom , class NodalIndexSet >
void Dumux::CCMpfa::addBoundaryVolVarsAtNode (std::vector< VolumeVariables > &volVars, std::vector< IndexType > &volVarIndices, const Problem &problem, const typename FVElemGeom::GridGeometry::GridView::template Codim< 0 >::Entity &element, const FVElemGeom &fvGeometry, const NodalIndexSet &nodalIndexSet)
 Adds the boundary volume variables found within a given nodal index set into the provided containers and stores the indices associated with them. More...
 
template<class VolumeVariables , class IndexType , class Problem , class FVElemGeom >
void Dumux::CCMpfa::addBoundaryVolVars (std::vector< VolumeVariables > &volVars, std::vector< IndexType > &volVarIndices, const Problem &problem, const typename FVElemGeom::GridGeometry::GridView::template Codim< 0 >::Entity &element, const FVElemGeom &fvGeometry)
 Adds the boundary volume variables found within the stencil to the provided containers and stores the indices associated with them. More...
 

Enumeration Type Documentation

◆ MpfaMethods

enum class Dumux::MpfaMethods : unsigned int
strong
Enumerator
oMethod 

Function Documentation

◆ addBoundaryVolVars()

template<class VolumeVariables , class IndexType , class Problem , class FVElemGeom >
void Dumux::CCMpfa::addBoundaryVolVars ( std::vector< VolumeVariables > &  volVars,
std::vector< IndexType > &  volVarIndices,
const Problem &  problem,
const typename FVElemGeom::GridGeometry::GridView::template Codim< 0 >::Entity &  element,
const FVElemGeom &  fvGeometry 
)
Parameters
volVarsThe container where the volume variables are stored
volVarIndicesThe container where the volume variable indices are stored
problemThe problem containing the Dirichlet boundary conditions
elementThe element to which the finite volume geometry was bound
fvGeometryThe element finite volume geometry

◆ addBoundaryVolVarsAtNode()

template<class VolumeVariables , class IndexType , class Problem , class FVElemGeom , class NodalIndexSet >
void Dumux::CCMpfa::addBoundaryVolVarsAtNode ( std::vector< VolumeVariables > &  volVars,
std::vector< IndexType > &  volVarIndices,
const Problem &  problem,
const typename FVElemGeom::GridGeometry::GridView::template Codim< 0 >::Entity &  element,
const FVElemGeom &  fvGeometry,
const NodalIndexSet &  nodalIndexSet 
)
Note
It only adds those boundary vol vars that do not live on scvfs that are inside the bound element. These have to be added separately
Parameters
volVarsThe container where the volume variables are stored
volVarIndicesThe container where the volume variable indices are stored
problemThe problem containing the Dirichlet boundary conditions
elementThe element to which the finite volume geometry is bound
fvGeometryThe element finite volume geometry
nodalIndexSetThe dual grid index set around a node

◆ computeMpfaTransmissibility() [1/2]

template<class EG , class IVSubControlVolume , class Tensor >
Dune::FieldVector< typename Tensor::field_type, IVSubControlVolume::myDimension > Dumux::computeMpfaTransmissibility ( const EG &  fvGeometry,
const IVSubControlVolume &  scv,
const typename EG::SubControlVolumeFace &  scvf,
const Tensor &  t,
typename IVSubControlVolume::ctype  extrusionFactor 
)
Parameters
fvGeometryThe element-centered control volume geometry
scvThe iv-local sub-control volume
scvfThe grid sub-control volume face
tThe tensor living in the scv
extrusionFactorThe extrusion factor of the scv

◆ computeMpfaTransmissibility() [2/2]

template<class EG , class IVSubControlVolume , class Tensor , std::enable_if_t< Dune::IsNumber< Tensor >::value, int > = 1>
Dune::FieldVector< Tensor, IVSubControlVolume::myDimension > Dumux::computeMpfaTransmissibility ( const EG &  fvGeometry,
const IVSubControlVolume &  scv,
const typename EG::SubControlVolumeFace &  scvf,
const Tensor &  t,
typename IVSubControlVolume::ctype  extrusionFactor 
)
Parameters
fvGeometryThe element-centered control volume geometry
scvThe iv-local sub-control volume
scvfThe grid sub-control volume face
tthe scalar quantity living in the scv
extrusionFactorThe extrusion factor of the scv

◆ maxNumBoundaryVolVars()

template<class FVElementGeometry >
std::size_t Dumux::CCMpfa::maxNumBoundaryVolVars ( const FVElementGeometry &  fvGeometry)
Parameters
fvGeometrythe element finite volume geometry