MFEM  v4.3.0
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-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 // 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 synching 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 values in a vector.
243 
244  @note any process can insert in any location
245  @note This is a collective operation, so all process needs to call it */
247 
248  /** @brief Add values in a vector.
249 
250  @note any process can add to any location
251  @note This is a collective operation, so all process needs to call it */
253 
254  /// Define operators for PETSc vectors.
260 
261  /** @brief Temporarily replace the data of the PETSc Vec object. To return to
262  the original data array, call ResetArray().
263 
264  @note This method calls PETSc's VecPlaceArray() function.
265  @note The inherited Vector::data pointer is not affected by this call. */
266  void PlaceArray(PetscScalar *temp_data);
267 
268  /** @brief Reset the PETSc Vec object to use its default data. Call this
269  method after the use of PlaceArray().
270 
271  @note This method calls PETSc's VecResetArray() function. */
272  void ResetArray();
273 
274  /** @brief This requests write access from where the memory is valid
275  and temporarily replaces the corresponding array used by the PETSc Vec
276  The bool parameter indicates read/write request */
277  void PlaceMemory(Memory<double>&,bool=false);
278 
279  /** @brief This requests read access from where the memory is valid
280  and temporarily replaces the corresponding array used by the PETSc Vec */
281  void PlaceMemory(const Memory<double>&);
282 
283  /** @brief Completes the operation started with PlaceMemory */
284  void ResetMemory();
285 
286  /** @brief Update PETSc's Vec after having accessed its data via GetMemory() */
287  void UpdateVecFromFlags();
288 
289  /// Set random values
290  void Randomize(PetscInt seed = 0);
291 
292  /// Prints the vector (to stdout if @a fname is NULL)
293  void Print(const char *fname = NULL, bool binary = false) const;
294 
295  const double *Read(bool=true) const override;
296  const double *HostRead() const override;
297  double *Write(bool=true) override;
298  double *HostWrite() override;
299  double *ReadWrite(bool=true) override;
300  double *HostReadWrite() override;
301  bool UseDevice() const override;
302  void UseDevice(bool) const override;
303 };
304 
305 
306 /// Wrapper for PETSc's matrix class
307 class PetscParMatrix : public Operator
308 {
309 protected:
310  /// The actual PETSc object
312 
313  /// Auxiliary vectors for typecasting
314  mutable PetscParVector *X, *Y;
315 
316  /// Initialize with defaults. Does not initialize inherited members.
317  void Init();
318 
319  /// Delete all owned data. Does not perform re-initialization with defaults.
320  void Destroy();
321 
322  /** @brief Creates a wrapper around a mfem::Operator @a op using PETSc's
323  MATSHELL object and returns the Mat in @a B.
324 
325  This does not take any reference to @a op, that should not be destroyed
326  until @a B is needed. */
327  void MakeWrapper(MPI_Comm comm, const Operator* op, petsc::Mat *B);
328 
329  /// Convert an mfem::Operator into a Mat @a B; @a op can be destroyed unless
330  /// tid == PETSC_MATSHELL or tid == PETSC_MATHYPRE
331  /// if op is a BlockOperator, the operator type is relevant to the individual
332  /// blocks
333  void ConvertOperator(MPI_Comm comm, const Operator& op, petsc::Mat *B,
334  Operator::Type tid);
335 
336  friend class PetscLinearSolver;
337  friend class PetscPreconditioner;
338 
339 private:
340  /// Constructs a block-diagonal Mat object
341  void BlockDiagonalConstructor(MPI_Comm comm, PetscInt *row_starts,
342  PetscInt *col_starts, SparseMatrix *diag,
343  bool assembled, petsc::Mat *A);
344 
345  void SetUpForDevice();
346 
347 public:
348  /// Create an empty matrix to be used as a reference to an existing matrix.
349  PetscParMatrix();
350 
351  /// Creates PetscParMatrix out of PETSc's Mat.
352  /** @param[in] a The PETSc Mat object.
353  @param[in] ref If true, we increase the reference count of @a a. */
354  PetscParMatrix(petsc::Mat a, bool ref=false);
355 
356  /** @brief Convert a PetscParMatrix @a pa with a new PETSc format @a tid.
357  Note that if @a pa is already a PetscParMatrix of the same type as
358  @a tid, the resulting PetscParMatrix will share the same Mat object */
359  explicit PetscParMatrix(const PetscParMatrix *pa, Operator::Type tid);
360 
361  /** @brief Creates a PetscParMatrix extracting the submatrix of @a A with
362  @a rows row indices and @a cols column indices */
363  PetscParMatrix(const PetscParMatrix& A, const Array<PetscInt>& rows,
364  const Array<PetscInt>& cols);
365 
366  /** @brief Convert a HypreParMatrix @a ha to a PetscParMatrix in the given
367  PETSc format @a tid. */
368  /** The supported type ids are: Operator::PETSC_MATAIJ,
369  Operator::PETSC_MATIS, Operator::PETSC_MATSHELL and
370  Operator::PETSC_MATHYPRE
371  @a ha can be destroyed unless tid == PETSC_MATSHELL or
372  tid == PETSC_MATHYPRE */
373  explicit PetscParMatrix(const HypreParMatrix *ha,
375 
376  /** @brief Convert a SparseMatrix @a ha to a PetscParMatrix in the given
377  PETSc format @a tid. */
378  explicit PetscParMatrix(const SparseMatrix *sa,
380 
381  /** @brief Convert an mfem::Operator into a PetscParMatrix in the given PETSc
382  format @a tid. */
383  /** If @a tid is Operator::PETSC_MATSHELL and @a op is not a PetscParMatrix,
384  it converts any mfem::Operator @a op implementing Operator::Mult() and
385  Operator::MultTranspose() into a PetscParMatrix. The Operator @a op
386  should not be deleted while the constructed PetscParMatrix is used.
387 
388  Otherwise, it tries to convert the operator in PETSc's classes.
389  @a op cannot be destroyed if tid == PETSC_MATHYPRE.
390 
391  In particular, if @a op is a BlockOperator, then a MATNEST Mat object is
392  created using @a tid as the type for the blocks.
393  Note that if @a op is already a PetscParMatrix of the same type as
394  @a tid, the resulting PetscParMatrix will share the same Mat object */
395  PetscParMatrix(MPI_Comm comm, const Operator *op,
397 
398  /// Creates block-diagonal square parallel matrix.
399  /** The block-diagonal is given by @a diag which must be in CSR format
400  (finalized). The new PetscParMatrix does not take ownership of any of the
401  input arrays. The type id @a tid can be either PETSC_MATAIJ (parallel
402  distributed CSR) or PETSC_MATIS. */
403  PetscParMatrix(MPI_Comm comm, PetscInt glob_size, PetscInt *row_starts,
404  SparseMatrix *diag, Operator::Type tid);
405 
406  /// Creates block-diagonal rectangular 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 global_num_rows,
412  PetscInt global_num_cols, PetscInt *row_starts,
413  PetscInt *col_starts, SparseMatrix *diag,
414  Operator::Type tid);
415 
416  /// Calls PETSc's destroy function.
417  virtual ~PetscParMatrix() { Destroy(); }
418 
419  /// Replace the inner Mat Object. The reference count of newA is increased
420  void SetMat(petsc::Mat newA);
421 
422  /// @name Assignment operators
423  ///@{
428  ///@}
429 
430  /// Matvec: @a y = @a a A @a x + @a b @a y.
431  void Mult(double a, const Vector &x, double b, Vector &y) const;
432 
433  /// Matvec transpose: @a y = @a a A^T @a x + @a b @a y.
434  void MultTranspose(double a, const Vector &x, double b, Vector &y) const;
435 
436  virtual void Mult(const Vector &x, Vector &y) const
437  { Mult(1.0, x, 0.0, y); }
438 
439  virtual void MultTranspose(const Vector &x, Vector &y) const
440  { MultTranspose(1.0, x, 0.0, y); }
441 
442  /// Get the associated MPI communicator
443  MPI_Comm GetComm() const;
444 
445  /// Typecasting to PETSc's Mat type
446  operator petsc::Mat() const { return A; }
447 
448  /// Typecasting to PETSc object
449  operator PetscObject() const { return (PetscObject)A; }
450 
451  /// Returns the global index of the first local row
452  PetscInt GetRowStart() const;
453 
454  /// Returns the global index of the first local column
455  PetscInt GetColStart() const;
456 
457  /// Returns the local number of rows
458  PetscInt GetNumRows() const;
459 
460  /// Returns the local number of columns
461  PetscInt GetNumCols() const;
462 
463  /// Returns the global number of rows
464  PetscInt M() const;
465 
466  /// Returns the global number of columns
467  PetscInt N() const;
468 
469  /// Returns the global number of rows
470  PetscInt GetGlobalNumRows() const { return M(); }
471 
472  /// Returns the global number of columns
473  PetscInt GetGlobalNumCols() const { return N(); }
474 
475  /// Returns the number of nonzeros.
476  /** Differently from HYPRE, this call is collective on the communicator,
477  as this number is not stored inside PETSc, but needs to be computed. */
478  PetscInt NNZ() const;
479 
480  /// Returns the inner vector in the domain of A (it creates it if needed)
481  PetscParVector* GetX() const;
482 
483  /// Returns the inner vector in the range of A (it creates it if needed)
484  PetscParVector* GetY() const;
485 
486  /// Returns the transpose of the PetscParMatrix.
487  /** If @a action is false, the new matrix is constructed with the PETSc
488  function MatTranspose().
489  If @a action is true, then the matrix is not actually transposed.
490  Instead, an object that behaves like the transpose is returned. */
491  PetscParMatrix* Transpose(bool action = false);
492 
493  /// Prints the matrix (to stdout if fname is NULL)
494  void Print(const char *fname = NULL, bool binary = false) const;
495 
496  /// Scale all entries by s: A_scaled = s*A.
497  void operator*=(double s);
498 
499  /** @brief Eliminate rows and columns from the matrix, and rows from the
500  vector @a B. Modify @a B with the BC values in @a X. Put @a diag
501  on the diagonal corresponding to eliminated entries */
502  void EliminateRowsCols(const Array<int> &rows_cols, const PetscParVector &X,
503  PetscParVector &B, double diag = 1.);
504  void EliminateRowsCols(const Array<int> &rows_cols, const HypreParVector &X,
505  HypreParVector &B, double diag = 1.);
506 
507  /** @brief Eliminate rows and columns from the matrix and store the
508  eliminated elements in a new matrix Ae (returned).
509 
510  The sum of the modified matrix and the returned matrix, Ae, is equal to
511  the original matrix. */
512  PetscParMatrix* EliminateRowsCols(const Array<int> &rows_cols);
513 
514  /// Scale the local row i by s(i).
515  void ScaleRows(const Vector & s);
516 
517  /// Scale the local col i by s(i).
518  void ScaleCols(const Vector & s);
519 
520  /// Shift diagonal by a constant
521  void Shift(double s);
522 
523  /// Shift diagonal by a vector
524  void Shift(const Vector & s);
525 
526  /** @brief Eliminate only the rows from the matrix */
527  void EliminateRows(const Array<int> &rows);
528 
529  /// Makes this object a reference to another PetscParMatrix
530  void MakeRef(const PetscParMatrix &master);
531 
532  /** @brief Release the PETSc Mat object. If @a dereference is true, decrement
533  the refcount of the Mat object. */
534  petsc::Mat ReleaseMat(bool dereference);
535 
536  Type GetType() const;
537 };
538 
539 /// Returns the matrix A * B
540 PetscParMatrix * ParMult(const PetscParMatrix *A, const PetscParMatrix *B);
541 
542 /// Returns the matrix Rt^t * A * P
543 PetscParMatrix * RAP(PetscParMatrix *Rt, PetscParMatrix *A, PetscParMatrix *P);
544 
545 /// Returns the matrix R * A * P
546 PetscParMatrix * TripleMatrixProduct(PetscParMatrix *R, PetscParMatrix *A,
547  PetscParMatrix *P);
548 
549 /// Returns the matrix P^t * A * P
550 PetscParMatrix * RAP(PetscParMatrix *A, PetscParMatrix *P);
551 
552 /// Returns the matrix P^t * A * P
553 PetscParMatrix * RAP(HypreParMatrix *A, PetscParMatrix *P);
554 
555 /** @brief Eliminate essential BC specified by @a ess_dof_list from the solution
556  @a X to the r.h.s. @a B.
557 
558  Here, @a A is a matrix with eliminated BC, while @a Ae is such that
559  (@a A + @a Ae) is the original (Neumann) matrix before elimination. */
560 void EliminateBC(PetscParMatrix &A, PetscParMatrix &Ae,
561  const Array<int> &ess_dof_list, const Vector &X, Vector &B);
562 
563 /// Helper class for handling essential boundary conditions.
565 {
566 public:
567  enum Type
568  {
570  CONSTANT, ///< Constant in time b.c.
572  };
573 
575  bctype(type_), setup(false), eval_t(0.0),
576  eval_t_cached(std::numeric_limits<double>::min()) {}
577  PetscBCHandler(Array<int>& ess_tdof_list, Type type_ = ZERO);
578 
579  virtual ~PetscBCHandler() {}
580 
581  /// Returns the type of boundary conditions
582  Type GetType() const { return bctype; }
583 
584  /// Sets the type of boundary conditions
585  void SetType(enum Type type_) { bctype = type_; setup = false; }
586 
587  /// Boundary conditions evaluation
588  /** In the result vector, @a g, only values at the essential dofs need to be
589  set. */
590  virtual void Eval(double t, Vector &g)
591  { mfem_error("PetscBCHandler::Eval method not overloaded"); }
592 
593  /// Sets essential dofs (local, per-process numbering)
594  void SetTDofs(Array<int>& list);
595 
596  /// Gets essential dofs (local, per-process numbering)
597  Array<int>& GetTDofs() { return ess_tdof_list; }
598 
599  /// Sets the current time
600  void SetTime(double t) { eval_t = t; }
601 
602  /// SetUp the helper object, where @a n is the size of the solution vector
603  void SetUp(PetscInt n);
604 
605  /// y = x on ess_tdof_list_c and y = g (internally evaluated) on ess_tdof_list
606  void ApplyBC(const Vector &x, Vector &y);
607 
608  /// Replace boundary dofs with the current value
609  void ApplyBC(Vector &x);
610 
611  /// y = x-g on ess_tdof_list, the rest of y is unchanged
612  void FixResidualBC(const Vector& x, Vector& y);
613 
614  /// Replace boundary dofs with 0
615  void Zero(Vector &x);
616 
617  /// y = x on ess_tdof_list_c and y = 0 on ess_tdof_list
618  void ZeroBC(const Vector &x, Vector &y);
619 
620 private:
621  enum Type bctype;
622  bool setup;
623 
624  double eval_t;
625  double eval_t_cached;
626  Vector eval_g;
627 
628  Array<int> ess_tdof_list; //Essential true dofs
629 };
630 
631 // Helper class for user-defined preconditioners that needs to be setup
633 {
634 private:
635  std::string name;
636 public:
637  PetscPreconditionerFactory(const std::string &name_ = "MFEM Factory")
638  : name(name_) { }
639  const char* GetName() { return name.c_str(); }
640  virtual Solver *NewPreconditioner(const OperatorHandle& oh) = 0;
642 };
643 
644 // Forward declarations of helper classes
645 class PetscSolverMonitor;
646 
647 /// Abstract class for PETSc's solvers.
649 {
650 protected:
651  /// Boolean to handle SetFromOptions calls.
652  mutable bool clcustom;
653 
654  /// The actual PETSc object (KSP, PC, SNES or TS).
656 
657  /// The class id of the actual PETSc object
659 
660  /// Right-hand side and solution vector
661  mutable PetscParVector *B, *X;
662 
663  /// Handler for boundary conditions
665 
666  /// Private context for solver
667  void *private_ctx;
668 
669  /// Boolean to handle SetOperator calls.
670  mutable bool operatorset;
671 
672 public:
673  /// Construct an empty PetscSolver. Initialize protected objects to NULL.
674  PetscSolver();
675 
676  /// Destroy the PetscParVectors allocated (if any).
677  virtual ~PetscSolver();
678 
679  /** @name Update of PETSc options.
680  The following Set methods can be used to update the internal PETSc
681  options.
682  @note They will be overwritten by the options in the input PETSc file. */
683  ///@{
684  void SetTol(double tol);
685  void SetRelTol(double tol);
686  void SetAbsTol(double tol);
687  void SetMaxIter(int max_iter);
688  void SetPrintLevel(int plev);
689  ///@}
690 
691  /// Customize object with options set
692  /** If @a customize is false, it disables any options customization. */
693  void Customize(bool customize = true) const;
694  int GetConverged();
695  int GetNumIterations();
696  double GetFinalNorm();
697 
698  /// Sets user-defined monitoring routine.
700 
701  /// Sets the object to handle essential boundary conditions
702  void SetBCHandler(PetscBCHandler *bch);
703 
704  /// Sets the object for the creation of the preconditioner
706 
707  /// Conversion function to PetscObject.
708  operator PetscObject() const { return obj; }
709 
710  /// Get the associated MPI communicator
711  MPI_Comm GetComm() const;
712 
713 protected:
714  /// These two methods handle creation and destructions of
715  /// private data for the Solver objects
716  void CreatePrivateContext();
717  void FreePrivateContext();
718 };
719 
720 
721 /// Abstract class for PETSc's linear solvers.
722 class PetscLinearSolver : public PetscSolver, public Solver
723 {
724 private:
725  /// Internal flag to handle HypreParMatrix conversion or not.
726  bool wrap;
727  void MultKernel(const Vector &b, Vector &x, bool trans) const;
728 
729 public:
730  PetscLinearSolver(MPI_Comm comm, const std::string &prefix = std::string(),
731  bool wrap = true);
733  const std::string &prefix = std::string());
734  /// Constructs a solver using a HypreParMatrix.
735  /** If @a wrap is true, then the MatMult ops of HypreParMatrix are wrapped.
736  No preconditioner can be automatically constructed from PETSc. If
737  @a wrap is false, the HypreParMatrix is converted into a the AIJ
738  PETSc format, which is suitable for most preconditioning methods. */
739  PetscLinearSolver(const HypreParMatrix &A, bool wrap = true,
740  const std::string &prefix = std::string());
741  virtual ~PetscLinearSolver();
742 
743  /// Sets the operator to be used for mat-vec operations and
744  /// for the construction of the preconditioner
745  virtual void SetOperator(const Operator &op);
746 
747  /// Allows to prescribe a different operator (@a pop) to construct
748  /// the preconditioner
749  void SetOperator(const Operator &op, const Operator &pop);
750 
751  /// Sets the solver to perform preconditioning
752  /// preserves the linear operator for the mat-vec
753  void SetPreconditioner(Solver &precond);
754 
755  /// Application of the solver.
756  virtual void Mult(const Vector &b, Vector &x) const;
757  virtual void MultTranspose(const Vector &b, Vector &x) const;
758 
759  /// Conversion function to PETSc's KSP type.
760  operator petsc::KSP() const { return (petsc::KSP)obj; }
761 };
762 
763 
765 {
766 public:
767  PetscPCGSolver(MPI_Comm comm, const std::string &prefix = std::string());
768  PetscPCGSolver(PetscParMatrix &A, const std::string &prefix = std::string());
769  PetscPCGSolver(HypreParMatrix &A,bool wrap=true,
770  const std::string &prefix = std::string());
771 };
772 
773 
774 /// Abstract class for PETSc's preconditioners.
775 class PetscPreconditioner : public PetscSolver, public Solver
776 {
777 private:
778  void MultKernel(const Vector &b, Vector &x, bool trans) const;
779 
780 public:
781  PetscPreconditioner(MPI_Comm comm,
782  const std::string &prefix = std::string());
784  const std::string &prefix = std::string());
785  PetscPreconditioner(MPI_Comm comm, Operator &op,
786  const std::string &prefix = std::string());
787  virtual ~PetscPreconditioner();
788 
789  virtual void SetOperator(const Operator &op);
790 
791  /// Application of the preconditioner.
792  virtual void Mult(const Vector &b, Vector &x) const;
793  virtual void MultTranspose(const Vector &b, Vector &x) const;
794 
795  /// Conversion function to PETSc's PC type.
796  operator petsc::PC() const { return (petsc::PC)obj; }
797 };
798 
799 
800 /// Auxiliary class for BDDC customization.
802 {
803 protected:
809  bool netflux;
810  friend class PetscBDDCSolver;
811 
812 public:
814  nat_dof(NULL), nat_dof_local(false), netflux(false)
815  {}
817 
818  /// Specify dofs on the essential boundary.
819  /** If @a loc is false, it is a list of true dofs in local ordering.
820  If @a loc is true, it is a marker for Vdofs in local ordering. */
821  void SetEssBdrDofs(const Array<int> *essdofs, bool loc = false)
822  {
823  ess_dof = essdofs;
824  ess_dof_local = loc;
825  }
826  /// Specify dofs on the natural boundary.
827  /** If @a loc is false, it is a list of true dofs in local ordering.
828  If @a loc is true, it is a marker for Vdofs in local ordering. */
829  void SetNatBdrDofs(const Array<int> *natdofs, bool loc = false)
830  {
831  nat_dof = natdofs;
832  nat_dof_local = loc;
833  }
834  /// Setup BDDC with no-net-flux local solvers. Needs a ParFiniteElementSpace attached
835  void SetComputeNetFlux(bool net = true)
836  {
837  netflux = net;
838  }
839 };
840 
841 
843 {
844 private:
845  void BDDCSolverConstructor(const PetscBDDCSolverParams &opts);
846 
847 public:
848  PetscBDDCSolver(MPI_Comm comm, Operator &op,
850  const std::string &prefix = std::string());
853  const std::string &prefix = std::string());
854 };
855 
856 
858 {
859 public:
860  PetscFieldSplitSolver(MPI_Comm comm, Operator &op,
861  const std::string &prefix = std::string());
862 };
863 
865 {
866 private:
867  void H2SolverConstructor(ParFiniteElementSpace *fes);
868 
869 public:
872  const std::string &prefix = std::string());
873 
874 };
875 
876 /// Abstract class for PETSc's nonlinear solvers.
878 {
879 public:
880  PetscNonlinearSolver(MPI_Comm comm,
881  const std::string &prefix = std::string());
882  PetscNonlinearSolver(MPI_Comm comm, Operator &op,
883  const std::string &prefix = std::string());
884  virtual ~PetscNonlinearSolver();
885 
886  /// Specification of the nonlinear operator.
887  virtual void SetOperator(const Operator &op);
888 
889  /// Specifies the desired format of the Jacobian in case a PetscParMatrix
890  /// is not returned by the GetGradient method.
891  void SetJacobianType(Operator::Type type);
892 
893  /// Application of the solver.
894  virtual void Mult(const Vector &b, Vector &x) const;
895 
896  /// Specification of an objective function to be used for line search.
897  void SetObjective(void (*obj)(Operator* op, const Vector &x, double *f));
898 
899  /// User-defined routine to be applied after a successful line search step.
900  /// The user can change the current direction Y and/or the updated solution W
901  /// (with W = X - lambda * Y) but not the previous solution X.
902  /// If Y or W have been changed, the corresponding booleans need to updated.
903  void SetPostCheck(void (*post)(Operator *op, const Vector &X, Vector &Y,
904  Vector &W, bool &changed_y, bool &changed_w));
905 
906  /// General purpose update function to be called at the beginning of each step
907  /// it is the current nonlinear iteration number
908  /// F is the current function value, X the current solution
909  /// D the previous step taken, and P the previous solution
910  void SetUpdate(void (*update)(Operator *op, int it,
911  const Vector& F, const Vector& X,
912  const Vector& D, const Vector& P));
913 
914  /// Conversion function to PETSc's SNES type.
915  operator petsc::SNES() const { return (petsc::SNES)obj; }
916 };
917 
918 
919 /// Abstract class for PETSc's ODE solvers.
920 class PetscODESolver : public PetscSolver, public ODESolver
921 {
922 public:
923  /// The type of the ODE. Use ODE_SOLVER_LINEAR if the Jacobians
924  /// are linear and independent of time.
925  enum Type
926  {
929  };
930 
931  PetscODESolver(MPI_Comm comm, const std::string &prefix = std::string());
932  virtual ~PetscODESolver();
933 
934  /// Initialize the ODE solver.
935  virtual void Init(TimeDependentOperator &f_,
936  enum PetscODESolver::Type type);
938 
941 
942  /// Specifies the desired format of the Jacobian in case a PetscParMatrix
943  /// is not returned by the GetGradient methods
944  void SetJacobianType(Operator::Type type);
945 
946  virtual void Step(Vector &x, double &t, double &dt);
947  virtual void Run(Vector &x, double &t, double &dt, double t_final);
948 
949  /// Conversion function to PETSc's TS type.
950  operator petsc::TS() const { return (petsc::TS)obj; }
951 };
952 
953 
954 /// Abstract class for monitoring PETSc's solvers.
956 {
957 public:
958  bool mon_sol;
959  bool mon_res;
960  PetscSolverMonitor(bool monitor_sol = false, bool monitor_res = true)
961  : mon_sol(monitor_sol), mon_res(monitor_res) {}
962  virtual ~PetscSolverMonitor() {}
963 
964  /// Monitor the solution vector x
965  virtual void MonitorSolution(PetscInt it, PetscReal norm, const Vector &x)
966  {
967  MFEM_ABORT("MonitorSolution() is not implemented!")
968  }
969 
970  /// Monitor the residual vector r
971  virtual void MonitorResidual(PetscInt it, PetscReal norm, const Vector &r)
972  {
973  MFEM_ABORT("MonitorResidual() is not implemented!")
974  }
975 
976  /// Generic monitor to take access to the solver
977  virtual void MonitorSolver(PetscSolver* solver) {}
978 };
979 
980 
981 } // namespace mfem
982 
983 #endif // MFEM_USE_MPI
984 #endif // MFEM_USE_PETSC
985 
986 #endif
void SetPreconditioner(Solver &precond)
Definition: petsc.cpp:3025
void SetPrintLevel(int plev)
Definition: petsc.cpp:2374
PetscParVector * B
Right-hand side and solution vector.
Definition: petsc.hpp:661
void CreatePrivateContext()
Definition: petsc.cpp:2673
void Print(const char *fname=NULL, bool binary=false) const
Prints the matrix (to stdout if fname is NULL)
Definition: petsc.cpp:1945
void trans(const Vector &u, Vector &x)
Definition: ex27.cpp:421
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:2857
void SetUpdate(void(*update)(Operator *op, int it, const Vector &F, const Vector &X, const Vector &D, const Vector &P))
Definition: petsc.cpp:4005
virtual void SetOperator(const Operator &op)
Specification of the nonlinear operator.
Definition: petsc.cpp:3912
Abstract class for PETSc&#39;s preconditioners.
Definition: petsc.hpp:775
PetscParMatrix & operator+=(const PetscParMatrix &B)
Definition: petsc.cpp:1078
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:2880
HypreParMatrix * RAP(const HypreParMatrix *A, const HypreParMatrix *P)
Returns the matrix P^t * A * P.
Definition: hypre.cpp:2485
bool clcustom
Boolean to handle SetFromOptions calls.
Definition: petsc.hpp:652
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:3980
HYPRE_Int PetscInt
Definition: petsc.hpp:47
double GetFinalNorm()
Definition: petsc.cpp:2648
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:3281
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:311
Wrapper for PETSc&#39;s matrix class.
Definition: petsc.hpp:307
PetscParVector & operator+=(const PetscParVector &y)
Definition: petsc.cpp:692
const Array< int > * nat_dof
Definition: petsc.hpp:807
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:2756
void SetHostValid() const
Definition: petsc.hpp:93
Abstract class for PETSc&#39;s linear solvers.
Definition: petsc.hpp:722
Abstract class for PETSc&#39;s solvers.
Definition: petsc.hpp:648
PetscODESolver::Type GetType() const
Definition: petsc.cpp:4163
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:4225
void SetPostCheck(void(*post)(Operator *op, const Vector &X, Vector &Y, Vector &W, bool &changed_y, bool &changed_w))
Definition: petsc.cpp:3991
PetscBCHandler(Type type_=ZERO)
Definition: petsc.hpp:574
virtual void Mult(const Vector &x, Vector &y) const
Operator application: y=A(x).
Definition: petsc.hpp:436
ParFiniteElementSpace * fespace
Definition: petsc.hpp:804
Type GetType() const
Returns the type of boundary conditions.
Definition: petsc.hpp:582
Base abstract class for first order time dependent operators.
Definition: operator.hpp:282
PetscInt GetNumCols() const
Returns the local number of columns.
Definition: petsc.cpp:921
void MakeAliasForSync(const Memory< double > &base_, int offset_, int size_, bool usedev_)
Definition: petsc.hpp:99
void SetType(PetscODESolver::Type)
Definition: petsc.cpp:4169
void ConvertOperator(MPI_Comm comm, const Operator &op, petsc::Mat *B, Operator::Type tid)
Definition: petsc.cpp:1294
virtual void MonitorSolver(PetscSolver *solver)
Generic monitor to take access to the solver.
Definition: petsc.hpp:977
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:1910
PetscBDDCSolver(MPI_Comm comm, Operator &op, const PetscBDDCSolverParams &opts=PetscBDDCSolverParams(), const std::string &prefix=std::string())
Definition: petsc.cpp:3783
PetscParVector & AddValues(const Array< PetscInt > &, const Array< PetscScalar > &)
Add values in a vector.
Definition: petsc.cpp:671
virtual ~PetscSolver()
Destroy the PetscParVectors allocated (if any).
Definition: petsc.cpp:2285
PetscNonlinearSolver(MPI_Comm comm, const std::string &prefix=std::string())
Definition: petsc.cpp:3872
bool DeviceRequested() const
Definition: petsc.hpp:140
PetscInt GetGlobalNumRows() const
Returns the global number of rows.
Definition: petsc.hpp:470
Type GetType() const
Definition: petsc.cpp:2227
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:4187
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:949
PetscPreconditioner(MPI_Comm comm, const std::string &prefix=std::string())
Definition: petsc.cpp:3155
PetscClassId cid
The class id of the actual PETSc object.
Definition: petsc.hpp:658
void SetDataAndSize_()
Definition: petsc.cpp:218
Abstract parallel finite element space.
Definition: pfespace.hpp:28
Constant in time b.c.
Definition: petsc.hpp:570
void SetRelTol(double tol)
Definition: petsc.cpp:2297
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:439
Wrapper for synching PETSc&#39;s vector memory.
Definition: petsc.hpp:84
void Randomize(PetscInt seed=0)
Set random values.
Definition: petsc.cpp:858
PetscPCGSolver(MPI_Comm comm, const std::string &prefix=std::string())
Definition: petsc.cpp:3125
void Destroy()
Delete all owned data. Does not perform re-initialization with defaults.
Definition: petsc.cpp:1690
void SyncBaseAndReset()
Definition: petsc.hpp:124
void SetMonitor(PetscSolverMonitor *ctx)
Sets user-defined monitoring routine.
Definition: petsc.cpp:2459
void SetMaxIter(int max_iter)
Definition: petsc.cpp:2347
Abstract class for PETSc&#39;s ODE solvers.
Definition: petsc.hpp:920
PetscParVector & operator*=(PetscScalar d)
Definition: petsc.cpp:706
HypreParMatrix * ParMult(const HypreParMatrix *A, const HypreParMatrix *B, bool own_matrix)
Definition: hypre.cpp:2464
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:2171
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:2549
PetscInt M() const
Returns the global number of rows.
Definition: petsc.cpp:928
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:4018
PetscParVector * Y
Definition: petsc.hpp:314
virtual ~PetscLinearSolver()
Definition: petsc.cpp:3115
PetscInt GetNumRows() const
Returns the local number of rows.
Definition: petsc.cpp:914
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:873
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:725
const double * GetHostPointer() const
Definition: petsc.cpp:198
PetscParVector & SetValues(const Array< PetscInt > &, const Array< PetscScalar > &)
Set values in a vector.
Definition: petsc.cpp:657
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:4074
PetscParMatrix * Transpose(bool action=false)
Returns the transpose of the PetscParMatrix.
Definition: petsc.cpp:1891
double f(const Vector &xvec)
Definition: lor_mms.hpp:32
PetscObject obj
The actual PETSc object (KSP, PC, SNES or TS).
Definition: petsc.hpp:655
Abstract class for monitoring PETSc&#39;s solvers.
Definition: petsc.hpp:955
Data type sparse matrix.
Definition: sparsemat.hpp:41
void SetDeviceValid() const
Definition: petsc.hpp:94
void SetJacobianType(Operator::Type type)
Definition: petsc.cpp:3974
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:153
double b
Definition: lissajous.cpp:42
Auxiliary class for BDDC customization.
Definition: petsc.hpp:801
virtual ~PetscPreconditioner()
Definition: petsc.cpp:3286
MPI_Comm GetComm() const
Get the associated MPI communicator.
Definition: petsc.cpp:1247
void FixResidualBC(const Vector &x, Vector &y)
y = x-g on ess_tdof_list, the rest of y is unchanged
Definition: petsc.cpp:2805
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:585
PetscParMatrix & operator=(const PetscParMatrix &B)
Definition: petsc.cpp:1062
bool ReadRequested() const
Definition: petsc.hpp:130
PetscInt GetColStart() const
Returns the global index of the first local column.
Definition: petsc.cpp:907
virtual ~PetscParVector()
Calls PETSc&#39;s destroy function.
Definition: petsc.cpp:518
PetscFieldSplitSolver(MPI_Comm comm, Operator &op, const std::string &prefix=std::string())
Definition: petsc.cpp:3792
PetscBCHandler * bchandler
Handler for boundary conditions.
Definition: petsc.hpp:664
PetscParMatrix & operator-=(const PetscParMatrix &B)
Definition: petsc.cpp:1093
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:1970
void Zero(Vector &x)
Replace boundary dofs with 0.
Definition: petsc.cpp:2824
MPI_Comm GetComm() const
Get the associated MPI communicator.
Definition: petsc.cpp:2454
bool WriteRequested() const
Definition: petsc.hpp:135
void Shift(double s)
Shift diagonal by a constant.
Definition: petsc.cpp:1992
double * HostReadWrite() override
Shortcut for mfem::ReadWrite(vec.GetMemory(), vec.Size(), false).
Definition: petsc.cpp:435
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:600
ID for class PetscParMatrix, MATSHELL format.
Definition: operator.hpp:259
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:1715
PetscParMatrix()
Create an empty matrix to be used as a reference to an existing matrix.
Definition: petsc.cpp:956
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:1260
Abstract class for PETSc&#39;s nonlinear solvers.
Definition: petsc.hpp:877
int GetNumIterations()
Definition: petsc.cpp:2615
void SetBCHandler(PetscBCHandler *bch)
Sets the object to handle essential boundary conditions.
Definition: petsc.cpp:2489
struct::_p_KSP * KSP
Definition: petsc.hpp:71
Wrapper for hypre&#39;s parallel vector class.
Definition: hypre.hpp:99
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:2833
virtual void Eval(double t, Vector &g)
Boundary conditions evaluation.
Definition: petsc.hpp:590
Type
Enumeration defining IDs for some classes derived from Operator.
Definition: operator.hpp:252
struct::_p_SNES * SNES
Definition: petsc.hpp:73
PetscInt NNZ() const
Returns the number of nonzeros.
Definition: petsc.cpp:942
PetscPreconditionerFactory(const std::string &name_="MFEM Factory")
Definition: petsc.hpp:637
PetscLinearSolver(MPI_Comm comm, const std::string &prefix=std::string(), bool wrap=true)
Definition: petsc.cpp:2845
virtual void Mult(const Vector &b, Vector &x) const
Application of the preconditioner.
Definition: petsc.cpp:3276
virtual void Init(TimeDependentOperator &f_)
Associate a TimeDependentOperator with the ODE solver.
Definition: petsc.hpp:937
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:935
PetscParVector * GetY() const
Returns the inner vector in the range of A (it creates it if needed)
Definition: petsc.cpp:1881
void SetJacobianType(Operator::Type type)
Definition: petsc.cpp:4157
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:670
PetscParVector * X
Auxiliary vectors for typecasting.
Definition: petsc.hpp:314
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:1927
bool IsAliasForSync() const
Definition: petsc.hpp:97
Pointer is an alias.
PetscParVector & operator-=(const PetscParVector &y)
Definition: petsc.cpp:699
virtual ~PetscSolverMonitor()
Definition: petsc.hpp:962
double a
Definition: lissajous.cpp:41
void MFEMFinalizePetsc()
Definition: petsc.cpp:192
virtual ~PetscNonlinearSolver()
Definition: petsc.cpp:3904
void SetComputeNetFlux(bool net=true)
Setup BDDC with no-net-flux local solvers. Needs a ParFiniteElementSpace attached.
Definition: petsc.hpp:835
virtual ~PetscODESolver()
Definition: petsc.cpp:4066
virtual ~PetscBCHandler()
Definition: petsc.hpp:579
int GetConverged()
Definition: petsc.cpp:2582
PetscSolver()
Construct an empty PetscSolver. Initialize protected objects to NULL.
Definition: petsc.cpp:2275
void SetPreconditionerFactory(PetscPreconditionerFactory *factory)
Sets the object for the creation of the preconditioner.
Definition: petsc.cpp:2508
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:2214
PetscParVector(MPI_Comm comm, PetscInt glob_size, PetscInt *col=NULL)
Creates vector with given global size and partitioning of the columns.
Definition: petsc.cpp:500
virtual void MonitorSolution(PetscInt it, PetscReal norm, const Vector &x)
Monitor the solution vector x.
Definition: petsc.hpp:965
ID for class PetscParMatrix, MATAIJ format.
Definition: operator.hpp:257
virtual void SetOperator(const Operator &op)
Set/update the solver for the given operator.
Definition: petsc.cpp:3190
void SetSpace(ParFiniteElementSpace *fe)
Definition: petsc.hpp:816
MPI_Comm GetComm() const
Get the associated MPI communicator.
Definition: petsc.cpp:620
PetscInt GlobalSize() const
Returns the global number of rows.
Definition: petsc.cpp:456
PetscParVector * X
Definition: petsc.hpp:661
const Array< int > * ess_dof
Definition: petsc.hpp:805
PetscInt GetRowStart() const
Returns the global index of the first local row.
Definition: petsc.cpp:900
struct _p_PetscObject * PetscObject
Definition: petsc.hpp:51
PetscH2Solver(Operator &op, ParFiniteElementSpace *fes, const std::string &prefix=std::string())
Definition: petsc.cpp:3827
void ResetMemory()
Completes the operation started with PlaceMemory.
Definition: petsc.cpp:814
void SetAbsTol(double tol)
Definition: petsc.cpp:2322
Helper class for handling essential boundary conditions.
Definition: petsc.hpp:564
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:1871
Array< int > & GetTDofs()
Gets essential dofs (local, per-process numbering)
Definition: petsc.hpp:597
void SetTol(double tol)
Definition: petsc.cpp:2292
PetscParVector & operator=(PetscScalar d)
Set constant values.
Definition: petsc.cpp:650
PetscSolverMonitor(bool monitor_sol=false, bool monitor_res=true)
Definition: petsc.hpp:960
RefCoord t[3]
virtual void MonitorResidual(PetscInt it, PetscReal norm, const Vector &r)
Monitor the residual vector r.
Definition: petsc.hpp:971
PetscParMatrix * TripleMatrixProduct(PetscParMatrix *R, PetscParMatrix *A, PetscParMatrix *P)
Returns the matrix R * A * P.
Definition: petsc.cpp:2011
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:720
void FreePrivateContext()
Definition: petsc.cpp:2705
RefCoord s[3]
Vector * GlobalVector() const
Returns the global vector in each processor.
Definition: petsc.cpp:625
PetscODESolver(MPI_Comm comm, const std::string &prefix=std::string())
Definition: petsc.cpp:4041
void EliminateRows(const Array< int > &rows)
Eliminate only the rows from the matrix.
Definition: petsc.cpp:2200
Base class for solvers.
Definition: operator.hpp:648
PetscInt GetGlobalNumCols() const
Returns the global number of columns.
Definition: petsc.hpp:473
void SetNatBdrDofs(const Array< int > *natdofs, bool loc=false)
Specify dofs on the natural boundary.
Definition: petsc.hpp:829
Abstract operator.
Definition: operator.hpp:24
Wrapper for hypre&#39;s ParCSR matrix class.
Definition: hypre.hpp:277
virtual void Mult(const Vector &b, Vector &x) const
Application of the solver.
Definition: petsc.cpp:3105
void operator*=(double s)
Scale all entries by s: A_scaled = s*A.
Definition: petsc.cpp:1905
void SetEssBdrDofs(const Array< int > *essdofs, bool loc=false)
Specify dofs on the essential boundary.
Definition: petsc.hpp:821
void SetTDofs(Array< int > &list)
Sets essential dofs (local, per-process numbering)
Definition: petsc.cpp:2733
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:3110
void PlaceMemory(Memory< double > &, bool=false)
This requests write access from where the memory is valid and temporarily replaces the corresponding ...
Definition: petsc.cpp:730
void * private_ctx
Private context for solver.
Definition: petsc.hpp:667
void SetUp(PetscInt n)
SetUp the helper object, where n is the size of the solution vector.
Definition: petsc.cpp:2740
const double * GetDevicePointer() const
Definition: petsc.cpp:207
void MakeRef(const PetscParMatrix &master)
Makes this object a reference to another PetscParMatrix.
Definition: petsc.cpp:1861
void ScaleCols(const Vector &s)
Scale the local col i by s(i).
Definition: petsc.cpp:1981
virtual ~PetscParMatrix()
Calls PETSc&#39;s destroy function.
Definition: petsc.hpp:417