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 SNESSolve_QN(SNES snes)
24: {
25: SNES_QN *qn = (SNES_QN *)snes->data;
26: Vec X, Xold;
27: Vec F, W;
28: Vec Y, D, Dold;
29: PetscInt i, i_r;
30: PetscReal fnorm, xnorm, ynorm;
31: SNESLineSearchReason lssucceed;
32: PetscBool badstep, powell, periodic;
33: PetscScalar DolddotD, DolddotDold;
34: SNESConvergedReason reason;
35: #if defined(PETSC_USE_INFO)
36: PetscReal gnorm;
37: #endif
39: /* basically just a regular newton's method except for the application of the Jacobian */
41: PetscFunctionBegin;
42: 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);
44: PetscCall(PetscCitationsRegister(SNESCitation, &SNEScite));
45: F = snes->vec_func; /* residual vector */
46: Y = snes->vec_sol_update; /* search direction generated by J^-1D*/
47: W = snes->work[3];
48: X = snes->vec_sol; /* solution vector */
49: Xold = snes->work[0];
51: /* directions generated by the preconditioned problem with F_pre = F or x - M(x, b) */
52: D = snes->work[1];
53: Dold = snes->work[2];
55: snes->reason = SNES_CONVERGED_ITERATING;
57: PetscCall(PetscObjectSAWsTakeAccess((PetscObject)snes));
58: snes->iter = 0;
59: snes->norm = 0.;
60: PetscCall(PetscObjectSAWsGrantAccess((PetscObject)snes));
62: if (snes->npc && snes->npcside == PC_LEFT && snes->functype == SNES_FUNCTION_PRECONDITIONED) {
63: PetscCall(SNESApplyNPC(snes, X, NULL, F));
64: PetscCall(SNESGetConvergedReason(snes->npc, &reason));
65: if (reason < 0 && reason != SNES_DIVERGED_MAX_IT) {
66: snes->reason = SNES_DIVERGED_INNER;
67: PetscFunctionReturn(PETSC_SUCCESS);
68: }
69: PetscCall(VecNorm(F, NORM_2, &fnorm));
70: } else {
71: if (!snes->vec_func_init_set) {
72: PetscCall(SNESComputeFunction(snes, X, F));
73: } else snes->vec_func_init_set = PETSC_FALSE;
75: PetscCall(VecNorm(F, NORM_2, &fnorm));
76: SNESCheckFunctionNorm(snes, fnorm);
77: }
78: if (snes->npc && snes->npcside == PC_LEFT && snes->functype == SNES_FUNCTION_UNPRECONDITIONED) {
79: PetscCall(SNESApplyNPC(snes, X, F, D));
80: PetscCall(SNESGetConvergedReason(snes->npc, &reason));
81: if (reason < 0 && reason != SNES_DIVERGED_MAX_IT) {
82: snes->reason = SNES_DIVERGED_INNER;
83: PetscFunctionReturn(PETSC_SUCCESS);
84: }
85: } else {
86: PetscCall(VecCopy(F, D));
87: }
89: PetscCall(PetscObjectSAWsTakeAccess((PetscObject)snes));
90: snes->norm = fnorm;
91: PetscCall(PetscObjectSAWsGrantAccess((PetscObject)snes));
92: PetscCall(SNESLogConvergenceHistory(snes, fnorm, 0));
93: PetscCall(SNESMonitor(snes, 0, fnorm));
95: /* test convergence */
96: PetscCall(SNESConverged(snes, 0, 0.0, 0.0, fnorm));
97: if (snes->reason) PetscFunctionReturn(PETSC_SUCCESS);
99: if (snes->npc && snes->npcside == PC_RIGHT) {
100: PetscCall(PetscLogEventBegin(SNES_NPCSolve, snes->npc, X, 0, 0));
101: PetscCall(SNESSolve(snes->npc, snes->vec_rhs, X));
102: PetscCall(PetscLogEventEnd(SNES_NPCSolve, snes->npc, X, 0, 0));
103: PetscCall(SNESGetConvergedReason(snes->npc, &reason));
104: if (reason < 0 && reason != SNES_DIVERGED_MAX_IT) {
105: snes->reason = SNES_DIVERGED_INNER;
106: PetscFunctionReturn(PETSC_SUCCESS);
107: }
108: PetscCall(SNESGetNPCFunction(snes, F, &fnorm));
109: PetscCall(VecCopy(F, D));
110: }
112: /* general purpose update */
113: PetscTryTypeMethod(snes, update, snes->iter);
115: /* scale the initial update */
116: if (qn->scale_type == SNES_QN_SCALE_JACOBIAN) {
117: PetscCall(SNESComputeJacobian(snes, X, snes->jacobian, snes->jacobian_pre));
118: SNESCheckJacobianDomainerror(snes);
119: PetscCall(KSPSetOperators(snes->ksp, snes->jacobian, snes->jacobian_pre));
120: PetscCall(MatLMVMSetJ0KSP(qn->B, snes->ksp));
121: }
123: for (i = 0, i_r = 0; i < snes->max_its; i++, i_r++) {
124: /* update QN approx and calculate step */
125: PetscCall(MatLMVMUpdate(qn->B, X, D));
126: PetscCall(MatSolve(qn->B, D, Y));
128: /* line search for lambda */
129: ynorm = 1;
130: #if defined(PETSC_USE_INFO)
131: gnorm = fnorm;
132: #endif
133: PetscCall(VecCopy(D, Dold));
134: PetscCall(VecCopy(X, Xold));
135: PetscCall(SNESLineSearchApply(snes->linesearch, X, F, &fnorm, Y));
136: if (snes->reason == SNES_DIVERGED_FUNCTION_COUNT) break;
137: PetscCall(SNESLineSearchGetReason(snes->linesearch, &lssucceed));
138: PetscCall(SNESLineSearchGetNorms(snes->linesearch, &xnorm, &fnorm, &ynorm));
139: badstep = PETSC_FALSE;
140: if (lssucceed) {
141: if (++snes->numFailures >= snes->maxFailures) {
142: snes->reason = SNES_DIVERGED_LINE_SEARCH;
143: break;
144: }
145: badstep = PETSC_TRUE;
146: }
148: /* convergence monitoring */
149: PetscCall(PetscInfo(snes, "fnorm=%18.16e, gnorm=%18.16e, ynorm=%18.16e, lssucceed=%d\n", (double)fnorm, (double)gnorm, (double)ynorm, (int)lssucceed));
151: if (snes->npc && snes->npcside == PC_RIGHT) {
152: PetscCall(PetscLogEventBegin(SNES_NPCSolve, snes->npc, X, 0, 0));
153: PetscCall(SNESSolve(snes->npc, snes->vec_rhs, X));
154: PetscCall(PetscLogEventEnd(SNES_NPCSolve, snes->npc, X, 0, 0));
155: PetscCall(SNESGetConvergedReason(snes->npc, &reason));
156: if (reason < 0 && reason != SNES_DIVERGED_MAX_IT) {
157: snes->reason = SNES_DIVERGED_INNER;
158: PetscFunctionReturn(PETSC_SUCCESS);
159: }
160: PetscCall(SNESGetNPCFunction(snes, F, &fnorm));
161: }
163: PetscCall(SNESSetIterationNumber(snes, i + 1));
164: snes->norm = fnorm;
165: snes->xnorm = xnorm;
166: snes->ynorm = ynorm;
168: PetscCall(SNESLogConvergenceHistory(snes, snes->norm, snes->iter));
170: /* set parameter for default relative tolerance convergence test */
171: PetscCall(SNESConverged(snes, snes->iter, xnorm, ynorm, fnorm));
172: PetscCall(SNESMonitor(snes, snes->iter, snes->norm));
173: if (snes->reason) PetscFunctionReturn(PETSC_SUCCESS);
174: if (snes->npc && snes->npcside == PC_LEFT && snes->functype == SNES_FUNCTION_UNPRECONDITIONED) {
175: PetscCall(SNESApplyNPC(snes, X, F, D));
176: PetscCall(SNESGetConvergedReason(snes->npc, &reason));
177: if (reason < 0 && reason != SNES_DIVERGED_MAX_IT) {
178: snes->reason = SNES_DIVERGED_INNER;
179: PetscFunctionReturn(PETSC_SUCCESS);
180: }
181: } else {
182: PetscCall(VecCopy(F, D));
183: }
185: /* general purpose update */
186: PetscTryTypeMethod(snes, update, snes->iter);
188: /* restart conditions */
189: powell = PETSC_FALSE;
190: if (qn->restart_type == SNES_QN_RESTART_POWELL && i_r > 1) {
191: /* check restart by Powell's Criterion: |F^T H_0 Fold| > powell_gamma * |Fold^T H_0 Fold| */
192: if (qn->scale_type == SNES_QN_SCALE_JACOBIAN) {
193: PetscCall(MatMult(snes->jacobian, Dold, W));
194: } else {
195: PetscCall(VecCopy(Dold, W));
196: }
197: PetscCall(VecDotBegin(W, Dold, &DolddotDold));
198: PetscCall(VecDotBegin(W, D, &DolddotD));
199: PetscCall(VecDotEnd(W, Dold, &DolddotDold));
200: PetscCall(VecDotEnd(W, D, &DolddotD));
201: if (PetscAbs(PetscRealPart(DolddotD)) > qn->powell_gamma * PetscAbs(PetscRealPart(DolddotDold))) powell = PETSC_TRUE;
202: }
203: periodic = PETSC_FALSE;
204: if (qn->restart_type == SNES_QN_RESTART_PERIODIC) {
205: if (i_r > qn->m - 1) periodic = PETSC_TRUE;
206: }
207: /* restart if either powell or periodic restart is satisfied. */
208: if (badstep || powell || periodic) {
209: if (qn->monflg) {
210: PetscCall(PetscViewerASCIIAddTab(qn->monitor, ((PetscObject)snes)->tablevel + 2));
211: if (powell) {
212: 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));
213: } else {
214: PetscCall(PetscViewerASCIIPrintf(qn->monitor, "Periodic restart! i_r = %" PetscInt_FMT "\n", i_r));
215: }
216: PetscCall(PetscViewerASCIISubtractTab(qn->monitor, ((PetscObject)snes)->tablevel + 2));
217: }
218: i_r = -1;
219: if (qn->scale_type == SNES_QN_SCALE_JACOBIAN) {
220: PetscCall(SNESComputeJacobian(snes, X, snes->jacobian, snes->jacobian_pre));
221: SNESCheckJacobianDomainerror(snes);
222: }
223: PetscCall(MatLMVMReset(qn->B, PETSC_FALSE));
224: }
225: }
226: PetscFunctionReturn(PETSC_SUCCESS);
227: }
229: static PetscErrorCode SNESSetUp_QN(SNES snes)
230: {
231: SNES_QN *qn = (SNES_QN *)snes->data;
232: DM dm;
233: PetscInt n, N;
235: PetscFunctionBegin;
237: if (!snes->vec_sol) {
238: PetscCall(SNESGetDM(snes, &dm));
239: PetscCall(DMCreateGlobalVector(dm, &snes->vec_sol));
240: }
241: PetscCall(SNESSetWorkVecs(snes, 4));
243: if (qn->scale_type == SNES_QN_SCALE_JACOBIAN) PetscCall(SNESSetUpMatrices(snes));
244: if (snes->npcside == PC_LEFT && snes->functype == SNES_FUNCTION_DEFAULT) snes->functype = SNES_FUNCTION_UNPRECONDITIONED;
246: /* set method defaults */
247: if (qn->scale_type == SNES_QN_SCALE_DEFAULT) {
248: if (qn->type == SNES_QN_BADBROYDEN) {
249: qn->scale_type = SNES_QN_SCALE_NONE;
250: } else {
251: qn->scale_type = SNES_QN_SCALE_SCALAR;
252: }
253: }
254: if (qn->restart_type == SNES_QN_RESTART_DEFAULT) {
255: if (qn->type == SNES_QN_LBFGS) {
256: qn->restart_type = SNES_QN_RESTART_POWELL;
257: } else {
258: qn->restart_type = SNES_QN_RESTART_PERIODIC;
259: }
260: }
262: /* Set up the LMVM matrix */
263: switch (qn->type) {
264: case SNES_QN_BROYDEN:
265: PetscCall(MatSetType(qn->B, MATLMVMBROYDEN));
266: qn->scale_type = SNES_QN_SCALE_NONE;
267: break;
268: case SNES_QN_BADBROYDEN:
269: PetscCall(MatSetType(qn->B, MATLMVMBADBROYDEN));
270: qn->scale_type = SNES_QN_SCALE_NONE;
271: break;
272: default:
273: PetscCall(MatSetType(qn->B, MATLMVMBFGS));
274: switch (qn->scale_type) {
275: case SNES_QN_SCALE_NONE:
276: PetscCall(MatLMVMSymBroydenSetScaleType(qn->B, MAT_LMVM_SYMBROYDEN_SCALE_NONE));
277: break;
278: case SNES_QN_SCALE_SCALAR:
279: PetscCall(MatLMVMSymBroydenSetScaleType(qn->B, MAT_LMVM_SYMBROYDEN_SCALE_SCALAR));
280: break;
281: case SNES_QN_SCALE_JACOBIAN:
282: PetscCall(MatLMVMSymBroydenSetScaleType(qn->B, MAT_LMVM_SYMBROYDEN_SCALE_USER));
283: break;
284: case SNES_QN_SCALE_DIAGONAL:
285: case SNES_QN_SCALE_DEFAULT:
286: default:
287: break;
288: }
289: break;
290: }
291: PetscCall(VecGetLocalSize(snes->vec_sol, &n));
292: PetscCall(VecGetSize(snes->vec_sol, &N));
293: PetscCall(MatSetSizes(qn->B, n, n, N, N));
294: PetscCall(MatSetUp(qn->B));
295: PetscCall(MatLMVMReset(qn->B, PETSC_TRUE));
296: PetscCall(MatLMVMSetHistorySize(qn->B, qn->m));
297: PetscCall(MatLMVMAllocate(qn->B, snes->vec_sol, snes->vec_func));
298: PetscFunctionReturn(PETSC_SUCCESS);
299: }
301: static PetscErrorCode SNESReset_QN(SNES snes)
302: {
303: SNES_QN *qn;
305: PetscFunctionBegin;
306: if (snes->data) {
307: qn = (SNES_QN *)snes->data;
308: PetscCall(MatDestroy(&qn->B));
309: }
310: PetscFunctionReturn(PETSC_SUCCESS);
311: }
313: static PetscErrorCode SNESDestroy_QN(SNES snes)
314: {
315: PetscFunctionBegin;
316: PetscCall(SNESReset_QN(snes));
317: PetscCall(PetscFree(snes->data));
318: PetscCall(PetscObjectComposeFunction((PetscObject)snes, "SNESQNSetScaleType_C", NULL));
319: PetscCall(PetscObjectComposeFunction((PetscObject)snes, "SNESQNSetRestartType_C", NULL));
320: PetscCall(PetscObjectComposeFunction((PetscObject)snes, "SNESQNSetType_C", NULL));
321: PetscFunctionReturn(PETSC_SUCCESS);
322: }
324: static PetscErrorCode SNESSetFromOptions_QN(SNES snes, PetscOptionItems *PetscOptionsObject)
325: {
326: SNES_QN *qn = (SNES_QN *)snes->data;
327: PetscBool flg;
328: SNESLineSearch linesearch;
329: SNESQNRestartType rtype = qn->restart_type;
330: SNESQNScaleType stype = qn->scale_type;
331: SNESQNType qtype = qn->type;
333: PetscFunctionBegin;
334: PetscOptionsHeadBegin(PetscOptionsObject, "SNES QN options");
335: PetscCall(PetscOptionsInt("-snes_qn_m", "Number of past states saved for L-BFGS methods", "SNESQN", qn->m, &qn->m, NULL));
336: PetscCall(PetscOptionsReal("-snes_qn_powell_gamma", "Powell angle tolerance", "SNESQN", qn->powell_gamma, &qn->powell_gamma, NULL));
337: PetscCall(PetscOptionsBool("-snes_qn_monitor", "Monitor for the QN methods", "SNESQN", qn->monflg, &qn->monflg, NULL));
338: PetscCall(PetscOptionsEnum("-snes_qn_scale_type", "Scaling type", "SNESQNSetScaleType", SNESQNScaleTypes, (PetscEnum)stype, (PetscEnum *)&stype, &flg));
339: if (flg) PetscCall(SNESQNSetScaleType(snes, stype));
341: PetscCall(PetscOptionsEnum("-snes_qn_restart_type", "Restart type", "SNESQNSetRestartType", SNESQNRestartTypes, (PetscEnum)rtype, (PetscEnum *)&rtype, &flg));
342: if (flg) PetscCall(SNESQNSetRestartType(snes, rtype));
344: PetscCall(PetscOptionsEnum("-snes_qn_type", "Quasi-Newton update type", "", SNESQNTypes, (PetscEnum)qtype, (PetscEnum *)&qtype, &flg));
345: if (flg) PetscCall(SNESQNSetType(snes, qtype));
346: PetscCall(MatSetFromOptions(qn->B));
347: PetscOptionsHeadEnd();
348: if (!snes->linesearch) {
349: PetscCall(SNESGetLineSearch(snes, &linesearch));
350: if (!((PetscObject)linesearch)->type_name) {
351: if (qn->type == SNES_QN_LBFGS) {
352: PetscCall(SNESLineSearchSetType(linesearch, SNESLINESEARCHCP));
353: } else if (qn->type == SNES_QN_BROYDEN) {
354: PetscCall(SNESLineSearchSetType(linesearch, SNESLINESEARCHBASIC));
355: } else {
356: PetscCall(SNESLineSearchSetType(linesearch, SNESLINESEARCHL2));
357: }
358: }
359: }
360: if (qn->monflg) PetscCall(PetscViewerASCIIGetStdout(PetscObjectComm((PetscObject)snes), &qn->monitor));
361: PetscFunctionReturn(PETSC_SUCCESS);
362: }
364: static PetscErrorCode SNESView_QN(SNES snes, PetscViewer viewer)
365: {
366: SNES_QN *qn = (SNES_QN *)snes->data;
367: PetscBool iascii;
369: PetscFunctionBegin;
370: PetscCall(PetscObjectTypeCompare((PetscObject)viewer, PETSCVIEWERASCII, &iascii));
371: if (iascii) {
372: 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]));
373: PetscCall(PetscViewerASCIIPrintf(viewer, " Stored subspace size: %" PetscInt_FMT "\n", qn->m));
374: }
375: PetscFunctionReturn(PETSC_SUCCESS);
376: }
378: /*@
379: SNESQNSetRestartType - Sets the restart type for `SNESQN`.
381: Logically Collective
383: Input Parameters:
384: + snes - the iterative context
385: - rtype - restart type, see `SNESQNRestartType`
387: Options Database Keys:
388: + -snes_qn_restart_type <powell,periodic,none> - set the restart type
389: - -snes_qn_m <m> - sets the number of stored updates and the restart period for periodic
391: Level: intermediate
393: .seealso: [](ch_snes), `SNES`, `SNESQN`, `SNESQNRestartType`, `SNES_QN_RESTART_NONE`, `SNES_QN_RESTART_POWELL`, `SNES_QN_RESTART_PERIODIC`,
394: `SNESQNType`, `SNESQNScaleType`
395: @*/
396: PetscErrorCode SNESQNSetRestartType(SNES snes, SNESQNRestartType rtype)
397: {
398: PetscFunctionBegin;
400: PetscTryMethod(snes, "SNESQNSetRestartType_C", (SNES, SNESQNRestartType), (snes, rtype));
401: PetscFunctionReturn(PETSC_SUCCESS);
402: }
404: /*@
405: SNESQNSetScaleType - Sets the scaling type for the inner inverse Jacobian in `SNESQN`.
407: Logically Collective
409: Input Parameters:
410: + snes - the nonlinear solver context
411: - stype - scale type, see `SNESQNScaleType`
413: Options Database Key:
414: . -snes_qn_scale_type <diagonal,none,scalar,jacobian> - Scaling type
416: Level: intermediate
418: .seealso: [](ch_snes), `SNES`, `SNESQN`, `SNESLineSearch`, `SNESQNScaleType`, `SNESSetJacobian()`, `SNESQNType`, `SNESQNRestartType`
419: @*/
420: PetscErrorCode SNESQNSetScaleType(SNES snes, SNESQNScaleType stype)
421: {
422: PetscFunctionBegin;
424: PetscTryMethod(snes, "SNESQNSetScaleType_C", (SNES, SNESQNScaleType), (snes, stype));
425: PetscFunctionReturn(PETSC_SUCCESS);
426: }
428: static PetscErrorCode SNESQNSetScaleType_QN(SNES snes, SNESQNScaleType stype)
429: {
430: SNES_QN *qn = (SNES_QN *)snes->data;
432: PetscFunctionBegin;
433: qn->scale_type = stype;
434: if (stype == SNES_QN_SCALE_JACOBIAN) snes->usesksp = PETSC_TRUE;
435: PetscFunctionReturn(PETSC_SUCCESS);
436: }
438: static PetscErrorCode SNESQNSetRestartType_QN(SNES snes, SNESQNRestartType rtype)
439: {
440: SNES_QN *qn = (SNES_QN *)snes->data;
442: PetscFunctionBegin;
443: qn->restart_type = rtype;
444: PetscFunctionReturn(PETSC_SUCCESS);
445: }
447: /*@
448: SNESQNSetType - Sets the quasi-Newton variant to be used in `SNESQN`.
450: Logically Collective
452: Input Parameters:
453: + snes - the iterative context
454: - qtype - variant type, see `SNESQNType`
456: Options Database Key:
457: . -snes_qn_type <lbfgs,broyden,badbroyden> - quasi-Newton type
459: Level: intermediate
461: .seealso: [](ch_snes), `SNESQN`, `SNES_QN_LBFGS`, `SNES_QN_BROYDEN`, `SNES_QN_BADBROYDEN`, `SNESQNType`, `SNESQNScaleType`, `TAOLMVM`, `TAOBLMVM`
462: @*/
463: PetscErrorCode SNESQNSetType(SNES snes, SNESQNType qtype)
464: {
465: PetscFunctionBegin;
467: PetscTryMethod(snes, "SNESQNSetType_C", (SNES, SNESQNType), (snes, qtype));
468: PetscFunctionReturn(PETSC_SUCCESS);
469: }
471: static PetscErrorCode SNESQNSetType_QN(SNES snes, SNESQNType qtype)
472: {
473: SNES_QN *qn = (SNES_QN *)snes->data;
475: PetscFunctionBegin;
476: qn->type = qtype;
477: PetscFunctionReturn(PETSC_SUCCESS);
478: }
480: /*MC
481: SNESQN - Limited-Memory Quasi-Newton methods for the solution of nonlinear systems.
483: Options Database Keys:
484: + -snes_qn_m <m> - Number of past states saved for the L-Broyden methods.
485: . -snes_qn_restart_type <powell,periodic,none> - set the restart type
486: . -snes_qn_powell_gamma - Angle condition for restart.
487: . -snes_qn_powell_descent - Descent condition for restart.
488: . -snes_qn_type <lbfgs,broyden,badbroyden> - QN type
489: . -snes_qn_scale_type <diagonal,none,scalar,jacobian> - scaling performed on inner Jacobian
490: . -snes_linesearch_type <cp, l2, basic> - Type of line search.
491: - -snes_qn_monitor - Monitors the quasi-newton Jacobian.
493: Level: beginner
495: Notes:
496: This implements the L-BFGS, Broyden, and "Bad" Broyden algorithms for the solution of $F(x) = b$ using
497: previous change in $F(x)$ and $x$ to form the approximate inverse Jacobian using a series of multiplicative rank-one
498: updates.
500: When using a nonlinear preconditioner, one has two options as to how the preconditioner is applied. The first of
501: these options, sequential, uses the preconditioner to generate a new solution and function and uses those at this
502: iteration as the current iteration's values when constructing the approximate Jacobian. The second, composed,
503: perturbs the problem the Jacobian represents to be $P(x, b) - x = 0 $, where $P(x, b)$ is the preconditioner.
505: Uses left nonlinear preconditioning by default.
507: See {cite}`byrd1994representations`, {cite}`kelley95`, {cite}`brown1985experiments`, {cite}`griewank2012broyden`, {cite}`gilbert1989some`,
508: {cite}`dener2019accelerating`, and {cite}`bruneknepleysmithtu15`
510: .seealso: [](ch_snes), `SNESQNRestartType`, `SNESQNSetRestartType()`, `SNESCreate()`, `SNES`, `SNESSetType()`, `SNESNEWTONLS`, `SNESNEWTONTR`,
511: `SNESQNScaleType`, `SNESQNSetScaleType()`, `SNESQNType`, `SNESQNSetType()`
512: M*/
513: PETSC_EXTERN PetscErrorCode SNESCreate_QN(SNES snes)
514: {
515: SNES_QN *qn;
516: const char *optionsprefix;
518: PetscFunctionBegin;
519: snes->ops->setup = SNESSetUp_QN;
520: snes->ops->solve = SNESSolve_QN;
521: snes->ops->destroy = SNESDestroy_QN;
522: snes->ops->setfromoptions = SNESSetFromOptions_QN;
523: snes->ops->view = SNESView_QN;
524: snes->ops->reset = SNESReset_QN;
526: snes->npcside = PC_LEFT;
528: snes->usesnpc = PETSC_TRUE;
529: snes->usesksp = PETSC_FALSE;
531: snes->alwayscomputesfinalresidual = PETSC_TRUE;
533: if (!snes->tolerancesset) {
534: snes->max_funcs = 30000;
535: snes->max_its = 10000;
536: }
538: PetscCall(PetscNew(&qn));
539: snes->data = (void *)qn;
540: qn->m = 10;
541: qn->scaling = 1.0;
542: qn->monitor = NULL;
543: qn->monflg = PETSC_FALSE;
544: qn->powell_gamma = 0.9999;
545: qn->scale_type = SNES_QN_SCALE_DEFAULT;
546: qn->restart_type = SNES_QN_RESTART_DEFAULT;
547: qn->type = SNES_QN_LBFGS;
549: PetscCall(MatCreate(PetscObjectComm((PetscObject)snes), &qn->B));
550: PetscCall(SNESGetOptionsPrefix(snes, &optionsprefix));
551: PetscCall(MatSetOptionsPrefix(qn->B, optionsprefix));
553: PetscCall(PetscObjectComposeFunction((PetscObject)snes, "SNESQNSetScaleType_C", SNESQNSetScaleType_QN));
554: PetscCall(PetscObjectComposeFunction((PetscObject)snes, "SNESQNSetRestartType_C", SNESQNSetRestartType_QN));
555: PetscCall(PetscObjectComposeFunction((PetscObject)snes, "SNESQNSetType_C", SNESQNSetType_QN));
556: PetscFunctionReturn(PETSC_SUCCESS);
557: }