Skip to content

Commit d837ef4

Browse files
ploskasnikosploskasnikos
authored
Add files via upload
1 parent d8da06f commit d837ef4

Some content is hidden

Large Commits have some content hidden by default. Use the searchbox below for content that may be hidden.

44 files changed

+5038
-0
lines changed

‎codes/appendix A/example1a.m

+34
Original file line numberDiff line numberDiff line change
@@ -0,0 +1,34 @@
1+
variables = {'x1', 'x2', 'x3', 'x4'}; % construct a cell
2+
% array with the names of the variables
3+
n = length(variables); % find the number of the variables
4+
for i = 1:n % create x for indexing
5+
eval([variables{i}, ' = ', num2str(i), ';']);
6+
end
7+
lb = zeros(1, n); % create a vector and add a lower bound
8+
% of zero to all variables
9+
lb([x1, x3]) = [1, 2]; % add lower bounds to variables x1,
10+
% x3
11+
ub = Inf(1, n); % create a vector and add an upper bound of
12+
% Inf to all variables
13+
ub(x3) = 10; % add an upper bound to variable x3
14+
A = zeros(2, n); % create the matrix A
15+
b = zeros(2, 1); % create the vector b
16+
% declare the 1st constraint
17+
A(1, [x1, x3]) = [-2, -3]; b(1) = -6;
18+
% declare the 2nd constraint
19+
A(2, [x1, x2, x4]) = [-3, 2, -4]; b(2) = -8;
20+
Aeq = zeros(2, n); % create the matrix Aeq
21+
beq = zeros(2, 1); % create the vector beq
22+
% declare the 1st equality constraint
23+
Aeq(1, [x1, x2, x3, x4]) = [4, -3, 8, -1]; beq(1) = 20;
24+
% declare the 2nd equality constraint
25+
Aeq(2, [x1, x3, x4]) = [4, -1, 4]; beq(2) = 18;
26+
c = zeros(n, 1); % create the objective function vector c
27+
c([x1 x2 x3 x4]) = [-2; 4; -2; 2];
28+
% call the linprog solver
29+
[x, objVal] = linprog(c, A, b, Aeq, beq, lb, ub);
30+
for i = 1:n % print results in formatted form
31+
fprintf('%s \t %20.4f\n', variables{i}, x(i))
32+
end
33+
fprintf(['The value of the objective function ' ...
34+
'is %.4f\n'], objVal)

‎codes/appendix A/example1b.m

+14
Original file line numberDiff line numberDiff line change
@@ -0,0 +1,14 @@
1+
A = [-2 0 -3 0; -3 2 0 -4]; % create the matrix A
2+
b = [-6; -8]; % create the vector b
3+
Aeq = [4 -3 8 -1; 4 0 -1 4]; %create the matrix Aeq
4+
beq = [20; 18]; % create the vector beq
5+
lb = [1 0 2 0]; % create the vector lb
6+
ub = [Inf Inf 10 Inf]; % create the vector ub
7+
c = [-2; 4; -2; 2]; %create the vector c
8+
% call the linprog solver
9+
[x, objVal] = linprog(c, A, b, Aeq, beq, lb, ub);
10+
for i = 1:4 % print results in formatted form
11+
fprintf('x%d \t %20.4f\n', i, x(i))
12+
end
13+
fprintf(['The value of the objective function ' ...
14+
'is %.4f\n'], objVal)

‎codes/appendix A/example2a.m

