13#ifndef DUMUX_STAGGERED_LOCAL_RESIDUAL_HH 
   14#define DUMUX_STAGGERED_LOCAL_RESIDUAL_HH 
   27template<
class TypeTag>
 
   35    using Element = 
typename GridView::template Codim<0>::Entity;
 
   40    using FVElementGeometry = 
typename GridGeometry::LocalView;
 
   41    using SubControlVolume = 
typename GridGeometry::SubControlVolume;
 
   42    using SubControlVolumeFace = 
typename GridGeometry::SubControlVolumeFace;
 
   43    using FaceSubControlVolume = 
typename GridGeometry::Traits::FaceSubControlVolume;
 
   67                                                           const FVElementGeometry& fvGeometry,
 
   68                                                           const ElementVolumeVariables& elemVolVars,
 
   69                                                           const ElementFaceVariables& elemFaceVars,
 
   70                                                           const ElementBoundaryTypes& bcTypes,
 
   71                                                           const ElementFluxVariablesCache& elemFluxVarsCache)
 const 
   76        for (
auto&& scv : scvs(fvGeometry))
 
   77            asImp().evalSourceForCellCenter(residual, this->
problem(), element, fvGeometry, elemVolVars, elemFaceVars, scv);
 
   80        for (
auto&& scvf : scvfs(fvGeometry))
 
   81            asImp().evalFluxForCellCenter(residual, this->
problem(), element, fvGeometry, elemVolVars, elemFaceVars, bcTypes, elemFluxVarsCache, scvf);
 
 
   89                               const Element& element,
 
   90                               const FVElementGeometry& fvGeometry,
 
   91                               const ElementVolumeVariables& elemVolVars,
 
   92                               const ElementFaceVariables& elemFaceVars,
 
   93                               const ElementBoundaryTypes& elemBcTypes,
 
   94                               const ElementFluxVariablesCache& elemFluxVarsCache,
 
   95                               const SubControlVolumeFace& scvf)
 const 
   98            residual += 
asImp_().computeFluxForCellCenter(
problem, element, fvGeometry, elemVolVars, elemFaceVars, scvf, elemFluxVarsCache);
 
  100            residual += 
asImp_().computeBoundaryFluxForCellCenter(
problem, element, fvGeometry, scvf, elemVolVars, elemFaceVars, elemBcTypes, elemFluxVarsCache);
 
 
  106                                 const Element& element,
 
  107                                 const FVElementGeometry& fvGeometry,
 
  108                                 const ElementVolumeVariables& curElemVolVars,
 
  109                                 const ElementFaceVariables& curElemFaceVars,
 
  110                                 const SubControlVolume& scv)
 const 
  112            const auto curExtrusionFactor = curElemVolVars[scv].extrusionFactor();
 
  115            auto source = 
asImp_().computeSourceForCellCenter(
problem, element, fvGeometry, curElemVolVars, curElemFaceVars, scv);
 
  116            source *= Extrusion::volume(fvGeometry, scv)*curExtrusionFactor;
 
 
  122                                                     const FVElementGeometry& fvGeometry,
 
  123                                                     const ElementVolumeVariables& prevElemVolVars,
 
  124                                                     const ElementVolumeVariables& curElemVolVars)
 const 
  126        assert(timeLoop_ && 
"no time loop set for storage term evaluation");
 
  129        for (
auto&& scv : scvs(fvGeometry))
 
  130            asImp().evalStorageForCellCenter(storage, 
problem(), element, fvGeometry, prevElemVolVars, curElemVolVars, scv);
 
 
  138                                  const Element &element,
 
  139                                  const FVElementGeometry& fvGeometry,
 
  140                                  const ElementVolumeVariables& prevElemVolVars,
 
  141                                  const ElementVolumeVariables& curElemVolVars,
 
  142                                  const SubControlVolume& scv)
 const 
  145        const auto& curVolVars = curElemVolVars[scv];
 
  146        const auto& prevVolVars = prevElemVolVars[scv];
 
  154        auto prevCCStorage = 
