矩阵对角化/约当标准形求解中的数值不准确性问题

• 你遇到的核心问题是矩阵存在重特征值或特征值极度接近，但通用数值计算工具误判其可对角化，本质上属于近不可对角化的病态矩阵——浮点精度下，矩阵的Jordan块结构被微小扰动掩盖，导致工具输出对角化结果而非真实的Jordan标准形。
• Numpy的np.linalg.eig、Mathematica的默认特征值算法均依赖QR迭代这类数值稳定但对重特征值敏感的方法，面对特征值分离度极小的矩阵时，会倾向于输出对角化结果，哪怕矩阵理论上不可对角化。

1. 符号计算（小规模精确矩阵首选） #

如果矩阵元素是精确有理数或符号表达式，用符号计算工具解析求解，彻底规避浮点误差：

• Mathematica中，确保矩阵以精确形式输入，调用JordanDecomposition[mat]（不要提前转成浮点型）；
• Python中使用sympy库的jordan_form函数，示例代码：

  import sympy as sp
  # 示例2阶Jordan块矩阵
  mat = sp.Matrix([[2, 1], [0, 2]])
  P, J = mat.jordan_form()
  print("Jordan标准形:\n", J)
  print("相似变换矩阵:\n", P)
  

  符号计算能精确识别Jordan块结构，不会被浮点精度干扰。

2. 数值Jordan标准形专用算法 #

针对大规模浮点矩阵，使用专门的Jordan分解算法：

• Python生态：用scipy.linalg.schur先做Schur分解，再通过检查Schur块的大小识别Jordan块；也可使用第三方库pyjordan（专注Jordan标准形数值计算）；
• Mathematica：尝试指定专用方法，如JordanDecomposition[mat, Method->"Jordan"]，或先通过**奇异值分解(SVD)**预处理降低矩阵病态程度后再分解。

3. 病态程度验证 #

先确认矩阵的病态性，判断误差来源：

• 用Numpy计算矩阵的Frobenius范数条件数：np.linalg.cond(mat, p='fro')，条件数越大，矩阵越病态，数值分解误差越大；
• 计算特征值分离度：用scipy.linalg.eigvals提取所有特征值，检查相邻特征值的差值，若差值远小于浮点精度（如<1e-12），说明特征值近似重合，通用算法会失效。

4. 手动构造广义特征向量（已知结构时） #

若通过理论分析已知Jordan块结构，可手动求解广义特征向量构建变换矩阵：

• 对重特征值λ，依次求解(mat - λ*I)^k x = 0且(mat - λ*I)^(k-1) x ≠ 0的向量，以此作为Jordan块的广义特征向量，组合成相似变换矩阵。

• 浮点数值计算中，严格不可对角化的矩阵几乎不存在——任何微小浮点扰动都可能让Jordan块“裂开”成对角矩阵，因此工具默认输出对角化结果是数值稳定的选择，但会掩盖原矩阵的理论结构；
• 只有通过高精度计算（如Python的mpmath库设置任意精度浮点）或符号计算，才能足够精确地捕捉到特征值重数和广义特征向量，得到准确的Jordan标准形。

内容的提问来源于stack exchange，提问作者PawatAka