1e7e93795SLois Curfman McInnes #ifndef lint 2*052efed2SBarry Smith static char vcid[] = "$Id: snesut.c,v 1.6 1995/11/01 23:20:55 bsmith Exp bsmith $"; 3e7e93795SLois Curfman McInnes #endif 4e7e93795SLois Curfman McInnes 5e7e93795SLois Curfman McInnes #include <math.h> 6e7e93795SLois Curfman McInnes #include "snesimpl.h" /*I "snes.h" I*/ 7e7e93795SLois Curfman McInnes 84b828684SBarry Smith /*@C 9e7e93795SLois Curfman McInnes SNESDefaultMonitor - Default SNES monitoring routine. 10e7e93795SLois Curfman McInnes 11e7e93795SLois Curfman McInnes Input Parameters: 12e7e93795SLois Curfman McInnes . snes - the SNES context 13e7e93795SLois Curfman McInnes . its - iteration number 14e7e93795SLois Curfman McInnes . fgnorm - 2-norm of residual (or gradient) 15e7e93795SLois Curfman McInnes . dummy - unused context 16e7e93795SLois Curfman McInnes 17e7e93795SLois Curfman McInnes Notes: 18e7e93795SLois Curfman McInnes For SNES_NONLINEAR_EQUATIONS methods the routine prints the 19e7e93795SLois Curfman McInnes residual norm at each iteration. 20e7e93795SLois Curfman McInnes 21e7e93795SLois Curfman McInnes For SNES_UNCONSTRAINED_MINIMIZATION methods the routine prints the 22e7e93795SLois Curfman McInnes function value and gradient norm at each iteration. 23e7e93795SLois Curfman McInnes 24e7e93795SLois Curfman McInnes .keywords: SNES, nonlinear, default, monitor, norm 25e7e93795SLois Curfman McInnes 26e7e93795SLois Curfman McInnes .seealso: SNESSetMonitor() 27e7e93795SLois Curfman McInnes @*/ 28e7e93795SLois Curfman McInnes int SNESDefaultMonitor(SNES snes,int its,double fgnorm,void *dummy) 29e7e93795SLois Curfman McInnes { 30e7e93795SLois Curfman McInnes if (snes->method_class == SNES_NONLINEAR_EQUATIONS) 31*052efed2SBarry Smith MPIU_printf(snes->comm, "iter = %d, SNES Function norm %g \n",its,fgnorm); 32e7e93795SLois Curfman McInnes else if (snes->method_class == SNES_UNCONSTRAINED_MINIMIZATION) 33e7e93795SLois Curfman McInnes MPIU_printf(snes->comm, 34e7e93795SLois Curfman McInnes "iter = %d, Function value %g, Gradient norm %g \n",its,snes->fc,fgnorm); 35e7e93795SLois Curfman McInnes else SETERRQ(1,"SNESDefaultMonitor:Unknown method class"); 36e7e93795SLois Curfman McInnes return 0; 37e7e93795SLois Curfman McInnes } 38e7e93795SLois Curfman McInnes /* ---------------------------------------------------------------- */ 39e7e93795SLois Curfman McInnes int SNESDefaultSMonitor(SNES snes,int its, double fgnorm,void *dummy) 40e7e93795SLois Curfman McInnes { 41e7e93795SLois Curfman McInnes if (snes->method_class == SNES_NONLINEAR_EQUATIONS) { 42e7e93795SLois Curfman McInnes if (fgnorm > 1.e-9 || fgnorm == 0.0) { 43e7e93795SLois Curfman McInnes MPIU_printf(snes->comm, "iter = %d, Function norm %g \n",its,fgnorm); 44e7e93795SLois Curfman McInnes } 45e7e93795SLois Curfman McInnes else if (fgnorm > 1.e-11){ 46e7e93795SLois Curfman McInnes MPIU_printf(snes->comm, "iter = %d, Function norm %5.3e \n",its,fgnorm); 47e7e93795SLois Curfman McInnes } 48e7e93795SLois Curfman McInnes else { 49e7e93795SLois Curfman McInnes MPIU_printf(snes->comm, "iter = %d, Function norm < 1.e-11\n",its); 50e7e93795SLois Curfman McInnes } 51e7e93795SLois Curfman McInnes } else if (snes->method_class == SNES_UNCONSTRAINED_MINIMIZATION) { 52e7e93795SLois Curfman McInnes if (fgnorm > 1.e-9 || fgnorm == 0.0) { 53e7e93795SLois Curfman McInnes MPIU_printf(snes->comm, 54e7e93795SLois Curfman McInnes "iter = %d, Function value %g, Gradient norm %g \n", 55e7e93795SLois Curfman McInnes its,snes->fc,fgnorm); 56e7e93795SLois Curfman McInnes } 57e7e93795SLois Curfman McInnes else if (fgnorm > 1.e-11) { 58e7e93795SLois Curfman McInnes MPIU_printf(snes->comm, 59e7e93795SLois Curfman McInnes "iter = %d, Function value %g, Gradient norm %5.3e \n", 60e7e93795SLois Curfman McInnes its,snes->fc,fgnorm); 61e7e93795SLois Curfman McInnes } 62e7e93795SLois Curfman McInnes else { 63e7e93795SLois Curfman McInnes MPIU_printf(snes->comm, 64e7e93795SLois Curfman McInnes "iter = %d, Function value %g, Gradient norm < 1.e-11\n", 65e7e93795SLois Curfman McInnes its,snes->fc); 66e7e93795SLois Curfman McInnes } 67e7e93795SLois Curfman McInnes } else SETERRQ(1,"SNESDefaultSMonitor:Unknown method class"); 68e7e93795SLois Curfman McInnes return 0; 69e7e93795SLois Curfman McInnes } 70e7e93795SLois Curfman McInnes /* ---------------------------------------------------------------- */ 714b828684SBarry Smith /*@C 72e7e93795SLois Curfman McInnes SNESDefaultConverged - Default test for monitoring the convergence 73e7e93795SLois Curfman McInnes of the solvers for systems of nonlinear equations. 74e7e93795SLois Curfman McInnes 75e7e93795SLois Curfman McInnes Input Parameters: 76e7e93795SLois Curfman McInnes . snes - the SNES context 77e7e93795SLois Curfman McInnes . xnorm - 2-norm of current iterate 78e7e93795SLois Curfman McInnes . pnorm - 2-norm of current step 79e7e93795SLois Curfman McInnes . fnorm - 2-norm of function 80e7e93795SLois Curfman McInnes . dummy - unused context 81e7e93795SLois Curfman McInnes 82e7e93795SLois Curfman McInnes Returns: 83e7e93795SLois Curfman McInnes $ 2 if ( fnorm < atol ), 84e7e93795SLois Curfman McInnes $ 3 if ( pnorm < xtol*xnorm ), 85e7e93795SLois Curfman McInnes $ -2 if ( nfct > maxf ), 86e7e93795SLois Curfman McInnes $ 0 otherwise, 87e7e93795SLois Curfman McInnes 88e7e93795SLois Curfman McInnes where 89e7e93795SLois Curfman McInnes $ maxf - maximum number of function evaluations, 90e7e93795SLois Curfman McInnes $ set with SNESSetMaxFunctionEvaluations() 91e7e93795SLois Curfman McInnes $ nfct - number of function evaluations, 92e7e93795SLois Curfman McInnes $ atol - absolute function norm tolerance, 93e7e93795SLois Curfman McInnes $ set with SNESSetAbsoluteTolerance() 94e7e93795SLois Curfman McInnes $ xtol - relative function norm tolerance, 95e7e93795SLois Curfman McInnes $ set with SNESSetRelativeTolerance() 96e7e93795SLois Curfman McInnes 97e7e93795SLois Curfman McInnes .keywords: SNES, nonlinear, default, converged, convergence 98e7e93795SLois Curfman McInnes 99e7e93795SLois Curfman McInnes .seealso: SNESSetConvergenceTest(), SNESEisenstatWalkerConverged() 100e7e93795SLois Curfman McInnes @*/ 101e7e93795SLois Curfman McInnes int SNESDefaultConverged(SNES snes,double xnorm,double pnorm,double fnorm, 102e7e93795SLois Curfman McInnes void *dummy) 103e7e93795SLois Curfman McInnes { 104e7e93795SLois Curfman McInnes if (snes->method_class != SNES_NONLINEAR_EQUATIONS) SETERRQ(1, 10548d91487SBarry Smith "SNESDefaultConverged:For SNES_NONLINEAR_EQUATIONS only"); 106e7e93795SLois Curfman McInnes /* Note: Reserve return code 1, -1 for compatibility with 107e7e93795SLois Curfman McInnes SNESTrustRegionDefaultConverged */ 108e7e93795SLois Curfman McInnes if (fnorm < snes->atol) { 109e7e93795SLois Curfman McInnes PLogInfo((PetscObject)snes, 1100de55854SLois Curfman McInnes "SNES: Converged due to function norm %g < %g\n",fnorm,snes->atol); 111e7e93795SLois Curfman McInnes return 2; 112e7e93795SLois Curfman McInnes } 113e7e93795SLois Curfman McInnes if (pnorm < snes->xtol*(xnorm)) { 114e7e93795SLois Curfman McInnes PLogInfo((PetscObject)snes, 115e7e93795SLois Curfman McInnes "SNES: Converged due to small update length: %g < %g * %g\n", 116e7e93795SLois Curfman McInnes pnorm,snes->xtol,xnorm); 117e7e93795SLois Curfman McInnes return 3; 118e7e93795SLois Curfman McInnes } 119e7e93795SLois Curfman McInnes if (snes->nfuncs > snes->max_funcs) { 120e7e93795SLois Curfman McInnes PLogInfo((PetscObject)snes, 121e7e93795SLois Curfman McInnes "SNES: Exceeded maximum number of function evaluations: %d > %d\n", 122e7e93795SLois Curfman McInnes snes->nfuncs, snes->max_funcs ); 123e7e93795SLois Curfman McInnes return -2; 124e7e93795SLois Curfman McInnes } 125e7e93795SLois Curfman McInnes return 0; 126e7e93795SLois Curfman McInnes } 127e7e93795SLois Curfman McInnes /* ------------------------------------------------------------ */ 128e7e93795SLois Curfman McInnes /*@ 129e7e93795SLois Curfman McInnes SNES_KSP_SetConvergenceTestEW - Sets alternative convergence test for 130e7e93795SLois Curfman McInnes for the linear solvers within an inexact Newton method. 131e7e93795SLois Curfman McInnes 132e7e93795SLois Curfman McInnes Input Parameter: 133e7e93795SLois Curfman McInnes . snes - SNES context 134e7e93795SLois Curfman McInnes 135e7e93795SLois Curfman McInnes Notes: 136e7e93795SLois Curfman McInnes Currently, the default is to use a constant relative tolerance for 137e7e93795SLois Curfman McInnes the inner linear solvers. Alternatively, one can use the 138e7e93795SLois Curfman McInnes Eisenstat-Walker method, where the relative convergence tolerance 139e7e93795SLois Curfman McInnes is reset at each Newton iteration according progress of the nonlinear 140e7e93795SLois Curfman McInnes solver. 141e7e93795SLois Curfman McInnes 142e7e93795SLois Curfman McInnes Reference: 143e7e93795SLois Curfman McInnes S. C. Eisenstat and H. F. Walker, "Choosing the forcing terms in an 144e7e93795SLois Curfman McInnes inexact Newton method", Utah State University Math. Stat. Dept. Res. 145e7e93795SLois Curfman McInnes Report 6/94/75, June, 1994, to appear in SIAM J. Sci. Comput. 146e7e93795SLois Curfman McInnes 147e7e93795SLois Curfman McInnes .keywords: SNES, KSP, Eisenstat, Walker, convergence, test, inexact, Newton 148e7e93795SLois Curfman McInnes @*/ 149e7e93795SLois Curfman McInnes int SNES_KSP_SetConvergenceTestEW(SNES snes) 150e7e93795SLois Curfman McInnes { 151e7e93795SLois Curfman McInnes snes->ksp_ewconv = 1; 152e7e93795SLois Curfman McInnes return 0; 153e7e93795SLois Curfman McInnes } 154e7e93795SLois Curfman McInnes 155e7e93795SLois Curfman McInnes /*@ 156e7e93795SLois Curfman McInnes SNES_KSP_SetParametersEW - Sets parameters for Eisenstat-Walker 157e7e93795SLois Curfman McInnes convergence criteria for the linear solvers within an inexact 158e7e93795SLois Curfman McInnes Newton method. 159e7e93795SLois Curfman McInnes 160e7e93795SLois Curfman McInnes Input Parameters: 161e7e93795SLois Curfman McInnes . snes - SNES context 162e7e93795SLois Curfman McInnes . version - version 1 or 2 (default is 2) 163e7e93795SLois Curfman McInnes . rtol_0 - initial relative tolerance 164e7e93795SLois Curfman McInnes $ (0 <= rtol_0 < 1) 165e7e93795SLois Curfman McInnes . rtol_max - maximum relative tolerance 166e7e93795SLois Curfman McInnes $ (0 <= rtol_max < 1) 167e7e93795SLois Curfman McInnes . alpha - power for version 2 rtol computation 168e7e93795SLois Curfman McInnes $ (1 < alpha <= 2) 169e7e93795SLois Curfman McInnes . alpha2 - power for safeguard 170e7e93795SLois Curfman McInnes . gamma2 - multiplicative factor for version 2 rtol computation 171e7e93795SLois Curfman McInnes $ (0 <= gamma2 <= 1) 172e7e93795SLois Curfman McInnes . threshold - threshold for imposing safeguard 173e7e93795SLois Curfman McInnes $ (0 < threshold < 1) 174e7e93795SLois Curfman McInnes 175e7e93795SLois Curfman McInnes Note: 176e7e93795SLois Curfman McInnes Use PETSC_DEFAULT to retain the default for any of the parameters. 177e7e93795SLois Curfman McInnes 178e7e93795SLois Curfman McInnes Reference: 179e7e93795SLois Curfman McInnes S. C. Eisenstat and H. F. Walker, "Choosing the forcing terms in an 180e7e93795SLois Curfman McInnes inexact Newton method", Utah State University Math. Stat. Dept. Res. 181e7e93795SLois Curfman McInnes Report 6/94/75, June, 1994, to appear in SIAM J. Sci. Comput. 182e7e93795SLois Curfman McInnes 183e7e93795SLois Curfman McInnes .keywords: SNES, KSP, Eisenstat, Walker, set, parameters 184e7e93795SLois Curfman McInnes 185e7e93795SLois Curfman McInnes .seealso: SNES_KSP_SetConvergenceTestEW() 186e7e93795SLois Curfman McInnes @*/ 187e7e93795SLois Curfman McInnes int SNES_KSP_SetParametersEW(SNES snes,int version,double rtol_0, 188e7e93795SLois Curfman McInnes double rtol_max,double gamma2,double alpha, 189e7e93795SLois Curfman McInnes double alpha2,double threshold) 190e7e93795SLois Curfman McInnes { 191e7e93795SLois Curfman McInnes SNES_KSP_EW_ConvCtx *kctx = (SNES_KSP_EW_ConvCtx*)snes->kspconvctx; 192e7e93795SLois Curfman McInnes if (!kctx) SETERRQ(1,"SNES_KSP_SetParametersEW:No context"); 193e7e93795SLois Curfman McInnes if (version != PETSC_DEFAULT) kctx->version = version; 194e7e93795SLois Curfman McInnes if (rtol_0 != PETSC_DEFAULT) kctx->rtol_0 = rtol_0; 195e7e93795SLois Curfman McInnes if (rtol_max != PETSC_DEFAULT) kctx->rtol_max = rtol_max; 196e7e93795SLois Curfman McInnes if (gamma2 != PETSC_DEFAULT) kctx->gamma = gamma2; 197e7e93795SLois Curfman McInnes if (alpha != PETSC_DEFAULT) kctx->alpha = alpha; 198e7e93795SLois Curfman McInnes if (alpha2 != PETSC_DEFAULT) kctx->alpha2 = alpha2; 199e7e93795SLois Curfman McInnes if (threshold != PETSC_DEFAULT) kctx->threshold = threshold; 200e7e93795SLois Curfman McInnes if (kctx->rtol_0 < 0.0 || kctx->rtol_0 >= 1.0) SETERRQ(1, 201e7e93795SLois Curfman McInnes "SNES_KSP_SetParametersEW: 0.0 <= rtol_0 < 1.0\n"); 202e7e93795SLois Curfman McInnes if (kctx->rtol_max < 0.0 || kctx->rtol_max >= 1.0) SETERRQ(1, 203e7e93795SLois Curfman McInnes "SNES_KSP_SetParametersEW: 0.0 <= rtol_max < 1.0\n"); 204e7e93795SLois Curfman McInnes if (kctx->threshold <= 0.0 || kctx->threshold >= 1.0) SETERRQ(1, 205e7e93795SLois Curfman McInnes "SNES_KSP_SetParametersEW: 0.0 < threshold < 1.0\n"); 206e7e93795SLois Curfman McInnes if (kctx->gamma < 0.0 || kctx->gamma > 1.0) SETERRQ(1, 207e7e93795SLois Curfman McInnes "SNES_KSP_SetParametersEW: 0.0 <= alpha <= 1.0\n"); 208e7e93795SLois Curfman McInnes if (kctx->alpha <= 1.0 || kctx->alpha > 2.0) SETERRQ(1, 209e7e93795SLois Curfman McInnes "SNES_KSP_SetParametersEW: 1.0 < alpha <= 2.0\n"); 210e7e93795SLois Curfman McInnes if (kctx->version != 1 && kctx->version !=2) SETERRQ(1, 211e7e93795SLois Curfman McInnes "SNES_KSP_SetParametersEW: Only versions 1 and 2 are supported"); 212e7e93795SLois Curfman McInnes return 0; 213e7e93795SLois Curfman McInnes } 214e7e93795SLois Curfman McInnes 215e7e93795SLois Curfman McInnes int SNES_KSP_EW_ComputeRelativeTolerance_Private(SNES snes,KSP ksp) 216e7e93795SLois Curfman McInnes { 217e7e93795SLois Curfman McInnes SNES_KSP_EW_ConvCtx *kctx = (SNES_KSP_EW_ConvCtx*)snes->kspconvctx; 218e7e93795SLois Curfman McInnes double rtol, stol; 219e7e93795SLois Curfman McInnes int ierr; 220e7e93795SLois Curfman McInnes if (!kctx) 221e7e93795SLois Curfman McInnes SETERRQ(1,"SNES_KSP_EW_ComputeRelativeTolerance_Private:No context"); 222e7e93795SLois Curfman McInnes if (snes->iter == 1) { 223e7e93795SLois Curfman McInnes rtol = kctx->rtol_0; 224e7e93795SLois Curfman McInnes } else { 225e7e93795SLois Curfman McInnes if (kctx->version == 1) { 226e7e93795SLois Curfman McInnes rtol = (snes->norm - kctx->lresid_last)/kctx->norm_last; 227e7e93795SLois Curfman McInnes if (rtol < 0.0) rtol = -rtol; 228e7e93795SLois Curfman McInnes stol = pow(kctx->rtol_last,kctx->alpha2); 2290452661fSBarry Smith if (stol > kctx->threshold) rtol = PetscMax(rtol,stol); 230e7e93795SLois Curfman McInnes } else if (kctx->version == 2) { 231e7e93795SLois Curfman McInnes rtol = kctx->gamma * pow(snes->norm/kctx->norm_last,kctx->alpha); 232e7e93795SLois Curfman McInnes stol = kctx->gamma * pow(kctx->rtol_last,kctx->alpha); 2330452661fSBarry Smith if (stol > kctx->threshold) rtol = PetscMax(rtol,stol); 234e7e93795SLois Curfman McInnes } else SETERRQ(1, 23548d91487SBarry Smith "SNES_KSP_EW_Converged_Private:Only versions 1 or 2 are supported"); 236e7e93795SLois Curfman McInnes } 2370452661fSBarry Smith rtol = PetscMin(rtol,kctx->rtol_max); 238e7e93795SLois Curfman McInnes kctx->rtol_last = rtol; 239e7e93795SLois Curfman McInnes PLogInfo((PetscObject)snes, 240e7e93795SLois Curfman McInnes "SNES: iter %d, Eisenstat-Walker (version %d) KSP rtol = %g\n", 241e7e93795SLois Curfman McInnes snes->iter,kctx->version,rtol); 242e7e93795SLois Curfman McInnes ierr = KSPSetTolerances(ksp,rtol,PETSC_DEFAULT,PETSC_DEFAULT,PETSC_DEFAULT); 243e7e93795SLois Curfman McInnes CHKERRQ(ierr); 244e7e93795SLois Curfman McInnes kctx->norm_last = snes->norm; 245e7e93795SLois Curfman McInnes return 0; 246e7e93795SLois Curfman McInnes } 247e7e93795SLois Curfman McInnes 248e7e93795SLois Curfman McInnes int SNES_KSP_EW_Converged_Private(KSP ksp,int n,double rnorm,void *ctx) 249e7e93795SLois Curfman McInnes { 250e7e93795SLois Curfman McInnes SNES snes = (SNES)ctx; 251e7e93795SLois Curfman McInnes SNES_KSP_EW_ConvCtx *kctx = (SNES_KSP_EW_ConvCtx*)snes->kspconvctx; 252e7e93795SLois Curfman McInnes int convinfo; 253e7e93795SLois Curfman McInnes 25448d91487SBarry Smith if (!kctx) SETERRQ(1,"SNES_KSP_EW_Converged_Private:No convergence context"); 255e7e93795SLois Curfman McInnes if (n == 0) SNES_KSP_EW_ComputeRelativeTolerance_Private(snes,ksp); 256e7e93795SLois Curfman McInnes convinfo = KSPDefaultConverged(ksp,n,rnorm,ctx); 257e7e93795SLois Curfman McInnes kctx->lresid_last = rnorm; 258e7e93795SLois Curfman McInnes if (convinfo) 259e7e93795SLois Curfman McInnes PLogInfo((PetscObject)snes,"SNES: KSP iterations=%d, rnorm=%g\n",n,rnorm); 260e7e93795SLois Curfman McInnes return convinfo; 261e7e93795SLois Curfman McInnes } 262e7e93795SLois Curfman McInnes 263e7e93795SLois Curfman McInnes 264