escript  Revision_Unversioneddirectory
Solver.h
Go to the documentation of this file.
1 
2 /*****************************************************************************
3 *
4 * Copyright (c) 2003-2016 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 by Centre for Geoscience Computing (GeoComp)
14 *
15 *****************************************************************************/
16 
17 
18 #ifndef __PASO_SOLVER_H__
19 #define __PASO_SOLVER_H__
20 
21 #include "SystemMatrix.h"
22 #include "performance.h"
23 #include "Functions.h"
24 
25 namespace paso {
26 
27 // error codes used in the solver
28 #define SOLVER_NO_ERROR 0
29 #define SOLVER_MAXITER_REACHED 1
30 #define SOLVER_INPUT_ERROR -1
31 #define SOLVER_MEMORY_ERROR -9
32 #define SOLVER_BREAKDOWN -10
33 #define SOLVER_NEGATIVE_NORM_ERROR -11
34 #define SOLVER_DIVERGENCE -12
35 
36 #define TOLERANCE_FOR_SCALARS (double)(0.)
37 
38 void solve(SystemMatrix_ptr A, double* out, double* in, Options* options);
39 
40 void solve_free(SystemMatrix* A);
41 
43 void Solver(SystemMatrix_ptr, double*, double*, Options*, Performance*);
44 
46 void Solver_free(SystemMatrix*);
47 
48 err_t Solver_BiCGStab(SystemMatrix_ptr A, double* B, double* X, dim_t* iter,
49  double* tolerance, Performance* pp);
50 
51 err_t Solver_PCG(SystemMatrix_ptr A, double* B, double* X, dim_t* iter,
52  double* tolerance, Performance* pp);
53 
54 err_t Solver_TFQMR(SystemMatrix_ptr A, double* B, double* X, dim_t* iter,
55  double* tolerance, Performance* pp);
56 
57 err_t Solver_MINRES(SystemMatrix_ptr A, double* B, double* X, dim_t* iter,
58  double* tolerance, Performance* pp);
59 
60 err_t Solver_GMRES(SystemMatrix_ptr A, double* r, double* x, dim_t* num_iter,
61  double* tolerance, dim_t length_of_recursion, dim_t restart,
62  Performance* pp);
63 
64 err_t Solver_GMRES2(Function* F, const double* f0, const double* x0, double* x,
65  dim_t* iter, double* tolerance, Performance* pp);
66 
67 err_t Solver_NewtonGMRES(Function* F, double* x, Options* options,
68  Performance* pp);
69 
70 } // namespace paso
71 
72 #endif // __PASO_SOLVER_H__
73 
err_t Solver_GMRES2(Function *F, const double *f0, const double *x0, double *dx, dim_t *iter, double *tolerance, Performance *pp)
Definition: GMRES2.cpp:23
#define PASO_DLL_API
Definition: Paso.h:41
void solve(SystemMatrix_ptr A, double *out, double *in, Options *options)
Definition: solve.cpp:38
err_t Solver_MINRES(SystemMatrix_ptr A, double *R, double *X, dim_t *iter, double *tolerance, Performance *pp)
Definition: MINRES.cpp:60
void Solver_free(SystemMatrix *A)
Definition: Solver.cpp:39
boost::shared_ptr< SystemMatrix > SystemMatrix_ptr
Definition: SystemMatrix.h:38
Definition: AMG.cpp:38
void Solver(SystemMatrix_ptr A, double *x, double *b, Options *options, Performance *pp)
calls the iterative solver
Definition: Solver.cpp:45
err_t Solver_PCG(SystemMatrix_ptr A, double *r, double *x, dim_t *iter, double *tolerance, Performance *pp)
Definition: PCG.cpp:63
int err_t
Definition: types.h:29
err_t Solver_NewtonGMRES(Function *F, double *x, Options *options, Performance *pp)
Definition: NewtonGMRES.cpp:40
err_t Solver_TFQMR(SystemMatrix_ptr A, double *B, double *X, dim_t *iter, double *tolerance, Performance *pp)
Definition: TFQMR.cpp:62
err_t Solver_GMRES(SystemMatrix_ptr A, double *r, double *x, dim_t *iter, double *tolerance, dim_t Length_of_recursion, dim_t restart, Performance *pp)
Definition: GMRES.cpp:67
index_t dim_t
Definition: types.h:27
void solve_free(SystemMatrix *in)
Definition: solve.cpp:116
err_t Solver_BiCGStab(SystemMatrix_ptr A, double *r, double *x, dim_t *iter, double *tolerance, Performance *pp)
Definition: BiCGStab.cpp:78