C++ OsiSolverInterface::getObjSense方法代码示例

// Create result
void OsiSolverResult::createResult(const OsiSolverInterface &solver, const double *lowerBefore,
  const double *upperBefore)
{
  delete[] primalSolution_;
  delete[] dualSolution_;
  if (solver.isProvenOptimal() && !solver.isDualObjectiveLimitReached()) {
    objectiveValue_ = solver.getObjValue() * solver.getObjSense();
    CoinWarmStartBasis *basis = dynamic_cast< CoinWarmStartBasis * >(solver.getWarmStart());
    assert(basis);
    basis_ = *basis;
    int numberRows = basis_.getNumArtificial();
    int numberColumns = basis_.getNumStructural();
    assert(numberColumns == solver.getNumCols());
    assert(numberRows == solver.getNumRows());
    primalSolution_ = CoinCopyOfArray(solver.getColSolution(), numberColumns);
    dualSolution_ = CoinCopyOfArray(solver.getRowPrice(), numberRows);
    fixed_.addBranch(-1, numberColumns, lowerBefore, solver.getColLower(),
      upperBefore, solver.getColUpper());
  } else {
    // infeasible
    objectiveValue_ = COIN_DBL_MAX;
    basis_ = CoinWarmStartBasis();
    ;
    primalSolution_ = NULL;
    dualSolution_ = NULL;
  }
}

BcpsBranchObject * 
BlisObjectInt::preferredNewFeasible(BcpsModel *m) const
{
    BlisModel *model = dynamic_cast<BlisModel* >(m);
    OsiSolverInterface * solver = model->solver();
    
    double value = solver->getColSolution()[columnIndex_];
    
    double nearest = floor(value + 0.5);
    double integerTolerance = model->BlisPar()->entry(BlisParams::integerTol);

    assert(fabs(value - nearest) <= integerTolerance);

    double dj = solver->getObjSense()*solver->getReducedCost()[columnIndex_];

    BlisBranchObjectInt * object = NULL;

    if (dj >= 0.0) {
	// Better go down
	if (nearest > originalLower_ + 0.5) {
	    // Has room to go down
	    object = new BlisBranchObjectInt(model,
                                             objectIndex_,
                                             -1,
                                             nearest - 1.0,
                                             nearest - 1.0);
	}
    } 
    else {
	// Better go up
	if (nearest < originalUpper_ - 0.5) {
	    // Has room to go up
	    object = new BlisBranchObjectInt(model, 
                                             objectIndex_, 
                                             -1,
                                             nearest + 1.0,
                                             nearest + 1.0);
	}
    }

    return object;
}

