1*1daac38eSTodd Munson /* 2*1daac38eSTodd Munson Context for a Newton line search method (unconstrained minimization) 3*1daac38eSTodd Munson */ 4*1daac38eSTodd Munson 5*1daac38eSTodd Munson #if !defined(__TAO_NLS_H) 6*1daac38eSTodd Munson #define __TAO_NLS_H 7*1daac38eSTodd Munson #include <petsc/private/taoimpl.h> 8*1daac38eSTodd Munson 9*1daac38eSTodd Munson typedef struct { 10*1daac38eSTodd Munson Mat M; 11*1daac38eSTodd Munson 12*1daac38eSTodd Munson Vec D; 13*1daac38eSTodd Munson Vec W; 14*1daac38eSTodd Munson 15*1daac38eSTodd Munson Vec Xold; 16*1daac38eSTodd Munson Vec Gold; 17*1daac38eSTodd Munson Vec Diag; 18*1daac38eSTodd Munson 19*1daac38eSTodd Munson /* Parameters when updating the perturbation added to the Hessian matrix 20*1daac38eSTodd Munson according to the following scheme: 21*1daac38eSTodd Munson 22*1daac38eSTodd Munson pert = sval; 23*1daac38eSTodd Munson 24*1daac38eSTodd Munson do until convergence 25*1daac38eSTodd Munson shift Hessian by pert 26*1daac38eSTodd Munson solve Newton system 27*1daac38eSTodd Munson 28*1daac38eSTodd Munson if (linear solver failed or did not compute a descent direction) 29*1daac38eSTodd Munson use steepest descent direction and increase perturbation 30*1daac38eSTodd Munson 31*1daac38eSTodd Munson if (0 == pert) 32*1daac38eSTodd Munson initialize perturbation 33*1daac38eSTodd Munson pert = min(imax, max(imin, imfac * norm(G))) 34*1daac38eSTodd Munson else 35*1daac38eSTodd Munson increase perturbation 36*1daac38eSTodd Munson pert = min(pmax, max(pgfac * pert, pmgfac * norm(G))) 37*1daac38eSTodd Munson fi 38*1daac38eSTodd Munson else 39*1daac38eSTodd Munson use linear solver direction and decrease perturbation 40*1daac38eSTodd Munson 41*1daac38eSTodd Munson pert = min(psfac * pert, pmsfac * norm(G)) 42*1daac38eSTodd Munson if (pert < pmin) 43*1daac38eSTodd Munson pert = 0 44*1daac38eSTodd Munson fi 45*1daac38eSTodd Munson fi 46*1daac38eSTodd Munson 47*1daac38eSTodd Munson perform line search 48*1daac38eSTodd Munson function and gradient evaluation 49*1daac38eSTodd Munson check convergence 50*1daac38eSTodd Munson od 51*1daac38eSTodd Munson */ 52*1daac38eSTodd Munson PetscReal sval; /* Starting perturbation value, default zero */ 53*1daac38eSTodd Munson 54*1daac38eSTodd Munson PetscReal imin; /* Minimum perturbation added during initialization */ 55*1daac38eSTodd Munson PetscReal imax; /* Maximum perturbation added during initialization */ 56*1daac38eSTodd Munson PetscReal imfac; /* Merit function factor during initialization */ 57*1daac38eSTodd Munson 58*1daac38eSTodd Munson PetscReal pmin; /* Minimim perturbation value */ 59*1daac38eSTodd Munson PetscReal pmax; /* Maximum perturbation value */ 60*1daac38eSTodd Munson PetscReal pgfac; /* Perturbation growth factor */ 61*1daac38eSTodd Munson PetscReal psfac; /* Perturbation shrink factor */ 62*1daac38eSTodd Munson PetscReal pmgfac; /* Merit function growth factor */ 63*1daac38eSTodd Munson PetscReal pmsfac; /* Merit function shrink factor */ 64*1daac38eSTodd Munson 65*1daac38eSTodd Munson /* Parameters when updating the trust-region radius based on steplength 66*1daac38eSTodd Munson if step < nu1 (very bad step) 67*1daac38eSTodd Munson radius = omega1 * min(norm(d), radius) 68*1daac38eSTodd Munson elif step < nu2 (bad step) 69*1daac38eSTodd Munson radius = omega2 * min(norm(d), radius) 70*1daac38eSTodd Munson elif step < nu3 (okay step) 71*1daac38eSTodd Munson radius = omega3 * radius; 72*1daac38eSTodd Munson elif step < nu4 (good step) 73*1daac38eSTodd Munson radius = max(omega4 * norm(d), radius) 74*1daac38eSTodd Munson else (very good step) 75*1daac38eSTodd Munson radius = max(omega5 * norm(d), radius) 76*1daac38eSTodd Munson fi 77*1daac38eSTodd Munson */ 78*1daac38eSTodd Munson PetscReal nu1; /* used to compute trust-region radius */ 79*1daac38eSTodd Munson PetscReal nu2; /* used to compute trust-region radius */ 80*1daac38eSTodd Munson PetscReal nu3; /* used to compute trust-region radius */ 81*1daac38eSTodd Munson PetscReal nu4; /* used to compute trust-region radius */ 82*1daac38eSTodd Munson 83*1daac38eSTodd Munson PetscReal omega1; /* factor used for trust-region update */ 84*1daac38eSTodd Munson PetscReal omega2; /* factor used for trust-region update */ 85*1daac38eSTodd Munson PetscReal omega3; /* factor used for trust-region update */ 86*1daac38eSTodd Munson PetscReal omega4; /* factor used for trust-region update */ 87*1daac38eSTodd Munson PetscReal omega5; /* factor used for trust-region update */ 88*1daac38eSTodd Munson 89*1daac38eSTodd Munson /* Parameters when updating the trust-region radius based on reduction 90*1daac38eSTodd Munson 91*1daac38eSTodd Munson kappa = ared / pred 92*1daac38eSTodd Munson if kappa < eta1 (very bad step) 93*1daac38eSTodd Munson radius = alpha1 * min(norm(d), radius) 94*1daac38eSTodd Munson elif kappa < eta2 (bad step) 95*1daac38eSTodd Munson radius = alpha2 * min(norm(d), radius) 96*1daac38eSTodd Munson elif kappa < eta3 (okay step) 97*1daac38eSTodd Munson radius = alpha3 * radius; 98*1daac38eSTodd Munson elif kappa < eta4 (good step) 99*1daac38eSTodd Munson radius = max(alpha4 * norm(d), radius) 100*1daac38eSTodd Munson else (very good step) 101*1daac38eSTodd Munson radius = max(alpha5 * norm(d), radius) 102*1daac38eSTodd Munson fi 103*1daac38eSTodd Munson */ 104*1daac38eSTodd Munson PetscReal eta1; /* used to compute trust-region radius */ 105*1daac38eSTodd Munson PetscReal eta2; /* used to compute trust-region radius */ 106*1daac38eSTodd Munson PetscReal eta3; /* used to compute trust-region radius */ 107*1daac38eSTodd Munson PetscReal eta4; /* used to compute trust-region radius */ 108*1daac38eSTodd Munson 109*1daac38eSTodd Munson PetscReal alpha1; /* factor used for trust-region update */ 110*1daac38eSTodd Munson PetscReal alpha2; /* factor used for trust-region update */ 111*1daac38eSTodd Munson PetscReal alpha3; /* factor used for trust-region update */ 112*1daac38eSTodd Munson PetscReal alpha4; /* factor used for trust-region update */ 113*1daac38eSTodd Munson PetscReal alpha5; /* factor used for trust-region update */ 114*1daac38eSTodd Munson 115*1daac38eSTodd Munson /* Parameters when updating the trust-region radius based on interpolation 116*1daac38eSTodd Munson 117*1daac38eSTodd Munson kappa = ared / pred 118*1daac38eSTodd Munson if kappa >= 1.0 - mu1 (very good step) 119*1daac38eSTodd Munson choose tau in [gamma3, gamma4] 120*1daac38eSTodd Munson radius = max(tau * norm(d), radius) 121*1daac38eSTodd Munson elif kappa >= 1.0 - mu2 (good step) 122*1daac38eSTodd Munson choose tau in [gamma2, gamma3] 123*1daac38eSTodd Munson if (tau >= 1.0) 124*1daac38eSTodd Munson radius = max(tau * norm(d), radius) 125*1daac38eSTodd Munson else 126*1daac38eSTodd Munson radius = tau * min(norm(d), radius) 127*1daac38eSTodd Munson fi 128*1daac38eSTodd Munson else (bad step) 129*1daac38eSTodd Munson choose tau in [gamma1, 1.0] 130*1daac38eSTodd Munson radius = tau * min(norm(d), radius) 131*1daac38eSTodd Munson fi 132*1daac38eSTodd Munson */ 133*1daac38eSTodd Munson PetscReal mu1; /* used for model agreement in interpolation */ 134*1daac38eSTodd Munson PetscReal mu2; /* used for model agreement in interpolation */ 135*1daac38eSTodd Munson 136*1daac38eSTodd Munson PetscReal gamma1; /* factor used for interpolation */ 137*1daac38eSTodd Munson PetscReal gamma2; /* factor used for interpolation */ 138*1daac38eSTodd Munson PetscReal gamma3; /* factor used for interpolation */ 139*1daac38eSTodd Munson PetscReal gamma4; /* factor used for interpolation */ 140*1daac38eSTodd Munson 141*1daac38eSTodd Munson PetscReal theta; /* factor used for interpolation */ 142*1daac38eSTodd Munson 143*1daac38eSTodd Munson /* Parameters when initializing trust-region radius based on interpolation */ 144*1daac38eSTodd Munson PetscReal mu1_i; /* used for model agreement in interpolation */ 145*1daac38eSTodd Munson PetscReal mu2_i; /* used for model agreement in interpolation */ 146*1daac38eSTodd Munson 147*1daac38eSTodd Munson PetscReal gamma1_i; /* factor used for interpolation */ 148*1daac38eSTodd Munson PetscReal gamma2_i; /* factor used for interpolation */ 149*1daac38eSTodd Munson PetscReal gamma3_i; /* factor used for interpolation */ 150*1daac38eSTodd Munson PetscReal gamma4_i; /* factor used for interpolation */ 151*1daac38eSTodd Munson 152*1daac38eSTodd Munson PetscReal theta_i; /* factor used for interpolation */ 153*1daac38eSTodd Munson 154*1daac38eSTodd Munson /* Other parameters */ 155*1daac38eSTodd Munson PetscReal min_radius; /* lower bound on initial radius value */ 156*1daac38eSTodd Munson PetscReal max_radius; /* upper bound on trust region radius */ 157*1daac38eSTodd Munson PetscReal epsilon; /* tolerance used when computing ared/pred */ 158*1daac38eSTodd Munson 159*1daac38eSTodd Munson PetscInt newt; /* Newton directions attempted */ 160*1daac38eSTodd Munson PetscInt bfgs; /* BFGS directions attempted */ 161*1daac38eSTodd Munson PetscInt sgrad; /* Scaled gradient directions attempted */ 162*1daac38eSTodd Munson PetscInt grad; /* Gradient directions attempted */ 163*1daac38eSTodd Munson 164*1daac38eSTodd Munson 165*1daac38eSTodd Munson PetscInt pc_type; /* Preconditioner for the code */ 166*1daac38eSTodd Munson PetscInt bfgs_scale_type; /* Scaling matrix to used for the bfgs preconditioner */ 167*1daac38eSTodd Munson PetscInt init_type; /* Trust-region initialization method */ 168*1daac38eSTodd Munson PetscInt update_type; /* Trust-region update method */ 169*1daac38eSTodd Munson 170*1daac38eSTodd Munson PetscInt ksp_atol; 171*1daac38eSTodd Munson PetscInt ksp_rtol; 172*1daac38eSTodd Munson PetscInt ksp_ctol; 173*1daac38eSTodd Munson PetscInt ksp_negc; 174*1daac38eSTodd Munson PetscInt ksp_dtol; 175*1daac38eSTodd Munson PetscInt ksp_iter; 176*1daac38eSTodd Munson PetscInt ksp_othr; 177*1daac38eSTodd Munson } TAO_NLS; 178*1daac38eSTodd Munson 179*1daac38eSTodd Munson #endif /* if !defined(__TAO_NLS_H) */ 180