17 #include <deal.II/base/geometry_info.h> 18 #include <deal.II/base/memory_consumption.h> 19 #include <deal.II/base/std_cxx14/memory.h> 21 #include <deal.II/fe/mapping_q1.h> 23 #include <deal.II/grid/grid_tools.h> 24 #include <deal.II/grid/magic_numbers.h> 25 #include <deal.II/grid/manifold.h> 26 #include <deal.II/grid/tria.h> 27 #include <deal.II/grid/tria_accessor.h> 28 #include <deal.II/grid/tria_faces.h> 29 #include <deal.II/grid/tria_iterator.h> 30 #include <deal.II/grid/tria_levels.h> 32 #include <deal.II/lac/full_matrix.h> 33 #include <deal.II/lac/vector.h> 44 DEAL_II_NAMESPACE_OPEN
62 namespace TriangulationImplementation
136 template <
int dim,
int spacedim>
138 cell_is_patch_level_1(
143 unsigned int n_active_children = 0;
144 for (
unsigned int i = 0; i < cell->n_children(); ++i)
145 if (cell->child(i)->active())
148 return (n_active_children == 0) ||
149 (n_active_children == cell->n_children());
159 template <
int dim,
int spacedim>
161 cell_will_be_coarsened(
167 if (cell->has_children())
169 unsigned int children_to_coarsen = 0;
170 const unsigned int n_children = cell->n_children();
172 for (
unsigned int c = 0; c < n_children; ++c)
173 if (cell->child(c)->active() && cell->child(c)->coarsen_flag_set())
174 ++children_to_coarsen;
175 if (children_to_coarsen == n_children)
178 for (
unsigned int c = 0; c < n_children; ++c)
179 if (cell->child(c)->active())
180 cell->child(c)->clear_coarsen_flag();
216 template <
int dim,
int spacedim>
218 face_will_be_refined_by_neighbor_internal(
220 const unsigned int face_no,
229 cell->neighbor(face_no);
236 if (neighbor->has_children())
241 if (cell_will_be_coarsened(neighbor))
250 expected_face_ref_case = cell->face(face_no)->refinement_case();
262 const unsigned int neighbor_neighbor = cell->neighbor_face_no(face_no);
270 neighbor->face_orientation(neighbor_neighbor),
271 neighbor->face_flip(neighbor_neighbor),
272 neighbor->face_rotation(neighbor_neighbor));
276 neighbor_face = neighbor->face(neighbor_neighbor);
277 const int this_face_index = cell->face_index(face_no);
283 if (neighbor_face->index() == this_face_index)
289 expected_face_ref_case = face_ref_case;
307 for (
unsigned int c = 0; c < neighbor_face->n_children(); ++c)
308 if (neighbor_face->child_index(c) == this_face_index)
317 if ((neighbor_face->refinement_case() | face_ref_case) ==
318 neighbor_face->refinement_case())
337 expected_face_ref_case =
338 ~neighbor_face->refinement_case();
365 template <
int dim,
int spacedim>
367 face_will_be_refined_by_neighbor(
369 const unsigned int face_no)
372 return face_will_be_refined_by_neighbor_internal(cell, face_no, dummy);
379 template <
int dim,
int spacedim>
381 face_will_be_refined_by_neighbor(
383 const unsigned int face_no,
386 return face_will_be_refined_by_neighbor_internal(cell,
388 expected_face_ref_case);
393 template <
int dim,
int spacedim>
395 satisfies_level1_at_vertex_rule(
398 std::vector<unsigned int> min_adjacent_cell_level(
399 triangulation.n_vertices(), triangulation.n_levels());
400 std::vector<unsigned int> max_adjacent_cell_level(
401 triangulation.n_vertices(), 0);
404 triangulation.begin_active();
405 cell != triangulation.end();
407 for (
unsigned int v = 0; v < GeometryInfo<dim>::vertices_per_cell; ++v)
409 min_adjacent_cell_level[cell->vertex_index(v)] =
410 std::min<unsigned int>(
411 min_adjacent_cell_level[cell->vertex_index(v)], cell->level());
412 max_adjacent_cell_level[cell->vertex_index(v)] =
413 std::max<unsigned int>(
414 min_adjacent_cell_level[cell->vertex_index(v)], cell->level());
417 for (
unsigned int k = 0; k < triangulation.n_vertices(); ++k)
418 if (triangulation.vertex_used(k))
419 if (max_adjacent_cell_level[k] - min_adjacent_cell_level[k] > 1)
432 template <
int dim,
int spacedim>
433 std::vector<unsigned int>
438 std::vector<unsigned int> line_cell_count(triangulation.n_raw_lines(),
441 cell = triangulation.begin(),
442 endc = triangulation.end();
443 for (; cell != endc; ++cell)
444 for (
unsigned int l = 0; l < GeometryInfo<dim>::lines_per_cell; ++
l)
445 ++line_cell_count[cell->line_index(l)];
446 return line_cell_count;
449 return std::vector<unsigned int>();
460 template <
int dim,
int spacedim>
461 std::vector<unsigned int>
466 std::vector<unsigned int> quad_cell_count(triangulation.n_raw_quads(),
469 cell = triangulation.begin(),
470 endc = triangulation.end();
471 for (; cell != endc; ++cell)
472 for (
unsigned int q = 0; q < GeometryInfo<dim>::faces_per_cell; ++q)
473 ++quad_cell_count[cell->quad_index(q)];
474 return quad_cell_count;
477 return std::vector<unsigned int>();
501 void reorder_compatibility(std::vector<
CellData<2>> &cells,
504 for (
auto &cell : cells)
505 std::swap(cell.vertices[2], cell.vertices[3]);
509 void reorder_compatibility(std::vector<
CellData<3>> &cells,
513 for (
auto &cell : cells)
515 for (
unsigned int i = 0; i < GeometryInfo<3>::vertices_per_cell; ++i)
516 tmp[i] = cell.vertices[i];
517 for (
unsigned int i = 0; i < GeometryInfo<3>::vertices_per_cell; ++i)
522 std::vector<CellData<2>>::iterator boundary_quad =
524 std::vector<CellData<2>>::iterator end_quad =
526 for (
unsigned int quad_no = 0; boundary_quad != end_quad;
527 ++boundary_quad, ++quad_no)
528 std::swap(boundary_quad->vertices[2], boundary_quad->vertices[3]);
550 template <
int dim,
int spacedim>
555 if (line->has_children())
556 return line->child(0)->vertex_index(1);
561 template <
int dim,
int spacedim>
566 switch (static_cast<unsigned char>(quad->refinement_case()))
569 return middle_vertex_index<dim, spacedim>(quad->child(0)->line(1));
572 return middle_vertex_index<dim, spacedim>(quad->child(0)->line(3));
575 return quad->child(0)->vertex_index(3);
584 template <
int dim,
int spacedim>
589 switch (static_cast<unsigned char>(hex->refinement_case()))
592 return middle_vertex_index<dim, spacedim>(hex->child(0)->quad(1));
595 return middle_vertex_index<dim, spacedim>(hex->child(0)->quad(3));
598 return middle_vertex_index<dim, spacedim>(hex->child(0)->quad(5));
601 return middle_vertex_index<dim, spacedim>(hex->child(0)->line(11));
604 return middle_vertex_index<dim, spacedim>(hex->child(0)->line(5));
607 return middle_vertex_index<dim, spacedim>(hex->child(0)->line(7));
610 return hex->child(0)->vertex_index(7);
631 template <
class TRIANGULATION>
632 inline typename TRIANGULATION::DistortedCellList
633 collect_distorted_coarse_cells(
const TRIANGULATION &)
635 return typename TRIANGULATION::DistortedCellList();
654 triangulation.begin(0);
655 cell != triangulation.end(0);
659 for (
unsigned int i = 0; i < GeometryInfo<dim>::vertices_per_cell; ++i)
660 vertices[i] = cell->vertex(i);
665 for (
unsigned int i = 0; i < GeometryInfo<dim>::vertices_per_cell; ++i)
666 if (determinants[i] <= 1e-9 * std::pow(cell->diameter(), 1. * dim))
668 distorted_cells.distorted_cells.push_back(cell);
673 return distorted_cells;
685 has_distorted_children(
687 std::integral_constant<int, dim>,
688 std::integral_constant<int, dim>)
692 for (
unsigned int c = 0; c < cell->n_children(); ++c)
695 for (
unsigned int i = 0; i < GeometryInfo<dim>::vertices_per_cell; ++i)
696 vertices[i] = cell->child(c)->vertex(i);
701 for (
unsigned int i = 0; i < GeometryInfo<dim>::vertices_per_cell; ++i)
702 if (determinants[i] <=
703 1e-9 * std::pow(cell->child(c)->diameter(), 1. * dim))
718 template <
int dim,
int spacedim>
720 has_distorted_children(
722 std::integral_constant<int, dim>,
723 std::integral_constant<int, spacedim>)
734 template <
int spacedim>
739 template <
int dim,
int spacedim>
794 static const unsigned int left_right_offset[2][6][2] = {
803 {{0, 1}, {1, 0}, {0, 1}, {1, 0}, {0, 1}, {1, 0}}};
814 std::vector<typename Triangulation<dim, spacedim>::cell_iterator>
815 adjacent_cells(2 * triangulation.n_raw_faces(), dummy);
821 for (; cell != endc; ++cell)
822 for (
unsigned int f = 0; f < GeometryInfo<dim>::faces_per_cell; ++f)
827 const unsigned int offset =
828 (cell->direction_flag() ?
829 left_right_offset[dim - 2][f][cell->face_orientation(f)] :
830 1 - left_right_offset[dim - 2][f][cell->face_orientation(f)]);
832 adjacent_cells[2 * face->index() + offset] = cell;
842 if (cell->active() && face->has_children())
844 adjacent_cells[2 * face->child(0)->index() + offset] = cell;
845 adjacent_cells[2 * face->child(1)->index() + offset] = cell;
870 if (face->has_children() &&
873 cell->refinement_case(), f) ==
876 for (
unsigned int c = 0; c < face->n_children(); ++c)
877 adjacent_cells[2 * face->child(c)->index() + offset] = cell;
878 if (face->child(0)->has_children())
880 adjacent_cells[2 * face->child(0)->child(0)->index() +
882 adjacent_cells[2 * face->child(0)->child(1)->index() +
885 if (face->child(1)->has_children())
887 adjacent_cells[2 * face->child(1)->child(0)->index() +
889 adjacent_cells[2 * face->child(1)->child(1)->index() +
901 for (cell = triangulation.begin(); cell != endc; ++cell)
902 for (
unsigned int f = 0; f < GeometryInfo<dim>::faces_per_cell; ++f)
904 const unsigned int offset =
905 (cell->direction_flag() ?
906 left_right_offset[dim - 2][f][cell->face_orientation(f)] :
907 1 - left_right_offset[dim - 2][f][cell->face_orientation(f)]);
909 f, adjacent_cells[2 * cell->face(f)->index() + 1 - offset]);
914 template <
int dim,
int spacedim>
916 update_periodic_face_map_recursively(
919 unsigned int n_face_1,
920 unsigned int n_face_2,
921 const std::bitset<3> & orientation,
927 std::bitset<3>>> &periodic_face_map)
930 const FaceIterator face_1 = cell_1->face(n_face_1);
931 const FaceIterator face_2 = cell_2->face(n_face_2);
933 const bool face_orientation = orientation[0];
934 const bool face_flip = orientation[1];
935 const bool face_rotation = orientation[2];
937 Assert((dim != 1) || (face_orientation ==
true && face_flip ==
false &&
938 face_rotation ==
false),
940 "(face_orientation, face_flip, face_rotation) " 941 "is invalid for 1D"));
943 Assert((dim != 2) || (face_orientation ==
true && face_rotation ==
false),
945 "(face_orientation, face_flip, face_rotation) " 946 "is invalid for 2D"));
950 Assert(face_1->at_boundary() && face_2->at_boundary(),
951 ExcMessage(
"Periodic faces must be on the boundary"));
957 std::pair<typename Triangulation<dim, spacedim>::cell_iterator,
959 const CellFace cell_face_1(cell_1, n_face_1);
960 const CellFace cell_face_2(cell_2, n_face_2);
961 const std::pair<CellFace, std::bitset<3>> cell_face_orientation_2(
962 cell_face_2, orientation);
964 const std::pair<CellFace, std::pair<CellFace, std::bitset<3>>>
965 periodic_faces(cell_face_1, cell_face_orientation_2);
969 periodic_face_map.insert(periodic_faces);
975 static const int lookup_table_2d[2][2] =
982 static const int lookup_table_3d[2][2][2][4] =
1001 if (cell_1->has_children())
1003 if (cell_2->has_children())
1009 Assert(face_1->n_children() ==
1015 for (
unsigned int i = 0;
1016 i < GeometryInfo<dim>::max_children_per_face;
1024 j = lookup_table_2d[face_flip][i];
1027 j = lookup_table_3d[face_orientation][face_flip]
1035 unsigned int child_cell_1 =
1037 cell_1->refinement_case(),
1040 cell_1->face_orientation(n_face_1),
1041 cell_1->face_flip(n_face_1),
1042 cell_1->face_rotation(n_face_1),
1043 face_1->refinement_case());
1044 unsigned int child_cell_2 =
1046 cell_2->refinement_case(),
1049 cell_2->face_orientation(n_face_2),
1050 cell_2->face_flip(n_face_2),
1051 cell_2->face_rotation(n_face_2),
1052 face_2->refinement_case());
1054 Assert(cell_1->child(child_cell_1)->face(n_face_1) ==
1057 Assert(cell_2->child(child_cell_2)->face(n_face_2) ==
1063 update_periodic_face_map_recursively<dim, spacedim>(
1064 cell_1->child(child_cell_1),
1065 cell_2->child(child_cell_2),
1074 for (
unsigned int i = 0;
1075 i < GeometryInfo<dim>::max_children_per_face;
1079 unsigned int child_cell_1 =
1081 cell_1->refinement_case(),
1084 cell_1->face_orientation(n_face_1),
1085 cell_1->face_flip(n_face_1),
1086 cell_1->face_rotation(n_face_1),
1087 face_1->refinement_case());
1090 update_periodic_face_map_recursively<dim, spacedim>(
1091 cell_1->child(child_cell_1),
1108 namespace TriangulationImplementation
1116 using ::Triangulation;
1124 <<
"Something went wrong when making cell " << arg1
1125 <<
". Read the docs and the source code " 1126 <<
"for more information.");
1133 <<
"Something went wrong upon construction of cell " 1147 <<
" has negative measure. This typically " 1148 <<
"indicates some distortion in the cell, or a mistakenly " 1149 <<
"swapped pair of vertices in the input to " 1150 <<
"Triangulation::create_triangulation().");
1162 <<
"Error while creating cell " << arg1
1163 <<
": the vertex index " << arg2 <<
" must be between 0 and " 1172 <<
"While trying to assign a boundary indicator to a line: " 1173 <<
"the line with end vertices " << arg1 <<
" and " << arg2
1174 <<
" does not exist.");
1184 <<
"While trying to assign a boundary indicator to a quad: " 1185 <<
"the quad with bounding lines " << arg1 <<
", " << arg2
1186 <<
", " << arg3 <<
", " << arg4 <<
" does not exist.");
1196 <<
"The input data for creating a triangulation contained " 1197 <<
"information about a line with indices " << arg1 <<
" and " << arg2
1198 <<
" that is described to have boundary indicator " 1199 << static_cast<int>(arg3)
1200 <<
". However, this is an internal line not located on the " 1201 <<
"boundary. You cannot assign a boundary indicator to it." << std::endl
1203 <<
"If this happened at a place where you call " 1204 <<
"Triangulation::create_triangulation() yourself, you need " 1205 <<
"to check the SubCellData object you pass to this function." 1208 <<
"If this happened in a place where you are reading a mesh " 1209 <<
"from a file, then you need to investigate why such a line " 1210 <<
"ended up in the input file. A typical case is a geometry " 1211 <<
"that consisted of multiple parts and for which the mesh " 1212 <<
"generator program assumes that the interface between " 1213 <<
"two parts is a boundary when that isn't supposed to be " 1214 <<
"the case, or where the mesh generator simply assigns " 1215 <<
"'geometry indicators' to lines at the perimeter of " 1216 <<
"a part that are not supposed to be interpreted as " 1217 <<
"'boundary indicators'.");
1229 <<
"The input data for creating a triangulation contained " 1230 <<
"information about a quad with indices " << arg1 <<
", " << arg2
1231 <<
", " << arg3 <<
", and " << arg4
1232 <<
" that is described to have boundary indicator " 1233 << static_cast<int>(arg5)
1234 <<
". However, this is an internal quad not located on the " 1235 <<
"boundary. You cannot assign a boundary indicator to it." << std::endl
1237 <<
"If this happened at a place where you call " 1238 <<
"Triangulation::create_triangulation() yourself, you need " 1239 <<
"to check the SubCellData object you pass to this function." 1242 <<
"If this happened in a place where you are reading a mesh " 1243 <<
"from a file, then you need to investigate why such a quad " 1244 <<
"ended up in the input file. A typical case is a geometry " 1245 <<
"that consisted of multiple parts and for which the mesh " 1246 <<
"generator program assumes that the interface between " 1247 <<
"two parts is a boundary when that isn't supposed to be " 1248 <<
"the case, or where the mesh generator simply assigns " 1249 <<
"'geometry indicators' to quads at the surface of " 1250 <<
"a part that are not supposed to be interpreted as " 1251 <<
"'boundary indicators'.");
1260 <<
"In SubCellData the line info of the line with vertex indices " << arg1
1261 <<
" and " << arg2 <<
" appears more than once. " 1262 <<
"This is not allowed.");
1272 <<
"In SubCellData the line info of the line with vertex indices " << arg1
1273 <<
" and " << arg2 <<
" appears multiple times with different (valid) " 1274 << arg3 <<
". This is not allowed.");
1286 <<
"In SubCellData the quad info of the quad with line indices " << arg1
1287 <<
", " << arg2 <<
", " << arg3 <<
" and " << arg4
1288 <<
" appears multiple times with different (valid) " << arg5
1289 <<
". This is not allowed.");
1399 template <
int dim,
int spacedim>
1403 const unsigned int level_objects,
1406 using line_iterator =
1410 if (level_objects > 0)
1412 for (
unsigned int level = 0; level < level_objects; ++level)
1413 if (triangulation.begin(level) != triangulation.end(level))
1432 for (
unsigned int level = 0; level < number_cache.
n_levels; ++level)
1438 line_iterator line = triangulation.begin_line(level),
1440 (level == number_cache.
n_levels - 1 ?
1441 line_iterator(triangulation.end_line()) :
1442 triangulation.begin_line(level + 1));
1443 for (; line != endc; ++line)
1446 if (line->has_children() ==
false)
1462 line_iterator line = triangulation.begin_line(),
1463 endc = triangulation.end_line();
1464 for (; line != endc; ++line)
1467 if (line->has_children() ==
false)
1487 template <
int dim,
int spacedim>
1491 const unsigned int level_objects,
1502 &compute_number_cache<dim, spacedim>),
1508 using quad_iterator =
1525 unsigned int n_levels = 0;
1526 if (level_objects > 0)
1528 for (
unsigned int level = 0; level < level_objects; ++level)
1529 if (triangulation.begin(level) != triangulation.end(level))
1530 n_levels = level + 1;
1535 for (
unsigned int level = 0; level < n_levels; ++level)
1541 quad_iterator quad = triangulation.begin_quad(level),
1543 (level == n_levels - 1 ?
1544 quad_iterator(triangulation.end_quad()) :
1545 triangulation.begin_quad(level + 1));
1546 for (; quad != endc; ++quad)
1549 if (quad->has_children() ==
false)
1565 quad_iterator quad = triangulation.begin_quad(),
1566 endc = triangulation.end_quad();
1567 for (; quad != endc; ++quad)
1570 if (quad->has_children() ==
false)
1576 update_lines.
join();
1594 template <
int dim,
int spacedim>
1598 const unsigned int level_objects,
1609 &compute_number_cache<dim, spacedim>),
1615 using hex_iterator =
1633 unsigned int n_levels = 0;
1634 if (level_objects > 0)
1636 for (
unsigned int level = 0; level < level_objects; ++level)
1637 if (triangulation.begin(level) != triangulation.end(level))
1638 n_levels = level + 1;
1643 for (
unsigned int level = 0; level < n_levels; ++level)
1649 hex_iterator hex = triangulation.begin_hex(level),
1650 endc = (level == n_levels - 1 ?
1651 hex_iterator(triangulation.end_hex()) :
1652 triangulation.begin_hex(level + 1));
1653 for (; hex != endc; ++hex)
1656 if (hex->has_children() ==
false)
1672 hex_iterator hex = triangulation.begin_hex(),
1673 endc = triangulation.end_hex();
1674 for (; hex != endc; ++hex)
1677 if (hex->has_children() ==
false)
1683 update_quads_and_lines.
join();
1694 template <
int spacedim>
1710 const unsigned int dim = 1;
1713 triangulation.vertices = v;
1714 triangulation.vertices_used = std::vector<bool>(v.size(),
true);
1722 if (dim == spacedim)
1724 for (
unsigned int cell_no = 0; cell_no < cells.size(); ++cell_no)
1730 if (!triangulation.check_for_distorted_cells)
1732 const double cell_measure =
1733 GridTools::cell_measure<1>(triangulation.vertices,
1734 cells[cell_no].vertices);
1746 std::vector<std::vector<int>> lines_at_vertex(v.size());
1749 triangulation.levels.push_back(
1750 std_cxx14::make_unique<
1752 triangulation.levels[0]->reserve_space(cells.size(), dim, spacedim);
1753 triangulation.levels[0]->cells.reserve_space(0, cells.size());
1756 typename Triangulation<dim, spacedim>::raw_line_iterator
1757 next_free_line = triangulation.begin_raw_line();
1758 for (
unsigned int cell = 0; cell < cells.size(); ++cell)
1760 while (next_free_line->used())
1763 next_free_line->set(
1764 internal::TriangulationImplementation ::TriaObject<1>(
1765 cells[cell].vertices[0], cells[cell].vertices[1]));
1766 next_free_line->set_used_flag();
1767 next_free_line->set_material_id(cells[cell].material_id);
1768 next_free_line->set_manifold_id(cells[cell].manifold_id);
1770 next_free_line->set_subdomain_id(0);
1774 lines_at_vertex[cells[cell].vertices[0]].push_back(cell);
1775 lines_at_vertex[cells[cell].vertices[1]].push_back(cell);
1781 unsigned int boundary_nodes = 0;
1782 for (
const auto &line : lines_at_vertex)
1783 switch (line.size())
1796 "You have a vertex in your triangulation " 1797 "at which more than two cells come together. " 1798 "(For one dimensional triangulation, cells are " 1801 "This is not currently supported because the " 1802 "Triangulation class makes the assumption that " 1803 "every cell has zero or one neighbors behind " 1804 "each face (here, behind each vertex), but in your " 1805 "situation there would be more than one." 1807 "Support for this is not currently implemented. " 1808 "If you need to work with triangulations where " 1809 "more than two cells come together at a vertex, " 1810 "duplicate the vertices once per cell (i.e., put " 1811 "multiple vertices at the same physical location, " 1812 "but using different vertex indices for each) " 1813 "and then ensure continuity of the solution by " 1814 "explicitly creating constraints that the degrees " 1815 "of freedom at these vertices have the same " 1816 "value, using the AffineConstraints class."));
1828 AssertThrow(((spacedim == 1) && (boundary_nodes == 2)) ||
1830 ExcMessage(
"The Triangulation has too many end points"));
1837 triangulation.begin_active_line();
1839 for (; line != triangulation.end(); ++line)
1841 for (
unsigned int vertex = 0;
1842 vertex < GeometryInfo<dim>::vertices_per_cell;
1849 if (lines_at_vertex[line->vertex_index(vertex)][0] == line->index())
1850 if (lines_at_vertex[line->vertex_index(vertex)].size() == 2)
1853 neighbor(&triangulation,
1855 lines_at_vertex[line->vertex_index(vertex)][1]);
1856 line->set_neighbor(vertex, neighbor);
1862 line->set_neighbor(vertex, triangulation.end());
1869 neighbor(&triangulation,
1871 lines_at_vertex[line->vertex_index(vertex)][0]);
1872 line->set_neighbor(vertex, neighbor);
1879 triangulation.vertex_to_boundary_id_map_1d->clear();
1880 triangulation.vertex_to_manifold_id_map_1d->clear();
1882 triangulation.begin_active();
1883 cell != triangulation.end();
1885 for (
unsigned int f = 0; f < GeometryInfo<dim>::faces_per_cell; ++f)
1887 (*triangulation.vertex_to_manifold_id_map_1d)
1890 if (cell->at_boundary(f))
1891 (*triangulation.vertex_to_boundary_id_map_1d)
1892 [cell->face(f)->vertex_index()] = f;
1904 template <
int spacedim>
1914 const unsigned int dim = 2;
1917 triangulation.vertices = v;
1918 triangulation.vertices_used = std::vector<bool>(v.size(),
true);
1926 if (dim == spacedim)
1928 for (
unsigned int cell_no = 0; cell_no < cells.size(); ++cell_no)
1931 if (!triangulation.check_for_distorted_cells)
1933 const double cell_measure =
1934 GridTools::cell_measure<2>(triangulation.vertices,
1935 cells[cell_no].vertices);
1956 std::map<std::pair<int, int>,
1959 for (
unsigned int cell = 0; cell < cells.size(); ++cell)
1961 for (
const auto vertex : cells[cell].vertices)
1962 AssertThrow(vertex < triangulation.vertices.size(),
1965 triangulation.vertices.size()));
1967 for (
unsigned int line = 0;
1968 line < GeometryInfo<dim>::faces_per_cell;
1974 std::pair<int, int> line_vertices(
1997 line_vertices.second, line_vertices.first)) ==
2005 needed_lines[line_vertices] = triangulation.end_line();
2013 std::vector<unsigned short int> vertex_touch_count(v.size(), 0);
2015 std::pair<int, int>,
2017 for (i = needed_lines.begin(); i != needed_lines.end(); ++i)
2021 ++vertex_touch_count[i->first.first];
2022 ++vertex_touch_count[i->first.second];
2029 AssertThrow(*(std::min_element(vertex_touch_count.begin(),
2030 vertex_touch_count.end())) >= 2,
2032 "During creation of a triangulation, a part of the " 2033 "algorithm encountered a vertex that is part of only " 2034 "a single adjacent line. However, in 2d, every vertex " 2035 "needs to be at least part of two lines."));
2039 triangulation.levels.push_back(
2040 std_cxx14::make_unique<
2042 triangulation.faces = std_cxx14::make_unique<
2044 triangulation.levels[0]->reserve_space(cells.size(), dim, spacedim);
2045 triangulation.faces->lines.reserve_space(0, needed_lines.size());
2046 triangulation.levels[0]->cells.reserve_space(0, cells.size());
2050 typename Triangulation<dim, spacedim>::raw_line_iterator line =
2051 triangulation.begin_raw_line();
2053 std::pair<int, int>,
2055 for (i = needed_lines.begin(); line != triangulation.end_line();
2059 i->first.first, i->first.second));
2060 line->set_used_flag();
2061 line->clear_user_flag();
2072 std::vector<typename Triangulation<dim, spacedim>::cell_iterator>>
2078 triangulation.begin_raw_quad();
2079 for (
unsigned int c = 0; c < cells.size(); ++c, ++cell)
2083 for (
unsigned int line = 0;
2084 line < GeometryInfo<dim>::lines_per_cell;
2086 lines[line] = needed_lines[std::make_pair(
2096 lines[3]->index()));
2098 cell->set_used_flag();
2099 cell->set_material_id(cells[c].material_id);
2100 cell->set_manifold_id(cells[c].manifold_id);
2101 cell->clear_user_data();
2102 cell->set_subdomain_id(0);
2107 for (
const auto &line : lines)
2108 adjacent_cells[line->index()].push_back(cell);
2114 triangulation.begin_line();
2115 line != triangulation.end_line();
2118 const unsigned int n_adj_cells =
2119 adjacent_cells[line->index()].size();
2127 AssertThrow((n_adj_cells >= 1) && (n_adj_cells <= 2),
2133 (n_adj_cells >= 1) && (n_adj_cells <= 2),
2134 ExcMessage(
"You have a line in your triangulation at which " 2135 "more than two cells come together." 2137 "This is not currently supported because the " 2138 "Triangulation class makes the assumption that " 2139 "every cell has zero or one neighbors behind each " 2140 "face (here, behind each line), but in your " 2141 "situation there would be more than one." 2143 "Support for this is not currently implemented. " 2144 "If you need to work with triangulations where " 2145 "more than two cells come together at a line, " 2146 "duplicate the vertices once per cell (i.e., put " 2147 "multiple vertices at the same physical location, " 2148 "but using different vertex indices for each) " 2149 "and then ensure continuity of the solution by " 2150 "explicitly creating constraints that the degrees " 2151 "of freedom at these lines have the same " 2152 "value, using the AffineConstraints class."));
2157 line->set_boundary_id_internal(
2166 std::pair<int, int> line_vertices(
2167 std::make_pair(subcell_line.vertices[0],
2168 subcell_line.vertices[1]));
2169 if (needed_lines.find(line_vertices) != needed_lines.end())
2171 line = needed_lines[line_vertices];
2175 std::swap(line_vertices.first, line_vertices.second);
2176 if (needed_lines.find(line_vertices) != needed_lines.end())
2177 line = needed_lines[line_vertices];
2182 line_vertices.second));
2187 line->boundary_id() !=
2190 line_vertices.second));
2195 line->manifold_id() == subcell_line.manifold_id,
2197 line_vertices.second,
2199 line->set_manifold_id(subcell_line.manifold_id);
2208 line->vertex_index(1),
2209 subcell_line.boundary_id));
2210 line->set_boundary_id_internal(subcell_line.boundary_id);
2217 triangulation.begin();
2218 cell != triangulation.end();
2220 for (
unsigned int side = 0; side < 4; ++side)
2221 if (adjacent_cells[cell->line(side)->index()][0] == cell)
2225 if (adjacent_cells[cell->line(side)->index()].size() == 2)
2229 side, adjacent_cells[cell->line(side)->index()][1]);
2235 cell->set_neighbor(side,
2236 adjacent_cells[cell->line(side)->index()][0]);
2284 template <
int spacedim>
2294 const unsigned int dim = 3;
2297 triangulation.vertices = v;
2298 triangulation.vertices_used = std::vector<bool>(v.size(),
true);
2304 for (
unsigned int cell_no = 0; cell_no < cells.size(); ++cell_no)
2307 if (!triangulation.check_for_distorted_cells)
2309 const double cell_measure =
2310 GridTools::cell_measure<3>(triangulation.vertices,
2311 cells[cell_no].vertices);
2336 typename std::map<std::pair<int, int>,
2339 for (
unsigned int cell = 0; cell < cells.size(); ++cell)
2343 for (
const auto vertex : cells[cell].vertices)
2344 AssertThrow(vertex < triangulation.vertices.size(),
2347 triangulation.vertices.size()));
2349 for (
unsigned int line = 0;
2350 line < GeometryInfo<dim>::lines_per_cell;
2359 std::pair<int, int> line_vertices(
2368 if ((needed_lines.find(std::make_pair(line_vertices.second,
2369 line_vertices.first)) ==
2370 needed_lines.end()))
2376 needed_lines[line_vertices] = triangulation.end_line();
2388 std::vector<unsigned short int> vertex_touch_count(v.size(), 0);
2390 std::pair<int, int>,
2392 for (i = needed_lines.begin(); i != needed_lines.end(); ++i)
2396 ++vertex_touch_count[i->first.first];
2397 ++vertex_touch_count[i->first.second];
2405 *(std::min_element(vertex_touch_count.begin(),
2406 vertex_touch_count.end())) >= 3,
2408 "During creation of a triangulation, a part of the " 2409 "algorithm encountered a vertex that is part of only " 2410 "one or two adjacent lines. However, in 3d, every vertex " 2411 "needs to be at least part of three lines."));
2419 triangulation.levels.push_back(
2420 std_cxx14::make_unique<
2422 triangulation.faces = std_cxx14::make_unique<
2424 triangulation.levels[0]->reserve_space(cells.size(), dim, spacedim);
2425 triangulation.faces->lines.reserve_space(0, needed_lines.size());
2429 typename Triangulation<dim, spacedim>::raw_line_iterator line =
2430 triangulation.begin_raw_line();
2432 std::pair<int, int>,
2434 for (i = needed_lines.begin(); line != triangulation.end_line();
2438 i->first.first, i->first.second));
2439 line->set_used_flag();
2440 line->clear_user_flag();
2464 std::map<internal::TriangulationImplementation::TriaObject<2>,
2465 std::pair<typename Triangulation<dim, spacedim>::quad_iterator,
2466 std::array<bool, GeometryInfo<dim>::lines_per_face>>,
2469 for (
const auto &cell : cells)
2497 std::array<bool, GeometryInfo<dim>::lines_per_face> orientation;
2499 for (
unsigned int line = 0;
2500 line < GeometryInfo<dim>::lines_per_cell;
2503 line_list[line] = std::pair<int, int>(
2508 inverse_line_list[line] = std::pair<int, int>(
2515 for (
unsigned int face = 0;
2516 face < GeometryInfo<dim>::faces_per_cell;
2527 for (
unsigned int l = 0; l < GeometryInfo<dim>::lines_per_face;
2529 if (needed_lines.find(
2531 face, l)]) == needed_lines.end())
2535 dim>::face_to_cell_lines(face, l)]]
2537 orientation[l] =
true;
2543 dim>::face_to_cell_lines(face, l)]]
2545 orientation[l] =
false;
2594 test_quad_1(quad.
face(2),
2600 test_quad_2(quad.
face(0),
2606 test_quad_3(quad.
face(3),
2612 test_quad_4(quad.
face(1),
2618 test_quad_5(quad.
face(2),
2624 test_quad_6(quad.
face(1),
2630 test_quad_7(quad.
face(3),
2636 if (needed_quads.find(test_quad_1) == needed_quads.end() &&
2637 needed_quads.find(test_quad_2) == needed_quads.end() &&
2638 needed_quads.find(test_quad_3) == needed_quads.end() &&
2639 needed_quads.find(test_quad_4) == needed_quads.end() &&
2640 needed_quads.find(test_quad_5) == needed_quads.end() &&
2641 needed_quads.find(test_quad_6) == needed_quads.end() &&
2642 needed_quads.find(test_quad_7) == needed_quads.end())
2643 needed_quads[quad] =
2644 std::make_pair(triangulation.end_quad(), orientation);
2654 triangulation.faces->quads.reserve_space(0, needed_quads.size());
2657 typename Triangulation<dim, spacedim>::raw_quad_iterator quad =
2658 triangulation.begin_raw_quad();
2661 std::pair<typename Triangulation<dim, spacedim>::quad_iterator,
2662 std::array<bool, GeometryInfo<dim>::lines_per_face>>,
2664 for (q = needed_quads.begin(); quad != triangulation.end_quad();
2667 quad->set(q->first);
2668 quad->set_used_flag();
2669 quad->clear_user_flag();
2672 quad->set_line_orientation(0, q->second.second[0]);
2673 quad->set_line_orientation(1, q->second.second[1]);
2674 quad->set_line_orientation(2, q->second.second[2]);
2675 quad->set_line_orientation(3, q->second.second[3]);
2680 q->second.first = quad;
2686 triangulation.levels[0]->cells.reserve_space(cells.size());
2692 std::vector<typename Triangulation<dim, spacedim>::cell_iterator>>
2698 triangulation.begin_raw_hex();
2699 for (
unsigned int c = 0; c < cells.size(); ++c, ++cell)
2714 unsigned int face_line_list[4];
2715 for (
unsigned int line = 0;
2716 line < GeometryInfo<dim>::lines_per_cell;
2719 line_list[line] = std::make_pair(
2724 inverse_line_list[line] = std::pair<int, int>(
2741 for (
unsigned int face = 0;
2742 face < GeometryInfo<dim>::faces_per_cell;
2745 for (
unsigned int l = 0;
2746 l < GeometryInfo<dim>::lines_per_face;
2749 dim>::face_to_cell_lines(face, l)]) ==
2753 dim>::face_to_cell_lines(face, l)]]
2758 dim>::face_to_cell_lines(face, l)]]
2767 if (needed_quads.find(quad) != needed_quads.end())
2778 face_iterator[face] = needed_quads[quad].first;
2779 face_orientation[face] =
true;
2780 face_flip[face] =
false;
2781 face_rotation[face] =
false;
2809 test_quad_4(quad.
face(1),
2827 test_quad_7(quad.
face(3),
2833 if (needed_quads.find(test_quad_1) != needed_quads.end())
2835 face_iterator[face] = needed_quads[test_quad_1].first;
2836 face_orientation[face] =
false;
2837 face_flip[face] =
false;
2838 face_rotation[face] =
false;
2840 else if (needed_quads.find(test_quad_2) !=
2843 face_iterator[face] = needed_quads[test_quad_2].first;
2844 face_orientation[face] =
false;
2845 face_flip[face] =
false;
2846 face_rotation[face] =
true;
2848 else if (needed_quads.find(test_quad_3) !=
2851 face_iterator[face] = needed_quads[test_quad_3].first;
2852 face_orientation[face] =
false;
2853 face_flip[face] =
true;
2854 face_rotation[face] =
false;
2856 else if (needed_quads.find(test_quad_4) !=
2859 face_iterator[face] = needed_quads[test_quad_4].first;
2860 face_orientation[face] =
false;
2861 face_flip[face] =
true;
2862 face_rotation[face] =
true;
2864 else if (needed_quads.find(test_quad_5) !=
2867 face_iterator[face] = needed_quads[test_quad_5].first;
2868 face_orientation[face] =
true;
2869 face_flip[face] =
false;
2870 face_rotation[face] =
true;
2872 else if (needed_quads.find(test_quad_6) !=
2875 face_iterator[face] = needed_quads[test_quad_6].first;
2876 face_orientation[face] =
true;
2877 face_flip[face] =
true;
2878 face_rotation[face] =
false;
2880 else if (needed_quads.find(test_quad_7) !=
2883 face_iterator[face] = needed_quads[test_quad_7].first;
2884 face_orientation[face] =
true;
2885 face_flip[face] =
true;
2886 face_rotation[face] =
true;
2900 cell->set(internal::TriangulationImplementation ::TriaObject<3>(
2901 face_iterator[0]->index(),
2902 face_iterator[1]->index(),
2903 face_iterator[2]->index(),
2904 face_iterator[3]->index(),
2905 face_iterator[4]->index(),
2906 face_iterator[5]->index()));
2908 cell->set_used_flag();
2909 cell->set_material_id(cells[c].material_id);
2910 cell->set_manifold_id(cells[c].manifold_id);
2911 cell->clear_user_flag();
2912 cell->clear_user_data();
2913 cell->set_subdomain_id(0);
2917 for (
unsigned int quad = 0;
2918 quad < GeometryInfo<dim>::faces_per_cell;
2921 cell->set_face_orientation(quad, face_orientation[quad]);
2922 cell->set_face_flip(quad, face_flip[quad]);
2923 cell->set_face_rotation(quad, face_rotation[quad]);
2930 for (
const auto &quad : face_iterator)
2931 adjacent_cells[quad->index()].push_back(cell);
2948 std::multimap<unsigned int, std::pair<unsigned int, unsigned int>>
2950 for (
unsigned int face = 0;
2951 face < GeometryInfo<dim>::faces_per_cell;
2953 for (
unsigned int line = 0;
2954 line < GeometryInfo<dim>::lines_per_face;
2956 cell_to_face_lines.insert(
2957 std::pair<
unsigned int,
2958 std::pair<unsigned int, unsigned int>>(
2960 std::pair<unsigned int, unsigned int>(face, line)));
2961 std::multimap<
unsigned int,
2962 std::pair<unsigned int, unsigned int>>::
2963 const_iterator map_iter = cell_to_face_lines.begin();
2965 for (; map_iter != cell_to_face_lines.end(); ++map_iter)
2967 const unsigned int cell_line = map_iter->first;
2968 const unsigned int face1 = map_iter->second.first;
2969 const unsigned int line1 = map_iter->second.second;
2971 Assert(map_iter != cell_to_face_lines.end(),
2973 Assert(map_iter->first == cell_line,
2975 const unsigned int face2 = map_iter->second.first;
2976 const unsigned int line2 = map_iter->second.second;
2983 Assert(face_iterator[face1]->line(
2986 face_orientation[face1],
2988 face_rotation[face1])) ==
2989 face_iterator[face2]->line(
2992 face_orientation[face2],
2994 face_rotation[face2])),
3006 triangulation.begin_quad();
3007 quad != triangulation.end_quad();
3010 const unsigned int n_adj_cells =
3011 adjacent_cells[quad->index()].size();
3014 AssertThrow((n_adj_cells >= 1) && (n_adj_cells <= 2),
3019 quad->set_boundary_id_internal(
3034 triangulation.begin_line();
3035 line != triangulation.end_line();
3061 triangulation.begin_quad();
3062 quad != triangulation.end_quad();
3064 if (quad->at_boundary())
3066 for (
unsigned int l = 0; l < 4; ++l)
3070 line->set_boundary_id_internal(0);
3082 std::pair<int, int> line_vertices(
3083 std::make_pair(subcell_line.vertices[0],
3084 subcell_line.vertices[1]));
3085 if (needed_lines.find(line_vertices) != needed_lines.end())
3088 line = needed_lines[line_vertices];
3094 std::swap(line_vertices.first, line_vertices.second);
3095 if (needed_lines.find(line_vertices) != needed_lines.end())
3096 line = needed_lines[line_vertices];
3101 line_vertices.second));
3104 if (line->at_boundary())
3109 line->boundary_id() == subcell_line.boundary_id,
3111 line_vertices.second,
3114 line->set_boundary_id_internal(subcell_line.boundary_id);
3118 line->manifold_id() == subcell_line.manifold_id,
3120 line_vertices.second,
3122 line->set_manifold_id(subcell_line.manifold_id);
3139 for (
unsigned int i = 0; i < 4; ++i)
3141 std::pair<int, int> line_vertices(
3151 if (needed_lines.find(line_vertices) != needed_lines.end())
3152 line[i] = needed_lines[line_vertices];
3157 std::swap(line_vertices.first, line_vertices.second);
3158 if (needed_lines.find(line_vertices) != needed_lines.end())
3159 line[i] = needed_lines[line_vertices];
3165 line_vertices.second));
3204 static const unsigned int lex2cclock[4] = {3, 1, 0, 2};
3211 line_counterclock[4];
3212 for (
unsigned int i = 0; i < 4; ++i)
3213 line_counterclock[lex2cclock[i]] = line[i];
3214 unsigned int n_rotations = 0;
3215 bool not_found_quad_1;
3216 while ((not_found_quad_1 = (needed_quads.find(quad_compare_1) ==
3217 needed_quads.end())) &&
3218 (needed_quads.find(quad_compare_2) == needed_quads.end()) &&
3223 std::rotate(line_counterclock,
3224 line_counterclock + 1,
3225 line_counterclock + 4);
3229 for (
unsigned int i = 0; i < 4; ++i)
3232 i, line_counterclock[lex2cclock[i]]->index());
3234 (i + 2) % 4, line_counterclock[lex2cclock[i]]->index());
3246 if (not_found_quad_1)
3247 quad = needed_quads[quad_compare_2].first;
3249 quad = needed_quads[quad_compare_1].first;
3253 if (quad->at_boundary())
3258 quad->boundary_id() == subcell_quad.boundary_id,
3265 quad->set_boundary_id_internal(subcell_quad.boundary_id);
3269 AssertThrow(quad->manifold_id() == subcell_quad.manifold_id,
3276 quad->set_manifold_id(subcell_quad.manifold_id);
3283 triangulation.begin();
3284 cell != triangulation.end();
3286 for (
unsigned int face = 0; face < 6; ++face)
3287 if (adjacent_cells[cell->quad(face)->index()][0] == cell)
3291 if (adjacent_cells[cell->quad(face)->index()].size() == 2)
3295 face, adjacent_cells[cell->quad(face)->index()][1]);
3301 cell->set_neighbor(face,
3302 adjacent_cells[cell->quad(face)->index()][0]);
3321 template <
int spacedim>
3325 std::vector<unsigned int> &,
3326 std::vector<unsigned int> &)
3328 const unsigned int dim = 1;
3340 Assert(!cell->child(0)->has_children() &&
3341 !cell->child(1)->has_children(),
3347 if (cell->neighbor(0)->has_children())
3354 neighbor = neighbor->child(1);
3357 Assert(neighbor->neighbor(1) == cell->child(0),
3359 neighbor->set_neighbor(1, cell);
3365 if (neighbor->has_children())
3366 neighbor = neighbor->child(1);
3375 if (cell->neighbor(1)->has_children())
3382 neighbor = neighbor->child(0);
3385 Assert(neighbor->neighbor(0) == cell->child(1),
3387 neighbor->set_neighbor(0, cell);
3393 if (neighbor->has_children())
3394 neighbor = neighbor->child(0);
3404 triangulation.vertices_used[cell->child(0)->vertex_index(1)] =
false;
3410 for (
unsigned int child = 0; child < cell->n_children(); ++child)
3413 cell->child(child)->clear_user_flag();
3414 cell->child(child)->clear_used_flag();
3419 cell->clear_children();
3420 cell->clear_user_flag();
3425 template <
int spacedim>
3429 std::vector<unsigned int> & line_cell_count,
3430 std::vector<unsigned int> &)
3432 const unsigned int dim = 2;
3435 Assert(line_cell_count.size() == triangulation.n_raw_lines(),
3440 std::vector<typename Triangulation<dim, spacedim>::line_iterator>
3443 lines_to_delete.reserve(4 * 2 + 4);
3448 for (
unsigned int c = 0; c < cell->n_children(); ++c)
3452 for (
unsigned int l = 0; l < GeometryInfo<dim>::lines_per_cell; ++l)
3453 --line_cell_count[child->line_index(l)];
3468 .vertices_used[cell->child(0)->line(1)->vertex_index(1)] =
false;
3470 lines_to_delete.push_back(cell->child(0)->line(1));
3471 lines_to_delete.push_back(cell->child(0)->line(3));
3472 lines_to_delete.push_back(cell->child(3)->line(0));
3473 lines_to_delete.push_back(cell->child(3)->line(2));
3477 unsigned int inner_face_no =
3482 lines_to_delete.push_back(cell->child(0)->line(inner_face_no));
3486 for (
unsigned int child = 0; child < cell->n_children(); ++child)
3488 cell->child(child)->clear_user_data();
3489 cell->child(child)->clear_user_flag();
3490 cell->child(child)->clear_used_flag();
3495 cell->clear_children();
3496 cell->clear_refinement_case();
3497 cell->clear_user_flag();
3503 for (
unsigned int line_no = 0;
3504 line_no < GeometryInfo<dim>::lines_per_cell;
3508 cell->line(line_no);
3510 if (line->has_children())
3515 Assert((line_cell_count[line->child_index(0)] == 0 &&
3516 line_cell_count[line->child_index(1)] == 0) ||
3517 (line_cell_count[line->child_index(0)] > 0 &&
3518 line_cell_count[line->child_index(1)] > 0),
3521 if (line_cell_count[line->child_index(0)] == 0)
3523 for (
unsigned int c = 0; c < 2; ++c)
3524 Assert(!line->child(c)->has_children(),
3532 .vertices_used[line->child(0)->vertex_index(1)] =
false;
3534 lines_to_delete.push_back(line->child(0));
3535 lines_to_delete.push_back(line->child(1));
3537 line->clear_children();
3549 typename std::vector<
3551 line = lines_to_delete.
begin(),
3552 endline = lines_to_delete.end();
3553 for (; line != endline; ++line)
3555 (*line)->clear_user_data();
3556 (*line)->clear_user_flag();
3557 (*line)->clear_used_flag();
3563 template <
int spacedim>
3567 std::vector<unsigned int> & line_cell_count,
3568 std::vector<unsigned int> & quad_cell_count)
3570 const unsigned int dim = 3;
3572 Assert(line_cell_count.size() == triangulation.n_raw_lines(),
3574 Assert(quad_cell_count.size() == triangulation.n_raw_quads(),
3582 std::vector<typename Triangulation<dim, spacedim>::line_iterator>
3584 std::vector<typename Triangulation<dim, spacedim>::quad_iterator>
3587 lines_to_delete.reserve(12 * 2 + 6 * 4 + 6);
3588 quads_to_delete.reserve(6 * 4 + 12);
3592 for (
unsigned int c = 0; c < cell->n_children(); ++c)
3596 for (
unsigned int l = 0; l < GeometryInfo<dim>::lines_per_cell; ++
l)
3597 --line_cell_count[child->line_index(l)];
3598 for (
unsigned int f = 0; f < GeometryInfo<dim>::faces_per_cell; ++f)
3599 --quad_cell_count[child->quad_index(f)];
3613 quads_to_delete.push_back(cell->child(0)->face(1));
3616 quads_to_delete.push_back(cell->child(0)->face(3));
3619 quads_to_delete.push_back(cell->child(0)->face(5));
3622 quads_to_delete.push_back(cell->child(0)->face(1));
3623 quads_to_delete.push_back(cell->child(0)->face(3));
3624 quads_to_delete.push_back(cell->child(3)->face(0));
3625 quads_to_delete.push_back(cell->child(3)->face(2));
3627 lines_to_delete.push_back(cell->child(0)->line(11));
3630 quads_to_delete.push_back(cell->child(0)->face(1));
3631 quads_to_delete.push_back(cell->child(0)->face(5));
3632 quads_to_delete.push_back(cell->child(3)->face(0));
3633 quads_to_delete.push_back(cell->child(3)->face(4));
3635 lines_to_delete.push_back(cell->child(0)->line(5));
3638 quads_to_delete.push_back(cell->child(0)->face(3));
3639 quads_to_delete.push_back(cell->child(0)->face(5));
3640 quads_to_delete.push_back(cell->child(3)->face(2));
3641 quads_to_delete.push_back(cell->child(3)->face(4));
3643 lines_to_delete.push_back(cell->child(0)->line(7));
3646 quads_to_delete.push_back(cell->child(0)->face(1));
3647 quads_to_delete.push_back(cell->child(2)->face(1));
3648 quads_to_delete.push_back(cell->child(4)->face(1));
3649 quads_to_delete.push_back(cell->child(6)->face(1));
3651 quads_to_delete.push_back(cell->child(0)->face(3));
3652 quads_to_delete.push_back(cell->child(1)->face(3));
3653 quads_to_delete.push_back(cell->child(4)->face(3));
3654 quads_to_delete.push_back(cell->child(5)->face(3));
3656 quads_to_delete.push_back(cell->child(0)->face(5));
3657 quads_to_delete.push_back(cell->child(1)->face(5));
3658 quads_to_delete.push_back(cell->child(2)->face(5));
3659 quads_to_delete.push_back(cell->child(3)->face(5));
3661 lines_to_delete.push_back(cell->child(0)->line(5));
3662 lines_to_delete.push_back(cell->child(0)->line(7));
3663 lines_to_delete.push_back(cell->child(0)->line(11));
3664 lines_to_delete.push_back(cell->child(7)->line(0));
3665 lines_to_delete.push_back(cell->child(7)->line(2));
3666 lines_to_delete.push_back(cell->child(7)->line(8));
3672 triangulation.vertices_used[cell->child(0)->vertex_index(7)] =
3684 for (
unsigned int child = 0; child < cell->n_children(); ++child)
3686 cell->child(child)->clear_user_data();
3687 cell->child(child)->clear_user_flag();
3689 for (
unsigned int f = 0; f < GeometryInfo<dim>::faces_per_cell; ++f)
3694 cell->child(child)->set_face_orientation(f,
true);
3695 cell->child(child)->set_face_flip(f,
false);
3696 cell->child(child)->set_face_rotation(f,
false);
3699 cell->child(child)->clear_used_flag();
3704 cell->clear_children();
3705 cell->clear_refinement_case();
3706 cell->clear_user_flag();
3716 for (
unsigned int quad_no = 0;
3717 quad_no < GeometryInfo<dim>::faces_per_cell;
3721 cell->face(quad_no);
3725 quad->has_children()) ||
3730 switch (quad->refinement_case())
3742 Assert((quad_cell_count[quad->child_index(0)] == 0 &&
3743 quad_cell_count[quad->child_index(1)] == 0) ||
3744 (quad_cell_count[quad->child_index(0)] > 0 &&
3745 quad_cell_count[quad->child_index(1)] > 0),
3751 unsigned int deleted_grandchildren = 0;
3752 unsigned int number_of_child_refinements = 0;
3754 for (
unsigned int c = 0; c < 2; ++c)
3755 if (quad->child(c)->has_children())
3757 ++number_of_child_refinements;
3762 (quad_cell_count[quad->child(c)->child_index(0)] ==
3764 quad_cell_count[quad->child(c)->child_index(1)] ==
3766 (quad_cell_count[quad->child(c)->child_index(0)] >
3768 quad_cell_count[quad->child(c)->child_index(1)] >
3771 if (quad_cell_count[quad->child(c)->child_index(0)] ==
3778 Assert(quad->refinement_case() +
3779 quad->child(c)->refinement_case() ==
3787 quads_to_delete.push_back(
3788 quad->child(c)->child(0));
3789 quads_to_delete.push_back(
3790 quad->child(c)->child(1));
3791 if (quad->child(c)->refinement_case() ==
3793 lines_to_delete.push_back(
3794 quad->child(c)->child(0)->line(1));
3796 lines_to_delete.push_back(
3797 quad->child(c)->child(0)->line(3));
3798 quad->child(c)->clear_children();
3799 quad->child(c)->clear_refinement_case();
3800 ++deleted_grandchildren;
3808 if (number_of_child_refinements > 0 &&
3809 deleted_grandchildren == number_of_child_refinements)
3814 middle_line = quad->child(0)->line(1);
3816 middle_line = quad->child(0)->line(3);
3818 lines_to_delete.push_back(middle_line->child(0));
3819 lines_to_delete.push_back(middle_line->child(1));
3821 .vertices_used[middle_vertex_index<dim, spacedim>(
3822 middle_line)] =
false;
3823 middle_line->clear_children();
3828 if (quad_cell_count[quad->child_index(0)] == 0)
3834 quads_to_delete.push_back(quad->child(0));
3835 quads_to_delete.push_back(quad->child(1));
3837 lines_to_delete.push_back(quad->child(0)->line(1));
3839 lines_to_delete.push_back(quad->child(0)->line(3));
3859 if (quad->child(0)->has_children())
3860 if (quad->refinement_case() ==
3894 quad_iterator switch_1 =
3895 quad->child(0)->child(1),
3897 quad->child(1)->child(0);
3899 Assert(!switch_1->has_children(),
3901 Assert(!switch_2->has_children(),
3904 const int switch_1_index = switch_1->index();
3905 const int switch_2_index = switch_2->index();
3906 for (
unsigned int l = 0;
3907 l < triangulation.levels.size();
3909 for (
unsigned int h = 0;
3910 h < triangulation.levels[
l]
3911 ->cells.cells.size();
3913 for (
unsigned int q = 0;
3914 q < GeometryInfo<dim>::faces_per_cell;
3917 const int index = triangulation.levels[
l]
3920 if (index == switch_1_index)
3921 triangulation.levels[
l]
3923 .set_face(q, switch_2_index);
3924 else if (index == switch_2_index)
3925 triangulation.levels[
l]
3927 .set_face(q, switch_1_index);
3932 const int switch_1_lines[4] = {
3933 static_cast<signed int>(
3934 switch_1->line_index(0)),
3935 static_cast<signed int>(
3936 switch_1->line_index(1)),
3937 static_cast<signed int>(
3938 switch_1->line_index(2)),
3939 static_cast<signed int>(
3940 switch_1->line_index(3))};
3941 const bool switch_1_line_orientations[4] = {
3942 switch_1->line_orientation(0),
3943 switch_1->line_orientation(1),
3944 switch_1->line_orientation(2),
3945 switch_1->line_orientation(3)};
3947 switch_1->boundary_id();
3948 const unsigned int switch_1_user_index =
3949 switch_1->user_index();
3950 const bool switch_1_user_flag =
3951 switch_1->user_flag_set();
3955 TriaObject<2>(switch_2->line_index(0),
3956 switch_2->line_index(1),
3957 switch_2->line_index(2),
3958 switch_2->line_index(3)));
3959 switch_1->set_line_orientation(
3960 0, switch_2->line_orientation(0));
3961 switch_1->set_line_orientation(
3962 1, switch_2->line_orientation(1));
3963 switch_1->set_line_orientation(
3964 2, switch_2->line_orientation(2));
3965 switch_1->set_line_orientation(
3966 3, switch_2->line_orientation(3));
3967 switch_1->set_boundary_id_internal(
3968 switch_2->boundary_id());
3969 switch_1->set_manifold_id(
3970 switch_2->manifold_id());
3971 switch_1->set_user_index(switch_2->user_index());
3972 if (switch_2->user_flag_set())
3973 switch_1->set_user_flag();
3975 switch_1->clear_user_flag();
3979 TriaObject<2>(switch_1_lines[0],
3982 switch_1_lines[3]));
3983 switch_2->set_line_orientation(
3984 0, switch_1_line_orientations[0]);
3985 switch_2->set_line_orientation(
3986 1, switch_1_line_orientations[1]);
3987 switch_2->set_line_orientation(
3988 2, switch_1_line_orientations[2]);
3989 switch_2->set_line_orientation(
3990 3, switch_1_line_orientations[3]);
3991 switch_2->set_boundary_id_internal(
3992 switch_1_boundary_id);
3993 switch_2->set_manifold_id(
3994 switch_1->manifold_id());
3995 switch_2->set_user_index(switch_1_user_index);
3996 if (switch_1_user_flag)
3997 switch_2->set_user_flag();
3999 switch_2->clear_user_flag();
4001 const unsigned int child_0 =
4002 quad->child(0)->child_index(0);
4003 const unsigned int child_2 =
4004 quad->child(1)->child_index(0);
4005 quad->clear_children();
4006 quad->clear_refinement_case();
4007 quad->set_refinement_case(
4009 quad->set_children(0, child_0);
4010 quad->set_children(2, child_2);
4012 quad_cell_count[child_2]);
4027 const unsigned int child_0 =
4028 quad->child(0)->child_index(0);
4029 const unsigned int child_2 =
4030 quad->child(1)->child_index(0);
4031 quad->clear_children();
4032 quad->clear_refinement_case();
4033 quad->set_refinement_case(
4035 quad->set_children(0, child_0);
4036 quad->set_children(2, child_2);
4040 quad->clear_children();
4041 quad->clear_refinement_case();
4052 Assert((quad_cell_count[quad->child_index(0)] == 0 &&
4053 quad_cell_count[quad->child_index(1)] == 0 &&
4054 quad_cell_count[quad->child_index(2)] == 0 &&
4055 quad_cell_count[quad->child_index(3)] == 0) ||
4056 (quad_cell_count[quad->child_index(0)] > 0 &&
4057 quad_cell_count[quad->child_index(1)] > 0 &&
4058 quad_cell_count[quad->child_index(2)] > 0 &&
4059 quad_cell_count[quad->child_index(3)] > 0),
4062 if (quad_cell_count[quad->child_index(0)] == 0)
4068 lines_to_delete.push_back(quad->child(0)->line(1));
4069 lines_to_delete.push_back(quad->child(3)->line(0));
4070 lines_to_delete.push_back(quad->child(0)->line(3));
4071 lines_to_delete.push_back(quad->child(3)->line(2));
4073 for (
unsigned int child = 0; child < quad->n_children();
4075 quads_to_delete.push_back(quad->child(child));
4078 .vertices_used[quad->child(0)->vertex_index(3)] =
4081 quad->clear_children();
4082 quad->clear_refinement_case();
4101 for (
unsigned int line_no = 0;
4102 line_no < GeometryInfo<dim>::lines_per_cell;
4106 cell->line(line_no);
4110 line->has_children()) ||
4115 if (line->has_children())
4120 Assert((line_cell_count[line->child_index(0)] == 0 &&
4121 line_cell_count[line->child_index(1)] == 0) ||
4122 (line_cell_count[line->child_index(0)] > 0 &&
4123 line_cell_count[line->child_index(1)] > 0),
4126 if (line_cell_count[line->child_index(0)] == 0)
4128 for (
unsigned int c = 0; c < 2; ++c)
4129 Assert(!line->child(c)->has_children(),
4137 .vertices_used[line->child(0)->vertex_index(1)] =
false;
4139 lines_to_delete.push_back(line->child(0));
4140 lines_to_delete.push_back(line->child(1));
4142 line->clear_children();
4154 typename std::vector<
4156 line = lines_to_delete.
begin(),
4157 endline = lines_to_delete.end();
4158 for (; line != endline; ++line)
4160 (*line)->clear_user_data();
4161 (*line)->clear_user_flag();
4162 (*line)->clear_used_flag();
4165 typename std::vector<
4167 quad = quads_to_delete.
begin(),
4168 endquad = quads_to_delete.end();
4169 for (; quad != endquad; ++quad)
4171 (*quad)->clear_user_data();
4172 (*quad)->clear_children();
4173 (*quad)->clear_refinement_case();
4174 (*quad)->clear_user_flag();
4175 (*quad)->clear_used_flag();
4197 template <
int spacedim>
4200 unsigned int & next_unused_vertex,
4201 typename Triangulation<2, spacedim>::raw_line_iterator
4207 const unsigned int dim = 2;
4210 cell->clear_refine_flag();
4307 int new_vertices[9];
4308 for (
unsigned int vertex_no = 0; vertex_no < 4; ++vertex_no)
4309 new_vertices[vertex_no] = cell->vertex_index(vertex_no);
4310 for (
unsigned int line_no = 0; line_no < 4; ++line_no)
4311 if (cell->line(line_no)->has_children())
4312 new_vertices[4 + line_no] =
4313 cell->line(line_no)->child(0)->vertex_index(1);
4322 while (triangulation.vertices_used[next_unused_vertex] ==
true)
4323 ++next_unused_vertex;
4325 next_unused_vertex < triangulation.vertices.size(),
4327 "Internal error: During refinement, the triangulation wants to access an element of the 'vertices' array but it turns out that the array is not large enough."));
4328 triangulation.vertices_used[next_unused_vertex] =
true;
4330 new_vertices[8] = next_unused_vertex;
4342 if (dim == spacedim)
4345 triangulation.vertices[next_unused_vertex] = cell->center(
true);
4352 if (cell->user_flag_set())
4355 cell->clear_user_flag();
4356 triangulation.vertices[next_unused_vertex] =
4357 cell->center(
true,
true);
4367 cell->clear_user_flag();
4371 triangulation.vertices[next_unused_vertex] = cell->center(
true);
4377 typename Triangulation<dim, spacedim>::raw_line_iterator new_lines[12];
4378 unsigned int lmin = 8;
4379 unsigned int lmax = 12;
4386 for (
unsigned int l = lmin; l < lmax; ++l)
4388 while (next_unused_line->used() ==
true)
4390 new_lines[l] = next_unused_line;
4394 new_lines[l]->used() ==
false,
4396 "Internal error: We want to use a cell during refinement that should be unused, but turns out not to be."));
4410 for (
unsigned int face_no = 0;
4411 face_no < GeometryInfo<dim>::faces_per_cell;
4413 for (
unsigned int c = 0; c < 2; ++c, ++l)
4414 new_lines[l] = cell->line(face_no)->child(c);
4419 new_vertices[6], new_vertices[8]));
4422 new_vertices[8], new_vertices[7]));
4425 new_vertices[4], new_vertices[8]));
4428 new_vertices[8], new_vertices[5]));
4437 new_lines[0] = cell->line(0);
4438 new_lines[1] = cell->line(1);
4439 new_lines[2] = cell->line(2)->child(0);
4440 new_lines[3] = cell->line(2)->child(1);
4441 new_lines[4] = cell->line(3)->child(0);
4442 new_lines[5] = cell->line(3)->child(1);
4445 new_vertices[6], new_vertices[7]));
4455 new_lines[0] = cell->line(0)->child(0);
4456 new_lines[1] = cell->line(0)->child(1);
4457 new_lines[2] = cell->line(1)->child(0);
4458 new_lines[3] = cell->line(1)->child(1);
4459 new_lines[4] = cell->line(2);
4460 new_lines[5] = cell->line(3);
4463 new_vertices[4], new_vertices[5]));
4466 for (
unsigned int l = lmin; l < lmax; ++l)
4468 new_lines[l]->set_used_flag();
4469 new_lines[l]->clear_user_flag();
4471 new_lines[l]->clear_children();
4473 new_lines[l]->set_boundary_id_internal(
4475 new_lines[l]->set_manifold_id(cell->manifold_id());
4482 while (next_unused_cell->used() ==
true)
4486 for (
unsigned int i = 0; i < n_children; ++i)
4489 next_unused_cell->used() ==
false,
4491 "Internal error: We want to use a cell during refinement that should be unused, but turns out not to be."));
4492 subcells[i] = next_unused_cell;
4494 if (i % 2 == 1 && i < n_children - 1)
4495 while (next_unused_cell->used() ==
true)
4515 new_lines[0]->index(),
4516 new_lines[8]->index(),
4517 new_lines[4]->index(),
4518 new_lines[10]->index()));
4521 new_lines[8]->index(),
4522 new_lines[2]->index(),
4523 new_lines[5]->index(),
4524 new_lines[11]->index()));
4527 new_lines[1]->index(),
4528 new_lines[9]->index(),
4529 new_lines[10]->index(),
4530 new_lines[6]->index()));
4533 new_lines[9]->index(),
4534 new_lines[3]->index(),
4535 new_lines[11]->index(),
4536 new_lines[7]->index()));
4554 new_lines[0]->index(),
4555 new_lines[6]->index(),
4556 new_lines[2]->index(),
4557 new_lines[4]->index()));
4560 new_lines[6]->index(),
4561 new_lines[1]->index(),
4562 new_lines[3]->index(),
4563 new_lines[5]->index()));
4582 new_lines[0]->index(),
4583 new_lines[2]->index(),
4584 new_lines[4]->index(),
4585 new_lines[6]->index()));
4588 new_lines[1]->index(),
4589 new_lines[3]->index(),
4590 new_lines[6]->index(),
4591 new_lines[5]->index()));
4596 for (
unsigned int i = 0; i < n_children; ++i)
4598 subcells[i]->set_used_flag();
4599 subcells[i]->clear_refine_flag();
4600 subcells[i]->clear_user_flag();
4601 subcells[i]->clear_user_data();
4602 subcells[i]->clear_children();
4605 subcells[i]->set_material_id(cell->material_id());
4606 subcells[i]->set_manifold_id(cell->manifold_id());
4607 subcells[i]->set_subdomain_id(subdomainid);
4610 subcells[i]->set_parent(cell->index());
4618 for (
unsigned int i = 0; i < n_children / 2; ++i)
4619 cell->set_children(2 * i, subcells[2 * i]->index());
4621 cell->set_refinement_case(ref_case);
4629 for (
unsigned int c = 0; c < n_children; ++c)
4630 cell->child(c)->set_direction_flag(cell->direction_flag());
4639 template <
int spacedim>
4644 const unsigned int dim = 1;
4651 cell = triangulation.begin_active(triangulation.levels.size() - 1),
4652 endc = triangulation.end();
4653 for (; cell != endc; ++cell)
4655 if (cell->refine_flag_set())
4657 triangulation.levels.push_back(
4658 std_cxx14::make_unique<
4669 unsigned int needed_vertices = 0;
4670 for (
int level = triangulation.levels.size() - 2; level >= 0; --level)
4674 unsigned int flagged_cells = 0;
4676 acell = triangulation.begin_active(level),
4677 aendc = triangulation.begin_active(level + 1);
4678 for (; acell != aendc; ++acell)
4679 if (acell->refine_flag_set())
4684 const unsigned int used_cells = std::count_if(
4685 triangulation.levels[level + 1]->cells.used.begin(),
4686 triangulation.levels[level + 1]->cells.used.end(),
4687 std::bind(std::equal_to<bool>(), std::placeholders::_1,
true));
4692 triangulation.levels[level + 1]->reserve_space(
4699 triangulation.levels[level + 1]->cells.reserve_space(
4702 needed_vertices += flagged_cells;
4707 needed_vertices += std::count_if(triangulation.vertices_used.begin(),
4708 triangulation.vertices_used.end(),
4709 std::bind(std::equal_to<bool>(),
4710 std::placeholders::_1,
4715 if (needed_vertices > triangulation.vertices.size())
4718 triangulation.vertices_used.resize(needed_vertices,
false);
4726 unsigned int next_unused_vertex = 0;
4728 for (
int level = triangulation.levels.size() - 2; level >= 0; --level)
4731 cell = triangulation.begin_active(level),
4732 endc = triangulation.begin_active(level + 1);
4735 next_unused_cell = triangulation.begin_raw(level + 1);
4737 for (; (cell != endc) && (cell->level() == level); ++cell)
4738 if (cell->refine_flag_set())
4741 cell->clear_refine_flag();
4745 while (triangulation.vertices_used[next_unused_vertex] ==
4747 ++next_unused_vertex;
4749 next_unused_vertex < triangulation.vertices.size(),
4751 "Internal error: During refinement, the triangulation wants to access an element of the 'vertices' array but it turns out that the array is not large enough."));
4756 triangulation.vertices[next_unused_vertex] =
4759 triangulation.vertices_used[next_unused_vertex] =
true;
4766 while (next_unused_cell->used() ==
true)
4768 first_child = next_unused_cell;
4769 first_child->set_used_flag();
4770 first_child->clear_user_data();
4773 next_unused_cell->used() ==
false,
4775 "Internal error: We want to use a cell during refinement that should be unused, but turns out not to be."));
4776 second_child = next_unused_cell;
4777 second_child->set_used_flag();
4778 second_child->clear_user_data();
4783 cell->set_children(0, first_child->index());
4784 first_child->clear_children();
4786 internal::TriangulationImplementation ::TriaObject<1>(
4787 cell->vertex_index(0), next_unused_vertex));
4788 first_child->set_material_id(cell->material_id());
4789 first_child->set_manifold_id(cell->manifold_id());
4790 first_child->set_subdomain_id(subdomainid);
4791 first_child->set_direction_flag(cell->direction_flag());
4793 first_child->set_parent(cell->index());
4797 first_child->face(1)->set_manifold_id(cell->manifold_id());
4802 first_child->set_neighbor(1, second_child);
4804 first_child->set_neighbor(0, cell->neighbor(0));
4805 else if (cell->neighbor(0)->active())
4811 Assert(cell->neighbor(0)->level() <= cell->level(),
4813 first_child->set_neighbor(0, cell->neighbor(0));
4819 const unsigned int nbnb = cell->neighbor_of_neighbor(0);
4820 first_child->set_neighbor(0,
4821 cell->neighbor(0)->child(nbnb));
4826 left_neighbor = cell->neighbor(0);
4827 while (left_neighbor->has_children())
4829 left_neighbor = left_neighbor->child(nbnb);
4830 left_neighbor->set_neighbor(nbnb, first_child);
4835 second_child->clear_children();
4837 internal::TriangulationImplementation ::TriaObject<1>(
4838 next_unused_vertex, cell->vertex_index(1)));
4839 second_child->set_neighbor(0, first_child);
4840 second_child->set_material_id(cell->material_id());
4841 second_child->set_manifold_id(cell->manifold_id());
4842 second_child->set_subdomain_id(subdomainid);
4843 second_child->set_direction_flag(cell->direction_flag());
4846 second_child->set_neighbor(1, cell->neighbor(1));
4847 else if (cell->neighbor(1)->active())
4849 Assert(cell->neighbor(1)->level() <= cell->level(),
4851 second_child->set_neighbor(1, cell->neighbor(1));
4856 const unsigned int nbnb = cell->neighbor_of_neighbor(1);
4857 second_child->set_neighbor(
4858 1, cell->neighbor(1)->child(nbnb));
4861 right_neighbor = cell->neighbor(1);
4862 while (right_neighbor->has_children())
4864 right_neighbor = right_neighbor->child(nbnb);
4865 right_neighbor->set_neighbor(nbnb, second_child);
4869 triangulation.signals.post_refinement_on_cell(cell);
4885 template <
int spacedim>
4888 const bool check_for_distorted_cells)
4890 const unsigned int dim = 2;
4897 cell = triangulation.begin_active(triangulation.levels.size() - 1),
4898 endc = triangulation.end();
4899 for (; cell != endc; ++cell)
4901 if (cell->refine_flag_set())
4903 triangulation.levels.push_back(
4904 std_cxx14::make_unique<
4913 triangulation.begin_line();
4914 line != triangulation.end_line();
4917 line->clear_user_flag();
4918 line->clear_user_data();
4923 unsigned int n_single_lines = 0;
4927 unsigned int n_lines_in_pairs = 0;
4933 unsigned int needed_vertices = 0;
4934 for (
int level = triangulation.levels.size() - 2; level >= 0; --level)
4938 unsigned int needed_cells = 0;
4941 cell = triangulation.begin_active(level),
4942 endc = triangulation.begin_active(level + 1);
4943 for (; cell != endc; ++cell)
4944 if (cell->refine_flag_set())
4955 n_single_lines += 4;
4961 triangulation.anisotropic_refinement =
true;
4967 n_single_lines += 1;
4975 for (
unsigned int line_no = 0;
4976 line_no < GeometryInfo<dim>::faces_per_cell;
4980 cell->refine_flag_set(), line_no) ==
4984 line = cell->line(line_no);
4985 if (line->has_children() ==
false)
4986 line->set_user_flag();
4993 const unsigned int used_cells = std::count_if(
4994 triangulation.levels[level + 1]->cells.used.begin(),
4995 triangulation.levels[level + 1]->cells.used.end(),
4996 std::bind(std::equal_to<bool>(), std::placeholders::_1,
true));
5002 triangulation.levels[level + 1]->reserve_space(
5003 used_cells + needed_cells, 2, spacedim);
5007 triangulation.levels[level + 1]->cells.reserve_space(needed_cells,
5013 triangulation.begin_line();
5014 line != triangulation.end_line();
5016 if (line->user_flag_set())
5019 n_lines_in_pairs += 2;
5020 needed_vertices += 1;
5029 triangulation.faces->lines.reserve_space(n_lines_in_pairs, 0);
5032 needed_vertices += std::count_if(triangulation.vertices_used.begin(),
5033 triangulation.vertices_used.end(),
5034 std::bind(std::equal_to<bool>(),
5035 std::placeholders::_1,
5040 if (needed_vertices > triangulation.vertices.size())
5043 triangulation.vertices_used.resize(needed_vertices,
false);
5051 unsigned int next_unused_vertex = 0;
5058 line = triangulation.begin_active_line(),
5059 endl = triangulation.end_line();
5060 typename Triangulation<dim, spacedim>::raw_line_iterator
5061 next_unused_line = triangulation.begin_raw_line();
5063 for (; line != endl; ++line)
5064 if (line->user_flag_set())
5070 while (triangulation.vertices_used[next_unused_vertex] ==
true)
5071 ++next_unused_vertex;
5073 next_unused_vertex < triangulation.vertices.size(),
5075 "Internal error: During refinement, the triangulation wants to access an element of the 'vertices' array but it turns out that the array is not large enough."));
5076 triangulation.vertices_used[next_unused_vertex] =
true;
5078 triangulation.vertices[next_unused_vertex] = line->center(
true);
5083 bool pair_found =
false;
5085 for (; next_unused_line != endl; ++next_unused_line)
5086 if (!next_unused_line->used() &&
5087 !(++next_unused_line)->used())
5100 line->set_children(0, next_unused_line->index());
5103 const typename Triangulation<dim, spacedim>::raw_line_iterator
5104 children[2] = {next_unused_line, ++next_unused_line};
5108 children[0]->used() ==
false,
5110 "Internal error: We want to use a cell during refinement that should be unused, but turns out not to be."));
5112 children[1]->used() ==
false,
5114 "Internal error: We want to use a cell during refinement that should be unused, but turns out not to be."));
5117 internal::TriangulationImplementation ::TriaObject<1>(
5118 line->vertex_index(0), next_unused_vertex));
5120 internal::TriangulationImplementation ::TriaObject<1>(
5121 next_unused_vertex, line->vertex_index(1)));
5123 children[0]->set_used_flag();
5124 children[1]->set_used_flag();
5125 children[0]->clear_children();
5126 children[1]->clear_children();
5129 children[0]->clear_user_flag();
5130 children[1]->clear_user_flag();
5133 children[0]->set_boundary_id_internal(line->boundary_id());
5134 children[1]->set_boundary_id_internal(line->boundary_id());
5136 children[0]->set_manifold_id(line->manifold_id());
5137 children[1]->set_manifold_id(line->manifold_id());
5141 line->clear_user_flag();
5150 triangulation.faces->lines.reserve_space(0, n_single_lines);
5153 cells_with_distorted_children;
5157 typename Triangulation<dim, spacedim>::raw_line_iterator
5158 next_unused_line = triangulation.begin_raw_line();
5161 level < static_cast<int>(triangulation.levels.size()) - 1;
5167 cell = triangulation.begin_active(level),
5168 endc = triangulation.begin_active(level + 1);
5171 next_unused_cell = triangulation.begin_raw(level + 1);
5173 for (; cell != endc; ++cell)
5174 if (cell->refine_flag_set())
5182 if (cell->at_boundary())
5183 cell->set_user_flag();
5193 if ((check_for_distorted_cells ==
true) &&
5194 has_distorted_children(
5196 std::integral_constant<int, dim>(),
5197 std::integral_constant<int, spacedim>()))
5201 triangulation.signals.post_refinement_on_cell(cell);
5205 return cells_with_distorted_children;
5213 template <
int spacedim>
5216 const bool check_for_distorted_cells)
5218 const unsigned int dim = 3;
5231 cell = triangulation.begin_active(triangulation.levels.size() - 1),
5232 endc = triangulation.end();
5233 for (; cell != endc; ++cell)
5235 if (cell->refine_flag_set())
5237 triangulation.levels.push_back(
5238 std_cxx14::make_unique<
5247 triangulation.faces->quads.clear_user_data();
5250 triangulation.begin_line();
5251 line != triangulation.end_line();
5253 line->clear_user_flag();
5255 triangulation.begin_quad();
5256 quad != triangulation.end_quad();
5258 quad->clear_user_flag();
5281 unsigned int needed_vertices = 0;
5282 unsigned int needed_lines_single = 0;
5283 unsigned int needed_quads_single = 0;
5284 unsigned int needed_lines_pair = 0;
5285 unsigned int needed_quads_pair = 0;
5286 for (
int level = triangulation.levels.size() - 2; level >= 0; --level)
5290 unsigned int new_cells = 0;
5293 acell = triangulation.begin_active(level),
5294 aendc = triangulation.begin_active(level + 1);
5295 for (; acell != aendc; ++acell)
5296 if (acell->refine_flag_set())
5306 ++needed_quads_single;
5308 triangulation.anisotropic_refinement =
true;
5314 ++needed_lines_single;
5315 needed_quads_single += 4;
5317 triangulation.anisotropic_refinement =
true;
5322 needed_lines_single += 6;
5323 needed_quads_single += 12;
5337 for (
unsigned int face = 0;
5338 face < GeometryInfo<dim>::faces_per_cell;
5342 aface = acell->face(face);
5349 acell->face_orientation(face),
5350 acell->face_flip(face),
5351 acell->face_rotation(face));
5356 if (face_ref_case ==
5359 if (aface->number_of_children() < 4)
5362 aface->set_user_flag();
5364 else if (aface->refinement_case() != face_ref_case)
5375 Assert(aface->refinement_case() ==
5377 dim - 1>::isotropic_refinement ||
5378 aface->refinement_case() ==
5381 aface->set_user_index(face_ref_case);
5387 for (
unsigned int line = 0;
5388 line < GeometryInfo<dim>::lines_per_cell;
5392 !acell->line(line)->has_children())
5393 acell->line(line)->set_user_flag();
5399 const unsigned int used_cells = std::count_if(
5400 triangulation.levels[level + 1]->cells.used.begin(),
5401 triangulation.levels[level + 1]->cells.used.end(),
5402 std::bind(std::equal_to<bool>(), std::placeholders::_1,
true));
5408 triangulation.levels[level + 1]->reserve_space(
5409 used_cells + new_cells, 3, spacedim);
5412 triangulation.levels[level + 1]->cells.reserve_space(new_cells);
5417 triangulation.begin_quad();
5418 quad != triangulation.end_quad();
5421 if (quad->user_flag_set())
5427 needed_quads_pair += 4;
5428 needed_lines_pair += 4;
5429 needed_vertices += 1;
5431 if (quad->user_index())
5435 needed_quads_pair += 2;
5436 needed_lines_single += 1;
5449 if (quad->has_children())
5451 Assert(quad->refinement_case() ==
5454 if ((face_refinement_cases[quad->user_index()] ==
5456 (quad->child(0)->line_index(1) + 1 !=
5457 quad->child(2)->line_index(1))) ||
5458 (face_refinement_cases[quad->user_index()] ==
5460 (quad->child(0)->line_index(3) + 1 !=
5461 quad->child(1)->line_index(3))))
5462 needed_lines_pair += 2;
5468 triangulation.begin_line();
5469 line != triangulation.end_line();
5471 if (line->user_flag_set())
5473 needed_lines_pair += 2;
5474 needed_vertices += 1;
5478 triangulation.faces->lines.reserve_space(needed_lines_pair,
5479 needed_lines_single);
5481 triangulation.faces->quads.reserve_space(needed_quads_pair,
5482 needed_quads_single);
5486 needed_vertices += std::count_if(triangulation.vertices_used.begin(),
5487 triangulation.vertices_used.end(),
5488 std::bind(std::equal_to<bool>(),
5489 std::placeholders::_1,
5494 if (needed_vertices > triangulation.vertices.size())
5497 triangulation.vertices_used.resize(needed_vertices,
false);
5513 triangulation.begin_active();
5514 cell != triangulation.end();
5516 if (!cell->refine_flag_set())
5517 for (
unsigned int line = 0;
5518 line < GeometryInfo<dim>::lines_per_cell;
5520 if (cell->line(line)->has_children())
5521 for (
unsigned int c = 0; c < 2; ++c)
5522 Assert(cell->line(line)->child(c)->user_flag_set() ==
false,
5534 unsigned int next_unused_vertex = 0;
5540 line = triangulation.begin_active_line(),
5541 endl = triangulation.end_line();
5542 typename Triangulation<dim, spacedim>::raw_line_iterator
5543 next_unused_line = triangulation.begin_raw_line();
5545 for (; line != endl; ++line)
5546 if (line->user_flag_set())
5552 while (triangulation.vertices_used[next_unused_vertex] ==
true)
5553 ++next_unused_vertex;
5555 next_unused_vertex < triangulation.vertices.size(),
5557 "Internal error: During refinement, the triangulation wants to access an element of the 'vertices' array but it turns out that the array is not large enough."));
5558 triangulation.vertices_used[next_unused_vertex] =
true;
5560 triangulation.vertices[next_unused_vertex] = line->center(
true);
5566 triangulation.faces->lines.next_free_pair_object(
5574 line->set_children(0, next_unused_line->index());
5577 const typename Triangulation<dim, spacedim>::raw_line_iterator
5578 children[2] = {next_unused_line, ++next_unused_line};
5583 children[0]->used() ==
false,
5585 "Internal error: We want to use a cell during refinement that should be unused, but turns out not to be."));
5587 children[1]->used() ==
false,
5589 "Internal error: We want to use a cell during refinement that should be unused, but turns out not to be."));
5592 internal::TriangulationImplementation ::TriaObject<1>(
5593 line->vertex_index(0), next_unused_vertex));
5595 internal::TriangulationImplementation ::TriaObject<1>(
5596 next_unused_vertex, line->vertex_index(1)));
5598 children[0]->set_used_flag();
5599 children[1]->set_used_flag();
5600 children[0]->clear_children();
5601 children[1]->clear_children();
5604 children[0]->clear_user_flag();
5605 children[1]->clear_user_flag();
5607 children[0]->set_boundary_id_internal(line->boundary_id());
5608 children[1]->set_boundary_id_internal(line->boundary_id());
5610 children[0]->set_manifold_id(line->manifold_id());
5611 children[1]->set_manifold_id(line->manifold_id());
5616 line->clear_user_flag();
5652 quad = triangulation.begin_quad(),
5653 endq = triangulation.end_quad();
5654 typename Triangulation<dim, spacedim>::raw_line_iterator
5655 next_unused_line = triangulation.begin_raw_line();
5656 typename Triangulation<dim, spacedim>::raw_quad_iterator
5657 next_unused_quad = triangulation.begin_raw_quad();
5659 for (; quad != endq; ++quad)
5661 if (quad->user_index())
5664 face_refinement_cases[quad->user_index()];
5673 if (aniso_quad_ref_case == quad->refinement_case())
5676 Assert(quad->refinement_case() ==
5678 quad->refinement_case() ==
5683 Assert(quad->user_index() ==
5685 quad->user_index() ==
5690 typename Triangulation<dim, spacedim>::raw_line_iterator
5694 triangulation.faces->lines.next_free_single_object(
5697 new_line->used() ==
false,
5699 "Internal error: We want to use a cell during refinement that should be unused, but turns out not to be."));
5714 unsigned int vertex_indices[2];
5718 quad->line(2)->child(0)->vertex_index(1);
5720 quad->line(3)->child(0)->vertex_index(1);
5725 quad->line(0)->child(0)->vertex_index(1);
5727 quad->line(1)->child(0)->vertex_index(1);
5732 vertex_indices[0], vertex_indices[1]));
5733 new_line->set_used_flag();
5734 new_line->clear_user_flag();
5736 new_line->clear_children();
5737 new_line->set_boundary_id_internal(quad->boundary_id());
5738 new_line->set_manifold_id(quad->manifold_id());
5746 const unsigned int index[2][2] = {
5752 typename Triangulation<dim, spacedim>::raw_quad_iterator
5756 triangulation.faces->quads.next_free_pair_object(
5758 new_quads[0] = next_unused_quad;
5760 new_quads[0]->used() ==
false,
5762 "Internal error: We want to use a cell during refinement that should be unused, but turns out not to be."));
5765 new_quads[1] = next_unused_quad;
5767 new_quads[1]->used() ==
false,
5769 "Internal error: We want to use a cell during refinement that should be unused, but turns out not to be."));
5775 internal::TriangulationImplementation ::TriaObject<2>(
5776 quad->line_index(0),
5779 ->child(index[0][quad->line_orientation(2)])
5782 ->child(index[0][quad->line_orientation(3)])
5785 internal::TriangulationImplementation ::TriaObject<2>(
5787 quad->line_index(1),
5789 ->child(index[1][quad->line_orientation(2)])
5792 ->child(index[1][quad->line_orientation(3)])
5798 internal::TriangulationImplementation ::TriaObject<2>(
5800 ->child(index[0][quad->line_orientation(0)])
5803 ->child(index[0][quad->line_orientation(1)])
5805 quad->line_index(2),
5806 new_line->index()));
5808 internal::TriangulationImplementation ::TriaObject<2>(
5810 ->child(index[1][quad->line_orientation(0)])
5813 ->child(index[1][quad->line_orientation(1)])
5816 quad->line_index(3)));
5819 for (
const auto &new_quad : new_quads)
5821 new_quad->set_used_flag();
5822 new_quad->clear_user_flag();
5824 new_quad->clear_children();
5825 new_quad->set_boundary_id_internal(quad->boundary_id());
5826 new_quad->set_manifold_id(quad->manifold_id());
5830 for (
unsigned int j = 0;
5831 j < GeometryInfo<dim>::lines_per_face;
5833 new_quad->set_line_orientation(j,
true);
5839 new_quads[0]->set_line_orientation(
5840 0, quad->line_orientation(0));
5841 new_quads[0]->set_line_orientation(
5842 2, quad->line_orientation(2));
5843 new_quads[1]->set_line_orientation(
5844 1, quad->line_orientation(1));
5845 new_quads[1]->set_line_orientation(
5846 3, quad->line_orientation(3));
5849 new_quads[0]->set_line_orientation(
5850 3, quad->line_orientation(3));
5851 new_quads[1]->set_line_orientation(
5852 2, quad->line_orientation(2));
5856 new_quads[0]->set_line_orientation(
5857 1, quad->line_orientation(1));
5858 new_quads[1]->set_line_orientation(
5859 0, quad->line_orientation(0));
5865 if (quad->refinement_case() ==
5888 if (aniso_quad_ref_case ==
5891 old_child[0] = quad->child(0)->line(1);
5892 old_child[1] = quad->child(2)->line(1);
5896 Assert(aniso_quad_ref_case ==
5900 old_child[0] = quad->child(0)->line(3);
5901 old_child[1] = quad->child(1)->line(3);
5904 if (old_child[0]->index() + 1 != old_child[1]->index())
5910 spacedim>::raw_line_iterator
5913 new_child[0] = new_child[1] =
5914 triangulation.faces->lines.next_free_pair_object(
5918 new_child[0]->set_used_flag();
5919 new_child[1]->set_used_flag();
5921 const int old_index_0 = old_child[0]->index(),
5922 old_index_1 = old_child[1]->index(),
5923 new_index_0 = new_child[0]->index(),
5924 new_index_1 = new_child[1]->index();
5928 for (
unsigned int q = 0;
5929 q < triangulation.faces->quads.cells.size();
5931 for (
unsigned int l = 0;
5932 l < GeometryInfo<dim>::lines_per_face;
5935 const int this_index =
5936 triangulation.faces->quads.cells[q].face(l);
5937 if (this_index == old_index_0)
5938 triangulation.faces->quads.cells[q]
5939 .set_face(l, new_index_0);
5940 else if (this_index == old_index_1)
5941 triangulation.faces->quads.cells[q]
5942 .set_face(l, new_index_1);
5946 for (
unsigned int i = 0; i < 2; ++i)
5948 Assert(!old_child[i]->has_children(),
5954 old_child[i]->vertex_index(
5956 new_child[i]->set_boundary_id_internal(
5957 old_child[i]->boundary_id());
5958 new_child[i]->set_manifold_id(
5959 old_child[i]->manifold_id());
5960 new_child[i]->set_user_index(
5961 old_child[i]->user_index());
5962 if (old_child[i]->user_flag_set())
5963 new_child[i]->set_user_flag();
5965 new_child[i]->clear_user_flag();
5967 new_child[i]->clear_children();
5969 old_child[i]->clear_user_flag();
5970 old_child[i]->clear_user_index();
5971 old_child[i]->clear_used_flag();
5977 if (aniso_quad_ref_case ==
5980 new_line->set_children(
5981 0, quad->child(0)->line_index(1));
5982 Assert(new_line->child(1) ==
5983 quad->child(2)->line(1),
6020 quad_iterator switch_1 = quad->child(1),
6021 switch_2 = quad->child(2);
6022 const int switch_1_index = switch_1->index();
6023 const int switch_2_index = switch_2->index();
6024 for (
unsigned int l = 0;
6025 l < triangulation.levels.size();
6027 for (
unsigned int h = 0;
6029 triangulation.levels[l]->cells.cells.size();
6031 for (
unsigned int q = 0;
6032 q < GeometryInfo<dim>::faces_per_cell;
6035 const int face_index =
6036 triangulation.levels[l]
6039 if (face_index == switch_1_index)
6040 triangulation.levels[l]
6042 .set_face(q, switch_2_index);
6043 else if (face_index == switch_2_index)
6044 triangulation.levels[l]
6046 .set_face(q, switch_1_index);
6050 const unsigned int switch_1_lines[4] = {
6051 switch_1->line_index(0),
6052 switch_1->line_index(1),
6053 switch_1->line_index(2),
6054 switch_1->line_index(3)};
6055 const bool switch_1_line_orientations[4] = {
6056 switch_1->line_orientation(0),
6057 switch_1->line_orientation(1),
6058 switch_1->line_orientation(2),
6059 switch_1->line_orientation(3)};
6061 switch_1->boundary_id();
6062 const unsigned int switch_1_user_index =
6063 switch_1->user_index();
6064 const bool switch_1_user_flag =
6065 switch_1->user_flag_set();
6067 switch_1_refinement_case =
6068 switch_1->refinement_case();
6069 const int switch_1_first_child_pair =
6070 (switch_1_refinement_case ?
6071 switch_1->child_index(0) :
6073 const int switch_1_second_child_pair =
6074 (switch_1_refinement_case ==
6076 switch_1->child_index(2) :
6081 2>(switch_2->line_index(0),
6082 switch_2->line_index(1),
6083 switch_2->line_index(2),
6084 switch_2->line_index(3)));
6085 switch_1->set_line_orientation(
6086 0, switch_2->line_orientation(0));
6087 switch_1->set_line_orientation(
6088 1, switch_2->line_orientation(1));
6089 switch_1->set_line_orientation(
6090 2, switch_2->line_orientation(2));
6091 switch_1->set_line_orientation(
6092 3, switch_2->line_orientation(3));
6093 switch_1->set_boundary_id_internal(
6094 switch_2->boundary_id());
6095 switch_1->set_manifold_id(switch_2->manifold_id());
6096 switch_1->set_user_index(switch_2->user_index());
6097 if (switch_2->user_flag_set())
6098 switch_1->set_user_flag();
6100 switch_1->clear_user_flag();
6101 switch_1->clear_refinement_case();
6102 switch_1->set_refinement_case(
6103 switch_2->refinement_case());
6104 switch_1->clear_children();
6105 if (switch_2->refinement_case())
6106 switch_1->set_children(0,
6107 switch_2->child_index(0));
6108 if (switch_2->refinement_case() ==
6110 switch_1->set_children(2,
6111 switch_2->child_index(2));
6115 2>(switch_1_lines[0],
6118 switch_1_lines[3]));
6119 switch_2->set_line_orientation(
6120 0, switch_1_line_orientations[0]);
6121 switch_2->set_line_orientation(
6122 1, switch_1_line_orientations[1]);
6123 switch_2->set_line_orientation(
6124 2, switch_1_line_orientations[2]);
6125 switch_2->set_line_orientation(
6126 3, switch_1_line_orientations[3]);
6127 switch_2->set_boundary_id_internal(
6128 switch_1_boundary_id);
6129 switch_2->set_manifold_id(switch_1->manifold_id());
6130 switch_2->set_user_index(switch_1_user_index);
6131 if (switch_1_user_flag)
6132 switch_2->set_user_flag();
6134 switch_2->clear_user_flag();
6135 switch_2->clear_refinement_case();
6136 switch_2->set_refinement_case(
6137 switch_1_refinement_case);
6138 switch_2->clear_children();
6139 switch_2->set_children(0,
6140 switch_1_first_child_pair);
6141 switch_2->set_children(2,
6142 switch_1_second_child_pair);
6144 new_quads[0]->set_refinement_case(
6146 new_quads[0]->set_children(0, quad->child_index(0));
6147 new_quads[1]->set_refinement_case(
6149 new_quads[1]->set_children(0, quad->child_index(2));
6153 new_quads[0]->set_refinement_case(
6155 new_quads[0]->set_children(0, quad->child_index(0));
6156 new_quads[1]->set_refinement_case(
6158 new_quads[1]->set_children(0, quad->child_index(2));
6159 new_line->set_children(
6160 0, quad->child(0)->line_index(3));
6161 Assert(new_line->child(1) ==
6162 quad->child(1)->line(3),
6165 quad->clear_children();
6169 quad->set_children(0, new_quads[0]->index());
6171 quad->set_refinement_case(aniso_quad_ref_case);
6178 if (quad->user_flag_set())
6188 while (triangulation.vertices_used[next_unused_vertex] ==
6190 ++next_unused_vertex;
6192 next_unused_vertex < triangulation.vertices.size(),
6194 "Internal error: During refinement, the triangulation wants to access an element of the 'vertices' array but it turns out that the array is not large enough."));
6202 quad->refinement_case();
6208 quad->child(0)->set_user_index(
6210 quad->child(1)->set_user_index(
6216 middle_line = quad->child(0)->line(1);
6218 middle_line = quad->child(0)->line(3);
6225 if (!middle_line->has_children())
6232 triangulation.vertices[next_unused_vertex] =
6233 middle_line->center(
true);
6234 triangulation.vertices_used[next_unused_vertex] =
6240 triangulation.faces->lines.next_free_pair_object(
6245 middle_line->set_children(
6246 0, next_unused_line->index());
6250 raw_line_iterator children[2] = {
6251 next_unused_line, ++next_unused_line};
6257 children[0]->used() ==
false,
6259 "Internal error: We want to use a cell during refinement that should be unused, but turns out not to be."));
6261 children[1]->used() ==
false,
6263 "Internal error: We want to use a cell during refinement that should be unused, but turns out not to be."));
6267 1>(middle_line->vertex_index(0),
6268 next_unused_vertex));
6271 1>(next_unused_vertex,
6272 middle_line->vertex_index(1)));
6274 children[0]->set_used_flag();
6275 children[1]->set_used_flag();
6276 children[0]->clear_children();
6277 children[1]->clear_children();
6280 children[0]->clear_user_flag();
6281 children[1]->clear_user_flag();
6283 children[0]->set_boundary_id_internal(
6284 middle_line->boundary_id());
6285 children[1]->set_boundary_id_internal(
6286 middle_line->boundary_id());
6288 children[0]->set_manifold_id(
6289 middle_line->manifold_id());
6290 children[1]->set_manifold_id(
6291 middle_line->manifold_id());
6297 quad->clear_user_flag();
6326 triangulation.vertices[next_unused_vertex] =
6327 quad->center(
true,
true);
6328 triangulation.vertices_used[next_unused_vertex] =
true;
6333 typename Triangulation<dim, spacedim>::raw_line_iterator
6336 for (
unsigned int i = 0; i < 4; ++i)
6346 triangulation.faces->lines.next_free_pair_object(
6349 new_lines[i] = next_unused_line;
6353 new_lines[i]->used() ==
false,
6355 "Internal error: We want to use a cell during refinement that should be unused, but turns out not to be."));
6375 const unsigned int vertex_indices[5] = {
6376 quad->line(0)->child(0)->vertex_index(1),
6377 quad->line(1)->child(0)->vertex_index(1),
6378 quad->line(2)->child(0)->vertex_index(1),
6379 quad->line(3)->child(0)->vertex_index(1),
6380 next_unused_vertex};
6384 vertex_indices[2], vertex_indices[4]));
6387 vertex_indices[4], vertex_indices[3]));
6390 vertex_indices[0], vertex_indices[4]));
6393 vertex_indices[4], vertex_indices[1]));
6395 for (
const auto &new_line : new_lines)
6397 new_line->set_used_flag();
6398 new_line->clear_user_flag();
6400 new_line->clear_children();
6401 new_line->set_boundary_id_internal(quad->boundary_id());
6402 new_line->set_manifold_id(quad->manifold_id());
6421 const unsigned int index[2][2] = {
6425 const int line_indices[12] = {
6427 ->child(index[0][quad->line_orientation(0)])
6430 ->child(index[1][quad->line_orientation(0)])
6433 ->child(index[0][quad->line_orientation(1)])
6436 ->child(index[1][quad->line_orientation(1)])
6439 ->child(index[0][quad->line_orientation(2)])
6442 ->child(index[1][quad->line_orientation(2)])
6445 ->child(index[0][quad->line_orientation(3)])
6448 ->child(index[1][quad->line_orientation(3)])
6450 new_lines[0]->index(),
6451 new_lines[1]->index(),
6452 new_lines[2]->index(),
6453 new_lines[3]->index()};
6458 typename Triangulation<dim, spacedim>::raw_quad_iterator
6462 triangulation.faces->quads.next_free_pair_object(
6465 new_quads[0] = next_unused_quad;
6467 new_quads[0]->used() ==
false,
6469 "Internal error: We want to use a cell during refinement that should be unused, but turns out not to be."));
6472 new_quads[1] = next_unused_quad;
6474 new_quads[1]->used() ==
false,
6476 "Internal error: We want to use a cell during refinement that should be unused, but turns out not to be."));
6479 triangulation.faces->quads.next_free_pair_object(
6481 new_quads[2] = next_unused_quad;
6483 new_quads[2]->used() ==
false,
6485 "Internal error: We want to use a cell during refinement that should be unused, but turns out not to be."));
6488 new_quads[3] = next_unused_quad;
6490 new_quads[3]->used() ==
false,
6492 "Internal error: We want to use a cell during refinement that should be unused, but turns out not to be."));
6495 quad->set_children(0, new_quads[0]->index());
6496 quad->set_children(2, new_quads[2]->index());
6498 internal::TriangulationImplementation ::TriaObject<2>(
6507 internal::TriangulationImplementation ::TriaObject<2>(
6513 internal::TriangulationImplementation ::TriaObject<2>(
6519 internal::TriangulationImplementation ::TriaObject<2>(
6525 internal::TriangulationImplementation ::TriaObject<2>(
6530 for (
const auto &new_quad : new_quads)
6532 new_quad->set_used_flag();
6533 new_quad->clear_user_flag();
6535 new_quad->clear_children();
6536 new_quad->set_boundary_id_internal(quad->boundary_id());
6537 new_quad->set_manifold_id(quad->manifold_id());
6541 for (
unsigned int j = 0;
6542 j < GeometryInfo<dim>::lines_per_face;
6544 new_quad->set_line_orientation(j,
true);
6550 new_quads[0]->set_line_orientation(
6551 0, quad->line_orientation(0));
6552 new_quads[0]->set_line_orientation(
6553 2, quad->line_orientation(2));
6554 new_quads[1]->set_line_orientation(
6555 1, quad->line_orientation(1));
6556 new_quads[1]->set_line_orientation(
6557 2, quad->line_orientation(2));
6558 new_quads[2]->set_line_orientation(
6559 0, quad->line_orientation(0));
6560 new_quads[2]->set_line_orientation(
6561 3, quad->line_orientation(3));
6562 new_quads[3]->set_line_orientation(
6563 1, quad->line_orientation(1));
6564 new_quads[3]->set_line_orientation(
6565 3, quad->line_orientation(3));
6569 quad->clear_user_flag();
6580 cells_with_distorted_children;
6582 for (
unsigned int level = 0; level != triangulation.levels.size() - 1;
6589 hex = triangulation.begin_active_hex(level),
6590 endh = triangulation.begin_active_hex(level + 1);
6591 typename Triangulation<dim, spacedim>::raw_hex_iterator
6592 next_unused_hex = triangulation.begin_raw_hex(level + 1);
6594 for (; hex != endh; ++hex)
6595 if (hex->refine_flag_set())
6603 hex->clear_refine_flag();
6604 hex->set_refinement_case(ref_case);
6615 unsigned int n_new_lines = 0;
6616 unsigned int n_new_quads = 0;
6617 unsigned int n_new_hexes = 0;
6647 typename Triangulation<dim, spacedim>::raw_line_iterator>
6648 new_lines(n_new_lines);
6649 for (
unsigned int i = 0; i < n_new_lines; ++i)
6652 triangulation.faces->lines.next_free_single_object(
6656 new_lines[i]->used() ==
false,
6658 "Internal error: We want to use a cell during refinement that should be unused, but turns out not to be."));
6659 new_lines[i]->set_used_flag();
6660 new_lines[i]->clear_user_flag();
6661 new_lines[i]->clear_user_data();
6662 new_lines[i]->clear_children();
6664 new_lines[i]->set_boundary_id_internal(
6668 new_lines[i]->set_manifold_id(hex->manifold_id());
6674 typename Triangulation<dim, spacedim>::raw_quad_iterator>
6675 new_quads(n_new_quads);
6676 for (
unsigned int i = 0; i < n_new_quads; ++i)
6679 triangulation.faces->quads.next_free_single_object(
6683 new_quads[i]->used() ==
false,
6685 "Internal error: We want to use a cell during refinement that should be unused, but turns out not to be."));
6686 new_quads[i]->set_used_flag();
6687 new_quads[i]->clear_user_flag();
6688 new_quads[i]->clear_user_data();
6689 new_quads[i]->clear_children();
6691 new_quads[i]->set_boundary_id_internal(
6695 new_quads[i]->set_manifold_id(hex->manifold_id());
6698 for (
unsigned int j = 0;
6699 j < GeometryInfo<dim>::lines_per_face;
6701 new_quads[i]->set_line_orientation(j,
true);
6709 typename Triangulation<dim, spacedim>::raw_hex_iterator>
6710 new_hexes(n_new_hexes);
6711 for (
unsigned int i = 0; i < n_new_hexes; ++i)
6715 triangulation.levels[level + 1]->cells.next_free_hex(
6716 triangulation, level + 1);
6720 new_hexes[i] = next_unused_hex;
6723 new_hexes[i]->used() ==
false,
6725 "Internal error: We want to use a cell during refinement that should be unused, but turns out not to be."));
6726 new_hexes[i]->set_used_flag();
6727 new_hexes[i]->clear_user_flag();
6728 new_hexes[i]->clear_user_data();
6729 new_hexes[i]->clear_children();
6732 new_hexes[i]->set_material_id(hex->material_id());
6733 new_hexes[i]->set_manifold_id(hex->manifold_id());
6734 new_hexes[i]->set_subdomain_id(subdomainid);
6737 new_hexes[i]->set_parent(hex->index());
6749 for (
unsigned int f = 0;
6750 f < GeometryInfo<dim>::faces_per_cell;
6753 new_hexes[i]->set_face_orientation(f,
true);
6754 new_hexes[i]->set_face_flip(f,
false);
6755 new_hexes[i]->set_face_rotation(f,
false);
6759 for (
unsigned int i = 0; i < n_new_hexes / 2; ++i)
6760 hex->set_children(2 * i, new_hexes[2 * i]->index());
6767 const bool f_or[6] = {hex->face_orientation(0),
6768 hex->face_orientation(1),
6769 hex->face_orientation(2),
6770 hex->face_orientation(3),
6771 hex->face_orientation(4),
6772 hex->face_orientation(5)};
6775 const bool f_fl[6] = {hex->face_flip(0),
6783 const bool f_ro[6] = {hex->face_rotation(0),
6784 hex->face_rotation(1),
6785 hex->face_rotation(2),
6786 hex->face_rotation(3),
6787 hex->face_rotation(4),
6788 hex->face_rotation(5)};
6798 const unsigned int child_at_origin[2][2][2] = {
6885 raw_line_iterator lines[4] = {
6886 hex->face(2)->child(0)->line(
6887 (hex->face(2)->refinement_case() ==
6891 hex->face(3)->child(0)->line(
6892 (hex->face(3)->refinement_case() ==
6896 hex->face(4)->child(0)->line(
6897 (hex->face(4)->refinement_case() ==
6901 hex->face(5)->child(0)->line(
6902 (hex->face(5)->refinement_case() ==
6908 unsigned int line_indices[4];
6909 for (
unsigned int i = 0; i < 4; ++i)
6910 line_indices[i] = lines[i]->index();
6917 bool line_orientation[4];
6923 const unsigned int middle_vertices[2] = {
6924 hex->line(2)->child(0)->vertex_index(1),
6925 hex->line(7)->child(0)->vertex_index(1)};
6927 for (
unsigned int i = 0; i < 4; ++i)
6928 if (lines[i]->vertex_index(i % 2) ==
6929 middle_vertices[i % 2])
6930 line_orientation[i] =
true;
6935 Assert(lines[i]->vertex_index((i + 1) % 2) ==
6936 middle_vertices[i % 2],
6938 line_orientation[i] =
false;
6944 internal::TriangulationImplementation ::TriaObject<
6950 new_quads[0]->set_line_orientation(
6951 0, line_orientation[0]);
6952 new_quads[0]->set_line_orientation(
6953 1, line_orientation[1]);
6954 new_quads[0]->set_line_orientation(
6955 2, line_orientation[2]);
6956 new_quads[0]->set_line_orientation(
6957 3, line_orientation[3]);
6990 const int quad_indices[11] = {
6991 new_quads[0]->index(),
6993 hex->face(0)->index(),
6995 hex->face(1)->index(),
6997 hex->face(2)->child_index(
6998 child_at_origin[hex->face(2)->refinement_case() -
6999 1][f_fl[2]][f_ro[2]]),
7000 hex->face(2)->child_index(
7002 child_at_origin[hex->face(2)->refinement_case() -
7003 1][f_fl[2]][f_ro[2]]),
7005 hex->face(3)->child_index(
7006 child_at_origin[hex->face(3)->refinement_case() -
7007 1][f_fl[3]][f_ro[3]]),
7008 hex->face(3)->child_index(
7010 child_at_origin[hex->face(3)->refinement_case() -
7011 1][f_fl[3]][f_ro[3]]),
7013 hex->face(4)->child_index(
7014 child_at_origin[hex->face(4)->refinement_case() -
7015 1][f_fl[4]][f_ro[4]]),
7016 hex->face(4)->child_index(
7018 child_at_origin[hex->face(4)->refinement_case() -
7019 1][f_fl[4]][f_ro[4]]),
7021 hex->face(5)->child_index(
7022 child_at_origin[hex->face(5)->refinement_case() -
7023 1][f_fl[5]][f_ro[5]]),
7024 hex->face(5)->child_index(
7026 child_at_origin[hex->face(5)->refinement_case() -
7027 1][f_fl[5]][f_ro[5]])
7032 internal::TriangulationImplementation ::TriaObject<
7040 internal::TriangulationImplementation ::TriaObject<
7115 raw_line_iterator lines[4] = {
7116 hex->face(0)->child(0)->line(
7117 (hex->face(0)->refinement_case() ==
7121 hex->face(1)->child(0)->line(
7122 (hex->face(1)->refinement_case() ==
7126 hex->face(4)->child(0)->line(
7127 (hex->face(4)->refinement_case() ==
7131 hex->face(5)->child(0)->line(
7132 (hex->face(5)->refinement_case() ==
7138 unsigned int line_indices[4];
7139 for (
unsigned int i = 0; i < 4; ++i)
7140 line_indices[i] = lines[i]->index();
7147 bool line_orientation[4];
7153 const unsigned int middle_vertices[2] = {
7154 hex->line(0)->child(0)->vertex_index(1),
7155 hex->line(5)->child(0)->vertex_index(1)};
7157 for (
unsigned int i = 0; i < 4; ++i)
7158 if (lines[i]->vertex_index(i % 2) ==
7159 middle_vertices[i % 2])
7160 line_orientation[i] =
true;
7164 Assert(lines[i]->vertex_index((i + 1) % 2) ==
7165 middle_vertices[i % 2],
7167 line_orientation[i] =
false;
7173 internal::TriangulationImplementation ::TriaObject<
7179 new_quads[0]->set_line_orientation(
7180 0, line_orientation[2]);
7181 new_quads[0]->set_line_orientation(
7182 1, line_orientation[3]);
7183 new_quads[0]->set_line_orientation(
7184 2, line_orientation[0]);
7185 new_quads[0]->set_line_orientation(
7186 3, line_orientation[1]);
7219 const int quad_indices[11] = {
7220 new_quads[0]->index(),
7222 hex->face(0)->child_index(
7223 child_at_origin[hex->face(0)->refinement_case() -
7224 1][f_fl[0]][f_ro[0]]),
7225 hex->face(0)->child_index(
7227 child_at_origin[hex->face(0)->refinement_case() -
7228 1][f_fl[0]][f_ro[0]]),
7230 hex->face(1)->child_index(
7231 child_at_origin[hex->face(1)->refinement_case() -
7232 1][f_fl[1]][f_ro[1]]),
7233 hex->face(1)->child_index(
7235 child_at_origin[hex->face(1)->refinement_case() -
7236 1][f_fl[1]][f_ro[1]]),
7238 hex->face(2)->index(),
7240 hex->face(3)->index(),
7242 hex->face(4)->child_index(
7243 child_at_origin[hex->face(4)->refinement_case() -
7244 1][f_fl[4]][f_ro[4]]),
7245 hex->face(4)->child_index(
7247 child_at_origin[hex->face(4)->refinement_case() -
7248 1][f_fl[4]][f_ro[4]]),
7250 hex->face(5)->child_index(
7251 child_at_origin[hex->face(5)->refinement_case() -
7252 1][f_fl[5]][f_ro[5]]),
7253 hex->face(5)->child_index(
7255 child_at_origin[hex->face(5)->refinement_case() -
7256 1][f_fl[5]][f_ro[5]])
7261 internal::TriangulationImplementation ::TriaObject<
7269 internal::TriangulationImplementation ::TriaObject<
7346 raw_line_iterator lines[4] = {
7347 hex->face(0)->child(0)->line(
7348 (hex->face(0)->refinement_case() ==
7352 hex->face(1)->child(0)->line(
7353 (hex->face(1)->refinement_case() ==
7357 hex->face(2)->child(0)->line(
7358 (hex->face(2)->refinement_case() ==
7362 hex->face(3)->child(0)->line(
7363 (hex->face(3)->refinement_case() ==
7369 unsigned int line_indices[4];
7370 for (
unsigned int i = 0; i < 4; ++i)
7371 line_indices[i] = lines[i]->index();
7378 bool line_orientation[4];
7384 const unsigned int middle_vertices[2] = {
7385 middle_vertex_index<dim, spacedim>(hex->line(8)),
7386 middle_vertex_index<dim, spacedim>(hex->line(11))};
7388 for (
unsigned int i = 0; i < 4; ++i)
7389 if (lines[i]->vertex_index(i % 2) ==
7390 middle_vertices[i % 2])
7391 line_orientation[i] =
true;
7395 Assert(lines[i]->vertex_index((i + 1) % 2) ==
7396 middle_vertices[i % 2],
7398 line_orientation[i] =
false;
7404 internal::TriangulationImplementation ::TriaObject<
7410 new_quads[0]->set_line_orientation(
7411 0, line_orientation[0]);
7412 new_quads[0]->set_line_orientation(
7413 1, line_orientation[1]);
7414 new_quads[0]->set_line_orientation(
7415 2, line_orientation[2]);
7416 new_quads[0]->set_line_orientation(
7417 3, line_orientation[3]);
7451 const int quad_indices[11] = {
7452 new_quads[0]->index(),
7454 hex->face(0)->child_index(
7455 child_at_origin[hex->face(0)->refinement_case() -
7456 1][f_fl[0]][f_ro[0]]),
7457 hex->face(0)->child_index(
7459 child_at_origin[hex->face(0)->refinement_case() -
7460 1][f_fl[0]][f_ro[0]]),
7462 hex->face(1)->child_index(
7463 child_at_origin[hex->face(1)->refinement_case() -
7464 1][f_fl[1]][f_ro[1]]),
7465 hex->face(1)->child_index(
7467 child_at_origin[hex->face(1)->refinement_case() -
7468 1][f_fl[1]][f_ro[1]]),
7470 hex->face(2)->child_index(
7471 child_at_origin[hex->face(2)->refinement_case() -
7472 1][f_fl[2]][f_ro[2]]),
7473 hex->face(2)->child_index(
7475 child_at_origin[hex->face(2)->refinement_case() -
7476 1][f_fl[2]][f_ro[2]]),
7478 hex->face(3)->child_index(
7479 child_at_origin[hex->face(3)->refinement_case() -
7480 1][f_fl[3]][f_ro[3]]),
7481 hex->face(3)->child_index(
7483 child_at_origin[hex->face(3)->refinement_case() -
7484 1][f_fl[3]][f_ro[3]]),
7486 hex->face(4)->index(),
7488 hex->face(5)->index()
7492 internal::TriangulationImplementation ::TriaObject<
7500 internal::TriangulationImplementation ::TriaObject<
7534 1>(middle_vertex_index<dim, spacedim>(
7536 middle_vertex_index<dim, spacedim>(
7605 spacedim>::raw_line_iterator lines[13] = {
7606 hex->face(0)->child(0)->line(
7607 (hex->face(0)->refinement_case() ==
7611 hex->face(1)->child(0)->line(
7612 (hex->face(1)->refinement_case() ==
7616 hex->face(2)->child(0)->line(
7617 (hex->face(2)->refinement_case() ==
7621 hex->face(3)->child(0)->line(
7622 (hex->face(3)->refinement_case() ==
7630 0, f_or[4], f_fl[4], f_ro[4]))
7633 1, f_or[4], f_fl[4], f_ro[4])),
7637 3, f_or[4], f_fl[4], f_ro[4]))
7640 0, f_or[4], f_fl[4], f_ro[4])),
7644 0, f_or[4], f_fl[4], f_ro[4]))
7647 3, f_or[4], f_fl[4], f_ro[4])),
7651 3, f_or[4], f_fl[4], f_ro[4]))
7654 2, f_or[4], f_fl[4], f_ro[4])),
7659 0, f_or[5], f_fl[5], f_ro[5]))
7662 1, f_or[5], f_fl[5], f_ro[5])),
7666 3, f_or[5], f_fl[5], f_ro[5]))
7669 0, f_or[5], f_fl[5], f_ro[5])),
7673 0, f_or[5], f_fl[5], f_ro[5]))
7676 3, f_or[5], f_fl[5], f_ro[5])),
7680 3, f_or[5], f_fl[5], f_ro[5]))
7683 2, f_or[5], f_fl[5], f_ro[5])),
7688 unsigned int line_indices[13];
7689 for (
unsigned int i = 0; i < 13; ++i)
7690 line_indices[i] = lines[i]->index();
7697 bool line_orientation[13];
7701 const unsigned int middle_vertices[4] = {
7702 hex->line(0)->child(0)->vertex_index(1),
7703 hex->line(1)->child(0)->vertex_index(1),
7704 hex->line(2)->child(0)->vertex_index(1),
7705 hex->line(3)->child(0)->vertex_index(1),
7711 for (
unsigned int i = 0; i < 4; ++i)
7712 if (lines[i]->vertex_index(0) == middle_vertices[i])
7713 line_orientation[i] =
true;
7717 Assert(lines[i]->vertex_index(1) ==
7720 line_orientation[i] =
false;
7729 for (
unsigned int i = 4; i < 12; ++i)
7730 if (lines[i]->vertex_index((i + 1) % 2) ==
7731 middle_vertex_index<dim, spacedim>(
7732 hex->face(3 + i / 4)))
7733 line_orientation[i] =
true;
7738 Assert(lines[i]->vertex_index(i % 2) ==
7739 (middle_vertex_index<dim, spacedim>(
7740 hex->face(3 + i / 4))),
7742 line_orientation[i] =
false;
7747 line_orientation[12] =
true;
7773 internal::TriangulationImplementation ::TriaObject<
7779 internal::TriangulationImplementation ::TriaObject<
7780 2>(line_indices[12],
7785 internal::TriangulationImplementation ::TriaObject<
7791 internal::TriangulationImplementation ::TriaObject<
7797 new_quads[0]->set_line_orientation(
7798 0, line_orientation[2]);
7799 new_quads[0]->set_line_orientation(
7800 2, line_orientation[4]);
7801 new_quads[0]->set_line_orientation(
7802 3, line_orientation[8]);
7804 new_quads[1]->set_line_orientation(
7805 1, line_orientation[3]);
7806 new_quads[1]->set_line_orientation(
7807 2, line_orientation[5]);
7808 new_quads[1]->set_line_orientation(
7809 3, line_orientation[9]);
7811 new_quads[2]->set_line_orientation(
7812 0, line_orientation[6]);
7813 new_quads[2]->set_line_orientation(
7814 1, line_orientation[10]);
7815 new_quads[2]->set_line_orientation(
7816 2, line_orientation[0]);
7818 new_quads[3]->set_line_orientation(
7819 0, line_orientation[7]);
7820 new_quads[3]->set_line_orientation(
7821 1, line_orientation[11]);
7822 new_quads[3]->set_line_orientation(
7823 3, line_orientation[1]);
7856 const int quad_indices[20] = {
7857 new_quads[0]->index(),
7858 new_quads[1]->index(),
7859 new_quads[2]->index(),
7860 new_quads[3]->index(),
7862 hex->face(0)->child_index(
7863 child_at_origin[hex->face(0)->refinement_case() -
7864 1][f_fl[0]][f_ro[0]]),
7865 hex->face(0)->child_index(
7867 child_at_origin[hex->face(0)->refinement_case() -
7868 1][f_fl[0]][f_ro[0]]),
7870 hex->face(1)->child_index(
7871 child_at_origin[hex->face(1)->refinement_case() -
7872 1][f_fl[1]][f_ro[1]]),
7873 hex->face(1)->child_index(
7875 child_at_origin[hex->face(1)->refinement_case() -
7876 1][f_fl[1]][f_ro[1]]),
7878 hex->face(2)->child_index(
7879 child_at_origin[hex->face(2)->refinement_case() -
7880 1][f_fl[2]][f_ro[2]]),
7881 hex->face(2)->child_index(
7883 child_at_origin[hex->face(2)->refinement_case() -
7884 1][f_fl[2]][f_ro[2]]),
7886 hex->face(3)->child_index(
7887 child_at_origin[hex->face(3)->refinement_case() -
7888 1][f_fl[3]][f_ro[3]]),
7889 hex->face(3)->child_index(
7891 child_at_origin[hex->face(3)->refinement_case() -
7892 1][f_fl[3]][f_ro[3]]),
7894 hex->face(4)->isotropic_child_index(
7896 0, f_or[4], f_fl[4], f_ro[4])),
7897 hex->face(4)->isotropic_child_index(
7899 1, f_or[4], f_fl[4], f_ro[4])),
7900 hex->face(4)->isotropic_child_index(
7902 2, f_or[4], f_fl[4], f_ro[4])),
7903 hex->face(4)->isotropic_child_index(
7905 3, f_or[4], f_fl[4], f_ro[4])),
7907 hex->face(5)->isotropic_child_index(
7909 0, f_or[5], f_fl[5], f_ro[5])),
7910 hex->face(5)->isotropic_child_index(
7912 1, f_or[5], f_fl[5], f_ro[5])),
7913 hex->face(5)->isotropic_child_index(
7915 2, f_or[5], f_fl[5], f_ro[5])),
7916 hex->face(5)->isotropic_child_index(
7918 3, f_or[5], f_fl[5], f_ro[5]))};
7921 internal::TriangulationImplementation ::TriaObject<
7929 internal::TriangulationImplementation ::TriaObject<
7937 internal::TriangulationImplementation ::TriaObject<
7945 internal::TriangulationImplementation ::TriaObject<
7979 1>(middle_vertex_index<dim, spacedim>(
7981 middle_vertex_index<dim, spacedim>(
8050 spacedim>::raw_line_iterator lines[13] = {
8051 hex->face(0)->child(0)->line(
8052 (hex->face(0)->refinement_case() ==
8056 hex->face(1)->child(0)->line(
8057 (hex->face(1)->refinement_case() ==
8061 hex->face(4)->child(0)->line(
8062 (hex->face(4)->refinement_case() ==
8066 hex->face(5)->child(0)->line(
8067 (hex->face(5)->refinement_case() ==
8075 0, f_or[2], f_fl[2], f_ro[2]))
8078 3, f_or[2], f_fl[2], f_ro[2])),
8082 3, f_or[2], f_fl[2], f_ro[2]))
8085 2, f_or[2], f_fl[2], f_ro[2])),
8089 0, f_or[2], f_fl[2], f_ro[2]))
8092 1, f_or[2], f_fl[2], f_ro[2])),
8096 3, f_or[2], f_fl[2], f_ro[2]))
8099 0, f_or[2], f_fl[2], f_ro[2])),
8104 0, f_or[3], f_fl[3], f_ro[3]))
8107 3, f_or[3], f_fl[3], f_ro[3])),
8111 3, f_or[3], f_fl[3], f_ro[3]))
8114 2, f_or[3], f_fl[3], f_ro[3])),
8118 0, f_or[3], f_fl[3], f_ro[3]))
8121 1, f_or[3], f_fl[3], f_ro[3])),
8125 3, f_or[3], f_fl[3], f_ro[3]))
8128 0, f_or[3], f_fl[3], f_ro[3])),
8133 unsigned int line_indices[13];
8134 for (
unsigned int i = 0; i < 13; ++i)
8135 line_indices[i] = lines[i]->index();
8142 bool line_orientation[13];
8146 const unsigned int middle_vertices[4] = {
8147 hex->line(8)->child(0)->vertex_index(1),
8148 hex->line(9)->child(0)->vertex_index(1),
8149 hex->line(2)->child(0)->vertex_index(1),
8150 hex->line(6)->child(0)->vertex_index(1),
8155 for (
unsigned int i = 0; i < 4; ++i)
8156 if (lines[i]->vertex_index(0) == middle_vertices[i])
8157 line_orientation[i] =
true;
8161 Assert(lines[i]->vertex_index(1) ==
8164 line_orientation[i] =
false;
8173 for (
unsigned int i = 4; i < 12; ++i)
8174 if (lines[i]->vertex_index((i + 1) % 2) ==
8175 middle_vertex_index<dim, spacedim>(
8176 hex->face(1 + i / 4)))
8177 line_orientation[i] =
true;
8182 Assert(lines[i]->vertex_index(i % 2) ==
8183 (middle_vertex_index<dim, spacedim>(
8184 hex->face(1 + i / 4))),
8186 line_orientation[i] =
false;
8191 line_orientation[12] =
true;
8218 internal::TriangulationImplementation ::TriaObject<
8224 internal::TriangulationImplementation ::TriaObject<
8225 2>(line_indices[12],
8230 internal::TriangulationImplementation ::TriaObject<
8236 internal::TriangulationImplementation ::TriaObject<
8242 new_quads[0]->set_line_orientation(
8243 0, line_orientation[0]);
8244 new_quads[0]->set_line_orientation(
8245 2, line_orientation[6]);
8246 new_quads[0]->set_line_orientation(
8247 3, line_orientation[10]);
8249 new_quads[1]->set_line_orientation(
8250 1, line_orientation[1]);
8251 new_quads[1]->set_line_orientation(
8252 2, line_orientation[7]);
8253 new_quads[1]->set_line_orientation(
8254 3, line_orientation[11]);
8256 new_quads[2]->set_line_orientation(
8257 0, line_orientation[4]);
8258 new_quads[2]->set_line_orientation(
8259 1, line_orientation[8]);
8260 new_quads[2]->set_line_orientation(
8261 2, line_orientation[2]);
8263 new_quads[3]->set_line_orientation(
8264 0, line_orientation[5]);
8265 new_quads[3]->set_line_orientation(
8266 1, line_orientation[9]);
8267 new_quads[3]->set_line_orientation(
8268 3, line_orientation[3]);
8302 const int quad_indices[20] = {
8303 new_quads[0]->index(),
8304 new_quads[1]->index(),
8305 new_quads[2]->index(),
8306 new_quads[3]->index(),
8308 hex->face(0)->child_index(
8309 child_at_origin[hex->face(0)->refinement_case() -
8310 1][f_fl[0]][f_ro[0]]),
8311 hex->face(0)->child_index(
8313 child_at_origin[hex->face(0)->refinement_case() -
8314 1][f_fl[0]][f_ro[0]]),
8316 hex->face(1)->child_index(
8317 child_at_origin[hex->face(1)->refinement_case() -
8318 1][f_fl[1]][f_ro[1]]),
8319 hex->face(1)->child_index(
8321 child_at_origin[hex->face(1)->refinement_case() -
8322 1][f_fl[1]][f_ro[1]]),
8324 hex->face(2)->isotropic_child_index(
8326 0, f_or[2], f_fl[2], f_ro[2])),
8327 hex->face(2)->isotropic_child_index(
8329 1, f_or[2], f_fl[2], f_ro[2])),
8330 hex->face(2)->isotropic_child_index(
8332 2, f_or[2], f_fl[2], f_ro[2])),
8333 hex->face(2)->isotropic_child_index(
8335 3, f_or[2], f_fl[2], f_ro[2])),
8337 hex->face(3)->isotropic_child_index(
8339 0, f_or[3], f_fl[3], f_ro[3])),
8340 hex->face(3)->isotropic_child_index(
8342 1, f_or[3], f_fl[3], f_ro[3])),
8343 hex->face(3)->isotropic_child_index(
8345 2, f_or[3], f_fl[3], f_ro[3])),
8346 hex->face(3)->isotropic_child_index(
8348 3, f_or[3], f_fl[3], f_ro[3])),
8350 hex->face(4)->child_index(
8351 child_at_origin[hex->face(4)->refinement_case() -
8352 1][f_fl[4]][f_ro[4]]),
8353 hex->face(4)->child_index(
8355 child_at_origin[hex->face(4)->refinement_case() -
8356 1][f_fl[4]][f_ro[4]]),
8358 hex->face(5)->child_index(
8359 child_at_origin[hex->face(5)->refinement_case() -
8360 1][f_fl[5]][f_ro[5]]),
8361 hex->face(5)->child_index(
8363 child_at_origin[hex->face(5)->refinement_case() -
8364 1][f_fl[5]][f_ro[5]])};
8376 internal::TriangulationImplementation ::TriaObject<
8384 internal::TriangulationImplementation ::TriaObject<
8392 internal::TriangulationImplementation ::TriaObject<
8400 internal::TriangulationImplementation ::TriaObject<
8435 1>(middle_vertex_index<dim, spacedim>(
8437 middle_vertex_index<dim, spacedim>(
8507 spacedim>::raw_line_iterator lines[13] = {
8508 hex->face(2)->child(0)->line(
8509 (hex->face(2)->refinement_case() ==
8513 hex->face(3)->child(0)->line(
8514 (hex->face(3)->refinement_case() ==
8518 hex->face(4)->child(0)->line(
8519 (hex->face(4)->refinement_case() ==
8523 hex->face(5)->child(0)->line(
8524 (hex->face(5)->refinement_case() ==
8532 0, f_or[0], f_fl[0], f_ro[0]))
8535 1, f_or[0], f_fl[0], f_ro[0])),
8539 3, f_or[0], f_fl[0], f_ro[0]))
8542 0, f_or[0], f_fl[0], f_ro[0])),
8546 0, f_or[0], f_fl[0], f_ro[0]))
8549 3, f_or[0], f_fl[0], f_ro[0])),
8553 3, f_or[0], f_fl[0], f_ro[0]))
8556 2, f_or[0], f_fl[0], f_ro[0])),
8561 0, f_or[1], f_fl[1], f_ro[1]))
8564 1, f_or[1], f_fl[1], f_ro[1])),
8568 3, f_or[1], f_fl[1], f_ro[1]))
8571 0, f_or[1], f_fl[1], f_ro[1])),
8575 0, f_or[1], f_fl[1], f_ro[1]))
8578 3, f_or[1], f_fl[1], f_ro[1])),
8582 3, f_or[1], f_fl[1], f_ro[1]))
8585 2, f_or[1], f_fl[1], f_ro[1])),
8590 unsigned int line_indices[13];
8592 for (
unsigned int i = 0; i < 13; ++i)
8593 line_indices[i] = lines[i]->index();
8600 bool line_orientation[13];
8604 const unsigned int middle_vertices[4] = {
8605 hex->line(8)->child(0)->vertex_index(1),
8606 hex->line(10)->child(0)->vertex_index(1),
8607 hex->line(0)->child(0)->vertex_index(1),
8608 hex->line(4)->child(0)->vertex_index(1),
8613 for (
unsigned int i = 0; i < 4; ++i)
8614 if (lines[i]->vertex_index(0) == middle_vertices[i])
8615 line_orientation[i] =
true;
8619 Assert(lines[i]->vertex_index(1) ==
8622 line_orientation[i] =
false;
8631 for (
unsigned int i = 4; i < 12; ++i)
8632 if (lines[i]->vertex_index((i + 1) % 2) ==
8633 middle_vertex_index<dim, spacedim>(
8634 hex->face(i / 4 - 1)))
8635 line_orientation[i] =
true;
8640 Assert(lines[i]->vertex_index(i % 2) ==
8641 (middle_vertex_index<dim, spacedim>(
8642 hex->face(i / 4 - 1))),
8644 line_orientation[i] =
false;
8649 line_orientation[12] =
true;
8670 internal::TriangulationImplementation ::TriaObject<
8676 internal::TriangulationImplementation ::TriaObject<
8682 internal::TriangulationImplementation ::TriaObject<
8688 internal::TriangulationImplementation ::TriaObject<
8689 2>(line_indices[12],
8694 new_quads[0]->set_line_orientation(
8695 0, line_orientation[6]);
8696 new_quads[0]->set_line_orientation(
8697 1, line_orientation[10]);
8698 new_quads[0]->set_line_orientation(
8699 2, line_orientation[0]);
8701 new_quads[1]->set_line_orientation(
8702 0, line_orientation[7]);
8703 new_quads[1]->set_line_orientation(
8704 1, line_orientation[11]);
8705 new_quads[1]->set_line_orientation(
8706 3, line_orientation[1]);
8708 new_quads[2]->set_line_orientation(
8709 0, line_orientation[2]);
8710 new_quads[2]->set_line_orientation(
8711 2, line_orientation[4]);
8712 new_quads[2]->set_line_orientation(
8713 3, line_orientation[8]);
8715 new_quads[3]->set_line_orientation(
8716 1, line_orientation[3]);
8717 new_quads[3]->set_line_orientation(
8718 2, line_orientation[5]);
8719 new_quads[3]->set_line_orientation(
8720 3, line_orientation[9]);
8754 const int quad_indices[20] = {
8755 new_quads[0]->index(),
8756 new_quads[1]->index(),
8757 new_quads[2]->index(),
8758 new_quads[3]->index(),
8760 hex->face(0)->isotropic_child_index(
8762 0, f_or[0], f_fl[0], f_ro[0])),
8763 hex->face(0)->isotropic_child_index(
8765 1, f_or[0], f_fl[0], f_ro[0])),
8766 hex->face(0)->isotropic_child_index(
8768 2, f_or[0], f_fl[0], f_ro[0])),
8769 hex->face(0)->isotropic_child_index(
8771 3, f_or[0], f_fl[0], f_ro[0])),
8773 hex->face(1)->isotropic_child_index(
8775 0, f_or[1], f_fl[1], f_ro[1])),
8776 hex->face(1)->isotropic_child_index(
8778 1, f_or[1], f_fl[1], f_ro[1])),
8779 hex->face(1)->isotropic_child_index(
8781 2, f_or[1], f_fl[1], f_ro[1])),
8782 hex->face(1)->isotropic_child_index(
8784 3, f_or[1], f_fl[1], f_ro[1])),
8786 hex->face(2)->child_index(
8787 child_at_origin[hex->face(2)->refinement_case() -
8788 1][f_fl[2]][f_ro[2]]),
8789 hex->face(2)->child_index(
8791 child_at_origin[hex->face(2)->refinement_case() -
8792 1][f_fl[2]][f_ro[2]]),
8794 hex->face(3)->child_index(
8795 child_at_origin[hex->face(3)->refinement_case() -
8796 1][f_fl[3]][f_ro[3]]),
8797 hex->face(3)->child_index(
8799 child_at_origin[hex->face(3)->refinement_case() -
8800 1][f_fl[3]][f_ro[3]]),
8802 hex->face(4)->child_index(
8803 child_at_origin[hex->face(4)->refinement_case() -
8804 1][f_fl[4]][f_ro[4]]),
8805 hex->face(4)->child_index(
8807 child_at_origin[hex->face(4)->refinement_case() -
8808 1][f_fl[4]][f_ro[4]]),
8810 hex->face(5)->child_index(
8811 child_at_origin[hex->face(5)->refinement_case() -
8812 1][f_fl[5]][f_ro[5]]),
8813 hex->face(5)->child_index(
8815 child_at_origin[hex->face(5)->refinement_case() -
8816 1][f_fl[5]][f_ro[5]])};
8819 internal::TriangulationImplementation ::TriaObject<
8827 internal::TriangulationImplementation ::TriaObject<
8835 internal::TriangulationImplementation ::TriaObject<
8843 internal::TriangulationImplementation ::TriaObject<
8879 triangulation.vertices_used[next_unused_vertex] ==
8881 ++next_unused_vertex;
8883 next_unused_vertex < triangulation.vertices.size(),
8885 "Internal error: During refinement, the triangulation wants to access an element of the 'vertices' array but it turns out that the array is not large enough."));
8886 triangulation.vertices_used[next_unused_vertex] =
8895 triangulation.vertices[next_unused_vertex] =
8896 hex->center(
true,
true);
8917 const unsigned int vertex_indices[7] = {
8918 middle_vertex_index<dim, spacedim>(hex->face(0)),
8919 middle_vertex_index<dim, spacedim>(hex->face(1)),
8920 middle_vertex_index<dim, spacedim>(hex->face(2)),
8921 middle_vertex_index<dim, spacedim>(hex->face(3)),
8922 middle_vertex_index<dim, spacedim>(hex->face(4)),
8923 middle_vertex_index<dim, spacedim>(hex->face(5)),
8924 next_unused_vertex};
8928 1>(vertex_indices[2], vertex_indices[6]));
8931 1>(vertex_indices[6], vertex_indices[3]));
8934 1>(vertex_indices[0], vertex_indices[6]));
8937 1>(vertex_indices[6], vertex_indices[1]));
8940 1>(vertex_indices[4], vertex_indices[6]));
8943 1>(vertex_indices[6], vertex_indices[5]));
9012 spacedim>::raw_line_iterator lines[30] = {
9016 0, f_or[0], f_fl[0], f_ro[0]))
9019 1, f_or[0], f_fl[0], f_ro[0])),
9023 3, f_or[0], f_fl[0], f_ro[0]))
9026 0, f_or[0], f_fl[0], f_ro[0])),
9030 0, f_or[0], f_fl[0], f_ro[0]))
9033 3, f_or[0], f_fl[0], f_ro[0])),
9037 3, f_or[0], f_fl[0], f_ro[0]))
9040 2, f_or[0], f_fl[0], f_ro[0])),
9045 0, f_or[1], f_fl[1], f_ro[1]))
9048 1, f_or[1], f_fl[1], f_ro[1])),
9052 3, f_or[1], f_fl[1], f_ro[1]))
9055 0, f_or[1], f_fl[1], f_ro[1])),
9059 0, f_or[1], f_fl[1], f_ro[1]))
9062 3, f_or[1], f_fl[1], f_ro[1])),
9066 3, f_or[1], f_fl[1], f_ro[1]))
9069 2, f_or[1], f_fl[1], f_ro[1])),
9074 0, f_or[2], f_fl[2], f_ro[2]))
9077 1, f_or[2], f_fl[2], f_ro[2])),
9081 3, f_or[2], f_fl[2], f_ro[2]))
9084 0, f_or[2], f_fl[2], f_ro[2])),
9088 0, f_or[2], f_fl[2], f_ro[2]))
9091 3, f_or[2], f_fl[2], f_ro[2])),
9095 3, f_or[2], f_fl[2], f_ro[2]))
9098 2, f_or[2], f_fl[2], f_ro[2])),
9103 0, f_or[3], f_fl[3], f_ro[3]))
9106 1, f_or[3], f_fl[3], f_ro[3])),
9110 3, f_or[3], f_fl[3], f_ro[3]))
9113 0, f_or[3], f_fl[3], f_ro[3])),
9117 0, f_or[3], f_fl[3], f_ro[3]))
9120 3, f_or[3], f_fl[3], f_ro[3])),
9124 3, f_or[3], f_fl[3], f_ro[3]))
9127 2, f_or[3], f_fl[3], f_ro[3])),
9132 0, f_or[4], f_fl[4], f_ro[4]))
9135 1, f_or[4], f_fl[4], f_ro[4])),
9139 3, f_or[4], f_fl[4], f_ro[4]))
9142 0, f_or[4], f_fl[4], f_ro[4])),
9146 0, f_or[4], f_fl[4], f_ro[4]))
9149 3, f_or[4], f_fl[4], f_ro[4])),
9153 3, f_or[4], f_fl[4], f_ro[4]))
9156 2, f_or[4], f_fl[4], f_ro[4])),
9161 0, f_or[5], f_fl[5], f_ro[5]))
9164 1, f_or[5], f_fl[5], f_ro[5])),
9168 3, f_or[5], f_fl[5], f_ro[5]))
9171 0, f_or[5], f_fl[5], f_ro[5])),
9175 0, f_or[5], f_fl[5], f_ro[5]))
9178 3, f_or[5], f_fl[5], f_ro[5])),
9182 3, f_or[5], f_fl[5], f_ro[5]))
9185 2, f_or[5], f_fl[5], f_ro[5])),
9195 unsigned int line_indices[30];
9196 for (
unsigned int i = 0; i < 30; ++i)
9197 line_indices[i] = lines[i]->index();
9204 bool line_orientation[30];
9212 for (
unsigned int i = 0; i < 24; ++i)
9213 if (lines[i]->vertex_index((i + 1) % 2) ==
9214 vertex_indices[i / 4])
9215 line_orientation[i] =
true;
9220 Assert(lines[i]->vertex_index(i % 2) ==
9221 vertex_indices[i / 4],
9223 line_orientation[i] =
false;
9228 for (
unsigned int i = 24; i < 30; ++i)
9229 line_orientation[i] =
true;
9262 internal::TriangulationImplementation ::TriaObject<
9263 2>(line_indices[10],
9268 internal::TriangulationImplementation ::TriaObject<
9269 2>(line_indices[28],
9274 internal::TriangulationImplementation ::TriaObject<
9275 2>(line_indices[11],
9280 internal::TriangulationImplementation ::TriaObject<
9281 2>(line_indices[29],
9286 internal::TriangulationImplementation ::TriaObject<
9287 2>(line_indices[18],
9292 internal::TriangulationImplementation ::TriaObject<
9293 2>(line_indices[26],
9298 internal::TriangulationImplementation ::TriaObject<
9299 2>(line_indices[19],
9304 internal::TriangulationImplementation ::TriaObject<
9305 2>(line_indices[27],
9310 internal::TriangulationImplementation ::TriaObject<
9316 internal::TriangulationImplementation ::TriaObject<
9317 2>(line_indices[24],
9322 internal::TriangulationImplementation ::TriaObject<
9328 internal::TriangulationImplementation ::TriaObject<
9329 2>(line_indices[25],
9337 new_quads[0]->set_line_orientation(
9338 0, line_orientation[10]);
9339 new_quads[0]->set_line_orientation(
9340 2, line_orientation[16]);
9342 new_quads[1]->set_line_orientation(
9343 1, line_orientation[14]);
9344 new_quads[1]->set_line_orientation(
9345 2, line_orientation[17]);
9347 new_quads[2]->set_line_orientation(
9348 0, line_orientation[11]);
9349 new_quads[2]->set_line_orientation(
9350 3, line_orientation[20]);
9352 new_quads[3]->set_line_orientation(
9353 1, line_orientation[15]);
9354 new_quads[3]->set_line_orientation(
9355 3, line_orientation[21]);
9357 new_quads[4]->set_line_orientation(
9358 0, line_orientation[18]);
9359 new_quads[4]->set_line_orientation(
9360 2, line_orientation[0]);
9362 new_quads[5]->set_line_orientation(
9363 1, line_orientation[22]);
9364 new_quads[5]->set_line_orientation(
9365 2, line_orientation[1]);
9367 new_quads[6]->set_line_orientation(
9368 0, line_orientation[19]);
9369 new_quads[6]->set_line_orientation(
9370 3, line_orientation[4]);
9372 new_quads[7]->set_line_orientation(
9373 1, line_orientation[23]);
9374 new_quads[7]->set_line_orientation(
9375 3, line_orientation[5]);
9377 new_quads[8]->set_line_orientation(
9378 0, line_orientation[2]);
9379 new_quads[8]->set_line_orientation(
9380 2, line_orientation[8]);
9382 new_quads[9]->set_line_orientation(
9383 1, line_orientation[6]);
9384 new_quads[9]->set_line_orientation(
9385 2, line_orientation[9]);
9387 new_quads[10]->set_line_orientation(
9388 0, line_orientation[3]);
9389 new_quads[10]->set_line_orientation(
9390 3, line_orientation[12]);
9392 new_quads[11]->set_line_orientation(
9393 1, line_orientation[7]);
9394 new_quads[11]->set_line_orientation(
9395 3, line_orientation[13]);
9436 const int quad_indices[36] = {
9437 new_quads[0]->index(),
9438 new_quads[1]->index(),
9439 new_quads[2]->index(),
9440 new_quads[3]->index(),
9441 new_quads[4]->index(),
9442 new_quads[5]->index(),
9443 new_quads[6]->index(),
9444 new_quads[7]->index(),
9445 new_quads[8]->index(),
9446 new_quads[9]->index(),
9447 new_quads[10]->index(),
9448 new_quads[11]->index(),
9450 hex->face(0)->isotropic_child_index(
9452 0, f_or[0], f_fl[0], f_ro[0])),
9453 hex->face(0)->isotropic_child_index(
9455 1, f_or[0], f_fl[0], f_ro[0])),
9456 hex->face(0)->isotropic_child_index(
9458 2, f_or[0], f_fl[0], f_ro[0])),
9459 hex->face(0)->isotropic_child_index(
9461 3, f_or[0], f_fl[0], f_ro[0])),
9463 hex->face(1)->isotropic_child_index(
9465 0, f_or[1], f_fl[1], f_ro[1])),
9466 hex->face(1)->isotropic_child_index(
9468 1, f_or[1], f_fl[1], f_ro[1])),
9469 hex->face(1)->isotropic_child_index(
9471 2, f_or[1], f_fl[1], f_ro[1])),
9472 hex->face(1)->isotropic_child_index(
9474 3, f_or[1], f_fl[1], f_ro[1])),
9476 hex->face(2)->isotropic_child_index(
9478 0, f_or[2], f_fl[2], f_ro[2])),
9479 hex->face(2)->isotropic_child_index(
9481 1, f_or[2], f_fl[2], f_ro[2])),
9482 hex->face(2)->isotropic_child_index(
9484 2, f_or[2], f_fl[2], f_ro[2])),
9485 hex->face(2)->isotropic_child_index(
9487 3, f_or[2], f_fl[2], f_ro[2])),
9489 hex->face(3)->isotropic_child_index(
9491 0, f_or[3], f_fl[3], f_ro[3])),
9492 hex->face(3)->isotropic_child_index(
9494 1, f_or[3], f_fl[3], f_ro[3])),
9495 hex->face(3)->isotropic_child_index(
9497 2, f_or[3], f_fl[3], f_ro[3])),
9498 hex->face(3)->isotropic_child_index(
9500 3, f_or[3], f_fl[3], f_ro[3])),
9502 hex->face(4)->isotropic_child_index(
9504 0, f_or[4], f_fl[4], f_ro[4])),
9505 hex->face(4)->isotropic_child_index(
9507 1, f_or[4], f_fl[4], f_ro[4])),
9508 hex->face(4)->isotropic_child_index(
9510 2, f_or[4], f_fl[4], f_ro[4])),
9511 hex->face(4)->isotropic_child_index(
9513 3, f_or[4], f_fl[4], f_ro[4])),
9515 hex->face(5)->isotropic_child_index(
9517 0, f_or[5], f_fl[5], f_ro[5])),
9518 hex->face(5)->isotropic_child_index(
9520 1, f_or[5], f_fl[5], f_ro[5])),
9521 hex->face(5)->isotropic_child_index(
9523 2, f_or[5], f_fl[5], f_ro[5])),
9524 hex->face(5)->isotropic_child_index(
9526 3, f_or[5], f_fl[5], f_ro[5]))};
9530 internal::TriangulationImplementation ::TriaObject<
9531 3>(quad_indices[12],
9538 internal::TriangulationImplementation ::TriaObject<
9546 internal::TriangulationImplementation ::TriaObject<
9547 3>(quad_indices[13],
9554 internal::TriangulationImplementation ::TriaObject<
9564 internal::TriangulationImplementation ::TriaObject<
9565 3>(quad_indices[14],
9572 internal::TriangulationImplementation ::TriaObject<
9580 internal::TriangulationImplementation ::TriaObject<
9581 3>(quad_indices[15],
9588 internal::TriangulationImplementation ::TriaObject<
9627 for (
unsigned int f = 0;
9628 f < GeometryInfo<dim>::faces_per_cell;
9630 for (
unsigned int s = 0;
9637 const unsigned int current_child =
9646 ref_case, f, f_or[f], f_fl[f], f_ro[f]));
9647 new_hexes[current_child]->set_face_orientation(f,
9649 new_hexes[current_child]->set_face_flip(f, f_fl[f]);
9650 new_hexes[current_child]->set_face_rotation(f, f_ro[f]);
9655 if ((check_for_distorted_cells ==
true) &&
9656 has_distorted_children(
9658 std::integral_constant<int, dim>(),
9659 std::integral_constant<int, spacedim>()))
9667 triangulation.signals.post_refinement_on_cell(hex);
9676 triangulation.faces->quads.clear_user_data();
9679 return cells_with_distorted_children;
9695 template <
int spacedim>
9700 template <
int dim,
int spacedim>
9711 triangulation.begin();
9712 cell != triangulation.end();
9714 if (cell->at_boundary() && cell->refine_flag_set() &&
9715 cell->refine_flag_set() !=
9722 for (
unsigned int face_no = 0;
9723 face_no < GeometryInfo<dim>::faces_per_cell;
9725 if (cell->face(face_no)->at_boundary())
9736 spacedim>::face_iterator
9737 face = cell->face(face_no);
9745 .transform_real_to_unit_cell(cell, new_bound);
9756 std::fabs(new_unit[face_no / 2] - face_no % 2);
9761 const double allowed = 0.25;
9764 cell->flag_for_face_refinement(face_no);
9783 template <
int dim,
int spacedim>
9788 ExcMessage(
"Wrong function called -- there should " 9789 "be a specialization."));
9793 template <
int spacedim>
9797 const unsigned int dim = 3;
9801 triangulation.clear_user_flags_line();
9805 triangulation.clear_user_flags_hex();
9809 bool mesh_changed =
false;
9813 mesh_changed =
false;
9818 triangulation.fix_coarsen_flags();
9824 triangulation.begin();
9825 cell != triangulation.end();
9827 if (cell->refine_flag_set())
9829 for (
unsigned int line = 0;
9830 line < GeometryInfo<dim>::lines_per_cell;
9833 cell->refine_flag_set(), line) ==
9837 cell->line(line)->set_user_flag();
9839 else if (cell->has_children() &&
9840 !cell->child(0)->coarsen_flag_set())
9842 for (
unsigned int line = 0;
9843 line < GeometryInfo<dim>::lines_per_cell;
9846 cell->refinement_case(), line) ==
9850 cell->line(line)->set_user_flag();
9852 else if (cell->has_children() &&
9853 cell->child(0)->coarsen_flag_set())
9854 cell->set_user_flag();
9862 cell = triangulation.last_active();
9863 cell != triangulation.end();
9865 for (
unsigned int line = 0;
9866 line < GeometryInfo<dim>::lines_per_cell;
9869 if (cell->line(line)->has_children())
9878 bool offending_line_found =
false;
9880 for (
unsigned int c = 0; c < 2; ++c)
9882 Assert(cell->line(line)->child(c)->has_children() ==
9886 if (cell->line(line)->child(c)->user_flag_set() &&
9888 cell->refine_flag_set(), line) ==
9892 cell->clear_coarsen_flag();
9897 if (triangulation.smooth_grid &
9899 allow_anisotropic_smoothing)
9900 cell->flag_for_line_refinement(line);
9902 cell->set_refine_flag();
9904 for (
unsigned int l = 0;
9905 l < GeometryInfo<dim>::lines_per_cell;
9908 cell->refine_flag_set(), line) ==
9911 cell->line(l)->set_user_flag();
9914 offending_line_found =
true;
9920 for (
unsigned int l = 0;
9921 l < GeometryInfo<dim>::lines_per_cell;
9923 if (!cell->line(l)->has_children() &&
9925 cell->refine_flag_set(), l) !=
9927 cell->line(l)->set_user_flag();
9933 if (offending_line_found)
9935 mesh_changed =
true;
9953 triangulation.last();
9954 cell != triangulation.end();
9957 if (cell->user_flag_set())
9958 for (
unsigned int line = 0;
9959 line < GeometryInfo<dim>::lines_per_cell;
9961 if (cell->line(line)->has_children() &&
9962 (cell->line(line)->child(0)->user_flag_set() ||
9963 cell->line(line)->child(1)->user_flag_set()))
9965 for (
unsigned int c = 0; c < cell->n_children(); ++c)
9966 cell->child(c)->clear_coarsen_flag();
9967 cell->clear_user_flag();
9968 for (
unsigned int l = 0;
9969 l < GeometryInfo<dim>::lines_per_cell;
9972 cell->refinement_case(),
l) ==
9976 cell->line(l)->set_user_flag();
9977 mesh_changed =
true;
9982 while (mesh_changed ==
true);
9993 template <
int dim,
int spacedim>
10004 for (
unsigned int n = 0; n < GeometryInfo<dim>::faces_per_cell; ++n)
10013 const unsigned int n_subfaces =
10016 if (n_subfaces == 0 || cell->at_boundary(n))
10018 for (
unsigned int c = 0; c < n_subfaces; ++c)
10021 child = cell->child(
10025 child_neighbor = child->neighbor(n);
10026 if (!child->neighbor_is_coarser(n))
10039 if ((child_neighbor->has_children() &&
10040 !child_neighbor->user_flag_set()) ||
10049 child_neighbor->refine_flag_set())
10059 template <
int dim,
int spacedim>
10061 get_default_flat_manifold()
10064 return flat_manifold;
10071 template <
int dim,
int spacedim>
10076 template <
int dim,
int spacedim>
10079 const bool check_for_distorted_cells)
10080 : smooth_grid(smooth_grid)
10081 , anisotropic_refinement(false)
10082 , check_for_distorted_cells(check_for_distorted_cells)
10087 std_cxx14::make_unique<std::map<unsigned int, types::boundary_id>>();
10089 std_cxx14::make_unique<std::map<unsigned int, types::manifold_id>>();
10101 template <
int dim,
int spacedim>
10105 , smooth_grid(tria.smooth_grid)
10106 , periodic_face_pairs_level_0(std::move(tria.periodic_face_pairs_level_0))
10107 , periodic_face_map(std::move(tria.periodic_face_map))
10108 , levels(std::move(tria.levels))
10109 , faces(std::move(tria.faces))
10110 , vertices(std::move(tria.vertices))
10111 , vertices_used(std::move(tria.vertices_used))
10112 , manifold(std::move(tria.manifold))
10113 , anisotropic_refinement(tria.anisotropic_refinement)
10114 , check_for_distorted_cells(tria.check_for_distorted_cells)
10115 , number_cache(std::move(tria.number_cache))
10116 , vertex_to_boundary_id_map_1d(std::move(tria.vertex_to_boundary_id_map_1d))
10117 , vertex_to_manifold_id_map_1d(std::move(tria.vertex_to_manifold_id_map_1d))
10123 template <
int dim,
int spacedim>
10130 smooth_grid = tria.smooth_grid;
10131 periodic_face_pairs_level_0 = std::move(tria.periodic_face_pairs_level_0);
10132 periodic_face_map = std::move(tria.periodic_face_map);
10133 levels = std::move(tria.levels);
10134 faces = std::move(tria.faces);
10135 vertices = std::move(tria.vertices);
10136 vertices_used = std::move(tria.vertices_used);
10137 manifold = std::move(tria.manifold);
10138 anisotropic_refinement = tria.anisotropic_refinement;
10139 number_cache = tria.number_cache;
10140 vertex_to_boundary_id_map_1d = std::move(tria.vertex_to_boundary_id_map_1d);
10141 vertex_to_manifold_id_map_1d = std::move(tria.vertex_to_manifold_id_map_1d);
10150 template <
int dim,
int spacedim>
10167 AssertNothrow((dim == 1) || (vertex_to_boundary_id_map_1d ==
nullptr),
10172 AssertNothrow((dim == 1) || (vertex_to_manifold_id_map_1d ==
nullptr),
10178 template <
int dim,
int spacedim>
10186 clear_despite_subscriptions();
10187 periodic_face_pairs_level_0.clear();
10188 periodic_face_map.clear();
10193 template <
int dim,
int spacedim>
10199 ExcTriangulationNotEmpty(vertices.size(), levels.size()));
10200 smooth_grid = mesh_smoothing;
10205 template <
int dim,
int spacedim>
10209 return smooth_grid;
10214 template <
int dim,
int spacedim>
10223 manifold[m_number] = manifold_object.
clone();
10228 template <
int dim,
int spacedim>
10232 reset_manifold(m_number);
10236 template <
int dim,
int spacedim>
10244 manifold.erase(m_number);
10248 template <
int dim,
int spacedim>
10256 template <
int dim,
int spacedim>
10264 "Error: set_all_manifold_ids() can not be called on an empty Triangulation."));
10267 cell = this->begin_active(),
10268 endc = this->end();
10270 for (; cell != endc; ++cell)
10271 cell->set_all_manifold_ids(m_number);
10275 template <
int dim,
int spacedim>
10283 "Error: set_all_manifold_ids_on_boundary() can not be called on an empty Triangulation."));
10286 cell = this->begin_active(),
10287 endc = this->end();
10289 for (; cell != endc; ++cell)
10290 for (
unsigned int f = 0; f < GeometryInfo<dim>::faces_per_cell; ++f)
10291 if (cell->face(f)->at_boundary())
10292 cell->face(f)->set_all_manifold_ids(m_number);
10296 template <
int dim,
int spacedim>
10305 "Error: set_all_manifold_ids_on_boundary() can not be called on an empty Triangulation."));
10307 bool boundary_found =
false;
10309 cell = this->begin_active(),
10310 endc = this->end();
10312 for (; cell != endc; ++cell)
10315 for (
unsigned int f = 0; f < GeometryInfo<dim>::faces_per_cell; ++f)
10316 if (cell->face(f)->at_boundary() &&
10317 cell->face(f)->boundary_id() == b_id)
10319 boundary_found =
true;
10320 cell->face(f)->set_manifold_id(m_number);
10325 for (
unsigned int e = 0; e < GeometryInfo<dim>::lines_per_cell; ++e)
10326 if (cell->line(e)->at_boundary() &&
10327 cell->line(e)->boundary_id() == b_id)
10329 boundary_found =
true;
10330 cell->line(e)->set_manifold_id(m_number);
10334 (void)boundary_found;
10335 Assert(boundary_found, ExcBoundaryIdNotFound(b_id));
10340 template <
int dim,
int spacedim>
10347 const auto it = manifold.find(m_number);
10349 if (it != manifold.end())
10352 return *(it->second);
10357 return internal::TriangulationImplementation::
10358 get_default_flat_manifold<dim, spacedim>();
10363 template <
int dim,
int spacedim>
10364 std::vector<types::boundary_id>
10371 std::vector<types::boundary_id> boundary_ids;
10372 for (std::map<unsigned int, types::boundary_id>::const_iterator p =
10373 vertex_to_boundary_id_map_1d->begin();
10374 p != vertex_to_boundary_id_map_1d->end();
10376 boundary_ids.push_back(p->second);
10378 return boundary_ids;
10382 std::set<types::boundary_id> b_ids;
10384 for (; cell != end(); ++cell)
10385 for (
unsigned int face = 0; face < GeometryInfo<dim>::faces_per_cell;
10387 if (cell->at_boundary(face))
10388 b_ids.insert(cell->face(face)->boundary_id());
10389 std::vector<types::boundary_id> boundary_ids(b_ids.begin(), b_ids.end());
10390 return boundary_ids;
10396 template <
int dim,
int spacedim>
10397 std::vector<types::manifold_id>
10400 std::set<types::manifold_id> m_ids;
10402 for (; cell != end(); ++cell)
10404 m_ids.insert(cell->manifold_id());
10406 for (
unsigned int face = 0; face < GeometryInfo<dim>::faces_per_cell;
10408 if (cell->at_boundary(face))
10409 m_ids.insert(cell->face(face)->manifold_id());
10411 std::vector<types::manifold_id> manifold_indicators(m_ids.begin(),
10413 return manifold_indicators;
10419 template <
int dim,
int spacedim>
10424 Assert((vertices.size() == 0) && (levels.size() == 0) && (faces ==
nullptr),
10425 ExcTriangulationNotEmpty(vertices.size(), levels.size()));
10427 (dim == 1 || other_tria.
faces !=
nullptr),
10429 "When calling Triangulation::copy_triangulation(), " 10430 "the target triangulation must be empty but the source " 10431 "triangulation (the argument to this function) must contain " 10432 "something. Here, it seems like the source does not " 10433 "contain anything at all."));
10443 faces = std_cxx14::make_unique<
10446 auto bdry_iterator = other_tria.
manifold.begin();
10447 for (; bdry_iterator != other_tria.
manifold.end(); ++bdry_iterator)
10448 manifold[bdry_iterator->first] = bdry_iterator->second->clone();
10451 levels.reserve(other_tria.
levels.size());
10452 for (
unsigned int level = 0; level < other_tria.
levels.size(); ++level)
10453 levels.push_back(std_cxx14::make_unique<
10455 *other_tria.
levels[level]));
10461 vertex_to_boundary_id_map_1d =
10462 std_cxx14::make_unique<std::map<unsigned int, types::boundary_id>>(
10465 vertex_to_manifold_id_map_1d =
10466 std_cxx14::make_unique<std::map<unsigned int, types::manifold_id>>(
10471 other_tria.
signals.copy(*
this);
10482 template <
int dim,
int spacedim>
10489 std::vector<CellData<dim>> reordered_cells(cells);
10493 reorder_compatibility(reordered_cells, reordered_subcelldata);
10497 create_triangulation(v, reordered_cells, reordered_subcelldata);
10502 template <
int dim,
int spacedim>
10509 Assert((vertices.size() == 0) && (levels.size() == 0) && (faces ==
nullptr),
10510 ExcTriangulationNotEmpty(vertices.size(), levels.size()));
10525 clear_despite_subscriptions();
10532 *
this, levels.size(), number_cache);
10533 reset_active_cell_indices();
10538 if (check_for_distorted_cells ==
true)
10540 DistortedCellList distorted_cells = collect_distorted_coarse_cells(*
this);
10544 AssertThrow(distorted_cells.distorted_cells.size() == 0, distorted_cells);
10577 if (dim < spacedim)
10585 bool values[][2] = {{
false,
true}, {
true,
false}};
10586 for (
unsigned int i = 0; i < GeometryInfo<dim>::faces_per_cell;
10588 for (
unsigned int j = 0; j < GeometryInfo<dim>::faces_per_cell;
10590 correct(i, j) = (values[i][j]);
10595 bool values[][4] = {{
false,
true,
true,
false},
10596 {
true,
false,
false,
true},
10597 {
true,
false,
false,
true},
10598 {
false,
true,
true,
false}};
10599 for (
unsigned int i = 0; i < GeometryInfo<dim>::faces_per_cell;
10601 for (
unsigned int j = 0; j < GeometryInfo<dim>::faces_per_cell;
10603 correct(i, j) = (values[i][j]);
10611 std::list<active_cell_iterator> this_round, next_round;
10614 this_round.push_back(begin_active());
10615 begin_active()->set_direction_flag(
true);
10616 begin_active()->set_user_flag();
10618 while (this_round.size() > 0)
10620 for (
typename std::list<active_cell_iterator>::iterator cell =
10621 this_round.begin();
10622 cell != this_round.end();
10625 for (
unsigned int i = 0; i < GeometryInfo<dim>::faces_per_cell;
10628 if (!((*cell)->face(i)->at_boundary()))
10630 neighbor = (*cell)->neighbor(i);
10632 unsigned int cf = (*cell)->face_index(i);
10633 unsigned int j = 0;
10634 while (neighbor->face_index(j) != cf)
10642 if (neighbor->user_flag_set())
10646 Assert(!(correct(i, j) ^
10647 (neighbor->direction_flag() ==
10648 (*cell)->direction_flag())),
10649 ExcNonOrientableTriangulation());
10653 next_round.push_back(neighbor);
10654 neighbor->set_user_flag();
10655 if ((correct(i, j) ^ (neighbor->direction_flag() ==
10656 (*cell)->direction_flag())))
10657 neighbor->set_direction_flag(
10658 !neighbor->direction_flag());
10668 if (next_round.size() == 0)
10671 if (cell->user_flag_set() ==
false)
10673 next_round.push_back(cell);
10674 cell->set_direction_flag(
true);
10675 cell->set_user_flag();
10679 this_round = next_round;
10680 next_round.clear();
10690 template <
int dim,
int spacedim>
10695 ExcMessage(
"Only works for dim == spacedim-1"));
10697 cell->set_direction_flag(!cell->direction_flag());
10702 template <
int dim,
int spacedim>
10707 ExcMessage(
"Error: An empty Triangulation can not be refined."));
10710 for (; cell != endc; ++cell)
10712 cell->clear_coarsen_flag();
10713 cell->set_refine_flag();
10719 template <
int dim,
int spacedim>
10723 for (
unsigned int i = 0; i < times; ++i)
10725 set_all_refine_flags();
10726 execute_coarsening_and_refinement();
10736 template <
int dim,
int spacedim>
10740 v.resize(dim * n_active_cells(),
false);
10741 std::vector<bool>::iterator i = v.begin();
10743 for (; cell != endc; ++cell)
10744 for (
unsigned int j = 0; j < dim; ++j, ++i)
10745 if (cell->refine_flag_set() & (1 << j))
10753 template <
int dim,
int spacedim>
10757 std::vector<bool> v;
10758 save_refine_flags(v);
10759 write_bool_vector(mn_tria_refine_flags_begin,
10761 mn_tria_refine_flags_end,
10767 template <
int dim,
int spacedim>
10771 std::vector<bool> v;
10772 read_bool_vector(mn_tria_refine_flags_begin, v, mn_tria_refine_flags_end, in);
10773 load_refine_flags(v);
10778 template <
int dim,
int spacedim>
10782 AssertThrow(v.size() == dim * n_active_cells(), ExcGridReadError());
10785 std::vector<bool>::const_iterator i = v.begin();
10786 for (; cell != endc; ++cell)
10788 unsigned int ref_case = 0;
10790 for (
unsigned int j = 0; j < dim; ++j, ++i)
10792 ref_case += 1 << j;
10794 ExcGridReadError());
10798 cell->clear_refine_flag();
10806 template <
int dim,
int spacedim>
10810 v.resize(n_active_cells(),
false);
10811 std::vector<bool>::iterator i = v.begin();
10813 for (; cell != endc; ++cell, ++i)
10814 *i = cell->coarsen_flag_set();
10821 template <
int dim,
int spacedim>
10825 std::vector<bool> v;
10826 save_coarsen_flags(v);
10827 write_bool_vector(mn_tria_coarsen_flags_begin,
10829 mn_tria_coarsen_flags_end,
10835 template <
int dim,
int spacedim>
10839 std::vector<bool> v;
10840 read_bool_vector(mn_tria_coarsen_flags_begin,
10842 mn_tria_coarsen_flags_end,
10844 load_coarsen_flags(v);
10849 template <
int dim,
int spacedim>
10853 Assert(v.size() == n_active_cells(), ExcGridReadError());
10856 std::vector<bool>::const_iterator i = v.begin();
10857 for (; cell != endc; ++cell, ++i)
10859 cell->set_coarsen_flag();
10861 cell->clear_coarsen_flag();
10867 template <
int dim,
int spacedim>
10871 return anisotropic_refinement;
10889 for (
unsigned int level = 0; level < levels.size(); ++level)
10890 levels[level]->cells.clear_user_data();
10904 faces->
lines.clear_user_data();
10911 faces->
lines.clear_user_data();
10917 template <
int dim,
int spacedim>
10922 ::clear_user_data(levels);
10923 ::clear_user_data(faces.get());
10930 void clear_user_flags_line(
10936 for (
const auto &level : levels)
10937 level->cells.clear_user_flags();
10942 clear_user_flags_line(
10947 faces->lines.clear_user_flags();
10952 template <
int dim,
int spacedim>
10956 ::clear_user_flags_line(levels, faces.get());
10963 void clear_user_flags_quad(
10971 void clear_user_flags_quad(
10977 for (
const auto &level : levels)
10978 level->cells.clear_user_flags();
10983 clear_user_flags_quad(
10988 faces->quads.clear_user_flags();
10993 template <
int dim,
int spacedim>
10997 ::clear_user_flags_quad(levels, faces.get());
11004 void clear_user_flags_hex(
11013 void clear_user_flags_hex(
11021 void clear_user_flags_hex(
11027 for (
const auto &level : levels)
11028 level->cells.clear_user_flags();
11033 template <
int dim,
int spacedim>
11037 ::clear_user_flags_hex(levels, faces.get());
11042 template <
int dim,
int spacedim>
11046 clear_user_flags_line();
11047 clear_user_flags_quad();
11048 clear_user_flags_hex();
11053 template <
int dim,
int spacedim>
11057 save_user_flags_line(out);
11060 save_user_flags_quad(out);
11063 save_user_flags_hex(out);
11071 template <
int dim,
int spacedim>
11079 std::vector<bool> tmp;
11081 save_user_flags_line(tmp);
11082 v.insert(v.end(), tmp.begin(), tmp.end());
11086 save_user_flags_quad(tmp);
11087 v.insert(v.end(), tmp.begin(), tmp.end());
11092 save_user_flags_hex(tmp);
11093 v.insert(v.end(), tmp.begin(), tmp.end());
11102 template <
int dim,
int spacedim>
11106 load_user_flags_line(in);
11109 load_user_flags_quad(in);
11112 load_user_flags_hex(in);
11120 template <
int dim,
int spacedim>
11125 std::vector<bool> tmp;
11129 tmp.insert(tmp.end(), v.begin(), v.begin() + n_lines());
11131 load_user_flags_line(tmp);
11136 tmp.insert(tmp.end(),
11137 v.begin() + n_lines(),
11138 v.begin() + n_lines() + n_quads());
11139 load_user_flags_quad(tmp);
11145 tmp.insert(tmp.end(),
11146 v.begin() + n_lines() + n_quads(),
11147 v.begin() + n_lines() + n_quads() + n_hexs());
11148 load_user_flags_hex(tmp);
11157 template <
int dim,
int spacedim>
11161 v.resize(n_lines(),
false);
11162 std::vector<bool>::iterator i = v.begin();
11164 for (; line != endl; ++line, ++i)
11165 *i = line->user_flag_set();
11172 template <
int dim,
int spacedim>
11176 std::vector<bool> v;
11177 save_user_flags_line(v);
11178 write_bool_vector(mn_tria_line_user_flags_begin,
11180 mn_tria_line_user_flags_end,
11186 template <
int dim,
int spacedim>
11190 std::vector<bool> v;
11191 read_bool_vector(mn_tria_line_user_flags_begin,
11193 mn_tria_line_user_flags_end,
11195 load_user_flags_line(v);
11200 template <
int dim,
int spacedim>
11204 Assert(v.size() == n_lines(), ExcGridReadError());
11207 std::vector<bool>::const_iterator i = v.begin();
11208 for (; line != endl; ++line, ++i)
11210 line->set_user_flag();
11212 line->clear_user_flag();
11220 template <
typename Iterator>
11222 get_user_flag(
const Iterator &i)
11224 return i->user_flag_set();
11229 template <
int structdim,
int dim,
int spacedim>
11239 template <
typename Iterator>
11241 set_user_flag(
const Iterator &i)
11243 i->set_user_flag();
11248 template <
int structdim,
int dim,
int spacedim>
11257 template <
typename Iterator>
11259 clear_user_flag(
const Iterator &i)
11261 i->clear_user_flag();
11266 template <
int structdim,
int dim,
int spacedim>
11276 template <
int dim,
int spacedim>
11280 v.resize(n_quads(),
false);
11284 std::vector<bool>::iterator i = v.begin();
11286 for (; quad != endq; ++quad, ++i)
11287 *i = get_user_flag(quad);
11295 template <
int dim,
int spacedim>
11299 std::vector<bool> v;
11300 save_user_flags_quad(v);
11301 write_bool_vector(mn_tria_quad_user_flags_begin,
11303 mn_tria_quad_user_flags_end,
11309 template <
int dim,
int spacedim>
11313 std::vector<bool> v;
11314 read_bool_vector(mn_tria_quad_user_flags_begin,
11316 mn_tria_quad_user_flags_end,
11318 load_user_flags_quad(v);
11323 template <
int dim,
int spacedim>
11327 Assert(v.size() == n_quads(), ExcGridReadError());
11332 std::vector<bool>::const_iterator i = v.begin();
11333 for (; quad != endq; ++quad, ++i)
11335 set_user_flag(quad);
11337 clear_user_flag(quad);
11345 template <
int dim,
int spacedim>
11349 v.resize(n_hexs(),
false);
11353 std::vector<bool>::iterator i = v.begin();
11355 for (; hex != endh; ++hex, ++i)
11356 *i = get_user_flag(hex);
11364 template <
int dim,
int spacedim>
11368 std::vector<bool> v;
11369 save_user_flags_hex(v);
11370 write_bool_vector(mn_tria_hex_user_flags_begin,
11372 mn_tria_hex_user_flags_end,
11378 template <
int dim,
int spacedim>
11382 std::vector<bool> v;
11383 read_bool_vector(mn_tria_hex_user_flags_begin,
11385 mn_tria_hex_user_flags_end,
11387 load_user_flags_hex(v);
11392 template <
int dim,
int spacedim>
11396 Assert(v.size() == n_hexs(), ExcGridReadError());
11401 std::vector<bool>::const_iterator i = v.begin();
11402 for (; hex != endh; ++hex, ++i)
11404 set_user_flag(hex);
11406 clear_user_flag(hex);
11414 template <
int dim,
int spacedim>
11417 std::vector<unsigned int> &v)
const 11423 std::vector<unsigned int> tmp;
11425 save_user_indices_line(tmp);
11426 v.insert(v.end(), tmp.begin(), tmp.end());
11430 save_user_indices_quad(tmp);
11431 v.insert(v.end(), tmp.begin(), tmp.end());
11436 save_user_indices_hex(tmp);
11437 v.insert(v.end(), tmp.begin(), tmp.end());
11446 template <
int dim,
int spacedim>
11449 const std::vector<unsigned int> &v)
11452 std::vector<unsigned int> tmp;
11456 tmp.insert(tmp.end(), v.begin(), v.begin() + n_lines());
11458 load_user_indices_line(tmp);
11463 tmp.insert(tmp.end(),
11464 v.begin() + n_lines(),
11465 v.begin() + n_lines() + n_quads());
11466 load_user_indices_quad(tmp);
11472 tmp.insert(tmp.end(),
11473 v.begin() + n_lines() + n_quads(),
11474 v.begin() + n_lines() + n_quads() + n_hexs());
11475 load_user_indices_hex(tmp);
11486 template <
typename Iterator>
11488 get_user_index(
const Iterator &i)
11490 return i->user_index();
11495 template <
int structdim,
int dim,
int spacedim>
11506 template <
typename Iterator>
11508 set_user_index(
const Iterator &i,
const unsigned int x)
11510 i->set_user_index(x);
11515 template <
int structdim,
int dim,
int spacedim>
11519 const unsigned int)
11526 template <
int dim,
int spacedim>
11529 std::vector<unsigned int> &v)
const 11531 v.resize(n_lines(), 0);
11532 std::vector<unsigned int>::iterator i = v.begin();
11534 for (; line != endl; ++line, ++i)
11535 *i = line->user_index();
11540 template <
int dim,
int spacedim>
11543 const std::vector<unsigned int> &v)
11545 Assert(v.size() == n_lines(), ExcGridReadError());
11548 std::vector<unsigned int>::const_iterator i = v.begin();
11549 for (; line != endl; ++line, ++i)
11550 line->set_user_index(*i);
11554 template <
int dim,
int spacedim>
11557 std::vector<unsigned int> &v)
const 11559 v.resize(n_quads(), 0);
11563 std::vector<unsigned int>::iterator i = v.begin();
11565 for (; quad != endq; ++quad, ++i)
11566 *i = get_user_index(quad);
11572 template <
int dim,
int spacedim>
11575 const std::vector<unsigned int> &v)
11577 Assert(v.size() == n_quads(), ExcGridReadError());
11582 std::vector<unsigned int>::const_iterator i = v.begin();
11583 for (; quad != endq; ++quad, ++i)
11584 set_user_index(quad, *i);
11589 template <
int dim,
int spacedim>
11592 std::vector<unsigned int> &v)
const 11594 v.resize(n_hexs(), 0);
11598 std::vector<unsigned int>::iterator i = v.begin();
11600 for (; hex != endh; ++hex, ++i)
11601 *i = get_user_index(hex);
11607 template <
int dim,
int spacedim>
11610 const std::vector<unsigned int> &v)
11612 Assert(v.size() == n_hexs(), ExcGridReadError());
11617 std::vector<unsigned int>::const_iterator i = v.begin();
11618 for (; hex != endh; ++hex, ++i)
11619 set_user_index(hex, *i);
11630 template <
typename Iterator>
11632 get_user_pointer(
const Iterator &i)
11634 return i->user_pointer();
11639 template <
int structdim,
int dim,
int spacedim>
11650 template <
typename Iterator>
11652 set_user_pointer(
const Iterator &i,
void *x)
11654 i->set_user_pointer(x);
11659 template <
int structdim,
int dim,
int spacedim>
11670 template <
int dim,
int spacedim>
11678 std::vector<void *> tmp;
11680 save_user_pointers_line(tmp);
11681 v.insert(v.end(), tmp.begin(), tmp.end());
11685 save_user_pointers_quad(tmp);
11686 v.insert(v.end(), tmp.begin(), tmp.end());
11691 save_user_pointers_hex(tmp);
11692 v.insert(v.end(), tmp.begin(), tmp.end());
11701 template <
int dim,
int spacedim>
11706 std::vector<void *> tmp;
11710 tmp.insert(tmp.end(), v.begin(), v.begin() + n_lines());
11712 load_user_pointers_line(tmp);
11717 tmp.insert(tmp.end(),
11718 v.begin() + n_lines(),
11719 v.begin() + n_lines() + n_quads());
11720 load_user_pointers_quad(tmp);
11726 tmp.insert(tmp.end(),
11727 v.begin() + n_lines() + n_quads(),
11728 v.begin() + n_lines() + n_quads() + n_hexs());
11729 load_user_pointers_hex(tmp);
11738 template <
int dim,
int spacedim>
11741 std::vector<void *> &v)
const 11743 v.resize(n_lines(),
nullptr);
11744 std::vector<void *>::iterator i = v.begin();
11746 for (; line != endl; ++line, ++i)
11747 *i = line->user_pointer();
11752 template <
int dim,
int spacedim>
11755 const std::vector<void *> &v)
11757 Assert(v.size() == n_lines(), ExcGridReadError());
11760 std::vector<void *>::const_iterator i = v.begin();
11761 for (; line != endl; ++line, ++i)
11762 line->set_user_pointer(*i);
11767 template <
int dim,
int spacedim>
11770 std::vector<void *> &v)
const 11772 v.resize(n_quads(),
nullptr);
11776 std::vector<void *>::iterator i = v.begin();
11778 for (; quad != endq; ++quad, ++i)
11779 *i = get_user_pointer(quad);
11785 template <
int dim,
int spacedim>
11788 const std::vector<void *> &v)
11790 Assert(v.size() == n_quads(), ExcGridReadError());
11795 std::vector<void *>::const_iterator i = v.begin();
11796 for (; quad != endq; ++quad, ++i)
11797 set_user_pointer(quad, *i);
11802 template <
int dim,
int spacedim>
11805 std::vector<void *> &v)
const 11807 v.resize(n_hexs(),
nullptr);
11811 std::vector<void *>::iterator i = v.begin();
11813 for (; hex != endh; ++hex, ++i)
11814 *i = get_user_pointer(hex);
11820 template <
int dim,
int spacedim>
11823 const std::vector<void *> &v)
11825 Assert(v.size() == n_hexs(), ExcGridReadError());
11830 std::vector<void *>::const_iterator i = v.begin();
11831 for (; hex != endh; ++hex, ++i)
11832 set_user_pointer(hex, *i);
11841 template <
int dim,
int spacedim>
11848 return begin_raw_line(level);
11850 return begin_raw_quad(level);
11852 return begin_raw_hex(level);
11861 template <
int dim,
int spacedim>
11868 return begin_line(level);
11870 return begin_quad(level);
11872 return begin_hex(level);
11881 template <
int dim,
int spacedim>
11888 return begin_active_line(level);
11890 return begin_active_quad(level);
11892 return begin_active_hex(level);
11901 template <
int dim,
int spacedim>
11905 const unsigned int level = levels.size() - 1;
11906 if (levels[level]->cells.cells.size() == 0)
11913 levels[level]->cells.cells.size() - 1);
11916 if (ri->used() ==
true)
11919 if (ri->used() ==
true)
11926 template <
int dim,
int spacedim>
11936 if (cell->active() ==
true)
11939 if (cell->active() ==
true)
11947 template <
int dim,
int spacedim>
11958 template <
int dim,
int spacedim>
11970 if (level >= levels.size())
11972 Assert(level < n_global_levels(),
11973 ExcInvalidLevel(level, n_global_levels()));
11979 Assert(level < levels.size(), ExcInvalidLevel(level, levels.size()));
11980 if (level < levels.size() - 1)
11981 return begin_raw(level + 1);
11987 template <
int dim,
int spacedim>
11999 if (level >= levels.size())
12001 Assert(level < n_global_levels(),
12002 ExcInvalidLevel(level, n_global_levels()));
12008 Assert(level < levels.size(), ExcInvalidLevel(level, levels.size()));
12009 if (level < levels.size() - 1)
12010 return begin(level + 1);
12016 template <
int dim,
int spacedim>
12028 if (level >= levels.size())
12030 Assert(level < n_global_levels(),
12031 ExcInvalidLevel(level, n_global_levels()));
12037 Assert(level < levels.size(), ExcInvalidLevel(level, levels.size()));
12039 begin_active(level + 1));
12044 template <
int dim,
int spacedim>
12053 template <
int dim,
int spacedim>
12064 template <
int dim,
int spacedim>
12067 const unsigned int level)
const 12070 begin(level), end(level));
12075 template <
int dim,
int spacedim>
12078 const unsigned int level)
const 12082 begin_active(level), end_active(level));
12089 template <
int dim,
int spacedim>
12099 return begin_line();
12101 return begin_quad();
12110 template <
int dim,
int spacedim>
12120 return begin_active_line();
12122 return begin_active_quad();
12131 template <
int dim,
int spacedim>
12152 template <
int dim,
int spacedim>
12158 begin_active_face(), end_face());
12164 template <
int dim,
int spacedim>
12173 while (i->used() ==
false)
12181 template <
int dim,
int spacedim>
12186 return begin_vertex();
12191 template <
int dim,
int spacedim>
12206 template <
int dim,
int spacedim>
12207 typename Triangulation<dim, spacedim>::raw_line_iterator
12218 if (level >= levels.size())
12220 Assert(level < n_global_levels(),
12221 ExcInvalidLevel(level, n_global_levels()));
12230 Assert(level < levels.size(), ExcInvalidLevel(level, levels.size()));
12232 if (level >= levels.size() || levels[level]->cells.cells.size() == 0)
12235 return raw_line_iterator(
12240 return raw_line_iterator(
12246 template <
int dim,
int spacedim>
12251 raw_line_iterator ri = begin_raw_line(level);
12254 while (ri->used() ==
false)
12262 template <
int dim,
int spacedim>
12270 while (i->has_children())
12278 template <
int dim,
int spacedim>
12292 template <
int dim,
int spacedim>
12293 typename Triangulation<dim, spacedim>::raw_quad_iterator
12304 if (level >= levels.size())
12306 Assert(level < n_global_levels(),
12307 ExcInvalidLevel(level, n_global_levels()));
12315 return raw_hex_iterator();
12320 Assert(level < levels.size(), ExcInvalidLevel(level, levels.size()));
12322 if (level >= levels.size() || levels[level]->cells.cells.size() == 0)
12325 return raw_quad_iterator(
12333 return raw_quad_iterator(
12340 return raw_hex_iterator();
12346 template <
int dim,
int spacedim>
12351 raw_quad_iterator ri = begin_raw_quad(level);
12354 while (ri->used() ==
false)
12362 template <
int dim,
int spacedim>
12370 while (i->has_children())
12378 template <
int dim,
int spacedim>
12391 template <
int dim,
int spacedim>
12392 typename Triangulation<dim, spacedim>::raw_hex_iterator
12403 if (level >= levels.size())
12405 Assert(level < n_global_levels(),
12406 ExcInvalidLevel(level, n_global_levels()));
12415 return raw_hex_iterator();
12420 Assert(level < levels.size(), ExcInvalidLevel(level, levels.size()));
12422 if (level >= levels.size() || levels[level]->cells.cells.size() == 0)
12425 return raw_hex_iterator(
12431 return raw_hex_iterator();
12437 template <
int dim,
int spacedim>
12442 raw_hex_iterator ri = begin_raw_hex(level);
12445 while (ri->used() ==
false)
12453 template <
int dim,
int spacedim>
12461 while (i->has_children())
12469 template <
int dim,
int spacedim>
12485 namespace TriangulationImplementation
12487 inline unsigned int 12494 inline unsigned int 12502 inline unsigned int 12509 inline unsigned int 12517 inline unsigned int 12524 inline unsigned int 12535 template <
int dim,
int spacedim>
12539 return internal::TriangulationImplementation::n_cells(number_cache);
12543 template <
int dim,
int spacedim>
12547 return internal::TriangulationImplementation::n_active_cells(number_cache);
12550 template <
int dim,
int spacedim>
12554 return n_active_cells();
12559 template <
int dim,
int spacedim>
12566 return n_used_vertices();
12578 template <
int dim,
int spacedim>
12585 return n_vertices();
12587 return n_raw_lines();
12589 return n_raw_quads();
12597 template <
int dim,
int spacedim>
12604 return n_used_vertices();
12606 return n_active_lines();
12608 return n_active_quads();
12616 template <
int dim,
int spacedim>
12623 return n_raw_lines(level);
12625 return n_raw_quads(level);
12627 return n_raw_hexs(level);
12636 template <
int dim,
int spacedim>
12643 return n_lines(level);
12645 return n_quads(level);
12647 return n_hexs(level);
12656 template <
int dim,
int spacedim>
12663 return n_active_lines(level);
12665 return n_active_quads(level);
12667 return n_active_hexs(level);
12675 template <
int dim,
int spacedim>
12679 for (
unsigned int lvl = 0; lvl < n_global_levels() - 1; lvl++)
12680 if (n_active_cells(lvl) != 0)
12687 template <
int dim,
int spacedim>
12691 return number_cache.n_lines;
12701 return levels[level]->cells.cells.size();
12720 return levels[level]->cells.cells.size();
12738 return levels[level]->cells.cells.size();
12751 template <
int dim,
int spacedim>
12760 template <
int dim,
int spacedim>
12764 return faces->
lines.cells.size();
12768 template <
int dim,
int spacedim>
12772 Assert(level < number_cache.n_lines_level.size(),
12773 ExcIndexRange(level, 0, number_cache.n_lines_level.size()));
12775 return number_cache.n_lines_level[level];
12779 template <
int dim,
int spacedim>
12783 return number_cache.n_active_lines;
12787 template <
int dim,
int spacedim>
12791 Assert(level < number_cache.n_lines_level.size(),
12792 ExcIndexRange(level, 0, number_cache.n_lines_level.size()));
12795 return number_cache.n_active_lines_level[level];
12945 template <
int dim,
int spacedim>
12949 return number_cache.n_quads;
12953 template <
int dim,
int spacedim>
12958 Assert(level < number_cache.n_quads_level.size(),
12959 ExcIndexRange(level, 0, number_cache.n_quads_level.size()));
12960 return number_cache.n_quads_level[level];
12970 return levels[level]->cells.cells.size();
12980 return levels[level]->cells.cells.size();
12994 template <
int dim,
int spacedim>
13013 template <
int dim,
int spacedim>
13017 return number_cache.n_active_quads;
13021 template <
int dim,
int spacedim>
13025 Assert(level < number_cache.n_quads_level.size(),
13026 ExcIndexRange(level, 0, number_cache.n_quads_level.size()));
13029 return number_cache.n_active_quads_level[level];
13033 template <
int dim,
int spacedim>
13042 template <
int dim,
int spacedim>
13051 template <
int dim,
int spacedim>
13059 template <
int dim,
int spacedim>
13068 template <
int dim,
int spacedim>
13080 return number_cache.n_hexes;
13089 Assert(level < number_cache.n_hexes_level.size(),
13090 ExcIndexRange(level, 0, number_cache.n_hexes_level.size()));
13092 return number_cache.n_hexes_level[level];
13102 return levels[level]->cells.cells.size();
13110 return number_cache.n_active_hexes;
13119 Assert(level < number_cache.n_hexes_level.size(),
13120 ExcIndexRange(level, 0, number_cache.n_hexes_level.size()));
13122 return number_cache.n_active_hexes_level[level];
13127 template <
int dim,
int spacedim>
13131 return std::count_if(vertices_used.begin(),
13132 vertices_used.end(),
13133 std::bind(std::equal_to<bool>(),
13134 std::placeholders::_1,
13140 template <
int dim,
int spacedim>
13141 const std::vector<bool> &
13144 return vertices_used;
13174 template <
int dim,
int spacedim>
13182 unsigned int max_vertex_index = 0;
13183 for (; cell != endc; ++cell)
13184 for (
unsigned int vertex = 0; vertex < GeometryInfo<dim>::vertices_per_cell;
13186 if (cell->vertex_index(vertex) > max_vertex_index)
13187 max_vertex_index = cell->vertex_index(vertex);
13193 std::vector<unsigned short int> usage_count(max_vertex_index + 1, 0);
13197 for (cell = begin(); cell != endc; ++cell)
13198 for (
unsigned int vertex = 0; vertex < GeometryInfo<dim>::vertices_per_cell;
13200 ++usage_count[cell->vertex_index(vertex)];
13203 static_cast<unsigned int>(
13204 *std::max_element(usage_count.begin(), usage_count.end())));
13209 template <
int dim,
int spacedim>
13218 template <
int dim,
int spacedim>
13227 template <
int dim,
int spacedim>
13236 template <
int dim,
int spacedim>
13240 &periodicity_vector)
13242 periodic_face_pairs_level_0.insert(periodic_face_pairs_level_0.end(),
13243 periodicity_vector.begin(),
13244 periodicity_vector.end());
13247 update_periodic_face_map();
13252 template <
int dim,
int spacedim>
13253 const typename std::map<
13254 std::pair<typename Triangulation<dim, spacedim>::cell_iterator,
unsigned int>,
13255 std::pair<std::pair<typename Triangulation<dim, spacedim>::cell_iterator,
13260 return periodic_face_map;
13265 template <
int dim,
int spacedim>
13269 prepare_coarsening_and_refinement();
13274 if (smooth_grid & limit_level_difference_at_vertices)
13278 signals.pre_refinement();
13280 execute_coarsening();
13282 const DistortedCellList cells_with_distorted_children = execute_refinement();
13287 if (smooth_grid & limit_level_difference_at_vertices)
13292 update_neighbors(*
this);
13293 reset_active_cell_indices();
13296 signals.post_refinement();
13298 AssertThrow(cells_with_distorted_children.distorted_cells.size() == 0,
13299 cells_with_distorted_children);
13301 update_periodic_face_map();
13306 template <
int dim,
int spacedim>
13310 unsigned int active_cell_index = 0;
13312 if ((cell->used() ==
false) || cell->has_children())
13316 cell->set_active_cell_index(active_cell_index);
13317 ++active_cell_index;
13324 template <
int dim,
int spacedim>
13329 periodic_face_map.clear();
13331 typename std::vector<
13333 for (it = periodic_face_pairs_level_0.begin();
13334 it != periodic_face_pairs_level_0.end();
13337 update_periodic_face_map_recursively<dim, spacedim>(it->cell[0],
13342 periodic_face_map);
13345 std::bitset<3> inverted_orientation;
13347 bool orientation, flip, rotation;
13348 orientation = it->orientation[0];
13349 rotation = it->orientation[2];
13350 flip = orientation ? rotation ^ it->orientation[1] : it->orientation[1];
13351 inverted_orientation[0] = orientation;
13352 inverted_orientation[1] = flip;
13353 inverted_orientation[2] = rotation;
13355 update_periodic_face_map_recursively<dim, spacedim>(it->cell[1],
13359 inverted_orientation,
13360 periodic_face_map);
13364 typename std::map<std::pair<cell_iterator, unsigned int>,
13365 std::pair<std::pair<cell_iterator, unsigned int>,
13366 std::bitset<3>>>::const_iterator it_test;
13367 for (it_test = periodic_face_map.begin(); it_test != periodic_face_map.end();
13371 it_test->first.first;
13373 it_test->second.first.first;
13374 if (cell_1->level() == cell_2->level())
13378 Assert(periodic_face_map[it_test->second.first].first ==
13387 template <
int dim,
int spacedim>
13395 vertices_used.clear();
13403 template <
int dim,
int spacedim>
13407 const DistortedCellList cells_with_distorted_children =
13409 *
this, check_for_distorted_cells);
13415 *
this, levels.size(), number_cache);
13418 for (
unsigned int level = 0; level < levels.size(); ++level)
13419 levels[level]->cells.monitor_memory(dim);
13426 while (cell != endc)
13430 return cells_with_distorted_children;
13435 template <
int dim,
int spacedim>
13443 std::vector<unsigned int> line_cell_count =
13444 count_cells_bounded_by_line(*
this);
13445 std::vector<unsigned int> quad_cell_count =
13446 count_cells_bounded_by_quad(*
this);
13458 clear_user_flags();
13461 for (; cell != endc; ++cell)
13462 if (!cell->active())
13463 if (cell->child(0)->coarsen_flag_set())
13465 cell->set_user_flag();
13466 for (
unsigned int child = 0; child < cell->n_children(); ++child)
13468 Assert(cell->child(child)->coarsen_flag_set(),
13470 cell->child(child)->clear_coarsen_flag();
13485 if (levels.size() >= 2)
13486 for (cell = last(); cell != endc; --cell)
13487 if (cell->level() <=
static_cast<int>(levels.size() - 2) &&
13488 cell->user_flag_set())
13491 signals.pre_coarsening_on_cell(cell);
13499 *
this, levels.size(), number_cache);
13503 for (cell = begin(); cell != endc; ++cell)
13510 template <
int dim,
int spacedim>
13527 std::vector<bool> previous_coarsen_flags(n_active_cells());
13528 save_coarsen_flags(previous_coarsen_flags);
13530 std::vector<int> vertex_level(vertices.size(), 0);
13532 bool continue_iterating =
true;
13536 if (smooth_grid & limit_level_difference_at_vertices)
13538 Assert(!anisotropic_refinement,
13539 ExcMessage(
"In case of anisotropic refinement the " 13540 "limit_level_difference_at_vertices flag for " 13541 "mesh smoothing must not be set!"));
13545 std::fill(vertex_level.begin(), vertex_level.end(), 0);
13547 for (; cell != endc; ++cell)
13549 if (cell->refine_flag_set())
13550 for (
unsigned int vertex = 0;
13551 vertex < GeometryInfo<dim>::vertices_per_cell;
13553 vertex_level[cell->vertex_index(vertex)] =
13554 std::max(vertex_level[cell->vertex_index(vertex)],
13555 cell->level() + 1);
13556 else if (!cell->coarsen_flag_set())
13557 for (
unsigned int vertex = 0;
13558 vertex < GeometryInfo<dim>::vertices_per_cell;
13560 vertex_level[cell->vertex_index(vertex)] =
13561 std::max(vertex_level[cell->vertex_index(vertex)],
13572 for (
unsigned int vertex = 0;
13573 vertex < GeometryInfo<dim>::vertices_per_cell;
13575 vertex_level[cell->vertex_index(vertex)] =
13576 std::max(vertex_level[cell->vertex_index(vertex)],
13577 cell->level() - 1);
13591 for (cell = last_active(); cell != endc; --cell)
13592 if (cell->refine_flag_set() ==
false)
13594 for (
unsigned int vertex = 0;
13595 vertex < GeometryInfo<dim>::vertices_per_cell;
13597 if (vertex_level[cell->vertex_index(vertex)] >=
13601 cell->clear_coarsen_flag();
13606 if (vertex_level[cell->vertex_index(vertex)] >
13609 cell->set_refine_flag();
13611 for (
unsigned int v = 0;
13612 v < GeometryInfo<dim>::vertices_per_cell;
13614 vertex_level[cell->vertex_index(v)] =
13615 std::max(vertex_level[cell->vertex_index(v)],
13616 cell->level() + 1);
13638 clear_user_flags();
13642 for (; acell != end_ac; ++acell)
13643 acell->clear_coarsen_flag();
13646 for (; cell != endc; ++cell)
13649 if (cell->active())
13652 const unsigned int n_children = cell->n_children();
13653 unsigned int flagged_children = 0;
13654 for (
unsigned int child = 0; child < n_children; ++child)
13655 if (cell->child(child)->active() &&
13656 cell->child(child)->coarsen_flag_set())
13658 ++flagged_children;
13660 cell->child(child)->clear_coarsen_flag();
13665 if (flagged_children == n_children)
13666 cell->set_user_flag();
13672 for (cell = begin(); cell != endc; ++cell)
13693 for (cell = last(); cell != endc; --cell)
13694 if (cell->user_flag_set())
13698 template coarsening_allowed<dim, spacedim>(cell))
13699 for (
unsigned int c = 0; c < cell->n_children(); ++c)
13701 Assert(cell->child(c)->refine_flag_set() ==
false,
13704 cell->child(c)->set_coarsen_flag();
13709 clear_user_flags();
13714 std::vector<bool> current_coarsen_flags(n_active_cells());
13715 save_coarsen_flags(current_coarsen_flags);
13717 continue_iterating = (current_coarsen_flags != previous_coarsen_flags);
13718 previous_coarsen_flags = current_coarsen_flags;
13720 while (continue_iterating ==
true);
13731 std::vector<bool> flags_before;
13732 save_coarsen_flags(flags_before);
13735 fix_coarsen_flags();
13737 std::vector<bool> flags_after;
13738 save_coarsen_flags(flags_after);
13740 return (flags_before != flags_after);
13750 std::vector<bool> flags_before;
13751 save_coarsen_flags(flags_before);
13754 fix_coarsen_flags();
13756 std::vector<bool> flags_after;
13757 save_coarsen_flags(flags_after);
13759 return (flags_before != flags_after);
13769 std::vector<bool> flags_before;
13770 save_coarsen_flags(flags_before);
13773 fix_coarsen_flags();
13775 std::vector<bool> flags_after;
13776 save_coarsen_flags(flags_after);
13778 return (flags_before != flags_after);
13789 template <
int dim,
int spacedim>
13791 possibly_do_not_produce_unrefined_islands(
13796 unsigned int n_neighbors = 0;
13799 unsigned int count = 0;
13800 for (
unsigned int n = 0; n < GeometryInfo<dim>::faces_per_cell; ++n)
13807 if (face_will_be_refined_by_neighbor(cell, n))
13814 if (count == n_neighbors ||
13815 (count >= n_neighbors - 1 &&
13818 for (
unsigned int c = 0; c < cell->n_children(); ++c)
13819 cell->child(c)->clear_coarsen_flag();
13821 for (
unsigned int face = 0; face < GeometryInfo<dim>::faces_per_cell;
13823 if (!cell->at_boundary(face) && (!cell->neighbor(face)->active()) &&
13824 (cell_will_be_coarsened(cell->neighbor(face))))
13825 possibly_do_not_produce_unrefined_islands<dim, spacedim>(
13826 cell->neighbor(face));
13841 template <
int dim,
int spacedim>
13843 possibly_refine_unrefined_island(
13845 const bool allow_anisotropic_smoothing)
13860 if (allow_anisotropic_smoothing ==
false)
13863 unsigned int refined_neighbors = 0, unrefined_neighbors = 0;
13864 for (
unsigned int face = 0; face < GeometryInfo<dim>::faces_per_cell;
13866 if (!cell->at_boundary(face))
13868 if (face_will_be_refined_by_neighbor(cell, face))
13869 ++refined_neighbors;
13871 ++unrefined_neighbors;
13874 if (unrefined_neighbors < refined_neighbors)
13876 cell->clear_coarsen_flag();
13877 cell->set_refine_flag();
13882 if (unrefined_neighbors > 0)
13883 for (
unsigned int face = 0;
13884 face < GeometryInfo<dim>::faces_per_cell;
13886 if (!cell->at_boundary(face) &&
13887 (face_will_be_refined_by_neighbor(cell, face) ==
false) &&
13888 (cell->neighbor(face)->has_children() ==
false) &&
13889 (cell->neighbor(face)->refine_flag_set() ==
false))
13890 possibly_refine_unrefined_island<dim, spacedim>(
13891 cell->neighbor(face), allow_anisotropic_smoothing);
13903 for (
unsigned int face_pair = 0;
13904 face_pair < GeometryInfo<dim>::faces_per_cell / 2;
13913 for (
unsigned int face_index = 0; face_index < 2; ++face_index)
13915 unsigned int face = 2 * face_pair + face_index;
13922 face_will_be_refined_by_neighbor<dim, spacedim>(
13923 cell, face, expected_face_ref_case);
13937 directional_cell_refinement_case =
13938 (directional_cell_refinement_case &
13941 expected_face_ref_case,
13943 cell->face_orientation(face),
13944 cell->face_flip(face),
13945 cell->face_rotation(face)));
13951 Assert(directional_cell_refinement_case <
13954 smoothing_cell_refinement_case =
13955 smoothing_cell_refinement_case | directional_cell_refinement_case;
13960 if (smoothing_cell_refinement_case)
13962 cell->clear_coarsen_flag();
13963 cell->set_refine_flag(cell->refine_flag_set() |
13964 smoothing_cell_refinement_case);
13971 template <
int dim,
int spacedim>
13977 std::vector<bool> flags_before[2];
13978 save_coarsen_flags(flags_before[0]);
13979 save_refine_flags(flags_before[1]);
13997 std::vector<bool> flags_before_loop[2] = {flags_before[0], flags_before[1]};
14074 if (((smooth_grid & coarsest_level_1) || (smooth_grid & patch_level_1)) &&
14077 for (
const auto &cell : active_cell_iterators_on_level(1))
14078 cell->clear_coarsen_flag();
14081 bool mesh_changed_in_this_loop =
false;
14090 if (smooth_grid & do_not_produce_unrefined_islands &&
14091 !(smooth_grid & patch_level_1))
14093 for (
const auto &cell : cell_iterators())
14097 if (!cell->active() && cell_will_be_coarsened(cell))
14098 possibly_do_not_produce_unrefined_islands<dim, spacedim>(cell);
14122 if (smooth_grid & (eliminate_refined_inner_islands |
14123 eliminate_refined_boundary_islands) &&
14124 !(smooth_grid & patch_level_1))
14126 for (
const auto &cell : cell_iterators())
14127 if (!cell->active() || (cell->active() && cell->refine_flag_set() &&
14128 cell->is_locally_owned()))
14135 bool all_children_active =
true;
14136 if (!cell->active())
14137 for (
unsigned int c = 0; c < cell->n_children(); ++c)
14138 if (!cell->child(c)->active() ||
14139 cell->child(c)->is_ghost() ||
14140 cell->child(c)->is_artificial())
14142 all_children_active =
false;
14146 if (all_children_active)
14153 unsigned int unrefined_neighbors = 0, total_neighbors = 0;
14155 for (
unsigned int n = 0;
14156 n < GeometryInfo<dim>::faces_per_cell;
14164 if (!face_will_be_refined_by_neighbor(cell, n))
14165 ++unrefined_neighbors;
14180 if ((unrefined_neighbors == total_neighbors) &&
14181 (((unrefined_neighbors ==
14183 (smooth_grid & eliminate_refined_inner_islands)) ||
14184 ((unrefined_neighbors <
14187 eliminate_refined_boundary_islands))) &&
14188 (total_neighbors != 0))
14190 if (!cell->active())
14191 for (
unsigned int c = 0; c < cell->n_children(); ++c)
14193 cell->child(c)->clear_refine_flag();
14194 cell->child(c)->set_coarsen_flag();
14197 cell->clear_refine_flag();
14212 if (smooth_grid & limit_level_difference_at_vertices)
14214 Assert(!anisotropic_refinement,
14215 ExcMessage(
"In case of anisotropic refinement the " 14216 "limit_level_difference_at_vertices flag for " 14217 "mesh smoothing must not be set!"));
14221 std::vector<int> vertex_level(vertices.size(), 0);
14222 for (
const auto &cell : active_cell_iterators())
14224 if (cell->refine_flag_set())
14225 for (
unsigned int vertex = 0;
14226 vertex < GeometryInfo<dim>::vertices_per_cell;
14228 vertex_level[cell->vertex_index(vertex)] =
14229 std::max(vertex_level[cell->vertex_index(vertex)],
14230 cell->level() + 1);
14231 else if (!cell->coarsen_flag_set())
14232 for (
unsigned int vertex = 0;
14233 vertex < GeometryInfo<dim>::vertices_per_cell;
14235 vertex_level[cell->vertex_index(vertex)] =
14236 std::max(vertex_level[cell->vertex_index(vertex)],
14245 for (
unsigned int vertex = 0;
14246 vertex < GeometryInfo<dim>::vertices_per_cell;
14248 vertex_level[cell->vertex_index(vertex)] =
14249 std::max(vertex_level[cell->vertex_index(vertex)],
14250 cell->level() - 1);
14265 if (cell->refine_flag_set() ==
false)
14267 for (
unsigned int vertex = 0;
14268 vertex < GeometryInfo<dim>::vertices_per_cell;
14270 if (vertex_level[cell->vertex_index(vertex)] >=
14274 cell->clear_coarsen_flag();
14279 if (vertex_level[cell->vertex_index(vertex)] >
14282 cell->set_refine_flag();
14284 for (
unsigned int v = 0;
14285 v < GeometryInfo<dim>::vertices_per_cell;
14287 vertex_level[cell->vertex_index(v)] =
14288 std::max(vertex_level[cell->vertex_index(v)],
14289 cell->level() + 1);
14310 if (smooth_grid & eliminate_unrefined_islands)
14315 if (cell->refine_flag_set() !=
14317 possibly_refine_unrefined_island<dim, spacedim>(
14318 cell, (smooth_grid & allow_anisotropic_smoothing) != 0);
14346 if (smooth_grid & patch_level_1)
14356 for (
const auto &cell : cell_iterators())
14357 if (!cell->active())
14363 if (cell->child(0)->has_children() ==
true)
14370 for (
unsigned int i = 0; i < cell->n_children(); ++i)
14371 combined_ref_case =
14372 combined_ref_case | cell->child(i)->refine_flag_set();
14374 for (
unsigned int i = 0; i < cell->n_children(); ++i)
14378 child->clear_coarsen_flag();
14379 child->set_refine_flag(combined_ref_case);
14392 for (
const auto &cell : cell_iterators())
14398 if (cell->active() || cell->child(0)->active())
14405 const unsigned int n_children = cell->n_children();
14406 bool has_active_grandchildren =
false;
14408 for (
unsigned int i = 0; i < n_children; ++i)
14409 if (cell->child(i)->child(0)->active())
14411 has_active_grandchildren =
true;
14415 if (has_active_grandchildren ==
false)
14421 unsigned int n_grandchildren = 0;
14424 unsigned int n_coarsen_flags = 0;
14429 for (
unsigned int c = 0; c < n_children; ++c)
14436 const unsigned int nn_children = child->n_children();
14437 n_grandchildren += nn_children;
14442 if (child->child(0)->active())
14443 for (
unsigned int cc = 0; cc < nn_children; ++cc)
14444 if (child->child(cc)->coarsen_flag_set())
14459 if ((n_coarsen_flags != n_grandchildren) && (n_coarsen_flags > 0))
14460 for (
unsigned int c = 0; c < n_children; ++c)
14463 if (child->child(0)->active())
14464 for (
unsigned int cc = 0; cc < child->n_children(); ++cc)
14465 child->child(cc)->clear_coarsen_flag();
14493 bool changed =
true;
14499 for (; cell != endc; --cell)
14500 if (cell->refine_flag_set())
14503 for (
unsigned int i = 0; i < GeometryInfo<dim>::faces_per_cell;
14509 const bool has_periodic_neighbor =
14510 cell->has_periodic_neighbor(i);
14511 const bool has_neighbor_or_periodic_neighbor =
14512 !cell->at_boundary(i) || has_periodic_neighbor;
14513 if (has_neighbor_or_periodic_neighbor &&
14515 cell->refine_flag_set(), i) !=
14528 if (cell->neighbor_or_periodic_neighbor(i)->active())
14530 if ((!has_periodic_neighbor &&
14531 cell->neighbor_is_coarser(i)) ||
14532 (has_periodic_neighbor &&
14533 cell->periodic_neighbor_is_coarser(i)))
14535 if (cell->neighbor_or_periodic_neighbor(i)
14536 ->coarsen_flag_set())
14537 cell->neighbor_or_periodic_neighbor(i)
14538 ->clear_coarsen_flag();
14555 allow_anisotropic_smoothing)
14557 has_periodic_neighbor ?
14558 cell->periodic_neighbor(i)
14559 ->flag_for_face_refinement(
14561 ->periodic_neighbor_of_coarser_periodic_neighbor(
14566 ->flag_for_face_refinement(
14568 ->neighbor_of_coarser_neighbor(
14575 ->neighbor_or_periodic_neighbor(
14577 ->refine_flag_set())
14579 cell->neighbor_or_periodic_neighbor(i)
14580 ->set_refine_flag();
14712 std::pair<unsigned int, unsigned int>
14714 has_periodic_neighbor ?
14716 ->periodic_neighbor_of_coarser_periodic_neighbor(
14718 cell->neighbor_of_coarser_neighbor(i);
14719 unsigned int refined_along_x = 0,
14720 refined_along_y = 0,
14721 to_be_refined_along_x = 0,
14722 to_be_refined_along_y = 0;
14724 const int this_face_index =
14725 cell->face_index(i);
14734 const auto parent_face = [&]() {
14735 if (has_periodic_neighbor)
14737 const auto neighbor =
14738 cell->periodic_neighbor(i);
14739 const auto parent_face_no =
14741 ->periodic_neighbor_of_periodic_neighbor(
14744 neighbor->periodic_neighbor(
14746 return parent->face(parent_face_no);
14749 return cell->neighbor(i)->face(
14753 if ((this_face_index ==
14754 parent_face->child_index(0)) ||
14755 (this_face_index ==
14756 parent_face->child_index(1)))
14765 parent_face->refinement_case();
14783 cell->refine_flag_set(),
14785 cell->face_orientation(i),
14786 cell->face_flip(i),
14787 cell->face_rotation(i));
14790 ++to_be_refined_along_x;
14793 ++to_be_refined_along_y;
14798 allow_anisotropic_smoothing) ||
14799 cell->neighbor_or_periodic_neighbor(i)
14800 ->refine_flag_set())
14802 if (refined_along_x +
14803 to_be_refined_along_x >
14807 ->neighbor_or_periodic_neighbor(i)
14808 ->flag_for_face_refinement(
14812 if (refined_along_y +
14813 to_be_refined_along_y >
14817 ->neighbor_or_periodic_neighbor(i)
14818 ->flag_for_face_refinement(
14826 ->neighbor_or_periodic_neighbor(i)
14827 ->refine_flag_set() !=
14829 dim>::isotropic_refinement)
14831 cell->neighbor_or_periodic_neighbor(i)
14832 ->set_refine_flag();
14838 cell->neighbor_or_periodic_neighbor(i);
14841 nb->refine_flag_set(),
14843 nb->face_orientation(nb_indices.first),
14844 nb->face_flip(nb_indices.first),
14845 nb->face_rotation(nb_indices.first));
14847 !(refined_along_x ||
14848 to_be_refined_along_x))
14849 changed |= cell->flag_for_face_refinement(
14853 !(refined_along_y ||
14854 to_be_refined_along_y))
14855 changed |= cell->flag_for_face_refinement(
14862 cell->neighbor_or_periodic_neighbor(i)
14863 ->clear_coarsen_flag();
14864 const unsigned int nb_nb =
14865 has_periodic_neighbor ?
14867 ->periodic_neighbor_of_periodic_neighbor(
14869 cell->neighbor_of_neighbor(i);
14871 cell->neighbor_or_periodic_neighbor(i);
14874 neighbor->refine_flag_set(),
14876 neighbor->face_orientation(nb_nb),
14877 neighbor->face_flip(nb_nb),
14878 neighbor->face_rotation(nb_nb));
14881 cell->refine_flag_set(),
14883 cell->face_orientation(i),
14884 cell->face_flip(i),
14885 cell->face_rotation(i));
14890 if ((face_ref_case ==
14892 needed_face_ref_case ==
14896 needed_face_ref_case ==
14899 changed = cell->flag_for_face_refinement(
14901 neighbor->flag_for_face_refinement(
14902 nb_nb, needed_face_ref_case);
14909 face_ref_case = cell->face(i)->refinement_case(),
14910 needed_face_ref_case =
14912 cell->refine_flag_set(),
14914 cell->face_orientation(i),
14915 cell->face_flip(i),
14916 cell->face_rotation(i));
14921 needed_face_ref_case ==
14924 needed_face_ref_case ==
14927 cell->flag_for_face_refinement(i,
14948 fix_coarsen_flags();
14951 std::vector<bool> flags_after_loop[2];
14952 save_coarsen_flags(flags_after_loop[0]);
14953 save_refine_flags(flags_after_loop[1]);
14957 mesh_changed_in_this_loop =
14958 ((flags_before_loop[0] != flags_after_loop[0]) ||
14959 (flags_before_loop[1] != flags_after_loop[1]));
14963 flags_before_loop[0].swap(flags_after_loop[0]);
14964 flags_before_loop[1].swap(flags_after_loop[1]);
14966 while (mesh_changed_in_this_loop);
14972 return ((flags_before[0] != flags_before_loop[0]) ||
14973 (flags_before[1] != flags_before_loop[1]));
14978 template <
int dim,
int spacedim>
14981 const unsigned int magic_number1,
14982 const std::vector<bool> &v,
14983 const unsigned int magic_number2,
14984 std::ostream & out)
14986 const unsigned int N = v.size();
14987 unsigned char * flags =
new unsigned char[N / 8 + 1];
14988 for (
unsigned int i = 0; i < N / 8 + 1; ++i)
14991 for (
unsigned int position = 0; position < N; ++position)
14992 flags[position / 8] |= (v[position] ? (1 << (position % 8)) : 0);
15001 out << magic_number1 <<
' ' << N << std::endl;
15002 for (
unsigned int i = 0; i < N / 8 + 1; ++i)
15003 out << static_cast<unsigned int>(flags[i]) <<
' ';
15005 out << std::endl << magic_number2 << std::endl;
15013 template <
int dim,
int spacedim>
15016 std::vector<bool> &v,
15017 const unsigned int magic_number2,
15022 unsigned int magic_number;
15023 in >> magic_number;
15024 AssertThrow(magic_number == magic_number1, ExcGridReadError());
15030 unsigned char * flags =
new unsigned char[N / 8 + 1];
15031 unsigned short int tmp;
15032 for (
unsigned int i = 0; i < N / 8 + 1; ++i)
15038 for (
unsigned int position = 0; position != N; ++position)
15039 v[position] = (flags[position / 8] & (1 << (position % 8)));
15041 in >> magic_number;
15042 AssertThrow(magic_number == magic_number2, ExcGridReadError());
15051 template <
int dim,
int spacedim>
15055 std::size_t mem = 0;
15057 for (
unsigned int i = 0; i < levels.size(); ++i)
15061 mem +=
sizeof(manifold);
15062 mem +=
sizeof(smooth_grid);
15064 mem +=
sizeof(faces);
15073 #include "tria.inst" 15075 DEAL_II_NAMESPACE_CLOSE
std::vector< CellData< 1 > > boundary_lines
static ::ExceptionBase & ExcInconsistentQuadInfoOfQuad(int arg1, int arg2, int arg3, int arg4, std::string arg5)
::internal::TriangulationImplementation::NumberCache< dim > number_cache
unsigned int n_active_cells() const
boost::signals2::signal< void()> post_refinement
virtual void copy_triangulation(const Triangulation< dim, spacedim > &other_tria)
void loop(ITERATOR begin, typename identity< ITERATOR >::type end, DOFINFO &dinfo, INFOBOX &info, const std::function< void(DOFINFO &, typename INFOBOX::CellInfo &)> &cell_worker, const std::function< void(DOFINFO &, typename INFOBOX::CellInfo &)> &boundary_worker, const std::function< void(DOFINFO &, DOFINFO &, typename INFOBOX::CellInfo &, typename INFOBOX::CellInfo &)> &face_worker, ASSEMBLER &assembler, const LoopControl &lctrl=LoopControl())
unsigned int n_used_vertices() const
std::vector< unsigned int > n_quads_level
const types::manifold_id flat_manifold_id
static const unsigned int invalid_unsigned_int
void load_user_flags_line(std::istream &in)
#define AssertNothrow(cond, exc)
#define DeclException2(Exception2, type1, type2, outsequence)
static ::ExceptionBase & ExcInternalErrorOnCell(int arg1)
hex_iterator begin_hex(const unsigned int level=0) const
const types::subdomain_id invalid_subdomain_id
active_face_iterator begin_active_face() const
typename IteratorSelector::line_iterator line_iterator
virtual void create_triangulation_compatibility(const std::vector< Point< spacedim >> &vertices, const std::vector< CellData< dim >> &cells, const SubCellData &subcelldata)
void save_user_flags_quad(std::ostream &out) const
virtual bool has_hanging_nodes() const
unsigned int n_raw_cells(const unsigned int level) const
cell_iterator last() const
void save_user_pointers_hex(std::vector< void *> &v) const
vertex_iterator begin_vertex() const
static void prepare_refinement_dim_dependent(const Triangulation< dim, spacedim > &)
static RefinementCase< 1 > line_refinement_case(const RefinementCase< dim > &cell_refinement_case, const unsigned int line_no)
unsigned int n_active_lines
unsigned int n_cells() const
static ::ExceptionBase & ExcIO()
Subscriptor & operator=(const Subscriptor &)
Task< RT > new_task(const std::function< RT()> &function)
std::vector< Point< spacedim > > vertices
std::map< types::manifold_id, std::unique_ptr< const Manifold< dim, spacedim > > > manifold
typename IteratorSelector::hex_iterator hex_iterator
static unsigned int line_to_cell_vertices(const unsigned int line, const unsigned int vertex)
IteratorRange< active_cell_iterator > active_cell_iterators() const
void set_all_refine_flags()
void load_user_flags(std::istream &in)
typename IteratorSelector::active_line_iterator active_line_iterator
virtual std::unique_ptr< Manifold< dim, spacedim > > clone() const =0
bool anisotropic_refinement
int face(const unsigned int i) const
static void read_bool_vector(const unsigned int magic_number1, std::vector< bool > &v, const unsigned int magic_number2, std::istream &in)
line_iterator begin_line(const unsigned int level=0) const
line_iterator end_line() const
void clear_user_flags_quad()
static void create_children(Triangulation< 2, spacedim > &triangulation, unsigned int &next_unused_vertex, typename Triangulation< 2, spacedim >::raw_line_iterator &next_unused_line, typename Triangulation< 2, spacedim >::raw_cell_iterator &next_unused_cell, typename Triangulation< 2, spacedim >::cell_iterator &cell)
unsigned int n_faces() const
static ::ExceptionBase & ExcCellHasNegativeMeasure(int arg1)
void flip_all_direction_flags()
active_cell_iterator begin_active(const unsigned int level=0) const
#define AssertThrow(cond, exc)
static ::ExceptionBase & ExcInvalidVertexIndex(int arg1, int arg2, int arg3)
static void compute_number_cache(const Triangulation< dim, spacedim > &triangulation, const unsigned int level_objects, internal::TriangulationImplementation::NumberCache< 1 > &number_cache)
static void write_bool_vector(const unsigned int magic_number1, const std::vector< bool > &v, const unsigned int magic_number2, std::ostream &out)
unsigned int n_active_hexs() const
void save_user_pointers_quad(std::vector< void *> &v) const
static ::ExceptionBase & ExcIndexRange(int arg1, int arg2, int arg3)
static ::ExceptionBase & ExcFacesHaveNoLevel()
void execute_coarsening()
static Triangulation< 2, spacedim >::DistortedCellList execute_refinement(Triangulation< 2, spacedim > &triangulation, const bool check_for_distorted_cells)
cell_iterator begin(const unsigned int level=0) const
void load_user_flags_quad(std::istream &in)
#define DeclException5( Exception5, type1, type2, type3, type4, type5, outsequence)
void set_manifold(const types::manifold_id number, const Manifold< dim, spacedim > &manifold_object)
raw_hex_iterator begin_raw_hex(const unsigned int level=0) const
void load_user_pointers_line(const std::vector< void *> &v)
unsigned int n_active_faces() const
raw_quad_iterator begin_raw_quad(const unsigned int level=0) const
std::unique_ptr<::internal::TriangulationImplementation::TriaFaces< dim > > faces
static ::ExceptionBase & ExcInteriorQuadCantBeBoundary(int arg1, int arg2, int arg3, int arg4, types::boundary_id arg5)
cell_iterator end() const
static ::ExceptionBase & ExcInconsistentLineInfoOfLine(int arg1, int arg2, std::string arg3)
virtual void add_periodicity(const std::vector< GridTools::PeriodicFacePair< cell_iterator >> &)
static Triangulation< 3, spacedim >::DistortedCellList execute_refinement(Triangulation< 3, spacedim > &triangulation, const bool check_for_distorted_cells)
void load_coarsen_flags(std::istream &out)
void load_user_flags_hex(std::istream &in)
active_quad_iterator begin_active_quad(const unsigned int level=0) const
virtual void execute_coarsening_and_refinement()
static ::ExceptionBase & ExcGridHasInvalidCell(int arg1)
unsigned int n_active_quads
unsigned int n_active_lines() const
IteratorRange< cell_iterator > cell_iterators_on_level(const unsigned int level) const
virtual bool prepare_coarsening_and_refinement()
static ::ExceptionBase & ExcMessage(std::string arg1)
bool check_consistency(const unsigned int dim) const
static ::ExceptionBase & ExcImpossibleInDim(int arg1)
void reset_active_cell_indices()
DistortedCellList execute_refinement()
unsigned int subdomain_id
void save_user_flags_line(std::ostream &out) const
void save_user_pointers_line(std::vector< void *> &v) const
Triangulation(const MeshSmoothing smooth_grid=none, const bool check_for_distorted_cells=false)
#define DeclException1(Exception1, type1, outsequence)
unsigned int n_active_hexes
void clear_user_flags_hex()
void load_user_indices_line(const std::vector< unsigned int > &v)
vertex_iterator end_vertex() const
unsigned int n_raw_quads() const
virtual void create_triangulation(const std::vector< Point< spacedim >> &vertices, const std::vector< CellData< dim >> &cells, const SubCellData &subcelldata)
#define Assert(cond, exc)
active_hex_iterator begin_active_hex(const unsigned int level=0) const
void save_user_indices_line(std::vector< unsigned int > &v) const
static void prevent_distorted_boundary_cells(const Triangulation< 1, spacedim > &)
static RefinementCase< dim - 1 > face_refinement_case(const RefinementCase< dim > &cell_refinement_case, const unsigned int face_no, const bool face_orientation=true, const bool face_flip=false, const bool face_rotation=false)
boost::signals2::signal< void()> clear
void set_all_manifold_ids(const types::manifold_id number)
static ::ExceptionBase & ExcInteriorLineCantBeBoundary(int arg1, int arg2, types::boundary_id arg3)
unsigned int n_raw_hexs(const unsigned int level) const
unsigned int max_adjacent_cells() const
bool get_anisotropic_refinement_flag() const
std::list< typename Triangulation< dim, spacedim >::cell_iterator > distorted_cells
unsigned int n_quads() const
std::unique_ptr< std::map< unsigned int, types::manifold_id > > vertex_to_manifold_id_map_1d
std::vector< bool > vertices_used
static void alternating_form_at_vertices(const Point< spacedim >(&vertices)[vertices_per_cell], Tensor< spacedim - dim, spacedim >(&forms)[vertices_per_cell])
void save_coarsen_flags(std::ostream &out) const
void reset_all_manifolds()
typename IteratorSelector::active_quad_iterator active_quad_iterator
void save_user_indices(std::vector< unsigned int > &v) const
boost::signals2::signal< void()> mesh_movement
void load_user_indices(const std::vector< unsigned int > &v)
void save_refine_flags(std::ostream &out) const
unsigned int n_lines() const
static ::ExceptionBase & ExcMultiplySetLineInfoOfLine(int arg1, int arg2)
void clear_despite_subscriptions()
raw_cell_iterator begin_raw(const unsigned int level=0) const
unsigned int n_raw_lines() const
virtual std::size_t memory_consumption() const
IteratorRange< active_face_iterator > active_face_iterators() const
typename IteratorSelector::quad_iterator quad_iterator
virtual void set_mesh_smoothing(const MeshSmoothing mesh_smoothing)
void load_user_pointers_hex(const std::vector< void *> &v)
static void create_triangulation(const std::vector< Point< spacedim >> &v, const std::vector< CellData< 2 >> &cells, const SubCellData &subcelldata, Triangulation< 2, spacedim > &triangulation)
void save_user_flags_hex(std::ostream &out) const
face_iterator begin_face() const
boost::signals2::signal< void()> create
std::unique_ptr< std::map< unsigned int, types::boundary_id > > vertex_to_boundary_id_map_1d
virtual types::global_dof_index n_global_active_cells() const
unsigned int n_hexs() const
std::vector< std::unique_ptr< ::internal::TriangulationImplementation::TriaLevel< dim > > > levels
unsigned int n_active_quads() const
active_line_iterator begin_active_line(const unsigned int level=0) const
raw_cell_iterator end_raw(const unsigned int level) const
unsigned int n_raw_faces() const
std::vector< unsigned int > n_lines_level
void load_user_indices_hex(const std::vector< unsigned int > &v)
static ::ExceptionBase & ExcLineInexistant(int arg1, int arg2)
static unsigned int n_children(const RefinementCase< dim > &refinement_case)
typename IteratorSelector::active_hex_iterator active_hex_iterator
void swap(Vector< Number > &u, Vector< Number > &v)
active_vertex_iterator begin_active_vertex() const
void reset_manifold(const types::manifold_id manifold_number)
unsigned int global_dof_index
static bool coarsening_allowed(const typename Triangulation< dim, spacedim >::cell_iterator &cell)
TriaObjects< TriaObject< 1 > > lines
void load_user_pointers_quad(const std::vector< void *> &v)
static void create_triangulation(const std::vector< Point< spacedim >> &v, const std::vector< CellData< 3 >> &cells, const SubCellData &subcelldata, Triangulation< 3, spacedim > &triangulation)
TriaObjects< TriaObject< 1 > > lines
static void create_triangulation(const std::vector< Point< spacedim >> &v, const std::vector< CellData< 1 >> &cells, const SubCellData &, Triangulation< 1, spacedim > &triangulation)
active_cell_iterator end_active(const unsigned int level) const
void update_periodic_face_map()
std::vector< unsigned int > n_active_hexes_level
static void compute_number_cache(const Triangulation< dim, spacedim > &triangulation, const unsigned int level_objects, internal::TriangulationImplementation::NumberCache< 3 > &number_cache)
void save_user_indices_quad(std::vector< unsigned int > &v) const
static unsigned int child_cell_on_face(const RefinementCase< dim > &ref_case, const unsigned int face, const unsigned int subface, const bool face_orientation=true, const bool face_flip=false, const bool face_rotation=false, const RefinementCase< dim - 1 > &face_refinement_case=RefinementCase< dim - 1 >::isotropic_refinement)
const std::vector< bool > & get_used_vertices() const
hex_iterator end_hex() const
std::vector< CellData< 2 > > boundary_quads
MeshSmoothing smooth_grid
virtual ~Triangulation() override
#define DeclException4(Exception4, type1, type2, type3, type4, outsequence)
IteratorRange< active_cell_iterator > active_cell_iterators_on_level(const unsigned int level) const
void load_user_pointers(const std::vector< void *> &v)
std::vector< unsigned int > n_active_lines_level
void refine_global(const unsigned int times=1)
void save_user_pointers(std::vector< void *> &v) const
void save_user_flags(std::ostream &out) const
static RefinementCase< dim > min_cell_refinement_case_for_face_refinement(const RefinementCase< dim - 1 > &face_refinement_case, const unsigned int face_no, const bool face_orientation=true, const bool face_flip=false, const bool face_rotation=false)
static ::ExceptionBase & ExcNotImplemented()
void clear_user_data(const unsigned int i)
Iterator points to a valid object.
quad_iterator end_quad() const
IteratorRange< cell_iterator > cell_iterators() const
#define DeclException3(Exception3, type1, type2, type3, outsequence)
face_iterator end_face() const
void set_all_manifold_ids_on_boundary(const types::manifold_id number)
static ::ExceptionBase & ExcQuadInexistant(int arg1, int arg2, int arg3, int arg4)
active_cell_iterator last_active() const
const types::boundary_id internal_face_boundary_id
std::vector< unsigned int > n_active_quads_level
void clear_user_flags_line()
IteratorState::IteratorStates state() const
quad_iterator begin_quad(const unsigned int level=0) const
static Triangulation< 1, spacedim >::DistortedCellList execute_refinement(Triangulation< 1, spacedim > &triangulation, const bool)
void set_face(const unsigned int i, const int index)
Triangulation & operator=(Triangulation< dim, spacedim > &&tria) noexcept
std::vector< unsigned int > n_hexes_level
std::vector< types::manifold_id > get_manifold_ids() const
static void delete_children(Triangulation< 1, spacedim > &triangulation, typename Triangulation< 1, spacedim >::cell_iterator &cell, std::vector< unsigned int > &, std::vector< unsigned int > &)
const Manifold< dim, spacedim > & get_manifold(const types::manifold_id number) const
void load_user_indices_quad(const std::vector< unsigned int > &v)
virtual types::subdomain_id locally_owned_subdomain() const
std::vector< types::boundary_id > get_boundary_ids() const
boost::signals2::signal< void()> any_change
raw_line_iterator begin_raw_line(const unsigned int level=0) const
virtual const MeshSmoothing & get_mesh_smoothing() const
std::enable_if< std::is_fundamental< T >::value, std::size_t >::type memory_consumption(const T &t)
const std::map< std::pair< cell_iterator, unsigned int >, std::pair< std::pair< cell_iterator, unsigned int >, std::bitset< 3 > > > & get_periodic_face_map() const
void load_refine_flags(std::istream &in)
Tensor< 2, dim, Number > l(const Tensor< 2, dim, Number > &F, const Tensor< 2, dim, Number > &dF_dt)
static ::ExceptionBase & ExcInternalError()
Triangulation< dim, spacedim > & get_triangulation()
static void compute_number_cache(const Triangulation< dim, spacedim > &triangulation, const unsigned int level_objects, internal::TriangulationImplementation::NumberCache< 2 > &number_cache)
void save_user_indices_hex(std::vector< unsigned int > &v) const