C++ CoinPackedMatrix::getNumCols方法代码示例

/* Modifies djs to allow for quadratic.
   returns quadratic offset */
CoinWorkDouble
ClpInterior::quadraticDjs(CoinWorkDouble * djRegion, const CoinWorkDouble * solution,
                          CoinWorkDouble scaleFactor)
{
     CoinWorkDouble quadraticOffset = 0.0;
#ifndef NO_RTTI
     ClpQuadraticObjective * quadraticObj = (dynamic_cast< ClpQuadraticObjective*>(objective_));
#else
     ClpQuadraticObjective * quadraticObj = NULL;
     if (objective_->type() == 2)
          quadraticObj = (static_cast< ClpQuadraticObjective*>(objective_));
#endif
     if (quadraticObj) {
          CoinPackedMatrix * quadratic = quadraticObj->quadraticObjective();
          const int * columnQuadratic = quadratic->getIndices();
          const CoinBigIndex * columnQuadraticStart = quadratic->getVectorStarts();
          const int * columnQuadraticLength = quadratic->getVectorLengths();
          double * quadraticElement = quadratic->getMutableElements();
          int numberColumns = quadratic->getNumCols();
          for (int iColumn = 0; iColumn < numberColumns; iColumn++) {
               CoinWorkDouble value = 0.0;
               for (CoinBigIndex j = columnQuadraticStart[iColumn];
                         j < columnQuadraticStart[iColumn] + columnQuadraticLength[iColumn]; j++) {
                    int jColumn = columnQuadratic[j];
                    CoinWorkDouble valueJ = solution[jColumn];
                    CoinWorkDouble elementValue = quadraticElement[j];
                    //value += valueI*valueJ*elementValue;
                    value += valueJ * elementValue;
                    quadraticOffset += solution[iColumn] * valueJ * elementValue;
               }
               djRegion[iColumn] += scaleFactor * value;
          }
     }
     return quadraticOffset;
}

/* Load in an problem by copying the arguments (the constraints on the
   rows are given by lower and upper bounds). If a pointer is NULL then the
   following values are the default:
   <ul>
   <li> <code>colub</code>: all columns have upper bound infinity
   <li> <code>collb</code>: all columns have lower bound 0 
   <li> <code>rowub</code>: all rows have upper bound infinity
   <li> <code>rowlb</code>: all rows have lower bound -infinity
   <li> <code>obj</code>: all variables have 0 objective coefficient
      </ul>
*/
void 
CoinSnapshot::loadProblem(const CoinPackedMatrix& matrix,
			  const double* collb, const double* colub,   
			  const double* obj,
			  const double* rowlb, const double* rowub,
			  bool makeRowCopy)
{
  // Keep scalars (apart from objective value etc)
  gutsOfDestructor(3+8);
  numRows_ = matrix.getNumRows();
  numCols_ = matrix.getNumCols();
  numElements_ = matrix.getNumElements();
  owned_.matrixByCol = 1;
  matrixByCol_ = new CoinPackedMatrix(matrix);
  if (makeRowCopy) {
    owned_.matrixByRow = 1;
    CoinPackedMatrix * matrixByRow = new CoinPackedMatrix(matrix);
    matrixByRow->reverseOrdering();
    matrixByRow_ = matrixByRow;
  }
  colLower_ = CoinCopyOfArray(collb,numCols_,0.0);
  colUpper_ = CoinCopyOfArray(colub,numCols_,infinity_);
  objCoefficients_ = CoinCopyOfArray(obj,numCols_,0.0);
  rowLower_ = CoinCopyOfArray(rowlb,numRows_,-infinity_);
  rowUpper_ = CoinCopyOfArray(rowub,numRows_,infinity_);
  // do rhs as well
  createRightHandSide();
}

lpHook::lpHook(const double* clb, const double* cub, const double* obj,
               const double* rhs, const char* sense,
               const CoinPackedMatrix& mat)
{
     const int colnum = mat.getNumCols();
     const int rownum = mat.getNumRows();

     colupper_ = new double[colnum];
     collower_ = new double[colnum];
     objcoeffs_ = new double[colnum];
     rhs_ = new double[rownum];
     sense_ = new char[rownum];

     std::copy(clb, clb + colnum, collower_);
     std::copy(cub, cub + colnum, colupper_);
     std::copy(obj, obj + colnum, objcoeffs_);
     std::copy(rhs, rhs + rownum, rhs_);
     std::copy(sense, sense + rownum, sense_);

     if (mat.isColOrdered()) {
          colMatrix_.copyOf(mat);
          rowMatrix_.reverseOrderedCopyOf(mat);
     } else {
          rowMatrix_.copyOf(mat);
          colMatrix_.reverseOrderedCopyOf(mat);
     }
}

