12#ifndef DUMUX_DISCRETIZATION_FACECENTERED_ELEMENTVOLUMEVARIABLES_HH 
   13#define DUMUX_DISCRETIZATION_FACECENTERED_ELEMENTVOLUMEVARIABLES_HH 
   28template<
class Gr
idVolumeVariables, 
bool cachingEnabled>
 
   39    using GridGeometry = std::decay_t<decltype(std::declval<GVV>().problem().gridGeometry())>;
 
   40    using FVElementGeometry = 
typename GridGeometry::LocalView;
 
   41    using SubControlVolume = 
typename GridGeometry::SubControlVolume;
 
   45    using GridVolumeVariables = GVV;
 
   48    using VolumeVariables = 
typename GridVolumeVariables::VolumeVariables;
 
   51    : gridVolumeVariablesPtr_(&gridVolumeVariables)
 
   52    , numScv_(gridVolumeVariables.problem().gridGeometry().numScv())
 
   56    const VolumeVariables& operator [](
const SubControlVolume& scv)
 const 
   58        if (scv.index() < numScv_)
 
   59            return gridVolVars().volVars(scv.index());
 
   61            return boundaryVolumeVariables_[getLocalIdx_(scv.index())];
 
   66    const VolumeVariables& operator [](
const std::size_t scvIdx)
 const 
   69            return gridVolVars().volVars(scvIdx);
 
   71            return boundaryVolumeVariables_[getLocalIdx_(scvIdx)];
 
   79    template<
class SolutionVector>
 
   80    ThisType bind(
const typename FVElementGeometry::Element& element,
 
   81                  const FVElementGeometry& fvGeometry,
 
   82                  const SolutionVector& sol) &&
 
   84        this->bind(element, fvGeometry, sol);
 
   85        return std::move(*
this);
 
   90    template<
