triangulation.vertices_used.size(),
ExcInternalError());
- {
- typename Triangulation<dim, spacedim>::raw_cell_iterator
- cell = triangulation.begin_active(triangulation.levels.size() - 1),
- endc = triangulation.end();
- for (; cell != endc; ++cell)
- if (cell->used())
- if (cell->refine_flag_set())
- {
- triangulation.levels.push_back(
- std::make_unique<
- internal::TriangulationImplementation::TriaLevel>(dim));
- break;
- }
- }
-
+ // Check whether a new level is needed. We have to check for
+ // this on the highest level only
+ for (const auto &cell : triangulation.active_cell_iterators_on_level(
+ triangulation.levels.size() - 1))
+ if (cell->refine_flag_set())
+ {
+ triangulation.levels.push_back(
+ std::make_unique<
+ internal::TriangulationImplementation::TriaLevel>(dim));
+ break;
+ }
+ // first clear user flags for quads and lines; we're going to
+ // use them to flag which lines and quads need refinement
triangulation.faces->quads.clear_user_data();
for (typename Triangulation<dim, spacedim>::line_iterator line =
++quad)
quad->clear_user_flag();
+ // check how much space is needed on every level we need not
+ // check the highest level since either
+ // - on the highest level no cells are flagged for refinement
+ // - there are, but prepare_refinement added another empty
+ // level which then is the highest level
+
+ // variables to hold the number of newly to be created
+ // vertices, lines and quads. as these are stored globally,
+ // declare them outside the loop over al levels. we need lines
+ // and quads in pairs for refinement of old ones and lines and
+ // quads, that can be stored as single ones, as they are newly
+ // created in the inside of an existing cell
unsigned int needed_vertices = 0;
unsigned int needed_lines_single = 0;
unsigned int needed_quads_single = 0;
{
unsigned int new_cells = 0;
- for (const auto &acell :
+ for (const auto &cell :
triangulation.active_cell_iterators_on_level(level))
- if (acell->refine_flag_set())
+ if (cell->refine_flag_set())
{
- Assert(acell->refine_flag_set() ==
+ Assert(cell->refine_flag_set() ==
RefinementCase<dim>::cut_xyz,
ExcInternalError());
- if (acell->reference_cell() == ReferenceCells::Hexahedron)
+ if (cell->reference_cell() == ReferenceCells::Hexahedron)
{
++needed_vertices;
needed_lines_single += 6;
needed_quads_single += 12;
new_cells += 8;
}
- else if (acell->reference_cell() ==
+ else if (cell->reference_cell() ==
ReferenceCells::Tetrahedron)
{
needed_lines_single += 1;
Assert(false, ExcInternalError());
}
- for (const auto face : acell->face_indices())
- if (acell->face(face)->number_of_children() < 4)
- acell->face(face)->set_user_flag();
+ for (const auto face : cell->face_indices())
+ if (cell->face(face)->number_of_children() < 4)
+ cell->face(face)->set_user_flag();
- for (const auto line : acell->line_indices())
- if (acell->line(line)->has_children() == false)
- acell->line(line)->set_user_flag();
+ for (const auto line : cell->line_indices())
+ if (cell->line(line)->has_children() == false)
+ cell->line(line)->set_user_flag();
}
const unsigned int used_cells =
reserve_space(triangulation.levels[level + 1]->cells, new_cells);
}
+ // now count the quads and lines which were flagged for
+ // refinement
for (typename Triangulation<dim, spacedim>::quad_iterator quad =
triangulation.begin_quad();
quad != triangulation.end_quad();
needed_quads_single);
+ // add to needed vertices how many vertices are already in use
needed_vertices += std::count(triangulation.vertices_used.begin(),
triangulation.vertices_used.end(),
true);
triangulation.vertices_used.resize(needed_vertices, false);
}
+ ///////////////////////////////////////////
+ // Before we start with the actual refinement, we do some
+ // sanity checks if in debug mode. especially, we try to catch
+ // the notorious problem with lines being twice refined,
+ // i.e. there are cells adjacent at one line ("around the
+ // edge", but not at a face), with two cells differing by more
+ // than one refinement level
+ //
+ // this check is very simple to implement here, since we have
+ // all lines flagged if they shall be refined
+#ifdef DEBUG
+ for (const auto &cell : triangulation.active_cell_iterators())
+ if (!cell->refine_flag_set())
+ for (unsigned int line_n = 0; line_n < cell->n_lines(); ++line_n)
+ if (cell->line(line_n)->has_children())
+ for (unsigned int c = 0; c < 2; ++c)
+ Assert(cell->line(line_n)->child(c)->user_flag_set() == false,
+ ExcInternalError());
+#endif
+
unsigned int current_vertex = 0;
// helper function - find the next available vertex number and mark it
Assert(next_unused_line.state() == IteratorState::valid,
ExcInternalError());
+ // now we found two consecutive unused lines, such
+ // that the children of a line will be consecutive.
+ // then set the child pointer of the present line
line->set_children(0, next_unused_line->index());
const typename Triangulation<dim, spacedim>::raw_line_iterator
}
}
- // 3) create new quads (properties are set below)
+ // 3) create new quads (properties are set below). Both triangles
+ // and quads are divided in four.
std::array<
typename Triangulation<dim, spacedim>::raw_quad_iterator,
4>
quad->set_refinement_case(RefinementCase<2>::cut_xy);
}
- std::array<unsigned int, 9> vertex_indices;
+ // Maximum of 9 vertices per refined quad (9 for Quadrilateral, 6
+ // for Triangle)
+ std::array<unsigned int, 9> vertex_indices = {};
{
- for (auto &i : vertex_indices)
- i = 0.0;
-
unsigned int k = 0;
for (const auto i : quad->vertex_indices())
vertex_indices[k++] = quad->vertex_index(i);
{
static constexpr std::array<std::array<unsigned int, 2>,
2>
- index = {
- {{{1, 0}}, // child 0, line_orientation=false and true
- {{0,
- 1}}}}; // child 1, line_orientation=false and true
+ index = {// child 0, line_orientation=false and true
+ {{{1, 0}},
+ // child 1, line_orientation=false and true
+ {{0, 1}}}};
lines[k++] = quad->line(l)->child(
index[c][quad->line_orientation(l)]);
hex->center(true, true);
}
- std::vector<
- typename Triangulation<dim, spacedim>::raw_line_iterator>
+ boost::container::small_vector<
+ typename Triangulation<dim, spacedim>::raw_line_iterator,
+ 6>
new_lines(n_new_lines);
for (unsigned int i = 0; i < n_new_lines; ++i)
{
new_lines[i]->set_manifold_id(hex->manifold_id());
}
- std::vector<
- typename Triangulation<dim, spacedim>::raw_quad_iterator>
+ boost::container::small_vector<
+ typename Triangulation<dim, spacedim>::raw_quad_iterator,
+ 12>
new_quads(n_new_quads);
for (unsigned int i = 0; i < n_new_quads; ++i)
{
new_quad->set_line_orientation(j, true);
}
- std::vector<
- typename Triangulation<dim, spacedim>::raw_hex_iterator>
- new_hexes(n_new_hexes);
+ // we always get 8 children per refined cell
+ std::array<
+ typename Triangulation<dim, spacedim>::raw_hex_iterator,
+ 8>
+ new_hexes;
{
for (unsigned int i = 0; i < n_new_hexes; ++i)
{
if (i % 2)
new_hex->set_parent(hex->index());
+ // set the face_orientation flag to true for all
+ // faces initially, as this is the default value
+ // which is true for all faces interior to the
+ // hex. later on go the other way round and
+ // reset faces that are at the boundary of the
+ // mother cube
+ //
+ // the same is true for the face_flip and
+ // face_rotation flags. however, the latter two
+ // are set to false by default as this is the
+ // standard value
for (const auto f : new_hex->face_indices())
{
new_hex->set_face_orientation(f, true);
new_hex->set_face_rotation(f, false);
}
}
- // note these hexes as children to the present cell
for (unsigned int i = 0; i < n_new_hexes / 2; ++i)
hex->set_children(2 * i, new_hexes[2 * i]->index());
}
{
// load vertex indices
- std::array<unsigned int, 27> vertex_indices;
- for (auto &i : vertex_indices)
- i = 0;
+ std::array<unsigned int, 27> vertex_indices = {};
{
unsigned int k = 0;
projects / dealii.git / commitdiff