C++ LinearFunctionPtr::eval方法代码示例

void ParQGHandler::addCut_(const double *nlpx, const double *lpx,
                        ConstraintPtr con, CutManager *cutman,
                        SeparationStatus *status)
{
  int error=0;
  ConstraintPtr newcon;
  std::stringstream sstm;
  double c, lpvio, act, cUb;
  FunctionPtr f = con->getFunction();
  LinearFunctionPtr lf = LinearFunctionPtr();

  act = con->getActivity(nlpx, &error);
  if (error == 0) {
    linearAt_(f, act, nlpx, &c, &lf, &error);
    if (error==0) {
      cUb = con->getUb();
      lpvio = std::max(lf->eval(lpx)-cUb+c, 0.0);
      if ((lpvio>solAbsTol_) &&
          ((cUb-c)==0 || (lpvio>fabs(cUb-c)*solRelTol_))) {
#if SPEW
        logger_->msgStream(LogDebug) << me_ << " linearization of constraint "
          << con->getName() << " violated at LP solution with violation = " <<
          lpvio << ". OA cut added." << std::endl;
#endif
        ++(stats_->cuts);
        sstm << "_OAcut_";
        sstm << stats_->cuts;
        *status = SepaResolve;
        f = (FunctionPtr) new Function(lf);
        newcon = rel_->newConstraint(f, -INFINITY, cUb-c, sstm.str());
        CutPtr cut = (CutPtr) new Cut(minlp_->getNumVars(),f, -INFINITY,
                                      cUb-c, false,false);
        cut->setCons(newcon);
        cutman->addCutToPool(cut);
        return;
      }
    }
  }	else {
    logger_->msgStream(LogError) << me_ << "Constraint not defined at"
      << " this point. "<<  std::endl;
#if SPEW
          logger_->msgStream(LogDebug) << me_ << "constraint " <<
            con->getName() << " not defined at this point." << std::endl;
#endif
  }
  return;
}

void ParQGHandler::oaCutEngLim_(const double *lpx, CutManager *,
                             SeparationStatus *status)
{
  int error=0;
  FunctionPtr f;
  LinearFunctionPtr lf;
  std::stringstream sstm;
  ConstraintPtr con, newcon;
  double c, lpvio, nlpact, cUb;

  for (CCIter it = nlCons_.begin(); it != nlCons_.end(); ++it) {
    con = *it;
    f = con->getFunction();
    lf = LinearFunctionPtr();
    nlpact =  con->getActivity(lpx, &error);
    if (error == 0) {
      cUb = con->getUb();
      if ((nlpact > cUb + solAbsTol_) &&
          (cUb == 0 || nlpact > cUb+fabs(cUb)*solRelTol_)) {
        linearAt_(f, nlpact, lpx, &c, &lf, &error);
        if (error==0) {
          lpvio = std::max(lf->eval(lpx)-cUb+c, 0.0);
          if ((lpvio>solAbsTol_) && ((cUb-c)==0 ||
                                   (lpvio>fabs(cUb-c)*solRelTol_))) {
            ++(stats_->cuts);
            sstm << "_OAcut_";
            sstm << stats_->cuts;
            *status = SepaResolve;
            f = (FunctionPtr) new Function(lf);
            newcon = rel_->newConstraint(f, -INFINITY, cUb-c, sstm.str());
            return;
          }
        }
      }
    }	else {
      logger_->msgStream(LogError) << me_ << " constraint not defined at" <<
        " this point. "<<  std::endl;
    }
  }
  return;
}

void ParQGHandler::oaCutToObj_(const double *nlpx, const double *lpx,
                            CutManager *, SeparationStatus *status)
{
  if (oNl_) {
    int error=0;
    FunctionPtr f;
    double c, vio, act;
    ConstraintPtr newcon;
    std::stringstream sstm;
    ObjectivePtr o = minlp_->getObjective();

    act = o->eval(lpx, &error);
    if (error == 0) {
      vio = std::max(act-relobj_, 0.0);
      if ((vio > solAbsTol_) &&
          (relobj_ == 0 || vio > fabs(relobj_)*solRelTol_)) {
#if SPEW
        logger_->msgStream(LogDebug) << me_ << " objective violated at LP "
          << " solution with violation = " << vio << std::endl;
#endif
        act = o->eval(nlpx, &error);
        if (error == 0) {
          f = o->getFunction();
          LinearFunctionPtr lf = LinearFunctionPtr();
          linearAt_(f, act, nlpx, &c, &lf, &error);
          if (error == 0) {
            vio = std::max(c+lf->eval(lpx)-relobj_, 0.0);
            if ((vio > solAbsTol_) && ((relobj_-c)==0
                                     || vio > fabs(relobj_-c)*solRelTol_)) {
#if SPEW
              logger_->msgStream(LogDebug) << me_ << "linearization of "
                "objective violated at LP solution with violation = " <<
                vio << ". OA cut added." << std::endl;
#endif
              ++(stats_->cuts);
              sstm << "_OAObjcut_";
              sstm << stats_->cuts;
              lf->addTerm(objVar_, -1.0);
              *status = SepaResolve;
              f = (FunctionPtr) new Function(lf);
              newcon = rel_->newConstraint(f, -INFINITY, -1.0*c, sstm.str());
            }
          }
        }
      }  else {
#if SPEW
        logger_->msgStream(LogDebug) << me_ << " objective feasible at LP "
          << " solution. No OA cut to be added." << std::endl;
#endif
      }
    }	else {
      logger_->msgStream(LogError) << me_
        << " objective not defined at this solution point." << std::endl;
#if SPEW
      logger_->msgStream(LogDebug) << me_ << " objective not defined at this "
        << " point." << std::endl;
#endif
    }
  }
  return;
}

