12#ifndef DUMUX_DISCRETIZATION_CC_MPFA_O_INTERACTIONVOLUME_HH 
   13#define DUMUX_DISCRETIZATION_CC_MPFA_O_INTERACTIONVOLUME_HH 
   17#include <dune/common/dynmatrix.hh> 
   18#include <dune/common/dynvector.hh> 
   19#include <dune/common/fvector.hh> 
   20#include <dune/common/reservedvector.hh> 
   46template< 
class NodalIndexSet, 
class Scalar >
 
   50    using GridIndexType = 
typename NodalIndexSet::GridIndexType;
 
   51    using LocalIndexType = 
typename NodalIndexSet::LocalIndexType;
 
   53    static constexpr int dim = NodalIndexSet::Traits::GridView::dimension;
 
   54    static constexpr int dimWorld = NodalIndexSet::Traits::GridView::dimensionworld;
 
   56    using DimVector = Dune::FieldVector<Scalar, dim>;
 
   57    using FaceOmegas = 
typename std::conditional< (dim<dimWorld),
 
   58                                                  std::vector<DimVector>,
 
   59                                                  Dune::ReservedVector<DimVector, 2> >::type;
 
   64        using OmegaStorage = std::vector< FaceOmegas >;
 
   66        using AMatrix = Dune::DynamicMatrix< Scalar >;
 
   67        using BMatrix = Dune::DynamicMatrix< Scalar >;
 
   68        using CMatrix = Dune::DynamicMatrix< Scalar >;
 
   69        using DMatrix = Dune::DynamicMatrix< Scalar >;
 
   70        using TMatrix = Dune::DynamicMatrix< Scalar >;
 
   71        using CellVector = Dune::DynamicVector< Scalar >;
 
   72        using FaceVector = Dune::DynamicVector< Scalar >;
 
   77    using GridView = 
typename NodalIndexSet::Traits::GridView;
 
   90    template<
class Problem, 
class FVElementGeometry, 
class ElemVolVars>
 
 
  100template< 
class Traits >
 
  104    using GridView = 
typename Traits::GridView;
 
  105    using Element = 
typename GridView::template Codim<0>::Entity;
 
  107    using IndexSet = 
typename Traits::IndexSet;
 
  108    using GridIndexType = 
typename IndexSet::GridIndexType;
 
  109    using LocalIndexType = 
typename IndexSet::LocalIndexType;
 
  110    using Stencil = 
typename IndexSet::NodalGridStencilType;
 
  112    using LocalScvType = 
typename Traits::LocalScvType;
 
  113    using LocalScvfType = 
typename Traits::LocalScvfType;
 
  114    using LocalFaceData = 
typename Traits::LocalFaceData;
 
  121        GridIndexType volVarIndex_;
 
 
  134    template< 
