MFEM  v3.3.2
Finite element discretization library
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Pages
strumpack.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.org.
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_STRUMPACK
13 #define MFEM_STRUMPACK
14 
15 #include "../config/config.hpp"
16 
17 #ifdef MFEM_USE_STRUMPACK
18 #ifdef MFEM_USE_MPI
19 #include "operator.hpp"
20 #include "hypre.hpp"
21 
22 #include <mpi.h>
23 
24 #include "StrumpackSparseSolverMPIDist.hpp"
25 
26 namespace mfem
27 {
28 
30 {
31 public:
32  /** Creates a general parallel matrix from a local CSR matrix on each
33  processor described by the I, J and data arrays. The local matrix should
34  be of size (local) nrows by (global) glob_ncols. The new parallel matrix
35  contains copies of all input arrays (so they can be deleted). */
36  STRUMPACKRowLocMatrix(MPI_Comm comm,
37  int num_loc_rows, int first_loc_row,
38  int glob_nrows, int glob_ncols,
39  int *I, int *J, double *data);
40 
41  /** Creates a copy of the parallel matrix hypParMat in STRUMPACK's RowLoc
42  format. All data is copied so the original matrix may be deleted. */
43  STRUMPACKRowLocMatrix(const HypreParMatrix & hypParMat);
44 
46 
47  void Mult(const Vector &x, Vector &y) const
48  {
49  mfem_error("STRUMPACKRowLocMatrix::Mult(...)\n"
50  " matrix vector products are not supported.");
51  }
52 
53  MPI_Comm GetComm() const { return comm_; }
54 
55  strumpack::CSRMatrixMPI<double,int>* getA() const { return A_; }
56 
57 private:
58  MPI_Comm comm_;
59  strumpack::CSRMatrixMPI<double,int>* A_;
60 
61 }; // mfem::STRUMPACKRowLocMatrix
62 
63 /** The MFEM STRUMPACK Direct Solver class.
64 
65  The mfem::STRUMPACKSolver class uses the STRUMPACK library to perform LU
66  factorization of a parallel sparse matrix. The solver is capable of handling
67  double precision types. See http://portal.nersc.gov/project/sparse/strumpack
68 */
70 {
71 public:
72  // Constructor with MPI_Comm parameter.
73  STRUMPACKSolver( int argc, char* argv[], MPI_Comm comm );
74 
75  // Constructor with STRUMPACK Matrix Object.
77 
78  // Default destructor.
79  ~STRUMPACKSolver( void );
80 
81  // Factor and solve the linear system y = Op^{-1} x.
82  void Mult( const Vector & x, Vector & y ) const;
83 
84  // Set the operator.
85  void SetOperator( const Operator & op );
86 
87  // Set various solver options. Refer to STRUMPACK documentation for
88  // details.
89  void SetFromCommandLine( );
90  void SetPrintFactorStatistics( bool print_stat );
91  void SetPrintSolveStatistics( bool print_stat );
92  void SetRelTol( double rtol );
93  void SetAbsTol( double atol );
94 
95  /**
96  * STRUMPACK is an (approximate) direct solver. It can be used as a direct
97  * solver or as a preconditioner. To use STRUMPACK as only a preconditioner,
98  * set the Krylov solver to DIRECT. STRUMPACK also provides iterative solvers
99  * which can use the preconditioner, and these iterative solvers can also be
100  * used without preconditioner.
101  *
102  * Supported values are:
103  * AUTO: Use iterative refinement if no HSS compression is used,
104  * otherwise use GMRes.
105  * DIRECT: No outer iterative solver, just a single application of
106  * the multifrontal solver.
107  * REFINE: Iterative refinement.
108  * PREC_GMRES: Preconditioned GMRes.
109  * The preconditioner is the (approx) multifrontal solver.
110  * GMRES: UN-preconditioned GMRes. (for testing mainly)
111  * PREC_BICGSTAB: Preconditioned BiCGStab.
112  * The preconditioner is the (approx) multifrontal solver.
113  * BICGSTAB: UN-preconditioned BiCGStab. (for testing mainly)
114  */
115  void SetKrylovSolver( strumpack::KrylovSolver method );
116 
117  /**
118  * Supported reorderings are:
119  * METIS, PARMETIS, SCOTCH, PTSCOTCH, RCM
120  */
121  void SetReorderingStrategy( strumpack::ReorderingStrategy method );
122 
123  /**
124  * MC64 performs (static) pivoting. Using a matching algorithm, it permutes
125  * the sparse input matrix in order to get nonzero elements on the
126  * diagonal. If the input matrix is already diagonally dominant, this
127  * reordering can be disabled.
128  * Possible values are:
129  * NONE: Don't do anything
130  * MAX_CARDINALITY: Maximum cardinality
131  * MAX_SMALLEST_DIAGONAL: Maximize smallest diagonal value
132  * MAX_SMALLEST_DIAGONAL_2: Same as MAX_SMALLEST_DIAGONAL, but
133  * different algorithm
134  * MAX_DIAGONAL_SUM: Maximize sum of diagonal values
135  * MAX_DIAGONAL_PRODUCT_SCALING: Maximize the product of the diagonal
136  * values and perform row & column scaling
137  */
138  void SetMC64Job( strumpack::MC64Job job );
139 
140 private:
141  void Init( int argc, char* argv[] );
142 
143 protected:
144 
145  MPI_Comm comm_;
147  int myid_;
148 
151 
153  strumpack::StrumpackSparseSolverMPIDist<double,int> * solver_;
154 
155 }; // mfem::STRUMPACKSolver class
156 
157 } // mfem namespace
158 
159 #endif // MFEM_USE_MPI
160 #endif // MFEM_USE_STRUMPACK
161 #endif // MFEM_STRUMPACK
MPI_Comm GetComm() const
Definition: strumpack.hpp:53
const STRUMPACKRowLocMatrix * APtr_
Definition: strumpack.hpp:152
void SetMC64Job(strumpack::MC64Job job)
Definition: strumpack.cpp:148
STRUMPACKSolver(int argc, char *argv[], MPI_Comm comm)
Definition: strumpack.cpp:86
strumpack::StrumpackSparseSolverMPIDist< double, int > * solver_
Definition: strumpack.hpp:153
strumpack::CSRMatrixMPI< double, int > * getA() const
Definition: strumpack.hpp:55
void SetPrintFactorStatistics(bool print_stat)
Definition: strumpack.cpp:127
STRUMPACKRowLocMatrix(MPI_Comm comm, int num_loc_rows, int first_loc_row, int glob_nrows, int glob_ncols, int *I, int *J, double *data)
Definition: strumpack.cpp:25
void SetPrintSolveStatistics(bool print_stat)
Definition: strumpack.cpp:132
void SetReorderingStrategy(strumpack::ReorderingStrategy method)
Definition: strumpack.cpp:142
void SetRelTol(double rtol)
Definition: strumpack.cpp:153
void mfem_error(const char *msg)
Definition: error.cpp:107
void Mult(const Vector &x, Vector &y) const
Operator application: y=A(x).
Definition: strumpack.cpp:164
void SetKrylovSolver(strumpack::KrylovSolver method)
Definition: strumpack.cpp:137
void SetAbsTol(double atol)
Definition: strumpack.cpp:158
void SetOperator(const Operator &op)
Set/update the solver for the given operator.
Definition: strumpack.cpp:198
Vector data type.
Definition: vector.hpp:41
void Mult(const Vector &x, Vector &y) const
Operator application: y=A(x).
Definition: strumpack.hpp:47
Base class for solvers.
Definition: operator.hpp:259
Abstract operator.
Definition: operator.hpp:21
Wrapper for hypre&#39;s ParCSR matrix class.
Definition: hypre.hpp:175