+33
Original file line numberDiff line numberDiff line change
@@ -0,0 +1,33 @@
1+
variables = {'x1', 'x2', 'x3', 'x4', 'x5'}; % construct a
2+
% cell array with the names of the variables
3+
n = length(variables); % find the number of the variables
4+
for i = 1:n % create x for indexing
5+
eval([variables{i}, ' = ', num2str(i), ';']);
6+
end
7+
lb = zeros(1, n); % create a vector and add a lower bound
8+
% of zero to all variables
9+
lb([x1, x2]) = [5, 1]; % add lower bounds to variables x1,
10+
% x2
11+
ub = Inf(1, n); % create a vector and add an upper bound
12+
% of Inf to all variables
13+
A = zeros(3, n); % create the matrix A
14+
b = zeros(3, 1); % create the vector b
15+
% declare the 1st inequality constraint
16+
A(1, [x1, x2, x3, x4]) = [1, 4, -2, -1]; b(1) = 8;
17+
% declare the 2nd inequality constraint
18+
A(2, [x1, x2, x3, x4]) = [1, 3, 2, -1]; b(2) = 10;
19+
% declare the 3rd inequality constraint
20+
A(3, [x1, x2, x3, x4, x5]) = [2, 1, 2, 3, -1]; b(3) = 20;
21+
Aeq = zeros(1, n); % create the matrix Aeq
22+
beq = zeros(1, 1); % create the vector beq
23+
% declare the equality constraint
24+
Aeq(1, [x1, x2, x3, x4, x5]) = [1, 3, -4, -1, 1]; beq(1) = 7;
25+
c = zeros(n, 1); % create the objective function vector c
26+
c([x1 x2 x3 x4 x5]) = [-2; -3; 1; 4; 1];
27+
% call the linprog solver
28+
[x, objVal] = linprog(c, A, b, Aeq, beq, lb, ub);
29+
for i = 1:n % print results in formatted form
30+
fprintf('%s \t %20.4f\n', variables{i}, x(i))
31+
end
32+
fprintf(['The value of the objective function ' ...
33+
'is %.4f\n'], objVal)

‎codes/appendix A/example2b.m

+14
Original file line numberDiff line numberDiff line change
@@ -0,0 +1,14 @@
1+
A = [1 4 -2 -1 0; 1 3 2 -1 0; 2 1 2 3 -1]; % create the matrix A
2+
b = [8; 10; 20]; % create the vector b
3+
Aeq = [1 3 -4 -1 1]; %create the matrix Aeq
4+
beq = [7]; % create the vector beq
5+
lb = [5 1 0 0 0]; % create the vector lb
6+
ub = [Inf Inf Inf Inf Inf]; % create the vector ub
7+
c = [-2; -3; 1; 4; 1]; %create the vector c
8+
% call the linprog solver
9+
[x, objVal] = linprog(c, A, b, Aeq, beq, lb, ub);
10+
for i = 1:5 % print results in formatted form
11+
fprintf('x%d \t %20.4f\n', i, x(i))
12+
end
13+
fprintf(['The value of the objective function ' ...
14+
'is %.4f\n'], objVal)

‎codes/appendix A/linprogSolver.m

