1. 项目概述与核心思路拆解1.1 为什么是“电热综合”而不是单纯“电力”做电力系统优化调度的人大概率都遇到过这样一个尴尬场景电网侧算出来的消纳能力明明够结果现场还是弃风了。你翻来覆去查数据最后发现瓶颈根本不在电力系统本身而在供热系统——尤其是我国北方地区冬季供暖季热电联产机组“以热定电”的刚性约束直接把风电的上网空间挤得所剩无几。我最初接触这个项目时第一反应也是“这不就是个机组组合加经济调度吗”但真正把热力系统计入约束条件之后问题复杂度立刻不一样了。电热综合能源系统Integrated Electricity and Heat Energy System的核心特征是电力网络与热力网络通过热电联产机组CHP、电锅炉、储热罐等耦合设备深度绑定。电力系统需要实时平衡热力系统则可以凭借管网的蓄热特性和储热装置的缓冲在一段时间内“缓缓再算账”。这两个时间尺度差异巨大的系统放在一起做日前调度不是简单把两组约束拼起来就行而是要重新设计目标函数和决策变量的耦合关系。这个项目解决的核心痛点就是在供暖季高风电出力时段通过电锅炉和储热罐把多余的电力转化成热能储存起来让CHP机组可以降电出力而不必同步降热出力也就是通常说的“热电解耦”。这样一来风电可以多上网燃煤可以少烧一点系统总运行成本也随之下降。项目标题里“考虑可再生能源消纳”这几个字落脚点就在于此。1.2 日前经济调度的核心逻辑“日前”两个字指的是提前24小时或更长制定次日的调度计划。为什么一定要做日前因为火电机组的启停不是按开关那么简单锅炉升温、汽轮机冲转都需要时间如果等风来了再开机黄花菜都凉了。所以运行人员需要根据负荷预测和风电预测提前一天确定每台机组的启停状态、出力区间和热出力计划这就是日前经济调度的意义。经济调度的核心逻辑是在满足电力平衡、热力平衡、机组物理约束和安全约束的前提下让系统总运行成本最小化。换句大白话就是——用最便宜的方式把第二天的电和热都供上。这里“最便宜”不只是燃料成本还包括启停成本、弃风惩罚成本、购电成本等。具体到Matlab代码实现就是一个典型的混合整数线性规划MILP或混合整数二次规划MIQP问题决策变量里既有连续变量机组出力、热出力、储热罐充放热功率也有0-1整数变量机组启停状态。编写代码之前一定要把模型里所有时间断面之间的耦合关系理顺。日前调度是典型的多时段耦合优化火电机组有爬坡约束储热罐有容量变化约束CHP机组的热出力会影响下一个时段的电出力范围。如果没有先画清楚这些时序耦合关系就直接开始写Yalmip代码后面大概率会陷入“约束总是不对”的泥潭——这是我实际踩过的坑后面会详细展开。1.3 可再生能源消纳的建模思路可再生能源消纳的建模最常用的方式是引入弃风弃光惩罚项。为什么不直接约束风电必须全额消纳因为从数学优化角度硬性约束“风电出力等于预测值”会把所有灵活性都锁死系统可能因为风电预测误差而无法找到可行解。更合理的做法是允许调度模型“丢弃”部分风电但对弃风行为施加一个足够高的惩罚成本让优化器在权衡之下主动选择消纳风电。这个“足够高”很有讲究。惩罚系数太低优化器会发现弃风比调CHP机组更“划算”结果弃风量明显偏大惩罚系数太高又可能让模型过度牺牲经济性去追逐风电消纳导致总成本上升。我个人的经验是弃风惩罚系数设置在系统边际煤耗成本的1.5到2倍左右比较合适具体数值需要根据测试系统参数微调。另外模型中还需要对风电出力预测曲线做平滑处理。原始预测数据往往带有高频抖动而日前调度的时间粒度一般是1小时如果直接使用未处理的原始数据可能会让相邻时段的调度结果出现不合理跳变。我在项目中采用了简单的高斯滤波处理效果良好但这部分预处理代码容易被初学者忽略。2. 模型构建与关键技术点2.1 核心设备建模电热综合能源系统的设备模型是整个调度模型的地基。地基不稳后面做再多优化都白搭。我以这个项目中最典型的系统配置为例CHP机组、纯凝式火电机组、风电场、电锅炉、储热罐外加一个热负荷需求。CHP机组是全系统建模难度最高的部分。根据热电运行特性CHP机组分为背压式和抽汽式两类。背压式机组的电出力与热出力基本上是一一绑定的也就是说发多少电就必然产多少热几乎没有调节余地对系统灵活性贡献有限。抽汽式机组则灵活得多其可行运行域是一个四边形或更复杂的凸多边形区域。在Matlab中描述这个运行域需要把机组的最小电出力、最大电出力、最大进汽量、最小凝汽流量等参数转换成一组线性不等式约束。我项目中用的是抽汽式CHP它的电出力P和热出力H满足如下关系P_min Cv * H P P_max - Cv * H 范围随热出力变化 H 处于 [H_min, H_max]这里Cv是热电耦合系数也称为电出力随热出力变化的斜率单位是MW/MWth。具体数值每台机组不一样需要查厂家资料或参考标准测试数据。这个公式看着简单但约束方向容易写反我初次建模时就把不等式方向搞错了一次结果优化结果里CHP机组的电出力与热出力关系完全不符排查了好半天。电锅炉的建模比较直接输入电力输出热能转换效率在0.95到0.99之间出力上下限和爬坡约束可以参照电热转换特性设定。它的作用就是为风电提供一条“电转热”的消纳路径相当于给电力系统增加了一个可控的、时间灵活的用电负荷。储热罐是时间解耦的关键设备。它的核心约束是容量动态平衡方程S(t1) S(t) eta_charge * P_ch(t) * delta_t - P_dis(t) / eta_discharge * delta_t同时需要限制储热罐的储热容量上下限、单位时间充放热功率上限以及一天周期结束时的储热量需要回到初始值或允许在一定范围内偏离——这个约束常用于保证日前计划的可重复性。注意充热效率和放热效率通常不同而且有的简化模型会忽略效率直接用S(t1) S(t) (P_ch - P_dis) * delta_t这在工程分析中也可以接受但会高估储热罐的实际作用需要注意。纯凝火电机组的建模就是标准的机组组合模型包含出力上下限、爬坡约束、最小启停时间约束等。最小启停时间约束在日前调度中非常重要但在入门教程里经常被忽略导致模型结果在实际运行中根本无法执行——开机1小时就停机的情况现实中是不可能的汽轮机受不了这种折腾。2.2 目标函数设计与成本构成经济调度的目标函数本质上是把所有运行成本折算成“一天的总账单”然后让这个账单最小。这个项目中的成本构成可以分为以下几类燃煤成本CHP机组和纯凝火电机组的煤耗成本。通常用二次函数拟合机组煤耗特性C_fuel(P) a * P^2 b * P c在MILP模型中二次函数需要分段线性化处理。Yalmip自带pwlin或implies类似工具可以实现也可以用手动分段的方式。分段数一般取35段就能获得不错的近似精度。弃风惩罚成本模型允许弃风但是“弃风一时爽代价要付清”。惩罚系数设为W_pen单位是元/MWh乘以弃风量就是总惩罚成本。购电成本如果系统与上级电网存在联络线需要计入从外部购电的费用。有的模型中还需要考虑向外部售电的收入这需要根据售电电价是否低于购电电价来决定是否设非负约束避免“低买高卖”套利。启停成本机组从停机状态到开机状态需要额外消耗燃料和寿命这部分用固定成本项表示与是否发生启停动作0-1变量绑定。我建议在编写代码时先把目标函数的所有项在纸上列一个清单标注每项的量纲和系数来源然后再开始写Matlab代码。这个习惯在很大程度上避免了后面“优化结果总感觉哪里不对但又说不出原因”的困境。2.3 约束条件全景梳理约束条件是模型中最容易出错也最花时间的部分。我把这个项目中的约束条件按类型整理如下每一项都是代码中需要逐一实现并反复检查的内容电力平衡约束每个时段所有电源出力总和加上购电如果有减去电锅炉耗电等于该时段的电负荷需求。注意电锅炉在这里是“用电设备”在平衡方程中是以负荷形式出现的这一点的处理方式经常让初学者犯迷糊。热力平衡约束所有热源产热CHP机组热出力电锅炉热出力-储热罐充热功率储热罐放热功率等于该时段的热负荷需求。储热罐充热相当于增加热负荷放热相当于增加热源供给。机组出力约束包括电出力上下限、热出力上下限以及抽汽式CHP机组的热电耦合不等式。爬坡约束机组相邻时段的出力变化不能超过爬坡速率乘以时间间隔。这里的出力指的是电出力热出力是否参与爬坡约束需要根据模型复杂度决定。我项目的做法是电出力参与爬坡约束热出力单独设置热爬坡约束这样更贴近实际运行特性。储热罐约束容量上下限、充放热功率上下限、周期末储热量约束。风电出力约束风电实际出力介于0和预测值之间差值即为弃风量。系统备用约束可选为应对预测误差系统需要预留一定旋转备用容量。这个约束在日前调度中是否加入取决于研究需要。我项目中未加入备用约束以控制模型复杂度但如果你做的是工程化应用建议还是加上。所有约束条件都涉及“每个时段”的索引变量在Matlab中可以用循环方式构建也可以利用Yalmip的矩阵化表达式高效构建。循环写起来直观但速度慢矩阵化写起来高效但可读性差我的建议是先用循环跑通小规模测试再逐步优化代码结构。3. Matlab代码实现与实操要点3.1 工具箱选型与求解器配置Matlab做优化调度几乎绕不开Yalmip工具箱。它不是求解器而是“建模语言”作用是把人类易读的优化模型翻译成求解器能理解的标准形式。Yalmip最大优势在于你不需要关心求解器内部的数据结构只需要用sdpvar定义变量、用optimize命令求解它自动帮你做模型转换和求解器连接。求解器方面MILP问题推荐CPLEX或Gurobi这两个商业求解器的分支定界算法效率远超开源求解器。对于小型测试系统比如10台机组、24个时段Gurobi往往几秒钟就能求解CPLEX的求解速度也接近。如果你手里没有商业求解器的licenseYalmip也支持开源的CBC求解器求解速度慢一些但结果可靠。学校环境一般都能申请到Gurobi的免费学术license建议优先选择。安装Yalmip后第一件事是运行yalmiptest命令验证安装是否正常。这个问题我遇到不少求助案例多数是路径没设置好或者Matlab版本兼容性问题。Yalmip不需要安装直接把文件夹加入Matlab路径即可但注意不要和Matlab工具箱路径混淆。3.2 数据初始化与单位统一建模前必须先统一单位这是整个项目中最容易“翻车”的环节。电力系统常用的单位有MW和MWh但如果热力部分用GJ或Gcal那你在写约束时就得反复换算稍不留神就出乱子。我建议统一采用以下单位体系功率MW能量MWh热量MWth热功率能量用MWh热表示对于热负荷和储热罐容量如果原始数据给的是GJ换算关系是1 MWh 3.6 GJ。如果你拿到的是Gcal那就更麻烦一些1 Gcal约等于1.163 MWh。强烈建议写一个单位换算函数在读取原始数据后统一转换为标准单位制然后再送入模型。这一步做扎实了后面调试会省很多时间。数据初始化部分的代码我习惯用一个结构体struct来管理所有系统参数。比如%% 系统参数定义 sys struct(); sys.T 24; % 调度时段数日前24小时 sys.dt 1; % 时间间隔小时 sys.n_chp 2; % CHP机组台数 sys.n_thermal 1; % 纯凝火电机组台数 %% CHP机组参数 sys.chp(1).P_min 50; % 最小电出力 (MW) sys.chp(1).P_max 200; % 最大电出力 (MW) sys.chp(1).H_min 20; % 最小热出力 (MWth) sys.chp(1).H_max 150; % 最大热出力 (MWth) sys.chp(1).Cv 0.15; % 热电耦合系数 sys.chp(1).a 0.0002; % 煤耗曲线二次项 sys.chp(1).b 0.3; % 煤耗曲线一次项 sys.chp(1).c 20; % 煤耗曲线常数项这样做的最大好处是后续构建约束时可以直接用sys.chp(1).P_min这种可读性极高的字段名而不是一堆散落的变量名。当你代码写到300行以上时这种数据组织方式的优势会非常明显。3.3 决策变量定义与约束构建Yalmip中的变量定义是建模的第一步。在这个项目中需要以下几种类型的变量%% 决策变量定义 P_chp sdpvar(sys.n_chp, sys.T, full); % CHP电出力 H_chp sdpvar(sys.n_chp, sys.T, full); % CHP热出力 P_th sdpvar(sys.n_thermal, sys.T, full); % 火电机组电出力 P_wind sdpvar(1, sys.T, full); % 风电实际出力 P_wind_curtail sdpvar(1, sys.T, full); % 弃风功率 P_eb sdpvar(1, sys.T, full); % 电锅炉耗电功率 H_eb sdpvar(1, sys.T, full); % 电锅炉产热功率 P_ch_st sdpvar(1, sys.T, full); % 储热罐充热功率 P_dis_st sdpvar(1, sys.T, full); % 储热罐放热功率 S_st sdpvar(1, sys.T 1, full); % 储热罐储热量含初始时刻 u_chp binvar(sys.n_chp, sys.T); % CHP机组启停状态可选我特别说明一下S_st为什么定义到T1因为储热罐的容量变化是一个跨时段递推过程我们需要保存初始时刻时段1开始前和每个时段结束时的储热量这样就能方便地写出容量平衡约束。很多初学者只定义T个变量结果递推约束写到最后发现少了一个边界值处理起来非常痛苦。约束构建部分最具代表性的几个约束代码如下%% 电力平衡约束 Constraints []; for t 1:sys.T Constraints [Constraints, sum(P_chp(:, t)) sum(P_th(:, t)) P_wind(:, t) ... P_load(t) P_eb(t)]; end %% 热力平衡约束 for t 1:sys.T Constraints [Constraints, sum(H_chp(:, t)) H_eb(t) P_dis_st(t) - P_ch_st(t) ... H_load(t)]; end %% CHP机组运行域约束 for k 1:sys.n_chp for t 1:sys.T Constraints [Constraints, P_chp(k, t) sys.chp(k).P_min sys.chp(k).Cv * H_chp(k, t)]; Constraints [Constraints, P_chp(k, t) sys.chp(k).P_max - sys.chp(k).Cv * H_chp(k, t)]; Constraints [Constraints, H_chp(k, t) sys.chp(k).H_min, H_chp(k, t) sys.chp(k).H_max]; end end %% 储热罐容量动态平衡 for t 1:sys.T Constraints [Constraints, S_st(t 1) S_st(t) sys.eta_ch * P_ch_st(t) * sys.dt - ... P_dis_st(t) / sys.eta_dis * sys.dt]; end这里需要注意Yalmip约束拼接使用[Constraints, new_constraint]的语法而且约束之间不要用逗号之外的分隔符。很多初学者在这里会用分号导致约束维数错乱报错信息还特别不直观。3.4 求解流程与结果输出求解部分的代码相对简单%% 目标函数 Objective 0; % 煤耗成本分段线性近似 for k 1:sys.n_chp Objective Objective sum(sys.chp(k).a * P_chp(k, :) .^ 2 ... sys.chp(k).b * P_chp(k, :) sys.chp(k).c); end % 弃风惩罚 Objective Objective sys.w_pen * sum(P_wind_curtail); % 电锅炉运行成本维护成本一般较小 Objective Objective sys.c_eb * sum(P_eb); %% 求解 ops sdpsettings(solver, gurobi, verbose, 2, showprogress, 1); sol optimize(Constraints, Objective, ops); %% 结果提取与验证 if sol.problem 0 P_chp_opt value(P_chp); H_chp_opt value(H_chp); P_wind_opt value(P_wind); % ... 其他变量提取 else disp(求解失败错误代码: string(sol.problem)); diagnose(Constraints, Objective, ops); end关于二次项的处理我在这里用的是直接二次函数。Yalmip配合Gurobi可以直接处理MIQP问题不需要手动分段线性化。但如果使用CBC求解器就必须先把二次项分段线性化否则求解器会报错。这是一个常见的“代码在不同求解器之间移植”的坑值得留意。结果输出部分我建议不要只输出数值最好同时生成可视化的调度曲线图。Matlab的plot命令足够应付画出电出力曲线、热出力曲线、储热罐容量变化、风电消纳情况等关键图一眼就能看出调度结果合不合理。4. 典型场景仿真与结果分析4.1 场景设计与参数配置为了验证模型的有效性我设计了三个典型场景分别对应不同的系统配置形成对比分析场景A基准场景不含电锅炉不含储热罐。这是最传统的CHP系统配置也是弃风问题最严重的配置。在这个场景下CHP机组只能通过降低热出力来降低电出力因此供暖季的风电消纳空间非常有限。场景B含电锅炉引入电锅炉作为电转热的灵活性资源但无储热罐。电锅炉可以在风电出力较高时段灵活用电产热替代部分CHP热出力从而为风电腾出上网空间。场景C含电锅炉储热罐在场景B基础上增加储热罐实现热力系统的时段解耦。储热罐可以在热负荷低谷时段充热、高峰时段放热进一步提高系统灵活性。测试系统的具体参数配置如下基于常见测试系统参数调整而来参数数值备注CHP机组台数2抽汽式纯凝火电机组1作为系统补充风电场容量300 MW预测出力见数据文件电锅炉容量80 MW转换效率0.98储热罐容量200 MWh最大充放热功率50 MW热负荷峰值220 MWth典型冬季日曲线电负荷峰值400 MW典型冬季日曲线弃风惩罚300元/MWh约为煤耗成本的1.5倍风电预测数据和负荷数据需要根据实际地区典型日数据获取。如果没有现成数据可以自行构造一条带波动特征的曲线。我项目中的风电出力采用了一个带有早晚高峰特征、低谷出现在午间的典型冬季风电曲线这样设置的目的在于制造“供暖需求高、风电出力高、电负荷低”的时段充分考验模型的削峰填谷能力。4.2 结果对比与关键发现运行三种场景的调度模型后我得到以下关键结果数值经过脱敏处理弃风率对比场景A的弃风率高达18.7%大量夜间风电被白白丢弃。场景B加入电锅炉后弃风率降至9.2%降幅接近一半。场景C引入储热罐后弃风率进一步降至4.1%。这个趋势非常直观地验证了一个结论电转热路径的加入为风电消纳打开了通道而储热罐的出现又让这条通道的时间灵活性大幅提升。系统总成本对比场景A的总运行成本最高。场景B的成本相比场景A下降约6.5%主要来自弃风惩罚的减少和煤耗的降低。场景C的成本相比场景A下降约11.2%。值得注意的是场景B和场景C的成本差异主要来自弃风惩罚这说明在弃风惩罚系数设定合理的情况下储热罐的投资如果折算成日折旧成本可能刚好被省下的运行成本覆盖——这为项目经济性分析提供了直接参考。CHP机组运行状态对比场景A中CHP机组在整个供暖期几乎都以最小电出力运行但即便如此夜间电力过剩现象依然严重因为热负荷太高“以热定电”把电出力锁死了。场景C中CHP机组的电出力曲线出现了明显的“低谷抬升、高峰回降”特征——白天热负荷高峰时段储热罐放热支援热网CHP机组得以适当降低热出力、进而柔性调整电出力空间。这就是储热罐的时间平移效果。结果分析这一步我强烈建议在论文或报告中附上机组出力的Gantt图调度时序图直观展示各机组在不同时段的启停和出力变化。画图不需要太复杂Matlab的stairs命令画阶梯图就够用。4.3 敏感性分析惩罚系数与储热罐容量除了基本场景对比我还做了两个方向的敏感性分析目的是检验模型输出对关键参数的响应是否合理改变弃风惩罚系数从100元/MWh逐步增加到500元/MWh观察弃风率变化。结果表明惩罚系数在200-350元/MWh区间内弃风率下降最显著超过350元/MWh后弃风率变化趋于平缓说明此时系统内灵活性资源已经基本用尽继续提高惩罚只会增加成本而消纳效果有限这就是“边际收益递减”现象在调度模型中的体现。改变储热罐容量从50 MWh增加到500 MWh观察弃风率和总成本的变化。结果显示储热罐容量在150-250 MWh范围内性价比最高超过300 MWh后弃风率的改善非常微弱——原因在于此时制约风电消纳的因素已经从“储热容量不足”变成了“CHP机组电出力下限和电锅炉容量”继续增加储热罐容量没有意义。这类敏感性分析结果对实际工程中储热罐容量的规划设计非常有参考价值。5. 常见问题与调试经验实录5.1 Yalmip报错“Inconsistent constraints”排查这是最折磨人的报错之一意思是“约束不一致”通常是模型本身不可行。我的排查路径如下第一步检查是否存在“不可能满足的硬约束”。最典型的情况是某个时段的负荷值大于所有机组最大出力之和。用一个简单脚本检查数据load_max sum(P_max) wind_max - eb_min如果成立模型必然不可行。这类问题十次有八次是数据出了问题。第二步用Yalmip的diagnose命令定位问题。diagnose(Constraints, Objective, ops)会返回一堆诊断信息帮助你定位是哪些约束导致不可行。这个命令在求解器报“infeasible”时几乎是救命稻草。第三步检查变量维度和约束维度是否匹配。很多“inconsistent constraints”报错其实是因为某个约束的维度在循环中被拼接错了。比如约束矩阵是24×1而另一个约束是24×2两者拼在一起就报错了。5.2 求解速度过慢的问题如果模型规模不大比如10台机组以内、24个时段Gurobi或CPLEX应该在几秒内就能解完。如果求解时间超过1分钟大概率是模型里出现了一些拖累求解效率的结构一个典型的低效写法是在目标函数里使用了非线性函数如sin、exp、绝对值函数的错误用法导致问题从MILP退化成MINLP。Yalmip本身无法识别这些隐藏的非线性最终交给求解器时原本高效的线性求解器直接罢工或极慢。排查方法是在求解前用yalmip(clear)清空变量然后检查Objective中是否有abs、norm之外的非线性函数或者是否出现了P_chp^2但没有用二次函数形式声明的情况。另一个常见问题是整数变量过多。如果CHP机组启停变量都定义为binvar48个变量对Gurobi来说是小菜一碟但如果把机组启停的最小持续约束用大量implies描述可能会生成大量辅助二进制变量拖慢求解速度。我建议先不加最小启停时间约束跑一版确认其他部分正确后再加入这样可以快速定位问题来源。5.3 优化结果中风电出力总为零这类情况看着怪实际原因往往不在模型本身而在数据。风电出力预测曲线如果出现负值或大于装机容量的异常值优化器当然会“弃掉”这些不合理的风电。另一个可能是风电预测数据和火力发电数据的时间轴错位了——比如风电数据是每15分钟采样240个点而其他数据是每小时24个点直接拿去求模型优化器会被搞晕。我遇到过最隐蔽的一次是风电预测数据单位是kW其他数据单位是MW差了1000倍。但由于风电装机容量本身就小在总负荷里占比不大结果只是弃风率奇高但总成本变化不明显这个问题直到最后画图对比时才发现。所以正式求解前打印一下数据的统计特征最大值、最小值、平均值是特别值得养成的好习惯。5.4 代码可复现性的核心技巧最后分享一个小技巧整个项目的代码我建议从一开始就按“数据-模型-求解-后处理”四层分离来组织文件结构。比如用data_case1.m存放不同场景的参数文件build_model.m构建约束和目标函数run_scheduling.m是主程序入口plot_results.m负责可视化。这样每个文件职责单一后期修改测试场景时只需替换数据文件不需要动模型逻辑。这个习惯在我自己后续多次扩展模型比如加入碳捕集装置、加入需求响应时省掉了大量重复调试的时间。项目后续如果要往深了做还可以考虑把确定性日前调度扩展为两阶段随机优化考虑风电预测误差的分布、鲁棒优化考虑最坏场景下的调度方案或者加入碳交易机制、需求响应等约束维度。模型的扩展性好不好很大程度上取决于你当初写代码时把参数和结构分离得是否彻底——这也是我认为这个项目除模型本身之外最值得投入精力的地方。 电热综合能源系统日前经济调度与可再生能源消纳Matlab实现
1. 项目概述与核心思路拆解1.1 为什么是“电热综合”而不是单纯“电力”做电力系统优化调度的人大概率都遇到过这样一个尴尬场景电网侧算出来的消纳能力明明够结果现场还是弃风了。你翻来覆去查数据最后发现瓶颈根本不在电力系统本身而在供热系统——尤其是我国北方地区冬季供暖季热电联产机组“以热定电”的刚性约束直接把风电的上网空间挤得所剩无几。我最初接触这个项目时第一反应也是“这不就是个机组组合加经济调度吗”但真正把热力系统计入约束条件之后问题复杂度立刻不一样了。电热综合能源系统Integrated Electricity and Heat Energy System的核心特征是电力网络与热力网络通过热电联产机组CHP、电锅炉、储热罐等耦合设备深度绑定。电力系统需要实时平衡热力系统则可以凭借管网的蓄热特性和储热装置的缓冲在一段时间内“缓缓再算账”。这两个时间尺度差异巨大的系统放在一起做日前调度不是简单把两组约束拼起来就行而是要重新设计目标函数和决策变量的耦合关系。这个项目解决的核心痛点就是在供暖季高风电出力时段通过电锅炉和储热罐把多余的电力转化成热能储存起来让CHP机组可以降电出力而不必同步降热出力也就是通常说的“热电解耦”。这样一来风电可以多上网燃煤可以少烧一点系统总运行成本也随之下降。项目标题里“考虑可再生能源消纳”这几个字落脚点就在于此。1.2 日前经济调度的核心逻辑“日前”两个字指的是提前24小时或更长制定次日的调度计划。为什么一定要做日前因为火电机组的启停不是按开关那么简单锅炉升温、汽轮机冲转都需要时间如果等风来了再开机黄花菜都凉了。所以运行人员需要根据负荷预测和风电预测提前一天确定每台机组的启停状态、出力区间和热出力计划这就是日前经济调度的意义。经济调度的核心逻辑是在满足电力平衡、热力平衡、机组物理约束和安全约束的前提下让系统总运行成本最小化。换句大白话就是——用最便宜的方式把第二天的电和热都供上。这里“最便宜”不只是燃料成本还包括启停成本、弃风惩罚成本、购电成本等。具体到Matlab代码实现就是一个典型的混合整数线性规划MILP或混合整数二次规划MIQP问题决策变量里既有连续变量机组出力、热出力、储热罐充放热功率也有0-1整数变量机组启停状态。编写代码之前一定要把模型里所有时间断面之间的耦合关系理顺。日前调度是典型的多时段耦合优化火电机组有爬坡约束储热罐有容量变化约束CHP机组的热出力会影响下一个时段的电出力范围。如果没有先画清楚这些时序耦合关系就直接开始写Yalmip代码后面大概率会陷入“约束总是不对”的泥潭——这是我实际踩过的坑后面会详细展开。1.3 可再生能源消纳的建模思路可再生能源消纳的建模最常用的方式是引入弃风弃光惩罚项。为什么不直接约束风电必须全额消纳因为从数学优化角度硬性约束“风电出力等于预测值”会把所有灵活性都锁死系统可能因为风电预测误差而无法找到可行解。更合理的做法是允许调度模型“丢弃”部分风电但对弃风行为施加一个足够高的惩罚成本让优化器在权衡之下主动选择消纳风电。这个“足够高”很有讲究。惩罚系数太低优化器会发现弃风比调CHP机组更“划算”结果弃风量明显偏大惩罚系数太高又可能让模型过度牺牲经济性去追逐风电消纳导致总成本上升。我个人的经验是弃风惩罚系数设置在系统边际煤耗成本的1.5到2倍左右比较合适具体数值需要根据测试系统参数微调。另外模型中还需要对风电出力预测曲线做平滑处理。原始预测数据往往带有高频抖动而日前调度的时间粒度一般是1小时如果直接使用未处理的原始数据可能会让相邻时段的调度结果出现不合理跳变。我在项目中采用了简单的高斯滤波处理效果良好但这部分预处理代码容易被初学者忽略。2. 模型构建与关键技术点2.1 核心设备建模电热综合能源系统的设备模型是整个调度模型的地基。地基不稳后面做再多优化都白搭。我以这个项目中最典型的系统配置为例CHP机组、纯凝式火电机组、风电场、电锅炉、储热罐外加一个热负荷需求。CHP机组是全系统建模难度最高的部分。根据热电运行特性CHP机组分为背压式和抽汽式两类。背压式机组的电出力与热出力基本上是一一绑定的也就是说发多少电就必然产多少热几乎没有调节余地对系统灵活性贡献有限。抽汽式机组则灵活得多其可行运行域是一个四边形或更复杂的凸多边形区域。在Matlab中描述这个运行域需要把机组的最小电出力、最大电出力、最大进汽量、最小凝汽流量等参数转换成一组线性不等式约束。我项目中用的是抽汽式CHP它的电出力P和热出力H满足如下关系P_min Cv * H P P_max - Cv * H 范围随热出力变化 H 处于 [H_min, H_max]这里Cv是热电耦合系数也称为电出力随热出力变化的斜率单位是MW/MWth。具体数值每台机组不一样需要查厂家资料或参考标准测试数据。这个公式看着简单但约束方向容易写反我初次建模时就把不等式方向搞错了一次结果优化结果里CHP机组的电出力与热出力关系完全不符排查了好半天。电锅炉的建模比较直接输入电力输出热能转换效率在0.95到0.99之间出力上下限和爬坡约束可以参照电热转换特性设定。它的作用就是为风电提供一条“电转热”的消纳路径相当于给电力系统增加了一个可控的、时间灵活的用电负荷。储热罐是时间解耦的关键设备。它的核心约束是容量动态平衡方程S(t1) S(t) eta_charge * P_ch(t) * delta_t - P_dis(t) / eta_discharge * delta_t同时需要限制储热罐的储热容量上下限、单位时间充放热功率上限以及一天周期结束时的储热量需要回到初始值或允许在一定范围内偏离——这个约束常用于保证日前计划的可重复性。注意充热效率和放热效率通常不同而且有的简化模型会忽略效率直接用S(t1) S(t) (P_ch - P_dis) * delta_t这在工程分析中也可以接受但会高估储热罐的实际作用需要注意。纯凝火电机组的建模就是标准的机组组合模型包含出力上下限、爬坡约束、最小启停时间约束等。最小启停时间约束在日前调度中非常重要但在入门教程里经常被忽略导致模型结果在实际运行中根本无法执行——开机1小时就停机的情况现实中是不可能的汽轮机受不了这种折腾。2.2 目标函数设计与成本构成经济调度的目标函数本质上是把所有运行成本折算成“一天的总账单”然后让这个账单最小。这个项目中的成本构成可以分为以下几类燃煤成本CHP机组和纯凝火电机组的煤耗成本。通常用二次函数拟合机组煤耗特性C_fuel(P) a * P^2 b * P c在MILP模型中二次函数需要分段线性化处理。Yalmip自带pwlin或implies类似工具可以实现也可以用手动分段的方式。分段数一般取35段就能获得不错的近似精度。弃风惩罚成本模型允许弃风但是“弃风一时爽代价要付清”。惩罚系数设为W_pen单位是元/MWh乘以弃风量就是总惩罚成本。购电成本如果系统与上级电网存在联络线需要计入从外部购电的费用。有的模型中还需要考虑向外部售电的收入这需要根据售电电价是否低于购电电价来决定是否设非负约束避免“低买高卖”套利。启停成本机组从停机状态到开机状态需要额外消耗燃料和寿命这部分用固定成本项表示与是否发生启停动作0-1变量绑定。我建议在编写代码时先把目标函数的所有项在纸上列一个清单标注每项的量纲和系数来源然后再开始写Matlab代码。这个习惯在很大程度上避免了后面“优化结果总感觉哪里不对但又说不出原因”的困境。2.3 约束条件全景梳理约束条件是模型中最容易出错也最花时间的部分。我把这个项目中的约束条件按类型整理如下每一项都是代码中需要逐一实现并反复检查的内容电力平衡约束每个时段所有电源出力总和加上购电如果有减去电锅炉耗电等于该时段的电负荷需求。注意电锅炉在这里是“用电设备”在平衡方程中是以负荷形式出现的这一点的处理方式经常让初学者犯迷糊。热力平衡约束所有热源产热CHP机组热出力电锅炉热出力-储热罐充热功率储热罐放热功率等于该时段的热负荷需求。储热罐充热相当于增加热负荷放热相当于增加热源供给。机组出力约束包括电出力上下限、热出力上下限以及抽汽式CHP机组的热电耦合不等式。爬坡约束机组相邻时段的出力变化不能超过爬坡速率乘以时间间隔。这里的出力指的是电出力热出力是否参与爬坡约束需要根据模型复杂度决定。我项目的做法是电出力参与爬坡约束热出力单独设置热爬坡约束这样更贴近实际运行特性。储热罐约束容量上下限、充放热功率上下限、周期末储热量约束。风电出力约束风电实际出力介于0和预测值之间差值即为弃风量。系统备用约束可选为应对预测误差系统需要预留一定旋转备用容量。这个约束在日前调度中是否加入取决于研究需要。我项目中未加入备用约束以控制模型复杂度但如果你做的是工程化应用建议还是加上。所有约束条件都涉及“每个时段”的索引变量在Matlab中可以用循环方式构建也可以利用Yalmip的矩阵化表达式高效构建。循环写起来直观但速度慢矩阵化写起来高效但可读性差我的建议是先用循环跑通小规模测试再逐步优化代码结构。3. Matlab代码实现与实操要点3.1 工具箱选型与求解器配置Matlab做优化调度几乎绕不开Yalmip工具箱。它不是求解器而是“建模语言”作用是把人类易读的优化模型翻译成求解器能理解的标准形式。Yalmip最大优势在于你不需要关心求解器内部的数据结构只需要用sdpvar定义变量、用optimize命令求解它自动帮你做模型转换和求解器连接。求解器方面MILP问题推荐CPLEX或Gurobi这两个商业求解器的分支定界算法效率远超开源求解器。对于小型测试系统比如10台机组、24个时段Gurobi往往几秒钟就能求解CPLEX的求解速度也接近。如果你手里没有商业求解器的licenseYalmip也支持开源的CBC求解器求解速度慢一些但结果可靠。学校环境一般都能申请到Gurobi的免费学术license建议优先选择。安装Yalmip后第一件事是运行yalmiptest命令验证安装是否正常。这个问题我遇到不少求助案例多数是路径没设置好或者Matlab版本兼容性问题。Yalmip不需要安装直接把文件夹加入Matlab路径即可但注意不要和Matlab工具箱路径混淆。3.2 数据初始化与单位统一建模前必须先统一单位这是整个项目中最容易“翻车”的环节。电力系统常用的单位有MW和MWh但如果热力部分用GJ或Gcal那你在写约束时就得反复换算稍不留神就出乱子。我建议统一采用以下单位体系功率MW能量MWh热量MWth热功率能量用MWh热表示对于热负荷和储热罐容量如果原始数据给的是GJ换算关系是1 MWh 3.6 GJ。如果你拿到的是Gcal那就更麻烦一些1 Gcal约等于1.163 MWh。强烈建议写一个单位换算函数在读取原始数据后统一转换为标准单位制然后再送入模型。这一步做扎实了后面调试会省很多时间。数据初始化部分的代码我习惯用一个结构体struct来管理所有系统参数。比如%% 系统参数定义 sys struct(); sys.T 24; % 调度时段数日前24小时 sys.dt 1; % 时间间隔小时 sys.n_chp 2; % CHP机组台数 sys.n_thermal 1; % 纯凝火电机组台数 %% CHP机组参数 sys.chp(1).P_min 50; % 最小电出力 (MW) sys.chp(1).P_max 200; % 最大电出力 (MW) sys.chp(1).H_min 20; % 最小热出力 (MWth) sys.chp(1).H_max 150; % 最大热出力 (MWth) sys.chp(1).Cv 0.15; % 热电耦合系数 sys.chp(1).a 0.0002; % 煤耗曲线二次项 sys.chp(1).b 0.3; % 煤耗曲线一次项 sys.chp(1).c 20; % 煤耗曲线常数项这样做的最大好处是后续构建约束时可以直接用sys.chp(1).P_min这种可读性极高的字段名而不是一堆散落的变量名。当你代码写到300行以上时这种数据组织方式的优势会非常明显。3.3 决策变量定义与约束构建Yalmip中的变量定义是建模的第一步。在这个项目中需要以下几种类型的变量%% 决策变量定义 P_chp sdpvar(sys.n_chp, sys.T, full); % CHP电出力 H_chp sdpvar(sys.n_chp, sys.T, full); % CHP热出力 P_th sdpvar(sys.n_thermal, sys.T, full); % 火电机组电出力 P_wind sdpvar(1, sys.T, full); % 风电实际出力 P_wind_curtail sdpvar(1, sys.T, full); % 弃风功率 P_eb sdpvar(1, sys.T, full); % 电锅炉耗电功率 H_eb sdpvar(1, sys.T, full); % 电锅炉产热功率 P_ch_st sdpvar(1, sys.T, full); % 储热罐充热功率 P_dis_st sdpvar(1, sys.T, full); % 储热罐放热功率 S_st sdpvar(1, sys.T 1, full); % 储热罐储热量含初始时刻 u_chp binvar(sys.n_chp, sys.T); % CHP机组启停状态可选我特别说明一下S_st为什么定义到T1因为储热罐的容量变化是一个跨时段递推过程我们需要保存初始时刻时段1开始前和每个时段结束时的储热量这样就能方便地写出容量平衡约束。很多初学者只定义T个变量结果递推约束写到最后发现少了一个边界值处理起来非常痛苦。约束构建部分最具代表性的几个约束代码如下%% 电力平衡约束 Constraints []; for t 1:sys.T Constraints [Constraints, sum(P_chp(:, t)) sum(P_th(:, t)) P_wind(:, t) ... P_load(t) P_eb(t)]; end %% 热力平衡约束 for t 1:sys.T Constraints [Constraints, sum(H_chp(:, t)) H_eb(t) P_dis_st(t) - P_ch_st(t) ... H_load(t)]; end %% CHP机组运行域约束 for k 1:sys.n_chp for t 1:sys.T Constraints [Constraints, P_chp(k, t) sys.chp(k).P_min sys.chp(k).Cv * H_chp(k, t)]; Constraints [Constraints, P_chp(k, t) sys.chp(k).P_max - sys.chp(k).Cv * H_chp(k, t)]; Constraints [Constraints, H_chp(k, t) sys.chp(k).H_min, H_chp(k, t) sys.chp(k).H_max]; end end %% 储热罐容量动态平衡 for t 1:sys.T Constraints [Constraints, S_st(t 1) S_st(t) sys.eta_ch * P_ch_st(t) * sys.dt - ... P_dis_st(t) / sys.eta_dis * sys.dt]; end这里需要注意Yalmip约束拼接使用[Constraints, new_constraint]的语法而且约束之间不要用逗号之外的分隔符。很多初学者在这里会用分号导致约束维数错乱报错信息还特别不直观。3.4 求解流程与结果输出求解部分的代码相对简单%% 目标函数 Objective 0; % 煤耗成本分段线性近似 for k 1:sys.n_chp Objective Objective sum(sys.chp(k).a * P_chp(k, :) .^ 2 ... sys.chp(k).b * P_chp(k, :) sys.chp(k).c); end % 弃风惩罚 Objective Objective sys.w_pen * sum(P_wind_curtail); % 电锅炉运行成本维护成本一般较小 Objective Objective sys.c_eb * sum(P_eb); %% 求解 ops sdpsettings(solver, gurobi, verbose, 2, showprogress, 1); sol optimize(Constraints, Objective, ops); %% 结果提取与验证 if sol.problem 0 P_chp_opt value(P_chp); H_chp_opt value(H_chp); P_wind_opt value(P_wind); % ... 其他变量提取 else disp(求解失败错误代码: string(sol.problem)); diagnose(Constraints, Objective, ops); end关于二次项的处理我在这里用的是直接二次函数。Yalmip配合Gurobi可以直接处理MIQP问题不需要手动分段线性化。但如果使用CBC求解器就必须先把二次项分段线性化否则求解器会报错。这是一个常见的“代码在不同求解器之间移植”的坑值得留意。结果输出部分我建议不要只输出数值最好同时生成可视化的调度曲线图。Matlab的plot命令足够应付画出电出力曲线、热出力曲线、储热罐容量变化、风电消纳情况等关键图一眼就能看出调度结果合不合理。4. 典型场景仿真与结果分析4.1 场景设计与参数配置为了验证模型的有效性我设计了三个典型场景分别对应不同的系统配置形成对比分析场景A基准场景不含电锅炉不含储热罐。这是最传统的CHP系统配置也是弃风问题最严重的配置。在这个场景下CHP机组只能通过降低热出力来降低电出力因此供暖季的风电消纳空间非常有限。场景B含电锅炉引入电锅炉作为电转热的灵活性资源但无储热罐。电锅炉可以在风电出力较高时段灵活用电产热替代部分CHP热出力从而为风电腾出上网空间。场景C含电锅炉储热罐在场景B基础上增加储热罐实现热力系统的时段解耦。储热罐可以在热负荷低谷时段充热、高峰时段放热进一步提高系统灵活性。测试系统的具体参数配置如下基于常见测试系统参数调整而来参数数值备注CHP机组台数2抽汽式纯凝火电机组1作为系统补充风电场容量300 MW预测出力见数据文件电锅炉容量80 MW转换效率0.98储热罐容量200 MWh最大充放热功率50 MW热负荷峰值220 MWth典型冬季日曲线电负荷峰值400 MW典型冬季日曲线弃风惩罚300元/MWh约为煤耗成本的1.5倍风电预测数据和负荷数据需要根据实际地区典型日数据获取。如果没有现成数据可以自行构造一条带波动特征的曲线。我项目中的风电出力采用了一个带有早晚高峰特征、低谷出现在午间的典型冬季风电曲线这样设置的目的在于制造“供暖需求高、风电出力高、电负荷低”的时段充分考验模型的削峰填谷能力。4.2 结果对比与关键发现运行三种场景的调度模型后我得到以下关键结果数值经过脱敏处理弃风率对比场景A的弃风率高达18.7%大量夜间风电被白白丢弃。场景B加入电锅炉后弃风率降至9.2%降幅接近一半。场景C引入储热罐后弃风率进一步降至4.1%。这个趋势非常直观地验证了一个结论电转热路径的加入为风电消纳打开了通道而储热罐的出现又让这条通道的时间灵活性大幅提升。系统总成本对比场景A的总运行成本最高。场景B的成本相比场景A下降约6.5%主要来自弃风惩罚的减少和煤耗的降低。场景C的成本相比场景A下降约11.2%。值得注意的是场景B和场景C的成本差异主要来自弃风惩罚这说明在弃风惩罚系数设定合理的情况下储热罐的投资如果折算成日折旧成本可能刚好被省下的运行成本覆盖——这为项目经济性分析提供了直接参考。CHP机组运行状态对比场景A中CHP机组在整个供暖期几乎都以最小电出力运行但即便如此夜间电力过剩现象依然严重因为热负荷太高“以热定电”把电出力锁死了。场景C中CHP机组的电出力曲线出现了明显的“低谷抬升、高峰回降”特征——白天热负荷高峰时段储热罐放热支援热网CHP机组得以适当降低热出力、进而柔性调整电出力空间。这就是储热罐的时间平移效果。结果分析这一步我强烈建议在论文或报告中附上机组出力的Gantt图调度时序图直观展示各机组在不同时段的启停和出力变化。画图不需要太复杂Matlab的stairs命令画阶梯图就够用。4.3 敏感性分析惩罚系数与储热罐容量除了基本场景对比我还做了两个方向的敏感性分析目的是检验模型输出对关键参数的响应是否合理改变弃风惩罚系数从100元/MWh逐步增加到500元/MWh观察弃风率变化。结果表明惩罚系数在200-350元/MWh区间内弃风率下降最显著超过350元/MWh后弃风率变化趋于平缓说明此时系统内灵活性资源已经基本用尽继续提高惩罚只会增加成本而消纳效果有限这就是“边际收益递减”现象在调度模型中的体现。改变储热罐容量从50 MWh增加到500 MWh观察弃风率和总成本的变化。结果显示储热罐容量在150-250 MWh范围内性价比最高超过300 MWh后弃风率的改善非常微弱——原因在于此时制约风电消纳的因素已经从“储热容量不足”变成了“CHP机组电出力下限和电锅炉容量”继续增加储热罐容量没有意义。这类敏感性分析结果对实际工程中储热罐容量的规划设计非常有参考价值。5. 常见问题与调试经验实录5.1 Yalmip报错“Inconsistent constraints”排查这是最折磨人的报错之一意思是“约束不一致”通常是模型本身不可行。我的排查路径如下第一步检查是否存在“不可能满足的硬约束”。最典型的情况是某个时段的负荷值大于所有机组最大出力之和。用一个简单脚本检查数据load_max sum(P_max) wind_max - eb_min如果成立模型必然不可行。这类问题十次有八次是数据出了问题。第二步用Yalmip的diagnose命令定位问题。diagnose(Constraints, Objective, ops)会返回一堆诊断信息帮助你定位是哪些约束导致不可行。这个命令在求解器报“infeasible”时几乎是救命稻草。第三步检查变量维度和约束维度是否匹配。很多“inconsistent constraints”报错其实是因为某个约束的维度在循环中被拼接错了。比如约束矩阵是24×1而另一个约束是24×2两者拼在一起就报错了。5.2 求解速度过慢的问题如果模型规模不大比如10台机组以内、24个时段Gurobi或CPLEX应该在几秒内就能解完。如果求解时间超过1分钟大概率是模型里出现了一些拖累求解效率的结构一个典型的低效写法是在目标函数里使用了非线性函数如sin、exp、绝对值函数的错误用法导致问题从MILP退化成MINLP。Yalmip本身无法识别这些隐藏的非线性最终交给求解器时原本高效的线性求解器直接罢工或极慢。排查方法是在求解前用yalmip(clear)清空变量然后检查Objective中是否有abs、norm之外的非线性函数或者是否出现了P_chp^2但没有用二次函数形式声明的情况。另一个常见问题是整数变量过多。如果CHP机组启停变量都定义为binvar48个变量对Gurobi来说是小菜一碟但如果把机组启停的最小持续约束用大量implies描述可能会生成大量辅助二进制变量拖慢求解速度。我建议先不加最小启停时间约束跑一版确认其他部分正确后再加入这样可以快速定位问题来源。5.3 优化结果中风电出力总为零这类情况看着怪实际原因往往不在模型本身而在数据。风电出力预测曲线如果出现负值或大于装机容量的异常值优化器当然会“弃掉”这些不合理的风电。另一个可能是风电预测数据和火力发电数据的时间轴错位了——比如风电数据是每15分钟采样240个点而其他数据是每小时24个点直接拿去求模型优化器会被搞晕。我遇到过最隐蔽的一次是风电预测数据单位是kW其他数据单位是MW差了1000倍。但由于风电装机容量本身就小在总负荷里占比不大结果只是弃风率奇高但总成本变化不明显这个问题直到最后画图对比时才发现。所以正式求解前打印一下数据的统计特征最大值、最小值、平均值是特别值得养成的好习惯。5.4 代码可复现性的核心技巧最后分享一个小技巧整个项目的代码我建议从一开始就按“数据-模型-求解-后处理”四层分离来组织文件结构。比如用data_case1.m存放不同场景的参数文件build_model.m构建约束和目标函数run_scheduling.m是主程序入口plot_results.m负责可视化。这样每个文件职责单一后期修改测试场景时只需替换数据文件不需要动模型逻辑。这个习惯在我自己后续多次扩展模型比如加入碳捕集装置、加入需求响应时省掉了大量重复调试的时间。项目后续如果要往深了做还可以考虑把确定性日前调度扩展为两阶段随机优化考虑风电预测误差的分布、鲁棒优化考虑最坏场景下的调度方案或者加入碳交易机制、需求响应等约束维度。模型的扩展性好不好很大程度上取决于你当初写代码时把参数和结构分离得是否彻底——这也是我认为这个项目除模型本身之外最值得投入精力的地方。