矩阵的特征值（eigenvalue）和特征向量（eigenvector）在很多应用中都具有重要的数学和物理意义。乘幂法（Power Iteration）是线性代数中一种重要的数值计算方法，用于估计矩阵的最大特征值及其对应特征向量的迭代算法，广泛应用于许多科学和工程领域。
本文将详细介绍乘幂法的基本原理和步骤，并给出其Python实现。
一、乘幂法
1. 天书
a. 乘幂法
本文仅考虑有唯一的主特征值情况，的主特征值不唯一情况不做介绍
b. 理论证明
c. 规范化的乘幂法
2. 人话
a. 数学原理
乘幂法（Power Iteration）是一种用于估计矩阵的最大特征值及其对应特征向量的迭代算法，基于以下的数学原理：
- 给定一个方阵
，如果
是
的最大特征值，而
是对应的特征向量
- 那么对于任意非零向量
，迭代
的结果会趋近于
。当
越来越大时，
的方向会趋近于
的方向，而其模（长度）则会趋近于
的绝对值。
- 理论证明详见上述天书部分
b. 基本步骤
- 选择初始向量
：通常选择一个非零向量作为初始向量，其选择可能影响到迭代的收敛速度。
- 迭代计算
：对于每一次迭代
，计算
。
- 规范化
：通过规范化
（即将其除以其最大分量），得到新的向量
。
- 判断收敛条件：检查
是否小于一个预定的收敛阈值，满足则停止迭代。
- 更新
：将
赋值给
，重复步骤2。
- 计算特征值：一旦迭代收敛，通过
的比值来估计矩阵
的最大特征值。
乘幂法的优点是它的简单性和易实现性。然而，它只能找到最大特征值，且其收敛速度可能相对较慢。对于一些特殊的矩阵，可能需要使用其他的迭代方法。
c. 注意事项
- 收敛性：
- 乘幂法只能估计最大特征值，并且其收敛速度取决于初始向量的选择以及特征值之间的差异。
- 对于对称正定矩阵，收敛是保证的。
- 复杂性：
- 乘幂法是一种简单且易于实现的方法，但对于某些情况下的矩阵，收敛速度可能较慢。
- 在某些情况下，可能需要使用其他迭代方法。
- 对称矩阵： 乘幂法在处理对称矩阵时效果更好，因为对称矩阵的特征向量是正交的。
- 扩展： 乘幂法的扩展形式包括反幂法、带有原点移位的乘幂法等。
3. 典例
4. 实现
代码语言：javascript复制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-4):
max_val0 = 0
for i in range(max_iter):
# 计算矩阵与向量的乘积
np.set_printoptions(precision=3, 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 x, eigenvalue
# 给定矩阵 A
A = np.array([[-4, 14, 0],
[-5, 13, 0],
[-1, 0, 2]])
A1 = np.array([[4, 2, 2],
[2, 5, 1],
[2, 1, 6]])
A2 = np.array([[7, 3, -2],
[3, 4, -1],
[-2, -1, 3]])
A3 = np.array([[-11, 11, 1],
[11, 9, -2],
[1, -2, 13]])
# 初始向量
x0 = np.array([1, 1, 1])
# 运行乘幂法迭代
# result_vector, eigenvalue = power_iteration(A, x0)
result_vector, eigenvalue = power_iteration(A1, x0)
# result_vector, eigenvalue = power_iteration(A2, x0)
# result_vector, eigenvalue = power_iteration(A3, x0)
print("估计的特征向量:", result_vector)
print("估计的特征值:", eigenvalue)-
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。 - 调用
power_iteration函数，分别传入不同的矩阵和初始向量进行乘幂法迭代。 - 打印估计得到的特征向量和特征值。
x1： [10 8 1]
第1次：
特征值max(x1)： 10
特征向量 y1： [1. 0.8 0.1]
x2： [ 7.2 5.4 -0.8]
第2次：
特征值max(x2)： 7.200000000000001
特征向量 y2： [ 1. 0.75 -0.1111]
x3： [ 6.5 4.75 -1.2222]
第3次：
特征值max(x3)： 6.499999999999998
特征向量 y3： [ 1. 0.7308 -0.188 ]
x4： [ 6.2308 4.5 -1.3761]
第4次：
特征值max(x4)： 6.23076923076923
特征向量 y4： [ 1. 0.7222 -0.2209]
x5： [ 6.1111 4.3889 -1.4417]
第5次：
特征值max(x5)： 6.1111111111111125
特征向量 y5： [ 1. 0.7182 -0.2359]
x6： [ 6.0545 4.3364 -1.4718]
第6次：
特征值max(x6)： 6.054545454545453
特征向量 y6： [ 1. 0.7162 -0.2431]
x7： [ 6.027 4.3108 -1.4862]
第7次：
特征值max(x7)： 6.0270270270270245
特征向量 y7： [ 1. 0.7152 -0.2466]
x8： [ 6.0135 4.2982 -1.4932]
第8次：
特征值max(x8)： 6.013452914798206
特征向量 y8： [ 1. 0.7148 -0.2483]
x9： [ 6.0067 4.2919 -1.4966]
第9次：
特征值max(x9)： 6.006711409395972
特征向量 y9： [ 1. 0.7145 -0.2492]
x10： [ 6.0034 4.2888 -1.4983]
第10次：
特征值max(x10)： 6.003351955307263
特征向量 y10： [ 1. 0.7144 -0.2496]
x11： [ 6.0017 4.2873 -1.4992]
第11次：
特征值max(x11)： 6.001675041876046
特征向量 y11： [ 1. 0.7143 -0.2498]
x12： [ 6.0008 4.2865 -1.4996]
第12次：
特征值max(x12)： 6.00083728718951
特征向量 y12： [ 1. 0.7143 -0.2499]
估计的特征向量: [ 1. 0.7143 -0.2499]
估计的特征值: 6.0016750418760475. 课后题
黄明游~《数值计算方法》 ~ P91
x1： [8 8 9]
第1次：
特征值max(x1)： 9
特征向量 y1： [0.8889 0.8889 1. ]
x2： [7.3333 7.2222 8.6667]
第2次：
特征值max(x2)： 8.666666666666666
特征向量 y2： [0.8462 0.8333 1. ]
x3： [7.0513 6.859 8.5256]
第3次：
特征值max(x3)： 8.525641025641026
特征向量 y3： [0.8271 0.8045 1. ]
x4： [6.9173 6.6767 8.4586]
第4次：
特征值max(x4)： 8.458646616541353
特征向量 y4： [0.8178 0.7893 1. ]
x5： [6.8498 6.5822 8.4249]
第5次：
特征值max(x5)： 8.424888888888889
特征向量 y5： [0.813 0.7813 1. ]
x6： [6.8147 6.5325 8.4074]
第6次：
特征值max(x6)： 8.407364422874025
特征向量 y6： [0.8106 0.777 1. ]
x7： [6.7963 6.5061 8.3981]
第7次：
特征值max(x7)： 8.398130137416075
特征向量 y7： [0.8093 0.7747 1. ]
x8： [6.7865 6.4921 8.3932]
第8次：
特征值max(x8)： 8.39322778520782
特征向量 y8： [0.8086 0.7735 1. ]
x9： [6.7812 6.4846 8.3906]
第9次：
特征值max(x9)： 8.390615458295574
特征向量 y9： [0.8082 0.7728 1. ]
x10： [6.7784 6.4806 8.3892]
第10次：
特征值max(x10)： 8.389220813597767
特征向量 y10： [0.808 0.7725 1. ]
x11： [6.777 6.4784 8.3885]
第11次：
特征值max(x11)： 8.388475552357708
特征向量 y11： [0.8079 0.7723 1. ]
估计的特征向量: [0.8079 0.7723 1. ]
估计的特征值: 8.389220813597769 x1： [8 6 0]
第1次：
特征值max(x1)： 8
特征向量 y1： [1. 0.75 0. ]
x2： [ 9.25 6. -2.75]
第2次：
特征值max(x2)： 9.25
特征向量 y2： [ 1. 0.6486 -0.2973]
x3： [ 9.5405 5.8919 -3.5405]
第3次：
特征值max(x3)： 9.54054054054054
特征向量 y3： [ 1. 0.6176 -0.3711]
x4： [ 9.5949 5.8414 -3.7309]
第4次：
特征值max(x4)： 9.594900849858357
特征向量 y4： [ 1. 0.6088 -0.3888]
x5： [ 9.6041 5.824 -3.7753]
第5次：
特征值max(x5)： 9.604074402125775
特征向量 y5： [ 1. 0.6064 -0.3931]
x6： [ 9.6054 5.8187 -3.7857]
第6次：
特征值max(x6)： 9.605429001813766
特征向量 y6： [ 1. 0.6058 -0.3941]
x7： [ 9.6056 5.8172 -3.7881]
第7次：
特征值max(x7)： 9.60557200236834
特征向量 y7： [ 1. 0.6056 -0.3944]
估计的特征向量: [ 1. 0.6056 -0.3944]
估计的特征值: 9.605429001813766
