如何用 python 代码复现 MATLAB simulink 的 PID

MATLAB
- 在 Simulink 里做以下设置
- MATLAB 脚本调用示例
python 实现
- 离散 PID 实现（并行形式）
Simulink 中两种 PID 结构（并联形式, I-形式）下连续/离散时域里积分增益 I 的表示
- 并联（Parallel） vs 理想（I-Form），use I*Ts 的勾选
- 什么时候用连续域 vs 离散域？

MATLAB

首先，MATLAB simulink 的 PID 长这样：
在这里插入图片描述

其中，PID 控制器以 1e-4 的步长运行，而物理模型用固定步长 1e-5 的求解器来集成。

这样，每仿真一次 StopTime = 1e-4，就得到了“10 步模型响应 + 1 次 PID 更新”之后的 $Z$ 和 $V\_IC1$ ，非常符合“每个控制电压需要响应 10 次”的需求。

在 Simulink 里做以下设置

Solver
- 打开 Model Settings (Ctrl+E) → Solver
- Type 选 Fixed-step
- Solver 选 ode4 (Runge–Kutta)（或者其他固定步长求解器）
- Fixed-step size 填 1e-5
PID Controller
- 双击 PID 块（或 PID Controller），在 Discrete 模式下
- Sample time (Ts) 填 1e-4
Rate Transition（可选，但推荐）
- 在 PID 输出和物理模型输入之间加一个 Rate Transition 块，保证两条不同速率信号安全切换。
- 默 1e-4 → 1e-5，不用改参数，Simulink 会自动用最近邻保持。
Initial Condition
- 把物理模型里初始化状态的那一支“Initial Condition”块的初值写成变量 x0（工作区定义）。
To Workspace
- 用 To Workspace 块把 outZ（物理模型的 Z）和 outV_IC1（PID 输出）导出。

MATLAB 脚本调用示例

%%—— 预先在 base workspace 定义好初始状态和参考值 ——%%
x0 = [...];             % 初始状态向量
Z_reference = 0;        % 参考值恒为 0assignin('base','x0',x0);
assignin('base','Z_reference',Z_reference);%%—— 设置仿真——%%
Ts_ctrl   = 1e-4;       % 控制器采样周期
Ts_solver = 1e-5;       % 求解器步长mdl = 'PID';
open_system(mdl);% 1) 全局求解器
set_param(mdl, ...'Solver',        'ode4', ...         % 固定步长 RK4'FixedStep',     num2str(Ts_solver), ...'StopTime',      num2str(Ts_ctrl)    % 仿真到一个控制步
);% 2) PID 块采样时间
set_param([mdl '/PID Controller'], ...'SampleTime', num2str(Ts_ctrl) ...
);%%—— 运行仿真 ——%%
simOut = sim(mdl, ...'SrcWorkspace','base', ...'SaveOutput','on', ...'ReturnWorkspaceOutputs','on');%%—— 提取最后一步 ——%%
Z_all     = simOut.get('outZ');        % time series
V_IC1_all = simOut.get('outV_IC1');Z     = Z_all(end);
V_IC1 = V_IC1_all(end);fprintf('t=%.0e 时刻： Z = %.6f,  V_IC1 = %.6f\n', Ts_ctrl, Z, V_IC1);

python 实现

对于相同的求解器初始值，如果 MATLAB PID 和 python PID 的 第一步输出相差一个数量级，最常见的坑就是：

在这里插入图片描述

MATLAB PID 块的参数含义可能和你认为的不一样（详见下文 Simulink 中两种 PID 结构中 I 增益的表示）

Simulink 默认的 PID 块有好几种形式（Parallel / I-form / P-only …），比如，连续时域中 它的 “I 增益” 在 Parallel 形式下 表示的是 $\frac{K_i}{s}$ 。如果直接把块里的 “I 增益” 当成 Python 里的 $K_i$ ，就会差好几倍。
检查：双击 Simulink PID 块，看它 是连续时间还是离散时间，是 Parallel 形式还是 I-form，然后 把参数转换成 Python 里的 $K_p,\,K_i,\,K_d,\,N)$ 平行形式。

对号入座，才能在 Python 里准确还原 MATLAB/Simulink 的 PID 积分行为。

离散 PID 实现（并行形式）

Simulink 的离散 PID Controller 块并不是把连续‐时间的

$u(t)=-\bigl(K_p e + K_i\!\int e\,dt + K_d\dot e\bigr)$

简单地用欧拉积分＋一阶滤波来近似。它用的是一整套 “ $z$ 域” 差分方程：

$C(z)\;=\;P \;+\; I\,T_s\frac{1}{z-1} \;+\; D\,\frac{N}{1+N\,T_s}\,\frac{z-1}{z}$

积分项

$u_I[k]=u_I[k-1]+I\,T_s\;e[k]$
微分项（带滤波）
令 $N_d=N\,T_s$ ，有

