单纯形法、线性规划实践报告

一、线性规划一一单纯形法程序设计1实验目的：使学生在程序设计方面得到进一步的训练;，掌握Matlab（C或VB）语言进行程序设计中一些常用方法。使学生对线性规划的单纯形法有更深的理解.2问题陈述本实验主要编写一般线性规划问题的计算程序:MinCTxAixSg;Azxb：；Asx=b3xNO,bib2b30引入松弛变量将其化为一般标准型线性规划问题：MinCTxAx=b;x0,b0A为m*n的矩阵，有m个约束，n个变量。求解上述线性规划采用单纯形算法，初始可行基由引入的m个人工变量对应的单位阵组成，并采用大M算法3算法描述Xi（1）将引入的人工变量对应的单位阵作为初始可行基，则原线性规划问题构造出下面的新线性规划问题：2内+%2丿/=1为勺+心+产0丿=1加St弋y=ix0,7=1”昇？+1刃+勿T厂一（2）通过判别数计算公式a=CCBA可求出n+m个变量的判别数，若全部判别数则得到一个最优基本可行解，运算结束；否则，转到下一步（3）找出判别数为负的最小判别数，其对应的变量为入基变量，记下标为k,然后看其对应的列向量Pk,若Pk中的所有元aik0,则原线性规划无最优解，否则，转下一步（4）计算则Xi为离基变量，然后对A进行初等变换，计算aik（z工0；$=/?,-化喩（心/）;（5）用入基变量与出基变量对应的列向量、判别、对应的系数均对换，则可用计算机编程循坏以上步骤直至得出结果4计算程序（matlab）fundionf=linpro(Aeq,b?cO)nijn=size(Aeq);E=eye(m):%引入人工变里对应的列向里M=10000;cl=ones(uJR).%人工变里对应的系数c=cOj.cl:A=AeqE;%构洼的新柜阵sigma=c-cl*A;策计真判别数whileleneth(sigma(sigma0%所有判别数全为正时为最优餡z3k=min(sigma):%找入基变里ifA(:k)=0fprintf此i可题无最优SS!rt):brea0：theta(j)=b(j)./A(k):.壬s=min(theta):%找出基变里endfori=l:m%前后两次迭代系数矩蹲的对应关系ifi=sb(s)=b(s)./A(ssk):A(s?:)=A(s:)./A(s?k):elseb(i)=b(i)-(Ak)/A(sk).*b(s);A(i,:)=A(iJ:)-(A(i5k)./A(s5k).*A(s,:):endendsigm.a=sigma-(sigma(k)/A(s,.k).*A(s,.:):%前后两次迭代判别救的关系sigma(k)=sigjRa(s+n);signva(s-n)=0:d=A(:,k):弼对换入基和出基变里对应的列向里A(:k)=A(:2s-t-n):A(:s+n)=d;r=c(s+n);滋对换其判别数c(s-Ha)=c(k):c(k)=r:fJR=O:pfori=l:m%计算最优值fm二fm十b(i)*cCn十i):enddisp(fm);滋输出最优值nd程序保存为文件5算例验证mniz=一015X-0.1兀2-0.08x3-0.12x4Xj兀？兀3兀4-x2-x3+x4O,J=1,2,.,4在窗口输入:Aeq=l-1-1-110;01J1011100;b=0;0;l;c0=00;linpro(Aeq,b,cO)1000010000说明通过三次迭代找到最优解为.用Matlab求解线性规划的命令linprog的计算结果:f=;A=1-1-1-10-1-11;b=0;0;Aeq=l111);beq=l;lb=zeros(4,l)；然后调用linprog函数：x,fval=linprog(f,AbAeq,beqb)：fval=最优值为，与上面的结果一致，说明程序正确。二、单变量在单峰区间上的极小点黄金分割法1.问题陈述设f(x)为区间a,b上的下单峰函数，在区间内取两点xl,x2,使每次迭代都把区间缩短率定为，若f(xl)f(x2)/则X气xl.b这样区间就减小了，直至找到最优解。2.算法描述(1)令x2=a+(b-a)/f2=f(x2);(2)令xl=a+(b-a)/fl=f(xl);若ba,则找到最优解X性(a+b)2否则转下一步;(4)若flf2,ma=xl/xl=x2,fl=f2,转第三步;(5)令x2=a+(b-a)/f2=f(x2)/转第三步。3程序代码(C+)以函数f(X)=x2-x+2为例编写:#includeusingnamespacestd;doublea,b,t;doublef(doublex)doubled;d=x*x-x+2;returnd;doubleminf(doublea,doubleb)doublexl,x2,fl,f2;xl=a+*(b-a);x2=a+*(b-a);fl=f(xl);f2=f(x2);while(b-ale-4)b=x2;x2=xl;f2=fl;xl=a+*(b-a);fl=f(xl)；elsea=xl;xl=x2;fl=f2;x2=a+*(b-a);f2=f(x2);elsea=xl;b=x2;xl=a+*(b-a);x2=a+*(b-a);fl=f(xl);f2=f(x2);t=(a+b)/;returnt;voidmain()doublek;cout请输入单峰区间的端点a和b：cinab;k=minf(a,b);cout最优值为:,k,M优解为：f(k)endl;3.算例验证x2-x+2在区间13上的最小值，程序运算如卞:请输入单峰区间的端点a和b；-13最优值为：0-500016最优解为：1-75Pressanykeytocontinue在matlab命令窗口输入如下:fun=xA2-x+2;x,fval=fminbnd(fun,-l,3)运算结果如下：fval=两者运行结果完全一致，说明程序正确。1. 三、运用非线性规划建模的实例问题描述：试设计一压缩圆柱螺旋弹簧，要求其质量最小。弹簧材料为65Mn,最人工作载荷Pmax=40N,最小工作载荷为0,载荷变化频率fr=25Hz,弹簧寿命为104h,弹簧钢丝直径d的取值范闱为l-4mm,中径D2的取值范围为10-30mm,工作圈数n不应小于圈，弹簧旋绕比C不应小于4,弹簧一端固定，一端自由，工作温度为50C,弹簧变形量不小于10mm。模型建立本题的优化目标是使弹簧质量最小，圆柱螺旋弹簧的质量可以表示为兀2M/(n+w2)2式中，7-弹簧材料的密度，对于钢材=X10-6kg/mm3；n工作圈数：n2死圈数，常取n2=，现取n2=2；D2弹簧中径，mm；d弹簧钢丝直径，mm。将已知参数代入公式，进行整理以后得到问题的目标函数为/(X)=M=0.192457x10-4(x2+2)xfx3根据弹簧性能和结构上的要求，可写出问题的约束条件：强度条件gi(X)=350-163.0兀严686二o刚度条件g2(X)=0.4xKF?xx2x-10.00稳定性条件g3(X)=3.7*3一(2+1-0.44xIO-2x4x2xl0不发生共振现彖，要求g4(X)=0.356x10抵1巧亍2375X0弹簧旋绕比的限制g5(X)=x3xf1-4.00对d,n,D2的限制1.0J4.0且应取标准值，即，“等。4.5w5010Z)230由上可知，该压缩圆柱螺旋弹簧的优化设计是一个三维的约束优化问题，其数学模型为:nwi/(X)=M=0.192457x10_4(x2+2)xfx3gl(X)=350-163xf286x860g3(X)=3.7兀3(兀2+15)兀1-0.44xIO-2xf4x2xj0g4(X)=0.356xlO62-3750&5(*)=兀3/兀1一40g9(X)=50-x20gioE-lOnOgH(X)=30-x30取初始设计参数为”。)=卩首先编写目标函数的M文件，返回x处的函数值f。functionf=myfun(x)f=*le-4*(x(2)+2)*x(l)A2*x(3);由于约束条件中有非线性约束，所以需要编写一个描述非线性约束条件的M文件：functionczceq=mycon(x)c(l)=350-163*x(l)A*x(3)A;c(2)=*x(l)A(-4)*x妆八3;c(3)=(x(2)+妆(l)+x八(4)妆(2)妆A*x(3);c(4)=*le6*x(l)*x(2)A(-l)*x(3)A(-2);c(5)=4-x(3)/x(l);然后设置线性约束的系数：A=-l001000-1001000-1001；b=-l;4;50;-10;30;卞一步给定初值，给定变量的下限约束，并调用优化过程(磁盘中M文件为xO=;；lb=zeros(3,l);xzfvalzexitflag,output,lambda=fmincon(opt25_3o,xO/A,b/Ib,opt25_3c)计算结果为：fval=exitflag=1output=iterations:3funcCount:16stepsize:1algorithm:medium-scale:SQRQuasi-Newton,line-search1firstorderopt:cgiterations:所以当弹簧钢丝的直径d、工作圈数n及中径D2分别取、和时弹簧质量最小,为克。考虑到实际情况，各参数可分别取、和x=；z=optim253(x)z=故此时弹簧质量仍为克。