运筹学线性规划二阶段法matlab矩阵描述

function Two_stage(zflag,b,A,Aflag,c) clear;clc format rat fprintf('\n'); %zflag=input('请输入所求目标函数的状态(求最大值为1，求最小值为-1):'); %b=input('请输入资源矩阵b:'); %A=input('请输入约束矩阵A:'); %Aflag=input('请输入约束方程的状态矩阵(小于号为-1，等于为0，大于号为1);'); %c=input('请输入目标函数价值系数矩阵c:'); zflag=1; b=[9 15 5]; A=[5 3 1;-5 6 15;2 1 1]; c=[10 15 12]; Aflag=[-1 -1 1]; %-------------------------------------------------标准化 if zflag==-1 c=-1*c; end nb=size(b); [mA,nA]=size(A); for i=1:nb if b(i)<0 b(i)=-1*b(i); A(i,:)=-1*A(i,:); end end num1=find(Aflag==-1); num2=find(Aflag==0); num3=find(Aflag==1); N=[]; c1=zeros(1,length(c)); %构造第一阶段目标函数系数矩阵 for i=1:length(num1) A(num1(i),nA+1)=1; c(nA+1)=0; c1(nA+1)=0; nA=nA+1; N(i)=nA; end for i=1:length(num3) A(num3(i),nA+1)=-1; A(num3(i),nA+2)=1; c(nA+1)=0; c1(nA+1)=0; c1(nA+2)=-1; nA=nA+2; N(end+1)=nA; end for i=1:length(num2) A(num2(i),nA+1)=1; c1(nA+1)=-1; nA=nA+1; N(end+1)=nA; end %以下部分初始化基变量的下标值矩阵N,基变量的系数矩阵CB %-------------------------------------------------------- CB=zeros(1,mA); [mA,nA]=size(A); for i=1:mA CB(i)=c1(N(i)) ; %确定基变量的系数 end if (length(num2)+length(num3))~=0 [A,b,c1,N,CB]=first_stage(A,b,c1,N,CB); [A,b,c,N,CB]=second_stage(A,b,c,c1,N,CB); else fprintf('标准化后初始单纯形表为：'); c p=[CB',N',b',A] disp('***********************************************************'); [A,b,c,N,CB,t1]=second_stage_diedai(A,b,c,N,CB);%只调用第二阶段迭代过程 end % 单纯形法第一阶段函数 %============================================================ function [A,b,c1,N,CB]=first_stage(A,b,c1,N,CB) k=1; % k是 times=0; Isend=0; temp1=' 第';temp2='次后单纯形表为：'; b=b'; while Isend==0 if (times==0) fprintf('标准化后初始单纯形表为：'); else fprintf('%s%d%s',temp1,times,temp2); end times=times+1; [mA,nA]=size(A); t=zeros(1,mA); for j=1:nA %start for t(j)=c1(j)-CB*A(:,j); end %End for c1 p=[CB',N',b,A] t [m1,j]=max(t); %确定导入基，导入基的列下标为col k=max(t); if m1<=0 X=zeros(1,nA); for i=1:length(N) X(N(i))=b(i); end break; %跳出while 循环 else times2=0; for i=1:mA if A(i,j)>0 f(i)=b(i)/A(i,j); else f(i)=10000; %如果相除是负值，把值设为10000 times2=times2+1; end end if times2==mA fprintf('此问题的无最优解\n') Isend=1; break; %跳出if m1<0 else [m2,l]=min(f); CB(l)=c1(j);% 确定导出基 导出基的行标为 l N(l)=j; B=[]; for i=1:length(N) B=[B,A(:,N(i))]; end b=B\b ; % end the for cycle,已改变了b 数组的值 A=B\A; end %endif end %endif end %end while Z=c1*(X'); if Z~=0 fprintf('人工变量仍在基变量中，该问题无可行解。\n') else fprintf('第一阶段最优解为：\n') X end disp('***********************************************************'); %*****************************************************第一阶段完毕 % 单纯形法第二阶段函数 % function [A,b,c,N,CB]=second_stage(A,b,c,c1,N,CB) if isempty(find(CB==-1))%判断最优解的基变量中是否含有非零的人工变量 n=length(find(c1==-1));%统计人工变量的个数 [mA,nA]=size(A); A=A(:,1:nA-n);%除去人工变量后的矩阵A [A,b,c,N,CB,t1]=second_stage_diedai(A,b,c,N,CB);%调用第二阶段迭代函数 else fprintf('此问题无解！') end % 单纯形法第二阶段迭代函数 % function [A,b,c,N,CB,t1]=second_stage_diedai(A,b,c,N,CB) [mA,nA]=size(A); for i=1:length(N) CB(i)=c(N(i)); end k1=1; times=0; Isend=0; while Isend==0 if (times==0) fprintf('初始单纯形表为：'); else fprintf('%s%d%s\n','第',times,'次迭代的单纯性表'); end times=times+1; for j=1:nA%(nA-n) t1(j)=c(j)-CB*A(:,j); %计算检验数 end p=[CB',N',b,A] c,t1 [m3,column]=max(t1);%找出检验数中最大的值及其所在的列 k1=m3; j=column; ftimes=0; if m3>0 for i=1:mA if A(i,j)*b(i)>0 f(i)=b(i)/A(i,j); else f(i)=10000; %如果相除是负值，把值设为10000 ftimes=ftimes+1; end end if ftimes==mA fprintf('此问题的无最优解\n') Isend=1; else [m4,row]=min(f); CB(row)=c(column); N(row)=column; B=[]; for i=1:length(N) B=[B,A(:,N(i))]; end b=B\b; % end the for cycle,已改变了b 数组的值 A=B\A; end else Isend=1; X=zeros(1,nA); for i=1:length(N) X(N(i))=b(i); fprintf('此问题的最优解为：\n') X fprintf('目标函数的最优值为：\n') Z=c*(X') disp('***********************************************************'); end end %end if end%end while %*******************************************************第二阶段迭代函数结束