精华内容
下载资源
问答
  • 共轭梯度法求解线性方程

    万次阅读 2017-04-09 20:09:14
    当线性方程组Ax= b的系数矩阵A是对称正定矩阵时,可以采用共轭梯度法对该方程组进行求解,可以证明,式(1)所示的n元二次函数  (1)  取得极小值点x*是方程Ax= b的解。共轭梯度法是把求解线性方程组的问题转化为...

    1.1算法原理及程序框图

    当线性方程组Ax= b的系数矩阵A是对称正定矩阵时,可以采用共轭梯度法对该方程组进行求解,可以证明,式(1)所示的n元二次函数

                                           (1)      

    取得极小值点x*是方程Ax= b的解。共轭梯度法是把求解线性方程组的问题转化为求解一个与之等价的二次函数极小值的问题。从任意给定的初始点出发,沿一组关于矩阵A的共轭方向进行线性搜索,在无舍入误差的假定下,最多迭代n次(其中n为矩阵A的阶数),就可求得二次函数的极小点,也就求得线性方程组Ax= b的解。其迭代格式为:

                                            (2)     

    经过一系列的证明和简化,最终可得共轭梯度法的计算过程如下:

    (1) 给定初始计算向量x(0)即精度ε>0;

    (2) 计算r(0)= b-Ax(0),取d(0) =r(0)

    (3) for k =0 to n-1 do

    (i);

    (ii);

    (iii);

    (iv)

    (v);

    (vi);

    End do

    计算程序如下:

    1.2程序使用说明

    共轭梯度法求解线性方程组的matlab已提交,该程序可以求解系数矩阵为对称正定矩阵的线性方程组。直接双击打开程序软件即可,在程序中根据问题需要更改系数矩阵A和矩阵b(直接输入矩阵或者用简单命令生成矩阵)。根据计算需要修改计算精度eps,点击运行即可得到结果。

    1.3算例计算

    选择课本第113页计算习题3.2计算。例题要求在N=100,N=200,N=400下的计算结果,以N=200为例先将N=100和N=400注释掉,通过diag()命令生成矩阵A和b,运行计算结果。当n=100,n=200,n=400时分别需要迭代的次数是50次,100次,200次,x的各元素为1。

    % ** 文件名:Conjugate Gradient.m
    % 
    % ** 日 期:2016.12.14
    % 
    % ** 描 述:共轭梯度法解线性方程组(Conjugate Gradient method)
    % 
    % ** 函数:113页计算实习3.2  
    % 
    % ** 参考教材:《数值分析》李乃成,梅立泉,科学出版社
    %
    clear;
    clc;
    %%
    %获得需要的N阶矩阵A
    N=100; %解向量的维数
    %N=200;
    %N=400;
    a0=eye(N);%a0为n阶的单位阵
    a1=eye(N-1);%a1为n-1阶的单位阵
    a11=diag(a1);%a11是a1的单位阵的元素
    a12=diag(a11,1);%a12是一个把a11的值放到其对角阵上一层的 矩阵
    A=a12'+a12+(-2)*a0;%A是对角线及其对角线上下(第n-1,n+1条对角线)有元素的矩阵
    %%
    %获得矩阵b
    b=zeros(N,1);
    b(1)=-1;
    b(N)=-1;
    b;
    %%  
    fprintf('库函数计算结果:');
    x=inv(A)*b      %库函数计算结果
    x=zeros(N,1);%迭代近似向量
    eps=0.0000001;%精度
    r=b-A*x;
    d=r;
    for k=0:N-1
        fprintf('第%d次迭代:',k+1);
        a=(norm(r)^2)/(d'*A*d)
        x=x+a*d;
        rr=b-A*x;    %rr=r(k+1)
        if (norm(rr)<=eps)||(k==N-1)
            break;
        end
        B=(norm(rr)^2)/(norm(r)^2);
        d=rr+B*d;
        r=rr;
    end
    x
        
    %由上述结果可知方程组由共轭梯度法求解所得的解与计算机直接计算所得的解相同,
    %故所得的的解释可靠的。
    %当n=100,n=200,n=400时分别需要迭代的次数是50次,100次,200次



    展开全文
  • 简单回顾一下,线性共轭梯度法是一种不需要矩阵求逆或矩阵分解,以迭代的方式求解线性方程的方法,而且可以保证在N迭代内收敛,其中N是变量的维度。对于维度较大的问题,辅以良好的预处理步骤,算法可以在很少的步...

    上篇文章介绍了线性共轭梯度法。简单回顾一下,线性共轭梯度法是一种不需要矩阵求逆或矩阵分解,以迭代的方式求解线性方程的方法,而且可以保证在N次迭代内收敛,其中N是变量的维度。对于维度较大的问题,辅以良好的预处理步骤,算法可以在很少的步数内收敛,得到一个比较好的近似解。

    在线性共轭梯度法提出约10年后,Fletcher和Reeves将其推广到非线性优化问题中,称为非线性共轭梯度法。这种新方法可以替代之前讲过的线搜索和信赖域方法,并且仍然在高维度问题中展现出非凡的性能。

    Fletcher-Reeves CG

    我们把Fletcher和Reeves最早提出的非线性共轭梯度法称为Fletcher-Reeves CG(FR-CG)算法。按照惯例,我们应该先介绍该方法提出的动机,并一步步推导其计算公式。但据我所知,FR-CG似乎不是从某个目标推导出来的,而是直接在CG(线性共轭梯度法)上改进而来的。直接替换CG算法中某些值的计算方法,如下所示

    与共轭梯度法(一):线性共轭梯度中的式(5.23)相比,只修改了两处。一是计算

    的方式,用线搜索(一):步长的选取中介绍的线搜索算法来搜索
    方向上可接受的步长。之所以这样做,是因为在非线性问题中,我们已经无从构造一组共轭向量了,而且也不可能保证N次迭代内收敛。FR-CG更像是一种线搜索方法,通过所谓共轭梯度的方式计算优化方向,然后用步长搜索算法计算一个可接受的步长。二是用梯度
    替代残差
    ,因为在CG中,残差恰好就是梯度,这里只不过将其推广而已。换句话说,式(5.40)不过是CG推广到普通非线性函数后的结果,CG则是其在二次函数下的特例。

    不过,强行推广到普通非线性函数就能奏效吗?那可未必。首先,我们要确认式(5.40)的每次迭代都能生成使得目标函数下降的方向。其次,我们还要探究该方法是否能保证全局收敛。

    第一个问题等价于证明

    ,我们代入式(5.40)的第4行,得到

    如果第

    次迭代求得了确切解,即
    方向上的最优步长,那么下一次迭代的梯度方向一定与上次迭代方向垂直,即
    ,于是
    成立。但是,确切解是不可能的,我们无法对普通的非线性函数求某个方向的极值,时间上也不允许,只能用前面提到的线搜索算法来搜索一个可接受的步长。所以式(5.41)中的
    很可能大于
    ,从而无法保证
    。不过,虽然线搜索无法找到确切最优解,我们仍可以把搜索的结果限制在可接受的范围内。在线搜索那篇文章中,我们讲解了Goldstein conditions和Wolfe conditions,它们用来限制搜索的区域。如果我们把后者的约束变得再严格一些,可以得到下面的strong Wolfe conditions

    它要求函数值有充分的下降,并要求导数的绝对值不要过大。可以证明,当常数因子满足

    时,
    一定成立。证明方式是数学归纳法,先验证
    成立。然后假设
    成立,再借助式(5.42)的第二行,可以推导出
    也成立。(证明过程略。)好在线搜索算法一般都在strong Wolfe conditions下进行,所以这个要求很容易满足。

    第二个问题,关于收敛性,则有更多值得讨论的地方。

    收敛性分析

    先从直观上感受FR-CG方法的收敛性能。想象非线性目标函数

    ,如果在非常接近极小值点的区域,
    恰好是一个严格的凸二次函数。那么我们就可以期待,在该区域内,算法能够像线性共轭梯度法一样在N次迭代内收敛。

    但是,注意到CG有效的一个前提是,初始方向必须是最速下降方向,这在FR-CG的迭代过程中是无法保证的。这时,我们需要一种称为重新启动(restart)的策略,让FR-CG每N次迭代后都从最速下降方向重新开始,具体操作是直接令

    。这样一来,一旦进入凸二次函数的区域,迟早会有一次重新启动,接下来收敛速度就非常快了。

    不过,由于FR-CG是用线搜索的方式计算步长

    ,因此即使进入凸二次函数区域,也不一定搜索出最优步长,这就违反了CG方法的要求。回顾线搜索那篇文章的最后一节,我们提到了用插值的方式加速线搜索。现在,我们会发现这一策略用在FR-CG上简直完美,一旦进入凸二次函数区域,二次插值的结果就会完全吻合原目标函数,直接就能求出最优步长。

    事实上,添加了重新启动策略的FR-CG方法,具有n-step意义上的二次收敛性质,即

    可以理解为,如果把每N次迭代合并为一次,得到的序列就是二次收敛的。

    此时,另一个问题出现了。每N次迭代进行一次重新启动并不一定是最优方案,当N特别大时,我们其实期望整个算法在远小于N次迭代时就收敛,所以重新启动可能永远不会执行。这种情况下,可以采用更灵活的重新启动策略,比如,当

    时,意味着连续两次迭代的梯度方向不垂直。但我们知道,线性共轭梯度法中,任意两次迭代的梯度都应该是垂直的,此时方可获得最快的收敛速度。因此,如果它们不垂直,就说明收敛速度慢下来了,需要重新启动。

    现在回归主题,既然是收敛性分析,就需要研究FR-CG是否具有全局收敛性。

    在线搜索的收敛性分析中提到过,当目标函数

    满足Lipschitz continuous性质,且线搜索步长满足Wolfe conditions时,下式成立

    我们知道,全局收敛的条件是

    ,要想满足该条件,需要
    。一个典型的例子是最速下降法中
    恒成立,因此最速下降法全局收敛。而在FR-CG中,每次重新启动也是满足
    的,所以所有重新启动对应的
    组成的序列满足
    。这意味着梯度的某个子序列是全局收敛的,即

    这里极限运算符中的

    表示取整个序列中极限值最小的那个子序列。这是一种打了折的全局收敛性,虽然不是严格意义上的全局收敛,但也不错。

    不过,如果没有重新启动呢,纯粹的FR-CG是否具有全局收敛性?答案是有的,可以通过反证法证明。假设

    ,利用strong Wolfe conditions,式(5.42)的第二行,可以推导出与假设相反的矛盾,详细过程本文就不讲了,需要的同学可以参考原书。

    缺陷与改进

    如果实际中使用FR-CG,会发现它收敛得很慢,而且有时会在很小的区域内停滞不前。这是因为其固有的一个缺陷,当某次迭代计算出的方向

    较差时,比如与梯度方向相差90°,线搜索的结果是一个极小的步长,这是显然的,因为
    恰好是函数
    的等高线方向。我们接下来看看此时会发生什么。

    在本文第一段的最后,我们简要提到了可以用数学归纳法证明

    。事实上,这部分证明的结果不只如此,完整的结论是

    其中

    。借助这个结果,我们继续证明刚才的问题。对上式左右同时乘上
    ,并带入
    的定义,得到

    其中,

    。回到最初的假设,一旦当前的
    相差90°,即
    ,就会导致式(5.52)的左右两端的值都是0,这意味着
    。另外,由于
    很小,迭代前后变量值基本不变,即
    ,于是迭代前后的梯度值也基本不变,即
    。由式(5.40)第三行可知,
    ,再由式(5.40)第四行可得,

    这一结果说明,从第

    次开始的每次迭代,都会前进一个相同且非常小的步长,相当于在做无用功。这也是导致FR-CG方法收敛慢的罪魁祸首。不过,既然找到了原因,自然可以想到改进方案。当出现这种情况的时候,让算法做一次重新启动,强行让迭代方向与梯度方向垂直,就可以避开这种无畏的消耗了。这一思想有个更直接的名称,PR-CG。

    Polak-Ribière CG

    PR-CG是FR-CG的一种改进,由Polak和Ribière提出,与FR-CG的唯一区别在于如何选取

    。选取方法为

    在凸二次函数情况下,

    是一样的,因为
    。但对于普通非线性函数,当出现上一节介绍的停滞不前的状况时,可以发现
    ,这恰好就相当于一次重新启动。所以PR-CG用改进
    计算公式的方式解决了FR-CG收敛慢的问题。

    但PR-CG也不是没有缺点,

    的定义会导致其不满足全局收敛性。更具体地,会使得某些迭代产生小于
    ,进而使得新的
    不再是下降方向。好在解决方法也很简单,用

    替代

    。这样做的意义是,当出现无法使得目标函数下降的方向时,执行重新启动策略。这种方法称为PR+-CG。

    总结

    本文介绍了非线性共轭梯度法中的FR-CG、PR-CG和PR+-CG三种方法。它们都源自于线性共轭梯度法,只不过修改了残差的计算方式和步长的确定方式,以适应任意的非线性目标函数。为了确保全局收敛性,选择在合适的时机使用重新启动策略非常重要,也是PR-CG和PR+-CG设计的初衷。

    至此,数值最优化系列文章已经推出了八篇,从宏观的层面上介绍了线搜索、信赖域和共轭梯度的整体思想和算法流程。从下篇文章开始,我们会进一步深入细节,探究各种方法的具体实现细节以及针对实际情况的改进。

    上一篇:共轭梯度法(一):线性共轭梯度
    下一篇:牛顿法进阶
    展开全文
  • 共轭梯度法是Hestenes和Stiefel(1952年)在求解大规模线性方程组时创立的一种迭代算法。1964年,Fletcher和Reeves将其引入非线性最优化问题,创立了非线性共轭梯度法。该方法基于前一迭代点的搜索方向对当前迭代点的...

    d9eb0839b34ab4ca9ea74ff67e1c935e.png

    共轭梯度法是Hestenes和Stiefel(1952年)在求解大规模线性方程组时创立的一种迭代算法。1964年,Fletcher和Reeves将其引入非线性最优化问题,创立了非线性共轭梯度法。

    该方法基于前一迭代点的搜索方向当前迭代点的负梯度方向进行修正来产生新的搜索方向,即搜索方向当前迭代点的负梯度方向与前一搜索方向的线性组合

    本期大纲

    1. 共轭方向

    2. 共轭方向法

    3. 正定二次函数的共轭梯度法

    4. 一般函数的共轭梯度法

    5. 算法实现

    6. 总结

    1共轭方向在介绍共轭梯度法之间,小编首先介绍一个重要的定义——共轭方向,并给出共轭方向的一些性质。

    定义(共轭方向):设Q是一个n×n的对称正定矩阵,若n维向量空间中的非零向量575e4fd959471d64d55025a769703233.png满足

    63017df2f41959e086ad54f4feaf9866.png  

    则称575e4fd959471d64d55025a769703233.png是Q的共轭向量或称向量575e4fd959471d64d55025a769703233.png关于Q共轭的(简称共轭),也称575e4fd959471d64d55025a769703233.png是Q的m个共轭方向。

    例如:

    20cf27d40a64c1dc9498d2737c7c8f7b.png

    注:当 Q=I(单位矩阵)时,公式(1)变为

    73c460ddf377f03922bb2efe5a4d0f65.png

    即向量 575e4fd959471d64d55025a769703233.png 互相正交。由此可见,正交是共轭的一种特殊情况共轭是正交的推广。 

    性质(线性无关性)设Q是一个n×n的对称正定矩阵,若575e4fd959471d64d55025a769703233.png是Q的m个非零共轭向量,则这m个向量是线性无关的。

    59c55c74680156d9c6f9a073b233b442.png

    注:在n维的向量空间中,非零的共轭向量个数不超过n

    2共轭方向法共轭方向法的算法步骤如下:

    58718a8e257a89470353be62194349f8.png

     在考虑普通函数之前,首先考虑如下的正定二次函数

    33a18ab04fbf99444c53a0f86159a2a6.png

    定理1:设有二次函数

    33a18ab04fbf99444c53a0f86159a2a6.png

    其中Q是一个n×n的对称正定矩阵,若575e4fd959471d64d55025a769703233.png是Q的m个非零共轭向量,从任意一点出发833324ee794f02d5903a691e627de3c4.png出发,依次沿这组向量进行一维搜索

    af3ee618e77898d745f5e8b68f63d4a1.png

    则有

    d5c49ec2a8a665d5ffdf14230a18423c.png

    当前点的梯度方向之前所有搜索方向正交

    定理2(算法收敛性——二次终止性):设有二次函数

    33a18ab04fbf99444c53a0f86159a2a6.png

    其中Q是一个n×n的对称正定矩阵,若575e4fd959471d64d55025a769703233.png是Q的m个非零共轭向量,从任意一点出发833324ee794f02d5903a691e627de3c4.png出发,依次沿这组向量进行一维搜索

    af3ee618e77898d745f5e8b68f63d4a1.png

    至多经过n步算法收敛,即2dd018a1af73fc6cce57a9e54709c8a4.png是f(x)在f955e2c9b46fcdcebe645fe8845a601f.png上的极小值点

    因此,对于任意n元正定二次函数,我们可以从任意点出发,然后沿着n个共轭方向最多做n次线搜索,就可以求的目标函数的最优点。3正定二次函数的共轭梯度法在共轭方向法中,我们需要预先给定共轭方向,这在实际应用中是很难实现的,下面我们给出一种新的算法——共轭梯度法,它可以随着迭代的进行不断产生共轭方向。共轭梯度法基本思想:在FR共轭梯度法中,初始点处的搜索方向取为该点的负梯度方向,即取  

    f98c6d1199e5cdf7fc1322207de8a740.png

    而之后点e8eb17ef239a9d5b7ec19c52f35d1106.png的搜索方向8875e068fb7fdd61b2944b6e3ca36ae7.png由该点的负梯度aa1c8233b23847799311192efec241a8.png与已得的共轭向量6db43fafb7daf92d3fb28106fbe122b0.png的线性组合来确定。

    命题1:设有二次函数

    33a18ab04fbf99444c53a0f86159a2a6.png

    其中Q是一个n×n的对称正定矩阵。若选取初始点为973e1ae3bbc734ef345c98043305d0a0.png,取faed79efddf779d1c25c8991669e5864.png。依次从d0acfee0bf7c1ce5d1191f80df7a5570.png出发,沿0f976331f00e90d098b12cc941c8e293.png在最优步长规则下做一维线搜索,可得:

    aa15cdb17c7ce7c8f5276e68fb19efb9.png

          共轭梯度法的算法步骤如下:

    Step 0. 选取初始点9529a2ef215ad8c0a145e886d9d7fcc5.png,允许误差ε> 0,令 k=1;

    Step 1. 计算1b79b7b79609f089694b71e5fb3397f5.png209b22685febe13d56b2bce10acab9ee.png,算法停止;.否则,转 Step 2;

    Step 2. 依次计算

    446e25d70d5f45c0286a25076f088088.png

    并置k=k+1,转 step 1.

    定理3:对于正定二次函数,最优步长规则下的共轭梯度法经过m≤n次一维搜索后终止,且对任意的1im,下列关系式成立:

    2d74ade002bc35aece58db83728a6565.png

    注:定理3中,(1)说明共轭梯度法生成的搜索方向是关于Q共轭的;(2)说明共轭梯度法生成的点列的梯度是相互正交的。

    命题2:

    1211b2146bec7aa728092fb4072cbb8c.png

    注:共轭梯度法中β的取值不同又可细分FR、PRP等方法,上述命题β的取值所对应的方法被称为FR共轭梯度法4一般函数的共轭梯度法

     上面所讲的迭代公式要想适用于非二次函数,存在以下问题:

    (1) 步长c1a1c7b87db477faca49228c4d08e391.png不能再用公式af4261ef87b5a656bdb98d69eeefc58f.png计算,必须使用其他的一维搜索方法来确定;(2) 凡用到矩阵Q之处,需要用现行点处的Hessian矩阵1e2391fdd0714647e76634d167d7e96f.png代替;(3) 有限步内迭代达不到极小点。


        迭代的延续:

    (1)  直接搜索:即总用公式71ef7613b901b21a0267aa2f76023586.png构造搜索方向。

    (2)  重新开始,把n步作为一轮,每搜索一轮之后,取一次最速下降方向,开始下一轮。

        一般函数的FR共轭梯度法的算法步骤如下:

    5f8520b43aa52692376d5d78a380b762.png

    3590da8cc19004fd580dc91938562781.png

    215512a6fe14a17782d461f1d364e664.png

    7f785b2ab5b4f557f03759a58b0936ba.png

    2ae5a607e5f41cdfc3a10c2eb77bbfca.png

    5算法实现

    例1:求解无约束优化问题

    53187c5e7b12a68cb5892803bfc9677d.png

    该问题有精确解c2efbca32cb09215829f86b3d6a705fd.png

    解:调用格式如下:

    function [X,Y,Error] = frcg(x0)

    %使用FR共轭梯度法求解无约束问题:min f(x)

    %使用非精确Armijo搜索步长规则

    %输入:x0为初始点,fun为目标函数,grad为梯度函数

    %输出:X,Y分别为每一步的迭代点及相应函数值,Error为误差

    tic

    %% 初始化参数

    f = @(x)100*(x(1)^2-x(2))^2+(x(1)-1)^2;                             %目标函数

    g = @(x)[400*x(1)*(x(1)^2-x(2))+2*(x(1)-1);-200*(x(1)^2-x(2))];     %一阶导数

    G = @(x)[1200*x(1)^2-400*x(2)+2, -400*x(1);-400*x(1), 200 ];        %二阶导数   

    %% 初始化算法参数

    k = 1;

    n = length(x0);

    e = norm(g(x0));

    X = x0;

    Y = f(x0);

    Error = e;

    rho=0.6;sigma=0.4;

    %% 设定终止条件

    N = 5000;                                                           %最大迭代步                

    E = 1e-6;                                                           %给定误差

    %% 迭代

    while k < N && e > E

        % 计算搜索方向

        if(k == 1)

            d = -g(x0);

        else       

            beta = (g(x0)'*g(x0))/(g0'*g0);

            d = -g(x0) + beta*d0;

            if(g(x0)'*d>=0)

                d = -g(x0);

            end

        end

        %使用Armijo法求步长

        m=0; mk=0;

        while(m<20)

            if(f(x0+rho^m*d)

                mk=m; break;

            end

            m=m+1;

        end

        g0 = g(x0);

        d0 = d;

        x0 = x0+rho^mk*d;                                       %迭代公式

        k = k+1;

        X(:,k) = x0;

        Y(k) = f(x0);

        e = norm(g(x0));

        Error(k) = e;

    end

    toc

    end

    调用格式如下:

    x0 = [0;0];                 %迭代初始点(可以根据需要修改)

    [X,Y,E]=frcg(x0);

    本问题求解结果如下:

    567d42b9eee8218145affc40efb85aa1.png

    ee9a6aaf4323a5919a78d2ae63bb6970.png

    e8aa2db5568e6ef6e900b3ff44275032.png

    6总结共轭方向法是介于最速下降法和Newton法之间的一种方法——它的收敛速度(二阶收敛)比最速下降法(线性收敛)快;它克服了最速下降法的锯齿现象,同时它的计算量和存储量又比牛顿法要小。共轭梯度法不需要预先给定Q的共轭方向,而是随着迭代的进行不断产生Q的共轭方向——在每次的迭代中,利用前一个搜索方向和当前迭代点的梯度的线性组合构造一个新的方向。7997d431feb806a0d999522a9dad3536.png几种下降算法的迭代过程

    往期回顾

    算法丨优化算法之牛顿法

    算法丨优化算法之最速下降法

    算法丨经典优化算法大合集(比赛备用,值得收藏)

    算法丨优化算法系列之遗传算法原理

    算法丨优化算法系列之遗传算法MATLAB程序设计及应用实例

    算法丨优化算法系列之模拟退火算法(1)

    算法丨优化算法系列之模拟退火算法(2)——0-1背包问题

    参考文献:马昌凤.最优化方法及其 Matlab 程序设计[M]. 科学出版社, 2010.

    c58b6b2fed7359f0976f3f19eef12f90.png

    展开全文
  • MATLAB - 共轭梯度法

    千次阅读 2019-10-07 10:25:29
    求解系数矩阵A是对称正定矩阵的线性方程组或求解二次函数的极小值点。 二、算法 共轭梯度法,步骤如下: 1 任意给定初始点及精度 2 3 对于,作 1) 2) 3),或 4) 若或,则输出,,取作为的解,...

    一、问题描述
           求解系数矩阵A是对称正定矩阵的线性方程组 Ax=b或求解二次函数f(x)=\frac{1}{2}x^{T}Ax-b^{T}x的极小值点。

    二、算法

          共轭梯度法,步骤如下:

         1 任意给定初始点x^{^{0}}及精度\varepsilon >0

         2  d^{(0)}=r^{(0)}=b-Ax

         3  对于k=0,1,...,n-1,作

              1)a_{(k)}=\frac{\left \| r^{(k)} \right \|\begin{matrix} 2 \\ 2 \end{matrix}}{<d^{(k)},Ad^{(k)}>}

              2) x^{(k+1)}=x^{(k)}+a_{(k)}d^{(k)}

              3) r^{^{(k+1)}}=b-Ax^{(k+1)} ,或 r^{(k+1)}=b-Ax^{(k+1)}=b-A(x^{(k)}+a_{(k)}d^{(k)})=r^{(k)}-a_{(k)}Ad^{(k)}

              4) 若\frac{\left \| r^{(k+1)} \right \|}{\left \| b \right \|}\leq \varepsilonk+1=n,则输出x^{(k+1)}r^{(k+1)},取r^{(k+1)}作为Ax=b的解,否则

              5)   

                      \beta _{(k+1),k}=\frac{\left \| r^{(k+1)} \right \|\begin{matrix} 2\\ 2 \end{matrix}}{\left \| r^{(k)}\right \|\begin{matrix} 2\\ 2 \end{matrix}}

              6)     d^{(k+1)}=r^{(k+1)}+\beta _{(k+1),k}d^{(k)}.

       MATLAB代码:

    function CGM_main()
    A = [10 1 2 3 4
        1 9 -1 2 -3
        2 -1 7 3 -5
        3 2 3 12 -1
        4 -3 -5 -1 15];
    b=[12 -27 14 -17 12]';
    x0=[0 0 0 0 0]';
    max_iter=10000;
    
    fprintf('\n');
    fprintf('共轭梯度法:\n');
    fprintf('==========================\n');
    [y,iter]=cgm (A,b,x0,max_iter);
    fprintf('\n');
    fprintf('迭代次数:\n   %d \n',iter);
    fprintf('方程的解: \n');
    fprintf('%10.6f',y);
    fprintf('\n\n=========================\n\n');
    end
    
    function [x,iter] = cgm (A,b,x0,max_iter)
    x=x0;
    epsilon=1.0e-6;
    fprintf('\n x0= ');
    fprintf('   %10.6f',x0);
    r=b-A*x;
    d=r;
    for k=0:max_iter
        alpha=(r'*r)/(d'*A*d);
        xx=x+alpha*d;
        rr=b-A*xx;
        if (norm(rr,2)/norm(b,2))<= epsilon
            fprintf('\n 找到啦~');
            iter = k+1;
            x=xx;
            r=rr;
            fprintf('\n x%d = ',k+1);
            fprintf('   %10.6f',x);
            fprintf('\n r%d = ',k+1);
            fprintf('   %10.6f',r);
            
            return
        end
        beta=(rr'*rr)/(r'*r);
        d=rr+beta*d;
        x=xx;
        r=rr;
        fprintf('\n x%d = ',k+1);
        fprintf('   %10.6f',x);
    end
    iter = max_iter;
    return 
    end
    

        运行结果:

    共轭梯度法:
    ==========================
    
     x0=      0.000000     0.000000     0.000000     0.000000     0.000000
     x1 =      1.073560    -2.415510     1.252487    -1.520877     1.073560
     x2 =      1.305605    -2.627981     2.146636    -1.694270     0.442393
     x3 =      1.446618    -2.225384     2.448048    -1.970691     0.620722
     x4 =      1.086550    -2.063574     2.792911    -2.101645     0.836386
     找到啦~
     x5 =      1.000000    -2.000000     3.000000    -2.000000     1.000000
     r5 =      0.000000     0.000000    -0.000000    -0.000000     0.000000
    迭代次数:
       5 
    方程的解: 
      1.000000 -2.000000  3.000000 -2.000000  1.000000
    
    ===================== 

    注:编写MATLAB程序时主函数一定请放在最前面...


    参考文献:
    [1] MATLAB R2014官方手册
    [2] 近代优化方法,徐成贤,陈志平,李乃成编著.科学出版社,2002.
    [3] 最优化方法(第二版).孙文瑜,徐成贤,朱德通编著.高等教育出版社,2010.
     

    展开全文
  • 最优化方法-共轭梯度法

    千次阅读 2019-03-31 17:16:23
    共轭梯度法最初由Hesteness和Stiefel于1952年为求解线性方程组而提出的。其基本思想是把共轭性与最速下降方法相结合,利用已知点处的梯度构造一组共轭方向,并沿这组方向进行搜素,求出目标函数的极小点。根据共轭...
  • 二次型 (Quadratic form)是一个形如f(x)=12xTAx−bTx+cf(x)=12xTAx−bTx+cf(x) = \frac{1}{2}x^TAx - b^Tx + c的标量二次方程,如果AAA是一个n∗nn∗nn*n对称正定阵,那么f(x)f(x)f(x)的最小值在Ax=bAx=bAx=b时取到...
  • (一)共轭梯度算法

    万次阅读 2019-06-19 16:37:26
    共轭梯度算法是干什么? 共轭梯度算法是一种迭代算法,在一次次的对待中最终...共轭梯度法是属于最小化类的迭代方法。为了求解Ax = b这样的一次函数,可以先构造一个二次齐函数 f(x)=12xTAx−bTxf(x) = \frac{1}{...
  • 共轭梯度法是Hestenes和Stiefel(1952年)在求解大规模线性方程组时创立的一种迭代算法。1964年,Fletcher和Reeves将其引入非线性最优化问题,创立了非线性共轭梯度法。该方法基于前一迭代点的搜索方向对当前迭代点的...
  •   此次实验要求采用共轭梯度下降来解决二次优化问题,在解决此问题前,首先分析下什么是共轭梯度法共轭梯度法能解决什么问题。   我们来看一个线性方程组Ax=b,求解方程组的过程可以看成是形如公式(1)的优化...
  • 上一节笔记:数值优化(4)——非线性共轭梯度法,信赖域法————————————————————————————————————大家好!俗话说得好,DDL是唯一生产力……在DDL的逼迫下,高产自然就来了2333。...
  • 多重网格法可能是目前为止解泊松方程最快的算法,n个格点需要n计算就可以收敛,而快速傅里叶变换的收敛速度是n*logn, 共轭梯度法是n^2.。多重网格法可以方便的应对各种边界条件,这一点比傅里叶变换之类的谱方法...
  • 最优化程序

    2015-11-13 12:43:24
    资源包含了最速下降法(精确步长搜索和wolfe条件非精确步长搜索),共轭梯度法,预处理共轭梯度法求解二次线性方程组的最优问题
  • 二次规划问题

    千次阅读 2013-05-21 10:56:22
    二次规划(Quadratic programming),在运筹学当中,...否则的话,常用的二次规划解法有:内点法(interior point)、有效集法(active set)和共轭梯度法等。凸集二次规划问题是凸优化问题的一个特例。 每个二次
  • 介绍的算法有:无约束问题的最速下降法、Newton法、拟Newton法、共轭梯度法、信赖域算法和直接法;非线性方程组和最小二乘问题的Newton法和拟Newton法;约束问题的罚函数法、乘子法、可行方向法、序列二次规划算法和...
  • .Taloy法求解常微分方程 1.重新学习Taloy公式: 基本思想:将复杂函数用多项式代替 p(x)=f(x)+f’(x)(x-x0) 低多项式表示复杂函数缺点: (1)精度不高,x0与x趋近时,误差为(x-x0)的高阶无穷小...
  • 赛德尔迭代法161 3.11 求解对称正定方程组的共轭梯度法165 3.12 求解线性最小二乘问题的豪斯荷尔德变换法169 3.13 求解线性最小二乘问题的广义逆法175 3.14 求解病态方程组189 第4章 非线性方程方程组的求解195 ...
  • 4.13 求解对称正定方程组的共轭梯度法 4.14 求解线性最小二乘问题的豪斯荷尔德变换法 4.]5 求解线性最小二乘问题的广义逆法 4.16 病态方程组的求解 第5章 非线性方程方程组的求解 5.1 非线性方程方程组类设计 ...
  • 科学与工程数值计算算法

    热门讨论 2010-12-21 12:45:03
    3.13 求解对称正定方程组的共轭梯度法 3.14 求解线性最小二乘问题的豪斯荷尔德变换法 3.15 求解线性最小二乘问题的广义逆法 3.16 病态方程组的求解 第4章 非线性方程方程组的求解 4.1 非线性方程方程...
  • C++常用算法程序集

    热门讨论 2011-03-24 09:15:32
    3.11 求解对称正定方程组的共轭梯度法165 3.12 求解线性最小二乘问题的豪斯荷尔德变换法169 3.13 求解线性最小二乘问题的广义逆法175 3.14 求解病态方程组189 第4章 非线性方程方程组的求解195 4.1 求非线性方程实...
  • C#数值计算算法编程 代码

    热门讨论 2009-12-11 11:31:58
    4.13 求解对称正定方程组的共轭梯度法 4.14 求解线性最小二乘问题的豪斯荷尔德变换法 4.]5 求解线性最小二乘问题的广义逆法 4.16 病态方程组的求解 第5章 非线性方程方程组的求解 5.1 非线性方程方程组类设计 ...
  • preconjgrad 预处理共轭梯度法求线性方程组Ax=b的解 BJ 块雅克比迭代法求线性方程组Ax=b的解 BGS 块高斯-赛德尔迭代法求线性方程组Ax=b的解 BSOR 块逐次超松弛迭代法求线性方程组Ax=b的解 第13章: 随机数生成 PFQZ...
  • MATLAB常用算法

    热门讨论 2010-04-05 10:34:28
    preconjgrad 预处理共轭梯度法求线性方程组Ax=b的解 BJ 块雅克比迭代法求线性方程组Ax=b的解 BGS 块高斯-赛德尔迭代法求线性方程组Ax=b的解 BSOR 块逐次超松弛迭代法求线性方程组Ax=b的解 第13章: 随机数生成 PFQZ...
  • 3.11 求解对称正定方程组的共轭梯度法165 3.12 求解线性最小二乘问题的豪斯荷尔德变换法169 3.13 求解线性最小二乘问题的广义逆法175 3.14 求解病态方程组189 第4章 非线性方程方程组的求解195 4.1 求非线性方程实...
  • 范例1-32 头插建立单链表 75 ∷相关函数:createlist函数 1.3.2 限制链表长度建立单链表 77 范例1-33 限制链表长度建立长单链表 77 ∷相关函数:createlist函数 1.3.3 尾插建立单链表 79 范例1-34 尾插...
  • C语言通用范例开发金典.part2.rar

    热门讨论 2012-08-31 14:18:18
    范例1-32 头插建立单链表 75 ∷相关函数:createlist函数 1.3.2 限制链表长度建立单链表 77 范例1-33 限制链表长度建立长单链表 77 ∷相关函数:createlist函数 1.3.3 尾插建立单链表 79 范例1-34 尾插...
  • C 开发金典

    2013-06-20 16:20:03
    范例1-32 头插建立单链表 75 ∷相关函数:createlist函数 1.3.2 限制链表长度建立单链表 77 范例1-33 限制链表长度建立长单链表 77 ∷相关函数:createlist函数 1.3.3 尾插建立单链表 79 范例1-34 尾插...
  • C语言通用范例开发金典.part1.rar

    热门讨论 2012-08-31 14:09:26
    范例1-32 头插建立单链表 75 ∷相关函数:createlist函数 1.3.2 限制链表长度建立单链表 77 范例1-33 限制链表长度建立长单链表 77 ∷相关函数:createlist函数 1.3.3 尾插建立单链表 79 范例1-34 尾插...

空空如也

空空如也

1 2
收藏数 27
精华内容 10
关键字:

共轭梯度法求解二次方程