/[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

revision 2155 by phornby, Thu Nov 20 16:10:10 2008 UTC revision 2156 by gross, Mon Dec 15 05:09:02 2008 UTC
# Line 56  void Paso_Solver(Paso_SystemMatrix* A,do Line 56  void Paso_Solver(Paso_SystemMatrix* A,do
56     dim_t numEqua = Paso_SystemMatrix_getTotalNumRows(A);     dim_t numEqua = Paso_SystemMatrix_getTotalNumRows(A);
57     double blocktimer_precond, blocktimer_start = blocktimer_time();     double blocktimer_precond, blocktimer_start = blocktimer_time();
58    
      tolerance=MAX(options->tolerance,EPSILON);  
59       Paso_resetError();       Paso_resetError();
60         tolerance=options->tolerance;
61         if (tolerance < 100.* EPSILON) {
62           Paso_setError(VALUE_ERROR,"Paso_Solver: Tolerance is too small.");
63         }
64         if (tolerance >1.) {
65           Paso_setError(VALUE_ERROR,"Paso_Solver: Tolerance mut be less than one.");
66         }
67       method=Paso_Options_getSolver(options->method,PASO_PASO,options->symmetric);       method=Paso_Options_getSolver(options->method,PASO_PASO,options->symmetric);
68       /* check matrix type */       /* check matrix type */
69       if ((A->type & MATRIX_FORMAT_CSC) || (A->type & MATRIX_FORMAT_OFFSET1) || (A->type & MATRIX_FORMAT_SYM) ) {       if ((A->type & MATRIX_FORMAT_CSC) || (A->type & MATRIX_FORMAT_OFFSET1) || (A->type & MATRIX_FORMAT_SYM) ) {
70         Paso_setError(TYPE_ERROR,"Iterative solver requires CSR format with unsymmetric storage scheme and index offset 0.");         Paso_setError(TYPE_ERROR,"Paso_Solver: Iterative solver requires CSR format with unsymmetric storage scheme and index offset 0.");
71       }       }
72       if (A->col_block_size != A->row_block_size) {       if (A->col_block_size != A->row_block_size) {
73          Paso_setError(TYPE_ERROR,"Iterative solver requires row and column block sizes to be equal.");          Paso_setError(TYPE_ERROR,"Paso_Solver: Iterative solver requires row and column block sizes to be equal.");
74       }       }
75       if (Paso_SystemMatrix_getGlobalNumCols(A) != Paso_SystemMatrix_getGlobalNumRows(A)) {       if (Paso_SystemMatrix_getGlobalNumCols(A) != Paso_SystemMatrix_getGlobalNumRows(A)) {
76          Paso_setError(TYPE_ERROR,"Iterative solver requires a square matrix.");          Paso_setError(TYPE_ERROR,"Paso_Solver: Iterative solver requires a square matrix.");
77          return;          return;
78       }       }
79       /* this for testing only */       /* this for testing only */
# Line 124  void Paso_Solver(Paso_SystemMatrix* A,do Line 130  void Paso_Solver(Paso_SystemMatrix* A,do
130              if (options->verbose) printf("right hand side is identical zero.\n");              if (options->verbose) printf("right hand side is identical zero.\n");
131           } else {           } else {
132              if (options->verbose) {              if (options->verbose) {
133                 printf("l2/lmax-norm of right hand side is  %e/%e.\n",norm2_of_b,norm_max_of_b);                 printf("Paso_Solver: l2/lmax-norm of right hand side is  %e/%e.\n",norm2_of_b,norm_max_of_b);
134                 printf("l2/lmax-stopping criterion is  %e/%e.\n",norm2_of_b*tolerance,norm_max_of_b*tolerance);                 printf("Paso_Solver: l2/lmax-stopping criterion is  %e/%e.\n",norm2_of_b*tolerance,norm_max_of_b*tolerance);
135                 switch (method) {                 switch (method) {
136                 case PASO_BICGSTAB:                 case PASO_BICGSTAB:
137                    printf("Iterative method is BiCGStab.\n");                    printf("Paso_Solver: Iterative method is BiCGStab.\n");
138                    break;                    break;
139                 case PASO_PCG:                 case PASO_PCG:
140                    printf("Iterative method is PCG.\n");                    printf("Paso_Solver: Iterative method is PCG.\n");
141                    break;                    break;
142                 case PASO_TFQMR:                 case PASO_TFQMR:
143                    printf("Iterative method is TFQMR.\n");                    printf("Paso_Solver: Iterative method is TFQMR.\n");
144                    break;                    break;
145                 case PASO_MINRES:                 case PASO_MINRES:
146                    printf("Iterative method is MINRES.\n");                    printf("Paso_Solver: Iterative method is MINRES.\n");
147                    break;                    break;
148                 case PASO_PRES20:                 case PASO_PRES20:
149                    printf("Iterative method is PRES20.\n");                    printf("Paso_Solver: Iterative method is PRES20.\n");
150                    break;                    break;
151                 case PASO_GMRES:                 case PASO_GMRES:
152                    if (options->restart>0) {                    if (options->restart>0) {
153                       printf("Iterative method is GMRES(%d,%d)\n",options->truncation,options->restart);                       printf("Paso_Solver: Iterative method is GMRES(%d,%d)\n",options->truncation,options->restart);
154                    } else {                    } else {
155                       printf("Iterative method is GMRES(%d)\n",options->truncation);                       printf("Paso_Solver: Iterative method is GMRES(%d)\n",options->truncation);
156                    }                    }
157                    break;                    break;
158                 }                 }
# Line 205  void Paso_Solver(Paso_SystemMatrix* A,do Line 211  void Paso_Solver(Paso_SystemMatrix* A,do
211                      #endif                      #endif
212                      norm2_of_residual =sqrt(norm2_of_residual);                      norm2_of_residual =sqrt(norm2_of_residual);
213                            
214                    if (options->verbose) printf("Step %5d: l2/lmax-norm of residual is  %e/%e",totIter,norm2_of_residual,norm_max_of_residual);                    if (options->verbose) printf("Paso_Solver: Step %5d: l2/lmax-norm of residual is  %e/%e",totIter,norm2_of_residual,norm_max_of_residual);
215                    if (totIter>0 && norm2_of_residual>=last_norm2_of_residual &&  norm_max_of_residual>=last_norm_max_of_residual) {                    if (totIter>0 && norm2_of_residual>=last_norm2_of_residual &&  norm_max_of_residual>=last_norm_max_of_residual) {
216                       if (options->verbose) printf(" divergence!\n");                       if (options->verbose) printf(" divergence!\n");
217                       Paso_setError(WARNING, "No improvement during iteration. Iterative solver gives up.");                       Paso_setError(WARNING, "Paso_Solver: No improvement during iteration. Iterative solver gives up.");
218                    } else {                    } else {
219                       /* */                       /* */
220                       if (norm2_of_residual>tolerance*norm2_of_b || norm_max_of_residual>tolerance*norm_max_of_b ) {                       if (norm2_of_residual>tolerance*norm2_of_b || norm_max_of_residual>tolerance*norm_max_of_b ) {
# Line 242  void Paso_Solver(Paso_SystemMatrix* A,do Line 248  void Paso_Solver(Paso_SystemMatrix* A,do
248                          if (errorCode==NO_ERROR) {                          if (errorCode==NO_ERROR) {
249                             finalizeIteration = FALSE;                             finalizeIteration = FALSE;
250                          } else if (errorCode==SOLVER_MAXITER_REACHED) {                          } else if (errorCode==SOLVER_MAXITER_REACHED) {
251                             Paso_setError(WARNING,"maximum number of iteration step reached.\nReturned solution does not fulfil stopping criterion.");                             Paso_setError(WARNING,"Paso_Solver: maximum number of iteration step reached.\nReturned solution does not fulfil stopping criterion.");
252                             if (options->verbose) printf("Maximum number of iterations reached.!\n");                             if (options->verbose) printf("Paso_Solver: Maximum number of iterations reached.!\n");
253                          } else if (errorCode == SOLVER_INPUT_ERROR ) {                          } else if (errorCode == SOLVER_INPUT_ERROR ) {
254                             Paso_setError(SYSTEM_ERROR,"illegal dimension in iterative solver.");                             Paso_setError(SYSTEM_ERROR,"Paso_Solver: illegal dimension in iterative solver.");
255                             if (options->verbose) printf("Internal error!\n");                             if (options->verbose) printf("Paso_Solver: Internal error!\n");
256                          } else if ( errorCode == SOLVER_BREAKDOWN ) {                          } else if ( errorCode == SOLVER_BREAKDOWN ) {
257                             if (cntIter <= 1) {                             if (cntIter <= 1) {
258                                Paso_setError(ZERO_DIVISION_ERROR, "fatal break down in iterative solver.");                                Paso_setError(ZERO_DIVISION_ERROR, "Paso_Solver: fatal break down in iterative solver.");
259                                if (options->verbose) printf("Uncurable break down!\n");                                if (options->verbose) printf("Paso_Solver: Uncurable break down!\n");
260                             } else {                             } else {
261                                if (options->verbose) printf("Breakdown at iter %d (residual = %e). Restarting ...\n", cntIter+totIter, tol);                                if (options->verbose) printf("Paso_Solver: Breakdown at iter %d (residual = %e). Restarting ...\n", cntIter+totIter, tol);
262                                finalizeIteration = FALSE;                                finalizeIteration = FALSE;
263                             }                             }
264                        } else {                        } else {
# Line 264  void Paso_Solver(Paso_SystemMatrix* A,do Line 270  void Paso_Solver(Paso_SystemMatrix* A,do
270                MEMFREE(r);                MEMFREE(r);
271                time_iter=Paso_timer()-time_iter;                time_iter=Paso_timer()-time_iter;
272                if (options->verbose)  {                if (options->verbose)  {
273                   printf("\ntiming: solver: %.4e sec\n",time_iter);                   printf("\ntiming: Paso_Solver:  %.4e sec\n",time_iter);
274                   if (totIter>0) printf("timing: per iteration step: %.4e sec\n",time_iter/totIter);                   if (totIter>0) printf("timing: per iteration step: %.4e sec\n",time_iter/totIter);
275                }                }
276             }             }

Legend:
Removed from v.2155  
changed lines
  Added in v.2156

  ViewVC Help
Powered by ViewVC 1.1.26