17 #include <deal.II/base/function.h> 18 #include <deal.II/base/logstream.h> 19 #include <deal.II/base/vectorization.h> 21 #include <deal.II/dofs/dof_accessor.h> 22 #include <deal.II/dofs/dof_tools.h> 24 #include <deal.II/fe/fe.h> 26 #include <deal.II/grid/tria.h> 27 #include <deal.II/grid/tria_iterator.h> 29 #include <deal.II/lac/la_parallel_vector.h> 31 #include <deal.II/matrix_free/evaluation_kernels.h> 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> 39 DEAL_II_NAMESPACE_OPEN
42 template <
int dim,
typename Number>
45 , element_is_continuous(false)
47 , n_child_cell_dofs(0)
52 template <
int dim,
typename Number>
65 template <
int dim,
typename Number>
75 template <
int dim,
typename Number>
95 template <
int dim,
typename Number>
101 std::vector<std::vector<Number>> weights_unvectorized;
105 internal::MGTransfer::setup_transfer<dim, Number>(
113 weights_unvectorized,
124 for (
unsigned int i = 0; i < elem_info.prolongation_matrix_1d.size(); i++)
129 const unsigned int n_levels = mg_dof.get_triangulation().n_global_levels();
131 const unsigned int n_weights_per_cell = Utilities::fixed_power<dim>(3);
133 for (
unsigned int level = 1; level < n_levels; ++level)
136 ((n_owned_level_cells[level - 1] + vec_size - 1) / vec_size) *
139 for (
unsigned int c = 0; c < n_owned_level_cells[level - 1]; ++c)
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)
146 weights_unvectorized[level - 1][c * n_weights_per_cell + i];
156 template <
int dim,
typename Number>
159 const unsigned int to_level,
179 do_prolongate_add<0>(to_level,
181 this->ghosted_level_vector[to_level - 1]);
183 do_prolongate_add<1>(to_level,
184 this->ghosted_level_vector[to_level],
185 this->ghosted_level_vector[to_level - 1]);
187 do_prolongate_add<2>(to_level,
188 this->ghosted_level_vector[to_level],
189 this->ghosted_level_vector[to_level - 1]);
191 do_prolongate_add<3>(to_level,
192 this->ghosted_level_vector[to_level],
193 this->ghosted_level_vector[to_level - 1]);
195 do_prolongate_add<4>(to_level,
196 this->ghosted_level_vector[to_level],
197 this->ghosted_level_vector[to_level - 1]);
199 do_prolongate_add<5>(to_level,
200 this->ghosted_level_vector[to_level],
201 this->ghosted_level_vector[to_level - 1]);
203 do_prolongate_add<6>(to_level,
204 this->ghosted_level_vector[to_level],
205 this->ghosted_level_vector[to_level - 1]);
207 do_prolongate_add<7>(to_level,
208 this->ghosted_level_vector[to_level],
209 this->ghosted_level_vector[to_level - 1]);
211 do_prolongate_add<8>(to_level,
212 this->ghosted_level_vector[to_level],
213 this->ghosted_level_vector[to_level - 1]);
215 do_prolongate_add<9>(to_level,
216 this->ghosted_level_vector[to_level],
217 this->ghosted_level_vector[to_level - 1]);
219 do_prolongate_add<10>(to_level,
220 this->ghosted_level_vector[to_level],
221 this->ghosted_level_vector[to_level - 1]);
224 this->ghosted_level_vector[to_level],
225 this->ghosted_level_vector[to_level - 1]);
233 template <
int dim,
typename Number>
236 const unsigned int from_level,
253 do_restrict_add<0>(from_level,
255 this->ghosted_level_vector[from_level]);
257 do_restrict_add<1>(from_level,
258 this->ghosted_level_vector[from_level - 1],
259 this->ghosted_level_vector[from_level]);
261 do_restrict_add<2>(from_level,
262 this->ghosted_level_vector[from_level - 1],
263 this->ghosted_level_vector[from_level]);
265 do_restrict_add<3>(from_level,
266 this->ghosted_level_vector[from_level - 1],
267 this->ghosted_level_vector[from_level]);
269 do_restrict_add<4>(from_level,
270 this->ghosted_level_vector[from_level - 1],
271 this->ghosted_level_vector[from_level]);
273 do_restrict_add<5>(from_level,
274 this->ghosted_level_vector[from_level - 1],
275 this->ghosted_level_vector[from_level]);
277 do_restrict_add<6>(from_level,
278 this->ghosted_level_vector[from_level - 1],
279 this->ghosted_level_vector[from_level]);
281 do_restrict_add<7>(from_level,
282 this->ghosted_level_vector[from_level - 1],
283 this->ghosted_level_vector[from_level]);
285 do_restrict_add<8>(from_level,
286 this->ghosted_level_vector[from_level - 1],
287 this->ghosted_level_vector[from_level]);
289 do_restrict_add<9>(from_level,
290 this->ghosted_level_vector[from_level - 1],
291 this->ghosted_level_vector[from_level]);
293 do_restrict_add<10>(from_level,
294 this->ghosted_level_vector[from_level - 1],
295 this->ghosted_level_vector[from_level]);
299 this->ghosted_level_vector[from_level - 1],
300 this->ghosted_level_vector[from_level]);
303 dst += this->ghosted_level_vector[from_level - 1];
310 template <
int dim,
int degree,
typename Number>
319 const int loop_length = degree != -1 ? 2 * degree + 1 : 2 * fe_degree + 1;
320 unsigned int degree_to_3[100];
322 for (
int i = 1; i < loop_length - 1; ++i)
324 degree_to_3[loop_length - 1] = 2;
326 for (
int k = 0; k < (dim > 2 ? loop_length : 1); ++k)
327 for (
int j = 0; j < (dim > 1 ? loop_length : 1); ++j)
329 const unsigned int shift = 9 * degree_to_3[k] + 3 * degree_to_3[j];
330 data[0] *= weights[shift];
333 for (
int i = 1; i < loop_length - 1; ++i)
334 data[i] *= weights[shift + 1];
335 data[loop_length - 1] *= weights[shift + 2];
343 template <
int dim,
typename Number>
344 template <
int degree>
347 const unsigned int to_level,
352 const unsigned int degree_size = (degree > -1 ? degree :
fe_degree) + 1;
354 const unsigned int n_scalar_cell_dofs =
355 Utilities::fixed_power<dim>(n_child_dofs_1d);
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 :
367 for (
unsigned int v = 0; v < n_lanes; ++v)
369 const unsigned int shift =
370 internal::MGTransfer::compute_shift_within_children<dim>(
374 const unsigned int *indices =
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)
386 indices[c * n_scalar_cell_dofs +
387 k * n_child_dofs_1d * n_child_dofs_1d +
388 j * n_child_dofs_1d + i]);
391 for (std::vector<unsigned short>::const_iterator i =
400 degree_size * n_child_dofs_1d);
402 if (element_is_continuous)
416 c * Utilities::fixed_power<dim>(degree_size),
420 weight_dofs_on_child<dim, degree, Number>(
438 c * Utilities::fixed_power<dim>(degree_size),
445 const unsigned int *indices =
447 for (
unsigned int v = 0; v < n_lanes; ++v)
458 template <
int dim,
typename Number>
459 template <
int degree>
462 const unsigned int from_level,
467 const unsigned int degree_size = (degree > -1 ? degree :
fe_degree) + 1;
469 const unsigned int n_scalar_cell_dofs =
470 Utilities::fixed_power<dim>(n_child_dofs_1d);
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 :
483 const unsigned int *indices =
485 for (
unsigned int v = 0; v < n_lanes; ++v)
494 degree_size * n_child_dofs_1d);
496 if (element_is_continuous)
498 weight_dofs_on_child<dim, degree, Number>(
500 [(cell / vec_size) * three_to_dim],
516 c * Utilities::fixed_power<dim>(degree_size),
534 c * Utilities::fixed_power<dim>(degree_size),
540 for (
unsigned int v = 0; v < n_lanes; ++v)
542 const unsigned int shift =
543 internal::MGTransfer::compute_shift_within_children<dim>(
552 const unsigned int *indices =
561 for (std::vector<unsigned short>::const_iterator i =
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)
571 indices[c * n_scalar_cell_dofs +
572 k * n_child_dofs_1d * n_child_dofs_1d +
573 j * n_child_dofs_1d + i]) +=
582 template <
int dim,
typename Number>
600 template <
int dim,
typename Number>
610 template <
int dim,
typename Number>
612 const std::vector<MGConstrainedDoFs> &mg_c)
615 for (
unsigned int i = 0; i < mg_c.size(); ++i)
621 template <
int dim,
typename Number>
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!"));
637 template <
int dim,
typename Number>
640 const std::vector<MGConstrainedDoFs> &mg_c)
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 " 649 for (
unsigned int i = 0; i < mg_c.size(); ++i)
655 template <
int dim,
typename Number>
664 template <
int dim,
typename Number>
675 template <
int dim,
typename Number>
681 for (
unsigned int i = 0; i < mg_dof.size(); ++i)
687 template <
int dim,
typename Number>
690 const unsigned int to_level,
694 const unsigned int n_blocks = src.
n_blocks();
700 for (
unsigned int b = 0; b < n_blocks; ++b)
711 template <
int dim,
typename Number>
714 const unsigned int from_level,
718 const unsigned int n_blocks = src.
n_blocks();
724 for (
unsigned int b = 0; b < n_blocks; ++b)
735 template <
int dim,
typename Number>
739 std::size_t total_memory_consumption = 0;
741 total_memory_consumption += el.memory_consumption();
742 return total_memory_consumption;
748 #include "mg_transfer_matrix_free.inst" 751 DEAL_II_NAMESPACE_CLOSE
void do_prolongate_add(const unsigned int to_level, LinearAlgebra::distributed::Vector< Number > &dst, const LinearAlgebra::distributed::Vector< Number > &src) const
#define AssertDimension(dim1, dim2)
Number local_element(const size_type local_index) const
std::vector< unsigned int > n_owned_level_cells
AlignedVector< VectorizedArray< Number > > evaluation_data
std::vector< std::vector< std::pair< unsigned int, unsigned int > > > copy_indices_global_mine
void copy_locally_owned_data_from(const Vector< Number2 > &src)
unsigned int n_components
std::vector< std::vector< std::vector< unsigned short > > > dirichlet_indices
void build(const DoFHandler< dim, dim > &mg_dof)
void initialize_constraints(const MGConstrainedDoFs &mg_constrained_dofs)
#define AssertIndexRange(index, range)
unsigned int n_blocks() const
std::size_t memory_consumption() const
static::ExceptionBase & ExcIndexRange(int arg1, int arg2, int arg3)
void resize(const size_type size_in)
size_type local_size() const
static::ExceptionBase & ExcMessage(std::string arg1)
AlignedVector< VectorizedArray< Number > > prolongation_matrix_1d
#define Assert(cond, exc)
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
std::size_t memory_consumption() const
std::vector< AlignedVector< VectorizedArray< Number > > > weights_on_refined
void fill_and_communicate_copy_indices(const DoFHandler< dim, spacedim > &mg_dof)
void do_restrict_add(const unsigned int from_level, LinearAlgebra::distributed::Vector< Number > &dst, const LinearAlgebra::distributed::Vector< Number > &src) const
bool element_is_continuous
std::vector< std::vector< unsigned int > > level_dof_indices
SmartPointer< const MGConstrainedDoFs, MGLevelGlobalTransfer< LinearAlgebra::distributed::Vector< Number > > > mg_constrained_dofs
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
std::vector< MGTransferMatrixFree< dim, Number > > matrix_free_transfer_vector
MGTransferBlockMatrixFree()=default
static::ExceptionBase & ExcNotImplemented()
unsigned int n_child_cell_dofs
unsigned int n_components(const DoFHandler< dim, spacedim > &dh)
void initialize_constraints(const MGConstrainedDoFs &mg_constrained_dofs)
BlockType & block(const unsigned int i)
MGLevelObject< LinearAlgebra::distributed::Vector< Number > > ghosted_level_vector
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)