Actual source code: tr.c

  1: #include <../src/snes/impls/tr/trimpl.h>

  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;

 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;

 19:   PetscFunctionBegin;
 20:   PetscCall((*ctx->convtest)(ksp, n, rnorm, reason, ctx->convctx));
 21:   if (*reason) PetscCall(PetscInfo(snes, "Default or user provided convergence test KSP iterations=%" PetscInt_FMT ", rnorm=%g\n", n, (double)rnorm));
 22:   /* Determine norm of solution */
 23:   PetscCall(KSPBuildSolution(ksp, NULL, &x));
 24:   PetscCall(VecNorm(x, NORM_2, &nrm));
 25:   if (nrm >= neP->delta) {
 26:     PetscCall(PetscInfo(snes, "Ending linear iteration early, delta=%g, length=%g\n", (double)neP->delta, (double)nrm));
 27:     *reason = KSP_CONVERGED_STEP_LENGTH;
 28:   }
 29:   PetscFunctionReturn(PETSC_SUCCESS);
 30: }

 32: static PetscErrorCode SNESTR_KSPConverged_Destroy(void *cctx)
 33: {
 34:   SNES_TR_KSPConverged_Ctx *ctx = (SNES_TR_KSPConverged_Ctx *)cctx;

 36:   PetscFunctionBegin;
 37:   PetscCall((*ctx->convdestroy)(ctx->convctx));
 38:   PetscCall(PetscFree(ctx));
 39:   PetscFunctionReturn(PETSC_SUCCESS);
 40: }

 42: /*
 43:    SNESTR_Converged_Private -test convergence JUST for
 44:    the trust region tolerance.
 45: */
 46: static PetscErrorCode SNESTR_Converged_Private(SNES snes, PetscInt it, PetscReal xnorm, PetscReal pnorm, PetscReal fnorm, SNESConvergedReason *reason, void *dummy)
 47: {
 48:   SNES_NEWTONTR *neP = (SNES_NEWTONTR *)snes->data;

 50:   PetscFunctionBegin;
 51:   *reason = SNES_CONVERGED_ITERATING;
 52:   if (neP->delta < xnorm * snes->deltatol) {
 53:     PetscCall(PetscInfo(snes, "Converged due to trust region param %g<%g*%g\n", (double)neP->delta, (double)xnorm, (double)snes->deltatol));
 54:     *reason = SNES_DIVERGED_TR_DELTA;
 55:   } else if (snes->nfuncs >= snes->max_funcs && snes->max_funcs >= 0) {
 56:     PetscCall(PetscInfo(snes, "Exceeded maximum number of function evaluations: %" PetscInt_FMT "\n", snes->max_funcs));
 57:     *reason = SNES_DIVERGED_FUNCTION_COUNT;
 58:   }
 59:   PetscFunctionReturn(PETSC_SUCCESS);
 60: }

 62: /*@C
 63:    SNESNewtonTRSetPreCheck - Sets a user function that is called before the search step has been determined.
 64:        Allows the user a chance to change or override the decision of the line search routine.

 66:    Deprecated use `SNESNEWTONDCTRDC`

 68:    Logically Collective

 70:    Input Parameters:
 71: +  snes - the nonlinear solver object
 72: .  func - [optional] function evaluation routine, see `SNESNewtonTRPreCheck()`  for the calling sequence
 73: -  ctx  - [optional] user-defined context for private data for the function evaluation routine (may be NULL)

 75:    Level: intermediate

 77:    Note:
 78:    This function is called BEFORE the function evaluation within the `SNESNEWTONTR` solver.

 80: .seealso: `SNESNEWTONDCTRDC`, `SNESNEWTONDCTR`, `SNESNewtonTRPreCheck()`, `SNESNewtonTRGetPreCheck()`, `SNESNewtonTRSetPostCheck()`, `SNESNewtonTRGetPostCheck()`
 81: @*/
 82: PetscErrorCode SNESNewtonTRSetPreCheck(SNES snes, PetscErrorCode (*func)(SNES, Vec, Vec, PetscBool *, void *), void *ctx)
 83: {
 84:   SNES_NEWTONTR *tr = (SNES_NEWTONTR *)snes->data;

 86:   PetscFunctionBegin;
 88:   if (func) tr->precheck = func;
 89:   if (ctx) tr->precheckctx = ctx;
 90:   PetscFunctionReturn(PETSC_SUCCESS);
 91: }

 93: /*@C
 94:    SNESNewtonTRGetPreCheck - Gets the pre-check function

 96:    Deprecated use `SNESNEWTONDCTRDC`

 98:    Not collective

100:    Input Parameter:
101: .  snes - the nonlinear solver context

103:    Output Parameters:
104: +  func - [optional] function evaluation routine, see for the calling sequence `SNESNewtonTRPreCheck()`
105: -  ctx  - [optional] user-defined context for private data for the function evaluation routine (may be NULL)

107:    Level: intermediate

109: .seealso: `SNESNEWTONDCTRDC`, `SNESNEWTONDCTR`, `SNESNewtonTRSetPreCheck()`, `SNESNewtonTRPreCheck()`
110: @*/
111: PetscErrorCode SNESNewtonTRGetPreCheck(SNES snes, PetscErrorCode (**func)(SNES, Vec, Vec, PetscBool *, void *), void **ctx)
112: {
113:   SNES_NEWTONTR *tr = (SNES_NEWTONTR *)snes->data;

115:   PetscFunctionBegin;
117:   if (func) *func = tr->precheck;
118:   if (ctx) *ctx = tr->precheckctx;
119:   PetscFunctionReturn(PETSC_SUCCESS);
120: }

