Вычислить матрицу перехода состояний в Python - PullRequest
Новая
       52

Вычислить матрицу перехода состояний в Python

0 голосов
/ 03 января 2019

Существует ли прямой способ вычисления матрицы перехода состояний (т. Е. E ^ (A * t), где A - матрица)?Я планировал рассчитать это следующим образом:

x(t) is the e^(A*t) and x(0)=1

, но не получилось:

failed when calculate inverse laplace tansform

И если я вычислю напрямую A * tсначала, а затем используйте expm(), он по-прежнему не может работать, поскольку в expm() не должно быть никаких переменных.

Надеюсь, я ясно проиллюстрирую мою проблему:)

РЕДАКТИРОВАТЬ: Вот кодЯ думаю, должно быть полезно решить мою проблему: 

import numpy as np
import sympy
import scipy
from scipy.integrate import quad


Ts=0.02
s=sympy.symbols('s')
t=sympy.symbols('t')
T0=np.matrix([[1,0,0],
              [0,1,0],
              [0,-1,1]])
M0=np.matrix([[1.735,0.15851,0.042262],
              [0.123728,0.07019322,0.02070838],
              [0.042262,0.0243628,0.014375212]])
F0=np.matrix([[-22.915,0,0],
              [0,-0.00969,0.00264],
              [0,0.00264,-0.00264]])
N0=np.matrix([[0,0,0],
              [0,1.553398,0],
              [0,0,0.4141676]])
G0=np.matrix([[11.887],[0],[0]])
Ky=np.matrix([1.0121,4.5728,6.3652,0.9117,1.5246,0.9989])
A21=T0*(M0.I)*N0*(T0.I)
A22=T0*(M0.I)*F0*(T0.I)
Z=np.zeros((3,3))
Y=(np.matrix([0,0,0])).T
by1=np.row_stack((Z,A21))
by2=np.row_stack((np.identity(3),A22))
A=np.column_stack((by1,by2))
G=scipy.linalg.expm(A*Ts)


B2=T0*(M0.I)*G0
B=np.row_stack((Y,B2))
S1=sympy.Matrix((s*np.identity(6))-A)
S2=S1.inv()
S=S2
for (i,j), orinm in scipy.ndenumerate(S2):
    S[i,j]=sympy.inverse_laplace_transform(orinm, s, t)

#integral
H=np.zeros(S2.shape, dtype=float)
for (i,j),func_sympy in scipy.ndenumerate(S2):
    func = sympy.lambdify( (t),func_sympy, 'math')
    H[i,j] = quad(func, 0, 0.02)[0]
print(H)

1 Ответ

0 голосов
/ 03 января 2019

Вы можете напрямую рассчитать экспоненциальную матрицу, используя scipy. 

import numpy as np
from scipy.linalg import expm

A = np.random.random((10, 10))
exp_A = expm(A)

Документация для этого здесь Используется приближение Паде.

Вот пример использования матрицы тождества 2x2. 

>>> expm(np.eye(2))
array([[2.71828183, 0.        ],
       [0.        , 2.71828183]])

Если вам нужна матричная экспонента символической матрицы (согласно вашему комментарию), то вы можете сделать это с помощью Sympy: 

t = sympy.symbols('t')
A = sympy.Matrix([[t, 0], [0, t]])
>>> sympy.exp(A)
Matrix([
[expt(t), 0],
[0,  exp(t)]])
...