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

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

Parent Directory Parent Directory | Revision Log Revision Log


Revision 2490 - (show annotations)
Tue Jun 23 06:02:26 2009 UTC (10 years, 3 months ago) by artak
File MIME type: text/plain
File size: 6979 byte(s)
adjusting to new Paso_Options.
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 coupling parameter for AMG]
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 rhs matrix - right hand side vector in CSR Matrix Market format.
23 coupling parameter for AMG - this is the threshold value used in AMG in courenening process. Default is 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: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 break;
91 case 'f':
92 filename = optarg;
93 A=MEMALLOC(1,Paso_SystemMatrix);
94 A=Paso_SystemMatrix_loadMM_toCSR(filename);
95 n=Paso_SystemMatrix_getTotalNumRows(A);
96 b=MEMALLOC(n,double);
97 x=MEMALLOC(n,double);
98 x_ref=MEMALLOC(n,double);
99 error=MEMALLOC(n,double);
100 for(i=0;i<n;i++) {
101 x_ref[i]=cos(i*PI/n);
102 }
103 Paso_SystemMatrix_MatrixVector_CSR_OFFSET0(DBLE(1),A,x_ref,DBLE(0),b);
104 break;
105 case 'r':
106 rhs=optarg;
107 if (A==NULL) {
108 fprintf(stderr,"Left hand side not loaded yet.\n");
109 break;
110 }
111 n=Paso_SystemMatrix_getTotalNumRows(A);
112 Paso_RHS_loadMM_toCSR(rhs,b,n);
113 break;
114 case 'l':
115 level=atoi(optarg);
116 break;
117 case 'c':
118 options.coarsening_threshold=atof(optarg);
119 break;
120 case '?':
121 printf("unknown arg %c\n", optopt);
122 break;
123 case 'h':
124 printf("Usage: PasoTests -f filename [-s solver] [-p preconditioner] [-l level] [-r rhs matrix] [-c coupling parameter for AMG]\n");
125 printf("\t filename - matrix to be loaded in CSR Matrix-Market format\n");
126 printf("\t solver - PCG, GMRES, PRES20, TFQMR and MINRES\n");
127 printf("\t preconditioner - ILU0, RILU, JACOBI, GS and AMG\n");
128 printf("\t level - options are 1,2 and 3\n");
129 printf("\t\t 0 - default option just solves with default of specified parameters\n");
130 printf("\t\t 1 - test all solvers with default preconditioner\n");
131 printf("\t\t 2 - test all preconditioners with default solver\n");
132 printf("\t\t 3 - compare solution obtained by using AMG and Jacobi precondioners\n");
133 printf("\trhs matrix - right hand side vector in CSR Matrix Market format.\n");
134 printf("\tcoupling parameter for AMG - this is the threshold value used in AMG in courenening process. Default is 0.05.\n");
135 break;
136 }
137 }
138
139 if (A==NULL) {
140 /*fprintf(stderr,"CSR Matrix not Loaded\n");*/
141 return 0;
142 }
143
144
145 if (level==0) {
146 Paso_solve(A,x,b,&options);
147 }
148 else if (level==3) {
149 options.method=PASO_PCG;
150 options.verbose=TRUE;
151 options.preconditioner=PASO_JACOBI;
152
153 Paso_solve(A,x,b,&options);
154
155 for(i=0;i<n;i++) {
156 error[i]=x[i]-x_ref[i];
157 }
158 Lsuperror=Lsup(error,n)/Lsup(x_ref,n);
159 fprintf(stdout,"Lsup error Jacobi %e\n",Lsuperror);
160
161 A->solver=NULL;
162 options.method=PASO_PCG;
163 options.preconditioner=PASO_AMG;
164
165 Paso_solve(A,x,b,&options);
166
167 for(i=0;i<n;i++) {
168 error[i]=x[i]-x_ref[i];
169 }
170 Lsuperror=Lsup(error,n)/Lsup(x_ref,n);
171 fprintf(stdout,"Lsup error AMG %e\n",Lsuperror);
172 }
173 else if (level==4) {
174 options.method=PASO_PCG;
175 options.verbose=TRUE;
176 options.preconditioner=PASO_AMG;
177
178 Paso_Solver_setPreconditioner(A,&options);
179 Paso_Solver_solvePreconditioner(A,x,b);
180
181
182 for(i=0;i<n;i++) {
183 error[i]=x[i]-x_ref[i];
184 }
185 Lsuperror=Lsup(error,n)/Lsup(x_ref,n);
186 fprintf(stdout,"Lsup error AMG as a solver %e\n\n",Lsuperror);
187
188 A->solver=NULL;
189 options.method=PASO_PCG;
190 options.preconditioner=PASO_AMG;
191
192 Paso_solve(A,x,b,&options);
193
194 for(i=0;i<n;i++) {
195 error[i]=x[i]-x_ref[i];
196 }
197 Lsuperror=Lsup(error,n)/Lsup(x_ref,n);
198 fprintf(stdout,"Lsup error AMG %e\n",Lsuperror);
199 }
200 else {
201 Paso_test_run(A,b,level);
202 }
203
204 if (A!=NULL) {
205 MEMFREE(b);
206 MEMFREE(x);
207 MEMFREE(x_ref);
208 MEMFREE(error);
209 Paso_SystemMatrix_free(A);
210 }
211
212 return 1;
213 }
214
215

  ViewVC Help
Powered by ViewVC 1.1.26