Actual source code: qn.c
1: #include <petsc/private/snesimpl.h>
2: #include <petscdm.h>
4: #define H(i, j) qn->dXdFmat[i * qn->m + j]
6: const char *const SNESQNScaleTypes[] = {"DEFAULT", "NONE", "SCALAR", "DIAGONAL", "JACOBIAN", "SNESQNScaleType", "SNES_QN_SCALING_", NULL};
7: const char *const SNESQNRestartTypes[] = {"DEFAULT", "NONE", "POWELL", "PERIODIC", "SNESQNRestartType", "SNES_QN_RESTART_", NULL};
8: const char *const SNESQNTypes[] = {"LBFGS", "BROYDEN", "BADBROYDEN", "SNESQNType", "SNES_QN_", NULL};
10: typedef struct {
11: Mat B; /* Quasi-Newton approximation Matrix (MATLMVM) */
12: PetscInt m; /* The number of kept previous steps */
13: PetscReal *lambda; /* The line search history of the method */
14: PetscBool monflg;
15: PetscViewer monitor;
16: PetscReal powell_gamma; /* Powell angle restart condition */
17: PetscReal scaling; /* scaling of H0 */
18: SNESQNType type; /* the type of quasi-newton method used */
19: SNESQNScaleType scale_type; /* the type of scaling used */
20: SNESQNRestartType restart_type; /* determine the frequency and type of restart conditions */
21: } SNES_QN;
23: static PetscErrorCode SNESQNGetMatrix_Private(SNES snes, Mat *B)
24: {
25: SNES_QN *qn = (SNES_QN *)snes->data;
26: const char *optionsprefix;
28: PetscFunctionBegin;
29: if (!qn->B) {
30: PetscCall(MatCreate(PetscObjectComm((PetscObject)snes), &qn->B));
31: PetscCall(SNESGetOptionsPrefix(snes, &optionsprefix));
32: PetscCall(MatSetOptionsPrefix(qn->B, optionsprefix));
33: PetscCall(MatAppendOptionsPrefix(qn->B, "qn_"));
34: switch (qn->type) {
35: case SNES_QN_BROYDEN:
36: PetscCall(MatSetType(qn->B, MATLMVMBROYDEN));
37: qn->scale_type = SNES_QN_SCALE_NONE;
38: break;
39: case SNES_QN_BADBROYDEN:
40: PetscCall(MatSetType(qn->B, MATLMVMBADBROYDEN));
41: qn->scale_type = SNES_QN_SCALE_NONE;
42: break;
43: default:
44: PetscCall(MatSetType(qn->B, MATLMVMBFGS));
45: switch (qn->scale_type) {
46: case SNES_QN_SCALE_NONE:
47: PetscCall(MatLMVMSymBroydenSetScaleType(qn->B, MAT_LMVM_SYMBROYDEN_SCALE_NONE));
48: break;
49: case SNES_QN_SCALE_SCALAR:
50: case SNES_QN_SCALE_DEFAULT:
51: PetscCall(MatLMVMSymBroydenSetScaleType(qn->B, MAT_LMVM_SYMBROYDEN_SCALE_SCALAR));
52: break;
53: case SNES_QN_SCALE_JACOBIAN:
54: PetscCall(MatLMVMSymBroydenSetScaleType(qn->B, MAT_LMVM_SYMBROYDEN_SCALE_USER));
55: break;
56: case SNES_QN_SCALE_DIAGONAL:
57: default:
58: break;
59: }
60: break;
61: }
62: }
63: *B = qn->B;
64: PetscFunctionReturn(PETSC_SUCCESS);
65: }
67: static PetscErrorCode SNESSolve_QN(SNES snes)
68: {
69: SNES_QN *qn = (SNES_QN *)snes->data;
70: Vec X, F, W;
71: Vec Y, D, Dold;
72: PetscInt i, i_r;
73: PetscReal fnorm, xnorm, ynorm;
74: SNESLineSearchReason lssucceed;
75: PetscBool badstep, powell, periodic, restart;
76: PetscScalar DolddotD, DolddotDold;
77: SNESConvergedReason reason;
78: #if defined(PETSC_USE_INFO)
79: PetscReal gnorm;
80: #endif
82: PetscFunctionBegin;
83: 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);
85: PetscCall(PetscCitationsRegister(SNESCitation, &SNEScite));
86: F = snes->vec_func; /* residual vector */
87: Y = snes->vec_sol_update; /* search direction generated by J^-1D*/
88: X = snes->vec_sol; /* solution vector */
90: /* directions */
91: W = snes->work[0];
92: D = snes->work[1];
93: Dold = snes->work[2];
95: snes->reason = SNES_CONVERGED_ITERATING;
97: PetscCall(PetscObjectSAWsTakeAccess((PetscObject)snes));
98: snes->iter = 0;
99: snes->norm = 0.;
100: PetscCall(PetscObjectSAWsGrantAccess((PetscObject)snes));
102: if (snes->npc && snes->npcside == PC_LEFT && snes->functype == SNES_FUNCTION_PRECONDITIONED) {
103: /* we need to sample the preconditioned function only once here,
104: since linesearch will always use the preconditioned function */
105: PetscCall(SNESApplyNPC(snes, X, NULL, F));
106: PetscCall(SNESGetConvergedReason(snes->npc, &reason));
107: if (reason < 0 && reason != SNES_DIVERGED_MAX_IT && reason != SNES_DIVERGED_TR_DELTA) {
108: snes->reason = SNES_DIVERGED_INNER;
109: PetscFunctionReturn(PETSC_SUCCESS);
110: }
111: } else {
112: if (!snes->vec_func_init_set) PetscCall(SNESComputeFunction(snes, X, F));
113: else snes->vec_func_init_set = PETSC_FALSE;
114: }
115: PetscCall(VecNorm(F, NORM_2, &fnorm));
116: SNESCheckFunctionNorm(snes, fnorm);
118: PetscCall(PetscObjectSAWsTakeAccess((PetscObject)snes));
119: snes->norm = fnorm;
120: PetscCall(PetscObjectSAWsGrantAccess((PetscObject)snes));
121: PetscCall(SNESLogConvergenceHistory(snes, fnorm, 0));
122: PetscCall(SNESMonitor(snes, 0, fnorm));
124: /* test convergence */
125: PetscCall(SNESConverged(snes, 0, 0.0, 0.0, fnorm));
126: if (snes->reason) PetscFunctionReturn(PETSC_SUCCESS);
128: restart = PETSC_TRUE;
129: for (i = 0, i_r = 0; i < snes->max_its; i++, i_r++) {
130: PetscScalar gdx;
132: /* general purpose update */
133: PetscTryTypeMethod(snes, update, snes->iter);
135: if (snes->npc) {
136: if (snes->npcside == PC_LEFT && snes->functype == SNES_FUNCTION_UNPRECONDITIONED) {
137: PetscCall(SNESApplyNPC(snes, X, F, D));
138: PetscCall(SNESGetConvergedReason(snes->npc, &reason));
139: if (reason < 0 && reason != SNES_DIVERGED_MAX_IT && reason != SNES_DIVERGED_TR_DELTA) {
140: snes->reason = SNES_DIVERGED_INNER;
141: PetscFunctionReturn(PETSC_SUCCESS);
142: }
143: } else if (snes->npcside == PC_RIGHT) {
144: PetscCall(PetscLogEventBegin(SNES_NPCSolve, snes->npc, X, 0, 0));
145: PetscCall(SNESSolve(snes->npc, snes->vec_rhs, X));
146: PetscCall(PetscLogEventEnd(SNES_NPCSolve, snes->npc, X, 0, 0));
147: PetscCall(SNESGetConvergedReason(snes->npc, &reason));
148: if (reason < 0 && reason != SNES_DIVERGED_MAX_IT && reason != SNES_DIVERGED_TR_DELTA) {
149: snes->reason = SNES_DIVERGED_INNER;
150: PetscFunctionReturn(PETSC_SUCCESS);
151: }
152: PetscCall(SNESGetNPCFunction(snes, F, &fnorm));
153: PetscCall(VecCopy(F, D));
154: } else {
155: PetscCall(VecCopy(F, D));
156: }
157: } else {
158: PetscCall(VecCopy(F, D));
159: }
161: /* scale the initial update */
162: if (qn->scale_type == SNES_QN_SCALE_JACOBIAN && restart) {
163: PetscCall(SNESComputeJacobian(snes, X, snes->jacobian, snes->jacobian_pre));
164: SNESCheckJacobianDomainerror(snes);
165: PetscCall(KSPSetOperators(snes->ksp, snes->jacobian, snes->jacobian_pre));
166: PetscCall(MatLMVMSetJ0KSP(qn->B, snes->ksp));
167: }
169: /* update QN approx and calculate step */
170: PetscCall(MatLMVMUpdate(qn->B, X, D));
171: PetscCall(MatSolve(qn->B, D, Y));
172: PetscCall(VecDot(F, Y, &gdx));
173: if (PetscAbsScalar(gdx) <= 0) {
174: /* Step is not descent or solve was not successful
175: Use steepest descent direction */
176: PetscCall(VecCopy(F, Y));
177: PetscCall(MatLMVMReset(qn->B, PETSC_FALSE));
178: }
180: /* line search for lambda */
181: ynorm = 1;
182: #if defined(PETSC_USE_INFO)
183: gnorm = fnorm;
184: #endif
185: PetscCall(VecCopy(D, Dold));
186: PetscCall(SNESLineSearchApply(snes->linesearch, X, F, &fnorm, Y));
187: if (snes->reason == SNES_DIVERGED_FUNCTION_COUNT) break;
188: PetscCall(SNESLineSearchGetReason(snes->linesearch, &lssucceed));
189: PetscCall(SNESLineSearchGetNorms(snes->linesearch, &xnorm, &fnorm, &ynorm));
190: badstep = PETSC_FALSE;
191: if (lssucceed) {
192: if (++snes->numFailures >= snes->maxFailures) {
193: snes->reason = SNES_DIVERGED_LINE_SEARCH;
194: break;
195: }
196: badstep = PETSC_TRUE;
197: }
199: /* convergence monitoring */
200: PetscCall(PetscInfo(snes, "fnorm=%18.16e, gnorm=%18.16e, ynorm=%18.16e, lssucceed=%d\n", (double)fnorm, (double)gnorm, (double)ynorm, (int)lssucceed));
202: PetscCall(PetscObjectSAWsTakeAccess((PetscObject)snes));
203: snes->iter = i + 1;
204: snes->norm = fnorm;
205: snes->xnorm = xnorm;
206: snes->ynorm = ynorm;
207: PetscCall(PetscObjectSAWsGrantAccess((PetscObject)snes));
209: PetscCall(SNESLogConvergenceHistory(snes, snes->norm, snes->iter));
211: /* set parameter for default relative tolerance convergence test */
212: PetscCall(SNESConverged(snes, snes->iter, xnorm, ynorm, fnorm));
213: PetscCall(SNESMonitor(snes, snes->iter, snes->norm));
214: if (snes->reason) PetscFunctionReturn(PETSC_SUCCESS);
216: if (snes->npc && snes->npcside == PC_LEFT && snes->functype == SNES_FUNCTION_UNPRECONDITIONED) {
217: PetscCall(SNESApplyNPC(snes, X, F, D));
218: PetscCall(SNESGetConvergedReason(snes->npc, &reason));
219: if (reason < 0 && reason != SNES_DIVERGED_MAX_IT && reason != SNES_DIVERGED_TR_DELTA) {
220: snes->reason = SNES_DIVERGED_INNER;
221: PetscFunctionReturn(PETSC_SUCCESS);
222: }
223: } else {
224: PetscCall(VecCopy(F, D));
225: }
227: /* restart conditions */
228: powell = PETSC_FALSE;
229: if (qn->restart_type == SNES_QN_RESTART_POWELL && i_r > 1) {
230: /* check restart by Powell's Criterion: |F^T H_0 Fold| > powell_gamma * |Fold^T H_0 Fold| */
231: if (qn->scale_type == SNES_QN_SCALE_JACOBIAN) {
232: PetscCall(MatMult(snes->jacobian, Dold, W));
233: } else {
234: PetscCall(VecCopy(Dold, W));
235: }
236: PetscCall(VecDotBegin(W, Dold, &DolddotDold));
237: PetscCall(VecDotBegin(W, D, &DolddotD));
238: PetscCall(VecDotEnd(W, Dold, &DolddotDold));
239: PetscCall(VecDotEnd(W, D, &DolddotD));
240: if (PetscAbs(PetscRealPart(DolddotD)) > qn->powell_gamma * PetscAbs(PetscRealPart(DolddotDold))) powell = PETSC_TRUE;
241: }
242: periodic = PETSC_FALSE;
243: if (qn->restart_type == SNES_QN_RESTART_PERIODIC) {
244: if (i_r > qn->m - 1) periodic = PETSC_TRUE;
245: }
247: /* restart if either powell or periodic restart is satisfied. */
248: restart = PETSC_FALSE;
249: if (badstep || powell || periodic) {
250: restart = PETSC_TRUE;
251: if (qn->monflg) {
252: PetscCall(PetscViewerASCIIAddTab(qn->monitor, ((PetscObject)snes)->tablevel + 2));
253: if (powell) {
254: PetscCall(PetscViewerASCIIPrintf(qn->monitor, "Powell restart! |%14.12e| > %6.4f*|%14.12e| i_r = %" PetscInt_FMT "\n", (double)PetscRealPart(DolddotD), (double)qn->powell_gamma, (double)PetscRealPart(DolddotDold), i_r));
255: } else {
256: PetscCall(PetscViewerASCIIPrintf(qn->monitor, "Periodic restart! i_r = %" PetscInt_FMT "\n", i_r));
257: }
258: PetscCall(PetscViewerASCIISubtractTab(qn->monitor, ((PetscObject)snes)->tablevel + 2));
259: }
260: i_r = -1;
261: PetscCall(MatLMVMReset(qn->B, PETSC_FALSE));
262: }
263: }
264: PetscFunctionReturn(PETSC_SUCCESS);
265: }
267: static PetscErrorCode SNESSetUp_QN(SNES snes)
268: {
269: SNES_QN *qn = (SNES_QN *)snes->data;
270: DM dm;
271: PetscInt n, N;
273: PetscFunctionBegin;
274: if (!snes->vec_sol) {
275: PetscCall(SNESGetDM(snes, &dm));
276: PetscCall(DMCreateGlobalVector(dm, &snes->vec_sol));
277: }
278: PetscCall(SNESSetWorkVecs(snes, 3));
279: PetscCall(SNESQNGetMatrix_Private(snes, &qn->B));
281: if (qn->scale_type == SNES_QN_SCALE_JACOBIAN) PetscCall(SNESSetUpMatrices(snes));
282: if (snes->npcside == PC_LEFT && snes->functype == SNES_FUNCTION_DEFAULT) snes->functype = SNES_FUNCTION_PRECONDITIONED;
284: /* set method defaults */
285: if (qn->scale_type == SNES_QN_SCALE_DEFAULT) {
286: if (qn->type == SNES_QN_BADBROYDEN) {
287: qn->scale_type = SNES_QN_SCALE_NONE;
288: } else {
289: qn->scale_type = SNES_QN_SCALE_SCALAR;
290: }
291: }
292: if (qn->restart_type == SNES_QN_RESTART_DEFAULT) {
293: if (qn->type == SNES_QN_LBFGS) {
294: qn->restart_type = SNES_QN_RESTART_POWELL;
295: } else {
296: qn->restart_type = SNES_QN_RESTART_PERIODIC;
297: }
298: }
299: /* Set up the LMVM matrix */
300: PetscCall(SNESQNGetMatrix_Private(snes, &qn->B));
301: PetscCall(VecGetLocalSize(snes->vec_sol, &n));
302: PetscCall(VecGetSize(snes->vec_sol, &N));
303: PetscCall(MatSetSizes(qn->B, n, n, N, N));
304: PetscCall(MatSetUp(qn->B));
305: PetscCall(MatLMVMReset(qn->B, PETSC_TRUE));
306: PetscCall(MatLMVMSetHistorySize(qn->B, qn->m));
307: PetscCall(MatLMVMAllocate(qn->B, snes->vec_sol, snes->vec_func));
308: PetscFunctionReturn(PETSC_SUCCESS);
309: }
311: static PetscErrorCode SNESReset_QN(SNES snes)
312: {
313: SNES_QN *qn = (SNES_QN *)snes->data;
315: PetscFunctionBegin;
316: PetscCall(MatDestroy(&qn->B));
317: PetscFunctionReturn(PETSC_SUCCESS);
318: }
320: static PetscErrorCode SNESDestroy_QN(SNES snes)
321: {
322: PetscFunctionBegin;
323: PetscCall(SNESReset_QN(snes));
324: PetscCall(PetscFree(snes->data));
325: PetscCall(PetscObjectComposeFunction((PetscObject)snes, "SNESQNSetScaleType_C", NULL));
326: PetscCall(PetscObjectComposeFunction((PetscObject)snes, "SNESQNSetRestartType_C", NULL));
327: PetscCall(PetscObjectComposeFunction((PetscObject)snes, "SNESQNSetType_C", NULL));
328: PetscFunctionReturn(PETSC_SUCCESS);
329: }
331: static PetscErrorCode SNESSetFromOptions_QN(SNES snes, PetscOptionItems *PetscOptionsObject)
332: {
333: SNES_QN *qn = (SNES_QN *)snes->data;
334: PetscBool flg;
335: SNESLineSearch linesearch;
336: SNESQNRestartType rtype = qn->restart_type;
337: SNESQNScaleType stype = qn->scale_type;
338: SNESQNType qtype = qn->type;
340: PetscFunctionBegin;
341: PetscOptionsHeadBegin(PetscOptionsObject, "SNES QN options");
342: PetscCall(PetscOptionsInt("-snes_qn_m", "Number of past states saved for L-BFGS methods", "SNESQN", qn->m, &qn->m, NULL));
343: PetscCall(PetscOptionsReal("-snes_qn_powell_gamma", "Powell angle tolerance", "SNESQN", qn->powell_gamma, &qn->powell_gamma, NULL));
344: PetscCall(PetscOptionsBool("-snes_qn_monitor", "Monitor for the QN methods", "SNESQN", qn->monflg, &qn->monflg, NULL));
345: PetscCall(PetscOptionsEnum("-snes_qn_scale_type", "Scaling type", "SNESQNSetScaleType", SNESQNScaleTypes, (PetscEnum)stype, (PetscEnum *)&stype, &flg));
346: if (flg) PetscCall(SNESQNSetScaleType(snes, stype));
348: PetscCall(PetscOptionsEnum("-snes_qn_restart_type", "Restart type", "SNESQNSetRestartType", SNESQNRestartTypes, (PetscEnum)rtype, (PetscEnum *)&rtype, &flg));
349: if (flg) PetscCall(SNESQNSetRestartType(snes, rtype));
351: PetscCall(PetscOptionsEnum("-snes_qn_type", "Quasi-Newton update type", "", SNESQNTypes, (PetscEnum)qtype, (PetscEnum *)&qtype, &flg));
352: if (flg) PetscCall(SNESQNSetType(snes, qtype));
353: PetscOptionsHeadEnd();
354: PetscCall(SNESQNGetMatrix_Private(snes, &qn->B));
355: PetscCall(MatSetFromOptions(qn->B));
356: if (!snes->linesearch) {
357: PetscCall(SNESGetLineSearch(snes, &linesearch));
358: if (!((PetscObject)linesearch)->type_name) {
359: if (qn->type == SNES_QN_LBFGS) {
360: PetscCall(SNESLineSearchSetType(linesearch, SNESLINESEARCHCP));
361: } else if (qn->type == SNES_QN_BROYDEN) {
362: PetscCall(SNESLineSearchSetType(linesearch, SNESLINESEARCHBASIC));
363: } else {
364: PetscCall(SNESLineSearchSetType(linesearch, SNESLINESEARCHL2));
365: }
366: }
367: }
368: if (qn->monflg) PetscCall(PetscViewerASCIIGetStdout(PetscObjectComm((PetscObject)snes), &qn->monitor));
369: PetscFunctionReturn(PETSC_SUCCESS);
370: }
372: static PetscErrorCode SNESView_QN(SNES snes, PetscViewer viewer)
373: {
374: SNES_QN *qn = (SNES_QN *)snes->data;
375: PetscBool iascii;
377: PetscFunctionBegin;
378: PetscCall(PetscObjectTypeCompare((PetscObject)viewer, PETSCVIEWERASCII, &iascii));
379: if (iascii) {
380: PetscCall(PetscViewerASCIIPrintf(viewer, " type is %s, restart type is %s, scale type is %s\n", SNESQNTypes[qn->type], SNESQNRestartTypes[qn->restart_type], SNESQNScaleTypes[qn->scale_type]));
381: PetscCall(PetscViewerASCIIPrintf(viewer, " Stored subspace size: %" PetscInt_FMT "\n", qn->m));
382: PetscCall(MatView(qn->B, viewer));
383: }
384: PetscFunctionReturn(PETSC_SUCCESS);
385: }
387: /*@
388: SNESQNSetRestartType - Sets the restart type for `SNESQN`.
390: Logically Collective
392: Input Parameters:
393: + snes - the iterative context
394: - rtype - restart type, see `SNESQNRestartType`
396: Options Database Keys:
397: + -snes_qn_restart_type <powell,periodic,none> - set the restart type
398: - -snes_qn_m <m> - sets the number of stored updates and the restart period for periodic
400: Level: intermediate
402: .seealso: [](ch_snes), `SNES`, `SNESQN`, `SNESQNRestartType`, `SNES_QN_RESTART_NONE`, `SNES_QN_RESTART_POWELL`, `SNES_QN_RESTART_PERIODIC`,
403: `SNESQNType`, `SNESQNScaleType`
404: @*/
405: PetscErrorCode SNESQNSetRestartType(SNES snes, SNESQNRestartType rtype)
406: {
407: PetscFunctionBegin;
409: PetscTryMethod(snes, "SNESQNSetRestartType_C", (SNES, SNESQNRestartType), (snes, rtype));
410: PetscFunctionReturn(PETSC_SUCCESS);
411: }
413: /*@
414: SNESQNSetScaleType - Sets the scaling type for the inner inverse Jacobian in `SNESQN`.
416: Logically Collective
418: Input Parameters:
419: + snes - the nonlinear solver context
420: - stype - scale type, see `SNESQNScaleType`
422: Options Database Key:
423: . -snes_qn_scale_type <diagonal,none,scalar,jacobian> - Scaling type
425: Level: intermediate
427: .seealso: [](ch_snes), `SNES`, `SNESQN`, `SNESLineSearch`, `SNESQNScaleType`, `SNESSetJacobian()`, `SNESQNType`, `SNESQNRestartType`
428: @*/
429: PetscErrorCode SNESQNSetScaleType(SNES snes, SNESQNScaleType stype)
430: {
431: PetscFunctionBegin;
433: PetscTryMethod(snes, "SNESQNSetScaleType_C", (SNES, SNESQNScaleType), (snes, stype));
434: PetscFunctionReturn(PETSC_SUCCESS);
435: }
437: static PetscErrorCode SNESQNSetScaleType_QN(SNES snes, SNESQNScaleType stype)
438: {
439: SNES_QN *qn = (SNES_QN *)snes->data;
441: PetscFunctionBegin;
442: qn->scale_type = stype;
443: if (stype == SNES_QN_SCALE_JACOBIAN) snes->usesksp = PETSC_TRUE;
444: PetscFunctionReturn(PETSC_SUCCESS);
445: }
447: static PetscErrorCode SNESQNSetRestartType_QN(SNES snes, SNESQNRestartType rtype)
448: {
449: SNES_QN *qn = (SNES_QN *)snes->data;
451: PetscFunctionBegin;
452: qn->restart_type = rtype;
453: PetscFunctionReturn(PETSC_SUCCESS);
454: }
456: /*@
457: SNESQNSetType - Sets the quasi-Newton variant to be used in `SNESQN`.
459: Logically Collective
461: Input Parameters:
462: + snes - the iterative context
463: - qtype - variant type, see `SNESQNType`
465: Options Database Key:
466: . -snes_qn_type <lbfgs,broyden,badbroyden> - quasi-Newton type
468: Level: intermediate
470: .seealso: [](ch_snes), `SNESQN`, `SNES_QN_LBFGS`, `SNES_QN_BROYDEN`, `SNES_QN_BADBROYDEN`, `SNESQNType`, `SNESQNScaleType`, `TAOLMVM`, `TAOBLMVM`
471: @*/
472: PetscErrorCode SNESQNSetType(SNES snes, SNESQNType qtype)
473: {
474: PetscFunctionBegin;
476: PetscTryMethod(snes, "SNESQNSetType_C", (SNES, SNESQNType), (snes, qtype));
477: PetscFunctionReturn(PETSC_SUCCESS);
478: }
480: static PetscErrorCode SNESQNSetType_QN(SNES snes, SNESQNType qtype)
481: {
482: SNES_QN *qn = (SNES_QN *)snes->data;
484: PetscFunctionBegin;
485: qn->type = qtype;
486: PetscFunctionReturn(PETSC_SUCCESS);
487: }
489: /*MC
490: SNESQN - Limited-Memory Quasi-Newton methods for the solution of nonlinear systems.
492: Options Database Keys:
493: + -snes_qn_m <m> - Number of past states saved for the L-Broyden methods.
494: . -snes_qn_restart_type <powell,periodic,none> - set the restart type
495: . -snes_qn_powell_gamma - Angle condition for restart.
496: . -snes_qn_powell_descent - Descent condition for restart.
497: . -snes_qn_type <lbfgs,broyden,badbroyden> - QN type
498: . -snes_qn_scale_type <diagonal,none,scalar,jacobian> - scaling performed on inner Jacobian
499: . -snes_linesearch_type <cp, l2, basic> - Type of line search.
500: - -snes_qn_monitor - Monitors the quasi-newton Jacobian.
502: Level: beginner
504: Notes:
505: This implements the L-BFGS, Broyden, and "Bad" Broyden algorithms for the solution of $F(x) = b$ using
506: previous change in $F(x)$ and $x$ to form the approximate inverse Jacobian using a series of multiplicative rank-one
507: updates.
509: When using a nonlinear preconditioner, one has two options as to how the preconditioner is applied. The first of
510: these options, sequential, uses the preconditioner to generate a new solution and function and uses those at this
511: iteration as the current iteration's values when constructing the approximate Jacobian. The second, composed,
512: perturbs the problem the Jacobian represents to be $P(x, b) - x = 0 $, where $P(x, b)$ is the preconditioner.
514: Uses left nonlinear preconditioning by default.
516: See {cite}`byrd1994representations`, {cite}`kelley95`, {cite}`brown1985experiments`, {cite}`griewank2012broyden`, {cite}`gilbert1989some`,
517: {cite}`dener2019accelerating`, and {cite}`bruneknepleysmithtu15`
519: .seealso: [](ch_snes), `SNESQNRestartType`, `SNESQNSetRestartType()`, `SNESCreate()`, `SNES`, `SNESSetType()`, `SNESNEWTONLS`, `SNESNEWTONTR`,
520: `SNESQNScaleType`, `SNESQNSetScaleType()`, `SNESQNType`, `SNESQNSetType()`
521: M*/
522: PETSC_EXTERN PetscErrorCode SNESCreate_QN(SNES snes)
523: {
524: SNES_QN *qn;
526: PetscFunctionBegin;
527: snes->ops->setup = SNESSetUp_QN;
528: snes->ops->solve = SNESSolve_QN;
529: snes->ops->destroy = SNESDestroy_QN;
530: snes->ops->setfromoptions = SNESSetFromOptions_QN;
531: snes->ops->view = SNESView_QN;
532: snes->ops->reset = SNESReset_QN;
534: snes->npcside = PC_LEFT;
536: snes->usesnpc = PETSC_TRUE;
537: snes->usesksp = PETSC_FALSE;
539: snes->alwayscomputesfinalresidual = PETSC_TRUE;
541: if (!snes->tolerancesset) {
542: snes->max_funcs = 30000;
543: snes->max_its = 10000;
544: }
546: PetscCall(PetscNew(&qn));
547: snes->data = (void *)qn;
548: qn->m = 10;
549: qn->scaling = 1.0;
550: qn->monitor = NULL;
551: qn->monflg = PETSC_FALSE;
552: qn->powell_gamma = 0.9999;
553: qn->scale_type = SNES_QN_SCALE_DEFAULT;
554: qn->restart_type = SNES_QN_RESTART_DEFAULT;
555: qn->type = SNES_QN_LBFGS;
557: PetscCall(PetscObjectComposeFunction((PetscObject)snes, "SNESQNSetScaleType_C", SNESQNSetScaleType_QN));
558: PetscCall(PetscObjectComposeFunction((PetscObject)snes, "SNESQNSetRestartType_C", SNESQNSetRestartType_QN));
559: PetscCall(PetscObjectComposeFunction((PetscObject)snes, "SNESQNSetType_C", SNESQNSetType_QN));
560: PetscFunctionReturn(PETSC_SUCCESS);
561: }