Reference documentation for deal.II version 9.1.1
\(\newcommand{\dealcoloneq}{\mathrel{\vcenter{:}}=}\)
mg_transfer_matrix_free.cc
1 // ---------------------------------------------------------------------
2 //
3 // Copyright (C) 2016 - 2019 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>
19 #include <deal.II/base/vectorization.h>
20 
21 #include <deal.II/dofs/dof_accessor.h>
22 #include <deal.II/dofs/dof_tools.h>
23 
24 #include <deal.II/fe/fe.h>
25 
26 #include <deal.II/grid/tria.h>
27 #include <deal.II/grid/tria_iterator.h>
28 
29 #include <deal.II/lac/la_parallel_vector.h>
30 
31 #include <deal.II/matrix_free/evaluation_kernels.h>
32 
33 #include <deal.II/multigrid/mg_tools.h>
34 #include <deal.II/multigrid/mg_transfer_internal.h>
35 #include <deal.II/multigrid/mg_transfer_matrix_free.h>
36 
37 #include <algorithm>
38 
39 DEAL_II_NAMESPACE_OPEN
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)
58  , n_child_cell_dofs(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;
84  n_child_cell_dofs = 0;
85  level_dof_indices.clear();
86  parent_child_connect.clear();
87  dirichlet_indices.clear();
88  n_owned_level_cells.clear();
89  prolongation_matrix_1d.clear();
90  evaluation_data.clear();
91  weights_on_refined.clear();
92 }
93 
94 
95 template <int dim, typename Number>
96 void
98 {
99  this->fill_and_communicate_copy_indices(mg_dof);
100 
101  std::vector<std::vector<Number>> weights_unvectorized;
102 
104 
105  internal::MGTransfer::setup_transfer<dim, Number>(
106  mg_dof,
107  this->mg_constrained_dofs,
108  elem_info,
109  level_dof_indices,
110  parent_child_connect,
111  n_owned_level_cells,
112  dirichlet_indices,
113  weights_unvectorized,
114  this->copy_indices_global_mine,
115  this->ghosted_level_vector);
116  // unpack element info data
117  fe_degree = elem_info.fe_degree;
118  element_is_continuous = elem_info.element_is_continuous;
119  n_components = elem_info.n_components;
120  n_child_cell_dofs = elem_info.n_child_cell_dofs;
121 
122  // duplicate and put into vectorized array
123  prolongation_matrix_1d.resize(elem_info.prolongation_matrix_1d.size());
124  for (unsigned int i = 0; i < elem_info.prolongation_matrix_1d.size(); i++)
125  prolongation_matrix_1d[i] = elem_info.prolongation_matrix_1d[i];
126 
127  // reshuffle into aligned vector of vectorized arrays
128  const unsigned int vec_size = VectorizedArray<Number>::n_array_elements;
129  const unsigned int n_levels = mg_dof.get_triangulation().n_global_levels();
130 
131  const unsigned int n_weights_per_cell = Utilities::fixed_power<dim>(3);
132  weights_on_refined.resize(n_levels - 1);
133  for (unsigned int level = 1; level < n_levels; ++level)
134  {
135  weights_on_refined[level - 1].resize(
136  ((n_owned_level_cells[level - 1] + vec_size - 1) / vec_size) *
137  n_weights_per_cell);
138 
139  for (unsigned int c = 0; c < n_owned_level_cells[level - 1]; ++c)
140  {
141  const unsigned int comp = c / vec_size;
142  const unsigned int v = c % vec_size;
143  for (unsigned int i = 0; i < n_weights_per_cell; ++i)
144  {
145  weights_on_refined[level - 1][comp * n_weights_per_cell + i][v] =
146  weights_unvectorized[level - 1][c * n_weights_per_cell + i];
147  }
148  }
149  }
150 
151  evaluation_data.resize(n_child_cell_dofs);
152 }
153 
154 
155 
156 template <int dim, typename Number>
157 void
159  const unsigned int to_level,
162 {
163  Assert((to_level >= 1) && (to_level <= level_dof_indices.size()),
164  ExcIndexRange(to_level, 1, level_dof_indices.size() + 1));
165 
166  AssertDimension(this->ghosted_level_vector[to_level].local_size(),
167  dst.local_size());
168  AssertDimension(this->ghosted_level_vector[to_level - 1].local_size(),
169  src.local_size());
170 
171  this->ghosted_level_vector[to_level - 1].copy_locally_owned_data_from(src);
172  this->ghosted_level_vector[to_level - 1].update_ghost_values();
173  this->ghosted_level_vector[to_level] = 0.;
174 
175  // the implementation in do_prolongate_add is templated in the degree of the
176  // element (for efficiency reasons), so we need to find the appropriate
177  // kernel here...
178  if (fe_degree == 0)
179  do_prolongate_add<0>(to_level,
180  this->ghosted_level_vector[to_level],
181  this->ghosted_level_vector[to_level - 1]);
182  else if (fe_degree == 1)
183  do_prolongate_add<1>(to_level,
184  this->ghosted_level_vector[to_level],
185  this->ghosted_level_vector[to_level - 1]);
186  else if (fe_degree == 2)
187  do_prolongate_add<2>(to_level,
188  this->ghosted_level_vector[to_level],
189  this->ghosted_level_vector[to_level - 1]);
190  else if (fe_degree == 3)
191  do_prolongate_add<3>(to_level,
192  this->ghosted_level_vector[to_level],
193  this->ghosted_level_vector[to_level - 1]);
194  else if (fe_degree == 4)
195  do_prolongate_add<4>(to_level,
196  this->ghosted_level_vector[to_level],
197  this->ghosted_level_vector[to_level - 1]);
198  else if (fe_degree == 5)
199  do_prolongate_add<5>(to_level,
200  this->ghosted_level_vector[to_level],
201  this->ghosted_level_vector[to_level - 1]);
202  else if (fe_degree == 6)
203  do_prolongate_add<6>(to_level,
204  this->ghosted_level_vector[to_level],
205  this->ghosted_level_vector[to_level - 1]);
206  else if (fe_degree == 7)
207  do_prolongate_add<7>(to_level,
208  this->ghosted_level_vector[to_level],
209  this->ghosted_level_vector[to_level - 1]);
210  else if (fe_degree == 8)
211  do_prolongate_add<8>(to_level,
212  this->ghosted_level_vector[to_level],
213  this->ghosted_level_vector[to_level - 1]);
214  else if (fe_degree == 9)
215  do_prolongate_add<9>(to_level,
216  this->ghosted_level_vector[to_level],
217  this->ghosted_level_vector[to_level - 1]);
218  else if (fe_degree == 10)
219  do_prolongate_add<10>(to_level,
220  this->ghosted_level_vector[to_level],
221  this->ghosted_level_vector[to_level - 1]);
222  else
223  do_prolongate_add<-1>(to_level,
224  this->ghosted_level_vector[to_level],
225  this->ghosted_level_vector[to_level - 1]);
226 
227  this->ghosted_level_vector[to_level].compress(VectorOperation::add);
228  dst.copy_locally_owned_data_from(this->ghosted_level_vector[to_level]);
229 }
230 
231 
232 
233 template <int dim, typename Number>
234 void
236  const unsigned int from_level,
239 {
240  Assert((from_level >= 1) && (from_level <= level_dof_indices.size()),
241  ExcIndexRange(from_level, 1, level_dof_indices.size() + 1));
242 
243  AssertDimension(this->ghosted_level_vector[from_level].local_size(),
244  src.local_size());
245  AssertDimension(this->ghosted_level_vector[from_level - 1].local_size(),
246  dst.local_size());
247 
248  this->ghosted_level_vector[from_level].copy_locally_owned_data_from(src);
249  this->ghosted_level_vector[from_level].update_ghost_values();
250  this->ghosted_level_vector[from_level - 1] = 0.;
251 
252  if (fe_degree == 0)
253  do_restrict_add<0>(from_level,
254  this->ghosted_level_vector[from_level - 1],
255  this->ghosted_level_vector[from_level]);
256  else if (fe_degree == 1)
257  do_restrict_add<1>(from_level,
258  this->ghosted_level_vector[from_level - 1],
259  this->ghosted_level_vector[from_level]);
260  else if (fe_degree == 2)
261  do_restrict_add<2>(from_level,
262  this->ghosted_level_vector[from_level - 1],
263  this->ghosted_level_vector[from_level]);
264  else if (fe_degree == 3)
265  do_restrict_add<3>(from_level,
266  this->ghosted_level_vector[from_level - 1],
267  this->ghosted_level_vector[from_level]);
268  else if (fe_degree == 4)
269  do_restrict_add<4>(from_level,
270  this->ghosted_level_vector[from_level - 1],
271  this->ghosted_level_vector[from_level]);
272  else if (fe_degree == 5)
273  do_restrict_add<5>(from_level,
274  this->ghosted_level_vector[from_level - 1],
275  this->ghosted_level_vector[from_level]);
276  else if (fe_degree == 6)
277  do_restrict_add<6>(from_level,
278  this->ghosted_level_vector[from_level - 1],
279  this->ghosted_level_vector[from_level]);
280  else if (fe_degree == 7)
281  do_restrict_add<7>(from_level,
282  this->ghosted_level_vector[from_level - 1],
283  this->ghosted_level_vector[from_level]);
284  else if (fe_degree == 8)
285  do_restrict_add<8>(from_level,
286  this->ghosted_level_vector[from_level - 1],
287  this->ghosted_level_vector[from_level]);
288  else if (fe_degree == 9)
289  do_restrict_add<9>(from_level,
290  this->ghosted_level_vector[from_level - 1],
291  this->ghosted_level_vector[from_level]);
292  else if (fe_degree == 10)
293  do_restrict_add<10>(from_level,
294  this->ghosted_level_vector[from_level - 1],
295  this->ghosted_level_vector[from_level]);
296  else
297  // go to the non-templated version of the evaluator
298  do_restrict_add<-1>(from_level,
299  this->ghosted_level_vector[from_level - 1],
300  this->ghosted_level_vector[from_level]);
301 
302  this->ghosted_level_vector[from_level - 1].compress(VectorOperation::add);
303  dst += this->ghosted_level_vector[from_level - 1];
304 }
305 
306 
307 
308 namespace
309 {
310  template <int dim, int degree, typename Number>
311  void
312  weight_dofs_on_child(const VectorizedArray<Number> *weights,
313  const unsigned int n_components,
314  const unsigned int fe_degree,
316  {
317  Assert(fe_degree > 0, ExcNotImplemented());
318  Assert(fe_degree < 100, ExcNotImplemented());
319  const int loop_length = degree != -1 ? 2 * degree + 1 : 2 * fe_degree + 1;
320  unsigned int degree_to_3[100];
321  degree_to_3[0] = 0;
322  for (int i = 1; i < loop_length - 1; ++i)
323  degree_to_3[i] = 1;
324  degree_to_3[loop_length - 1] = 2;
325  for (unsigned int c = 0; c < n_components; ++c)
326  for (int k = 0; k < (dim > 2 ? loop_length : 1); ++k)
327  for (int j = 0; j < (dim > 1 ? loop_length : 1); ++j)
328  {
329  const unsigned int shift = 9 * degree_to_3[k] + 3 * degree_to_3[j];
330  data[0] *= weights[shift];
331  // loop bound as int avoids compiler warnings in case loop_length
332  // == 1 (polynomial degree 0)
333  for (int i = 1; i < loop_length - 1; ++i)
334  data[i] *= weights[shift + 1];
335  data[loop_length - 1] *= weights[shift + 2];
336  data += loop_length;
337  }
338  }
339 } // namespace
340 
341 
342 
343 template <int dim, typename Number>
344 template <int degree>
345 void
347  const unsigned int to_level,
350 {
351  const unsigned int vec_size = VectorizedArray<Number>::n_array_elements;
352  const unsigned int degree_size = (degree > -1 ? degree : fe_degree) + 1;
353  const unsigned int n_child_dofs_1d = 2 * degree_size - element_is_continuous;
354  const unsigned int n_scalar_cell_dofs =
355  Utilities::fixed_power<dim>(n_child_dofs_1d);
356  constexpr unsigned int three_to_dim = Utilities::pow(3, dim);
357 
358  for (unsigned int cell = 0; cell < n_owned_level_cells[to_level - 1];
359  cell += vec_size)
360  {
361  const unsigned int n_lanes =
362  cell + vec_size > n_owned_level_cells[to_level - 1] ?
363  n_owned_level_cells[to_level - 1] - cell :
364  vec_size;
365 
366  // read from source vector
367  for (unsigned int v = 0; v < n_lanes; ++v)
368  {
369  const unsigned int shift =
370  internal::MGTransfer::compute_shift_within_children<dim>(
371  parent_child_connect[to_level - 1][cell + v].second,
372  fe_degree + 1 - element_is_continuous,
373  fe_degree);
374  const unsigned int *indices =
375  &level_dof_indices[to_level - 1]
376  [parent_child_connect[to_level - 1][cell + v]
377  .first *
378  n_child_cell_dofs +
379  shift];
380  for (unsigned int c = 0, m = 0; c < n_components; ++c)
381  {
382  for (unsigned int k = 0; k < (dim > 2 ? degree_size : 1); ++k)
383  for (unsigned int j = 0; j < (dim > 1 ? degree_size : 1); ++j)
384  for (unsigned int i = 0; i < degree_size; ++i, ++m)
385  evaluation_data[m][v] = src.local_element(
386  indices[c * n_scalar_cell_dofs +
387  k * n_child_dofs_1d * n_child_dofs_1d +
388  j * n_child_dofs_1d + i]);
389 
390  // apply Dirichlet boundary conditions on parent cell
391  for (std::vector<unsigned short>::const_iterator i =
392  dirichlet_indices[to_level - 1][cell + v].begin();
393  i != dirichlet_indices[to_level - 1][cell + v].end();
394  ++i)
395  evaluation_data[*i][v] = 0.;
396  }
397  }
398 
399  AssertDimension(prolongation_matrix_1d.size(),
400  degree_size * n_child_dofs_1d);
401  // perform tensorized operation
402  if (element_is_continuous)
403  {
404  // must go through the components backwards because we want to write
405  // the output to the same array as the input
406  for (int c = n_components - 1; c >= 0; --c)
408  dim,
409  degree + 1,
410  2 * degree + 1,
411  1,
414  do_forward(prolongation_matrix_1d,
415  evaluation_data.begin() +
416  c * Utilities::fixed_power<dim>(degree_size),
417  evaluation_data.begin() + c * n_scalar_cell_dofs,
418  fe_degree + 1,
419  2 * fe_degree + 1);
420  weight_dofs_on_child<dim, degree, Number>(
421  &weights_on_refined[to_level - 1][(cell / vec_size) * three_to_dim],
422  n_components,
423  fe_degree,
424  evaluation_data.begin());
425  }
426  else
427  {
428  for (int c = n_components - 1; c >= 0; --c)
430  dim,
431  degree + 1,
432  2 * degree + 2,
433  1,
436  do_forward(prolongation_matrix_1d,
437  evaluation_data.begin() +
438  c * Utilities::fixed_power<dim>(degree_size),
439  evaluation_data.begin() + c * n_scalar_cell_dofs,
440  fe_degree + 1,
441  2 * fe_degree + 2);
442  }
443 
444  // write into dst vector
445  const unsigned int *indices =
446  &level_dof_indices[to_level][cell * n_child_cell_dofs];
447  for (unsigned int v = 0; v < n_lanes; ++v)
448  {
449  for (unsigned int i = 0; i < n_child_cell_dofs; ++i)
450  dst.local_element(indices[i]) += evaluation_data[i][v];
451  indices += n_child_cell_dofs;
452  }
453  }
454 }
455 
456 
457 
458 template <int dim, typename Number>
459 template <int degree>
460 void
462  const unsigned int from_level,
465 {
466  const unsigned int vec_size = VectorizedArray<Number>::n_array_elements;
467  const unsigned int degree_size = (degree > -1 ? degree : fe_degree) + 1;
468  const unsigned int n_child_dofs_1d = 2 * degree_size - element_is_continuous;
469  const unsigned int n_scalar_cell_dofs =
470  Utilities::fixed_power<dim>(n_child_dofs_1d);
471  constexpr unsigned int three_to_dim = Utilities::pow(3, dim);
472 
473  for (unsigned int cell = 0; cell < n_owned_level_cells[from_level - 1];
474  cell += vec_size)
475  {
476  const unsigned int n_lanes =
477  cell + vec_size > n_owned_level_cells[from_level - 1] ?
478  n_owned_level_cells[from_level - 1] - cell :
479  vec_size;
480 
481  // read from source vector
482  {
483  const unsigned int *indices =
484  &level_dof_indices[from_level][cell * n_child_cell_dofs];
485  for (unsigned int v = 0; v < n_lanes; ++v)
486  {
487  for (unsigned int i = 0; i < n_child_cell_dofs; ++i)
488  evaluation_data[i][v] = src.local_element(indices[i]);
489  indices += n_child_cell_dofs;
490  }
491  }
492 
493  AssertDimension(prolongation_matrix_1d.size(),
494  degree_size * n_child_dofs_1d);
495  // perform tensorized operation
496  if (element_is_continuous)
497  {
498  weight_dofs_on_child<dim, degree, Number>(
499  &weights_on_refined[from_level - 1]
500  [(cell / vec_size) * three_to_dim],
501  n_components,
502  fe_degree,
503  evaluation_data.data());
504  for (unsigned int c = 0; c < n_components; ++c)
506  dim,
507  degree + 1,
508  2 * degree + 1,
509  1,
512  do_backward(prolongation_matrix_1d,
513  false,
514  evaluation_data.begin() + c * n_scalar_cell_dofs,
515  evaluation_data.begin() +
516  c * Utilities::fixed_power<dim>(degree_size),
517  fe_degree + 1,
518  2 * fe_degree + 1);
519  }
520  else
521  {
522  for (unsigned int c = 0; c < n_components; ++c)
524  dim,
525  degree + 1,
526  2 * degree + 2,
527  1,
530  do_backward(prolongation_matrix_1d,
531  false,
532  evaluation_data.begin() + c * n_scalar_cell_dofs,
533  evaluation_data.begin() +
534  c * Utilities::fixed_power<dim>(degree_size),
535  fe_degree + 1,
536  2 * fe_degree + 2);
537  }
538 
539  // write into dst vector
540  for (unsigned int v = 0; v < n_lanes; ++v)
541  {
542  const unsigned int shift =
543  internal::MGTransfer::compute_shift_within_children<dim>(
544  parent_child_connect[from_level - 1][cell + v].second,
545  fe_degree + 1 - element_is_continuous,
546  fe_degree);
548  parent_child_connect[from_level - 1][cell + v].first *
549  n_child_cell_dofs +
550  n_child_cell_dofs - 1,
551  level_dof_indices[from_level - 1].size());
552  const unsigned int *indices =
553  &level_dof_indices[from_level - 1]
554  [parent_child_connect[from_level - 1][cell + v]
555  .first *
556  n_child_cell_dofs +
557  shift];
558  for (unsigned int c = 0, m = 0; c < n_components; ++c)
559  {
560  // apply Dirichlet boundary conditions on parent cell
561  for (std::vector<unsigned short>::const_iterator i =
562  dirichlet_indices[from_level - 1][cell + v].begin();
563  i != dirichlet_indices[from_level - 1][cell + v].end();
564  ++i)
565  evaluation_data[*i][v] = 0.;
566 
567  for (unsigned int k = 0; k < (dim > 2 ? degree_size : 1); ++k)
568  for (unsigned int j = 0; j < (dim > 1 ? degree_size : 1); ++j)
569  for (unsigned int i = 0; i < degree_size; ++i, ++m)
570  dst.local_element(
571  indices[c * n_scalar_cell_dofs +
572  k * n_child_dofs_1d * n_child_dofs_1d +
573  j * n_child_dofs_1d + i]) +=
574  evaluation_data[m][v];
575  }
576  }
577  }
578 }
579 
580 
581 
582 template <int dim, typename Number>
583 std::size_t
585 {
586  std::size_t memory = MGLevelGlobalTransfer<
587  LinearAlgebra::distributed::Vector<Number>>::memory_consumption();
588  memory += MemoryConsumption::memory_consumption(level_dof_indices);
589  memory += MemoryConsumption::memory_consumption(parent_child_connect);
590  memory += MemoryConsumption::memory_consumption(n_owned_level_cells);
591  memory += MemoryConsumption::memory_consumption(prolongation_matrix_1d);
592  memory += MemoryConsumption::memory_consumption(evaluation_data);
593  memory += MemoryConsumption::memory_consumption(weights_on_refined);
594  memory += MemoryConsumption::memory_consumption(dirichlet_indices);
595  return memory;
596 }
597 
598 
599 
600 template <int dim, typename Number>
602  const MGConstrainedDoFs &mg_c)
603  : same_for_all(true)
604 {
605  matrix_free_transfer_vector.emplace_back(mg_c);
606 }
607 
608 
609 
610 template <int dim, typename Number>
612  const std::vector<MGConstrainedDoFs> &mg_c)
613  : same_for_all(false)
614 {
615  for (const auto &constrained_block_dofs : mg_c)
616  matrix_free_transfer_vector.emplace_back(constrained_block_dofs);
617 }
618 
619 
620 
621 template <int dim, typename Number>
622 void
624  const MGConstrainedDoFs &mg_c)
625 {
626  Assert(same_for_all,
627  ExcMessage("This object was initialized with support for usage with "
628  "one DoFHandler for each block, but this method assumes "
629  "that the same DoFHandler is used for all the blocks!"));
630  AssertDimension(matrix_free_transfer_vector.size(), 1);
631 
632  matrix_free_transfer_vector[0].initialize_constraints(mg_c);
633 }
634 
635 
636 
637 template <int dim, typename Number>
638 void
640  const std::vector<MGConstrainedDoFs> &mg_c)
641 {
642  Assert(!same_for_all,
643  ExcMessage("This object was initialized with support for using "
644  "the same DoFHandler for all the blocks, but this "
645  "method assumes that there is a separate DoFHandler "
646  "for each block!"));
647  AssertDimension(matrix_free_transfer_vector.size(), mg_c.size());
648 
649  for (unsigned int i = 0; i < mg_c.size(); ++i)
650  matrix_free_transfer_vector[i].initialize_constraints(mg_c[i]);
651 }
652 
653 
654 
655 template <int dim, typename Number>
656 void
658 {
659  matrix_free_transfer_vector.clear();
660 }
661 
662 
663 
664 template <int dim, typename Number>
665 void
667  const DoFHandler<dim, dim> &mg_dof)
668 {
669  AssertDimension(matrix_free_transfer_vector.size(), 1);
670  matrix_free_transfer_vector[0].build(mg_dof);
671 }
672 
673 
674 
675 template <int dim, typename Number>
676 void
678  const std::vector<const DoFHandler<dim, dim> *> &mg_dof)
679 {
680  AssertDimension(matrix_free_transfer_vector.size(), mg_dof.size());
681  for (unsigned int i = 0; i < mg_dof.size(); ++i)
682  matrix_free_transfer_vector[i].build(*mg_dof[i]);
683 }
684 
685 
686 
687 template <int dim, typename Number>
688 void
690  const unsigned int to_level,
693 {
694  const unsigned int n_blocks = src.n_blocks();
695  AssertDimension(dst.n_blocks(), n_blocks);
696 
697  if (!same_for_all)
698  AssertDimension(matrix_free_transfer_vector.size(), n_blocks);
699 
700  for (unsigned int b = 0; b < n_blocks; ++b)
701  {
702  const unsigned int data_block = same_for_all ? 0 : b;
703  matrix_free_transfer_vector[data_block].prolongate(to_level,
704  dst.block(b),
705  src.block(b));
706  }
707 }
708 
709 
710 
711 template <int dim, typename Number>
712 void
714  const unsigned int from_level,
717 {
718  const unsigned int n_blocks = src.n_blocks();
719  AssertDimension(dst.n_blocks(), n_blocks);
720 
721  if (!same_for_all)
722  AssertDimension(matrix_free_transfer_vector.size(), n_blocks);
723 
724  for (unsigned int b = 0; b < n_blocks; ++b)
725  {
726  const unsigned int data_block = same_for_all ? 0 : b;
727  matrix_free_transfer_vector[data_block].restrict_and_add(from_level,
728  dst.block(b),
729  src.block(b));
730  }
731 }
732 
733 
734 
735 template <int dim, typename Number>
736 std::size_t
738 {
739  std::size_t total_memory_consumption = 0;
740  for (const auto &el : matrix_free_transfer_vector)
741  total_memory_consumption += el.memory_consumption();
742  return total_memory_consumption;
743 }
744 
745 
746 
747 // explicit instantiation
748 #include "mg_transfer_matrix_free.inst"
749 
750 
751 DEAL_II_NAMESPACE_CLOSE
void copy_locally_owned_data_from(const Vector< Number2, MemorySpace > &src)
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:1567
std::size_t memory_consumption() const
void build(const DoFHandler< dim, dim > &mg_dof)
void initialize_constraints(const MGConstrainedDoFs &mg_constrained_dofs)
#define AssertIndexRange(index, range)
Definition: exceptions.h:1637
static ::ExceptionBase & ExcIndexRange(int arg1, int arg2, int arg3)
std::size_t memory_consumption() const
Number local_element(const size_type local_index) const
static ::ExceptionBase & ExcMessage(std::string arg1)
#define Assert(cond, exc)
Definition: exceptions.h:1407
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
constexpr unsigned int pow(const unsigned int base, const int iexp)
Definition: utilities.h:428
virtual void prolongate(const unsigned int to_level, LinearAlgebra::distributed::Vector< Number > &dst, const LinearAlgebra::distributed::Vector< Number > &src) const override
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:616
virtual void prolongate(const unsigned int to_level, LinearAlgebra::distributed::BlockVector< Number > &dst, const LinearAlgebra::distributed::BlockVector< Number > &src) const override
std::vector< MGTransferMatrixFree< dim, Number > > matrix_free_transfer_vector
MGTransferBlockMatrixFree()=default
static ::ExceptionBase & ExcNotImplemented()
unsigned int n_components(const DoFHandler< dim, spacedim > &dh)
void initialize_constraints(const MGConstrainedDoFs &mg_constrained_dofs)
BlockType & block(const unsigned int i)
void build(const DoFHandler< dim, dim > &mg_dof)
std::enable_if< std::is_fundamental< T >::value, std::size_t >::type memory_consumption(const T &t)