Actual source code: ntr.h

  1: /*$Id$*/

  3: // Context for a Newton trust region method (unconstrained minimization)

  5: #ifndef __TAO_NTR_H
  7: #include "tao_solver.h"
  8: #include "src/matrix/lmvmmat.h"

 10: typedef struct {
 11:   TaoLMVMMat *M;

 13:   TaoVec *G;
 14:   TaoVec *D;
 15:   TaoVec *W;

 17:   TaoVec *Diag;

 19:   // Parameters when updating the trust-region radius based on reduction
 20:   double eta1;          // used to compute trust-region radius
 21:   double eta2;          // used to compute trust-region radius
 22:   double eta3;          // used to compute trust-region radius
 23:   double eta4;          // used to compute trust-region radius

 25:   double alpha1;        // factor used for trust-region update
 26:   double alpha2;        // factor used for trust-region update
 27:   double alpha3;        // factor used for trust-region update
 28:   double alpha4;        // factor used for trust-region update
 29:   double alpha5;        // factor used for trust-region update

 31:   // kappa = ared / pred
 32:   // if   kappa < eta1          (very bad step)
 33:   //   radius = alpha1 * min(norm(d), radius)
 34:   // elif kappa < eta2          (bad step)
 35:   //   radius = alpha2 * min(norm(d), radius)
 36:   // elif kappa < eta3          (okay step)
 37:   //   radius = alpha3 * radius;
 38:   // elif kappa < eta4          (good step)
 39:   //   radius = max(alpha4 * norm(d), radius)
 40:   // else                       (very good step)
 41:   //   radius = max(alpha5 * norm(d), radius)
 42:   // fi

 44:   // Parameters when updating the trust-region radius based on interpolation
 45:   double mu1;                // used for model agreement in radius update
 46:   double mu2;                // used for model agreement in radius update

 48:   double gamma1;        // factor used for radius update
 49:   double gamma2;        // factor used for radius update
 50:   double gamma3;        // factor used for radius update
 51:   double gamma4;        // factor used for radius update

 53:   double theta;                // factor used for radius update

 55:   // kappa = ared / pred
 56:   // if   kappa >= 1.0 - mu1    (very good step)
 57:   //   choose tau in [gamma3, gamma4]
 58:   //   radius = max(tau * norm(d), radius)
 59:   // elif kappa >= 1.0 - mu2    (good step)
 60:   //   choose tau in [gamma2, gamma3]
 61:   //   if (tau >= 1.0)
 62:   //     radius = max(tau * norm(d), radius)
 63:   //   else
 64:   //     radius = tau * min(norm(d), radius)
 65:   //   fi
 66:   // else                       (bad step)
 67:   //   choose tau in [gamma1, 1.0]
 68:   //   radius = tau * min(norm(d), radius)
 69:   // fi

 71:   // Parameters when initializing trust-region radius based on interpolation
 72:   double mu1_i;         // used for model agreement in interpolation
 73:   double mu2_i;         // used for model agreement in interpolation

 75:   double gamma1_i;      // factor used for interpolation
 76:   double gamma2_i;      // factor used for interpolation
 77:   double gamma3_i;      // factor used for interpolation
 78:   double gamma4_i;      // factor used for interpolation

 80:   double theta_i;       // factor used for interpolation

 82:   double min_radius;        // lower bound on initial radius value
 83:   double max_radius;        // upper bound on trust region radius
 84:   double epsilon;        // tolerance used when computing actred/prered

 86:   int ksp_type;         // KSP method for the code
 87:   int pc_type;          // Preconditioner for the code
 88:   int bfgs_scale_type;  // Scaling matrix for the bfgs preconditioner
 89:   int init_type;        // Trust-region initialization method
 90:   int update_type;        // Trust-region update method
 91: } TAO_NTR;

 93: #endif