简介面向电气工程与自动化专业毕业设计这套MATLAB源程序围绕基于氢储能的热电联供型微电网优化调度课题提供了完整的两阶段调度实现方案。程序针对质子交换膜燃料电池与电解槽的热电联供特性搭建了含光伏、风机、燃料电池、蓄电池、电锅炉及燃气锅炉的微电网系统重点体现电氢转换时的余热回收与热电氢耦合机制可支撑电、热、氢三种能量的协调利用与相互转换。第一阶段调度基于日前风光发电与负荷预测以运行成本最小为目标采用混合整数线性规划完成全局优化第二阶段调度利用超短期预测结果以模型预测控制嵌入混合整数二次规划算法有效弱化预测误差带来的经济性影响。包体共4个文件含2个M源文件及2张PNG结果图总大小仅377KB轻量易读。已有142人学习浏览包含可运行的两阶段调度MATLAB代码与结果示意图适合正在开展微电网优化调度、氢储能或综合能源系统方向毕业设计的学生参考与二次开发。1. 氢储能热电联供微电网优化调度论文能跑通不等于你能复现我第一次把氢储能放进热电联供微电网模型里花了三个晚上才跑出第一个可行解。问题不在算法而在热负荷曲线和燃料电池余热的对接只要把燃料电池的冷却热从热平衡里漏掉整个调度结果就会畸形。后来我意识到这个方向的关键不是追求更高级的求解器而是先把电、热、氢三条能量流之间的约束和参数设定搞清楚才谈得上优化调度。这篇文章面向想复现类似论文源程序的从业者和研究生把建模、代码、参数和踩坑一步步讲清楚。标题对应的论文可以在知网检索到但代码逻辑更多得靠自己去搭这正好是本文想帮你省下的时间。2. 建模先立住氢储能热电联供微电网里到底有哪些能量流和设备约束2.1 热电联供微电网的拓扑与四条能量母线先说清楚“氢储能热电联供型微电网”到底长什么样。常见拓扑是把光伏、风电、CHP机组、电解槽、储氢罐、燃料电池、储热罐、电负荷和热负荷接到同一个系统里。实际建模时我习惯按四条母线拆开看电母线、热母线、氢母线和气母线。电母线连接光伏、风电、CHP、燃料电池、电解槽、购电点和电负荷热母线连接CHP的余热、燃料电池的冷却热、储热罐和热负荷氢母线连接电解槽、储氢罐、燃料电池气母线连接外部天然气和CHP。这个拆法决定了优化调度模型的基本结构。四条母线通过设备存在耦合CHP烧气产电产热电解槽耗电产氢燃料电池耗氢产电同时产热。只写电平衡而忽略热平衡模型根本不能叫热电联供。很多复现失败的代码问题就出在少了热母线或者把热负荷当成固定值结果电出力算得很大热出力却对不上现场。在这个拓扑里氢储能更像一个跨母线的缓冲器光伏大发或电价低时电解槽制氢把电能转成化学能峰时燃料电池放电同时回收余热给热母线。这套逻辑不需要在模型里刻意指定运行模式只要把设备参数和能量平衡写对目标函数自然会引导求解器选到合理的运行区间。固定规则做不到这种跨时段套利这也是为什么这个方向一定要走优化调度。2.2 氢储能单元的关键参数制氢效率、燃料电池余热和储氢容量氢储能单元由三部分组成电解槽、储氢罐、燃料电池。在优化调度模型里这三个部件的参数直接决定解的质量和可行性。电解槽的核心参数是电转氢效率通常按氢的低热值折算数值在60%到75%之间。实际运行中效率会随负载率变化但论文里的源程序为了求解方便普遍用恒定效率处理。我复现时也先用常数算通以后再考虑变效率。储氢罐的关键参数有三个容量上限、初始储氢量、最低储量。建模时有两种思路一种直接用氢气的能量单位kWh另一种用压力或质量。用kWh最省事因为可以直接和电解槽、燃料电池的效率做乘除。如果论文给的是压力单位需要先完成热值换算第5章会专门讲这个坑。燃料电池的参数除了发电效率还有一个常被忽略的余热系数。燃料电池发电效率一般在45%到55%其余能量大部分变成热这部分热参与热电联供后系统总能效能明显提上去。模型里我一般加一个等式燃料电池余热功率等于发电功率乘余热回收系数系数取0.8到1.0。如果设成0热负荷缺口会被CHP强行顶上电出力虚高成本自然难看。储热罐参数也要先定好容量、每小时热损失率、充放热功率上限。热损失率一般取1%到3%每小时不要设成0。设成0会让储热罐变成无损耗的理想电池调度结果会过度依赖储热和现实差距很大。下表是我常用的参数起点具体数值按论文和实际设备调整。参数符号单位典型范围电解槽效率eta_elz%60-75燃料电池发电效率eta_fc%45-55燃料电池余热系数k_fc_heat无量纲0.8-1.0CHP发电效率eta_chp%35-45CHP热电比k_chp无量纲1.0-1.5储热罐热损失率heat_loss%/h1-32.3 为什么说固定规则调度在分时电价下一定不是最优的有人觉得微电网规模小直接设定优先顺序也能用有光伏先光伏不够就CHP再不够就买电。但一旦加入氢储能和储热罐固定规则的组合空间就爆炸了。同一个时段电解槽可以开也可以不开燃料电池可以放电也可以等峰时再放CHP可以多发电多产热也可以为了少出热而降出力。这些选择互相耦合固定规则只能拍到一种策略而且往往不是成本最低的。最典型的例子是夜间谷电。低价时段买电开电解槽制氢白天峰时燃料电池放电并回收余热这套操作在固定规则里不会出现因为固定规则只会在负荷不够时启动发电不会主动囤氢。但优化调度会算这笔账夜间买电成本低制氢存储后再放电虽然有损耗但峰谷价差如果超过损耗就是划算的。有没有这种策略取决于分时电价差和系统效率不是一句“氢储能不经济”能盖过去的。所以微电网优化调度不是单纯“哪个便宜开哪个”而是在时间维度上把各种储能设备的充放策略一起找出来。这套东西用数学语言描述就是下一章的目标函数和约束条件。读懂这一层比拿到一段能跑的代码更重要。3. 优化调度数学框架目标函数、三层平衡约束和求解器选型3.1 目标函数购电成本、天然气成本和启停成本怎么落进一个式子优化调度的目标函数绝大多数是最小化运行总成本典型构成是购电成本加天然气成本加启停成本。购电成本等于每个时段购电功率乘对应电价再乘时间步长。天然气成本主要是CHP消耗的燃料CHP发电功率除以发电效率得到燃料输入功率再乘天然气单价。燃料电池在这个模型里默认烧电解槽产的氢所以氢气成本通过氢平衡约束和购电成本间接体现不需要在目标函数里额外写一条“氢气单价”。启停成本在最小可运行模型里可以先不加因为加启停就要引入二进制变量模型从线性规划变成混合整数线性规划。但论文里的模型基本都会加因为不加剧离实际仿真里设备可以随便开停现场设备启停一次有寿命损耗和维护成本。我一般会在目标函数里加“启停罚项”用相邻两个时段的设备启停状态变量之差来判断是否有启停动作。目标函数的单位一致性常被忽略。购电成本用元天然气成本也用元如果天然气单价用了元/立方米必须先把立方米换算成千瓦时。天然气低热值大约每标准立方米9.5到10.5千瓦时具体以本地气源参数为准。量纲错了求解器算出来的调度策略一定偏而且很难排查。3.2 约束条件电、热、氢三条平衡约束的最小集合最小可行模型我把约束分成四组设备上下限约束、能量平衡约束、储能动态约束、初始和终端约束。设备上下限就是每个设备出力的最高最低限制。能量平衡是核心电平衡写“购电加CHP加燃料电池等于电负荷加电解槽”热平衡写“CHP余热加燃料电池余热加储热罐净充放等于热负荷”氢平衡写“储氢量变化等于电解槽制氢减燃料电池耗氢”。储能动态约束把上一时段的储量、本时段的充放和损耗串起来。储热罐必须加每小时损失率否则模型会把热量存上十几个小时再放现实中热量早散光了。储氢如果假设无泄漏可以不加强制损失项但容量上下限一定要有。约束的表达式写成文本是这样的# 电平衡功率单位统一为 kW pbuy[t] pchp[t] pfc[t] pload[t] pelz[t] # 热平衡qtes 为正表示向储热罐充热为负表示放热 k_chp * pchp[t] k_fc_heat * pfc[t] - qtes[t] qload[t] # 氢平衡储氢量单位 kWh het[t] het[t-1] eta_elz * pelz[t] * dt - (pfc[t] / eta_fc) * dt # 储热动态 tes[t] tes[t-1] - qtes[t] * dt - heat_loss * tes[t-1] * dt这里有一个容易懵的地方qtes 的正负。热源富余时 qtes 为正多余的热存进储热罐热源不足时 qtes 为负储热罐放热补缺口。使用净变量可以省掉“充放不能同时进行”的二进制约束代价是充放效率变成对称的。第4章的代码沿用了这个简化复制到工程项目时再考虑用两个变量。初始约束把第一个时段的储能初值固定终端约束要求调度结束时储能不低于初始值。很多复现代码只写初始不写终端结果最后两个时段储能被放空成本光鲜但第二天没法接着跑。这一点第5章还要展开。3.3 求解器选型与问题规模Gurobi、CBC、SCIP 怎么选模型写好后求解器选择直接决定能不能在可接受时间内拿到解。小规模调度比如一台CHP、一台电解槽、一台燃料电池、几十个时段用开源的CBC就行装个pip包就能跑。一旦加二进制变量CBC求解速度会明显下降几分钟甚至半小时都出不来。这时可以换商业求解器比如Gurobi或CPLEX。我的习惯是先用免费求解器把模型逻辑跑通验证约束没错再换商业求解器做敏感性分析。模型有bug时用再贵的求解器也是白搭。如果目标是复现论文数据最好看论文里写了什么求解器尽量保持一致。不同求解器对MIP的容差和割平面策略不一样同一组约束可能得到略有差异的解。论文没写就用开源的结果趋势对不用纠结完全复现成一样的数字。求解器对问题规模的敏感点主要在二进制变量数量和Big-M约束质量。Big-M的M值取得太大数值稳定性会崩取太小又会误伤可行域。常见做法是把M设成该设备功率上限的1.2到1.5倍不要拍脑袋给个十万八这是新手最容易踩的隐性坑。4. 用 Python-mip 复现最小调度模型完整代码与参数解读4.1 最小可运行代码在本地跑通 24 小时调度这一章直接给一份能跑的Python代码用开源的python-mip库实现氢储能热电联供微电网的24小时优化调度。它把上一章的约束全部落进去不依赖商业求解器和MATLAB。先安装依赖pip install mip。# 最小可运行氢储能热电联供微电网 24 小时调度 # 依赖pip install mip from mip import Model, xsum, minimize, CONTINUOUS T 24 # 时段数步长 1 小时 DT 1.0 # 时间步长h # 输入分时电价元/kWh23点后到次日8点为谷段 c_buy [0.32 if (t 8 or t 23) else 0.72 for t in range(T)] c_gas 0.40 # 天然气折价元/kWh 天然气 # 输入电负荷、热负荷kW用典型日形状代替实测数据 pload [180 90 * (1 if 8 t 20 else 0) for t in range(T)] qload [220 40 * (1 if (6 t 10 or 16 t 20) else 0) for t in range(T)] # 设备参数 eta_elz 0.7 # 电解槽电转氢效率 eta_fc 0.5 # 燃料电池氢转电效率 k_fc_heat 0.8 # 燃料电池余热系数 eta_chp 0.4 # CHP 发电效率 k_chp 1.2 # CHP 热电比 heat_loss 0.01 # 储热罐每小时热损失率 h_init 300.0 # 初始储氢量kWh tes_init 200.0 # 初始储热量kWh h_cap 1000.0 # 储氢罐容量上限kWh tes_cap 500.0 # 储热罐容量上限kWh m Model(hydrogen_cchp_scheduling) # 变量每个时段的设备功率和储能状态 pbuy [m.add_var(var_typeCONTINUOUS, lb0, ub1000, namefpbuy_{t}) for t in range(T)] pchp [m.add_var(var_typeCONTINUOUS, lb0, ub200, namefpchp_{t}) for t in range(T)] pfc [m.add_var(var_typeCONTINUOUS, lb0, ub100, namefpfc_{t}) for t in range(T)] pelz [m.add_var(var_typeCONTINUOUS, lb0, ub100, namefpelz_{t}) for t in range(T)] qtes [m.add_var(var_typeCONTINUOUS, lb-200, ub200, namefqtes_{t}) for t in range(T)] het [m.add_var(var_typeCONTINUOUS, lb0, ubh_cap, namefhet_{t}) for t in range(T)] tes [m.add_var(var_typeCONTINUOUS, lb0, ubtes_cap, nameftes_{t}) for t in range(T)] # 约束电、热、氢、储热 for t in range(T): # 电平衡购电 CHP 燃料电池 电负荷 电解槽 m pbuy[t] pchp[t] pfc[t] - pelz[t] pload[t] # 热平衡qtes 为正表示向储热罐充热为负表示放热 m k_chp * pchp[t] k_fc_heat * pfc[t] - qtes[t] qload[t] # 氢平衡储氢量 上时段储氢量 制氢 - 耗氢 h_prev h_init if t 0 else het[t-1] m het[t] h_prev eta_elz * pelz[t] * DT - (pfc[t] / eta_fc) * DT # 储热动态损失与当前储热量成正比 tes_prev tes_init if t 0 else tes[t-1] m tes[t] tes_prev - qtes[t] * DT - heat_loss * tes_prev * DT # 终端约束调度结束时储能不低于初始值保证方案可持续 m het[T-1] h_init m tes[T-1] tes_init # 目标购电成本 CHP 天然气成本 m.objective minimize( xsum((c_buy[t] * pbuy[t] c_gas * (pchp[t] / eta_chp)) * DT for t in range(T)) ) status m.optimize() print(solver status:, status.name) if m.num_solutions: print(total cost:, round(m.objective_value, 2)) for t in range(T): print(t, buy, round(pbuy[t].x, 1), chp, round(pchp[t].x, 1), fc, round(pfc[t].x, 1), elz, round(pelz[t].x, 1), het, round(het[t].x, 1), tes, round(tes[t].x, 1))代码逻辑上每个时段的电平衡把购电、CHP、燃料电池和电解槽放在一个等式里这等于强制系统每一小时发电和用电必须相等。热平衡里的 qtes 是一个净变量正负号决定储热罐是充电还是放电省去了充放同时存在的非线性判断。氢平衡用上一时段的储氢量累加把电解槽产氢和燃料电池耗氢绑定。t0 时上一时段用的是初始储量这是储能模型的固定写法。参数说明c_buy 设置的是峰谷电价0.32元和0.72元分别对应谷段和峰段实际项目直接替换成当地分时电价表即可。pload 和 qload 用了一个简单的日负荷形状白天电负荷抬升、早晚热负荷抬升你可以换成自己工程里的实测曲线。设备上下限在 add_var 的 lb 和 ub 里约束这里把CHP上限设200kW、燃料电池100kW、电解槽100kW是因为示例系统容量就按这个规模配。跑通以后把 ub 改成你的真实设备铭牌参数模型不用动。4.2 从最小模型到论文源程序还需要补哪些模块第4.1节的代码是一个可复现的骨架。论文里的源程序通常会多几个模块第一个是启停状态变量。给每个可启停设备加一个二进制变量z[t]表示该时段是否运行然后让出力上下限乘以z[t]目标函数里加启停成本。这会增加几十个二进制变量但模型更接近现实。第二个是爬坡约束相邻时段出力变化量限制在设备爬坡率以内。CHP和电解槽都有爬坡限制不做这一步调度结果可能在现场执行不了。第三个是储氢罐的SOC曲线。论文如果研究压力变化会把储氢量换算成压力再加一个压力和功率上下限的关系如果只是经济调度用kWh能量量纲就够了。第四个是热回收的上限约束。燃料电池余热不是无上限的除了固定比例还应加一个最大热回收功率防止模型在热负荷大时让燃料电池超发。这些模块我在工程里会逐步加上每次只加一类改完就跑一遍对比成本而不是一次性全塞进去。4.3 MATLABYalmip 写法与 Python 写法的对应关系论文源程序用MATLABYalmip写很常见如果你习惯MATLAB核心框架和Python几乎一一对应。变量声明从sdpvar开始约束用Constraints [Constraints, ...]拼接求解调用optimize。下面是氢平衡这一段在YALMIP里的对照写法%% YALMIP 与 Python-mip 的核心对应关系片段 T 24; dt 1; H0 300; % 初始储氢量 Pelz sdpvar(T,1); % 电解槽功率 Pfc sdpvar(T,1); % 燃料电池功率 H sdpvar(T,1); % 储氢量状态 Constraints []; % t 1 使用初始值 Constraints [Constraints, H(1) H0 eta_elz * Pelz(1) * dt - Pfc(1)/eta_fc * dt]; % t 2..T 使用上一时刻状态 for t 2:T Constraints [Constraints, H(t) H(t-1) eta_elz * Pelz(t) * dt - Pfc(t)/eta_fc * dt]; end % 其余平衡约束同理最后调用 % optimize(Constraints, objective, sdpsettings(solver,gurobi,verbose,0));对应关系很直观python-mip 用m ...加约束YALMIP 用Constraints [Constraints, ...]拼接变量下标从0开始还是从1开始是两个工具最明显的差异。这个片段不是一份完整源码而是帮你把第4.1节的模型翻译到MATLAB环境的对照。真正落地时把电平衡、热平衡、储热动态和终端约束都按同样方式补上即可。5. 氢储能调度复现的 4 个典型坑现象、原因与解决办法5.1 储氢罐单位换算kWh、kg、标方混用直接导致无解现象模型看起来没问题但求解器报不可行或者氢平衡约束总对不上。检查约束时发现储氢量有时候是kg有时候是kWh甚至状态变量初始值和约束表达式里用的是两套单位。原因论文数据来源不统一。氢气能量密度按低热值算约33.3 kWh/kg按标准状况下的体积算约3 kWh/Nm³。有人把kg直接当kWh用一个数量级的偏差就足以让约束矛盾。解决建模一开始就规定所有氢相关变量统一用kWhkg和Nm³只在参数输入时换算。我给每条氢相关的约束都加注释写明“此变量单位是kWh不是kg”。代码里加一行强制检查初始化时把换算后的数值范围和物理常识对比超过100倍直接报错。提示换算储氢量时低压储氢罐的密度随压力变化很大经济调度阶段不必模拟压力动态直接用能量量纲更稳。5.2 燃料电池余热不参与热平衡热负荷越大成本越离谱现象调度结果里热负荷峰值时段CHP出力很大电出力跟着超标甚至需要向电网卖电整体成本高得离谱。对比论文里的运行曲线发现燃料电池原本在那些时段有出力。原因复现代码只把CHP的余热打进热母线燃料电池只建模成发电设备k_fc_heat这个参数被忽略或设成了0。热平衡缺了这部分热源系统只能让CHP多发电来满足供热电侧被迫多出力。解决在热平衡约束里加一项k_fc_heat * pfc[t]并把余热系数设为0.8到1.0。同时给燃料电池余热功率单独设一个上限值比如qfc_heat_max 80 kW因为受换热器面积限制余热回收不可能无限大。加上这项以后热负荷缺口明显减少整个调度成本也会降下来。5.3 终端储能约束缺失最后两小时储能被免费放空现象24小时调度结果很漂亮总成本很低。但把这24小时结果作为第二天初始状态输入模型时第二天直接不可行或者开机后储氢罐见底。原因模型没有写终端约束优化器会认为调度结束后储能没有价值于是把储氢罐和储热罐在最后几个时段全部放空用来省成本。这在单时段看起来是“最优”但破坏了日与日之间的连续性。解决加两个终端约束het[T-1] h_init和tes[T-1] tes_init。如果希望更严格可改成等于初始值。我一般还会把终端约束写成储能状态直接关联比如“最后4小时储能量不低于初始值的80%”给调度留一点富余量防止数值误差造成第二天初值越界。5.4 量纲与数值缩放MW 和 kW 混用让求解器报 Numeric issues现象模型在几台设备参数规模不同时求解器报数值问题或者收敛很慢。同样的模型把参数统一改写一遍后求解速度明显提升。原因有的数据从论文里直接抄个别功率用了MW负荷曲线用kW成本用万元电价用元每千瓦时。约束矩阵里的系数跨了三个数量级对求解器的数值稳定性是很不友好的。解决所有功率统一到kW电价统一到元/kWh成本结果统一到元。设备上下限都写成显式常数不要在约束里做隐式换算。我通常在数据加载后加一行断言检查所有数组的最大值不超同一量纲范围比如功率数据应该在0到1e4之间如果出现1e6就说明单位混了。这一步能省下大量排查时间。6. 把日调度升级成在线调度滚动时域、边界验证和一点习惯6.1 滚动时域调度把一次性离线优化改成四小时一轮日调度在实际运行中很难直接用因为负荷和电价都有预测误差。常见做法是把24小时窗口按滚动方式执行每次求一个24小时的最优解但只执行前4小时的指令4小时后再拿最新预测重新求解。这样既保留了全局优化的视野又能吸收新的预测信息。实现上只需要在代码外面套一层循环核心求解函数不需要大改。滚动窗口长度一般取24小时或更长的36小时执行步长取2到4小时。执行步长越长对预测误差的抵抗力越弱步长越短求解频率越高对算力要求也越高。我一般先做4小时一轮够用再尝试缩短。6.2 两个验证习惯影子价格检查和容量拉满测试调度结果出来以后不要急着布置现场。我做两个验证。第一个是影子价格检查。求解器会给每条约束返回一个对偶值表示该约束每放宽1个单位目标函数能降多少。如果某时段储氢容量约束的影子价格很高说明系统在那些时段特别需要储能增加储氢容量会直接带来收益。如果整个周期内所有储能容量的影子价格都是0说明储能设备在当前负荷和电价下没有边际价值白配了。第二个是容量拉满测试。把燃料电池上限从100kW改成1000kW重新求解。如果总成本下降几乎为零说明燃料电池在当前模型里没有发挥关键作用要么参数设置有问题要么氢储能参与调度的空间本来就不大。这个测试比看一堆运行曲线更能说明问题。我的习惯是每次改动参数都跑一轮对比把成本、储能终值、峰时购电三个数字记下来形成一个小表时间长了就能看出模型边界在哪。这套流程坚持下来复现或自建微电网优化调度会越来越顺手。很多问题看起来像求解器玄学实际都是单位、终端约束和余热参数在捣鬼。希望这篇笔记能帮你少加点班把精力留在真正值得研究的调度策略上。本文还有配套的精品资源点击获取