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