冷热电微网双层优化这几年在综合能源方向几乎成了标配选题MATLABYALMIP也成了最常见的实现组合。但说实话光看论文摘要觉得很简单真到自己写代码时才发现一堆坑上层和下层怎么嵌套储能SOC怎么递推才不会出现同时充放电YALMIP报Infeasible之后从哪查起这篇内容我按自己做模拟项目X时踩过的路把双层优化的建模思路、求解路线、代码骨架和调参心得一次讲清楚。适合正在做微网储能优化、课程设计或者准备论文的同行也适合刚接触YALMIP、想快速搭起可运行模型的同学。先说清楚一个根本问题双层优化为什么非要用“两层”而不是一个大模型搞定这直接决定了后面所有代码怎么写。1. 先想清楚双层优化到底在优化什么1.1 时间尺度差异让问题天然分层冷热电微网里面有容量规划和运行调度两类决策它们的决策频率完全不同。容量规划决定“装多大的燃气轮机、装多少储能、配多少光伏”这个决策一做就是十几年属于年尺度甚至项目周期尺度的长期投资问题。运行调度决定“某一天某一时刻燃气轮机发多少电、储能充多少放多少、从电网买多少电”这个决策是小时尺度的短期运行问题。如果把两类变量塞进一个模型目标函数里既要算二十年投资回报又要算8760个小时的逐时购电成本求解器面对的量级和语义完全不一样很难在同一个优化框架下做权衡。而且容量和运行本身就存在主从关系定了容量之后系统的运行受容量上限约束反过来运行的经济性又决定了大容量值不值。这种“你先定、我后跑”的主从关系用双层模型表达最自然。我在做模拟项目X时最大的体会是先别急着写代码先把上层变量、下层变量、各自目标、各自约束列成一张表你就能发现上下层交换的信息到底是什么。上面传下来的是容量上限下面传回去的是年运行成本或收益。这不只是概念清楚的问题它直接决定你外层适应度函数怎么写。1.2 上层与下层的数学模型上层优化通常是最小化全生命周期总成本或者最大化项目净收益。一般写成min C_total C_inv C_om C_run其中C_inv是设备投资的等年值C_om是年维护费用C_run是内层优化返回的年运行费用购电、购气成本减去售电收益等。上层决策变量是容量比如燃气轮机额定功率、储能额定容量和功率、电制冷机容量、光伏装机容量。下层优化在给定容量下求解典型日或全年时段的运行方案min C_run sum(购电费用 购气费用 - 售电收益)下层决策变量是各时刻的设备出力、储能充放电功率、购售电功率、启停状态等。约束包括功率平衡、设备出力上下限、储能SOC递推、爬坡约束、备用约束等。注意上层的C_inv需要折算成年值公式是C_inv C_cap_total * r * (1r)^n / ((1r)^n - 1)这里的r是折现率n是项目年限。很多新手直接把总投资加进目标函数内层返回的却是一年费用这样量纲不一致优化出来的容量必然偏小。我一开始也在这个地方迷糊过后来统一把上层所有成本都折算到年才得到合理结果。1.3 储能是链接上下层的关键枢纽储能装置在这类模型中非常特殊它的容量属于上层变量它的充放电功率属于下层变量而且它的运行状态不能孤立确定必须满足跨时段递推关系。换句话说储能是上下层信息交换的核心载体。上层决定了电池装多大下层只能在“这块电池”的范围内动态调度。实话说如果模型里没有储能双层优化就会退化成单层配置问题加一堆运行约束复杂度降一个量级。正因为储能引入了SOC的状态变量和时序耦合下层才真正变成一个动态优化问题上层也才需要通过反复迭代来寻找容量最优值。这也是为什么绝大多数论文和项目都把储能作为优化重点的原因。2. 微网模型搭建把冷热电和储能写成数学表达式2.1 能量流与设备库先搭能量流框架。冷热电微网系统通常包含三条母线电母线、热母线和冷母线。设备按“源—转换—存储—负荷”四类划分。常见组成源侧光伏、风电、电网购电转换设备燃气轮机发电兼余热回收、燃气锅炉补热、吸收式制冷机用热制冷、电制冷机用电制冷存储设备蓄电池存电、蓄热罐存热、蓄冷罐存冷负荷电负荷、热负荷、冷负荷这三条母线之间通过设备耦合比如燃气轮机发电的同时产生余热余热进入热母线一部分直接供热一部分进吸收式制冷机转成冷。你可以想象成电力、热水、冷水三套管路在一个园区里交织优化要做的就是决定每一时刻能量往哪条管路走、走多少。2.2 关键设备数学模型设备建模不用搞得太复杂工程上常用简化线性模型保主要特性即可。燃气轮机或内燃机的发电模型P_gt eta_gt * LHV * F_gas其中eta_gt是发电效率LHV是天然气低位热值F_gas是燃料消耗速率。它的余热回收量可以写成Q_hr alpha_hr * P_gtalpha_hr是热电比系数不同机型的温差、排烟热量不同取值通常在0.5到1.5之间。我曾经见过有人把Q_hr写成(1-eta_gt)*P_gt这其实是错的因为排烟损失、散热损失都没算进去回收不到那么多热量。燃气锅炉模型Q_boiler eta_boiler * LHV * F_boiler吸收式制冷机Q_cold_abs COP_abs * Q_heat_in电制冷机Q_cold_ec COP_ec * P_ec光伏模型可以简化成P_pv(t) P_pv_cap * irradiance_norm(t)irradiance_norm是归一化辐照序列如果有温度系数可以再加一个温度修正项但对容量规划问题来说归一化辐照已经够用。把这些设备约束写成线性不等式组交给YALMIP处理后交给求解器整个过程会很顺畅。关键是别在建模阶段就引入非线性比如燃料成本函数做成分段线性、设备效率做成常数线性模型跑出来的结果已经足够支撑容量决策趋势分析。2.3 储能建模是双层耦合的核心蓄电池SOC递推公式是这类模型最容易出错的地方SOC(t1) SOC(t) eta_ch * P_ch(t) * dt / E_cap - P_dis(t) * dt / (eta_dis * E_cap)注意充电效率eta_ch是乘在充电功率上的放电效率eta_dis是放在分母上的。含义很简单充进去1度电电池实际存下的只有eta_ch度放出来1度电电池内部需要消耗1/eta_dis度的能量。效率位置一旦写反SOC递推就会出现能量不守恒运行成本也会虚高。重要约束0 P_ch(t) P_ch_max 0 P_dis(t) P_dis_max SOC_min SOC(t) SOC_max还要处理一个工程细节防止同一时刻既充电又放电。YALMIP里可以直接加互斥约束也可以引入二进制变量P_ch(t) M * u_ch(t) P_dis(t) M * (1 - u_ch(t))其中M是足够大的常数u_ch是0-1变量。有时候电池模型本身就允许轻微同时充放但为了竞赛和论文评审的严谨性建议还是加上互斥约束。蓄热罐、蓄冷罐的模型与电池类似只是把电功率换成热功率/冷功率效率和容量约束同理。储能SOC的初值和终值也要写成等式或不等式比如要求一个调度周期结束时SOC回到初始值这样才能让日调度方案可重复执行。3. 求解路线选择主流的三条打法3.1 三种路线对比双层优化的求解方法有不少我实际用下来主要分三条路线求解路线基本思路优点缺点适用范围智能算法内层数学规划外层用遗传算法/粒子群搜容量内层用YALMIP求解运行调度实现灵活、通用性强、容易调试计算量大结果有一定随机性绝大多数工程和论文场景KKT条件转单层把下层问题用KKT条件等价替换合并成单层约束数学严格解是确定的推导繁琐要求下层问题凸且强对偶成立下层为线性规划的简化场景MPEC/强对偶变换将对偶变量纳入模型形成含均衡约束的数学规划理论完整学术价值高建模复杂求解器支持有限小规模理论证明类题目我最推荐也最常用的是第一条路线——遗传算法YALMIP。理由很现实冷热电微网下层通常是混合整数线性规划KKT条件写起来极其痛苦MPEC对求解器和建模技巧要求又高而智能算法套YALMIP几乎不需要额外推导只要把上下层接口写好就能跑。计算量确实大但通过典型日加权、并行计算完全能接受。3.2 遗传算法作为外层框架的核心流程外层GA的框架大概是这样的初始化种群每个个体是一组容量配置比如[P_gt_cap, E_sto_cap, P_sto_ch_max, P_ec_cap, P_pv_cap]对每个个体调用内层运行优化得到典型日运行成本用运行成本加上投资年值计算适应度选择、交叉、变异生成下一代迭代直到收敛或达到最大代数这里有个很关键的设计点GA的适应度函数中运行成本必须是从内层优化器返回的不能自己凭经验估。因为内层调度结果直接反映容量配置是否合理容量太小会导致无法满足负荷或大量购电容量太大会造成投资浪费。GA搜索的过程本质上就是不断让内层调度去验证各种容量组合的运行效果。我在编码时用实数编码而不是二进制编码原因是容量变量连续且边界明确实数编码直接对应实际工程单位避免编解码误差也方便后期加约束。搜索范围设置成设备可能容量的1.2倍左右不要设太大否则收敛很慢。3.3 为什么内层必须用YALMIP而不是手写求解内层运行调度是多时段、多设备、多约束的优化问题手写求解器或者用单纯形法自己编码实现工作量巨大且极易出错。YALMIP的核心价值在于它让你用声明式的方式描述优化问题而不需要关心标准型转换、求解器接口这些底层细节。你只需要定义sdpvar决策变量用Constraints [Constraints, ...]累加约束再调用optimize剩下的交给求解器。还有一个好处YALMIP天然支持混合整数变量。当你需要处理机组启停、储能互斥、阶梯电价档位选择时只需要把一个变量声明为binvar求解器会自动按混合整数规划处理。如果手写求解器光是分支定界就得写一堆代码。4. MATLABYALMIP代码骨架与实现要点4.1 工程目录结构我习惯把项目拆成清晰的功能模块模拟项目X的目录如下/project_x/ data/load_curve.mat % 负荷、辐照、电价等时间序列 scripts/main.m % 主程序入口 scripts/upper_ga.m % 外层遗传算法 scripts/inner_schedule.m % 内层运行调度 utils/calc_investment.m % 投资年值计算 output/results.mat % 结果输出主程序只需要读数据、设置优化参数、调用上层GA、输出结果。内层调度函数接收容量参数和基础数据返回运行成本。这样写的好处是可以单独调试内层确认无误后再接外层避免两边同时出错时无从排查。4.2 内层运行调度代码骨架内层函数核心是YALMIP建模。以一个24小时典型日为例function run_cost inner_schedule(cap, data) T 24; % 决策变量 P_gt sdpvar(T, 1); % 燃气轮机发电功率 P_ch sdpvar(T, 1); % 储能充电功率 P_dis sdpvar(T, 1); % 储能放电功率 SOC sdpvar(T 1, 1); % 储能SOC P_grid sdpvar(T, 1); % 购电功率 P_ec sdpvar(T, 1); % 电制冷机功率 % 基础数据 P_load_e data.P_load_e; % 电负荷 Q_load_h data.Q_load_h; % 热负荷 Q_load_c data.Q_load_c; % 冷负荷 P_pv data.P_pv; % 光伏出力序列 price data.price; % 分时电价 Constraints []; for t 1:T % 设备出力上下限 Constraints [Constraints, 0 P_gt(t) cap.P_gt_cap]; Constraints [Constraints, 0 P_ch(t) cap.P_ch_max]; Constraints [Constraints, 0 P_dis(t) cap.P_dis_max]; Constraints [Constraints, 0 P_ec(t) cap.P_ec_cap]; % 储能SOC递推 Constraints [Constraints, ... SOC(t1) SOC(t) 0.95 * P_ch(t) * 1 / cap.E_sto ... - P_dis(t) * 1 / (0.95 * cap.E_sto)]; Constraints [Constraints, 0.1 SOC(t1) 0.9]; % 电功率平衡 Constraints [Constraints, ... P_gt(t) P_dis(t) P_pv(t) P_grid(t) ... P_load_e(t) P_ch(t) P_ec(t)]; % 冷功率平衡吸收式制冷用热、电制冷用电 % 这里简化吸收式制冷量 alpha_cool * Q_hr(t) % 详细模型视具体设备配置而定 end % 储能SOC首尾约束 Constraints [Constraints, SOC(1) SOC(T 1)]; % 目标函数购电费用 - 售电收益 购气费用 gas_cost sum(P_gt) * data.gas_price / data.eta_gt; elec_cost sum(P_grid .* price); Objective elec_cost gas_cost; % 求解 ops sdpsettings(solver, gurobi, verbose, 0); optimize(Constraints, Objective, ops); % 返回目标值 run_cost value(Objective); end注意我在SOC递推里用的充电效率0.95放电效率0.95单位是小时所以dt1。实际项目中如果时间步长是15分钟记得把dt带入公式否则SOC变化速度会明显失真。还有一个容易被忽略的点电功率平衡中的P_grid(t)如果允许负数代表向电网售电。这时目标函数里要区分购电价格和售电价格否则求解器会把“先多买再卖”当成套利手段模糊化处理会让结果失真。4.3 外层遗传算法主循环与并行加速外层GA的核心循环大致如下function [best_sol, best_fit] upper_ga(data, ga_options) pop_size ga_options.pop_size; max_gen ga_options.max_gen; bounds ga_options.bounds; % 初始化种群均匀随机分布在边界内 pop zeros(pop_size, n_vars); for i 1:pop_size for j 1:n_vars pop(i, j) bounds(j, 1) rand * (bounds(j, 2) - bounds(j, 1)); end end for gen 1:max_gen fitness zeros(pop_size, 1); parfor i 1:pop_size cap decode(pop(i, :)); C_inv calc_investment(cap); % 投资年值 C_run inner_schedule(cap, data); % 内层运行成本 fitness(i) C_inv C_run; % 总适应度 end % 选择、交叉、变异 % ... end end这里有几个实战经验parfor并行很有用但要确保每个worker都能访问YALMIP和求解器。如果并行池初始化后报找不到求解器就在并行循环内加上一行yalmip(clear)或重新sdpsettings。先串行跑一遍小规模用例确认内层全部没问题再开并行。并行调试一旦报错错误提示往往非常隐晦。GA种群一般设30到60迭代次数50到100。太小收敛不好太大计算时间成倍增加性价比不高。4.4 求解器配置要点YALMIP本身不是求解器它只是建模层。如果没有安装求解器optimize会报错。常用的组合是商业求解器如Gurobi或CPLEX它们对线性规划和混合整数线性规划的求解效率非常高。个人或者学术用途可以申请免费许可证。配置方式ops sdpsettings(solver, gurobi, verbose, 1);verbose设为1时会输出求解过程中的迭代信息方便查看问题规模和解的状态。跑成功后检查value(Objective)是否合理同时确认YALMIP返回的problem字段为00代表求解成功1代表不可行其他值代表数值问题。如果实在装不上商业求解器MATLAB自带的linprog和intlinprog也能救急。但要注意纯纯的linprog不能处理0-1变量遇到启停约束就得用intlinprog或另想办法。我的建议是尽早把商业求解器配好后面能省掉无数折腾时间。5. 实测复盘模拟项目X的完整参数与结果5.1 项目参数与数据设置模拟项目X是一个园区级冷热电微网基本参数如下参数数值说明光伏装机范围50 ~ 300 kW上层待优化燃气轮机容量范围0 ~ 500 kW上层待优化储能容量范围50 ~ 500 kWh上层待优化储能功率容量/4按2小时充放倍率设置折现率6%投资年值计算项目年限10年等年值折算负荷数据用三个典型日加权冬季典型日权重90天、夏季典型日权重90天、过渡季典型日权重185天。每个典型日24个点既有冷热电负荷曲线也有归一化光伏辐照曲线。这样做比直接跑8760小时快得多同时仍能反映季节差异。分时电价采用峰谷平三段式峰时9点到12点、17点到20点电价高平时6点到9点、12点到17点、20点到22点谷时22点到次日6点。这个设置很常见也是储能套利空间的主要来源。5.2 双层迭代的收敛情况我试了两组参数第一组GA种群30代、迭代40代第二组种群60代、迭代80代。第二组最后得到的容量配置比第一组总成本低3%左右但计算时间从40分钟涨到近3小时。这说明GA收敛性对参数敏感建议先小规模跑通再按需增加精度。内层单次调用耗时约0.3秒到1秒主要取决于YALMIP的约束个数和求解器的效率。上层一次完整迭代需要调用种群规模次内层这也是总计算时间的大头。收敛曲线通常表现为前10代快速下降后面进入缓慢优化阶段。如果前10代没有明显下降检查一下适应度函数有没有问题比如投资年值折算错误、内层返回值是负数等。5.3 储能容量对经济性的影响我单独固定其他容量扫描储能容量从50 kWh到500 kWh得到年总成本的变化趋势。结果并不是储能越大越好容量从50增加到200 kWh时年运行成本下降明显购电费用减少储能套利收益增加但超过200 kWh后投资年值的增加开始超过运行费用的节省总成本反而回升。这背后逻辑很好理解储能越大能转移的电量越多但边际套利收益是递减的——每天就那么多峰谷价差电池一天只能充放一次两次装得再大也消化不了而投资成本是线性的所以在某个拐点之后储能就成了“赔钱货”。拐点位置取决于峰谷价差、电池单位投资成本和循环寿命。这就解释了为什么储能容量必须靠优化而不是拍脑袋定。输出结果中我还对比了是否考虑SOC首尾约束的差异。不加首尾约束时模型会把有效容量放大实际运行时会出现第一天从满电开始、最后一天放空的“作弊”现象年成本偏低结论失真。加上SOC回位约束之后结果才真正可落地。6. 常见问题与排查心得6.1 问题速查表实际跑模型时遇到的问题我整理了一张速查表症状可能原因排查思路YALMIP报Infeasible约束冲突、变量无界、系数单位错误逐条注释约束缩小变量范围定位求解成功但结果明显不合理SOC递推效率写反、电价单位不一致检查公式推导和单位换算运行成本为负售电收益率远高于购电成本模型在套利区分购电/售电价格加数值上限GA长时间不收敛种群太小、搜索范围过大、适应度噪声大增大种群、缩小边界、固定随机种子parfor报求解器缺失并行池worker没有加载YALMIP和求解器在并行循环内重新配置求解器路径优化速度极慢约束太多、变量规模大、求解器配置不合理用典型日代替8760h简化设备模型6.2 不可行问题的排查流程遇到Infeasible千万别慌它不是“代码坏了”而是数学上约束自相矛盾。我的排查顺序是先求解一个“只有运行调度、不包含投资决策”的内层问题如果这个都不可行说明负荷平衡或SOC约束本身有问题。把SOC的上下界从0.1~0.9放宽到0~1排除SOC边界太紧导致的无解。把购电功率上限从1000 kW放宽到1e6 kW排除购电上限限制导致的不可行。检查冷热电平衡公式尤其是单位kW和kWh不能混热量单位要用一致。最后用YALMIP的diagnostic函数或者把Constraints拆成两部分分别求解很快就能定位冲突的那个约束组。还有一个容易被忽略的原因储能SOC初值和终值相等加上充放电互斥可能在高负荷日没有可行解。这时可以把首尾约束从等式改成不等式比如允许终值不低于初值的90%既保住可重复性又给求解器更多自由度。6.3 计算时间优化心得双层优化最烦人的就是慢。8760小时全时段内层优化会让YALMIP建模变得很重变量数量动辄十几万以上求解时间让人崩溃。我的经验是按季节选典型日每个典型日24个点乘上对应的天数权重。这样内层问题规模小两个量级结果却能反映全年的经济性走向。如果一定要跑全年可以考虑把储能SOC的跨天耦合拆开或者用滚动窗口近似。不过对容量规划来说典型日法完全够用没必要一开始就追求全年8760小时。还有一个提速技巧内层问题尽量保持线性避免分段二次或非线性目标。YALMIP遇到非线性会尝试调用非线性求解器速度慢且容易不稳定。把设备效率设为常数、燃料成本线性化、阶梯电价用0-1变量表达虽然约束多了但求解起来反而快。最后再说一个小技巧。写外层GA时给每个容量变量设一个合理的离散精度比如储能容量按10 kWh为步长取整光伏按5 kW取整。这样既压缩了搜索空间也让结果更贴近实际工程选型不至于优化出一个“203.7 kWh”的诡异容量。实话说这个项目反复调试下来我最深的感受是双层优化的难点不在算法本身而在于把物理问题严谨地转成数学问题再严谨地转成代码。先把单层运行模型调通再套外层迭代一步步来比任何捷径都稳。