☰
基于Matlab的热电联产机组联合优化控制消纳风电实践
2026/10/11 7:53:39 网站建设 项目流程

这两天整理风电消纳的东西,又翻出之前用Matlab写热电联产机组联合优化控制的代码。说个调度台常见的画面:冬季大风天的深夜,风电大发,电负荷却在低谷,而供热机组因为“背着”居民供暖任务,电出力被压在一个高位降不下来。电网想多消纳风电,可全系统的下调空间已经见底,结果只能眼睁睁看着风机限功率运行。这个问题的根源,就是热电联产机组的电出力被热负荷“绑架”了。而今天要聊的这套联合优化控制,正是通过优化热电联产机组的“电-热”工作点,把供热系统里藏着的柔性空间挖出来,用于最大化消纳风电。我会结合Matlab代码实现,把从建模到求解、再到踩坑的经验一次说清楚。适合正在做新能源调度、热电联产机组建模,或者想入电力系统优化这个方向的同学参考。

1. 热电联产机组为什么总跟风电“打架”:先搞懂以热定电

1.1 机组类型决定了调节自由度:背压式与抽汽式

热电联产机组分两大类。背压式机组,汽轮机排汽全部进入供热管网,发电和供热是同一个蒸汽流量决定的,电出力完全跟着热负荷走,从调度角度看等于没有调节余地。抽汽式机组相对灵活一点,中压缸排汽一部分抽出去供热,另一部分继续进入低压缸发电,电热之间存在耦合,但耦合区间仍然非常“紧”。

实际工程里,抽汽式机组在供热工况下的电出力区间,并不是固定不变的。热负荷越高,为了维持供热抽汽压力,主蒸汽流量就得越大,机组最小电出力也跟着被迫抬高。冬季供热高峰时,一台机组的电出力可能被热负荷限定在额定容量的70%甚至80%以上运行,剩下能调节的“抽屉”只有最顶部那一条窄缝。这就是为什么供热期系统调峰能力骤降。

1.2 冬季弃风的真正推手:系统下调节能力不足

风电出力本身带有强烈的反调峰特性——白天负荷高的时候风往往小,夜里负荷低的时候风反而大。冬季更极端,夜间热负荷处于一天峰值,热电机组全部在“满供”状态,电出力下限被抬到很高。

这时候看整个电网的调节资源:常规火电已经压到锅炉稳燃极限,大约30%到40%额定出力;水电站受来水或库容限制也未必能扛;剩下的负荷变化空间本来就小,风电一上来,系统没有地方“腾位置”。从机组个体的角度,它觉得自己只是“按供热情务完成发电”,但从系统整体看,它占用了大量的消纳空间,最终风电被逼到极限边缘,只能用弃风这种方式来维持功率平衡。这就是业内常说的“以热定电”死结。

2. 联合优化控制的破题思路:把供热系统的柔性挖出来

2.1 供热系统原本藏着大量可调节能力

很多人以为热电联产机组的灵活性只能靠机组本体改造,比如低压缸切缸、旁路供热之类的。其实供热系统本身就有大量柔性资源被浪费掉了。热网管道里有存水量,整个管网的蓄热能力相当于一个巨大的“热电池”;建筑物本身也有热惯性,室温在一定范围内波动,人在屋里是基本感知不到的;蓄热罐更是直接可以充热放热的缓冲装置。

联合优化控制的核心思路,就是把机组的电调峰压力和供热系统的热惯性放在同一个优化问题里考虑。让热电机组的一部分供热量在某一时段由蓄热罐或热网存量来承担,机组的热出力就能降下来,电出力下限随之下降,风电就能多发一些。反过来,在风电小或者电价高的时候,机组多供热,把蓄热罐充起来。这一放一充,就是调峰空间。

2.2 本项目里联合优化控制的控制架构

我代码里采用的架构是这样的:输入侧是未来一个调度周期(通常24小时,时间分辨率15分钟或1小时)的风电预测出力、热负荷预测、电负荷预测,以及机组的运行参数。优化层输出的是每台机组各时段的电出力计划、热出力计划、启停状态,以及蓄热罐的充放热策略。输出后直接下发到机组集控和供热集控执行。

