Actual source code: taoimpl.h

  1: #pragma once

  3: #include <petsctao.h>
  4: #include <petsctaolinesearch.h>
  5: #include <petsc/private/petscimpl.h>

  7: PETSC_EXTERN PetscBool      TaoRegisterAllCalled;
  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 (*computeresidual)(Tao, Vec, Vec, void *);
 19:   PetscErrorCode (*computeresidualjacobian)(Tao, Vec, Mat, Mat, 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 *);
 29:   PetscErrorCode (*update)(Tao, PetscInt, 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)(Tao, PetscOptionItems *);
 39:   PetscErrorCode (*destroy)(Tao);
 40: };

 42: #define MAXTAOMONITORS 10

 44: struct _p_Tao {
 45:   PETSCHEADER(struct _TaoOps);
 46:   void *ctx; /* user provided context */
 47:   void *user_objP;
 48:   void *user_objgradP;
 49:   void *user_gradP;
 50:   void *user_hessP;
 51:   void *user_lsresP;
 52:   void *user_lsjacP;
 53:   void *user_conP;
 54:   void *user_con_equalityP;
 55:   void *user_con_inequalityP;
 56:   void *user_jacP;
 57:   void *user_jac_equalityP;
 58:   void *user_jac_inequalityP;
 59:   void *user_jac_stateP;
 60:   void *user_jac_designP;
 61:   void *user_boundsP;
 62:   void *user_update;

 64:   PetscErrorCode (*monitor[MAXTAOMONITORS])(Tao, void *);
 65:   PetscCtxDestroyFn *monitordestroy[MAXTAOMONITORS];
 66:   void              *monitorcontext[MAXTAOMONITORS];
 67:   PetscInt           numbermonitors;
 68:   void              *cnvP;
 69:   TaoConvergedReason reason;

 71:   PetscBool setupcalled;
 72:   void     *data;

 74:   Vec        solution;
 75:   Vec        gradient;
 76:   Vec        stepdirection;
 77:   Vec        XL;
 78:   Vec        XU;
 79:   Vec        IL;
 80:   Vec        IU;
 81:   Vec        DI;
 82:   Vec        DE;
 83:   Mat        hessian;
 84:   Mat        hessian_pre;
 85:   Mat        gradient_norm;
 86:   Vec        gradient_norm_tmp;
 87:   Vec        ls_res;
 88:   Mat        ls_jac;
 89:   Mat        ls_jac_pre;
 90:   Vec        res_weights_v;
 91:   PetscInt   res_weights_n;
 92:   PetscInt  *res_weights_rows;
 93:   PetscInt  *res_weights_cols;
 94:   PetscReal *res_weights_w;
 95:   Vec        constraints;
 96:   Vec        constraints_equality;
 97:   Vec        constraints_inequality;
 98:   Mat        jacobian;
 99:   Mat        jacobian_pre;
100:   Mat        jacobian_inequality;
101:   Mat        jacobian_inequality_pre;
102:   Mat        jacobian_equality;
103:   Mat        jacobian_equality_pre;
104:   Mat        jacobian_state;
105:   Mat        jacobian_state_inv;
106:   Mat        jacobian_design;
107:   Mat        jacobian_state_pre;
108:   Mat        jacobian_design_pre;
109:   IS         state_is;
110:   IS         design_is;
111:   PetscReal  step;
112:   PetscReal  residual;
113:   PetscReal  gnorm0;
114:   PetscReal  cnorm;
115:   PetscReal  cnorm0;
116:   PetscReal  fc;

118:   PetscInt max_constraints;
119:   PetscInt nfuncs;
120:   PetscInt ngrads;
121:   PetscInt nfuncgrads;
122:   PetscInt nhess;
123:   PetscInt niter;
124:   PetscInt ntotalits;
125:   PetscInt nconstraints;
126:   PetscInt niconstraints;
127:   PetscInt neconstraints;
128:   PetscInt njac;
129:   PetscInt njac_equality;
130:   PetscInt njac_inequality;
131:   PetscInt njac_state;
132:   PetscInt njac_design;

134:   PetscInt ksp_its;     /* KSP iterations for this solver iteration */
135:   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     trust; /* Current trust region */

142:   /* EW type forcing term */
143:   PetscBool ksp_ewconv;
144:   SNES      snes_ewdummy;

146:   PetscObjectParameterDeclare(PetscReal, gatol);
147:   PetscObjectParameterDeclare(PetscReal, grtol);
148:   PetscObjectParameterDeclare(PetscReal, gttol);
149:   PetscObjectParameterDeclare(PetscReal, catol);
150:   PetscObjectParameterDeclare(PetscReal, crtol);
151:   PetscObjectParameterDeclare(PetscReal, steptol);
152:   PetscObjectParameterDeclare(PetscReal, fmin);
153:   PetscObjectParameterDeclare(PetscInt, max_it);
154:   PetscObjectParameterDeclare(PetscInt, max_funcs);
155:   PetscObjectParameterDeclare(PetscReal, trust0); /* initial trust region radius */

157:   PetscBool printreason;
158:   PetscBool viewsolution;
159:   PetscBool viewgradient;
160:   PetscBool viewconstraints;
161:   PetscBool viewhessian;
162:   PetscBool viewjacobian;
163:   PetscBool bounded;
164:   PetscBool constrained;
165:   PetscBool eq_constrained;
166:   PetscBool ineq_constrained;
167:   PetscBool ineq_doublesided;
168:   PetscBool header_printed;
169:   PetscBool recycle;

171:   TaoSubsetType subset_type;
172:   PetscInt      hist_max;   /* Number of iteration histories to keep */
173:   PetscReal    *hist_obj;   /* obj value at each iteration */
174:   PetscReal    *hist_resid; /* residual at each iteration */
175:   PetscReal    *hist_cnorm; /* constraint norm at each iteration */
176:   PetscInt     *hist_lits;  /* number of ksp its at each TAO iteration */
177:   PetscInt      hist_len;
178:   PetscBool     hist_reset;
179:   PetscBool     hist_malloc;
180: };

