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

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

Parent Directory Parent Directory | Revision Log Revision Log | View Patch Patch

revision 2168 by artak, Wed Dec 17 00:50:37 2008 UTC revision 2275 by artak, Tue Feb 17 04:25:08 2009 UTC
# Line 3  Line 3 
3  #include "paso/Solver.h"  #include "paso/Solver.h"
4  #include "paso/SystemMatrix.h"  #include "paso/SystemMatrix.h"
5  #include "Paso_tests.h"  #include "Paso_tests.h"
6    #include <math.h>
7    
8    #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  int main (int argc, char *argv[]) {  int main (int argc, char *argv[]) {
22      Paso_SystemMatrix *A = NULL;      Paso_SystemMatrix *A = NULL;
23      double *b,*x;      double *b,*x,*x_ref;
24      dim_t i,n,level=1;      dim_t i,n,level=1;
25        double *error;
26        double Lsuperror;
27      
28        Paso_Options options;
29        Paso_Options_setDefaults(&options);
30    
31        
32      if (argc<2) {      if (argc<2) {
33          fprintf(stderr,"Please enter the filename\n");          fprintf(stderr,"Please enter the filename\n");
# Line 21  int main (int argc, char *argv[]) { Line 41  int main (int argc, char *argv[]) {
41      A=MEMALLOC(1,Paso_SystemMatrix);      A=MEMALLOC(1,Paso_SystemMatrix);
42        
43      A=Paso_SystemMatrix_loadMM_toCSR(argv[1]);      A=Paso_SystemMatrix_loadMM_toCSR(argv[1]);
44    
45      if (A==NULL) {      if (A==NULL) {
46        fprintf(stderr,"CSR Matrix not Loaded\n");        fprintf(stderr,"CSR Matrix not Loaded\n");
47        return 0;        return 0;
# Line 28  int main (int argc, char *argv[]) { Line 49  int main (int argc, char *argv[]) {
49      n=Paso_SystemMatrix_getTotalNumRows(A);      n=Paso_SystemMatrix_getTotalNumRows(A);
50      b=MEMALLOC(n,double);      b=MEMALLOC(n,double);
51      x=MEMALLOC(n,double);      x=MEMALLOC(n,double);
52        x_ref=MEMALLOC(n,double);
53        error=MEMALLOC(n,double);
54            
55      if(argc==4) {      if(argc==4) {
56          Paso_RHS_loadMM_toCSR(argv[3],b,n);          Paso_RHS_loadMM_toCSR(argv[3],b,n);
57            Paso_test_run(A,b,level);
58      }      }
59      else {      else {
60          for(i=0;i<n;i++) {          for(i=0;i<n;i++) {
61            b[i]=1;            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        }
99    
100    
101            
     Paso_test_run(A,b,level);  
102      MEMFREE(b);      MEMFREE(b);
103      MEMFREE(x);      MEMFREE(x);
104        MEMFREE(x_ref);
105        MEMFREE(error);
106      Paso_SystemMatrix_free(A);      Paso_SystemMatrix_free(A);
107  return 1;  return 1;
108  }                                                                                                                                                                                                      }
109    
110    

Legend:
Removed from v.2168  
changed lines
  Added in v.2275

  ViewVC Help
Powered by ViewVC 1.1.26