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: }