12#ifndef DUMUX_DISCRETIZATION_STAGGERED_GRID_FLUXVARSCACHE_HH 
   13#define DUMUX_DISCRETIZATION_STAGGERED_GRID_FLUXVARSCACHE_HH 
   29template<
class P, 
class FVC, 
class FVCF, 
int upwOrder>
 
   36    template<
class Gr
idFluxVariablesCache, 
bool cachingEnabled>
 
 
   45template<
class Problem,
 
   46         class FluxVariablesCache,
 
   47         class FluxVariablesCacheFiller,
 
   48         bool EnableGridFluxVariablesCache = 
false,
 
   49         int upwindSchemeOrder = 1,
 
   58template<
class P, 
class FVC, 
class FVCF, 
int upwindSchemeOrder, 
class TheTraits>
 
   61    using Problem = 
typename TheTraits::Problem;
 
   65    using FluxVariablesCacheFiller = 
typename TheTraits::FluxVariablesCacheFiller;
 
   68    using Traits = TheTraits;
 
   71    using FluxVariablesCache = 
typename Traits::FluxVariablesCache;
 
   72    using Scalar = 
typename FluxVariablesCache::Scalar;
 
   74    static constexpr bool useHigherOrder = upwindSchemeOrder > 1;
 
   80    static constexpr bool cachingEnabled = 
true;
 
   83    using LocalView = 
typename Traits::template LocalView<ThisType, cachingEnabled>;
 
   86    : problemPtr_(&problem)
 
   87    , staggeredUpwindMethods_(problem.paramGroup())
 
   91    template<
class Gr
idGeometry, 
class Gr
idVolumeVariables, 
class SolutionVector>
 
   92    void update(
const GridGeometry& gridGeometry,
 
   93                const GridVolumeVariables& gridVolVars,
 
   94                const SolutionVector& sol,
 
   95                bool forceUpdate = 
false)
 
   98        if (FluxVariablesCacheFiller::isSolDependent || forceUpdate)
 
  102            FluxVariablesCacheFiller filler(problem());
 
  104            fluxVarsCache_.resize(gridGeometry.numScvf());
 
  105            auto fvGeometry = 
localView(gridGeometry);
 
  106            auto elemVolVars = 
localView(gridVolVars);
 
  107            for (
const auto& element : elements(gridGeometry.gridView()))
 
  110                fvGeometry.bind(element);
 
  111                elemVolVars.bind(element, fvGeometry, sol);
 
  113                for (
auto&& scvf : scvfs(fvGeometry))
 
  115                    filler.fill(*
this, fluxVarsCache_[scvf.index()], element, fvGeometry, elemVolVars, scvf, forceUpdate);
 
  124        return staggeredUpwindMethods_;
 
  127    const Problem& problem()
 const 
  128    { 
return *problemPtr_; }
 
  131    const FluxVariablesCache& operator [](std::size_t scvfIdx)
 const 
  132    { 
return fluxVarsCache_[scvfIdx]; }
 
  134    FluxVariablesCache& operator [](std::size_t scvfIdx)
 
  135    { 
return fluxVarsCache_[scvfIdx]; }
 
  138    const Problem* problemPtr_;
 
  141    std::vector<FluxVariablesCache> fluxVarsCache_;
 
  142    std::vector<std::size_t> globalScvfIndices_;
 
  150template<
class P, 
class FVC, 
class FVCF, 
int upwindSchemeOrder, 
class TheTraits>
 
  153    using Problem = 
typename TheTraits::Problem;
 
  154    using ThisType = StaggeredGridFluxVariablesCache<P, FVC, FVCF, false, upwindSchemeOrder, TheTraits>;
 
  157    using FluxVariablesCacheFiller = 
typename TheTraits::FluxVariablesCacheFiller;
 
  160    using Traits = TheTraits;
 
  163    using FluxVariablesCache = 
typename Traits::FluxVariablesCache;
 
  166    using Scalar = 
typename FluxVariablesCache::Scalar;
 
  169    static constexpr bool cachingEnabled = 
false;
 
  172    using LocalView = 
typename Traits::template LocalView<ThisType, cachingEnabled>;
 
  175    using UpwindScheme = StaggeredUpwindMethods<Scalar, upwindSchemeOrder>;
 
  177    StaggeredGridFluxVariablesCache(
const Problem& problem)
 
  178    : problemPtr_(&problem)
 
  179    , staggeredUpwindMethods_(problem.paramGroup())
 
  183    template<
class Gr
idGeometry, 
class Gr
idVolumeVariables, 
class SolutionVector>
 
  184    void update(
const GridGeometry& gridGeometry,
 
  185                const GridVolumeVariables& gridVolVars,
 
  186                const SolutionVector& sol,
 
  187                bool forceUpdate = 
false) {}
 
  189    const Problem& problem()
 const 
  190    { 
return *problemPtr_; }
 
  195        return staggeredUpwindMethods_;
 
  199    const Problem* problemPtr_;
 
Base class for the stencil local flux variables cache for the staggered model.
Definition discretization/staggered/elementfluxvariablescache.hh:30
Flux variables cache class for staggered models.
Definition discretization/staggered/gridfluxvariablescache.hh:51
This file contains different higher order methods for approximating the velocity.
Definition staggeredupwindmethods.hh:50
Base class for the stencil local flux variables cache for the staggered model.
GridCache::LocalView localView(const GridCache &gridCache)
Free function to get the local view of a grid cache object.
Definition localview.hh:26
UpwindSchemeImpl< GridGeometry, typename GridGeometry::DiscretizationMethod > UpwindScheme
The upwind scheme used for the advective fluxes. This depends on the chosen discretization method.
Definition flux/upwindscheme.hh:30
Free function to get the local view of a grid cache object.
This file contains different higher order methods for approximating the velocity.
Traits class to be used for the StaggeredGridVFluxVariablesCache.
Definition discretization/staggered/gridfluxvariablescache.hh:31
StaggeredElementFluxVariablesCache< GridFluxVariablesCache, cachingEnabled > LocalView
Definition discretization/staggered/gridfluxvariablescache.hh:37
static constexpr int upwindSchemeOrder
Definition discretization/staggered/gridfluxvariablescache.hh:38
P Problem
Definition discretization/staggered/gridfluxvariablescache.hh:32
FVCF FluxVariablesCacheFiller
Definition discretization/staggered/gridfluxvariablescache.hh:34
FVC FluxVariablesCache
Definition discretization/staggered/gridfluxvariablescache.hh:33