12#ifndef DUMUX_DISCRETIZATION_STAGGERED_GRID_VOLUMEVARIABLES_HH 
   13#define DUMUX_DISCRETIZATION_STAGGERED_GRID_VOLUMEVARIABLES_HH 
   15#include <dune/common/exceptions.hh> 
   16#include <dune/common/rangeutilities.hh> 
   25template<
class P, 
class VV>
 
   32    template<
class Gr
idVolumeVariables, 
bool cachingEnabled>
 
   37    template<
class Problem, 
class SolutionVector, 
class Element, 
class SubControlVolumeFace>
 
   39                                               const SolutionVector& sol,
 
   40                                               const Element& element,
 
   41                                               const SubControlVolumeFace& scvf)
 
   43        using CellCenterPrimaryVariables = 
typename SolutionVector::value_type;
 
   44        using Indices = 
typename VolumeVariables::Indices;
 
   45        static constexpr auto dim = PrimaryVariables::dimension - CellCenterPrimaryVariables::dimension;
 
   46        static constexpr auto offset = dim;
 
   48        const auto bcTypes = problem.boundaryTypes(element, scvf);
 
   52        for(
int i = 0; i < dim; ++i)
 
   54            if(bcTypes.isOutflow(Indices::velocity(i)))
 
   55                DUNE_THROW(Dune::InvalidStateException, 
"Outflow condition cannot be used for velocity. Set only a Dirichlet value for pressure instead.");
 
   58        if(bcTypes.isOutflow(Indices::pressureIdx))
 
   59            DUNE_THROW(Dune::InvalidStateException, 
"Outflow condition cannot be used for pressure. Set only a Dirichlet value for velocity instead.");
 
   63        if(bcTypes.isDirichlet(Indices::velocity(scvf.directionIndex())))
 
   65            if(bcTypes.isDirichlet(Indices::pressureIdx))
 
   66                DUNE_THROW(Dune::InvalidStateException, 
"A Dirichlet condition for velocity must not be combined with a Dirichlet condition for pressure");
 
   68                boundaryPriVars[Indices::pressureIdx] = sol[scvf.insideScvIdx()][Indices::pressureIdx - offset];
 
   74        if(bcTypes.isDirichlet(Indices::pressureIdx))
 
   76            if(bcTypes.isDirichlet(Indices::velocity(scvf.directionIndex())))
 
   77                DUNE_THROW(Dune::InvalidStateException, 
"A Dirichlet condition for velocity must not be combined with a Dirichlet condition for pressure");
 
   79                boundaryPriVars[Indices::pressureIdx] = problem.dirichlet(element, scvf)[Indices::pressureIdx];
 
   83        if(CellCenterPrimaryVariables::dimension == 1)
 
   84            return boundaryPriVars;
 
   87        for(
int eqIdx = offset; eqIdx < PrimaryVariables::dimension; ++eqIdx)
 
   89            if(eqIdx == Indices::pressureIdx)
 
   92            if(bcTypes.isDirichlet(eqIdx))
 
   93                boundaryPriVars[eqIdx] = problem.dirichlet(element, scvf)[eqIdx];
 
   94            else if(bcTypes.isOutflow(eqIdx) || bcTypes.isSymmetry() || bcTypes.isNeumann(eqIdx))
 
   95                boundaryPriVars[eqIdx] = sol[scvf.insideScvIdx()][eqIdx - offset];
 
   99        std::array<bool, VolumeVariables::numFluidComponents() - 1> isComponentOutflow;
 
  100        for(
int compIdx = 1; compIdx < VolumeVariables::numFluidComponents(); ++compIdx)
 
  102            const auto eqIdx = VolumeVariables::Indices::conti0EqIdx + compIdx;
 
  103            isComponentOutflow[compIdx -1] = bcTypes.isOutflow(eqIdx);
 
  106        if(Dune::any_true(isComponentOutflow) && !Dune::all_true(isComponentOutflow))
 
  107            DUNE_THROW(Dune::InvalidStateException, 
"Outflow condition must be set for all components!");
 
  109        return boundaryPriVars;
 
 
 
  117template<
class Traits, 
bool cachingEnabled>
 
  125template<
class Traits>
 
  129    using PrimaryVariables = 
typename Traits::VolumeVariables::PrimaryVariables;
 
  133    using Problem = 
typename Traits::Problem;
 
  136    using Indices = 
typename Traits::VolumeVariables::Indices;
 
  139    using VolumeVariables = 
typename Traits::VolumeVariables;
 
  142    static constexpr bool cachingEnabled = 
true;
 
  145    using LocalView = 
typename Traits::template LocalView<ThisType, cachingEnabled>;
 
  150    template<
class Gr
idGeometry, 
class SolutionVector>
 
  151    void update(
const GridGeometry& gridGeometry, 
const SolutionVector& sol)
 
  153        if (sol.size() != gridGeometry.numScv())
 
  154            DUNE_THROW(Dune::InvalidStateException, 
"The solution vector passed to the GridVolumeVariables has the wrong size.\n" 
  155                                                     << 
"Make sure to initialize the gridVariables correctly: \n\n" 
  156                                                     << 
"auto ffSol = partial(sol, ffFaceIdx, ffCellCenterIdx); \n" 
  157                                                     << 
"ffGridVariables->init(ffSol);\n\n");
 
  159        volumeVariables_.resize(gridGeometry.numScv());
 
  160        auto fvGeometry = 
