Actual source code: petsc-snesimpl.h
petsc-3.4.5 2014-06-29
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*,MatStructure*,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 gridsequence; /* number of grid sequence steps to take; defaults to zero */
100: PetscBool tolerancesset; /* SNESSetTolerances() called and tolerances should persist through SNESCreate_XXX()*/
102: PetscReal norm_init; /* the initial norm value */
103: PetscBool norm_init_set; /* the initial norm has been set */
104: PetscBool vec_func_init_set; /* the initial function has been set */
106: SNESNormType normtype; /* Norm computation type for SNES instance */
108: /* ------------------------ Default work-area management ---------------------- */
110: PetscInt nwork;
111: Vec *work;
113: /* ------------------------- Miscellaneous Information ------------------------ */
115: PetscReal *conv_hist; /* If !0, stores function norm (or
116: gradient norm) at each iteration */
117: PetscInt *conv_hist_its; /* linear iterations for each Newton step */
118: PetscInt conv_hist_len; /* size of convergence history array */
119: PetscInt conv_hist_max; /* actual amount of data in conv_history */
120: PetscBool conv_hist_reset; /* reset counter for each new SNES solve */
121: PetscBool conv_malloc;
123: /* the next two are used for failures in the line search; they should be put elsewhere */
124: PetscInt numFailures; /* number of unsuccessful step attempts */
125: PetscInt maxFailures; /* maximum number of unsuccessful step attempts */
127: PetscInt numLinearSolveFailures;
128: PetscInt maxLinearSolveFailures;
130: PetscBool domainerror; /* set with SNESSetFunctionDomainError() */
132: PetscBool ksp_ewconv; /* flag indicating use of Eisenstat-Walker KSP convergence criteria */
133: void *kspconvctx; /* Eisenstat-Walker KSP convergence context */
135: PetscReal ttol; /* used by default convergence test routine */
137: Vec *vwork; /* more work vectors for Jacobian approx */
138: PetscInt nvwork;
140: PetscBool mf; /* -snes_mf was used on this snes */
141: PetscBool mf_operator; /* -snes_mf_operator was used on this snes */
142: PetscInt mf_version; /* The version of snes_mf used */
144: Vec xl,xu; /* upper and lower bounds for box constrained VI problems */
145: PetscInt ntruebounds; /* number of non-infinite bounds set for VI box constraints */
146: PetscBool usersetbounds; /* bounds have been set via SNESVISetVariableBounds(), rather than via computevariablebounds() callback. */
147: };
149: typedef struct _p_DMSNES *DMSNES;
150: typedef struct _DMSNESOps *DMSNESOps;
151: struct _DMSNESOps {
152: PetscErrorCode (*computefunction)(SNES,Vec,Vec,void*);
153: PetscErrorCode (*computejacobian)(SNES,Vec,Mat*,Mat*,MatStructure*,void*);
155: /* objective */
156: PetscErrorCode (*computeobjective)(SNES,Vec,PetscReal*,void*);
158: /* Picard iteration functions */
159: PetscErrorCode (*computepfunction)(SNES,Vec,Vec,void*);
160: PetscErrorCode (*computepjacobian)(SNES,Vec,Mat*,Mat*,MatStructure*,void*);
162: /* User-defined smoother */
163: PetscErrorCode (*computegs)(SNES,Vec,Vec,void*);
165: PetscErrorCode (*destroy)(DMSNES);
166: PetscErrorCode (*duplicate)(DMSNES,DMSNES);
167: };
169: struct _p_DMSNES {
170: PETSCHEADER(struct _DMSNESOps);
171: void *functionctx;
172: void *gsctx;
173: void *pctx;
174: void *jacobianctx;
175: void *objectivectx;
177: void *data;
179: /* This is NOT reference counted. The DM on which this context was first created is cached here to implement one-way
180: * copy-on-write. When DMGetDMSNESWrite() sees a request using a different DM, it makes a copy. Thus, if a user
181: * only interacts directly with one level, e.g., using SNESSetFunction(), then SNESSetUp_FAS() is called to build
182: * coarse levels, then the user changes the routine with another call to SNESSetFunction(), it automatically
183: * propagates to all the levels. If instead, they get out a specific level and set the function on that level,
184: * subsequent changes to the original level will no longer propagate to that level.
185: */
186: DM originaldm;
187: };
188: PETSC_EXTERN PetscErrorCode DMGetDMSNES(DM,DMSNES*);
189: PETSC_EXTERN PetscErrorCode DMSNESView(DMSNES,PetscViewer);
190: PETSC_EXTERN PetscErrorCode DMSNESLoad(DMSNES,PetscViewer);
191: PETSC_EXTERN PetscErrorCode DMGetDMSNESWrite(DM,DMSNES*);
192: PETSC_EXTERN PetscErrorCode DMCopyDMSNES(DM,DM);
195: /* Context for Eisenstat-Walker convergence criteria for KSP solvers */
196: typedef struct {
197: PetscInt version; /* flag indicating version 1 or 2 of test */
198: PetscReal rtol_0; /* initial rtol */
199: PetscReal rtol_last; /* last rtol */
200: PetscReal rtol_max; /* maximum rtol */
201: PetscReal gamma; /* mult. factor for version 2 rtol computation */
202: PetscReal alpha; /* power for version 2 rtol computation */
203: PetscReal alpha2; /* power for safeguard */
204: PetscReal threshold; /* threshold for imposing safeguard */
205: PetscReal lresid_last; /* linear residual from last iteration */
206: PetscReal norm_last; /* function norm from last iteration */
207: } SNESKSPEW;
211: PETSC_STATIC_INLINE PetscErrorCode SNESLogConvergenceHistory(SNES snes,PetscReal res,PetscInt its)
212: {
216: PetscObjectAMSTakeAccess((PetscObject)snes);
217: if (snes->conv_hist && snes->conv_hist_max > snes->conv_hist_len) {
218: if (snes->conv_hist) snes->conv_hist[snes->conv_hist_len] = res;
219: if (snes->conv_hist_its) snes->conv_hist_its[snes->conv_hist_len] = its;
220: snes->conv_hist_len++;
221: }
222: PetscObjectAMSGrantAccess((PetscObject)snes);
223: return(0);
224: }
226: PETSC_EXTERN PetscErrorCode SNESVIProjectOntoBounds(SNES,Vec);
227: PETSC_INTERN PetscErrorCode SNESVICheckLocalMin_Private(SNES,Mat,Vec,Vec,PetscReal,PetscBool*);
228: PETSC_INTERN PetscErrorCode SNESReset_VI(SNES);
229: PETSC_INTERN PetscErrorCode SNESDestroy_VI(SNES);
230: PETSC_INTERN PetscErrorCode SNESView_VI(SNES,PetscViewer);
231: PETSC_INTERN PetscErrorCode SNESSetFromOptions_VI(SNES);
232: PETSC_INTERN PetscErrorCode SNESSetUp_VI(SNES);
233: PETSC_EXTERN_TYPEDEF typedef PetscErrorCode (*SNESVIComputeVariableBoundsFunction)(SNES,Vec,Vec);
234: PETSC_INTERN PetscErrorCode SNESVISetComputeVariableBounds_VI(SNES,SNESVIComputeVariableBoundsFunction);
235: PETSC_INTERN PetscErrorCode SNESVISetVariableBounds_VI(SNES,Vec,Vec);
236: PETSC_INTERN PetscErrorCode SNESConvergedDefault_VI(SNES,PetscInt,PetscReal,PetscReal,PetscReal,SNESConvergedReason*,void*);
238: PetscErrorCode SNESScaleStep_Private(SNES,Vec,PetscReal*,PetscReal*,PetscReal*,PetscReal*);
240: PETSC_EXTERN PetscLogEvent SNES_Solve, SNES_LineSearch, SNES_FunctionEval, SNES_JacobianEval, SNES_GSEval, SNES_NPCSolve;
242: #endif