dune-grid  2.9.0
geometrygrid/intersection.hh
Go to the documentation of this file.
1 // SPDX-FileCopyrightText: Copyright (C) DUNE Project contributors, see file LICENSE.md in module root
2 // SPDX-License-Identifier: LicenseRef-GPL-2.0-only-with-DUNE-exception
3 // -*- tab-width: 4; indent-tabs-mode: nil; c-basic-offset: 2 -*-
4 // vi: set et ts=4 sw=2 sts=2:
5 #ifndef DUNE_GEOGRID_INTERSECTION_HH
6 #define DUNE_GEOGRID_INTERSECTION_HH
7 
10 
11 namespace Dune
12 {
13 
14  namespace GeoGrid
15  {
16 
17  // Intersection
18  // ------------
19 
20  template< class Grid, class HostIntersection >
22  {
23  typedef typename HostIntersection::Geometry HostGeometry;
24  typedef typename HostIntersection::LocalGeometry HostLocalGeometry;
25 
26  typedef typename std::remove_const< Grid >::type::Traits Traits;
27 
28  public:
29  typedef typename Traits::ctype ctype;
30 
31  static const int dimension = Traits::dimension;
32  static const int dimensionworld = Traits::dimensionworld;
33 
34  typedef typename Traits::template Codim< 0 >::Entity Entity;
35  typedef typename Traits::template Codim< 1 >::Geometry Geometry;
36  typedef typename Traits::template Codim< 1 >::LocalGeometry LocalGeometry;
37 
38  typedef typename Traits::template Codim< 0 >::Geometry ElementGeometry;
39 
40  private:
42 
43  typedef typename Traits::template Codim< 0 >::EntityImpl EntityImpl;
44 
45  typedef typename Traits::template Codim< 1 >::GeometryImpl GeometryImpl;
46  typedef typename Traits::template Codim< 0 >::GeometryImpl ElementGeometryImpl;
47 
48  public:
49 
51  {}
52 
53  explicit Intersection ( const HostIntersection &hostIntersection, const ElementGeometryImpl &insideGeo )
54  : hostIntersection_( hostIntersection )
55  , insideGeo_ ( insideGeo )
56  , geo_( grid() )
57  {}
58 
59  explicit Intersection ( HostIntersection&& hostIntersection, const ElementGeometryImpl &insideGeo )
60  : hostIntersection_( std::move( hostIntersection ) )
61  , insideGeo_ ( insideGeo )
62  , geo_( grid() )
63  {}
64 
65  bool equals ( const Intersection &other) const
66  {
67  return hostIntersection_ == other.hostIntersection_;
68  }
69 
70  explicit operator bool () const { return bool( hostIntersection_ ); }
71 
72  Entity inside () const
73  {
74  return EntityImpl( insideGeo_, hostIntersection().inside() );
75  }
76 
77  Entity outside () const
78  {
79  return EntityImpl( grid(), hostIntersection().outside() );
80  }
81 
82  bool boundary () const { return hostIntersection().boundary(); }
83 
84  bool conforming () const { return hostIntersection().conforming(); }
85 
86  bool neighbor () const { return hostIntersection().neighbor(); }
87 
88  size_t boundarySegmentIndex () const
89  {
90  return hostIntersection().boundarySegmentIndex();
91  }
92 
94  {
95  return hostIntersection().geometryInInside();
96  }
97 
99  {
100  return hostIntersection().geometryInOutside();
101  }
102 
104  {
105  if( !geo_ )
106  {
107  CoordVector coords( insideGeo_, geometryInInside() );
108  geo_ = GeometryImpl( grid(), type(), coords );
109  }
110  return Geometry( geo_ );
111  }
112 
113  GeometryType type () const { return hostIntersection().type(); }
114 
115  int indexInInside () const
116  {
117  return hostIntersection().indexInInside();
118  }
119 
120  int indexInOutside () const
121  {
122  return hostIntersection().indexInOutside();
123  }
124 
125  FieldVector< ctype, dimensionworld >
126  integrationOuterNormal ( const FieldVector< ctype, dimension-1 > &local ) const
127  {
128  const LocalGeometry geoInInside = geometryInInside();
129  const int idxInInside = indexInInside();
130 
131  auto refElement = referenceElement< ctype, dimension >( insideGeo_.type() );
132 
133  FieldVector< ctype, dimension > x( geoInInside.global( local ) );
134  const typename ElementGeometryImpl::JacobianInverseTransposed &jit = insideGeo_.jacobianInverseTransposed( x );
135  FieldVector< ctype, dimension > refNormal = refElement.integrationOuterNormal( idxInInside );
136 
137  FieldVector< ctype, dimensionworld > normal;
138  jit.mv( refNormal, normal );
139  if( !conforming() )
140  normal *= geoInInside.volume() / refElement.template geometry< 1 >( idxInInside ).volume();
141  normal *= jit.detInv();
142  //normal *= insideGeo_.integrationElement( x );
143  return normal;
144  }
145 
146  FieldVector< ctype, dimensionworld >
147  outerNormal ( const FieldVector< ctype, dimension-1 > &local ) const
148  {
149  auto refElement = referenceElement< ctype, dimension >( insideGeo_.type() );
150 
151  FieldVector< ctype, dimension > x( geometryInInside().global( local ) );
152  const typename ElementGeometryImpl::JacobianInverseTransposed &jit = insideGeo_.jacobianInverseTransposed( x );
153  FieldVector< ctype, dimension > refNormal = refElement.integrationOuterNormal( indexInInside() );
154 
155  FieldVector< ctype, dimensionworld > normal;
156  jit.mv( refNormal, normal );
157  return normal;
158  }
159 
160  FieldVector< ctype, dimensionworld >
161  unitOuterNormal ( const FieldVector< ctype, dimension-1 > &local ) const
162  {
163  FieldVector< ctype, dimensionworld > normal = outerNormal( local );
164  normal *= (ctype( 1 ) / normal.two_norm());
165  return normal;
166  }
167 
168  FieldVector< ctype, dimensionworld > centerUnitOuterNormal () const
169  {
170  auto refFace = referenceElement< ctype, dimension-1 >( type() );
171  return unitOuterNormal( refFace.position( 0, 0 ) );
172  }
173 
174  const HostIntersection &hostIntersection () const
175  {
176  return hostIntersection_;
177  }
178 
179  const Grid &grid () const { return insideGeo_.grid(); }
180 
181  private:
182  HostIntersection hostIntersection_;
183  ElementGeometryImpl insideGeo_;
184  mutable GeometryImpl geo_;
185  };
186 
187  } // namespace GeoGrid
188 
189 } // namespace Dune
190 
191 #endif // #ifndef DUNE_GEOGRID_INTERSECTION_HH
Include standard header files.
Definition: agrid.hh:60
auto referenceElement(const Geometry< mydim, cdim, GridImp, GeometryImp > &geo) -> decltype(referenceElement(geo, geo.impl()))
Definition: common/geometry.hh:558
GeometryType
Type representing VTK's entity geometry types.
Definition: common.hh:132
Grid abstract base class.
Definition: common/grid.hh:375
Definition: cornerstorage.hh:123
Definition: geometrygrid/intersection.hh:22
FieldVector< ctype, dimensionworld > integrationOuterNormal(const FieldVector< ctype, dimension-1 > &local) const
Definition: geometrygrid/intersection.hh:126
Traits::template Codim< 0 >::Entity Entity
Definition: geometrygrid/intersection.hh:34
Traits::ctype ctype
Definition: geometrygrid/intersection.hh:29
const HostIntersection & hostIntersection() const
Definition: geometrygrid/intersection.hh:174
GeometryType type() const
Definition: geometrygrid/intersection.hh:113
LocalGeometry geometryInInside() const
Definition: geometrygrid/intersection.hh:93
Geometry geometry() const
Definition: geometrygrid/intersection.hh:103
bool boundary() const
Definition: geometrygrid/intersection.hh:82
Traits::template Codim< 1 >::Geometry Geometry
Definition: geometrygrid/intersection.hh:35
Intersection()
Definition: geometrygrid/intersection.hh:50
bool neighbor() const
Definition: geometrygrid/intersection.hh:86
size_t boundarySegmentIndex() const
Definition: geometrygrid/intersection.hh:88
FieldVector< ctype, dimensionworld > outerNormal(const FieldVector< ctype, dimension-1 > &local) const
Definition: geometrygrid/intersection.hh:147
FieldVector< ctype, dimensionworld > unitOuterNormal(const FieldVector< ctype, dimension-1 > &local) const
Definition: geometrygrid/intersection.hh:161
Traits::template Codim< 1 >::LocalGeometry LocalGeometry
Definition: geometrygrid/intersection.hh:36
Intersection(HostIntersection &&hostIntersection, const ElementGeometryImpl &insideGeo)
Definition: geometrygrid/intersection.hh:59
Traits::template Codim< 0 >::Geometry ElementGeometry
Definition: geometrygrid/intersection.hh:38
const Grid & grid() const
Definition: geometrygrid/intersection.hh:179
Intersection(const HostIntersection &hostIntersection, const ElementGeometryImpl &insideGeo)
Definition: geometrygrid/intersection.hh:53
static const int dimensionworld
Definition: geometrygrid/intersection.hh:32
int indexInOutside() const
Definition: geometrygrid/intersection.hh:120
bool equals(const Intersection &other) const
Definition: geometrygrid/intersection.hh:65
bool conforming() const
Definition: geometrygrid/intersection.hh:84
Entity outside() const
Definition: geometrygrid/intersection.hh:77
int indexInInside() const
Definition: geometrygrid/intersection.hh:115
FieldVector< ctype, dimensionworld > centerUnitOuterNormal() const
Definition: geometrygrid/intersection.hh:168
static const int dimension
Definition: geometrygrid/intersection.hh:31
Entity inside() const
Definition: geometrygrid/intersection.hh:72
LocalGeometry geometryInOutside() const
Definition: geometrygrid/intersection.hh:98