A cell-centered finite volume scheme with multi-point flux approximation. More...
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}
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... | |
|
strong |
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 | ||
) |
volVars | The container where the volume variables are stored |
volVarIndices | The container where the volume variable indices are stored |
problem | The problem containing the Dirichlet boundary conditions |
element | The element to which the finite volume geometry was bound |
fvGeometry | The element finite volume geometry |
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 | ||
) |
volVars | The container where the volume variables are stored |
volVarIndices | The container where the volume variable indices are stored |
problem | The problem containing the Dirichlet boundary conditions |
element | The element to which the finite volume geometry is bound |
fvGeometry | The element finite volume geometry |
nodalIndexSet | The dual grid index set around a node |
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 | ||
) |
fvGeometry | The element-centered control volume geometry |
scv | The iv-local sub-control volume |
scvf | The grid sub-control volume face |
t | The tensor living in the scv |
extrusionFactor | The extrusion factor of the scv |
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 | ||
) |
fvGeometry | The element-centered control volume geometry |
scv | The iv-local sub-control volume |
scvf | The grid sub-control volume face |
t | the scalar quantity living in the scv |
extrusionFactor | The extrusion factor of the scv |
std::size_t Dumux::CCMpfa::maxNumBoundaryVolVars | ( | const FVElementGeometry & | fvGeometry | ) |
fvGeometry | the element finite volume geometry |