Solve(Opti​misationPr​oblem), intlinprog and cplex disagree

I am trying to run a MILP. I made an OptimisationProblem which I can run with Solve(OptimisationProblem). Teh answer this gives me is correct. However when I run the same probelm in cplex (for calculation speed inprovements) or just intlinprog , I get a different (incorrect) answer.

%% using solve and timing
tic
[sol,fval,exitflag,output]=solve(scheduleprob);
toc

%% converting the problem
SP=prob2struct(scheduleprob);

%% using intlinprog and timing
tic
[sol2,fval2, exitflag2, output2] = intlinprog(SP.f,SP.intcon,SP.Aineq,SP.bineq,...
    SP.Aeq,SP.beq,SP.lb,SP.ub,SP.x0,SP.options);
toc


%% creating the ctype array
ctype(1:numel(SP.f))='C';
ctype(SP.intcon)='I';

%$ using cplex and timing
tic
[sol3,fval3, exitflag3, output3] = cplexmilp(SP.f, SP.Aineq, SP.bineq, ...
    SP.Aeq, SP.beq,[], [], [], SP.lb, SP.ub, ctype, [], SP.options);
toc

fval,fval2,fval3

when I run this I get the following output:

LP:                Optimal objective value is -33898.197234.                                        

Cut Generation:    Applied 19 implication cuts, 5 mir cuts,                                         
                   13 Gomory cuts, 4 cover cuts,                                                    
                   and 4 flow cover cuts.                                                           
                   Lower bound is -31811.836555.                                                    

Heuristics:        Found 1 solution using diving.                                                   
                   Upper bound is -25000.000000.                                                    
                   Relative gap is 27.25%.                                                         

Cut Generation:    Applied 5 Gomory cuts,                                                           
                   8 implication cuts, and 1 flow cover cut.                                        
                   Lower bound is -31789.993780.                                                    
                   Relative gap is 27.16%.                                                         

Branch and Bound:

   nodes     total   num int        integer       relative                                          
explored  time (s)  solution           fval        gap (%)                                         
     161      9.87         2  -2.773075e+04   1.277874e+01                                          
     441     11.88         3  -2.849219e+04   6.454425e+00                                          
     651     13.01         4  -2.859219e+04   1.618574e+00                                          
     860     13.52         4  -2.859219e+04   0.000000e+00                                          

Optimal solution found.

Intlinprog stopped because the objective value is within a gap tolerance of the
optimal value, options.AbsoluteGapTolerance = 0 (the default value). The intcon
variables are integer within tolerance, options.IntegerTolerance = 1e-05 (the
default value).

Elapsed time is 13.978858 seconds.
LP:                Optimal objective value is -33898.197234.                                        

Cut Generation:    Applied 19 implication cuts, 5 mir cuts,                                         
                   13 Gomory cuts, 4 cover cuts,                                                    
                   and 4 flow cover cuts.                                                           
                   Lower bound is -31811.836555.                                                    

Heuristics:        Found 1 solution using diving.                                                   
                   Upper bound is -25000.000000.                                                    
                   Relative gap is 27.25%.                                                         

Cut Generation:    Applied 5 Gomory cuts,                                                           
                   8 implication cuts, and 1 flow cover cut.                                        
                   Lower bound is -31789.993780.                                                    
                   Relative gap is 27.16%.                                                         

Branch and Bound:

   nodes     total   num int        integer       relative                                          
explored  time (s)  solution           fval        gap (%)                                         
     161      9.45         2  -2.773075e+04   1.277874e+01                                          
     441     11.40         3  -2.849219e+04   6.454425e+00                                          
     651     12.51         4  -2.859219e+04   1.618574e+00                                          
     860     13.05         4  -2.859219e+04   0.000000e+00                                          

Optimal solution found.

Intlinprog stopped because the objective value is within a gap tolerance of the
optimal value, options.AbsoluteGapTolerance = 0 (the default value). The intcon
variables are integer within tolerance, options.IntegerTolerance = 1e-05 (the
default value).

Elapsed time is 13.083694 seconds.
Elapsed time is 0.740988 seconds.
fval =
   3.5922e+03
fval2 =
  -2.8592e+04
fval3 =
  -2.8592e+04
  

the optimal solution for the problem should be 3592.2 so only the "solve" command got the right answer. Remarkable is that during calculation the solve command does show the value reached by the other two commands. Is there a step I am missing?

[Pepijn Baart]

When I use

[sol4,fval4,exitflag4,output4]=solve(scheduleprob,'Solver','intlinprog');

I also get the correct answer, I think that means it is not the solver that is the problem.

[Torsten]

And does the solution from "solve" satisfy all constraints you imposed ?

