MFEM  v4.5.1
Finite element discretization library
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Pages
petsc.hpp
Go to the documentation of this file.
1 // Copyright (c) 2010-2022, 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 // Author: Stefano Zampini <stefano.zampini@gmail.com>
13 
14 #ifndef MFEM_PETSC
15 #define MFEM_PETSC
16 
17 #include "../config/config.hpp"
18 
19 #ifdef MFEM_USE_PETSC
20 #ifdef MFEM_USE_MPI
21 
22 #include <limits>
23 
24 #include "handle.hpp"
25 #include "hypre.hpp"
26 #include "ode.hpp"
27 #include "../general/mem_manager.hpp"
28 
29 #include "petscconf.h"
30 #if !defined(PETSC_USE_REAL_DOUBLE)
31 #error "MFEM does not work with PETSc compiled without double precision"
32 #endif
33 #if defined(PETSC_USE_COMPLEX)
34 #error "MFEM does not work with PETSc compiled with complex numbers support"
35 #endif
36 #if defined(PETSC_USE_64BIT_INDICES) && !defined(HYPRE_BIGINT)
37 #error "Mismatch between HYPRE (32bit) and PETSc (64bit) integer types"
38 #endif
39 #if !defined(PETSC_USE_64BIT_INDICES) && defined(HYPRE_BIGINT)
40 #error "Mismatch between HYPRE (64bit) and PETSc (32bit) integer types"
41 #endif
42 
43 #include "petscversion.h"
44 #if PETSC_VERSION_GE(3,12,0)
45 #include "petscsystypes.h"
46 #else
47 typedef HYPRE_Int PetscInt;
48 typedef double PetscScalar;
49 typedef double PetscReal;
50 typedef int PetscClassId;
51 typedef struct _p_PetscObject *PetscObject;
52 #endif
53 
54 // forward declarations of PETSc internal structs
55 struct _p_Vec;
56 struct _p_Mat;
57 struct _p_KSP;
58 struct _p_PC;
59 struct _p_SNES;
60 struct _p_TS;
61 
62 
63 namespace mfem
64 {
65 
66 // Declare aliases of PETSc's types inside the namespace mfem::petsc:
67 namespace petsc
68 {
69 typedef struct ::_p_Vec *Vec;
70 typedef struct ::_p_Mat *Mat;
71 typedef struct ::_p_KSP *KSP;
72 typedef struct ::_p_PC *PC;
73 typedef struct ::_p_SNES *SNES;
74 typedef struct ::_p_TS *TS;
75 }
76 
77 /// Convenience functions to initialize/finalize PETSc
78 void MFEMInitializePetsc();
79 void MFEMInitializePetsc(int*,char***);
80 void MFEMInitializePetsc(int*,char***,const char[],const char[]);
81 void MFEMFinalizePetsc();
82 
83 /// Wrapper for syncing PETSc's vector memory
84 class PetscMemory : public Memory<double>
85 {
86 private:
87  Memory<double> *base;
88  bool read;
89  bool write;
90  bool usedev;
91 public:
92  PetscMemory() { Reset(); base = nullptr; }
93  void SetHostValid() const { flags |= VALID_HOST; }
94  void SetDeviceValid() const { flags |= VALID_DEVICE; }
95  void SetHostInvalid() const { flags &= ~VALID_HOST; }
96  void SetDeviceInvalid() const { flags &= ~VALID_DEVICE; }
97  inline bool IsAliasForSync() const { return base && (flags & ALIAS); }
98 
99  inline void MakeAliasForSync(const Memory<double> &base_, int offset_,
100  int size_, bool usedev_)
101  {
102  MFEM_VERIFY(!IsAliasForSync(),"Already alias");
103  base = (Memory<double>*)&base_;
104  read = true;
105  write = false;
106  usedev = usedev_;
107  MakeAlias(base_,offset_,size_);
108  }
109  inline void MakeAliasForSync(Memory<double> &base_, int offset_, int size_,
110  bool read_, bool write_, bool usedev_)
111  {
112  MFEM_VERIFY(!IsAliasForSync(),"Already alias");
113  base = (Memory<double>*)&base_;
114  read = read_;
115  write = write_;
116  usedev = usedev_;
117  MakeAlias(base_,offset_,size_);
118  }
119  inline void SyncBase()
120  {
121  MFEM_VERIFY(IsAliasForSync(),"MakeAliasForSynch not called");
122  base->Sync(*this);
123  }
124  inline void SyncBaseAndReset()
125  {
126  SyncBase();
127  base = nullptr;
128  Reset();
129  }
130  inline bool ReadRequested() const
131  {
132  MFEM_VERIFY(IsAliasForSync(),"MakeAliasForSynch not called");
133  return read;
134  }
135  inline bool WriteRequested() const
136  {
137  MFEM_VERIFY(IsAliasForSync(),"MakeAliasForSynch not called");
138  return write;
139  }
140  inline bool DeviceRequested() const
141  {
142  MFEM_VERIFY(IsAliasForSync(),"MakeAliasForSynch not called");
143  return usedev;
144  }
145  const double *GetHostPointer() const;
146  const double *GetDevicePointer() const;
147 };
148 
149 /// Wrapper for PETSc's vector class
150 class ParFiniteElementSpace;
151 class PetscParMatrix;
152 
153 class PetscParVector : public Vector
154 {
155 protected:
156  /// The actual PETSc object
158 
160 
161  friend class PetscParMatrix;
162  friend class PetscODESolver;
163  friend class PetscLinearSolver;
164  friend class PetscPreconditioner;
165  friend class PetscNonlinearSolver;
166  friend class PetscBDDCSolver;
167 
168  // Set Vector::data and Vector::size from x
169  void SetDataAndSize_();
170 
171  // Set Vec type from Device type
172  void SetVecType_();
173 
174  // Update Memory flags from PETSc offloadmask
175  void SetFlagsFromMask_() const;
176 
177 public:
178  /// Creates vector with given global size and partitioning of the columns.
179  /** If @a col is provided, processor P owns columns [col[P],col[P+1]).
180  Otherwise, PETSc decides the partitioning */
181  PetscParVector(MPI_Comm comm, PetscInt glob_size, PetscInt *col = NULL);
182 
183  /** @brief Creates vector with given global size, partitioning of the
184  columns, and data.
185 
186  The data must be allocated and destroyed outside. If @a data_ is NULL, a
187  dummy vector without a valid data array will be created. */
188  PetscParVector(MPI_Comm comm, PetscInt glob_size, PetscScalar *data_,
189  PetscInt *col);
190 
191  /// Creates vector compatible with @a y
192  PetscParVector(const PetscParVector &y);
193 
194  /** @brief Creates a PetscParVector from a Vector
195  @param[in] comm MPI communicator on which the new object lives
196  @param[in] x_ The mfem Vector (data is not shared)
197  @param[in] copy Whether to copy the data in x_ or not */
198  PetscParVector(MPI_Comm comm, const Vector &x_, bool copy = false);
199 
200  /** @brief Creates vector compatible with the Operator (i.e. in the domain
201  of) @a op or its adjoint. */
202  /** The argument @a allocate determines if the memory is actually allocated
203  to store the data. */
204  explicit PetscParVector(MPI_Comm comm, const Operator &op,
205  bool transpose = false, bool allocate = true);
206 
207  /// Creates vector compatible with (i.e. in the domain of) @a A or @a A^T
208  /** The argument @a allocate determines if the memory is actually allocated
209  to store the data. */
210  explicit PetscParVector(const PetscParMatrix &A, bool transpose = false,
211  bool allocate = true);
212 
213  /// Creates PetscParVector out of PETSc Vec object.
214  /** @param[in] y The PETSc Vec object.
215  @param[in] ref If true, we increase the reference count of @a y. */
216  explicit PetscParVector(petsc::Vec y, bool ref=false);
217 
218  /// Create a true dof parallel vector on a given ParFiniteElementSpace
219  explicit PetscParVector(ParFiniteElementSpace *pfes);
220 
221  /// Calls PETSc's destroy function
222  virtual ~PetscParVector();
223 
224  /// Get the associated MPI communicator
225  MPI_Comm GetComm() const;
226 
227  /// Returns the global number of rows
228  PetscInt GlobalSize() const;
229 
230  /// Typecasting to PETSc's Vec type
231  operator petsc::Vec() const { return x; }
232 
233  /// Typecasting to PETSc object
234  operator PetscObject() const { return (PetscObject)x; }
235 
236  /// Returns the global vector in each processor
237  Vector* GlobalVector() const;
238 
239  /// Set constant values
241 
242  /** @brief Set block size of a vector.
243 
244  @note This will error if the local size of the vector is not a multiple
245  of the block size @a bs.
246  @note This is a logically collective operation, so all processes need
247  to call it. */
248  void SetBlockSize(PetscInt bs);
249 
250  /** @brief Set values in a vector.
251 
252  @note Any process can insert in any location.
253  @note This is a collective operation, so all processes need to call it. */
255 
256  /** @brief Add values in a vector.
257 
258  @note Any process can add to any location.
259  @note This is a collective operation, so all processes need to call it. */
261 
262  /// Define operators for PETSc vectors.
268 
269  /** @brief Temporarily replace the data of the PETSc Vec object. To return to
270  the original data array, call ResetArray().
271 
272  @note This method calls PETSc's VecPlaceArray() function.
273  @note The inherited Vector::data pointer is not affected by this call. */
274  void PlaceArray(PetscScalar *temp_data);
275 
276  /** @brief Reset the PETSc Vec object to use its default data. Call this
277  method after the use of PlaceArray().
278 
279  @note This method calls PETSc's VecResetArray() function. */
280  void ResetArray();
281 
282  /** @brief This requests write access from where the memory is valid
283  and temporarily replaces the corresponding array used by the PETSc Vec
284  The bool parameter indicates read/write request */
285  void PlaceMemory(Memory<double>&,bool=false);
286 
287  /** @brief This requests read access from where the memory is valid
288  and temporarily replaces the corresponding array used by the PETSc Vec */
289  void PlaceMemory(const Memory<double>&);
290 
291  /** @brief Completes the operation started with PlaceMemory */
292  void ResetMemory();
293 
294  /** @brief Update PETSc's Vec after having accessed its data via GetMemory() */
295  void UpdateVecFromFlags();
296 
297  /// Set random values
298  void Randomize(PetscInt seed = 0);
299 
300  /// Prints the vector (to stdout if @a fname is NULL)
301  void Print(const char *fname = NULL, bool binary = false) const;
302 
303  const double *Read(bool=true) const override;
304  const double *HostRead() const override;
305  double *Write(bool=true) override;
306  double *HostWrite() override;
307  double *ReadWrite(bool=true) override;
308  double *HostReadWrite() override;
309  bool UseDevice() const override;
310  void UseDevice(bool) const override;
311 };
312 
313 
314 /// Wrapper for PETSc's matrix class
315 class PetscParMatrix : public Operator
316 {
317 protected:
318  /// The actual PETSc object
320 
321  /// Auxiliary vectors for typecasting
322  mutable PetscParVector *X, *Y;
323 
324  /// Initialize with defaults. Does not initialize inherited members.
325  void Init();
326 
327  /// Delete all owned data. Does not perform re-initialization with defaults.
328  void Destroy();
329 
330  /** @brief Creates a wrapper around a mfem::Operator @a op using PETSc's
331  MATSHELL object and returns the Mat in @a B.
332 
333  This does not take any reference to @a op, that should not be destroyed
334  until @a B is needed. */
335  void MakeWrapper(MPI_Comm comm, const Operator* op, petsc::Mat *B);
336 
337  /// Convert an mfem::Operator into a Mat @a B; @a op can be destroyed unless
338  /// tid == PETSC_MATSHELL or tid == PETSC_MATHYPRE
339  /// if op is a BlockOperator, the operator type is relevant to the individual
340  /// blocks
341  void ConvertOperator(MPI_Comm comm, const Operator& op, petsc::Mat *B,
342  Operator::Type tid);
343 
344  friend class PetscLinearSolver;
345  friend class PetscPreconditioner;
346 
347 private:
348  /// Constructs a block-diagonal Mat object
349  void BlockDiagonalConstructor(MPI_Comm comm, PetscInt *row_starts,
350  PetscInt *col_starts, SparseMatrix *diag,
351  bool assembled, petsc::Mat *A);
352 
353  void SetUpForDevice();
354 
355 public:
356  /// Create an empty matrix to be used as a reference to an existing matrix.
357  PetscParMatrix();
358 
359  /// Creates PetscParMatrix out of PETSc's Mat.
360  /** @param[in] a The PETSc Mat object.
361  @param[in] ref If true, we increase the reference count of @a a. */
362  PetscParMatrix(petsc::Mat a, bool ref=false);
363 
364  /** @brief Convert a PetscParMatrix @a pa with a new PETSc format @a tid.
365  Note that if @a pa is already a PetscParMatrix of the same type as
366  @a tid, the resulting PetscParMatrix will share the same Mat object */
367  explicit PetscParMatrix(const PetscParMatrix *pa, Operator::Type tid);
368 
369  /** @brief Creates a PetscParMatrix extracting the submatrix of @a A with
370  @a rows row indices and @a cols column indices */
371  PetscParMatrix(const PetscParMatrix& A, const Array<PetscInt>& rows,
372  const Array<PetscInt>& cols);
373 
374  /** @brief Convert a HypreParMatrix @a ha to a PetscParMatrix in the given
375  PETSc format @a tid. */
376  /** The supported type ids are: Operator::PETSC_MATAIJ,
377  Operator::PETSC_MATIS, Operator::PETSC_MATSHELL and
378  Operator::PETSC_MATHYPRE
379  @a ha can be destroyed unless tid == PETSC_MATSHELL or
380  tid == PETSC_MATHYPRE */
381  explicit PetscParMatrix(const HypreParMatrix *ha,
383 
384  /** @brief Convert a SparseMatrix @a ha to a PetscParMatrix in the given
385  PETSc format @a tid. */
386  explicit PetscParMatrix(const SparseMatrix *sa,
388 
389  /** @brief Convert an mfem::Operator into a PetscParMatrix in the given PETSc
390  format @a tid. */
391  /** If @a tid is Operator::PETSC_MATSHELL and @a op is not a PetscParMatrix,
392  it converts any mfem::Operator @a op implementing Operator::Mult() and
393  Operator::MultTranspose() into a PetscParMatrix. The Operator @a op
394  should not be deleted while the constructed PetscParMatrix is used.
395 
396  Otherwise, it tries to convert the operator in PETSc's classes.
397  @a op cannot be destroyed if tid == PETSC_MATHYPRE.
398 
399  In particular, if @a op is a BlockOperator, then a MATNEST Mat object is
400  created using @a tid as the type for the blocks.
401  Note that if @a op is already a PetscParMatrix of the same type as
402  @a tid, the resulting PetscParMatrix will share the same Mat object */
403  PetscParMatrix(MPI_Comm comm, const Operator *op,
405 
406  /// Creates block-diagonal square parallel matrix.
407  /** The block-diagonal is given by @a diag which must be in CSR format
408  (finalized). The new PetscParMatrix does not take ownership of any of the
409  input arrays. The type id @a tid can be either PETSC_MATAIJ (parallel
410  distributed CSR) or PETSC_MATIS. */
411  PetscParMatrix(MPI_Comm comm, PetscInt glob_size, PetscInt *row_starts,
412  SparseMatrix *diag, Operator::Type tid);
413 
414  /// Creates block-diagonal rectangular parallel matrix.
415  /** The block-diagonal is given by @a diag which must be in CSR format
416  (finalized). The new PetscParMatrix does not take ownership of any of the
417  input arrays. The type id @a tid can be either PETSC_MATAIJ (parallel
418  distributed CSR) or PETSC_MATIS. */
419  PetscParMatrix(MPI_Comm comm, PetscInt global_num_rows,
420  PetscInt global_num_cols, PetscInt *row_starts,
421  PetscInt *col_starts, SparseMatrix *diag,
422  Operator::Type tid);
423 
424  /// Calls PETSc's destroy function.
425  virtual ~PetscParMatrix() { Destroy(); }
426 
427  /// Replace the inner Mat Object. The reference count of newA is increased
428  void SetMat(petsc::Mat newA);
429 
430  /// @name Assignment operators
431  ///@{
436  ///@}
437 
438  /// Matvec: @a y = @a a A @a x + @a b @a y.
439  void Mult(double a, const Vector &x, double b, Vector &y) const;
440 
441  /// Matvec transpose: @a y = @a a A^T @a x + @a b @a y.
442  void MultTranspose(double a, const Vector &x, double b, Vector &y) const;
443 
444  virtual void Mult(const Vector &x, Vector &y) const
445  { Mult(1.0, x, 0.0, y); }
446 
447  virtual void MultTranspose(const Vector &x, Vector &y) const
448  { MultTranspose(1.0, x, 0.0, y); }
449 
450  /// Get the associated MPI communicator
451  MPI_Comm GetComm() const;
452 
453  /// Typecasting to PETSc's Mat type
454  operator petsc::Mat() const { return A; }
455 
456  /// Typecasting to PETSc object
457  operator PetscObject() const { return (PetscObject)A; }
458 
459  /// Returns the global index of the first local row
460  PetscInt GetRowStart() const;
461 
462  /// Returns the global index of the first local column
463  PetscInt GetColStart() const;
464 
465  /// Returns the local number of rows
466  PetscInt GetNumRows() const;
467 
468  /// Returns the local number of columns
469  PetscInt GetNumCols() const;
470 
471  /// Returns the global number of rows
472  PetscInt M() const;
473 
474  /// Returns the global number of columns
475  PetscInt N() const;
476 
477  /// Returns the global number of rows
478  PetscInt GetGlobalNumRows() const { return M(); }
479 
480  /// Returns the global number of columns
481  PetscInt GetGlobalNumCols() const { return N(); }
482 
483  /// Returns the number of nonzeros.
484  /** Differently from HYPRE, this call is collective on the communicator,
485  as this number is not stored inside PETSc, but needs to be computed. */
486  PetscInt NNZ() const;
487 
488  /// Returns the inner vector in the domain of A (it creates it if needed)
489  PetscParVector* GetX() const;
490 
491  /// Returns the inner vector in the range of A (it creates it if needed)
492  PetscParVector* GetY() const;
493 
494  /// Returns the transpose of the PetscParMatrix.
495  /** If @a action is false, the new matrix is constructed with the PETSc
496  function MatTranspose().
497  If @a action is true, then the matrix is not actually transposed.
498  Instead, an object that behaves like the transpose is returned. */
499  PetscParMatrix* Transpose(bool action = false);
500 
501  /// Prints the matrix (to stdout if fname is NULL)
502  void Print(const char *fname = NULL, bool binary = false) const;
503 
504  /// Scale all entries by s: A_scaled = s*A.
505  void operator*=(double s);
506 
507  /** @brief Eliminate rows and columns from the matrix, and rows from the
508  vector @a B. Modify @a B with the BC values in @a X. Put @a diag
509  on the diagonal corresponding to eliminated entries */
510  void EliminateRowsCols(const Array<int> &rows_cols, const PetscParVector &X,
511  PetscParVector &B, double diag = 1.);
512  void EliminateRowsCols(const Array<int> &rows_cols, const HypreParVector &X,
513  HypreParVector &B, double diag = 1.);
514 
515  /** @brief Eliminate rows and columns from the matrix and store the
516  eliminated elements in a new matrix Ae (returned).
517 
518  The sum of the modified matrix and the returned matrix, Ae, is equal to
519  the original matrix. */
520  PetscParMatrix* EliminateRowsCols(const Array<int> &rows_cols);
521 
522  /// Scale the local row i by s(i).
523  void ScaleRows(const Vector & s);
524 
525  /// Scale the local col i by s(i).
526  void ScaleCols(const Vector & s);
527 
528  /// Shift diagonal by a constant
529  void Shift(double s);
530 
531  /// Shift diagonal by a vector
532  void Shift(const Vector & s);
533 
534  /** @brief Eliminate only the rows from the matrix */
535  void EliminateRows(const Array<int> &rows);
536 
537  /** @brief Set row and column block sizes of a matrix.
538 
539  @note This will error if the local sizes of the matrix are not a
540  multiple of the block sizes.
541  @note This is a logically collective operation, so all processes need
542  to call it. */
543  void SetBlockSize(PetscInt rbs,PetscInt cbs=-1);
544 
545  /// Makes this object a reference to another PetscParMatrix
546  void MakeRef(const PetscParMatrix &master);
547 
548  /** @brief Release the PETSc Mat object. If @a dereference is true, decrement
549  the refcount of the Mat object. */
550  petsc::Mat ReleaseMat(bool dereference);
551 
552  Type GetType() const;
553 };
554 
555 /// Returns the matrix A * B
556 PetscParMatrix * ParMult(const PetscParMatrix *A, const PetscParMatrix *B);
557 
558 /// Returns the matrix Rt^t * A * P
559 PetscParMatrix * RAP(PetscParMatrix *Rt, PetscParMatrix *A, PetscParMatrix *P);
560 
561 /// Returns the matrix R * A * P
562 PetscParMatrix * TripleMatrixProduct(PetscParMatrix *R, PetscParMatrix *A,
563  PetscParMatrix *P);
564 
565 /// Returns the matrix P^t * A * P
566 PetscParMatrix * RAP(PetscParMatrix *A, PetscParMatrix *P);
567 
568 /// Returns the matrix P^t * A * P
569 PetscParMatrix * RAP(HypreParMatrix *A, PetscParMatrix *P);
570 
571 /** @brief Eliminate essential BC specified by @a ess_dof_list from the solution
572  @a X to the r.h.s. @a B.
573 
574  Here, @a A is a matrix with eliminated BC, while @a Ae is such that
575  (@a A + @a Ae) is the original (Neumann) matrix before elimination. */
576 void EliminateBC(PetscParMatrix &A, PetscParMatrix &Ae,
577  const Array<int> &ess_dof_list, const Vector &X, Vector &B);
578 
579 /// Helper class for handling essential boundary conditions.
581 {
582 public:
583  enum Type
584  {
586  CONSTANT, ///< Constant in time b.c.
588  };
589 
591  bctype(type_), setup(false), eval_t(0.0),
592  eval_t_cached(std::numeric_limits<double>::min()) {}
593  PetscBCHandler(Array<int>& ess_tdof_list, Type type_ = ZERO);
594 
595  virtual ~PetscBCHandler() {}
596 
597  /// Returns the type of boundary conditions
598  Type GetType() const { return bctype; }
599 
600  /// Sets the type of boundary conditions
601  void SetType(enum Type type_) { bctype = type_; setup = false; }
602 
603  /// Boundary conditions evaluation
604  /** In the result vector, @a g, only values at the essential dofs need to be
605  set. */
606  virtual void Eval(double t, Vector &g)
607  { mfem_error("PetscBCHandler::Eval method not overloaded"); }
608 
609  /// Sets essential dofs (local, per-process numbering)
610  void SetTDofs(Array<int>& list);
611 
612  /// Gets essential dofs (local, per-process numbering)
613  Array<int>& GetTDofs() { return ess_tdof_list; }
614 
615  /// Sets the current time
616  void SetTime(double t) { eval_t = t; }
617 
618  /// SetUp the helper object, where @a n is the size of the solution vector
619  void SetUp(PetscInt n);
620 
621  /// y = x on ess_tdof_list_c and y = g (internally evaluated) on ess_tdof_list
622  void ApplyBC(const Vector &x, Vector &y);
623 
624  /// Replace boundary dofs with the current value
625  void ApplyBC(Vector &x);
626 
627  /// y = x-g on ess_tdof_list, the rest of y is unchanged
628  void FixResidualBC(const Vector& x, Vector& y);
629 
630  /// Replace boundary dofs with 0
631  void Zero(Vector &x);
632 
633  /// y = x on ess_tdof_list_c and y = 0 on ess_tdof_list
634  void ZeroBC(const Vector &x, Vector &y);
635 
636 private:
637  enum Type bctype;
638  bool setup;
639 
640  double eval_t;
641  double eval_t_cached;
642  Vector eval_g;
643 
644  Array<int> ess_tdof_list; //Essential true dofs
645 };
646 
647 // Helper class for user-defined preconditioners that needs to be setup
649 {
650 private:
651  std::string name;
652 public:
653  PetscPreconditionerFactory(const std::string &name_ = "MFEM Factory")
654  : name(name_) { }
655  const char* GetName() { return name.c_str(); }
656  virtual Solver *NewPreconditioner(const OperatorHandle& oh) = 0;
658 };
659 
660 // Forward declarations of helper classes
661 class PetscSolverMonitor;
662 
663 /// Abstract class for PETSc's solvers.
665 {
666 protected:
667  /// Boolean to handle SetFromOptions calls.
668  mutable bool clcustom;
669 
670  /// The actual PETSc object (KSP, PC, SNES or TS).
672 
673  /// The class id of the actual PETSc object
675 
676  /// Right-hand side and solution vector
677  mutable PetscParVector *B, *X;
678 
679  /// Handler for boundary conditions
681 
682  /// Private context for solver
683  void *private_ctx;
684 
685  /// Boolean to handle SetOperator calls.
686  mutable bool operatorset;
687 
688 public:
689  /// Construct an empty PetscSolver. Initialize protected objects to NULL.
690  PetscSolver();
691 
692  /// Destroy the PetscParVectors allocated (if any).
693  virtual ~PetscSolver();
694 
695  /** @name Update of PETSc options.
696  The following Set methods can be used to update the internal PETSc
697  options.
698  @note They will be overwritten by the options in the input PETSc file. */
699  ///@{
700  void SetTol(double tol);
701  void SetRelTol(double tol);
702  void SetAbsTol(double tol);
703  void SetMaxIter(int max_iter);
704  void SetPrintLevel(int plev);
705  ///@}
706 
707  /// Customize object with options set
708  /** If @a customize is false, it disables any options customization. */
709  void Customize(bool customize = true) const;
710  int GetConverged();
711  int GetNumIterations();
712  double GetFinalNorm();
713 
714  /// Sets user-defined monitoring routine.
716 
717  /// Sets the object to handle essential boundary conditions
718  void SetBCHandler(PetscBCHandler *bch);
719 
720  /// Sets the object for the creation of the preconditioner
722 
723  /// Conversion function to PetscObject.
724  operator PetscObject() const { return obj; }
725 
726  /// Get the associated MPI communicator
727  MPI_Comm GetComm() const;
728 
729 protected:
730  /// These two methods handle creation and destructions of
731  /// private data for the Solver objects
732  void CreatePrivateContext();
733  void FreePrivateContext();
734 };
735 
736 
737 /// Abstract class for PETSc's linear solvers.
738 class PetscLinearSolver : public PetscSolver, public Solver
739 {
740 private:
741  /// Internal flag to handle HypreParMatrix conversion or not.
742  bool wrap;
743  void MultKernel(const Vector &b, Vector &x, bool trans) const;
744 
745 public:
746  PetscLinearSolver(MPI_Comm comm, const std::string &prefix = std::string(),
747  bool wrap = true, bool iter_mode = false);
749  const std::string &prefix = std::string(), bool iter_mode = false);
750  /// Constructs a solver using a HypreParMatrix.
751  /** If @a wrap is true, then the MatMult ops of HypreParMatrix are wrapped.
752  No preconditioner can be automatically constructed from PETSc. If
753  @a wrap is false, the HypreParMatrix is converted into a the AIJ
754  PETSc format, which is suitable for most preconditioning methods. */
755  PetscLinearSolver(const HypreParMatrix &A, bool wrap = true,
756  const std::string &prefix = std::string(), bool iter_mode = false);
757  virtual ~PetscLinearSolver();
758 
759  /// Sets the operator to be used for mat-vec operations and
760  /// for the construction of the preconditioner
761  virtual void SetOperator(const Operator &op);
762 
763  /// Allows to prescribe a different operator (@a pop) to construct
764  /// the preconditioner
765  void SetOperator(const Operator &op, const Operator &pop);
766 
767  /// Sets the solver to perform preconditioning
768  /// preserves the linear operator for the mat-vec
769  void SetPreconditioner(Solver &precond);
770 
771  /// Application of the solver.
772  virtual void Mult(const Vector &b, Vector &x) const;
773  virtual void MultTranspose(const Vector &b, Vector &x) const;
774 
775  /// Conversion function to PETSc's KSP type.
776  operator petsc::KSP() const { return (petsc::KSP)obj; }
777 };
778 
779 
781 {
782 public:
783  PetscPCGSolver(MPI_Comm comm, const std::string &prefix = std::string(),
784  bool iter_mode = false);
785  PetscPCGSolver(PetscParMatrix &A, const std::string &prefix = std::string(),
786  bool iter_mode = false);
787  PetscPCGSolver(HypreParMatrix &A, bool wrap = true,
788  const std::string &prefix = std::string(), bool iter_mode = false);
789 };
790 
791 
792 /// Abstract class for PETSc's preconditioners.
793 class PetscPreconditioner : public PetscSolver, public Solver
794 {
795 private:
796  void MultKernel(const Vector &b, Vector &x, bool trans) const;
797 
798 public:
799  PetscPreconditioner(MPI_Comm comm,
800  const std::string &prefix = std::string());
802  const std::string &prefix = std::string());
803  PetscPreconditioner(MPI_Comm comm, Operator &op,
804  const std::string &prefix = std::string());
805  virtual ~PetscPreconditioner();
806 
807  virtual void SetOperator(const Operator &op);
808 
809  /// Application of the preconditioner.
810  virtual void Mult(const Vector &b, Vector &x) const;
811  virtual void MultTranspose(const Vector &b, Vector &x) const;
812 
813  /// Conversion function to PETSc's PC type.
814  operator petsc::PC() const { return (petsc::PC)obj; }
815 };
816 
817 
818 /// Auxiliary class for BDDC customization.
820 {
821 protected:
827  bool netflux;
828  friend class PetscBDDCSolver;
829 
830 public:
832  nat_dof(NULL), nat_dof_local(false), netflux(false)
833  {}
835 
836  /// Specify dofs on the essential boundary.
837  /** If @a loc is false, it is a list of true dofs in local ordering.
838  If @a loc is true, it is a marker for Vdofs in local ordering. */
839  void SetEssBdrDofs(const Array<int> *essdofs, bool loc = false)
840  {
841  ess_dof = essdofs;
842  ess_dof_local = loc;
843  }
844  /// Specify dofs on the natural boundary.
845  /** If @a loc is false, it is a list of true dofs in local ordering.
846  If @a loc is true, it is a marker for Vdofs in local ordering. */
847  void SetNatBdrDofs(const Array<int> *natdofs, bool loc = false)
848  {
849  nat_dof = natdofs;
850  nat_dof_local = loc;
851  }
852  /// Setup BDDC with no-net-flux local solvers. Needs a ParFiniteElementSpace attached
853  void SetComputeNetFlux(bool net = true)
854  {
855  netflux = net;
856  }
857 };
858 
859 
861 {
862 private:
863  void BDDCSolverConstructor(const PetscBDDCSolverParams &opts);
864 
865 public:
866  PetscBDDCSolver(MPI_Comm comm, Operator &op,
868  const std::string &prefix = std::string());
871  const std::string &prefix = std::string());
872 };
873 
874 
876 {
877 public:
878  PetscFieldSplitSolver(MPI_Comm comm, Operator &op,
879  const std::string &prefix = std::string());
880 };
881 
883 {
884 private:
885  void H2SolverConstructor(ParFiniteElementSpace *fes);
886 
887 public:
890  const std::string &prefix = std::string());
891 
892 };
893 
894 /// Abstract class for PETSc's nonlinear solvers.
896 {
897 public:
898  PetscNonlinearSolver(MPI_Comm comm,
899  const std::string &prefix = std::string());
900  PetscNonlinearSolver(MPI_Comm comm, Operator &op,
901  const std::string &prefix = std::string());
902  virtual ~PetscNonlinearSolver();
903 
904  /// Specification of the nonlinear operator.
905  virtual void SetOperator(const Operator &op);
906 
907  /// Specifies the desired format of the Jacobian in case a PetscParMatrix
908  /// is not returned by the GetGradient method.
909  void SetJacobianType(Operator::Type type);
910 
911  /// Application of the solver.
912  virtual void Mult(const Vector &b, Vector &x) const;
913 
914  /// Specification of an objective function to be used for line search.
915  void SetObjective(void (*obj)(Operator* op, const Vector &x, double *f));
916 
917  /// User-defined routine to be applied after a successful line search step.
918  /// The user can change the current direction Y and/or the updated solution W
919  /// (with W = X - lambda * Y) but not the previous solution X.
920  /// If Y or W have been changed, the corresponding booleans need to updated.
921  void SetPostCheck(void (*post)(Operator *op, const Vector &X, Vector &Y,
922  Vector &W, bool &changed_y, bool &changed_w));
923 
924  /// General purpose update function to be called at the beginning of each step
925  /// it is the current nonlinear iteration number
926  /// F is the current function value, X the current solution
927  /// D the previous step taken, and P the previous solution
928  void SetUpdate(void (*update)(Operator *op, int it,
929  const Vector& F, const Vector& X,
930  const Vector& D, const Vector& P));
931 
932  /// Conversion function to PETSc's SNES type.
933  operator petsc::SNES() const { return (petsc::SNES)obj; }
934 };
935 
936 
937 /// Abstract class for PETSc's ODE solvers.
938 class PetscODESolver : public PetscSolver, public ODESolver
939 {
940 public:
941  /// The type of the ODE. Use ODE_SOLVER_LINEAR if the Jacobians
942  /// are linear and independent of time.
943  enum Type
944  {
947  };
948 
949  PetscODESolver(MPI_Comm comm, const std::string &prefix = std::string());
950  virtual ~PetscODESolver();
951 
952  /// Initialize the ODE solver.
953  virtual void Init(TimeDependentOperator &f_,
954  enum PetscODESolver::Type type);
956 
959 
960  /// Specifies the desired format of the Jacobian in case a PetscParMatrix
961  /// is not returned by the GetGradient methods
962  void SetJacobianType(Operator::Type type);
963 
964  virtual void Step(Vector &x, double &t, double &dt);
965  virtual void Run(Vector &x, double &t, double &dt, double t_final);
966 
967  /// Conversion function to PETSc's TS type.
968  operator petsc::TS() const { return (petsc::TS)obj; }
969 };
970 
971 
972 /// Abstract class for monitoring PETSc's solvers.
974 {
975 public:
976  bool mon_sol;
977  bool mon_res;
978  PetscSolverMonitor(bool monitor_sol = false, bool monitor_res = true)
979  : mon_sol(monitor_sol), mon_res(monitor_res) {}
980  virtual ~PetscSolverMonitor() {}
981 
982  /// Monitor the solution vector x
983  virtual void MonitorSolution(PetscInt it, PetscReal norm, const Vector &x)
984  {
985  MFEM_ABORT("MonitorSolution() is not implemented!")
986  }
987 
988  /// Monitor the residual vector r
989  virtual void MonitorResidual(PetscInt it, PetscReal norm, const Vector &r)
990  {
991  MFEM_ABORT("MonitorResidual() is not implemented!")
992  }
993 
994  /// Generic monitor to take access to the solver
995  virtual void MonitorSolver(PetscSolver* solver) {}
996 };
997 
998 
999 } // namespace mfem
1000 
1001 #endif // MFEM_USE_MPI
1002 #endif // MFEM_USE_PETSC
1003 
1004 #endif
void SetPreconditioner(Solver &precond)
Definition: petsc.cpp:3043
void SetPrintLevel(int plev)
Definition: petsc.cpp:2386
PetscParVector * B
Right-hand side and solution vector.
Definition: petsc.hpp:677
void CreatePrivateContext()
Definition: petsc.cpp:2685
void Print(const char *fname=NULL, bool binary=false) const
Prints the matrix (to stdout if fname is NULL)
Definition: petsc.cpp:1957
void trans(const Vector &u, Vector &x)
Definition: ex27.cpp:412
void EliminateBC(const HypreParMatrix &A, const HypreParMatrix &Ae, const Array< int > &ess_dof_list, const Vector &X, Vector &B)
Eliminate essential BC specified by ess_dof_list from the solution X to the r.h.s. B.
Definition: hypre.cpp:3263
void SetUpdate(void(*update)(Operator *op, int it, const Vector &F, const Vector &X, const Vector &D, const Vector &P))
Definition: petsc.cpp:4034
virtual void SetOperator(const Operator &op)
Specification of the nonlinear operator.
Definition: petsc.cpp:3941
Abstract class for PETSc&#39;s preconditioners.
Definition: petsc.hpp:793
PetscParMatrix & operator+=(const PetscParMatrix &B)
Definition: petsc.cpp:1085
double PetscScalar
Definition: petsc.hpp:48
PetscMemory pdata
Definition: petsc.hpp:159
petsc::Vec x
The actual PETSc object.
Definition: petsc.hpp:157
virtual void SetOperator(const Operator &op)
Definition: petsc.cpp:2898
bool clcustom
Boolean to handle SetFromOptions calls.
Definition: petsc.hpp:668
void SetObjective(void(*obj)(Operator *op, const Vector &x, double *f))
Specification of an objective function to be used for line search.
Definition: petsc.cpp:4009
HYPRE_Int PetscInt
Definition: petsc.hpp:47
double GetFinalNorm()
Definition: petsc.cpp:2660
virtual void MultTranspose(const Vector &b, Vector &x) const
Action of the transpose operator: y=A^t(x). The default behavior in class Operator is to generate an ...
Definition: petsc.cpp:3301
void MakeAliasForSync(Memory< double > &base_, int offset_, int size_, bool read_, bool write_, bool usedev_)
Definition: petsc.hpp:109
petsc::Mat A
The actual PETSc object.
Definition: petsc.hpp:319
Wrapper for PETSc&#39;s matrix class.
Definition: petsc.hpp:315
PetscParVector & operator+=(const PetscParVector &y)
Definition: petsc.cpp:693
const Array< int > * nat_dof
Definition: petsc.hpp:825
void ApplyBC(const Vector &x, Vector &y)
y = x on ess_tdof_list_c and y = g (internally evaluated) on ess_tdof_list
Definition: petsc.cpp:2768
void SetHostValid() const
Definition: petsc.hpp:93
Abstract class for PETSc&#39;s linear solvers.
Definition: petsc.hpp:738
Abstract class for PETSc&#39;s solvers.
Definition: petsc.hpp:664
PetscODESolver::Type GetType() const
Definition: petsc.cpp:4192
virtual void Run(Vector &x, double &t, double &dt, double t_final)
Perform time integration from time t [in] to time tf [in].
Definition: petsc.cpp:4254
void SetPostCheck(void(*post)(Operator *op, const Vector &X, Vector &Y, Vector &W, bool &changed_y, bool &changed_w))
Definition: petsc.cpp:4020
PetscBCHandler(Type type_=ZERO)
Definition: petsc.hpp:590
virtual void Mult(const Vector &x, Vector &y) const
Operator application: y=A(x).
Definition: petsc.hpp:444
ParFiniteElementSpace * fespace
Definition: petsc.hpp:822
Type GetType() const
Returns the type of boundary conditions.
Definition: petsc.hpp:598
Base abstract class for first order time dependent operators.
Definition: operator.hpp:289
PetscInt GetNumCols() const
Returns the local number of columns.
Definition: petsc.cpp:922
void MakeAliasForSync(const Memory< double > &base_, int offset_, int size_, bool usedev_)
Definition: petsc.hpp:99
void SetType(PetscODESolver::Type)
Definition: petsc.cpp:4198
void ConvertOperator(MPI_Comm comm, const Operator &op, petsc::Mat *B, Operator::Type tid)
Definition: petsc.cpp:1301
virtual void MonitorSolver(PetscSolver *solver)
Generic monitor to take access to the solver.
Definition: petsc.hpp:995
Pointer to an Operator of a specified type.
Definition: handle.hpp:33
void Mult(double a, const Vector &x, double b, Vector &y) const
Matvec: y = a A x + b y.
Definition: petsc.cpp:1922
PetscBDDCSolver(MPI_Comm comm, Operator &op, const PetscBDDCSolverParams &opts=PetscBDDCSolverParams(), const std::string &prefix=std::string())
Definition: petsc.cpp:3812
PetscParVector & AddValues(const Array< PetscInt > &, const Array< PetscScalar > &)
Add values in a vector.
Definition: petsc.cpp:672
virtual ~PetscSolver()
Destroy the PetscParVectors allocated (if any).
Definition: petsc.cpp:2297
PetscNonlinearSolver(MPI_Comm comm, const std::string &prefix=std::string())
Definition: petsc.cpp:3901
bool DeviceRequested() const
Definition: petsc.hpp:140
PetscInt GetGlobalNumRows() const
Returns the global number of rows.
Definition: petsc.hpp:478
Type GetType() const
Definition: petsc.cpp:2239
void SetBlockSize(PetscInt rbs, PetscInt cbs=-1)
Set row and column block sizes of a matrix.
Definition: petsc.cpp:950
virtual void Step(Vector &x, double &t, double &dt)
Perform a time step from time t [in] to time t [out] based on the requested step size dt [in]...
Definition: petsc.cpp:4216
struct::_p_Mat * Mat
Definition: petsc.hpp:70
Abstract class for solving systems of ODEs: dx/dt = f(x,t)
Definition: ode.hpp:22
void Init()
Initialize with defaults. Does not initialize inherited members.
Definition: petsc.cpp:956
PetscPreconditioner(MPI_Comm comm, const std::string &prefix=std::string())
Definition: petsc.cpp:3175
PetscClassId cid
The class id of the actual PETSc object.
Definition: petsc.hpp:674
void SetDataAndSize_()
Definition: petsc.cpp:218
Abstract parallel finite element space.
Definition: pfespace.hpp:28
Constant in time b.c.
Definition: petsc.hpp:586
void SetRelTol(double tol)
Definition: petsc.cpp:2309
virtual void MultTranspose(const Vector &x, Vector &y) const
Action of the transpose operator: y=A^t(x). The default behavior in class Operator is to generate an ...
Definition: petsc.hpp:447
Wrapper for syncing PETSc&#39;s vector memory.
Definition: petsc.hpp:84
void Randomize(PetscInt seed=0)
Set random values.
Definition: petsc.cpp:859
void Destroy()
Delete all owned data. Does not perform re-initialization with defaults.
Definition: petsc.cpp:1702
void SyncBaseAndReset()
Definition: petsc.hpp:124
void SetMonitor(PetscSolverMonitor *ctx)
Sets user-defined monitoring routine.
Definition: petsc.cpp:2471
void SetMaxIter(int max_iter)
Definition: petsc.cpp:2359
Abstract class for PETSc&#39;s ODE solvers.
Definition: petsc.hpp:938
PetscParVector & operator*=(PetscScalar d)
Definition: petsc.cpp:707
HypreParMatrix * ParMult(const HypreParMatrix *A, const HypreParMatrix *B, bool own_matrix)
Definition: hypre.cpp:2849
void EliminateRowsCols(const Array< int > &rows_cols, const PetscParVector &X, PetscParVector &B, double diag=1.)
Eliminate rows and columns from the matrix, and rows from the vector B. Modify B with the BC values i...
Definition: petsc.cpp:2183
double PetscReal
Definition: petsc.hpp:49
const double * HostRead() const override
Shortcut for mfem::Read(vec.GetMemory(), vec.Size(), false).
Definition: petsc.cpp:377
void Customize(bool customize=true) const
Customize object with options set.
Definition: petsc.cpp:2561
PetscInt M() const
Returns the global number of rows.
Definition: petsc.cpp:929
void SyncBase()
Definition: petsc.hpp:119
unsigned flags
Bit flags defined from the FlagMask enum.
int PetscClassId
Definition: petsc.hpp:50
virtual void Mult(const Vector &b, Vector &x) const
Application of the solver.
Definition: petsc.cpp:4047
PetscParVector * Y
Definition: petsc.hpp:322
virtual ~PetscLinearSolver()
Definition: petsc.cpp:3133
PetscInt GetNumRows() const
Returns the local number of rows.
Definition: petsc.cpp:915
void SetHostInvalid() const
Definition: petsc.hpp:95
struct::_p_PC * PC
Definition: petsc.hpp:72
void Print(const char *fname=NULL, bool binary=false) const
Prints the vector (to stdout if fname is NULL)
Definition: petsc.cpp:874
void SetDeviceInvalid() const
Definition: petsc.hpp:96
void ResetArray()
Reset the PETSc Vec object to use its default data. Call this method after the use of PlaceArray()...
Definition: petsc.cpp:726
const double * GetHostPointer() const
Definition: petsc.cpp:198
PetscParVector & SetValues(const Array< PetscInt > &, const Array< PetscScalar > &)
Set values in a vector.
Definition: petsc.cpp:658
double * Write(bool=true) override
Shortcut for mfem::Write(vec.GetMemory(), vec.Size(), on_dev).
Definition: petsc.cpp:382
virtual void Init(TimeDependentOperator &f_, enum PetscODESolver::Type type)
Initialize the ODE solver.
Definition: petsc.cpp:4103
PetscParMatrix * Transpose(bool action=false)
Returns the transpose of the PetscParMatrix.
Definition: petsc.cpp:1903
double f(const Vector &xvec)
Definition: lor_mms.hpp:32
PetscObject obj
The actual PETSc object (KSP, PC, SNES or TS).
Definition: petsc.hpp:671
Abstract class for monitoring PETSc&#39;s solvers.
Definition: petsc.hpp:973
Data type sparse matrix.
Definition: sparsemat.hpp:50
void SetDeviceValid() const
Definition: petsc.hpp:94
PetscLinearSolver(MPI_Comm comm, const std::string &prefix=std::string(), bool wrap=true, bool iter_mode=false)
Definition: petsc.cpp:2857
void SetJacobianType(Operator::Type type)
Definition: petsc.cpp:4003
struct::_p_TS * TS
Definition: petsc.hpp:74
void mfem_error(const char *msg)
Function called when an error is encountered. Used by the macros MFEM_ABORT, MFEM_ASSERT, MFEM_VERIFY.
Definition: error.cpp:154
double b
Definition: lissajous.cpp:42
Auxiliary class for BDDC customization.
Definition: petsc.hpp:819
virtual ~PetscPreconditioner()
Definition: petsc.cpp:3306
MPI_Comm GetComm() const
Get the associated MPI communicator.
Definition: petsc.cpp:1254
void FixResidualBC(const Vector &x, Vector &y)
y = x-g on ess_tdof_list, the rest of y is unchanged
Definition: petsc.cpp:2817
bool UseDevice() const override
Return the device flag of the Memory object used by the Vector.
Definition: petsc.cpp:448
void SetFlagsFromMask_() const
Definition: petsc.cpp:268
void SetType(enum Type type_)
Sets the type of boundary conditions.
Definition: petsc.hpp:601
PetscParMatrix & operator=(const PetscParMatrix &B)
Definition: petsc.cpp:1069
bool ReadRequested() const
Definition: petsc.hpp:130
PetscInt GetColStart() const
Returns the global index of the first local column.
Definition: petsc.cpp:908
virtual ~PetscParVector()
Calls PETSc&#39;s destroy function.
Definition: petsc.cpp:523
PetscFieldSplitSolver(MPI_Comm comm, Operator &op, const std::string &prefix=std::string())
Definition: petsc.cpp:3821
PetscBCHandler * bchandler
Handler for boundary conditions.
Definition: petsc.hpp:680
PetscParMatrix & operator-=(const PetscParMatrix &B)
Definition: petsc.cpp:1100
struct::_p_Vec * Vec
Definition: petsc.hpp:69
virtual Solver * NewPreconditioner(const OperatorHandle &oh)=0
double * ReadWrite(bool=true) override
Shortcut for mfem::ReadWrite(vec.GetMemory(), vec.Size(), on_dev).
Definition: petsc.cpp:411
void ScaleRows(const Vector &s)
Scale the local row i by s(i).
Definition: petsc.cpp:1982
void Zero(Vector &x)
Replace boundary dofs with 0.
Definition: petsc.cpp:2836
MPI_Comm GetComm() const
Get the associated MPI communicator.
Definition: petsc.cpp:2466
bool WriteRequested() const
Definition: petsc.hpp:135
void Shift(double s)
Shift diagonal by a constant.
Definition: petsc.cpp:2004
double * HostReadWrite() override
Shortcut for mfem::ReadWrite(vec.GetMemory(), vec.Size(), false).
Definition: petsc.cpp:435
void RAP(const DenseMatrix &A, const DenseMatrix &P, DenseMatrix &RAP)
Definition: densemat.cpp:3213
void UpdateVecFromFlags()
Update PETSc&#39;s Vec after having accessed its data via GetMemory()
Definition: petsc.cpp:301
void SetTime(double t)
Sets the current time.
Definition: petsc.hpp:616
ID for class PetscParMatrix, MATSHELL format.
Definition: operator.hpp:263
void Reset()
Reset the memory to be empty, ensuring that Delete() will be a no-op.
void SetMat(petsc::Mat newA)
Replace the inner Mat Object. The reference count of newA is increased.
Definition: petsc.cpp:1727
PetscParMatrix()
Create an empty matrix to be used as a reference to an existing matrix.
Definition: petsc.cpp:963
void MakeWrapper(MPI_Comm comm, const Operator *op, petsc::Mat *B)
Creates a wrapper around a mfem::Operator op using PETSc&#39;s MATSHELL object and returns the Mat in B...
Definition: petsc.cpp:1267
Abstract class for PETSc&#39;s nonlinear solvers.
Definition: petsc.hpp:895
int GetNumIterations()
Definition: petsc.cpp:2627
void SetBCHandler(PetscBCHandler *bch)
Sets the object to handle essential boundary conditions.
Definition: petsc.cpp:2501
struct::_p_KSP * KSP
Definition: petsc.hpp:71
Wrapper for hypre&#39;s parallel vector class.
Definition: hypre.hpp:161
void SetBlockSize(PetscInt bs)
Set block size of a vector.
Definition: petsc.cpp:463
const double * Read(bool=true) const override
Shortcut for mfem::Read(vec.GetMemory(), vec.Size(), on_dev).
Definition: petsc.cpp:354
void ZeroBC(const Vector &x, Vector &y)
y = x on ess_tdof_list_c and y = 0 on ess_tdof_list
Definition: petsc.cpp:2845
virtual void Eval(double t, Vector &g)
Boundary conditions evaluation.
Definition: petsc.hpp:606
Type
Enumeration defining IDs for some classes derived from Operator.
Definition: operator.hpp:256
struct::_p_SNES * SNES
Definition: petsc.hpp:73
PetscInt NNZ() const
Returns the number of nonzeros.
Definition: petsc.cpp:943
PetscPreconditionerFactory(const std::string &name_="MFEM Factory")
Definition: petsc.hpp:653
virtual void Mult(const Vector &b, Vector &x) const
Application of the preconditioner.
Definition: petsc.cpp:3296
virtual void Init(TimeDependentOperator &f_)
Associate a TimeDependentOperator with the ODE solver.
Definition: petsc.hpp:955
void Sync(const Memory &other) const
Copy the host/device pointer validity flags from other to *this.
PetscInt N() const
Returns the global number of columns.
Definition: petsc.cpp:936
PetscParVector * GetY() const
Returns the inner vector in the range of A (it creates it if needed)
Definition: petsc.cpp:1893
void SetJacobianType(Operator::Type type)
Definition: petsc.cpp:4186
double * HostWrite() override
Shortcut for mfem::Write(vec.GetMemory(), vec.Size(), false).
Definition: petsc.cpp:406
bool operatorset
Boolean to handle SetOperator calls.
Definition: petsc.hpp:686
PetscParVector * X
Auxiliary vectors for typecasting.
Definition: petsc.hpp:322
void MakeAlias(const Memory &base, int offset, int size)
Create a memory object that points inside the memory object base.
void MultTranspose(double a, const Vector &x, double b, Vector &y) const
Matvec transpose: y = a A^T x + b y.
Definition: petsc.cpp:1939
bool IsAliasForSync() const
Definition: petsc.hpp:97
Pointer is an alias.
PetscParVector & operator-=(const PetscParVector &y)
Definition: petsc.cpp:700
virtual ~PetscSolverMonitor()
Definition: petsc.hpp:980
double a
Definition: lissajous.cpp:41
void MFEMFinalizePetsc()
Definition: petsc.cpp:192
virtual ~PetscNonlinearSolver()
Definition: petsc.cpp:3933
void SetComputeNetFlux(bool net=true)
Setup BDDC with no-net-flux local solvers. Needs a ParFiniteElementSpace attached.
Definition: petsc.hpp:853
virtual ~PetscODESolver()
Definition: petsc.cpp:4095
virtual ~PetscBCHandler()
Definition: petsc.hpp:595
int GetConverged()
Definition: petsc.cpp:2594
PetscSolver()
Construct an empty PetscSolver. Initialize protected objects to NULL.
Definition: petsc.cpp:2287
void SetPreconditionerFactory(PetscPreconditionerFactory *factory)
Sets the object for the creation of the preconditioner.
Definition: petsc.cpp:2520
void MFEMInitializePetsc()
Convenience functions to initialize/finalize PETSc.
Definition: petsc.cpp:168
petsc::Mat ReleaseMat(bool dereference)
Release the PETSc Mat object. If dereference is true, decrement the refcount of the Mat object...
Definition: petsc.cpp:2226
PetscParVector(MPI_Comm comm, PetscInt glob_size, PetscInt *col=NULL)
Creates vector with given global size and partitioning of the columns.
Definition: petsc.cpp:505
virtual void MonitorSolution(PetscInt it, PetscReal norm, const Vector &x)
Monitor the solution vector x.
Definition: petsc.hpp:983
ID for class PetscParMatrix, MATAIJ format.
Definition: operator.hpp:261
virtual void SetOperator(const Operator &op)
Set/update the solver for the given operator.
Definition: petsc.cpp:3210
void SetSpace(ParFiniteElementSpace *fe)
Definition: petsc.hpp:834
MPI_Comm GetComm() const
Get the associated MPI communicator.
Definition: petsc.cpp:621
PetscInt GlobalSize() const
Returns the global number of rows.
Definition: petsc.cpp:456
PetscParVector * X
Definition: petsc.hpp:677
const Array< int > * ess_dof
Definition: petsc.hpp:823
PetscInt GetRowStart() const
Returns the global index of the first local row.
Definition: petsc.cpp:901
struct _p_PetscObject * PetscObject
Definition: petsc.hpp:51
PetscH2Solver(Operator &op, ParFiniteElementSpace *fes, const std::string &prefix=std::string())
Definition: petsc.cpp:3856
void ResetMemory()
Completes the operation started with PlaceMemory.
Definition: petsc.cpp:815
void SetAbsTol(double tol)
Definition: petsc.cpp:2334
Helper class for handling essential boundary conditions.
Definition: petsc.hpp:580
Class used by MFEM to store pointers to host and/or device memory.
PetscParVector * GetX() const
Returns the inner vector in the domain of A (it creates it if needed)
Definition: petsc.cpp:1883
Array< int > & GetTDofs()
Gets essential dofs (local, per-process numbering)
Definition: petsc.hpp:613
void SetTol(double tol)
Definition: petsc.cpp:2304
PetscParVector & operator=(PetscScalar d)
Set constant values.
Definition: petsc.cpp:651
PetscSolverMonitor(bool monitor_sol=false, bool monitor_res=true)
Definition: petsc.hpp:978
RefCoord t[3]
virtual void MonitorResidual(PetscInt it, PetscReal norm, const Vector &r)
Monitor the residual vector r.
Definition: petsc.hpp:989
PetscParMatrix * TripleMatrixProduct(PetscParMatrix *R, PetscParMatrix *A, PetscParMatrix *P)
Returns the matrix R * A * P.
Definition: petsc.cpp:2023
Vector data type.
Definition: vector.hpp:60
void PlaceArray(PetscScalar *temp_data)
Temporarily replace the data of the PETSc Vec object. To return to the original data array...
Definition: petsc.cpp:721
void FreePrivateContext()
Definition: petsc.cpp:2717
RefCoord s[3]
Vector * GlobalVector() const
Returns the global vector in each processor.
Definition: petsc.cpp:626
PetscODESolver(MPI_Comm comm, const std::string &prefix=std::string())
Definition: petsc.cpp:4070
void EliminateRows(const Array< int > &rows)
Eliminate only the rows from the matrix.
Definition: petsc.cpp:2212
Base class for solvers.
Definition: operator.hpp:655
PetscInt GetGlobalNumCols() const
Returns the global number of columns.
Definition: petsc.hpp:481
PetscPCGSolver(MPI_Comm comm, const std::string &prefix=std::string(), bool iter_mode=false)
Definition: petsc.cpp:3143
void SetNatBdrDofs(const Array< int > *natdofs, bool loc=false)
Specify dofs on the natural boundary.
Definition: petsc.hpp:847
Abstract operator.
Definition: operator.hpp:24
Wrapper for hypre&#39;s ParCSR matrix class.
Definition: hypre.hpp:343
virtual void Mult(const Vector &b, Vector &x) const
Application of the solver.
Definition: petsc.cpp:3123
void operator*=(double s)
Scale all entries by s: A_scaled = s*A.
Definition: petsc.cpp:1917
void SetEssBdrDofs(const Array< int > *essdofs, bool loc=false)
Specify dofs on the essential boundary.
Definition: petsc.hpp:839
void SetTDofs(Array< int > &list)
Sets essential dofs (local, per-process numbering)
Definition: petsc.cpp:2745
virtual void MultTranspose(const Vector &b, Vector &x) const
Action of the transpose operator: y=A^t(x). The default behavior in class Operator is to generate an ...
Definition: petsc.cpp:3128
void PlaceMemory(Memory< double > &, bool=false)
This requests write access from where the memory is valid and temporarily replaces the corresponding ...
Definition: petsc.cpp:731
void * private_ctx
Private context for solver.
Definition: petsc.hpp:683
void SetUp(PetscInt n)
SetUp the helper object, where n is the size of the solution vector.
Definition: petsc.cpp:2752
const double * GetDevicePointer() const
Definition: petsc.cpp:207
void MakeRef(const PetscParMatrix &master)
Makes this object a reference to another PetscParMatrix.
Definition: petsc.cpp:1873
void ScaleCols(const Vector &s)
Scale the local col i by s(i).
Definition: petsc.cpp:1993
virtual ~PetscParMatrix()
Calls PETSc&#39;s destroy function.
Definition: petsc.hpp:425