76 template <
int dim,
typename Number =
double>
88 double downstream_size = 0;
90 for (
unsigned int d = 0; d < dim; ++d)
92 downstream_size += (rhs[d] - lhs[d]) * weight;
95 if (downstream_size < 0)
97 else if (downstream_size > 0)
101 for (
unsigned int d = 0; d < dim; ++d)
103 if (lhs[d] == rhs[d])
105 return lhs[d] < rhs[d];
125 template <
int dim,
int spacedim>
126 std::vector<unsigned char>
130 std::vector<unsigned char> local_component_association(
139 local_component_association[i] =
150 const unsigned int first_comp =
155 local_component_association[i] = first_comp;
162 for (
unsigned int c = first_comp; c < fe.
n_components(); ++c)
163 if (component_mask[c] ==
true)
165 local_component_association[i] = c;
170 Assert(std::find(local_component_association.begin(),
171 local_component_association.end(),
172 static_cast<unsigned char>(-1)) ==
173 local_component_association.end(),
176 return local_component_association;
196 template <
int dim,
int spacedim>
200 std::vector<unsigned char> &dofs_by_component)
202 const ::hp::FECollection<dim, spacedim> &fe_collection =
216 std::vector<std::vector<unsigned char>> local_component_association(
217 fe_collection.size());
218 for (
unsigned int f = 0; f < fe_collection.size(); ++f)
221 local_component_association[f] =
226 std::vector<types::global_dof_index> indices;
230 const unsigned int fe_index = c->active_fe_index();
231 const unsigned int dofs_per_cell = c->get_fe().n_dofs_per_cell();
232 indices.resize(dofs_per_cell);
233 c->get_dof_indices(indices);
234 for (
unsigned int i = 0; i < dofs_per_cell; ++i)
237 indices[i])] = local_component_association[fe_index][i];
249 template <
int dim,
int spacedim>
252 std::vector<unsigned char> & dofs_by_block)
254 const ::hp::FECollection<dim, spacedim> &fe_collection =
268 std::vector<std::vector<unsigned char>> local_block_association(
269 fe_collection.size());
270 for (
unsigned int f = 0; f < fe_collection.size(); ++f)
274 static_cast<unsigned char>(-1));
278 Assert(std::find(local_block_association[f].begin(),
279 local_block_association[f].end(),
280 static_cast<unsigned char>(-1)) ==
281 local_block_association[f].end(),
286 std::vector<types::global_dof_index> indices;
288 if (cell->is_locally_owned())
290 const unsigned int fe_index = cell->active_fe_index();
291 const unsigned int dofs_per_cell = cell->get_fe().n_dofs_per_cell();
292 indices.resize(dofs_per_cell);
293 cell->get_dof_indices(indices);
294 for (
unsigned int i = 0; i < dofs_per_cell; ++i)
297 indices[i])] = local_block_association[fe_index][i];
304 template <
int dim,
int spacedim,
typename Number>
309 const unsigned int component)
318 for (
unsigned int i = 0; i < fe_collection.size(); ++i)
320 Assert(fe_collection[i].is_primitive() ==
true,
327 const bool consider_components =
331 if (consider_components ==
false)
335 std::vector<unsigned char> component_dofs(
342 for (
unsigned int i = 0; i < dof_data.
size(); ++i)
343 if (component_dofs[i] ==
static_cast<unsigned char>(component))
348 std::vector<unsigned char> touch_count(dof_handler.
n_dofs(), 0);
350 std::vector<types::global_dof_index> dof_indices;
351 dof_indices.reserve(fe_collection.max_dofs_per_cell());
355 const unsigned int dofs_per_cell = cell->get_fe().n_dofs_per_cell();
356 dof_indices.resize(dofs_per_cell);
357 cell->get_dof_indices(dof_indices);
359 for (
unsigned int i = 0; i < dofs_per_cell; ++i)
362 if (!consider_components ||
363 (cell->get_fe().system_to_component_index(i).first == component))
366 dof_data(dof_indices[i]) += cell_data(cell->active_cell_index());
369 ++touch_count[dof_indices[i]];
379 Assert(consider_components || (touch_count[i] != 0),
381 if (touch_count[i] != 0)
382 dof_data(i) /= touch_count[i];
388 template <
int dim,
int spacedim>
396 "The given component mask is not sized correctly to represent the "
397 "components of the given finite element."));
415 std::vector<types::global_dof_index> selected_dofs;
418 if (component_mask[dofs_by_component[i]] ==
true)
423 result.
add_indices(selected_dofs.begin(), selected_dofs.end());
429 template <
int dim,
int spacedim>
435 return extract_dofs<dim, spacedim>(
441 template <
int dim,
int spacedim>
442 std::vector<IndexSet>
449 "The given component mask is not sized correctly to represent the "
450 "components of the given finite element."));
459 std::vector<IndexSet> index_per_comp(n_comps,
IndexSet(dof.
n_dofs()));
463 const auto &comp_i = dofs_by_component[i];
464 if (component_mask[comp_i])
465 index_per_comp[comp_i].add_index(
466 locally_owned_dofs.nth_index_in_set(i));
468 for (
auto &c : index_per_comp)
470 return index_per_comp;
475 template <
int dim,
int spacedim>
480 std::vector<bool> & selected_dofs)
487 "The given component mask is not sized correctly to represent the "
488 "components of the given finite element."));
497 std::fill_n(selected_dofs.begin(), dof.
n_dofs(
level),
false);
504 std::fill_n(selected_dofs.begin(), dof.
n_dofs(
level),
true);
509 std::fill_n(selected_dofs.begin(), dof.
n_dofs(
level),
false);
513 std::vector<unsigned char> local_component_asssociation =
517 local_selected_dofs[i] = component_mask[local_component_asssociation[i]];
523 cell->get_mg_dof_indices(indices);
525 selected_dofs[indices[i]] = local_selected_dofs[i];
531 template <
int dim,
int spacedim>
536 std::vector<bool> & selected_dofs)
547 template <
int dim,
int spacedim>
551 std::vector<bool> & selected_dofs,
552 const std::set<types::boundary_id> &boundary_ids)
558 "This function can not be used with distributed triangulations. "
559 "See the documentation for more information."));
565 selected_dofs.clear();
566 selected_dofs.resize(dof_handler.
n_dofs(),
false);
574 template <
int dim,
int spacedim>
579 const std::set<types::boundary_id> &boundary_ids)
588 template <
int dim,
int spacedim>
592 const std::set<types::boundary_id> &boundary_ids)
596 ExcMessage(
"Component mask has invalid size."));
605 const bool check_boundary_id = (boundary_ids.size() != 0);
609 const bool check_vector_component =
615 std::vector<types::global_dof_index> face_dof_indices;
616 face_dof_indices.reserve(
628 if (cell->is_artificial() ==
false)
629 for (
const unsigned int face : cell->face_indices())
630 if (cell->at_boundary(face))
631 if (!check_boundary_id ||
632 (boundary_ids.find(cell->face(face)->boundary_id()) !=
637 const auto reference_cell = cell->reference_cell();
639 const unsigned int n_vertices_per_cell =
640 reference_cell.n_vertices();
641 const unsigned int n_lines_per_cell = reference_cell.n_lines();
642 const unsigned int n_vertices_per_face =
643 reference_cell.face_reference_cell(face).n_vertices();
644 const unsigned int n_lines_per_face =
645 reference_cell.face_reference_cell(face).n_lines();
648 face_dof_indices.resize(dofs_per_face);
649 cell->face(face)->get_dof_indices(face_dof_indices,
650 cell->active_fe_index());
653 if (!check_vector_component)
654 selected_dofs.
add_index(face_dof_indices[i]);
669 (dim == 3 ? (i < n_vertices_per_face *
672 (i < n_vertices_per_face *
676 (i - n_vertices_per_face *
678 n_vertices_per_cell *
681 n_vertices_per_face *
685 n_vertices_per_cell *
695 selected_dofs.
add_index(face_dof_indices[i]);
699 const unsigned int first_nonzero_comp =
705 if (component_mask[first_nonzero_comp] ==
true)
706 selected_dofs.
add_index(face_dof_indices[i]);
711 return selected_dofs;
716 template <
int dim,
int spacedim>
721 std::vector<bool> & selected_dofs,
722 const std::set<types::boundary_id> &boundary_ids)
726 ExcMessage(
"This component mask has the wrong size."));
733 const bool check_boundary_id = (boundary_ids.size() != 0);
737 const bool check_vector_component =
741 selected_dofs.clear();
742 selected_dofs.resize(dof_handler.
n_dofs(),
false);
743 std::vector<types::global_dof_index> cell_dof_indices;
744 cell_dof_indices.reserve(
754 for (
const unsigned int face : cell->face_indices())
755 if (cell->at_boundary(face))
756 if (!check_boundary_id ||
757 (boundary_ids.find(cell->face(face)->boundary_id()) !=
763 cell_dof_indices.resize(dofs_per_cell);
764 cell->get_dof_indices(cell_dof_indices);
769 if (!check_vector_component)
770 selected_dofs[cell_dof_indices[i]] =
true;
777 selected_dofs[cell_dof_indices[i]] =
782 const unsigned int first_nonzero_comp =
788 selected_dofs[cell_dof_indices[i]] =
789 (component_mask[first_nonzero_comp] ==
true);
798 template <
int dim,
int spacedim,
typename number>
807 const std::function<
bool(
812 ->
bool {
return cell->is_locally_owned() && predicate(cell); };
814 std::vector<types::global_dof_index> local_dof_indices;
815 local_dof_indices.reserve(
819 std::set<types::global_dof_index> predicate_dofs;
822 if (!cell->is_artificial() && predicate(cell))
824 local_dof_indices.resize(cell->get_fe().n_dofs_per_cell());
825 cell->get_dof_indices(local_dof_indices);
826 predicate_dofs.insert(local_dof_indices.begin(),
827 local_dof_indices.end());
831 std::set<types::global_dof_index> dofs_with_support_on_halo_cells;
833 const std::vector<typename DoFHandler<dim, spacedim>::active_cell_iterator>
837 active_cell_iterator>::const_iterator it =
839 it != halo_cells.end();
851 const unsigned int dofs_per_cell = (*it)->get_fe().n_dofs_per_cell();
852 local_dof_indices.resize(dofs_per_cell);
853 (*it)->get_dof_indices(local_dof_indices);
854 dofs_with_support_on_halo_cells.insert(local_dof_indices.begin(),
855 local_dof_indices.end());
867 std::set<types::global_dof_index> extra_halo;
868 for (
const auto dof : dofs_with_support_on_halo_cells)
874 const unsigned int line_size = line_ptr->size();
875 for (
unsigned int j = 0; j < line_size; ++j)
876 extra_halo.insert((*line_ptr)[j].first);
879 dofs_with_support_on_halo_cells.insert(extra_halo.begin(),
886 support_set.
add_indices(predicate_dofs.begin(), predicate_dofs.end());
890 halo_set.
add_indices(dofs_with_support_on_halo_cells.begin(),
891 dofs_with_support_on_halo_cells.end());
906 template <
int spacedim>
909 const ::DoFHandler<1, spacedim> &dof_handler)
912 return IndexSet(dof_handler.n_dofs());
916 template <
int spacedim>
919 const ::DoFHandler<2, spacedim> &dof_handler)
921 const unsigned int dim = 2;
923 IndexSet selected_dofs(dof_handler.n_dofs());
929 for (
const auto &cell : dof_handler.active_cell_iterators())
930 if (!cell->is_artificial())
932 for (
const unsigned int face : cell->face_indices())
933 if (cell->face(face)->has_children())
935 const typename ::DoFHandler<dim,
936 spacedim>::line_iterator
937 line = cell->face(face);
941 selected_dofs.add_index(
942 line->child(0)->vertex_dof_index(1, dof));
944 for (
unsigned int child = 0; child < 2; ++child)
946 if (cell->neighbor_child_on_subface(face, child)
951 selected_dofs.add_index(
952 line->child(child)->dof_index(dof));
957 selected_dofs.compress();
958 return selected_dofs;
962 template <
int spacedim>
965 const ::DoFHandler<3, spacedim> &dof_handler)
967 const unsigned int dim = 3;
969 IndexSet selected_dofs(dof_handler.n_dofs());
970 IndexSet unconstrained_dofs(dof_handler.n_dofs());
974 for (
const auto &cell : dof_handler.active_cell_iterators())
975 if (!cell->is_artificial())
976 for (
auto f : cell->face_indices())
978 const typename ::DoFHandler<dim, spacedim>::face_iterator
979 face = cell->face(f);
980 if (cell->face(f)->has_children())
982 for (
unsigned int child = 0; child < 4; ++child)
983 if (!cell->neighbor_child_on_subface(f, child)
987 std::vector<types::global_dof_index> ldi(
989 face->child(child)->get_dof_indices(ldi);
990 selected_dofs.add_indices(ldi.begin(), ldi.end());
995 for (
unsigned int vertex = 0; vertex < 4; ++vertex)
998 unconstrained_dofs.add_index(
999 face->vertex_dof_index(vertex, dof));
1002 selected_dofs.subtract_set(unconstrained_dofs);
1003 return selected_dofs;
1010 template <
int dim,
int spacedim>
1014 return internal::extract_hanging_node_dofs(dof_handler);
1019 template <
int dim,
int spacedim>
1023 std::vector<bool> & selected_dofs)
1029 std::fill_n(selected_dofs.begin(), dof_handler.
n_dofs(),
false);
1031 std::vector<types::global_dof_index> local_dof_indices;
1032 local_dof_indices.reserve(
1038 if (cell->subdomain_id() == subdomain_id)
1040 const unsigned int dofs_per_cell = cell->get_fe().n_dofs_per_cell();
1041 local_dof_indices.resize(dofs_per_cell);
1042 cell->get_dof_indices(local_dof_indices);
1043 for (
unsigned int i = 0; i < dofs_per_cell; ++i)
1044 selected_dofs[local_dof_indices[i]] =
true;
1050 template <
int dim,
int spacedim>
1060 std::vector<types::global_dof_index> dof_indices;
1061 std::set<types::global_dof_index> global_dof_indices;
1066 dof_indices.resize(cell->get_fe().n_dofs_per_cell());
1067 cell->get_dof_indices(dof_indices);
1071 global_dof_indices.insert(dof_index);
1074 dof_set.
add_indices(global_dof_indices.begin(), global_dof_indices.end());
1083 template <
int dim,
int spacedim>
1093 template <
int dim,
int spacedim>
1097 const unsigned int level)
1105 std::vector<types::global_dof_index> dof_indices;
1106 std::set<types::global_dof_index> global_dof_indices;
1108 const auto filtered_iterators_range =
1111 for (
const auto &cell : filtered_iterators_range)
1113 dof_indices.resize(cell->get_fe().n_dofs_per_cell());
1114 cell->get_mg_dof_indices(dof_indices);
1118 global_dof_indices.insert(dof_index);
1121 dof_set.
add_indices(global_dof_indices.begin(), global_dof_indices.end());
1130 template <
int dim,
int spacedim>
1135 const unsigned int level)
1142 template <
int dim,
int spacedim>
1157 std::vector<types::global_dof_index> dof_indices;
1158 std::vector<types::global_dof_index> dofs_on_ghosts;
1161 if (cell->is_ghost())
1163 dof_indices.resize(cell->get_fe().n_dofs_per_cell());
1164 cell->get_dof_indices(dof_indices);
1165 for (
const auto dof_index : dof_indices)
1167 dofs_on_ghosts.push_back(dof_index);
1171 std::sort(dofs_on_ghosts.begin(), dofs_on_ghosts.end());
1173 std::unique(dofs_on_ghosts.begin(),
1174 dofs_on_ghosts.end()));
1182 template <
int dim,
int spacedim>
1192 template <
int dim,
int spacedim>
1196 const unsigned int level)
1209 std::vector<types::global_dof_index> dof_indices;
1210 std::vector<types::global_dof_index> dofs_on_ghosts;
1221 dof_indices.resize(cell->get_fe().n_dofs_per_cell());
1222 cell->get_mg_dof_indices(dof_indices);
1223 for (
const auto dof_index : dof_indices)
1225 dofs_on_ghosts.push_back(dof_index);
1229 std::sort(dofs_on_ghosts.begin(), dofs_on_ghosts.end());
1231 std::unique(dofs_on_ghosts.begin(),
1232 dofs_on_ghosts.end()));
1241 template <
int dim,
int spacedim>
1245 const unsigned int level,
1253 template <
int dim,
int spacedim>
1257 std::vector<std::vector<bool>> & constant_modes)
1263 constant_modes = std::vector<std::vector<bool>>(0);
1271 std::vector<unsigned char> dofs_by_component(
1276 unsigned int n_selected_dofs = 0;
1277 for (
unsigned int i = 0; i < n_components; ++i)
1278 if (component_mask[i] ==
true)
1280 std::count(dofs_by_component.begin(), dofs_by_component.end(), i);
1284 std::vector<unsigned int> component_numbering(
1286 for (
unsigned int i = 0, count = 0; i < locally_owned_dofs.
n_elements();
1288 if (component_mask[dofs_by_component[i]])
1289 component_numbering[i] = count++;
1296 const ::hp::FECollection<dim, spacedim> &fe_collection =
1298 std::vector<Table<2, bool>> element_constant_modes;
1299 std::vector<std::vector<std::pair<unsigned int, unsigned int>>>
1300 constant_mode_to_component_translation(n_components);
1302 unsigned int n_constant_modes = 0;
1303 int first_non_empty_constant_mode = -1;
1304 for (
unsigned int f = 0; f < fe_collection.size(); ++f)
1306 std::pair<Table<2, bool>, std::vector<unsigned int>> data =
1307 fe_collection[f].get_constant_modes();
1311 if (first_non_empty_constant_mode < 0 && data.first.n_rows() > 0)
1313 first_non_empty_constant_mode = f;
1317 for (
unsigned int i = 0; i < data.second.size(); ++i)
1318 if (component_mask[data.second[i]])
1319 constant_mode_to_component_translation[data.second[i]]
1320 .emplace_back(n_constant_modes++, i);
1326 element_constant_modes.push_back(data.first);
1328 element_constant_modes.back().n_rows() == 0 ||
1329 element_constant_modes.back().n_rows() ==
1330 element_constant_modes[first_non_empty_constant_mode].n_rows(),
1337 constant_modes.clear();
1338 constant_modes.resize(n_constant_modes,
1339 std::vector<bool>(n_selected_dofs,
false));
1343 std::vector<types::global_dof_index> dof_indices;
1347 dof_indices.resize(cell->get_fe().n_dofs_per_cell());
1348 cell->get_dof_indices(dof_indices);
1350 for (
unsigned int i = 0; i < dof_indices.size(); ++i)
1351 if (locally_owned_dofs.
is_element(dof_indices[i]))
1353 const unsigned int loc_index =
1355 const unsigned int comp = dofs_by_component[loc_index];
1356 if (component_mask[comp])
1357 for (
auto &indices :
1358 constant_mode_to_component_translation[comp])
1359 constant_modes[indices
1360 .first][component_numbering[loc_index]] =
1361 element_constant_modes[cell->active_fe_index()](
1369 template <
int dim,
int spacedim>
1372 std::vector<unsigned int> & active_fe_indices)
1378 active_fe_indices[cell->active_cell_index()] = cell->active_fe_index();
1381 template <
int dim,
int spacedim>
1382 std::vector<IndexSet>
1386 ExcMessage(
"The given DoFHandler has no DoFs."));
1394 "For parallel::distributed::Triangulation objects and "
1395 "associated DoF handler objects, asking for any information "
1396 "related to a subdomain other than the locally owned one does "
1397 "not make sense."));
1401 std::vector<::types::subdomain_id> subdomain_association(
1404 subdomain_association);
1418 const unsigned int n_subdomains =
1422 unsigned int max_subdomain_id = 0;
1425 std::max(max_subdomain_id, cell->subdomain_id());
1426 return max_subdomain_id + 1;
1431 ->get_communicator()));
1432 Assert(n_subdomains > *std::max_element(subdomain_association.begin(),
1433 subdomain_association.end()),
1436 std::vector<::IndexSet> index_sets(
1445 index < subdomain_association.size();
1449 if (subdomain_association[index] != this_subdomain)
1451 index_sets[this_subdomain].add_range(i_min, index);
1453 this_subdomain = subdomain_association[index];
1458 if (i_min == subdomain_association.size() - 1)
1460 index_sets[this_subdomain].add_index(i_min);
1466 index_sets[this_subdomain].add_range(i_min,
1467 subdomain_association.size());
1470 for (
unsigned int i = 0; i < n_subdomains; ++i)
1471 index_sets[i].compress();
1476 template <
int dim,
int spacedim>
1477 std::vector<IndexSet>
1487 "For parallel::distributed::Triangulation objects and "
1488 "associated DoF handler objects, asking for any information "
1489 "related to a subdomain other than the locally owned one does "
1490 "not make sense."));
1498 std::vector<IndexSet> dof_set =
1500 const ::types::subdomain_id n_subdomains = dof_set.size();
1506 subdomain_id < n_subdomains;
1519 std::vector<types::global_dof_index> local_dof_indices;
1520 std::set<types::global_dof_index> subdomain_halo_global_dof_indices;
1521 for (
typename std::vector<
1523 const_iterator it_cell = active_halo_layer.begin();
1524 it_cell != active_halo_layer.end();
1530 cell->subdomain_id() != subdomain_id,
1532 "The subdomain ID of the halo cell should not match that of the vector entry."));
1535 cell->get_dof_indices(local_dof_indices);
1539 subdomain_halo_global_dof_indices.insert(local_dof_index);
1542 dof_set[subdomain_id].add_indices(
1543 subdomain_halo_global_dof_indices.begin(),
1544 subdomain_halo_global_dof_indices.end());
1546 dof_set[subdomain_id].compress();
1552 template <
int dim,
int spacedim>
1556 std::vector<types::subdomain_id> &subdomain_association)
1564 "For parallel::distributed::Triangulation objects and "
1565 "associated DoF handler objects, asking for any subdomain other "
1566 "than the locally owned one does not make sense."));
1568 Assert(subdomain_association.size() == dof_handler.
n_dofs(),
1582 std::vector<types::subdomain_id> cell_owners(
1588 cell_owners = tr->get_true_subdomain_ids_of_cells();
1589 Assert(tr->get_true_subdomain_ids_of_cells().size() ==
1590 tr->n_active_cells(),
1597 cell_owners[cell->active_cell_index()] = cell->subdomain_id();
1601 std::fill_n(subdomain_association.begin(),
1605 std::vector<types::global_dof_index> local_dof_indices;
1606 local_dof_indices.reserve(
1614 cell_owners[cell->active_cell_index()];
1615 const unsigned int dofs_per_cell = cell->get_fe().n_dofs_per_cell();
1616 local_dof_indices.resize(dofs_per_cell);
1617 cell->get_dof_indices(local_dof_indices);
1623 for (
unsigned int i = 0; i < dofs_per_cell; ++i)
1624 if (subdomain_association[local_dof_indices[i]] ==
1626 subdomain_association[local_dof_indices[i]] = subdomain_id;
1627 else if (subdomain_association[local_dof_indices[i]] > subdomain_id)
1629 subdomain_association[local_dof_indices[i]] = subdomain_id;
1633 Assert(std::find(subdomain_association.begin(),
1634 subdomain_association.end(),
1636 subdomain_association.end(),
1642 template <
int dim,
int spacedim>
1648 std::vector<types::subdomain_id> subdomain_association(
1652 return std::count(subdomain_association.begin(),
1653 subdomain_association.end(),
1659 template <
int dim,
int spacedim>
1672 "For parallel::distributed::Triangulation objects and "
1673 "associated DoF handler objects, asking for any subdomain other "
1674 "than the locally owned one does not make sense."));
1678 std::vector<types::global_dof_index> local_dof_indices;
1679 local_dof_indices.reserve(
1687 std::vector<types::global_dof_index> subdomain_indices;
1690 if ((cell->is_artificial() ==
false) &&
1691 (cell->subdomain_id() == subdomain))
1693 const unsigned int dofs_per_cell = cell->get_fe().n_dofs_per_cell();
1694 local_dof_indices.resize(dofs_per_cell);
1695 cell->get_dof_indices(local_dof_indices);
1696 subdomain_indices.insert(subdomain_indices.end(),
1697 local_dof_indices.begin(),
1698 local_dof_indices.end());
1701 std::sort(subdomain_indices.begin(), subdomain_indices.end());
1702 subdomain_indices.erase(std::unique(subdomain_indices.begin(),
1703 subdomain_indices.end()),
1704 subdomain_indices.end());
1707 index_set.
add_indices(subdomain_indices.begin(), subdomain_indices.end());
1715 template <
int dim,
int spacedim>
1720 std::vector<unsigned int> & n_dofs_on_subdomain)
1725 std::fill(n_dofs_on_subdomain.begin(), n_dofs_on_subdomain.end(), 0);
1734 return cell.subdomain_id() == subdomain;
1736 ExcMessage(
"There are no cells for the given subdomain!"));
1738 std::vector<types::subdomain_id> subdomain_association(
1742 std::vector<unsigned char> component_association(dof_handler.
n_dofs());
1744 std::vector<bool>(),
1745 component_association);
1747 for (
unsigned int c = 0; c < dof_handler.
get_fe(0).n_components(); ++c)
1750 if ((subdomain_association[i] == subdomain) &&
1751 (component_association[i] ==
static_cast<unsigned char>(c)))
1752 ++n_dofs_on_subdomain[c];
1764 template <
int dim,
int spacedim>
1767 const std::vector<unsigned char> & dofs_by_component,
1768 const std::vector<unsigned int> & target_component,
1769 const bool only_once,
1770 std::vector<types::global_dof_index> &dofs_per_component,
1771 unsigned int & component)
1790 for (
unsigned int dd = 0; dd < d; ++dd, ++component)
1791 dofs_per_component[target_component[component]] +=
1792 std::count(dofs_by_component.begin(),
1793 dofs_by_component.end(),
1800 for (
unsigned int dd = 1; dd < d; ++dd)
1801 dofs_per_component[target_component[component - d + dd]] =
1802 dofs_per_component[target_component[component - d]];
1809 template <
int dim,
int spacedim>
1812 const std::vector<unsigned char> & dofs_by_component,
1813 const std::vector<unsigned int> & target_component,
1814 const bool only_once,
1815 std::vector<types::global_dof_index> &dofs_per_component,
1816 unsigned int & component)
1821 for (
unsigned int fe = 1; fe < fe_collection.
size(); ++fe)
1823 Assert(fe_collection[fe].n_components() ==
1824 fe_collection[0].n_components(),
1826 Assert(fe_collection[fe].n_base_elements() ==
1827 fe_collection[0].n_base_elements(),
1829 for (
unsigned int b = 0; b < fe_collection[0].n_base_elements(); ++b)
1831 Assert(fe_collection[fe].base_element(b).n_components() ==
1832 fe_collection[0].base_element(b).n_components(),
1834 Assert(fe_collection[fe].base_element(b).n_base_elements() ==
1835 fe_collection[0].base_element(b).n_base_elements(),
1858 template <
int dim,
int spacedim>
1870 template <
int dim,
int spacedim>
1872 all_elements_are_primitive(
1873 const ::hp::FECollection<dim, spacedim> &fe_collection)
1875 for (
unsigned int i = 0; i < fe_collection.size(); ++i)
1876 if (fe_collection[i].is_primitive() ==
false)
1886 template <
int dim,
int spacedim>
1887 std::vector<types::global_dof_index>
1890 const bool only_once,
1891 const std::vector<unsigned int> &target_component_)
1897 std::vector<unsigned int> target_component = target_component_;
1898 if (target_component.size() == 0)
1900 target_component.resize(n_components);
1901 for (
unsigned int i = 0; i < n_components; ++i)
1902 target_component[i] = i;
1905 Assert(target_component.size() == n_components,
1909 const unsigned int max_component =
1910 *std::max_element(target_component.begin(), target_component.end());
1911 const unsigned int n_target_components = max_component + 1;
1913 std::vector<types::global_dof_index> dofs_per_component(
1918 if (n_components == 1)
1921 return dofs_per_component;
1927 std::vector<unsigned char> dofs_by_component(
1934 unsigned int component = 0;
1945 Assert((internal::all_elements_are_primitive(
1947 (std::accumulate(dofs_per_component.begin(),
1948 dofs_per_component.end(),
1954#ifdef DEAL_II_WITH_MPI
1960 std::vector<types::global_dof_index> local_dof_count =
1963 const int ierr = MPI_Allreduce(local_dof_count.data(),
1964 dofs_per_component.data(),
1965 n_target_components,
1973 return dofs_per_component;
1978 template <
int dim,
int spacedim>
1979 std::vector<types::global_dof_index>
1981 const std::vector<unsigned int> &target_block_)
1983 const ::hp::FECollection<dim, spacedim> &fe_collection =
1992 const unsigned int n_blocks = fe_collection[0].n_blocks();
1994 std::vector<unsigned int> target_block = target_block_;
1995 if (target_block.size() == 0)
1997 target_block.resize(fe_collection[0].n_blocks());
1998 for (
unsigned int i = 0; i < n_blocks; ++i)
1999 target_block[i] = i;
2002 Assert(target_block.size() == n_blocks,
2004 for (
unsigned int f = 1; f < fe_collection.size(); ++f)
2005 Assert(fe_collection[0].n_blocks() == fe_collection[f].n_blocks(),
2006 ExcMessage(
"This function can only work if all elements in a "
2007 "collection have the same number of blocks."));
2013 std::vector<types::global_dof_index> dofs_per_block(1);
2014 dofs_per_block[0] = dof_handler.
n_dofs();
2015 return dofs_per_block;
2019 const unsigned int max_block =
2020 *std::max_element(target_block.begin(), target_block.end());
2021 const unsigned int n_target_blocks = max_block + 1;
2023 std::vector<types::global_dof_index> dofs_per_block(n_target_blocks);
2027 for (
unsigned int this_fe = fe_collection.size() - 1;
2028 this_fe < fe_collection.size();
2033 std::vector<unsigned char> dofs_by_block(
2038 for (
unsigned int block = 0; block < fe.
n_blocks(); ++block)
2039 dofs_per_block[target_block[block]] +=
2040 std::count(dofs_by_block.begin(), dofs_by_block.end(), block);
2042#ifdef DEAL_II_WITH_MPI
2049 std::vector<types::global_dof_index> local_dof_count =
2051 const int ierr = MPI_Allreduce(local_dof_count.data(),
2052 dofs_per_block.data(),
2062 return dofs_per_block;
2067 template <
int dim,
int spacedim>
2070 std::vector<types::global_dof_index> &mapping)
2073 mapping.insert(mapping.end(),
2077 std::vector<types::global_dof_index> dofs_on_face;
2088 for (
const unsigned int f : cell->face_indices())
2089 if (cell->at_boundary(f))
2091 const unsigned int dofs_per_face =
2092 cell->get_fe().n_dofs_per_face(f);
2093 dofs_on_face.resize(dofs_per_face);
2094 cell->face(f)->get_dof_indices(dofs_on_face,
2095 cell->active_fe_index());
2096 for (
unsigned int i = 0; i < dofs_per_face; ++i)
2098 mapping[dofs_on_face[i]] = next_boundary_index++;
2106 template <
int dim,
int spacedim>
2109 const std::set<types::boundary_id> &boundary_ids,
2110 std::vector<types::global_dof_index> &mapping)
2117 mapping.insert(mapping.end(),
2122 if (boundary_ids.size() == 0)
2125 std::vector<types::global_dof_index> dofs_on_face;
2130 for (
const unsigned int f : cell->face_indices())
2131 if (boundary_ids.find(cell->face(f)->boundary_id()) !=
2134 const unsigned int dofs_per_face =
2135 cell->get_fe().n_dofs_per_face(f);
2136 dofs_on_face.resize(dofs_per_face);
2137 cell->face(f)->get_dof_indices(dofs_on_face,
2138 cell->active_fe_index());
2139 for (
unsigned int i = 0; i < dofs_per_face; ++i)
2141 mapping[dofs_on_face[i]] = next_boundary_index++;
2152 template <
int dim,
int spacedim>
2164 for (
unsigned int fe_index = 0; fe_index < fe_collection.
size();
2168 Assert((fe_collection[fe_index].n_dofs_per_cell() == 0) ||
2169 (fe_collection[fe_index].has_support_points()),
2172 fe_collection[fe_index].get_unit_support_points()));
2177 (in_mask.
size() == 0 ?
2193 std::vector<types::global_dof_index> local_dof_indices;
2196 if (cell->is_artificial() ==
false)
2198 hp_fe_values.reinit(cell);
2202 local_dof_indices.resize(cell->get_fe().n_dofs_per_cell());
2203 cell->get_dof_indices(local_dof_indices);
2205 const std::vector<Point<spacedim>> &points =
2207 for (
unsigned int i = 0; i < cell->get_fe().n_dofs_per_cell();
2210 const unsigned int dof_comp =
2211 cell->get_fe().system_to_component_index(i).first;
2215 support_points[local_dof_indices[i]] = points[i];
2221 template <
int dim,
int spacedim>
2230 std::map<types::global_dof_index, Point<spacedim>> x_support_points;
2240 Assert(x_support_points.find(i) != x_support_points.end(),
2243 support_points[i] = x_support_points[i];
2249 template <
int dim,
int spacedim>
2261 "This function can not be used with distributed triangulations. "
2262 "See the documentation for more information."));
2268 internal::map_dofs_to_support_points(mapping_collection,
2275 template <
int dim,
int spacedim>
2288 "This function can not be used with distributed triangulations. "
2289 "See the documentation for more information."));
2293 internal::map_dofs_to_support_points(mapping,
2300 template <
int dim,
int spacedim>
2308 support_points.clear();
2314 internal::map_dofs_to_support_points(mapping_collection,
2321 template <
int dim,
int spacedim>
2329 support_points.clear();
2333 internal::map_dofs_to_support_points(mapping,
2339 template <
int spacedim>
2347 std::map<Point<spacedim>,
2348 std::vector<types::global_dof_index>,
2353 for (
const auto &it : support_points)
2355 std::vector<types::global_dof_index> &v = point_map[it.second];
2356 v.push_back(it.first);
2360 for (
const auto &it : point_map)
2362 out << it.first <<
" \"";
2363 const std::vector<types::global_dof_index> &v = it.second;
2364 for (
unsigned int i = 0; i < v.size(); ++i)
2378 template <
int dim,
int spacedim>
2387 const unsigned int nb = fe.
n_blocks();
2389 tables_by_block.resize(1);
2390 tables_by_block[0].reinit(nb, nb);
2391 tables_by_block[0].fill(
none);
2399 tables_by_block[0](ib, jb) |= table(i, j);
2407 tables_by_block.resize(fe_collection.
size());
2409 for (
unsigned int f = 0; f < fe_collection.
size(); ++f)
2413 const unsigned int nb = fe.
n_blocks();
2414 tables_by_block[f].reinit(nb, nb);
2415 tables_by_block[f].fill(
none);
2422 tables_by_block[f](ib, jb) |= table(i, j);
2431 template <
int dim,
int spacedim>
2435 const unsigned int level,
2436 const std::vector<bool> & selected_dofs,
2439 std::vector<types::global_dof_index> indices;
2444 if (cell->is_locally_owned_on_level())
2451 if (cell->is_locally_owned_on_level())
2453 indices.resize(cell->get_fe().n_dofs_per_cell());
2454 cell->get_mg_dof_indices(indices);
2456 if (selected_dofs.size() != 0)
2461 if (selected_dofs.size() == 0)
2462 block_list.
add(i, indices[j] - offset);
2465 if (selected_dofs[j])
2466 block_list.
add(i, indices[j] - offset);
2474 template <
int dim,
int spacedim>
2478 const unsigned int level,
2479 const bool interior_only)
2484 std::vector<types::global_dof_index> indices;
2485 std::vector<bool> exclude;
2489 indices.resize(cell->get_fe().n_dofs_per_cell());
2490 cell->get_mg_dof_indices(indices);
2496 std::fill(exclude.begin(), exclude.end(),
false);
2498 for (
const unsigned int face : cell->face_indices())
2499 if (cell->at_boundary(face) ||
2500 cell->neighbor(face)->level() != cell->level())
2505 block_list.
add(0, indices[j]);
2509 for (
const auto index : indices)
2510 block_list.
add(0, index);
2516 template <
int dim,
int spacedim>
2520 const unsigned int level,
2521 const bool interior_dofs_only,
2522 const bool boundary_dofs)
2527 std::vector<types::global_dof_index> indices;
2528 std::vector<bool> exclude;
2530 unsigned int block = 0;
2533 if (pcell->is_active())
2536 for (
unsigned int child = 0; child < pcell->n_children(); ++child)
2538 const auto cell = pcell->child(child);
2543 indices.resize(n_dofs);
2544 exclude.resize(n_dofs);
2545 std::fill(exclude.begin(), exclude.end(),
false);
2546 cell->get_mg_dof_indices(indices);
2548 if (interior_dofs_only)
2552 for (
unsigned int d = 0; d < dim; ++d)
2554 const unsigned int face =
2563 for (
const unsigned int face :
2565 if (cell->at_boundary(face))
2571 for (
unsigned int i = 0; i < n_dofs; ++i)
2573 block_list.
add(block, indices[i]);
2579 template <
int dim,
int spacedim>
2580 std::vector<unsigned int>
2583 const unsigned int level,
2584 const bool interior_only,
2585 const bool boundary_patches,
2586 const bool level_boundary_patches,
2587 const bool single_cell_patches,
2588 const bool invert_vertex_mapping)
2595 exclude_boundary_dofs,
2597 level_boundary_patches,
2598 single_cell_patches,
2599 invert_vertex_mapping);
2602 template <
int dim,
int spacedim>
2603 std::vector<unsigned int>
2606 const unsigned int level,
2607 const BlockMask & exclude_boundary_dofs,
2608 const bool boundary_patches,
2609 const bool level_boundary_patches,
2610 const bool single_cell_patches,
2611 const bool invert_vertex_mapping)
2615 std::vector<unsigned int> vertex_cell_count(
2619 std::vector<bool> vertex_boundary(
2622 std::vector<unsigned int> vertex_mapping(
2627 std::vector<unsigned int> vertex_dof_count(
2633 for (
const unsigned int v : cell->vertex_indices())
2635 const unsigned int vg = cell->vertex_index(v);
2636 vertex_dof_count[vg] += cell->get_fe().n_dofs_per_cell();
2637 ++vertex_cell_count[vg];
2638 for (
unsigned int d = 0; d < dim; ++d)
2641 if (cell->at_boundary(face))
2642 vertex_boundary[vg] =
true;
2643 else if ((!level_boundary_patches) &&
2644 (cell->neighbor(face)->level() !=
2645 static_cast<int>(
level)))
2646 vertex_boundary[vg] =
true;
2652 for (
unsigned int vg = 0; vg < vertex_dof_count.size(); ++vg)
2653 if ((!single_cell_patches && vertex_cell_count[vg] < 2) ||
2654 (!boundary_patches && vertex_boundary[vg]))
2655 vertex_dof_count[vg] = 0;
2658 unsigned int n_vertex_count = 0;
2659 for (
unsigned int vg = 0; vg < vertex_mapping.size(); ++vg)
2660 if (vertex_dof_count[vg] != 0)
2661 vertex_mapping[vg] = n_vertex_count++;
2664 for (
unsigned int vg = 0; vg < vertex_mapping.size(); ++vg)
2665 if (vertex_dof_count[vg] != 0)
2666 vertex_dof_count[vertex_mapping[vg]] = vertex_dof_count[vg];
2670 vertex_dof_count.resize(n_vertex_count);
2674 block_list.
reinit(vertex_dof_count.size(),
2678 std::vector<types::global_dof_index> indices;
2679 std::vector<bool> exclude;
2685 cell->get_mg_dof_indices(indices);
2687 for (
const unsigned int v : cell->vertex_indices())
2689 const unsigned int vg = cell->vertex_index(v);
2690 const unsigned int block = vertex_mapping[vg];
2696 if (exclude_boundary_dofs.
size() == 0 ||
2702 std::fill(exclude.begin(), exclude.end(),
false);
2704 for (
unsigned int d = 0; d < dim; ++d)
2706 const unsigned int a_face =
2708 const unsigned int face =
2721 for (
unsigned int j = 0; j < indices.size(); ++j)
2723 block_list.
add(block, indices[j]);
2727 for (
const auto index : indices)
2728 block_list.
add(block, index);
2733 if (invert_vertex_mapping)
2736 unsigned int n_vertex_count = 0;
2737 for (
unsigned int vg = 0; vg < vertex_mapping.size(); ++vg)
2739 vertex_mapping[n_vertex_count++] = vg;
2742 vertex_mapping.resize(n_vertex_count);
2745 return vertex_mapping;
2749 template <
int dim,
int spacedim>
2755 std::set<types::global_dof_index> dofs_on_patch;
2756 std::vector<types::global_dof_index> local_dof_indices;
2761 for (
unsigned int i = 0; i < patch.size(); ++i)
2765 Assert(cell->is_artificial() ==
false,
2766 ExcMessage(
"This function can not be called with cells that are "
2767 "not either locally owned or ghost cells."));
2769 cell->get_dof_indices(local_dof_indices);
2770 dofs_on_patch.insert(local_dof_indices.begin(),
2771 local_dof_indices.end());
2775 return dofs_on_patch.size();
2780 template <
int dim,
int spacedim>
2781 std::vector<types::global_dof_index>
2786 std::set<types::global_dof_index> dofs_on_patch;
2787 std::vector<types::global_dof_index> local_dof_indices;
2792 for (
unsigned int i = 0; i < patch.size(); ++i)
2796 Assert(cell->is_artificial() ==
false,
2797 ExcMessage(
"This function can not be called with cells that are "
2798 "not either locally owned or ghost cells."));
2800 cell->get_dof_indices(local_dof_indices);
2801 dofs_on_patch.insert(local_dof_indices.begin(),
2802 local_dof_indices.end());
2805 Assert((dofs_on_patch.size() == count_dofs_on_patch<dim, spacedim>(patch)),
2812 return std::vector<types::global_dof_index>(dofs_on_patch.begin(),
2813 dofs_on_patch.end());
2823#include "dof_tools.inst"
const std::vector< std::pair< size_type, number > > * get_constraint_entries(const size_type line_n) const
size_type n_constraints() const
unsigned int size() const
unsigned int n_selected_blocks(const unsigned int overall_number_of_blocks=numbers::invalid_unsigned_int) const
bool represents_n_components(const unsigned int n) const
bool represents_the_all_selected_mask() const
unsigned int size() const
unsigned int n_selected_components(const unsigned int overall_number_of_components=numbers::invalid_unsigned_int) const
unsigned int first_selected_component(const unsigned int overall_number_of_components=numbers::invalid_unsigned_int) const
const FiniteElement< dim, spacedim > & get_fe(const unsigned int index=0) const
const hp::FECollection< dim, spacedim > & get_fe_collection() const
const IndexSet & locally_owned_mg_dofs(const unsigned int level) const
types::global_dof_index n_boundary_dofs() const
const Triangulation< dim, spacedim > & get_triangulation() const
const IndexSet & locally_owned_dofs() const
active_cell_iterator begin_active(const unsigned int level=0) const
bool has_hp_capabilities() const
types::global_dof_index n_dofs() const
cell_iterator begin(const unsigned int level=0) const
types::global_dof_index n_locally_owned_dofs() const
const std::vector< Point< spacedim > > & get_quadrature_points() const
const FEValues< dim, spacedim > & get_present_fe_values() const
unsigned int n_dofs_per_vertex() const
unsigned int n_dofs_per_cell() const
unsigned int n_dofs_per_line() const
unsigned int n_dofs_per_face(unsigned int face_no=0, unsigned int child=0) const
unsigned int n_blocks() const
unsigned int n_components() 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
virtual const FiniteElement< dim, spacedim > & base_element(const unsigned int index) const
ComponentMask component_mask(const FEValuesExtractors::Scalar &scalar) const
bool is_primitive() const
virtual bool has_support_on_face(const unsigned int shape_index, const unsigned int face_index) const
virtual unsigned int face_to_cell_index(const unsigned int face_dof_index, const unsigned int face, const bool face_orientation=true, const bool face_flip=false, const bool face_rotation=false) const
unsigned int component_to_block_index(const unsigned int component) const
std::pair< unsigned int, unsigned int > system_to_component_index(const unsigned int index) const
unsigned int element_multiplicity(const unsigned int index) const
unsigned int n_base_elements() const
std::pair< unsigned int, unsigned int > face_system_to_component_index(const unsigned int index, const unsigned int face_no=0) 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
void add_index(const size_type index)
void fill_binary_vector(VectorType &vector) 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)
Abstract base class for mapping classes.
void add(const size_type i, const size_type j)
virtual void reinit(const size_type m, const size_type n, const ArrayView< const unsigned int > &row_lengths) override
virtual MPI_Comm get_communicator() const
virtual types::subdomain_id locally_owned_subdomain() const
unsigned int n_active_cells() const
unsigned int n_levels() const
unsigned int n_vertices() const
unsigned int size() const
unsigned int max_dofs_per_face() const
ComponentMask component_mask(const FEValuesExtractors::Scalar &scalar) const
unsigned int n_components() const
unsigned int max_dofs_per_cell() const
void push_back(const Quadrature< dim_in > &new_quadrature)
#define DEAL_II_NAMESPACE_OPEN
#define DEAL_II_NAMESPACE_CLOSE
@ update_quadrature_points
Transformed quadrature points.
IteratorRange< active_cell_iterator > active_cell_iterators() const
IteratorRange< FilteredIterator< BaseIterator > > filter_iterators(IteratorRange< BaseIterator > i, const Predicate &p)
IteratorRange< cell_iterator > cell_iterators_on_level(const unsigned int level) const
static ::ExceptionBase & ExcIO()
static ::ExceptionBase & ExcInvalidBoundaryIndicator()
static ::ExceptionBase & ExcNotImplemented()
#define Assert(cond, exc)
#define AssertDimension(dim1, dim2)
#define AssertThrowMPI(error_code)
#define AssertIndexRange(index, range)
static ::ExceptionBase & ExcInternalError()
static ::ExceptionBase & ExcIndexRange(std::size_t arg1, std::size_t arg2, std::size_t arg3)
static ::ExceptionBase & ExcDimensionMismatch(std::size_t arg1, std::size_t arg2)
static ::ExceptionBase & ExcMessage(std::string arg1)
#define AssertThrow(cond, exc)
typename ActiveSelector::active_cell_iterator active_cell_iterator
typename ActiveSelector::CellAccessor cell_accessor
unsigned int n_mpi_processes(const MPI_Comm &mpi_communicator)
const types::boundary_id internal_face_boundary_id
const types::subdomain_id artificial_subdomain_id
const types::subdomain_id invalid_subdomain_id
static const unsigned int invalid_unsigned_int
const types::global_dof_index invalid_dof_index
::VectorizedArray< Number, width > max(const ::VectorizedArray< Number, width > &, const ::VectorizedArray< Number, width > &)
const ::Triangulation< dim, spacedim > & tria
#define DEAL_II_DOF_INDEX_MPI_TYPE