+137
Original file line numberDiff line numberDiff line change
@@ -0,0 +1,137 @@
1+
function [xsol, fval, exitflag, iterations] = ...
2+
linprogSolver(A, c, b, Eqin, MinMaxLP, c0, ...
3+
algorithm)
4+
% Filename: linprogSolver.m
5+
% Description: the function is a MATLAB code to solve LPs
6+
% using the linprog solver
7+
% Authors: Ploskas, N., & Samaras, N.
8+
%
9+
% Syntax: [x, fval, exitflag, iterations] = ...
10+
% linprogSolver(A, c, b, Eqin, MinMaxLP, c0, ...
11+
% algorithm)
12+
%
13+
% Input:
14+
% -- A: matrix of coefficients of the constraints
15+
% (size m x n)
16+
% -- c: vector of coefficients of the objective function
17+
% (size n x 1)
18+
% -- b: vector of the right-hand side of the constraints
19+
% (size m x 1)
20+
% -- Eqin: vector of the type of the constraints
21+
% (size m x 1)
22+
% -- MinMaxLP: the type of optimization (optional:
23+
% default value -1 - minimization)
24+
% -- c0: constant term of the objective function
25+
% (optional: default value 0)
26+
% -- algorithm: the LP algorithm that will be used
27+
% (possible values: 'interior-point', 'dual-simplex',
28+
% 'simplex', 'active-set')
29+
%
30+
% Output:
31+
% -- xsol: the solution found by the solver (size m x 1)
32+
% -- fval: the value of the objective function at the
33+
% solution xsol
34+
% -- exitflag: the reason that the algorithm terminated
35+
% (1: the solver converged to a solution x, 0: the
36+
% number of iterations exceeded the MaxIter option,
37+
% -1: the input data is not logically or numerically
38+
% correct, -2: no feasible point was found, -3: the LP
39+
% problem is unbounded, -4: NaN value was encountered
40+
% during the execution of the algorithm, -5: both primal
41+
% and dual problems are infeasible, -7: the search
42+
% direction became too small and no further progress
43+
% could be made)
44+
% -- iterations: the number of iterations
45+
46+
% set default values to missing inputs
47+
if ~exist('MinMaxLP')
48+
MinMaxLP = -1;
49+
end
50+
if ~exist('c0')
51+
c0 = 0;
52+
end
53+
[m, n] = size(A); % find the size of matrix A
54+
[m2, n2] = size(c); % find the size of vector c
55+
[m3, n3] = size(Eqin); % find the size of vector Eqin
56+
[m4, n4] = size(b); % find the size of vector b
57+
% check if input data is logically correct
58+
if n2 ~= 1
59+
disp('Vector c is not a column vector.')
60+
exitflag = -1;
61+
return
62+
end
63+
if n ~= m2
64+
disp(['The number of columns in matrix A and ' ...
65+
'the number of rows in vector c do not match.'])
66+
exitflag = -1;
67+
return
68+
end
69+
if m4 ~= m
70+
disp(['The number of the right-hand side values ' ...
71+
'is not equal to the number of constraints.'])
72+
exitflag = -1;
73+
return
74+
end
75+
if n3 ~= 1
76+
disp('Vector Eqin is not a column vector')
77+
exitflag = -1;
78+
return
79+
end
80+
if n4 ~= 1
81+
disp('Vector b is not a column vector')
82+
exitflag = -1;
83+
return
84+
end
85+
if m4 ~= m3
86+
disp('The size of vectors Eqin and b does not match')
87+
exitflag = -1;
88+
return
89+
end
90+
% if the problem is a maximization problem, transform it to
91+
% a minimization problem
92+
if MinMaxLP == 1
93+
c = -c;
94+
end
95+
Aeq = []; % matrix constraint of equality constraints
96+
beq = []; % right-hand side of equality constraints
97+
% check if all constraints are equalities
98+
flag = isequal(zeros(m3, 1), Eqin);
99+
if flag == 0 % some or all constraints are inequalities
100+
indices = []; % indices of the equality constraints
101+
for i = 1:m3
102+
if Eqin(i, 1) == 0 % find the equality constraints
103+
indices = [indices i];
104+
% convert 'greater than or equal to' type of
105+
% constraints to 'less than or equal to' type
106+
% of constraints
107+
elseif Eqin(i, 1) == 1
108+
A(i, :) = -A(i, :);
109+
b(i) = -b(i);
110+
end
111+
end
112+
% create the matrix constraint of equality constraints
113+
Aeq = A(indices, :);
114+
A(indices, :) = []; % delete equality constraints
115+
% create the right-hand side of equality constraints
116+
beq = b(indices);
117+
b(indices) = []; % delete equality constraints
118+
else % all constraints are equalities
119+
Aeq = A;
120+
beq = b;
121+
A = [];
122+
b = [];
123+
end
124+
lb = zeros(n, 1); % create zero lower bounds
125+
% choose LP algorithm
126+
options = optimoptions(@linprog, 'Algorithm', algorithm);
127+
% call the linprog solver
128+
[xsol, fval, exitflag, output] = linprog(c, A, b, ...
129+
Aeq, beq, lb, [], [], options);
130+
iterations = output.iterations;
131+
% calculate the value of the objective function
132+
if MinMaxLP == 1 % maximization
133+
fval = -(fval + c0);
134+
else % minimization
135+
fval = fval + c0;
136+
end
137+
end

‎codes/appendix B/clpOptimizer.m

