C++ SmartPtr::Nrm2方法代码示例

  bool CGPenaltyLSAcceptor::IsAcceptableToPiecewisePenalty(Number alpha_primal_test)
  {
    DBG_START_METH("CGPenaltyLSAcceptor::IsAcceptableToPiecewisePenalty",
                   dbg_verbosity);

    // If the current infeasibility is small, we require the barrier to be decreased.
    bool accept = false;
    Number infeasibility = IpCq().curr_primal_infeasibility(NORM_MAX);
    SmartPtr<const Vector> dx = IpData().delta()->x();
    SmartPtr<const Vector> ds = IpData().delta()->s();
    Number curr_barr = IpCq().curr_barrier_obj();
    Number trial_barr = IpCq().trial_barrier_obj();
    Number nrm_dx_ds = pow(dx->Nrm2(),2.) + pow(ds->Nrm2(),2.);
    if (infeasibility < theta_min_) {
      Number biggest_barr = PiecewisePenalty_.BiggestBarr();
      accept = Compare_le(trial_barr-biggest_barr, -alpha_primal_test*
                          piecewisepenalty_gamma_obj_*nrm_dx_ds, curr_barr);
      if (!accept) {
        return false;
      }
    }
    Number Fzconst, Fzlin;
    Fzconst = IpCq().trial_barrier_obj() + alpha_primal_test *
              piecewisepenalty_gamma_obj_ * nrm_dx_ds;
    Fzlin = IpCq().trial_constraint_violation() + alpha_primal_test *
            piecewisepenalty_gamma_infeasi_ * IpCq().curr_constraint_violation();
    accept  = PiecewisePenalty_.Acceptable(Fzconst, Fzlin);
    return accept;
  }

Number InexactLSAcceptor::ComputeAlphaForY(
   Number                    alpha_primal,
   Number                    alpha_dual,
   SmartPtr<IteratesVector>& delta
   )
{
   DBG_START_METH("InexactLSAcceptor::ComputeAlphaForY",
      dbg_verbosity);

   // Here, we choose as stepsize for y either alpha_primal, if the
   // conditions from the ineqaxt paper is satisfied for it, or we
   // compute the step size closest to alpha_primal but great than
   // it, that does give the same progress as the full step would
   // give.

   Number alpha_y = alpha_primal;

   SmartPtr<Vector> gx = IpCq().curr_grad_barrier_obj_x()->MakeNewCopy();
   gx->AddTwoVectors(1., *IpCq().curr_jac_cT_times_curr_y_c(), 1., *IpCq().curr_jac_dT_times_curr_y_d(), 1.);
   SmartPtr<Vector> Jxy = gx->MakeNew();
   IpCq().curr_jac_c()->TransMultVector(1., *delta->y_c(), 0., *Jxy);
   IpCq().curr_jac_d()->TransMultVector(1., *delta->y_d(), 1., *Jxy);

   SmartPtr<const Vector> curr_scaling_slacks = InexCq().curr_scaling_slacks();
   SmartPtr<Vector> gs = curr_scaling_slacks->MakeNew();
   gs->AddTwoVectors(1., *IpCq().curr_grad_barrier_obj_s(), -1., *IpData().curr()->y_d(), 0.);
   gs->ElementWiseMultiply(*curr_scaling_slacks);

   SmartPtr<Vector> Sdy = delta->y_d()->MakeNewCopy();
   Sdy->ElementWiseMultiply(*curr_scaling_slacks);

   // using the magic formula in my notebook
   Number a = pow(Jxy->Nrm2(), 2) + pow(Sdy->Nrm2(), 2);
   Number b = 2 * (gx->Dot(*Jxy) - gs->Dot(*Sdy));
   Number c = pow(gx->Nrm2(), 2) + pow(gs->Nrm2(), 2);

   // First we check if the primal step size is good enough:
   Number val_ap = alpha_primal * alpha_primal * a + alpha_primal * b + c;
   Number val_1 = a + b + c;

   if( val_ap <= val_1 )
   {
      Jnlst().Printf(J_DETAILED, J_LINE_SEARCH, "  Step size for y: using alpha_primal\n.");
   }
   else
   {
      Number alpha_2 = -b / a - 1.;
      Jnlst().Printf(J_DETAILED, J_LINE_SEARCH, "  Step size for y candidate: %8.2e - ", alpha_2);
      if( alpha_2 > alpha_primal && alpha_2 < 1. )
      {
         alpha_y = alpha_2;
         Jnlst().Printf(J_DETAILED, J_LINE_SEARCH, "using that one\n.");
      }
      else
      {
         alpha_y = 1.;
         Jnlst().Printf(J_DETAILED, J_LINE_SEARCH, "using 1 instead\n");
      }
   }

   return alpha_y;
}