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

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

Parent Directory Parent Directory | Revision Log Revision Log


Revision 4325 - (show annotations)
Wed Mar 20 02:28:59 2013 UTC (6 years, 8 months ago) by jfenwick
File size: 5467 byte(s)
all MEMALLOCs removed from PASO. further restructuring required
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 /************************************************************************************/
18
19 /* Paso: preconditioner set up */
20
21 /************************************************************************************/
22
23 /* Author: Lutz Gross, l.gross@uq.edu.au */
24
25 /************************************************************************************/
26
27 #include "Paso.h"
28 #include "SystemMatrix.h"
29 #include "PasoUtil.h"
30 #include "Preconditioner.h"
31
32 /*********************************************************************************************************/
33
34 /* free space */
35
36 void Paso_Preconditioner_free(Paso_Preconditioner* in) {
37 if (in!=NULL) {
38 Paso_Preconditioner_Smoother_free(in->jacobi);
39 Paso_Preconditioner_Smoother_free(in->gs);
40 Paso_Preconditioner_AMG_Root_free(in->amg);
41 /*********************************/
42 Paso_Solver_ILU_free(in->ilu);
43 Paso_Solver_RILU_free(in->rilu);
44 /*********************************/
45
46 delete 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=new 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->amg=NULL;
64
65 /*********************************/
66 prec->rilu=NULL;
67 prec->ilu=NULL;
68 /*********************************/
69
70 if (options->verbose && options->use_local_preconditioner) printf("Paso: Apply preconditioner locally only.\n");
71
72 switch (options->preconditioner) {
73 default:
74 case PASO_JACOBI:
75 if (options->verbose) {
76 if (options->sweeps >0 ) {
77 printf("Paso_Preconditioner: Jacobi(%d) preconditioner is used.\n",options->sweeps);
78 } else {
79 printf("Paso_Preconditioner: Jacobi preconditioner is used.\n");
80 }
81 }
82 prec->jacobi=Paso_Preconditioner_Smoother_alloc(A, TRUE, options->use_local_preconditioner, options->verbose);
83 prec->type=PASO_JACOBI;
84 prec->sweeps=options->sweeps;
85 break;
86 case PASO_GS:
87 if (options->verbose) {
88 if (options->sweeps >0 ) {
89 printf("Paso_Preconditioner: Gauss-Seidel(%d) preconditioner is used.\n",options->sweeps);
90 } else {
91 printf("Paso_Preconditioner: Gauss-Seidel preconditioner is used.\n");
92 }
93 }
94 prec->gs=Paso_Preconditioner_Smoother_alloc(A, FALSE, options->use_local_preconditioner, options->verbose);
95 prec->type=PASO_GS;
96 prec->sweeps=options->sweeps;
97 break;
98 case PASO_BOOMERAMG:
99 case PASO_AMLI:
100 case PASO_AMG:
101 prec->amg=Paso_Preconditioner_AMG_Root_alloc(A, options);
102 prec->type=PASO_AMG;
103 break;
104
105 /*************************************************************************************************************/
106 case PASO_ILU0:
107 if (options->verbose) printf("Paso_Preconditioner: ILU preconditioner is used.\n");
108 prec->ilu=Paso_Solver_getILU(A->mainBlock,options->verbose);
109 prec->type=PASO_ILU0;
110 Esys_MPIInfo_noError(A->mpi_info);
111 break;
112 case PASO_RILU:
113 if (options->verbose) printf("Paso_Preconditioner: RILU preconditioner is used.\n");
114 prec->rilu=Paso_Solver_getRILU(A->mainBlock,options->verbose);
115 Esys_MPIInfo_noError(A->mpi_info);
116 prec->type=PASO_RILU;
117 break;
118 case PASO_NO_PRECONDITIONER:
119 if (options->verbose) printf("Paso_Preconditioner: no preconditioner is applied.\n");
120 prec->type=PASO_NO_PRECONDITIONER;
121 break;
122 }
123 }
124 if (! Esys_noError() ){
125 Paso_Preconditioner_free(prec);
126 return NULL;
127 } else {
128 return prec;
129 }
130 }
131
132 /* Applies the preconditioner. */
133 /* Has to be called within a parallel region. */
134 /* Barrier synchronization is performed before the evaluation to make sure that the input vector is available */
135 void Paso_Preconditioner_solve(Paso_Preconditioner* prec, Paso_SystemMatrix* A,double* x,double* b){
136 dim_t n=0;
137
138 switch (prec->type) {
139 default:
140 case PASO_JACOBI:
141 Paso_Preconditioner_Smoother_solve(A, prec->jacobi,x,b,prec->sweeps, FALSE);
142 break;
143 case PASO_GS:
144 Paso_Preconditioner_Smoother_solve(A, prec->gs,x,b,prec->sweeps, FALSE);
145 break;
146 case PASO_AMG:
147 Paso_Preconditioner_AMG_Root_solve(A, prec->amg ,x,b);
148 break;
149
150 /*=========================================================*/
151 case PASO_ILU0:
152 Paso_Solver_solveILU(A->mainBlock, prec->ilu,x,b);
153 break;
154 case PASO_RILU:
155 Paso_Solver_solveRILU(prec->rilu,x,b);
156 break;
157 case PASO_NO_PRECONDITIONER:
158 n = MIN(Paso_SystemMatrix_getTotalNumCols(A),Paso_SystemMatrix_getTotalNumRows(A));
159 Paso_Copy(n,x,b);
160 break;
161
162 }
163 }
164

Properties

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

  ViewVC Help
Powered by ViewVC 1.1.26