Skip to content

Commit

Permalink
Merge pull request opencobra#967 from shjchan/fix_solveCobraMILP_solv…
Browse files Browse the repository at this point in the history
…erParam

Fix Gurobi-specific parameter structure in solveCobraMILP.m
  • Loading branch information
syarra authored Oct 29, 2017
2 parents 71cf13c + 2c04ef2 commit 6aeef2a
Show file tree
Hide file tree
Showing 3 changed files with 113 additions and 16 deletions.
Original file line number Diff line number Diff line change
Expand Up @@ -48,7 +48,7 @@
P.osense = 1;
P.csense = [repmat('E', 2*size(Np,2) + length(c), 1); 'L'; 'G'; 'G'];
P.vartype = repmat('I', size(P.A,2), 1);
P.x0 = suma*[c; c];
P.x0 = suma' .* [c; c];

% Run MILP
solution = solveCobraMILP(P);
Expand Down
48 changes: 34 additions & 14 deletions src/base/solvers/solveCobraMILP.m
Original file line number Diff line number Diff line change
Expand Up @@ -88,7 +88,7 @@
'feasTol', 'optTol', 'absMipGapTol', 'NUMERICALEMPHASIS', 'solver'};

parameters = [];

parametersStructureFlag = false;
% First input can be 'default' or a solver-specific parameter structure
if ~isempty(varargin)
isdone = false(size(varargin));
Expand All @@ -99,8 +99,8 @@
varargin = varargin(~isdone);

elseif isstruct(varargin{1}) % solver-specific parameter structure
solverParams = varargin{1};

[solverParams, directParamStruct] = deal(varargin{1});
parametersStructureFlag = true;
isdone(1) = true;
varargin = varargin(~isdone);
end
Expand All @@ -111,8 +111,8 @@
isdone = false(size(varargin));

if isstruct(varargin{end})
solverParams = varargin{end};

[solverParams, directParamStruct] = deal(varargin{end});
parametersStructureFlag = true;
isdone(end) = true;
varargin = varargin(~isdone);
end
Expand All @@ -122,7 +122,6 @@
if mod(length(varargin), 2) == 0
try
parameters = struct(varargin{:});
parametersStructureFlag = 0;
catch
error('solveCobraLP: Invalid parameter name-value pairs.')
end
Expand Down Expand Up @@ -156,12 +155,16 @@
else
error('solveCobraMILP: Invalid number of parameters/values')
end
% [minNorm, printLevel, primalOnlyFlag, saveInput, feasTol, optTol] = ...
% getCobraSolverParams('LP', optParamNames(1:6), parameters);
[minNorm, printLevel, primalOnlyFlag, saveInput, feasTol, optTol] = ...
getCobraSolverParams('LP', optParamNames(1:6), parameters);
getCobraSolverParams('LP', {'minNorm', 'printLevel', 'primalOnlyFlag', 'saveInput', 'feasTol', 'optTol'}, parameters);
else
parametersStructureFlag = 0;
% [minNorm, printLevel, primalOnlyFlag, saveInput, feasTol, optTol] = ...
% getCobraSolverParams('LP', optParamNames(1:6));
[minNorm, printLevel, primalOnlyFlag, saveInput, feasTol, optTol] = ...
getCobraSolverParams('LP', optParamNames(1:6));
getCobraSolverParams('LP', {'minNorm', 'printLevel', 'primalOnlyFlag', 'saveInput', 'feasTol', 'optTol'}, parameters);
% parameters will later be accessed and should be initialized.
parameters = '';
end
Expand Down Expand Up @@ -190,6 +193,10 @@
xCont = [];
f = [];

if ~isfield(MILPproblem, 'x0')
MILPproblem.x0 = [];
end

[A, b, c, lb, ub, csense, osense, vartype, x0] = ...
deal(MILPproblem.A, MILPproblem.b, MILPproblem.c, MILPproblem.lb, MILPproblem.ub, ...
MILPproblem.csense, MILPproblem.osense, MILPproblem.vartype, MILPproblem.x0);
Expand Down Expand Up @@ -316,8 +323,8 @@
end

% minimum intTol for gurobi = 1e-9
if solverParams.intTol<1e-9,
solverParams.intTol=1e-9
if solverParams.intTol<1e-9
solverParams.intTol = 1e-9;
end

opts.TimeLimit=solverParams.timeLimit;
Expand Down Expand Up @@ -387,8 +394,16 @@
cplexlp.Param.mip.tolerances.integrality.Cur = solverParams.intTol;
cplexlp.Param.timelimit.Cur = solverParams.timeLimit;
cplexlp.Param.output.writelevel.Cur = solverParams.printLevel;

outputfile = fopen(solverParams.logFile,'a');


if isscalar(solverParams.logFile) && solverParams.logFile == 1
% allow print to command window by setting solverParams.logFile == 1
outputfile = 1;
logToConsole = true;
else
outputfile = fopen(solverParams.logFile,'a');
logToConsole = false;
end
cplexlp.DisplayFunc = @redirect;

cplexlp.Param.simplex.tolerances.optimality.Cur = solverParams.optTol;
Expand All @@ -411,8 +426,10 @@
% Solve problem
Result = cplexlp.solve();

% Close the output file
fclose(outputfile);
if ~logToConsole
% Close the output file
fclose(outputfile);
end

