// $Id$
// Version: $Name$
//
-// Copyright (C) 1998, 1999, 2000, 2001, 2002, 2003, 2004, 2005, 2006, 2007, 2008, 2009 by the deal.II authors
+// Copyright (C) 1998, 1999, 2000, 2001, 2002, 2003, 2004, 2005, 2006, 2007, 2008, 2009, 2010 by the deal.II authors
//
// This file is subject to QPL and may not be distributed
// without copyright and license information. Please refer
*/
namespace TriaAccessorExceptions
{
+//TODO: Write documentation!
/**
* @ingroup Exceptions
*/
DeclException0 (ExcCellNotUsed);
/**
+ * The cell is not an @ref
+ * GlossActive "active" cell, but
+ * it already has children. Some
+ * operations, like setting
+ * refinement flags or accessing
+ * degrees of freedom are only
+ * possible on active cells.
+ *
* @ingroup Exceptions
*/
DeclException0 (ExcCellNotActive);
/**
+ * Trying to access the children of
+ * a cell which is in fact active.
+ *
* @ingroup Exceptions
*/
DeclException0 (ExcCellHasNoChildren);
+//TODO: Write documentation!
/**
* @ingroup Exceptions
*/
DeclException0 (ExcUnusedCellAsChild);
+//TODO: Write documentation!
/**
* @ingroup Exceptions
*/
<< "You can only set the child index if the cell has no "
<< "children, or clear it. The given "
<< "index was " << arg1 << " (-1 means: clear children)");
+//TODO: Write documentation!
/**
* @ingroup Exceptions
*/
DeclException0 (ExcUnusedCellAsNeighbor);
+//TODO: Write documentation!
/**
* @ingroup Exceptions
*/
DeclException0 (ExcUncaughtCase);
+//TODO: Write documentation!
/**
* @ingroup Exceptions
*/
DeclException0 (ExcDereferenceInvalidObject);
+//TODO: Write documentation!
/**
* @ingroup Exceptions
*/
DeclException0 (ExcCantCompareIterators);
+//TODO: Write documentation!
/**
* @ingroup Exceptions
*/
DeclException0 (ExcNeighborIsCoarser);
+//TODO: Write documentation!
/**
* @ingroup Exceptions
*/
DeclException0 (ExcNeighborIsNotCoarser);
/**
+ * You are trying to access the
+ * level of a face, but faces have
+ * no inherent level. The level of
+ * a face can only be determined by
+ * the level of an adjacent face,
+ * which in turn implies that a
+ * face can have several levels.
+ *
* @ingroup Exceptions
*/
DeclException0 (ExcFacesHaveNoLevel);
+//TODO: Write documentation!
/**
* @ingroup Exceptions
*/
* operators and the like), but has no functionality to actually
* dereference data. This is done in the derived classes.
*
- * This template is used for faces and edges, which have no level in
- * the triangulation hierarchy associated with them. There exists a
- * partial specialization of the current class template where
- * <tt>structdim</tt> (the dimensionality of the object represented,
- * for example 1 for a line) equals <tt>dim</tt> (the dimensionality
- * of the space the object lives in, for example 3 if we solve
- * three-dimensional problems).
+ * In the implementation, the behavior of this class differs between
+ * the cases where <tt>structdim==dim</tt> (cells of a mesh) and
+ * <tt>structdim<dim</tt> (faces and edges). For the latter,
+ * #present_level is always equal to zero and the constructors may not
+ * receive a positive value there. For cells, any level is
+ * possible, but only those within the range of the levels of the
+ * Triangulation are reasonable. Furthermore, the function objects()
+ * returns either the container with all cells on the same level or
+ * the container with all objects of this dimension (<tt>structdim<dim</tt>).
*
* Some internals of this class are discussed in @ref IteratorAccessorInternals .
*
* @ingroup grid
* @ingroup Accessors
- * @author Wolfgang Bangerth, 1998
+ * @author Wolfgang Bangerth, Guido Kanschat, 1998, 2010
*/
template <int structdim, int dim, int spacedim=dim>
class TriaAccessorBase
* pointed to belongs to.
* This is only valid for cells.
*/
- static int level ();
-
- /**
- * Return the index of the
- * element presently pointed to
- * on the present level.
- */
- int index () const;
-
- /**
- * Return the state of the
- * iterator. For the different
- * states an accessor can be in,
- * refer to the
- * TriaRawIterator
- * documentation.
- */
- IteratorState::IteratorStates state () const;
-
- /**
- * Return a pointer to the
- * triangulation which the object
- * pointed to by this class
- * belongs to.
- */
- const Triangulation<dim,spacedim> & get_triangulation () const;
-
- /**
- * @}
- */
- protected:
-
- /**
- * Used to store the index of
- * the element presently pointed
- * to on the level presentl
- * used.
- */
- int present_index;
-
- /**
- * Pointer to the triangulation
- * which we act on.
- */
- const Triangulation<dim,spacedim> *tria;
-
- private:
-
- template <typename Accessor> friend class TriaRawIterator;
- template <typename Accessor> friend class TriaIterator;
- template <typename Accessor> friend class TriaActiveIterator;
-};
-
-
-
-/**
- * A base class for the accessor classes used by TriaRawIterator and
- * derived classes.
- *
- * This class offers only the basic functionality required by the
- * iterators (stores the necessary data members, offers comparison
- * operators and the like), but has no functionality to actually
- * dereference data. This is done in the derived classes.
- *
- * This template is a partial specialization of the general template
- * used for cells only, i.e. for the case where <tt>structdim</tt>
- * (the dimensionality of the object represented, for example 1 for a
- * line) equals <tt>dim</tt> (the dimensionality of the triangulation
- * the object lives in, for example 2 if we solve problems on a
- * two-dimensional manifold -- whether this manifold is itself the
- * space $R^2$ or embedded in a higher-dimensional space whose
- * dimensionality is determined by the <code>spacedim</code> template
- * argument). The reason for using a partial specialization rather
- * than the general template is that cells have a level in the
- * triangulation hierarchy associated with them, whereas faces and
- * edges do not.
- *
- * Some internals of this class are discussed in @ref IteratorAccessorInternals .
- *
- * @ingroup grid
- * @ingroup Accessors
- * @author Wolfgang Bangerth, 1998
- */
-template <int dim, int spacedim>
-class TriaAccessorBase<dim,dim,spacedim>
-{
- public:
- /**
- * Dimension of the space the
- * object represented by this
- * accessor lives in. For
- * example, if this accessor
- * represents a quad that is
- * part of a two-dimensional
- * surface in four-dimensional
- * space, then this value is
- * four.
- */
- static const unsigned int space_dimension = spacedim;
-
- /**
- * Dimensionality of the object
- * that the thing represented by
- * this accessopr is part of. For
- * example, if this accessor
- * represents a line that is part
- * of a hexahedron, then this
- * value will be three.
- */
- static const unsigned int dimension = dim;
-
- /**
- * Dimensionality of the current
- * object represented by this
- * accessor. For example, if it
- * is line (irrespective of
- * whether it is part of a quad
- * or hex, and what dimension we
- * are in), then this value
- * equals 1.
- */
- static const unsigned int structure_dimension = dim;
-
- protected:
- /**
- * Declare the data type that
- * this accessor class expects to
- * get passed from the iterator
- * classes. Since the pure
- * triangulation iterators need
- * no additional data, this data
- * type is @p void.
- */
- typedef void AccessorData;
-
- /**
- * Constructor. Protected, thus
- * only callable from friend
- * classes.
- */
- TriaAccessorBase (const Triangulation<dim,spacedim> *parent = 0,
- const int level = -1,
- const int index = -1,
- const AccessorData * = 0);
-
- /**
- * Copy operator. Since this is
- * only called from iterators,
- * do not return anything, since
- * the iterator will return
- * itself.
- *
- * This method is protected,
- * since it is only to be called
- * from the iterator class.
- */
- void copy_from (const TriaAccessorBase &);
-
- /**
- * Copy operator. This is normally
- * used in a context like
- * <tt>iterator a,b; *a=*b;</tt>. Since
- * the meaning is to copy the object
- * pointed to by @p b to the object
- * pointed to by @p a and since
- * accessors are not real but
- * virtual objects, this operation
- * is not useful for iterators on
- * triangulations. We declare this
- * function here private, thus it may
- * not be used from outside.
- * Furthermore it is not implemented
- * and will give a linker error if
- * used anyway.
- */
- void operator = (const TriaAccessorBase *);
-
- /**
- * Same as above.
- */
- TriaAccessorBase & operator = (const TriaAccessorBase &);
-
- /**
- * @name Advancement of iterators
- */
- /**
- * @{
- */
-
- /**
- * This operator advances the
- * iterator to the next element.
- *
- * For cells only:
- * The next element is next on
- * this level if there are
- * more. If the present element
- * is the last on this level,
- * the first on the next level
- * is accessed.
- */
- void operator ++ ();
-
- /**
- * This operator moves the
- * iterator to the previous
- * element.
- *
- * For cells only:
- * The previous element is
- * previous on this level if
- * <tt>index>0</tt>. If the present
- * element is the first on this
- * level, the last on the
- * previous level is accessed.
- */
- void operator -- ();
-
- /**
- * Access to the other objects of
- * a Triangulation with same
- * dimension.
- */
- internal::Triangulation::TriaObjects<internal::Triangulation::TriaObject<dim> > &
- objects () const;
-
- protected:
-
- /**
- * Compare for equality.
- */
- bool operator == (const TriaAccessorBase &) const;
-
- /**
- * Compare for inequality.
- */
- bool operator != (const TriaAccessorBase &) const;
-
- public:
- /**
- * Data type to be used for passing
- * parameters from iterators to the
- * accessor classes in a unified
- * way, no matter what the type of
- * number of these parameters is.
- */
- typedef void * LocalData;
-
- /** @name Iterator address and state
- */
- /**
- * @{
- */
- /**
- * Return the level the element
- * pointed to belongs to.
- * This is only valid for cells.
- */
int level () const;
/**
*/
protected:
/**
- * Used to store the level
- * presently pointed to. For
- * faces this should always be
- * 0.
+ * The level if this is a cell
+ * (<tt>structdim==dim</tt>). Else,
+ * contains zero.
*/
int present_level;
* which we act on.
*/
const Triangulation<dim,spacedim> *tria;
-
+
private:
template <typename Accessor> friend class TriaRawIterator;
template <typename Accessor> friend class TriaIterator;
template <typename Accessor> friend class TriaActiveIterator;
-
- template <int, int, int> friend class TriaAccessor;
};
+
/**
* A class that represents accessor objects to iterators that don't
* make sense such as quad iterators in on 1d meshes. This class can
template <int structdim, int dim, int spacedim>
inline
-TriaAccessorBase<structdim,dim,spacedim>::TriaAccessorBase (const Triangulation<dim,spacedim> *parent,
- const int level,
- const int index,
- const AccessorData *)
+TriaAccessorBase<structdim,dim,spacedim>::TriaAccessorBase (
+ const Triangulation<dim,spacedim>* tria,
+ const int level,
+ const int index,
+ const AccessorData*)
:
+ present_level((structdim==dim) ? level : 0),
present_index (index),
- tria (parent)
+ tria (tria)
{
-
+
// non-cells have no level, so a 0
// should have been passed, or a -1
// for an end-iterator, or -2 for
// an invalid (default constructed)
// iterator
- Assert ((level == 0) || (level == -1) || (level == -2),
- ExcInternalError());
+ if (structdim != dim)
+ {
+ Assert ((level == 0) || (level == -1) || (level == -2),
+ ExcInternalError());
+ }
}
void
TriaAccessorBase<structdim,dim,spacedim>::copy_from (const TriaAccessorBase<structdim,dim,spacedim> &a)
{
+ present_level = a.present_level;
present_index = a.present_index;
tria = a.tria;
+
+ if (structdim != dim)
+ {
+ Assert ((present_level == 0) || (present_level == -1) || (present_level == -2),
+ ExcInternalError());
+ }
+}
+
+
+
+template <int structdim, int dim, int spacedim>
+inline
+TriaAccessorBase<structdim,dim,spacedim>&
+TriaAccessorBase<structdim,dim,spacedim>::operator= (const TriaAccessorBase<structdim,dim,spacedim> &a)
+{
+ present_level = a.present_level;
+ present_index = a.present_index;
+ tria = a.tria;
+
+ if (structdim != dim)
+ {
+ Assert ((present_level == 0) || (present_level == -1) || (present_level == -2),
+ ExcInternalError());
+ }
+ return *this;
}
TriaAccessorBase<structdim,dim,spacedim>::operator == (const TriaAccessorBase<structdim,dim,spacedim> &a) const
{
Assert (tria == a.tria, TriaAccessorExceptions::ExcCantCompareIterators());
- return (present_index == a.present_index);
+ return ((present_level == a.present_level) &&
+ (present_index == a.present_index));
}
TriaAccessorBase<structdim,dim,spacedim>::operator != (const TriaAccessorBase<structdim,dim,spacedim> &a) const
{
Assert (tria == a.tria, TriaAccessorExceptions::ExcCantCompareIterators());
- return (present_index != a.present_index);
+ return ((present_level != a.present_level) ||
+ (present_index != a.present_index));
}
template <int structdim, int dim, int spacedim>
inline
int
-TriaAccessorBase<structdim,dim,spacedim>::level ()
+TriaAccessorBase<structdim,dim,spacedim>::level () const
{
- return 0;
+ // This is always zero or invalid
+ // if the object is not a cell
+ return present_level;
}
IteratorState::IteratorStates
TriaAccessorBase<structdim,dim,spacedim>::state () const
{
- if (present_index>=0)
+ if ((present_level>=0) && (present_index>=0))
return IteratorState::valid;
else
if (present_index==-1)
// this iterator is used for
// objects without level
++this->present_index;
- // is index still in the range of
- // the vector? (note that we don't
- // have to set the level, since
- // dim!=1 and the object therefore
- // has no level)
- if (this->present_index
- >=
- static_cast<int>(objects().cells.size()))
- this->present_index = -1;
+
+ if (structdim != dim)
+ {
+ // is index still in the range of
+ // the vector? (note that we don't
+ // have to set the level, since
+ // dim!=1 and the object therefore
+ // has no level)
+ if (this->present_index
+ >=
+ static_cast<int>(objects().cells.size()))
+ this->present_index = -1;
+ }
+ else
+ {
+ while (this->present_index
+ >=
+ static_cast<int>(this->tria->levels[this->present_level]->cells.cells.size()))
+ {
+ // no -> go one level up until we find
+ // one with more than zero cells
+ ++this->present_level;
+ this->present_index = 0;
+ // highest level reached?
+ if (this->present_level >= static_cast<int>(this->tria->levels.size()))
+ {
+ // return with past the end pointer
+ this->present_level = this->present_index = -1;
+ return;
+ }
+ }
+ }
}
// same as operator++
--this->present_index;
- if (this->present_index < 0)
- this->present_index = -1;
+ if (structdim != dim)
+ {
+ if (this->present_index < 0)
+ this->present_index = -1;
+ }
+ else
+ {
+ while (this->present_index < 0)
+ {
+ // no -> go one level down
+ --this->present_level;
+ // lowest level reached?
+ if (this->present_level == -1)
+ {
+ // return with past the end pointer
+ this->present_level = this->present_index = -1;
+ return;
+ }
+ // else
+ this->present_index = this->tria->levels[this->present_level]->cells.cells.size()-1;
+ }
+ }
}
*/
template <int dim>
inline
- internal::Triangulation::TriaObjects<internal::Triangulation::TriaObject<1> > &
+ internal::Triangulation::TriaObjects<internal::Triangulation::TriaObject<1> >*
get_objects (internal::Triangulation::TriaFaces<dim> *faces,
const internal::int2type<1>)
{
- return faces->lines;
+ return &faces->lines;
}
template <int dim>
inline
- internal::Triangulation::TriaObjects<internal::Triangulation::TriaObject<2> > &
+ internal::Triangulation::TriaObjects<internal::Triangulation::TriaObject<2> >*
get_objects (internal::Triangulation::TriaFaces<dim> *faces,
const internal::int2type<2>)
{
- return faces->quads;
+ return &faces->quads;
}
- }
-}
-
-
-
-
-template <int structdim, int dim, int spacedim>
-inline
-internal::Triangulation::TriaObjects<internal::Triangulation::TriaObject<structdim> > &
-TriaAccessorBase<structdim,dim,spacedim>::objects() const
-{
- // get sub-objects. note that the
- // current class is only used for
- // objects that are *not* cells
- return internal::TriaAccessorBase::get_objects (this->tria->faces,
- internal::int2type<structdim> ());
-}
-
-
-
-
-/*------------------------ Functions: TriaAccessorBase<dim,dim> ---------------------------*/
-
-template <int dim, int spacedim>
-inline
-TriaAccessorBase<dim,dim,spacedim>::TriaAccessorBase (const Triangulation<dim,spacedim> *parent,
- const int level,
- const int index,
- const AccessorData *)
- :
- present_level (level),
- present_index (index),
- tria (parent)
-{}
-
-
-
-template <int dim, int spacedim>
-inline
-void
-TriaAccessorBase<dim,dim,spacedim>::copy_from (const TriaAccessorBase<dim,dim,spacedim> &a)
-{
- present_level = a.present_level;
- present_index = a.present_index;
- tria = a.tria;
-}
-
-
-
-template <int dim, int spacedim>
-inline
-bool
-TriaAccessorBase<dim,dim,spacedim>::operator == (const TriaAccessorBase<dim,dim,spacedim> &a) const
-{
- Assert (tria == a.tria, TriaAccessorExceptions::ExcCantCompareIterators());
- return ((present_index == a.present_index) &&
- (present_level == a.present_level));
-}
-
-
-
-template <int dim, int spacedim>
-inline
-bool
-TriaAccessorBase<dim,dim,spacedim>::operator != (const TriaAccessorBase<dim,dim,spacedim> &a) const
-{
- Assert (tria == a.tria, TriaAccessorExceptions::ExcCantCompareIterators());
- return ((present_index != a.present_index) ||
- (present_level != a.present_level));
-}
-
-
-
-template <int dim, int spacedim>
-inline
-int
-TriaAccessorBase<dim,dim,spacedim>::level () const
-{
- return present_level;
-}
-
-
-
-template <int dim, int spacedim>
-inline
-int
-TriaAccessorBase<dim,dim,spacedim>::index () const
-{
- return present_index;
-}
-
-
-
-template <int dim, int spacedim>
-inline
-IteratorState::IteratorStates
-TriaAccessorBase<dim,dim,spacedim>::state () const
-{
- if ((present_level>=0) && (present_index>=0))
- return IteratorState::valid;
- else
- if ((present_index==-1) && (present_index==-1))
- return IteratorState::past_the_end;
- else
- return IteratorState::invalid;
-}
-
-
-
-template <int dim, int spacedim>
-inline
-const Triangulation<dim,spacedim> &
-TriaAccessorBase<dim,dim,spacedim>::get_triangulation () const
-{
- return *tria;
-}
-
-
-
-template <int dim, int spacedim>
-inline
-void
-TriaAccessorBase<dim,dim,spacedim>::operator ++ ()
-{
- // in structdim == dim, so we have
- // cells here. these have levels
- ++this->present_index;
- // is index still in the range of
- // the vector?
- while (this->present_index
- >=
- static_cast<int>(this->tria->levels[this->present_level]->cells.cells.size()))
+
+ inline
+ internal::Triangulation::TriaObjects<internal::Triangulation::TriaObject<1> >*
+ get_objects (internal::Triangulation::TriaFaces<1>*,
+ const internal::int2type<1>)
{
- // no -> go one level up until we find
- // one with more than zero cells
- ++this->present_level;
- this->present_index = 0;
- // highest level reached?
- if (this->present_level >= static_cast<int>(this->tria->levels.size()))
- {
- // return with past the end pointer
- this->present_level = this->present_index = -1;
- return;
- }
+ Assert (false, ExcInternalError());
+ return 0;
+ }
+
+ inline
+ internal::Triangulation::TriaObjects<internal::Triangulation::TriaObject<2> >*
+ get_objects (internal::Triangulation::TriaFaces<2>*,
+ const internal::int2type<2>)
+ {
+ Assert (false, ExcInternalError());
+ return 0;
+ }
+
+ inline
+ internal::Triangulation::TriaObjects<internal::Triangulation::TriaObject<3> >*
+ get_objects (internal::Triangulation::TriaFaces<3>*,
+ const internal::int2type<3>)
+ {
+ Assert (false, ExcInternalError());
+ return 0;
+ }
+
+ /**
+ * This function should never be
+ * used, but we need it for the
+ * template instantiation of TriaAccessorBase<dim,dim,spacedim>::objects() const
+ */
+ template <int dim>
+ inline
+ internal::Triangulation::TriaObjects<internal::Triangulation::TriaObject<3> >*
+ get_objects (internal::Triangulation::TriaFaces<dim> *faces,
+ const internal::int2type<3>)
+ {
+ Assert (false, ExcInternalError());
+ return 0;
}
-}
+ /**
+ * Copy the above functions for
+ * cell objects.
+ */
+ template <int structdim, int dim>
+ inline
+ internal::Triangulation::TriaObjects<internal::Triangulation::TriaObject<structdim> >*
+ get_objects (internal::Triangulation::TriaObjects<internal::Triangulation::TriaObject<dim> >*,
+ const internal::int2type<structdim>)
+ {
+ Assert (false, ExcInternalError());
+ return 0;
+ }
-template <int dim, int spacedim>
-inline
-void
-TriaAccessorBase<dim,dim,spacedim>::operator -- ()
-{
- // same as operator++
- --this->present_index;
- // is index still in the range of
- // the vector?
- while (this->present_index < 0)
+ template <int dim>
+ inline
+ internal::Triangulation::TriaObjects<internal::Triangulation::TriaObject<dim> >*
+ get_objects (internal::Triangulation::TriaObjects<internal::Triangulation::TriaObject<dim> >* cells,
+ const internal::int2type<dim>)
{
- // no -> go one level down
- --this->present_level;
- // lowest level reached?
- if (this->present_level == -1)
- {
- // return with past the end pointer
- this->present_level = this->present_index = -1;
- return;
- }
- // else
- this->present_index = this->tria->levels[this->present_level]->cells.cells.size()-1;
+ return cells;
}
+ }
}
-template <int dim, int spacedim>
+template <int structdim, int dim, int spacedim>
inline
-internal::Triangulation::TriaObjects<internal::Triangulation::TriaObject<dim> > &
-TriaAccessorBase<dim,dim,spacedim>::objects() const
+internal::Triangulation::TriaObjects<internal::Triangulation::TriaObject<structdim> > &
+TriaAccessorBase<structdim,dim,spacedim>::objects() const
{
- return this->tria->levels[this->present_level]->cells;
+ if (structdim != dim)
+ // get sub-objects. note that the
+ // current class is only used for
+ // objects that are *not* cells
+ return *internal::TriaAccessorBase::get_objects (this->tria->faces,
+ internal::int2type<structdim> ());
+ else
+ return *internal::TriaAccessorBase::get_objects (&this->tria->levels[this->present_level]->cells,
+ internal::int2type<structdim> ());
}
-
-
/*------------------------ Functions: InvalidAccessor ---------------------------*/
template <int structdim, int dim, int spacedim>
// $Id$
// Version: $Name$
//
-// Copyright (C) 2003, 2004, 2005, 2006, 2007, 2008 by the deal.II authors
+// Copyright (C) 2003, 2004, 2005, 2006, 2007, 2008, 2010 by the deal.II authors
//
// This file is subject to QPL and may not be distributed
// without copyright and license information. Please refer
template <int dim, int spacedim> class CellAccessor;
template <int, int, int> class TriaAccessorBase;
-template <int dim, int spacedim> class TriaAccessorBase<dim,dim,spacedim>;
template <int, int, int> class InvalidAccessor;
template <int, int, int> class TriaAccessor;
template <int dim, int spacedim> class TriaAccessor<0, dim, spacedim>;