Reference documentation for deal.II version GIT 741a3088b8 2023-04-01 07:05: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.h
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 #ifndef dealii_mg_transfer_matrix_free_h
17 #define dealii_mg_transfer_matrix_free_h
18 
19 #include <deal.II/base/config.h>
20 
23 
25 
28 
33 
34 
36 
37 
55 template <int dim, typename Number>
57  : public MGLevelGlobalTransfer<LinearAlgebra::distributed::Vector<Number>>
58 {
59 public:
65 
71 
75  virtual ~MGTransferMatrixFree() override = default;
76 
80  void
82 
86  void
87  clear();
88 
104  void
105  build(const DoFHandler<dim, dim> &dof_handler,
106  const std::vector<std::shared_ptr<const Utilities::MPI::Partitioner>>
107  &external_partitioners =
108  std::vector<std::shared_ptr<const Utilities::MPI::Partitioner>>());
109 
114  void
115  build(const DoFHandler<dim, dim> &dof_handler,
116  const std::function<void(const unsigned int,
119 
134  virtual void
135  prolongate(
136  const unsigned int to_level,
138  const LinearAlgebra::distributed::Vector<Number> &src) const override;
139 
140  virtual void
142  const unsigned int to_level,
144  const LinearAlgebra::distributed::Vector<Number> &src) const override;
145 
164  virtual void
166  const unsigned int from_level,
168  const LinearAlgebra::distributed::Vector<Number> &src) const override;
169 
184  template <typename BlockVectorType2, int spacedim>
185  void
187  const DoFHandler<dim, spacedim> & dof_handler,
189  const BlockVectorType2 & src) const;
190 
195 
199  std::size_t
200  memory_consumption() const;
201 
202 private:
208  unsigned int fe_degree;
209 
215 
220  unsigned int n_components;
221 
227  unsigned int n_child_cell_dofs;
228 
239  std::vector<std::vector<unsigned int>> level_dof_indices;
240 
245  std::vector<std::vector<std::pair<unsigned int, unsigned int>>>
247 
252  std::vector<unsigned int> n_owned_level_cells;
253 
259 
264 
276  std::vector<AlignedVector<VectorizedArray<Number>>> weights_on_refined;
277 
283  std::vector<std::vector<std::vector<unsigned short>>> dirichlet_indices;
284 
293 
297  template <int degree>
298  void
300  const unsigned int to_level,
303 
307  template <int degree>
308  void
309  do_restrict_add(const unsigned int from_level,
312 };
313 
314 
315 
321 template <int dim, typename Number, typename TransferType>
323  : public MGTransferBase<LinearAlgebra::distributed::BlockVector<Number>>
324 {
325 public:
328  {}
329 
344  virtual void
346  const unsigned int to_level,
348  const LinearAlgebra::distributed::BlockVector<Number> &src) const override;
349 
350  virtual void
352  const unsigned int to_level,
354  const LinearAlgebra::distributed::BlockVector<Number> &src) const override;
355 
374  virtual void
376  const unsigned int from_level,
378  const LinearAlgebra::distributed::BlockVector<Number> &src) const override;
379 
390  template <typename BlockVectorType2, int spacedim>
391  void
393  const DoFHandler<dim, spacedim> & dof_handler,
395  const BlockVectorType2 & src) const;
396 
400  template <typename BlockVectorType2, int spacedim>
401  void
403  const std::vector<const DoFHandler<dim, spacedim> *> & dof_handler,
405  const BlockVectorType2 & src) const;
406 
410  template <typename BlockVectorType2, int spacedim>
411  void
413  const DoFHandler<dim, spacedim> &dof_handler,
414  BlockVectorType2 & dst,
416  const;
417 
421  template <typename BlockVectorType2, int spacedim>
422  void
424  const std::vector<const DoFHandler<dim, spacedim> *> &dof_handler,
425  BlockVectorType2 & dst,
427  const;
428 
433  static const bool supports_dof_handler_vector = true;
434 
435 protected:
440  virtual const TransferType &
441  get_matrix_free_transfer(const unsigned int b) const = 0;
442 
447  const bool same_for_all;
448 };
449 
450 
451 
465 template <int dim, typename Number>
467  : public MGTransferBlockMatrixFreeBase<dim,
468  Number,
469  MGTransferMatrixFree<dim, Number>>
470 {
471 public:
477 
482  MGTransferBlockMatrixFree(const MGConstrainedDoFs &mg_constrained_dofs);
483 
488  const std::vector<MGConstrainedDoFs> &mg_constrained_dofs);
489 
493  virtual ~MGTransferBlockMatrixFree() override = default;
494 
498  void
499  initialize_constraints(const MGConstrainedDoFs &mg_constrained_dofs);
500 
504  void
506  const std::vector<MGConstrainedDoFs> &mg_constrained_dofs);
507 
511  void
512  clear();
513 
517  void
518  build(const DoFHandler<dim, dim> &dof_handler);
519 
523  void
524  build(const std::vector<const DoFHandler<dim, dim> *> &dof_handler);
525 
529  std::size_t
530  memory_consumption() const;
531 
532 protected:
534  get_matrix_free_transfer(const unsigned int b) const override;
535 
536 private:
540  std::vector<MGTransferMatrixFree<dim, Number>> matrix_free_transfer_vector;
541 };
542 
543 
547 //------------------------ templated functions -------------------------
548 #ifndef DOXYGEN
549 
550 
551 template <int dim, typename Number>
552 template <typename BlockVectorType2, int spacedim>
553 void
555  const DoFHandler<dim, spacedim> & dof_handler,
557  const BlockVectorType2 & src) const
558 {
559  const unsigned int min_level = dst.min_level();
560  const unsigned int max_level = dst.max_level();
561 
562  Assert(max_level == dof_handler.get_triangulation().n_global_levels() - 1,
564  max_level, dof_handler.get_triangulation().n_global_levels() - 1));
565 
566  const FiniteElement<dim, spacedim> &fe = dof_handler.get_fe();
567 
568  for (unsigned int level = min_level; level <= max_level; ++level)
569  if (dst[level].size() != dof_handler.n_dofs(level) ||
570  dst[level].locally_owned_size() !=
571  dof_handler.locally_owned_mg_dofs(level).n_elements())
572  dst[level].reinit(this->vector_partitioners[level]);
573 
574  // copy fine level vector to active cells in MG hierarchy
575  this->copy_to_mg(dof_handler, dst, src, true);
576 
577  // FIXME: maybe need to store hanging nodes constraints per level?
578  // MGConstrainedDoFs does NOT keep this info right now, only periodicity
579  // constraints...
580 
581  // do the transfer from level to level-1:
582  dst[max_level].update_ghost_values();
583 
584  for (unsigned int level = max_level; level > min_level; --level)
585  {
586  // auxiliary vector which always has ghost elements
587  const LinearAlgebra::distributed::Vector<Number> *input = nullptr;
589  if (dst[level].get_partitioner().get() ==
590  this->vector_partitioners[level].get())
591  input = &dst[level];
592  else
593  {
594  ghosted_fine.reinit(this->vector_partitioners[level]);
595  ghosted_fine.copy_locally_owned_data_from(dst[level]);
596  ghosted_fine.update_ghost_values();
597  input = &ghosted_fine;
598  }
599 
600  std::vector<Number> dof_values_coarse(fe.n_dofs_per_cell());
601  Vector<Number> dof_values_fine(fe.n_dofs_per_cell());
603  std::vector<types::global_dof_index> dof_indices(fe.n_dofs_per_cell());
604  for (const auto &cell : dof_handler.cell_iterators_on_level(level - 1))
605  if (cell->is_locally_owned_on_level())
606  {
607  // if we get to a cell without children (== active), we can
608  // skip it as there values should be already set by the
609  // equivalent of copy_to_mg()
610  if (cell->is_active())
611  continue;
612 
613  std::fill(dof_values_coarse.begin(), dof_values_coarse.end(), 0.);
614  for (unsigned int child = 0; child < cell->n_children(); ++child)
615  {
616  cell->child(child)->get_mg_dof_indices(dof_indices);
617 
619  this->mg_constrained_dofs, level, dof_indices);
620 
621  for (unsigned int i = 0; i < fe.n_dofs_per_cell(); ++i)
622  dof_values_fine(i) = (*input)(dof_indices[i]);
623  fe.get_restriction_matrix(child, cell->refinement_case())
624  .vmult(tmp, dof_values_fine);
625  for (unsigned int i = 0; i < fe.n_dofs_per_cell(); ++i)
626  if (fe.restriction_is_additive(i))
627  dof_values_coarse[i] += tmp[i];
628  else if (tmp(i) != 0.)
629  dof_values_coarse[i] = tmp[i];
630  }
631  cell->get_mg_dof_indices(dof_indices);
632  for (unsigned int i = 0; i < fe.n_dofs_per_cell(); ++i)
633  if (dof_handler.locally_owned_mg_dofs(level - 1).is_element(
634  dof_indices[i]))
635  dst[level - 1](dof_indices[i]) = dof_values_coarse[i];
636  }
637 
638  dst[level - 1].update_ghost_values();
639  }
640 }
641 
642 
643 
644 template <int dim, typename Number, typename TransferType>
645 template <typename BlockVectorType2, int spacedim>
646 void
648  const DoFHandler<dim, spacedim> & dof_handler,
650  const BlockVectorType2 & src) const
651 {
652  Assert(same_for_all,
653  ExcMessage(
654  "This object was initialized with support for usage with one "
655  "DoFHandler for each block, but this method assumes that "
656  "the same DoFHandler is used for all the blocks!"));
657  const std::vector<const DoFHandler<dim, spacedim> *> mg_dofs(src.n_blocks(),
658  &dof_handler);
659 
660  copy_to_mg(mg_dofs, dst, src);
661 }
662 
663 
664 
665 template <int dim, typename Number, typename TransferType>
666 template <typename BlockVectorType2, int spacedim>
667 void
669  const std::vector<const DoFHandler<dim, spacedim> *> & dof_handler,
671  const BlockVectorType2 & src) const
672 {
673  const unsigned int n_blocks = src.n_blocks();
674  AssertDimension(dof_handler.size(), n_blocks);
675 
676  if (n_blocks == 0)
677  return;
678 
679  const unsigned int min_level = dst.min_level();
680  const unsigned int max_level = dst.max_level();
681 
682  for (unsigned int level = min_level; level <= max_level; ++level)
683  if (dst[level].n_blocks() != n_blocks)
684  dst[level].reinit(n_blocks);
685 
686  // FIXME: this a quite ugly as we need a temporary object:
688  min_level, max_level);
689 
690  for (unsigned int b = 0; b < n_blocks; ++b)
691  {
692  const unsigned int data_block = same_for_all ? 0 : b;
693  get_matrix_free_transfer(data_block)
694  .copy_to_mg(*dof_handler[b], dst_non_block, src.block(b));
695 
696  for (unsigned int l = min_level; l <= max_level; ++l)
697  dst[l].block(b) = dst_non_block[l];
698  }
699 
700  for (unsigned int level = min_level; level <= max_level; ++level)
701  dst[level].collect_sizes();
702 }
703 
704 template <int dim, typename Number, typename TransferType>
705 template <typename BlockVectorType2, int spacedim>
706 void
708  const DoFHandler<dim, spacedim> &dof_handler,
709  BlockVectorType2 & dst,
711  const
712 {
713  Assert(same_for_all,
714  ExcMessage(
715  "This object was initialized with support for usage with one "
716  "DoFHandler for each block, but this method assumes that "
717  "the same DoFHandler is used for all the blocks!"));
718  const std::vector<const DoFHandler<dim, spacedim> *> mg_dofs(dst.n_blocks(),
719  &dof_handler);
720 
721  copy_from_mg(mg_dofs, dst, src);
722 }
723 
724 template <int dim, typename Number, typename TransferType>
725 template <typename BlockVectorType2, int spacedim>
726 void
728  const std::vector<const DoFHandler<dim, spacedim> *> &dof_handler,
729  BlockVectorType2 & dst,
731  const
732 {
733  const unsigned int n_blocks = dst.n_blocks();
734  AssertDimension(dof_handler.size(), n_blocks);
735 
736  if (n_blocks == 0)
737  return;
738 
739  const unsigned int min_level = src.min_level();
740  const unsigned int max_level = src.max_level();
741 
742  for (unsigned int l = min_level; l <= max_level; ++l)
743  AssertDimension(src[l].n_blocks(), dst.n_blocks());
744 
745  // FIXME: this a quite ugly as we need a temporary object:
747  min_level, max_level);
748 
749  for (unsigned int b = 0; b < n_blocks; ++b)
750  {
751  for (unsigned int l = min_level; l <= max_level; ++l)
752  {
753  src_non_block[l].reinit(src[l].block(b));
754  src_non_block[l] = src[l].block(b);
755  }
756  const unsigned int data_block = same_for_all ? 0 : b;
757  get_matrix_free_transfer(data_block)
758  .copy_from_mg(*dof_handler[b], dst.block(b), src_non_block);
759  }
760 }
761 
762 
763 
764 #endif // DOXYGEN
765 
766 
768 
769 #endif
const Triangulation< dim, spacedim > & get_triangulation() const
types::global_dof_index n_dofs() const
const IndexSet & locally_owned_mg_dofs(const unsigned int level) const
const FiniteElement< dim, spacedim > & get_fe(const types::fe_index index=0) const
unsigned int n_dofs_per_cell() const
bool restriction_is_additive(const unsigned int index) const
virtual const FullMatrix< double > & get_restriction_matrix(const unsigned int child, const RefinementCase< dim > &refinement_case=RefinementCase< dim >::isotropic_refinement) const
size_type n_elements() const
Definition: index_set.h:1796
bool is_element(const size_type index) const
Definition: index_set.h:1756
void copy_locally_owned_data_from(const Vector< Number2, MemorySpace > &src)
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::function< void(const unsigned int, LinearAlgebra::distributed::Vector< Number > &)> initialize_dof_vector
Definition: mg_transfer.h:619
virtual void prolongate(const unsigned int to_level, LinearAlgebra::distributed::BlockVector< Number > &dst, const LinearAlgebra::distributed::BlockVector< Number > &src) const override
virtual const TransferType & get_matrix_free_transfer(const unsigned int b) const =0
void copy_to_mg(const DoFHandler< dim, spacedim > &dof_handler, MGLevelObject< LinearAlgebra::distributed::BlockVector< Number >> &dst, const BlockVectorType2 &src) const
virtual void prolongate_and_add(const unsigned int to_level, LinearAlgebra::distributed::BlockVector< Number > &dst, const LinearAlgebra::distributed::BlockVector< Number > &src) const override
void copy_from_mg(const std::vector< const DoFHandler< dim, spacedim > * > &dof_handler, BlockVectorType2 &dst, const MGLevelObject< LinearAlgebra::distributed::BlockVector< Number >> &src) const
virtual void restrict_and_add(const unsigned int from_level, LinearAlgebra::distributed::BlockVector< Number > &dst, const LinearAlgebra::distributed::BlockVector< Number > &src) const override
MGTransferBlockMatrixFreeBase(const bool same_for_all)
void copy_to_mg(const std::vector< const DoFHandler< dim, spacedim > * > &dof_handler, MGLevelObject< LinearAlgebra::distributed::BlockVector< Number >> &dst, const BlockVectorType2 &src) const
void copy_from_mg(const DoFHandler< dim, spacedim > &dof_handler, BlockVectorType2 &dst, const MGLevelObject< LinearAlgebra::distributed::BlockVector< Number >> &src) const
std::vector< MGTransferMatrixFree< dim, Number > > matrix_free_transfer_vector
std::size_t memory_consumption() const
MGTransferBlockMatrixFree()=default
virtual ~MGTransferBlockMatrixFree() override=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 ~MGTransferMatrixFree() override=default
std::vector< std::vector< unsigned int > > level_dof_indices
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)
std::vector< std::vector< std::vector< unsigned short > > > dirichlet_indices
MGLevelObject< std::shared_ptr< const Utilities::MPI::Partitioner > > vector_partitioners
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::vector< std::vector< std::pair< unsigned int, unsigned int > > > parent_child_connect
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
std::vector< AlignedVector< VectorizedArray< Number > > > weights_on_refined
void interpolate_to_mg(const DoFHandler< dim, spacedim > &dof_handler, MGLevelObject< LinearAlgebra::distributed::Vector< Number >> &dst, const BlockVectorType2 &src) const
std::vector< unsigned int > n_owned_level_cells
AlignedVector< VectorizedArray< Number > > prolongation_matrix_1d
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
Definition: vector.h:109
void update_ghost_values() const
#define DEAL_II_NAMESPACE_OPEN
Definition: config.h:474
#define DEAL_II_NAMESPACE_CLOSE
Definition: config.h:475
unsigned int level
Definition: grid_out.cc:4618
IteratorRange< cell_iterator > cell_iterators_on_level(const unsigned int level) const
#define DeclException0(Exception0)
Definition: exceptions.h:465
static ::ExceptionBase & ExcNoProlongation()
static ::ExceptionBase & ExcDimensionMismatch(std::size_t arg1, std::size_t arg2)
#define Assert(cond, exc)
Definition: exceptions.h:1586
#define AssertDimension(dim1, dim2)
Definition: exceptions.h:1759
static ::ExceptionBase & ExcMessage(std::string arg1)
std::enable_if_t< IsBlockVector< VectorType >::value, unsigned int > n_blocks(const VectorType &vector)
Definition: operators.h:50
std::enable_if_t< IsBlockVector< VectorType >::value, void > collect_sizes(VectorType &vector)
Definition: operators.h:96
Tensor< 2, dim, Number > l(const Tensor< 2, dim, Number > &F, const Tensor< 2, dim, Number > &dF_dt)
SymmetricTensor< 2, dim, Number > b(const Tensor< 2, dim, Number > &F)
void resolve_identity_constraints(const MGConstrainedDoFs *mg_constrained_dofs, const unsigned int level, std::vector< types::global_dof_index > &dof_indices)