13#ifndef DUMUX_CVFE_LOCAL_RESIDUAL_HH 
   14#define DUMUX_CVFE_LOCAL_RESIDUAL_HH 
   16#include <dune/common/std/type_traits.hh> 
   17#include <dune/geometry/type.hh> 
   18#include <dune/istl/matrix.hh> 
   27template<
class Imp, 
class P, 
class G, 
class S, 
class V>
 
   29    std::declval<Imp>().computeStorage(
 
   30        std::declval<P>(), std::declval<G>(), std::declval<S>(), std::declval<V>(), 
true 
   34template<
class Imp, 
class P, 
class G, 
class S, 
class V>
 
   36{ 
return Dune::Std::is_detected<TimeInfoInterfaceCVFEDetector, Imp, P, G, S, V>::value; }
 
 
   40    std::declval<Imp>().isOverlapping()
 
   45{ 
return Dune::Std::is_detected<SCVFIsOverlappingDetector, Imp>::value; }
 
 
 
   58template<
class TypeTag>
 
   66    using GridView = 
typename GridGeometry::GridView;
 
   67    using Element = 
typename GridView::template Codim<0>::Entity;
 
   69    using FVElementGeometry = 
typename GridGeometry::LocalView;
 
   71    using ElementVolumeVariables = 
typename GridVolumeVariables::LocalView;
 
   72    using VolumeVariables = 
typename GridVolumeVariables::VolumeVariables;
 
   73    using SubControlVolume = 
typename FVElementGeometry::SubControlVolume;
 
   74    using SubControlVolumeFace = 
typename FVElementGeometry::SubControlVolumeFace;
 
   82    using ParentType::ParentType;
 
   87                  const Element& element,
 
   88                  const FVElementGeometry& fvGeometry,
 
   89                  const ElementVolumeVariables& elemVolVars,
 
   90                  const ElementBoundaryTypes& elemBcTypes,
 
   91                  const ElementFluxVariablesCache& elemFluxVarsCache,
 
   92                  const SubControlVolumeFace& scvf)
 const 
   94        const auto flux = 
evalFlux(
problem, element, fvGeometry, elemVolVars, elemBcTypes, elemFluxVarsCache, scvf);
 
   97            const auto& insideScv = fvGeometry.scv(scvf.insideScvIdx());
 
   98            const auto& outsideScv = fvGeometry.scv(scvf.outsideScvIdx());
 
   99            residual[insideScv.localDofIndex()] += flux;
 
  104                if (!scvf.isOverlapping())
 
  105                    residual[outsideScv.localDofIndex()] -= flux;
 
  108                residual[outsideScv.localDofIndex()] -= flux;
 
  112            const auto& insideScv = fvGeometry.scv(scvf.insideScvIdx());
 
  113            residual[insideScv.localDofIndex()] += flux;
 
 
  119                         const Element& element,
 
  120                         const FVElementGeometry& fvGeometry,
 
  121                         const ElementVolumeVariables& elemVolVars,
 
  122                         const ElementBoundaryTypes& elemBcTypes,
 
  123                         const ElementFluxVariablesCache& elemFluxVarsCache,
 
  124                         const SubControlVolumeFace& scvf)
 const 
  126        NumEqVector flux(0.0);
 
  129        if (!scvf.boundary())
 
  130            flux += this->
asImp().computeFlux(problem, element, fvGeometry, elemVolVars, scvf, elemFluxVarsCache);
 
  135            const auto& scv = fvGeometry.scv(scvf.insideScvIdx());
 
  136            const auto& bcTypes = elemBcTypes.get(fvGeometry, scv);
 
  141            if (bcTypes.hasNeumann())
 
  143                auto neumannFluxes = 
