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

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

Parent Directory Parent Directory | Revision Log Revision Log


Revision 4829 - (show annotations)
Thu Apr 3 04:02:53 2014 UTC (5 years, 5 months ago) by caltinay
File MIME type: text/plain
File size: 13092 byte(s)
checkpointing some SparseMatrix cleanup.

1
2 /*****************************************************************************
3 *
4 * Copyright (c) 2003-2014 by University of Queensland
5 * http://www.uq.edu.au
6 *
7 * Primary Business: Queensland, Australia
8 * Licensed under the Open Software License version 3.0
9 * http://www.opensource.org/licenses/osl-3.0.php
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 #ifndef __PASO_PRECONDITIONER_H__
18 #define __PASO_PRECONDITIONER_H__
19
20 #include "SystemMatrix.h"
21 #include "performance.h"
22 #include "BOOMERAMG.h"
23 #include "MergedSolver.h"
24
25 #define PRECONDITIONER_NO_ERROR 0
26 #define PRECONDITIONER_MAXITER_REACHED 1
27 #define PRECONDITIONER_INPUT_ERROR -1
28 #define PRECONDITIONER_MEMORY_ERROR -9
29 #define PRECONDITIONER_BREAKDOWN -10
30 #define PRECONDITIONER_NEGATIVE_NORM_ERROR -11
31 #define PRECONDITIONER_DIVERGENCE -12
32
33
34 /*
35 #define PASO_AMG_UNDECIDED -1
36 #define PASO_AMG_IN_F 0
37 #define PASO_AMG_IN_C 1
38 */
39
40
41 /* GAUSS SEIDEL & Jacobi */
42 typedef struct Paso_Preconditioner_LocalSmoother {
43 bool Jacobi;
44 double* diag;
45 double* buffer;
46 index_t* pivot;
47 } Paso_Preconditioner_LocalSmoother;
48
49 typedef struct Paso_Preconditioner_Smoother {
50 Paso_Preconditioner_LocalSmoother* localSmoother;
51 bool is_local;
52 } Paso_Preconditioner_Smoother;
53
54 void Paso_Preconditioner_Smoother_free(Paso_Preconditioner_Smoother * in);
55 void Paso_Preconditioner_LocalSmoother_free(Paso_Preconditioner_LocalSmoother * in);
56
57
58 Paso_Preconditioner_Smoother* Paso_Preconditioner_Smoother_alloc(
59 Paso_SystemMatrix* A_p, bool jacobi, bool is_local, bool verbose);
60
61 Paso_Preconditioner_LocalSmoother* Paso_Preconditioner_LocalSmoother_alloc(
62 paso::SparseMatrix_ptr A, bool jacobi, bool verbose);
63
64 void Paso_Preconditioner_Smoother_solve(Paso_SystemMatrix* A,
65 Paso_Preconditioner_Smoother* gs, double* x, const double* b,
66 dim_t sweeps, bool x_is_initial);
67
68 void Paso_Preconditioner_LocalSmoother_solve(paso::SparseMatrix_ptr A,
69 Paso_Preconditioner_LocalSmoother* gs, double* x, const double* b,
70 dim_t sweeps, bool x_is_initial);
71
72 err_t Paso_Preconditioner_Smoother_solve_byTolerance(Paso_SystemMatrix* A,
73 Paso_Preconditioner_Smoother* gs, double* x, const double* b,
74 double atol, dim_t* sweeps, bool x_is_initial);
75
76 void Paso_Preconditioner_LocalSmoother_Sweep(paso::SparseMatrix_ptr A,
77 Paso_Preconditioner_LocalSmoother* gs, double* x);
78
79 void Paso_Preconditioner_LocalSmoother_Sweep_sequential(
80 paso::SparseMatrix_ptr A, Paso_Preconditioner_LocalSmoother* gs,
81 double* x);
82
83 void Paso_Preconditioner_LocalSmoother_Sweep_tiled(paso::SparseMatrix_ptr A,
84 Paso_Preconditioner_LocalSmoother* gs, double* x);
85
86 void Paso_Preconditioner_LocalSmoother_Sweep_colored(paso::SparseMatrix_ptr A,
87 Paso_Preconditioner_LocalSmoother* gs, double* x);
88
89
90 /* Local preconditioner */
91 struct Paso_Preconditioner_AMG {
92 dim_t level;
93 Paso_SystemMatrix * A_C; /* coarse level matrix */
94 Paso_SystemMatrix * P; /* prolongation n x n_C*/
95 Paso_SystemMatrix * R; /* restriction n_C x n */
96
97 Paso_Preconditioner_Smoother* Smoother;
98 dim_t post_sweeps;
99 dim_t pre_sweeps;
100 dim_t options_smoother; /* used in direct solver */
101 bool verbose; /* used in direct solver */
102 index_t reordering; /* applied reordering in direct solver */
103 dim_t refinements; /* number of refinements in direct solver (typically =0) */
104 double* r; /* buffer for residual */
105 double* x_C; /* solution of coarse level system */
106 double* b_C; /* right hand side of coarse level system */
107 Paso_MergedSolver* merged_solver; /* used on the coarsest level */
108 struct Paso_Preconditioner_AMG * AMG_C;
109 };
110 typedef struct Paso_Preconditioner_AMG Paso_Preconditioner_AMG;
111
112 void Paso_Preconditioner_AMG_free(Paso_Preconditioner_AMG * in);
113 Paso_Preconditioner_AMG* Paso_Preconditioner_AMG_alloc(Paso_SystemMatrix * A_p,dim_t level,Paso_Options* options);
114 void Paso_Preconditioner_AMG_solve(Paso_SystemMatrix* A, Paso_Preconditioner_AMG * amg, double * x, double * b);
115 void Paso_Preconditioner_AMG_setStrongConnections(Paso_SystemMatrix* A, dim_t *degree_S, index_t* offset_S, index_t *S, const double theta, const double tau);
116 void Paso_Preconditioner_AMG_setStrongConnections_Block(Paso_SystemMatrix* A, dim_t *degree_S, index_t* offset_S, index_t *S, const double theta, const double tau);
117 Paso_SystemMatrix* Paso_Preconditioner_AMG_getProlongation(Paso_SystemMatrix* A_p, const index_t* offset_S, const dim_t* degree_S, const index_t* S, const dim_t n_C, index_t* counter_C, const index_t interpolation_method);
118 void Paso_Preconditioner_AMG_setClassicProlongation(Paso_SystemMatrix* P, Paso_SystemMatrix* A, const index_t* offset_S, const dim_t* degree_S, const index_t* S,const index_t *counter_C);
119 void Paso_Preconditioner_AMG_setClassicProlongation_Block(Paso_SystemMatrix* P, Paso_SystemMatrix* A, const index_t* offset_S, const dim_t* degree_S, const index_t* S,const index_t *counter_C);
120 void Paso_Preconditioner_AMG_setDirectProlongation(Paso_SystemMatrix* P, Paso_SystemMatrix* A, const index_t* offset_S, const dim_t* degree_S, const index_t* S,const index_t *counter_C);
121 void Paso_Preconditioner_AMG_setDirectProlongation_Block(Paso_SystemMatrix* P, Paso_SystemMatrix* A, const index_t* offset_S, const dim_t* degree_S, const index_t* S,const index_t *counter_C);
122 double Paso_Preconditioner_AMG_getCoarseLevelSparsity(const Paso_Preconditioner_AMG * in);
123 dim_t Paso_Preconditioner_AMG_getNumCoarseUnknwons(const Paso_Preconditioner_AMG * in);
124 index_t Paso_Preconditioner_AMG_getMaxLevel(const Paso_Preconditioner_AMG * in);
125 void Paso_Preconditioner_AMG_transposeStrongConnections(const dim_t n, const dim_t* degree_S, const index_t* offset_S, const index_t* S, const dim_t nT, dim_t *degree_ST, index_t* offset_ST,index_t* ST);
126 void Paso_Preconditioner_AMG_CIJPCoarsening(const dim_t n, const dim_t my_n, AMGBlockSelect*split_marker,
127 const dim_t* degree_S, const index_t* offset_S, const index_t* S,
128 const dim_t* degree_ST, const index_t* offset_ST, const index_t* ST,
129 paso::Connector_ptr col_connector, paso::const_Distribution_ptr col_dist);
130 Paso_SystemMatrix* Paso_Preconditioner_AMG_getRestriction(Paso_SystemMatrix* P);
131 Paso_SystemMatrix* Paso_Preconditioner_AMG_buildInterpolationOperator(
132 Paso_SystemMatrix* A, Paso_SystemMatrix* P, Paso_SystemMatrix* R);
133
134 Paso_SystemMatrix* Paso_Preconditioner_AMG_buildInterpolationOperatorBlock(
135 Paso_SystemMatrix* A, Paso_SystemMatrix* P, Paso_SystemMatrix* R);
136
137 paso::SparseMatrix_ptr Paso_Preconditioner_AMG_mergeSystemMatrix(Paso_SystemMatrix* A);
138
139 void Paso_Preconditioner_AMG_mergeSolve(Paso_Preconditioner_AMG* amg);
140
141 /* Local AMG preconditioner */
142 struct Paso_Preconditioner_LocalAMG {
143 dim_t level;
144 paso::SparseMatrix_ptr A_C; /* coarse level matrix */
145 paso::SparseMatrix_ptr P; /* prolongation n x n_C*/
146 paso::SparseMatrix_ptr R; /* restriction n_C x n */
147
148 Paso_Preconditioner_LocalSmoother* Smoother;
149 dim_t post_sweeps;
150 dim_t pre_sweeps;
151 index_t reordering; /* applied reordering in direct solver */
152 dim_t refinements; /* number of refinements in direct solver (typically =0) */
153 double* r; /* buffer for residual */
154 double* x_C; /* solution of coarse level system */
155 double* b_C; /* right hand side of coarse level system */
156 struct Paso_Preconditioner_LocalAMG * AMG_C;
157 };
158 typedef struct Paso_Preconditioner_LocalAMG Paso_Preconditioner_LocalAMG;
159
160 void Paso_Preconditioner_LocalAMG_free(Paso_Preconditioner_LocalAMG * in);
161 Paso_Preconditioner_LocalAMG* Paso_Preconditioner_LocalAMG_alloc(paso::SparseMatrix_ptr A_p, dim_t level, Paso_Options* options);
162 void Paso_Preconditioner_LocalAMG_solve(paso::SparseMatrix_ptr A, Paso_Preconditioner_LocalAMG * amg, double * x, double * b);
163
164 void Paso_Preconditioner_LocalAMG_RungeStuebenSearch(const dim_t n, const index_t* offset, const dim_t* degree, const index_t* S, AMGBlockSelect*split_marker, const bool usePanel);
165 void Paso_Preconditioner_LocalAMG_setStrongConnections_Block(paso::SparseMatrix_ptr A, dim_t *degree, index_t *S, const double theta, const double tau);
166 void Paso_Preconditioner_LocalAMG_setStrongConnections(paso::SparseMatrix_ptr A, dim_t *degree, index_t *S, const double theta, const double tau);
167
168 paso::SparseMatrix_ptr Paso_Preconditioner_LocalAMG_getProlongation(
169 paso::SparseMatrix_ptr A_p, const index_t* offset_S,
170 const dim_t* degree_S, const index_t* S, dim_t n_C,
171 const index_t* counter_C, index_t interpolation_method);
172
173 void Paso_Preconditioner_LocalAMG_setDirectProlongation_Block(paso::SparseMatrix_ptr P_p, paso::const_SparseMatrix_ptr A_p, const index_t *counter_C);
174
175 void Paso_Preconditioner_LocalAMG_setDirectProlongation(paso::SparseMatrix_ptr P_p, paso::const_SparseMatrix_ptr A_p, const index_t *counter_C);
176 void Paso_Preconditioner_LocalAMG_setClassicProlongation(paso::SparseMatrix_ptr P_p, paso::SparseMatrix_ptr A_p, const index_t* offset_S, const dim_t* degree_S, const index_t* S, const index_t *counter_C);
177 void Paso_Preconditioner_LocalAMG_setClassicProlongation_Block(paso::SparseMatrix_ptr P_p, paso::SparseMatrix_ptr A_p, const index_t* offset_S, const dim_t* degree_S, const index_t* S, const index_t *counter_C);
178 index_t Paso_Preconditioner_LocalAMG_getMaxLevel(const Paso_Preconditioner_LocalAMG * in);
179 double Paso_Preconditioner_LocalAMG_getCoarseLevelSparsity(const Paso_Preconditioner_LocalAMG * in);
180 dim_t Paso_Preconditioner_LocalAMG_getNumCoarseUnknwons(const Paso_Preconditioner_LocalAMG * in);
181 void Paso_Preconditioner_LocalAMG_enforceFFConnectivity(const dim_t n, const index_t* offset_S, const dim_t* degree_S, const index_t* S, AMGBlockSelect*split_marker);
182
183
184 struct Paso_Preconditioner_BoomerAMG
185 {
186 Paso_BOOMERAMG_Handler* pt;
187 };
188 typedef struct Paso_Preconditioner_BoomerAMG Paso_Preconditioner_BoomerAMG;
189 void Paso_Preconditioner_BoomerAMG_free(Paso_Preconditioner_BoomerAMG * in);
190 Paso_Preconditioner_BoomerAMG* Paso_Preconditioner_BoomerAMG_alloc(Paso_SystemMatrix * A_p,Paso_Options* options);
191 void Paso_Preconditioner_BoomerAMG_solve(Paso_SystemMatrix* A, Paso_Preconditioner_BoomerAMG * amg, double * x, double * b);
192
193
194 struct Paso_Preconditioner_AMG_Root
195 {
196 bool is_local;
197 Paso_Preconditioner_AMG* amg;
198 Paso_Preconditioner_LocalAMG* localamg;
199 Paso_Preconditioner_BoomerAMG* boomeramg;
200 dim_t sweeps;
201 Paso_Preconditioner_Smoother* amgsubstitute;
202 };
203 typedef struct Paso_Preconditioner_AMG_Root Paso_Preconditioner_AMG_Root;
204
205 Paso_Preconditioner_AMG_Root* Paso_Preconditioner_AMG_Root_alloc(Paso_SystemMatrix *A, Paso_Options* options);
206 void Paso_Preconditioner_AMG_Root_free(Paso_Preconditioner_AMG_Root * in);
207 void Paso_Preconditioner_AMG_Root_solve(Paso_SystemMatrix* A, Paso_Preconditioner_AMG_Root * amg, double * x, double * b);
208
209 /*===============================================*/
210 /* ILU preconditioner */
211 struct Paso_Solver_ILU {
212 double* factors;
213 };
214 typedef struct Paso_Solver_ILU Paso_Solver_ILU;
215
216
217
218 /* RILU preconditioner */
219 struct Paso_Solver_RILU {
220 dim_t n;
221 dim_t n_block;
222 dim_t n_F;
223 dim_t n_C;
224 double* inv_A_FF;
225 index_t* A_FF_pivot;
226 paso::SparseMatrix_ptr A_FC;
227 paso::SparseMatrix_ptr A_CF;
228 index_t* rows_in_F;
229 index_t* rows_in_C;
230 index_t* mask_F;
231 index_t* mask_C;
232 double* x_F;
233 double* b_F;
234 double* x_C;
235 double* b_C;
236 struct Paso_Solver_RILU * RILU_of_Schur;
237 };
238 typedef struct Paso_Solver_RILU Paso_Solver_RILU;
239
240
241
242 /* general preconditioner interface */
243
244 typedef struct Paso_Preconditioner {
245 dim_t type;
246 dim_t sweeps;
247 /* jacobi preconditioner */
248 Paso_Preconditioner_Smoother* jacobi;
249 /* Gauss-Seidel preconditioner */
250 Paso_Preconditioner_Smoother* gs;
251 /* amg preconditioner */
252 Paso_Preconditioner_AMG_Root *amg;
253
254 /* ilu preconditioner */
255 Paso_Solver_ILU* ilu;
256 /* rilu preconditioner */
257 Paso_Solver_RILU* rilu;
258
259 } Paso_Preconditioner;
260
261 void Paso_Preconditioner_free(Paso_Preconditioner*);
262 Paso_Preconditioner* Paso_Preconditioner_alloc(Paso_SystemMatrix* A,Paso_Options* options);
263 void Paso_Preconditioner_solve(Paso_Preconditioner* prec, Paso_SystemMatrix* A,double*,double*);
264
265
266 /*******************************************/
267 void Paso_Solver_ILU_free(Paso_Solver_ILU * in);
268 Paso_Solver_ILU* Paso_Solver_getILU(paso::SparseMatrix_ptr A, bool verbose);
269 void Paso_Solver_solveILU(paso::SparseMatrix_ptr A, Paso_Solver_ILU* ilu, double* x, const double* b);
270
271 void Paso_Solver_RILU_free(Paso_Solver_RILU* in);
272 Paso_Solver_RILU* Paso_Solver_getRILU(paso::SparseMatrix_ptr A, bool verbose);
273 void Paso_Solver_solveRILU(Paso_Solver_RILU* rilu, double* x, double* b);
274
275 void Paso_Solver_updateIncompleteSchurComplement(paso::SparseMatrix_ptr A_CC,
276 paso::SparseMatrix_ptr A_CF, double* invA_FF, index_t* A_FF_pivot,
277 paso::SparseMatrix_ptr A_FC);
278
279
280 #endif // __PASO_PRECONDITIONER_H__
281

  ViewVC Help
Powered by ViewVC 1.1.26