/[escript]/trunk/paso/src/Solver.h
ViewVC logotype

Contents of /trunk/paso/src/Solver.h

Parent Directory Parent Directory | Revision Log Revision Log


Revision 2828 - (show annotations)
Tue Dec 22 01:24:40 2009 UTC (9 years, 8 months ago) by artak
File MIME type: text/plain
File size: 8058 byte(s)
Smoother for AMG now can be selected from python. Now only Jacobi and Gauss-Seidel are available as smoothers.
1
2 /*******************************************************
3 *
4 * Copyright (c) 2003-2009 by University of Queensland
5 * Earth Systems Science Computational Center (ESSCC)
6 * http://www.uq.edu.au/esscc
7 *
8 * Primary Business: Queensland, Australia
9 * Licensed under the Open Software License version 3.0
10 * http://www.opensource.org/licenses/osl-3.0.php
11 *
12 *******************************************************/
13
14
15 #ifndef INC_SOLVER
16 #define INC_SOLVER
17
18 #include "SystemMatrix.h"
19 #include "performance.h"
20 #include "Functions.h"
21
22 #define PASO_TRACE
23 /* error codes used in the solver */
24 #define SOLVER_NO_ERROR 0
25 #define SOLVER_MAXITER_REACHED 1
26 #define SOLVER_INPUT_ERROR -1
27 #define SOLVER_MEMORY_ERROR -9
28 #define SOLVER_BREAKDOWN -10
29 #define SOLVER_NEGATIVE_NORM_ERROR -11
30
31 #define TOLERANCE_FOR_SCALARS (double)(0.)
32 #define PASO_ONE (double)(1.0)
33 #define PASO_ZERO (double)(0.0)
34
35 #define MAX_BLOCK_SIZE 3
36
37 /* static double ONE=1.; */
38 /* static double ZERO=0.;*/
39 /*static double TOLERANCE_FOR_SCALARS=0.;*/
40
41 /* jacobi preconditioner */
42
43 typedef struct Paso_Solver_Jacobi {
44 dim_t n_block;
45 dim_t n;
46 double* values;
47 index_t* pivot;
48 } Paso_Solver_Jacobi;
49
50
51 /* ILU preconditioner */
52 struct Paso_Solver_ILU {
53 dim_t n_block;
54 dim_t n;
55 index_t num_colors;
56 index_t* colorOf;
57 index_t* main_iptr;
58 double* factors;
59 Paso_Pattern* pattern;
60 };
61 typedef struct Paso_Solver_ILU Paso_Solver_ILU;
62
63 /* GS preconditioner */
64 struct Paso_Solver_GS {
65 dim_t n_block;
66 dim_t n;
67 index_t num_colors;
68 index_t* colorOf;
69 index_t* main_iptr;
70 double* diag;
71 Paso_SparseMatrix * factors;
72 Paso_Pattern* pattern;
73 dim_t sweeps;
74 double* x_old;
75 };
76 typedef struct Paso_Solver_GS Paso_Solver_GS;
77
78 /* RILU preconditioner */
79 struct Paso_Solver_RILU {
80 dim_t n;
81 dim_t n_block;
82 dim_t n_F;
83 dim_t n_C;
84 double* inv_A_FF;
85 index_t* A_FF_pivot;
86 Paso_SparseMatrix * A_FC;
87 Paso_SparseMatrix * A_CF;
88 index_t* rows_in_F;
89 index_t* rows_in_C;
90 index_t* mask_F;
91 index_t* mask_C;
92 double* x_F;
93 double* b_F;
94 double* x_C;
95 double* b_C;
96 struct Paso_Solver_RILU * RILU_of_Schur;
97 };
98 typedef struct Paso_Solver_RILU Paso_Solver_RILU;
99
100 struct Paso_Solver_Smoother {
101 dim_t ID;
102 Paso_Solver_Jacobi* Jacobi;
103 Paso_Solver_GS* GS;
104 };
105 typedef struct Paso_Solver_Smoother Paso_Solver_Smoother;
106
107 /* AMG preconditioner */
108 struct Paso_Solver_AMG {
109 dim_t n;
110 dim_t level;
111 bool_t coarsest_level;
112 dim_t n_block;
113 dim_t n_F;
114 dim_t n_C;
115
116 Paso_SparseMatrix * A_FF;
117 Paso_SparseMatrix * A_FC;
118 Paso_SparseMatrix * A_CF;
119 Paso_SparseMatrix * W_FC;
120
121 Paso_SparseMatrix * P;
122 Paso_SparseMatrix * R;
123
124 index_t* rows_in_F;
125 index_t* rows_in_C;
126 index_t* mask_F;
127 index_t* mask_C;
128 double* x_F;
129 double* b_F;
130 double* x_C;
131 double* b_C;
132
133 dim_t post_sweeps;
134 dim_t pre_sweeps;
135
136 Paso_SparseMatrix * A;
137 Paso_SparseMatrix * AOffset1;
138 void* solver;
139 Paso_Solver_Smoother* Smoother;
140 struct Paso_Solver_AMG * AMG_of_Coarse;
141 };
142 typedef struct Paso_Solver_AMG Paso_Solver_AMG;
143
144
145 /* AMLI preconditioner */
146 struct Paso_Solver_AMLI {
147 dim_t n;
148 dim_t level;
149 bool_t coarsest_level;
150 dim_t n_block;
151 dim_t n_F;
152 dim_t n_C;
153 double* inv_A_FF;
154 index_t* A_FF_pivot;
155 Paso_SparseMatrix * A_FC;
156 Paso_SparseMatrix * A_CF;
157 index_t* rows_in_F;
158 index_t* rows_in_C;
159 index_t* mask_F;
160 index_t* mask_C;
161 double* x_F;
162 double* b_F;
163 double* x_C;
164 double* b_C;
165
166 dim_t post_sweeps;
167 dim_t pre_sweeps;
168
169 Paso_SparseMatrix * A;
170 Paso_SparseMatrix * AOffset1;
171 void* solver;
172 Paso_Solver_Jacobi* GS;
173 struct Paso_Solver_AMLI * AMLI_of_Schur;
174 };
175 typedef struct Paso_Solver_AMLI Paso_Solver_AMLI;
176
177
178 /* AMLI preconditioner on blocks*/
179 struct Paso_Solver_AMLI_System {
180 dim_t block_size;
181 Paso_SparseMatrix *block[MAX_BLOCK_SIZE];
182 Paso_Solver_AMLI *amliblock[MAX_BLOCK_SIZE];
183 };
184 typedef struct Paso_Solver_AMLI_System Paso_Solver_AMLI_System;
185
186
187 /* AMG preconditioner on blocks*/
188 struct Paso_Solver_AMG_System {
189 dim_t block_size;
190 Paso_SparseMatrix *block[MAX_BLOCK_SIZE];
191 Paso_Solver_AMG *amgblock[MAX_BLOCK_SIZE];
192 };
193 typedef struct Paso_Solver_AMG_System Paso_Solver_AMG_System;
194
195 /* general preconditioner interface */
196
197 typedef struct Paso_Solver_Preconditioner {
198 dim_t type;
199 /* jacobi preconditioner */
200 Paso_Solver_Jacobi* jacobi;
201 /* ilu preconditioner */
202 Paso_Solver_ILU* ilu;
203 /* rilu preconditioner */
204 Paso_Solver_RILU* rilu;
205 /* Gauss-Seidel preconditioner */
206 Paso_Solver_GS* gs;
207 /* amg preconditioner */
208 Paso_Solver_AMG* amg;
209 /* amg on System */
210 Paso_Solver_AMG_System* amgSystem;
211 /* amg preconditioner */
212 Paso_Solver_AMLI* amli;
213 /* amg on System */
214 Paso_Solver_AMLI_System* amliSystem;
215
216 } Paso_Solver_Preconditioner;
217
218 void Paso_Solver(Paso_SystemMatrix*,double*,double*,Paso_Options*,Paso_Performance* pp);
219 void Paso_Solver_free(Paso_SystemMatrix*);
220 err_t Paso_Solver_BiCGStab( Paso_SystemMatrix * A, double* B, double * X, dim_t *iter, double * tolerance, Paso_Performance* pp);
221 err_t Paso_Solver_PCG( Paso_SystemMatrix * A, double* B, double * X, dim_t *iter, double * tolerance, Paso_Performance* pp);
222 err_t Paso_Solver_TFQMR( Paso_SystemMatrix * A, double* B, double * X, dim_t *iter, double * tolerance, Paso_Performance* pp);
223 err_t Paso_Solver_MINRES( Paso_SystemMatrix * A, double* B, double * X, dim_t *iter, double * tolerance, Paso_Performance* pp);
224 err_t Paso_Solver_GMRES(Paso_SystemMatrix * A, double * r, double * x, dim_t *num_iter, double * tolerance,dim_t length_of_recursion,dim_t restart, Paso_Performance* pp);
225 void Paso_Preconditioner_free(Paso_Solver_Preconditioner*);
226 void Paso_Solver_setPreconditioner(Paso_SystemMatrix* A,Paso_Options* options);
227 void Paso_Solver_solvePreconditioner(Paso_SystemMatrix* A,double*,double*);
228 void Paso_Solver_applyBlockDiagonalMatrix(dim_t n_block,dim_t n,double* D,index_t* pivot,double* x,double* b);
229
230 void Paso_Solver_ILU_free(Paso_Solver_ILU * in);
231 Paso_Solver_ILU* Paso_Solver_getILU(Paso_SparseMatrix * A_p,bool_t verbose);
232 void Paso_Solver_solveILU(Paso_Solver_ILU * ilu, double * x, double * b);
233
234 void Paso_Solver_GS_free(Paso_Solver_GS * in);
235 Paso_Solver_GS* Paso_Solver_getGS(Paso_SparseMatrix * A_p,bool_t verbose);
236 void Paso_Solver_solveGS(Paso_Solver_GS * gs, double * x, double * b);
237
238 void Paso_Solver_RILU_free(Paso_Solver_RILU * in);
239 Paso_Solver_RILU* Paso_Solver_getRILU(Paso_SparseMatrix * A_p,bool_t verbose);
240 void Paso_Solver_solveRILU(Paso_Solver_RILU * rilu, double * x, double * b);
241
242 void Paso_Solver_AMG_System_free(Paso_Solver_AMG_System * in);
243 void Paso_Solver_AMG_free(Paso_Solver_AMG * in);
244 Paso_Solver_AMG* Paso_Solver_getAMG(Paso_SparseMatrix * A_p,dim_t level,Paso_Options* options);
245 void Paso_Solver_solveAMG(Paso_Solver_AMG * amg, double * x, double * b);
246
247 void Paso_Solver_AMLI_System_free(Paso_Solver_AMLI_System * in);
248 void Paso_Solver_AMLI_free(Paso_Solver_AMLI * in);
249 Paso_Solver_AMLI* Paso_Solver_getAMLI(Paso_SparseMatrix * A_p,dim_t level,Paso_Options* options);
250 void Paso_Solver_solveAMLI(Paso_Solver_AMLI * amli, double * x, double * b);
251
252 void Paso_Solver_updateIncompleteSchurComplement(Paso_SparseMatrix* A_CC, Paso_SparseMatrix *A_CF,double* invA_FF,index_t* A_FF_pivot, Paso_SparseMatrix *A_FC);
253 Paso_Solver_Jacobi* Paso_Solver_getJacobi(Paso_SparseMatrix * A_p);
254 void Paso_Solver_solveJacobi(Paso_Solver_Jacobi * prec, double * x, double * b);
255 void Paso_Solver_Jacobi_free(Paso_Solver_Jacobi * in);
256
257 err_t Paso_Solver_GMRES2(Paso_Function * F, const double* f0, const double* x0, double * x, dim_t *iter, double* tolerance, Paso_Performance* pp);
258 err_t Paso_Solver_NewtonGMRES(Paso_Function *F, double *x, Paso_Options* options, Paso_Performance* pp);
259
260 Paso_Function * Paso_Function_LinearSystem_alloc(Paso_SystemMatrix* A, double* b, Paso_Options* options);
261 err_t Paso_Function_LinearSystem_call(Paso_Function * F,double* value, const double* arg, Paso_Performance *pp);
262 void Paso_Function_LinearSystem_free(Paso_Function * F);
263 err_t Paso_Function_LinearSystem_setInitialGuess(Paso_SystemMatrix* A, double* x, Paso_Performance *pp);
264
265 #endif /* #ifndef INC_SOLVER */

Properties

Name Value
svn:eol-style native
svn:keywords Author Date Id Revision

  ViewVC Help
Powered by ViewVC 1.1.26