problem.neumann(element, fvGeometry, elemVolVars, elemFluxVarsCache, scvf);
 
  146                neumannFluxes *= Extrusion::area(fvGeometry, scvf)*elemVolVars[scv].extrusionFactor();
 
  149                for (
int eqIdx = 0; eqIdx < NumEqVector::dimension; ++eqIdx)
 
  150                    if (bcTypes.isNeumann(eqIdx))
 
  151                        flux[eqIdx] += neumannFluxes[eqIdx];
 
 
  176                     const Element& element,
 
  177                     const FVElementGeometry& fvGeometry,
 
  178                     const ElementVolumeVariables& prevElemVolVars,
 
  179                     const ElementVolumeVariables& curElemVolVars,
 
  180                     const SubControlVolume& scv)
 const 
  182        const auto& curVolVars = curElemVolVars[scv];
 
  183        const auto& prevVolVars = prevElemVolVars[scv];
 
  188        NumEqVector prevStorage = computeStorageImpl_(
problem, fvGeometry, scv, prevVolVars, 
true);
 
  189        NumEqVector storage = computeStorageImpl_(
problem, fvGeometry, scv, curVolVars, 
false);
 
  191        prevStorage *= prevVolVars.extrusionFactor();
 
  192        storage *= curVolVars.extrusionFactor();
 
  194        storage -= prevStorage;
 
  195        storage *= Extrusion::volume(fvGeometry, scv);
 
  198        residual[scv.localDofIndex()] += storage;
 
 
  203                                    const FVElementGeometry& fvGeometry,
 
  204                                    const SubControlVolume& scv,
 
  205                                    const VolumeVariables& volVars,
 
  206                                    [[maybe_unused]] 
bool isPreviousTimeStep)
 const 
  209            return this->
asImp().computeStorage(problem, fvGeometry, scv, volVars, isPreviousTimeStep);
 
  211            return this->
asImp().computeStorage(problem, scv, volVars);
 
 
The element-wise residual for control-volume finite element schemes.
Definition cvfelocalresidual.hh:60
NumEqVector evalFlux(const Problem &problem, const Element &element, const FVElementGeometry &fvGeometry, const ElementVolumeVariables &elemVolVars, const ElementBoundaryTypes &elemBcTypes, const ElementFluxVariablesCache &elemFluxVarsCache, const SubControlVolumeFace &scvf) const
evaluate flux residuals for one sub control volume face
Definition cvfelocalresidual.hh:118
void evalStorage(ElementResidualVector &residual, const Problem &problem, const Element &element, const FVElementGeometry &fvGeometry, const ElementVolumeVariables &prevElemVolVars, const ElementVolumeVariables &curElemVolVars, const SubControlVolume &scv) const
Compute the storage local residual, i.e. the deviation of the storage term from zero for instationary...
Definition cvfelocalresidual.hh:174
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 flux residuals for one sub control volume face and add to residual
Definition cvfelocalresidual.hh:85
typename ParentType::ElementResidualVector ElementResidualVector
Definition cvfelocalresidual.hh:81
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
const TimeLoop & timeLoop() const
Definition fvlocalresidual.hh:478
ElementResidualVector evalStorage(const Problem &problem, const Element &element, const GridGeometry &gridGeometry, const GridVariables &gridVariables, const SolutionVector &sol) const
Compute the storage term for the current solution.
Definition fvlocalresidual.hh:86
FVLocalResidual(const Problem *problem, const TimeLoop *timeLoop=nullptr)
the constructor
Definition fvlocalresidual.hh:60
virtual Scalar timeStepSize() const =0
Returns the suggested time step length .
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
Distance implementation details.
Definition cvfelocalresidual.hh:25
constexpr bool hasScvfIsOverlapping()
Definition cvfelocalresidual.hh:44
constexpr bool hasTimeInfoInterfaceCVFE()
Definition cvfelocalresidual.hh:35
decltype( std::declval< Imp >().computeStorage( std::declval< P >(), std::declval< G >(), std::declval< S >(), std::declval< V >(), true)) TimeInfoInterfaceCVFEDetector
Definition cvfelocalresidual.hh:28
decltype( std::declval< Imp >().isOverlapping()) SCVFIsOverlappingDetector
Definition cvfelocalresidual.hh:39
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.