$y_D[k] =\frac{1}{1+N_d}\,y_D[k-1] \;+\;\frac{K_d\,N_d}{1+N_d}\bigl(e[k]-e[k-1]\bigr)$
总输出

$u[k] = K_p\,e[k] + u_I[k] + y_D[k]$

而如果在 Python 里用的是连续时间那套，

new_integral = integral + e*dt
derivative  = (N*dt*raw_derivative + prev_filtered)/(1+N*dt)
u =  Kp*e + Ki*integral + Kd*derivative

它们并不等价。

下面这段代码给出了一个跟 Simulink 离散 PID Controller 块 一模一样的差分实现。它对应于块上方显示的 并行形式，

$\;+\; I\,T_s\frac{1}{z-1} \;+\; D\frac{N}{1+N\,T_s}\frac{z-1}{z},$

其中 “积分(I)”框里显示的就是 $I\,T_s$ （注意勾选），所以在代码里直接用它来做积分增量。积分输出 和 最终的控制量 都被钳位在 $[-V_{\rm lim},\,V_{\rm lim}]$ 以内。

$P$ 项没有啥好说的，首先根据上面说的实现 $I$ 项，

# 内部状态
self.uI     = 0.0   # 上一步积分项 u_I[k-1]
self.yD     = 0.0   # 上一步滤波微分 y_D[k-1]
self.e_prev = 0.0   # 上一次误差 e[k-1]

# 积分项（差分形式）
uI_candidate  = self.uI_prev + self.Ki * error   # use I*Ts in simulink
# Anti-windup: limit integral term 将 integral 限在 [-V_lim, +V_lim]
uI_candidate  = torch.clamp(uI_candidate, -self.V_lim, self.V_lim)

注意，在 Simulink 离散 PID Controller 块中，积分项的更新总是滞后一拍，它遵循“先用上一步（或初始）的状态算输出，再更新状态” 的时序。

可以在 MATLAB 中只打开 I 控制器实验，会发现第一步的控制器输出为 0。

初始时刻 $t = 0$ 的输出

块会先把“初始积累值” $u_I[0]$ （默认由 Integrator Initial Condition 参数决定，缺省为 0）和“初始滤波值” $y_D[0]$ （同样缺省为 0）带入控制律中。
此时即便误差 $e[0]\neq0$ ，积分项 还没有 加上 $\,K_iT_s\cdot e[0]$ 的那部分增量，因此积分输出项是 0；

在第一个采样周期结束（ $t=T_s$ ）时

Simulink 才会用 $e [0]$ 来更新积分器：
$u_I[1] = u_I[0] + (K_iT_s)\,e[0],$
并在此刻把新的 $u_I[1]$ 参与下一次输出计算。

把 forward() 实现改成这样，就能和 Simulink 保持一致：

# 先用 P项、旧滤波、微分项 计算输出
u_unsat = self.Kp*e + self.uI_prev + self.yD# 再更新积分状态，为下一步准备
uI_candidate  = self.uI_prev + self.Ki * error   # use I*Ts in simulink

接下来是微分项 $D$ 的实现，

在这里插入图片描述

 # 更新 D 项 —— 前向欧拉离散滤波, 差分方程来源于： yD[k] = yD[k-1] + Ts*( Kd*N*(e[k]-e[k-1])/Ts - N*yD[k-1] )
raw_deriv = (error - self.e_prev) / self.dt
self.yD = self.yD_prev + self.dt * ( self.Kd * self.N * raw_deriv - self.N * self.yD_prev )

Simulink 中两种 PID 结构（并联形式, I-形式）下连续/离散时域里积分增益 I 的表示

先把几个概念理清：

连续域的 $s$ （拉普拉斯算子）
- 把时间域信号 $f (t)$ 变换到复频域后，微分 $\tfrac{d}{dt}$ 对应乘以 $s$ ：
  
  $\mathcal{L}\{\,\dot x(t)\} = s X(s).$
- 因此， $\dfrac{K_i}{s}$ 就是“对误差做积分”—— $K_i$ 是积分增益，单位大致是“控制量／（误差×秒）”。
离散域的 $\tfrac{1}{z-1}$ （差分算子）
- 把时间序列 $x [k]$ 做 $z$ 变换后，
  
  $\mathcal{Z}\{\,x[k]-x[k-1]\} = (1 - z^{-1})X(z) \quad\Longrightarrow\quad \mathcal{Z}\Bigl\{\sum_{i=0}^k x[i]\Bigr\} = \frac{1}{1-z^{-1}}X(z) = \frac{z^{-1}}{1-z^{-1}}\,X(z)\,.$
- 换一种常见写法 $\tfrac{1}{z-1} = -\tfrac{z^{-1}}{1-z^{-1}}$ ，它就对应“累加求和”那一块。离散积分项一般写成
  
  $K_i\,T_s\;\frac{1}{z-1} \quad\Longleftrightarrow\quad u_I[k] = u_I[k-1] + K_i\,T_s\,e[k].$
