单纯形法MATLAB程序(线性与非线性规划大作业)
单纯形法matlab
数学软件与实验数学与信息科学学院信息与计算科学单纯形法的Matlab程序如下：function [xx,fm]=myprgmh(m,n,A,b,c)B0=A(:,1:m);cb=c(:,1:m);xx=1:n;sgm=c-cb*B0^-1*A;h=-1;sta=ones(m,1);for i=m+1:nif sgm(i)>0h=1;endendwhile h>0[msg,mk]=max(sgm);for i=1:msta(i)=b(i)/A(i,mk);end[mst,mr]=min(sta);zy=A(mr,mk);for i=1:mif i==mrfor j=1:nA(i,j)=A(i,j)/zy;endb(i)=b(i)/zy;endendfor i=1:mif i~=mrfor j=1:nA(i,j)=A(i,j)-A(i,mk)*A(mr,j);endb(i)=b(i)-A(i,mk)*b(mr);endendB1=A(:,1:m);cb(mr)=c(mk);xx(mr)=mk;sgm=c-cb*B1*A;for i=m+1:nif sgm(i)>0h=1;endendendfm=c*xx;例题：编写下列求解如下线性规划问题的单纯形法函数min f'xs.t ax<=b(其中b>=0)函数形式function [x,fval,it,op]=singl(f,a,b) 输出中x为最优解fval为最优值it为迭代次数无最优解op=0有最优解op=1编写程序如下：function [x,fval,it,op]=singl(f,a,b)[m,n]=size(a);c=[a eye(m) b;f' zeros(1,m+1)];fval=0;x=zeros(m+n,1);op=1;it=0;e=zeros(1,m);lie=find(f<0);l=length(lie);while(l>0)for j=1:ld=find(c(:,lie(j)));d_l=length(d);if d_l>0for i=1:mif c(i,lie(j))>0e(i)=c(i,end)/c(i,lie(j));elsee(i)=inf;endend[g,h]=min(e);for w=1:m+1if w==hc(w,:)=c(w,:)/c(h,lie(j));elsec(w,:)=c(w,:)-c(h,:)*c(w,lie(j))/c(h,lie(j));endendit=it+1;elseop=0;endendlie=find(c(end,:)<0);l=length(lie);endfor i=1:(m+n)ix=find(c(:,i));if(length(ix)==1)&(ix<=m)&(c(ix,i)==1) x(i)=c(ix,end)elsex(i)=0endendfval=-c(end,end);。
matlab单纯形法
%求解标准型线性规划:max c*x;s.t. A*x=b;x>=0%本函数中的A是单纯初始表，包括:最后一行是初始的检验数，最后一列是资源向量b %N是初始的基变量的下标%输出变量sol是最优解%输出变量val是最优值，kk是迭代次数function [sol,val,kk]=ssimplex(A,N)[mA,nA]=size(A);kk=0; %迭代次数flag=1;while flagkk=kk+1;if A(mA,:)<=0 %已找到最优解flag=0;sol=zeros(1,nA-1);%给每个变量赋初值0for i=1:mA-1sol(N(i))=A(i,nA);%给基变量赋新值(替换0)end %给出最优解val=-A(mA,nA);elsefor i=1:nA-1if A(mA,i)>0&A(1:mA-1,i)<=0 %问题有无界解disp('have infinite solution!');flag=0;break;endendif flag %还不是最优表，进行转轴运算temp=0;for i=1:nA-1if A(mA,i)>temptemp=A(mA,i);inb=i; % 进基变量的下标endend %选择最大检验数纵向对应的变量为进基变量sita=zeros(1,mA-1);for i=1:mA-1if A(i,inb)>0sita(i)=A(i,nA)/A(i,inb);endendtemp=inf;for i=1:mA-1if sita(i)>0&sita(i)<temptemp=sita(i);outb=i; %出基变量下标endend %选择最小的sita横向对应的变量为出基变量%以下更新Nfor i=1:mA-1if i==outbN(i)=inb;%以进基变量的下标替代出基变量的下标endend%以下进行转轴运算A(outb,:)=A(outb,:)/A(outb,inb);%将主元化为1for i=1:mAif i~=outbA(i,:)=A(i,:)-A(outb,:)*A(i,inb);%将进基变量所在列除主元外的其余元素化为0endendendendend。
单纯形法MATLAB程序(线性与非线性规划大作业)
单纯形法程序整理.txt [hB,lB]=size(B);%测算B的行、列数 for j=1:hB if(A(B(j,1),B(j,2))~=1)%如果主元不为1，将对应行除以主元化为1 A(B(j,1),:)=A(B(j,1),:)/A(B(j,1),B(j,2)); end for i=1:hA%将其他元化为零 if(i==B(j,1)) continue; end A(i,:)=A(i,:)-A(B(j,1),:)*A(i,B(j,2)); end end end %————————————————% %% 子程序-选主元 function [r,s]=SubFc_XuanZhuYuan(TableK,h,l) %选主元素 % 选出表矩阵TableK中的主元素 %% 选列主元 t=0; for i=1:l if(TableK(h+1,i)<=t) t=TableK(h+1,i);%记录判别数最小值 s=i;%记录下标 end end if(t==0) error('判别数全为正'); end %% 选行主元 for i=1:h%检查元素是否为负，并赋予t初值 if(TableK(i,s)>0) t=TableK(i,l+1)/TableK(i,s); break; end if(i==h)%若元素全为负，则报错 error('主列元素全为非正，规划问题有无界解！'); end end %k=1; for i=1:h %若未知数出现负值，则是初始迭代标准性表产生的，应更换初始主元； if(TableK(i,l+1)<0) error('b向量出现负值，更换初始主元'); end if(TableK(i,s)<=0)%排除A表负元素 continue; end if(TableK(i,l+1)/TableK(i,s)<=t) t=TableK(i,l+1)/TableK(i,s);%记录判别数最小值 r=i;%记录下标 end end %补充判据 k=1; 第 2 页
单纯形法matlab代码
单纯形法（Simplex Method）——优化问题的强有力工具一、引言优化问题在数学和工程领域扮演着重要的角色，而单纯形法则是一种常用的解决优化问题的方法。
它可以用于线性规划问题的求解，通过逐步迭代，不断优化目标函数的值。
本文将介绍单纯形法的原理和基本步骤，并使用MATLAB代码展示其在实际问题中的应用。
二、单纯形法原理单纯形法是一种基于几何直觉的算法，它通过多次迭代寻找可行解空间中的顶点，并在每次迭代中逐渐改进目标函数的值，直至找到最优解或确定问题无解。
该方法的基本思想是从初始可行解出发，通过交换基变量和非基变量，不断向较优的顶点移动，直至达到最优解。
三、单纯形法步骤单纯形法的求解过程可以分为以下几个步骤：3.1 构建初始单纯形表构建初始单纯形表包括将优化问题转化为标准型，并引入松弛变量将约束条件转化为等式，同时引入人工变量以确保可行解存在。
3.2 选择入基变量和出基变量在每一次迭代中，需要选择一个入基变量和一个出基变量。
入基变量是指由非基变量变为基变量的变量，而出基变量是指由基变量变为非基变量的变量。
3.3 计算回代数通过计算回代数，可以确定迭代的方向和距离。
回代数是指基变量要离开基础解所需要沿其可行方向运动的最大距离。
3.4 更新单纯形表通过更新单纯形表，可以得到下一次迭代的基变量和非基变量，并计算相应的解。
更新单纯形表的方法一般有高斯型和乘子分析型两种。
3.5 判断是否达到终止条件在每一次迭代后，需要判断是否满足终止条件。
终止条件可以是目标函数的值不再改变或约束条件不再发生变化等。
3.6 迭代直至达到最优解如果未达到终止条件，则继续进行下一次迭代，直至达到最优解或确定问题无解。
四、单纯形法在MATLAB中的应用在MATLAB中，可以使用线性规划工具箱（Linear Programming Toolbox）来实现单纯形法的运算。
以下是一个简单的例子，以详细介绍如何使用MATLAB代码解决线性规划问题。
实验二：MATLAB编程单纯形法求解
北京联合大学实验报告项目名称：运筹学专题实验报告学院：自动化专业：物流工程班级： 1201B 学号：2012100358081 姓名：管水城成绩：2015 年 5 月 6 日实验二：MATLAB 编程单纯形法求解一、实验目的：(1)使学生在程序设计方面得到进一步的训练;，掌握Matlab (C 或VB)语言进行程序设计中一些常用方法。
(2)使学生对线性规划的单纯形法有更深的理解. 二、实验用仪器设备、器材或软件环境 计算机, Matlab R2006三、算法步骤、计算框图、计算程序等本实验主要编写如下线性规划问题的计算程序:⎩⎨⎧≥≥≤0,0..min b x b Ax t s cx 其中初始可行基为松弛变量对应的列组成. 对于一般标准线性规划问题:⎩⎨⎧≥≥=0,0..min b x b Ax t s cx１．求解上述一般标准线性规划的单纯形算法（修正）步骤如下：对于一般的标准形式线性规划问题(求极小问题),首先给定一个初始基本可行解。
设初始基为B,然后执行如下步骤： (1).解B Bx b=,求得1B x B b -=，0,N B B x f c x ==令计算目标函数值 1(1,2,...,)i m B b i -=i 以b 记的第个分量(2).计算单纯形乘子w,BwB C =,得到1B wC B -=,对于非基变量，计算判别数1i i i B i i z c c B p c σ-=-=-,可直接计算σ=1B A c c B --令max{}k i Rσσ∈=,R 为非基变量集合若判别数0k σ≤ ,则得到一个最优基本可行解，运算结束；否则，转到下一步 (3).解k kBy p =,得到1k k y B p -=；若0k y ≤,即k y 的每个分量均非正数， 则停止计算，问题不存在有限最优解，否则，进行步骤(4).确定下标r,使{}:0min ,0t rrktktk b b tk y y t y y >=>且rB x 为离基变量,,r k B x p k 为进基变量,用p 替换得到新的基矩阵B,还回步骤(1);２、计算框图为：图1 3．计算程序(Matlab)：A=input('A=');b=input('b=');c=input('c=');format rat%可以让结果用分数输出[m,n]=size(A);E=1:m;E=E';F=n-m+1:n;F=F';D=[E,F]; %创建一个一一映射，为了结果能够标准输出X=zeros(1,n); %初始化Xif(n<m) %判断是否为标准型fprintf('不符合要求需引入松弛变量')flag=0;elseflag=1;B=A(:,n-m+1:n); %找基矩阵cB=c(n-m+1:n); %基矩阵对应目标值的cwhile flagw=cB/B; %计算单纯形乘子，cB/B=cB*inv(B),用cB/B的目的是，为了提高运行速度。
单纯形法MATLAB程序
单纯形法（Mat lab程序）%%单纯形法（Mat lab程序）a= input (' input the major matrix A '); b=input (' input the matrix b '); n=input C input the judgement ');%%为计数器(确定循环次数)萨0;while g<40%%确定非负alength=max(size(n));blength二max(size(b));m=0;for i=l:alength辻n(i)〉=0m二m+1;endend;if m==alengthx=b;breakend;%%找Ks二min(n);for i=l:alengthif n(i) ==sk二i;breakend;end;%%a[i,k]的非负性m=0;for i=l:blengthif a(i, k)<0m二m+1;end;end;if m==blengthdisp('x does not exit');judge二1；breakend;%%找L确定主元cc=100000;for i=l:blengthif a (i, k) >0if(b(i)/a(i, k))<cccc=b(i)/a(i, k);endend end; for i=l:blengthif a(i, k)~=0if (b(i)/a(i, k))==cc1二i；breakendend end; %%计算,a 标准化zu=a(l, k); aa=a; for i=l:1-1 for j=l:alength aa(i, j)=a(i, j)-a(l, j)*a(i, k)/a(l, k);end end; for i=l+l:blengthfor j=l :alength aa(i, j)=a(i, j)-a(l, j)*a(i, k)/a(l, k);end end; for j=l:alengthaa(l, j)=a(l, j)/zu; end;%%b 勺判别bb=b; bb(l)=b(l)/zu;for i=l: 1~1 bb(i)=b(i)~b⑴*a(i, k)/a(l, k);end;for i=l+l:blength bb(i)二b(i)-b(l)*a(i, k)/a(l, k);end;b二bb; %%确定判别数tt 二n;for j=l:alength11 (j) =n(j)-a(1, j)*n(k)/a(1, k) ; end; n=tt;a=aa;%%显示单纯形表sa sa二[b' aa;0 n];dispC单纯表示例’)；disp(g+1);disp(sa);g二g+l;judge=2;end;if judge==2q二0; result=zeros (alength, 2); for j=l+q:alengthif n(j)=0 t=a(:, j) ; zu=find( t) ; resu lt( j, l)=j ; result (j, 2)=x(zu) ; q 二q+1 ；endif n(j)>0 result(j,l)=q+l； q=q+l；endend;dispC最优解’)；disp (result);dispC循环次数')；end。
实验二MATLAB编程单纯形法求解
北京联合大学实验报告项目名称: 运筹学专题实验报告学院: 自动化专业: 物流工程班级: 1201B 学号:21姓名: 管水城成绩: 2015 年 5 月 6 日实验二:MATLAB 编程单纯形法求解一、实验目的:(1)使学生在程序设计方面得到进一步的训练;,掌握Matlab (C 或VB)语言进行程序设计中一些常用方法。
(2)使学生对线性规划的单纯形法有更深的理解、二、实验用仪器设备、器材或软件环境计算机, Matlab R2006三、算法步骤、计算框图、计算程序等本实验主要编写如下线性规划问题的计算程序:⎩⎨⎧≥≥≤0,0..min b x b Ax t s cx其中初始可行基为松弛变量对应的列组成、对于一般标准线性规划问题:⎩⎨⎧≥≥=0,0..min b x b Ax t s cx１.求解上述一般标准线性规划的单纯形算法(修正)步骤如下:对于一般的标准形式线性规划问题(求极小问题),首先给定一个初始基本可行解。
设初始基为B,然后执行如下步骤:(1)、解B Bx b =,求得1B x B b -=,0,N B B x f c x ==令计算目标函数值 1(1,2,...,)i m B b i -=i 以b 记的第个分量(2)、计算单纯形乘子w, BwB C =,得到1B w C B -=,对于非基变量,计算判别数1i i i B i i z c c B p c σ-=-=-,可直接计算σ=1B A c c B --令 max{}k i R σσ∈=,R 为非基变量集合若判别数0k σ≤ ,则得到一个最优基本可行解,运算结束;否则,转到下一步(3)、解k k By p =,得到1k k y B p -=;若0k y ≤,即k y 的每个分量均非正数, 则停止计算,问题不存在有限最优解,否则,进行步骤(4)、确定下标r,使 {}:0min ,0t rrk tk tk b b tk y y t y y >=>且r B x 为离基变量,,r k B x p k 为进基变量,用p 替换得到新的基矩阵B,还回步骤(1);２、计算框图为:图1 3.计算程序(Matlab):A=input('A=');b=input('b=');c=input('c=');format rat%可以让结果用分数输出[m,n]=size(A);E=1:m;E=E';F=n-m+1:n;F=F';D=[E,F]; %创建一个一一映射,为了结果能够标准输出X=zeros(1,n); %初始化Xif(n<m) %判断就是否为标准型fprintf('不符合要求需引入松弛变量')flag=0;elseflag=1;B=A(:,n-m+1:n); %找基矩阵cB=c(n-m+1:n); %基矩阵对应目标值的cwhile flagw=cB/B; %计算单纯形乘子,cB/B=cB*inv(B),用cB/B的目的就是,为了提高运行速度。
MATLAB线性规划非线性规划
解 编写M文件xxgh1.m如下： c=[-0.4 -0.28 -0.32 -0.72 -0.64 -0.6]; A=[0.01 0.01 0.01 0.03 0.03 0.03;0.02 0 0 0.05 0 0;0 0.02 0 0 0.05 0;0 0 0.03 0 0 0.08]; b=[850;700;100;900];
返 回
解答
线性规划模型的一般形式
目标函数和所有的约束条件都是设计变量 的线性函数.
min u ci xi
i 1
n
n aik xk bi , i 1, 2,..., n. s.t. k 1 x 0, i 1, 2,..., n. i
矩阵形式： min u cx Ax b s.t. vlb x vub
Aeq=[]; beq=[];
vlb=[0;0;0;0;0;0]; vub=[];
To MATLAB (xxgh1)
[x,fval]=linprog(c,A,b,Aeq,beq,vlb,vub)
例2
min z 6x1 3x2 4x3 s.t. x1 x2 x3 120 x1 30 0 x2 50 x3 20
编写M文件xxgh3.m如下: f = [13 9 10 11 12 8]; A = [0.4 1.1 1 0 0 0 0 0 0 0.5 1.2 1.3]; b = [800; 900]; Aeq=[1 0 0 1 0 0 0 1 0 0 1 0 0 0 1 0 0 1]; To MATLAB (xxgh3) beq=[400 600 500]; vlb = zeros(6,1); vub=[]; [x,fval] = linprog(f,A,b,Aeq,beq,vlb,vub)
单纯形法的matlab编程
单纯形法的matlab实现首先输入三个值系数矩阵A目标函数系数行向量C列向量b根据大M法进行扩列A,C,b.使得行数不变，列数增加M 进行的到基向量的坐标，非基变量的坐Cb,Cn,Xb,Xn,此时的值便是典式，不在需要进行进一步化简，只需求解检验变量delta的值迭代过程输入上一步得到A,C,b,Cb,Cn,Xb,Xn,输出值为最优解为X，得到目标函数的最优解Z的值迭代循化用while循环当找到解时结束循环break或者当发现循化结果没有最优解时跳出循环，这里涉及两个判断，两个判断量初始值都可以写在循环外，两者的值共同决定循环的执行与否循化最开始进行判断初始可行解是否为最最优解，若是直接跳出循化，若上面的判断不成立，接下来进行下一个判断，若不符合进行下面入基和出基变量的选值入基和出基变量的循化是两次循化，第一次找到k的值，第二次根据上一次的k找r的值注意因为值有约束，而且是找函数最小值，需要对这个列向量进行变换一下将小于等于0的都变成无穷大，接下来进形下一次的循化，进而找到转轴元将A,b,delta合成一个新的矩阵，进行旋转变化，得到值后反变回相应的值，接下来需要对Xb，Xn的值进行交换这个步骤要两个循环，第一个循化对Ark的所在行进行变化，接下来进行对整个矩阵进行行变换，包括两种情况，两次循化嵌套分别是r==1时和r~=1的时候建立总体X的坐标列向量发生交换时出基变量找Xb,入基变量从X中找有先后顺序先解决Xn的变化。
在解决Xb的值直接解决基变量其他为0A=input('输入系数矩阵\n');b=input('输入列向量b\n');C=input('输入目标函数行向量\n');M=5200;global m;global n;global X;[m,n]=size(A);I=eye(m);A=[A,I];Xb=[];Xn=[];for i=1:mC(i+n)=-M;Xb(i)=n+i;endXb=Xb';Cb=C(1,n+1:n+m);for i=1:nXn(i)=i;endXn=Xn';X=[Xn;Xb];[m,n]=size(A);diedai(A,C,b,Cb,Xb);function[Z]=diedai(A,C,b,Cb,Xb)delta=C-Cb*A;global m;global n;global X;while1s2=0;s1=0;for j=1:nif delta(j)>0s1=1;for i=1:mif A(i,j)>0s2=1;endendendendif s1==0disp('目标函数最优解')Z=Cb*b;disp(Z)disp('基变量为');[Xb,index]=sort(Xb);disp(Xb)b=b(index);disp('基可行解为');disp(b)break;endif s2==0disp('目标函数无界，无最优解');break;end[~,k]=max(delta);p=A(:,k);zhuan=[];for i=1:mzhuan(i)=b(i)/p(i);if zhuan(i)<=0zhuan(i)=inf;endend[~,r]=min(zhuan);b(m+1)=0;Z=[A;delta];Z=[Z,b];z=Z;ark=A(r,k);for j=1:n+1Z(r,j)=Z(r,j)/ark;endif r==1for i=2:m+1for j=1:n+1Z(i,j)=Z(i,j)-z(i,k)*Z(r,j);endendelse for i=[1:r-1,r+1:m+1]for j=1:n+1Z(i,j)=Z(i,j)-z(i,k)*Z(r,j);endendendA=Z(1:m,1:n);delta=Z(m+1,1:n);b=Z(1:m,n+1);Cb(r)=C(k);Xb(r)=X(k);endend。
线性规划单纯形法matlab解法
%单纯形法matlab程序-ssimplex% 求解标准型线性规划:max c*x; . A*x=b; x>=0% 本函数中的A是单纯初始表，包括:最后一行是初始的检验数，最后一列是资源向量b% N是初始的基变量的下标% 输出变量sol是最优解, 其中松弛变量（或剩余变量）可能不为0% 输出变量val是最优目标值，kk是迭代次数% 例：max 2*x1+3*x2% . x1+2*x2<=8% 4*x1<=16% 4*x2<=12% x1,x2>=0% 加入松驰变量，化为标准型，得到% A=[1 2 1 0 0 8;% 4 0 0 1 0 16;% 0 4 0 0 1 12;% 2 3 0 0 0 0];% N=[3 4 5];% [sol,val,kk]=ssimplex(A,N)% 然后执行 [sol,val,kk]=ssimplex(A,N)就可以了。
function [sol,val,kk]=ssimplex(A,N)[mA,nA]=size(A);kk=0; % 迭代次数flag=1;while flagkk=kk+1;if A(mA,:)<=0 % 已找到最优解flag=0;sol=zeros(1,nA-1);for i=1:mA-1sol(N(i))=A(i,nA);endval=-A(mA,nA);elsefor i=1:nA-1if A(mA,i)>0&A(1:mA-1,i)<=0 % 问题有无界解disp('have infinite solution!');flag=0;break;endendif flag % 还不是最优表，进行转轴运算temp=0;for i=1:nA-1if A(mA,i)>temptemp=A(mA,i);inb=i; % 进基变量的下标endendsita=zeros(1,mA-1);for i=1:mA-1if A(i,inb)>0sita(i)=A(i,nA)/A(i,inb);endendtemp=inf;for i=1:mA-1if sita(i)>0&sita(i)<temptemp=sita(i);outb=i; % 出基变量下标endend% 以下更新Nfor i=1:mA-1if i==outbN(i)=inb;endend% 以下进行转轴运算A(outb,:)=A(outb,:)/A(outb,inb);for i=1:mAif i~=outbA(i,:)=A(i,:)-A(outb,:)*A( i,inb);EndEndEndEndend。
单纯形法的MATLAB代码
单纯形法的MATLAB代码% 求解标准型线性规划:max c*x; s.t. A*x=b;x>=0%A1是标准系数矩阵及最后一列是资源向量，C是目标函数的系数向量% N是(初始的)基变量的下标%M=10000 人工变量系数% 本函数中的A是单纯形表，包括:最后一行是初始的检验数，最后一列是资源向量b%c1是基变量系数%输出变量sol是最优解%输出变量val是最优值，k是迭代次数%flag1的值代表有无最优解，0无界解，1无可行解，2无穷多解，3唯一最优解function [sol,val,k,flag1]=ssimplex(A1,C,N)M=10000;[mA1,nA1]=size(A1);C1=[C,0];val=zeros(1,length(C));for i=1:length(N)c1(i)=C1(N(i));endfor i=1:nA1a(i)=C1(i)-c1*A1(:,i);%计算初始检验数endA=[A1;a]; %构造初始单纯形表[mA,nA]=size(A);k=0; % 迭代次数flag=1;while flagfor i=1:(nA-1)if A(mA,i)<=0flag=0;elseflag=1;break;endendif flag==0 % 已找到最优解val1=A(1:(mA-1),nA)';for i=1:length(N)if (val1(i)~=0&&abs(C(N(i)))==M)disp('无可行解');sol=inf;val=inf;flag3=0;flag1=1;break;elseflag3=1;endendif flag3if length(find(A(mA,1:(nA-1))==0))>length(N) disp('存在无穷多最优解');flag1=2;elsedisp('存在最优解');flag1=3;endsol=c1*val1';endelseif flag==1for j=1:(mA-1)if A(j,i)<=0flag2=1;elseflag2=0;break;endendif flag==1&&flag2==1disp('此线性规划问题存在无界解');sol=inf;val=inf;flag1=0;flag=0; %跳出while循环break;endmaxq=max(A(mA,1:(nA-1)));[m,nb]=find(A(mA,:)==maxq); %确定入基变量的纵坐标for s=1:(mA-1)if A(s,nb)>0temp(s)=A(s,nA)/A(s,nb);elsetemp(s)=10000;endendk=k+1;mino=min(temp);[n,mb]=find(temp==mino); %确定入基变量的横坐标if length(mb)>1mb=mb(1);endab=A(mb,nb);A2=A;for i=1:(mA-1)for j=1:nAif i==mbA(mb,j)=A2(mb,j)/ab;elseA(i,j)=A2(i,j)-A2(i,nb)*(A2(mb,j)/ab); endendendfor i=1:length(N)if i==mbN(i)=nb;endendfor i=1:length(N)c1(i)=C(N(i));endfor i=1:nAA(mA,i)=C1(i)-c1*A(1:(mA-1),i); endendendif sol~=inffor i=1:length(C)for j=1:length(N)if i==N(j) val(i)=val1(j); endendendend。
单纯形法的matlab实现(极小化问题)
实验报告实验题目 :单纯形法的matlab 实现学生姓名 :学号:实验时间 :2013-4-15一．实验名称 :单纯形法的MATLAB实现二．实验目的及要求 :1.了解单纯形算法的原理及其matlab 实现 .2.运用 MATLAB 编辑单纯形法程序解决线性规划的极小化问题, 求出最优解及目标函数值 .三．实验内容 :1.单纯形方法原理 :单纯形方法的基本思想 , 是从一个基本可行解出发 , 求一个使目标函数值有所改善的基本可行解 ; 通过不断改进基本可行解 , 力图达到最优基本可行解 . 对问题min f def cxs.t.Ax b,x0.其中 A 是一个 m× n 矩阵 , 且秩为 m, c 为n维行向量,x 为n维列向量,b为m维非负列向量 . 符号“def”表示右端的表达式是左端的定义式,即目标函数f 的具体形式就是cx .记A ( p1, p2 ,..., p n )令 A =(B,N), B为基矩阵, N为非基矩阵,设x( 0) B -1b是基本可行解 , 在x(0 )处的目标函数值f 0 cx( 0)(c B,c N )B-1b-1c B B b , 0其中 c B是c中与基变量对应的分量组成的m 维行向量 ;c N是c中与非基变量对应的分量组成的 n-m 维行向量 .现由基本可行解 x(0)出发求解一个改进的基本可行解.设 x x B是任一可行解 , 则由Ax b得到x Nx B B-1b B-1N x N,在点 x 处的目标函数值f cx (c B , c N )x Bf 0(z j c j ) x j, x N j R其中 R 是非基变量下标集,z j c B B-1 p j.2.单纯形方法计算步骤 :首先给定一个初始基本可行解, 设初始基为 B, 然后执行下列主要步骤:）解 Bx B b ,求得 x B B 1b _（1 b ,令 x N0 ,计算目标函数值 f c B x B.（2）求单纯形乘子w ,解wB c B,得到 w c B B 1. 对于所有非基变量,计算判别数z j - c j w j p j c j.令z k - c k max{z j- c j }j R.若 z k - c k0 ,则对于所有非基变量z j- c j0, 对应基变量的判别数总是为零, 因此停止计算 , 现行基本可行解是最优解. 否则 , 进行下一步 .（3）解 By k pk ,得到y k B 1 p k,若 y k0, 即y k的每个分量均非正数,则停止计算 ,问题不存在有限最优解. 否则进行步骤（4） .（4）确定下标r, 使x k =b r min bi y ik 0 ,y rk y ikxB r 为离基变量 ,x k 为进基变量 . 用 p k 替换 p B r , 得到新的基矩阵 B, 返回步骤（ 1） .3. 单纯形方法表格形式 :表 3.1.1fx Bx N 右端x BI mB -1NB -1bfc B B -1 N - c Nc B B -1b1表 3.1.2（ 3.1.1 略去左端列后的详表）xB 1xBrxBmxjxkx B 110 0y 1 jy 1k_b 1x B r1 0y rjy rk_b rxB m1ymjymk_b mf0 0 z j - c jz k - c kc B b假设 bB -1b0 , 由上表得 x B b, x N 0 .若 c B B -1 N - c N 0 , 则现行基本可行解是最优解 .若 c B B -1N - c N0, 则用主元消去法求改进的基本可行解.先 根 据z k - c k max{z j - c j }b rb iy ik 0 找主行 , 主元为 y rk , 然后j R选择主列 , 再根据 minyrkyik进行主元消去 , 得到新单纯形表. 表的最后一行是判别数和函数目标值.四．实验流程图及其 MATLAB 实现 :1.流程图 :开始初始基本可行解B_解 Bx B b ,求得 x B B 1b b ,令 x N0 ,计算目标函数值f c B x B 求单纯形乘子 w ,解wB c B,得到 w c B B 1. 对于所有非基变量 , 计算判别数z j- c j w j p j c j.令z k - c k max{z j - c j}j RNz k - c k0Y解Bykpk,得到y k B1p k现行基本可行解是最优解NYy k0确定下标 r, 使 x k =b r min bi y ik 0赋 b i以正的大值 Ny rk yikyikN Ymin=N问题不存在有限最优解x B r为离基变量,x k为进基变量.用p k替换 p B r,得到新的基矩阵B2.代码及数值算例 :(1) 程序源代码 :function[x,f]=DCmin(c,A,b,AR,y0,d)%x: 最优解%f: 目标函数最优值%c: 目标函数系数向量%A: 系数矩阵%b: m 维列向量%AR: 松弛变量系数矩阵%y0: 基矩阵初始向量%d: 补充向量（非目标系数向量, 为一零向量）N=10000;B=[A,AR,b];[m,n]=size(B);C=[c,d];y=y0;x=zeros(1,length(c));for k=1:Nk;z=B(:,end);%右端for j=1:n-1t(j)=y*B(:,j)-C(j);%检验数endt;f=y*z;%%======== 选取主元 ==========%%%---------选取主列---------% [alpha,q]=max(t);q;W(k)=q;%x下标矩阵%-------------------------%%--------选取主元----------%for p=1:mif B(p,q)<=0r(p)=N;else r(p)=z(p)/B(p,q);endend[beta,p]=min(r);p;y(p)=C(q);%-------------------------%%%==========================%% B(p,:)=B(p,:)/B(p,q);for i=1:mif i~=pB(i,:)=B(i,:)-B(p,:)*B(i,q);endendif max(t)<=0break ;end B; end%++++++++++++++++++++++++++++++++++++++% Z=B(:,end);if length(x(W))~=length(Z)x=char(' NONE' ); f=char( ' NONE');disp( '不存在有限最优解 ' );else x(W)=Z'; end(2) 数值算例 :例 3.1.2 用单纯形方法解下列问题min x 1 - 2x 2 x 3s.t.x 1 x 2 -2x 3 x 4 ，102x 1 - x 2 4x 3，8- x 1 2x 2 - x 3 ，4x j ， ，，，0 j 1 2 3 4.引进松弛变量 x 5 , x 6 , 问题标准化 :min x 1 - 2x 2 x 3s.t.x 1 x 2-2x 3x 4，102x 1 - x 2 4x 3 x 5 ，8 - x 1 2x 2 - x 3 x 6 4.x j ， ，，，，，0 j 123456.(i) 输出命令 :>>c=[1 -2 1];A=[1 1 -2 1;2 -1 4 0;-1 2 -4 0];b=[10;8;4];AR=[0 0;1 0;0 1];y0=[0 0 0];d=[0 0 0]; >>[x,f]=DCmin(c,A,b,AR,y0,d)(ii)运行结果 :B =11-2100102-140108-12-40014k =1t =-12-1000f =B =1.500000 1.00000-0.50008.00001.500002.00000 1.00000.500010.0000-0.5000 1.0000-2.0000000.5000 2.0000k =2t =00300-1f =-4B =1.500000 1.00000-0.50008.00000.75000 1.000000.50000.2500 5.00001.0000 1.000000 1.0000 1.000012.0000k =3t =-2.2500000-1.5000-1.7500f =-19x =0125f =-19五．总结 :在单纯形法求解过程中 , 每一个基本可行解 x 都以某个经过初等行变换的约束方程组中的单位矩阵为可行基 . 对于极大化的线性规划问题 , 先标准化 , 即将极大化问题变换为极小化问题 :max cx min - cx然后利用单纯形方法求解.六．参考文献 :陈宝林编著《最优化理论与算法》清华大学出版社2005 年 10 月第 2 版。
单纯形法的matlab实现
结果分析：
可以看出只需要90根原料，其中，方案1需要10根，方案2需要50根，方案4需要30根，即可达到要求，此时总剩余废料最小，为16m。
附：
方案
合计
余料
1
2
0
1
2
1
2
0
3
1
1
1
4
1
0
3
0
5
0
3
0
6
0
22Biblioteka 7013
8
0
0
4
6
实验的启示：
通过本次实验加深了我对单纯形法的进一步理解，利用matlab程序求解线性规划问题，在生活中有不少问题都是线性规划问题。例如本次实验中的钢材下料问题，学习用已学习的知识解决实际问题是我们的最终目标。
在上述算法中，当存在不止一个简约价值系数 时，选取最负的 的指标为 ，并以 作为入基变量。
Matlab计算程序：
Function[x,f]=zuiyouhua(A,b,c)
Size(A)=[m,n];
i=n+1:n+m;
N=1:n;
B=eye(m,m);
xb=b’;
xn=zeros(m,1);
f1=0;
在matlab的输入区域输入：
A=[2,1,1,1,0,0,0,0;0,2,1,0,3,2,1,0;1,0,1,3,0,2,3,4];
b=[100,100,100]’;
c=[,,,0,,,,];
[x,f]=zuiyouhua(A,b,c)
Matlab输出内容：
x=10 50 0 30 0 0 0 0
w=zeros(1,m);
z=-c;
MATLAB优化工具箱--线性规划-非线性规划
fmincon函数求解形如下面的有约束非线性规 划模型
一般形式：
min f ( X ) s.t. AX b
Aeq X beq l X u c(X ) 0 ceq ( X ) 0
Matlab求解有约束非线性最小化 1.约束中可以有等式约束 2.可以含线性、非线性约束均可
数学实验
输入参数语法：
靠着楼房建有一个温室，温室高10英尺，延伸进 花园7英尺。清洁工要打扫温室上方的楼房的窗户。 他只有借助于梯子，一头放在花园中，一头靠在 楼房的墙上，攀援上去进行工作。他只有一架20 米长的梯子，你认为他能否成功？能满足要求的 梯子的最小长度是多少？请就以上问题建立数学 模型，并编程求解。
18
提示：
19
算失败
7
例题的求解程序
模型： max 6x1+4x2 s.t. 2x1+5x2 ≤100
4x1+2x2 ≤120
Matlab求解程序：
A=[2 5;4 2]; b=[100 120]; f=-[6 4]; [optx ,funvalue,exitflag]=linprog(f,A,b,[],[],[0
lb,ub,nonlcon)
16
学习小结
最优化问题建模的关键是先要确定三要素，再转
化为数学表达式（数学模型）。
学习中既要初步掌握最优化问题的建模步骤，也
要善于运用Matlab的优化工具箱求解优化模型。
有些模型可以采用多个Matlab函数求解，可以比
较结果，加深认识。
17
思考题
一幢楼房的后面是一个很大的花园。在花园中紧
3
Matlab求解线性规划模型 函数linprog
求解下列形式的线性规划模型：
- 1、下载文档前请自行甄别文档内容的完整性，平台不提供额外的编辑、内容补充、找答案等附加服务。
- 2、"仅部分预览"的文档,不可在线预览部分如存在完整性等问题,可反馈申请退款(可完整预览的文档不适用该条件!)。
- 3、如文档侵犯您的权益，请联系客服反馈,我们会尽快为您处理(人工客服工作时间：9:00-18:30)。
第 4 页
单纯形法程序整理.txt for i=1:h%判断有无相等的最小值及最小值的重数 if(TableK(i,l+1)/TableK(i,s)==t) R(k)=i;%记录下标 k=k+1; end end %移动比较列 flag=0;%标记跳出循环位置,判断标记跳出多重循环 for j=1:l if(k>2)%有两个及以上最小值时，移动分子的列位 k=k-1; for i=1:k%检查比值是否为负 if(TableK(R(i),j)/TableK(R(i),s)>=0) t=TableK(R(i),j)/TableK(R(i),s); break; end if(i==k)%若比值全为负，则寻找下一列 flag=1; break; end end if(j==6) asdfa=1; end if(flag==1)%判断是否需要移动列位 flag=0;%重置标记 k=k+1;%恢复k continue; end for i=1:k if(TableK(R(i),j)/TableK(R(i),s)<0)%排除负元素 continue; end if(TableK(R(i),j)/TableK(R(i),s)<=t) t=TableK(R(i),j)/TableK(R(i),s);%记录判别数最小值 r=i;%记录下标 end end end k=1; for i=1:h%判断有无相等的最小值及最小值的重数 if(TableK(i,s)==0) continue;%跳过为0的项 end if(TableK(i,j)/TableK(i,s)==t) R(k)=i;%记录下标 k=k+1; end end if(k==2)%若已无重复最小值，退出循环 break; end end end %————————————————% 第 3 页
单纯形法程序整理.txt [hB,lB]=size(B);%测算B的行、列数 for j=1:hB if(A(B(j,1),B(j,2))~=1)%如果主元不为1，将对应行除以主元化为1 A(B(j,1),:)=A(B(j,1),:)/A(B(j,1),B(j,2)); end for i=1:hA%将其他元化为零 if(i==B(j,1)) continue; end A(i,:)=A(i,:)-A(B(j,1),:)*A(i,B(j,2)); end end end %————————————————% %% 子程序-选主元 function [r,s]=SubFc_XuanZhuYuan(TableK,h,l) %选主元素 % 选出表矩阵TableK中的主元素 %% 选列主元 t=0; for i=1:l if(TableK(h+1,i)<=t) t=TableK(h+1,i);%记录判别数最小值 s=i;%记录下标 end end if(t==0) error('判别数全为正'); end %% 选行主元 for i=1:h%检查元素是否为负，并赋予t初值 if(TableK(i,s)>0) t=TableK(i,l+1)/TableK(i,s); break; end if(i==h)%若元素全为负，则报错 error('主列元素全为非正，规划问题有无界解！'); end end %k=1; for i=1:h %若未知数出现负值，则是初始迭代标准性表产生的，应更换初始主元； if(TableK(i,l+1)<0) error('b向量出现负值，更换初始主元'); end if(TableK(i,s)<=0)%排除A表负元素 continue; end if(TableK(i,l+1)/TableK(i,s)<=t) t=TableK(i,l+1)/TableK(i,s);%记录判别数最小值 r=i;%记录下标 end end %补充判据 k=1; 第 2 页
单纯形法程序整理.txt %% 子程序-迭代 function [ TableK,BK ] = SubFc_DieDai( TableK,BK,h,l ) %迭代一步 % 给出上一步得到的TableK，和这一步的主元位置BK，得到这一步的TableK和下一步的BK %% 初等行变换 将主元对应列化成单位列向量 [TableK]=SubFc_ChuDengHangBianHuan(TableK,BK); disp('当前标准性表：'); disp(TableK); %检查判别数是否全为正 t=0; for i=1:l if(TableK(h+1,i)<=t) t=TableK(h+1,i);%记录判别数最小值 s=i;%记录下标 end end if(t==0) return; end % TableK=function[Table,B]; % Tபைடு நூலகம்bleK=Table; %% 选主元素 [r,s]= SubFc_XuanZhuYuan(TableK,h,l); disp('更换主元位置：行/列'); disp([r,s]); %% 更换主元 BK(r,:)=[r,s];%更换主元 end
单纯形法程序整理.txt %% 单纯形法 %% 主程序 %————————————————% %% 单纯形法 % 求解线性规划问题 % 编程：范中允 %% 标准矩阵A % 输入标准形式的约束条件系数矩阵A、右端列向量b、约束条件未知数对应目标函数的系数 c； % 有需要可以指定迭代的初始单位阵（一组基本可行解）； A=[1,1,2,1,0; 1,4,-1,0,1; ]; c=[-2,-1,1,0,0]; b=[6;4]; Table=[A,b;c,0];%初始表 [h,l]=size(A); %% 给出初始基本主元 B=[1,1;2,4]; BK=B; TableK=Table; %% 迭代 for j=1:100 [TableK,BK]=SubFc_DieDai(TableK,BK,h,l); t=0; disp(['第',num2str(j),'次迭代']); for i=1:l if(TableK(h+1,i)<=t) t=TableK(h+1,i);%记录判别数最小值 s=i;%记录下标 end end if(t==0) disp('迭代结束'); break; end end if(j==100) sprintf('迭代100次未收敛'); end disp('最终的标准型表：'); disp(TableK); disp('最优解的目标函数：'); disp(-TableK(h+1,l+1)); disp('基本可行解下标：'); disp(BK(:,2)); disp('基本可行解系数：'); disp(TableK(1:h,l+1)); %————————————————% %% 子程序-初等行变换 function [ A ] = SubFc_ChuDengHangBianHuan( A,B ) %初等行变换 % 将矩阵指定的主元对应列，由初等行变换，化为主元为1其他元为0 % A为需要变换的矩阵，B为用于指定主元位置的矩阵 [hA,lA]=size(A);%测算A的行、列数 第 1 页