122: /*@C
123:    SNESNewtonTRSetPostCheck - Sets a user function that is called after the search step has been determined but before the next
124:        function evaluation. Allows the user a chance to change or override the decision of the line search routine

126:    Deprecated use `SNESNEWTONDCTRDC`

128:    Logically Collective

130:    Input Parameters:
131: +  snes - the nonlinear solver object
132: .  func - [optional] function evaluation routine, see `SNESNewtonTRPostCheck()`  for the calling sequence
133: -  ctx  - [optional] user-defined context for private data for the function evaluation routine (may be NULL)

135:    Level: intermediate

137:    Note:
138:    This function is called BEFORE the function evaluation within the `SNESNEWTONTR` solver while the function set in
139:    `SNESLineSearchSetPostCheck()` is called AFTER the function evaluation.

141: .seealso: `SNESNEWTONDCTRDC`, `SNESNEWTONDCTR`, `SNESNewtonTRPostCheck()`, `SNESNewtonTRGetPostCheck()`
142: @*/
143: PetscErrorCode SNESNewtonTRSetPostCheck(SNES snes, PetscErrorCode (*func)(SNES, Vec, Vec, Vec, PetscBool *, PetscBool *, void *), void *ctx)
144: {
145:   SNES_NEWTONTR *tr = (SNES_NEWTONTR *)snes->data;

147:   PetscFunctionBegin;
149:   if (func) tr->postcheck = func;
150:   if (ctx) tr->postcheckctx = ctx;
151:   PetscFunctionReturn(PETSC_SUCCESS);
152: }

154: /*@C
155:    SNESNewtonTRGetPostCheck - Gets the post-check function

157:    Deprecated use `SNESNEWTONDCTRDC`

159:    Not collective

161:    Input Parameter:
162: .  snes - the nonlinear solver context

164:    Output Parameters:
165: +  func - [optional] function evaluation routine, see for the calling sequence `SNESNewtonTRPostCheck()`
166: -  ctx  - [optional] user-defined context for private data for the function evaluation routine (may be NULL)

168:    Level: intermediate

170: .seealso: `SNESNEWTONDCTRDC`, `SNESNEWTONDCTR`, `SNESNewtonTRSetPostCheck()`, `SNESNewtonTRPostCheck()`
171: @*/
172: PetscErrorCode SNESNewtonTRGetPostCheck(SNES snes, PetscErrorCode (**func)(SNES, Vec, Vec, Vec, PetscBool *, PetscBool *, void *), void **ctx)
173: {
174:   SNES_NEWTONTR *tr = (SNES_NEWTONTR *)snes->data;

176:   PetscFunctionBegin;
178:   if (func) *func = tr->postcheck;
179:   if (ctx) *ctx = tr->postcheckctx;
180:   PetscFunctionReturn(PETSC_SUCCESS);
181: }

