Actual source code: petsc-snesimpl.h
petsc-3.5.4 2015-05-23
2: #ifndef __SNESIMPL_H
5: #include <petscsnes.h>
6: #include <petsc-private/petscimpl.h>
8: typedef struct _SNESOps *SNESOps;
10: struct _SNESOps {
11: PetscErrorCode (*computeinitialguess)(SNES,Vec,void*);
12: PetscErrorCode (*computescaling)(Vec,Vec,void*);
13: PetscErrorCode (*update)(SNES, PetscInt); /* General purpose function for update */
14: PetscErrorCode (*converged)(SNES,PetscInt,PetscReal,PetscReal,PetscReal,SNESConvergedReason*,void*);
15: PetscErrorCode (*convergeddestroy)(void*);
16: PetscErrorCode (*setup)(SNES); /* routine to set up the nonlinear solver */
17: PetscErrorCode (*solve)(SNES); /* actual nonlinear solver */
18: PetscErrorCode (*view)(SNES,PetscViewer);
19: PetscErrorCode (*setfromoptions)(SNES); /* sets options from database */
20: PetscErrorCode (*destroy)(SNES);
21: PetscErrorCode (*reset)(SNES);
22: PetscErrorCode (*usercompute)(SNES,void**);
23: PetscErrorCode (*userdestroy)(void**);
24: PetscErrorCode (*computevariablebounds)(SNES,Vec,Vec); /* user provided routine to set box constrained variable bounds */
25: PetscErrorCode (*computepfunction)(SNES,Vec,Vec,void*);
26: PetscErrorCode (*computepjacobian)(SNES,Vec,Mat,Mat,void*);
27: PetscErrorCode (*load)(SNES,PetscViewer);
28: };
30: /*
31: Nonlinear solver context
32: */
33: #define MAXSNESMONITORS 5
35: struct _p_SNES {
36: PETSCHEADER(struct _SNESOps);
37: DM dm;
38: PetscBool dmAuto; /* SNES created currently used DM automatically */
39: SNES pc;
40: PCSide pcside;
41: PetscBool usespc;
43: /* ------------------------ User-provided stuff -------------------------------*/
44: void *user; /* user-defined context */
46: Vec vec_rhs; /* If non-null, solve F(x) = rhs */
47: Vec vec_sol; /* pointer to solution */
49: Vec vec_func; /* pointer to function */
51: Mat jacobian; /* Jacobian matrix */
52: Mat jacobian_pre; /* preconditioner matrix */
53: void *initialguessP; /* user-defined initial guess context */
54: KSP ksp; /* linear solver context */
55: SNESLineSearch linesearch; /* line search context */
56: PetscBool usesksp;
57: MatStructure matstruct; /* Used by Picard solver */
59: Vec vec_sol_update; /* pointer to solution update */
61: Vec scaling; /* scaling vector */
62: void *scaP; /* scaling context */
64: PetscReal precheck_picard_angle; /* For use with SNESLineSearchPreCheckPicard */
66: /* ------------------------Time stepping hooks-----------------------------------*/
68: /* ---------------- PETSc-provided (or user-provided) stuff ---------------------*/
70: PetscErrorCode (*monitor[MAXSNESMONITORS])(SNES,PetscInt,PetscReal,void*); /* monitor routine */
71: PetscErrorCode (*monitordestroy[MAXSNESMONITORS])(void**); /* monitor context destroy routine */
72: void *monitorcontext[MAXSNESMONITORS]; /* monitor context */
73: PetscInt numbermonitors; /* number of monitors */
74: void *cnvP; /* convergence context */
75: SNESConvergedReason reason;
76: PetscBool errorifnotconverged;
78: /* --- Routines and data that are unique to each particular solver --- */
80: PetscBool setupcalled; /* true if setup has been called */
81: void *data; /* implementation-specific data */
83: /* -------------------------- Parameters -------------------------------------- */
85: PetscInt max_its; /* max number of iterations */
86: PetscInt max_funcs; /* max number of function evals */
87: PetscInt nfuncs; /* number of function evaluations */
88: PetscInt iter; /* global iteration number */
89: PetscInt linear_its; /* total number of linear solver iterations */
90: PetscReal norm; /* residual norm of current iterate */
91: PetscReal rtol; /* relative tolerance */
92: PetscReal abstol; /* absolute tolerance */
93: PetscReal stol; /* step length tolerance*/
94: PetscReal deltatol; /* trust region convergence tolerance */
95: PetscBool printreason; /* print reason for convergence/divergence after each solve */
96: PetscInt lagpreconditioner; /* SNESSetLagPreconditioner() */
97: PetscInt lagjacobian; /* SNESSetLagJacobian() */
98: PetscInt jac_iter; /* The present iteration of the Jacobian lagging */
99: PetscBool lagjac_persist; /* The jac_iter persists until reset */
100: PetscInt pre_iter; /* The present iteration of the Preconditioner lagging */
101: PetscBool lagpre_persist; /* The pre_iter persists until reset */
102: PetscInt gridsequence; /* number of grid sequence steps to take; defaults to zero */
104: PetscBool tolerancesset; /* SNESSetTolerances() called and tolerances should persist through SNESCreate_XXX()*/
106: PetscBool vec_func_init_set; /* the initial function has been set */
108: SNESNormSchedule normschedule; /* Norm computation type for SNES instance */
109: SNESFunctionType functype; /* Function type for the SNES instance */
111: /* ------------------------ Default work-area management ---------------------- */
113: PetscInt nwork;
114: Vec *work;
116: /* ------------------------- Miscellaneous Information ------------------------ */
118: PetscReal *conv_hist; /* If !0, stores function norm (or
119: gradient norm) at each iteration */
120: PetscInt *conv_hist_its; /* linear iterations for each Newton step */
121: PetscInt conv_hist_len; /* size of convergence history array */
122: PetscInt conv_hist_max; /* actual amount of data in conv_history */
123: PetscBool conv_hist_reset; /* reset counter for each new SNES solve */
124: PetscBool conv_malloc;
126: PetscBool counters_reset; /* reset counter for each new SNES solve */
128: /* the next two are used for failures in the line search; they should be put elsewhere */
129: PetscInt numFailures; /* number of unsuccessful step attempts */
130: PetscInt maxFailures; /* maximum number of unsuccessful step attempts */
132: PetscInt numLinearSolveFailures;
133: PetscInt maxLinearSolveFailures;
135: PetscBool domainerror; /* set with SNESSetFunctionDomainError() */
137: PetscBool ksp_ewconv; /* flag indicating use of Eisenstat-Walker KSP convergence criteria */
138: void *kspconvctx; /* Eisenstat-Walker KSP convergence context */
140: /* SNESConvergedDefault context: split it off into a separate var/struct to be passed as context to SNESConvergedDefault? */
141: PetscReal ttol; /* rtol*initial_residual_norm */
143: Vec *vwork; /* more work vectors for Jacobian approx */
144: PetscInt nvwork;
146: PetscBool mf; /* -snes_mf was used on this snes */
147: PetscBool mf_operator; /* -snes_mf_operator was used on this snes */
148: PetscInt mf_version; /* The version of snes_mf used */
150: Vec xl,xu; /* upper and lower bounds for box constrained VI problems */
151: PetscInt ntruebounds; /* number of non-infinite bounds set for VI box constraints */
152: PetscBool usersetbounds; /* bounds have been set via SNESVISetVariableBounds(), rather than via computevariablebounds() callback. */
154: };
156: typedef struct _p_DMSNES *DMSNES;
157: typedef struct _DMSNESOps *DMSNESOps;
158: struct _DMSNESOps {
159: PetscErrorCode (*computefunction)(SNES,Vec,Vec,void*);
160: PetscErrorCode (*computejacobian)(SNES,Vec,Mat,Mat,void*);
162: /* objective */
163: PetscErrorCode (*computeobjective)(SNES,Vec,PetscReal*,void*);
165: /* Picard iteration functions */
166: PetscErrorCode (*computepfunction)(SNES,Vec,Vec,void*);
167: PetscErrorCode (*computepjacobian)(SNES,Vec,Mat,Mat,void*);
169: /* User-defined smoother */
170: PetscErrorCode (*computegs)(SNES,Vec,Vec,void*);
172: PetscErrorCode (*destroy)(DMSNES);
173: PetscErrorCode (*duplicate)(DMSNES,DMSNES);
174: };
176: struct _p_DMSNES {
177: PETSCHEADER(struct _DMSNESOps);
178: void *functionctx;
179: void *gsctx;
180: void *pctx;
181: void *jacobianctx;
182: void *objectivectx;
184: void *data;
186: /* This is NOT reference counted. The DM on which this context was first created is cached here to implement one-way
187: * copy-on-write. When DMGetDMSNESWrite() sees a request using a different DM, it makes a copy. Thus, if a user
188: * only interacts directly with one level, e.g., using SNESSetFunction(), then SNESSetUp_FAS() is called to build
189: * coarse levels, then the user changes the routine with another call to SNESSetFunction(), it automatically
190: * propagates to all the levels. If instead, they get out a specific level and set the function on that level,
191: * subsequent changes to the original level will no longer propagate to that level.
192: */
193: DM originaldm;
194: };
195: PETSC_EXTERN PetscErrorCode DMGetDMSNES(DM,DMSNES*);
196: PETSC_EXTERN PetscErrorCode DMSNESView(DMSNES,PetscViewer);
197: PETSC_EXTERN PetscErrorCode DMSNESLoad(DMSNES,PetscViewer);
198: PETSC_EXTERN PetscErrorCode DMGetDMSNESWrite(DM,DMSNES*);
199: PETSC_EXTERN PetscErrorCode DMCopyDMSNES(DM,DM);
202: /* Context for Eisenstat-Walker convergence criteria for KSP solvers */
203: typedef struct {
204: PetscInt version; /* flag indicating version 1 or 2 of test */
205: PetscReal rtol_0; /* initial rtol */
206: PetscReal rtol_last; /* last rtol */
207: PetscReal rtol_max; /* maximum rtol */
208: PetscReal gamma; /* mult. factor for version 2 rtol computation */
209: PetscReal alpha; /* power for version 2 rtol computation */
210: PetscReal alpha2; /* power for safeguard */
211: PetscReal threshold; /* threshold for imposing safeguard */
212: PetscReal lresid_last; /* linear residual from last iteration */
213: PetscReal norm_last; /* function norm from last iteration */
214: PetscReal norm_first; /* function norm from the beginning of the first iteration. */
215: } SNESKSPEW;
219: PETSC_STATIC_INLINE PetscErrorCode SNESLogConvergenceHistory(SNES snes,PetscReal res,PetscInt its)
220: {
224: PetscObjectSAWsTakeAccess((PetscObject)snes);
225: if (snes->conv_hist && snes->conv_hist_max > snes->conv_hist_len) {
226: if (snes->conv_hist) snes->conv_hist[snes->conv_hist_len] = res;
227: if (snes->conv_hist_its) snes->conv_hist_its[snes->conv_hist_len] = its;
228: snes->conv_hist_len++;
229: }
230: PetscObjectSAWsGrantAccess((PetscObject)snes);
231: return(0);
232: }
234: PETSC_EXTERN PetscErrorCode SNESVIProjectOntoBounds(SNES,Vec);
235: PETSC_INTERN PetscErrorCode SNESVICheckLocalMin_Private(SNES,Mat,Vec,Vec,PetscReal,PetscBool*);
236: PETSC_INTERN PetscErrorCode SNESReset_VI(SNES);
237: PETSC_INTERN PetscErrorCode SNESDestroy_VI(SNES);
238: PETSC_INTERN PetscErrorCode SNESView_VI(SNES,PetscViewer);
239: PETSC_INTERN PetscErrorCode SNESSetFromOptions_VI(SNES);
240: PETSC_INTERN PetscErrorCode SNESSetUp_VI(SNES);
241: PETSC_EXTERN_TYPEDEF typedef PetscErrorCode (*SNESVIComputeVariableBoundsFunction)(SNES,Vec,Vec);
242: PETSC_INTERN PetscErrorCode SNESVISetComputeVariableBounds_VI(SNES,SNESVIComputeVariableBoundsFunction);
243: PETSC_INTERN PetscErrorCode SNESVISetVariableBounds_VI(SNES,Vec,Vec);
244: PETSC_INTERN PetscErrorCode SNESConvergedDefault_VI(SNES,PetscInt,PetscReal,PetscReal,PetscReal,SNESConvergedReason*,void*);
246: PetscErrorCode SNESScaleStep_Private(SNES,Vec,PetscReal*,PetscReal*,PetscReal*,PetscReal*);
248: PETSC_EXTERN PetscLogEvent SNES_Solve, SNES_LineSearch, SNES_FunctionEval, SNES_JacobianEval, SNES_NGSEval, SNES_NGSFuncEval, SNES_NPCSolve;
250: extern PetscBool SNEScite;
251: extern const char SNESCitation[];
253: #endif