class SolutionVector>
 
   91    void bind(
const typename FVElementGeometry::Element& element,
 
   92              const FVElementGeometry& fvGeometry,
 
   93              const SolutionVector& sol) &
 
   95        if (!fvGeometry.hasBoundaryScvf())
 
  100        boundaryVolVarIndices_.reserve(fvGeometry.numScvf()-
element.subEntities(1));
 
  101        boundaryVolumeVariables_.reserve(fvGeometry.numScvf()-
element.subEntities(1));
 
  103        for (
const auto& scvf : scvfs(fvGeometry))
 
  105            if (!scvf.boundary() || scvf.isFrontal() || scvf.processorBoundary())
 
  109            const auto& problem = gridVolVars().problem();
 
  110            const auto bcTypes = problem.boundaryTypes(element, scvf);
 
  114                const auto& scvI = fvGeometry.scv(scvFace.insideScvIdx());
 
  115                typename VolumeVariables::PrimaryVariables pv(
 
  116                    problem.dirichlet(element, scvFace)[scvI.dofAxis()]
 
  120                VolumeVariables volVars;
 
  121                volVars.update(dirichletPriVars, problem, element, scvI);
 
  123                boundaryVolumeVariables_.emplace_back(std::move(volVars));
 
  124                boundaryVolVarIndices_.push_back(scvFace.outsideScvIdx());
 
  127            if (bcTypes.hasDirichlet())
 
  134            if (
const auto& orthogonalScvf = fvGeometry.lateralOrthogonalScvf(scvf); orthogonalScvf.boundary())
 
  135                if (problem.boundaryTypes(element, orthogonalScvf).hasDirichlet())
 
  140        assert(boundaryVolumeVariables_.size() == boundaryVolVarIndices_.size());
 
  148    template<
class SolutionVector>
 
  149    ThisType bindElement(
const typename FVElementGeometry::Element& element,
 
  150                         const FVElementGeometry& fvGeometry,
 
  151                         const SolutionVector& sol) &&
 
  153        this->bindElement(element, fvGeometry, sol);
 
  154        return std::move(*
this);
 
  158    template<
class SolutionVector>
 
  159    void bindElement(
const typename FVElementGeometry::Element& element,
 
  160                     const FVElementGeometry& fvGeometry,
 
  161                     const SolutionVector& sol) &
 
  166    const GridVolumeVariables& gridVolVars()
 const 
  167    { 
return *gridVolumeVariablesPtr_; }
 
  170    bool hasVolVars(
const std::size_t scvIdx)
 const 
  172        if (scvIdx < numScv_)
 
  176            const auto it = std::find(boundaryVolVarIndices_.begin(), boundaryVolVarIndices_.end(), scvIdx);
 
  177            return it != boundaryVolVarIndices_.end();
 
  185        boundaryVolVarIndices_.clear();
 
  186        boundaryVolumeVariables_.clear();
 
  190    int getLocalIdx_(
const std::size_t volVarIdx)
 const 
  192        const auto it = std::find(boundaryVolVarIndices_.begin(), boundaryVolVarIndices_.end(), volVarIdx);
 
  193        assert(it != boundaryVolVarIndices_.end() && 
"Could not find the current volume variables for volVarIdx!");
 
  194        return std::distance(boundaryVolVarIndices_.begin(), it);
 
  197    const GridVolumeVariables* gridVolumeVariablesPtr_;
 
  198    const std::size_t numScv_;
 
  199    std::vector<std::size_t> boundaryVolVarIndices_;
 
  200    std::vector<VolumeVariables> boundaryVolumeVariables_;
 
  210    using ThisType = FaceCenteredStaggeredElementVolumeVariables<GVV, 
false>;
 
  211    using GridGeometry = std::decay_t<decltype(std::declval<GVV>().problem().gridGeometry())>;
 
  212    using FVElementGeometry = 
typename GridGeometry::LocalView;
 
  213    using SubControlVolume = 
typename GridGeometry::SubControlVolume;
 
  215    static constexpr auto dim = GridGeometry::GridView::dimension;
 
  216    static constexpr auto numInsideVolVars = dim * 2;
 
  217    static constexpr auto numOutsideVolVars = numInsideVolVars * 2 * (dim - 1);
 
  221    using GridVolumeVariables = GVV;
 
  224    using VolumeVariables = 
typename GridVolumeVariables::VolumeVariables;
 
  226    FaceCenteredStaggeredElementVolumeVariables(
const GridVolumeVariables& globalFacesVars)
 
  227    : gridVolumeVariablesPtr_(&globalFacesVars) {}
 
  230    const VolumeVariables& operator [](
const SubControlVolume& scv)
 const 
  231    { 
return volumeVariables_[getLocalIdx_(scv.index())]; }
 
  234    const VolumeVariables& operator [](
const std::size_t scvIdx)
 const 
  235    { 
return volumeVariables_[getLocalIdx_(scvIdx)]; }
 
  238    VolumeVariables& operator [](
const SubControlVolume& scv)
 
  239    { 
return volumeVariables_[getLocalIdx_(scv.index())]; }
 
  242    VolumeVariables& operator [](
const std::size_t scvIdx)
 
  243    { 
return volumeVariables_[getLocalIdx_(scvIdx)]; }
 
  250    template<
class SolutionVector>
 
  251    ThisType bind(
const typename FVElementGeometry::Element& element,
 
  252                  const FVElementGeometry& fvGeometry,
 
  253                  const SolutionVector& sol) &&
 
  255        this->bind_(element, fvGeometry, sol);
 
  256        return std::move(*
this);
 
  259    template<
class SolutionVector>
 
  260    void bind(
const typename FVElementGeometry::Element& element,
 
  261              const FVElementGeometry& fvGeometry,
 
  262              const SolutionVector& sol) &
 
  264        this->bind_(element, fvGeometry, sol);
 
  272    template<
class SolutionVector>
 
  273    ThisType bindElement(
const typename FVElementGeometry::Element& element,
 
  274                         const FVElementGeometry& fvGeometry,
 
  275                         const SolutionVector& sol) &&
 
  277        this->bindElement_(element, fvGeometry, sol);
 
  278        return std::move(*
this);
 
  281    template<
class SolutionVector>
 
  282    void bindElement(
const typename FVElementGeometry::Element& element,
 
  283                     const FVElementGeometry& fvGeometry,
 
  284                     const SolutionVector& sol) &
 
  285    { this->bindElement_(element, fvGeometry, sol); }
 
  288    const GridVolumeVariables& gridVolVars()
 const 
  289    { 
return *gridVolumeVariablesPtr_; }
 
  292    bool hasVolVars(
const std::size_t scvIdx)
 const 
  293    { 
return volVarsInserted_(scvIdx); }
 
  298    template<
class SolutionVector>
 
  299    void bind_(
const typename FVElementGeometry::Element& element,
 
  300               const FVElementGeometry& fvGeometry,
 
  301               const SolutionVector& sol)
 
  305        const auto& problem = gridVolVars().problem();
 
  306        const auto& gridGeometry = fvGeometry.gridGeometry();
 
  308        volVarIndices_.reserve(numInsideVolVars + numInsideVolVars);
 
  309        volumeVariables_.reserve(numInsideVolVars + numInsideVolVars);
 
  311        for (
const auto& scv : scvs(fvGeometry))
 
  313            for (
const auto otherScvIdx : gridGeometry.connectivityMap()[scv.index()])
 
  315                if (!volVarsInserted_(otherScvIdx))
 
  317                    const auto& otherScv = fvGeometry.scv(otherScvIdx);
 
  318                    volVarIndices_.push_back(otherScvIdx);
 
  319                    volumeVariables_.emplace_back();
 
  320                    const auto& otherElement = gridGeometry.element(otherScv.elementIndex());
 
  321                    volumeVariables_.back().update(
 
  323                        problem, otherElement, otherScv
 
  329        if (fvGeometry.hasBoundaryScvf())
 
  331            for (
const auto& scvf : scvfs(fvGeometry))
 
  333                if (!scvf.boundary() || scvf.isFrontal())
 
  337                const auto& problem = gridVolVars().problem();
 
  338                const auto bcTypes = problem.boundaryTypes(element, scvf);
 
  342                    const auto& scvI = fvGeometry.scv(scvFace.insideScvIdx());
 
  343                    typename VolumeVariables::PrimaryVariables pv(
 
  344                        problem.dirichlet(element, scvFace)[scvI.dofAxis()]
 
  348                    VolumeVariables volVars;
 
  349                    volVars.update(dirichletPriVars,
 
  354                    volumeVariables_.emplace_back(std::move(volVars));
 
  355                    volVarIndices_.push_back(scvFace.outsideScvIdx());
 
  358                if (bcTypes.hasDirichlet())
 
  365                if (
const auto& orthogonalScvf = fvGeometry.lateralOrthogonalScvf(scvf); orthogonalScvf.boundary())
 
  366                    if (problem.boundaryTypes(element, orthogonalScvf).hasDirichlet())
 
  375    template<
class SolutionVector>
 
  376    void bindElement_(
const typename FVElementGeometry::Element& element,
 
  377                      const FVElementGeometry& fvGeometry,
 
  378                      const SolutionVector& sol)
 
  381        const auto& problem = gridVolVars().problem();
 
  382        const auto& gridGeometry = fvGeometry.gridGeometry();
 
  383        volVarIndices_.reserve(numInsideVolVars);
 
  385        for (
const auto& scv : scvs(fvGeometry))
 
  387            volVarIndices_.push_back(scv.index());
 
  388            volumeVariables_.emplace_back();
 
  389            volumeVariables_.back().update(
 
  391                problem, element, scv
 
  399        volVarIndices_.clear();
 
  400        volumeVariables_.clear();
 
  403    bool volVarsInserted_(
const std::size_t scvIdx)
 const 
  405        return std::find(volVarIndices_.begin(), volVarIndices_.end(), scvIdx) != volVarIndices_.end();
 
  408    int getLocalIdx_(
const int scvfIdx)
 const 
  410        const auto it = std::find(volVarIndices_.begin(), volVarIndices_.end(), scvfIdx);
 
  411        assert(it != volVarIndices_.end() && 
"Could not find the current face variables for scvfIdx!");
 
  412        return std::distance(volVarIndices_.begin(), it);
 
  415    const GridVolumeVariables* gridVolumeVariablesPtr_;
 
  416    std::vector<std::size_t> volVarIndices_;
 
  417    std::vector<VolumeVariables> volumeVariables_;
 
Base class for the face variables vector.
Definition facecentered/staggered/elementvolumevariables.hh:29
Element solution classes and factory functions.
void 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...
Definition cellcentered/mpfa/elementvolumevariables.hh:124
auto elementSolution(const Element &element, const SolutionVector &sol, const GridGeometry &gg) -> std::enable_if_t< GridGeometry::discMethod==DiscretizationMethods::cctpfa||GridGeometry::discMethod==DiscretizationMethods::ccmpfa, CCElementSolution< typename GridGeometry::LocalView, std::decay_t< decltype(std::declval< SolutionVector >()[0])> > >
Make an element solution for cell-centered schemes.
Definition cellcentered/elementsolution.hh:101