Reference documentation for deal.II version GIT 6a72d26406 2023-06-07 13:05:01+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.h
Go to the documentation of this file.
1 // ---------------------------------------------------------------------
2 //
3 // Copyright (C) 2001 - 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_h
17 #define dealii_mg_transfer_h
18 
19 #include <deal.II/base/config.h>
20 
22 
24 
26 
36 
39 
40 #include <memory>
41 
42 
44 
45 
46 namespace internal
47 {
48  template <typename VectorType>
50  {
53 
54  static const bool requires_distributed_sparsity_pattern = false;
55 
56  template <typename SparsityPatternType, int dim, int spacedim>
57  static void
59  Sparsity & sparsity,
60  int level,
61  const SparsityPatternType &sp,
63  {
64  sparsity.copy_from(sp);
65  (void)level;
66  matrix.reinit(sparsity);
67  }
68  };
69 
70 #ifdef DEAL_II_WITH_TRILINOS
71  template <typename Number>
73  {
76 
77  static const bool requires_distributed_sparsity_pattern = false;
78 
79  template <typename SparsityPatternType, int dim, int spacedim>
80  static void
82  Sparsity &,
83  int level,
84  const SparsityPatternType & sp,
85  const DoFHandler<dim, spacedim> &dh)
86  {
87  const MPI_Comm communicator = dh.get_communicator();
88 
89  matrix.reinit(dh.locally_owned_mg_dofs(level + 1),
91  sp,
92  communicator,
93  true);
94  }
95  };
96 
97  template <>
99  {
102 
103  static const bool requires_distributed_sparsity_pattern = false;
104 
105  template <typename SparsityPatternType, int dim, int spacedim>
106  static void
108  Sparsity &,
109  int level,
110  const SparsityPatternType & sp,
111  const DoFHandler<dim, spacedim> &dh)
112  {
113  const MPI_Comm communicator = dh.get_communicator();
114 
115  matrix.reinit(dh.locally_owned_mg_dofs(level + 1),
117  sp,
118  communicator,
119  true);
120  }
121  };
122 
123 # ifdef DEAL_II_WITH_MPI
124 # ifdef DEAL_II_TRILINOS_WITH_TPETRA
125  template <typename Number>
127  {
130 
131  static const bool requires_distributed_sparsity_pattern = false;
132 
133  template <typename SparsityPatternType, int dim, int spacedim>
134  static void
136  Sparsity &,
137  int level,
138  const SparsityPatternType & sp,
139  const DoFHandler<dim, spacedim> &dh)
140  {
141  const MPI_Comm communicator = dh.get_communicator();
142 
143  matrix.reinit(dh.locally_owned_mg_dofs(level + 1),
145  sp,
146  communicator,
147  true);
148  }
149  };
150 # endif
151 
152  template <>
154  {
157 
158  static const bool requires_distributed_sparsity_pattern = false;
159 
160  template <typename SparsityPatternType, int dim, int spacedim>
161  static void
163  Sparsity &,
164  int level,
165  const SparsityPatternType & sp,
166  const DoFHandler<dim, spacedim> &dh)
167  {
168  const MPI_Comm communicator = dh.get_communicator();
169 
170  matrix.reinit(dh.locally_owned_mg_dofs(level + 1),
172  sp,
173  communicator,
174  true);
175  }
176  };
177 # endif
178 
179 #else
180  // ! DEAL_II_WITH_TRILINOS
181  template <typename Number>
182  struct MatrixSelector<LinearAlgebra::distributed::Vector<Number>>
183  {
184  using Sparsity = ::SparsityPattern;
186 
187  static const bool requires_distributed_sparsity_pattern = false;
188 
189  template <typename SparsityPatternType, int dim, int spacedim>
190  static void
191  reinit(Matrix &,
192  Sparsity &,
193  int,
194  const SparsityPatternType &,
196  {
197  AssertThrow(
198  false,
200  "ERROR: MGTransferPrebuilt with LinearAlgebra::distributed::Vector currently "
201  "needs deal.II to be configured with Trilinos."));
202  }
203  };
204 
205 #endif
206 
207 #ifdef DEAL_II_WITH_PETSC
208  template <>
210  {
213 
214  static const bool requires_distributed_sparsity_pattern = true;
215 
216  template <typename SparsityPatternType, int dim, int spacedim>
217  static void
219  Sparsity &,
220  int level,
221  const SparsityPatternType & sp,
222  const DoFHandler<dim, spacedim> &dh)
223  {
224  const MPI_Comm communicator = dh.get_communicator();
225 
226  // Reinit PETSc matrix
227  matrix.reinit(dh.locally_owned_mg_dofs(level + 1),
229  sp,
230  communicator);
231  }
232  };
233 #endif
234 } // namespace internal
235 
236 /*
237  * MGTransferBase is defined in mg_base.h
238  */
239 
249 template <typename VectorType>
250 class MGLevelGlobalTransfer : public MGTransferBase<VectorType>
251 {
252 public:
256  void
257  clear();
258 
265  template <int dim, class InVector, int spacedim>
266  void
269  const InVector & src) const;
270 
278  template <int dim, class OutVector, int spacedim>
279  void
281  OutVector & dst,
282  const MGLevelObject<VectorType> &src) const;
283 
289  template <int dim, class OutVector, int spacedim>
290  void
292  OutVector & dst,
293  const MGLevelObject<VectorType> &src) const;
294 
310  void
311  set_component_to_block_map(const std::vector<unsigned int> &map);
312 
316  std::size_t
317  memory_consumption() const;
318 
322  void
323  print_indices(std::ostream &os) const;
324 
325 protected:
329  template <int dim, int spacedim>
330  void
332  const DoFHandler<dim, spacedim> &dof_handler);
333 
337  std::vector<types::global_dof_index> sizes;
338 
346  std::vector<
347  std::vector<std::pair<types::global_dof_index, types::global_dof_index>>>
349 
357  std::vector<
358  std::vector<std::pair<types::global_dof_index, types::global_dof_index>>>
360 
368  std::vector<
369  std::vector<std::pair<types::global_dof_index, types::global_dof_index>>>
371 
379 
384  std::vector<unsigned int> component_to_block_map;
385 
391 
392 private:
396  template <int dim, int spacedim>
397  void
398  assert_built(const DoFHandler<dim, spacedim> &dof_handler) const;
399 };
400 
401 
402 
411 template <typename Number>
413  : public MGTransferBase<LinearAlgebra::distributed::Vector<Number>>
414 {
415 public:
419  void
420  clear();
421 
428  template <int dim, typename Number2, int spacedim>
429  void
433 
441  template <int dim, typename Number2, int spacedim>
442  void
444  const DoFHandler<dim, spacedim> & dof_handler,
447 
453  template <int dim, typename Number2, int spacedim>
454  void
456  const DoFHandler<dim, spacedim> & dof_handler,
459 
475  void
476  set_component_to_block_map(const std::vector<unsigned int> &map);
477 
481  std::size_t
482  memory_consumption() const;
483 
487  void
488  print_indices(std::ostream &os) const;
489 
490 protected:
495  template <int dim, typename Number2, int spacedim>
496  void
500  const bool solution_transfer) const;
501 
505  template <int dim, int spacedim>
506  void
508  const DoFHandler<dim, spacedim> &dof_handler);
509 
513  std::vector<types::global_dof_index> sizes;
514 
523  std::vector<Table<2, unsigned int>> copy_indices;
524 
528  std::vector<Table<2, unsigned int>> solution_copy_indices;
529 
537  std::vector<Table<2, unsigned int>> copy_indices_global_mine;
538 
542  std::vector<Table<2, unsigned int>> solution_copy_indices_global_mine;
543 
551  std::vector<Table<2, unsigned int>> copy_indices_level_mine;
552 
556  std::vector<Table<2, unsigned int>> solution_copy_indices_level_mine;
557 
565 
573 
578  std::vector<unsigned int> component_to_block_map;
579 
583  SmartPointer<
584  const MGConstrainedDoFs,
587 
594 
600 
607 
613 
617  std::function<void(const unsigned int,
620 
621 private:
625  template <int dim, int spacedim>
626  void
627  assert_built(const DoFHandler<dim, spacedim> &dof_handler) const;
628 };
629 
630 
631 
642 template <typename VectorType>
643 class MGTransferPrebuilt : public MGLevelGlobalTransfer<VectorType>
644 {
645 public:
650  MGTransferPrebuilt() = default;
651 
657 
661  virtual ~MGTransferPrebuilt() override = default;
662 
666  void
668 
672  void
673  clear();
674 
679  template <int dim, int spacedim>
680  void
681  build(const DoFHandler<dim, spacedim> &dof_handler);
682 
694  virtual void
695  prolongate(const unsigned int to_level,
696  VectorType & dst,
697  const VectorType & src) const override;
698 
714  virtual void
715  restrict_and_add(const unsigned int from_level,
716  VectorType & dst,
717  const VectorType & src) const override;
718 
723 
728 
732  std::size_t
733  memory_consumption() const;
734 
738  void
739  print_matrices(std::ostream &os) const;
740 
741 private:
745  std::vector<
746  std::shared_ptr<typename internal::MatrixSelector<VectorType>::Sparsity>>
748 
754  std::vector<
755  std::shared_ptr<typename internal::MatrixSelector<VectorType>::Matrix>>
757 
762  std::vector<std::vector<bool>> interface_dofs;
763 };
764 
765 
770 
771 #endif
const IndexSet & locally_owned_mg_dofs(const unsigned int level) const
MPI_Comm get_communicator() const
SmartPointer< const MGConstrainedDoFs, MGLevelGlobalTransfer< LinearAlgebra::distributed::Vector< Number > > > mg_constrained_dofs
Definition: mg_transfer.h:586
std::vector< Table< 2, unsigned int > > solution_copy_indices_level_mine
Definition: mg_transfer.h:556
void copy_from_mg(const DoFHandler< dim, spacedim > &dof_handler, LinearAlgebra::distributed::Vector< Number2 > &dst, const MGLevelObject< LinearAlgebra::distributed::Vector< Number >> &src) const
void copy_to_mg(const DoFHandler< dim, spacedim > &dof_handler, MGLevelObject< LinearAlgebra::distributed::Vector< Number >> &dst, const LinearAlgebra::distributed::Vector< Number2 > &src, const bool solution_transfer) const
void copy_to_mg(const DoFHandler< dim, spacedim > &dof_handler, MGLevelObject< LinearAlgebra::distributed::Vector< Number >> &dst, const LinearAlgebra::distributed::Vector< Number2 > &src) const
MGLevelObject< LinearAlgebra::distributed::Vector< Number > > ghosted_level_vector
Definition: mg_transfer.h:606
void set_component_to_block_map(const std::vector< unsigned int > &map)
std::vector< Table< 2, unsigned int > > solution_copy_indices
Definition: mg_transfer.h:528
std::vector< Table< 2, unsigned int > > solution_copy_indices_global_mine
Definition: mg_transfer.h:542
LinearAlgebra::distributed::Vector< Number > solution_ghosted_global_vector
Definition: mg_transfer.h:599
void assert_built(const DoFHandler< dim, spacedim > &dof_handler) const
std::vector< Table< 2, unsigned int > > copy_indices_global_mine
Definition: mg_transfer.h:537
MGLevelObject< LinearAlgebra::distributed::Vector< Number > > solution_ghosted_level_vector
Definition: mg_transfer.h:612
LinearAlgebra::distributed::Vector< Number > ghosted_global_vector
Definition: mg_transfer.h:593
std::function< void(const unsigned int, LinearAlgebra::distributed::Vector< Number > &)> initialize_dof_vector
Definition: mg_transfer.h:619
std::vector< Table< 2, unsigned int > > copy_indices_level_mine
Definition: mg_transfer.h:551
void copy_from_mg_add(const DoFHandler< dim, spacedim > &dof_handler, LinearAlgebra::distributed::Vector< Number2 > &dst, const MGLevelObject< LinearAlgebra::distributed::Vector< Number >> &src) const
void print_indices(std::ostream &os) const
void assert_built(const DoFHandler< dim, spacedim > &dof_handler) const
std::vector< std::vector< std::pair< types::global_dof_index, types::global_dof_index > > > copy_indices_level_mine
Definition: mg_transfer.h:370
std::size_t memory_consumption() const
void set_component_to_block_map(const std::vector< unsigned int > &map)
std::vector< std::vector< std::pair< types::global_dof_index, types::global_dof_index > > > copy_indices
Definition: mg_transfer.h:348
void fill_and_communicate_copy_indices(const DoFHandler< dim, spacedim > &dof_handler)
std::vector< unsigned int > component_to_block_map
Definition: mg_transfer.h:384
void copy_from_mg_add(const DoFHandler< dim, spacedim > &dof_handler, OutVector &dst, const MGLevelObject< VectorType > &src) const
std::vector< types::global_dof_index > sizes
Definition: mg_transfer.h:337
std::vector< std::vector< std::pair< types::global_dof_index, types::global_dof_index > > > copy_indices_global_mine
Definition: mg_transfer.h:359
void copy_from_mg(const DoFHandler< dim, spacedim > &dof_handler, OutVector &dst, const MGLevelObject< VectorType > &src) const
void copy_to_mg(const DoFHandler< dim, spacedim > &dof_handler, MGLevelObject< VectorType > &dst, const InVector &src) const
SmartPointer< const MGConstrainedDoFs, MGLevelGlobalTransfer< VectorType > > mg_constrained_dofs
Definition: mg_transfer.h:390
virtual ~MGTransferPrebuilt() override=default
std::vector< std::vector< bool > > interface_dofs
Definition: mg_transfer.h:762
virtual void restrict_and_add(const unsigned int from_level, VectorType &dst, const VectorType &src) const override
void build(const DoFHandler< dim, spacedim > &dof_handler)
std::vector< std::shared_ptr< typename internal::MatrixSelector< VectorType >::Sparsity > > prolongation_sparsities
Definition: mg_transfer.h:747
std::size_t memory_consumption() const
MGTransferPrebuilt()=default
virtual void prolongate(const unsigned int to_level, VectorType &dst, const VectorType &src) const override
std::vector< std::shared_ptr< typename internal::MatrixSelector< VectorType >::Matrix > > prolongation_matrices
Definition: mg_transfer.h:756
void print_matrices(std::ostream &os) const
void initialize_constraints(const MGConstrainedDoFs &mg_constrained_dofs)
void copy_from(const size_type n_rows, const size_type n_cols, const ForwardIterator begin, const ForwardIterator end)
#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
#define DeclException0(Exception0)
Definition: exceptions.h:465
static ::ExceptionBase & ExcMatricesNotBuilt()
static ::ExceptionBase & ExcNotImplemented()
static ::ExceptionBase & ExcNoProlongation()
#define AssertThrow(cond, exc)
Definition: exceptions.h:1703
@ matrix
Contents is actually a matrix.
SparseMatrix< double > SparseMatrix
static void reinit(Matrix &matrix, Sparsity &, int level, const SparsityPatternType &sp, const DoFHandler< dim, spacedim > &dh)
Definition: mg_transfer.h:81
static void reinit(Matrix &matrix, Sparsity &, int level, const SparsityPatternType &sp, const DoFHandler< dim, spacedim > &dh)
Definition: mg_transfer.h:162
static void reinit(Matrix &matrix, Sparsity &, int level, const SparsityPatternType &sp, const DoFHandler< dim, spacedim > &dh)
Definition: mg_transfer.h:135
static void reinit(Matrix &matrix, Sparsity &, int level, const SparsityPatternType &sp, const DoFHandler< dim, spacedim > &dh)
Definition: mg_transfer.h:218
static void reinit(Matrix &matrix, Sparsity &, int level, const SparsityPatternType &sp, const DoFHandler< dim, spacedim > &dh)
Definition: mg_transfer.h:107
static const bool requires_distributed_sparsity_pattern
Definition: mg_transfer.h:54
static void reinit(Matrix &matrix, Sparsity &sparsity, int level, const SparsityPatternType &sp, const DoFHandler< dim, spacedim > &)
Definition: mg_transfer.h:58
::SparseMatrix< typename VectorType::value_type > Matrix
Definition: mg_transfer.h:52
::SparsityPattern Sparsity
Definition: mg_transfer.h:51