前几个月在做一个综合能源系统调度方向的课题,核心就是标题里这个模型:考虑可再生能源消纳的电热综合能源系统日前经济调度模型,并且用Matlab写了一套可跑的代码。这个课题的典型场景是冬季供暖期,热电联产机组为了保供热必须压着最小电出力运行,导致夜间大风时段的风电被大量弃掉。单纯在电网侧做调度很难解决这个问题,必须把热网的灵活性也拉进来一起优化。这篇博文就把我踩过的坑、模型里每个约束到底在干什么、以及Matlab里怎么落地写代码,完整梳理一遍。适合正在做综合能源优化调度、电力系统经济运行或者想用Matlab复现类似模型的学生和工程师参考。
1. 项目背景与问题定义
1.1 为什么是“电热综合能源系统”
传统电力调度只考虑电网,热力系统归热力公司管,两者基本是各调各的。但城市里大量供暖来自热电联产机组(CHP),这种机组同时产电和产热,电出力会受到热出力的硬约束。冬季为了满足热负荷,CHP往往被迫维持较高的电出力,给风电留出的消纳空间被严重压缩。如果只做电力系统的经济调度,CHP的最小电出力就得当作固定边界来处理,算出来的结果当然会牺牲风电。
电热综合能源系统的核心思路,是把电网、热网以及电锅炉、储热罐、电储能这类耦合设备放进同一个优化问题里。调度中心在安排次日各机组出力时,同时决定电功率和热功率怎么分配、储热罐什么时候充放热、电锅炉要不要启动。这样一来,热力系统就变成了电力系统的灵活性资源,而不是僵硬的负荷。这也是近年来“热电解耦”概念被反复提起的原因。
1.2 可再生能源消纳的核心矛盾
可再生能源消纳问题最突出的时间点通常不是用电高峰,反而是负荷低谷加上大风、光照强的时段。电网没办法凭空消耗那么多电,火电和CHP又都有最小技术出力限制,如果热负荷把CHP的出力“焊死”在一个高值上,那风电就只能弃掉。从调度模型的角度看,消纳能力取决于系统向下调节能力有多强:常规火电机组能不能少发一点,CHP能不能通过储热来替代一部分供热出力,储能能不能多充电,这些调节手段的数学描述,最终都会变成约束条件和目标函数里的成本项。
这个模型要做的,就是在保证电热平衡的前提下,把总运行成本降到最低。可再生能源发电的边际成本接近零,所以目标函数里只要给弃风弃光加一个惩罚成本,优化器就会自动倾向于减少弃电。但这个惩罚成本不能随便设,设得太高相当于逼着系统用其他手段强行消纳,设得太低又起不到效果,实际调试的时候需要结合弃电率目标来标定。
1.3 日前经济调度的目标与场景
日前调度的意思是提前一天制定未来24小时的调度计划,时间分辨率通常取1小时,也有取15分钟的。它和实时调度不一样,日前计划更侧重经济性,为的是给次日发电机组、热电联产设备、储能设备一个可执行的基准。实际运行中,风电出力预测、负荷预测都会有误差,但这些误差可以放到日内滚动优化里去修正,日前模型就把预测值当成已知输入来处理。
我复现的场景里,系统包括一个风电场、一个常规火电机组、一个热电联产机组、一个电锅炉、一个储热罐和一个电储能装置。电网与大电网相连,可以购电也可以向大电网售电。热网侧只考虑供热功率平衡,不做复杂的管网水力计算,这样既能抓住问题的本质,又不会把建模难度抬得太高。对于刚刚接触这个方向的人来说,这是一个既能说明原理又足够上手的切入点。
2. 日前经济调度模型的核心设计
2.1 系统结构:电网、热网、耦合设备
整个系统可以看成两个能量母线:电母线和热母线。电母线上连接着风电、火电、热电联产机组的电输出、电储能、电锅炉、固定电负荷以及与大电网的交换功率。热母线上连接着热电联产机组的热输出、电锅炉的热输出、储热罐的充放热以及固定热负荷。电锅炉和热电联产机组是电热耦合的核心,储热罐则是热侧的时间转移装置。
在建模时,电网和热网之间的耦合主要靠CHP的电热运行可行域来体现。CHP机组的电出力P和热出力H之间不是独立的,通常被限制在由多个顶点构成的多边形区域内。我在代码里用的是线性不等式约束来表达这个多边形,这在混合整数线性规划(MILP)里是很标准的做法。储热罐的储能状态用热量的单位(比如MWh)来描述,充放热效率、容量上限和最大充放热功率都要考虑。
2.2 目标函数:经济性指标怎么量化
目标函数要覆盖一整天的总成本,每一项都有一个明确的物理或经济含义。我的模型里把总成本拆成了五个部分:常规火电机组的煤耗成本、热电联产机组的燃料成本、与大电网交换功率的成本、弃风弃光的惩罚成本、储能设备的充放电损耗成本。煤耗成本用一个二次函数近似,也可以用分段线性函数来逼近,然后转成MILP。
具体来说,火电成本是C_C = a_C * P_C(t)^2 + b_C * P_C(t) + c_C,CHP成本类似,但需要用机组当前的电、热出力共同计算。购电成本是电价乘以从电网输入的功率,如果是向电网倒送电,则按上网电价计算收益。弃风惩罚成本是弃风功率乘以单位惩罚因子,弃光同理。目标函数写成min sum_t (所有成本项),这就是日前经济调度的“经济”二字所在。
我在代码里没有把储能损耗单独列项,而是通过充放电效率来隐式表示,因为电池充进去1MWh能放出来不到1MWh,这个能量损失已经体现在功率平衡等式里了。如果你想让模型更精细,可以考虑增加一个设备运行维护成本项,但初版模型不必过度复杂,先把主框架跑通再逐步加细节。
2.3 约束条件:功率平衡、设备出力、爬坡、储能
约束条件按性质分层。第一层是系统平衡约束,也就是每个时段电功率平衡和热功率平衡。电平衡:风电出力+光伏出力+火电出力+CHP电出力+电储能放电-电储能充电+购电功率-售电功率=电负荷+电锅炉耗电。热平衡同理。这是所有调度模型必须满足的硬约束。
第二层是设备自身约束。火电机组和CHP机组有出力上下限,CHP的电出力上限和热出力上限还要满足可行域约束。常规火电还有爬坡约束:相邻两个时段出力的变化量不能超过爬坡速率乘以时段长度。电锅炉的耗电功率有上限,储热罐有容量上限和最大充放热功率限制。第三层是储能的动态约束,也就是荷电状态(SOC)的递推关系:下一时段SOC等于当前SOC加上充热/充电量再减去放热/放电量,还要考虑效率。最后还可以加一个调度周期始末SOC相等的约束,保证储能设备的日循环性,这样第二天调度不会把昨天存下的能量都偷走。
3. 可再生能源消纳的关键建模细节
3.1 风电/光伏出力建模与预测误差
日前调度中风电出力和光伏出力通常作为确定性预测值输入,但为了研究消纳问题,我们往往会在预测值基础上乘以不同渗透率或者叠加一个预测误差场景。最简单的做法是给定一个24小时的风电预测曲线,风电实际出力等于预测值减去弃风量。这样弃风量作为一个非负决策变量进入模型,优化器会权衡弃风惩罚和其他调节成本。
光伏的特征是白天有出力、夜间为零,和负荷曲线有一定的相关性。在北方冬季,光伏对消纳的贡献有限,但模型里可以一并考虑。如果想进一步做不确定性研究,可以把风电出力的多个随机场景作为输入,采用随机规划;如果不想增加太多求解负担,可以只做确定性模型,再在结果分析中对预测误差做灵敏度测试。我在复现时先用确定性模型跑通,后续再改成场景法,这样调试起来逻辑更清晰。
3.2 弃风弃光惩罚与消纳率约束
很多文献里会用弃风率作为评价指标,但这个指标在优化模型里不是一个天然的量,需要通过惩罚项或约束来体现。我给弃风功率设定了一个较高的单位惩罚成本,这样优化目标会让各机组出力组合尽量给风电让路。惩罚成本的数值通常设为基准电价的1到2倍,但具体需要试算:如果太小,结果里弃风量依然很大;如果太大,可能会让储能设备出现不必要的深充深放,反而增加成本。
另一种做法是显式加入可再生能源消纳率约束,比如要求一天总弃风率不超过5%。这个约束可以用弃风总量除以风电预测总量来写,但因为这是一个带有分数线形式的表达式,线性化时要把分母乘过去:弃风总量 <= 5% * 风电预测总量。这个约束本身就是线性的,非常好处理。实际工程里两种方式都有用,惩罚项更平滑,消纳率约束更硬性。
3.3 热电联产机组“以热定电”的破解思路
传统热电联产机组遵循“以热定电”,就是热负荷确定了,电出力就被限制在一个比较窄的范围内。为了打破这个限制,工程上常见做法是加装储热罐或者电锅炉。储热罐可以在热负荷低谷时吸收CHP的富余热出力,让CHP在电负荷低谷时段降低电出力;电锅炉则是在弃风严重时直接用电来产热,替代一部分CHP供热。
在模型里,储热罐和电锅炉的加入改变了热平衡约束。原本CHP热出力必须等于热负荷,现在变成CHP热出力+电锅炉热出力+储热罐放热-储热罐充热=热负荷。这样一来,CHP的热出力可以低于热负荷,多出的热负荷由电锅炉在弃风时段承担,等于用风电替代了CHP的一部分燃煤供热。这就是热电解耦在数学模型中的具体体现。我复现的结果显示,加入合适的储热容量并配置电锅炉后,弃风率可以下降不少,虽然总成本未必大幅下降,但可再生能源利用率明显改善。
4. Matlab代码实现与求解方案
4.1 整体代码框架与文件结构
我习惯把Matlab代码拆成多个脚本和函数,避免所有内容堆在同一个文件里。主程序叫run_day_ahead.m,负责设置参数、读取数据、调用模型构建函数、求解并展示结果。数据文件用Excel或.mat格式,保存24小时的负荷、风电预测出力、电价等序列。核心的模型构建部分放在一个函数里,比如build_model(params, data),返回一个包含决策变量、约束和目标函数的Yalmip结构体。
整个文件结构大致是:
run_day_ahead.m:主程序,设置路径、调用流程。load_data.m:读取基础数据,做单位换算。build_model.m:用Yalmip定义决策变量、约束和目标函数。solve_model.m:调用优化求解器,检查求解状态。plot_results.m:绘制结果图表。calc_metrics.m:计算总成本、弃风率等指标。
这样的好处是,改数据不需要动模型,改模型不需要动绘图。特别是当你想对比不同场景(有储热和无储热)时,只需要循环调用模型构建函数并传不同的参数。我第一次做的时候把所有代码堆在一个脚本里,改一个约束要来回滚动,非常痛苦,后来才重构。
4.2 数据准备与参数设置
数据准备是Matlab代码里最容易被忽略但最容易出错的部分。首先是时间序列数据,包括24个小时的电负荷、热负荷、风电出力系数和光伏出力系数。风电出力系数乘以风电装机容量就是预测出力。电负荷单位统一用MW,热负荷单位我用MWth,两者虽然物理维度都是功率,但一定要在变量名里标注清楚,避免混淆。
参数设置我建议写成结构体params,比如params.gen_Pmax、params.chp_area、params.ess_eta这些字段,每个字段都有注释。这样做的好处是,调用params的时候代码可读性会高很多,而且不容易写错数值。参数的来源要注明,比如从参考文献、实际机组手册或者典型日数据里来。我复现时用的数据一部分参考了某北方园区微网的数据,CHP参数采用典型的30MW抽凝式机组参数,电锅炉容量设为10MW,储热罐容量设为50MWh。
4.3 基于Yalmip建模与求解器配置
Matlab里建优化模型,强烈建议用Yalmip工具箱再加Gurobi或CPLEX求解器。Yalmip写约束非常直观,接近数学表达式;Gurobi和CPLEX在求解大规模线性规划/MILP时性能远好于Matlab自带的linprog和intlinprog。Yalmip的设置很简单,下载后添加到Matlab路径即可,求解器需要单独安装。
下面是一段核心建模代码的骨架,主要展示决策变量和约束的写法:
%% 决策变量定义 P_gen = sdpvar(T, 1); % 火电出力 P_chp = sdpvar(T, 1); % CHP电出力 H_chp = sdpvar(T, 1); % CHP热出力 P_wind = sdpvar(T, 1); % 风电上网出力 P_wind_curt = sdpvar(T, 1); % 弃风功率 P_eb = sdpvar(T, 1); % 电锅炉耗电 H_eb = sdpvar(T, 1); % 电锅炉产热 soc_hs = sdpvar(T+1, 1); % 储热罐SOC P_buy = sdpvar(T, 1); % 购电功率 P_sell = sdpvar(T, 1); % 售电功率 %% 变量定义结束后,写约束 Constraints = []; % 电功率平衡 for t = 1:T Constraints = [Constraints, ... P_gen(t) + P_chp(t) + P_wind(t) + ess_dis(t) - ess_chg(t) + P_buy(t) - P_sell(t) ... == P_load(t) + P_eb(t)]; end % CHP可行域约束(简单多边形顶点组合) for t = 1:T Constraints = [Constraints, ... P_chp(t) >= params.chp_Pmin_whenH0 + params.chp_k * H_chp(t)]; Constraints = [Constraints, ... P_chp(t) <= params.chp_Pmax - params.chp_k2 * H_chp(t)]; end % 储热罐动态 for t = 1:T Constraints = [Constraints, ... soc_hs(t+1) == soc_hs(t) + params.hs_eta_ch * H_chp_heat_hs(t) - H_hs_dis(t)]; end % 最终调用求解 optimize(Constraints, Objective, sdpsettings('solver','gurobi','verbose',2));注意上面的ess_dis、ess_chg、H_chp_heat_hs等变量我还没有定义完整,实际代码里都需要补全。这里想强调的是,构造约束时尽量用循环而不是一堆重复代码,虽然性能差别不大,但可维护性好很多。
4.4 结果可视化与指标计算
求解完成后,第一件事是检查optimize的返回状态和求解器输出。Yalmip中可以用problem = optimize(...),如果problem ~= 0,说明模型出了问题。结果可视化方面,我一般画四张图:第一张是电功率平衡堆叠图,第二张是热功率平衡堆叠图,第三张是储热罐SOC和电锅炉出力曲线,第四张是机组出力与弃风曲线。Matlab里用area函数画堆叠图最方便,比如:
figure; area(1:T, [P_gen P_chp P_wind P_buy], 'LineWidth', 1.5); hold on; plot(1:T, P_load + P_eb, 'k--', 'LineWidth', 2); legend('火电','CHP','风电','购电','负荷+电锅炉');指标计算包括总燃料成本、购电成本、弃风率、消纳率、CHP热出力占比等。这些指标可以用来横向比较不同方案,比如“有储热”和“无储热”两种结果。我计算弃风率用的是sum(P_wind_curt)/sum(P_wind_forecast),注意分子分母求和,不要用逐时比例再平均,后者会跟总弃风率不一致。
5. 实操过程中的常见问题与调试经验
5.1 求解器不可用或许可证报错
Matlab提示找不到求解器是新手最常遇到的问题。Yalmip内部自带的求解器只能处理一些小规模线性规划,像我们这种几十个变量、上百个约束的模型,效率很受影响。推荐装Gurobi或CPLEX。Gurobi给学术界提供免费许可证,申请的时候用学校邮箱即可,下载安装包后要把gurobi的Matlab接口路径添加到Matlab环境里。
如果不想用商业求解器,Matlab自带的intlinprog也能解MILP,但Yalmip对它的兼容性不如Gurobi。一个偷懒的办法是把模型里所有整数变量去掉,改成纯线性规划,用linprog就能跑。不过大多数调度模型里设备启停状态是整数变量,纯LP会失去启停逻辑,结果不一定合理。所以我个人建议还是花点时间装Gurobi,后面调到更复杂的模型也顺手。
5.2 模型不可行:怎么定位冲突的约束
模型不可行是所有优化建模者都会遇到的坎。前几天刚调一个带储热罐的模型,怎么求解都返回不可行,后来发现是储热罐的SOC递推关系里,同一个时段既规定了初始容量又要满足最大充放热功率,两个约束加在一起超出了物理极限。定位不可行约束有个通用方法:先把模型简化,例如只保留电平衡和机组上下限,跑通了再加热平衡,然后加储能,逐步引入新约束。这样一旦某步不可行,就知道问题出在刚加的那组约束里。
Yalmip本身不直接给出不可行约束列表,但Gurobi求解器在调试模式下有IIS(Irreducible Inconsistent Subsystem)功能,可以返回一组导致不可行的最小约束集合。在Yalmip里可以通过sdpsettings('debug',1)进行诊断,或者直接让Gurobi输出IIS,不过这个功能需要额外调用Gurobi的Matlab命令。更实用的笨办法是逐时段检查数据:比如某个时段热负荷大于所有供热设备的最大供热能力,那这个时段必不可行,检查数据合理性往往能快速发现问题。
5.3 非线性项与求解慢的线性化技巧
经济调度里最典型的非线性项是火电煤耗二次函数。二次规划(QP)用Gurobi能直接求解,但如果你的模型里还有其他整数变量,就变成MIQP,求解速度明显变慢。工程上常用分段线性化把二次项近似成多段线性函数,这样模型变为MILP,求解速度和稳定性都好很多。分段线性化在Yalmip里可以手写,用一个整数变量和一个连续变量表示每一段的状态,也可以用binvar加约束实现,但手写时要注意不要引入不等式方向错误。
CHP机组的电热可行域通常是非凸多边形,线性化思路是多边形顶点。如果引入0/1变量表示CHP运行状态,可行域可以拆成多个三角形,每个三角形对应一组线性不等式,再用大M法连接。这个技巧在文献里很常见,但实现起来需要细心。如果你的模型允许CHP始终在线,不启停,就可以直接把整个可行域用外部凸包线性不等式表示,省掉整数变量,速度大幅提升。
5.4 数据量纲和时段对齐问题
最后说一个特别容易犯的错:数据量纲。风电出力如果给了千瓦,负荷用的是兆瓦,目标函数里数值差了1000倍,优化器会把大数项优先优化掉,结果完全失真。我处理的方法是在load_data里统一除以一个基准值变成MW,热功率统一用MWth,时间统一用小时。这样所有的成本系数、功率、能量都在一个量纲体系里,结果才有意义。
还有时段对齐问题:储热罐SOC有T+1个状态,充电功率有T个值。如果你给SOC初始值赋值到了soc_hs(1),但约束里写的是soc_hs(2) == soc_hs(1) + ...,没问题;但如果循环里把时间向量和SOC索引搞错一位,模型就会在第一个时段把初始SOC当作一个大变量来处理,结果完全不符合直觉。我每次都会用size()检查所有变量的维度,尽量在建模前就打印出来核对一遍。
6. 案例复现与效果分析
6.1 典型日数据与场景设置
为了验证模型能真实反映消纳问题,我构造了一个冬季典型日场景:夜间0点到6点是电负荷低谷,风电出力较高,热负荷也处于全天最高时段。全天电负荷峰值出现在白天10点到12点以及晚上19点到21点。火电机组最大出力50MW,CHP最大电出力30MW,风电场容量40MW。热负荷峰值在夜间约30MWth。
在这个场景下,如果不做电热协调,CHP在夜间必须保持较高热出力,电出力至少在15MW以上,风电预测出力可能达到30MW,电负荷只有20MW左右,除非把多余电功率外送否则就会弃风。但外送通道受限,购售电价也有峰谷差。正好可以测试模型能不能通过储热罐在夜间多充热、让CHP降热出力,从而多给风电让路。
6.2 调度结果与经济性分析
跑完模型后,我最关心的三个数字是:总成本、弃风率、CHP热电比变化。无储热、无电锅炉的基准场景里,弃风率在15%左右,总成本较高。加入电锅炉和储热罐后,模型会在夜间增加电锅炉耗电,同时让CHP降低热出力,储热罐储存多余热量供早高峰使用。这样夜间的风电上网空间增大,白天的电锅炉关闭,储热罐放热满足热负荷。最终弃风率降到5%以内,总成本因为惩罚成本下降而小幅下降,虽然电锅炉增加了设备损耗,但同时减少CHP燃料消耗。
具体数字会随数据浮动,但规律非常稳定:储热容量越大,弃风率下降越明显,但超过一定容量后边际收益锐减。电锅炉功率增加也有类似效果。这说明调度模型不仅算出了最优值,还能帮助我们做容量配置分析,这是这个项目最值钱的地方。
6.3 可再生能源消纳提升效果对比
我做了三组对比:基准无耦合、仅加储热、加储热电锅炉。三组模型的目标函数形式和约束不完全一样,但框架完全复用。对比后可以清楚看到,仅加储热就能把弃风率从15%降到9%,但电锅炉投入后进一步降到4%。这个对比过程用代码实现没有多少难度,难的是怎么解读结果:储热罐通过时间平移消纳弃风,电锅炉通过能量转换消纳弃风,两种手段的成本结构不同。比如储热罐主要在夜间充热、白天放热,燃料替代效果有限;电锅炉直接消耗多余风电产生热,替代的是CHP的燃气消耗,燃料成本下降更直接。
在做论文或项目报告时,这类对比分析是很有说服力的。你甚至可以加一个参数扫描循环,把电锅炉容量从0扫到20MW,画出弃风率随容量变化的曲线,这样一眼就能看出最佳配置区间。代码层面只需要在4.2节的基础上增加一个for循环,每次修改params.eb_Pmax,重新求解并记录结果。
最后分享一点体会:这种电热联合调度的模型,难点通常不在数学表述,而在工程约束的取舍。你可以在模型里把热网建模得非常精细,比如加入供热管道温度动态和节点流量平衡,但那样会让模型规模和求解难度都急剧上升;对于一个从零开始的项目,先用简化的热平衡把概念验证打通,比一上来就追求全面细致要重要得多。把日前经济调度模型跑通之后,后面扩展碳交易成本、需求响应、风光储联合参与现货市场都是顺路的事。模型框架本身,尤其是Matlab里的代码结构,是不用推翻重来的。