INMOST
Mathematical Modelling Toolkit
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Pages
inmost_solver.h
Go to the documentation of this file.
1 
2 #ifndef INMOST_SOLVER_INCLUDED
3 #define INMOST_SOLVER_INCLUDED
4 
5 
6 #include "inmost_common.h"
7 #include "inmost_sparse.h"
8 //#include "solver_prototypes.hpp"
9 
10 
11 #define DEFAULT_ADDITIVE_SCHWARTZ_OVERLAP 1
12 #define DEFAULT_ABSOLUTE_TOLERANCE 1.0e-5
13 #define DEFAULT_RELATIVE_TOLERANCE 1.0e-12
14 #define DEFAULT_DIVERGENCE_TOLERANCE 1.0e+100
15 #define DEFAULT_MAXIMUM_ITERATIONS 2500
16 #define DEFAULT_SOLVER_GMRES_SUBSTEPS 2
17 #define DEFAULT_PRECONDITIONER_DROP_TOLERANCE 0.005
18 #define DEFAULT_PRECONDITIONER_REUSE_TOLERANCE 0.00005
19 #define DEFAULT_PRECONDITIONER_FILL_LEVEL 3
20 #define DEFAULT_PRECONDITIONER_DDPQ_TOLERANCE 0.75
21 #define DEFAULT_PRECONDITIONER_REORDER_NONZEROS 1
22 #define DEFAULT_PRECONDITIONER_RESCALE_ITERS 6
23 #define DEFAULT_PRECONDITIONER_CONDITION_ESTIMATION 1
24 #define DEFAULT_PRECONDITIONER_ADAPT_DDPQ_TOLERANCE 1
25 
26 #if defined(USE_SOLVER)
27 namespace INMOST
28 {
38  class Solver
39  {
40  private:
41  static INMOST_MPI_Type RowEntryType; //prepared in Initialize
42  public:
44  enum Type
45  {
55  ANI,
58  };
59 
60  static std::string TypeName(Type t);
61 
62 
63  //solver.cpp::::::::::::::::::::::::::::::::::::::::::::::::::::
64  public:
65 
67  class OrderInfo
68  {
69  private:
70  typedef std::vector<INMOST_DATA_ENUM_TYPE> storage_type;
71  storage_type global_to_proc; //stores ends of all non-overlapping intervals of elements, owned by this processor
72  storage_type global_overlap; //stores pairs: [begin,end) of overlapping intervals of rows
73  std::vector<INMOST_DATA_ENUM_TYPE> vector_exchange_recv, vector_exchange_send;
74  std::vector<INMOST_DATA_REAL_TYPE> send_storage, recv_storage;
75  std::vector<INMOST_MPI_Request> send_requests, recv_requests;
76  std::vector<INMOST_DATA_ENUM_TYPE> extended_indexes;
77 
78  //remote indexes
79  INMOST_DATA_ENUM_TYPE local_vector_begin, local_vector_end;
80  INMOST_DATA_ENUM_TYPE initial_matrix_begin, initial_matrix_end; //local interval of matrix
81  INMOST_DATA_ENUM_TYPE local_matrix_begin, local_matrix_end; //local interval of matrix
82 
83  bool have_matrix;
84  INMOST_MPI_Comm comm;
85  int rank,size;
86  public:
87  void Clear();
89  bool & HaveMatrix() { return have_matrix; }
90  OrderInfo();
91  OrderInfo(const OrderInfo & other);
92  OrderInfo & operator =(OrderInfo const & other);
93  ~OrderInfo();
101  void RestoreMatrix(Sparse::Matrix & m);
103  void PrepareVector(Sparse::Vector & v) const;
105  void RestoreVector(Sparse::Vector & v) const;
107  INMOST_DATA_ENUM_TYPE GetProcessor(INMOST_DATA_ENUM_TYPE gind) const; //retrieve processor by binary search in global_to_proc
112  void GetVectorRegion(INMOST_DATA_ENUM_TYPE & mbeg, INMOST_DATA_ENUM_TYPE & mend) const {mbeg = local_vector_begin; mend = local_vector_end;}
114  INMOST_DATA_ENUM_TYPE GetRank() const {return rank;}
116  INMOST_DATA_ENUM_TYPE GetSize() const {return size;}
118  void Update (Sparse::Vector & x); // update parallel vector
120  void Accumulate(Sparse::Vector & x); // sum shared values in parallel vector
122  void Integrate(INMOST_DATA_REAL_TYPE * inout, INMOST_DATA_ENUM_TYPE num) const;
124  INMOST_MPI_Comm GetComm() const {return comm;}
125  // Access to arrays below allows to organize manual exchange
126  INMOST_MPI_Request * GetSendRequests() {assert(!send_requests.empty()); return &send_requests[0];}
127  INMOST_MPI_Request * GetRecvRequests() {assert(!recv_requests.empty()); return &recv_requests[0];}
128  INMOST_DATA_ENUM_TYPE GetSendRequestsSize() {return static_cast<INMOST_DATA_ENUM_TYPE>(send_requests.size());}
129  INMOST_DATA_ENUM_TYPE GetRecvRequestsSize() {return static_cast<INMOST_DATA_ENUM_TYPE>(recv_requests.size());}
130  INMOST_DATA_ENUM_TYPE * GetSendExchangeArray() {assert(!vector_exchange_send.empty()); return &vector_exchange_send[0];}
131  INMOST_DATA_ENUM_TYPE GetSendExchangeSize() {return static_cast<INMOST_DATA_ENUM_TYPE>(send_storage.size());}
132  INMOST_DATA_ENUM_TYPE * GetRecvExchangeArray() {assert(!vector_exchange_recv.empty()); return &vector_exchange_recv[0];}
133  INMOST_DATA_ENUM_TYPE GetRecvExchangeSize() {return static_cast<INMOST_DATA_ENUM_TYPE>(recv_storage.size());}
134  //for debug
135  //~ void BeginSequentialCode() {for(int i = 0; i < rank; i++) MPI_Barrier(comm);}
136  //~ void EndSequentialCode() {for(int i = rank; i < size; i++) MPI_Barrier(comm);}
137 
138  // Get the scalar product of the specified interval of the distributed vector.
139  // Conflicts with OpenMP, should not be used in future
140  //void ScalarProd(Vector const & left, Vector const & right, INMOST_DATA_ENUM_TYPE index_begin, INMOST_DATA_ENUM_TYPE index_end, INMOST_DATA_REAL_TYPE & sum) const;
141  };
142 
143 
144 
145  private:
146  static bool is_initialized, is_finalized;
147  INMOST_MPI_Comm comm;
148  std::string name;
149  INMOST_DATA_ENUM_TYPE local_size, global_size;
150  INMOST_DATA_ENUM_TYPE last_it;
151  INMOST_DATA_REAL_TYPE last_resid;
152  OrderInfo info;
153 
154  INMOST_DATA_ENUM_TYPE additive_schwartz_overlap;
155 
156  INMOST_DATA_ENUM_TYPE maximum_iterations;
157  INMOST_DATA_REAL_TYPE absolute_tolerance;
158  INMOST_DATA_REAL_TYPE relative_tolerance;
159  INMOST_DATA_REAL_TYPE divergence_tolerance;
160 
161  INMOST_DATA_REAL_TYPE preconditioner_drop_tolerance;
162  INMOST_DATA_REAL_TYPE preconditioner_reuse_tolerance;
163  INMOST_DATA_REAL_TYPE preconditioner_ddpq_tolerance;
164  INMOST_DATA_ENUM_TYPE preconditioner_reorder_nonzero;
165  INMOST_DATA_REAL_TYPE preconditioner_fill_level;
166  INMOST_DATA_ENUM_TYPE preconditioner_rescale_iterations;
167  INMOST_DATA_ENUM_TYPE preconditioner_condition_estimation;
168  INMOST_DATA_ENUM_TYPE preconditioner_adapt_ddpq_tolerance;
169 
170  INMOST_DATA_ENUM_TYPE solver_gmres_substeps;
171 
172  std::string return_reason;
173 
174  void * solver_data;
175  void * matrix_data;
176  void * precond_data;
177 
178  void * rhs_data;
179  void * solution_data;
180 
181  Type _pack;
182  Solver(const Solver & other);// prohibit copy
183  Solver & operator =(Solver const & other); //prohibit assignment
184  public:
218  void SetParameterEnum(std::string name, INMOST_DATA_ENUM_TYPE value);
253  void SetParameterReal(std::string name, INMOST_DATA_REAL_TYPE value);
255  std::string GetName() {return name;}
256 
258  Type GetPackage() const {return _pack;}
259 
272  void SetMatrix(Sparse::Matrix & A, bool OldPreconditioner = false);
283  bool Solve(Sparse::Vector & RHS, Sparse::Vector & SOL);
292  std::string GetReason();
316  Solver(Type pack, std::string _name = "", INMOST_MPI_Comm comm = INMOST_MPI_COMM_WORLD);
317  ~Solver();
341  static void Initialize(int * argc, char *** argv, const char * database = "");
347  static void Finalize();
348  static bool isInitialized() {return is_initialized;}
349  static bool isFinalized() {return is_finalized;}
351  void Clear();
352  };
353 }
354 
355 #endif // USE_SOLVER
356 
357 #endif // INMOST_SOLVER_INCLUDED
Base class for low level operations with objects of Solver class.
Definition: inmost_solver.h:67
INMOST_DATA_ENUM_TYPE GetRank() const
Get the rank of the current communicator, i.e. the current process index.
void SetParameterEnum(std::string name, INMOST_DATA_ENUM_TYPE value)
inner Solver based on BiCGStab(L) solver with second order Crout-ILU with inversed-based condition es...
Definition: inmost_solver.h:48
INMOST_MPI_Comm GetComm() const
Get the communicator which the solver is associated with.
static bool isFinalized()
INMOST_DATA_REAL_TYPE GetConditionNumberL()
void Integrate(INMOST_DATA_REAL_TYPE *inout, INMOST_DATA_ENUM_TYPE num) const
Get the sum of num elements of real array on all processes.
INMOST_DATA_ENUM_TYPE GetSize() const
Get the size of the current communicator, i.e. the total number of processes used.
void RestoreMatrix(Sparse::Matrix &m)
Restore initial nonparallel state of the Matrix with no overlap.
external Solver PETSc,
Definition: inmost_solver.h:54
INMOST_DATA_ENUM_TYPE * GetRecvExchangeArray()
inner Solver based on BiCGStab(L) solver with second order Crout-ILU with inversed-based condition es...
Definition: inmost_solver.h:47
Type GetPackage() const
Get the package Type.
external Solver AztecOO with ML preconditioner.
Definition: inmost_solver.h:52
#define INMOST_MPI_Comm
#define INMOST_MPI_COMM_WORLD
bool & HaveMatrix()
Return true if Matrix data have already been specified.
Definition: inmost_solver.h:89
#define INMOST_DATA_REAL_TYPE
void Clear()
Clear all internal data of the current solver including matrix, preconditioner etc.
void PrepareVector(Sparse::Vector &v) const
Prepare parallel state of the Vector.
INMOST_MPI_Request * GetSendRequests()
INMOST_DATA_ENUM_TYPE Iterations()
static void Initialize(int *argc, char ***argv, const char *database="")
void PrepareMatrix(Sparse::Matrix &m, INMOST_DATA_ENUM_TYPE overlap)
INMOST_DATA_ENUM_TYPE GetRecvRequestsSize()
void RestoreVector(Sparse::Vector &v) const
Restore initial nonparallel state of the Vector.
#define INMOST_MPI_Type
std::string GetName()
Get the used defined name of the Solver.
inner K3BIILU2 Solver (BIILU2 parallel version).
Definition: inmost_solver.h:57
void GetLocalRegion(INMOST_DATA_ENUM_TYPE proc, INMOST_DATA_ENUM_TYPE &mbeg, INMOST_DATA_ENUM_TYPE &mend) const
Get the local index region for the specified process.
void Accumulate(Sparse::Vector &x)
Sum shared values in parallel vector.
void Update(Sparse::Vector &x)
Update the shared data in parallel vector.
Type
Type of the Solver can be currently used in this version of INMOST.
Definition: inmost_solver.h:44
OrderInfo & operator=(OrderInfo const &other)
bool Solve(Sparse::Vector &RHS, Sparse::Vector &SOL)
static bool isInitialized()
INMOST_MPI_Request * GetRecvRequests()
external FCBIILU2 Solver (BIILU2 parallel F2C version).
Definition: inmost_solver.h:56
external Solver AztecOO from Trilinos package.
Definition: inmost_solver.h:50
INMOST_DATA_ENUM_TYPE * GetSendExchangeArray()
INMOST_DATA_ENUM_TYPE GetSendExchangeSize()
INMOST_DATA_REAL_TYPE Residual()
external Solver AztecOO with Ifpack preconditioner.
Definition: inmost_solver.h:53
void SetMatrix(Sparse::Matrix &A, bool OldPreconditioner=false)
#define INMOST_DATA_ENUM_TYPE
INMOST_DATA_REAL_TYPE GetConditionNumberU()
static void Finalize()
void SetParameterReal(std::string name, INMOST_DATA_REAL_TYPE value)
void GetVectorRegion(INMOST_DATA_ENUM_TYPE &mbeg, INMOST_DATA_ENUM_TYPE &mend) const
Get the local index region for the current process.
std::string GetReason()
inner Solver based on BiCGStab(L) solver with second order ILU and maximum product transversal reorde...
Definition: inmost_solver.h:49
external Solver Belos from Trilinos package, currently without preconditioner.
Definition: inmost_solver.h:51
#define INMOST_MPI_Request
static std::string TypeName(Type t)
INMOST_DATA_ENUM_TYPE GetRecvExchangeSize()
INMOST_DATA_REAL_TYPE Condest(INMOST_DATA_REAL_TYPE tol, INMOST_DATA_ENUM_TYPE maxits=100)
external Solver from ANI3D based on ILU2 (sequential Fortran version),
Definition: inmost_solver.h:55
INMOST_DATA_ENUM_TYPE GetSendRequestsSize()
INMOST_DATA_ENUM_TYPE GetProcessor(INMOST_DATA_ENUM_TYPE gind) const
Retrieve the processor number by binary search for the specified global index.
inner Solver based on BiCGStab(L) solver with second order ILU factorization as preconditioner.
Definition: inmost_solver.h:46
void GetOverlapRegion(INMOST_DATA_ENUM_TYPE proc, INMOST_DATA_ENUM_TYPE &mbeg, INMOST_DATA_ENUM_TYPE &mend) const