Reference documentation for deal.II version 9.1.1
\(\newcommand{\dealcoloneq}{\mathrel{\vcenter{:}}=}\)
mg_block_smoother.h
1 // ---------------------------------------------------------------------
2 //
3 // Copyright (C) 2005 - 2018 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 
22 #include <deal.II/base/mg_level_object.h>
23 #include <deal.II/base/smartpointer.h>
24 
25 #include <deal.II/lac/block_vector.h>
26 #include <deal.II/lac/linear_operator.h>
27 #include <deal.II/lac/vector_memory.h>
28 
29 #include <deal.II/multigrid/mg_smoother.h>
30 
31 #include <vector>
32 
33 DEAL_II_NAMESPACE_OPEN
34 
35 /*
36  * MGSmootherBase is defined in mg_base.h
37  */
38 
41 
50 template <typename MatrixType, class RelaxationType, typename number>
51 class MGSmootherBlock : public MGSmoother<BlockVector<number>>
52 {
53 public:
60  DEAL_II_DEPRECATED
62  const unsigned int steps = 1,
63  const bool variable = false,
64  const bool symmetric = false,
65  const bool transpose = false,
66  const bool reverse = false);
67 
71  MGSmootherBlock(const unsigned int steps = 1,
72  const bool variable = false,
73  const bool symmetric = false,
74  const bool transpose = false,
75  const bool reverse = false);
76 
90  template <class MGMatrixType, class MGRelaxationType>
91  void
92  initialize(const MGMatrixType &matrices, const MGRelaxationType &smoothers);
93 
97  void
98  clear();
99 
103  void
104  set_reverse(const bool);
105 
111  virtual void
112  smooth(const unsigned int level,
114  const BlockVector<number> &rhs) const;
115 
119  std::size_t
120  memory_consumption() const;
121 
122 private:
127 
132 
136  bool reverse;
137 
144 };
145 
148 //---------------------------------------------------------------------------
149 
150 #ifndef DOXYGEN
151 
152 template <typename MatrixType, class RelaxationType, typename number>
155  const unsigned int steps,
156  const bool variable,
157  const bool symmetric,
158  const bool transpose,
159  const bool reverse)
160  : MGSmoother<BlockVector<number>>(steps, variable, symmetric, transpose)
161  , reverse(reverse)
162  , mem(&mem)
163 {}
164 
165 template <typename MatrixType, class RelaxationType, typename number>
167  const unsigned int steps,
168  const bool variable,
169  const bool symmetric,
170  const bool transpose,
171  const bool reverse)
172  : MGSmoother<BlockVector<number>>(steps, variable, symmetric, transpose)
173  , reverse(reverse)
174  , mem(&this->vector_memory)
175 {}
176 
177 
178 template <typename MatrixType, class RelaxationType, typename number>
179 inline void
181 {
182  unsigned int i = matrices.min_level(), max_level = matrices.max_level();
183  for (; i <= max_level; ++i)
184  {
185  smoothers[i] = LinearOperator<BlockVector<number>>();
186  matrices[i] = LinearOperator<BlockVector<number>>();
187  }
188 }
189 
190 
191 template <typename MatrixType, class RelaxationType, typename number>
192 template <class MGMatrixType, class MGRelaxationType>
193 inline void
195  const MGMatrixType & m,
196  const MGRelaxationType &s)
197 {
198  const unsigned int min = m.min_level();
199  const unsigned int max = m.max_level();
200 
201  matrices.resize(min, max);
202  smoothers.resize(min, max);
203 
204  for (unsigned int i = min; i <= max; ++i)
205  {
206  // Workaround: Unfortunately, not every "m[i]" object has a
207  // rich enough interface to populate reinit_(domain|range)_vector.
208  // Thus, apply an empty LinearOperator exemplar.
209  matrices[i] = linear_operator<BlockVector<number>>(
211  smoothers[i] = linear_operator<BlockVector<number>>(matrices[i], s[i]);
212  }
213 }
214 
215 
216 template <typename MatrixType, class RelaxationType, typename number>
217 inline void
219  const bool flag)
220 {
221  reverse = flag;
222 }
223 
224 
225 template <typename MatrixType, class RelaxationType, typename number>
226 inline std::size_t
228 {
229  return sizeof(*this) + matrices.memory_consumption() +
230  smoothers.memory_consumption() +
231  this->vector_memory.memory_consumption();
232 }
233 
234 
235 template <typename MatrixType, class RelaxationType, typename number>
236 inline void
238  const unsigned int level,
240  const BlockVector<number> &rhs) const
241 {
242  LogStream::Prefix prefix("Smooth");
243 
244  unsigned int maxlevel = matrices.max_level();
245  unsigned int steps2 = this->steps;
246 
247  if (this->variable)
248  steps2 *= (1 << (maxlevel - level));
249 
250  typename VectorMemory<BlockVector<number>>::Pointer r(*this->mem);
251  typename VectorMemory<BlockVector<number>>::Pointer d(*this->mem);
252  r->reinit(u);
253  d->reinit(u);
254 
255  bool T = this->transpose;
256  if (this->symmetric && (steps2 % 2 == 0))
257  T = false;
258 
259  for (unsigned int i = 0; i < steps2; ++i)
260  {
261  if (T)
262  {
263  matrices[level].vmult(*r, u);
264  r->sadd(-1., 1., rhs);
265  smoothers[level].Tvmult(*d, *r);
266  }
267  else
268  {
269  matrices[level].vmult(*r, u);
270  r->sadd(-1., 1., rhs);
271  smoothers[level].vmult(*d, *r);
272  }
273  u += *d;
274  if (this->symmetric)
275  T = !T;
276  }
277 }
278 
279 #endif // DOXYGEN
280 
281 DEAL_II_NAMESPACE_CLOSE
282 
283 #endif
Matrix is symmetric.
void initialize(const MGMatrixType &matrices, const MGRelaxationType &smoothers)
MGLevelObject< LinearOperator< BlockVector< number > > > smoothers
SmartPointer< VectorMemory< BlockVector< number > >, MGSmootherBlock< MatrixType, RelaxationType, number > > mem
void set_reverse(const bool)
virtual void smooth(const unsigned int level, BlockVector< number > &u, const BlockVector< number > &rhs) const
MGSmootherBlock(VectorMemory< BlockVector< number >> &mem, const unsigned int steps=1, const bool variable=false, const bool symmetric=false, const bool transpose=false, const bool reverse=false)
VectorizedArray< Number > min(const ::VectorizedArray< Number > &x, const ::VectorizedArray< Number > &y)
SymmetricTensor< 2, dim, Number > d(const Tensor< 2, dim, Number > &F, const Tensor< 2, dim, Number > &dF_dt)
MGLevelObject< LinearOperator< BlockVector< number > > > matrices
std::size_t memory_consumption() const
DerivativeForm< 1, spacedim, dim, Number > transpose(const DerivativeForm< 1, dim, spacedim, Number > &DF)
VectorizedArray< Number > max(const ::VectorizedArray< Number > &x, const ::VectorizedArray< Number > &y)