Actual source code: taoimpl.h
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 (*computeresidual)(Tao,Vec,Vec,void*);
20: PetscErrorCode (*computeresidualjacobian)(Tao,Vec,Mat,Mat,void*);
21: PetscErrorCode (*computeconstraints)(Tao,Vec,Vec,void*);
22: PetscErrorCode (*computeinequalityconstraints)(Tao,Vec,Vec,void*);
23: PetscErrorCode (*computeequalityconstraints)(Tao,Vec,Vec,void*);
24: PetscErrorCode (*computejacobian)(Tao,Vec,Mat,Mat,void*);
25: PetscErrorCode (*computejacobianstate)(Tao,Vec,Mat,Mat,Mat,void*);
26: PetscErrorCode (*computejacobiandesign)(Tao,Vec,Mat,void*);
27: PetscErrorCode (*computejacobianinequality)(Tao,Vec,Mat,Mat,void*);
28: PetscErrorCode (*computejacobianequality)(Tao,Vec,Mat,Mat,void*);
29: PetscErrorCode (*computebounds)(Tao,Vec,Vec,void*);
30: PetscErrorCode (*update)(Tao,PetscInt,void*);
31: PetscErrorCode (*convergencetest)(Tao,void*);
32: PetscErrorCode (*convergencedestroy)(void*);
34: /* Methods set by solver */
35: PetscErrorCode (*computedual)(Tao,Vec,Vec);
36: PetscErrorCode (*setup)(Tao);
37: PetscErrorCode (*solve)(Tao);
38: PetscErrorCode (*view)(Tao,PetscViewer);
39: PetscErrorCode (*setfromoptions)(PetscOptionItems*,Tao);
40: PetscErrorCode (*destroy)(Tao);
41: };
43: #define MAXTAOMONITORS 10
45: struct _p_Tao {
46: PETSCHEADER(struct _TaoOps);
47: void *user;
48: void *user_objP;
49: void *user_objgradP;
50: void *user_gradP;
51: void *user_hessP;
52: void *user_lsresP;
53: void *user_lsjacP;
54: void *user_conP;
55: void *user_con_equalityP;
56: void *user_con_inequalityP;
57: void *user_jacP;
58: void *user_jac_equalityP;
59: void *user_jac_inequalityP;
60: void *user_jac_stateP;
61: void *user_jac_designP;
62: void *user_boundsP;
63: void *user_update;
65: PetscErrorCode (*monitor[MAXTAOMONITORS])(Tao,void*);
66: PetscErrorCode (*monitordestroy[MAXTAOMONITORS])(void**);
67: void *monitorcontext[MAXTAOMONITORS];
68: PetscInt numbermonitors;
69: void *cnvP;
70: TaoConvergedReason reason;
72: PetscBool setupcalled;
73: void *data;
75: Vec solution;
76: Vec gradient;
77: Vec stepdirection;
78: Vec XL;
79: Vec XU;
80: Vec IL;
81: Vec IU;
82: Vec DI;
83: Vec DE;
84: Mat hessian;
85: Mat hessian_pre;
86: Mat gradient_norm;
87: Vec gradient_norm_tmp;
88: Vec ls_res;
89: Mat ls_jac;
90: Mat ls_jac_pre;
91: Vec res_weights_v;
92: PetscInt res_weights_n;
93: PetscInt *res_weights_rows;
94: PetscInt *res_weights_cols;
95: PetscReal *res_weights_w;
96: Vec constraints;
97: Vec constraints_equality;
98: Vec constraints_inequality;
99: Mat jacobian;
100: Mat jacobian_pre;
101: Mat jacobian_inequality;
102: Mat jacobian_inequality_pre;
103: Mat jacobian_equality;
104: Mat jacobian_equality_pre;
105: Mat jacobian_state;
106: Mat jacobian_state_inv;
107: Mat jacobian_design;
108: Mat jacobian_state_pre;
109: Mat jacobian_design_pre;
110: IS state_is;
111: IS design_is;
112: PetscReal step;
113: PetscReal residual;
114: PetscReal gnorm0;
115: PetscReal cnorm;
116: PetscReal cnorm0;
117: PetscReal fc;
119: PetscInt max_it;
120: PetscInt max_funcs;
121: PetscInt max_constraints;
122: PetscInt nfuncs;
123: PetscInt ngrads;
124: PetscInt nfuncgrads;
125: PetscInt nhess;
126: PetscInt niter;
127: PetscInt ntotalits;
128: PetscInt nconstraints;
129: PetscInt niconstraints;
130: PetscInt neconstraints;
131: PetscInt njac;
132: PetscInt njac_equality;
133: PetscInt njac_inequality;
134: PetscInt njac_state;
135: PetscInt njac_design;
137: PetscInt ksp_its; /* KSP iterations for this solver iteration */
138: PetscInt ksp_tot_its; /* Total (cumulative) KSP iterations */
140: TaoLineSearch linesearch;
141: PetscBool lsflag; /* goes up when line search fails */
142: KSP ksp;
143: PetscReal trust0; /* initial trust region radius */
144: PetscReal trust; /* Current trust region */
146: PetscReal gatol;
147: PetscReal grtol;
148: PetscReal gttol;
149: PetscReal catol;
150: PetscReal crtol;
151: PetscReal steptol;
152: PetscReal fmin;
153: PetscBool max_funcs_changed;
154: PetscBool max_it_changed;
155: PetscBool gatol_changed;
156: PetscBool grtol_changed;
157: PetscBool gttol_changed;
158: PetscBool fmin_changed;
159: PetscBool catol_changed;
160: PetscBool crtol_changed;
161: PetscBool steptol_changed;
162: PetscBool trust0_changed;
163: PetscBool printreason;
164: PetscBool viewsolution;
165: PetscBool viewgradient;
166: PetscBool viewconstraints;
167: PetscBool viewhessian;
168: PetscBool viewjacobian;
169: PetscBool bounded;
170: PetscBool constrained;
171: PetscBool eq_constrained;
172: PetscBool ineq_constrained;
173: PetscBool ineq_doublesided;
174: PetscBool header_printed;
175: PetscBool recycle;
177: TaoSubsetType subset_type;
178: PetscInt hist_max; /* Number of iteration histories to keep */
179: PetscReal *hist_obj; /* obj value at each iteration */
180: PetscReal *hist_resid; /* residual at each iteration */
181: PetscReal *hist_cnorm; /* constraint norm at each iteration */
182: PetscInt *hist_lits; /* number of ksp its at each TAO iteration */
183: PetscInt hist_len;
184: PetscBool hist_reset;
185: PetscBool hist_malloc;
186: };
188: PETSC_EXTERN PetscLogEvent TAO_Solve;
189: PETSC_EXTERN PetscLogEvent TAO_ObjectiveEval;
190: PETSC_EXTERN PetscLogEvent TAO_GradientEval;
191: PETSC_EXTERN PetscLogEvent TAO_ObjGradEval;
192: PETSC_EXTERN PetscLogEvent TAO_HessianEval;
193: PETSC_EXTERN PetscLogEvent TAO_ConstraintsEval;
194: PETSC_EXTERN PetscLogEvent TAO_JacobianEval;
196: static inline PetscErrorCode TaoLogConvergenceHistory(Tao tao, PetscReal obj, PetscReal resid, PetscReal cnorm, PetscInt totits)
197: {
198: if (tao->hist_max > tao->hist_len) {
199: if (tao->hist_obj) tao->hist_obj[tao->hist_len] = obj;
200: if (tao->hist_resid) tao->hist_resid[tao->hist_len] = resid;
201: if (tao->hist_cnorm) tao->hist_cnorm[tao->hist_len] = cnorm;
202: if (tao->hist_lits) {
203: if (tao->hist_len <= 0) tao->hist_lits[0] = totits;
204: else tao->hist_lits[tao->hist_len] = totits - tao->hist_lits[tao->hist_len-1];
205: }
206: tao->hist_len++;
207: }
208: return 0;
209: }
211: PETSC_INTERN PetscErrorCode TaoVecGetSubVec(Vec,IS,TaoSubsetType,PetscReal,Vec*);
212: PETSC_INTERN PetscErrorCode TaoMatGetSubMat(Mat,IS,Vec,TaoSubsetType,Mat*);
213: PETSC_INTERN PetscErrorCode TaoGradientNorm(Tao,Vec,NormType,PetscReal*);
214: PETSC_INTERN PetscErrorCode TaoEstimateActiveBounds(Vec,Vec,Vec,Vec,Vec,Vec,PetscReal,PetscReal*,IS*,IS*,IS*,IS*,IS*);
215: PETSC_INTERN PetscErrorCode TaoBoundStep(Vec,Vec,Vec,IS,IS,IS,PetscReal,Vec);
216: PETSC_INTERN PetscErrorCode TaoBoundSolution(Vec,Vec,Vec,PetscReal,PetscInt*,Vec);
218: #endif