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

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

Parent Directory Parent Directory | Revision Log Revision Log


Revision 3402 - (show annotations)
Tue Dec 7 07:36:12 2010 UTC (8 years, 9 months ago) by gross
File MIME type: text/plain
File size: 6409 byte(s)
AMG support now block size > 3 (requires clapack or MKL).

additional diagnostics for AMG 


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 /**************************************************************/
16
17 /* Paso: SystemMatrix: sets-up the preconditioner */
18
19 /**************************************************************/
20
21 /* Copyrights by ACcESS Australia 2003/04 */
22 /* Author: Lutz Gross, l.gross@uq.edu.au */
23
24 /**************************************************************/
25
26 #include "Paso.h"
27 #include "SystemMatrix.h"
28 #include "PasoUtil.h"
29 #include "Preconditioner.h"
30
31 /***********************************************************************************/
32
33 /* free space */
34
35 void Paso_Preconditioner_free(Paso_Preconditioner* in) {
36 if (in!=NULL) {
37 Paso_Preconditioner_Smoother_free(in->jacobi);
38 Paso_Preconditioner_Smoother_free(in->gs);
39 Paso_Preconditioner_LocalAMG_free(in->localamg);
40 Paso_Preconditioner_LocalSmoother_free(in->localamgsubstitute);
41 /*********************************/
42 Paso_Solver_ILU_free(in->ilu);
43 Paso_Solver_RILU_free(in->rilu);
44 /*********************************/
45
46 MEMFREE(in);
47 }
48 }
49
50 Paso_Preconditioner* Paso_Preconditioner_alloc(Paso_SystemMatrix* A,Paso_Options* options) {
51
52 Paso_Preconditioner* prec=NULL;
53
54 prec=MEMALLOC(1,Paso_Preconditioner);
55
56 if (! Esys_checkPtr(prec)) {
57
58 prec->type=UNKNOWN;
59
60
61 prec->jacobi=NULL;
62 prec->gs=NULL;
63 prec->localamg=NULL;
64 prec->localamgsubstitute=NULL;
65
66 /*********************************/
67 prec->rilu=NULL;
68 prec->ilu=NULL;
69 /*********************************/
70
71 if (options->verbose && options->use_local_preconditioner) printf("Paso: Apply preconditioner locally only.\n");
72
73 switch (options->preconditioner) {
74 default:
75 case PASO_JACOBI:
76 if (options->verbose) printf("Paso_Preconditioner: Jacobi(%d) preconditioner is used.\n",options->sweeps);
77 prec->jacobi=Paso_Preconditioner_Smoother_alloc(A, TRUE, options->use_local_preconditioner, options->verbose);
78 prec->type=PASO_JACOBI;
79 prec->sweeps=options->sweeps;
80 break;
81 case PASO_GS:
82 if (options->verbose) printf("Paso_Preconditioner: Gauss-Seidel(%d) preconditioner is used.\n",options->sweeps);
83 prec->gs=Paso_Preconditioner_Smoother_alloc(A, FALSE, options->use_local_preconditioner, options->verbose);
84 prec->type=PASO_GS;
85 prec->sweeps=options->sweeps;
86 break;
87 case PASO_AMLI:
88 case PASO_AMG:
89 if (options->verbose) {
90 printf("Paso_Preconditioner: AMG preconditioner is used. Smoother is ");
91 if (options->smoother == PASO_JACOBI) {
92 printf("Jacobi");
93 } else {
94 printf("Gauss-Seidel");
95 }
96 printf(" with %d/%d pre/post sweeps.\n",options->pre_sweeps, options->post_sweeps);
97 }
98 prec->localamg=Paso_Preconditioner_LocalAMG_alloc(A->mainBlock,1,options);
99 prec->sweeps=options->sweeps;
100 /* if NULL is returned (and no error) no AMG has been constructed because the system is too small or not big enough
101 we now use the Smoother as a preconditioner */
102 if ( Esys_noError() ) {
103 if (prec->localamg == NULL) {
104 if (options->verbose) {
105 if (options->smoother == PASO_JACOBI) {
106 printf("Paso_Preconditioner: Jacobi(%d) preconditioner is used.\n",prec->sweeps);
107 } else {
108 printf("PPaso_Preconditioner: Gauss-Seidel(%d) preconditioner is used.\n",prec->sweeps);
109 }
110 }
111 options->num_level=0;
112 prec->localamgsubstitute=Paso_Preconditioner_LocalSmoother_alloc(A->mainBlock, (options->smoother == PASO_JACOBI), options->verbose);
113 } else {
114 options->num_level=Paso_Preconditioner_LocalAMG_getMaxLevel(prec->localamg);
115 options->coarse_level_sparsity=Paso_Preconditioner_LocalAMG_getCoarseLevelSparsity(prec->localamg);
116 options->num_coarse_unknowns=Paso_Preconditioner_LocalAMG_getNumCoarseUnknwons(prec->localamg);
117 }
118 }
119 prec->type=PASO_AMG;
120 Esys_MPIInfo_noError(A->mpi_info);
121 break;
122
123 /***************************************************************************************/
124 case PASO_ILU0:
125 if (options->verbose) printf("Paso_Preconditioner: ILU preconditioner is used.\n");
126 prec->ilu=Paso_Solver_getILU(A->mainBlock,options->verbose);
127 prec->type=PASO_ILU0;
128 Esys_MPIInfo_noError(A->mpi_info);
129 break;
130 case PASO_RILU:
131 if (options->verbose) printf("Paso_Preconditioner: RILU preconditioner is used.\n");
132 prec->rilu=Paso_Solver_getRILU(A->mainBlock,options->verbose);
133 Esys_MPIInfo_noError(A->mpi_info);
134 prec->type=PASO_RILU;
135 break;
136 }
137 }
138 if (! Esys_MPIInfo_noError(A->mpi_info ) ){
139 Paso_Preconditioner_free(prec);
140 return NULL;
141 } else {
142 return prec;
143 }
144 }
145
146 /* applies the preconditioner */
147 /* has to be called within a parallel reqion */
148 /* barrier synchronization is performed before the evaluation to make sure that the input vector is available */
149 void Paso_Preconditioner_solve(Paso_Preconditioner* prec, Paso_SystemMatrix* A,double* x,double* b){
150
151 switch (prec->type) {
152 default:
153 case PASO_JACOBI:
154 Paso_Preconditioner_Smoother_solve(A, prec->jacobi,x,b,prec->sweeps, FALSE);
155 break;
156 case PASO_GS:
157 Paso_Preconditioner_Smoother_solve(A, prec->gs,x,b,prec->sweeps, FALSE);
158 break;
159 case PASO_AMG:
160 if (prec->localamg == NULL) {
161 Paso_Preconditioner_LocalSmoother_solve(A->mainBlock, prec->localamgsubstitute,x,b,prec->sweeps, FALSE);
162 } else {
163 Paso_Preconditioner_LocalAMG_solve(A->mainBlock, prec->localamg,x,b);
164 }
165 break;
166
167 /*=========================================================*/
168 case PASO_ILU0:
169 Paso_Solver_solveILU(A->mainBlock, prec->ilu,x,b);
170 break;
171 case PASO_RILU:
172 Paso_Solver_solveRILU(prec->rilu,x,b);
173 break;
174
175 }
176
177 }

Properties

Name Value
svn:eol-style native
svn:keywords Author Date Id Revision

  ViewVC Help
Powered by ViewVC 1.1.26