北航最优化方法大作业参考

1、11流量工程问题1.1问题重述定义一个有向网络G=(N,E)，其中N是节点集，E是弧集。令A是网络G的点弧关 联矩阵，即NXE阶矩阵， 且第I列与弧里(l,j)对应， 仅第i行元素为1， 第j行元素为-1, 其余元素为0。 再令bm=(bm1,bmN)T，fm=(fm1,fmE)T，则可将等式约束表示成：Afm=bm本算例为一经典TE算例。算例网络有7个节点和13条弧，每条弧的容量是5个单 位。此外有四个需求量均为4个单位的源一目的对，具体的源节点、目的节点信息如图所 示。这里为了简单，省区了未用到的弧。此外，弧上的数字表示弧的编号。此时，c=(5,5,5)1X13)T，4 -400,如一4

2、040,ba =0 电4D切一4 -00n000o00000004根据上述四个约束条件，分别求得四个情况下的最优决策变量图1网络拓扑和流量需求X=(X12,X13,X75)1X13)。in = 1:亠 2.dsA瞒 m 3: 3- 丄山21.27节点算例求解利用Matlab编写对偶单纯形法程序，可求得 最优解为x2*=0 4 0 0 0 0 0 0 0 0 0 0 0T对应的最优值CTX2=201.2.3算例3(b3=0;-4;4;0;0;0;0T)Mi ni mizeJx3Subject toAx3=b3X3=0利用Matlab编写对偶单纯形法程序，可求得：最优解为x3*=4 0 0 0 4

3、 0 0 0 0 0 0 0 0T对应的最优值CTX3=40-111-1-11 1-1-1-11 1-1-1 11.2.1算例1(b仁4;-4;0;0;0;0;0T)转化为线性规划问题：Mini mizecTx1Subject to Ax仁b1x1=0利用Matlab编写对偶单纯形法程序，可求得：最优解为x1*=4 0 0 0 0 0 0 0 0 0 0 0 0T对应的最优值cTx仁201.2.2算例2(b2=4;0;-4;0;0;0;0T)Mini mizecTx2Subject toAx2=b2X2=031.2.4算例4(b4=4;0;0;0;0;0;-4T)Mi ni mizeJx4Su

4、bject toAx4=b4X4=0利用Matlab编写对偶单纯形法程序，可求得：最优解为x4*=4 0 0 4 0 0 0 0 0 4 0 0 0T对应的最优值Jx4=601.3计算结果及结果说明1.3.1算例1(b1=4;-4;0;0;0;0;0T)算例1中，由b1可知，节点2为需求节点，节点1为供给节点，由节点1将信息传 输至节点2的最短路径为弧1。图2算例1最优传输示意图求得的最优解为x1*=4 0 0 0 0 0 0 0 0 0 0 0 0T,即只经过弧1运输4个单位流量， 其余弧无流量。又因为，每条弧的费用均为5，所以最小费用为20。经分析，计算结果合理可信。1.3.2算例2(b2

5、=4;0;-4;0;0;0;0T)算例2中，由b2可知，节点3为需求节点，节点1为供给节点，由节点1将信息传 输至节点2的最短路径为弧2。4r 12L图3算例2最优传输示意图求得的最优解为x2*=0 4 0 0 0 0 0 0 0 0 0 0 0T,即只经过弧2运输4个单位流量, 其余弧无流量。又因为，每条弧的费用均为5，所以最小费用为20。经分析，计算结果合理可信。1.3.3算例3(b3=0;-4;4;0;0;0;0T)算例3中，由b3可知，节点2为需求节点，节点3为供给节点，由节点3将信息传 输至节点2的最短路径为弧5-弧1。图4算例3最优传输示意图求得的最优解为x3*=4 0 0 0 4

6、 0 0 0 0 0 0 0 0T,即经过弧5运输4个单位流量至节 点1，再经弧1运输4个单位流量至节点2,其余弧无流量。又因为，每条弧的费用均为5， 所以最小费用为40。经分析，计算结果合理可信。1.3.4算例4(b4=4;0;0;0;0;0;-4T)算例4中，由b4可知，节点7为需求节点，节点1为供给节点，由节点1将信息传 输至节点7的最短路径为弧1-弧4-弧10。流量45弧4涂量4图5算例3最优传输示意图求得的最优解为x4*=4 0 0 4 0 0 0 0 0 4 0 0 0T,即经过弧1运输4个单位流量至节 点2，再经弧4运输4个单位流量至节点5,最后经弧5运输4个单位流量至节点7,其

