如何使用python-controlIVP问题模拟系统传递函数的时间响应？

问题描述

我正在尝试演示如何使用系统传递函数和python-control模块的定义来“求解”（模拟求解）微分方程初值问题（IVP）。事实是我真的是控制方面的新手。

我有一个简单的微分作为示例：y'' - 4y' + 13y = 0，具有以下初始条件：y(0) = 1和y'(0) = 0。

我可以手动实现此传递函数： Y(s) = (s - 4)/(s^2 - 4*s + 13)。

因此，在Python中，我正在编写这段代码（请注意，y_ans是该差分IVP的答案为seen here）：

import numpy as np
import control as ctl
import matplotlib.pyplot as plt

t = np.linspace(0.,1.5,100)
sys = ctl.tf([1.,-4.],[1.,-4.,13.])
T,yout,_ = ctl.forced_response(sys,T=t,X0=[1,0])

y_ans = lambda x: 1/3*np.exp(2*x)*(3*np.cos(3*x) - 2*np.sin(3*x))

plt.plot(t,y_ans(t),'-.',color='gray',alpha=0.5,linewidth=3,label='correct answer')
plt.plot(T,'r',label='simulated')
plt.legend()

此代码为我提供了这张图：

但是当我在yout前面插入一个负号时，我得到了一个匹配，因为它是这样的：

plt.plot(T,-yout,label='simulated') ### yout with negative sign

我在做什么错？ Python控制文档对我来说不是很清楚。另外，我不知道我对X0的{{1}}参数的解释是否正确。可以按照我的意愿去做吗？

欢迎任何有此想法的人贡献力量。

编辑

设置control.forced_response可以得到以下图形：

解决方法

我认为这里最好的做法是将您的系统转换为状态空间并看看发生了什么（有许多可能的状态空间表示）：

sys_tf = ctl.tf([1.,-4.],[1.,-4.,13.])
sys_ss = ctl.tf2ss(sys_tf)
print(sys_ss)

输出：

A = [[ 4.00000000e+00  1.30000000e+00]
     [-1.00000000e+01  1.33226763e-15]]

B = [[-1.]
     [ 0.]]

C = [[-1.  -0.4]]

D = [[0.]]

我们想要找到 x(0) 使得 y(0) = Cx(0) = 1 和 y'(0) = CAx(0) = 0。

我们可以写出这些方程并手动求解，也可以使用线性代数：

A = np.vstack([sys_ss.C,sys_ss.C @ sys_ss.A])
b = np.array([[1],[0]])
x0 = np.linalg.solve(A,b)
print(x0)

给出：

[[-1.00000000e+00]
 [-1.36642834e-15]]

因此，这应该有效：

T,yout = ctl.forced_response(sys_ss,T=t,X0=[-1,0])

此外，由于您只对初始条件（即 u(t)=0）的瞬态响应感兴趣，您可以使用 initial_response 函数：

T,yout = ctl.initial_response(sys_ss,0])

由于@LutzLehmann的评论，我一直在思考“两次编码”的含义。因此，回到平方，我意识到此传递函数包含了输入（时间或斜坡）和初始条件。它实际上是一个输出。我需要某种拉普拉斯逆变换，或者，正如我开始思考的那样，我只需要按原样对其进行仿真，而无需进一步的信息。

因此，我设法使用了脉冲输入（拉普拉斯变换等于1），并且能够获得恰好是我的tf实时模拟的输出。

import numpy as np
import control as ctl
import matplotlib.pyplot as plt

t = np.linspace(0.,1.5,100)
sys = ctl.tf([1.,13.])
T,yout = ctl.impulse_response(sys,T=t) # HERE is what I wanted

y_ans = lambda x: 1/3*np.exp(2*x)*(3*np.cos(3*x) - 2*np.sin(3*x))

plt.plot(t,y_ans(t),'-.',color='gray',alpha=0.5,linewidth=3,label='correct answer')
plt.plot(T,yout,'r',label='simulated')
plt.legend()

现在我想我可以展示如何使用python-control间接模拟微分方程的答案。 :-D

differential-equations python python-control transfer-function