跳至主要內容

Matlab 模型仿真

大约 7 分钟

Matlab 模型仿真

系统数学模型

系统数学模型函数用于创建系统的数学模型 (传递参数), 函数参数即模型的参数, 返回一个模型对象

tf 函数

  • 功能: 创建传递函数模型
  • 函数原型 sys = tf(num, den)
  • 参数
    • num 实数数组
      传递函数分子多项式系数, 从高次项到低次项
    • den 实数数组
      传递函数分母多项式系数, 从高次项到低次项
  • 返回值
    返回创建的传递函数模型
  • 注意
    传递函数的分母次数必须大于分子, 即 size(den)>size(num)
  • 示例
num = [1 2];
den = [1 2 3];
sys = tf(num, den)

示例代码将创建传递函数

G(s)=1s+2s2+2s+3 G(s)=\frac{1s+2}{s^2+2s+3}

zpk 函数

  • 功能: 创建零极点的传递函数模型
  • 函数原型 sys = zpk(z, p, k)
  • 参数
    • z 复数数组
      传递函数的所有零点
    • p 复数数组
      传递函数的所有极点
    • K 实数
      传递函数的增益
  • 返回值
    返回创建的传递函数模型
  • 注意
    • 允许传入复数极点或零点, 但必须同时传入两个共轭极点/零点
    • 对于多重根重复输入即可
    • 注意, 如下方公式所示, zpk 得到的是 G2G_2, 与教材上的传递函数 G1G_1 形式不同, 特别注意 K2K1K_2\neq K_1

G1(s)=K1(Tis+1)(Tis+1)G2(s)=K2(szi)(Tispi) \begin{split}G_1(s)=\frac{K_1\prod(T_i s+1)}{\prod(T_i' s+1)}\\ G_2(s)=\frac{K_2\prod(s-z_i)}{\prod(T_i' s-p_i)} \end{split}

  • 示例
z = [4 5];
p = [1 2 3];
sys = zpk(z, p, 6);

示例代码将创建传递函数

G(s)=6(s4)(s5)(s1)(s2)(s3) G(s)=\frac{6(s-4)(s-5)}{(s-1)(s-2)(s-3)}

ss 函数

  • 功能: 创建系统状态空间模型
  • 函数原型 ss_sys = ss(A, B, C, D)
  • 参数
    • A 状态矩阵 nx×nxn_x\times n_x 的矩阵
    • B 输入矩阵 nx×nun_x\times n_u 的矩阵
    • C 输出矩阵 ny×nxn_y\times n_x 的矩阵
    • D 直接传输矩阵 ny×nun_y\times n_u 的矩阵
  • 返回值
    返回值为系统的状态空间模型对象
  • 注意
    • 通过调用 ss_sys(a, b) 可以获取第 a 个输出变量到第 b 个输入变量的单输入单输出的传递函数对象
    • 除了根轨迹, 状态模型可直接使用模型仿真模型分析函数分析各个输入与输出间的关系
    • 通过对象的成员变量 A, B, C, D 可访问空间状态方程的四个矩阵

模型连接

series 函数

  • 功能: 串联两个数学模型
  • 函数原型 sys = serise(sys1, sys2)
  • 参数
    • sys1, sys2 数学模型
      参与串联的两个模型
  • 返回值
    返回串联后的函数模型
  • 示例
sys1 = tf([2 9], [4 7 2]);
sys2 = tf([1 6], [1 7 1]);
sys = series(sys1, sys2);

示例将串联两个传递函数, 得到

G(s)=2s+94s2+7s+2×s+61s2+7s+1 G(s)=\frac{2s+9}{4s^2+7s+2}\times\frac{s+6}{1s^2+7s+1}

parallel 函数

  • 功能: 并联两个数学模型
  • 函数原型 sys = parallel(sys1, sys2)
  • 参数
    • sys1, sys2 数学模型
      参与并联的两个模型
  • 返回值
    返回并联后的函数模型
  • 示例
sys1 = tf([2 9], [4 7 2]);
sys2 = tf([1 6], [1 7 1]);
sys = parallel(sys1, sys2);

示例将并联两个传递函数, 得到

G(s)=2s+94s2+7s+2+s+61s2+7s+1 G(s)=\frac{2s+9}{4s^2+7s+2}+\frac{s+6}{1s^2+7s+1}

feedback 函数

  • 功能: 反馈连接两个数学模型
  • 函数原型 sys = parallel(sys1, sys2, sign=-1)
  • 参数
    • sys1 数学模型
      前向通道上的数学模型
    • sys2 数学模型
      反馈通道上的数学模型
    • sign 1-1
      反馈类型, 取 1 时为正反馈, 取 -1 时为负反馈
  • 返回值
    返回反馈连接后的函数模型
  • 示例
sys1 = tf([2 9], [2 6 5]);
sys2 = tf(1, 1);
sysp = feedback(sys1, sys2, 1);

示例将传递函数 Gk=2s+92s2+6s+5G_k=\frac{2s+9}{2s^2+6s+5} 进行单位正反馈, 得到

G(s)=Gk(s)1Gk(s) G(s)=\frac{G_k(s)}{1-G_k(s)}

模型仿真

对于模型仿真图像, 可使用右键菜单, 以显示系统的响应特征

step 函数

  • 功能: 获取模型的单位阶跃响应
  • 函数原型 [Y T] = step(sys, t)
  • 参数
    • sys 数学模型
      用于仿真的模型
    • t 实数数组
      采样响应结果的时间轴
  • 返回值
    • Y 实数数组
      响应结果
    • T 实数数组
      响应结果的时间轴
  • 示例