localView(gridGeometry);
 
  161        for (
const auto& element : elements(gridGeometry.gridView()))
 
  163            fvGeometry.bindElement(element);
 
  164            for (
auto&& scv : scvs(fvGeometry))
 
  167                const auto& cellCenterPriVars = sol[scv.dofIndex()];
 
  171                volumeVariables_[scv.dofIndex()].update(elemSol, problem(), element, scv);
 
  176    const VolumeVariables& volVars(
const std::size_t scvIdx)
 const 
  177    { 
return volumeVariables_[scvIdx]; }
 
  179    VolumeVariables& volVars(
const std::size_t scvIdx)
 
  180    { 
return volumeVariables_[scvIdx]; }
 
  182    template<class SubControlVolume, typename std::enable_if_t<!std::is_integral<SubControlVolume>::value, 
int> = 0>
 
  183    const VolumeVariables& volVars(
const SubControlVolume& scv)
 const 
  184    { 
return volumeVariables_[scv.dofIndex()]; }
 
  186    template<class SubControlVolume, typename std::enable_if_t<!std::is_integral<SubControlVolume>::value, 
int> = 0>
 
  187    VolumeVariables& volVars(
const SubControlVolume& scv)
 
  188    { 
return volumeVariables_[scv.dofIndex()]; }
 
  190    const Problem& problem()
 const 
  191    { 
return *problemPtr_; }
 
  195    template<
class... Args>
 
  196    PrimaryVariables getBoundaryPriVars(Args&&... args)
 const 
  198        return Traits::getBoundaryPriVars(std::forward<Args>(args)...);
 
  202    const Problem* problemPtr_;
 
  203    std::vector<VolumeVariables> volumeVariables_;
 
  212template<
class Traits>
 
  215    using ThisType = StaggeredGridVolumeVariables<Traits, false>;
 
  216    using PrimaryVariables = 
typename Traits::VolumeVariables::PrimaryVariables;
 
  220    using Problem = 
typename Traits::Problem;
 
  223    using VolumeVariables = 
typename Traits::VolumeVariables;
 
  226    static constexpr bool cachingEnabled = 
false;
 
  229    using LocalView = 
typename Traits::template LocalView<ThisType, cachingEnabled>;
 
  231    StaggeredGridVolumeVariables(
const Problem& problem) : problemPtr_(&problem) {}
 
  233    template<
class Gr
idGeometry, 
class SolutionVector>
 
  234    void update(
const GridGeometry& gridGeometry, 
const SolutionVector& sol)
 
  236        if (sol.size() != gridGeometry.numScv())
 
  237            DUNE_THROW(Dune::InvalidStateException, 
"The solution vector passed to the GridVolumeVariables has the wrong size.\n" 
  238                                                     << 
"Make sure to initialize the gridVariables correctly: \n\n" 
  239                                                     << 
"auto ffSol = partial(sol, ffFaceIdx, ffCellCenterIdx); \n" 
  240                                                     << 
"ffGridVariables->init(ffSol);\n\n");
 
  243    const Problem& problem()
 const 
  244    { 
return *problemPtr_;}
 
  248    template<
class... Args>
 
  249    PrimaryVariables getBoundaryPriVars(Args&&... args)
 const 
  251        return Traits::getBoundaryPriVars(std::forward<Args>(args)...);
 
  256    const Problem* problemPtr_;
 
Base class for the element volume variables vector for the staggered model.
Definition staggered/freeflow/elementvolumevariables.hh:32
Grid volume variables class for staggered models.
Definition staggered/freeflow/gridvolumevariables.hh:118
GridCache::LocalView localView(const GridCache &gridCache)
Free function to get the local view of a grid cache object.
Definition localview.hh:26
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
PrimaryVariables makePriVarsFromCellCenterPriVars(const CellCenterPrimaryVariables &cellCenterPriVars)
Helper function to create a PrimaryVariables object from CellCenterPrimaryVariables.
Definition staggered/elementsolution.hh:29
Free function to get the local view of a grid cache object.
The local element solution class for staggered methods.
Base class for the element volume variables vector for the staggered model.
Definition staggered/freeflow/gridvolumevariables.hh:27
VV VolumeVariables
Definition staggered/freeflow/gridvolumevariables.hh:29
StaggeredElementVolumeVariables< GridVolumeVariables, cachingEnabled > LocalView
Definition staggered/freeflow/gridvolumevariables.hh:33
P Problem
Definition staggered/freeflow/gridvolumevariables.hh:28
static PrimaryVariables getBoundaryPriVars(const Problem &problem, const SolutionVector &sol, const Element &element, const SubControlVolumeFace &scvf)
Definition staggered/freeflow/gridvolumevariables.hh:38
typename VV::PrimaryVariables PrimaryVariables
Definition staggered/freeflow/gridvolumevariables.hh:30