void
OsiVolSolverInterface::loadProblem(const CoinPackedMatrix& matrix,
				   const double* collb, const double* colub,
				   const double* obj,
				   const char* rowsen, const double* rowrhs,   
				   const double* rowrng)
{
   gutsOfDestructor_();
   const int rownum = matrix.getNumRows();
   const int colnum = matrix.getNumCols();

   if (matrix.isColOrdered()) {
      colMatrix_ = matrix;
      colMatrixCurrent_ = true;
      rowMatrixCurrent_ = false;
      maxNumcols_ = colMatrix_.getMaxMajorDim();
      maxNumrows_ = static_cast<int>((1+colMatrix_.getExtraGap()) *
				     colMatrix_.getMinorDim());
   } else {
      rowMatrix_ = matrix;
      rowMatrixCurrent_ = true;
      colMatrixCurrent_ = false;
      maxNumcols_ = static_cast<int>((1+rowMatrix_.getExtraGap()) *
				     rowMatrix_.getMinorDim());
      maxNumrows_ = rowMatrix_.getMaxMajorDim();
   }

   initFromRhsSenseRange(rownum, rowsen, rowrhs, rowrng);
   initFromClbCubObj(colnum, collb, colub, obj);
}

   /** Equivalence.
       Two matrices are equivalent if they are both by rows or both by columns,
       they have the same dimensions, and each vector is equivalent. 
       In this method the FloatEqual function operator can be specified. 
   */
   template <class FloatEqual> bool 
   isEquivalent(const CoinPackedMatrix& rhs, const FloatEqual& eq) const
   {
      // Both must be column order or both row ordered and must be of same size
      if ((isColOrdered() ^ rhs.isColOrdered()) ||
	  (getNumCols() != rhs.getNumCols()) ||
	  (getNumRows() != rhs.getNumRows()) ||
	  (getNumElements() != rhs.getNumElements()))
	 return false;
     
      for (int i=getMajorDim()-1; i >= 0; --i) {
        CoinShallowPackedVector pv = getVector(i);
        CoinShallowPackedVector rhsPv = rhs.getVector(i);
        if ( !pv.isEquivalent(rhsPv,eq) )
          return false;
      }
      return true;
   }

 //TODO: ugh force both row/col format?
 void loadDataAndSolution(const CoinPackedMatrix & rowMatrix,
                          const CoinPackedMatrix & colMatrix,
                          const double           * collb, 
                          const double           * colub,   
                          const double           * obj,
                          const double           * rowlb, 
                          const double           * rowub,
                          const char             * colType,
                          const double           * primalSol,
                          const double             infinity){
    data_->setMatrixByRow(&rowMatrix);
    data_->setMatrixByCol(&colMatrix);      
    data_->setNrow(rowMatrix.getNumRows());
    data_->setNcol(rowMatrix.getNumCols());
    data_->setColLower(collb);
    data_->setColUpper(colub);
    data_->setRowLower(rowlb);
    data_->setRowUpper(rowub);
    data_->setObj(obj);
    data_->setColType(colType);
    data_->setPrimalSol(primalSol);
    data_->setInfinity(infinity);
    data_->initializeOtherData();
 }

