diff --git a/grid/yaspgrid.hh b/grid/yaspgrid.hh
index 7c1226dd7e0f307a5f5fb333e9506115915fe3b0..4d817ade586591c1308468374d5d33fefb43fa7a 100644
--- a/grid/yaspgrid.hh
+++ b/grid/yaspgrid.hh
@@ -50,6 +50,7 @@ namespace Dune {
template<int mydim, int cdim, class GridImp> class YaspGeometry;
template<int codim, int dim, class GridImp> class YaspEntity;
template<int codim, PartitionIteratorType pitype, class GridImp> class YaspLevelIterator;
+ template<class GridImp> class YaspLeafIterator;
template<class GridImp> class YaspIntersectionIterator;
template<class GridImp> class YaspHierarchicIterator;
template<class GridImp> class YaspBoundaryEntity;
@@ -286,7 +287,7 @@ namespace Dune {
s << " " << extension[i];
s << " missing is " << missing;
}
-
+ private:
const YaspGeometry<mydim,cdim,GridImp>&
operator = (const YaspGeometry<mydim,cdim,GridImp>& g)
{
@@ -306,8 +307,8 @@ namespace Dune {
// IMPORTANT midpoint and extension can't be references,
// because they must stay the same when the iterator changes
- FieldVector<ctype, cdim> midpoint; // the midpoint
- FieldVector<ctype, cdim> extension; // the extension
+ const FieldVector<ctype, cdim> & midpoint; // the midpoint
+ const FieldVector<ctype, cdim> & extension; // the extension
int& missing; // the missing, i.e. constant direction
// In addition we need memory in order to return references.
@@ -429,7 +430,7 @@ namespace Dune {
for (int i=0; i<mydim; i++)
s << " " << extension[i];
}
-
+ private:
const YaspGeometry<mydim,mydim,GridImp>&
operator = (const YaspGeometry<mydim,mydim,GridImp>& g)
{
@@ -448,8 +449,8 @@ namespace Dune {
// IMPORTANT midpoint and extension can't be references,
// because they must stay the same when the iterator changes
- FieldVector<ctype, mydim> midpoint; // the midpoint
- FieldVector<ctype, mydim> extension; // the extension
+ const FieldVector<ctype, mydim> & midpoint; // the midpoint
+ const FieldVector<ctype, mydim> & extension; // the extension
// In addition we need memory in order to return references.
// Possibly we should change this in the interface ...
@@ -495,7 +496,7 @@ namespace Dune {
s << "YaspGeometry<"<<0<<","<<cdim<< "> ";
s << "position " << position;
}
-
+ private:
const YaspGeometry<0,cdim,GridImp>&
operator = (const YaspGeometry<0,cdim,GridImp>& g)
{
@@ -506,7 +507,7 @@ namespace Dune {
private:
// IMPORTANT position can't be references,
// because they must stay the same when the iterator changes
- FieldVector<ctype, cdim> position; //!< where the vertex is
+ const FieldVector<ctype, cdim> & position; //!< where the vertex is
};
// operator<< for all YaspGeometrys
@@ -546,11 +547,11 @@ namespace Dune {
{};
const TSI& transformingsubiterator () const
{
- return this->realEntity._it;
+ return this->realEntity.transformingsubiterator();
}
const YGLI& gridlevel () const
{
- return this->realEntity._g;
+ return this->realEntity.gridlevel();
}
};
@@ -629,7 +630,7 @@ namespace Dune {
int level () const {return _g.level();}
//! index is unique and consecutive per level
- int index () const {return _it.superindex();} // superindex works also for iteration over subgrids
+ int index () const { return _it.superindex();} // superindex works also for iteration over subgrids
//! return partition type attribute
PartitionType partitionType () const
@@ -673,6 +674,10 @@ namespace Dune {
return YaspLevelIterator<cc,All_Partition,GridImp>(_g,_g.vertex_overlapfront().tsubbegin(coord));
}
+ if (cc==0)
+ {
+ return *this;
+ }
DUNE_THROW(GridError, "codim not (yet) implemented");
}
@@ -728,6 +733,14 @@ namespace Dune {
return _g;
}
+ bool isLeaf() const
+ {
+ std::cout << _g.level() << " ... " << _g.mg()->maxlevel() << std::endl;
+ if (_g.level() == _g.mg()->maxlevel())
+ std::cout << "is Leaf\n";
+ return (_g.level() == _g.mg()->maxlevel());
+ }
+
//! returns intersection iterator for first intersection
IntersectionIterator ibegin () const
{
@@ -906,6 +919,7 @@ namespace Dune {
enum { dim=GridImp::dimension };
enum { dimworld=GridImp::dimensionworld };
typedef typename GridImp::ctype ctype;
+ YaspIntersectionIterator();
public:
// types used from grids
typedef typename MultiYGrid<dim,ctype>::YGridLevelIterator YGLI;
@@ -1009,19 +1023,20 @@ namespace Dune {
return _nb;
}
-
//! return unit outer normal, this should be dependent on local coordinates for higher order boundary
- FieldVector<ctype, dimworld>& unitOuterNormal (FieldVector<ctype, dim-1>& local) const
+ FieldVector<ctype, dimworld> outerNormal (const FieldVector<ctype, dim-1>& local) const
{
return _normal;
}
- //! return unit outer normal, if you know it is constant use this function instead
- FieldVector<ctype, dimworld>& unitOuterNormal () const
+
+ //! return unit outer normal, this should be dependent on local coordinates for higher order boundary
+ FieldVector<ctype, dimworld> unitOuterNormal (const FieldVector<ctype, dim-1>& local) const
{
return _normal;
}
+
/*! 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.
*/
@@ -1069,7 +1084,7 @@ namespace Dune {
_ext_local(1.0),
_is_self_local(_pos_self_local,_ext_local,_dir),
_is_nb_local(_pos_nb_local,_ext_local,_dir),
- _is_global(_pos_world,myself.transformingsubiterator().meshsize(),_dir),
+ _is_global(_pos_world,_myself.transformingsubiterator().meshsize(),_dir),
_normal(0.0)
{
// making an end iterator?
@@ -1079,7 +1094,6 @@ namespace Dune {
_count = 2*dim;
return;
}
-
// initialize to first neighbor
_count = 0;
_dir = 0;
@@ -1097,6 +1111,26 @@ namespace Dune {
_normal[0] = -1.0;
}
+ //! copy constructor
+ YaspIntersectionIterator (const YaspIntersectionIterator<GridImp>& it)
+ // important! _itnb and _nb must recreated not copied
+ : _count(it._count),
+ _dir(it._dir),
+ _face(it._face),
+ _itnb(it._itnb),
+ _myself(it._myself),
+ _nb(_myself.gridlevel(),_itnb),
+ _pos_self_local(it._pos_self_local),
+ _pos_nb_local(it._pos_nb_local),
+ _pos_world(it._pos_world),
+ _ext_local(it._ext_local),
+ // important! _is_* must recreated not copied
+ _is_self_local(_pos_self_local,_ext_local,_dir),
+ _is_nb_local(_pos_nb_local,_ext_local,_dir),
+ _is_global(_pos_world,_myself.transformingsubiterator().meshsize(),_dir),
+ _normal(it._normal)
+ {}
+
private:
int _count; //!< valid neighbor count in 0 .. 2*dim-1
int _dir; //!< count/2
@@ -1107,7 +1141,7 @@ namespace Dune {
mutable SpecialEntity _nb; //!< virtual neighbor entity, built on the fly
FieldVector<ctype, dim> _pos_self_local; //!< center of face in own local coordinates
FieldVector<ctype, dim> _pos_nb_local; //!< center of face in neighbors local coordinates
- FieldVector<ctype, dim> _pos_world; //!< center of face in world coordinates
+ FieldVector<ctype, dimworld>_pos_world; //!< center of face in world coordinates
FieldVector<ctype, dim> _ext_local; //!< extension of face in local coordinates
SpecialLocalGeometry _is_self_local; //!< intersection in own local coordinates
SpecialLocalGeometry _is_nb_local; //!< intersection in neighbors local coordinates
@@ -1160,6 +1194,12 @@ namespace Dune {
pop_tos();
}
+ //! constructor
+ YaspHierarchicIterator (const YaspHierarchicIterator& it) :
+ _g(it._g), _it(it._it), _entity(_g,_it),
+ _maxlevel(it._maxlevel), stack(it.stack)
+ {}
+
//! increment
void increment ()
{
@@ -1332,6 +1372,20 @@ namespace Dune {
mutable SpecialEntity _entity; //!< virtual entity
};
+ template <class GridImp>
+ class YaspLeafIterator : public YaspLevelIterator<0,All_Partition,GridImp>
+ {
+ typedef typename GridImp::ctype ctype;
+ public:
+ typedef typename MultiYGrid<0,ctype>::YGridLevelIterator YGLI;
+ typedef typename SubYGrid<0,ctype>::TransformingSubIterator TSI;
+ YaspLeafIterator(const YGLI & g, const TSI & it) :
+ YaspLevelIterator<0,All_Partition,GridImp>(g,it)
+ {}
+ YaspLeafIterator(const YaspLevelIterator<0,All_Partition,GridImp>& i) :
+ YaspLevelIterator<0,All_Partition,GridImp>(i)
+ {}
+ };
//************************************************************************
diff --git a/grid/yaspgrid/grids.hh b/grid/yaspgrid/grids.hh
index 05ecc808c342a046daf10e7977c078aa8168ec94..966bc8596ea67cb152d758e3d5c2d8a756f83135 100644
--- a/grid/yaspgrid/grids.hh
+++ b/grid/yaspgrid/grids.hh
@@ -1719,6 +1719,11 @@ namespace Dune {
i=start; l=level;
}
+ //! make iterator pointing to level k (no check made)
+ YGridLevelIterator (const YGridLevelIterator & it)
+ : i(it.i), l(it.l)
+ {}
+
//! return number of this grid level
int level () const
{