asImp_().computeStorageForCellCenter(
problem, scv, prevVolVars);
 
  155        auto curCCStorage = 
asImp_().computeStorageForCellCenter(
problem, scv, curVolVars);
 
  157        prevCCStorage *= prevVolVars.extrusionFactor();
 
  158        curCCStorage *= curVolVars.extrusionFactor();
 
  160        storage = std::move(curCCStorage);
 
  161        storage -= std::move(prevCCStorage);
 
  162        storage *= Extrusion::volume(fvGeometry, scv);
 
  163        storage /= timeLoop_->timeStepSize();
 
 
  169    template<
class... Args>
 
  176    template<
class... Args>
 
  191                                               const FVElementGeometry& fvGeometry,
 
  192                                               const ElementVolumeVariables& elemVolVars,
 
  193                                               const ElementFaceVariables& elemFaceVars,
 
  194                                               const ElementBoundaryTypes& bcTypes,
 
  195                                               const ElementFluxVariablesCache& elemFluxVarsCache,
 
  196                                               const SubControlVolumeFace& scvf)
 const 
  199        asImp().evalSourceForFace(residual, this->
problem(), element, fvGeometry, elemVolVars, elemFaceVars, scvf);
 
  200        asImp().evalFluxForFace(residual, this->
problem(), element, fvGeometry, elemVolVars, elemFaceVars, bcTypes, elemFluxVarsCache, scvf);
 
 
  208                         const Element& element,
 
  209                         const FVElementGeometry& fvGeometry,
 
  210                         const ElementVolumeVariables& elemVolVars,
 
  211                         const ElementFaceVariables& elemFaceVars,
 
  212                         const ElementBoundaryTypes& elemBcTypes,
 
  213                         const ElementFluxVariablesCache& elemFluxVarsCache,
 
  214                         const SubControlVolumeFace& scvf)
 const 
  216        if (!scvf.boundary())
 
  217            residual += 
asImp_().computeFluxForFace(
problem, element, scvf, fvGeometry, elemVolVars, elemFaceVars, elemFluxVarsCache);
 
  219            residual += 
asImp_().computeBoundaryFluxForFace(
problem, element, fvGeometry, scvf, elemVolVars, elemFaceVars, elemBcTypes, elemFluxVarsCache);
 
 
  225                           const Element& element,
 
  226                           const FVElementGeometry& fvGeometry,
 
  227                           const ElementVolumeVariables& elemVolVars,
 
  228                           const ElementFaceVariables& elemFaceVars,
 
  229                           const SubControlVolumeFace& scvf)
 const 
  232        auto source = 
asImp_().computeSourceForFace(
problem, element, fvGeometry, scvf, elemVolVars, elemFaceVars);
 
  233        const auto& scv = fvGeometry.scv(scvf.insideScvIdx());
 
  234        const auto extrusionFactor = elemVolVars[scv].extrusionFactor();
 
  237        auto faceScvCenter = scvf.center() + scv.center();
 
  238        faceScvCenter *= 0.5;
 
  239        FaceSubControlVolume faceScv(faceScvCenter, 0.5*scv.volume());
 
  241        source *= Extrusion::volume(fvGeometry, faceScv)*extrusionFactor;
 
 
  247                                         const FVElementGeometry& fvGeometry,
 
  248                                         const ElementVolumeVariables& prevElemVolVars,
 
  249                                         const ElementVolumeVariables& curElemVolVars,
 
  250                                         const ElementFaceVariables& prevElemFaceVars,
 
  251                                         const ElementFaceVariables& curElemFaceVars,
 
  252                                         const SubControlVolumeFace& scvf)
 const 
  254        assert(timeLoop_ && 
"no time loop set for storage term evaluation");
 
  256        asImp().evalStorageForFace(storage, 
problem(), element, fvGeometry, prevElemVolVars, curElemVolVars, prevElemFaceVars, curElemFaceVars, scvf);
 
 
  263                            const Element& element,
 
  264                            const FVElementGeometry& fvGeometry,
 
  265                            const ElementVolumeVariables& prevElemVolVars,
 
  266                            const ElementVolumeVariables& curElemVolVars,
 
  267                            const ElementFaceVariables& prevElemFaceVars,
 
  268                            const ElementFaceVariables& curElemFaceVars,
 
  269                            const SubControlVolumeFace& scvf)
 const 
  271        const auto& scv = fvGeometry.scv(scvf.insideScvIdx());
 
  273        auto storage = 
