Actual source code: qn.c

  1: #include <petsc/private/snesimpl.h>
  2: #include <petscdm.h>

  4: const char *const SNESQNScaleTypes[]   = {"DEFAULT", "NONE", "SCALAR", "DIAGONAL", "JACOBIAN", "SNESQNScaleType", "SNES_QN_SCALING_", NULL};
  5: const char *const SNESQNRestartTypes[] = {"DEFAULT", "NONE", "POWELL", "PERIODIC", "SNESQNRestartType", "SNES_QN_RESTART_", NULL};
  6: const char *const SNESQNTypes[]        = {"LBFGS", "BROYDEN", "BADBROYDEN", "SNESQNType", "SNES_QN_", NULL};

  8: typedef struct {
  9:   Mat               B;      /* Quasi-Newton approximation Matrix (MATLMVM) */
 10:   PetscInt          m;      /* The number of kept previous steps */
 11:   PetscReal        *lambda; /* The line search history of the method */
 12:   PetscBool         monflg;
 13:   PetscViewer       monitor;
 14:   PetscReal         powell_gamma; /* Powell angle restart condition */
 15:   PetscReal         scaling;      /* scaling of H0 */
 16:   SNESQNType        type;         /* the type of quasi-newton method used */
 17:   SNESQNScaleType   scale_type;   /* the type of scaling used */
 18:   SNESQNRestartType restart_type; /* determine the frequency and type of restart conditions */
 19: } SNES_QN;

 21: static PetscErrorCode SNESQNGetMatrix_Private(SNES snes, Mat *B)
 22: {
 23:   SNES_QN    *qn = (SNES_QN *)snes->data;
 24:   const char *optionsprefix;

 26:   PetscFunctionBegin;
 27:   if (!qn->B) {
 28:     PetscCall(MatCreate(PetscObjectComm((PetscObject)snes), &qn->B));
 29:     PetscCall(SNESGetOptionsPrefix(snes, &optionsprefix));
 30:     PetscCall(MatSetOptionsPrefix(qn->B, optionsprefix));
 31:     PetscCall(MatAppendOptionsPrefix(qn->B, "qn_"));
 32:     switch (qn->type) {
 33:     case SNES_QN_BROYDEN:
 34:       PetscCall(MatSetType(qn->B, MATLMVMBROYDEN));
 35:       qn->scale_type = SNES_QN_SCALE_NONE;
 36:       break;
 37:     case SNES_QN_BADBROYDEN:
 38:       PetscCall(MatSetType(qn->B, MATLMVMBADBROYDEN));
 39:       qn->scale_type = SNES_QN_SCALE_NONE;
 40:       break;
 41:     default:
 42:       PetscCall(MatSetType(qn->B, MATLMVMBFGS));
 43:       switch (qn->scale_type) {
 44:       case SNES_QN_SCALE_NONE:
 45:       case SNES_QN_SCALE_JACOBIAN:
 46:         PetscCall(MatLMVMSymBroydenSetScaleType(qn->B, MAT_LMVM_SYMBROYDEN_SCALE_NONE));
 47:         break;
 48:       case SNES_QN_SCALE_SCALAR:
 49:       case SNES_QN_SCALE_DEFAULT:
 50:         PetscCall(MatLMVMSymBroydenSetScaleType(qn->B, MAT_LMVM_SYMBROYDEN_SCALE_SCALAR));
 51:         break;
 52:       default:
 53:         break;
 54:       }
 55:       break;
 56:     }
 57:   }
 58:   *B = qn->B;
 59:   PetscFunctionReturn(PETSC_SUCCESS);
 60: }

 62: static PetscErrorCode SNESSolve_QN(SNES snes)
 63: {
 64:   SNES_QN             *qn = (SNES_QN *)snes->data;
 65:   Vec                  X, F, W;
 66:   Vec                  Y, D, Dold;
 67:   PetscInt             i, i_r;
 68:   PetscReal            fnorm, xnorm, ynorm;
 69:   SNESLineSearchReason lsreason;
 70:   PetscBool            powell, periodic, restart;
 71:   PetscScalar          DolddotD, DolddotDold;
 72:   SNESConvergedReason  reason;
 73: #if defined(PETSC_USE_INFO)
 74:   PetscReal gnorm;
 75: #endif

 77:   PetscFunctionBegin;
 78:   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);

 80:   PetscCall(PetscCitationsRegister(SNESCitation, &SNEScite));
 81:   F = snes->vec_func;       /* residual vector */
 82:   Y = snes->vec_sol_update; /* search direction generated by J^-1D*/
 83:   X = snes->vec_sol;        /* solution vector */

 85:   /* directions */
 86:   W    = snes->work[0];
 87:   D    = snes->work[1];
 88:   Dold = snes->work[2];

 90:   snes->reason = SNES_CONVERGED_ITERATING;

 92:   PetscCall(PetscObjectSAWsTakeAccess((PetscObject)snes));
 93:   snes->iter = 0;
 94:   snes->norm = 0.;
 95:   PetscCall(PetscObjectSAWsGrantAccess((PetscObject)snes));

 97:   if (snes->npc && snes->npcside == PC_LEFT && snes->functype == SNES_FUNCTION_PRECONDITIONED) {
 98:     /* we need to sample the preconditioned function only once here,
 99:        since linesearch will always use the preconditioned function */
100:     PetscCall(SNESApplyNPC(snes, X, NULL, F));
101:     PetscCall(SNESGetConvergedReason(snes->npc, &reason));
102:     if (reason < 0 && reason != SNES_DIVERGED_MAX_IT && reason != SNES_DIVERGED_TR_DELTA) {
103:       snes->reason = SNES_DIVERGED_INNER;
104:       PetscFunctionReturn(PETSC_SUCCESS);
105:     }
106:   } else {
107:     if (!snes->vec_func_init_set) PetscCall(SNESComputeFunction(snes, X, F));
108:     else snes->vec_func_init_set = PETSC_FALSE;
109:   }
110:   PetscCall(VecNorm(F, NORM_2, &fnorm));
111:   SNESCheckFunctionDomainError(snes, fnorm);

