Dune::IntersectionIterator< GridImp, IntersectionIteratorImp > Class Template Reference
[IntersectionIterator]

#include <intersectioniterator.hh>

List of all members.


Detailed Description

template<class GridImp, template< class > class IntersectionIteratorImp>
class Dune::IntersectionIterator< GridImp, IntersectionIteratorImp >

Mesh entities of codimension 0 ("elements") allow to visit all intersections with "neighboring" elements and with the domain boundary.

Template parameters are:

Warning:
The name IntersectionIterator may be somewhat misleading. This class has neither an operator* nor an operator->. It iterates over codimension 1 intersections with other entities and the (sub-)domain boundary.

The number of neigbors may be different from the number of faces/edges of an element!

Overview

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.

Engine Concept

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.

Intersections, leaf grid and level grid

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.

islocalref.png

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.

Methods neighbor and boundary

The following holds for both the level and the leaf intersection iterator: The intersection iterator is started on a codimension 0 entity of the grid. If this entity belongs to the interior or the overlap region (see. ???) then the union of all intersections is identical to the boundary of the entity. On ghost elements the iterator only stops on the border of the domain, i.e., only on the intersections with entities in the interior region. Depending on the boolean values returned by the methods boundary() and neighbor() one can detect the position of an intersection relative to the boundary. In general the method boundary() returns true if and only if the intersection is part of the physical boundary of the domain. The method neighbor() returns true only if the method outside() has a well defined return value.

The following cases are possible if the intersection iterator is started on an entity in the interior or overlap region. More details are given below:

intersectionneighbor()boundary()outside()
1with inner, overlap
or ghost entity
truefalse the neighbor entity
2on domain boundary falsetrueundefined
3on periodic boundary truetrueGhost-/Overlap cell (with transformed geometry)
4on processor boundary false if grid has no ghosts
true otherwise
false ghost entity (if it exists)

  1. Inner Intersections:
    The type of the neighboring entity can be determined through methods defined on the outside entity.
  2. Handling physical boundaries: Different types of physical boundaries can be modeled using either the global coordinates of the intersection or by using the boundaryID method. On some grids (AluGrid, AlbertaGrid) this method returns an integer value which can be individually assigned to each boundary intersection of the macro grid and which is prolonged to higher levels during grid refinement. A more general concept will be included in latter releases along the following guidelines:
    • We require differently constructed geometries outside the domain
    • The kind of construction depends on the discrete problem
    • Therefor these constructions can't be part of the Grid interface
    • Utility classes are required to do this construction
    • The utility classes must be parameterized with the intersection (in our case the IntersectionIterator)
    • The utility classes return a suitable transformation of the inner() entitys geometry (with respect to the intersection), e.g., reflection at the intersection point reflection reflection combined with translation...
  3. Handling periodic boundaries:
    • The IntersectionIterator stops at periodic boundaries
    • periodic grids are handled in correspondence to parallel grids
    • At the periodic boundary one can adjust an overlap- or ghost-layer.
    • outer() returns a ghost or overlap cell (for ghost and overlap look into the documentation of the parallel grid interface)
    • outer() cell has a periodically transformed geometry (so that one does not see a jump or something like that)
    • outer() cell has its own index
    • outer() cell has the same id as the corresponding "original" cell
  4. Processor boundaries:
    At processor boundaries, i.e. when an element has an intersection with another element in the sequential grid but this element is only stored in other processors the intersection iterator stops but neither neighbor() nor boundary() are true.

    Geometry of an intersection

The method intersectionGlobal returns a geometry mapping the intersection as a codim one structure to global coordinates. The methods intersectionSelfLocal and intersectionNeighborLocal return geometries mapping the intersection into the reference elements of the originating entity and the neighboring entity, respectively. The numberInSelf and numberInNeighbor methods return the codim one subentities which contain the intersection.


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 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 IntersectionIterator
< GridImp,
IntersectionIteratorImp > 
Intersection
typedef GridImp::ctype ctype
 define type used for coordinates in grid module

Public Member Functions

const Intersectionoperator* () const
 Dereferencing operator.
