【数值计算方法（黄明游）】矩阵特征值与特征向量的计算（二）：乘幂法的加速（带有原点移位的乘幂法）【理论到程序】

Qomolangma

发布于 2024-07-30 10:40:34

1000

发布于 2024-07-30 10:40:34

文章被收录于专栏：深度学习

矩阵的特征值（eigenvalue）和特征向量（eigenvector）在很多应用中都具有重要的数学和物理意义。乘幂法（Power Iteration）是线性代数中一种重要的数值计算方法，用于估计矩阵的最大特征值及其对应特征向量的迭代算法，广泛应用于许多科学和工程领域。

本文将详细介绍带有原点移位的乘幂法，并给出其Python实现。

一、乘幂法

给定一个方阵

，如果

\lambda

是

的最大特征值，而

是对应的特征向量

那么对于任意非零向量

，迭代

A^n v

的结果会趋近于

\lambda^n x

。当

越来越大时，

A^n v

的方向会趋近于

的方向，而其模（长度）则会趋近于

\lambda^n

的绝对值。

【数值计算方法（黄明游）】矩阵特征值与特征向量的计算（一）：乘幂法【理论到程序】

二、乘幂法的加速

1. 天书

2. 人话

乘幂法（Power Iteration）是一种用于寻找矩阵的最大特征值及其对应特征向量的迭代算法。它通过迭代计算矩阵与向量的乘积，并规范化得到新的向量，最终收敛到矩阵的最大特征值和对应的特征向量。然而，对于某些矩阵，乘幂法的收敛速度可能相对较慢。为了加速乘幂法的收敛，一种常见的做法是通过平移（Shift）矩阵的方式。

a. 基本乘幂法

选择初始向量

x_0

：通常选择一个非零向量作为初始向量，其选择可能影响到迭代的收敛速度。

迭代计算

A x_k

：对于每一次迭代

，计算

A x_k

。

规范化

A x_k

：通过规范化

A x_k

（即将其除以其最大分量），得到新的向量

x_{k+1}

。

判断收敛条件：检查

\frac{\|x_{k+1} - x_k\|}{\|x_{k+1}\|}

是否小于一个预定的收敛阈值，满足则停止迭代。

更新

x_k

：将

x_{k+1}

赋值给

x_k

，重复步骤2。

计算特征值：一旦迭代收敛，通过

\frac{A x_k}{x_k}

的比值来估计矩阵

的最大特征值。

b. 加速乘幂法的思想

加速乘幂法的思想是通过平移矩阵

，将其转化为

B = A - \lambda I

的形式。这个平移操作使得新的矩阵

的最大特征值与原矩阵

的最大特征值之间的差距变大，从而加速收敛。

c. 具体步骤

平移操作： 对给定的矩阵

，选择一个数值

\lambda

，构造新矩阵

B = A - \lambda I

。

基本乘幂法： 使用

进行上述乘幂法迭代计算，直至收敛

计算特征值： 一旦迭代收敛，通过

\frac{(Bx_k)^T x_k}{(x_k)^T x_k} + \lambda

来估计矩阵

的最大特征值。

3. 典例

4. 实现

a.代码

import numpy as np

def normalize_vector(vec):
    # 计算向量的最大分量
    # max_val = np.max(np.abs(vec))
    max_index = np.argmax(np.abs(vec))
    max = vec[max_index]
    # 对向量进行规范化
    return vec / vec[max_index], max


def power_iteration(A, x, max_iter=1000, tol=1e-6):
    max_val0 = 0
    for i in range(max_iter):
        # 计算矩阵与向量的乘积
        np.set_printoptions(precision=4, suppress=True)
        Ax = np.dot(A, x)
        # 对乘积结果进行规范化
        print("          x{}：".format(i + 1), Ax)
        x, max_val = normalize_vector(Ax)
        print("第{}次：".format(i + 1))
        print("特征值max(x{})：".format(i + 1), max_val)
        print("特征向量   y{}：".format(i + 1), x)
        # 如果最大分量小于收敛阈值，提前退出迭代
        if abs(max_val0 - max_val) < tol:
            break
        max_val0 = max_val
        # 计算对应的特征值
        eigenvalue = np.dot(Ax, x) / np.dot(x, x)

    return eigenvalue


