☰
微能源网优化调度:冷热电联供系统的MILP建模与求解
2026/10/3 9:31:29 网站建设 项目流程

先说结论:冷热电气多能互补的微能源网优化调度,本质上是一个混合整数线性规划(MILP)问题,而且是一个特别适合用Matlab+YALMIP+Gurobi组合求解的问题。你可能看过不少论文用粒子群、遗传算法做这个,但我在实际对比之后,还是坚定地选了MILP这条路——原因后面展开讲。

这篇内容适合正在做综合能源系统、微电网、园区级冷热电联供系统研究的同学,也适合刚接触优化调度、想快速搭出一个能跑的Matlab模型的工程师。我会把设备建模、耦合约束、代码骨架、典型工况结果分析、还有我实际调试中踩过的坑都摊开讲,尽量做到看完就能照着搭一版。

1. 为什么我最终选择了YALMIP+Gurobi这条技术路线

1.1 先说我走过的那条弯路:智能算法不是万能的

一开始我也跟风用过粒子群优化(PSO)和遗传算法(GA)来做冷热电联供系统的经济调度。当时觉得"智能算法不需要求导、不需要凸性假设,多灵活"。实际跑起来才发现问题相当大。

调度问题里大量存在设备启停、储能充放这样的离散决策——买不买电、启不启机、充不充能,这些都是0-1变量。连续型的智能算法处理这类变量特别别扭,最常见的做法是把0-1变量放松成连续变量再阈值化,这一放松,可行域就变了,解得的最优性根本站不住脚。约束一多更是灾难,我最多试过28个设备、96个时段,PSO的罚函数法调惩罚系数调到怀疑人生,罚轻了不可行,罚重了目标被带偏。

还有个更现实的问题:智能算法每次运行结果都不一样。写论文的时候,审稿人一句"请提供全局最优解证明"直接卡壳。不是说智能算法没用,它在非线性强、模型不精确、可行域模糊的场合确实有优势。但微能源网调度模型本质上是线性的,完全没有必要用启发式算法去碰运气。

1.2 MILP建模的优势:为什么它是微能源网调度的主流范式

微能源网里最常见的设备——燃气轮机、燃气锅炉、电制冷机、吸收式制冷机、电锅炉、储电/储热/储冷——在稳态运行下,能耗关系基本都能用线性或分段线性的效率模型表达。这意味着目标函数(购电购气成本、设备运行成本、碳排放成本)、能量平衡约束(电热冷气四条母线)、设备运行区间和爬坡约束、储能容量约束,统统可以写成线性数学表达式;设备启停决策由二进制变量表达。整个模型毫无悬念地落入标准MILP框架。

对MILP问题,Gurobi、CPLEX这类商业求解器可以在给定精度下证明最优性(最优间隙可以收敛到0),这是启发式算法做不到的。对科研写论文、做工程方案比选,这个"可复现、可验证最优"的属性太重要了。而且MILP求解器内部的branch and cut、presolve这些技术,已经成熟到可以处理几千个变量的模型,求解速度远快于你手写任何启发式算法。

1.3 为什么用Matlab而不是Python

说实话,用Python+Gurobi也能做,我身边不少同事就是这么干的。我之所以继续用Matlab,主要是三个理由。

一是数据前后处理和可视化方便。Matlab的plot家族画设备出力堆叠图、负荷曲线对比图,代码量最小,尤其area函数画堆叠面积图,几乎是为调度结果展示量身定做的。二是YALMIP这个建模工具箱在Matlab里生态最成熟,网上资料多,遇到问题搜一下就有答案。三是我所在的团队前期潮流计算、设备建模代码都是Matlab写的,整个项目在一套环境里跑通最省事。

当然,如果你所在的团队统一用Python,那用Pyomo或直接Gurobi Python API也很顺畅。工具只是手段,建模思想才是核心。下面所有代码示例我用的是Matlab+YALMIP,但约束表达方式换成Python只是语法差异。

2. 微能源网设备建模:每一个方程背后的物理含义

2.1 系统的四条母线:电、热、冷、气

我做的微能源网结构是典型的园区级配置。理解这个结构是建模的第一步,后面所有约束都是围绕母线平衡展开的。

电母线连接上级电网、燃气轮机发电、光伏、电储能、电负荷、电制冷机、电锅炉;气母线连接上级气网、燃气轮机、燃气锅炉;热母线连接燃气轮机余热、燃气锅炉、电锅炉、储热罐、热负荷、吸收式制冷机(消耗热来制冷);冷母线连接电制冷机、吸收式制冷机、储冷罐、冷负荷。每个时段t(通常取24个一小时时段或96个15分钟时段),四条母线必须各自满足功率平衡。

