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

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

Parent Directory Parent Directory | Revision Log Revision Log


Revision 3125 - (hide annotations)
Tue Aug 31 05:37:30 2010 UTC (9 years ago) by gross
File MIME type: text/plain
File size: 7378 byte(s)
more code cleaning
1 gross 3120
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     /* jacobi preconditioner */
24    
25 gross 3125 typedef struct Paso_Preconditioner_Jacobi {
26 gross 3120 double* values;
27     index_t* pivot;
28 gross 3125 } Paso_Preconditioner_Jacobi;
29 gross 3120
30     /* GS preconditioner */
31 gross 3125 typedef struct Paso_Preconditioner_LocalGS {
32 gross 3120 double* diag;
33     double* buffer;
34     index_t* pivot;
35     dim_t sweeps;
36 gross 3125 } Paso_Preconditioner_LocalGS;
37 gross 3120
38 gross 3125 typedef struct Paso_Preconditioner_GS {
39     Paso_Preconditioner_LocalGS* localGS;
40 gross 3120 bool_t is_local;
41 gross 3125 } Paso_Preconditioner_GS;
42 gross 3120
43    
44    
45     /*===============================================*/
46     /* ILU preconditioner */
47     struct Paso_Solver_ILU {
48     double* factors;
49     };
50     typedef struct Paso_Solver_ILU Paso_Solver_ILU;
51    
52    
53    
54    
55     /* RILU preconditioner */
56     struct Paso_Solver_RILU {
57     dim_t n;
58     dim_t n_block;
59     dim_t n_F;
60     dim_t n_C;
61     double* inv_A_FF;
62     index_t* A_FF_pivot;
63     Paso_SparseMatrix * A_FC;
64     Paso_SparseMatrix * A_CF;
65     index_t* rows_in_F;
66     index_t* rows_in_C;
67     index_t* mask_F;
68     index_t* mask_C;
69     double* x_F;
70     double* b_F;
71     double* x_C;
72     double* b_C;
73     struct Paso_Solver_RILU * RILU_of_Schur;
74     };
75     typedef struct Paso_Solver_RILU Paso_Solver_RILU;
76    
77     struct Paso_Solver_Smoother {
78     dim_t ID;
79 gross 3125 Paso_Preconditioner_Jacobi* Jacobi;
80     Paso_Preconditioner_LocalGS* GS;
81 gross 3120 };
82     typedef struct Paso_Solver_Smoother Paso_Solver_Smoother;
83    
84     /* AMG preconditioner */
85     struct Paso_Solver_AMG {
86     dim_t n;
87     dim_t level;
88     bool_t coarsest_level;
89     dim_t n_block;
90     dim_t n_F;
91     dim_t n_C;
92    
93     Paso_SparseMatrix * A_FF;
94     Paso_SparseMatrix * A_FC;
95     Paso_SparseMatrix * A_CF;
96     Paso_SparseMatrix * W_FC;
97    
98     Paso_SparseMatrix * P;
99     Paso_SparseMatrix * R;
100    
101     index_t* rows_in_F;
102     index_t* rows_in_C;
103     index_t* mask_F;
104     index_t* mask_C;
105     double* x_F;
106     double* b_F;
107     double* x_C;
108     double* b_C;
109    
110     dim_t post_sweeps;
111     dim_t pre_sweeps;
112    
113     Paso_SparseMatrix * A;
114     Paso_SparseMatrix * AOffset1;
115     Paso_SparseMatrix * AUnrolled;
116     void* solver;
117     Paso_Solver_Smoother* Smoother;
118     struct Paso_Solver_AMG * AMG_of_Coarse;
119     };
120     typedef struct Paso_Solver_AMG Paso_Solver_AMG;
121    
122    
123     /* AMLI preconditioner */
124     struct Paso_Solver_AMLI {
125     dim_t n;
126     dim_t level;
127     bool_t coarsest_level;
128     dim_t n_block;
129     dim_t n_F;
130     dim_t n_C;
131     double* inv_A_FF;
132     index_t* A_FF_pivot;
133     Paso_SparseMatrix * A_FC;
134     Paso_SparseMatrix * A_CF;
135     index_t* rows_in_F;
136     index_t* rows_in_C;
137     index_t* mask_F;
138     index_t* mask_C;
139     double* x_F;
140     double* b_F;
141     double* x_C;
142     double* b_C;
143    
144     dim_t post_sweeps;
145     dim_t pre_sweeps;
146    
147     Paso_SparseMatrix * A;
148     Paso_SparseMatrix * AOffset1;
149     void* solver;
150 gross 3125 Paso_Preconditioner_Jacobi* GS;
151 gross 3120 struct Paso_Solver_AMLI * AMLI_of_Schur;
152     };
153     typedef struct Paso_Solver_AMLI Paso_Solver_AMLI;
154    
155    
156     /* AMLI preconditioner on blocks*/
157     struct Paso_Solver_AMLI_System {
158     dim_t block_size;
159     Paso_SparseMatrix *block[MAX_BLOCK_SIZE];
160     Paso_Solver_AMLI *amliblock[MAX_BLOCK_SIZE];
161     };
162     typedef struct Paso_Solver_AMLI_System Paso_Solver_AMLI_System;
163    
164    
165     /* AMG preconditioner on blocks*/
166     struct Paso_Solver_AMG_System {
167     dim_t block_size;
168     Paso_SparseMatrix *block[MAX_BLOCK_SIZE];
169     Paso_Solver_AMG *amgblock[MAX_BLOCK_SIZE];
170     };
171     typedef struct Paso_Solver_AMG_System Paso_Solver_AMG_System;
172    
173     /* general preconditioner interface */
174    
175     typedef struct Paso_Solver_Preconditioner {
176     dim_t type;
177    
178     /* jacobi preconditioner */
179 gross 3125 Paso_Preconditioner_Jacobi* jacobi;
180 gross 3120 /* Gauss-Seidel preconditioner */
181 gross 3125 Paso_Preconditioner_GS* gs;
182 gross 3120
183     /* ilu preconditioner */
184     Paso_Solver_ILU* ilu;
185     /* rilu preconditioner */
186     Paso_Solver_RILU* rilu;
187     /* amg preconditioner */
188     Paso_Solver_AMG* amg;
189     /* amg on System */
190     Paso_Solver_AMG_System* amgSystem;
191     /* amg preconditioner */
192     Paso_Solver_AMLI* amli;
193     /* amg on System */
194     Paso_Solver_AMLI_System* amliSystem;
195    
196     } Paso_Solver_Preconditioner;
197    
198     void Paso_Preconditioner_free(Paso_Solver_Preconditioner*);
199     Paso_Solver_Preconditioner* Paso_Preconditioner_alloc(Paso_SystemMatrix* A,Paso_Options* options);
200     void Paso_Preconditioner_solve(Paso_Solver_Preconditioner* prec, Paso_SystemMatrix* A,double*,double*);
201    
202 gross 3125 /* JACOBI */
203     Paso_Preconditioner_Jacobi* Paso_Preconditioner_Jacobi_alloc(Paso_SystemMatrix * A_p);
204     Paso_Preconditioner_Jacobi* Paso_Preconditioner_LocalJacobi_alloc(Paso_SparseMatrix * A_p);
205     void Paso_Preconditioner_Jacobi_free(Paso_Preconditioner_Jacobi * in);
206 gross 3120
207 gross 3125 void Paso_Preconditioner_Jacobi_solve(Paso_SystemMatrix * A_p, Paso_Preconditioner_Jacobi * prec, double * x, double * b);
208     void Paso_Preconditioner_LocalJacobi_solve(Paso_SparseMatrix * A_p, Paso_Preconditioner_Jacobi * prec, double * x, double * b);
209 gross 3120
210 gross 3125 /* GAUSS SEIDEL */
211     void Paso_Preconditioner_GS_free(Paso_Preconditioner_GS * in);
212     void Paso_Preconditioner_LocalGS_free(Paso_Preconditioner_LocalGS * in);
213     Paso_Preconditioner_GS* Paso_Preconditioner_GS_alloc(Paso_SystemMatrix * A_p, dim_t sweeps, bool_t is_local, bool_t verbose);
214 gross 3120
215 gross 3125 Paso_Preconditioner_LocalGS* Paso_Preconditioner_LocalGS_alloc(Paso_SparseMatrix * A_p, dim_t sweeps, bool_t verbose);
216     void Paso_Preconditioner_GS_solve(Paso_SystemMatrix* A, Paso_Preconditioner_GS * gs, double * x, const double * b);
217     void Paso_Preconditioner_LocalGS_solve(Paso_SparseMatrix* A, Paso_Preconditioner_LocalGS * gs, double * x, const double * b);
218    
219     void Paso_Preconditioner_LocalGS_Sweep(Paso_SparseMatrix* A, Paso_Preconditioner_LocalGS * gs, double * x);
220     void Paso_Preconditioner_LocalGS_Sweep_sequential(Paso_SparseMatrix* A, Paso_Preconditioner_LocalGS * gs, double * x);
221     void Paso_Preconditioner_LocalGS_Sweep_tiled(Paso_SparseMatrix* A, Paso_Preconditioner_LocalGS * gs, double * x);
222     void Paso_Preconditioner_LocalGS_Sweep_colored(Paso_SparseMatrix* A, Paso_Preconditioner_LocalGS * gs, double * x);
223    
224 gross 3120 /*******************************************/
225     void Paso_Solver_ILU_free(Paso_Solver_ILU * in);
226     Paso_Solver_ILU* Paso_Solver_getILU(Paso_SparseMatrix * A_p,bool_t verbose);
227     void Paso_Solver_solveILU(Paso_SparseMatrix * A, Paso_Solver_ILU * ilu, double * x, const double * b);
228    
229     void Paso_Solver_RILU_free(Paso_Solver_RILU * in);
230     Paso_Solver_RILU* Paso_Solver_getRILU(Paso_SparseMatrix * A_p,bool_t verbose);
231     void Paso_Solver_solveRILU(Paso_Solver_RILU * rilu, double * x, double * b);
232    
233     void Paso_Solver_AMG_System_free(Paso_Solver_AMG_System * in);
234     void Paso_Solver_AMG_free(Paso_Solver_AMG * in);
235     Paso_Solver_AMG* Paso_Solver_getAMG(Paso_SparseMatrix * A_p,dim_t level,Paso_Options* options);
236     void Paso_Solver_solveAMG(Paso_Solver_AMG * amg, double * x, double * b);
237    
238     void Paso_Solver_AMLI_System_free(Paso_Solver_AMLI_System * in);
239     void Paso_Solver_AMLI_free(Paso_Solver_AMLI * in);
240     Paso_Solver_AMLI* Paso_Solver_getAMLI(Paso_SparseMatrix * A_p,dim_t level,Paso_Options* options);
241     void Paso_Solver_solveAMLI(Paso_Solver_AMLI * amli, 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