/[escript]/branches/doubleplusgood/paso/src/Preconditioner.h
ViewVC logotype

Contents of /branches/doubleplusgood/paso/src/Preconditioner.h

Parent Directory Parent Directory | Revision Log Revision Log


Revision 4257 - (show annotations)
Wed Feb 27 03:42:40 2013 UTC (6 years, 9 months ago) by jfenwick
File MIME type: text/plain
File size: 12780 byte(s)
Some simple experiments for c++ Finley

1
2 /*****************************************************************************
3 *
4 * Copyright (c) 2003-2013 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 since 2012 by School of Earth Sciences
13 *
14 *****************************************************************************/
15
16
17 #ifndef INC_PRECS
18 #define INC_PRECS
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 #define PASO_AMG_UNDECIDED -1
35 #define PASO_AMG_IN_F 0
36 #define PASO_AMG_IN_C 1
37
38 /* GAUSS SEIDEL & Jacobi */
39 typedef struct Paso_Preconditioner_LocalSmoother {
40 bool_t Jacobi;
41 double* diag;
42 double* buffer;
43 index_t* pivot;
44 } Paso_Preconditioner_LocalSmoother;
45
46 typedef struct Paso_Preconditioner_Smoother {
47 Paso_Preconditioner_LocalSmoother* localSmoother;
48 bool_t is_local;
49 } Paso_Preconditioner_Smoother;
50
51 void Paso_Preconditioner_Smoother_free(Paso_Preconditioner_Smoother * in);
52 void Paso_Preconditioner_LocalSmoother_free(Paso_Preconditioner_LocalSmoother * in);
53
54
55 Paso_Preconditioner_Smoother* Paso_Preconditioner_Smoother_alloc(Paso_SystemMatrix * A_p, const bool_t jacobi, const bool_t is_local, const bool_t verbose);
56 Paso_Preconditioner_LocalSmoother* Paso_Preconditioner_LocalSmoother_alloc(Paso_SparseMatrix * A_p,const bool_t jacobi, const bool_t verbose);
57
58 void Paso_Preconditioner_Smoother_solve(Paso_SystemMatrix* A, Paso_Preconditioner_Smoother * gs, double * x, const double * b, const dim_t sweeps, const bool_t x_is_initial);
59 void Paso_Preconditioner_LocalSmoother_solve(Paso_SparseMatrix* A, Paso_Preconditioner_LocalSmoother * gs, double * x, const double * b, const dim_t sweeps, const bool_t x_is_initial);
60 err_t Paso_Preconditioner_Smoother_solve_byTolerance(Paso_SystemMatrix* A, Paso_Preconditioner_Smoother * gs, double * x, const double * b, const double atol, dim_t *sweeps, const bool_t x_is_initial);
61
62 void Paso_Preconditioner_LocalSmoother_Sweep(Paso_SparseMatrix* A, Paso_Preconditioner_LocalSmoother * gs, double * x);
63 void Paso_Preconditioner_LocalSmoother_Sweep_sequential(Paso_SparseMatrix* A, Paso_Preconditioner_LocalSmoother * gs, double * x);
64 void Paso_Preconditioner_LocalSmoother_Sweep_tiled(Paso_SparseMatrix* A, Paso_Preconditioner_LocalSmoother * gs, double * x);
65 void Paso_Preconditioner_LocalSmoother_Sweep_colored(Paso_SparseMatrix* A, Paso_Preconditioner_LocalSmoother * gs, double * x);
66
67
68
69 /* Local preconditioner */
70 struct Paso_Preconditioner_AMG {
71 dim_t level;
72 Paso_SystemMatrix * A_C; /* coarse level matrix */
73 Paso_SystemMatrix * P; /* prolongation n x n_C*/
74 Paso_SystemMatrix * R; /* restriction n_C x n */
75
76 Paso_Preconditioner_Smoother* Smoother;
77 dim_t post_sweeps;
78 dim_t pre_sweeps;
79 dim_t options_smoother; /* used in direct solver */
80 bool_t verbose; /* used in direct solver */
81 index_t reordering; /* applied reordering in direct solver */
82 dim_t refinements; /* number of refinements in direct solver (typically =0) */
83 double* r; /* buffer for residual */
84 double* x_C; /* solution of coarse level system */
85 double* b_C; /* right hand side of coarse level system */
86 Paso_MergedSolver* merged_solver; /* used on the coarsest level */
87 struct Paso_Preconditioner_AMG * AMG_C;
88 };
89 typedef struct Paso_Preconditioner_AMG Paso_Preconditioner_AMG;
90
91 void Paso_Preconditioner_AMG_free(Paso_Preconditioner_AMG * in);
92 Paso_Preconditioner_AMG* Paso_Preconditioner_AMG_alloc(Paso_SystemMatrix * A_p,dim_t level,Paso_Options* options);
93 void Paso_Preconditioner_AMG_solve(Paso_SystemMatrix* A, Paso_Preconditioner_AMG * amg, double * x, double * b);
94 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);
95 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);
96 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);
97 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);
98 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);
99 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);
100 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);
101 double Paso_Preconditioner_AMG_getCoarseLevelSparsity(const Paso_Preconditioner_AMG * in);
102 dim_t Paso_Preconditioner_AMG_getNumCoarseUnknwons(const Paso_Preconditioner_AMG * in);
103 index_t Paso_Preconditioner_AMG_getMaxLevel(const Paso_Preconditioner_AMG * in);
104 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);
105 void Paso_Preconditioner_AMG_CIJPCoarsening(const dim_t n, const dim_t my_n, index_t*split_marker,
106 const dim_t* degree_S, const index_t* offset_S, const index_t* S,
107 const dim_t* degree_ST, const index_t* offset_ST, const index_t* ST,
108 Paso_Connector* col_connector, Paso_Distribution* col_dist);
109 Paso_SystemMatrix* Paso_Preconditioner_AMG_getRestriction(Paso_SystemMatrix* P);
110 Paso_SystemMatrix* Paso_Preconditioner_AMG_buildInterpolationOperator(Paso_SystemMatrix* A, Paso_SystemMatrix* P, Paso_SystemMatrix* R);
111 Paso_SystemMatrix* Paso_Preconditioner_AMG_buildInterpolationOperatorBlock(Paso_SystemMatrix* A, Paso_SystemMatrix* P, Paso_SystemMatrix* R);
112 Paso_SparseMatrix* Paso_Preconditioner_AMG_mergeSystemMatrix(Paso_SystemMatrix* A);
113 void Paso_Preconditioner_AMG_mergeSolve(Paso_Preconditioner_AMG* amg);
114
115 /* Local AMG preconditioner */
116 struct Paso_Preconditioner_LocalAMG {
117 dim_t level;
118 Paso_SparseMatrix * A_C; /* coarse level matrix */
119 Paso_SparseMatrix * P; /* prolongation n x n_C*/
120 Paso_SparseMatrix * R; /* restriction n_C x n */
121
122 Paso_Preconditioner_LocalSmoother* Smoother;
123 dim_t post_sweeps;
124 dim_t pre_sweeps;
125 index_t reordering; /* applied reordering in direct solver */
126 dim_t refinements; /* number of refinements in direct solver (typically =0) */
127 double* r; /* buffer for residual */
128 double* x_C; /* solution of coarse level system */
129 double* b_C; /* right hand side of coarse level system */
130 struct Paso_Preconditioner_LocalAMG * AMG_C;
131 };
132 typedef struct Paso_Preconditioner_LocalAMG Paso_Preconditioner_LocalAMG;
133
134 void Paso_Preconditioner_LocalAMG_free(Paso_Preconditioner_LocalAMG * in);
135 Paso_Preconditioner_LocalAMG* Paso_Preconditioner_LocalAMG_alloc(Paso_SparseMatrix * A_p,dim_t level,Paso_Options* options);
136 void Paso_Preconditioner_LocalAMG_solve(Paso_SparseMatrix* A, Paso_Preconditioner_LocalAMG * amg, double * x, double * b);
137
138 void Paso_Preconditioner_LocalAMG_RungeStuebenSearch(const dim_t n, const index_t* offset, const dim_t* degree, const index_t* S, index_t*split_marker, const bool_t usePanel);
139 void Paso_Preconditioner_LocalAMG_setStrongConnections_Block(Paso_SparseMatrix* A, dim_t *degree, index_t *S, const double theta, const double tau);
140 void Paso_Preconditioner_LocalAMG_setStrongConnections(Paso_SparseMatrix* A, dim_t *degree, index_t *S, const double theta, const double tau);
141 Paso_SparseMatrix* Paso_Preconditioner_LocalAMG_getProlongation(Paso_SparseMatrix* A_p, const index_t* offset_S, const dim_t* degree_S, const index_t* S, const dim_t n_C, const index_t* counter_C, const index_t interpolation_method);
142 void Paso_Preconditioner_LocalAMG_setDirectProlongation_Block(Paso_SparseMatrix* P_p, const Paso_SparseMatrix* A_p, const index_t *counter_C);
143 void Paso_Preconditioner_LocalAMG_setDirectProlongation(Paso_SparseMatrix* P_p, const Paso_SparseMatrix* A_p, const index_t *counter_C);
144 void Paso_Preconditioner_LocalAMG_setClassicProlongation(Paso_SparseMatrix* P_p, Paso_SparseMatrix* A_p, const index_t* offset_S, const dim_t* degree_S, const index_t* S, const index_t *counter_C);
145 void Paso_Preconditioner_LocalAMG_setClassicProlongation_Block(Paso_SparseMatrix* P_p, Paso_SparseMatrix* A_p, const index_t* offset_S, const dim_t* degree_S, const index_t* S, const index_t *counter_C);
146 index_t Paso_Preconditioner_LocalAMG_getMaxLevel(const Paso_Preconditioner_LocalAMG * in);
147 double Paso_Preconditioner_LocalAMG_getCoarseLevelSparsity(const Paso_Preconditioner_LocalAMG * in);
148 dim_t Paso_Preconditioner_LocalAMG_getNumCoarseUnknwons(const Paso_Preconditioner_LocalAMG * in);
149 void Paso_Preconditioner_LocalAMG_enforceFFConnectivity(const dim_t n, const index_t* offset_S, const dim_t* degree_S, const index_t* S, index_t*split_marker);
150
151
152 struct Paso_Preconditioner_BoomerAMG
153 {
154 Paso_BOOMERAMG_Handler* pt;
155 };
156 typedef struct Paso_Preconditioner_BoomerAMG Paso_Preconditioner_BoomerAMG;
157 void Paso_Preconditioner_BoomerAMG_free(Paso_Preconditioner_BoomerAMG * in);
158 Paso_Preconditioner_BoomerAMG* Paso_Preconditioner_BoomerAMG_alloc(Paso_SystemMatrix * A_p,Paso_Options* options);
159 void Paso_Preconditioner_BoomerAMG_solve(Paso_SystemMatrix* A, Paso_Preconditioner_BoomerAMG * amg, double * x, double * b);
160
161
162 struct Paso_Preconditioner_AMG_Root
163 {
164 bool_t is_local;
165 Paso_Preconditioner_AMG* amg;
166 Paso_Preconditioner_LocalAMG* localamg;
167 Paso_Preconditioner_BoomerAMG* boomeramg;
168 dim_t sweeps;
169 Paso_Preconditioner_Smoother* amgsubstitute;
170 };
171 typedef struct Paso_Preconditioner_AMG_Root Paso_Preconditioner_AMG_Root;
172
173 Paso_Preconditioner_AMG_Root* Paso_Preconditioner_AMG_Root_alloc(Paso_SystemMatrix *A, Paso_Options* options);
174 void Paso_Preconditioner_AMG_Root_free(Paso_Preconditioner_AMG_Root * in);
175 void Paso_Preconditioner_AMG_Root_solve(Paso_SystemMatrix* A, Paso_Preconditioner_AMG_Root * amg, double * x, double * b);
176
177 /*===============================================*/
178 /* ILU preconditioner */
179 struct Paso_Solver_ILU {
180 double* factors;
181 };
182 typedef struct Paso_Solver_ILU Paso_Solver_ILU;
183
184
185
186 /* RILU preconditioner */
187 struct Paso_Solver_RILU {
188 dim_t n;
189 dim_t n_block;
190 dim_t n_F;
191 dim_t n_C;
192 double* inv_A_FF;
193 index_t* A_FF_pivot;
194 Paso_SparseMatrix * A_FC;
195 Paso_SparseMatrix * A_CF;
196 index_t* rows_in_F;
197 index_t* rows_in_C;
198 index_t* mask_F;
199 index_t* mask_C;
200 double* x_F;
201 double* b_F;
202 double* x_C;
203 double* b_C;
204 struct Paso_Solver_RILU * RILU_of_Schur;
205 };
206 typedef struct Paso_Solver_RILU Paso_Solver_RILU;
207
208
209
210 /* general preconditioner interface */
211
212 typedef struct Paso_Preconditioner {
213 dim_t type;
214 dim_t sweeps;
215 /* jacobi preconditioner */
216 Paso_Preconditioner_Smoother* jacobi;
217 /* Gauss-Seidel preconditioner */
218 Paso_Preconditioner_Smoother* gs;
219 /* amg preconditioner */
220 Paso_Preconditioner_AMG_Root *amg;
221
222 /* ilu preconditioner */
223 Paso_Solver_ILU* ilu;
224 /* rilu preconditioner */
225 Paso_Solver_RILU* rilu;
226
227 } Paso_Preconditioner;
228
229 void Paso_Preconditioner_free(Paso_Preconditioner*);
230 Paso_Preconditioner* Paso_Preconditioner_alloc(Paso_SystemMatrix* A,Paso_Options* options);
231 void Paso_Preconditioner_solve(Paso_Preconditioner* prec, Paso_SystemMatrix* A,double*,double*);
232
233
234 /*******************************************/
235 void Paso_Solver_ILU_free(Paso_Solver_ILU * in);
236 Paso_Solver_ILU* Paso_Solver_getILU(Paso_SparseMatrix * A_p,bool_t verbose);
237 void Paso_Solver_solveILU(Paso_SparseMatrix * A, Paso_Solver_ILU * ilu, double * x, const double * b);
238
239 void Paso_Solver_RILU_free(Paso_Solver_RILU * in);
240 Paso_Solver_RILU* Paso_Solver_getRILU(Paso_SparseMatrix * A_p,bool_t verbose);
241 void Paso_Solver_solveRILU(Paso_Solver_RILU * rilu, double * x, double * b);
242
243 void Paso_Solver_updateIncompleteSchurComplement(Paso_SparseMatrix* A_CC, Paso_SparseMatrix *A_CF,double* invA_FF,index_t* A_FF_pivot, Paso_SparseMatrix *A_FC);
244
245
246 #endif /* #ifndef INC_PRECS */

  ViewVC Help
Powered by ViewVC 1.1.26