xref: /petsc/src/snes/impls/tr/tr.c (revision ebead697dbf761eb322f829370bbe90b3bd93fa3)
1 #include <../src/snes/impls/tr/trimpl.h>                /*I   "petscsnes.h"   I*/
2 
3 typedef struct {
4   SNES           snes;
5   /*  Information on the regular SNES convergence test; which may have been user provided */
6   PetscErrorCode (*convtest)(KSP,PetscInt,PetscReal,KSPConvergedReason*,void*);
7   PetscErrorCode (*convdestroy)(void*);
8   void           *convctx;
9 } SNES_TR_KSPConverged_Ctx;
10 
11 static PetscErrorCode SNESTR_KSPConverged_Private(KSP ksp,PetscInt n,PetscReal rnorm,KSPConvergedReason *reason,void *cctx)
12 {
13   SNES_TR_KSPConverged_Ctx *ctx = (SNES_TR_KSPConverged_Ctx*)cctx;
14   SNES                     snes = ctx->snes;
15   SNES_NEWTONTR            *neP = (SNES_NEWTONTR*)snes->data;
16   Vec                      x;
17   PetscReal                nrm;
18 
19   PetscFunctionBegin;
20   PetscCall((*ctx->convtest)(ksp,n,rnorm,reason,ctx->convctx));
21   if (*reason) {
22     PetscCall(PetscInfo(snes,"Default or user provided convergence test KSP iterations=%" PetscInt_FMT ", rnorm=%g\n",n,(double)rnorm));
23   }
24   /* Determine norm of solution */
25   PetscCall(KSPBuildSolution(ksp,NULL,&x));
26   PetscCall(VecNorm(x,NORM_2,&nrm));
27   if (nrm >= neP->delta) {
28     PetscCall(PetscInfo(snes,"Ending linear iteration early, delta=%g, length=%g\n",(double)neP->delta,(double)nrm));
29     *reason = KSP_CONVERGED_STEP_LENGTH;
30   }
31   PetscFunctionReturn(0);
32 }
33 
34 static PetscErrorCode SNESTR_KSPConverged_Destroy(void *cctx)
35 {
36   SNES_TR_KSPConverged_Ctx *ctx = (SNES_TR_KSPConverged_Ctx*)cctx;
37 
38   PetscFunctionBegin;
39   PetscCall((*ctx->convdestroy)(ctx->convctx));
40   PetscCall(PetscFree(ctx));
41   PetscFunctionReturn(0);
42 }
43 
44 /* ---------------------------------------------------------------- */
45 /*
46    SNESTR_Converged_Private -test convergence JUST for
47    the trust region tolerance.
48 
49 */
50 static PetscErrorCode SNESTR_Converged_Private(SNES snes,PetscInt it,PetscReal xnorm,PetscReal pnorm,PetscReal fnorm,SNESConvergedReason *reason,void *dummy)
51 {
52   SNES_NEWTONTR  *neP = (SNES_NEWTONTR*)snes->data;
53 
54   PetscFunctionBegin;
55   *reason = SNES_CONVERGED_ITERATING;
56   if (neP->delta < xnorm * snes->deltatol) {
57     PetscCall(PetscInfo(snes,"Converged due to trust region param %g<%g*%g\n",(double)neP->delta,(double)xnorm,(double)snes->deltatol));
58     *reason = SNES_DIVERGED_TR_DELTA;
59   } else if (snes->nfuncs >= snes->max_funcs && snes->max_funcs >= 0) {
60     PetscCall(PetscInfo(snes,"Exceeded maximum number of function evaluations: %" PetscInt_FMT "\n",snes->max_funcs));
61     *reason = SNES_DIVERGED_FUNCTION_COUNT;
62   }
63   PetscFunctionReturn(0);
64 }
65 
66 /*@C
67    SNESNewtonTRSetPreCheck - Sets a user function that is called before the search step has been determined.
68        Allows the user a chance to change or override the decision of the line search routine.
69 
70    Logically Collective on snes
71 
72    Input Parameters:
73 +  snes - the nonlinear solver object
74 .  func - [optional] function evaluation routine, see SNESNewtonTRPreCheck()  for the calling sequence
75 -  ctx  - [optional] user-defined context for private data for the function evaluation routine (may be NULL)
76 
77    Level: intermediate
78 
79    Note: This function is called BEFORE the function evaluation within the SNESNEWTONTR solver.
80 
81 .seealso: `SNESNewtonTRPreCheck()`, `SNESNewtonTRGetPreCheck()`, `SNESNewtonTRSetPostCheck()`, `SNESNewtonTRGetPostCheck()`
82 @*/
83 PetscErrorCode  SNESNewtonTRSetPreCheck(SNES snes, PetscErrorCode (*func)(SNES,Vec,Vec,PetscBool*,void*),void *ctx)
84 {
85   SNES_NEWTONTR  *tr = (SNES_NEWTONTR*)snes->data;
86 
87   PetscFunctionBegin;
88   PetscValidHeaderSpecific(snes,SNES_CLASSID,1);
89   if (func) tr->precheck    = func;
90   if (ctx)  tr->precheckctx = ctx;
91   PetscFunctionReturn(0);
92 }
93 
94 /*@C
95    SNESNewtonTRGetPreCheck - Gets the pre-check function
96 
97    Not collective
98 
99    Input Parameter:
100 .  snes - the nonlinear solver context
101 
102    Output Parameters:
103 +  func - [optional] function evaluation routine, see for the calling sequence SNESNewtonTRPreCheck()
104 -  ctx  - [optional] user-defined context for private data for the function evaluation routine (may be NULL)
105 
106    Level: intermediate
107 
108 .seealso: `SNESNewtonTRSetPreCheck()`, `SNESNewtonTRPreCheck()`
109 @*/
110 PetscErrorCode  SNESNewtonTRGetPreCheck(SNES snes, PetscErrorCode (**func)(SNES,Vec,Vec,PetscBool*,void*),void **ctx)
111 {
112   SNES_NEWTONTR  *tr = (SNES_NEWTONTR*)snes->data;
113 
114   PetscFunctionBegin;
115   PetscValidHeaderSpecific(snes,SNES_CLASSID,1);
116   if (func) *func = tr->precheck;
117   if (ctx)  *ctx  = tr->precheckctx;
118   PetscFunctionReturn(0);
119 }
120 
121 /*@C
122    SNESNewtonTRSetPostCheck - Sets a user function that is called after the search step has been determined but before the next
123        function evaluation. Allows the user a chance to change or override the decision of the line search routine
124 
125    Logically Collective on snes
126 
127    Input Parameters:
128 +  snes - the nonlinear solver object
129 .  func - [optional] function evaluation routine, see SNESNewtonTRPostCheck()  for the calling sequence
130 -  ctx  - [optional] user-defined context for private data for the function evaluation routine (may be NULL)
131 
132    Level: intermediate
133 
134    Note: This function is called BEFORE the function evaluation within the SNESNEWTONTR solver while the function set in
135    SNESLineSearchSetPostCheck() is called AFTER the function evaluation.
136 
137 .seealso: `SNESNewtonTRPostCheck()`, `SNESNewtonTRGetPostCheck()`
138 @*/
139 PetscErrorCode  SNESNewtonTRSetPostCheck(SNES snes,PetscErrorCode (*func)(SNES,Vec,Vec,Vec,PetscBool*,PetscBool*,void*),void *ctx)
140 {
141   SNES_NEWTONTR  *tr = (SNES_NEWTONTR*)snes->data;
142 
143   PetscFunctionBegin;
144   PetscValidHeaderSpecific(snes,SNES_CLASSID,1);
145   if (func) tr->postcheck    = func;
146   if (ctx)  tr->postcheckctx = ctx;
147   PetscFunctionReturn(0);
148 }
149 
150 /*@C
151    SNESNewtonTRGetPostCheck - Gets the post-check function
152 
153    Not collective
154 
155    Input Parameter:
156 .  snes - the nonlinear solver context
157 
158    Output Parameters:
159 +  func - [optional] function evaluation routine, see for the calling sequence SNESNewtonTRPostCheck()
160 -  ctx  - [optional] user-defined context for private data for the function evaluation routine (may be NULL)
161 
162    Level: intermediate
163 
164 .seealso: `SNESNewtonTRSetPostCheck()`, `SNESNewtonTRPostCheck()`
165 @*/
166 PetscErrorCode  SNESNewtonTRGetPostCheck(SNES snes,PetscErrorCode (**func)(SNES,Vec,Vec,Vec,PetscBool*,PetscBool*,void*),void **ctx)
167 {
168   SNES_NEWTONTR  *tr = (SNES_NEWTONTR*)snes->data;
169 
170   PetscFunctionBegin;
171   PetscValidHeaderSpecific(snes,SNES_CLASSID,1);
172   if (func) *func = tr->postcheck;
173   if (ctx)  *ctx  = tr->postcheckctx;
174   PetscFunctionReturn(0);
175 }
176 
177 /*@C
178    SNESNewtonTRPreCheck - Called before the step has been determined in SNESNEWTONTR
179 
180    Logically Collective on snes
181 
182    Input Parameters:
183 +  snes - the solver
184 .  X - The last solution
185 -  Y - The step direction
186 
187    Output Parameters:
188 .  changed_Y - Indicator that the step direction Y has been changed.
189 
190    Level: developer
191 
192 .seealso: `SNESNewtonTRSetPreCheck()`, `SNESNewtonTRGetPreCheck()`
193 @*/
194 static PetscErrorCode SNESNewtonTRPreCheck(SNES snes,Vec X,Vec Y,PetscBool *changed_Y)
195 {
196   SNES_NEWTONTR  *tr = (SNES_NEWTONTR*)snes->data;
197 
198   PetscFunctionBegin;
199   *changed_Y = PETSC_FALSE;
200   if (tr->precheck) {
201     PetscCall((*tr->precheck)(snes,X,Y,changed_Y,tr->precheckctx));
202     PetscValidLogicalCollectiveBool(snes,*changed_Y,4);
203   }
204   PetscFunctionReturn(0);
205 }
206 
207 /*@C
208    SNESNewtonTRPostCheck - Called after the step has been determined in SNESNEWTONTR but before the function evaluation
209 
210    Logically Collective on snes
211 
212    Input Parameters:
213 +  snes - the solver
214 .  X - The last solution
215 .  Y - The full step direction
216 -  W - The updated solution, W = X - Y
217 
218    Output Parameters:
219 +  changed_Y - indicator if step has been changed
220 -  changed_W - Indicator if the new candidate solution W has been changed.
221 
222    Notes:
223      If Y is changed then W is recomputed as X - Y
224 
225    Level: developer
226 
227 .seealso: `SNESNewtonTRSetPostCheck()`, `SNESNewtonTRGetPostCheck()`
228 @*/
229 static PetscErrorCode SNESNewtonTRPostCheck(SNES snes,Vec X,Vec Y,Vec W,PetscBool *changed_Y,PetscBool *changed_W)
230 {
231   SNES_NEWTONTR  *tr = (SNES_NEWTONTR*)snes->data;
232 
233   PetscFunctionBegin;
234   *changed_Y = PETSC_FALSE;
235   *changed_W = PETSC_FALSE;
236   if (tr->postcheck) {
237     PetscCall((*tr->postcheck)(snes,X,Y,W,changed_Y,changed_W,tr->postcheckctx));
238     PetscValidLogicalCollectiveBool(snes,*changed_Y,5);
239     PetscValidLogicalCollectiveBool(snes,*changed_W,6);
240   }
241   PetscFunctionReturn(0);
242 }
243 
244 /*
245    SNESSolve_NEWTONTR - Implements Newton's Method with a very simple trust
246    region approach for solving systems of nonlinear equations.
247 
248 */
249 static PetscErrorCode SNESSolve_NEWTONTR(SNES snes)
250 {
251   SNES_NEWTONTR            *neP = (SNES_NEWTONTR*)snes->data;
252   Vec                      X,F,Y,G,Ytmp,W;
253   PetscInt                 maxits,i,lits;
254   PetscReal                rho,fnorm,gnorm,gpnorm,xnorm=0,delta,nrm,ynorm,norm1;
255   PetscScalar              cnorm;
256   KSP                      ksp;
257   SNESConvergedReason      reason = SNES_CONVERGED_ITERATING;
258   PetscBool                breakout = PETSC_FALSE;
259   SNES_TR_KSPConverged_Ctx *ctx;
260   PetscErrorCode           (*convtest)(KSP,PetscInt,PetscReal,KSPConvergedReason*,void*),(*convdestroy)(void*);
261   void                     *convctx;
262 
263   PetscFunctionBegin;
264   PetscCheck(!snes->xl && !snes->xu && !snes->ops->computevariablebounds,PetscObjectComm((PetscObject)snes),PETSC_ERR_ARG_WRONGSTATE, "SNES solver %s does not support bounds", ((PetscObject)snes)->type_name);
265 
266   maxits = snes->max_its;               /* maximum number of iterations */
267   X      = snes->vec_sol;               /* solution vector */
268   F      = snes->vec_func;              /* residual vector */
269   Y      = snes->work[0];               /* work vectors */
270   G      = snes->work[1];
271   Ytmp   = snes->work[2];
272   W      = snes->work[3];
273 
274   PetscCall(PetscObjectSAWsTakeAccess((PetscObject)snes));
275   snes->iter = 0;
276   PetscCall(PetscObjectSAWsGrantAccess((PetscObject)snes));
277 
278   /* Set the linear stopping criteria to use the More' trick. */
279   PetscCall(SNESGetKSP(snes,&ksp));
280   PetscCall(KSPGetConvergenceTest(ksp,&convtest,&convctx,&convdestroy));
281   if (convtest != SNESTR_KSPConverged_Private) {
282     PetscCall(PetscNew(&ctx));
283     ctx->snes             = snes;
284     PetscCall(KSPGetAndClearConvergenceTest(ksp,&ctx->convtest,&ctx->convctx,&ctx->convdestroy));
285     PetscCall(KSPSetConvergenceTest(ksp,SNESTR_KSPConverged_Private,ctx,SNESTR_KSPConverged_Destroy));
286     PetscCall(PetscInfo(snes,"Using Krylov convergence test SNESTR_KSPConverged_Private\n"));
287   }
288 
289   if (!snes->vec_func_init_set) {
290     PetscCall(SNESComputeFunction(snes,X,F));          /* F(X) */
291   } else snes->vec_func_init_set = PETSC_FALSE;
292 
293   PetscCall(VecNorm(F,NORM_2,&fnorm));             /* fnorm <- || F || */
294   SNESCheckFunctionNorm(snes,fnorm);
295   PetscCall(VecNorm(X,NORM_2,&xnorm));             /* xnorm <- || X || */
296   PetscCall(PetscObjectSAWsTakeAccess((PetscObject)snes));
297   snes->norm = fnorm;
298   PetscCall(PetscObjectSAWsGrantAccess((PetscObject)snes));
299   delta      = xnorm ? neP->delta0*xnorm : neP->delta0;
300   neP->delta = delta;
301   PetscCall(SNESLogConvergenceHistory(snes,fnorm,0));
302   PetscCall(SNESMonitor(snes,0,fnorm));
303 
304   /* test convergence */
305   PetscCall((*snes->ops->converged)(snes,snes->iter,0.0,0.0,fnorm,&snes->reason,snes->cnvP));
306   if (snes->reason) PetscFunctionReturn(0);
307 
308   for (i=0; i<maxits; i++) {
309 
310     /* Call general purpose update function */
311     if (snes->ops->update) PetscCall((*snes->ops->update)(snes, snes->iter));
312 
313     /* Solve J Y = F, where J is Jacobian matrix */
314     PetscCall(SNESComputeJacobian(snes,X,snes->jacobian,snes->jacobian_pre));
315     SNESCheckJacobianDomainerror(snes);
316     PetscCall(KSPSetOperators(snes->ksp,snes->jacobian,snes->jacobian_pre));
317     PetscCall(KSPSolve(snes->ksp,F,Ytmp));
318     PetscCall(KSPGetIterationNumber(snes->ksp,&lits));
319     snes->linear_its += lits;
320 
321     PetscCall(PetscInfo(snes,"iter=%" PetscInt_FMT ", linear solve iterations=%" PetscInt_FMT "\n",snes->iter,lits));
322     PetscCall(VecNorm(Ytmp,NORM_2,&nrm));
323     norm1 = nrm;
324 
325     while (1) {
326       PetscBool changed_y;
327       PetscBool changed_w;
328       PetscCall(VecCopy(Ytmp,Y));
329       nrm  = norm1;
330 
331       /* Scale Y if need be and predict new value of F norm */
332       if (nrm >= delta) {
333         nrm    = delta/nrm;
334         gpnorm = (1.0 - nrm)*fnorm;
335         cnorm  = nrm;
336         PetscCall(PetscInfo(snes,"Scaling direction by %g\n",(double)nrm));
337         PetscCall(VecScale(Y,cnorm));
338         nrm    = gpnorm;
339         ynorm  = delta;
340       } else {
341         gpnorm = 0.0;
342         PetscCall(PetscInfo(snes,"Direction is in Trust Region\n"));
343         ynorm  = nrm;
344       }
345       /* PreCheck() allows for updates to Y prior to W <- X - Y */
346       PetscCall(SNESNewtonTRPreCheck(snes,X,Y,&changed_y));
347       PetscCall(VecWAXPY(W,-1.0,Y,X));         /* W <- X - Y */
348       PetscCall(SNESNewtonTRPostCheck(snes,X,Y,W,&changed_y,&changed_w));
349       if (changed_y) PetscCall(VecWAXPY(W,-1.0,Y,X));
350       PetscCall(VecCopy(Y,snes->vec_sol_update));
351       PetscCall(SNESComputeFunction(snes,W,G)); /*  F(X-Y) = G */
352       PetscCall(VecNorm(G,NORM_2,&gnorm));      /* gnorm <- || g || */
353       SNESCheckFunctionNorm(snes,gnorm);
354       if (fnorm == gpnorm) rho = 0.0;
355       else rho = (fnorm*fnorm - gnorm*gnorm)/(fnorm*fnorm - gpnorm*gpnorm);
356 
357       /* Update size of trust region */
358       if      (rho < neP->mu)  delta *= neP->delta1;
359       else if (rho < neP->eta) delta *= neP->delta2;
360       else                     delta *= neP->delta3;
361       PetscCall(PetscInfo(snes,"fnorm=%g, gnorm=%g, ynorm=%g\n",(double)fnorm,(double)gnorm,(double)ynorm));
362       PetscCall(PetscInfo(snes,"gpred=%g, rho=%g, delta=%g\n",(double)gpnorm,(double)rho,(double)delta));
363 
364       neP->delta = delta;
365       if (rho > neP->sigma) break;
366       PetscCall(PetscInfo(snes,"Trying again in smaller region\n"));
367 
368       /* check to see if progress is hopeless */
369       neP->itflag = PETSC_FALSE;
370       PetscCall(SNESTR_Converged_Private(snes,snes->iter,xnorm,ynorm,fnorm,&reason,snes->cnvP));
371       if (!reason) PetscCall((*snes->ops->converged)(snes,snes->iter,xnorm,ynorm,fnorm,&reason,snes->cnvP));
372       if (reason == SNES_CONVERGED_SNORM_RELATIVE) reason = SNES_DIVERGED_INNER;
373       if (reason) {
374         /* We're not progressing, so return with the current iterate */
375         PetscCall(SNESMonitor(snes,i+1,fnorm));
376         breakout = PETSC_TRUE;
377         break;
378       }
379       snes->numFailures++;
380     }
381     if (!breakout) {
382       /* Update function and solution vectors */
383       fnorm = gnorm;
384       PetscCall(VecCopy(G,F));
385       PetscCall(VecCopy(W,X));
386       /* Monitor convergence */
387       PetscCall(PetscObjectSAWsTakeAccess((PetscObject)snes));
388       snes->iter = i+1;
389       snes->norm = fnorm;
390       snes->xnorm = xnorm;
391       snes->ynorm = ynorm;
392       PetscCall(PetscObjectSAWsGrantAccess((PetscObject)snes));
393       PetscCall(SNESLogConvergenceHistory(snes,snes->norm,lits));
394       PetscCall(SNESMonitor(snes,snes->iter,snes->norm));
395       /* Test for convergence, xnorm = || X || */
396       neP->itflag = PETSC_TRUE;
397       if (snes->ops->converged != SNESConvergedSkip) PetscCall(VecNorm(X,NORM_2,&xnorm));
398       PetscCall((*snes->ops->converged)(snes,snes->iter,xnorm,ynorm,fnorm,&reason,snes->cnvP));
399       if (reason) break;
400     } else break;
401   }
402 
403   if (i == maxits) {
404     PetscCall(PetscInfo(snes,"Maximum number of iterations has been reached: %" PetscInt_FMT "\n",maxits));
405     if (!reason) reason = SNES_DIVERGED_MAX_IT;
406   }
407   PetscCall(PetscObjectSAWsTakeAccess((PetscObject)snes));
408   snes->reason = reason;
409   PetscCall(PetscObjectSAWsGrantAccess((PetscObject)snes));
410   if (convtest != SNESTR_KSPConverged_Private) {
411     PetscCall(KSPGetAndClearConvergenceTest(ksp,&ctx->convtest,&ctx->convctx,&ctx->convdestroy));
412     PetscCall(PetscFree(ctx));
413     PetscCall(KSPSetConvergenceTest(ksp,convtest,convctx,convdestroy));
414   }
415   PetscFunctionReturn(0);
416 }
417 
418 /*------------------------------------------------------------*/
419 static PetscErrorCode SNESSetUp_NEWTONTR(SNES snes)
420 {
421   PetscFunctionBegin;
422   PetscCall(SNESSetWorkVecs(snes,4));
423   PetscCall(SNESSetUpMatrices(snes));
424   PetscFunctionReturn(0);
425 }
426 
427 PetscErrorCode SNESReset_NEWTONTR(SNES snes)
428 {
429   PetscFunctionBegin;
430   PetscFunctionReturn(0);
431 }
432 
433 static PetscErrorCode SNESDestroy_NEWTONTR(SNES snes)
434 {
435   PetscFunctionBegin;
436   PetscCall(SNESReset_NEWTONTR(snes));
437   PetscCall(PetscFree(snes->data));
438   PetscFunctionReturn(0);
439 }
440 /*------------------------------------------------------------*/
441 
442 static PetscErrorCode SNESSetFromOptions_NEWTONTR(PetscOptionItems *PetscOptionsObject,SNES snes)
443 {
444   SNES_NEWTONTR  *ctx = (SNES_NEWTONTR*)snes->data;
445 
446   PetscFunctionBegin;
447   PetscOptionsHeadBegin(PetscOptionsObject,"SNES trust region options for nonlinear equations");
448   PetscCall(PetscOptionsReal("-snes_trtol","Trust region tolerance","SNESSetTrustRegionTolerance",snes->deltatol,&snes->deltatol,NULL));
449   PetscCall(PetscOptionsReal("-snes_tr_mu","mu","None",ctx->mu,&ctx->mu,NULL));
450   PetscCall(PetscOptionsReal("-snes_tr_eta","eta","None",ctx->eta,&ctx->eta,NULL));
451   PetscCall(PetscOptionsReal("-snes_tr_sigma","sigma","None",ctx->sigma,&ctx->sigma,NULL));
452   PetscCall(PetscOptionsReal("-snes_tr_delta0","delta0","None",ctx->delta0,&ctx->delta0,NULL));
453   PetscCall(PetscOptionsReal("-snes_tr_delta1","delta1","None",ctx->delta1,&ctx->delta1,NULL));
454   PetscCall(PetscOptionsReal("-snes_tr_delta2","delta2","None",ctx->delta2,&ctx->delta2,NULL));
455   PetscCall(PetscOptionsReal("-snes_tr_delta3","delta3","None",ctx->delta3,&ctx->delta3,NULL));
456   PetscOptionsHeadEnd();
457   PetscFunctionReturn(0);
458 }
459 
460 static PetscErrorCode SNESView_NEWTONTR(SNES snes,PetscViewer viewer)
461 {
462   SNES_NEWTONTR  *tr = (SNES_NEWTONTR*)snes->data;
463   PetscBool      iascii;
464 
465   PetscFunctionBegin;
466   PetscCall(PetscObjectTypeCompare((PetscObject)viewer,PETSCVIEWERASCII,&iascii));
467   if (iascii) {
468     PetscCall(PetscViewerASCIIPrintf(viewer,"  Trust region tolerance %g (-snes_trtol)\n",(double)snes->deltatol));
469     PetscCall(PetscViewerASCIIPrintf(viewer,"  mu=%g, eta=%g, sigma=%g\n",(double)tr->mu,(double)tr->eta,(double)tr->sigma));
470     PetscCall(PetscViewerASCIIPrintf(viewer,"  delta0=%g, delta1=%g, delta2=%g, delta3=%g\n",(double)tr->delta0,(double)tr->delta1,(double)tr->delta2,(double)tr->delta3));
471   }
472   PetscFunctionReturn(0);
473 }
474 /* ------------------------------------------------------------ */
475 /*MC
476       SNESNEWTONTR - Newton based nonlinear solver that uses a trust region
477 
478    Options Database:
479 +    -snes_trtol <tol> - trust region tolerance
480 .    -snes_tr_mu <mu> - trust region parameter
481 .    -snes_tr_eta <eta> - trust region parameter
482 .    -snes_tr_sigma <sigma> - trust region parameter
483 .    -snes_tr_delta0 <delta0> -  initial size of the trust region is delta0*norm2(x)
484 .    -snes_tr_delta1 <delta1> - trust region parameter
485 .    -snes_tr_delta2 <delta2> - trust region parameter
486 -    -snes_tr_delta3 <delta3> - trust region parameter
487 
488    The basic algorithm is taken from "The Minpack Project", by More',
489    Sorensen, Garbow, Hillstrom, pages 88-111 of "Sources and Development
490    of Mathematical Software", Wayne Cowell, editor.
491 
492    Level: intermediate
493 
494 .seealso: `SNESCreate()`, `SNES`, `SNESSetType()`, `SNESNEWTONLS`, `SNESSetTrustRegionTolerance()`
495 
496 M*/
497 PETSC_EXTERN PetscErrorCode SNESCreate_NEWTONTR(SNES snes)
498 {
499   SNES_NEWTONTR  *neP;
500 
501   PetscFunctionBegin;
502   snes->ops->setup          = SNESSetUp_NEWTONTR;
503   snes->ops->solve          = SNESSolve_NEWTONTR;
504   snes->ops->destroy        = SNESDestroy_NEWTONTR;
505   snes->ops->setfromoptions = SNESSetFromOptions_NEWTONTR;
506   snes->ops->view           = SNESView_NEWTONTR;
507   snes->ops->reset          = SNESReset_NEWTONTR;
508 
509   snes->usesksp = PETSC_TRUE;
510   snes->usesnpc = PETSC_FALSE;
511 
512   snes->alwayscomputesfinalresidual = PETSC_TRUE;
513 
514   PetscCall(PetscNewLog(snes,&neP));
515   snes->data  = (void*)neP;
516   neP->mu     = 0.25;
517   neP->eta    = 0.75;
518   neP->delta  = 0.0;
519   neP->delta0 = 0.2;
520   neP->delta1 = 0.3;
521   neP->delta2 = 0.75;
522   neP->delta3 = 2.0;
523   neP->sigma  = 0.0001;
524   neP->itflag = PETSC_FALSE;
525   neP->rnorm0 = 0.0;
526   neP->ttol   = 0.0;
527   PetscFunctionReturn(0);
528 }
529