113:   PetscCall(PetscObjectSAWsTakeAccess((PetscObject)snes));
114:   snes->norm = fnorm;
115:   PetscCall(PetscObjectSAWsGrantAccess((PetscObject)snes));
116:   PetscCall(SNESLogConvergenceHistory(snes, fnorm, 0));
117:   PetscCall(SNESMonitor(snes, 0, fnorm));

119:   /* test convergence */
120:   PetscCall(SNESConverged(snes, 0, 0.0, 0.0, fnorm));
121:   if (snes->reason) PetscFunctionReturn(PETSC_SUCCESS);

123:   restart = PETSC_TRUE;
124:   for (i = 0, i_r = 0; i < snes->max_its; i++, i_r++) {
125:     PetscScalar gdx;

127:     /* general purpose update */
128:     PetscTryTypeMethod(snes, update, snes->iter);

130:     if (snes->npc) {
131:       if (snes->npcside == PC_LEFT && snes->functype == SNES_FUNCTION_UNPRECONDITIONED) {
132:         PetscCall(SNESApplyNPC(snes, X, F, D));
133:         PetscCall(SNESGetConvergedReason(snes->npc, &reason));
134:         if (reason < 0 && reason != SNES_DIVERGED_MAX_IT && reason != SNES_DIVERGED_TR_DELTA) {
135:           snes->reason = SNES_DIVERGED_INNER;
136:           PetscFunctionReturn(PETSC_SUCCESS);
137:         }
138:       } else if (snes->npcside == PC_RIGHT) {
139:         PetscCall(PetscLogEventBegin(SNES_NPCSolve, snes->npc, X, 0, 0));
140:         PetscCall(SNESSolve(snes->npc, snes->vec_rhs, X));
141:         PetscCall(PetscLogEventEnd(SNES_NPCSolve, snes->npc, X, 0, 0));
142:         PetscCall(SNESGetConvergedReason(snes->npc, &reason));
143:         if (reason < 0 && reason != SNES_DIVERGED_MAX_IT && reason != SNES_DIVERGED_TR_DELTA) {
144:           snes->reason = SNES_DIVERGED_INNER;
145:           PetscFunctionReturn(PETSC_SUCCESS);
146:         }
147:         PetscCall(SNESGetNPCFunction(snes, F, &fnorm));
148:         PetscCall(VecCopy(F, D));
149:       } else {
150:         PetscCall(VecCopy(F, D));
151:       }
152:     } else {
153:       PetscCall(VecCopy(F, D));
154:     }

156:     /* scale the initial update */
157:     if (qn->scale_type == SNES_QN_SCALE_JACOBIAN && restart) {
158:       PetscCall(SNESComputeJacobian(snes, X, snes->jacobian, snes->jacobian_pre));
159:       SNESCheckJacobianDomainError(snes);
160:       PetscCall(KSPSetOperators(snes->ksp, snes->jacobian, snes->jacobian_pre));
161:       PetscCall(MatLMVMSetJ0KSP(qn->B, snes->ksp));
162:     }

164:     /* update QN approx and calculate step */
165:     PetscCall(MatLMVMUpdate(qn->B, X, D));
166:     PetscCall(MatSolve(qn->B, D, Y));
167:     PetscCall(VecDot(F, Y, &gdx));
168:     if (PetscAbsScalar(gdx) <= 0) {
169:       /* Step is not descent or solve was not successful
170:          Use steepest descent direction */
171:       PetscCall(VecCopy(F, Y));
172:       PetscCall(MatLMVMReset(qn->B, PETSC_FALSE));
173:     }

175:     /* line search for lambda */
176:     ynorm = 1;
177: #if defined(PETSC_USE_INFO)
178:     gnorm = fnorm;
179: #endif
180:     PetscCall(VecCopy(D, Dold));
181:     PetscCall(SNESLineSearchApply(snes->linesearch, X, F, &fnorm, Y));
182:     if (snes->reason) break;
183:     SNESCheckLineSearchFailure(snes);
184:     PetscCall(SNESLineSearchGetReason(snes->linesearch, &lsreason));

186:     PetscCall(SNESLineSearchGetNorms(snes->linesearch, &xnorm, &fnorm, &ynorm));
187:     /* convergence monitoring */
188:     PetscCall(PetscInfo(snes, "fnorm=%18.16e, gnorm=%18.16e, ynorm=%18.16e, lssucceed=%d\n", (double)fnorm, (double)gnorm, (double)ynorm, (int)lsreason));

190:     PetscCall(PetscObjectSAWsTakeAccess((PetscObject)snes));
191:     snes->iter  = i + 1;
192:     snes->norm  = fnorm;
193:     snes->xnorm = xnorm;
194:     snes->ynorm = ynorm;
195:     PetscCall(PetscObjectSAWsGrantAccess((PetscObject)snes));

197:     PetscCall(SNESLogConvergenceHistory(snes, snes->norm, snes->iter));

199:     /* set parameter for default relative tolerance convergence test */
200:     PetscCall(SNESConverged(snes, snes->iter, xnorm, ynorm, fnorm));
201:     PetscCall(SNESMonitor(snes, snes->iter, snes->norm));
202:     if (snes->reason) PetscFunctionReturn(PETSC_SUCCESS);

204:     if (snes->npc && snes->npcside == PC_LEFT && snes->functype == SNES_FUNCTION_UNPRECONDITIONED) {
205:       PetscCall(SNESApplyNPC(snes, X, F, D));
206:       PetscCall(SNESGetConvergedReason(snes->npc, &reason));
207:       if (reason < 0 && reason != SNES_DIVERGED_MAX_IT && reason != SNES_DIVERGED_TR_DELTA) {
208:         snes->reason = SNES_DIVERGED_INNER;
209:         PetscFunctionReturn(PETSC_SUCCESS);
210:       }
211:     } else {
212:       PetscCall(VecCopy(F, D));
213:     }

215:     /* restart conditions */
216:     powell = PETSC_FALSE;
217:     if (qn->restart_type == SNES_QN_RESTART_POWELL && i_r > 1) {
218:       /* check restart by Powell's Criterion: |F^T H_0 Fold| > powell_gamma * |Fold^T H_0 Fold| */
219:       if (qn->scale_type == SNES_QN_SCALE_JACOBIAN) PetscCall(MatMult(snes->jacobian, Dold, W));
220:       else PetscCall(VecCopy(Dold, W));
221:       PetscCall(VecDotBegin(W, Dold, &DolddotDold));
222:       PetscCall(VecDotBegin(W, D, &DolddotD));
223:       PetscCall(VecDotEnd(W, Dold, &DolddotDold));
224:       PetscCall(VecDotEnd(W, D, &DolddotD));
225:       if (PetscAbs(PetscRealPart(DolddotD)) > qn->powell_gamma * PetscAbs(PetscRealPart(DolddotDold))) powell = PETSC_TRUE;
226:     }
227:     periodic = PETSC_FALSE;
228:     if (qn->restart_type == SNES_QN_RESTART_PERIODIC) {
229:       if (i_r > qn->m - 1) periodic = PETSC_TRUE;
230:     }

232:     /* restart if either powell or periodic restart is satisfied. */
233:     restart = PETSC_FALSE;
234:     if (lsreason || powell || periodic) {
235:       restart = PETSC_TRUE;
236:       if (qn->monflg) {
237:         PetscCall(PetscViewerASCIIAddTab(qn->monitor, ((PetscObject)snes)->tablevel + 2));
238:         if (powell) {
239:           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));
240:         } else {
241:           PetscCall(PetscViewerASCIIPrintf(qn->monitor, "Periodic restart! i_r = %" PetscInt_FMT "\n", i_r));
242:         }
243:         PetscCall(PetscViewerASCIISubtractTab(qn->monitor, ((PetscObject)snes)->tablevel + 2));
244:       }
245:       i_r = -1;
246:       PetscCall(MatLMVMReset(qn->B, PETSC_FALSE));
247:     }
248:   }
249:   PetscFunctionReturn(PETSC_SUCCESS);
250: }