这个架构里有一个关键点:热负荷预测和电负荷预测同步纳入模型。很多同学做调度只看电平衡,把热负荷当成一个固定的“外挂条件”硬塞进去,结果优化出来的机组工作点根本不在可行域内。真正的热电联合优化,必须把热系统变量当成决策变量的一部分,让模型自己决定“这一小时热由机组烧、由蓄热罐放、还是依靠管网惯性扛”。

2.3 调度层优化与底层控制的边界

需要说明,这里的联合优化控制解决的是调度层的机组工作点优化问题,时间尺度是分钟级到小时级。底层执行还是靠机炉协调控制系统和供热调节阀的秒级调节。调度层给出工作点,底层负责跟踪,两层之间靠指令对接。

我在做这个项目的时候没有重写底层控制逻辑,也不需要。只要调度指令是机组可行域内的、爬坡速率可达的,底层就能稳定跟踪。很多文章把问题描述成“实时优化控制”,实际算例里跑的还是日前调度模型,时间尺度对不上。做代码时先把边界划清楚,后续扩滚动调度才不容易乱。

3. 目标函数与约束条件的取舍:风电最大化消纳怎么落到数学上

3.1 目标函数:弃风惩罚系数是唯一的关键旋钮

优化目标我选的是系统运行总成本最小,包含四块:机组燃料成本、启停成本、弃风惩罚、备用成本。燃料成本用线性或分段线性函数逼近,启停成本只在机组状态变化时计一次,备用成本按旋转备用量计。

弃风惩罚项是“风电最大化消纳”从口号变成数学语言的关键。目标函数写成弃风量乘以惩罚系数,模型为了降低总成本,会自动减少弃风。但惩罚系数不能拍脑袋定:定得太低,模型觉得弃风无所谓,消纳效果差;定得太高,模型会让机组疯了一样爬坡去消纳最后1兆瓦风电,结果煤耗暴涨、爬坡越限,算出来的方案根本不现实。按我调参的经验,弃风惩罚系数取机组边际运行成本的1.2到1.5倍比较稳妥,既能让模型优先消纳风电,又不至于为了消纳而过度牺牲机组经济性。

3.2 抽汽机组电热可行域约束:别用简单的上下限代替

这是建模里最容易出错的地方。抽汽式机组在电热平面上不是一个矩形可行域,而是一个凸多边形。它的边界由汽轮机进汽量上限、低压缸最小冷却流量、抽汽压力、锅炉最大蒸发量等共同决定。常见做法是用一组线性不等式来包络这个多边形。

我在Matlab里采用的方法是顶点法:先根据机组设计参数和典型工况表,提取出电热可行域的若干个顶点,然后通过顶点生成半空间不等式约束。以一台典型的300兆瓦抽汽机组为例,不供热时电出力区间约为150到300兆瓦;热出力50兆瓦时,电出力下限抬到约190兆瓦;热出力100兆瓦时,电出力下限接近230兆瓦。这些工况点连起来,就形成了一个向右上方倾斜的下边界。直接用固定上下限建模,会高估甚至低估机组的调峰能力,导致优化结果不可行。

3.3 功率平衡、爬坡与备用约束:把系统规则逐条写进去

电功率平衡约束要求各机组电出力之和加上实际上网风电,等于电负荷;热功率平衡约束要求各机组热出力加上蓄热罐放热量,等于热负荷。爬坡约束限定机组相邻时段电出力变化量不超过上限;旋转备用约束要求机组在某一方向上的可调容量满足系统备用需求。

这里有个细节:热出力变化本身也受锅炉燃烧调整速度限制,不能只看电爬坡。很多模型只约束电出力爬坡,忽略热出力变化率,结果算出来热负荷曲线剧烈抖动,现场执行不了。我实际建模时会给热出力也加上变化率约束,数值上通常取电爬坡上限的1.5到2倍,因为热负荷调节相对缓慢。

