Reference documentation for deal.II version GIT 29f9da0a34 2023-12-07 10:00: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\}}\)
precondition_selector.h
Go to the documentation of this file.
1 // ---------------------------------------------------------------------
2 //
3 // Copyright (C) 1999 - 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_precondition_selector_h
17 #define dealii_precondition_selector_h
18 
19 
20 #include <deal.II/base/config.h>
21 
24 
25 #include <string>
26 
28 
29 // Forward declarations
30 #ifndef DOXYGEN
31 template <class number>
32 class Vector;
33 template <class number>
34 class SparseMatrix;
35 #endif
36 
37 
99 template <typename MatrixType = SparseMatrix<double>,
100  typename VectorType = ::Vector<double>>
102 {
103 public:
108 
113  PreconditionSelector(const std::string &preconditioning,
114  const typename VectorType::value_type &omega = 1.);
115 
119  virtual ~PreconditionSelector() override;
120 
125  void
126  use_matrix(const MatrixType &M);
127 
132  size_type
133  m() const;
134 
139  size_type
140  n() const;
141 
146  virtual void
147  vmult(VectorType &dst, const VectorType &src) const;
148 
153  virtual void
154  Tvmult(VectorType &dst, const VectorType &src) const;
155 
166  static std::string
168 
179 
181 protected:
185  std::string preconditioning;
186 
187 private:
193  A;
194 
198  const typename VectorType::value_type omega;
199 };
200 
202 /* --------------------- Inline and template functions ------------------- */
203 
204 
205 template <typename MatrixType, typename VectorType>
207  const std::string &preconditioning,
208  const typename VectorType::value_type &omega)
209  : preconditioning(preconditioning)
210  , omega(omega)
211 {}
212 
213 
214 template <typename MatrixType, typename VectorType>
216 {
217  // release the matrix A
218  A = nullptr;
219 }
220 
221 
222 template <typename MatrixType, typename VectorType>
223 void
225 {
226  A = &M;
227 }
228 
229 
230 template <typename MatrixType, typename VectorType>
233 {
234  Assert(A != nullptr, ExcNoMatrixGivenToUse());
235  return A->m();
236 }
237 
238 
239 template <typename MatrixType, typename VectorType>
242 {
243  Assert(A != nullptr, ExcNoMatrixGivenToUse());
244  return A->n();
245 }
246 
247 
248 
249 template <typename MatrixType, typename VectorType>
250 void
252  const VectorType &src) const
253 {
254  if (preconditioning == "none")
255  {
256  dst = src;
257  }
258  else
259  {
260  Assert(A != nullptr, ExcNoMatrixGivenToUse());
261 
262  if (preconditioning == "jacobi")
263  {
264  A->precondition_Jacobi(dst, src, omega);
265  }
266  else if (preconditioning == "sor")
267  {
268  A->precondition_SOR(dst, src, omega);
269  }
270  else if (preconditioning == "ssor")
271  {
272  A->precondition_SSOR(dst, src, omega);
273  }
274  else
275  Assert(false, ExcNotImplemented());
276  }
277 }
278 
279 
280 template <typename MatrixType, typename VectorType>
281 void
283  VectorType &dst,
284  const VectorType &src) const
285 {
286  if (preconditioning == "none")
287  {
288  dst = src;
289  }
290  else
291  {
292  Assert(A != nullptr, ExcNoMatrixGivenToUse());
293 
294  if (preconditioning == "jacobi")
295  {
296  A->precondition_Jacobi(dst, src, omega); // Symmetric operation
297  }
298  else if (preconditioning == "sor")
299  {
300  A->precondition_TSOR(dst, src, omega);
301  }
302  else if (preconditioning == "ssor")
303  {
304  A->precondition_SSOR(dst, src, omega); // Symmetric operation
305  }
306  else
307  Assert(false, ExcNotImplemented());
308  }
309 }
310 
311 
312 template <typename MatrixType, typename VectorType>
313 std::string
315 {
316  return "none|jacobi|sor|ssor";
317 }
318 
319 
321 
322 #endif
virtual void vmult(VectorType &dst, const VectorType &src) const
virtual void Tvmult(VectorType &dst, const VectorType &src) const
typename MatrixType::size_type size_type
const VectorType::value_type omega
static std::string get_precondition_names()
void use_matrix(const MatrixType &M)
virtual ~PreconditionSelector() override
PreconditionSelector(const std::string &preconditioning, const typename VectorType::value_type &omega=1.)
SmartPointer< const MatrixType, PreconditionSelector< MatrixType, VectorType > > A
Definition: vector.h:110
#define DEAL_II_NAMESPACE_OPEN
Definition: config.h:477
#define DEAL_II_NAMESPACE_CLOSE
Definition: config.h:478
#define DeclException0(Exception0)
Definition: exceptions.h:472
#define Assert(cond, exc)
Definition: exceptions.h:1631
static ::ExceptionBase & ExcNotImplemented()
static ::ExceptionBase & ExcNoMatrixGivenToUse()
static const char A
types::global_dof_index size_type
Definition: cuda_kernels.h:45