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

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

Parent Directory Parent Directory | Revision Log Revision Log


Revision 3510 - (show annotations)
Fri May 13 06:09:49 2011 UTC (8 years, 4 months ago) by gross
File MIME type: text/plain
File size: 10570 byte(s)
some fixes for the compilation without BOOMERAMG
1
2 /*******************************************************
3 *
4 * Copyright (c) 2003-2010 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_PRECS
16 #define INC_PRECS
17
18 #include "SystemMatrix.h"
19 #include "performance.h"
20 #include "BOOMERAMG.h"
21
22
23 #define PASO_AMG_UNDECIDED -1
24 #define PASO_AMG_IN_F 0
25 #define PASO_AMG_IN_C 1
26
27 /* GAUSS SEIDEL & Jacobi */
28 typedef struct Paso_Preconditioner_LocalSmoother {
29 bool_t Jacobi;
30 double* diag;
31 double* buffer;
32 index_t* pivot;
33 } Paso_Preconditioner_LocalSmoother;
34
35 typedef struct Paso_Preconditioner_Smoother {
36 Paso_Preconditioner_LocalSmoother* localSmoother;
37 bool_t is_local;
38 } Paso_Preconditioner_Smoother;
39
40 void Paso_Preconditioner_Smoother_free(Paso_Preconditioner_Smoother * in);
41 void Paso_Preconditioner_LocalSmoother_free(Paso_Preconditioner_LocalSmoother * in);
42
43 Paso_Preconditioner_Smoother* Paso_Preconditioner_Smoother_alloc(Paso_SystemMatrix * A_p, const bool_t jacobi, const bool_t is_local, const bool_t verbose);
44 Paso_Preconditioner_LocalSmoother* Paso_Preconditioner_LocalSmoother_alloc(Paso_SparseMatrix * A_p,const bool_t jacobi, const bool_t verbose);
45
46 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);
47 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);
48
49 void Paso_Preconditioner_LocalSmoother_Sweep(Paso_SparseMatrix* A, Paso_Preconditioner_LocalSmoother * gs, double * x);
50 void Paso_Preconditioner_LocalSmoother_Sweep_sequential(Paso_SparseMatrix* A, Paso_Preconditioner_LocalSmoother * gs, double * x);
51 void Paso_Preconditioner_LocalSmoother_Sweep_tiled(Paso_SparseMatrix* A, Paso_Preconditioner_LocalSmoother * gs, double * x);
52 void Paso_Preconditioner_LocalSmoother_Sweep_colored(Paso_SparseMatrix* A, Paso_Preconditioner_LocalSmoother * gs, double * x);
53
54
55 /* Local preconditioner */
56 struct Paso_Preconditioner_AMG {
57 dim_t level;
58 dim_t n;
59 dim_t n_F;
60 dim_t n_block;
61 Paso_SystemMatrix * A_C; /* coarse level matrix */
62 Paso_SystemMatrix * P; /* prolongation n x n_C*/
63 Paso_SystemMatrix * R; /* restriction n_C x n */
64
65 Paso_Preconditioner_Smoother* Smoother;
66 dim_t post_sweeps;
67 dim_t pre_sweeps;
68 index_t reordering; /* applied reordering in direct solver */
69 dim_t refinements; /* number of refinements in direct solver (typically =0) */
70 double* r; /* buffer for residual */
71 double* x_C; /* solution of coarse level system */
72 double* b_C; /* right hand side of coarse level system */
73 struct Paso_Preconditioner_AMG * AMG_C;
74 };
75 typedef struct Paso_Preconditioner_AMG Paso_Preconditioner_AMG;
76
77 void Paso_Preconditioner_AMG_free(Paso_Preconditioner_AMG * in);
78 Paso_Preconditioner_AMG* Paso_Preconditioner_AMG_alloc(Paso_SystemMatrix * A_p,dim_t level,Paso_Options* options);
79 void Paso_Preconditioner_AMG_solve(Paso_SystemMatrix* A, Paso_Preconditioner_AMG * amg, double * x, double * b);
80 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);
81 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);
82 double Paso_Preconditioner_AMG_getCoarseLevelSparsity(const Paso_Preconditioner_AMG * in);
83 dim_t Paso_Preconditioner_AMG_getNumCoarseUnknwons(const Paso_Preconditioner_AMG * in);
84 index_t Paso_Preconditioner_AMG_getMaxLevel(const Paso_Preconditioner_AMG * in);
85 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);
86 void Paso_Preconditioner_AMG_CIJPCoarsening(const dim_t n, const dim_t my_n, index_t*split_marker,
87 const dim_t* degree_S, const index_t* offset_S, const index_t* S,
88 const dim_t* degree_ST, const index_t* offset_ST, const index_t* ST,
89 Paso_Connector* col_connector, Paso_Distribution* col_dist);
90 /* Local AMG preconditioner */
91 struct Paso_Preconditioner_LocalAMG {
92 dim_t level;
93 dim_t n;
94 dim_t n_F;
95 dim_t n_block;
96
97 Paso_SparseMatrix * A_C; /* coarse level matrix */
98 Paso_SparseMatrix * P; /* prolongation n x n_C*/
99 Paso_SparseMatrix * R; /* restriction n_C x n */
100
101 Paso_Preconditioner_LocalSmoother* Smoother;
102 dim_t post_sweeps;
103 dim_t pre_sweeps;
104 index_t reordering; /* applied reordering in direct solver */
105 dim_t refinements; /* number of refinements in direct solver (typically =0) */
106 double* r; /* buffer for residual */
107 double* x_C; /* solution of coarse level system */
108 double* b_C; /* right hand side of coarse level system */
109 struct Paso_Preconditioner_LocalAMG * AMG_C;
110 };
111 typedef struct Paso_Preconditioner_LocalAMG Paso_Preconditioner_LocalAMG;
112
113 void Paso_Preconditioner_LocalAMG_free(Paso_Preconditioner_LocalAMG * in);
114 Paso_Preconditioner_LocalAMG* Paso_Preconditioner_LocalAMG_alloc(Paso_SparseMatrix * A_p,dim_t level,Paso_Options* options);
115 void Paso_Preconditioner_LocalAMG_solve(Paso_SparseMatrix* A, Paso_Preconditioner_LocalAMG * amg, double * x, double * b);
116
117 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);
118 void Paso_Preconditioner_LocalAMG_setStrongConnections_Block(Paso_SparseMatrix* A, dim_t *degree, index_t *S, const double theta, const double tau);
119 void Paso_Preconditioner_LocalAMG_setStrongConnections(Paso_SparseMatrix* A, dim_t *degree, index_t *S, const double theta, const double tau);
120 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);
121 void Paso_Preconditioner_LocalAMG_setDirectProlongation_Block(Paso_SparseMatrix* P_p, const Paso_SparseMatrix* A_p, const index_t *counter_C);
122 void Paso_Preconditioner_LocalAMG_setDirectProlongation(Paso_SparseMatrix* P_p, const Paso_SparseMatrix* A_p, const index_t *counter_C);
123 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);
124 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);
125 index_t Paso_Preconditioner_LocalAMG_getMaxLevel(const Paso_Preconditioner_LocalAMG * in);
126 double Paso_Preconditioner_LocalAMG_getCoarseLevelSparsity(const Paso_Preconditioner_LocalAMG * in);
127 dim_t Paso_Preconditioner_LocalAMG_getNumCoarseUnknwons(const Paso_Preconditioner_LocalAMG * in);
128 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);
129
130
131 struct Paso_Preconditioner_BoomerAMG
132 {
133 Paso_BOOMERAMG_Handler* pt;
134 };
135 typedef struct Paso_Preconditioner_BoomerAMG Paso_Preconditioner_BoomerAMG;
136 void Paso_Preconditioner_BoomerAMG_free(Paso_Preconditioner_BoomerAMG * in);
137 Paso_Preconditioner_BoomerAMG* Paso_Preconditioner_BoomerAMG_alloc(Paso_SystemMatrix * A_p,Paso_Options* options);
138 void Paso_Preconditioner_BoomerAMG_solve(Paso_SystemMatrix* A, Paso_Preconditioner_BoomerAMG * amg, double * x, double * b);
139
140
141 struct Paso_Preconditioner_AMG_Root
142 {
143 bool_t is_local;
144 Paso_Preconditioner_AMG* amg;
145 Paso_Preconditioner_LocalAMG* localamg;
146 Paso_Preconditioner_BoomerAMG* boomeramg;
147 dim_t sweeps;
148 Paso_Preconditioner_Smoother* amgsubstitute;
149 };
150 typedef struct Paso_Preconditioner_AMG_Root Paso_Preconditioner_AMG_Root;
151
152 Paso_Preconditioner_AMG_Root* Paso_Preconditioner_AMG_Root_alloc(Paso_SystemMatrix *A, Paso_Options* options);
153 void Paso_Preconditioner_AMG_Root_free(Paso_Preconditioner_AMG_Root * in);
154 void Paso_Preconditioner_AMG_Root_solve(Paso_SystemMatrix* A, Paso_Preconditioner_AMG_Root * amg, double * x, double * b);
155
156 /*===============================================*/
157 /* ILU preconditioner */
158 struct Paso_Solver_ILU {
159 double* factors;
160 };
161 typedef struct Paso_Solver_ILU Paso_Solver_ILU;
162
163
164
165 /* RILU preconditioner */
166 struct Paso_Solver_RILU {
167 dim_t n;
168 dim_t n_block;
169 dim_t n_F;
170 dim_t n_C;
171 double* inv_A_FF;
172 index_t* A_FF_pivot;
173 Paso_SparseMatrix * A_FC;
174 Paso_SparseMatrix * A_CF;
175 index_t* rows_in_F;
176 index_t* rows_in_C;
177 index_t* mask_F;
178 index_t* mask_C;
179 double* x_F;
180 double* b_F;
181 double* x_C;
182 double* b_C;
183 struct Paso_Solver_RILU * RILU_of_Schur;
184 };
185 typedef struct Paso_Solver_RILU Paso_Solver_RILU;
186
187
188
189 /* general preconditioner interface */
190
191 typedef struct Paso_Preconditioner {
192 dim_t type;
193 dim_t sweeps;
194 /* jacobi preconditioner */
195 Paso_Preconditioner_Smoother* jacobi;
196 /* Gauss-Seidel preconditioner */
197 Paso_Preconditioner_Smoother* gs;
198 /* amg preconditioner */
199 Paso_Preconditioner_AMG_Root *amg;
200
201 /* ilu preconditioner */
202 Paso_Solver_ILU* ilu;
203 /* rilu preconditioner */
204 Paso_Solver_RILU* rilu;
205
206 } Paso_Preconditioner;
207
208 void Paso_Preconditioner_free(Paso_Preconditioner*);
209 Paso_Preconditioner* Paso_Preconditioner_alloc(Paso_SystemMatrix* A,Paso_Options* options);
210 void Paso_Preconditioner_solve(Paso_Preconditioner* prec, Paso_SystemMatrix* A,double*,double*);
211
212
213 /*******************************************/
214 void Paso_Solver_ILU_free(Paso_Solver_ILU * in);
215 Paso_Solver_ILU* Paso_Solver_getILU(Paso_SparseMatrix * A_p,bool_t verbose);
216 void Paso_Solver_solveILU(Paso_SparseMatrix * A, Paso_Solver_ILU * ilu, double * x, const double * b);
217
218 void Paso_Solver_RILU_free(Paso_Solver_RILU * in);
219 Paso_Solver_RILU* Paso_Solver_getRILU(Paso_SparseMatrix * A_p,bool_t verbose);
220 void Paso_Solver_solveRILU(Paso_Solver_RILU * rilu, double * x, double * b);
221
222 void Paso_Solver_updateIncompleteSchurComplement(Paso_SparseMatrix* A_CC, Paso_SparseMatrix *A_CF,double* invA_FF,index_t* A_FF_pivot, Paso_SparseMatrix *A_FC);
223
224
225 #endif /* #ifndef INC_PRECS */

  ViewVC Help
Powered by ViewVC 1.1.26