Reference documentation for deal.II version 9.3.3
\(\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
48template <typename MatrixType, class RelaxationType, typename number>
49class MGSmootherBlock : public MGSmoother<BlockVector<number>>
50{
51public:
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
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
106private:
111
116
121
128};
129
132//---------------------------------------------------------------------------
133
134#ifndef DOXYGEN
135
136template <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
149template <typename MatrixType, class RelaxationType, typename number>
150inline 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
162template <typename MatrixType, class RelaxationType, typename number>
163template <class MGMatrixType, class MGRelaxationType>
164inline 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
187template <typename MatrixType, class RelaxationType, typename number>
188inline void
190 const bool flag)
191{
192 reverse = flag;
193}
194
195
196template <typename MatrixType, class RelaxationType, typename number>
197inline std::size_t
199{
200 return sizeof(*this) + matrices.memory_consumption() +
201 smoothers.memory_consumption() +
202 this->vector_memory.memory_consumption();
203}
204
205
206template <typename MatrixType, class RelaxationType, typename number>
207inline 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:402
#define DEAL_II_NAMESPACE_CLOSE
Definition: config.h:403
unsigned int level
Definition: grid_out.cc:4590
@ 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)