用matlab实现共轭梯度法求解实例.doc

用MATLAB 实现共轭梯度法求解实例 康福 201103710031 1．无约束优化方法 1.1 无约束优化方法的必要性 一般机械优化设计问题，都是在一定的限制条件下追求某一指标为最小，它 们都属于约束优化问题。但是为什么要研究无约束优化问题?(1)有些实际问题，其数学模型本身就是一个无约束优化问题。(2)通过熟悉它的解法可以为研究约束优化问题打下良好的基础。(3)约束优化问题的求解可以通过一系列无约束优化方法来达到。所以无约束 优化问题的解法是优化设计方法的基本组成部分，也是优化方法的基础。 (4)对于多维无约束问题来说，古典极值理论中令一阶导数为零，但要求二阶 可微，且要判断海赛矩阵为正定才能求得极小点，这种方法有理论意义， 但无实用价值。和一维问题一样，若多元函数F(X)不可微，亦无法求解。 但古典极值理论是无约束优化方法发展的基础。 1.2共轭梯度法 目前已研究出很多种无约束优化方法，它们的主要不同点在于构造搜索方向 上的差别。 (1)间接法——要使用导数，如梯度法、 (阻尼)牛顿法、变尺度法、共轭梯度 法等。 (2)直接法——不使用导数信息，如坐标轮换法、鲍威尔法单纯形法等。 用直接法寻找极小点时，不必求函数的导数，只要计算目标函数值。这类方 法较适用于解决变量个数较少的(n ≤20)问题，一般情况下比间接法效率低。 间接法除要计算目标函数值外，还要计算目标函数的梯度，有的还要计算其海赛 矩阵。 搜索方向的构成问题乃是无约束优化方法的关键。 共轭梯度法是沿着共轭方向进行搜索，属于共轭方向法中的一种，该方法中 每一个共轭向量都是依赖于迭代点处的负梯度而构造出来。共轭梯度法作为一种 实用的迭代法，它主要有下面的优点： (1)算法中，系数矩阵Ａ的作用仅仅是用来由已知向量P产生向量W=AP，这不仅 可充分利用Ａ的稀疏性，而且对某些提供矩阵Ａ较为困难而由已知向量P产 生向量W=AP又十分方便的应用问题是很有益的。 (2)不需要预先估计任何参数就可以计算，这一点不像SOR等； (3)每次迭代所需的计算，主要是向量之间的运算，便于并行化。 共轭梯度法原理的知识较多，请详见《机械优化设计》第四章的第四、五节。 图1为共轭梯度法的程度框图 1 ( 0,1,2, ) k k k k s k       x x图1为共轭梯度法的程度框图 2．设计题目及要求 2.1设计题目 用共轭梯度法求二次函数 2 2 1 2 1 2 1 1 2 ( , ) 2 4 2 f x x x x x x x    的极小点及极小值。 2.2设计要求 (1)使用matlab编写程序，熟练撑握matlab编程方法。(2)学习并撑握共轭梯度法的原理、方法及应用，并了解不同无约束优化方法的 区别、优缺点及特殊要求。 (3)编写程序，计算出二次函数的极小点及极小值，并适当选取不同的初始点及 迭代精度精度，分析比较结果。 三．计算步骤 3.1计算求解 解：已知初始点[1,1] T迭代精度 0.001   1)第一次沿负梯度方向搜寻 计算初始点处的梯度： 为一维搜索最佳步长，应满足 得：2)第二次迭代代入目标函数 由 得 从而有： 因 收敛。 0 1 2 0 2 1 2 2 4 4 ( ) 4 2 2 x x f x x                    x x 0 1 0 0 0 0 0 1 4 1 4 1 2 1 2                               x x d 1 0 0 2 ( ) min ( ) min(40 20 3) f f           x x d 0 0.25   1 2 0.5        x 1 1 ( ) 2 f           x 2 1 2 0 0 ( ) 5 0.25 20 ( ) f f       x x 1 1 0 0 2 ( ) 1.5 f             d x d 2 1 1 2 2 2 2 0.5 1.5 0.5 1.5                              x x d 2 2 ( ) (2 2 ) 2(0.5 1.5 ) 2(2 2 )(0.5 1.5 ) 4(2 2 ) ( ) f x                  ( ) 0     1   2 2 2 4 0 , ( ) 8, ( ) 2 0 f f                  x x x 2 ( ) 0 f     x3.2运行与程序 运行：打开matlab,确定conjugate_grad_2d.m文件夹为当前目录。在命令窗中输入：f=conjugate_grad_2d([1,1],0.001)选择不同的初始点坐标[0,0],[0,1],[1,0],和迭代精度0.01，0.0001,进 行运行时，需要多次调用conjugate_grad_2d函数。 程序及说明： function f=conjugate_grad_2d(x0,t) %用共轭梯度法求已知函数f(x1,x2)=x1^2+2*x2^2-4*x1-2*x1*x2的极值点 %已知初始点坐标：x0 %已知收敛精度：t %求得已知函数的极值：f x=x0; syms xi yi a; %定义自变量，步长为符号变量 f=xi^2+2*yi^2-4*xi-2*xi*yi; %创建符号表达式f fx=diff(f,xi); %求表达式f对xi的一阶求导 fy=diff(f,yi); %求表达式f对yi的一阶求导 fx=subs(fx,{xi,yi},x0); %代入初始点坐标计算对xi的一阶求导实值 fy=subs(fy,{xi,yi},x0); %代入初始点坐标计算对yi的一阶求导实值 fi=[fx,fy]; %初始点梯度向量 count=0; %搜索次数初始为0 while double(sqrt(fx^2+fy^2))>t %搜索精度不满足已知条件s=-fi; %第一次搜索的方向为负梯度方向if count<=0s=-fi;elses=s1;endx=x+a*s; %进行一次搜索后的点坐标f=subs(f,{xi,yi},x); %构造一元搜索的一元函数φ(a)f1=diff(f); %对函数φ(a)进行求导f1=solve(f1); %得到最佳步长aif f1~=0ai=double(f1); %强制转换数据类型为双精度数值elsebreak %若a=0，则直接跳出循环，此点即为极值点endx=subs(x,a,ai); %得到一次搜索后的点坐标值f=xi^2+2*yi^2-4*xi-2*xi*yi;fxi=diff(f,xi);fyi=diff(f,yi);fxi=subs(fxi,{xi,