28#ifndef EWOMS_FV_BASE_LOCAL_RESIDUAL_HH
29#define EWOMS_FV_BASE_LOCAL_RESIDUAL_HH
36#include <opm/material/common/Valgrind.hpp>
38#include <dune/istl/bvector.hh>
39#include <dune/grid/common/geometry.hh>
41#include <dune/common/fvector.hh>
43#include <dune/common/classname.hh>
56template<
class TypeTag>
63 using Element =
typename GridView::template Codim<0>::Entity;
81 using EvalVector = Dune::FieldVector<Evaluation, numEq>;
88 using LocalEvalBlockVector = Dune::BlockVector<EvalVector,
aligned_allocator<EvalVector,
alignof(EvalVector)> >;
107 {
return internalResidual_; }
116 {
return internalResidual_[
dofIdx]; }
128 void eval(
const Problem& problem,
const Element& element)
130 ElementContext
elemCtx(problem);
146 size_t numDof =
elemCtx.numDof(0);
147 internalResidual_.resize(numDof);
148 asImp_().eval(internalResidual_,
elemCtx);
174 if (useVolumetricResidual) {
177 size_t numDof =
elemCtx.numDof(0);
183 assert(std::isfinite(dofVolume));
184 Valgrind::CheckDefined(dofVolume);
214 size_t numPrimaryDof =
elemCtx.numPrimaryDof(0);
237 Dune::FieldVector<Scalar, numEq> tmp;
238 asImp_().computeStorage(tmp,
251 if (
elemCtx.enableStorageCache()) {
263 Dune::FieldVector<Scalar, numEq> tmp;
264 size_t numPrimaryDof =
elemCtx.numPrimaryDof(0);
267 asImp_().computeStorage(tmp,
282 size_t numPrimaryDof =
elemCtx.numPrimaryDof(0);
310 unsigned i =
face.interiorIndex();
311 unsigned j =
face.exteriorIndex();
313 Valgrind::SetUndefined(
flux);
315 Valgrind::CheckDefined(
flux);
322 alpha *=
face.area();
323 Valgrind::CheckDefined(alpha);
353 for (
unsigned i=0; i < numDof; i++) {
354 for (
unsigned j = 0; j < numEq; ++ j) {
356 Valgrind::CheckDefined(
residual[i][j]);
376 const ElementContext&,
380 throw std::logic_error(
"Not implemented: The local residual "+Dune::className<Implementation>()
381 +
" does not implement the required method 'computeStorage()'");
392 const ElementContext&,
396 throw std::logic_error(
"Not implemented: The local residual "+Dune::className<Implementation>()
397 +
" does not implement the required method 'computeFlux()'");
407 const ElementContext&,
411 throw std::logic_error(
"Not implemented: The local residual "+Dune::className<Implementation>()
412 +
" does not implement the required method 'computeSource()'");
432 size_t numBoundaryFaces =
boundaryCtx.numBoundaryFaces(0);
444 size_t numDof =
elemCtx.numDof(0);
445 for (
unsigned i=0; i < numDof; i++) {
446 for (
unsigned j = 0; j < numEq; ++ j) {
448 Valgrind::CheckDefined(
residual[i][j]);
464 BoundaryRateVector values;
466 Valgrind::SetUndefined(values);
468 Valgrind::CheckDefined(values);
478 Valgrind::CheckDefined(values[
eqIdx]);
502 size_t numPrimaryDof =
elemCtx.numPrimaryDof(0);
504 Scalar extrusionFactor =
506 Valgrind::CheckDefined(extrusionFactor);
508 assert(extrusionFactor > 0.0);
510 elemCtx.stencil(0).subControlVolume(
dofIdx).volume() * extrusionFactor;
518 if (!extensiveStorageTerm &&
519 !std::is_same<Scalar, Evaluation>::value &&
530 Valgrind::CheckDefined(tmp);
535 if (
elemCtx.enableStorageCache()) {
536 const auto& model =
elemCtx.model();
538 if (model.newtonMethod().numIterations() == 0 &&
539 !
elemCtx.haveStashedIntensiveQuantities())
541 if (!
elemCtx.problem().recycleFirstIterationStorage()) {
546 elemCtx.updatePrimaryIntensiveQuantities(1);
561 Valgrind::CheckDefined(
tmp2);
570 Valgrind::CheckDefined(
tmp2);
578 Valgrind::CheckDefined(
tmp2);
583 double dt =
elemCtx.simulator().timeStepSize();
599 if (!extensiveStorageTerm &&
600 !std::is_same<Scalar, Evaluation>::value &&
618 size_t numDof =
elemCtx.numDof(0);
619 for (
unsigned i=0; i < numDof; i++) {
620 for (
unsigned j = 0; j < numEq; ++ j) {
622 Valgrind::CheckDefined(
residual[i][j]);
630 Implementation& asImp_()
631 {
return *
static_cast<Implementation*
>(
this); }
633 const Implementation& asImp_()
const
634 {
return *
static_cast<const Implementation*
>(
this); }
636 LocalEvalBlockVector internalResidual_;
This is a stand-alone version of boost::alignment::aligned_allocator from Boost 1....
Element-wise caculation of the residual matrix for models based on a finite volume spatial discretiza...
Definition fvbaselocalresidual.hh:58
void evalVolumeTerms_(LocalEvalBlockVector &residual, ElementContext &elemCtx) const
Add the change in the storage terms and the source term to the local residual of all sub-control volu...
Definition fvbaselocalresidual.hh:491
static void registerParameters()
Register all run-time parameters for the local residual.
Definition fvbaselocalresidual.hh:99
void computeSource(RateVector &, const ElementContext &, unsigned, unsigned) const
Calculate the source term of the equation.
Definition fvbaselocalresidual.hh:406
void evalFluxes(LocalEvalBlockVector &residual, const ElementContext &elemCtx, unsigned timeIdx) const
Add the flux term to a local residual.
Definition fvbaselocalresidual.hh:299
void computeStorage(EqVector &, const ElementContext &, unsigned, unsigned) const
Evaluate the amount all conservation quantities (e.g.
Definition fvbaselocalresidual.hh:375
void evalStorage(LocalEvalBlockVector &storage, const ElementContext &elemCtx, unsigned timeIdx) const
Calculate the amount of all conservation quantities stored in all element's sub-control volumes for a...
Definition fvbaselocalresidual.hh:204
void eval(const Problem &problem, const Element &element)
Compute the local residual, i.e.
Definition fvbaselocalresidual.hh:128
void evalBoundarySegment_(LocalEvalBlockVector &residual, const BoundaryContext &boundaryCtx, unsigned boundaryFaceIdx, unsigned timeIdx) const
Evaluate all boundary conditions for a single sub-control volume face to the local residual.
Definition fvbaselocalresidual.hh:459
void eval(LocalEvalBlockVector &residual, ElementContext &elemCtx) const
Compute the local residual, i.e.
Definition fvbaselocalresidual.hh:158
void computeFlux(RateVector &, const ElementContext &, unsigned, unsigned) const
Evaluates the total mass flux of all conservation quantities over a face of a sub-control volume.
Definition fvbaselocalresidual.hh:391
void evalBoundary_(LocalEvalBlockVector &residual, const ElementContext &elemCtx, unsigned timeIdx) const
Evaluate the boundary conditions of an element.
Definition fvbaselocalresidual.hh:419
const LocalEvalBlockVector & residual() const
Return the result of the eval() call using internal storage.
Definition fvbaselocalresidual.hh:106
void eval(ElementContext &elemCtx)
Compute the local residual, i.e.
Definition fvbaselocalresidual.hh:144
const EvalVector & residual(unsigned dofIdx) const
Return the result of the eval() call using internal storage.
Definition fvbaselocalresidual.hh:115
Definition alignedallocator.hh:94
Declare the properties used by the infrastructure code of the finite volume discretizations.
This file contains a set of helper functions used by VFPProd / VFPInj.
Definition blackoilboundaryratevector.hh:37
constexpr auto getPropValue()
get the value data member of a property
Definition propertysystem.hh:242
typename Properties::Detail::GetPropImpl< TypeTag, Property >::type::type GetPropType
get the type alias defined in the property (equivalent to old macro GET_PROP_TYPE(....
Definition propertysystem.hh:235
This file provides the infrastructure to retrieve run-time parameters.