12#ifndef DUMUX_FV_LOCAL_RESIDUAL_HH 
   13#define DUMUX_FV_LOCAL_RESIDUAL_HH 
   15#include <dune/common/exceptions.hh> 
   16#include <dune/istl/bvector.hh> 
   33template<
class TypeTag>
 
   40    using Element = 
typename GridView::template Codim<0>::Entity;
 
   44    using SubControlVolume = 
typename GridGeometry::SubControlVolume;
 
   45    using SubControlVolumeFace = 
typename GridGeometry::SubControlVolumeFace;
 
   87                                      const Element &element,
 
   88                                      const GridGeometry& gridGeometry,
 
   89                                      const GridVariables& gridVariables,
 
   90                                      const SolutionVector& sol)
 const 
   93        const auto fvGeometry = 
localView(gridGeometry).bind(element);
 
   94        const auto elemVolVars = 
localView(gridVariables.curGridVolVars()).bind(element, fvGeometry, sol);
 
  100        for (
auto&& scv : scvs(fvGeometry))
 
  102            auto localScvIdx = scv.localDofIndex();
 
  103            const auto& volVars = elemVolVars[scv];
 
  104            storage[localScvIdx] = 
asImp().computeStorage(
problem, scv, volVars);
 
  105            storage[localScvIdx] *= Extrusion::volume(fvGeometry, scv) * volVars.extrusionFactor();
 
 
  132                                      const FVElementGeometry& fvGeometry,
 
  133                                      const ElementVolumeVariables& prevElemVolVars,
 
  134                                      const ElementVolumeVariables& curElemVolVars)
 const 
  136        assert(timeLoop_ && 
"no time loop set for storage term evaluation");
 
  143        for (
auto&& scv : scvs(fvGeometry))
 
  144            asImp().evalStorage(residual, this->
problem(), element, fvGeometry, prevElemVolVars, curElemVolVars, scv);
 
 
  150                                            const FVElementGeometry& fvGeometry,
 
  151                                            const ElementVolumeVariables& elemVolVars,
 
  152                                            const ElementFluxVariablesCache& elemFluxVarsCache,
 
  153                                            const ElementBoundaryTypes &bcTypes)
 const 
  160        for (
auto&& scv : scvs(fvGeometry))
 
  161            asImp().evalSource(residual, this->
problem(), element, fvGeometry, elemVolVars, scv);
 
  164        for (
auto&& scvf : scvfs(fvGeometry))
 
  165            asImp().evalFlux(residual, this->
problem(), element, fvGeometry, elemVolVars, bcTypes, elemFluxVarsCache, scvf);
 
 
  189                               const SubControlVolume& scv,
 
  190                               const VolumeVariables& volVars)
 const 
  192        DUNE_THROW(Dune::NotImplemented, 
"This model does not implement a storage method!");
 
 
  209                              const Element& element,
 
  210                              const FVElementGeometry& fvGeometry,
 
  211                              const ElementVolumeVariables& elemVolVars,
 
  212                              const SubControlVolume &scv)
 const 
  214        NumEqVector source(0.0);
 
  217        source += 
problem.source(element, fvGeometry, elemVolVars, scv);
 
  220        if (!
problem.pointSourceMap().empty())
 
  221            source += 
problem.scvPointSources(element, fvGeometry, elemVolVars, scv);
 
 
  241                            const Element& element,
 
  242                            const FVElementGeometry& fvGeometry,
 
  243                            const ElementVolumeVariables& elemVolVars,
 
  244                            const SubControlVolumeFace& scvf,
 
  245                            const ElementFluxVariablesCache& elemFluxVarsCache)
 const 
  247        DUNE_THROW(Dune::NotImplemented, 
"This model does not implement a flux method!");
 
 
  275                     const Element& element,
 
  276                     const FVElementGeometry& fvGeometry,
 
  277                     const ElementVolumeVariables& prevElemVolVars,
 
  278                     const ElementVolumeVariables& curElemVolVars,
 
  279                     const SubControlVolume& scv)
 const 
  281        const auto& curVolVars = curElemVolVars[scv];
 
  282        const auto& prevVolVars = prevElemVolVars[scv];
 
  292        NumEqVector prevStorage = 
asImp().computeStorage(
problem, scv, prevVolVars);
 
  293        NumEqVector storage = 
asImp().computeStorage(
problem, scv, curVolVars);
 
  295        prevStorage *= prevVolVars.extrusionFactor();
 
  296        storage *= curVolVars.extrusionFactor();
 
  298        storage -= prevStorage;
 
  299        storage *= Extrusion::volume(fvGeometry, scv);
 
  300        storage /= timeLoop_->timeStepSize();
 
  302        residual[scv.localDofIndex()] += storage;
 
 
  320                    const Element& element,
 
  321                    const FVElementGeometry& fvGeometry,
 
  322                    const ElementVolumeVariables& curElemVolVars,
 
  323                    const SubControlVolume& scv)
 const 
  326        const auto& curVolVars = curElemVolVars[scv];
 
  327        NumEqVector source = 
