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>
42 #endif // HAVE_DUNE_LOCALFUNCTIONS
44 #include <dune/common/fvector.hh>
48 #ifdef HAVE_DUNE_LOCALFUNCTIONS
49 #define EWOMS_NO_LOCALFUNCTIONS_UNUSED
51 #define EWOMS_NO_LOCALFUNCTIONS_UNUSED OPM_UNUSED
64 template<
class TypeTag>
68 typedef typename GET_PROP_TYPE(TypeTag, Scalar) Scalar;
69 typedef typename GET_PROP_TYPE(TypeTag, Evaluation) Evaluation;
70 typedef typename GET_PROP_TYPE(TypeTag, GridView) GridView;
71 typedef typename GET_PROP_TYPE(TypeTag, ElementContext) ElementContext;
73 enum { dim = GridView::dimension };
78 enum { maxDof = (2 << dim) };
79 enum { maxFap = maxDof };
81 typedef typename GridView::ctype CoordScalar;
82 typedef Dune::FieldVector<Scalar, dim> DimVector;
84 #if HAVE_DUNE_LOCALFUNCTIONS
85 typedef Dune::PQkLocalFiniteElementCache<CoordScalar, Scalar, dim, 1> LocalFiniteElementCache;
86 typedef typename LocalFiniteElementCache::FiniteElementType LocalFiniteElement;
87 typedef typename LocalFiniteElement::Traits::LocalBasisType::Traits LocalBasisTraits;
88 typedef typename LocalBasisTraits::JacobianType ShapeJacobian;
89 #endif // HAVE_DUNE_LOCALFUNCTIONS
98 template <
bool prepareValues = true,
bool prepareGradients = true>
99 void prepare(
const ElementContext& elemCtx EWOMS_NO_LOCALFUNCTIONS_UNUSED,
100 unsigned timeIdx EWOMS_NO_LOCALFUNCTIONS_UNUSED)
103 #if !HAVE_DUNE_LOCALFUNCTIONS
105 OPM_THROW(std::logic_error,
"The dune-localfunctions module is required in oder to use"
106 " finite element gradients");
108 const auto& stencil = elemCtx.stencil(timeIdx);
110 const LocalFiniteElement& localFE = feCache_.get(elemCtx.element().type());
111 localFiniteElement_ = &localFE;
114 for (
unsigned faceIdx = 0; faceIdx < stencil.numInteriorFaces(); ++faceIdx) {
115 const auto& localFacePos = stencil.interiorFace(faceIdx).localPos();
120 localFE.localBasis().evaluateFunction(localFacePos, p1Value_[faceIdx]);
122 if (prepareGradients) {
124 std::vector<ShapeJacobian> localGradient;
125 localFE.localBasis().evaluateJacobian(localFacePos, localGradient);
130 const auto& geom = elemCtx.element().geometry();
131 const auto& jacInvT =
132 geom.jacobianInverseTransposed(localFacePos);
134 size_t numVertices = elemCtx.numDof(timeIdx);
135 for (
unsigned vertIdx = 0; vertIdx < numVertices; vertIdx++) {
136 jacInvT.mv(localGradient[vertIdx][0],
137 p1Gradient_[faceIdx][vertIdx]);
144 ParentType::template prepare<prepareValues, prepareGradients>(elemCtx, timeIdx);
158 template <
class QuantityCallback>
160 unsigned fapIdx EWOMS_NO_LOCALFUNCTIONS_UNUSED,
161 const QuantityCallback& quantityCallback EWOMS_NO_LOCALFUNCTIONS_UNUSED)
const
162 ->
typename std::remove_reference<typename QuantityCallback::ResultType>::type
165 #if !HAVE_DUNE_LOCALFUNCTIONS
167 OPM_THROW(std::logic_error,
"The dune-localfunctions module is required in oder to use"
168 " finite element gradients");
170 typedef typename std::remove_reference<typename QuantityCallback::ResultType>::type QuantityConstType;
171 typedef typename std::remove_const<QuantityConstType>::type QuantityType;
172 typedef Opm::MathToolbox<QuantityType> Toolbox;
176 QuantityType value(0.0);
177 for (
unsigned vertIdx = 0; vertIdx < elemCtx.numDof(0); ++vertIdx) {
178 if (std::is_same<QuantityType, Scalar>::value ||
179 elemCtx.focusDofIndex() == vertIdx)
180 value += quantityCallback(vertIdx)*p1Value_[fapIdx][vertIdx];
182 value += Toolbox::value(quantityCallback(vertIdx))*p1Value_[fapIdx][vertIdx];
203 template <
class QuantityCallback>
205 unsigned fapIdx EWOMS_NO_LOCALFUNCTIONS_UNUSED,
206 const QuantityCallback& quantityCallback EWOMS_NO_LOCALFUNCTIONS_UNUSED)
const
207 ->
typename std::remove_reference<typename QuantityCallback::ResultType>::type
210 #if !HAVE_DUNE_LOCALFUNCTIONS
212 OPM_THROW(std::logic_error,
"The dune-localfunctions module is required in oder to use"
213 " finite element gradients");
215 typedef typename std::remove_reference<typename QuantityCallback::ResultType>::type QuantityConstType;
216 typedef typename std::remove_const<QuantityConstType>::type QuantityType;
218 typedef decltype(std::declval<QuantityType>()[0]) RawFieldType;
219 typedef typename std::remove_const<typename std::remove_reference<RawFieldType>::type>::type FieldType;
221 typedef Opm::MathToolbox<FieldType> Toolbox;
225 QuantityType value(0.0);
226 for (
unsigned vertIdx = 0; vertIdx < elemCtx.numDof(0); ++vertIdx) {
227 if (std::is_same<QuantityType, Scalar>::value ||
228 elemCtx.focusDofIndex() == vertIdx)
230 const auto& tmp = quantityCallback(vertIdx);
231 for (
unsigned k = 0; k < tmp.size(); ++k)
232 value[k] += tmp[k]*p1Value_[fapIdx][vertIdx];
235 const auto& tmp = quantityCallback(vertIdx);
236 for (
unsigned k = 0; k < tmp.size(); ++k)
237 value[k] += Toolbox::value(tmp[k])*p1Value_[fapIdx][vertIdx];
259 template <
class QuantityCallback,
class EvalDimVector>
261 const ElementContext& elemCtx EWOMS_NO_LOCALFUNCTIONS_UNUSED,
262 unsigned fapIdx EWOMS_NO_LOCALFUNCTIONS_UNUSED,
263 const QuantityCallback& quantityCallback EWOMS_NO_LOCALFUNCTIONS_UNUSED)
const
266 #if !HAVE_DUNE_LOCALFUNCTIONS
268 OPM_THROW(std::logic_error,
"The dune-localfunctions module is required in oder to use"
269 " finite element gradients");
271 typedef typename std::remove_reference<typename QuantityCallback::ResultType>::type QuantityConstType;
272 typedef typename std::remove_const<QuantityConstType>::type QuantityType;
277 for (
unsigned vertIdx = 0; vertIdx < elemCtx.numDof(0); ++vertIdx) {
278 if (std::is_same<QuantityType, Scalar>::value ||
279 elemCtx.focusDofIndex() == vertIdx)
281 const auto& dofVal = quantityCallback(vertIdx);
282 const auto& tmp = p1Gradient_[fapIdx][vertIdx];
283 for (
int dimIdx = 0; dimIdx < dim; ++ dimIdx)
284 quantityGrad[dimIdx] += dofVal*tmp[dimIdx];
287 const auto& dofVal = quantityCallback(vertIdx);
288 const auto& tmp = p1Gradient_[fapIdx][vertIdx];
289 for (
int dimIdx = 0; dimIdx < dim; ++ dimIdx)
290 quantityGrad[dimIdx] += Opm::scalarValue(dofVal)*tmp[dimIdx];
313 template <
class QuantityCallback>
316 const QuantityCallback& quantityCallback)
334 template <
class QuantityCallback,
class EvalDimVector>
336 const ElementContext& elemCtx,
338 const QuantityCallback& quantityCallback)
const
341 #if HAVE_DUNE_LOCALFUNCTIONS
342 static LocalFiniteElementCache& localFiniteElementCache()
347 #if HAVE_DUNE_LOCALFUNCTIONS
348 static LocalFiniteElementCache feCache_;
350 const LocalFiniteElement* localFiniteElement_;
351 std::vector<Dune::FieldVector<Scalar, 1>> p1Value_[maxFap];
352 DimVector p1Gradient_[maxFap][maxDof];
353 #endif // HAVE_DUNE_LOCALFUNCTIONS
356 #if HAVE_DUNE_LOCALFUNCTIONS
357 template<
class TypeTag>
358 typename P1FeGradientCalculator<TypeTag>::LocalFiniteElementCache
359 P1FeGradientCalculator<TypeTag>::feCache_;
363 #undef EWOMS_NO_LOCALFUNCTIONS_UNUSED
auto calculateVectorValue(const ElementContext &elemCtx EWOMS_NO_LOCALFUNCTIONS_UNUSED, unsigned fapIdx EWOMS_NO_LOCALFUNCTIONS_UNUSED, const QuantityCallback &quantityCallback EWOMS_NO_LOCALFUNCTIONS_UNUSED) 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:204
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:294
#define GET_PROP_VALUE(TypeTag, PropTagName)
Access the value attribute of a property for a type tag.
Definition: propertysystem.hh:469
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:314
void prepare(const ElementContext &elemCtx EWOMS_NO_LOCALFUNCTIONS_UNUSED, unsigned timeIdx EWOMS_NO_LOCALFUNCTIONS_UNUSED)
Precomputes the common values to calculate gradients and values of quantities at any flux approximati...
Definition: p1fegradientcalculator.hh:99
This class calculates gradients of arbitrary quantities at flux integration points using the two-poin...
This class calculates gradients of arbitrary quantities at flux integration points using the two-poin...
Definition: fvbasegradientcalculator.hh:48
This class calculates gradients of arbitrary quantities at flux integration points using first order ...
Definition: p1fegradientcalculator.hh:65
void calculateGradient(EvalDimVector &quantityGrad EWOMS_NO_LOCALFUNCTIONS_UNUSED, const ElementContext &elemCtx EWOMS_NO_LOCALFUNCTIONS_UNUSED, unsigned fapIdx EWOMS_NO_LOCALFUNCTIONS_UNUSED, const QuantityCallback &quantityCallback EWOMS_NO_LOCALFUNCTIONS_UNUSED) const
Calculates the gradient of an arbitrary quantity at any flux approximation point. ...
Definition: p1fegradientcalculator.hh:260
Declares the basic properties used by the common infrastructure of the vertex-centered finite volume ...
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:211
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:98
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:335
auto calculateBoundaryValue(const ElementContext &elemCtx OPM_UNUSED, unsigned fapIdx OPM_UNUSED, 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:273
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:146
auto calculateScalarValue(const ElementContext &elemCtx EWOMS_NO_LOCALFUNCTIONS_UNUSED, unsigned fapIdx EWOMS_NO_LOCALFUNCTIONS_UNUSED, const QuantityCallback &quantityCallback EWOMS_NO_LOCALFUNCTIONS_UNUSED) 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:159