escript  Revision_Unversioneddirectory
Preconditioner.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 #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;
35 typedef boost::shared_ptr<Preconditioner> Preconditioner_ptr;
36 typedef boost::shared_ptr<const Preconditioner> const_Preconditioner_ptr;
37 
40 struct Solver_ILU;
41 struct Solver_RILU;
42 
43 // general preconditioner interface
45 {
58 };
59 
62 void Preconditioner_solve(Preconditioner* prec, SystemMatrix_ptr A, double*, double*);
63 
64 
65 // GAUSS SEIDEL & Jacobi
67 {
68  bool Jacobi;
69  double* diag;
70  double* buffer;
72 };
73 
75 {
77  bool is_local;
78 };
79 
82 
84  SystemMatrix_ptr A, bool jacobi, bool is_local, bool verbose);
85 
87  SparseMatrix_ptr A, bool jacobi, bool verbose);
88 
90  Preconditioner_Smoother* gs, double* x, const double* b,
91  dim_t sweeps, bool x_is_initial);
92 
94  Preconditioner_LocalSmoother* gs, double* x, const double* b,
95  dim_t sweeps, bool x_is_initial);
96 
98  Preconditioner_Smoother* gs, double* x, const double* b,
99  double atol, dim_t* sweeps, bool x_is_initial);
100 
102  Preconditioner_LocalSmoother* gs, double* x);
103 
106  double* x);
107 
109  Preconditioner_LocalSmoother* gs, double* x);
110 
112  Preconditioner_LocalSmoother* gs, double* x);
113 
114 
115 typedef enum
116 {
121 
124 {
125  int level;
132 
139  bool verbose;
145  double* r;
147  double* x_C;
149  double* b_C;
153 };
154 
156  Options* options);
157 
159 
161  double* x, double* b);
162 
164  dim_t* degree_S, index_t* offset_S, index_t* S,
165  double theta, double tau);
166 
168  dim_t* degree_S, index_t* offset_S, index_t* S,
169  double theta, double tau);
170 
172  const index_t* offset_S, const dim_t* degree_S,
173  const index_t* S, dim_t n_C, index_t* counter_C,
174  index_t interpolation_method);
175 
177  SystemMatrix_ptr A, const index_t* offset_S,
178  const dim_t* degree_S, const index_t* S,
179  const index_t* counter_C);
180 
182  SystemMatrix_ptr A, const index_t* offset_S,
183  const dim_t* degree_S, const index_t* S,
184  const index_t* counter_C);
185 
187  SystemMatrix_ptr A, const index_t* offset_S,
188  const dim_t* degree_S, const index_t* S,
189  const index_t* counter_C);
190 
192  SystemMatrix_ptr A, const index_t* offset_S,
193  const dim_t* degree_S, const index_t* S,
194  const index_t* counter_C);
196 
198 
200 
202  const dim_t* degree_S, const index_t* offset_S,
203  const index_t* S, dim_t nT, dim_t* degree_ST,
204  index_t* offset_ST, index_t* ST);
205 
207  AMGBlockSelect* split_marker, const dim_t* degree_S,
208  const index_t* offset_S, const index_t* S,
209  const dim_t* degree_ST, const index_t* offset_ST,
210  const index_t* ST, const_Connector_ptr col_connector,
211  const_Distribution_ptr col_dist);
212 
214 
217 
220 
222 
224 
227 {
229  SparseMatrix_ptr A_C; // coarse level matrix
230  SparseMatrix_ptr P; // prolongation n x n_C
231  SparseMatrix_ptr R; // restriction n_C x n
232 
236  index_t reordering; // applied reordering in direct solver
237  int refinements; // number of refinements in direct solver (typically=0)
238  double* r; // buffer for residual
239  double* x_C; // solution of coarse level system
240  double* b_C; // right hand side of coarse level system
242 };
243 
245  int level, Options* options);
248  Preconditioner_LocalAMG* amg, double* x, double* b);
249 
251  const dim_t* degree, const index_t* S,
252  AMGBlockSelect* split_marker, bool usePanel);
253 
255  dim_t* degree, index_t* S, double theta, double tau);
256 
258  dim_t* degree, index_t* S, double theta, double tau);
259 
261  const index_t* offset_S, const dim_t* degree_S,
262  const index_t* S, dim_t n_C, const index_t* counter_C,
263  index_t interpolation_method);
264 
266  const_SparseMatrix_ptr A, const index_t* counter_C);
267 
269  const_SparseMatrix_ptr A, const index_t* counter_C);
270 
272  SparseMatrix_ptr A, const index_t* offset_S,
273  const dim_t* degree_S, const index_t* S,
274  const index_t* counter_C);
275 
277  SparseMatrix_ptr A, const index_t* offset_S,
278  const dim_t* degree_S, const index_t* S,
279  const index_t* counter_C);
280 
285  const index_t* offset_S, const dim_t* degree_S,
286  const index_t* S, AMGBlockSelect* split_marker);
287 
288 
290 {
291  bool is_local;
295  int sweeps;
297 };
298 
300  Options* options);
303  Preconditioner_AMG_Root* amg, double* x, double* b);
304 
307 {
308  double* factors;
309 };
310 
313 {
318  double* inv_A_FF;
326  double* x_F;
327  double* b_F;
328  double* x_C;
329  double* b_C;
331 };
332 
333 void Solver_ILU_free(Solver_ILU * in);
334 Solver_ILU* Solver_getILU(SparseMatrix_ptr A, bool verbose);
335 void Solver_solveILU(SparseMatrix_ptr A, Solver_ILU* ilu, double* x, const double* b);
336 
337 void Solver_RILU_free(Solver_RILU* in);
339 void Solver_solveRILU(Solver_RILU* rilu, double* x, double* b);
340 
342  SparseMatrix_ptr A_CF, double* invA_FF, index_t* A_FF_pivot,
343  SparseMatrix_ptr A_FC);
344 
345 } // namespace paso
346 
347 #endif // __PASO_PRECONDITIONER_H__
348 
Solver_RILU * RILU_of_Schur
Definition: Preconditioner.h:330
bool is_local
Definition: Preconditioner.h:291
Solver_RILU * Solver_getRILU(SparseMatrix_ptr A, bool verbose)
Definition: RILU.cpp:71
void Preconditioner_Smoother_free(Preconditioner_Smoother *in)
Definition: Smoother.cpp:34
void Preconditioner_Smoother_solve(SystemMatrix_ptr A, Preconditioner_Smoother *gs, double *x, const double *b, dim_t sweeps, bool x_is_initial)
Definition: Smoother.cpp:111
dim_t type
Definition: Preconditioner.h:46
Preconditioner_AMG * AMG_C
Definition: Preconditioner.h:152
Preconditioner_BoomerAMG * boomeramg
Definition: Preconditioner.h:294
SystemMatrix_ptr Preconditioner_AMG_getRestriction(SystemMatrix_ptr P)
Definition: AMG_Restriction.cpp:47
double * inv_A_FF
Definition: Preconditioner.h:318
SparseMatrix_ptr P
Definition: Preconditioner.h:230
void Preconditioner_LocalSmoother_free(Preconditioner_LocalSmoother *in)
Definition: Smoother.cpp:42
AMGBlockSelect
Definition: Preconditioner.h:115
boost::shared_ptr< const SparseMatrix > const_SparseMatrix_ptr
Definition: SparseMatrix.h:37
void Preconditioner_LocalAMG_solve(SparseMatrix_ptr A, Preconditioner_LocalAMG *amg, double *x, double *b)
Definition: LocalAMG.cpp:327
double Preconditioner_LocalAMG_getCoarseLevelSparsity(const Preconditioner_LocalAMG *in)
Definition: LocalAMG.cpp:62
int post_sweeps
Definition: Preconditioner.h:234
Definition: Preconditioner.h:117
SystemMatrix_ptr A_C
coarse level matrix
Definition: Preconditioner.h:127
index_t * rows_in_C
Definition: Preconditioner.h:323
Preconditioner_Smoother * amgsubstitute
Definition: Preconditioner.h:296
double * b_C
Definition: Preconditioner.h:240
err_t Preconditioner_Smoother_solve_byTolerance(SystemMatrix_ptr A, Preconditioner_Smoother *gs, double *x, const double *b, double atol, dim_t *sweeps, bool x_is_initial)
Definition: Smoother.cpp:137
Definition: Preconditioner.h:289
SparseMatrix_ptr Preconditioner_AMG_mergeSystemMatrix(SystemMatrix_ptr A)
double * x_C
Definition: Preconditioner.h:328
void Preconditioner_LocalAMG_setStrongConnections(SparseMatrix_ptr A, dim_t *degree_S, index_t *S, const double theta, const double tau)
Definition: LocalAMG.cpp:394
void Preconditioner_AMG_solve(SystemMatrix_ptr A, Preconditioner_AMG *amg, double *x, double *b)
Definition: AMG.cpp:311
void Preconditioner_LocalSmoother_Sweep_colored(SparseMatrix_ptr A, Preconditioner_LocalSmoother *gs, double *x)
Definition: Smoother.cpp:340
double * x_C
solution of coarse level system
Definition: Preconditioner.h:147
Definition: Options.h:90
Definition: BOOMERAMG.h:44
Preconditioner_LocalAMG * localamg
Definition: Preconditioner.h:293
void Solver_ILU_free(Solver_ILU *in)
Definition: ILU.cpp:35
Preconditioner * Preconditioner_alloc(SystemMatrix_ptr A, Options *options)
Definition: Preconditioner.cpp:47
dim_t Preconditioner_AMG_getNumCoarseUnknowns(const Preconditioner_AMG *in)
Definition: AMG.cpp:74
index_t * pivot
Definition: Preconditioner.h:71
void Preconditioner_AMG_Root_solve(SystemMatrix_ptr A, Preconditioner_AMG_Root *prec, double *x, double *b)
Definition: AMG_Root.cpp:125
double * x_C
Definition: Preconditioner.h:239
Preconditioner_Smoother * jacobi
Jacobi preconditioner.
Definition: Preconditioner.h:49
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)
Definition: LocalAMG_Prolongation.cpp:369
void Preconditioner_AMG_Root_free(Preconditioner_AMG_Root *in)
Definition: AMG_Root.cpp:35
Preconditioner_AMG_Root * amg
AMG preconditioner.
Definition: Preconditioner.h:53
void Preconditioner_LocalAMG_setStrongConnections_Block(SparseMatrix_ptr A, dim_t *degree_S, index_t *S, double theta, double tau)
Definition: LocalAMG.cpp:442
double * b_F
Definition: Preconditioner.h:327
void Preconditioner_AMG_setStrongConnections_Block(SystemMatrix_ptr A, dim_t *degree_S, index_t *offset_S, index_t *S, double theta, double tau)
Definition: AMG.cpp:484
void Preconditioner_LocalSmoother_Sweep(SparseMatrix_ptr A, Preconditioner_LocalSmoother *gs, double *x)
Definition: Smoother.cpp:210
void Solver_updateIncompleteSchurComplement(SparseMatrix_ptr A_CC, SparseMatrix_ptr A_CF, double *invA_FF, index_t *A_FF_pivot, SparseMatrix_ptr A_FC)
Definition: SchurComplement.cpp:39
bool is_local
Definition: Preconditioner.h:77
int refinements
Definition: Preconditioner.h:237
Preconditioner_LocalAMG * Preconditioner_LocalAMG_alloc(SparseMatrix_ptr A_p, int level, Options *options)
Definition: LocalAMG.cpp:91
double * r
buffer for residual
Definition: Preconditioner.h:145
Definition: Preconditioner.h:119
Preconditioner_LocalSmoother * Preconditioner_LocalSmoother_alloc(SparseMatrix_ptr A, bool jacobi, bool verbose)
Definition: Smoother.cpp:69
void Preconditioner_LocalAMG_enforceFFConnectivity(dim_t n, const index_t *offset_S, const dim_t *degree_S, const index_t *S, AMGBlockSelect *split_marker)
ensures that two F nodes are connected via a C node
Definition: LocalAMG.cpp:723
int pre_sweeps
Definition: Preconditioner.h:235
int post_sweeps
Definition: Preconditioner.h:134
dim_t level
Definition: Preconditioner.h:228
void Preconditioner_LocalSmoother_solve(SparseMatrix_ptr A, Preconditioner_LocalSmoother *gs, double *x, const double *b, dim_t sweeps, bool x_is_initial)
Definition: Smoother.cpp:170
boost::shared_ptr< SparseMatrix > SparseMatrix_ptr
Definition: SparseMatrix.h:35
void Preconditioner_LocalAMG_setDirectProlongation_Block(SparseMatrix_ptr P_p, const_SparseMatrix_ptr A_p, const index_t *counter_C)
Definition: LocalAMG_Prolongation.cpp:241
boost::shared_ptr< Preconditioner > Preconditioner_ptr
Definition: Preconditioner.h:34
struct Preconditioner_LocalAMG * AMG_C
Definition: Preconditioner.h:241
index_t * mask_F
Definition: Preconditioner.h:324
void Solver_solveRILU(Solver_RILU *rilu, double *x, double *b)
Definition: RILU.cpp:306
Solver_RILU * rilu
RILU preconditioner.
Definition: Preconditioner.h:57
void Preconditioner_LocalSmoother_Sweep_sequential(SparseMatrix_ptr A, Preconditioner_LocalSmoother *gs, double *x)
inplace Gauss-Seidel sweep in sequential mode
Definition: Smoother.cpp:226
SparseMatrix_ptr R
Definition: Preconditioner.h:231
boost::shared_ptr< SystemMatrix > SystemMatrix_ptr
Definition: SystemMatrix.h:38
index_t reordering
Definition: Preconditioner.h:236
Solver_ILU * Solver_getILU(SparseMatrix_ptr A, bool verbose)
constructs the incomplete block factorization
Definition: ILU.cpp:44
double Preconditioner_AMG_getCoarseLevelSparsity(const Preconditioner_AMG *in)
Definition: AMG.cpp:62
Preconditioner_Smoother * gs
Gauss-Seidel preconditioner.
Definition: Preconditioner.h:51
Preconditioner_AMG_Root * Preconditioner_AMG_Root_alloc(SystemMatrix_ptr A, Options *options)
Definition: AMG_Root.cpp:46
Definition: AMG.cpp:38
dim_t Preconditioner_LocalAMG_getNumCoarseUnknowns(const Preconditioner_LocalAMG *in)
Definition: LocalAMG.cpp:74
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)
Definition: LocalAMG_Prolongation.cpp:465
dim_t n_block
Definition: Preconditioner.h:315
Preconditioner_Smoother * Preconditioner_Smoother_alloc(SystemMatrix_ptr A, bool jacobi, bool is_local, bool verbose)
constructs the symmetric Gauss-Seidel preconditioner
Definition: Smoother.cpp:54
void Preconditioner_LocalSmoother_Sweep_tiled(SparseMatrix_ptr A, Preconditioner_LocalSmoother *gs, double *x)
SparseMatrix_ptr A_CF
Definition: Preconditioner.h:321
Local AMG preconditioner.
Definition: Preconditioner.h:226
ILU preconditioner.
Definition: Preconditioner.h:306
double * buffer
Definition: Preconditioner.h:70
Preconditioner_LocalSmoother * localSmoother
Definition: Preconditioner.h:76
Preconditioner_Smoother * Smoother
Definition: Preconditioner.h:133
int refinements
number of refinements in direct solver (typically =0)
Definition: Preconditioner.h:143
void Preconditioner_free(Preconditioner *in)
Definition: Preconditioner.cpp:35
Definition: Preconditioner.h:118
void Preconditioner_AMG_mergeSolve(Preconditioner_AMG *amg)
Preconditioner_AMG * Preconditioner_AMG_alloc(SystemMatrix_ptr A, int level, Options *options)
Definition: AMG.cpp:91
Definition: Preconditioner.h:44
int Preconditioner_LocalAMG_getMaxLevel(const Preconditioner_LocalAMG *in)
Definition: LocalAMG.cpp:54
SystemMatrix_ptr R
restriction n_C x n
Definition: Preconditioner.h:131
int pre_sweeps
Definition: Preconditioner.h:135
void Solver_solveILU(SparseMatrix_ptr A, Solver_ILU *ilu, double *x, const double *b)
Definition: ILU.cpp:321
void Solver_RILU_free(Solver_RILU *in)
Definition: RILU.cpp:35
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)
Definition: AMG_Prolongation.cpp:699
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)
Definition: AMG_Prolongation.cpp:504
int err_t
Definition: types.h:29
Preconditioner_LocalSmoother * Smoother
Definition: Preconditioner.h:233
RILU preconditioner.
Definition: Preconditioner.h:312
index_t * A_FF_pivot
Definition: Preconditioner.h:319
Definition: Preconditioner.h:74
int level
Definition: Preconditioner.h:125
dim_t n
Definition: Preconditioner.h:314
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)
Definition: AMG_Prolongation.cpp:929
int Preconditioner_AMG_getMaxLevel(const Preconditioner_AMG *in)
Definition: AMG.cpp:53
SystemMatrix_ptr P
prolongation n x n_C
Definition: Preconditioner.h:129
Definition: MergedSolver.h:35
index_t reordering
applied reordering in direct solver
Definition: Preconditioner.h:141
int index_t
Definition: types.h:24
Local preconditioner.
Definition: Preconditioner.h:123
void Preconditioner_AMG_setStrongConnections(SystemMatrix_ptr A, dim_t *degree_S, index_t *offset_S, index_t *S, double theta, double tau)
Definition: AMG.cpp:368
bool verbose
used in direct solver
Definition: Preconditioner.h:139
MergedSolver * merged_solver
used on the coarsest level
Definition: Preconditioner.h:151
index_t * rows_in_F
Definition: Preconditioner.h:322
int sweeps
Definition: Preconditioner.h:295
void Preconditioner_LocalAMG_free(Preconditioner_LocalAMG *in)
Definition: LocalAMG.cpp:42
void Preconditioner_solve(Preconditioner *prec, SystemMatrix_ptr A, double *x, double *b)
Definition: Preconditioner.cpp:127
double * b_C
right hand side of coarse level system
Definition: Preconditioner.h:149
boost::shared_ptr< const Distribution > const_Distribution_ptr
Definition: Distribution.h:38
SparseMatrix_ptr A_FC
Definition: Preconditioner.h:320
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)
Definition: AMG_Prolongation.cpp:367
boost::shared_ptr< const Preconditioner > const_Preconditioner_ptr
Definition: Preconditioner.h:36
SparseMatrix_ptr A_C
Definition: Preconditioner.h:229
index_t dim_t
Definition: types.h:27
dim_t sweeps
Definition: Preconditioner.h:47
double * b_C
Definition: Preconditioner.h:329
double * x_F
Definition: Preconditioner.h:326
boost::shared_ptr< const Connector > const_Connector_ptr
Definition: Coupler.h:37
SystemMatrix_ptr Preconditioner_AMG_getProlongation(SystemMatrix_ptr 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)
Definition: AMG_Prolongation.cpp:56
void Preconditioner_LocalAMG_setDirectProlongation(SparseMatrix_ptr P_p, const_SparseMatrix_ptr A_p, const index_t *counter_C)
Definition: LocalAMG_Prolongation.cpp:160
SystemMatrix_ptr Preconditioner_AMG_buildInterpolationOperatorBlock(SystemMatrix_ptr A, SystemMatrix_ptr P, SystemMatrix_ptr R)
Definition: AMG_Interpolation.cpp:1989
dim_t n_C
Definition: Preconditioner.h:317
dim_t options_smoother
used in direct solver
Definition: Preconditioner.h:137
Definition: Preconditioner.h:66
SystemMatrix_ptr Preconditioner_AMG_buildInterpolationOperator(SystemMatrix_ptr A, SystemMatrix_ptr P, SystemMatrix_ptr R)
Definition: AMG_Interpolation.cpp:570
void Preconditioner_AMG_CIJPCoarsening(dim_t n, dim_t my_n, AMGBlockSelect *split_marker, const dim_t *degree_S, const index_t *offset_S, const index_t *S, const dim_t *degree_ST, const index_t *offset_ST, const index_t *ST, const_Connector_ptr col_connector, const_Distribution_ptr col_dist)
Definition: AMG.cpp:667
void Preconditioner_AMG_transposeStrongConnections(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)
Definition: AMG.cpp:640
SparseMatrix_ptr Preconditioner_LocalAMG_getProlongation(SparseMatrix_ptr A_p, const index_t *offset_S, const dim_t *degree_S, const index_t *S, dim_t n_C, const index_t *counter_C, index_t interpolation_method)
Definition: LocalAMG_Prolongation.cpp:54
#define S(_J_, _I_)
Definition: ShapeFunctions.cpp:121
index_t * mask_C
Definition: Preconditioner.h:325
dim_t n_F
Definition: Preconditioner.h:316
double * diag
Definition: Preconditioner.h:69
double * factors
Definition: Preconditioner.h:308
double * r
Definition: Preconditioner.h:238
void Preconditioner_LocalAMG_RungeStuebenSearch(dim_t n, const index_t *offset_S, const dim_t *degree_S, const index_t *S, AMGBlockSelect *split_marker, bool usePanel)
Definition: LocalAMG.cpp:508
Preconditioner_AMG * amg
Definition: Preconditioner.h:292
void Preconditioner_AMG_free(Preconditioner_AMG *in)
Definition: AMG.cpp:40
Solver_ILU * ilu
ILU preconditioner.
Definition: Preconditioner.h:55
bool Jacobi
Definition: Preconditioner.h:68