183: /*@C
184:    SNESNewtonTRPreCheck - Called before the step has been determined in `SNESNEWTONTR`

186:    Deprecated use `SNESNEWTONDCTRDC`

188:    Logically Collective

190:    Input Parameters:
191: +  snes - the solver
192: .  X - The last solution
193: -  Y - The step direction

195:    Output Parameters:
196: .  changed_Y - Indicator that the step direction Y has been changed.

198:    Level: developer

200: .seealso: `SNESNEWTONDCTRDC`, `SNESNEWTONDCTR`, `SNESNewtonTRSetPreCheck()`, `SNESNewtonTRGetPreCheck()`
201: @*/
202: static PetscErrorCode SNESNewtonTRPreCheck(SNES snes, Vec X, Vec Y, PetscBool *changed_Y)
203: {
204:   SNES_NEWTONTR *tr = (SNES_NEWTONTR *)snes->data;

206:   PetscFunctionBegin;
207:   *changed_Y = PETSC_FALSE;
208:   if (tr->precheck) {
209:     PetscCall((*tr->precheck)(snes, X, Y, changed_Y, tr->precheckctx));
211:   }
212:   PetscFunctionReturn(PETSC_SUCCESS);
213: }

215: /*@C
216:    SNESNewtonTRPostCheck - Called after the step has been determined in `SNESNEWTONTR` but before the function evaluation

218:    Deprecated use `SNESNEWTONDCTRDC`

220:    Logically Collective

222:    Input Parameters:
223: +  snes - the solver
224: .  X - The last solution
225: .  Y - The full step direction
226: -  W - The updated solution, W = X - Y

228:    Output Parameters:
229: +  changed_Y - indicator if step has been changed
230: -  changed_W - Indicator if the new candidate solution W has been changed.

232:    Note:
233:      If Y is changed then W is recomputed as X - Y

235:    Level: developer

237: .seealso: `SNESNEWTONDCTRDC`, `SNESNEWTONDCTR`, `SNESNewtonTRSetPostCheck()`, `SNESNewtonTRGetPostCheck()`
238: @*/
239: static PetscErrorCode SNESNewtonTRPostCheck(SNES snes, Vec X, Vec Y, Vec W, PetscBool *changed_Y, PetscBool *changed_W)
240: {
241:   SNES_NEWTONTR *tr = (SNES_NEWTONTR *)snes->data;

243:   PetscFunctionBegin;
244:   *changed_Y = PETSC_FALSE;
245:   *changed_W = PETSC_FALSE;
246:   if (tr->postcheck) {
247:     PetscCall((*tr->postcheck)(snes, X, Y, W, changed_Y, changed_W, tr->postcheckctx));
250:   }
251:   PetscFunctionReturn(PETSC_SUCCESS);
252: }

254: /*
255:    SNESSolve_NEWTONTR - Implements Newton's Method with a very simple trust
256:    region approach for solving systems of nonlinear equations.

258: */
259: static PetscErrorCode SNESSolve_NEWTONTR(SNES snes)
260: {
261:   SNES_NEWTONTR            *neP = (SNES_NEWTONTR *)snes->data;
262:   Vec                       X, F, Y, G, Ytmp, W;
263:   PetscInt                  maxits, i, lits;
264:   PetscReal                 rho, fnorm, gnorm, gpnorm, xnorm = 0, delta, nrm, ynorm, norm1;
265:   PetscScalar               cnorm;
266:   KSP                       ksp;
267:   SNESConvergedReason       reason   = SNES_CONVERGED_ITERATING;
268:   PetscBool                 breakout = PETSC_FALSE;
269:   SNES_TR_KSPConverged_Ctx *ctx;
270:   PetscErrorCode (*convtest)(KSP, PetscInt, PetscReal, KSPConvergedReason *, void *), (*convdestroy)(void *);
271:   void *convctx;

273:   PetscFunctionBegin;
274:   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);

