传递函数是动态系统的输入到输出的表示。频域(相对于时域)的一个好处是将微分方程转换为代数方程。这些代数方程可以重新排列并转换回时域以获得解决方案,或者进一步与其他传递函数相结合以创建更复杂的系统。创建传递函数的第一步是用拉普拉斯变换来转换微分方程。
传递函数G(s)将输入U(s)与输出Y(s)联系起来:
使用Python进行拉普拉斯变换
要得到传递函数,首先要对系统的输出与输入量进行拉普拉斯变换。
这里要使用到Sympy,Sympy是一个支持符号计算的Python库,符号计算是指处理数学对象的计算,与数值计算相比,在符号计算中,数学对象是精确表示的,而不是近似的,未计算的数学表达式会以符号形式保留。
典型的符号计算有:
Sympy同样支持拉普拉斯变换(laplace transform)、拉普拉斯逆变换(inverse laplace transform)、部分分数展开( apart)、多项式展开(expand)和多项式根(roots)等系统建模分析中常用的运算。
下面以常见的拉普拉斯变换为例进行介绍:
import sympy as sym
from sympy.abc import s,t,x,y,z,a
from sympy.integrals import laplace_transform
from sympy.integrals import inverse_laplace_transform
单位阶跃
U = laplace_transform(1,t,s)
U[0]
sin(at)
计算案例
这里假设一包含阶跃和斜坡函数的输入,其图像为:
系统传递函数是一个具有两个极点(分母根)和一个零点(分子根)的稳定系统:
基于Sympy的符号解
首先将阶跃和斜坡信号分解为三个单独的函数。分别计算系统对这三个输入的响应,然后对信号求和。
这里利用到了延迟定理(t域平移定理):
求解代码:
import sympy as sym
from sympy.abc import s,t,x,y,z
import numpy as np
from sympy.integrals import inverse_laplace_transform
import matplotlib.pyplot as plt
# 定义输入
# 阶跃
U1 = 2/s*sym.exp(-s)
# 3s处的斜坡
U2 = -1/s**2*sym.exp(-3*s)
# 5s处斜坡
U3 = 1/s**2*sym.exp(-5*s)
# 定义传递函数
G = 5*(s+1)/(s+3)**2
# 计算输出
Y1 = G * U1
Y2 = G * U2
Y3 = G * U3
# 拉普拉斯逆变换
u1 = inverse_laplace_transform(U1,s,t)
u2 = inverse_laplace_transform(U2,s,t)
u3 = inverse_laplace_transform(U3,s,t)
y1 = inverse_laplace_transform(Y1,s,t)
y2 = inverse_laplace_transform(Y2,s,t)
y3 = inverse_laplace_transform(Y3,s,t)
print('y1')
print(y1)
# 创建绘图数据
tm = np.linspace(0,8,100)
us = np.zeros(len(tm))
ys = np.zeros(len(tm))
# 数值化u和y
for u in [u1,u2,u3]:
for i in range(len(tm)):
us[i] += u.subs(t,tm[i])
for y in [y1,y2,y3]:
for i in range(len(tm)):
ys[i] += y.subs(t,tm[i])
# 结果图
plt.figure()
plt.plot(tm,us,label='u(t)')
plt.plot(tm,ys,'--r',label='y(t)')
plt.legend()
plt.xlabel('Time')
plt.show()
基于GEKKO的数值解法
符号解法的一个替代方法是在时域中用数值法计算响应。传递函数必须首先被转化为微分方程。
from gekko import GEKKO
import numpy as np
import matplotlib.pyplot as plt
m = GEKKO()
# 创建数据点
nt = 81
m.time = np.linspace(0,8,nt)
# 定义输入
ut = np.zeros(nt)
# 阶跃
ut[11:31] = 2.0
# 斜坡
for i in range(31,51):
ut[i] = ut[i-1] - 0.1
# 定义模型
u = m.Param(value=ut)
ud = m.Var()
y = m.Var()
dydt = m.Var()
m.Equation(ud==u)
m.Equation(dydt==y.dt())
m.Equation(dydt.dt() + 6*y.dt() + 9*y==5*ud.dt()+5*u)
# 计算设置
m.options.IMODE=7
m.options.NODES=4
m.solve(disp=False)
#结果图
plt.figure()
plt.plot(m.time,u.value,label='u(t)')
plt.plot(m.time,y.value,'--g',label='y(t)')
plt.legend()
plt.xlabel('Time')
plt.show()
符号(解析)解的优点是高度准确,不依赖数值方法来近似解。此外,解决方案是一个紧凑的形式,可以用于进一步的分析。符号解仅限于输入函数和系统传递函数可以用拉普拉斯变换的情况。拉普拉斯变换的符号解法也无法用于非线性或复杂的系统
数值解可以处理成千上万或数百万的非线性关系的方程。数值解法的缺点是它是对真实解法的一种近似,可能有不准确的地方。另一个缺点是,求解器可能会无法收敛。
联系客服