const Intersectionoperator-> () const
 Pointer operator.
IntersectionIteratoroperator++ ()
 Preincrement operator. Proceed to next intersection.
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 LocalGeometryintersectionSelfLocal () const
 geometrical information about this intersection in local coordinates of the inside() entity.
const LocalGeometryintersectionNeighborLocal () const
 geometrical information about this intersection in local coordinates of the outside() entity.
const GeometryintersectionGlobal () 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).
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.
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

Member Typedef Documentation

template<class GridImp, template< class > class IntersectionIteratorImp>
typedef IntersectionIterator<GridImp, IntersectionIteratorImp> Dune::IntersectionIterator< GridImp, IntersectionIteratorImp >::Intersection

forward compatibility to new Dune Intersection's


Member Enumeration Documentation

template<class GridImp, template< class > class IntersectionIteratorImp>
anonymous enum

Export grid dimension.

Enumerator:
dimension  grid dimension

template<class GridImp, template< class > class IntersectionIteratorImp>
anonymous enum

Export dimension of world.

Enumerator:
dimensionworld  dimension of world


Constructor & Destructor Documentation

template<class GridImp, template< class > class IntersectionIteratorImp>
Dune::IntersectionIterator< GridImp, IntersectionIteratorImp >::IntersectionIterator ( const IntersectionIteratorImp< const GridImp > &  i  )  [inline]

Copy Constructor from IntersectionIteratorImp

template<class GridImp, template< class > class IntersectionIteratorImp>
Dune::IntersectionIterator< GridImp, IntersectionIteratorImp >::IntersectionIterator ( const IntersectionIterator< GridImp, IntersectionIteratorImp > &  i  )  [inline]

Copy constructor


Member Function Documentation

template<class GridImp, template< class > class IntersectionIteratorImp>
int Dune::IntersectionIterator< GridImp, IntersectionIteratorImp >::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:

  • Id==0 for all intersections without boundary()==false
  • Id>=0 for all intersections without boundary()==true

The way the Identifiers are attached to the grid may differ between the different grid implementations.

template<class GridImp, template< class > class IntersectionIteratorImp>
EntityPointer Dune::IntersectionIterator< GridImp, IntersectionIteratorImp >::outside (  )  const [inline]

return EntityPointer to the Entity on the outside of this intersection. That is the neighboring Entity.

Warning:
Don't call this method if there is no neighboring Entity (neighbor() returns false). In this case the result is undefined.

template<class GridImp, template< class > class IntersectionIteratorImp>
const LocalGeometry& Dune::IntersectionIterator< GridImp, IntersectionIteratorImp >::intersectionSelfLocal (  )  const [inline]

geometrical information about this intersection in local coordinates of the inside() entity.

This method returns a Geometry object that provides a mapping from local coordinates of the intersection to local coordinates of the inside() entity.

template<class GridImp, template< class > class IntersectionIteratorImp>
const LocalGeometry& Dune::IntersectionIterator< GridImp, IntersectionIteratorImp >::intersectionNeighborLocal (  )  const [inline]

geometrical information about this intersection in local coordinates of the outside() entity.

This method returns a Geometry object that provides a mapping from local coordinates of the intersection to local coordinates of the outside() entity.

template<class GridImp, template< class > class IntersectionIteratorImp>
const Geometry& Dune::IntersectionIterator< GridImp, IntersectionIteratorImp >::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.

template<class GridImp, template< class > class IntersectionIteratorImp>
FieldVector<ctype, dimworld> Dune::IntersectionIterator< GridImp, IntersectionIteratorImp >::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.

template<class GridImp, template< class > class IntersectionIteratorImp>
FieldVector<ctype, dimworld> Dune::IntersectionIterator< GridImp, IntersectionIteratorImp >::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().

template<class GridImp, template< class > class IntersectionIteratorImp>
FieldVector<ctype, dimworld> Dune::IntersectionIterator< GridImp, IntersectionIteratorImp >::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.


The documentation for this class was generated from the following file:

Generated on 6 Nov 2008 with Doxygen (ver 1.5.6) [logfile].