Reference documentation for deal.II version 8.1.0
mg_block_smoother.h
1 // ---------------------------------------------------------------------
2 // @f$Id: mg_block_smoother.h 30036 2013-07-18 16:55:32Z maier @f$
3 //
4 // Copyright (C) 2005 - 2013 by the deal.II authors
5 //
6 // This file is part of the deal.II library.
7 //
8 // The deal.II library is free software; you can use it, redistribute
9 // it, and/or modify it under the terms of the GNU Lesser General
10 // Public License as published by the Free Software Foundation; either
11 // version 2.1 of the License, or (at your option) any later version.
12 // The full text of the license can be found in the file LICENSE at
13 // the top level of the deal.II distribution.
14 //
15 // ---------------------------------------------------------------------
16 
17 #ifndef __deal2__mg_block_smoother_h
18 #define __deal2__mg_block_smoother_h
19 
20 
21 #include <deal.II/base/config.h>
22 #include <deal.II/base/smartpointer.h>
23 #include <deal.II/lac/pointer_matrix.h>
24 #include <deal.II/lac/vector_memory.h>
25 #include <deal.II/lac/block_vector.h>
26 #include <deal.II/multigrid/mg_base.h>
27 #include <deal.II/base/mg_level_object.h>
28 #include <vector>
29 
31 
32 /*
33  * MGSmootherBase is defined in mg_base.h
34  */
35 
38 
47 template <class MATRIX, class RELAX, typename number>
49  : public MGSmootherBase<BlockVector<number> >
50 {
51 public:
57  const unsigned int steps = 1,
58  const bool variable = false,
59  const bool symmetric = false,
60  const bool transpose = false,
61  const bool reverse = false);
62 
85  template <class MGMATRIX, class MGRELAX>
86  void initialize (const MGMATRIX &matrices,
87  const MGRELAX &smoothers);
88 
92  void clear ();
93 
98  void set_steps (const unsigned int);
99 
104  void set_variable (const bool);
105 
110  void set_symmetric (const bool);
111 
117  void set_transpose (const bool);
118 
124  void set_reverse (const bool);
125 
136  virtual void smooth (const unsigned int level,
138  const BlockVector<number> &rhs) const;
139 private:
144 
149 
153  unsigned int steps;
154 
158  bool variable;
159 
163  bool symmetric;
164 
165  /*
166  * Transposed?
167  */
168  bool transpose;
169 
173  bool reverse;
174 
179 
180 };
181 
184 //---------------------------------------------------------------------------
185 
186 #ifndef DOXYGEN
187 
188 template <class MATRIX, class RELAX, typename number>
189 inline
192  const unsigned int steps,
193  const bool variable,
194  const bool symmetric,
195  const bool transpose,
196  const bool reverse)
197  :
198  steps(steps),
199  variable(variable),
200  symmetric(symmetric),
201  transpose(transpose),
202  reverse(reverse),
203  mem(mem)
204 {}
205 
206 
207 template <class MATRIX, class RELAX, typename number>
208 inline void
210 {
211  unsigned int i=matrices.min_level(),
212  max_level=matrices.max_level();
213  for (; i<=max_level; ++i)
214  {
215  smoothers[i] = 0;
216  matrices[i] = 0;
217  }
218 }
219 
220 
221 template <class MATRIX, class RELAX, typename number>
222 template <class MGMATRIX, class MGRELAX>
223 inline void
225  const MGMATRIX &m,
226  const MGRELAX &s)
227 {
228  const unsigned int min = m.min_level();
229  const unsigned int max = m.max_level();
230 
231  matrices.resize(min, max);
232  smoothers.resize(min, max);
233 
234  for (unsigned int i=min; i<=max; ++i)
235  {
236  matrices[i] = &m[i];
237  smoothers[i] = &s[i];
238  }
239 }
240 
241 template <class MATRIX, class RELAX, typename number>
242 inline void
244 set_steps (const unsigned int s)
245 {
246  steps = s;
247 }
248 
249 
250 template <class MATRIX, class RELAX, typename number>
251 inline void
253 set_variable (const bool flag)
254 {
255  variable = flag;
256 }
257 
258 
259 template <class MATRIX, class RELAX, typename number>
260 inline void
262 set_symmetric (const bool flag)
263 {
264  symmetric = flag;
265 }
266 
267 
268 template <class MATRIX, class RELAX, typename number>
269 inline void
271 set_transpose (const bool flag)
272 {
273  transpose = flag;
274 }
275 
276 
277 template <class MATRIX, class RELAX, typename number>
278 inline void
280 set_reverse (const bool flag)
281 {
282  reverse = flag;
283 }
284 
285 
286 template <class MATRIX, class RELAX, typename number>
287 inline void
289  const unsigned int level,
291  const BlockVector<number> &rhs) const
292 {
293  deallog.push("Smooth");
294 
295  unsigned int maxlevel = matrices.max_level();
296  unsigned int steps2 = steps;
297 
298  if (variable)
299  steps2 *= (1<<(maxlevel-level));
300 
301  BlockVector<number> *r = mem.alloc();
302  BlockVector<number> *d = mem.alloc();
303  r->reinit(u);
304  d->reinit(u);
305 
306  bool T = transpose;
307  if (symmetric && (steps2 % 2 == 0))
308  T = false;
309 // cerr << 'S' << level;
310 // cerr << '(' << matrices[level]->m() << ',' << matrices[level]->n() << ')';
311 
312  for (unsigned int i=0; i<steps2; ++i)
313  {
314  if (T)
315  {
316 // cerr << 'T';
317  matrices[level].vmult(*r,u);
318  r->sadd(-1.,1.,rhs);
319  smoothers[level].Tvmult(*d, *r);
320  }
321  else
322  {
323 // cerr << 'N';
324  matrices[level].vmult(*r,u);
325  r->sadd(-1.,1.,rhs);
326  smoothers[level].vmult(*d, *r);
327  }
328 // cerr << '{' << r->l2_norm() << '}';
329  u += *d;
330  if (symmetric)
331  T = !T;
332  }
333 
334  mem.free(r);
335  mem.free(d);
336  deallog.pop();
337 }
338 
339 #endif // DOXYGEN
340 
341 DEAL_II_NAMESPACE_CLOSE
342 
343 #endif
void pop()
void set_transpose(const bool)
void set_steps(const unsigned int)
void reinit(const unsigned int num_blocks, const size_type block_size=0, const bool fast=false)
VectorizedArray< Number > min(const ::VectorizedArray< Number > &x, const ::VectorizedArray< Number > &y)
SymmetricTensor< rank, dim, Number > transpose(const SymmetricTensor< rank, dim, Number > &t)
VectorMemory< BlockVector< number > > & mem
void initialize(const MGMATRIX &matrices, const MGRELAX &smoothers)
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)
void set_reverse(const bool)
void sadd(const value_type s, const BlockVectorBase &V)
void set_variable(const bool)
void push(const std::string &text)
virtual void smooth(const unsigned int level, BlockVector< number > &u, const BlockVector< number > &rhs) const
unsigned int steps
MGLevelObject< PointerMatrix< RELAX, BlockVector< number > > > smoothers
void set_symmetric(const bool)
VectorizedArray< Number > max(const ::VectorizedArray< Number > &x, const ::VectorizedArray< Number > &y)
MGLevelObject< PointerMatrix< MATRIX, BlockVector< number > > > matrices