* @ingroup Exceptions
*/
DeclException0 (ExcCellHasNoChildren);
+ /**
+ * Trying to access the parent of
+ * a cell which is in the coarsest
+ * level of the triangulation.
+ *
+ * @ingroup Exceptions
+ */
+ DeclException0 (ExcCellHasNoParent);
//TODO: Write documentation!
/**
* @ingroup Exceptions
* in the library.
*/
void clear_refinement_case () const;
+
+ /**
+ * Set the parent of a cell.
+ */
+ void set_parent (const unsigned int parent_index);
/**
* Set the index of the ith
* cell.
*/
void set_subdomain_id (const unsigned int new_subdomain_id) const;
-
+
+ /**
+ * Return an iterator to the
+ * parent.
+ */
+ TriaIterator<CellAccessor<dim,spacedim> >
+ parent () const;
+
/**
* @}
*/
+template <int structdim, int dim, int spacedim>
+void
+TriaAccessor<structdim, dim, spacedim>::set_parent (const unsigned int parent_index)
+{
+ Assert (this->used(), TriaAccessorExceptions::ExcCellNotUsed());
+ Assert (this->present_level > 0, TriaAccessorExceptions::ExcCellHasNoParent ());
+ this->tria->levels[this->present_level]->parents[this->present_index / 2]
+ = parent_index;
+}
+
+
+
template <int structdim, int dim, int spacedim>
void
TriaAccessor<structdim, dim, spacedim>::clear_children () const
* number.
*/
std::vector<unsigned int> subdomain_ids;
+
+ /**
+ * One integer for every consecutive
+ * pair of cells to store which
+ * index their parent has.
+ */
+ std::vector<int> parents;
/**
* Reserve enough space to accomodate
std::vector<unsigned char> refine_flags;
std::vector<bool> coarsen_flags;
std::vector<std::pair<int,int> > neighbors;
+ std::vector<int> parents;
std::vector<unsigned int> subdomain_ids;
void reserve_space (const unsigned int total_cells,
const unsigned int dimension);
// properties
subcells[i]->set_material_id (cell->material_id());
subcells[i]->set_subdomain_id (cell->subdomain_id());
+
+ if (i%2==0)
+ subcells[i]->set_parent (cell->index ());
}
// set child index for
next_unused_vertex));
first_child->set_material_id (cell->material_id());
first_child->set_subdomain_id (cell->subdomain_id());
+ first_child->set_parent (cell->index ());
// reset neighborship info (refer
// to
{
if (i%2==0)
next_unused_hex=triangulation.levels[level+1]->cells.next_free_hex(triangulation,level+1);
+
else
++next_unused_hex;
// properties
new_hexes[i]->set_material_id (hex->material_id());
new_hexes[i]->set_subdomain_id (hex->subdomain_id());
+
+ if (i%2)
+ new_hexes[i]->set_parent (hex->index ());
// set the face_orientation
// flag to true for all faces
// initially, as this is the
}
+template <int dim, int spacedim>
+TriaIterator<CellAccessor<dim,spacedim> >
+CellAccessor<dim, spacedim>::parent () const
+{
+ Assert (this->used(), TriaAccessorExceptions::ExcCellNotUsed());
+ Assert (this->present_level > 0, TriaAccessorExceptions::ExcCellHasNoParent ());
+ TriaIterator<CellAccessor<dim,spacedim> >
+ q (this->tria, this->present_level-1,
+ this->tria->levels[this->present_level]->parents[this->present_index / 2]);
+
+ return q;
+}
+
template <int dim, int spacedim>
void CellAccessor<dim, spacedim>::set_neighbor (const unsigned int i,
subdomain_ids.insert (subdomain_ids.end(),
total_cells - subdomain_ids.size(),
0);
+
+ parents.reserve ((int) (total_cells + 1) / 2);
+ parents.insert (parents.end (),
+ (total_cells + 1) / 2 - parents.size (),
+ -1);
neighbors.reserve (total_cells*(2*dimension));
neighbors.insert (neighbors.end(),
subdomain_ids.insert (subdomain_ids.end(),
total_cells - subdomain_ids.size(),
0);
+
+ parents.reserve ((int) (total_cells + 1) / 2);
+ parents.insert (parents.end (),
+ (total_cells + 1) / 2 - parents.size (),
+ -1);
neighbors.reserve (total_cells*(2*dimension));
neighbors.insert (neighbors.end(),
<ol>
+ <li>
+ <p>
+ New: Data structures and a function parent to get the parent of a cell.
+ <br>
+ (Markus Buerg 2010/08/26)
+ </p>
+
<li>
<p>
Improved: DoFHandler has a default constructor, so that it can be used in containers.