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

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

Parent Directory Parent Directory | Revision Log Revision Log


Revision 3158 - (show annotations)
Mon Sep 6 06:09:11 2010 UTC (9 years ago) by gross
File MIME type: text/plain
File size: 6622 byte(s)
GS and Jacobi are now used through the same interface.
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
21 #define MAX_BLOCK_SIZE 3
22
23 /* GAUSS SEIDEL & Jacobi */
24 typedef struct Paso_Preconditioner_LocalSmoother {
25 bool_t Jacobi;
26 double* diag;
27 double* buffer;
28 index_t* pivot;
29 } Paso_Preconditioner_LocalSmoother;
30
31 typedef struct Paso_Preconditioner_Smoother {
32 Paso_Preconditioner_LocalSmoother* localSmoother;
33 bool_t is_local;
34 } Paso_Preconditioner_Smoother;
35
36 void Paso_Preconditioner_Smoother_free(Paso_Preconditioner_Smoother * in);
37 void Paso_Preconditioner_LocalSmoother_free(Paso_Preconditioner_LocalSmoother * in);
38
39 Paso_Preconditioner_Smoother* Paso_Preconditioner_Smoother_alloc(Paso_SystemMatrix * A_p, const bool_t jacobi, const bool_t is_local, const bool_t verbose);
40 Paso_Preconditioner_LocalSmoother* Paso_Preconditioner_LocalSmoother_alloc(Paso_SparseMatrix * A_p,const bool_t jacobi, const bool_t verbose);
41
42 void Paso_Preconditioner_Smoother_solve(Paso_SystemMatrix* A, Paso_Preconditioner_Smoother * gs, double * x, const double * b, const dim_t sweeps);
43 void Paso_Preconditioner_LocalSmoother_solve(Paso_SparseMatrix* A, Paso_Preconditioner_LocalSmoother * gs, double * x, const double * b, const dim_t sweeps);
44
45 void Paso_Preconditioner_LocalSmoother_Sweep(Paso_SparseMatrix* A, Paso_Preconditioner_LocalSmoother * gs, double * x);
46 void Paso_Preconditioner_LocalSmoother_Sweep_sequential(Paso_SparseMatrix* A, Paso_Preconditioner_LocalSmoother * gs, double * x);
47 void Paso_Preconditioner_LocalSmoother_Sweep_tiled(Paso_SparseMatrix* A, Paso_Preconditioner_LocalSmoother * gs, double * x);
48 void Paso_Preconditioner_LocalSmoother_Sweep_colored(Paso_SparseMatrix* A, Paso_Preconditioner_LocalSmoother * gs, double * x);
49
50 /*===============================================*/
51 /* ILU preconditioner */
52 struct Paso_Solver_ILU {
53 double* factors;
54 };
55 typedef struct Paso_Solver_ILU Paso_Solver_ILU;
56
57
58
59
60 /* RILU preconditioner */
61 struct Paso_Solver_RILU {
62 dim_t n;
63 dim_t n_block;
64 dim_t n_F;
65 dim_t n_C;
66 double* inv_A_FF;
67 index_t* A_FF_pivot;
68 Paso_SparseMatrix * A_FC;
69 Paso_SparseMatrix * A_CF;
70 index_t* rows_in_F;
71 index_t* rows_in_C;
72 index_t* mask_F;
73 index_t* mask_C;
74 double* x_F;
75 double* b_F;
76 double* x_C;
77 double* b_C;
78 struct Paso_Solver_RILU * RILU_of_Schur;
79 };
80 typedef struct Paso_Solver_RILU Paso_Solver_RILU;
81
82 struct Paso_Solver_Smoother {
83 dim_t ID;
84 Paso_Preconditioner_LocalSmoother* Jacobi;
85 Paso_Preconditioner_LocalSmoother* GS;
86 };
87 typedef struct Paso_Solver_Smoother Paso_Solver_Smoother;
88
89 /* AMG preconditioner */
90 struct Paso_Solver_AMG {
91 dim_t n;
92 dim_t level;
93 bool_t coarsest_level;
94 dim_t n_block;
95 dim_t n_F;
96 dim_t n_C;
97
98 Paso_SparseMatrix * A_FF;
99 Paso_SparseMatrix * A_FC;
100 Paso_SparseMatrix * A_CF;
101 Paso_SparseMatrix * W_FC;
102
103 Paso_SparseMatrix * P;
104 Paso_SparseMatrix * R;
105
106 index_t* rows_in_F;
107 index_t* rows_in_C;
108 index_t* mask_F;
109 index_t* mask_C;
110 double* x_F;
111 double* b_F;
112 double* x_C;
113 double* b_C;
114
115 dim_t post_sweeps;
116 dim_t pre_sweeps;
117
118 Paso_SparseMatrix * A;
119 Paso_SparseMatrix * AOffset1;
120 Paso_SparseMatrix * AUnrolled;
121 void* solver;
122 Paso_Solver_Smoother* Smoother;
123 struct Paso_Solver_AMG * AMG_of_Coarse;
124 };
125 typedef struct Paso_Solver_AMG Paso_Solver_AMG;
126
127
128 /* AMLI preconditioner */
129 struct Paso_Solver_AMLI {
130 dim_t n;
131 dim_t level;
132 bool_t coarsest_level;
133 dim_t n_block;
134 dim_t n_F;
135 dim_t n_C;
136 double* inv_A_FF;
137 index_t* A_FF_pivot;
138 Paso_SparseMatrix * A_FC;
139 Paso_SparseMatrix * A_CF;
140 index_t* rows_in_F;
141 index_t* rows_in_C;
142 index_t* mask_F;
143 index_t* mask_C;
144 double* x_F;
145 double* b_F;
146 double* x_C;
147 double* b_C;
148
149 dim_t post_sweeps;
150 dim_t pre_sweeps;
151
152 Paso_SparseMatrix * A;
153 Paso_SparseMatrix * AOffset1;
154 void* solver;
155 Paso_Preconditioner_LocalSmoother* GS;
156 struct Paso_Solver_AMLI * AMLI_of_Schur;
157 };
158 typedef struct Paso_Solver_AMLI Paso_Solver_AMLI;
159
160 /* AMLI preconditioner on blocks*/
161 struct Paso_Solver_AMLI_System {
162 dim_t block_size;
163 Paso_SparseMatrix *block[MAX_BLOCK_SIZE];
164 Paso_Solver_AMLI *amliblock[MAX_BLOCK_SIZE];
165 };
166 typedef struct Paso_Solver_AMLI_System Paso_Solver_AMLI_System;
167
168
169 /* general preconditioner interface */
170
171 typedef struct Paso_Solver_Preconditioner {
172 dim_t type;
173 dim_t sweeps;
174 /* jacobi preconditioner */
175 Paso_Preconditioner_Smoother* jacobi;
176 /* Gauss-Seidel preconditioner */
177 Paso_Preconditioner_Smoother* gs;
178
179 /* ilu preconditioner */
180 Paso_Solver_ILU* ilu;
181 /* rilu preconditioner */
182 Paso_Solver_RILU* rilu;
183 /* amg preconditioner */
184 Paso_Solver_AMG* amg;
185 /* amg preconditioner */
186 Paso_Solver_AMLI* amli;
187 /* amg on System */
188 Paso_Solver_AMLI_System* amliSystem;
189
190 } Paso_Solver_Preconditioner;
191
192 void Paso_Preconditioner_free(Paso_Solver_Preconditioner*);
193 Paso_Solver_Preconditioner* Paso_Preconditioner_alloc(Paso_SystemMatrix* A,Paso_Options* options);
194 void Paso_Preconditioner_solve(Paso_Solver_Preconditioner* prec, Paso_SystemMatrix* A,double*,double*);
195
196
197
198
199 /* AMG: */
200 void Paso_Solver_AMG_free(Paso_Solver_AMG * in);
201 Paso_Solver_AMG* Paso_Solver_getAMG(Paso_SparseMatrix * A_p,dim_t level,Paso_Options* options);
202 void Paso_Solver_solveAMG(Paso_Solver_AMG * amg, double * x, double * b);
203
204 /*******************************************/
205 void Paso_Solver_ILU_free(Paso_Solver_ILU * in);
206 Paso_Solver_ILU* Paso_Solver_getILU(Paso_SparseMatrix * A_p,bool_t verbose);
207 void Paso_Solver_solveILU(Paso_SparseMatrix * A, Paso_Solver_ILU * ilu, double * x, const double * b);
208
209 void Paso_Solver_RILU_free(Paso_Solver_RILU * in);
210 Paso_Solver_RILU* Paso_Solver_getRILU(Paso_SparseMatrix * A_p,bool_t verbose);
211 void Paso_Solver_solveRILU(Paso_Solver_RILU * rilu, double * x, double * b);
212
213
214
215 void Paso_Solver_AMLI_System_free(Paso_Solver_AMLI_System * in);
216 void Paso_Solver_AMLI_free(Paso_Solver_AMLI * in);
217 Paso_Solver_AMLI* Paso_Solver_getAMLI(Paso_SparseMatrix * A_p,dim_t level,Paso_Options* options);
218 void Paso_Solver_solveAMLI(Paso_Solver_AMLI * amli, double * x, double * b);
219
220 void Paso_Solver_updateIncompleteSchurComplement(Paso_SparseMatrix* A_CC, Paso_SparseMatrix *A_CF,double* invA_FF,index_t* A_FF_pivot, Paso_SparseMatrix *A_FC);
221
222
223 #endif /* #ifndef INC_PRECS */

  ViewVC Help
Powered by ViewVC 1.1.26