python-3.x - 求解特征值问题的 Sympy 符号矩阵会返回一个空列表

Question

我是 Python 的新手。我会很感激你能提供的任何帮助。我正在尝试解决特征值问题。我有一个矩阵，我将其称为“data3”，其中包含未知变量“omega”的元素。矩阵的行列式将为我提供以“欧米茄”表示的特征多项式，从中我可以确定我需要的根或特征值。

以下是我的代码：

import sympy as sp
import numpy as np
import scipy as sc

omega=sp.Symbol('omega');
data1=sp.Matrix(sc.genfromtxt('./Z_real.csv',dtype=float,delimiter=',')); 
data2=sp.Matrix(sc.genfromtxt('./Z_imag.csv',dtype=float,delimiter=','));
data2a=data2*1j*omega;
data3=sp.Matrix(data1+data2a);

eqn=sp.Eq(sp.det(data3),0);
print(sp.solve(eqn));

如上所示，我导入了 data1 和 data2，然后将它们组合成一个复数矩阵（需要解决的 'omega'）。例如，当输出到 iPython 控制台时，我的 data3 矩阵如下所示：

Matrix([
[             0.536,  -1.119e-8*I*omega, -1.3558e-8*I*omega, -3.8778e-9*I*omega],
[ -1.119e-8*I*omega,              0.536, -3.8778e-9*I*omega, -1.3558e-8*I*omega],
[-1.3558e-8*I*omega, -3.8778e-9*I*omega,              0.536,  -1.119e-8*I*omega],
[-3.8778e-9*I*omega, -1.3558e-8*I*omega,  -1.119e-8*I*omega,              0.536]])

但是，代码的 sp.solve(eqn) 行为我返回了一个空列表。我期望看到的是矩阵多项式的根（来自行列式）。有人可以告诉我我做错了什么吗？此外，如果有人可以向我展示解决“欧米茄”的替代方法，那也很棒。如果您需要有关 data1 和 data2 的信息，以下是我用于测试的矩阵：

数据1

Matrix([
[0.536,   0.0,   0.0,   0.0],
[  0.0, 0.536,   0.0,   0.0],
[  0.0,   0.0, 0.536,   0.0],
[  0.0,   0.0,   0.0, 0.536]])

和数据2

Matrix([
[       0.0,  -1.119e-8, -1.3558e-8, -3.8778e-9],
[ -1.119e-8,        0.0, -3.8778e-9, -1.3558e-8],
[-1.3558e-8, -3.8778e-9,        0.0,  -1.119e-8],
[-3.8778e-9, -1.3558e-8,  -1.119e-8,        0.0]])

非常感谢您的宝贵时间。

标签： python-3.xsympydeterminants

Answer 1

感谢所有阅读帖子的人。我现在已经解决了这个问题。我发现不是直接求解矩阵的特征多项式，而是以多项式方程形式获得行列式（现在使用 Berkowitz 算法 - 没有特殊原因）要容易得多，使用 sympy Poly 提取多项式的系数，然后插入numpy 根函数中的根。奇迹般有效！请在下面找到代码（我导入的数据无关紧要，该方法适用于任何 sympy 矩阵）：

import numpy as np
import scipy as sc
import sympy as sp
import csv

omega=sp.Symbol('omega');
data_Za=sp.Matrix(sc.genfromtxt('./Z_matrix.csv',dtype=float,delimiter=','));
L_mata=sp.Matrix(sc.genfromtxt('./L_mat.csv',dtype=float,delimiter=','));
C_mata=sp.Matrix(sc.genfromtxt('./C_mat.csv',dtype=float,delimiter=','));
N_total=int(sc.genfromtxt('./N_total.csv',dtype=float,delimiter=','));

data_Z=data_Za*omega;
L_mat=L_mata*omega;
C_mat=C_mata*(1/omega);

data_matrix=data_Z+L_mat+C_mat;
det_eqn=data_matrix.berkowitz_det();

这就是我导入矩阵的方式，将它们相加以获得我想要的最终矩阵并计算行列式 - 直到这里没有什么新鲜事。

det_poly_coeff=sp.Poly(det_eqn*(omega**N_total),omega);
#print(det_poly_coeff.coeffs());

我使用的矩阵有点复杂，因为正如您可以看到我对 C_mat 所做的那样，我最终会在行列式（特征多项式）中得到omega即使在分母中也有变量的项。Sympy 的 Poly 函数不喜欢这样。因此，我最终将行列式乘以omega**N_total其中 N_total 是多项式方程的次数（以及分母中出现的最高 omega 次数）。

root_find_sq=np.roots(det_poly_coeff.coeffs());
#print(root_find_sq);

将系数插入 numpy 根给了我方程的根，这是我想要解决的值。

root_find=np.sqrt(abs(root_find_sq));

最后是平方根，因为求解的根是 foromega**2而不是直接 for omega。

希望其他可能遇到此问题的人发现这篇文章有帮助。干杯。:)