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_block_smoother.h
Go to the documentation of this file.
1// ---------------------------------------------------------------------
2//
3// Copyright (C) 2005 - 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_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
50template <typename MatrixType, class RelaxationType, typename number>
51class MGSmootherBlock : public MGSmoother<BlockVector<number>>
52{
53public:
57 MGSmootherBlock(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
76 template <class MGMatrixType, class MGRelaxationType>
77 void
78 initialize(const MGMatrixType &matrices, const MGRelaxationType &smoothers);
79
83 void
85
89 void
90 set_reverse(const bool);
91
97 virtual void
98 smooth(const unsigned int level,
100 const BlockVector<number> &rhs) const;
101
105 std::size_t
107
108private:
113
118
123
130};
131
134//---------------------------------------------------------------------------
135
136#ifndef DOXYGEN
137
138template <typename MatrixType, class RelaxationType, typename number>
140 const unsigned int steps,
141 const bool variable,
142 const bool symmetric,
143 const bool transpose,
144 const bool reverse)
145 : MGSmoother<BlockVector<number>>(steps, variable, symmetric, transpose)
146 , reverse(reverse)
147 , mem(&this->vector_memory)
148{}
149
150
151template <typename MatrixType, class RelaxationType, typename number>
152inline void
154{
155 unsigned int i = matrices.min_level(), max_level = matrices.max_level();
156 for (; i <= max_level; ++i)
157 {
158 smoothers[i] = LinearOperator<BlockVector<number>>();
159 matrices[i] = LinearOperator<BlockVector<number>>();
160 }
161}
162
163
164template <typename MatrixType, class RelaxationType, typename number>
165template <class MGMatrixType, class MGRelaxationType>
166inline void
168 const MGMatrixType & m,
169 const MGRelaxationType &s)
170{
171 const unsigned int min = m.min_level();
172 const unsigned int max = m.max_level();
173
174 matrices.resize(min, max);
175 smoothers.resize(min, max);
176
177 for (unsigned int i = min; i <= max; ++i)
178 {
179 // Workaround: Unfortunately, not every "m[i]" object has a
180 // rich enough interface to populate reinit_(domain|range)_vector.
181 // Thus, apply an empty LinearOperator exemplar.
182 matrices[i] = linear_operator<BlockVector<number>>(
184 smoothers[i] = linear_operator<BlockVector<number>>(matrices[i], s[i]);
185 }
186}
187
188
189template <typename MatrixType, class RelaxationType, typename number>
190inline void
192 const bool flag)
193{
194 reverse = flag;
195}
196
197
198template <typename MatrixType, class RelaxationType, typename number>
199inline std::size_t
201{
202 return sizeof(*this) + matrices.memory_consumption() +
203 smoothers.memory_consumption() +
204 this->vector_memory.memory_consumption();
205}
206
207
208template <typename MatrixType, class RelaxationType, typename number>
209inline void
211 const unsigned int level,
213 const BlockVector<number> &rhs) const
214{
215 LogStream::Prefix prefix("Smooth");
216
217 unsigned int maxlevel = matrices.max_level();
218 unsigned int steps2 = this->steps;
219
220 if (this->variable)
221 steps2 *= (1 << (maxlevel - level));
222
223 typename VectorMemory<BlockVector<number>>::Pointer r(*this->mem);
224 typename VectorMemory<BlockVector<number>>::Pointer d(*this->mem);
225 r->reinit(u);
226 d->reinit(u);
227
228 bool T = this->transpose;
229 if (this->symmetric && (steps2 % 2 == 0))
230 T = false;
231
232 for (unsigned int i = 0; i < steps2; ++i)
233 {
234 if (T)
235 {
236 matrices[level].vmult(*r, u);
237 r->sadd(-1., 1., rhs);
238 smoothers[level].Tvmult(*d, *r);
239 }
240 else
241 {
242 matrices[level].vmult(*r, u);
243 r->sadd(-1., 1., rhs);
244 smoothers[level].vmult(*d, *r);
245 }
246 u += *d;
247 if (this->symmetric)
248 T = !T;
249 }
250}
251
252#endif // DOXYGEN
253
255
256#endif
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)
#define DEAL_II_NAMESPACE_OPEN
Definition config.h:472
#define DEAL_II_NAMESPACE_CLOSE
Definition config.h:473
DerivativeForm< 1, spacedim, dim, Number > transpose(const DerivativeForm< 1, dim, spacedim, Number > &DF)
unsigned int level
Definition grid_out.cc:4618
SymmetricTensor< 2, dim, Number > d(const Tensor< 2, dim, Number > &F, const Tensor< 2, dim, Number > &dF_dt)