252: static PetscErrorCode SNESSetUp_QN(SNES snes)
253: {
254:   SNES_QN *qn = (SNES_QN *)snes->data;
255:   DM       dm;
256:   PetscInt n, N;

258:   PetscFunctionBegin;
259:   if (!snes->vec_sol) {
260:     PetscCall(SNESGetDM(snes, &dm));
261:     PetscCall(DMCreateGlobalVector(dm, &snes->vec_sol));
262:   }
263:   PetscCall(SNESSetWorkVecs(snes, 3));
264:   PetscCall(SNESQNGetMatrix_Private(snes, &qn->B));

266:   if (qn->scale_type == SNES_QN_SCALE_JACOBIAN) PetscCall(SNESSetUpMatrices(snes));
267:   if (snes->npcside == PC_LEFT && snes->functype == SNES_FUNCTION_DEFAULT) snes->functype = SNES_FUNCTION_PRECONDITIONED;

269:   /* set method defaults */
270:   if (qn->scale_type == SNES_QN_SCALE_DEFAULT) {
271:     if (qn->type == SNES_QN_BADBROYDEN) {
272:       qn->scale_type = SNES_QN_SCALE_NONE;
273:     } else {
274:       qn->scale_type = SNES_QN_SCALE_SCALAR;
275:     }
276:   }
277:   if (qn->restart_type == SNES_QN_RESTART_DEFAULT) {
278:     if (qn->type == SNES_QN_LBFGS) {
279:       qn->restart_type = SNES_QN_RESTART_POWELL;
280:     } else {
281:       qn->restart_type = SNES_QN_RESTART_PERIODIC;
282:     }
283:   }
284:   /* Set up the LMVM matrix */
285:   PetscCall(SNESQNGetMatrix_Private(snes, &qn->B));
286:   PetscCall(VecGetLocalSize(snes->vec_sol, &n));
287:   PetscCall(VecGetSize(snes->vec_sol, &N));
288:   PetscCall(MatSetSizes(qn->B, n, n, N, N));
289:   PetscCall(MatSetUp(qn->B));
290:   PetscCall(MatLMVMReset(qn->B, PETSC_TRUE));
291:   PetscCall(MatLMVMSetHistorySize(qn->B, qn->m));
292:   PetscCall(MatLMVMAllocate(qn->B, snes->vec_sol, snes->vec_func));
293:   PetscFunctionReturn(PETSC_SUCCESS);
294: }