# 给定矩阵 A
A = np.array([[-4, 14, 0],
              [-5, 13, 0],
              [-1, 0, 2]])

A0 = np.array([[1, 1, 0.5],
               [1, 1, 0.25],
               [0.5, 0.25, 2]])
# lambda0 = 0.75
# B = A0-lambda0*np.eye(3)

A3 = np.array([[-11, 11, 1],
               [11, 9, -2],
               [1, -2, 13]])

lambda0 = 15
B = A3 - lambda0 * np.eye(3)
# 初始向量
x0 = np.array([1, 1, 1])

# 运行乘幂法迭代
# result_vector, eigenvalue = power_iteration(A, x0)
# result_vector, eigenvalue = power_iteration(A0, x0)
# result_vector, eigenvalue = power_iteration(A3, x0)

eigenvalue0 = power_iteration(B, x0)
eigenvalue = lambda0 + eigenvalue0
print("估计的特征值:", eigenvalue)
# np.dot(A, result_vector)
# np.dot(eigenvalue, result_vector)

b.解析

normalize_vector 函数：
- 输入：一个向量 vec。
- 功能：计算向量的最大分量，并将向量规范化。
- 输出：规范化后的向量和最大分量。
power_iteration 函数：
- 输入：
  - A：一个方阵（矩阵）。
  - x：初始向量。
  - max_iter：最大迭代次数，默认为1000。
  - tol：收敛阈值，默认为1e-4。
- 功能：使用乘幂法迭代来估计矩阵的最大特征值及其对应的特征向量。
  - 计算矩阵 A 与向量 x 的乘积，得到 Ax。
  - 调用 normalize_vector 函数对 Ax 进行规范化，得到规范化后的向量和最大分量。
  - 打印每次迭代的结果，即特征值、特征向量。
  - 判断是否满足收敛条件，如果最大分量变化小于阈值 tol，则提前退出迭代。
  - 计算对应的特征值，更新最大分量，并继续迭代。
- 输出：估计得到的特征值。
主程序部分：
- 教材例题及课后题的矩阵 A、A1、A2、A3。
- 定义了初始向量 x0、系数
\lambda
- 构造新矩阵
B = A - \lambda I
。
- 调用 power_iteration 函数，传入矩阵
B
和初始向量进行乘幂法迭代。
B
的特征值+
\lambda
=
A
的最大特征值。

c.应用

          x1： [1.75 1.5  2.  ]
第1次：
特征值max(x1)： 2.0
特征向量   y1： [0.875 0.75  1.   ]
          x2： [1.4688 1.3125 1.875 ]
第2次：
特征值max(x2)： 1.875
特征向量   y2： [0.7833 0.7    1.    ]
          x3： [1.3958 1.2083 1.8167]
第3次：
特征值max(x3)： 1.8166666666666667
特征向量   y3： [0.7683 0.6651 1.    ]
          x4： [1.3572 1.1846 1.8005]
第4次：
特征值max(x4)： 1.8004587155963303
特征向量   y4： [0.7538 0.658  1.    ]
          x5： [1.3464 1.1683 1.7914]
第5次：
特征值max(x5)： 1.7914012738853502
特征向量   y5： [0.7516 0.6522 1.    ]
          x6： [1.3401 1.1646 1.7888]
第6次：
特征值max(x6)： 1.7888444444444445
特征向量   y6： [0.7491 0.6511 1.    ]
          x7： [1.3383 1.1619 1.7873]
第7次：
特征值max(x7)： 1.7873301200029814
特征向量   y7： [0.7488 0.6501 1.    ]
          x8： [1.3373 1.1613 1.7869]
第8次：
特征值max(x8)： 1.7869153405872398
特征向量   y8： [0.7484 0.6499 1.    ]
          x9： [1.337  1.1608 1.7867]