- 这里 $T_s$ 是采样周期，把每步的积分增益拉成了和连续域里 $K_i/s$ 数值等价的 “每步累加” 形式。

并联（Parallel） vs 理想（I-Form），use I*Ts 的勾选

特性	并联（Parallel）形式	I-Form（Ideal）形式
连续时域公式	$C(s)=K_p + \dfrac{K_i}{s} + K_d\,\frac{N}{1+N \frac{1}{s}}$	$C(s)=K_p + \dfrac{K_p}{T_i}\,\dfrac{1}{s} + K_p K_d \frac{N}{1+N \frac{1}{s}}$
离散时域公式	$C(z)=K_p + K_i T_s \dfrac{1}{z-1} + K_d \frac{N}{1 + N T_s\frac{1}{z-1}}$	$C(z)=K_p + \dfrac{K_p T_s}{T_i}\,\dfrac{1}{z-1} + K_p K_d \frac{N}{1 + N T_s\frac{1}{z-1}}$
离散时域差分实现	$u_I[k]=u_I[k-1]+K_i\,T_s\,e[k]$	$\;u_I[k]=u_I[k-1]+\tfrac{K_p\,T_s}{T_i}\,e[k]$

在 离散模式 下，
- 勾选 “Use I*Ts” 时，PID Controller 块中的“积分 (I)”文本框显示的就是 $K_i\,T_s$ ，这个数值直接用于差分方程 $u_I[k]=u_I[k-1]+(K_iT_s)\,e[k]$ 。
- 若 不勾选 “Use I*Ts”，文本框则显示纯粹的离散积分增益 $K_i$ ，块会在后台自动再乘以 $T_s$ 来进行积分运算，即等效地使用 $K_iT_s$ 。
在 连续模式 下，无论是否有“Use I*Ts”，文本框显示的都是连续域的积分增益
$\dfrac{K_i}{s}$ 中的 $K_i$ ，即单位为“控制量/（误差·秒）”的参数。

结构类型	文本框显示（离散）	差分方程增量	Simulink 帮你做的事
并联（Parallel）	$K_iT_s$ （勾选 Use I*Ts） $K_i$ （不勾选）	$u_I[k]=u_I[k-1]+(K_iT_s)\,e[k]$	若不勾选，块会自动乘以 $T_s$ ；勾选则用文本框数值直接累加。
理想（I-Form）	$T_i$ （始终）	$\;u_I[k]=u_I[k-1]+\tfrac{K_p\,T_s}{T_i}\,e[k]$	勾选后内部把 $T_i$ 换算为 $\tfrac{K_pT_s}{T_i}$ 。

根据 MATLAB 官方文档，

在这里插入图片描述

【Tip】

如果你已经在并联下直接获得了离散 $K_iT_s$ （比如 Simulink 给你的那一栏），就直接拿来用，不要再二次换算；
在并联（Parallel）结构下，Simulink 参数栏 “积分(I)”框里显示的正是 $K_i\,T_s$ （离散）或 $K_i$ （连续），直接把 Simulink 中读到的 “积分(I)” 值当做 $K_i\,T_s$ （离散）或 $K_i$ （连续）使用即可。
如果是在 理想（I-Form） 下得到一个积分时间 $T_i$ ，“积分( I )” 填的是 $T_i$ ，它会自动在内部换算成 $K_pT_s/T_i$ ；python 实现需要再用
$K_{i,\mathrm{discrete}} = \frac{K_p\,T_s}{T_i}$
换算成“每步累加增量”。

什么时候用连续域 vs 离散域？

如果把整个闭环模型当成一个离散‐时间系统 —— 比如用

$x_{k+1} = e^{A\,T_s}\,x_k \;+\; A^{-1}(e^{A\,T_s}-I)\,B\,u_k$

这样的精确离散更新公式（基于矩阵指数）来推进状态，那么整个控制器的设计、仿真与实现，就更自然地放在离散时域里：
- 在每个采样点 $k$ 计算误差 $e [k]$ ，更新离散 PID 差分方程，输出 $u_k$ ，然后用上面的公式算 $x_{k+1}$ 。
- 这时把 Simulink 里“并联形式”的离散 PID 差分方程直接拿过来用，就保证“数值上一致”。
- 连续域的 $K_i/s$ 并不直接作用于这个仿真步长里，它只是一个设计思想——真正跑的时候都要离散化。
如果你打算在连续‐时间下设计 PID（比如先用频域调参工具、根轨迹、Bode 曲线，得到 $K_i,K_p,K_d$ ），再把它离散化（如向后差分、双线性变换 Tustin 法、零阶保持 ZOH），那么就先在连续时域下选好参数，然后再用某种离散化方法转换成离散差分形式。