+111
Original file line numberDiff line numberDiff line change
@@ -0,0 +1,111 @@
1+
function [xsol, fval, exitflag, iterations] = ...
2+
clpOptimizer(mpsFullPath, algorithm, primalTol, ...
3+
dualTol, maxIter, maxTime, displayLevel, objBias, ...
4+
numPresolvePasses, factorFreq, numberRefinements, ...
5+
primalObjLim, dualObjLim, numThreads, abcState)
6+
% Filename: clpOptimizer.m
7+
% Description: the function is a MATLAB code to solve LPs
8+
% using CLP (via OPTI Toolbox)
9+
% Authors: Ploskas, N., & Samaras, N.
10+
%
11+
% Syntax: [xsol, fval, exitflag, iterations] = ...
12+
% clpOptimizer(mpsFullPath, algorithm, primalTol, ...
13+
% dualTol, maxIter, maxTime, display, objBias, ...
14+
% numPresolvePasses, factorFreq, numberRefinements, ...
15+
% primalObjLim, dualObjLim, numThreads, abcState)
16+
%
17+
% Input:
18+
% -- mpsFullPath: the full path of the MPS file
19+
% -- algorithm: the CLP algorithm that will be used
20+
% (optional, possible values: 'DualSimplex',
21+
% 'PrimalSimplex', 'PrimalSimplexOrSprint',
22+
% 'Barrier', 'BarrierNoCross', 'Automatic')
23+
% -- primalTol: the primal tolerance (optional)
24+
% -- dualTol: the dual tolerance (optional)
25+
% -- maxIter: the maximum number of iterations
26+
% -- (optional)
27+
% -- maxTime: the maximum execution time in seconds
28+
% -- (optional)
29+
% -- displayLevel: the display level (optional,
30+
% possible values: 0, 1 -> 100 increasing)
31+
% -- objBias: the objective bias term (optional)
32+
% -- numPresolvePasses: the number of presolver passes
33+
% (optional)
34+
% -- factorFreq: every factorFreq number of iterations,
35+
% the basis inverse is re-computed from scratch
36+
% (optional)
37+
% -- numberRefinements: the number of iterative simplex
38+
% refinements (optional)
39+
% -- primalObjLim: the primal objective limit (optional)
40+
% -- dualObjLim: the dual objective limit (optional)
41+
% -- numThreads: the number of Cilk worker threads
42+
% (optional, only with Aboca CLP build)
43+
% -- abcState: Aboca's partition size (optional)
44+
%
45+
% Output:
46+
% -- xsol: the solution found by the solver (size m x 1)
47+
% -- fval: the value of the objective function at the
48+
% solution xsol
49+
% -- exitflag: the reason that the algorithm terminated
50+
% (1: the solver converged to a solution x, 0: the
51+
% number of iterations exceeded the maxIter option or
52+
% time reached, -1: the LP problem is infeasible or
53+
% unbounded)
54+
% -- iterations: the number of iterations
55+
56+
% set user defined values to options
57+
opts = optiset('solver', 'clp');
58+
if exist('algorithm')
59+
opts.solverOpts.algorithm = algorithm;
60+
end
61+
if exist('primalTol')
62+
opts.solverOpts.primalTol = primalTol;
63+
opts.tolrfun = primalTol;
64+
end
65+
if exist('dualTol')
66+
opts.solverOpts.dualTol = dualTol;
67+
end
68+
if exist('maxIter')
69+
opts.maxiter = maxIter;
70+
else
71+
opts.maxiter = 1000000;
72+
end
73+
if exist('maxTime')
74+
opts.maxtime = maxTime;
75+
else
76+
opts.maxtime = 1000000;
77+
end
78+
if exist('displayLevel')
79+
opts.display = displayLevel;
80+
end
81+
if exist('objBias')
82+
opts.solverOpts.objbias = objBias;
83+
end
84+
if exist('numPresolvePasses')
85+
opts.solverOpts.numPresolvePasses = numPresolvePasses;
86+
end
87+
if exist('factorFreq')
88+
opts.solverOpts.factorFreq = factorFreq;
89+
end
90+
if exist('numberRefinements')
91+
opts.solverOpts.numberRefinements = numberRefinements;
92+
end
93+
if exist('primalObjLim')
94+
opts.solverOpts.primalObjLim = primalObjLim;
95+
end
96+
if exist('dualObjLim')
97+
opts.solverOpts.dualObjLim = dualObjLim;
98+
end
99+
if exist('numThreads')
100+
opts.solverOpts.numThreads = numThreads;
101+
end
102+
if exist('abcState')
103+
opts.solverOpts.abcState = abcState;
104+
end
105+
% read the MPS file
106+
prob = coinRead(mpsFullPath);
107+
% call CLP solver
108+
[xsol, fval, exitflag, info] = opti_clp([], prob.f, ...
109+
prob.A, prob.rl, prob.ru, prob.lb, prob.ub, opts);
110+
iterations = info.Iterations;
111+
end

0 commit comments

Comments
 (0)