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)(Tao, PetscOptionItems *);
 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:   /* EW type forcing term */
147:   PetscBool ksp_ewconv;
148:   SNES      snes_ewdummy;

150:   PetscReal gatol;
151:   PetscReal grtol;
152:   PetscReal gttol;
153:   PetscReal catol;
154:   PetscReal crtol;
155:   PetscReal steptol;
156:   PetscReal fmin;
157:   PetscBool max_funcs_changed;
158:   PetscBool max_it_changed;
159:   PetscBool gatol_changed;
160:   PetscBool grtol_changed;
161:   PetscBool gttol_changed;
162:   PetscBool fmin_changed;
163:   PetscBool catol_changed;
164:   PetscBool crtol_changed;
165:   PetscBool steptol_changed;
166:   PetscBool trust0_changed;
167:   PetscBool printreason;
168:   PetscBool viewsolution;
169:   PetscBool viewgradient;
170:   PetscBool viewconstraints;
171:   PetscBool viewhessian;
172:   PetscBool viewjacobian;
173:   PetscBool bounded;
174:   PetscBool constrained;
175:   PetscBool eq_constrained;
176:   PetscBool ineq_constrained;
177:   PetscBool ineq_doublesided;
178:   PetscBool header_printed;
179:   PetscBool recycle;

181:   TaoSubsetType subset_type;
182:   PetscInt      hist_max;   /* Number of iteration histories to keep */
183:   PetscReal    *hist_obj;   /* obj value at each iteration */
184:   PetscReal    *hist_resid; /* residual at each iteration */
185:   PetscReal    *hist_cnorm; /* constraint norm at each iteration */
186:   PetscInt     *hist_lits;  /* number of ksp its at each TAO iteration */
187:   PetscInt      hist_len;
188:   PetscBool     hist_reset;
189:   PetscBool     hist_malloc;
190: };

192: PETSC_EXTERN PetscLogEvent TAO_Solve;
193: PETSC_EXTERN PetscLogEvent TAO_ObjectiveEval;
194: PETSC_EXTERN PetscLogEvent TAO_GradientEval;
195: PETSC_EXTERN PetscLogEvent TAO_ObjGradEval;
196: PETSC_EXTERN PetscLogEvent TAO_HessianEval;
197: PETSC_EXTERN PetscLogEvent TAO_ConstraintsEval;
198: PETSC_EXTERN PetscLogEvent TAO_JacobianEval;

200: static inline PetscErrorCode TaoLogConvergenceHistory(Tao tao, PetscReal obj, PetscReal resid, PetscReal cnorm, PetscInt totits)
201: {
202:   PetscFunctionBegin;
203:   if (tao->hist_max > tao->hist_len) {
204:     if (tao->hist_obj) tao->hist_obj[tao->hist_len] = obj;
205:     if (tao->hist_resid) tao->hist_resid[tao->hist_len] = resid;
206:     if (tao->hist_cnorm) tao->hist_cnorm[tao->hist_len] = cnorm;
207:     if (tao->hist_lits) {
208:       PetscInt sits = totits;
209:       PetscCheck(tao->hist_len >= 0, PetscObjectComm((PetscObject)tao), PETSC_ERR_ARG_WRONGSTATE, "History length cannot be negative");
210:       for (PetscInt i = 0; i < tao->hist_len; i++) sits -= tao->hist_lits[i];
211:       tao->hist_lits[tao->hist_len] = sits;
212:     }
213:     tao->hist_len++;
214:   }
215:   PetscFunctionReturn(PETSC_SUCCESS);
216: }

218: PETSC_INTERN PetscErrorCode TaoVecGetSubVec(Vec, IS, TaoSubsetType, PetscReal, Vec *);
219: PETSC_INTERN PetscErrorCode TaoMatGetSubMat(Mat, IS, Vec, TaoSubsetType, Mat *);
220: PETSC_INTERN PetscErrorCode TaoGradientNorm(Tao, Vec, NormType, PetscReal *);
221: PETSC_INTERN PetscErrorCode TaoEstimateActiveBounds(Vec, Vec, Vec, Vec, Vec, Vec, PetscReal, PetscReal *, IS *, IS *, IS *, IS *, IS *);
222: PETSC_INTERN PetscErrorCode TaoBoundStep(Vec, Vec, Vec, IS, IS, IS, PetscReal, Vec);
223: PETSC_INTERN PetscErrorCode TaoBoundSolution(Vec, Vec, Vec, PetscReal, PetscInt *, Vec);

225: #endif