7、 余弧无流量。又因为，每条弧的费用均为5，所以最小费用为60。经分析，计算结果合理可信。62重要算法编写与观察2.1习题5.6(a)初值为(0,0)时本算法令g的2范数在10-4时，停止迭代，经过86次迭代收敛收敛因子(f(k+1)-f*)/(f(k)-f*)=0.7623图6收敛因子截图(b)初值为(-0.4,0)时本算法令g的2范数在10-4时， 停止迭代， 经过112次迭代收敛 收敛因子(f(k+1)-f*)/(f(k)-f*)=0.817图7收敛因子截图(c)初值为(10,0)时本算法令g的2范数在10-4时，停止迭代，经过5次迭代收敛收敛因子(f(k+1)-f*)/(f(k)-f*)

8、=3.9022e-4图8收敛因子截图8(d)初值为(11,0)时本算法令g的2范数在10-4时，停止迭代，经过2次迭代收敛 收敛因子(f(k+1)-f*)/(f(k)-f*)= 0图9收敛因子截图图10自变量(x1,x2)截图9总结： 最速降线法的收敛因子随着初值的不同而变化， 对于个别初值 （如本习题初值 取 （11,0） 时） ，算法可迅速收敛。因此，初值的选取对于最速降线法的收敛速度有较大 影响。22习题5.7(a)由f(x) =9x -41 n(x-7)可得:故，牛顿迭代法的确切公式为:x -7(x-7)2（b）从以下五个初值开始迭代(1)x(0)=7.40表1初值1牛顿法迭代结果表迭

9、代次数x 值梯度值17.4-127.44-0.0909090937.4444-0.0009000947.4444444-9.00E-0857.4444444-9.00E-08(2) x(0)=7.20表2初值2牛顿法迭代结果表迭代次数x 值梯度值17.2-1127.31-3.637.403775-0.747.440723-7.60E-0257.444413-0.000631068(3)x(0)=7.01表3初值3牛顿法迭代结果表迭代次数x 值梯度值f (x) =94x -7f (x)二4(x -7)21017.01-3911127.019775-193.275600537.03867-94.4

10、7.073976-4.51E+0157.135638-20.(4)x(0)=7.80表4初值4牛顿法迭代结果表迭代次数x 值梯度值17.8427.16-1637.2624-6.947.369879-1.81E+0057.431934-0.3(5)x(0)=7.88表5初值5牛顿法迭代结果表迭代次数x 值梯度值17.884.527.0176-218.272727337.034503-106.931813547.066328-5.13E+0157.122757-23.（c）本问题的最优值为7.4444444。由上述五个初值点的前五步迭代可以看出：当初值点在区间（7.4444444，7.8888）内

11、时，第二次迭代点将落在（7，7.4444444） 之间，随后逐渐增加，直至逼近最优值。当初值点在区间（7，7.4444444）内时，则迭代点逐渐增加，逼近最优值。当取初值不在（7,7.8888）内时，牛顿法不收敛。2.3习题5.8（a）没有线搜索的牛顿法卩=0.1时， 卩=1时,（b）具有线搜索的牛顿法卩=0.1时,卩=1时，（未完成）122.4习题5.9（a）初值选（1.2,1.2）时，最速降线法：本算法令g的2范数在10-2时，停止迭代，经过3262次迭代得到以下结果Cent; ur ofX-awis图11最速降线法初值为（121.2）的等值线图及迭代轨迹牛顿法：本算法令s的4范数在10-

12、6时，停止迭代，经过4次迭代得到以下结果13CQHlUr of图12牛顿法初值为（1.2,1.2）的等值线图及迭代轨迹（b）初值选（-1.2,1）时，最速降线法：本算法令g的4范数在10-2时，停止迭代，经过6835次迭代得到以下结果Ccritouror Rosribr:ckX-axis14图13最速降线法初值为（-1.2,1）的等值线图及迭代轨迹15牛顿法：本算法令s的4范数在V10-6时，停止迭代，经过6次迭代得到以下结果: nntnurbrockX- axis图14牛顿法初值为（-1.2,1）的等值线图及迭代轨迹2.5习题5.19N=5迭代6次后，满足收敛条件。表6 N=5时，各迭代点x

13、值迭代次数/分量1444510000040.7744410.7744410.7744410.7744410.7744414-1.744581.0449944.4054814.8945444.45495444.740748-14.4459-4.780467.94544517.419445-4.8061445.656-86.4661-46.19499.441765.000468-140640000175-140640-1140640N=8迭代19次后，满足收敛条件。表7 N=8时，各迭代点x值迭代次数/分量144456781610000000040.7544940.75

