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

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

Parent Directory Parent Directory | Revision Log Revision Log


Revision 5148 - (show annotations)
Mon Sep 15 01:25:23 2014 UTC (4 years, 11 months ago) by caltinay
File MIME type: text/plain
File size: 11467 byte(s)
Merging ripley diagonal storage + CUDA support into trunk.
Options file version has been incremented due to new options
'cuda' and 'nvccflags'.

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 "BOOMERAMG.h"
22
23 #define PRECONDITIONER_NO_ERROR 0
24 #define PRECONDITIONER_MAXITER_REACHED 1
25 #define PRECONDITIONER_INPUT_ERROR -1
26 #define PRECONDITIONER_MEMORY_ERROR -9
27 #define PRECONDITIONER_BREAKDOWN -10
28 #define PRECONDITIONER_NEGATIVE_NORM_ERROR -11
29 #define PRECONDITIONER_DIVERGENCE -12
30
31 namespace paso {
32
33 struct MergedSolver;
34 struct Preconditioner;
35 typedef boost::shared_ptr<Preconditioner> Preconditioner_ptr;
36 typedef boost::shared_ptr<const Preconditioner> const_Preconditioner_ptr;
37
38 struct Preconditioner_Smoother;
39 struct Preconditioner_AMG_Root;
40 struct Solver_ILU;
41 struct Solver_RILU;
42
43 // general preconditioner interface
44 struct Preconditioner
45 {
46 dim_t type;
47 dim_t sweeps;
48 /// Jacobi preconditioner
49 Preconditioner_Smoother* jacobi;
50 /// Gauss-Seidel preconditioner
51 Preconditioner_Smoother* gs;
52 /// AMG preconditioner
53 Preconditioner_AMG_Root* amg;
54 /// ILU preconditioner
55 Solver_ILU* ilu;
56 /// RILU preconditioner
57 Solver_RILU* rilu;
58 };
59
60 void Preconditioner_free(Preconditioner*);
61 Preconditioner* Preconditioner_alloc(SystemMatrix_ptr A, Options* options);
62 void Preconditioner_solve(Preconditioner* prec, SystemMatrix_ptr A, double*, double*);
63
64
65 // GAUSS SEIDEL & Jacobi
66 struct Preconditioner_LocalSmoother
67 {
68 bool Jacobi;
69 double* diag;
70 double* buffer;
71 index_t* pivot;
72 };
73
74 struct Preconditioner_Smoother
75 {
76 Preconditioner_LocalSmoother* localSmoother;
77 bool is_local;
78 };
79
80 void Preconditioner_Smoother_free(Preconditioner_Smoother * in);
81 void Preconditioner_LocalSmoother_free(Preconditioner_LocalSmoother * in);
82
83 Preconditioner_Smoother* Preconditioner_Smoother_alloc(
84 SystemMatrix_ptr A, bool jacobi, bool is_local, bool verbose);
85
86 Preconditioner_LocalSmoother* Preconditioner_LocalSmoother_alloc(
87 SparseMatrix_ptr A, bool jacobi, bool verbose);
88
89 void Preconditioner_Smoother_solve(SystemMatrix_ptr A,
90 Preconditioner_Smoother* gs, double* x, const double* b,
91 dim_t sweeps, bool x_is_initial);
92
93 void Preconditioner_LocalSmoother_solve(SparseMatrix_ptr A,
94 Preconditioner_LocalSmoother* gs, double* x, const double* b,
95 dim_t sweeps, bool x_is_initial);
96
97 err_t Preconditioner_Smoother_solve_byTolerance(SystemMatrix_ptr A,
98 Preconditioner_Smoother* gs, double* x, const double* b,
99 double atol, dim_t* sweeps, bool x_is_initial);
100
101 void Preconditioner_LocalSmoother_Sweep(SparseMatrix_ptr A,
102 Preconditioner_LocalSmoother* gs, double* x);
103
104 void Preconditioner_LocalSmoother_Sweep_sequential(
105 SparseMatrix_ptr A, Preconditioner_LocalSmoother* gs,
106 double* x);
107
108 void Preconditioner_LocalSmoother_Sweep_tiled(SparseMatrix_ptr A,
109 Preconditioner_LocalSmoother* gs, double* x);
110
111 void Preconditioner_LocalSmoother_Sweep_colored(SparseMatrix_ptr A,
112 Preconditioner_LocalSmoother* gs, double* x);
113
114
115 typedef enum
116 {
117 PASO_AMG_UNDECIDED=-1,
118 PASO_AMG_IN_F=0,
119 PASO_AMG_IN_C=1
120 } AMGBlockSelect;
121
122 /// Local preconditioner
123 struct Preconditioner_AMG
124 {
125 dim_t level;
126 /// coarse level matrix
127 SystemMatrix_ptr A_C;
128 /// prolongation n x n_C
129 SystemMatrix_ptr P;
130 /// restriction n_C x n
131 SystemMatrix_ptr R;
132
133 Preconditioner_Smoother* Smoother;
134 dim_t post_sweeps;
135 dim_t pre_sweeps;
136 /// used in direct solver
137 dim_t options_smoother;
138 /// used in direct solver
139 bool verbose;
140 /// applied reordering in direct solver
141 index_t reordering;
142 /// number of refinements in direct solver (typically =0)
143 dim_t refinements;
144 /// buffer for residual
145 double* r;
146 /// solution of coarse level system
147 double* x_C;
148 /// right hand side of coarse level system
149 double* b_C;
150 /// used on the coarsest level
151 MergedSolver* merged_solver;
152 Preconditioner_AMG* AMG_C;
153 };
154
155 void Preconditioner_AMG_free(Preconditioner_AMG* in);
156 Preconditioner_AMG* Preconditioner_AMG_alloc(SystemMatrix_ptr A, dim_t level, Options* options);
157 void Preconditioner_AMG_solve(SystemMatrix_ptr A, Preconditioner_AMG* amg, double* x, double* b);
158 void Preconditioner_AMG_setStrongConnections(SystemMatrix_ptr A, dim_t* degree_S, index_t* offset_S, index_t* S, double theta, double tau);
159 void Preconditioner_AMG_setStrongConnections_Block(SystemMatrix_ptr A, dim_t* degree_S, index_t* offset_S, index_t* S, double theta, double tau);
160 SystemMatrix_ptr Preconditioner_AMG_getProlongation(SystemMatrix_ptr A_p, const index_t* offset_S, const dim_t* degree_S, const index_t* S, dim_t n_C, index_t* counter_C, index_t interpolation_method);
161 void Preconditioner_AMG_setClassicProlongation(SystemMatrix_ptr P, SystemMatrix_ptr A, const index_t* offset_S, const dim_t* degree_S, const index_t* S,const index_t *counter_C);
162 void Preconditioner_AMG_setClassicProlongation_Block(SystemMatrix_ptr P, SystemMatrix_ptr A, const index_t* offset_S, const dim_t* degree_S, const index_t* S, const index_t* counter_C);
163 void Preconditioner_AMG_setDirectProlongation(SystemMatrix_ptr P, SystemMatrix_ptr A, const index_t* offset_S, const dim_t* degree_S, const index_t* S, const index_t* counter_C);
164 void Preconditioner_AMG_setDirectProlongation_Block(SystemMatrix_ptr P, SystemMatrix_ptr A, const index_t* offset_S, const dim_t* degree_S, const index_t* S, const index_t* counter_C);
165 double Preconditioner_AMG_getCoarseLevelSparsity(const Preconditioner_AMG* in);
166 dim_t Preconditioner_AMG_getNumCoarseUnknowns(const Preconditioner_AMG* in);
167 index_t Preconditioner_AMG_getMaxLevel(const Preconditioner_AMG* in);
168 void Preconditioner_AMG_transposeStrongConnections(dim_t n, const dim_t* degree_S, const index_t* offset_S, const index_t* S, dim_t nT, dim_t* degree_ST, index_t* offset_ST, index_t* ST);
169 void Preconditioner_AMG_CIJPCoarsening(dim_t n, dim_t my_n,
170 AMGBlockSelect* split_marker, const dim_t* degree_S,
171 const index_t* offset_S, const index_t* S, const dim_t* degree_ST,
172 const index_t* offset_ST, const index_t* ST,
173 const_Connector_ptr col_connector, const_Distribution_ptr col_dist);
174
175 SystemMatrix_ptr Preconditioner_AMG_getRestriction(SystemMatrix_ptr P);
176
177 SystemMatrix_ptr Preconditioner_AMG_buildInterpolationOperator(
178 SystemMatrix_ptr A, SystemMatrix_ptr P, SystemMatrix_ptr R);
179
180 SystemMatrix_ptr Preconditioner_AMG_buildInterpolationOperatorBlock(
181 SystemMatrix_ptr A, SystemMatrix_ptr P, SystemMatrix_ptr R);
182
183 SparseMatrix_ptr Preconditioner_AMG_mergeSystemMatrix(SystemMatrix_ptr A);
184
185 void Preconditioner_AMG_mergeSolve(Preconditioner_AMG* amg);
186
187 /// Local AMG preconditioner
188 struct Preconditioner_LocalAMG
189 {
190 dim_t level;
191 SparseMatrix_ptr A_C; /* coarse level matrix */
192 SparseMatrix_ptr P; /* prolongation n x n_C*/
193 SparseMatrix_ptr R; /* restriction n_C x n */
194
195 Preconditioner_LocalSmoother* Smoother;
196 dim_t post_sweeps;
197 dim_t pre_sweeps;
198 index_t reordering; /* applied reordering in direct solver */
199 dim_t refinements; /* number of refinements in direct solver (typically =0) */
200 double* r; /* buffer for residual */
201 double* x_C; /* solution of coarse level system */
202 double* b_C; /* right hand side of coarse level system */
203 struct Preconditioner_LocalAMG * AMG_C;
204 };
205
206 void Preconditioner_LocalAMG_free(Preconditioner_LocalAMG * in);
207 Preconditioner_LocalAMG* Preconditioner_LocalAMG_alloc(SparseMatrix_ptr A, dim_t level, Options* options);
208 void Preconditioner_LocalAMG_solve(SparseMatrix_ptr A, Preconditioner_LocalAMG * amg, double * x, double * b);
209
210 void 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);
211 void Preconditioner_LocalAMG_setStrongConnections_Block(SparseMatrix_ptr A, dim_t *degree, index_t *S, const double theta, const double tau);
212 void Preconditioner_LocalAMG_setStrongConnections(SparseMatrix_ptr A, dim_t *degree, index_t *S, const double theta, const double tau);
213
214 SparseMatrix_ptr Preconditioner_LocalAMG_getProlongation(
215 SparseMatrix_ptr A_p, const index_t* offset_S,
216 const dim_t* degree_S, const index_t* S, dim_t n_C,
217 const index_t* counter_C, index_t interpolation_method);
218
219 void Preconditioner_LocalAMG_setDirectProlongation_Block(SparseMatrix_ptr P_p, const_SparseMatrix_ptr A_p, const index_t *counter_C);
220
221 void Preconditioner_LocalAMG_setDirectProlongation(SparseMatrix_ptr P_p, const_SparseMatrix_ptr A_p, const index_t *counter_C);
222 void Preconditioner_LocalAMG_setClassicProlongation(SparseMatrix_ptr P_p, SparseMatrix_ptr A_p, const index_t* offset_S, const dim_t* degree_S, const index_t* S, const index_t *counter_C);
223 void Preconditioner_LocalAMG_setClassicProlongation_Block(SparseMatrix_ptr P_p, SparseMatrix_ptr A_p, const index_t* offset_S, const dim_t* degree_S, const index_t* S, const index_t *counter_C);
224 index_t Preconditioner_LocalAMG_getMaxLevel(const Preconditioner_LocalAMG * in);
225 double Preconditioner_LocalAMG_getCoarseLevelSparsity(const Preconditioner_LocalAMG * in);
226 dim_t Preconditioner_LocalAMG_getNumCoarseUnknowns(const Preconditioner_LocalAMG * in);
227 void Preconditioner_LocalAMG_enforceFFConnectivity(const dim_t n, const index_t* offset_S, const dim_t* degree_S, const index_t* S, AMGBlockSelect*split_marker);
228
229
230 struct Preconditioner_AMG_Root
231 {
232 bool is_local;
233 Preconditioner_AMG* amg;
234 Preconditioner_LocalAMG* localamg;
235 Preconditioner_BoomerAMG* boomeramg;
236 dim_t sweeps;
237 Preconditioner_Smoother* amgsubstitute;
238 };
239
240 Preconditioner_AMG_Root* Preconditioner_AMG_Root_alloc(SystemMatrix_ptr A, Options* options);
241 void Preconditioner_AMG_Root_free(Preconditioner_AMG_Root * in);
242 void Preconditioner_AMG_Root_solve(SystemMatrix_ptr A, Preconditioner_AMG_Root * amg, double * x, double * b);
243
244 /// ILU preconditioner
245 struct Solver_ILU
246 {
247 double* factors;
248 };
249
250 /// RILU preconditioner
251 struct Solver_RILU
252 {
253 dim_t n;
254 dim_t n_block;
255 dim_t n_F;
256 dim_t n_C;
257 double* inv_A_FF;
258 index_t* A_FF_pivot;
259 SparseMatrix_ptr A_FC;
260 SparseMatrix_ptr A_CF;
261 index_t* rows_in_F;
262 index_t* rows_in_C;
263 index_t* mask_F;
264 index_t* mask_C;
265 double* x_F;
266 double* b_F;
267 double* x_C;
268 double* b_C;
269 Solver_RILU* RILU_of_Schur;
270 };
271
272 void Solver_ILU_free(Solver_ILU * in);
273 Solver_ILU* Solver_getILU(SparseMatrix_ptr A, bool verbose);
274 void Solver_solveILU(SparseMatrix_ptr A, Solver_ILU* ilu, double* x, const double* b);
275
276 void Solver_RILU_free(Solver_RILU* in);
277 Solver_RILU* Solver_getRILU(SparseMatrix_ptr A, bool verbose);
278 void Solver_solveRILU(Solver_RILU* rilu, double* x, double* b);
279
280 void Solver_updateIncompleteSchurComplement(SparseMatrix_ptr A_CC,
281 SparseMatrix_ptr A_CF, double* invA_FF, index_t* A_FF_pivot,
282 SparseMatrix_ptr A_FC);
283
284 } // namespace paso
285
286 #endif // __PASO_PRECONDITIONER_H__
287

  ViewVC Help
Powered by ViewVC 1.1.26