% Get results
x = Result.x;
Expand Down Expand Up @@ -494,6 +511,9 @@
MILPproblem.vtype = vartype;
MILPproblem.modelsense = MILPproblem.osense;
[MILPproblem.A,MILPproblem.rhs,MILPproblem.obj,MILPproblem.sense] = deal(sparse(MILPproblem.A),MILPproblem.b,double(MILPproblem.c),MILPproblem.csense);
if ~isempty(x0)
MILPproblem.start = x0;
end
resultgurobi = gurobi(MILPproblem,params);

stat = resultgurobi.status;
Expand Down
79 changes: 78 additions & 1 deletion test/verifiedTests/base/testSolvers/testSolveCobraMILP.m
Original file line number Diff line number Diff line change
Expand Up @@ -17,14 +17,18 @@
fileDir = fileparts(which('testSolveCobraMILP'));
cd(fileDir);

% save original solve
global CBT_MILP_SOLVER;
orig_solver = CBT_MILP_SOLVER;

% test solver packages
solverPkgs = {'cplex_direct', 'ibm_cplex', 'tomlab_cplex', 'gurobi6', 'glpk'};

% set the tolerance
tol = 1e-8;

for k = 1:length(solverPkgs)
fprintf(' Running solveCobraLPCPLEX using %s ... ', solverPkgs{k});
fprintf(' Running solveCobraMILP using %s ... ', solverPkgs{k});

% change the COBRA solver (LP)
solverOK = changeCobraSolver(solverPkgs{k}, 'MILP', 0);
Expand All @@ -51,13 +55,18 @@
parameters.relMipGapTol = 1e-12;
parameters.intTol = 1e-12;
MILPsolution = solveCobraMILP(MILPproblem, parameters);
% check if MILP can be solved without x0 supplied
MILPsolution2 = solveCobraMILP(rmfield(MILPproblem, 'x0'), parameters);
else
MILPsolution = solveCobraMILP(MILPproblem);
% check if MILP can be solved without x0 supplied
MILPsolution2 = solveCobraMILP(rmfield(MILPproblem, 'x0'));
end

% check results with expected answer.
assert(all(abs(MILPsolution.int - [0; 31; 46]) < tol))
assert(abs(MILPsolution.obj - 554) < tol)
assert(abs(MILPsolution2.obj - 554) < tol)

if strcmp(solverPkgs{k}, 'ibm_cplex')
% test IBM-Cplex-specific parameters. Solve with the below parameters changed
Expand Down Expand Up @@ -102,17 +111,85 @@
delete(['testIBMcplexMILPparam' num2str(jTest) '.log']);
end
end

if strcmp(solverPkgs{k}, 'gurobi6')
% check additional parameters for Gurobi
% temporarily shut down warning
warning_stat = warning;
warning off
% test TimeLimit as a gurobi-specific parameter
sol = solveCobraMILP(MILPproblem, struct('TimeLimit',0));
assert(strcmp(sol.origStat, 'TIME_LIMIT'))
% restore previous warning state
warning(warning_stat)

% check user-supplied x0
MILPproblem.A = rand(10,20);
MILPproblem.b = 1000 * ones(10,1);
MILPproblem.c = zeros(20,1);
MILPproblem.lb = zeros(20,1);
MILPproblem.ub = ones(20,1);
MILPproblem.vartype = char(['C'*ones(1,10), 'B'*ones(1,10)]);
MILPproblem.csense = char('L' * ones(1, 10));

% no objective function. The supplied should be the returned
% (if not everything becomes zero after presolve)
MILPproblem.x0 = zeros(20, 1);
sol = solveCobraMILP(MILPproblem);
assert(isequal(sol.full, MILPproblem.x0));

MILPproblem.x0 = ones(20, 1);
sol = solveCobraMILP(MILPproblem);
assert(isequal(sol.full, MILPproblem.x0));

end
end

% output a success message
fprintf('Done.\n');
end

% test ibm_cplex output to command window
solverOK = changeCobraSolver('ibm_cplex', 'MILP', 0);
if solverOK
fprintf('Test ibm_cplex output to command window ...\n')
% solve without logToFile = 1
diary test_ibm_cplex_output_to_console1.txt
sol = solveCobraMILP(MILPproblem);
diary off
% read the diary, which should be empty
f = fopen('test_ibm_cplex_output_to_console1.txt', 'r');
l = fgets(f);
assert(isequal(l, -1))
fclose(f);
delete('test_ibm_cplex_output_to_console1.txt')

% solve wit logToFile = 1
diary test_ibm_cplex_output_to_console2.txt
sol = solveCobraMILP(MILPproblem, 'logFile', 1);
diary off
% read the diary, which should be non-empty
f = fopen('test_ibm_cplex_output_to_console2.txt', 'r');
l = fgets(f);
line = 0;
while ~isequal(l, -1)
line = line + 1;
l = fgets(f);
end
fclose(f);
assert(line > 3)
delete('test_ibm_cplex_output_to_console2.txt')
fprintf('Test ibm_cplex output to command window ... Done\n')
end

% remove the generated file
fullFileNamePath = [fileparts(which(mfilename)), filesep, 'MILPProblem.mat'];
if exist(fullFileNamePath, 'file') == 2
delete(fullFileNamePath);
end

% change back to the original solver
changeCobraSolver(orig_solver, 'MILP', 0);

% change the directory
cd(currentDir)

0 comments on commit 6aeef2a

Please sign in to comment.