Reference documentation for deal.II version Git 3f1f337db3 2021-10-23 13:19:02 -0600
\(\newcommand{\dealvcentcolon}{\mathrel{\mathop{:}}}\) \(\newcommand{\dealcoloneq}{\dealvcentcolon\mathrel{\mkern-1.2mu}=}\) \(\newcommand{\jump}[1]{\left[\!\left[ #1 \right]\!\right]}\) \(\newcommand{\average}[1]{\left\{\!\left\{ #1 \right\}\!\right\}}\)
mg_transfer_matrix_free.cc
Go to the documentation of this file.
1 // ---------------------------------------------------------------------
2 //
3 // Copyright (C) 2016 - 2021 by the deal.II authors
4 //
5 // This file is part of the deal.II library.
6 //
7 // The deal.II library is free software; you can use it, redistribute
8 // it, and/or modify it under the terms of the GNU Lesser General
9 // Public License as published by the Free Software Foundation; either
10 // version 2.1 of the License, or (at your option) any later version.
11 // The full text of the license can be found in the file LICENSE.md at
12 // the top level directory of deal.II.
13 //
14 // ---------------------------------------------------------------------
15 
16 
17 #include <deal.II/base/function.h>
18 #include <deal.II/base/logstream.h>
20 
22 #include <deal.II/dofs/dof_tools.h>
23 
24 #include <deal.II/fe/fe.h>
25 
26 #include <deal.II/grid/tria.h>
28 
30 
32 
36 
37 #include <algorithm>
38 
40 
41 
42 template <int dim, typename Number>
44  : fe_degree(0)
45  , element_is_continuous(false)
46  , n_components(0)
47  , n_child_cell_dofs(0)
48 {}
49 
50 
51 
52 template <int dim, typename Number>
54  const MGConstrainedDoFs &mg_c)
55  : fe_degree(0)
56  , element_is_continuous(false)
57  , n_components(0)
59 {
60  this->mg_constrained_dofs = &mg_c;
61 }
62 
63 
64 
65 template <int dim, typename Number>
66 void
68  const MGConstrainedDoFs &mg_c)
69 {
70  this->mg_constrained_dofs = &mg_c;
71 }
72 
73 
74 
75 template <int dim, typename Number>
76 void
78 {
81  fe_degree = 0;
82  element_is_continuous = false;
83  n_components = 0;
85  level_dof_indices.clear();
86  parent_child_connect.clear();
87  dirichlet_indices.clear();
88  n_owned_level_cells.clear();
91  weights_on_refined.clear();
92 }
93 
94 
95 template <int dim, typename Number>
96 void
98  const DoFHandler<dim, dim> &dof_handler,
99  const std::vector<std::shared_ptr<const Utilities::MPI::Partitioner>>
100  &external_partitioners)
101 {
102  Assert(dof_handler.has_level_dofs(),
103  ExcMessage(
104  "The underlying DoFHandler object has not had its "
105  "distribute_mg_dofs() function called, but this is a prerequisite "
106  "for multigrid transfers. You will need to call this function, "
107  "probably close to where you already call distribute_dofs()."));
108 
109  this->fill_and_communicate_copy_indices(dof_handler);
110 
112  dof_handler.get_triangulation().n_global_levels() -
113  1);
114  for (unsigned int level = 0; level <= this->ghosted_level_vector.max_level();
115  ++level)
117  this->ghosted_level_vector[level].get_partitioner();
118 
119  std::vector<std::vector<Number>> weights_unvectorized;
120 
122 
123  internal::MGTransfer::setup_transfer<dim, Number>(
124  dof_handler,
125  this->mg_constrained_dofs,
126  external_partitioners,
127  elem_info,
132  weights_unvectorized,
135 
136  // reinit the ghosted level vector in case its partitioner differs from the
137  // one the user intends to use externally
139  for (unsigned int level = 0; level <= vector_partitioners.max_level();
140  ++level)
141  if (external_partitioners.size() == vector_partitioners.max_level() + 1 &&
142  external_partitioners[level].get() == vector_partitioners[level].get())
144  else
146 
147  // unpack element info data
148  fe_degree = elem_info.fe_degree;
149  element_is_continuous = elem_info.element_is_continuous;
150  n_components = elem_info.n_components;
151  n_child_cell_dofs = elem_info.n_child_cell_dofs;
152 
153  // duplicate and put into vectorized array
154  prolongation_matrix_1d.resize(elem_info.prolongation_matrix_1d.size());
155  for (unsigned int i = 0; i < elem_info.prolongation_matrix_1d.size(); ++i)
156  prolongation_matrix_1d[i] = elem_info.prolongation_matrix_1d[i];
157 
158  // reshuffle into aligned vector of vectorized arrays
159  const unsigned int vec_size = VectorizedArray<Number>::size();
160  const unsigned int n_levels =
161  dof_handler.get_triangulation().n_global_levels();
162 
163  const unsigned int n_weights_per_cell = Utilities::fixed_power<dim>(3);
164  weights_on_refined.resize(n_levels - 1);
165  for (unsigned int level = 1; level < n_levels; ++level)
166  {
167  weights_on_refined[level - 1].resize(
168  ((n_owned_level_cells[level - 1] + vec_size - 1) / vec_size) *
169  n_weights_per_cell);
170 
171  for (unsigned int c = 0; c < n_owned_level_cells[level - 1]; ++c)
172  {
173  const unsigned int comp = c / vec_size;
174  const unsigned int v = c % vec_size;
175  for (unsigned int i = 0; i < n_weights_per_cell; ++i)
176  {
177  weights_on_refined[level - 1][comp * n_weights_per_cell + i][v] =
178  weights_unvectorized[level - 1][c * n_weights_per_cell + i];
179  }
180  }
181  }
182 
184 }
185 
186 
187 
188 template <int dim, typename Number>
189 void
191  const unsigned int to_level,
194 {
195  dst = 0;
196  prolongate_and_add(to_level, dst, src);
197 }
198 
199 
200 
201 template <int dim, typename Number>
202 void
204  const unsigned int to_level,
207 {
208  Assert((to_level >= 1) && (to_level <= level_dof_indices.size()),
209  ExcIndexRange(to_level, 1, level_dof_indices.size() + 1));
210 
211  const bool src_inplace = src.get_partitioner().get() ==
212  this->vector_partitioners[to_level - 1].get();
213  if (src_inplace == false)
214  {
215  if (this->ghosted_level_vector[to_level - 1].get_partitioner().get() !=
216  this->vector_partitioners[to_level - 1].get())
217  this->ghosted_level_vector[to_level - 1].reinit(
218  this->vector_partitioners[to_level - 1]);
219  this->ghosted_level_vector[to_level - 1].copy_locally_owned_data_from(
220  src);
221  }
222 
223  const bool dst_inplace =
224  dst.get_partitioner().get() == this->vector_partitioners[to_level].get();
225  if (dst_inplace == false)
226  {
227  if (this->ghosted_level_vector[to_level].get_partitioner().get() !=
228  this->vector_partitioners[to_level].get())
229  this->ghosted_level_vector[to_level].reinit(
230  this->vector_partitioners[to_level]);
231  AssertDimension(this->ghosted_level_vector[to_level].locally_owned_size(),
232  dst.locally_owned_size());
233  this->ghosted_level_vector[to_level] = 0.;
234  }
235 
237  src_inplace ? src : this->ghosted_level_vector[to_level - 1];
239  dst_inplace ? dst : this->ghosted_level_vector[to_level];
240 
241  src_vec.update_ghost_values();
242  // the implementation in do_prolongate_add is templated in the degree of the
243  // element (for efficiency reasons), so we need to find the appropriate
244  // kernel here...
245  if (fe_degree == 0)
246  do_prolongate_add<0>(to_level, dst_vec, src_vec);
247  else if (fe_degree == 1)
248  do_prolongate_add<1>(to_level, dst_vec, src_vec);
249  else if (fe_degree == 2)
250  do_prolongate_add<2>(to_level, dst_vec, src_vec);
251  else if (fe_degree == 3)
252  do_prolongate_add<3>(to_level, dst_vec, src_vec);
253  else if (fe_degree == 4)
254  do_prolongate_add<4>(to_level, dst_vec, src_vec);
255  else if (fe_degree == 5)
256  do_prolongate_add<5>(to_level, dst_vec, src_vec);
257  else if (fe_degree == 6)
258  do_prolongate_add<6>(to_level, dst_vec, src_vec);
259  else if (fe_degree == 7)
260  do_prolongate_add<7>(to_level, dst_vec, src_vec);
261  else if (fe_degree == 8)
262  do_prolongate_add<8>(to_level, dst_vec, src_vec);
263  else if (fe_degree == 9)
264  do_prolongate_add<9>(to_level, dst_vec, src_vec);
265  else if (fe_degree == 10)
266  do_prolongate_add<10>(to_level, dst_vec, src_vec);
267  else
268  do_prolongate_add<-1>(to_level, dst_vec, src_vec);
269 
270  dst_vec.compress(VectorOperation::add);
271  if (dst_inplace == false)
272  dst += dst_vec;
273 
274  if (src_inplace == true)
275  src.zero_out_ghost_values();
276 }
277 
278 
279 
280 template <int dim, typename Number>
281 void
283  const unsigned int from_level,
286 {
287  Assert((from_level >= 1) && (from_level <= level_dof_indices.size()),
288  ExcIndexRange(from_level, 1, level_dof_indices.size() + 1));
289 
290  const bool src_inplace =
291  src.get_partitioner().get() == this->vector_partitioners[from_level].get();
292  if (src_inplace == false)
293  {
294  if (this->ghosted_level_vector[from_level].get_partitioner().get() !=
295  this->vector_partitioners[from_level].get())
296  this->ghosted_level_vector[from_level].reinit(
297  this->vector_partitioners[from_level]);
298  this->ghosted_level_vector[from_level].copy_locally_owned_data_from(src);
299  }
300 
301  const bool dst_inplace = dst.get_partitioner().get() ==
302  this->vector_partitioners[from_level - 1].get();
303  if (dst_inplace == false)
304  {
305  if (this->ghosted_level_vector[from_level - 1].get_partitioner().get() !=
306  this->vector_partitioners[from_level - 1].get())
307  this->ghosted_level_vector[from_level - 1].reinit(
308  this->vector_partitioners[from_level - 1]);
310  this->ghosted_level_vector[from_level - 1].locally_owned_size(),
311  dst.locally_owned_size());
312  this->ghosted_level_vector[from_level - 1] = 0.;
313  }
314 
316  src_inplace ? src : this->ghosted_level_vector[from_level];
318  dst_inplace ? dst : this->ghosted_level_vector[from_level - 1];
319 
320  src_vec.update_ghost_values();
321 
322  if (fe_degree == 0)
323  do_restrict_add<0>(from_level, dst_vec, src_vec);
324  else if (fe_degree == 1)
325  do_restrict_add<1>(from_level, dst_vec, src_vec);
326  else if (fe_degree == 2)
327  do_restrict_add<2>(from_level, dst_vec, src_vec);
328  else if (fe_degree == 3)
329  do_restrict_add<3>(from_level, dst_vec, src_vec);
330  else if (fe_degree == 4)
331  do_restrict_add<4>(from_level, dst_vec, src_vec);
332  else if (fe_degree == 5)
333  do_restrict_add<5>(from_level, dst_vec, src_vec);
334  else if (fe_degree == 6)
335  do_restrict_add<6>(from_level, dst_vec, src_vec);
336  else if (fe_degree == 7)
337  do_restrict_add<7>(from_level, dst_vec, src_vec);
338  else if (fe_degree == 8)
339  do_restrict_add<8>(from_level, dst_vec, src_vec);
340  else if (fe_degree == 9)
341  do_restrict_add<9>(from_level, dst_vec, src_vec);
342  else if (fe_degree == 10)
343  do_restrict_add<10>(from_level, dst_vec, src_vec);
344  else
345  // go to the non-templated version of the evaluator
346  do_restrict_add<-1>(from_level, dst_vec, src_vec);
347 
348  dst_vec.compress(VectorOperation::add);
349  if (dst_inplace == false)
350  dst += dst_vec;
351 
352  if (src_inplace == true)
353  src.zero_out_ghost_values();
354 }
355 
356 
357 
358 namespace
359 {
360  template <int dim, int degree, typename Number>
361  void
362  weight_dofs_on_child(const VectorizedArray<Number> *weights,
363  const unsigned int n_components,
364  const unsigned int fe_degree,
366  {
367  Assert(fe_degree > 0, ExcNotImplemented());
368  Assert(fe_degree < 100, ExcNotImplemented());
369  const int loop_length = degree != -1 ? 2 * degree + 1 : 2 * fe_degree + 1;
370  unsigned int degree_to_3[100];
371  degree_to_3[0] = 0;
372  for (int i = 1; i < loop_length - 1; ++i)
373  degree_to_3[i] = 1;
374  degree_to_3[loop_length - 1] = 2;
375  for (unsigned int c = 0; c < n_components; ++c)
376  for (int k = 0; k < (dim > 2 ? loop_length : 1); ++k)
377  for (int j = 0; j < (dim > 1 ? loop_length : 1); ++j)
378  {
379  const unsigned int shift = 9 * degree_to_3[k] + 3 * degree_to_3[j];
380  data[0] *= weights[shift];
381  // loop bound as int avoids compiler warnings in case loop_length
382  // == 1 (polynomial degree 0)
383  for (int i = 1; i < loop_length - 1; ++i)
384  data[i] *= weights[shift + 1];
385  data[loop_length - 1] *= weights[shift + 2];
386  data += loop_length;
387  }
388  }
389 } // namespace
390 
391 
392 
393 template <int dim, typename Number>
394 template <int degree>
395 void
397  const unsigned int to_level,
400 {
401  const unsigned int vec_size = VectorizedArray<Number>::size();
402  const unsigned int degree_size = (degree > -1 ? degree : fe_degree) + 1;
403  const unsigned int n_child_dofs_1d = 2 * degree_size - element_is_continuous;
404  const unsigned int n_scalar_cell_dofs =
405  Utilities::fixed_power<dim>(n_child_dofs_1d);
406  constexpr unsigned int three_to_dim = Utilities::pow(3, dim);
407 
408  // If we have user defined MG constraints, we must create
409  // a non-const, ghosted version of the source vector to distribute
410  // constraints.
411  const LinearAlgebra::distributed::Vector<Number> *to_use = &src;
413  if (this->mg_constrained_dofs != nullptr &&
415  .get_local_lines()
416  .size() > 0)
417  {
419 
420  // Distribute any user defined constraints
422  .distribute(copy_src);
423 
424  // Re-initialize new ghosted vector with correct constraints
425  new_src.reinit(copy_src);
426  new_src = copy_src;
427  new_src.update_ghost_values();
428 
429  to_use = &new_src;
430  }
431 
432  for (unsigned int cell = 0; cell < n_owned_level_cells[to_level - 1];
433  cell += vec_size)
434  {
435  const unsigned int n_lanes =
436  (cell + vec_size > n_owned_level_cells[to_level - 1]) ?
437  (n_owned_level_cells[to_level - 1] - cell) :
438  vec_size;
439 
440  // read from source vector
441  for (unsigned int v = 0; v < n_lanes; ++v)
442  {
443  const unsigned int shift =
444  internal::MGTransfer::compute_shift_within_children<dim>(
445  parent_child_connect[to_level - 1][cell + v].second,
447  fe_degree);
448  const unsigned int *indices =
449  &level_dof_indices[to_level - 1]
450  [parent_child_connect[to_level - 1][cell + v]
451  .first *
453  shift];
454  for (unsigned int c = 0, m = 0; c < n_components; ++c)
455  {
456  for (unsigned int k = 0; k < (dim > 2 ? degree_size : 1); ++k)
457  for (unsigned int j = 0; j < (dim > 1 ? degree_size : 1); ++j)
458  for (unsigned int i = 0; i < degree_size; ++i, ++m)
459  evaluation_data[m][v] = to_use->local_element(
460  indices[c * n_scalar_cell_dofs +
461  k * n_child_dofs_1d * n_child_dofs_1d +
462  j * n_child_dofs_1d + i]);
463 
464  // apply Dirichlet boundary conditions on parent cell
465  for (std::vector<unsigned short>::const_iterator i =
466  dirichlet_indices[to_level - 1][cell + v].begin();
467  i != dirichlet_indices[to_level - 1][cell + v].end();
468  ++i)
469  evaluation_data[*i][v] = 0.;
470  }
471  }
472 
474  degree_size * n_child_dofs_1d);
475  // perform tensorized operation
476  if (element_is_continuous)
477  {
478  // must go through the components backwards because we want to write
479  // the output to the same array as the input
480  for (int c = n_components - 1; c >= 0; --c)
484  dim,
485  degree + 1,
486  2 * degree + 1,
488  VectorizedArray<Number>>::do_forward(1,
492  dim>(degree_size),
494  c * n_scalar_cell_dofs,
495  fe_degree + 1,
496  2 * fe_degree + 1);
497  weight_dofs_on_child<dim, degree, Number>(
498  &weights_on_refined[to_level - 1][(cell / vec_size) * three_to_dim],
499  n_components,
500  fe_degree,
502  }
503  else
504  {
505  for (int c = n_components - 1; c >= 0; --c)
509  dim,
510  degree + 1,
511  2 * degree + 2,
513  VectorizedArray<Number>>::do_forward(1,
517  dim>(degree_size),
519  c * n_scalar_cell_dofs,
520  fe_degree + 1,
521  2 * fe_degree + 2);
522  }
523 
524  // write into dst vector
525  const unsigned int *indices =
526  &level_dof_indices[to_level][cell * n_child_cell_dofs];
527  for (unsigned int v = 0; v < n_lanes; ++v)
528  {
529  for (unsigned int i = 0; i < n_child_cell_dofs; ++i)
530  dst.local_element(indices[i]) += evaluation_data[i][v];
531  indices += n_child_cell_dofs;
532  }
533  }
534 }
535 
536 
537 
538 template <int dim, typename Number>
539 template <int degree>
540 void
542  const unsigned int from_level,
545 {
546  const unsigned int vec_size = VectorizedArray<Number>::size();
547  const unsigned int degree_size = (degree > -1 ? degree : fe_degree) + 1;
548  const unsigned int n_child_dofs_1d = 2 * degree_size - element_is_continuous;
549  const unsigned int n_scalar_cell_dofs =
550  Utilities::fixed_power<dim>(n_child_dofs_1d);
551  constexpr unsigned int three_to_dim = Utilities::pow(3, dim);
552 
553  for (unsigned int cell = 0; cell < n_owned_level_cells[from_level - 1];
554  cell += vec_size)
555  {
556  const unsigned int n_lanes =
557  (cell + vec_size > n_owned_level_cells[from_level - 1]) ?
558  (n_owned_level_cells[from_level - 1] - cell) :
559  vec_size;
560 
561  // read from source vector
562  {
563  const unsigned int *indices =
564  &level_dof_indices[from_level][cell * n_child_cell_dofs];
565  for (unsigned int v = 0; v < n_lanes; ++v)
566  {
567  for (unsigned int i = 0; i < n_child_cell_dofs; ++i)
568  evaluation_data[i][v] = src.local_element(indices[i]);
569  indices += n_child_cell_dofs;
570  }
571  }
572 
574  degree_size * n_child_dofs_1d);
575  // perform tensorized operation
576  if (element_is_continuous)
577  {
578  weight_dofs_on_child<dim, degree, Number>(
579  &weights_on_refined[from_level - 1]
580  [(cell / vec_size) * three_to_dim],
581  n_components,
582  fe_degree,
584  for (unsigned int c = 0; c < n_components; ++c)
588  dim,
589  degree + 1,
590  2 * degree + 1,
592  VectorizedArray<Number>>::do_backward(1,
594  false,
596  c * n_scalar_cell_dofs,
598  c *
600  dim>(degree_size),
601  fe_degree + 1,
602  2 * fe_degree + 1);
603  }
604  else
605  {
606  for (unsigned int c = 0; c < n_components; ++c)
610  dim,
611  degree + 1,
612  2 * degree + 2,
614  VectorizedArray<Number>>::do_backward(1,
616  false,
618  c * n_scalar_cell_dofs,
620  c *
622  dim>(degree_size),
623  fe_degree + 1,
624  2 * fe_degree + 2);
625  }
626 
627  // write into dst vector
628  for (unsigned int v = 0; v < n_lanes; ++v)
629  {
630  const unsigned int shift =
631  internal::MGTransfer::compute_shift_within_children<dim>(
632  parent_child_connect[from_level - 1][cell + v].second,
634  fe_degree);
636  parent_child_connect[from_level - 1][cell + v].first *
638  n_child_cell_dofs - 1,
639  level_dof_indices[from_level - 1].size());
640  const unsigned int *indices =
641  &level_dof_indices[from_level - 1]
642  [parent_child_connect[from_level - 1][cell + v]
643  .first *
645  shift];
646  for (unsigned int c = 0, m = 0; c < n_components; ++c)
647  {
648  // apply Dirichlet boundary conditions on parent cell
649  for (std::vector<unsigned short>::const_iterator i =
650  dirichlet_indices[from_level - 1][cell + v].begin();
651  i != dirichlet_indices[from_level - 1][cell + v].end();
652  ++i)
653  evaluation_data[*i][v] = 0.;
654 
655  for (unsigned int k = 0; k < (dim > 2 ? degree_size : 1); ++k)
656  for (unsigned int j = 0; j < (dim > 1 ? degree_size : 1); ++j)
657  for (unsigned int i = 0; i < degree_size; ++i, ++m)
658  dst.local_element(
659  indices[c * n_scalar_cell_dofs +
660  k * n_child_dofs_1d * n_child_dofs_1d +
661  j * n_child_dofs_1d + i]) +=
662  evaluation_data[m][v];
663  }
664  }
665  }
666 }
667 
668 
669 
670 template <int dim, typename Number>
671 std::size_t
673 {
674  std::size_t memory = MGLevelGlobalTransfer<
683  return memory;
684 }
685 
686 
687 
688 template <int dim, typename Number>
690  const MGConstrainedDoFs &mg_c)
691  : same_for_all(true)
692 {
693  matrix_free_transfer_vector.emplace_back(mg_c);
694 }
695 
696 
697 
698 template <int dim, typename Number>
700  const std::vector<MGConstrainedDoFs> &mg_c)
701  : same_for_all(false)
702 {
703  for (const auto &constrained_block_dofs : mg_c)
704  matrix_free_transfer_vector.emplace_back(constrained_block_dofs);
705 }
706 
707 
708 
709 template <int dim, typename Number>
710 void
712  const MGConstrainedDoFs &mg_c)
713 {
715  ExcMessage("This object was initialized with support for usage with "
716  "one DoFHandler for each block, but this method assumes "
717  "that the same DoFHandler is used for all the blocks!"));
719 
720  matrix_free_transfer_vector[0].initialize_constraints(mg_c);
721 }
722 
723 
724 
725 template <int dim, typename Number>
726 void
728  const std::vector<MGConstrainedDoFs> &mg_c)
729 {
731  ExcMessage("This object was initialized with support for using "
732  "the same DoFHandler for all the blocks, but this "
733  "method assumes that there is a separate DoFHandler "
734  "for each block!"));
735  AssertDimension(matrix_free_transfer_vector.size(), mg_c.size());
736 
737  for (unsigned int i = 0; i < mg_c.size(); ++i)
739 }
740 
741 
742 
743 template <int dim, typename Number>
744 void
746 {
748 }
749 
750 
751 
752 template <int dim, typename Number>
753 void
755  const DoFHandler<dim, dim> &dof_handler)
756 {
758  matrix_free_transfer_vector[0].build(dof_handler);
759 }
760 
761 
762 
763 template <int dim, typename Number>
764 void
766  const std::vector<const DoFHandler<dim, dim> *> &dof_handler)
767 {
768  AssertDimension(matrix_free_transfer_vector.size(), dof_handler.size());
769  for (unsigned int i = 0; i < dof_handler.size(); ++i)
770  matrix_free_transfer_vector[i].build(*dof_handler[i]);
771 }
772 
773 
774 
775 template <int dim, typename Number>
776 void
778  const unsigned int to_level,
781 {
782  const unsigned int n_blocks = src.n_blocks();
784 
785  if (!same_for_all)
787 
788  for (unsigned int b = 0; b < n_blocks; ++b)
789  {
790  const unsigned int data_block = same_for_all ? 0 : b;
791  matrix_free_transfer_vector[data_block].prolongate(to_level,
792  dst.block(b),
793  src.block(b));
794  }
795 }
796 
797 
798 
799 template <int dim, typename Number>
800 void
802  const unsigned int to_level,
805 {
806  const unsigned int n_blocks = src.n_blocks();
808 
809  if (!same_for_all)
811 
812  for (unsigned int b = 0; b < n_blocks; ++b)
813  {
814  const unsigned int data_block = same_for_all ? 0 : b;
815  matrix_free_transfer_vector[data_block].prolongate_and_add(to_level,
816  dst.block(b),
817  src.block(b));
818  }
819 }
820 
821 
822 
823 template <int dim, typename Number>
824 void
826  const unsigned int from_level,
829 {
830  const unsigned int n_blocks = src.n_blocks();
832 
833  if (!same_for_all)
835 
836  for (unsigned int b = 0; b < n_blocks; ++b)
837  {
838  const unsigned int data_block = same_for_all ? 0 : b;
839  matrix_free_transfer_vector[data_block].restrict_and_add(from_level,
840  dst.block(b),
841  src.block(b));
842  }
843 }
844 
845 
846 
847 template <int dim, typename Number>
848 std::size_t
850 {
851  std::size_t total_memory_consumption = 0;
852  for (const auto &el : matrix_free_transfer_vector)
853  total_memory_consumption += el.memory_consumption();
854  return total_memory_consumption;
855 }
856 
857 
858 
859 // explicit instantiation
860 #include "mg_transfer_matrix_free.inst"
861 
862 
unsigned int max_level() const
void resize(const size_type new_size)
const Triangulation< dim, spacedim > & get_triangulation() const
virtual void prolongate_and_add(const unsigned int to_level, LinearAlgebra::distributed::Vector< Number > &dst, const LinearAlgebra::distributed::Vector< Number > &src) const override
std::vector< Table< 2, unsigned int > > copy_indices_global_mine
Definition: mg_transfer.h:537
void reinit(MatrixBlock< MatrixType > &v, const BlockSparsityPattern &p)
Definition: matrix_block.h:618
void do_restrict_add(const unsigned int from_level, LinearAlgebra::distributed::Vector< Number > &dst, const LinearAlgebra::distributed::Vector< Number > &src) const
#define AssertDimension(dim1, dim2)
Definition: exceptions.h:1655
std::vector< unsigned int > n_owned_level_cells
void build(const DoFHandler< dim, dim > &dof_handler, const std::vector< std::shared_ptr< const Utilities::MPI::Partitioner >> &external_partitioners=std::vector< std::shared_ptr< const Utilities::MPI::Partitioner >>())
AlignedVector< VectorizedArray< Number > > evaluation_data
std::size_t memory_consumption() const
void build(const DoFHandler< dim, dim > &dof_handler)
std::vector< std::vector< std::vector< unsigned short > > > dirichlet_indices
void initialize_constraints(const MGConstrainedDoFs &mg_constrained_dofs)
#define AssertIndexRange(index, range)
Definition: exceptions.h:1720
pointer data()
void reinit(const size_type size, const bool omit_zeroing_entries=false)
void resize(const unsigned int new_minlevel, const unsigned int new_maxlevel, Args &&...args)
size_type locally_owned_size() const
std::size_t memory_consumption() const
constexpr T pow(const T base, const int iexp)
Definition: utilities.h:461
size_type size() const
Definition: index_set.h:1639
T fixed_power(const T t)
Definition: utilities.h:1082
virtual void prolongate_and_add(const unsigned int to_level, LinearAlgebra::distributed::BlockVector< Number > &dst, const LinearAlgebra::distributed::BlockVector< Number > &src) const override
Number local_element(const size_type local_index) const
static ::ExceptionBase & ExcMessage(std::string arg1)
AlignedVector< VectorizedArray< Number > > prolongation_matrix_1d
#define Assert(cond, exc)
Definition: exceptions.h:1461
virtual void restrict_and_add(const unsigned int from_level, LinearAlgebra::distributed::Vector< Number > &dst, const LinearAlgebra::distributed::Vector< Number > &src) const override
virtual void restrict_and_add(const unsigned int from_level, LinearAlgebra::distributed::BlockVector< Number > &dst, const LinearAlgebra::distributed::BlockVector< Number > &src) const override
virtual void prolongate(const unsigned int to_level, LinearAlgebra::distributed::Vector< Number > &dst, const LinearAlgebra::distributed::Vector< Number > &src) const override
#define DEAL_II_NAMESPACE_CLOSE
Definition: config.h:401
unsigned int level
Definition: grid_out.cc:4589
std::vector< AlignedVector< VectorizedArray< Number > > > weights_on_refined
std::enable_if< IsBlockVector< VectorType >::value, unsigned int >::type n_blocks(const VectorType &vector)
Definition: operators.h:50
void fill_and_communicate_copy_indices(const DoFHandler< dim, spacedim > &dof_handler)
std::vector< std::vector< unsigned int > > level_dof_indices
const IndexSet & get_local_lines() const
SymmetricTensor< 2, dim, Number > b(const Tensor< 2, dim, Number > &F)
Point< 2 > first
Definition: grid_out.cc:4586
void do_prolongate_add(const unsigned int to_level, LinearAlgebra::distributed::Vector< Number > &dst, const LinearAlgebra::distributed::Vector< Number > &src) const
SmartPointer< const MGConstrainedDoFs, MGLevelGlobalTransfer< LinearAlgebra::distributed::Vector< Number > > > mg_constrained_dofs
Definition: mg_transfer.h:586
void distribute(VectorType &vec) const
std::vector< std::vector< std::pair< unsigned int, unsigned int > > > parent_child_connect
virtual void prolongate(const unsigned int to_level, LinearAlgebra::distributed::BlockVector< Number > &dst, const LinearAlgebra::distributed::BlockVector< Number > &src) const override
#define DEAL_II_NAMESPACE_OPEN
Definition: config.h:400
static ::ExceptionBase & ExcIndexRange(std::size_t arg1, std::size_t arg2, std::size_t arg3)
VectorType::value_type * begin(VectorType &V)
std::vector< MGTransferMatrixFree< dim, Number > > matrix_free_transfer_vector
const std::shared_ptr< const Utilities::MPI::Partitioner > & get_partitioner() const
MGTransferBlockMatrixFree()=default
iterator begin()
size_type size() const
static ::ExceptionBase & ExcNotImplemented()
MGLevelObject< std::shared_ptr< const Utilities::MPI::Partitioner > > vector_partitioners
void initialize_constraints(const MGConstrainedDoFs &mg_constrained_dofs)
const AffineConstraints< double > & get_user_constraint_matrix(const unsigned int level) const
bool has_level_dofs() const
BlockType & block(const unsigned int i)
MGLevelObject< LinearAlgebra::distributed::Vector< Number > > ghosted_level_vector
Definition: mg_transfer.h:606
void shift(const Tensor< 1, spacedim > &shift_vector, Triangulation< dim, spacedim > &triangulation)
Definition: grid_tools.cc:2013
std::enable_if< std::is_fundamental< T >::value, std::size_t >::type memory_consumption(const T &t)