/[escript]/trunk/paso/src/Preconditioner.c
ViewVC logotype

Diff of /trunk/paso/src/Preconditioner.c

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

revision 3000 by lgao, Tue Mar 30 03:10:13 2010 UTC revision 3005 by gross, Thu Apr 22 05:59:31 2010 UTC
# Line 44  void Paso_Preconditioner_free(Paso_Solve Line 44  void Paso_Preconditioner_free(Paso_Solve
44        MEMFREE(in);        MEMFREE(in);
45      }      }
46  }  }
47    
48    
49  /*  call the iterative solver: */  /*  call the iterative solver: */
50    
51  void Paso_Solver_setPreconditioner(Paso_SystemMatrix* A,Paso_Options* options) {  void Paso_Solver_setPreconditioner(Paso_SystemMatrix* A,Paso_Options* options) {
# Line 54  void Paso_Solver_setPreconditioner(Paso_ Line 56  void Paso_Solver_setPreconditioner(Paso_
56      if (A->solver==NULL) {      if (A->solver==NULL) {
57          /* allocate structure to hold preconditioner */          /* allocate structure to hold preconditioner */
58          prec=MEMALLOC(1,Paso_Solver_Preconditioner);          prec=MEMALLOC(1,Paso_Solver_Preconditioner);
59          prec->amgSystem=MEMALLOC(1,Paso_Solver_AMG_System);  
         prec->amliSystem=MEMALLOC(1,Paso_Solver_AMLI_System);  
60          if (Paso_checkPtr(prec)) return;          if (Paso_checkPtr(prec)) return;
61                    
         if (Paso_checkPtr(prec->amliSystem)) return;  
           
         if (Paso_checkPtr(prec->amgSystem)) return;  
           
62          prec->type=UNKNOWN;          prec->type=UNKNOWN;
63          prec->rilu=NULL;          prec->rilu=NULL;
64          prec->ilu=NULL;          prec->ilu=NULL;
# Line 69  void Paso_Solver_setPreconditioner(Paso_ Line 66  void Paso_Solver_setPreconditioner(Paso_
66          prec->gs=NULL;          prec->gs=NULL;
67          prec->amg=NULL;          prec->amg=NULL;
68          prec->amli=NULL;          prec->amli=NULL;
69                    prec->amgSystem=NULL;
70          prec->amgSystem->block_size=A->row_block_size;          prec->amliSystem=NULL;
         for (i=0;i<A->row_block_size;++i) {  
           prec->amgSystem->amgblock[i]=NULL;  
           prec->amgSystem->block[i]=NULL;  
         }  
           
         prec->amliSystem->block_size=A->row_block_size;  
         for (i=0;i<A->row_block_size;++i) {  
           prec->amliSystem->amliblock[i]=NULL;  
           prec->amliSystem->block[i]=NULL;  
         }  
71    
72          A->solver=prec;          A->solver=prec;
73          switch (options->preconditioner) {          switch (options->preconditioner) {
# Line 107  void Paso_Solver_setPreconditioner(Paso_ Line 94  void Paso_Solver_setPreconditioner(Paso_
94                prec->type=PASO_GS;                prec->type=PASO_GS;
95                break;                break;
96              case PASO_AMG:              case PASO_AMG:
97              prec->amgSystem=MEMALLOC(1,Paso_Solver_AMG_System);
98              if (Paso_checkPtr(prec->amgSystem)) return;
99              prec->amgSystem->block_size=A->row_block_size;
100              for (i=0;i<A->row_block_size;++i) {
101            prec->amgSystem->amgblock[i]=NULL;
102            prec->amgSystem->block[i]=NULL;
103              }
104                                    
105                if (options->verbose) printf("AMG preconditioner is used.\n");                if (options->verbose) printf("AMG preconditioner is used.\n");            
                 
106                /*For performace reasons we check if block_size is one. If yes, then we do not need to separate blocks.*/                /*For performace reasons we check if block_size is one. If yes, then we do not need to separate blocks.*/
107                if (A->row_block_size==1) {                if (A->row_block_size==1) {
108                  prec->amg=Paso_Solver_getAMG(A->mainBlock,options->level_max,options);                    prec->amg=Paso_Solver_getAMG(A->mainBlock,options->level_max,options);  
# Line 129  void Paso_Solver_setPreconditioner(Paso_ Line 122  void Paso_Solver_setPreconditioner(Paso_
122                break;                break;
123    
124              case PASO_AMLI:              case PASO_AMLI:
125                if (options->verbose) printf("AMLI preconditioner is used.\n");            prec->amliSystem->block_size=A->row_block_size;
126              if (Paso_checkPtr(prec->amliSystem)) return;
127    
128              for (i=0;i<A->row_block_size;++i) {
129              prec->amliSystem->amliblock[i]=NULL;
130              prec->amliSystem->block[i]=NULL;
131              }
132    
133              if (options->verbose) printf("AMLI preconditioner is used.\n");
134                                
135                /*For performace reasons we check if block_size is one. If yes, then we do not need to separate blocks.*/                /*For performace reasons we check if block_size is one. If yes, then we do not need to separate blocks.*/
136                if (A->row_block_size==1) {                if (A->row_block_size==1) {
# Line 164  void Paso_Solver_solvePreconditioner(Pas Line 165  void Paso_Solver_solvePreconditioner(Pas
165            
166      xx=MEMALLOC(A->row_block_size,double*);      xx=MEMALLOC(A->row_block_size,double*);
167      bb=MEMALLOC(A->row_block_size,double*);      bb=MEMALLOC(A->row_block_size,double*);
168      if (Paso_checkPtr(xx) && Paso_checkPtr(bb)) return;      if (Paso_checkPtr(xx) || Paso_checkPtr(bb)) return;
169    
170      switch (prec->type) {      switch (prec->type) {
171          default:          default:
# Line 213  void Paso_Solver_solvePreconditioner(Pas Line 214  void Paso_Solver_solvePreconditioner(Pas
214                     Paso_SparseMatrix_MatrixVector_CSR_OFFSET0(DBLE(-1), A->mainBlock, x, DBLE(1), bnew);                     Paso_SparseMatrix_MatrixVector_CSR_OFFSET0(DBLE(-1), A->mainBlock, x, DBLE(1), bnew);
215                     /* Go round again*/                     /* Go round again*/
216                     Paso_Solver_solveGS(prec->gs,x,bnew);                     Paso_Solver_solveGS(prec->gs,x,bnew);
217                     sweeps=sweeps-1;                     sweeps--;
218                  }                  }
219                  MEMFREE(bnew);                  MEMFREE(bnew);
220             }             }

Legend:
Removed from v.3000  
changed lines
  Added in v.3005

  ViewVC Help
Powered by ViewVC 1.1.26