#include <intersectioniterator.hh>
Template parameters are:
GridImp
Type that is a model of Dune::GridIntersectionIteratorImp
Class template that is a model of Dune::IntersectionIterator
Intersections are codimension 1 objects. These intersections are accessed via an IntersectionIterator. This allows the implementation of non-matching grids, as one face can now consist of several intersections. In a conforming mesh such an intersection corresponds to an entity of codimension 1 but in the general non-conforming case there will be no entity in the mesh that directly corresponds to the intersection. Thus, the IntersectionIterator describes these intersections implicitly.
The IntersectionIterator class template wraps an object of type IntersectionIteratorImp and forwards all member function calls to corresponding members of this class. In that sense IntersectionIterator defines the interface and IntersectionIteratorImp supplies the implementation.
On an entity e of codimension zero there are two ways to create IntersectionIterators by either using ilevelbegin() / ilevelend() or ileafbegin()/ileafend(). In the first case intersections with neighboring entities having the same level as e are traversed; in the second case ileafbegin()==ileafend() if e is not a leaf otherwise all intersections with neighboring leaf entities are traversed.
Consider a situation where two elements a and b have a common intersection. Element b has been refined into an element c and d, while a has not been refined. In one space dimension this situation is depicted in the figure below.
IntersectionIterator in a locally refined mesh.
Here the rule is the following: The LevelIntersectionIterator delivers all intersections with elements on the same level, the LeafIntersectionIterator delivers the intersections with all leaf elements if it has been started on a leaf element. Both iterators also stop at intersections with the grid boundary. According to this rule the level intersection iterator started at element a in the example above delivers an intersection with b and the left grid boundary, whereas the leaf intersection iterator returns c instead of b. Starting on entity c the level intersection iterator returns d and the intersection with the left boundary of the level 1 grid, but the leaf intersection iterator returns both d and a. Finally, starting on b the level intersection iterator returns a and the right boundary, but the leaf intersection iterator is empty since b is not a leaf entity of the grid. Starting on d both the level and the leaf intersection iterators will return the element c together with the right grid boundary.
Public Types | |
enum | { dimension = dim } |
Export grid dimension. More... | |
enum | { dimensionworld = dimworld } |
Export dimension of world. More... | |
typedef GridImp::template Codim< 0 >::Entity | Entity |
Type of entity that this IntersectionIterator belongs to. | |
typedef GridImp::template Codim< 0 >::EntityPointer | EntityPointer |
Pointer to the type of entities that this IntersectionIterator belongs to. | |
typedef Dune::Intersection < const GridImp, IntersectionImp > | Intersection |
Type of Intersection this IntersectionIterator points to. | |
typedef GridImp::template Codim< 1 >::Geometry | Geometry |
Codim 1 geometry returned by intersectionGlobal(). | |
typedef GridImp::template Codim< 1 >::LocalGeometry | LocalGeometry |
Codim 1 geometry returned by intersectionSelfLocal() and intersectionNeighborLocal(). | |
typedef GridImp::ctype | ctype |
define type used for coordinates in grid module | |
Public Member Functions | |
IntersectionIterator & | operator++ () |
Preincrement operator. Proceed to next intersection. | |
Dereferencing | |
const Intersection & | operator* () const |
Dereferencing operator. | |
const Intersection * | operator-> () const |
Pointer operator. | |
Compare methods | |
bool | operator== (const IntersectionIterator &rhs) const |
Checks for equality. Only Iterators pointing to the same intersection from the same Entity are equal. Pointing to the same intersection from neighbor is unequal as inside and outside are permuted. | |
bool | operator!= (const IntersectionIterator &rhs) const |
Checks for inequality. Only Iterators pointing to the same intersection from the same Entity are equal. Pointing to the same intersection from neighbor is unequal as inside and outside are permuted. | |
Query methods | |
bool | boundary () const |
return true if intersection is with interior or exterior boundary (see the cases above) | |
int | boundaryId () const |
Identifier for boundary segment from macro grid. | |
bool | neighbor () const |
return true if intersection is shared with another element. | |
EntityPointer | inside () const |
return EntityPointer to the Entity on the inside of this intersection. That is the Entity where we started this Iterator. | |
EntityPointer | outside () const |
return EntityPointer to the Entity on the outside of this intersection. That is the neighboring Entity. | |
const LocalGeometry & | intersectionSelfLocal () const |
geometrical information about this intersection in local coordinates of the inside() entity. | |
const LocalGeometry & | intersectionNeighborLocal () const |
geometrical information about this intersection in local coordinates of the outside() entity. | |
const Geometry & | intersectionGlobal () const |
geometrical information about this intersection in global coordinates. | |
int | numberInSelf () const |
Local number of codim 1 entity in the inside() Entity where intersection is contained in. | |
int | numberInNeighbor () const |
Local number of codim 1 entity in outside() Entity where intersection is contained in. | |
FieldVector< ctype, dimworld > | outerNormal (const FieldVector< ctype, dim-1 > &local) const |
Return an outer normal (length not necessarily 1). | |
FieldVector< ctype, dimworld > | integrationOuterNormal (const FieldVector< ctype, dim-1 > &local) const |
return outer normal scaled with the integration element | |
FieldVector< ctype, dimworld > | unitOuterNormal (const FieldVector< ctype, dim-1 > &local) const |
Return unit outer normal (length == 1). | |
Implementor interface | |
bool | equals (const IntersectionIterator &rhs) const |
forward equality check to realIterator | |
IntersectionIterator (const IntersectionIteratorImp< const GridImp > &i) | |
IntersectionIterator (const IntersectionIterator &i) | |
Protected Member Functions | |
ImplementationType & | getRealImp () |
return reference to the real implementation | |
const ImplementationType & | getRealImp () const |
return reference to the real implementation |
anonymous enum |
anonymous enum |
Dune::IntersectionIterator< GridImp, IntersectionIteratorImp, IntersectionImp >::IntersectionIterator | ( | const IntersectionIteratorImp< const GridImp > & | i | ) | [inline] |
Copy Constructor from IntersectionIteratorImp
Dune::IntersectionIterator< GridImp, IntersectionIteratorImp, IntersectionImp >::IntersectionIterator | ( | const IntersectionIterator< GridImp, IntersectionIteratorImp, IntersectionImp > & | i | ) | [inline] |
Copy constructor
bool Dune::IntersectionIterator< GridImp, IntersectionIteratorImp, IntersectionImp >::boundary | ( | ) | const [inline] |
int Dune::IntersectionIterator< GridImp, IntersectionIteratorImp, IntersectionImp >::boundaryId | ( | ) | const [inline] |
Identifier for boundary segment from macro grid.
One can attach a boundary Id to a boundary segment on the macro grid. This Id will also be used for all fragments of these boundary segments.
The numbering is defined as:
The way the Identifiers are attached to the grid may differ between the different grid implementations.
bool Dune::IntersectionIterator< GridImp, IntersectionIteratorImp, IntersectionImp >::neighbor | ( | ) | const [inline] |
EntityPointer Dune::IntersectionIterator< GridImp, IntersectionIteratorImp, IntersectionImp >::inside | ( | ) | const [inline] |
return EntityPointer to the Entity on the inside of this intersection. That is the Entity where we started this Iterator.
EntityPointer Dune::IntersectionIterator< GridImp, IntersectionIteratorImp, IntersectionImp >::outside | ( | ) | const [inline] |
return EntityPointer to the Entity on the outside of this intersection. That is the neighboring Entity.
const LocalGeometry& Dune::IntersectionIterator< GridImp, IntersectionIteratorImp, IntersectionImp >::intersectionSelfLocal | ( | ) | const [inline] |
const LocalGeometry& Dune::IntersectionIterator< GridImp, IntersectionIteratorImp, IntersectionImp >::intersectionNeighborLocal | ( | ) | const [inline] |
const Geometry& Dune::IntersectionIterator< GridImp, IntersectionIteratorImp, IntersectionImp >::intersectionGlobal | ( | ) | const [inline] |
geometrical information about this intersection in global coordinates.
This method returns a Geometry object that provides a mapping from local coordinates of the intersection to global (world) coordinates.
int Dune::IntersectionIterator< GridImp, IntersectionIteratorImp, IntersectionImp >::numberInSelf | ( | ) | const [inline] |
int Dune::IntersectionIterator< GridImp, IntersectionIteratorImp, IntersectionImp >::numberInNeighbor | ( | ) | const [inline] |
FieldVector<ctype, dimworld> Dune::IntersectionIterator< GridImp, IntersectionIteratorImp, IntersectionImp >::outerNormal | ( | const FieldVector< ctype, dim-1 > & | local | ) | const [inline] |
Return an outer normal (length not necessarily 1).
The returned vector may depend on local position within the intersection.
FieldVector<ctype, dimworld> Dune::IntersectionIterator< GridImp, IntersectionIteratorImp, IntersectionImp >::integrationOuterNormal | ( | const FieldVector< ctype, dim-1 > & | local | ) | const [inline] |
return outer normal scaled with the integration element
Return an outer normal (length not necessarily 1).
The returned vector may depend on local position within the intersection.
The normal is scaled with the integration element of the intersection. This method is redundant but it may be more efficent to use this function rather than computing the integration element via intersectionGlobal().FieldVector<ctype, dimworld> Dune::IntersectionIterator< GridImp, IntersectionIteratorImp, IntersectionImp >::unitOuterNormal | ( | const FieldVector< ctype, dim-1 > & | local | ) | const [inline] |
Return unit outer normal (length == 1).
The returned vector may depend on the local position within the intersection. It is scaled to have unit length.