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

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

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

trunk/esys2/paso/src/Solvers/Solver.c revision 150 by jgs, Thu Sep 15 03:44:45 2005 UTC trunk/paso/src/Solvers/Solver.c revision 425 by gross, Tue Jan 10 04:10:39 2006 UTC
# Line 15  Line 15 
15  #include "SystemMatrix.h"  #include "SystemMatrix.h"
16  #include "Solver.h"  #include "Solver.h"
17    
 #if PTR_OFFSET !=0 || INDEX_OFFSET!=0  
 #error Paso library usage requires PTR_OFFSET=0 and INDEX_OFFSET=0  
 #endif  
   
   
18  /***********************************************************************************/  /***********************************************************************************/
19    
20  /*  free space */  /*  free space */
21    
22  void Paso_Solver_free(Paso_SystemMatrix* A) {  void Paso_Solver_free(Paso_SystemMatrix* A) {
23      Paso_Preconditioner_free(A->iterative);      Paso_Preconditioner_free(A->solver);
24      A->iterative=NULL;      A->solver=NULL;
25  }  }
26  /*  call the iterative solver: */  /*  call the iterative solver: */
27    
# Line 34  void Paso_Solver(Paso_SystemMatrix* A,do Line 29  void Paso_Solver(Paso_SystemMatrix* A,do
29      double norm2_of_b,tol,tolerance,time_iter,time_prec,*r=NULL,norm2_of_residual,last_norm2_of_residual,norm_max_of_b,      double norm2_of_b,tol,tolerance,time_iter,time_prec,*r=NULL,norm2_of_residual,last_norm2_of_residual,norm_max_of_b,
30             norm2_of_b_local,norm_max_of_b_local,norm2_of_residual_local,norm_max_of_residual_local,norm_max_of_residual,             norm2_of_b_local,norm_max_of_b_local,norm2_of_residual_local,norm_max_of_residual_local,norm_max_of_residual,
31             last_norm_max_of_residual,*scaling;             last_norm_max_of_residual,*scaling;
32      dim_t i,totIter,cntIter,method;       dim_t i,totIter,cntIter,method;
33      bool_t finalizeIteration;       bool_t finalizeIteration;
34      err_t errorCode;       err_t errorCode;
35      dim_t n_col = A->num_cols * A-> col_block_size;       dim_t n_col = A->num_cols * A-> col_block_size;
36      dim_t n_row = A->num_rows * A-> row_block_size;       dim_t n_row = A->num_rows * A-> row_block_size;
37    
38      tolerance=MAX(options->tolerance,EPSILON);       tolerance=MAX(options->tolerance,EPSILON);
39      Paso_resetError();       Paso_resetError();
40      method=Paso_Options_getSolver(options->method,PASO_PASO,options->symmetric);       method=Paso_Options_getSolver(options->method,PASO_PASO,options->symmetric);
41      /* check matrix type */       /* check matrix type */
42      if (A->type!=CSR) {       if ((A->type & MATRIX_FORMAT_CSC) || (A->type & MATRIX_FORMAT_OFFSET1) || (A->type & MATRIX_FORMAT_SYM) ) {
43        Paso_setError(TYPE_ERROR,"Iterative solver requires CSR format.");         Paso_setError(TYPE_ERROR,"Iterative solver requires CSR format with unsymmetric storage scheme and index offset 0.");
44      }       }
45      if (A->col_block_size != A->row_block_size) {       if (A->col_block_size != A->row_block_size) {
46         Paso_setError(TYPE_ERROR,"Iterative solver requires row and column block sizes to be equal.");          Paso_setError(TYPE_ERROR,"Iterative solver requires row and column block sizes to be equal.");
47       }       }
48       if (A->num_cols != A->num_rows) {       if (A->num_cols != A->num_rows) {
49         Paso_setError(TYPE_ERROR,"Iterative solver requires requires a square matrix.");          Paso_setError(TYPE_ERROR,"Iterative solver requires requires a square matrix.");
50         return;          return;
51       }       }
52       if (Paso_noError()) {       if (Paso_noError()) {
53          /* get normalization */          /* get normalization */
# Line 137  void Paso_Solver(Paso_SystemMatrix* A,do Line 132  void Paso_Solver(Paso_SystemMatrix* A,do
132                #pragma omp for private(i) schedule(static)                #pragma omp for private(i) schedule(static)
133                for (i = 0; i < n_row; i++) r[i]=b[i];                for (i = 0; i < n_row; i++) r[i]=b[i];
134            
135                Paso_SystemMatrix_MatrixVector(DBLE(-1), A, x, DBLE(1), r);                Paso_SystemMatrix_MatrixVector_CSR_OFFSET0(DBLE(-1), A, x, DBLE(1), r);
136            
137                norm2_of_residual_local = 0;                norm2_of_residual_local = 0;
138                norm_max_of_residual_local = 0;                norm_max_of_residual_local = 0;
# Line 187  void Paso_Solver(Paso_SystemMatrix* A,do Line 182  void Paso_Solver(Paso_SystemMatrix* A,do
182                          finalizeIteration = FALSE;                          finalizeIteration = FALSE;
183                    } else if (errorCode==SOLVER_MAXITER_REACHED) {                    } else if (errorCode==SOLVER_MAXITER_REACHED) {
184                           Paso_setError(WARNING,"maximum number of iteration step reached.\nReturned solution does not fulfil stopping criterion.");                           Paso_setError(WARNING,"maximum number of iteration step reached.\nReturned solution does not fulfil stopping criterion.");
185                             if (options->verbose) printf("Maximum number of iterations reached.!\n");
186                    } else if (errorCode == SOLVER_INPUT_ERROR ) {                    } else if (errorCode == SOLVER_INPUT_ERROR ) {
187                           Paso_setError(SYSTEM_ERROR,"illegal dimension in iterative solver.");                           Paso_setError(SYSTEM_ERROR,"illegal dimension in iterative solver.");
188                             if (options->verbose) printf("Internal error!\n");
189                    } else if ( errorCode == SOLVER_BREAKDOWN ) {                    } else if ( errorCode == SOLVER_BREAKDOWN ) {
190                    if (cntIter <= 1) {                    if (cntIter <= 1) {
191                             Paso_setError(ZERO_DIVISION_ERROR, "fatal break down in iterative solver.");                             Paso_setError(ZERO_DIVISION_ERROR, "fatal break down in iterative solver.");
192                               if (options->verbose) printf("Uncurable break down!\n");
193                        } else {                        } else {
194                        if (options->verbose) printf("Breakdown at iter %d (residual = %e). Restarting ...\n", cntIter+totIter, tol);                        if (options->verbose) printf("Breakdown at iter %d (residual = %e). Restarting ...\n", cntIter+totIter, tol);
195                            finalizeIteration = FALSE;                            finalizeIteration = FALSE;
196                        }                        }
197                    } else if (errorCode == SOLVER_MEMORY_ERROR) {                    } else if (errorCode == SOLVER_MEMORY_ERROR) {
198                         Paso_setError(MEMORY_ERROR,"memory allocation failed.");                     Paso_setError(MEMORY_ERROR,"memory allocation failed.");
199                           if (options->verbose) printf("Memory allocation failed!\n");
200                    } else if (errorCode !=SOLVER_NO_ERROR ) {                    } else if (errorCode !=SOLVER_NO_ERROR ) {
201                         Paso_setError(SYSTEM_ERROR,"unidentified error in iterative solver.");                     Paso_setError(SYSTEM_ERROR,"unidentified error in iterative solver.");
202                           if (options->verbose) printf("Unidentified error!\n");
203                    }                    }
204               } else {               } else {
205                  if (options->verbose) printf(". convergence! \n");                  if (options->verbose) printf(". convergence! \n");
# Line 215  void Paso_Solver(Paso_SystemMatrix* A,do Line 215  void Paso_Solver(Paso_SystemMatrix* A,do
215          }          }
216     }     }
217  }  }
   
 /*  
  * $Log$  
  * Revision 1.2  2005/09/15 03:44:40  jgs  
  * Merge of development branch dev-02 back to main trunk on 2005-09-15  
  *  
  * Revision 1.1.2.2  2005/09/07 00:59:09  gross  
  * some inconsistent renaming fixed to make the linking work.  
  *  
  * Revision 1.1.2.1  2005/09/05 06:29:49  gross  
  * These files have been extracted from finley to define a stand alone libray for iterative  
  * linear solvers on the ALTIX. main entry through Paso_solve. this version compiles but  
  * has not been tested yet.  
  *  
  *  
  */  

Legend:
Removed from v.150  
changed lines
  Added in v.425

  ViewVC Help
Powered by ViewVC 1.1.26