当前位置:文档之家› (完整版)运筹学实验报告

(完整版)运筹学实验报告

运筹学实验报告班级:数电四班姓名:刘文搏学号:一、实验目的运用MATLAB程序设计语言完成单纯性算法求解线性规划问题。

二、实验内容编写一个MATLAB的函数文件:linp.m用于求解标准形的线性规划问题:min f=c*x subject to :A*x=b ; x>=0;1、函数基本调用形式:[x,minf,optmatrx,flag]=linp(A,b,c)2、参数介绍:A:线性规划问题的约束A*x=b且x>=0中变量的系数组成的矩阵,是一个m*n的矩阵。

c :线性规划问题的目标函数f=c*x中各变量的系数向量,是一个n 维的行向量。

b :线性规划问题的约束A*x=b且x>=0中的常数向量,是一个m维的列向量。

x :输出线性规划问题的最优解,当线性规划问题没有可行解或有可行解无最优解时x=[].minf :输出线性规划问题的最优值,当线性规划问题没有可行解时minf=[],当线性规划问题有可行解无最优解时minf=-Inf。

flag :线性规划问题的求解结果标志值,当线性规划问题有最优解时flag=1,当线性规划问题有可行解无最优解时flag=0,当线性规划问题没有可行解时flag=-1.cpt:输出最优解对应的单纯性表,当线性规划问题没有可行解或有可行解无最优解时cpt=[].三、Linp函数%此函数是使用两阶段算法求解线性规划问题function [x,minf,flag,cpt]=linp(A,b,c);for i=1:p %判断b是否<=0;将b转换成大于0;if b(i)<0A(i,:)=-1*A(i,:);b(i)=-1*b(i);endend%返回值:x,第一张单纯形表,基,标志参数 A,c,b%********第一张单纯形表的初始化[m,n]=size(A);%获得矩阵A的维数[p,q]=size(b);dcxb=zeros(m+2,m+n+1);%确定第一张单纯形表的大小dcxb(1,:)=[-c,zeros(1,m+1)];%¸给表的第一行赋值dcxb(2,:)=[zeros(1,n),-1*ones(1,m),0];%¸给表的第二行赋值dcxb([3:m+2],:)=[A,eye(m,m),b];%添A和b到表中jxl=[n+1:n+m];for i=3:m+2dcxb(2,:)=dcxb(2,:)+dcxb(i,:);for i=3:m+2dxcb(2,:)=dxcb1(2,:)+dxcb1(i,:);enddxcb;%************辅助问题换基迭代**********************dyl=find(dcxb(2,[1:m+n])>0);while ~isempty(dyl)firstnum=dyl(1);dll=dcxb([3:m+2],firstnum);youduanb=dcxb([3:m+2],m+n+1);look=find(dll>0);if isempty(look)dcxb(2,firstnum)=0;elsemin=Inf;for i=3:m+2if dll(i-2)>0&youduanb(i-2)/dll(i-2)<min min=youduanb(i-2)/dll(i-2);line1=i;endenddcxb(line1,:)=dcxb(line1,:)/dcxb(line1,firstnum);for i=1:m+2if i~=line1dcxb(i,:)=dcxb(i,:)+(-1*dcxb(i,firstnum)*dcxb(line1,: ));endendjxl(line1-2)=firstnum;enddyl=find(dcxb(2,[1:m+n])>0);dcxbendif dcxb(2,m+n+1)>0fprintf('g>0,´此问题没有可行解');x=[];minf=inf;cpt=[];flag=-1;returnendlook1=find(jxl>n);if dcxb(2,m+n+1)==0 %等于0,判断基变量中是否有人工变量;if ~isempty(look1)%´存在时进行处理while ~isempty(look1)line2=look1(1)+2;chdy0=find(dcxb(line2,[1:n])~=0);if isempty(chdy0)%´存在人工变量都为零的那一行,去掉该行dcxb(line2,:)=[];look1(1)=[];jxl(line2-2)=[];else%否则进行换基迭代secondnum=chdy0(1);dcxb(line2,:)=dcxb(line2,:)/dcxb(line2,secondnum);jxl(line2-2)=secondnum;for i=1:m+2if i~=line2dcxb(i,:)=dcxb(i,:)+(-1*dcxb(i,secondnum)*dcxb(line2, :));end;endlook1(1)=[];end;endendend%去掉人工变量,得到单纯性的第一张表dcxb(2,:)=[];dcxb(:,[n+1:n+m])=[];%有可行解,判断zdcxb2=dcxb;look2=find(dcxb2(1,[1:n])>0);while ~isempty(look2)thirdnum=look2(1);duilie=dcxb2([2:m+1],thirdnum);youduanb1=dcxb2([2:m+1],n+1);look3=find(duilie>0);if isempty(look3)fprintf('´此问题有可行解,但没有最优解'); x=zeros(n,1);[mi,n1]=size(jxl);for i=1:n1x(jxl(i))=dcxb2(i+1,n+1);endfprintf(' 可行解为');xminf=-Infcpt=[]flag=0returnendmin1=Inf;for i=1:mif duilie(i)>0&youduanb1(i)/duilie(i)<min1%找最小比值min1=youduanb1(i)/duilie(i);line=i+1;%记录行数endenddcxb2(line,:)=dcxb2(line,:)/dcxb2(line,thirdnum);for i=1:m+1if i~=linedcxb2(i,:)=dcxb2(i,:)+(-1*dcxb2(i,thirdnum)*dcxb2(lin e,:));endendjxl(line-1)=thirdnum;dcxb2look2=find(dcxb2(1,[1:n])>0);endminf=dcxb2(1,n+1);x=zeros(n,1);[p,q]=size(jxl);fprintf('\最优解已找到n');for i=1:qx(jxl(i))=dcxb2(i+1,n+1); endfprintf('最优可行解为:');xfprintf('最优值为:');minfcpt=dcxb2;fprintf('最优解对应的单纯形表为:');cptflag=1return例题1.A=[1/2 1 1/2 -2/3; 3/2 0 3/4 0];b=[2; 3];c=[4 0 3 0];运行结果:>> [x,minf,flag,cpt]=linp(A,b,c)请一次输入系数矩阵A;输入右端向量b; 输入所求问题的向量cA=[1/2 1 1/2 -2/3; 3/2 0 3/4 0];b=[2; 3];c=[4 0 3 0];dcxb =-4.0000 0 -3.0000 0 0 0 02.0000 1.0000 1.2500 -0.6667 0 05.00000.5000 1.0000 0.5000 -0.6667 1.0000 02.00001.5000 0 0.7500 0 0 1.00003.0000dcxb =0 0 -1.0000 0 0 2.66678.00000 1.0000 0.2500 -0.6667 0 -1.33331.00000 1.0000 0.2500 -0.6667 1.0000 -0.33331.00001.0000 0 0.5000 0 0 0.66672.0000dcxb =0 0 -1.0000 0 0 2.66678.00000 0 0 0 -1.0000 -1.0000 00 1.0000 0.2500 -0.6667 1.0000 -0.33331.00001.0000 0 0.5000 0 0 0.66672.0000最优解已找到!最优可行解为:x =21最优值为:minf =8最优解对应的单纯形表为:cpt =0 0 -1.0000 0 8.00000 1.0000 0.2500 -0.6667 1.00001.0000 0 0.5000 02.0000flag =1ans =21例题2A=[1/2 1 1/2 -2/3; 3/2 0 3/4 0;3 -6 0 4];b=[2; 3;0];c=[4 0 3 0];运行结果>> [x,minf,flag,cpt]=linp(A,b,c)请一次输入系数矩阵A;输入右端向量b; 输入所求问题的向量cA=[1/2 1 1/2 -2/3; 3/2 0 3/4 0;3 -6 0 4];b=[2; 3;0];c=[4 0 3 0];dcxb =Columns 1 through 6-4.0000 0 -3.0000 0 0 0 5.0000 -5.0000 1.2500 3.3333 0 00.5000 1.0000 0.5000 -0.6667 1.0000 01.5000 0 0.7500 0 0 1.0000 3.0000 -6.0000 0 4.0000 0 0 Columns 7 through 80 00 5.00000 2.00000 3.0000dcxb =Columns 1 through 60 -8.0000 -3.0000 5.3333 0 0 0 5.0000 1.2500 -3.3333 0 0 0 2.0000 0.5000 -1.3333 1.0000 00 3.0000 0.7500 -2.0000 0 1.00001.0000 -2.0000 0 1.3333 0 0 Columns 7 through 81.3333 0-1.6667 5.0000-0.1667 2.0000-0.5000 3.00000.3333 0dcxb =Columns 1 through 60 0 -1.0000 0 4.0000 0 0 0 0 0.0000 -2.5000 0 0 1.0000 0.2500 -0.6667 0.5000 00 0 0 0 -1.5000 1.00001.0000 0 0.5000 0 1.0000 0 Columns 7 through 80.6667 8.0000-1.2500 0-0.0833 1.00000.1667 2.0000dcxb =Columns 1 through 60 0 -1.0000 0 4.0000 0 0 0 0 0 -2.5000 0 0 1.0000 0.2500 -0.6667 0.5000 00 0 0 0 -1.5000 1.00001.0000 0 0.5000 0 1.0000 0 Columns 7 through 80.6667 8.0000-1.2500 0-0.0833 1.0000-0.2500 00.1667 2.0000最优解已找到!最优可行解为:x =21最优值为:minf =8最优解对应的单纯形表为:cpt =0 0 -1.0000 0 8.00000 1.0000 0.2500 -0.6667 1.00001.0000 0 0.5000 02.0000flag =1ans =21例题3A= [1/2 1 1/2 -2/3; 3/2 0 1/2 0;3 -6 0 4];b=[2; 3;0];c=[4 0 3 0];运行结果:>> [x,minf,flag,cpt]=linp(A,b,c)请一次输入系数矩阵A;输入右端向量b; 输入所求问题的向量cA=[1/2 1 1/2 -2/3; 3/2 0 1/2 0;3 -6 0 4];b=[2; 3;0];c=[4 0 3 0];dcxb =Columns 1 through 6-4.0000 0 -3.0000 0 0 05.0000 -5.0000 1.0000 3.3333 0 00.5000 1.0000 0.5000 -0.6667 1.0000 01.5000 0 0.5000 0 0 1.0000 3.0000 -6.0000 0 4.0000 0 0 Columns 7 through 80 00 5.00000 2.00000 3.00001.0000 0dcxb =Columns 1 through 60 -8.0000 -3.0000 5.3333 0 0 0 5.0000 1.0000 -3.3333 0 0 0 2.0000 0.5000 -1.3333 1.0000 00 3.0000 0.5000 -2.0000 0 1.00001.0000 -2.0000 0 1.3333 0 0 Columns 7 through 81.3333 0-1.6667 5.0000-0.1667 2.0000-0.5000 3.00000.3333 0dcxb =Columns 1 through 60 0 -1.0000 0 4.0000 0 0 0 -0.2500 0.0000 -2.5000 0 0 1.0000 0.2500 -0.6667 0.5000 00 0 -0.2500 0 -1.5000 1.00001.0000 0 0.5000 0 1.0000 0 Columns 7 through 80.6667 8.0000-1.2500 0-0.0833 1.0000-0.2500 00.1667 2.0000dcxb =Columns 1 through 60 0 -1.0000 0 4.0000 0 0 0 -0.2500 0 -2.5000 0 0 1.0000 0.2500 -0.6667 0.5000 00 0 -0.2500 0 -1.5000 1.00001.0000 0 0.5000 0 1.0000 0 Columns 7 through 80.6667 8.0000-1.2500 0-0.0833 1.0000-0.2500 00.1667 2.0000最优解已找到!最优可行解为:x =21最优值为:minf =8最优解对应的单纯形表为:cpt =0 0 0 0 8.0000 0 1.0000 0 -0.6667 1.00000 0 1.0000 0 01.0000 0 0 02.0000 flag =1ans =21例题4A=[-1 1 -1 0;-1 -1 0 -1];b=[1;2];c=[2 2 0 0];运行结果>>[x,minf,flag,cpt]=linp(A,b,c)请一次输入系数矩阵A;输入右端向量b; 输入所求问题的向量c A=[-1 1 -1 0;-1 -1 0 -1];b=[1;2];c=[2 2 0 0];dcxb =-2 -2 0 0 0 0 0-2 0 -1 -1 0 0 3-1 1 -1 0 1 0 1-1 -1 0 -1 0 1 2g>0,此问题没有可行解ans =[]例题5A=[-1 0 1 0 0;0 1 0 1 0;-1 2 0 0 1];b=[4; 3 ; 8];c=[-2 -5 0 0 0];运行结果:>> [x,minf,flag,cpt]=linp(A,b,c)请一次输入系数矩阵A;输入右端向量b; 输入所求问题的向量c A=[-1 0 1 0 0;0 1 0 1 0;-1 2 0 0 1];b=[4; 3 ; 8];c=[-2 -5 0 0 0];dcxb =2 5 0 0 0 0 0 0 0 -23 1 1 1 0 0 0 15 -1 0 1 0 0 1 0 04 0 1 0 1 0 0 1 0 3 -1 2 0 0 1 0 0 1 8 dcxb =2 0 0 -5 0 0 -5 0 -15 -2 0 1 -2 1 0 -3 0 6 -1 0 1 0 0 1 0 04 0 1 0 1 0 0 1 0 3 -1 0 0 -2 1 0 -2 1 2 dcxb =2 0 0 -5 0 0 -5 0 -15 -1 0 0 -2 1 -1 -3 0 2 -1 0 1 0 0 1 0 04 0 1 0 1 0 0 1 0 3 -1 0 0 -2 1 0 -2 1 2 dcxb =2 0 0 -5 0 0 -5 0 -15 0 0 0 0 0 -1 -1 -1 0 -1 0 1 0 0 1 0 0 4 0 1 0 1 0 0 1 03 -1 0 0 -2 1 0 -2 1 2此问题有可行解,但是没有最优解可行解为:x =342minf =-Infcpt =[]flag =ans =342。

相关主题