运筹单纯型算法.docx
《运筹单纯型算法.docx》由会员分享,可在线阅读,更多相关《运筹单纯型算法.docx(23页珍藏版)》请在冰豆网上搜索。
运筹单纯型算法
运筹学实验报告
题目:
用单纯形法求解线性规划问题
姓名王赛赛
学号2012436138
年级专业12数电类2
指导教师苏珂
2013年11月15日
一、实验目的:
运用MATLAB科学计算软件完成单纯型算法求解线性规划问题。
二、实验内容:
编写一个MATLAB的函数文件:
matrix.m,用于求解标准型的线性规划问题:
Ax=b
minf=c’*xsubjectto
x≥0
1.函数调用基本形式:
[x,minf,optmatrx,flag]=matrix(A,c,b)
2.参数介绍:
Ax=b
A:
线性规划问题的约束中各变量的系数组成的矩阵,是一个m×n的矩阵。
x≥0
c:
线性规划问题的目标函数f=c’*x中各变量的系数变量,是一个n维的行向量。
Ax=b
b:
线性规划问题的约束中的常数向量,是一个mn维德列向量。
x≥0
X:
输出线性规划问题的最优解,当线性规划问题没有可行解或者有可行解但没有优解X=[]。
minf:
输出线性规划问题的最优解,当线性规划问题没有可行解minf[],当线性规划问题有可行解但没有最优解时minf=-inf。
optmartx:
输出最优解对应的单纯型表,当线性规划问题没有可行解或者有可行解但没有最优解时optmartx=[]。
flag:
相性规划问题的求解结果标志值,当线性规划问题有最优解时flag=1,当线性规划问题有可行解但没有最优解时flag=0,当线性规划问题没有可行解时flag=-1.
三、MATLAB函数linprog的使用
Function[x,fval,exitflag,output,lambda]=linprog(f,A,B,Aeq,Beq,lb,ub,x0,options)
这个函数的功能是解决线性规划问题,它可以有变化的函数传递个数,函数内部根据传值的多少进行不同的功能操作。
一、函数的传值
1、X=LINPROG(f,A,b)试图解决线性问题:
求出f’*x的最小值,且f’x服从A*x<=b,其中x为自变量
2、X=LINPROG(f,A,b,Aeq,beq)解决满足添加约束等式:
Aeq*x=beq的上述问题
3、X=LINPROG(f,A,b,Aeq,beq,LB,UB)在以上的条件加下,与此同时定义了自变量的上下线,LB<=X<=UB,其中X为问题中的自变量。
如果问题中没有上下线(即没有UB和LB),则用空矩阵来代替。
如果问题中没有下线,则LB(i)=-Inf;如果没有上线,则UB(i)=Inf。
4、X=LINPROG(f,A,b,Aeq,beq,x0),在以上条件下,设置起点x0。
这一项是唯一可用的有效集算法。
默认的内点算法将忽略任何非空的起点。
5.X=LINPROG(F,A,B,AEQ,BEQ,LB,UBX0,OPTIONS),最大限度地减少使用默认的值在股权结构。
二、函数的返回值
1、[X,FVAL]=LINPROG(f,A,b),返回值FVAL是目标函数的最优解,X是相应于最优解的为指数的向量
2、[X,FVAL,EXITFLAG]=LINPROG(f,A,b),前两个与上述相同,最后一个EXITFLAG返回z值描述退出函数的不同情况,其中,可能返回的EXITFLAG的值及其对应的情况如下:
1——LINPROG收敛于结果X。
0——达到最大迭代次数
-2——找不到可行解
-3——问题无界
-4——当执行算法的时候出现不确定解
-5——原问题及其对偶问题都无可行解
-7——搜索方向的数量级过小;没有办法进行更加深入的搜索。
问题条件有错或者变量有错误。
3、[X,FVAL,EXITFLAG,OUTPUT]=LINPROG(f,A,b),前面与上相同,最后返回一个OUIPUI的结构,此结构是带有从OUOPUT中得到的iterations的数值的最大值。
参数Iterations是OUIPUI中约束问题的最大值。
参数Constrviolation是一个在OUIPUI算法中使用的算法。
Algorithm是在OUIPUI中共轭梯度的iterations的值。
4、[X,FVAL,EXITFLAG,OUTPUT,LAMBDA]=LINPROG(f,A,b),前面与上相同,最后返回拉格朗日乘法算子LAMBDA,在这个解中:
LAMBDA.ineqlin代表先行不等式A,LAMBDA.eqlin代表线性等式Aeq,LAMBDA.lower代表LB,LAMBDA.upper代表UB。
四、matrix函数
%两阶段算法求解线性规划问题
%子函数function使用格式;
%function[输出参数列表]=函数名(输入参数列表)
%函数体
%输入参数列表中的c,A,b的含义
%c目标函数中的系数向量
%A约束中的系数矩阵
%b约束中的常数向量
%输出参数列表中的x,minf,optmatrx,flag的含义:
%x最优解
%minf最优值
%optmatrx最优值对应的单纯形表
%flag线性规划求解结果的标志值
function[x,minf,optmatrx,flag]=linpa(A,b,c)
A=input('请输入约束中的系数矩阵A=')
b=input('请输入约束中的常数向量b=')%b为单列矩阵
c=input('请输入目标函数中的系数向量c=')%c为单行矩阵,注意目标函数中的系数向量c=-c
%判断c是否为单行矩阵,b是否为单列矩阵
%
[i,j]=size(c)
%如果输入的c为单列矩阵
ifj==1
c=c'
disp('输入的目标函数系数向量c为单位矩阵,已将其转化为单行矩阵')
fprintf('\nc=\n')
disp(c)
end
[i,j]=size(b)
%如果输入的b为单行矩阵
ifi==1
b=b'
disp('输入的约束中的常数向量b为单行矩阵,已将其转化为单列矩阵')
fprintf('\nb=\n')
disp(b)
end
%
%将c=-c处理
c=-c
%
%size()用于获取数组的行数和列数
[m,n]=size(A);%第一个输出变量m用于存储数组的行数,第二个输出变量用于存储数组的列数
m1=m;%保存下原来的行数
s=eye(m);%建立大小行数m的单位矩阵,eye()用于获取单位矩阵
A=[As];
A=[Ab];
g1=zeros(1,n);
x=ones(1,m);
g1=[g1-x];
g=[0];
g1=[g1g];%人工变量的检验向量
s1=n+m+1;%目前列数之和
s2=zeros(1,m+1);
c=[cs2];
A1=zeros(m,1);
fori1=1:
m
A1(i1,1)=i1+n;%基变量的数值存储区
end
fori=1:
m
g1(1,:
)=g1(1,:
)+A(i,:
);
end
%%***********************第一阶段*****************************
decide=find(g1(1,1:
m+n)>0);%寻找g1中大于零的数值列数存于dicide数组中
while~isempty(decide)%当存在大于零的值时
i=decide
(1);%当列数赋给i
text=find(A(1:
m,i)>0);%text存放该列中所有数值大于零的行数
ifisempty(text)%当没有大于零的数值时无解
flag=0;
break;
end
min=inf;
fori1=1:
m%当该列存在大于零的数值时
ifA(i1,i)>0&A(i1,s1)/A(i1,i)min=A(i1,s1)/A(i1,i);
x1=i1;%将比值最小的行数赋给x1
end
end
A(x1,:
)=A(x1,:
)/A(x1,i);%单位化
c(1,:
)=c(1,:
)+(-1*c(1,i)*A(x1,:
));
g1(1,:
)=g1(1,:
)+(-1*g1(1,i)*A(x1,:
));
fori1=1:
m
ifi1~=x1
A(i1,:
)=A(i1,:
)+(-1*A(i1,i)*A(x1,:
));%进行换基迭代
end
end
A1(x1,1)=i;%将列数既对应的非基变量转换为基变量
decide=find(g1(1,[1:
m+n])>0);%再进行查找
end
%%*************************************************************
ifg1(1,s1)>0%此时对应的为无解情况
flag=-1;
end
%%*******************判断人工变量**************************
ifg1(1,s1)==0
g1(1,:
)=[];
fori6=1:
m
A(:
n+1)=[];
end
fori6=1:
m
c(n+1)=[];
end
decide=find(A1(1:
m,1)>n);%查找是否有人工变量
if~isempty(decide)%存在人工变量
while~isempty(decide)
x1=decide
(1);%存放的是人工变量的行数
text=find(A(x1,1:
n)~=0);%该行的所有元素都不为零的列坐标
ifisempty(text)%都为零,则该行为多余的行,删除
decide
(1)=[];
A1(x1,:
)=[];
A(x1,:
)=[];
m=m-1;
else
i=text
(1);%i保存的是列数
A(x1,:
)=A(x1,:
)/A(x1,i);
A1(x1,1)=i;
c(1,:
)=c(1,:
)+(-1*c(1,i)*A(x1,:
));
fori1=1:
m
ifi1~=x1
A(i1,:
)=A(i1,:
)+(-1*A(i1,i)*A(x1,:
));
end
end
decide
(1)=[];
text
(1)=[];
end
end
end
%%**************************第二阶段**************************
decide=find(c(1,1:
n)>0);%查找检验向量是否存在大于零的数,将列数保存
while~isempty(decide)
i=decide
(1);%将列数赋给i
text=find(A(:
i)>0);%保存大于零的项的行数
ifisempty(text)%为空的时候,则无解
disp('有可行解但没有最优解')
flag=0
c=[c;A];
optmatrx=c;
x=zeros(n,1);
foro=1:
n
foro1=1:
m
ifo==A1(o1,1)
x(o,1)=A(o1,n+1);
end
end
end
minf=-inf
optmatrx
x
return;
end
min=inf;
fori1=1:
m
ifA(i1,i)>0&A(i1,n+1)/A(i1,i)min=A(i1,n+1)/A(i1,i);
x1=i1;%保存最小项的行数
end
end
A(x1,:
)=A(x1,:
)/A(x1,i);%单位化
j=c(1,i);
forj1=1:
n+1
c(1,j1)=c(1,j1)+A(x1,j1)*(-1)*j;%换基迭代
end
fori1=1:
m
ifi1~=x1
A(i1,:
)=A(i1,:
)+(-1*A(x1,:
)*A(i1,i));
end
end
A1(x1,1)=i;
decide=find(c(1,1:
n)>0);
end
%%*************************结果赋值结算***********************
%%*
c=[c;A];
optmatrx=c;
x=zeros(n,1);
foro=1:
n
foro1=1:
m
ifo==A1(o1,1)
x(o,1)=A(o1,n+1);
end
end
end
minf=c(1,n+1);
flag=1;
end
switchflag
case1%有可行解
disp('输出最优解对应的单纯形表')
optmatrx
disp('输出最优解')
x
disp('输出线性规划问题的最优值')
minf
case0%有可行解但没有最优值
disp('输出可行解对应的单纯形表')
optmatrx
disp('输出可行解')
x
disp('输出线性规划问题的最优值')
minf
case-1%没有可行解
disp('没有可行解â')
disp('输出最优解对应的单纯形表')
optmatrx=[]
disp('输出最优解')
x=[]
disp('输出线性规划问题的最优值')
minf=-inf
end
五、例题求解
例题1(有最优解)
题目:
minz=4x1+3x3
s.t.1/2x1+x2+1/2x3–2/3x4=2
3/2x1+1/2x3=3
3x1-6x2+4x4=0
xj≥0,j=1,…,4
输入及函数的调用:
测试结果:
经过测试当线性规划问题存在最优解时,程序运行没有问题。
例题2(有可行解但无最优解)
题目:
minz=-2x1-5x2
s.t.-1x1+x3=4
x2+x4=3
-x1+2x2+x5=8
xj≥0,j=1,…,5
输入及函数的调用:
测试结果:
经过测试当线性规划问题有可行解但无最优解时,程序运行没有问题。
例题3(没有可行解)
题目:
minz=2x1+2x2
s.t.-1x1+x2-x3=1
-x1-x2-x4=2
xj≥0,j=1,…,4
输入及函数的调用:
测试结果:
经过测试当线性规划问题没有可行解时,程序运行没有问题。
>>linpa
请输入约束中的系数矩阵A=[100001-100000;1100000-10000;01100000-1000;001100000-100;0001100000-10;00001100000-1]
A=
Columns1through11
100001-10000
1100000-1000
01100000-100
001100000-10
0001100000-1
00001100000
Column12
0
0
0
0
0
-1
请输入约束中的常数向量b=[101525201812]
b=
101525201812
请输入目标函数中的系数向量c=[111111000000]
c=
Columns1through11
11111100000
Column12
0
i=
1
j=
12
i=
1
j=
6
b=
10
15
25
20
18
12
输入的约束中的常数向量b为单行矩阵,已将其转化为单列矩阵
b=
10
15
25
20
18
12
c=
Columns1through11
-1-1-1-1-1-100000
Column12
0
输出最优解对应的单纯形表
optmatrx=
Columns1through11
000000-10-10-1
100001-10000
01000-100-11-1
001001000-11
00010-10000-1
000000-11-11-1
00001100000
Columns12through13
053
010
111
-114
16
16
-112
输出最优解
x=
10
11
14
6
12
0
0
6
0
0
0
0
输出线性规划问题的最优值
minf=
53
ans=
10
11
14
6
12
0
0
6
0
0
0
0