Solve an integer programming problem with LinearProgramming?

Questions
, ,
1

I am trying to solve an Integer programming problem of modest size in Mathematica (the same problem is solved optimally in CPLEX in a few seconds), however I am getting strange results:

  1. Mathematica gives a solution very quickly but one of the restrictions is not satisfied (even though it does not issue an error or warning with this respect).

  2. Even stranger, when I relax some restrictions LinearProgramming finds no feasible solution.

     c = ReadList["matrizesjjC.txt"][[1]];
 M = ReadList["matrizesjjM.txt"][[1]];
 b = ReadList["matrizesjjb.txt"][[1]];
 sol = LinearProgramming[c, M, b, Table[{0, 1}, {i, 1, Length[c]}], Integers]

Any ideas? I attach the datafiles, just in case it may be useful. Unless something very strange is going on, this makes LinearProgramming totally unreliable!

P.S. The problem has in fact no admissible solution.

2

this is what I get using version 11

c = ReadList["matrizesjjC.txt"][[1]];
M = ReadList["matrizesjjM.txt"][[1]];
b = ReadList["matrizesjjb.txt"][[1]];


sol = LinearProgramming[c, M, b, Table[{0, 1}, {i, 1, Length[c]}], 
  Integers]

LinearProgramming::lpip: Warning: integer linear programming will use a machine-precision approximation of the inputs.

LinearProgramming::lpsnf: No solution can be found that satisfies the constraints.

LinearProgramming::lpip: Warning: integer linear programming will use a machine-precision approximation of the inputs.

LinearProgramming::lpsnf: No solution can be found that satisfies the constraints.
3

First of all, thank you for taking the time to test this example. I am also using version 11.0 but the result is not the same. This is what I get.

In[1]:= c = ReadList["~/Desktop/matrizesjjC.txt"][[1]];
M = ReadList["~/Desktop/matrizesjjM.txt"][[1]];
b = ReadList["~/Desktop/matrizesjjb.txt"][[1]];

In[4]:= LinearProgramming[c, M, b, 
Table[{0, 1}, {i, 1, Length[c]}], Integers]

During evaluation of In[4]:= LinearProgramming::lpip: Warning: integer linear programming will use a machine-precision approximation of the inputs.

Out[4]= {0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, \
0, 0, 1, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, \
0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, \
0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, \
0, 0, 0, 0, 0, 0, 0}

LinearProgramming returns a solution that is unfeasible, as the last restriction states that the sum of all variables must be 3.

I am using the Linux version… I always assumed that versions would be equivalent across platforms but this is obviously not the case. I’ll try to use other platforms/versions to identify the affected systems.

4

I used the Windows version

5

What is the actual solution? I have tried this a few ways and gotten no feasible values each time.

6

There is no feasible solution… The problem here is that with the Linux version of mathematica 11.0 LinearProgramming[] returns a solution that does not satisfy the restrictions. This makes me very suspicious of results returned for real instances where we typically just want to rely on the solver without further checks.

7

The Macintosh does the exact same as the Linux version. Why are you trying to replace CPLEX (it is probably the best MIP Solver available)? How are you interfacing to CPLEX now?

8

Ah, no solution. When you state “the same problem is solved optimally in CPLEX”, that rather strongly implies that there is a solution. So people might actually spend time searching for it, wondering why there seems not to be any solution.

As for the return of an infeasible solution, I believe that may be from an uncaught exception deep in the version of the COIN-CLP library used in Mathematica. I also believe we are going to upgrade this at some point but I neither know when nor whether that will fix this particular issue.

When doing ILPs and wanting to avoid this sort of issue I often just code a branch-and-prune loop using NMinimize along with some heuristic cut-off value for deciding when an approximate number is or is not an integer. For this example it could be done as below.

AbsoluteTiming[
 Module[{vars, obj, x, vecs, constraints, program, vals, stack, soln, 
   solns = {}, badvar, varvals, val, counter = 1, var, tmp, 
   min = Infinity, numsols = 0, mylist},
  vars = Array[x, Length[c]];
  constraints = 
   Join[Thread[Most[mat].vars <= 1], 
    Thread[0 <= vars <= 1], {Total[vars] == 3}];
  obj = c.vars;
  program = {obj, constraints, 0}; stack = {program, {}};
  While[stack =!= {},
   program = First[stack]; stack = Last[stack];
   If[Last[program] > min, Continue[]];
   counter++; program = Drop[program, -1];
   Quiet[soln = NMinimize[program, vars], NMinimize::"nsol"];
   If[! FreeQ[soln, Indeterminate], Continue[]];
   {tmp, soln} = soln; If[tmp > min, Continue[]];
   soln = Chop[vars /. soln];
   badvar = 
    Position[soln, (a_ /; Chop[a - Round[a], 10^(-5)] =!= 0), {1}, 1, 
     Heads -> False];
   If[badvar == {},
    If[tmp < min, solns = {}; min = tmp; numsols = 0];
    numsols++;
    solns = {Apply[mylist, soln], solns};
    ,
    badvar = badvar[[1, 1]]; constraints = Last[program];
    stack = {{obj, Append[constraints, vars[[badvar]] == 0], tmp}, 
      stack};
    stack = {{obj, Append[constraints, vars[[badvar]] == 1], tmp}, 
      stack};];];
  {numsols, counter, Flatten[solns, Infinity] /. mylist -> List}
  ]]

(* Out[194]= {63.196984, {0, 270, {}}} *)

Upshot: 63 seconds, no solutions, 270 attempts on relaxed subproblems.

9

You are totally right! Sorry for that… I was dealing with several instances and the one I posted was one with no solution. Also thank you for your contribution with the “branch-and-prune”.

10

Well, I’m not trying to replace CPLEX. Mathematica is a possibility for an Operations research course for undergraduate management students. They have natural aversion to any reasonable interface to CPLEX and the excel solver is too limited in terms of number of variables/restrictions.

11

You might want to look at

https://courses.edx.org/courses/course-v1:MITx+15.053x+3T2016/info

It is the archive of an online course on “Optimization Methods in Business Analytics”
using the Excel solver and the open-source Julia language.