asImp_().computeStorageForFace(
problem, scvf, curElemVolVars[scv], curElemFaceVars);
 
  274        storage -= 
asImp_().computeStorageForFace(
problem, scvf, prevElemVolVars[scv], prevElemFaceVars);
 
  276        const auto extrusionFactor = curElemVolVars[scv].extrusionFactor();
 
  279        auto faceScvCenter = scvf.center() + scv.center();
 
  280        faceScvCenter *= 0.5;
 
  281        FaceSubControlVolume faceScv(faceScvCenter, 0.5*scv.volume());
 
  283        storage *= Extrusion::volume(fvGeometry, faceScv)*extrusionFactor;
 
  284        storage /= timeLoop_->timeStepSize();
 
 
  291    { 
return !timeLoop_; }
 
 
  295    { 
return *problem_; }
 
 
  300    { 
return *
static_cast<Implementation*
>(
this); }
 
 
  303    { 
return *
static_cast<const Implementation*
>(
this); }
 
 
  307    { 
return *timeLoop_; }
 
 
  310    { 
return *timeLoop_; }
 
 
  313    { 
return *
static_cast<Implementation*
>(
this); }
 
 
  316    { 
return *
static_cast<const Implementation*
>(
this); }
 
 
  319    const Problem* problem_; 
 
 
