Answered step by step
Verified Expert Solution
Link Copied!

Question

1 Approved Answer

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

blur-text-image

Get Instant Access to Expert-Tailored Solutions

See step-by-step solutions with expert insights and AI powered tools for academic success

Step: 2

blur-text-image

Step: 3

blur-text-image

Ace Your Homework with AI

Get the answers you need in no time with our AI-driven, step-by-step assistance

Get Started

Recommended Textbook for

Graph Databases

Authors: Ian Robinson, Jim Webber, Emil Eifrem

1st Edition

1449356265, 978-1449356262

More Books

Students also viewed these Databases questions

Question

Does the series converge as N ? (3.90)

Answered: 1 week ago