第9次：
特征值max(x9)： 1.7866588534107755
特征向量   y9： [0.7483 0.6497 1.    ]
          x10： [1.3368 1.1608 1.7866]
第10次：
特征值max(x10)： 1.786591427357151
特征向量   y10： [0.7482 0.6497 1.    ]
          x11： [1.3368 1.1607 1.7865]
第11次：
特征值max(x11)： 1.7865478412280433
特征向量   y11： [0.7482 0.6497 1.    ]
          x12： [1.3367 1.1607 1.7865]
第12次：
特征值max(x12)： 1.7865369093758368
特征向量   y12： [0.7482 0.6497 1.    ]
          x13： [1.3367 1.1606 1.7865]
第13次：
特征值max(x13)： 1.786529489372963
特征向量   y13： [0.7482 0.6497 1.    ]
          x14： [1.3367 1.1606 1.7865]
第14次：
特征值max(x14)： 1.7865277239374533
特征向量   y14： [0.7482 0.6497 1.    ]
          x15： [1.3367 1.1606 1.7865]
第15次：
特征值max(x15)： 1.7865264587497982
特征向量   y15： [0.7482 0.6497 1.    ]
          x16： [1.3367 1.1606 1.7865]
第16次：
特征值max(x16)： 1.786526175001885
特征向量   y16： [0.7482 0.6497 1.    ]
估计的特征值: 2.536526458749798

5. 课后题

a. 乘幂法

        x1： [ 1 18 12]
第1次：
特征值max(x1)： 18
特征向量   y1： [0.0556 1.     0.6667]
          x2： [11.0556  8.2778  6.7222]
第2次：
特征值max(x2)： 11.055555555555555
特征向量   y2： [1.     0.7487 0.608 ]
          x3： [-2.1558 16.5226  7.407 ]
第3次：
特征值max(x3)： 16.52261306532663
特征向量   y3： [-0.1305  1.      0.4483]
          x4： [12.8835  6.6682  3.6974]
第4次：
特征值max(x4)： 12.88351581508516
特征向量   y4： [1.     0.5176 0.287 ]
          x5： [-5.0197 15.0842  3.6957]
第5次：
特征值max(x5)： 15.084204811029013
特征向量   y5： [-0.3328  1.      0.245 ]
          x6： [14.9056  4.8494  0.8523]
第6次：
特征值max(x6)： 14.90555759004166
特征向量   y6： [1.     0.3253 0.0572]
          x7： [-7.364  13.8137  1.0926]
第7次：
特征值max(x7)： 13.813746465256594
特征向量   y7： [-0.5331  1.      0.0791]
          x8： [16.9431  2.9778 -1.5049]
第8次：
特征值max(x8)： 16.94313702089086
特征向量   y8： [ 1.      0.1758 -0.0888]
          x9： [-9.1556 12.7594 -0.5061]
第9次：
特征值max(x9)： 12.759391141602638
特征向量   y9： [-0.7176  1.     -0.0397]
          x10： [18.8534  1.1862 -3.2332]
第10次：
特征值max(x10)： 18.853433942481207
特征向量   y10： [ 1.      0.0629 -0.1715]
          x11： [-10.4794  11.9093  -1.3552]
第11次：
特征值max(x11)： 11.909254019299365
特征向量   y11： [-0.8799  1.     -0.1138]
          x12： [20.5655 -0.4517 -4.3593]
第12次：
特征值max(x12)： 20.565503848485033
特征向量   y12： [ 1.    -0.022 -0.212]
          x13： [-11.4536  11.2263  -1.7117]
第13次：
特征值max(x13)： -11.453578681607251
特征向量   y13： [ 1.     -0.9802  0.1494]
          x14： [-21.6322   1.8797   4.9031]
第14次：
特征值max(x14)： -21.632241268767345
特征向量   y14： [ 1.     -0.0869 -0.2267]
          x15： [-12.1825  10.6713  -1.7728]