276:   maxits = snes->max_its;  /* maximum number of iterations */
277:   X      = snes->vec_sol;  /* solution vector */
278:   F      = snes->vec_func; /* residual vector */
279:   Y      = snes->work[0];  /* work vectors */
280:   G      = snes->work[1];
281:   Ytmp   = snes->work[2];
282:   W      = snes->work[3];

284:   PetscCall(PetscObjectSAWsTakeAccess((PetscObject)snes));
285:   snes->iter = 0;
286:   PetscCall(PetscObjectSAWsGrantAccess((PetscObject)snes));

288:   /* Set the linear stopping criteria to use the More' trick. */
289:   PetscCall(SNESGetKSP(snes, &ksp));
290:   PetscCall(KSPGetConvergenceTest(ksp, &convtest, &convctx, &convdestroy));
291:   if (convtest != SNESTR_KSPConverged_Private) {
292:     PetscCall(PetscNew(&ctx));
293:     ctx->snes = snes;
294:     PetscCall(KSPGetAndClearConvergenceTest(ksp, &ctx->convtest, &ctx->convctx, &ctx->convdestroy));
295:     PetscCall(KSPSetConvergenceTest(ksp, SNESTR_KSPConverged_Private, ctx, SNESTR_KSPConverged_Destroy));
296:     PetscCall(PetscInfo(snes, "Using Krylov convergence test SNESTR_KSPConverged_Private\n"));
297:   }

299:   if (!snes->vec_func_init_set) {
300:     PetscCall(SNESComputeFunction(snes, X, F)); /* F(X) */
301:   } else snes->vec_func_init_set = PETSC_FALSE;

303:   PetscCall(VecNorm(F, NORM_2, &fnorm)); /* fnorm <- || F || */
304:   SNESCheckFunctionNorm(snes, fnorm);
305:   PetscCall(VecNorm(X, NORM_2, &xnorm)); /* xnorm <- || X || */
306:   PetscCall(PetscObjectSAWsTakeAccess((PetscObject)snes));
307:   snes->norm = fnorm;
308:   PetscCall(PetscObjectSAWsGrantAccess((PetscObject)snes));
309:   delta      = xnorm ? neP->delta0 * xnorm : neP->delta0;
310:   neP->delta = delta;
311:   PetscCall(SNESLogConvergenceHistory(snes, fnorm, 0));
312:   PetscCall(SNESMonitor(snes, 0, fnorm));

314:   /* test convergence */
315:   PetscUseTypeMethod(snes, converged, snes->iter, 0.0, 0.0, fnorm, &snes->reason, snes->cnvP);
316:   if (snes->reason) PetscFunctionReturn(PETSC_SUCCESS);

