17 #include <deal.II/base/std_cxx14/memory.h> 18 #include <deal.II/base/memory_consumption.h> 19 #include <deal.II/base/table.h> 20 #include <deal.II/base/geometry_info.h> 22 #include <deal.II/grid/tria.h> 23 #include <deal.II/grid/tria_levels.h> 24 #include <deal.II/grid/tria_faces.h> 25 #include <deal.II/grid/manifold.h> 26 #include <deal.II/grid/tria_boundary.h> 27 #include <deal.II/grid/tria_accessor.h> 28 #include <deal.II/grid/tria_iterator.h> 29 #include <deal.II/grid/grid_tools.h> 30 #include <deal.II/grid/magic_numbers.h> 31 #include <deal.II/fe/mapping_q1.h> 32 #include <deal.II/lac/vector.h> 33 #include <deal.II/lac/full_matrix.h> 44 DEAL_II_NAMESPACE_OPEN
63 namespace TriangulationImplementation
139 template <
int dim,
int spacedim>
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) || (n_active_children == cell->n_children());
159 template <
int dim,
int spacedim>
165 if (cell->has_children())
167 unsigned int children_to_coarsen=0;
168 const unsigned int n_children=cell->n_children();
170 for (
unsigned int c=0; c<n_children; ++c)
171 if (cell->child(c)->active() &&
172 cell->child(c)->coarsen_flag_set())
173 ++children_to_coarsen;
174 if (children_to_coarsen==n_children)
177 for (
unsigned int c=0; c<n_children; ++c)
178 if (cell->child(c)->active())
179 cell->child(c)->clear_coarsen_flag();
215 template <
int dim,
int spacedim>
218 const unsigned int face_no,
233 if (neighbor->has_children())
238 if (cell_will_be_coarsened(neighbor))
247 expected_face_ref_case=cell->face(face_no)->refinement_case();
259 const unsigned int neighbor_neighbor=cell->neighbor_face_no(face_no);
266 neighbor->face_orientation(neighbor_neighbor),
267 neighbor->face_flip(neighbor_neighbor),
268 neighbor->face_rotation(neighbor_neighbor));
272 const int this_face_index=cell->face_index(face_no);
278 if (neighbor_face->index()==this_face_index)
284 expected_face_ref_case = face_ref_case;
303 for (
unsigned int c=0; c<neighbor_face->n_children(); ++c)
304 if (neighbor_face->child_index(c)==this_face_index)
313 if ((neighbor_face->refinement_case() | face_ref_case)
314 == neighbor_face->refinement_case())
330 expected_face_ref_case = ~neighbor_face->refinement_case();
357 template <
int dim,
int spacedim>
360 const unsigned int face_no)
363 return face_will_be_refined_by_neighbor_internal(cell, face_no, dummy);
370 template <
int dim,
int spacedim>
373 const unsigned int face_no,
376 return face_will_be_refined_by_neighbor_internal(cell, face_no,
377 expected_face_ref_case);
382 template <
int dim,
int spacedim>
386 std::vector<unsigned int> min_adjacent_cell_level (triangulation.n_vertices(),
387 triangulation.n_levels());
388 std::vector<unsigned int> max_adjacent_cell_level (triangulation.n_vertices(),
392 cell = triangulation.begin_active();
393 cell != triangulation.end(); ++cell)
394 for (
unsigned int v=0; v<GeometryInfo<dim>::vertices_per_cell; ++v)
396 min_adjacent_cell_level[cell->vertex_index(v)]
397 = std::min<unsigned int>
398 (min_adjacent_cell_level[cell->vertex_index(v)],
400 max_adjacent_cell_level[cell->vertex_index(v)]
401 = std::max<unsigned int> (min_adjacent_cell_level[cell->vertex_index(v)],
405 for (
unsigned int k=0; k<triangulation.n_vertices(); ++k)
406 if (triangulation.vertex_used(k))
407 if (max_adjacent_cell_level[k] -
408 min_adjacent_cell_level[k] > 1)
421 template <
int dim,
int spacedim>
422 std::vector<unsigned int>
427 std::vector<unsigned int> line_cell_count(triangulation.n_raw_lines(),0);
429 cell=triangulation.begin(),
430 endc=triangulation.end();
431 for (; cell!=endc; ++cell)
432 for (
unsigned int l=0; l<GeometryInfo<dim>::lines_per_cell; ++
l)
433 ++line_cell_count[cell->line_index(l)];
434 return line_cell_count;
437 return std::vector<unsigned int>();
449 template <
int dim,
int spacedim>
450 std::vector<unsigned int>
455 std::vector<unsigned int> quad_cell_count (triangulation.n_raw_quads(),0);
457 cell=triangulation.begin(),
458 endc=triangulation.end();
459 for (; cell!=endc; ++cell)
460 for (
unsigned int q=0; q<GeometryInfo<dim>::faces_per_cell; ++q)
461 ++quad_cell_count[cell->quad_index(q)];
462 return quad_cell_count;
465 return std::vector<unsigned int>();
482 reorder_compatibility (
const std::vector<
CellData<1> > &,
491 reorder_compatibility (std::vector<
CellData<2> > &cells,
494 for (
unsigned int cell=0; cell<cells.size(); ++cell)
495 std::swap(cells[cell].vertices[2],cells[cell].vertices[3]);
500 reorder_compatibility (std::vector<
CellData<3> > &cells,
504 for (
unsigned int cell=0; cell<cells.size(); ++cell)
506 for (
unsigned int i=0; i<GeometryInfo<3>::vertices_per_cell; ++i)
507 tmp[i] = cells[cell].vertices[i];
508 for (
unsigned int i=0; i<GeometryInfo<3>::vertices_per_cell; ++i)
513 std::vector<CellData<2> >::iterator boundary_quad
515 std::vector<CellData<2> >::iterator end_quad
517 for (
unsigned int quad_no=0; boundary_quad!=end_quad; ++boundary_quad, ++quad_no)
518 std::swap(boundary_quad->vertices[2], boundary_quad->vertices[3]);
540 template <
int dim,
int spacedim>
542 middle_vertex_index(
const typename Triangulation<dim,spacedim>::line_iterator &line)
544 if (line->has_children())
545 return line->child(0)->vertex_index(1);
550 template <
int dim,
int spacedim>
552 middle_vertex_index(
const typename Triangulation<dim,spacedim>::quad_iterator &quad)
554 switch (static_cast<unsigned char> (quad->refinement_case()))
557 return middle_vertex_index<dim,spacedim>(quad->child(0)->line(1));
560 return middle_vertex_index<dim,spacedim>(quad->child(0)->line(3));
563 return quad->child(0)->vertex_index(3);
572 template <
int dim,
int spacedim>
574 middle_vertex_index(
const typename Triangulation<dim,spacedim>::hex_iterator &hex)
576 switch (static_cast<unsigned char> (hex->refinement_case()))
579 return middle_vertex_index<dim,spacedim>(hex->child(0)->quad(1));
582 return middle_vertex_index<dim,spacedim>(hex->child(0)->quad(3));
585 return middle_vertex_index<dim,spacedim>(hex->child(0)->quad(5));
588 return middle_vertex_index<dim,spacedim>(hex->child(0)->line(11));
591 return middle_vertex_index<dim,spacedim>(hex->child(0)->line(5));
594 return middle_vertex_index<dim,spacedim>(hex->child(0)->line(7));
597 return hex->child(0)->vertex_index(7);
618 template <
class TRIANGULATION>
620 typename TRIANGULATION::DistortedCellList
621 collect_distorted_coarse_cells (
const TRIANGULATION &)
623 return typename TRIANGULATION::DistortedCellList();
643 cell = triangulation.begin(0); cell != triangulation.end(0); ++cell)
646 for (
unsigned int i=0; i<GeometryInfo<dim>::vertices_per_cell; ++i)
647 vertices[i] = cell->vertex(i);
653 for (
unsigned int i=0; i<GeometryInfo<dim>::vertices_per_cell; ++i)
654 if (determinants[i] <= 1e-9 * std::pow (cell->diameter(),
657 distorted_cells.distorted_cells.push_back (cell);
662 return distorted_cells;
675 std::integral_constant<int, dim>,
676 std::integral_constant<int, dim>)
680 for (
unsigned int c=0; c<cell->n_children(); ++c)
683 for (
unsigned int i=0; i<GeometryInfo<dim>::vertices_per_cell; ++i)
684 vertices[i] = cell->child(c)->vertex(i);
690 for (
unsigned int i=0; i<GeometryInfo<dim>::vertices_per_cell; ++i)
691 if (determinants[i] <= 1e-9 * std::pow (cell->child(c)->diameter(),
707 template <
int dim,
int spacedim>
710 std::integral_constant<int, dim>,
711 std::integral_constant<int, spacedim>)
722 template <
int spacedim>
729 template <
int dim,
int spacedim>
784 static const unsigned int left_right_offset[2][6][2] =
813 std::vector<typename Triangulation<dim,spacedim>::cell_iterator>
814 adjacent_cells(2*triangulation.n_raw_faces(), dummy);
817 cell = triangulation.begin(),
818 endc = triangulation.end();
819 for (; cell != endc; ++cell)
820 for (
unsigned int f=0; f<GeometryInfo<dim>::faces_per_cell; ++f)
826 offset = (cell->direction_flag()
828 left_right_offset[dim-2][f][cell->face_orientation(f)]
830 1-left_right_offset[dim-2][f][cell->face_orientation(f)]);
832 adjacent_cells[2*face->index() + offset] = cell;
842 if (cell->active() && face->has_children())
844 adjacent_cells[2*face->child(0)->index() + offset] = cell;
845 adjacent_cells[2*face->child(1)->index() + offset] = cell;
870 if (face->has_children() &&
875 for (
unsigned int c=0; c<face->n_children(); ++c)
876 adjacent_cells[2*face->child(c)->index() + offset] = cell;
877 if (face->child(0)->has_children())
879 adjacent_cells[2*face->child(0)->child(0)->index() + offset] = cell;
880 adjacent_cells[2*face->child(0)->child(1)->index() + offset] = cell;
882 if (face->child(1)->has_children())
884 adjacent_cells[2*face->child(1)->child(0)->index() + offset] = cell;
885 adjacent_cells[2*face->child(1)->child(1)->index() + offset] = cell;
896 for (cell=triangulation.begin(); cell != endc; ++cell)
897 for (
unsigned int f=0; f<GeometryInfo<dim>::faces_per_cell; ++f)
900 offset = (cell->direction_flag()
902 left_right_offset[dim-2][f][cell->face_orientation(f)]
904 1-left_right_offset[dim-2][f][cell->face_orientation(f)]);
905 cell->set_neighbor(f,
906 adjacent_cells[2*cell->face(f)->index() + 1 - offset]);
911 template <
int dim,
int spacedim>
912 void update_periodic_face_map_recursively
915 unsigned int n_face_1,
unsigned int n_face_2,
916 const std::bitset<3> &orientation,
921 const FaceIterator face_1 = cell_1->face(n_face_1);
922 const FaceIterator face_2 = cell_2->face(n_face_2);
924 const bool face_orientation = orientation[0];
925 const bool face_flip = orientation[1];
926 const bool face_rotation = orientation[2];
928 Assert((dim != 1) || (face_orientation ==
true && face_flip ==
false && face_rotation ==
false),
929 ExcMessage (
"The supplied orientation " "(face_orientation, face_flip, face_rotation) " "is invalid for 1D"));
931 Assert((dim != 2) || (face_orientation ==
true && face_rotation ==
false),
932 ExcMessage (
"The supplied orientation " "(face_orientation, face_flip, face_rotation) " "is invalid for 2D"));
937 Assert(face_1->at_boundary() &&face_2->at_boundary(),
938 ExcMessage (
"Periodic faces must be on the boundary"));
941 typedef std::pair<typename Triangulation<dim,spacedim>::cell_iterator,
unsigned int> CellFace;
942 const CellFace cell_face_1 (cell_1, n_face_1);
943 const CellFace cell_face_2 (cell_2, n_face_2);
944 const std::pair<CellFace, std::bitset<3> > cell_face_orientation_2 (cell_face_2, orientation);
946 const std::pair<CellFace, std::pair<CellFace, std::bitset<3> > > periodic_faces (cell_face_1, cell_face_orientation_2);
950 periodic_face_map.insert(periodic_faces);
956 static const int lookup_table_2d[2][2] =
963 static const int lookup_table_3d[2][2][2][4] =
982 if (cell_1->has_children())
984 if (cell_2->has_children())
993 for (
unsigned int i = 0; i < GeometryInfo<dim>::max_children_per_face; ++i)
1000 j = lookup_table_2d[face_flip][i];
1003 j = lookup_table_3d[face_orientation][face_flip][face_rotation][i];
1010 unsigned int child_cell_1
1012 (cell_1->refinement_case(), n_face_1, i, cell_1->face_orientation(n_face_1),
1013 cell_1->face_flip(n_face_1), cell_1->face_rotation(n_face_1), face_1->refinement_case());
1014 unsigned int child_cell_2
1016 (cell_2->refinement_case(), n_face_2, j, cell_2->face_orientation(n_face_2),
1017 cell_2->face_flip(n_face_2), cell_2->face_rotation(n_face_2), face_2->refinement_case());
1024 update_periodic_face_map_recursively<dim, spacedim>
1025 (cell_1->child(child_cell_1), cell_2->child(child_cell_2),
1026 n_face_1, n_face_2, orientation, periodic_face_map);
1031 for (
unsigned int i = 0; i < GeometryInfo<dim>::max_children_per_face; ++i)
1034 unsigned int child_cell_1
1036 cell_1->face_flip(n_face_1), cell_1->face_rotation(n_face_1), face_1->refinement_case());
1039 update_periodic_face_map_recursively<dim, spacedim>
1040 (cell_1->child(child_cell_1), cell_2,
1041 n_face_1, n_face_2, orientation, periodic_face_map);
1054 namespace TriangulationImplementation
1062 using ::Triangulation;
1070 <<
"Something went wrong when making cell " << arg1
1071 <<
". Read the docs and the source code " 1072 <<
"for more information.");
1079 <<
"Something went wrong upon construction of cell " 1092 <<
"Cell " << arg1 <<
" has negative measure. This typically " 1093 <<
"indicates some distortion in the cell, or a mistakenly " 1094 <<
"swapped pair of vertices in the input to " 1095 <<
"Triangulation::create_triangulation().");
1105 <<
"Error while creating cell " << arg1
1106 <<
": the vertex index " << arg2 <<
" must be between 0 and " 1114 <<
"While trying to assign a boundary indicator to a line: " 1115 <<
"the line with end vertices " << arg1 <<
" and " 1116 << arg2 <<
" does not exist.");
1123 <<
"While trying to assign a boundary indicator to a quad: " 1124 <<
"the quad with bounding lines " << arg1 <<
", " << arg2
1125 <<
", " << arg3 <<
", " << arg4 <<
" does not exist.");
1133 <<
"The input data for creating a triangulation contained " 1134 <<
"information about a line with indices " 1135 << arg1 <<
" and " << arg2
1136 <<
" that is described to have boundary indicator " 1138 <<
". However, this is an internal line not located on the " 1139 <<
"boundary. You cannot assign a boundary indicator to it." 1142 <<
"If this happened at a place where you call " 1143 <<
"Triangulation::create_triangulation() yourself, you need " 1144 <<
"to check the SubCellData object you pass to this function." 1147 <<
"If this happened in a place where you are reading a mesh " 1148 <<
"from a file, then you need to investigate why such a line " 1149 <<
"ended up in the input file. A typical case is a geometry " 1150 <<
"that consisted of multiple parts and for which the mesh " 1151 <<
"generator program assumes that the interface between " 1152 <<
"two parts is a boundary when that isn't supposed to be " 1153 <<
"the case, or where the mesh generator simply assigns " 1154 <<
"'geometry indicators' to lines at the perimeter of " 1155 <<
"a part that are not supposed to be interpreted as " 1156 <<
"'boundary indicators'.");
1164 <<
"The input data for creating a triangulation contained " 1165 <<
"information about a quad with indices " 1166 << arg1 <<
", " << arg2 <<
", " << arg3 <<
", and " << arg4
1167 <<
" that is described to have boundary indicator " 1169 <<
". However, this is an internal quad not located on the " 1170 <<
"boundary. You cannot assign a boundary indicator to it." 1173 <<
"If this happened at a place where you call " 1174 <<
"Triangulation::create_triangulation() yourself, you need " 1175 <<
"to check the SubCellData object you pass to this function." 1178 <<
"If this happened in a place where you are reading a mesh " 1179 <<
"from a file, then you need to investigate why such a quad " 1180 <<
"ended up in the input file. A typical case is a geometry " 1181 <<
"that consisted of multiple parts and for which the mesh " 1182 <<
"generator program assumes that the interface between " 1183 <<
"two parts is a boundary when that isn't supposed to be " 1184 <<
"the case, or where the mesh generator simply assigns " 1185 <<
"'geometry indicators' to quads at the surface of " 1186 <<
"a part that are not supposed to be interpreted as " 1187 <<
"'boundary indicators'.");
1194 <<
"In SubCellData the line info of the line with vertex indices " 1195 << arg1 <<
" and " << arg2 <<
" appears more than once. " 1196 <<
"This is not allowed.");
1307 template <
int dim,
int spacedim>
1310 const unsigned int level_objects,
1314 typename Triangulation<dim,spacedim>::line_iterator line_iterator;
1317 if (level_objects > 0)
1319 for (
unsigned int level=0; level<level_objects; ++level)
1320 if (triangulation.begin(level) !=
1321 triangulation.end(level))
1340 for (
unsigned int level=0; level<number_cache.
n_levels; ++level)
1346 line_iterator line = triangulation.begin_line (level),
1347 endc = (level == number_cache.
n_levels-1 ?
1348 line_iterator(triangulation.end_line()) :
1349 triangulation.begin_line (level+1));
1350 for (; line!=endc; ++line)
1353 if (line->has_children() ==
false)
1368 line_iterator line = triangulation.begin_line (),
1369 endc = triangulation.end_line();
1370 for (; line!=endc; ++line)
1373 if (line->has_children() ==
false)
1393 template <
int dim,
int spacedim>
1396 const unsigned int level_objects,
1406 (&compute_number_cache<dim,spacedim>),
1413 typename Triangulation<dim,spacedim>::quad_iterator quad_iterator;
1429 unsigned int n_levels = 0;
1430 if (level_objects > 0)
1432 for (
unsigned int level=0; level<level_objects; ++level)
1433 if (triangulation.begin(level) !=
1434 triangulation.end(level))
1440 for (
unsigned int level=0; level<n_levels; ++level)
1446 quad_iterator quad = triangulation.begin_quad (level),
1447 endc = (level == n_levels-1 ?
1448 quad_iterator(triangulation.end_quad()) :
1449 triangulation.begin_quad (level+1));
1450 for (; quad!=endc; ++quad)
1453 if (quad->has_children() ==
false)
1468 quad_iterator quad = triangulation.begin_quad (),
1469 endc = triangulation.end_quad();
1470 for (; quad!=endc; ++quad)
1473 if (quad->has_children() ==
false)
1479 update_lines.
join ();
1497 template <
int dim,
int spacedim>
1500 const unsigned int level_objects,
1510 (&compute_number_cache<dim,spacedim>),
1517 typename Triangulation<dim,spacedim>::hex_iterator hex_iterator;
1534 unsigned int n_levels = 0;
1535 if (level_objects > 0)
1537 for (
unsigned int level=0; level<level_objects; ++level)
1538 if (triangulation.begin(level) !=
1539 triangulation.end(level))
1545 for (
unsigned int level=0; level<n_levels; ++level)
1551 hex_iterator hex = triangulation.begin_hex (level),
1552 endc = (level == n_levels-1 ?
1553 hex_iterator(triangulation.end_hex()) :
1554 triangulation.begin_hex (level+1));
1555 for (; hex!=endc; ++hex)
1558 if (hex->has_children() ==
false)
1573 hex_iterator hex = triangulation.begin_hex (),
1574 endc = triangulation.end_hex();
1575 for (; hex!=endc; ++hex)
1578 if (hex->has_children() ==
false)
1584 update_quads_and_lines.
join ();
1595 template <
int spacedim>
1612 const unsigned int dim=1;
1615 triangulation.vertices = v;
1616 triangulation.vertices_used = std::vector<bool> (v.size(),
true);
1623 if (dim == spacedim)
1625 for (
unsigned int cell_no = 0; cell_no<cells.size(); ++cell_no)
1631 if (!triangulation.check_for_distorted_cells)
1633 const double cell_measure = GridTools::cell_measure<1>
1634 (triangulation.vertices, cells[cell_no].vertices);
1645 std::vector<std::vector<int> > lines_at_vertex (v.size());
1649 triangulation.levels[0]->reserve_space (cells.size(), dim, spacedim);
1650 triangulation.levels[0]->cells.reserve_space (0, cells.size());
1653 typename Triangulation<dim,spacedim>::raw_line_iterator
1654 next_free_line = triangulation.begin_raw_line ();
1655 for (
unsigned int cell=0; cell<cells.size(); ++cell)
1657 while (next_free_line->used())
1662 cells[cell].vertices[1]));
1663 next_free_line->set_used_flag ();
1664 next_free_line->set_material_id (cells[cell].material_id);
1665 next_free_line->set_manifold_id (cells[cell].manifold_id);
1667 next_free_line->set_subdomain_id (0);
1671 lines_at_vertex[cells[cell].vertices[0]].push_back (cell);
1672 lines_at_vertex[cells[cell].vertices[1]].push_back (cell);
1678 unsigned int boundary_nodes = 0;
1679 for (
unsigned int i=0; i<lines_at_vertex.size(); ++i)
1680 switch (lines_at_vertex[i].size())
1691 ExcMessage (
"You have a vertex in your triangulation " 1692 "at which more than two cells come together. " 1693 "(For one dimensional triangulation, cells are " 1696 "This is not currently supported because the " 1697 "Triangulation class makes the assumption that " 1698 "every cell has zero or one neighbors behind " 1699 "each face (here, behind each vertex), but in your " 1700 "situation there would be more than one." 1702 "Support for this is not currently implemented. " 1703 "If you need to work with triangulations where " 1704 "more than two cells come together at a vertex, " 1705 "duplicate the vertices once per cell (i.e., put " 1706 "multiple vertices at the same physical location, " 1707 "but using different vertex indices for each) " 1708 "and then ensure continuity of the solution by " 1709 "explicitly creating constraints that the degrees " 1710 "of freedom at these vertices have the same " 1711 "value, using the ConstraintMatrix class."));
1723 AssertThrow (((spacedim == 1) && (boundary_nodes == 2))
1726 ExcMessage(
"The Triangulation has too many end points"));
1732 typename Triangulation<dim,spacedim>::active_line_iterator
1733 line = triangulation.begin_active_line ();
1735 for (; line!=triangulation.end(); ++line)
1737 for (
unsigned int vertex=0; vertex<GeometryInfo<dim>::vertices_per_cell; ++vertex)
1743 if (lines_at_vertex[line->vertex_index(vertex)][0] == line->index())
1744 if (lines_at_vertex[line->vertex_index(vertex)].size() == 2)
1747 neighbor (&triangulation,
1749 lines_at_vertex[line->vertex_index(vertex)][1]);
1750 line->set_neighbor (vertex, neighbor);
1756 line->set_neighbor (vertex, triangulation.end());
1763 neighbor (&triangulation,
1765 lines_at_vertex[line->vertex_index(vertex)][0]);
1766 line->set_neighbor (vertex, neighbor);
1773 triangulation.vertex_to_boundary_id_map_1d->clear();
1774 triangulation.vertex_to_manifold_id_map_1d->clear();
1776 cell = triangulation.begin_active();
1777 cell != triangulation.end(); ++cell)
1778 for (
unsigned int f=0; f<GeometryInfo<dim>::faces_per_cell; ++f)
1781 .vertex_to_manifold_id_map_1d)[cell->face(f)->vertex_index()]
1784 if (cell->at_boundary(f))
1786 .vertex_to_boundary_id_map_1d)[cell->face(f)->vertex_index()]
1799 template <
int spacedim>
1810 const unsigned int dim=2;
1813 triangulation.vertices = v;
1814 triangulation.vertices_used = std::vector<bool> (v.size(),
true);
1821 if (dim == spacedim)
1823 for (
unsigned int cell_no = 0; cell_no<cells.size(); ++cell_no)
1826 if (!triangulation.check_for_distorted_cells)
1828 const double cell_measure = GridTools::cell_measure<2>
1829 (triangulation.vertices, cells[cell_no].vertices);
1849 std::map<std::pair<int,int>,
1850 typename Triangulation<dim,spacedim>::line_iterator> needed_lines;
1851 for (
unsigned int cell=0; cell<cells.size(); ++cell)
1853 for (
unsigned int vertex=0; vertex<4; ++vertex)
1854 AssertThrow (cells[cell].vertices[vertex] < triangulation.vertices.size(),
1856 triangulation.vertices.size()));
1858 for (
unsigned int line=0; line<GeometryInfo<dim>::faces_per_cell; ++line)
1863 std::pair<int,int> line_vertices(
1883 AssertThrow (needed_lines.find(std::make_pair(line_vertices.second,
1884 line_vertices.first))
1893 needed_lines[line_vertices] = triangulation.end_line();
1901 std::vector<unsigned short int> vertex_touch_count (v.size(), 0);
1902 typename std::map<std::pair<int,int>,
1903 typename Triangulation<dim,spacedim>::line_iterator>::iterator i;
1904 for (i=needed_lines.begin(); i!=needed_lines.end(); ++i)
1908 ++vertex_touch_count[i->first.first];
1909 ++vertex_touch_count[i->first.second];
1916 AssertThrow (* (std::min_element(vertex_touch_count.begin(),
1917 vertex_touch_count.end())) >= 2,
1918 ExcMessage(
"During creation of a triangulation, a part of the " 1919 "algorithm encountered a vertex that is part of only " 1920 "a single adjacent line. However, in 2d, every vertex " 1921 "needs to be at least part of two lines."));
1926 triangulation.faces = std_cxx14::make_unique<internal::TriangulationImplementation::TriaFaces<dim>> ();
1927 triangulation.levels[0]->reserve_space (cells.size(), dim, spacedim);
1928 triangulation.faces->lines.reserve_space (0,needed_lines.size());
1929 triangulation.levels[0]->cells.reserve_space (0,cells.size());
1933 typename Triangulation<dim,spacedim>::raw_line_iterator
1934 line = triangulation.begin_raw_line();
1935 typename std::map<std::pair<int,int>,
1936 typename Triangulation<dim,spacedim>::line_iterator>::iterator i;
1937 for (i = needed_lines.begin();
1938 line!=triangulation.end_line(); ++line, ++i)
1942 line->set_used_flag ();
1943 line->clear_user_flag ();
1952 std::map<int,std::vector<typename Triangulation<dim,spacedim>::cell_iterator> >
1958 cell = triangulation.begin_raw_quad();
1959 for (
unsigned int c=0; c<cells.size(); ++c, ++cell)
1961 typename Triangulation<dim,spacedim>::line_iterator
1963 for (
unsigned int line=0; line<GeometryInfo<dim>::lines_per_cell; ++line)
1964 lines[line]=needed_lines[std::make_pair(
1971 lines[3]->index()));
1973 cell->set_used_flag ();
1974 cell->set_material_id (cells[c].material_id);
1975 cell->set_manifold_id (cells[c].manifold_id);
1976 cell->clear_user_data ();
1977 cell->set_subdomain_id (0);
1982 for (
unsigned int line=0; line<GeometryInfo<dim>::lines_per_cell; ++line)
1983 adjacent_cells[lines[line]->index()].push_back (cell);
1988 for (
typename Triangulation<dim,spacedim>::line_iterator
1989 line=triangulation.begin_line();
1990 line!=triangulation.end_line(); ++line)
1992 const unsigned int n_adj_cells = adjacent_cells[line->index()].size();
2005 ExcMessage (
"You have a line in your triangulation " 2006 "at which more than two cells come together. " 2008 "This is not currently supported because the " 2009 "Triangulation class makes the assumption that " 2010 "every cell has zero or one neighbors behind " 2011 "each face (here, behind each line), but in your " 2012 "situation there would be more than one." 2014 "Support for this is not currently implemented. " 2015 "If you need to work with triangulations where " 2016 "more than two cells come together at a line, " 2017 "duplicate the vertices once per cell (i.e., put " 2018 "multiple vertices at the same physical location, " 2019 "but using different vertex indices for each) " 2020 "and then ensure continuity of the solution by " 2021 "explicitly creating constraints that the degrees " 2022 "of freedom at these lines have the same " 2023 "value, using the ConstraintMatrix class."));
2032 std::vector<CellData<1> >::const_iterator boundary_line
2034 std::vector<CellData<1> >::const_iterator end_boundary_line
2036 for (; boundary_line!=end_boundary_line; ++boundary_line)
2038 typename Triangulation<dim,spacedim>::line_iterator line;
2039 std::pair<int,int> line_vertices(std::make_pair(boundary_line->vertices[0],
2040 boundary_line->vertices[1]));
2041 if (needed_lines.find(line_vertices) != needed_lines.end())
2043 line = needed_lines[line_vertices];
2047 std::swap (line_vertices.first, line_vertices.second);
2048 if (needed_lines.find(line_vertices) != needed_lines.end())
2049 line = needed_lines[line_vertices];
2053 line_vertices.second));
2060 line_vertices.second));
2084 line->vertex_index(1),
2085 boundary_line->boundary_id));
2089 line->set_manifold_id (boundary_line->manifold_id);
2093 line->set_boundary_id_internal(boundary_line->boundary_id);
2095 line->set_manifold_id (boundary_line->manifold_id);
2101 cell=triangulation.begin(); cell!=triangulation.end(); ++cell)
2102 for (
unsigned int side=0; side<4; ++side)
2103 if (adjacent_cells[cell->line(side)->index()][0] == cell)
2107 if (adjacent_cells[cell->line(side)->index()].size() == 2)
2110 cell->set_neighbor (side,
2111 adjacent_cells[cell->line(side)->index()][1]);
2117 cell->set_neighbor (side,
2118 adjacent_cells[cell->line(side)->index()][0]);
2167 template <
int spacedim>
2178 const unsigned int dim=3;
2181 triangulation.vertices = v;
2182 triangulation.vertices_used = std::vector<bool> (v.size(),
true);
2187 for (
unsigned int cell_no = 0; cell_no<cells.size(); ++cell_no)
2190 if (!triangulation.check_for_distorted_cells)
2192 const double cell_measure = GridTools::cell_measure<3>
2193 (triangulation.vertices, cells[cell_no].vertices);
2218 typename std::map<std::pair<int,int>,
2219 typename Triangulation<dim,spacedim>::line_iterator> needed_lines;
2220 for (
unsigned int cell=0; cell<cells.size(); ++cell)
2224 for (
unsigned int vertex=0; vertex<GeometryInfo<dim>::vertices_per_cell; ++vertex)
2225 AssertThrow (cells[cell].vertices[vertex] < triangulation.vertices.size(),
2227 triangulation.vertices.size()));
2229 for (
unsigned int line=0; line<GeometryInfo<dim>::lines_per_cell; ++line)
2237 std::pair<int,int> line_vertices(
2244 if ( (needed_lines.find(std::make_pair(line_vertices.second,
2245 line_vertices.first))
2247 needed_lines.end()))
2253 needed_lines[line_vertices] = triangulation.end_line();
2265 std::vector<unsigned short int> vertex_touch_count (v.size(), 0);
2266 typename std::map<std::pair<int,int>,
2267 typename Triangulation<dim,spacedim>::line_iterator>::iterator i;
2268 for (i=needed_lines.begin(); i!=needed_lines.end(); ++i)
2272 ++vertex_touch_count[i->first.first];
2273 ++vertex_touch_count[i->first.second];
2280 AssertThrow (* (std::min_element(vertex_touch_count.begin(),
2281 vertex_touch_count.end())) >= 3,
2282 ExcMessage(
"During creation of a triangulation, a part of the " 2283 "algorithm encountered a vertex that is part of only " 2284 "one or two adjacent lines. However, in 3d, every vertex " 2285 "needs to be at least part of three lines."));
2294 triangulation.faces = std_cxx14::make_unique<internal::TriangulationImplementation::TriaFaces<dim>> ();
2295 triangulation.levels[0]->reserve_space (cells.size(), dim, spacedim);
2296 triangulation.faces->lines.reserve_space (0,needed_lines.size());
2300 typename Triangulation<dim,spacedim>::raw_line_iterator
2301 line = triangulation.begin_raw_line();
2302 typename std::map<std::pair<int,int>,
2303 typename Triangulation<dim,spacedim>::line_iterator>::iterator i;
2304 for (i = needed_lines.begin(); line!=triangulation.end_line(); ++line, ++i)
2308 line->set_used_flag ();
2309 line->clear_user_flag ();
2333 std::map<internal::TriangulationImplementation::TriaObject<2>,
2334 std::pair<typename Triangulation<dim,spacedim>::quad_iterator,
2335 std::array<bool,GeometryInfo<dim>::lines_per_face> >,
2338 for (
unsigned int cell=0; cell<cells.size(); ++cell)
2366 std::array<bool,GeometryInfo<dim>::lines_per_face> orientation;
2368 for (
unsigned int line=0; line<GeometryInfo<dim>::lines_per_cell; ++line)
2370 line_list[line]=std::pair<int,int> (
2373 inverse_line_list[line]=std::pair<int,int> (
2378 for (
unsigned int face=0; face<GeometryInfo<dim>::faces_per_cell; ++face)
2388 for (
unsigned int l=0; l<GeometryInfo<dim>::lines_per_face; ++l)
2390 face_to_cell_lines(face,l)]) == needed_lines.end())
2394 orientation[l]=
true;
2400 orientation[l]=
false;
2405 quad(face_line_list[0],
2449 test_quad_1(quad.
face(2), quad.
face(3),
2451 test_quad_2(quad.
face(0), quad.
face(1),
2453 test_quad_3(quad.
face(3), quad.
face(2),
2455 test_quad_4(quad.
face(1), quad.
face(0),
2457 test_quad_5(quad.
face(2), quad.
face(3),
2459 test_quad_6(quad.
face(1), quad.
face(0),
2461 test_quad_7(quad.
face(3), quad.
face(2),
2463 if (needed_quads.find (test_quad_1) == needed_quads.end() &&
2464 needed_quads.find (test_quad_2) == needed_quads.end() &&
2465 needed_quads.find (test_quad_3) == needed_quads.end() &&
2466 needed_quads.find (test_quad_4) == needed_quads.end() &&
2467 needed_quads.find (test_quad_5) == needed_quads.end() &&
2468 needed_quads.find (test_quad_6) == needed_quads.end() &&
2469 needed_quads.find (test_quad_7) == needed_quads.end())
2470 needed_quads[quad] = std::make_pair(triangulation.end_quad(),orientation);
2480 triangulation.faces->quads.reserve_space (0,needed_quads.size());
2483 typename Triangulation<dim,spacedim>::raw_quad_iterator
2484 quad = triangulation.begin_raw_quad();
2485 typename std::map<internal::TriangulationImplementation::TriaObject<2>,
2486 std::pair<typename Triangulation<dim,spacedim>::quad_iterator,
2487 std::array<bool,GeometryInfo<dim>::lines_per_face> >,
2490 for (q = needed_quads.begin(); quad!=triangulation.end_quad(); ++quad, ++q)
2492 quad->set (q->first);
2493 quad->set_used_flag ();
2494 quad->clear_user_flag ();
2497 quad->set_line_orientation(0,q->second.second[0]);
2498 quad->set_line_orientation(1,q->second.second[1]);
2499 quad->set_line_orientation(2,q->second.second[2]);
2500 quad->set_line_orientation(3,q->second.second[3]);
2505 q->second.first = quad;
2511 triangulation.levels[0]->cells.reserve_space (cells.size());
2515 std::map<int,std::vector<typename Triangulation<dim,spacedim>::cell_iterator> >
2521 cell = triangulation.begin_raw_hex();
2522 for (
unsigned int c=0; c<cells.size(); ++c, ++cell)
2537 unsigned int face_line_list[4];
2538 for (
unsigned int line=0; line<GeometryInfo<dim>::lines_per_cell; ++line)
2540 line_list[line]=std::make_pair(
2543 inverse_line_list[line]=std::pair<int,int> (
2553 typename Triangulation<dim,spacedim>::quad_iterator
2558 for (
unsigned int face=0; face<GeometryInfo<dim>::faces_per_cell; ++face)
2560 for (
unsigned int l=0; l<GeometryInfo<dim>::lines_per_face; ++l)
2562 face_to_cell_lines(face,l)]) == needed_lines.end())
2564 face_to_cell_lines(face,l)]]->index();
2570 quad(face_line_list[0],
2575 if (needed_quads.find (quad) != needed_quads.end())
2586 face_iterator[face] = needed_quads[quad].first;
2587 face_orientation[face] =
true;
2588 face_flip[face]=
false;
2589 face_rotation[face]=
false;
2599 test_quad_1(quad.
face(2), quad.
face(3),
2601 test_quad_2(quad.
face(0), quad.
face(1),
2603 test_quad_3(quad.
face(3), quad.
face(2),
2605 test_quad_4(quad.
face(1), quad.
face(0),
2607 test_quad_5(quad.
face(2), quad.
face(3),
2609 test_quad_6(quad.
face(1), quad.
face(0),
2611 test_quad_7(quad.
face(3), quad.
face(2),
2613 if (needed_quads.find (test_quad_1) != needed_quads.end())
2615 face_iterator[face] = needed_quads[test_quad_1].first;
2616 face_orientation[face] =
false;
2617 face_flip[face]=
false;
2618 face_rotation[face]=
false;
2620 else if (needed_quads.find (test_quad_2) != needed_quads.end())
2622 face_iterator[face] = needed_quads[test_quad_2].first;
2623 face_orientation[face] =
false;
2624 face_flip[face]=
false;
2625 face_rotation[face]=
true;
2627 else if (needed_quads.find (test_quad_3) != needed_quads.end())
2629 face_iterator[face] = needed_quads[test_quad_3].first;
2630 face_orientation[face] =
false;
2631 face_flip[face]=
true;
2632 face_rotation[face]=
false;
2634 else if (needed_quads.find (test_quad_4) != needed_quads.end())
2636 face_iterator[face] = needed_quads[test_quad_4].first;
2637 face_orientation[face] =
false;
2638 face_flip[face]=
true;
2639 face_rotation[face]=
true;
2641 else if (needed_quads.find (test_quad_5) != needed_quads.end())
2643 face_iterator[face] = needed_quads[test_quad_5].first;
2644 face_orientation[face] =
true;
2645 face_flip[face]=
false;
2646 face_rotation[face]=
true;
2648 else if (needed_quads.find (test_quad_6) != needed_quads.end())
2650 face_iterator[face] = needed_quads[test_quad_6].first;
2651 face_orientation[face] =
true;
2652 face_flip[face]=
true;
2653 face_rotation[face]=
false;
2655 else if (needed_quads.find (test_quad_7) != needed_quads.end())
2657 face_iterator[face] = needed_quads[test_quad_7].first;
2658 face_orientation[face] =
true;
2659 face_flip[face]=
true;
2660 face_rotation[face]=
true;
2677 face_iterator[1]->index(),
2678 face_iterator[2]->index(),
2679 face_iterator[3]->index(),
2680 face_iterator[4]->index(),
2681 face_iterator[5]->index()));
2683 cell->set_used_flag ();
2684 cell->set_material_id (cells[c].material_id);
2685 cell->set_manifold_id (cells[c].manifold_id);
2686 cell->clear_user_flag ();
2687 cell->clear_user_data ();
2688 cell->set_subdomain_id (0);
2692 for (
unsigned int quad=0; quad<GeometryInfo<dim>::faces_per_cell; ++quad)
2694 cell->set_face_orientation (quad, face_orientation[quad]);
2695 cell->set_face_flip (quad, face_flip[quad]);
2696 cell->set_face_rotation (quad, face_rotation[quad]);
2703 for (
unsigned int quad=0; quad<GeometryInfo<dim>::faces_per_cell; ++quad)
2704 adjacent_cells[face_iterator[quad]->index()].push_back (cell);
2721 std::multimap<unsigned int, std::pair<unsigned int, unsigned int> >
2723 for (
unsigned int face=0; face<GeometryInfo<dim>::faces_per_cell; ++face)
2724 for (
unsigned int line=0; line<GeometryInfo<dim>::lines_per_face; ++line)
2725 cell_to_face_lines.insert(
2726 std::pair<
unsigned int, std::pair<unsigned int, unsigned int> > (
2728 std::pair<unsigned int, unsigned int> (face,line)));
2729 std::multimap<unsigned int, std::pair<unsigned int, unsigned int> >::const_iterator
2730 map_iter=cell_to_face_lines.begin();
2732 for (; map_iter!=cell_to_face_lines.end(); ++map_iter)
2734 const unsigned int cell_line=map_iter->first;
2735 const unsigned int face1=map_iter->second.first;
2736 const unsigned int line1=map_iter->second.second;
2740 const unsigned int face2=map_iter->second.first;
2741 const unsigned int line2=map_iter->second.second;
2750 face_orientation[face1],
2752 face_rotation[face1])) ==
2755 face_orientation[face2],
2757 face_rotation[face2])),
2768 for (
typename Triangulation<dim,spacedim>::quad_iterator
2769 quad=triangulation.begin_quad(); quad!=triangulation.end_quad(); ++quad)
2771 const unsigned int n_adj_cells = adjacent_cells[quad->index()].size();
2793 for (
typename Triangulation<dim,spacedim>::line_iterator
2794 line=triangulation.begin_line(); line!=triangulation.end_line(); ++line)
2818 for (
typename Triangulation<dim,spacedim>::quad_iterator
2819 quad=triangulation.begin_quad(); quad!=triangulation.end_quad(); ++quad)
2820 if (quad->at_boundary())
2822 for (
unsigned int l=0; l<4; ++l)
2824 typename Triangulation<dim,spacedim>::line_iterator line = quad->line(l);
2825 line->set_boundary_id_internal(0);
2834 std::vector<CellData<1> >::const_iterator boundary_line
2836 std::vector<CellData<1> >::const_iterator end_boundary_line
2838 for (; boundary_line!=end_boundary_line; ++boundary_line)
2840 typename Triangulation<dim,spacedim>::line_iterator line;
2841 std::pair <int, int> line_vertices(std::make_pair(boundary_line->vertices[0],
2842 boundary_line->vertices[1]));
2843 if (needed_lines.find(line_vertices) != needed_lines.end())
2846 line = needed_lines[line_vertices];
2852 std::swap (line_vertices.first, line_vertices.second);
2853 if (needed_lines.find(line_vertices) != needed_lines.end())
2854 line = needed_lines[line_vertices];
2858 line_vertices.second));
2865 line->vertex_index(1),
2866 boundary_line->boundary_id));
2873 if (line->boundary_id() != 0)
2874 AssertThrow (line->boundary_id() == boundary_line->boundary_id,
2875 ExcMessage (
"Duplicate boundary lines are only allowed " 2876 "if they carry the same boundary indicator."));
2878 line->set_boundary_id_internal(boundary_line->boundary_id);
2880 line->set_manifold_id(boundary_line->manifold_id);
2885 std::vector<CellData<2> >::const_iterator boundary_quad
2887 std::vector<CellData<2> >::const_iterator end_boundary_quad
2889 for (; boundary_quad!=end_boundary_quad; ++boundary_quad)
2891 typename Triangulation<dim,spacedim>::quad_iterator quad;
2892 typename Triangulation<dim,spacedim>::line_iterator line[4];
2901 for (
unsigned int i=0; i<4; ++i)
2903 std::pair<int, int> line_vertices(
2909 if (needed_lines.find(line_vertices) != needed_lines.end())
2910 line[i] = needed_lines[line_vertices];
2915 std::swap (line_vertices.first, line_vertices.second);
2916 if (needed_lines.find(line_vertices) != needed_lines.end())
2917 line[i] = needed_lines[line_vertices];
2922 line_vertices.second));
2940 quad_compare_1(line[0]->index(), line[1]->index(),
2941 line[2]->index(), line[3]->index());
2943 quad_compare_2(line[2]->index(), line[3]->index(),
2944 line[0]->index(), line[1]->index());
2957 static const unsigned int lex2cclock[4]= {3,1,0,2};
2963 typename Triangulation<dim,spacedim>::line_iterator
2964 line_counterclock[4];
2965 for (
unsigned int i=0; i<4; ++i)
2966 line_counterclock[lex2cclock[i]]=line[i];
2967 unsigned int n_rotations=0;
2968 bool not_found_quad_1;
2969 while ( (not_found_quad_1=(needed_quads.find(quad_compare_1) == needed_quads.end())) &&
2970 ( needed_quads.find(quad_compare_2) == needed_quads.end()) &&
2975 std::rotate(line_counterclock, line_counterclock+1, line_counterclock+4);
2979 for (
unsigned int i=0; i<4; ++i)
2981 quad_compare_1.
set_face(i, line_counterclock[lex2cclock[i]]->index());
2982 quad_compare_2.
set_face((i+2)%4, line_counterclock[lex2cclock[i]]->index());
2990 line[2]->index(), line[3]->index()));
2992 if (not_found_quad_1)
2993 quad = needed_quads[quad_compare_2].first;
2995 quad = needed_quads[quad_compare_1].first;
3001 quad->vertex_index(1),
3002 quad->vertex_index(2),
3003 quad->vertex_index(3),
3004 boundary_quad->boundary_id));
3011 if (quad->boundary_id() != 0)
3012 AssertThrow (quad->boundary_id() == boundary_quad->boundary_id,
3013 ExcMessage (
"Duplicate boundary quads are only allowed " 3014 "if they carry the same boundary indicator."));
3016 quad->set_boundary_id_internal (boundary_quad->boundary_id);
3017 quad->set_manifold_id (boundary_quad->manifold_id);
3024 cell=triangulation.begin(); cell!=triangulation.end(); ++cell)
3025 for (
unsigned int face=0; face<6; ++face)
3026 if (adjacent_cells[cell->quad(face)->index()][0] == cell)
3030 if (adjacent_cells[cell->quad(face)->index()].size() == 2)
3033 cell->set_neighbor (face,
3034 adjacent_cells[cell->quad(face)->index()][1]);
3040 cell->set_neighbor (face,
3041 adjacent_cells[cell->quad(face)->index()][0]);
3060 template <
int spacedim>
3065 std::vector<unsigned int> &,
3066 std::vector<unsigned int> &)
3068 const unsigned int dim = 1;
3080 Assert (!cell->child(0)->has_children() && !cell->child(1)->has_children(),
3086 if (cell->neighbor(0)->has_children())
3089 neighbor = cell->neighbor(0);
3093 neighbor = neighbor->child(1);
3096 Assert (neighbor->neighbor(1) == cell->child(0),
3098 neighbor->set_neighbor (1, cell);
3104 if (neighbor->has_children())
3105 neighbor = neighbor->child(1);
3114 if (cell->neighbor(1)->has_children())
3117 neighbor = cell->neighbor(1);
3121 neighbor = neighbor->child(0);
3124 Assert (neighbor->neighbor(0) == cell->child(1),
3126 neighbor->set_neighbor (0, cell);
3132 if (neighbor->has_children())
3133 neighbor = neighbor->child(0);
3143 triangulation.vertices_used[cell->child(0)->vertex_index(1)] =
false;
3149 for (
unsigned int child=0; child<cell->n_children(); ++child)
3152 cell->child(child)->clear_user_flag();
3153 cell->child(child)->clear_used_flag();
3158 cell->clear_children ();
3159 cell->clear_user_flag();
3164 template <
int spacedim>
3169 std::vector<unsigned int> &line_cell_count,
3170 std::vector<unsigned int> &)
3172 const unsigned int dim=2;
3179 std::vector<typename Triangulation<dim,spacedim>::line_iterator>
3182 lines_to_delete.reserve(4*2+4);
3187 for (
unsigned int c=0; c<cell->n_children(); ++c)
3190 child=cell->child(c);
3191 for (
unsigned int l=0; l<GeometryInfo<dim>::lines_per_cell; ++l)
3192 --line_cell_count[child->line_index(l)];
3206 triangulation.vertices_used[cell->child(0)->line(1)->vertex_index(1)] =
false;
3208 lines_to_delete.push_back(cell->child(0)->line(1));
3209 lines_to_delete.push_back(cell->child(0)->line(3));
3210 lines_to_delete.push_back(cell->child(3)->line(0));
3211 lines_to_delete.push_back(cell->child(3)->line(2));
3219 lines_to_delete.push_back(cell->child(0)->line(inner_face_no));
3223 for (
unsigned int child=0; child<cell->n_children(); ++child)
3225 cell->child(child)->clear_user_data();
3226 cell->child(child)->clear_user_flag();
3227 cell->child(child)->clear_used_flag();
3232 cell->clear_children ();
3233 cell->clear_refinement_case();
3234 cell->clear_user_flag();
3240 for (
unsigned int line_no=0; line_no<GeometryInfo<dim>::lines_per_cell; ++line_no)
3242 typename Triangulation<dim,spacedim>::line_iterator
3243 line=cell->line(line_no);
3245 if (line->has_children())
3250 Assert((line_cell_count[line->child_index(0)] == 0 &&
3251 line_cell_count[line->child_index(1)] == 0) ||
3252 (line_cell_count[line->child_index(0)] > 0 &&
3253 line_cell_count[line->child_index(1)] > 0),
3256 if (line_cell_count[line->child_index(0)]==0)
3258 for (
unsigned int c=0; c<2; ++c)
3259 Assert (!line->child(c)->has_children(),
3266 triangulation.vertices_used[line->child(0)->vertex_index(1)] =
false;
3268 lines_to_delete.push_back(line->child(0));
3269 lines_to_delete.push_back(line->child(1));
3271 line->clear_children();
3283 typename std::vector<typename Triangulation<dim,spacedim>::line_iterator>::iterator
3284 line=lines_to_delete.
begin(),
3285 endline=lines_to_delete.end();
3286 for (; line!=endline; ++line)
3288 (*line)->clear_user_data();
3289 (*line)->clear_user_flag();
3290 (*line)->clear_used_flag();
3296 template <
int spacedim>
3301 std::vector<unsigned int> &line_cell_count,
3302 std::vector<unsigned int> &quad_cell_count)
3304 const unsigned int dim=3;
3314 std::vector<typename Triangulation<dim,spacedim>::line_iterator>
3316 std::vector<typename Triangulation<dim,spacedim>::quad_iterator>
3319 lines_to_delete.reserve(12*2+6*4+6);
3320 quads_to_delete.reserve(6*4+12);
3324 for (
unsigned int c=0; c<cell->n_children(); ++c)
3327 child=cell->child(c);
3328 for (
unsigned int l=0; l<GeometryInfo<dim>::lines_per_cell; ++
l)
3329 --line_cell_count[child->line_index(l)];
3330 for (
unsigned int f=0; f<GeometryInfo<dim>::faces_per_cell; ++f)
3331 --quad_cell_count[child->quad_index(f)];
3345 quads_to_delete.push_back(cell->child(0)->face(1));
3348 quads_to_delete.push_back(cell->child(0)->face(3));
3351 quads_to_delete.push_back(cell->child(0)->face(5));
3354 quads_to_delete.push_back(cell->child(0)->face(1));
3355 quads_to_delete.push_back(cell->child(0)->face(3));
3356 quads_to_delete.push_back(cell->child(3)->face(0));
3357 quads_to_delete.push_back(cell->child(3)->face(2));
3359 lines_to_delete.push_back(cell->child(0)->line(11));
3362 quads_to_delete.push_back(cell->child(0)->face(1));
3363 quads_to_delete.push_back(cell->child(0)->face(5));
3364 quads_to_delete.push_back(cell->child(3)->face(0));
3365 quads_to_delete.push_back(cell->child(3)->face(4));
3367 lines_to_delete.push_back(cell->child(0)->line(5));
3370 quads_to_delete.push_back(cell->child(0)->face(3));
3371 quads_to_delete.push_back(cell->child(0)->face(5));
3372 quads_to_delete.push_back(cell->child(3)->face(2));
3373 quads_to_delete.push_back(cell->child(3)->face(4));
3375 lines_to_delete.push_back(cell->child(0)->line(7));
3378 quads_to_delete.push_back(cell->child(0)->face(1));
3379 quads_to_delete.push_back(cell->child(2)->face(1));
3380 quads_to_delete.push_back(cell->child(4)->face(1));
3381 quads_to_delete.push_back(cell->child(6)->face(1));
3383 quads_to_delete.push_back(cell->child(0)->face(3));
3384 quads_to_delete.push_back(cell->child(1)->face(3));
3385 quads_to_delete.push_back(cell->child(4)->face(3));
3386 quads_to_delete.push_back(cell->child(5)->face(3));
3388 quads_to_delete.push_back(cell->child(0)->face(5));
3389 quads_to_delete.push_back(cell->child(1)->face(5));
3390 quads_to_delete.push_back(cell->child(2)->face(5));
3391 quads_to_delete.push_back(cell->child(3)->face(5));
3393 lines_to_delete.push_back(cell->child(0)->line(5));
3394 lines_to_delete.push_back(cell->child(0)->line(7));
3395 lines_to_delete.push_back(cell->child(0)->line(11));
3396 lines_to_delete.push_back(cell->child(7)->line(0));
3397 lines_to_delete.push_back(cell->child(7)->line(2));
3398 lines_to_delete.push_back(cell->child(7)->line(8));
3404 triangulation.vertices_used[cell->child(0)->vertex_index(7)] =
false;
3415 for (
unsigned int child=0; child<cell->n_children(); ++child)
3417 cell->child(child)->clear_user_data();
3418 cell->child(child)->clear_user_flag();
3420 for (
unsigned int f=0; f<GeometryInfo<dim>::faces_per_cell; ++f)
3425 cell->child(child)->set_face_orientation (f,
true);
3426 cell->child(child)->set_face_flip(f,
false);
3427 cell->child(child)->set_face_rotation(f,
false);
3430 cell->child(child)->clear_used_flag();
3435 cell->clear_children ();
3436 cell->clear_refinement_case ();
3437 cell->clear_user_flag();
3447 for (
unsigned int quad_no=0; quad_no<GeometryInfo<dim>::faces_per_cell; ++quad_no)
3449 typename Triangulation<dim,spacedim>::quad_iterator
3450 quad=cell->face(quad_no);
3456 switch (quad->refinement_case())
3468 Assert((quad_cell_count[quad->child_index(0)] == 0 &&
3469 quad_cell_count[quad->child_index(1)] == 0) ||
3470 (quad_cell_count[quad->child_index(0)] > 0 &&
3471 quad_cell_count[quad->child_index(1)] > 0),
3477 unsigned int deleted_grandchildren=0;
3478 unsigned int number_of_child_refinements=0;
3480 for (
unsigned int c=0; c<2; ++c)
3481 if (quad->child(c)->has_children())
3483 ++number_of_child_refinements;
3487 Assert((quad_cell_count[quad->child(c)->child_index(0)] == 0 &&
3488 quad_cell_count[quad->child(c)->child_index(1)] == 0) ||
3489 (quad_cell_count[quad->child(c)->child_index(0)] > 0 &&
3490 quad_cell_count[quad->child(c)->child_index(1)] > 0),
3492 if (quad_cell_count[quad->child(c)->child_index(0)]==0)
3505 quads_to_delete.push_back(quad->child(c)->child(0));
3506 quads_to_delete.push_back(quad->child(c)->child(1));
3508 lines_to_delete.push_back(quad->child(c)->child(0)->line(1));
3510 lines_to_delete.push_back(quad->child(c)->child(0)->line(3));
3511 quad->child(c)->clear_children();
3512 quad->child(c)->clear_refinement_case();
3513 ++deleted_grandchildren;
3521 if (number_of_child_refinements>0 &&
3522 deleted_grandchildren==number_of_child_refinements)
3524 typename Triangulation<dim,spacedim>::line_iterator
3527 middle_line=quad->child(0)->line(1);
3529 middle_line=quad->child(0)->line(3);
3531 lines_to_delete.push_back(middle_line->child(0));
3532 lines_to_delete.push_back(middle_line->child(1));
3533 triangulation.vertices_used[middle_vertex_index<dim,spacedim>(middle_line)]
3535 middle_line->clear_children();
3540 if (quad_cell_count[quad->child_index(0)]==0)
3546 quads_to_delete.push_back(quad->child(0));
3547 quads_to_delete.push_back(quad->child(1));
3549 lines_to_delete.push_back(quad->child(0)->line(1));
3551 lines_to_delete.push_back(quad->child(0)->line(3));
3571 if (quad->child(0)->has_children())
3604 const typename Triangulation<dim,spacedim>::quad_iterator
3605 switch_1=quad->child(0)->child(1),
3606 switch_2=quad->child(1)->child(0);
3611 const int switch_1_index=switch_1->index();
3612 const int switch_2_index=switch_2->index();
3613 for (
unsigned int l=0;
l<triangulation.levels.size(); ++
l)
3614 for (
unsigned int h=0; h<triangulation.levels[
l]->cells.cells.size(); ++h)
3615 for (
unsigned int q=0; q<GeometryInfo<dim>::faces_per_cell; ++q)
3617 const int index=triangulation.levels[
l]->cells.cells[h].face(q);
3618 if (index==switch_1_index)
3619 triangulation.levels[
l]->cells.cells[h].set_face(q,switch_2_index);
3620 else if (index==switch_2_index)
3621 triangulation.levels[
l]->cells.cells[h].set_face(q,switch_1_index);
3626 const int switch_1_lines[4]=
3628 static_cast<signed int>(switch_1->line_index(0)),
3629 static_cast<signed int>(switch_1->line_index(1)),
3630 static_cast<signed int>(switch_1->line_index(2)),
3631 static_cast<signed int>(switch_1->line_index(3))
3633 const bool switch_1_line_orientations[4]=
3635 switch_1->line_orientation(0),
3636 switch_1->line_orientation(1),
3637 switch_1->line_orientation(2),
3638 switch_1->line_orientation(3)
3641 const unsigned int switch_1_user_index=switch_1->user_index();
3642 const bool switch_1_user_flag=switch_1->user_flag_set();
3645 switch_2->line_index(1),
3646 switch_2->line_index(2),
3647 switch_2->line_index(3)));
3648 switch_1->set_line_orientation(0, switch_2->line_orientation(0));
3649 switch_1->set_line_orientation(1, switch_2->line_orientation(1));
3650 switch_1->set_line_orientation(2, switch_2->line_orientation(2));
3651 switch_1->set_line_orientation(3, switch_2->line_orientation(3));
3652 switch_1->set_boundary_id_internal(switch_2->boundary_id());
3653 switch_1->set_manifold_id(switch_2->manifold_id());
3654 switch_1->set_user_index(switch_2->user_index());
3655 if (switch_2->user_flag_set())
3656 switch_1->set_user_flag();
3658 switch_1->clear_user_flag();
3663 switch_1_lines[3]));
3664 switch_2->set_line_orientation(0, switch_1_line_orientations[0]);
3665 switch_2->set_line_orientation(1, switch_1_line_orientations[1]);
3666 switch_2->set_line_orientation(2, switch_1_line_orientations[2]);
3667 switch_2->set_line_orientation(3, switch_1_line_orientations[3]);
3668 switch_2->set_boundary_id_internal(switch_1_boundary_id);
3669 switch_2->set_manifold_id(switch_1->manifold_id());
3670 switch_2->set_user_index(switch_1_user_index);
3671 if (switch_1_user_flag)
3672 switch_2->set_user_flag();
3674 switch_2->clear_user_flag();
3676 const unsigned int child_0=quad->child(0)->child_index(0);
3677 const unsigned int child_2=quad->child(1)->child_index(0);
3678 quad->clear_children();
3679 quad->clear_refinement_case();
3681 quad->set_children(0,child_0);
3682 quad->set_children(2,child_2);
3683 std::swap(quad_cell_count[child_0+1],quad_cell_count[child_2]);
3698 const unsigned int child_0=quad->child(0)->child_index(0);
3699 const unsigned int child_2=quad->child(1)->child_index(0);
3700 quad->clear_children();
3701 quad->clear_refinement_case();
3703 quad->set_children(0,child_0);
3704 quad->set_children(2,child_2);
3708 quad->clear_children();
3709 quad->clear_refinement_case();
3722 Assert((quad_cell_count[quad->child_index(0)] == 0 &&
3723 quad_cell_count[quad->child_index(1)] == 0 &&
3724 quad_cell_count[quad->child_index(2)] == 0 &&
3725 quad_cell_count[quad->child_index(3)] == 0) ||
3726 (quad_cell_count[quad->child_index(0)] > 0 &&
3727 quad_cell_count[quad->child_index(1)] > 0 &&
3728 quad_cell_count[quad->child_index(2)] > 0 &&
3729 quad_cell_count[quad->child_index(3)] > 0),
3732 if (quad_cell_count[quad->child_index(0)]==0)
3738 lines_to_delete.push_back(quad->child(0)->line(1));
3739 lines_to_delete.push_back(quad->child(3)->line(0));
3740 lines_to_delete.push_back(quad->child(0)->line(3));
3741 lines_to_delete.push_back(quad->child(3)->line(2));
3743 for (
unsigned int child=0; child<quad->n_children(); ++child)
3744 quads_to_delete.push_back(quad->child(child));
3746 triangulation.vertices_used[quad->child(0)->vertex_index(3)] =
false;
3748 quad->clear_children();
3749 quad->clear_refinement_case();
3769 for (
unsigned int line_no=0; line_no<GeometryInfo<dim>::lines_per_cell; ++line_no)
3771 typename Triangulation<dim,spacedim>::line_iterator
3772 line=cell->line(line_no);
3778 if (line->has_children())
3783 Assert((line_cell_count[line->child_index(0)] == 0 &&
3784 line_cell_count[line->child_index(1)] == 0) ||
3785 (line_cell_count[line->child_index(0)] > 0 &&
3786 line_cell_count[line->child_index(1)] > 0),
3789 if (line_cell_count[line->child_index(0)]==0)
3791 for (
unsigned int c=0; c<2; ++c)
3792 Assert (!line->child(c)->has_children(),
3799 triangulation.vertices_used[line->child(0)->vertex_index(1)] =
false;
3801 lines_to_delete.push_back(line->child(0));
3802 lines_to_delete.push_back(line->child(1));
3804 line->clear_children();
3816 typename std::vector<typename Triangulation<dim,spacedim>::line_iterator>::iterator
3817 line=lines_to_delete.
begin(),
3818 endline=lines_to_delete.end();
3819 for (; line!=endline; ++line)
3821 (*line)->clear_user_data();
3822 (*line)->clear_user_flag();
3823 (*line)->clear_used_flag();
3826 typename std::vector<typename Triangulation<dim,spacedim>::quad_iterator>::iterator
3827 quad=quads_to_delete.
begin(),
3828 endquad=quads_to_delete.end();
3829 for (; quad!=endquad; ++quad)
3831 (*quad)->clear_user_data();
3832 (*quad)->clear_children();
3833 (*quad)->clear_refinement_case();
3834 (*quad)->clear_user_flag();
3835 (*quad)->clear_used_flag();
3857 template <
int spacedim>
3861 unsigned int &next_unused_vertex,
3862 typename Triangulation<2,spacedim>::raw_line_iterator &next_unused_line,
3866 const unsigned int dim=2;
3869 cell->clear_refine_flag ();
3965 int new_vertices[9];
3966 for (
unsigned int vertex_no=0; vertex_no<4; ++vertex_no)
3967 new_vertices[vertex_no]=cell->vertex_index(vertex_no);
3968 for (
unsigned int line_no=0; line_no<4; ++line_no)
3969 if (cell->line(line_no)->has_children())
3970 new_vertices[4+line_no]=cell->line(line_no)->child(0)->vertex_index(1);
3980 while (triangulation.vertices_used[next_unused_vertex] ==
true)
3981 ++next_unused_vertex;
3982 Assert (next_unused_vertex < triangulation.vertices.size(),
3983 ExcMessage(
"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."));
3984 triangulation.vertices_used[next_unused_vertex] =
true;
3986 new_vertices[8] = next_unused_vertex;
3998 if (dim == spacedim)
4001 triangulation.vertices[next_unused_vertex] = cell->center(
true);
4008 if (cell->user_flag_set())
4011 cell->clear_user_flag();
4012 triangulation.vertices[next_unused_vertex]
4013 =cell->center(
true,
true);
4023 cell->clear_user_flag();
4027 triangulation.vertices[next_unused_vertex] =
4034 typename Triangulation<dim,spacedim>::raw_line_iterator new_lines[12];
4035 unsigned int lmin=8;
4036 unsigned int lmax=12;
4043 for (
unsigned int l=lmin; l<lmax; ++l)
4045 while (next_unused_line->used() ==
true)
4047 new_lines[l] = next_unused_line;
4050 Assert (new_lines[l]->used() ==
false,
4051 ExcMessage(
"Internal error: We want to use a cell during refinement that should be unused, but turns out not to be."));
4065 for (
unsigned int face_no=0; face_no<GeometryInfo<dim>::faces_per_cell; ++face_no)
4066 for (
unsigned int c=0; c<2; ++c, ++l)
4067 new_lines[l]=cell->line(face_no)->child(c);
4086 new_lines[0]=cell->line(0);
4087 new_lines[1]=cell->line(1);
4088 new_lines[2]=cell->line(2)->child(0);
4089 new_lines[3]=cell->line(2)->child(1);
4090 new_lines[4]=cell->line(3)->child(0);
4091 new_lines[5]=cell->line(3)->child(1);
4103 new_lines[0]=cell->line(0)->child(0);
4104 new_lines[1]=cell->line(0)->child(1);
4105 new_lines[2]=cell->line(1)->child(0);
4106 new_lines[3]=cell->line(1)->child(1);
4107 new_lines[4]=cell->line(2);
4108 new_lines[5]=cell->line(3);
4113 for (
unsigned int l=lmin; l<lmax; ++l)
4115 new_lines[l]->set_used_flag();
4116 new_lines[l]->clear_user_flag();
4118 new_lines[l]->clear_children();
4121 new_lines[l]->set_manifold_id(cell->manifold_id());
4128 while (next_unused_cell->used() ==
true)
4131 const unsigned int n_children=
4133 for (
unsigned int i=0; i<n_children; ++i)
4135 Assert (next_unused_cell->used() ==
false,
4136 ExcMessage(
"Internal error: We want to use a cell during refinement that should be unused, but turns out not to be."));
4137 subcells[i] = next_unused_cell;
4139 if (i%2==1 && i<n_children-1)
4140 while (next_unused_cell->used() ==
true)
4160 new_lines[8]->index(),
4161 new_lines[4]->index(),
4162 new_lines[10]->index()));
4165 new_lines[2]->index(),
4166 new_lines[5]->index(),
4167 new_lines[11]->index()));
4170 new_lines[9]->index(),
4171 new_lines[10]->index(),
4172 new_lines[6]->index()));
4175 new_lines[3]->index(),
4176 new_lines[11]->index(),
4177 new_lines[7]->index()));
4195 new_lines[6]->index(),
4196 new_lines[2]->index(),
4197 new_lines[4]->index()));
4200 new_lines[1]->index(),
4201 new_lines[3]->index(),
4202 new_lines[5]->index()));
4221 new_lines[2]->index(),
4222 new_lines[4]->index(),
4223 new_lines[6]->index()));
4226 new_lines[3]->index(),
4227 new_lines[6]->index(),
4228 new_lines[5]->index()));
4233 for (
unsigned int i=0; i<n_children; ++i)
4235 subcells[i]->set_used_flag();
4236 subcells[i]->clear_refine_flag();
4237 subcells[i]->clear_user_flag();
4238 subcells[i]->clear_user_data();
4239 subcells[i]->clear_children();
4242 subcells[i]->set_material_id (cell->material_id());
4243 subcells[i]->set_manifold_id (cell->manifold_id());
4244 subcells[i]->set_subdomain_id (subdomainid);
4247 subcells[i]->set_parent (cell->index ());
4255 for (
unsigned int i=0; i<n_children/2; ++i)
4256 cell->set_children (2*i, subcells[2*i]->index());
4258 cell->set_refinement_case(ref_case);
4266 for (
unsigned int c=0; c<n_children; ++c)
4267 cell->child(c)->set_direction_flag (cell->direction_flag());
4277 template <
int spacedim>
4283 const unsigned int dim = 1;
4290 cell = triangulation.begin_active (triangulation.levels.size()-1),
4291 endc = triangulation.end();
4292 for (; cell != endc; ++cell)
4294 if (cell->refine_flag_set())
4296 triangulation.levels
4307 unsigned int needed_vertices = 0;
4308 for (
int level=triangulation.levels.size()-2; level>=0; --level)
4312 unsigned int flagged_cells = 0;
4314 acell = triangulation.begin_active(level),
4315 aendc = triangulation.begin_active(level+1);
4316 for (; acell!=aendc; ++acell)
4317 if (acell->refine_flag_set())
4322 const unsigned int used_cells
4323 = std::count_if (triangulation.levels[level+1]->cells.used.begin(),
4324 triangulation.levels[level+1]->cells.used.end(),
4325 std::bind (std::equal_to<bool>(), std::placeholders::_1,
true));
4330 triangulation.levels[level+1]
4331 ->reserve_space(used_cells+
4338 triangulation.levels[level+1]->cells
4343 needed_vertices += flagged_cells;
4348 needed_vertices += std::count_if (triangulation.vertices_used.begin(),
4349 triangulation.vertices_used.end(),
4350 std::bind (std::equal_to<bool>(),
4351 std::placeholders::_1,
4356 if (needed_vertices > triangulation.vertices.size())
4358 triangulation.vertices.resize (needed_vertices,
4360 triangulation.vertices_used.resize (needed_vertices,
false);
4368 unsigned int next_unused_vertex = 0;
4370 for (
int level=triangulation.levels.size()-2; level>=0; --level)
4373 cell = triangulation.begin_active(level),
4374 endc = triangulation.begin_active(level+1);
4377 next_unused_cell = triangulation.begin_raw (level+1);
4379 for (; (cell!=endc) && (cell->level()==level); ++cell)
4380 if (cell->refine_flag_set())
4383 cell->clear_refine_flag ();
4387 while (triangulation.vertices_used[next_unused_vertex] ==
true)
4388 ++next_unused_vertex;
4389 Assert (next_unused_vertex < triangulation.vertices.size(),
4390 ExcMessage(
"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."));
4395 triangulation.vertices[next_unused_vertex] =
4398 triangulation.vertices_used[next_unused_vertex] =
true;
4405 while (next_unused_cell->used() ==
true)
4407 first_child = next_unused_cell;
4408 first_child->set_used_flag ();
4409 first_child->clear_user_data ();
4411 Assert (next_unused_cell->used() ==
false,
4412 ExcMessage(
"Internal error: We want to use a cell during refinement that should be unused, but turns out not to be."));
4413 second_child = next_unused_cell;
4414 second_child->set_used_flag ();
4415 second_child->clear_user_data ();
4420 cell->set_children (0, first_child->index());
4421 first_child->clear_children ();
4424 next_unused_vertex));
4425 first_child->set_material_id (cell->material_id());
4426 first_child->set_manifold_id (cell->manifold_id());
4427 first_child->set_subdomain_id (subdomainid);
4428 first_child->set_direction_flag (cell->direction_flag());
4430 first_child->set_parent (cell->index ());
4434 first_child->face(1)->set_manifold_id(cell->manifold_id());
4439 first_child->set_neighbor (1, second_child);
4441 first_child->set_neighbor (0, cell->neighbor(0));
4442 else if (cell->neighbor(0)->active())
4448 Assert (cell->neighbor (0)->level () <= cell->level (),
4450 first_child->set_neighbor (0, cell->neighbor(0));
4456 const unsigned int nbnb = cell->neighbor_of_neighbor (0);
4457 first_child->set_neighbor (0, cell->neighbor(0)->child(nbnb));
4462 left_neighbor = cell->neighbor(0);
4463 while (left_neighbor->has_children())
4465 left_neighbor = left_neighbor->child(nbnb);
4466 left_neighbor->set_neighbor (nbnb, first_child);
4471 second_child->clear_children ();
4474 cell->vertex_index(1)));
4475 second_child->set_neighbor (0, first_child);
4476 second_child->set_material_id (cell->material_id());
4477 second_child->set_manifold_id (cell->manifold_id());
4478 second_child->set_subdomain_id (subdomainid);
4479 second_child->set_direction_flag (cell->direction_flag());
4482 second_child->set_neighbor (1, cell->neighbor(1));
4483 else if (cell->neighbor(1)->active())
4485 Assert (cell->neighbor (1)->level () <= cell->level (),
4487 second_child->set_neighbor (1, cell->neighbor(1));
4492 const unsigned int nbnb = cell->neighbor_of_neighbor (1);
4493 second_child->set_neighbor (1, cell->neighbor(1)->child(nbnb));
4496 right_neighbor = cell->neighbor(1);
4497 while (right_neighbor->has_children())
4499 right_neighbor = right_neighbor->child(nbnb);
4500 right_neighbor->set_neighbor (nbnb, second_child);
4504 triangulation.signals.post_refinement_on_cell(cell);
4520 template <
int spacedim>
4524 const bool check_for_distorted_cells)
4526 const unsigned int dim = 2;
4534 cell = triangulation.begin_active (triangulation.levels.size()-1),
4535 endc = triangulation.end();
4536 for (; cell != endc; ++cell)
4538 if (cell->refine_flag_set())
4548 for (
typename Triangulation<dim,spacedim>::line_iterator
4549 line=triangulation.begin_line(); line!=triangulation.end_line(); ++line)
4551 line->clear_user_flag();
4552 line->clear_user_data();
4557 unsigned int n_single_lines=0;
4561 unsigned int n_lines_in_pairs = 0;
4567 unsigned int needed_vertices = 0;
4568 for (
int level=triangulation.levels.size()-2; level>=0; --level)
4572 unsigned int needed_cells = 0;
4575 cell = triangulation.begin_active(level),
4576 endc = triangulation.begin_active(level+1);
4577 for (; cell!=endc; ++cell)
4578 if (cell->refine_flag_set())
4589 n_single_lines += 4;
4595 triangulation.anisotropic_refinement =
true;
4601 n_single_lines += 1;
4610 for (
unsigned int line_no=0; line_no<GeometryInfo<dim>::faces_per_cell;
4616 typename Triangulation<dim,spacedim>::line_iterator
4617 line = cell->line(line_no);
4618 if (line->has_children() ==
false)
4620 line->set_user_flag ();
4628 if (line->at_boundary())
4631 line->set_user_index(line->boundary_id());
4636 line->set_user_index(cell->material_id());
4645 const unsigned int used_cells
4646 = std::count_if (triangulation.levels[level+1]->cells.used.begin(),
4647 triangulation.levels[level+1]->cells.used.end(),
4648 std::bind (std::equal_to<bool>(), std::placeholders::_1,
true));
4654 triangulation.levels[level+1]
4655 ->reserve_space (used_cells+needed_cells, 2, spacedim);
4659 triangulation.levels[level+1]->cells.
4660 reserve_space (needed_cells,0);
4664 for (
typename Triangulation<dim,spacedim>::line_iterator
4665 line=triangulation.begin_line(); line!=triangulation.end_line(); ++line)
4666 if (line->user_flag_set())
4669 n_lines_in_pairs += 2;
4670 needed_vertices += 1;
4679 triangulation.faces->lines.
4680 reserve_space (n_lines_in_pairs, 0);
4683 needed_vertices += std::count_if (triangulation.vertices_used.begin(), triangulation.vertices_used.end(),
4684 std::bind (std::equal_to<bool>(), std::placeholders::_1,
true));
4688 if (needed_vertices > triangulation.vertices.size())
4691 triangulation.vertices_used.resize (needed_vertices,
false);
4699 unsigned int next_unused_vertex = 0;
4706 typename Triangulation<dim,spacedim>::active_line_iterator
4707 line = triangulation.begin_active_line(),
4708 endl = triangulation.end_line();
4709 typename Triangulation<dim,spacedim>::raw_line_iterator
4710 next_unused_line = triangulation.begin_raw_line ();
4712 for (; line!=endl; ++line)
4713 if (line->user_flag_set())
4719 while (triangulation.vertices_used[next_unused_vertex] ==
true)
4720 ++next_unused_vertex;
4721 Assert (next_unused_vertex < triangulation.vertices.size(),
4722 ExcMessage(
"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."));
4723 triangulation.vertices_used[next_unused_vertex] =
true;
4725 if (spacedim == dim)
4732 triangulation.vertices[next_unused_vertex]
4733 = line->center(
true);
4742 triangulation.vertices[next_unused_vertex]
4743 = triangulation.get_manifold(line->user_index()).get_new_point_on_line (line);
4745 triangulation.vertices[next_unused_vertex]
4746 = line->center(
true);
4751 bool pair_found=
false;
4753 for (; next_unused_line!=endl; ++next_unused_line)
4754 if (!next_unused_line->used() &&
4755 !(++next_unused_line)->used())
4768 line->set_children (0, next_unused_line->index());
4771 const typename Triangulation<dim,spacedim>::raw_line_iterator
4772 children[2] = { next_unused_line,
4777 Assert (children[0]->used() ==
false,
ExcMessage(
"Internal error: We want to use a cell during refinement that should be unused, but turns out not to be."));
4778 Assert (children[1]->used() ==
false,
ExcMessage(
"Internal error: We want to use a cell during refinement that should be unused, but turns out not to be."));
4782 next_unused_vertex));
4785 line->vertex_index(1)));
4787 children[0]->set_used_flag();
4788 children[1]->set_used_flag();
4789 children[0]->clear_children();
4790 children[1]->clear_children();
4793 children[0]->clear_user_flag();
4794 children[1]->clear_user_flag();
4797 children[0]->set_boundary_id_internal(line->boundary_id());
4798 children[1]->set_boundary_id_internal(line->boundary_id());
4800 children[0]->set_manifold_id (line->manifold_id());
4801 children[1]->set_manifold_id (line->manifold_id());
4805 line->clear_user_flag ();
4814 triangulation.faces->lines.
4815 reserve_space (0,n_single_lines);
4818 cells_with_distorted_children;
4822 typename Triangulation<dim,spacedim>::raw_line_iterator
4823 next_unused_line = triangulation.begin_raw_line ();
4825 for (
int level=0; level<static_cast<int>(triangulation.levels.size())-1; ++level)
4831 cell = triangulation.begin_active(level),
4832 endc = triangulation.begin_active(level+1);
4835 next_unused_cell = triangulation.begin_raw (level+1);
4837 for (; cell!=endc; ++cell)
4838 if (cell->refine_flag_set())
4846 if (cell->at_boundary())
4847 cell->set_user_flag();
4857 if ((check_for_distorted_cells ==
true)
4859 has_distorted_children (cell,
4860 std::integral_constant<int, dim>(),
4861 std::integral_constant<int, spacedim>()))
4864 triangulation.signals.post_refinement_on_cell(cell);
4868 return cells_with_distorted_children;
4876 template <
int spacedim>
4880 const bool check_for_distorted_cells)
4882 const unsigned int dim = 3;
4896 cell = triangulation.begin_active (triangulation.levels.size()-1),
4897 endc = triangulation.end();
4898 for (; cell != endc; ++cell)
4900 if (cell->refine_flag_set())
4910 triangulation.faces->quads.clear_user_data();
4912 for (
typename Triangulation<dim,spacedim>::line_iterator
4913 line=triangulation.begin_line(); line!=triangulation.end_line(); ++line)
4914 line->clear_user_flag();
4915 for (
typename Triangulation<dim,spacedim>::quad_iterator
4916 quad=triangulation.begin_quad(); quad!=triangulation.end_quad(); ++quad)
4917 quad->clear_user_flag();
4942 unsigned int needed_vertices = 0;
4943 unsigned int needed_lines_single = 0;
4944 unsigned int needed_quads_single = 0;
4945 unsigned int needed_lines_pair = 0;
4946 unsigned int needed_quads_pair = 0;
4947 for (
int level=triangulation.levels.size()-2; level>=0; --level)
4951 unsigned int new_cells = 0;
4954 acell = triangulation.begin_active(level),
4955 aendc = triangulation.begin_active(level+1);
4956 for (; acell!=aendc; ++acell)
4957 if (acell->refine_flag_set())
4967 ++needed_quads_single;
4969 triangulation.anisotropic_refinement=
true;
4975 ++needed_lines_single;
4976 needed_quads_single += 4;
4978 triangulation.anisotropic_refinement=
true;
4983 needed_lines_single += 6;
4984 needed_quads_single += 12;
4998 for (
unsigned int face=0; face<GeometryInfo<dim>::faces_per_cell;
5002 aface = acell->face(face);
5008 acell->face_orientation(face),
5009 acell->face_flip(face),
5010 acell->face_rotation(face));
5017 if (aface->number_of_children()<4)
5020 aface->set_user_flag();
5022 else if (aface->refinement_case()!=face_ref_case)
5036 aface->set_user_index(face_ref_case);
5042 for (
unsigned int line=0; line<GeometryInfo<dim>::lines_per_cell; ++line)
5044 !acell->line(line)->has_children())
5045 acell->line(line)->set_user_flag();
5051 const unsigned int used_cells
5052 = std::count_if (triangulation.levels[level+1]->cells.used.begin(),
5053 triangulation.levels[level+1]->cells.used.end(),
5054 std::bind (std::equal_to<bool>(), std::placeholders::_1,
true));
5060 triangulation.levels[level+1]
5061 ->reserve_space (used_cells+new_cells, 3, spacedim);
5064 triangulation.levels[level+1]->cells.reserve_space (new_cells);
5068 for (
typename Triangulation<dim,spacedim>::quad_iterator
5069 quad=triangulation.begin_quad(); quad!=triangulation.end_quad(); ++quad)
5071 if (quad->user_flag_set())
5077 needed_quads_pair += 4;
5078 needed_lines_pair += 4;
5079 needed_vertices += 1;
5081 if (quad->user_index())
5085 needed_quads_pair += 2;
5086 needed_lines_single += 1;
5099 if (quad->has_children())
5102 if ((face_refinement_cases[quad->user_index()]==
RefinementCase<dim-1>::cut_x
5103 && (quad->child(0)->line_index(1)+1!=quad->child(2)->line_index(1))) ||
5104 (face_refinement_cases[quad->user_index()]==
RefinementCase<dim-1>::cut_y
5105 && (quad->child(0)->line_index(3)+1!=quad->child(1)->line_index(3))))
5106 needed_lines_pair +=2;
5111 for (
typename Triangulation<dim,spacedim>::line_iterator
5112 line=triangulation.begin_line(); line!=triangulation.end_line(); ++line)
5113 if (line->user_flag_set())
5115 needed_lines_pair += 2;
5116 needed_vertices += 1;
5120 triangulation.faces->lines.
5121 reserve_space (needed_lines_pair,needed_lines_single);
5123 triangulation.faces->quads.
5124 reserve_space (needed_quads_pair,needed_quads_single);
5128 needed_vertices += std::count_if (triangulation.vertices_used.begin(), triangulation.vertices_used.end(),
5129 std::bind (std::equal_to<bool>(), std::placeholders::_1,
true));
5133 if (needed_vertices > triangulation.vertices.size())
5136 triangulation.vertices_used.resize (needed_vertices,
false);
5152 cell=triangulation.begin_active(); cell!=triangulation.end(); ++cell)
5153 if (!cell->refine_flag_set())
5154 for (
unsigned int line=0; line<GeometryInfo<dim>::lines_per_cell; ++line)
5155 if (cell->line(line)->has_children())
5156 for (
unsigned int c=0; c<2; ++c)
5157 Assert (cell->line(line)->child(c)->user_flag_set() ==
false,
5169 unsigned int next_unused_vertex = 0;
5175 typename Triangulation<dim,spacedim>::active_line_iterator
5176 line = triangulation.begin_active_line(),
5177 endl = triangulation.end_line();
5178 typename Triangulation<dim,spacedim>::raw_line_iterator
5179 next_unused_line = triangulation.begin_raw_line ();
5181 for (; line!=endl; ++line)
5182 if (line->user_flag_set())
5188 while (triangulation.vertices_used[next_unused_vertex] ==
true)
5189 ++next_unused_vertex;
5190 Assert (next_unused_vertex < triangulation.vertices.size(),
5191 ExcMessage(
"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."));
5192 triangulation.vertices_used[next_unused_vertex] =
true;
5194 triangulation.vertices[next_unused_vertex]
5195 = line->center(
true);
5200 next_unused_line=triangulation.faces->lines.next_free_pair_object(triangulation);
5207 line->set_children (0, next_unused_line->index());
5210 const typename Triangulation<dim,spacedim>::raw_line_iterator
5211 children[2] = { next_unused_line,
5217 Assert (children[0]->used() ==
false,
ExcMessage(
"Internal error: We want to use a cell during refinement that should be unused, but turns out not to be."));
5218 Assert (children[1]->used() ==
false,
ExcMessage(
"Internal error: We want to use a cell during refinement that should be unused, but turns out not to be."));
5222 next_unused_vertex));
5225 line->vertex_index(1)));
5227 children[0]->set_used_flag();
5228 children[1]->set_used_flag();
5229 children[0]->clear_children();
5230 children[1]->clear_children();
5233 children[0]->clear_user_flag();
5234 children[1]->clear_user_flag();
5236 children[0]->set_boundary_id_internal(line->boundary_id());
5237 children[1]->set_boundary_id_internal(line->boundary_id());
5239 children[0]->set_manifold_id (line->manifold_id());
5240 children[1]->set_manifold_id (line->manifold_id());
5245 line->clear_user_flag ();
5280 typename Triangulation<dim,spacedim>::quad_iterator
5281 quad = triangulation.begin_quad(),
5282 endq = triangulation.end_quad();
5283 typename Triangulation<dim,spacedim>::raw_line_iterator
5284 next_unused_line = triangulation.begin_raw_line ();
5285 typename Triangulation<dim,spacedim>::raw_quad_iterator
5286 next_unused_quad = triangulation.begin_raw_quad ();
5288 for (; quad!=endq; ++quad)
5290 if (quad->user_index())
5292 RefinementCase<dim-1> aniso_quad_ref_case=face_refinement_cases[quad->user_index()];
5301 if (aniso_quad_ref_case == quad->refinement_case())
5314 typename Triangulation<dim,spacedim>::raw_line_iterator new_line;
5316 new_line=triangulation.faces->lines.next_free_single_object(triangulation);
5317 Assert (new_line->used() ==
false,
5318 ExcMessage(
"Internal error: We want to use a cell during refinement that should be unused, but turns out not to be."));
5333 unsigned int vertex_indices[2];
5336 vertex_indices[0]=quad->line(2)->child(0)->vertex_index(1);
5337 vertex_indices[1]=quad->line(3)->child(0)->vertex_index(1);
5341 vertex_indices[0]=quad->line(0)->child(0)->vertex_index(1);
5342 vertex_indices[1]=quad->line(1)->child(0)->vertex_index(1);
5347 new_line->set_used_flag();
5348 new_line->clear_user_flag();
5350 new_line->clear_children();
5351 new_line->set_boundary_id_internal(quad->boundary_id());
5352 new_line->set_manifold_id(quad->manifold_id());
5360 const unsigned int index[2][2]=
5368 typename Triangulation<dim,spacedim>::raw_quad_iterator new_quads[2];
5370 next_unused_quad=triangulation.faces->quads.next_free_pair_object(triangulation);
5371 new_quads[0] = next_unused_quad;
5372 Assert (new_quads[0]->used() ==
false,
ExcMessage(
"Internal error: We want to use a cell during refinement that should be unused, but turns out not to be."));
5375 new_quads[1] = next_unused_quad;
5376 Assert (new_quads[1]->used() ==
false,
ExcMessage(
"Internal error: We want to use a cell during refinement that should be unused, but turns out not to be."));
5384 quad->line(2)->child(index[0][quad->line_orientation(2)])->index(),
5385 quad->line(3)->child(index[0][quad->line_orientation(3)])->index()));
5388 quad->line_index(1),
5389 quad->line(2)->child(index[1][quad->line_orientation(2)])->index(),
5390 quad->line(3)->child(index[1][quad->line_orientation(3)])->index()));
5395 ::
TriaObject<2>(quad->line(0)->child(index[0][quad->line_orientation(0)])->index(),
5396 quad->line(1)->child(index[0][quad->line_orientation(1)])->index(),
5397 quad->line_index(2),
5398 new_line->index()));
5400 ::
TriaObject<2>(quad->line(0)->child(index[1][quad->line_orientation(0)])->index(),
5401 quad->line(1)->child(index[1][quad->line_orientation(1)])->index(),
5403 quad->line_index(3)));
5406 for (
unsigned int i=0; i<2; ++i)
5408 new_quads[i]->set_used_flag();
5409 new_quads[i]->clear_user_flag();
5411 new_quads[i]->clear_children();
5412 new_quads[i]->set_boundary_id_internal(quad->boundary_id());
5413 new_quads[i]->set_manifold_id (quad->manifold_id());
5417 for (
unsigned int j=0; j<GeometryInfo<dim>::lines_per_face; ++j)
5418 new_quads[i]->set_line_orientation(j,
true);
5424 new_quads[0]->set_line_orientation(0,quad->line_orientation(0));
5425 new_quads[0]->set_line_orientation(2,quad->line_orientation(2));
5426 new_quads[1]->set_line_orientation(1,quad->line_orientation(1));
5427 new_quads[1]->set_line_orientation(3,quad->line_orientation(3));
5430 new_quads[0]->set_line_orientation(3,quad->line_orientation(3));
5431 new_quads[1]->set_line_orientation(2,quad->line_orientation(2));
5435 new_quads[0]->set_line_orientation(1,quad->line_orientation(1));
5436 new_quads[1]->set_line_orientation(0,quad->line_orientation(0));
5462 typename Triangulation<dim,spacedim>::line_iterator old_child[2];
5465 old_child[0]=quad->child(0)->line(1);
5466 old_child[1]=quad->child(2)->line(1);
5472 old_child[0]=quad->child(0)->line(3);
5473 old_child[1]=quad->child(1)->line(3);
5476 if (old_child[0]->index()+1 != old_child[1]->index())
5481 typename Triangulation<dim,spacedim>::raw_line_iterator new_child[2];
5483 new_child[0]=new_child[1]=triangulation.faces->lines.next_free_pair_object(triangulation);
5486 new_child[0]->set_used_flag();
5487 new_child[1]->set_used_flag();
5489 const int old_index_0=old_child[0]->index(),
5490 old_index_1=old_child[1]->index(),
5491 new_index_0=new_child[0]->index(),
5492 new_index_1=new_child[1]->index();
5496 for (
unsigned int q=0; q<triangulation.faces->quads.cells.size(); ++q)
5497 for (
unsigned int l=0; l<GeometryInfo<dim>::lines_per_face; ++l)
5499 const int this_index=triangulation.faces->quads.cells[q].face(l);
5500 if (this_index==old_index_0)
5501 triangulation.faces->quads.cells[q].set_face(l,new_index_0);
5502 else if (this_index==old_index_1)
5503 triangulation.faces->quads.cells[q].set_face(l,new_index_1);
5507 for (
unsigned int i=0; i<2; ++i)
5512 old_child[i]->vertex_index(1)));
5513 new_child[i]->set_boundary_id_internal(old_child[i]->boundary_id());
5514 new_child[i]->set_manifold_id(old_child[i]->manifold_id());
5515 new_child[i]->set_user_index(old_child[i]->user_index());
5516 if (old_child[i]->user_flag_set())
5517 new_child[i]->set_user_flag();
5519 new_child[i]->clear_user_flag();
5521 new_child[i]->clear_children();
5523 old_child[i]->clear_user_flag();
5524 old_child[i]->clear_user_index();
5525 old_child[i]->clear_used_flag();
5533 new_line->set_children(0, quad->child(0)->line_index(1));
5534 Assert(new_line->child(1)==quad->child(2)->line(1),
5570 const typename Triangulation<dim,spacedim>::quad_iterator
5571 switch_1=quad->child(1),
5572 switch_2=quad->child(2);
5573 const int switch_1_index=switch_1->index();
5574 const int switch_2_index=switch_2->index();
5575 for (
unsigned int l=0; l<triangulation.levels.size(); ++l)
5576 for (
unsigned int h=0; h<triangulation.levels[l]->cells.cells.size(); ++h)
5577 for (
unsigned int q=0; q<GeometryInfo<dim>::faces_per_cell; ++q)
5579 const int face_index=triangulation.levels[l]->cells.cells[h].face(q);
5580 if (face_index==switch_1_index)
5581 triangulation.levels[l]->cells.cells[h].set_face(q,switch_2_index);
5582 else if (face_index==switch_2_index)
5583 triangulation.levels[l]->cells.cells[h].set_face(q,switch_1_index);
5587 const unsigned int switch_1_lines[4]=
5589 switch_1->line_index(0),
5590 switch_1->line_index(1),
5591 switch_1->line_index(2),
5592 switch_1->line_index(3)
5594 const bool switch_1_line_orientations[4]=
5596 switch_1->line_orientation(0),
5597 switch_1->line_orientation(1),
5598 switch_1->line_orientation(2),
5599 switch_1->line_orientation(3)
5602 const unsigned int switch_1_user_index=switch_1->user_index();
5603 const bool switch_1_user_flag=switch_1->user_flag_set();
5604 const RefinementCase<dim-1> switch_1_refinement_case=switch_1->refinement_case();
5605 const int switch_1_first_child_pair=(switch_1_refinement_case ? switch_1->child_index(0) : -1);
5606 const int switch_1_second_child_pair=(switch_1_refinement_case==
RefinementCase<dim-1>::cut_xy ? switch_1->child_index(2) : -1);
5609 switch_2->line_index(1),
5610 switch_2->line_index(2),
5611 switch_2->line_index(3)));
5612 switch_1->set_line_orientation(0, switch_2->line_orientation(0));
5613 switch_1->set_line_orientation(1, switch_2->line_orientation(1));
5614 switch_1->set_line_orientation(2, switch_2->line_orientation(2));
5615 switch_1->set_line_orientation(3, switch_2->line_orientation(3));
5616 switch_1->set_boundary_id_internal(switch_2->boundary_id());
5617 switch_1->set_manifold_id(switch_2->manifold_id());
5618 switch_1->set_user_index(switch_2->user_index());
5619 if (switch_2->user_flag_set())
5620 switch_1->set_user_flag();
5622 switch_1->clear_user_flag();
5623 switch_1->clear_refinement_case();
5624 switch_1->set_refinement_case(switch_2->refinement_case());
5625 switch_1->clear_children();
5626 if (switch_2->refinement_case())
5627 switch_1->set_children(0, switch_2->child_index(0));
5629 switch_1->set_children(2, switch_2->child_index(2));
5634 switch_1_lines[3]));
5635 switch_2->set_line_orientation(0, switch_1_line_orientations[0]);
5636 switch_2->set_line_orientation(1, switch_1_line_orientations[1]);
5637 switch_2->set_line_orientation(2, switch_1_line_orientations[2]);
5638 switch_2->set_line_orientation(3, switch_1_line_orientations[3]);
5639 switch_2->set_boundary_id_internal(switch_1_boundary_id);
5640 switch_2->set_manifold_id(switch_1->manifold_id());
5641 switch_2->set_user_index(switch_1_user_index);
5642 if (switch_1_user_flag)
5643 switch_2->set_user_flag();
5645 switch_2->clear_user_flag();
5646 switch_2->clear_refinement_case();
5647 switch_2->set_refinement_case(switch_1_refinement_case);
5648 switch_2->clear_children();
5649 switch_2->set_children(0, switch_1_first_child_pair);
5650 switch_2->set_children(2, switch_1_second_child_pair);
5653 new_quads[0]->set_children(0, quad->child_index(0));
5655 new_quads[1]->set_children(0, quad->child_index(2));
5660 new_quads[0]->set_children(0, quad->child_index(0));
5662 new_quads[1]->set_children(0, quad->child_index(2));
5663 new_line->set_children(0, quad->child(0)->line_index(3));
5664 Assert(new_line->child(1)==quad->child(1)->line(3),
5667 quad->clear_children();
5671 quad->set_children (0, new_quads[0]->index());
5673 quad->set_refinement_case(aniso_quad_ref_case);
5680 if (quad->user_flag_set())
5690 while (triangulation.vertices_used[next_unused_vertex] ==
true)
5691 ++next_unused_vertex;
5692 Assert (next_unused_vertex < triangulation.vertices.size(),
5693 ExcMessage(
"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."));
5700 const RefinementCase<dim-1> quad_ref_case=quad->refinement_case();
5709 typename Triangulation<dim,spacedim>::line_iterator middle_line;
5711 middle_line=quad->child(0)->line(1);
5713 middle_line=quad->child(0)->line(3);
5720 if (!middle_line->has_children())
5727 triangulation.vertices[next_unused_vertex]
5728 = middle_line->center(
true);
5729 triangulation.vertices_used[next_unused_vertex] =
true;
5733 next_unused_line=triangulation.faces->lines.next_free_pair_object(triangulation);
5737 middle_line->set_children (0, next_unused_line->index());
5740 const typename Triangulation<dim,spacedim>::raw_line_iterator
5741 children[2] = { next_unused_line,
5748 Assert (children[0]->used() ==
false,
ExcMessage(
"Internal error: We want to use a cell during refinement that should be unused, but turns out not to be."));
5749 Assert (children[1]->used() ==
false,
ExcMessage(
"Internal error: We want to use a cell during refinement that should be unused, but turns out not to be."));
5753 next_unused_vertex));
5756 middle_line->vertex_index(1)));
5758 children[0]->set_used_flag();
5759 children[1]->set_used_flag();
5760 children[0]->clear_children();
5761 children[1]->clear_children();
5764 children[0]->clear_user_flag();
5765 children[1]->clear_user_flag();
5767 children[0]->set_boundary_id_internal (middle_line->boundary_id());
5768 children[1]->set_boundary_id_internal (middle_line->boundary_id());
5770 children[0]->set_manifold_id (middle_line->manifold_id());
5771 children[1]->set_manifold_id (middle_line->manifold_id());
5777 quad->clear_user_flag();
5804 triangulation.vertices[next_unused_vertex] =
5805 quad->center(
true,
true);
5806 triangulation.vertices_used[next_unused_vertex] =
true;
5811 typename Triangulation<dim,spacedim>::raw_line_iterator new_lines[4];
5813 for (
unsigned int i=0; i<4; ++i)
5822 next_unused_line=triangulation.faces->lines.next_free_pair_object(triangulation);
5824 new_lines[i] = next_unused_line;
5827 Assert (new_lines[i]->used() ==
false,
5828 ExcMessage(
"Internal error: We want to use a cell during refinement that should be unused, but turns out not to be."));
5848 const unsigned int vertex_indices[5]
5849 = { quad->line(0)->child(0)->vertex_index(1),
5850 quad->line(1)->child(0)->vertex_index(1),
5851 quad->line(2)->child(0)->vertex_index(1),
5852 quad->line(3)->child(0)->vertex_index(1),
5865 for (
unsigned int i=0; i<4; ++i)
5867 new_lines[i]->set_used_flag();
5868 new_lines[i]->clear_user_flag();
5870 new_lines[i]->clear_children();
5871 new_lines[i]->set_boundary_id_internal(quad->boundary_id());
5872 new_lines[i]->set_manifold_id(quad->manifold_id());
5891 const unsigned int index[2][2]=
5897 const int line_indices[12]
5898 = { quad->line(0)->child(index[0][quad->line_orientation(0)])->index(),
5899 quad->line(0)->child(index[1][quad->line_orientation(0)])->index(),
5900 quad->line(1)->child(index[0][quad->line_orientation(1)])->index(),
5901 quad->line(1)->child(index[1][quad->line_orientation(1)])->index(),
5902 quad->line(2)->child(index[0][quad->line_orientation(2)])->index(),
5903 quad->line(2)->child(index[1][quad->line_orientation(2)])->index(),
5904 quad->line(3)->child(index[0][quad->line_orientation(3)])->index(),
5905 quad->line(3)->child(index[1][quad->line_orientation(3)])->index(),
5906 new_lines[0]->index(),
5907 new_lines[1]->index(),
5908 new_lines[2]->index(),
5909 new_lines[3]->index()
5915 typename Triangulation<dim,spacedim>::raw_quad_iterator new_quads[4];
5917 next_unused_quad=triangulation.faces->quads.next_free_pair_object(triangulation);
5919 new_quads[0] = next_unused_quad;
5920 Assert (new_quads[0]->used() ==
false,
ExcMessage(
"Internal error: We want to use a cell during refinement that should be unused, but turns out not to be."));
5923 new_quads[1] = next_unused_quad;
5924 Assert (new_quads[1]->used() ==
false,
ExcMessage(
"Internal error: We want to use a cell during refinement that should be unused, but turns out not to be."));
5926 next_unused_quad=triangulation.faces->quads.next_free_pair_object(triangulation);
5927 new_quads[2] = next_unused_quad;
5928 Assert (new_quads[2]->used() ==
false,
ExcMessage(
"Internal error: We want to use a cell during refinement that should be unused, but turns out not to be."));
5931 new_quads[3] = next_unused_quad;
5932 Assert (new_quads[3]->used() ==
false,
ExcMessage(
"Internal error: We want to use a cell during refinement that should be unused, but turns out not to be."));
5935 quad->set_children (0, new_quads[0]->index());
5936 quad->set_children (2, new_quads[2]->index());
5965 for (
unsigned int i=0; i<4; ++i)
5967 new_quads[i]->set_used_flag();
5968 new_quads[i]->clear_user_flag();
5970 new_quads[i]->clear_children();
5971 new_quads[i]->set_boundary_id_internal (quad->boundary_id());
5972 new_quads[i]->set_manifold_id (quad->manifold_id());
5976 for (
unsigned int j=0; j<GeometryInfo<dim>::lines_per_face; ++j)
5977 new_quads[i]->set_line_orientation(j,
true);
5983 new_quads[0]->set_line_orientation(0,quad->line_orientation(0));
5984 new_quads[0]->set_line_orientation(2,quad->line_orientation(2));
5985 new_quads[1]->set_line_orientation(1,quad->line_orientation(1));
5986 new_quads[1]->set_line_orientation(2,quad->line_orientation(2));
5987 new_quads[2]->set_line_orientation(0,quad->line_orientation(0));
5988 new_quads[2]->set_line_orientation(3,quad->line_orientation(3));
5989 new_quads[3]->set_line_orientation(1,quad->line_orientation(1));
5990 new_quads[3]->set_line_orientation(3,quad->line_orientation(3));
5994 quad->clear_user_flag ();
6005 cells_with_distorted_children;
6007 for (
unsigned int level=0; level!=triangulation.levels.size()-1; ++level)
6012 typename Triangulation<dim,spacedim>::active_hex_iterator
6013 hex = triangulation.begin_active_hex(level),
6014 endh = triangulation.begin_active_hex(level+1);
6015 typename Triangulation<dim,spacedim>::raw_hex_iterator
6016 next_unused_hex = triangulation.begin_raw_hex (level+1);
6018 for (; hex!=endh; ++hex)
6019 if (hex->refine_flag_set())
6027 hex->clear_refine_flag ();
6028 hex->set_refinement_case(ref_case);
6039 unsigned int n_new_lines=0;
6040 unsigned int n_new_quads=0;
6041 unsigned int n_new_hexes=0;
6070 std::vector<typename Triangulation<dim,spacedim>::raw_line_iterator>
6071 new_lines(n_new_lines);
6072 for (
unsigned int i=0; i<n_new_lines; ++i)
6074 new_lines[i] = triangulation.faces->lines.next_free_single_object(triangulation);
6076 Assert (new_lines[i]->used() ==
false,
6077 ExcMessage(
"Internal error: We want to use a cell during refinement that should be unused, but turns out not to be."));
6078 new_lines[i]->set_used_flag();
6079 new_lines[i]->clear_user_flag();
6080 new_lines[i]->clear_user_data();
6081 new_lines[i]->clear_children();
6085 new_lines[i]->set_manifold_id(hex->manifold_id());
6090 std::vector<typename Triangulation<dim,spacedim>::raw_quad_iterator>
6091 new_quads(n_new_quads);
6092 for (
unsigned int i=0; i<n_new_quads; ++i)
6094 new_quads[i] = triangulation.faces->quads.next_free_single_object(triangulation);
6096 Assert (new_quads[i]->used() ==
false,
6097 ExcMessage(
"Internal error: We want to use a cell during refinement that should be unused, but turns out not to be."));
6098 new_quads[i]->set_used_flag();
6099 new_quads[i]->clear_user_flag();
6100 new_quads[i]->clear_user_data();
6101 new_quads[i]->clear_children();
6105 new_quads[i]->set_manifold_id (hex->manifold_id());
6108 for (
unsigned int j=0; j<GeometryInfo<dim>::lines_per_face; ++j)
6109 new_quads[i]->set_line_orientation(j,
true);
6116 std::vector<typename Triangulation<dim,spacedim>::raw_hex_iterator>
6117 new_hexes(n_new_hexes);
6118 for (
unsigned int i=0; i<n_new_hexes; ++i)
6121 next_unused_hex=triangulation.levels[level+1]->cells.next_free_hex(triangulation,level+1);
6125 new_hexes[i]=next_unused_hex;
6127 Assert (new_hexes[i]->used() ==
false,
6128 ExcMessage(
"Internal error: We want to use a cell during refinement that should be unused, but turns out not to be."));
6129 new_hexes[i]->set_used_flag();
6130 new_hexes[i]->clear_user_flag();
6131 new_hexes[i]->clear_user_data();
6132 new_hexes[i]->clear_children();
6135 new_hexes[i]->set_material_id (hex->material_id());
6136 new_hexes[i]->set_manifold_id (hex->manifold_id());
6137 new_hexes[i]->set_subdomain_id (subdomainid);
6140 new_hexes[i]->set_parent (hex->index ());
6152 for (
unsigned int f=0; f<GeometryInfo<dim>::faces_per_cell; ++f)
6154 new_hexes[i]->set_face_orientation(f,
true);
6155 new_hexes[i]->set_face_flip(f,
false);
6156 new_hexes[i]->set_face_rotation(f,
false);
6160 for (
unsigned int i=0; i<n_new_hexes/2; ++i)
6161 hex->set_children (2*i, new_hexes[2*i]->index());
6169 = { hex->face_orientation (0),
6170 hex->face_orientation (1),
6171 hex->face_orientation (2),
6172 hex->face_orientation (3),
6173 hex->face_orientation (4),
6174 hex->face_orientation (5)
6179 = { hex->face_flip (0),
6189 = { hex->face_rotation (0),
6190 hex->face_rotation (1),
6191 hex->face_rotation (2),
6192 hex->face_rotation (3),
6193 hex->face_rotation (4),
6194 hex->face_rotation (5)
6205 const unsigned int child_at_origin[2][2][2]=
6291 const typename Triangulation<dim,spacedim>::raw_line_iterator
6295 hex->face(2)->child(0)
6297 hex->face(3)->child(0)
6299 hex->face(4)->child(0)
6301 hex->face(5)->child(0)
6305 unsigned int line_indices[4];
6306 for (
unsigned int i=0; i<4; ++i)
6307 line_indices[i] = lines[i]->index();
6314 bool line_orientation[4];
6320 const unsigned int middle_vertices[2]=
6322 hex->line(2)->child(0)->vertex_index(1),
6323 hex->line(7)->child(0)->vertex_index(1)
6326 for (
unsigned int i=0; i<4; ++i)
6327 if (lines[i]->vertex_index(i%2)==middle_vertices[i%2])
6328 line_orientation[i]=
true;
6333 Assert(lines[i]->vertex_index((i+1)%2)==middle_vertices[i%2],
6335 line_orientation[i]=
false;
6346 new_quads[0]->set_line_orientation(0,line_orientation[0]);
6347 new_quads[0]->set_line_orientation(1,line_orientation[1]);
6348 new_quads[0]->set_line_orientation(2,line_orientation[2]);
6349 new_quads[0]->set_line_orientation(3,line_orientation[3]);
6382 const int quad_indices[11]
6385 new_quads[0]->index(),
6387 hex->face(0)->index(),
6389 hex->face(1)->index(),
6391 hex->face(2)->child_index( child_at_origin[hex->face(2)->refinement_case()-1][f_fl[2]][f_ro[2]]),
6392 hex->face(2)->child_index(1-child_at_origin[hex->face(2)->refinement_case()-1][f_fl[2]][f_ro[2]]),
6394 hex->face(3)->child_index( child_at_origin[hex->face(3)->refinement_case()-1][f_fl[3]][f_ro[3]]),
6395 hex->face(3)->child_index(1-child_at_origin[hex->face(3)->refinement_case()-1][f_fl[3]][f_ro[3]]),
6397 hex->face(4)->child_index( child_at_origin[hex->face(4)->refinement_case()-1][f_fl[4]][f_ro[4]]),
6398 hex->face(4)->child_index(1-child_at_origin[hex->face(4)->refinement_case()-1][f_fl[4]][f_ro[4]]),
6400 hex->face(5)->child_index( child_at_origin[hex->face(5)->refinement_case()-1][f_fl[5]][f_ro[5]]),
6401 hex->face(5)->child_index(1-child_at_origin[hex->face(5)->refinement_case()-1][f_fl[5]][f_ro[5]])
6486 const typename Triangulation<dim,spacedim>::raw_line_iterator
6490 hex->face(0)->child(0)
6492 hex->face(1)->child(0)
6494 hex->face(4)->child(0)
6496 hex->face(5)->child(0)
6500 unsigned int line_indices[4];
6501 for (
unsigned int i=0; i<4; ++i)
6502 line_indices[i]=lines[i]->index();
6509 bool line_orientation[4];
6515 const unsigned int middle_vertices[2]=
6517 hex->line(0)->child(0)->vertex_index(1),
6518 hex->line(5)->child(0)->vertex_index(1)
6521 for (
unsigned int i=0; i<4; ++i)
6522 if (lines[i]->vertex_index(i%2)==middle_vertices[i%2])
6523 line_orientation[i]=
true;
6527 Assert(lines[i]->vertex_index((i+1)%2)==middle_vertices[i%2],
6529 line_orientation[i]=
false;
6540 new_quads[0]->set_line_orientation(0,line_orientation[2]);
6541 new_quads[0]->set_line_orientation(1,line_orientation[3]);
6542 new_quads[0]->set_line_orientation(2,line_orientation[0]);
6543 new_quads[0]->set_line_orientation(3,line_orientation[1]);
6576 const int quad_indices[11]
6579 new_quads[0]->index(),
6581 hex->face(0)->child_index( child_at_origin[hex->face(0)->refinement_case()-1][f_fl[0]][f_ro[0]]),
6582 hex->face(0)->child_index(1-child_at_origin[hex->face(0)->refinement_case()-1][f_fl[0]][f_ro[0]]),
6584 hex->face(1)->child_index( child_at_origin[hex->face(1)->refinement_case()-1][f_fl[1]][f_ro[1]]),
6585 hex->face(1)->child_index(1-child_at_origin[hex->face(1)->refinement_case()-1][f_fl[1]][f_ro[1]]),
6587 hex->face(2)->index(),
6589 hex->face(3)->index(),
6591 hex->face(4)->child_index( child_at_origin[hex->face(4)->refinement_case()-1][f_fl[4]][f_ro[4]]),
6592 hex->face(4)->child_index(1-child_at_origin[hex->face(4)->refinement_case()-1][f_fl[4]][f_ro[4]]),
6594 hex->face(5)->child_index( child_at_origin[hex->face(5)->refinement_case()-1][f_fl[5]][f_ro[5]]),
6595 hex->face(5)->child_index(1-child_at_origin[hex->face(5)->refinement_case()-1][f_fl[5]][f_ro[5]])
6682 const typename Triangulation<dim,spacedim>::raw_line_iterator
6686 hex->face(0)->child(0)
6688 hex->face(1)->child(0)
6690 hex->face(2)->child(0)
6692 hex->face(3)->child(0)
6696 unsigned int line_indices[4];
6697 for (
unsigned int i=0; i<4; ++i)
6698 line_indices[i]=lines[i]->index();
6705 bool line_orientation[4];
6711 const unsigned int middle_vertices[2]=
6713 middle_vertex_index<dim,spacedim>(hex->line(8)),
6714 middle_vertex_index<dim,spacedim>(hex->line(11))
6717 for (
unsigned int i=0; i<4; ++i)
6718 if (lines[i]->vertex_index(i%2)==middle_vertices[i%2])
6719 line_orientation[i]=
true;
6723 Assert(lines[i]->vertex_index((i+1)%2)==middle_vertices[i%2],
6725 line_orientation[i]=
false;
6736 new_quads[0]->set_line_orientation(0,line_orientation[0]);
6737 new_quads[0]->set_line_orientation(1,line_orientation[1]);
6738 new_quads[0]->set_line_orientation(2,line_orientation[2]);
6739 new_quads[0]->set_line_orientation(3,line_orientation[3]);
6773 const int quad_indices[11]
6776 new_quads[0]->index(),
6778 hex->face(0)->child_index( child_at_origin[hex->face(0)->refinement_case()-1][f_fl[0]][f_ro[0]]),
6779 hex->face(0)->child_index(1-child_at_origin[hex->face(0)->refinement_case()-1][f_fl[0]][f_ro[0]]),
6781 hex->face(1)->child_index( child_at_origin[hex->face(1)->refinement_case()-1][f_fl[1]][f_ro[1]]),
6782 hex->face(1)->child_index(1-child_at_origin[hex->face(1)->refinement_case()-1][f_fl[1]][f_ro[1]]),
6784 hex->face(2)->child_index( child_at_origin[hex->face(2)->refinement_case()-1][f_fl[2]][f_ro[2]]),
6785 hex->face(2)->child_index(1-child_at_origin[hex->face(2)->refinement_case()-1][f_fl[2]][f_ro[2]]),
6787 hex->face(3)->child_index( child_at_origin[hex->face(3)->refinement_case()-1][f_fl[3]][f_ro[3]]),
6788 hex->face(3)->child_index(1-child_at_origin[hex->face(3)->refinement_case()-1][f_fl[3]][f_ro[3]]),
6790 hex->face(4)->index(),
6792 hex->face(5)->index()
6835 TriaObject<1>(middle_vertex_index<dim,spacedim>(hex->face(4)),
6836 middle_vertex_index<dim,spacedim>(hex->face(5))));
6902 const typename Triangulation<dim,spacedim>::raw_line_iterator
6906 hex->face(0)->child(0)
6908 hex->face(1)->child(0)
6910 hex->face(2)->child(0)
6912 hex->face(3)->child(0)
6936 unsigned int line_indices[13];
6937 for (
unsigned int i=0; i<13; ++i)
6938 line_indices[i]=lines[i]->index();
6945 bool line_orientation[13];
6949 const unsigned int middle_vertices[4]=
6951 hex->line(0)->child(0)->vertex_index(1),
6952 hex->line(1)->child(0)->vertex_index(1),
6953 hex->line(2)->child(0)->vertex_index(1),
6954 hex->line(3)->child(0)->vertex_index(1),
6960 for (
unsigned int i=0; i<4; ++i)
6961 if (lines[i]->vertex_index(0)==middle_vertices[i])
6962 line_orientation[i]=
true;
6966 Assert(lines[i]->vertex_index(1)==middle_vertices[i],
6968 line_orientation[i]=
false;
6977 for (
unsigned int i=4; i<12; ++i)
6978 if (lines[i]->vertex_index((i+1)%2) ==
6979 middle_vertex_index<dim,spacedim>(hex->face(3+i/4)))
6980 line_orientation[i]=
true;
6985 Assert(lines[i]->vertex_index(i%2) ==
6986 (middle_vertex_index<dim,spacedim>(hex->face(3+i/4))),
6988 line_orientation[i]=
false;
6993 line_orientation[12]=
true;
7039 new_quads[0]->set_line_orientation(0,line_orientation[2]);
7040 new_quads[0]->set_line_orientation(2,line_orientation[4]);
7041 new_quads[0]->set_line_orientation(3,line_orientation[8]);
7043 new_quads[1]->set_line_orientation(1,line_orientation[3]);
7044 new_quads[1]->set_line_orientation(2,line_orientation[5]);
7045 new_quads[1]->set_line_orientation(3,line_orientation[9]);
7047 new_quads[2]->set_line_orientation(0,line_orientation[6]);
7048 new_quads[2]->set_line_orientation(1,line_orientation[10]);
7049 new_quads[2]->set_line_orientation(2,line_orientation[0]);
7051 new_quads[3]->set_line_orientation(0,line_orientation[7]);
7052 new_quads[3]->set_line_orientation(1,line_orientation[11]);
7053 new_quads[3]->set_line_orientation(3,line_orientation[1]);
7086 const int quad_indices[20]
7089 new_quads[0]->index(),
7090 new_quads[1]->index(),
7091 new_quads[2]->index(),
7092 new_quads[3]->index(),
7094 hex->face(0)->child_index( child_at_origin[hex->face(0)->refinement_case()-1][f_fl[0]][f_ro[0]]),
7095 hex->face(0)->child_index(1-child_at_origin[hex->face(0)->refinement_case()-1][f_fl[0]][f_ro[0]]),
7097 hex->face(1)->child_index( child_at_origin[hex->face(1)->refinement_case()-1][f_fl[1]][f_ro[1]]),
7098 hex->face(1)->child_index(1-child_at_origin[hex->face(1)->refinement_case()-1][f_fl[1]][f_ro[1]]),
7100 hex->face(2)->child_index( child_at_origin[hex->face(2)->refinement_case()-1][f_fl[2]][f_ro[2]]),
7101 hex->face(2)->child_index(1-child_at_origin[hex->face(2)->refinement_case()-1][f_fl[2]][f_ro[2]]),
7103 hex->face(3)->child_index( child_at_origin[hex->face(3)->refinement_case()-1][f_fl[3]][f_ro[3]]),
7104 hex->face(3)->child_index(1-child_at_origin[hex->face(3)->refinement_case()-1][f_fl[3]][f_ro[3]]),
7171 TriaObject<1>(middle_vertex_index<dim,spacedim>(hex->face(2)),
7172 middle_vertex_index<dim,spacedim>(hex->face(3))));
7238 const typename Triangulation<dim,spacedim>::raw_line_iterator
7242 hex->face(0)->child(0)
7244 hex->face(1)->child(0)
7246 hex->face(4)->child(0)
7248 hex->face(5)->child(0)
7272 unsigned int line_indices[13];
7273 for (
unsigned int i=0; i<13; ++i)
7274 line_indices[i]=lines[i]->index();
7281 bool line_orientation[13];
7285 const unsigned int middle_vertices[4]=
7287 hex->line(8)->child(0)->vertex_index(1),
7288 hex->line(9)->child(0)->vertex_index(1),
7289 hex->line(2)->child(0)->vertex_index(1),
7290 hex->line(6)->child(0)->vertex_index(1),
7295 for (
unsigned int i=0; i<4; ++i)
7296 if (lines[i]->vertex_index(0)==middle_vertices[i])
7297 line_orientation[i]=
true;
7301 Assert(lines[i]->vertex_index(1)==middle_vertices[i],
7303 line_orientation[i]=
false;
7312 for (
unsigned int i=4; i<12; ++i)
7313 if (lines[i]->vertex_index((i+1)%2) ==
7314 middle_vertex_index<dim,spacedim>(hex->face(1+i/4)))
7315 line_orientation[i]=
true;
7320 Assert(lines[i]->vertex_index(i%2) ==
7321 (middle_vertex_index<dim,spacedim>(hex->face(1+i/4))),
7323 line_orientation[i]=
false;
7328 line_orientation[12]=
true;
7375 new_quads[0]->set_line_orientation(0,line_orientation[0]);
7376 new_quads[0]->set_line_orientation(2,line_orientation[6]);
7377 new_quads[0]->set_line_orientation(3,line_orientation[10]);
7379 new_quads[1]->set_line_orientation(1,line_orientation[1]);
7380 new_quads[1]->set_line_orientation(2,line_orientation[7]);
7381 new_quads[1]->set_line_orientation(3,line_orientation[11]);
7383 new_quads[2]->set_line_orientation(0,line_orientation[4]);
7384 new_quads[2]->set_line_orientation(1,line_orientation[8]);
7385 new_quads[2]->set_line_orientation(2,line_orientation[2]);
7387 new_quads[3]->set_line_orientation(0,line_orientation[5]);
7388 new_quads[3]->set_line_orientation(1,line_orientation[9]);
7389 new_quads[3]->set_line_orientation(3,line_orientation[3]);
7423 const int quad_indices[20]
7426 new_quads[0]->index(),
7427 new_quads[1]->index(),
7428 new_quads[2]->index(),
7429 new_quads[3]->index(),
7431 hex->face(0)->child_index( child_at_origin[hex->face(0)->refinement_case()-1][f_fl[0]][f_ro[0]]),
7432 hex->face(0)->child_index(1-child_at_origin[hex->face(0)->refinement_case()-1][f_fl[0]][f_ro[0]]),
7434 hex->face(1)->child_index( child_at_origin[hex->face(1)->refinement_case()-1][f_fl[1]][f_ro[1]]),
7435 hex->face(1)->child_index(1-child_at_origin[hex->face(1)->refinement_case()-1][f_fl[1]][f_ro[1]]),
7447 hex->face(4)->child_index( child_at_origin[hex->face(4)->refinement_case()-1][f_fl[4]][f_ro[4]]),
7448 hex->face(4)->child_index(1-child_at_origin[hex->face(4)->refinement_case()-1][f_fl[4]][f_ro[4]]),
7450 hex->face(5)->child_index( child_at_origin[hex->face(5)->refinement_case()-1][f_fl[5]][f_ro[5]]),
7451 hex->face(5)->child_index(1-child_at_origin[hex->face(5)->refinement_case()-1][f_fl[5]][f_ro[5]])
7518 TriaObject<1>(middle_vertex_index<dim,spacedim>(hex->face(0)),
7519 middle_vertex_index<dim,spacedim>(hex->face(1))));
7586 const typename Triangulation<dim,spacedim>::raw_line_iterator
7590 hex->face(2)->child(0)
7592 hex->face(3)->child(0)
7594 hex->face(4)->child(0)
7596 hex->face(5)->child(0)
7620 unsigned int line_indices[13];
7622 for (
unsigned int i=0; i<13; ++i)
7623 line_indices[i]=lines[i]->index();
7630 bool line_orientation[13];
7634 const unsigned int middle_vertices[4]=
7636 hex->line(8)->child(0)->vertex_index(1),
7637 hex->line(10)->child(0)->vertex_index(1),
7638 hex->line(0)->child(0)->vertex_index(1),
7639 hex->line(4)->child(0)->vertex_index(1),
7644 for (
unsigned int i=0; i<4; ++i)
7645 if (lines[i]->vertex_index(0)==middle_vertices[i])
7646 line_orientation[i]=
true;
7650 Assert(lines[i]->vertex_index(1)==middle_vertices[i],
7652 line_orientation[i]=
false;
7661 for (
unsigned int i=4; i<12; ++i)
7662 if (lines[i]->vertex_index((i+1)%2) ==
7663 middle_vertex_index<dim,spacedim>(hex->face(i/4-1)))
7664 line_orientation[i]=
true;
7669 Assert(lines[i]->vertex_index(i%2) ==
7670 (middle_vertex_index<dim,spacedim>(hex->face(i/4-1))),
7672 line_orientation[i]=
false;
7677 line_orientation[12]=
true;
7718 new_quads[0]->set_line_orientation(0,line_orientation[6]);
7719 new_quads[0]->set_line_orientation(1,line_orientation[10]);
7720 new_quads[0]->set_line_orientation(2,line_orientation[0]);
7722 new_quads[1]->set_line_orientation(0,line_orientation[7]);
7723 new_quads[1]->set_line_orientation(1,line_orientation[11]);
7724 new_quads[1]->set_line_orientation(3,line_orientation[1]);
7726 new_quads[2]->set_line_orientation(0,line_orientation[2]);
7727 new_quads[2]->set_line_orientation(2,line_orientation[4]);
7728 new_quads[2]->set_line_orientation(3,line_orientation[8]);
7730 new_quads[3]->set_line_orientation(1,line_orientation[3]);
7731 new_quads[3]->set_line_orientation(2,line_orientation[5]);
7732 new_quads[3]->set_line_orientation(3,line_orientation[9]);
7766 const int quad_indices[20]
7769 new_quads[0]->index(),
7770 new_quads[1]->index(),
7771 new_quads[2]->index(),
7772 new_quads[3]->index(),
7784 hex->face(2)->child_index( child_at_origin[hex->face(2)->refinement_case()-1][f_fl[2]][f_ro[2]]),
7785 hex->face(2)->child_index(1-child_at_origin[hex->face(2)->refinement_case()-1][f_fl[2]][f_ro[2]]),
7787 hex->face(3)->child_index( child_at_origin[hex->face(3)->refinement_case()-1][f_fl[3]][f_ro[3]]),
7788 hex->face(3)->child_index(1-child_at_origin[hex->face(3)->refinement_case()-1][f_fl[3]][f_ro[3]]),
7790 hex->face(4)->child_index( child_at_origin[hex->face(4)->refinement_case()-1][f_fl[4]][f_ro[4]]),
7791 hex->face(4)->child_index(1-child_at_origin[hex->face(4)->refinement_case()-1][f_fl[4]][f_ro[4]]),
7793 hex->face(5)->child_index( child_at_origin[hex->face(5)->refinement_case()-1][f_fl[5]][f_ro[5]]),
7794 hex->face(5)->child_index(1-child_at_origin[hex->face(5)->refinement_case()-1][f_fl[5]][f_ro[5]])
7853 while (triangulation.vertices_used[next_unused_vertex] ==
true)
7854 ++next_unused_vertex;
7855 Assert (next_unused_vertex < triangulation.vertices.size(),
7856 ExcMessage(
"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."));
7857 triangulation.vertices_used[next_unused_vertex] =
true;
7865 triangulation.vertices[next_unused_vertex] =
7866 hex->center(
true,
true);
7887 const unsigned int vertex_indices[7]
7888 = { middle_vertex_index<dim,spacedim>(hex->face(0)),
7889 middle_vertex_index<dim,spacedim>(hex->face(1)),
7890 middle_vertex_index<dim,spacedim>(hex->face(2)),
7891 middle_vertex_index<dim,spacedim>(hex->face(3)),
7892 middle_vertex_index<dim,spacedim>(hex->face(4)),
7893 middle_vertex_index<dim,spacedim>(hex->face(5)),
7975 const typename Triangulation<dim,spacedim>::raw_line_iterator
8041 unsigned int line_indices[30];
8042 for (
unsigned int i=0; i<30; ++i)
8043 line_indices[i]=lines[i]->index();
8050 bool line_orientation[30];
8058 for (
unsigned int i=0; i<24; ++i)
8059 if (lines[i]->vertex_index((i+1)%2)==vertex_indices[i/4])
8060 line_orientation[i]=
true;
8065 Assert(lines[i]->vertex_index(i%2)==vertex_indices[i/4],
8067 line_orientation[i]=
false;
8072 for (
unsigned int i=24; i<30; ++i)
8073 line_orientation[i]=
true;
8169 new_quads[0]->set_line_orientation(0,line_orientation[10]);
8170 new_quads[0]->set_line_orientation(2,line_orientation[16]);
8172 new_quads[1]->set_line_orientation(1,line_orientation[14]);
8173 new_quads[1]->set_line_orientation(2,line_orientation[17]);
8175 new_quads[2]->set_line_orientation(0,line_orientation[11]);
8176 new_quads[2]->set_line_orientation(3,line_orientation[20]);
8178 new_quads[3]->set_line_orientation(1,line_orientation[15]);
8179 new_quads[3]->set_line_orientation(3,line_orientation[21]);
8181 new_quads[4]->set_line_orientation(0,line_orientation[18]);
8182 new_quads[4]->set_line_orientation(2,line_orientation[0]);
8184 new_quads[5]->set_line_orientation(1,line_orientation[22]);
8185 new_quads[5]->set_line_orientation(2,line_orientation[1]);
8187 new_quads[6]->set_line_orientation(0,line_orientation[19]);
8188 new_quads[6]->set_line_orientation(3,line_orientation[4]);
8190 new_quads[7]->set_line_orientation(1,line_orientation[23]);
8191 new_quads[7]->set_line_orientation(3,line_orientation[5]);
8193 new_quads[8]->set_line_orientation(0,line_orientation[2]);
8194 new_quads[8]->set_line_orientation(2,line_orientation[8]);
8196 new_quads[9]->set_line_orientation(1,line_orientation[6]);
8197 new_quads[9]->set_line_orientation(2,line_orientation[9]);
8199 new_quads[10]->set_line_orientation(0,line_orientation[3]);
8200 new_quads[10]->set_line_orientation(3,line_orientation[12]);
8202 new_quads[11]->set_line_orientation(1,line_orientation[7]);
8203 new_quads[11]->set_line_orientation(3,line_orientation[13]);
8244 const int quad_indices[36]
8247 new_quads[0]->index(),
8248 new_quads[1]->index(),
8249 new_quads[2]->index(),
8250 new_quads[3]->index(),
8251 new_quads[4]->index(),
8252 new_quads[5]->index(),
8253 new_quads[6]->index(),
8254 new_quads[7]->index(),
8255 new_quads[8]->index(),
8256 new_quads[9]->index(),
8257 new_quads[10]->index(),
8258 new_quads[11]->index(),
8382 for (
unsigned int f=0; f<GeometryInfo<dim>::faces_per_cell; ++f)
8383 for (
unsigned int s=0;
8388 const unsigned int current_child
8400 new_hexes[current_child]->set_face_orientation (f, f_or[f]);
8401 new_hexes[current_child]->set_face_flip (f, f_fl[f]);
8402 new_hexes[current_child]->set_face_rotation (f, f_ro[f]);
8407 if ((check_for_distorted_cells ==
true)
8409 has_distorted_children (hex,
8410 std::integral_constant<int, dim>(),
8411 std::integral_constant<int, spacedim>()))
8418 triangulation.signals.post_refinement_on_cell(hex);
8427 triangulation.faces->quads.clear_user_data();
8430 return cells_with_distorted_children;
8446 template <
int spacedim>
8452 template <
int dim,
int spacedim>
8459 if (spacedim>dim)
return;
8462 cell=triangulation.begin(); cell!=triangulation.end(); ++cell)
8463 if (cell->at_boundary() &&
8464 cell->refine_flag_set() &&
8471 for (
unsigned int face_no=0;
8472 face_no<GeometryInfo<dim>::faces_per_cell;
8474 if (cell->face(face_no)->at_boundary())
8484 face = cell->face(face_no);
8488 = face->center(
true);
8493 transform_real_to_unit_cell(cell,
8504 const double dist = std::fabs(new_unit[face_no/2] - face_no%2);
8509 const double allowed=0.25;
8512 cell->flag_for_face_refinement(face_no);
8531 template <
int dim,
int spacedim>
8537 ExcMessage (
"Wrong function called -- there should " 8538 "be a specialization."));
8542 template <
int spacedim>
8547 const unsigned int dim = 3;
8551 triangulation.clear_user_flags_line();
8555 triangulation.clear_user_flags_hex();
8559 bool mesh_changed =
false;
8563 mesh_changed =
false;
8568 triangulation.fix_coarsen_flags ();
8574 cell=triangulation.begin(); cell!=triangulation.end(); ++cell)
8575 if (cell->refine_flag_set())
8577 for (
unsigned int line=0; line<GeometryInfo<dim>::lines_per_cell; ++line)
8582 cell->line(line)->set_user_flag();
8584 else if (cell->has_children() && !cell->child(0)->coarsen_flag_set())
8586 for (
unsigned int line=0; line<GeometryInfo<dim>::lines_per_cell; ++line)
8591 cell->line(line)->set_user_flag();
8593 else if (cell->has_children() && cell->child(0)->coarsen_flag_set())
8594 cell->set_user_flag();
8602 cell=triangulation.last_active(); cell!=triangulation.end(); --cell)
8603 for (
unsigned int line=0; line<GeometryInfo<dim>::lines_per_cell; ++line)
8605 if (cell->line(line)->has_children())
8614 bool offending_line_found =
false;
8616 for (
unsigned int c=0; c<2; ++c)
8618 Assert (cell->line(line)->child(c)->has_children() ==
false,
8621 if (cell->line(line)->child(c)->user_flag_set () &&
8627 cell->clear_coarsen_flag ();
8632 if (triangulation.smooth_grid &
8634 cell->flag_for_line_refinement(line);
8636 cell->set_refine_flag();
8638 for (
unsigned int l=0; l<GeometryInfo<dim>::lines_per_cell; ++l)
8642 cell->line(l)->set_user_flag();
8645 offending_line_found =
true;
8651 for (
unsigned int l=0;
8652 l<GeometryInfo<dim>::lines_per_cell; ++l)
8653 if (!cell->line(l)->has_children() &&
8657 cell->line(l)->set_user_flag();
8663 if (offending_line_found)
8665 mesh_changed =
true;
8683 cell=triangulation.last(); cell!=triangulation.end(); --cell)
8685 if (cell->user_flag_set())
8686 for (
unsigned int line=0; line<GeometryInfo<dim>::lines_per_cell; ++line)
8687 if (cell->line(line)->has_children() &&
8688 (cell->line(line)->child(0)->user_flag_set() ||
8689 cell->line(line)->child(1)->user_flag_set()))
8691 for (
unsigned int c=0; c<cell->n_children(); ++c)
8692 cell->child(c)->clear_coarsen_flag ();
8693 cell->clear_user_flag();
8694 for (
unsigned int l=0; l<GeometryInfo<dim>::lines_per_cell; ++
l)
8699 cell->line(l)->set_user_flag();
8700 mesh_changed =
true;
8705 while (mesh_changed ==
true);
8716 template <
int dim,
int spacedim>
8727 for (
unsigned int n=0; n<GeometryInfo<dim>::faces_per_cell; ++n)
8736 const unsigned int n_subfaces
8739 if (n_subfaces == 0 || cell->at_boundary(n))
8741 for (
unsigned int c=0; c<n_subfaces; ++c)
8745 child_cell_on_face(ref_case,
8749 child_neighbor = child->neighbor(n);
8750 if (!child->neighbor_is_coarser(n))
8763 if ((child_neighbor->has_children() &&
8764 !child_neighbor->user_flag_set())||
8773 child_neighbor->refine_flag_set())
8783 template <
int dim,
int spacedim>
8785 get_default_flat_manifold()
8788 return flat_manifold;
8794 template <
int dim,
int spacedim>
8800 template <
int dim,
int spacedim>
8806 template <
int dim,
int spacedim>
8809 const bool check_for_distorted_cells)
8811 smooth_grid(smooth_grid),
8812 anisotropic_refinement(false),
8813 check_for_distorted_cells(check_for_distorted_cells)
8830 template <
int dim,
int spacedim>
8835 smooth_grid(tria.smooth_grid),
8836 periodic_face_pairs_level_0(std::move(tria.periodic_face_pairs_level_0)),
8837 periodic_face_map(std::move(tria.periodic_face_map)),
8838 levels(std::move(tria.levels)),
8839 faces(std::move(tria.faces)),
8840 vertices(std::move(tria.vertices)),
8841 vertices_used(std::move(tria.vertices_used)),
8842 manifold(std::move(tria.manifold)),
8843 anisotropic_refinement(tria.anisotropic_refinement),
8844 check_for_distorted_cells(tria.check_for_distorted_cells),
8845 number_cache(std::move(tria.number_cache)),
8846 vertex_to_boundary_id_map_1d(std::move(tria.vertex_to_boundary_id_map_1d)),
8847 vertex_to_manifold_id_map_1d(std::move(tria.vertex_to_manifold_id_map_1d))
8853 template <
int dim,
int spacedim>
8859 smooth_grid = tria.smooth_grid;
8860 periodic_face_pairs_level_0 = std::move(tria.periodic_face_pairs_level_0);
8861 periodic_face_map = std::move(tria.periodic_face_map);
8862 levels = std::move(tria.levels);
8863 faces = std::move(tria.faces);
8864 vertices = std::move(tria.vertices);
8865 vertices_used = std::move(tria.vertices_used);
8866 manifold = std::move(tria.manifold);
8867 anisotropic_refinement = tria.anisotropic_refinement;
8868 number_cache = tria.number_cache;
8869 vertex_to_boundary_id_map_1d = std::move(tria.vertex_to_boundary_id_map_1d);
8870 vertex_to_manifold_id_map_1d = std::move(tria.vertex_to_manifold_id_map_1d);
8879 template <
int dim,
int spacedim>
8898 (vertex_to_boundary_id_map_1d ==
nullptr),
8905 (vertex_to_manifold_id_map_1d ==
nullptr),
8911 template <
int dim,
int spacedim>
8918 clear_despite_subscriptions();
8919 periodic_face_pairs_level_0.clear();
8920 periodic_face_map.clear();
8925 template <
int dim,
int spacedim>
8930 ExcTriangulationNotEmpty (vertices.size(), levels.size()));
8931 smooth_grid=mesh_smoothing;
8936 template <
int dim,
int spacedim>
8945 template <
int dim,
int spacedim>
8950 set_manifold(m_number, boundary_object);
8953 template <
int dim,
int spacedim>
8961 manifold[m_number] = manifold_object.
clone();
8965 template <
int dim,
int spacedim>
8969 set_manifold(m_number);
8973 template <
int dim,
int spacedim>
8977 reset_manifold(m_number);
8981 template<
int dim,
int spacedim>
8988 manifold.erase(m_number);
8992 template<
int dim,
int spacedim>
8999 template <
int dim,
int spacedim>
9004 ExcMessage(
"Error: set_all_manifold_ids() can not be called on an empty Triangulation."));
9007 cell=this->begin_active(), endc=this->end();
9009 for (; cell != endc; ++cell)
9010 cell->set_all_manifold_ids(m_number);
9014 template <
int dim,
int spacedim>
9019 ExcMessage(
"Error: set_all_manifold_ids_on_boundary() can not be called on an empty Triangulation."));
9022 cell=this->begin_active(), endc=this->end();
9024 for (; cell != endc; ++cell)
9025 for (
unsigned int f=0; f<GeometryInfo<dim>::faces_per_cell; ++f)
9026 if (cell->face(f)->at_boundary())
9027 cell->face(f)->set_all_manifold_ids(m_number);
9031 template <
int dim,
int spacedim>
9037 ExcMessage(
"Error: set_all_manifold_ids_on_boundary() can not be called on an empty Triangulation."));
9039 bool boundary_found =
false;
9041 cell=this->begin_active(), endc=this->end();
9043 for (; cell != endc; ++cell)
9046 for (
unsigned int f=0; f<GeometryInfo<dim>::faces_per_cell; ++f)
9047 if (cell->face(f)->at_boundary() && cell->face(f)->boundary_id()==b_id)
9049 boundary_found =
true;
9050 cell->face(f)->set_manifold_id(m_number);
9055 for (
unsigned int e=0; e<GeometryInfo<dim>::lines_per_cell; ++e)
9056 if (cell->line(e)->at_boundary() && cell->line(e)->boundary_id()==b_id)
9058 boundary_found =
true;
9059 cell->line(e)->set_manifold_id(m_number);
9063 (void)boundary_found;
9064 Assert(boundary_found, ExcBoundaryIdNotFound(b_id));
9068 template <
int dim,
int spacedim>
9075 ExcMessage(
"You tried to get a Boundary, but I only have a Manifold."));
9081 template <
int dim,
int spacedim>
9087 const auto it = manifold.find(m_number);
9089 if (it != manifold.end())
9092 return *(it->second);
9097 return internal::TriangulationImplementation::get_default_flat_manifold<dim,spacedim>();
9103 template <
int dim,
int spacedim>
9104 std::vector<types::boundary_id>
9111 std::vector<types::boundary_id> boundary_ids;
9112 for (std::map<unsigned int, types::boundary_id>::const_iterator
9113 p = vertex_to_boundary_id_map_1d->begin();
9114 p != vertex_to_boundary_id_map_1d->end();
9116 boundary_ids.push_back (p->second);
9118 return boundary_ids;
9122 std::set<types::boundary_id> b_ids;
9124 for (; cell!=end(); ++cell)
9125 for (
unsigned int face=0; face<GeometryInfo<dim>::faces_per_cell; ++face)
9126 if (cell->at_boundary(face))
9127 b_ids.insert(cell->face(face)->boundary_id());
9128 std::vector<types::boundary_id> boundary_ids(b_ids.begin(), b_ids.end());
9129 return boundary_ids;
9135 template <
int dim,
int spacedim>
9136 std::vector<types::manifold_id>
9139 std::set<types::manifold_id> m_ids;
9141 for (; cell!=end(); ++cell)
9143 m_ids.insert(cell->manifold_id());
9145 for (
unsigned int face=0; face<GeometryInfo<dim>::faces_per_cell; ++face)
9146 if (cell->at_boundary(face))
9147 m_ids.insert(cell->face(face)->manifold_id());
9149 std::vector<types::manifold_id> manifold_indicators(m_ids.begin(), m_ids.end());
9150 return manifold_indicators;
9156 template <
int dim,
int spacedim>
9161 Assert ((vertices.size() == 0) &&
9162 (levels.size () == 0) &&
9164 ExcTriangulationNotEmpty(vertices.size(), levels.size()));
9166 (other_tria.
vertices.size() != 0) &&
9167 (dim == 1 || other_tria.
faces !=
nullptr),
9168 ExcMessage(
"When calling Triangulation::copy_triangulation(), " 9169 "the target triangulation must be empty but the source " 9170 "triangulation (the argument to this function) must contain " 9171 "something. Here, it seems like the source does not " 9172 "contain anything at all."));
9182 faces = std_cxx14::make_unique<internal::TriangulationImplementation::TriaFaces<dim>> (*other_tria.
faces);
9184 auto bdry_iterator = other_tria.
manifold.begin();
9185 for (; bdry_iterator != other_tria.
manifold.end() ; ++bdry_iterator)
9186 manifold[bdry_iterator->first] = bdry_iterator->second->clone();
9189 levels.reserve (other_tria.
levels.size());
9190 for (
unsigned int level=0; level<other_tria.
levels.size(); ++level)
9192 (*other_tria.
levels[level]));
9198 vertex_to_boundary_id_map_1d = std_cxx14::make_unique<std::map<unsigned int, types::boundary_id>>
9201 vertex_to_manifold_id_map_1d = std_cxx14::make_unique<std::map<unsigned int, types::manifold_id>>
9206 other_tria.
signals.copy (*
this);
9217 template <
int dim,
int spacedim>
9224 std::vector<CellData<dim> > reordered_cells (cells);
9228 reorder_compatibility (reordered_cells, reordered_subcelldata);
9232 create_triangulation(v, reordered_cells, reordered_subcelldata);
9237 template <
int dim,
int spacedim>
9244 Assert ((vertices.size() == 0) &&
9245 (levels.size () == 0) &&
9247 ExcTriangulationNotEmpty(vertices.size(), levels.size()));
9261 clear_despite_subscriptions();
9269 reset_active_cell_indices ();
9274 if (check_for_distorted_cells ==
true)
9276 DistortedCellList distorted_cells = collect_distorted_coarse_cells (*
this);
9280 AssertThrow (distorted_cells.distorted_cells.size() == 0,
9321 bool values [][2] = {{
false,
true},
9324 for (
unsigned int i=0; i< GeometryInfo< dim >::faces_per_cell; ++i)
9325 for (
unsigned int j=0; j< GeometryInfo< dim >::faces_per_cell; ++j)
9326 correct(i,j) = ( values[i][j]);
9331 bool values [][4]= {{
false,
true,
true,
false},
9332 {
true,
false,
false,
true },
9333 {
true,
false,
false,
true },
9334 {
false,
true,
true,
false}
9336 for (
unsigned int i=0; i< GeometryInfo< dim >::faces_per_cell; ++i)
9337 for (
unsigned int j=0; j< GeometryInfo< dim >::faces_per_cell; ++j)
9338 correct(i,j) = ( values[i][j]);
9346 std::list<active_cell_iterator> this_round, next_round;
9349 this_round.push_back (begin_active());
9350 begin_active()->set_direction_flag (
true);
9351 begin_active()->set_user_flag ();
9353 while (this_round.size() > 0)
9355 for (
typename std::list<active_cell_iterator>::iterator cell = this_round.begin();
9356 cell != this_round.end(); ++cell)
9358 for (
unsigned int i = 0; i < GeometryInfo< dim >::faces_per_cell; ++i)
9360 if ( !((*cell)->face(i)->at_boundary()) )
9362 neighbor = (*cell)->neighbor(i);
9364 unsigned int cf = (*cell)->face_index(i);
9366 while (neighbor->face_index(j) != cf)
9373 if (neighbor->user_flag_set())
9377 Assert(!(correct(i,j) ^ (neighbor->direction_flag() == (*cell)->direction_flag() )),
9378 ExcNonOrientableTriangulation());
9382 next_round.push_back (neighbor);
9383 neighbor->set_user_flag ();
9384 if ( (correct(i,j) ^ ( neighbor->direction_flag() == (*cell)->direction_flag() ) ) )
9385 neighbor->set_direction_flag(!neighbor->direction_flag());
9395 if (next_round.size() == 0)
9397 cell != end(); ++cell)
9398 if (cell->user_flag_set() ==
false)
9400 next_round.push_back (cell);
9401 cell->set_direction_flag (
true);
9402 cell->set_user_flag ();
9406 this_round = next_round;
9418 template <
int dim,
int spacedim>
9425 cell != end(); ++cell)
9426 cell->set_direction_flag (!cell->direction_flag());
9431 template <
int dim,
int spacedim>
9434 Assert(n_cells()>0,
ExcMessage(
"Error: An empty Triangulation can not be refined."));
9438 for (; cell != endc; ++cell)
9440 cell->clear_coarsen_flag();
9441 cell->set_refine_flag ();
9447 template <
int dim,
int spacedim>
9450 for (
unsigned int i=0; i<times; ++i)
9452 set_all_refine_flags();
9453 execute_coarsening_and_refinement ();
9463 template <
int dim,
int spacedim>
9466 v.resize (dim*n_active_cells(),
false);
9467 std::vector<bool>::iterator i = v.begin();
9470 for (; cell!=endc; ++cell)
9471 for (
unsigned int j=0; j<dim; ++j,++i)
9472 if (cell->refine_flag_set() & (1<<j) )
9480 template <
int dim,
int spacedim>
9483 std::vector<bool> v;
9484 save_refine_flags (v);
9485 write_bool_vector (mn_tria_refine_flags_begin, v, mn_tria_refine_flags_end,
9491 template <
int dim,
int spacedim>
9494 std::vector<bool> v;
9495 read_bool_vector (mn_tria_refine_flags_begin, v, mn_tria_refine_flags_end,
9497 load_refine_flags (v);
9502 template <
int dim,
int spacedim>
9505 AssertThrow (v.size() == dim*n_active_cells(), ExcGridReadError());
9509 std::vector<bool>::const_iterator i = v.begin();
9510 for (; cell!=endc; ++cell)
9512 unsigned int ref_case=0;
9514 for (
unsigned int j=0; j<dim; ++j, ++i)
9518 ExcGridReadError());
9522 cell->clear_refine_flag();
9530 template <
int dim,
int spacedim>
9533 v.resize (n_active_cells(),
false);
9534 std::vector<bool>::iterator i = v.begin();
9537 for (; cell!=endc; ++cell, ++i)
9538 *i = cell->coarsen_flag_set();
9545 template <
int dim,
int spacedim>
9548 std::vector<bool> v;
9549 save_coarsen_flags (v);
9550 write_bool_vector (mn_tria_coarsen_flags_begin, v, mn_tria_coarsen_flags_end,
9556 template <
int dim,
int spacedim>
9559 std::vector<bool> v;
9560 read_bool_vector (mn_tria_coarsen_flags_begin, v, mn_tria_coarsen_flags_end,
9562 load_coarsen_flags (v);
9567 template <
int dim,
int spacedim>
9570 Assert (v.size() == n_active_cells(), ExcGridReadError());
9574 std::vector<bool>::const_iterator i = v.begin();
9575 for (; cell!=endc; ++cell, ++i)
9577 cell->set_coarsen_flag();
9579 cell->clear_coarsen_flag();
9585 template <
int dim,
int spacedim>
9588 return anisotropic_refinement;
9602 for (
unsigned int level=0; level<levels.size(); ++level)
9603 levels[level]->cells.clear_user_data();
9616 faces->
lines.clear_user_data();
9622 faces->
lines.clear_user_data();
9628 template <
int dim,
int spacedim>
9632 ::clear_user_data (levels);
9633 ::clear_user_data (faces.get());
9643 for (
unsigned int level=0; level<levels.size(); ++level)
9644 levels[level]->cells.clear_user_flags();
9651 faces->lines.clear_user_flags();
9656 template <
int dim,
int spacedim>
9659 ::clear_user_flags_line (levels, faces.get());
9675 for (
unsigned int level=0; level<levels.size(); ++level)
9676 levels[level]->cells.clear_user_flags();
9683 faces->quads.clear_user_flags();
9688 template <
int dim,
int spacedim>
9691 ::clear_user_flags_quad (levels, faces.get());
9714 for (
unsigned int level=0; level<levels.size(); ++level)
9715 levels[level]->cells.clear_user_flags();
9720 template <
int dim,
int spacedim>
9723 ::clear_user_flags_hex (levels, faces.get());
9728 template <
int dim,
int spacedim>
9731 clear_user_flags_line ();
9732 clear_user_flags_quad ();
9733 clear_user_flags_hex ();
9738 template <
int dim,
int spacedim>
9741 save_user_flags_line (out);
9744 save_user_flags_quad (out);
9747 save_user_flags_hex (out);
9755 template <
int dim,
int spacedim>
9762 std::vector<bool> tmp;
9764 save_user_flags_line (tmp);
9765 v.insert (v.end(), tmp.begin(), tmp.end());
9769 save_user_flags_quad (tmp);
9770 v.insert (v.end(), tmp.begin(), tmp.end());
9775 save_user_flags_hex (tmp);
9776 v.insert (v.end(), tmp.begin(), tmp.end());
9785 template <
int dim,
int spacedim>
9788 load_user_flags_line (in);
9791 load_user_flags_quad (in);
9794 load_user_flags_hex (in);
9802 template <
int dim,
int spacedim>
9806 std::vector<bool> tmp;
9810 tmp.insert (tmp.end(),
9811 v.begin(), v.begin()+n_lines());
9813 load_user_flags_line (tmp);
9818 tmp.insert (tmp.end(),
9819 v.begin()+n_lines(), v.begin()+n_lines()+n_quads());
9820 load_user_flags_quad (tmp);
9826 tmp.insert (tmp.end(),
9827 v.begin()+n_lines()+n_quads(), v.begin()+n_lines()+n_quads()+n_hexs());
9828 load_user_flags_hex (tmp);
9837 template <
int dim,
int spacedim>
9840 v.resize (n_lines(),
false);
9841 std::vector<bool>::iterator i = v.begin();
9842 line_iterator line = begin_line(),
9844 for (; line!=endl; ++line, ++i)
9845 *i = line->user_flag_set();
9852 template <
int dim,
int spacedim>
9855 std::vector<bool> v;
9856 save_user_flags_line (v);
9857 write_bool_vector (mn_tria_line_user_flags_begin, v, mn_tria_line_user_flags_end,
9863 template <
int dim,
int spacedim>
9866 std::vector<bool> v;
9867 read_bool_vector (mn_tria_line_user_flags_begin, v, mn_tria_line_user_flags_end,
9869 load_user_flags_line (v);
9874 template <
int dim,
int spacedim>
9877 Assert (v.size() == n_lines(), ExcGridReadError());
9879 line_iterator line = begin_line(),
9881 std::vector<bool>::const_iterator i = v.begin();
9882 for (; line!=endl; ++line, ++i)
9884 line->set_user_flag();
9886 line->clear_user_flag();
9894 template <
typename Iterator>
9895 bool get_user_flag (
const Iterator &i)
9897 return i->user_flag_set();
9902 template <
int structdim,
int dim,
int spacedim>
9911 template <
typename Iterator>
9912 void set_user_flag (
const Iterator &i)
9919 template <
int structdim,
int dim,
int spacedim>
9927 template <
typename Iterator>
9928 void clear_user_flag (
const Iterator &i)
9930 i->clear_user_flag();
9935 template <
int structdim,
int dim,
int spacedim>
9943 template <
int dim,
int spacedim>
9946 v.resize (n_quads(),
false);
9950 std::vector<bool>::iterator i = v.begin();
9951 quad_iterator quad = begin_quad(),
9953 for (; quad!=endq; ++quad, ++i)
9954 *i = get_user_flag (quad);
9962 template <
int dim,
int spacedim>
9965 std::vector<bool> v;
9966 save_user_flags_quad (v);
9967 write_bool_vector (mn_tria_quad_user_flags_begin, v, mn_tria_quad_user_flags_end,
9973 template <
int dim,
int spacedim>
9976 std::vector<bool> v;
9977 read_bool_vector (mn_tria_quad_user_flags_begin, v, mn_tria_quad_user_flags_end,
9979 load_user_flags_quad (v);
9984 template <
int dim,
int spacedim>
9987 Assert (v.size() == n_quads(), ExcGridReadError());
9991 quad_iterator quad = begin_quad(),
9993 std::vector<bool>::const_iterator i = v.begin();
9994 for (; quad!=endq; ++quad, ++i)
9996 set_user_flag(quad);
9998 clear_user_flag(quad);
10006 template <
int dim,
int spacedim>
10009 v.resize (n_hexs(),
false);
10013 std::vector<bool>::iterator i = v.begin();
10014 hex_iterator hex = begin_hex(),
10016 for (; hex!=endh; ++hex, ++i)
10017 *i = get_user_flag (hex);
10025 template <
int dim,
int spacedim>
10028 std::vector<bool> v;
10029 save_user_flags_hex (v);
10030 write_bool_vector (mn_tria_hex_user_flags_begin, v, mn_tria_hex_user_flags_end,
10036 template <
int dim,
int spacedim>
10039 std::vector<bool> v;
10040 read_bool_vector (mn_tria_hex_user_flags_begin, v, mn_tria_hex_user_flags_end,
10042 load_user_flags_hex (v);
10047 template <
int dim,
int spacedim>
10050 Assert (v.size() == n_hexs(), ExcGridReadError());
10054 hex_iterator hex = begin_hex(),
10056 std::vector<bool>::const_iterator i = v.begin();
10057 for (; hex!=endh; ++hex, ++i)
10059 set_user_flag(hex);
10061 clear_user_flag(hex);
10069 template <
int dim,
int spacedim>
10076 std::vector<unsigned int> tmp;
10078 save_user_indices_line (tmp);
10079 v.insert (v.end(), tmp.begin(), tmp.end());
10083 save_user_indices_quad (tmp);
10084 v.insert (v.end(), tmp.begin(), tmp.end());
10089 save_user_indices_hex (tmp);
10090 v.insert (v.end(), tmp.begin(), tmp.end());
10099 template <
int dim,
int spacedim>
10103 std::vector<unsigned int> tmp;
10107 tmp.insert (tmp.end(),
10108 v.begin(), v.begin()+n_lines());
10110 load_user_indices_line (tmp);
10115 tmp.insert (tmp.end(),
10116 v.begin()+n_lines(), v.begin()+n_lines()+n_quads());
10117 load_user_indices_quad (tmp);
10123 tmp.insert (tmp.end(),
10124 v.begin()+n_lines()+n_quads(), v.begin()+n_lines()+n_quads()+n_hexs());
10125 load_user_indices_hex (tmp);
10136 template <
typename Iterator>
10137 unsigned int get_user_index (
const Iterator &i)
10139 return i->user_index();
10144 template <
int structdim,
int dim,
int spacedim>
10153 template <
typename Iterator>
10154 void set_user_index (
const Iterator &i,
10155 const unsigned int x)
10157 i->set_user_index(x);
10162 template <
int structdim,
int dim,
int spacedim>
10164 const unsigned int)
10171 template <
int dim,
int spacedim>
10174 v.resize (n_lines(), 0);
10175 std::vector<unsigned int>::iterator i = v.begin();
10176 line_iterator line = begin_line(),
10178 for (; line!=endl; ++line, ++i)
10179 *i = line->user_index();
10184 template <
int dim,
int spacedim>
10187 Assert (v.size() == n_lines(), ExcGridReadError());
10189 line_iterator line = begin_line(),
10191 std::vector<unsigned int>::const_iterator i = v.begin();
10192 for (; line!=endl; ++line, ++i)
10193 line->set_user_index(*i);
10197 template <
int dim,
int spacedim>
10200 v.resize (n_quads(), 0);
10204 std::vector<unsigned int>::iterator i = v.begin();
10205 quad_iterator quad = begin_quad(),
10207 for (; quad!=endq; ++quad, ++i)
10208 *i = get_user_index(quad);
10214 template <
int dim,
int spacedim>
10217 Assert (v.size() == n_quads(), ExcGridReadError());
10221 quad_iterator quad = begin_quad(),
10223 std::vector<unsigned int>::const_iterator i = v.begin();
10224 for (; quad!=endq; ++quad, ++i)
10225 set_user_index(quad, *i);
10230 template <
int dim,
int spacedim>
10233 v.resize (n_hexs(), 0);
10237 std::vector<unsigned int>::iterator i = v.begin();
10238 hex_iterator hex = begin_hex(),
10240 for (; hex!=endh; ++hex, ++i)
10241 *i = get_user_index(hex);
10247 template <
int dim,
int spacedim>
10250 Assert (v.size() == n_hexs(), ExcGridReadError());
10254 hex_iterator hex = begin_hex(),
10256 std::vector<unsigned int>::const_iterator i = v.begin();
10257 for (; hex!=endh; ++hex, ++i)
10258 set_user_index(hex, *i);
10269 template <
typename Iterator>
10270 void *get_user_pointer (
const Iterator &i)
10272 return i->user_pointer();
10277 template <
int structdim,
int dim,
int spacedim>
10286 template <
typename Iterator>
10287 void set_user_pointer (
const Iterator &i,
10290 i->set_user_pointer(x);
10295 template <
int structdim,
int dim,
int spacedim>
10304 template <
int dim,
int spacedim>
10311 std::vector<void *> tmp;
10313 save_user_pointers_line (tmp);
10314 v.insert (v.end(), tmp.begin(), tmp.end());
10318 save_user_pointers_quad (tmp);
10319 v.insert (v.end(), tmp.begin(), tmp.end());
10324 save_user_pointers_hex (tmp);
10325 v.insert (v.end(), tmp.begin(), tmp.end());
10334 template <
int dim,
int spacedim>
10338 std::vector<void *> tmp;
10342 tmp.insert (tmp.end(),
10343 v.begin(), v.begin()+n_lines());
10345 load_user_pointers_line (tmp);
10350 tmp.insert (tmp.end(),
10351 v.begin()+n_lines(), v.begin()+n_lines()+n_quads());
10352 load_user_pointers_quad (tmp);
10358 tmp.insert (tmp.end(),
10359 v.begin()+n_lines()+n_quads(), v.begin()+n_lines()+n_quads()+n_hexs());
10360 load_user_pointers_hex (tmp);
10369 template <
int dim,
int spacedim>
10372 v.resize (n_lines(),
nullptr);
10373 std::vector<void *>::iterator i = v.begin();
10374 line_iterator line = begin_line(),
10376 for (; line!=endl; ++line, ++i)
10377 *i = line->user_pointer();
10382 template <
int dim,
int spacedim>
10385 Assert (v.size() == n_lines(), ExcGridReadError());
10387 line_iterator line = begin_line(),
10389 std::vector<void *>::const_iterator i = v.begin();
10390 for (; line!=endl; ++line, ++i)
10391 line->set_user_pointer(*i);
10396 template <
int dim,
int spacedim>
10399 v.resize (n_quads(),
nullptr);
10403 std::vector<void *>::iterator i = v.begin();
10404 quad_iterator quad = begin_quad(),
10406 for (; quad!=endq; ++quad, ++i)
10407 *i = get_user_pointer(quad);
10413 template <
int dim,
int spacedim>
10416 Assert (v.size() == n_quads(), ExcGridReadError());
10420 quad_iterator quad = begin_quad(),
10422 std::vector<void *>::const_iterator i = v.begin();
10423 for (; quad!=endq; ++quad, ++i)
10424 set_user_pointer(quad, *i);
10429 template <
int dim,
int spacedim>
10432 v.resize (n_hexs(),
nullptr);
10436 std::vector<void *>::iterator i = v.begin();
10437 hex_iterator hex = begin_hex(),
10439 for (; hex!=endh; ++hex, ++i)
10440 *i = get_user_pointer(hex);
10446 template <
int dim,
int spacedim>
10449 Assert (v.size() == n_hexs(), ExcGridReadError());
10453 hex_iterator hex = begin_hex(),
10455 std::vector<void *>::const_iterator i = v.begin();
10456 for (; hex!=endh; ++hex, ++i)
10457 set_user_pointer(hex, *i);
10466 template <
int dim,
int spacedim>
10473 return begin_raw_line (level);
10475 return begin_raw_quad (level);
10477 return begin_raw_hex (level);
10486 template <
int dim,
int spacedim>
10493 return begin_line (level);
10495 return begin_quad (level);
10497 return begin_hex (level);
10506 template <
int dim,
int spacedim>
10513 return begin_active_line (level);
10515 return begin_active_quad (level);
10517 return begin_active_hex (level);
10526 template <
int dim,
int spacedim>
10530 const unsigned int level = levels.size()-1;
10532 Assert (level<n_global_levels() || level<levels.size(), ExcInvalidLevel(level));
10533 if (levels[level]->cells.cells.size() ==0)
10540 levels[level]->cells.cells.size()-1);
10543 if (ri->used()==
true)
10546 if (ri->used()==
true)
10553 template <
int dim,
int spacedim>
10563 if (cell->active()==
true)
10566 if (cell->active()==
true)
10574 template <
int dim,
int spacedim>
10585 template <
int dim,
int spacedim>
10589 Assert (level<n_global_levels(), ExcInvalidLevel(level));
10590 if (level < levels.size()-1)
10591 return begin_raw (level+1);
10597 template <
int dim,
int spacedim>
10601 if (level < levels.size()-1)
10602 return begin (level+1);
10603 Assert (level<n_global_levels() || level<levels.size(), ExcInvalidLevel(level));
10608 template <
int dim,
int spacedim>
10612 Assert (level<n_global_levels() || level < levels.size(), ExcInvalidLevel(level));
10613 return (level >= levels.size()-1 ?
10615 begin_active (level+1));
10620 template <
int dim,
int spacedim>
10630 template <
int dim,
int spacedim>
10636 (begin_active(), end());
10641 template <
int dim,
int spacedim>
10647 (begin(level), end(level));
10652 template <
int dim,
int spacedim>
10658 (begin_active(level), end_active(level));
10665 template <
int dim,
int spacedim>
10675 return begin_line ();
10677 return begin_quad ();
10686 template <
int dim,
int spacedim>
10696 return begin_active_line ();
10698 return begin_active_quad ();
10707 template <
int dim,
int spacedim>
10717 return end_line ();
10719 return end_quad ();
10730 template <
int dim,
int spacedim>
10749 while (i->used() ==
false)
10758 template <
int dim,
int spacedim>
10762 return begin_vertex();
10767 template <
int dim,
int spacedim>
10789 template <
int dim,
int spacedim>
10790 typename Triangulation<dim, spacedim>::raw_line_iterator
10796 Assert (level<n_global_levels() || level<levels.size(), ExcInvalidLevel(level));
10798 if (level >= levels.size() || levels[level]->cells.cells.size() == 0)
10814 template <
int dim,
int spacedim>
10815 typename Triangulation<dim, spacedim>::line_iterator
10819 raw_line_iterator ri = begin_raw_line (level);
10822 while (ri->used() ==
false)
10830 template <
int dim,
int spacedim>
10831 typename Triangulation<dim, spacedim>::active_line_iterator
10835 line_iterator i = begin_line (level);
10838 while (i->has_children())
10846 template <
int dim,
int spacedim>
10847 typename Triangulation<dim, spacedim>::line_iterator
10860 template <
int dim,
int spacedim>
10861 typename Triangulation<dim,spacedim>::raw_quad_iterator
10868 return raw_hex_iterator();
10871 Assert (level<n_global_levels() || level<levels.size(), ExcInvalidLevel(level));
10873 if (level >= levels.size() || levels[level]->cells.cells.size() == 0)
10893 return raw_hex_iterator();
10899 template <
int dim,
int spacedim>
10900 typename Triangulation<dim,spacedim>::quad_iterator
10904 raw_quad_iterator ri = begin_raw_quad (level);
10907 while (ri->used() ==
false)
10915 template <
int dim,
int spacedim>
10916 typename Triangulation<dim,spacedim>::active_quad_iterator
10920 quad_iterator i = begin_quad (level);
10923 while (i->has_children())
10931 template <
int dim,
int spacedim>
10932 typename Triangulation<dim,spacedim>::quad_iterator
10944 template <
int dim,
int spacedim>
10945 typename Triangulation<dim,spacedim>::raw_hex_iterator
10953 return raw_hex_iterator();
10956 Assert (level<n_global_levels() || level<levels.size(), ExcInvalidLevel(level));
10958 if (level >= levels.size() || levels[level]->cells.cells.size() == 0)
10968 return raw_hex_iterator();
10974 template <
int dim,
int spacedim>
10975 typename Triangulation<dim,spacedim>::hex_iterator
10979 raw_hex_iterator ri = begin_raw_hex (level);
10982 while (ri->used() ==
false)
10990 template <
int dim,
int spacedim>
10991 typename Triangulation<dim, spacedim>::active_hex_iterator
10995 hex_iterator i = begin_hex (level);
10998 while (i->has_children())
11006 template <
int dim,
int spacedim>
11007 typename Triangulation<dim, spacedim>::hex_iterator
11023 namespace TriangulationImplementation
11076 template <
int dim,
int spacedim>
11079 return internal::TriangulationImplementation::n_cells (number_cache);
11083 template <
int dim,
int spacedim>
11086 return internal::TriangulationImplementation::n_active_cells (number_cache);
11089 template <
int dim,
int spacedim>
11092 return n_active_cells();
11097 template <
int dim,
int spacedim>
11115 template <
int dim,
int spacedim>
11121 return n_raw_lines();
11123 return n_raw_quads();
11131 template <
int dim,
int spacedim>
11139 return n_active_lines();
11141 return n_active_quads();
11149 template <
int dim,
int spacedim>
11155 return n_raw_lines(level);
11157 return n_raw_quads(level);
11159 return n_raw_hexs(level);
11168 template <
int dim,
int spacedim>
11174 return n_lines(level);
11176 return n_quads(level);
11178 return n_hexs(level);
11187 template <
int dim,
int spacedim>
11193 return n_active_lines(level);
11195 return n_active_quads(level);
11197 return n_active_hexs(level);
11205 template <
int dim,
int spacedim>
11208 for (
unsigned int lvl = 0; lvl<n_global_levels()-1; lvl++)
11209 if (n_active_cells(lvl) != 0)
11216 template <
int dim,
int spacedim>
11219 return number_cache.n_lines;
11228 return levels[level]->cells.cells.size();
11245 return levels[level]->cells.cells.size();
11261 return levels[level]->cells.cells.size();
11273 template <
int dim,
int spacedim>
11281 template <
int dim,
int spacedim>
11284 return faces->
lines.cells.size();
11288 template <
int dim,
int spacedim>
11291 Assert (level < number_cache.n_lines_level.size(),
11292 ExcIndexRange (level, 0, number_cache.n_lines_level.size()));
11294 return number_cache.n_lines_level[level];
11298 template <
int dim,
int spacedim>
11301 return number_cache.n_active_lines;
11305 template <
int dim,
int spacedim>
11308 Assert (level < number_cache.n_lines_level.size(),
11309 ExcIndexRange (level, 0, number_cache.n_lines_level.size()));
11312 return number_cache.n_active_lines_level[level];
11445 template <
int dim,
int spacedim>
11448 return number_cache.n_quads;
11452 template <
int dim,
int spacedim>
11456 Assert (level < number_cache.n_quads_level.size(),
11457 ExcIndexRange (level, 0, number_cache.n_quads_level.size()));
11458 return number_cache.n_quads_level[level];
11467 return levels[level]->cells.cells.size();
11476 return levels[level]->cells.cells.size();
11491 template <
int dim,
int spacedim>
11508 template <
int dim,
int spacedim>
11511 return number_cache.n_active_quads;
11515 template <
int dim,
int spacedim>
11518 Assert (level < number_cache.n_quads_level.size(),
11519 ExcIndexRange (level, 0, number_cache.n_quads_level.size()));
11522 return number_cache.n_active_quads_level[level];
11526 template <
int dim,
int spacedim>
11534 template <
int dim,
int spacedim>
11542 template <
int dim,
int spacedim>
11549 template <
int dim,
int spacedim>
11557 template <
int dim,
int spacedim>
11567 return number_cache.n_hexes;
11575 Assert (level < number_cache.n_hexes_level.size(),
11576 ExcIndexRange (level, 0, number_cache.n_hexes_level.size()));
11578 return number_cache.n_hexes_level[level];
11587 return levels[level]->cells.cells.size();
11594 return number_cache.n_active_hexes;
11602 Assert (level < number_cache.n_hexes_level.size(),
11603 ExcIndexRange (level, 0, number_cache.n_hexes_level.size()));
11605 return number_cache.n_active_hexes_level[level];
11610 template <
int dim,
int spacedim>
11614 return std::count_if (vertices_used.begin(), vertices_used.end(),
11615 std::bind (std::equal_to<bool>(), std::placeholders::_1,
true));
11620 template <
int dim,
int spacedim>
11621 const std::vector<bool> &
11624 return vertices_used;
11652 template <
int dim,
int spacedim>
11659 unsigned int max_vertex_index = 0;
11660 for (; cell!=endc; ++cell)
11661 for (
unsigned int vertex=0; vertex<GeometryInfo<dim>::vertices_per_cell; ++vertex)
11662 if (cell->vertex_index(vertex) > max_vertex_index)
11663 max_vertex_index = cell->vertex_index(vertex);
11669 std::vector<unsigned short int> usage_count (max_vertex_index+1, 0);
11673 for (cell=begin(); cell!=endc; ++cell)
11674 for (
unsigned int vertex=0; vertex<GeometryInfo<dim>::vertices_per_cell; ++vertex)
11675 ++usage_count[cell->vertex_index(vertex)];
11678 static_cast<unsigned int>(*std::max_element (usage_count.begin(),
11679 usage_count.end())));
11684 template <
int dim,
int spacedim>
11693 template <
int dim,
int spacedim>
11702 template <
int dim,
int spacedim>
11711 template <
int dim,
int spacedim>
11715 periodicity_vector)
11717 periodic_face_pairs_level_0.insert(periodic_face_pairs_level_0.end(),
11718 periodicity_vector.begin(),
11719 periodicity_vector.end());
11722 update_periodic_face_map();
11727 template <
int dim,
int spacedim>
11728 const typename std::map<std::pair<typename Triangulation<dim, spacedim>::cell_iterator,
unsigned int>,
11729 std::pair<std::pair<typename Triangulation<dim, spacedim>::cell_iterator,
unsigned int>, std::bitset<3> > > &
11732 return periodic_face_map;
11737 template <
int dim,
int spacedim>
11741 prepare_coarsening_and_refinement ();
11746 if (smooth_grid & limit_level_difference_at_vertices)
11747 Assert (satisfies_level1_at_vertex_rule (*
this) ==
true,
11751 signals.pre_refinement();
11753 execute_coarsening();
11755 const DistortedCellList
11756 cells_with_distorted_children = execute_refinement();
11761 if (smooth_grid & limit_level_difference_at_vertices)
11762 Assert (satisfies_level1_at_vertex_rule (*
this) ==
true,
11767 update_neighbors(*
this);
11768 reset_active_cell_indices ();
11771 signals.post_refinement();
11773 AssertThrow (cells_with_distorted_children.distorted_cells.size() == 0,
11774 cells_with_distorted_children);
11776 update_periodic_face_map();
11781 template <
int dim,
int spacedim>
11785 unsigned int active_cell_index = 0;
11787 if ((cell->used() ==
false) || cell->has_children())
11791 cell->set_active_cell_index (active_cell_index);
11792 ++active_cell_index;
11799 template <
int dim,
int spacedim>
11804 periodic_face_map.clear();
11806 typename std::vector<GridTools::PeriodicFacePair<cell_iterator> >::const_iterator it;
11807 for (it=periodic_face_pairs_level_0.begin(); it!=periodic_face_pairs_level_0.end(); ++it)
11809 update_periodic_face_map_recursively<dim, spacedim>
11810 (it->cell[0], it->cell[1], it->face_idx[0], it->face_idx[1],
11811 it->orientation, periodic_face_map);
11814 std::bitset<3> inverted_orientation;
11816 bool orientation, flip, rotation;
11817 orientation = it->orientation[0];
11818 rotation = it->orientation[2];
11819 flip = orientation ? rotation ^ it->orientation[1] : it->orientation[1];
11820 inverted_orientation[0] = orientation;
11821 inverted_orientation[1] = flip;
11822 inverted_orientation[2] = rotation;
11824 update_periodic_face_map_recursively<dim, spacedim>
11825 (it->cell[1], it->cell[0], it->face_idx[1], it->face_idx[0],
11826 inverted_orientation, periodic_face_map);
11830 typename std::map<std::pair<cell_iterator, unsigned int>,
11831 std::pair<std::pair<cell_iterator, unsigned int>, std::bitset<3> > >::const_iterator it_test;
11832 for (it_test=periodic_face_map.begin(); it_test!=periodic_face_map.end(); ++it_test)
11836 if (cell_1->level() == cell_2->level())
11840 Assert (periodic_face_map[it_test->second.first].first == it_test->first,
11848 template <
int dim,
int spacedim>
11856 vertices_used.clear ();
11864 template <
int dim,
int spacedim>
11868 const DistortedCellList
11869 cells_with_distorted_children
11881 for (
unsigned int level=0; level<levels.size(); ++level)
11882 levels[level]->cells.monitor_memory (dim);
11890 while (cell != endc)
11894 return cells_with_distorted_children;
11899 template <
int dim,
int spacedim>
11906 std::vector<unsigned int> line_cell_count = count_cells_bounded_by_line (*
this);
11907 std::vector<unsigned int> quad_cell_count = count_cells_bounded_by_quad (*
this);
11919 clear_user_flags ();
11923 for (; cell!=endc; ++cell)
11924 if (!cell->active())
11925 if (cell->child(0)->coarsen_flag_set())
11927 cell->set_user_flag();
11928 for (
unsigned int child=0; child<cell->n_children(); ++child)
11930 Assert (cell->child(child)->coarsen_flag_set(),
11932 cell->child(child)->clear_coarsen_flag();
11947 if (levels.size() >= 2)
11948 for (cell = last(); cell!=endc; --cell)
11949 if (cell->level()<=
static_cast<int>(levels.size()-2) && cell->user_flag_set())
11952 signals.pre_coarsening_on_cell(cell);
11964 for (cell=begin(); cell!=endc; ++cell)
11971 template <
int dim,
int spacedim>
11987 std::vector<bool> previous_coarsen_flags (n_active_cells());
11988 save_coarsen_flags (previous_coarsen_flags);
11990 std::vector<int> vertex_level (vertices.size(), 0);
11992 bool continue_iterating =
true;
11996 if (smooth_grid & limit_level_difference_at_vertices)
11998 Assert(!anisotropic_refinement,
11999 ExcMessage(
"In case of anisotropic refinement the " 12000 "limit_level_difference_at_vertices flag for " 12001 "mesh smoothing must not be set!"));
12005 std::fill (vertex_level.begin(), vertex_level.end(), 0);
12008 for (; cell!=endc; ++cell)
12010 if (cell->refine_flag_set())
12011 for (
unsigned int vertex=0; vertex<GeometryInfo<dim>::vertices_per_cell;
12013 vertex_level[cell->vertex_index(vertex)]
12014 = std::max (vertex_level[cell->vertex_index(vertex)],
12016 else if (!cell->coarsen_flag_set())
12017 for (
unsigned int vertex=0; vertex<GeometryInfo<dim>::vertices_per_cell;
12019 vertex_level[cell->vertex_index(vertex)]
12020 = std::max (vertex_level[cell->vertex_index(vertex)],
12031 for (
unsigned int vertex=0; vertex<GeometryInfo<dim>::vertices_per_cell;
12033 vertex_level[cell->vertex_index(vertex)]
12034 = std::max (vertex_level[cell->vertex_index(vertex)],
12049 for (cell=last_active(); cell != endc; --cell)
12050 if (cell->refine_flag_set() ==
false)
12052 for (
unsigned int vertex=0;
12053 vertex<GeometryInfo<dim>::vertices_per_cell; ++vertex)
12054 if (vertex_level[cell->vertex_index(vertex)] >=
12058 cell->clear_coarsen_flag();
12063 if (vertex_level[cell->vertex_index(vertex)] >
12066 cell->set_refine_flag();
12068 for (
unsigned int v=0; v<GeometryInfo<dim>::vertices_per_cell;
12070 vertex_level[cell->vertex_index(v)]
12071 = std::max (vertex_level[cell->vertex_index(v)],
12094 clear_user_flags ();
12098 end_ac = end_active(0);
12099 for (; acell!=end_ac; ++acell)
12100 acell->clear_coarsen_flag();
12104 for (; cell!=endc; ++cell)
12107 if (cell->active())
12110 const unsigned int n_children=cell->n_children();
12111 unsigned int flagged_children=0;
12112 for (
unsigned int child=0; child<n_children; ++child)
12113 if (cell->child(child)->active() &&
12114 cell->child(child)->coarsen_flag_set())
12116 ++flagged_children;
12118 cell->child(child)->clear_coarsen_flag();
12123 if (flagged_children == n_children)
12124 cell->set_user_flag();
12130 for (cell=begin(); cell!=endc; ++cell)
12151 for (cell=last(); cell!=endc; --cell)
12152 if (cell->user_flag_set())
12155 if (internal::TriangulationImplementation::Implementation::template coarsening_allowed<dim,spacedim>(cell))
12156 for (
unsigned int c=0; c<cell->n_children(); ++c)
12158 Assert (cell->child(c)->refine_flag_set()==
false,
12161 cell->child(c)->set_coarsen_flag();
12166 clear_user_flags ();
12171 std::vector<bool> current_coarsen_flags (n_active_cells());
12172 save_coarsen_flags (current_coarsen_flags);
12174 continue_iterating = (current_coarsen_flags != previous_coarsen_flags);
12175 previous_coarsen_flags = current_coarsen_flags;
12177 while (continue_iterating ==
true);
12187 std::vector<bool> flags_before;
12188 save_coarsen_flags (flags_before);
12191 fix_coarsen_flags ();
12193 std::vector<bool> flags_after;
12194 save_coarsen_flags (flags_after);
12196 return (flags_before != flags_after);
12205 std::vector<bool> flags_before;
12206 save_coarsen_flags (flags_before);
12209 fix_coarsen_flags ();
12211 std::vector<bool> flags_after;
12212 save_coarsen_flags (flags_after);
12214 return (flags_before != flags_after);
12223 std::vector<bool> flags_before;
12224 save_coarsen_flags (flags_before);
12227 fix_coarsen_flags ();
12229 std::vector<bool> flags_after;
12230 save_coarsen_flags (flags_after);
12232 return (flags_before != flags_after);
12245 template <
int dim,
int spacedim>
12247 possibly_do_not_produce_unrefined_islands(
12252 unsigned int n_neighbors=0;
12255 unsigned int count=0;
12256 for (
unsigned int n=0; n<GeometryInfo<dim>::faces_per_cell; ++n)
12262 if (face_will_be_refined_by_neighbor(cell,n))
12269 if (count==n_neighbors ||
12270 (count>=n_neighbors-1 &&
12273 for (
unsigned int c=0; c<cell->n_children(); ++c)
12274 cell->child(c)->clear_coarsen_flag();
12276 for (
unsigned int face=0; face<GeometryInfo<dim>::faces_per_cell; ++face)
12277 if (!cell->at_boundary(face)
12279 ( !cell->neighbor(face)->active() )
12280 && (cell_will_be_coarsened(cell->neighbor(face))) )
12281 possibly_do_not_produce_unrefined_islands<dim,spacedim>( cell->neighbor(face) );
12296 template <
int dim,
int spacedim>
12298 possibly_refine_unrefined_island
12300 const bool allow_anisotropic_smoothing)
12315 if (allow_anisotropic_smoothing ==
false)
12318 unsigned int refined_neighbors = 0,
12319 unrefined_neighbors = 0;
12320 for (
unsigned int face=0; face<GeometryInfo<dim>::faces_per_cell; ++face)
12321 if (!cell->at_boundary(face))
12323 if (face_will_be_refined_by_neighbor(cell,face))
12324 ++refined_neighbors;
12326 ++unrefined_neighbors;
12329 if (unrefined_neighbors < refined_neighbors)
12331 cell->clear_coarsen_flag();
12332 cell->set_refine_flag ();
12337 if (unrefined_neighbors > 0)
12338 for (
unsigned int face=0; face<GeometryInfo<dim>::faces_per_cell; ++face)
12339 if (!cell->at_boundary(face)
12341 (face_will_be_refined_by_neighbor(cell,face) ==
false)
12343 (cell->neighbor(face)->has_children() ==
false)
12345 (cell->neighbor(face)->refine_flag_set() ==
false))
12346 possibly_refine_unrefined_island<dim,spacedim>
12347 (cell->neighbor(face),
12348 allow_anisotropic_smoothing);
12360 for (
unsigned int face_pair=0;
12361 face_pair<GeometryInfo<dim>::faces_per_cell/2; ++face_pair)
12369 for (
unsigned int face_index=0; face_index<2; ++face_index)
12371 unsigned int face=2*face_pair+face_index;
12378 face_will_be_refined_by_neighbor<dim,spacedim>(cell,face,expected_face_ref_case);
12392 directional_cell_refinement_case
12393 = (directional_cell_refinement_case &
12395 expected_face_ref_case,
12397 cell->face_orientation(face),
12398 cell->face_flip(face),
12399 cell->face_rotation(face)));
12405 Assert(directional_cell_refinement_case <
12408 smoothing_cell_refinement_case = smoothing_cell_refinement_case |
12409 directional_cell_refinement_case;
12414 if (smoothing_cell_refinement_case)
12416 cell->clear_coarsen_flag();
12417 cell->set_refine_flag(cell->refine_flag_set() |
12418 smoothing_cell_refinement_case);
12425 template <
int dim,
int spacedim>
12430 std::vector<bool> flags_before[2];
12431 save_coarsen_flags (flags_before[0]);
12432 save_refine_flags (flags_before[1]);
12450 std::vector<bool> flags_before_loop[2] = {flags_before[0],
12529 if (((smooth_grid & coarsest_level_1) ||
12530 (smooth_grid & patch_level_1)) && n_levels()>=2)
12533 cell=begin_active(1),
12534 endc=end_active(1);
12536 for (; cell!=endc; ++cell)
12537 cell->clear_coarsen_flag();
12540 bool mesh_changed_in_this_loop =
false;
12549 if (smooth_grid & do_not_produce_unrefined_islands &&
12550 !(smooth_grid & patch_level_1))
12555 for (cell=begin(); cell!=endc; ++cell)
12559 if (!cell->active() && cell_will_be_coarsened(cell))
12560 possibly_do_not_produce_unrefined_islands<dim,spacedim>(cell);
12584 if (smooth_grid & (eliminate_refined_inner_islands |
12585 eliminate_refined_boundary_islands) &&
12586 !(smooth_grid & patch_level_1))
12591 for (cell=begin(); cell!=endc; ++cell)
12592 if (!cell->active() ||
12594 cell->refine_flag_set() &&
12595 cell->is_locally_owned()))
12602 bool all_children_active =
true;
12603 if (!cell->active())
12604 for (
unsigned int c=0; c<cell->n_children(); ++c)
12605 if (!cell->child(c)->active() ||
12606 cell->child(c)->is_ghost() ||
12607 cell->child(c)->is_artificial())
12609 all_children_active =
false;
12613 if (all_children_active)
12620 unsigned int unrefined_neighbors = 0,
12621 total_neighbors = 0;
12623 for (
unsigned int n=0; n<GeometryInfo<dim>::faces_per_cell; ++n)
12630 if (!face_will_be_refined_by_neighbor(cell,n))
12631 ++unrefined_neighbors;
12647 if ((unrefined_neighbors == total_neighbors)
12650 (smooth_grid & eliminate_refined_inner_islands)) ||
12652 (smooth_grid & eliminate_refined_boundary_islands)) )
12654 (total_neighbors != 0))
12656 if (!cell->active())
12657 for (
unsigned int c=0; c<cell->n_children(); ++c)
12659 cell->child(c)->clear_refine_flag ();
12660 cell->child(c)->set_coarsen_flag ();
12663 cell->clear_refine_flag();
12678 if (smooth_grid & limit_level_difference_at_vertices)
12680 Assert(!anisotropic_refinement,
12681 ExcMessage(
"In case of anisotropic refinement the " 12682 "limit_level_difference_at_vertices flag for " 12683 "mesh smoothing must not be set!"));
12687 std::vector<int> vertex_level (vertices.size(), 0);
12690 for (; cell!=endc; ++cell)
12692 if (cell->refine_flag_set())
12693 for (
unsigned int vertex=0; vertex<GeometryInfo<dim>::vertices_per_cell;
12695 vertex_level[cell->vertex_index(vertex)]
12696 = std::max (vertex_level[cell->vertex_index(vertex)],
12698 else if (!cell->coarsen_flag_set())
12699 for (
unsigned int vertex=0; vertex<GeometryInfo<dim>::vertices_per_cell;
12701 vertex_level[cell->vertex_index(vertex)]
12702 = std::max (vertex_level[cell->vertex_index(vertex)],
12711 for (
unsigned int vertex=0; vertex<GeometryInfo<dim>::vertices_per_cell;
12713 vertex_level[cell->vertex_index(vertex)]
12714 = std::max (vertex_level[cell->vertex_index(vertex)],
12729 for (cell=last_active(); cell != endc; --cell)
12730 if (cell->refine_flag_set() ==
false)
12732 for (
unsigned int vertex=0;
12733 vertex<GeometryInfo<dim>::vertices_per_cell; ++vertex)
12734 if (vertex_level[cell->vertex_index(vertex)] >=
12738 cell->clear_coarsen_flag();
12743 if (vertex_level[cell->vertex_index(vertex)] >
12746 cell->set_refine_flag();
12748 for (
unsigned int v=0; v<GeometryInfo<dim>::vertices_per_cell;
12750 vertex_level[cell->vertex_index(v)]
12751 = std::max (vertex_level[cell->vertex_index(v)],
12773 if (smooth_grid & eliminate_unrefined_islands)
12778 for (; cell != endc; --cell)
12782 possibly_refine_unrefined_island<dim,spacedim>
12784 (smooth_grid & allow_anisotropic_smoothing) != 0);
12812 if (smooth_grid & patch_level_1)
12824 if (!cell->active())
12829 Assert (cell_is_patch_level_1(cell),
12831 if (cell->child(0)->has_children() ==
true)
12837 for (
unsigned int i=0; i<cell->n_children(); ++i)
12838 combined_ref_case = combined_ref_case |
12839 cell->child(i)->refine_flag_set();
12841 for (
unsigned int i=0; i<cell->n_children(); ++i)
12845 child->clear_coarsen_flag();
12846 child->set_refine_flag(combined_ref_case);
12867 cell->child(0)->active())
12874 const unsigned int n_children=cell->n_children();
12875 bool has_active_grandchildren =
false;
12877 for (
unsigned int i=0; i<n_children; ++i)
12878 if (cell->child(i)->child(0)->active())
12880 has_active_grandchildren =
true;
12884 if (has_active_grandchildren ==
false)
12890 unsigned int n_grandchildren=0;
12893 unsigned int n_coarsen_flags=0;
12898 for (
unsigned int c=0; c<n_children; ++c)
12905 const unsigned int nn_children=child->n_children();
12906 n_grandchildren += nn_children;
12911 if (child->child(0)->active())
12912 for (
unsigned int cc=0; cc<nn_children; ++cc)
12913 if (child->child(cc)->coarsen_flag_set())
12928 if ((n_coarsen_flags != n_grandchildren)
12930 (n_coarsen_flags > 0))
12931 for (
unsigned int c=0; c<n_children; ++c)
12934 if (child->child(0)->active())
12935 for (
unsigned int cc=0; cc<child->n_children(); ++cc)
12936 child->child(cc)->clear_coarsen_flag();
12963 bool changed =
true;
12970 for (; cell != endc; --cell)
12971 if (cell->refine_flag_set())
12974 for (
unsigned int i=0; i<GeometryInfo<dim>::faces_per_cell; ++i)
12994 if (cell->neighbor(i)->active())
12996 if (cell->neighbor_is_coarser(i))
12998 if (cell->neighbor(i)->coarsen_flag_set())
12999 cell->neighbor(i)->clear_coarsen_flag();
13015 if (smooth_grid & allow_anisotropic_smoothing)
13016 changed=cell->neighbor(i)->flag_for_face_refinement(cell->neighbor_of_coarser_neighbor(i).first,
13020 if (!cell->neighbor(i)->refine_flag_set())
13022 cell->neighbor(i)->set_refine_flag();
13145 std::pair<unsigned int, unsigned int> nb_indices
13146 =cell->neighbor_of_coarser_neighbor(i);
13147 unsigned int refined_along_x=0,
13149 to_be_refined_along_x=0,
13150 to_be_refined_along_y=0;
13152 const int this_face_index=cell->face_index(i);
13155 if ((this_face_index
13156 == cell->neighbor(i)->face(nb_indices.first)->child_index(0)) ||
13158 == cell->neighbor(i)->face(nb_indices.first)->child_index(1)))
13166 RefinementCase<dim-1> frc=cell->neighbor(i)->face(nb_indices.first)->refinement_case();
13184 cell->face_orientation(i),
13185 cell->face_flip(i),
13186 cell->face_rotation(i));
13188 ++to_be_refined_along_x;
13190 ++to_be_refined_along_y;
13193 if ((smooth_grid & allow_anisotropic_smoothing) ||
13194 cell->neighbor(i)->refine_flag_set())
13196 if (refined_along_x + to_be_refined_along_x > 1)
13197 changed |= cell->neighbor(i)->flag_for_face_refinement(nb_indices.first,
13199 if (refined_along_y + to_be_refined_along_y > 1)
13200 changed |= cell->neighbor(i)->flag_for_face_refinement(nb_indices.first,
13207 cell->neighbor(i)->set_refine_flag();
13215 nb->face_orientation(nb_indices.first),
13216 nb->face_flip(nb_indices.first),
13217 nb->face_rotation(nb_indices.first));
13219 !(refined_along_x || to_be_refined_along_x))
13222 !(refined_along_y || to_be_refined_along_y))
13228 cell->neighbor(i)->clear_coarsen_flag();
13229 const unsigned int nb_nb=cell->neighbor_of_neighbor(i);
13234 neighbor->face_orientation(nb_nb),
13235 neighbor->face_flip(nb_nb),
13236 neighbor->face_rotation(nb_nb));
13240 cell->face_orientation(i),
13241 cell->face_flip(i),
13242 cell->face_rotation(i));
13250 changed=cell->flag_for_face_refinement(i, face_ref_case);
13251 neighbor->flag_for_face_refinement(nb_nb, needed_face_ref_case);
13257 RefinementCase<dim-1> face_ref_case = cell->face(i)->refinement_case(),
13260 cell->face_orientation(i),
13261 cell->face_flip(i),
13262 cell->face_rotation(i));
13268 changed=cell->flag_for_face_refinement(i, face_ref_case);
13287 fix_coarsen_flags ();
13290 std::vector<bool> flags_after_loop[2];
13291 save_coarsen_flags (flags_after_loop[0]);
13292 save_refine_flags (flags_after_loop[1]);
13296 mesh_changed_in_this_loop
13297 = ((flags_before_loop[0] != flags_after_loop[0]) ||
13298 (flags_before_loop[1] != flags_after_loop[1]));
13302 flags_before_loop[0].swap(flags_after_loop[0]);
13303 flags_before_loop[1].swap(flags_after_loop[1]);
13305 while (mesh_changed_in_this_loop);
13311 return ((flags_before[0] != flags_before_loop[0]) ||
13312 (flags_before[1] != flags_before_loop[1]));
13318 template <
int dim,
int spacedim>
13320 const std::vector<bool> &v,
13321 const unsigned int magic_number2,
13324 const unsigned int N = v.size();
13325 unsigned char *flags =
new unsigned char[N/8+1];
13326 for (
unsigned int i=0; i<N/8+1; ++i) flags[i]=0;
13328 for (
unsigned int position=0; position<N; ++position)
13329 flags[position/8] |= (v[position] ? (1<<(position%8)) : 0);
13338 out << magic_number1 <<
' ' << N << std::endl;
13339 for (
unsigned int i=0; i<N/8+1; ++i)
13340 out << static_cast<unsigned int>(flags[i]) <<
' ';
13342 out << std::endl << magic_number2 << std::endl;
13350 template <
int dim,
int spacedim>
13352 std::vector<bool> &v,
13353 const unsigned int magic_number2,
13358 unsigned int magic_number;
13359 in >> magic_number;
13360 AssertThrow (magic_number==magic_number1, ExcGridReadError());
13366 unsigned char *flags =
new unsigned char[N/8+1];
13367 unsigned short int tmp;
13368 for (
unsigned int i=0; i<N/8+1; ++i)
13374 for (
unsigned int position=0; position!=N; ++position)
13375 v[position] = (flags[position/8] & (1<<(position%8)));
13377 in >> magic_number;
13378 AssertThrow (magic_number==magic_number2, ExcGridReadError());
13387 template <
int dim,
int spacedim>
13391 std::size_t mem = 0;
13393 for (
unsigned int i=0; i<levels.size(); ++i)
13397 mem +=
sizeof(manifold);
13398 mem +=
sizeof(smooth_grid);
13400 mem +=
sizeof (faces);
13410 template <
int dim,
int spacedim>
13421 #include "tria.inst" 13423 DEAL_II_NAMESPACE_CLOSE
std::vector< CellData< 1 > > boundary_lines
boost::signals2::signal< void()> any_change
::internal::TriangulationImplementation::NumberCache< dim > number_cache
unsigned int n_active_cells() const
void set_boundary(const types::manifold_id number, const Boundary< dim, spacedim > &boundary_object)
virtual void copy_triangulation(const Triangulation< dim, spacedim > &other_tria)
void loop(ITERATOR begin, typename identity< ITERATOR >::type end, DOFINFO &dinfo, INFOBOX &info, const std::function< void(DOFINFO &, typename INFOBOX::CellInfo &)> &cell_worker, const std::function< void(DOFINFO &, typename INFOBOX::CellInfo &)> &boundary_worker, const std::function< void(DOFINFO &, DOFINFO &, typename INFOBOX::CellInfo &, typename INFOBOX::CellInfo &)> &face_worker, ASSEMBLER &assembler, const LoopControl &lctrl=LoopControl())
unsigned int n_used_vertices() const
std::vector< unsigned int > n_quads_level
const types::manifold_id flat_manifold_id
static const unsigned int invalid_unsigned_int
void load_user_flags_line(std::istream &in)
#define AssertNothrow(cond, exc)
#define DeclException2(Exception2, type1, type2, outsequence)
static ::ExceptionBase & ExcInternalErrorOnCell(int arg1)
hex_iterator begin_hex(const unsigned int level=0) const
const types::subdomain_id invalid_subdomain_id
virtual ~DistortedCellList() noexcept
active_face_iterator begin_active_face() const
void save_user_flags_quad(std::ostream &out) const
virtual bool has_hanging_nodes() const
unsigned int n_raw_cells(const unsigned int level) const
cell_iterator last() const
void save_user_pointers_hex(std::vector< void *> &v) const
vertex_iterator begin_vertex() const
static void prepare_refinement_dim_dependent(const Triangulation< dim, spacedim > &)
static RefinementCase< 1 > line_refinement_case(const RefinementCase< dim > &cell_refinement_case, const unsigned int line_no)
unsigned int n_active_lines
unsigned int n_cells() const
static unsigned int face_to_cell_lines(const unsigned int face, const unsigned int line, const bool face_orientation=true, const bool face_flip=false, const bool face_rotation=false)
static ::ExceptionBase & ExcIO()
Subscriptor & operator=(const Subscriptor &)
Task< RT > new_task(const std::function< RT()> &function)
std::list< typename Triangulation< dim, spacedim >::cell_iterator > distorted_cells
static Triangulation< 1, spacedim >::DistortedCellList execute_refinement(Triangulation< 1, spacedim > &triangulation, const bool)
std::vector< Point< spacedim > > vertices
static void create_triangulation(const std::vector< Point< spacedim > > &v, const std::vector< CellData< 1 > > &cells, const SubCellData &, Triangulation< 1, spacedim > &triangulation)
boost::signals2::signal< void()> clear
static unsigned int line_to_cell_vertices(const unsigned int line, const unsigned int vertex)
IteratorRange< active_cell_iterator > active_cell_iterators() const
void set_all_refine_flags()
void load_user_flags(std::istream &in)
bool anisotropic_refinement
int face(const unsigned int i) const
static void read_bool_vector(const unsigned int magic_number1, std::vector< bool > &v, const unsigned int magic_number2, std::istream &in)
line_iterator begin_line(const unsigned int level=0) const
line_iterator end_line() const
boost::signals2::signal< void()> mesh_movement
void clear_user_flags_quad()
static void create_children(Triangulation< 2, spacedim > &triangulation, unsigned int &next_unused_vertex, typename Triangulation< 2, spacedim >::raw_line_iterator &next_unused_line, typename Triangulation< 2, spacedim >::raw_cell_iterator &next_unused_cell, typename Triangulation< 2, spacedim >::cell_iterator &cell)
unsigned int n_faces() const
static ::ExceptionBase & ExcCellHasNegativeMeasure(int arg1)
void flip_all_direction_flags()
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)
static Triangulation< 2, spacedim >::DistortedCellList execute_refinement(Triangulation< 2, spacedim > &triangulation, const bool check_for_distorted_cells)
active_cell_iterator begin_active(const unsigned int level=0) const
#define AssertThrow(cond, exc)
static ::ExceptionBase & ExcInvalidVertexIndex(int arg1, int arg2, int arg3)
static void compute_number_cache(const Triangulation< dim, spacedim > &triangulation, const unsigned int level_objects, internal::TriangulationImplementation::NumberCache< 1 > &number_cache)
static void write_bool_vector(const unsigned int magic_number1, const std::vector< bool > &v, const unsigned int magic_number2, std::ostream &out)
unsigned int n_active_hexs() const
void save_user_pointers_quad(std::vector< void *> &v) const
static ::ExceptionBase & ExcIndexRange(int arg1, int arg2, int arg3)
static ::ExceptionBase & ExcFacesHaveNoLevel()
void execute_coarsening()
const Boundary< dim, spacedim > & get_boundary(const types::manifold_id number) const
cell_iterator begin(const unsigned int level=0) const
void load_user_flags_quad(std::istream &in)
void set_manifold(const types::manifold_id number, const Manifold< dim, spacedim > &manifold_object)
static Triangulation< 3, spacedim >::DistortedCellList execute_refinement(Triangulation< 3, spacedim > &triangulation, const bool check_for_distorted_cells)
raw_hex_iterator begin_raw_hex(const unsigned int level=0) const
std::vector< std::unique_ptr<::internal::TriangulationImplementation::TriaLevel< dim > > > levels
void load_user_pointers_line(const std::vector< void *> &v)
unsigned int n_active_faces() const
raw_quad_iterator begin_raw_quad(const unsigned int level=0) const
std::unique_ptr<::internal::TriangulationImplementation::TriaFaces< dim > > faces
static ::ExceptionBase & ExcInteriorQuadCantBeBoundary(int arg1, int arg2, int arg3, int arg4, types::boundary_id arg5)
cell_iterator end() const
void load_coarsen_flags(std::istream &out)
void load_user_flags_hex(std::istream &in)
active_quad_iterator begin_active_quad(const unsigned int level=0) const
virtual void execute_coarsening_and_refinement()
static ::ExceptionBase & ExcGridHasInvalidCell(int arg1)
unsigned int n_active_quads
unsigned int n_active_lines() const
IteratorRange< cell_iterator > cell_iterators_on_level(const unsigned int level) const
virtual bool prepare_coarsening_and_refinement()
static ::ExceptionBase & ExcMessage(std::string arg1)
virtual std::unique_ptr< Manifold< dim, spacedim > > clone() const =0
bool check_consistency(const unsigned int dim) const
static ::ExceptionBase & ExcImpossibleInDim(int arg1)
void reset_active_cell_indices()
DistortedCellList execute_refinement()
void save_user_flags_line(std::ostream &out) const
void save_user_pointers_line(std::vector< void *> &v) const
Triangulation(const MeshSmoothing smooth_grid=none, const bool check_for_distorted_cells=false)
#define DeclException1(Exception1, type1, outsequence)
unsigned int n_active_hexes
void clear_user_flags_hex()
void load_user_indices_line(const std::vector< unsigned int > &v)
unsigned int global_dof_index
vertex_iterator end_vertex() const
unsigned int n_raw_quads() const
static unsigned int child_cell_on_face(const RefinementCase< dim > &ref_case, const unsigned int face, const unsigned int subface, const bool face_orientation=true, const bool face_flip=false, const bool face_rotation=false, const RefinementCase< dim-1 > &face_refinement_case=RefinementCase< dim-1 >::isotropic_refinement)
#define Assert(cond, exc)
active_hex_iterator begin_active_hex(const unsigned int level=0) const
void save_user_indices_line(std::vector< unsigned int > &v) const
static void prevent_distorted_boundary_cells(const Triangulation< 1, spacedim > &)
void set_all_manifold_ids(const types::manifold_id number)
static ::ExceptionBase & ExcInteriorLineCantBeBoundary(int arg1, int arg2, types::boundary_id arg3)
boost::signals2::signal< void()> create
unsigned int n_raw_hexs(const unsigned int level) const
unsigned int max_adjacent_cells() const
bool get_anisotropic_refinement_flag() const
unsigned int n_quads() const
std::unique_ptr< std::map< unsigned int, types::manifold_id > > vertex_to_manifold_id_map_1d
std::vector< bool > vertices_used
virtual void create_triangulation(const std::vector< Point< spacedim > > &vertices, const std::vector< CellData< dim > > &cells, const SubCellData &subcelldata)
void save_coarsen_flags(std::ostream &out) const
void reset_all_manifolds()
void save_user_indices(std::vector< unsigned int > &v) const
void load_user_indices(const std::vector< unsigned int > &v)
void save_refine_flags(std::ostream &out) const
unsigned int n_lines() const
static ::ExceptionBase & ExcMultiplySetLineInfoOfLine(int arg1, int arg2)
void clear_despite_subscriptions()
raw_cell_iterator begin_raw(const unsigned int level=0) const
unsigned int n_raw_lines() const
virtual std::size_t memory_consumption() const
#define DeclException5(Exception5, type1, type2, type3, type4, type5, outsequence)
virtual void set_mesh_smoothing(const MeshSmoothing mesh_smoothing)
void load_user_pointers_hex(const std::vector< void *> &v)
void save_user_flags_hex(std::ostream &out) const
face_iterator begin_face() const
unsigned int subdomain_id
boost::signals2::signal< void()> post_refinement
std::unique_ptr< std::map< unsigned int, types::boundary_id > > vertex_to_boundary_id_map_1d
virtual types::global_dof_index n_global_active_cells() const
unsigned int n_hexs() const
static void create_triangulation(const std::vector< Point< spacedim > > &v, const std::vector< CellData< 2 > > &cells, const SubCellData &subcelldata, Triangulation< 2, spacedim > &triangulation)
unsigned int n_active_quads() const
active_line_iterator begin_active_line(const unsigned int level=0) const
raw_cell_iterator end_raw(const unsigned int level) const
unsigned int n_raw_faces() const
std::vector< unsigned int > n_lines_level
void load_user_indices_hex(const std::vector< unsigned int > &v)
static ::ExceptionBase & ExcLineInexistant(int arg1, int arg2)
static unsigned int n_children(const RefinementCase< dim > &refinement_case)
void swap(Vector< Number > &u, Vector< Number > &v)
active_vertex_iterator begin_active_vertex() const
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])
static bool coarsening_allowed(const typename Triangulation< dim, spacedim >::cell_iterator &cell)
std::vector< types::boundary_id > get_boundary_ids() const
const types::manifold_id invalid_manifold_id
TriaObjects< TriaObject< 1 > > lines
void load_user_pointers_quad(const std::vector< void *> &v)
TriaObjects< TriaObject< 1 > > lines
active_cell_iterator end_active(const unsigned int level) const
void update_periodic_face_map()
std::vector< unsigned int > n_active_hexes_level
static void compute_number_cache(const Triangulation< dim, spacedim > &triangulation, const unsigned int level_objects, internal::TriangulationImplementation::NumberCache< 3 > &number_cache)
void save_user_indices_quad(std::vector< unsigned int > &v) const
virtual void add_periodicity(const std::vector< GridTools::PeriodicFacePair< cell_iterator > > &)
const std::vector< bool > & get_used_vertices() const
hex_iterator end_hex() const
std::vector< CellData< 2 > > boundary_quads
MeshSmoothing smooth_grid
#define DeclException4(Exception4, type1, type2, type3, type4, outsequence)
IteratorRange< active_cell_iterator > active_cell_iterators_on_level(const unsigned int level) const
void load_user_pointers(const std::vector< void *> &v)
std::vector< unsigned int > n_active_lines_level
void refine_global(const unsigned int times=1)
void save_user_pointers(std::vector< void *> &v) const
void save_user_flags(std::ostream &out) const
static ::ExceptionBase & ExcNotImplemented()
void clear_user_data(const unsigned int i)
Iterator points to a valid object.
quad_iterator end_quad() const
IteratorRange< cell_iterator > cell_iterators() const
#define DeclException3(Exception3, type1, type2, type3, outsequence)
face_iterator end_face() const
void set_all_manifold_ids_on_boundary(const types::manifold_id number)
static ::ExceptionBase & ExcQuadInexistant(int arg1, int arg2, int arg3, int arg4)
active_cell_iterator last_active() const
const types::boundary_id internal_face_boundary_id
std::map< types::manifold_id, std::unique_ptr< const Manifold< dim, spacedim > > > manifold
std::vector< unsigned int > n_active_quads_level
void clear_user_flags_line()
IteratorState::IteratorStates state() const
quad_iterator begin_quad(const unsigned int level=0) const
void set_face(const unsigned int i, const int index)
Triangulation & operator=(Triangulation< dim, spacedim > &&tria) noexcept
std::vector< unsigned int > n_hexes_level
std::vector< types::manifold_id > get_manifold_ids() const
static void delete_children(Triangulation< 1, spacedim > &triangulation, typename Triangulation< 1, spacedim >::cell_iterator &cell, std::vector< unsigned int > &, std::vector< unsigned int > &)
const Manifold< dim, spacedim > & get_manifold(const types::manifold_id number) const
virtual void create_triangulation_compatibility(const std::vector< Point< spacedim > > &vertices, const std::vector< CellData< dim > > &cells, const SubCellData &subcelldata)
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)
void load_user_indices_quad(const std::vector< unsigned int > &v)
virtual types::subdomain_id locally_owned_subdomain() const
raw_line_iterator begin_raw_line(const unsigned int level=0) const
virtual const MeshSmoothing & get_mesh_smoothing() const
std::enable_if< std::is_fundamental< T >::value, std::size_t >::type memory_consumption(const T &t)
const std::map< std::pair< cell_iterator, unsigned int >, std::pair< std::pair< cell_iterator, unsigned int >, std::bitset< 3 > > > & get_periodic_face_map() const
void load_refine_flags(std::istream &in)
Tensor< 2, dim, Number > l(const Tensor< 2, dim, Number > &F, const Tensor< 2, dim, Number > &dF_dt)
static ::ExceptionBase & ExcInternalError()
Triangulation< dim, spacedim > & get_triangulation()
static void create_triangulation(const std::vector< Point< spacedim > > &v, const std::vector< CellData< 3 > > &cells, const SubCellData &subcelldata, Triangulation< 3, spacedim > &triangulation)
static void compute_number_cache(const Triangulation< dim, spacedim > &triangulation, const unsigned int level_objects, internal::TriangulationImplementation::NumberCache< 2 > &number_cache)
void save_user_indices_hex(std::vector< unsigned int > &v) const