296: static PetscErrorCode SNESReset_QN(SNES snes)
297: {
298:   SNES_QN *qn = (SNES_QN *)snes->data;

300:   PetscFunctionBegin;
301:   PetscCall(MatDestroy(&qn->B));
302:   PetscFunctionReturn(PETSC_SUCCESS);
303: }

305: static PetscErrorCode SNESDestroy_QN(SNES snes)
306: {
307:   PetscFunctionBegin;
308:   PetscCall(SNESReset_QN(snes));
309:   PetscCall(PetscFree(snes->data));
310:   PetscCall(PetscObjectComposeFunction((PetscObject)snes, "SNESQNSetScaleType_C", NULL));
311:   PetscCall(PetscObjectComposeFunction((PetscObject)snes, "SNESQNSetRestartType_C", NULL));
312:   PetscCall(PetscObjectComposeFunction((PetscObject)snes, "SNESQNSetType_C", NULL));
313:   PetscFunctionReturn(PETSC_SUCCESS);
314: }

316: static PetscErrorCode SNESSetFromOptions_QN(SNES snes, PetscOptionItems PetscOptionsObject)
317: {
318:   SNES_QN          *qn = (SNES_QN *)snes->data;
319:   PetscBool         flg;
320:   SNESLineSearch    linesearch;
321:   SNESQNRestartType rtype = qn->restart_type;
322:   SNESQNScaleType   stype = qn->scale_type;
323:   SNESQNType        qtype = qn->type;

325:   PetscFunctionBegin;
326:   PetscOptionsHeadBegin(PetscOptionsObject, "SNES QN options");
327:   PetscCall(PetscOptionsInt("-snes_qn_m", "Number of past states saved for L-BFGS methods", "SNESQN", qn->m, &qn->m, NULL));
328:   PetscCall(PetscOptionsReal("-snes_qn_powell_gamma", "Powell angle tolerance", "SNESQN", qn->powell_gamma, &qn->powell_gamma, NULL));
329:   PetscCall(PetscOptionsBool("-snes_qn_monitor", "Monitor for the QN methods", "SNESQN", qn->monflg, &qn->monflg, NULL));
330:   PetscCall(PetscOptionsEnum("-snes_qn_scale_type", "Scaling type", "SNESQNSetScaleType", SNESQNScaleTypes, (PetscEnum)stype, (PetscEnum *)&stype, &flg));
331:   if (flg) PetscCall(SNESQNSetScaleType(snes, stype));

333:   PetscCall(PetscOptionsEnum("-snes_qn_restart_type", "Restart type", "SNESQNSetRestartType", SNESQNRestartTypes, (PetscEnum)rtype, (PetscEnum *)&rtype, &flg));
334:   if (flg) PetscCall(SNESQNSetRestartType(snes, rtype));

336:   PetscCall(PetscOptionsEnum("-snes_qn_type", "Quasi-Newton update type", "", SNESQNTypes, (PetscEnum)qtype, (PetscEnum *)&qtype, &flg));
337:   if (flg) PetscCall(SNESQNSetType(snes, qtype));
338:   PetscOptionsHeadEnd();
339:   PetscCall(SNESQNGetMatrix_Private(snes, &qn->B));
340:   PetscCall(MatSetFromOptions(qn->B));
341:   if (!snes->linesearch) {
342:     PetscCall(SNESGetLineSearch(snes, &linesearch));
343:     if (!((PetscObject)linesearch)->type_name) {
344:       if (qn->type == SNES_QN_LBFGS) {
345:         PetscCall(SNESLineSearchSetType(linesearch, SNESLINESEARCHCP));
346:       } else if (qn->type == SNES_QN_BROYDEN) {
347:         PetscCall(SNESLineSearchSetType(linesearch, SNESLINESEARCHBASIC));
348:       } else {
349:         PetscCall(SNESLineSearchSetType(linesearch, SNESLINESEARCHSECANT));
350:       }
351:     }
352:   }
353:   if (qn->monflg) PetscCall(PetscViewerASCIIGetStdout(PetscObjectComm((PetscObject)snes), &qn->monitor));
354:   PetscFunctionReturn(PETSC_SUCCESS);
355: }

