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

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

Parent Directory Parent Directory | Revision Log Revision Log


Revision 2275 - (hide annotations)
Tue Feb 17 04:25:08 2009 UTC (12 years, 5 months ago) by artak
File MIME type: text/plain
File size: 2415 byte(s)
Extra tests are added for checking the quality of the solution from Jacobi and AMG.
1 artak 2160 #include <stdio.h>
2     #include "paso/Common.h"
3     #include "paso/Solver.h"
4     #include "paso/SystemMatrix.h"
5 artak 2163 #include "Paso_tests.h"
6 artak 2275 #include <math.h>
7 artak 2160
8 artak 2275 #define PI (3.141592653589793)
9    
10     double Lsup(double* x, int n) {
11     double max=0;
12     int i;
13    
14     for (i=0;i<n;i++) {
15     max=MAX(ABS(x[i]),max);
16     }
17    
18     return max;
19     }
20    
21 artak 2160 int main (int argc, char *argv[]) {
22     Paso_SystemMatrix *A = NULL;
23 artak 2275 double *b,*x,*x_ref;
24 artak 2165 dim_t i,n,level=1;
25 artak 2275 double *error;
26     double Lsuperror;
27 artak 2160
28 artak 2275 Paso_Options options;
29     Paso_Options_setDefaults(&options);
30    
31    
32 artak 2160 if (argc<2) {
33     fprintf(stderr,"Please enter the filename\n");
34     return -1;
35     }
36 artak 2165
37 artak 2168 if (argc>2) {
38 artak 2165 level=atoi(argv[2]);
39     }
40 artak 2160
41     A=MEMALLOC(1,Paso_SystemMatrix);
42    
43     A=Paso_SystemMatrix_loadMM_toCSR(argv[1]);
44 artak 2275
45 artak 2160 if (A==NULL) {
46 artak 2165 fprintf(stderr,"CSR Matrix not Loaded\n");
47 artak 2160 return 0;
48     }
49     n=Paso_SystemMatrix_getTotalNumRows(A);
50     b=MEMALLOC(n,double);
51     x=MEMALLOC(n,double);
52 artak 2275 x_ref=MEMALLOC(n,double);
53     error=MEMALLOC(n,double);
54 artak 2168
55     if(argc==4) {
56     Paso_RHS_loadMM_toCSR(argv[3],b,n);
57 artak 2275 Paso_test_run(A,b,level);
58 artak 2160 }
59 artak 2168 else {
60     for(i=0;i<n;i++) {
61 artak 2275 x_ref[i]=cos(i*PI/n);
62     }
63     /* raw scaled vector update operation: out = alpha * A * in + beta * out */
64     Paso_SystemMatrix_MatrixVector_CSR_OFFSET0(DBLE(1),A,x_ref,DBLE(0),b);
65    
66     if (level==3) {
67     options.method=PASO_PCG;
68     options.verbose=FALSE;
69     options.preconditioner=PASO_JACOBI;
70    
71    
72    
73     Paso_solve(A,x,b,&options);
74    
75    
76     for(i=0;i<n;i++) {
77     error[i]=x[i]-x_ref[i];
78     }
79     Lsuperror=Lsup(error,n)/Lsup(x_ref,n);
80     fprintf(stdout,"Lsup error Jacobi %e\n",Lsuperror);
81    
82     A->solver=NULL;
83     options.method=PASO_PCG;
84     options.preconditioner=PASO_AMG;
85    
86     Paso_solve(A,x,b,&options);
87    
88     for(i=0;i<n;i++) {
89     error[i]=x[i]-x_ref[i];
90     }
91     Lsuperror=Lsup(error,n)/Lsup(x_ref,n);
92     fprintf(stdout,"Lsup error AMG %e\n",Lsuperror);
93     }
94     else {
95     Paso_test_run(A,b,level);
96     }
97    
98 artak 2168 }
99 artak 2275
100    
101 artak 2160
102     MEMFREE(b);
103     MEMFREE(x);
104 artak 2275 MEMFREE(x_ref);
105     MEMFREE(error);
106 artak 2160 Paso_SystemMatrix_free(A);
107     return 1;
108 artak 2275 }
109    
110    

  ViewVC Help
Powered by ViewVC 1.1.26