简介这份资源面向电力系统研究、储能技术开发领域的工程师与研究人员聚焦可再生能源高占比下的一次调频难题复现并优化了基于自适应SOC的电池-飞轮混合储能控制策略。内容结合正/负虚拟惯性控制与虚拟下垂控制通过SOC自适应调整控制参数权重将系统频率波动限制在±0.2Hz内并有效防止电池过充过放。资源包为1个PDF文件约736KB内含完整可运行的Python代码及逐段解释涵盖类初始化、电力系统频率响应模型、自适应控制策略、SOC更新机制与仿真流程并借助Matplotlib完成结果可视化。读者可据此理解权重分配策略、SOC自适应调节、混合控制模式切换及SOC安全保护机制等关键技术点评估不同策略下的频率响应与SOC变化进而优化调频性能、延长储能设备寿命。目前已有155人学习适合在掌握相关理论后逐步实践代码并对照仿真结果深入分析。1. 从一次 0.2pu 阶跃扰动说起这套混合储能调频代码到底在算什么电网里可再生能源渗透率一高一次调频的压力就从调得动变成调得准还得护设备。这份复现资源盯的就是这个矛盾电池能量密度高但怕过充过放飞轮功率密度高但储能量小单靠一种储能做一次调频要么 SOC 很快撞上下限要么频率压不住。它用 Python 把一篇论文里的自适应 SOC 控制策略完整跑了一遍核心是正/负虚拟惯性和虚拟下垂两路控制按频率状态动态分配权重再乘上随 SOC 变化的出力系数最后把频率偏差摁在 ±0.2Hz 内同时让电池 SOC 别掉出 0.2~1.0 的安全区。适合做电力系统仿真、储能控制策略验证的工程师尤其是想拿一份能直接跑、能改参数、能看曲线的代码当起点的人。2. 电力系统频率响应模型与 lfilter 离散化实现2.1 为什么用二阶传递函数而不是完整机组模型一次调频仿真里最容易被过度建模的就是同步机组。完整模型要写调速器、汽轮机、再热器一串环节参数多到没法从论文里对齐。这份代码走的是经典简化路线把系统等效成一个惯性加阻尼的二阶环节传递函数写成Δf(s)/ΔP(s) 1/(2Hs D)其中 H 是惯性常数、D 是阻尼系数。这个式子的物理含义很直白——功率缺额先被转子动能吸收表现为频率变化率然后阻尼项慢慢把频率拉向稳态偏差。对一次调频这种秒级过程二阶模型足够反映频率跌多深、跌多快而且参数只有两个复现时不容易跑偏。代码里power_system_model用scipy.signal.lfilter做离散化而不是手写差分方程。lfilter对num[1, D]、den[2H, D]这种低阶传递函数等价于双线性变换后的 IIR 滤波省去了自己推 Z 变换的麻烦。要注意的是lfilter默认从零初始条件开始所以第一次调用时输出会有一个瞬态这也是为什么仿真里负荷扰动设在第 5 秒而不是第 0 秒——给系统留出稳定时间。2.2 参数取值与物理量纲对齐参数代码取值物理含义调整影响H5.0 s系统惯性常数越大频率变化越慢初始 RoCoF 越小D1.0阻尼系数越大稳态频差越小但响应变钝R0.05调差系数决定一次调频稳态出力比例Tg0.2 s调速器时间常数影响机组出力爬升速度Tch0.3 s汽轮机时间常数影响机械功率滞后dt0.01 s仿真步长需小于最小时间常数的 1/10这张表里 R、Tg、Tch 在当前简化模型里其实没直接进传递函数但它们是后续如果要接完整机组模型时的接口参数。我一般会先把 H 和 D 调到位让空载阶跃下的稳态频差和理论值ΔP/D对上再去动其他参数。比如 0.2pu 扰动、D1.0 时理论稳态频差约 0.2Hz正好卡在论文说的 ±0.2Hz 边界上这说明储能必须介入才能压下去。2.3 离散化步长与数值稳定性import numpy as np from scipy.signal import lfilter # 系统参数 H 5.0 # 惯性常数(s) D 1.0 # 阻尼系数 dt 0.01 # 仿真步长(s) total_time 30 time np.arange(0, total_time, dt) # 负荷扰动第5秒起0.2pu阶跃 delta_PL np.zeros_like(time) delta_PL[int(5/dt):] 0.2 # 传递函数离散化Δf(s)/ΔP(s) 1/(2Hs D) num [1, D] den [2*H, D] delta_f lfilter(num, den, delta_PL) # 检查稳态频差 print(稳态频差理论值:, 0.2/D, Hz) print(仿真末值:, delta_f[-1], Hz)这段代码先构造扰动序列再用lfilter一次性算出全时段频率响应。num和den的系数顺序是降幂排列[1, D]对应s D[2H, D]对应2Hs D。运行后如果仿真末值接近 0.2Hz说明模型量纲对上了如果差一个数量级多半是 H 或 D 的单位搞混了——H 用秒D 用标幺值别把 D 当成有名值。步长 dt0.01s 对 Tg0.2s 来说有 20 个采样点够用如果后面要接更快的飞轮动态建议降到 0.001s否则会看到数值振荡。提示lfilter的初始条件默认全零等价于系统从静止启动。如果想让仿真从稳态开始可以用lfilter_zi算初始状态再传入zi参数。3. 自适应 SOC 权重分配与虚拟惯性/下垂控制实现3.1 频率恶化与恢复阶段的判据这套策略最核心的判断就一个d_delta_f * delta_f 0说明频率正在远离额定值属于恶化阶段反之是恢复阶段。恶化阶段要猛出力所以虚拟惯性权重alpha1用exp(-10*|Δf|)频率偏差越大权重反而越小——这看起来反直觉其实是因为偏差大时下垂控制更有效惯性控制主要负责抑制变化率。恢复阶段则用 Sigmoid 函数1/(1exp(5000*|Δf|))把惯性权重压到接近零避免频率往回走时储能还在猛冲导致超调。代码里delta_f_db 0.015Hz是死区小于这个值直接不动作防止储能频繁微调磨损设备。delta_f_set 0.2Hz是额定偏差界限超过它权重函数切换成另一套 Sigmoid本质是进入保频率优先模式暂时牺牲 SOC 均衡。3.2 SOC 自适应系数与 Logistic 衰减def soc_adjustment(soc, soc_min, soc_low, soc_high, soc_max, K_max, n7): SOC自适应调节系数安全区满出力临界区非线性衰减 if soc soc_min: return 0 elif soc_min soc soc_low: # 低SOC临界区指数衰减 term np.exp(n * (soc_low - soc) / 0.1) return (1 0.2) * K_max / 2 * term / (term 1) elif soc_low soc soc_high: # 安全区满出力 return K_max elif soc_high soc soc_max: # 高SOC临界区指数衰减 term np.exp(n * (soc_max - soc) / 0.1) return K_max * 0.01 * term / (term 1) else: return 0这个函数把 SOC 分成五段低于下限直接闭锁下限到低警戒线之间用指数函数平滑爬升安全区给满系数高警戒线到上限之间衰减到 1% 左右超过上限闭锁。n7控制衰减陡度n 越大过渡越急。我一般会把soc_low和soc_high设成 0.4 和 0.6给安全区留 20% 的缓冲这样电池在大多数扰动下都能满出力只有连续同向扰动才会触发衰减。飞轮的soc_low/high可以放宽到 0.3/0.7因为它循环寿命长不怕深充深放。3.3 混合出力合成与 SOC 更新class HybridESS: def __init__(self): self.battery_SOC 0.5 self.flywheel_SOC 0.5 self.battery_min, self.battery_max 0.2, 1.0 self.flywheel_min, self.flywheel_max 0.05, 0.95 self.Km, self.Kk 0.8, 0.6 self.delta_f_db, self.delta_f_set 0.015, 0.2 self.dt 0.01 def adaptive_control(self, delta_f, d_delta_f): if abs(delta_f) self.delta_f_db: alpha1, alpha2 0, 0 elif abs(delta_f) self.delta_f_set: if d_delta_f * delta_f 0: alpha1 np.exp(-10 * abs(delta_f)) alpha2 1 - alpha1 else: alpha1 1 / (1 np.exp(5000 * abs(delta_f))) alpha2 1 - alpha1 else: alpha1 1 / (1 np.exp(5000 * abs(delta_f))) alpha2 1 - alpha1 K_battery 1 - 2 * abs(self.battery_SOC - 0.5) K_flywheel 1 - 2 * abs(self.flywheel_SOC - 0.5) delta_P_MH -self.Km * d_delta_f * alpha1 delta_P_KH -self.Kk * delta_f * alpha2 return K_battery * delta_P_MH K_flywheel * delta_P_KH def update_SOC(self, delta_P_H): self.battery_SOC - delta_P_H * self.dt / 1000 self.battery_SOC np.clip(self.battery_SOC, self.battery_min, self.battery_max) self.flywheel_SOC - delta_P_H * self.dt / 100 self.flywheel_SOC np.clip(self.flywheel_SOC, self.flywheel_min, self.flywheel_max)adaptive_control先算权重再算 SOC 系数最后把惯性项和下垂项按权重合成。注意K_battery和K_flywheel用的是1 - 2*|SOC-0.5|SOC 越偏离 0.5 系数越小这是最简单的线性 SOC 均衡。update_SOC里电池容量按 1000kWh、飞轮按 100kWh 折算所以同样出力下飞轮 SOC 变化快 10 倍——这符合功率型储能的特性但也意味着飞轮更容易撞限实际工程里会给飞轮加更宽的 SOC 区间或者更快的恢复策略。注意update_SOC里电池和飞轮用的是同一个delta_P_H没有按各自出力拆分。如果要精确模拟应该把K_battery*delta_P_MH和K_flywheel*delta_P_KH分开更新各自的 SOC否则飞轮的 SOC 变化会被电池的出力污染。4. 仿真主循环、结果可视化与常见复现坑4.1 主循环里那个容易写错的频率更新原始代码里delta_f[i] self.power_system_model(delta_PL[:i1] - delta_P_H[:i1])[-1]这一行每次循环都重新跑一遍全时段lfilter计算量是 O(n²)。30 秒 3000 步还能忍但如果仿真时长拉到 300 秒耗时会明显上去。更麻烦的是这种写法隐含了净功率序列已知的假设而实际上delta_P_H是逐步算出来的所以每次重跑lfilter时前面的储能出力已经固定只有最后一步是新的。正确做法是用状态空间递推或者维护lfilter的zi状态每步只算一次。from scipy.signal import lfilter, lfilter_zi # 初始化滤波器状态 zi lfilter_zi(num, den) * 0 # 零初始条件 delta_f np.zeros_like(time) delta_P_H np.zeros_like(time) for i in range(1, len(time)): d_delta_f (delta_f[i-1] - delta_f[i-2]) / dt if i 1 else 0 delta_P_H[i] ess.adaptive_control(delta_f[i-1], d_delta_f) net_power delta_PL[i] - delta_P_H[i] # 单步递推只算当前时刻 delta_f[i], zi lfilter(num, den, [net_power], zizi) delta_f[i] delta_f[i][0] ess.update_SOC(delta_P_H[i])这样每步只做一次一阶 IIR 计算复杂度降到 O(n)。zi是滤波器内部状态必须跨步保留否则每步都从零初始条件重启频率响应会碎成一段段瞬态。这个坑我在复现其他调频代码时踩过不止一次表现是频率曲线在扰动后出现高频毛刺其实就是状态没接上。4.2 可视化三件套与判读要点import matplotlib.pyplot as plt fig, axes plt.subplots(3, 1, figsize(12, 8), sharexTrue) axes[0].plot(time, delta_f) axes[0].axhline(y0.2, colorr, linestyle--) axes[0].axhline(y-0.2, colorr, linestyle--) axes[0].set_ylabel(Frequency Deviation (Hz)) axes[0].set_title(System Frequency Response) axes[1].plot(time, delta_P_H) axes[1].set_ylabel(Power Output (pu)) axes[1].set_title(Hybrid ESS Power Output) axes[2].plot(time, battery_SOC, labelBattery SOC) axes[2].plot(time, flywheel_SOC, labelFlywheel SOC) axes[2].legend() axes[2].set_ylabel(SOC) axes[2].set_xlabel(Time (s)) axes[2].set_title(State of Charge) plt.tight_layout() plt.show()三张图分别看频率、出力、SOC。频率图重点看两个指标最大频差有没有超过 0.2Hz 红线以及恢复时间——从扰动发生到频率回到死区内用了多久。出力图看储能是先猛后缓还是全程平推前者说明惯性控制起作用了后者说明下垂占主导。SOC 图看两条曲线的分离程度如果电池 SOC 几乎不动而飞轮 SOC 大幅波动说明 SOC 系数把出力更多分给了飞轮这是符合设计意图的反过来则要检查K_battery和K_flywheel的初始值是不是设反了。4.3 复现时最常遇到的三个问题第一个是频率曲线不收敛或者发散。九成是lfilter的den系数符号写反了2Hs D对应[2H, D]如果写成[2H, -D]系统就变成不稳定极点。第二个是 SOC 曲线一条直线不动。检查update_SOC里的容量除数1000 和 100 差一个数量级如果两个都写成 1000飞轮 SOC 变化会小到看不出来。第三个是权重alpha1和alpha2加起来不等于 1。原始代码在恢复阶段用alpha2 1 - alpha1是对的但扩展版里出现过alpha2 1 1/(1exp(...))的写法那会导致总出力翻倍频率过冲。复现时以原始代码为准扩展版当参考。提示如果跑出来的频率稳态值不是 0.2Hz 而是接近 0说明储能出力把扰动完全抵消了这通常意味着Km或Kk设得过大实际工程里一次调频允许留有稳态偏差不会完全归零。5. 从仿真到工程SOC 安全监控与参数整定技巧5.1 给 SOC 加一层独立的安全监控仿真里 SOC 靠np.clip硬限幅工程上不能这么干因为硬限幅意味着功率突变对变流器不友好。常见做法是加一层独立的安全监控在 SOC 接近临界值之前就主动降出力而不是等到撞限才切。def soc_safety_monitor(soc, device_type): SOC安全监控返回建议的出力缩放系数 thresholds { battery: {warning: 0.25, critical: 0.15}, flywheel: {warning: 0.15, critical: 0.05} } th thresholds[device_type] if soc th[critical]: return 0.0 # 紧急闭锁 elif soc th[warning]: # 线性降出力从warning到critical系数从1降到0 return (soc - th[critical]) / (th[warning] - th[critical]) else: return 1.0这个函数返回一个 0~1 的缩放系数乘到 SOC 自适应系数上。好处是降出力是连续的不会出现功率阶跃。warning和critical的取值要比soc_min高一点留出反应时间。电池的 warning 设 0.25、critical 设 0.15对应soc_min0.2这样在 SOC 掉到 0.2 之前就已经开始降出力了。5.2 参数整定的三个经验值Km和Kk的比例决定惯性控制和下垂控制的相对强度。我一般从Km/Kk 1.5开始试因为惯性控制抑制的是频率变化率对一次调频的初始阶段更重要。如果频率最大偏差压不下去先加Kk如果频率变化率太大导致 RoCoF 保护动作先加Km。delta_f_db死区别设太小0.01~0.02Hz 比较合适太小会让储能频繁动作太大则调频死区超标。delta_f_set按并网标准来多数场景是 0.2Hz但如果是弱电网可以放宽到 0.5Hz。5.3 用阶跃加随机扰动验证鲁棒性单一阶跃扰动只能验证稳态性能实际电网的扰动是连续的。把阶跃和随机波动叠加能看出控制策略在噪声下的表现np.random.seed(42) step_disturbance 0.2 * (time 5) random_disturbance 0.02 * np.random.randn(len(time)) disturbance step_disturbance random_disturbance # 用带安全监控的控制器跑一遍 # 重点观察SOC是否在随机扰动下持续漂移随机扰动幅值取阶跃的 10% 左右比较合理太大就变成连续大扰动不是一次调频的典型场景。跑完后重点看 SOC 曲线的长期趋势——如果电池 SOC 在 30 秒内单调下降超过 0.05说明 SOC 均衡策略没起作用需要调大K_battery的 SOC 偏离惩罚项或者给电池加一个慢速恢复回路。飞轮 SOC 波动大是正常的但均值应该稳定在 0.5 附近如果持续偏向一侧说明K_flywheel的符号或者权重分配有问题。本文还有配套的精品资源点击获取
混合储能一次调频：自适应SOC控制与虚拟惯性/下垂策略Python复现
简介这份资源面向电力系统研究、储能技术开发领域的工程师与研究人员聚焦可再生能源高占比下的一次调频难题复现并优化了基于自适应SOC的电池-飞轮混合储能控制策略。内容结合正/负虚拟惯性控制与虚拟下垂控制通过SOC自适应调整控制参数权重将系统频率波动限制在±0.2Hz内并有效防止电池过充过放。资源包为1个PDF文件约736KB内含完整可运行的Python代码及逐段解释涵盖类初始化、电力系统频率响应模型、自适应控制策略、SOC更新机制与仿真流程并借助Matplotlib完成结果可视化。读者可据此理解权重分配策略、SOC自适应调节、混合控制模式切换及SOC安全保护机制等关键技术点评估不同策略下的频率响应与SOC变化进而优化调频性能、延长储能设备寿命。目前已有155人学习适合在掌握相关理论后逐步实践代码并对照仿真结果深入分析。1. 从一次 0.2pu 阶跃扰动说起这套混合储能调频代码到底在算什么电网里可再生能源渗透率一高一次调频的压力就从调得动变成调得准还得护设备。这份复现资源盯的就是这个矛盾电池能量密度高但怕过充过放飞轮功率密度高但储能量小单靠一种储能做一次调频要么 SOC 很快撞上下限要么频率压不住。它用 Python 把一篇论文里的自适应 SOC 控制策略完整跑了一遍核心是正/负虚拟惯性和虚拟下垂两路控制按频率状态动态分配权重再乘上随 SOC 变化的出力系数最后把频率偏差摁在 ±0.2Hz 内同时让电池 SOC 别掉出 0.2~1.0 的安全区。适合做电力系统仿真、储能控制策略验证的工程师尤其是想拿一份能直接跑、能改参数、能看曲线的代码当起点的人。2. 电力系统频率响应模型与 lfilter 离散化实现2.1 为什么用二阶传递函数而不是完整机组模型一次调频仿真里最容易被过度建模的就是同步机组。完整模型要写调速器、汽轮机、再热器一串环节参数多到没法从论文里对齐。这份代码走的是经典简化路线把系统等效成一个惯性加阻尼的二阶环节传递函数写成Δf(s)/ΔP(s) 1/(2Hs D)其中 H 是惯性常数、D 是阻尼系数。这个式子的物理含义很直白——功率缺额先被转子动能吸收表现为频率变化率然后阻尼项慢慢把频率拉向稳态偏差。对一次调频这种秒级过程二阶模型足够反映频率跌多深、跌多快而且参数只有两个复现时不容易跑偏。代码里power_system_model用scipy.signal.lfilter做离散化而不是手写差分方程。lfilter对num[1, D]、den[2H, D]这种低阶传递函数等价于双线性变换后的 IIR 滤波省去了自己推 Z 变换的麻烦。要注意的是lfilter默认从零初始条件开始所以第一次调用时输出会有一个瞬态这也是为什么仿真里负荷扰动设在第 5 秒而不是第 0 秒——给系统留出稳定时间。2.2 参数取值与物理量纲对齐参数代码取值物理含义调整影响H5.0 s系统惯性常数越大频率变化越慢初始 RoCoF 越小D1.0阻尼系数越大稳态频差越小但响应变钝R0.05调差系数决定一次调频稳态出力比例Tg0.2 s调速器时间常数影响机组出力爬升速度Tch0.3 s汽轮机时间常数影响机械功率滞后dt0.01 s仿真步长需小于最小时间常数的 1/10这张表里 R、Tg、Tch 在当前简化模型里其实没直接进传递函数但它们是后续如果要接完整机组模型时的接口参数。我一般会先把 H 和 D 调到位让空载阶跃下的稳态频差和理论值ΔP/D对上再去动其他参数。比如 0.2pu 扰动、D1.0 时理论稳态频差约 0.2Hz正好卡在论文说的 ±0.2Hz 边界上这说明储能必须介入才能压下去。2.3 离散化步长与数值稳定性import numpy as np from scipy.signal import lfilter # 系统参数 H 5.0 # 惯性常数(s) D 1.0 # 阻尼系数 dt 0.01 # 仿真步长(s) total_time 30 time np.arange(0, total_time, dt) # 负荷扰动第5秒起0.2pu阶跃 delta_PL np.zeros_like(time) delta_PL[int(5/dt):] 0.2 # 传递函数离散化Δf(s)/ΔP(s) 1/(2Hs D) num [1, D] den [2*H, D] delta_f lfilter(num, den, delta_PL) # 检查稳态频差 print(稳态频差理论值:, 0.2/D, Hz) print(仿真末值:, delta_f[-1], Hz)这段代码先构造扰动序列再用lfilter一次性算出全时段频率响应。num和den的系数顺序是降幂排列[1, D]对应s D[2H, D]对应2Hs D。运行后如果仿真末值接近 0.2Hz说明模型量纲对上了如果差一个数量级多半是 H 或 D 的单位搞混了——H 用秒D 用标幺值别把 D 当成有名值。步长 dt0.01s 对 Tg0.2s 来说有 20 个采样点够用如果后面要接更快的飞轮动态建议降到 0.001s否则会看到数值振荡。提示lfilter的初始条件默认全零等价于系统从静止启动。如果想让仿真从稳态开始可以用lfilter_zi算初始状态再传入zi参数。3. 自适应 SOC 权重分配与虚拟惯性/下垂控制实现3.1 频率恶化与恢复阶段的判据这套策略最核心的判断就一个d_delta_f * delta_f 0说明频率正在远离额定值属于恶化阶段反之是恢复阶段。恶化阶段要猛出力所以虚拟惯性权重alpha1用exp(-10*|Δf|)频率偏差越大权重反而越小——这看起来反直觉其实是因为偏差大时下垂控制更有效惯性控制主要负责抑制变化率。恢复阶段则用 Sigmoid 函数1/(1exp(5000*|Δf|))把惯性权重压到接近零避免频率往回走时储能还在猛冲导致超调。代码里delta_f_db 0.015Hz是死区小于这个值直接不动作防止储能频繁微调磨损设备。delta_f_set 0.2Hz是额定偏差界限超过它权重函数切换成另一套 Sigmoid本质是进入保频率优先模式暂时牺牲 SOC 均衡。3.2 SOC 自适应系数与 Logistic 衰减def soc_adjustment(soc, soc_min, soc_low, soc_high, soc_max, K_max, n7): SOC自适应调节系数安全区满出力临界区非线性衰减 if soc soc_min: return 0 elif soc_min soc soc_low: # 低SOC临界区指数衰减 term np.exp(n * (soc_low - soc) / 0.1) return (1 0.2) * K_max / 2 * term / (term 1) elif soc_low soc soc_high: # 安全区满出力 return K_max elif soc_high soc soc_max: # 高SOC临界区指数衰减 term np.exp(n * (soc_max - soc) / 0.1) return K_max * 0.01 * term / (term 1) else: return 0这个函数把 SOC 分成五段低于下限直接闭锁下限到低警戒线之间用指数函数平滑爬升安全区给满系数高警戒线到上限之间衰减到 1% 左右超过上限闭锁。n7控制衰减陡度n 越大过渡越急。我一般会把soc_low和soc_high设成 0.4 和 0.6给安全区留 20% 的缓冲这样电池在大多数扰动下都能满出力只有连续同向扰动才会触发衰减。飞轮的soc_low/high可以放宽到 0.3/0.7因为它循环寿命长不怕深充深放。3.3 混合出力合成与 SOC 更新class HybridESS: def __init__(self): self.battery_SOC 0.5 self.flywheel_SOC 0.5 self.battery_min, self.battery_max 0.2, 1.0 self.flywheel_min, self.flywheel_max 0.05, 0.95 self.Km, self.Kk 0.8, 0.6 self.delta_f_db, self.delta_f_set 0.015, 0.2 self.dt 0.01 def adaptive_control(self, delta_f, d_delta_f): if abs(delta_f) self.delta_f_db: alpha1, alpha2 0, 0 elif abs(delta_f) self.delta_f_set: if d_delta_f * delta_f 0: alpha1 np.exp(-10 * abs(delta_f)) alpha2 1 - alpha1 else: alpha1 1 / (1 np.exp(5000 * abs(delta_f))) alpha2 1 - alpha1 else: alpha1 1 / (1 np.exp(5000 * abs(delta_f))) alpha2 1 - alpha1 K_battery 1 - 2 * abs(self.battery_SOC - 0.5) K_flywheel 1 - 2 * abs(self.flywheel_SOC - 0.5) delta_P_MH -self.Km * d_delta_f * alpha1 delta_P_KH -self.Kk * delta_f * alpha2 return K_battery * delta_P_MH K_flywheel * delta_P_KH def update_SOC(self, delta_P_H): self.battery_SOC - delta_P_H * self.dt / 1000 self.battery_SOC np.clip(self.battery_SOC, self.battery_min, self.battery_max) self.flywheel_SOC - delta_P_H * self.dt / 100 self.flywheel_SOC np.clip(self.flywheel_SOC, self.flywheel_min, self.flywheel_max)adaptive_control先算权重再算 SOC 系数最后把惯性项和下垂项按权重合成。注意K_battery和K_flywheel用的是1 - 2*|SOC-0.5|SOC 越偏离 0.5 系数越小这是最简单的线性 SOC 均衡。update_SOC里电池容量按 1000kWh、飞轮按 100kWh 折算所以同样出力下飞轮 SOC 变化快 10 倍——这符合功率型储能的特性但也意味着飞轮更容易撞限实际工程里会给飞轮加更宽的 SOC 区间或者更快的恢复策略。注意update_SOC里电池和飞轮用的是同一个delta_P_H没有按各自出力拆分。如果要精确模拟应该把K_battery*delta_P_MH和K_flywheel*delta_P_KH分开更新各自的 SOC否则飞轮的 SOC 变化会被电池的出力污染。4. 仿真主循环、结果可视化与常见复现坑4.1 主循环里那个容易写错的频率更新原始代码里delta_f[i] self.power_system_model(delta_PL[:i1] - delta_P_H[:i1])[-1]这一行每次循环都重新跑一遍全时段lfilter计算量是 O(n²)。30 秒 3000 步还能忍但如果仿真时长拉到 300 秒耗时会明显上去。更麻烦的是这种写法隐含了净功率序列已知的假设而实际上delta_P_H是逐步算出来的所以每次重跑lfilter时前面的储能出力已经固定只有最后一步是新的。正确做法是用状态空间递推或者维护lfilter的zi状态每步只算一次。from scipy.signal import lfilter, lfilter_zi # 初始化滤波器状态 zi lfilter_zi(num, den) * 0 # 零初始条件 delta_f np.zeros_like(time) delta_P_H np.zeros_like(time) for i in range(1, len(time)): d_delta_f (delta_f[i-1] - delta_f[i-2]) / dt if i 1 else 0 delta_P_H[i] ess.adaptive_control(delta_f[i-1], d_delta_f) net_power delta_PL[i] - delta_P_H[i] # 单步递推只算当前时刻 delta_f[i], zi lfilter(num, den, [net_power], zizi) delta_f[i] delta_f[i][0] ess.update_SOC(delta_P_H[i])这样每步只做一次一阶 IIR 计算复杂度降到 O(n)。zi是滤波器内部状态必须跨步保留否则每步都从零初始条件重启频率响应会碎成一段段瞬态。这个坑我在复现其他调频代码时踩过不止一次表现是频率曲线在扰动后出现高频毛刺其实就是状态没接上。4.2 可视化三件套与判读要点import matplotlib.pyplot as plt fig, axes plt.subplots(3, 1, figsize(12, 8), sharexTrue) axes[0].plot(time, delta_f) axes[0].axhline(y0.2, colorr, linestyle--) axes[0].axhline(y-0.2, colorr, linestyle--) axes[0].set_ylabel(Frequency Deviation (Hz)) axes[0].set_title(System Frequency Response) axes[1].plot(time, delta_P_H) axes[1].set_ylabel(Power Output (pu)) axes[1].set_title(Hybrid ESS Power Output) axes[2].plot(time, battery_SOC, labelBattery SOC) axes[2].plot(time, flywheel_SOC, labelFlywheel SOC) axes[2].legend() axes[2].set_ylabel(SOC) axes[2].set_xlabel(Time (s)) axes[2].set_title(State of Charge) plt.tight_layout() plt.show()三张图分别看频率、出力、SOC。频率图重点看两个指标最大频差有没有超过 0.2Hz 红线以及恢复时间——从扰动发生到频率回到死区内用了多久。出力图看储能是先猛后缓还是全程平推前者说明惯性控制起作用了后者说明下垂占主导。SOC 图看两条曲线的分离程度如果电池 SOC 几乎不动而飞轮 SOC 大幅波动说明 SOC 系数把出力更多分给了飞轮这是符合设计意图的反过来则要检查K_battery和K_flywheel的初始值是不是设反了。4.3 复现时最常遇到的三个问题第一个是频率曲线不收敛或者发散。九成是lfilter的den系数符号写反了2Hs D对应[2H, D]如果写成[2H, -D]系统就变成不稳定极点。第二个是 SOC 曲线一条直线不动。检查update_SOC里的容量除数1000 和 100 差一个数量级如果两个都写成 1000飞轮 SOC 变化会小到看不出来。第三个是权重alpha1和alpha2加起来不等于 1。原始代码在恢复阶段用alpha2 1 - alpha1是对的但扩展版里出现过alpha2 1 1/(1exp(...))的写法那会导致总出力翻倍频率过冲。复现时以原始代码为准扩展版当参考。提示如果跑出来的频率稳态值不是 0.2Hz 而是接近 0说明储能出力把扰动完全抵消了这通常意味着Km或Kk设得过大实际工程里一次调频允许留有稳态偏差不会完全归零。5. 从仿真到工程SOC 安全监控与参数整定技巧5.1 给 SOC 加一层独立的安全监控仿真里 SOC 靠np.clip硬限幅工程上不能这么干因为硬限幅意味着功率突变对变流器不友好。常见做法是加一层独立的安全监控在 SOC 接近临界值之前就主动降出力而不是等到撞限才切。def soc_safety_monitor(soc, device_type): SOC安全监控返回建议的出力缩放系数 thresholds { battery: {warning: 0.25, critical: 0.15}, flywheel: {warning: 0.15, critical: 0.05} } th thresholds[device_type] if soc th[critical]: return 0.0 # 紧急闭锁 elif soc th[warning]: # 线性降出力从warning到critical系数从1降到0 return (soc - th[critical]) / (th[warning] - th[critical]) else: return 1.0这个函数返回一个 0~1 的缩放系数乘到 SOC 自适应系数上。好处是降出力是连续的不会出现功率阶跃。warning和critical的取值要比soc_min高一点留出反应时间。电池的 warning 设 0.25、critical 设 0.15对应soc_min0.2这样在 SOC 掉到 0.2 之前就已经开始降出力了。5.2 参数整定的三个经验值Km和Kk的比例决定惯性控制和下垂控制的相对强度。我一般从Km/Kk 1.5开始试因为惯性控制抑制的是频率变化率对一次调频的初始阶段更重要。如果频率最大偏差压不下去先加Kk如果频率变化率太大导致 RoCoF 保护动作先加Km。delta_f_db死区别设太小0.01~0.02Hz 比较合适太小会让储能频繁动作太大则调频死区超标。delta_f_set按并网标准来多数场景是 0.2Hz但如果是弱电网可以放宽到 0.5Hz。5.3 用阶跃加随机扰动验证鲁棒性单一阶跃扰动只能验证稳态性能实际电网的扰动是连续的。把阶跃和随机波动叠加能看出控制策略在噪声下的表现np.random.seed(42) step_disturbance 0.2 * (time 5) random_disturbance 0.02 * np.random.randn(len(time)) disturbance step_disturbance random_disturbance # 用带安全监控的控制器跑一遍 # 重点观察SOC是否在随机扰动下持续漂移随机扰动幅值取阶跃的 10% 左右比较合理太大就变成连续大扰动不是一次调频的典型场景。跑完后重点看 SOC 曲线的长期趋势——如果电池 SOC 在 30 秒内单调下降超过 0.05说明 SOC 均衡策略没起作用需要调大K_battery的 SOC 偏离惩罚项或者给电池加一个慢速恢复回路。飞轮 SOC 波动大是正常的但均值应该稳定在 0.5 附近如果持续偏向一侧说明K_flywheel的符号或者权重分配有问题。本文还有配套的精品资源点击获取