3 #ifndef DUNE_FUNCTIONS_GRIDFUNCTIONS_FACENORMALGRIDFUNCTION_HH
4 #define DUNE_FUNCTIONS_GRIDFUNCTIONS_FACENORMALGRIDFUNCTION_HH
9 #include <dune/common/exceptions.hh>
10 #include <dune/common/typeutilities.hh>
11 #include <dune/common/rangeutilities.hh>
12 #include <dune/geometry/referenceelements.hh>
14 #include <dune/functions/gridfunctions/gridfunction.hh>
15 #include <dune/functions/gridfunctions/gridviewentityset.hh>
18 namespace Dune::Functions {
23 template<
class ReferenceElement,
class Coordinate>
24 auto closestFaceIndex(
const ReferenceElement& re,
const Coordinate& x)
26 auto closestFaceIndex = decltype(re.subEntity(0,1,0,1)){};
27 double closestFaceDistance = std::numeric_limits<double>::max();
28 for(
auto&& faceIndex : Dune::range(re.size(1)))
33 auto normal = re.integrationOuterNormal(faceIndex);
34 normal /= normal.two_norm();
35 auto c = re.position(faceIndex,1);
37 auto faceDistance = (c*normal);
38 if (faceDistance<closestFaceDistance)
40 closestFaceDistance = faceDistance;
41 closestFaceIndex = faceIndex;
44 return closestFaceIndex;
73 using Domain =
typename EntitySet::GlobalCoordinate;
74 using Range =
typename EntitySet::GlobalCoordinate;
82 using Geometry =
typename Element::Geometry;
83 static const int dimension = GV::dimension;
97 void bind(
const Element& element)
100 geometry_.emplace(element_.geometry());
112 return static_cast<bool>(geometry_);
126 auto&& re = Dune::referenceElement(*geometry_);
128 auto face = Impl::closestFaceIndex(re, x);
129 auto localNormal = re.integrationOuterNormal(face);
133 auto normal = Range{};
134 geometry_->jacobianInverseTransposed(x).mv(localNormal, normal);
135 normal /= normal.two_norm();
140 const Element& localContext()
const
146 friend typename Traits::LocalFunctionTraits::DerivativeInterface
derivative(
const LocalFunction& t)
148 DUNE_THROW(NotImplemented,
"not implemented");
152 std::optional<Geometry> geometry_;
165 DUNE_THROW(NotImplemented,
"not implemented");
171 DUNE_THROW(NotImplemented,
"not implemented");
177 return LocalFunction{};
187 EntitySet entitySet_;
Grid function implementing the piecewise element face normal.
Definition: facenormalgridfunction.hh:66
friend Traits::DerivativeInterface derivative(const FaceNormalGridFunction &t)
Not implemented.
Definition: facenormalgridfunction.hh:169
friend LocalFunction localFunction(const FaceNormalGridFunction &t)
Return a local-function associated to FaceNormalGridFunction.
Definition: facenormalgridfunction.hh:175
FaceNormalGridFunction(const GridView &gridView)
Construct the FaceNormalGridFunction.
Definition: facenormalgridfunction.hh:158
Range operator()(const Domain &x) const
Not implemented.
Definition: facenormalgridfunction.hh:163
const EntitySet & entitySet() const
Return the stored GridViewEntitySet.
Definition: facenormalgridfunction.hh:181
GridView::template Codim< codim >::Entity Element
Type of Elements contained in this EntitySet.
Definition: gridviewentityset.hh:32
Element::Geometry::LocalCoordinate LocalCoordinate
Type of local coordinates with respect to the Element.
Definition: gridviewentityset.hh:35
Default implementation for derivative traits.
Definition: defaultderivativetraits.hh:37