CellCenterResidualValue evalStorageForCellCenter(const Element &element, const FVElementGeometry &fvGeometry, const ElementVolumeVariables &prevElemVolVars, const ElementVolumeVariables &curElemVolVars) const
Evaluate the storage terms for a cell center residual.
Definition staggeredlocalresidual.hh:121
GetPropType< TypeTag, Properties::CellCenterPrimaryVariables > CellCenterResidualValue
Definition staggeredlocalresidual.hh:54
Implementation & asImp()
Definition staggeredlocalresidual.hh:312
void evalStorageForCellCenter(CellCenterResidualValue &residual, const Problem &problem, const Element &element, const FVElementGeometry &fvGeometry, const ElementVolumeVariables &prevElemVolVars, const ElementVolumeVariables &curElemVolVars, const SubControlVolume &scv) const
Evaluate the storage terms for a cell center residual.
Definition staggeredlocalresidual.hh:136
Implementation & asImp_()
Definition staggeredlocalresidual.hh:299
const TimeLoop & timeLoop() const
Definition staggeredlocalresidual.hh:309
FaceResidualValue evalFluxAndSourceForFace(const Element &element, const FVElementGeometry &fvGeometry, const ElementVolumeVariables &elemVolVars, const ElementFaceVariables &elemFaceVars, const ElementBoundaryTypes &bcTypes, const ElementFluxVariablesCache &elemFluxVarsCache, const SubControlVolumeFace &scvf) const
Convenience function to evaluate the flux and source terms for the face residual.
Definition staggeredlocalresidual.hh:190
GetPropType< TypeTag, Properties::FacePrimaryVariables > FaceResidualValue
Definition staggeredlocalresidual.hh:55
const Problem & problem() const
the problem
Definition staggeredlocalresidual.hh:294
const Implementation & asImp_() const
Definition staggeredlocalresidual.hh:302
CellCenterResidualValue evalStorage(Args &&... args) const
for compatibility with FVLocalAssemblerBase
Definition staggeredlocalresidual.hh:177
void evalSourceForFace(FaceResidualValue &residual, const Problem &problem, const Element &element, const FVElementGeometry &fvGeometry, const ElementVolumeVariables &elemVolVars, const ElementFaceVariables &elemFaceVars, const SubControlVolumeFace &scvf) const
Evaluate the source terms for a face residual.
Definition staggeredlocalresidual.hh:223
bool isStationary() const
If no solution has been set, we treat the problem as stationary.
Definition staggeredlocalresidual.hh:290
TimeLoop & timeLoop()
Definition staggeredlocalresidual.hh:306
StaggeredLocalResidual(const Problem *problem, const TimeLoop *timeLoop=nullptr)
the constructor
Definition staggeredlocalresidual.hh:59
CellCenterResidualValue evalFluxAndSource(Args &&... args) const
for compatibility with FVLocalAssemblerBase
Definition staggeredlocalresidual.hh:170
void evalStorageForFace(FaceResidualValue &residual, const Problem &problem, const Element &element, const FVElementGeometry &fvGeometry, const ElementVolumeVariables &prevElemVolVars, const ElementVolumeVariables &curElemVolVars, const ElementFaceVariables &prevElemFaceVars, const ElementFaceVariables &curElemFaceVars, const SubControlVolumeFace &scvf) const
Evaluate the storage terms for a face residual.
Definition staggeredlocalresidual.hh:261
CellCenterResidualValue evalFluxAndSourceForCellCenter(const Element &element, const FVElementGeometry &fvGeometry, const ElementVolumeVariables &elemVolVars, const ElementFaceVariables &elemFaceVars, const ElementBoundaryTypes &bcTypes, const ElementFluxVariablesCache &elemFluxVarsCache) const
Convenience function to evaluate the flux and source terms for the cell center residual.
Definition staggeredlocalresidual.hh:66
void evalFluxForFace(FaceResidualValue &residual, const Problem &problem, const Element &element, const FVElementGeometry &fvGeometry, const ElementVolumeVariables &elemVolVars, const ElementFaceVariables &elemFaceVars, const ElementBoundaryTypes &elemBcTypes, const ElementFluxVariablesCache &elemFluxVarsCache, const SubControlVolumeFace &scvf) const
Evaluate the flux terms for a face residual.
Definition staggeredlocalresidual.hh:206
void evalSourceForCellCenter(CellCenterResidualValue &residual, const Problem &problem, const Element &element, const FVElementGeometry &fvGeometry, const ElementVolumeVariables &curElemVolVars, const ElementFaceVariables &curElemFaceVars, const SubControlVolume &scv) const
Evaluate the source terms for a cell center residual.
Definition staggeredlocalresidual.hh:104
CellCenterResidualValue ElementResidualVector
Definition staggeredlocalresidual.hh:56
const Implementation & asImp() const
Definition staggeredlocalresidual.hh:315
void evalFluxForCellCenter(CellCenterResidualValue &residual, const Problem &problem, const Element &element, const FVElementGeometry &fvGeometry, const ElementVolumeVariables &elemVolVars, const ElementFaceVariables &elemFaceVars, const ElementBoundaryTypes &elemBcTypes, const ElementFluxVariablesCache &elemFluxVarsCache, const SubControlVolumeFace &scvf) const
Evaluate the flux terms for a cell center residual.
Definition staggeredlocalresidual.hh:87
FaceResidualValue evalStorageForFace(const Element &element, const FVElementGeometry &fvGeometry, const ElementVolumeVariables &prevElemVolVars, const ElementVolumeVariables &curElemVolVars, const ElementFaceVariables &prevElemFaceVars, const ElementFaceVariables &curElemFaceVars, const SubControlVolumeFace &scvf) const
Evaluate the storage terms for a face residual.
Definition staggeredlocalresidual.hh:246
Manages the handling of time dependent problems.
Definition common/timeloop.hh:84
The default time loop for instationary simulations.
Definition common/timeloop.hh:139
Defines all properties used in Dumux.
Manages the handling of time dependent problems.
Helper classes to compute the integration elements.
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