FR共轭梯度法与拟牛顿法计算机实现及仿真.DOC

1、FR 共轭梯度法与拟牛顿法计算机实现及仿真信计 0701 作者：翟铮 学号：0120714410105摘要：最速下降法和牛顿法是最基本的无约束最优化方法，但由于最速下降法收敛速度慢而牛顿法虽然收敛速度比较快，但需要计算目标函数的 Hesse 矩阵及其逆矩阵，计算量比较大。我们这里采用采用一种无需计算二阶导数，但收敛速度比较快的一种算法，即 FR 共轭梯度。关键词:FR 共轭梯度法、拟牛顿法、黄金分割法、最优一维线性搜索、收敛性1、 FR 共轭梯度法与拟牛顿法1、共轭梯度法在共轭梯度方向法中，如果取初始的搜索方向 ，而以以下各共轭方向 由第 k 次迭代点的负)(00xfdd梯度 与已经得到的共轭

2、方向 的线性组合来确定，这样就构造了一种具体的共轭方向法。因为每一个)(kxf1kd共轭方向都依赖于迭代点处的负梯度，所以称之为共轭梯度法（conjugate gradient method）给定初始点 ，取初始搜索方向 ，从 出发，沿 进行一维搜索，得到迭代点 ，一nR0 )(00xf00d1x下按 2,1,)(11 nkdfdkk 来构造搜索方向， 的选取应该使得 和 是共轭的。因为kk kTkTkkT QdxfQd)(11且要使 和 是共轭的，应有 ，所以1kd0k 2,10,)(1nkdxfTkk 于是得到 n 个搜索方向 如下：110,ndkTkkkkQdxffd)()(11100F

3、letcher_Reeves 公式为 21|)(|kkxf2、拟牛顿法：牛顿法成功的关键是利用了 Hesse 矩阵提供的曲率信息，但计算 Hesse 矩阵工作量大，并且有的目标函数的 Hesse 矩阵很难计算，甚至不好求出。针对这一问题，拟牛顿法比牛顿法更为有效。这类算法仅利用目标函数值和一阶导数的信息，构造出目标函数的曲率近似，使方法具有类似牛顿法的收敛速度快的优点。一般的拟牛顿法算法为：给出初始点 .0,)(,00 kRHBRxnn或如果 停 止 。,kg解 kkd得搜索方向 kd（或计算 .kkgH）由线性搜索得步长因子 .,1kkkdx并 令校正 ，产 生或 校 正产 生 )H(B1k

4、kkk 使得拟牛顿条件成立。本次实验采用 DFP 校正，校正公式为：kTkTkk yys12、 计算机实现步骤由于该问题设计的问题比较复杂，所用变量多，为防止变量之间的冲突，在下面的程序中全部用函数来实现。函数一：计算函数在任意点处的梯度值。函数名：gradient_my(f,y,num)，参数： f:待求梯度函数 y：待求梯度点 num:函数变量数函数二：黄金分割法求最优步长。函数名：golden_search(f,a,b) ， 参数：f: 待求梯度函数 a:搜索区间左端点 b:搜索区间右端点函数三：共轭梯度法。函数名：conjugate_gradient（f,x0,error） ，参数：f

5、:待求梯度函数 x0:初始点 error：允许误差函数四：拟牛顿法函数名：quasi_Newton(f,x0,error),参数：f:待求梯度函数 x0:初始点 error：允许误差三、实验结果及仿真第一题、Rosenbrock 函数 ,)1()(10)( 22xxf 0,.*0 fxTT设定黄金分割法的步长区间为 ，搜索精度为 0.001，设定共轭梯度法误差为：0.0011.,调用函数：conjugate_gradient(f,-1.2,1,0.001);-1.5 -1 -0.5 0 0.5 1 1.5-0.200.20.40.60.811.2XY与与与与与与与 70与与与与与与与与收 敛

6、点初 始 点迭 代 中 间 过 程结论：从图上可以看到利用共轭梯度法求解过程中，从初始点向收敛点迭代的步长是阶段性的变化的，期间有的阶段步长很小，阶段间的步长比较大。此方法在搜索在收敛中比较快，有效的避免了最速下降法的锯齿效应。第三题：Wood 函数 )1(8.19)()1(.0 .9042242 2433 xxxf ,3, *fTT解：我们用上述的共轭梯度法求解，导致搜索精度低，达不到函数停止迭代的下限，当搜索次数增大时，导致了迭代不收敛，因此我们试图改用拟牛顿法来求解，观察拟牛顿法在处理高维函数时的优异性。我们用拟牛顿法设定搜索步长区间为-15,15，误差精度为：0.001 得到的结果如下

7、：迭代次数 X Y Z W1 -3.00000 -1.00000 -3.00000 -1.00000 2 -2.19396 -0.86038 -2.27451 -0.87381 3 -2.44735 -0.57344 -2.05230 -0.54747 4 -3.08965 5.91018 -2.94903 5.31308 5 1.83003 6.50315 0.98570 5.77727 6 -2.58836 7.31350 -2.55296 5.74866 7 -2.07720 4.07977 -2.88577 8.54865 8 -2.08072 3.93519 -2.93265 8.6

8、3671 9 -1.42740 1.08110 -3.22100 10.46009 10 -0.80527 0.27410 -3.06216 8.88119 76 1.07314 1.14360 0.91394 0.81642 77 1.06164 1.12478 0.92729 0.84102 78 1.04445 1.09663 0.94733 0.87837 79 1.01963 1.04288 0.96666 0.91740 80 1.01630 1.03885 0.96776 0.92561 81 1.01066 1.03129 0.97072 0.94075 82 0.97768

9、0.96009 1.01125 1.02000 83 0.96822 0.93931 1.02478 1.04712 84 0.98839 0.97846 1.00975 1.02145 85 1.00025 1.00087 1.00022 1.00027 86 1.00020 1.00042 0.99978 0.99957 87 0.99999 0.99998 1.00000 1.00001 0 10 20 30 40 50 60 70 80 90-4-2024681012 与与与与与与与与与与与xyzw-4 -3-2 -10 12-202468-2024681012X与与与与Wood与与与

10、与87与与与与与与与与YW初 始 点 （ -3， -1， -1）收 敛 点 （ 1， 1， 1）-50510-4-3-2-1012-2024681012Y与与与与与与 87与与与与与与与与ZW初 始 点 （ -1， -3， -1）收 敛 点 （ 1， 1， 1）第四题：Powell 奇异函数： 414324321 )(0)()(5)0() xxxxxf ,0, *fTT解：设定黄金分割法的步长区间为：-0.1,0.1，搜索精度为：0.0001，共轭梯度法的终止点误差范围为。 迭代终止点为：0.0673 -0.0067 0.0332 0.033401.|)(|2*xf调用函数： conjugat

11、e_gradient(f,3,-1,0,1,0.0001); 可得结果为：由于函数为四维，我们以下下分别用三个图形来表示收敛效果，分别列出了图 1、X-Y-Z 图 2、X-Y-W 图三 Y-Z-W0 0.51 1.52 2.53-1-0.5000.10.20.30.40.50.60.7X与与与与与与与与与与 0.0001与与与与与YZ收 敛 点 初 始 点图 10 0.51 1.52 2.53-1-0.8-0.6-0.4-0.2000.511.522.5X与与与与与与与与与与与与与0.0001与与与与与与YW初 始 点收 敛 点图 2-1 -0.8-0.6 -0.4-0.2 000.20.40

12、.60.800.511.522.5Y与与与与与与与与与与与与与0.0001与与与与与ZW初 始 点 收 敛 点图 3结论：从上图可以看到，共轭梯度法搜索过程中，初始步长较大，进行搜索，随着距离收敛点的距离面小，步长逐步变小，向最优解逼近。附程序：共轭梯度法主程序：function A=conjugate_gradient(f,x0,error)a,b=size(x0);initial_gradient=gradient_my(f,x0,b);norm=0;norm0=0;syms step_zzh;A=;for i=1:bnorm0=norm0+(initial_gradient(i)2;en

13、dsearch_direction=-initial_gradient;x=x0+step_zzh*search_direction;f_step=subs(f,findsym(f),x);best_step=golden_search(f_step,-0.5,0.5);x_1=x0+best_step*search_direction;norm=norm0;k=1;while normerror norm=0;for i=1:bnorm=norm+(gradient(i)2;endnew_direction=-gradient+norm/norm0*search_direction;x=x_

14、1+step_zzh*new_direction;f_step=subs(f,findsym(f),x);best_step=golden_search(f_step,-0.1,0.1)x_2=x_1+best_step*new_directionA=A;x_2;norm0=norm;search_direction=new_direction;x_1=x_2;k=k+1;endw,z=size(A);if z=2plot(A);enda_return=x;拟牛顿法主程序：function A=quasi_Newton(f,x0,error)a,b=size(x0);G0=eye(b);initial_gradient=gradient_my(f,x0,b);norm0=0;norm0=initial_gradient*initial_gradient;syms step_zzh;A=x0;search_direction=-initial_gradient;x=x0+step_zzh*search_direction;f_step=subs(f,findsym(f),x);