28#ifndef EWOMS_FV_BASE_GRADIENT_CALCULATOR_HH
29#define EWOMS_FV_BASE_GRADIENT_CALCULATOR_HH
33#include <dune/common/fvector.hh>
36template<
class TypeTag>
37class EcfvDiscretization;
45template<
class TypeTag>
54 enum { dim = GridView::dimension };
55 enum { dimWorld = GridView::dimensionworld };
59 enum { maxFap = 2 << dim };
61 using DimVector = Dune::FieldVector<Scalar, dimWorld>;
62 using EvalDimVector = Dune::FieldVector<Evaluation, dimWorld>;
79 template <
bool prepareValues = true,
bool prepareGradients = true>
80 void prepare(
const ElementContext&,
unsigned)
94 template <
class QuantityCallback>
98 ->
typename std::remove_reference<
decltype(
quantityCallback.operator()(0))>::type
101 using ReturnType =
typename std::remove_const<typename std::remove_reference<RawReturnType>::type>::type;
107 const auto&
face = elemCtx.stencil(0).interiorFace(
fapIdx);
108 auto i =
face.interiorIndex();
109 auto j =
face.exteriorIndex();
140 template <
class QuantityCallback>
144 ->
typename std::remove_reference<
decltype(
quantityCallback.operator()(0))>::type
147 using ReturnType =
typename std::remove_const<typename std::remove_reference<RawReturnType>::type>::type;
153 const auto&
face = elemCtx.stencil(0).interiorFace(
fapIdx);
154 auto i =
face.interiorIndex();
155 auto j =
face.exteriorIndex();
162 for (
int k = 0;
k < value.size(); ++
k)
183 for (
int k = 0;
k < value.size(); ++
k)
200 template <
class QuantityCallback>
202 const ElementContext& elemCtx,
206 const auto& stencil = elemCtx.stencil(0);
207 const auto&
face = stencil.interiorFace(
fapIdx);
209 auto i =
face.interiorIndex();
210 auto j =
face.exteriorIndex();
211 auto focusIdx = elemCtx.focusDofIndex();
213 const auto&
interiorPos = stencil.subControlVolume(i).globalPos();
214 const auto&
exteriorPos = stencil.subControlVolume(j).globalPos();
262 template <
class QuantityCallback>
283 template <
class QuantityCallback>
285 const ElementContext& elemCtx,
289 const auto& stencil = elemCtx.stencil(0);
290 const auto&
face = stencil.boundaryFace(faceIdx);
293 if (
face.interiorIndex() == elemCtx.focusDofIndex())
301 const auto&
interiorPos = stencil.subControlVolume(
face.interiorIndex()).center();
323 const ElementContext& elemCtx,
326 const auto& stencil = elemCtx.stencil(0);
327 const auto&
face = stencil.interiorFace(
fapIdx);
331 const auto& normal =
face.normal();
332 auto i =
face.interiorIndex();
333 auto j =
face.exteriorIndex();
334 const auto&
interiorPos = stencil.subControlVolume(i).globalPos();
335 const auto&
exteriorPos = stencil.subControlVolume(j).globalPos();
336 const auto& integrationPos =
face.integrationPos();
This class calculates gradients of arbitrary quantities at flux integration points using the two-poin...
Definition fvbasegradientcalculator.hh:47
void calculateBoundaryGradient(EvalDimVector &quantityGrad, const ElementContext &elemCtx, unsigned faceIdx, const QuantityCallback &quantityCallback) const
Calculates the gradient of an arbitrary quantity at any flux approximation point on the boundary.
Definition fvbasegradientcalculator.hh:284
void calculateGradient(EvalDimVector &quantityGrad, const ElementContext &elemCtx, unsigned fapIdx, const QuantityCallback &quantityCallback) const
Calculates the gradient of an arbitrary quantity at any flux approximation point.
Definition fvbasegradientcalculator.hh:201
void prepare(const ElementContext &, unsigned)
Precomputes the common values to calculate gradients and values of quantities at every interior flux ...
Definition fvbasegradientcalculator.hh:80
auto calculateScalarValue(const ElementContext &elemCtx, unsigned fapIdx, const QuantityCallback &quantityCallback) const -> typename std::remove_reference< decltype(quantityCallback.operator()(0))>::type
Calculates the value of an arbitrary scalar quantity at any interior flux approximation point.
Definition fvbasegradientcalculator.hh:95
static void registerParameters()
Register all run-time parameters for the gradient calculator of the base class of the discretization.
Definition fvbasegradientcalculator.hh:69
auto calculateBoundaryValue(const ElementContext &, unsigned, const QuantityCallback &quantityCallback) -> decltype(quantityCallback.boundaryValue())
Calculates the value of an arbitrary quantity at any flux approximation point on the grid boundary.
Definition fvbasegradientcalculator.hh:263
auto calculateVectorValue(const ElementContext &elemCtx, unsigned fapIdx, const QuantityCallback &quantityCallback) const -> typename std::remove_reference< decltype(quantityCallback.operator()(0))>::type
Calculates the value of an arbitrary vectorial quantity at any interior flux approximation point.
Definition fvbasegradientcalculator.hh:141
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