Loading [MathJax]/extensions/TeX/newcommand.js
 deal.II version GIT relicensing-3065-g80d88f1806 2025-04-14 16:50:00+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\}}
All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Modules Pages Concepts
precondition_selector.h
Go to the documentation of this file.
1// ------------------------------------------------------------------------
2//
3// SPDX-License-Identifier: LGPL-2.1-or-later
4// Copyright (C) 1999 - 2024 by the deal.II authors
5//
6// This file is part of the deal.II library.
7//
8// Part of the source code is dual licensed under Apache-2.0 WITH
9// LLVM-exception OR LGPL-2.1-or-later. Detailed license information
10// governing the source code and code contributions can be found in
11// LICENSE.md and CONTRIBUTING.md at the top level directory of deal.II.
12//
13// ------------------------------------------------------------------------
14
15#ifndef dealii_precondition_selector_h
16#define dealii_precondition_selector_h
17
18
19#include <deal.II/base/config.h>
20
24
25#include <string>
26
28
29// Forward declarations
30#ifndef DOXYGEN
31template <class number>
32class Vector;
33template <class number>
34class SparseMatrix;
35#endif
36
37
99template <typename MatrixType = SparseMatrix<double>,
100 typename VectorType = ::Vector<double>>
102{
103public:
107 using size_type = typename MatrixType::size_type;
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
133 m() const;
134
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
181protected:
185 std::string preconditioning;
186
187private:
192 ObserverPointer<const MatrixType,
195
199 const typename VectorType::value_type omega;
200};
201
203/* --------------------- Inline and template functions ------------------- */
204
205
206template <typename MatrixType, typename VectorType>
208 const std::string &preconditioning,
209 const typename VectorType::value_type &omega)
210 : preconditioning(preconditioning)
211 , omega(omega)
212{}
213
214
215template <typename MatrixType, typename VectorType>
217{
218 // release the matrix A
219 A = nullptr;
220}
221
222
223template <typename MatrixType, typename VectorType>
224void
226{
227 A = &M;
228}
229
230
231template <typename MatrixType, typename VectorType>
234{
235 Assert(A != nullptr, ExcNoMatrixGivenToUse());
236 return A->m();
237}
238
239
240template <typename MatrixType, typename VectorType>
243{
244 Assert(A != nullptr, ExcNoMatrixGivenToUse());
245 return A->n();
246}
247
248
249
250template <typename MatrixType, typename VectorType>
251void
253 const VectorType &src) const
254{
255 if (preconditioning == "none")
256 {
257 dst = src;
258 }
259 else
260 {
261 Assert(A != nullptr, ExcNoMatrixGivenToUse());
262
263 if (preconditioning == "jacobi")
264 {
265 A->precondition_Jacobi(dst, src, omega);
266 }
267 else if (preconditioning == "sor")
268 {
269 A->precondition_SOR(dst, src, omega);
270 }
271 else if (preconditioning == "ssor")
272 {
273 A->precondition_SSOR(dst, src, omega);
274 }
275 else
277 }
278}
279
280
281template <typename MatrixType, typename VectorType>
282void
284 VectorType &dst,
285 const VectorType &src) const
286{
287 if (preconditioning == "none")
288 {
289 dst = src;
290 }
291 else
292 {
293 Assert(A != nullptr, ExcNoMatrixGivenToUse());
294
295 if (preconditioning == "jacobi")
296 {
297 A->precondition_Jacobi(dst, src, omega); // Symmetric operation
298 }
299 else if (preconditioning == "sor")
300 {
301 A->precondition_TSOR(dst, src, omega);
302 }
303 else if (preconditioning == "ssor")
304 {
305 A->precondition_SSOR(dst, src, omega); // Symmetric operation
306 }
307 else
309 }
310}
311
312
313template <typename MatrixType, typename VectorType>
314std::string
316{
317 return "none|jacobi|sor|ssor";
318}
319
320
322
323#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.)
ObserverPointer< const MatrixType, PreconditionSelector< MatrixType, VectorType > > A
#define DEAL_II_NAMESPACE_OPEN
Definition config.h:35
#define DEAL_II_NAMESPACE_CLOSE
Definition config.h:36
#define DEAL_II_NOT_IMPLEMENTED()
static ::ExceptionBase & ExcNoMatrixGivenToUse()
#define DeclException0(Exception0)
#define Assert(cond, exc)