// inner part of dive
int 
CbcHeuristicDive::solution(double & solutionValue, int & numberNodes,
			   int & numberCuts, OsiRowCut ** cuts,
			   CbcSubProblem ** & nodes,
			   double * newSolution)
{
#ifdef DIVE_DEBUG
    int nRoundInfeasible = 0;
    int nRoundFeasible = 0;
#endif
    int reasonToStop = 0;
    double time1 = CoinCpuTime();
    int numberSimplexIterations = 0;
    int maxSimplexIterations = (model_->getNodeCount()) ? maxSimplexIterations_
                               : maxSimplexIterationsAtRoot_;
    // but can't be exactly coin_int_max
    maxSimplexIterations = CoinMin(maxSimplexIterations,COIN_INT_MAX>>3);
    OsiSolverInterface * solver = cloneBut(6); // was model_->solver()->clone();
# ifdef COIN_HAS_CLP
    OsiClpSolverInterface * clpSolver
    = dynamic_cast<OsiClpSolverInterface *> (solver);
    if (clpSolver) {
      ClpSimplex * clpSimplex = clpSolver->getModelPtr();
      int oneSolveIts = clpSimplex->maximumIterations();
      oneSolveIts = CoinMin(1000+2*(clpSimplex->numberRows()+clpSimplex->numberColumns()),oneSolveIts);
      clpSimplex->setMaximumIterations(oneSolveIts);
      if (!nodes) {
        // say give up easily
        clpSimplex->setMoreSpecialOptions(clpSimplex->moreSpecialOptions() | 64);
      } else {
	// get ray
	int specialOptions = clpSimplex->specialOptions();
	specialOptions &= ~0x3100000;
	specialOptions |= 32;
        clpSimplex->setSpecialOptions(specialOptions);
        clpSolver->setSpecialOptions(clpSolver->specialOptions() | 1048576);
	if ((model_->moreSpecialOptions()&16777216)!=0) {
	  // cutoff is constraint
	  clpSolver->setDblParam(OsiDualObjectiveLimit, COIN_DBL_MAX);
	}
      }
    }
# endif
    const double * lower = solver->getColLower();
    const double * upper = solver->getColUpper();
    const double * rowLower = solver->getRowLower();
    const double * rowUpper = solver->getRowUpper();
    const double * solution = solver->getColSolution();
    const double * objective = solver->getObjCoefficients();
    double integerTolerance = model_->getDblParam(CbcModel::CbcIntegerTolerance);
    double primalTolerance;
    solver->getDblParam(OsiPrimalTolerance, primalTolerance);

    int numberRows = matrix_.getNumRows();
    assert (numberRows <= solver->getNumRows());
    int numberIntegers = model_->numberIntegers();
    const int * integerVariable = model_->integerVariable();
    double direction = solver->getObjSense(); // 1 for min, -1 for max
    double newSolutionValue = direction * solver->getObjValue();
    int returnCode = 0;
    // Column copy
    const double * element = matrix_.getElements();
    const int * row = matrix_.getIndices();
    const CoinBigIndex * columnStart = matrix_.getVectorStarts();
    const int * columnLength = matrix_.getVectorLengths();
#ifdef DIVE_FIX_BINARY_VARIABLES
    // Row copy
    const double * elementByRow = matrixByRow_.getElements();
    const int * column = matrixByRow_.getIndices();
    const CoinBigIndex * rowStart = matrixByRow_.getVectorStarts();
    const int * rowLength = matrixByRow_.getVectorLengths();
#endif

    // Get solution array for heuristic solution
    int numberColumns = solver->getNumCols();
    memcpy(newSolution, solution, numberColumns*sizeof(double));

    // vectors to store the latest variables fixed at their bounds
    int* columnFixed = new int [numberIntegers];
    double* originalBound = new double [numberIntegers+2*numberColumns];
    double * lowerBefore = originalBound+numberIntegers;
    double * upperBefore = lowerBefore+numberColumns;
    memcpy(lowerBefore,lower,numberColumns*sizeof(double));
    memcpy(upperBefore,upper,numberColumns*sizeof(double));
    double * lastDjs=newSolution+numberColumns;
    bool * fixedAtLowerBound = new bool [numberIntegers];
    PseudoReducedCost * candidate = new PseudoReducedCost [numberIntegers];
    double * random = new double [numberIntegers];

    int maxNumberAtBoundToFix = static_cast<int> (floor(percentageToFix_ * numberIntegers));
    assert (!maxNumberAtBoundToFix||!nodes);

    // count how many fractional variables
    int numberFractionalVariables = 0;
    for (int i = 0; i < numberIntegers; i++) {
        random[i] = randomNumberGenerator_.randomDouble() + 0.3;
        int iColumn = integerVariable[i];
        double value = newSolution[iColumn];
        if (fabs(floor(value + 0.5) - value) > integerTolerance) {
//.........这里部分代码省略.........

//.........这里部分代码省略.........
    else{
       dynamic_cast<OsiSymSolverInterface *> 
	  (lSolver)->setSymParam("prep_level", -1);
       
       dynamic_cast<OsiSymSolverInterface *> 
	  (lSolver)->setSymParam("verbosity", -2);
       
       dynamic_cast<OsiSymSolverInterface *> 
	  (lSolver)->setSymParam("max_active_nodes", 1);
    }

    lSolver->branchAndBound();

    if (lSolver->isProvenOptimal()){
       const double * sol = hSolver->getColSolution();
       double objVal(lSolver->getObjValue() * objSense);
       double etol(etol_);
       double lowerObj = getLowerObj(sol, objSense);  
       
       double * optUpperSolutionOrd = new double[uCols];
       double * optLowerSolutionOrd = new double[lCols];
       
       CoinZeroN(optUpperSolutionOrd, uCols);
       CoinZeroN(optLowerSolutionOrd, lCols);
       
       if(fabs(objVal - lowerObj) < etol){
	  
	  /** Current solution is bilevel feasible **/
	  
	  for(i = 0; i < tCols; i++)
	     upperObjVal += 
		hSolver->getColSolution()[i] * oSolver->getObjCoefficients()[i];
	  
	  mibSol = new MibSSolution(hSolver->getNumCols(),
				    hSolver->getColSolution(),
				    upperObjVal,
				    model);
	  
	  model->storeSolution(BlisSolutionTypeHeuristic, mibSol);
	  mibSol = NULL;
       }
       else{
	  
	  /* solution is not bilevel feasible, create one that is */
	  
	  const double * uSol = hSolver->getColSolution();
	  const double * lSol = lSolver->getColSolution();
	  int numElements(hSolver->getNumCols());
	  int i(0), pos(0), index(0);
	  double * lpSolution = new double[numElements];
	  double upperObj(0.0);
	  
	  //FIXME: problem is still here.  indices may be wrong.  
	  //also is all this necessary, or can we just paste together uSol and lSol?
	  //this may be an old comment
	  
	  for(i = 0; i < numElements; i++){
	     pos = model->bS_->binarySearch(0, lCols - 1, i, lColIndices);
	     if(pos < 0){
		pos = model->bS_->binarySearch(0, uCols - 1, i, uColIndices);
		if (pos >= 0){
		   optUpperSolutionOrd[pos] = uSol[i];
		}
	     }
	     else{
		optLowerSolutionOrd[pos] = lSol[pos];
	     }
	  }
	  
	  for(i = 0; i < uCols; i++){
	     index = uColIndices[i];
	     lpSolution[index] = optUpperSolutionOrd[i];
	     upperObj += 
		optUpperSolutionOrd[i] * oSolver->getObjCoefficients()[index];
	  }
	  
	  for(i = 0; i < lCols; i++){
	     index = lColIndices[i];
	     lpSolution[index] = optLowerSolutionOrd[i];
	     upperObj += 
		optLowerSolutionOrd[i] * oSolver->getObjCoefficients()[index];
	  }
	  
	  if(model->checkUpperFeasibility(lpSolution)){
	     mibSol = new MibSSolution(hSolver->getNumCols(),
				       lpSolution,
				       upperObj * oSolver->getObjSense(),
				       model);
	     
	     model->storeSolution(BlisSolutionTypeHeuristic, mibSol);
	     mibSol = NULL;
	  }
	  delete [] lpSolution;
       }
    }
    delete lSolver;
  }
  delete hSolver;

}

//#############################################################################
void 
MibSHeuristic::greedyHeuristic()
{

  MibSModel * model = MibSModel_;
  //OsiSolverInterface * oSolver = model->getSolver();
  OsiSolverInterface * oSolver = model->solver();
  
  double uObjSense(oSolver->getObjSense());
  double lObjSense(model->getLowerObjSense());  
  int lCols(model->getLowerDim());
  int uCols(model->getUpperDim());
  int * uColIndices = model->getUpperColInd();
  int * lColIndices = model->getLowerColInd();
  double * lObjCoeffs = model->getLowerObjCoeffs();
  double * intCost = model->getInterdictCost();
  double intBudget = model->getInterdictBudget();

  int tCols(uCols + lCols); 

  assert(tCols == oSolver->getNumCols());

  int i(0), ind_min_wt(0);
  double usedBudget(0.0); 
  double * fixedVars = new double[lCols];
  double * testsol = new double[tCols];
  CoinZeroN(fixedVars, lCols);
  CoinZeroN(testsol, tCols);

  std::multimap<double, int> lObjCoeffsOrd;

  for(i = 0; i < lCols; i++)
    lObjCoeffsOrd.insert(std::pair<double, int>(lObjCoeffs[i] * lObjSense, i));
  
  if(!bestSol_)
    bestSol_ = new double[tCols];

  //initialize the best solution information
  //bestObjVal_ = model->getSolver()->getInfinity() * uObjSense;
  //CoinZeroN(bestSol_, tCols);

  std::multimap<double, int>::iterator iter;
  //std::multimap<double, int>::iterator first;
  //std::multimap<double, int>::iterator last;
  //int dist = std::distance(first, last);
  srandom((unsigned) time(NULL));

  int randchoice(0); 
  if(0)
    std::cout << "randchoice " << randchoice << std::endl;


  double cost(0.0);

  //starting from the largest, fix corr upper-level variables
  //then, with these fixed, solve the lower-level problem
  //this yields a feasible solution

  iter = lObjCoeffsOrd.begin(); 
  
  while((usedBudget < intBudget) && (iter != lObjCoeffsOrd.end())){
    
    ind_min_wt = iter->second;
    cost = intCost[ind_min_wt];
    testsol[uColIndices[ind_min_wt]] = 1.0;
    double min_wt = iter->first;
    
    if(0){
      std::cout << "upper: " << ind_min_wt << " " 
		<< uColIndices[ind_min_wt] << " "  
		<< oSolver->getColUpper()[uColIndices[ind_min_wt]] << " " 
		<< oSolver->getColLower()[uColIndices[ind_min_wt]] << std::endl;
      
      std::cout << "lower: " << ind_min_wt << " " 
		<< lColIndices[ind_min_wt] << " "  
		<< oSolver->getColUpper()[lColIndices[ind_min_wt]] << std::endl;
    }

    //if((oSolver->getColUpper()[uColIndices[ind_min_wt]] == 1.0) 
       //&& (oSolver->getColUpper()[lColIndices[ind_min_wt]] > 0)){
    if(oSolver->getColUpper()[uColIndices[ind_min_wt]] > etol_){ 
      
      //if(((usedBudget + cost) <= intBudget) 
      // && checkLowerFeasibility(oSolver, testsol)){
      if((usedBudget + cost) <= intBudget){
	
	//FIXME: SHOULD BE CHECKING FOR CURRENT BOUNDS HERE  
	//fix the corresponding upper-level variable to 1
	randchoice = random() % 2;
	if(0)
	  std::cout << "randchoice " << random << std::endl;
	if(randchoice){
	  fixedVars[ind_min_wt] = 1.0;
	  usedBudget += intCost[ind_min_wt];
	}
      }
    }
    else{
      
//.........这里部分代码省略.........

/** Create a set of candidate branching objects. */
int 
BlisBranchStrategyPseudo::createCandBranchObjects(int numPassesLeft,
						  double ub)
{
    int bStatus = 0;
    int i, pass, colInd;

    int preferDir, saveLimit;
    int numFirsts  = 0;
    int numInfs = 0;
    int minCount = 0;
    int numLowerTightens = 0;
    int numUpperTightens = 0;
    double lpX, score, infeasibility, downDeg, upDeg, sumDeg = 0.0; 
    
    bool roundAgain, downKeep, downGood, upKeep, upGood;


    int *lbInd = NULL;
    int *ubInd = NULL;
    double *newLB = NULL;
    double *newUB = NULL;

    double *saveUpper = NULL;
    double *saveLower = NULL;
    double *saveSolution = NULL;

    BlisModel *model = dynamic_cast<BlisModel *>(model_);
    OsiSolverInterface *solver = model->solver();
    
    int numCols = model->getNumCols();
    int numObjects = model->numObjects();
    int aveIterations = model->getAveIterations();


    //std::cout <<  "aveIterations = " <<  aveIterations << std::endl;

     //------------------------------------------------------
    // Check if max time is reached or no pass is left.
    //------------------------------------------------------
    
    double timeLimit = model->AlpsPar()->entry(AlpsParams::timeLimit);
    AlpsKnowledgeBroker *broker = model->getKnowledgeBroker();
    bool maxTimeReached = (broker->timer().getTime() > timeLimit);
    bool selectNow = false;
    
    if (maxTimeReached || !numPassesLeft) {
        selectNow = true;
#ifdef BLIS_DEBUG
        printf("PSEUDO: CREATE: maxTimeReached %d, numPassesLeft %d\n", 
               maxTimeReached, numPassesLeft);
#endif
    }
    
    // Store first time objects.
    std::vector<BlisObjectInt *> firstObjects;

    // Store infeasible objects.
    std::vector<BlisObjectInt *> infObjects;

    // TODO: check if sorting is expensive.
    std::multimap<double, BcpsBranchObject*, BlisPseuoGreater> candObjects;

    double objValue = solver->getObjSense() * solver->getObjValue();

    const double * lower = solver->getColLower();
    const double * upper = solver->getColUpper();
    saveSolution = new double[numCols];
    memcpy(saveSolution, solver->getColSolution(), numCols*sizeof(double));

    //--------------------------------------------------
    // Find the infeasible objects.
    // NOTE: we might go round this loop twice if we are feed in
    //       a "feasible" solution.
    //--------------------------------------------------
    
    for (pass = 0; pass < 2; ++pass) {
	
        numInfs = 0;

        BcpsObject * object = NULL;
        BlisObjectInt * intObject = NULL;
            
        infObjects.clear();
        firstObjects.clear();
        
        for (i = 0; i < numObjects; ++i) {
                
            object = model->objects(i);
            infeasibility = object->infeasibility(model, preferDir);
            
            if (infeasibility) {
                
                ++numInfs;
                intObject = dynamic_cast<BlisObjectInt *>(object);
                
                if (intObject) {
                    infObjects.push_back(intObject);
                    
//.........这里部分代码省略.........

/** Create a set of candidate branching objects. */
int
BlisBranchStrategyRel::createCandBranchObjects(int numPassesLeft)
{
    int bStatus = 0;
    int i, pass, colInd;

    int preferDir, saveLimit;
    int numFirsts  = 0;
    int numInfs = 0;
    int minCount = 0;
    int numLowerTightens = 0;
    int numUpperTightens = 0;

    double lpX, score, infeasibility, downDeg, upDeg, sumDeg = 0.0;

    bool roundAgain, downKeep, downGood, upKeep, upGood;


    int *lbInd = NULL;
    int *ubInd = NULL;
    double *newLB = NULL;
    double *newUB = NULL;

    double * saveUpper = NULL;
    double * saveLower = NULL;
    double * saveSolution = NULL;


    BlisModel *model = dynamic_cast<BlisModel *>(model_);
    OsiSolverInterface * solver = model->solver();

    int numCols = model->getNumCols();
    int numObjects = model->numObjects();

    //int lookAhead = dynamic_cast<BlisParams*>
    //  (model->blisPar())->entry(BlisParams::lookAhead);

    //------------------------------------------------------
    // Check if max time is reached or no pass is left.
    //------------------------------------------------------

    double timeLimit = model->AlpsPar()->entry(AlpsParams::timeLimit);
    bool maxTimeReached = (CoinCpuTime() - model->startTime_  > timeLimit);
    bool selectNow = false;

    if (maxTimeReached || !numPassesLeft) {
        selectNow = true;
#ifdef BLIS_DEBUG
        printf("REL: CREATE: maxTimeReached %d, numPassesLeft %d\n",
               maxTimeReached, numPassesLeft);
#endif
    }


    // Store first time objects.
    std::vector<BlisObjectInt *> firstObjects;

    // Store infeasible objects.
    std::vector<BlisObjectInt *> infObjects;

    // TODO: check if sorting is expensive.
    std::multimap<double, BlisObjectInt*, BlisPseuoGreater> sortedObjects;

    double objValue = solver->getObjSense() * solver->getObjValue();

    const double * lower = solver->getColLower();
    const double * upper = solver->getColUpper();

    int lookAhead = dynamic_cast<BlisParams*>
                    (model->BlisPar())->entry(BlisParams::lookAhead);

    BlisObjectInt * intObject = NULL;

    //------------------------------------------------------
    // Backup solver status and mark hot start.
    //-----------------------------------------------------

    saveSolution = new double[numCols];
    memcpy(saveSolution, solver->getColSolution(), numCols*sizeof(double));
    saveLower = new double[numCols];
    saveUpper = new double[numCols];
    memcpy(saveLower, lower, numCols * sizeof(double));
    memcpy(saveUpper, upper, numCols * sizeof(double));

    //------------------------------------------------------
    // Find the infeasible objects.
    // NOTE: we might go round this loop twice if we are feed in
    //       a "feasible" solution.
    //------------------------------------------------------

    for (pass = 0; pass < 2; ++pass) {

        numInfs = 0;

        BcpsObject * object = NULL;


        infObjects.clear();
        firstObjects.clear();
//.........这里部分代码省略.........

int * analyze(OsiClpSolverInterface * solverMod, int & numberChanged,
		     double & increment, bool changeInt,
		     CoinMessageHandler * generalMessageHandler, bool noPrinting)
{
    bool noPrinting_ = noPrinting;
    OsiSolverInterface * solver = solverMod->clone();
    char generalPrint[200];
    if (0) {
        // just get increment
        CbcModel model(*solver);
        model.analyzeObjective();
        double increment2 = model.getCutoffIncrement();
        printf("initial cutoff increment %g\n", increment2);
    }
    const double *objective = solver->getObjCoefficients() ;
    const double *lower = solver->getColLower() ;
    const double *upper = solver->getColUpper() ;
    int numberColumns = solver->getNumCols() ;
    int numberRows = solver->getNumRows();
    double direction = solver->getObjSense();
    int iRow, iColumn;

    // Row copy
    CoinPackedMatrix matrixByRow(*solver->getMatrixByRow());
    const double * elementByRow = matrixByRow.getElements();
    const int * column = matrixByRow.getIndices();
    const CoinBigIndex * rowStart = matrixByRow.getVectorStarts();
    const int * rowLength = matrixByRow.getVectorLengths();

    // Column copy
    CoinPackedMatrix  matrixByCol(*solver->getMatrixByCol());
    const double * element = matrixByCol.getElements();
    const int * row = matrixByCol.getIndices();
    const CoinBigIndex * columnStart = matrixByCol.getVectorStarts();
    const int * columnLength = matrixByCol.getVectorLengths();

    const double * rowLower = solver->getRowLower();
    const double * rowUpper = solver->getRowUpper();

    char * ignore = new char [numberRows];
    int * changed = new int[numberColumns];
    int * which = new int[numberRows];
    double * changeRhs = new double[numberRows];
    memset(changeRhs, 0, numberRows*sizeof(double));
    memset(ignore, 0, numberRows);
    numberChanged = 0;
    int numberInteger = 0;
    for (iColumn = 0; iColumn < numberColumns; iColumn++) {
        if (upper[iColumn] > lower[iColumn] + 1.0e-8 && solver->isInteger(iColumn))
            numberInteger++;
    }
    bool finished = false;
    while (!finished) {
        int saveNumberChanged = numberChanged;
        for (iRow = 0; iRow < numberRows; iRow++) {
            int numberContinuous = 0;
            double value1 = 0.0, value2 = 0.0;
            bool allIntegerCoeff = true;
            double sumFixed = 0.0;
            int jColumn1 = -1, jColumn2 = -1;
            for (CoinBigIndex j = rowStart[iRow]; j < rowStart[iRow] + rowLength[iRow]; j++) {
                int jColumn = column[j];
                double value = elementByRow[j];
                if (upper[jColumn] > lower[jColumn] + 1.0e-8) {
                    if (!solver->isInteger(jColumn)) {
                        if (numberContinuous == 0) {
                            jColumn1 = jColumn;
                            value1 = value;
                        } else {
                            jColumn2 = jColumn;
                            value2 = value;
                        }
                        numberContinuous++;
                    } else {
                        if (fabs(value - floor(value + 0.5)) > 1.0e-12)
                            allIntegerCoeff = false;
                    }
                } else {
                    sumFixed += lower[jColumn] * value;
                }
            }
            double low = rowLower[iRow];
            if (low > -1.0e20) {
                low -= sumFixed;
                if (fabs(low - floor(low + 0.5)) > 1.0e-12)
                    allIntegerCoeff = false;
            }
            double up = rowUpper[iRow];
            if (up < 1.0e20) {
                up -= sumFixed;
                if (fabs(up - floor(up + 0.5)) > 1.0e-12)
                    allIntegerCoeff = false;
            }
            if (!allIntegerCoeff)
                continue; // can't do
            if (numberContinuous == 1) {
                // see if really integer
                // This does not allow for complicated cases
                if (low == up) {
                    if (fabs(value1) > 1.0e-3) {
//.........这里部分代码省略.........

int main (int argc, const char *argv[])
{

  OsiClpSolverInterface solver1;
  //#define USE_OSI_NAMES
#ifdef USE_OSI_NAMES
  // Say we are keeping names (a bit slower this way)
  solver1.setIntParam(OsiNameDiscipline,1);
#endif
  // Read in model using argv[1]
  // and assert that it is a clean model
  std::string mpsFileName;
#if defined(SAMPLEDIR)
  mpsFileName = SAMPLEDIR "/p0033.mps";
#else
  if (argc < 2) {
    fprintf(stderr, "Do not know where to find sample MPS files.\n");
    exit(1);
  }
#endif
  if (argc>=2) mpsFileName = argv[1];
  int numMpsReadErrors = solver1.readMps(mpsFileName.c_str(),"");
  assert(numMpsReadErrors==0);

  // Strip off integer information and save
  int numberColumns = solver1.getNumCols();
  char * integer = new char[numberColumns];
  int i;
  for (i=0;i<numberColumns;i++) {
    if (solver1.isInteger(i)) {
      integer[i]=1;
      solver1.setContinuous(i);
    } else {
      integer[i]=0;
    }
  }
  // Pass to Cbc initialize defaults 
  CbcModel model(solver1);    
  CbcMain0(model);

  // Solve just to show there are no integers
  model.branchAndBound();
  // Set cutoff etc back in model and solver
  model.resetModel();
  // Solver was cloned so get it
  OsiSolverInterface * solver = model.solver();
  // Put back integers.  Here the user could do anything really
#define ADD_DIRECTLY
#ifndef ADD_DIRECTLY
  for (i=0;i<numberColumns;i++) {
    if (integer[i])
      solver->setInteger(i);
  }
#else
  CbcObject ** objects = new CbcObject * [ numberColumns];
  int n=0;
  for (i=0;i<numberColumns;i++) {
    if (integer[i]) {
      CbcSimpleIntegerDynamicPseudoCost * newObject =
	new CbcSimpleIntegerDynamicPseudoCost(&model,i);
      objects[n++]=newObject;
    }
  }
  model.addObjects(n,objects);
  for (i=0;i<n;i++)
    delete objects[i];
  delete [] objects;
#endif
  delete [] integer;
  /* Now go into code for standalone solver
     Could copy arguments and add -quit at end to be safe
     but this will do
  */
  if (argc>2) {
    CbcMain1(argc-1,argv+1,model);
  } else {
    const char * argv2[]={"driver3","-solve","-quit"};
    CbcMain1(3,argv2,model);
  }

  // Print solution if finished (could get from model.bestSolution() as well

  if (solver->getObjValue()*solver->getObjSense()<1.0e50) {
    
    const double * solution = solver->getColSolution();
    
    int iColumn;
    std::cout<<std::setiosflags(std::ios::fixed|std::ios::showpoint)<<std::setw(14);
    
    std::cout<<"--------------------------------------"<<std::endl;
#ifdef USE_OSI_NAMES
    
    for (iColumn=0;iColumn<numberColumns;iColumn++) {
      double value=solution[iColumn];
      if (fabs(value)>1.0e-7&&solver->isInteger(iColumn)) 
	std::cout<<std::setw(6)<<iColumn<<" "<<std::setw(8)<<setiosflags(std::ios::left)<<solver->getColName(iColumn)
		 <<resetiosflags(std::ios::adjustfield)<<std::setw(14)<<" "<<value<<std::endl;
    }
#else
    // names may not be in current solver - use original
//.........这里部分代码省略.........

// See if rounding will give solution
// Sets value of solution
// Assumes rhs for original matrix still okay
// At present only works with integers 
// Fix values if asked for
// Returns 1 if solution, 0 if not
int
AbcRounding::solution(double & solutionValue,
		      double * betterSolution)
{
    // Get a copy of original matrix (and by row for rounding);
    matrix_ = *(model_->solver()->getMatrixByCol());
    matrixByRow_ = *(model_->solver()->getMatrixByRow());
    seed_=1;

    OsiSolverInterface * solver = model_->solver();
    const double * lower = solver->getColLower();
    const double * upper = solver->getColUpper();
    const double * rowLower = solver->getRowLower();
    const double * rowUpper = solver->getRowUpper();
    const double * solution = solver->getColSolution();
    const double * objective = solver->getObjCoefficients();
    double integerTolerance = 1.0e-5;
    //model_->getDblParam(AbcModel::AbcIntegerTolerance);
    double primalTolerance;
    solver->getDblParam(OsiPrimalTolerance, primalTolerance);

    int numberRows = matrix_.getNumRows();

    int numberIntegers = model_->numberIntegers();
    const int * integerVariable = model_->integerVariable();
    int i;
    double direction = solver->getObjSense();
    double newSolutionValue = direction * solver->getObjValue();
    int returnCode = 0;

    // Column copy
    const double * element = matrix_.getElements();
    const int * row = matrix_.getIndices();
    const int * columnStart = matrix_.getVectorStarts();
    const int * columnLength = matrix_.getVectorLengths();
    // Row copy
    const double * elementByRow = matrixByRow_.getElements();
    const int * column = matrixByRow_.getIndices();
    const int * rowStart = matrixByRow_.getVectorStarts();
    const int * rowLength = matrixByRow_.getVectorLengths();

    // Get solution array for heuristic solution
    int numberColumns = solver->getNumCols();
    double * newSolution = new double [numberColumns];
    memcpy(newSolution, solution, numberColumns * sizeof(double));

    double * rowActivity = new double[numberRows];
    memset(rowActivity, 0, numberRows*sizeof(double));
    for (i = 0; i < numberColumns; i++) {
	int j;
	double value = newSolution[i];
	if (value) {
	    for (j = columnStart[i];
		 j < columnStart[i] + columnLength[i]; j++) {
		int iRow = row[j];
		rowActivity[iRow] += value*element[j];
	    }
	}
    }
    // check was feasible - if not adjust (cleaning may move)
    for (i = 0; i < numberRows; i++) {
	if(rowActivity[i] < rowLower[i]) {
	    //assert (rowActivity[i]>rowLower[i]-1000.0*primalTolerance);
	    rowActivity[i] = rowLower[i];
	} else if(rowActivity[i] > rowUpper[i]) {
	    //assert (rowActivity[i]<rowUpper[i]+1000.0*primalTolerance);
	    rowActivity[i] = rowUpper[i];
	}
    }
    for (i = 0; i < numberIntegers; i++) {
	int iColumn = integerVariable[i];
	double value = newSolution[iColumn];
	if (fabs(floor(value + 0.5) - value) > integerTolerance) {
	    double below = floor(value);
	    double newValue = newSolution[iColumn];
	    double cost = direction * objective[iColumn];
	    double move;
	    if (cost > 0.0) {
		// try up
		move = 1.0 - (value - below);
	    } else if (cost < 0.0) {
		// try down
		move = below - value;
	    } else {
		// won't be able to move unless we can grab another variable
		// just for now go down
		move = below-value;
	    }
	    newValue += move;
	    newSolution[iColumn] = newValue;
	    newSolutionValue += move * cost;
	    int j;
	    for (j = columnStart[iColumn];
		 j < columnStart[iColumn] + columnLength[iColumn]; j++) {
//.........这里部分代码省略.........

BlisReturnStatus
BlisStrongBranch(BlisModel *model, double objValue, int colInd, double x,
		 const double *saveLower, const double *saveUpper,
		 bool &downKeep, bool &downFinished, double &downDeg,
		 bool &upKeep, bool &upFinished, double &upDeg)
{
    BlisReturnStatus status = BlisReturnStatusOk;
    int lpStatus = 0;

    int j, numIntInfDown, numObjInfDown;

    double newObjValue;
    
    OsiSolverInterface * solver = model->solver();
    
    int numCols = solver->getNumCols();
    const double * lower = solver->getColLower();
    const double * upper = solver->getColUpper();

    // Restore bounds
    int numDiff = 0;

    BlisSolution* ksol = NULL;

    int ind = model->getIntObjIndices()[colInd];
    BlisObjectInt *intObj = dynamic_cast<BlisObjectInt *>(model->objects(ind));
    
#ifdef BLIS_DEBUG_MORE
    for (j = 0; j < numCols; ++j) {
	if (saveLower[j] != lower[j]) {
	    //solver->setColLower(j, saveLower[j]);
            ++numDiff;
	}
	if (saveUpper[j] != upper[j]) {
	    //solver->setColUpper(j, saveUpper[j]);
            ++numDiff;
	}
    }
    std::cout << "BEFORE: numDiff = " << numDiff << std::endl;
#endif	 
   
    //------------------------------------------------------
    // Branching down.
    //------------------------------------------------------

    solver->setColUpper(colInd, floor(x));
    solver->solveFromHotStart();
    
    newObjValue = solver->getObjSense() * solver->getObjValue();
    downDeg = newObjValue - objValue;
    
    if (solver->isProvenOptimal()) {
	lpStatus = 0; // optimal
#ifdef BLIS_DEBUG_MORE
        printf("STRONG: COL[%d]: downDeg=%g, x=%g\n", colInd, downDeg, x);
#endif
        // Update pseudocost
        intObj->pseudocost().update(-1, downDeg, x);
        model->setSharedObjectMark(ind);        

        // Check if ip feasible
        ksol = model->feasibleSolution(numIntInfDown, numObjInfDown);
        if (ksol) {
#ifdef BLIS_DEBUG_MORE
            printf("STRONG:Down:found a feasible solution\n");
#endif
            
            model->storeSolution(BlisSolutionTypeStrong, ksol);
	    downKeep = false;
        }
	else {
	    downKeep = true;
	}
	downFinished = true;
    }
    else if (solver->isIterationLimitReached() && 
	     !solver->isDualObjectiveLimitReached()) {
	lpStatus = 2;      // unknown 
	downKeep = true;
	downFinished = false;
    }
    else {
        downDeg = 1.0e20;
	lpStatus = 1; // infeasible
	downKeep = false;
	downFinished = false;
    }       
            
#ifdef BLIS_DEBUG_MORE
    std::cout << "Down: lpStatus = " << lpStatus << std::endl;
#endif
    
    // restore bounds
    numDiff = 0;
    for (j = 0; j < numCols; ++j) {
	if (saveLower[j] != lower[j]) {
	    solver->setColLower(j, saveLower[j]);
            ++numDiff;
	}
	if (saveUpper[j] != upper[j]) {
//.........这里部分代码省略.........

/*
  First tries setting a variable to better value.  If feasible then
  tries setting others.  If not feasible then tries swaps
  Returns 1 if solution, 0 if not 
  The main body of this routine implements an O((q^2)/2) brute force search
  around the current solution, for q = number of integer variables. Call this
  the inc/dec heuristic.  For each integer variable x<i>, first decrement the
  value. Then, for integer variables x<i+1>, ..., x<q-1>, try increment and
  decrement. If one of these permutations produces a better solution,
  remember it.  Then repeat, with x<i> incremented. If we find a better
  solution, update our notion of current solution and continue.

  The net effect is a greedy walk: As each improving pair is found, the
  current solution is updated and the search continues from this updated
  solution.

  Way down at the end, we call solutionFix, which will create a drastically
  restricted problem based on variables marked as used, then do mini-BaC on
  the restricted problem. This can occur even if we don't try the inc/dec
  heuristic. This would be more obvious if the inc/dec heuristic were broken
  out as a separate routine and solutionFix had a name that reflected where
  it was headed.

  The return code of 0 is grossly overloaded, because it maps to a return
  code of 0 from solutionFix, which is itself grossly overloaded. See
  comments in solutionFix and in CbcHeuristic::smallBranchAndBound.
  */
int
CbcHeuristicLocal::solution(double & solutionValue,
                            double * betterSolution)
{
/*
  Execute only if a new solution has been discovered since the last time we
  were called.
*/

    numCouldRun_++;
    // See if frequency kills off idea
    int swap = swap_%100;
    int skip = swap_/100;
    int nodeCount = model_->getNodeCount();
    if (nodeCount<lastRunDeep_+skip && nodeCount != lastRunDeep_+1) 
      return 0;
    if (numberSolutions_ == model_->getSolutionCount() &&
	(numberSolutions_ == howOftenShallow_ ||
	 nodeCount < lastRunDeep_+2*skip))
        return 0;
    howOftenShallow_ = numberSolutions_;
    numberSolutions_ = model_->getSolutionCount();
    if (nodeCount<lastRunDeep_+skip ) 
      return 0;
    lastRunDeep_ = nodeCount;
    howOftenShallow_ = numberSolutions_;

    if ((swap%10) == 2) {
        // try merge
        return solutionFix( solutionValue, betterSolution, NULL);
    }
/*
  Exclude long (column), thin (row) systems.

  Given the n^2 nature of the search, more than 100,000 columns could get
  expensive. But I don't yet see the rationale for the second part of the
  condition (cols > 10*rows). And cost is proportional to number of integer
  variables --- shouldn't we use that?

  Why wait until we have more than one solution?
*/
    if ((model_->getNumCols() > 100000 && model_->getNumCols() >
            10*model_->getNumRows()) || numberSolutions_ <= 1)
        return 0; // probably not worth it
    // worth trying

    OsiSolverInterface * solver = model_->solver();
    const double * rowLower = solver->getRowLower();
    const double * rowUpper = solver->getRowUpper();
    const double * solution = model_->bestSolution();
/*
  Shouldn't this test be redundant if we've already checked that
  numberSolutions_ > 1? Stronger: shouldn't this be an assertion?
*/
    if (!solution)
        return 0; // No solution found yet
    const double * objective = solver->getObjCoefficients();
    double primalTolerance;
    solver->getDblParam(OsiPrimalTolerance, primalTolerance);

    int numberRows = matrix_.getNumRows();

    int numberIntegers = model_->numberIntegers();
    const int * integerVariable = model_->integerVariable();

    int i;
    double direction = solver->getObjSense();
    double newSolutionValue = model_->getObjValue() * direction;
    int returnCode = 0;
    numRuns_++;
    // Column copy
    const double * element = matrix_.getElements();
    const int * row = matrix_.getIndices();
//.........这里部分代码省略.........

//.........这里部分代码省略.........
        for (i = 0; i < n; i++) {
            int iColumn = integerVariable[i];
            if (used_[iColumn] <= allow) {
                newSolver->setColUpper(iColumn, colLower[iColumn]);
                nFix2++;
            } else {
                break;
            }
        }
        delete [] which;
        nFix += nFix2;
#ifdef CLP_INVESTIGATE2
        printf("Number fixed increased from %d to %d\n",
               nFix - nFix2, nFix);
#endif
    }
    if (nFix*10 > numberIntegers) {
        returnCode = smallBranchAndBound(newSolver, numberNodes_, newSolution, objectiveValue,
                                         objectiveValue, "CbcHeuristicLocal");
 /*
  -2 is return due to user event, and -1 is overloaded with what look to be
  two contradictory meanings.
*/
       if (returnCode < 0) {
            returnCode = 0; // returned on size
            int numberColumns = newSolver->getNumCols();
            int numberContinuous = numberColumns - numberIntegers;
            if (numberContinuous > 2*numberIntegers &&
                    nFix*10 < numberColumns) {
#define LOCAL_FIX_CONTINUOUS
#ifdef LOCAL_FIX_CONTINUOUS
                //const double * colUpper = newSolver->getColUpper();
                const double * colLower = newSolver->getColLower();
                int nAtLb = 0;
                //double sumDj=0.0;
                const double * dj = newSolver->getReducedCost();
                double direction = newSolver->getObjSense();
                for (int iColumn = 0; iColumn < numberColumns; iColumn++) {
                    if (!newSolver->isInteger(iColumn)) {
                        if (!used_[iColumn]) {
                            //double djValue = dj[iColumn]*direction;
                            nAtLb++;
                            //sumDj += djValue;
                        }
                    }
                }
                if (nAtLb) {
                    // fix some continuous
                    double * sort = new double[nAtLb];
                    int * which = new int [nAtLb];
                    //double threshold = CoinMax((0.01*sumDj)/static_cast<double>(nAtLb),1.0e-6);
                    int nFix2 = 0;
                    for (int iColumn = 0; iColumn < numberColumns; iColumn++) {
                        if (!newSolver->isInteger(iColumn)) {
                            if (!used_[iColumn]) {
                                double djValue = dj[iColumn] * direction;
                                if (djValue > 1.0e-6) {
                                    sort[nFix2] = -djValue;
                                    which[nFix2++] = iColumn;
                                }
                            }
                        }
                    }
                    CoinSort_2(sort, sort + nFix2, which);
                    int divisor = 2;
                    nFix2 = CoinMin(nFix2, (numberColumns - nFix) / divisor);
                    for (int i = 0; i < nFix2; i++) {
                        int iColumn = which[i];
                        newSolver->setColUpper(iColumn, colLower[iColumn]);
                    }
                    delete [] sort;
                    delete [] which;
#ifdef CLP_INVESTIGATE2
                    printf("%d integers have zero value, and %d continuous fixed at lb\n",
                           nFix, nFix2);
#endif
                    returnCode = smallBranchAndBound(newSolver,
                                                     numberNodes_, newSolution,
                                                     objectiveValue,
                                                     objectiveValue, "CbcHeuristicLocal");
                    if (returnCode < 0)
                        returnCode = 0; // returned on size
                }
#endif
            }
        }
    }
/*
  If the result is complete exploration with a solution (3) or proven
  infeasibility (2), we could generate a cut (the AI folks would call it a
  nogood) to prevent us from going down this route in the future.
*/
    if ((returnCode&2) != 0) {
        // could add cut
        returnCode &= ~2;
    }

    delete newSolver;
    return returnCode;
}