100 const std::vector<std::shared_ptr<const Utilities::MPI::Partitioner>>
101 &external_partitioners)
105 "The underlying DoFHandler object has not had its "
106 "distribute_mg_dofs() function called, but this is a prerequisite "
107 "for multigrid transfers. You will need to call this function, "
108 "probably close to where you already call distribute_dofs()."));
109 if (external_partitioners.size() > 0)
111 Assert(this->initialize_dof_vector ==
nullptr,
112 ExcMessage(
"An initialize_dof_vector function has already "
113 "been registered in the constructor!"));
115 this->initialize_dof_vector =
116 [external_partitioners](
117 const unsigned int level,
119 vec.reinit(external_partitioners[
level]);
123 this->fill_and_communicate_copy_indices(dof_handler);
125 vector_partitioners.resize(0,
128 for (
unsigned int level = 0;
level <= this->ghosted_level_vector.max_level();
130 vector_partitioners[
level] =
131 this->ghosted_level_vector[
level].get_partitioner();
133 std::vector<std::vector<Number>> weights_unvectorized;
137 internal::MGTransfer::setup_transfer<dim, Number>(
139 this->mg_constrained_dofs,
140 external_partitioners,
143 parent_child_connect,
146 weights_unvectorized,
147 this->copy_indices_global_mine,
148 vector_partitioners);
152 this->ghosted_level_vector.resize(0, vector_partitioners.max_level());
153 for (
unsigned int level = 0;
level <= vector_partitioners.max_level();
155 if (external_partitioners.size() == vector_partitioners.max_level() + 1 &&
156 external_partitioners[
level].get() == vector_partitioners[
level].get())
157 this->ghosted_level_vector[
level].reinit(0);
159 this->ghosted_level_vector[
level].reinit(vector_partitioners[
level]);
174 const unsigned int n_levels =
177 const unsigned int n_weights_per_cell = Utilities::fixed_power<dim>(3);
178 weights_on_refined.resize(n_levels - 1);
181 weights_on_refined[
level - 1].resize(
182 ((n_owned_level_cells[
level - 1] + vec_size - 1) / vec_size) *
185 for (
unsigned int c = 0; c < n_owned_level_cells[
level - 1]; ++c)
187 const unsigned int comp = c / vec_size;
188 const unsigned int v = c % vec_size;
189 for (
unsigned int i = 0; i < n_weights_per_cell; ++i)
191 weights_on_refined[
level - 1][comp * n_weights_per_cell + i][v] =
192 weights_unvectorized[
level - 1][c * n_weights_per_cell + i];
197 evaluation_data.resize(n_child_cell_dofs);
251 const unsigned int to_level,
255 Assert((to_level >= 1) && (to_level <= level_dof_indices.size()),
259 this->vector_partitioners[to_level - 1].get();
260 if (src_inplace ==
false)
262 if (this->ghosted_level_vector[to_level - 1].get_partitioner().get() !=
263 this->vector_partitioners[to_level - 1].get())
264 this->ghosted_level_vector[to_level - 1].reinit(
265 this->vector_partitioners[to_level - 1]);
266 this->ghosted_level_vector[to_level - 1].copy_locally_owned_data_from(
270 const bool dst_inplace =
271 dst.
get_partitioner().get() == this->vector_partitioners[to_level].get();
272 if (dst_inplace ==
false)
274 if (this->ghosted_level_vector[to_level].get_partitioner().get() !=
275 this->vector_partitioners[to_level].get())
276 this->ghosted_level_vector[to_level].reinit(
277 this->vector_partitioners[to_level]);
280 this->ghosted_level_vector[to_level] = 0.;
284 src_inplace ? src : this->ghosted_level_vector[to_level - 1];
286 dst_inplace ? dst : this->ghosted_level_vector[to_level];
293 do_prolongate_add<0>(to_level, dst_vec, src_vec);
294 else if (fe_degree == 1)
295 do_prolongate_add<1>(to_level, dst_vec, src_vec);
296 else if (fe_degree == 2)
297 do_prolongate_add<2>(to_level, dst_vec, src_vec);
298 else if (fe_degree == 3)
299 do_prolongate_add<3>(to_level, dst_vec, src_vec);
300 else if (fe_degree == 4)
301 do_prolongate_add<4>(to_level, dst_vec, src_vec);
302 else if (fe_degree == 5)
303 do_prolongate_add<5>(to_level, dst_vec, src_vec);
304 else if (fe_degree == 6)
305 do_prolongate_add<6>(to_level, dst_vec, src_vec);
306 else if (fe_degree == 7)
307 do_prolongate_add<7>(to_level, dst_vec, src_vec);
308 else if (fe_degree == 8)
309 do_prolongate_add<8>(to_level, dst_vec, src_vec);
310 else if (fe_degree == 9)
311 do_prolongate_add<9>(to_level, dst_vec, src_vec);
312 else if (fe_degree == 10)
313 do_prolongate_add<10>(to_level, dst_vec, src_vec);
315 do_prolongate_add<-1>(to_level, dst_vec, src_vec);
318 if (dst_inplace ==
false)
321 if (src_inplace ==
true)
330 const unsigned int from_level,
334 Assert((from_level >= 1) && (from_level <= level_dof_indices.size()),
337 const bool src_inplace =
338 src.
get_partitioner().get() == this->vector_partitioners[from_level].get();
339 if (src_inplace ==
false)
341 if (this->ghosted_level_vector[from_level].get_partitioner().get() !=
342 this->vector_partitioners[from_level].get())
343 this->ghosted_level_vector[from_level].reinit(
344 this->vector_partitioners[from_level]);
345 this->ghosted_level_vector[from_level].copy_locally_owned_data_from(src);
349 this->vector_partitioners[from_level - 1].get();
350 if (dst_inplace ==
false)
352 if (this->ghosted_level_vector[from_level - 1].get_partitioner().get() !=
353 this->vector_partitioners[from_level - 1].get())
354 this->ghosted_level_vector[from_level - 1].reinit(
355 this->vector_partitioners[from_level - 1]);
359 this->ghosted_level_vector[from_level - 1] = 0.;
363 src_inplace ? src : this->ghosted_level_vector[from_level];
365 dst_inplace ? dst : this->ghosted_level_vector[from_level - 1];
370 do_restrict_add<0>(from_level, dst_vec, src_vec);
371 else if (fe_degree == 1)
372 do_restrict_add<1>(from_level, dst_vec, src_vec);
373 else if (fe_degree == 2)
374 do_restrict_add<2>(from_level, dst_vec, src_vec);
375 else if (fe_degree == 3)
376 do_restrict_add<3>(from_level, dst_vec, src_vec);
377 else if (fe_degree == 4)
378 do_restrict_add<4>(from_level, dst_vec, src_vec);
379 else if (fe_degree == 5)
380 do_restrict_add<5>(from_level, dst_vec, src_vec);
381 else if (fe_degree == 6)
382 do_restrict_add<6>(from_level, dst_vec, src_vec);
383 else if (fe_degree == 7)
384 do_restrict_add<7>(from_level, dst_vec, src_vec);
385 else if (fe_degree == 8)
386 do_restrict_add<8>(from_level, dst_vec, src_vec);
387 else if (fe_degree == 9)
388 do_restrict_add<9>(from_level, dst_vec, src_vec);
389 else if (fe_degree == 10)
390 do_restrict_add<10>(from_level, dst_vec, src_vec);
393 do_restrict_add<-1>(from_level, dst_vec, src_vec);
396 if (dst_inplace ==
false)
399 if (src_inplace ==
true)
409 const unsigned int to_level,
414 const unsigned int degree_size = (degree > -1 ? degree : fe_degree) + 1;
415 const unsigned int n_child_dofs_1d = 2 * degree_size - element_is_continuous;
416 const unsigned int n_scalar_cell_dofs =
417 Utilities::fixed_power<dim>(n_child_dofs_1d);
425 if (this->mg_constrained_dofs !=
nullptr &&
426 this->mg_constrained_dofs->get_user_constraint_matrix(to_level - 1)
433 this->mg_constrained_dofs->get_user_constraint_matrix(to_level - 1)
434 .distribute(copy_src);
444 for (
unsigned int cell = 0; cell < n_owned_level_cells[to_level - 1];
447 const unsigned int n_lanes =
448 (cell + vec_size > n_owned_level_cells[to_level - 1]) ?
449 (n_owned_level_cells[to_level - 1] - cell) :
453 for (
unsigned int v = 0; v < n_lanes; ++v)
455 const unsigned int shift =
456 internal::MGTransfer::compute_shift_within_children<dim>(
457 parent_child_connect[to_level - 1][cell + v].
second,
458 fe_degree + 1 - element_is_continuous,
460 const unsigned int *indices =
461 &level_dof_indices[to_level - 1]
462 [parent_child_connect[to_level - 1][cell + v]
466 for (
unsigned int c = 0, m = 0; c < n_components; ++c)
468 for (
unsigned int k = 0; k < (dim > 2 ? degree_size : 1); ++k)
469 for (
unsigned int j = 0; j < (dim > 1 ? degree_size : 1); ++j)
470 for (
unsigned int i = 0; i < degree_size; ++i, ++m)
472 indices[c * n_scalar_cell_dofs +
473 k * n_child_dofs_1d * n_child_dofs_1d +
474 j * n_child_dofs_1d + i]);
477 for (std::vector<unsigned short>::const_iterator i =
478 dirichlet_indices[to_level - 1][cell + v].begin();
479 i != dirichlet_indices[to_level - 1][cell + v].end();
481 evaluation_data[*i][v] = 0.;
486 degree_size * n_child_dofs_1d);
488 if (element_is_continuous)
492 for (
int c = n_components - 1; c >= 0; --c)
498 2 * degree + 1>::do_forward(1,
499 prolongation_matrix_1d,
500 evaluation_data.begin() +
501 c * Utilities::fixed_power<dim>(
503 evaluation_data.begin() +
504 c * n_scalar_cell_dofs,
508 degree != -1 ? 2 * degree + 1 :
511 &weights_on_refined[to_level - 1][(cell / vec_size) * three_to_dim],
514 evaluation_data.begin());
518 for (
int c = n_components - 1; c >= 0; --c)
524 2 * degree + 2>::do_forward(1,
525 prolongation_matrix_1d,
526 evaluation_data.begin() +
527 c * Utilities::fixed_power<dim>(
529 evaluation_data.begin() +
530 c * n_scalar_cell_dofs,
536 const unsigned int *indices =
537 &level_dof_indices[to_level][cell * n_child_cell_dofs];
538 for (
unsigned int v = 0; v < n_lanes; ++v)
540 for (
unsigned int i = 0; i < n_child_cell_dofs; ++i)
542 indices += n_child_cell_dofs;
553 const unsigned int from_level,
558 const unsigned int degree_size = (degree > -1 ? degree : fe_degree) + 1;
559 const unsigned int n_child_dofs_1d = 2 * degree_size - element_is_continuous;
560 const unsigned int n_scalar_cell_dofs =
561 Utilities::fixed_power<dim>(n_child_dofs_1d);
564 for (
unsigned int cell = 0; cell < n_owned_level_cells[from_level - 1];
567 const unsigned int n_lanes =
568 (cell + vec_size > n_owned_level_cells[from_level - 1]) ?
569 (n_owned_level_cells[from_level - 1] - cell) :
574 const unsigned int *indices =
575 &level_dof_indices[from_level][cell * n_child_cell_dofs];
576 for (
unsigned int v = 0; v < n_lanes; ++v)
578 for (
unsigned int i = 0; i < n_child_cell_dofs; ++i)
580 indices += n_child_cell_dofs;
585 degree_size * n_child_dofs_1d);
587 if (element_is_continuous)
590 degree != -1 ? 2 * degree + 1 :
593 &weights_on_refined[from_level - 1]
594 [(cell / vec_size) * three_to_dim],
597 evaluation_data.data());
598 for (
unsigned int c = 0; c < n_components; ++c)
604 2 * degree + 1>::do_backward(1,
605 prolongation_matrix_1d,
607 evaluation_data.begin() +
608 c * n_scalar_cell_dofs,
609 evaluation_data.begin() +
610 c * Utilities::fixed_power<dim>(
617 for (
unsigned int c = 0; c < n_components; ++c)
623 2 * degree + 2>::do_backward(1,
624 prolongation_matrix_1d,
626 evaluation_data.begin() +
627 c * n_scalar_cell_dofs,
628 evaluation_data.begin() +
629 c * Utilities::fixed_power<dim>(
636 for (
unsigned int v = 0; v < n_lanes; ++v)
638 const unsigned int shift =
639 internal::MGTransfer::compute_shift_within_children<dim>(
640 parent_child_connect[from_level - 1][cell + v].
second,
641 fe_degree + 1 - element_is_continuous,
644 parent_child_connect[from_level - 1][cell + v].
first *
646 n_child_cell_dofs - 1,
647 level_dof_indices[from_level - 1].
size());
648 const unsigned int *indices =
649 &level_dof_indices[from_level - 1]
650 [parent_child_connect[from_level - 1][cell + v]
654 for (
unsigned int c = 0, m = 0; c < n_components; ++c)
657 for (std::vector<unsigned short>::const_iterator i =
658 dirichlet_indices[from_level - 1][cell + v].begin();
659 i != dirichlet_indices[from_level - 1][cell + v].end();
661 evaluation_data[*i][v] = 0.;
663 for (
unsigned int k = 0; k < (dim > 2 ? degree_size : 1); ++k)
664 for (
unsigned int j = 0; j < (dim > 1 ? degree_size : 1); ++j)
665 for (
unsigned int i = 0; i < degree_size; ++i, ++m)
667 indices[c * n_scalar_cell_dofs +
668 k * n_child_dofs_1d * n_child_dofs_1d +
669 j * n_child_dofs_1d + i]) +=
670 evaluation_data[m][v];