asImp().computeSource(
problem, element, fvGeometry, curElemVolVars, scv);
 
  328        source *= Extrusion::volume(fvGeometry, scv)*curVolVars.extrusionFactor();
 
  331        residual[scv.localDofIndex()] -= source;
 
 
  350                  const Element& element,
 
  351                  const FVElementGeometry& fvGeometry,
 
  352                  const ElementVolumeVariables& elemVolVars,
 
  353                  const ElementBoundaryTypes& elemBcTypes,
 
  354                  const ElementFluxVariablesCache& elemFluxVarsCache,
 
  355                  const SubControlVolumeFace& scvf)
 const {}
 
 
  371                         const Element& element,
 
  372                         const FVElementGeometry& fvGeometry,
 
  373                         const ElementVolumeVariables& elemVolVars,
 
  374                         const ElementFluxVariablesCache& elemFluxVarsCache,
 
  375                         const SubControlVolumeFace& scvf)
 const 
  377        return asImp().evalFlux(
problem, element, fvGeometry, elemVolVars, elemFluxVarsCache, scvf);
 
 
  388    template<
class PartialDerivativeMatrix>
 
  391                               const Element& element,
 
  392                               const FVElementGeometry& fvGeometry,
 
  393                               const VolumeVariables& curVolVars,
 
  394                               const SubControlVolume& scv)
 const 
  396        DUNE_THROW(Dune::NotImplemented, 
"analytic storage derivative");
 
 
  400    template<
class PartialDerivativeMatrix>
 
  403                              const Element& element,
 
  404                              const FVElementGeometry& fvGeometry,
 
  405                              const VolumeVariables& curVolVars,
 
  406                              const SubControlVolume& scv)
 const 
  408        DUNE_THROW(Dune::NotImplemented, 
"analytic source derivative");
 
 
  412    template<
class PartialDerivativeMatrices, 
class T = TypeTag>
 
  416                            const Element& element,
 
  417                            const FVElementGeometry& fvGeometry,
 
  418                            const ElementVolumeVariables& curElemVolVars,
 
  419                            const ElementFluxVariablesCache& elemFluxVarsCache,
 
  420                            const SubControlVolumeFace& scvf)
 const 
  422        DUNE_THROW(Dune::NotImplemented, 
"analytic flux derivative for cell-centered models");
 
 
  426    template<
class JacobianMatrix, 
class T = TypeTag>
 
  430                            const Element& element,
 
  431                            const FVElementGeometry& fvGeometry,
 
  432                            const ElementVolumeVariables& curElemVolVars,
 
  433                            const ElementFluxVariablesCache& elemFluxVarsCache,
 
  434                            const SubControlVolumeFace& scvf)
 const 
  436        DUNE_THROW(Dune::NotImplemented, 
"analytic flux derivative for box models");
 
 
  440    template<
class PartialDerivativeMatrices>
 
  443                                     const Element& element,
 
  444                                     const FVElementGeometry& fvGeometry,
 
  445                                     const ElementVolumeVariables& curElemVolVars,
 
  446                                     const ElementFluxVariablesCache& elemFluxVarsCache,
 
  447                                     const SubControlVolumeFace& scvf)
 const 
  449        DUNE_THROW(Dune::NotImplemented, 
"analytic Dirichlet flux derivative");
 
 
  453    template<
class PartialDerivativeMatrices>
 
  456                                 const Element& element,
 
  457                                 const FVElementGeometry& fvGeometry,
 
  458                                 const ElementVolumeVariables& curElemVolVars,
 
  459                                 const ElementFluxVariablesCache& elemFluxVarsCache,
 
  460                                 const SubControlVolumeFace& scvf)
 const 
  462        DUNE_THROW(Dune::NotImplemented, 
"analytic Robin flux derivative");
 
 
  474    { 
return *problem_; }
 
 
  479    { 
return *timeLoop_; }
 
 
  483    { 
return !timeLoop_; }
 
 
  489    { 
return *
static_cast<Implementation*
>(
this); }
 
 
  492    { 
return *
static_cast<const Implementation*
>(
this); }
 
 
  495    const Problem* problem_; 
 
 
