1 /* 2 Context for a Newton trust-region, line-search method for unconstrained 3 minimization 4 */ 5 6 #pragma once 7 #include <petsc/private/taoimpl.h> 8 9 typedef struct { 10 Mat M; 11 PC bfgs_pre; 12 13 Vec W; 14 Vec Xold; 15 Vec Gold; 16 17 /* Parameters when updating the trust-region radius based on steplength 18 19 if step < nu1 (very bad step) 20 radius = omega1 * min(norm(d), radius) 21 elif step < nu2 (bad step) 22 radius = omega2 * min(norm(d), radius) 23 elif step < nu3 (okay step) 24 radius = omega3 * radius; 25 elif step < nu4 (good step) 26 radius = max(omega4 * norm(d), radius) 27 else (very good step) 28 radius = max(omega5 * norm(d), radius) 29 fi 30 */ 31 32 PetscReal nu1; /* used to compute trust-region radius */ 33 PetscReal nu2; /* used to compute trust-region radius */ 34 PetscReal nu3; /* used to compute trust-region radius */ 35 PetscReal nu4; /* used to compute trust-region radius */ 36 37 PetscReal omega1; /* factor used for trust-region update */ 38 PetscReal omega2; /* factor used for trust-region update */ 39 PetscReal omega3; /* factor used for trust-region update */ 40 PetscReal omega4; /* factor used for trust-region update */ 41 PetscReal omega5; /* factor used for trust-region update */ 42 43 /* Parameters when updating the trust-region radius based on reduction 44 45 kappa = ared / pred 46 if kappa < eta1 (very bad step) 47 radius = alpha1 * min(norm(d), radius) 48 elif kappa < eta2 (bad step) 49 radius = alpha2 * min(norm(d), radius) 50 elif kappa < eta3 (okay step) 51 radius = alpha3 * radius; 52 elif kappa < eta4 (good step) 53 radius = max(alpha4 * norm(d), radius) 54 else (very good step) 55 radius = max(alpha5 * norm(d), radius) 56 fi 57 */ 58 59 PetscReal eta1; /* used to compute trust-region radius */ 60 PetscReal eta2; /* used to compute trust-region radius */ 61 PetscReal eta3; /* used to compute trust-region radius */ 62 PetscReal eta4; /* used to compute trust-region radius */ 63 64 PetscReal alpha1; /* factor used for trust-region update */ 65 PetscReal alpha2; /* factor used for trust-region update */ 66 PetscReal alpha3; /* factor used for trust-region update */ 67 PetscReal alpha4; /* factor used for trust-region update */ 68 PetscReal alpha5; /* factor used for trust-region update */ 69 70 /* Parameters when updating the trust-region radius based on interpolation 71 kappa = ared / pred 72 if kappa >= 1.0 - mu1 (very good step) 73 choose tau in [gamma3, gamma4] 74 radius = max(tau * norm(d), radius) 75 elif kappa >= 1.0 - mu2 (good step) 76 choose tau in [gamma2, gamma3] 77 if (tau >= 1.0) 78 radius = max(tau * norm(d), radius) 79 else 80 radius = tau * min(norm(d), radius) 81 fi 82 else (bad step) 83 choose tau in [gamma1, 1.0] 84 radius = tau * min(norm(d), radius) 85 fi 86 */ 87 88 PetscReal mu1; /* used for model agreement in interpolation */ 89 PetscReal mu2; /* used for model agreement in interpolation */ 90 91 PetscReal gamma1; /* factor used for interpolation */ 92 PetscReal gamma2; /* factor used for interpolation */ 93 PetscReal gamma3; /* factor used for interpolation */ 94 PetscReal gamma4; /* factor used for interpolation */ 95 96 PetscReal theta; /* factor used for interpolation */ 97 98 /* Parameters when initializing trust-region radius based on interpolation */ 99 PetscReal mu1_i; /* used for model agreement in interpolation */ 100 PetscReal mu2_i; /* used for model agreement in interpolation */ 101 102 PetscReal gamma1_i; /* factor used for interpolation */ 103 PetscReal gamma2_i; /* factor used for interpolation */ 104 PetscReal gamma3_i; /* factor used for interpolation */ 105 PetscReal gamma4_i; /* factor used for interpolation */ 106 107 PetscReal theta_i; /* factor used for interpolation */ 108 109 /* Other parameters */ 110 PetscReal min_radius; /* lower bound on initial radius value */ 111 PetscReal max_radius; /* upper bound on trust region radius */ 112 PetscReal epsilon; /* tolerance used when computing ared/pred */ 113 114 PetscInt ntrust; /* Trust-region steps accepted */ 115 PetscInt newt; /* Newton directions attempted */ 116 PetscInt bfgs; /* BFGS directions attempted */ 117 PetscInt grad; /* Gradient directions attempted */ 118 119 PetscInt init_type; /* Trust-region initialization method */ 120 PetscInt update_type; /* Trust-region update method */ 121 } TAO_NTL; 122