【MATLAB例程12】永磁同步电机基于增广无迹卡尔曼滤波的未知负载转矩在线观测

原创代码包运行成功文章目录程序简介算法原理与实现仿真场景与结果分析程序特点与适用方向运行结果MATLAB源代码程序简介本程序基于MATLAB搭建永磁同步电机Permanent Magnet Synchronous MotorPMSM仿真模型采用增广无迹卡尔曼滤波Unscented Kalman FilterUKF实现电机电流、机械转速、转子电角度及未知负载转矩的联合估计。在永磁同步电机运行过程中外部负载转矩变化会影响电机的动态响应和转速稳定性。如果负载转矩无法直接测量就需要借助状态观测算法根据电机数学模型和已有测量信息对其进行估计。本程序将未知负载转矩作为增广状态引入UKF利用电压和带噪电流测量信息进行递推估计并通过负载阶跃变化工况研究UKF对未知负载的跟踪能力。算法原理与实现程序采用五维增广状态模型包含d轴电流、q轴电流、机械角速度、转子电角度和未知负载转矩。与传统扩展卡尔曼滤波EKF不同UKF通过Sigma点传播状态均值和协方差无需对非线性模型求取雅可比矩阵。程序主要包括以下功能增广状态估计将未知负载转矩纳入状态向量采用随机游走模型描述其变化并通过电机动力学方程建立负载与其他状态之间的联系。Sigma点无迹变换通过生成和传播Sigma点完成非线性状态预测及电流量测更新。电角度环绕处理针对转子电角度的周期性采用角度环绕方法计算均值和估计误差。数值稳定性检查对无迹变换、协方差正定性及状态估计结果进行检查。仿真场景与结果分析程序构建磁场定向控制FOC下的永磁同步电机仿真环境采用恒定转速给定并模拟外部负载转矩发生两次阶跃变化的工况。仿真总时长为1.2秒采样周期为0.1 ms电流测量噪声标准差为0.04 A。负载转矩设置如下00.45秒负载转矩为0.35 N·m。0.450.85秒负载转矩增加至0.90 N·m。0.851.20秒负载转矩降低至0.45 N·m。程序主要输出以下结果未知负载转矩观测对比负载转矩真实值与UKF估计值分析负载阶跃变化时的跟踪情况。负载转矩估计误差展示负载估计误差随时间的变化并计算指定稳定区间内的均方根误差RMSE。转速与电角度估计对比机械转速、转子电角度的真实值与估计结果分析状态观测精度。电流滤波结果展示α-β轴带噪电流测量、真实电流及滤波估计电流。滤波性能评价输出转速、角度及负载转矩的RMSE以及归一化新息平方NIS和协方差数值检查结果。程序特点与适用方向本程序采用基础MATLAB函数实现包含完整的PMSM数学模型、FOC控制仿真、UKF状态预测、非线性量测更新、未知负载估计及结果可视化过程。电机参数、负载阶跃时刻、负载变化幅度、电流测量噪声、采样周期及UKF参数等均可自行调整便于开展不同负载工况下的状态观测研究。本例适用于永磁同步电机未知负载转矩观测、增广状态估计、无迹卡尔曼滤波及非线性状态估计算法的学习与仿真研究也可以作为后续研究负载扰动补偿和电机抗扰控制的基础例程。运行结果状态估计误差电压与新息诊断电流观测与估计负载转矩观测转速与角度观测命令行输出MATLAB源代码部分代码% 未知负载转矩观测独立MATLAB例程% 2026-10-16/Ver1作者matlabfilterV同号可获取代码定制、讲解等% 只使用基础MATLAB。参数为普通标量信号为数值数组不使用结构体。% 点击运行即可末尾含全部局部函数。Figure无水印配套PNG由维护工具导出。clear;clc;close all;rng(0);%% 共用教学参数与单位不将本例解释为Word中15 kW电机的实机验证Rs0.45;Ld3e-3;Lq3.6e-3;psi_f0.075;pole_pairs4;J2e-3;B5e-4;Ts1e-4;simulation_time1.2;dc_voltage120;voltage_limitdc_voltage/sqrt(3);current_limit10;Kp_speed0.25;Ki_speed4.0;current_bandwidth1200;Kp_dLd*current_bandwidth;Kp_qLq*current_bandwidth;Ki_dRs*current_bandwidth;Ki_qRs*current_bandwidth;time0:Ts:simulation_time;sample_countnumel(time);Kt1.5*pole_pairs*psi_f;Kr1.5*pole_pairs*(Ld-Lq);% 维护用汇总也是数值向量各指标仍以独立、有意义的变量保留。speed_rmseNaN;angle_rmse_degNaN;load_rmseNaN;parameter_relative_rmse_pctNaN;secondary_parameter_rmse_pctNaN;improvement_pctNaN;current_rmseNaN;peak_speed_errorNaN;gap_speed_rmseNaN;recovery_speed_rmseNaN;final_parameterNaN;final_secondary_parameterNaN;disturbance_speed_peakNaN;control_current_peakNaN;takeover_speed_peakNaN;recovery_timeNaN;covariance_minInf;covariance_chol_failures0;ut_error0;state_is_finite1;%% 编码器FOC台架仅负责生成实际电压、电流和物理一致的真值speed_reference8020*(1tanh((time-0.25)/0.035))/2...-10*(1tanh((time-0.85)/0.035))/2;load_true0.35*ones(1,sample_count);speed_reference(:)80;load_true(time0.45)0.90;load_true(time0.85)0.45;state_truezeros(4,sample_count);voltage_abzeros(2,sample_count);plant[0;(0.35B*80)/Kt;80;0.30];state_true(:,1)plant;controller_integralszeros(3,1);iq_referencezeros(1,sample_count);voltage_saturatedfalse(1,sample_count);iq_bias(0.35B*80)/Kt;fork2:sample_count feedback_abcurrent_measurement(plant);[uab,controller_integrals,iq_reference(k-1),voltage_saturated(k-1)]foc_controller(...feedback_ab,plant(4),plant(3),speed_reference(k-1),controller_integrals,Ts,...Rs,Ld,Lq,psi_f,pole_pairs,voltage_limit,current_limit,...Kp_speed,Ki_speed,Kp_d,Ki_d,Kp_q,Ki_q,iq_bias);voltage_ab(:,k-1)uab;plantmotor_truth_step(plant,uab,Ts,Rs,Ld,Lq,psi_f,pole_pairs,J,B,load_true(k-1));state_true(:,k)plant;endvoltage_ab(:,end)voltage_ab(:,end-1);iq_reference(end)iq_reference(end-1);current_ab_true[cos(state_true(4,:)).*state_true(1,:)-sin(state_true(4,:)).*state_true(2,:);...sin(state_true(4,:)).*state_true(1,:)cos(state_true(4,:)).*state_true(2,:)];noise_std0.04*ones(2,sample_count);current_ab_meascurrent_ab_truenoise_std.*randn(2,sample_count);measurement_validtrue(1,sample_count);current_ab_meas(:,~measurement_valid)NaN;%% 状态观测不使用台架的真实转速、真实角度或未知负载state_dimension5;ut_errorcheck_sigma_transform_ukf();assert(abs(wrap_angle((-pi0.01)-(pi-0.01))-0.02)1e-12);filter_count1;state_historyzeros(state_dimension,sample_count,filter_count);nis_historynan(filter_count,sample_count);R_historyzeros(2,sample_count,filter_count);weight_historyones(2,sample_count);innovation_historynan(2,sample_count,filter_count);sigma_speed_historyzeros(filter_count,sample_count);forfilter_index1:filter_count x[0.15;(0.35B*80)/Kt-0.1;85;0.42];Pdiag([0.3^2,0.3^2,12^2,0.18^2]);Qdiag([2e-6,2e-6,2e-3,2e-8]);ifstate_dimension5x[x;0.20];Pblkdiag(P,0.25^2);Qblkdiag(Q,3e-5);endR0.04^2*eye(2);innovation_covarianceR;state_history(:,1,filter_index)x;R_history(:,1,filter_index)diag(R);sigma_speed_history(filter_index,1)sqrt(P(3,3));fork2:sample_count[xp,Pp]motor_predict_ukf(x,P,Q,voltage_ab(:,k-1),Ts,Rs,Ld,Lq,psi_f,pole_pairs,J,B,0.35);ifmeasurement_valid(k)[predicted,Czz,Pxz]current_moments_ukf(xp,Pp);residualcurrent_ab_meas(:,k)-predicted;effective_RR;[x,P,residual,nis_value]state_update_ukf(xp,Pp,current_ab_meas(:,k),effective_R);innovation_history(:,k,filter_index)residual;nis_history(filter_index,k)nis_value;elsexxp;x(4)wrap_angle(x(4));P(PpPp)/2;end[~,chol_flag]chol(P);covariance_chol_failurescovariance_chol_failures(chol_flag~0);covariance_minmin(covariance_min,min(eig(P)));assert(chol_flag0all(isfinite(x)),状态或协方差异常);state_history(:,k,filter_index)x;R_history(:,k,filter_index)diag(R);sigma_speed_history(filter_index,k)sqrt(P(3,3));endendstate_eststate_history(:,:,end);current_est[cos(state_est(4,:)).*state_est(1,:)-sin(state_est(4,:)).*state_est(2,:);...sin(state_est(4,:)).*state_est(1,:)cos(state_est(4,:)).*state_est(2,:)];innovationinnovation_history(:,:,end);nis_valuesnis_history(end,:);speed_truestate_true(3,:);speed_eststate_est(3,:);angle_truestate_true(4,:);angle_eststate_est(4,:);speed_errorspeed_est-speed_true;angle_errorwrap_angle(angle_est-angle_true);evaluationtime0.24;speed_rmsesqrt(mean(speed_error(evaluation).^2));angle_rmse_degsqrt(mean(angle_error(evaluation).^2))*180/pi;peak_speed_errormax(abs(speed_error(evaluation)));valid_evaluationevaluationmeasurement_valid;current_rmsesqrt(mean((current_est(:,valid_evaluation)-current_ab_true(:,valid_evaluation)).^2,all));control_current_peakmax(abs(state_true(2,:)));load_eststate_est(5,:);plateau(time0.65time0.80)|time1.05;load_rmsesqrt(mean((load_est(plateau)-load_true(plateau)).^2));%% 绘图输出的Figure保持无水印完整代码和函数如需帮助或有滤波相关的MATLAB代码定制需求可从个人主页左侧联系我