Reference documentation for deal.II version GIT 3d869ba6cd 2022-05-16 18:15: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_block_smoother.h
Go to the documentation of this file.
1 // ---------------------------------------------------------------------
2 //
3 // Copyright (C) 2005 - 2020 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_block_smoother_h
17 #define dealii_mg_block_smoother_h
18 
19 
20 #include <deal.II/base/config.h>
21 
24 
28 
30 
31 #include <vector>
32 
34 
35 /*
36  * MGSmootherBase is defined in mg_base.h
37  */
38 
41 
48 template <typename MatrixType, class RelaxationType, typename number>
49 class MGSmootherBlock : public MGSmoother<BlockVector<number>>
50 {
51 public:
55  MGSmootherBlock(const unsigned int steps = 1,
56  const bool variable = false,
57  const bool symmetric = false,
58  const bool transpose = false,
59  const bool reverse = false);
60 
74  template <class MGMatrixType, class MGRelaxationType>
75  void
76  initialize(const MGMatrixType &matrices, const MGRelaxationType &smoothers);
77 
81  void
82  clear();
83 
87  void
88  set_reverse(const bool);
89 
95  virtual void
96  smooth(const unsigned int level,
98  const BlockVector<number> &rhs) const;
99 
103  std::size_t
105 
106 private:
111 
116 
120  bool reverse;
121 
128 };
129 
132 //---------------------------------------------------------------------------
133 
134 #ifndef DOXYGEN
135 
136 template <typename MatrixType, class RelaxationType, typename number>
138  const unsigned int steps,
139  const bool variable,
140  const bool symmetric,
141  const bool transpose,
142  const bool reverse)
143  : MGSmoother<BlockVector<number>>(steps, variable, symmetric, transpose)
144  , reverse(reverse)
145  , mem(&this->vector_memory)
146 {}
147 
148 
149 template <typename MatrixType, class RelaxationType, typename number>
150 inline void
152 {
153  unsigned int i = matrices.min_level(), max_level = matrices.max_level();
154  for (; i <= max_level; ++i)
155  {
156  smoothers[i] = LinearOperator<BlockVector<number>>();
157  matrices[i] = LinearOperator<BlockVector<number>>();
158  }
159 }
160 
161 
162 template <typename MatrixType, class RelaxationType, typename number>
163 template <class MGMatrixType, class MGRelaxationType>
164 inline void
166  const MGMatrixType & m,
167  const MGRelaxationType &s)
168 {
169  const unsigned int min = m.min_level();
170  const unsigned int max = m.max_level();
171 
172  matrices.resize(min, max);
173  smoothers.resize(min, max);
174 
175  for (unsigned int i = min; i <= max; ++i)
176  {
177  // Workaround: Unfortunately, not every "m[i]" object has a
178  // rich enough interface to populate reinit_(domain|range)_vector.
179  // Thus, apply an empty LinearOperator exemplar.
180  matrices[i] = linear_operator<BlockVector<number>>(
182  smoothers[i] = linear_operator<BlockVector<number>>(matrices[i], s[i]);
183  }
184 }
185 
186 
187 template <typename MatrixType, class RelaxationType, typename number>
188 inline void
190  const bool flag)
191 {
192  reverse = flag;
193 }
194 
195 
196 template <typename MatrixType, class RelaxationType, typename number>
197 inline std::size_t
199 {
200  return sizeof(*this) + matrices.memory_consumption() +
201  smoothers.memory_consumption() +
202  this->vector_memory.memory_consumption();
203 }
204 
205 
206 template <typename MatrixType, class RelaxationType, typename number>
207 inline void
209  const unsigned int level,
211  const BlockVector<number> &rhs) const
212 {
213  LogStream::Prefix prefix("Smooth");
214 
215  unsigned int maxlevel = matrices.max_level();
216  unsigned int steps2 = this->steps;
217 
218  if (this->variable)
219  steps2 *= (1 << (maxlevel - level));
220 
221  typename VectorMemory<BlockVector<number>>::Pointer r(*this->mem);
222  typename VectorMemory<BlockVector<number>>::Pointer d(*this->mem);
223  r->reinit(u);
224  d->reinit(u);
225 
226  bool T = this->transpose;
227  if (this->symmetric && (steps2 % 2 == 0))
228  T = false;
229 
230  for (unsigned int i = 0; i < steps2; ++i)
231  {
232  if (T)
233  {
234  matrices[level].vmult(*r, u);
235  r->sadd(-1., 1., rhs);
236  smoothers[level].Tvmult(*d, *r);
237  }
238  else
239  {
240  matrices[level].vmult(*r, u);
241  r->sadd(-1., 1., rhs);
242  smoothers[level].vmult(*d, *r);
243  }
244  u += *d;
245  if (this->symmetric)
246  T = !T;
247  }
248 }
249 
250 #endif // DOXYGEN
251 
253 
254 #endif
DerivativeForm< 1, spacedim, dim, Number > transpose(const DerivativeForm< 1, dim, spacedim, Number > &DF)
virtual void smooth(const unsigned int level, BlockVector< number > &u, const BlockVector< number > &rhs) const
MGLevelObject< LinearOperator< BlockVector< number > > > smoothers
std::size_t memory_consumption() const
SmartPointer< VectorMemory< BlockVector< number > >, MGSmootherBlock< MatrixType, RelaxationType, number > > mem
void initialize(const MGMatrixType &matrices, const MGRelaxationType &smoothers)
void set_reverse(const bool)
MGLevelObject< LinearOperator< BlockVector< number > > > matrices
MGSmootherBlock(const unsigned int steps=1, const bool variable=false, const bool symmetric=false, const bool transpose=false, const bool reverse=false)
VectorizedArray< Number, width > max(const ::VectorizedArray< Number, width > &x, const ::VectorizedArray< Number, width > &y)
VectorizedArray< Number, width > min(const ::VectorizedArray< Number, width > &x, const ::VectorizedArray< Number, width > &y)
#define DEAL_II_NAMESPACE_OPEN
Definition: config.h:416
#define DEAL_II_NAMESPACE_CLOSE
Definition: config.h:417
unsigned int level
Definition: grid_out.cc:4606
@ symmetric
Matrix is symmetric.
static const char T
SymmetricTensor< 2, dim, Number > d(const Tensor< 2, dim, Number > &F, const Tensor< 2, dim, Number > &dF_dt)