Reference documentation for deal.II version GIT c545eda070 2023-01-27 00:25:02+00:00
\(\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 - 2022 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 
33 
36 #include <deal.II/multigrid/mg_transfer_matrix_free.templates.h>
37 
38 #include <algorithm>
39 
41 
42 
43 template <int dim, typename Number>
45  : fe_degree(0)
46  , element_is_continuous(false)
47  , n_components(0)
48  , n_child_cell_dofs(0)
49 {}
50 
51 
52 
53 template <int dim, typename Number>
55  const MGConstrainedDoFs &mg_c)
56  : fe_degree(0)
57  , element_is_continuous(false)
58  , n_components(0)
59  , n_child_cell_dofs(0)
60 {
61  this->mg_constrained_dofs = &mg_c;
62 }
63 
64 
65 
66 template <int dim, typename Number>
67 void
69  const MGConstrainedDoFs &mg_c)
70 {
71  this->mg_constrained_dofs = &mg_c;
72 }
73 
74 
75 
76 template <int dim, typename Number>
77 void
79 {
82  fe_degree = 0;
83  element_is_continuous = false;
84  n_components = 0;
85  n_child_cell_dofs = 0;
86  level_dof_indices.clear();
87  parent_child_connect.clear();
88  dirichlet_indices.clear();
89  n_owned_level_cells.clear();
90  prolongation_matrix_1d.clear();
91  evaluation_data.clear();
92  weights_on_refined.clear();
93 }
94 
95 
96 
97 template <int dim, typename Number>
98 void
100  const DoFHandler<dim, dim> &dof_handler,
101  const std::vector<std::shared_ptr<const Utilities::MPI::Partitioner>>
102  &external_partitioners)
103 {
104  Assert(dof_handler.has_level_dofs(),
105  ExcMessage(
106  "The underlying DoFHandler object has not had its "
107  "distribute_mg_dofs() function called, but this is a prerequisite "
108  "for multigrid transfers. You will need to call this function, "
109  "probably close to where you already call distribute_dofs()."));
110  if (external_partitioners.size() > 0)
111  {
112  Assert(
113  this->initialize_dof_vector == nullptr,
114  ExcMessage(
115  "A initialize_dof_vector function has already been registered in the constructor!"));
116 
117  this->initialize_dof_vector =
118  [external_partitioners](
119  const unsigned int level,
121  vec.reinit(external_partitioners[level]);
122  };
123  }
124 
125  this->fill_and_communicate_copy_indices(dof_handler);
126 
127  vector_partitioners.resize(0,
128  dof_handler.get_triangulation().n_global_levels() -
129  1);
130  for (unsigned int level = 0; level <= this->ghosted_level_vector.max_level();
131  ++level)
132  vector_partitioners[level] =
133  this->ghosted_level_vector[level].get_partitioner();
134 
135  std::vector<std::vector<Number>> weights_unvectorized;
136 
138 
139  internal::MGTransfer::setup_transfer<dim, Number>(
140  dof_handler,
141  this->mg_constrained_dofs,
142  external_partitioners,
143  elem_info,
144  level_dof_indices,
145  parent_child_connect,
146  n_owned_level_cells,
147  dirichlet_indices,
148  weights_unvectorized,
149  this->copy_indices_global_mine,
150  vector_partitioners);
151 
152  // reinit the ghosted level vector in case its partitioner differs from the
153  // one the user intends to use externally
154  this->ghosted_level_vector.resize(0, vector_partitioners.max_level());
155  for (unsigned int level = 0; level <= vector_partitioners.max_level();
156  ++level)
157  if (external_partitioners.size() == vector_partitioners.max_level() + 1 &&
158  external_partitioners[level].get() == vector_partitioners[level].get())
159  this->ghosted_level_vector[level].reinit(0);
160  else
161  this->ghosted_level_vector[level].reinit(vector_partitioners[level]);
162 
163  // unpack element info data
164  fe_degree = elem_info.fe_degree;
165  element_is_continuous = elem_info.element_is_continuous;
166  n_components = elem_info.n_components;
167  n_child_cell_dofs = elem_info.n_child_cell_dofs;
168 
169  // duplicate and put into vectorized array
170  prolongation_matrix_1d.resize(elem_info.prolongation_matrix_1d.size());
171  for (unsigned int i = 0; i < elem_info.prolongation_matrix_1d.size(); ++i)
172  prolongation_matrix_1d[i] = elem_info.prolongation_matrix_1d[i];
173 
174  // reshuffle into aligned vector of vectorized arrays
175  const unsigned int vec_size = VectorizedArray<Number>::size();
176  const unsigned int n_levels =
177  dof_handler.get_triangulation().n_global_levels();
178 
179  const unsigned int n_weights_per_cell = Utilities::fixed_power<dim>(3);
180  weights_on_refined.resize(n_levels - 1);
181  for (unsigned int level = 1; level < n_levels; ++level)
182  {
183  weights_on_refined[level - 1].resize(
184  ((n_owned_level_cells[level - 1] + vec_size - 1) / vec_size) *
185  n_weights_per_cell);
186 
187  for (unsigned int c = 0; c < n_owned_level_cells[level - 1]; ++c)
188  {
189  const unsigned int comp = c / vec_size;
190  const unsigned int v = c % vec_size;
191  for (unsigned int i = 0; i < n_weights_per_cell; ++i)
192  {
193  weights_on_refined[level - 1][comp * n_weights_per_cell + i][v] =
194  weights_unvectorized[level - 1][c * n_weights_per_cell + i];
195  }
196  }
197  }
198 
199  evaluation_data.resize(n_child_cell_dofs);
200 }
201 
202 
203 
204 template <int dim, typename Number>
205 void
207  const DoFHandler<dim, dim> &dof_handler,
208  const std::function<void(const unsigned int,
210  &initialize_dof_vector)
211 {
212  if (initialize_dof_vector)
213  {
214  const unsigned int n_levels =
215  dof_handler.get_triangulation().n_global_levels();
216 
217  std::vector<std::shared_ptr<const Utilities::MPI::Partitioner>>
218  external_partitioners(n_levels);
219 
220  for (unsigned int level = 0; level < n_levels; ++level)
221  {
223  initialize_dof_vector(level, vector);
224  external_partitioners[level] = vector.get_partitioner();
225  }
226 
227  build(dof_handler, external_partitioners);
228  }
229  else
230  {
231  build(dof_handler);
232  }
233 }
234 
235 
236 
237 template <int dim, typename Number>
238 void
240  const unsigned int to_level,
243 {
244  dst = 0;
245  prolongate_and_add(to_level, dst, src);
246 }
247 
248 
249 
250 template <int dim, typename Number>
251 void
253  const unsigned int to_level,
256 {
257  Assert((to_level >= 1) && (to_level <= level_dof_indices.size()),
258  ExcIndexRange(to_level, 1, level_dof_indices.size() + 1));
259 
260  const bool src_inplace = src.get_partitioner().get() ==
261  this->vector_partitioners[to_level - 1].get();
262  if (src_inplace == false)
263  {
264  if (this->ghosted_level_vector[to_level - 1].get_partitioner().get() !=
265  this->vector_partitioners[to_level - 1].get())
266  this->ghosted_level_vector[to_level - 1].reinit(
267  this->vector_partitioners[to_level - 1]);
268  this->ghosted_level_vector[to_level - 1].copy_locally_owned_data_from(
269  src);
270  }
271 
272  const bool dst_inplace =
273  dst.get_partitioner().get() == this->vector_partitioners[to_level].get();
274  if (dst_inplace == false)
275  {
276  if (this->ghosted_level_vector[to_level].get_partitioner().get() !=
277  this->vector_partitioners[to_level].get())
278  this->ghosted_level_vector[to_level].reinit(
279  this->vector_partitioners[to_level]);
280  AssertDimension(this->ghosted_level_vector[to_level].locally_owned_size(),
281  dst.locally_owned_size());
282  this->ghosted_level_vector[to_level] = 0.;
283  }
284 
286  src_inplace ? src : this->ghosted_level_vector[to_level - 1];
288  dst_inplace ? dst : this->ghosted_level_vector[to_level];
289 
290  src_vec.update_ghost_values();
291  // the implementation in do_prolongate_add is templated in the degree of the
292  // element (for efficiency reasons), so we need to find the appropriate
293  // kernel here...
294  if (fe_degree == 0)
295  do_prolongate_add<0>(to_level, dst_vec, src_vec);
296  else if (fe_degree == 1)
297  do_prolongate_add<1>(to_level, dst_vec, src_vec);
298  else if (fe_degree == 2)
299  do_prolongate_add<2>(to_level, dst_vec, src_vec);
300  else if (fe_degree == 3)
301  do_prolongate_add<3>(to_level, dst_vec, src_vec);
302  else if (fe_degree == 4)
303  do_prolongate_add<4>(to_level, dst_vec, src_vec);
304  else if (fe_degree == 5)
305  do_prolongate_add<5>(to_level, dst_vec, src_vec);
306  else if (fe_degree == 6)
307  do_prolongate_add<6>(to_level, dst_vec, src_vec);
308  else if (fe_degree == 7)
309  do_prolongate_add<7>(to_level, dst_vec, src_vec);
310  else if (fe_degree == 8)
311  do_prolongate_add<8>(to_level, dst_vec, src_vec);
312  else if (fe_degree == 9)
313  do_prolongate_add<9>(to_level, dst_vec, src_vec);
314  else if (fe_degree == 10)
315  do_prolongate_add<10>(to_level, dst_vec, src_vec);
316  else
317  do_prolongate_add<-1>(to_level, dst_vec, src_vec);
318 
320  if (dst_inplace == false)
321  dst += dst_vec;
322 
323  if (src_inplace == true)
324  src.zero_out_ghost_values();
325 }
326 
327 
328 
329 template <int dim, typename Number>
330 void
332  const unsigned int from_level,
335 {
336  Assert((from_level >= 1) && (from_level <= level_dof_indices.size()),
337  ExcIndexRange(from_level, 1, level_dof_indices.size() + 1));
338 
339  const bool src_inplace =
340  src.get_partitioner().get() == this->vector_partitioners[from_level].get();
341  if (src_inplace == false)
342  {
343  if (this->ghosted_level_vector[from_level].get_partitioner().get() !=
344  this->vector_partitioners[from_level].get())
345  this->ghosted_level_vector[from_level].reinit(
346  this->vector_partitioners[from_level]);
347  this->ghosted_level_vector[from_level].copy_locally_owned_data_from(src);
348  }
349 
350  const bool dst_inplace = dst.get_partitioner().get() ==
351  this->vector_partitioners[from_level - 1].get();
352  if (dst_inplace == false)
353  {
354  if (this->ghosted_level_vector[from_level - 1].get_partitioner().get() !=
355  this->vector_partitioners[from_level - 1].get())
356  this->ghosted_level_vector[from_level - 1].reinit(
357  this->vector_partitioners[from_level - 1]);
359  this->ghosted_level_vector[from_level - 1].locally_owned_size(),
360  dst.locally_owned_size());
361  this->ghosted_level_vector[from_level - 1] = 0.;
362  }
363 
365  src_inplace ? src : this->ghosted_level_vector[from_level];
367  dst_inplace ? dst : this->ghosted_level_vector[from_level - 1];
368 
369  src_vec.update_ghost_values();
370 
371  if (fe_degree == 0)
372  do_restrict_add<0>(from_level, dst_vec, src_vec);
373  else if (fe_degree == 1)
374  do_restrict_add<1>(from_level, dst_vec, src_vec);
375  else if (fe_degree == 2)
376  do_restrict_add<2>(from_level, dst_vec, src_vec);
377  else if (fe_degree == 3)
378  do_restrict_add<3>(from_level, dst_vec, src_vec);
379  else if (fe_degree == 4)
380  do_restrict_add<4>(from_level, dst_vec, src_vec);
381  else if (fe_degree == 5)
382  do_restrict_add<5>(from_level, dst_vec, src_vec);
383  else if (fe_degree == 6)
384  do_restrict_add<6>(from_level, dst_vec, src_vec);
385  else if (fe_degree == 7)
386  do_restrict_add<7>(from_level, dst_vec, src_vec);
387  else if (fe_degree == 8)
388  do_restrict_add<8>(from_level, dst_vec, src_vec);
389  else if (fe_degree == 9)
390  do_restrict_add<9>(from_level, dst_vec, src_vec);
391  else if (fe_degree == 10)
392  do_restrict_add<10>(from_level, dst_vec, src_vec);
393  else
394  // go to the non-templated version of the evaluator
395  do_restrict_add<-1>(from_level, dst_vec, src_vec);
396 
398  if (dst_inplace == false)
399  dst += dst_vec;
400 
401  if (src_inplace == true)
402  src.zero_out_ghost_values();
403 }
404 
405 
406 
407 template <int dim, typename Number>
408 template <int degree>
409 void
411  const unsigned int to_level,
414 {
415  const unsigned int vec_size = VectorizedArray<Number>::size();
416  const unsigned int degree_size = (degree > -1 ? degree : fe_degree) + 1;
417  const unsigned int n_child_dofs_1d = 2 * degree_size - element_is_continuous;
418  const unsigned int n_scalar_cell_dofs =
419  Utilities::fixed_power<dim>(n_child_dofs_1d);
420  constexpr unsigned int three_to_dim = Utilities::pow(3, dim);
421 
422  // If we have user defined MG constraints, we must create
423  // a non-const, ghosted version of the source vector to distribute
424  // constraints.
425  const LinearAlgebra::distributed::Vector<Number> *to_use = &src;
427  if (this->mg_constrained_dofs != nullptr &&
428  this->mg_constrained_dofs->get_user_constraint_matrix(to_level - 1)
429  .get_local_lines()
430  .size() > 0)
431  {
433 
434  // Distribute any user defined constraints
435  this->mg_constrained_dofs->get_user_constraint_matrix(to_level - 1)
436  .distribute(copy_src);
437 
438  // Re-initialize new ghosted vector with correct constraints
439  new_src.reinit(copy_src);
440  new_src = copy_src;
441  new_src.update_ghost_values();
442 
443  to_use = &new_src;
444  }
445 
446  for (unsigned int cell = 0; cell < n_owned_level_cells[to_level - 1];
447  cell += vec_size)
448  {
449  const unsigned int n_lanes =
450  (cell + vec_size > n_owned_level_cells[to_level - 1]) ?
451  (n_owned_level_cells[to_level - 1] - cell) :
452  vec_size;
453 
454  // read from source vector
455  for (unsigned int v = 0; v < n_lanes; ++v)
456  {
457  const unsigned int shift =
458  internal::MGTransfer::compute_shift_within_children<dim>(
459  parent_child_connect[to_level - 1][cell + v].second,
460  fe_degree + 1 - element_is_continuous,
461  fe_degree);
462  const unsigned int *indices =
463  &level_dof_indices[to_level - 1]
464  [parent_child_connect[to_level - 1][cell + v]
465  .first *
466  n_child_cell_dofs +
467  shift];
468  for (unsigned int c = 0, m = 0; c < n_components; ++c)
469  {
470  for (unsigned int k = 0; k < (dim > 2 ? degree_size : 1); ++k)
471  for (unsigned int j = 0; j < (dim > 1 ? degree_size : 1); ++j)
472  for (unsigned int i = 0; i < degree_size; ++i, ++m)
473  evaluation_data[m][v] = to_use->local_element(
474  indices[c * n_scalar_cell_dofs +
475  k * n_child_dofs_1d * n_child_dofs_1d +
476  j * n_child_dofs_1d + i]);
477 
478  // apply Dirichlet boundary conditions on parent cell
479  for (std::vector<unsigned short>::const_iterator i =
480  dirichlet_indices[to_level - 1][cell + v].begin();
481  i != dirichlet_indices[to_level - 1][cell + v].end();
482  ++i)
483  evaluation_data[*i][v] = 0.;
484  }
485  }
486 
487  AssertDimension(prolongation_matrix_1d.size(),
488  degree_size * n_child_dofs_1d);
489  // perform tensorized operation
490  if (element_is_continuous)
491  {
492  // must go through the components backwards because we want to write
493  // the output to the same array as the input
494  for (int c = n_components - 1; c >= 0; --c)
498  dim,
499  degree + 1,
500  2 * degree + 1,
502  VectorizedArray<Number>>::do_forward(1,
503  prolongation_matrix_1d,
504  evaluation_data.begin() +
506  dim>(degree_size),
507  evaluation_data.begin() +
508  c * n_scalar_cell_dofs,
509  fe_degree + 1,
510  2 * fe_degree + 1);
512  degree != -1 ? 2 * degree + 1 :
513  -1,
515  &weights_on_refined[to_level - 1][(cell / vec_size) * three_to_dim],
516  n_components,
517  2 * fe_degree + 1,
518  evaluation_data.begin());
519  }
520  else
521  {
522  for (int c = n_components - 1; c >= 0; --c)
526  dim,
527  degree + 1,
528  2 * degree + 2,
530  VectorizedArray<Number>>::do_forward(1,
531  prolongation_matrix_1d,
532  evaluation_data.begin() +
534  dim>(degree_size),
535  evaluation_data.begin() +
536  c * n_scalar_cell_dofs,
537  fe_degree + 1,
538  2 * fe_degree + 2);
539  }
540 
541  // write into dst vector
542  const unsigned int *indices =
543  &level_dof_indices[to_level][cell * n_child_cell_dofs];
544  for (unsigned int v = 0; v < n_lanes; ++v)
545  {
546  for (unsigned int i = 0; i < n_child_cell_dofs; ++i)
547  dst.local_element(indices[i]) += evaluation_data[i][v];
548  indices += n_child_cell_dofs;
549  }
550  }
551 }
552 
553 
554 
555 template <int dim, typename Number>
556 template <int degree>
557 void
559  const unsigned int from_level,
562 {
563  const unsigned int vec_size = VectorizedArray<Number>::size();
564  const unsigned int degree_size = (degree > -1 ? degree : fe_degree) + 1;
565  const unsigned int n_child_dofs_1d = 2 * degree_size - element_is_continuous;
566  const unsigned int n_scalar_cell_dofs =
567  Utilities::fixed_power<dim>(n_child_dofs_1d);
568  constexpr unsigned int three_to_dim = Utilities::pow(3, dim);
569 
570  for (unsigned int cell = 0; cell < n_owned_level_cells[from_level - 1];
571  cell += vec_size)
572  {
573  const unsigned int n_lanes =
574  (cell + vec_size > n_owned_level_cells[from_level - 1]) ?
575  (n_owned_level_cells[from_level - 1] - cell) :
576  vec_size;
577 
578  // read from source vector
579  {
580  const unsigned int *indices =
581  &level_dof_indices[from_level][cell * n_child_cell_dofs];
582  for (unsigned int v = 0; v < n_lanes; ++v)
583  {
584  for (unsigned int i = 0; i < n_child_cell_dofs; ++i)
585  evaluation_data[i][v] = src.local_element(indices[i]);
586  indices += n_child_cell_dofs;
587  }
588  }
589 
590  AssertDimension(prolongation_matrix_1d.size(),
591  degree_size * n_child_dofs_1d);
592  // perform tensorized operation
593  if (element_is_continuous)
594  {
596  degree != -1 ? 2 * degree + 1 :
597  -1,
599  &weights_on_refined[from_level - 1]
600  [(cell / vec_size) * three_to_dim],
601  n_components,
602  2 * fe_degree + 1,
603  evaluation_data.data());
604  for (unsigned int c = 0; c < n_components; ++c)
608  dim,
609  degree + 1,
610  2 * degree + 1,
612  VectorizedArray<Number>>::do_backward(1,
613  prolongation_matrix_1d,
614  false,
615  evaluation_data.begin() +
616  c * n_scalar_cell_dofs,
617  evaluation_data.begin() +
618  c *
620  dim>(degree_size),
621  fe_degree + 1,
622  2 * fe_degree + 1);
623  }
624  else
625  {
626  for (unsigned int c = 0; c < n_components; ++c)
630  dim,
631  degree + 1,
632  2 * degree + 2,
634  VectorizedArray<Number>>::do_backward(1,
635  prolongation_matrix_1d,
636  false,
637  evaluation_data.begin() +
638  c * n_scalar_cell_dofs,
639  evaluation_data.begin() +
640  c *
642  dim>(degree_size),
643  fe_degree + 1,
644  2 * fe_degree + 2);
645  }
646 
647  // write into dst vector
648  for (unsigned int v = 0; v < n_lanes; ++v)
649  {
650  const unsigned int shift =
651  internal::MGTransfer::compute_shift_within_children<dim>(
652  parent_child_connect[from_level - 1][cell + v].second,
653  fe_degree + 1 - element_is_continuous,
654  fe_degree);
656  parent_child_connect[from_level - 1][cell + v].first *
657  n_child_cell_dofs +
658  n_child_cell_dofs - 1,
659  level_dof_indices[from_level - 1].size());
660  const unsigned int *indices =
661  &level_dof_indices[from_level - 1]
662  [parent_child_connect[from_level - 1][cell + v]
663  .first *
664  n_child_cell_dofs +
665  shift];
666  for (unsigned int c = 0, m = 0; c < n_components; ++c)
667  {
668  // apply Dirichlet boundary conditions on parent cell
669  for (std::vector<unsigned short>::const_iterator i =
670  dirichlet_indices[from_level - 1][cell + v].begin();
671  i != dirichlet_indices[from_level - 1][cell + v].end();
672  ++i)
673  evaluation_data[*i][v] = 0.;
674 
675  for (unsigned int k = 0; k < (dim > 2 ? degree_size : 1); ++k)
676  for (unsigned int j = 0; j < (dim > 1 ? degree_size : 1); ++j)
677  for (unsigned int i = 0; i < degree_size; ++i, ++m)
678  dst.local_element(
679  indices[c * n_scalar_cell_dofs +
680  k * n_child_dofs_1d * n_child_dofs_1d +
681  j * n_child_dofs_1d + i]) +=
682  evaluation_data[m][v];
683  }
684  }
685  }
686 }
687 
688 
689 
690 template <int dim, typename Number>
691 std::size_t
693 {
694  std::size_t memory = MGLevelGlobalTransfer<
696  memory += MemoryConsumption::memory_consumption(level_dof_indices);
697  memory += MemoryConsumption::memory_consumption(parent_child_connect);
698  memory += MemoryConsumption::memory_consumption(n_owned_level_cells);
699  memory += MemoryConsumption::memory_consumption(prolongation_matrix_1d);
700  memory += MemoryConsumption::memory_consumption(evaluation_data);
701  memory += MemoryConsumption::memory_consumption(weights_on_refined);
702  memory += MemoryConsumption::memory_consumption(dirichlet_indices);
703  return memory;
704 }
705 
706 
707 
708 template <int dim, typename Number>
710  const MGConstrainedDoFs &mg_c)
712  Number,
713  MGTransferMatrixFree<dim, Number>>(true)
714 {
715  this->matrix_free_transfer_vector.emplace_back(mg_c);
716 }
717 
718 
719 
720 template <int dim, typename Number>
722  const std::vector<MGConstrainedDoFs> &mg_c)
724  Number,
725  MGTransferMatrixFree<dim, Number>>(false)
726 {
727  for (const auto &constrained_block_dofs : mg_c)
728  this->matrix_free_transfer_vector.emplace_back(constrained_block_dofs);
729 }
730 
731 
732 
733 template <int dim, typename Number>
734 void
736  const MGConstrainedDoFs &mg_c)
737 {
738  Assert(this->same_for_all,
739  ExcMessage("This object was initialized with support for usage with "
740  "one DoFHandler for each block, but this method assumes "
741  "that the same DoFHandler is used for all the blocks!"));
742  AssertDimension(this->matrix_free_transfer_vector.size(), 1);
743 
744  this->matrix_free_transfer_vector[0].initialize_constraints(mg_c);
745 }
746 
747 
748 
749 template <int dim, typename Number>
750 void
752  const std::vector<MGConstrainedDoFs> &mg_c)
753 {
754  Assert(!this->same_for_all,
755  ExcMessage("This object was initialized with support for using "
756  "the same DoFHandler for all the blocks, but this "
757  "method assumes that there is a separate DoFHandler "
758  "for each block!"));
759  AssertDimension(this->matrix_free_transfer_vector.size(), mg_c.size());
760 
761  for (unsigned int i = 0; i < mg_c.size(); ++i)
762  this->matrix_free_transfer_vector[i].initialize_constraints(mg_c[i]);
763 }
764 
765 
766 
767 template <int dim, typename Number>
768 void
770  const DoFHandler<dim, dim> &dof_handler)
771 {
772  AssertDimension(this->matrix_free_transfer_vector.size(), 1);
773  this->matrix_free_transfer_vector[0].build(dof_handler);
774 }
775 
776 
777 
778 template <int dim, typename Number>
779 void
781  const std::vector<const DoFHandler<dim, dim> *> &dof_handler)
782 {
783  AssertDimension(this->matrix_free_transfer_vector.size(), dof_handler.size());
784  for (unsigned int i = 0; i < dof_handler.size(); ++i)
785  this->matrix_free_transfer_vector[i].build(*dof_handler[i]);
786 }
787 
788 
789 
790 template <int dim, typename Number>
791 std::size_t
793 {
794  std::size_t total_memory_consumption = 0;
795  for (const auto &el : this->matrix_free_transfer_vector)
796  total_memory_consumption += el.memory_consumption();
797  return total_memory_consumption;
798 }
799 
800 
801 
802 template <int dim, typename Number>
805  const unsigned int b) const
806 {
807  return matrix_free_transfer_vector[b];
808 }
809 
810 
811 
812 template <int dim, typename Number>
813 void
815 {
816  this->matrix_free_transfer_vector.clear();
817 }
818 
819 
820 
821 // explicit instantiation
822 #include "mg_transfer_matrix_free.inst"
823 
824 
const Triangulation< dim, spacedim > & get_triangulation() const
bool has_level_dofs() const
Number local_element(const size_type local_index) const
size_type locally_owned_size() const
const std::shared_ptr< const Utilities::MPI::Partitioner > & get_partitioner() const
virtual void compress(::VectorOperation::values operation) override
void reinit(const size_type size, const bool omit_zeroing_entries=false)
SmartPointer< const MGConstrainedDoFs, MGLevelGlobalTransfer< LinearAlgebra::distributed::Vector< Number > > > mg_constrained_dofs
Definition: mg_transfer.h:586
std::vector< MGTransferMatrixFree< dim, Number > > matrix_free_transfer_vector
std::size_t memory_consumption() const
MGTransferBlockMatrixFree()=default
const MGTransferMatrixFree< dim, Number > & get_matrix_free_transfer(const unsigned int b) const override
void build(const DoFHandler< dim, dim > &dof_handler)
void initialize_constraints(const MGConstrainedDoFs &mg_constrained_dofs)
void do_restrict_add(const unsigned int from_level, LinearAlgebra::distributed::Vector< Number > &dst, const LinearAlgebra::distributed::Vector< Number > &src) const
void do_prolongate_add(const unsigned int to_level, LinearAlgebra::distributed::Vector< Number > &dst, const LinearAlgebra::distributed::Vector< Number > &src) const
virtual void prolongate(const unsigned int to_level, LinearAlgebra::distributed::Vector< Number > &dst, const LinearAlgebra::distributed::Vector< Number > &src) const override
void initialize_constraints(const MGConstrainedDoFs &mg_constrained_dofs)
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 >>())
virtual void prolongate_and_add(const unsigned int to_level, LinearAlgebra::distributed::Vector< Number > &dst, const LinearAlgebra::distributed::Vector< Number > &src) const override
std::size_t memory_consumption() const
virtual void restrict_and_add(const unsigned int from_level, LinearAlgebra::distributed::Vector< Number > &dst, const LinearAlgebra::distributed::Vector< Number > &src) const override
virtual unsigned int n_global_levels() const
#define DEAL_II_NAMESPACE_OPEN
Definition: config.h:458
#define DEAL_II_NAMESPACE_CLOSE
Definition: config.h:459
Point< 2 > second
Definition: grid_out.cc:4606
Point< 2 > first
Definition: grid_out.cc:4605
unsigned int level
Definition: grid_out.cc:4608
#define Assert(cond, exc)
Definition: exceptions.h:1509
#define AssertDimension(dim1, dim2)
Definition: exceptions.h:1703
#define AssertIndexRange(index, range)
Definition: exceptions.h:1768
static ::ExceptionBase & ExcIndexRange(std::size_t arg1, std::size_t arg2, std::size_t arg3)
static ::ExceptionBase & ExcMessage(std::string arg1)
void shift(const Tensor< 1, spacedim > &shift_vector, Triangulation< dim, spacedim > &triangulation)
Definition: grid_tools.cc:2117
std::enable_if_t< std::is_fundamental< T >::value, std::size_t > memory_consumption(const T &t)
SymmetricTensor< 2, dim, Number > b(const Tensor< 2, dim, Number > &F)
VectorType::value_type * begin(VectorType &V)
constexpr T pow(const T base, const int iexp)
Definition: utilities.h:450
T fixed_power(const T t)
Definition: utilities.h:969
void weight_fe_q_dofs_by_entity(const Number *weights, const unsigned int n_components, const int n_points_1d_non_template, Number *data)
std::vector< Number > prolongation_matrix_1d