Triangulation< dim, spacedim > Class Template Reference
[Grid classes]

Inheritance diagram for Triangulation< dim, spacedim >:
Inheritance graph
[legend]

List of all members.

Classes

class  ExcEmptyLevel
class  ExcFacesHaveNoLevel
class  ExcGridReadError
class  ExcInvalidLevel
class  ExcTriangulationNotEmpty
class  RefinementListener

Public Types

enum  MeshSmoothing {
  none = 0x0, limit_level_difference_at_vertices = 0x1, eliminate_unrefined_islands = 0x2, patch_level_1 = 0x4,
  coarsest_level_1 = 0x8, allow_anisotropic_smoothing = 0x10, eliminate_refined_inner_islands = 0x100, eliminate_refined_boundary_islands = 0x200,
  do_not_produce_unrefined_islands = 0x400, smoothing_on_refinement, smoothing_on_coarsening, maximum_smoothing = 0xffff ^ allow_anisotropic_smoothing
}
typedef
IteratorSelector::raw_line_iterator 
raw_line_iterator
typedef
IteratorSelector::line_iterator 
line_iterator
typedef
IteratorSelector::active_line_iterator 
active_line_iterator
typedef
IteratorSelector::raw_quad_iterator 
raw_quad_iterator
typedef
IteratorSelector::quad_iterator 
quad_iterator
typedef
IteratorSelector::active_quad_iterator 
active_quad_iterator
typedef
IteratorSelector::raw_hex_iterator 
raw_hex_iterator
typedef
IteratorSelector::hex_iterator 
hex_iterator
typedef
IteratorSelector::active_hex_iterator 
active_hex_iterator
typedef
IteratorSelector::raw_cell_iterator 
raw_cell_iterator
typedef
IteratorSelector::cell_iterator 
cell_iterator
typedef
IteratorSelector::active_cell_iterator 
active_cell_iterator
typedef
IteratorSelector::raw_face_iterator 
raw_face_iterator
typedef
IteratorSelector::face_iterator 
face_iterator
typedef
IteratorSelector::active_face_iterator 
active_face_iterator

Public Member Functions

 Triangulation (const MeshSmoothing smooth_grid=none)
 Triangulation (const Triangulation< dim, spacedim > &t)
virtual ~Triangulation ()
void clear ()
void set_boundary (const unsigned int number, const Boundary< dim, spacedim > &boundary_object)
void set_boundary (const unsigned int number)
const Boundary< dim, spacedim > & get_boundary (const unsigned int number) const
std::vector< unsigned char > get_boundary_indicators () const
virtual void copy_triangulation (const Triangulation< dim, spacedim > &old_tria)
virtual void create_triangulation (const std::vector< Point< spacedim > > &vertices, const std::vector< CellData< dim > > &cells, const SubCellData &subcelldata)
virtual void create_triangulation_compatibility (const std::vector< Point< spacedim > > &vertices, const std::vector< CellData< dim > > &cells, const SubCellData &subcelldata)
void distort_random (const double factor, const bool keep_boundary=true)
virtual unsigned int memory_consumption () const
Mesh refinement



void set_all_refine_flags ()
void refine_global (const unsigned int times)
virtual void execute_coarsening_and_refinement ()
bool prepare_coarsening_and_refinement ()
void add_refinement_listener (RefinementListener &listener) const
void remove_refinement_listener (RefinementListener &listener) const
History of a triangulation



void save_refine_flags (std::ostream &out) const
void save_refine_flags (std::vector< bool > &v) const
void load_refine_flags (std::istream &in)
void load_refine_flags (const std::vector< bool > &v)
void save_coarsen_flags (std::ostream &out) const
void save_coarsen_flags (std::vector< bool > &v) const
void load_coarsen_flags (std::istream &out)
void load_coarsen_flags (const std::vector< bool > &v)
bool get_anisotropic_refinement_flag () const
User data



void clear_user_flags ()
void save_user_flags (std::ostream &out) const
void save_user_flags (std::vector< bool > &v) const
void load_user_flags (std::istream &in)
void load_user_flags (const std::vector< bool > &v)
void clear_user_flags_line ()
void save_user_flags_line (std::ostream &out) const
void save_user_flags_line (std::vector< bool > &v) const
void load_user_flags_line (std::istream &in)
void load_user_flags_line (const std::vector< bool > &v)
void clear_user_flags_quad ()
void save_user_flags_quad (std::ostream &out) const
void save_user_flags_quad (std::vector< bool > &v) const
void load_user_flags_quad (std::istream &in)
void load_user_flags_quad (const std::vector< bool > &v)
void clear_user_flags_hex ()
void save_user_flags_hex (std::ostream &out) const
void save_user_flags_hex (std::vector< bool > &v) const
void load_user_flags_hex (std::istream &in)
void load_user_flags_hex (const std::vector< bool > &v)
void clear_user_data ()
void clear_user_pointers ()
void save_user_indices (std::vector< unsigned int > &v) const
void load_user_indices (const std::vector< unsigned int > &v)
void save_user_pointers (std::vector< void * > &v) const
void load_user_pointers (const std::vector< void * > &v)
void save_user_indices_line (std::vector< unsigned int > &v) const
void load_user_indices_line (const std::vector< unsigned int > &v)
void save_user_indices_quad (std::vector< unsigned int > &v) const
void load_user_indices_quad (const std::vector< unsigned int > &v)
void save_user_indices_hex (std::vector< unsigned int > &v) const
void load_user_indices_hex (const std::vector< unsigned int > &v)
void save_user_pointers_line (std::vector< void * > &v) const
void load_user_pointers_line (const std::vector< void * > &v)
void save_user_pointers_quad (std::vector< void * > &v) const
void load_user_pointers_quad (const std::vector< void * > &v)
void save_user_pointers_hex (std::vector< void * > &v) const
void load_user_pointers_hex (const std::vector< void * > &v)
Cell iterator functions



raw_cell_iterator begin_raw (const unsigned int level=0) const
cell_iterator begin (const unsigned int level=0) const
active_cell_iterator begin_active (const unsigned int level=0) const
raw_cell_iterator end () const
cell_iterator end (const unsigned int level) const
raw_cell_iterator end_raw (const unsigned int level) const
active_cell_iterator end_active (const unsigned int level) const
raw_cell_iterator last_raw () const
raw_cell_iterator last_raw (const unsigned int level) const
cell_iterator last () const
cell_iterator last (const unsigned int level) const
active_cell_iterator last_active () const
active_cell_iterator last_active (const unsigned int level) const
Face iterator functions



raw_face_iterator begin_raw_face () const
face_iterator begin_face () const
active_face_iterator begin_active_face () const
raw_face_iterator end_face () const
raw_face_iterator end_raw_face () const
active_face_iterator end_active_face () const
raw_face_iterator last_raw_face () const
face_iterator last_face () const
active_face_iterator last_active_face () const
Line iterator functions



raw_line_iterator begin_raw_line (const unsigned int level=0) const
line_iterator begin_line (const unsigned int level=0) const
active_line_iterator begin_active_line (const unsigned int level=0) const
raw_line_iterator end_line () const
line_iterator end_line (const unsigned int level) const
raw_line_iterator end_raw_line (const unsigned int level) const
active_line_iterator end_active_line (const unsigned int level) const
raw_line_iterator last_raw_line () const
raw_line_iterator last_raw_line (const unsigned int level) const
line_iterator last_line () const
line_iterator last_line (const unsigned int level) const
active_line_iterator last_active_line () const
active_line_iterator last_active_line (const unsigned int level) const
Quad iterator functions



raw_quad_iterator begin_raw_quad (const unsigned int level=0) const
quad_iterator begin_quad (const unsigned int level=0) const
active_quad_iterator begin_active_quad (const unsigned int level=0) const
raw_quad_iterator end_quad () const
quad_iterator end_quad (const unsigned int level) const
raw_quad_iterator end_raw_quad (const unsigned int level) const
active_quad_iterator end_active_quad (const unsigned int level) const
raw_quad_iterator last_raw_quad () const
raw_quad_iterator last_raw_quad (const unsigned int level) const
quad_iterator last_quad () const
quad_iterator last_quad (const unsigned int level) const
active_quad_iterator last_active_quad () const
active_quad_iterator last_active_quad (const unsigned int level) const
Hex iterator functions



raw_hex_iterator begin_raw_hex (const unsigned int level=0) const
hex_iterator begin_hex (const unsigned int level=0) const
active_hex_iterator begin_active_hex (const unsigned int level=0) const
raw_hex_iterator end_hex () const
hex_iterator end_hex (const unsigned int level) const
raw_hex_iterator end_raw_hex (const unsigned int level) const
active_hex_iterator end_active_hex (const unsigned int level) const
raw_hex_iterator last_raw_hex () const
raw_hex_iterator last_raw_hex (const unsigned int level) const
hex_iterator last_hex () const
hex_iterator last_hex (const unsigned int level) const
active_hex_iterator last_active_hex () const
active_hex_iterator last_active_hex (const unsigned int level) const
Information about the triangulation



