[Dune] Dune::Intersection
Christian Engwer
christi at uni-hd.de
Fri Mar 28 15:03:16 CET 2008
Dear all,
as promised on the last developer meeting, I implemented a Dune::Intersection.
The Dune::IntersectionIteratorDefault offers an implementation of
dereference, so that the existing implementations for
IntersectionIterator can also be used for both IntersectionIterator
and Intersection.
The Intersection can parametrized with what was the
IntersectionIteratorImp. This allows an easy transition. The Grid
Traits now also export the types LeafIntersection and
LevelIntersection and take two new parameters for the
LeafIntersectionImp and LevelIntersectionImp. I updated all grids
using the Dune::GridTraits class, Alberta and Alugrid use their own
Traits. I didn't want to break anything. This means that someone from
Freiburg has to add these typedefs to the Traits.
I was also wandering what the IntersectionIteratorWrapper is. I have a
rough idea of what it is used for, but there is no documentation on
how to use it. Therefor I'm not sure that nothing got broken for this
wrapper.
You can find the upcoming patch attached. I will commit later today.
Christian
-------------- next part --------------
Index: grid/test/checkintersectionit.cc
===================================================================
--- grid/test/checkintersectionit.cc (Revision 4067)
+++ grid/test/checkintersectionit.cc (Arbeitskopie)
@@ -98,16 +98,22 @@
// /////////////////////////////////////////////////////////
// Check the types defined by the iterator
// /////////////////////////////////////////////////////////
- IsTrue< is_same<
+ dune_static_assert((is_same<
typename IntersectionIterator::ctype,
- typename GridType::ctype>::value == true >::yes();
+ typename GridType::ctype>::value),
+ "IntersectionIterator has wrong ctype");
- IsTrue<static_cast<int>(IntersectionIterator::dimension)
- == static_cast<int>(GridType::dimension)>::yes();
+ dune_static_assert((is_same<
+ typename IntersectionIterator::Intersection,
+ typename GridPartType::IntersectionType>::value),
+ "IntersectionIterator has wrong Intersection type");
- IsTrue<static_cast<int>(IntersectionIterator::dimensionworld)
- == static_cast<int>(GridType::dimensionworld)>::yes();
+ dune_static_assert((static_cast<int>(IntersectionIterator::dimension)
+ == static_cast<int>(GridType::dimension)),"IntersectionIterator has wrong dimension");
+ dune_static_assert((static_cast<int>(IntersectionIterator::dimensionworld)
+ == static_cast<int>(GridType::dimensionworld)),"IntersectionIterator has wrong dimensionworld");
+
IntersectionIterator iIt = gridPart.ibegin(*eIt);
IntersectionIterator iEndIt = gridPart.iend(*eIt);
@@ -121,25 +127,25 @@
// //////////////////////////////////////////////////////////////////////
const int interDim = IntersectionIterator::LocalGeometry::mydimension;
const QuadratureRule<double, interDim>& quad
- = QuadratureRules<double, interDim>::rule(iIt.intersectionSelfLocal().type(), interDim);
+ = QuadratureRules<double, interDim>::rule(iIt->intersectionSelfLocal().type(), interDim);
for (size_t i=0; i<quad.size(); i++)
- sumNormal.axpy(quad[i].weight(), iIt.integrationOuterNormal(quad[i].position()));
+ sumNormal.axpy(quad[i].weight(), iIt->integrationOuterNormal(quad[i].position()));
typedef typename IntersectionIterator::Entity EntityType;
typedef typename EntityType::EntityPointer EntityPointer;
- assert(eIt == iIt.inside());
+ assert(eIt == iIt->inside());
// check that boundary id has positive value
- if( iIt.boundary() )
+ if( iIt->boundary() )
{
// entity has boundary intersections
hasBoundaryIntersection = true;
- if( iIt.boundaryId() <= 0 )
+ if( iIt->boundaryId() <= 0 )
{
- DUNE_THROW(GridError, "boundary id has non-positive value (" << iIt.boundaryId() << ") !");
+ DUNE_THROW(GridError, "boundary id has non-positive value (" << iIt->boundaryId() << ") !");
}
}
@@ -147,9 +153,9 @@
// Check whether the 'has-intersection-with'-relation is symmetric
// //////////////////////////////////////////////////////////////////////
- if (iIt.neighbor() && checkOutside )
+ if (iIt->neighbor() && checkOutside )
{
- EntityPointer outside = iIt.outside();
+ EntityPointer outside = iIt->outside();
bool insideFound = false;
IntersectionIterator outsideIIt = gridPart.ibegin(*outside);
@@ -157,9 +163,9 @@
for (; outsideIIt!=outsideIEndIt; ++outsideIIt) {
- if (outsideIIt.neighbor() && outsideIIt.outside() == iIt.inside()) {
+ if (outsideIIt->neighbor() && outsideIIt->outside() == iIt->inside()) {
- if (outsideIIt.numberInSelf() != iIt.numberInNeighbor())
+ if (outsideIIt->numberInSelf() != iIt->numberInNeighbor())
DUNE_THROW(GridError, "outside()->outside() == inside(), but with incorrect numbering!");
else
insideFound = true;
@@ -186,11 +192,11 @@
// Check the consistency of numberInSelf, numberInNeighbor
// and the indices of the subface between.
// /////////////////////////////////////////////////////////////
- if ( GridPartType::conforming && iIt.neighbor() )
+ if ( GridPartType::conforming && iIt->neighbor() )
{
- EntityPointer outside = iIt.outside();
- int numberInSelf = iIt.numberInSelf();
- int numberInNeighbor = iIt.numberInNeighbor();
+ EntityPointer outside = iIt->outside();
+ int numberInSelf = iIt->numberInSelf();
+ int numberInNeighbor = iIt->numberInNeighbor();
assert(indexSet.template subIndex<1>(*eIt, numberInSelf)
== indexSet.template subIndex<1>(*outside, numberInNeighbor));
@@ -207,7 +213,7 @@
// Check the geometry returned by intersectionGlobal()
// //////////////////////////////////////////////////////////
typedef typename IntersectionIterator::Geometry Geometry;
- const Geometry& intersectionGlobal = iIt.intersectionGlobal();
+ const Geometry& intersectionGlobal = iIt->intersectionGlobal();
checkGeometry(intersectionGlobal);
@@ -215,7 +221,7 @@
// Check the geometry returned by intersectionSelfLocal()
// //////////////////////////////////////////////////////////
- const typename IntersectionIterator::LocalGeometry& intersectionSelfLocal = iIt.intersectionSelfLocal();
+ const typename IntersectionIterator::LocalGeometry& intersectionSelfLocal = iIt->intersectionSelfLocal();
checkGeometry(intersectionSelfLocal);
// Check the consistency of intersectionSelfLocal() and intersectionGlobal
@@ -228,7 +234,7 @@
{
// check integrationOuterNormal
double det = intersectionGlobal.integrationElement(quad[i].position());
- det -= iIt.integrationOuterNormal(quad[i].position()).two_norm();
+ det -= iIt->integrationOuterNormal(quad[i].position()).two_norm();
if( std::abs( det ) > 1e-8 )
{
DUNE_THROW(GridError, "integrationElement and length of integrationOuterNormal do no match!");
@@ -246,10 +252,10 @@
// Check the geometry returned by intersectionNeighborLocal()
// ////////////////////////////////////////////////////////////////
- if (iIt.neighbor() )
+ if (iIt->neighbor() )
{
- const typename IntersectionIterator::LocalGeometry& intersectionNeighborLocal = iIt.intersectionNeighborLocal();
+ const typename IntersectionIterator::LocalGeometry& intersectionNeighborLocal = iIt->intersectionNeighborLocal();
checkGeometry(intersectionNeighborLocal);
@@ -265,7 +271,7 @@
{
FieldVector<double,dimworld> globalPos = intersectionGlobal.global(quad[i].position());
- FieldVector<double,dimworld> localPos = iIt.outside()->geometry().global(intersectionNeighborLocal.global(quad[i].position()));
+ FieldVector<double,dimworld> localPos = iIt->outside()->geometry().global(intersectionNeighborLocal.global(quad[i].position()));
if ( (globalPos - localPos).infinity_norm() > 1e-6)
DUNE_THROW(GridError, "global( intersectionNeighborLocal(global() ) is not the same as intersectionGlobal.global() at " << quad[i].position() << "!");
Index: grid/test/checkgeometryinfather.cc
===================================================================
--- grid/test/checkgeometryinfather.cc (Revision 4067)
+++ grid/test/checkgeometryinfather.cc (Arbeitskopie)
@@ -105,21 +105,21 @@
// Check for types and constants
// //////////////////////////////////////////////////////
- IsTrue< is_same<
+ dune_static_assert((is_same<
typename Geometry::ctype,
- typename GridType::ctype>::value == true >::yes();
+ typename GridType::ctype>::value == true),"Geometry has wrong ctype");
- IsTrue<static_cast<int>(Geometry::dimension)
- == static_cast<int>(GridType::dimension)>::yes();
+ dune_static_assert((static_cast<int>(Geometry::dimension)
+ == static_cast<int>(GridType::dimension)),"Geometry has wrong dimension");
- IsTrue<static_cast<int>(Geometry::mydimension)
- == static_cast<int>(GridType::dimension)>::yes();
+ dune_static_assert((static_cast<int>(Geometry::mydimension)
+ == static_cast<int>(GridType::dimension)),"Geometry has wrong mydimension");
- IsTrue<static_cast<int>(Geometry::coorddimension)
- == static_cast<int>(GridType::dimensionworld)>::yes();
+ dune_static_assert((static_cast<int>(Geometry::coorddimension)
+ == static_cast<int>(GridType::dimensionworld)),"Geometry has wrong coorddimension");
- IsTrue<static_cast<int>(Geometry::dimensionworld)
- == static_cast<int>(GridType::dimensionworld)>::yes();
+ dune_static_assert((static_cast<int>(Geometry::dimensionworld)
+ == static_cast<int>(GridType::dimensionworld)),"Geometry has wrong dimensionworld");
// ///////////////////////////////////////////////////////
// Check the different methods
Index: grid/test/checkindexset.cc
===================================================================
--- grid/test/checkindexset.cc (Revision 4067)
+++ grid/test/checkindexset.cc (Arbeitskopie)
@@ -466,10 +466,10 @@
IntersectionIterator endnit = it->ilevelend();
for(IntersectionIterator nit = it->ilevelbegin(); nit != endnit; ++nit)
{
- if(nit.neighbor())
+ if(nit->neighbor())
{
typedef typename GridType :: template Codim<0> :: EntityPointer EnPointer;
- EnPointer ep = nit.outside();
+ EnPointer ep = nit->outside();
checkSubEntity<codim> (grid, *ep, lset, sout,
subEntities, vertices, vertexCoordsMap);
@@ -493,10 +493,10 @@
IntersectionIterator endnit = it->ileafend();
for(IntersectionIterator nit = it->ileafbegin(); nit != endnit; ++nit)
{
- if(nit.neighbor())
+ if(nit->neighbor())
{
typedef typename GridType :: template Codim<0> :: EntityPointer EnPointer;
- EnPointer ep = nit.outside();
+ EnPointer ep = nit->outside();
checkSubEntity<codim> (grid, *ep, lset, sout,
subEntities, vertices, vertexCoordsMap);
Index: grid/test/gridcheck.cc
===================================================================
--- grid/test/gridcheck.cc (Revision 4067)
+++ grid/test/gridcheck.cc (Arbeitskopie)
@@ -519,8 +519,8 @@
}
// recursively check entity-interface
// ... we only allow grids with codim 0 zero entites
- IsTrue<Dune::Capabilities::hasEntity<Grid, 0>::v>::yes();
- IsTrue<Dune::Capabilities::hasEntity<const Grid, 0>::v>::yes();
+ dune_static_assert((Dune::Capabilities::hasEntity<Grid, 0>::v),"Grid must have codim 0 entities");
+ dune_static_assert((Dune::Capabilities::hasEntity<const Grid, 0>::v),"Grid must have codim 0 entities");
//EntityInterface<Grid, 0, Grid::dimension,
// Dune::Capabilities::hasEntity<Grid, 0>::v >();
@@ -691,10 +691,10 @@
IntersectionIterator endit = e->ilevelend();
IntersectionIterator it = e->ilevelbegin();
// state
- it.boundary();
- it.neighbor();
+ it->boundary();
+ it->neighbor();
// id of boundary segment
- it.boundaryId();
+ it->boundaryId();
// check id
//assert(globalid.id(*e) >= 0);
assert(it != endit);
@@ -714,46 +714,46 @@
for(; it != endit; ++it)
{
// mark visited face
- visited[it.numberInSelf()] = true;
+ visited[it->numberInSelf()] = true;
// check id
- assert(globalid.id(*(it.inside())) ==
+ assert(globalid.id(*(it->inside())) ==
globalid.id(*e));
// numbering
- int num = it.numberInSelf();
+ int num = it->numberInSelf();
assert( num >= 0 && num < e->template count<1> () );
- if(it.neighbor())
+ if(it->neighbor())
{
// geometry
- it.intersectionNeighborLocal();
+ it->intersectionNeighborLocal();
// numbering
- num = it.numberInNeighbor();
- assert( num >= 0 && num < it.outside()->template count<1> () );
+ num = it->numberInNeighbor();
+ assert( num >= 0 && num < it->outside()->template count<1> () );
}
// geometry
- it.intersectionSelfLocal();
- it.intersectionGlobal();
+ it->intersectionSelfLocal();
+ it->intersectionGlobal();
// normal vectors
Dune::FieldVector<ct, dim-1> v(0);
- it.outerNormal(v);
- it.integrationOuterNormal(v);
- it.unitOuterNormal(v);
+ it->outerNormal(v);
+ it->integrationOuterNormal(v);
+ it->unitOuterNormal(v);
// search neighbouring cell
- if (it.neighbor())
+ if (it->neighbor() && it->outside()->partitionType() != Dune::InteriorEntity)
{
- //assert(globalid.id(*(it.outside())) >= 0);
- assert(globalid.id(*(it.outside())) !=
+ //assert(globalid.id(*(it->outside())) >= 0);
+ assert(globalid.id(*(it->outside())) !=
globalid.id(*e));
LevelIterator n = g.template lbegin<0>(e->level());
LevelIterator nend = g.template lend<0> (e->level());
- while (n != it.outside() && n != nend)
+ while (n != it->outside() && n != nend)
{
- assert(globalid.id(*(it.outside())) !=
+ assert(globalid.id(*(it->outside())) !=
globalid.id(*n));
++n;
}
@@ -922,15 +922,15 @@
i == l2; \
i == h2; \
i == L2; \
- if (i2 != l2->ileafend()) i == i2.inside(); \
- if (i2 != l2->ileafend() && i2.neighbor()) i == i2.outside(); \
+ if (i2 != l2->ileafend()) i == i2->inside(); \
+ if (i2 != l2->ileafend() && i2->neighbor()) i == i2->outside(); \
}
TestEquals(e1);
TestEquals(l1);
TestEquals(h1);
TestEquals(L1);
- if (i1 != l1->ileafend()) TestEquals(i1.inside());
- if (i1 != l1->ileafend() && i1.neighbor()) TestEquals(i1.outside());
+ if (i1 != l1->ileafend()) TestEquals(i1->inside());
+ if (i1 != l1->ileafend() && i1->neighbor()) TestEquals(i1->outside());
}
template <class Grid>
Index: grid/uggrid.hh
===================================================================
--- grid/uggrid.hh (Revision 4067)
+++ grid/uggrid.hh (Arbeitskopie)
@@ -102,6 +102,8 @@
UGGridEntity,
UGGridEntityPointer,
UGGridLevelIterator,
+ UGGridLeafIntersectionIterator, // leaf intersection
+ UGGridLevelIntersectionIterator, // level intersection
UGGridLeafIntersectionIterator, // leaf intersection iterartor
UGGridLevelIntersectionIterator, // level intersection iterartor
UGGridHierarchicIterator,
@@ -512,7 +514,7 @@
const FieldVector<double, dim>& pos);
/** \brief For a point on the grid boundary return its position on the domain boundary */
- FieldVector<ctype,dim> getBoundaryPosition(const IntersectionIterator<const UGGrid<dim>, UGGridLevelIntersectionIterator>& iIt,
+ FieldVector<ctype,dim> getBoundaryPosition(const typename Traits::LevelIntersectionIterator& iIt,
const FieldVector<ctype,dim-1>& localPos) const;
/** \brief Does uniform refinement
Index: grid/yaspgrid.hh
===================================================================
--- grid/yaspgrid.hh (Revision 4070)
+++ grid/yaspgrid.hh (Arbeitskopie)
@@ -2102,6 +2102,8 @@
typedef GridTraits<dim,dimworld,Dune::YaspGrid<dim,dimworld>,
YaspGeometry,YaspEntity,
YaspEntityPointer,YaspLevelIterator,
+ YaspIntersectionIterator, // leaf intersection
+ YaspIntersectionIterator, // level intersection
YaspIntersectionIterator, // leaf intersection iter
YaspIntersectionIterator, // level intersection iter
YaspHierarchicIterator,
Index: grid/uggrid/ugintersectionit.hh
===================================================================
--- grid/uggrid/ugintersectionit.hh (Revision 4067)
+++ grid/uggrid/ugintersectionit.hh (Arbeitskopie)
@@ -21,7 +21,8 @@
of an element!
*/
template<class GridImp>
-class UGGridLevelIntersectionIterator
+class UGGridLevelIntersectionIterator :
+ public IntersectionIteratorDefaultImplementation<GridImp,UGGridLevelIntersectionIterator>
{
enum {dim=GridImp::dimension};
@@ -164,7 +165,8 @@
template<class GridImp>
-class UGGridLeafIntersectionIterator
+class UGGridLeafIntersectionIterator :
+ public IntersectionIteratorDefaultImplementation<GridImp,UGGridLeafIntersectionIterator>
{
enum {dim=GridImp::dimension};
Index: grid/uggrid/uggrid.cc
===================================================================
--- grid/uggrid/uggrid.cc (Revision 4067)
+++ grid/uggrid/uggrid.cc (Arbeitskopie)
@@ -1212,7 +1212,7 @@
}
template <int dim>
-Dune::FieldVector<typename UGGrid<dim>::ctype,dim> Dune::UGGrid<dim>::getBoundaryPosition(const IntersectionIterator<const UGGrid<dim>, UGGridLevelIntersectionIterator>& iIt,
+Dune::FieldVector<typename UGGrid<dim>::ctype,dim> Dune::UGGrid<dim>::getBoundaryPosition(const typename Traits::LevelIntersectionIterator& iIt,
const FieldVector<ctype,dim-1>& localPos) const
{
if (!iIt.boundary())
Index: grid/common/intersection.hh
===================================================================
--- grid/common/intersection.hh (Revision 0)
+++ grid/common/intersection.hh (Revision 0)
@@ -0,0 +1,796 @@
+#ifndef DUNE_GRID_INTERSECTION_HH
+#define DUNE_GRID_INTERSECTION_HH
+
+namespace Dune
+{
+
+/** \brief %Intersection of a mesh entities of codimension 0 ("elements")
+ with a "neighboring" elements and with the domain
+ boundary.
+
+ Template parameters are:
+
+ - <tt>GridImp</tt> Type that is a model of Dune::Grid
+ - <tt>IntersectionImp</tt> Class template that is a model of
+ Dune::Intersection
+
+ @warning The name Intersection 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.
+
+ @warning The number of neigbors may be different from the number of
+ faces/edges of an element!
+
+ <h2>Overview</h2>
+
+ Intersections are codimension 1 objects. These
+ intersections are accessed via an Intersection. 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
+ Intersection describes these intersections implicitly.
+
+ <H2>Engine Concept</H2>
+
+ The Intersection class template wraps an object of type IntersectionImp
+ and forwards all member
+ function calls to corresponding members of this class. In that sense Intersection
+ defines the interface and IntersectionImp supplies the implementation.
+
+ <h2>Methods neighbor and boundary </h2>
+
+ The following holds for both the level and the leaf intersection
+ :
+ The %intersection 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 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 is
+ started on an entity in the interior or overlap region. More
+ details are given below:
+
+ <table>
+ <tr>
+ <td></td><td>intersection</td><td>neighbor()</td><td>boundary()</td><td>outside()</td></tr>
+ <tr>
+ <td>1</td><td>with inner, overlap <br>
+ or ghost entity</td>
+ <td>true</td><td>false</td>
+ <td>the neighbor entity</td></tr>
+ <tr>
+ <td>2</td><td>on domain boundary</td>
+ <td>false</td><td>true</td><td><em>undefined</em></td></tr>
+ <tr>
+ <td>3</td><td>on periodic boundary</td>
+ <td>true</td><td>true</td><td>Ghost-/Overlap cell at br (with transformed geometry)</td></tr>
+ <tr>
+ <td>4</td><td>on processor boundary</td>
+ <td>false <em>if grid has no ghosts</em><br>true <em>otherwise</em></td><td>false </td>
+ <td>ghost entity <em>(if it exists)</em></td></tr>
+ </table>
+
+ -# <b> Inner Intersections: </b> \n
+ The type of the neighboring entity can be determined through
+ methods defined on the outside entity.
+ -# <b> Handling physical boundaries: </b>
+ 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. \br
+ 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 Intersection)
+ - 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...
+ .
+ -# <b> Handling periodic boundaries: </b>
+ - The Intersection 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
+ -# <b> Processor boundaries: </b> \n
+ 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 stops but neither
+ neighbor() nor boundary()
+ are true.
+ .
+
+
+ <h2>Geometry of an intersection</h2>
+
+ 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.
+
+
+ @ingroup GIIntersectionIterator
+ */
+template<class GridImp, template<class> class IntersectionImp>
+class Intersection
+{
+ IntersectionImp<const GridImp> real;
+
+ enum { dim=GridImp::dimension };
+ enum { dimworld=GridImp::dimensionworld };
+
+public:
+
+ // type of real implementation
+ typedef IntersectionImp<const GridImp> ImplementationType;
+
+ /** \brief Type of entity that this Intersection belongs to */
+ typedef typename GridImp::template Codim<0>::Entity Entity;
+
+ /** \brief Pointer to the type of entities that this Intersection belongs to */
+ typedef typename GridImp::template Codim<0>::EntityPointer EntityPointer;
+
+ /** \brief Codim 1 geometry returned by intersectionGlobal() */
+ typedef typename GridImp::template Codim<1>::Geometry Geometry;
+
+ /** \brief Codim 1 geometry returned by intersectionSelfLocal()
+ and intersectionNeighborLocal() */
+ typedef typename GridImp::template Codim<1>::LocalGeometry LocalGeometry;
+
+ //! @brief Export grid dimension
+ enum { dimension=dim /*!< grid dimension */ };
+
+ //! @brief Export dimension of world
+ enum { dimensionworld=dimworld /*!< dimension of world */ };
+
+ //! define type used for coordinates in grid module
+ typedef typename GridImp::ctype ctype;
+
+ //! return true if intersection is with interior or exterior boundary (see the cases above)
+ bool boundary () const
+ {
+ return this->real.boundary();
+ }
+
+ /**
+ \brief 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.
+
+ \deprecated */
+ int boundaryId () const
+ {
+ return this->real.boundaryId();
+ }
+
+ //! @brief return true if intersection is shared with another element.
+ bool neighbor () const
+ {
+ return this->real.neighbor();
+ }
+
+ /*! @brief return EntityPointer to the Entity on the inside of this
+ intersection. That is the Entity where we started this .
+ */
+ EntityPointer inside() const
+ {
+ return this->real.inside();
+ }
+
+ /*! @brief 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.
+ */
+ EntityPointer outside() const
+ {
+ return this->real.outside();
+ }
+
+ /*! @brief 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.
+ */
+ const LocalGeometry& intersectionSelfLocal () const
+ {
+ return this->real.intersectionSelfLocal();
+ }
+ /*! @brief 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.
+ */
+ const LocalGeometry& intersectionNeighborLocal () const
+ {
+ return this->real.intersectionNeighborLocal();
+ }
+
+ /*! @brief 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.
+ */
+ const Geometry& intersectionGlobal () const
+ {
+ return this->real.intersectionGlobal();
+ }
+
+ //! Local number of codim 1 entity in the inside() Entity where intersection is contained in
+ int numberInSelf () const
+ {
+ return this->real.numberInSelf ();
+ }
+
+ //! Local number of codim 1 entity in outside() Entity where intersection is contained in
+ int numberInNeighbor () const
+ {
+ return this->real.numberInNeighbor ();
+ }
+
+ /*! @brief Return an outer normal (length not necessarily 1)
+
+ The returned vector may depend on local position within the intersection.
+ */
+ FieldVector<ctype, dimworld> outerNormal (const FieldVector<ctype, dim-1>& local) const
+ {
+ return this->real.outerNormal(local);
+ }
+
+ /*! @brief return outer normal scaled with the integration element
+ @copydoc outerNormal
+ 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> integrationOuterNormal (const FieldVector<ctype, dim-1>& local) const
+ {
+ return this->real.integrationOuterNormal(local);
+ }
+
+ /*! @brief 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.
+ */
+ FieldVector<ctype, dimworld> unitOuterNormal (const FieldVector<ctype, dim-1>& local) const
+ {
+ return this->real.unitOuterNormal(local);
+ }
+
+ /** @brief Checks for equality.
+ Only s 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 Intersection& rhs) const
+ {
+ return rhs.equals(*this);
+ }
+
+ /** @brief Checks for inequality.
+ Only s 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 Intersection& rhs) const
+ {
+ return ! rhs.equals(*this);
+ }
+
+
+ //===========================================================
+ /** @name Implementor interface
+ */
+ //@{
+ //===========================================================
+
+ /** @brief forward equality check to real */
+ bool equals(const Intersection& rhs) const
+ {
+ return this->real.equals(rhs.real);
+ }
+
+ /** Copy Constructor from IntersectionImp */
+ Intersection(const IntersectionImp<const GridImp> & i) :
+ real(i) {};
+
+ /** Copy constructor */
+ Intersection(const Intersection& i) :
+ real(i.real) {}
+ //@}
+
+ typedef typename remove_const<GridImp>::type mutableGridImp;
+protected:
+ // give the GridDefaultImplementation class access to the realImp
+ friend class GridDefaultImplementation<
+ GridImp::dimension, GridImp::dimensionworld,
+ typename GridImp::ctype,
+ typename GridImp::GridFamily> ;
+
+ //! return reference to the real implementation
+ ImplementationType & getRealImp() { return real; }
+ //! return reference to the real implementation
+ const ImplementationType & getRealImp() const { return real; }
+
+};
+
+//**********************************************************************
+/**
+ @brief Default Implementations of integrationOuterNormal and unitOuterNormal for IntersectionImp
+
+ @ingroup GridDevel
+*/
+template<class GridImp, template<class> class IntersectionImp>
+class IntersectionDefaultNormalVectors
+{
+ enum { dim=GridImp::dimension };
+ enum { dimworld=GridImp::dimensionworld };
+ typedef typename GridImp::ctype ct;
+public:
+ //! return unit outer normal, this should be dependent on
+ //! local coordinates for higher order boundary
+ //! the normal is scaled with the integration element of the intersection.
+ FieldVector<ct, dimworld> integrationOuterNormal (const FieldVector<ct, dim-1>& local) const
+ {
+ FieldVector<ct, dimworld> n = unitOuterNormal(local);
+ n *= asImp().intersectionGlobal().integrationElement(local);
+ return n;
+ }
+ //! return unit outer normal
+ FieldVector<ct, dimworld> unitOuterNormal (const FieldVector<ct, dim-1>& local) const
+ {
+ FieldVector<ct, dimworld> n = asImp().outerNormal(local);
+ n /= n.two_norm();
+ return n;
+ }
+
+private:
+ //! Barton-Nackman trick
+ IntersectionImp<GridImp>& asImp ()
+ {return static_cast<IntersectionImp<GridImp>&>(*this);}
+ const IntersectionImp<GridImp>& asImp () const
+ {return static_cast<const IntersectionImp<GridImp>&>(*this);}
+};
+
+}
+
+#endif // DUNE_GRID_INTERSECTION_HH
+#ifndef DUNE_GRID_INTERSECTION_HH
+#define DUNE_GRID_INTERSECTION_HH
+
+namespace Dune
+{
+
+/** \brief %Intersection of a mesh entities of codimension 0 ("elements")
+ with a "neighboring" elements and with the domain
+ boundary.
+
+ Template parameters are:
+
+ - <tt>GridImp</tt> Type that is a model of Dune::Grid
+ - <tt>IntersectionImp</tt> Class template that is a model of
+ Dune::Intersection
+
+ @warning The name Intersection 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.
+
+ @warning The number of neigbors may be different from the number of
+ faces/edges of an element!
+
+ <h2>Overview</h2>
+
+ Intersections are codimension 1 objects. These
+ intersections are accessed via an Intersection. 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
+ Intersection describes these intersections implicitly.
+
+ <H2>Engine Concept</H2>
+
+ The Intersection class template wraps an object of type IntersectionImp
+ and forwards all member
+ function calls to corresponding members of this class. In that sense Intersection
+ defines the interface and IntersectionImp supplies the implementation.
+
+ <h2>Methods neighbor and boundary </h2>
+
+ The following holds for both the level and the leaf intersection
+ :
+ The %intersection 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 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 is
+ started on an entity in the interior or overlap region. More
+ details are given below:
+
+ <table>
+ <tr>
+ <td></td><td>intersection</td><td>neighbor()</td><td>boundary()</td><td>outside()</td></tr>
+ <tr>
+ <td>1</td><td>with inner, overlap <br>
+ or ghost entity</td>
+ <td>true</td><td>false</td>
+ <td>the neighbor entity</td></tr>
+ <tr>
+ <td>2</td><td>on domain boundary</td>
+ <td>false</td><td>true</td><td><em>undefined</em></td></tr>
+ <tr>
+ <td>3</td><td>on periodic boundary</td>
+ <td>true</td><td>true</td><td>Ghost-/Overlap cell at br (with transformed geometry)</td></tr>
+ <tr>
+ <td>4</td><td>on processor boundary</td>
+ <td>false <em>if grid has no ghosts</em><br>true <em>otherwise</em></td><td>false </td>
+ <td>ghost entity <em>(if it exists)</em></td></tr>
+ </table>
+
+ -# <b> Inner Intersections: </b> \n
+ The type of the neighboring entity can be determined through
+ methods defined on the outside entity.
+ -# <b> Handling physical boundaries: </b>
+ 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. \br
+ 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 Intersection)
+ - 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...
+ .
+ -# <b> Handling periodic boundaries: </b>
+ - The Intersection 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
+ -# <b> Processor boundaries: </b> \n
+ 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 stops but neither
+ neighbor() nor boundary()
+ are true.
+ .
+
+
+ <h2>Geometry of an intersection</h2>
+
+ 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.
+
+
+ @ingroup GIIntersectionIterator
+ */
+template<class GridImp, template<class> class IntersectionImp>
+class Intersection
+{
+ IntersectionImp<const GridImp> real;
+
+ enum { dim=GridImp::dimension };
+ enum { dimworld=GridImp::dimensionworld };
+
+public:
+
+ // type of real implementation
+ typedef IntersectionImp<const GridImp> ImplementationType;
+
+ /** \brief Type of entity that this Intersection belongs to */
+ typedef typename GridImp::template Codim<0>::Entity Entity;
+
+ /** \brief Pointer to the type of entities that this Intersection belongs to */
+ typedef typename GridImp::template Codim<0>::EntityPointer EntityPointer;
+
+ /** \brief Codim 1 geometry returned by intersectionGlobal() */
+ typedef typename GridImp::template Codim<1>::Geometry Geometry;
+
+ /** \brief Codim 1 geometry returned by intersectionSelfLocal()
+ and intersectionNeighborLocal() */
+ typedef typename GridImp::template Codim<1>::LocalGeometry LocalGeometry;
+
+ //! @brief Export grid dimension
+ enum { dimension=dim /*!< grid dimension */ };
+
+ //! @brief Export dimension of world
+ enum { dimensionworld=dimworld /*!< dimension of world */ };
+
+ //! define type used for coordinates in grid module
+ typedef typename GridImp::ctype ctype;
+
+ //! return true if intersection is with interior or exterior boundary (see the cases above)
+ bool boundary () const
+ {
+ return this->real.boundary();
+ }
+
+ /**
+ \brief 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.
+
+ \deprecated */
+ int boundaryId () const
+ {
+ return this->real.boundaryId();
+ }
+
+ //! @brief return true if intersection is shared with another element.
+ bool neighbor () const
+ {
+ return this->real.neighbor();
+ }
+
+ /*! @brief return EntityPointer to the Entity on the inside of this
+ intersection. That is the Entity where we started this .
+ */
+ EntityPointer inside() const
+ {
+ return this->real.inside();
+ }
+
+ /*! @brief 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.
+ */
+ EntityPointer outside() const
+ {
+ return this->real.outside();
+ }
+
+ /*! @brief 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.
+ */
+ const LocalGeometry& intersectionSelfLocal () const
+ {
+ return this->real.intersectionSelfLocal();
+ }
+ /*! @brief 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.
+ */
+ const LocalGeometry& intersectionNeighborLocal () const
+ {
+ return this->real.intersectionNeighborLocal();
+ }
+
+ /*! @brief 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.
+ */
+ const Geometry& intersectionGlobal () const
+ {
+ return this->real.intersectionGlobal();
+ }
+
+ //! Local number of codim 1 entity in the inside() Entity where intersection is contained in
+ int numberInSelf () const
+ {
+ return this->real.numberInSelf ();
+ }
+
+ //! Local number of codim 1 entity in outside() Entity where intersection is contained in
+ int numberInNeighbor () const
+ {
+ return this->real.numberInNeighbor ();
+ }
+
+ /*! @brief Return an outer normal (length not necessarily 1)
+
+ The returned vector may depend on local position within the intersection.
+ */
+ FieldVector<ctype, dimworld> outerNormal (const FieldVector<ctype, dim-1>& local) const
+ {
+ return this->real.outerNormal(local);
+ }
+
+ /*! @brief return outer normal scaled with the integration element
+ @copydoc outerNormal
+ 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> integrationOuterNormal (const FieldVector<ctype, dim-1>& local) const
+ {
+ return this->real.integrationOuterNormal(local);
+ }
+
+ /*! @brief 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.
+ */
+ FieldVector<ctype, dimworld> unitOuterNormal (const FieldVector<ctype, dim-1>& local) const
+ {
+ return this->real.unitOuterNormal(local);
+ }
+
+ /** @brief Checks for equality.
+ Only s 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 Intersection& rhs) const
+ {
+ return rhs.equals(*this);
+ }
+
+ /** @brief Checks for inequality.
+ Only s 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 Intersection& rhs) const
+ {
+ return ! rhs.equals(*this);
+ }
+
+
+ //===========================================================
+ /** @name Implementor interface
+ */
+ //@{
+ //===========================================================
+
+ /** @brief forward equality check to real */
+ bool equals(const Intersection& rhs) const
+ {
+ return this->real.equals(rhs.real);
+ }
+
+ /** Copy Constructor from IntersectionImp */
+ Intersection(const IntersectionImp<const GridImp> & i) :
+ real(i) {};
+
+ /** Copy constructor */
+ Intersection(const Intersection& i) :
+ real(i.real) {}
+ //@}
+
+ typedef typename remove_const<GridImp>::type mutableGridImp;
+protected:
+ // give the GridDefaultImplementation class access to the realImp
+ friend class GridDefaultImplementation<
+ GridImp::dimension, GridImp::dimensionworld,
+ typename GridImp::ctype,
+ typename GridImp::GridFamily> ;
+
+ //! return reference to the real implementation
+ ImplementationType & getRealImp() { return real; }
+ //! return reference to the real implementation
+ const ImplementationType & getRealImp() const { return real; }
+
+};
+
+//**********************************************************************
+/**
+ @brief Default Implementations of integrationOuterNormal and unitOuterNormal for IntersectionImp
+
+ @ingroup GridDevel
+*/
+template<class GridImp, template<class> class IntersectionImp>
+class IntersectionDefaultNormalVectors
+{
+ enum { dim=GridImp::dimension };
+ enum { dimworld=GridImp::dimensionworld };
+ typedef typename GridImp::ctype ct;
+public:
+ //! return unit outer normal, this should be dependent on
+ //! local coordinates for higher order boundary
+ //! the normal is scaled with the integration element of the intersection.
+ FieldVector<ct, dimworld> integrationOuterNormal (const FieldVector<ct, dim-1>& local) const
+ {
+ FieldVector<ct, dimworld> n = unitOuterNormal(local);
+ n *= asImp().intersectionGlobal().integrationElement(local);
+ return n;
+ }
+ //! return unit outer normal
+ FieldVector<ct, dimworld> unitOuterNormal (const FieldVector<ct, dim-1>& local) const
+ {
+ FieldVector<ct, dimworld> n = asImp().outerNormal(local);
+ n /= n.two_norm();
+ return n;
+ }
+
+private:
+ //! Barton-Nackman trick
+ IntersectionImp<GridImp>& asImp ()
+ {return static_cast<IntersectionImp<GridImp>&>(*this);}
+ const IntersectionImp<GridImp>& asImp () const
+ {return static_cast<const IntersectionImp<GridImp>&>(*this);}
+};
+
+}
+
+#endif // DUNE_GRID_INTERSECTION_HH
Index: grid/common/entity.hh
===================================================================
--- grid/common/entity.hh (Revision 4067)
+++ grid/common/entity.hh (Arbeitskopie)
@@ -605,7 +605,7 @@
IntersectionIterator end = asImp().ilevelend();
for(IntersectionIterator it = asImp().ilevelbegin(); it != end; ++it)
{
- if( it.boundary() ) return true;
+ if( it->boundary() ) return true;
}
}
@@ -615,7 +615,7 @@
IntersectionIterator end = asImp().ileafend();
for(IntersectionIterator it = asImp().ileafbegin(); it != end; ++it)
{
- if( it.boundary() ) return true;
+ if( it->boundary() ) return true;
}
}
Index: grid/common/grid.hh
===================================================================
--- grid/common/grid.hh (Revision 4067)
+++ grid/common/grid.hh (Arbeitskopie)
@@ -215,14 +215,17 @@
of a given entity of codimension 0. %EntityPointer is copy-constructible from a
%HierarchicIterator.
- - %IntersectionIterator which is a model of Dune::IntersectionIterator
- provides access to all entities of codimension 0 that have an intersection of
- codimension 1 with a given entity of codimension 0. In a conforming mesh these
- are the face neighbors of an element. For two entities with a common intersection
- the %IntersectionIterator also provides information about the geometric location
+ - %Intersection which is a model of Dune::Intersection
+ provides access an intersection of codimension 1 of two entity of codimension 0
+ or one entity and the boundary. In a conforming mesh this
+ is a face of an element. For two entities with a common intersection
+ the %Intersection also provides information about the geometric location
of the intersection. Furthermore it also provides information about intersections
of an entity with the internal or external boundaries.
+ - %IntersectionIterator which is a model of Dune::IntersectionIterator
+ provides access to all intersections of a given entity of codimension 0.
+
- %LevelIndexSet and %LeafIndexSet which are both models of Dune::IndexSet are
used to attach any kind of user-defined data to (subsets of) entities of the grid.
This data is supposed to be stored in one-dimensional arrays for reasons
@@ -300,6 +303,13 @@
<TD>no</TD>
</TR>
<TR>
+ <TD>Intersection</TD>
+ <TD>yes</TD>
+ <TD>no</TD>
+ <TD>yes</TD>
+ <TD>no</TD>
+ </TR>
+ <TR>
<TD>IntersectionIterator</TD>
<TD>yes</TD>
<TD>no</TD>
@@ -411,7 +421,8 @@
template<class GridImp, class EntityPointerImp> class EntityPointer;
template<int codim, PartitionIteratorType pitype, class GridImp,
template<int,PartitionIteratorType,class> class LevelIteratorImp> class LevelIterator;
-template<class GridImp, template<class> class IntersectionIteratorImp> class IntersectionIterator;
+template<class GridImp, template<class> class IntersectionImp> class Intersection;
+template<class GridImp, template<class> class IntersectionIteratorImp, template<class> class IntersectionImp> class IntersectionIterator;
template<class GridImp, template<class> class HierarchicIteratorImp> class HierarchicIterator;
template<int codim, PartitionIteratorType pitype, class GridImp,
template<int,PartitionIteratorType,class> class LeafIteratorImp> class LeafIterator;
@@ -516,6 +527,16 @@
typedef typename GridFamily::Traits::template Codim<cd>::template Partition<pitype>::LeafIterator LeafIterator;
};
+ /*! \brief A type that is a model of Dune::LeafIntersection, an
+ intersections of two codimension 1 of two codimension 0 entities in the leaf view.
+ */
+ typedef typename GridFamily::Traits::LeafIntersection LeafIntersection;
+
+ /*! \brief A type that is a model of Dune::Intersection, an
+ intersections of two codimension 1 of two codimension 0 entities in a level view.
+ */
+ typedef typename GridFamily::Traits::LevelIntersection LevelIntersection;
+
/*! \brief A type that is a model of Dune::IntersectionIterator
which is an iterator that allows to examine, but not to modify, the
intersections of codimension 1 of an leaf element (entity of codimension 0)
@@ -1129,6 +1150,8 @@
template<int,int,class> class EntityImp,
template<int,class> class EntityPointerImp,
template<int,PartitionIteratorType,class> class LevelIteratorImp,
+ template<class> class LeafIntersectionImp,
+ template<class> class LevelIntersectionImp,
template<class> class LeafIntersectionIteratorImp,
template<class> class LevelIntersectionIteratorImp,
template<class> class HierarchicIteratorImp,
@@ -1139,8 +1162,10 @@
{
typedef GridImp Grid;
- typedef Dune::IntersectionIterator<const GridImp, LeafIntersectionIteratorImp> LeafIntersectionIterator;
- typedef Dune::IntersectionIterator<const GridImp, LevelIntersectionIteratorImp> LevelIntersectionIterator;
+ typedef Dune::Intersection<const GridImp, LeafIntersectionImp> LeafIntersection;
+ typedef Dune::Intersection<const GridImp, LevelIntersectionImp> LevelIntersection;
+ typedef Dune::IntersectionIterator<const GridImp, LeafIntersectionIteratorImp, LeafIntersectionImp> LeafIntersectionIterator;
+ typedef Dune::IntersectionIterator<const GridImp, LevelIntersectionIteratorImp, LevelIntersectionImp> LevelIntersectionIterator;
typedef Dune::HierarchicIterator<const GridImp, HierarchicIteratorImp> HierarchicIterator;
@@ -1208,6 +1233,7 @@
#include "entity.hh"
#include "entitypointer.hh"
#include "leveliterator.hh"
+#include "intersection.hh"
#include "intersectioniterator.hh"
#include "hierarchiciterator.hh"
#include "leafiterator.hh"
Index: grid/common/gridpart.hh
===================================================================
--- grid/common/gridpart.hh (Revision 4067)
+++ grid/common/gridpart.hh (Arbeitskopie)
@@ -235,6 +235,9 @@
//! Level index set that corresponds to the grid
typedef typename Traits::IndexSetType IndexSetType;
+ //! The corresponding Intersection
+ typedef typename Traits::IntersectionType IntersectionType ;
+
//! The corresponding IntersectionIterator
typedef typename Traits::IntersectionIteratorType IntersectionIteratorType ;
@@ -319,6 +322,10 @@
/** \brief The appropriate index set */
typedef WrappedLevelIndexSet<GridType> IndexSetType;
+ /** \brief The appropriate intersection */
+ typedef typename GridType::Traits::
+ LevelIntersection IntersectionType;
+
/** \brief The appropriate intersection iterator */
typedef typename GridType::template Codim<0>::Entity::
LevelIntersectionIterator IntersectionIteratorType;
@@ -347,6 +354,9 @@
//! The leaf index set of the grid implementation
typedef typename Traits::IndexSetType IndexSetType;
+ //! The corresponding Intersection
+ typedef typename Traits::IntersectionType IntersectionType ;
+
//! The corresponding IntersectionIterator
typedef typename Traits::IntersectionIteratorType IntersectionIteratorType ;
@@ -422,6 +432,10 @@
/** \brief The appropriate index set */
typedef WrappedLeafIndexSet<GridType> IndexSetType;
+ /** \brief The appropriate intersection */
+ typedef typename GridType::Traits::
+ LeafIntersection IntersectionType;
+
/** \brief The appropriate intersection iterator */
typedef typename GridType::template Codim<0>::Entity::
LeafIntersectionIterator IntersectionIteratorType;
Index: grid/common/intersectioniterator.hh
===================================================================
--- grid/common/intersectioniterator.hh (Revision 4067)
+++ grid/common/intersectioniterator.hh (Arbeitskopie)
@@ -15,11 +15,15 @@
- <tt>GridImp</tt> Type that is a model of Dune::Grid
- <tt>IntersectionIteratorImp</tt> Class template that is a model of
Dune::IntersectionIterator
-
- @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.
+ @warning the IntersectionIterator used to be both, Intersection and IntersectionIterator,
+ at the same time. The two concepts are now properly seperated. The IntersectionIterator
+ still offers the old methods, but these are forwarded to the Intersection. All these methods
+ are now marked deprecated.
+
+ \deprecated All Intersection methods are deprecated,
+ dereference IntersectionIterator to get the Intersection and call methods there.
+
@warning The number of neigbors may be different from the number of
faces/edges of an element!
@@ -76,105 +80,9 @@
level and the leaf intersection iterators will return the element \b c
together with the right grid boundary.
- <h2>Methods neighbor and boundary </h2>
-
- 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:
-
- <table>
- <tr>
- <td></td><td>intersection</td><td>neighbor()</td><td>boundary()</td><td>outside()</td></tr>
- <tr>
- <td>1</td><td>with inner, overlap <br>
- or ghost entity</td>
- <td>true</td><td>false</td>
- <td>the neighbor entity</td></tr>
- <tr>
- <td>2</td><td>on domain boundary</td>
- <td>false</td><td>true</td><td><em>undefined</em></td></tr>
- <tr>
- <td>3</td><td>on periodic boundary</td>
- <td>true</td><td>true</td><td>Ghost-/Overlap cell at br (with transformed geometry)</td></tr>
- <tr>
- <td>4</td><td>on processor boundary</td>
- <td>false <em>if grid has no ghosts</em><br>true <em>otherwise</em></td><td>false </td>
- <td>ghost entity <em>(if it exists)</em></td></tr>
- </table>
-
- -# <b> Inner Intersections: </b> \n
- The type of the neighboring entity can be determined through
- methods defined on the outside entity.
- -# <b> Handling physical boundaries: </b>
- 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. \br
- 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...
- .
- -# <b> Handling periodic boundaries: </b>
- - 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
- -# <b> Processor boundaries: </b> \n
- 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.
- .
-
-
- <h2>Geometry of an intersection</h2>
-
- 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.
-
-
@ingroup GIIntersectionIterator
*/
-template<class GridImp, template<class> class IntersectionIteratorImp>
+template<class GridImp, template<class> class IntersectionIteratorImp, template<class> class IntersectionImp>
class IntersectionIterator
{
IntersectionIteratorImp<const GridImp> realIterator;
@@ -193,6 +101,9 @@
/** \brief Pointer to the type of entities that this IntersectionIterator belongs to */
typedef typename GridImp::template Codim<0>::EntityPointer EntityPointer;
+ /** \brief Type of Intersection this IntersectionIterator points to */
+ typedef Dune::Intersection< const GridImp, IntersectionImp > Intersection;
+
/** \brief Codim 1 geometry returned by intersectionGlobal() */
typedef typename GridImp::template Codim<1>::Geometry Geometry;
@@ -209,6 +120,53 @@
//! define type used for coordinates in grid module
typedef typename GridImp::ctype ctype;
+ //===========================================================
+ /** @name Dereferencing
+ */
+ //@{
+ //===========================================================
+
+ /** \brief Dereferencing operator. */
+ const Intersection & operator*() const
+ {
+ return this->realIterator.dereference();
+ }
+
+ /** \brief Pointer operator. */
+ const Intersection * operator->() const
+ {
+ return & this->realIterator.dereference();
+ }
+ //@}
+
+
+ //===========================================================
+ /** @name Compare methods
+ */
+ //@{
+ //===========================================================
+
+ /** @brief 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
+ {
+ return rhs.equals(*this);
+ }
+
+ /** @brief 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.
+ */
+ bool operator!=(const IntersectionIterator& rhs) const
+ {
+ return ! rhs.equals(*this);
+ }
+ //@}
+
/** @brief Preincrement operator. Proceed to next intersection.*/
IntersectionIterator& operator++()
{
@@ -216,8 +174,14 @@
return *this;
}
+ //===========================================================
+ /** @name Query methods
+ */
+ //@{
+ //===========================================================
+
//! return true if intersection is with interior or exterior boundary (see the cases above)
- bool boundary () const
+ bool boundary () const DUNE_DEPRECATED
{
return this->realIterator.boundary();
}
@@ -236,13 +200,13 @@
The way the Identifiers are attached to the grid may differ
between the different grid implementations.
*/
- int boundaryId () const
+ int boundaryId () const DUNE_DEPRECATED
{
return this->realIterator.boundaryId();
- }
+ }
//! @brief return true if intersection is shared with another element.
- bool neighbor () const
+ bool neighbor () const DUNE_DEPRECATED
{
return this->realIterator.neighbor();
}
@@ -250,7 +214,7 @@
/*! @brief return EntityPointer to the Entity on the inside of this
intersection. That is the Entity where we started this Iterator.
*/
- EntityPointer inside() const
+ EntityPointer inside() const DUNE_DEPRECATED
{
return this->realIterator.inside();
}
@@ -261,7 +225,7 @@
@warning Don't call this method if there is no neighboring Entity
(neighbor() returns false). In this case the result is undefined.
*/
- EntityPointer outside() const
+ EntityPointer outside() const DUNE_DEPRECATED
{
return this->realIterator.outside();
}
@@ -273,7 +237,7 @@
local coordinates of the intersection to local coordinates of the
inside() entity.
*/
- const LocalGeometry& intersectionSelfLocal () const
+ const LocalGeometry& intersectionSelfLocal () const DUNE_DEPRECATED
{
return this->realIterator.intersectionSelfLocal();
}
@@ -284,7 +248,7 @@
local coordinates of the intersection to local coordinates of the
outside() entity.
*/
- const LocalGeometry& intersectionNeighborLocal () const
+ const LocalGeometry& intersectionNeighborLocal () const DUNE_DEPRECATED
{
return this->realIterator.intersectionNeighborLocal();
}
@@ -294,19 +258,19 @@
This method returns a Geometry object that provides a mapping from
local coordinates of the intersection to global (world) coordinates.
*/
- const Geometry& intersectionGlobal () const
+ const Geometry& intersectionGlobal () const DUNE_DEPRECATED
{
return this->realIterator.intersectionGlobal();
}
//! Local number of codim 1 entity in the inside() Entity where intersection is contained in
- int numberInSelf () const
+ int numberInSelf () const DUNE_DEPRECATED
{
return this->realIterator.numberInSelf ();
}
//! Local number of codim 1 entity in outside() Entity where intersection is contained in
- int numberInNeighbor () const
+ int numberInNeighbor () const DUNE_DEPRECATED
{
return this->realIterator.numberInNeighbor ();
}
@@ -315,7 +279,7 @@
The returned vector may depend on local position within the intersection.
*/
- FieldVector<ctype, dimworld> outerNormal (const FieldVector<ctype, dim-1>& local) const
+ FieldVector<ctype, dimworld> outerNormal (const FieldVector<ctype, dim-1>& local) const DUNE_DEPRECATED
{
return this->realIterator.outerNormal(local);
}
@@ -326,7 +290,7 @@
method is redundant but it may be more efficent to use this function
rather than computing the integration element via intersectionGlobal().
*/
- FieldVector<ctype, dimworld> integrationOuterNormal (const FieldVector<ctype, dim-1>& local) const
+ FieldVector<ctype, dimworld> integrationOuterNormal (const FieldVector<ctype, dim-1>& local) const DUNE_DEPRECATED
{
return this->realIterator.integrationOuterNormal(local);
}
@@ -336,32 +300,12 @@
The returned vector may depend on the local position within the intersection.
It is scaled to have unit length.
*/
- FieldVector<ctype, dimworld> unitOuterNormal (const FieldVector<ctype, dim-1>& local) const
+ FieldVector<ctype, dimworld> unitOuterNormal (const FieldVector<ctype, dim-1>& local) const DUNE_DEPRECATED
{
return this->realIterator.unitOuterNormal(local);
}
+ //@}
- /** @brief 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
- {
- return rhs.equals(*this);
- }
-
- /** @brief 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.
- */
- bool operator!=(const IntersectionIterator& rhs) const
- {
- return ! rhs.equals(*this);
- }
-
-
//===========================================================
/** @name Implementor interface
*/
@@ -398,19 +342,21 @@
};
-//**********************************************************************
-/**
- @brief Default Implementations for IntersectionIteratorImp
+} // namespace Dune
- @ingroup GridDevel
-*/
-template<class GridImp, template<class> class IntersectionIteratorImp>
+#include "intersection.hh"
+
+namespace Dune {
+
+template<class GridImp, template<class> class IntersectionAndIteratorImp>
class IntersectionIteratorDefaultImplementation
{
enum { dim=GridImp::dimension };
enum { dimworld=GridImp::dimensionworld };
typedef typename GridImp::ctype ct;
public:
+ typedef Dune::Intersection<const GridImp, IntersectionAndIteratorImp> Intersection;
+ typedef IntersectionAndIteratorImp<const GridImp> ImplementationType;
//! return unit outer normal, this should be dependent on
//! local coordinates for higher order boundary
//! the normal is scaled with the integration element of the intersection.
@@ -427,15 +373,19 @@
n /= n.two_norm();
return n;
}
-
+ //! \brief dereferencing
+ const Intersection & dereference() const
+ {
+ return reinterpret_cast<const Intersection&>(*this);
+ }
private:
//! Barton-Nackman trick
- IntersectionIteratorImp<GridImp>& asImp ()
- {return static_cast<IntersectionIteratorImp<GridImp>&>(*this);}
- const IntersectionIteratorImp<GridImp>& asImp () const
- {return static_cast<const IntersectionIteratorImp<GridImp>&>(*this);}
-};
+ ImplementationType& asImp ()
+ {return static_cast<ImplementationType&>(*this);}
+ const ImplementationType& asImp () const
+ {return static_cast<const ImplementationType&>(*this);}
+} DUNE_DEPRECATED;
-}
+} // namespace Dune
#endif // DUNE_GRID_INTERSECTIONITERATOR_HH
Index: grid/onedgrid.hh
===================================================================
--- grid/onedgrid.hh (Revision 4067)
+++ grid/onedgrid.hh (Arbeitskopie)
@@ -156,19 +156,21 @@
OneDGridEntity,
OneDGridEntityPointer,
OneDGridLevelIterator,
+ OneDGridLeafIntersectionIterator, // leaf intersection
+ OneDGridLevelIntersectionIterator, // level intersection
OneDGridLeafIntersectionIterator, // leaf intersection iter
OneDGridLevelIntersectionIterator, // level intersection iter
OneDGridHierarchicIterator,
OneDGridLeafIterator,
OneDGridLevelIndexSet<const OneDGrid>,
- OneDGridLevelIndexSetTypes<const OneDGrid>,
+ OneDGridLevelIndexSetTypes<const OneDGrid>,
OneDGridLeafIndexSet<const OneDGrid>,
- OneDGridLeafIndexSetTypes<const OneDGrid>,
+ OneDGridLeafIndexSetTypes<const OneDGrid>,
OneDGridIdSet<const OneDGrid>,
unsigned int,
OneDGridIdSet<const OneDGrid>,
unsigned int,
- CollectiveCommunication<Dune::OneDGrid> >
+ CollectiveCommunication<Dune::OneDGrid> >
Traits;
};
Index: grid/sgrid.hh
===================================================================
--- grid/sgrid.hh (Revision 4067)
+++ grid/sgrid.hh (Arbeitskopie)
@@ -803,6 +803,13 @@
return normal;
}
+ //! return integration outer normal
+ FieldVector<ctype, GridImp::dimensionworld> integrationOuterNormal (const FieldVector<ctype, GridImp::dimension-1>& local) const
+ {
+ FieldVector<sgrid_ctype, dimworld> n = unitOuterNormal(local);
+ n *= is_global.integrationElement(local);
+ return n;
+ }
/*! intersection of codimension 1 of this neighbor with element where iteration started.
Here returned element is in LOCAL coordinates of the element where iteration started.
@@ -836,7 +843,6 @@
if (partition != it.partition)
DUNE_THROW(GridError, "assignment of SIntersectionIterator "
<< "with different partition");
- // (zred == it.zred) if (self == it.self)
/* Assign current position and update ne */
ne = it.ne;
@@ -1163,6 +1169,8 @@
typedef GridTraits<dim,dimworld,Dune::SGrid<dim,dimworld>,
SGeometry,SEntity,
SEntityPointer,SLevelIterator,
+ SIntersectionIterator, // leaf intersection
+ SIntersectionIterator, // level intersection
SIntersectionIterator, // leaf intersection iter
SIntersectionIterator, // level intersection iter
SHierarchicIterator,
Index: grid/sgrid/sgrid.cc
===================================================================
--- grid/sgrid/sgrid.cc (Revision 4067)
+++ grid/sgrid/sgrid.cc (Arbeitskopie)
@@ -765,7 +765,7 @@
inline void SGrid<dim,dimworld>::makeSGrid (const int* N_,
const sgrid_ctype* L_, const sgrid_ctype* H_)
{
- IsTrue< dimworld <= std::numeric_limits<int>::digits >::yes();
+ dune_static_assert(dimworld <= std::numeric_limits<int>::digits,"world dimension too high, must be <= # of bits of int");
N = new array<int,dim>[MAXL];
h = new FieldVector<sgrid_ctype, dim>[MAXL];
@@ -793,7 +793,7 @@
template<int dim, int dimworld>
inline SGrid<dim,dimworld>::SGrid (const int* N_, const sgrid_ctype* H_)
{
- IsTrue< dimworld <= std::numeric_limits<int>::digits >::yes();
+ dune_static_assert(dimworld <= std::numeric_limits<int>::digits,"world dimension too high, must be <= # of bits of int");
sgrid_ctype L_[dim];
for (int i=0; i<dim; i++)
Index: grid/modules
===================================================================
--- grid/modules (Revision 4067)
+++ grid/modules (Arbeitskopie)
@@ -23,7 +23,7 @@
@ingroup GridInterface */
/** @defgroup GIEntityPointer EntityPointer and Iterators
@ingroup GridInterface */
- /** @defgroup GIIntersectionIterator IntersectionIterator
+ /** @defgroup GIIntersectionIterator Intersection and IntersectionIterator
@ingroup GridInterface */
/** @defgroup IndexIdSets IndexSet and IdSet
@ingroup GridInterface */
More information about the Dune
mailing list