用Sympy实现特征向量的“手动”求解

Question

设置

from sympy import Matrix, symbols, eye, factor, solve, simplify, solveset
from fractions import Fraction
_lambda = symbols("lambda")
x1, x2 = symbols("x_1, x_2")

这是我正在处理的矩阵

S = Matrix([[5,4],[4,5]])
M = Matrix([[1,0],[0,4]])
H = M**Fraction(-1,2)*S*M**Fraction(-1,2)

我想做的下一步是这样的

eigen = []
for eigenvalue in solve((S-_lambda*M).det()):
    eigen.append((
                   eigenvalue,
                   solve((S-_lambda*M).subs(_lambda,eigenvalue) * Matrix([x1,x2]))
                ))

但同情找不到答案..。它可能只是不明白我要求什么，但我不知道如何更清楚地告诉它。

print(eigen)
> [(25/8 - sqrt(481)/8, {x_1: nan, x_2: nan}),
 (sqrt(481)/8 + 25/8, {x_1: nan, x_2: nan})]

Answer 1

我不太明白你想做什么。矩阵H未使用。看起来你想要矩阵S的特征值和特征向量，但是我不知道M矩阵在做什么。我假设你想要找到S的特征值和特征向量。

我们有矩阵S

In [60]: S = Matrix([[5,4],[4,5]])                                                                                                             

In [61]: S                                                                                                                                     
Out[61]: 
⎡5  4⎤
⎢    ⎥
⎣4  5⎦

我们首先找到它的特征值：

In [62]: l = Symbol('lambda')                                                                                                                  

In [63]: det(S - l*eye(2))                                                                                                                     
Out[63]: 
       2     
(5 - λ)  - 16

In [64]: l1, l2 = roots(_, l)                                                                                                                  

In [65]: l1                                                                                                                                    
Out[65]: 9

In [66]: l2                                                                                                                                    
Out[66]: 1

现在，这些特征值中的每一个都有一个相应的特征向量。实际上，它不是一个单一的特征向量，而是一个一维向量族。如果向量是[x1, x2]，那么我们就不能对x1和x2进行唯一的求解，因为有一个无穷大的解族。因此，使用solve给出了x1的x2参数化解。

In [67]: x1, x2 = symbols('x1, x2')                                                                                                            

In [68]: solve((S - l1*eye(2)) * Matrix([x1, x2]), [x1, x2])                                                                                   
Out[68]: {x₁: x₂}

这告诉我们，向量[x1, x2]将是提供x1 == x2的解决方案。因此，我们得到了特征向量[a, a]族，其特征值为l1 (9)：

In [69]: S*Matrix([1, 1])                                                                                                                      
Out[69]: 
⎡9⎤
⎢ ⎥
⎣9⎦

更直接地说，我们真正想要的是S - l*I的空空间，我们可以将其计算为

In [70]: (S - l1*eye(2)).nullspace()                                                                                                           
Out[70]: 
⎡⎡1⎤⎤
⎢⎢ ⎥⎥
⎣⎣1⎦⎦

In [71]: (S - l2*eye(2)).nullspace()                                                                                                           
Out[71]: 
⎡⎡-1⎤⎤
⎢⎢  ⎥⎥
⎣⎣1 ⎦⎦

nullspace方法返回作为空空间基础的向量列表。在这种情况下，这意味着它们是给定特征值的特征空间的基础。由于我们只有一个向量，所以它是一个一维特征空间，它由[1, 1] (或l2的[-1, 1] )的所有标量倍数组成。如果有重复的特征值，那么空空间就有可能有更大的维数。