第三章 线性代数方程组的共轭梯度法

华北电力大学数理学院
School of mathematics & physics
§2 最速下降法 近似解的误差估计
x
(k )
−x
* A
⎛ λn − λ1 ⎞ ( 0 ) ⎟ x − x* =⎜ ⎜λ +λ ⎟ ⎝ n 1⎠
k
A
注：λ1、λ n分别是A的最小和最大特征值, 收敛速度 λ n − λ1 由q = 决定。q越小收敛越快。 λ n + λ1
df = ( p ( k ) )T Ap ( k ) ⋅ t + ( p ( k ) )T ( Ax ( k ) − b) = 0 dt
x ( k +1) ： 导出
华北电力大学数理学院
School of mathematics & physics
§1 下降法的一般理论
r
(k )
= b − Ax
(k ) T
A的特征值：0.0440
4.0000
7.9560
A的真解：x=[-0.51428571428573 6.85714285714288 -4.11428571428571] q=0.9890 最速下降法迭代568次，得到满足精度0.0005的解： X=[ -0.5088 6.8513 -4.1159]
例2：用最速下降法解方程组
3.5 0 ⎤ ⎡ x1 ⎤ ⎡ 24 ⎤ ⎡4 ⎢3.5 4 − 1.5⎥ ⎢ x 2 ⎥ = ⎢ 30 ⎥ ⎢ ⎥⎢ ⎥ ⎢ ⎥ ⎢ 0 − 1.5 4 ⎥ ⎢ x 3 ⎥ ⎢− 24⎥ ⎣ ⎦⎣ ⎦ ⎣ ⎦
A的特征值：0.1921 A的真解：x=[6 0 -6] q=0.9520
∵ ∇ f ( x ) = Ax − b , ∇ 2 f ( x ) = A 注： ∴ 当 A对称正定时， ∇ f ( x ) = Ax − b = 0 ⇔ min f ( x )。 n
x∈ R
如何寻找f(x)的极小点呢？
华北电力大学数理学院
School of mathematics & physics
k =1
n
(n)
= b − Ax
(n)
− ∑ tk Ap ( k )
k =1
∴ (r ( n ) , p ( j ) ) = (r (0) , p ( j ) ) − (t j Ap ( j ) , p ( j ) ) = (r (0) , p ( j ) ) − (r ( j -1) − r ( j ) , p ( j ) ) = (r (0) , p ( j ) ) − (r ( j -1) , p ( j ) ) = (r (0) , p ( j ) ) − (r (0) − ∑ tk Ap ( k ) , p ( j ) )
华北电力大学数理学院
School of mathematics & physics
§3 共轭梯度法
对于 k = 1, 2 ,
现有x ( k )及共轭方向p ( k −1)，则 r ( k ) = b − Ax ( k ) 在r ( k )和p ( k −1) 确定的超平面上 找共轭方向 p ( k ) = r ( k ) + s k p ( k −1) 其中 (p ( k −1) ) T Ar ( k ) s k = − ( k −1) T ( k −1) (p ) Ap
( k +1)
=x
(k )
+tkr
(k )
华北电力大学数理学院
School of mathematics & physics
§2 最速下降法 function x=zuisu(a,b,x,e) n=length(b); r=b-a*x; q=r’*r; while q>e w=a*r; t=q/(r'*w); x=x+t*r; r=r-t*w; q=r’*r; end
(p ( k ) ) T Ar ( k +1) sk = − (k) T (k) (p ) Ap p ( k +1) = r ( k +1) + s k p ( k )
华北电力大学数理学院
School of mathematics & physics
§3 共轭梯度法
可以证明：CG产生的序列有 1、 p , p ,
华北电力大学数理学院
School of mathematics & physics
§1 下降法的一般理论 具体过程： 已求得x
(k )
1：在 x ( k )处确定一个使 f ( x )下降的方向 p ( k )；
2：在射线 x = x ( k ) + tp ( k ) 上求 f ( x )的极小点 x ( k +1) , 即 f ( x ( k +1) ) = min f ( x ( k ) + tp ( k ) )
§1 下降法的一般理论 基本思想就是下降法:
从某一点x (0)出发, 逐步产生一串点： x (0) , x (1) ,
使 f ( x (0) > f ( x (1) ) > f ( x(k ) ) > 并以“最快的速度”下降到f ( x)的极小值。
, x(k ) ,
完成任务的关键就是确定每步的下降方向！
华北电力大学数理学院
School of mathematics & physics
§3 共轭梯度法
对于 k = 1 , 2 ,
现有x 及共轭方向p , 则
(k) (k )
(r ( k ) ) T p ( k ) t k = (k) T (k) (p ) Ap x ( k +1) = x ( k ) + t k p ( k ) r ( k +1) = b − Ax ( k +1)
华北电力大学数理学院
School of mathematics & physics
第三章 线性代数方程组的共轭梯度法(CG) 问题：解Ax=b，其中，A为对称正定矩阵。
原理： Ax = b ⇔ min f ( x)。A是对称正定矩阵。 n
x∈R
其中， = 1 xT Ax − bT x。 f(x) 2
4.0000 7.8079
最速下降法迭代141次，得到满足精度0.0005的解： X=[5.9987 0.0015 -5.9994]
华北电力大学数理学院
School of mathematics & physics
§2 最速下降法
例3：用最速下降法解方程组
3.8 0 ⎤ ⎡ x1 ⎤ ⎡ 24 ⎤ ⎡4 ⎢3.8 4 − 1.1⎥ ⎢ x 2 ⎥ = ⎢ 30 ⎥ ⎢ ⎥⎢ ⎥ ⎢ ⎥ ⎢ 0 − 1.1 4 ⎥ ⎢ x 3 ⎥ ⎢− 24⎥ ⎣ ⎦⎣ ⎦ ⎣ ⎦
华北电力大学数理学院
School of mathematics & physics
§2 最速下降法 小 结
最速下降法有简单易行，保稀疏性等特点，但 当A的最大特征值远远大于最小特征值时收敛速度变 得非常慢。最速下降法并非最速！ 共轭梯度法可使这一问题得到一定改善！
华北电力大学数理学院
School of mathematics & physics
§3 共轭梯度法
定义：向量p、q称之为关于对称正定矩阵A共轭，如 p Τ Aq = 0 。 果满足
华北电力大学数理学院
School of mathematics & physics
是关于n阶对称正定 p ( k )为下降方向的算 矩阵A共轭的向量组，则以 法：
x (0) ∈ R n ; r (0) = b − Ax (0) ; k = 1, 2, (r ( k −1) )T p ( k ) ； t k = (k) T (k) (p ) A p x ( k ) = x ( k −1) + t k p ( k )； r ( k ) = b − Ax ( k ) = r ( k −1) − t k Ap ( k )
定理1：设p 1 , p 2 ,
pn
有Ax
(n)
= b,即r
(n)
= 0.
华北电力大学数理学院
School of mathematics & physics
§3 共轭梯度法
证明： Ax ∵ r
(n)
= Ax
( n −1)
+ tn Ap =r
(0)
(n)
= Ax
n
(0)
+ ∑ tk Ap ( k )
(k ) (k )
(r ) p t k = (k) T (k) (p ) Ap x
( k +1)
=x
(k )
+t k p
(k )
t k 称之为x (k)到x (k +1)的步长。
华北电力大学数理学院
School of mathematics & physics
§2 最速下降法
p ( k ) 的确定！ 下降方向
t >0
3： 判 断
∇f (x
( k + 1)
) ≤ ε, 或
x
(k + 1 )
−x
(k)
≤ε ！
是 ！ x * ≈ x ( k +1)， 停 止 ; 否 ！ k = k + 1, 转 1.
华北电力大学数理学院
School of mathematics & physics
§1 下降法的一般理论 问题 f ( x ( k +1) ) = min f ( x ( k ) + tp ( k ) ) 的解决！
k =1 j -1
=0 ∴ r ( n )与p (1)，p (2)， ，p ( n )正交，必有r ( n ) = 0
华北电力大学数理学院
School of mathematics & physics
§3 共轭梯度法
定理1说明：取A共轭的方向的下降算法至多n步 得到n阶线性代数方程组的精确解。把这种方法 称之为共轭梯度法（CG）。
如何寻找A共轭的方向呢？
华北电力大学数理学院
School of mathematics & physics
§3 共轭梯度法
∀x ( 0 ) ∈ R n , r ( 0) = b − Ax ( 0 ) , p ( 0) = r ( 0) ; 开始：
(p ( 0) ) T p ( 0) t 0 = (k) T , x (1) = x ( 0 ) + t 0 p ( 0 ) . (p ) A p (k)
A的特征值：.8377,4.0000,7.1623 0
A的真解：x=[3 q=0.7906.
4
-5]
最速下降法迭代29次，得到满足精度0.0005的解： X=[2.9998 4.0003 -4.9999]
华北电力大学数理学院
School of mathematics & physics
§2 最速下降法
(0) (1)
,p ,
(k)
是两两关于A共轭的； , k − 1;
2、 p ) r (
(i ) T ( k )
= 0, i = 0,1, , r (k) ,
3、 r ( 0 ) , r (1) ,
是两两正交的；
函数的负梯度方向是函数值下降最快的方向，因 此，取 p ( k ) = − f ′( x ( k ) ) = b − Ax ( k ) = r ( k )为下降方向。 最速下降法：
r
(k )
= b − Ax
(k ) T
(k ) (k )
(r ) r tk = ( k ) T ( k ) (r ) Ar x
华北电力大学数理学院
School of mathematics & physics
wk.baidu.com
§2 最速下降法
例1：用最速下降法解方程组
⎡4 3 0 ⎤ ⎡ x1 ⎤ ⎡ 24 ⎤ ⎢3 4 − 1⎥ ⎢ x ⎥ = ⎢ 30 ⎥ ⎢ ⎥⎢ 2 ⎥ ⎢ ⎥ ⎢0 − 1 4 ⎥ ⎢ x 3 ⎥ ⎢− 24⎥ ⎣ ⎦⎣ ⎦ ⎣ ⎦
t >0
f ( x ( k ) + tp ( k ) ) = 1 ( x ( k ) + tp ( k ) )T A( x ( k ) + tp ( k ) ) − bT ( x ( k ) + tp ( k ) ) 2 = 1 ( p ( k ) )T Ap ( k ) ⋅ t 2 + ( p ( k ) )T ( Ax ( k ) − b) ⋅ t 2 + 1 ( x ( k ) )T Ax ( k ) − bT x ( k ) 2
华北电力大学数理学院
School of mathematics & physics
§3 共轭梯度法 综合上述，有下面的共轭梯度法：
∀x ( 0 ) ∈ R n , r ( 0 ) = b − Ax ( 0 ) , p ( 0 ) = r ( 0 ) ;
对于 k = 0 ,1 , 2 ,
(r ( k ) ) T p ( k ) t k = (k ) T (k ) (p ) Ap x ( k +1) = x ( k ) + t k p ( k ) r ( k +1) = b − Ax ( k +1) = r ( k ) − t k Ap( k )