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, gnorm;
31: SNESLineSearchReason lssucceed;
32: PetscBool badstep, powell, periodic;
33: PetscScalar DolddotD, DolddotDold;
34: SNESConvergedReason reason;
36: /* basically just a regular newton's method except for the application of the Jacobian */
40: PetscCitationsRegister(SNESCitation, &SNEScite);
41: F = snes->vec_func; /* residual vector */
42: Y = snes->vec_sol_update; /* search direction generated by J^-1D*/
43: W = snes->work[3];
44: X = snes->vec_sol; /* solution vector */
45: Xold = snes->work[0];
47: /* directions generated by the preconditioned problem with F_pre = F or x - M(x, b) */
48: D = snes->work[1];
49: Dold = snes->work[2];
51: snes->reason = SNES_CONVERGED_ITERATING;
53: PetscObjectSAWsTakeAccess((PetscObject)snes);
54: snes->iter = 0;
55: snes->norm = 0.;
56: PetscObjectSAWsGrantAccess((PetscObject)snes);
58: if (snes->npc && snes->npcside == PC_LEFT && snes->functype == SNES_FUNCTION_PRECONDITIONED) {
59: SNESApplyNPC(snes, X, NULL, F);
60: SNESGetConvergedReason(snes->npc, &reason);
61: if (reason < 0 && reason != SNES_DIVERGED_MAX_IT) {
62: snes->reason = SNES_DIVERGED_INNER;
63: return 0;
64: }
65: VecNorm(F, NORM_2, &fnorm);
66: } else {
67: if (!snes->vec_func_init_set) {
68: SNESComputeFunction(snes, X, F);
69: } else snes->vec_func_init_set = PETSC_FALSE;
71: VecNorm(F, NORM_2, &fnorm);
72: SNESCheckFunctionNorm(snes, fnorm);
73: }
74: if (snes->npc && snes->npcside == PC_LEFT && snes->functype == SNES_FUNCTION_UNPRECONDITIONED) {
75: SNESApplyNPC(snes, X, F, D);
76: SNESGetConvergedReason(snes->npc, &reason);
77: if (reason < 0 && reason != SNES_DIVERGED_MAX_IT) {
78: snes->reason = SNES_DIVERGED_INNER;
79: return 0;
80: }
81: } else {
82: VecCopy(F, D);
83: }
85: PetscObjectSAWsTakeAccess((PetscObject)snes);
86: snes->norm = fnorm;
87: PetscObjectSAWsGrantAccess((PetscObject)snes);
88: SNESLogConvergenceHistory(snes, fnorm, 0);
89: SNESMonitor(snes, 0, fnorm);
91: /* test convergence */
92: PetscUseTypeMethod(snes, converged, 0, 0.0, 0.0, fnorm, &snes->reason, snes->cnvP);
93: if (snes->reason) return 0;
95: if (snes->npc && snes->npcside == PC_RIGHT) {
96: PetscLogEventBegin(SNES_NPCSolve, snes->npc, X, 0, 0);
97: SNESSolve(snes->npc, snes->vec_rhs, X);
98: PetscLogEventEnd(SNES_NPCSolve, snes->npc, X, 0, 0);
99: SNESGetConvergedReason(snes->npc, &reason);
100: if (reason < 0 && reason != SNES_DIVERGED_MAX_IT) {
101: snes->reason = SNES_DIVERGED_INNER;
102: return 0;
103: }
104: SNESGetNPCFunction(snes, F, &fnorm);
105: VecCopy(F, D);
106: }
108: /* general purpose update */
109: PetscTryTypeMethod(snes, update, snes->iter);
111: /* scale the initial update */
112: if (qn->scale_type == SNES_QN_SCALE_JACOBIAN) {
113: SNESComputeJacobian(snes, X, snes->jacobian, snes->jacobian_pre);
114: SNESCheckJacobianDomainerror(snes);
115: KSPSetOperators(snes->ksp, snes->jacobian, snes->jacobian_pre);
116: MatLMVMSetJ0KSP(qn->B, snes->ksp);
117: }
119: for (i = 0, i_r = 0; i < snes->max_its; i++, i_r++) {
120: /* update QN approx and calculate step */
121: MatLMVMUpdate(qn->B, X, D);
122: MatSolve(qn->B, D, Y);
124: /* line search for lambda */
125: ynorm = 1;
126: gnorm = fnorm;
127: VecCopy(D, Dold);
128: VecCopy(X, Xold);
129: SNESLineSearchApply(snes->linesearch, X, F, &fnorm, Y);
130: if (snes->reason == SNES_DIVERGED_FUNCTION_COUNT) break;
131: SNESLineSearchGetReason(snes->linesearch, &lssucceed);
132: SNESLineSearchGetNorms(snes->linesearch, &xnorm, &fnorm, &ynorm);
133: badstep = PETSC_FALSE;
134: if (lssucceed) {
135: if (++snes->numFailures >= snes->maxFailures) {
136: snes->reason = SNES_DIVERGED_LINE_SEARCH;
137: break;
138: }
139: badstep = PETSC_TRUE;
140: }
142: /* convergence monitoring */
143: PetscInfo(snes, "fnorm=%18.16e, gnorm=%18.16e, ynorm=%18.16e, lssucceed=%d\n", (double)fnorm, (double)gnorm, (double)ynorm, (int)lssucceed);
145: if (snes->npc && snes->npcside == PC_RIGHT) {
146: PetscLogEventBegin(SNES_NPCSolve, snes->npc, X, 0, 0);
147: SNESSolve(snes->npc, snes->vec_rhs, X);
148: PetscLogEventEnd(SNES_NPCSolve, snes->npc, X, 0, 0);
149: SNESGetConvergedReason(snes->npc, &reason);
150: if (reason < 0 && reason != SNES_DIVERGED_MAX_IT) {
151: snes->reason = SNES_DIVERGED_INNER;
152: return 0;
153: }
154: SNESGetNPCFunction(snes, F, &fnorm);
155: }
157: SNESSetIterationNumber(snes, i + 1);
158: snes->norm = fnorm;
159: snes->xnorm = xnorm;
160: snes->ynorm = ynorm;
162: SNESLogConvergenceHistory(snes, snes->norm, snes->iter);
163: SNESMonitor(snes, snes->iter, snes->norm);
165: /* set parameter for default relative tolerance convergence test */
166: PetscUseTypeMethod(snes, converged, snes->iter, xnorm, ynorm, fnorm, &snes->reason, snes->cnvP);
167: if (snes->reason) return 0;
168: if (snes->npc && snes->npcside == PC_LEFT && snes->functype == SNES_FUNCTION_UNPRECONDITIONED) {
169: SNESApplyNPC(snes, X, F, D);
170: SNESGetConvergedReason(snes->npc, &reason);
171: if (reason < 0 && reason != SNES_DIVERGED_MAX_IT) {
172: snes->reason = SNES_DIVERGED_INNER;
173: return 0;
174: }
175: } else {
176: VecCopy(F, D);
177: }
179: /* general purpose update */
180: PetscTryTypeMethod(snes, update, snes->iter);
182: /* restart conditions */
183: powell = PETSC_FALSE;
184: if (qn->restart_type == SNES_QN_RESTART_POWELL && i_r > 1) {
185: /* check restart by Powell's Criterion: |F^T H_0 Fold| > powell_gamma * |Fold^T H_0 Fold| */
186: if (qn->scale_type == SNES_QN_SCALE_JACOBIAN) {
187: MatMult(snes->jacobian, Dold, W);
188: } else {
189: VecCopy(Dold, W);
190: }
191: VecDotBegin(W, Dold, &DolddotDold);
192: VecDotBegin(W, D, &DolddotD);
193: VecDotEnd(W, Dold, &DolddotDold);
194: VecDotEnd(W, D, &DolddotD);
195: if (PetscAbs(PetscRealPart(DolddotD)) > qn->powell_gamma * PetscAbs(PetscRealPart(DolddotDold))) powell = PETSC_TRUE;
196: }
197: periodic = PETSC_FALSE;
198: if (qn->restart_type == SNES_QN_RESTART_PERIODIC) {
199: if (i_r > qn->m - 1) periodic = PETSC_TRUE;
200: }
201: /* restart if either powell or periodic restart is satisfied. */
202: if (badstep || powell || periodic) {
203: if (qn->monflg) {
204: PetscViewerASCIIAddTab(qn->monitor, ((PetscObject)snes)->tablevel + 2);
205: if (powell) {
206: 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);
207: } else {
208: PetscViewerASCIIPrintf(qn->monitor, "Periodic restart! i_r = %" PetscInt_FMT "\n", i_r);
209: }
210: PetscViewerASCIISubtractTab(qn->monitor, ((PetscObject)snes)->tablevel + 2);
211: }
212: i_r = -1;
213: if (qn->scale_type == SNES_QN_SCALE_JACOBIAN) {
214: SNESComputeJacobian(snes, X, snes->jacobian, snes->jacobian_pre);
215: SNESCheckJacobianDomainerror(snes);
216: }
217: MatLMVMReset(qn->B, PETSC_FALSE);
218: }
219: }
220: if (i == snes->max_its) {
221: PetscInfo(snes, "Maximum number of iterations has been reached: %" PetscInt_FMT "\n", snes->max_its);
222: if (!snes->reason) snes->reason = SNES_DIVERGED_MAX_IT;
223: }
224: return 0;
225: }
227: static PetscErrorCode SNESSetUp_QN(SNES snes)
228: {
229: SNES_QN *qn = (SNES_QN *)snes->data;
230: DM dm;
231: PetscInt n, N;
234: if (!snes->vec_sol) {
235: SNESGetDM(snes, &dm);
236: DMCreateGlobalVector(dm, &snes->vec_sol);
237: }
238: SNESSetWorkVecs(snes, 4);
240: if (qn->scale_type == SNES_QN_SCALE_JACOBIAN) SNESSetUpMatrices(snes);
241: if (snes->npcside == PC_LEFT && snes->functype == SNES_FUNCTION_DEFAULT) snes->functype = SNES_FUNCTION_UNPRECONDITIONED;
243: /* set method defaults */
244: if (qn->scale_type == SNES_QN_SCALE_DEFAULT) {
245: if (qn->type == SNES_QN_BADBROYDEN) {
246: qn->scale_type = SNES_QN_SCALE_NONE;
247: } else {
248: qn->scale_type = SNES_QN_SCALE_SCALAR;
249: }
250: }
251: if (qn->restart_type == SNES_QN_RESTART_DEFAULT) {
252: if (qn->type == SNES_QN_LBFGS) {
253: qn->restart_type = SNES_QN_RESTART_POWELL;
254: } else {
255: qn->restart_type = SNES_QN_RESTART_PERIODIC;
256: }
257: }
259: /* Set up the LMVM matrix */
260: switch (qn->type) {
261: case SNES_QN_BROYDEN:
262: MatSetType(qn->B, MATLMVMBROYDEN);
263: qn->scale_type = SNES_QN_SCALE_NONE;
264: break;
265: case SNES_QN_BADBROYDEN:
266: MatSetType(qn->B, MATLMVMBADBROYDEN);
267: qn->scale_type = SNES_QN_SCALE_NONE;
268: break;
269: default:
270: MatSetType(qn->B, MATLMVMBFGS);
271: switch (qn->scale_type) {
272: case SNES_QN_SCALE_NONE:
273: MatLMVMSymBroydenSetScaleType(qn->B, MAT_LMVM_SYMBROYDEN_SCALE_NONE);
274: break;
275: case SNES_QN_SCALE_SCALAR:
276: MatLMVMSymBroydenSetScaleType(qn->B, MAT_LMVM_SYMBROYDEN_SCALE_SCALAR);
277: break;
278: case SNES_QN_SCALE_JACOBIAN:
279: MatLMVMSymBroydenSetScaleType(qn->B, MAT_LMVM_SYMBROYDEN_SCALE_USER);
280: break;
281: case SNES_QN_SCALE_DIAGONAL:
282: case SNES_QN_SCALE_DEFAULT:
283: default:
284: break;
285: }
286: break;
287: }
288: VecGetLocalSize(snes->vec_sol, &n);
289: VecGetSize(snes->vec_sol, &N);
290: MatSetSizes(qn->B, n, n, N, N);
291: MatSetUp(qn->B);
292: MatLMVMReset(qn->B, PETSC_TRUE);
293: MatLMVMSetHistorySize(qn->B, qn->m);
294: MatLMVMAllocate(qn->B, snes->vec_sol, snes->vec_func);
295: return 0;
296: }
298: static PetscErrorCode SNESReset_QN(SNES snes)
299: {
300: SNES_QN *qn;
302: if (snes->data) {
303: qn = (SNES_QN *)snes->data;
304: MatDestroy(&qn->B);
305: }
306: return 0;
307: }
309: static PetscErrorCode SNESDestroy_QN(SNES snes)
310: {
311: SNESReset_QN(snes);
312: PetscFree(snes->data);
313: PetscObjectComposeFunction((PetscObject)snes, "SNESQNSetScaleType_C", NULL);
314: PetscObjectComposeFunction((PetscObject)snes, "SNESQNSetRestartType_C", NULL);
315: PetscObjectComposeFunction((PetscObject)snes, "SNESQNSetType_C", NULL);
316: return 0;
317: }
319: static PetscErrorCode SNESSetFromOptions_QN(SNES snes, PetscOptionItems *PetscOptionsObject)
320: {
321: SNES_QN *qn = (SNES_QN *)snes->data;
322: PetscBool flg;
323: SNESLineSearch linesearch;
324: SNESQNRestartType rtype = qn->restart_type;
325: SNESQNScaleType stype = qn->scale_type;
326: SNESQNType qtype = qn->type;
328: PetscOptionsHeadBegin(PetscOptionsObject, "SNES QN options");
329: PetscOptionsInt("-snes_qn_m", "Number of past states saved for L-BFGS methods", "SNESQN", qn->m, &qn->m, NULL);
330: PetscOptionsReal("-snes_qn_powell_gamma", "Powell angle tolerance", "SNESQN", qn->powell_gamma, &qn->powell_gamma, NULL);
331: PetscOptionsBool("-snes_qn_monitor", "Monitor for the QN methods", "SNESQN", qn->monflg, &qn->monflg, NULL);
332: PetscOptionsEnum("-snes_qn_scale_type", "Scaling type", "SNESQNSetScaleType", SNESQNScaleTypes, (PetscEnum)stype, (PetscEnum *)&stype, &flg);
333: if (flg) SNESQNSetScaleType(snes, stype);
335: PetscOptionsEnum("-snes_qn_restart_type", "Restart type", "SNESQNSetRestartType", SNESQNRestartTypes, (PetscEnum)rtype, (PetscEnum *)&rtype, &flg);
336: if (flg) SNESQNSetRestartType(snes, rtype);
338: PetscOptionsEnum("-snes_qn_type", "Quasi-Newton update type", "", SNESQNTypes, (PetscEnum)qtype, (PetscEnum *)&qtype, &flg);
339: if (flg) SNESQNSetType(snes, qtype);
340: MatSetFromOptions(qn->B);
341: PetscOptionsHeadEnd();
342: if (!snes->linesearch) {
343: SNESGetLineSearch(snes, &linesearch);
344: if (!((PetscObject)linesearch)->type_name) {
345: if (qn->type == SNES_QN_LBFGS) {
346: SNESLineSearchSetType(linesearch, SNESLINESEARCHCP);
347: } else if (qn->type == SNES_QN_BROYDEN) {
348: SNESLineSearchSetType(linesearch, SNESLINESEARCHBASIC);
349: } else {
350: SNESLineSearchSetType(linesearch, SNESLINESEARCHL2);
351: }
352: }
353: }
354: if (qn->monflg) PetscViewerASCIIGetStdout(PetscObjectComm((PetscObject)snes), &qn->monitor);
355: return 0;
356: }
358: static PetscErrorCode SNESView_QN(SNES snes, PetscViewer viewer)
359: {
360: SNES_QN *qn = (SNES_QN *)snes->data;
361: PetscBool iascii;
363: PetscObjectTypeCompare((PetscObject)viewer, PETSCVIEWERASCII, &iascii);
364: if (iascii) {
365: 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]);
366: PetscViewerASCIIPrintf(viewer, " Stored subspace size: %" PetscInt_FMT "\n", qn->m);
367: }
368: return 0;
369: }
371: /*@
372: SNESQNSetRestartType - Sets the restart type for `SNESQN`.
374: Logically Collective
376: Input Parameters:
377: + snes - the iterative context
378: - rtype - restart type
380: Options Database Keys:
381: + -snes_qn_restart_type <powell,periodic,none> - set the restart type
382: - -snes_qn_m <m> - sets the number of stored updates and the restart period for periodic
384: Level: intermediate
386: `SNESQNRestartType`s:
387: + `SNES_QN_RESTART_NONE` - never restart
388: . `SNES_QN_RESTART_POWELL` - restart based upon descent criteria
389: - `SNES_QN_RESTART_PERIODIC` - restart after a fixed number of iterations
391: .seealso: `SNESQN`, `SNESQNRestartType`, `SNES_QN_RESTART_NONE`, `SNES_QN_RESTART_POWELL`, `SNES_QN_RESTART_PERIODIC`
392: @*/
393: PetscErrorCode SNESQNSetRestartType(SNES snes, SNESQNRestartType rtype)
394: {
396: PetscTryMethod(snes, "SNESQNSetRestartType_C", (SNES, SNESQNRestartType), (snes, rtype));
397: return 0;
398: }
400: /*@
401: SNESQNSetScaleType - Sets the scaling type for the inner inverse Jacobian in `SNESQN`.
403: Logically Collective
405: Input Parameters:
406: + snes - the nonlinear solver context
407: - stype - scale type
409: Options Database Key:
410: . -snes_qn_scale_type <diagonal,none,scalar,jacobian> - Scaling type
412: Level: intermediate
414: `SNESQNScaleType`s:
415: + `SNES_QN_SCALE_NONE` - don't scale the problem
416: . `SNES_QN_SCALE_SCALAR` - use Shanno scaling
417: . `SNES_QN_SCALE_DIAGONAL` - scale with a diagonalized BFGS formula (see Gilbert and Lemarechal 1989), available
418: - `SNES_QN_SCALE_JACOBIAN` - scale by solving a linear system coming from the Jacobian you provided with SNESSetJacobian() computed at the first iteration
419: of QN and at ever restart.
421: .seealso: `SNES`, `SNESQN`, `SNESLineSearch`, `SNESQNScaleType`, `SNESSetJacobian()`
422: @*/
424: PetscErrorCode SNESQNSetScaleType(SNES snes, SNESQNScaleType stype)
425: {
427: PetscTryMethod(snes, "SNESQNSetScaleType_C", (SNES, SNESQNScaleType), (snes, stype));
428: return 0;
429: }
431: PetscErrorCode SNESQNSetScaleType_QN(SNES snes, SNESQNScaleType stype)
432: {
433: SNES_QN *qn = (SNES_QN *)snes->data;
435: qn->scale_type = stype;
436: if (stype == SNES_QN_SCALE_JACOBIAN) snes->usesksp = PETSC_TRUE;
437: return 0;
438: }
440: PetscErrorCode SNESQNSetRestartType_QN(SNES snes, SNESQNRestartType rtype)
441: {
442: SNES_QN *qn = (SNES_QN *)snes->data;
444: qn->restart_type = rtype;
445: return 0;
446: }
448: /*@
449: SNESQNSetType - Sets the quasi-Newton variant to be used in `SNESQN`.
451: Logically Collective
453: Input Parameters:
454: + snes - the iterative context
455: - qtype - variant type
457: Options Database Key:
458: . -snes_qn_type <lbfgs,broyden,badbroyden> - quasi-Newton type
460: Level: beginner
462: `SNESQNType`s:
463: + `SNES_QN_LBFGS` - LBFGS variant
464: . `SNES_QN_BROYDEN` - Broyden variant
465: - `SNES_QN_BADBROYDEN` - Bad Broyden variant
467: .seealso: `SNESQN`, `SNES_QN_LBFGS`, `SNES_QN_BROYDEN`, `SNES_QN_BADBROYDEN`, `SNESQNType`, `TAOLMVM`, `TAOBLMVM`
468: @*/
470: PetscErrorCode SNESQNSetType(SNES snes, SNESQNType qtype)
471: {
473: PetscTryMethod(snes, "SNESQNSetType_C", (SNES, SNESQNType), (snes, qtype));
474: return 0;
475: }
477: PetscErrorCode SNESQNSetType_QN(SNES snes, SNESQNType qtype)
478: {
479: SNES_QN *qn = (SNES_QN *)snes->data;
481: qn->type = qtype;
482: return 0;
483: }
485: /*MC
486: SNESQN - Limited-Memory Quasi-Newton methods for the solution of nonlinear systems.
488: Options Database Keys:
489: + -snes_qn_m <m> - Number of past states saved for the L-Broyden methods.
490: . -snes_qn_restart_type <powell,periodic,none> - set the restart type
491: . -snes_qn_powell_gamma - Angle condition for restart.
492: . -snes_qn_powell_descent - Descent condition for restart.
493: . -snes_qn_type <lbfgs,broyden,badbroyden> - QN type
494: . -snes_qn_scale_type <diagonal,none,scalar,jacobian> - scaling performed on inner Jacobian
495: . -snes_linesearch_type <cp, l2, basic> - Type of line search.
496: - -snes_qn_monitor - Monitors the quasi-newton Jacobian.
498: References:
499: + * - Kelley, C.T., Iterative Methods for Linear and Nonlinear Equations, Chapter 8, SIAM, 1995.
500: . * - R. Byrd, J. Nocedal, R. Schnabel, Representations of Quasi Newton Matrices and their use in Limited Memory Methods,
501: Technical Report, Northwestern University, June 1992.
502: . * - Peter N. Brown, Alan C. Hindmarsh, Homer F. Walker, Experiments with Quasi-Newton Methods in Solving Stiff ODE
503: Systems, SIAM J. Sci. Stat. Comput. Vol 6(2), April 1985.
504: . * - Peter R. Brune, Matthew G. Knepley, Barry F. Smith, and Xuemin Tu, "Composing Scalable Nonlinear Algebraic Solvers",
505: SIAM Review, 57(4), 2015
506: . * - Griewank, Andreas. "Broyden updating, the good and the bad!." Doc. Math (2012): 301-315.
507: . * - Gilbert, Jean Charles, and Claude Lemar{\'e}chal. "Some numerical experiments with variable-storage quasi-Newton algorithms."
508: Mathematical programming 45.1-3 (1989): 407-435.
509: - * - Dener A., Munson T. "Accelerating Limited-Memory Quasi-Newton Convergence for Large-Scale Optimization"
510: Computational Science - ICCS 2019. ICCS 2019. Lecture Notes in Computer Science, vol 11538. Springer, Cham
512: Level: beginner
514: Notes:
515: This implements the L-BFGS, Broyden, and "Bad" Broyden algorithms for the solution of F(x) = b using
516: previous change in F(x) and x to form the approximate inverse Jacobian using a series of multiplicative rank-one
517: updates.
519: When using a nonlinear preconditioner, one has two options as to how the preconditioner is applied. The first of
520: these options, sequential, uses the preconditioner to generate a new solution and function and uses those at this
521: iteration as the current iteration's values when constructing the approximate Jacobian. The second, composed,
522: perturbs the problem the Jacobian represents to be P(x, b) - x = 0, where P(x, b) is the preconditioner.
524: Uses left nonlinear preconditioning by default.
526: .seealso: `SNESQNRestartType`, `SNESQNSetRestartType()`, `SNESCreate()`, `SNES`, `SNESSetType()`, `SNESNEWTONLS`, `SNESNEWTONTR`,
527: `SNESQNScaleType`, `SNESQNSetScaleType()`, `SNESQNSetType`, `SNESQNSetType ()`
528: M*/
529: PETSC_EXTERN PetscErrorCode SNESCreate_QN(SNES snes)
530: {
531: SNES_QN *qn;
532: const char *optionsprefix;
534: snes->ops->setup = SNESSetUp_QN;
535: snes->ops->solve = SNESSolve_QN;
536: snes->ops->destroy = SNESDestroy_QN;
537: snes->ops->setfromoptions = SNESSetFromOptions_QN;
538: snes->ops->view = SNESView_QN;
539: snes->ops->reset = SNESReset_QN;
541: snes->npcside = PC_LEFT;
543: snes->usesnpc = PETSC_TRUE;
544: snes->usesksp = PETSC_FALSE;
546: snes->alwayscomputesfinalresidual = PETSC_TRUE;
548: if (!snes->tolerancesset) {
549: snes->max_funcs = 30000;
550: snes->max_its = 10000;
551: }
553: PetscNew(&qn);
554: snes->data = (void *)qn;
555: qn->m = 10;
556: qn->scaling = 1.0;
557: qn->monitor = NULL;
558: qn->monflg = PETSC_FALSE;
559: qn->powell_gamma = 0.9999;
560: qn->scale_type = SNES_QN_SCALE_DEFAULT;
561: qn->restart_type = SNES_QN_RESTART_DEFAULT;
562: qn->type = SNES_QN_LBFGS;
564: MatCreate(PetscObjectComm((PetscObject)snes), &qn->B);
565: SNESGetOptionsPrefix(snes, &optionsprefix);
566: MatSetOptionsPrefix(qn->B, optionsprefix);
568: PetscObjectComposeFunction((PetscObject)snes, "SNESQNSetScaleType_C", SNESQNSetScaleType_QN);
569: PetscObjectComposeFunction((PetscObject)snes, "SNESQNSetRestartType_C", SNESQNSetRestartType_QN);
570: PetscObjectComposeFunction((PetscObject)snes, "SNESQNSetType_C", SNESQNSetType_QN);
571: return 0;
572: }