2.2 关键设备的线性化建模

燃气轮机组(CHP机组)是整个系统的核心,它同时产出电和热,是"多能互补"的标志性设备。核心关系是热电联产方程:发电功率等于天然气消耗量乘以发电效率再乘以低位热值,回收余热功率等于天然气消耗量乘以热回收效率乘以低位热值。这里的天然气消耗量单位用kW(按热值折算),低位热值取9.78 kWh/Nm³,也就是每立方米天然气完全燃烧产生的可用热量。

实际建模时有几个容易漏掉的细节:运行区间约束必须绑定启停变量,即出力下限乘启停变量小于等于出力,小于等于出力上限乘启停变量。这个约束的物理含义是"机组没启动时出力必须是0,启动后出力落在区间内"。爬坡约束限制相邻时段出力变化速率,即上一时段出力减去本时段出力小于等于向下爬坡速率,本时段出力减去上一时段出力小于等于向上爬坡速率。最小连续运行/停机时间在简单日前调度研究里可以忽略,但如果做日内滚动调度,最好加上,否则求解器会让机组做出"每小时启停一次"的极端策略,物理上不可行。

燃气锅炉就简单多了:热出力等于耗气量乘以锅炉效率乘以低位热值,再配上出力上下限。燃气锅炉一般不考虑启停二进制变量,因为启停成本低、动作灵活,这个假设在园区级规划里是被广泛接受的。

电制冷机和吸收式制冷机是冷母线的两个互补来源。电制冷机的制冷功率等于耗电功率乘以能效比,能效比一般在3到5之间,含义是"消耗1份电可以搬移3到5份热量"。吸收式制冷机的制冷功率等于输入热功率乘以性能系数,性能系数一般在0.7到1.2之间。这个数值差距很关键——电制冷的能效远高于吸收式制冷,吸收式制冷的价值不在效率,而在它消耗的是燃气轮机余热这种"原本可能被浪费掉"的热量。精细模型里,吸收式制冷机的性能系数在低负荷时会明显下降,简单模型保持恒定值即可,需要更高精度时按分段线性函数处理。

储能设备的数学模型在电、热、冷三种形式上是完全一样的,区别只在自放率、充放效率和容量尺度。荷电状态等于上一时段荷电状态乘以自放率系数,加上充电功率乘以充电效率,减去放电功率除以放电效率。约束包括容量上下限、充电功率上限乘充电状态变量、放电功率上限乘放电状态变量、充电状态加放电状态小于等于1。

最后这条互补约束极其关键,如果不写,求解器会在同一个时段同时充电和放电,白白浪费能量还做出一个虚假的"最优"结果。我见过不少初学者栽在这里,看起来成本很低,实际是储能系统在"左手倒右手"。

2.3 目标函数的设计:经济成本还是碳排放

我做过的两版目标函数差别挺大,值得单独说一下。

经济调度版是最常见的:最小化各时段购电成本加购气成本加设备运维成本之和。购电成本用的是分时电价乘以购电功率,购气成本是天然气单价乘以耗气量再除以低位热值(把气量从kWh换算回Nm³),运维成本按各设备出力乘以单位运维成本。这里有个关键细节:如果允许向电网售电,必须把电网交互功率拆成购电和售电两个非负变量,并加上二者不同时为正的互补约束。这个约束用大M法线性化:二者之和小于等于一个大数乘以一个0-1变量,同时二者之和小于等于大数乘以一减这个变量。否则求解器会同时买入卖出赚差价,得到一个完全虚假的最优解。

低碳调度版是在目标里加入碳排放成本,即碳配额价格乘以电网购电的碳排放因子加天然气消耗的碳排放因子。两个版本可以加权组合,调节系数alpha和beta,alpha加beta等于1,实现经济性和低碳性的平衡。这个改动在代码层面只是目标函数多几行,但对结果的影响非常大,我在后面的工况分析里会具体说。

3. 冷热电多能耦合的关键:能源集线器与矩阵化约束表达

3.1 能源集线器模型怎么用

能源集线器是苏黎世联邦理工学院提出的经典概念,核心思想是把多能源输入映射到多能源输出,写成矩阵形式:负荷向量等于耦合矩阵乘以输入向量。冷热电联供系统用这个框架表达特别清晰。

