Reference documentation for deal.II version GIT relicensing-233-g802318d791 2024-03-28 20:20:02+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\}}\)
Loading...
Searching...
No Matches
solver_relaxation.h
Go to the documentation of this file.
1// ------------------------------------------------------------------------
2//
3// SPDX-License-Identifier: LGPL-2.1-or-later
4// Copyright (C) 2010 - 2023 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_solver_relaxation_h
16#define dealii_solver_relaxation_h
17
18
19#include <deal.II/base/config.h>
20
23
24#include <deal.II/lac/solver.h>
26
28
55template <typename VectorType = Vector<double>>
56class SolverRelaxation : public SolverBase<VectorType>
57{
58public:
64 {};
65
70 const AdditionalData &data = AdditionalData());
71
77 template <typename MatrixType, typename RelaxationType>
78 void
79 solve(const MatrixType &A,
80 VectorType &x,
81 const VectorType &b,
82 const RelaxationType &R);
83};
84
85//----------------------------------------------------------------------//
86
87template <typename VectorType>
92
93
94
95template <typename VectorType>
96template <typename MatrixType, typename RelaxationType>
97void
99 VectorType &x,
100 const VectorType &b,
101 const RelaxationType &R)
102{
105
106 // Memory allocation
107 typename VectorMemory<VectorType>::Pointer Vr(mem);
108 VectorType &r = *Vr;
109 r.reinit(x);
110 typename VectorMemory<VectorType>::Pointer Vd(mem);
111 VectorType &d = *Vd;
112 d.reinit(x);
113
114 LogStream::Prefix prefix("Relaxation");
115
116 int iter = 0;
117 // Main loop
118 for (; conv == SolverControl::iterate; ++iter)
119 {
120 // Compute residual
121 A.vmult(r, x);
122 r.sadd(-1., 1., b);
123
124 // The required norm of the
125 // (preconditioned)
126 // residual is computed in
127 // criterion() and stored
128 // in res.
129 conv = this->iteration_status(iter, r.l2_norm(), x);
130 if (conv != SolverControl::iterate)
131 break;
132 R.step(x, b);
133 }
134
135 // in case of failure: throw exception
137 SolverControl::NoConvergence(iter, r.l2_norm()));
138 // otherwise exit as normal
139}
140
141
143
144#endif
@ iterate
Continue iteration.
@ success
Stop iteration, goal reached.
SolverRelaxation(SolverControl &cn, const AdditionalData &data=AdditionalData())
void solve(const MatrixType &A, VectorType &x, const VectorType &b, const RelaxationType &R)
#define DEAL_II_NAMESPACE_OPEN
Definition config.h:502
#define DEAL_II_NAMESPACE_CLOSE
Definition config.h:503
#define AssertThrow(cond, exc)