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> 20 #include <deal.II/base/table.h> 22 #include <deal.II/fe/mapping_q1.h> 24 #include <deal.II/grid/grid_tools.h> 25 #include <deal.II/grid/magic_numbers.h> 26 #include <deal.II/grid/manifold.h> 27 #include <deal.II/grid/tria.h> 28 #include <deal.II/grid/tria_accessor.h> 29 #include <deal.II/grid/tria_faces.h> 30 #include <deal.II/grid/tria_iterator.h> 31 #include <deal.II/grid/tria_levels.h> 33 #include <deal.II/lac/full_matrix.h> 34 #include <deal.II/lac/vector.h> 45 DEAL_II_NAMESPACE_OPEN
63 namespace TriangulationImplementation
137 template <
int dim,
int spacedim>
139 cell_is_patch_level_1(
144 unsigned int n_active_children = 0;
145 for (
unsigned int i = 0; i < cell->n_children(); ++i)
146 if (cell->child(i)->active())
149 return (n_active_children == 0) ||
150 (n_active_children == cell->n_children());
160 template <
int dim,
int spacedim>
162 cell_will_be_coarsened(
168 if (cell->has_children())
170 unsigned int children_to_coarsen = 0;
171 const unsigned int n_children = cell->n_children();
173 for (
unsigned int c = 0; c < n_children; ++c)
174 if (cell->child(c)->active() && cell->child(c)->coarsen_flag_set())
175 ++children_to_coarsen;
176 if (children_to_coarsen == n_children)
179 for (
unsigned int c = 0; c < n_children; ++c)
180 if (cell->child(c)->active())
181 cell->child(c)->clear_coarsen_flag();
217 template <
int dim,
int spacedim>
219 face_will_be_refined_by_neighbor_internal(
221 const unsigned int face_no,
230 cell->neighbor(face_no);
237 if (neighbor->has_children())
242 if (cell_will_be_coarsened(neighbor))
251 expected_face_ref_case = cell->face(face_no)->refinement_case();
263 const unsigned int neighbor_neighbor = cell->neighbor_face_no(face_no);
271 neighbor->face_orientation(neighbor_neighbor),
272 neighbor->face_flip(neighbor_neighbor),
273 neighbor->face_rotation(neighbor_neighbor));
277 neighbor_face = neighbor->face(neighbor_neighbor);
278 const int this_face_index = cell->face_index(face_no);
284 if (neighbor_face->index() == this_face_index)
290 expected_face_ref_case = face_ref_case;
308 for (
unsigned int c = 0; c < neighbor_face->n_children(); ++c)
309 if (neighbor_face->child_index(c) == this_face_index)
318 if ((neighbor_face->refinement_case() | face_ref_case) ==
319 neighbor_face->refinement_case())
338 expected_face_ref_case =
339 ~neighbor_face->refinement_case();
366 template <
int dim,
int spacedim>
368 face_will_be_refined_by_neighbor(
370 const unsigned int face_no)
373 return face_will_be_refined_by_neighbor_internal(cell, face_no, dummy);
380 template <
int dim,
int spacedim>
382 face_will_be_refined_by_neighbor(
384 const unsigned int face_no,
387 return face_will_be_refined_by_neighbor_internal(cell,
389 expected_face_ref_case);
394 template <
int dim,
int spacedim>
396 satisfies_level1_at_vertex_rule(
399 std::vector<unsigned int> min_adjacent_cell_level(
401 std::vector<unsigned int> max_adjacent_cell_level(
406 cell != triangulation.
end();
408 for (
unsigned int v = 0; v < GeometryInfo<dim>::vertices_per_cell; ++v)
410 min_adjacent_cell_level[cell->vertex_index(v)] =
411 std::min<unsigned int>(
412 min_adjacent_cell_level[cell->vertex_index(v)], cell->level());
413 max_adjacent_cell_level[cell->vertex_index(v)] =
414 std::max<unsigned int>(
415 min_adjacent_cell_level[cell->vertex_index(v)], cell->level());
418 for (
unsigned int k = 0; k < triangulation.
n_vertices(); ++k)
420 if (max_adjacent_cell_level[k] - min_adjacent_cell_level[k] > 1)
433 template <
int dim,
int spacedim>
434 std::vector<unsigned int>
439 std::vector<unsigned int> line_cell_count(triangulation.
n_raw_lines(),
442 cell = triangulation.
begin(),
443 endc = triangulation.
end();
444 for (; cell != endc; ++cell)
445 for (
unsigned int l = 0; l < GeometryInfo<dim>::lines_per_cell; ++l)
446 ++line_cell_count[cell->line_index(l)];
447 return line_cell_count;
450 return std::vector<unsigned int>();
461 template <
int dim,
int spacedim>
462 std::vector<unsigned int>
467 std::vector<unsigned int> quad_cell_count(triangulation.
n_raw_quads(),
470 cell = triangulation.
begin(),
471 endc = triangulation.
end();
472 for (; cell != endc; ++cell)
473 for (
unsigned int q = 0; q < GeometryInfo<dim>::faces_per_cell; ++q)
474 ++quad_cell_count[cell->quad_index(q)];
475 return quad_cell_count;
478 return std::vector<unsigned int>();
502 void reorder_compatibility(std::vector<
CellData<2>> &cells,
505 for (
unsigned int cell = 0; cell < cells.size(); ++cell)
506 std::swap(cells[cell].vertices[2], cells[cell].vertices[3]);
510 void reorder_compatibility(std::vector<
CellData<3>> &cells,
514 for (
unsigned int cell = 0; cell < cells.size(); ++cell)
516 for (
unsigned int i = 0; i < GeometryInfo<3>::vertices_per_cell; ++i)
517 tmp[i] = cells[cell].vertices[i];
518 for (
unsigned int i = 0; i < GeometryInfo<3>::vertices_per_cell; ++i)
523 std::vector<CellData<2>>::iterator boundary_quad =
525 std::vector<CellData<2>>::iterator end_quad =
527 for (
unsigned int quad_no = 0; boundary_quad != end_quad;
528 ++boundary_quad, ++quad_no)
529 std::swap(boundary_quad->vertices[2], boundary_quad->vertices[3]);
551 template <
int dim,
int spacedim>
554 const typename Triangulation<dim, spacedim>::line_iterator &line)
556 if (line->has_children())
557 return line->child(0)->vertex_index(1);
562 template <
int dim,
int spacedim>
565 const typename Triangulation<dim, spacedim>::quad_iterator &quad)
567 switch (static_cast<unsigned char>(quad->refinement_case()))
570 return middle_vertex_index<dim, spacedim>(quad->child(0)->line(1));
573 return middle_vertex_index<dim, spacedim>(quad->child(0)->line(3));
576 return quad->child(0)->vertex_index(3);
585 template <
int dim,
int spacedim>
588 const typename Triangulation<dim, spacedim>::hex_iterator &hex)
590 switch (static_cast<unsigned char>(hex->refinement_case()))
593 return middle_vertex_index<dim, spacedim>(hex->child(0)->quad(1));
596 return middle_vertex_index<dim, spacedim>(hex->child(0)->quad(3));
599 return middle_vertex_index<dim, spacedim>(hex->child(0)->quad(5));
602 return middle_vertex_index<dim, spacedim>(hex->child(0)->line(11));
605 return middle_vertex_index<dim, spacedim>(hex->child(0)->line(5));
608 return middle_vertex_index<dim, spacedim>(hex->child(0)->line(7));
611 return hex->child(0)->vertex_index(7);
632 template <
class TRIANGULATION>
633 inline typename TRIANGULATION::DistortedCellList
634 collect_distorted_coarse_cells(
const TRIANGULATION &)
636 return typename TRIANGULATION::DistortedCellList();
655 triangulation.
begin(0);
656 cell != triangulation.
end(0);
660 for (
unsigned int i = 0; i < GeometryInfo<dim>::vertices_per_cell; ++i)
661 vertices[i] = cell->vertex(i);
666 for (
unsigned int i = 0; i < GeometryInfo<dim>::vertices_per_cell; ++i)
667 if (determinants[i] <= 1e-9 * std::pow(cell->diameter(), 1. * dim))
669 distorted_cells.distorted_cells.push_back(cell);
674 return distorted_cells;
686 has_distorted_children(
688 std::integral_constant<int, dim>,
689 std::integral_constant<int, dim>)
693 for (
unsigned int c = 0; c < cell->n_children(); ++c)
696 for (
unsigned int i = 0; i < GeometryInfo<dim>::vertices_per_cell; ++i)
697 vertices[i] = cell->child(c)->vertex(i);
702 for (
unsigned int i = 0; i < GeometryInfo<dim>::vertices_per_cell; ++i)
703 if (determinants[i] <=
704 1e-9 * std::pow(cell->child(c)->diameter(), 1. * dim))
719 template <
int dim,
int spacedim>
721 has_distorted_children(
723 std::integral_constant<int, dim>,
724 std::integral_constant<int, spacedim>)
735 template <
int spacedim>
740 template <
int dim,
int spacedim>
795 static const unsigned int left_right_offset[2][6][2] = {
804 {{0, 1}, {1, 0}, {0, 1}, {1, 0}, {0, 1}, {1, 0}}};
815 std::vector<typename Triangulation<dim, spacedim>::cell_iterator>
816 adjacent_cells(2 * triangulation.
n_raw_faces(), dummy);
822 for (; cell != endc; ++cell)
823 for (
unsigned int f = 0; f < GeometryInfo<dim>::faces_per_cell; ++f)
828 const unsigned int offset =
829 (cell->direction_flag() ?
830 left_right_offset[dim - 2][f][cell->face_orientation(f)] :
831 1 - left_right_offset[dim - 2][f][cell->face_orientation(f)]);
833 adjacent_cells[2 * face->index() + offset] = cell;
843 if (cell->active() && face->has_children())
845 adjacent_cells[2 * face->child(0)->index() + offset] = cell;
846 adjacent_cells[2 * face->child(1)->index() + offset] = cell;
871 if (face->has_children() &&
874 cell->refinement_case(), f) ==
877 for (
unsigned int c = 0; c < face->n_children(); ++c)
878 adjacent_cells[2 * face->child(c)->index() + offset] = cell;
879 if (face->child(0)->has_children())
881 adjacent_cells[2 * face->child(0)->child(0)->index() +
883 adjacent_cells[2 * face->child(0)->child(1)->index() +
886 if (face->child(1)->has_children())
888 adjacent_cells[2 * face->child(1)->child(0)->index() +
890 adjacent_cells[2 * face->child(1)->child(1)->index() +
902 for (cell = triangulation.
begin(); cell != endc; ++cell)
903 for (
unsigned int f = 0; f < GeometryInfo<dim>::faces_per_cell; ++f)
905 const unsigned int offset =
906 (cell->direction_flag() ?
907 left_right_offset[dim - 2][f][cell->face_orientation(f)] :
908 1 - left_right_offset[dim - 2][f][cell->face_orientation(f)]);
910 f, adjacent_cells[2 * cell->face(f)->index() + 1 - offset]);
915 template <
int dim,
int spacedim>
917 update_periodic_face_map_recursively(
920 unsigned int n_face_1,
921 unsigned int n_face_2,
922 const std::bitset<3> & orientation,
928 std::bitset<3>>> &periodic_face_map)
931 const FaceIterator face_1 = cell_1->face(n_face_1);
932 const FaceIterator face_2 = cell_2->face(n_face_2);
934 const bool face_orientation = orientation[0];
935 const bool face_flip = orientation[1];
936 const bool face_rotation = orientation[2];
938 Assert((dim != 1) || (face_orientation ==
true && face_flip ==
false &&
939 face_rotation ==
false),
941 "(face_orientation, face_flip, face_rotation) " 942 "is invalid for 1D"));
944 Assert((dim != 2) || (face_orientation ==
true && face_rotation ==
false),
946 "(face_orientation, face_flip, face_rotation) " 947 "is invalid for 2D"));
951 Assert(face_1->at_boundary() && face_2->at_boundary(),
952 ExcMessage(
"Periodic faces must be on the boundary"));
956 std::pair<typename Triangulation<dim, spacedim>::cell_iterator,
958 const CellFace cell_face_1(cell_1, n_face_1);
959 const CellFace cell_face_2(cell_2, n_face_2);
960 const std::pair<CellFace, std::bitset<3>> cell_face_orientation_2(
961 cell_face_2, orientation);
963 const std::pair<CellFace, std::pair<CellFace, std::bitset<3>>>
964 periodic_faces(cell_face_1, cell_face_orientation_2);
968 periodic_face_map.insert(periodic_faces);
974 static const int lookup_table_2d[2][2] =
981 static const int lookup_table_3d[2][2][2][4] =
1000 if (cell_1->has_children())
1002 if (cell_2->has_children())
1008 Assert(face_1->n_children() ==
1014 for (
unsigned int i = 0;
1015 i < GeometryInfo<dim>::max_children_per_face;
1023 j = lookup_table_2d[face_flip][i];
1026 j = lookup_table_3d[face_orientation][face_flip]
1034 unsigned int child_cell_1 =
1036 cell_1->refinement_case(),
1039 cell_1->face_orientation(n_face_1),
1040 cell_1->face_flip(n_face_1),
1041 cell_1->face_rotation(n_face_1),
1042 face_1->refinement_case());
1043 unsigned int child_cell_2 =
1045 cell_2->refinement_case(),
1048 cell_2->face_orientation(n_face_2),
1049 cell_2->face_flip(n_face_2),
1050 cell_2->face_rotation(n_face_2),
1051 face_2->refinement_case());
1053 Assert(cell_1->child(child_cell_1)->face(n_face_1) ==
1056 Assert(cell_2->child(child_cell_2)->face(n_face_2) ==
1062 update_periodic_face_map_recursively<dim, spacedim>(
1063 cell_1->child(child_cell_1),
1064 cell_2->child(child_cell_2),
1073 for (
unsigned int i = 0;
1074 i < GeometryInfo<dim>::max_children_per_face;
1078 unsigned int child_cell_1 =
1080 cell_1->refinement_case(),
1083 cell_1->face_orientation(n_face_1),
1084 cell_1->face_flip(n_face_1),
1085 cell_1->face_rotation(n_face_1),
1086 face_1->refinement_case());
1089 update_periodic_face_map_recursively<dim, spacedim>(
1090 cell_1->child(child_cell_1),
1107 namespace TriangulationImplementation
1115 using ::Triangulation;
1123 <<
"Something went wrong when making cell " << arg1
1124 <<
". Read the docs and the source code " 1125 <<
"for more information.");
1132 <<
"Something went wrong upon construction of cell " 1146 <<
" has negative measure. This typically " 1147 <<
"indicates some distortion in the cell, or a mistakenly " 1148 <<
"swapped pair of vertices in the input to " 1149 <<
"Triangulation::create_triangulation().");
1161 <<
"Error while creating cell " << arg1
1162 <<
": the vertex index " << arg2 <<
" must be between 0 and " 1171 <<
"While trying to assign a boundary indicator to a line: " 1172 <<
"the line with end vertices " << arg1 <<
" and " << arg2
1173 <<
" does not exist.");
1183 <<
"While trying to assign a boundary indicator to a quad: " 1184 <<
"the quad with bounding lines " << arg1 <<
", " << arg2
1185 <<
", " << arg3 <<
", " << arg4 <<
" does not exist.");
1195 <<
"The input data for creating a triangulation contained " 1196 <<
"information about a line with indices " << arg1 <<
" and " << arg2
1197 <<
" that is described to have boundary indicator " << (
int)arg3
1198 <<
". However, this is an internal line not located on the " 1199 <<
"boundary. You cannot assign a boundary indicator to it." << std::endl
1201 <<
"If this happened at a place where you call " 1202 <<
"Triangulation::create_triangulation() yourself, you need " 1203 <<
"to check the SubCellData object you pass to this function." 1206 <<
"If this happened in a place where you are reading a mesh " 1207 <<
"from a file, then you need to investigate why such a line " 1208 <<
"ended up in the input file. A typical case is a geometry " 1209 <<
"that consisted of multiple parts and for which the mesh " 1210 <<
"generator program assumes that the interface between " 1211 <<
"two parts is a boundary when that isn't supposed to be " 1212 <<
"the case, or where the mesh generator simply assigns " 1213 <<
"'geometry indicators' to lines at the perimeter of " 1214 <<
"a part that are not supposed to be interpreted as " 1215 <<
"'boundary indicators'.");
1227 <<
"The input data for creating a triangulation contained " 1228 <<
"information about a quad with indices " << arg1 <<
", " << arg2
1229 <<
", " << arg3 <<
", and " << arg4
1230 <<
" that is described to have boundary indicator " << (
int)arg5
1231 <<
". However, this is an internal quad not located on the " 1232 <<
"boundary. You cannot assign a boundary indicator to it." << std::endl
1234 <<
"If this happened at a place where you call " 1235 <<
"Triangulation::create_triangulation() yourself, you need " 1236 <<
"to check the SubCellData object you pass to this function." 1239 <<
"If this happened in a place where you are reading a mesh " 1240 <<
"from a file, then you need to investigate why such a quad " 1241 <<
"ended up in the input file. A typical case is a geometry " 1242 <<
"that consisted of multiple parts and for which the mesh " 1243 <<
"generator program assumes that the interface between " 1244 <<
"two parts is a boundary when that isn't supposed to be " 1245 <<
"the case, or where the mesh generator simply assigns " 1246 <<
"'geometry indicators' to quads at the surface of " 1247 <<
"a part that are not supposed to be interpreted as " 1248 <<
"'boundary indicators'.");
1257 <<
"In SubCellData the line info of the line with vertex indices " << arg1
1258 <<
" and " << arg2 <<
" appears more than once. " 1259 <<
"This is not allowed.");
1370 template <
int dim,
int spacedim>
1374 const unsigned int level_objects,
1377 using line_iterator =
1378 typename Triangulation<dim, spacedim>::line_iterator;
1381 if (level_objects > 0)
1383 for (
unsigned int level = 0; level < level_objects; ++level)
1384 if (triangulation.
begin(level) != triangulation.
end(level))
1403 for (
unsigned int level = 0; level < number_cache.
n_levels; ++level)
1409 line_iterator line = triangulation.
begin_line(level),
1411 (level == number_cache.
n_levels - 1 ?
1412 line_iterator(triangulation.
end_line()) :
1414 for (; line != endc; ++line)
1417 if (line->has_children() ==
false)
1433 line_iterator line = triangulation.
begin_line(),
1435 for (; line != endc; ++line)
1438 if (line->has_children() ==
false)
1458 template <
int dim,
int spacedim>
1462 const unsigned int level_objects,
1472 &compute_number_cache<dim, spacedim>),
1478 using quad_iterator =
1479 typename Triangulation<dim, spacedim>::quad_iterator;
1495 unsigned int n_levels = 0;
1496 if (level_objects > 0)
1498 for (
unsigned int level = 0; level < level_objects; ++level)
1499 if (triangulation.
begin(level) != triangulation.
end(level))
1500 n_levels = level + 1;
1505 for (
unsigned int level = 0; level < n_levels; ++level)
1511 quad_iterator quad = triangulation.
begin_quad(level),
1513 (level == n_levels - 1 ?
1514 quad_iterator(triangulation.
end_quad()) :
1516 for (; quad != endc; ++quad)
1519 if (quad->has_children() ==
false)
1535 quad_iterator quad = triangulation.
begin_quad(),
1537 for (; quad != endc; ++quad)
1540 if (quad->has_children() ==
false)
1546 update_lines.
join();
1564 template <
int dim,
int spacedim>
1568 const unsigned int level_objects,
1578 &compute_number_cache<dim, spacedim>),
1584 using hex_iterator =
1585 typename Triangulation<dim, spacedim>::hex_iterator;
1602 unsigned int n_levels = 0;
1603 if (level_objects > 0)
1605 for (
unsigned int level = 0; level < level_objects; ++level)
1606 if (triangulation.
begin(level) != triangulation.
end(level))
1607 n_levels = level + 1;
1612 for (
unsigned int level = 0; level < n_levels; ++level)
1618 hex_iterator hex = triangulation.
begin_hex(level),
1619 endc = (level == n_levels - 1 ?
1620 hex_iterator(triangulation.
end_hex()) :
1622 for (; hex != endc; ++hex)
1625 if (hex->has_children() ==
false)
1641 hex_iterator hex = triangulation.
begin_hex(),
1642 endc = triangulation.
end_hex();
1643 for (; hex != endc; ++hex)
1646 if (hex->has_children() ==
false)
1652 update_quads_and_lines.
join();
1663 template <
int spacedim>
1679 const unsigned int dim = 1;
1683 triangulation.
vertices_used = std::vector<bool>(v.size(),
true);
1691 if (dim == spacedim)
1693 for (
unsigned int cell_no = 0; cell_no < cells.size(); ++cell_no)
1701 const double cell_measure =
1702 GridTools::cell_measure<1>(triangulation.
vertices,
1703 cells[cell_no].vertices);
1715 std::vector<std::vector<int>> lines_at_vertex(v.size());
1718 triangulation.
levels.push_back(
1719 std_cxx14::make_unique<
1721 triangulation.
levels[0]->reserve_space(cells.size(), dim, spacedim);
1722 triangulation.
levels[0]->cells.reserve_space(0, cells.size());
1725 typename Triangulation<dim, spacedim>::raw_line_iterator
1727 for (
unsigned int cell = 0; cell < cells.size(); ++cell)
1729 while (next_free_line->used())
1732 next_free_line->set(
1734 cells[cell].vertices[0], cells[cell].vertices[1]));
1735 next_free_line->set_used_flag();
1736 next_free_line->set_material_id(cells[cell].material_id);
1737 next_free_line->set_manifold_id(cells[cell].manifold_id);
1739 next_free_line->set_subdomain_id(0);
1743 lines_at_vertex[cells[cell].vertices[0]].push_back(cell);
1744 lines_at_vertex[cells[cell].vertices[1]].push_back(cell);
1750 unsigned int boundary_nodes = 0;
1751 for (
unsigned int i = 0; i < lines_at_vertex.size(); ++i)
1752 switch (lines_at_vertex[i].size())
1765 "You have a vertex in your triangulation " 1766 "at which more than two cells come together. " 1767 "(For one dimensional triangulation, cells are " 1770 "This is not currently supported because the " 1771 "Triangulation class makes the assumption that " 1772 "every cell has zero or one neighbors behind " 1773 "each face (here, behind each vertex), but in your " 1774 "situation there would be more than one." 1776 "Support for this is not currently implemented. " 1777 "If you need to work with triangulations where " 1778 "more than two cells come together at a vertex, " 1779 "duplicate the vertices once per cell (i.e., put " 1780 "multiple vertices at the same physical location, " 1781 "but using different vertex indices for each) " 1782 "and then ensure continuity of the solution by " 1783 "explicitly creating constraints that the degrees " 1784 "of freedom at these vertices have the same " 1785 "value, using the AffineConstraints class."));
1797 AssertThrow(((spacedim == 1) && (boundary_nodes == 2)) ||
1799 ExcMessage(
"The Triangulation has too many end points"));
1805 typename Triangulation<dim, spacedim>::active_line_iterator line =
1808 for (; line != triangulation.
end(); ++line)
1810 for (
unsigned int vertex = 0;
1811 vertex < GeometryInfo<dim>::vertices_per_cell;
1818 if (lines_at_vertex[line->vertex_index(vertex)][0] == line->index())
1819 if (lines_at_vertex[line->vertex_index(vertex)].size() == 2)
1822 neighbor(&triangulation,
1824 lines_at_vertex[line->vertex_index(vertex)][1]);
1825 line->set_neighbor(vertex, neighbor);
1831 line->set_neighbor(vertex, triangulation.
end());
1838 neighbor(&triangulation,
1840 lines_at_vertex[line->vertex_index(vertex)][0]);
1841 line->set_neighbor(vertex, neighbor);
1852 cell != triangulation.
end();
1854 for (
unsigned int f = 0; f < GeometryInfo<dim>::faces_per_cell; ++f)
1859 if (cell->at_boundary(f))
1861 [cell->face(f)->vertex_index()] = f;
1873 template <
int spacedim>
1883 const unsigned int dim = 2;
1887 triangulation.
vertices_used = std::vector<bool>(v.size(),
true);
1895 if (dim == spacedim)
1897 for (
unsigned int cell_no = 0; cell_no < cells.size(); ++cell_no)
1902 const double cell_measure =
1903 GridTools::cell_measure<2>(triangulation.
vertices,
1904 cells[cell_no].vertices);
1925 std::map<std::pair<int, int>,
1926 typename Triangulation<dim, spacedim>::line_iterator>
1928 for (
unsigned int cell = 0; cell < cells.size(); ++cell)
1930 for (
unsigned int vertex = 0; vertex < 4; ++vertex)
1934 cells[cell].vertices[vertex],
1937 for (
unsigned int line = 0;
1938 line < GeometryInfo<dim>::faces_per_cell;
1944 std::pair<int, int> line_vertices(
1967 line_vertices.second, line_vertices.first)) ==
1975 needed_lines[line_vertices] = triangulation.
end_line();
1983 std::vector<unsigned short int> vertex_touch_count(v.size(), 0);
1985 std::pair<int, int>,
1986 typename Triangulation<dim, spacedim>::line_iterator>::iterator i;
1987 for (i = needed_lines.begin(); i != needed_lines.end(); ++i)
1991 ++vertex_touch_count[i->first.first];
1992 ++vertex_touch_count[i->first.second];
1999 AssertThrow(*(std::min_element(vertex_touch_count.begin(),
2000 vertex_touch_count.end())) >= 2,
2002 "During creation of a triangulation, a part of the " 2003 "algorithm encountered a vertex that is part of only " 2004 "a single adjacent line. However, in 2d, every vertex " 2005 "needs to be at least part of two lines."));
2009 triangulation.
levels.push_back(
2010 std_cxx14::make_unique<
2012 triangulation.
faces = std_cxx14::make_unique<
2014 triangulation.
levels[0]->reserve_space(cells.size(), dim, spacedim);
2015 triangulation.
faces->lines.reserve_space(0, needed_lines.size());
2016 triangulation.
levels[0]->cells.reserve_space(0, cells.size());
2020 typename Triangulation<dim, spacedim>::raw_line_iterator line =
2023 std::pair<int, int>,
2024 typename Triangulation<dim, spacedim>::line_iterator>::iterator i;
2025 for (i = needed_lines.begin(); line != triangulation.
end_line();
2029 i->first.first, i->first.second));
2030 line->set_used_flag();
2031 line->clear_user_flag();
2042 std::vector<typename Triangulation<dim, spacedim>::cell_iterator>>
2049 for (
unsigned int c = 0; c < cells.size(); ++c, ++cell)
2051 typename Triangulation<dim, spacedim>::line_iterator
2053 for (
unsigned int line = 0;
2054 line < GeometryInfo<dim>::lines_per_cell;
2056 lines[line] = needed_lines[std::make_pair(
2066 lines[3]->index()));
2068 cell->set_used_flag();
2069 cell->set_material_id(cells[c].material_id);
2070 cell->set_manifold_id(cells[c].manifold_id);
2071 cell->clear_user_data();
2072 cell->set_subdomain_id(0);
2077 for (
unsigned int line = 0;
2078 line < GeometryInfo<dim>::lines_per_cell;
2080 adjacent_cells[lines[line]->index()].push_back(cell);
2085 for (
typename Triangulation<dim, spacedim>::line_iterator line =
2090 const unsigned int n_adj_cells =
2091 adjacent_cells[line->index()].size();
2099 AssertThrow((n_adj_cells >= 1) && (n_adj_cells <= 2),
2105 (n_adj_cells >= 1) && (n_adj_cells <= 2),
2106 ExcMessage(
"You have a line in your triangulation at which " 2107 "more than two cells come together." 2109 "This is not currently supported because the " 2110 "Triangulation class makes the assumption that " 2111 "every cell has zero or one neighbors behind each " 2112 "face (here, behind each line), but in your " 2113 "situation there would be more than one." 2115 "Support for this is not currently implemented. " 2116 "If you need to work with triangulations where " 2117 "more than two cells come together at a line, " 2118 "duplicate the vertices once per cell (i.e., put " 2119 "multiple vertices at the same physical location, " 2120 "but using different vertex indices for each) " 2121 "and then ensure continuity of the solution by " 2122 "explicitly creating constraints that the degrees " 2123 "of freedom at these lines have the same " 2124 "value, using the AffineConstraints class."));
2129 line->set_boundary_id_internal(
2137 typename Triangulation<dim, spacedim>::line_iterator line;
2138 std::pair<int, int> line_vertices(
2139 std::make_pair(subcell_line.vertices[0],
2140 subcell_line.vertices[1]));
2141 if (needed_lines.find(line_vertices) != needed_lines.end())
2143 line = needed_lines[line_vertices];
2147 std::swap(line_vertices.first, line_vertices.second);
2148 if (needed_lines.find(line_vertices) != needed_lines.end())
2149 line = needed_lines[line_vertices];
2154 line_vertices.second));
2159 line->boundary_id() !=
2162 line_vertices.second));
2191 line->vertex_index(0),
2192 line->vertex_index(1),
2193 subcell_line.boundary_id));
2197 line->set_manifold_id(subcell_line.manifold_id);
2201 line->set_boundary_id_internal(subcell_line.boundary_id);
2203 line->set_manifold_id(subcell_line.manifold_id);
2209 triangulation.
begin();
2210 cell != triangulation.
end();
2212 for (
unsigned int side = 0; side < 4; ++side)
2213 if (adjacent_cells[cell->line(side)->index()][0] == cell)
2217 if (adjacent_cells[cell->line(side)->index()].size() == 2)
2221 side, adjacent_cells[cell->line(side)->index()][1]);
2227 cell->set_neighbor(side,
2228 adjacent_cells[cell->line(side)->index()][0]);
2276 template <
int spacedim>
2286 const unsigned int dim = 3;
2290 triangulation.
vertices_used = std::vector<bool>(v.size(),
true);
2296 for (
unsigned int cell_no = 0; cell_no < cells.size(); ++cell_no)
2301 const double cell_measure =
2302 GridTools::cell_measure<3>(triangulation.
vertices,
2303 cells[cell_no].vertices);
2328 typename std::map<std::pair<int, int>,
2329 typename Triangulation<dim, spacedim>::line_iterator>
2331 for (
unsigned int cell = 0; cell < cells.size(); ++cell)
2335 for (
unsigned int vertex = 0;
2336 vertex < GeometryInfo<dim>::vertices_per_cell;
2341 cells[cell].vertices[vertex],
2344 for (
unsigned int line = 0;
2345 line < GeometryInfo<dim>::lines_per_cell;
2354 std::pair<int, int> line_vertices(
2363 if ((needed_lines.find(std::make_pair(line_vertices.second,
2364 line_vertices.first)) ==
2365 needed_lines.end()))
2371 needed_lines[line_vertices] = triangulation.
end_line();
2383 std::vector<unsigned short int> vertex_touch_count(v.size(), 0);
2385 std::pair<int, int>,
2386 typename Triangulation<dim, spacedim>::line_iterator>::iterator i;
2387 for (i = needed_lines.begin(); i != needed_lines.end(); ++i)
2391 ++vertex_touch_count[i->first.first];
2392 ++vertex_touch_count[i->first.second];
2400 *(std::min_element(vertex_touch_count.begin(),
2401 vertex_touch_count.end())) >= 3,
2403 "During creation of a triangulation, a part of the " 2404 "algorithm encountered a vertex that is part of only " 2405 "one or two adjacent lines. However, in 3d, every vertex " 2406 "needs to be at least part of three lines."));
2414 triangulation.
levels.push_back(
2415 std_cxx14::make_unique<
2417 triangulation.
faces = std_cxx14::make_unique<
2419 triangulation.
levels[0]->reserve_space(cells.size(), dim, spacedim);
2420 triangulation.
faces->lines.reserve_space(0, needed_lines.size());
2424 typename Triangulation<dim, spacedim>::raw_line_iterator line =
2427 std::pair<int, int>,
2428 typename Triangulation<dim, spacedim>::line_iterator>::iterator i;
2429 for (i = needed_lines.begin(); line != triangulation.
end_line();
2433 i->first.first, i->first.second));
2434 line->set_used_flag();
2435 line->clear_user_flag();
2459 std::map<internal::TriangulationImplementation::TriaObject<2>,
2460 std::pair<typename Triangulation<dim, spacedim>::quad_iterator,
2461 std::array<bool, GeometryInfo<dim>::lines_per_face>>,
2464 for (
unsigned int cell = 0; cell < cells.size(); ++cell)
2492 std::array<bool, GeometryInfo<dim>::lines_per_face> orientation;
2494 for (
unsigned int line = 0;
2495 line < GeometryInfo<dim>::lines_per_cell;
2498 line_list[line] = std::pair<int, int>(
2503 inverse_line_list[line] = std::pair<int, int>(
2510 for (
unsigned int face = 0;
2511 face < GeometryInfo<dim>::faces_per_cell;
2522 for (
unsigned int l = 0; l < GeometryInfo<dim>::lines_per_face;
2524 if (needed_lines.find(
2526 face, l)]) == needed_lines.end())
2530 dim>::face_to_cell_lines(face, l)]]
2532 orientation[l] =
true;
2538 dim>::face_to_cell_lines(face, l)]]
2540 orientation[l] =
false;
2589 test_quad_1(quad.
face(2),
2595 test_quad_2(quad.
face(0),
2601 test_quad_3(quad.
face(3),
2607 test_quad_4(quad.
face(1),
2613 test_quad_5(quad.
face(2),
2619 test_quad_6(quad.
face(1),
2625 test_quad_7(quad.
face(3),
2631 if (needed_quads.find(test_quad_1) == needed_quads.end() &&
2632 needed_quads.find(test_quad_2) == needed_quads.end() &&
2633 needed_quads.find(test_quad_3) == needed_quads.end() &&
2634 needed_quads.find(test_quad_4) == needed_quads.end() &&
2635 needed_quads.find(test_quad_5) == needed_quads.end() &&
2636 needed_quads.find(test_quad_6) == needed_quads.end() &&
2637 needed_quads.find(test_quad_7) == needed_quads.end())
2638 needed_quads[quad] =
2639 std::make_pair(triangulation.
end_quad(), orientation);
2649 triangulation.
faces->quads.reserve_space(0, needed_quads.size());
2652 typename Triangulation<dim, spacedim>::raw_quad_iterator quad =
2656 std::pair<typename Triangulation<dim, spacedim>::quad_iterator,
2657 std::array<bool, GeometryInfo<dim>::lines_per_face>>,
2659 for (q = needed_quads.begin(); quad != triangulation.
end_quad();
2662 quad->set(q->first);
2663 quad->set_used_flag();
2664 quad->clear_user_flag();
2667 quad->set_line_orientation(0, q->second.second[0]);
2668 quad->set_line_orientation(1, q->second.second[1]);
2669 quad->set_line_orientation(2, q->second.second[2]);
2670 quad->set_line_orientation(3, q->second.second[3]);
2675 q->second.first = quad;
2681 triangulation.
levels[0]->cells.reserve_space(cells.size());
2687 std::vector<typename Triangulation<dim, spacedim>::cell_iterator>>
2694 for (
unsigned int c = 0; c < cells.size(); ++c, ++cell)
2709 unsigned int face_line_list[4];
2710 for (
unsigned int line = 0;
2711 line < GeometryInfo<dim>::lines_per_cell;
2714 line_list[line] = std::make_pair(
2719 inverse_line_list[line] = std::pair<int, int>(
2731 typename Triangulation<dim, spacedim>::quad_iterator
2736 for (
unsigned int face = 0;
2737 face < GeometryInfo<dim>::faces_per_cell;
2740 for (
unsigned int l = 0;
2741 l < GeometryInfo<dim>::lines_per_face;
2744 dim>::face_to_cell_lines(face, l)]) ==
2748 dim>::face_to_cell_lines(face, l)]]
2753 dim>::face_to_cell_lines(face, l)]]
2762 if (needed_quads.find(quad) != needed_quads.end())
2773 face_iterator[face] = needed_quads[quad].first;
2774 face_orientation[face] =
true;
2775 face_flip[face] =
false;
2776 face_rotation[face] =
false;
2804 test_quad_4(quad.
face(1),
2822 test_quad_7(quad.
face(3),
2828 if (needed_quads.find(test_quad_1) != needed_quads.end())
2830 face_iterator[face] = needed_quads[test_quad_1].first;
2831 face_orientation[face] =
false;
2832 face_flip[face] =
false;
2833 face_rotation[face] =
false;
2835 else if (needed_quads.find(test_quad_2) !=
2838 face_iterator[face] = needed_quads[test_quad_2].first;
2839 face_orientation[face] =
false;
2840 face_flip[face] =
false;
2841 face_rotation[face] =
true;
2843 else if (needed_quads.find(test_quad_3) !=
2846 face_iterator[face] = needed_quads[test_quad_3].first;
2847 face_orientation[face] =
false;
2848 face_flip[face] =
true;
2849 face_rotation[face] =
false;
2851 else if (needed_quads.find(test_quad_4) !=
2854 face_iterator[face] = needed_quads[test_quad_4].first;
2855 face_orientation[face] =
false;
2856 face_flip[face] =
true;
2857 face_rotation[face] =
true;
2859 else if (needed_quads.find(test_quad_5) !=
2862 face_iterator[face] = needed_quads[test_quad_5].first;
2863 face_orientation[face] =
true;
2864 face_flip[face] =
false;
2865 face_rotation[face] =
true;
2867 else if (needed_quads.find(test_quad_6) !=
2870 face_iterator[face] = needed_quads[test_quad_6].first;
2871 face_orientation[face] =
true;
2872 face_flip[face] =
true;
2873 face_rotation[face] =
false;
2875 else if (needed_quads.find(test_quad_7) !=
2878 face_iterator[face] = needed_quads[test_quad_7].first;
2879 face_orientation[face] =
true;
2880 face_flip[face] =
true;
2881 face_rotation[face] =
true;
2896 face_iterator[0]->index(),
2897 face_iterator[1]->index(),
2898 face_iterator[2]->index(),
2899 face_iterator[3]->index(),
2900 face_iterator[4]->index(),
2901 face_iterator[5]->index()));
2903 cell->set_used_flag();
2904 cell->set_material_id(cells[c].material_id);
2905 cell->set_manifold_id(cells[c].manifold_id);
2906 cell->clear_user_flag();
2907 cell->clear_user_data();
2908 cell->set_subdomain_id(0);
2912 for (
unsigned int quad = 0;
2913 quad < GeometryInfo<dim>::faces_per_cell;
2916 cell->set_face_orientation(quad, face_orientation[quad]);
2917 cell->set_face_flip(quad, face_flip[quad]);
2918 cell->set_face_rotation(quad, face_rotation[quad]);
2925 for (
unsigned int quad = 0;
2926 quad < GeometryInfo<dim>::faces_per_cell;
2928 adjacent_cells[face_iterator[quad]->index()].push_back(cell);
2945 std::multimap<unsigned int, std::pair<unsigned int, unsigned int>>
2947 for (
unsigned int face = 0;
2948 face < GeometryInfo<dim>::faces_per_cell;
2950 for (
unsigned int line = 0;
2951 line < GeometryInfo<dim>::lines_per_face;
2953 cell_to_face_lines.insert(
2954 std::pair<
unsigned int,
2955 std::pair<unsigned int, unsigned int>>(
2957 std::pair<unsigned int, unsigned int>(face, line)));
2958 std::multimap<
unsigned int,
2959 std::pair<unsigned int, unsigned int>>::
2960 const_iterator map_iter = cell_to_face_lines.begin();
2962 for (; map_iter != cell_to_face_lines.end(); ++map_iter)
2964 const unsigned int cell_line = map_iter->first;
2965 const unsigned int face1 = map_iter->second.first;
2966 const unsigned int line1 = map_iter->second.second;
2968 Assert(map_iter != cell_to_face_lines.end(),
2970 Assert(map_iter->first == cell_line,
2972 const unsigned int face2 = map_iter->second.first;
2973 const unsigned int line2 = map_iter->second.second;
2980 Assert(face_iterator[face1]->line(
2983 face_orientation[face1],
2985 face_rotation[face1])) ==
2986 face_iterator[face2]->line(
2989 face_orientation[face2],
2991 face_rotation[face2])),
3002 for (
typename Triangulation<dim, spacedim>::quad_iterator quad =
3007 const unsigned int n_adj_cells =
3008 adjacent_cells[quad->index()].size();
3011 AssertThrow((n_adj_cells >= 1) && (n_adj_cells <= 2),
3016 quad->set_boundary_id_internal(
3030 for (
typename Triangulation<dim, spacedim>::line_iterator line =
3057 for (
typename Triangulation<dim, spacedim>::quad_iterator quad =
3061 if (quad->at_boundary())
3063 for (
unsigned int l = 0; l < 4; ++l)
3065 typename Triangulation<dim, spacedim>::line_iterator line =
3067 line->set_boundary_id_internal(0);
3078 typename Triangulation<dim, spacedim>::line_iterator line;
3079 std::pair<int, int> line_vertices(
3080 std::make_pair(subcell_line.vertices[0],
3081 subcell_line.vertices[1]));
3082 if (needed_lines.find(line_vertices) != needed_lines.end())
3085 line = needed_lines[line_vertices];
3091 std::swap(line_vertices.first, line_vertices.second);
3092 if (needed_lines.find(line_vertices) != needed_lines.end())
3093 line = needed_lines[line_vertices];
3098 line_vertices.second));
3101 if (line->at_boundary())
3105 if (line->boundary_id() != 0)
3106 AssertThrow(line->boundary_id() == subcell_line.boundary_id,
3108 "Duplicate boundary lines are only allowed " 3109 "if they carry the same boundary indicator."));
3111 line->set_boundary_id_internal(subcell_line.boundary_id);
3115 AssertThrow(line->manifold_id() == subcell_line.manifold_id,
3117 "Duplicate lines are only allowed " 3118 "if they carry the same manifold indicator."));
3119 line->set_manifold_id(subcell_line.manifold_id);
3126 typename Triangulation<dim, spacedim>::quad_iterator quad;
3127 typename Triangulation<dim, spacedim>::line_iterator line[4];
3136 for (
unsigned int i = 0; i < 4; ++i)
3138 std::pair<int, int> line_vertices(
3148 if (needed_lines.find(line_vertices) != needed_lines.end())
3149 line[i] = needed_lines[line_vertices];
3154 std::swap(line_vertices.first, line_vertices.second);
3155 if (needed_lines.find(line_vertices) != needed_lines.end())
3156 line[i] = needed_lines[line_vertices];
3162 line_vertices.second));
3201 static const unsigned int lex2cclock[4] = {3, 1, 0, 2};
3207 typename Triangulation<dim, spacedim>::line_iterator
3208 line_counterclock[4];
3209 for (
unsigned int i = 0; i < 4; ++i)
3210 line_counterclock[lex2cclock[i]] = line[i];
3211 unsigned int n_rotations = 0;
3212 bool not_found_quad_1;
3213 while ((not_found_quad_1 = (needed_quads.find(quad_compare_1) ==
3214 needed_quads.end())) &&
3215 (needed_quads.find(quad_compare_2) == needed_quads.end()) &&
3220 std::rotate(line_counterclock,
3221 line_counterclock + 1,
3222 line_counterclock + 4);
3226 for (
unsigned int i = 0; i < 4; ++i)
3229 i, line_counterclock[lex2cclock[i]]->index());
3231 (i + 2) % 4, line_counterclock[lex2cclock[i]]->index());
3243 if (not_found_quad_1)
3244 quad = needed_quads[quad_compare_2].first;
3246 quad = needed_quads[quad_compare_1].first;
3250 if (quad->at_boundary())
3254 if (quad->boundary_id() != 0)
3255 AssertThrow(quad->boundary_id() == subcell_quad.boundary_id,
3257 "Duplicate boundary quads are only allowed " 3258 "if they carry the same boundary indicator."));
3260 quad->set_boundary_id_internal(subcell_quad.boundary_id);
3264 AssertThrow(quad->manifold_id() == subcell_quad.manifold_id,
3266 "Duplicate boundary quads are only allowed " 3267 "if they carry the same manifold indicator."));
3269 quad->set_manifold_id(subcell_quad.manifold_id);
3276 triangulation.
begin();
3277 cell != triangulation.
end();
3279 for (
unsigned int face = 0; face < 6; ++face)
3280 if (adjacent_cells[cell->quad(face)->index()][0] == cell)
3284 if (adjacent_cells[cell->quad(face)->index()].size() == 2)
3288 face, adjacent_cells[cell->quad(face)->index()][1]);
3294 cell->set_neighbor(face,
3295 adjacent_cells[cell->quad(face)->index()][0]);
3314 template <
int spacedim>
3318 std::vector<unsigned int> &,
3319 std::vector<unsigned int> &)
3321 const unsigned int dim = 1;
3333 Assert(!cell->child(0)->has_children() &&
3334 !cell->child(1)->has_children(),
3340 if (cell->neighbor(0)->has_children())
3347 neighbor = neighbor->child(1);
3350 Assert(neighbor->neighbor(1) == cell->child(0),
3352 neighbor->set_neighbor(1, cell);
3358 if (neighbor->has_children())
3359 neighbor = neighbor->child(1);
3368 if (cell->neighbor(1)->has_children())
3375 neighbor = neighbor->child(0);
3378 Assert(neighbor->neighbor(0) == cell->child(1),
3380 neighbor->set_neighbor(0, cell);
3386 if (neighbor->has_children())
3387 neighbor = neighbor->child(0);
3397 triangulation.
vertices_used[cell->child(0)->vertex_index(1)] =
false;
3403 for (
unsigned int child = 0; child < cell->n_children(); ++child)
3406 cell->child(child)->clear_user_flag();
3407 cell->child(child)->clear_used_flag();
3412 cell->clear_children();
3413 cell->clear_user_flag();
3418 template <
int spacedim>
3419 static void delete_children(
3422 std::vector<unsigned int> & line_cell_count,
3423 std::vector<unsigned int> &)
3425 const unsigned int dim = 2;
3433 std::vector<typename Triangulation<dim, spacedim>::line_iterator>
3436 lines_to_delete.reserve(4 * 2 + 4);
3441 for (
unsigned int c = 0; c < cell->n_children(); ++c)
3445 for (
unsigned int l = 0; l < GeometryInfo<dim>::lines_per_cell; ++l)
3446 --line_cell_count[child->line_index(l)];
3461 .
vertices_used[cell->child(0)->line(1)->vertex_index(1)] =
false;
3463 lines_to_delete.push_back(cell->child(0)->line(1));
3464 lines_to_delete.push_back(cell->child(0)->line(3));
3465 lines_to_delete.push_back(cell->child(3)->line(0));
3466 lines_to_delete.push_back(cell->child(3)->line(2));
3470 unsigned int inner_face_no =
3475 lines_to_delete.push_back(cell->child(0)->line(inner_face_no));
3479 for (
unsigned int child = 0; child < cell->n_children(); ++child)
3481 cell->child(child)->clear_user_data();
3482 cell->child(child)->clear_user_flag();
3483 cell->child(child)->clear_used_flag();
3488 cell->clear_children();
3489 cell->clear_refinement_case();
3490 cell->clear_user_flag();
3496 for (
unsigned int line_no = 0;
3497 line_no < GeometryInfo<dim>::lines_per_cell;
3500 typename Triangulation<dim, spacedim>::line_iterator line =
3501 cell->line(line_no);
3503 if (line->has_children())
3508 Assert((line_cell_count[line->child_index(0)] == 0 &&
3509 line_cell_count[line->child_index(1)] == 0) ||
3510 (line_cell_count[line->child_index(0)] > 0 &&
3511 line_cell_count[line->child_index(1)] > 0),
3514 if (line_cell_count[line->child_index(0)] == 0)
3516 for (
unsigned int c = 0; c < 2; ++c)
3517 Assert(!line->child(c)->has_children(),
3527 lines_to_delete.push_back(line->child(0));
3528 lines_to_delete.push_back(line->child(1));
3530 line->clear_children();
3542 typename std::vector<
3543 typename Triangulation<dim, spacedim>::line_iterator>::iterator
3544 line = lines_to_delete.
begin(),
3545 endline = lines_to_delete.end();
3546 for (; line != endline; ++line)
3548 (*line)->clear_user_data();
3549 (*line)->clear_user_flag();
3550 (*line)->clear_used_flag();
3556 template <
int spacedim>
3557 static void delete_children(
3560 std::vector<unsigned int> & line_cell_count,
3561 std::vector<unsigned int> & quad_cell_count)
3563 const unsigned int dim = 3;
3575 std::vector<typename Triangulation<dim, spacedim>::line_iterator>
3577 std::vector<typename Triangulation<dim, spacedim>::quad_iterator>
3580 lines_to_delete.reserve(12 * 2 + 6 * 4 + 6);
3581 quads_to_delete.reserve(6 * 4 + 12);
3585 for (
unsigned int c = 0; c < cell->n_children(); ++c)
3589 for (
unsigned int l = 0; l < GeometryInfo<dim>::lines_per_cell; ++l)
3590 --line_cell_count[child->line_index(l)];
3591 for (
unsigned int f = 0; f < GeometryInfo<dim>::faces_per_cell; ++f)
3592 --quad_cell_count[child->quad_index(f)];
3606 quads_to_delete.push_back(cell->child(0)->face(1));
3609 quads_to_delete.push_back(cell->child(0)->face(3));
3612 quads_to_delete.push_back(cell->child(0)->face(5));
3615 quads_to_delete.push_back(cell->child(0)->face(1));
3616 quads_to_delete.push_back(cell->child(0)->face(3));
3617 quads_to_delete.push_back(cell->child(3)->face(0));
3618 quads_to_delete.push_back(cell->child(3)->face(2));
3620 lines_to_delete.push_back(cell->child(0)->line(11));
3623 quads_to_delete.push_back(cell->child(0)->face(1));
3624 quads_to_delete.push_back(cell->child(0)->face(5));
3625 quads_to_delete.push_back(cell->child(3)->face(0));
3626 quads_to_delete.push_back(cell->child(3)->face(4));
3628 lines_to_delete.push_back(cell->child(0)->line(5));
3631 quads_to_delete.push_back(cell->child(0)->face(3));
3632 quads_to_delete.push_back(cell->child(0)->face(5));
3633 quads_to_delete.push_back(cell->child(3)->face(2));
3634 quads_to_delete.push_back(cell->child(3)->face(4));
3636 lines_to_delete.push_back(cell->child(0)->line(7));
3639 quads_to_delete.push_back(cell->child(0)->face(1));
3640 quads_to_delete.push_back(cell->child(2)->face(1));
3641 quads_to_delete.push_back(cell->child(4)->face(1));
3642 quads_to_delete.push_back(cell->child(6)->face(1));
3644 quads_to_delete.push_back(cell->child(0)->face(3));
3645 quads_to_delete.push_back(cell->child(1)->face(3));
3646 quads_to_delete.push_back(cell->child(4)->face(3));
3647 quads_to_delete.push_back(cell->child(5)->face(3));
3649 quads_to_delete.push_back(cell->child(0)->face(5));
3650 quads_to_delete.push_back(cell->child(1)->face(5));
3651 quads_to_delete.push_back(cell->child(2)->face(5));
3652 quads_to_delete.push_back(cell->child(3)->face(5));
3654 lines_to_delete.push_back(cell->child(0)->line(5));
3655 lines_to_delete.push_back(cell->child(0)->line(7));
3656 lines_to_delete.push_back(cell->child(0)->line(11));
3657 lines_to_delete.push_back(cell->child(7)->line(0));
3658 lines_to_delete.push_back(cell->child(7)->line(2));
3659 lines_to_delete.push_back(cell->child(7)->line(8));
3665 triangulation.
vertices_used[cell->child(0)->vertex_index(7)] =
3677 for (
unsigned int child = 0; child < cell->n_children(); ++child)
3679 cell->child(child)->clear_user_data();
3680 cell->child(child)->clear_user_flag();
3682 for (
unsigned int f = 0; f < GeometryInfo<dim>::faces_per_cell; ++f)
3687 cell->child(child)->set_face_orientation(f,
true);
3688 cell->child(child)->set_face_flip(f,
false);
3689 cell->child(child)->set_face_rotation(f,
false);
3692 cell->child(child)->clear_used_flag();
3697 cell->clear_children();
3698 cell->clear_refinement_case();
3699 cell->clear_user_flag();
3709 for (
unsigned int quad_no = 0;
3710 quad_no < GeometryInfo<dim>::faces_per_cell;
3713 typename Triangulation<dim, spacedim>::quad_iterator quad =
3714 cell->face(quad_no);
3718 quad->has_children()) ||
3723 switch (quad->refinement_case())
3735 Assert((quad_cell_count[quad->child_index(0)] == 0 &&
3736 quad_cell_count[quad->child_index(1)] == 0) ||
3737 (quad_cell_count[quad->child_index(0)] > 0 &&
3738 quad_cell_count[quad->child_index(1)] > 0),
3744 unsigned int deleted_grandchildren = 0;
3745 unsigned int number_of_child_refinements = 0;
3747 for (
unsigned int c = 0; c < 2; ++c)
3748 if (quad->child(c)->has_children())
3750 ++number_of_child_refinements;
3755 (quad_cell_count[quad->child(c)->child_index(0)] ==
3757 quad_cell_count[quad->child(c)->child_index(1)] ==
3759 (quad_cell_count[quad->child(c)->child_index(0)] >
3761 quad_cell_count[quad->child(c)->child_index(1)] >
3764 if (quad_cell_count[quad->child(c)->child_index(0)] ==
3771 Assert(quad->refinement_case() +
3772 quad->child(c)->refinement_case() ==
3780 quads_to_delete.push_back(
3781 quad->child(c)->child(0));
3782 quads_to_delete.push_back(
3783 quad->child(c)->child(1));
3784 if (quad->child(c)->refinement_case() ==
3786 lines_to_delete.push_back(
3787 quad->child(c)->child(0)->line(1));
3789 lines_to_delete.push_back(
3790 quad->child(c)->child(0)->line(3));
3791 quad->child(c)->clear_children();
3792 quad->child(c)->clear_refinement_case();
3793 ++deleted_grandchildren;
3801 if (number_of_child_refinements > 0 &&
3802 deleted_grandchildren == number_of_child_refinements)
3804 typename Triangulation<dim, spacedim>::line_iterator
3807 middle_line = quad->child(0)->line(1);
3809 middle_line = quad->child(0)->line(3);
3811 lines_to_delete.push_back(middle_line->child(0));
3812 lines_to_delete.push_back(middle_line->child(1));
3815 middle_line)] =
false;
3816 middle_line->clear_children();
3821 if (quad_cell_count[quad->child_index(0)] == 0)
3827 quads_to_delete.push_back(quad->child(0));
3828 quads_to_delete.push_back(quad->child(1));
3830 lines_to_delete.push_back(quad->child(0)->line(1));
3832 lines_to_delete.push_back(quad->child(0)->line(3));
3852 if (quad->child(0)->has_children())
3853 if (quad->refinement_case() ==
3887 quad_iterator switch_1 =
3888 quad->child(0)->child(1),
3890 quad->child(1)->child(0);
3892 Assert(!switch_1->has_children(),
3894 Assert(!switch_2->has_children(),
3897 const int switch_1_index = switch_1->index();
3898 const int switch_2_index = switch_2->index();
3899 for (
unsigned int l = 0;
3900 l < triangulation.
levels.size();
3902 for (
unsigned int h = 0;
3903 h < triangulation.
levels[l]
3904 ->cells.cells.size();
3906 for (
unsigned int q = 0;
3907 q < GeometryInfo<dim>::faces_per_cell;
3910 const int index = triangulation.
levels[l]
3913 if (index == switch_1_index)
3916 .set_face(q, switch_2_index);
3917 else if (index == switch_2_index)
3920 .set_face(q, switch_1_index);
3925 const int switch_1_lines[4] = {
3926 static_cast<signed int>(
3927 switch_1->line_index(0)),
3928 static_cast<signed int>(
3929 switch_1->line_index(1)),
3930 static_cast<signed int>(
3931 switch_1->line_index(2)),
3932 static_cast<signed int>(
3933 switch_1->line_index(3))};
3934 const bool switch_1_line_orientations[4] = {
3935 switch_1->line_orientation(0),
3936 switch_1->line_orientation(1),
3937 switch_1->line_orientation(2),
3938 switch_1->line_orientation(3)};
3940 switch_1->boundary_id();
3941 const unsigned int switch_1_user_index =
3942 switch_1->user_index();
3943 const bool switch_1_user_flag =
3944 switch_1->user_flag_set();
3949 switch_2->line_index(1),
3950 switch_2->line_index(2),
3951 switch_2->line_index(3)));
3952 switch_1->set_line_orientation(
3953 0, switch_2->line_orientation(0));
3954 switch_1->set_line_orientation(
3955 1, switch_2->line_orientation(1));
3956 switch_1->set_line_orientation(
3957 2, switch_2->line_orientation(2));
3958 switch_1->set_line_orientation(
3959 3, switch_2->line_orientation(3));
3960 switch_1->set_boundary_id_internal(
3961 switch_2->boundary_id());
3962 switch_1->set_manifold_id(
3963 switch_2->manifold_id());
3964 switch_1->set_user_index(switch_2->user_index());
3965 if (switch_2->user_flag_set())
3966 switch_1->set_user_flag();
3968 switch_1->clear_user_flag();
3975 switch_1_lines[3]));
3976 switch_2->set_line_orientation(
3977 0, switch_1_line_orientations[0]);
3978 switch_2->set_line_orientation(
3979 1, switch_1_line_orientations[1]);
3980 switch_2->set_line_orientation(
3981 2, switch_1_line_orientations[2]);
3982 switch_2->set_line_orientation(
3983 3, switch_1_line_orientations[3]);
3984 switch_2->set_boundary_id_internal(
3985 switch_1_boundary_id);
3986 switch_2->set_manifold_id(
3987 switch_1->manifold_id());
3988 switch_2->set_user_index(switch_1_user_index);
3989 if (switch_1_user_flag)
3990 switch_2->set_user_flag();
3992 switch_2->clear_user_flag();
3994 const unsigned int child_0 =
3995 quad->child(0)->child_index(0);
3996 const unsigned int child_2 =
3997 quad->child(1)->child_index(0);
3998 quad->clear_children();
3999 quad->clear_refinement_case();
4000 quad->set_refinement_case(
4002 quad->set_children(0, child_0);
4003 quad->set_children(2, child_2);
4005 quad_cell_count[child_2]);
4020 const unsigned int child_0 =
4021 quad->child(0)->child_index(0);
4022 const unsigned int child_2 =
4023 quad->child(1)->child_index(0);
4024 quad->clear_children();
4025 quad->clear_refinement_case();
4026 quad->set_refinement_case(
4028 quad->set_children(0, child_0);
4029 quad->set_children(2, child_2);
4033 quad->clear_children();
4034 quad->clear_refinement_case();
4045 Assert((quad_cell_count[quad->child_index(0)] == 0 &&
4046 quad_cell_count[quad->child_index(1)] == 0 &&
4047 quad_cell_count[quad->child_index(2)] == 0 &&
4048 quad_cell_count[quad->child_index(3)] == 0) ||
4049 (quad_cell_count[quad->child_index(0)] > 0 &&
4050 quad_cell_count[quad->child_index(1)] > 0 &&
4051 quad_cell_count[quad->child_index(2)] > 0 &&
4052 quad_cell_count[quad->child_index(3)] > 0),
4055 if (quad_cell_count[quad->child_index(0)] == 0)
4061 lines_to_delete.push_back(quad->child(0)->line(1));
4062 lines_to_delete.push_back(quad->child(3)->line(0));
4063 lines_to_delete.push_back(quad->child(0)->line(3));
4064 lines_to_delete.push_back(quad->child(3)->line(2));
4066 for (
unsigned int child = 0; child < quad->n_children();
4068 quads_to_delete.push_back(quad->child(child));
4074 quad->clear_children();
4075 quad->clear_refinement_case();
4094 for (
unsigned int line_no = 0;
4095 line_no < GeometryInfo<dim>::lines_per_cell;
4098 typename Triangulation<dim, spacedim>::line_iterator line =
4099 cell->line(line_no);
4103 line->has_children()) ||
4108 if (line->has_children())
4113 Assert((line_cell_count[line->child_index(0)] == 0 &&
4114 line_cell_count[line->child_index(1)] == 0) ||
4115 (line_cell_count[line->child_index(0)] > 0 &&
4116 line_cell_count[line->child_index(1)] > 0),
4119 if (line_cell_count[line->child_index(0)] == 0)
4121 for (
unsigned int c = 0; c < 2; ++c)
4122 Assert(!line->child(c)->has_children(),
4132 lines_to_delete.push_back(line->child(0));
4133 lines_to_delete.push_back(line->child(1));
4135 line->clear_children();
4147 typename std::vector<
4148 typename Triangulation<dim, spacedim>::line_iterator>::iterator
4149 line = lines_to_delete.begin(),
4150 endline = lines_to_delete.end();
4151 for (; line != endline; ++line)
4153 (*line)->clear_user_data();
4154 (*line)->clear_user_flag();
4155 (*line)->clear_used_flag();
4158 typename std::vector<
4159 typename Triangulation<dim, spacedim>::quad_iterator>::iterator
4160 quad = quads_to_delete.begin(),
4161 endquad = quads_to_delete.end();
4162 for (; quad != endquad; ++quad)
4164 (*quad)->clear_user_data();
4165 (*quad)->clear_children();
4166 (*quad)->clear_refinement_case();
4167 (*quad)->clear_user_flag();
4168 (*quad)->clear_used_flag();
4190 template <
int spacedim>
4193 unsigned int & next_unused_vertex,
4194 typename Triangulation<2, spacedim>::raw_line_iterator
4200 const unsigned int dim = 2;
4203 cell->clear_refine_flag();
4300 int new_vertices[9];
4301 for (
unsigned int vertex_no = 0; vertex_no < 4; ++vertex_no)
4302 new_vertices[vertex_no] = cell->vertex_index(vertex_no);
4303 for (
unsigned int line_no = 0; line_no < 4; ++line_no)
4304 if (cell->line(line_no)->has_children())
4305 new_vertices[4 + line_no] =
4306 cell->line(line_no)->child(0)->vertex_index(1);
4315 while (triangulation.
vertices_used[next_unused_vertex] ==
true)
4316 ++next_unused_vertex;
4318 next_unused_vertex < triangulation.
vertices.size(),
4320 "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."));
4323 new_vertices[8] = next_unused_vertex;
4335 if (dim == spacedim)
4338 triangulation.
vertices[next_unused_vertex] = cell->center(
true);
4345 if (cell->user_flag_set())
4348 cell->clear_user_flag();
4349 triangulation.
vertices[next_unused_vertex] =
4350 cell->center(
true,
true);
4360 cell->clear_user_flag();
4364 triangulation.
vertices[next_unused_vertex] = cell->center(
true);
4370 typename Triangulation<dim, spacedim>::raw_line_iterator new_lines[12];
4371 unsigned int lmin = 8;
4372 unsigned int lmax = 12;
4379 for (
unsigned int l = lmin; l < lmax; ++l)
4381 while (next_unused_line->used() ==
true)
4383 new_lines[l] = next_unused_line;
4387 new_lines[l]->used() ==
false,
4389 "Internal error: We want to use a cell during refinement that should be unused, but turns out not to be."));
4403 for (
unsigned int face_no = 0;
4404 face_no < GeometryInfo<dim>::faces_per_cell;
4406 for (
unsigned int c = 0; c < 2; ++c, ++l)
4407 new_lines[l] = cell->line(face_no)->child(c);
4412 new_vertices[6], new_vertices[8]));
4415 new_vertices[8], new_vertices[7]));
4418 new_vertices[4], new_vertices[8]));
4421 new_vertices[8], new_vertices[5]));
4430 new_lines[0] = cell->line(0);
4431 new_lines[1] = cell->line(1);
4432 new_lines[2] = cell->line(2)->child(0);
4433 new_lines[3] = cell->line(2)->child(1);
4434 new_lines[4] = cell->line(3)->child(0);
4435 new_lines[5] = cell->line(3)->child(1);
4438 new_vertices[6], new_vertices[7]));
4448 new_lines[0] = cell->line(0)->child(0);
4449 new_lines[1] = cell->line(0)->child(1);
4450 new_lines[2] = cell->line(1)->child(0);
4451 new_lines[3] = cell->line(1)->child(1);
4452 new_lines[4] = cell->line(2);
4453 new_lines[5] = cell->line(3);
4456 new_vertices[4], new_vertices[5]));
4459 for (
unsigned int l = lmin; l < lmax; ++l)
4461 new_lines[l]->set_used_flag();
4462 new_lines[l]->clear_user_flag();
4464 new_lines[l]->clear_children();
4466 new_lines[l]->set_boundary_id_internal(
4468 new_lines[l]->set_manifold_id(cell->manifold_id());
4475 while (next_unused_cell->used() ==
true)
4479 for (
unsigned int i = 0; i < n_children; ++i)
4482 next_unused_cell->used() ==
false,
4484 "Internal error: We want to use a cell during refinement that should be unused, but turns out not to be."));
4485 subcells[i] = next_unused_cell;
4487 if (i % 2 == 1 && i < n_children - 1)
4488 while (next_unused_cell->used() ==
true)
4508 new_lines[0]->index(),
4509 new_lines[8]->index(),
4510 new_lines[4]->index(),
4511 new_lines[10]->index()));
4514 new_lines[8]->index(),
4515 new_lines[2]->index(),
4516 new_lines[5]->index(),
4517 new_lines[11]->index()));
4520 new_lines[1]->index(),
4521 new_lines[9]->index(),
4522 new_lines[10]->index(),
4523 new_lines[6]->index()));
4526 new_lines[9]->index(),
4527 new_lines[3]->index(),
4528 new_lines[11]->index(),
4529 new_lines[7]->index()));
4547 new_lines[0]->index(),
4548 new_lines[6]->index(),
4549 new_lines[2]->index(),
4550 new_lines[4]->index()));
4553 new_lines[6]->index(),
4554 new_lines[1]->index(),
4555 new_lines[3]->index(),
4556 new_lines[5]->index()));
4575 new_lines[0]->index(),
4576 new_lines[2]->index(),
4577 new_lines[4]->index(),
4578 new_lines[6]->index()));
4581 new_lines[1]->index(),
4582 new_lines[3]->index(),
4583 new_lines[6]->index(),
4584 new_lines[5]->index()));
4589 for (
unsigned int i = 0; i < n_children; ++i)
4591 subcells[i]->set_used_flag();
4592 subcells[i]->clear_refine_flag();
4593 subcells[i]->clear_user_flag();
4594 subcells[i]->clear_user_data();
4595 subcells[i]->clear_children();
4598 subcells[i]->set_material_id(cell->material_id());
4599 subcells[i]->set_manifold_id(cell->manifold_id());
4600 subcells[i]->set_subdomain_id(subdomainid);
4603 subcells[i]->set_parent(cell->index());
4611 for (
unsigned int i = 0; i < n_children / 2; ++i)
4612 cell->set_children(2 * i, subcells[2 * i]->index());
4614 cell->set_refinement_case(ref_case);
4622 for (
unsigned int c = 0; c < n_children; ++c)
4623 cell->child(c)->set_direction_flag(cell->direction_flag());
4632 template <
int spacedim>
4637 const unsigned int dim = 1;
4645 endc = triangulation.
end();
4646 for (; cell != endc; ++cell)
4648 if (cell->refine_flag_set())
4650 triangulation.
levels.push_back(
4651 std_cxx14::make_unique<
4662 unsigned int needed_vertices = 0;
4663 for (
int level = triangulation.
levels.size() - 2; level >= 0; --level)
4667 unsigned int flagged_cells = 0;
4671 for (; acell != aendc; ++acell)
4672 if (acell->refine_flag_set())
4677 const unsigned int used_cells = std::count_if(
4678 triangulation.
levels[level + 1]->cells.used.begin(),
4679 triangulation.
levels[level + 1]->cells.used.end(),
4680 std::bind(std::equal_to<bool>(), std::placeholders::_1,
true));
4685 triangulation.
levels[level + 1]->reserve_space(
4692 triangulation.
levels[level + 1]->cells.reserve_space(
4695 needed_vertices += flagged_cells;
4700 needed_vertices += std::count_if(triangulation.
vertices_used.begin(),
4702 std::bind(std::equal_to<bool>(),
4703 std::placeholders::_1,
4708 if (needed_vertices > triangulation.
vertices.size())
4719 unsigned int next_unused_vertex = 0;
4721 for (
int level = triangulation.
levels.size() - 2; level >= 0; --level)
4728 next_unused_cell = triangulation.
begin_raw(level + 1);
4730 for (; (cell != endc) && (cell->level() == level); ++cell)
4731 if (cell->refine_flag_set())
4734 cell->clear_refine_flag();
4740 ++next_unused_vertex;
4742 next_unused_vertex < triangulation.
vertices.size(),
4744 "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."));
4749 triangulation.
vertices[next_unused_vertex] =
4759 while (next_unused_cell->used() ==
true)
4761 first_child = next_unused_cell;
4762 first_child->set_used_flag();
4763 first_child->clear_user_data();
4766 next_unused_cell->used() ==
false,
4768 "Internal error: We want to use a cell during refinement that should be unused, but turns out not to be."));
4769 second_child = next_unused_cell;
4770 second_child->set_used_flag();
4771 second_child->clear_user_data();
4776 cell->set_children(0, first_child->index());
4777 first_child->clear_children();
4780 cell->vertex_index(0), next_unused_vertex));
4781 first_child->set_material_id(cell->material_id());
4782 first_child->set_manifold_id(cell->manifold_id());
4783 first_child->set_subdomain_id(subdomainid);
4784 first_child->set_direction_flag(cell->direction_flag());
4786 first_child->set_parent(cell->index());
4790 first_child->face(1)->set_manifold_id(cell->manifold_id());
4795 first_child->set_neighbor(1, second_child);
4797 first_child->set_neighbor(0, cell->neighbor(0));
4798 else if (cell->neighbor(0)->active())
4804 Assert(cell->neighbor(0)->level() <= cell->level(),
4806 first_child->set_neighbor(0, cell->neighbor(0));
4812 const unsigned int nbnb = cell->neighbor_of_neighbor(0);
4813 first_child->set_neighbor(0,
4814 cell->neighbor(0)->child(nbnb));
4819 left_neighbor = cell->neighbor(0);
4820 while (left_neighbor->has_children())
4822 left_neighbor = left_neighbor->child(nbnb);
4823 left_neighbor->set_neighbor(nbnb, first_child);
4828 second_child->clear_children();
4831 next_unused_vertex, cell->vertex_index(1)));
4832 second_child->set_neighbor(0, first_child);
4833 second_child->set_material_id(cell->material_id());
4834 second_child->set_manifold_id(cell->manifold_id());
4835 second_child->set_subdomain_id(subdomainid);
4836 second_child->set_direction_flag(cell->direction_flag());
4839 second_child->set_neighbor(1, cell->neighbor(1));
4840 else if (cell->neighbor(1)->active())
4842 Assert(cell->neighbor(1)->level() <= cell->level(),
4844 second_child->set_neighbor(1, cell->neighbor(1));
4849 const unsigned int nbnb = cell->neighbor_of_neighbor(1);
4850 second_child->set_neighbor(
4851 1, cell->neighbor(1)->child(nbnb));
4854 right_neighbor = cell->neighbor(1);
4855 while (right_neighbor->has_children())
4857 right_neighbor = right_neighbor->child(nbnb);
4858 right_neighbor->set_neighbor(nbnb, second_child);
4862 triangulation.
signals.post_refinement_on_cell(cell);
4878 template <
int spacedim>
4881 const bool check_for_distorted_cells)
4883 const unsigned int dim = 2;
4893 endc = triangulation.
end();
4894 for (; cell != endc; ++cell)
4896 if (cell->refine_flag_set())
4898 triangulation.
levels.push_back(
4899 std_cxx14::make_unique<
4909 for (
typename Triangulation<dim, spacedim>::line_iterator line =
4914 line->clear_user_flag();
4915 line->clear_user_data();
4920 unsigned int n_single_lines = 0;
4924 unsigned int n_lines_in_pairs = 0;
4930 unsigned int needed_vertices = 0;
4931 for (
int level = triangulation.
levels.size() - 2; level >= 0; --level)
4935 unsigned int needed_cells = 0;
4940 for (; cell != endc; ++cell)
4941 if (cell->refine_flag_set())
4952 n_single_lines += 4;
4964 n_single_lines += 1;
4972 for (
unsigned int line_no = 0;
4973 line_no < GeometryInfo<dim>::faces_per_cell;
4977 cell->refine_flag_set(), line_no) ==
4980 typename Triangulation<dim, spacedim>::line_iterator
4981 line = cell->line(line_no);
4982 if (line->has_children() ==
false)
4984 line->set_user_flag();
4995 if (line->at_boundary())
4998 line->set_user_index(line->boundary_id());
5003 line->set_user_index(cell->material_id());
5012 const unsigned int used_cells = std::count_if(
5013 triangulation.
levels[level + 1]->cells.used.begin(),
5014 triangulation.
levels[level + 1]->cells.used.end(),
5015 std::bind(std::equal_to<bool>(), std::placeholders::_1,
true));
5021 triangulation.
levels[level + 1]->reserve_space(
5022 used_cells + needed_cells, 2, spacedim);
5026 triangulation.
levels[level + 1]->cells.reserve_space(needed_cells,
5031 for (
typename Triangulation<dim, spacedim>::line_iterator line =
5035 if (line->user_flag_set())
5038 n_lines_in_pairs += 2;
5039 needed_vertices += 1;
5048 triangulation.
faces->lines.reserve_space(n_lines_in_pairs, 0);
5051 needed_vertices += std::count_if(triangulation.
vertices_used.begin(),
5053 std::bind(std::equal_to<bool>(),
5054 std::placeholders::_1,
5059 if (needed_vertices > triangulation.
vertices.size())
5070 unsigned int next_unused_vertex = 0;
5077 typename Triangulation<dim, spacedim>::active_line_iterator
5080 typename Triangulation<dim, spacedim>::raw_line_iterator
5083 for (; line != endl; ++line)
5084 if (line->user_flag_set())
5092 ++next_unused_vertex;
5094 next_unused_vertex < triangulation.
vertices.size(),
5096 "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."));
5099 if (spacedim == dim)
5106 triangulation.
vertices[next_unused_vertex] =
5116 triangulation.
vertices[next_unused_vertex] =
5118 .get_new_point_on_line(line);
5120 triangulation.
vertices[next_unused_vertex] =
5126 bool pair_found =
false;
5128 for (; next_unused_line != endl; ++next_unused_line)
5129 if (!next_unused_line->used() &&
5130 !(++next_unused_line)->used())
5143 line->set_children(0, next_unused_line->index());
5146 const typename Triangulation<dim, spacedim>::raw_line_iterator
5147 children[2] = {next_unused_line, ++next_unused_line};
5151 children[0]->used() ==
false,
5153 "Internal error: We want to use a cell during refinement that should be unused, but turns out not to be."));
5155 children[1]->used() ==
false,
5157 "Internal error: We want to use a cell during refinement that should be unused, but turns out not to be."));
5161 line->vertex_index(0), next_unused_vertex));
5164 next_unused_vertex, line->vertex_index(1)));
5166 children[0]->set_used_flag();
5167 children[1]->set_used_flag();
5168 children[0]->clear_children();
5169 children[1]->clear_children();
5172 children[0]->clear_user_flag();
5173 children[1]->clear_user_flag();
5176 children[0]->set_boundary_id_internal(line->boundary_id());
5177 children[1]->set_boundary_id_internal(line->boundary_id());
5179 children[0]->set_manifold_id(line->manifold_id());
5180 children[1]->set_manifold_id(line->manifold_id());
5184 line->clear_user_flag();
5193 triangulation.
faces->lines.reserve_space(0, n_single_lines);
5196 cells_with_distorted_children;
5200 typename Triangulation<dim, spacedim>::raw_line_iterator
5204 level < static_cast<int>(triangulation.
levels.size()) - 1;
5214 next_unused_cell = triangulation.
begin_raw(level + 1);
5216 for (; cell != endc; ++cell)
5217 if (cell->refine_flag_set())
5225 if (cell->at_boundary())
5226 cell->set_user_flag();
5230 create_children(triangulation,
5236 if ((check_for_distorted_cells ==
true) &&
5237 has_distorted_children(
5239 std::integral_constant<int, dim>(),
5240 std::integral_constant<int, spacedim>()))
5244 triangulation.
signals.post_refinement_on_cell(cell);
5248 return cells_with_distorted_children;
5256 template <
int spacedim>
5259 const bool check_for_distorted_cells)
5261 const unsigned int dim = 3;
5277 endc = triangulation.
end();
5278 for (; cell != endc; ++cell)
5280 if (cell->refine_flag_set())
5282 triangulation.
levels.push_back(
5283 std_cxx14::make_unique<
5293 triangulation.
faces->quads.clear_user_data();
5295 for (
typename Triangulation<dim, spacedim>::line_iterator line =
5299 line->clear_user_flag();
5300 for (
typename Triangulation<dim, spacedim>::quad_iterator quad =
5304 quad->clear_user_flag();
5327 unsigned int needed_vertices = 0;
5328 unsigned int needed_lines_single = 0;
5329 unsigned int needed_quads_single = 0;
5330 unsigned int needed_lines_pair = 0;
5331 unsigned int needed_quads_pair = 0;
5332 for (
int level = triangulation.
levels.size() - 2; level >= 0; --level)
5336 unsigned int new_cells = 0;
5341 for (; acell != aendc; ++acell)
5342 if (acell->refine_flag_set())
5352 ++needed_quads_single;
5360 ++needed_lines_single;
5361 needed_quads_single += 4;
5368 needed_lines_single += 6;
5369 needed_quads_single += 12;
5383 for (
unsigned int face = 0;
5384 face < GeometryInfo<dim>::faces_per_cell;
5388 aface = acell->face(face);
5395 acell->face_orientation(face),
5396 acell->face_flip(face),
5397 acell->face_rotation(face));
5402 if (face_ref_case ==
5405 if (aface->number_of_children() < 4)
5408 aface->set_user_flag();
5410 else if (aface->refinement_case() != face_ref_case)
5421 Assert(aface->refinement_case() ==
5423 dim - 1>::isotropic_refinement ||
5424 aface->refinement_case() ==
5427 aface->set_user_index(face_ref_case);
5433 for (
unsigned int line = 0;
5434 line < GeometryInfo<dim>::lines_per_cell;
5438 !acell->line(line)->has_children())
5439 acell->line(line)->set_user_flag();
5445 const unsigned int used_cells = std::count_if(
5446 triangulation.
levels[level + 1]->cells.used.begin(),
5447 triangulation.
levels[level + 1]->cells.used.end(),
5448 std::bind(std::equal_to<bool>(), std::placeholders::_1,
true));
5454 triangulation.
levels[level + 1]->reserve_space(
5455 used_cells + new_cells, 3, spacedim);
5458 triangulation.
levels[level + 1]->cells.reserve_space(new_cells);
5462 for (
typename Triangulation<dim, spacedim>::quad_iterator quad =
5467 if (quad->user_flag_set())
5473 needed_quads_pair += 4;
5474 needed_lines_pair += 4;
5475 needed_vertices += 1;
5477 if (quad->user_index())
5481 needed_quads_pair += 2;
5482 needed_lines_single += 1;
5495 if (quad->has_children())
5497 Assert(quad->refinement_case() ==
5500 if ((face_refinement_cases[quad->user_index()] ==
5502 (quad->child(0)->line_index(1) + 1 !=
5503 quad->child(2)->line_index(1))) ||
5504 (face_refinement_cases[quad->user_index()] ==
5506 (quad->child(0)->line_index(3) + 1 !=
5507 quad->child(1)->line_index(3))))
5508 needed_lines_pair += 2;
5513 for (
typename Triangulation<dim, spacedim>::line_iterator line =
5517 if (line->user_flag_set())
5519 needed_lines_pair += 2;
5520 needed_vertices += 1;
5524 triangulation.
faces->lines.reserve_space(needed_lines_pair,
5525 needed_lines_single);
5527 triangulation.
faces->quads.reserve_space(needed_quads_pair,
5528 needed_quads_single);
5532 needed_vertices += std::count_if(triangulation.
vertices_used.begin(),
5534 std::bind(std::equal_to<bool>(),
5535 std::placeholders::_1,
5540 if (needed_vertices > triangulation.
vertices.size())
5560 cell != triangulation.
end();
5562 if (!cell->refine_flag_set())
5563 for (
unsigned int line = 0;
5564 line < GeometryInfo<dim>::lines_per_cell;
5566 if (cell->line(line)->has_children())
5567 for (
unsigned int c = 0; c < 2; ++c)
5568 Assert(cell->line(line)->child(c)->user_flag_set() ==
false,
5580 unsigned int next_unused_vertex = 0;
5586 typename Triangulation<dim, spacedim>::active_line_iterator
5589 typename Triangulation<dim, spacedim>::raw_line_iterator
5592 for (; line != endl; ++line)
5593 if (line->user_flag_set())
5601 ++next_unused_vertex;
5603 next_unused_vertex < triangulation.
vertices.size(),
5605 "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."));
5608 triangulation.
vertices[next_unused_vertex] =
5615 triangulation.
faces->lines.next_free_pair_object(
5623 line->set_children(0, next_unused_line->index());
5626 const typename Triangulation<dim, spacedim>::raw_line_iterator
5627 children[2] = {next_unused_line, ++next_unused_line};
5632 children[0]->used() ==
false,
5634 "Internal error: We want to use a cell during refinement that should be unused, but turns out not to be."));
5636 children[1]->used() ==
false,
5638 "Internal error: We want to use a cell during refinement that should be unused, but turns out not to be."));
5642 line->vertex_index(0), next_unused_vertex));
5645 next_unused_vertex, line->vertex_index(1)));
5647 children[0]->set_used_flag();
5648 children[1]->set_used_flag();
5649 children[0]->clear_children();
5650 children[1]->clear_children();
5653 children[0]->clear_user_flag();
5654 children[1]->clear_user_flag();
5656 children[0]->set_boundary_id_internal(line->boundary_id());
5657 children[1]->set_boundary_id_internal(line->boundary_id());
5659 children[0]->set_manifold_id(line->manifold_id());
5660 children[1]->set_manifold_id(line->manifold_id());
5665 line->clear_user_flag();
5700 typename Triangulation<dim, spacedim>::quad_iterator
5703 typename Triangulation<dim, spacedim>::raw_line_iterator
5705 typename Triangulation<dim, spacedim>::raw_quad_iterator
5708 for (; quad != endq; ++quad)
5710 if (quad->user_index())
5713 face_refinement_cases[quad->user_index()];
5722 if (aniso_quad_ref_case == quad->refinement_case())
5725 Assert(quad->refinement_case() ==
5727 quad->refinement_case() ==
5732 Assert(quad->user_index() ==
5734 quad->user_index() ==
5739 typename Triangulation<dim, spacedim>::raw_line_iterator
5743 triangulation.
faces->lines.next_free_single_object(
5746 new_line->used() ==
false,
5748 "Internal error: We want to use a cell during refinement that should be unused, but turns out not to be."));
5763 unsigned int vertex_indices[2];
5767 quad->line(2)->child(0)->vertex_index(1);
5769 quad->line(3)->child(0)->vertex_index(1);
5774 quad->line(0)->child(0)->vertex_index(1);
5776 quad->line(1)->child(0)->vertex_index(1);
5781 vertex_indices[0], vertex_indices[1]));
5782 new_line->set_used_flag();
5783 new_line->clear_user_flag();
5785 new_line->clear_children();
5786 new_line->set_boundary_id_internal(quad->boundary_id());
5787 new_line->set_manifold_id(quad->manifold_id());
5795 const unsigned int index[2][2] = {
5801 typename Triangulation<dim, spacedim>::raw_quad_iterator
5805 triangulation.
faces->quads.next_free_pair_object(
5807 new_quads[0] = next_unused_quad;
5809 new_quads[0]->used() ==
false,
5811 "Internal error: We want to use a cell during refinement that should be unused, but turns out not to be."));
5814 new_quads[1] = next_unused_quad;
5816 new_quads[1]->used() ==
false,
5818 "Internal error: We want to use a cell during refinement that should be unused, but turns out not to be."));
5825 quad->line_index(0),
5828 ->child(index[0][quad->line_orientation(2)])
5831 ->child(index[0][quad->line_orientation(3)])
5836 quad->line_index(1),
5838 ->child(index[1][quad->line_orientation(2)])
5841 ->child(index[1][quad->line_orientation(3)])
5849 ->child(index[0][quad->line_orientation(0)])
5852 ->child(index[0][quad->line_orientation(1)])
5854 quad->line_index(2),
5855 new_line->index()));
5859 ->child(index[1][quad->line_orientation(0)])
5862 ->child(index[1][quad->line_orientation(1)])
5865 quad->line_index(3)));
5868 for (
unsigned int i = 0; i < 2; ++i)
5870 new_quads[i]->set_used_flag();
5871 new_quads[i]->clear_user_flag();
5873 new_quads[i]->clear_children();
5874 new_quads[i]->set_boundary_id_internal(
5875 quad->boundary_id());
5876 new_quads[i]->set_manifold_id(quad->manifold_id());
5880 for (
unsigned int j = 0;
5881 j < GeometryInfo<dim>::lines_per_face;
5883 new_quads[i]->set_line_orientation(j,
true);
5889 new_quads[0]->set_line_orientation(
5890 0, quad->line_orientation(0));
5891 new_quads[0]->set_line_orientation(
5892 2, quad->line_orientation(2));
5893 new_quads[1]->set_line_orientation(
5894 1, quad->line_orientation(1));
5895 new_quads[1]->set_line_orientation(
5896 3, quad->line_orientation(3));
5899 new_quads[0]->set_line_orientation(
5900 3, quad->line_orientation(3));
5901 new_quads[1]->set_line_orientation(
5902 2, quad->line_orientation(2));
5906 new_quads[0]->set_line_orientation(
5907 1, quad->line_orientation(1));
5908 new_quads[1]->set_line_orientation(
5909 0, quad->line_orientation(0));
5915 if (quad->refinement_case() ==
5936 typename Triangulation<dim, spacedim>::line_iterator
5938 if (aniso_quad_ref_case ==
5941 old_child[0] = quad->child(0)->line(1);
5942 old_child[1] = quad->child(2)->line(1);
5946 Assert(aniso_quad_ref_case ==
5950 old_child[0] = quad->child(0)->line(3);
5951 old_child[1] = quad->child(1)->line(3);
5954 if (old_child[0]->index() + 1 != old_child[1]->index())
5960 spacedim>::raw_line_iterator
5963 new_child[0] = new_child[1] =
5964 triangulation.
faces->lines.next_free_pair_object(
5968 new_child[0]->set_used_flag();
5969 new_child[1]->set_used_flag();
5971 const int old_index_0 = old_child[0]->index(),
5972 old_index_1 = old_child[1]->index(),
5973 new_index_0 = new_child[0]->index(),
5974 new_index_1 = new_child[1]->index();
5978 for (
unsigned int q = 0;
5979 q < triangulation.
faces->quads.cells.size();
5981 for (
unsigned int l = 0;
5982 l < GeometryInfo<dim>::lines_per_face;
5985 const int this_index =
5986 triangulation.
faces->quads.cells[q].face(l);
5987 if (this_index == old_index_0)
5988 triangulation.
faces->quads.cells[q]
5989 .set_face(l, new_index_0);
5990 else if (this_index == old_index_1)
5991 triangulation.
faces->quads.cells[q]
5992 .set_face(l, new_index_1);
5996 for (
unsigned int i = 0; i < 2; ++i)
5998 Assert(!old_child[i]->has_children(),
6004 old_child[i]->vertex_index(
6006 new_child[i]->set_boundary_id_internal(
6007 old_child[i]->boundary_id());
6008 new_child[i]->set_manifold_id(
6009 old_child[i]->manifold_id());
6010 new_child[i]->set_user_index(
6011 old_child[i]->user_index());
6012 if (old_child[i]->user_flag_set())
6013 new_child[i]->set_user_flag();
6015 new_child[i]->clear_user_flag();
6017 new_child[i]->clear_children();
6019 old_child[i]->clear_user_flag();
6020 old_child[i]->clear_user_index();
6021 old_child[i]->clear_used_flag();
6027 if (aniso_quad_ref_case ==
6030 new_line->set_children(
6031 0, quad->child(0)->line_index(1));
6032 Assert(new_line->child(1) ==
6033 quad->child(2)->line(1),
6070 quad_iterator switch_1 = quad->child(1),
6071 switch_2 = quad->child(2);
6072 const int switch_1_index = switch_1->index();
6073 const int switch_2_index = switch_2->index();
6074 for (
unsigned int l = 0;
6075 l < triangulation.
levels.size();
6077 for (
unsigned int h = 0;
6079 triangulation.
levels[l]->cells.cells.size();
6081 for (
unsigned int q = 0;
6082 q < GeometryInfo<dim>::faces_per_cell;
6085 const int face_index =
6089 if (face_index == switch_1_index)
6092 .set_face(q, switch_2_index);
6093 else if (face_index == switch_2_index)
6096 .set_face(q, switch_1_index);
6100 const unsigned int switch_1_lines[4] = {
6101 switch_1->line_index(0),
6102 switch_1->line_index(1),
6103 switch_1->line_index(2),
6104 switch_1->line_index(3)};
6105 const bool switch_1_line_orientations[4] = {
6106 switch_1->line_orientation(0),
6107 switch_1->line_orientation(1),
6108 switch_1->line_orientation(2),
6109 switch_1->line_orientation(3)};
6111 switch_1->boundary_id();
6112 const unsigned int switch_1_user_index =
6113 switch_1->user_index();
6114 const bool switch_1_user_flag =
6115 switch_1->user_flag_set();
6117 switch_1_refinement_case =
6118 switch_1->refinement_case();
6119 const int switch_1_first_child_pair =
6120 (switch_1_refinement_case ?
6121 switch_1->child_index(0) :
6123 const int switch_1_second_child_pair =
6124 (switch_1_refinement_case ==
6126 switch_1->child_index(2) :
6131 2>(switch_2->line_index(0),
6132 switch_2->line_index(1),
6133 switch_2->line_index(2),
6134 switch_2->line_index(3)));
6135 switch_1->set_line_orientation(
6136 0, switch_2->line_orientation(0));
6137 switch_1->set_line_orientation(
6138 1, switch_2->line_orientation(1));
6139 switch_1->set_line_orientation(
6140 2, switch_2->line_orientation(2));
6141 switch_1->set_line_orientation(
6142 3, switch_2->line_orientation(3));
6143 switch_1->set_boundary_id_internal(
6144 switch_2->boundary_id());
6145 switch_1->set_manifold_id(switch_2->manifold_id());
6146 switch_1->set_user_index(switch_2->user_index());
6147 if (switch_2->user_flag_set())
6148 switch_1->set_user_flag();
6150 switch_1->clear_user_flag();
6151 switch_1->clear_refinement_case();
6152 switch_1->set_refinement_case(
6153 switch_2->refinement_case());
6154 switch_1->clear_children();
6155 if (switch_2->refinement_case())
6156 switch_1->set_children(0,
6157 switch_2->child_index(0));
6158 if (switch_2->refinement_case() ==
6160 switch_1->set_children(2,
6161 switch_2->child_index(2));
6165 2>(switch_1_lines[0],
6168 switch_1_lines[3]));
6169 switch_2->set_line_orientation(
6170 0, switch_1_line_orientations[0]);
6171 switch_2->set_line_orientation(
6172 1, switch_1_line_orientations[1]);
6173 switch_2->set_line_orientation(
6174 2, switch_1_line_orientations[2]);
6175 switch_2->set_line_orientation(
6176 3, switch_1_line_orientations[3]);
6177 switch_2->set_boundary_id_internal(
6178 switch_1_boundary_id);
6179 switch_2->set_manifold_id(switch_1->manifold_id());
6180 switch_2->set_user_index(switch_1_user_index);
6181 if (switch_1_user_flag)
6182 switch_2->set_user_flag();
6184 switch_2->clear_user_flag();
6185 switch_2->clear_refinement_case();
6186 switch_2->set_refinement_case(
6187 switch_1_refinement_case);
6188 switch_2->clear_children();
6189 switch_2->set_children(0,
6190 switch_1_first_child_pair);
6191 switch_2->set_children(2,
6192 switch_1_second_child_pair);
6194 new_quads[0]->set_refinement_case(
6196 new_quads[0]->set_children(0, quad->child_index(0));
6197 new_quads[1]->set_refinement_case(
6199 new_quads[1]->set_children(0, quad->child_index(2));
6203 new_quads[0]->set_refinement_case(
6205 new_quads[0]->set_children(0, quad->child_index(0));
6206 new_quads[1]->set_refinement_case(
6208 new_quads[1]->set_children(0, quad->child_index(2));
6209 new_line->set_children(
6210 0, quad->child(0)->line_index(3));
6211 Assert(new_line->child(1) ==
6212 quad->child(1)->line(3),
6215 quad->clear_children();
6219 quad->set_children(0, new_quads[0]->index());
6221 quad->set_refinement_case(aniso_quad_ref_case);
6228 if (quad->user_flag_set())
6240 ++next_unused_vertex;
6242 next_unused_vertex < triangulation.
vertices.size(),
6244 "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."));
6252 quad->refinement_case();
6258 quad->child(0)->set_user_index(
6260 quad->child(1)->set_user_index(
6263 typename Triangulation<dim, spacedim>::line_iterator
6266 middle_line = quad->child(0)->line(1);
6268 middle_line = quad->child(0)->line(3);
6275 if (!middle_line->has_children())
6282 triangulation.
vertices[next_unused_vertex] =
6283 middle_line->center(
true);
6290 triangulation.
faces->lines.next_free_pair_object(
6295 middle_line->set_children(
6296 0, next_unused_line->index());
6300 raw_line_iterator children[2] = {
6301 next_unused_line, ++next_unused_line};
6307 children[0]->used() ==
false,
6309 "Internal error: We want to use a cell during refinement that should be unused, but turns out not to be."));
6311 children[1]->used() ==
false,
6313 "Internal error: We want to use a cell during refinement that should be unused, but turns out not to be."));
6317 1>(middle_line->vertex_index(0),
6318 next_unused_vertex));
6321 1>(next_unused_vertex,
6322 middle_line->vertex_index(1)));
6324 children[0]->set_used_flag();
6325 children[1]->set_used_flag();
6326 children[0]->clear_children();
6327 children[1]->clear_children();
6330 children[0]->clear_user_flag();
6331 children[1]->clear_user_flag();
6333 children[0]->set_boundary_id_internal(
6334 middle_line->boundary_id());
6335 children[1]->set_boundary_id_internal(
6336 middle_line->boundary_id());
6338 children[0]->set_manifold_id(
6339 middle_line->manifold_id());
6340 children[1]->set_manifold_id(
6341 middle_line->manifold_id());
6347 quad->clear_user_flag();
6376 triangulation.
vertices[next_unused_vertex] =
6377 quad->center(
true,
true);
6383 typename Triangulation<dim, spacedim>::raw_line_iterator
6386 for (
unsigned int i = 0; i < 4; ++i)
6396 triangulation.
faces->lines.next_free_pair_object(
6399 new_lines[i] = next_unused_line;
6403 new_lines[i]->used() ==
false,
6405 "Internal error: We want to use a cell during refinement that should be unused, but turns out not to be."));
6425 const unsigned int vertex_indices[5] = {
6426 quad->line(0)->child(0)->vertex_index(1),
6427 quad->line(1)->child(0)->vertex_index(1),
6428 quad->line(2)->child(0)->vertex_index(1),
6429 quad->line(3)->child(0)->vertex_index(1),
6430 next_unused_vertex};
6434 vertex_indices[2], vertex_indices[4]));
6437 vertex_indices[4], vertex_indices[3]));
6440 vertex_indices[0], vertex_indices[4]));
6443 vertex_indices[4], vertex_indices[1]));
6445 for (
unsigned int i = 0; i < 4; ++i)
6447 new_lines[i]->set_used_flag();
6448 new_lines[i]->clear_user_flag();
6450 new_lines[i]->clear_children();
6451 new_lines[i]->set_boundary_id_internal(
6452 quad->boundary_id());
6453 new_lines[i]->set_manifold_id(quad->manifold_id());
6472 const unsigned int index[2][2] = {
6476 const int line_indices[12] = {
6478 ->child(index[0][quad->line_orientation(0)])
6481 ->child(index[1][quad->line_orientation(0)])
6484 ->child(index[0][quad->line_orientation(1)])
6487 ->child(index[1][quad->line_orientation(1)])
6490 ->child(index[0][quad->line_orientation(2)])
6493 ->child(index[1][quad->line_orientation(2)])
6496 ->child(index[0][quad->line_orientation(3)])
6499 ->child(index[1][quad->line_orientation(3)])
6501 new_lines[0]->index(),
6502 new_lines[1]->index(),
6503 new_lines[2]->index(),
6504 new_lines[3]->index()};
6509 typename Triangulation<dim, spacedim>::raw_quad_iterator
6513 triangulation.
faces->quads.next_free_pair_object(
6516 new_quads[0] = next_unused_quad;
6518 new_quads[0]->used() ==
false,
6520 "Internal error: We want to use a cell during refinement that should be unused, but turns out not to be."));
6523 new_quads[1] = next_unused_quad;
6525 new_quads[1]->used() ==
false,
6527 "Internal error: We want to use a cell during refinement that should be unused, but turns out not to be."));
6530 triangulation.
faces->quads.next_free_pair_object(
6532 new_quads[2] = next_unused_quad;
6534 new_quads[2]->used() ==
false,
6536 "Internal error: We want to use a cell during refinement that should be unused, but turns out not to be."));
6539 new_quads[3] = next_unused_quad;
6541 new_quads[3]->used() ==
false,
6543 "Internal error: We want to use a cell during refinement that should be unused, but turns out not to be."));
6546 quad->set_children(0, new_quads[0]->index());
6547 quad->set_children(2, new_quads[2]->index());
6581 for (
unsigned int i = 0; i < 4; ++i)
6583 new_quads[i]->set_used_flag();
6584 new_quads[i]->clear_user_flag();
6586 new_quads[i]->clear_children();
6587 new_quads[i]->set_boundary_id_internal(
6588 quad->boundary_id());
6589 new_quads[i]->set_manifold_id(quad->manifold_id());
6593 for (
unsigned int j = 0;
6594 j < GeometryInfo<dim>::lines_per_face;
6596 new_quads[i]->set_line_orientation(j,
true);
6602 new_quads[0]->set_line_orientation(
6603 0, quad->line_orientation(0));
6604 new_quads[0]->set_line_orientation(
6605 2, quad->line_orientation(2));
6606 new_quads[1]->set_line_orientation(
6607 1, quad->line_orientation(1));
6608 new_quads[1]->set_line_orientation(
6609 2, quad->line_orientation(2));
6610 new_quads[2]->set_line_orientation(
6611 0, quad->line_orientation(0));
6612 new_quads[2]->set_line_orientation(
6613 3, quad->line_orientation(3));
6614 new_quads[3]->set_line_orientation(
6615 1, quad->line_orientation(1));
6616 new_quads[3]->set_line_orientation(
6617 3, quad->line_orientation(3));
6621 quad->clear_user_flag();
6632 cells_with_distorted_children;
6634 for (
unsigned int level = 0; level != triangulation.
levels.size() - 1;
6640 typename Triangulation<dim, spacedim>::active_hex_iterator
6643 typename Triangulation<dim, spacedim>::raw_hex_iterator
6646 for (; hex != endh; ++hex)
6647 if (hex->refine_flag_set())
6655 hex->clear_refine_flag();
6656 hex->set_refinement_case(ref_case);
6667 unsigned int n_new_lines = 0;
6668 unsigned int n_new_quads = 0;
6669 unsigned int n_new_hexes = 0;
6699 typename Triangulation<dim, spacedim>::raw_line_iterator>
6700 new_lines(n_new_lines);
6701 for (
unsigned int i = 0; i < n_new_lines; ++i)
6704 triangulation.
faces->lines.next_free_single_object(
6708 new_lines[i]->used() ==
false,
6710 "Internal error: We want to use a cell during refinement that should be unused, but turns out not to be."));
6711 new_lines[i]->set_used_flag();
6712 new_lines[i]->clear_user_flag();
6713 new_lines[i]->clear_user_data();
6714 new_lines[i]->clear_children();
6716 new_lines[i]->set_boundary_id_internal(
6720 new_lines[i]->set_manifold_id(hex->manifold_id());
6726 typename Triangulation<dim, spacedim>::raw_quad_iterator>
6727 new_quads(n_new_quads);
6728 for (
unsigned int i = 0; i < n_new_quads; ++i)
6731 triangulation.
faces->quads.next_free_single_object(
6735 new_quads[i]->used() ==
false,
6737 "Internal error: We want to use a cell during refinement that should be unused, but turns out not to be."));
6738 new_quads[i]->set_used_flag();
6739 new_quads[i]->clear_user_flag();
6740 new_quads[i]->clear_user_data();
6741 new_quads[i]->clear_children();
6743 new_quads[i]->set_boundary_id_internal(
6747 new_quads[i]->set_manifold_id(hex->manifold_id());
6750 for (
unsigned int j = 0;
6751 j < GeometryInfo<dim>::lines_per_face;
6753 new_quads[i]->set_line_orientation(j,
true);
6761 typename Triangulation<dim, spacedim>::raw_hex_iterator>
6762 new_hexes(n_new_hexes);
6763 for (
unsigned int i = 0; i < n_new_hexes; ++i)
6767 triangulation.
levels[level + 1]->cells.next_free_hex(
6768 triangulation, level + 1);
6772 new_hexes[i] = next_unused_hex;
6775 new_hexes[i]->used() ==
false,
6777 "Internal error: We want to use a cell during refinement that should be unused, but turns out not to be."));
6778 new_hexes[i]->set_used_flag();
6779 new_hexes[i]->clear_user_flag();
6780 new_hexes[i]->clear_user_data();
6781 new_hexes[i]->clear_children();
6784 new_hexes[i]->set_material_id(hex->material_id());
6785 new_hexes[i]->set_manifold_id(hex->manifold_id());
6786 new_hexes[i]->set_subdomain_id(subdomainid);
6789 new_hexes[i]->set_parent(hex->index());
6801 for (
unsigned int f = 0;
6802 f < GeometryInfo<dim>::faces_per_cell;
6805 new_hexes[i]->set_face_orientation(f,
true);
6806 new_hexes[i]->set_face_flip(f,
false);
6807 new_hexes[i]->set_face_rotation(f,
false);
6811 for (
unsigned int i = 0; i < n_new_hexes / 2; ++i)
6812 hex->set_children(2 * i, new_hexes[2 * i]->index());
6819 const bool f_or[6] = {hex->face_orientation(0),
6820 hex->face_orientation(1),
6821 hex->face_orientation(2),
6822 hex->face_orientation(3),
6823 hex->face_orientation(4),
6824 hex->face_orientation(5)};
6827 const bool f_fl[6] = {hex->face_flip(0),
6835 const bool f_ro[6] = {hex->face_rotation(0),
6836 hex->face_rotation(1),
6837 hex->face_rotation(2),
6838 hex->face_rotation(3),
6839 hex->face_rotation(4),
6840 hex->face_rotation(5)};
6850 const unsigned int child_at_origin[2][2][2] = {
6937 raw_line_iterator lines[4] = {
6938 hex->face(2)->child(0)->line(
6939 (hex->face(2)->refinement_case() ==
6943 hex->face(3)->child(0)->line(
6944 (hex->face(3)->refinement_case() ==
6948 hex->face(4)->child(0)->line(
6949 (hex->face(4)->refinement_case() ==
6953 hex->face(5)->child(0)->line(
6954 (hex->face(5)->refinement_case() ==
6960 unsigned int line_indices[4];
6961 for (
unsigned int i = 0; i < 4; ++i)
6962 line_indices[i] = lines[i]->index();
6969 bool line_orientation[4];
6975 const unsigned int middle_vertices[2] = {
6976 hex->line(2)->child(0)->vertex_index(1),
6977 hex->line(7)->child(0)->vertex_index(1)};
6979 for (
unsigned int i = 0; i < 4; ++i)
6980 if (lines[i]->vertex_index(i % 2) ==
6981 middle_vertices[i % 2])
6982 line_orientation[i] =
true;
6987 Assert(lines[i]->vertex_index((i + 1) % 2) ==
6988 middle_vertices[i % 2],
6990 line_orientation[i] =
false;
7002 new_quads[0]->set_line_orientation(
7003 0, line_orientation[0]);
7004 new_quads[0]->set_line_orientation(
7005 1, line_orientation[1]);
7006 new_quads[0]->set_line_orientation(
7007 2, line_orientation[2]);
7008 new_quads[0]->set_line_orientation(
7009 3, line_orientation[3]);
7042 const int quad_indices[11] = {
7043 new_quads[0]->index(),
7045 hex->face(0)->index(),
7047 hex->face(1)->index(),
7049 hex->face(2)->child_index(
7050 child_at_origin[hex->face(2)->refinement_case() -
7051 1][f_fl[2]][f_ro[2]]),
7052 hex->face(2)->child_index(
7054 child_at_origin[hex->face(2)->refinement_case() -
7055 1][f_fl[2]][f_ro[2]]),
7057 hex->face(3)->child_index(
7058 child_at_origin[hex->face(3)->refinement_case() -
7059 1][f_fl[3]][f_ro[3]]),
7060 hex->face(3)->child_index(
7062 child_at_origin[hex->face(3)->refinement_case() -
7063 1][f_fl[3]][f_ro[3]]),
7065 hex->face(4)->child_index(
7066 child_at_origin[hex->face(4)->refinement_case() -
7067 1][f_fl[4]][f_ro[4]]),
7068 hex->face(4)->child_index(
7070 child_at_origin[hex->face(4)->refinement_case() -
7071 1][f_fl[4]][f_ro[4]]),
7073 hex->face(5)->child_index(
7074 child_at_origin[hex->face(5)->refinement_case() -
7075 1][f_fl[5]][f_ro[5]]),
7076 hex->face(5)->child_index(
7078 child_at_origin[hex->face(5)->refinement_case() -
7079 1][f_fl[5]][f_ro[5]])
7167 raw_line_iterator lines[4] = {
7168 hex->face(0)->child(0)->line(
7169 (hex->face(0)->refinement_case() ==
7173 hex->face(1)->child(0)->line(
7174 (hex->face(1)->refinement_case() ==
7178 hex->face(4)->child(0)->line(
7179 (hex->face(4)->refinement_case() ==
7183 hex->face(5)->child(0)->line(
7184 (hex->face(5)->refinement_case() ==
7190 unsigned int line_indices[4];
7191 for (
unsigned int i = 0; i < 4; ++i)
7192 line_indices[i] = lines[i]->index();
7199 bool line_orientation[4];
7205 const unsigned int middle_vertices[2] = {
7206 hex->line(0)->child(0)->vertex_index(1),
7207 hex->line(5)->child(0)->vertex_index(1)};
7209 for (
unsigned int i = 0; i < 4; ++i)
7210 if (lines[i]->vertex_index(i % 2) ==
7211 middle_vertices[i % 2])
7212 line_orientation[i] =
true;
7216 Assert(lines[i]->vertex_index((i + 1) % 2) ==
7217 middle_vertices[i % 2],
7219 line_orientation[i] =
false;
7231 new_quads[0]->set_line_orientation(
7232 0, line_orientation[2]);
7233 new_quads[0]->set_line_orientation(
7234 1, line_orientation[3]);
7235 new_quads[0]->set_line_orientation(
7236 2, line_orientation[0]);
7237 new_quads[0]->set_line_orientation(
7238 3, line_orientation[1]);
7271 const int quad_indices[11] = {
7272 new_quads[0]->index(),
7274 hex->face(0)->child_index(
7275 child_at_origin[hex->face(0)->refinement_case() -
7276 1][f_fl[0]][f_ro[0]]),
7277 hex->face(0)->child_index(
7279 child_at_origin[hex->face(0)->refinement_case() -
7280 1][f_fl[0]][f_ro[0]]),
7282 hex->face(1)->child_index(
7283 child_at_origin[hex->face(1)->refinement_case() -
7284 1][f_fl[1]][f_ro[1]]),
7285 hex->face(1)->child_index(
7287 child_at_origin[hex->face(1)->refinement_case() -
7288 1][f_fl[1]][f_ro[1]]),
7290 hex->face(2)->index(),
7292 hex->face(3)->index(),
7294 hex->face(4)->child_index(
7295 child_at_origin[hex->face(4)->refinement_case() -
7296 1][f_fl[4]][f_ro[4]]),
7297 hex->face(4)->child_index(
7299 child_at_origin[hex->face(4)->refinement_case() -
7300 1][f_fl[4]][f_ro[4]]),
7302 hex->face(5)->child_index(
7303 child_at_origin[hex->face(5)->refinement_case() -
7304 1][f_fl[5]][f_ro[5]]),
7305 hex->face(5)->child_index(
7307 child_at_origin[hex->face(5)->refinement_case() -
7308 1][f_fl[5]][f_ro[5]])
7398 raw_line_iterator lines[4] = {
7399 hex->face(0)->child(0)->line(
7400 (hex->face(0)->refinement_case() ==
7404 hex->face(1)->child(0)->line(
7405 (hex->face(1)->refinement_case() ==
7409 hex->face(2)->child(0)->line(
7410 (hex->face(2)->refinement_case() ==
7414 hex->face(3)->child(0)->line(
7415 (hex->face(3)->refinement_case() ==
7421 unsigned int line_indices[4];
7422 for (
unsigned int i = 0; i < 4; ++i)
7423 line_indices[i] = lines[i]->index();
7430 bool line_orientation[4];
7436 const unsigned int middle_vertices[2] = {
7437 middle_vertex_index<dim, spacedim>(hex->line(8)),
7438 middle_vertex_index<dim, spacedim>(hex->line(11))};
7440 for (
unsigned int i = 0; i < 4; ++i)
7441 if (lines[i]->vertex_index(i % 2) ==
7442 middle_vertices[i % 2])
7443 line_orientation[i] =
true;
7447 Assert(lines[i]->vertex_index((i + 1) % 2) ==
7448 middle_vertices[i % 2],
7450 line_orientation[i] =
false;
7462 new_quads[0]->set_line_orientation(
7463 0, line_orientation[0]);
7464 new_quads[0]->set_line_orientation(
7465 1, line_orientation[1]);
7466 new_quads[0]->set_line_orientation(
7467 2, line_orientation[2]);
7468 new_quads[0]->set_line_orientation(
7469 3, line_orientation[3]);
7503 const int quad_indices[11] = {
7504 new_quads[0]->index(),
7506 hex->face(0)->child_index(
7507 child_at_origin[hex->face(0)->refinement_case() -
7508 1][f_fl[0]][f_ro[0]]),
7509 hex->face(0)->child_index(
7511 child_at_origin[hex->face(0)->refinement_case() -
7512 1][f_fl[0]][f_ro[0]]),
7514 hex->face(1)->child_index(
7515 child_at_origin[hex->face(1)->refinement_case() -
7516 1][f_fl[1]][f_ro[1]]),
7517 hex->face(1)->child_index(
7519 child_at_origin[hex->face(1)->refinement_case() -
7520 1][f_fl[1]][f_ro[1]]),
7522 hex->face(2)->child_index(
7523 child_at_origin[hex->face(2)->refinement_case() -
7524 1][f_fl[2]][f_ro[2]]),
7525 hex->face(2)->child_index(
7527 child_at_origin[hex->face(2)->refinement_case() -
7528 1][f_fl[2]][f_ro[2]]),
7530 hex->face(3)->child_index(
7531 child_at_origin[hex->face(3)->refinement_case() -
7532 1][f_fl[3]][f_ro[3]]),
7533 hex->face(3)->child_index(
7535 child_at_origin[hex->face(3)->refinement_case() -
7536 1][f_fl[3]][f_ro[3]]),
7538 hex->face(4)->index(),
7540 hex->face(5)->index()
7586 1>(middle_vertex_index<dim, spacedim>(
7588 middle_vertex_index<dim, spacedim>(
7657 spacedim>::raw_line_iterator lines[13] = {
7658 hex->face(0)->child(0)->line(
7659 (hex->face(0)->refinement_case() ==
7663 hex->face(1)->child(0)->line(
7664 (hex->face(1)->refinement_case() ==
7668 hex->face(2)->child(0)->line(
7669 (hex->face(2)->refinement_case() ==
7673 hex->face(3)->child(0)->line(
7674 (hex->face(3)->refinement_case() ==
7682 0, f_or[4], f_fl[4], f_ro[4]))
7685 1, f_or[4], f_fl[4], f_ro[4])),
7689 3, f_or[4], f_fl[4], f_ro[4]))
7692 0, f_or[4], f_fl[4], f_ro[4])),
7696 0, f_or[4], f_fl[4], f_ro[4]))
7699 3, f_or[4], f_fl[4], f_ro[4])),
7703 3, f_or[4], f_fl[4], f_ro[4]))
7706 2, f_or[4], f_fl[4], f_ro[4])),
7711 0, f_or[5], f_fl[5], f_ro[5]))
7714 1, f_or[5], f_fl[5], f_ro[5])),
7718 3, f_or[5], f_fl[5], f_ro[5]))
7721 0, f_or[5], f_fl[5], f_ro[5])),
7725 0, f_or[5], f_fl[5], f_ro[5]))
7728 3, f_or[5], f_fl[5], f_ro[5])),
7732 3, f_or[5], f_fl[5], f_ro[5]))
7735 2, f_or[5], f_fl[5], f_ro[5])),
7740 unsigned int line_indices[13];
7741 for (
unsigned int i = 0; i < 13; ++i)
7742 line_indices[i] = lines[i]->index();
7749 bool line_orientation[13];
7753 const unsigned int middle_vertices[4] = {
7754 hex->line(0)->child(0)->vertex_index(1),
7755 hex->line(1)->child(0)->vertex_index(1),
7756 hex->line(2)->child(0)->vertex_index(1),
7757 hex->line(3)->child(0)->vertex_index(1),
7763 for (
unsigned int i = 0; i < 4; ++i)
7764 if (lines[i]->vertex_index(0) == middle_vertices[i])
7765 line_orientation[i] =
true;
7769 Assert(lines[i]->vertex_index(1) ==
7772 line_orientation[i] =
false;
7781 for (
unsigned int i = 4; i < 12; ++i)
7782 if (lines[i]->vertex_index((i + 1) % 2) ==
7783 middle_vertex_index<dim, spacedim>(
7784 hex->face(3 + i / 4)))
7785 line_orientation[i] =
true;
7790 Assert(lines[i]->vertex_index(i % 2) ==
7791 (middle_vertex_index<dim, spacedim>(
7792 hex->face(3 + i / 4))),
7794 line_orientation[i] =
false;
7799 line_orientation[12] =
true;
7832 2>(line_indices[12],
7849 new_quads[0]->set_line_orientation(
7850 0, line_orientation[2]);
7851 new_quads[0]->set_line_orientation(
7852 2, line_orientation[4]);
7853 new_quads[0]->set_line_orientation(
7854 3, line_orientation[8]);
7856 new_quads[1]->set_line_orientation(
7857 1, line_orientation[3]);
7858 new_quads[1]->set_line_orientation(
7859 2, line_orientation[5]);
7860 new_quads[1]->set_line_orientation(
7861 3, line_orientation[9]);
7863 new_quads[2]->set_line_orientation(
7864 0, line_orientation[6]);
7865 new_quads[2]->set_line_orientation(
7866 1, line_orientation[10]);
7867 new_quads[2]->set_line_orientation(
7868 2, line_orientation[0]);
7870 new_quads[3]->set_line_orientation(
7871 0, line_orientation[7]);
7872 new_quads[3]->set_line_orientation(
7873 1, line_orientation[11]);
7874 new_quads[3]->set_line_orientation(
7875 3, line_orientation[1]);
7908 const int quad_indices[20] = {
7909 new_quads[0]->index(),
7910 new_quads[1]->index(),
7911 new_quads[2]->index(),
7912 new_quads[3]->index(),
7914 hex->face(0)->child_index(
7915 child_at_origin[hex->face(0)->refinement_case() -
7916 1][f_fl[0]][f_ro[0]]),
7917 hex->face(0)->child_index(
7919 child_at_origin[hex->face(0)->refinement_case() -
7920 1][f_fl[0]][f_ro[0]]),
7922 hex->face(1)->child_index(
7923 child_at_origin[hex->face(1)->refinement_case() -
7924 1][f_fl[1]][f_ro[1]]),
7925 hex->face(1)->child_index(
7927 child_at_origin[hex->face(1)->refinement_case() -
7928 1][f_fl[1]][f_ro[1]]),
7930 hex->face(2)->child_index(
7931 child_at_origin[hex->face(2)->refinement_case() -
7932 1][f_fl[2]][f_ro[2]]),
7933 hex->face(2)->child_index(
7935 child_at_origin[hex->face(2)->refinement_case() -
7936 1][f_fl[2]][f_ro[2]]),
7938 hex->face(3)->child_index(
7939 child_at_origin[hex->face(3)->refinement_case() -
7940 1][f_fl[3]][f_ro[3]]),
7941 hex->face(3)->child_index(
7943 child_at_origin[hex->face(3)->refinement_case() -
7944 1][f_fl[3]][f_ro[3]]),
7946 hex->face(4)->isotropic_child_index(
7948 0, f_or[4], f_fl[4], f_ro[4])),
7949 hex->face(4)->isotropic_child_index(
7951 1, f_or[4], f_fl[4], f_ro[4])),
7952 hex->face(4)->isotropic_child_index(
7954 2, f_or[4], f_fl[4], f_ro[4])),
7955 hex->face(4)->isotropic_child_index(
7957 3, f_or[4], f_fl[4], f_ro[4])),
7959 hex->face(5)->isotropic_child_index(
7961 0, f_or[5], f_fl[5], f_ro[5])),
7962 hex->face(5)->isotropic_child_index(
7964 1, f_or[5], f_fl[5], f_ro[5])),
7965 hex->face(5)->isotropic_child_index(
7967 2, f_or[5], f_fl[5], f_ro[5])),
7968 hex->face(5)->isotropic_child_index(
7970 3, f_or[5], f_fl[5], f_ro[5]))};
8031 1>(middle_vertex_index<dim, spacedim>(
8033 middle_vertex_index<dim, spacedim>(
8102 spacedim>::raw_line_iterator lines[13] = {
8103 hex->face(0)->child(0)->line(
8104 (hex->face(0)->refinement_case() ==
8108 hex->face(1)->child(0)->line(
8109 (hex->face(1)->refinement_case() ==
8113 hex->face(4)->child(0)->line(
8114 (hex->face(4)->refinement_case() ==
8118 hex->face(5)->child(0)->line(
8119 (hex->face(5)->refinement_case() ==
8127 0, f_or[2], f_fl[2], f_ro[2]))
8130 3, f_or[2], f_fl[2], f_ro[2])),
8134 3, f_or[2], f_fl[2], f_ro[2]))
8137 2, f_or[2], f_fl[2], f_ro[2])),
8141 0, f_or[2], f_fl[2], f_ro[2]))
8144 1, f_or[2], f_fl[2], f_ro[2])),
8148 3, f_or[2], f_fl[2], f_ro[2]))
8151 0, f_or[2], f_fl[2], f_ro[2])),
8156 0, f_or[3], f_fl[3], f_ro[3]))
8159 3, f_or[3], f_fl[3], f_ro[3])),
8163 3, f_or[3], f_fl[3], f_ro[3]))
8166 2, f_or[3], f_fl[3], f_ro[3])),
8170 0, f_or[3], f_fl[3], f_ro[3]))
8173 1, f_or[3], f_fl[3], f_ro[3])),
8177 3, f_or[3], f_fl[3], f_ro[3]))
8180 0, f_or[3], f_fl[3], f_ro[3])),
8185 unsigned int line_indices[13];
8186 for (
unsigned int i = 0; i < 13; ++i)
8187 line_indices[i] = lines[i]->index();
8194 bool line_orientation[13];
8198 const unsigned int middle_vertices[4] = {
8199 hex->line(8)->child(0)->vertex_index(1),
8200 hex->line(9)->child(0)->vertex_index(1),
8201 hex->line(2)->child(0)->vertex_index(1),
8202 hex->line(6)->child(0)->vertex_index(1),
8207 for (
unsigned int i = 0; i < 4; ++i)
8208 if (lines[i]->vertex_index(0) == middle_vertices[i])
8209 line_orientation[i] =
true;
8213 Assert(lines[i]->vertex_index(1) ==
8216 line_orientation[i] =
false;
8225 for (
unsigned int i = 4; i < 12; ++i)
8226 if (lines[i]->vertex_index((i + 1) % 2) ==
8227 middle_vertex_index<dim, spacedim>(
8228 hex->face(1 + i / 4)))
8229 line_orientation[i] =
true;
8234 Assert(lines[i]->vertex_index(i % 2) ==
8235 (middle_vertex_index<dim, spacedim>(
8236 hex->face(1 + i / 4))),
8238 line_orientation[i] =
false;
8243 line_orientation[12] =
true;
8277 2>(line_indices[12],
8294 new_quads[0]->set_line_orientation(
8295 0, line_orientation[0]);
8296 new_quads[0]->set_line_orientation(
8297 2, line_orientation[6]);
8298 new_quads[0]->set_line_orientation(
8299 3, line_orientation[10]);
8301 new_quads[1]->set_line_orientation(
8302 1, line_orientation[1]);
8303 new_quads[1]->set_line_orientation(
8304 2, line_orientation[7]);
8305 new_quads[1]->set_line_orientation(
8306 3, line_orientation[11]);
8308 new_quads[2]->set_line_orientation(
8309 0, line_orientation[4]);
8310 new_quads[2]->set_line_orientation(
8311 1, line_orientation[8]);
8312 new_quads[2]->set_line_orientation(
8313 2, line_orientation[2]);
8315 new_quads[3]->set_line_orientation(
8316 0, line_orientation[5]);
8317 new_quads[3]->set_line_orientation(
8318 1, line_orientation[9]);
8319 new_quads[3]->set_line_orientation(
8320 3, line_orientation[3]);
8354 const int quad_indices[20] = {
8355 new_quads[0]->index(),
8356 new_quads[1]->index(),
8357 new_quads[2]->index(),
8358 new_quads[3]->index(),
8360 hex->face(0)->child_index(
8361 child_at_origin[hex->face(0)->refinement_case() -
8362 1][f_fl[0]][f_ro[0]]),
8363 hex->face(0)->child_index(
8365 child_at_origin[hex->face(0)->refinement_case() -
8366 1][f_fl[0]][f_ro[0]]),
8368 hex->face(1)->child_index(
8369 child_at_origin[hex->face(1)->refinement_case() -
8370 1][f_fl[1]][f_ro[1]]),
8371 hex->face(1)->child_index(
8373 child_at_origin[hex->face(1)->refinement_case() -
8374 1][f_fl[1]][f_ro[1]]),
8376 hex->face(2)->isotropic_child_index(
8378 0, f_or[2], f_fl[2], f_ro[2])),
8379 hex->face(2)->isotropic_child_index(
8381 1, f_or[2], f_fl[2], f_ro[2])),
8382 hex->face(2)->isotropic_child_index(
8384 2, f_or[2], f_fl[2], f_ro[2])),
8385 hex->face(2)->isotropic_child_index(
8387 3, f_or[2], f_fl[2], f_ro[2])),
8389 hex->face(3)->isotropic_child_index(
8391 0, f_or[3], f_fl[3], f_ro[3])),
8392 hex->face(3)->isotropic_child_index(
8394 1, f_or[3], f_fl[3], f_ro[3])),
8395 hex->face(3)->isotropic_child_index(
8397 2, f_or[3], f_fl[3], f_ro[3])),
8398 hex->face(3)->isotropic_child_index(
8400 3, f_or[3], f_fl[3], f_ro[3])),
8402 hex->face(4)->child_index(
8403 child_at_origin[hex->face(4)->refinement_case() -
8404 1][f_fl[4]][f_ro[4]]),
8405 hex->face(4)->child_index(
8407 child_at_origin[hex->face(4)->refinement_case() -
8408 1][f_fl[4]][f_ro[4]]),
8410 hex->face(5)->child_index(
8411 child_at_origin[hex->face(5)->refinement_case() -
8412 1][f_fl[5]][f_ro[5]]),
8413 hex->face(5)->child_index(
8415 child_at_origin[hex->face(5)->refinement_case() -
8416 1][f_fl[5]][f_ro[5]])};
8487 1>(middle_vertex_index<dim, spacedim>(
8489 middle_vertex_index<dim, spacedim>(
8559 spacedim>::raw_line_iterator lines[13] = {
8560 hex->face(2)->child(0)->line(
8561 (hex->face(2)->refinement_case() ==
8565 hex->face(3)->child(0)->line(
8566 (hex->face(3)->refinement_case() ==
8570 hex->face(4)->child(0)->line(
8571 (hex->face(4)->refinement_case() ==
8575 hex->face(5)->child(0)->line(
8576 (hex->face(5)->refinement_case() ==
8584 0, f_or[0], f_fl[0], f_ro[0]))
8587 1, f_or[0], f_fl[0], f_ro[0])),
8591 3, f_or[0], f_fl[0], f_ro[0]))
8594 0, f_or[0], f_fl[0], f_ro[0])),
8598 0, f_or[0], f_fl[0], f_ro[0]))
8601 3, f_or[0], f_fl[0], f_ro[0])),
8605 3, f_or[0], f_fl[0], f_ro[0]))
8608 2, f_or[0], f_fl[0], f_ro[0])),
8613 0, f_or[1], f_fl[1], f_ro[1]))
8616 1, f_or[1], f_fl[1], f_ro[1])),
8620 3, f_or[1], f_fl[1], f_ro[1]))
8623 0, f_or[1], f_fl[1], f_ro[1])),
8627 0, f_or[1], f_fl[1], f_ro[1]))
8630 3, f_or[1], f_fl[1], f_ro[1])),
8634 3, f_or[1], f_fl[1], f_ro[1]))
8637 2, f_or[1], f_fl[1], f_ro[1])),
8642 unsigned int line_indices[13];
8644 for (
unsigned int i = 0; i < 13; ++i)
8645 line_indices[i] = lines[i]->index();
8652 bool line_orientation[13];
8656 const unsigned int middle_vertices[4] = {
8657 hex->line(8)->child(0)->vertex_index(1),
8658 hex->line(10)->child(0)->vertex_index(1),
8659 hex->line(0)->child(0)->vertex_index(1),
8660 hex->line(4)->child(0)->vertex_index(1),
8665 for (
unsigned int i = 0; i < 4; ++i)
8666 if (lines[i]->vertex_index(0) == middle_vertices[i])
8667 line_orientation[i] =
true;
8671 Assert(lines[i]->vertex_index(1) ==
8674 line_orientation[i] =
false;
8683 for (
unsigned int i = 4; i < 12; ++i)
8684 if (lines[i]->vertex_index((i + 1) % 2) ==
8685 middle_vertex_index<dim, spacedim>(
8686 hex->face(i / 4 - 1)))
8687 line_orientation[i] =
true;
8692 Assert(lines[i]->vertex_index(i % 2) ==
8693 (middle_vertex_index<dim, spacedim>(
8694 hex->face(i / 4 - 1))),
8696 line_orientation[i] =
false;
8701 line_orientation[12] =
true;
8741 2>(line_indices[12],
8746 new_quads[0]->set_line_orientation(
8747 0, line_orientation[6]);
8748 new_quads[0]->set_line_orientation(
8749 1, line_orientation[10]);
8750 new_quads[0]->set_line_orientation(
8751 2, line_orientation[0]);
8753 new_quads[1]->set_line_orientation(
8754 0, line_orientation[7]);
8755 new_quads[1]->set_line_orientation(
8756 1, line_orientation[11]);
8757 new_quads[1]->set_line_orientation(
8758 3, line_orientation[1]);
8760 new_quads[2]->set_line_orientation(
8761 0, line_orientation[2]);
8762 new_quads[2]->set_line_orientation(
8763 2, line_orientation[4]);
8764 new_quads[2]->set_line_orientation(
8765 3, line_orientation[8]);
8767 new_quads[3]->set_line_orientation(
8768 1, line_orientation[3]);
8769 new_quads[3]->set_line_orientation(
8770 2, line_orientation[5]);
8771 new_quads[3]->set_line_orientation(
8772 3, line_orientation[9]);
8806 const int quad_indices[20] = {
8807 new_quads[0]->index(),
8808 new_quads[1]->index(),
8809 new_quads[2]->index(),
8810 new_quads[3]->index(),
8812 hex->face(0)->isotropic_child_index(
8814 0, f_or[0], f_fl[0], f_ro[0])),
8815 hex->face(0)->isotropic_child_index(
8817 1, f_or[0], f_fl[0], f_ro[0])),
8818 hex->face(0)->isotropic_child_index(
8820 2, f_or[0], f_fl[0], f_ro[0])),
8821 hex->face(0)->isotropic_child_index(
8823 3, f_or[0], f_fl[0], f_ro[0])),
8825 hex->face(1)->isotropic_child_index(
8827 0, f_or[1], f_fl[1], f_ro[1])),
8828 hex->face(1)->isotropic_child_index(
8830 1, f_or[1], f_fl[1], f_ro[1])),
8831 hex->face(1)->isotropic_child_index(
8833 2, f_or[1], f_fl[1], f_ro[1])),
8834 hex->face(1)->isotropic_child_index(
8836 3, f_or[1], f_fl[1], f_ro[1])),
8838 hex->face(2)->child_index(
8839 child_at_origin[hex->face(2)->refinement_case() -
8840 1][f_fl[2]][f_ro[2]]),
8841 hex->face(2)->child_index(
8843 child_at_origin[hex->face(2)->refinement_case() -
8844 1][f_fl[2]][f_ro[2]]),
8846 hex->face(3)->child_index(
8847 child_at_origin[hex->face(3)->refinement_case() -
8848 1][f_fl[3]][f_ro[3]]),
8849 hex->face(3)->child_index(
8851 child_at_origin[hex->face(3)->refinement_case() -
8852 1][f_fl[3]][f_ro[3]]),
8854 hex->face(4)->child_index(
8855 child_at_origin[hex->face(4)->refinement_case() -
8856 1][f_fl[4]][f_ro[4]]),
8857 hex->face(4)->child_index(
8859 child_at_origin[hex->face(4)->refinement_case() -
8860 1][f_fl[4]][f_ro[4]]),
8862 hex->face(5)->child_index(
8863 child_at_origin[hex->face(5)->refinement_case() -
8864 1][f_fl[5]][f_ro[5]]),
8865 hex->face(5)->child_index(
8867 child_at_origin[hex->face(5)->refinement_case() -
8868 1][f_fl[5]][f_ro[5]])};
8933 ++next_unused_vertex;
8935 next_unused_vertex < triangulation.
vertices.size(),
8937 "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."));
8947 triangulation.
vertices[next_unused_vertex] =
8948 hex->center(
true,
true);
8969 const unsigned int vertex_indices[7] = {
8970 middle_vertex_index<dim, spacedim>(hex->face(0)),
8971 middle_vertex_index<dim, spacedim>(hex->face(1)),
8972 middle_vertex_index<dim, spacedim>(hex->face(2)),
8973 middle_vertex_index<dim, spacedim>(hex->face(3)),
8974 middle_vertex_index<dim, spacedim>(hex->face(4)),
8975 middle_vertex_index<dim, spacedim>(hex->face(5)),
8976 next_unused_vertex};
8980 1>(vertex_indices[2], vertex_indices[6]));
8983 1>(vertex_indices[6], vertex_indices[3]));
8986 1>(vertex_indices[0], vertex_indices[6]));
8989 1>(vertex_indices[6], vertex_indices[1]));
8992 1>(vertex_indices[4], vertex_indices[6]));
8995 1>(vertex_indices[6], vertex_indices[5]));
9064 spacedim>::raw_line_iterator lines[30] = {
9068 0, f_or[0], f_fl[0], f_ro[0]))
9071 1, f_or[0], f_fl[0], f_ro[0])),
9075 3, f_or[0], f_fl[0], f_ro[0]))
9078 0, f_or[0], f_fl[0], f_ro[0])),
9082 0, f_or[0], f_fl[0], f_ro[0]))
9085 3, f_or[0], f_fl[0], f_ro[0])),
9089 3, f_or[0], f_fl[0], f_ro[0]))
9092 2, f_or[0], f_fl[0], f_ro[0])),
9097 0, f_or[1], f_fl[1], f_ro[1]))
9100 1, f_or[1], f_fl[1], f_ro[1])),
9104 3, f_or[1], f_fl[1], f_ro[1]))
9107 0, f_or[1], f_fl[1], f_ro[1])),
9111 0, f_or[1], f_fl[1], f_ro[1]))
9114 3, f_or[1], f_fl[1], f_ro[1])),
9118 3, f_or[1], f_fl[1], f_ro[1]))
9121 2, f_or[1], f_fl[1], f_ro[1])),
9126 0, f_or[2], f_fl[2], f_ro[2]))
9129 1, f_or[2], f_fl[2], f_ro[2])),
9133 3, f_or[2], f_fl[2], f_ro[2]))
9136 0, f_or[2], f_fl[2], f_ro[2])),
9140 0, f_or[2], f_fl[2], f_ro[2]))
9143 3, f_or[2], f_fl[2], f_ro[2])),
9147 3, f_or[2], f_fl[2], f_ro[2]))
9150 2, f_or[2], f_fl[2], f_ro[2])),
9155 0, f_or[3], f_fl[3], f_ro[3]))
9158 1, f_or[3], f_fl[3], f_ro[3])),
9162 3, f_or[3], f_fl[3], f_ro[3]))
9165 0, f_or[3], f_fl[3], f_ro[3])),
9169 0, f_or[3], f_fl[3], f_ro[3]))
9172 3, f_or[3], f_fl[3], f_ro[3])),
9176 3, f_or[3], f_fl[3], f_ro[3]))
9179 2, f_or[3], f_fl[3], f_ro[3])),
9184 0, f_or[4], f_fl[4], f_ro[4]))
9187 1, f_or[4], f_fl[4], f_ro[4])),
9191 3, f_or[4], f_fl[4], f_ro[4]))
9194 0, f_or[4], f_fl[4], f_ro[4])),
9198 0, f_or[4], f_fl[4], f_ro[4]))
9201 3, f_or[4], f_fl[4], f_ro[4])),
9205 3, f_or[4], f_fl[4], f_ro[4]))
9208 2, f_or[4], f_fl[4], f_ro[4])),
9213 0, f_or[5], f_fl[5], f_ro[5]))
9216 1, f_or[5], f_fl[5], f_ro[5])),
9220 3, f_or[5], f_fl[5], f_ro[5]))
9223 0, f_or[5], f_fl[5], f_ro[5])),
9227 0, f_or[5], f_fl[5], f_ro[5]))
9230 3, f_or[5], f_fl[5], f_ro[5])),
9234 3, f_or[5], f_fl[5], f_ro[5]))
9237 2, f_or[5], f_fl[5], f_ro[5])),
9247 unsigned int line_indices[30];
9248 for (
unsigned int i = 0; i < 30; ++i)
9249 line_indices[i] = lines[i]->index();
9256 bool line_orientation[30];
9264 for (
unsigned int i = 0; i < 24; ++i)
9265 if (lines[i]->vertex_index((i + 1) % 2) ==
9266 vertex_indices[i / 4])
9267 line_orientation[i] =
true;
9272 Assert(lines[i]->vertex_index(i % 2) ==
9273 vertex_indices[i / 4],
9275 line_orientation[i] =
false;
9280 for (
unsigned int i = 24; i < 30; ++i)
9281 line_orientation[i] =
true;
9315 2>(line_indices[10],
9321 2>(line_indices[28],
9327 2>(line_indices[11],
9333 2>(line_indices[29],
9339 2>(line_indices[18],
9345 2>(line_indices[26],
9351 2>(line_indices[19],
9357 2>(line_indices[27],
9369 2>(line_indices[24],
9381 2>(line_indices[25],
9389 new_quads[0]->set_line_orientation(
9390 0, line_orientation[10]);
9391 new_quads[0]->set_line_orientation(
9392 2, line_orientation[16]);
9394 new_quads[1]->set_line_orientation(
9395 1, line_orientation[14]);
9396 new_quads[1]->set_line_orientation(
9397 2, line_orientation[17]);
9399 new_quads[2]->set_line_orientation(
9400 0, line_orientation[11]);
9401 new_quads[2]->set_line_orientation(
9402 3, line_orientation[20]);
9404 new_quads[3]->set_line_orientation(
9405 1, line_orientation[15]);
9406 new_quads[3]->set_line_orientation(
9407 3, line_orientation[21]);
9409 new_quads[4]->set_line_orientation(
9410 0, line_orientation[18]);
9411 new_quads[4]->set_line_orientation(
9412 2, line_orientation[0]);
9414 new_quads[5]->set_line_orientation(
9415 1, line_orientation[22]);
9416 new_quads[5]->set_line_orientation(
9417 2, line_orientation[1]);
9419 new_quads[6]->set_line_orientation(
9420 0, line_orientation[19]);
9421 new_quads[6]->set_line_orientation(
9422 3, line_orientation[4]);
9424 new_quads[7]->set_line_orientation(
9425 1, line_orientation[23]);
9426 new_quads[7]->set_line_orientation(
9427 3, line_orientation[5]);
9429 new_quads[8]->set_line_orientation(
9430 0, line_orientation[2]);
9431 new_quads[8]->set_line_orientation(
9432 2, line_orientation[8]);
9434 new_quads[9]->set_line_orientation(
9435 1, line_orientation[6]);
9436 new_quads[9]->set_line_orientation(
9437 2, line_orientation[9]);
9439 new_quads[10]->set_line_orientation(
9440 0, line_orientation[3]);
9441 new_quads[10]->set_line_orientation(
9442 3, line_orientation[12]);
9444 new_quads[11]->set_line_orientation(
9445 1, line_orientation[7]);
9446 new_quads[11]->set_line_orientation(
9447 3, line_orientation[13]);
9488 const int quad_indices[36] = {
9489 new_quads[0]->index(),
9490 new_quads[1]->index(),
9491 new_quads[2]->index(),
9492 new_quads[3]->index(),
9493 new_quads[4]->index(),
9494 new_quads[5]->index(),
9495 new_quads[6]->index(),
9496 new_quads[7]->index(),
9497 new_quads[8]->index(),
9498 new_quads[9]->index(),
9499 new_quads[10]->index(),
9500 new_quads[11]->index(),
9502 hex->face(0)->isotropic_child_index(
9504 0, f_or[0], f_fl[0], f_ro[0])),
9505 hex->face(0)->isotropic_child_index(
9507 1, f_or[0], f_fl[0], f_ro[0])),
9508 hex->face(0)->isotropic_child_index(
9510 2, f_or[0], f_fl[0], f_ro[0])),
9511 hex->face(0)->isotropic_child_index(
9513 3, f_or[0], f_fl[0], f_ro[0])),
9515 hex->face(1)->isotropic_child_index(
9517 0, f_or[1], f_fl[1], f_ro[1])),
9518 hex->face(1)->isotropic_child_index(
9520 1, f_or[1], f_fl[1], f_ro[1])),
9521 hex->face(1)->isotropic_child_index(
9523 2, f_or[1], f_fl[1], f_ro[1])),
9524 hex->face(1)->isotropic_child_index(
9526 3, f_or[1], f_fl[1], f_ro[1])),
9528 hex->face(2)->isotropic_child_index(
9530 0, f_or[2], f_fl[2], f_ro[2])),
9531 hex->face(2)->isotropic_child_index(
9533 1, f_or[2], f_fl[2], f_ro[2])),
9534 hex->face(2)->isotropic_child_index(
9536 2, f_or[2], f_fl[2], f_ro[2])),
9537 hex->face(2)->isotropic_child_index(
9539 3, f_or[2], f_fl[2], f_ro[2])),
9541 hex->face(3)->isotropic_child_index(
9543 0, f_or[3], f_fl[3], f_ro[3])),
9544 hex->face(3)->isotropic_child_index(
9546 1, f_or[3], f_fl[3], f_ro[3])),
9547 hex->face(3)->isotropic_child_index(
9549 2, f_or[3], f_fl[3], f_ro[3])),
9550 hex->face(3)->isotropic_child_index(
9552 3, f_or[3], f_fl[3], f_ro[3])),
9554 hex->face(4)->isotropic_child_index(
9556 0, f_or[4], f_fl[4], f_ro[4])),
9557 hex->face(4)->isotropic_child_index(
9559 1, f_or[4], f_fl[4], f_ro[4])),
9560 hex->face(4)->isotropic_child_index(
9562 2, f_or[4], f_fl[4], f_ro[4])),
9563 hex->face(4)->isotropic_child_index(
9565 3, f_or[4], f_fl[4], f_ro[4])),
9567 hex->face(5)->isotropic_child_index(
9569 0, f_or[5], f_fl[5], f_ro[5])),
9570 hex->face(5)->isotropic_child_index(
9572 1, f_or[5], f_fl[5], f_ro[5])),
9573 hex->face(5)->isotropic_child_index(
9575 2, f_or[5], f_fl[5], f_ro[5])),
9576 hex->face(5)->isotropic_child_index(
9578 3, f_or[5], f_fl[5], f_ro[5]))};
9583 3>(quad_indices[12],
9599 3>(quad_indices[13],
9617 3>(quad_indices[14],
9633 3>(quad_indices[15],
9679 for (
unsigned int f = 0;
9680 f < GeometryInfo<dim>::faces_per_cell;
9682 for (
unsigned int s = 0;
9689 const unsigned int current_child =
9698 ref_case, f, f_or[f], f_fl[f], f_ro[f]));
9699 new_hexes[current_child]->set_face_orientation(f,
9701 new_hexes[current_child]->set_face_flip(f, f_fl[f]);
9702 new_hexes[current_child]->set_face_rotation(f, f_ro[f]);
9707 if ((check_for_distorted_cells ==
true) &&
9708 has_distorted_children(
9710 std::integral_constant<int, dim>(),
9711 std::integral_constant<int, spacedim>()))
9719 triangulation.
signals.post_refinement_on_cell(hex);
9728 triangulation.
faces->quads.clear_user_data();
9731 return cells_with_distorted_children;
9747 template <
int spacedim>
9752 template <
int dim,
int spacedim>
9754 prevent_distorted_boundary_cells(
9763 triangulation.
begin();
9764 cell != triangulation.
end();
9766 if (cell->at_boundary() && cell->refine_flag_set() &&
9767 cell->refine_flag_set() !=
9774 for (
unsigned int face_no = 0;
9775 face_no < GeometryInfo<dim>::faces_per_cell;
9777 if (cell->face(face_no)->at_boundary())
9788 spacedim>::face_iterator
9789 face = cell->face(face_no);
9797 .transform_real_to_unit_cell(cell, new_bound);
9808 std::fabs(new_unit[face_no / 2] - face_no % 2);
9813 const double allowed = 0.25;
9816 cell->flag_for_face_refinement(face_no);
9835 template <
int dim,
int spacedim>
9840 ExcMessage(
"Wrong function called -- there should " 9841 "be a specialization."));
9845 template <
int spacedim>
9846 static void prepare_refinement_dim_dependent(
9849 const unsigned int dim = 3;
9861 bool mesh_changed =
false;
9865 mesh_changed =
false;
9876 triangulation.
begin();
9877 cell != triangulation.
end();
9879 if (cell->refine_flag_set())
9881 for (
unsigned int line = 0;
9882 line < GeometryInfo<dim>::lines_per_cell;
9885 cell->refine_flag_set(), line) ==
9889 cell->line(line)->set_user_flag();
9891 else if (cell->has_children() &&
9892 !cell->child(0)->coarsen_flag_set())
9894 for (
unsigned int line = 0;
9895 line < GeometryInfo<dim>::lines_per_cell;
9898 cell->refinement_case(), line) ==
9902 cell->line(line)->set_user_flag();
9904 else if (cell->has_children() &&
9905 cell->child(0)->coarsen_flag_set())
9906 cell->set_user_flag();
9915 cell != triangulation.
end();
9917 for (
unsigned int line = 0;
9918 line < GeometryInfo<dim>::lines_per_cell;
9921 if (cell->line(line)->has_children())
9930 bool offending_line_found =
false;
9932 for (
unsigned int c = 0; c < 2; ++c)
9934 Assert(cell->line(line)->child(c)->has_children() ==
9938 if (cell->line(line)->child(c)->user_flag_set() &&
9940 cell->refine_flag_set(), line) ==
9944 cell->clear_coarsen_flag();
9951 allow_anisotropic_smoothing)
9952 cell->flag_for_line_refinement(line);
9954 cell->set_refine_flag();
9956 for (
unsigned int l = 0;
9957 l < GeometryInfo<dim>::lines_per_cell;
9960 cell->refine_flag_set(), line) ==
9963 cell->line(l)->set_user_flag();
9966 offending_line_found =
true;
9972 for (
unsigned int l = 0;
9973 l < GeometryInfo<dim>::lines_per_cell;
9975 if (!cell->line(l)->has_children() &&
9977 cell->refine_flag_set(), l) !=
9979 cell->line(l)->set_user_flag();
9985 if (offending_line_found)
9987 mesh_changed =
true;
10005 triangulation.
last();
10006 cell != triangulation.
end();
10009 if (cell->user_flag_set())
10010 for (
unsigned int line = 0;
10011 line < GeometryInfo<dim>::lines_per_cell;
10013 if (cell->line(line)->has_children() &&
10014 (cell->line(line)->child(0)->user_flag_set() ||
10015 cell->line(line)->child(1)->user_flag_set()))
10017 for (
unsigned int c = 0; c < cell->n_children(); ++c)
10018 cell->child(c)->clear_coarsen_flag();
10019 cell->clear_user_flag();
10020 for (
unsigned int l = 0;
10021 l < GeometryInfo<dim>::lines_per_cell;
10024 cell->refinement_case(), l) ==
10028 cell->line(l)->set_user_flag();
10029 mesh_changed =
true;
10034 while (mesh_changed ==
true);
10045 template <
int dim,
int spacedim>
10056 for (
unsigned int n = 0; n < GeometryInfo<dim>::faces_per_cell; ++n)
10065 const unsigned int n_subfaces =
10068 if (n_subfaces == 0 || cell->at_boundary(n))
10070 for (
unsigned int c = 0; c < n_subfaces; ++c)
10073 child = cell->child(
10077 child_neighbor = child->neighbor(n);
10078 if (!child->neighbor_is_coarser(n))
10091 if ((child_neighbor->has_children() &&
10092 !child_neighbor->user_flag_set()) ||
10101 child_neighbor->refine_flag_set())
10111 template <
int dim,
int spacedim>
10113 get_default_flat_manifold()
10116 return flat_manifold;
10123 template <
int dim,
int spacedim>
10128 template <
int dim,
int spacedim>
10131 const bool check_for_distorted_cells)
10132 : smooth_grid(smooth_grid)
10133 , anisotropic_refinement(false)
10134 , check_for_distorted_cells(check_for_distorted_cells)
10139 std_cxx14::make_unique<std::map<unsigned int, types::boundary_id>>();
10141 std_cxx14::make_unique<std::map<unsigned int, types::manifold_id>>();
10153 template <
int dim,
int spacedim>
10175 template <
int dim,
int spacedim>
10185 levels = std::move(tria.levels);
10186 faces = std::move(tria.faces);
10187 vertices = std::move(tria.vertices);
10189 manifold = std::move(tria.manifold);
10202 template <
int dim,
int spacedim>
10230 template <
int dim,
int spacedim>
10245 template <
int dim,
int spacedim>
10257 template <
int dim,
int spacedim>
10266 template <
int dim,
int spacedim>
10280 template <
int dim,
int spacedim>
10288 template <
int dim,
int spacedim>
10300 template <
int dim,
int spacedim>
10308 template <
int dim,
int spacedim>
10316 "Error: set_all_manifold_ids() can not be called on an empty Triangulation."));
10320 endc = this->
end();
10322 for (; cell != endc; ++cell)
10323 cell->set_all_manifold_ids(m_number);
10327 template <
int dim,
int spacedim>
10335 "Error: set_all_manifold_ids_on_boundary() can not be called on an empty Triangulation."));
10339 endc = this->
end();
10341 for (; cell != endc; ++cell)
10342 for (
unsigned int f = 0; f < GeometryInfo<dim>::faces_per_cell; ++f)
10343 if (cell->face(f)->at_boundary())
10344 cell->face(f)->set_all_manifold_ids(m_number);
10348 template <
int dim,
int spacedim>
10357 "Error: set_all_manifold_ids_on_boundary() can not be called on an empty Triangulation."));
10359 bool boundary_found =
false;
10362 endc = this->
end();
10364 for (; cell != endc; ++cell)
10367 for (
unsigned int f = 0; f < GeometryInfo<dim>::faces_per_cell; ++f)
10368 if (cell->face(f)->at_boundary() &&
10369 cell->face(f)->boundary_id() == b_id)
10371 boundary_found =
true;
10372 cell->face(f)->set_manifold_id(m_number);
10377 for (
unsigned int e = 0; e < GeometryInfo<dim>::lines_per_cell; ++e)
10378 if (cell->line(e)->at_boundary() &&
10379 cell->line(e)->boundary_id() == b_id)
10381 boundary_found =
true;
10382 cell->line(e)->set_manifold_id(m_number);
10386 (void)boundary_found;
10392 template <
int dim,
int spacedim>
10399 const auto it =
manifold.find(m_number);
10404 return *(it->second);
10409 return internal::TriangulationImplementation::
10410 get_default_flat_manifold<dim, spacedim>();
10415 template <
int dim,
int spacedim>
10416 std::vector<types::boundary_id>
10423 std::vector<types::boundary_id> boundary_ids;
10424 for (std::map<unsigned int, types::boundary_id>::const_iterator p =
10428 boundary_ids.push_back(p->second);
10430 return boundary_ids;
10434 std::set<types::boundary_id> b_ids;
10436 for (; cell !=
end(); ++cell)
10437 for (
unsigned int face = 0; face < GeometryInfo<dim>::faces_per_cell;
10439 if (cell->at_boundary(face))
10440 b_ids.insert(cell->face(face)->boundary_id());
10441 std::vector<types::boundary_id> boundary_ids(b_ids.begin(), b_ids.end());
10442 return boundary_ids;
10448 template <
int dim,
int spacedim>
10449 std::vector<types::manifold_id>
10452 std::set<types::manifold_id> m_ids;
10454 for (; cell !=
end(); ++cell)
10456 m_ids.insert(cell->manifold_id());
10458 for (
unsigned int face = 0; face < GeometryInfo<dim>::faces_per_cell;
10460 if (cell->at_boundary(face))
10461 m_ids.insert(cell->face(face)->manifold_id());
10463 std::vector<types::manifold_id> manifold_indicators(m_ids.begin(),
10465 return manifold_indicators;
10471 template <
int dim,
int spacedim>
10479 (dim == 1 || other_tria.
faces !=
nullptr),
10481 "When calling Triangulation::copy_triangulation(), " 10482 "the target triangulation must be empty but the source " 10483 "triangulation (the argument to this function) must contain " 10484 "something. Here, it seems like the source does not " 10485 "contain anything at all."));
10495 faces = std_cxx14::make_unique<
10498 auto bdry_iterator = other_tria.
manifold.begin();
10499 for (; bdry_iterator != other_tria.
manifold.end(); ++bdry_iterator)
10500 manifold[bdry_iterator->first] = bdry_iterator->second->clone();
10504 for (
unsigned int level = 0; level < other_tria.
levels.size(); ++level)
10505 levels.push_back(std_cxx14::make_unique<
10507 *other_tria.
levels[level]));
10514 std_cxx14::make_unique<std::map<unsigned int, types::boundary_id>>(
10518 std_cxx14::make_unique<std::map<unsigned int, types::manifold_id>>(
10534 template <
int dim,
int spacedim>
10541 std::vector<CellData<dim>> reordered_cells(cells);
10545 reorder_compatibility(reordered_cells, reordered_subcelldata);
10554 template <
int dim,
int spacedim>
10629 if (dim < spacedim)
10637 bool values[][2] = {{
false,
true}, {
true,
false}};
10638 for (
unsigned int i = 0; i < GeometryInfo<dim>::faces_per_cell;
10640 for (
unsigned int j = 0; j < GeometryInfo<dim>::faces_per_cell;
10642 correct(i, j) = (values[i][j]);
10647 bool values[][4] = {{
false,
true,
true,
false},
10648 {
true,
false,
false,
true},
10649 {
true,
false,
false,
true},
10650 {
false,
true,
true,
false}};
10651 for (
unsigned int i = 0; i < GeometryInfo<dim>::faces_per_cell;
10653 for (
unsigned int j = 0; j < GeometryInfo<dim>::faces_per_cell;
10655 correct(i, j) = (values[i][j]);
10663 std::list<active_cell_iterator> this_round, next_round;
10670 while (this_round.size() > 0)
10672 for (
typename std::list<active_cell_iterator>::iterator cell =
10673 this_round.begin();
10674 cell != this_round.end();
10677 for (
unsigned int i = 0; i < GeometryInfo<dim>::faces_per_cell;
10680 if (!((*cell)->face(i)->at_boundary()))
10682 neighbor = (*cell)->neighbor(i);
10684 unsigned int cf = (*cell)->face_index(i);
10685 unsigned int j = 0;
10686 while (neighbor->face_index(j) != cf)
10694 if (neighbor->user_flag_set())
10698 Assert(!(correct(i, j) ^
10699 (neighbor->direction_flag() ==
10700 (*cell)->direction_flag())),
10705 next_round.push_back(neighbor);
10706 neighbor->set_user_flag();
10707 if ((correct(i, j) ^ (neighbor->direction_flag() ==
10708 (*cell)->direction_flag())))
10709 neighbor->set_direction_flag(
10710 !neighbor->direction_flag());
10720 if (next_round.size() == 0)
10723 if (cell->user_flag_set() ==
false)
10725 next_round.push_back(cell);
10726 cell->set_direction_flag(
true);
10727 cell->set_user_flag();
10731 this_round = next_round;
10732 next_round.clear();
10742 template <
int dim,
int spacedim>
10747 ExcMessage(
"Only works for dim == spacedim-1"));
10749 cell->set_direction_flag(!cell->direction_flag());
10754 template <
int dim,
int spacedim>
10759 ExcMessage(
"Error: An empty Triangulation can not be refined."));
10762 for (; cell != endc; ++cell)
10764 cell->clear_coarsen_flag();
10765 cell->set_refine_flag();
10771 template <
int dim,
int spacedim>
10775 for (
unsigned int i = 0; i < times; ++i)
10788 template <
int dim,
int spacedim>
10793 std::vector<bool>::iterator i = v.begin();
10795 for (; cell != endc; ++cell)
10796 for (
unsigned int j = 0; j < dim; ++j, ++i)
10797 if (cell->refine_flag_set() & (1 << j))
10805 template <
int dim,
int spacedim>
10809 std::vector<bool> v;
10813 mn_tria_refine_flags_end,
10819 template <
int dim,
int spacedim>
10823 std::vector<bool> v;
10824 read_bool_vector(mn_tria_refine_flags_begin, v, mn_tria_refine_flags_end, in);
10830 template <
int dim,
int spacedim>
10837 std::vector<bool>::const_iterator i = v.begin();
10838 for (; cell != endc; ++cell)
10840 unsigned int ref_case = 0;
10842 for (
unsigned int j = 0; j < dim; ++j, ++i)
10844 ref_case += 1 << j;
10850 cell->clear_refine_flag();
10858 template <
int dim,
int spacedim>
10863 std::vector<bool>::iterator i = v.begin();
10865 for (; cell != endc; ++cell, ++i)
10866 *i = cell->coarsen_flag_set();
10873 template <
int dim,
int spacedim>
10877 std::vector<bool> v;
10881 mn_tria_coarsen_flags_end,
10887 template <
int dim,
int spacedim>
10891 std::vector<bool> v;
10894 mn_tria_coarsen_flags_end,
10901 template <
int dim,
int spacedim>
10908 std::vector<bool>::const_iterator i = v.begin();
10909 for (; cell != endc; ++cell, ++i)
10911 cell->set_coarsen_flag();
10913 cell->clear_coarsen_flag();
10919 template <
int dim,
int spacedim>
10941 for (
unsigned int level = 0; level <
levels.size(); ++level)
10942 levels[level]->cells.clear_user_data();
10956 faces->
lines.clear_user_data();
10963 faces->
lines.clear_user_data();
10969 template <
int dim,
int spacedim>
10988 for (
unsigned int level = 0; level <
levels.size(); ++level)
10989 levels[level]->cells.clear_user_flags();
10999 faces->lines.clear_user_flags();
11004 template <
int dim,
int spacedim>
11029 for (
unsigned int level = 0; level <
levels.size(); ++level)
11030 levels[level]->cells.clear_user_flags();
11040 faces->quads.clear_user_flags();
11045 template <
int dim,
int spacedim>
11079 for (
unsigned int level = 0; level <
levels.size(); ++level)
11080 levels[level]->cells.clear_user_flags();
11085 template <
int dim,
int spacedim>
11094 template <
int dim,
int spacedim>
11105 template <
int dim,
int spacedim>
11123 template <
int dim,
int spacedim>
11131 std::vector<bool> tmp;
11134 v.insert(v.end(), tmp.begin(), tmp.end());
11139 v.insert(v.end(), tmp.begin(), tmp.end());
11145 v.insert(v.end(), tmp.begin(), tmp.end());
11154 template <
int dim,
int spacedim>
11172 template <
int dim,
int spacedim>
11177 std::vector<bool> tmp;
11181 tmp.insert(tmp.end(), v.begin(), v.begin() +
n_lines());
11188 tmp.insert(tmp.end(),
11197 tmp.insert(tmp.end(),
11209 template <
int dim,
int spacedim>
11214 std::vector<bool>::iterator i = v.begin();
11216 for (; line != endl; ++line, ++i)
11217 *i = line->user_flag_set();
11224 template <
int dim,
int spacedim>
11228 std::vector<bool> v;
11232 mn_tria_line_user_flags_end,
11238 template <
int dim,
int spacedim>
11242 std::vector<bool> v;
11245 mn_tria_line_user_flags_end,
11252 template <
int dim,
int spacedim>
11259 std::vector<bool>::const_iterator i = v.begin();
11260 for (; line != endl; ++line, ++i)
11262 line->set_user_flag();
11264 line->clear_user_flag();
11272 template <
typename Iterator>
11274 get_user_flag(
const Iterator &i)
11276 return i->user_flag_set();
11281 template <
int structdim,
int dim,
int spacedim>
11291 template <
typename Iterator>
11293 set_user_flag(
const Iterator &i)
11295 i->set_user_flag();
11300 template <
int structdim,
int dim,
int spacedim>
11309 template <
typename Iterator>
11311 clear_user_flag(
const Iterator &i)
11313 i->clear_user_flag();
11318 template <
int structdim,
int dim,
int spacedim>
11328 template <
int dim,
int spacedim>
11336 std::vector<bool>::iterator i = v.begin();
11338 for (; quad != endq; ++quad, ++i)
11339 *i = get_user_flag(quad);
11347 template <
int dim,
int spacedim>
11351 std::vector<bool> v;
11355 mn_tria_quad_user_flags_end,
11361 template <
int dim,
int spacedim>
11365 std::vector<bool> v;
11368 mn_tria_quad_user_flags_end,
11375 template <
int dim,
int spacedim>
11384 std::vector<bool>::const_iterator i = v.begin();
11385 for (; quad != endq; ++quad, ++i)
11387 set_user_flag(quad);
11389 clear_user_flag(quad);
11397 template <
int dim,
int spacedim>
11401 v.resize(
n_hexs(),
false);
11405 std::vector<bool>::iterator i = v.begin();
11407 for (; hex != endh; ++hex, ++i)
11408 *i = get_user_flag(hex);
11416 template <
int dim,
int spacedim>
11420 std::vector<bool> v;
11424 mn_tria_hex_user_flags_end,
11430 template <
int dim,
int spacedim>
11434 std::vector<bool> v;
11437 mn_tria_hex_user_flags_end,
11444 template <
int dim,
int spacedim>
11453 std::vector<bool>::const_iterator i = v.begin();
11454 for (; hex != endh; ++hex, ++i)
11456 set_user_flag(hex);
11458 clear_user_flag(hex);
11466 template <
int dim,
int spacedim>
11469 std::vector<unsigned int> &v)
const 11475 std::vector<unsigned int> tmp;
11478 v.insert(v.end(), tmp.begin(), tmp.end());
11483 v.insert(v.end(), tmp.begin(), tmp.end());
11489 v.insert(v.end(), tmp.begin(), tmp.end());
11498 template <
int dim,
int spacedim>
11501 const std::vector<unsigned int> &v)
11504 std::vector<unsigned int> tmp;
11508 tmp.insert(tmp.end(), v.begin(), v.begin() +
n_lines());
11515 tmp.insert(tmp.end(),
11524 tmp.insert(tmp.end(),
11538 template <
typename Iterator>
11540 get_user_index(
const Iterator &i)
11542 return i->user_index();
11547 template <
int structdim,
int dim,
int spacedim>
11558 template <
typename Iterator>
11560 set_user_index(
const Iterator &i,
const unsigned int x)
11562 i->set_user_index(x);
11567 template <
int structdim,
int dim,
int spacedim>
11571 const unsigned int)
11578 template <
int dim,
int spacedim>
11581 std::vector<unsigned int> &v)
const 11584 std::vector<unsigned int>::iterator i = v.begin();
11586 for (; line != endl; ++line, ++i)
11587 *i = line->user_index();
11592 template <
int dim,
int spacedim>
11595 const std::vector<unsigned int> &v)
11600 std::vector<unsigned int>::const_iterator i = v.begin();
11601 for (; line != endl; ++line, ++i)
11602 line->set_user_index(*i);
11606 template <
int dim,
int spacedim>
11609 std::vector<unsigned int> &v)
const 11615 std::vector<unsigned int>::iterator i = v.begin();
11617 for (; quad != endq; ++quad, ++i)
11618 *i = get_user_index(quad);
11624 template <
int dim,
int spacedim>
11627 const std::vector<unsigned int> &v)
11634 std::vector<unsigned int>::const_iterator i = v.begin();
11635 for (; quad != endq; ++quad, ++i)
11636 set_user_index(quad, *i);
11641 template <
int dim,
int spacedim>
11644 std::vector<unsigned int> &v)
const 11650 std::vector<unsigned int>::iterator i = v.begin();
11652 for (; hex != endh; ++hex, ++i)
11653 *i = get_user_index(hex);
11659 template <
int dim,
int spacedim>
11662 const std::vector<unsigned int> &v)
11669 std::vector<unsigned int>::const_iterator i = v.begin();
11670 for (; hex != endh; ++hex, ++i)
11671 set_user_index(hex, *i);
11682 template <
typename Iterator>
11684 get_user_pointer(
const Iterator &i)
11686 return i->user_pointer();
11691 template <
int structdim,
int dim,
int spacedim>
11702 template <
typename Iterator>
11704 set_user_pointer(
const Iterator &i,
void *x)
11706 i->set_user_pointer(x);
11711 template <
int structdim,
int dim,
int spacedim>
11722 template <
int dim,
int spacedim>
11730 std::vector<void *> tmp;
11733 v.insert(v.end(), tmp.begin(), tmp.end());
11738 v.insert(v.end(), tmp.begin(), tmp.end());
11744 v.insert(v.end(), tmp.begin(), tmp.end());
11753 template <
int dim,
int spacedim>
11758 std::vector<void *> tmp;
11762 tmp.insert(tmp.end(), v.begin(), v.begin() +
n_lines());
11769 tmp.insert(tmp.end(),
11778 tmp.insert(tmp.end(),
11790 template <
int dim,
int spacedim>
11793 std::vector<void *> &v)
const 11795 v.resize(
n_lines(),
nullptr);
11796 std::vector<void *>::iterator i = v.begin();
11798 for (; line != endl; ++line, ++i)
11799 *i = line->user_pointer();
11804 template <
int dim,
int spacedim>
11807 const std::vector<void *> &v)
11812 std::vector<void *>::const_iterator i = v.begin();
11813 for (; line != endl; ++line, ++i)
11814 line->set_user_pointer(*i);
11819 template <
int dim,
int spacedim>
11822 std::vector<void *> &v)
const 11824 v.resize(
n_quads(),
nullptr);
11828 std::vector<void *>::iterator i = v.begin();
11830 for (; quad != endq; ++quad, ++i)
11831 *i = get_user_pointer(quad);
11837 template <
int dim,
int spacedim>
11840 const std::vector<void *> &v)
11847 std::vector<void *>::const_iterator i = v.begin();
11848 for (; quad != endq; ++quad, ++i)
11849 set_user_pointer(quad, *i);
11854 template <
int dim,
int spacedim>
11857 std::vector<void *> &v)
const 11859 v.resize(
n_hexs(),
nullptr);
11863 std::vector<void *>::iterator i = v.begin();
11865 for (; hex != endh; ++hex, ++i)
11866 *i = get_user_pointer(hex);
11872 template <
int dim,
int spacedim>
11875 const std::vector<void *> &v)
11882 std::vector<void *>::const_iterator i = v.begin();
11883 for (; hex != endh; ++hex, ++i)
11884 set_user_pointer(hex, *i);
11893 template <
int dim,
int spacedim>
11913 template <
int dim,
int spacedim>
11933 template <
int dim,
int spacedim>
11953 template <
int dim,
int spacedim>
11957 const unsigned int level =
levels.size() - 1;
11961 if (
levels[level]->cells.cells.size() == 0)
11968 levels[level]->cells.cells.size() - 1);
11971 if (ri->used() ==
true)
11974 if (ri->used() ==
true)
11981 template <
int dim,
int spacedim>
11991 if (cell->active() ==
true)
11994 if (cell->active() ==
true)
12002 template <
int dim,
int spacedim>
12013 template <
int dim,
int spacedim>
12018 if (level <
levels.size() - 1)
12025 template <
int dim,
int spacedim>
12029 if (level <
levels.size() - 1)
12030 return begin(level + 1);
12037 template <
int dim,
int spacedim>
12049 template <
int dim,
int spacedim>
12058 template <
int dim,
int spacedim>
12069 template <
int dim,
int spacedim>
12072 const unsigned int level)
const 12080 template <
int dim,
int spacedim>
12083 const unsigned int level)
const 12094 template <
int dim,
int spacedim>
12115 template <
int dim,
int spacedim>
12136 template <
int dim,
int spacedim>
12159 template <
int dim,
int spacedim>
12179 while (i->used() ==
false)
12188 template <
int dim,
int spacedim>
12197 template <
int dim,
int spacedim>
12218 template <
int dim,
int spacedim>
12219 typename Triangulation<dim, spacedim>::raw_line_iterator
12228 if (level >=
levels.size() ||
levels[level]->cells.cells.size() == 0)
12231 return raw_line_iterator(
12236 return raw_line_iterator(
12242 template <
int dim,
int spacedim>
12243 typename Triangulation<dim, spacedim>::line_iterator
12250 while (ri->used() ==
false)
12258 template <
int dim,
int spacedim>
12259 typename Triangulation<dim, spacedim>::active_line_iterator
12266 while (i->has_children())
12274 template <
int dim,
int spacedim>
12275 typename Triangulation<dim, spacedim>::line_iterator
12288 template <
int dim,
int spacedim>
12289 typename Triangulation<dim, spacedim>::raw_quad_iterator
12296 return raw_hex_iterator();
12302 if (level >=
levels.size() ||
levels[level]->cells.cells.size() == 0)
12305 return raw_quad_iterator(
12313 return raw_quad_iterator(
12320 return raw_hex_iterator();
12326 template <
int dim,
int spacedim>
12327 typename Triangulation<dim, spacedim>::quad_iterator
12334 while (ri->used() ==
false)
12342 template <
int dim,
int spacedim>
12343 typename Triangulation<dim, spacedim>::active_quad_iterator
12350 while (i->has_children())
12358 template <
int dim,
int spacedim>
12359 typename Triangulation<dim, spacedim>::quad_iterator
12371 template <
int dim,
int spacedim>
12372 typename Triangulation<dim, spacedim>::raw_hex_iterator
12380 return raw_hex_iterator();
12386 if (level >=
levels.size() ||
levels[level]->cells.cells.size() == 0)
12389 return raw_hex_iterator(
12395 return raw_hex_iterator();
12401 template <
int dim,
int spacedim>
12402 typename Triangulation<dim, spacedim>::hex_iterator
12409 while (ri->used() ==
false)
12417 template <
int dim,
int spacedim>
12418 typename Triangulation<dim, spacedim>::active_hex_iterator
12425 while (i->has_children())
12433 template <
int dim,
int spacedim>
12434 typename Triangulation<dim, spacedim>::hex_iterator
12449 namespace TriangulationImplementation
12451 inline unsigned int 12458 inline unsigned int 12466 inline unsigned int 12473 inline unsigned int 12481 inline unsigned int 12488 inline unsigned int 12499 template <
int dim,
int spacedim>
12503 return internal::TriangulationImplementation::n_cells(
number_cache);
12507 template <
int dim,
int spacedim>
12511 return internal::TriangulationImplementation::n_active_cells(
number_cache);
12514 template <
int dim,
int spacedim>
12523 template <
int dim,
int spacedim>
12542 template <
int dim,
int spacedim>
12559 template <
int dim,
int spacedim>
12578 template <
int dim,
int spacedim>
12598 template <
int dim,
int spacedim>
12618 template <
int dim,
int spacedim>
12637 template <
int dim,
int spacedim>
12649 template <
int dim,
int spacedim>
12663 return levels[level]->cells.cells.size();
12682 return levels[level]->cells.cells.size();
12700 return levels[level]->cells.cells.size();
12713 template <
int dim,
int spacedim>
12722 template <
int dim,
int spacedim>
12726 return faces->lines.cells.size();
12730 template <
int dim,
int spacedim>
12741 template <
int dim,
int spacedim>
12749 template <
int dim,
int spacedim>
12907 template <
int dim,
int spacedim>
12915 template <
int dim,
int spacedim>
12932 return levels[level]->cells.cells.size();
12942 return levels[level]->cells.cells.size();
12956 template <
int dim,
int spacedim>
12970 return faces->quads.cells.size();
12975 template <
int dim,
int spacedim>
12983 template <
int dim,
int spacedim>
12995 template <
int dim,
int spacedim>
13004 template <
int dim,
int spacedim>
13013 template <
int dim,
int spacedim>
13021 template <
int dim,
int spacedim>
13030 template <
int dim,
int spacedim>
13064 return levels[level]->cells.cells.size();
13089 template <
int dim,
int spacedim>
13095 std::bind(std::equal_to<bool>(),
13096 std::placeholders::_1,
13102 template <
int dim,
int spacedim>
13103 const std::vector<bool> &
13136 template <
int dim,
int spacedim>
13144 unsigned int max_vertex_index = 0;
13145 for (; cell != endc; ++cell)
13146 for (
unsigned int vertex = 0; vertex < GeometryInfo<dim>::vertices_per_cell;
13148 if (cell->vertex_index(vertex) > max_vertex_index)
13149 max_vertex_index = cell->vertex_index(vertex);
13155 std::vector<unsigned short int> usage_count(max_vertex_index + 1, 0);
13159 for (cell =
begin(); cell != endc; ++cell)
13160 for (
unsigned int vertex = 0; vertex < GeometryInfo<dim>::vertices_per_cell;
13162 ++usage_count[cell->vertex_index(vertex)];
13165 static_cast<unsigned int>(
13166 *std::max_element(usage_count.begin(), usage_count.end())));
13171 template <
int dim,
int spacedim>
13180 template <
int dim,
int spacedim>
13189 template <
int dim,
int spacedim>
13198 template <
int dim,
int spacedim>
13202 &periodicity_vector)
13205 periodicity_vector.begin(),
13206 periodicity_vector.end());
13214 template <
int dim,
int spacedim>
13215 const typename std::map<
13216 std::pair<typename Triangulation<dim, spacedim>::cell_iterator,
unsigned int>,
13217 std::pair<std::pair<typename Triangulation<dim, spacedim>::cell_iterator,
13227 template <
int dim,
int spacedim>
13249 if (
smooth_grid & limit_level_difference_at_vertices)
13254 update_neighbors(*
this);
13261 cells_with_distorted_children);
13268 template <
int dim,
int spacedim>
13272 unsigned int active_cell_index = 0;
13274 if ((cell->used() ==
false) || cell->has_children())
13278 cell->set_active_cell_index(active_cell_index);
13279 ++active_cell_index;
13286 template <
int dim,
int spacedim>
13293 typename std::vector<
13299 update_periodic_face_map_recursively<dim, spacedim>(it->cell[0],
13307 std::bitset<3> inverted_orientation;
13309 bool orientation, flip, rotation;
13311 rotation = it->orientation[2];
13312 flip = orientation ? rotation ^ it->orientation[1] : it->orientation[1];
13313 inverted_orientation[0] = orientation;
13314 inverted_orientation[1] = flip;
13315 inverted_orientation[2] = rotation;
13317 update_periodic_face_map_recursively<dim, spacedim>(it->cell[1],
13321 inverted_orientation,
13326 typename std::map<std::pair<cell_iterator, unsigned int>,
13327 std::pair<std::pair<cell_iterator, unsigned int>,
13328 std::bitset<3>>>::const_iterator it_test;
13333 it_test->first.first;
13335 it_test->second.first.first;
13336 if (cell_1->level() == cell_2->level())
13349 template <
int dim,
int spacedim>
13365 template <
int dim,
int spacedim>
13380 for (
unsigned int level = 0; level <
levels.size(); ++level)
13381 levels[level]->cells.monitor_memory(dim);
13388 while (cell != endc)
13392 return cells_with_distorted_children;
13397 template <
int dim,
int spacedim>
13405 std::vector<unsigned int> line_cell_count =
13406 count_cells_bounded_by_line(*
this);
13407 std::vector<unsigned int> quad_cell_count =
13408 count_cells_bounded_by_quad(*
this);
13423 for (; cell != endc; ++cell)
13424 if (!cell->active())
13425 if (cell->child(0)->coarsen_flag_set())
13427 cell->set_user_flag();
13428 for (
unsigned int child = 0; child < cell->n_children(); ++child)
13430 Assert(cell->child(child)->coarsen_flag_set(),
13432 cell->child(child)->clear_coarsen_flag();
13448 for (cell =
last(); cell != endc; --cell)
13449 if (cell->level() <=
static_cast<int>(
levels.size() - 2) &&
13450 cell->user_flag_set())
13465 for (cell =
begin(); cell != endc; ++cell)
13472 template <
int dim,
int spacedim>
13492 std::vector<int> vertex_level(
vertices.size(), 0);
13494 bool continue_iterating =
true;
13501 ExcMessage(
"In case of anisotropic refinement the " 13502 "limit_level_difference_at_vertices flag for " 13503 "mesh smoothing must not be set!"));
13507 std::fill(vertex_level.begin(), vertex_level.end(), 0);
13509 for (; cell != endc; ++cell)
13511 if (cell->refine_flag_set())
13512 for (
unsigned int vertex = 0;
13513 vertex < GeometryInfo<dim>::vertices_per_cell;
13515 vertex_level[cell->vertex_index(vertex)] =
13516 std::max(vertex_level[cell->vertex_index(vertex)],
13517 cell->level() + 1);
13518 else if (!cell->coarsen_flag_set())
13519 for (
unsigned int vertex = 0;
13520 vertex < GeometryInfo<dim>::vertices_per_cell;
13522 vertex_level[cell->vertex_index(vertex)] =
13523 std::max(vertex_level[cell->vertex_index(vertex)],
13534 for (
unsigned int vertex = 0;
13535 vertex < GeometryInfo<dim>::vertices_per_cell;
13537 vertex_level[cell->vertex_index(vertex)] =
13538 std::max(vertex_level[cell->vertex_index(vertex)],
13539 cell->level() - 1);
13554 if (cell->refine_flag_set() ==
false)
13556 for (
unsigned int vertex = 0;
13557 vertex < GeometryInfo<dim>::vertices_per_cell;
13559 if (vertex_level[cell->vertex_index(vertex)] >=
13563 cell->clear_coarsen_flag();
13568 if (vertex_level[cell->vertex_index(vertex)] >
13571 cell->set_refine_flag();
13573 for (
unsigned int v = 0;
13574 v < GeometryInfo<dim>::vertices_per_cell;
13576 vertex_level[cell->vertex_index(v)] =
13577 std::max(vertex_level[cell->vertex_index(v)],
13578 cell->level() + 1);
13604 for (; acell != end_ac; ++acell)
13605 acell->clear_coarsen_flag();
13608 for (; cell != endc; ++cell)
13611 if (cell->active())
13614 const unsigned int n_children = cell->n_children();
13615 unsigned int flagged_children = 0;
13616 for (
unsigned int child = 0; child < n_children; ++child)
13617 if (cell->child(child)->active() &&
13618 cell->child(child)->coarsen_flag_set())
13620 ++flagged_children;
13622 cell->child(child)->clear_coarsen_flag();
13627 if (flagged_children == n_children)
13628 cell->set_user_flag();
13634 for (cell =
begin(); cell != endc; ++cell)
13655 for (cell =
last(); cell != endc; --cell)
13656 if (cell->user_flag_set())
13660 template coarsening_allowed<dim, spacedim>(cell))
13661 for (
unsigned int c = 0; c < cell->n_children(); ++c)
13663 Assert(cell->child(c)->refine_flag_set() ==
false,
13666 cell->child(c)->set_coarsen_flag();
13679 continue_iterating = (current_coarsen_flags != previous_coarsen_flags);
13680 previous_coarsen_flags = current_coarsen_flags;
13682 while (continue_iterating ==
true);
13693 std::vector<bool> flags_before;
13699 std::vector<bool> flags_after;
13702 return (flags_before != flags_after);
13712 std::vector<bool> flags_before;
13718 std::vector<bool> flags_after;
13721 return (flags_before != flags_after);
13731 std::vector<bool> flags_before;
13737 std::vector<bool> flags_after;
13740 return (flags_before != flags_after);
13751 template <
int dim,
int spacedim>
13753 possibly_do_not_produce_unrefined_islands(
13758 unsigned int n_neighbors = 0;
13761 unsigned int count = 0;
13762 for (
unsigned int n = 0; n < GeometryInfo<dim>::faces_per_cell; ++n)
13769 if (face_will_be_refined_by_neighbor(cell, n))
13776 if (count == n_neighbors ||
13777 (count >= n_neighbors - 1 &&
13780 for (
unsigned int c = 0; c < cell->n_children(); ++c)
13781 cell->child(c)->clear_coarsen_flag();
13783 for (
unsigned int face = 0; face < GeometryInfo<dim>::faces_per_cell;
13785 if (!cell->at_boundary(face) && (!cell->neighbor(face)->active()) &&
13786 (cell_will_be_coarsened(cell->neighbor(face))))
13787 possibly_do_not_produce_unrefined_islands<dim, spacedim>(
13788 cell->neighbor(face));
13803 template <
int dim,
int spacedim>
13805 possibly_refine_unrefined_island(
13822 if (allow_anisotropic_smoothing ==
false)
13825 unsigned int refined_neighbors = 0, unrefined_neighbors = 0;
13826 for (
unsigned int face = 0; face < GeometryInfo<dim>::faces_per_cell;
13828 if (!cell->at_boundary(face))
13830 if (face_will_be_refined_by_neighbor(cell, face))
13831 ++refined_neighbors;
13833 ++unrefined_neighbors;
13836 if (unrefined_neighbors < refined_neighbors)
13838 cell->clear_coarsen_flag();
13839 cell->set_refine_flag();
13844 if (unrefined_neighbors > 0)
13845 for (
unsigned int face = 0;
13846 face < GeometryInfo<dim>::faces_per_cell;
13848 if (!cell->at_boundary(face) &&
13849 (face_will_be_refined_by_neighbor(cell, face) ==
false) &&
13850 (cell->neighbor(face)->has_children() ==
false) &&
13851 (cell->neighbor(face)->refine_flag_set() ==
false))
13852 possibly_refine_unrefined_island<dim, spacedim>(
13865 for (
unsigned int face_pair = 0;
13866 face_pair < GeometryInfo<dim>::faces_per_cell / 2;
13875 for (
unsigned int face_index = 0; face_index < 2; ++face_index)
13877 unsigned int face = 2 * face_pair + face_index;
13884 face_will_be_refined_by_neighbor<dim, spacedim>(
13885 cell, face, expected_face_ref_case);
13899 directional_cell_refinement_case =
13900 (directional_cell_refinement_case &
13903 expected_face_ref_case,
13905 cell->face_orientation(face),
13906 cell->face_flip(face),
13907 cell->face_rotation(face)));
13913 Assert(directional_cell_refinement_case <
13916 smoothing_cell_refinement_case =
13917 smoothing_cell_refinement_case | directional_cell_refinement_case;
13922 if (smoothing_cell_refinement_case)
13924 cell->clear_coarsen_flag();
13925 cell->set_refine_flag(cell->refine_flag_set() |
13926 smoothing_cell_refinement_case);
13933 template <
int dim,
int spacedim>
13939 std::vector<bool> flags_before[2];
13959 std::vector<bool> flags_before_loop[2] = {flags_before[0], flags_before[1]};
14041 for (; cell != endc; ++cell)
14042 cell->clear_coarsen_flag();
14045 bool mesh_changed_in_this_loop =
false;
14060 for (cell =
begin(); cell != endc; ++cell)
14064 if (!cell->active() && cell_will_be_coarsened(cell))
14065 possibly_do_not_produce_unrefined_islands<dim, spacedim>(cell);
14096 for (cell =
begin(); cell != endc; ++cell)
14097 if (!cell->active() || (cell->active() && cell->refine_flag_set() &&
14098 cell->is_locally_owned()))
14105 bool all_children_active =
true;
14106 if (!cell->active())
14107 for (
unsigned int c = 0; c < cell->n_children(); ++c)
14108 if (!cell->child(c)->active() ||
14109 cell->child(c)->is_ghost() ||
14110 cell->child(c)->is_artificial())
14112 all_children_active =
false;
14116 if (all_children_active)
14123 unsigned int unrefined_neighbors = 0, total_neighbors = 0;
14125 for (
unsigned int n = 0;
14126 n < GeometryInfo<dim>::faces_per_cell;
14134 if (!face_will_be_refined_by_neighbor(cell, n))
14135 ++unrefined_neighbors;
14150 if ((unrefined_neighbors == total_neighbors) &&
14151 (((unrefined_neighbors ==
14154 ((unrefined_neighbors <
14158 (total_neighbors != 0))
14160 if (!cell->active())
14161 for (
unsigned int c = 0; c < cell->n_children(); ++c)
14163 cell->child(c)->clear_refine_flag();
14164 cell->child(c)->set_coarsen_flag();
14167 cell->clear_refine_flag();
14185 ExcMessage(
"In case of anisotropic refinement the " 14186 "limit_level_difference_at_vertices flag for " 14187 "mesh smoothing must not be set!"));
14191 std::vector<int> vertex_level(
vertices.size(), 0);
14193 for (; cell != endc; ++cell)
14195 if (cell->refine_flag_set())
14196 for (
unsigned int vertex = 0;
14197 vertex < GeometryInfo<dim>::vertices_per_cell;
14199 vertex_level[cell->vertex_index(vertex)] =
14200 std::max(vertex_level[cell->vertex_index(vertex)],
14201 cell->level() + 1);
14202 else if (!cell->coarsen_flag_set())
14203 for (
unsigned int vertex = 0;
14204 vertex < GeometryInfo<dim>::vertices_per_cell;
14206 vertex_level[cell->vertex_index(vertex)] =
14207 std::max(vertex_level[cell->vertex_index(vertex)],
14216 for (
unsigned int vertex = 0;
14217 vertex < GeometryInfo<dim>::vertices_per_cell;
14219 vertex_level[cell->vertex_index(vertex)] =
14220 std::max(vertex_level[cell->vertex_index(vertex)],
14221 cell->level() - 1);
14236 if (cell->refine_flag_set() ==
false)
14238 for (
unsigned int vertex = 0;
14239 vertex < GeometryInfo<dim>::vertices_per_cell;
14241 if (vertex_level[cell->vertex_index(vertex)] >=
14245 cell->clear_coarsen_flag();
14250 if (vertex_level[cell->vertex_index(vertex)] >
14253 cell->set_refine_flag();
14255 for (
unsigned int v = 0;
14256 v < GeometryInfo<dim>::vertices_per_cell;
14258 vertex_level[cell->vertex_index(v)] =
14259 std::max(vertex_level[cell->vertex_index(v)],
14260 cell->level() + 1);
14285 for (; cell != endc; --cell)
14288 if (cell->refine_flag_set() !=
14290 possibly_refine_unrefined_island<dim, spacedim>(
14330 if (!cell->active())
14336 if (cell->child(0)->has_children() ==
true)
14343 for (
unsigned int i = 0; i < cell->n_children(); ++i)
14344 combined_ref_case =
14345 combined_ref_case | cell->child(i)->refine_flag_set();
14347 for (
unsigned int i = 0; i < cell->n_children(); ++i)
14351 child->clear_coarsen_flag();
14352 child->set_refine_flag(combined_ref_case);
14371 if (cell->active() || cell->child(0)->active())
14378 const unsigned int n_children = cell->n_children();
14379 bool has_active_grandchildren =
false;
14381 for (
unsigned int i = 0; i < n_children; ++i)
14382 if (cell->child(i)->child(0)->active())
14384 has_active_grandchildren =
true;
14388 if (has_active_grandchildren ==
false)
14394 unsigned int n_grandchildren = 0;
14397 unsigned int n_coarsen_flags = 0;
14402 for (
unsigned int c = 0; c < n_children; ++c)
14409 const unsigned int nn_children = child->n_children();
14410 n_grandchildren += nn_children;
14415 if (child->child(0)->active())
14416 for (
unsigned int cc = 0; cc < nn_children; ++cc)
14417 if (child->child(cc)->coarsen_flag_set())
14432 if ((n_coarsen_flags != n_grandchildren) && (n_coarsen_flags > 0))
14433 for (
unsigned int c = 0; c < n_children; ++c)
14436 if (child->child(0)->active())
14437 for (
unsigned int cc = 0; cc < child->n_children(); ++cc)
14438 child->child(cc)->clear_coarsen_flag();
14466 bool changed =
true;
14472 for (; cell != endc; --cell)
14473 if (cell->refine_flag_set())
14476 for (
unsigned int i = 0; i < GeometryInfo<dim>::faces_per_cell;
14484 cell->refine_flag_set(), i) !=
14497 if (cell->neighbor(i)->active())
14499 if (cell->neighbor_is_coarser(i))
14501 if (cell->neighbor(i)->coarsen_flag_set())
14502 cell->neighbor(i)->clear_coarsen_flag();
14522 ->flag_for_face_refinement(
14524 ->neighbor_of_coarser_neighbor(i)
14529 if (!cell->neighbor(i)
14530 ->refine_flag_set())
14532 cell->neighbor(i)->set_refine_flag();
14664 std::pair<unsigned int, unsigned int>
14666 cell->neighbor_of_coarser_neighbor(i);
14667 unsigned int refined_along_x = 0,
14668 refined_along_y = 0,
14669 to_be_refined_along_x = 0,
14670 to_be_refined_along_y = 0;
14672 const int this_face_index =
14673 cell->face_index(i);
14677 if ((this_face_index ==
14679 ->face(nb_indices.first)
14680 ->child_index(0)) ||
14681 (this_face_index ==
14683 ->face(nb_indices.first)
14694 ->face(nb_indices.first)
14695 ->refinement_case();
14713 cell->refine_flag_set(),
14715 cell->face_orientation(i),
14716 cell->face_flip(i),
14717 cell->face_rotation(i));
14720 ++to_be_refined_along_x;
14723 ++to_be_refined_along_y;
14729 cell->neighbor(i)->refine_flag_set())
14731 if (refined_along_x +
14732 to_be_refined_along_x >
14736 ->flag_for_face_refinement(
14740 if (refined_along_y +
14741 to_be_refined_along_y >
14745 ->flag_for_face_refinement(
14752 if (cell->neighbor(i)
14753 ->refine_flag_set() !=
14755 dim>::isotropic_refinement)
14757 cell->neighbor(i)->set_refine_flag();
14765 nb->refine_flag_set(),
14767 nb->face_orientation(nb_indices.first),
14768 nb->face_flip(nb_indices.first),
14769 nb->face_rotation(nb_indices.first));
14770 if ((nb_frc & RefinementCase<dim>::cut_x) &&
14771 !(refined_along_x ||
14772 to_be_refined_along_x))
14773 changed |= cell->flag_for_face_refinement(
14776 if ((nb_frc & RefinementCase<dim>::cut_y) &&
14777 !(refined_along_y ||
14778 to_be_refined_along_y))
14779 changed |= cell->flag_for_face_refinement(
14786 cell->neighbor(i)->clear_coarsen_flag();
14787 const unsigned int nb_nb =
14788 cell->neighbor_of_neighbor(i);
14793 neighbor->refine_flag_set(),
14795 neighbor->face_orientation(nb_nb),
14796 neighbor->face_flip(nb_nb),
14797 neighbor->face_rotation(nb_nb));
14800 cell->refine_flag_set(),
14802 cell->face_orientation(i),
14803 cell->face_flip(i),
14804 cell->face_rotation(i));
14809 if ((face_ref_case ==
14811 needed_face_ref_case ==
14815 needed_face_ref_case ==
14818 changed = cell->flag_for_face_refinement(
14820 neighbor->flag_for_face_refinement(
14821 nb_nb, needed_face_ref_case);
14828 face_ref_case = cell->face(i)->refinement_case(),
14829 needed_face_ref_case =
14831 cell->refine_flag_set(),
14833 cell->face_orientation(i),
14834 cell->face_flip(i),
14835 cell->face_rotation(i));
14840 needed_face_ref_case ==
14843 needed_face_ref_case ==
14846 cell->flag_for_face_refinement(i,
14870 std::vector<bool> flags_after_loop[2];
14876 mesh_changed_in_this_loop =
14877 ((flags_before_loop[0] != flags_after_loop[0]) ||
14878 (flags_before_loop[1] != flags_after_loop[1]));
14882 flags_before_loop[0].swap(flags_after_loop[0]);
14883 flags_before_loop[1].swap(flags_after_loop[1]);
14885 while (mesh_changed_in_this_loop);
14891 return ((flags_before[0] != flags_before_loop[0]) ||
14892 (flags_before[1] != flags_before_loop[1]));
14897 template <
int dim,
int spacedim>
14900 const unsigned int magic_number1,
14901 const std::vector<bool> &v,
14902 const unsigned int magic_number2,
14903 std::ostream & out)
14905 const unsigned int N = v.size();
14906 unsigned char * flags =
new unsigned char[N / 8 + 1];
14907 for (
unsigned int i = 0; i < N / 8 + 1; ++i)
14910 for (
unsigned int position = 0; position < N; ++position)
14911 flags[position / 8] |= (v[position] ? (1 << (position % 8)) : 0);
14920 out << magic_number1 <<
' ' << N << std::endl;
14921 for (
unsigned int i = 0; i < N / 8 + 1; ++i)
14922 out << static_cast<unsigned int>(flags[i]) <<
' ';
14924 out << std::endl << magic_number2 << std::endl;
14932 template <
int dim,
int spacedim>
14935 std::vector<bool> &v,
14936 const unsigned int magic_number2,
14941 unsigned int magic_number;
14942 in >> magic_number;
14949 unsigned char * flags =
new unsigned char[N / 8 + 1];
14950 unsigned short int tmp;
14951 for (
unsigned int i = 0; i < N / 8 + 1; ++i)
14957 for (
unsigned int position = 0; position != N; ++position)
14958 v[position] = (flags[position / 8] & (1 << (position % 8)));
14960 in >> magic_number;
14970 template <
int dim,
int spacedim>
14974 std::size_t mem = 0;
14976 for (
unsigned int i = 0; i <
levels.size(); ++i)
14983 mem +=
sizeof(
faces);
14992 template <
int dim,
int spacedim>
15003 #include "tria.inst" 15005 DEAL_II_NAMESPACE_CLOSE
std::vector< CellData< 1 > > boundary_lines
boost::signals2::signal< void(const Triangulation< dim, spacedim > &destination_tria)> copy
TriaIterator< TriaAccessor< dim-1, dim, spacedim >> face_iterator
::internal::TriangulationImplementation::NumberCache< dim > number_cache
unsigned int n_raw_faces() const
boost::signals2::signal< void()> post_refinement
active_line_iterator begin_active_line(const unsigned int level=0) const
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())
std::vector< unsigned int > n_quads_level
const types::manifold_id flat_manifold_id
static const unsigned int invalid_unsigned_int
static::ExceptionBase & ExcInvalidLevel(int arg1)
void load_user_flags_line(std::istream &in)
boost::signals2::signal< void(const typename Triangulation< dim, spacedim >::cell_iterator &cell)> pre_coarsening_on_cell
#define AssertNothrow(cond, exc)
#define DeclException2(Exception2, type1, type2, outsequence)
static::ExceptionBase & ExcInternalErrorOnCell(int arg1)
const types::subdomain_id invalid_subdomain_id
line_iterator begin_line(const unsigned int level=0) const
virtual void create_triangulation_compatibility(const std::vector< Point< spacedim >> &vertices, const std::vector< CellData< dim >> &cells, const SubCellData &subcelldata)
cell_iterator end() const
virtual unsigned int n_global_levels() 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
void save_user_flags_hex(std::ostream &out) const
unsigned int n_active_lines() const
active_face_iterator begin_active_face() const
unsigned int n_used_vertices() const
static::ExceptionBase & ExcIO()
Subscriptor & operator=(const Subscriptor &)
std::vector< GridTools::PeriodicFacePair< cell_iterator > > periodic_face_pairs_level_0
Task< RT > new_task(const std::function< RT()> &function)
TriaActiveIterator< TriaAccessor< dim-1, dim, spacedim >> active_face_iterator
bool check_consistency(const unsigned int dim) const
unsigned int n_faces() const
IteratorRange< active_cell_iterator > active_cell_iterators() const
std::vector< Point< spacedim > > vertices
std::map< types::manifold_id, std::unique_ptr< const Manifold< dim, spacedim > > > manifold
static unsigned int line_to_cell_vertices(const unsigned int line, const unsigned int vertex)
void set_all_refine_flags()
void load_user_flags(std::istream &in)
virtual std::unique_ptr< Manifold< dim, spacedim > > clone() const =0
bool anisotropic_refinement
static void read_bool_vector(const unsigned int magic_number1, std::vector< bool > &v, const unsigned int magic_number2, std::istream &in)
static::ExceptionBase & ExcGridReadError()
void clear_user_flags_quad()
face_iterator begin_face() const
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)
void save_user_pointers_quad(std::vector< void * > &v) const
static::ExceptionBase & ExcCellHasNegativeMeasure(int arg1)
void flip_all_direction_flags()
raw_hex_iterator begin_raw_hex(const unsigned int level=0) const
const std::vector< bool > & get_used_vertices() const
cell_iterator last() const
line_iterator end_line() 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)
static::ExceptionBase & ExcIndexRange(int arg1, int arg2, int arg3)
cell_iterator begin(const unsigned int level=0) const
void execute_coarsening()
static Triangulation< 2, spacedim >::DistortedCellList execute_refinement(Triangulation< 2, spacedim > &triangulation, const bool check_for_distorted_cells)
void load_user_flags_quad(std::istream &in)
void save_user_flags_line(std::ostream &out) const
unsigned int n_active_cells() const
unsigned long long int global_dof_index
void load_user_pointers_line(const std::vector< void * > &v)
unsigned int n_levels() const
void set_manifold(const types::manifold_id number, const Manifold< dim, spacedim > &manifold_object)
IteratorRange< cell_iterator > cell_iterators() const
virtual ~DistortedCellList() noexceptoverride
static::ExceptionBase & ExcFacesHaveNoLevel()
std::unique_ptr<::internal::TriangulationImplementation::TriaFaces< dim > > faces
bool get_anisotropic_refinement_flag() const
static::ExceptionBase & ExcInteriorQuadCantBeBoundary(int arg1, int arg2, int arg3, int arg4, types::boundary_id arg5)
vertex_iterator end_vertex() const
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)
virtual void execute_coarsening_and_refinement()
static::ExceptionBase & ExcGridHasInvalidCell(int arg1)
void save_user_indices_hex(std::vector< unsigned int > &v) const
unsigned int n_active_quads
unsigned int n_active_quads() const
void save_coarsen_flags(std::ostream &out) const
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 save_user_indices(std::vector< unsigned int > &v) const
virtual bool prepare_coarsening_and_refinement()
static::ExceptionBase & ExcMessage(std::string arg1)
unsigned int n_raw_hexs(const unsigned int level) const
boost::signals2::signal< void()> pre_refinement
static::ExceptionBase & ExcImpossibleInDim(int arg1)
unsigned int n_raw_lines() const
static::ExceptionBase & ExcTriangulationNotEmpty(int arg1, int arg2)
void reset_active_cell_indices()
DistortedCellList execute_refinement()
const bool check_for_distorted_cells
unsigned int subdomain_id
quad_iterator end_quad() const
static::ExceptionBase & ExcNonOrientableTriangulation()
Triangulation(const MeshSmoothing smooth_grid=none, const bool check_for_distorted_cells=false)
unsigned int n_cells() const
#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)
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)
virtual void create_triangulation(const std::vector< Point< spacedim >> &vertices, const std::vector< CellData< dim >> &cells, const SubCellData &subcelldata)
active_cell_iterator end_active(const unsigned int level) const
IteratorRange< cell_iterator > cell_iterators_on_level(const unsigned int level) const
#define Assert(cond, exc)
static void prevent_distorted_boundary_cells(const Triangulation< 1, spacedim > &)
std::vector< types::manifold_id > get_manifold_ids() const
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)
void save_user_flags(std::ostream &out) const
hex_iterator end_hex() const
unsigned int n_raw_quads() const
boost::signals2::signal< void()> clear
unsigned int n_hexs() const
void set_all_manifold_ids(const types::manifold_id number)
static::ExceptionBase & ExcInteriorLineCantBeBoundary(int arg1, int arg2, types::boundary_id arg3)
void load_user_pointers_quad(const std::vector< void * > &v)
active_vertex_iterator begin_active_vertex() const
std::list< typename Triangulation< dim, spacedim >::cell_iterator > distorted_cells
std::unique_ptr< std::map< unsigned int, types::manifold_id > > vertex_to_manifold_id_map_1d
std::vector< bool > vertices_used
virtual bool has_hanging_nodes() const
void reset_all_manifolds()
boost::signals2::signal< void()> mesh_movement
void load_user_indices(const std::vector< unsigned int > &v)
void save_user_pointers(std::vector< void * > &v) const
raw_line_iterator begin_raw_line(const unsigned int level=0) const
static::ExceptionBase & ExcMultiplySetLineInfoOfLine(int arg1, int arg2)
void clear_despite_subscriptions()
const Manifold< dim, spacedim > & get_manifold(const types::manifold_id number) const
TriaRawIterator< CellAccessor< dim, spacedim >> raw_cell_iterator
int face(const unsigned int i) const
std::vector< types::boundary_id > get_boundary_ids() const
hex_iterator begin_hex(const unsigned int level=0) const
raw_cell_iterator begin_raw(const unsigned int level=0) const
void save_refine_flags(std::ostream &out) const
unsigned int n_lines() const
std::map< std::pair< cell_iterator, unsigned int >, std::pair< std::pair< cell_iterator, unsigned int >, std::bitset< 3 > > > periodic_face_map
#define DeclException5(Exception5, type1, type2, type3, type4, type5, outsequence)
vertex_iterator begin_vertex() const
IteratorState::IteratorStates state() const
void save_user_indices_line(std::vector< unsigned int > &v) const
virtual void set_mesh_smoothing(const MeshSmoothing mesh_smoothing)
static void create_triangulation(const std::vector< Point< spacedim >> &v, const std::vector< CellData< 2 >> &cells, const SubCellData &subcelldata, Triangulation< 2, spacedim > &triangulation)
virtual types::subdomain_id locally_owned_subdomain() const
boost::signals2::signal< void()> create
active_cell_iterator last_active() const
std::unique_ptr< std::map< unsigned int, types::boundary_id > > vertex_to_boundary_id_map_1d
std::vector< std::unique_ptr< ::internal::TriangulationImplementation::TriaLevel< dim > > > levels
bool vertex_used(const unsigned int index) const
IteratorRange< active_cell_iterator > active_cell_iterators_on_level(const unsigned int level) const
void save_user_pointers_line(std::vector< void * > &v) const
active_quad_iterator begin_active_quad(const unsigned int level=0) const
virtual types::global_dof_index n_global_active_cells() const
virtual const MeshSmoothing & get_mesh_smoothing() 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)
void swap(Vector< Number > &u, Vector< Number > &v)
void reset_manifold(const types::manifold_id manifold_number)
static void alternating_form_at_vertices(const Point< spacedim >(&vertices)[vertices_per_cell], Tensor< spacedim-dim, spacedim >(&forms)[vertices_per_cell])
void save_user_indices_quad(std::vector< unsigned int > &v) const
unsigned int n_raw_cells(const unsigned int level) const
static bool coarsening_allowed(const typename Triangulation< dim, spacedim >::cell_iterator &cell)
const types::manifold_id invalid_manifold_id
TriaObjects< TriaObject< 1 > > lines
unsigned int n_active_faces() const
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
void save_user_pointers_hex(std::vector< void * > &v) const
static void create_triangulation(const std::vector< Point< spacedim >> &v, const std::vector< CellData< 1 >> &cells, const SubCellData &, Triangulation< 1, spacedim > &triangulation)
void update_periodic_face_map()
std::vector< unsigned int > n_active_hexes_level
face_iterator end_face() const
static void compute_number_cache(const Triangulation< dim, spacedim > &triangulation, const unsigned int level_objects, internal::TriangulationImplementation::NumberCache< 3 > &number_cache)
raw_cell_iterator end_raw(const unsigned int level) const
std::vector< CellData< 2 > > boundary_quads
MeshSmoothing smooth_grid
virtual ~Triangulation() override
#define DeclException4(Exception4, type1, type2, type3, type4, outsequence)
quad_iterator begin_quad(const unsigned int level=0) const
unsigned int n_quads() const
std::vector< unsigned int > n_active_lines_level
active_hex_iterator begin_active_hex(const unsigned int level=0) const
void refine_global(const unsigned int times=1)
static::ExceptionBase & ExcNotImplemented()
void clear_user_data(const unsigned int i)
Iterator points to a valid object.
#define DeclException3(Exception3, type1, type2, type3, outsequence)
void load_user_pointers_hex(const std::vector< void * > &v)
void set_all_manifold_ids_on_boundary(const types::manifold_id number)
static::ExceptionBase & ExcQuadInexistant(int arg1, int arg2, int arg3, int arg4)
const types::boundary_id internal_face_boundary_id
void load_user_pointers(const std::vector< void * > &v)
active_cell_iterator begin_active(const unsigned int level=0) const
std::vector< unsigned int > n_active_quads_level
void clear_user_flags_line()
static::ExceptionBase & ExcBoundaryIdNotFound(types::boundary_id arg1)
unsigned int max_adjacent_cells() const
static Triangulation< 1, spacedim >::DistortedCellList execute_refinement(Triangulation< 1, spacedim > &triangulation, const bool)
void set_face(const unsigned int i, const int index)
unsigned int n_vertices() const
Triangulation & operator=(Triangulation< dim, spacedim > &&tria) noexcept
std::vector< unsigned int > n_hexes_level
static void delete_children(Triangulation< 1, spacedim > &triangulation, typename Triangulation< 1, spacedim >::cell_iterator &cell, std::vector< unsigned int > &, std::vector< unsigned int > &)
unsigned int n_active_hexs() const
raw_quad_iterator begin_raw_quad(const unsigned int level=0) 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)
TriaIterator< CellAccessor< dim, spacedim >> cell_iterator
TriaActiveIterator< CellAccessor< dim, spacedim >> active_cell_iterator
void load_user_indices_quad(const std::vector< unsigned int > &v)
boost::signals2::signal< void()> any_change
std::enable_if< std::is_fundamental< T >::value, std::size_t >::type memory_consumption(const T &t)
void save_user_flags_quad(std::ostream &out) const
void load_refine_flags(std::istream &in)
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)
virtual std::size_t memory_consumption() const