3.4 风电出力相关约束:弃风量的上下界要卡住

风电实际上网功率等于预测可用功率减去弃风量。弃风量是非负变量,且不能超过预测功率。这里还要注意,风电出力可能存在最大技术出力限制,通常按场站装机容量乘一个同时率来设定。我见过有人把风电当成“来了就能全额消纳”的刚性电源,只写一个上下限,结果模型为了降成本,让风电在预测值以内随便削减,这其实也是一种变相弃风,必须通过目标函数中的惩罚项和弃风变量共同约束住。

4. 那个最麻烦的地方:非线性项的线性化处理

4.1 为什么不直接上非线性求解器

热电联产机组运行优化问题,如果煤耗曲线用二次函数描述,可行域用线性不等式描述,加上启停状态的整数变量,就是一个混合整数二次规划。小规模算例可以直接丢给非线性求解器跑,但规模化之后求解时间不可控,而且容易陷入局部最优。

工程上更稳妥的做法是把燃料成本函数分段线性化,把整个问题转成混合整数线性规划(MILP)。MILP有一大优势:商用求解器能在有限时间内给出带最优性间隙(gap)的全局解,调度人员能清楚知道当前解离最优解有多远。这对实际运行非常重要——你不能告诉调度员“这个解可能好,也可能不好”,你得告诉他“当前解与最优解的差距不超过2%”。

4.2 煤耗曲面的分段线性化思路

电热联产机组的煤耗随电出力和热出力同时变化,形成一个三维曲面。线性化处理时,按热出力区间分段,每个区间内煤耗近似为电出力的线性函数。分段的边界点来自于机组热力试验报告,通常取额定工况、纯凝工况、最大抽汽工况等几个典型工况点。

具体代码实现时,我会为每台机组预先算好各分段点的煤耗率,然后用凸组合的方式把煤耗成本表达成线性约束。这里的细节是:如果用分段线性函数逼近凸函数,需要确保分段点是凸组合的,否则模型会“钻空子”选择非物理的分段权重。

4.3 整数变量的三个典型用途

这个模型里整数变量有三个来源。第一是机组启停状态,0表示停机1表示运行,启停成本通过状态跳变的0-1约束表达。第二是蓄热罐的充放状态,同一时刻只能充或者只能放,需要一个二进制变量切换。第三是分段线性化的区间选择,每个分段对应一个0-1变量。

整数变量一多,求解时间会指数上涨。我的经验是:24时段、2台机组、1个蓄热罐的小规模问题,大概有200个左右的0-1变量,商用求解器几十秒就能到1%的gap;如果机组数量到10台、时段到96个,变量数上千,这时候就要靠设置求解时间上限和gap容忍度来控制计算时间。

4.4 热网的时间延迟怎么处理才合理

严格意义的热网动态需要偏微分方程描述,直接放进优化模型不现实。我在代码里做了两层简化:第一,把建筑物和管网的热惯性等效成一个一阶惯性环节,在热平衡方程里加入一个热储能状态变量;第二,通过蓄热罐的充放策略来吸收热负荷预测误差和机组调节延迟。

更细一点的做法是给热负荷预测做平滑处理。实际热负荷曲线没有电负荷那么“毛糙”,居民供暖的室温变化周期很长,半小时以内的波动完全可以被管网和建筑本体吸收。所以模型里热平衡约束用的是15分钟或1小时平均热负荷,而不是瞬时值。这套处理方式在工程上足够用,而且大幅降低了建模复杂度。

5. Matlab代码实现:从变量声明到求解器调用的完整框架

5.1 代码文件组织与数据准备

整个Matlab工程我按模块拆了四个文件:

文件职责
input_data.m定义机组参数、负荷曲线、风电预测
build_model.m构建决策变量、约束和目标函数
solve_model.m调用求解器并输出结果
plot_result.m绘制电热平衡图和弃风分析图

