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

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

在这里插入图片描述

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

一、乘幂法

  • 给定一个方阵 A A A,如果 λ \lambda λ A A A 的最大特征值,而 x x x 是对应的特征向量
  • 那么对于任意非零向量 v v v,迭代 A n v A^n v Anv 的结果会趋近于 λ n x \lambda^n x λnx。当 n n n 越来越大时, A n v A^n v Anv 的方向会趋近于 x x x 的方向,而其模(长度)则会趋近于 λ n \lambda^n λn 的绝对值。
    在这里插入图片描述

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

二、乘幂法的加速

1. 天书

在这里插入图片描述
在这里插入图片描述

2. 人话

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

a. 基本乘幂法

  1. 选择初始向量 x 0 x_0 x0:通常选择一个非零向量作为初始向量,其选择可能影响到迭代的收敛速度。

  2. 迭代计算 A x k A x_k Axk:对于每一次迭代 k k k,计算 A x k A x_k Axk

  3. 规范化 A x k A x_k Axk:通过规范化 A x k A x_k Axk(即将其除以其最大分量),得到新的向量 x k + 1 x_{k+1} xk+1

  4. 判断收敛条件:检查 ∥ x k + 1 − x k ∥ ∥ x k + 1 ∥ \frac{\|x_{k+1} - x_k\|}{\|x_{k+1}\|} xk+1xk+1xk 是否小于一个预定的收敛阈值,满足则停止迭代。

  5. 更新 x k x_k xk:将 x k + 1 x_{k+1} xk+1 赋值给 x k x_k xk,重复步骤2。

  6. 计算特征值:一旦迭代收敛,通过 A x k x k \frac{A x_k}{x_k} xkAxk 的比值来估计矩阵 A A A 的最大特征值。

b. 加速乘幂法的思想

  加速乘幂法的思想是通过平移矩阵 A A A,将其转化为 B = A − λ I B = A - \lambda I B=AλI 的形式。这个平移操作使得新的矩阵 B B B 的最大特征值与原矩阵 A A A 的最大特征值之间的差距变大,从而加速收敛。

c. 具体步骤

  • 平移操作: 对给定的矩阵 A A A,选择一个数值 λ \lambda λ,构造新矩阵 B = A − λ I B = A - \lambda I B=AλI

  • 基本乘幂法: 使用 B B B 进行上述乘幂法迭代计算,直至收敛

  • 计算特征值: 一旦迭代收敛,通过 ( B x k ) T x k ( x k ) T x k + λ \frac{(Bx_k)^T x_k}{(x_k)^T x_k} + \lambda (xk)Txk(Bxk)Txk+λ 来估计矩阵 A A A的最大特征值。

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.解析

  1. normalize_vector 函数:

    • 输入:一个向量 vec
    • 功能:计算向量的最大分量,并将向量规范化。
    • 输出:规范化后的向量和最大分量。
  2. power_iteration 函数:

    • 输入:
      • A:一个方阵(矩阵)。
      • x:初始向量。
      • max_iter:最大迭代次数,默认为1000。
      • tol:收敛阈值,默认为1e-4。
    • 功能:使用乘幂法迭代来估计矩阵的最大特征值及其对应的特征向量。
      • 计算矩阵 A 与向量 x 的乘积,得到 Ax。
      • 调用 normalize_vector 函数对 Ax 进行规范化,得到规范化后的向量和最大分量。
      • 打印每次迭代的结果,即特征值、特征向量。
      • 判断是否满足收敛条件,如果最大分量变化小于阈值 tol,则提前退出迭代。
      • 计算对应的特征值,更新最大分量,并继续迭代。
    • 输出:估计得到的特征值。
  3. 主程序部分:

    • 教材例题及课后题的矩阵 A、A1、A2、A3
    • 定义了初始向量 x0、系数 λ \lambda λ
    • 构造新矩阵 B = A − λ I B = A - \lambda I B=AλI
    • 调用 power_iteration 函数,传入矩阵 B B B 和初始向量进行乘幂法迭代。
    • B B B的特征值+ λ \lambda λ = A A 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
评论
添加红包

请填写红包祝福语或标题

红包个数最小为10个

红包金额最低5元

当前余额3.43前往充值 >
需支付:10.00
成就一亿技术人!
领取后你会自动成为博主和红包主的粉丝 规则
hope_wisdom
发出的红包

打赏作者

QomolangmaH

你的鼓励将是我创作的最大动力

¥1 ¥2 ¥4 ¥6 ¥10 ¥20
扫码支付:¥1
获取中
扫码支付

您的余额不足,请更换扫码支付或充值

打赏作者

实付
使用余额支付
点击重新获取
扫码支付
钱包余额 0

抵扣说明:

1.余额是钱包充值的虚拟货币,按照1:1的比例进行支付金额的抵扣。
2.余额无法直接购买下载,可以购买VIP、付费专栏及课程。

余额充值