28#ifndef EWOMS_P1FE_GRADIENT_CALCULATOR_HH
29#define EWOMS_P1FE_GRADIENT_CALCULATOR_HH
35#include <dune/common/fvector.hh>
36#include <dune/common/version.hh>
38#include <dune/geometry/type.hh>
40#if HAVE_DUNE_LOCALFUNCTIONS
41#include <dune/localfunctions/lagrange/pqkfactory.hh>
44#include <dune/common/fvector.hh>
58template<
class TypeTag>
67 enum { dim = GridView::dimension };
72 enum { maxDof = (2 << dim) };
73 enum { maxFap = maxDof };
75 using CoordScalar =
typename GridView::ctype;
76 using DimVector = Dune::FieldVector<Scalar, dim>;
78#if HAVE_DUNE_LOCALFUNCTIONS
81 using LocalBasisTraits =
typename LocalFiniteElement::Traits::LocalBasisType::Traits;
82 using ShapeJacobian =
typename LocalBasisTraits::JacobianType;
92 template <
bool prepareValues = true,
bool prepareGradients = true>
97#if !HAVE_DUNE_LOCALFUNCTIONS
99 throw std::logic_error(
"The dune-localfunctions module is required in oder to use"
100 " finite element gradients");
102 const auto& stencil = elemCtx.stencil(
timeIdx);
108 for (
unsigned faceIdx = 0; faceIdx < stencil.numInteriorFaces(); ++faceIdx) {
109 const auto&
localFacePos = stencil.interiorFace(faceIdx).localPos();
124 const auto&
geom = elemCtx.element().geometry();
128 size_t numVertices = elemCtx.numDof(
timeIdx);
152 template <
class QuantityCallback>
156 ->
typename std::remove_reference<typename QuantityCallback::ResultType>::type
159#if !HAVE_DUNE_LOCALFUNCTIONS
161 throw std::logic_error(
"The dune-localfunctions module is required in oder to use"
162 " finite element gradients");
164 using QuantityConstType =
typename std::remove_reference<typename QuantityCallback::ResultType>::type;
165 using QuantityType =
typename std::remove_const<QuantityConstType>::type;
172 if (std::is_same<QuantityType, Scalar>::value ||
173 elemCtx.focusDofIndex() ==
vertIdx)
197 template <
class QuantityCallback>
201 ->
typename std::remove_reference<typename QuantityCallback::ResultType>::type
204#if !HAVE_DUNE_LOCALFUNCTIONS
206 throw std::logic_error(
"The dune-localfunctions module is required in oder to use"
207 " finite element gradients");
209 using QuantityConstType =
typename std::remove_reference<typename QuantityCallback::ResultType>::type;
210 using QuantityType =
typename std::remove_const<QuantityConstType>::type;
212 using RawFieldType =
decltype(std::declval<QuantityType>()[0]);
213 using FieldType =
typename std::remove_const<typename std::remove_reference<RawFieldType>::type>::type;
221 if (std::is_same<QuantityType, Scalar>::value ||
222 elemCtx.focusDofIndex() ==
vertIdx)
225 for (
unsigned k = 0;
k < tmp.size(); ++
k)
230 for (
unsigned k = 0;
k < tmp.size(); ++
k)
253 template <
class QuantityCallback,
class EvalDimVector>
260#if !HAVE_DUNE_LOCALFUNCTIONS
262 throw std::logic_error(
"The dune-localfunctions module is required in oder to use"
263 " finite element gradients");
265 using QuantityConstType =
typename std::remove_reference<typename QuantityCallback::ResultType>::type;
266 using QuantityType =
typename std::remove_const<QuantityConstType>::type;
272 if (std::is_same<QuantityType, Scalar>::value ||
273 elemCtx.focusDofIndex() ==
vertIdx)
307 template <
class QuantityCallback>
328 template <
class QuantityCallback,
class EvalDimVector>
330 const ElementContext& elemCtx,
335#if HAVE_DUNE_LOCALFUNCTIONS
341#if HAVE_DUNE_LOCALFUNCTIONS
345 std::vector<Dune::FieldVector<Scalar, 1>>
p1Value_[maxFap]{};
350#if HAVE_DUNE_LOCALFUNCTIONS
351template<
class TypeTag>
352typename P1FeGradientCalculator<TypeTag>::LocalFiniteElementCache
353P1FeGradientCalculator<TypeTag>::feCache_;
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
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
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
This class calculates gradients of arbitrary quantities at flux integration points using first order ...
Definition p1fegradientcalculator.hh:60
auto calculateScalarValue(const ElementContext &elemCtx, unsigned fapIdx, const QuantityCallback &quantityCallback) const -> typename std::remove_reference< typename QuantityCallback::ResultType >::type
Calculates the value of an arbitrary quantity at any interior flux approximation point.
Definition p1fegradientcalculator.hh:153
auto calculateVectorValue(const ElementContext &elemCtx, unsigned fapIdx, const QuantityCallback &quantityCallback) const -> typename std::remove_reference< typename QuantityCallback::ResultType >::type
Calculates the value of an arbitrary quantity at any interior flux approximation point.
Definition p1fegradientcalculator.hh:198
void calculateBoundaryGradient(EvalDimVector &quantityGrad, const ElementContext &elemCtx, unsigned fapIdx, const QuantityCallback &quantityCallback) const
Calculates the gradient of an arbitrary quantity at any flux approximation point on the boundary.
Definition p1fegradientcalculator.hh:329
auto calculateBoundaryValue(const ElementContext &elemCtx, unsigned fapIdx, const QuantityCallback &quantityCallback) -> decltype(ParentType::calculateBoundaryValue(elemCtx, fapIdx, quantityCallback))
Calculates the value of an arbitrary quantity at any flux approximation point on the grid boundary.
Definition p1fegradientcalculator.hh:308
void prepare(const ElementContext &elemCtx, unsigned timeIdx)
Precomputes the common values to calculate gradients and values of quantities at any flux approximati...
Definition p1fegradientcalculator.hh:93
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 p1fegradientcalculator.hh:254
This class calculates gradients of arbitrary quantities at flux integration points using the two-poin...
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
Declares the basic properties used by the common infrastructure of the vertex-centered finite volume ...