I have a matlab script for diet linear programming problem. I am trying to compare results of linearprog or(intlinprog) with results of cplexlp(or cplexmilp) functions provided by toolbox CPLEX for matlab by IBM here is the script
%% Defining Variables
clear;clc
Pnames = ["BEEF";
"CHK";
"FISH";
"HAM";
"MCH";
"MTL";
"SPG";
"TUR" ];
Packs = optimvar('Packs',Pnames,'Type','integer');
Packs.LowerBound = 2*ones(length(Pnames),1);
Packs.UpperBound = 10*ones(length(Pnames),1);
%% Setting the problem data
CostPerPack = [3.19;2.59;2.29;2.89;1.89;1.99;1.99;2.49];
VitA = [60;8;8;40;15;70;25;60];
VitC = [20;0;10;40;35;30;50;20];
VitB1 = [10;20;15;35;15;15;25;15];
VitB2 = [15;20;10;10;15;15;15;10];
NA = [938;2180;945;278;1182;896;1329;1397];
CAL = [295;770;440;430;315;400;370;450];
% Amount Per package table
TAMT = table(VitA,VitC,VitB1,VitB2,NA,CAL,...
'VariableNames',{'A','C','B1','B2','NA','CAL'},...
'RowNames',Pnames);
%% Objective
TotCost = sum(CostPerPack .* Packs); % Objective
obj = TotCost;
%% Constraints
prob = optimproblem('Objective',obj,'ObjectiveSense','min');
prob.Constraints.c1 = sum(VitA.*Packs) >= 700;
prob.Constraints.c1a = sum(VitA.*Packs) <= 20000;
prob.Constraints.c2 = sum(VitC.*Packs) >= 700;
prob.Constraints.c2a = sum(VitC.*Packs) <= 20000;
prob.Constraints.c3 = sum(VitB1.*Packs) >= 700;
prob.Constraints.c3a = sum(VitB1.*Packs) <= 20000;
prob.Constraints.c4 = sum(VitB2.*Packs) >= 700;
prob.Constraints.c4a = sum(VitB2.*Packs) <= 20000;
prob.Constraints.c5 = sum(NA.*Packs) >= 0;
prob.Constraints.c5a = sum(NA.*Packs) <= 40000;
prob.Constraints.c5 = sum(CAL.*Packs) >= 16000;
prob.Constraints.c5a = sum(CAL.*Packs) <= 24000;
%% Putting the problem together and solving
problem = prob2struct(prob);
% LP
% [sol,fval,exitflag,output1] = linprog(problem);
% [sol2,fval2,exitflag2,output2] = cplexlp(problem);
% MILP problem
[sol,fval,exitflag,output1] = intlinprog(problem);
[sol2,fval2,exitflag2,output2] =cplexmilp(problem);
%% Display
disp('Linear Prog function results')
if (~isempty(sol) )
T1 = table(sol,sol.*VitA,sol.*VitC,sol.*VitB1,sol.*VitB2,sol.*NA,sol.*CAL,sol.*CostPerPack,...
'VariableNames',{'NuofPacks','PerVitA','PerVitC','PerVitB1','PerVitB2','NA','CAL','Cost'},...
'RowNames',Pnames);
sumrow = array2table(sum(T1.Variables),...
'VariableNames',...
{'NuofPacks','PerVitA','PerVitC','PerVitB1','PerVitB2','NA','CAL','Cost'},...
'RowNames',"Sum");
T1 = [T1;sumrow];
disp(T1)
else
disp(['No feasible Solution with exit flag = ' ,num2str(exitflag)])
end
%% Display Cplex
disp('Cplex function results')
if (~isempty(sol2) )
T2 = table(sol2,sol2.*VitA,sol2.*VitC,sol2.*VitB1,sol2.*VitB2,sol2.*NA,sol2.*CAL,sol2.*CostPerPack,...
'VariableNames',{'NuofPacks','PerVitA','PerVitC','PerVitB1','PerVitB2','NA','CAL','Cost'},...
'RowNames',Pnames);
sumrow2 = array2table(sum(T2.Variables),...
'VariableNames',...
{'NuofPacks','PerVitA','PerVitC','PerVitB1','PerVitB2','NA','CAL','Cost'},...
'RowNames',"Sum");
T2 = [T2;sumrow2];
disp(T2)
else
disp(['No feasible Solution with exit flag = ' ,num2str(exitflag2)])
end
and here Cplex function results
NuofPacks PerVitA PerVitC PerVitB1 PerVitB2 NA CAL Cost
_________ _______ _______ ________ ________ _____ ______ ______
BEEF 2 120 40 20 30 1876 590 6.38
CHK 10 80 0 200 200 21800 7700 25.9
FISH 2 16 20 30 20 1890 880 4.58
HAM 2 80 80 70 20 556 860 5.78
MCH 10 150 350 150 150 11820 3150 18.9
MTL 10 700 300 150 150 8960 4000 19.9
SPG 7.3333 183.33 366.67 183.33 110 9746 2713.3 14.593
TUR 2 120 40 30 20 2794 900 4.98
Sum 45.333 1449.3 1196.7 833.33 700 59442 20793 101.01
The reults for this MILP from intlinprog seems good but from cplexmilp I get non-ineger value for SPG package. Can anyone help me to know the problem here ? Thanks in advance
The problem is the following: when you call problem = prob2struct(prob); the integrality information from prob is carried over into problem.intcon, but unfortunately it is not carried over into cplex.
This issue came up before, but was not addressed directly in cplex since there is a trivial workaround: replace the call
[sol2,fval2,exitflag2,output2] =cplexmilp(problem);
with
[sol2,fval2,exitflag2,output2] = ...
cplexmilp(problem.f, problem.Aineq, problem.bineq, problem.Aeq, ...
problem.beq, [], [], [], problem.lb, problem.ub, ...
intcon_ctype(length(problem.lb),problem.intcon));
where the intcon_ctype function is:
function ctype = intcon_ctype(numvars, intcon)
ctype = char(ones(1, numvars) * 'C');
for i = 1:length(intcon)
ctype(intcon(i)) = 'I'
end
end
Note that this is a very simplistic intcon_ctype conversion function. You may want to check the lb/ub values and set the type to 'B' if the variable is binary, you can specify semi-continuous variables, etc.