48 #define BOOST_BIND_GLOBAL_PLACEHOLDERS
49 #include <boost/config.hpp>
50 #include <boost/graph/adjacency_list.hpp>
51 #include <boost/graph/bandwidth.hpp>
52 #include <boost/graph/cuthill_mckee_ordering.hpp>
53 #include <boost/graph/king_ordering.hpp>
54 #include <boost/graph/minimum_degree_ordering.hpp>
55 #include <boost/graph/properties.hpp>
56 #include <boost/random.hpp>
57 #include <boost/random/uniform_int_distribution.hpp>
59 #undef BOOST_BIND_GLOBAL_PLACEHOLDERS
77 using namespace ::
boost;
79 using Graph = adjacency_list<vecS,
82 property<vertex_color_t,
84 property<vertex_degree_t, int>>>;
85 using Vertex = graph_traits<Graph>::vertex_descriptor;
86 using size_type = graph_traits<Graph>::vertices_size_type;
88 using Pair = std::pair<size_type, size_type>;
94 template <
int dim,
int spacedim>
97 const bool use_constraints,
100 boosttypes::vertex_degree_t>::type
115 for (
unsigned int row = 0; row < dsp.
n_rows(); ++row)
116 for (
unsigned int col = 0; col < dsp.
row_length(row); ++col)
120 boosttypes::graph_traits<boosttypes::Graph>::vertex_iterator ui, ui_end;
122 graph_degree = get(::boost::vertex_degree, graph);
123 for (::boost::tie(ui, ui_end) =
vertices(graph); ui != ui_end; ++ui)
124 graph_degree[*ui] = degree(*ui, graph);
129 template <
int dim,
int spacedim>
132 const bool reversed_numbering,
133 const bool use_constraints)
135 std::vector<types::global_dof_index> renumbering(
149 template <
int dim,
int spacedim>
153 const bool reversed_numbering,
154 const bool use_constraints)
158 boosttypes::vertex_degree_t>::type graph_degree;
163 boosttypes::vertex_index_t>::type index_map =
164 get(::boost::vertex_index, graph);
167 std::vector<boosttypes::Vertex> inv_perm(num_vertices(graph));
169 if (reversed_numbering ==
false)
170 ::boost::cuthill_mckee_ordering(graph,
172 get(::boost::vertex_color, graph),
173 make_degree_map(graph));
175 ::boost::cuthill_mckee_ordering(graph,
177 get(::boost::vertex_color, graph),
178 make_degree_map(graph));
181 new_dof_indices[index_map[inv_perm[c]]] = c;
183 Assert(std::find(new_dof_indices.begin(),
184 new_dof_indices.end(),
191 template <
int dim,
int spacedim>
194 const bool reversed_numbering,
195 const bool use_constraints)
197 std::vector<types::global_dof_index> renumbering(
211 template <
int dim,
int spacedim>
215 const bool reversed_numbering,
216 const bool use_constraints)
220 boosttypes::vertex_degree_t>::type graph_degree;
225 boosttypes::vertex_index_t>::type index_map =
226 get(::boost::vertex_index, graph);
229 std::vector<boosttypes::Vertex> inv_perm(num_vertices(graph));
231 if (reversed_numbering ==
false)
237 new_dof_indices[index_map[inv_perm[c]]] = c;
239 Assert(std::find(new_dof_indices.begin(),
240 new_dof_indices.end(),
247 template <
int dim,
int spacedim>
250 const bool reversed_numbering,
251 const bool use_constraints)
253 std::vector<types::global_dof_index> renumbering(
267 template <
int dim,
int spacedim>
270 std::vector<types::global_dof_index> &new_dof_indices,
272 const bool reversed_numbering,
273 const bool use_constraints)
275 (void)use_constraints;
284 using namespace ::
boost;
289 using Graph = adjacency_list<vecS, vecS, directedS>;
291 int n = dof_handler.
n_dofs();
295 std::vector<::types::global_dof_index> dofs_on_this_cell;
299 const unsigned int dofs_per_cell = cell->get_fe().n_dofs_per_cell();
301 dofs_on_this_cell.resize(dofs_per_cell);
303 cell->get_active_or_mg_dof_indices(dofs_on_this_cell);
304 for (
unsigned int i = 0; i < dofs_per_cell; ++i)
305 for (
unsigned int j = 0; j < dofs_per_cell; ++j)
306 if (dofs_on_this_cell[i] > dofs_on_this_cell[j])
308 add_edge(dofs_on_this_cell[i], dofs_on_this_cell[j], G);
309 add_edge(dofs_on_this_cell[j], dofs_on_this_cell[i], G);
314 using Vector = std::vector<int>;
317 Vector inverse_perm(n, 0);
322 Vector supernode_sizes(n, 1);
325 ::boost::property_map<Graph, vertex_index_t>::type
id =
326 get(vertex_index, G);
332 minimum_degree_ordering(
334 make_iterator_property_map(degree.
begin(),
id, degree[0]),
337 make_iterator_property_map(supernode_sizes.
begin(),
344 for (
int i = 0; i < n; ++i)
354 if (reversed_numbering ==
true)
359 new_dof_indices.begin());
366 template <
int dim,
int spacedim>
369 const bool reversed_numbering,
370 const bool use_constraints,
371 const std::vector<types::global_dof_index> &starting_indices)
373 std::vector<types::global_dof_index> renumbering(
390 template <
int dim,
int spacedim>
393 std::vector<types::global_dof_index> & new_indices,
395 const bool reversed_numbering,
396 const bool use_constraints,
397 const std::vector<types::global_dof_index> &starting_indices,
398 const unsigned int level)
400 const bool reorder_level_dofs =
417 if (reorder_level_dofs ==
false)
420 locally_relevant_dofs);
429 locally_relevant_dofs);
440 constraints.
reinit(locally_relevant_dofs);
451 locally_owned_dofs.
size());
452 if (reorder_level_dofs ==
false)
464 if (reversed_numbering)
478 if (reorder_level_dofs ==
false)
481 locally_active_dofs);
490 bool needs_locally_active =
false;
491 for (
const auto starting_index : starting_indices)
493 if ((needs_locally_active ==
495 (locally_owned_dofs.
is_element(starting_index) ==
false))
498 locally_active_dofs.
is_element(starting_index),
500 "You specified global degree of freedom " +
502 " as a starting index, but this index is not among the "
503 "locally active ones on this processor, as required "
504 "for this function."));
505 needs_locally_active =
true;
510 (needs_locally_active ? locally_active_dofs : locally_owned_dofs);
521 index_set_to_use.
size(),
523 if (reorder_level_dofs ==
false)
534 std::vector<types::global_dof_index> row_entries;
535 for (
unsigned int i = 0; i < index_set_to_use.
n_elements(); ++i)
539 const unsigned int row_length = dsp.
row_length(row);
541 for (
unsigned int j = 0; j < row_length; ++j)
544 if (col != row && index_set_to_use.
is_element(col))
554 std::vector<types::global_dof_index> local_starting_indices(
555 starting_indices.size());
556 for (
unsigned int i = 0; i < starting_indices.size(); ++i)
557 local_starting_indices[i] =
562 std::vector<types::global_dof_index> my_new_indices(
566 local_starting_indices);
567 if (reversed_numbering)
576 if (needs_locally_active ==
true)
579 IndexSet active_but_not_owned_dofs = locally_active_dofs;
580 active_but_not_owned_dofs.
subtract_set(locally_owned_dofs);
582 std::set<types::global_dof_index> erase_these_indices;
583 for (
const auto p : active_but_not_owned_dofs)
588 erase_these_indices.insert(my_new_indices[index]);
591 Assert(erase_these_indices.size() ==
594 Assert(
static_cast<unsigned int>(
595 std::count(my_new_indices.begin(),
596 my_new_indices.end(),
602 std::vector<types::global_dof_index> translate_indices(
603 my_new_indices.size());
605 std::set<types::global_dof_index>::const_iterator
606 next_erased_index = erase_these_indices.begin();
608 for (
unsigned int i = 0; i < translate_indices.size(); ++i)
609 if ((next_erased_index != erase_these_indices.end()) &&
610 (*next_erased_index == i))
617 translate_indices[i] = next_new_index;
627 new_indices.reserve(locally_owned_dofs.
n_elements());
628 for (
const auto &p : my_new_indices)
633 new_indices.push_back(translate_indices[p]);
639 new_indices = std::move(my_new_indices);
652 template <
int dim,
int spacedim>
655 const unsigned int level,
656 const bool reversed_numbering,
657 const std::vector<types::global_dof_index> &starting_indices)
662 std::vector<types::global_dof_index> new_indices(
681 template <
int dim,
int spacedim>
684 const std::vector<unsigned int> &component_order_arg)
686 std::vector<types::global_dof_index> renumbering(
690 compute_component_wise<dim, spacedim>(renumbering,
712 (result <= dof_handler.
n_dofs())),
720 template <
int dim,
int spacedim>
723 const unsigned int level,
724 const std::vector<unsigned int> &component_order_arg)
729 std::vector<types::global_dof_index> renumbering(
739 compute_component_wise<dim, spacedim>(
740 renumbering, start,
end, component_order_arg,
true);
753 template <
int dim,
int spacedim,
typename CellIterator>
756 const CellIterator & start,
758 const std::vector<unsigned int> &component_order_arg,
759 const bool is_level_operation)
762 start->get_dof_handler().get_fe_collection());
768 new_indices.resize(0);
775 const IndexSet &locally_owned_dofs =
777 start->get_dof_handler().locally_owned_mg_dofs(start->level()) :
778 start->get_dof_handler().locally_owned_dofs();
782 std::vector<unsigned int> component_order(component_order_arg);
787 if (component_order.size() == 0)
788 for (
unsigned int i = 0; i < fe_collection.
n_components(); ++i)
789 component_order.push_back(i);
795 for (
const unsigned int component : component_order)
803 std::vector<types::global_dof_index> local_dof_indices;
815 std::vector<std::vector<unsigned int>> component_list(fe_collection.
size());
816 for (
unsigned int f = 0; f < fe_collection.
size(); ++f)
820 component_list[f].resize(dofs_per_cell);
821 for (
unsigned int i = 0; i < dofs_per_cell; ++i)
823 component_list[f][i] =
827 const unsigned int comp =
833 component_list[f][i] = component_order[comp];
848 std::vector<std::vector<types::global_dof_index>> component_to_dof_map(
850 for (CellIterator cell = start; cell !=
end; ++cell)
852 if (is_level_operation)
855 if (!cell->is_locally_owned_on_level())
862 if (!cell->is_locally_owned())
866 if (is_level_operation)
868 cell->level() == start->level(),
870 "Multigrid renumbering in compute_component_wise() needs to be applied to a single level!"));
875 const unsigned int fe_index = cell->active_fe_index();
876 const unsigned int dofs_per_cell =
877 fe_collection[fe_index].n_dofs_per_cell();
878 local_dof_indices.resize(dofs_per_cell);
879 cell->get_active_or_mg_dof_indices(local_dof_indices);
881 for (
unsigned int i = 0; i < dofs_per_cell; ++i)
882 if (locally_owned_dofs.
is_element(local_dof_indices[i]))
883 component_to_dof_map[component_list[fe_index][i]].
push_back(
884 local_dof_indices[i]);
910 for (
unsigned int component = 0; component < fe_collection.
n_components();
913 std::sort(component_to_dof_map[component].
begin(),
914 component_to_dof_map[component].
end());
915 component_to_dof_map[component].erase(
916 std::unique(component_to_dof_map[component].
begin(),
917 component_to_dof_map[component].
end()),
918 component_to_dof_map[component].
end());
923 const unsigned int n_buckets = fe_collection.
n_components();
924 std::vector<types::global_dof_index> shifts(n_buckets);
928 &start->get_dof_handler().get_triangulation())))
930 #ifdef DEAL_II_WITH_MPI
931 std::vector<types::global_dof_index> local_dof_count(n_buckets);
933 for (
unsigned int c = 0; c < n_buckets; ++c)
934 local_dof_count[c] = component_to_dof_map[c].size();
936 std::vector<types::global_dof_index> prefix_dof_count(n_buckets);
937 const int ierr = MPI_Exscan(local_dof_count.data(),
938 prefix_dof_count.data(),
945 std::vector<types::global_dof_index> global_dof_count(n_buckets);
952 for (
unsigned int c = 0; c < n_buckets; ++c)
954 shifts[c] = prefix_dof_count[c] + cumulated;
955 cumulated += global_dof_count[c];
965 for (
unsigned int c = 1; c < fe_collection.
n_components(); ++c)
966 shifts[c] = shifts[c - 1] + component_to_dof_map[c - 1].size();
975 for (
unsigned int component = 0; component < fe_collection.
n_components();
978 next_free_index = shifts[component];
981 component_to_dof_map[component])
993 return next_free_index;
998 template <
int dim,
int spacedim>
1002 std::vector<types::global_dof_index> renumbering(
1023 (result <= dof_handler.
n_dofs())),
1031 template <
int dim,
int spacedim>
1038 std::vector<types::global_dof_index> renumbering(
1067 template <
int dim,
int spacedim,
class ITERATOR,
class ENDITERATOR>
1070 const ITERATOR & start,
1071 const ENDITERATOR &
end,
1072 const bool is_level_operation)
1075 start->get_dof_handler().get_fe_collection());
1081 new_indices.resize(0);
1088 const IndexSet &locally_owned_dofs =
1089 is_level_operation ?
1090 start->get_dof_handler().locally_owned_mg_dofs(start->level()) :
1091 start->get_dof_handler().locally_owned_dofs();
1095 std::vector<types::global_dof_index> local_dof_indices;
1100 std::vector<std::vector<types::global_dof_index>> block_list(
1101 fe_collection.
size());
1102 for (
unsigned int f = 0; f < fe_collection.
size(); ++f)
1121 std::vector<std::vector<types::global_dof_index>> block_to_dof_map(
1123 for (ITERATOR cell = start; cell !=
end; ++cell)
1125 if (is_level_operation)
1128 if (!cell->is_locally_owned_on_level())
1135 if (!cell->is_locally_owned())
1139 if (is_level_operation)
1141 cell->level() == start->level(),
1143 "Multigrid renumbering in compute_block_wise() needs to be applied to a single level!"));
1148 const unsigned int fe_index = cell->active_fe_index();
1149 const unsigned int dofs_per_cell =
1150 fe_collection[fe_index].n_dofs_per_cell();
1151 local_dof_indices.resize(dofs_per_cell);
1152 cell->get_active_or_mg_dof_indices(local_dof_indices);
1154 for (
unsigned int i = 0; i < dofs_per_cell; ++i)
1155 if (locally_owned_dofs.
is_element(local_dof_indices[i]))
1156 block_to_dof_map[block_list[fe_index][i]].
push_back(
1157 local_dof_indices[i]);
1172 for (
unsigned int block = 0; block < fe_collection.
n_blocks(); ++block)
1174 std::sort(block_to_dof_map[block].
begin(),
1175 block_to_dof_map[block].
end());
1176 block_to_dof_map[block].erase(
1177 std::unique(block_to_dof_map[block].
begin(),
1178 block_to_dof_map[block].
end()),
1179 block_to_dof_map[block].
end());
1184 const unsigned int n_buckets = fe_collection.
n_blocks();
1185 std::vector<types::global_dof_index> shifts(n_buckets);
1189 &start->get_dof_handler().get_triangulation())))
1191 #ifdef DEAL_II_WITH_MPI
1192 std::vector<types::global_dof_index> local_dof_count(n_buckets);
1194 for (
unsigned int c = 0; c < n_buckets; ++c)
1195 local_dof_count[c] = block_to_dof_map[c].size();
1197 std::vector<types::global_dof_index> prefix_dof_count(n_buckets);
1198 const int ierr = MPI_Exscan(local_dof_count.data(),
1199 prefix_dof_count.data(),
1206 std::vector<types::global_dof_index> global_dof_count(n_buckets);
1213 for (
unsigned int c = 0; c < n_buckets; ++c)
1215 shifts[c] = prefix_dof_count[c] + cumulated;
1216 cumulated += global_dof_count[c];
1226 for (
unsigned int c = 1; c < fe_collection.
n_blocks(); ++c)
1227 shifts[c] = shifts[c - 1] + block_to_dof_map[c - 1].size();
1236 for (
unsigned int block = 0; block < fe_collection.
n_blocks(); ++block)
1238 const typename std::vector<types::global_dof_index>::const_iterator
1239 begin_of_component = block_to_dof_map[block].begin(),
1240 end_of_component = block_to_dof_map[block].end();
1242 next_free_index = shifts[block];
1244 for (
typename std::vector<types::global_dof_index>::const_iterator
1245 dof_index = begin_of_component;
1246 dof_index != end_of_component;
1257 return next_free_index;
1269 template <
int dim,
class CellIteratorType>
1271 compute_hierarchical_recursive(
1274 const CellIteratorType & cell,
1275 const IndexSet & locally_owned_dof_indices,
1276 std::vector<types::global_dof_index> &new_indices)
1279 next_free_dof_offset;
1281 if (cell->has_children())
1284 for (
unsigned int c = 0; c < GeometryInfo<dim>::max_children_per_cell;
1286 current_next_free_dof_offset =
1287 compute_hierarchical_recursive<dim>(current_next_free_dof_offset,
1290 locally_owned_dof_indices,
1310 if (cell->is_locally_owned())
1313 const unsigned int dofs_per_cell =
1314 cell->get_fe().n_dofs_per_cell();
1315 std::vector<types::global_dof_index> local_dof_indices(
1317 cell->get_dof_indices(local_dof_indices);
1336 for (
unsigned int i = 0; i < dofs_per_cell; ++i)
1337 if (locally_owned_dof_indices.
is_element(local_dof_indices[i]))
1341 const unsigned int idx =
1343 local_dof_indices[i]);
1347 my_starting_index + current_next_free_dof_offset;
1348 ++current_next_free_dof_offset;
1354 return current_next_free_dof_offset;
1360 template <
int dim,
int spacedim>
1364 std::vector<types::global_dof_index> renumbering(
1387 #ifdef DEAL_II_WITH_MPI
1390 const int ierr = MPI_Exscan(&locally_owned_size,
1404 #ifdef DEAL_II_WITH_P4EST
1409 for (
unsigned int c = 0; c <
tria->
n_cells(0); ++c)
1411 const unsigned int coarse_cell_index =
1412 tria->get_p4est_tree_to_coarse_cell_permutation()[c];
1415 this_cell(
tria, 0, coarse_cell_index, &dof_handler);
1417 next_free_dof_offset =
1418 compute_hierarchical_recursive<dim>(next_free_dof_offset,
1433 dof_handler.
begin(0);
1434 cell != dof_handler.
end(0);
1436 next_free_dof_offset =
1437 compute_hierarchical_recursive<dim>(next_free_dof_offset,
1450 (next_free_dof_offset <= dof_handler.
n_dofs())),
1454 Assert(std::find(renumbering.begin(),
1464 template <
int dim,
int spacedim>
1467 const std::vector<bool> & selected_dofs)
1469 std::vector<types::global_dof_index> renumbering(
1478 template <
int dim,
int spacedim>
1481 const std::vector<bool> & selected_dofs,
1482 const unsigned int level)
1487 std::vector<types::global_dof_index> renumbering(
1499 template <
int dim,
int spacedim>
1502 std::vector<types::global_dof_index> &new_indices,
1504 const std::vector<bool> & selected_dofs)
1507 Assert(selected_dofs.size() == n_dofs,
1512 Assert(new_indices.size() == n_dofs,
1516 std::count(selected_dofs.begin(), selected_dofs.end(),
false);
1521 if (selected_dofs[i] ==
false)
1523 new_indices[i] = next_unselected;
1528 new_indices[i] = next_selected;
1537 template <
int dim,
int spacedim>
1540 std::vector<types::global_dof_index> &new_indices,
1542 const std::vector<bool> & selected_dofs,
1543 const unsigned int level)
1548 const unsigned int n_dofs = dof_handler.
n_dofs(
level);
1549 Assert(selected_dofs.size() == n_dofs,
1554 Assert(new_indices.size() == n_dofs,
1557 const unsigned int n_selected_dofs =
1558 std::count(selected_dofs.begin(), selected_dofs.end(),
false);
1560 unsigned int next_unselected = 0;
1561 unsigned int next_selected = n_selected_dofs;
1562 for (
unsigned int i = 0; i < n_dofs; ++i)
1563 if (selected_dofs[i] ==
false)
1565 new_indices[i] = next_unselected;
1570 new_indices[i] = next_selected;
1579 template <
int dim,
int spacedim>
1586 std::vector<types::global_dof_index> renumbering(
1595 template <
int dim,
int spacedim>
1598 std::vector<types::global_dof_index> &new_indices,
1599 std::vector<types::global_dof_index> &reverse,
1601 const typename std::vector<
1626 std::vector<bool> already_sorted(n_owned_dofs,
false);
1627 std::vector<types::global_dof_index> cell_dofs;
1631 unsigned int index = 0;
1633 for (
const auto &cell : cells)
1637 const unsigned int n_cell_dofs = cell->get_fe().n_dofs_per_cell();
1638 cell_dofs.resize(n_cell_dofs);
1640 cell->get_active_or_mg_dof_indices(cell_dofs);
1644 std::sort(cell_dofs.begin(), cell_dofs.end());
1646 for (
const auto dof : cell_dofs)
1648 const auto local_dof = owned_dofs.index_within_set(dof);
1650 !already_sorted[local_dof])
1652 already_sorted[local_dof] =
true;
1653 reverse[index++] = local_dof;
1657 Assert(index == n_owned_dofs,
1659 "Traversing over the given set of cells did not cover all "
1660 "degrees of freedom in the DoFHandler. Does the set of cells "
1661 "not include all active cells?"));
1664 new_indices[reverse[i]] = owned_dofs.nth_index_in_set(i);
1669 template <
int dim,
int spacedim>
1672 const unsigned int level,
1673 const typename std::vector<
1679 std::vector<types::global_dof_index> renumbering(dof.
n_dofs(
level));
1680 std::vector<types::global_dof_index> reverse(dof.
n_dofs(
level));
1688 template <
int dim,
int spacedim>
1691 std::vector<types::global_dof_index> &new_order,
1692 std::vector<types::global_dof_index> &reverse,
1694 const unsigned int level,
1695 const typename std::vector<
1709 std::vector<bool> already_sorted(n_global_dofs,
false);
1710 std::vector<types::global_dof_index> cell_dofs(n_cell_dofs);
1714 for (
const auto &cell : cells)
1718 cell->get_active_or_mg_dof_indices(cell_dofs);
1719 std::sort(cell_dofs.begin(), cell_dofs.end());
1721 for (
unsigned int i = 0; i < n_cell_dofs; ++i)
1723 if (!already_sorted[cell_dofs[i]])
1725 already_sorted[cell_dofs[i]] =
true;
1732 "Traversing over the given set of cells did not cover all "
1733 "degrees of freedom in the DoFHandler. Does the set of cells "
1734 "not include all cells of the specified level?"));
1737 new_order[reverse[i]] = i;
1742 template <
int dim,
int spacedim>
1746 const bool dof_wise_renumbering)
1748 std::vector<types::global_dof_index> renumbering(dof.
n_dofs());
1749 std::vector<types::global_dof_index> reverse(dof.
n_dofs());
1751 renumbering, reverse, dof, direction, dof_wise_renumbering);
1758 template <
int dim,
int spacedim>
1761 std::vector<types::global_dof_index> &reverse,
1764 const bool dof_wise_renumbering)
1770 if (dof_wise_renumbering ==
false)
1772 std::vector<typename DoFHandler<dim, spacedim>::active_cell_iterator>
1778 comparator(direction);
1781 ordered_cells.push_back(cell);
1783 std::sort(ordered_cells.begin(), ordered_cells.end(), comparator);
1795 const unsigned int n_dofs = dof.
n_dofs();
1796 std::vector<std::pair<Point<spacedim>,
unsigned int>>
1797 support_point_list(n_dofs);
1800 Assert(fe_collection[0].has_support_points(),
1803 for (
unsigned int comp = 0; comp < fe_collection.
size(); ++comp)
1805 Assert(fe_collection[comp].has_support_points(),
1811 quadrature_collection,
1814 std::vector<bool> already_touched(n_dofs,
false);
1816 std::vector<types::global_dof_index> local_dof_indices;
1820 const unsigned int dofs_per_cell = cell->get_fe().n_dofs_per_cell();
1821 local_dof_indices.resize(dofs_per_cell);
1822 hp_fe_values.
reinit(cell);
1825 cell->get_active_or_mg_dof_indices(local_dof_indices);
1826 const std::vector<Point<spacedim>> &points =
1828 for (
unsigned int i = 0; i < dofs_per_cell; ++i)
1829 if (!already_touched[local_dof_indices[i]])
1831 support_point_list[local_dof_indices[i]].first = points[i];
1832 support_point_list[local_dof_indices[i]].second =
1833 local_dof_indices[i];
1834 already_touched[local_dof_indices[i]] =
true;
1839 std::sort(support_point_list.begin(),
1840 support_point_list.end(),
1843 new_indices[support_point_list[i].
second] = i;
1849 template <
int dim,
int spacedim>
1852 const unsigned int level,
1854 const bool dof_wise_renumbering)
1856 std::vector<types::global_dof_index> renumbering(dof.
n_dofs(
level));
1857 std::vector<types::global_dof_index> reverse(dof.
n_dofs(
level));
1859 renumbering, reverse, dof,
level, direction, dof_wise_renumbering);
1866 template <
int dim,
int spacedim>
1869 std::vector<types::global_dof_index> &reverse,
1871 const unsigned int level,
1873 const bool dof_wise_renumbering)
1875 if (dof_wise_renumbering ==
false)
1877 std::vector<typename DoFHandler<dim, spacedim>::level_cell_iterator>
1883 comparator(direction);
1892 ordered_cells.push_back(p);
1895 std::sort(ordered_cells.begin(), ordered_cells.end(), comparator);
1904 std::vector<std::pair<Point<spacedim>,
unsigned int>>
1905 support_point_list(n_dofs);
1912 std::vector<bool> already_touched(dof.
n_dofs(),
false);
1915 std::vector<types::global_dof_index> local_dof_indices(dofs_per_cell);
1923 &begin_tria =
begin;
1924 begin->get_active_or_mg_dof_indices(local_dof_indices);
1925 fe_values.
reinit(begin_tria);
1926 const std::vector<Point<spacedim>> &points =
1928 for (
unsigned int i = 0; i < dofs_per_cell; ++i)
1929 if (!already_touched[local_dof_indices[i]])
1931 support_point_list[local_dof_indices[i]].first = points[i];
1932 support_point_list[local_dof_indices[i]].second =
1933 local_dof_indices[i];
1934 already_touched[local_dof_indices[i]] =
true;
1939 std::sort(support_point_list.begin(),
1940 support_point_list.end(),
1943 new_indices[support_point_list[i].
second] = i;
1977 template <
class DHCellIterator>
1979 operator()(
const DHCellIterator &c1,
const DHCellIterator &c2)
const
1983 return compare(c1, c2, std::integral_constant<int, dim>());
1990 template <
class DHCellIterator,
int xdim>
1993 const DHCellIterator &c2,
1994 std::integral_constant<int, xdim>)
const
2000 return (
counter ? (s1 > s2) : (s2 > s1));
2008 template <
class DHCellIterator>
2011 const DHCellIterator &,
2012 std::integral_constant<int, 1>)
const
2015 ExcMessage(
"This operation only makes sense for dim>=2."));
2023 template <
int dim,
int spacedim>
2029 std::vector<types::global_dof_index> renumbering(dof.
n_dofs());
2037 template <
int dim,
int spacedim>
2044 std::vector<typename DoFHandler<dim, spacedim>::active_cell_iterator>
2050 ordered_cells.push_back(cell);
2052 std::sort(ordered_cells.begin(), ordered_cells.end(), comparator);
2054 std::vector<types::global_dof_index> reverse(new_indices.size());
2060 template <
int dim,
int spacedim>
2063 const unsigned int level,
2067 std::vector<typename DoFHandler<dim, spacedim>::level_cell_iterator>
2079 ordered_cells.push_back(p);
2082 std::sort(ordered_cells.begin(), ordered_cells.end(), comparator);
2089 template <
int dim,
int spacedim>
2093 std::vector<types::global_dof_index> renumbering(
2102 template <
int dim,
int spacedim>
2109 std::vector<types::global_dof_index> renumbering(
2120 template <
int dim,
int spacedim>
2126 Assert(new_indices.size() == n_dofs,
2129 std::iota(new_indices.begin(),
2137 ::boost::mt19937 random_number_generator;
2138 for (
unsigned int i = 1; i < n_dofs; ++i)
2141 const unsigned int j =
2142 ::boost::random::uniform_int_distribution<>(0, i)(
2143 random_number_generator);
2147 std::swap(new_indices[i], new_indices[j]);
2153 template <
int dim,
int spacedim>
2157 const unsigned int level)
2160 Assert(new_indices.size() == n_dofs,
2163 std::iota(new_indices.begin(),
2171 ::boost::mt19937 random_number_generator;
2172 for (
unsigned int i = 1; i < n_dofs; ++i)
2175 const unsigned int j =
2176 ::boost::random::uniform_int_distribution<>(0, i)(
2177 random_number_generator);
2181 std::swap(new_indices[i], new_indices[j]);
2187 template <
int dim,
int spacedim>
2195 "Parallel triangulations are already enumerated according to their MPI process id."));
2197 std::vector<types::global_dof_index> renumbering(
2206 template <
int dim,
int spacedim>
2212 Assert(new_dof_indices.size() == n_dofs,
2218 std::vector<types::subdomain_id> subdomain_association(n_dofs);
2220 const unsigned int n_subdomains =
2221 *std::max_element(subdomain_association.begin(),
2222 subdomain_association.end()) +
2232 std::fill(new_dof_indices.begin(),
2233 new_dof_indices.end(),
2239 if (subdomain_association[i] == subdomain)
2243 new_dof_indices[i] = next_free_index;
2249 Assert(std::find(new_dof_indices.begin(),
2250 new_dof_indices.end(),
2257 template <
int dim,
int spacedim>
2261 std::vector<types::global_dof_index> renumbering(
2274 template <
int dim,
int spacedim>
2277 std::vector<types::global_dof_index> &new_dof_indices,
2281 Assert(new_dof_indices.size() == n_dofs,
2293 std::vector<types::global_dof_index> component_renumbering(
2295 compute_component_wise<dim, spacedim>(component_renumbering,
2298 std::vector<unsigned int>(),
2303 const std::vector<types::global_dof_index> dofs_per_component =
2305 for (
const auto &dpc : dofs_per_component)
2309 const unsigned int n_components =
2313 dof_handler.
n_dofs() / n_components;
2319 if (component_renumbering.size() == 0)
2321 new_dof_indices.resize(0);
2324 std::fill(new_dof_indices.begin(),
2325 new_dof_indices.end(),
2333 std::vector<types::global_dof_index> component_dofs(
2334 local_dofs_per_component);
2335 for (
unsigned int component = 0; component < n_components; ++component)
2337 for (std::size_t i = 0; i < local_dofs_per_component; ++i)
2339 component_renumbering[n_components * i + component];
2340 component_renumbered_dofs.
add_indices(component_dofs.begin(),
2341 component_dofs.end());
2343 component_renumbered_dofs.
compress();
2347 component_renumbered_dofs2.
add_indices(component_renumbering.begin(),
2348 component_renumbering.end());
2349 Assert(component_renumbered_dofs2 == component_renumbered_dofs,
2356 AssertThrow(fe.dofs_per_cell == 0 || fe.has_support_points(),
2358 for (
unsigned int i = 0; i < fe.n_base_elements(); ++i)
2360 fe.base_element(0).get_unit_support_points() ==
2361 fe.base_element(i).get_unit_support_points(),
2363 "All base elements should have the same support points."));
2366 std::vector<types::global_dof_index> component_to_nodal(
2370 std::vector<types::global_dof_index> cell_dofs;
2371 std::vector<types::global_dof_index> component_renumbered_cell_dofs;
2375 auto next_dof_it = locally_owned_dofs.
begin();
2377 if (cell->is_locally_owned())
2382 cell->get_dof_indices(cell_dofs);
2389 if (locally_owned_dofs.
is_element(cell_dofs[i]))
2391 const auto local_index =
2393 component_renumbered_cell_dofs[i] =
2394 component_renumbering[local_index];
2398 component_renumbered_cell_dofs[i] =
2407 component_renumbered_cell_dofs[i]))
2409 for (
unsigned int component = 0;
2416 const auto local_index =
2418 component_renumbered_cell_dofs[i] +
2419 dofs_per_component * component);
2421 if (component_to_nodal[local_index] ==
2424 component_to_nodal[local_index] = *next_dof_it;
2435 const auto local_index =
2437 new_dof_indices[i] = component_to_nodal[local_index];
2446 typename VectorizedArrayType>
2452 const std::vector<types::global_dof_index> new_global_numbers =
2462 template <
int dim,
int spacedim,
typename Number,
typename AdditionalDataType>
2466 const AdditionalDataType & matrix_free_data)
2468 const std::vector<types::global_dof_index> new_global_numbers =
2475 dof_handler.
renumber_dofs(matrix_free_data.mg_level, new_global_numbers);
2480 template <
int dim,
int spacedim,
typename Number,
typename AdditionalDataType>
2481 std::vector<types::global_dof_index>
2485 const AdditionalDataType & matrix_free_data)
2487 AdditionalDataType my_mf_data = matrix_free_data;
2488 my_mf_data.initialize_mapping =
false;
2491 typename AdditionalDataType::MatrixFreeType separate_matrix_free;
2510 std::vector<std::vector<unsigned int>>
2511 group_dofs_by_rank_access(
2516 std::vector<unsigned int> touch_count(partitioner.locally_owned_size());
2517 for (
const auto &p : partitioner.import_indices())
2518 for (
unsigned int i = p.first; i < p.second; ++i)
2522 std::vector<std::vector<unsigned int>> result(1);
2523 for (
unsigned int i = 0; i < touch_count.size(); ++i)
2524 if (touch_count[i] == 0)
2525 result.back().push_back(i);
2530 std::map<unsigned int, std::vector<unsigned int>>
2531 multiple_ranks_access_dof;
2532 const std::vector<std::pair<unsigned int, unsigned int>> &import_targets =
2533 partitioner.import_targets();
2534 auto it = partitioner.import_indices().begin();
2535 for (
const std::pair<unsigned int, unsigned int> &proc : import_targets)
2537 result.emplace_back();
2538 unsigned int count_dofs = 0;
2539 while (count_dofs < proc.second)
2541 for (
unsigned int i = it->first; i < it->
second;
2544 if (touch_count[i] == 1)
2545 result.back().push_back(i);
2547 multiple_ranks_access_dof[i].push_back(proc.first);
2557 std::map<std::vector<unsigned int>,
2558 std::vector<unsigned int>,
2559 std::function<
bool(
const std::vector<unsigned int> &,
2560 const std::vector<unsigned int> &)>>
2561 dofs_by_rank{[](
const std::vector<unsigned int> &a,
2562 const std::vector<unsigned int> &
b) {
2563 if (a.size() <
b.size())
2565 if (a.size() ==
b.size())
2567 for (
unsigned int i = 0; i < a.size(); ++i)
2570 else if (a[i] >
b[i])
2575 for (
const auto &entry : multiple_ranks_access_dof)
2576 dofs_by_rank[entry.second].push_back(entry.first);
2578 for (
const auto &procs : dofs_by_rank)
2579 result.push_back(procs.second);
2590 template <
int dim,
typename Number,
typename VectorizedArrayType>
2591 std::pair<std::vector<unsigned int>, std::vector<unsigned char>>
2592 compute_mf_numbering(
2594 const unsigned int component)
2598 const unsigned int n_comp =
2601 matrix_free.
get_dof_handler(component).get_fe().n_base_elements() == 1,
2609 const unsigned int fe_degree =
2611 const unsigned int nn = fe_degree - 1;
2620 std::array<std::pair<unsigned int, unsigned int>,
2621 ::Utilities::pow(3, dim)>
2625 dofs_on_objects[0] = std::make_pair(0
U, 1U);
2626 dofs_on_objects[1] = std::make_pair(2 * n_comp, nn);
2627 dofs_on_objects[2] = std::make_pair(n_comp, 1U);
2631 dofs_on_objects[0] = std::make_pair(0
U, 1U);
2632 dofs_on_objects[1] = std::make_pair(n_comp * (4 + 2 * nn), nn);
2633 dofs_on_objects[2] = std::make_pair(n_comp, 1U);
2634 dofs_on_objects[3] = std::make_pair(n_comp * 4, nn);
2635 dofs_on_objects[4] = std::make_pair(n_comp * (4 + 4 * nn), nn * nn);
2636 dofs_on_objects[5] = std::make_pair(n_comp * (4 + 1 * nn), nn);
2637 dofs_on_objects[6] = std::make_pair(2 * n_comp, 1U);
2638 dofs_on_objects[7] = std::make_pair(n_comp * (4 + 3 * nn), nn);
2639 dofs_on_objects[8] = std::make_pair(3 * n_comp, 1U);
2643 dofs_on_objects[0] = std::make_pair(0
U, 1U);
2644 dofs_on_objects[1] = std::make_pair(n_comp * (8 + 2 * nn), nn);
2645 dofs_on_objects[2] = std::make_pair(n_comp, 1U);
2646 dofs_on_objects[3] = std::make_pair(n_comp * 8, nn);
2647 dofs_on_objects[4] =
2648 std::make_pair(n_comp * (8 + 12 * nn + 4 * nn * nn), nn * nn);
2649 dofs_on_objects[5] = std::make_pair(n_comp * (8 + 1 * nn), nn);
2650 dofs_on_objects[6] = std::make_pair(n_comp * 2, 1U);
2651 dofs_on_objects[7] = std::make_pair(n_comp * (8 + 3 * nn), nn);
2652 dofs_on_objects[8] = std::make_pair(n_comp * 3, 1U);
2653 dofs_on_objects[9] = std::make_pair(n_comp * (8 + 8 * nn), nn);
2654 dofs_on_objects[10] =
2655 std::make_pair(n_comp * (8 + 12 * nn + 2 * nn * nn), nn * nn);
2656 dofs_on_objects[11] = std::make_pair(n_comp * (8 + 9 * nn), nn);
2657 dofs_on_objects[12] = std::make_pair(n_comp * (8 + 12 * nn), nn * nn);
2658 dofs_on_objects[13] =
2659 std::make_pair(n_comp * (8 + 12 * nn + 6 * nn * nn), nn * nn * nn);
2660 dofs_on_objects[14] =
2661 std::make_pair(n_comp * (8 + 12 * nn + 1 * nn * nn), nn * nn);
2662 dofs_on_objects[15] = std::make_pair(n_comp * (8 + 10 * nn), nn);
2663 dofs_on_objects[16] =
2664 std::make_pair(n_comp * (8 + 12 * nn + 3 * nn * nn), nn * nn);
2665 dofs_on_objects[17] = std::make_pair(n_comp * (8 + 11 * nn), nn);
2666 dofs_on_objects[18] = std::make_pair(n_comp * 4, 1U);
2667 dofs_on_objects[19] = std::make_pair(n_comp * (8 + 6 * nn), nn);
2668 dofs_on_objects[20] = std::make_pair(n_comp * 5, 1U);
2669 dofs_on_objects[21] = std::make_pair(n_comp * (8 + 4 * nn), nn);
2670 dofs_on_objects[22] =
2671 std::make_pair(n_comp * (8 + 12 * nn + 5 * nn * nn), nn * nn);
2672 dofs_on_objects[23] = std::make_pair(n_comp * (8 + 5 * nn), nn);
2673 dofs_on_objects[24] = std::make_pair(n_comp * 6, 1U);
2674 dofs_on_objects[25] = std::make_pair(n_comp * (8 + 7 * nn), nn);
2675 dofs_on_objects[26] = std::make_pair(n_comp * 7, 1U);
2680 std::vector<unsigned int> & result,
2681 unsigned int &counter_dof_numbers) {
2682 if (owned_dofs.is_element(dof_index))
2684 const unsigned int local_dof_index =
2685 owned_dofs.index_within_set(dof_index);
2688 result[local_dof_index] = counter_dof_numbers++;
2692 unsigned int counter_dof_numbers = 0;
2693 std::vector<unsigned int> local_dofs, dofs_extracted;
2694 std::vector<types::global_dof_index> dof_indices(
2711 std::vector<unsigned int> dof_numbers_mf_order(
2713 std::vector<unsigned char> touch_count(owned_dofs.n_elements());
2715 const auto operation_on_cell_range =
2718 const unsigned int &,
2719 const std::pair<unsigned int, unsigned int> &cell_range) {
2721 for (
unsigned int cell = cell_range.first; cell < cell_range.second;
2725 for (
unsigned int v = 0;
2732 ->get_dof_indices(dof_indices);
2735 ->get_mg_dof_indices(dof_indices);
2737 for (
unsigned int a = 0; a < dofs_on_objects.size(); ++a)
2739 const auto &r = dofs_on_objects[a];
2740 if (a == 10 || a == 16)
2743 for (
unsigned int i1 = 0; i1 < nn; ++i1)
2744 for (
unsigned int i0 = 0; i0 < nn; ++i0)
2745 for (
unsigned int c = 0; c < n_comp; ++c)
2746 renumber_func(dof_indices[r.first + r.second * c +
2749 dof_numbers_mf_order,
2750 counter_dof_numbers);
2752 for (
unsigned int i = 0; i < r.second; ++i)
2753 for (
unsigned int c = 0; c < n_comp; ++c)
2755 dof_indices[r.first + r.second * c + i],
2757 dof_numbers_mf_order,
2758 counter_dof_numbers);
2765 dofs_extracted, cell);
2766 local_dofs.insert(local_dofs.end(),
2767 dofs_extracted.begin(),
2768 dofs_extracted.end());
2771 std::sort(local_dofs.begin(), local_dofs.end());
2772 local_dofs.erase(std::unique(local_dofs.begin(), local_dofs.end()),
2775 for (
unsigned int i : local_dofs)
2776 if (i < touch_count.size())
2784 "version of MatrixFree::cell_loop"));
2786 matrix_free.template cell_loop<unsigned int, unsigned int>(
2787 operation_on_cell_range, counter_dof_numbers, counter_dof_numbers);
2791 return std::make_pair(dof_numbers_mf_order, touch_count);
2801 typename VectorizedArrayType>
2802 std::vector<types::global_dof_index>
2808 ExcMessage(
"You need to set up indices in MatrixFree "
2809 "to be able to compute a renumbering!"));
2812 unsigned int component = 0;
2813 for (; component < matrix_free.
n_components(); ++component)
2818 ExcMessage(
"Could not locate the given DoFHandler in MatrixFree"));
2831 const std::vector<std::vector<unsigned int>> dofs_by_rank_access =
2832 group_dofs_by_rank_access(
2835 const std::pair<std::vector<unsigned int>, std::vector<unsigned char>>
2836 local_numbering = compute_mf_numbering(matrix_free, component);
2841 std::vector<unsigned int> new_numbers;
2842 new_numbers.reserve(owned_dofs.
n_elements());
2848 const std::vector<unsigned int> &purely_local_dofs = dofs_by_rank_access[0];
2849 for (
unsigned int i : purely_local_dofs)
2850 if (local_numbering.second[i] == 1)
2851 new_numbers.push_back(i);
2853 const auto comparator = [&](
const unsigned int a,
const unsigned int b) {
2854 return (local_numbering.first[a] < local_numbering.first[
b]);
2857 std::sort(new_numbers.begin(), new_numbers.end(), comparator);
2858 unsigned int sorted_size = new_numbers.size();
2862 for (
auto i : purely_local_dofs)
2863 if (local_numbering.second[i] > 1)
2864 new_numbers.push_back(i);
2865 std::sort(new_numbers.begin() + sorted_size, new_numbers.end(), comparator);
2866 sorted_size = new_numbers.size();
2869 for (
unsigned int chunk = 1; chunk < dofs_by_rank_access.size(); ++chunk)
2870 for (
auto i : dofs_by_rank_access[chunk])
2871 new_numbers.push_back(i);
2872 std::sort(new_numbers.begin() + sorted_size, new_numbers.end(), comparator);
2873 sorted_size = new_numbers.size();
2876 for (
auto i : purely_local_dofs)
2877 if (local_numbering.second[i] == 0)
2878 new_numbers.push_back(i);
2879 std::sort(new_numbers.begin() + sorted_size, new_numbers.end(), comparator);
2883 std::vector<::types::global_dof_index> new_global_numbers(
2885 for (
unsigned int i = 0; i < new_numbers.size(); ++i)
2888 return new_global_numbers;
2896 #include "dof_renumbering.inst"
void reinit(const IndexSet &local_constraints=IndexSet())
unsigned int first_selected_component(const unsigned int overall_number_of_components=numbers::invalid_unsigned_int) const
cell_iterator end() const
void renumber_dofs(const std::vector< types::global_dof_index > &new_numbers)
const Triangulation< dim, spacedim > & get_triangulation() const
active_cell_iterator begin_active(const unsigned int level=0) const
types::global_dof_index n_dofs() const
const FiniteElement< dim, spacedim > & get_fe(const unsigned int index=0) const
const IndexSet & locally_owned_mg_dofs(const unsigned int level) const
cell_iterator begin(const unsigned int level=0) const
MPI_Comm get_communicator() const
const IndexSet & locally_owned_dofs() const
const hp::FECollection< dim, spacedim > & get_fe_collection() const
types::global_dof_index n_locally_owned_dofs() const
const std::vector< Point< spacedim > > & get_quadrature_points() const
void reinit(const TriaIterator< DoFCellAccessor< dim, spacedim, level_dof_access >> &cell)
unsigned int n_dofs_per_cell() const
unsigned int n_components() const
const unsigned int dofs_per_cell
std::pair< unsigned int, unsigned int > system_to_component_index(const unsigned int index) const
const std::vector< Point< dim > > & get_unit_support_points() const
bool has_support_points() const
std::pair< unsigned int, types::global_dof_index > system_to_block_index(const unsigned int component) const
const ComponentMask & get_nonzero_components(const unsigned int i) const
bool is_primitive() const
size_type index_within_set(const size_type global_index) const
size_type n_elements() const
bool is_element(const size_type index) const
ElementIterator begin() const
void subtract_set(const IndexSet &other)
size_type nth_index_in_set(const size_type local_index) const
void add_indices(const ForwardIterator &begin, const ForwardIterator &end)
unsigned int n_active_entries_per_cell_batch(const unsigned int cell_batch_index) const
unsigned int get_mg_level() const
const internal::MatrixFreeFunctions::DoFInfo & get_dof_info(const unsigned int dof_handler_index_component=0) const
bool indices_initialized() const
DoFHandler< dim >::cell_iterator get_cell_iterator(const unsigned int cell_batch_index, const unsigned int lane_index, const unsigned int dof_handler_index=0) const
unsigned int n_components() const
const DoFHandler< dim > & get_dof_handler(const unsigned int dof_handler_index=0) const
const internal::MatrixFreeFunctions::TaskInfo & get_task_info() const
virtual MPI_Comm get_communicator() const
unsigned int n_active_cells() const
unsigned int n_cells() const
unsigned int size() const
void push_back(const FiniteElement< dim, spacedim > &new_fe)
unsigned int n_blocks() const
unsigned int n_components() const
const FEValuesType & get_present_fe_values() const
void reinit(const TriaIterator< DoFCellAccessor< dim, spacedim, lda >> &cell, const unsigned int q_index=numbers::invalid_unsigned_int, const unsigned int mapping_index=numbers::invalid_unsigned_int, const unsigned int fe_index=numbers::invalid_unsigned_int)
void push_back(const Quadrature< dim_in > &new_quadrature)
#define DEAL_II_NAMESPACE_OPEN
#define DEAL_II_DISABLE_EXTRA_DIAGNOSTICS
#define DEAL_II_NAMESPACE_CLOSE
#define DEAL_II_ENABLE_EXTRA_DIAGNOSTICS
@ update_quadrature_points
Transformed quadrature points.
IteratorRange< active_cell_iterator > active_cell_iterators() const
static ::ExceptionBase & ExcInternalError()
static ::ExceptionBase & ExcDimensionMismatch(std::size_t arg1, std::size_t arg2)
#define Assert(cond, exc)
std::string to_string(const T &t)
static ::ExceptionBase & ExcNotImplemented()
static ::ExceptionBase & ExcDoFHandlerNotInitialized()
#define AssertDimension(dim1, dim2)
#define AssertThrowMPI(error_code)
#define AssertIndexRange(index, range)
static ::ExceptionBase & ExcMessage(std::string arg1)
#define AssertThrow(cond, exc)
void add_entries(const size_type row, ForwardIterator begin, ForwardIterator end, const bool indices_are_unique_and_sorted=false)
size_type row_length(const size_type row) const
size_type column_number(const size_type row, const size_type index) const
void make_hanging_node_constraints(const DoFHandler< dim, spacedim > &dof_handler, AffineConstraints< number > &constraints)
void make_sparsity_pattern(const DoFHandler< dim, spacedim > &dof_handler, SparsityPatternType &sparsity_pattern, const AffineConstraints< number > &constraints=AffineConstraints< number >(), const bool keep_constrained_dofs=true, const types::subdomain_id subdomain_id=numbers::invalid_subdomain_id)
Expression atan2(const Expression &y, const Expression &x)
graph_traits< Graph >::vertices_size_type size_type
graph_traits< Graph >::vertex_descriptor Vertex
adjacency_list< vecS, vecS, undirectedS, property< vertex_color_t, default_color_type, property< vertex_degree_t, int > >> Graph
std::pair< size_type, size_type > Pair
void create_graph(const DoFHandler< dim, spacedim > &dof_handler, const bool use_constraints, boosttypes::Graph &graph, boosttypes::property_map< boosttypes::Graph, boosttypes::vertex_degree_t >::type &graph_degree)
void compute_Cuthill_McKee(std::vector< types::global_dof_index > &new_dof_indices, const DoFHandler< dim, spacedim > &, const bool reversed_numbering=false, const bool use_constraints=false)
void Cuthill_McKee(DoFHandler< dim, spacedim > &dof_handler, const bool reversed_numbering=false, const bool use_constraints=false)
void king_ordering(DoFHandler< dim, spacedim > &dof_handler, const bool reversed_numbering=false, const bool use_constraints=false)
void compute_king_ordering(std::vector< types::global_dof_index > &new_dof_indices, const DoFHandler< dim, spacedim > &, const bool reversed_numbering=false, const bool use_constraints=false)
void compute_minimum_degree(std::vector< types::global_dof_index > &new_dof_indices, const DoFHandler< dim, spacedim > &, const bool reversed_numbering=false, const bool use_constraints=false)
void minimum_degree(DoFHandler< dim, spacedim > &dof_handler, const bool reversed_numbering=false, const bool use_constraints=false)
void compute_support_point_wise(std::vector< types::global_dof_index > &new_dof_indices, const DoFHandler< dim, spacedim > &dof_handler)
void compute_subdomain_wise(std::vector< types::global_dof_index > &new_dof_indices, const DoFHandler< dim, spacedim > &dof_handler)
void matrix_free_data_locality(DoFHandler< dim, spacedim > &dof_handler, const MatrixFree< dim, Number, VectorizedArrayType > &matrix_free)
void subdomain_wise(DoFHandler< dim, spacedim > &dof_handler)
void hierarchical(DoFHandler< dim, spacedim > &dof_handler)
void compute_random(std::vector< types::global_dof_index > &new_dof_indices, const DoFHandler< dim, spacedim > &dof_handler)
void component_wise(DoFHandler< dim, spacedim > &dof_handler, const std::vector< unsigned int > &target_component=std::vector< unsigned int >())
void downstream(DoFHandler< dim, spacedim > &dof_handler, const Tensor< 1, spacedim > &direction, const bool dof_wise_renumbering=false)
types::global_dof_index compute_component_wise(std::vector< types::global_dof_index > &new_dof_indices, const CellIterator &start, const typename identity< CellIterator >::type &end, const std::vector< unsigned int > &target_component, const bool is_level_operation)
void compute_Cuthill_McKee(std::vector< types::global_dof_index > &new_dof_indices, const DoFHandler< dim, spacedim > &, const bool reversed_numbering=false, const bool use_constraints=false, const std::vector< types::global_dof_index > &starting_indices=std::vector< types::global_dof_index >(), const unsigned int level=numbers::invalid_unsigned_int)
void block_wise(DoFHandler< dim, spacedim > &dof_handler)
void Cuthill_McKee(DoFHandler< dim, spacedim > &dof_handler, const bool reversed_numbering=false, const bool use_constraints=false, const std::vector< types::global_dof_index > &starting_indices=std::vector< types::global_dof_index >())
void compute_sort_selected_dofs_back(std::vector< types::global_dof_index > &new_dof_indices, const DoFHandler< dim, spacedim > &dof_handler, const std::vector< bool > &selected_dofs)
void sort_selected_dofs_back(DoFHandler< dim, spacedim > &dof_handler, const std::vector< bool > &selected_dofs)
void support_point_wise(DoFHandler< dim, spacedim > &dof_handler)
void compute_cell_wise(std::vector< types::global_dof_index > &renumbering, std::vector< types::global_dof_index > &inverse_renumbering, const DoFHandler< dim, spacedim > &dof_handler, const std::vector< typename DoFHandler< dim, spacedim >::active_cell_iterator > &cell_order)
std::vector< types::global_dof_index > compute_matrix_free_data_locality(const DoFHandler< dim, spacedim > &dof_handler, const MatrixFree< dim, Number, VectorizedArrayType > &matrix_free)
void clockwise_dg(DoFHandler< dim, spacedim > &dof_handler, const Point< spacedim > ¢er, const bool counter=false)
void random(DoFHandler< dim, spacedim > &dof_handler)
void compute_downstream(std::vector< types::global_dof_index > &new_dof_indices, std::vector< types::global_dof_index > &reverse, const DoFHandler< dim, spacedim > &dof_handler, const Tensor< 1, spacedim > &direction, const bool dof_wise_renumbering)
void cell_wise(DoFHandler< dim, spacedim > &dof_handler, const std::vector< typename DoFHandler< dim, spacedim >::active_cell_iterator > &cell_order)
types::global_dof_index compute_block_wise(std::vector< types::global_dof_index > &new_dof_indices, const ITERATOR &start, const ENDITERATOR &end, bool is_level_operation)
void compute_clockwise_dg(std::vector< types::global_dof_index > &new_dof_indices, const DoFHandler< dim, spacedim > &dof_handler, const Point< spacedim > ¢er, const bool counter)
void swap(MemorySpaceData< Number, MemorySpace > &, MemorySpaceData< Number, MemorySpace > &)
SymmetricTensor< 2, dim, Number > b(const Tensor< 2, dim, Number > &F)
VectorType::value_type * begin(VectorType &V)
VectorType::value_type * end(VectorType &V)
TrilinosWrappers::types::int_type global_index(const Epetra_BlockMap &map, const ::types::global_dof_index i)
T sum(const T &t, const MPI_Comm &mpi_communicator)
T max(const T &t, const MPI_Comm &mpi_communicator)
std::vector< Integer > reverse_permutation(const std::vector< Integer > &permutation)
void copy(const T *begin, const T *end, U *dest)
static const unsigned int invalid_unsigned_int
const types::global_dof_index invalid_dof_index
unsigned int global_dof_index
unsigned int subdomain_id
bool compare(const DHCellIterator &c1, const DHCellIterator &c2, std::integral_constant< int, xdim >) const
bool compare(const DHCellIterator &, const DHCellIterator &, std::integral_constant< int, 1 >) const
ClockCells(const Point< dim > ¢er, bool counter)
const Point< dim > & center
bool operator()(const DHCellIterator &c1, const DHCellIterator &c2) const
void get_dof_indices_on_cell_batch(std::vector< unsigned int > &local_indices, const unsigned int cell_batch, const bool with_constraints=true) const
std::shared_ptr< const Utilities::MPI::Partitioner > vector_partitioner
TasksParallelScheme scheme
const ::Triangulation< dim, spacedim > & tria
#define DEAL_II_DOF_INDEX_MPI_TYPE