Решение символьной матрицы Sympy для проблемы собственных значений возвращает пустой список

Я новичок в Python. Буду признателен за любую помощь, которую вы можете оказать. Я пытаюсь решить проблему с собственными значениями. У меня есть матрица, которую я называю «данные3» с элементами в терминах неизвестной переменной «омега». Определитель матрицы даст мне характеристический полином в терминах «омега», по которому я могу определить корни или собственные значения, которые мне нужны.

Вот мой код:

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, которые затем объединяю, чтобы сформировать матрицу комплексных чисел (с «омегой», которую мне нужно решить). Например, моя матрица data3 при выводе на консоль iPython выглядит так:

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, ниже приведены матрицы, которые я использую для своего теста:

data1

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]])

и data2

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.x sympy determinants
person Sathya    schedule 23.11.2018    source источник

Ответы (1)


arrow_upward
1
arrow_downward

Спасибо всем, кто прочитал этот пост. Теперь я исправил эту проблему. Вместо того, чтобы напрямую решать характеристический полином матрицы, я обнаружил, что было намного проще получить определитель (теперь с использованием алгоритма Берковица - без особой причины) в форме полиномиального уравнения, извлечь коэффициенты полинома с помощью sympy Poly, а затем подключить корни внутри функции numpy root. Работает как шарм! Пожалуйста, найдите приведенный ниже код (данные, которые я импортирую, не имеют значения, метод работает для любой симпозиумной матрицы):

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's Poly не любит этого. Поэтому я умножаю определитель на omega**N_total, где N_total - это степень полиномиального уравнения (а также наивысшая степень омеги, которая встречается в знаменателе).

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

Включение коэффициентов в несколько корней дает мне корни уравнения, которые являются значениями, которые я хотел решить.

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

Наконец, квадратный корень, потому что корни, которые были решены, были для omega**2, а не напрямую для omega.

Надеюсь, этот пост окажется полезным для всех, у кого есть эта проблема. Ваше здоровье. :)

person Sathya    schedule 04.12.2018