escript  Revision_
Solver.h
Go to the documentation of this file.
1 
2 /*****************************************************************************
3 *
4 * Copyright (c) 2003-2020 by The University of Queensland
5 * http://www.uq.edu.au
6 *
7 * Primary Business: Queensland, Australia
8 * Licensed under the Apache License, version 2.0
9 * http://www.apache.org/licenses/LICENSE-2.0
10 *
11 * Development until 2012 by Earth Systems Science Computational Center (ESSCC)
12 * Development 2012-2013 by School of Earth Sciences
13 * Development from 2014-2017 by Centre for Geoscience Computing (GeoComp)
14 * Development from 2019 by School of Earth and Environmental Sciences
15 **
16 *****************************************************************************/
17 
18 
19 #ifndef __PASO_SOLVER_H__
20 #define __PASO_SOLVER_H__
21 
22 #include "Paso.h"
23 #include "Functions.h"
24 #include "performance.h"
25 #include "SystemMatrix.h"
26 
27 namespace paso {
28 
29 #define TOLERANCE_FOR_SCALARS (double)(0.)
30 
31 struct Function;
32 
33 template <typename T>
34 void solve_free(SystemMatrix<T>* A);
35 
36 SolverResult Solver(SystemMatrix_ptr<double>, double*, double*, Options*, Performance*);
37 void PASO_DLL_API Solver_free(SystemMatrix<double>*);
38 
39 SolverResult Solver(SystemMatrix_ptr<cplx_t>, cplx_t*, cplx_t*, Options*, Performance*);
40 void PASO_DLL_API Solver_free(SystemMatrix<cplx_t>*);
41 
42 SolverResult Solver_BiCGStab(SystemMatrix_ptr<double> A, double* B, double* X,
43  dim_t* iter, double* tolerance, Performance* pp);
44 
45 SolverResult Solver_PCG(SystemMatrix_ptr<double> A, double* B, double* X, dim_t* iter,
46  double* tolerance, Performance* pp);
47 
48 SolverResult Solver_TFQMR(SystemMatrix_ptr<double> A, double* B, double* X, dim_t* iter,
49  double* tolerance, Performance* pp);
50 
51 SolverResult Solver_MINRES(SystemMatrix_ptr<double> A, double* B, double* X,
52  dim_t* iter, double* tolerance, Performance* pp);
53 
54 SolverResult Solver_GMRES(SystemMatrix_ptr<double> A, double* r, double* x,
55  dim_t* num_iter, double* tolerance,
56  dim_t length_of_recursion, dim_t restart,
57  Performance* pp);
58 
59 SolverResult Solver_GMRES2(Function* F, const double* f0, const double* x0,
60  double* x, dim_t* iter, double* tolerance,
61  Performance* pp);
62 
63 SolverResult Solver_NewtonGMRES(Function* F, double* x, Options* options,
64  Performance* pp);
65 
66 } // namespace paso
67 
68 #include "Preconditioner.h"
69 #include "MKL.h"
70 #include "UMFPACK.h"
71 #include "MUMPS.h"
72 
73 namespace paso {
74 
75 struct Preconditioner_Smoother;
76 void PASO_DLL_API Preconditioner_Smoother_free(Preconditioner_Smoother * in);
77 
78 template <typename T>
80 {
81  if (!in) return;
82 
83  switch(in->solver_package) {
84  case PASO_PASO:
85  Solver_free(in);
86  break;
87 
88  case PASO_SMOOTHER:
90  break;
91 
92  case PASO_MKL:
93  MKL_free(in->mainBlock.get());
94  break;
95 
96  case PASO_UMFPACK:
97  UMFPACK_free(in->mainBlock.get());
98  break;
99 
100  case PASO_MUMPS:
101  MUMPS_free(in->mainBlock.get());
102  break;
103  }
104 }
105 
106 } // namespace paso
107 
108 #endif // __PASO_SOLVER_H__
109 
#define PASO_SMOOTHER
Definition: Options.h:75
#define PASO_MUMPS
Definition: Options.h:57
#define PASO_UMFPACK
Definition: Options.h:51
#define PASO_PASO
Definition: Options.h:56
#define PASO_MKL
Definition: Options.h:50
this class holds a (distributed) stiffness matrix
Definition: SystemMatrix.h:50
void * solver_p
pointer to data needed by a solver
Definition: SystemMatrix.h:353
index_t solver_package
package code controlling the solver pointer
Definition: SystemMatrix.h:350
SparseMatrix_ptr< T > mainBlock
main block
Definition: SystemMatrix.h:329
std::complex< real_t > cplx_t
complex data type
Definition: DataTypes.h:55
index_t dim_t
Definition: DataTypes.h:66
Definition: BiCGStab.cpp:25
SolverResult Solver_PCG(SystemMatrix_ptr< double > A, double *r, double *x, dim_t *iter, double *tolerance, Performance *pp)
Definition: PCG.cpp:62
void Solver_free(SystemMatrix< double > *A)
Definition: Solver.cpp:40
SolverResult Solver_TFQMR(SystemMatrix_ptr< double > A, double *B, double *X, dim_t *iter, double *tolerance, Performance *pp)
Definition: TFQMR.cpp:62
SolverResult Solver(SystemMatrix_ptr< double > A, double *x, double *b, Options *options, Performance *pp)
calls the iterative solver
Definition: Solver.cpp:46
void Preconditioner_Smoother_free(Preconditioner_Smoother *in)
Definition: Smoother.cpp:34
SolverResult Solver_GMRES(SystemMatrix_ptr< double > A, double *r, double *x, dim_t *iter, double *tolerance, dim_t Length_of_recursion, dim_t restart, Performance *pp)
Definition: GMRES.cpp:68
void UMFPACK_free(SparseMatrix< double > *A)
frees any UMFPACK related data from the matrix
Definition: UMFPACK.cpp:35
void solve_free(SystemMatrix< T > *A)
Definition: Solver.h:79
SolverResult
Definition: Paso.h:44
SolverResult Solver_MINRES(SystemMatrix_ptr< double > A, double *R, double *X, dim_t *iter, double *tolerance, Performance *pp)
Definition: MINRES.cpp:59
SolverResult Solver_BiCGStab(SystemMatrix_ptr< double > A, double *r, double *x, dim_t *iter, double *tolerance, Performance *pp)
Definition: BiCGStab.cpp:77
void MUMPS_free(SparseMatrix< T > *A)
frees any MUMPS related data from the matrix
Definition: MUMPS.h:118
void MKL_free(SparseMatrix< double > *A)
Definition: MKL.cpp:35
SolverResult Solver_GMRES2(Function *F, const double *f0, const double *x0, double *dx, dim_t *iter, double *tolerance, Performance *pp)
Definition: GMRES2.cpp:24
SolverResult Solver_NewtonGMRES(Function *F, double *x, Options *options, Performance *pp)
Definition: NewtonGMRES.cpp:43
#define PASO_DLL_API
Definition: paso/src/system_dep.h:29
Definition: Preconditioner.h:65