357: static PetscErrorCode SNESView_QN(SNES snes, PetscViewer viewer)
358: {
359:   SNES_QN  *qn = (SNES_QN *)snes->data;
360:   PetscBool isascii;

362:   PetscFunctionBegin;
363:   PetscCall(PetscObjectTypeCompare((PetscObject)viewer, PETSCVIEWERASCII, &isascii));
364:   if (isascii) {
365:     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]));
366:     PetscCall(PetscViewerASCIIPrintf(viewer, "  Stored subspace size: %" PetscInt_FMT "\n", qn->m));
367:     PetscCall(MatView(qn->B, viewer));
368:   }
369:   PetscFunctionReturn(PETSC_SUCCESS);
370: }

372: /*@
373:   SNESQNSetRestartType - Sets the restart type for `SNESQN`.

375:   Logically Collective

377:   Input Parameters:
378: + snes  - the iterative context
379: - rtype - restart type, see `SNESQNRestartType`

381:   Options Database Keys:
382: + -snes_qn_restart_type (powell|periodic|none) - set the restart type
383: - -snes_qn_m m                                 - sets the number of stored updates and the restart period for periodic

385:   Level: intermediate

387: .seealso: [](ch_snes), `SNES`, `SNESQN`, `SNESQNRestartType`, `SNES_QN_RESTART_NONE`, `SNES_QN_RESTART_POWELL`, `SNES_QN_RESTART_PERIODIC`,
388:           `SNESQNType`, `SNESQNScaleType`
389: @*/
390: PetscErrorCode SNESQNSetRestartType(SNES snes, SNESQNRestartType rtype)
391: {
392:   PetscFunctionBegin;
394:   PetscTryMethod(snes, "SNESQNSetRestartType_C", (SNES, SNESQNRestartType), (snes, rtype));
395:   PetscFunctionReturn(PETSC_SUCCESS);
396: }

