用Python验证Cayley-Hamilton定理：从理论到代码的完整实现

最新推荐文章于 2026-06-30 20:27:17 发布

原创

最新推荐文章于 2026-06-30 20:27:17 发布 · 852 阅读

用Python验证Cayley-Hamilton定理：从理论到代码的完整实现

你是否曾觉得线性代数中的那些优美定理，比如Cayley-Hamilton定理，虽然理论上令人着迷，但总感觉离实际的编程世界有些遥远？作为一个喜欢用代码解决问题的开发者，我也有过同样的困惑。直到有一次，我在优化一个涉及矩阵幂运算的算法时，偶然发现如果能利用这个定理，计算量可以大幅减少。这让我意识到，将抽象的数学定理转化为具体的代码，不仅能加深对理论的理解，更能解锁解决实际工程问题的新思路。这篇文章就是为你——一位对数学有好奇心、更热衷于动手实践的Python开发者——准备的。我们将完全避开繁琐的纯符号证明，直接深入NumPy的世界，用代码亲手“触摸”这个定理，并探索它在计算中的实际威力。

1. 定理的直观理解与Python环境搭建

在开始敲代码之前，我们得先搞清楚Cayley-Hamilton定理到底在说什么。用最直白的话讲：任何一个方阵，都是它自己特征方程的“根”。这听起来有点反直觉，矩阵怎么能成为方程的“根”呢？这里的“根”不是指代数值，而是指当你把矩阵本身代入它的特征多项式后，结果会神奇地得到一个零矩阵。

假设我们有一个 n×n 的矩阵 A。它的特征多项式通常写作： χ(λ) = det(λI - A) = λ^n + c_{n-1}λ^{n-1} + ... + c_1λ + c_0 其中 λ 是特征值，I 是单位阵，c_i 是系数。定理告诉我们，如果把多项式里的变量 λ 替换成矩阵 A 本身，那么计算结果： χ(A) = A^n + c_{n-1}A^{n-1} + ... + c_1A + c_0I 一定等于零矩阵 O。

注意：这里 c_0I 中的 c_0 是常数项，它需要乘以单位矩阵 I 才能与前面的矩阵项相加，这是矩阵多项式运算的基本规则。

为了验证它，我们需要一个强大的数值计算环境。Python的SciPy生态是绝佳选择。首先确保你的环境已经就绪：

# 使用pip安装必要的库
pip install numpy scipy matplotlib ipython

我们将主要依赖 numpy 进行核心的矩阵运算，scipy.linalg 提供一些更高级的线性代数函数。下面是一个快速的环境检查脚本，你可以新建一个Jupyter Notebook或Python文件来跟随操作：

import numpy as np
import scipy.linalg as la
import sys

print(f"Python 版本: {sys.version}")
print(f"NumPy 版本: {np.__version__}")
print(f"SciPy 版本: {la.__version__}")

# 快速测试矩阵创建
A = np.array([[2, 1], [1, 2]])
print("测试矩阵 A:\n", A)

如果一切正常，你会看到相关库的版本信息和一个小矩阵。准备工作就绪，我们可以进入实战了。

2. 核心实现：计算特征多项式与矩阵代入

验证定理的第一步，是获取给定矩阵的特征多项式系数。在NumPy中，我们可以通过 np.poly 函数轻松得到特征多项式的系数，但这里我们更倾向于从原理出发，分步实现，以加深理解。

2.1 手动计算特征多项式系数

虽然 np.poly(A) 能直接给出结果，但了解其背后的计算过程很有价值。特征多项式 det(λI - A) 的系数可以通过矩阵的迹（trace）和更高级的幂和关系来推导，但对于代码验证，最直接的方法是构造特征矩阵并符号化（或数值化）计算行列式。不过，对于中小型矩阵，我们可以利用 scipy.linalg.eigvals 先求出特征值，因为特征值 λ_i 就是特征方程的根，多项式系数可以通过这些根来重构。

def char_poly_coeffs_from_eigvals(A):
    """
    通过特征值计算矩阵A的特征多项式系数。
    返回的系数列表coeffs满足：
    λ^n + coeffs[0]*λ^{n-1} + ... + coeffs[n-2]*λ + coeffs[n-1] = 0
    """
    # 计算特征值
    eigvals = la.eigvals(A)
    n = A.shape[0]
    # 根据根构造多项式系数：多项式 = (λ - λ1)(λ - λ2)...(λ - λn)
    # np.poly函数接收根序列，返回最高次项系数为1的多项式系数
    coeffs_full = np.poly(eigvals)  # 返回 [1, c_{n-1}, ..., c_0]
    # 去掉最高次的1，返回 [c_{n-1}, ..., c_0]
    return coeffs_full[1:]

# 测试一个2x2矩阵
A_2x2 = np.array([[4, -2], [1, 1]])
coeffs = char_poly_coeffs_from_eigvals(A_2x2)
print("矩阵 A_2x2 的特征多项式系数 (从c_{n-1}到c_0):", coeffs)
# 验证：特征多项式应为 λ^2 - 5λ + 6 = 0
# 系数应为 [-5, 6]

这种方法直观，但依赖于特征值计算的数值精度。对于更大的矩阵或病态矩阵，特征值求解可能引入误差，这是我们后续验证时需要留意的地方。