[Clp] How to get integral solution for totally unimodular matrix?

[Lasse Kliemann]

* Message by -Michael Hennebry- from Tue 2009-06-23:
> On Tue, 23 Jun 2009, Lasse Kliemann wrote:
> 
> >I tried using CLP via 'OsiClpSolverInterface.initialSolve()' to
> >solve a bipartite matching problem, which is known to have a
> >totally unimodular matrix. So, Simplex method should give only
> >integral solutions AFAIK. Indeed, this works well with GLPK's
> >'glp_simplex(...)'. However, CLP on several occasions delivers
> >clearly non-integral solutions. What could be the reason for
> >this? Is there a configuration parameter to change that?
> 
> Every basic solution, feasible or not, is integral.
> Any other answer is wrong.
> You are doing smeting wrong or Clp is doing something wrong.

Let's find out. It is not too easy to reproduce the problem. 
Attached is a self-contained example that triggers it. A 
bipartite graph is constructed "randomly", then an LP expressing 
maximum-cardinality matching is solved. When a non-integral 
solution is encountered, it is printed.

I tested it on 2 machines. Pretty soon after starting, I see 
numbers printed to the screen, like:

0.0133177
0.986682
0.0133177
0.986682
0.0133177
0.986682
0.986682

The choice of parameter 'p', the edge probability, seems crucial.  
For instance, using 'p=0.5' instead of 'p=0.4' does not trigger 
the problem.
-------------- next part --------------
#include <OsiSolverInterface.hpp>
#include <OsiClpSolverInterface.hpp>
#include <CoinPackedMatrix.hpp>
#include <stdlib.h>

void doit() {
	const unsigned int n = 250; // number of vertices
	const unsigned int n0 = 104; // number of vertices in first partition
	CoinPackedMatrix *matrix = new CoinPackedMatrix(true, 0.5, 0.0);
	matrix->setDimensions(n, 0);
	const double ones[2] = {1, 1};
	int idx[2];
	unsigned int m = 0; // counter for the number of edges
	const double p = 0.4; // edge probability
	for (unsigned int i=0; i<n0; ++i)
		for (unsigned int j=n0; j<n; ++j)
			if (drand48() < p) {
				idx[0] = i;
				idx[1] = j;
				matrix->appendCol(2, &idx[0], &ones[0]);
				++m;
			}
	double *row_lb = new double[n];
	double *row_ub = new double[n];
	double *col_lb = new double[m];
	double *col_ub = new double[m];
	double *objective = new double[m];
	for (unsigned int i=0; i<n; ++i) {
		row_lb[i] = 0.0;
		row_ub[i] = 1.0;
	}
	for (unsigned int i=0; i<m; ++i) {
		col_lb[i] = 0.0;
		col_ub[i] = 1.0;
		objective[i] = 1.0;
	}
	OsiClpSolverInterface *const si = new OsiClpSolverInterface;
	si->messageHandler()->setLogLevel(0);
	si->assignProblem(matrix, col_lb, col_ub, objective, row_lb, row_ub);
	si->setObjSense(-1); // 1 min, -1 max
	si->initialSolve();
	assert(si->isProvenOptimal());
	const double *t = si->getColSolution();
	for (unsigned int i=0; i<m; ++i)
		if ( !((t[i]<1e-2) || (t[i]>1-(1e-2))) )
			std::cout << t[i] << std::endl;
	delete si;
	delete matrix;
	delete[] row_lb; delete[] row_ub; delete[] col_lb; delete[] col_ub; delete[] objective;
}

int main() {
	srand48(0);
	for (unsigned int i=0; i<9999; ++i) doit();
	return 0;
}
-------------- next part --------------
A non-text attachment was scrubbed...
Name: not available
Type: application/pgp-signature
Size: 835 bytes
Desc: not available
URL: <http://list.coin-or.org/pipermail/clp/attachments/20090624/919f6dbd/attachment.sig>