1. 先把问题说清楚CCHP调度到底在优化什么最近我接手了一个园区综合能源系统的优化调度模拟项目系统里包含燃气轮机、余热回收、吸收式制冷、电制冷、燃气锅炉外加蓄电池和蓄热槽典型的一套冷热电联供CCHP配置。目标很明确在满足电、热、冷三类负荷的前提下把一天的运行成本压到最低。一开始我按习惯上了遗传算法后来换成粒子群最后真正跑顺的反而是很多人不太在意的POA算法——渐进优化算法英文全称 Progressive Optimality Algorithm。这篇博文把POA怎么用来解CCHP调度、MATLAB代码框架怎么搭、中间踩了哪些坑一次讲清楚。先说结论POA这套算法在电力系统里其实早就有应用尤其是梯级水电调度几乎属于经典方法。但把它搬到冷热电联供系统里的人不算多。我试过之后发现CCHP这种“多设备、强耦合、含储能、分时电价”的场景恰恰是POA的舒适区。它不需要像动态规划那样把状态变量离散到爆炸也不像遗传算法那样跑一次一个样、让人没法跟领导解释结果为什么每次都不同。只要模型写对、初值给得合理POA的收敛结果是稳定可复现的这一点在工程交付里非常重要。1.1 冷热电联供系统里到底有哪些决策要做CCHP系统的核心逻辑是“能量梯级利用”燃气轮机烧天然气发电发电后的高温烟气不扔掉进余热锅炉回收变成蒸汽或热水这部分热能再去驱动吸收式制冷机吹出冷风或者直接供热。这样一来一套系统同时产出电、热、冷三种能量综合能源利用率可以做到70%到85%比单独发电高出一大截。但问题也出在这个“同时产出”上。电负荷、热负荷、冷负荷的变化规律不一样设备之间又有强耦合比如燃气轮机发多少电就会连带产生多少可利用的余热吸收式制冷机用多少热又会反过来挤占供热能力。再加上分时电价的影响白天电价高你可能希望多开燃气轮机、少买电晚上电价低你又想让蓄电池充电、把燃气轮机压到最低。这里面的决策变量包括每个时段的燃气轮机出力、燃气锅炉出力、吸收式制冷量、电制冷量、蓄电池充放电功率、蓄热槽充放热功率林林总总算下来一个24时段的调度问题决策变量轻松超过100个。而且这个系统不是拍脑袋定的它必须满足三条硬约束电功率平衡、热功率平衡、冷功率平衡。储能设备还额外带来时序耦合——这一小时充了多少电直接影响下一小时能放多少电蓄电池的荷电状态SOC是一条贯穿全天的链条不能只看单时段。这就是为什么这类问题不能简单地“每个时段独立优化”必须放在一个完整的时间轴上整体求解。1.2 为什么我没选动态规划、也没选遗传算法最开始我试图用动态规划DP来做。DP在理论上是最优的但它的痛点非常现实状态变量每增加一个维度状态数量就指数级膨胀。CCHP系统里至少有蓄电池SOC、蓄热槽SOC两个储能状态每个状态离散成100个点组合起来就是100乘100等于一万个状态再乘以24个时段每个状态还要做一次设备组合优化算到后面内存和时间都吃不消。想要精度更高离散点数加多直接算不动。遗传算法和粒子群这类智能算法我也跑了。它们的好处是能处理非凸、非线性目标不用求梯度但问题是结果有随机性每次运行成本可能差1%到2%调参还极其玄学种群大小、交叉概率、变异概率随便动一下收敛曲线就换一副面孔。更麻烦的是工程上你很难判断算法到底是收敛到了好的解还是碰巧撞了一个好结果。我在模拟项目X里跑了20次遗传算法最好和最差之间能差出将近5%的成本这对一个以“降本”为卖点的系统来说决策风险太大了。POA的思路则完全不同。它把一整个N时段优化问题拆成N减1个“两阶段子问题”每次只优化相邻两个时段的决策其他时段全部固定。听着简单但它有理论依据基于Bellman最优性原理一条最优决策轨迹的任何子轨迹相对其端点来说也是最优的。所以固定前后端点、反复扫描相邻时段就能逐步逼近全局最优。这个算法不需要离散状态、不需要求梯度、不依赖随机数种子收敛稳定内存开销也小。1.3 POA在CCHP场景下的整体求解思路具体到CCHP调度我把每个时段的设备出力组合定义为决策变量把蓄电池SOC和蓄热槽SOC作为跨时段状态量。POA的流程分三层外层是迭代循环从初始轨迹出发反复执行“从第1时段扫到第N时段”的扫描过程中层是两阶段子问题求解。对任意相邻时段t和t1固定其他时段的决策只调整这两个时段的设备出力和储能功率用非线性规划求解器做局部优化内层是状态更新子问题解完后更新该时段的SOC和下一时段的起始SOC再进入下一对相邻时段。一个关键细节是每次子问题里两端的状态边界要取“上一轮扫描得到的值”不是当前轮已经更新过的值。如果顺序更新时不加区分误差会从第一对时段一路传导到最后一对导致目标函数不降反升。这个问题我一开始没注意结果POA迭代了好几十轮还在振荡最后把状态传递逻辑改成“先全部读旧值、解完再统一写新值”才稳定下来。2. MATLAB代码框架与参数准备从数据到求解器代码不是一上来就写求解器的。我个人的习惯是先把系统参数、负荷曲线、能源价格整理成结构体再写目标函数和约束最后才碰POA循环。顺序反了的话后面改一个参数就要翻半天代码非常痛苦。2.1 设备参数、负荷曲线与分时电价准备我用的是24个时段每时段1小时。如果要跑96点每15分钟一个点模型结构不用大改只要把数组维度从24改成96再把爬坡约束的时间尺度换算一下就行。下面这张表是我在模拟项目X里用到的参数单位统一为kW和kW·h供参考设备关键参数数值燃气轮机额定电功率800 kW燃气轮机发电效率0.33燃气轮机热电比1.4余热锅炉余热回收效率0.8吸收式制冷机COP1.2电制冷机COP4.5燃气锅炉热效率0.9蓄电池容量600 kW·h蓄电池最大充放电功率200 kW蓄电池充放电效率0.92蓄热槽容量800 kW·h蓄热槽最大充放热功率250 kW负荷数据我用的是夏季典型日的曲线电负荷白天高、晚上低冷负荷在下午两三点冲到峰值热负荷主要是生活热水相对平稳。分时电价采用常见的峰平谷三段峰时段8点到11点、18点到21点电价1.2元/kW·h平时段11点到18点、21点到23点电价0.7元/kW·h谷时段23点到次日8点电价0.4元/kW·h。天然气价格固定为2.6元/立方米折算成单位热值成本后参与目标函数计算。这里有个容易踩的坑蓄电池和蓄热槽的SOC初值由我手动指定但一天结束后的终值不能随意。工程上通常要求一天结束回到初始SOC比如0.5保证第二天调度可以连续运行。这个“终端SOC约束”必须在模型的约束函数里显式写出来否则优化结果会让储能设备在最后几个时段疯狂放电把SOC榨干第二天直接没法用。2.2 代码模块划分与关键数据结构我的MATLAB工程目录不复杂坑都在每个文件内部。整体分四个模块data_par.m定义设备参数、负荷、电价、气价统一放到一个结构体par里obj_fun.m目标函数输入是决策变量矩阵和参数结构体输出是全天运行总成本cons_fun.m约束函数包括等式平衡约束和不等式容量约束输出给非线性规划求解器poa_main.mPOA主循环负责初始化轨迹、调用子问题求解、判断收敛、输出结果。决策变量我用一个二维矩阵X表示行数等于时段数24列数等于每个时段的决策变量个数。具体列的顺序我固定为燃气轮机出力、燃气锅炉出力、吸收式制冷功率、电制冷功率、蓄电池充电功率、蓄电池放电功率、蓄热槽充热功率、蓄热槽放热功率。列顺序一旦定了后面所有函数都按这个约定索引别换来换去。实际编码中我强烈建议把所有功率统一归一到基准值比如以100kW为基准或者以燃气轮机额定功率为基准。我自己最开始没归一化燃气轮机出力是几百kW、成本是几万元、效率是0.3几个数量级完全不一样的东西一起进优化器结果fmincon一直在原地打转目标函数下降极慢。归一化之后所有变量基本都在0到1的范围内求解器的数值稳定性好很多。2.3 初始轨迹的生成方式POA属于局部寻优算法初始轨迹会给定的好收敛就快结果也更优给得差虽然也能收敛但有可能落在一个不太好的局部解。所以我不建议随机初始化更推荐用工程启发式方法生成两个候选初始解跑一遍取成本低的那条。第一种是“电跟随”策略燃气轮机按电负荷的某个比例出力余热全部回收优先满足热负荷不够的部分用燃气锅炉补冷负荷先让吸收式制冷机吃余热不足再用电制冷。第二种是“热跟随”策略燃气轮机按热负荷反推出力发出的电不够就用购电补冷负荷仍然先吸收式后电制冷。两种策略都有道理但成本差异不小。在模拟项目X的夏季典型日里电跟随初始解对应的全天运行成本大约是7900元热跟随大约是8100元最后POA用更好的那条初始轨迹起跑迭代次数明显更少。代码里就是把这两种策略各写一个30行的循环生成两条初始轨迹选目标函数值更小的作为启动轨迹。3. POA求解CCHP的核心代码逐段拆解这一部分是文章的重点。我不会贴整份工程代码那样篇幅太长也不利于理解但会把目标函数、约束、子问题求解、主循环这几块最核心的代码框架拿出来逐行解释。你能看明白这个框架在自己的项目里改参数就能直接复现。3.1 目标函数与约束的MATLAB实现目标函数的物理含义是全天运行总成本我把它拆成四块燃气费用、购电费用减去售电收入、设备维护费用、储能折旧折算。总表达式写成function J obj_fun(X, par) N par.N; % 时段数 Pgt X(:,1); % 燃气轮机出力 Pgb X(:,2); % 燃气锅炉热出力 Pac X(:,3); % 吸收式制冷功率 Pec X(:,4); % 电制冷功率 Pch X(:,5); % 蓄电池充电功率 Pdis X(:,6); % 蓄电池放电功率 Qch X(:,7); % 蓄热槽充热功率 Qdis X(:,8); % 蓄热槽放热功率 % 燃气轮机耗气量与电出力和热电比有关 gas_gt Pgt ./ par.eta_gt; % 单位kW 天然气热量 % 燃气锅炉耗气量 gas_gb Pgb ./ par.eta_gb; % 燃料成本 fuel_cost sum(gas_gt gas_gb) * par.gas_price; % 购电成本假设不能向电网售电简化处理 Pbuy max(0, par.Pload - Pgt - Pdis Pec Pch); buy_cost sum(Pbuy .* par.elec_price); % 维护成本按出力线性折算 om_cost sum(Pgt) * par.om_gt sum(Pgb) * par.om_gb ... sum(Pac) * par.om_ac sum(Pec) * par.om_ec; % 储能折旧按充放电量折算延长寿命 bat_cost sum(Pch Pdis) * par.om_bat; tank_cost sum(Qch Qdis) * par.om_tank; J fuel_cost buy_cost om_cost bat_cost tank_cost; end这里有几个要点。第一购电功率Pbuy我写成max(0, ...)的形式相当于默认系统不向电网倒送电。如果你的项目里允许余电上网就要额外加一个售电变量和售电收益项去点max改成等式约束。第二维护成本和储能折旧不在物理平衡约束里体现直接体现在目标函数中让优化器自动权衡“多用设备省钱”和“过度损耗设备多掏维护费”这两件事。约束函数重点处理三类约束。第一类是功率平衡等式约束典型写法如下function [c, ceq] cons_fun(X, par) N par.N; Pgt X(:,1); Pgb X(:,2); Pac X(:,3); Pec X(:,4); Pch X(:,5); Pdis X(:,6); Qch X(:,7); Qdis X(:,8); % 电平衡发电购电放电 负荷电制冷充电 Pbuy max(0, par.Pload - Pgt - Pdis Pec Pch); ceq1 Pgt Pdis Pbuy - par.Pload - Pec - Pch; % 热平衡余热锅炉蓄放热 热负荷供热型吸收式耗热 q_gt par.hpr * Pgt; % 燃气轮机余热 q_hrsg q_gt * par.eta_hrsg; ceq2 q_hrsg Pgb Qdis - par.Qload - Qch - Pac / par.cop_ac; % 冷平衡吸收式电制冷 冷负荷 ceq3 Pac Pec - par.Cload; % SOC递推约束蓄电池 soc_bat par.soc_bat_init; for t 1:N soc_bat soc_bat (Pch(t)*par.eta_ch - Pdis(t)/par.eta_dis) / par.bat_cap; end ceq4 soc_bat - par.soc_bat_init; % 终端回到初值 % 不等式约束设备容量边界、SOC上下限、爬坡约束 c [...]; ceq [ceq1; ceq2; ceq3; ceq4]; end热平衡等式里吸收式制冷消耗的热量是Pac / par.cop_ac这个单位换算是新手最容易漏的地方。很多人把吸收式制冷功率Pac直接当成消耗的热功率实际不对吸收式制冷机的COP通常是1.2左右意思是消耗1份热可以产出1.2份冷所以热消耗是制冷量除以COP。3.2 POA主循环和两阶段子问题怎么搭POA主循环的核心思想是一次只优化两个相邻时段。我把子问题封装在一个函数里输入是当前轨迹X、时段t、参数par输出是更新后的第t和第t1时段决策。function X poa_main(X0, par) X X0; N par.N; maxIter par.maxIter; tol par.tol; J_prev inf; J_hist zeros(maxIter, 1); for iter 1:maxIter for t 1:N-1 X solve_subproblem(X, t, par); end J obj_fun(X, par); J_hist(iter) J; if abs(J - J_prev) tol break; end J_prev J; end end子问题的实现是重头戏。相邻时段t和t1的决策变量总共有16个每个时段8个变量我把它们从整个轨迹里提出来调用fmincon求解解完再放回去。关键点在于第t时段之前的SOCl状态和第t2时段之后的状态都要固定为上一轮的值不能动。这样才符合POA“只动两个端点之间”的原则。function X solve_subproblem(X, t, par) n par.nVar; % 每个时段的决策变量个数 % 提取当前两个时段决策作为初值 x0 [X(t,:), X(t1,:)]; lb [par.lb_t, par.lb_t]; % 两个时段的下界 ub [par.ub_t, par.ub_t]; % 等式约束函数第三个参数是固定边界信息 ceq (x) sub_cons(x, X, t, par); options optimoptions(fmincon, Algorithm, sqp, ... MaxIterations, 200, MaxFunctionEvaluations, 3000, ... OptimalityTolerance, 1e-6, Display, off); [xopt, ~] fmincon((x) sub_obj(x, par, t), x0, [], [], [], [], ... lb, ub, ceq, options); % 解回填到轨迹 X(t,:) xopt(1:n); X(t1,:) xopt(n1:2*n); end关于sub_obj和sub_cons的实现我再补充一点相邻两个时段的子问题目标函数不用重新写一个完整版直接把两个时段各自的成本项加总就行。约束条件里除了两个时段的功率平衡还要把储能SOC递推关系写进去并且固定边界SOC第t时段的初始SOC等于当前轨迹中第t-1时段算出来的值第t1时段的末端SOC等于当前轨迹中第t2时段起始的SOC。这里最隐蔽的一个坑是固定边界的“当前值”到底取哪一轮。我在调试过程中发现如果子问题里用的边界SOC是本轮已经更新过的整个迭代过程就会出现“追尾巴”效应目标函数降几轮之后开始震荡。解决方案是在主循环里做一个“影子轨迹”X_old X子问题约束里统一读X_old的边界值解完再更新X。虽然多占了一点内存但稳定性提升非常明显。3.3 收敛判据与加速收敛的实用手段POA的外层迭代不需要跑很多轮我在24时段模型上通常设定最大迭代次数为30实际上15到20轮就收敛了。判据用目标函数值相对变化量小于1e-4就停。有个细节是第一次迭代的J_prev不能设成inf吗可以但建议第一次迭代完成后直接把J_prev J赋上同时记录一下初值轨迹的成本方便对比POA到底降了多少成本。加速收敛还有一个非常有效的手段双向扫描。标准POA是从t1扫到tN-1但单向扫描有时会在目标函数比较平坦的区域“卡住”几个轮次。我在主循环里改成一次正向扫描加一次反向扫描也就是t从1到N-1再从N-1到1各来一遍一正一反作为一个外层迭代。实测在同样的容差下总迭代次数能少个4到6轮。如果觉得fmincon求解子问题比较慢可以换个思路把两个时段的变量用“决策变量离散化枚举”的方式直接搜索。因为相邻时段只有16个决策变量你可以在初值附近构造一个网格枚举几十个候选方案选目标函数最小的那个。这种方式的好处是稳定、不用调求解器参数缺点是精度受网格密度限制。我建议先用精细网格跑一遍拿到结果再用这个结果作为fmincon的初值做一次微调等于“粗筛精修”两段式求解速度和精度都能兼顾。4. 算例结果、常见问题与调试避坑模型搭好、代码跑通只是开始真正花时间的是看结果、调参数、处理各种“说不上来哪里出错但就是不对”的问题。这一节我把典型结果和一个经常要翻的排查清单都整理出来。4.1 典型日结果长什么样在模拟项目X的夏季典型日里负荷峰值是电负荷820kW、冷负荷600kW、热负荷280kW。分别用规则策略、遗传算法、POA求解全天总成本如下模拟数据仅供对比趋势求解方式全天运行成本/元备注电跟随规则策略7963无优化纯启发式遗传算法7428种群100迭代200标准差约110POA7369确定性结果多次运行一致POA跑出来的设备运行规律也符合直觉白天电价高、冷负荷大燃气轮机基本顶着额定功率运行发电满足电负荷后余热优先喂给吸收式制冷机下午冷负荷峰值时段电制冷机会补上一部分缺口因为此时燃气轮机的余热已经用满了。蓄电池在凌晨谷电时段充电上午峰电价时段放电下午再充再放做一次峰谷套利。蓄热槽则更像一个“热缓冲池”燃气轮机余热多的时候先把热存在槽里晚上热负荷高的时候再放出来。收敛曲线方面POA前5轮成本下降最快大概能降掉总降幅的80%之后是缓慢修正从第二轮到第18轮累计只降了不到1%。这就提示我们如果只做工程估算跑个8轮就停也能接受但如果做正式方案比选建议还是跑到完整收敛反正也就几秒钟的事。4.2 常见问题排查速查表我把自己调试过程中遇到的高频问题整理成了一张表按症状、原因、解决办法三列排列。症状可能原因解法目标函数前几轮下降后面开始振荡子问题边界读的是本轮更新值改用影子轨迹统一读上轮值蓄电池SOC被优化到边界终值不归位缺终端SOC约束约束函数显式加ceq4终端回到初值电平衡/热平衡等式相对误差大罚函数系数太小或fmincon容差太宽提高罚系数把OptimalityTolerance收严到1e-6检查单位是否统一迭代很多轮但成本几乎不动初始轨迹太差落进局部解换电跟随/热跟随两条初值取较优或加双向扫描96点模型跑起来非常慢单次fmincon调用次数太多先用24点调通逻辑子问题改用粗网格枚举fmincon精修三套平衡约束总有一两个时段不满足冷、热、电量纲混了COP换算漏除逐时段检查三个平衡残差定位是哪类设备项写错其中“量纲混了”是我见过最隐蔽的问题。冷负荷单位是kW但吸收式制冷消耗的热量是Pac/cop_ac电制冷消耗的电量是Pec/cop_ec如果不除以COP冷平衡表面看起来凑得上实际在热平衡和电平衡里就会多出一块说不清的能量。排查的时候我习惯把约束残差单独打印出来对每个时段看电、热、冷三条残差曲线哪条不平整就去查哪个平衡式。4.3 几条独家调试经验最后分享几条从实战里总结出来的经验都不是能直接搜到的内容。第一fmincon的算法设置建议用sqp而不是interior-point。内点法在变量边界接近约束边界时偶尔会卡在可行域内部导致SOC进不了高位区间SQP处理这类边界约束更直接收敛更干净。当然这不是绝对的但我在两个不同配置的CCHP模型上都观察到类似现象。第二罚函数的系数不是越大越好但要保证等式约束的残差在最终解里落到一个很小的量级。我自己习惯的做法是先在目标函数里不带罚项跑一次看约束残差水平再把罚系数从1e3起步翻倍往上试找到“约束残差降低速度明显放缓”的那个拐点。这个拐点附近的系数通常就是比较合适的值。第三把24点扩展到96点时不要直接改维度就开跑。我吃过这个亏第一次扩到96点运行时间直接翻了不止4倍因为fmincon在每个子问题里的函数评估次数没有跟着时段长度缩放。后来我把子问题里的MaxFunctionEvaluations从3000降到1500外层POA容差从1e-4放宽到5e-4求解时间降了一半多结果质量差别小于0.3%完全可以接受。说白了工程上不是每个环节都要追求数学意义上的最优够用就行。第四结果出来后一定要画一张“全设备出力堆叠图”而不是只看成本数字。堆叠图把燃气轮机、购电、蓄电池放电、蓄热槽放热按顺序叠起来和负荷曲线放在一起看能在五秒钟内发现“电负荷曲线有个尖峰但蓄电池在那个时段恰好放不出电”这类潜在问题。我曾经就对着一张堆叠图发现了蓄热槽容量差200kW·h导致凌晨热负荷缺口改完参数成本又降了1.2%。我个人做完这个项目最大的体会是POA算法在冷热电联供调度里属于“上限不算惊艳但下限非常高”的方案。它不像遗传算法那样能给你一种“正在全局搜索”的幻觉也不像动态规划那样让你盯着状态网格数量发愁。它的优势就是稳定、可复现、可解释而且对工程落地的人来说这三条比“理论最优”重要得多。最后再分享一个小技巧如果你后续想把单目标扩展成“成本加碳排放”的双目标不需要推翻POA框架只要把目标函数里的成本项换成加权和就行外层迭代逻辑完全不用动。这个算法框架的扩展性是我当初选择它的时候没预料到的惊喜。