Ginkgo Generated from branch based on master. Ginkgo version 1.7.0
A numerical linear algebra library targeting many-core architectures
Loading...
Searching...
No Matches
gmres.hpp
1/*******************************<GINKGO LICENSE>******************************
2Copyright (c) 2017-2023, the Ginkgo authors
3All rights reserved.
4
5Redistribution and use in source and binary forms, with or without
6modification, are permitted provided that the following conditions
7are met:
8
91. Redistributions of source code must retain the above copyright
10notice, this list of conditions and the following disclaimer.
11
122. Redistributions in binary form must reproduce the above copyright
13notice, this list of conditions and the following disclaimer in the
14documentation and/or other materials provided with the distribution.
15
163. Neither the name of the copyright holder nor the names of its
17contributors may be used to endorse or promote products derived from
18this software without specific prior written permission.
19
20THIS SOFTWARE IS PROVIDED BY THE COPYRIGHT HOLDERS AND CONTRIBUTORS "AS
21IS" AND ANY EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT LIMITED
22TO, THE IMPLIED WARRANTIES OF MERCHANTABILITY AND FITNESS FOR A
23PARTICULAR PURPOSE ARE DISCLAIMED. IN NO EVENT SHALL THE COPYRIGHT
24HOLDER OR CONTRIBUTORS BE LIABLE FOR ANY DIRECT, INDIRECT, INCIDENTAL,
25SPECIAL, EXEMPLARY, OR CONSEQUENTIAL DAMAGES (INCLUDING, BUT NOT
26LIMITED TO, PROCUREMENT OF SUBSTITUTE GOODS OR SERVICES; LOSS OF USE,
27DATA, OR PROFITS; OR BUSINESS INTERRUPTION) HOWEVER CAUSED AND ON ANY
28THEORY OF LIABILITY, WHETHER IN CONTRACT, STRICT LIABILITY, OR TORT
29(INCLUDING NEGLIGENCE OR OTHERWISE) ARISING IN ANY WAY OUT OF THE USE
30OF THIS SOFTWARE, EVEN IF ADVISED OF THE POSSIBILITY OF SUCH DAMAGE.
31******************************<GINKGO LICENSE>*******************************/
32
33#ifndef GKO_PUBLIC_CORE_SOLVER_GMRES_HPP_
34#define GKO_PUBLIC_CORE_SOLVER_GMRES_HPP_
35
36
37#include <vector>
38
39
40#include <ginkgo/core/base/array.hpp>
41#include <ginkgo/core/base/exception_helpers.hpp>
42#include <ginkgo/core/base/lin_op.hpp>
43#include <ginkgo/core/base/math.hpp>
44#include <ginkgo/core/base/types.hpp>
45#include <ginkgo/core/log/logger.hpp>
46#include <ginkgo/core/matrix/dense.hpp>
47#include <ginkgo/core/matrix/identity.hpp>
48#include <ginkgo/core/solver/solver_base.hpp>
49#include <ginkgo/core/stop/combined.hpp>
50#include <ginkgo/core/stop/criterion.hpp>
51
52
53namespace gko {
54namespace solver {
55
56
57[[deprecated]] constexpr size_type default_krylov_dim = 100u;
58
59constexpr size_type gmres_default_krylov_dim = 100u;
60
61
75template <typename ValueType = default_precision>
76class Gmres
77 : public EnableLinOp<Gmres<ValueType>>,
78 public EnablePreconditionedIterativeSolver<ValueType, Gmres<ValueType>>,
79 public Transposable {
80 friend class EnableLinOp<Gmres>;
81 friend class EnablePolymorphicObject<Gmres, LinOp>;
82
83public:
84 using value_type = ValueType;
86
87 std::unique_ptr<LinOp> transpose() const override;
88
89 std::unique_ptr<LinOp> conj_transpose() const override;
90
96 bool apply_uses_initial_guess() const override { return true; }
97
103 size_type get_krylov_dim() const { return parameters_.krylov_dim; }
104
111
112
113 class Factory;
114
126
127protected:
128 void apply_impl(const LinOp* b, LinOp* x) const override;
129
130 template <typename VectorType>
131 void apply_dense_impl(const VectorType* b, VectorType* x) const;
132
133 void apply_impl(const LinOp* alpha, const LinOp* b, const LinOp* beta,
134 LinOp* x) const override;
135
136 explicit Gmres(std::shared_ptr<const Executor> exec)
137 : EnableLinOp<Gmres>(std::move(exec))
138 {}
139
140 explicit Gmres(const Factory* factory,
141 std::shared_ptr<const LinOp> system_matrix)
142 : EnableLinOp<Gmres>(factory->get_executor(),
143 gko::transpose(system_matrix->get_size())),
144 EnablePreconditionedIterativeSolver<ValueType, Gmres<ValueType>>{
145 std::move(system_matrix), factory->get_parameters()},
146 parameters_{factory->get_parameters()}
147 {
148 if (!parameters_.krylov_dim) {
149 parameters_.krylov_dim = gmres_default_krylov_dim;
150 }
151 }
152};
153
154
155template <typename ValueType>
156struct workspace_traits<Gmres<ValueType>> {
157 using Solver = Gmres<ValueType>;
158 // number of vectors used by this workspace
159 static int num_vectors(const Solver&);
160 // number of arrays used by this workspace
161 static int num_arrays(const Solver&);
162 // array containing the num_vectors names for the workspace vectors
163 static std::vector<std::string> op_names(const Solver&);
164 // array containing the num_arrays names for the workspace vectors
165 static std::vector<std::string> array_names(const Solver&);
166 // array containing all varying scalar vectors (independent of problem size)
167 static std::vector<int> scalars(const Solver&);
168 // array containing all varying vectors (dependent on problem size)
169 static std::vector<int> vectors(const Solver&);
170
171 // residual vector
172 constexpr static int residual = 0;
173 // preconditioned vector
174 constexpr static int preconditioned_vector = 1;
175 // krylov basis multivector
176 constexpr static int krylov_bases = 2;
177 // hessenberg matrix
178 constexpr static int hessenberg = 3;
179 // givens sin parameters
180 constexpr static int givens_sin = 4;
181 // givens cos parameters
182 constexpr static int givens_cos = 5;
183 // coefficients of the residual in Krylov space
184 constexpr static int residual_norm_collection = 6;
185 // residual norm scalar
186 constexpr static int residual_norm = 7;
187 // solution of the least-squares problem in Krylov space
188 constexpr static int y = 8;
189 // solution of the least-squares problem mapped to the full space
190 constexpr static int before_preconditioner = 9;
191 // preconditioned solution of the least-squares problem
192 constexpr static int after_preconditioner = 10;
193 // constant 1.0 scalar
194 constexpr static int one = 11;
195 // constant -1.0 scalar
196 constexpr static int minus_one = 12;
197 // temporary norm vector of next_krylov to copy into hessenberg matrix
198 constexpr static int next_krylov_norm_tmp = 13;
199 // preconditioned krylov basis multivector
200 constexpr static int preconditioned_krylov_bases = 14;
201
202 // stopping status array
203 constexpr static int stop = 0;
204 // reduction tmp array
205 constexpr static int tmp = 1;
206 // final iteration number array
207 constexpr static int final_iter_nums = 2;
208};
209
210
211} // namespace solver
212} // namespace gko
213
214
215#endif // GKO_PUBLIC_CORE_SOLVER_GMRES_HPP_
The EnableLinOp mixin can be used to provide sensible default implementations of the majority of the ...
Definition lin_op.hpp:908
This mixin inherits from (a subclass of) PolymorphicObject and provides a base implementation of a ne...
Definition polymorphic_object.hpp:691
Definition lin_op.hpp:146
std::shared_ptr< const Executor > get_executor() const noexcept
Returns the Executor of the object.
Definition polymorphic_object.hpp:263
Linear operators which support transposition should implement the Transposable interface.
Definition lin_op.hpp:462
A LinOp implementing this interface stores a system matrix and stopping criterion factory.
Definition solver_base.hpp:816
Definition gmres.hpp:124
GMRES or the generalized minimal residual method is an iterative type Krylov subspace method which is...
Definition gmres.hpp:79
std::unique_ptr< LinOp > conj_transpose() const override
Returns a LinOp representing the conjugate transpose of the Transposable object.
std::unique_ptr< LinOp > transpose() const override
Returns a LinOp representing the transpose of the Transposable object.
size_type get_krylov_dim() const
Gets the Krylov dimension of the solver.
Definition gmres.hpp:103
bool apply_uses_initial_guess() const override
Return true as iterative solvers use the data in x as an initial guess.
Definition gmres.hpp:96
void set_krylov_dim(size_type other)
Sets the Krylov dimension.
Definition gmres.hpp:110
#define GKO_FACTORY_PARAMETER_SCALAR(_name, _default)
Creates a scalar factory parameter in the factory parameters structure.
Definition abstract_factory.hpp:473
#define GKO_ENABLE_BUILD_METHOD(_factory_name)
Defines a build method for the factory, simplifying its construction by removing the repetitive typin...
Definition abstract_factory.hpp:422
#define GKO_ENABLE_LIN_OP_FACTORY(_lin_op, _parameters_name, _factory_name)
This macro will generate a default implementation of a LinOpFactory for the LinOp subclass it is defi...
Definition lin_op.hpp:1046
The Ginkgo namespace.
Definition abstract_factory.hpp:48
constexpr T one()
Returns the multiplicative identity for T.
Definition math.hpp:803
std::size_t size_type
Integral type used for allocation quantities.
Definition types.hpp:120
Definition gmres.hpp:117
size_type krylov_dim
Krylov subspace dimension/restart value.
Definition gmres.hpp:119
bool flexible
Flexible GMRES.
Definition gmres.hpp:122
Traits class providing information on the type and location of workspace vectors inside a solver.
Definition solver_base.hpp:267