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