unsigned int n_raw_lines () const
unsigned int n_raw_lines (const unsigned int level) const
unsigned int n_lines () const
unsigned int n_lines (const unsigned int level) const
unsigned int n_active_lines () const
unsigned int n_active_lines (const unsigned int level) const
unsigned int n_raw_quads () const
unsigned int n_raw_quads (const unsigned int level) const
unsigned int n_quads () const
unsigned int n_quads (const unsigned int level) const
unsigned int n_active_quads () const
unsigned int n_active_quads (const unsigned int level) const
unsigned int n_raw_hexs () const
unsigned int n_raw_hexs (const unsigned int level) const
unsigned int n_hexs () const
unsigned int n_hexs (const unsigned int level) const
unsigned int n_active_hexs () const
unsigned int n_active_hexs (const unsigned int level) const
unsigned int n_raw_cells (const unsigned int level) const
unsigned int n_cells () const
unsigned int n_cells (const unsigned int level) const
unsigned int n_active_cells () const
unsigned int n_active_cells (const unsigned int level) const
unsigned int n_faces () const
unsigned int n_active_faces () const
unsigned int n_levels () const
unsigned int n_vertices () const
const std::vector< Point
< spacedim > > & 
get_vertices () const
unsigned int n_used_vertices () const
bool vertex_used (const unsigned int index) const
const std::vector< bool > & get_used_vertices () const
unsigned int max_adjacent_cells () const

Static Public Attributes

static const StraightBoundary
< dim, spacedim > 
straight_boundary
static const StraightBoundary
< spacedim, spacedim > 
straight_manifold_description
static const unsigned int dimension = dim
static const unsigned int space_dimension = spacedim

Static Protected Member Functions

static void write_bool_vector (const unsigned int magic_number1, const std::vector< bool > &v, const unsigned int magic_number2, std::ostream &out)
static void read_bool_vector (const unsigned int magic_number1, std::vector< bool > &v, const unsigned int magic_number2, std::istream &in)

Private Types

typedef
internal::Triangulation::Iterators
< dim, spacedim > 
IteratorSelector

Private Member Functions

void clear_despite_subscriptions ()
void execute_refinement ()
void execute_coarsening ()
void fix_coarsen_flags ()

Private Attributes

std::vector
< internal::Triangulation::TriaLevel
< dim > * > 
levels
internal::Triangulation::TriaFaces
< dim > * 
faces
std::vector< Point< spacedim > > vertices
std::vector< boolvertices_used
SmartPointer< const Boundary
< dim, spacedim > > 
boundary [255]
SmartPointer< const Boundary
< spacedim, spacedim > > 
manifold_description [255]
bool anisotropic_refinement
MeshSmoothing smooth_grid
internal::Triangulation::NumberCache
< dim > 
number_cache
std::list< RefinementListener * > refinement_listeners

Friends

class TriaAccessorBase
class TriaAccessor
class CellAccessor< dim, spacedim >
class internal::TriaAccessor::Implementation
class hp::DoFHandler< dim, spacedim >
class internal::Triangulation::Implementation

Detailed Description

template<int dim, int spacedim = dim>
class Triangulation< dim, spacedim >

Triangulations denote a hierarchy of levels of elements which together form a dim -dimensional manifold in spacedim spatial dimensions (if spacedim is not specified it takes the default value @ spacedim=dim).

Thus, for example, an object of type Triangulation<1,1> (or simply Triangulation<1> since spacedim==dim by default) is used to represent and handle the usual one-dimensional triangulation used in the finite element method (so, segments on a straight line). On the other hand, objects such as Triangulation<1,2> or Triangulation<2,3> (that are associated with curves in 2D or surfaces in 3D) are the ones one wants to use in the boundary element method.

This class is written to be as independent of the dimension as possible (thus the complex construction of the internal::Triangulation::TriaLevel classes) to allow code-sharing, to allow reducing the need to mirror changes in the code for one dimension to the code for other dimensions. Nonetheless, some of the functions are dependent of the dimension and there only exist specialized versions for distinct dimensions.

Structure and iterators

The actual data structure of a Triangulation object is rather complex and quite inconvenient if one attempted to operate on it directly, since data is spread over quite a lot of arrays and other places. However, there are ways powerful enough to work on these data structures without knowing their exact relations. This is done through the concept of iterators (see the STL documentation and TriaRawIterator). In order to make things as easy and dimension independent as possible, use of class local typedefs is made, see below.

The Triangulation class provides iterator which enable looping over all lines, cells, etc without knowing the exact representation used to describe them. Their names are typedefs imported from the Iterators class (thus making them local types to this class) and are as follows:

Additionaly, for dim==1, the following identities hold:

 *    typedef raw_line_iterator    raw_cell_iterator;
 *    typedef line_iterator        cell_iterator;
 *    typedef active_line_iterator active_cell_iterator;
 *  

while for dim==2

 *    typedef quad_line_iterator   raw_cell_iterator;    
 *    typedef quad_iterator        cell_iterator;
 *    typedef active_quad_iterator active_cell_iterator;
 *
 *    typedef raw_line_iterator    raw_face_iterator;
 *    typedef line_iterator        face_iterator;
 *    typedef active_line_iterator active_face_iterator;    
 *  

By using the cell iterators, you can write code nearly independent of the spatial dimension. The same applies for substructure iterators, where a substructure is defined as a face of a cell. The face of a cell is a vertex in 1D and a line in 2D; however, vertices are handled in a different way and therefore lines have no faces.

The Triangulation class offers functions like begin_active which gives you an iterator to the first active cell. There are quite a lot of functions returning iterators. Take a look at the class doc to get an overview.

Usage of these iterators works mostly like with the STL iterators. Some examples taken from the Triangulation source code follow (notice that in the last two examples the template parameter spacedim has been omitted, so it takes the default value @ dim).

Usage

Usage of a Triangulation is mainly done through the use of iterators. An example probably shows best how to use it:

 *  void main () {
 *    Triangulation<2> tria;
 *
 *    // read in a coarse grid file
 *
 *                                     // we want to log the
 *                                     // refinement history
 *    ofstream history ("mesh.history");
 *    
 *                                     // refine first cell
 *    tria.begin_active()->set_refine_flag();
 *    tria.save_refine_flags (history);
 *    tria.execute_coarsening_and_refinement ();
 *
 *                                     // refine first active cell
 *                                     // on coarsest level
 *    tria.begin_active()->set_refine_flag ();
 *    tria.save_refine_flags (history);
 *    tria.execute_coarsening_and_refinement ();
 *    
 *    Triangulation<2>::active_cell_iterator cell;
 *    for (int i=0; i<17; ++i) 
 *      {
 *                                         // refine the presently
 *                                         // second last cell 17
 *                                         // times
 *        cell = tria.last_active(tria.n_levels()-1);
 *        --cell;
 *        cell->set_refine_flag ();
 *        tria.save_refine_flags (history);
 *        tria.execute_coarsening_and_refinement ();
 *      };
 *                                       // output the grid
 *    ofstream out("grid.1");
 *    GridOut::write_gnuplot (tria, out);
 *  };  
 *  

Creating a triangulation

There are several possibilities to create a triangulation:

Finally, there is a special function for folks who like bad grids: Triangulation<dim,spacedim>::distort_random. It moves all the vertices in the grid a bit around by a random value, leaving behind a distorted mesh. Note that you should apply this function to the final mesh, since refinement smoothes the mesh a bit.

