12#ifndef DUMUX_FV_LOCAL_ASSEMBLER_BASE_HH 
   13#define DUMUX_FV_LOCAL_ASSEMBLER_BASE_HH 
   15#include <dune/common/reservedvector.hh> 
   16#include <dune/grid/common/gridenums.hh>  
   17#include <dune/istl/matrixindexset.hh> 
   34template<
class TypeTag, 
class Assembler, 
class Implementation, 
bool useImplicitAssembly>
 
   42    using SolutionVector = 
typename Assembler::SolutionVector;
 
   45    using SubControlVolume = 
typename FVElementGeometry::SubControlVolume;
 
   46    using SubControlVolumeFace = 
typename FVElementGeometry::SubControlVolumeFace;
 
   51    using Element = 
typename GridView::template Codim<0>::Entity;
 
   55    using LocalResidual = std::decay_t<decltype(std::declval<Assembler>().localResidual())>;
 
   63                                  const SolutionVector& 
curSol)
 
 
   80                                  const SolutionVector& 
curSol,
 
 
  102    { 
return useImplicitAssembly; }
 
 
  111            if (this->
assembler().isStationaryProblem())
 
  112                DUNE_THROW(Dune::InvalidStateException, 
"Using explicit jacobian assembler with stationary local residual");
 
 
  156        return localResidual_.evalFluxAndSource(element_, fvGeometry_, elemVolVars, elemFluxVarsCache_, elemBcTypes_);
 
 
  166        return localResidual_.evalStorage(element_, fvGeometry_, prevElemVolVars_, curElemVolVars_);
 
 
  178        const auto& prevSol = this->
assembler().prevSol();
 
  191            if (!this->
assembler().isStationaryProblem())
 
 
  205    template<
typename ApplyFunction, 
class P = Problem, 
typename std::enable_if_t<P::enableInternalDirichletConstra
ints(), 
int> = 0>
 
  210        for (
const auto& scvI : scvs(this->
fvGeometry()))
 
  212            const auto internalDirichletConstraints = 
asImp_().problem().hasInternalDirichletConstraint(this->
element(), scvI);
 
  213            if (internalDirichletConstraints.any())
 
  215                const auto dirichletValues = 
asImp_().problem().internalDirichlet(this->
element(), scvI);
 
  217                for (
int eqIdx = 0; eqIdx < internalDirichletConstraints.size(); ++eqIdx)
 
  218                    if (internalDirichletConstraints[eqIdx])
 
  219                        applyDirichlet(scvI, dirichletValues, eqIdx, eqIdx);
 
 
  224    template<
typename ApplyFunction, 
class P = Problem, 
typename std::enable_if_t<!P::enableInternalDirichletConstra
ints(), 
int> = 0>
 
  230    { 
return assembler_.problem(); }
 
 
  234    { 
return assembler_; }
 
 
  242    { 
return elementIsGhost_; }
 
 
  250    { 
return fvGeometry_; }
 
 
  254    { 
return curElemVolVars_; }
 
 
  258    { 
return prevElemVolVars_; }
 
 
  262    { 
return elemFluxVarsCache_; }
 
 
  266    { 
return localResidual_; }
 
 
  270    { 
return elemBcTypes_; }
 
 
  274    { 
return fvGeometry_; }
 
 
  278    { 
return curElemVolVars_; }
 
 
  282    { 
return prevElemVolVars_; }
 
 
  286    { 
return elemFluxVarsCache_; }
 
 
  290    { 
return elemBcTypes_; }
 
 
  294    { 
return localResidual_; }
 
 
  298    { 
return *
static_cast<Implementation*
>(
this); }
 
 
  301    { 
return *
static_cast<const Implementation*
>(
this); }
 
 
  303    template<class T = TypeTag, typename std::enable_if_t<!GetPropType<T, Properties::GridVariables>::GridVolumeVariables::cachingEnabled, 
int> = 0>
 
  304    VolumeVariables& 
getVolVarAccess(GridVolumeVariables& gridVolVars, ElementVolumeVariables& elemVolVars, 
const SubControlVolume& scv)
 
  305    { 
return elemVolVars[scv]; }
 
 
  307    template<class T = TypeTag, typename std::enable_if_t<GetPropType<T, Properties::GridVariables>::GridVolumeVariables::cachingEnabled, 
int> = 0>
 
  308    VolumeVariables& 
getVolVarAccess(GridVolumeVariables& gridVolVars, ElementVolumeVariables& elemVolVars, 
const SubControlVolume& scv)
 
  309    { 
return gridVolVars.volVars(scv); }
 
 
  313    const Assembler& assembler_; 
 
  314    const Element& element_; 
 
  315    const SolutionVector& curSol_; 
 
  317    FVElementGeometry fvGeometry_;
 
  318    ElementVolumeVariables curElemVolVars_;
 
  319    ElementVolumeVariables prevElemVolVars_;
 
  320    ElementFluxVariablesCache elemFluxVarsCache_;
 
  321    ElementBoundaryTypes elemBcTypes_;
 
  324    bool elementIsGhost_; 
 
 
