#include <base/forward-declarations.h>
#include <base/point.h>
#include <grid/geometry_info.h>
+//#include <grid/tria_boundary.h>
#include <base/subscriptor.h>
-
-
/**
* Declare some symbolic names for mesh smoothing algorithms. The meaning of
* these flags is documented in the #Triangulation# class.
*
* \subsection{Boundary approximation}
*
+ * To be updated!
+ *
* You can specify a boundary function: if a new vertex is created on a
* side or face at the boundary, this function is used to compute where
* it will be placed. See \Ref{Boundary} for the details. Usage with
* @see TriaRawIterator
* @author Wolfgang Bangerth, 1998 */
template <int dim>
-class Triangulation : public TriaDimensionInfo<dim>, public Subscriptor {
+class Triangulation
+ : public TriaDimensionInfo<dim>,
+ public Subscriptor
+{
+ private:
+ /**
+ * Default boundary object. This declaration is used
+ * for the default argument in #set_boundary#.
+ */
+ static const StraightBoundary<dim>& straight_boundary;
+
public:
typedef typename TriaDimensionInfo<dim>::raw_line_iterator raw_line_iterator;
typedef typename TriaDimensionInfo<dim>::line_iterator line_iterator;
/**
* Assign a boundary object to the
- * triangulation, which is used to find
- * the point where to place new vertices
- * on the boundary. Ownership of the
- * object remains with the caller of this
- * function, i.e. you have to make sure
- * that it is not destroyed before
- * #Triangulation<>::execute_refinement_and_refinement()#
- * is called the last time. Checking this
- * is mostly done by the library by
- * using a #Subscriptor# object as base
- * class for boundary objects.
- *
- * If you copy this triangulation to
- * another one using the
- * #copy_triangulation# function, you must
- * make sure that the lifetime of boundary
- * object extends to the lifetime of the
- * new triangulation as well.
- *
- * Because of the use of the #Subscriptor#
- * base class, you must wait for the
- * #Triangulation# object to release the
- * lock to the boundary object before
- * destroying that. This usually happens
- * when the destructor of the
- * triangulation is executed, but you may
- * manually do so by calling
- * #set_boundary(0)#, passing a Null
- * pointer as new boundary object. This
- * is the interpreted as "take the
- * default boundary object", which
- * actually is a #StraightBoundary#
- * object.
- *
- * You are allowed to remove or replace
+ * triangulation. If a face with boundary number #number#
+ * is refined, this object is used to find
+ * the location of new vertices
+ * on the boundary. This will also be true for non-linear
+ * transformations of cells to the unit cell in shape
+ * function calculations. Multiple calls to this function are
+ * allowed to store different boundary curves or surfaces.
+ *
+ * Numbers of boundary object correspond to material numbers of
+ * faces at the boundary, for instance the material id in a UCD
+ * input file. They are not necessarily consecutive but must be
+ * in the range 0-254.
+ *
+ * The #boundary_object# is not copied and MUST persist
+ * until the triangulation is destroyed. Otherwise,
+ * the #Subscriptor# class will issue #ExcObjectInUse#.
+ * This is also true for triangulations generated from this
+ * one by #copy_triangulation#.
+ *
+ * It is possible to remove or replace
* the boundary object during the lifetime
- * of a non-empty triangulation, but
- * you should be sure about what you do.
- */
- void set_boundary (const Boundary<dim> *boundary_object);
-
- /**
- * Return a constant reference to the boundary
- * object used for this triangulation. This
- * is necessary for a host of applications
- * that need to use the exact shape of the
- * boundary. Since the boundary used by all
- * these applications, it is useful to ask
- * the triangulation rather than asking
- * for another parameter for each of these
- * functions.
- *
- * If you need to store a pointer to the
- * boundary somewhere, for example
- * in the constructor of an object, for
- * later use in the member functions,
- * you should consider calling #subscribe#
- * on the boundary to guarantee that the
- * boundary object is still alive during
- * the lifetime of your object.
- */
- const Boundary<dim> & get_boundary () const;
+ * of a non-empty triangulation. Usually, this is done
+ * before the first refinement and is dangerous afterwards.
+ * Removal of a boundary object is done by
+ * #set_boundary(number)#.
+ */
+ void set_boundary (unsigned int number,
+ const Boundary<dim>& boundary_object = straight_boundary);
+
+ /**
+ * Return a constant reference to a boundary
+ * object used for this triangulation.
+ * Number is the same as in #set_boundary#
+ */
+ const Boundary<dim> & get_boundary (unsigned int number) const;
/**
* Copy a triangulation. This operation is
vector<bool> vertices_used;
/**
- * Pointer to a boundary object.
+ * Collection of boundary objects.
*/
- const Boundary<dim> *boundary;
+ const Boundary<dim>* boundary[255];
/**
* Do some smoothing in the process
+template <int dim>
+const StraightBoundary<dim>& Triangulation<dim>::straight_boundary
+= StraightBoundary<dim>();
+
+
template <int dim>
Triangulation<dim>::Triangulation (const MeshSmoothing smooth_grid) :
Subscriptor (),
- boundary (0),
smooth_grid(smooth_grid)
{
- // set default boundary
- set_boundary (0);
+ // set default boundary for all possible components
+ for (unsigned int i=0;i<255;++i)
+ {
+ boundary[i] = &straight_boundary;
+ straight_boundary.subscribe();
+ }
levels.push_back (new TriangulationLevel<dim>);
};
template <int dim>
-Triangulation<dim>::Triangulation (const Triangulation<dim> &) :
+Triangulation<dim>::Triangulation (const Triangulation<dim> &)
+ :
Subscriptor () // do not set any subscriptors; anyway,
// calling this constructor is an error!
{
template <int dim>
-Triangulation<dim>::~Triangulation () {
+Triangulation<dim>::~Triangulation ()
+{
for (unsigned int i=0; i<levels.size(); ++i)
delete levels[i];
- boundary->unsubscribe ();
+ for (unsigned int i=0;i<255;++i)
+ boundary[i]->unsubscribe ();
levels.clear ();
};
-
template <int dim>
-void Triangulation<dim>::set_boundary (const Boundary<dim> *boundary_object) {
- static StraightBoundary<dim> default_boundary;
-
- if (boundary != 0)
- boundary->unsubscribe ();
-
- if (boundary_object != 0)
- boundary = boundary_object;
- else
- boundary = &default_boundary;
-
- boundary->subscribe ();
+void
+Triangulation<dim>::set_boundary (unsigned int number,
+ const Boundary<dim>& boundary_object)
+{
+ Assert(number<255, ExcIndexRange(number,0,255));
+
+ boundary[number]->unsubscribe ();
+ boundary[number] = &boundary_object;
+ boundary_object.subscribe();
};
template <int dim>
const Boundary<dim> &
-Triangulation<dim>::get_boundary () const
+Triangulation<dim>::get_boundary (unsigned int number) const
{
- return *boundary;
+ Assert(number<255, ExcIndexRange(number,0,255));
+
+ return *(boundary[number]);
};
template <int dim>
-void Triangulation<dim>::copy_triangulation (const Triangulation<dim> &old_tria) {
+void Triangulation<dim>::copy_triangulation (const Triangulation<dim> &old_tria)
+{
Assert (vertices.size() == 0, ExcTriangulationNotEmpty());
Assert (n_cells () == 0, ExcTriangulationNotEmpty());
Assert (levels.size () == 1, ExcTriangulationNotEmpty());
vertices_used = old_tria.vertices_used;
smooth_grid = old_tria.smooth_grid;
- if (boundary)
- boundary->unsubscribe ();
- boundary = old_tria.boundary;
- boundary->subscribe ();
-
+ for (unsigned i=0;i<255;++i)
+ {
+ boundary[i]->unsubscribe ();
+ boundary[i] = old_tria.boundary[i];
+ boundary[i]->subscribe ();
+ }
// delete only level previously existing,
// reserve space for new and copy
delete levels[0];
// where shall we put the new
// vertex?
Point<2> new_point;
- if (cell->at_boundary(nb))
+
+ face_iterator face=cell->line(nb);
+
+ if ( face->boundary_indicator() != 255 )
+ {
// boundary vertex
- new_point = boundary->get_new_point_on_line (cell->line(nb));
- else {
+ new_point = boundary[face->boundary_indicator()]->
+ get_new_point_on_line (face);
+ } else {
// vertex between two
// normal cells
new_point = vertices[new_vertices[2*nb]];
if (line->at_boundary())
vertices[next_unused_vertex]
- = boundary->get_new_point_on_line (line);
+ = boundary[line->boundary_indicator()]->get_new_point_on_line (line);
else
vertices[next_unused_vertex]
= (line->vertex(0) + line->vertex(1)) / 2;
if (quad->at_boundary())
vertices[next_unused_vertex]
- = boundary->get_new_point_on_quad (quad);
+ = boundary[quad->boundary_indicator()]->get_new_point_on_quad (quad);
else
vertices[next_unused_vertex]
= (quad->vertex(0) + quad->vertex(1) +