dune-grid 3.0-git
intersection.cc
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_ALBERTA_INTERSECTION_CC
4#define DUNE_ALBERTA_INTERSECTION_CC
5
7
8namespace Dune
9{
10
11 // AlbertaGridIntersectionBase
12 // ---------------------------
13
14 template< class Grid >
15 inline AlbertaGridIntersectionBase< Grid >
16 ::AlbertaGridIntersectionBase ()
17 : grid_( nullptr ),
18 elementInfo_(),
19 oppVertex_( -1 ) // mark invalid intersection
20 {}
21
22 template< class Grid >
24 ::AlbertaGridIntersectionBase ( const EntityImp &entity, const int oppVertex )
25 : grid_( &entity.grid() ),
26 elementInfo_( entity.elementInfo() ),
27 oppVertex_( oppVertex )
28 {}
29
30
31 template< class Grid >
32 inline typename Grid::template Codim< 0 >::Entity
34 {
36 return EntityImp( grid(), elementInfo(), 0 );
37 }
38
39
40 template< class Grid >
42 {
43 return elementInfo().isBoundary( oppVertex_ );
44 }
45
46
47 template< class Grid >
49 {
50 if( boundary() )
51 {
52 const int id = elementInfo().boundaryId( oppVertex_ );
53 assert( id != 0 );
54 return id;
55 }
56 else
57 return 0;
58 }
59
61 template< class Grid >
63 {
64 assert( boundary() );
65 const Alberta::BasicNodeProjection *projection = elementInfo().boundaryProjection( oppVertex_ );
66 assert( projection );
67 return projection->boundaryIndex();
68 }
69
71 template< class Grid >
73 {
74 const int face = (dimension > 1 ? oppVertex_ : 1-oppVertex_);
75 return grid().alberta2generic( 1, face );
76 }
77
79 template< class Grid >
80 inline GeometryType AlbertaGridIntersectionBase< Grid >::type () const
81 {
82 typedef typename GenericGeometry::SimplexTopology< dimension-1 >::type Topology;
83 return GeometryType( Topology() );
84 }
85
87 template< class Grid >
90 {
91 const typename Entity::Geometry geoInside = inside().geometry();
92
93 const int face = indexInInside();
94 const ReferenceElement< ctype, dimension > &refSimplex = ReferenceElements< ctype, dimension >::simplex();
95 const FieldVector< ctype, dimension > &refNormal = refSimplex.integrationOuterNormal( face );
96
97 const typename Entity::Geometry::JacobianInverseTransposed &jInvT
98 = Grid::getRealImplementation( geoInside ).jacobianInverseTransposed();
99
100 NormalVector normal;
101 jInvT.mv( refNormal, normal );
102 normal *= Grid::getRealImplementation( geoInside ).integrationElement();
103
104 return normal;
105 }
106
107 template<>
109 AlbertaGridIntersectionBase< const AlbertaGrid< 1, 1 > >::centerIntegrationOuterNormal () const
110 {
111 const Alberta::GlobalVector &oppCoord = grid().getCoord( elementInfo(), oppVertex_ );
112 const Alberta::GlobalVector &myCoord = grid().getCoord( elementInfo(), 1-oppVertex_ );
113 NormalVector n;
114 n[ 0 ] = (myCoord[ 0 ] > oppCoord[ 0 ] ? ctype( 1 ) : -ctype( 1 ));
115 return n;
116 }
117
118#if defined GRIDDIM && GRIDDIM > 1
119 template<>
121 AlbertaGridIntersectionBase< const AlbertaGrid< 2, 2 > >::centerIntegrationOuterNormal () const
122 {
123 const Alberta::GlobalVector &coordOne = grid().getCoord( elementInfo(), (oppVertex_+1)%3 );
124 const Alberta::GlobalVector &coordTwo = grid().getCoord( elementInfo(), (oppVertex_+2)%3 );
125
126 NormalVector n;
127 n[ 0 ] = -(coordOne[ 1 ] - coordTwo[ 1 ]);
128 n[ 1 ] = coordOne[ 0 ] - coordTwo[ 0 ];
129 return n;
130 }
131#endif // defined GRIDDIM && GRIDDIM > 1
132
133 template<>
134 inline AlbertaGridIntersectionBase< const AlbertaGrid< 3, 3 > >::NormalVector
135 AlbertaGridIntersectionBase< const AlbertaGrid< 3, 3 > >::centerIntegrationOuterNormal () const
136 {
137 // in this case the orientation is negative, multiply by -1
138 const ALBERTA EL_INFO &elInfo = elementInfo().elInfo();
139 const ctype val = (elInfo.orientation > 0) ? 1.0 : -1.0;
140
141 static const int faceVertices[ 4 ][ 3 ]
142 = { {1,3,2}, {0,2,3}, {0,3,1}, {0,1,2} };
143 const int *localFaces = faceVertices[ oppVertex_ ];
144
145 const Alberta::GlobalVector &coord0 = grid().getCoord( elementInfo(), localFaces[ 0 ] );
146 const Alberta::GlobalVector &coord1 = grid().getCoord( elementInfo(), localFaces[ 1 ] );
147 const Alberta::GlobalVector &coord2 = grid().getCoord( elementInfo(), localFaces[ 2 ] );
148
149 FieldVector< ctype, dimensionworld > u;
150 FieldVector< ctype, dimensionworld > v;
151 for( int i = 0; i < dimension; ++i )
152 {
153 v[ i ] = coord1[ i ] - coord0[ i ];
154 u[ i ] = coord2[ i ] - coord1[ i ];
155 }
156
157 NormalVector n;
158 for( int i = 0; i < dimension; ++i )
159 {
160 const int j = (i+1)%dimension;
161 const int k = (i+2)%dimension;
162 n[ i ] = val * (u[ j ] * v[ k ] - u[ k ] * v[ j ]);
163 }
164 return n;
165 }
166
167
168 template< class Grid >
171 {
172 return centerIntegrationOuterNormal();
173 }
174
175
176 template< class Grid >
179 {
180 NormalVector normal = centerOuterNormal();
181 normal *= (1.0 / normal.two_norm());
182 return normal;
183 }
184
185
186 template< class Grid >
189 {
190 return centerIntegrationOuterNormal();
191 }
192
193
194 template< class Grid >
197 {
198 return centerOuterNormal();
199 }
200
201
202 template< class Grid >
205 {
206 return centerUnitOuterNormal();
207 }
208
209
210 template< class Grid >
213 {
214 return AlbertaTransformation( elementInfo().transformation( oppVertex_ ) );
215 }
216
217
218 template< class Grid >
220 {
221 return *grid_;
222 }
223
224
225 template< class Grid >
228 {
229 assert( !!elementInfo_ );
230 return elementInfo_;
231 }
232
233
234
235 // AlbertaGridIntersectionBase::GlobalCoordReader
236 // ----------------------------------------------
237
238 template< class GridImp >
240 {
241 typedef typename std::remove_const< GridImp >::type Grid;
242
243 static const int dimension = Grid::dimension;
244 static const int codimension = 1;
245 static const int mydimension = dimension - codimension;
246 static const int coorddimension = Grid::dimensionworld;
247
249
251 typedef FieldVector< ctype, coorddimension > Coordinate;
252
253 private:
254 const Grid &grid_;
256 const int subEntity_;
257 const int twist_;
258
259 public:
260 GlobalCoordReader ( const GridImp &grid,
262 int subEntity )
263 : grid_( grid ),
265 subEntity_( subEntity ),
266 twist_( elementInfo.template twist< codimension >( subEntity ) )
267 {}
268
269 void coordinate ( int i, Coordinate &x ) const
270 {
271 assert( !elementInfo_ == false );
272 assert( (i >= 0) && (i <= mydimension) );
273
274 const int ti = Alberta::applyInverseTwist< mydimension >( twist_, i );
275 const int k = mapVertices( subEntity_, ti );
276 const Alberta::GlobalVector &coord = grid_.getCoord( elementInfo_, k );
277 for( int j = 0; j < coorddimension; ++j )
278 x[ j ] = coord[ j ];
279 }
280
281 bool hasDeterminant () const
282 {
283 return false;
284 }
285
287 {
288 assert( false );
289 return ctype( 0 );
290 }
291
292 private:
293 static int mapVertices ( int subEntity, int i )
294 {
296 }
297 };
298
299
300
301
302 // AlbertaGridIntersectionBase::LocalCoordReader
303 // ---------------------------------------------
304
305 template< class GridImp >
307 {
308 typedef typename std::remove_const< GridImp >::type Grid;
309
310 static const int dimension = Grid::dimension;
311 static const int codimension = 1;
312 static const int mydimension = dimension - codimension;
313 static const int coorddimension = dimension;
314
316
317 typedef FieldVector< ctype, coorddimension > Coordinate;
318
319 typedef typename Grid::template Codim< 0 >::Geometry ElementGeometry;
320 typedef typename Grid::template Codim< 1 >::Geometry FaceGeometry;
321
322 private:
323 const ElementGeometry &elementGeometry_;
324 const FaceGeometry &faceGeometry_;
325
326 public:
327 LocalCoordReader ( const ElementGeometry &elementGeometry,
328 const FaceGeometry &faceGeometry )
329 : elementGeometry_( elementGeometry ),
330 faceGeometry_( faceGeometry )
331 {}
332
333 void coordinate ( int i, Coordinate &x ) const
334 {
335 x = elementGeometry_.local( faceGeometry_.corner( i ) );
336 }
337
338 bool hasDeterminant () const
339 {
340 return false;
341 }
342
344 {
345 return ctype( 0 );
346 }
347 };
348
349
350
351 // AlbertaGridLeafIntersection
352 // ---------------------------
353
354 template< class GridImp >
356 ::AlbertaGridLeafIntersection ( const EntityImp &entity, const int n )
357 : Base( entity, n ),
358 neighborInfo_()
359 {}
360
361
362 template< class GridImp >
365 : Base( other ),
366 neighborInfo_()
367 {}
368
369
370 template< class GridImp >
373 {
374 *((Base *)this) = other;
375 neighborInfo_ = ElementInfo();
376 return *this;
377 }
378
379
380 template< class GridImp >
381 inline bool
383 {
384 return ((elementInfo() == other.elementInfo()) && (oppVertex_ == other.oppVertex_));
385 }
386
387 template< class GridImp >
389 {
390 assert( oppVertex_ <= dimension );
391 ++oppVertex_;
392 neighborInfo_ = ElementInfo();
393 }
394
395 template< class GridImp >
396 inline typename GridImp::template Codim< 0 >::Entity
398 {
400
401 if( !neighborInfo_ )
402 {
403 assert( neighbor() );
404
405 neighborInfo_ = elementInfo().leafNeighbor( oppVertex_ );
406 }
407
408 assert( !neighborInfo_ == false );
409 assert( neighborInfo_.el() != NULL );
410 return EntityImp( grid(), neighborInfo_, 0 );
411 }
412
413 template< class GridImp >
415 {
416 return true;
417 }
418
419 template< class GridImp >
421 {
422 assert( oppVertex_ <= dimension );
423 return elementInfo().hasLeafNeighbor( oppVertex_ );
424 }
425
426
427 template< class GridImp >
430 {
431 typedef AlbertaGridLocalGeometryProvider< GridImp > LocalGeoProvider;
432 const int twist = elementInfo().template twist< 1 >( oppVertex_ );
433 const int face = (dimension > 1 ? oppVertex_ : 1-oppVertex_);
434 return LocalGeometry( LocalGeoProvider::instance().faceGeometry( face, twist ) );
435 }
436
437
438 template< class GridImp >
441 {
442 assert( neighbor() );
443
444 typedef AlbertaGridLocalGeometryProvider< GridImp > LocalGeoProvider;
445 const ALBERTA EL_INFO &elInfo = elementInfo().elInfo();
446 const int oppVertex = elInfo.opp_vertex[ oppVertex_ ];
447 const int twist = elementInfo().twistInNeighbor( oppVertex_ );
448 const int face = (dimension > 1 ? oppVertex : 1-oppVertex);
449 return LocalGeometry( LocalGeoProvider::instance().faceGeometry( face, twist ) );
450 }
451
452
453 template< class GridImp >
456 {
457 const int face = (dimension > 1 ? oppVertex_ : 1-oppVertex_);
458 const GlobalCoordReader coordReader( grid(), elementInfo(), face );
459 return Geometry( GeometryImpl( coordReader ) );
460 }
461
462
463 template< class GridImp >
465 {
466 const ALBERTA EL_INFO &elInfo = elementInfo().elInfo();
467 const int oppVertex = elInfo.opp_vertex[ oppVertex_ ];
468 const int face = (dimension > 1 ? oppVertex : 1-oppVertex);
469 return grid().alberta2generic( 1, face );
470 }
471
472
473 template< class GridImp >
474 inline int
476 {
477 return elementInfo().template twist< 1 >( oppVertex_ );
478 }
479
480
481 template< class GridImp >
482 inline int
484 {
485 return elementInfo().twistInNeighbor( oppVertex_ );
486 }
487
488} // namespace Dune
489
490#endif // #ifndef DUNE_ALBERTA_INTERSECTION_CC
#define ALBERTA
Definition albertaheader.hh:27
Include standard header files.
Definition agrid.hh:60
ALBERTA REAL Real
Definition misc.hh:45
ALBERTA REAL_D GlobalVector
Definition misc.hh:47
Definition albertagrid/intersection.hh:103
LocalGeometry geometryInOutside() const
Definition intersection.cc:440
GridImp::template Codim< 0 >::Entity outside() const
Definition intersection.cc:397
int twistInInside() const
Definition intersection.cc:475
Base::EntityImp EntityImp
Definition albertagrid/intersection.hh:126
void next()
Definition intersection.cc:388
bool operator==(const This &other) const
Definition intersection.cc:382
Base::GeometryImpl GeometryImpl
Definition albertagrid/intersection.hh:128
int indexInOutside() const
Definition intersection.cc:464
LocalGeometry geometryInInside() const
Definition intersection.cc:429
Base::GlobalCoordReader GlobalCoordReader
Definition albertagrid/intersection.hh:131
Base::Geometry Geometry
Definition albertagrid/intersection.hh:120
This & operator=(const This &other)
Definition intersection.cc:372
Geometry geometry() const
Definition intersection.cc:455
bool neighbor() const
Definition intersection.cc:420
bool conforming() const
Definition intersection.cc:414
int twistInOutside() const
Definition intersection.cc:483
Base::LocalGeometry LocalGeometry
Definition albertagrid/intersection.hh:121
Definition albertagrid/entity.hh:47
Definition albertagrid/geometry.hh:474
FieldVector< ctype, coorddimension > Coordinate
Definition intersection.cc:251
Alberta::Real ctype
Definition intersection.cc:248
std::remove_const< GridImp >::type Grid
Definition intersection.cc:241
Alberta::ElementInfo< dimension > ElementInfo
Definition intersection.cc:250
ctype determinant() const
Definition intersection.cc:286
bool hasDeterminant() const
Definition intersection.cc:281
GlobalCoordReader(const GridImp &grid, const ElementInfo &elementInfo, int subEntity)
Definition intersection.cc:260
void coordinate(int i, Coordinate &x) const
Definition intersection.cc:269
Grid::template Codim< 1 >::Geometry FaceGeometry
Definition intersection.cc:320
Alberta::Real ctype
Definition intersection.cc:315
LocalCoordReader(const ElementGeometry &elementGeometry, const FaceGeometry &faceGeometry)
Definition intersection.cc:327
std::remove_const< GridImp >::type Grid
Definition intersection.cc:308
void coordinate(int i, Coordinate &x) const
Definition intersection.cc:333
bool hasDeterminant() const
Definition intersection.cc:338
ctype determinant() const
Definition intersection.cc:343
FieldVector< ctype, coorddimension > Coordinate
Definition intersection.cc:317
Grid::template Codim< 0 >::Geometry ElementGeometry
Definition intersection.cc:319
Definition albertagrid/intersection.hh:30
NormalVector centerOuterNormal() const
Definition intersection.cc:170
const ElementInfo & elementInfo() const
Definition intersection.cc:227
FieldVector< ctype, dimensionworld > NormalVector
Definition albertagrid/intersection.hh:39
NormalVector centerIntegrationOuterNormal() const
Definition intersection.cc:89
Grid::ctype ctype
Definition albertagrid/intersection.hh:34
NormalVector integrationOuterNormal(const LocalCoordType &local) const
Definition intersection.cc:188
ElementInfo elementInfo_
Definition albertagrid/intersection.hh:91
FieldVector< ctype, dimension-1 > LocalCoordType
Definition albertagrid/intersection.hh:40
NormalVector centerUnitOuterNormal() const
Definition intersection.cc:178
int indexInInside() const
Definition intersection.cc:72
Entity inside() const
Definition intersection.cc:33
NormalVector unitOuterNormal(const LocalCoordType &local) const
Definition intersection.cc:204
const Grid & grid() const
Definition intersection.cc:219
int oppVertex_
Definition albertagrid/intersection.hh:92
static const int dimension
Definition albertagrid/intersection.hh:36
bool boundary() const
Definition intersection.cc:41
size_t boundarySegmentIndex() const
Definition intersection.cc:62
GeometryType type() const
Definition intersection.cc:80
AlbertaTransformation transformation() const
Definition intersection.cc:212
NormalVector outerNormal(const LocalCoordType &local) const
Definition intersection.cc:196
int boundaryId() const
Definition intersection.cc:48
const Grid * grid_
Definition albertagrid/intersection.hh:90
Definition misc.hh:440
Definition albertagrid/projection.hh:206
unsigned int boundaryIndex() const
Definition albertagrid/projection.hh:216
Definition transformation.hh:16
GridImp::template Codim< cd >::Geometry Geometry
The corresponding geometry type.
Definition common/entity.hh:97
Grid abstract base class.
Definition common/grid.hh:373
@ dimension
The dimension of the grid.
Definition common/grid.hh:387
@ dimensionworld
The dimension of the world the grid lives in.
Definition common/grid.hh:393