MFEM v4.8.0
Finite element discretization library
Loading...
Searching...
No Matches
ginkgo.hpp
Go to the documentation of this file.
1// Copyright (c) 2010-2025, Lawrence Livermore National Security, LLC. Produced
2// at the Lawrence Livermore National Laboratory. All Rights reserved. See files
3// LICENSE and NOTICE for details. LLNL-CODE-806117.
4//
5// This file is part of the MFEM library. For more information and source code
6// availability visit https://mfem.org.
7//
8// MFEM is free software; you can redistribute it and/or modify it under the
9// terms of the BSD-3 license. We welcome feedback and contributions, see file
10// CONTRIBUTING.md for details.
11
12#ifndef MFEM_GINKGO
13#define MFEM_GINKGO
14
15#include "../config/config.hpp"
16
17#ifdef MFEM_USE_GINKGO
18
19#include "operator.hpp"
20#include "sparsemat.hpp"
21#include "solvers.hpp"
22
23#include <ginkgo/ginkgo.hpp>
24
25#include <iomanip>
26#include <ios>
27#include <string>
28#include <vector>
29#include <fstream>
30#include <iostream>
31
32#define MFEM_GINKGO_VERSION \
33 ((GKO_VERSION_MAJOR*100 + GKO_VERSION_MINOR)*100 + GKO_VERSION_PATCH)
34
35namespace mfem
36{
37namespace Ginkgo
38{
39
40template <typename T> using gko_array = gko::array<T>;
41
42/**
43* Helper class for a case where a wrapped MFEM Vector
44* should be owned by Ginkgo, and deleted when the wrapper
45* object goes out of scope.
46*/
47template <typename T>
49{
50public:
51 using pointer = T *;
52
53 // Destroys an MFEM object. Requires object to have a Destroy() method.
54 void operator()(pointer ptr) const noexcept { ptr->Destroy(); }
55};
56
57/**
58* This class wraps an MFEM vector object for Ginkgo's use. It
59* allows Ginkgo and MFEM to operate directly on the same
60* data, and is necessary to use MFEM Operators with Ginkgo
61* solvers.
62*
63* @ingroup Ginkgo
64*/
65
66class VectorWrapper : public gko::matrix::Dense<real_t>
67{
68public:
69 VectorWrapper(std::shared_ptr<const gko::Executor> exec,
70 gko::size_type size, Vector *mfem_vec,
71 bool ownership = false)
72 : gko::matrix::Dense<real_t>(
73 exec,
74 gko::dim<2> {size, 1},
76 size,
77 mfem_vec->ReadWrite(
78 exec != exec->get_master() ? true : false)),
79 1)
80 {
81 // This controls whether or not we want Ginkgo to own its MFEM Vector.
82 // Normally, when we are wrapping an MFEM Vector created outside
83 // Ginkgo, we do not want ownership to be true. However, Ginkgo
84 // creates its own temporary vectors as part of its solvers, and
85 // these will be owned (and deleted) by Ginkgo.
86 if (ownership)
87 {
88 using deleter = gko_mfem_destroy<Vector>;
89 wrapped_vec = std::unique_ptr<Vector,
90 std::function<void(Vector *)>>(
91 mfem_vec, deleter{});
92 }
93 else
94 {
95 using deleter = gko::null_deleter<Vector>;
96 wrapped_vec = std::unique_ptr<Vector,
97 std::function<void(Vector *)>>(
98 mfem_vec, deleter{});
99 }
100 }
101
102 static std::unique_ptr<VectorWrapper> create(
103 std::shared_ptr<const gko::Executor> exec,
104 gko::size_type size,
105 Vector *mfem_vec,
106 bool ownership = false)
107 {
108 return std::unique_ptr<VectorWrapper>(
109 new VectorWrapper(exec, size, mfem_vec, ownership));
110 }
111
112 // Return reference to MFEM Vector object
113 Vector &get_mfem_vec_ref() { return *(this->wrapped_vec.get()); }
114
115 // Return const reference to MFEM Vector object
116 const Vector &get_mfem_vec_const_ref() const { return *(this->wrapped_vec.get()); }
117
118 // Override base Dense class implementation for creating new vectors
119 // with same executor and size as self
120 std::unique_ptr<gko::matrix::Dense<real_t>>
121 create_with_same_config() const override
122 {
123 Vector *mfem_vec = new Vector(
124 this->get_size()[0],
125 this->wrapped_vec.get()->GetMemory().GetMemoryType());
126
127 mfem_vec->UseDevice(this->wrapped_vec.get()->UseDevice());
128
129 // If this function is called, Ginkgo is creating this
130 // object and should control the memory, so ownership is
131 // set to true
132 return VectorWrapper::create(this->get_executor(),
133 this->get_size()[0],
134 mfem_vec,
135 true);
136 }
137
138 // Override base Dense class implementation for creating new vectors
139 // with same executor and type as self, but with a different size.
140 // This function will create "one large VectorWrapper" of size
141 // size[0] * size[1], since MFEM Vectors only have one dimension.
142 std::unique_ptr<gko::matrix::Dense<real_t>> create_with_type_of_impl(
143 std::shared_ptr<const gko::Executor> exec,
144 const gko::dim<2> &size,
145 gko::size_type stride) const override
146 {
147 // Only stride of 1 is allowed for VectorWrapper type
148 if (stride > 1)
149 {
150 throw gko::Error(
151 __FILE__, __LINE__,
152 "VectorWrapper cannot be created with stride > 1");
153 }
154 // Compute total size of new Vector
155 gko::size_type total_size = size[0]*size[1];
156 Vector *mfem_vec = new Vector(
157 total_size,
158 this->wrapped_vec.get()->GetMemory().GetMemoryType());
159
160 mfem_vec->UseDevice(this->wrapped_vec.get()->UseDevice());
161
162 // If this function is called, Ginkgo is creating this
163 // object and should control the memory, so ownership is
164 // set to true
166 this->get_executor(), total_size, mfem_vec,
167 true);
168 }
169
170 // Override base Dense class implementation for creating new sub-vectors
171 // from a larger vector.
172 std::unique_ptr<gko::matrix::Dense<real_t>> create_submatrix_impl(
173 const gko::span &rows,
174 const gko::span &columns,
175 const gko::size_type stride) override
176 {
177
178 gko::size_type num_rows = rows.end - rows.begin;
179 gko::size_type num_cols = columns.end - columns.begin;
180 // Data in the Dense matrix will be stored in row-major format.
181 // Check that we only have one column, and that the stride = 1
182 // (only allowed value for VectorWrappers).
183 if (num_cols > 1 || stride > 1)
184 {
185 throw gko::BadDimension(
186 __FILE__, __LINE__, __func__, "new_submatrix", num_rows,
187 num_cols,
188 "VectorWrapper submatrix must have one column and stride = 1");
189 }
190 int data_size = static_cast<int>(num_rows * num_cols);
191 int start = static_cast<int>(rows.begin);
192 // Create a new MFEM Vector pointing to this starting point in the data
193 Vector *mfem_vec = new Vector();
194 mfem_vec->MakeRef(*(this->wrapped_vec.get()), start, data_size);
195 mfem_vec->UseDevice(this->wrapped_vec.get()->UseDevice());
196
197 // If this function is called, Ginkgo is creating this
198 // object and should control the memory, so ownership is
199 // set to true (but MFEM doesn't own and won't delete
200 // the data, at it's only a reference to the parent Vector)
202 this->get_executor(), data_size, mfem_vec,
203 true);
204 }
205
206private:
207 std::unique_ptr<Vector, std::function<void(Vector *)>> wrapped_vec;
208};
209
210/**
211* This class wraps an MFEM Operator for Ginkgo, to make its Mult()
212* function available to Ginkgo, provided the input and output vectors
213* are of the VectorWrapper type.
214* Note that this class does NOT take ownership of the MFEM Operator.
215*
216* @ingroup Ginkgo
217*/
219 : public gko::EnableLinOp<OperatorWrapper>,
220 public gko::EnableCreateMethod<OperatorWrapper>
221{
222public:
223 OperatorWrapper(std::shared_ptr<const gko::Executor> exec,
224 gko::size_type size = 0,
225 const Operator *oper = NULL)
226 : gko::EnableLinOp<OperatorWrapper>(exec, gko::dim<2> {size, size}),
227 gko::EnableCreateMethod<OperatorWrapper>()
228 {
229 this->wrapped_oper = oper;
230 }
231
232protected:
233 void apply_impl(const gko::LinOp *b, gko::LinOp *x) const override;
234 void apply_impl(const gko::LinOp *alpha, const gko::LinOp *b,
235 const gko::LinOp *beta, gko::LinOp *x) const override;
236
237private:
238 const Operator *wrapped_oper;
239};
240
241// Utility function which gets the scalar value of a Ginkgo gko::matrix::Dense
242// matrix representing the norm of a vector.
243template <typename ValueType=real_t>
244real_t get_norm(const gko::matrix::Dense<ValueType> *norm)
245{
246 // Put the value on CPU thanks to the master executor
247 auto cpu_norm = clone(norm->get_executor()->get_master(), norm);
248 // Return the scalar value contained at position (0, 0)
249 return cpu_norm->at(0, 0);
250}
251
252// Utility function which computes the norm of a Ginkgo gko::matrix::Dense
253// vector.
254template <typename ValueType=real_t>
255real_t compute_norm(const gko::matrix::Dense<ValueType> *b)
256{
257 // Get the executor of the vector
258 auto exec = b->get_executor();
259 // Initialize a result scalar containing the value 0.0.
260 auto b_norm = gko::initialize<gko::matrix::Dense<ValueType>>({0.0}, exec);
261 // Use the dense `compute_norm2` function to compute the norm.
262 b->compute_norm2(b_norm);
263 // Use the other utility function to return the norm contained in `b_norm``
264 return std::pow(get_norm(b_norm.get()),2);
265}
266
267/**
268 * Custom logger class which intercepts the residual norm scalar and solution
269 * vector in order to print a table of real vs recurrent (internal to the
270 * solvers) residual norms.
271 *
272 * This has been taken from the custom-logger example of Ginkgo. See the
273 * custom-logger example to understand how to write and modify your own loggers
274 * with Ginkgo.
275 *
276 * @ingroup Ginkgo
277 */
278template <typename ValueType=real_t>
279struct ResidualLogger : gko::log::Logger
280{
281 // Output the logger's data in a table format
282 void write() const
283 {
284 // Print a header for the table
285 if (compute_real_residual)
286 {
287 mfem::out << "Iteration log with real residual norms:" << std::endl;
288 }
289 else
290 {
291 mfem::out << "Iteration log with residual norms:" << std::endl;
292 }
293 mfem::out << '|' << std::setw(10) << "Iteration" << '|' << std::setw(25)
294 << "Residual Norm" << '|' << std::endl;
295 // Print a separation line. Note that for creating `10` characters
296 // `std::setw()` should be set to `11`.
297 mfem::out << '|' << std::setfill('-') << std::setw(11) << '|' <<
298 std::setw(26) << '|' << std::setfill(' ') << std::endl;
299 // Print the data one by one in the form
300 mfem::out << std::scientific;
301 for (std::size_t i = 0; i < iterations.size(); i++)
302 {
303 mfem::out << '|' << std::setw(10) << iterations[i] << '|'
304 << std::setw(25) << residual_norms[i] << '|' << std::endl;
305 }
306 // std::defaultfloat could be used here but some compilers do not support
307 // it properly, e.g. the Intel compiler
308 mfem::out.unsetf(std::ios_base::floatfield);
309 // Print a separation line
310 mfem::out << '|' << std::setfill('-') << std::setw(11) << '|' <<
311 std::setw(26) << '|' << std::setfill(' ') << std::endl;
312 }
313
314 using gko_dense = gko::matrix::Dense<ValueType>;
315
316 // Ginkgo 1.5 and older: version for solver that doesn't log implicit res norm
317 void on_iteration_complete(const gko::LinOp *op,
318 const gko::size_type &iteration,
319 const gko::LinOp *residual,
320 const gko::LinOp *solution,
321 const gko::LinOp *residual_norm) const override
322 {
323 iteration_complete_core(iteration, residual, solution, residual_norm,
324 nullptr);
325 }
326 // Ginkgo 1.5 and older: version with implicit residual norm
327 void on_iteration_complete(const gko::LinOp *op,
328 const gko::size_type &iteration,
329 const gko::LinOp *residual,
330 const gko::LinOp *solution,
331 const gko::LinOp *residual_norm,
332 const gko::LinOp *implicit_sq_residual_norm) const override
333 {
334 iteration_complete_core(iteration, residual, solution, residual_norm,
335 implicit_sq_residual_norm);
336 }
337 // Ginkgo 1.6 and newer
338 void on_iteration_complete(const gko::LinOp *op,
339 const gko::LinOp *rhs,
340 const gko::LinOp *solution,
341 const gko::size_type &iteration,
342 const gko::LinOp *residual,
343 const gko::LinOp *residual_norm,
344 const gko::LinOp *implicit_sq_residual_norm,
345 const gko::array<gko::stopping_status>* status,
346 bool stopped) const override
347 {
348 iteration_complete_core(iteration, residual, solution, residual_norm,
349 implicit_sq_residual_norm);
350 }
351
352 // Construct the logger and store the system matrix and b vectors
353 ResidualLogger(std::shared_ptr<const gko::Executor> exec,
354 const gko::LinOp *matrix, const gko_dense *b,
355 bool compute_real_residual=false)
356 :
357 gko::log::Logger(gko::log::Logger::iteration_complete_mask),
358 matrix {matrix},
359 b{b},
360 compute_real_residual{compute_real_residual}
361 {
362 if (compute_real_residual == true)
363 {
364 if (dynamic_cast<const VectorWrapper*>(b))
365 {
366 const VectorWrapper *b_cast = gko::as<const VectorWrapper>(b);
367 res = std::move(gko_dense::create_with_config_of(b_cast).release());
368 }
369 else
370 {
371 res = std::move(gko::clone(b).release());
372 }
373 }
374 }
375
376private:
377 // Customize the logging hook which is called every time an iteration is
378 // completed.
379 void iteration_complete_core(const gko::size_type &iteration,
380 const gko::LinOp *residual,
381 const gko::LinOp *solution,
382 const gko::LinOp *residual_norm,
383 const gko::LinOp *implicit_sq_residual_norm) const
384 {
385 // If the solver shares the current solution vector and we want to
386 // compute the residual from that
387 if (solution && compute_real_residual)
388 {
389 // Store the matrix's executor
390 auto exec = matrix->get_executor();
391 // Compute the real residual vector by calling apply on the system
392 // First, compute res = A * x
393 matrix->apply(solution, res);
394 // Now do res = res - b, depending on which vector/oper type
395 // Check if b is a Ginkgo vector or wrapped MFEM Vector
396 if (dynamic_cast<const VectorWrapper*>(b))
397 {
398 const VectorWrapper *b_cast = gko::as<const VectorWrapper>(b);
399 // Copy the MFEM Vector stored in b
400 VectorWrapper *res_cast = gko::as<VectorWrapper>(res);
401 res_cast->get_mfem_vec_ref() -= b_cast->get_mfem_vec_const_ref();
402 }
403 else
404 {
405 // Create a scalar containing the value -1.0
406 auto neg_one = gko::initialize<gko_dense>({-1.0}, exec);
407 res->add_scaled(neg_one, b);
408 }
409
410 // Compute the norm of the residual vector and add it to the
411 // `residual_norms` vector
412 residual_norms.push_back(compute_norm(res));
413 }
414 else
415 {
416 // If the solver shares an implicit or recurrent residual norm, log its value
417 if (implicit_sq_residual_norm)
418 {
419 auto dense_norm = gko::as<gko_dense>(implicit_sq_residual_norm);
420 // Add the norm to the `residual_norms` vector
421 residual_norms.push_back(get_norm(dense_norm));
422 // Otherwise, use the recurrent residual vector
423 }
424 else if (residual_norm)
425 {
426 auto dense_norm = gko::as<gko_dense>(residual_norm);
427 // Add the norm to the `residual_norms` vector
428 residual_norms.push_back(get_norm(dense_norm));
429 // Otherwise, use the recurrent residual vector
430 }
431 else
432 {
433 auto dense_residual = gko::as<gko_dense>(residual);
434 // Compute the residual vector's norm
435 auto norm = compute_norm(dense_residual);
436 // Add the computed norm to the `residual_norms` vector
437 residual_norms.push_back(norm);
438 }
439 }
440 // Add the current iteration number to the `iterations` vector
441 iterations.push_back(iteration);
442 }
443
444 // Pointer to the system matrix
445 const gko::LinOp *matrix;
446 // Pointer to the right hand sides
447 const gko_dense *b;
448 // Pointer to the residual workspace vector
449 gko_dense *res;
450 // Vector which stores all the residual norms
451 mutable std::vector<ValueType> residual_norms{};
452 // Vector which stores all the iteration numbers
453 mutable std::vector<std::size_t> iterations{};
454 // Whether or not to compute the residual at every iteration,
455 // rather than using the recurrent norm
456 const bool compute_real_residual;
457};
458
459/**
460* This class wraps a Ginkgo Executor for use in MFEM.
461* Note that objects in the Ginkgo namespace intended to work
462* together, e.g. a Ginkgo solver and preconditioner, should use the same
463* GinkgoExecutor object. In general, most users will want to create
464* one GinkgoExecutor object for use with all Ginkgo-related objects.
465* The wrapper can be created to match MFEM's device configuration.
466*/
468{
469public:
470 // Types of Ginkgo Executors.
472 {
473 /// Reference CPU Executor.
475 /// OpenMP CPU Executor.
476 OMP = 1,
477 /// CUDA GPU Executor.
478 CUDA = 2,
479 /// HIP GPU Executor.
480 HIP = 3
481 };
482 /**
483 * Constructor.
484 * Takes an @p GinkgoExecType argument and creates an Executor.
485 * In Ginkgo, GPU Executors must have an associated host Executor.
486 * This routine will select a CPU Executor based on the OpenMP support
487 * for Ginkgo.
488 */
489 GinkgoExecutor(ExecType exec_type);
490
491 /**
492 * Constructor.
493 * Takes an @p GinkgoExecType argument and creates an Executor.
494 * In Ginkgo, GPU Executors must have an associated host Executor.
495 * This routine allows for explicite setting of the CPU Executor
496 * for GPU backends.
497 */
498 GinkgoExecutor(ExecType exec_type, ExecType host_exec_type);
499
500 /**
501 * Constructor.
502 * Takes an MFEM @p Device object and creates an Executor
503 * that "matches" (e.g., if MFEM is using the CPU, Ginkgo
504 * will choose the Reference or OmpExecutor based on MFEM's
505 * configuration and Ginkgo's capabilities; if MFEM is using
506 * CUDA, Ginkgo will choose the CudaExecutor with a default
507 * CPU Executor based on Ginkgo's OpenMP support).
508 */
509 GinkgoExecutor(Device &mfem_device);
510
511 /**
512 * Constructor.
513 * Takes an MFEM @p Device object and creates an Executor
514 * that "matches", but allows the user to specify the host
515 * Executor for GPU backends.
516 */
517 GinkgoExecutor(Device &mfem_device, ExecType host_exec_type);
518
519 /**
520 * Destructor.
521 */
522 virtual ~GinkgoExecutor() = default;
523
524 std::shared_ptr<gko::Executor> GetExecutor() const
525 {
526 return this->executor;
527 };
528
529private:
530 std::shared_ptr<gko::Executor> executor;
531
532};
533
534/**
535* This class forms the base class for all of Ginkgo's preconditioners. The
536* various derived classes only take the additional data that is specific to them.
537* The entire collection of preconditioners that Ginkgo implements is available
538* at the Ginkgo documentation and manual pages,
539* https://ginkgo-project.github.io/ginkgo/doc/develop.
540*
541* @ingroup Ginkgo
542*/
544{
545public:
546 /**
547 * Constructor.
548 *
549 * The @p exec defines the paradigm where the solution is computed.
550 * Ginkgo currently supports four different executor types:
551 *
552 * + OmpExecutor specifies that the data should be stored and the
553 * associated operations executed on an OpenMP-supporting device (e.g.
554 * host CPU);
555 * + CudaExecutor specifies that the data should be stored and the
556 * operations executed on the NVIDIA GPU accelerator;
557 * + HipExecutor specifies that the data should be stored and the
558 * operations executed on the GPU accelerator using HIP;
559 * + ReferenceExecutor executes a non-optimized reference implementation,
560 * which can be used to debug the library.
561 */
563
564 /**
565 * Destructor.
566 */
567 virtual ~GinkgoPreconditioner() = default;
568
569 /**
570 * Generate the preconditioner for the given matrix @p op,
571 * which must be of MFEM SparseMatrix type.
572 * Calling this function is only required when creating a
573 * preconditioner for use with another MFEM solver; to use with
574 * a Ginkgo solver, get the LinOpFactory pointer through @p GetFactory()
575 * and pass to the Ginkgo solver constructor.
576 */
577 void SetOperator(const Operator &op) override;
578
579 /**
580 * Apply the preconditioner to input vector @p x, with out @p y.
581 */
582 void Mult(const Vector &x, Vector &y) const override;
583
584 /**
585 * Return a pointer to the LinOpFactory that will generate the preconditioner
586 * with the parameters set through the specific constructor.
587 */
588 const std::shared_ptr<gko::LinOpFactory> GetFactory() const
589 {
590 return this->precond_gen;
591 };
592
593 /**
594 * Return a pointer to the generated preconditioner for a specific matrix
595 * (that has previously been set with @p SetOperator).
596 */
597 const std::shared_ptr<gko::LinOp> GetGeneratedPreconditioner() const
598 {
599 return this->generated_precond;
600 };
601
602 /**
603 * Return whether this GinkgoPreconditioner object has an explicitly-
604 * generated preconditioner, built for a specific matrix.
605 */
607 {
608 return this->has_generated_precond;
609 };
610
611protected:
612 /**
613 * The Ginkgo generated solver factory object.
614 */
615 std::shared_ptr<gko::LinOpFactory> precond_gen;
616
617 /**
618 * Generated Ginkgo preconditioner for a specific matrix, created through
619 * @p SetOperator(), or a wrapped MFEM preconditioner.
620 * Must exist to use @p Mult().
621 */
622 std::shared_ptr<gko::LinOp> generated_precond;
623
624 /**
625 * The execution paradigm in Ginkgo. The choices are between
626 * `gko::OmpExecutor`, `gko::CudaExecutor` and `gko::ReferenceExecutor`
627 * and more details can be found in Ginkgo's documentation.
628 */
629 std::shared_ptr<gko::Executor> executor;
630
631 /**
632 * Whether or not we have generated a specific preconditioner for
633 * a matrix.
634 */
636
637};
638
639/**
640* This class forms the base class for all of Ginkgo's iterative solvers.
641* It is not intended to be used directly by MFEM applications. The various
642* derived classes only take the additional data that is specific to them
643* and solve the given linear system. The entire collection of solvers that
644* Ginkgo implements is available at the Ginkgo documentation and manual pages,
645* https://ginkgo-project.github.io/ginkgo/doc/develop.
646*
647* @ingroup Ginkgo
648*/
650{
651public:
652 /**
653 * Return a pointer to the LinOpFactory that will generate the solver
654 * with the parameters set through the specific constructor.
655 */
656 const std::shared_ptr<gko::LinOpFactory> GetFactory() const
657 {
658 return this->solver_gen;
659 };
660
661 void SetPrintLevel(int print_lvl) { print_level = print_lvl; }
662
663 int GetNumIterations() const { return final_iter; }
664 int GetConverged() const { return converged; }
665 real_t GetFinalNorm() const { return final_norm; }
666
667 /**
668 * If the Operator is a SparseMatrix, set up a Ginkgo Csr matrix
669 * to use its data directly. If the Operator is not a matrix,
670 * create an OperatorWrapper for it and store.
671 */
672 void SetOperator(const Operator &op) override;
673
674 /**
675 * Solve the linear system <tt>Ax=y</tt>. Dependent on the information
676 * provided by derived classes one of Ginkgo's linear solvers is chosen.
677 */
678 void Mult(const Vector &x, Vector &y) const override;
679
680 /**
681 * Return whether this GinkgoIterativeSolver object will use
682 * VectorWrapper types for input and output vectors.
683 * Note that Mult() will automatically create these wrappers if needed.
684 */
686 {
687 return this->needs_wrapped_vecs;
688 };
689
690 /**
691 * Destructor.
692 */
693 virtual ~GinkgoIterativeSolver() = default;
694
695protected:
696 /**
697 * Constructor.
698 *
699 * The @p exec defines the paradigm where the solution is computed.
700 * @p use_implicit_res_norm is for internal use by the derived classes
701 * for specific Ginkgo solvers; it indicates whether the solver makes
702 * an implicit residual norm estimate available for convergence checking.
703 * Each derived class automatically sets the correct value when calling this
704 * base class constructor.
705 *
706 */
709
716 mutable int final_iter;
717 mutable int converged;
718
719 /**
720 * The Ginkgo solver factory object, to generate specific solvers.
721 */
722 std::shared_ptr<gko::LinOpFactory> solver_gen;
723
724 /**
725 * The Ginkgo solver object, generated for a specific operator.
726 */
727 std::shared_ptr<gko::LinOp> solver;
728
729 /**
730 * The residual criterion object that controls the reduction of the residual
731 * relative to the initial residual.
732 */
733 std::shared_ptr<gko::stop::ResidualNorm<>::Factory>
735
736 /**
737 * The residual criterion object that controls the reduction of the residual
738 * based on an absolute tolerance.
739 */
740 std::shared_ptr<gko::stop::ResidualNorm<>::Factory>
742
743 /**
744 * The implicit residual criterion object that controls the reduction of the residual
745 * relative to the initial residual, based on an implicit residual norm value.
746 */
747 std::shared_ptr<gko::stop::ImplicitResidualNorm<>::Factory>
749
750 /**
751 * The implicit residual criterion object that controls the reduction of the residual
752 * based on an absolute tolerance, based on an implicit residual norm value.
753 */
754 std::shared_ptr<gko::stop::ImplicitResidualNorm<>::Factory>
756
757 /**
758 * The Ginkgo convergence logger used to check for convergence and other
759 * solver data if needed.
760 */
761 mutable std::shared_ptr<gko::log::Convergence<>> convergence_logger;
762
763 /**
764 * The residual logger object used to check for convergence and other solver
765 * data if needed.
766 */
767 mutable std::shared_ptr<ResidualLogger<>> residual_logger;
768
769 /**
770 * The Ginkgo combined factory object is used to create a combined stopping
771 * criterion to be passed to the solver.
772 */
773 std::shared_ptr<gko::stop::Combined::Factory> combined_factory;
774
775 /**
776 * The execution paradigm in Ginkgo. The choices are between
777 * `gko::OmpExecutor`, `gko::CudaExecutor` and `gko::ReferenceExecutor`
778 * and more details can be found in Ginkgo's documentation.
779 */
780 std::shared_ptr<gko::Executor> executor;
781
782 /**
783 * Whether or not we need to use VectorWrapper types with this solver.
784 */
786
787 /**
788 * Whether or not we need to use VectorWrapper types with the preconditioner
789 * or an inner solver. This value is set upon creation of the
790 * GinkgoIterativeSolver object and should never change.
791 */
793
794 /** Rebuild the Ginkgo stopping criterion factory with the latest values
795 * of rel_tol, abs_tol, and max_iter.
796 */
797 void update_stop_factory();
798
799private:
800 /**
801 * Initialize the Ginkgo logger object with event masks. Refer to the logging
802 * event masks in Ginkgo's .../include/ginkgo/core/log/logger.hpp.
803 */
804 void
805 initialize_ginkgo_log(gko::matrix::Dense<real_t>* b) const;
806
807 /**
808 * Pointer to either a Ginkgo CSR matrix or an OperatorWrapper wrapping
809 * an MFEM Operator (for matrix-free evaluation).
810 */
811 std::shared_ptr<gko::LinOp> system_oper;
812
813};
814
815/**
816 * This class adds helper functions for updating Ginkgo factories
817 * and solvers, when the full class type is needed. The derived classes
818 * should inherit from this class, rather than from GinkgoIterativeSolver
819 * directly.
820 */
821template<typename SolverType>
823{
824public:
827
828 void SetRelTol(real_t rtol)
829 {
830 rel_tol = rtol;
831 this->update_stop_factory();
832 auto current_params = gko::as<typename SolverType::Factory>
833 (solver_gen)->get_parameters();
834 this->solver_gen = current_params.with_criteria(this->combined_factory)
835 .on(this->executor);
836 if (solver)
837 {
838 gko::as<SolverType>(solver)->set_stop_criterion_factory(combined_factory);
839 }
840 }
841
842 void SetAbsTol(real_t atol)
843 {
844 abs_tol = atol;
845 this->update_stop_factory();
846 auto current_params = gko::as<typename SolverType::Factory>
847 (solver_gen)->get_parameters();
848 this->solver_gen = current_params.with_criteria(this->combined_factory)
849 .on(this->executor);
850 if (solver)
851 {
852 gko::as<SolverType>(solver)->set_stop_criterion_factory(combined_factory);
853 }
854 }
855
856 void SetMaxIter(int max_it)
857 {
858 max_iter = max_it;
859 this->update_stop_factory();
860 auto current_params = gko::as<typename SolverType::Factory>
861 (solver_gen)->get_parameters();
862 this->solver_gen = current_params.with_criteria(this->combined_factory)
863 .on(this->executor);
864 if (solver)
865 {
866 gko::as<SolverType>(solver)->set_stop_criterion_factory(combined_factory);
867 }
868 }
869};
870
871
872/**
873 * An implementation of the solver interface using the Ginkgo CG solver.
874 *
875 * @ingroup Ginkgo
876 */
877class CGSolver : public EnableGinkgoSolver<gko::solver::Cg<real_t>>
878{
879public:
880 /**
881 * Constructor.
882 *
883 * @param[in] exec The execution paradigm for the solver.
884 */
886
887 /**
888 * Constructor.
889 *
890 * @param[in] exec The execution paradigm for the solver.
891 * @param[in] preconditioner The preconditioner for the solver.
892 */
894 const GinkgoPreconditioner &preconditioner);
895};
896
897
898/**
899 * An implementation of the solver interface using the Ginkgo BiCGStab solver.
900 *
901 * @ingroup Ginkgo
902 */
903class BICGSTABSolver : public EnableGinkgoSolver<gko::solver::Bicgstab<real_t>>
904{
905public:
906 /**
907 * Constructor.
908 *
909 * @param[in] exec The execution paradigm for the solver.
910 */
912
913 /**
914 * Constructor.
915 *
916 * @param[in] exec The execution paradigm for the solver.
917 * @param[in] preconditioner The preconditioner for the solver.
918 */
920 const GinkgoPreconditioner &preconditioner);
921
922};
923
924/**
925 * An implementation of the solver interface using the Ginkgo CGS solver.
926 *
927 * CGS or the conjugate gradient square method is an iterative type Krylov
928 * subspace method which is suitable for general systems.
929 *
930 * @ingroup Ginkgo
931 */
932class CGSSolver : public EnableGinkgoSolver<gko::solver::Cgs<real_t>>
933{
934public:
935 /**
936 * Constructor.
937 *
938 * @param[in] exec The execution paradigm for the solver.
939 */
941
942 /**
943 * Constructor.
944 *
945 * @param[in] exec The execution paradigm for the solver.
946 * @param[in] preconditioner The preconditioner for the solver.
947 */
949 const GinkgoPreconditioner &preconditioner);
950};
951
952/**
953 * An implementation of the solver interface using the Ginkgo FCG solver.
954 *
955 * FCG or the flexible conjugate gradient method is an iterative type Krylov
956 * subspace method which is suitable for symmetric positive definite methods.
957 *
958 * Though this method performs very well for symmetric positive definite
959 * matrices, it is in general not suitable for general matrices.
960 *
961 * In contrast to the standard CG based on the Polack-Ribiere formula, the
962 * flexible CG uses the Fletcher-Reeves formula for creating the orthonormal
963 * vectors spanning the Krylov subspace. This increases the computational cost
964 * of every Krylov solver iteration but allows for non-constant preconditioners.
965 *
966 * @ingroup Ginkgo
967 */
968class FCGSolver : public EnableGinkgoSolver<gko::solver::Fcg<real_t>>
969{
970public:
971 /**
972 * Constructor.
973 *
974 * @param[in] exec The execution paradigm for the solver.
975 */
977
978 /**
979 * Constructor.
980 *
981 * @param[in] exec The execution paradigm for the solver.
982 * @param[in] preconditioner The preconditioner for the solver.
983 */
985 const GinkgoPreconditioner &preconditioner);
986};
987
988/**
989 * An implementation of the solver interface using the Ginkgo GMRES solver.
990 *
991 * @ingroup Ginkgo
992 */
993class GMRESSolver : public EnableGinkgoSolver<gko::solver::Gmres<real_t>>
994{
995public:
996 /**
997 * Constructor.
998 *
999 * @param[in] exec The execution paradigm for the solver.
1000 * @param[in] dim The Krylov dimension of the solver. Value of 0 will
1001 * let Ginkgo use its own internal default value.
1002 */
1003 GMRESSolver(GinkgoExecutor &exec, int dim = 0);
1004
1005 /**
1006 * Constructor.
1007 *
1008 * @param[in] exec The execution paradigm for the solver.
1009 * @param[in] preconditioner The preconditioner for the solver.
1010 * @param[in] dim The Krylov dimension of the solver. Value of 0 will
1011 * let Ginkgo use its own internal default value.
1012 */
1014 const GinkgoPreconditioner &preconditioner,
1015 int dim = 0);
1016
1017 /**
1018 * Change the Krylov dimension of the solver.
1019 */
1020 void SetKDim(int dim);
1021
1022protected:
1023 int m; // Dimension of Krylov subspace
1024};
1025
1026using gko::solver::cb_gmres::storage_precision;
1027/**
1028 * An implementation of the solver interface using the Ginkgo
1029 * Compressed Basis GMRES solver. With CB-GMRES, the Krylov basis
1030 * is "compressed" by storing in a lower precision. Currently, computations
1031 * are always performed in the MFEM-defined `real_t` precision, when using this
1032 * MFEM integration.
1033 * The Ginkgo storage precision options are accessed
1034 * through Ginkgo::storage_precision::*. The default choice
1035 * is Ginkgo::storage_precision::reduce1, i.e., store in float
1036 * instead of double or half instead of float.
1037 *
1038 * @ingroup Ginkgo
1039 */
1040class CBGMRESSolver : public EnableGinkgoSolver<gko::solver::CbGmres<real_t>>
1041{
1042public:
1043 /**
1044 * Constructor.
1045 *
1046 * @param[in] exec The execution paradigm for the solver.
1047 * @param[in] dim The Krylov dimension of the solver. Value of 0 will
1048 * let Ginkgo use its own internal default value.
1049 * @param[in] prec The storage precision used in the CB-GMRES. Options
1050 * are: keep (keep `real_t` precision), reduce1 (double -> float
1051 * or float -> half), reduce2 (double -> half or float -> half),
1052 * integer (`real_t` -> int64), ireduce1 (double -> int32 or
1053 * float -> int16), ireduce2 (double -> int16 or float -> int16).
1054 * See Ginkgo documentation for more about CB-GMRES.
1055 */
1056 CBGMRESSolver(GinkgoExecutor &exec, int dim = 0,
1057 storage_precision prec = storage_precision::reduce1);
1058
1059 /**
1060 * Constructor.
1061 *
1062 * @param[in] exec The execution paradigm for the solver.
1063 * @param[in] preconditioner The preconditioner for the solver.
1064 * @param[in] dim The Krylov dimension of the solver. Value of 0 will
1065 * let Ginkgo use its own internal default value.
1066 * @param[in] prec The storage precision used in the CB-GMRES. Options
1067 * are: keep (keep `real_t` precision), reduce1 (double -> float
1068 * or float -> half), reduce2 (double -> half or float -> half),
1069 * integer (`real_t` -> int64), ireduce1 (double -> int32 or
1070 * float -> int16), ireduce2 (double -> int16 or float -> int16).
1071 * See Ginkgo documentation for more about CB-GMRES.
1072 */
1074 const GinkgoPreconditioner &preconditioner,
1075 int dim = 0,
1076 storage_precision prec = storage_precision::reduce1);
1077
1078 /**
1079 * Change the Krylov dimension of the solver.
1080 */
1081 void SetKDim(int dim);
1082
1083protected:
1084 int m; // Dimension of Krylov subspace
1085};
1086
1087/**
1088 * An implementation of the solver interface using the Ginkgo IR solver.
1089 *
1090 * Iterative refinement (IR) is an iterative method that uses another coarse
1091 * method to approximate the error of the current solution via the current
1092 * residual.
1093 *
1094 * @ingroup Ginkgo
1095 */
1096class IRSolver : public EnableGinkgoSolver<gko::solver::Ir<real_t>>
1097{
1098public:
1099 /**
1100 * Constructor.
1101 *
1102 * @param[in] exec The execution paradigm for the solver.
1103 */
1104 IRSolver(GinkgoExecutor &exec);
1105
1106 /**
1107 * Constructor.
1108 *
1109 * @param[in] exec The execution paradigm for the solver.
1110 * @param[in] inner_solver The inner solver for the main solver.
1111 */
1113 const GinkgoIterativeSolver &inner_solver);
1114
1115};
1116
1117/**
1118 * An implementation of the preconditioner interface using the Ginkgo Jacobi
1119 * preconditioner.
1120 *
1121 * @ingroup Ginkgo
1122 */
1124{
1125public:
1126 /**
1127 * Constructor.
1128 *
1129 * @param[in] exec The execution paradigm for the preconditioner.
1130 * @param[in] storage_opt The storage optimization parameter.
1131 * @param[in] accuracy The relative accuracy for the adaptive version.
1132 * @param[in] max_block_size Maximum block size.
1133 * See the Ginkgo documentation for more information on these parameters.
1134 */
1136 GinkgoExecutor &exec,
1137 const std::string &storage_opt = "none",
1138 const real_t accuracy = 1.e-1,
1139 const int max_block_size = 32
1140 );
1141};
1142
1143/**
1144 * An implementation of the preconditioner interface using the Ginkgo
1145 * Incomplete LU preconditioner (ILU(0)).
1146 *
1147 * @ingroup Ginkgo
1148 */
1150{
1151public:
1152 /**
1153 * Constructor.
1154 *
1155 * @param[in] exec The execution paradigm for the preconditioner.
1156 * @param[in] factorization_type The factorization type: "exact" or
1157 * "parilu".
1158 * @param[in] sweeps The number of sweeps to do in the ParIlu
1159 * factorization algorithm. A value of 0 tells Ginkgo to use its
1160 * internal default value. This parameter is ignored in the case
1161 * of an exact factorization.
1162 * @param[in] skip_sort Only set this to true if the input matrix
1163 * that will be used to generate this preconditioner is guaranteed
1164 * to be sorted by column.
1165 *
1166 * Note: The use of this preconditioner will sort any input matrix
1167 * given to it, potentially changing the order of the stored values.
1168 */
1170 GinkgoExecutor &exec,
1171 const std::string &factorization_type = "exact",
1172 const int sweeps = 0,
1173 const bool skip_sort = false
1174 );
1175};
1176
1177/**
1178 * An implementation of the preconditioner interface using the Ginkgo
1179 * Incomplete LU-Incomplete Sparse Approximate Inverse preconditioner.
1180 * The Ilu-ISAI preconditioner differs from the Ilu preconditioner in
1181 * that Incomplete Sparse Approximate Inverses (ISAIs) are formed
1182 * to approximate solving the triangular systems defined by L and U.
1183 * When the preconditioner is applied, these ISAI matrices are applied
1184 * through matrix-vector multiplication.
1185 *
1186 * @ingroup Ginkgo
1187 */
1189{
1190public:
1191 /**
1192 * Constructor.
1193 *
1194 * @param[in] exec The execution paradigm for the preconditioner.
1195 * @param[in] factorization_type The factorization type: "exact" or
1196 * "parilu".
1197 * @param[in] sweeps The number of sweeps to do in the ParIlu
1198 * factorization algorithm. A value of 0 tells Ginkgo to use its
1199 * internal default value. This parameter is ignored in the case
1200 * of an exact factorization.
1201 * @param[in] sparsity_power Parameter determining the sparsity pattern of
1202 * the ISAI approximations.
1203 * @param[in] skip_sort Only set this to true if the input matrix
1204 * that will be used to generate this preconditioner is guaranteed
1205 * to be sorted by column.
1206 * See the Ginkgo documentation for more information on these parameters.
1207 *
1208 * Note: The use of this preconditioner will sort any input matrix
1209 * given to it, potentially changing the order of the stored values.
1210 */
1212 GinkgoExecutor &exec,
1213 const std::string &factorization_type = "exact",
1214 const int sweeps = 0,
1215 const int sparsity_power = 1,
1216 const bool skip_sort = false
1217 );
1218};
1219
1220/**
1221 * An implementation of the preconditioner interface using the Ginkgo
1222 * Incomplete Cholesky preconditioner (IC(0)).
1223 *
1224 * @ingroup Ginkgo
1225 */
1227{
1228public:
1229 /**
1230 * Constructor.
1231 *
1232 * @param[in] exec The execution paradigm for the preconditioner.
1233 * @param[in] factorization_type The factorization type: "exact" or
1234 * "paric".
1235 * @param[in] sweeps The number of sweeps to do in the ParIc
1236 * factorization algorithm. A value of 0 tells Ginkgo to use its
1237 * internal default value. This parameter is ignored in the case
1238 * of an exact factorization.
1239 * @param[in] skip_sort Only set this to true if the input matrix
1240 * that will be used to generate this preconditioner is guaranteed
1241 * to be sorted by column.
1242 *
1243 * Note: The use of this preconditioner will sort any input matrix
1244 * given to it, potentially changing the order of the stored values.
1245 */
1247 GinkgoExecutor &exec,
1248 const std::string &factorization_type = "exact",
1249 const int sweeps = 0,
1250 const bool skip_sort = false
1251 );
1252};
1253
1254/**
1255 * An implementation of the preconditioner interface using the Ginkgo
1256 * Incomplete Cholesky-Incomplete Sparse Approximate Inverse preconditioner.
1257 * The Ic-ISAI preconditioner differs from the Ic preconditioner in
1258 * that Incomplete Sparse Approximate Inverses (ISAIs) are formed
1259 * to approximate solving the triangular systems defined by L and L^T.
1260 * When the preconditioner is applied, these ISAI matrices are applied
1261 * through matrix-vector multiplication.
1262 *
1263 * @ingroup Ginkgo
1264 */
1266{
1267public:
1268 /**
1269 * Constructor.
1270 *
1271 * @param[in] exec The execution paradigm for the preconditioner.
1272 * @param[in] factorization_type The factorization type: "exact" or
1273 * "paric".
1274 * @param[in] sweeps The number of sweeps to do in the ParIc
1275 * factorization algorithm. A value of 0 tells Ginkgo to use its
1276 * internal default value. This parameter is ignored in the case
1277 * of an exact factorization.
1278 * @param[in] sparsity_power Parameter determining the sparsity pattern of
1279 * the ISAI approximations.
1280 * @param[in] skip_sort Only set this to true if the input matrix
1281 * that will be used to generate this preconditioner is guaranteed
1282 * to be sorted by column.
1283 * See the Ginkgo documentation for more information on these parameters.
1284 *
1285 * Note: The use of this preconditioner will sort any input matrix
1286 * given to it, potentially changing the order of the stored values.
1287 */
1289 GinkgoExecutor &exec,
1290 const std::string &factorization_type = "exact",
1291 const int sweeps = 0,
1292 const int sparsity_power = 1,
1293 const bool skip_sort = false
1294 );
1295};
1296
1297/**
1298 * A wrapper that allows Ginkgo to use MFEM preconditioners.
1299 *
1300 * @ingroup Ginkgo
1301 */
1303{
1304public:
1305 /**
1306 * Constructor.
1307 *
1308 * @param[in] exec The execution paradigm for the preconditioner.
1309 * @param[in] mfem_precond The MFEM Preconditioner to wrap.
1310 */
1312 GinkgoExecutor &exec,
1313 const Solver &mfem_precond
1314 );
1315
1316 /**
1317 * SetOperator is not allowed for this type of preconditioner;
1318 * this function overrides the base class in order to give an
1319 * error if SetOperator() is called for this class.
1320 */
1321 void SetOperator(const Operator &op) override
1322 {
1323 MFEM_ABORT("Ginkgo::MFEMPreconditioner must be constructed "
1324 "with the MFEM Operator that it will wrap as an argument;\n"
1325 "calling SetOperator() is not allowed.");
1326 };
1327};
1328} // namespace Ginkgo
1329
1330} // namespace mfem
1331
1332#endif // MFEM_USE_GINKGO
1333
1334#endif // MFEM_GINKGO
The MFEM Device class abstracts hardware devices such as GPUs, as well as programming models such as ...
Definition device.hpp:123
BICGSTABSolver(GinkgoExecutor &exec)
Definition ginkgo.cpp:734
CBGMRESSolver(GinkgoExecutor &exec, int dim=0, storage_precision prec=storage_precision::reduce1)
Definition ginkgo.cpp:943
CGSSolver(GinkgoExecutor &exec)
Definition ginkgo.cpp:773
CGSolver(GinkgoExecutor &exec)
Definition ginkgo.cpp:695
void SetRelTol(real_t rtol)
Definition ginkgo.hpp:828
void SetAbsTol(real_t atol)
Definition ginkgo.hpp:842
EnableGinkgoSolver(GinkgoExecutor &exec, bool use_implicit_res_norm)
Definition ginkgo.hpp:825
FCGSolver(GinkgoExecutor &exec)
Definition ginkgo.cpp:811
GMRESSolver(GinkgoExecutor &exec, int dim=0)
Definition ginkgo.cpp:849
void SetKDim(int dim)
Definition ginkgo.cpp:928
std::shared_ptr< gko::Executor > GetExecutor() const
Definition ginkgo.hpp:524
GinkgoExecutor(ExecType exec_type)
Definition ginkgo.cpp:33
virtual ~GinkgoExecutor()=default
@ CUDA
CUDA GPU Executor.
Definition ginkgo.hpp:478
@ HIP
HIP GPU Executor.
Definition ginkgo.hpp:480
@ REFERENCE
Reference CPU Executor.
Definition ginkgo.hpp:474
@ OMP
OpenMP CPU Executor.
Definition ginkgo.hpp:476
void SetOperator(const Operator &op) override
Definition ginkgo.cpp:631
std::shared_ptr< gko::stop::ImplicitResidualNorm<>::Factory > imp_abs_criterion
Definition ginkgo.hpp:755
void Mult(const Vector &x, Vector &y) const override
Definition ginkgo.cpp:503
std::shared_ptr< ResidualLogger<> > residual_logger
Definition ginkgo.hpp:767
std::shared_ptr< gko::stop::ResidualNorm<>::Factory > abs_criterion
Definition ginkgo.hpp:741
std::shared_ptr< gko::stop::Combined::Factory > combined_factory
Definition ginkgo.hpp:773
const std::shared_ptr< gko::LinOpFactory > GetFactory() const
Definition ginkgo.hpp:656
std::shared_ptr< gko::log::Convergence<> > convergence_logger
Definition ginkgo.hpp:761
GinkgoIterativeSolver(GinkgoExecutor &exec, bool use_implicit_res_norm)
Definition ginkgo.cpp:342
void SetPrintLevel(int print_lvl)
Definition ginkgo.hpp:661
std::shared_ptr< gko::LinOpFactory > solver_gen
Definition ginkgo.hpp:722
std::shared_ptr< gko::stop::ResidualNorm<>::Factory > rel_criterion
Definition ginkgo.hpp:734
virtual ~GinkgoIterativeSolver()=default
std::shared_ptr< gko::stop::ImplicitResidualNorm<>::Factory > imp_rel_criterion
Definition ginkgo.hpp:748
std::shared_ptr< gko::Executor > executor
Definition ginkgo.hpp:780
std::shared_ptr< gko::LinOp > solver
Definition ginkgo.hpp:727
void SetOperator(const Operator &op) override
Definition ginkgo.cpp:1110
std::shared_ptr< gko::Executor > executor
Definition ginkgo.hpp:629
const std::shared_ptr< gko::LinOpFactory > GetFactory() const
Definition ginkgo.hpp:588
void Mult(const Vector &x, Vector &y) const override
Definition ginkgo.cpp:1080
virtual ~GinkgoPreconditioner()=default
std::shared_ptr< gko::LinOp > generated_precond
Definition ginkgo.hpp:622
GinkgoPreconditioner(GinkgoExecutor &exec)
Definition ginkgo.cpp:1071
std::shared_ptr< gko::LinOpFactory > precond_gen
Definition ginkgo.hpp:615
const std::shared_ptr< gko::LinOp > GetGeneratedPreconditioner() const
Definition ginkgo.hpp:597
IRSolver(GinkgoExecutor &exec)
Definition ginkgo.cpp:1045
IcIsaiPreconditioner(GinkgoExecutor &exec, const std::string &factorization_type="exact", const int sweeps=0, const int sparsity_power=1, const bool skip_sort=false)
Definition ginkgo.cpp:1312
IcPreconditioner(GinkgoExecutor &exec, const std::string &factorization_type="exact", const int sweeps=0, const bool skip_sort=false)
Definition ginkgo.cpp:1276
IluIsaiPreconditioner(GinkgoExecutor &exec, const std::string &factorization_type="exact", const int sweeps=0, const int sparsity_power=1, const bool skip_sort=false)
Definition ginkgo.cpp:1219
IluPreconditioner(GinkgoExecutor &exec, const std::string &factorization_type="exact", const int sweeps=0, const bool skip_sort=false)
Definition ginkgo.cpp:1185
JacobiPreconditioner(GinkgoExecutor &exec, const std::string &storage_opt="none", const real_t accuracy=1.e-1, const int max_block_size=32)
Definition ginkgo.cpp:1154
MFEMPreconditioner(GinkgoExecutor &exec, const Solver &mfem_precond)
Definition ginkgo.cpp:1357
void SetOperator(const Operator &op) override
Definition ginkgo.hpp:1321
OperatorWrapper(std::shared_ptr< const gko::Executor > exec, gko::size_type size=0, const Operator *oper=NULL)
Definition ginkgo.hpp:223
void apply_impl(const gko::LinOp *b, gko::LinOp *x) const override
Definition ginkgo.cpp:418
static std::unique_ptr< VectorWrapper > create(std::shared_ptr< const gko::Executor > exec, gko::size_type size, Vector *mfem_vec, bool ownership=false)
Definition ginkgo.hpp:102
std::unique_ptr< gko::matrix::Dense< real_t > > create_with_same_config() const override
Definition ginkgo.hpp:121
std::unique_ptr< gko::matrix::Dense< real_t > > create_submatrix_impl(const gko::span &rows, const gko::span &columns, const gko::size_type stride) override
Definition ginkgo.hpp:172
std::unique_ptr< gko::matrix::Dense< real_t > > create_with_type_of_impl(std::shared_ptr< const gko::Executor > exec, const gko::dim< 2 > &size, gko::size_type stride) const override
Definition ginkgo.hpp:142
VectorWrapper(std::shared_ptr< const gko::Executor > exec, gko::size_type size, Vector *mfem_vec, bool ownership=false)
Definition ginkgo.hpp:69
const Vector & get_mfem_vec_const_ref() const
Definition ginkgo.hpp:116
void operator()(pointer ptr) const noexcept
Definition ginkgo.hpp:54
Abstract operator.
Definition operator.hpp:25
Base class for solvers.
Definition operator.hpp:780
Vector data type.
Definition vector.hpp:82
virtual void UseDevice(bool use_dev) const
Enable execution of Vector operations using the mfem::Device.
Definition vector.hpp:144
void MakeRef(Vector &base, int offset, int size)
Reset the Vector to be a reference to a sub-vector of base.
Definition vector.hpp:622
Vector beta
const real_t alpha
Definition ex15.cpp:369
int dim
Definition ex24.cpp:53
real_t b
Definition lissajous.cpp:42
real_t get_norm(const gko::matrix::Dense< ValueType > *norm)
Definition ginkgo.hpp:244
gko::array< T > gko_array
Definition ginkgo.hpp:40
real_t compute_norm(const gko::matrix::Dense< ValueType > *b)
Definition ginkgo.hpp:255
OutStream out(std::cout)
Global stream used by the library for standard output. Initially it uses the same std::streambuf as s...
Definition globals.hpp:66
float real_t
Definition config.hpp:43
void on_iteration_complete(const gko::LinOp *op, const gko::size_type &iteration, const gko::LinOp *residual, const gko::LinOp *solution, const gko::LinOp *residual_norm) const override
Definition ginkgo.hpp:317
void on_iteration_complete(const gko::LinOp *op, const gko::LinOp *rhs, const gko::LinOp *solution, const gko::size_type &iteration, const gko::LinOp *residual, const gko::LinOp *residual_norm, const gko::LinOp *implicit_sq_residual_norm, const gko::array< gko::stopping_status > *status, bool stopped) const override
Definition ginkgo.hpp:338
ResidualLogger(std::shared_ptr< const gko::Executor > exec, const gko::LinOp *matrix, const gko_dense *b, bool compute_real_residual=false)
Definition ginkgo.hpp:353
gko::matrix::Dense< ValueType > gko_dense
Definition ginkgo.hpp:314
void on_iteration_complete(const gko::LinOp *op, const gko::size_type &iteration, const gko::LinOp *residual, const gko::LinOp *solution, const gko::LinOp *residual_norm, const gko::LinOp *implicit_sq_residual_norm) const override
Definition ginkgo.hpp:327