尧图网站建设 
简介电力系统次同步谐振的数字仿真PDF资料面向电力系统运行分析、电气工程相关专业的读者聚焦串联补偿线路引发的次同步谐振机理及数字仿真方法。内容包含同步发电机电气与轴系数学模型、隐式梯形积分法下的差分模型推导并基于MATLAB编程再现次同步谐振过程系统分析网络参数、串联补偿度对轴系扭振响应的影响为抑制次同步谐振提供理论依据与仿真参考。压缩包内提供单份PDF格式的学术论文大小仅1.1MB适合深入研读与对照推导。已有88人学习下载。读者可获得完整的机理分析、模型构建思路与仿真试验结论有助于在电力系统稳定分析、串补输电工程中快速建立相关认知并指导参数设计可作为相关课题研究的基础参考资料。1. 次同步谐振的数字仿真先弄清楚它在“谐振”什么做电力系统数字仿真的工程师迟早会遇到次同步谐振这个词。常见于装有串联补偿电容的输电线路上或风电变流器与弱网交互时标题里“数字仿真”通常不是指教科书例题而是要在时域里把发电机轴系的扭振过程复现出来。很多人按机电暂态模型搭出来发现频率对不上或波形看不到振荡问题大多出在机械侧自由度不足、电气网络等效过粗、时间步长取太大。这篇文章以电力系统次同步谐振数字仿真为主线从六质量块轴系建模到串联补偿度扫参再到结果可信度验证给出可直接改参数运行的最小化方案。2. 仿真之前先把“次同步谐振”拧成可计算的模型2.1 为什么是“多质量块轴系”而不是单刚体发电机模型次同步谐振的本质是机械侧和电气侧在低于工频的频率上交换能量。机械侧的能量载体不是发电机整机而是旋转轴系各段之间的扭转变形。如果仿真里只保留一个刚体转子就只能看到整个轴系整体摇摆的一个模态看不到高压缸、中压缸、低压缸和发电机之间相对扭转形成的扭振模态也就无从判断“哪个频率被激起来了”。工程上处理这个问题普遍把轴系离散成若干质量块相邻质量块之间用扭转弹簧连接每个质量块带一个惯性常数和一个阻尼系数。发电机模型可以用两质量块先跑通现象但要用于次同步谐振风险分析至少要做到四到六质量块。通常的划分方式是高压缸、中压缸、两个低压缸、发电机和励磁机一共六个质量块。因为励磁机质量小、刚度大很多资料里也会把励磁机去掉但对于数字仿真多保留一两个质量块的代价只是矩阵多几行没必要牺牲这一自由度。为什么机械阻尼在仿真里特别关键因为轴系本身的机械阻尼很小模态阻尼的主要来源反而在电气侧。当电气网络在某个次同步频率上呈现等效负阻尼时扭振模态就会从电气侧吸收能量幅值越振越大。数字仿真若只给一个常数正阻尼等于提前把最危险的失稳条件抹掉了。2.2 电气侧频率怎么算串补度把谐振点抬到哪电气侧的次同步谐振频率由串联补偿电容容抗和线路等值感抗共同决定。以最简单的单机无穷大系统为例补偿度为 k工频为 f0则电气谐振频率约为f_er f0 * sqrt(k)如果系统频率是 60Hz串补度 10% 时f_er 大约 19Hz串补度 20% 时大约 26.8Hz。这个区间和大型汽轮发电机组轴系的一阶到三阶扭振模态高度重叠所以“串补度越高越接近工频”并不意味着越安全而是要逐点看它落在哪个扭振模态附近。不过这只是粗略估算。真正决定是否失稳的是这个频率下电气回路的相位关系也就是电气转矩相对转子角变化的响应。严格做法是用小信号扫频得到次同步频段的等效阻抗再折算成复转矩系数数字仿真的价值恰恰就是在时域里把这个相位滞后出现的后果完整地演出来。2.3 给轴系定参数IEEE FBM 与模态估算代码做次同步谐振仿真最常用的是电力系统工程师熟知的 First Benchmark Model即 IEEE 次同步谐振标准模型。它有公开的机械参数适合用来验证算法和对比结果。下面这组参数是六质量块配置单位采用标幺值坐标统一到电弧度质块惯性常数 Hs相邻轴段扭转刚度 Kpu/radHP0.0929HP-IP19.303IP0.1556IP-LP134.929LP10.8582LP1-LP252.038LP20.8842LP2-GEN70.858GEN0.8685GEN-EXC2.822EXC0.0342拿到这组参数先别急着进时域仿真我一般会先用广义特征值估算一次轴系固有扭振频率确认模型自由度没有丢。下面的 MATLAB 代码可以完成这一步% freq_ssr.m —— 由机械参数估算轴系固有扭振频率 f0 60; % 系统额定频率 omega0 2*pi*f0; H [0.0929; 0.1556; 0.8582; 0.8842; 0.8685; 0.0342]; K [19.303; 34.929; 52.038; 70.858; 2.822]; M diag(2*H); Kmat zeros(6); for i 1:5 Kmat(i,i) Kmat(i,i) K(i); Kmat(i1,i1) Kmat(i1,i1) K(i); Kmat(i,i1) -K(i); Kmat(i1,i) -K(i); end [~, D] eig(omega0 * Kmat, M); % 广义特征值omega0 统一角度与角速度量纲 ev real(diag(D)); fn sort(sqrt(-ev(ev -1e-8)))/(2*pi); disp(轴系固有扭振频率(Hz)); disp(fn);代码里的 M 是惯性常数对角阵2H 是运动方程中需要的时间常数Kmat 是标准的刚性对称刚度矩阵每一段的刚度只有相邻两个自由度有耦合。eig 返回的特征值包含一个零特征值对应整根轴同步旋转的刚体模态对扭振没有意义所以过滤掉。负特征值开方后得到的就是各阶扭振频率单位是 Hz。这套计算里最容易出错的是量纲。如果 H 用秒、K 用标幺值、角度用电弧度那么必须乘上工频角速度 omega0否则算出来的频率会整体偏小或偏大。工程上拿到机组实测扭振频率后通常会用这个脚本先对一遍模型参数。3. 用 MATLAB/Simulink 搭一个 SSR 最小可复现仿真3.1 时域仿真和特征值分析谁先谁后次同步谐振数字仿真有两条并行路径一是线性化状态空间模型求特征值看各模态阻尼二是时域积分看扰动后扭振波形。特征值法快能直接判断某个模态是否失稳时域法直观能把非线性元件、控制保护和波形叠加都装进去。我的顺序通常是先做特征值分析把危险工况筛出来再用时域仿真去验证危险工况的波形和阻尼。这样既避免在安全工况上浪费积分时间也避免了只看特征值而漏掉饱和、死区或保护动作导致的非线性行为。实际操作中纯 Simulink 框图搭六质量块轴系并不方便我建议把机械部分写成 MATLAB 状态空间函数再在 Simulink 里用 State-Space 模块对接电气模型。3.2 六质量块状态空间搭建ode45 能跑的最小脚本下面这个脚本是“电力系统次同步谐振数字仿真”里最精简的可运行版本。它把六质量块轴系做成线性状态空间在发电机质块上施加短时制动转矩用来模拟故障切除后电气转矩对轴系的冲击。状态变量取为各质量块的转子电角度偏差和标幺转差外部输入是发电机质块上的电磁转矩增量。% ssr_torque_sim.m —— 六质量块轴系次同步扭振时域仿真 function ssr_torque_sim() f0 60; omega0 2*pi*f0; H [0.0929; 0.1556; 0.8582; 0.8842; 0.8685; 0.0342]; K [19.303; 34.929; 52.038; 70.858; 2.822]; M diag(2*H); Kmat build_Kmat(K); Dm diag([0.002; 0.003; 0.005; 0.005; 0.004; 0.001]); % 机械阻尼需按实测细化 Ks 1.2; % 电气同步转矩系数 De -0.6; % 次同步等效阻尼负值代表负阻尼 A [zeros(6), omega0*eye(6); -M\Kmat, -M\Dm]; % 电气转矩增量落在发电机质块索引5ΔTe Ks*Δδ De*Δω A(10,5) A(10,5) - Ks / M(5,5); A(10,10) A(10,10) - De / M(5,5); tStart 0.02; tEnd 0.08; amp -0.03; % 0.04s冲击转矩 Text (t) event_torque(6, 5, amp, tStart, tEnd, t); f_ode (t,x) A*x [zeros(6,1); M\Text(t)]; [t, x] ode45(f_ode, [0 3.0], zeros(12,1), odeset(MaxStep,1e-4)); delta x(:,1:6); T_shaft Kmat * delta; % 各轴段扭矩增量 plot(t, T_shaft(2,:)*100, LineWidth, 1); xlabel(时间 / s); ylabel(轴扭矩增量 / %); title(IP-LP1 轴段扭振响应); end function Kmat build_Kmat(K) Kmat zeros(6); for i 1:5 Kmat(i,i) Kmat(i,i) K(i); Kmat(i1,i1) Kmat(i1,i1) K(i); Kmat(i,i1) -K(i); Kmat(i1,i) -K(i); end end function Tq event_torque(n, idx, amp, t0, t1, t) Tq zeros(n,1); if t t0 t t1 Tq(idx) amp; end end这段代码最需要关注的是 A 矩阵最后两行的修正。发电机质块上的电磁转矩增量同时包含角度项和转速项角度项相当于给轴系增加了一个电气刚度转速项则是电气阻尼。把 De 设成负值后发电机质块会在转速偏离同步点时获得一个同向的转矩增量这就是次同步振荡能够自激的条件。改动 De 的符号和大小可以快速观察同一套机械参数在“正阻尼”和“负阻尼”两种场景下的区别。3.3 电气负阻尼是怎么进来的一个可调的 De常数 De 看起来太简单但它抓住了次同步谐振失稳的核心在特定次同步频率附近电气网络呈现的阻尼可能是负的。产生负阻尼的物理原因是串联补偿电容把系统谐振频率推到工频以下转子扭振产生的次同步电流经过网络阻抗后形成的电磁转矩和转子速度分量之间的相位关系不再构成制动转矩。严格的 De 不是常数而是频率的函数。工程上更常见的方法是扫频求复转矩把不同频率下的电气阻尼曲线拿到再和轴系机械模态叠加判断稳定性。上面脚本里的 De 相当于在某个危险频率上取了一个代表值适合做机理验证如果要把仿真结果用于工程结论必须在目标频率点上做网络阻抗扫频来更新 De。3.4 求解器与步长的配置表用 ode45 跑线性轴系模型很顺但一旦串上限流、饱和或控制坏就要换成刚性问题求解器。下面这组配置是我常用的起点参数推荐值说明仿真总时长3~5 s扭振模态衰减周期长短于 2s 看不到趋势MaxStep1e-4 s最高关心 50Hz每周期至少 200 个点RelTol1e-6默认容差容易丢掉弱阻尼模态AbsTol1e-8角度和转差量级差异大需要分开设置扰动时长0.04~0.08 s模仿故障切除后的电气冲击这里 MaxStep 不是随便压到 1e-4。如果只关心 20Hz 的扭振模态1ms 步长已经够要同时捕捉 60Hz 工频和 50Hz 扭转耦合步长必须低于 0.2ms。判断标准是最高关心频率的每个周期至少有 20 个采样点否则 FFT 结果会失真阻尼比拟合也不可靠。4. 扰动与串补度扫参把危险模式逼出来4.1 扰动怎么加用短时制动转矩而不是三相短路有些人会把故障事件直接设成三相短路但这在小信号轴系仿真里不是好选择。三相短路瞬间电气转矩几乎消失同时伴随复杂的电磁暂态对轴系模型的激励成分太杂很难定位到底激发了哪个扭振模态。我更常用的是短时制动转矩在发电机质块上施加一个幅值固定、时间很短的负向转矩相当于故障切除瞬间电气转矩突变。这样所有扭振模态都会被冲击起来随后按各自的阻尼衰减或增长。步骤可以固定成三条把扰动转矩幅值设为基础转矩的 1%5%避免进入非线性区。扰动持续时间取 4080ms覆盖次同步频段最高频率的多个周期。仿真启动后先让系统稳定 20ms 再加入扰动方便后续截取基线。如果必须模拟短路后串补重新投入需要把模型换到 Simulink/SimPowerSystems 里加真实线路和电容开关机械侧再用 S-Function 或状态空间模块对接不能继续用上面的线性 A 矩阵。4.2 把串补度从 5% 扫到 75%看电气谐振频率落到哪个模态上次同步谐振数字仿真的核心实操是扫串联补偿度。每改变一次补偿度电气谐振频率就移动一次可能跨过某个扭振模态。最简单的一轮扫参可以这样跑% scan_subsync.m —— 串补度与电气谐振频率对应关系 f0 60; k_range 0.05:0.05:0.75; for i 1:length(k_range) f_er f0 * sqrt(k_range(i)); fprintf(串补度 %.2f - 电气谐振频率 %.2f Hz\n, k_range(i), f_er); end这个循环只能定位“可能出现问题”的频率范围真正的风险判断要回到上一章的特征值脚本对每个串补度重算 De 和模态阻尼。常见做法是把这段扫频结果和轴系固有频率表放在一起做成一张风险对照表串补度电气谐振频率Hz潜在风险0.1523.2可能落入中频扭振模态0.3032.9与较高频扭振模态接近0.4540.2进入高段扭振区0.6046.5接近同步频率伴随感应发电机效应0.7552.0高于工频风险转向另一类问题工程上要注意表格里的 f_er 只是不计电阻的近似值。线路电阻和发电机次暂态电抗都会把实际谐振频率拉低几个 Hz所以扫参时建议把补偿度步长设到 1%2%否则可能跳过阻尼穿越点。4.3 扫参时的参数表步长、时长、扫频步长扫参不是把补偿度设置成数组然后循环那么简单。串补度扫参有两个细节一是每个工况要跑到足够的时长二是相邻工况之间的补偿度增量要足够细否则特征值阻尼的“由正变负”会被漏掉。参数取值说明串补度扫描步长0.01~0.02阻尼穿越曲线通常很尖复转矩扫频步长0.2~0.5 Hz太粗会把相位穿越点平滑掉单工况仿真时长2~5 s提取阻尼比需要足够的衰减数据记录数据截取长度至少 2 个最低模态周期FFT 分辨率不够时无法分离邻近模态有一点需要提醒扫参过程中如果发现某一档补偿度下模态阻尼快速变化不要直接用相邻两档线性插值要回到这一档重新细化补偿度。次同步区域内的模态耦合有时会是窄带强耦合线性插值出来的“临界串补度”并不可信。5. 结果可信度验证FFT 与阻尼比拟合时域仿真结束以后第一步不是看彩色曲线而是先做 FFT确认波形里确实存在轴系扭振频率。取 IP-LP1 轴段扭矩数据截取扰动结束后的 1~2s 窗口去掉直流分量然后做频谱分析fs 10000; % 如果 ode45 输出不规则先重采样到固定频率 y T_ip_lp1(2000:20000); % 只取扰动结束后的稳定段 N length(y); Y fft(y - mean(y)); f (0:N-1)*fs/N; Y_mag abs(Y(2:N/2)); [pks, loc] findpeaks(Y_mag, MinPeakHeight, 0.05*max(Y_mag), NPeaks, 4); disp(f(loc));这里最容易犯的错是整段数据一起做 FFT。如果信号同时包含衰减分量和增长分量一个 FFT 窗口会把两种趋势混在一起峰值频率看起来好像稳定实际模态已经在变化。正确做法是把时域数据切成多个窗口比较每个窗口里同一频率的幅值是在上升还是下降。确认频率后再做阻尼比拟合。对单一模态做带通滤波然后用最小二乘拟合衰减正弦信号model (p,t) p(1) * exp(p(2)*t) .* sin(2*pi*p(3)*t p(4)); p0 [max(y), -0.5, f_peak, 0]; p lsqnonlin((p) model(p,t)-y, p0); damping_ratio -p(2) / (2*pi*p(3));当 p(2) 为正时阻尼比是负值表示这个模态正在从电气侧吸收能量。把时域拟合出的频率和阻尼比和第二章特征值计算的结果放在一起对比如果差异超过 2%优先回查机械参数和单位是否统一而不是急着换积分算法。这个验证流程才是“数字仿真”真正能够支撑工程结论的最后一环。本文还有配套的精品资源点击获取 次同步谐振数字仿真：六质量块轴系建模与串补度扫参
简介电力系统次同步谐振的数字仿真PDF资料面向电力系统运行分析、电气工程相关专业的读者聚焦串联补偿线路引发的次同步谐振机理及数字仿真方法。内容包含同步发电机电气与轴系数学模型、隐式梯形积分法下的差分模型推导并基于MATLAB编程再现次同步谐振过程系统分析网络参数、串联补偿度对轴系扭振响应的影响为抑制次同步谐振提供理论依据与仿真参考。压缩包内提供单份PDF格式的学术论文大小仅1.1MB适合深入研读与对照推导。已有88人学习下载。读者可获得完整的机理分析、模型构建思路与仿真试验结论有助于在电力系统稳定分析、串补输电工程中快速建立相关认知并指导参数设计可作为相关课题研究的基础参考资料。1. 次同步谐振的数字仿真先弄清楚它在“谐振”什么做电力系统数字仿真的工程师迟早会遇到次同步谐振这个词。常见于装有串联补偿电容的输电线路上或风电变流器与弱网交互时标题里“数字仿真”通常不是指教科书例题而是要在时域里把发电机轴系的扭振过程复现出来。很多人按机电暂态模型搭出来发现频率对不上或波形看不到振荡问题大多出在机械侧自由度不足、电气网络等效过粗、时间步长取太大。这篇文章以电力系统次同步谐振数字仿真为主线从六质量块轴系建模到串联补偿度扫参再到结果可信度验证给出可直接改参数运行的最小化方案。2. 仿真之前先把“次同步谐振”拧成可计算的模型2.1 为什么是“多质量块轴系”而不是单刚体发电机模型次同步谐振的本质是机械侧和电气侧在低于工频的频率上交换能量。机械侧的能量载体不是发电机整机而是旋转轴系各段之间的扭转变形。如果仿真里只保留一个刚体转子就只能看到整个轴系整体摇摆的一个模态看不到高压缸、中压缸、低压缸和发电机之间相对扭转形成的扭振模态也就无从判断“哪个频率被激起来了”。工程上处理这个问题普遍把轴系离散成若干质量块相邻质量块之间用扭转弹簧连接每个质量块带一个惯性常数和一个阻尼系数。发电机模型可以用两质量块先跑通现象但要用于次同步谐振风险分析至少要做到四到六质量块。通常的划分方式是高压缸、中压缸、两个低压缸、发电机和励磁机一共六个质量块。因为励磁机质量小、刚度大很多资料里也会把励磁机去掉但对于数字仿真多保留一两个质量块的代价只是矩阵多几行没必要牺牲这一自由度。为什么机械阻尼在仿真里特别关键因为轴系本身的机械阻尼很小模态阻尼的主要来源反而在电气侧。当电气网络在某个次同步频率上呈现等效负阻尼时扭振模态就会从电气侧吸收能量幅值越振越大。数字仿真若只给一个常数正阻尼等于提前把最危险的失稳条件抹掉了。2.2 电气侧频率怎么算串补度把谐振点抬到哪电气侧的次同步谐振频率由串联补偿电容容抗和线路等值感抗共同决定。以最简单的单机无穷大系统为例补偿度为 k工频为 f0则电气谐振频率约为f_er f0 * sqrt(k)如果系统频率是 60Hz串补度 10% 时f_er 大约 19Hz串补度 20% 时大约 26.8Hz。这个区间和大型汽轮发电机组轴系的一阶到三阶扭振模态高度重叠所以“串补度越高越接近工频”并不意味着越安全而是要逐点看它落在哪个扭振模态附近。不过这只是粗略估算。真正决定是否失稳的是这个频率下电气回路的相位关系也就是电气转矩相对转子角变化的响应。严格做法是用小信号扫频得到次同步频段的等效阻抗再折算成复转矩系数数字仿真的价值恰恰就是在时域里把这个相位滞后出现的后果完整地演出来。2.3 给轴系定参数IEEE FBM 与模态估算代码做次同步谐振仿真最常用的是电力系统工程师熟知的 First Benchmark Model即 IEEE 次同步谐振标准模型。它有公开的机械参数适合用来验证算法和对比结果。下面这组参数是六质量块配置单位采用标幺值坐标统一到电弧度质块惯性常数 Hs相邻轴段扭转刚度 Kpu/radHP0.0929HP-IP19.303IP0.1556IP-LP134.929LP10.8582LP1-LP252.038LP20.8842LP2-GEN70.858GEN0.8685GEN-EXC2.822EXC0.0342拿到这组参数先别急着进时域仿真我一般会先用广义特征值估算一次轴系固有扭振频率确认模型自由度没有丢。下面的 MATLAB 代码可以完成这一步% freq_ssr.m —— 由机械参数估算轴系固有扭振频率 f0 60; % 系统额定频率 omega0 2*pi*f0; H [0.0929; 0.1556; 0.8582; 0.8842; 0.8685; 0.0342]; K [19.303; 34.929; 52.038; 70.858; 2.822]; M diag(2*H); Kmat zeros(6); for i 1:5 Kmat(i,i) Kmat(i,i) K(i); Kmat(i1,i1) Kmat(i1,i1) K(i); Kmat(i,i1) -K(i); Kmat(i1,i) -K(i); end [~, D] eig(omega0 * Kmat, M); % 广义特征值omega0 统一角度与角速度量纲 ev real(diag(D)); fn sort(sqrt(-ev(ev -1e-8)))/(2*pi); disp(轴系固有扭振频率(Hz)); disp(fn);代码里的 M 是惯性常数对角阵2H 是运动方程中需要的时间常数Kmat 是标准的刚性对称刚度矩阵每一段的刚度只有相邻两个自由度有耦合。eig 返回的特征值包含一个零特征值对应整根轴同步旋转的刚体模态对扭振没有意义所以过滤掉。负特征值开方后得到的就是各阶扭振频率单位是 Hz。这套计算里最容易出错的是量纲。如果 H 用秒、K 用标幺值、角度用电弧度那么必须乘上工频角速度 omega0否则算出来的频率会整体偏小或偏大。工程上拿到机组实测扭振频率后通常会用这个脚本先对一遍模型参数。3. 用 MATLAB/Simulink 搭一个 SSR 最小可复现仿真3.1 时域仿真和特征值分析谁先谁后次同步谐振数字仿真有两条并行路径一是线性化状态空间模型求特征值看各模态阻尼二是时域积分看扰动后扭振波形。特征值法快能直接判断某个模态是否失稳时域法直观能把非线性元件、控制保护和波形叠加都装进去。我的顺序通常是先做特征值分析把危险工况筛出来再用时域仿真去验证危险工况的波形和阻尼。这样既避免在安全工况上浪费积分时间也避免了只看特征值而漏掉饱和、死区或保护动作导致的非线性行为。实际操作中纯 Simulink 框图搭六质量块轴系并不方便我建议把机械部分写成 MATLAB 状态空间函数再在 Simulink 里用 State-Space 模块对接电气模型。3.2 六质量块状态空间搭建ode45 能跑的最小脚本下面这个脚本是“电力系统次同步谐振数字仿真”里最精简的可运行版本。它把六质量块轴系做成线性状态空间在发电机质块上施加短时制动转矩用来模拟故障切除后电气转矩对轴系的冲击。状态变量取为各质量块的转子电角度偏差和标幺转差外部输入是发电机质块上的电磁转矩增量。% ssr_torque_sim.m —— 六质量块轴系次同步扭振时域仿真 function ssr_torque_sim() f0 60; omega0 2*pi*f0; H [0.0929; 0.1556; 0.8582; 0.8842; 0.8685; 0.0342]; K [19.303; 34.929; 52.038; 70.858; 2.822]; M diag(2*H); Kmat build_Kmat(K); Dm diag([0.002; 0.003; 0.005; 0.005; 0.004; 0.001]); % 机械阻尼需按实测细化 Ks 1.2; % 电气同步转矩系数 De -0.6; % 次同步等效阻尼负值代表负阻尼 A [zeros(6), omega0*eye(6); -M\Kmat, -M\Dm]; % 电气转矩增量落在发电机质块索引5ΔTe Ks*Δδ De*Δω A(10,5) A(10,5) - Ks / M(5,5); A(10,10) A(10,10) - De / M(5,5); tStart 0.02; tEnd 0.08; amp -0.03; % 0.04s冲击转矩 Text (t) event_torque(6, 5, amp, tStart, tEnd, t); f_ode (t,x) A*x [zeros(6,1); M\Text(t)]; [t, x] ode45(f_ode, [0 3.0], zeros(12,1), odeset(MaxStep,1e-4)); delta x(:,1:6); T_shaft Kmat * delta; % 各轴段扭矩增量 plot(t, T_shaft(2,:)*100, LineWidth, 1); xlabel(时间 / s); ylabel(轴扭矩增量 / %); title(IP-LP1 轴段扭振响应); end function Kmat build_Kmat(K) Kmat zeros(6); for i 1:5 Kmat(i,i) Kmat(i,i) K(i); Kmat(i1,i1) Kmat(i1,i1) K(i); Kmat(i,i1) -K(i); Kmat(i1,i) -K(i); end end function Tq event_torque(n, idx, amp, t0, t1, t) Tq zeros(n,1); if t t0 t t1 Tq(idx) amp; end end这段代码最需要关注的是 A 矩阵最后两行的修正。发电机质块上的电磁转矩增量同时包含角度项和转速项角度项相当于给轴系增加了一个电气刚度转速项则是电气阻尼。把 De 设成负值后发电机质块会在转速偏离同步点时获得一个同向的转矩增量这就是次同步振荡能够自激的条件。改动 De 的符号和大小可以快速观察同一套机械参数在“正阻尼”和“负阻尼”两种场景下的区别。3.3 电气负阻尼是怎么进来的一个可调的 De常数 De 看起来太简单但它抓住了次同步谐振失稳的核心在特定次同步频率附近电气网络呈现的阻尼可能是负的。产生负阻尼的物理原因是串联补偿电容把系统谐振频率推到工频以下转子扭振产生的次同步电流经过网络阻抗后形成的电磁转矩和转子速度分量之间的相位关系不再构成制动转矩。严格的 De 不是常数而是频率的函数。工程上更常见的方法是扫频求复转矩把不同频率下的电气阻尼曲线拿到再和轴系机械模态叠加判断稳定性。上面脚本里的 De 相当于在某个危险频率上取了一个代表值适合做机理验证如果要把仿真结果用于工程结论必须在目标频率点上做网络阻抗扫频来更新 De。3.4 求解器与步长的配置表用 ode45 跑线性轴系模型很顺但一旦串上限流、饱和或控制坏就要换成刚性问题求解器。下面这组配置是我常用的起点参数推荐值说明仿真总时长3~5 s扭振模态衰减周期长短于 2s 看不到趋势MaxStep1e-4 s最高关心 50Hz每周期至少 200 个点RelTol1e-6默认容差容易丢掉弱阻尼模态AbsTol1e-8角度和转差量级差异大需要分开设置扰动时长0.04~0.08 s模仿故障切除后的电气冲击这里 MaxStep 不是随便压到 1e-4。如果只关心 20Hz 的扭振模态1ms 步长已经够要同时捕捉 60Hz 工频和 50Hz 扭转耦合步长必须低于 0.2ms。判断标准是最高关心频率的每个周期至少有 20 个采样点否则 FFT 结果会失真阻尼比拟合也不可靠。4. 扰动与串补度扫参把危险模式逼出来4.1 扰动怎么加用短时制动转矩而不是三相短路有些人会把故障事件直接设成三相短路但这在小信号轴系仿真里不是好选择。三相短路瞬间电气转矩几乎消失同时伴随复杂的电磁暂态对轴系模型的激励成分太杂很难定位到底激发了哪个扭振模态。我更常用的是短时制动转矩在发电机质块上施加一个幅值固定、时间很短的负向转矩相当于故障切除瞬间电气转矩突变。这样所有扭振模态都会被冲击起来随后按各自的阻尼衰减或增长。步骤可以固定成三条把扰动转矩幅值设为基础转矩的 1%5%避免进入非线性区。扰动持续时间取 4080ms覆盖次同步频段最高频率的多个周期。仿真启动后先让系统稳定 20ms 再加入扰动方便后续截取基线。如果必须模拟短路后串补重新投入需要把模型换到 Simulink/SimPowerSystems 里加真实线路和电容开关机械侧再用 S-Function 或状态空间模块对接不能继续用上面的线性 A 矩阵。4.2 把串补度从 5% 扫到 75%看电气谐振频率落到哪个模态上次同步谐振数字仿真的核心实操是扫串联补偿度。每改变一次补偿度电气谐振频率就移动一次可能跨过某个扭振模态。最简单的一轮扫参可以这样跑% scan_subsync.m —— 串补度与电气谐振频率对应关系 f0 60; k_range 0.05:0.05:0.75; for i 1:length(k_range) f_er f0 * sqrt(k_range(i)); fprintf(串补度 %.2f - 电气谐振频率 %.2f Hz\n, k_range(i), f_er); end这个循环只能定位“可能出现问题”的频率范围真正的风险判断要回到上一章的特征值脚本对每个串补度重算 De 和模态阻尼。常见做法是把这段扫频结果和轴系固有频率表放在一起做成一张风险对照表串补度电气谐振频率Hz潜在风险0.1523.2可能落入中频扭振模态0.3032.9与较高频扭振模态接近0.4540.2进入高段扭振区0.6046.5接近同步频率伴随感应发电机效应0.7552.0高于工频风险转向另一类问题工程上要注意表格里的 f_er 只是不计电阻的近似值。线路电阻和发电机次暂态电抗都会把实际谐振频率拉低几个 Hz所以扫参时建议把补偿度步长设到 1%2%否则可能跳过阻尼穿越点。4.3 扫参时的参数表步长、时长、扫频步长扫参不是把补偿度设置成数组然后循环那么简单。串补度扫参有两个细节一是每个工况要跑到足够的时长二是相邻工况之间的补偿度增量要足够细否则特征值阻尼的“由正变负”会被漏掉。参数取值说明串补度扫描步长0.01~0.02阻尼穿越曲线通常很尖复转矩扫频步长0.2~0.5 Hz太粗会把相位穿越点平滑掉单工况仿真时长2~5 s提取阻尼比需要足够的衰减数据记录数据截取长度至少 2 个最低模态周期FFT 分辨率不够时无法分离邻近模态有一点需要提醒扫参过程中如果发现某一档补偿度下模态阻尼快速变化不要直接用相邻两档线性插值要回到这一档重新细化补偿度。次同步区域内的模态耦合有时会是窄带强耦合线性插值出来的“临界串补度”并不可信。5. 结果可信度验证FFT 与阻尼比拟合时域仿真结束以后第一步不是看彩色曲线而是先做 FFT确认波形里确实存在轴系扭振频率。取 IP-LP1 轴段扭矩数据截取扰动结束后的 1~2s 窗口去掉直流分量然后做频谱分析fs 10000; % 如果 ode45 输出不规则先重采样到固定频率 y T_ip_lp1(2000:20000); % 只取扰动结束后的稳定段 N length(y); Y fft(y - mean(y)); f (0:N-1)*fs/N; Y_mag abs(Y(2:N/2)); [pks, loc] findpeaks(Y_mag, MinPeakHeight, 0.05*max(Y_mag), NPeaks, 4); disp(f(loc));这里最容易犯的错是整段数据一起做 FFT。如果信号同时包含衰减分量和增长分量一个 FFT 窗口会把两种趋势混在一起峰值频率看起来好像稳定实际模态已经在变化。正确做法是把时域数据切成多个窗口比较每个窗口里同一频率的幅值是在上升还是下降。确认频率后再做阻尼比拟合。对单一模态做带通滤波然后用最小二乘拟合衰减正弦信号model (p,t) p(1) * exp(p(2)*t) .* sin(2*pi*p(3)*t p(4)); p0 [max(y), -0.5, f_peak, 0]; p lsqnonlin((p) model(p,t)-y, p0); damping_ratio -p(2) / (2*pi*p(3));当 p(2) 为正时阻尼比是负值表示这个模态正在从电气侧吸收能量。把时域拟合出的频率和阻尼比和第二章特征值计算的结果放在一起对比如果差异超过 2%优先回查机械参数和单位是否统一而不是急着换积分算法。这个验证流程才是“数字仿真”真正能够支撑工程结论的最后一环。本文还有配套的精品资源点击获取