1 #pragma once
2
3 #include <petscsnes.h>
4 #include <petsc/private/petscimpl.h>
5
6 /* SUBMANSEC = SNES */
7
8 PETSC_EXTERN PetscBool SNESRegisterAllCalled;
9 PETSC_EXTERN PetscErrorCode SNESRegisterAll(void);
10
11 typedef struct _SNESOps *SNESOps;
12
13 struct _SNESOps {
14 PetscErrorCode (*computeinitialguess)(SNES, Vec, void *);
15 PetscErrorCode (*computescaling)(Vec, Vec, void *);
16 PetscErrorCode (*update)(SNES, PetscInt); /* General purpose function for update */
17 PetscErrorCode (*converged)(SNES, PetscInt, PetscReal, PetscReal, PetscReal, SNESConvergedReason *, void *);
18 PetscCtxDestroyFn *convergeddestroy;
19 PetscErrorCode (*setup)(SNES); /* routine to set up the nonlinear solver */
20 PetscErrorCode (*solve)(SNES); /* actual nonlinear solver */
21 PetscErrorCode (*view)(SNES, PetscViewer);
22 PetscErrorCode (*setfromoptions)(SNES, PetscOptionItems); /* sets options from database */
23 PetscErrorCode (*destroy)(SNES);
24 PetscErrorCode (*reset)(SNES);
25 PetscErrorCode (*ctxcompute)(SNES, PetscCtxRt);
26 PetscCtxDestroyFn *ctxdestroy;
27 PetscErrorCode (*computevariablebounds)(SNES, Vec, Vec); /* user provided routine to set box constrained variable bounds */
28 PetscErrorCode (*computepfunction)(SNES, Vec, Vec, void *);
29 PetscErrorCode (*computepjacobian)(SNES, Vec, Mat, Mat, void *);
30 PetscErrorCode (*load)(SNES, PetscViewer);
31 };
32
33 /*
34 Nonlinear solver context
35 */
36 #define MAXSNESMONITORS 5
37 #define MAXSNESREASONVIEWS 5
38
39 struct _p_SNES {
40 PETSCHEADER(struct _SNESOps);
41 DM dm;
42 PetscBool dmAuto; /* SNES created currently used DM automatically */
43 SNES npc;
44 PCSide npcside;
45 PetscBool usesnpc; /* type can use a nonlinear preconditioner */
46
47 /* ------------------------ User-provided stuff -------------------------------*/
48 PetscCtx ctx; /* user-defined context */
49
50 Vec vec_rhs; /* If non-null, solve F(x) = rhs */
51 Vec vec_sol; /* pointer to solution */
52
53 Vec vec_func; /* pointer to function */
54
55 Mat jacobian; /* Jacobian matrix */
56 Mat jacobian_pre; /* matrix used to construct the preconditioner of the Jacobian */
57 Mat picard; /* copy of jacobian_pre needed for Picard with -snes_mf_operator */
58 void *initialguessP; /* user-defined initial guess context */
59 KSP ksp; /* linear solver context */
60 SNESLineSearch linesearch; /* line search context */
61 PetscBool usesksp;
62 MatStructure matstruct; /* Used by Picard solver */
63
64 Vec vec_sol_update; /* pointer to solution update */
65
66 Vec scaling; /* scaling vector */
67 void *scaP; /* scaling context */
68
69 PetscReal precheck_picard_angle; /* For use with SNESLineSearchPreCheckPicard */
70
71 /* ------------------------Time stepping hooks-----------------------------------*/
72
73 /* ---------------- PETSc-provided (or user-provided) stuff ---------------------*/
74
75 PetscErrorCode (*monitor[MAXSNESMONITORS])(SNES, PetscInt, PetscReal, void *); /* monitor routine */
76 PetscCtxDestroyFn *monitordestroy[MAXSNESMONITORS]; /* monitor context destroy routine */
77 void *monitorcontext[MAXSNESMONITORS]; /* monitor context */
78 PetscInt numbermonitors; /* number of monitors */
79 PetscBool pauseFinal; /* pause all drawing monitor at the final iterate */
80 void *cnvP; /* convergence context */
81 SNESConvergedReason reason; /* converged reason */
82
83 PetscViewer convergedreasonviewer;
84 PetscViewerFormat convergedreasonformat;
85 PetscErrorCode (*reasonview[MAXSNESREASONVIEWS])(SNES, void *); /* snes converged reason view */
86 PetscCtxDestroyFn *reasonviewdestroy[MAXSNESREASONVIEWS]; /* reason view context destroy routine */
87 void *reasonviewcontext[MAXSNESREASONVIEWS]; /* reason view context */
88 PetscInt numberreasonviews; /* number of reason views */
89 PetscBool errorifnotconverged;
90
91 /* --- Routines and data that are unique to each particular solver --- */
92
93 PetscBool setupcalled; /* true if setup has been called */
94 void *data; /* implementation-specific data */
95
96 /* -------------------------- Parameters -------------------------------------- */
97
98 PetscInt nfuncs; /* number of function evaluations */
99 PetscInt iter; /* global iteration number */
100 PetscInt linear_its; /* total number of linear solver iterations */
101 PetscReal norm; /* residual norm of current iterate */
102 PetscReal ynorm; /* update norm of current iterate */
103 PetscReal xnorm; /* solution norm of current iterate */
104 PetscBool forceiteration; /* Force SNES to take at least one iteration regardless of the initial residual norm */
105 PetscInt lagpreconditioner; /* SNESSetLagPreconditioner() */
106 PetscInt lagjacobian; /* SNESSetLagJacobian() */
107 PetscInt jac_iter; /* The present iteration of the Jacobian lagging */
108 PetscBool lagjac_persist; /* The jac_iter persists until reset */
109 PetscInt pre_iter; /* The present iteration of the Preconditioner lagging */
110 PetscBool lagpre_persist; /* The pre_iter persists until reset */
111 PetscInt gridsequence; /* number of grid sequence steps to take; defaults to zero */
112 PetscObjectParameterDeclare(PetscInt, max_its); /* max number of iterations */
113 PetscObjectParameterDeclare(PetscInt, max_funcs); /* max number of function evals */
114 PetscObjectParameterDeclare(PetscReal, rtol); /* relative tolerance */
115 PetscObjectParameterDeclare(PetscReal, divtol); /* relative divergence tolerance */
116 PetscObjectParameterDeclare(PetscReal, abstol); /* absolute tolerance */
117 PetscObjectParameterDeclare(PetscReal, stol); /* step length tolerance*/
118
119 PetscBool vec_func_init_set; /* the initial function has been set */
120
121 SNESNormSchedule normschedule; /* Norm computation type for SNES instance */
122 SNESFunctionType functype; /* Function type for the SNES instance */
123
124 /* ------------------------ Default work-area management ---------------------- */
125
126 PetscInt nwork;
127 Vec *work;
128
129 /* ---------------------------------- Testing --------------------------------- */
130 PetscBool testFunc; // Test the function routine
131 PetscBool testJac; // Test the Jacobian routine
132
133 /* ------------------------- Miscellaneous Information ------------------------ */
134
135 PetscInt setfromoptionscalled;
136 PetscReal *conv_hist; /* If !0, stores function norm (or
137 gradient norm) at each iteration */
138 PetscInt *conv_hist_its; /* linear iterations for each Newton step */
139 size_t conv_hist_len; /* size of convergence history array */
140 size_t conv_hist_max; /* actual amount of data in conv_history */
141 PetscBool conv_hist_reset; /* reset counter for each new SNES solve */
142 PetscBool conv_hist_alloc;
143 PetscBool counters_reset; /* reset counter for each new SNES solve */
144
145 /* the next two are used for failures in the line search; they should be put elsewhere */
146 PetscInt numFailures; /* number of unsuccessful step attempts */
147 PetscInt maxFailures; /* maximum number of unsuccessful step attempts */
148
149 PetscInt numLinearSolveFailures;
150 PetscInt maxLinearSolveFailures;
151
152 PetscBool functiondomainerror; /* set with SNESSetFunctionDomainError() */
153 PetscBool objectivedomainerror; /* set with SNESSetObjectiveDomainError() */
154 PetscBool jacobiandomainerror; /* set with SNESSetJacobianDomainError() */
155 PetscBool checkjacdomainerror; /* if SNESCheckJacobianDomainError() is called after Jacobian evaluations */
156
157 PetscBool ksp_ewconv; /* flag indicating use of Eisenstat-Walker KSP convergence criteria */
158 void *kspconvctx; /* Eisenstat-Walker KSP convergence context */
159
160 /* SNESConvergedDefault context: split it off into a separate var/struct to be passed as context to SNESConvergedDefault? */
161 PetscReal ttol; /* rtol*initial_residual_norm */
162 PetscReal rnorm0; /* initial residual norm (used for divergence testing) */
163
164 Vec *vwork; /* more work vectors for Jacobian approx */
165 PetscInt nvwork;
166
167 PetscBool mf; /* -snes_mf OR -snes_mf_operator was used on this snes */
168 PetscBool mf_operator; /* -snes_mf_operator was used on this snes */
169 PetscInt mf_version; /* The version of snes_mf used */
170
171 PetscReal vizerotolerance; /* tolerance for considering an x[] value to be on the bound */
172 Vec xl, xu; /* upper and lower bounds for box constrained VI problems */
173 PetscInt ntruebounds; /* number of non-infinite bounds set for VI box constraints */
174 PetscBool usersetbounds; /* bounds have been set via SNESVISetVariableBounds(), rather than via computevariablebounds() callback. */
175
176 PetscBool alwayscomputesfinalresidual; /* Does SNESSolve_XXX always compute the value of the residual at the final
177 * solution and put it in vec_func? Used inside SNESSolve_FAS to determine
178 * if the final residual must be computed before restricting or prolonging
179 * it. */
180 };
181
182 typedef struct _p_DMSNES *DMSNES;
183 typedef struct _DMSNESOps *DMSNESOps;
184 struct _DMSNESOps {
185 SNESFunctionFn *computefunction;
186 SNESFunctionFn *computemffunction;
187 SNESJacobianFn *computejacobian;
188
189 /* objective */
190 SNESObjectiveFn *computeobjective;
191
192 /* Picard iteration functions */
193 SNESFunctionFn *computepfunction;
194 SNESJacobianFn *computepjacobian;
195
196 /* User-defined smoother */
197 SNESNGSFn *computegs;
198
199 PetscErrorCode (*destroy)(DMSNES);
200 PetscErrorCode (*duplicate)(DMSNES, DMSNES);
201 };
202
203 /*S
204 DMSNES - Object held by a `DM` that contains all the callback functions and their contexts needed by a `SNES`
205
206 Level: developer
207
208 Notes:
209 Users provides callback functions and their contexts to `SNES` using, for example, `SNESSetFunction()`. These values are stored
210 in a `DMSNES` that is contained in the `DM` associated with the `SNES`. If no `DM` was provided by
211 the user with `SNESSetDM()` it is automatically created by `SNESGetDM()` with `DMShellCreate()`.
212
213 Users very rarely need to worked directly with the `DMSNES` object, rather they work with the `SNES` and the `DM` they created
214
215 Multiple `DM` can share a single `DMSNES`, often each `DM` is associated with
216 a grid refinement level. `DMGetDMSNES()` returns the `DMSNES` associated with a `DM`. `DMGetDMSNESWrite()` returns a unique
217 `DMSNES` that is only associated with the current `DM`, making a copy of the shared `DMSNES` if needed (copy-on-write).
218
219 See `DMKSP` for details on why there is a needed for `DMSNES` instead of simply storing the user callbacks directly in the `DM` or the `TS`
220
221 Developer Note:
222 The `originaldm` inside the `DMSNES` is NOT reference counted (to prevent a reference count loop between a `DM` and a `DMSNES`).
223 The `DM` on which this context was first created is cached here to implement one-way
224 copy-on-write. When `DMGetDMSNESWrite()` sees a request using a different `DM`, it makes a copy of the `TSDM`. Thus, if a user
225 only interacts directly with one level, e.g., using `TSSetIFunction()`, then coarse levels of a multilevel item
226 integrator are built, then the user changes the routine with another call to `TSSetIFunction()`, it automatically
227 propagates to all the levels. If instead, they get out a specific level and set the function on that level,
228 subsequent changes to the original level will no longer propagate to that level.
229
230 .seealso: [](ch_snes), `SNES`, `SNESCreate()`, `DM`, `DMGetDMSNESWrite()`, `DMGetDMSNES()`, `DMKSP`, `DMTS`, `DMSNESSetFunction()`, `DMSNESGetFunction()`,
231 `DMSNESSetFunctionContextDestroy()`, `DMSNESSetMFFunction()`, `DMSNESSetNGS()`, `DMSNESGetNGS()`, `DMSNESSetJacobian()`, `DMSNESGetJacobian()`,
232 `DMSNESSetJacobianContextDestroy()`, `DMSNESSetPicard()`, `DMSNESGetPicard()`, `DMSNESSetObjective()`, `DMSNESGetObjective()`, `DMCopyDMSNES()`
233 S*/
234 struct _p_DMSNES {
235 PETSCHEADER(struct _DMSNESOps);
236 PetscContainer functionctxcontainer;
237 PetscContainer jacobianctxcontainer;
238 void *mffunctionctx;
239 void *gsctx;
240 void *pctx;
241 void *objectivectx;
242
243 void *data;
244
245 /* See developer note for DMSNES above */
246 DM originaldm;
247 };
248 PETSC_EXTERN PetscErrorCode DMGetDMSNES(DM, DMSNES *);
249 PETSC_EXTERN PetscErrorCode DMSNESView(DMSNES, PetscViewer);
250 PETSC_EXTERN PetscErrorCode DMSNESLoad(DMSNES, PetscViewer);
251 PETSC_EXTERN PetscErrorCode DMGetDMSNESWrite(DM, DMSNES *);
252
253 /* Context for Eisenstat-Walker convergence criteria for KSP solvers */
254 typedef struct {
255 PetscInt version; /* flag indicating version (1,2,3 or 4) */
256 PetscReal rtol_0; /* initial rtol */
257 PetscReal rtol_last; /* last rtol */
258 PetscReal rtol_max; /* maximum rtol */
259 PetscReal gamma; /* mult. factor for version 2 rtol computation */
260 PetscReal alpha; /* power for version 2 rtol computation */
261 PetscReal alpha2; /* power for safeguard */
262 PetscReal threshold; /* threshold for imposing safeguard */
263 PetscReal lresid_last; /* linear residual from last iteration */
264 PetscReal norm_last; /* function norm from last iteration */
265 PetscReal norm_first; /* function norm from the beginning of the first iteration. */
266 PetscReal rtol_last_2, rk_last, rk_last_2;
267 PetscReal v4_p1, v4_p2, v4_p3, v4_m1, v4_m2, v4_m3, v4_m4;
268 } SNESKSPEW;
269
SNESLogConvergenceHistory(SNES snes,PetscReal res,PetscInt its)270 static inline PetscErrorCode SNESLogConvergenceHistory(SNES snes, PetscReal res, PetscInt its)
271 {
272 PetscFunctionBegin;
273 PetscCall(PetscObjectSAWsTakeAccess((PetscObject)snes));
274 if (snes->conv_hist && snes->conv_hist_max > snes->conv_hist_len) {
275 if (snes->conv_hist) snes->conv_hist[snes->conv_hist_len] = res;
276 if (snes->conv_hist_its) snes->conv_hist_its[snes->conv_hist_len] = its;
277 snes->conv_hist_len++;
278 }
279 PetscCall(PetscObjectSAWsGrantAccess((PetscObject)snes));
280 PetscFunctionReturn(PETSC_SUCCESS);
281 }
282
283 PETSC_EXTERN PetscErrorCode SNESVIProjectOntoBounds(SNES, Vec);
284 PETSC_INTERN PetscErrorCode SNESVICheckLocalMin_Private(SNES, Mat, Vec, Vec, PetscReal, PetscBool *);
285 PETSC_INTERN PetscErrorCode SNESReset_VI(SNES);
286 PETSC_INTERN PetscErrorCode SNESDestroy_VI(SNES);
287 PETSC_INTERN PetscErrorCode SNESView_VI(SNES, PetscViewer);
288 PETSC_INTERN PetscErrorCode SNESSetFromOptions_VI(SNES, PetscOptionItems);
289 PETSC_INTERN PetscErrorCode SNESSetUp_VI(SNES);
290 PETSC_EXTERN_TYPEDEF typedef PetscErrorCode SNESVIComputeVariableBoundsFn(SNES, Vec, Vec);
291 PETSC_EXTERN_TYPEDEF typedef SNESVIComputeVariableBoundsFn *SNESVIComputeVariableBoundsFunction; // deprecated version
292 PETSC_INTERN PetscErrorCode SNESVISetComputeVariableBounds_VI(SNES, SNESVIComputeVariableBoundsFn);
293 PETSC_INTERN PetscErrorCode SNESVISetVariableBounds_VI(SNES, Vec, Vec);
294 PETSC_INTERN PetscErrorCode SNESConvergedDefault_VI(SNES, PetscInt, PetscReal, PetscReal, PetscReal, SNESConvergedReason *, void *);
295
296 PETSC_INTERN PetscErrorCode DMSNESUnsetFunctionContext_Internal(DM);
297 PETSC_EXTERN PetscErrorCode DMSNESUnsetJacobianContext_Internal(DM);
298
299 PETSC_EXTERN PetscLogEvent SNES_Solve;
300 PETSC_EXTERN PetscLogEvent SNES_SetUp;
301 PETSC_EXTERN PetscLogEvent SNES_FunctionEval;
302 PETSC_EXTERN PetscLogEvent SNES_JacobianEval;
303 PETSC_EXTERN PetscLogEvent SNES_NGSEval;
304 PETSC_EXTERN PetscLogEvent SNES_NGSFuncEval;
305 PETSC_EXTERN PetscLogEvent SNES_NewtonALEval;
306 PETSC_EXTERN PetscLogEvent SNES_NPCSolve;
307 PETSC_EXTERN PetscLogEvent SNES_ObjectiveEval;
308
309 PETSC_INTERN PetscBool SNEScite;
310 PETSC_INTERN const char SNESCitation[];
311
312 /* Used by TAOBNK solvers */
313 PETSC_SINGLE_LIBRARY_VISIBILITY_INTERNAL PetscErrorCode KSPPostSolve_SNESEW(KSP, Vec, Vec, void *);
314 PETSC_SINGLE_LIBRARY_VISIBILITY_INTERNAL PetscErrorCode KSPPreSolve_SNESEW(KSP, Vec, Vec, void *);
315 PETSC_SINGLE_LIBRARY_VISIBILITY_INTERNAL PetscErrorCode SNESEWSetFromOptions_Private(SNESKSPEW *, PetscBool, MPI_Comm, const char *);
316
317 /*MC
318 SNESCheckFunctionDomainError - Called after a `SNESComputeFunction()` and `VecNorm()` in a `SNES` solver to check if the function norm is infinity or NaN and
319 if the function callback set with `SNESSetFunction()` called `SNESSetFunctionDomainError()`.
320
321 Synopsis:
322 #include <snesimpl.h>
323 void SNESCheckFunctionDomainError(SNES snes, PetscReal fnorm)
324
325 Collective
326
327 Input Parameters:
328 + snes - the `SNES` solver object
329 - fnorm - the value of the norm
330
331 Level: developer
332
333 Notes:
334 If `fnorm` is infinity or NaN and `SNESSetErrorIfNotConverged()` was set, this immediately generates a `PETSC_ERR_CONV_FAILED`.
335
336 If `fnorm` is infinity or NaN and `SNESSetFunctionDomainError()` was called, this sets the `SNESConvergedReason` to `SNES_DIVERGED_FUNCTION_DOMAIN`
337 and exits the solver
338
339 Otherwise if `fnorm` is infinity or NaN, this sets the `SNESConvergedReason` to `SNES_DIVERGED_FUNCTION_NANORINF` and exits the solver
340
341 Developer Note:
342 This function exists so that `SNESSetFunctionDomainError()` does not need to be a collective operation since making it collective
343 would be cumbersome in most applications and require extra communication. Instead, `SNESSetFunctionDomainError()` sets the `functiondomainerror`
344 flag in the `SNES` object to true, `SNESComputeFunction()` checks that flag and sets a NaN into its local part of the vector if the flag has been set.
345 Then, when `VecNorm()` is called on the vector containing the computed function value, any NaN is propagated to all MPI processes without
346 any additional communication. Virtually all nonlinear solvers need to compute the function norm at some point so no extra communication needs to take place.
347
348 .seealso: [](ch_snes), `SNESSetFunctionDomainError()`, `SNESCheckObjectiveDomainError()`, `PETSC_ERR_CONV_FAILED`, `SNESSetErrorIfNotConverged()`, `SNES_DIVERGED_FUNCTION_DOMAIN`,
349 `SNESConvergedReason`, `SNES_DIVERGED_FUNCTION_NAN`
350 MC*/
351 #define SNESCheckFunctionDomainError(snes, fnorm) \
352 do { \
353 if (PetscIsInfOrNanReal(fnorm)) { \
354 PetscCheck(!snes->errorifnotconverged, PetscObjectComm((PetscObject)snes), PETSC_ERR_NOT_CONVERGED, "SNESSolve has not converged due to infinity or NaN norm"); \
355 { \
356 PetscBool domainerror; \
357 PetscCallMPI(MPIU_Allreduce(&snes->functiondomainerror, &domainerror, 1, MPI_C_BOOL, MPI_LOR, PetscObjectComm((PetscObject)snes))); \
358 if (domainerror) snes->reason = SNES_DIVERGED_FUNCTION_DOMAIN; \
359 else snes->reason = SNES_DIVERGED_FUNCTION_NANORINF; \
360 PetscFunctionReturn(PETSC_SUCCESS); \
361 } \
362 } \
363 } while (0)
364
365 /*MC
366 SNESCheckObjectiveDomainError - Called after a `SNESComputeObjective()` in a `SNES` solver to check if the objective value is infinity or NaN and/or
367 if the function callback set with `SNESSetObjective()` called `SNESSetObjectiveDomainError()`.
368
369 Synopsis:
370 #include <snesimpl.h>
371 void SNESCheckObjectiveDomainError(SNES snes, PetscReal fobj)
372
373 Collective
374
375 Input Parameters:
376 + snes - the `SNES` solver object
377 - fobj - the value of the objective function
378
379 Level: developer
380
381 Notes:
382 If `fobj` is infinity or NaN and `SNESSetErrorIfNotConverged()` was set, this immediately generates a `PETSC_ERR_CONV_FAILED`.
383
384 If `SNESSetObjectiveDomainError()` was called, this sets the `SNESConvergedReason` to `SNES_DIVERGED_OBJECTIVE_DOMAIN`
385 and exits the solver
386
387 Otherwise if `fobj` is infinity or NaN, this sets the `SNESConvergedReason` to `SNES_DIVERGED_OBJECTIVE_NANORINF` and exits the solver
388
389 .seealso: [](ch_snes), `SNESSetObjectiveDomainError()`, `PETSC_ERR_CONV_FAILED`, `SNESSetErrorIfNotConverged()`, `SNES_DIVERGED_OBJECTIVE_DOMAIN`, `SNES_DIVERGED_FUNCTION_DOMAIN`,
390 `SNESSetFunctionDomainError()`, `SNESConvergedReason`, `SNES_DIVERGED_OBJECTIVE_NANORINF`, `SNES_DIVERGED_FUNCTION_NAN`, `SNESLineSearchCheckObjectiveDomainError()`
391 MC*/
392 #define SNESCheckObjectiveDomainError(snes, fobj) \
393 do { \
394 if (snes->errorifnotconverged) { \
395 PetscCheck(!snes->objectivedomainerror, PetscObjectComm((PetscObject)snes), PETSC_ERR_NOT_CONVERGED, "SNESSolve has not converged due objective domain error"); \
396 PetscCheck(!PetscIsInfOrNanReal(fobj), PetscObjectComm((PetscObject)snes), PETSC_ERR_NOT_CONVERGED, "SNESSolve has not converged due to infinity or NaN norm"); \
397 } \
398 if (snes->objectivedomainerror) { \
399 snes->reason = SNES_DIVERGED_OBJECTIVE_DOMAIN; \
400 PetscFunctionReturn(PETSC_SUCCESS); \
401 } else if (PetscIsInfOrNanReal(fobj)) { \
402 snes->reason = SNES_DIVERGED_OBJECTIVE_NANORINF; \
403 PetscFunctionReturn(PETSC_SUCCESS); \
404 } \
405 } while (0)
406
407 /*MC
408 SNESCheckJacobianDomainError - Called after a `SNESComputeJacobian()` in a SNES solver to check if `SNESSetJacobianDomainError()` has been called.
409
410 Synopsis:
411 #include <snesimpl.h>
412 void SNESCheckJacobianDomainError(SNES snes)
413
414 Collective
415
416 Input Parameters:
417 . snes - the `SNES` solver object
418
419 Level: developer
420
421 Notes:
422 This turns the non-collective `SNESSetJacobianDomainError()` into a collective operation
423
424 This check is done in debug mode or if `SNESSetCheckJacobianDomainError()` has been called
425
426 .seealso: [](ch_snes), `SNESSetCheckJacobianDomainError()`, `SNESCheckObjectiveDomainError()`, `SNESSetFunctionDomainError()`, `PETSC_ERR_CONV_FAILED`, `SNESSetErrorIfNotConverged()`, `SNES_DIVERGED_FUNCTION_DOMAIN`,
427 `SNESConvergedReason`, `SNES_DIVERGED_FUNCTION_NAN`
428 MC*/
429 #define SNESCheckJacobianDomainError(snes) \
430 do { \
431 if (snes->checkjacdomainerror) { \
432 PetscBool domainerror; \
433 PetscCallMPI(MPIU_Allreduce(&snes->jacobiandomainerror, &domainerror, 1, MPI_C_BOOL, MPI_LOR, PetscObjectComm((PetscObject)snes))); \
434 if (domainerror) { \
435 snes->reason = SNES_DIVERGED_JACOBIAN_DOMAIN; \
436 PetscCheck(!snes->errorifnotconverged, PetscObjectComm((PetscObject)snes), PETSC_ERR_NOT_CONVERGED, "SNESSolve has not converged due to Jacobian domain error"); \
437 PetscFunctionReturn(PETSC_SUCCESS); \
438 } \
439 } \
440 } while (0)
441
442 /*MC
443 SNESCheckLineSearchFailure - Checks if a `SNESLineSearchApply()` has failed and possibly ends the current `SNESSolve()` if so
444
445 Synopsis:
446 #include <snesimpl.h>
447 void SNESCheckLineSearchFailure(SNES snes)
448
449 Collective
450
451 Input Parameters:
452 . snes - the `SNES` solver object
453
454 Level: developer
455
456 Notes:
457 If `SNESLineSearchApply()` produces a `SNES_LINESEARCH_FAILED_NANORINF` or `SNES_LINESEARCH_FAILED_NANORINF` the `SNESSolve()` is ended.
458
459 If the `SNESLineSearchApply()` produces any other failure reason and the number of failed steps is greater than the number set with
460 `SNESSetMaxNonlinearStepFailures()` the `SNESSolve()` is ended
461
462 .seealso: [](ch_snes), `SNESLineSearchApply()`, `SNESSetFunctionDomainError()`, `PETSC_ERR_CONV_FAILED`, `SNESSetErrorIfNotConverged()`, `SNES_DIVERGED_FUNCTION_DOMAIN`,
463 `SNESConvergedReason`, `SNES_DIVERGED_FUNCTION_NAN`, `SNESSolve()`, `SNESSetMaxNonlinearStepFailures()`
464 MC*/
465 #define SNESCheckLineSearchFailure(snes) \
466 do { \
467 SNESLineSearchReason lsreason; \
468 PetscCall(SNESLineSearchGetReason(snes->linesearch, &lsreason)); \
469 if (lsreason) { \
470 if (lsreason == SNES_LINESEARCH_FAILED_FUNCTION_DOMAIN) { \
471 PetscCheck(!snes->errorifnotconverged, PetscObjectComm((PetscObject)snes), PETSC_ERR_NOT_CONVERGED, "SNESLineSearchApply() has produced failure with function domain"); \
472 snes->reason = SNES_DIVERGED_FUNCTION_DOMAIN; \
473 PetscFunctionReturn(PETSC_SUCCESS); \
474 } \
475 if (lsreason == SNES_LINESEARCH_FAILED_NANORINF) { \
476 PetscCheck(!snes->errorifnotconverged, PetscObjectComm((PetscObject)snes), PETSC_ERR_NOT_CONVERGED, "SNESLineSearchApply() has produced failure with infinity or NaN"); \
477 snes->reason = SNES_DIVERGED_FUNCTION_NANORINF; \
478 PetscFunctionReturn(PETSC_SUCCESS); \
479 } \
480 if (lsreason == SNES_LINESEARCH_FAILED_OBJECTIVE_DOMAIN) { \
481 PetscCheck(!snes->errorifnotconverged, PetscObjectComm((PetscObject)snes), PETSC_ERR_NOT_CONVERGED, "SNESLineSearchApply() has produced failure with objective function domain"); \
482 snes->reason = SNES_DIVERGED_FUNCTION_DOMAIN; \
483 PetscFunctionReturn(PETSC_SUCCESS); \
484 } \
485 if (lsreason == SNES_LINESEARCH_FAILED_JACOBIAN_DOMAIN) { \
486 PetscCheck(!snes->errorifnotconverged, PetscObjectComm((PetscObject)snes), PETSC_ERR_NOT_CONVERGED, "SNESLineSearchApply() has produced failure with Jacobian domain"); \
487 snes->reason = SNES_DIVERGED_JACOBIAN_DOMAIN; \
488 PetscFunctionReturn(PETSC_SUCCESS); \
489 } \
490 if (++snes->numFailures >= snes->maxFailures) { \
491 PetscCheck(!snes->errorifnotconverged, PetscObjectComm((PetscObject)snes), PETSC_ERR_NOT_CONVERGED, "SNESLineSearchApply() has produced failure"); \
492 snes->reason = SNES_DIVERGED_LINE_SEARCH; \
493 PetscFunctionReturn(PETSC_SUCCESS); \
494 } \
495 } \
496 } while (0)
497
498 #define SNESCheckKSPSolve(snes) \
499 do { \
500 KSPConvergedReason kspreason; \
501 PetscInt lits; \
502 PetscCall(KSPGetIterationNumber(snes->ksp, &lits)); \
503 snes->linear_its += lits; \
504 PetscCall(KSPGetConvergedReason(snes->ksp, &kspreason)); \
505 if (kspreason < 0) { \
506 if (kspreason == KSP_DIVERGED_NANORINF) { \
507 PetscBool domainerror; \
508 PetscCallMPI(MPIU_Allreduce(&snes->functiondomainerror, &domainerror, 1, MPI_C_BOOL, MPI_LOR, PetscObjectComm((PetscObject)snes))); \
509 if (domainerror) { \
510 snes->reason = SNES_DIVERGED_FUNCTION_DOMAIN; \
511 snes->functiondomainerror = PETSC_FALSE; \
512 } else snes->reason = SNES_DIVERGED_LINEAR_SOLVE; \
513 PetscFunctionReturn(PETSC_SUCCESS); \
514 } else { \
515 if (++snes->numLinearSolveFailures >= snes->maxLinearSolveFailures) { \
516 PetscCall(PetscInfo(snes, "iter=%" PetscInt_FMT ", number linear solve failures %" PetscInt_FMT " greater than current SNES allowed %" PetscInt_FMT ", stopping solve\n", snes->iter, snes->numLinearSolveFailures, snes->maxLinearSolveFailures)); \
517 snes->reason = SNES_DIVERGED_LINEAR_SOLVE; \
518 PetscFunctionReturn(PETSC_SUCCESS); \
519 } \
520 } \
521 } \
522 } while (0)
523
524 #define SNESNeedNorm_Private(snes, iter) \
525 (((iter) == (snes)->max_its && ((snes)->normschedule == SNES_NORM_FINAL_ONLY || (snes)->normschedule == SNES_NORM_INITIAL_FINAL_ONLY)) || ((iter) == 0 && ((snes)->normschedule == SNES_NORM_INITIAL_ONLY || (snes)->normschedule == SNES_NORM_INITIAL_FINAL_ONLY)) || \
526 (snes)->normschedule == SNES_NORM_ALWAYS)
527