318:   for (i = 0; i < maxits; i++) {
319:     /* Call general purpose update function */
320:     PetscTryTypeMethod(snes, update, snes->iter);

322:     /* Solve J Y = F, where J is Jacobian matrix */
323:     PetscCall(SNESComputeJacobian(snes, X, snes->jacobian, snes->jacobian_pre));
324:     SNESCheckJacobianDomainerror(snes);
325:     PetscCall(KSPSetOperators(snes->ksp, snes->jacobian, snes->jacobian_pre));
326:     PetscCall(KSPSolve(snes->ksp, F, Ytmp));
327:     PetscCall(KSPGetIterationNumber(snes->ksp, &lits));
328:     snes->linear_its += lits;

330:     PetscCall(PetscInfo(snes, "iter=%" PetscInt_FMT ", linear solve iterations=%" PetscInt_FMT "\n", snes->iter, lits));
331:     PetscCall(VecNorm(Ytmp, NORM_2, &nrm));
332:     norm1 = nrm;

334:     while (1) {
335:       PetscBool changed_y;
336:       PetscBool changed_w;
337:       PetscCall(VecCopy(Ytmp, Y));
338:       nrm = norm1;

340:       /* Scale Y if need be and predict new value of F norm */
341:       if (nrm >= delta) {
342:         nrm    = delta / nrm;
343:         gpnorm = (1.0 - nrm) * fnorm;
344:         cnorm  = nrm;
345:         PetscCall(PetscInfo(snes, "Scaling direction by %g\n", (double)nrm));
346:         PetscCall(VecScale(Y, cnorm));
347:         nrm   = gpnorm;
348:         ynorm = delta;
349:       } else {
350:         gpnorm = 0.0;
351:         PetscCall(PetscInfo(snes, "Direction is in Trust Region\n"));
352:         ynorm = nrm;
353:       }
354:       /* PreCheck() allows for updates to Y prior to W <- X - Y */
355:       PetscCall(SNESNewtonTRPreCheck(snes, X, Y, &changed_y));
356:       PetscCall(VecWAXPY(W, -1.0, Y, X)); /* W <- X - Y */
357:       PetscCall(SNESNewtonTRPostCheck(snes, X, Y, W, &changed_y, &changed_w));
358:       if (changed_y) PetscCall(VecWAXPY(W, -1.0, Y, X));
359:       PetscCall(VecCopy(Y, snes->vec_sol_update));
360:       PetscCall(SNESComputeFunction(snes, W, G)); /*  F(X-Y) = G */
361:       PetscCall(VecNorm(G, NORM_2, &gnorm));      /* gnorm <- || g || */
362:       SNESCheckFunctionNorm(snes, gnorm);
363:       if (fnorm == gpnorm) rho = 0.0;
364:       else rho = (fnorm * fnorm - gnorm * gnorm) / (fnorm * fnorm - gpnorm * gpnorm);

366:       /* Update size of trust region */
367:       if (rho < neP->mu) delta *= neP->delta1;
368:       else if (rho < neP->eta) delta *= neP->delta2;
369:       else delta *= neP->delta3;
370:       PetscCall(PetscInfo(snes, "fnorm=%g, gnorm=%g, ynorm=%g\n", (double)fnorm, (double)gnorm, (double)ynorm));
371:       PetscCall(PetscInfo(snes, "gpred=%g, rho=%g, delta=%g\n", (double)gpnorm, (double)rho, (double)delta));

373:       neP->delta = delta;
374:       if (rho > neP->sigma) break;
375:       PetscCall(PetscInfo(snes, "Trying again in smaller region\n"));

377:       /* check to see if progress is hopeless */
378:       neP->itflag = PETSC_FALSE;
379:       PetscCall(SNESTR_Converged_Private(snes, snes->iter, xnorm, ynorm, fnorm, &reason, snes->cnvP));
380:       if (!reason) PetscUseTypeMethod(snes, converged, snes->iter, xnorm, ynorm, fnorm, &reason, snes->cnvP);
381:       if (reason == SNES_CONVERGED_SNORM_RELATIVE) reason = SNES_DIVERGED_INNER;
382:       if (reason) {
383:         /* We're not progressing, so return with the current iterate */
384:         PetscCall(SNESMonitor(snes, i + 1, fnorm));
385:         breakout = PETSC_TRUE;
386:         break;
387:       }
388:       snes->numFailures++;
389:     }
390:     if (!breakout) {
391:       /* Update function and solution vectors */
392:       fnorm = gnorm;
393:       PetscCall(VecCopy(G, F));
394:       PetscCall(VecCopy(W, X));
395:       /* Monitor convergence */
396:       PetscCall(PetscObjectSAWsTakeAccess((PetscObject)snes));
397:       snes->iter  = i + 1;
398:       snes->norm  = fnorm;
399:       snes->xnorm = xnorm;
400:       snes->ynorm = ynorm;
401:       PetscCall(PetscObjectSAWsGrantAccess((PetscObject)snes));
402:       PetscCall(SNESLogConvergenceHistory(snes, snes->norm, lits));
403:       PetscCall(SNESMonitor(snes, snes->iter, snes->norm));
404:       /* Test for convergence, xnorm = || X || */
405:       neP->itflag = PETSC_TRUE;
406:       if (snes->ops->converged != SNESConvergedSkip) PetscCall(VecNorm(X, NORM_2, &xnorm));
407:       PetscUseTypeMethod(snes, converged, snes->iter, xnorm, ynorm, fnorm, &reason, snes->cnvP);
408:       if (reason) break;
409:     } else break;
410:   }

