MFEM  v3.0
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Pages
hypre.hpp
Go to the documentation of this file.
1 // Copyright (c) 2010, Lawrence Livermore National Security, LLC. Produced at
2 // the Lawrence Livermore National Laboratory. LLNL-CODE-443211. All Rights
3 // reserved. See file COPYRIGHT for details.
4 //
5 // This file is part of the MFEM library. For more information and source code
6 // availability see http://mfem.googlecode.com.
7 //
8 // MFEM is free software; you can redistribute it and/or modify it under the
9 // terms of the GNU Lesser General Public License (as published by the Free
10 // Software Foundation) version 2.1 dated February 1999.
11 
12 #ifndef MFEM_HYPRE
13 #define MFEM_HYPRE
14 
15 #include "../config/config.hpp"
16 
17 #ifdef MFEM_USE_MPI
18 
19 #include <mpi.h>
20 
21 // Enable internal hypre timing routines
22 #define HYPRE_TIMING
23 
24 // hypre header files
25 #include "seq_mv.h"
26 #include "_hypre_parcsr_mv.h"
27 #include "_hypre_parcsr_ls.h"
28 
29 #include "sparsemat.hpp"
30 
31 namespace mfem
32 {
33 
34 class ParFiniteElementSpace;
35 class HypreParMatrix;
36 
38 class HypreParVector : public Vector
39 {
40 private:
41  int own_ParVector;
42 
44  hypre_ParVector *x;
45 
46  friend class HypreParMatrix;
47 
48 public:
51  HypreParVector(MPI_Comm comm, int glob_size, int *col);
54  HypreParVector(MPI_Comm comm, int glob_size, double *_data, int *col);
58  HypreParVector(HypreParMatrix &A, int tr = 0);
60  HypreParVector(HYPRE_ParVector y);
63 
65  operator hypre_ParVector*() const;
66 #ifndef HYPRE_PAR_VECTOR_STRUCT
67  operator HYPRE_ParVector() const;
69 #endif
70  hypre_ParVector *StealParVector() { own_ParVector = 0; return x; }
72 
75 
77  HypreParVector& operator= (double d);
80 
85  void SetData(double *_data);
86 
88  int Randomize(int seed);
89 
91  void Print(const char *fname);
92 
95 };
96 
98 double InnerProduct(HypreParVector &x, HypreParVector &y);
99 double InnerProduct(HypreParVector *x, HypreParVector *y);
100 
101 
103 class HypreParMatrix : public Operator
104 {
105 private:
107  hypre_ParCSRMatrix *A;
108 
110  hypre_ParCSRCommPkg *CommPkg;
111 
113  mutable HypreParVector *X, *Y;
114 
115 public:
117  HypreParMatrix(hypre_ParCSRMatrix *a) : A(a)
118  { height = GetNumRows(); width = GetNumCols(); X = Y = 0; CommPkg = 0; }
120  HypreParMatrix(MPI_Comm comm, int glob_size, int *row_starts,
121  SparseMatrix *diag);
124  HypreParMatrix(MPI_Comm comm, int global_num_rows, int global_num_cols,
125  int *row_starts, int *col_starts, SparseMatrix *diag);
127  HypreParMatrix(MPI_Comm comm, int global_num_rows, int global_num_cols,
128  int *row_starts, int *col_starts,
129  SparseMatrix *diag, SparseMatrix *offd, int *cmap);
130 
132  HypreParMatrix(MPI_Comm comm, int *row_starts, int *col_starts,
133  SparseMatrix *a);
134 
136  HypreParMatrix(MPI_Comm comm, int global_num_rows, int global_num_cols,
137  int *row_starts, int *col_starts, Table *diag);
139  HypreParMatrix(MPI_Comm comm, int id, int np, int *row, int *col,
140  int *i_diag, int *j_diag, int *i_offd, int *j_offd,
141  int *cmap, int cmap_size);
142 
147  HypreParMatrix(MPI_Comm comm, int nrows, int glob_nrows, int glob_ncols,
148  int *I, int *J, double *data, int *rows, int *cols);
149 
150  // hypre's communication package object
151  void SetCommPkg(hypre_ParCSRCommPkg *comm_pkg);
152  void CheckCommPkg();
153  void DestroyCommPkg();
154 
156  MPI_Comm GetComm() { return A->comm; }
157 
159  operator hypre_ParCSRMatrix*();
160 #ifndef HYPRE_PAR_CSR_MATRIX_STRUCT
161  operator HYPRE_ParCSRMatrix();
163 #endif
164  hypre_ParCSRMatrix* StealData();
166 
168  inline int NNZ() { return A->num_nonzeros; }
170  inline int * RowPart() { return A->row_starts; }
172  inline int * ColPart() { return A->col_starts; }
174  inline int M() { return A -> global_num_rows; }
176  inline int N() { return A -> global_num_cols; }
177 
179  void GetDiag(Vector &diag);
182 
184  int GetNumRows() const
185  { return hypre_CSRMatrixNumRows(hypre_ParCSRMatrixDiag(A)); }
186 
188  int GetNumCols() const
189  { return hypre_CSRMatrixNumCols(hypre_ParCSRMatrixDiag(A)); }
190 
191  int GetGlobalNumRows() const { return hypre_ParCSRMatrixGlobalNumRows(A); }
192 
193  int GetGlobalNumCols() const { return hypre_ParCSRMatrixGlobalNumCols(A); }
194 
195  int *GetRowStarts() const { return hypre_ParCSRMatrixRowStarts(A); }
196 
197  int *GetColStarts() const { return hypre_ParCSRMatrixColStarts(A); }
198 
201  double alpha = 1.0, double beta = 0.0);
203  int Mult(HYPRE_ParVector x, HYPRE_ParVector y,
204  double alpha = 1.0, double beta = 0.0);
207  double alpha = 1.0, double beta = 0.0);
208 
209  void Mult(double a, const Vector &x, double b, Vector &y) const;
210  void MultTranspose(double a, const Vector &x, double b, Vector &y) const;
211 
212  virtual void Mult(const Vector &x, Vector &y) const
213  { Mult(1.0, x, 0.0, y); }
214  virtual void MultTranspose(const Vector &x, Vector &y) const
215  { MultTranspose(1.0, x, 0.0, y); }
216 
218  void ScaleRows(const Vector & s);
220  void InvScaleRows(const Vector & s);
222  void operator*=(double s);
223 
225  void Print(const char *fname, int offi = 0, int offj = 0);
227  void Read(MPI_Comm comm, const char *fname);
228 
230  virtual ~HypreParMatrix();
231 };
232 
234 HypreParMatrix * ParMult(HypreParMatrix *A, HypreParMatrix *B);
235 
237 HypreParMatrix * RAP(HypreParMatrix *A, HypreParMatrix *P);
239 HypreParMatrix * RAP(HypreParMatrix * Rt, HypreParMatrix *A, HypreParMatrix *P);
240 
244 void EliminateBC(HypreParMatrix &A, HypreParMatrix &Ae,
245  Array<int> &ess_dof_list,
246  HypreParVector &x, HypreParVector &b);
247 
248 
250 class HypreSmoother : public Solver
251 {
252 protected:
256  mutable HypreParVector *B, *X;
258  mutable HypreParVector *V, *Z;
260  mutable HypreParVector *X0, *X1;
261 
264  int type;
268  double relax_weight;
270  double omega;
277 
279  double lambda;
280  double mu;
282 
284  double *l1_norms;
286  double max_eig_est;
288  double min_eig_est;
290  double window_params[3];
291 
293  double* fir_coeffs;
294 
295 public:
306  enum Type { Jacobi = 0, l1Jacobi = 1, l1GS = 2, l1GStr = 4, lumpedJacobi = 5,
307  GS = 6, Chebyshev = 16, Taubin = 1001, FIR = 1002 };
308 
309  HypreSmoother();
310 
312  int relax_times = 1, double relax_weight = 1.0,
313  double omega = 1.0, int poly_order = 2,
314  double poly_fraction = .3);
315 
319  void SetSOROptions(double relax_weight, double omega);
321  void SetPolyOptions(int poly_order, double poly_fraction);
323  void SetTaubinOptions(double lambda, double mu, int iter);
324 
326  void SetWindowByName(const char* window_name);
328  void SetWindowParameters(double a, double b, double c);
330  void SetFIRCoefficients(double max_eig);
331 
334  virtual void SetOperator(const Operator &op);
335 
337  virtual void Mult(const HypreParVector &b, HypreParVector &x) const;
338  virtual void Mult(const Vector &b, Vector &x) const;
339 
340  virtual ~HypreSmoother();
341 };
342 
343 
345 class HypreSolver : public Solver
346 {
347 protected:
350 
352  mutable HypreParVector *B, *X;
353 
355  mutable int setup_called;
356 
357 public:
358  HypreSolver();
359 
361 
363  virtual operator HYPRE_Solver() const = 0;
364 
366  virtual HYPRE_PtrToParSolverFcn SetupFcn() const = 0;
368  virtual HYPRE_PtrToParSolverFcn SolveFcn() const = 0;
369 
370  virtual void SetOperator(const Operator &op)
371  { mfem_error("HypreSolvers do not support SetOperator!"); }
372 
374  virtual void Mult(const HypreParVector &b, HypreParVector &x) const;
375  virtual void Mult(const Vector &b, Vector &x) const;
376 
377  virtual ~HypreSolver();
378 };
379 
381 class HyprePCG : public HypreSolver
382 {
383 private:
384  int print_level;
385  HYPRE_Solver pcg_solver;
386 
387 public:
389 
390  void SetTol(double tol);
391  void SetMaxIter(int max_iter);
392  void SetLogging(int logging);
393  void SetPrintLevel(int print_lvl);
394 
396  void SetPreconditioner(HypreSolver &precond);
397 
401  void SetResidualConvergenceOptions(int res_frequency=-1, double rtol=0.0);
402 
405 
406  void GetNumIterations(int &num_iterations)
407  { HYPRE_ParCSRPCGGetNumIterations(pcg_solver, &num_iterations); }
408 
410  virtual operator HYPRE_Solver() const { return pcg_solver; }
411 
413  virtual HYPRE_PtrToParSolverFcn SetupFcn() const
414  { return (HYPRE_PtrToParSolverFcn) HYPRE_ParCSRPCGSetup; }
416  virtual HYPRE_PtrToParSolverFcn SolveFcn() const
417  { return (HYPRE_PtrToParSolverFcn) HYPRE_ParCSRPCGSolve; }
418 
420  virtual void Mult(const HypreParVector &b, HypreParVector &x) const;
421  using HypreSolver::Mult;
422 
423  virtual ~HyprePCG();
424 };
425 
427 class HypreGMRES : public HypreSolver
428 {
429 private:
430  int print_level;
431  HYPRE_Solver gmres_solver;
432 
433 public:
435 
436  void SetTol(double tol);
437  void SetMaxIter(int max_iter);
438  void SetKDim(int dim);
439  void SetLogging(int logging);
440  void SetPrintLevel(int print_lvl);
441 
443  void SetPreconditioner(HypreSolver &precond);
444 
447 
449  virtual operator HYPRE_Solver() const { return gmres_solver; }
450 
452  virtual HYPRE_PtrToParSolverFcn SetupFcn() const
453  { return (HYPRE_PtrToParSolverFcn) HYPRE_ParCSRGMRESSetup; }
455  virtual HYPRE_PtrToParSolverFcn SolveFcn() const
456  { return (HYPRE_PtrToParSolverFcn) HYPRE_ParCSRGMRESSolve; }
457 
459  virtual void Mult (const HypreParVector &b, HypreParVector &x) const;
460  using HypreSolver::Mult;
461 
462  virtual ~HypreGMRES();
463 };
464 
467 {
468 public:
469  virtual operator HYPRE_Solver() const { return NULL; }
470 
471  virtual HYPRE_PtrToParSolverFcn SetupFcn() const
472  { return (HYPRE_PtrToParSolverFcn) hypre_ParKrylovIdentitySetup; }
473  virtual HYPRE_PtrToParSolverFcn SolveFcn() const
474  { return (HYPRE_PtrToParSolverFcn) hypre_ParKrylovIdentity; }
475 
476  virtual ~HypreIdentity() { }
477 };
478 
481 {
482 public:
485  virtual operator HYPRE_Solver() const { return NULL; }
486 
487  virtual HYPRE_PtrToParSolverFcn SetupFcn() const
488  { return (HYPRE_PtrToParSolverFcn) HYPRE_ParCSRDiagScaleSetup; }
489  virtual HYPRE_PtrToParSolverFcn SolveFcn() const
490  { return (HYPRE_PtrToParSolverFcn) HYPRE_ParCSRDiagScale; }
491 
492  HypreParMatrix* GetData() { return A; }
493  virtual ~HypreDiagScale() { }
494 };
495 
498 {
499 private:
500  HYPRE_Solver sai_precond;
501 
502 public:
504 
505  void SetSymmetry(int sym);
506 
508  virtual operator HYPRE_Solver() const { return sai_precond; }
509 
510  virtual HYPRE_PtrToParSolverFcn SetupFcn() const
511  { return (HYPRE_PtrToParSolverFcn) HYPRE_ParaSailsSetup; }
512  virtual HYPRE_PtrToParSolverFcn SolveFcn() const
513  { return (HYPRE_PtrToParSolverFcn) HYPRE_ParaSailsSolve; }
514 
515  virtual ~HypreParaSails();
516 };
517 
520 {
521 private:
522  HYPRE_Solver amg_precond;
523 
524 public:
526 
530  void SetSystemsOptions(int dim);
531 
532  void SetPrintLevel(int print_level)
533  { HYPRE_BoomerAMGSetPrintLevel(amg_precond, print_level); }
534 
536  virtual operator HYPRE_Solver() const { return amg_precond; }
537 
538  virtual HYPRE_PtrToParSolverFcn SetupFcn() const
539  { return (HYPRE_PtrToParSolverFcn) HYPRE_BoomerAMGSetup; }
540  virtual HYPRE_PtrToParSolverFcn SolveFcn() const
541  { return (HYPRE_PtrToParSolverFcn) HYPRE_BoomerAMGSolve; }
542 
543  virtual ~HypreBoomerAMG();
544 };
545 
547 HypreParMatrix* DiscreteGrad(ParFiniteElementSpace *edge_fespace,
548  ParFiniteElementSpace *vert_fespace);
550 HypreParMatrix* DiscreteCurl(ParFiniteElementSpace *face_fespace,
551  ParFiniteElementSpace *edge_fespace);
552 
554 class HypreAMS : public HypreSolver
555 {
556 private:
557  HYPRE_Solver ams;
558 
560  HypreParVector *x, *y, *z;
562  HypreParMatrix *G;
564  HypreParMatrix *Pi, *Pix, *Piy, *Piz;
565 
566 public:
568 
570  virtual operator HYPRE_Solver() const { return ams; }
571 
572  virtual HYPRE_PtrToParSolverFcn SetupFcn() const
573  { return (HYPRE_PtrToParSolverFcn) HYPRE_AMSSetup; }
574  virtual HYPRE_PtrToParSolverFcn SolveFcn() const
575  { return (HYPRE_PtrToParSolverFcn) HYPRE_AMSSolve; }
576 
577  virtual ~HypreAMS();
578 };
579 
581 class HypreADS : public HypreSolver
582 {
583 private:
584  HYPRE_Solver ads;
585 
587  HypreParVector *x, *y, *z;
589  HypreParMatrix *G;
591  HypreParMatrix *C;
593  HypreParMatrix *ND_Pi, *ND_Pix, *ND_Piy, *ND_Piz;
595  HypreParMatrix *RT_Pi, *RT_Pix, *RT_Piy, *RT_Piz;
596 
597 public:
599 
601  virtual operator HYPRE_Solver() const { return ads; }
602 
603  virtual HYPRE_PtrToParSolverFcn SetupFcn() const
604  { return (HYPRE_PtrToParSolverFcn) HYPRE_ADSSetup; }
605  virtual HYPRE_PtrToParSolverFcn SolveFcn() const
606  { return (HYPRE_PtrToParSolverFcn) HYPRE_ADSSolve; }
607 
608  virtual ~HypreADS();
609 };
610 
611 }
612 
613 #endif // MFEM_USE_MPI
614 
615 #endif
virtual ~HypreBoomerAMG()
Definition: hypre.cpp:1637
void SetTol(double tol)
Definition: hypre.cpp:1339
void SetPreconditioner(HypreSolver &precond)
Set the hypre solver to be used as a preconditioner.
Definition: hypre.cpp:1491
void DestroyCommPkg()
Definition: hypre.cpp:508
The Auxiliary-space Maxwell Solver in hypre.
Definition: hypre.hpp:554
HypreParVector * X0
FIR Filter Temporary Vectors.
Definition: hypre.hpp:260
double min_eig_est
Minimal eigenvalue estimate for polynomial smoothing.
Definition: hypre.hpp:288
virtual HYPRE_PtrToParSolverFcn SolveFcn() const
hypre&#39;s internal Solve function
Definition: hypre.hpp:473
The Auxiliary-space Divergence Solver in hypre.
Definition: hypre.hpp:581
int setup_called
Was hypre&#39;s Setup function called already?
Definition: hypre.hpp:355
HypreParVector * B
Right-hand side and solution vector.
Definition: hypre.hpp:352
HypreDiagScale(HypreParMatrix &A)
Definition: hypre.hpp:484
double window_params[3]
Paramters for windowing function of FIR filter.
Definition: hypre.hpp:290
virtual void Mult(const Vector &x, Vector &y) const
Operator application.
Definition: hypre.hpp:212
virtual HYPRE_PtrToParSolverFcn SolveFcn() const
hypre&#39;s internal Solve function
Definition: hypre.hpp:574
int * ColPart()
Returns the column partitioning.
Definition: hypre.hpp:172
HypreADS(HypreParMatrix &A, ParFiniteElementSpace *face_fespace)
Definition: hypre.cpp:1786
void SetWindowByName(const char *window_name)
Convenience function for setting canonical windowing parameters.
Definition: hypre.cpp:1037
HypreParMatrix(hypre_ParCSRMatrix *a)
Converts hypre&#39;s format to HypreParMatrix.
Definition: hypre.hpp:117
virtual ~HypreGMRES()
Definition: hypre.cpp:1562
Abstract parallel finite element space.
Definition: pfespace.hpp:28
bool iterative_mode
If true, use the second argument of Mult as an initial guess.
Definition: operator.hpp:106
int M()
Returns the global number of rows.
Definition: hypre.hpp:174
virtual ~HypreAMS()
Definition: hypre.cpp:1771
void SetTol(double tol)
Definition: hypre.cpp:1465
virtual HYPRE_PtrToParSolverFcn SolveFcn() const
GMRES Solve function.
Definition: hypre.hpp:455
void Print(const char *fname, int offi=0, int offj=0)
Prints the locally owned rows in parallel.
Definition: hypre.cpp:716
void SetPrintLevel(int print_lvl)
Definition: hypre.cpp:1485
void ScaleRows(const Vector &s)
Scale the local row i by s(i).
Definition: hypre.cpp:637
int Mult(HypreParVector &x, HypreParVector &y, double alpha=1.0, double beta=0.0)
Computes y = alpha * A * x + beta * y.
Definition: hypre.cpp:572
double poly_fraction
Fraction of spectrum to smooth for polynomial relaxation.
Definition: hypre.hpp:274
void SetPrintLevel(int print_lvl)
Definition: hypre.cpp:1354
HypreParMatrix * DiscreteGrad(ParFiniteElementSpace *edge_fespace, ParFiniteElementSpace *vert_fespace)
Compute the discrete gradient matrix between the nodal linear and ND1 spaces.
int poly_scale
Apply the polynomial smoother to A or D^{-1/2} A D^{-1/2}.
Definition: hypre.hpp:276
HypreParVector(MPI_Comm comm, int glob_size, int *col)
Definition: hypre.cpp:29
HypreParVector * X1
Definition: hypre.hpp:260
The identity operator as a hypre solver.
Definition: hypre.hpp:466
void SetWindowParameters(double a, double b, double c)
Set parameters for windowing function for FIR smoother.
Definition: hypre.cpp:1050
virtual void Mult(const HypreParVector &b, HypreParVector &x) const
Solve Ax=b with hypre&#39;s GMRES.
Definition: hypre.cpp:1499
void SetSystemsOptions(int dim)
Definition: hypre.cpp:1628
virtual void Mult(const HypreParVector &b, HypreParVector &x) const
Relax the linear system Ax=b.
Definition: hypre.cpp:1159
hypre_ParVector * StealParVector()
Changes the ownership of the the vector.
Definition: hypre.hpp:71
The BoomerAMG solver in hypre.
Definition: hypre.hpp:519
void SetSymmetry(int sym)
Definition: hypre.cpp:1591
virtual void MultTranspose(const Vector &x, Vector &y) const
Action of the transpose operator.
Definition: hypre.hpp:214
virtual ~HypreParaSails()
Definition: hypre.cpp:1596
The ParaSails preconditioner in hypre.
Definition: hypre.hpp:497
void Print(const char *fname)
Prints the locally owned rows in parallel.
Definition: hypre.cpp:155
int MultTranspose(HypreParVector &x, HypreParVector &y, double alpha=1.0, double beta=0.0)
Computes y = alpha * A^t * x + beta * y.
Definition: hypre.cpp:631
MPI_Comm GetComm()
MPI communicator.
Definition: hypre.hpp:156
double * l1_norms
l1 norms of the rows of A
Definition: hypre.hpp:284
Data type sparse matrix.
Definition: sparsemat.hpp:38
int * RowPart()
Returns the row partitioning.
Definition: hypre.hpp:170
void SetLogging(int logging)
Definition: hypre.cpp:1349
void SetZeroInintialIterate()
non-hypre setting
Definition: hypre.hpp:404
Jacobi preconditioner in hypre.
Definition: hypre.hpp:480
virtual ~HypreSolver()
Definition: hypre.cpp:1320
HyprePCG(HypreParMatrix &_A)
Definition: hypre.cpp:1327
virtual HYPRE_PtrToParSolverFcn SolveFcn() const =0
hypre&#39;s internal Solve function
void SetKDim(int dim)
Definition: hypre.cpp:1475
void SetResidualConvergenceOptions(int res_frequency=-1, double rtol=0.0)
Definition: hypre.cpp:1368
void SetPrintLevel(int print_level)
Definition: hypre.hpp:532
virtual void Mult(const HypreParVector &b, HypreParVector &x) const
Solve the linear system Ax=b.
Definition: hypre.cpp:1275
HypreParMatrix * RAP(HypreParMatrix *A, HypreParMatrix *P)
Returns the matrix P^t * A * P.
Definition: hypre.cpp:781
virtual HYPRE_PtrToParSolverFcn SolveFcn() const
PCG Solve function.
Definition: hypre.hpp:416
int Randomize(int seed)
Set random values.
Definition: hypre.cpp:150
HypreParMatrix * DiscreteCurl(ParFiniteElementSpace *face_fespace, ParFiniteElementSpace *edge_fespace)
Compute the discrete curl matrix between the ND1 and RT0 spaces.
double relax_weight
Damping coefficient (usually &lt;= 1)
Definition: hypre.hpp:268
Parallel smoothers in hypre.
Definition: hypre.hpp:250
virtual HYPRE_PtrToParSolverFcn SetupFcn() const
hypre&#39;s internal Setup function
Definition: hypre.hpp:538
void EliminateBC(HypreParMatrix &A, HypreParMatrix &Ae, Array< int > &ess_dof_list, HypreParVector &x, HypreParVector &b)
Definition: hypre.cpp:827
virtual HYPRE_PtrToParSolverFcn SetupFcn() const
hypre&#39;s internal Setup function
Definition: hypre.hpp:471
HypreParMatrix * A
The linear system matrix.
Definition: hypre.hpp:254
HypreAMS(HypreParMatrix &A, ParFiniteElementSpace *edge_fespace)
Definition: hypre.cpp:1643
virtual HYPRE_PtrToParSolverFcn SolveFcn() const
hypre&#39;s internal Solve function
Definition: hypre.hpp:489
virtual void SetOperator(const Operator &op)
Definition: hypre.cpp:1057
virtual HYPRE_PtrToParSolverFcn SetupFcn() const =0
hypre&#39;s internal Setup function
void SetLogging(int logging)
Definition: hypre.cpp:1480
virtual ~HypreADS()
Definition: hypre.cpp:1955
virtual HYPRE_PtrToParSolverFcn SetupFcn() const
hypre&#39;s internal Setup function
Definition: hypre.hpp:487
void SetMaxIter(int max_iter)
Definition: hypre.cpp:1344
HypreParMatrix * GetData()
Definition: hypre.hpp:492
void SetSOROptions(double relax_weight, double omega)
Set SOR-related parameters.
Definition: hypre.cpp:1017
Wrapper for hypre&#39;s parallel vector class.
Definition: hypre.hpp:38
HypreParVector * X
Definition: hypre.hpp:352
virtual HYPRE_PtrToParSolverFcn SetupFcn() const
hypre&#39;s internal Setup function
Definition: hypre.hpp:603
int N()
Returns the global number of columns.
Definition: hypre.hpp:176
int * GetColStarts() const
Definition: hypre.hpp:197
PCG solver in hypre.
Definition: hypre.hpp:381
GMRES solver in hypre.
Definition: hypre.hpp:427
Vector * GlobalVector()
Returns the global vector in each processor.
Definition: hypre.cpp:117
HypreParVector * X
Definition: hypre.hpp:256
void mfem_error(const char *msg)
Definition: error.cpp:23
virtual HYPRE_PtrToParSolverFcn SolveFcn() const
hypre&#39;s internal Solve function
Definition: hypre.hpp:512
int GetGlobalNumCols() const
Definition: hypre.hpp:193
virtual ~HypreDiagScale()
Definition: hypre.hpp:493
int GetNumRows() const
Returns the number of rows in the diagonal block of the ParCSRMatrix.
Definition: hypre.hpp:184
virtual void SetOperator(const Operator &op)
Set/update the solver for the given operator.
Definition: hypre.hpp:370
double max_eig_est
Maximal eigenvalue estimate for polynomial smoothing.
Definition: hypre.hpp:286
HypreParMatrix * Transpose()
Returns the transpose of *this.
Definition: hypre.cpp:561
virtual HYPRE_PtrToParSolverFcn SetupFcn() const
hypre&#39;s internal Setup function
Definition: hypre.hpp:510
virtual ~HypreIdentity()
Definition: hypre.hpp:476
int NNZ()
Returns the number of nonzeros.
Definition: hypre.hpp:168
void SetFIRCoefficients(double max_eig)
Compute window and Chebyshev coefficients for given polynomial order.
Definition: hypre.cpp:1121
double * fir_coeffs
Combined coefficients for windowing and Chebyshev polynomials.
Definition: hypre.hpp:293
double InnerProduct(HypreParVector *x, HypreParVector *y)
Definition: hypre.cpp:167
virtual ~HypreSmoother()
Definition: hypre.cpp:1245
int * GetRowStarts() const
Definition: hypre.hpp:195
void SetData(double *_data)
Definition: hypre.cpp:145
void InvScaleRows(const Vector &s)
Scale the local row i by 1./s(i)
Definition: hypre.cpp:665
void SetTaubinOptions(double lambda, double mu, int iter)
Set parameters for Taubin&#39;s lambda-mu method.
Definition: hypre.cpp:1029
HypreParVector * B
Right-hand side and solution vectors.
Definition: hypre.hpp:256
void SetCommPkg(hypre_ParCSRCommPkg *comm_pkg)
Definition: hypre.cpp:490
void SetZeroInintialIterate()
non-hypre setting
Definition: hypre.hpp:446
virtual ~HyprePCG()
Definition: hypre.cpp:1440
int GetNumCols() const
Returns the number of columns in the diagonal block of the ParCSRMatrix.
Definition: hypre.hpp:188
void SetPreconditioner(HypreSolver &precond)
Set the hypre solver to be used as a preconditioner.
Definition: hypre.cpp:1360
void SetPolyOptions(int poly_order, double poly_fraction)
Set parameters for polynomial smoothing.
Definition: hypre.cpp:1023
void GetNumIterations(int &num_iterations)
Definition: hypre.hpp:406
int relax_times
Number of relaxation sweeps.
Definition: hypre.hpp:266
void operator*=(double s)
Scale all entries by s: A_scaled = s*A.
Definition: hypre.cpp:697
Abstract class for hypre&#39;s solvers and preconditioners.
Definition: hypre.hpp:345
int GetGlobalNumRows() const
Definition: hypre.hpp:191
HypreParVector & operator=(double d)
Set constant values.
Definition: hypre.cpp:127
Vector data type.
Definition: vector.hpp:29
void GetDiag(Vector &diag)
Get the diagonal of the matrix.
Definition: hypre.cpp:546
HypreGMRES(HypreParMatrix &_A)
Definition: hypre.cpp:1446
virtual ~HypreParMatrix()
Calls hypre&#39;s destroy function.
Definition: hypre.cpp:731
hypre_ParCSRMatrix * StealData()
Changes the ownership of the the matrix.
Definition: hypre.cpp:539
HypreParMatrix * A
The linear system matrix.
Definition: hypre.hpp:349
void SetMaxIter(int max_iter)
Definition: hypre.cpp:1470
double omega
SOR parameter (usually in (0,2))
Definition: hypre.hpp:270
Base class for solvers.
Definition: operator.hpp:102
virtual HYPRE_PtrToParSolverFcn SetupFcn() const
hypre&#39;s internal Setup function
Definition: hypre.hpp:572
virtual HYPRE_PtrToParSolverFcn SetupFcn() const
PCG Setup function.
Definition: hypre.hpp:413
Abstract operator.
Definition: operator.hpp:21
Wrapper for hypre&#39;s ParCSR matrix class.
Definition: hypre.hpp:103
virtual HYPRE_PtrToParSolverFcn SolveFcn() const
hypre&#39;s internal Solve function
Definition: hypre.hpp:540
~HypreParVector()
Calls hypre&#39;s destroy function.
Definition: hypre.cpp:160
virtual void Mult(const HypreParVector &b, HypreParVector &x) const
Solve Ax=b with hypre&#39;s PCG.
Definition: hypre.cpp:1377
HypreParVector * Z
Definition: hypre.hpp:258
void Read(MPI_Comm comm, const char *fname)
Reads the matrix from a file.
Definition: hypre.cpp:721
HypreParVector * V
Temporary vectors.
Definition: hypre.hpp:258
void SetType(HypreSmoother::Type type, int relax_times=1)
Set the relaxation type and number of sweeps.
Definition: hypre.cpp:1011
HypreParaSails(HypreParMatrix &A)
Definition: hypre.cpp:1568
virtual HYPRE_PtrToParSolverFcn SetupFcn() const
GMRES Setup function.
Definition: hypre.hpp:452
int poly_order
Order of the smoothing polynomial.
Definition: hypre.hpp:272
double lambda
Taubin&#39;s lambda-mu method parameters.
Definition: hypre.hpp:279
HypreParMatrix * ParMult(HypreParMatrix *A, HypreParMatrix *B)
Returns the matrix A * B.
Definition: hypre.cpp:771
HypreBoomerAMG(HypreParMatrix &A)
Definition: hypre.cpp:1602
virtual HYPRE_PtrToParSolverFcn SolveFcn() const
hypre&#39;s internal Solve function
Definition: hypre.hpp:605