数据准备阶段我会把热负荷曲线、电负荷曲线、风电预测曲线都整理成矩阵,行对应时段,列对应节点或机组。这里建议做一次单位统一检查:功率用兆瓦,能量用兆瓦时,热值换算系数单独定义成常量。我吃过一次亏,热负荷单位用了吉焦每小时,电功率用了兆瓦,等式两边硬是差了一个换算系数,查了很久才发现。

5.2 核心代码:变量定义与约束组装

这是利用YALMIP工具箱构建优化模型的核心片段:

%% 决策变量定义 Pg = sdpvar(nG, T); % 机组电出力, MW Hg = sdpvar(nG, T); % 机组热出力, MWth Pw = sdpvar(nW, T); % 风电实际上网, MW Pc = sdpvar(nW, T); % 弃风量, MW u_on = binvar(nG, T); % 机组启停状态, 0/1 Hst = sdpvar(1, T); % 蓄热罐放热功率, 正为放热负为充热

约束组装的写法有一些固定的套路。电功率平衡、热功率平衡放在最前面,然后是机组可行域约束、爬坡约束、蓄热罐约束,最后才是目标函数中的各项成本:

Constraints = []; %% 电功率平衡: 各机组出力之和 + 实际上网风电 = 电负荷 Constraints = [Constraints, sum(Pg, 1) + sum(Pw, 1) == P_load']; %% 热功率平衡: 机组供热量 + 蓄热罐放热 = 热负荷 Constraints = [Constraints, sum(Hg, 1) + Hst == H_load']; %% 风电电量约束: 实际上网 = 预测值 - 弃风, 弃风非负 Constraints = [Constraints, Pw == P_forecast - Pc, Pc >= 0]; %% 蓄热罐容量与充放速率 Constraints = [Constraints, Hst >= -Hst_max, Hst <= Hst_max]; Constraints = [Constraints, SOC_min <= SOC_init + cumsum(Hst) * dt / Cap <= SOC_max];

目标函数按成本项分别累加:

%% 目标函数: 燃料成本 + 启停成本 + 弃风惩罚 fuel_cost = sum(sum(fuel_cost_coef .* Pg)); start_cost = sum(sum(start_cost_coef .* max(0, diff(u_on, 1, 2)))); curtail_penalty = lambda_curtail * sum(sum(Pc)); Objective = fuel_cost + start_cost + curtail_penalty;

这里燃料成本我用了线性系数,实际工程中如果煤耗曲线弯曲比较明显,建议按4.2节的分段线性方法处理,效果会更好。

5.3 求解器选择与参数设置

求解这块,YALMIP的好处是可以无缝切换底层求解器。学习和小规模算例,直接用Matlab自带的intlinprog就够了,不需要额外的license。算例规模变大、机组多于4台或者时段多于48个,建议换Cplex或Gurobi,求解速度差距非常明显。

ops = sdpsettings('solver', 'gurobi', 'mipgap', 0.01, 'verbose', 2); optimize(Constraints, Objective, ops);

注意到我把mipgap设成了1%,不是默认的0%。这不是偷懒,是工程上的必要取舍——对于多机组联合优化这种大规模MILP问题,1%的gap已经能满足调度精度要求,而求解时间可以缩短一个数量级。

5.4 结果可视化与输出

结果分析阶段,我最看重两张图:一张是电力平衡堆叠图,能直观看到风电出力被压缩在哪个时段;一张是热平衡图,能看到蓄热罐的充放节奏是否和风电消纳高峰错开。可视化代码用Matlab的area函数就能实现:

figure; area([Pg' Pw']); hold on; plot(P_load, 'k-', 'LineWidth', 2); legend('机组1', '机组2', '风电', '电负荷'); xlabel('时段'); ylabel('功率/MW');

还有个细节:输出结果时一定要保存弃风量曲线和蓄热罐SOC曲线,这对后续做灵敏度分析非常有用。我从一开始就把所有决策变量存成.mat文件,后面写论文、回查数据都靠它。

6. 算例验证与跑完代码后的三个大坑

