Actual source code: snesimpl.h
1: #pragma once
3: #include <petscsnes.h>
4: #include <petsc/private/petscimpl.h>
6: PETSC_EXTERN PetscBool SNESRegisterAllCalled;
7: PETSC_EXTERN PetscErrorCode SNESRegisterAll(void);
9: typedef struct _SNESOps *SNESOps;
11: struct _SNESOps {
12: PetscErrorCode (*computeinitialguess)(SNES, Vec, void *);
13: PetscErrorCode (*computescaling)(Vec, Vec, void *);
14: PetscErrorCode (*update)(SNES, PetscInt); /* General purpose function for update */
15: PetscErrorCode (*converged)(SNES, PetscInt, PetscReal, PetscReal, PetscReal, SNESConvergedReason *, void *);
16: PetscErrorCode (*convergeddestroy)(void *);
17: PetscErrorCode (*setup)(SNES); /* routine to set up the nonlinear solver */
18: PetscErrorCode (*solve)(SNES); /* actual nonlinear solver */
19: PetscErrorCode (*view)(SNES, PetscViewer);
20: PetscErrorCode (*setfromoptions)(SNES, PetscOptionItems *); /* sets options from database */
21: PetscErrorCode (*destroy)(SNES);
22: PetscErrorCode (*reset)(SNES);
23: PetscErrorCode (*usercompute)(SNES, void **);
24: PetscErrorCode (*userdestroy)(void **);
25: PetscErrorCode (*computevariablebounds)(SNES, Vec, Vec); /* user provided routine to set box constrained variable bounds */
26: PetscErrorCode (*computepfunction)(SNES, Vec, Vec, void *);
27: PetscErrorCode (*computepjacobian)(SNES, Vec, Mat, Mat, void *);
28: PetscErrorCode (*load)(SNES, PetscViewer);
29: };
31: /*
32: Nonlinear solver context
33: */
34: #define MAXSNESMONITORS 5
35: #define MAXSNESREASONVIEWS 5
37: struct _p_SNES {
38: PETSCHEADER(struct _SNESOps);
39: DM dm;
40: PetscBool dmAuto; /* SNES created currently used DM automatically */
41: SNES npc;
42: PCSide npcside;
43: PetscBool usesnpc; /* type can use a nonlinear preconditioner */
45: /* ------------------------ User-provided stuff -------------------------------*/
46: void *user; /* user-defined context */
48: Vec vec_rhs; /* If non-null, solve F(x) = rhs */
49: Vec vec_sol; /* pointer to solution */
51: Vec vec_func; /* pointer to function */
53: Mat jacobian; /* Jacobian matrix */
54: Mat jacobian_pre; /* preconditioner matrix */
55: Mat picard; /* copy of preconditioner matrix needed for Picard with -snes_mf_operator */
56: void *initialguessP; /* user-defined initial guess context */
57: KSP ksp; /* linear solver context */
58: SNESLineSearch linesearch; /* line search context */
59: PetscBool usesksp;
60: MatStructure matstruct; /* Used by Picard solver */
62: Vec vec_sol_update; /* pointer to solution update */
64: Vec scaling; /* scaling vector */
65: void *scaP; /* scaling context */
67: PetscReal precheck_picard_angle; /* For use with SNESLineSearchPreCheckPicard */
69: /* ------------------------Time stepping hooks-----------------------------------*/
71: /* ---------------- PETSc-provided (or user-provided) stuff ---------------------*/
73: PetscErrorCode (*monitor[MAXSNESMONITORS])(SNES, PetscInt, PetscReal, void *); /* monitor routine */
74: PetscErrorCode (*monitordestroy[MAXSNESMONITORS])(void **); /* monitor context destroy routine */
75: void *monitorcontext[MAXSNESMONITORS]; /* monitor context */
76: PetscInt numbermonitors; /* number of monitors */
77: PetscBool pauseFinal; /* pause all drawing monitor at the final iterate */
78: void *cnvP; /* convergence context */
79: SNESConvergedReason reason; /* converged reason */
80: PetscErrorCode (*reasonview[MAXSNESREASONVIEWS])(SNES, void *); /* snes converged reason view */
81: PetscErrorCode (*reasonviewdestroy[MAXSNESREASONVIEWS])(void **); /* reason view context destroy routine */
82: void *reasonviewcontext[MAXSNESREASONVIEWS]; /* reason view context */
83: PetscInt numberreasonviews; /* number of reason views */
84: PetscBool errorifnotconverged;
86: /* --- Routines and data that are unique to each particular solver --- */
88: PetscBool setupcalled; /* true if setup has been called */
89: void *data; /* implementation-specific data */
91: /* -------------------------- Parameters -------------------------------------- */
93: PetscInt max_its; /* max number of iterations */
94: PetscInt max_funcs; /* max number of function evals */
95: PetscInt nfuncs; /* number of function evaluations */
96: PetscInt iter; /* global iteration number */
97: PetscInt linear_its; /* total number of linear solver iterations */
98: PetscReal norm; /* residual norm of current iterate */
99: PetscReal ynorm; /* update norm of current iterate */
100: PetscReal xnorm; /* solution norm of current iterate */
101: PetscReal rtol; /* relative tolerance */
102: PetscReal divtol; /* relative divergence tolerance */
103: PetscReal abstol; /* absolute tolerance */
104: PetscReal stol; /* step length tolerance*/
105: PetscReal deltatol; /* trust region convergence tolerance */
106: PetscBool forceiteration; /* Force SNES to take at least one iteration regardless of the initial residual norm */
107: PetscInt lagpreconditioner; /* SNESSetLagPreconditioner() */
108: PetscInt lagjacobian; /* SNESSetLagJacobian() */
109: PetscInt jac_iter; /* The present iteration of the Jacobian lagging */
110: PetscBool lagjac_persist; /* The jac_iter persists until reset */
111: PetscInt pre_iter; /* The present iteration of the Preconditioner lagging */
112: PetscBool lagpre_persist; /* The pre_iter persists until reset */
113: PetscInt gridsequence; /* number of grid sequence steps to take; defaults to zero */
115: PetscBool tolerancesset; /* SNESSetTolerances() called and tolerances should persist through SNESCreate_XXX()*/
117: PetscBool vec_func_init_set; /* the initial function has been set */
119: SNESNormSchedule normschedule; /* Norm computation type for SNES instance */
120: SNESFunctionType functype; /* Function type for the SNES instance */
122: /* ------------------------ Default work-area management ---------------------- */
124: PetscInt nwork;
125: Vec *work;
127: /* ------------------------- Miscellaneous Information ------------------------ */
129: PetscInt setfromoptionscalled;
130: PetscReal *conv_hist; /* If !0, stores function norm (or
131: gradient norm) at each iteration */
132: PetscInt *conv_hist_its; /* linear iterations for each Newton step */
133: size_t conv_hist_len; /* size of convergence history array */
134: size_t conv_hist_max; /* actual amount of data in conv_history */
135: PetscBool conv_hist_reset; /* reset counter for each new SNES solve */
136: PetscBool conv_hist_alloc;
137: PetscBool counters_reset; /* reset counter for each new SNES solve */
139: /* the next two are used for failures in the line search; they should be put elsewhere */
140: PetscInt numFailures; /* number of unsuccessful step attempts */
141: PetscInt maxFailures; /* maximum number of unsuccessful step attempts */
143: PetscInt numLinearSolveFailures;
144: PetscInt maxLinearSolveFailures;
146: PetscBool domainerror; /* set with SNESSetFunctionDomainError() */
147: PetscBool jacobiandomainerror; /* set with SNESSetJacobianDomainError() */
148: PetscBool checkjacdomainerror; /* if or not check Jacobian domain error after Jacobian evaluations */
150: PetscBool ksp_ewconv; /* flag indicating use of Eisenstat-Walker KSP convergence criteria */
151: void *kspconvctx; /* Eisenstat-Walker KSP convergence context */
153: /* SNESConvergedDefault context: split it off into a separate var/struct to be passed as context to SNESConvergedDefault? */
154: PetscReal ttol; /* rtol*initial_residual_norm */
155: PetscReal rnorm0; /* initial residual norm (used for divergence testing) */
157: Vec *vwork; /* more work vectors for Jacobian approx */
158: PetscInt nvwork;
160: PetscBool mf; /* -snes_mf was used on this snes */
161: PetscBool mf_operator; /* -snes_mf_operator was used on this snes */
162: PetscInt mf_version; /* The version of snes_mf used */
164: PetscReal vizerotolerance; /* tolerance for considering an x[] value to be on the bound */
165: Vec xl, xu; /* upper and lower bounds for box constrained VI problems */
166: PetscInt ntruebounds; /* number of non-infinite bounds set for VI box constraints */
167: PetscBool usersetbounds; /* bounds have been set via SNESVISetVariableBounds(), rather than via computevariablebounds() callback. */
169: PetscBool alwayscomputesfinalresidual; /* Does SNESSolve_XXX always compute the value of the residual at the final
170: * solution and put it in vec_func? Used inside SNESSolve_FAS to determine
171: * if the final residual must be computed before restricting or prolonging
172: * it. */
173: };
175: typedef struct _p_DMSNES *DMSNES;
176: typedef struct _DMSNESOps *DMSNESOps;
177: struct _DMSNESOps {
178: PetscErrorCode (*computefunction)(SNES, Vec, Vec, void *);
179: PetscErrorCode (*computemffunction)(SNES, Vec, Vec, void *);
180: PetscErrorCode (*computejacobian)(SNES, Vec, Mat, Mat, void *);
182: /* objective */
183: PetscErrorCode (*computeobjective)(SNES, Vec, PetscReal *, void *);
185: /* Picard iteration functions */
186: PetscErrorCode (*computepfunction)(SNES, Vec, Vec, void *);
187: PetscErrorCode (*computepjacobian)(SNES, Vec, Mat, Mat, void *);
189: /* User-defined smoother */
190: PetscErrorCode (*computegs)(SNES, Vec, Vec, void *);
192: PetscErrorCode (*destroy)(DMSNES);
193: PetscErrorCode (*duplicate)(DMSNES, DMSNES);
194: };
196: struct _p_DMSNES {
197: PETSCHEADER(struct _DMSNESOps);
198: PetscContainer functionctxcontainer;
199: PetscContainer jacobianctxcontainer;
200: void *mffunctionctx;
201: void *gsctx;
202: void *pctx;
203: void *objectivectx;
205: void *data;
207: /* This is NOT reference counted. The DM on which this context was first created is cached here to implement one-way
208: * copy-on-write. When DMGetDMSNESWrite() sees a request using a different DM, it makes a copy. Thus, if a user
209: * only interacts directly with one level, e.g., using SNESSetFunction(), then SNESSetUp_FAS() is called to build
210: * coarse levels, then the user changes the routine with another call to SNESSetFunction(), it automatically
211: * propagates to all the levels. If instead, they get out a specific level and set the function on that level,
212: * subsequent changes to the original level will no longer propagate to that level.
213: */
214: DM originaldm;
215: };
216: PETSC_EXTERN PetscErrorCode DMGetDMSNES(DM, DMSNES *);
217: PETSC_EXTERN PetscErrorCode DMSNESView(DMSNES, PetscViewer);
218: PETSC_EXTERN PetscErrorCode DMSNESLoad(DMSNES, PetscViewer);
219: PETSC_EXTERN PetscErrorCode DMGetDMSNESWrite(DM, DMSNES *);
221: /* Context for Eisenstat-Walker convergence criteria for KSP solvers */
222: typedef struct {
223: PetscInt version; /* flag indicating version (1,2,3 or 4) */
224: PetscReal rtol_0; /* initial rtol */
225: PetscReal rtol_last; /* last rtol */
226: PetscReal rtol_max; /* maximum rtol */
227: PetscReal gamma; /* mult. factor for version 2 rtol computation */
228: PetscReal alpha; /* power for version 2 rtol computation */
229: PetscReal alpha2; /* power for safeguard */
230: PetscReal threshold; /* threshold for imposing safeguard */
231: PetscReal lresid_last; /* linear residual from last iteration */
232: PetscReal norm_last; /* function norm from last iteration */
233: PetscReal norm_first; /* function norm from the beginning of the first iteration. */
234: PetscReal rtol_last_2, rk_last, rk_last_2;
235: PetscReal v4_p1, v4_p2, v4_p3, v4_m1, v4_m2, v4_m3, v4_m4;
236: } SNESKSPEW;
238: static inline PetscErrorCode SNESLogConvergenceHistory(SNES snes, PetscReal res, PetscInt its)
239: {
240: PetscFunctionBegin;
241: PetscCall(PetscObjectSAWsTakeAccess((PetscObject)snes));
242: if (snes->conv_hist && snes->conv_hist_max > snes->conv_hist_len) {
243: if (snes->conv_hist) snes->conv_hist[snes->conv_hist_len] = res;
244: if (snes->conv_hist_its) snes->conv_hist_its[snes->conv_hist_len] = its;
245: snes->conv_hist_len++;
246: }
247: PetscCall(PetscObjectSAWsGrantAccess((PetscObject)snes));
248: PetscFunctionReturn(PETSC_SUCCESS);
249: }
251: PETSC_EXTERN PetscErrorCode SNESVIProjectOntoBounds(SNES, Vec);
252: PETSC_INTERN PetscErrorCode SNESVICheckLocalMin_Private(SNES, Mat, Vec, Vec, PetscReal, PetscBool *);
253: PETSC_INTERN PetscErrorCode SNESReset_VI(SNES);
254: PETSC_INTERN PetscErrorCode SNESDestroy_VI(SNES);
255: PETSC_INTERN PetscErrorCode SNESView_VI(SNES, PetscViewer);
256: PETSC_INTERN PetscErrorCode SNESSetFromOptions_VI(SNES, PetscOptionItems *);
257: PETSC_INTERN PetscErrorCode SNESSetUp_VI(SNES);
258: PETSC_EXTERN_TYPEDEF typedef PetscErrorCode (*SNESVIComputeVariableBoundsFunction)(SNES, Vec, Vec);
259: PETSC_INTERN PetscErrorCode SNESVISetComputeVariableBounds_VI(SNES, SNESVIComputeVariableBoundsFunction);
260: PETSC_INTERN PetscErrorCode SNESVISetVariableBounds_VI(SNES, Vec, Vec);
261: PETSC_INTERN PetscErrorCode SNESConvergedDefault_VI(SNES, PetscInt, PetscReal, PetscReal, PetscReal, SNESConvergedReason *, void *);
263: PETSC_EXTERN PetscErrorCode DMSNESUnsetFunctionContext_Internal(DM);
264: PETSC_EXTERN PetscErrorCode DMSNESUnsetJacobianContext_Internal(DM);
265: PETSC_EXTERN PetscErrorCode DMSNESCheck_Internal(SNES, DM, Vec);
267: PETSC_EXTERN PetscLogEvent SNES_Solve;
268: PETSC_EXTERN PetscLogEvent SNES_SetUp;
269: PETSC_EXTERN PetscLogEvent SNES_FunctionEval;
270: PETSC_EXTERN PetscLogEvent SNES_JacobianEval;
271: PETSC_EXTERN PetscLogEvent SNES_NGSEval;
272: PETSC_EXTERN PetscLogEvent SNES_NGSFuncEval;
273: PETSC_EXTERN PetscLogEvent SNES_NPCSolve;
274: PETSC_EXTERN PetscLogEvent SNES_ObjectiveEval;
276: PETSC_INTERN PetscBool SNEScite;
277: PETSC_INTERN const char SNESCitation[];
279: /* Used by TAOBNK solvers */
280: PETSC_EXTERN PetscErrorCode KSPPostSolve_SNESEW(KSP, Vec, Vec, void *);
281: PETSC_EXTERN PetscErrorCode KSPPreSolve_SNESEW(KSP, Vec, Vec, void *);
282: PETSC_EXTERN PetscErrorCode SNESEWSetFromOptions_Private(SNESKSPEW *, PetscBool, MPI_Comm, const char *);
284: /*
285: Either generate an error or mark as diverged when a real from a SNES function norm is Nan or Inf.
286: domainerror is reset here, once reason is set, to allow subsequent iterations to be feasible (e.g. line search).
287: */
288: #define SNESCheckFunctionNorm(snes, beta) \
289: do { \
290: if (PetscIsInfOrNanReal(beta)) { \
291: PetscCheck(!snes->errorifnotconverged, PetscObjectComm((PetscObject)snes), PETSC_ERR_NOT_CONVERGED, "SNESSolve has not converged due to Nan or Inf norm"); \
292: { \
293: PetscBool domainerror; \
294: PetscCall(MPIU_Allreduce(&snes->domainerror, &domainerror, 1, MPIU_BOOL, MPI_LOR, PetscObjectComm((PetscObject)snes))); \
295: if (domainerror) { \
296: snes->reason = SNES_DIVERGED_FUNCTION_DOMAIN; \
297: snes->domainerror = PETSC_FALSE; \
298: } else snes->reason = SNES_DIVERGED_FNORM_NAN; \
299: PetscFunctionReturn(PETSC_SUCCESS); \
300: } \
301: } \
302: } while (0)
304: #define SNESCheckJacobianDomainerror(snes) \
305: do { \
306: if (snes->checkjacdomainerror) { \
307: PetscBool domainerror; \
308: PetscCall(MPIU_Allreduce(&snes->jacobiandomainerror, &domainerror, 1, MPIU_BOOL, MPI_LOR, PetscObjectComm((PetscObject)snes))); \
309: if (domainerror) { \
310: snes->reason = SNES_DIVERGED_JACOBIAN_DOMAIN; \
311: snes->jacobiandomainerror = PETSC_FALSE; \
312: PetscCheck(!snes->errorifnotconverged, PetscObjectComm((PetscObject)snes), PETSC_ERR_NOT_CONVERGED, "SNESSolve has not converged due to Jacobian domain error"); \
313: PetscFunctionReturn(PETSC_SUCCESS); \
314: } \
315: } \
316: } while (0)
318: #define SNESCheckKSPSolve(snes) \
319: do { \
320: KSPConvergedReason kspreason; \
321: PetscInt lits; \
322: PetscCall(KSPGetIterationNumber(snes->ksp, &lits)); \
323: snes->linear_its += lits; \
324: PetscCall(KSPGetConvergedReason(snes->ksp, &kspreason)); \
325: if (kspreason < 0) { \
326: if (kspreason == KSP_DIVERGED_NANORINF) { \
327: PetscBool domainerror; \
328: PetscCall(MPIU_Allreduce(&snes->domainerror, &domainerror, 1, MPIU_BOOL, MPI_LOR, PetscObjectComm((PetscObject)snes))); \
329: if (domainerror) snes->reason = SNES_DIVERGED_FUNCTION_DOMAIN; \
330: else snes->reason = SNES_DIVERGED_LINEAR_SOLVE; \
331: PetscFunctionReturn(PETSC_SUCCESS); \
332: } else { \
333: if (++snes->numLinearSolveFailures >= snes->maxLinearSolveFailures) { \
334: 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)); \
335: snes->reason = SNES_DIVERGED_LINEAR_SOLVE; \
336: PetscFunctionReturn(PETSC_SUCCESS); \
337: } \
338: } \
339: } \
340: } while (0)
342: #define SNESNeedNorm_Private(snes, iter) \
343: (((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)) || \
344: (snes)->normschedule == SNES_NORM_ALWAYS)