12#ifndef DUMUX_PNM_2P_GRID_FLUXVARSCACHE_HH 
   13#define DUMUX_PNM_2P_GRID_FLUXVARSCACHE_HH 
   27template<
class P, 
class FVC, 
class IS = TwoPInvasionState<P>>
 
   34    template<
class Gr
idFluxVariablesCache, 
bool cachingEnabled>
 
 
   43template<
class Problem,
 
   44         class FluxVariablesCache,
 
   54template<
class P, 
class FVC, 
class Traits>
 
   57    using Problem = 
typename Traits::Problem;
 
   59    using InvasionState = 
typename Traits::InvasionState;
 
   63    using FluxVariablesCache = 
typename Traits::FluxVariablesCache;
 
   66    static constexpr bool cachingEnabled = 
true;
 
   69    using LocalView = 
typename Traits::template LocalView<ThisType, cachingEnabled>;
 
   72    : problemPtr_(&problem)
 
   73    , invasionState_(problem) {}
 
   75    template<
class Gr
idGeometry, 
class Gr
idVolumeVariables, 
class SolutionVector>
 
   77                const GridVolumeVariables& gridVolVars,
 
   78                const SolutionVector& sol,
 
   79                bool forceUpdate = 
true)
 
   81        fluxVarsCache_.resize(gridGeometry.gridView().size(0));
 
   82        auto fvGeometry = 
localView(gridGeometry);
 
   83        auto elemVolVars = 
localView(gridVolVars);
 
   84        for (
const auto& element : elements(gridGeometry.gridView()))
 
   86            auto eIdx = gridGeometry.elementMapper().index(element);
 
   89            fvGeometry.bind(element);
 
   90            elemVolVars.bind(element, fvGeometry, sol);
 
   92            for (
auto&& scvf : scvfs(fvGeometry))
 
   93                cache(eIdx, scvf.index()).update(problem(), element, fvGeometry, elemVolVars, scvf, invasionState().invaded(element));
 
   97    template<
class FVElementGeometry, 
class ElementVolumeVariables>
 
   98    void updateElement(
const typename FVElementGeometry::Element& element,
 
   99                       const FVElementGeometry& fvGeometry,
 
  100                       const ElementVolumeVariables& elemVolVars)
 
  102        if constexpr (FluxVariablesCache::isSolDependent)
 
  104            const auto eIdx = fvGeometry.gridGeometry().elementMapper().index(element);
 
  105            for (
const auto& scvf : scvfs(fvGeometry))
 
  106                cache(eIdx, scvf.index()).update(problem(), element, fvGeometry, elemVolVars, scvf, invasionState().invaded(element));
 
  110    const Problem& problem()
 const 
  111    { 
return *problemPtr_; }
 
  114    const FluxVariablesCache& cache(std::size_t eIdx, [[maybe_unused]] std::size_t scvfIdx)
 const 
  115    { 
return fluxVarsCache_[eIdx]; }
 
  118    FluxVariablesCache& cache(std::size_t eIdx, [[maybe_unused]] std::size_t scvfIdx)
 
  119    { 
return fluxVarsCache_[eIdx]; }
 
  121    const InvasionState& invasionState()
 const 
  122    { 
return invasionState_; }
 
  124    InvasionState& invasionState()
 
  125    { 
return invasionState_; }
 
  128    const Problem* problemPtr_;
 
  129    std::vector<FluxVariablesCache> fluxVarsCache_;
 
  130    InvasionState invasionState_;
 
  138template<
class P, 
class FVC, 
class Traits>
 
  141    using Problem = 
typename Traits::Problem;
 
  142    using ThisType = PNMTwoPGridFluxVariablesCache<P, FVC, false, Traits>;
 
  143    using InvasionState = 
typename Traits::InvasionState;
 
  147    using FluxVariablesCache = 
typename Traits::FluxVariablesCache;
 
  150    static constexpr bool cachingEnabled = 
false;
 
  153    using LocalView = 
typename Traits::template LocalView<ThisType, cachingEnabled>;
 
  155    PNMTwoPGridFluxVariablesCache(
const Problem& problem)
 
  156    : problemPtr_(&problem)
 
  157    , invasionState_(problem) {}
 
  159    template<
class Gr
idGeometry, 
class Gr
idVolumeVariables, 
class SolutionVector>
 
  160    void update(
const GridGeometry& gridGeometry,
 
  161                const GridVolumeVariables& gridVolVars,
 
  162                const SolutionVector& sol,
 
  163                bool forceUpdate = 
true) {}
 
  165    const Problem& problem()
 const 
  166    { 
return *problemPtr_; }
 
  168    const InvasionState& invasionState()
 const 
  169    { 
return invasionState_; }
 
  171    InvasionState& invasionState()
 
  172    { 
return invasionState_; }
 
  175    const Problem* problemPtr_;
 
  176    InvasionState invasionState_;
 
Base class for the finite volume geometry for porenetwork models.
Definition discretization/porenetwork/gridgeometry.hh:477
The flux variables caches for an element.
Definition porenetwork/2p/elementfluxvariablescache.hh:28
Flux variable caches on a gridview.
Definition porenetwork/2p/gridfluxvariablescache.hh:47
Global flux variable cache.
GridCache::LocalView localView(const GridCache &gridCache)
Free function to get the local view of a grid cache object.
Definition localview.hh:26
Invasion state class for the two-phase PNM.
Free function to get the local view of a grid cache object.
Definition discretization/porenetwork/fvelementgeometry.hh:24
The infrastructure to retrieve run-time parameters from Dune::ParameterTrees.
Element flux variable cache.
Flux variable caches traits.
Definition porenetwork/2p/gridfluxvariablescache.hh:29
FVC FluxVariablesCache
Definition porenetwork/2p/gridfluxvariablescache.hh:31
IS InvasionState
Definition porenetwork/2p/gridfluxvariablescache.hh:32
PNMTwoPElementFluxVariablesCache< GridFluxVariablesCache, cachingEnabled > LocalView
Definition porenetwork/2p/gridfluxvariablescache.hh:35
P Problem
Definition porenetwork/2p/gridfluxvariablescache.hh:30