现在科学计算方法.docxVIP

实验报告一实验时间：2012.10.23地点：理学院机房 报告人：王振宇题目共轭梯度法解线性方程组实验目标运用matlab实现共轭梯度法解线性方程组数学原理共轭梯度法又称共轭斜量法，是解线性代数方程组和非线性方程组的一种数值方法，例如对线性代数方程组 Ax=?, (1)式中A为n阶矩阵，x和?为n维列向量,当A对称正定时,可以证明求(1)的解X*和求二次泛函 　　 (2)的极小值问题是等价的。此处(x,у)表示向量x和у的内积。由此，给定了初始向量x（0），按某一方向去求(2)式取极小值的点x(1),就得到下一个迭代值x(2),再由x(2)出发，求x(3)等等，这样来逼近x*。若取求极小值的方向为F在 x(k=1,2,…)处的负梯度方向就是所谓最速下降法，然而理论和实际计算表明这个方法的收敛速度较慢，共轭梯度法则是在 x(k-1)处的梯度方向r(k-1)和这一步的修正方向p(k-1)所构成的二维平面内，寻找使F减小最快的方向作为下一步的修正方向p(k),即求极小值的方向p（其第一步仍取负梯度方向）。计算公式为 　　 公式再逐次计算 公式(k=1,2,…)。可以证明当i≠j时， 　　 公式从而平p(1),p(2)形成一共轭向量组;r(0),r(1),…形成一正交向量组。后者说明若没有舍入误差的话，至多 n次迭代就可得到(1)的精确解。然而在实际计算中,一般都有舍入误差，所以r(0),r(1),…并不真正互相正交，而r(0)，r(1),…等也只是逐步逼近(1)的真解,故一般将共轭梯度法作为迭代法来使用。 　　近来在解方程组(1)时,常将共轭梯度法同其他一些迭代法结合作用。特别是对病态方程组这种方法往往能收到比较显著的效果。其方法是选取一对称正定矩阵B并进行三角分解，得B=LLT。将方程组(1)化为 hу=b, (3) 　　此处y=lTx,b=l-1?,h=l-1Al-T,而 　　 公式再对(3)用共轭梯度法，计算公式为 　　 公式k=0,1,2,…)适当选取B，当B很接近A时，h的条件数较之A大大减小,从而可使共轭梯度法的收敛速度大为加快，由一些迭代法的矩阵分裂A=M -N，可选取M 为这里的B，例如对称超松弛迭代(SSOR)，强隐式迭代(SIP)等,这类方法常称为广义共轭梯度法或预条件共轭梯度法，它也可用于解代数特征值问题。四．程序设计(请将Matlab代码打印页附在该报告扉页)原方法function x = cgm(A,b)n=length(b);x=zeros(n,1);r=b-(A*x);p=r;while (norm(r)~=0) s=r; a=(r*r)/(p*A*p); x=x+a*p; r=r-a*A*p; c=(r*r)/(s*s); p=r+c*p;endend A=[0 1 -1 2；1 0 -1 1；-2 1 0 1];b=[1 0 1]; cgm(A,b)ans =0 0 11 A=[23.6 -6.8 4.5 6.3 7.2 8.6;4.2 7.62 0.8 -3.1 -4.5 5.3;9.8 -6.8 45.2 3.07 -4.37 7.82;8.4 -9.05 24.4 2.65 6.53 -3.64;-9.5 4.54 15.8 4.89 4.8 -7.84;26.32 3.08 4.6 -3.5 -7.8 25.6];b=[-9.5 8.12 42.3 -21.3 4.7 3.47]; cgm(A,b) 进入死循环修改后function x = rcgm(A,b,tol) x = b; r = b - A*x; p = r; for k = 1:numel(b); z = A*p; a= (r*r)/(p*z); x = x + a*p; s = r*r; r = r - a*z; if( norm(r) tol ) return; end B = (r*r)/s; p = r + B*p; endA=[23.6 -6.8 4.5 6.3 7.2 8.6;4.2 7.62 0.8 -3.1 -4.5 5.3;9.8 -6.8 45.2 3.07 -4.37 7.82;8.4 -9.05 24.4 2.65 6.53 -3.64;-9.5 4.54 15.8 4.89 4.8 -7.84;26.32 3.08 4.6 -3.5 -7.8 25.6];b=[-9.5 8.12 42.3 -21.3 4.7 3.47];tol=1e-6 rcgm(A,b,to