1 #include <../src/snes/impls/ngmres/snesngmres.h> /*I "petscsnes.h" I*/
2 #include <petscblaslapack.h>
3 #include <petscdm.h>
4
5 const char *const SNESNGMRESRestartTypes[] = {"NONE", "PERIODIC", "DIFFERENCE", "SNESNGMRESRestartType", "SNES_NGMRES_RESTART_", NULL};
6 const char *const SNESNGMRESSelectTypes[] = {"NONE", "DIFFERENCE", "LINESEARCH", "SNESNGMRESSelectType", "SNES_NGMRES_SELECT_", NULL};
7
SNESReset_NGMRES(SNES snes)8 PetscErrorCode SNESReset_NGMRES(SNES snes)
9 {
10 SNES_NGMRES *ngmres = (SNES_NGMRES *)snes->data;
11
12 PetscFunctionBegin;
13 PetscCall(VecDestroyVecs(ngmres->msize, &ngmres->Fdot));
14 PetscCall(VecDestroyVecs(ngmres->msize, &ngmres->Xdot));
15 PetscCall(SNESLineSearchDestroy(&ngmres->additive_linesearch));
16 PetscFunctionReturn(PETSC_SUCCESS);
17 }
18
SNESDestroy_NGMRES(SNES snes)19 PetscErrorCode SNESDestroy_NGMRES(SNES snes)
20 {
21 SNES_NGMRES *ngmres = (SNES_NGMRES *)snes->data;
22
23 PetscFunctionBegin;
24 PetscCall(SNESReset_NGMRES(snes));
25 PetscCall(PetscFree4(ngmres->h, ngmres->beta, ngmres->xi, ngmres->q));
26 PetscCall(PetscFree3(ngmres->xnorms, ngmres->fnorms, ngmres->s));
27 #if defined(PETSC_USE_COMPLEX)
28 PetscCall(PetscFree(ngmres->rwork));
29 #endif
30 PetscCall(PetscFree(ngmres->work));
31 PetscCall(PetscObjectComposeFunction((PetscObject)snes, "SNESNGMRESSetSelectType_C", NULL));
32 PetscCall(PetscObjectComposeFunction((PetscObject)snes, "SNESNGMRESSetRestartType_C", NULL));
33 PetscCall(PetscObjectComposeFunction((PetscObject)snes, "SNESNGMRESSetRestartFmRise_C", NULL));
34 PetscCall(PetscObjectComposeFunction((PetscObject)snes, "SNESNGMRESGetRestartFmRise_C", NULL));
35 PetscCall(PetscFree(snes->data));
36 PetscFunctionReturn(PETSC_SUCCESS);
37 }
38
SNESSetUp_NGMRES(SNES snes)39 PetscErrorCode SNESSetUp_NGMRES(SNES snes)
40 {
41 SNES_NGMRES *ngmres = (SNES_NGMRES *)snes->data;
42 PetscInt msize, hsize;
43 DM dm;
44
45 PetscFunctionBegin;
46 PetscCheck(!snes->npc || snes->npcside != PC_LEFT || snes->functype != SNES_FUNCTION_UNPRECONDITIONED, PetscObjectComm((PetscObject)snes), PETSC_ERR_ARG_WRONGSTATE, "SNESNGMRES does not support left preconditioning with unpreconditioned function");
47 if (snes->npcside == PC_LEFT && snes->functype == SNES_FUNCTION_DEFAULT) snes->functype = SNES_FUNCTION_PRECONDITIONED;
48 PetscCall(SNESSetWorkVecs(snes, 5));
49
50 if (!snes->vec_sol) {
51 PetscCall(SNESGetDM(snes, &dm));
52 PetscCall(DMCreateGlobalVector(dm, &snes->vec_sol));
53 }
54
55 if (!ngmres->Xdot) PetscCall(VecDuplicateVecs(snes->vec_sol, ngmres->msize, &ngmres->Xdot));
56 if (!ngmres->Fdot) PetscCall(VecDuplicateVecs(snes->vec_sol, ngmres->msize, &ngmres->Fdot));
57 if (!ngmres->setup_called) {
58 msize = ngmres->msize; /* restart size */
59 hsize = msize * msize;
60
61 /* explicit least squares minimization solve */
62 PetscCall(PetscCalloc4(hsize, &ngmres->h, msize, &ngmres->beta, msize, &ngmres->xi, hsize, &ngmres->q));
63 PetscCall(PetscMalloc3(msize, &ngmres->xnorms, msize, &ngmres->fnorms, msize, &ngmres->s));
64 ngmres->nrhs = 1;
65 PetscCall(PetscBLASIntCast(msize, &ngmres->lda));
66 PetscCall(PetscBLASIntCast(msize, &ngmres->ldb));
67 PetscCall(PetscBLASIntCast(12 * msize, &ngmres->lwork));
68 #if defined(PETSC_USE_COMPLEX)
69 PetscCall(PetscMalloc1(ngmres->lwork, &ngmres->rwork));
70 #endif
71 PetscCall(PetscMalloc1(ngmres->lwork, &ngmres->work));
72 }
73
74 ngmres->setup_called = PETSC_TRUE;
75 PetscFunctionReturn(PETSC_SUCCESS);
76 }
77
SNESSetFromOptions_NGMRES(SNES snes,PetscOptionItems PetscOptionsObject)78 static PetscErrorCode SNESSetFromOptions_NGMRES(SNES snes, PetscOptionItems PetscOptionsObject)
79 {
80 SNES_NGMRES *ngmres = (SNES_NGMRES *)snes->data;
81 PetscBool debug = PETSC_FALSE;
82
83 PetscFunctionBegin;
84 PetscOptionsHeadBegin(PetscOptionsObject, "SNES NGMRES options");
85 PetscCall(PetscOptionsEnum("-snes_ngmres_select_type", "Select type", "SNESNGMRESSetSelectType", SNESNGMRESSelectTypes, (PetscEnum)ngmres->select_type, (PetscEnum *)&ngmres->select_type, NULL));
86 PetscCall(PetscOptionsEnum("-snes_ngmres_restart_type", "Restart type", "SNESNGMRESSetRestartType", SNESNGMRESRestartTypes, (PetscEnum)ngmres->restart_type, (PetscEnum *)&ngmres->restart_type, NULL));
87 PetscCall(PetscOptionsBool("-snes_ngmres_candidate", "Use candidate storage", "SNES", ngmres->candidate, &ngmres->candidate, NULL));
88 PetscCall(PetscOptionsBool("-snes_ngmres_approxfunc", "Linearly approximate the function", "SNES", ngmres->approxfunc, &ngmres->approxfunc, NULL));
89 PetscCall(PetscOptionsInt("-snes_ngmres_m", "Number of directions", "SNES", ngmres->msize, &ngmres->msize, NULL));
90 PetscCall(PetscOptionsInt("-snes_ngmres_restart", "Iterations before forced restart", "SNES", ngmres->restart_periodic, &ngmres->restart_periodic, NULL));
91 PetscCall(PetscOptionsInt("-snes_ngmres_restart_it", "Tolerance iterations before restart", "SNES", ngmres->restart_it, &ngmres->restart_it, NULL));
92 PetscCall(PetscOptionsBool("-snes_ngmres_monitor", "Monitor actions of NGMRES", "SNES", ngmres->monitor ? PETSC_TRUE : PETSC_FALSE, &debug, NULL));
93 if (debug) ngmres->monitor = PETSC_VIEWER_STDOUT_(PetscObjectComm((PetscObject)snes));
94 PetscCall(PetscOptionsReal("-snes_ngmres_gammaA", "Residual selection constant", "SNES", ngmres->gammaA, &ngmres->gammaA, NULL));
95 PetscCall(PetscOptionsReal("-snes_ngmres_gammaC", "Residual restart constant", "SNES", ngmres->gammaC, &ngmres->gammaC, NULL));
96 PetscCall(PetscOptionsReal("-snes_ngmres_epsilonB", "Difference selection constant", "SNES", ngmres->epsilonB, &ngmres->epsilonB, NULL));
97 PetscCall(PetscOptionsReal("-snes_ngmres_deltaB", "Difference residual selection constant", "SNES", ngmres->deltaB, &ngmres->deltaB, NULL));
98 PetscCall(PetscOptionsBool("-snes_ngmres_restart_fm_rise", "Restart on F_M residual rise", "SNESNGMRESSetRestartFmRise", ngmres->restart_fm_rise, &ngmres->restart_fm_rise, NULL));
99 PetscOptionsHeadEnd();
100 if (ngmres->gammaA > ngmres->gammaC && ngmres->gammaC > 2.) ngmres->gammaC = ngmres->gammaA;
101 if (ngmres->select_type == SNES_NGMRES_SELECT_LINESEARCH) {
102 PetscCall(SNESNGMRESGetAdditiveLineSearch_Private(snes, &ngmres->additive_linesearch));
103 PetscCall(SNESLineSearchSetFromOptions(ngmres->additive_linesearch));
104 }
105 PetscFunctionReturn(PETSC_SUCCESS);
106 }
107
SNESView_NGMRES(SNES snes,PetscViewer viewer)108 PetscErrorCode SNESView_NGMRES(SNES snes, PetscViewer viewer)
109 {
110 SNES_NGMRES *ngmres = (SNES_NGMRES *)snes->data;
111 PetscBool isascii;
112
113 PetscFunctionBegin;
114 PetscCall(PetscObjectTypeCompare((PetscObject)viewer, PETSCVIEWERASCII, &isascii));
115 if (isascii) {
116 PetscCall(PetscViewerASCIIPrintf(viewer, " Number of stored past updates: %" PetscInt_FMT "\n", ngmres->msize));
117 if (ngmres->select_type == SNES_NGMRES_SELECT_DIFFERENCE) {
118 PetscCall(PetscViewerASCIIPrintf(viewer, " Residual selection: gammaA=%1.0e, gammaC=%1.0e\n", (double)ngmres->gammaA, (double)ngmres->gammaC));
119 PetscCall(PetscViewerASCIIPrintf(viewer, " Difference restart: epsilonB=%1.0e, deltaB=%1.0e\n", (double)ngmres->epsilonB, (double)ngmres->deltaB));
120 PetscCall(PetscViewerASCIIPrintf(viewer, " Restart on F_M residual increase: %s\n", PetscBools[ngmres->restart_fm_rise]));
121 }
122 if (ngmres->additive_linesearch) {
123 PetscCall(PetscViewerASCIIPrintf(viewer, " Additive line-search details:\n"));
124 PetscCall(SNESLineSearchView(ngmres->additive_linesearch, viewer));
125 }
126 }
127 PetscFunctionReturn(PETSC_SUCCESS);
128 }
129
SNESSolve_NGMRES(SNES snes)130 static PetscErrorCode SNESSolve_NGMRES(SNES snes)
131 {
132 SNES_NGMRES *ngmres = (SNES_NGMRES *)snes->data;
133 Vec X, F, B, D, Y; /* present solution, residual, and preconditioned residual */
134 Vec XA, FA, XM, FM; /* candidate linear combination answers */
135 PetscReal fnorm, fMnorm, fAnorm; /* coefficients and RHS to the minimization problem */
136 PetscReal xnorm, xMnorm, xAnorm;
137 PetscReal ynorm, yMnorm, yAnorm;
138 PetscInt k, k_restart, l, ivec, restart_count = 0;
139 PetscReal objmin, objM, objA, obj; /* support for objective functions minimization */
140 PetscBool selectRestart; /* solution selection data */
141 SNESConvergedReason reason;
142 SNESObjectiveFn *objective;
143 /*
144 These two variables are initialized to prevent compilers/analyzers from producing false warnings about these variables being passed
145 to SNESNGMRESSelect_Private() without being set when SNES_NGMRES_RESTART_DIFFERENCE, the values are not used in the subroutines in that case
146 so the code is correct as written.
147 */
148 PetscReal dnorm = 0.0, dminnorm = 0.0;
149
150 PetscFunctionBegin;
151 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);
152
153 PetscCall(PetscCitationsRegister(SNESCitation, &SNEScite));
154 /* variable initialization */
155 snes->reason = SNES_CONVERGED_ITERATING;
156 X = snes->vec_sol;
157 F = snes->vec_func;
158 B = snes->vec_rhs;
159 XA = snes->work[2];
160 FA = snes->work[0];
161 D = snes->work[1];
162
163 /* work for the line search */
164 Y = snes->vec_sol_update;
165 XM = snes->work[3];
166 FM = snes->work[4];
167
168 PetscCall(PetscObjectSAWsTakeAccess((PetscObject)snes));
169 snes->iter = 0;
170 snes->norm = 0.;
171 PetscCall(PetscObjectSAWsGrantAccess((PetscObject)snes));
172
173 /* initialization */
174
175 if (snes->npc && snes->npcside == PC_LEFT) {
176 PetscCall(SNESApplyNPC(snes, X, NULL, F));
177 PetscCall(SNESGetConvergedReason(snes->npc, &reason));
178 if (reason < 0 && reason != SNES_DIVERGED_MAX_IT) {
179 snes->reason = SNES_DIVERGED_INNER;
180 PetscFunctionReturn(PETSC_SUCCESS);
181 }
182 } else {
183 if (!snes->vec_func_init_set) PetscCall(SNESComputeFunction(snes, X, F));
184 else snes->vec_func_init_set = PETSC_FALSE;
185 }
186 PetscCall(VecNorm(F, NORM_2, &fnorm));
187 SNESCheckFunctionDomainError(snes, fnorm);
188 PetscCall(SNESGetObjective(snes, &objective, NULL));
189 objmin = fnorm;
190 if (objective) {
191 PetscCall(SNESComputeObjective(snes, X, &objmin));
192 SNESCheckObjectiveDomainError(snes, objmin);
193 }
194 obj = objmin;
195
196 PetscCall(PetscObjectSAWsTakeAccess((PetscObject)snes));
197 snes->norm = fnorm;
198 PetscCall(PetscObjectSAWsGrantAccess((PetscObject)snes));
199 PetscCall(SNESLogConvergenceHistory(snes, fnorm, 0));
200 PetscCall(SNESConverged(snes, 0, 0.0, 0.0, fnorm));
201 PetscCall(SNESMonitor(snes, 0, fnorm));
202 if (snes->reason) PetscFunctionReturn(PETSC_SUCCESS);
203 PetscCall(SNESNGMRESUpdateSubspace_Private(snes, 0, 0, F, fnorm, X));
204
205 k_restart = 1;
206 l = 1;
207 ivec = 0;
208 for (k = 1; k < snes->max_its + 1; k++) {
209 /* Call general purpose update function */
210 PetscTryTypeMethod(snes, update, snes->iter);
211
212 /* Computation of x^M */
213 if (snes->npc && snes->npcside == PC_RIGHT) {
214 PetscCall(VecCopy(X, XM));
215 PetscCall(SNESSetInitialFunction(snes->npc, F));
216
217 PetscCall(PetscLogEventBegin(SNES_NPCSolve, snes->npc, XM, B, 0));
218 PetscCall(SNESSolve(snes->npc, B, XM));
219 PetscCall(PetscLogEventEnd(SNES_NPCSolve, snes->npc, XM, B, 0));
220
221 PetscCall(SNESGetConvergedReason(snes->npc, &reason));
222 if (reason < 0 && reason != SNES_DIVERGED_MAX_IT) {
223 snes->reason = SNES_DIVERGED_INNER;
224 PetscFunctionReturn(PETSC_SUCCESS);
225 }
226 PetscCall(SNESGetNPCFunction(snes, FM, &fMnorm));
227 } else {
228 /* no preconditioner -- just take gradient descent with line search */
229 PetscCall(VecCopy(F, Y));
230 PetscCall(VecCopy(F, FM));
231 PetscCall(VecCopy(X, XM));
232
233 fMnorm = fnorm;
234
235 PetscCall(SNESLineSearchApply(snes->linesearch, XM, FM, &fMnorm, Y));
236 if (snes->reason) PetscFunctionReturn(PETSC_SUCCESS);
237 }
238 if (objective) {
239 PetscCall(SNESComputeObjective(snes, XM, &objM));
240 SNESCheckObjectiveDomainError(snes, objM);
241 } else objM = fMnorm;
242 objmin = PetscMin(objmin, objM);
243
244 PetscCall(SNESNGMRESFormCombinedSolution_Private(snes, ivec, l, XM, FM, fMnorm, X, XA, FA));
245
246 /* differences for selection and restart */
247 if (ngmres->restart_type == SNES_NGMRES_RESTART_DIFFERENCE || ngmres->select_type == SNES_NGMRES_SELECT_DIFFERENCE) {
248 PetscCall(SNESNGMRESNorms_Private(snes, l, X, F, XM, FM, XA, FA, D, &dnorm, &dminnorm, &xMnorm, NULL, &yMnorm, &xAnorm, &fAnorm, &yAnorm));
249 } else {
250 PetscCall(SNESNGMRESNorms_Private(snes, l, X, F, XM, FM, XA, FA, D, NULL, NULL, &xMnorm, NULL, &yMnorm, &xAnorm, &fAnorm, &yAnorm));
251 }
252 SNESCheckFunctionDomainError(snes, fnorm);
253 if (objective) {
254 PetscCall(SNESComputeObjective(snes, XA, &objA));
255 SNESCheckObjectiveDomainError(snes, objA);
256 } else objA = fAnorm;
257
258 /* combination (additive) or selection (multiplicative) of the N-GMRES solution */
259 PetscCall(SNESNGMRESSelect_Private(snes, k_restart, XM, FM, xMnorm, fMnorm, yMnorm, objM, XA, FA, xAnorm, fAnorm, yAnorm, objA, dnorm, objmin, dminnorm, X, F, Y, &xnorm, &fnorm, &ynorm));
260 SNESCheckFunctionDomainError(snes, fnorm);
261 if (objective) {
262 PetscCall(SNESComputeObjective(snes, X, &obj));
263 SNESCheckObjectiveDomainError(snes, obj);
264 } else obj = fnorm;
265 selectRestart = PETSC_FALSE;
266
267 if (ngmres->restart_type == SNES_NGMRES_RESTART_DIFFERENCE) {
268 PetscCall(SNESNGMRESSelectRestart_Private(snes, l, obj, objM, objA, dnorm, objmin, dminnorm, &selectRestart));
269
270 /* if the restart conditions persist for more than restart_it iterations, restart. */
271 if (selectRestart) restart_count++;
272 else restart_count = 0;
273 } else if (ngmres->restart_type == SNES_NGMRES_RESTART_PERIODIC) {
274 if (k_restart > ngmres->restart_periodic) {
275 if (ngmres->monitor) PetscCall(PetscViewerASCIIPrintf(ngmres->monitor, "periodic restart after %" PetscInt_FMT " iterations\n", k_restart));
276 restart_count = ngmres->restart_it;
277 }
278 }
279
280 ivec = k_restart % ngmres->msize; /* replace the last used part of the subspace */
281
282 /* restart after restart conditions have persisted for a fixed number of iterations */
283 if (restart_count >= ngmres->restart_it) {
284 if (ngmres->monitor) PetscCall(PetscViewerASCIIPrintf(ngmres->monitor, "Restarted at iteration %" PetscInt_FMT "\n", k_restart));
285 restart_count = 0;
286 k_restart = 1;
287 l = 1;
288 ivec = 0;
289 /* q_{00} = nu */
290 PetscCall(SNESNGMRESUpdateSubspace_Private(snes, 0, 0, FM, fMnorm, XM));
291 } else {
292 /* select the current size of the subspace */
293 if (l < ngmres->msize) l++;
294 k_restart++;
295 /* place the current entry in the list of previous entries */
296 if (ngmres->candidate) {
297 objmin = PetscMin(objmin, objM);
298 PetscCall(SNESNGMRESUpdateSubspace_Private(snes, ivec, l, FM, fMnorm, XM));
299 } else {
300 objmin = PetscMin(objmin, obj);
301 PetscCall(SNESNGMRESUpdateSubspace_Private(snes, ivec, l, F, fnorm, X));
302 }
303 }
304
305 PetscCall(PetscObjectSAWsTakeAccess((PetscObject)snes));
306 snes->iter = k;
307 snes->norm = fnorm;
308 snes->ynorm = ynorm;
309 snes->xnorm = xnorm;
310 PetscCall(PetscObjectSAWsGrantAccess((PetscObject)snes));
311 PetscCall(SNESLogConvergenceHistory(snes, snes->norm, snes->iter));
312 PetscCall(SNESConverged(snes, snes->iter, 0, 0, fnorm));
313 PetscCall(SNESMonitor(snes, snes->iter, snes->norm));
314 if (snes->reason) PetscFunctionReturn(PETSC_SUCCESS);
315 }
316 snes->reason = SNES_DIVERGED_MAX_IT;
317 PetscFunctionReturn(PETSC_SUCCESS);
318 }
319
320 /*@
321 SNESNGMRESSetRestartFmRise - Increase the restart count if the step x_M increases the residual F_M inside a `SNESNGMRES` solve
322
323 Input Parameters:
324 + snes - the `SNES` context.
325 - flg - boolean value deciding whether to use the option or not, default is `PETSC_FALSE`
326
327 Options Database Key:
328 . -snes_ngmres_restart_fm_rise - Increase the restart count if the step x_M increases the residual F_M
329
330 Level: advanced
331
332 Notes:
333 If the proposed step x_M increases the residual F_M, it might be trying to get out of a stagnation area.
334 To help the solver do that, remove the current stored solutions and residuals whenever F_M increases.
335
336 This option must be used with the `SNESNGMRES` `SNESNGMRESRestartType` of `SNES_NGMRES_RESTART_DIFFERENCE`
337
338 .seealso: [](ch_snes), `SNES`, `SNES_NGMRES_RESTART_DIFFERENCE`, `SNESNGMRES`, `SNESNGMRESRestartType`, `SNESNGMRESSetRestartType()`
339 @*/
SNESNGMRESSetRestartFmRise(SNES snes,PetscBool flg)340 PetscErrorCode SNESNGMRESSetRestartFmRise(SNES snes, PetscBool flg)
341 {
342 PetscErrorCode (*f)(SNES, PetscBool);
343
344 PetscFunctionBegin;
345 PetscCall(PetscObjectQueryFunction((PetscObject)snes, "SNESNGMRESSetRestartFmRise_C", &f));
346 if (f) PetscCall((f)(snes, flg));
347 PetscFunctionReturn(PETSC_SUCCESS);
348 }
349
SNESNGMRESSetRestartFmRise_NGMRES(SNES snes,PetscBool flg)350 static PetscErrorCode SNESNGMRESSetRestartFmRise_NGMRES(SNES snes, PetscBool flg)
351 {
352 SNES_NGMRES *ngmres = (SNES_NGMRES *)snes->data;
353
354 PetscFunctionBegin;
355 ngmres->restart_fm_rise = flg;
356 PetscFunctionReturn(PETSC_SUCCESS);
357 }
358
SNESNGMRESGetRestartFmRise(SNES snes,PetscBool * flg)359 PetscErrorCode SNESNGMRESGetRestartFmRise(SNES snes, PetscBool *flg)
360 {
361 PetscErrorCode (*f)(SNES, PetscBool *);
362
363 PetscFunctionBegin;
364 PetscCall(PetscObjectQueryFunction((PetscObject)snes, "SNESNGMRESGetRestartFmRise_C", &f));
365 if (f) PetscCall((f)(snes, flg));
366 PetscFunctionReturn(PETSC_SUCCESS);
367 }
368
SNESNGMRESGetRestartFmRise_NGMRES(SNES snes,PetscBool * flg)369 static PetscErrorCode SNESNGMRESGetRestartFmRise_NGMRES(SNES snes, PetscBool *flg)
370 {
371 SNES_NGMRES *ngmres = (SNES_NGMRES *)snes->data;
372
373 PetscFunctionBegin;
374 *flg = ngmres->restart_fm_rise;
375 PetscFunctionReturn(PETSC_SUCCESS);
376 }
377
378 /*@
379 SNESNGMRESSetRestartType - Sets the restart type for `SNESNGMRES`.
380
381 Logically Collective
382
383 Input Parameters:
384 + snes - the iterative context
385 - rtype - restart type, see `SNESNGMRESRestartType`
386
387 Options Database Keys:
388 + -snes_ngmres_restart_type<difference,periodic,none> - set the restart type
389 - -snes_ngmres_restart <30> - sets the number of iterations before restart for periodic
390
391 Level: intermediate
392
393 .seealso: [](ch_snes), `SNES`, `SNES_NGMRES_RESTART_DIFFERENCE`, `SNESNGMRES`, `SNESNGMRESRestartType`, `SNESNGMRESSetRestartFmRise()`,
394 `SNESNGMRESSetSelectType()`
395 @*/
SNESNGMRESSetRestartType(SNES snes,SNESNGMRESRestartType rtype)396 PetscErrorCode SNESNGMRESSetRestartType(SNES snes, SNESNGMRESRestartType rtype)
397 {
398 PetscFunctionBegin;
399 PetscValidHeaderSpecific(snes, SNES_CLASSID, 1);
400 PetscTryMethod(snes, "SNESNGMRESSetRestartType_C", (SNES, SNESNGMRESRestartType), (snes, rtype));
401 PetscFunctionReturn(PETSC_SUCCESS);
402 }
403
404 /*@
405 SNESNGMRESSetSelectType - Sets the selection type for `SNESNGMRES`. This determines how the candidate solution and
406 combined solution are used to create the next iterate.
407
408 Logically Collective
409
410 Input Parameters:
411 + snes - the iterative context
412 - stype - selection type, see `SNESNGMRESSelectType`
413
414 Options Database Key:
415 . -snes_ngmres_select_type<difference,none,linesearch> - select type
416
417 Level: intermediate
418
419 Note:
420 The default line search used is the `SNESLINESEARCHSECANT` line search and it requires two additional function evaluations.
421
422 .seealso: [](ch_snes), `SNES`, `SNESNGMRES`, `SNESNGMRESSelectType`, `SNES_NGMRES_SELECT_NONE`, `SNES_NGMRES_SELECT_DIFFERENCE`, `SNES_NGMRES_SELECT_LINESEARCH`,
423 `SNESNGMRESSetRestartType()`
424 @*/
SNESNGMRESSetSelectType(SNES snes,SNESNGMRESSelectType stype)425 PetscErrorCode SNESNGMRESSetSelectType(SNES snes, SNESNGMRESSelectType stype)
426 {
427 PetscFunctionBegin;
428 PetscValidHeaderSpecific(snes, SNES_CLASSID, 1);
429 PetscTryMethod(snes, "SNESNGMRESSetSelectType_C", (SNES, SNESNGMRESSelectType), (snes, stype));
430 PetscFunctionReturn(PETSC_SUCCESS);
431 }
432
SNESNGMRESSetSelectType_NGMRES(SNES snes,SNESNGMRESSelectType stype)433 static PetscErrorCode SNESNGMRESSetSelectType_NGMRES(SNES snes, SNESNGMRESSelectType stype)
434 {
435 SNES_NGMRES *ngmres = (SNES_NGMRES *)snes->data;
436
437 PetscFunctionBegin;
438 ngmres->select_type = stype;
439 PetscFunctionReturn(PETSC_SUCCESS);
440 }
441
SNESNGMRESSetRestartType_NGMRES(SNES snes,SNESNGMRESRestartType rtype)442 static PetscErrorCode SNESNGMRESSetRestartType_NGMRES(SNES snes, SNESNGMRESRestartType rtype)
443 {
444 SNES_NGMRES *ngmres = (SNES_NGMRES *)snes->data;
445
446 PetscFunctionBegin;
447 ngmres->restart_type = rtype;
448 PetscFunctionReturn(PETSC_SUCCESS);
449 }
450
451 /*MC
452 SNESNGMRES - An implementation of the Nonlinear Generalized Minimum Residual method, Nonlinear GMRES, or N-GMRES {cite}`ow1`, {cite}`bruneknepleysmithtu15` for solving
453 nonlinear systems with `SNES`.
454
455 Level: beginner
456
457 Options Database Keys:
458 + -snes_ngmres_select_type<difference,none,linesearch> - choose the select between candidate and combined solution
459 . -snes_ngmres_restart_type<difference,none,periodic> - choose the restart conditions
460 . -snes_ngmres_candidate - Use `SNESNGMRES` variant which combines candidate solutions instead of actual solutions
461 . -snes_ngmres_m - Number of stored previous solutions and residuals
462 . -snes_ngmres_restart_it - Number of iterations the restart conditions hold before restart
463 . -snes_ngmres_gammaA - Residual tolerance for solution select between the candidate and combination
464 . -snes_ngmres_gammaC - Residual tolerance for restart
465 . -snes_ngmres_epsilonB - Difference tolerance between subsequent solutions triggering restart
466 . -snes_ngmres_deltaB - Difference tolerance between residuals triggering restart
467 . -snes_ngmres_restart_fm_rise - Restart on residual rise from $x_M$ step
468 . -snes_ngmres_monitor - Prints relevant information about the nonlinear GNMRES iterations
469 . -snes_linesearch_type <basic,l2,cp> - Line search type used for the default smoother
470 - -snes_ngmres_additive_snes_linesearch_type - line search type used to select between the candidate and combined solution with additive select type
471
472 Notes:
473 The N-GMRES method combines m previous solutions into a minimum-residual solution by solving a small linearized
474 optimization problem at each iteration.
475
476 Very similar to the `SNESANDERSON` algorithm.
477
478 Unlike the linear GMRES algorithm this algorithm does not compute a Krylov subspace using the Arnoldi process. Instead it stores a
479 collection of previous solutions and the residuals $ F(x) - b $ at those solutions.
480
481 This algorithm ignores any Jacobian provided with `SNESSetJacobian()`
482
483 Only supports left non-linear preconditioning.
484
485 .seealso: [](ch_snes), `SNESCreate()`, `SNES`, `SNESSetType()`, `SNESType`, `SNESANDERSON`, `SNESNGMRESSetSelectType()`, `SNESNGMRESSetRestartType()`,
486 `SNESNGMRESSetRestartFmRise()`, `SNESNGMRESSelectType`, `SNESNGMRESRestartType`
487 M*/
488
SNESCreate_NGMRES(SNES snes)489 PETSC_EXTERN PetscErrorCode SNESCreate_NGMRES(SNES snes)
490 {
491 SNES_NGMRES *ngmres;
492 SNESLineSearch linesearch;
493
494 PetscFunctionBegin;
495 snes->ops->destroy = SNESDestroy_NGMRES;
496 snes->ops->setup = SNESSetUp_NGMRES;
497 snes->ops->setfromoptions = SNESSetFromOptions_NGMRES;
498 snes->ops->view = SNESView_NGMRES;
499 snes->ops->solve = SNESSolve_NGMRES;
500 snes->ops->reset = SNESReset_NGMRES;
501
502 snes->usesnpc = PETSC_TRUE;
503 snes->usesksp = PETSC_FALSE;
504 snes->npcside = PC_RIGHT;
505
506 snes->alwayscomputesfinalresidual = PETSC_TRUE;
507
508 PetscCall(PetscNew(&ngmres));
509 snes->data = (void *)ngmres;
510 ngmres->msize = 30;
511
512 PetscCall(SNESParametersInitialize(snes));
513 PetscObjectParameterSetDefault(snes, max_funcs, 30000);
514 PetscObjectParameterSetDefault(snes, max_its, 10000);
515
516 ngmres->candidate = PETSC_FALSE;
517
518 PetscCall(SNESGetLineSearch(snes, &linesearch));
519 if (!((PetscObject)linesearch)->type_name) PetscCall(SNESLineSearchSetType(linesearch, SNESLINESEARCHBASIC));
520
521 ngmres->additive_linesearch = NULL;
522 ngmres->approxfunc = PETSC_FALSE;
523 ngmres->restart_it = 2;
524 ngmres->restart_periodic = 30;
525 ngmres->gammaA = 2.0;
526 ngmres->gammaC = 2.0;
527 ngmres->deltaB = 0.9;
528 ngmres->epsilonB = 0.1;
529 ngmres->restart_fm_rise = PETSC_FALSE;
530
531 ngmres->restart_type = SNES_NGMRES_RESTART_DIFFERENCE;
532 ngmres->select_type = SNES_NGMRES_SELECT_DIFFERENCE;
533
534 PetscCall(PetscObjectComposeFunction((PetscObject)snes, "SNESNGMRESSetSelectType_C", SNESNGMRESSetSelectType_NGMRES));
535 PetscCall(PetscObjectComposeFunction((PetscObject)snes, "SNESNGMRESSetRestartType_C", SNESNGMRESSetRestartType_NGMRES));
536 PetscCall(PetscObjectComposeFunction((PetscObject)snes, "SNESNGMRESSetRestartFmRise_C", SNESNGMRESSetRestartFmRise_NGMRES));
537 PetscCall(PetscObjectComposeFunction((PetscObject)snes, "SNESNGMRESGetRestartFmRise_C", SNESNGMRESGetRestartFmRise_NGMRES));
538 PetscFunctionReturn(PETSC_SUCCESS);
539 }
540