13#ifndef DUMUX_CC_LOCAL_RESIDUAL_HH 
   14#define DUMUX_CC_LOCAL_RESIDUAL_HH 
   29template<
class TypeTag>
 
   40    using FVElementGeometry = 
typename GridGeometry::LocalView;
 
   42    using SubControlVolumeFace = 
typename FVElementGeometry::SubControlVolumeFace;
 
   46    using ParentType::ParentType;
 
   51                  const Element& element,
 
   52                  const FVElementGeometry& fvGeometry,
 
   53                  const ElementVolumeVariables& elemVolVars,
 
   54                  const ElementBoundaryTypes& elemBcTypes,
 
   55                  const ElementFluxVariablesCache& elemFluxVarsCache,
 
   56                  const SubControlVolumeFace& scvf)
 const 
   58        const auto& scv = fvGeometry.scv(scvf.insideScvIdx());
 
   59        const auto localScvIdx = scv.localDofIndex();
 
   60        residual[localScvIdx] += this->
asImp().evalFlux(problem, element, fvGeometry, elemVolVars, elemFluxVarsCache, scvf);
 
 
   65                         const Element& element,
 
   66                         const FVElementGeometry& fvGeometry,
 
   67                         const ElementVolumeVariables& elemVolVars,
 
   68                         const ElementFluxVariablesCache& elemFluxVarsCache,
 
   69                         const SubControlVolumeFace& scvf)
 const 
   71        NumEqVector flux(0.0);
 
   76            flux += this->
asImp().computeFlux(problem, element, fvGeometry, elemVolVars, scvf, elemFluxVarsCache);
 
   82            const auto& bcTypes = 
problem.boundaryTypes(element, scvf);
 
   85            if (bcTypes.hasDirichlet() && !bcTypes.hasNeumann())
 
   86                flux += this->
asImp().computeFlux(problem, element, fvGeometry, elemVolVars, scvf, elemFluxVarsCache);
 
   89            else if (bcTypes.hasNeumann() && !bcTypes.hasDirichlet())
 
   91                auto neumannFluxes = 
problem.neumann(element, fvGeometry, elemVolVars, elemFluxVarsCache, scvf);
 
   94                const auto& scv = fvGeometry.scv(scvf.insideScvIdx());
 
   95                neumannFluxes *= Extrusion::area(fvGeometry, scvf)*elemVolVars[scv].extrusionFactor();
 
   97                flux += neumannFluxes;
 
  101                DUNE_THROW(Dune::NotImplemented, 
"Mixed boundary conditions. Use pure boundary conditions by converting Dirichlet BCs to Robin BCs");
 
 
 
Calculates the element-wise residual for the cell-centered discretization schemes.
Definition cclocalresidual.hh:31
NumEqVector evalFlux(const Problem &problem, const Element &element, const FVElementGeometry &fvGeometry, const ElementVolumeVariables &elemVolVars, const ElementFluxVariablesCache &elemFluxVarsCache, const SubControlVolumeFace &scvf) const
evaluate the flux residual for a sub control volume face
Definition cclocalresidual.hh:64
void evalFlux(ElementResidualVector &residual, const Problem &problem, const Element &element, const FVElementGeometry &fvGeometry, const ElementVolumeVariables &elemVolVars, const ElementBoundaryTypes &elemBcTypes, const ElementFluxVariablesCache &elemFluxVarsCache, const SubControlVolumeFace &scvf) const
evaluate the flux residual for a sub control volume face and add to residual
Definition cclocalresidual.hh:49
typename ParentType::ElementResidualVector ElementResidualVector
Definition cclocalresidual.hh:45
Implementation & asImp()
Definition fvlocalresidual.hh:488
const Problem & problem() const
the problem
Definition fvlocalresidual.hh:473
ReservedBlockVector< NumEqVector, FVElementGeometry::maxNumElementScvs > ElementResidualVector
the container storing all element residuals
Definition fvlocalresidual.hh:57
FVLocalResidual(const Problem *problem, const TimeLoop *timeLoop=nullptr)
the constructor
Definition fvlocalresidual.hh:60
Defines all properties used in Dumux.
Helper classes to compute the integration elements.
The element-wise residual for finite volume schemes.
typename NumEqVectorTraits< PrimaryVariables >::type NumEqVector
A vector with the same size as numbers of equations This is the default implementation and has to be ...
Definition numeqvector.hh:34
typename GetProp< TypeTag, Property >::type GetPropType
get the type alias defined in the property
Definition propertysystem.hh:296
typename Extrusion< T >::type Extrusion_t
Convenience alias for obtaining the extrusion type.
Definition extrusion.hh:166
A helper to deduce a vector with the same size as numbers of equations.
A arithmetic block vector type based on DUNE's reserved vector.