/[escript]/trunk/paso/profiling/Test.c
ViewVC logotype

Contents of /trunk/paso/profiling/Test.c

Parent Directory Parent Directory | Revision Log Revision Log


Revision 2820 - (show annotations)
Thu Dec 10 05:03:11 2009 UTC (11 years, 6 months ago) by artak
File MIME type: text/plain
File size: 7788 byte(s)
Paso Testing tools now can select coarseing.The builld_PasoTests target is also added to scons. 
1 #include <stdio.h>
2 #include <unistd.h>
3 #include "paso/Common.h"
4 #include "paso/Solver.h"
5 #include "paso/SystemMatrix.h"
6 #include "Paso_tests.h"
7 #include "getopt.h"
8 #include <math.h>
9
10 #define PI (3.141592653589793)
11
12 /*
13 Usage: PasoTests -f filename [-s solver] [-p preconditioner] [-l level] [-r rhs matrix] [-c coarsening method] [-t threshold parameter for AMG coarsening]
14 filename - matrix to be loaded in CSR Matrix-Market format
15 solver - PCG, GMRES, PRES20, TFQMR and MINRES
16 preconditioner - ILU0, RILU, JACOBI, GS and AMG
17 level - options are 1,2 and 3
18 0 - default option just solves with default of specified parameters
19 1 - test all solvers with default preconditioner
20 2 - test all preconditioners with default solver
21 3 - compare solution obtained by using AMG and Jacobi precondioners
22 coarsening method - YS, RS, AGG and STD.
23 threshold parameter for AMG - this is the threshold value used in AMG in courenening process. Default is 0.25. For YS and AGG, please, use 0.05.
24 */
25
26 double Lsup(double* x, int n) {
27 double max=0;
28 int i;
29
30 for (i=0;i<n;i++) {
31 max=MAX(ABS(x[i]),max);
32 }
33
34 return max;
35 }
36
37 int main (int argc, char *argv[]) {
38 Paso_SystemMatrix *A = NULL;
39 double *b,*x,*x_ref;
40 dim_t i,n,level=0;
41 double *error;
42 double Lsuperror;
43
44 int c;
45 char *filename,*solver,*prec,*rhs;
46 extern char *optarg;
47 extern int optopt;
48
49 Paso_Options options;
50 Paso_Options_setDefaults(&options);
51
52
53
54 options.verbose=TRUE;
55
56 while ((c = getopt(argc, argv, "s:p:f:r:l:t:c:h")) != -1) {
57 switch(c) {
58 case 's':
59 solver=optarg;
60 if (strcmp(solver,"PCG")==0)
61 options.method=PASO_PCG;
62 else if (strcmp(solver,"GMRES")==0)
63 options.method=PASO_GMRES;
64 else if (strcmp(solver,"PRES20")==0)
65 options.method=PASO_PRES20;
66 else if (strcmp(solver,"BICGSTAB")==0)
67 options.method=PASO_BICGSTAB;
68 else if (strcmp(solver,"TFQMR")==0)
69 options.method=PASO_TFQMR;
70 else if (strcmp(solver,"MINRES")==0)
71 options.method=PASO_MINRES;
72 else if (strcmp(solver,"DIRECT")==0) {
73 options.package=PASO_MKL;
74 options.verbose=1;
75 }
76 break;
77 case 'p':
78 prec=optarg;
79 if (strcmp(prec,"JACOBI")==0)
80 options.preconditioner=PASO_JACOBI;
81 else if (strcmp(prec,"RILU")==0)
82 options.preconditioner=PASO_RILU;
83 else if (strcmp(prec,"ILU0")==0)
84 options.preconditioner=PASO_ILU0;
85 else if (strcmp(prec,"GS")==0)
86 options.preconditioner=PASO_GS;
87 else if (strcmp(prec,"AMG")==0) {
88 options.preconditioner=PASO_AMG;
89 }
90 else if (strcmp(prec,"AMLI")==0) {
91 options.preconditioner=PASO_AMG;
92 }
93 break;
94 case 'f':
95 filename = optarg;
96 A=MEMALLOC(1,Paso_SystemMatrix);
97 A=Paso_SystemMatrix_loadMM_toCSR(filename);
98 n=Paso_SystemMatrix_getTotalNumRows(A);
99 b=MEMALLOC(n,double);
100 x=MEMALLOC(n,double);
101 x_ref=MEMALLOC(n,double);
102 error=MEMALLOC(n,double);
103 for(i=0;i<n;i++) {
104 x_ref[i]=cos(i*PI/n);
105 }
106 Paso_SystemMatrix_MatrixVector_CSR_OFFSET0(DBLE(1),A,x_ref,DBLE(0),b);
107 break;
108 case 'r':
109 rhs=optarg;
110 if (A==NULL) {
111 fprintf(stderr,"System matrix is not loaded yet.\n");
112 break;
113 }
114 n=Paso_SystemMatrix_getTotalNumRows(A);
115 Paso_RHS_loadMM_toCSR(rhs,b,n);
116 break;
117 case 'l':
118 level=atoi(optarg);
119 break;
120 case 't':
121 options.coarsening_threshold=atof(optarg);
122 case 'c':
123 prec=optarg;
124 if (strcmp(prec,"RS")==0)
125 options.coarsening_method=PASO_RUGE_STUEBEN_COARSENING;
126 else if (strcmp(prec,"YS")==0)
127 options.coarsening_method=PASO_YAIR_SHAPIRA_COARSENING;
128 else if (strcmp(prec,"AGG")==0)
129 options.coarsening_method=PASO_AGGREGATION_COARSENING;
130 else if (strcmp(prec,"STD")==0)
131 options.coarsening_method=PASO_STANDARD_COARSENING;
132 break;
133 case '?':
134 printf("unknown arg %c\n", optopt);
135 break;
136 case 'h':
137 printf("Usage: PasoTests -f filename [-s solver] [-p preconditioner] [-l level] [-r rhs vector] [-c coarsening method] [-t threshold parameter for AMG coarsening] \n");
138 printf("\t filename - matrix to be loaded in CSR Matrix-Market format\n");
139 printf("\t solver - PCG, GMRES, PRES20, TFQMR and MINRES\n");
140 printf("\t preconditioner - ILU0, RILU, JACOBI, GS and AMG\n");
141 printf("\t level - options are 1,2 and 3\n");
142 printf("\t\t 0 - default option just solves with default of specified parameters\n");
143 printf("\t\t 1 - test all solvers with default preconditioner\n");
144 printf("\t\t 2 - test all preconditioners with default solver\n");
145 printf("\t\t 3 - compare solution obtained by using AMG and Jacobi precondioners\n");
146 printf("\trhs vector - right hand side vector in CSR Matrix Market format.\n");
147 printf("\tcoarsening method - YS, RS, AGG and STD.\n");
148 printf("\tthreshold parameter for AMG - this is the threshold value used in AMG in courenening process. Default is 0.25. For YS and AGG, please, use 0.05.\n");
149 break;
150 }
151 }
152
153 if (A==NULL) {
154 /*fprintf(stderr,"CSR Matrix not Loaded\n");*/
155 return 0;
156 }
157
158
159 if (level==0) {
160 Paso_solve(A,x,b,&options);
161 }
162 else if (level==3) {
163 options.method=PASO_PCG;
164 options.verbose=TRUE;
165 options.preconditioner=PASO_JACOBI;
166
167 Paso_solve(A,x,b,&options);
168
169 for(i=0;i<n;i++) {
170 error[i]=x[i]-x_ref[i];
171 }
172 Lsuperror=Lsup(error,n)/Lsup(x_ref,n);
173 fprintf(stdout,"Lsup error Jacobi %e\n",Lsuperror);
174
175 A->solver=NULL;
176 options.method=PASO_PCG;
177 options.preconditioner=PASO_AMG;
178
179 Paso_solve(A,x,b,&options);
180
181 for(i=0;i<n;i++) {
182 error[i]=x[i]-x_ref[i];
183 }
184 Lsuperror=Lsup(error,n)/Lsup(x_ref,n);
185 fprintf(stdout,"Lsup error AMG %e\n",Lsuperror);
186 }
187 else if (level==4) {
188 options.method=PASO_PCG;
189 options.verbose=TRUE;
190 options.preconditioner=PASO_AMG;
191
192 Paso_Solver_setPreconditioner(A,&options);
193 Paso_Solver_solvePreconditioner(A,x,b);
194
195
196 for(i=0;i<n;i++) {
197 error[i]=x[i]-x_ref[i];
198 }
199 Lsuperror=Lsup(error,n)/Lsup(x_ref,n);
200 fprintf(stdout,"Lsup error AMG as a solver %e\n\n",Lsuperror);
201
202 A->solver=NULL;
203 options.method=PASO_PCG;
204 options.preconditioner=PASO_AMG;
205
206 Paso_solve(A,x,b,&options);
207
208 for(i=0;i<n;i++) {
209 error[i]=x[i]-x_ref[i];
210 }
211 Lsuperror=Lsup(error,n)/Lsup(x_ref,n);
212 fprintf(stdout,"Lsup error AMG %e\n",Lsuperror);
213 }
214 else {
215 Paso_test_run(A,b,level);
216 }
217
218 if (A!=NULL) {
219 MEMFREE(b);
220 MEMFREE(x);
221 MEMFREE(x_ref);
222 MEMFREE(error);
223 Paso_SystemMatrix_free(A);
224 }
225
226 return 1;
227 }
228
229

  ViewVC Help
Powered by ViewVC 1.1.26