Implementation & asImp()
Definition fvlocalresidual.hh:488
NumEqVector computeSource(const Problem &problem, const Element &element, const FVElementGeometry &fvGeometry, const ElementVolumeVariables &elemVolVars, const SubControlVolume &scv) const
Calculate the source term of the equation.
Definition fvlocalresidual.hh:208
void addRobinFluxDerivatives(PartialDerivativeMatrices &derivativeMatrices, const Problem &problem, const Element &element, const FVElementGeometry &fvGeometry, const ElementVolumeVariables &curElemVolVars, const ElementFluxVariablesCache &elemFluxVarsCache, const SubControlVolumeFace &scvf) const
Compute the derivative of Robin type boundary conditions ("solution dependent Neumann")
Definition fvlocalresidual.hh:454
NumEqVector computeStorage(const Problem &problem, const SubControlVolume &scv, const VolumeVariables &volVars) const
Calculate the source term of the equation.
Definition fvlocalresidual.hh:188
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
Compute the flux local residual, i.e. the deviation of the flux term from zero.
Definition fvlocalresidual.hh:348
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
void addSourceDerivatives(PartialDerivativeMatrix &partialDerivatives, const Problem &problem, const Element &element, const FVElementGeometry &fvGeometry, const VolumeVariables &curVolVars, const SubControlVolume &scv) const
Compute the derivative of the source residual.
Definition fvlocalresidual.hh:401
const TimeLoop & timeLoop() const
Definition fvlocalresidual.hh:478
void addStorageDerivatives(PartialDerivativeMatrix &partialDerivatives, const Problem &problem, const Element &element, const FVElementGeometry &fvGeometry, const VolumeVariables &curVolVars, const SubControlVolume &scv) const
Compute the derivative of the storage residual.
Definition fvlocalresidual.hh:389
bool isStationary() const
returns true if the residual is stationary
Definition fvlocalresidual.hh:482
NumEqVector evalFlux(const Problem &problem, const Element &element, const FVElementGeometry &fvGeometry, const ElementVolumeVariables &elemVolVars, const ElementFluxVariablesCache &elemFluxVarsCache, const SubControlVolumeFace &scvf) const
Compute the flux local residual, i.e. the deviation of the flux term from zero.
Definition fvlocalresidual.hh:370
const Implementation & asImp() const
Definition fvlocalresidual.hh:491
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
ElementResidualVector evalStorage(const Element &element, const FVElementGeometry &fvGeometry, const ElementVolumeVariables &prevElemVolVars, const ElementVolumeVariables &curElemVolVars) const
Compute the storage local residual, i.e. the deviation of the storage term from zero for instationary...
Definition fvlocalresidual.hh:131
std::enable_if_t< GetPropType< T, Properties::GridGeometry >::discMethod !=DiscretizationMethods::box, void > addFluxDerivatives(PartialDerivativeMatrices &derivativeMatrices, const Problem &problem, const Element &element, const FVElementGeometry &fvGeometry, const ElementVolumeVariables &curElemVolVars, const ElementFluxVariablesCache &elemFluxVarsCache, const SubControlVolumeFace &scvf) const
Compute the derivative of the flux residual.
Definition fvlocalresidual.hh:414
void evalSource(ElementResidualVector &residual, const Problem &problem, const Element &element, const FVElementGeometry &fvGeometry, const ElementVolumeVariables &curElemVolVars, const SubControlVolume &scv) const
Compute the source local residual, i.e. the deviation of the source term from zero.
Definition fvlocalresidual.hh:318
std::enable_if_t< GetPropType< T, Properties::GridGeometry >::discMethod==DiscretizationMethods::box, void > addFluxDerivatives(JacobianMatrix &A, const Problem &problem, const Element &element, const FVElementGeometry &fvGeometry, const ElementVolumeVariables &curElemVolVars, const ElementFluxVariablesCache &elemFluxVarsCache, const SubControlVolumeFace &scvf) const
Compute the derivative of the flux residual for the box method.
Definition fvlocalresidual.hh:428
ElementResidualVector evalFluxAndSource(const Element &element, const FVElementGeometry &fvGeometry, const ElementVolumeVariables &elemVolVars, const ElementFluxVariablesCache &elemFluxVarsCache, const ElementBoundaryTypes &bcTypes) const
Definition fvlocalresidual.hh:149
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 fvlocalresidual.hh:273
void addCCDirichletFluxDerivatives(PartialDerivativeMatrices &derivativeMatrices, const Problem &problem, const Element &element, const FVElementGeometry &fvGeometry, const ElementVolumeVariables &curElemVolVars, const ElementFluxVariablesCache &elemFluxVarsCache, const SubControlVolumeFace &scvf) const
Compute the derivative of the Dirichlet flux residual for cell-centered schemes.
Definition fvlocalresidual.hh:441
NumEqVector computeFlux(const Problem &problem, const Element &element, const FVElementGeometry &fvGeometry, const ElementVolumeVariables &elemVolVars, const SubControlVolumeFace &scvf, const ElementFluxVariablesCache &elemFluxVarsCache) const
Calculate the flux term of the equation.
Definition fvlocalresidual.hh:240
A arithmetic block vector type based on DUNE's reserved vector.
Definition reservedblockvector.hh:26
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 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
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
The available discretization methods in Dumux.
constexpr Box box
Definition method.hh:147
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.