Reference documentation for deal.II version 9.5.0
\(\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\}}\)
Loading...
Searching...
No Matches
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
46namespace 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
58 reinit(Matrix & matrix,
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>
72 struct MatrixSelector<LinearAlgebra::distributed::Vector<Number>>
73 {
76
77 static const bool requires_distributed_sparsity_pattern = false;
78
79 template <typename SparsityPatternType, int dim, int spacedim>
80 static void
81 reinit(Matrix &matrix,
82 Sparsity &,
83 int level,
84 const SparsityPatternType & sp,
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
107 reinit(Matrix &matrix,
108 Sparsity &,
109 int level,
110 const SparsityPatternType & sp,
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
135 reinit(Matrix &matrix,
136 Sparsity &,
137 int level,
138 const SparsityPatternType & sp,
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
162 reinit(Matrix &matrix,
163 Sparsity &,
164 int level,
165 const SparsityPatternType & sp,
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 {
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 {
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
215
216 template <typename SparsityPatternType, int dim, int spacedim>
217 static void
218 reinit(Matrix &matrix,
219 Sparsity &,
220 int level,
221 const SparsityPatternType & sp,
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
249template <typename VectorType>
250class MGLevelGlobalTransfer : public MGTransferBase<VectorType>
251{
252public:
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
325protected:
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
392private:
396 template <int dim, int spacedim>
397 void
398 assert_built(const DoFHandler<dim, spacedim> &dof_handler) const;
399};
400
401
402
411template <typename Number>
412class MGLevelGlobalTransfer<LinearAlgebra::distributed::Vector<Number>>
413 : public MGTransferBase<LinearAlgebra::distributed::Vector<Number>>
414{
415public:
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
490protected:
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
584 const MGConstrainedDoFs,
587
594
600
607
613
617 std::function<void(const unsigned int,
620
621private:
625 template <int dim, int spacedim>
626 void
627 assert_built(const DoFHandler<dim, spacedim> &dof_handler) const;
628};
629
630
631
642template <typename VectorType>
644{
645public:
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
741private:
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
void copy_from_mg(const DoFHandler< dim, spacedim > &dof_handler, LinearAlgebra::distributed::Vector< Number2 > &dst, const MGLevelObject< LinearAlgebra::distributed::Vector< Number > > &src) const
SmartPointer< const MGConstrainedDoFs, MGLevelGlobalTransfer< LinearAlgebra::distributed::Vector< Number > > > mg_constrained_dofs
std::vector< Table< 2, unsigned int > > solution_copy_indices_level_mine
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
void set_component_to_block_map(const std::vector< unsigned int > &map)
void copy_from_mg_add(const DoFHandler< dim, spacedim > &dof_handler, LinearAlgebra::distributed::Vector< Number2 > &dst, const MGLevelObject< LinearAlgebra::distributed::Vector< Number > > &src) const
std::vector< Table< 2, unsigned int > > solution_copy_indices_global_mine
LinearAlgebra::distributed::Vector< Number > solution_ghosted_global_vector
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 assert_built(const DoFHandler< dim, spacedim > &dof_handler) const
std::vector< Table< 2, unsigned int > > copy_indices_global_mine
MGLevelObject< LinearAlgebra::distributed::Vector< Number > > solution_ghosted_level_vector
LinearAlgebra::distributed::Vector< Number > ghosted_global_vector
std::function< void(const unsigned int, LinearAlgebra::distributed::Vector< Number > &)> initialize_dof_vector
std::vector< Table< 2, unsigned int > > copy_indices_level_mine
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
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
void fill_and_communicate_copy_indices(const DoFHandler< dim, spacedim > &dof_handler)
std::vector< unsigned int > component_to_block_map
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
std::vector< std::vector< std::pair< types::global_dof_index, types::global_dof_index > > > copy_indices_global_mine
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
virtual ~MGTransferPrebuilt() override=default
std::vector< std::vector< bool > > interface_dofs
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
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
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:472
#define DEAL_II_NAMESPACE_CLOSE
Definition config.h:473
unsigned int level
Definition grid_out.cc:4618
#define DeclException0(Exception0)
Definition exceptions.h:465
static ::ExceptionBase & ExcMatricesNotBuilt()
static ::ExceptionBase & ExcNoProlongation()
static ::ExceptionBase & ExcNotImplemented()
#define AssertThrow(cond, exc)
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)
static void reinit(Matrix &matrix, Sparsity &, int level, const SparsityPatternType &sp, const DoFHandler< dim, spacedim > &dh)
static void reinit(Matrix &matrix, Sparsity &, int level, const SparsityPatternType &sp, const DoFHandler< dim, spacedim > &dh)
static void reinit(Matrix &matrix, Sparsity &, int level, const SparsityPatternType &sp, const DoFHandler< dim, spacedim > &dh)
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