Actual source code: taoimpl.h
petsc-3.7.6 2017-04-24
1: #ifndef __TAO_IMPL_H
4: #include <petsctaolinesearch.h>
5: #include <petsc/private/petscimpl.h>
6: #include <petscksp.h>
8: PETSC_EXTERN PetscErrorCode TaoRegisterAll(void);
10: typedef struct _TaoOps *TaoOps;
12: struct _TaoOps {
13: /* Methods set by application */
14: PetscErrorCode (*computeobjective)(Tao, Vec, PetscReal*, void*);
15: PetscErrorCode (*computeobjectiveandgradient)(Tao, Vec, PetscReal*, Vec, void*);
16: PetscErrorCode (*computegradient)(Tao, Vec, Vec, void*);
17: PetscErrorCode (*computehessian)(Tao, Vec, Mat, Mat, void*);
18: PetscErrorCode (*computeseparableobjective)(Tao, Vec, Vec, void*);
19: PetscErrorCode (*computeconstraints)(Tao, Vec, Vec, void*);
20: PetscErrorCode (*computeinequalityconstraints)(Tao, Vec, Vec, void*);
21: PetscErrorCode (*computeequalityconstraints)(Tao, Vec, Vec, void*);
22: PetscErrorCode (*computejacobian)(Tao, Vec, Mat, Mat, void*);
23: PetscErrorCode (*computejacobianstate)(Tao, Vec, Mat, Mat, Mat, void*);
24: PetscErrorCode (*computejacobiandesign)(Tao, Vec, Mat, void*);
25: PetscErrorCode (*computejacobianinequality)(Tao, Vec, Mat, Mat, void*);
26: PetscErrorCode (*computejacobianequality)(Tao, Vec, Mat, Mat, void*);
27: PetscErrorCode (*computebounds)(Tao, Vec, Vec, void*);
29: PetscErrorCode (*convergencetest)(Tao,void*);
30: PetscErrorCode (*convergencedestroy)(void*);
32: /* Methods set by solver */
33: PetscErrorCode (*computedual)(Tao, Vec, Vec);
34: PetscErrorCode (*setup)(Tao);
35: PetscErrorCode (*solve)(Tao);
36: PetscErrorCode (*view)(Tao, PetscViewer);
37: PetscErrorCode (*setfromoptions)(PetscOptionItems*,Tao);
38: PetscErrorCode (*destroy)(Tao);
39: };
41: #define MAXTAOMONITORS 10
43: struct _p_Tao {
44: PETSCHEADER(struct _TaoOps);
45: void *user;
46: void *user_objP;
47: void *user_objgradP;
48: void *user_gradP;
49: void *user_hessP;
50: void *user_sepobjP;
51: void *user_conP;
52: void *user_con_equalityP;
53: void *user_con_inequalityP;
54: void *user_jacP;
55: void *user_jac_equalityP;
56: void *user_jac_inequalityP;
57: void *user_jac_stateP;
58: void *user_jac_designP;
59: void *user_boundsP;
61: PetscErrorCode (*monitor[MAXTAOMONITORS])(Tao,void*);
62: PetscErrorCode (*monitordestroy[MAXTAOMONITORS])(void**);
63: void *monitorcontext[MAXTAOMONITORS];
64: PetscInt numbermonitors;
65: void *cnvP;
66: TaoConvergedReason reason;
68: PetscBool setupcalled;
69: void *data;
71: Vec solution;
72: Vec gradient;
73: Vec stepdirection;
74: Vec XL;
75: Vec XU;
76: Vec IL;
77: Vec IU;
78: Vec DI;
79: Vec DE;
80: Mat hessian;
81: Mat hessian_pre;
82: Mat gradient_norm;
83: Vec gradient_norm_tmp;
84: Vec sep_objective;
85: Vec sep_weights_v;
86: PetscInt sep_weights_n;
87: PetscInt *sep_weights_rows;
88: PetscInt *sep_weights_cols;
89: PetscReal *sep_weights_w;
90: Vec constraints;
91: Vec constraints_equality;
92: Vec constraints_inequality;
93: Mat jacobian;
94: Mat jacobian_pre;
95: Mat jacobian_inequality;
96: Mat jacobian_inequality_pre;
97: Mat jacobian_equality;
98: Mat jacobian_equality_pre;
99: Mat jacobian_state;
100: Mat jacobian_state_inv;
101: Mat jacobian_design;
102: Mat jacobian_state_pre;
103: Mat jacobian_design_pre;
104: IS state_is;
105: IS design_is;
106: PetscReal step;
107: PetscReal residual;
108: PetscReal gnorm0;
109: PetscReal cnorm;
110: PetscReal cnorm0;
111: PetscReal fc;
114: PetscInt max_it;
115: PetscInt max_funcs;
116: PetscInt max_constraints;
117: PetscInt nfuncs;
118: PetscInt ngrads;
119: PetscInt nfuncgrads;
120: PetscInt nhess;
121: PetscInt niter;
122: PetscInt ntotalits;
123: PetscInt nconstraints;
124: PetscInt niconstraints;
125: PetscInt neconstraints;
126: PetscInt njac;
127: PetscInt njac_equality;
128: PetscInt njac_inequality;
129: PetscInt njac_state;
130: PetscInt njac_design;
132: PetscInt ksp_its; /* KSP iterations for this solver iteration */
133: PetscInt ksp_tot_its; /* Total (cumulative) KSP iterations */
136: TaoLineSearch linesearch;
137: PetscBool lsflag; /* goes up when line search fails */
138: KSP ksp;
139: PetscReal trust0; /* initial trust region radius */
140: PetscReal trust; /* Current trust region */
142: PetscReal gatol;
143: PetscReal grtol;
144: PetscReal gttol;
145: PetscReal catol;
146: PetscReal crtol;
147: PetscReal steptol;
148: PetscReal fmin;
149: PetscBool max_funcs_changed;
150: PetscBool max_it_changed;
151: PetscBool gatol_changed;
152: PetscBool grtol_changed;
153: PetscBool gttol_changed;
154: PetscBool fmin_changed;
155: PetscBool catol_changed;
156: PetscBool crtol_changed;
157: PetscBool steptol_changed;
158: PetscBool trust0_changed;
159: PetscBool printreason;
160: PetscBool viewsolution;
161: PetscBool viewgradient;
162: PetscBool viewconstraints;
163: PetscBool viewhessian;
164: PetscBool viewjacobian;
166: TaoSubsetType subset_type;
167: PetscInt hist_max;/* Number of iteration histories to keep */
168: PetscReal *hist_obj; /* obj value at each iteration */
169: PetscReal *hist_resid; /* residual at each iteration */
170: PetscReal *hist_cnorm; /* constraint norm at each iteration */
171: PetscInt *hist_lits; /* number of ksp its at each TAO iteration */
172: PetscInt hist_len;
173: PetscBool hist_reset;
174: PetscBool hist_malloc;
175: };
177: extern PetscLogEvent Tao_Solve, Tao_ObjectiveEval, Tao_ObjGradientEval, Tao_GradientEval, Tao_HessianEval, Tao_ConstraintsEval, Tao_JacobianEval;
181: PETSC_STATIC_INLINE PetscErrorCode TaoLogConvergenceHistory(Tao tao, PetscReal obj, PetscReal resid, PetscReal cnorm, PetscInt totits)
182: {
184: if (tao->hist_max > tao->hist_len) {
185: if (tao->hist_obj) tao->hist_obj[tao->hist_len]=obj;
186: if (tao->hist_resid) tao->hist_resid[tao->hist_len]=resid;
187: if (tao->hist_cnorm) tao->hist_cnorm[tao->hist_len]=cnorm;
188: if (tao->hist_lits) {
189: if (tao->hist_len <= 0) {
190: tao->hist_lits[0] = totits;
191: } else {
192: tao->hist_lits[tao->hist_len]=totits - tao->hist_lits[tao->hist_len-1];
193: }
194: }
195: tao->hist_len++;
196: }
197: return(0);
198: }
200: PETSC_INTERN PetscErrorCode TaoVecGetSubVec(Vec, IS, TaoSubsetType, PetscReal, Vec*);
201: PETSC_INTERN PetscErrorCode TaoMatGetSubMat(Mat, IS, Vec, TaoSubsetType, Mat*);
202: PETSC_INTERN PetscErrorCode TaoGradientNorm(Tao, Vec, NormType, PetscReal*);
204: #endif