//.........这里部分代码省略.........
    vector<vector<double> > V;
    allExtreme(nlvars, lb, ub, V);
    for(UInt j = 0; j < V.size(); j++) {
      VariablePtr vnew = VariablePtr(new Variable(nv,0.0,1.0,Continuous));
      variables_.push_back(vnew);
      lambdavars[i].push_back(vnew);
      nv++;
    }

    for(UInt k = 0; k < nlvars.size(); k++) {
      // Create 'x' portion of this constraint
      LinearFunctionPtr lf = LinearFunctionPtr(new LinearFunction());
      lf->addTerm(xvars[nlvars[k]], -1.0);

      for(UInt j = 0; j < V.size(); j++) {
        lf->addTerm(lambdavars[i][j], V[j][k]);
      }                
      
      FunctionPtr f = (FunctionPtr) new Function(lf);
      ConstraintPtr c = (ConstraintPtr) new Constraint(f, 0.0, 0.0);
      constraints_.push_back(c);  

#if defined(DEBUG_ROW_BY_ROW)
      c->display();
#endif
    }
          
    
    // Now add the 'z' var definition
    LinearFunctionPtr tlf = LinearFunctionPtr(new LinearFunction());
    tlf->addTerm(zvar, -1.0);
    const FunctionPtr mlf = mlc->getFunction();
    const LinearFunctionPtr mlf_lf = mlf->getLinearFunction();
    

    for(int j = 0; j < V.size(); j++) {
      //XXX Kludgy, but just create a vector for the (full) point
      int onv = originalInstance_->getNumVars();
      vector<double> xe(onv,0.0);
      for(UInt k = 0; k < nlvars.size(); k++) {
        xe[nlvars[k]] = V[j][k];
      }
      double zval = mlf->eval(xe);
      // Need to subtract off linear part, since yuo keep those variables in
      zval -= mlf_lf->eval(xe);

      tlf->addTerm(lambdavars[i][j], zval);
    }
    FunctionPtr tf = (FunctionPtr) new Function(tlf);
    ConstraintPtr tc = (ConstraintPtr) new Constraint(tf, 0.0, 0.0);
    constraints_.push_back(tc);  

#if defined(DEBUG_ROW_BY_ROW)
      tc->display();
#endif
    
    // Now add the linear constraint involving linear x and t
    LinearFunctionPtr xtlf = LinearFunctionPtr(new LinearFunction());
    for(ConstVariableGroupIterator it = mlf_lf->varsBegin(); it != mlf_lf->varsEnd(); ++it) {
      xtlf->addTerm(xvars[it->first->getId()], it->second);
    }
    xtlf->addTerm(zvar,1.0);
    FunctionPtr xtf = (FunctionPtr) new Function(xtlf);
    ConstraintPtr xtc = (ConstraintPtr) new Constraint(xtf, mlc->getLb(), mlc->getUb());
    constraints_.push_back(xtc);  

#if defined(DEBUG_ROW_BY_ROW)
    xtc->display();
#endif


    // Also add sum (lambda) = 1
    LinearFunctionPtr convex_lf =  LinearFunctionPtr(new LinearFunction());
    for(int j = 0; j < V.size(); j++) {
      convex_lf->addTerm(lambdavars[i][j], 1.0);
    }
    FunctionPtr convex_f = (FunctionPtr) new Function(convex_lf);
    ConstraintPtr convex_c = (ConstraintPtr) new Constraint(convex_f, 1.0, 1.0);
    constraints_.push_back(convex_c);      

#if defined(DEBUG_ROW_BY_ROW)
    convex_c->display();
#endif

    
  }



  LinearFunctionPtr olf = LinearFunctionPtr(new LinearFunction());
  const LinearFunctionPtr originalInstance_olf = originalInstance_->getObjective()->getLinearFunction();
  for (ConstVariableGroupIterator it = originalInstance_olf->varsBegin(); 
       it != originalInstance_olf->varsEnd(); ++it) {
    olf->addTerm(xvars[it->first->getId()], it->second);    
  }
  FunctionPtr of = (FunctionPtr) new Function(olf);
  objective_ = ObjectivePtr(new Objective(of, 0));    
  

}