const ElementVolumeVariables & prevElemVolVars() const
The element volume variables of the provious time step.
Definition assembly/fvlocalassemblerbase.hh:281
void bindLocalViews()
Convenience function bind and prepare all relevant variables required for the evaluation of the local...
Definition assembly/fvlocalassemblerbase.hh:173
ElementResidualVector evalLocalResidual(const ElementVolumeVariables &elemVolVars) const
Evaluates the complete local residual for the current element.
Definition assembly/fvlocalassemblerbase.hh:125
ElementVolumeVariables & curElemVolVars()
Definition assembly/fvlocalassemblerbase.hh:253
typename LocalResidual::ElementResidualVector ElementResidualVector
Definition assembly/fvlocalassemblerbase.hh:56
ElementBoundaryTypes & elemBcTypes()
The element's boundary types.
Definition assembly/fvlocalassemblerbase.hh:269
Implementation & asImp_()
Definition assembly/fvlocalassemblerbase.hh:297
ElementVolumeVariables & prevElemVolVars()
Definition assembly/fvlocalassemblerbase.hh:257
ElementResidualVector evalLocalResidual() const
Convenience function to evaluate the complete local residual for the current element....
Definition assembly/fvlocalassemblerbase.hh:108
const Problem & problem() const
The problem.
Definition assembly/fvlocalassemblerbase.hh:229
const ElementBoundaryTypes & elemBcTypes() const
The element's boundary types.
Definition assembly/fvlocalassemblerbase.hh:289
const ElementFluxVariablesCache & elemFluxVarsCache() const
The element flux variables cache.
Definition assembly/fvlocalassemblerbase.hh:285
const Implementation & asImp_() const
Definition assembly/fvlocalassemblerbase.hh:300
FVLocalAssemblerBase(const Assembler &assembler, const Element &element, const SolutionVector &curSol, const FVElementGeometry &fvGeometry, const ElementVolumeVariables &curElemVolVars, const ElementVolumeVariables &prevElemVolVars, const ElementFluxVariablesCache &elemFluxVarsCache, const LocalResidual &localResidual, const bool elementIsGhost)
The constructor. General version explicitly expecting each argument.
Definition assembly/fvlocalassemblerbase.hh:78
const ElementVolumeVariables & curElemVolVars() const
The current element volume variables.
Definition assembly/fvlocalassemblerbase.hh:277
FVLocalAssemblerBase(const Assembler &assembler, const Element &element, const SolutionVector &curSol)
The constructor. Delegates to the general constructor.
Definition assembly/fvlocalassemblerbase.hh:61
const LocalResidual & localResidual() const
The local residual for the current element.
Definition assembly/fvlocalassemblerbase.hh:293
FVElementGeometry & fvGeometry()
Definition assembly/fvlocalassemblerbase.hh:249
ElementResidualVector evalLocalFluxAndSourceResidual(const ElementVolumeVariables &elemVolVars) const
Evaluates the flux and source terms (i.e, the terms without a time derivative) of the local residual ...
Definition assembly/fvlocalassemblerbase.hh:154
const Assembler & assembler() const
Definition assembly/fvlocalassemblerbase.hh:233
ElementFluxVariablesCache & elemFluxVarsCache()
Definition assembly/fvlocalassemblerbase.hh:261
bool elementIsGhost() const
Definition assembly/fvlocalassemblerbase.hh:241
std::decay_t< decltype(std::declval< Assembler >().localResidual())> LocalResidual
Definition assembly/fvlocalassemblerbase.hh:55
LocalResidual & localResidual()
Definition assembly/fvlocalassemblerbase.hh:265
static constexpr bool isImplicit()
Returns true if the assembler considers implicit assembly.
Definition assembly/fvlocalassemblerbase.hh:101
ElementResidualVector evalLocalFluxAndSourceResidual() const
Convenience function to evaluate the flux and source terms (i.e, the terms without a time derivative)...
Definition assembly/fvlocalassemblerbase.hh:142
void enforceInternalDirichletConstraints(const ApplyFunction &applyDirichlet)
Enforces Dirichlet constraints if enabled in the problem.
Definition assembly/fvlocalassemblerbase.hh:206
const Element & element() const
Definition assembly/fvlocalassemblerbase.hh:237
VolumeVariables & getVolVarAccess(GridVolumeVariables &gridVolVars, ElementVolumeVariables &elemVolVars, const SubControlVolume &scv)
Definition assembly/fvlocalassemblerbase.hh:304
const FVElementGeometry & fvGeometry() const
The finite volume geometry.
Definition assembly/fvlocalassemblerbase.hh:273
ElementResidualVector evalLocalStorageResidual() const
Convenience function to evaluate storage term (i.e, the term with a time derivative) of the local res...
Definition assembly/fvlocalassemblerbase.hh:164
const SolutionVector & curSol() const
Definition assembly/fvlocalassemblerbase.hh:245
Defines all properties used in Dumux.
An enum class to define various differentiation methods available in order to compute the derivatives...
GridCache::LocalView localView(const GridCache &gridCache)
Free function to get the local view of a grid cache object.
Definition localview.hh:26
typename GetProp< TypeTag, Property >::type GetPropType
get the type alias defined in the property
Definition propertysystem.hh:296
Definition common/pdesolver.hh:24
The infrastructure to retrieve run-time parameters from Dune::ParameterTrees.
A arithmetic block vector type based on DUNE's reserved vector.