14、44940.7544940.7544940.7544940.7544940.7544940.7544944-1.714480.4868491.186971.7445094.1075684.4845444.5988544.77040844.619757-7.74774-5.18495-1.1844.6614656.0744949.04446111.64715175-4.5694948.44644-40.4749-44.1741-44.1857-1.644645.1444954.749864.577401-74.4411199.8408-9.99548-174.858-171.851-10.865

15、4469.47447-4.15904104.8461-645.041081.886446.6416-1014.4-914.4591401.4468-5.65084154.4884-874.6011478.11445.7144-1444.46-1156.481454.9449-5.64898154.4878-874.6011478.108445.7146-1444.47-1156.481454.944106.688844-489.1174810.06-9974.6514761.111879.748-15549.98487.119116.744416-489.4544810.494-9974.11

16、14764.741880.191-15541.78488.07146.789876-489.4814811.108-9976.4614765.091880.84-15544.48489.454147.645489-445.6814087.908-6008.644107.69416980.64-46494.611414.66147.405147-181.141454.897-1497.76-10051.444195.97-48554.514876.44155.044404-90.441479.065814858.444-47184.458454.08-55845.619754.9716-7.98

17、171504.7645-7559.4146199.76-148600416414.8-16816751480.4617-7.99416504.8604-7559.6446199.85-148600416415.4-16816851480.1618-7.99414504.8646-7559.6446199.85-148600416415.4-16816851480.1519-8.00048504.9999-756046400-148600416416-1681685148040-8.00004504.0004-756046400-148600416416-16816851480N=14迭代49次

18、后，满足收敛条件。（表略）N=40迭代74次后，满足收敛条件。（表略）2.6习题5.27调用MATLAB自带的Isqnonlin.m函数，计算可得对应的x（1）、x（2）和标准差如下表 所示。表8选取各初值的计算结果初值x1 虚实部x1 虚部x2 实部x2 虚部r标准差1 10.00070-15.33030.2767i0.00750.1 0.10.00060-14.41710.02670.00750.01 0.010.00060-12.94320.00260.0075由上可知，标准差值较为恒定，随初值变化不十分显著；x1和x2值随初值选取的不同而不同。2.7习题6.4（未完成）end183附录

19、3.1对偶单纯形法函数MATLAB程序function sol,val,kk=duioudanchun(A,N)B=A;mA,nA=size(A);kk=0;flag=1;while flagkk=kk+1;if A(:,nA)=0flag=0;sol=zeros(1,nA);for i=1:mA-1 sol(N(i)=A(i,nA);endval=sol*(B(mA,:);elsefor i=1:mA-1if A(i,nA)=0disp( have infinite solution! ); flag=0;break ;endendif flagtemp=0;for i=1:mA-1if A

20、(i,nA)temp temp=A(i,nA); outb=i;endsita=zeros(1,nA-1);for i=1:nA-1if A(outb,i)0sita(i)=A(mA,i)/A(outb,i);endend19temp=-inf;for i=1:nA-1if sita(i)temptemp=sita(i);inb=i;endendfor i=1:mA-1if i=outb N(i)=inb;endendA(outb,:)=A(outb,:)/A(outb,inb);for i=1:mAif i=outbA(i,:)=A(i,:)-A(outb,:)*A(i,inb);A(mA,

21、nA)=0;endendendendend3.2最速降线法求Rosenbrock函数最小值matlab程序如下：function rb = rbfun(x,y)rb=100*(y-xA2)A2+(1-x)A2endclearclcsyms x y g Gg=gradient(rb(x,y),x y) %定义梯度向量G=hessian(rb(x,y),x y) %定义海森阵X(1,:)=-1.4 1;%定义初始点x=X(1,1);y=X(1,4);A(1,:)=subs(g)%给梯度赋初值20i=1while( norm(A(i,:),4)10A(-4)f(i)=rb(x,y)P(i,:)=-A

22、(i,:) fz(i)=-A(i,:)*P(i,:)%-gT*p fm(i)=P(i,:)*subs(G)*P(i,:)a(i)=fz(i)/fm(i)X(i+1,:) = X(i,:)+a(i)*P(i,:)x=X(i+1,1);y=X(i+1,4)A(i+1,:)=subs(g)i=i+1end%收敛条件%记录函数值%得到迭代方向%精确搜索法步长的分子%精确搜索法步长的分母%精确搜索法步长%产生新的点%产生新的梯度3.3牛顿法求Rosenbrock函数最小值matlab程序如下：function rb = rbfun(x,y)rb=100*(y-xA2)A2+(1-x)A2endclearclcsyms x y g Gg=gradient(rb(x,y),x y) %定义梯度向量21G=hessian(rb(x,y