1流量工程问题1.1问题重述定义一个有向网络G=(N,E)，其中N是节点集，E是弧集。

令A是网络G的点弧关联矩阵，即N×E阶矩阵，且第l列与弧里(I,j)对应，仅第i行元素为1，第j行元素为-1，其余元素为0。

再令bm =(bm1,…,bmN)T，fm=(fm1,…,fmE)T，则可将等式约束表示成：Af m=b m本算例为一经典TE算例。

算例网络有7个节点和13条弧，每条弧的容量是5个单位。

此外有四个需求量均为4个单位的源一目的对，具体的源节点、目的节点信息如图所示。

这里为了简单，省区了未用到的弧。

此外，弧上的数字表示弧的编号。

此时，c=((5,5 (5)1×13)T，根据上述四个约束条件，分别求得四个情况下的最优决策变量x=((x12,x13,…,x75)1×13)。

图 1 网络拓扑和流量需求1.27节点算例求解1.2.1\T）1.2.2算例1（b1=[4;-4;0;0;0;0;0]转化为线性规划问题：Minimize c T x1Subject to Ax1=b1x1>=0利用Matlab编写对偶单纯形法程序，可求得:最优解为x1*=[4 0 0 0 0 0 0 0 0 0 0 0 0]T对应的最优值c T x1=201.2.3算例2（b2=[4;0;-4;0;0;0;0]T）Minimize c T x2Subject to Ax2=b2\X2>=0利用Matlab编写对偶单纯形法程序，可求得:最优解为x2*=[0 4 0 0 0 0 0 0 0 0 0 0 0]T对应的最优值c T x2=201.2.4算例3（b3=[0;-4;4;0;0;0;0]T）Minimize c T x3Subject to Ax3=b3X3>=0利用Matlab编写对偶单纯形法程序，可求得:最优解为x3*=[4 0 0 0 4 0 0 0 0 0 0 0 0]T@对应的最优值c T x3=401.2.5算例4（b4=[4;0;0;0;0;0;-4]T）Minimize c T x4Subject to Ax4=b4X4>=0利用Matlab编写对偶单纯形法程序，可求得:最优解为x4*=[4 0 0 4 0 0 0 0 0 4 0 0 0]T对应的最优值c T x4=601.3计算结果及结果说明1.3.1算例1（b1=[4;-4;0;0;0;0;0]T）\算例1中，由b1可知，节点2为需求节点，节点1为供给节点，由节点1将信息传输至节点2的最短路径为弧1。

图 2 算例1最优传输示意图求得的最优解为x1*=[4 0 0 0 0 0 0 0 0 0 0 0 0]T，即只经过弧1运输4个单位流量，其余弧无流量。

又因为，每条弧的费用均为5，所以最小费用为20。

经分析，计算结果合理可信。

1.3.2算例2（b2=[4;0;-4;0;0;0;0]T）算例2中，由b2可知，节点3为需求节点，节点1为供给节点，由节点1将信息传输至节点2的最短路径为弧2。

图 3 算例2最优传输示意图求得的最优解为x2*=[0 4 0 0 0 0 0 0 0 0 0 0 0]T，即只经过弧2运输4个单位流量，其余弧无流量。

又因为，每条弧的费用均为5，所以最小费用为20。

经分析，计算结果合理可信。

1.3.3，T）1.3.4算例3（b3=[0;-4;4;0;0;0;0]算例3中，由b3可知，节点2为需求节点，节点3为供给节点，由节点3将信息传输至节点2的最短路径为弧5->弧1。

图 4 算例3最优传输示意图求得的最优解为x3*=[4 0 0 0 4 0 0 0 0 0 0 0 0]T，即经过弧5运输4个单位流量至节点1，再经弧1运输4个单位流量至节点2，其余弧无流量。

又因为，每条弧的费用均为5，所以最小费用为40。

经分析，计算结果合理可信。

1.3.5算例4（b4=[4;0;0;0;0;0;-4]T）算例4中，由b4可知，节点7为需求节点，节点1为供给节点，由节点1将信息传输至节点7的最短路径为弧1->弧4->弧10。

图 5 算例3最优传输示意图求得的最优解为x4*=[4 0 0 4 0 0 0 0 0 4 0 0 0]T，即经过弧1运输4个单位流量至节点2，再经弧4运输4个单位流量至节点5，最后经弧5运输4个单位流量至节点7，其余弧无流量。

又因为，每条弧的费用均为5，所以最小费用为60。

;经分析，计算结果合理可信。

2重要算法编写与观察2.1习题（a）初值为（0,0）时本算法令g的2范数在<10-4时，停止迭代，经过86次迭代收敛。

收敛因子（f(k+1)-f*）/（f(k)-f*）=图 6 收敛因子截图（b）初值为（,0）时本算法令g的2范数在<10-4时，停止迭代，经过112次迭代收敛。

·收敛因子（f(k+1)-f*）/（f(k)-f*）=图 7 收敛因子截图（c）初值为（10,0）时本算法令g的2范数在<10-4时，停止迭代，经过5次迭代收敛。

收敛因子（f(k+1)-f*）/（f(k)-f*）=图 8 收敛因子截图（d）初值为（11,0）时本算法令g的2范数在<10-4时，停止迭代，经过2次迭代收敛。

—收敛因子（f(k+1)-f*）/（f(k)-f*）= 0图 9 收敛因子截图图 10 自变量（x1,x2）截图总结：最速降线法的收敛因子随着初值的不同而变化，对于个别初值（如本习题初值取（11,0）时），算法可迅速收敛。

因此，初值的选取对于最速降线法的收敛速度有较大影响。

2.2 习题（a ） 由()94ln(7)f x x x =--可得：4'()97f x x =-- 24"()(7)f x x =-故，牛顿迭代法的确切公式为：—24974(7)x s x --=--（b ）从以下五个初值开始迭代 （1）x(0)=表、（2）x(0)=表（3）x(0)=￥表（4）x(0)=表（5）x(0)=表（c）本问题的最优值为。

由上述五个初值点的前五步迭代可以看出：当初值点在区间（，）内时，第二次迭代点将落在（7，）之间，随后逐渐增加，直至逼近最优值。

当初值点在区间（7，）内时，则迭代点逐渐增加，逼近最优值。

当取初值不在（7,）内时，牛顿法不收敛。

2.3习题（a）{（b）没有线搜索的牛顿法μ=时，μ=1时，（c）具有线搜索的牛顿法μ=时，μ=1时，，（未完成）2.4习题（a）初值选（,）时，最速降线法：本算法令g的2范数在<10-2时，停止迭代，经过3262次迭代得到以下结果。

图 11 最速降线法初值为（,）的等值线图及迭代轨迹牛顿法：本算法令s的4范数在<10-6时，停止迭代，经过4次迭代得到以下结果。

图 12 牛顿法初值为（,）的等值线图及迭代轨迹）（b）初值选（,1）时，最速降线法：本算法令g的4范数在<10-2时，停止迭代，经过6835次迭代得到以下结果。

图 13 最速降线法初值为（,1）的等值线图及迭代轨迹牛顿法：本算法令s的4范数在<10-6时，停止迭代，经过6次迭代得到以下结果。

图 14 牛顿法初值为（,1）的等值线图及迭代轨迹2.5[2.6习题N=5迭代6次后，满足收敛条件。

N=8迭代19次后，满足收敛条件。

N=14迭代49次后，满足收敛条件。

（表略）N=40迭代74次后，满足收敛条件。

（表略）2.7习题调用MATLAB自带的函数，计算可得对应的x(1)、x(2)和标准差如下表所示。

由上可知，标准差值较为恒定，随初值变化不十分显著；x1和x2值随初值选取的不同而不同。

2.8习题（未完成）3附录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-1sol(N(i))=A(i,nA);?endval=sol*(B(mA,:))';elsefor i=1:mA-1if A(i,nA)<0&A(i,1:nA-1)>=0 disp('have infinite solution!');flag=0;break;endendif flag^temp=0;for i=1:mA-1if A(i,nA)<temptemp=A(i,nA);outb=i;endendsita=zeros(1,nA-1);for i=1:nA-1if A(outb,i)<0sita(i)=A(mA,i)/A(outb,i);`endendtemp=-inf;for i=1:nA-1if sita(i)<0&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~=outb A(i,:)=A(i,:)-A(outb,:)*A(i,inb);A(mA,nA)=0;endendend(endend3.2最速降线法求Rosenbrock函数最小值matlab程序如下：function rb = rbfun(x,y)rb=100*(y-x^2)^2+(1-x)^2endclearclcsyms x y g G》g=gradient(rb(x,y),[x y]) %定义梯度向量G=hessian(rb(x,y),[x y]) %定义海森阵X(1,:)=[ 1]; %定义初始点x=X(1,1);y=X(1,4);A(1,:)=subs(g) %给梯度赋初值i=1while(norm(A(i,:),4)>10^(-4)) %收敛条件f(i)=rb(x,y) %记录函数值P(i,:)=-A(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+1end3.3牛顿法求Rosenbrock函数最小值matlab程序如下：function rb = rbfun(x,y)rb=100*(y-x^2)^2+(1-x)^2endclearclcsyms x y g Gg=gradient(rb(x,y),[x y]) %定义梯度向量G=hessian(rb(x,y),[x y]) %定义海森阵X(1,:)=[ 1]; %定义初值x=X(1,1);y=X(1,4);A(1,:)=subs(g) %给梯度赋初值H=subs(inv(G)) %得到海森阵初值S(1,:)=-A(1,:)*H %得到s初值i=1while (norm(S(i,:),4)>10^(-6)) %收敛条件f(i)=rb(x,y) %定义函数值X(i+1,:) = X(i,:)+S(i,:) %得到下一迭代点x=X(i+1,1);y=X(i+1,4) %给x,y分别赋值A(i+1,:)=subs(g) %得到新的梯度值H=subs(inv(G)) %得到新的海森阵S(i+1,:)=-A(i+1,:)*H %得到新的增量si=i+1end3.4共轭梯度法求解习题程序如下：clearclcK=40G=zeros(K,K)for m = 1: Kfor n = 1: KG(m,n)=1/(m+n-1)endendX(1,:)=zeros(1,K)b=ones(1,K)A(1,:)=X(1,:)*G-bP(1,:)=-A(1,:)i=1while (norm(A(i,:),4)>10^(-6)) %收敛条件d=P(i,:)*Gfz(i)=A(i,:)*A(i,:)' %精确搜索法步长的分子fm(i)=P(i,:)*d' %精确搜索法步长的分母a(i)=fz(i)/fm(i) %精确搜索法步长X(i+1,:) = X(i,:)+a(i)*P(i,:) %产生新的点A(i+1,:)=A(i,:)+a(i)*dbeta(i+1)=(A(i+1,:)*A(i+1,:)')/(A(i,:)*A(i,:)')P(i+1,:)=-A(i+1,:)+beta(i+1)*P(i,:)i=i+1end。