412:   if (i == maxits) {
413:     PetscCall(PetscInfo(snes, "Maximum number of iterations has been reached: %" PetscInt_FMT "\n", maxits));
414:     if (!reason) reason = SNES_DIVERGED_MAX_IT;
415:   }
416:   PetscCall(PetscObjectSAWsTakeAccess((PetscObject)snes));
417:   snes->reason = reason;
418:   PetscCall(PetscObjectSAWsGrantAccess((PetscObject)snes));
419:   if (convtest != SNESTR_KSPConverged_Private) {
420:     PetscCall(KSPGetAndClearConvergenceTest(ksp, &ctx->convtest, &ctx->convctx, &ctx->convdestroy));
421:     PetscCall(PetscFree(ctx));
422:     PetscCall(KSPSetConvergenceTest(ksp, convtest, convctx, convdestroy));
423:   }
424:   PetscFunctionReturn(PETSC_SUCCESS);
425: }

427: static PetscErrorCode SNESSetUp_NEWTONTR(SNES snes)
428: {
429:   PetscFunctionBegin;
430:   PetscCall(SNESSetWorkVecs(snes, 4));
431:   PetscCall(SNESSetUpMatrices(snes));
432:   PetscFunctionReturn(PETSC_SUCCESS);
433: }

435: PetscErrorCode SNESReset_NEWTONTR(SNES snes)
436: {
437:   PetscFunctionBegin;
438:   PetscFunctionReturn(PETSC_SUCCESS);
439: }

441: static PetscErrorCode SNESDestroy_NEWTONTR(SNES snes)
442: {
443:   PetscFunctionBegin;
444:   PetscCall(SNESReset_NEWTONTR(snes));
445:   PetscCall(PetscFree(snes->data));
446:   PetscFunctionReturn(PETSC_SUCCESS);
447: }

449: static PetscErrorCode SNESSetFromOptions_NEWTONTR(SNES snes, PetscOptionItems *PetscOptionsObject)
450: {
451:   SNES_NEWTONTR *ctx = (SNES_NEWTONTR *)snes->data;

453:   PetscFunctionBegin;
454:   PetscOptionsHeadBegin(PetscOptionsObject, "SNES trust region options for nonlinear equations");
455:   PetscCall(PetscOptionsReal("-snes_trtol", "Trust region tolerance", "SNESSetTrustRegionTolerance", snes->deltatol, &snes->deltatol, NULL));
456:   PetscCall(PetscOptionsReal("-snes_tr_mu", "mu", "None", ctx->mu, &ctx->mu, NULL));
457:   PetscCall(PetscOptionsReal("-snes_tr_eta", "eta", "None", ctx->eta, &ctx->eta, NULL));
458:   PetscCall(PetscOptionsReal("-snes_tr_sigma", "sigma", "None", ctx->sigma, &ctx->sigma, NULL));
459:   PetscCall(PetscOptionsReal("-snes_tr_delta0", "delta0", "None", ctx->delta0, &ctx->delta0, NULL));
460:   PetscCall(PetscOptionsReal("-snes_tr_delta1", "delta1", "None", ctx->delta1, &ctx->delta1, NULL));
461:   PetscCall(PetscOptionsReal("-snes_tr_delta2", "delta2", "None", ctx->delta2, &ctx->delta2, NULL));
462:   PetscCall(PetscOptionsReal("-snes_tr_delta3", "delta3", "None", ctx->delta3, &ctx->delta3, NULL));
463:   PetscOptionsHeadEnd();
464:   PetscFunctionReturn(PETSC_SUCCESS);
465: }