void CglOddHole::generateCuts(const OsiRowCutDebugger * /*debugger*/,
			      const CoinPackedMatrix & rowCopy, 
				 const double * solution, 
			      const double * dj, OsiCuts & cs,
				 const int * suitableRow,
			      const int * fixedColumn,
			      const CglTreeInfo info,
			      bool packed)
{
  CoinPackedMatrix columnCopy = rowCopy;
  columnCopy.reverseOrdering();

  // Get basic problem information
  int nRows=columnCopy.getNumRows(); 
  int nCols=columnCopy.getNumCols(); 
  
  const int * column = rowCopy.getIndices();
  const CoinBigIndex * rowStart = rowCopy.getVectorStarts();
  const int * rowLength = rowCopy.getVectorLengths(); 
  
  const int * row = columnCopy.getIndices();
  const CoinBigIndex * columnStart = columnCopy.getVectorStarts();
  const int * columnLength = columnCopy.getVectorLengths(); 

  // we need only look at suitable rows and variables with unsatisfied 0-1
  // lookup from true row to compressed matrix
  int * mrow = new int[nRows];
  // lookup from true column to compressed
  int * lookup = new int[nCols];
  // number of columns in compressed matrix
  int nSmall=0;
  int i;
  //do lookup from true sequence to compressed
  int n=0;
  for (i=0;i<nRows;i++) {
    if (suitableRow[i]>0) {
      mrow[i]=n++;
    } else {
      mrow[i]=-1;
    }
  }
  for (i=0;i<nCols;i++) {
    if (!fixedColumn[i]) {
      lookup[i]=nSmall++;
    } else {
      lookup[i]=-1;
    }
  }
  int nSmall2=2*nSmall;
  // we don't know how big matrix will be
#define MAXELS 50000
  int maxels=MAXELS;
  //How do I do reallocs in C++?
  // 1.0 - value x(i) - value x(j) for each node pair (or reverse if cover) 
  double * cost = reinterpret_cast<double *> (malloc(maxels*sizeof(double)));
  // arc i.e. j which can be reached from i
  int * to= reinterpret_cast<int *> (malloc(maxels*sizeof(int)));
  //original row for each arc
  int * rowfound=reinterpret_cast<int *> (malloc(maxels*sizeof(int)));
  // start of each column
  int * starts=new int[2*nSmall+1];
  starts[0]=0;
  // useful array for marking if already connected
  int * mark =new int[nSmall2];
  memset(mark,0,nSmall2*sizeof(int));
  n=0; //number of elements in matrix
  for (i=0;i<nCols;i++) {
    int icol=lookup[i];
    if (icol>=0) {
      // column in compressed matrix
      int k;
      double dd=1.0000001-solution[i];
      mark[icol]=1;
      // reallocate if matrix reached size limit
      if (n+nCols>maxels) {
	maxels*=2;
	cost=reinterpret_cast<double *> (realloc(cost,maxels*sizeof(double)));
	to=reinterpret_cast<int *> (realloc(to,maxels*sizeof(int)));
	rowfound=reinterpret_cast<int *> (realloc(rowfound,maxels*sizeof(int)));
      }
      // get all other connected variables
      for (k=columnStart[i];k<columnStart[i]+columnLength[i];k++) {
	int irow=row[k];
	int jrow=mrow[irow];
	// but only if row in compressed matrix
	if (jrow>=0) {
	  int j;
	  for (j=rowStart[irow];j<rowStart[irow]+rowLength[irow];j++) {
	    int jcol=column[j];
	    int kcol=lookup[jcol];
	    if (kcol>=0&&!mark[kcol]) {
	      cost[n]=dd-solution[jcol];
	      to[n]=kcol;
	      rowfound[n++]=irow;//original row
	      mark[kcol]=1;
	    }
	  }
	}
      }
      starts[icol+1]=n;
//.........这里部分代码省略.........

//.........这里部分代码省略.........
        }
    }
    int * markKnapsack = NULL;
    double * coefficient = NULL;
    double * linear = NULL;
    int * whichRow = NULL;
    int * lookupRow = NULL;
    badModel = (canDo == 0);
    if (numberKnapsack && canDo) {
        /* double check - OK if
           no nonlinear
           nonlinear only on columns in knapsack
           nonlinear only on columns in knapsack * ONE other - same for all in knapsack
           AND that is only row connected to knapsack
           (theoretically could split knapsack if two other and small numbers)
           also ONE could be ONE expression - not just a variable
        */
        int iKnapsack;
        markKnapsack = new int [numberKnapsack];
        coefficient = new double [numberKnapsack];
        linear = new double [numberColumns];
        for (iKnapsack = 0; iKnapsack < numberKnapsack; iKnapsack++)
            markKnapsack[iKnapsack] = -1;
        if (canDo == 2) {
            for (iRow = -1; iRow < numberRows; iRow++) {
                int numberOdd;
                CoinPackedMatrix * row = coinModel.quadraticRow(iRow, linear, numberOdd);
                if (row) {
                    // see if valid
                    const double * element = row->getElements();
                    const int * column = row->getIndices();
                    const CoinBigIndex * columnStart = row->getVectorStarts();
                    const int * columnLength = row->getVectorLengths();
                    int numberLook = row->getNumCols();
                    for (int i = 0; i < numberLook; i++) {
                        int iKnapsack = whichKnapsack[i];
                        if (iKnapsack < 0) {
                            // might be able to swap - but for now can't have knapsack in
                            for (int j = columnStart[i]; j < columnStart[i] + columnLength[i]; j++) {
                                int iColumn = column[j];
                                if (whichKnapsack[iColumn] >= 0) {
                                    canDo = 0; // no good
                                    badModel = true;
                                    break;
                                }
                            }
                        } else {
                            // OK if in same knapsack - or maybe just one
                            int marked = markKnapsack[iKnapsack];
                            for (int j = columnStart[i]; j < columnStart[i] + columnLength[i]; j++) {
                                int iColumn = column[j];
                                if (whichKnapsack[iColumn] != iKnapsack && whichKnapsack[iColumn] >= 0) {
                                    canDo = 0; // no good
                                    badModel = true;
                                    break;
                                } else if (marked == -1) {
                                    markKnapsack[iKnapsack] = iColumn;
                                    marked = iColumn;
                                    coefficient[iKnapsack] = element[j];
                                    coinModel.associateElement(coinModel.columnName(iColumn), 1.0);
                                } else if (marked != iColumn) {
                                    badModel = true;
                                    canDo = 0; // no good
                                    break;
                                } else {
                                    // could manage with different coefficients - but for now ...