电负荷的组成部分是燃气轮机发电加上电网购电减去售电减去电制冷机耗电减去电锅炉耗电加上储能放电减去储能充电;热负荷的组成部分是燃气轮机余热加上燃气锅炉供热减去吸收式制冷机消耗的热量加上储热放热减去储热充热(如果储热存在);冷负荷的组成部分是电制冷机产冷加上吸收式制冷机产冷加上储冷放出减去储冷充入。在Matlab里直接用矩阵写,本质就是把四条母线的平衡方程整理成线性等式组的形式。

这样做的好处有两个。系统扩展方便,后面加光伏、加储能、加氢能设备,只需要在矩阵对应位置填系数,不需要重写求解逻辑。与YALMIP的约束风格天然契合,YALMIP本身就是面向等式不等式约束的建模环境,直接传入矩阵约束就能高效处理。

3.2 为什么耦合约束这么容易写错

我调试时发现,最容易出错的不是单个设备的模型,而是设备之间的耦合链条。举两个我实际踩过的例子。

第一个例子,燃气轮机的余热不是想用多少就用多少。如果系统没有烟气旁通装置,燃气轮机只要开机,余热就必须被热负荷或吸收式制冷机消纳,多出来的热量只能通过散热器排掉。这时热母线平衡方程里必须有一个弃热变量,否则容易造成两种后果:要么问题不可行(因为热量无路可去),要么求解器通过降低燃气轮机出力来"规避"余热,结果完全偏离了以热定电的实际运行逻辑。

第二个例子,电制冷机和吸收式制冷机是典型的并联互补关系。夏季供冷时,二者会争抢"电"和"热"两种不同来源的能量。单看冷母线平衡方程,看不出这种竞争关系;只有把目标函数里电价和天然气价的比值放进去,才能解释为什么某个时段特定机组多出力。有一次我把气价从2.5元调到3.5元,吸收式制冷机的出力占比直接从60%跌到20%,电制冷机顶上来了。这就是多能互补系统"整体最优而非单机最优"的直观体现。

3.3 代码层面的矩阵化做法

我的实际做法是把每个时段t的所有决策变量排成一个列向量,包含燃气轮机发电功率、耗气量、余热功率、燃气锅炉热出力、耗气量、电制冷机耗电、吸收式制冷机输入热功率、电锅炉耗电、购电、售电、储能充放电、荷电状态、弃热、各启停变量等。

然后用YALMIP的变量函数定义全时段的变量矩阵,用循环或向量化操作构建约束。这里有一个非常实际的性能经验:不要用循环一条一条写96个时段的标量约束,要基于向量化操作构建约束矩阵。我做过对比,YALMIP处理一个96行向量的约束比处理96条独立标量约束快一个数量级以上。构建完成后,用竖线连接所有约束,形成一个约束集合,交给优化函数求解。

4. 调度代码的骨架:从数据准备到结果导出

4.1 环境配置里最容易忽略的细节

我的环境是Matlab R2022a加YALMIP R20210430加Gurobi 10.0.1学术授权。Gurobi在Matlab里的路径设置其实很简单,安装后添加Matlab接口路径到MATLAB搜索路径,保存路径后重启Matlab,再运行求解器设置函数即可。

新手最容易被卡住的是license文件。Gurobi的license文件放在用户主目录下,路径不能有中文,必须和hostid匹配。我见过不少报license expired的,最后发现是系统时间不对,不是license真的过期。另外YALMIP的安装就简单多了,把下载的文件夹添加到路径里就行,没有编译过程。

4.2 核心代码结构:一个能跑通的主干

下面给出一段简化但能跑通的主干代码。为了可读性,我省略了部分储能约束和全部爬坡约束,但骨架逻辑是完整的。实际运行前必须把每个变量的上下限、储能荷电状态递推式补全。

