1 /* 2 Context for bounded Newton-Krylov type optimization algorithms 3 */ 4 5 #if !defined(__TAO_BNK_H) 6 #define __TAO_BNK_H 7 #include <petsc/private/taoimpl.h> 8 #include <../src/tao/bound/impls/bncg/bncg.h> 9 10 typedef struct { 11 /* Embedded TAOBNCG */ 12 Tao bncg; 13 TAO_BNCG *bncg_ctx; 14 PetscInt max_cg_its, tot_cg_its; 15 Vec bncg_sol; 16 17 /* Allocated vectors */ 18 Vec W, Xwork, Gwork, Xold, Gold; 19 Vec unprojected_gradient, unprojected_gradient_old; 20 Vec Diag_min, Diag_max; 21 22 /* Unallocated matrices and vectors */ 23 Mat H_inactive, Hpre_inactive, M, M_inactive; 24 Vec Diag, Diag_red, X_inactive, G_inactive, inactive_work, active_work; 25 IS inactive_idx, active_idx, active_lower, active_upper, active_fixed; 26 27 /* Scalar values for the solution and step */ 28 PetscReal fold, f, gnorm, dnorm; 29 30 /* Parameters for active set estimation */ 31 PetscReal as_tol; 32 PetscReal as_step; 33 34 /* Parameters when updating the perturbation added to the Hessian matrix 35 according to the following scheme: 36 37 pert = sval; 38 39 do until convergence 40 shift Hessian by pert 41 solve Newton system 42 43 if (linear solver failed or did not compute a descent direction) 44 use steepest descent direction and increase perturbation 45 46 if (0 == pert) 47 initialize perturbation 48 pert = min(imax, max(imin, imfac * norm(G))) 49 else 50 increase perturbation 51 pert = min(pmax, max(pgfac * pert, pmgfac * norm(G))) 52 fi 53 else 54 use linear solver direction and decrease perturbation 55 56 pert = min(psfac * pert, pmsfac * norm(G)) 57 if (pert < pmin) 58 pert = 0 59 fi 60 fi 61 62 perform line search 63 function and gradient evaluation 64 check convergence 65 od 66 */ 67 PetscReal sval; /* Starting perturbation value, default zero */ 68 69 PetscReal imin; /* Minimum perturbation added during initialization */ 70 PetscReal imax; /* Maximum perturbation added during initialization */ 71 PetscReal imfac; /* Merit function factor during initialization */ 72 73 PetscReal pert; /* Current perturbation value */ 74 PetscReal pmin; /* Minimim perturbation value */ 75 PetscReal pmax; /* Maximum perturbation value */ 76 PetscReal pgfac; /* Perturbation growth factor */ 77 PetscReal psfac; /* Perturbation shrink factor */ 78 PetscReal pmgfac; /* Merit function growth factor */ 79 PetscReal pmsfac; /* Merit function shrink factor */ 80 81 /* Parameters when updating the trust-region radius based on steplength 82 if step < nu1 (very bad step) 83 radius = omega1 * min(norm(d), radius) 84 elif step < nu2 (bad step) 85 radius = omega2 * min(norm(d), radius) 86 elif step < nu3 (okay step) 87 radius = omega3 * radius; 88 elif step < nu4 (good step) 89 radius = max(omega4 * norm(d), radius) 90 else (very good step) 91 radius = max(omega5 * norm(d), radius) 92 fi 93 */ 94 PetscReal nu1; /* used to compute trust-region radius */ 95 PetscReal nu2; /* used to compute trust-region radius */ 96 PetscReal nu3; /* used to compute trust-region radius */ 97 PetscReal nu4; /* used to compute trust-region radius */ 98 99 PetscReal omega1; /* factor used for trust-region update */ 100 PetscReal omega2; /* factor used for trust-region update */ 101 PetscReal omega3; /* factor used for trust-region update */ 102 PetscReal omega4; /* factor used for trust-region update */ 103 PetscReal omega5; /* factor used for trust-region update */ 104 105 /* Parameters when updating the trust-region radius based on reduction 106 107 kappa = ared / pred 108 if kappa < eta1 (very bad step) 109 radius = alpha1 * min(norm(d), radius) 110 elif kappa < eta2 (bad step) 111 radius = alpha2 * min(norm(d), radius) 112 elif kappa < eta3 (okay step) 113 radius = alpha3 * radius; 114 elif kappa < eta4 (good step) 115 radius = max(alpha4 * norm(d), radius) 116 else (very good step) 117 radius = max(alpha5 * norm(d), radius) 118 fi 119 */ 120 PetscReal eta1; /* used to compute trust-region radius */ 121 PetscReal eta2; /* used to compute trust-region radius */ 122 PetscReal eta3; /* used to compute trust-region radius */ 123 PetscReal eta4; /* used to compute trust-region radius */ 124 125 PetscReal alpha1; /* factor used for trust-region update */ 126 PetscReal alpha2; /* factor used for trust-region update */ 127 PetscReal alpha3; /* factor used for trust-region update */ 128 PetscReal alpha4; /* factor used for trust-region update */ 129 PetscReal alpha5; /* factor used for trust-region update */ 130 131 /* Parameters when updating the trust-region radius based on interpolation 132 133 kappa = ared / pred 134 if kappa >= 1.0 - mu1 (very good step) 135 choose tau in [gamma3, gamma4] 136 radius = max(tau * norm(d), radius) 137 elif kappa >= 1.0 - mu2 (good step) 138 choose tau in [gamma2, gamma3] 139 if (tau >= 1.0) 140 radius = max(tau * norm(d), radius) 141 else 142 radius = tau * min(norm(d), radius) 143 fi 144 else (bad step) 145 choose tau in [gamma1, 1.0] 146 radius = tau * min(norm(d), radius) 147 fi 148 */ 149 PetscReal mu1; /* used for model agreement in interpolation */ 150 PetscReal mu2; /* used for model agreement in interpolation */ 151 152 PetscReal gamma1; /* factor used for interpolation */ 153 PetscReal gamma2; /* factor used for interpolation */ 154 PetscReal gamma3; /* factor used for interpolation */ 155 PetscReal gamma4; /* factor used for interpolation */ 156 157 PetscReal theta; /* factor used for interpolation */ 158 159 /* Parameters when initializing trust-region radius based on interpolation */ 160 PetscReal mu1_i; /* used for model agreement in interpolation */ 161 PetscReal mu2_i; /* used for model agreement in interpolation */ 162 163 PetscReal gamma1_i; /* factor used for interpolation */ 164 PetscReal gamma2_i; /* factor used for interpolation */ 165 PetscReal gamma3_i; /* factor used for interpolation */ 166 PetscReal gamma4_i; /* factor used for interpolation */ 167 168 PetscReal theta_i; /* factor used for interpolation */ 169 170 /* Other parameters */ 171 PetscReal min_radius; /* lower bound on initial radius value */ 172 PetscReal max_radius; /* upper bound on trust region radius */ 173 PetscReal epsilon; /* tolerance used when computing ared/pred */ 174 PetscReal dmin, dmax; /* upper and lower bounds for the Hessian diagonal vector */ 175 176 PetscInt newt; /* Newton directions attempted */ 177 PetscInt bfgs; /* BFGS directions attempted */ 178 PetscInt sgrad; /* Scaled gradient directions attempted */ 179 PetscInt grad; /* Gradient directions attempted */ 180 181 PetscInt as_type; /* Active set estimation method */ 182 PetscInt pc_type; /* Preconditioner for the code */ 183 PetscInt bfgs_scale_type; /* Scaling matrix to used for the bfgs preconditioner */ 184 PetscInt init_type; /* Trust-region initialization method */ 185 PetscInt update_type; /* Trust-region update method */ 186 187 /* Trackers for KSP solution type and convergence reasons */ 188 PetscInt ksp_atol; 189 PetscInt ksp_rtol; 190 PetscInt ksp_ctol; 191 PetscInt ksp_negc; 192 PetscInt ksp_dtol; 193 PetscInt ksp_iter; 194 PetscInt ksp_othr; 195 PetscBool is_nash, is_stcg, is_gltr; 196 } TAO_BNK; 197 198 #endif /* if !defined(__TAO_BNK_H) */ 199 200 #define BNK_NEWTON 0 201 #define BNK_BFGS 1 202 #define BNK_SCALED_GRADIENT 2 203 #define BNK_GRADIENT 3 204 205 #define BNK_PC_NONE 0 206 #define BNK_PC_AHESS 1 207 #define BNK_PC_BFGS 2 208 #define BNK_PC_PETSC 3 209 #define BNK_PC_TYPES 4 210 211 #define BFGS_SCALE_AHESS 0 212 #define BFGS_SCALE_PHESS 1 213 #define BFGS_SCALE_BFGS 2 214 #define BFGS_SCALE_TYPES 3 215 216 #define BNK_INIT_CONSTANT 0 217 #define BNK_INIT_DIRECTION 1 218 #define BNK_INIT_INTERPOLATION 2 219 #define BNK_INIT_TYPES 3 220 221 #define BNK_UPDATE_STEP 0 222 #define BNK_UPDATE_REDUCTION 1 223 #define BNK_UPDATE_INTERPOLATION 2 224 #define BNK_UPDATE_TYPES 3 225 226 #define BNK_AS_NONE 0 227 #define BNK_AS_BERTSEKAS 1 228 #define BNK_AS_TYPES 2 229 230 PETSC_INTERN PetscErrorCode TaoCreate_BNK(Tao); 231 PETSC_INTERN PetscErrorCode TaoSetUp_BNK(Tao); 232 233 PETSC_INTERN PetscErrorCode TaoBNKPreconBFGS(PC, Vec, Vec); 234 PETSC_INTERN PetscErrorCode TaoBNKInitialize(Tao, PetscInt, PetscBool*); 235 PETSC_INTERN PetscErrorCode TaoBNKEstimateActiveSet(Tao, PetscInt); 236 PETSC_INTERN PetscErrorCode TaoBNKComputeHessian(Tao); 237 PETSC_INTERN PetscErrorCode TaoBNKBoundStep(Tao, PetscInt, Vec); 238 PETSC_INTERN PetscErrorCode TaoBNKTakeCGSteps(Tao, PetscBool*); 239 PETSC_INTERN PetscErrorCode TaoBNKComputeStep(Tao, PetscBool, KSPConvergedReason*); 240 PETSC_INTERN PetscErrorCode TaoBNKRecomputePred(Tao, Vec, PetscReal*); 241 PETSC_INTERN PetscErrorCode TaoBNKSafeguardStep(Tao, KSPConvergedReason, PetscInt*); 242 PETSC_INTERN PetscErrorCode TaoBNKPerformLineSearch(Tao, PetscInt*, PetscReal*, TaoLineSearchConvergedReason*); 243 PETSC_INTERN PetscErrorCode TaoBNKUpdateTrustRadius(Tao, PetscReal, PetscReal, PetscInt, PetscInt, PetscBool*); 244 PETSC_INTERN PetscErrorCode TaoBNKAddStepCounts(Tao, PetscInt); 245