398: /*@
399:   SNESQNSetScaleType - Sets the scaling type for the inner inverse Jacobian in `SNESQN`.

401:   Logically Collective

403:   Input Parameters:
404: + snes  - the nonlinear solver context
405: - stype - scale type, see `SNESQNScaleType`

407:   Options Database Key:
408: . -snes_qn_scale_type (diagonal|none|scalar|jacobian) - Scaling type

410:   Level: intermediate

412: .seealso: [](ch_snes), `SNES`, `SNESQN`, `SNESLineSearch`, `SNESQNScaleType`, `SNESSetJacobian()`, `SNESQNType`, `SNESQNRestartType`
413: @*/
414: PetscErrorCode SNESQNSetScaleType(SNES snes, SNESQNScaleType stype)
415: {
416:   PetscFunctionBegin;
418:   PetscTryMethod(snes, "SNESQNSetScaleType_C", (SNES, SNESQNScaleType), (snes, stype));
419:   PetscFunctionReturn(PETSC_SUCCESS);
420: }

422: static PetscErrorCode SNESQNSetScaleType_QN(SNES snes, SNESQNScaleType stype)
423: {
424:   SNES_QN *qn = (SNES_QN *)snes->data;

426:   PetscFunctionBegin;
427:   qn->scale_type = stype;
428:   if (stype == SNES_QN_SCALE_JACOBIAN) snes->usesksp = PETSC_TRUE;
429:   PetscFunctionReturn(PETSC_SUCCESS);
430: }

432: static PetscErrorCode SNESQNSetRestartType_QN(SNES snes, SNESQNRestartType rtype)
433: {
434:   SNES_QN *qn = (SNES_QN *)snes->data;

436:   PetscFunctionBegin;
437:   qn->restart_type = rtype;
438:   PetscFunctionReturn(PETSC_SUCCESS);
439: }

441: /*@
442:   SNESQNSetType - Sets the quasi-Newton variant to be used in `SNESQN`.

444:   Logically Collective

446:   Input Parameters:
447: + snes  - the iterative context
448: - qtype - variant type, see `SNESQNType`

450:   Options Database Key:
451: . -snes_qn_type (lbfgs|broyden|badbroyden) - quasi-Newton type

453:   Level: intermediate

455: .seealso: [](ch_snes), `SNESQN`, `SNES_QN_LBFGS`, `SNES_QN_BROYDEN`, `SNES_QN_BADBROYDEN`, `SNESQNType`, `SNESQNScaleType`, `TAOLMVM`, `TAOBLMVM`
456: @*/
457: PetscErrorCode SNESQNSetType(SNES snes, SNESQNType qtype)
458: {
459:   PetscFunctionBegin;
461:   PetscTryMethod(snes, "SNESQNSetType_C", (SNES, SNESQNType), (snes, qtype));
462:   PetscFunctionReturn(PETSC_SUCCESS);
463: }

465: static PetscErrorCode SNESQNSetType_QN(SNES snes, SNESQNType qtype)
466: {
467:   SNES_QN *qn = (SNES_QN *)snes->data;

469:   PetscFunctionBegin;
470:   qn->type = qtype;
471:   PetscFunctionReturn(PETSC_SUCCESS);
472: }

