自学教程：C++ CoinMin函数代码示例

//-------------------------------------------------------------------// Generate Stored cuts//------------------------------------------------------------------- void CglStoredUser::generateCuts(const OsiSolverInterface & si, OsiCuts & cs,			     const CglTreeInfo info) const{  // Get basic problem information  const double * solution = si.getColSolution();  if (info.inTree&&info.pass>numberPasses_) {    // only continue if integer feasible    int numberColumns=si.getNumCols();     int i;    const double * colUpper = si.getColUpper();    const double * colLower = si.getColLower();    int numberAway=0;    for (i=0;i<numberColumns;i++) {      double value = solution[i];      // In case slightly away from bounds      value = CoinMax(colLower[i],value);      value = CoinMin(colUpper[i],value);      if (si.isInteger(i)&&fabs(value-fabs(value+0.5))>1.0e-5) 	numberAway++;    }    if (numberAway)      return; // let code branch  }  int numberRowCuts = cuts_.sizeRowCuts();  for (int i=0;i<numberRowCuts;i++) {    const OsiRowCut * rowCutPointer = cuts_.rowCutPtr(i);    double violation = rowCutPointer->violated(solution);    if (violation>=requiredViolation_)      cs.insert(*rowCutPointer);  }}

voidCglClique::generateCuts(const OsiSolverInterface& si, OsiCuts & cs,			const CglTreeInfo info){   int i;   bool has_petol_set = petol != -1.0;   if (! has_petol_set)      si.getDblParam(OsiPrimalTolerance, petol);   int numberOriginalRows = si.getNumRows();   if (info.inTree&&justOriginalRows_)     numberOriginalRows = info.formulation_rows;   int numberRowCutsBefore = cs.sizeRowCuts();   // First select which rows/columns we are interested in.   if (!setPacking_) {      selectFractionalBinaries(si);      if (!sp_orig_row_ind) {	 selectRowCliques(si,numberOriginalRows);      }   } else {      selectFractionals(si);      delete[] sp_orig_row_ind;      sp_numrows = numberOriginalRows;      //sp_numcols = si.getNumCols();      sp_orig_row_ind = new int[sp_numrows];      for (i = 0; i < sp_numrows; ++i)	 sp_orig_row_ind[i] = i;   }   // Just original rows   if (justOriginalRows_&&info.inTree)      sp_numrows = CoinMin(info.formulation_rows,sp_numrows);        createSetPackingSubMatrix(si);   fgraph.edgenum = createNodeNode();   createFractionalGraph();   cl_indices = new int[sp_numcols];   cl_del_indices = new int[sp_numcols];   if (do_row_clique)      find_rcl(cs);   if (do_star_clique)      find_scl(cs);   if (!info.inTree&&((info.options&4)==4||((info.options&8)&&!info.pass))) {     int numberRowCutsAfter = cs.sizeRowCuts();     for (int i=numberRowCutsBefore;i<numberRowCutsAfter;i++)       cs.rowCutPtr(i)->setGloballyValid();   }   delete[] cl_indices;     cl_indices = 0;   delete[] cl_del_indices; cl_del_indices = 0;   deleteFractionalGraph();   delete[] node_node;      node_node = 0;   deleteSetPackingSubMatrix();   if (! has_petol_set)      petol = -1;}

CbcBranchingObject *CbcNWay::createCbcBranch(OsiSolverInterface * solver, const OsiBranchingInformation * /*info*/, int /*way*/){    int numberFree = 0;    int j;    //OsiSolverInterface * solver = model_->solver();    const double * solution = model_->testSolution();    const double * lower = solver->getColLower();    const double * upper = solver->getColUpper();    int * list = new int[numberMembers_];    double * sort = new double[numberMembers_];    for (j = 0; j < numberMembers_; j++) {        int iColumn = members_[j];        double value = solution[iColumn];        value = CoinMax(value, lower[iColumn]);        value = CoinMin(value, upper[iColumn]);        if (upper[iColumn] > lower[iColumn]) {            double distance = upper[iColumn] - value;            list[numberFree] = j;            sort[numberFree++] = distance;        }    }    assert (numberFree);    // sort    CoinSort_2(sort, sort + numberFree, list);    // create object    CbcBranchingObject * branch;    branch = new CbcNWayBranchingObject(model_, this, numberFree, list);    branch->setOriginalObject(this);    delete [] list;    delete [] sort;    return branch;}

// Return "down" estimatedoubleCbcSimpleIntegerPseudoCost::downEstimate() const{    OsiSolverInterface * solver = model_->solver();    const double * solution = model_->testSolution();    const double * lower = solver->getColLower();    const double * upper = solver->getColUpper();    double value = solution[columnNumber_];    value = CoinMax(value, lower[columnNumber_]);    value = CoinMin(value, upper[columnNumber_]);    if (upper[columnNumber_] == lower[columnNumber_]) {        // fixed        return 0.0;    }    double integerTolerance =        model_->getDblParam(CbcModel::CbcIntegerTolerance);    double below = floor(value + integerTolerance);    double above = below + 1.0;    if (above > upper[columnNumber_]) {        above = below;        below = above - 1;    }    double downCost = CoinMax((value - below) * downPseudoCost_, 0.0);    return downCost;}

// Creates a branching objectCbcBranchingObject * CbcSimpleIntegerFixed::createBranch(OsiSolverInterface * solver,					    const OsiBranchingInformation * info, int way)  {  const double * solution = model_->testSolution();  const double * lower = solver->getColLower();  const double * upper = solver->getColUpper();  double value = solution[columnNumber_];  value = CoinMax(value, lower[columnNumber_]);  value = CoinMin(value, upper[columnNumber_]);  assert (upper[columnNumber_]>lower[columnNumber_]);  if (!model_->hotstartSolution()) {    double nearest = floor(value+0.5);    double integerTolerance =     model_->getDblParam(CbcModel::CbcIntegerTolerance);    if (fabs(value-nearest)<integerTolerance) {      // adjust value      if (nearest!=upper[columnNumber_])	value = nearest+2.0*integerTolerance;      else	value = nearest-2.0*integerTolerance;    }  } else {    const double * hotstartSolution = model_->hotstartSolution();    double targetValue = hotstartSolution[columnNumber_];    if (way>0)      value = targetValue-0.1;    else      value = targetValue+0.1;  }  CbcBranchingObject * branch = new CbcIntegerBranchingObject(model_,columnNumber_,way,					     value);  branch->setOriginalObject(this);  return branch;}

// Fill most valuesvoidClpDualRowSteepest::fill(const ClpDualRowSteepest& rhs){     state_ = rhs.state_;     mode_ = rhs.mode_;     persistence_ = rhs.persistence_;     assert (model_->numberRows() == rhs.model_->numberRows());     model_ = rhs.model_;     assert(model_);     int number = model_->numberRows();     if (rhs.savedWeights_)          number = CoinMin(number, rhs.savedWeights_->capacity());     if (rhs.infeasible_ != NULL) {          if (!infeasible_)               infeasible_ = new CoinIndexedVector(rhs.infeasible_);          else               *infeasible_ = *rhs.infeasible_;     } else {          delete infeasible_;          infeasible_ = NULL;     }     if (rhs.weights_ != NULL) {          if (!weights_)               weights_ = new double[number];          ClpDisjointCopyN(rhs.weights_, number, weights_);     } else {          delete [] weights_;          weights_ = NULL;     }     if (rhs.alternateWeights_ != NULL) {          if (!alternateWeights_)               alternateWeights_ = new CoinIndexedVector(rhs.alternateWeights_);          else               *alternateWeights_ = *rhs.alternateWeights_;     } else {          delete alternateWeights_;          alternateWeights_ = NULL;     }     if (rhs.savedWeights_ != NULL) {          if (!savedWeights_)               savedWeights_ = new CoinIndexedVector(rhs.savedWeights_);          else               *savedWeights_ = *rhs.savedWeights_;     } else {          delete savedWeights_;          savedWeights_ = NULL;     }     if (rhs.dubiousWeights_) {          assert(model_);          int number = model_->numberRows();          if (!dubiousWeights_)               dubiousWeights_ = new int[number];          ClpDisjointCopyN(rhs.dubiousWeights_, number, dubiousWeights_);     } else {          delete [] dubiousWeights_;          dubiousWeights_ = NULL;     }}

/* Fills in array with problem name  */COINLIBAPI void COINLINKAGEClp_problemName(Clp_Simplex * model, int maxNumberCharacters, char * array){     std::string name = model->model_->problemName();     maxNumberCharacters = CoinMin(maxNumberCharacters,     				   ((int) strlen(name.c_str()))+1) ;     strncpy(array, name.c_str(), maxNumberCharacters - 1);     array[maxNumberCharacters-1] = '/0';}

int main(int argc, const char *argv[]){     // Get model in some way     ClpInterior model;     // Open graph and parameter files     //FILE *fpin = fopen("./g.graph","r");     //FILE *fpp = fopen("./gparm","r");     FILE *fpin = fopen("./g.tiny", "r");     FILE *fpp = fopen("./gparm.tiny", "r");     assert(fpin);     assert(fpp);     myPdco stuff(model, fpin, fpp);     Info info;     Outfo outfo;     Options options;     /*      *     Set the input parameters for LSQR.      */     options.gamma = stuff.getD1();     options.delta = stuff.getD2();     options.MaxIter = 40;     options.FeaTol = 5.0e-4;     options.OptTol = 5.0e-4;     options.StepTol = 0.99;     //  options.x0min = 10.0/num_cols;     options.x0min = 0.01;     options.z0min = 0.01;     options.mu0 = 1.0e-6;     options.LSmethod = 3;   // 1=Cholesky    2=QR    3=LSQR     options.LSproblem = 1;  // See below     options.LSQRMaxIter = 999;     options.LSQRatol1 = 1.0e-3; // Initial  atol     options.LSQRatol2 = 1.0e-6; // Smallest atol (unless atol1 is smaller)     options.LSQRconlim = 1.0e12;     info.atolmin = options.LSQRatol2;     info.LSdamp = 0.0;     // These are already set?     model.xsize_ = 50.0 / (model.numberColumns());     model.xsize_ = CoinMin(1.0, model.xsize_);     /*      *     Solve the test problem      */     model.pdco(&stuff, options, info, outfo);     /*      *     Examine the results.      *     Print the residual norms RNORM and ARNORM given by LSQR, and then compute      */     return 0;}

  /// heuristic method  int  GuessHeuristic::solution(double &solutionValue, double *betterSolution)  {    // Get pointer to pseudo costs object    const BonChooseVariable* chooseMethod = dynamic_cast<BonChooseVariable*>(model_->branchingMethod()->chooseMethod());    if (!chooseMethod) {      (*model_->messageHandler()) << "Can't get pseudo costs!!!/n";      solutionValue = model_->getCurrentMinimizationObjValue();      return -1;    }    const OsiPseudoCosts& pseudoCosts = chooseMethod->pseudoCosts();    int numberObjects = pseudoCosts.numberObjects();    assert(numberObjects == model_->numberObjects());    const double* upTotalChange = pseudoCosts.upTotalChange();    const double* downTotalChange = pseudoCosts.downTotalChange();    const int* upNumber = pseudoCosts.upNumber();    const int* downNumber = pseudoCosts.downNumber();    double sumUpTot = 0.;    int numberUpTot = 0;    double sumDownTot = 0.;    int numberDownTot = 0;    for (int i=0;i<numberObjects;i++) {      sumUpTot += upTotalChange[i];      numberUpTot += upNumber[i];      sumDownTot += downTotalChange[i];      numberDownTot += downNumber[i];    }    if (!numberUpTot || !numberDownTot) {      // don't have ANY pseudo-costs information yet      solutionValue = COIN_DBL_MAX;      return -1;    }    double upAvrg=sumUpTot/numberUpTot;    double downAvrg=sumDownTot/numberDownTot;    OsiObject** object =  model_->objects();    solutionValue = model_->getCurrentMinimizationObjValue();    for (int iObj = 0; iObj < numberObjects; iObj++) {      //printf("%3d upest=%e uptot=%e upnum=%d downest=%e downtot=%e downnum=%d  ", iObj, object[iObj]->upEstimate(), upTotalChange[iObj], upNumber[iObj], object[iObj]->downEstimate(), downTotalChange[iObj], downNumber[iObj]);      double upEstimate = upNumber[iObj] ? object[iObj]->upEstimate()*upTotalChange[iObj]/upNumber[iObj] : object[iObj]->upEstimate()*upAvrg;      double downEstimate = downNumber[iObj] ? object[iObj]->downEstimate()*downTotalChange[iObj]/downNumber[iObj] : object[iObj]->downEstimate()*downAvrg;      //printf("up=%e down=%e/n", upEstimate, downEstimate);      solutionValue += CoinMin(upEstimate,downEstimate);    }    //printf("solutionValue = %e/n", solutionValue);    return -1;  }

// Resize objectivevoid ClpQuadraticObjective::resize(int newNumberColumns){  if (numberColumns_ != newNumberColumns) {    int newExtended = newNumberColumns + (numberExtendedColumns_ - numberColumns_);    int i;    double *newArray = new double[newExtended];    if (objective_)      CoinMemcpyN(objective_, CoinMin(newExtended, numberExtendedColumns_), newArray);    delete[] objective_;    objective_ = newArray;    for (i = numberColumns_; i < newNumberColumns; i++)      objective_[i] = 0.0;    if (gradient_) {      newArray = new double[newExtended];      if (gradient_)        CoinMemcpyN(gradient_, CoinMin(newExtended, numberExtendedColumns_), newArray);      delete[] gradient_;      gradient_ = newArray;      for (i = numberColumns_; i < newNumberColumns; i++)        gradient_[i] = 0.0;    }    if (quadraticObjective_) {      if (newNumberColumns < numberColumns_) {        int *which = new int[numberColumns_ - newNumberColumns];        int i;        for (i = newNumberColumns; i < numberColumns_; i++)          which[i - newNumberColumns] = i;        quadraticObjective_->deleteRows(numberColumns_ - newNumberColumns, which);        quadraticObjective_->deleteCols(numberColumns_ - newNumberColumns, which);        delete[] which;      } else {        quadraticObjective_->setDimensions(newNumberColumns, newNumberColumns);      }    }    numberColumns_ = newNumberColumns;    numberExtendedColumns_ = newExtended;  }}

DSP_RTN_CODE BdMaster::init(){	BGN_TRY_CATCH	/** create problem */	createProblem();	/** set node limit */	si_->setNodeLimit(par_->getIntParam("NODE_LIM"));	/** set print level */	si_->setPrintLevel(CoinMin(par_->getIntParam("LOG_LEVEL") + 2, 5));	END_TRY_CATCH_RTN(;,DSP_RTN_ERR)

template <typename T> voidCoinDenseVector<T>::resize(int newsize, T value){  if (newsize != nElements_){    assert(newsize > 0);    T *newarray = new T[newsize];    int cpysize = CoinMin(newsize, nElements_);    CoinMemcpyN( elements_,cpysize,newarray);    delete[] elements_;    elements_ = newarray;    nElements_ = newsize;    for(int i=cpysize; i<newsize; i++)      elements_[i] = value;  }}

//-------------------------------------------------------------------// Copy constructor//-------------------------------------------------------------------ClpDualRowSteepest::ClpDualRowSteepest (const ClpDualRowSteepest & rhs)     : ClpDualRowPivot(rhs){     state_ = rhs.state_;     mode_ = rhs.mode_;     persistence_ = rhs.persistence_;     model_ = rhs.model_;     if ((model_ && model_->whatsChanged() & 1) != 0) {          int number = model_->numberRows();          if (rhs.savedWeights_)               number = CoinMin(number, rhs.savedWeights_->capacity());          if (rhs.infeasible_) {               infeasible_ = new CoinIndexedVector(rhs.infeasible_);          } else {               infeasible_ = NULL;          }          if (rhs.weights_) {               weights_ = new double[number];               ClpDisjointCopyN(rhs.weights_, number, weights_);          } else {               weights_ = NULL;          }          if (rhs.alternateWeights_) {               alternateWeights_ = new CoinIndexedVector(rhs.alternateWeights_);          } else {               alternateWeights_ = NULL;          }          if (rhs.savedWeights_) {               savedWeights_ = new CoinIndexedVector(rhs.savedWeights_);          } else {               savedWeights_ = NULL;          }          if (rhs.dubiousWeights_) {               assert(model_);               int number = model_->numberRows();               dubiousWeights_ = new int[number];               ClpDisjointCopyN(rhs.dubiousWeights_, number, dubiousWeights_);          } else {               dubiousWeights_ = NULL;          }     } else {          infeasible_ = NULL;          weights_ = NULL;          alternateWeights_ = NULL;          savedWeights_ = NULL;          dubiousWeights_ = NULL;     }}

// Resize objectivevoidClpLinearObjective::resize(int newNumberColumns){     if (numberColumns_ != newNumberColumns) {          int i;          double * newArray = new double[newNumberColumns];          if (objective_)               CoinMemcpyN(objective_, CoinMin(newNumberColumns, numberColumns_), newArray);          delete [] objective_;          objective_ = newArray;          for (i = numberColumns_; i < newNumberColumns; i++)               objective_[i] = 0.0;          numberColumns_ = newNumberColumns;     }}

// Creates a branching object from this infeasible object.BcpsBranchObject * BlisObjectInt::createBranchObject(BcpsModel *m, int direction) const{    BlisModel *model = dynamic_cast<BlisModel* >(m);    OsiSolverInterface * solver = model->solver();        double integerTolerance = model->BlisPar()->entry(BlisParams::integerTol);        const double * solution = solver->getColSolution();    const double * lower = solver->getColLower();    const double * upper = solver->getColUpper();        double value = solution[columnIndex_];    //std::cout << "COL"<< columnIndex_ << ": x = " << value << std::endl;        // Force value in bounds.    value = CoinMax(value, lower[columnIndex_]);    value = CoinMin(value, upper[columnIndex_]);        double nearest = floor(value + 0.5);        assert (upper[columnIndex_] > lower[columnIndex_]);        int hotstartStrategy = model->getHotstartStrategy();        if (hotstartStrategy <= 0) {        if (fabs(value - nearest) < integerTolerance) {            // Already integeral.            std::cout << "ERROR: COL" << columnIndex_ << ": x=" << value                       << ", nearest=" << nearest                       << ", intTol=" << integerTolerance << std::endl;            assert(0);        }    }     else {	const double * incumbent = model->incumbent();	double targetValue = incumbent[columnIndex_];	if (direction > 0) {	    value = targetValue - 0.1;	}	else {	    value = targetValue + 0.1;	}    }        return new BlisBranchObjectInt(model, objectIndex_, direction, value);}

voidOsiTestSolverInterface::applyColCut(const OsiColCut& cc){   int i;   const double* lb_elem = cc.lbs().getElements();   const int* lb_ind = cc.lbs().getIndices();   for (i = cc.lbs().getNumElements() - 1; i >= 0; --i) {      collower_[lb_ind[i]] = CoinMax(collower_[lb_ind[i]], lb_elem[i]);   }      const double* ub_elem = cc.ubs().getElements();   const int* ub_ind = cc.ubs().getIndices();   for (i = cc.ubs().getNumElements() - 1; i >= 0; --i) {      colupper_[ub_ind[i]] = CoinMin(colupper_[ub_ind[i]], ub_elem[i]);   }}

doubleCbcTreeArray::getBestPossibleObjective(){    double bestPossibleObjective = 1e100;    for (int i = 0 ; i < static_cast<int> (nodes_.size()) ; i++) {        if (nodes_[i] && nodes_[i]->objectiveValue() < bestPossibleObjective) {            bestPossibleObjective = nodes_[i]->objectiveValue();        }    }    if (lastNode_) {        bestPossibleObjective = CoinMin(bestPossibleObjective, lastNode_->objectiveValue());    }    CbcCompareDefault * compareDefault    = dynamic_cast<CbcCompareDefault *> (comparison_.test_);    assert (compareDefault);    compareDefault->setBestPossible(bestPossibleObjective);    return bestPossibleObjective;}

// Returns true if current solution satsifies one side of branchbool OsiSolverBranch::feasibleOneWay(const OsiSolverInterface & solver) const{  bool feasible = false;  int numberColumns = solver.getNumCols();  const double * columnLower = solver.getColLower();  const double * columnUpper = solver.getColUpper();  const double * columnSolution = solver.getColSolution();  double primalTolerance;  solver.getDblParam(OsiPrimalTolerance,primalTolerance);  for (int base = 0; base<4; base +=2) {    feasible=true;    int i;    for (i=start_[base];i<start_[base+1];i++) {      int iColumn = indices_[i];      if (iColumn<numberColumns) {        double value = CoinMax(bound_[i],columnLower[iColumn]);        if (columnSolution[iColumn]<value-primalTolerance) {          feasible=false;          break;        }      } else {        abort(); // do later (other stuff messed up anyway - e.g. CBC)      }    }    if (!feasible)      break;    for (i=start_[base+1];i<start_[base+2];i++) {      int iColumn = indices_[i];      if (iColumn<numberColumns) {        double value = CoinMin(bound_[i],columnUpper[iColumn]);        if (columnSolution[iColumn]>value+primalTolerance) {          feasible=false;          break;        }      } else {        abort(); // do later (other stuff messed up anyway - e.g. CBC)      }    }    if (feasible)      break; // OK this way  }  return feasible;}

OsiRowCut *BlisConstraint::createOsiRowCut(){    double lower = CoinMax(getLbHard(), getLbSoft());    double upper = CoinMin(getUbHard(), getUbSoft());    OsiRowCut * cut = new OsiRowCut;    if (!cut) {        /* Out of memory. */        throw CoinError("Out of Memory", "Blis_constraintToOsiCut", "NONE");    }    assert(size_ > 0);    cut->setLb(lower);    cut->setUb(upper);    cut->setRow(size_, indices_, values_);    return cut;}

bool CglGMI::checkDynamism(const double* cutElem, const int* cutIndex,			   int cutNz) {  double min = param.getINFINIT();  double max = 0.0;  double val = 0.0;  for (int i = 0; i < cutNz; ++i) {    if (!isZero(cutElem[i])) {      val = fabs(cutElem[i]);      min = CoinMin(min, val);      max = CoinMax(max, val);    }  }  if (max > min*param.getMAXDYN()) {#if defined GMI_TRACE_CLEAN    printf("Max elem %g, min elem %g, dyn %g; cut discarded/n", max, min, max/min);#endif    return false;  }  else{    return true;  }  } /* checkDynamism */

// Useful constructor (which are integer indices)CbcGeneralDepth::CbcGeneralDepth (CbcModel * model, int maximumDepth)        : CbcGeneral(model),        maximumDepth_(maximumDepth),        maximumNodes_(0),        whichSolution_(-1),        numberNodes_(0),        nodeInfo_(NULL){    assert(maximumDepth_ < 1000000);    if (maximumDepth_ > 0)        maximumNodes_ = (1 << maximumDepth_) + 1 + maximumDepth_;    else if (maximumDepth_ < 0)        maximumNodes_ = 1 + 1 - maximumDepth_;    else        maximumNodes_ = 0;#define MAX_NODES 100    maximumNodes_ = CoinMin(maximumNodes_, 1 + maximumDepth_ + MAX_NODES);    if (maximumNodes_) {        nodeInfo_ = new ClpNodeStuff();        nodeInfo_->maximumNodes_ = maximumNodes_;        ClpNodeStuff * info = nodeInfo_;        // for reduced costs and duals        info->solverOptions_ |= 7;        if (maximumDepth_ > 0) {            info->nDepth_ = maximumDepth_;        } else {            info->nDepth_ = - maximumDepth_;            info->solverOptions_ |= 32;        }        ClpNode ** nodeInfo = new ClpNode * [maximumNodes_];        for (int i = 0; i < maximumNodes_; i++)            nodeInfo[i] = NULL;        info->nodeInfo_ = nodeInfo;    } else {        nodeInfo_ = NULL;    }}

//.........这里部分代码省略.........     }     double * lower = new double[numberRows];     double * upper = new double[numberRows];     const double * rowLower = model.rowLower();     const double * rowUpper = model.rowUpper();     for (iColumn = 0; iColumn < numberColumns; iColumn++)          mark[iColumn] = -1;     CoinPackedMatrix * matrix = model.matrix();     // get row copy     CoinPackedMatrix rowCopy = *matrix;     rowCopy.reverseOrdering();     const int * column = rowCopy.getIndices();     const int * rowLength = rowCopy.getVectorLengths();     const CoinBigIndex * rowStart = rowCopy.getVectorStarts();     const double * element = rowCopy.getElements();     int putGub = numberRows;     int putNonGub = numberRows;     int * rowIsGub = new int [numberRows];     for (iRow = numberRows - 1; iRow >= 0; iRow--) {          bool gubRow = true;          int first = numberColumns + 1;          int last = -1;          for (int j = rowStart[iRow]; j < rowStart[iRow] + rowLength[iRow]; j++) {               if (element[j] != 1.0) {                    gubRow = false;                    break;               } else {                    int iColumn = column[j];                    if (mark[iColumn] >= 0) {                         gubRow = false;                         break;                    } else {                         last = CoinMax(last, iColumn);                         first = CoinMin(first, iColumn);                    }               }          }          if (last - first + 1 != rowLength[iRow] || !gubRow) {               which[--putNonGub] = iRow;               rowIsGub[iRow] = 0;          } else {               for (int j = rowStart[iRow]; j < rowStart[iRow] + rowLength[iRow]; j++) {                    int iColumn = column[j];                    mark[iColumn] = iRow;               }               rowIsGub[iRow] = -1;               putGub--;               gubStart[putGub] = first;               gubEnd[putGub] = last + 1;               lower[putGub] = rowLower[iRow];               upper[putGub] = rowUpper[iRow];               whichGub[putGub] = iRow;          }     }     int numberNonGub = numberRows - putNonGub;     int numberGub = numberRows - putGub;     if (numberGub > 0) {          printf("** %d gub rows/n", numberGub);          int numberNormal = 0;          const int * row = matrix->getIndices();          const int * columnLength = matrix->getVectorLengths();          const CoinBigIndex * columnStart = matrix->getVectorStarts();          const double * elementByColumn = matrix->getElements();          int numberElements = 0;          bool doLower = false;          bool doUpper = false;

// inner part of diveint 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) {//.........这里部分代码省略.........

// Fix other variables at boundsintCbcHeuristicDive::fixOtherVariables(OsiSolverInterface * solver,                                    const double * solution,                                    PseudoReducedCost * candidate,                                    const double * random){    const double * lower = solver->getColLower();    const double * upper = solver->getColUpper();    double integerTolerance = model_->getDblParam(CbcModel::CbcIntegerTolerance);    double primalTolerance;    solver->getDblParam(OsiPrimalTolerance, primalTolerance);    int numberIntegers = model_->numberIntegers();    const int * integerVariable = model_->integerVariable();    const double* reducedCost = solver->getReducedCost();    // fix other integer variables that are at their bounds    int cnt = 0;#ifdef GAP    double direction = solver->getObjSense(); // 1 for min, -1 for max    double gap = 1.0e30;#endif#ifdef GAP    double cutoff = model_->getCutoff() ;    if (cutoff < 1.0e20 && false) {        double direction = solver->getObjSense() ;        gap = cutoff - solver->getObjValue() * direction ;        gap *= 0.1; // Fix more if plausible        double tolerance;        solver->getDblParam(OsiDualTolerance, tolerance) ;        if (gap <= 0.0)            gap = tolerance;        gap += 100.0 * tolerance;    }    int nOverGap = 0;#endif    int numberFree = 0;    int numberFixedAlready = 0;    for (int i = 0; i < numberIntegers; i++) {        int iColumn = integerVariable[i];        if (upper[iColumn] > lower[iColumn]) {            numberFree++;            double value = solution[iColumn];            if (fabs(floor(value + 0.5) - value) <= integerTolerance) {                candidate[cnt].var = iColumn;                candidate[cnt++].pseudoRedCost =                    fabs(reducedCost[iColumn] * random[i]);#ifdef GAP                if (fabs(reducedCost[iColumn]) > gap)                    nOverGap++;#endif            }        } else {            numberFixedAlready++;        }    }#ifdef GAP    int nLeft = maxNumberToFix - numberFixedAlready;#ifdef CLP_INVESTIGATE4    printf("cutoff %g obj %g nover %d - %d free, %d fixed/n",           cutoff, solver->getObjValue(), nOverGap, numberFree,           numberFixedAlready);#endif    if (nOverGap > nLeft && true) {        nOverGap = CoinMin(nOverGap, nLeft + maxNumberToFix / 2);        maxNumberToFix += nOverGap - nLeft;    }#else#ifdef CLP_INVESTIGATE4    printf("cutoff %g obj %g - %d free, %d fixed/n",           model_->getCutoff(), solver->getObjValue(), numberFree,           numberFixedAlready);#endif#endif    return cnt;}