sys = tf([2 9], [2 6 5]);
t = 0:0.05:30;
y = step(sys, t);
plot(t, y);

impluse 函数

  • 功能: 获取模型的单位冲激响应
  • 函数原型 [Y T] = impulse(sys, t)
  • 参数
    • sys 数学模型
      用于仿真的模型
    • t 实数数组
      采样响应结果的时间轴
  • 返回值
    • Y 实数数组
      响应结果
    • T 实数数组
      响应结果的时间轴
  • 示例
    step 函数 类似

lsim 函数

  • 功能: 获取模型对于给定信号的响应
  • 函数原型 [Y T] = lsim(sys, u, t)
  • 参数
    • sys 数学模型
      用于仿真的模型
    • u 实数数组
      用于仿真的信号, 采样时间与 t 对应
    • t 实数数组
      采样响应结果的时间轴
  • 返回值
    • Y 实数数组
      响应结果
    • T 实数数组
      响应结果的时间轴
  • 示例
sys = tf([2 9], [2 6 5]);
t = 0:0.05:4 * pi;
u = sin(t);
y = lsim(sys, u, t);
plot(t, y);

模型分析

对于模型分析图像, 可使用右键菜单, 以显示系统特征图中的特殊点

nyqusit 函数

  • 功能: 绘制指定系统的 Nyquist 图, 或获取其 Nyquist 图上的点
  • 函数原型 [re im w] = nyquist(sys)
  • 参数
    • sys 数学模型
      用于分析的数学模型
  • 返回值
    • re, im 实数数组
      Nyquist 图的实部与虚部
    • w 实数数组
      绘制 Nyquist 图时自动选择的频率
  • 注意
    • nyquist 函数没有输出被接收时, 将会自动绘制 Nyquist 图
    • 此函数通常配合 Nyquist 判据使用, 自动输出时将标记曲线方向与 (0,1)(0,-1)
  • 示例
sys = tf([1 2], [1 6 -5]);
nyquist(sys);

示例将绘制给定传递函数 G(s)=s+2s2+6s5G(s)=\frac{s+2}{s^2+6s-5} 的 Nyquist 图

bode 函数

  • 功能: 绘制指定系统的 Bode 图, 或获取其 Bode 图的幅频与相频曲线上的点
  • 函数原型 [mag phase w] = bode(sys)
  • 参数
    • sys 数学模型
      用于分析的数学模型
  • 返回值
    • mag 三维实数数组
      Bode 图的幅频曲线
    • phase 三维实数数组
      Bode 图的相频曲线
    • w 一维实数数组
      Bode 图曲线的频率轴
  • 注意
    • bode 函数没有输出被接收时, 将会自动绘制 Bode 图
    • magphase 为 1x1xn 的三维数组, 可使用 phase = phase(:) 的方法转换为一维
    • 直接接收得到的 mag 为原始振幅, 还需要 mag = 20 * log10(mag(:)) 将其转化为分贝
  • 示例
sys = tf([1 2], [1 6 -5]);
bode(sys);

示例将绘制给定传递函数 G(s)=s+2s2+6s5G(s)=\frac{s+2}{s^2+6s-5} 的 Bode 图

margin 函数

  • 功能: 求出给定开环系统的幅值裕度与相位裕度
  • 函数原型 [Gm Pm Wcg Wcp] = margin(sys)
  • 参数
    • sys 数学模型
      用于分析的数学模型
  • 返回值
    • Gm 实数
      幅值裕度 KgK_g
    • Wcg 实数
      相位交界频率 ωg\omega_g, 即幅值裕度取值对应的频率 (相频曲线第一次穿越 φ=180°\varphi=-180\degree)
    • Pm 实数
      相位裕度 γ\gamma
    • Wcp 实数
      剪切频率 ωc\omega_c, 即相位裕度取值对应的频率 (幅频曲线第一次穿越 20lgGH=0dB20\lg|GH|=0dB)
  • 注意
    • margin 函数没有输出被接收时, 将会自动绘制 Bode 图, 并且在图上标记出 GmPm
    • Gm 为原始振幅, 还需要 Gm = 20 * log10(Gm) 将其转化为分贝
  • 示例
sys = tf(3.5, [1, 2, 3, 2]);
margin(sys)

示例将绘制给定开环传递函数 G(s)=3.5s3+2s2+3s+2G(s)=\frac{3.5}{s^3+2s^2+3s+2} 的 Bode 图

rlocus 函数

建模技巧

非零初始状态的系统响应

以弹簧阻尼质量块系统为例, 系统有微分方程

mx¨+bx˙+kx=0 m\ddot{x}+b\dot{x}+kx=0

假设质量块相对平衡位置有偏置, 则有 x(0)=x0,x˙(0)=0x(0)=x_0,\dot{x}(0)=0, 根据微分特性定理, 此时方程的 Laplace 变换为

m[s2Xsx(0)x˙(0)]+b[sXx(0)]+kX=0mx(0)s+mx˙(0)+bx(0)ms2+bs+k=X(s) \begin{split} m[s^2X-sx(0)-\dot{x}(0)]+b[sX-x(0)]+kX&=0\\ \frac{mx(0)s+m\dot{x}(0)+bx(0)}{ms^2+bs+k}&=X(s) \end{split}

X(s)X(s) 视为一个传递函数, 而单位冲激函数的 Laplace 变换为 11, 因此通过求传递函数 X(s)X(s)单位冲激响应即可得到系统在非零状态下的响应

斜坡响应的仿真

单位斜坡响应可理解为单位阶跃响应的积分

  • 传统方法
  • 改变传递函数
  • 增加状态方程的输出变量