[Pepijn Baart]

Yes it does. I have the same problem when I use other optimisation problems. It is an maximisation problem is there maybe something in the prob2struct command that does something weird for that? making the optimisation function (SP.f) negative does not fix this.

[Pepijn Baart]

switching the problem to a minimisation problem (and inverting the objective) gives the negative version of the answer for the solve function (-3.5922e+03), as expected, but the same answer as before for the other two options (-2.8592e+04) which is wierd.

Answer 1

[Matt J]

If scheduleprob contains an additive offset value f0 in the objective function, then solve() will account for it, but intlinprog will not,

https://www.mathworks.com/help/optim/ug/optim.problemdef.optimizationproblem.prob2struct.html#mw_d550a523-5099-44f9-8457-187518647a55

and cplex may not (though that is probably a question for the cplex developers). Were the optimized variables different among the solvers as well? I expect not.

[Pepijn Baart]

You are my hero, there was indeed an offset of 25000 in the objective function. Do you perhaps know how to transform the answer given by intlinprog (which is an array of 1xAmountOfVariables) to the format of the anser of solve (which is a nice struct with the variable names)?

[Matt J]

This is already being done for you by solve, which merely calls intlinprog internally.

[Pepijn Baart]

Yes that what I meant, I like the way solve presents that answers, and I want cplexmilp to do the same. Since the anser cplexmilp gives me is the same as the answer intlingprog gives me, and solve can convert this to a nice answer, is there a way to do this with the answer of cplexmilp?

[Matt J]

Assuming cplexmilp gives the variables in the same order as intlinprog, you could use findindex to re-organize cplexmilp output into problem-based output.

[Pepijn Baart]

I tried to do as you said, but since I use multiple named variabels it does not seem to work:

SI  = optimvar('SI', 1,       1,       J,N,'Type','integer','Lowerbound',0,'Upperbound',1);
SO  = optimvar('SO', 1,       1,       J,N,'Type','integer','Lowerbound',0,'Upperbound',1);
SD  = optimvar('SD', 1,       1,       KD+J,N,'Type','integer','Lowerbound',0,'Upperbound',1);
X   = optimvar('X',  1,       numel(I),J,N,'Type','integer','Lowerbound',0,'Upperbound',1);
Y   = optimvar('Y',  1,       numel(I),J,N,'Type','integer','Lowerbound',0,'Upperbound',1);
test= optimvar('test',1,      numel(I),J,N,'Type','integer','Lowerbound',0,'Upperbound',1);

Z   = optimvar('Z',  1,       1,       J,N,'Lowerbound',0,'Upperbound',1);
E   = optimvar('E',  1,       1,       J,N,'Lowerbound',0,'Upperbound',1);
W   = optimvar('W',  1,       1,       J,N,'Lowerbound',0,'Upperbound',1);
T   = optimvar('T',  1,       1,       1,N,'Lowerbound',0,'Upperbound',H);
TLB = optimvar('TLB',1,       1,       J,N,'Lowerbound',0);
TEE = optimvar('TEE',1,       1,       J,N,'Lowerbound',0);
TS  = optimvar('TS', 1,       1,       J,N,'Lowerbound',0);
TW  = optimvar('TW', 1,       1,       J,N,'Lowerbound',0);
BS  = optimvar('BS', 1,       numel(I),J,N,'Lowerbound',0);
BE  = optimvar('BE', 1,       numel(I),J,N,'Lowerbound',0);
BP  = optimvar('BP', 1,       numel(I),J,N,'Lowerbound',0);
II  = optimvar('II', numel(M),1,       J,N,'Lowerbound',0);
IO  = optimvar('IO', numel(M),1,       J,N,'Lowerbound',0);
IV  = optimvar('IV', numel(M),1,       K+J,N,'Lowerbound',0);
FVU = optimvar('FVU',  numel(M),K+J,J,N,'Lowerbound',0);
FUV = optimvar('FUV',  numel(M),J,K+J,N,'Lowerbound',0);
FUU = optimvar('FUU',  numel(M),J,J,N,'Lowerbound',0);
FVV = optimvar('FVV',  numel(M),K+J,K+J,N,'Lowerbound',0);
Q   = optimvar('Q',  1,1,numel(R),N,'Lowerbound',0);

The solve() command gives me a sol that look slike this:

which is very nice and usable. with the findindex command each variable starts counting at 1 which does not help recreating this structure from the initial array.

[Pepijn Baart]

Or should i make a new question for this?

[Matt J]

I think it would be worth posting a new question, because it looks like a worthwhile topic for wider discussion. It is strange to me that a command for this doesn't already exist.

[Pepijn Baart]

Thank you for your help. I have created a New question.