474: /*MC
475:       SNESQN - Limited-Memory Quasi-Newton methods for the solution of nonlinear systems.

477:       Options Database Keys:
478: +     -snes_qn_m m                                        - Number of past states saved for the L-Broyden methods.
479: .     -snes_qn_restart_type (powell|periodic|none)        - set the restart type
480: .     -snes_qn_powell_gamma gamma                         - Angle condition for restart.
481: .     -snes_qn_powell_descent                             - Descent condition for restart.
482: .     -snes_qn_type (lbfgs|broyden|badbroyden)            - QN type
483: .     -snes_qn_scale_type (diagonal|none|scalar|jacobian) - scaling performed on inner Jacobian
484: .     -snes_linesearch_type (cp|l2|basic)                 - Type of line search.
485: -     -snes_qn_monitor                                    - Monitors the quasi-newton Jacobian.

487:       Level: beginner

489:       Notes:
490:       This implements the L-BFGS, Broyden, and "Bad" Broyden algorithms for the solution of $F(x) = b$ using
491:       previous change in $F(x)$ and $x$ to form the approximate inverse Jacobian using a series of multiplicative rank-one
492:       updates.

494:       When using a nonlinear preconditioner, one has two options as to how the preconditioner is applied.  The first of
495:       these options, sequential, uses the preconditioner to generate a new solution and function and uses those at this
496:       iteration as the current iteration's values when constructing the approximate Jacobian.  The second, composed,
497:       perturbs the problem the Jacobian represents to be $P(x, b) - x = 0 $, where $P(x, b)$ is the preconditioner.

499:       Uses left nonlinear preconditioning by default.

501:       See {cite}`byrd1994representations`, {cite}`kelley95`, {cite}`brown1985experiments`, {cite}`griewank2012broyden`, {cite}`gilbert1989some`,
502:       {cite}`dener2019accelerating`, and {cite}`bruneknepleysmithtu15`

504: .seealso: [](ch_snes), `SNESQNRestartType`, `SNESQNSetRestartType()`, `SNESCreate()`, `SNES`, `SNESSetType()`, `SNESNEWTONLS`, `SNESNEWTONTR`,
505:           `SNESQNScaleType`, `SNESQNSetScaleType()`, `SNESQNType`, `SNESQNSetType()`
506: M*/
507: PETSC_EXTERN PetscErrorCode SNESCreate_QN(SNES snes)
508: {
509:   SNES_QN *qn;

511:   PetscFunctionBegin;
512:   snes->ops->setup          = SNESSetUp_QN;
513:   snes->ops->solve          = SNESSolve_QN;
514:   snes->ops->destroy        = SNESDestroy_QN;
515:   snes->ops->setfromoptions = SNESSetFromOptions_QN;
516:   snes->ops->view           = SNESView_QN;
517:   snes->ops->reset          = SNESReset_QN;

519:   snes->npcside = PC_LEFT;

521:   snes->usesnpc = PETSC_TRUE;
522:   snes->usesksp = PETSC_FALSE;

524:   snes->alwayscomputesfinalresidual = PETSC_TRUE;

526:   PetscCall(SNESParametersInitialize(snes));
527:   PetscObjectParameterSetDefault(snes, max_funcs, 30000);
528:   PetscObjectParameterSetDefault(snes, max_its, 10000);

530:   PetscCall(PetscNew(&qn));
531:   snes->data       = (void *)qn;
532:   qn->m            = 10;
533:   qn->scaling      = 1.0;
534:   qn->monitor      = NULL;
535:   qn->monflg       = PETSC_FALSE;
536:   qn->powell_gamma = 0.9999;
537:   qn->scale_type   = SNES_QN_SCALE_DEFAULT;
538:   qn->restart_type = SNES_QN_RESTART_DEFAULT;
539:   qn->type         = SNES_QN_LBFGS;

541:   PetscCall(PetscObjectComposeFunction((PetscObject)snes, "SNESQNSetScaleType_C", SNESQNSetScaleType_QN));
542:   PetscCall(PetscObjectComposeFunction((PetscObject)snes, "SNESQNSetRestartType_C", SNESQNSetRestartType_QN));
543:   PetscCall(PetscObjectComposeFunction((PetscObject)snes, "SNESQNSetType_C", SNESQNSetType_QN));
544:   PetscFunctionReturn(PETSC_SUCCESS);
545: }