dune-grid 2.8.0
Loading...
Searching...
No Matches
yaspgridintersection.hh
Go to the documentation of this file.
1// -*- tab-width: 4; indent-tabs-mode: nil; c-basic-offset: 2 -*-
2// vi: set et ts=4 sw=2 sts=2:
3#ifndef DUNE_GRID_YASPGRIDINTERSECTION_HH
4#define DUNE_GRID_YASPGRIDINTERSECTION_HH
5
13namespace Dune {
14
18 template<class GridImp>
20 {
21 enum { dim=GridImp::dimension };
22 enum { dimworld=GridImp::dimensionworld };
23 typedef typename GridImp::ctype ctype;
24
25 typedef typename GridImp::Traits::template Codim< 1 >::GeometryImpl GeometryImpl;
26 typedef typename GridImp::Traits::template Codim< 1 >::LocalGeometryImpl LocalGeometryImpl;
27
28 friend class YaspIntersectionIterator<GridImp>;
29
30 public:
31 // types used from grids
32 typedef typename GridImp::YGridLevelIterator YGLI;
33 typedef typename GridImp::YGrid::Iterator I;
34 typedef typename GridImp::template Codim<0>::Entity Entity;
35 typedef typename GridImp::template Codim<1>::Geometry Geometry;
36 typedef typename GridImp::template Codim<1>::LocalGeometry LocalGeometry;
37
38 void update() {
39
40 // vector with per-direction movements
41 std::array<int,dim> dist{{0}};
42
43 // first move: back to center
44 dist[_dir] = 1 - 2*_face;
45
46 // update face info
47 _dir = _count / 2;
48 _face = _count % 2;
49
50 // second move: to new neighbor
51 dist[_dir] += -1 + 2*_face;
52
53 // move transforming iterator
54 _outside.transformingsubiterator().move(dist);
55 }
56
60 bool boundary () const
61 {
62 // Coordinate of intersection in its direction
63 int coord = _inside.transformingsubiterator().coord(_dir) + _face;
64 if (_inside.gridlevel()->mg->isPeriodic(_dir))
65 return false;
66 else
67 return coord == 0
68 ||
69 coord == _inside.gridlevel()->mg->levelSize(_inside.gridlevel()->level(),_dir);
70 }
71
73 bool neighbor () const
74 {
75 // Coordinate of intersection in its direction
76 int coord = _inside.transformingsubiterator().coord(_dir) + _face;
77 return coord > _inside.gridlevel()->overlap[0].dataBegin()->min(_dir)
78 &&
79 coord <= _inside.gridlevel()->overlap[0].dataBegin()->max(_dir);
80 }
81
83 bool conforming () const
84 {
85 return true;
86 }
87
90 Entity inside() const
91 {
92 return Entity(_inside);
93 }
94
97 {
98 return Entity(_outside);
99 }
100
104 {
105 if(! boundary())
106 DUNE_THROW(GridError, "called boundarySegmentIndex while boundary() == false");
107 // size of local macro grid
108 const std::array<int, dim> & size = _inside.gridlevel()->mg->begin()->overlap[0].dataBegin()->size();
109 const std::array<int, dim> & origin = _inside.gridlevel()->mg->begin()->overlap[0].dataBegin()->origin();
110 std::array<int, dim> sides;
111 {
112 for (int i=0; i<dim; i++)
113 {
114 sides[i] =
115 ((_inside.gridlevel()->mg->begin()->overlap[0].dataBegin()->origin(i)
116 == 0)+
117 (_inside.gridlevel()->mg->begin()->overlap[0].dataBegin()->origin(i) +
118 _inside.gridlevel()->mg->begin()->overlap[0].dataBegin()->size(i)
119 ==
120 _inside.gridlevel()->mg->levelSize(0,i)));
121
122 }
123 }
124 // global position of the cell on macro grid
125 std::array<int, dim> pos = _inside.transformingsubiterator().coord();
126 for(int i=0; i<dim; i++)
127 {
128 pos[i] = pos[i] / (1<<_inside.level());
129 pos[i] = pos[i] - origin[i];
130 }
131 // compute unit-cube-face-sizes
132 std::array<int, dim> fsize;
133 {
134 int vol = 1;
135 for (int k=0; k<dim; k++)
136 vol *= size[k];
137 for (int k=0; k<dim; k++)
138 fsize[k] = vol/size[k];
139 }
140 // compute index in the unit-cube-face
141 int index = 0;
142 {
143 int localoffset = 1;
144 for (int k=dim-1; k>=0; k--)
145 {
146 if (k == _dir) continue;
147 index += (pos[k]) * localoffset;
148 localoffset *= size[k];
149 }
150 }
151 // add unit-cube-face-offsets
152 {
153 for (int k=0; k<_dir; k++)
154 index += sides[k] * fsize[k];
155 // add fsize if we are on the right face and there is a left-face-boundary on this processor
156 index += _face * (sides[_dir]>1) * fsize[_dir];
157 }
158
159 return index;
160 }
161
163 FieldVector<ctype, dimworld> outerNormal (const FieldVector<ctype, dim-1>& /* local */) const
164 {
165 return centerUnitOuterNormal();
166 }
167
169 FieldVector<ctype, dimworld> unitOuterNormal (const FieldVector<ctype, dim-1>& /* local */) const
170 {
171 return centerUnitOuterNormal();
172 }
173
175 FieldVector<ctype, dimworld> centerUnitOuterNormal () const
176 {
177 FieldVector<ctype, dimworld> normal(0);
178 normal[_dir] = (_face==0) ? -1.0 : 1.0;
179 return normal;
180 }
181
185 FieldVector<ctype, dimworld> integrationOuterNormal ([[maybe_unused]] const FieldVector<ctype, dim-1>& local) const
186 {
187 return geometry().volume() * centerUnitOuterNormal();
188 }
189
194 {
195 // set of dimensions that span the intersection
196 std::bitset<dim> s;
197 s.set();
198 s[_dir] = false;
199
200 // lower-left and upper-right corners
201 Dune::FieldVector<ctype, dim> ll(0.0);
202 Dune::FieldVector<ctype, dim> ur(1.0);
203
204 ll[_dir] = ur[_dir] = (_face==0) ? 0.0 : 1.0;
205
206 return LocalGeometry(LocalGeometryImpl(ll,ur,s));
207 }
208
213 {
214 // set of dimensions that span the intersection
215 std::bitset<dim> s;
216 s.set();
217 s[_dir] = false;
218
219 // lower-left and upper-right corners
220 Dune::FieldVector<ctype, dim> ll(0.0);
221 Dune::FieldVector<ctype, dim> ur(1.0);
222
223 ll[_dir] = ur[_dir] = (_face==1) ? 0.0 : 1.0;
224
225 return LocalGeometry(LocalGeometryImpl(ll,ur,s));
226 }
227
231 {
232
233 std::bitset<dim> shift;
234 shift.set();
235 shift[_dir] = false;
236
237 Dune::FieldVector<ctype,dimworld> ll, ur;
238 for (int i=0; i<dimworld; i++)
239 {
240 int coord = _inside.transformingsubiterator().coord(i);
241
242 if ((i == _dir) and (_face))
243 coord++;
244
245 ll[i] = _inside.transformingsubiterator().coordCont()->coordinate(i,coord);
246 if (i != _dir)
247 coord++;
248 ur[i] = _inside.transformingsubiterator().coordCont()->coordinate(i,coord);
249
250 // If on periodic overlap, transform coordinates by domain size
251 if (_inside.gridlevel()->mg->isPeriodic(i)) {
252 int coordPeriodic = _inside.transformingsubiterator().coord(i);
253 if (coordPeriodic < 0) {
254 auto size = _inside.gridlevel()->mg->domainSize()[i];
255 ll[i] += size;
256 ur[i] += size;
257 } else if (coordPeriodic + 1 > _inside.gridlevel()->mg->levelSize(_inside.gridlevel()->level(),i)) {
258 auto size = _inside.gridlevel()->mg->domainSize()[i];
259 ll[i] -= size;
260 ur[i] -= size;
261 }
262 }
263 }
264
265 GeometryImpl _is_global(ll,ur,shift);
266 return Geometry( _is_global );
267 }
268
270 GeometryType type () const
271 {
272 return GeometryTypes::cube(dim-1);
273 }
274
276 int indexInInside () const
277 {
278 return _count;
279 }
280
282 int indexInOutside () const
283 {
284 // flip the last bit
285 return _count^1;
286 }
287
289 : _count(~std::uint8_t(0)) // Use as marker for invalid intersection
290 , _dir(0)
291 , _face(0)
292 {}
293
295 YaspIntersection (const YaspEntity<0,dim,GridImp>& myself, bool toend) :
296 _inside(myself.gridlevel(),
297 myself.transformingsubiterator()),
298 _outside(myself.gridlevel(),
299 myself.transformingsubiterator()),
300 // initialize to first neighbor
301 _count(0),
302 _dir(0),
303 _face(0)
304 {
305 if (toend)
306 {
307 // initialize end iterator
308 _count = 2*dim;
309 return;
310 }
311 _count = 0;
312
313 // move transforming iterator
314 _outside.transformingsubiterator().move(_dir,-1);
315 }
316
318
320 void assign (const YaspIntersection& it)
321 {
322 *this = it;
323 }
324
325 bool equals(const YaspIntersection& other) const
326 {
327 // compare counts first -- that's cheaper if the test fails
328 return _count == other._count && _inside.equals(other._inside);
329 }
330
331 private:
332 /* The two entities that make up the intersection */
335 /* current position */
336 std::uint8_t _count;
337 std::uint8_t _dir;
338 std::uint8_t _face;
339 };
340} // namespace Dune
341
342#endif // DUNE_GRID_YASPGRIDINTERSECTION_HH
STL namespace.
Include standard header files.
Definition: agrid.hh:58
Base class for exceptions in Dune grid modules.
Definition: exceptions.hh:18
Definition: yaspgridentity.hh:262
int level() const
level of this element
Definition: yaspgridentity.hh:276
const YGLI & gridlevel() const
Definition: yaspgridentity.hh:406
const I & transformingsubiterator() const
Definition: yaspgridentity.hh:405
bool equals(const YaspEntity &e) const
Return true when two iterators over the same grid are equal (!).
Definition: yaspgridentity.hh:333
YaspIntersectionIterator enables iteration over intersections with neighboring codim 0 entities.
Definition: yaspgridintersectioniterator.hh:20
YaspIntersection provides data about intersection with neighboring codim 0 entities.
Definition: yaspgridintersection.hh:20
FieldVector< ctype, dimworld > centerUnitOuterNormal() const
return unit outer normal at center of intersection geometry
Definition: yaspgridintersection.hh:175
bool equals(const YaspIntersection &other) const
Definition: yaspgridintersection.hh:325
Entity inside() const
Definition: yaspgridintersection.hh:90
Geometry geometry() const
Definition: yaspgridintersection.hh:230
GridImp::template Codim< 1 >::LocalGeometry LocalGeometry
Definition: yaspgridintersection.hh:36
YaspIntersection()
Definition: yaspgridintersection.hh:288
int indexInOutside() const
local index of codim 1 entity in neighbor where intersection is contained in
Definition: yaspgridintersection.hh:282
LocalGeometry geometryInInside() const
Definition: yaspgridintersection.hh:193
int boundarySegmentIndex() const
Definition: yaspgridintersection.hh:103
GridImp::template Codim< 1 >::Geometry Geometry
Definition: yaspgridintersection.hh:35
bool conforming() const
Yasp is always conform.
Definition: yaspgridintersection.hh:83
GridImp::template Codim< 0 >::Entity Entity
Definition: yaspgridintersection.hh:34
bool neighbor() const
return true if neighbor across intersection exists in this processor
Definition: yaspgridintersection.hh:73
FieldVector< ctype, dimworld > outerNormal(const FieldVector< ctype, dim-1 > &) const
return unit outer normal, this should be dependent on local coordinates for higher order boundary
Definition: yaspgridintersection.hh:163
YaspIntersection(const YaspEntity< 0, dim, GridImp > &myself, bool toend)
make intersection iterator from entity, initialize to first neighbor
Definition: yaspgridintersection.hh:295
GeometryType type() const
obtain the type of reference element for this intersection
Definition: yaspgridintersection.hh:270
void assign(const YaspIntersection &it)
copy constructor – use default
Definition: yaspgridintersection.hh:320
int indexInInside() const
local index of codim 1 entity in self where intersection is contained in
Definition: yaspgridintersection.hh:276
FieldVector< ctype, dimworld > integrationOuterNormal(const FieldVector< ctype, dim-1 > &local) const
Definition: yaspgridintersection.hh:185
GridImp::YGrid::Iterator I
Definition: yaspgridintersection.hh:33
LocalGeometry geometryInOutside() const
Definition: yaspgridintersection.hh:212
GridImp::YGridLevelIterator YGLI
Definition: yaspgridintersection.hh:32
Entity outside() const
return Entity on the outside of this intersection
Definition: yaspgridintersection.hh:96
bool boundary() const
Definition: yaspgridintersection.hh:60
FieldVector< ctype, dimworld > unitOuterNormal(const FieldVector< ctype, dim-1 > &) const
return unit outer normal, this should be dependent on local coordinates for higher order boundary
Definition: yaspgridintersection.hh:169
void update()
Definition: yaspgridintersection.hh:38