第15次：
特征值max(x15)： -12.182498255541232
特征向量   y15： [ 1.     -0.876   0.1455]
          x16： [-20.4899   2.8254   4.6436]
第16次：
特征值max(x16)： -20.48993757175785
特征向量   y16： [ 1.     -0.1379 -0.2266]
          x17： [-12.7434  10.2122  -1.6704]
第17次：
特征值max(x17)： -12.743448525598359
特征向量   y17： [ 1.     -0.8014  0.1311]
          x18： [-19.684    3.5255   4.3068]
第18次：
特征值max(x18)： -19.683997941414024
特征向量   y18： [ 1.     -0.1791 -0.2188]
          x19： [-13.189    9.8256  -1.4861]
第19次：
特征值max(x19)： -13.1889525203576
特征向量   y19： [ 1.     -0.745   0.1127]
          x20： [-19.0822   4.0697   3.9548]
第20次：
特征值max(x20)： -19.08221542867637
特征向量   y20： [ 1.     -0.2133 -0.2073]
          x21： [-13.5533   9.495   -1.2677]

………………

第148次：
特征值max(x148)： -15.96947870451083
特征向量   y148： [ 1.     -0.4457 -0.0653]
          x149： [-15.9683   7.1192   1.0421]
第149次：
特征值max(x149)： -15.96828831423285
特征向量   y149： [ 1.     -0.4458 -0.0653]
          x150： [-15.9694   7.118    1.0433]
第150次：
特征值max(x150)： -15.969403255924263
特征向量   y150： [ 1.     -0.4457 -0.0653]
          x151： [-15.9684   7.1191   1.0422]
第151次：
特征值max(x151)： -15.96835897584962
特征向量   y151： [ 1.     -0.4458 -0.0653]
          x152： [-15.9693   7.1181   1.0432]
第152次：
特征值max(x152)： -15.969337068297888
特征向量   y152： [ 1.     -0.4457 -0.0653]
估计的特征向量: [ 1.     -0.4457 -0.0653]
估计的特征值: -15.96835897584962

b. 带原点移位的乘幂法

          x1： [-14.   3.  -3.]
第1次：
特征值max(x1)： -14.0
特征向量   y1： [ 1.     -0.2143  0.2143]
          x2： [-28.1429  11.8571   1.    ]
第2次：
特征值max(x2)： -28.142857142857142
特征向量   y2： [ 1.     -0.4213 -0.0355]
          x3： [-30.6701  13.599    1.9137]
第3次：
特征值max(x3)： -30.67005076142132
特征向量   y3： [ 1.     -0.4434 -0.0624]
          x4： [-30.9398  13.7852   2.0116]
第4次：
特征值max(x4)： -30.93975504799735
特征向量   y4： [ 1.     -0.4455 -0.065 ]
          x5： [-30.9661  13.8033   2.0211]
第5次：
特征值max(x5)： -30.966052915940043
特征向量   y5： [ 1.     -0.4458 -0.0653]
          x6： [-30.9686  13.8051   2.0221]
第6次：
特征值max(x6)： -30.968592776449142
特征向量   y6： [ 1.     -0.4458 -0.0653]
          x7： [-30.9688  13.8052   2.0221]
第7次：
特征值max(x7)： -30.968837849183018
特征向量   y7： [ 1.     -0.4458 -0.0653]
          x8： [-30.9689  13.8053   2.0221]
第8次：
特征值max(x8)： -30.968861494278652
特征向量   y8： [ 1.     -0.4458 -0.0653]
          x9： [-30.9689  13.8053   2.0222]
第9次：
特征值max(x9)： -30.968863775583827
特征向量   y9： [ 1.     -0.4458 -0.0653]
          x10： [-30.9689  13.8053   2.0222]
第10次：
特征值max(x10)： -30.968863995686508
特征向量   y10： [ 1.     -0.4458 -0.0653]
估计的特征值: -15.96886377558382

本文参与腾讯云自媒体同步曝光计划，分享自作者个人站点/博客。

原始发表：2024-01-08，如有侵权请联系 cloudcommunity@tencent.com 删除

程序