57 template <
int dim,
int spacedim>
60 const ::DoFHandler<dim, spacedim> &dof_handler,
62 std::vector<std::vector<
63 std::pair<types::global_dof_index, types::global_dof_index>>>
65 std::vector<std::vector<
66 std::pair<types::global_dof_index, types::global_dof_index>>>
67 ©_indices_global_mine,
68 std::vector<std::vector<
69 std::pair<types::global_dof_index, types::global_dof_index>>>
70 & copy_indices_level_mine,
71 const bool skip_interface_dofs)
84 std::vector<DoFPair> send_data_temp;
86 const unsigned int n_levels =
87 dof_handler.get_triangulation().n_global_levels();
88 copy_indices.resize(n_levels);
89 copy_indices_global_mine.resize(n_levels);
90 copy_indices_level_mine.resize(n_levels);
94 const unsigned int dofs_per_cell = dof_handler.get_fe().dofs_per_cell;
95 std::vector<types::global_dof_index> global_dof_indices(dofs_per_cell);
96 std::vector<types::global_dof_index> level_dof_indices(dofs_per_cell);
100 std::vector<bool> dof_touched(globally_relevant.
n_elements(),
false);
106 std::vector<types::global_dof_index> unrolled_copy_indices;
108 copy_indices_level_mine[
level].clear();
109 copy_indices_global_mine[
level].clear();
111 for (
const auto &level_cell :
112 dof_handler.active_cell_iterators_on_level(
level))
114 if (dof_handler.get_triangulation().locally_owned_subdomain() !=
116 (level_cell->level_subdomain_id() ==
118 level_cell->subdomain_id() ==
122 unrolled_copy_indices.resize(
123 dof_handler.locally_owned_dofs().n_elements(),
128 level_cell->get_dof_indices(global_dof_indices);
129 level_cell->get_mg_dof_indices(level_dof_indices);
131 for (
unsigned int i = 0; i < dofs_per_cell; ++i)
135 if (skip_interface_dofs && mg_constrained_dofs !=
nullptr &&
137 level, level_dof_indices[i]))
145 dof_handler.locally_owned_dofs().is_element(
146 global_dof_indices[i]);
148 dof_handler.locally_owned_mg_dofs(
level).is_element(
149 level_dof_indices[i]);
151 if (global_mine && level_mine)
156 unrolled_copy_indices[dof_handler.locally_owned_dofs()
158 global_dof_indices[i])] =
159 level_dof_indices[i];
168 global_dof_indices[i]);
172 if (dof_touched[relevant_idx] ==
false)
176 copy_indices_global_mine[
level].emplace_back(
177 global_dof_indices[i], level_dof_indices[i]);
180 send_data_temp.emplace_back(
level,
181 global_dof_indices[i],
182 level_dof_indices[i]);
189 dof_touched[relevant_idx] =
true;
197 if (!unrolled_copy_indices.empty())
199 copy_indices[
level].clear();
203 if (copy_indices_global_mine[
level].empty())
204 copy_indices[
level].reserve(unrolled_copy_indices.size());
209 for (
unsigned int i = 0; i < unrolled_copy_indices.size(); ++i)
211 copy_indices[
level].emplace_back(
212 dof_handler.locally_owned_dofs().nth_index_in_set(i),
213 unrolled_copy_indices[i]);
217 const ::parallel::TriangulationBase<dim, spacedim> *tria =
218 (
dynamic_cast<const ::parallel::TriangulationBase<dim, spacedim>
219 *
>(&dof_handler.get_triangulation()));
221 send_data_temp.size() == 0 || tria !=
nullptr,
223 "We should only be sending information with a parallel Triangulation!"));
225 #ifdef DEAL_II_WITH_MPI
227 tria->get_communicator()) > 0)
229 const std::set<types::subdomain_id> &neighbors =
230 tria->level_ghost_owners();
231 std::map<int, std::vector<DoFPair>> send_data;
233 std::sort(send_data_temp.begin(),
234 send_data_temp.end(),
236 if (lhs.level < rhs.level)
238 if (lhs.level > rhs.level)
241 if (lhs.level_dof_index < rhs.level_dof_index)
243 if (lhs.level_dof_index > rhs.level_dof_index)
246 if (lhs.global_dof_index < rhs.global_dof_index)
251 send_data_temp.erase(
252 std::unique(send_data_temp.begin(),
253 send_data_temp.end(),
255 return (lhs.level == rhs.level) &&
256 (lhs.level_dof_index == rhs.level_dof_index) &&
257 (lhs.global_dof_index == rhs.global_dof_index);
259 send_data_temp.end());
264 dof_handler.locally_owned_mg_dofs(
level);
266 std::vector<types::global_dof_index> level_dof_indices;
267 std::vector<types::global_dof_index> global_dof_indices;
268 for (
const auto &dofpair : send_data_temp)
269 if (dofpair.level ==
level)
271 level_dof_indices.push_back(dofpair.level_dof_index);
272 global_dof_indices.push_back(dofpair.global_dof_index);
277 level_dof_indices.end());
282 const auto index_owner =
285 tria->get_communicator());
287 AssertThrow(level_dof_indices.size() == index_owner.size(),
290 for (
unsigned int i = 0; i < index_owner.size(); i++)
291 send_data[index_owner[i]].emplace_back(
level,
292 global_dof_indices[i],
293 level_dof_indices[i]);
300 mutex, tria->get_communicator());
306 std::vector<MPI_Request> requests;
308 for (
const auto dest : neighbors)
311 std::vector<DoFPair> &data = send_data[dest];
317 const int ierr = MPI_Isend(data.data(),
318 data.size() *
sizeof(data[0]),
322 tria->get_communicator(),
323 &*requests.rbegin());
328 const int ierr = MPI_Isend(
nullptr,
333 tria->get_communicator(),
334 &*requests.rbegin());
343 std::vector<DoFPair> receive_buffer;
344 for (
unsigned int counter = 0; counter < neighbors.size();
348 int ierr = MPI_Probe(MPI_ANY_SOURCE,
350 tria->get_communicator(),
354 ierr = MPI_Get_count(&status, MPI_BYTE, &len);
359 ierr = MPI_Recv(
nullptr,
364 tria->get_communicator(),
370 int count = len /
sizeof(
DoFPair);
373 receive_buffer.resize(count);
375 void *ptr = receive_buffer.data();
381 tria->get_communicator(),
385 for (
const auto &dof_pair : receive_buffer)
387 copy_indices_level_mine[dof_pair.level].emplace_back(
388 dof_pair.global_dof_index, dof_pair.level_dof_index);
394 if (requests.size() > 0)
396 const int ierr = MPI_Waitall(requests.size(),
398 MPI_STATUSES_IGNORE);
406 const int ierr = MPI_Barrier(tria->get_communicator());
416 std::less<std::pair<types::global_dof_index, types::global_dof_index>>
418 for (
auto &level_indices : copy_indices_level_mine)
419 std::sort(level_indices.begin(), level_indices.end(), compare);
420 for (
auto &level_indices : copy_indices_global_mine)
421 std::sort(level_indices.begin(), level_indices.end(), compare);
431 std::vector<types::global_dof_index> &ghosted_level_dofs,
432 const std::shared_ptr<const Utilities::MPI::Partitioner>
433 & external_partitioner,
435 std::shared_ptr<const Utilities::MPI::Partitioner> &target_partitioner,
438 std::sort(ghosted_level_dofs.begin(), ghosted_level_dofs.end());
440 ghosted_dofs.
add_indices(ghosted_level_dofs.begin(),
441 std::unique(ghosted_level_dofs.begin(),
442 ghosted_level_dofs.end()));
446 if (target_partitioner.get() !=
nullptr &&
447 target_partitioner->
size() == locally_owned.
size())
454 const int ghosts_locally_contained =
455 (external_partitioner.get() !=
nullptr &&
460 if (external_partitioner.get() !=
nullptr &&
466 if (target_partitioner.get() !=
nullptr &&
467 target_partitioner->
size() == locally_owned.
size())
468 for (
unsigned int i = 0; i < copy_indices_global_mine.n_cols(); ++i)
469 copy_indices_global_mine(1, i) =
472 copy_indices_global_mine(1, i)));
473 target_partitioner = external_partitioner;
477 if (target_partitioner.get() !=
nullptr &&
478 target_partitioner->
size() == locally_owned.
size())
479 for (
unsigned int i = 0; i < copy_indices_global_mine.n_cols(); ++i)
480 copy_indices_global_mine(1, i) =
484 copy_indices_global_mine(1, i)));
486 locally_owned, ghosted_dofs, communicator));
496 const std::vector<types::global_dof_index> &mine,
497 const std::vector<types::global_dof_index> &remote,
498 std::vector<unsigned int> & localized_indices)
500 localized_indices.resize(mine.size() + remote.size(),
502 for (
unsigned int i = 0; i < mine.size(); ++i)
506 for (
unsigned int i = 0; i < remote.size(); ++i)
518 const unsigned int fe_shift_1d,
519 const unsigned int fe_degree)
523 unsigned int c_tensor_index[dim];
524 unsigned int tmp = child;
525 for (
unsigned int d = 0;
d < dim; ++
d)
527 c_tensor_index[
d] = tmp % 2;
530 const unsigned int n_child_dofs_1d = fe_degree + 1 + fe_shift_1d;
531 unsigned int factor = 1;
532 unsigned int shift = fe_shift_1d * c_tensor_index[0];
533 for (
unsigned int d = 1;
d < dim; ++
d)
535 factor *= n_child_dofs_1d;
536 shift =
shift + factor * fe_shift_1d * c_tensor_index[
d];
548 const unsigned int child,
549 const unsigned int fe_shift_1d,
550 const unsigned int fe_degree,
551 const std::vector<unsigned int> & lexicographic_numbering,
552 const std::vector<types::global_dof_index> &local_dof_indices,
555 const unsigned int n_child_dofs_1d = fe_degree + 1 + fe_shift_1d;
556 const unsigned int shift =
557 compute_shift_within_children<dim>(child, fe_shift_1d, fe_degree);
559 local_dof_indices.size() / Utilities::fixed_power<dim>(fe_degree + 1);
561 const unsigned int n_scalar_cell_dofs =
562 Utilities::fixed_power<dim>(n_child_dofs_1d);
564 for (
unsigned int k = 0; k < (dim > 2 ? (fe_degree + 1) : 1); ++k)
565 for (
unsigned int j = 0; j < (dim > 1 ? (fe_degree + 1) : 1); ++j)
566 for (
unsigned int i = 0; i < (fe_degree + 1); ++i, ++m)
568 const unsigned int index =
569 c * n_scalar_cell_dofs +
570 k * n_child_dofs_1d * n_child_dofs_1d + j * n_child_dofs_1d +
574 local_dof_indices[lexicographic_numbering[m]],
576 indices[index] = local_dof_indices[lexicographic_numbering[m]];
582 template <
int dim,
typename Number>
586 const ::DoFHandler<dim> &dof_handler)
589 elem_info.
n_components = dof_handler.get_fe().element_multiplicity(0);
592 dof_handler.get_fe().dofs_per_cell);
618 const unsigned int n_child_dofs_1d =
623 elem_info.
n_components * Utilities::fixed_power<dim>(n_child_dofs_1d);
627 shape_info.
reinit(dummy_quadrature, dof_handler.get_fe(), 0);
634 for (
unsigned int c = 0; c < GeometryInfo<1>::max_children_per_cell; ++c)
638 .prolongation_matrix_1d[i * n_child_dofs_1d + j + c *
shift] =
653 const unsigned int level,
654 std::vector<types::global_dof_index> &dof_indices)
656 if (mg_constrained_dofs !=
nullptr &&
659 for (
auto &ind : dof_indices)
679 template <
int dim,
typename Number>
682 const ::DoFHandler<dim> &dof_handler,
684 const std::vector<std::shared_ptr<const Utilities::MPI::Partitioner>>
685 & external_partitioners,
687 std::vector<std::vector<unsigned int>> &level_dof_indices,
688 std::vector<std::vector<std::pair<unsigned int, unsigned int>>>
689 & parent_child_connect,
690 std::vector<unsigned int> &n_owned_level_cells,
691 std::vector<std::vector<std::vector<unsigned short>>> &dirichlet_indices,
692 std::vector<std::vector<Number>> & weights_on_refined,
694 MGLevelObject<std::shared_ptr<const Utilities::MPI::Partitioner>>
695 &target_partitioners)
697 level_dof_indices.clear();
698 parent_child_connect.clear();
699 n_owned_level_cells.clear();
700 dirichlet_indices.clear();
701 weights_on_refined.clear();
707 const ::Triangulation<dim> &tria = dof_handler.get_triangulation();
713 std::string fe_name = dof_handler.get_fe().base_element(0).get_name();
715 const std::size_t template_starts = fe_name.find_first_of(
'<');
716 Assert(fe_name[template_starts + 1] ==
717 (dim == 1 ?
'1' : (dim == 2 ?
'2' :
'3')),
719 fe_name[template_starts + 1] =
'1';
721 const std::unique_ptr<FiniteElement<1>> fe(
722 FETools::get_fe_by_name<1, 1>(fe_name));
729 const unsigned int n_levels = tria.n_global_levels();
730 level_dof_indices.resize(n_levels);
731 parent_child_connect.resize(n_levels - 1);
732 n_owned_level_cells.resize(n_levels - 1);
733 std::vector<std::vector<unsigned int>> coarse_level_indices(n_levels - 1);
734 for (
unsigned int level = 0;
737 coarse_level_indices[
level].resize(tria.n_raw_cells(
level),
739 std::vector<types::global_dof_index> local_dof_indices(
740 dof_handler.get_fe().dofs_per_cell);
741 dirichlet_indices.resize(n_levels - 1);
744 Assert(external_partitioners.empty() ||
745 external_partitioners.size() == n_levels,
750 unsigned int counter = 0;
751 std::vector<types::global_dof_index> global_level_dof_indices;
752 std::vector<types::global_dof_index> global_level_dof_indices_remote;
753 std::vector<types::global_dof_index> ghosted_level_dofs;
754 std::vector<types::global_dof_index> global_level_dof_indices_l0;
755 std::vector<types::global_dof_index> ghosted_level_dofs_l0;
758 typename ::DoFHandler<dim>::cell_iterator cell,
759 endc = dof_handler.end(
level - 1);
760 for (cell = dof_handler.begin(
level - 1); cell != endc; ++cell)
764 if (!cell->has_children())
768 (tria.locally_owned_subdomain() ==
770 cell->level_subdomain_id() == tria.locally_owned_subdomain());
775 const bool cell_is_remote = !consider_cell;
776 for (
unsigned int c = 0;
777 c < GeometryInfo<dim>::max_children_per_cell;
779 if (cell->child(c)->level_subdomain_id() ==
780 tria.locally_owned_subdomain())
782 consider_cell =
true;
797 std::vector<types::global_dof_index> &next_indices =
798 cell_is_remote ? global_level_dof_indices_remote :
799 global_level_dof_indices;
800 const std::size_t start_index = next_indices.size();
803 for (
unsigned int c = 0;
804 c < GeometryInfo<dim>::max_children_per_cell;
807 if (cell_is_remote && cell->child(c)->level_subdomain_id() !=
808 tria.locally_owned_subdomain())
810 cell->child(c)->get_mg_dof_indices(local_dof_indices);
812 replace(mg_constrained_dofs,
level, local_dof_indices);
815 dof_handler.locally_owned_mg_dofs(
level);
816 for (
const auto local_dof_index : local_dof_indices)
817 if (!owned_level_dofs.
is_element(local_dof_index))
818 ghosted_level_dofs.push_back(local_dof_index);
820 add_child_indices<dim>(c,
826 &next_indices[start_index]);
829 if (cell->child(c)->has_children() &&
830 (tria.locally_owned_subdomain() ==
832 cell->child(c)->level_subdomain_id() ==
833 tria.locally_owned_subdomain()))
835 const unsigned int child_index =
836 coarse_level_indices[
level][cell->child(c)->index()];
838 parent_child_connect[
level].size());
839 unsigned int parent_index = counter;
850 parent_child_connect[
level][child_index] =
851 std::make_pair(parent_index, c);
853 static_cast<unsigned short>(-1));
857 if (mg_constrained_dofs !=
nullptr)
858 for (
unsigned int i = 0;
859 i < dof_handler.get_fe().dofs_per_cell;
865 dirichlet_indices[
level][child_index].push_back(i);
871 coarse_level_indices[
level - 1].size());
872 coarse_level_indices[
level - 1][cell->index()] = counter++;
879 if (
level == 1 && !cell_is_remote)
881 cell->get_mg_dof_indices(local_dof_indices);
883 replace(mg_constrained_dofs,
level - 1, local_dof_indices);
885 const IndexSet &owned_level_dofs_l0 =
886 dof_handler.locally_owned_mg_dofs(0);
887 for (
const auto local_dof_index : local_dof_indices)
888 if (!owned_level_dofs_l0.
is_element(local_dof_index))
889 ghosted_level_dofs_l0.push_back(local_dof_index);
891 const std::size_t start_index =
892 global_level_dof_indices_l0.size();
893 global_level_dof_indices_l0.resize(
896 add_child_indices<dim>(
898 fe->dofs_per_cell - fe->dofs_per_vertex,
902 &global_level_dof_indices_l0[start_index]);
904 dirichlet_indices[0].emplace_back();
905 if (mg_constrained_dofs !=
nullptr)
906 for (
unsigned int i = 0;
907 i < dof_handler.get_fe().dofs_per_cell;
911 local_dof_indices[elem_info
913 dirichlet_indices[0].back().push_back(i);
921 global_level_dof_indices.size());
922 n_owned_level_cells[
level - 1] = counter;
923 dirichlet_indices[
level - 1].resize(counter);
924 parent_child_connect[
level - 1].resize(
932 if (
level < n_levels - 1)
933 for (std::vector<std::pair<unsigned int, unsigned int>>::iterator
935 i != parent_child_connect[
level].end();
937 if (i->first >= tria.n_cells(
level))
939 i->first -= tria.n_cells(
level);
949 const ::parallel::TriangulationBase<dim, dim> *ptria =
951 const ::parallel::TriangulationBase<dim, dim> *
>(&tria));
953 ptria !=
nullptr ? ptria->get_communicator() : MPI_COMM_SELF;
957 external_partitioners.empty() ?
959 external_partitioners[
level],
961 target_partitioners[
level],
962 copy_indices_global_mine[
level]);
965 global_level_dof_indices,
966 global_level_dof_indices_remote,
967 level_dof_indices[
level]);
971 for (
unsigned int i = 0; i < parent_child_connect[0].size(); ++i)
972 parent_child_connect[0][i] = std::make_pair(i, 0
U);
975 ghosted_level_dofs_l0,
976 external_partitioners.empty() ?
978 external_partitioners[0],
980 target_partitioners[0],
981 copy_indices_global_mine[0]);
984 *target_partitioners[0],
985 global_level_dof_indices_l0,
986 std::vector<types::global_dof_index>(),
987 level_dof_indices[0]);
993 const unsigned int n_child_dofs_1d =
994 fe->degree + 1 + fe->dofs_per_cell - fe->dofs_per_vertex;
998 weights_on_refined.resize(n_levels - 1);
1002 target_partitioners[
level]);
1003 for (
unsigned int c = 0; c < n_owned_level_cells[
level - 1]; ++c)
1011 std::vector<unsigned int> degree_to_3(n_child_dofs_1d);
1013 for (
unsigned int i = 1; i < n_child_dofs_1d - 1; ++i)
1015 degree_to_3.back() = 2;
1019 weights_on_refined[
level - 1].resize(n_owned_level_cells[
level - 1] *
1020 Utilities::fixed_power<dim>(3));
1021 for (
unsigned int c = 0; c < n_owned_level_cells[
level - 1]; ++c)
1022 for (
unsigned int k = 0, m = 0; k < (dim > 2 ? n_child_dofs_1d : 1);
1024 for (
unsigned int j = 0; j < (dim > 1 ? n_child_dofs_1d : 1); ++j)
1026 unsigned int shift = 9 * degree_to_3[k] + 3 * degree_to_3[j];
1027 for (
unsigned int i = 0; i < n_child_dofs_1d; ++i, ++m)
1028 weights_on_refined[
level -
1029 1][c * Utilities::fixed_power<dim>(3) +
1030 shift + degree_to_3[i]] =
1033 level_dof_indices[
level]
1044 #include "mg_transfer_internal.inst"