Question
With this code I'm trying to solve revised simplex method and 2 phase method on mat lab. However it is giving me errors. Please help
With this code I'm trying to solve revised simplex method and 2 phase method on mat lab. However it is giving me errors. Please help me fix it.
%SIMPLE 2p method
function simplex2p(type, c, A, rel, b)
% The Two-Phase Method for solving the LP problem
% min(or max) z = c*x
% Subject to Ax rel b
% x >= 0
% The input parameter type holds information about type of the LP
% problem to be solved. For the minimization problem type = 'min' and
% for the maximization problem type = 'max'.
% The input parameter rel is a string holding the relation signs.
% For instance, if rel = '<=>', then the constraint system consists
% of one inequality <=, one equation, and one inequality >=.
clc
if (type == 'min')
mm = 0;
else
mm = 1;
c = -c;
end
b = b(:);
c = c(:)';
[m, n] = size(A);
n1 = n;
les = 0;
neq = 0;
red = 0;
if length(c) < n
c = [c zeros(1,n-length(c))];
end
for i=1:m
if(rel(i) == '<')
A = [A vr(m,i)];
les = les + 1;
elseif(rel(i) == '>')
A = [A -vr(m,i)];
else
neq = neq + 1;
end
end
ncol = length(A);
if les == m
c = [c zeros(1,ncol-length(c))];
A = [A;c];
A = [A [b;0]];
[subs, A, z, p1] = loop(A, n1+1:ncol, mm, 1, 1);
disp(' End of Phase 1')
disp(' *********************************')
else
A = [A eye(m) b];
if m > 1
w = -sum(A(1:m,1:ncol));
else
w = -A(1,1:ncol);
end
14
c = [c zeros(1,length(A)-length(c))];
A = [A;c];
A = [A;[w zeros(1,m) -sum(b)]];
subs = ncol+1:ncol+m;
av = subs;
[subs, A, z, p1] = loop(A, subs, mm, 2, 1);
if p1 == 'y'
disp(' End of Phase 1')
disp(' *********************************')
end
nc = ncol + m + 1;
x = zeros(nc-1,1);
x(subs) = A(1:m,nc);
xa = x(av);
com = intersect(subs,av);
if (any(xa) ~= 0)
disp(sprintf(' Empty feasible region '))
return
else
if ~isempty(com)
red = 1;
end
end
A = A(1:m+1,1:nc);
A =[A(1:m+1,1:ncol) A(1:m+1,nc)];
[subs, A, z, p1] = loop(A, subs, mm, 1, 2);
if p1 == 'y'
disp(' End of Phase 2')
disp(' *********************************')
end
end
if (z == inf | z == -inf)
return
end
[m, n] = size(A);
x = zeros(n,1);
x(subs) = A(1:m-1,n);
x = x(1:n1);
if mm == 0
z = -A(m,n);
else
z = A(m,n);
end
disp(sprintf(' Problem has a finite optimal solution '))
disp(sprintf(' Values of the legitimate variables: '))
for i=1:n1
disp(sprintf(' x(%d)= %f ',i,x(i)))
end
disp(sprintf(' Objective value at the optimal point: '))
disp(sprintf(' z= %f',z))
t = find(A(m,1:n-1) == 0);
if length(t) > m-1
str = 'Problem has infinitely many solutions';
msgbox(str,'Warning Window','warn')
end
if red == 1
disp(sprintf(' Constraint system is redundant '))
15
end
function [subs, A, z, p1]= loop(A, subs, mm, k, ph)
% Main loop of the simplex primal algorithm.
% Bland's rule to prevente cycling is used.
tbn = 0;
str1 = 'Would you like to monitor the progress of Phase 1?';
str2 = 'Would you like to monitor the progress of Phase 2?';
if ph == 1
str = str1;
else
str = str2;
end
question_ans = questdlg(str,'Make a choice Window','Yes','No','No');
if strcmp(question_ans,'Yes')
p1 = 'y';
end
if p1 == 'y' & ph == 1
disp(sprintf(' Initial tableau'))
A
disp(sprintf(' Press any key to continue ... '))
pause
end
if p1 == 'y' & ph == 2
tbn = 1;
disp(sprintf(' Tableau %g',tbn))
A
disp(sprintf(' Press any key to continue ... '))
pause
end
[m, n] = size(A);
[mi, col] = Br(A(m,1:n-1));
while ~isempty(mi) & mi < 0 & abs(mi) > eps
t = A(1:m-k,col);
if all(t <= 0)
if mm == 0
z = -inf;
else
z = inf;
end
disp(sprintf(' Unbounded optimal solution with z= ',z))
return
end
[row, small] = MRT(A(1:m-k,n),A(1:m-k,col));
if ~isempty(row)
if abs(small) <= 100*eps & k == 1
[s,col] = Br(A(m,1:n-1));
end
if p1 == 'y'
fprintf(' pivot row-> %g pivot column->)%g ',...row,col; end
16
A(row,:)= A(row,:)/A(row,col);
subs(row) = col;
for i = 1:m
if i ~= row
A(i,:)= A(i,:)-A(i,col)*A(row,:);
end
end
[mi, col] = Br(A(m,1:n-1));
end
tbn = tbn + 1;
if p1 == 'y'
disp(sprintf(' Tableau %g',tbn))
A
disp(sprintf(' Press any key to continue ... '))
pause
end
end
z = A(m,n);
MATLAB subfunction loop is included in the m-file simplex2p. It implements the main loop of
the simplex algorithm.
SIMPLEX METHOD
% Create the initial Basis matrix, compute its inverse and
B=A(:,d);
invB = inv(B);
x(d) = invB*b;
%Iteration phase
n_max = 10; % At most n_max iterations
for n = 1:n_max
% Compute the vector of reduced costs c_r
c_B = c(d); % Basic variable costs
p = (c_B'*invB)'; % Dual variables
c_r = c' - p'*A; % Vector of reduced costs
% Check if the solution is optimal. If optimal, use
% 'return' to break from the function, e.g.
J = find(c_r < 0); % Find indices with negative reduced costs
if (isempty(J))
f_opt = c'*x;
x_opt = x;
return;
end
% Choose the entering variable
j_in = J(1);
% Compute the vector u (i.e., the reverse of the basic directions)
u = invB*A(:,j_in);
I = find(u > 0);
if (isempty(I))
f_opt = -inf; % Optimal objective function cost = -inf
x_opt = []; % Produce empty vector []
return % Break from the function
end
% Compute the optimal step length delta
delta = min(x(C(I))./u(I));
L = find(x(C)./u == delta); % Find all indices with ratio delta
% Select the exiting variable
l_out = L(1);
% Move to the adjacent solution
x(d) = x(d) - delta*u;
% Value of the incoming changeable is delta
x(j_in) = delta;
% Update the set of basic indices d
d(l_out) = j_in;
% calculate the original contrary basis B^-1 by the stage plain row
% operation on [B^-1 u] (pivot row index is l_out). The vector u is trans-
% shaped into a unit vector with u(l_out) = 1 and u(i) = 0 for
% other i.
M=horzcat(invB,u);
[f g]=size(M);
R(l_out,:)=M(l_out,:)/M(l_out,j_in); % Copy row l_out, normalizing M(l_out,j_in) to 1
u(l_out)=1;
for k = 1:f % For all matrix rows
if (k ~= l_out) % Other then l_out
u(k)=0;
R(k,:)=M(k,:)-M(k,j_in)*R(l_out,:); % Set them equal to the original matrix Minus a multiple of normalized row l_out, making R(k,j_in)=0
end
end
invM=horzcat(u,invB);
% make sure if too several iterations are performed (increase n_max to
% allow more iterations)
if(n == n_max)
fprintf('Max number of iterations performed! ');
return
end
end % End for (the main iteration loop)
end % End function
Step by Step Solution
There are 3 Steps involved in it
Step: 1
Get Instant Access to Expert-Tailored Solutions
See step-by-step solutions with expert insights and AI powered tools for academic success
Step: 2
Step: 3
Ace Your Homework with AI
Get the answers you need in no time with our AI-driven, step-by-step assistance
Get Started