The function will make sure that vertices on restricted faces (hanging nodes) will result in the correct place, i.e. in the middle of the two other vertices of the mother line, and the analogue in higher space dimensions (vertices on the boundary are not corrected, so don't distort boundary vertices in more than two space dimension, i.e. in dimensions where boundary vertices can be hanging nodes). Applying the algorithm has another drawback related to the placement of cells, however: the children of a cell will not occupy the same region of the domain as the mother cell does. While this is the usual behavior with cells at the boundary, here you may get into trouble when using multigrid algorithms or when transferring solutions from coarse to fine grids and back. In general, the use of this function is only safe if you only use the most refined level of the triangulation for computations.

Refinement and coarsening of a triangulation

Refinement of a triangulation may be done through several ways. The most low-level way is directly through iterators: let i be an iterator to an active cell (i.e. the cell pointed to has no children), then the function call i->set_refine_flag() marks the respective cell for refinement. Marking non-active cells results in an error.

After all the cells you wanted to mark for refinement, call the execute_coarsening_and_refinement function to actually perform the refinement. This function itself first calls the prepare_coarsening_and_refinement function to regularize the resulting triangulation: since a face between two adjacent cells may only be subdivided once (i.e. the levels of two adjacent cells may differ by one at most; it is not possible to have a cell refined twice while the neighboring one is not refined), some additional cells are flagged for refinement to smooth the grid. This enlarges the number of resulting cells but makes the grid more regular, thus leading to better approximation properties and, above all, making the handling of data structures and algorithms much much easier. To be honest, this is mostly an algorithmic step than one needed by the finite element method.

To coarsen a grid, the same way as above is possible by using i->set_coarsen_flag and calling execute_coarsening_and_refinement.

The reason for first coarsening, then refining is that the refinement usually adds some additional cells to keep the triangulation regular and thus satisfies all refinement requests, while the coarsening does not delete cells not requested for; therefore the refinement will often revert some effects of coarsening while the opposite is not true. The stated order of coarsening before refinement will thus normally lead to a result closer to the intended one.

Marking cells for refinement 'by hand' through iterators is one way to produce a new grid, especially if you know what kind of grid you are looking for, e.g. if you want to have a grid successively refined towards the boundary or always at the center (see the example programs, they do exactly these things). There are more advanced functions, however, which are more suitable for automatic generation of hierarchical grids in the context of a-posteriori error estimation and adaptive finite elements. These functions can be found in the GridRefinement class.

Smoothing of a triangulation

Some degradation of approximation properties has been observed for grids which are too unstructured. Therefore, the prepare_coarsening_and_refinement function which is automatically called by the execute_coarsening_and_refinement function can do some smoothing of the triangulation. Note that mesh smoothing is only done for two or more space dimensions, no smoothing is available at present for one spatial dimension. In the sequel, let execute_* stand for execute_coarsening_and_refinement.

For the purpose of smoothing, the Triangulation constructor takes an argument specifying whether a smoothing step shall be performed on the grid each time execute_* is called. The default is that such a step not be done, since this results in additional cells being produced, which may not be necessary in all cases. If switched on, calling execute_* results in flagging additional cells for refinement to avoid vertices as the ones mentioned. The algorithms for both regularization and smoothing of triangulations are described below in the section on technical issues. The reason why this parameter must be given to the constructor rather than to execute_* is that it would result in algorithmic problems if you called execute_* once without and once with smoothing, since then in some refinement steps would need to be refined twice.

The parameter taken by the constructor is an integer which may be composed bitwise by the constants defined in the enum MeshSmoothing. The meaning of these constants is explained in the following:

Material and boundary information

Each line, quad, etc stores one byte of information denoting the material or the part of the boundary that an object belongs to. The material of a cell may be used during matrix generation in order to implement different coefficients in different parts of the domain. It is not used by functions of the grid and dof handling libraries.

This material_id may be set upon construction of a triangulation (through the CellData data structure), or later through use of cell iterators. For a typical use of this functionality, see the step-28 tutorial program. The functions of the GridGenerator namespace typically set the material ID of all cells to zero. When reading a triangulation, the material id must be specified in the input file (UCD format) or is otherwise set to zero. Material IDs are inherited by child cells from their parent upon mesh refinement.

Boundary indicators on lower dimensional objects (these have no material id) indicate the number of a boundary component. These are used for two purposes: First, they specify a boundary curve. When a cell is refined, a function can be used to place new vertices on this curve. See the section on boundary approximation below. Furthermore, the weak formulation of the partial differential equation may have different boundary conditions on different parts of the boundary. The boundary indicator can be used in creating the matrix or the right hand side vector to indicate these different parts of the model (this use is like the material id of cells).

Material and boundary indicators may be in the range from zero to 254. The value 255 is reserved to denote interior lines (in 2D) and interior lines and quads (in 3D), which do not have a boundary or material indicator. This way, a program can easily determine, whether such an object is at the boundary or not.

Lines in two dimensions and quads in three dimensions inherit their boundary indicator to their children upon refinement. You should therefore make sure that if you have different boundary parts, the different parts are separated by a vertex (in 2D) or a line (in 3D) such that each boundary line or quad has a unique boundary indicator.

Since in one dimension, no substructures of lower dimension exist to cells (of course apart from vertices, but these are handled in another way than the structures and substructures with dimension one or greater), there is no way to denote boundary indicators to boundary vertices (the endpoints). This is not a big thing, however, since you will normally not want to do a loop over all vertices, but rather work on the cells, which in turn have a possibility to find out whether they are at one of the endpoints. Only the handling of boundary values gets a bit more complicated, but this seems to be the price to be paid for the different handling of vertices from lines and quads.

History of a triangulation

It is possible to reconstruct a grid from its refinement history, which can be stored and loaded through the save_refine_flags and load_refine_flags functions. Normally, the code will look like this:

 *                                 // open output file
 *     ofstream history("mesh.history");
 *                                 // do 10 refinement steps
 *     for (int step=0; step<10; ++step) {
 *       ...;
 *       // flag cells according to some criterion
 *       ...;
 *       tria.save_refine_flags (history);
 *       tria.execute_coarsening_and_refinement ();
 *     };        
 *   

If you want to re-create the grid from the stored information, you write:

 *                                 // open input file
 *     ifstream history("mesh.history");
 *                                 // do 10 refinement steps
 *     for (int step=0; step<10; ++step) {
 *       tria.load_refine_flags (history);
 *       tria.execute_coarsening_and_refinement ();
 *     };        
 *   

The same scheme is employed for coarsening and the coarsening flags.

You may write other information to the output file between different sets of refinement information, as long as you read it upon re-creation of the grid. You should make sure that the other information in the new triangulation which is to be created from the saved flags, matches that of the old triangulation, for example the smoothing level; if not, the cells actually created from the flags may be other ones, since smoothing adds additional cells, but their number may be depending on the smoothing level.

There actually are two sets of save_*_flags and load_*_flags functions. One takes a stream as argument and reads/writes the information from/to the stream, thus enabling storing flags to files. The other set takes an argument of type vector<bool>. This enables the user to temporarily store some flags, e.g. if another function needs them, and restore them afterwards.

User flags and data

A triangulation offers one bit per line, quad, etc for user data. This field can be accessed as all other data using iterators. Normally, this user flag is used if an algorithm walks over all cells and needs information whether another cell, e.g. a neighbor, has already been processed. It can also be used to flag the lines subject to constraints in 2D, as for example the functions in the DoFHandler classes do.

There are two functions, save_user_flags and load_user_flags which write and read these flags to and from a stream. Unlike save_refine_flags and load_refine_flags, these two functions store and read the flags of all used lines, quads, etc, not only of the active ones (well, activity is a concept which really only applies to cells, not for example to lines in 2D, so the abovementioned generalization to all lines, quads, etc seems plausible).

If you want to store more specific user flags, you can use the functions save_user_flags_line and load_user_flags_line and the generalizations for quads, etc.

As for the refinement and coarsening flags, there exist two versions of these functions, one which reads/writes from a stream and one which does so from a vector<bool>. The latter is used to store flags temporarily, while the first is used to store them in a file.

It is convention to clear the user flags using the Triangulation<>::clear_user_flags() function before usage, since it is often necessary to use the flags in more than one function consecutively and is then error prone to dedicate one of these to clear the flags.

It is recommended that a functions using the flags states so in its documentation. For example, the execute_coarsening_and_refinement() function uses some of the user flags and will therefore destroy any content stored in them.

There is another set of user data, which can be either an unsigned int or a void *, for each line, quad, etc. You can access these through the functions listed under User data in the accessor classes. These pointers are not used nor changed in many places of the library, and those classes and functions that do use them should document this clearly; the most prominent user of these pointers is the SolutionTransfer class which uses the cell->user_pointers.

The value of these user indices or pointers is NULL by default. Note that the pointers are not inherited to children upon refinement. Still, after a remeshing they are available on all cells, where they were set on the previous mesh (unless, of course, overwritten by the SolutionTransfer class).

The usual warning about the missing type safety of void pointers are obviously in place here; responsibility for correctness of types etc lies entirely with the user of the pointer.

Just like the user flags, this field is not available for vertices, which does no harm since the vertices have a unique and continuous number unlike the structured objects lines and quads.

Boundary approximation

You can specify a boundary function for each boundary component. 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. The boundary indicator of the face will be used to determine the proper component. See Boundary for the details. Usage with the Triangulation object is then like this (let Ball be a class derived from Boundary<2>):

 *     void main () {
 *       Triangulation<2> tria;
 *                                        // set the boundary function
 *                                        // for all boundaries with
 *                                        // boundary indicator 0
 *       Ball ball;
 *       tria.set_boundary (0, ball);
 *
 *       // read some coarse grid
 *
 *
 *       Triangulation<2>::active_cell_iterator cell, endc;
 *       for (int i=0; i<8; ++i) 
 *         {
 *           cell = tria.begin_active();
 *           endc = tria.end();
 *
 *                                            // refine all
 *                                            // boundary cells
 *           for (; cell!=endc; ++cell)
 *             if (cell->at_boundary())
 *               cell->set_refine_flag();
 *
 *           tria.execute_coarsening_and_refinement();
 *         };
 *     };            
 *   

You should take note of one caveat: if you have concave boundaries, you must make sure that a new boundary vertex does not lie too much inside the cell which is to be refined. The reason is that the center vertex is placed at the point which is the arithmetic mean of the vertices of the original cell. Therefore if your new boundary vertex is too near the center of the old quadrilateral or hexahedron, the distance to the midpoint vertex will become too small, thus generating distorted cells. Remedy: take care of such situations when defining the coarse grid.

Technical details

Algorithms for mesh regularization and smoothing upon refinement

We chose an inductive point of view: since upon creation of the triangulation all cells are on the same level, all regularity assumptions regarding the maximum difference in level of cells sharing a common face, edge or vertex hold. Since we use the regularization and smoothing in each step of the mesh history, when coming to the point of refining it further the assumptions also hold.

The regularization and smoothing is done in the prepare_coarsening_and_refinement function, which is called by execute_coarsening_and_refinement at the very beginning. It decides which additional cells to flag for refinement by looking at the old grid and the refinement flags for each cell.

Regularization and smoothing are a bit complementary in that we check whether we need to set additional refinement flags when being on a cell flagged for refinement (regularization) or on a cell not flagged for refinement. This makes readable programming easier.

All the described algorithms apply only for more than one space dimension, since for one dimension no restrictions apply. It may be necessary to apply some smoothing for multigrid algorithms, but this has to be decided upon later.

Warning

It seems impossible to preserve constness of a triangulation through iterator usage. Thus, if you declare pointers to a const triangulation object, you should be well aware that you might involuntarily alter the data stored in the triangulation.

Author:
Wolfgang Bangerth, 1998; Ralf Hartmann, 2005

Member Typedef Documentation

template<int dim, int spacedim = dim>
typedef internal::Triangulation::Iterators<dim, spacedim> Triangulation< dim, spacedim >::IteratorSelector [private]

An internal typedef to make the definition of the iterator classes simpler.

template<int dim, int spacedim = dim>
typedef IteratorSelector::raw_line_iterator Triangulation< dim, spacedim >::raw_line_iterator
template<int dim, int spacedim = dim>
typedef IteratorSelector::line_iterator Triangulation< dim, spacedim >::line_iterator
template<int dim, int spacedim = dim>
typedef IteratorSelector::active_line_iterator Triangulation< dim, spacedim >::active_line_iterator
template<int dim, int spacedim = dim>
typedef IteratorSelector::raw_quad_iterator Triangulation< dim, spacedim >::raw_quad_iterator
template<int dim, int spacedim = dim>
typedef IteratorSelector::quad_iterator Triangulation< dim, spacedim >::quad_iterator
template<int dim, int spacedim = dim>
typedef IteratorSelector::active_quad_iterator Triangulation< dim, spacedim >::active_quad_iterator
template<int dim, int spacedim = dim>
typedef IteratorSelector::raw_hex_iterator Triangulation< dim, spacedim >::raw_hex_iterator
template<int dim, int spacedim = dim>
typedef IteratorSelector::hex_iterator Triangulation< dim, spacedim >::hex_iterator
template<int dim, int spacedim = dim>
typedef IteratorSelector::active_hex_iterator Triangulation< dim, spacedim >::active_hex_iterator
template<int dim, int spacedim = dim>
typedef IteratorSelector::raw_cell_iterator Triangulation< dim, spacedim >::raw_cell_iterator
template<int dim, int spacedim = dim>
typedef IteratorSelector::cell_iterator Triangulation< dim, spacedim >::cell_iterator
template<int dim, int spacedim = dim>
typedef IteratorSelector::active_cell_iterator Triangulation< dim, spacedim >::active_cell_iterator
template<int dim, int spacedim = dim>
typedef IteratorSelector::raw_face_iterator Triangulation< dim, spacedim >::raw_face_iterator
template<int dim, int spacedim = dim>
typedef IteratorSelector::face_iterator Triangulation< dim, spacedim >::face_iterator
template<int dim, int spacedim = dim>
typedef IteratorSelector::active_face_iterator Triangulation< dim, spacedim >::active_face_iterator

Member Enumeration Documentation

template<int dim, int spacedim = dim>
enum Triangulation::MeshSmoothing

Declare some symbolic names for mesh smoothing algorithms. The meaning of these flags is documented in the Triangulation class.

Enumerator:
none 
limit_level_difference_at_vertices 
eliminate_unrefined_islands 
patch_level_1 
coarsest_level_1 
allow_anisotropic_smoothing 
eliminate_refined_inner_islands 
eliminate_refined_boundary_islands 
do_not_produce_unrefined_islands 
smoothing_on_refinement 
smoothing_on_coarsening 
maximum_smoothing 

Constructor & Destructor Documentation

template<int dim, int spacedim = dim>
Triangulation< dim, spacedim >::Triangulation ( const MeshSmoothing  smooth_grid = none  ) 

Create an empty triangulation. Do not create any cells.

template<int dim, int spacedim = dim>
Triangulation< dim, spacedim >::Triangulation ( const Triangulation< dim, spacedim > &  t  ) 

Copy constructor.

You should really use the copy_triangulation function, so we declare this function but let it throw an internal error. The reason for this is that we may want to use triangulation objects in collections. However, C++ containers require that the objects stored in them are copyable, so we need to provide a copy constructor. On the other hand, copying triangulations is so expensive that we do not want such objects copied by accident, for example in compiler-generated temporary objects. By defining a copy constructor but throwing an error, we satisfy the formal requirements of containers, but at the same time disallow actual copies. Finally, through the exception, one easily finds the places where code has to be changed to avoid copies.

template<int dim, int spacedim = dim>
virtual Triangulation< dim, spacedim >::~Triangulation (  )  [virtual]

Delete the object and all levels of the hierarchy.


Member Function Documentation

template<int dim, int spacedim = dim>
void Triangulation< dim, spacedim >::clear (  ) 

Reset this triangulation into a virgin state by deleting all data.

Note that this operation is only allowed if no subscriptions to this object exist any more, such as DoFHandler objects using it.

template<int dim, int spacedim = dim>
void Triangulation< dim, spacedim >::set_boundary ( const unsigned int  number,
const Boundary< dim, spacedim > &  boundary_object 
)

If dim==spacedim, assign a boundary object to a certain part of the boundary of a the triangulation. If a face with boundary number number is refined, this object is used to find the location of new vertices on the boundary. It is also used for for non-linear (i.e.: non-Q1) transformations of cells to the unit cell in shape function calculations.

If dim!=spacedim the boundary object is in fact the exact manifold that the triangulation is approximating (for example a circle approximated by a polygon triangulation). As above, the refinement is made in such a way that the new points are located on the exact manifold.

Numbers of boundary objects 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. Material IDs on boundaries are also called boundary indicators and are accessed with accessor functions of that name.

The boundary_object is not copied and MUST persist until the triangulation is destroyed. 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. Usually, this is done before the first refinement and is dangerous afterwards. Removal of a boundary object is done by set_boundary(number), i.e. the function of same name but only one argument. This operation then replaces the boundary object given before by a straight boundary approximation.

template<int dim, int spacedim = dim>
void Triangulation< dim, spacedim >::set_boundary ( const unsigned int  number  ) 

Reset those parts of the boundary with the given number to use a straight boundary approximation. This is the default state of a triangulation, and undoes assignment of a different boundary object by the function of same name and two arguments.

template<int dim, int spacedim = dim>
const Boundary<dim,spacedim>& Triangulation< dim, spacedim >::get_boundary ( const unsigned int  number  )  const

Return a constant reference to a boundary object used for this triangulation. Number is the same as in set_boundary

template<int dim, int spacedim = dim>
std::vector<unsigned char> Triangulation< dim, spacedim >::get_boundary_indicators (  )  const

Returns a vector containing all boundary indicators assigned to boundary faces of this Triangulation object. Note, that each boundary indicator is reported only once. The size of the return vector will represent the number of different indicators (which is greater or equal one).

template<int dim, int spacedim = dim>
virtual void Triangulation< dim, spacedim >::copy_triangulation ( const Triangulation< dim, spacedim > &  old_tria  )  [virtual]

Copy a triangulation. This operation is not cheap, so you should be careful with using this. We do not implement this function as a copy constructor, since it makes it easier to maintain collections of triangulations if you can assign them values later on.

Keep in mind that this function also copies the pointer to the boundary descriptor previously set by the set_boundary function. You must therefore also guarantee that the boundary objects has a lifetime at least as long as the copied triangulation.

This triangulation must be empty beforehand.

The function is made virtual since some derived classes might want to disable or extend the functionality of this function.

Note, that the list of RefinementListeners is not copied. However, each RefinementListener is notified of the copy operation through the RefinementListener::copy_notification() function, so it might subscribe to the new Triangulation as well, if that is desired.

template<int dim, int spacedim = dim>
virtual void Triangulation< dim, spacedim >::create_triangulation ( const std::vector< Point< spacedim > > &  vertices,
const std::vector< CellData< dim > > &  cells,
const SubCellData subcelldata 
) [virtual]

Create a triangulation from a list of vertices and a list of cells, each of the latter being a list of 1<<dim vertex indices. The triangulation must be empty upon calling this function and the cell list should be useful (connected domain, etc.).

Material data for the cells is given within the cells array, while boundary information is given in the subcelldata field.

The numbering of vertices within the cells array is subject to some constraints; see the general class documentation for this.

This function is made virtual to allow derived classes to set up some data structures as well.

For conditions when this function can generate a valid triangulation, see the documentation of this class, and the GridIn and GridReordering class.

Reimplemented in PersistentTriangulation< dim >.

template<int dim, int spacedim = dim>
virtual void Triangulation< dim, spacedim >::create_triangulation_compatibility ( const std::vector< Point< spacedim > > &  vertices,
const std::vector< CellData< dim > > &  cells,
const SubCellData subcelldata 
) [virtual]

For backward compatibility, only. This function takes the cell data in the ordering as requested by deal.II versions up to 5.2, converts it to the new (lexicographic) ordering and calls create_triangulation().

Reimplemented in PersistentTriangulation< dim >.

template<int dim, int spacedim = dim>
void Triangulation< dim, spacedim >::distort_random ( const double  factor,
const bool  keep_boundary = true 
)

Distort the grid by randomly moving around all the vertices of the grid. The direction of moving is random, while the length of the shift vector has a value of factor times the minimal length of the active lines adjacent to this vertex. Note that factor should obviously be well below 0.5.

If keep_boundary is set to true (which is the default), then boundary vertices are not moved.

template<int dim, int spacedim = dim>
void Triangulation< dim, spacedim >::set_all_refine_flags (  ) 

Flag all active cells for refinement. This will refine all cells of all levels which are not already refined (i.e. only cells are refined which do not yet have children). The cells are only flagged, not refined, thus you have the chance to save the refinement flags.

template<int dim, int spacedim = dim>
void Triangulation< dim, spacedim >::refine_global ( const unsigned int  times  ) 

Refine all cells times times, by alternatingly calling set_all_refine_flags() and execute_coarsening_and_refinement(). This function actually starts the refinement process, so you have no way to store the refinement flags unless you overload the execute_coarsening_and_refinement function.

template<int dim, int spacedim = dim>
virtual void Triangulation< dim, spacedim >::execute_coarsening_and_refinement (  )  [virtual]

Execute both refinement and coarsening of the triangulation.

The function resets all refinement and coarsening flags to false. It uses the user flags for internal purposes. They will therefore be overwritten by undefined content.

See the general docs for more information.

Note that this function is virtual to allow derived classes to insert hooks, such as saving refinement flags and the like (see e.g. the PersistentTriangulation class).

Reimplemented in PersistentTriangulation< dim >.

template<int dim, int spacedim = dim>
bool Triangulation< dim, spacedim >::prepare_coarsening_and_refinement (  ) 

Do both preparation for refinement and coarsening as well as mesh smoothing.

Regarding the refinement process it fixes the closure of the refinement in dim>=2 (make sure that no two cells are adjacent with a refinement level differing with more than one), etc. It performs some mesh smoothing if the according flag was given to the constructor of this class. The function returns whether additional cells have been flagged for refinement.

See the general doc of this class for more information on smoothing upon refinement.

Regarding the coarsening part, flagging and deflagging cells in preparation of the actual coarsening step are done. This includes deleting coarsen flags from cells which may not be deleted (e.g. because one neighbor is more refined than the cell), doing some smoothing, etc.

The effect is that only those cells are flagged for coarsening which will actually be coarsened. This includes the fact that all flagged cells belong to parent cells of which all children are flagged.

The function returns whether some cells' flagging has been changed in the process.

This function uses the user flags, so store them if you still need them afterwards.

template<int dim, int spacedim = dim>
void Triangulation< dim, spacedim >::add_refinement_listener ( RefinementListener listener  )  const

Add a RefinementListener. Adding listeners to the Triangulation allows other classes to be informed, when the Triangulation is refined.

template<int dim, int spacedim = dim>
void Triangulation< dim, spacedim >::remove_refinement_listener ( RefinementListener listener  )  const

Removes a RefinementListener. When some class needs no longer to be informed about refinements, the listener should be removed from the Triangulation.

template<int dim, int spacedim = dim>
void Triangulation< dim, spacedim >::save_refine_flags ( std::ostream &  out  )  const

Save the addresses of the cells which are flagged for refinement to out. For usage, read the general documentation for this class.

template<int dim, int spacedim = dim>
void Triangulation< dim, spacedim >::save_refine_flags ( std::vector< bool > &  v  )  const

Same as above, but store the flags to a bitvector rather than to a file.

template<int dim, int spacedim = dim>
void Triangulation< dim, spacedim >::load_refine_flags ( std::istream &  in  ) 

Read the information stored by save_refine_flags.

template<int dim, int spacedim = dim>
void Triangulation< dim, spacedim >::load_refine_flags ( const std::vector< bool > &  v  ) 

Read the information stored by save_refine_flags.

template<int dim, int spacedim = dim>
void Triangulation< dim, spacedim >::save_coarsen_flags ( std::ostream &  out  )  const

Analogue to save_refine_flags.

template<int dim, int spacedim = dim>
void Triangulation< dim, spacedim >::save_coarsen_flags ( std::vector< bool > &  v  )  const

Same as above, but store the flags to a bitvector rather than to a file.

template<int dim, int spacedim = dim>
void Triangulation< dim, spacedim >::load_coarsen_flags ( std::istream &  out  ) 

Analogue to load_refine_flags.

template<int dim, int spacedim = dim>
void Triangulation< dim, spacedim >::load_coarsen_flags ( const std::vector< bool > &  v  ) 

Analogue to load_refine_flags.

template<int dim, int spacedim = dim>
bool Triangulation< dim, spacedim >::get_anisotropic_refinement_flag (  )  const

Return whether this triangulation has ever undergone anisotropic (as opposed to only isotropic) refinement.

template<int dim, int spacedim = dim>
void Triangulation< dim, spacedim >::clear_user_flags (  ) 

Clear all user flags.

template<int dim, int spacedim = dim>
void Triangulation< dim, spacedim >::save_user_flags ( std::ostream &  out  )  const

Save all user flags. See the general documentation for this class and the documentation for the save_refine_flags for more details.

template<int dim, int spacedim = dim>
void Triangulation< dim, spacedim >::save_user_flags ( std::vector< bool > &  v  )  const

Same as above, but store the flags to a bitvector rather than to a file. The output vector is resized if necessary.

template<int dim, int spacedim = dim>
void Triangulation< dim, spacedim >::load_user_flags ( std::istream &  in  ) 

Read the information stored by save_user_flags.

template<int dim, int spacedim = dim>
void Triangulation< dim, spacedim >::load_user_flags ( const std::vector< bool > &  v  ) 

Read the information stored by save_user_flags.

template<int dim, int spacedim = dim>
void Triangulation< dim, spacedim >::clear_user_flags_line (  ) 

Clear all user flags on lines.

template<int dim, int spacedim = dim>
void Triangulation< dim, spacedim >::save_user_flags_line ( std::ostream &  out  )  const

Save the user flags on lines.

template<int dim, int spacedim = dim>
void Triangulation< dim, spacedim >::save_user_flags_line ( std::vector< bool > &  v  )  const

Same as above, but store the flags to a bitvector rather than to a file. The output vector is resized if necessary.

template<int dim, int spacedim = dim>
void Triangulation< dim, spacedim >::load_user_flags_line ( std::istream &  in  ) 

Load the user flags located on lines.

template<int dim, int spacedim = dim>
void Triangulation< dim, spacedim >::load_user_flags_line ( const std::vector< bool > &  v  ) 

Load the user flags located on lines.

template<int dim, int spacedim = dim>
void Triangulation< dim, spacedim >::clear_user_flags_quad (  ) 

Clear all user flags on quads.

template<int dim, int spacedim = dim>
void Triangulation< dim, spacedim >::save_user_flags_quad ( std::ostream &  out  )  const

Save the user flags on quads.

template<int dim, int spacedim = dim>
void Triangulation< dim, spacedim >::save_user_flags_quad ( std::vector< bool > &  v  )  const

Same as above, but store the flags to a bitvector rather than to a file. The output vector is resized if necessary.

template<int dim, int spacedim = dim>
void Triangulation< dim, spacedim >::load_user_flags_quad ( std::istream &  in  ) 

Load the user flags located on quads.

template<int dim, int spacedim = dim>
void Triangulation< dim, spacedim >::load_user_flags_quad ( const std::vector< bool > &  v  ) 

Load the user flags located on quads.

template<int dim, int spacedim = dim>
void Triangulation< dim, spacedim >::clear_user_flags_hex (  ) 

Clear all user flags on quads.

template<int dim, int spacedim = dim>
void Triangulation< dim, spacedim >::save_user_flags_hex ( std::ostream &  out  )  const

Save the user flags on hexs.

template<int dim, int spacedim = dim>
void Triangulation< dim, spacedim >::save_user_flags_hex ( std::vector< bool > &  v  )  const

Same as above, but store the flags to a bitvector rather than to a file. The output vector is resized if necessary.

template<int dim, int spacedim = dim>
void Triangulation< dim, spacedim >::load_user_flags_hex ( std::istream &  in  ) 

Load the user flags located on hexs.

template<int dim, int spacedim = dim>
void Triangulation< dim, spacedim >::load_user_flags_hex ( const std::vector< bool > &  v  ) 

Load the user flags located on hexs.

template<int dim, int spacedim = dim>
void Triangulation< dim, spacedim >::clear_user_data (  ) 

Clear all user pointers and indices and allow the use of both for next access.

template<int dim, int spacedim = dim>
void Triangulation< dim, spacedim >::clear_user_pointers (  ) 
Deprecated:
User clear_user_data() instead.

Clear all user pointers.

template<int dim, int spacedim = dim>
void Triangulation< dim, spacedim >::save_user_indices ( std::vector< unsigned int > &  v  )  const

Save all user indices. The output vector is resized if necessary.

template<int dim, int spacedim = dim>
void Triangulation< dim, spacedim >::load_user_indices ( const std::vector< unsigned int > &  v  ) 

Read the information stored by save_user_indices().

template<int dim, int spacedim = dim>
void Triangulation< dim, spacedim >::save_user_pointers ( std::vector< void * > &  v  )  const

Save all user pointers. The output vector is resized if necessary.

template<int dim, int spacedim = dim>
void Triangulation< dim, spacedim >::load_user_pointers ( const std::vector< void * > &  v  ) 

Read the information stored by save_user_pointers().

template<int dim, int spacedim = dim>
void Triangulation< dim, spacedim >::save_user_indices_line ( std::vector< unsigned int > &  v  )  const

Save the user indices on lines. The output vector is resized if necessary.

template<int dim, int spacedim = dim>
void Triangulation< dim, spacedim >::load_user_indices_line ( const std::vector< unsigned int > &  v  ) 

Load the user indices located on lines.

template<int dim, int spacedim = dim>
void Triangulation< dim, spacedim >::save_user_indices_quad ( std::vector< unsigned int > &  v  )  const

Save the user indices on quads. The output vector is resized if necessary.

template<int dim, int spacedim = dim>
void Triangulation< dim, spacedim >::load_user_indices_quad ( const std::vector< unsigned int > &  v  ) 

Load the user indices located on quads.

template<int dim, int spacedim = dim>
void Triangulation< dim, spacedim >::save_user_indices_hex ( std::vector< unsigned int > &  v  )  const

Save the user indices on hexes. The output vector is resized if necessary.

template<int dim, int spacedim = dim>
void Triangulation< dim, spacedim >::load_user_indices_hex ( const std::vector< unsigned int > &  v  ) 

Load the user indices located on hexs.

template<int dim, int spacedim = dim>
void Triangulation< dim, spacedim >::save_user_pointers_line ( std::vector< void * > &  v  )  const

Save the user indices on lines. The output vector is resized if necessary.

template<int dim, int spacedim = dim>
void Triangulation< dim, spacedim >::load_user_pointers_line ( const std::vector< void * > &  v  ) 

Load the user pointers located on lines.

template<int dim, int spacedim = dim>
void Triangulation< dim, spacedim >::save_user_pointers_quad ( std::vector< void * > &  v  )  const

Save the user pointers on quads. The output vector is resized if necessary.

template<int dim, int spacedim = dim>
void Triangulation< dim, spacedim >::load_user_pointers_quad ( const std::vector< void * > &  v  ) 

Load the user pointers located on quads.

template<int dim, int spacedim = dim>
void Triangulation< dim, spacedim >::save_user_pointers_hex ( std::vector< void * > &  v  )  const

Save the user pointers on hexes. The output vector is resized if necessary.

template<int dim, int spacedim = dim>
void Triangulation< dim, spacedim >::load_user_pointers_hex ( const std::vector< void * > &  v  ) 

Load the user pointers located on hexs.

template<int dim, int spacedim = dim>
raw_cell_iterator Triangulation< dim, spacedim >::begin_raw ( const unsigned int  level = 0  )  const

Iterator to the first cell, used or not, on level level. If a level has no cells, a past-the-end iterator is returned.

This function calls begin_raw_line() in 1D and begin_raw_quad() in 2D.

template<int dim, int spacedim = dim>
cell_iterator Triangulation< dim, spacedim >::begin ( const unsigned int  level = 0  )  const

Iterator to the first used cell on level level.

This function calls begin_line in 1D and begin_quad in 2D.

template<int dim, int spacedim = dim>
active_cell_iterator Triangulation< dim, spacedim >::begin_active ( const unsigned int  level = 0  )  const

Iterator to the first active cell on level level.

This function calls begin_active_line in 1D and begin_active_quad in 2D.

Referenced by GridTools::transform().

template<int dim, int spacedim = dim>
raw_cell_iterator Triangulation< dim, spacedim >::end (  )  const

Iterator past the end; this iterator serves for comparisons of iterators with past-the-end or before-the-beginning states.

This function calls end_line in 1D and end_quad in 2D.

Referenced by GridTools::transform().

template<int dim, int spacedim = dim>
cell_iterator Triangulation< dim, spacedim >::end ( const unsigned int  level  )  const

Return an iterator which is the first iterator not on level. If level is the last level, then this returns end().

template<int dim, int spacedim = dim>
raw_cell_iterator Triangulation< dim, spacedim >::end_raw ( const unsigned int  level  )  const

Return a raw iterator which is the first iterator not on level. If level is the last level, then this returns end().

template<int dim, int spacedim = dim>
active_cell_iterator Triangulation< dim, spacedim >::end_active ( const unsigned int  level  )  const

Return an active iterator which is the first iterator not on level. If level is the last level, then this returns end().

template<int dim, int spacedim = dim>
raw_cell_iterator Triangulation< dim, spacedim >::last_raw (  )  const

Return an iterator pointing to the last cell, used or not.

This function calls last_raw_line in 1D and last_raw_quad in 2D.

template<int dim, int spacedim = dim>
raw_cell_iterator Triangulation< dim, spacedim >::last_raw ( const unsigned int  level  )  const

Return an iterator pointing to the last cell of the level level, used or not.

This function calls last_raw_line in 1D and last_raw_quad in 2D.

template<int dim, int spacedim = dim>
cell_iterator Triangulation< dim, spacedim >::last (  )  const

Return an iterator pointing to the last used cell.

This function calls last_line in 1D and last_quad in 2D.

template<int dim, int spacedim = dim>
cell_iterator Triangulation< dim, spacedim >::last ( const unsigned int  level  )  const

Return an iterator pointing to the last used cell on level level.

This function calls last_line in 1D and last_quad in 2D.

template<int dim, int spacedim = dim>
active_cell_iterator Triangulation< dim, spacedim >::last_active (  )  const

Return an iterator pointing to the last active cell.

This function calls last_active_line in 1D and last_active_quad in 2D.

template<int dim, int spacedim = dim>
active_cell_iterator Triangulation< dim, spacedim >::last_active ( const unsigned int  level  )  const

Return an iterator pointing to the last active cell on level level.

This function calls last_active_line in 1D and last_active_quad in 2D.

template<int dim, int spacedim = dim>
raw_face_iterator Triangulation< dim, spacedim >::begin_raw_face (  )  const

Iterator to the first face, used or not. As faces have no level, no argument can be given.

This function calls begin_raw_line in 2D and begin_raw_quad in 3D.

template<int dim, int spacedim = dim>
face_iterator Triangulation< dim, spacedim >::begin_face (  )  const

Iterator to the first used face.

This function calls begin_line in 2D and begin_quad in 3D.

template<int dim, int spacedim = dim>
active_face_iterator Triangulation< dim, spacedim >::begin_active_face (  )  const

Iterator to the first active face.

This function calls begin_active_line in 2D and begin_active_quad in 3D.

template<int dim, int spacedim = dim>
raw_face_iterator Triangulation< dim, spacedim >::end_face (  )  const

Iterator past the end; this iterator serves for comparisons of iterators with past-the-end or before-the-beginning states.

This function calls end_line in 2D and end_quad in 3D.

template<int dim, int spacedim = dim>
raw_face_iterator Triangulation< dim, spacedim >::end_raw_face (  )  const

Return a raw iterator which is past the end. This is the same as end() and is only for combatibility with older versions.

template<int dim, int spacedim = dim>
active_face_iterator Triangulation< dim, spacedim >::end_active_face (  )  const

Return an active iterator which is past the end. This is the same as end() and is only for combatibility with older versions.

template<int dim, int spacedim = dim>
raw_face_iterator Triangulation< dim, spacedim >::last_raw_face (  )  const

Return an iterator pointing to the last face, used or not.

This function calls last_raw_line in 2D and last_raw_quad in 3D.

template<int dim, int spacedim = dim>
face_iterator Triangulation< dim, spacedim >::last_face (  )  const

Return an iterator pointing to the last used face.

This function calls last_line in 2D and last_quad in 3D.

template<int dim, int spacedim = dim>
active_face_iterator Triangulation< dim, spacedim >::last_active_face (  )  const

Return an iterator pointing to the last active face.

This function calls last_active_line in 2D and last_active_quad in 3D.

template<int dim, int spacedim = dim>
raw_line_iterator Triangulation< dim, spacedim >::begin_raw_line ( const unsigned int  level = 0  )  const

Iterator to the first line, used or not, on level level. If a level has no lines, a past-the-end iterator is returned. If lines are no cells, i.e. for dim>1 no level argument must be given. The same applies for all the other functions above, of course.

template<int dim, int spacedim = dim>
line_iterator Triangulation< dim, spacedim >::begin_line ( const unsigned int  level = 0  )  const

Iterator to the first used line on level level.

template<int dim, int spacedim = dim>
active_line_iterator Triangulation< dim, spacedim >::begin_active_line ( const unsigned int  level = 0  )  const

Iterator to the first active line on level level.

template<int dim, int spacedim = dim>
raw_line_iterator Triangulation< dim, spacedim >::end_line (  )  const

Iterator past the end; this iterator serves for comparisons of iterators with past-the-end or before-the-beginning states.

template<int dim, int spacedim = dim>
line_iterator Triangulation< dim, spacedim >::end_line ( const unsigned int  level  )  const

Return an iterator which is the first iterator not on level. If level is the last level, then this returns end().

template<int dim, int spacedim = dim>
raw_line_iterator Triangulation< dim, spacedim >::end_raw_line ( const unsigned int  level  )  const

Return a raw iterator which is the first iterator not on level. If level is the last level, then this returns end().

template<int dim, int spacedim = dim>
active_line_iterator Triangulation< dim, spacedim >::end_active_line ( const unsigned int  level  )  const

Return an active iterator which is the first iterator not on level. If level is the last level, then this returns end().

template<int dim, int spacedim = dim>
raw_line_iterator Triangulation< dim, spacedim >::last_raw_line (  )  const

Return an iterator pointing to the last line, used or not.

template<int dim, int spacedim = dim>
raw_line_iterator Triangulation< dim, spacedim >::last_raw_line ( const unsigned int  level  )  const

Return an iterator pointing to the last line of the level level, used or not.

template<int dim, int spacedim = dim>
line_iterator Triangulation< dim, spacedim >::last_line (  )  const

Return an iterator pointing to the last used line.

template<int dim, int spacedim = dim>
line_iterator Triangulation< dim, spacedim >::last_line ( const unsigned int  level  )  const

Return an iterator pointing to the last used line on level level.

template<int dim, int spacedim = dim>
active_line_iterator Triangulation< dim, spacedim >::last_active_line (  )  const

Return an iterator pointing to the last active line.

template<int dim, int spacedim = dim>
active_line_iterator Triangulation< dim, spacedim >::last_active_line ( const unsigned int  level  )  const

Return an iterator pointing to the last active line on level level.

template<int dim, int spacedim = dim>
raw_quad_iterator Triangulation< dim, spacedim >::begin_raw_quad ( const unsigned int  level = 0  )  const

Iterator to the first quad, used or not, on the given level. If a level has no quads, a past-the-end iterator is returned. If quads are no cells, i.e. for $dim>2$ no level argument must be given.

template<int dim, int spacedim = dim>
quad_iterator Triangulation< dim, spacedim >::begin_quad ( const unsigned int  level = 0  )  const

Iterator to the first used quad on level level.

template<int dim, int spacedim = dim>
active_quad_iterator Triangulation< dim, spacedim >::begin_active_quad ( const unsigned int  level = 0  )  const

Iterator to the first active quad on level level.

template<int dim, int spacedim = dim>
raw_quad_iterator Triangulation< dim, spacedim >::end_quad (  )  const

Iterator past the end; this iterator serves for comparisons of iterators with past-the-end or before-the-beginning states.

template<int dim, int spacedim = dim>
quad_iterator Triangulation< dim, spacedim >::end_quad ( const unsigned int  level  )  const

Return an iterator which is the first iterator not on level. If level is the last level, then this returns end().

template<int dim, int spacedim = dim>
raw_quad_iterator Triangulation< dim, spacedim >::end_raw_quad ( const unsigned int  level  )  const

Return a raw iterator which is the first iterator not on level. If level is the last level, then this returns end().

template<int dim, int spacedim = dim>
active_quad_iterator Triangulation< dim, spacedim >::end_active_quad ( const unsigned int  level  )  const

Return an active iterator which is the first iterator not on level. If level is the last level, then this returns end().

template<int dim, int spacedim = dim>
raw_quad_iterator Triangulation< dim, spacedim >::last_raw_quad (  )  const

Return an iterator pointing to the last quad, used or not.

template<int dim, int spacedim = dim>
raw_quad_iterator Triangulation< dim, spacedim >::last_raw_quad ( const unsigned int  level  )  const

Return an iterator pointing to the last quad of the level level, used or not.

template<int dim, int spacedim = dim>
quad_iterator Triangulation< dim, spacedim >::last_quad (  )  const

Return an iterator pointing to the last used quad.

template<int dim, int spacedim = dim>
quad_iterator Triangulation< dim, spacedim >::last_quad ( const unsigned int  level  )  const

Return an iterator pointing to the last used quad on level level.

template<int dim, int spacedim = dim>
active_quad_iterator Triangulation< dim, spacedim >::last_active_quad (  )  const

Return an iterator pointing to the last active quad.

template<int dim, int spacedim = dim>
active_quad_iterator Triangulation< dim, spacedim >::last_active_quad ( const unsigned int  level  )  const

Return an iterator pointing to the last active quad on level level.

template<int dim, int spacedim = dim>
raw_hex_iterator Triangulation< dim, spacedim >::begin_raw_hex ( const unsigned int  level = 0  )  const

Iterator to the first hex, used or not, on level level. If a level has no hexs, a past-the-end iterator is returned.

template<int dim, int spacedim = dim>
hex_iterator Triangulation< dim, spacedim >::begin_hex ( const unsigned int  level = 0  )  const

Iterator to the first used hex on level level.

template<int dim, int spacedim = dim>
active_hex_iterator Triangulation< dim, spacedim >::begin_active_hex ( const unsigned int  level = 0  )  const

Iterator to the first active hex on level level.

template<int dim, int spacedim = dim>
raw_hex_iterator Triangulation< dim, spacedim >::end_hex (  )  const

Iterator past the end; this iterator serves for comparisons of iterators with past-the-end or before-the-beginning states.

template<int dim, int spacedim = dim>
hex_iterator Triangulation< dim, spacedim >::end_hex ( const unsigned int  level  )  const

Return an iterator which is the first iterator not on level. If level is the last level, then this returns end().

template<int dim, int spacedim = dim>
raw_hex_iterator Triangulation< dim, spacedim >::end_raw_hex ( const unsigned int  level  )  const

Return a raw iterator which is the first iterator not on level. If level is the last level, then this returns end().

template<int dim, int spacedim = dim>
active_hex_iterator Triangulation< dim, spacedim >::end_active_hex ( const unsigned int  level  )  const

Return an active iterator which is the first iterator not on level. If level is the last level, then this returns end().

template<int dim, int spacedim = dim>
raw_hex_iterator Triangulation< dim, spacedim >::last_raw_hex (  )  const

Return an iterator pointing to the last hex, used or not.

template<int dim, int spacedim = dim>
raw_hex_iterator Triangulation< dim, spacedim >::last_raw_hex ( const unsigned int  level  )  const

Return an iterator pointing to the last hex of the level level, used or not.

template<int dim, int spacedim = dim>
hex_iterator Triangulation< dim, spacedim >::last_hex (  )  const

Return an iterator pointing to the last used hex.

template<int dim, int spacedim = dim>
hex_iterator Triangulation< dim, spacedim >::last_hex ( const unsigned int  level  )  const

Return an iterator pointing to the last used hex on level level.

template<int dim, int spacedim = dim>
active_hex_iterator Triangulation< dim, spacedim >::last_active_hex (  )  const

Return an iterator pointing to the last active hex.

template<int dim, int spacedim = dim>
active_hex_iterator Triangulation< dim, spacedim >::last_active_hex ( const unsigned int  level  )  const

Return an iterator pointing to the last active hex on level level.

template<int dim, int spacedim = dim>
unsigned int Triangulation< dim, spacedim >::n_raw_lines (  )  const

In the following, most functions are provided in two versions, with and without an argument describing the level. The versions with this argument are only applicable for objects descibing the cells of the present triangulation. For example: in 2D n_lines(level) cannot be called, only n_lines(), as lines are faces in 2D and therefore have no level. Total Number of lines, used or unused.

template<int dim, int spacedim = dim>
unsigned int Triangulation< dim, spacedim >::n_raw_lines ( const unsigned int  level  )  const

Number of lines, used or unused, on the given level.

template<int dim, int spacedim = dim>
unsigned int Triangulation< dim, spacedim >::n_lines (  )  const

Return total number of used lines, active or not.

template<int dim, int spacedim = dim>
unsigned int Triangulation< dim, spacedim >::n_lines ( const unsigned int  level  )  const

Return total number of used lines, active or not on level level.

template<int dim, int spacedim = dim>
unsigned int Triangulation< dim, spacedim >::n_active_lines (  )  const

Return total number of active lines.

template<int dim, int spacedim = dim>
unsigned int Triangulation< dim, spacedim >::n_active_lines ( const unsigned int  level  )  const

Return total number of active lines, on level level.

template<int dim, int spacedim = dim>
unsigned int Triangulation< dim, spacedim >::n_raw_quads (  )  const

Total number of quads, used or unused.

template<int dim, int spacedim = dim>
unsigned int Triangulation< dim, spacedim >::n_raw_quads ( const unsigned int  level  )  const

Number of quads, used or unused, on the given level.

template<int dim, int spacedim = dim>
unsigned int Triangulation< dim, spacedim >::n_quads (  )  const

Return total number of used quads, active or not.

template<int dim, int spacedim = dim>
unsigned int Triangulation< dim, spacedim >::n_quads ( const unsigned int  level  )  const

Return total number of used quads, active or not on level level.

template<int dim, int spacedim = dim>
unsigned int Triangulation< dim, spacedim >::n_active_quads (  )  const

Return total number of active quads, active or not.

template<int dim, int spacedim = dim>
unsigned int Triangulation< dim, spacedim >::n_active_quads ( const unsigned int  level  )  const

Return total number of active quads, active or not on level level.

template<int dim, int spacedim = dim>
unsigned int Triangulation< dim, spacedim >::n_raw_hexs (  )  const

Total number of hexs, used or unused.

template<int dim, int spacedim = dim>
unsigned int Triangulation< dim, spacedim >::n_raw_hexs ( const unsigned int  level  )  const

Number of hexs, used or unused, on the given level.

template<int dim, int spacedim = dim>
unsigned int Triangulation< dim, spacedim >::n_hexs (  )  const

Return total number of used hexahedra, active or not.

template<int dim, int spacedim = dim>
unsigned int Triangulation< dim, spacedim >::n_hexs ( const unsigned int  level  )  const

Return total number of used hexahedra, active or not on level level.

template<int dim, int spacedim = dim>
unsigned int Triangulation< dim, spacedim >::n_active_hexs (  )  const

Return total number of active hexahedra, active or not.

template<int dim, int spacedim = dim>
unsigned int Triangulation< dim, spacedim >::n_active_hexs ( const unsigned int  level  )  const

Return total number of active hexahedra, active or not on level level.

template<int dim, int spacedim = dim>
unsigned int Triangulation< dim, spacedim >::n_raw_cells ( const unsigned int  level  )  const

Number of cells, used or unused, on the given level.

template<int dim, int spacedim = dim>
unsigned int Triangulation< dim, spacedim >::n_cells (  )  const

Return total number of used cells, active or not. Maps to n_lines() in one space dimension and so on.

template<int dim, int spacedim = dim>
unsigned int Triangulation< dim, spacedim >::n_cells ( const unsigned int  level  )  const

Return total number of used cells, active or not, on level level. Maps to n_lines(level) in one space dimension and so on.

template<int dim, int spacedim = dim>
unsigned int Triangulation< dim, spacedim >::n_active_cells (  )  const

Return total number of active cells. Maps to n_active_lines() in one space dimension and so on.

template<int dim, int spacedim = dim>
unsigned int Triangulation< dim, spacedim >::n_active_cells ( const unsigned int  level  )  const

Return total number of active cells on level level. Maps to n_active_lines(level) in one space dimension and so on.

template<int dim, int spacedim = dim>
unsigned int Triangulation< dim, spacedim >::n_faces (  )  const

Return total number of used faces, active or not. In 2D, the result equals n_lines(), while in 3D it equals n_quads(). Since there are no face objects in 1d, the function returns zero in 1d.

template<int dim, int spacedim = dim>
unsigned int Triangulation< dim, spacedim >::n_active_faces (  )  const

Return total number of active faces, active or not. In 2D, the result equals n_active_lines(), while in 3D it equals n_active_quads(). Since there are no face objects in 1d, the function returns zero in 1d.

template<int dim, int spacedim = dim>
unsigned int Triangulation< dim, spacedim >::n_levels (  )  const

Return number of levels in use. This may be less than the number of levels existing in the triangulation if by coarsening the highest level was completely depopulated. That level is not removed, since it will most likely be repopulated soon by the next refinement process.

Referenced by GridTools::transform().

template<int dim, int spacedim = dim>
unsigned int Triangulation< dim, spacedim >::n_vertices (  )  const

Return the total number of vertices. Some of them may not be used, which usually happens upon coarsening of a triangulation when some vertices are discarded, but we do not want to renumber the remaining ones, leading to holes in the numbers of used vertices. You can get the number of used vertices using n_used_vertices function.

Referenced by GridTools::transform().

template<int dim, int spacedim = dim>
const std::vector<Point<spacedim> >& Triangulation< dim, spacedim >::get_vertices (  )  const

Return a constant reference to all the vertices present in this triangulation. Note that not necessarily all vertices in this array are actually used; for example, if you coarsen a mesh, then some vertices are deleted, but their positions in this array are unchanged as the indices of vertices are only allocated once. You can find out about which vertices are actually used by the function get_used_vertices().

template<int dim, int spacedim = dim>
unsigned int Triangulation< dim, spacedim >::n_used_vertices (  )  const

Return the number of vertices that are presently in use, i.e. belong to at least one used element.

template<int dim, int spacedim = dim>
bool Triangulation< dim, spacedim >::vertex_used ( const unsigned int  index  )  const

Return true if the vertex with this index is used.

template<int dim, int spacedim = dim>
const std::vector<bool>& Triangulation< dim, spacedim >::get_used_vertices (  )  const

Return a constant reference to the array of bools indicating whether an entry in the vertex array is used or not.

template<int dim, int spacedim = dim>
unsigned int Triangulation< dim, spacedim >::max_adjacent_cells (  )  const

Return the maximum number of cells meeting at a common vertex. Since this number is an invariant under refinement, only the cells on the coarsest level are considered. The operation is thus reasonably fast. The invariance is only true for sufficiently many cells in the coarsest triangulation (e.g. for a single cell one would be returned), so a minimum of four is returned in two dimensions, 8 in three dimensions, etc, which is how many cells meet if the triangulation is refined.

In one space dimension, two is returned.

template<int dim, int spacedim = dim>
virtual unsigned int Triangulation< dim, spacedim >::memory_consumption (  )  const [virtual]

Determine an estimate for the memory consumption (in bytes) of this object.

This function is made virtual, since a triangulation object might be accessed through a pointer to this base class, even if the actual object is a derived class.

Reimplemented in PersistentTriangulation< dim >.

template<int dim, int spacedim = dim>
static void Triangulation< dim, spacedim >::write_bool_vector ( const unsigned int  magic_number1,
const std::vector< bool > &  v,
const unsigned int  magic_number2,
std::ostream &  out 
) [static, protected]

Write a bool vector to the given stream, writing a pre- and a postfix magic number. The vector is written in an almost binary format, i.e. the bool flags are packed but the data is written as ASCII text.

The flags are stored in a binary format: for each true, a 1 bit is stored, a 0 bit otherwise. The bits are stored as unsigned char, thus avoiding endianess. They are written to out in plain text, thus amounting to 3.6 bits in the output per bits in the input on the average. Other information (magic numbers and number of elements of the input vector) is stored as plain text as well. The format should therefore be interplatform compatible.

template<int dim, int spacedim = dim>
static void Triangulation< dim, spacedim >::read_bool_vector ( const unsigned int  magic_number1,
std::vector< bool > &  v,
const unsigned int  magic_number2,
std::istream &  in 
) [static, protected]

Re-read a vector of bools previously written by write_bool_vector and compare with the magic numbers.

template<int dim, int spacedim = dim>
void Triangulation< dim, spacedim >::clear_despite_subscriptions (  )  [private]

The (public) clear() will only work when the triangulation is not subscriped to by other users. The clear_despite_subscriptions() function now allows the triangulation being cleared even when there are subscriptions.

Make sure, you know what you do, when calling this function, as its use is reasonable in very rare cases, only. For example, when the subscriptions were for the initially empty Triangulation and the Triangulation object wants to release its memory before throwing an assertion due to input errors (e.g. in the create_triangulation() function).

template<int dim, int spacedim = dim>
void Triangulation< dim, spacedim >::execute_refinement (  )  [private]

Refine all cells on all levels which were previously flagged for refinement.

Note, that this function uses the line->user_flags for dim=2,3 and the quad->user_flags for dim=3.

template<int dim, int spacedim = dim>
void Triangulation< dim, spacedim >::execute_coarsening (  )  [private]

Coarsen all cells which were flagged for coarsening, or rather: delete all children of those cells of which all child cells are flagged for coarsening and several other constraints hold (see the general doc of this class).

template<int dim, int spacedim = dim>
void Triangulation< dim, spacedim >::fix_coarsen_flags (  )  [private]

Make sure that either all or none of the children of a cell are tagged for coarsening.


Friends And Related Function Documentation

template<int dim, int spacedim = dim>
friend class TriaAccessorBase [friend]
template<int dim, int spacedim = dim>
friend class TriaAccessor [friend]
template<int dim, int spacedim = dim>
friend class CellAccessor< dim, spacedim > [friend]
template<int dim, int spacedim = dim>
friend class internal::TriaAccessor::Implementation [friend]
template<int dim, int spacedim = dim>
friend class hp::DoFHandler< dim, spacedim > [friend]
template<int dim, int spacedim = dim>
friend class internal::Triangulation::Implementation [friend]

Member Data Documentation

template<int dim, int spacedim = dim>
const StraightBoundary<dim,spacedim> Triangulation< dim, spacedim >::straight_boundary [static]

Default boundary object. This is used for those boundaries for which no boundary object has been explicitly set using set_boundary().

template<int dim, int spacedim = dim>
const StraightBoundary<spacedim,spacedim> Triangulation< dim, spacedim >::straight_manifold_description [static]

Default boundary object to be used in the codimension 1 case. This is used for those boundaries for which no boundary object has been explicitly set using set_boundary().

template<int dim, int spacedim = dim>
const unsigned int Triangulation< dim, spacedim >::dimension = dim [static]

Make the dimension available in function templates.

Reimplemented in PersistentTriangulation< dim >.

template<int dim, int spacedim = dim>
const unsigned int Triangulation< dim, spacedim >::space_dimension = spacedim [static]

Make the space-dimension available in function templates.

template<int dim, int spacedim = dim>
std::vector<internal::Triangulation::TriaLevel<dim>*> Triangulation< dim, spacedim >::levels [private]

Array of pointers pointing to the objects storing the cell data on the different levels.

template<int dim, int spacedim = dim>
internal::Triangulation::TriaFaces<dim>* Triangulation< dim, spacedim >::faces [private]

Pointer to the faces of the triangulation. In 1d this contains nothing, in 2D it contains data concerning lines and in 3D quads and lines. All of these have no level and are therefore treated seperately.

template<int dim, int spacedim = dim>
std::vector<Point<spacedim> > Triangulation< dim, spacedim >::vertices [private]

Array of the vertices of this triangulation.

template<int dim, int spacedim = dim>
std::vector<bool> Triangulation< dim, spacedim >::vertices_used [private]

Array storing a bit-pattern which vertices are used.

template<int dim, int spacedim = dim>
SmartPointer<const Boundary<dim,spacedim> > Triangulation< dim, spacedim >::boundary[255] [private]

Collection of boundary objects. We only need 255 objects rather than 256, since the indicator 255 is reserved for interior faces and can thus never be associated with a boundary.

template<int dim, int spacedim = dim>
SmartPointer<const Boundary<spacedim,spacedim> > Triangulation< dim, spacedim >::manifold_description[255] [private]

Collection of ly objects. We only need 255 objects rather than 256, since the indicator 255 is reserved for interior faces and can thus never be associated with a boundary.

template<int dim, int spacedim = dim>
bool Triangulation< dim, spacedim >::anisotropic_refinement [private]

Flag indicating whether anisotropic refinement took place.

template<int dim, int spacedim = dim>
MeshSmoothing Triangulation< dim, spacedim >::smooth_grid [private]

Do some smoothing in the process of refining the triangulation. See the general doc of this class for more information about this.

template<int dim, int spacedim = dim>
internal::Triangulation::NumberCache<dim> Triangulation< dim, spacedim >::number_cache [private]

Cache to hold the numbers of lines, quads, hexes, etc. These numbers are set at the end of the refinement and coarsening functions and enable faster access later on. In the old days, whenever one wanted to access one of these numbers, one had to perform a loop over all lines, e.g., and count the elements until we hit the end iterator. This is time consuming and since access to the number of lines etc is a rather frequent operation, this was not an optimal solution.

template<int dim, int spacedim = dim>
std::list<RefinementListener *> Triangulation< dim, spacedim >::refinement_listeners [mutable, private]

List of RefinementListeners, which want to be informed if the Triangulation is refined.


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

deal.II documentation generated on Mon Nov 23 22:58:15 2009 by doxygen 1.6.1