Actual source code: snesimpl.h
2: #ifndef __SNESIMPL_H
5: #include petscsnes.h
7: typedef struct _SNESOps *SNESOps;
9: struct _SNESOps {
10: PetscErrorCode (*computefunction)(SNES,Vec,Vec,void*);
11: PetscErrorCode (*computejacobian)(SNES,Vec,Mat*,Mat*,MatStructure*,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: };
23: /*
24: Nonlinear solver context
25: */
26: #define MAXSNESMONITORS 5
28: struct _p_SNES {
29: PETSCHEADER(struct _SNESOps);
31: /* ------------------------ User-provided stuff -------------------------------*/
32: void *user; /* user-defined context */
34: Vec vec_rhs; /* If non-null, solve F(x) = rhs */
35: Vec vec_sol; /* pointer to solution */
37: Vec vec_func; /* pointer to function */
38: void *funP; /* user-defined function context */
39:
40: Mat jacobian; /* Jacobian matrix */
41: Mat jacobian_pre; /* preconditioner matrix */
42: void *jacP; /* user-defined Jacobian context */
43: KSP ksp; /* linear solver context */
45: Vec vec_sol_update; /* pointer to solution update */
46:
47: Vec scaling; /* scaling vector */
48: void *scaP; /* scaling context */
50: /* ------------------------Time stepping hooks-----------------------------------*/
52: /* ---------------- PETSc-provided (or user-provided) stuff ---------------------*/
54: PetscErrorCode (*monitor[MAXSNESMONITORS])(SNES,PetscInt,PetscReal,void*); /* monitor routine */
55: PetscErrorCode (*monitordestroy[MAXSNESMONITORS])(void*); /* monitor context destroy routine */
56: void *monitorcontext[MAXSNESMONITORS]; /* monitor context */
57: PetscInt numbermonitors; /* number of monitors */
58: void *cnvP; /* convergence context */
59: SNESConvergedReason reason;
61: /* --- Routines and data that are unique to each particular solver --- */
63: PetscTruth setupcalled; /* true if setup has been called */
64: void *data; /* implementation-specific data */
66: /* -------------------------- Parameters -------------------------------------- */
68: PetscInt max_its; /* max number of iterations */
69: PetscInt max_funcs; /* max number of function evals */
70: PetscInt nfuncs; /* number of function evaluations */
71: PetscInt iter; /* global iteration number */
72: PetscInt linear_its; /* total number of linear solver iterations */
73: PetscReal norm; /* residual norm of current iterate */
74: PetscReal rtol; /* relative tolerance */
75: PetscReal abstol; /* absolute tolerance */
76: PetscReal xtol; /* relative tolerance in solution */
77: PetscReal deltatol; /* trust region convergence tolerance */
78: PetscTruth printreason; /* print reason for convergence/divergence after each solve */
79: PetscInt lagpreconditioner; /* SNESSetLagPreconditioner() */
80: PetscInt lagjacobian; /* SNESSetLagJacobian() */
81: /* ------------------------ Default work-area management ---------------------- */
83: PetscInt nwork;
84: Vec *work;
86: /* ------------------------- Miscellaneous Information ------------------------ */
88: PetscReal *conv_hist; /* If !0, stores function norm (or
89: gradient norm) at each iteration */
90: PetscInt *conv_hist_its; /* linear iterations for each Newton step */
91: PetscInt conv_hist_len; /* size of convergence history array */
92: PetscInt conv_hist_max; /* actual amount of data in conv_history */
93: PetscTruth conv_hist_reset; /* reset counter for each new SNES solve */
95: /* the next two are used for failures in the line search; they should be put into the LS struct */
96: PetscInt numFailures; /* number of unsuccessful step attempts */
97: PetscInt maxFailures; /* maximum number of unsuccessful step attempts */
99: PetscInt numLinearSolveFailures;
100: PetscInt maxLinearSolveFailures;
102: PetscTruth domainerror; /* set with SNESSetFunctionDomainError() */
104: PetscTruth ksp_ewconv; /* flag indicating use of Eisenstat-Walker KSP convergence criteria */
105: void *kspconvctx; /* Eisenstat-Walker KSP convergence context */
107: PetscReal ttol; /* used by default convergence test routine */
109: Vec *vwork; /* more work vectors for Jacobian approx */
110: PetscInt nvwork;
111: };
113: /* Context for Eisenstat-Walker convergence criteria for KSP solvers */
114: typedef struct {
115: PetscInt version; /* flag indicating version 1 or 2 of test */
116: PetscReal rtol_0; /* initial rtol */
117: PetscReal rtol_last; /* last rtol */
118: PetscReal rtol_max; /* maximum rtol */
119: PetscReal gamma; /* mult. factor for version 2 rtol computation */
120: PetscReal alpha; /* power for version 2 rtol computation */
121: PetscReal alpha2; /* power for safeguard */
122: PetscReal threshold; /* threshold for imposing safeguard */
123: PetscReal lresid_last; /* linear residual from last iteration */
124: PetscReal norm_last; /* function norm from last iteration */
125: } SNESKSPEW;
127: #define SNESLogConvHistory(snes,res,its) \
128: { if (snes->conv_hist && snes->conv_hist_max > snes->conv_hist_len) \
129: { if (snes->conv_hist) snes->conv_hist[snes->conv_hist_len] = res; \
130: if (snes->conv_hist_its) snes->conv_hist_its[snes->conv_hist_len] = its; \
131: snes->conv_hist_len++;\
132: }}
134: #define SNESMonitor(snes,it,rnorm) \
135: { PetscErrorCode _ierr; PetscInt _i,_im = snes->numbermonitors; \
136: for (_i=0; _i<_im; _i++) {\
137: _(*snes->monitor[_i])(snes,it,rnorm,snes->monitorcontext[_i]);CHKERRQ(_ierr); \
138: } \
139: }
141: PetscErrorCode SNES_KSPSolve(SNES,KSP,Vec,Vec);
142: PetscErrorCode SNESScaleStep_Private(SNES,Vec,PetscReal*,PetscReal*,PetscReal*,PetscReal*);
149: #endif