%% 微能源网优化调度主程序骨架 % 数据准备 T = 24; % 调度时段数 load_ele = [...]; % 电负荷 1xT load_heat = [...]; % 热负荷 1xT load_cool = [...]; % 冷负荷 1xT price_ele = [...]; % 分时电价 1xT price_gas = 2.5; % 天然气单价 元/Nm3 LHV = 9.78; % 天然气低位热值 kWh/Nm3 % 决策变量定义 u_gt = binvar(1, T); % 燃气轮机启停 P_gt = sdpvar(1, T); % 燃气轮机发电功率 F_gt = sdpvar(1, T); % 燃气轮机耗气量 P_gb = sdpvar(1, T); % 燃气锅炉热出力 F_gb = sdpvar(1, T); % 燃气锅炉耗气量 P_ec = sdpvar(1, T); % 电制冷机耗电功率 P_eb = sdpvar(1, T); % 电锅炉耗电功率 P_buy = sdpvar(1, T); % 电网购电功率 P_sell = sdpvar(1, T); % 电网售电功率 SOC_ess = sdpvar(1, T); % 储能荷电状态 P_ch = sdpvar(1, T); % 储能充电功率 P_dis = sdpvar(1, T); % 储能放电功率 Q_curt = sdpvar(1, T); % 弃热功率 % 约束构建 C = []; % 电母线平衡 C = [C, 0.35*P_gt + P_buy - P_sell - P_ec - P_eb + P_dis - P_ch == load_ele]; % 热母线平衡 C = [C, 0.45*P_gt + P_gb == load_heat + Q_curt]; % 冷母线平衡 C = [C, 3.5*P_ec == load_cool]; % 燃气轮机运行区间 C = [C, 100*u_gt <= P_gt <= 1000*u_gt]; % 燃气轮机热电关系 C = [C, F_gt == P_gt/0.35]; % 储能荷电状态递推 C = [C, SOC_ess(2:T) == SOC_ess(1:T-1)*0.98 + 0.95*P_ch(2:T) - P_dis(2:T)/0.95]; C = [C, SOC_ess >= 0.2, SOC_ess <= 0.9]; % 目标函数:购电成本 + 购气成本 + 储能折旧成本 Objective = sum(price_ele.*(P_buy - P_sell)) ... + sum(price_gas*(F_gt + F_gb)/LHV) ... + 0.05*sum(P_ch + P_dis); % 求解 ops = sdpsettings('solver','gurobi','verbose',2); result = optimize(C, Objective, ops); % 结果提取 P_gt_opt = value(P_gt);

这段代码故意省略了很多边界约束,实际跑之前必须把每个变量的上下限、储能SOC的上下限都补全。YALMIP有个好处,约束没写全时它会提示问题无界之类的警告,提醒你先补约束再求解。

4.3 数据前处理:单位统一是第一位的

我在做第一个案例时犯过一个特别低级的错误:光伏出力给了kW,电负荷给了MW,天然气热值用了kJ/Nm³,结果目标函数量纲完全乱了,Gurobi直接报数值错误。后来统一成功率一律用kW、能量用kWh、气流量用Nm³/h、热值用kWh/Nm³,这才顺利求解。单位统一之后,建议用检查函数对每个变量做一轮合理性检查,YALMIP的check命令可以逐条检查约束是否满足,比你自己写断言方便得多。

5. 调度结果怎么解读:三种典型工况的对照

5.1 冬季采暖工况下的运行特征

冬季热负荷很大,燃气轮机的余热基本全部被热负荷消纳。此时系统倾向于以热定电——燃气轮机为了供热而发电,多余的电可以卖给电网(如果允许售电),或者存到电池,或者干脆让电锅炉把电转成热。

时段电价水平燃气轮机电锅炉储能
低谷0.35元/kWh低出力高出力充电
高峰1.2元/kWh全力运行低出力放电
平段0.7元/kWh中等出力中等持平

这个结果逻辑很清晰:低谷期电价便宜,电网购电划算,燃气轮机少发电,热负荷缺口由电锅炉补;高峰期电价贵,燃气轮机满发,余热供采暖,多余的电要么自用要么储能放电应对晚高峰。如果售电价很低,多余的电宁可给电锅炉弃电转热,也不愿意低价卖出。这个现象在结果里看得特别清楚。

5.2 夏季供冷工况下的运行特征

夏季冷负荷主导,这里会出现一个反直觉的结论——从能效比看,吸收式制冷机只有0.7到1.2,远低于电制冷机的3到5,按理说应该多用电动制冷,但吸收式制冷机消耗的热量来自燃气轮机余热,而燃气轮机发电本身能抵消一部分购电成本。所以在气价便宜、电价贵的组合下,发电加吸收式制冷的组合策略反而更经济。

气电价格比吸收式制冷占比电制冷占比系统总成本
气价低/电价高60%+40%以下最低
气价高/电价低20%左右80%左右最低

这就是多能互补价值最直观的体现——不是某个设备单独最优,而是整体最优。我调整气价和电价做了几组敏感性分析,这个规律相当稳定。

5.3 过渡季节的综合指标

过渡季节冷热负荷都不大,系统会处于低出力多储能的状态。我一般会重点看三个指标:24小时总购电成本、天然气总消耗量、弃热弃光总量。这三个指标放在一起,就能快速判断调度结果的合理性。

