# Annotation of /trunk/paso/src/NewtonGMRES.c

Revision 1637 - (hide annotations)
Mon Jul 14 05:34:59 2008 UTC (14 years, 4 months ago) by ksteube
File MIME type: text/plain
File size: 3895 byte(s)
```Resolved some compiler warnings
Changed blocktimer to not use strdup, intead malloc and strcpy

```
 1 gross 1476 /* \$Id:\$ */ 2 3 /******************************************************* 4 * 5 * Copyright 2008 by University of Queensland 6 * 7 * http://esscc.uq.edu.au 8 * Primary Business: Queensland, Australia 9 * Licensed under the Open Software License version 3.0 10 * http://www.opensource.org/licenses/osl-3.0.php 11 * 12 *******************************************************/ 13 14 /* 15 * Purpose 16 * ======= 17 * 18 * NewtonGMRES solves the non-linear system F(x)=0 using the 19 * restarted GMRES method. 20 * 21 * Convergence test: norm(F(x))< rtol*norm(F(x0))+atol 22 * 23 * x input initial guess 24 * output solution approximation 25 * 26 * based on code by C. T. Kelley, July 1, 1994. 27 * 28 */ 29 30 #include "Common.h" 31 #include "Solver.h" 32 artak 1480 #include "Util.h" 33 gross 1476 #ifdef _OPENMP 34 #include 35 #endif 36 err_t Paso_Solver_NewtonGMRES( 37 Paso_Function *F, 38 double *x, 39 Paso_Options* options, 40 Paso_Performance* pp) 41 42 { 43 const double gamma=.9; 44 dim_t gmres_iter; 45 double stop_tol, fnrm,fnrmo, rat,etaold,etanew, gmres_tol; 46 bool_t convergeFlag=FALSE, maxIterFlag=FALSE, breakFlag=FALSE; 47 double *f=NULL, *step=NULL; 48 err_t Status=SOLVER_NO_ERROR; 49 bool_t debug=options->verbose; 50 dim_t n = n=F->local_n; 51 dim_t itc=0; 52 double atol=options->absolute_tolerance; 53 double rtol=options->tolerance; 54 dim_t maxit=options->iter_max; 55 dim_t lmaxit=options->inner_iter_max; 56 bool_t adapt_eta=options->adapt_inner_tolerance; 57 double etamax=options->inner_tolerance; 58 double eta=etamax; 59 /* 60 * etamax = Maximum error tolerance for residual in inner 61 * iteration. The inner iteration terminates 62 * when the relative linear residual is 63 * smaller than eta*| F(x_c) |. eta is determined 64 * by the modified Eisenstat-Walker formula if etamax > 0. 65 * If etamax < 0, then eta = |etamax| for the entire 66 * iteration. 67 */ 68 69 f=TMPMEMALLOC(n,double); 70 step=TMPMEMALLOC(n,double); 71 if (f==NULL || step ==NULL) { 72 Status=SOLVER_MEMORY_ERROR; 73 } else { 74 /* 75 * initial evaluation of F 76 */ 77 Paso_FunctionCall(F,f,x); 78 fnrm=Paso_l2(n,f,F->mpi_info); 79 /* 80 * stoping criterium: 81 */ 82 stop_tol=atol + rtol*fnrm; 83 itc=0; 84 if (debug) printf("Start Jacobi-free Newton scheme:\n\ttolerance = %e\n\tstoping tolerance = %e\n",rtol,stop_tol); 85 /* 86 * main iteration loop 87 */ 88 while (! (convergeFlag || maxIterFlag || breakFlag)) { 89 /* 90 * keep track of the ratio (rat = fnrm/frnmo) of successive residual norms and 91 * the iteration counter (itc) 92 */ 93 itc++; 94 ksteube 1637 if (debug) printf("iteration step %d: norm of F =%e\n",itc,fnrm); 95 gross 1476 /* 96 * call GMRES to get increment 97 */ 98 gmres_iter=lmaxit; 99 gmres_tol=eta; 100 ksteube 1637 if (debug) printf("GMRES called with tolerance = %e\n",eta); 101 gross 1476 Status=Paso_Solver_NLGMRES(F,f,x,step,&gmres_iter,&gmres_tol,pp); 102 itc+=gmres_iter; 103 if ((Status==SOLVER_NO_ERROR) || (Status==SOLVER_MAXITER_REACHED)) { 104 gross 1511 Status=SOLVER_NO_ERROR; 105 gross 1476 /* 106 * update x: 107 */ 108 Paso_Update(n,1.,x,1.,step); 109 Paso_FunctionCall(F,f,x); 110 fnrm=Paso_l2(n,f,F->mpi_info); 111 fnrmo=fnrm; 112 rat=fnrm/fnrmo; 113 /* 114 * adjust eta 115 */ 116 if (adapt_eta) { 117 etaold=eta; 118 etanew=gamma*rat*rat; 119 if (gamma*etaold*etaold >1.) etanew=MAX(etanew,gamma*etaold*etaold); 120 eta=MAX(MIN(etanew,etamax), .5*stop_tol/fnrm); 121 } 122 } else { 123 breakFlag=TRUE; 124 } 125 maxIterFlag = (itc > maxit); 126 convergeFlag = (fnrm <= stop_tol); 127 } 128 } 129 TMPMEMFREE(f); 130 TMPMEMFREE(step); 131 return Status; 132 }