做微网优化运行这个方向的同学大概率绕不开Matlab。特别是“基于多能互补的热电联供型微网优化运行”这类题目可以说是目前能源系统领域里最典型、也最容易出成果的一类项目。很多硕博论文、期刊文章、竞赛作品都围绕这个主题展开一方面是因为它紧扣双碳背景下的综合能源系统趋势另一方面是因为它的物理模型、数学规划、代码实现都有比较成熟的套路可循。这篇文章我打算抛开教科书式的表述把自己实际建模和写码过程中踩过的坑、验证过有效的思路、以及Matlab实现里的关键细节整理出来。内容涉及整体方案设计、热电联供系统的能量流拆分、目标函数与约束条件的数学化、YalmipGurobi或Cplex求解流程、典型日算例做结果分析以及常见代码报错的排查思路。无论你是刚接触微网优化、正在复现某篇论文还是准备自己搭一套完整的调度程序这篇文章都能帮你少走弯路。1. 整体框架与运行逻辑拆解1.1 热电联供型微网的物理架构与能量流先理清楚我们在优化什么。典型的热电联供型微网CHP-Microgrid由分布式电源、储能设备、热电转换设备以及负荷构成。最常出现的设备包括微型燃气轮机MT、燃气锅炉GB、电锅炉EB、蓄电池ESS、蓄热罐TST、光伏PV、风电WT有时候还带制冷机组变成冷热电三联供CCHP。能量流是“多能互补”的核心电网买入的电、光伏发出的电、风电发出的电汇入母线供给电负荷微型燃气轮机燃烧天然气发电同时产生高温烟气通过余热回收装置转换成热能供向热母线热母线同时接受燃气锅炉产热、电锅炉产热、蓄热罐放热满足热负荷蓄电池和蓄热罐分别作为电、热两个时间维度上的缓冲装置起到“削峰填谷”的作用。实操中我习惯先把画一个能流拓扑图把设备节点和母线关系固定下来再写代码。这一步省不了因为后面所有的约束、变量、平衡方程都必须对应到这张图上。1.2 多能互补与热电联供在调度中的价值多能互补的本质是把原本各自独立的电力系统、热力系统、燃气系统在微网层面耦合在一起通过能源品种之间的替代和转换来提高整体效率、降低成本。举一个最简单的例子电价高的时段微型燃气轮机满发不仅供电还供热替代了燃气锅炉的供热份额提升了天然气的综合利用效率电价低的时段则多从电网购电同时让蓄热罐蓄热、蓄电池充电把低成本的能源存储下来留到高峰时段使用。电和热在时间上的“耦合解耦”就是热电联供微网优化运行最大的发挥空间。所以目标并不是单纯让某一台设备效率最高而是在整个调度周期内所有设备的出力组合使系统总运行成本最低。这里才会出现“电跟热走”还是“热跟电走”两种运行策略的选择代码实现中通常通过约束条件来体现。1.3 为什么选择Matlab作为实现工具虽然Python在数据分析领域势头很猛但微网优化运行这种“建模求解”任务里Matlab仍然是最顺手的选择。矩阵运算天然高效多时段、多设备的优化变量转成向量/矩阵构造约束很方便Yalmip工具箱把优化建模做得极其顺手写约束条件几乎和数学公式一一对应直接调用Gurobi、Cplex等商用求解器读参数、调间隙、查解状态都封装得很好绘图功能成熟结果可视化和论文出图一体化尤其是画负荷曲线、设备出力堆叠图Matlab一套代码全搞定。如果你手头有Matlab和Yalmip环境这篇文章里所有思路都可以直接落地如果没有装Yalmip也可以改用linprog线性规划或者intlinprog混合整数线性规划实现但建模体验会差不少后期改约束会很痛苦。2. 数学模型构建与关键约束梳理2.1 目标函数经济运行成本最小化绝大多数热电联供微网优化运行的论文目标函数都写成“日运行总成本最小”组成项一般包括以下几部分从上级电网购电费用C_grid sum(price_e(t) * P_buy(t))分时电价是关键天然气燃料费用包括给微型燃气轮机和燃气锅炉供气的费用设备运行维护费用一般按各设备的输出功率乘以一个运维成本系数来估算弃风弃光惩罚费用有新能源时建议加碳排放成本有些论文会加看题目要求。我给出的建议是第一版代码不要把所有成本项都塞进去先做“购电费用燃气费用运维费用”三项主体部分求解跑通了再逐步增加碳排放约束、需求响应或不确定性处理。这样做的好处是每一步都能定位到问题否则模型一出错根本不知道是目标函数写错了还是约束条件起冲突了。目标函数写出来是离散时间求和的形式优化变量为每个时段的各设备出力值总调度周期T通常取24小时步长1小时。2.2 约束条件设备运行边界与功率平衡约束条件是最容易翻车的地方。我把实际建模中最常遇到的约束梳理一遍电功率平衡约束电网购电功率光伏出力风电出力燃气轮机发电功率蓄电池放电功率 电负荷电锅炉耗电功率蓄电池充电功率。注意电锅炉是耗电设备必须在等式右侧热功率平衡约束燃气轮机余热回收功率燃气锅炉产热功率电锅炉产热功率蓄热罐放热功率 热负荷蓄热罐吸热功率设备出力上下限约束每台设备都有最大、最小出力爬坡约束燃气轮机和燃气锅炉由于物理惯性相邻时段出力变化量有上限蓄电池SOC约束荷电状态在允许范围内同时考虑充放电效率不能同时充放电蓄热罐容量约束相当于热力侧的“储能电池”同样有蓄热/放热功率限制和容量限制。对于初次接触这些约束的读者不要被这么多条件吓到。在Yalmip里它们基本就是一行一行的Constraints [Constraints, ...]。把每条约束当成一张“电网物理规则清单”逐条翻译出错了也能很好地追踪。2.3 不确定性处理从确定性模型到鲁棒/随机优化基础版本做的是确定性优化——光伏出力、风电出力、负荷大小都是固定数值。但真实的微网运行中新能源出力和负荷天然带随机性。论文里常见的处理手段有三种场景法随机优化假设光伏、风电、负荷服从某种概率分布抽样生成大量场景将单场景扩展为多场景期望成本最小化鲁棒优化用不确定区间描述新能源出力求解最坏情况下的最优方案模型预测控制MPC把日前计划与日内滚动修正结合不断用最新测量数据更新调度指令。在Matlab里场景法和MPC都比较好实现鲁棒优化则需要引入对偶变换对数学功底要求高一点。我建议先掌握确定性版本再慢慢扩展。3. 求解方法与Matlab代码实现要点3.1 求解器选型YalmipGurobi是黄金组合关于求解器直接说结论能用Gurobi就用Gurobi没有许可license的用Cplex也行再没有就用Matlab自带的intlinprog注意是intlinprog不是linprog因为有蓄电池充放电状态这种0-1变量时模型是混合整数规划。Gurobi在学术圈的普及度非常高很多学校都有免费学术license申请渠道。在Matlab里只需要% 添加路径 addpath(genpath(D:\yalmip)); addpath(genpath(D:\gurobi)); % 定义优化变量 P_mt sdpvar(1, 24); % 燃气轮机发电功率 P_gb sdpvar(1, 24); % 燃气锅炉产热功率 P_eb sdpvar(1, 24); % 电锅炉耗电功率 P_buy sdpvar(1, 24); % 购电功率 P_ch sdpvar(1, 24); % 蓄电池充电功率 P_dis sdpvar(1, 24); % 蓄电池放电功率 u_ch binvar(1, 24); % 充电状态0-1变量 u_dis binvar(1, 24); % 放电状态0-1变量 SOC sdpvar(1, 24); % 蓄电池荷电状态Yalmip最舒服的一点是连续变量用sdpvar0-1整数变量用binvar优化变量定义清楚后后面不管写什么都非常直观。这里强调一点点充电和放电状态建议用两个独立的0-1变量来表示再加互斥约束避免求解器出现同时充放电这种不合理的解。3.2 初始化与sdpvar变量定义规范写代码第一步不是直接敲约束而是先“拆时段、拆设备、拆变量”。我会先定义三个列表时段列表T 24小时制设备列表燃气轮机、燃气锅炉、电锅炉、蓄电池、蓄热罐、光伏、风电每个设备对应一组变量名参数列表设备容量、效率、爬坡速率、储能容量、分时电价、气价。变量全部用向量定义不要写成24个标量。一个sdpvar(1,24)一次搞定后面约束和计算都靠矩阵、向量运算速度提升几个数量级代码也短很多。这里给一个初学者非常容易犯的错变量维度和时段数量不一致。比如你定义了P_mt sdpvar(1, 24)后面约束却写循环for t 1:24下P_mt(t)这没问题。但有些人会把优化变量定义成sdpvar(24, 1)即列向量后面索引时P_mt(t)还是对的但写矩阵乘积时就容易维度错乱。我的习惯是统一用行向量所有约束都默认1x24维度保持一致性。3.3 约束条件写入与求解指令约束写入非常流水线化下面是我实际项目中的代码片段摘取了核心部分保证可直接理解Constraints []; % 电功率平衡约束每个时段都满足 % 购电光伏风电燃气轮机发电蓄电池放电 电负荷电锅炉耗电蓄电池充电 Constraints [Constraints, ... P_buy P_pv P_wt P_mt P_dis P_load P_eb P_ch]; % 热功率平衡约束 % 燃气轮机余热燃气锅炉产热蓄热罐放热 热负荷电锅炉产热蓄热罐吸热 Constraints [Constraints, ... eta_hr * P_mt P_gb P_tst_out P_hload eta_eb * P_eb P_tst_in]; % 燃气轮机出力上下限 Constraints [Constraints, ... 0 P_mt P_mt_max]; % 燃气轮机爬坡约束 Constraints [Constraints, ... -delta_mt diff(P_mt) delta_mt]; % 蓄电池充放电互斥 Constraints [Constraints, ... P_ch u_ch * P_ch_max, ... P_dis u_dis * P_dis_max, ... u_ch u_dis 1]; % SOC递推约束初始SOC给定 Constraints [Constraints, ... SOC(1) SOC_init, ... SOC(2:end) SOC(1:end-1) (eta_ch * P_ch(1:end-1) - P_dis(1:end-1) / eta_dis) / Cap_ess, ... SOC_min SOC SOC_max];求解部分则非常简单ops sdpsettings(solver, gurobi, verbose, 2); optimize(Constraints, objective, ops);求解完成后直接从value()提取各个变量的结果即可P_mt_opt value(P_mt); P_buy_opt value(P_buy); SOC_opt value(SOC);如果求解器报告infeasible不可行不要慌先检查约束是否冲突尤其是SOC递推约束初值和容量限制、以及储能状态与系统平衡约束之间是否矛盾。后面章节我会专门列常见问题。3.4 智能优化算法在微网问题中的补充用法额外提一句“改进的麻雀搜索算法优化长短期记忆神经网络MISSA-LSTM”这类的热搜词这几年一直很火。在热电联供微网这个项目里神经网络和启发式算法并不是用来直接求解调度模型的而是出现在两处用LSTM做光伏出力或负荷的预测预测结果作为确定性模型的输入参数用麻雀搜索算法SSA的改进版优化LSTM的初始学习率、隐含层节点数、正则化系数等超参数提升预测精度。如果在Matlab里做负荷预测流程一般是构造输入特征序列划分训练集测试集把LSTM网络定义好然后用MISSA在超参数空间搜索一组最优配置再重训LSTM完成预测。预测精度提高了输入到优化模型的新能源/负荷曲线就更可靠整个调度结果的可信度自然也就上去了。这里提醒一句LSTM预测和微网优化是两个相对独立的模块。不要指望用SSA直接去搜微网调度方案的0-1组合——启发式算法求解这种高维混合整数问题很容易陷入局部最优也不稳定。正确姿势还是预测模块用智能算法LSTM调度模块用混合整数规划商用求解器。3.5 绘图与结果展示细节结果可视化的基本功是出“设备出力堆叠图”和“电/热平衡曲线”。有的论文需要把光伏出力的数据点标记为“空心黑圆点”在Matlab里实现非常简单plot(time, P_pv_opt, ko, MarkerFaceColor, none, MarkerEdgeColor, k, MarkerSize, 6);其中的ko表示黑色圆圈MarkerFaceColor, none就是让圆点中间空心MarkerEdgeColor, k设置边框为黑色。有时候画微网能量流图也需要这种空心点的标记用来标注光伏、风电这类间歇性电源的出力比实心点看起来清爽很多尤其适合彩色打印的论文插图。另外我推荐把电平衡的各个分量堆叠成面积图area再叠加一条总负荷曲线这样最直观。类似的热平衡单独一张图蓄热罐的充放热状态单独画不要挤在一张图里。4. 典型算例设计与运行结果分析4.1 算例数据构造在代码跑起来之前需要先准备24小时的输入数据。典型场景如下分时电价峰时段10:00-15:0018:00-21:001.2元/kWh平时段7:00-10:0015:00-18:0021:00-23:000.8元/kWh谷时段23:00-次日7:00还有午间11:00-13:00视地区政策可设定0.4元/kWh天然气价格2.4元/m³天然气热值按9.7kWh/m³折算微型燃气轮机额定功率500kW电效率35%热电比1.2燃气锅炉额定功率400kW热效率90%电锅炉额定功率200kW热效率95%蓄电池容量1000kWh最大充放电功率200kW充电效率0.95放电效率0.95SOC范围20%-90%光伏出力曲线和电热负荷曲线用典型日数据代入。这里我们可以简单验算一下微型燃气轮机的效率表达式P_mt eta_mt * F_gas * LHV也就是燃气轮机发电功率等于天然气消耗量乘以热值再乘发电效率。从这一项同时推导出燃气费用和余热回收量所以代码中P_mt既出现在目标函数燃气成本里也出现在电平衡和热平衡中是一个“枢纽型”变量。4.2 运行结果的经济性分析跑完优化后把24小时的储能SOC、购电功率、燃气轮机出力、燃气锅炉出力画在一起观察通常会看到几个典型现象深夜谷时段蓄电池从电网大量购电充电热负荷由蓄热罐放热和燃气锅炉供给微型燃气轮机停机或低载上午峰时段微型燃气轮机启动电功率爬坡热功率同步上升余热优先满足热负荷不足部分由燃气锅炉补充傍晚高峰蓄电池放电支撑电负荷微型燃气轮机继续运行蓄热罐此时通常已经放空。从成本构成来看购电费用大致占40%-55%燃气费用占35%-45%运维费用占剩余部分。加入蓄热罐之后热系统的运行灵活性会明显提升燃气锅炉的启停次数减少、平均负载率提高热力侧的运行成本通常能下降5%-15%。这也是论文里最有价值的对比结论之一。4.3 灵敏度分析与方案对比为了证明模型有效通常还需要做几组对比方案1不含储能、不含蓄热罐纯燃气轮机锅炉电网直接供电供热方案2只加蓄电池不加蓄热罐方案3只加蓄热罐不加蓄电池方案4完整的热电联供微网含蓄电池蓄热罐。每一组跑一次优化把总成本、购电费用、燃气费用列成一张表。通常方案4的总成本最低因为储能装置跨时段搬移了能量蓄热罐利用了热力系统的柔性两者叠加的效果大于单独运行之和。我实测的感受是蓄热罐对热负荷波动大尤其北方地区冬季供暖场景的项目增益更明显蓄电池则对分时电价差大的场景增益更明显。如果你的算例里热负荷比较平稳、电价波动也不大那储能带来的成本下降空间有限这会直接影响你对模型价值的判断。写论文时一定要根据自身数据来设计对比方案不要上来就主观断定“储能越多越好”。5. 常见问题与排查技巧实录5.1 求解报错“Infeasible problem”怎么办这是新手最容易卡住的地方。模型不可行通常的原因有几个SOC递推约束写得不对比如初始SOC设0.5容量1000kWh储能约束又要求每时段SOC不低于0.2如果某时段放电过大递推值会跌破下限导致整体不可行功率平衡约束缺项电平衡中忘了加电锅炉的耗电项或热平衡中忘了加蓄热罐的吸热项设备上下限约束与平衡约束矛盾比如燃气轮机最小出力设为600kW但额定功率只有500kW这种低级错误也见过不少0-1变量与连续变量关系不匹配充放电互斥约束中如果P_ch上限与状态变量u_ch的逻辑关系错了比如充电功率上限写成0.2而不是200就会导致模型过于保守甚至无解。排查方法很简单把约束一条一条注释掉看哪条约束去掉后问题变可行基本就能锁定冲突位置。我用这个方法排查过无数次效率非常高。5.2 求解器返回“Numerical issues”或收敛慢如果模型带了大数量级差异的参数比如储能容量1000kWh、气价2.4元/m³、电功率因子1e5Gurobi在数值上会很难受。缓解方式包括统一单位功率用kW、能量用kWh、价格用元/kWh避免出现1e6数量级的中间量给决策变量设置合理的边界不要用[-Inf, Inf]的默认边界用sdpsettings(gurobi, MIPGap, 0.01)控制相对间隙不要盲目追求0绝对最优如果模型是纯线性可以设置gurobi, Presolve, 1让预求解把模型压缩一下。5.3 优化结果中蓄电池出现同时充放电这个问题非常经典。即使你写了互斥约束仍然可能出现“充电功率和放电功率都为正”的情况原因是充放电效率不对称目标函数借助充放电循环套利。处理办法要么像前面那样增加互斥0-1变量要么把充放电效率统一成相同的数值。如果你确认用了互斥约束还是出现同值现象检查一下辅助变量是否被约束准确了P_ch应当同时满足P_ch u_ch*P_ch_max和P_ch 0两个缺一不可。5.4 常见问题速查表现象可能原因处理建议模型不可行SOC递推与容量约束冲突逐步注释约束定位冲突来源求解时间过长整数变量过多或MIPGap过紧增加MIPGap容差减少储能分段线性化段数蓄电池同时充放互斥约束缺失或效率不对称增加0-1互斥变量统一效率热平衡不满足余热回收系数乘错仔细核对热电比与换热效率的定义结果中燃气轮机全天不出力气价太高或热电比设置过低检查气价和电价的相对关系调整算例参数画图时用了空心圆但显示为实心MarkerFaceColor未设置加MarkerFaceColor, none5.5 代码调试的5条实战心得先跑24小时单时段简化模型确保程序能跑通再扩展到完整时段每次只加一个设备/约束出问题能快速回溯所有参数集中放到脚本顶部不要散落在代码中间否则调参改到崩溃用value()查看结果后先人工验算一下某个时段是否满足功率平衡不要直接相信求解器输出把每次跑出来的最优成本、求解时间、间隙记录下来形成日志方便前后对比。6. 从复现到创新的扩展方向代码跑通、结果合理之后如果想往更高水平走可以考虑几个方向。第一个是引入需求响应。把电负荷和热负荷从固定值变成部分可调变量设置可平移、可削减负荷的比例目标函数里加入补偿成本项本质是在用户舒适度和系统经济性之间做权衡。第二个是加入碳交易机制。给碳排放配额、碳价建模目标函数变成“运行成本碳交易成本”这样高碳排放的燃气设备出力会受到约束低碳的电气化供热设备会更受青睐。第三个是设备退化建模。蓄电池和蓄热罐的容量衰减、健康状态SOH变化可以提升模型的长期可靠性但代码复杂度会上升不少。第四个是用MPC框架做日内滚动优化把日前计划、日内修正、实时反馈串联起来让调度从“开环”变成“闭环”更贴近实际工程需求。我个人在实际操作中的体会是优化运行项目没有太多玄学核心就是把物理过程量化准确、把求解问题建模清楚、把结果分析做透。代码实现上先在确定性模型上站稳脚跟再逐步加入不确定性和多时间尺度这才是最稳的路线。如果做这个项目是为了论文发文建议对比实验多做几组、数据表格做得美观灵敏度分析一定不要省——审稿人很喜欢看这个。