467: static PetscErrorCode SNESView_NEWTONTR(SNES snes, PetscViewer viewer)
468: {
469:   SNES_NEWTONTR *tr = (SNES_NEWTONTR *)snes->data;
470:   PetscBool      iascii;

472:   PetscFunctionBegin;
473:   PetscCall(PetscObjectTypeCompare((PetscObject)viewer, PETSCVIEWERASCII, &iascii));
474:   if (iascii) {
475:     PetscCall(PetscViewerASCIIPrintf(viewer, "  Trust region tolerance %g (-snes_trtol)\n", (double)snes->deltatol));
476:     PetscCall(PetscViewerASCIIPrintf(viewer, "  mu=%g, eta=%g, sigma=%g\n", (double)tr->mu, (double)tr->eta, (double)tr->sigma));
477:     PetscCall(PetscViewerASCIIPrintf(viewer, "  delta0=%g, delta1=%g, delta2=%g, delta3=%g\n", (double)tr->delta0, (double)tr->delta1, (double)tr->delta2, (double)tr->delta3));
478:   }
479:   PetscFunctionReturn(PETSC_SUCCESS);
480: }

482: /*MC
483:       SNESNEWTONTR - Newton based nonlinear solver that uses a trust region

485:       Deprecated use `SNESNEWTONTRDC`

487:    Options Database Keys:
488: +    -snes_trtol <tol> - trust region tolerance
489: .    -snes_tr_mu <mu> - trust region parameter
490: .    -snes_tr_eta <eta> - trust region parameter
491: .    -snes_tr_sigma <sigma> - trust region parameter
492: .    -snes_tr_delta0 <delta0> -  initial size of the trust region is delta0*norm2(x)
493: .    -snes_tr_delta1 <delta1> - trust region parameter
494: .    -snes_tr_delta2 <delta2> - trust region parameter
495: -    -snes_tr_delta3 <delta3> - trust region parameter

497:    Reference:
498: .  - *  "The Minpack Project", by More', Sorensen, Garbow, Hillstrom, pages 88-111 of "Sources and Development
499:    of Mathematical Software", Wayne Cowell, editor.

501:    Level: intermediate

503: .seealso: `SNESNEWTONTRDC`, `SNESCreate()`, `SNES`, `SNESSetType()`, `SNESNEWTONLS`, `SNESSetTrustRegionTolerance()`
504: M*/
505: PETSC_EXTERN PetscErrorCode SNESCreate_NEWTONTR(SNES snes)
506: {
507:   SNES_NEWTONTR *neP;

509:   PetscFunctionBegin;
510:   snes->ops->setup          = SNESSetUp_NEWTONTR;
511:   snes->ops->solve          = SNESSolve_NEWTONTR;
512:   snes->ops->destroy        = SNESDestroy_NEWTONTR;
513:   snes->ops->setfromoptions = SNESSetFromOptions_NEWTONTR;
514:   snes->ops->view           = SNESView_NEWTONTR;
515:   snes->ops->reset          = SNESReset_NEWTONTR;

517:   snes->usesksp = PETSC_TRUE;
518:   snes->usesnpc = PETSC_FALSE;

520:   snes->alwayscomputesfinalresidual = PETSC_TRUE;

522:   PetscCall(PetscNew(&neP));
523:   snes->data  = (void *)neP;
524:   neP->mu     = 0.25;
525:   neP->eta    = 0.75;
526:   neP->delta  = 0.0;
527:   neP->delta0 = 0.2;
528:   neP->delta1 = 0.3;
529:   neP->delta2 = 0.75;
530:   neP->delta3 = 2.0;
531:   neP->sigma  = 0.0001;
532:   neP->itflag = PETSC_FALSE;
533:   neP->rnorm0 = 0.0;
534:   neP->ttol   = 0.0;
535:   PetscFunctionReturn(PETSC_SUCCESS);
536: }