class Problem, 
class FVElementGeometry >
 
  135    void bind(
const IndexSet& indexSet,
 
  136              const Problem& problem,
 
  137              const FVElementGeometry& fvGeometry)
 
  141        stencil_ = &indexSet.nodalIndexSet().gridScvIndices();
 
  144        numFaces_ = indexSet.numFaces();
 
  145        const auto numLocalScvs = indexSet.numScvs();
 
  146        const auto numGlobalScvfs = indexSet.nodalIndexSet().numScvfs();
 
  149        elements_.clear();      elements_.reserve(numLocalScvs);
 
  150        scvs_.clear();          scvs_.reserve(numLocalScvs);
 
  151        scvfs_.clear();         scvfs_.reserve(numFaces_);
 
  152        localFaceData_.clear(); localFaceData_.reserve(numGlobalScvfs);
 
  153        dirichletData_.clear(); dirichletData_.reserve(numFaces_);
 
  156        for (LocalIndexType scvIdxLocal = 0; scvIdxLocal < numLocalScvs; scvIdxLocal++)
 
  158            elements_.emplace_back(fvGeometry.gridGeometry().element( 
stencil()[scvIdxLocal] ));
 
  159            scvs_.emplace_back(fvGeometry.gridGeometry().mpfaHelper(),
 
  161                               fvGeometry.scv( 
stencil()[scvIdxLocal] ),
 
  168        numKnowns_ = numLocalScvs;
 
  171        for (LocalIndexType faceIdxLocal = 0; faceIdxLocal < numFaces_; ++faceIdxLocal)
 
  173            const auto& scvf = fvGeometry.scvf(indexSet.gridScvfIndex(faceIdxLocal));
 
  176            const auto& neighborScvIndicesLocal = indexSet.neighboringLocalScvIndices(faceIdxLocal);
 
  177            const auto numNeighborScvs = neighborScvIndicesLocal.size();
 
  178            localFaceData_.emplace_back(faceIdxLocal, neighborScvIndicesLocal[0], scvf.index());
 
  183                if (problem.boundaryTypes(elements_[neighborScvIndicesLocal[0]], scvf).hasOnlyDirichlet())
 
  185                    scvfs_.emplace_back(scvf, neighborScvIndicesLocal, numKnowns_++, 
true);
 
  186                    dirichletData_.emplace_back(scvf.outsideScvIdx());
 
  189                    scvfs_.emplace_back(scvf, neighborScvIndicesLocal, numUnknowns_++, 
false);
 
  193                scvfs_.emplace_back(scvf, neighborScvIndicesLocal, numUnknowns_++, 
false);
 
  196                for (LocalIndexType i = 1; i < numNeighborScvs; ++i)
 
  199                    const auto outsideLocalScvIdx = neighborScvIndicesLocal[i];
 
  200                    const auto& flipScvfIndex = fvGeometry.gridGeometry().flipScvfIndexSet()[scvf.index()][i-1];
 
  201                    const auto& flipScvf = fvGeometry.scvf(flipScvfIndex);
 
  202                    localFaceData_.emplace_back(faceIdxLocal,       
 
 
  213    { 
return numFaces_; }
 
 
  217    { 
return numUnknowns_; }
 
 
  221    { 
return numKnowns_; }
 
 
  225    { 
return scvs_.size(); }
 
 
  229    { 
return *stencil_; }
 
 
  232    const Element& 
element(LocalIndexType ivLocalScvIdx)
 const 
  233    { 
return elements_[ivLocalScvIdx]; }
 
 
  236    const LocalScvfType& 
localScvf(LocalIndexType ivLocalScvfIdx)
 const 
  237    { 
return scvfs_[ivLocalScvfIdx]; }
 
 
  240    const LocalScvType& 
localScv(LocalIndexType ivLocalScvIdx)
 const 
  241    { 
return scvs_[ivLocalScvIdx]; }
 
 
  245    { 
return localFaceData_; }
 
 
  249    { 
return dirichletData_; }
 
 
  252    template< 
class FVElementGeometry >
 
  253    auto getScvGeometry(LocalIndexType ivLocalScvIdx, 
const FVElementGeometry& fvGeometry)
 const 
 
  263    template< 
class IvIndexSetContainer,
 
  266              class FlipScvfIndexSet >
 
  268                               ScvfIndexMap& scvfIndexMap,
 
  269                               const NodalIndexSet& nodalIndexSet,
 
  270                               const FlipScvfIndexSet& flipScvfIndexSet)
 
  273        const auto curGlobalIndex = ivIndexSetContainer.size();
 
  276        ivIndexSetContainer.emplace_back(nodalIndexSet, flipScvfIndexSet);
 
  279        for (
const auto scvfIdx : nodalIndexSet.gridScvfIndices())
 
  280            scvfIndexMap[scvfIdx] = curGlobalIndex;
 
 
  285    const Stencil* stencil_;
 
  288    std::vector<Element> elements_;
 
  289    std::vector<LocalScvType> scvs_;
 
  290    std::vector<LocalScvfType> scvfs_;
 
  291    std::vector<LocalFaceData> localFaceData_;
 
  292    std::vector<DirichletData> dirichletData_;
 
  295    std::size_t numFaces_;
 
  296    std::size_t numUnknowns_;
 
  297    std::size_t numKnowns_;
 
 
Base class for the interaction volumes of mpfa methods. It defines the interface and actual implement...
Definition interactionvolumebase.hh:56
GridIndexType volVarIndex() const
Return corresponding vol var index.
Definition discretization/cellcentered/mpfa/omethod/interactionvolume.hh:127
DirichletData(const GridIndexType index)
Constructor.
Definition discretization/cellcentered/mpfa/omethod/interactionvolume.hh:124
Forward declaration of the o-method's interaction volume.
Definition discretization/cellcentered/mpfa/omethod/interactionvolume.hh:103
const LocalScvfType & localScvf(LocalIndexType ivLocalScvfIdx) const
returns the local scvf entity corresponding to a given iv-local scvf idx
Definition discretization/cellcentered/mpfa/omethod/interactionvolume.hh:236
std::size_t numUnknowns() const
returns the number of intermediate unknowns within this interaction volume
Definition discretization/cellcentered/mpfa/omethod/interactionvolume.hh:216
static constexpr std::size_t numIVAtVertex(const NI &nodalIndexSet)
returns the number of interaction volumes living around a vertex
Definition discretization/cellcentered/mpfa/omethod/interactionvolume.hh:258
void bind(const IndexSet &indexSet, const Problem &problem, const FVElementGeometry &fvGeometry)
Sets up the local scope for a given iv index set.
Definition discretization/cellcentered/mpfa/omethod/interactionvolume.hh:135
static constexpr MpfaMethods MpfaMethod
publicly state the mpfa-scheme this interaction volume is associated with
Definition discretization/cellcentered/mpfa/omethod/interactionvolume.hh:131
static void addIVIndexSets(IvIndexSetContainer &ivIndexSetContainer, ScvfIndexMap &scvfIndexMap, const NodalIndexSet &nodalIndexSet, const FlipScvfIndexSet &flipScvfIndexSet)
Definition discretization/cellcentered/mpfa/omethod/interactionvolume.hh:267
std::size_t numKnowns() const
returns the number of (in this context) known solution values within this interaction volume
Definition discretization/cellcentered/mpfa/omethod/interactionvolume.hh:220
const Element & element(LocalIndexType ivLocalScvIdx) const
returns the grid element corresponding to a given iv-local scv idx
Definition discretization/cellcentered/mpfa/omethod/interactionvolume.hh:232
const LocalScvType & localScv(LocalIndexType ivLocalScvIdx) const
returns the local scv entity corresponding to a given iv-local scv idx
Definition discretization/cellcentered/mpfa/omethod/interactionvolume.hh:240
const std::vector< LocalFaceData > & localFaceData() const
returns a reference to the container with the local face data
Definition discretization/cellcentered/mpfa/omethod/interactionvolume.hh:244
const Stencil & stencil() const
returns the cell-stencil of this interaction volume
Definition discretization/cellcentered/mpfa/omethod/interactionvolume.hh:228
const std::vector< DirichletData > & dirichletData() const
returns a reference to the information container on Dirichlet BCs within this iv
Definition discretization/cellcentered/mpfa/omethod/interactionvolume.hh:248
std::size_t numScvs() const
returns the number of scvs embedded in this interaction volume
Definition discretization/cellcentered/mpfa/omethod/interactionvolume.hh:224
std::size_t numFaces() const
returns the number of primary scvfs of this interaction volume
Definition discretization/cellcentered/mpfa/omethod/interactionvolume.hh:212
auto getScvGeometry(LocalIndexType ivLocalScvIdx, const FVElementGeometry &fvGeometry) const
returns the geometry of the i-th local scv
Definition discretization/cellcentered/mpfa/omethod/interactionvolume.hh:253
static ScvGeometry computeScvGeometry(LocalIndexType ivLocalScvIdx, const InteractionVolume &iv, const FVElementGeometry &fvGeometry)
returns the geometry of the i-th local scv
Definition scvgeometryhelper.hh:61
General implementation of a data structure holding interaction volume-local information for a grid su...
Definition localfacedata.hh:30
Specialization of the interaction volume-local assembler class for the schemes using an mpfa-o type a...
Definition discretization/cellcentered/mpfa/omethod/localassembler.hh:36
Class for the assembly of the local systems of equations involved in the transmissibility computation...
Classes for sub control entities of the mpfa-o method.
MpfaMethods
The available mpfa schemes in Dumux.
Definition methods.hh:22
@ oMethod
Definition methods.hh:23
Base class for interaction volumes of mpfa methods.
Class for the index set within an interaction volume of the mpfa-o scheme.
Data structure holding interaction volume-local information for a grid subb-control volume face embed...
The available mpfa schemes in Dumux.
Helper class providing functionality to compute the geometry of the interaction-volume local sub-cont...
The default interaction volume traits class for the mpfa-o method. This uses dynamic types types for ...
Definition discretization/cellcentered/mpfa/omethod/interactionvolume.hh:48
CCMpfaOInteractionVolumeIndexSet< NodalIndexSet > IndexSet
export the type for the interaction volume index set
Definition discretization/cellcentered/mpfa/omethod/interactionvolume.hh:79
MVTraits MatVecTraits
export the matrix/vector traits to be used by the iv
Definition discretization/cellcentered/mpfa/omethod/interactionvolume.hh:87
typename NodalIndexSet::Traits::GridView GridView
export the type of grid view
Definition discretization/cellcentered/mpfa/omethod/interactionvolume.hh:77
MpfaOInteractionVolumeAssembler< Problem, FVElementGeometry, ElemVolVars > LocalAssembler
the type of assembler used for the o-method's iv-local eq systems
Definition discretization/cellcentered/mpfa/omethod/interactionvolume.hh:91
CCMpfaOInteractionVolumeLocalScv< IndexSet, Scalar, dim, dimWorld > LocalScvType
export the type of interaction-volume local scvs
Definition discretization/cellcentered/mpfa/omethod/interactionvolume.hh:81
CCMpfaOInteractionVolumeLocalScvf< IndexSet > LocalScvfType
export the type of interaction-volume local scvfs
Definition discretization/cellcentered/mpfa/omethod/interactionvolume.hh:83
InteractionVolumeLocalFaceData< GridIndexType, LocalIndexType > LocalFaceData
export the type of used for the iv-local face data
Definition discretization/cellcentered/mpfa/omethod/interactionvolume.hh:85