成本低了但弃热严重,说明模型里可能缺少热负荷的灵活性,或者储能容量配置偏小;天然气消耗大了但电费省了很多,说明气电比价策略下的选择还需要仔细核对。我还会把日前的调度结果和实际负荷曲线叠在一起看,如果某些时段设备出力波动过于剧烈,往往是爬坡约束设置不合理,或者负荷预测数据本身有问题。

6. 我踩过的坑:求解器报错、数值尺度与调试技巧

6.1 模型不可行到底是谁的锅

这是所有初学者最容易崩溃的瞬间。我总结了一个排查顺序。

先看优化返回的诊断信息,如果显示问题不可行,可以用赋值函数手动给变量赋一组预期合理值,然后用检查函数逐条确认哪个约束超出范围。我遇到过的两种情况特别典型。

第一种,储能荷电状态递推关系的首末衔接问题。如果要求调度周期结束时荷电状态回到初始值,而对初始值和末值的差范围限制得太紧,很容易不可行。比如初始荷电状态是0.5,要求结束时回到0.5,但储能容量小、负荷大,根本充不回去,问题就无解了。解决方法是放宽末值范围,或者允许一定的偏离并加惩罚项。

第二种,负荷数据本身超出了系统供给能力。比如夏季冷负荷峰值超过电制冷机和吸收式制冷机的总装机容量,这就是设计容量选小了。这不是约束写错,而是系统配置本身不合理。遇到这种情况,要么加大设备容量参数,要么允许切负荷——在目标函数里加一个高惩罚系数的切负荷变量。

6.2 数值尺度问题:大数和小数混在一起

Gurobi对数值尺度非常敏感。举个例子,电负荷最大5000 kW,储能荷电状态是0到1的小数,天然气热值是9.78,目标函数里的成本系数跨度可以达到10的6次方以上。这种情况下求解器内部处理会变得很不稳定,甚至出现误判不可行。

解决办法有三个。把所有量纲统一到同一基准,功率用kW、能量用kWh、成本用元,不要混入MW和MWh的表述。对储能荷电状态约束,把荷电状态乘容量作为一个整体变量来建模,避免0到1的小数和MW级的功率直接混乘。对目标函数,把所有成本统一除以一个常数,比如除以10000,不会影响最优解,但能让求解器的数值条件好很多。我在一次仿真里把所有数值除以1000之后,Gurobi的求解时间从47秒降到3秒,这是实打实的收益。

6.3 一个容易忽略的现实约束:启动成本和阶梯电价

很多教材里的简单模型不包含机组启动成本,但实际调度里,启停频繁会显著影响设备寿命,而启动成本的存在会直接改变最优解的时序结构。加上启动成本很简单,定义一个启动标志变量大于等于当前启停状态减去上一时段启停状态,同时大于等于0,再把启动标志变量乘启动成本加进目标函数。

另外,如果有阶梯电价或者容量电价,购电成本就变成了分段线性函数。分段线性函数在MILP里可以用标准的分段线性建模方法处理,也可以直接引入多个购电变量,每个变量对应一个价格区间,用二进制变量保证只落在其中一个区间。这些现实约束加进去之后,模型复杂度上一个台阶,但调度结果的可信度也会明显提升。

6.4 两个亲测有效的调试习惯

第一个,先用少时段调试,再放到全时段。6个时段的MILP几乎是秒出结果,你可以在几分钟内反复调整约束定义,等到逻辑全部正确再扩展到24或96时段。我见过太多人一上来就上96时段,结果连约束手滑写错都发现不了,求解器报个不可行就傻眼了,排查难度直接翻好几倍。

第二个,把每个设备的出力和对应负荷曲线一起画出来。解除变量值之后直接用plot和area函数,面积图画各设备的贡献叠加,一眼就能看出某个时段总有设备出力异常。可视化不只是写论文时才用的工具,调试阶段它是最直观的"照妖镜",比盯着数字矩阵高效得多。

回到最初的问题。我踩过不少坑之后才意识到,做微能源网优化调度,最难的部分往往不是最后那个求解器调用,而是建模时对物理过程的理解、对耦合约束的表达、对数据尺度的把控。先把一个能跑的MILP骨架搭起来,再逐步加设备、加约束、加不确定性场景,这套方法论让我从第一个案例到后来做园区级综合能源规划项目都走得很顺。如果你正准备开始做类似方向,建议别在最开始追求模型有多复杂,先把代码跑通,再迭代升级,这是效率最高的路径。

需要专业的网站建设服务?

联系我们获取免费的网站建设咨询和方案报价,让我们帮助您实现业务目标

立即咨询