CglSimpleRounding Class Reference

#include <CglSimpleRounding.hpp>

Inheritance diagram for CglSimpleRounding:

List of all members.

Public Member Functions
Generate Cuts
virtual void	generateCuts (const OsiSolverInterface &si, OsiCuts &cs) const
Constructors and destructors
	CglSimpleRounding ()
	Default constructor.
	CglSimpleRounding (const CglSimpleRounding &)
	Copy constructor.
CglSimpleRounding &	operator= (const CglSimpleRounding &rhs)
	Assignment operator.
virtual	~CglSimpleRounding ()
	Destructor.
Private Member Functions
Private methods
bool	deriveAnIntegerRow (const OsiSolverInterface &si, int rowIndex, const CoinShallowPackedVector &matrixRow, CoinPackedVector &irow, double &b, bool *negative) const
	Derive a <= inequality in integer variables from the rowIndex-th constraint.
int	power10ToMakeDoubleAnInt (int size, const double *x, double dataTol) const
Greatest common denominators methods
int	gcd (int a, int b) const
	Returns the greatest common denominator of two positive integers, a and b.
int	gcdv (int n, const int *const vi) const
Private Attributes
Private member data
double	epsilon_
	A value within an epsilon_ neighborhood of 0 is considered to be 0.
Friends
void	CglSimpleRoundingUnitTest (const OsiSolverInterface *siP, const std::string mpdDir)

---

Detailed Description

Simple Rounding Cut Generator Class

This class generates simple rounding cuts via the following method: For each contraint, attempt to derive a <= inequality in all integer variables by netting out any continuous variables. Divide the resulting integer inequality through by the greatest common denomimator (gcd) of the lhs coefficients. Round down the rhs.

Warning: Use with careful attention to data precision.

(Reference: Nemhauser and Wolsey, Integer and Combinatorial Optimization, 1988, pg 211.)

Definition at line 26 of file CglSimpleRounding.hpp.

---

Member Function Documentation

int CglSimpleRounding::gcdv (  int  n, const int *const  vi )  const [inline, private]

Returns the greatest common denominator of a vector of positive integers, vi, of length n. Definition at line 141 of file CglSimpleRounding.hpp. References gcd(). Referenced by generateCuts(). { if (n==0) abort(); if (n==1) return vi[0]; int retval=gcd(vi[0], vi[1]); for (int i=2; i<n; i++){ retval=gcd(retval,vi[i]); } return retval; }

void CglSimpleRounding::generateCuts (  const OsiSolverInterface &  si, OsiCuts &  cs )  const [virtual]