6.1 一个能复现的小系统算例结果

为了验证代码逻辑,我搭了一个两机组加一个风电场加一个蓄热罐的小系统做测试。供热机组参数参照典型抽汽机组,风电预测数据取了冬季大风天的典型曲线,时间分辨率为1小时,调度的周期为24小时。

场景弃风电量/MWh弃风率系统煤耗/tce
固定热电比(不参与优化)142.612.8%321.5
联合优化(机组可调+蓄热罐)32.52.9%323.6
联合优化+热网惯性21.81.9%324.2

从结果看,单纯把机组电热工作点放开,弃风率就能从12.8%降到2.9%,代价是煤耗增加了约2吨标准煤。这说明风电消纳不是免费的——消纳风电需要机组在低负荷工况运行,锅炉效率会下降,煤耗必然微涨。实际做工程决策时,这就是一个“经济性换清洁性”的取舍问题。

6.2 坑一:弃风惩罚系数没有做灵敏度分析

我最早跑算例时,惩罚系数直接拍了0.8元/千瓦时,结果模型把2号机组在夜间来回启动停机,就为了多消纳最后几兆瓦风电。从弃风率看数字很漂亮,但看机组动作曲线完全没法用。

后来我补了一组灵敏度分析:惩罚系数从0.2到1.0元/千瓦时,步长0.05,每个系数各跑一遍,把弃风率和机组燃料成本的变化画成曲线。当系数超过一定阈值后,弃风率下降曲线明显变平缓,说明再增加惩罚只是让机组做无谓的出力调整。这个“拐点”位置就是该算例的最优惩罚系数。

6.3 坑二:忽视热网延时导致调度指令超前

联合优化算出来的方案,如果直接下发执行,热负荷侧的响应是滞后的——管网有传输延时,建筑物吸收热量也需要时间。我在最早一版代码里没考虑这层,热平衡约束用的是即时热负荷,结果蓄热罐白天放热,晚上充热,现场反馈说热用户端温度波动超标。

解决办法就是4.4节说的热惯性等效和蓄热罐SOC约束。我在热平衡方程里加了一个等效热储能项,并用SOC的上下限和充放速率来限制热计划的剧烈波动。修改之后,优化结果里蓄热罐的动作节奏明显平缓了,热网实际运行时室温波动控制在1摄氏度以内。

6.4 坑三:求解gap的目标设得太死,反而让方案失去时效性

调度计划是有时效性的——你现在算出来的结果是给下一个小时用的,如果求解器跑了两个小时才收敛,方案已经过期了。我一开始把mipgap设成0.001%,9台机组96时段的问题跑了一晚上没停机,第二天起来发现目标值已经2小时没有明显下降了。

后来我把gap放宽到1%,求解时间从接近2小时缩短到3分钟,目标函数值只多花了0.7%的成本。这个误差在工程调度完全可接受。建议大家在跑大型算例时,先设一个较宽的gap拿到可行解,再逐步收紧,同时设定求解时间上限,比如600秒,到时间就取当前最优可行解。

7. 写在最后:这套代码后续还能往哪走

现在这套Matlab代码还在持续更新,我最近在做两个方向的扩展。一个是最小化弃风和最大化供热的滚动协调问题,调度计划从“一次性算出24小时”改成“每小时滚动刷新,锁定前3小时计划”,这样能利用更新的风电预测信息,显著降低预测误差带来的弃风。另一个方向是加入储电和蓄热罐的联合协调,在算例里同时考虑电储能和热储能,让系统既有电的灵活性又有热的灵活性。

在跑大规模算例之前,建议先在小系统上把代码通一遍,确认所有约束的单位、可行域边界、惩罚系数都没有问题,再逐步扩大规模。从我个人的使用体验来说,这套代码框架最大的价值不在于某一个公式,而在于它把“电热耦合优化”这件事从理论推导变成了可以反复试错的工具——调参、加约束、换场景,都是几分钟的事。做优化的乐趣也恰恰在这里。

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

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

立即咨询