182: PETSC_EXTERN PetscLogEvent TAO_Solve;
183: PETSC_EXTERN PetscLogEvent TAO_ObjectiveEval;
184: PETSC_EXTERN PetscLogEvent TAO_GradientEval;
185: PETSC_EXTERN PetscLogEvent TAO_ObjGradEval;
186: PETSC_EXTERN PetscLogEvent TAO_HessianEval;
187: PETSC_EXTERN PetscLogEvent TAO_ConstraintsEval;
188: PETSC_EXTERN PetscLogEvent TAO_JacobianEval;

190: static inline PetscErrorCode TaoLogConvergenceHistory(Tao tao, PetscReal obj, PetscReal resid, PetscReal cnorm, PetscInt totits)
191: {
192:   PetscFunctionBegin;
193:   if (tao->hist_max > tao->hist_len) {
194:     if (tao->hist_obj) tao->hist_obj[tao->hist_len] = obj;
195:     if (tao->hist_resid) tao->hist_resid[tao->hist_len] = resid;
196:     if (tao->hist_cnorm) tao->hist_cnorm[tao->hist_len] = cnorm;
197:     if (tao->hist_lits) {
198:       PetscInt sits = totits;
199:       PetscCheck(tao->hist_len >= 0, PetscObjectComm((PetscObject)tao), PETSC_ERR_ARG_WRONGSTATE, "History length cannot be negative");
200:       for (PetscInt i = 0; i < tao->hist_len; i++) sits -= tao->hist_lits[i];
201:       tao->hist_lits[tao->hist_len] = sits;
202:     }
203:     tao->hist_len++;
204:   }
205:   PetscFunctionReturn(PETSC_SUCCESS);
206: }