Generate simple rounding cuts for the model accessed through the solver interface. Insert generated cuts into the cut set cs. Implements CglCutGenerator. Definition at line 20 of file CglSimpleRounding.cpp. References deriveAnIntegerRow(), epsilon_, gcd(), gcdv(), CoinPackedVector::getElements(), CoinPackedVector::getIndices(), CoinPackedVector::getNumElements(), CoinPackedVector::insert(), power10ToMakeDoubleAnInt(), and CoinPackedVector::setVector(). { int nRows=si.getNumRows(); // number of rows in the coefficient matrix int nCols=si.getNumCols(); // number of columns in the coefficient matrix int rowIndex; // index into the constraint matrix stored in row // order CoinPackedVector irow; // "integer row": working space to hold the integer // <= inequality derived from the rowIndex-th // constraint double b=0; // working space for the rhs of integer <= inequality bool * negative= new bool[nCols]; // negative[i]= true if coefficient of the // ith variable is negative and false // otherwise int k; // dummy iterator variable for ( k=0; k<nCols; k++ ) negative[k] = false; const CoinPackedMatrix * rowCopy = si.getMatrixByRow(); // row copy: matrix stored in row order // Main loop: // // For every row in the matrix, // // if we can derive a valid <= inequality in integer variables, then // // try to construct a simple rounding cut from the integer inequality. // // Add the resulting cut to the set of cuts. // for (rowIndex=0; rowIndex<nRows; rowIndex++){ // Only look at tight rows // double * pi=ekk_rowduals(model); // if (fabs(pi[row]) < epsilon_){ // continue; // } // Try to derive an <= inequality in integer variables from the row // by netting out the continuous variables. // Store the value and the sign of the coefficients separately: // irow.getElements() contains the absolute values of the coefficients. // negative is a boolean vector indicating the sign of the coeffcients. // b is the rhs of the <= integer inequality if (!deriveAnIntegerRow( si, rowIndex, rowCopy->getVector(rowIndex), irow, b, negative)) { // Reset local data for the next iteration of the rowIndex-loop for(k=0; k<irow.getNumElements(); k++) negative[irow.getIndices()[k]]=false; irow.setVector(0,NULL,NULL); continue; } // Euclid's greatest common divisor (gcd) algorithm applies to positive // INTEGERS. // Determine the power of 10 needed, so that multipylying the integer // inequality through by 10**power makes all coefficients essentially // integral. int power = power10ToMakeDoubleAnInt(irow.getNumElements(),irow.getElements(),epsilon_*1.0e-4); // Now a vector to store the integer-ized values. For instance, // if x[i] is .66 and power is 1000 then xInt[i] will be 660 int * xInt = NULL; if (power >=0) { xInt = new int[irow.getNumElements()]; double dxInt; // a double version of xInt for error trapping #ifdef CGL_DEBUG printf("The (double) coefficients and their integer-ized counterparts:\n"); #endif for (k=0; k<irow.getNumElements(); k++){ dxInt = irow.getElements()[k]*pow(10.0,power); xInt[k]= (int) (dxInt+0.5); // Need to add the 0.5 // so that a dxInt=9.999 will give a xInt=1 #ifdef CGL_DEBUG printf("%g %g \n",irow.getElements()[k],dxInt); #endif } } else { // If overflow is detected, one warning message is printed and // the row is skipped. #ifdef CGL_DEBUG printf("SimpleRounding: Warning: Overflow detected \n"); printf(" on %i of vars in processing row %i. Row skipped.\n", -power, rowIndex); #endif // reset local data for next iteration for(k=0; k<irow.getNumElements(); k++) negative[irow.getIndices()[k]]=false; irow.setVector(0,NULL,NULL); continue; } // find greatest common divisor of the irow.elements int gcd = gcdv(irow.getNumElements(), xInt); #ifdef CGL_DEBUG printf("The gcd of xInt is %i\n",gcd); #endif // construct new cut by dividing through by gcd and // rounding down rhs and accounting for negatives CoinPackedVector cut; for (k=0; k<irow.getNumElements(); k++){ cut.insert(irow.getIndices()[k],xInt[k]/gcd); } double cutRhs = floor((b*pow(10.0,power))/gcd); // un-negate the negated variables in the cut { const int s = cut.getNumElements(); const int * indices = cut.getIndices(); double* elements = cut.getElements(); for (k=0; k<s; k++){ int column=indices[k]; if (negative[column]) { elements[k] *= -1; } } } // Create the row cut and add it to the set of cuts // It may not be violated if (fabs(cutRhs*gcd-b)> epsilon_){ // if the cut and row are different. OsiRowCut rc; rc.setRow(cut.getNumElements(),cut.getIndices(),cut.getElements()); rc.setLb(-DBL_MAX); rc.setUb(cutRhs); cs.insert(rc); #ifdef CGL_DEBUG printf("Row %i had a simple rounding cut:\n",rowIndex); printf("Cut size: %i Cut rhs: %g Index Element \n", cut.getNumElements(), cutRhs); for (k=0; k<cut.getNumElements(); k++){ printf("%i %g\n",cut.getIndices()[k], cut.getElements()[k]); } printf("\n"); #endif } // Reset local data for the next iteration of the rowIndex-loop for(k=0; k<irow.getNumElements(); k++) negative[irow.getIndices()[k]]=false; irow.setVector(0,NULL,NULL); delete [] xInt; } delete [] negative; }

int CglSimpleRounding::power10ToMakeDoubleAnInt (  int  size, const double *  x, double  dataTol )  const [private]

Given a vector of doubles, x, with size elements and a positive tolerance, dataTol, this method returns the smallest power of 10 needed so that x[i]*10**power "is integer" for all i=0,...,size-1. change of definition of dataTol so that it refers to original data, not to scaled data as that seems to lead to problems. So if xScaled is x[i]*10**power and xInt is rounded(xScaled) then fabs(xScaled-xInt) <= dataTol*10**power. This means that dataTol should be smaller - say 1.0e-12 rather tahn 1.0e-8 Returns -number of times overflowed if the power is so big that it will cause overflow (i.e. integer stored will be bigger than 2**31). Test in cut generator. Definition at line 330 of file CglSimpleRounding.cpp. Referenced by generateCuts(). { // Assumption: data precision is positive assert( dataTol > 0 ); int i; // loop iterator int maxPower=0; // maximum power of 10 used to convert any x[i] to an // integer // this is the number we are after. int power = 0; // power of 10 used to convert a particular x[i] to an // integer #ifdef OLD_MULT double intPart; // the integer part of the number #endif double fracPart; // the fractional part of the number // we keep multiplying by 10 until the fractional part is 0 // (well, really just until the factional part is less than // dataTol) // JJF - code seems to fail sometimes as multiplying by 10 - so // definition of dataTol changed - see header file const double multiplier[16]={1.0,1.0e1,1.0e2,1.0e3,1.0e4,1.0e5, 1.0e6,1.0e7,1.0e8,1.0e9,1.0e10,1.0e11, 1.0e12,1.0e13,1.0e14,1.0e15}; // Loop through every element in the array in x for (i=0; i<size; i++){ power = 0; #ifdef OLD_MULT // look at the fractional part of x[i] // FYI: if you want to modify this member function to take an input // vector x of arbitary sign, change this line below to // fracPart = modf(fabs(x[i]),&intPart); fracPart = modf(x[i],&intPart); // if the fractional part is close enough to 0 or 1, we're done with this // value while(!(fracPart < dataTol || 1-fracPart < dataTol )) { // otherwise, multiply by 10 and look at the fractional part of the // result. ++power; fracPart = fracPart*10.0; fracPart = modf(fracPart,&intPart); } #else // use fabs as safer and does no harm double value = fabs(x[i]); double scaledValue; // Do loop - always using original value to stop round off error. // If we don't find in 15 goes give up for (power=0;power<16;power++) { double tolerance = dataTol*multiplier[power]; scaledValue = value*multiplier[power]; fracPart = scaledValue-floor(scaledValue); if(fracPart < tolerance || 1.0-fracPart < tolerance ) { break; } } if (power==16||scaledValue>2147483647) { #ifdef CGL_DEBUG printf("Overflow %g => %g, power %d\n",x[i],scaledValue,power); #endif return -1; } #endif #ifdef CGL_DEBUG printf("The smallest power of 10 to make %g integral = %i\n",x[i],power); #endif // keep track of the largest power needed so that at the end of the for // loop // x[i]*10**maxPower will be integral for all i if (maxPower < power) maxPower=power; } return maxPower; }

---

Friends And Related Function Documentation

void CglSimpleRoundingUnitTest (  const OsiSolverInterface *  siP, const std::string  mpdDir )  [friend]

A function that tests the methods in the CglSimpleRounding class. The only reason for it not to be a member method is that this way it doesn't have to be compiled into the library. And that's a gain, because the library should be compiled with optimization on, but this method should be compiled with debugging. Definition at line 20 of file CglSimpleRoundingTest.cpp. { // Test default constructor { CglSimpleRounding cg; } // Test copy & assignment { CglSimpleRounding rhs; { CglSimpleRounding cg; CglSimpleRounding cgC(cg); rhs=cg; } } // Test gcd and gcdn { CglSimpleRounding cg; int v = cg.gcd(122,356); assert(v==2); v=cg.gcd(356,122); assert(v==2); v=cg.gcd(54,67); assert(v==1); v=cg.gcd(67,54); assert(v==1); v=cg.gcd(485,485); assert(v==485); v=cg.gcd(17*13,17*23); assert( v==17); v=cg.gcd(17*13*5,17*23); assert( v==17); v=cg.gcd(17*13*23,17*23); assert(v==17*23); int a[4] = {12, 20, 32, 400}; v= cg.gcdv(4,a); assert(v== 4); int b[4] = {782, 4692, 51, 2754}; v= cg.gcdv(4,b); assert(v== 17); int c[4] = {50, 40, 30, 10}; v= cg.gcdv(4,c); assert(v== 10); } // Test generate cuts method on exmip1.5.mps { CglSimpleRounding cg; OsiSolverInterface * siP = baseSiP->clone(); std::string fn = mpsDir+"exmip1.5"; siP->readMps(fn.c_str(),"mps"); OsiCuts cuts; cg.generateCuts(*siP,cuts); // there should be 3 cuts int nRowCuts = cuts.sizeRowCuts(); assert(nRowCuts==3); // get the last "sr"=simple rounding cut that was derived OsiRowCut srRowCut2 = cuts.rowCut(2); CoinPackedVector srRowCutPV2 = srRowCut2.row(); // this is what the last cut should look like: i.e. the "solution" const int solSize = 2; int solCols[solSize]={2,3}; double solCoefs[solSize]={5.0, 4.0}; OsiRowCut solRowCut; solRowCut.setRow(solSize,solCols,solCoefs); solRowCut.setLb(-DBL_MAX); solRowCut.setUb(2.0); // Test for equality between the derived cut and the solution cut // Note: testing two OsiRowCuts are equal invokes testing two // CoinPackedVectors are equal which invokes testing two doubles // are equal. Usually not a good idea to test that two doubles are equal, // but in this cut the "doubles" represent integer values. assert(srRowCut2 == solRowCut); delete siP; } // Test generate cuts method on p0033 { CglSimpleRounding cg; OsiSolverInterface * siP = baseSiP->clone(); std::string fn = mpsDir+"p0033"; siP->readMps(fn.c_str(),"mps"); OsiCuts cuts; cg.generateCuts(*siP,cuts); // p0033 is the optimal solution to p0033 int objIndices[14] = { 0, 6, 7, 9, 13, 17, 18, 22, 24, 25, 26, 27, 28, 29 }; CoinPackedVector p0033(14,objIndices,1.0); // test that none of the generated cuts // chops off the optimal solution int nRowCuts = cuts.sizeRowCuts(); OsiRowCut rcut; CoinPackedVector rpv; int i; for (i=0; i<nRowCuts; i++){ rcut = cuts.rowCut(i); rpv = rcut.row(); double p0033Sum = (rpv*p0033).sum(); double rcutub = rcut.ub(); assert (p0033Sum <= rcutub); } // test that the cuts improve the // lp objective function value siP->initialSolve(); double lpRelaxBefore=siP->getObjValue(); OsiSolverInterface::ApplyCutsReturnCode rc = siP->applyCuts(cuts); siP->resolve(); double lpRelaxAfter=siP->getObjValue(); #ifdef CGL_DEBUG printf("\n\nOrig LP min=%f\n",lpRelaxBefore); printf("Final LP min=%f\n\n",lpRelaxAfter); #endif assert( lpRelaxBefore < lpRelaxAfter ); delete siP; } }

---

The documentation for this class was generated from the following files:

• CglSimpleRounding.hpp
• CglSimpleRounding.cpp

---