手把手教你学Simulink——基于扩展卡尔曼滤波（EKF）的列车牵引电机转子温度在线估计仿真

目录手把手教你学Simulink——基于扩展卡尔曼滤波EKF的列车牵引电机转子温度在线估计仿真一、为什么要用 EKF 估计电机转子温度二、电机热‑电增广模型EKF 用2.1 定子电阻温度关系2.2 一阶热近似简化2.3 PMSM dq 方程连续2.4 EKF 状态向量三、关键参数四、Simulink 建模手把手4.1 Step 1️⃣ —— PMSM FOC 电流环真值4.2 Step 2️⃣ —— EKF离散推荐 MATLAB Function 或 Extended Kalman Filter块 R2019b■ 若用 Control System Toolbox → extendedKalmanFilterobject推荐■ 若旧版无 EKF 块手写 5 方程预测 更新用 MATLAB Function可提供但较冗4.3 Step 3️⃣ —— 热真值校验4.4 Step 4️⃣ —— 运行 Scenario五、结果解读典型✅ 稳态 空载 (0~2s)✅ 加载 (2s~)✅ 电阻修正六、工程注意点**七、结论**手把手教你学Simulink——基于扩展卡尔曼滤波EKF的列车牵引电机转子温度在线估计仿真一、为什么要用EKF 估计电机转子温度列车牵引电机异步 / PMSM运行中不能直接测转子 / 绕组温度无滑环、成本高、安装点少温升影响永磁体磁链ψ_f(T) 下降PMSM→ 转矩偏差铜损 Rs​(T)Rs0​[1α(T−T0​)]→ 影响电流闭环、效率EKF 优势利用易测信号id​,iq​,ud​,uq​,ωm​​ 在线估计不可测状态Tr​,ψf​或 Rs​温变融合电机模型 过程 / 量测噪声 → 给出估计值 协方差可信度典型状态增广x[id​iq​ψf​Tr​]Torx[id​iq​Rs​]T本例用PMSM 增广 Rs​(T)状态一阶热模型耦合目标借鉴前几篇牵引 PMSMPMSMp3,ψf0​0.175Wb25℃,Rs25​0.05Ω,αCu​0.00393/K热模型T˙r​(Pcu​Pfe​−hA(Tr​−Ta​))/Cth​, 简化 T˙r​−(Tr​−Ta​)/τth​Kth​⋅i2Rs​EKF 观测量id​,iq​量测输入 u[ud​,uq​,ωm​]工况0~2s 加速 nref​1500rpm2s 加负载 TL​80NmTa​25℃观测Tr_est​,Rs_est​,id​,iq​,Tr_true​真值仅校验对比直接用 Rs25​与 EKF 修正 Rs​对 iq​跟踪影响基于 SimulinkPMSM FOC EKF 离散块 / MATLAB Function 热子模型架构是破局关键。无论你是轨道交通牵引传动 / 电机热监控工程师这篇硬核指南都成为你手中“无传感器转子温度 EKF 模板”。二、电机热‑电增广模型EKF 用2.1 定子电阻温度关系Rs​(Tr​)Rs0​⋅[1αCu​⋅(Tr​−Tref​)],Tref​25℃2.2 一阶热近似简化T˙r​−τth​1​(Tr​−Ta​)Cth​Irms2​⋅Rs​(Tr​)​τth​600s热时间常数 10minCth​Rth​⋅τth​,Rth​1K/W⇒ Cth​600J/KIrms​≈(id2​iq2​)/2​近似2.3 PMSM dq 方程连续dtdid​​dtdiq​​ψf​​Ld​1​(ud​−Rs​(Tr​)id​ωe​Lq​iq​)Lq​1​(uq​−Rs​(Tr​)iq​−ωe​Ld​id​−ωe​ψf​) 视为 const (或可增广慢变)​2.4 EKF 状态向量x​id​iq​Tr​​​,u​ud​uq​ωe​​​,y[id​iq​​]Process noise Q小Meas noise R电流传感 noise三、关键参数参数值极对数 p3ψf0​25℃0.175 WbRs25​0.05 ΩαCu​0.00393 /KLd​Lq​1.2 mHJ0.08 kg·m²牵引惯量小示B0.002 N·m·s/radTa​环境温度25 ℃τth​600 sCth​600 J/KRth​1 K/WEKF Ts1e‑4 s同控制采样Imax​±250 Anref​1500 rpmTL​0 → 80 Nm 2s仿真 Ts5e‑6电控制/ EKF Ts1e‑4四、Simulink 建模手把手4.1 Step 1️⃣ —— PMSM FOC 电流环真值PMSM 参数如上Rs​Rs25​初值热模型独立给 Tr​真值三相全桥 400VClarke/Park → 电流 PIKp​0.5,Ki​200→ SVPWM转速 PI → iqref​sat ±Imaxidref​0机械负载Step 0→80Nm 2s量测供 EKFid_meas​,iq_meas​ small white noiseBand‑Limited White Noiseσ0.05Aud​,uq​电压指令 * V_dc/√3 或 反 Park 输出 before PWM 映射ωe​ωm​⋅pfrom machineTrue Tr​校验用一阶热积分用真实 Rs​(Tr​)→ Tr​(t)4.2 Step 2️⃣ —— EKF离散推荐 MATLAB Function 或Extended Kalman Filter块 R2019b■ 若用Control System Toolbox →extendedKalmanFilterobject推荐在 ModelInitFcn%% EKF for PMSM Rotor Temp Ts_ekf 1e-4; % State x [id; iq; Tr] x0 [0; 0; 25]; % start at ambient % Process measurement functions fcnStateTransition (x,u,~) pmms_temp_fcn(x,u); fcnMeas (x,~) [x(1); x(2)]; % we measure id,iq ekf_obj extendedKalmanFilter(... StateTransitionFcn, fcnStateTransition,... MeasurementFcn, fcnMeas,... State, x0,... StateCovariance, diag([1e-3 1e-3 4]),... ProcessNoise, diag([1e-6 1e-6 1e-4]),... MeasurementNoise, diag([0.05^2 0.05^2])); assignin(base,ekf_obj,ekf_obj); assignin(base,Ts_ekf,Ts_ekf);pmms_temp_fcn.m离散预测function x_next pmms_temp_fcn(x,u) % x [id; iq; Tr] % u [ud; uq; we (rad/s)] persistent Ld Lq Rs0 alpha Rth Cth Ta if isempty(Ld) Ld 1.2e-3; Lq 1.2e-3; Rs0 0.05; alpha 0.00393; Rth 1; Cth 600; Ta 25; end id x(1); iq x(2); Tr x(3); Rs Rs0 * (1 alpha*(Tr - Ta)); ud u(1); uq u(2); we u(3); % current derivatives didt ( ud - Rs*id we*Lq*iq ) / Ld; diqdt ( uq - Rs*iq - we*Ld*id - we*0.175 ) / Lq; % psi_f0.175Wb % thermal (1st order) I2 (id^2 iq^2)/2; % rough rms^2 dTrdt -(Tr - Ta)/600 I2*Rs / Cth; % Euler Ts 1e-4; x_next x Ts * [didt; diqdt; dTrdt]; endSimulink 中用Extended Kalman Filter​ 块ekf_obj名Input 1MeasuredOutputs [i_d_meas; i_q_meas]Input 2Inputs [u_d; u_q; w_e]OutputEstimated States [i_d_est; i_q_est; T_r_est]OptionalState Covariance■ 若旧版无 EKF 块手写 5 方程预测 更新用 MATLAB Function可提供但较冗4.3 Step 3️⃣ —— 热真值校验同 Rs​(Tr​) 一阶热模型IntegratorT˙r_true​−(Tr_true​−Ta​)/τth​(id2​iq2​)/2⋅Rs​(Tr_true​)/Cth​初值 25℃Scope 对比 Tr_est​4.4 Step 4️⃣ —— 运行 Scenario时间事件0~2snref​1500rpm, TL​02sTL​80Nm→ iq​↑→ 铜损 ↑ → Tr​↑Scopenrpm​,iq​,Tr_true​,Tr_est​,Rs_est​Rs0​(1α(Tr_est​−25))✅ 期望Tr_est​跟 Tr_true​稳态 err 1~2℃取决于 Q/R 调2s 负载 ↑ → Tr​缓升τ600s ⇒ 几分内到稳EKF 给平滑 Tr_est​无跳✔若直接用 Rs25​算 i_q 参考 → 微偏差高温 Ψ_f 也可同估五、结果解读典型✅ 稳态 空载 (0~2s)iq​≈0→ Tr​微变 (~25℃)EKF Tr_est​≈25℃✔✅ 加载 (2s~)iq​≈TL​/(1.5pψf​)80/(1.5⋅3⋅0.175)≈101.6A铜损 Pcu​Irms2​Rs​≈((02101.62)/2)⋅0.05≈258WTr​按 τ600s 缓升例 10min → 4~5KEKF Tr_est​跟踪 true Tr​滞后仅 EKF din ≈ model mismatch✔✅ 电阻修正Rs_est​0.05⋅[10.00393⋅(Tr_est​−25)]温升 10K → Rs​≈0.05197Ω→ 若电流 PI 用 Rs_est​替代固定 0.05 ⇒ 微小 but 系统一致 ✔六、工程注意点**实因素Simulink 处理Ψ_f(T) 老化可增状态 ψf​慢变EKF 增维铁损 Pfe​加 Pfe​Kh​fBm2​Ke​f2Bm2​→ 热 inflow初值 Tr0​未知用 last park / default 25℃ 大 P0(3,3)10传感器 biasR 矩阵调EKF 自校正慢漂HIL / 实车EKF 定点化 (fixed‑point), 周期 100µs~1ms七、结论**✅ 你掌握了基于 EKF 的牵引电机转子 / 定子温度在线估计 Simulink 原型PMSM FOC 电流环p3,ψf​0.175Wb,Rs​250.05Ω,L1.2mH增广状态 x[id​,iq​,Tr​]热模型 T˙r​−(Tr​−Ta​)/τIrms2​Rs​(Tr​)/Cth​extendedKalmanFilter块或 obj→ 输入 u[ud​,uq​,ωe​], meas [id​,iq​]noisy​输出 Tr_est​,Rs_est​(Tr_est​)加载 → Tr​↑, EKF 跟 true Tr​±1~2K✔ 此法是列车 / EV 牵引电机热监控、过载保护、磁链补偿基础IEC 61800‑9 能效 热模型​可直接扩展Ψ_f(T) 同估4‑state EKF 或 dual EKF异步电机转子电阻 Rr​(T)估Reduced‑order Thermal EKF 融合低算耗HIL 定点 C 生成CMSIS‑DSP / IQmath如果下一步你想要手写 EKF 预测更新5公式MATLAB Function无 toolbox​PMSM Ψ_f(T) 慢变增广状态联合估​异步电机转子 Rr​(T)EKF 估计感应牵引电机​EKF 协方差收敛性 / 噪声 Q,R 调参指南