%% Math 464 Linear Optimization (Spring 2023) %% Lecture 27, 04/18/2023 %% Bland's rule (min-index) for pivoting to avoid cycling %% We go back to the instance from Lecture 18 that cycled when using %% default pivot rules. >> T = [ 0 -3/4 20 -1/2 6 0 0 0; 0 1/4 -8 -1 9 1 0 0; 0 1/2 -12 -1/2 3 0 1 0; 1 0 0 1 6 0 0 1] T = 0 -3/4 20 -1/2 6 0 0 0 0 1/4 -8 -1 9 1 0 0 0 1/2 -12 -1/2 3 0 1 0 1 0 0 1 6 0 0 1 >> [m,n] = size(T); >> Bind = [5 6 7]; %% We have to maintain and update the basis (as specified in the %% vectore Bind) - we need this info to select the leaving variable %% with the smallest indiex in each case >> % Iteration 1 >> j=2; %% This is the smallest j with c'_j < 0. It's just a coincidence that %% this c'_j is also the most negative in this case. %% We have to first find all candidates for the min-ratio test, then %% find all among them that give the min-ratio (theta). From among all %% these candidates, we need to pick the one with the smallest index %% of the corresponding x_l. >> j=2; inds = find(T(2:m,j)>0) inds = 1 2 %% We add 1 to the inds in order to make them apply directly to the %% tableau T (whose first rwo is actually Row-. >> j=2; inds = find(T(2:m,j)>0)+1 inds = 2 3 >> j=2; inds = find(T(2:m,j)>0)+1;ratios = T(inds,1)./T(inds,j) ratios = 0 0 %% Note that theta = 0 here, and that both the ratios are equal to %% this min value. But we need to do these calculations in the general %% setting. >> j=2; inds = find(T(2:m,j)>0)+1;ratios = T(inds,1)./T(inds,j);theta=min(ratios);candxjs = find(ratios==theta) candxjs = 1 2 >> j=2; inds = find(T(2:m,j)>0)+1;ratios = T(inds,1)./T(inds,j);theta=min(ratios);candxjs = Bind(find(ratios==theta)) candxjs = 5 6 %% That was a first attempt at indexing back to the x_j's. While it %% appears to work in this case, we actually need to do the %% back-indexing properly >> j=2;inds = find(T(2:m,j) > 0)+1;ratios = T(inds,1)./T(inds,j);theta=min(ratios);candxjs = Bind(inds(find(ratios==theta))-1) candxjs = 5 6 %% We actually need to a bit more - see commands below. Note how we %% store the min index among the candidate xj's in xl. >> j=2;inds = find(T(2:m,j) > 0)+1;ratios = T(inds,1)./T(inds,j);theta=min(ratios);candxjs = Bind(inds(find(ratios==theta))-1);[xl,l]=min(candxjs) xl = 5 l = 1 %% We need to do one last re-indexing for the argmin l of the min %% command to go back to the tableau indexing. >> j=2;inds = find(T(2:m,j) > 0)+1;ratios = T(inds,1)./T(inds,j);theta=min(ratios);candxjs = Bind(inds(find(ratios==theta))-1);[xl,l]=min(candxjs);l=find(Bind==xl)+1 l = 2 %% Let's look at T again... >> T T = 0 -3/4 20 -1/2 6 0 0 0 0 1/4 -8 -1 9 1 0 0 0 1/2 -12 -1/2 3 0 1 0 1 0 0 1 6 0 0 1 %% and do the pivot steps >> T(l,:)=T(l,:)/T(l,j); for i=setdiff(1:m,l);T(i,:) = T(i,:) - T(l,:)*T(i,j);end;T T = 0 0 -4 -7/2 33 3 0 0 0 1 -32 -4 36 4 0 0 0 0 4 3/2 -15 -2 1 0 1 0 0 1 6 0 0 1 %% Finally, we update Bind >> Bind Bind = 5 6 7 >> Bind(find(Bind==xl))=j-1 Bind = 1 6 7 >> % Iteratio 2 %% x2 enters because it's the leftmost x_j with c'_j < 0 >> j=3;inds = find(T(2:m,j) > 0)+1;ratios = T(inds,1)./T(inds,j);theta=min(ratios);candxjs = Bind(inds(find(ratios==theta))-1);[xl,l]=min(candxjs);l=find(Bind==xl)+1 l = 3 >> T(l,:)=T(l,:)/T(l,j); for i=setdiff(1:m,l);T(i,:) = T(i,:) - T(l,:)*T(i,j);end;T T = 0 0 0 -2 18 1 1 0 0 1 0 8 -84 -12 8 0 0 0 1 3/8 -15/4 -1/2 1/4 0 1 0 0 1 6 0 0 1 >> Bind(find(Bind==xl))=j-1 Bind = 1 2 7 >> % Iteration 3 >> j=4;inds = find(T(2:m,j) > 0)+1;ratios = T(inds,1)./T(inds,j);theta=min(ratios);candxjs = Bind(inds(find(ratios==theta))-1);[xl,l]=min(candxjs);l=find(Bind==xl)+1 l = 2 >> T(l,:)=T(l,:)/T(l,j); for i=setdiff(1:m,l);T(i,:) = T(i,:) - T(l,:)*T(i,j);end;T T = 0 1/4 0 0 -3 -2 3 0 0 1/8 0 1 -21/2 -3/2 1 0 0 -3/64 1 0 3/16 1/16 -1/8 0 1 -1/8 0 0 33/2 3/2 -1 1 >> Bind(find(Bind==xl))=j-1 Bind = 3 2 7 >> % Iteration 4 >> j=5;inds = find(T(2:m,j) > 0)+1;ratios = T(inds,1)./T(inds,j);theta=min(ratios);candxjs = Bind(inds(find(ratios==theta))-1);[xl,l]=min(candxjs);l=find(Bind==xl)+1 l = 3 >> T(l,:)=T(l,:)/T(l,j); for i=setdiff(1:m,l);T(i,:) = T(i,:) - T(l,:)*T(i,j);end;T T = 0 -1/2 16 0 0 -1 1 0 0 -5/2 56 1 0 2 -6 0 0 -1/4 16/3 0 1 1/3 -2/3 0 1 4 -88 0 0 -4 10 1 >> Bind(find(Bind==xl))=j-1 Bind = 3 4 7 >> % Iteration 5 >> j=2;inds = find(T(2:m,j) > 0)+1;ratios = T(inds,1)./T(inds,j);theta=min(ratios);candxjs = Bind(inds(find(ratios==theta))-1);[xl,l]=min(candxjs);l=find(Bind==xl)+1 l = 4 >> T(l,:)=T(l,:)/T(l,j); for i=setdiff(1:m,l);T(i,:) = T(i,:) - T(l,:)*T(i,j);end;T T = 1/8 0 5 0 0 -3/2 9/4 1/8 5/8 0 1 1 0 -1/2 1/4 5/8 1/16 0 -1/6 0 1 1/12 -1/24 1/16 1/4 1 -22 0 0 -1 5/2 1/4 >> Bind(find(Bind==xl))=j-1 Bind = 3 4 1 >> % Iteration 6 >> j=6;inds = find(T(2:m,j) > 0)+1;ratios = T(inds,1)./T(inds,j);theta=min(ratios);candxjs = Bind(inds(find(ratios==theta))-1);[xl,l]=min(candxjs);l=find(Bind==xl)+1 l = 3 >> T(l,:)=T(l,:)/T(l,j); for i=setdiff(1:m,l);T(i,:) = T(i,:) - T(l,:)*T(i,j);end;T T = 5/4 0 2 0 18 0 3/2 5/4 1 0 * 1 6 0 * 1 3/4 0 -2 0 12 1 -1/2 3/4 1 1 -24 0 12 0 2 1 >> Bind(find(Bind==xl))=j-1 Bind = 3 5 1 %% We have the optimal solution now!