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