前几天有个师妹拿着一版IEEE 33节点配电网的日前调度结果来找我眉头皱得很紧模型跑出来明明电压无功都合格、购电成本也压到了最低可拿同一组数据回放之后发现好几个时段节点电压已经掉到0.93实际电费比模型输出高了差不多15%。我扫了一眼她的代码典型的单阶段确定性优化光伏和风电出力按预测曲线写死储能SOC也按预测值一路推过去所有约束都在“预测完全正确”这个大前提下自洽运行。问题恰恰出在这里——当分布式电源大量接入配电网后预测误差不再是锦上添花时可以忽略的细节而是一个能直接改写运行方式的变量。这篇就围绕“含分布式电源的配电网日前两阶段优化调度模型”展开把为什么单阶段会失效、两阶段模型怎么建模、Matlab代码怎么落地、调试时最容易踩的坑一次说清楚。适合正准备做配电网优化调度的研究生、刚入门的调度工程师以及想用YALMIP快速搭建算例的读者。1. 分布式电源接入后单阶段调度为什么会失效1.1 三个被打破的传统假设传统配电网调度能够用一套确定性优化模型跑到底背后其实有三个隐含前提。第一个前提是电源完全可控。传统配电网的电源主要来自上级电网变电站出口的功率和电压可以主动调节调度员说买多少就买多少说发多少就发多少。分布式光伏和风电接入之后情况变了。光伏出力跟着光照强度走云层一飘过出力可能几分钟内掉一半风电出力跟着风速走日夜波动毫无规律。这些电源没有燃料库存可以调节本质上属于“看天吃饭”可控性很弱。第二个前提是潮流单向流动。传统配电网的功率从变电站母线向负荷末端单向输送电压沿着馈线逐渐降低调压手段集中在变电站侧。分布式电源接入后局部的功率流向可能反转某些时段节点电压甚至可能高于变电站出口电压。这时再用“电压沿馈线单调下降”的经验去判断很容易得出错误结论。第三个前提是负荷和发电的预测误差足够小。传统配电网的负荷曲线有比较强的规律性预测误差通常只有几个百分点上级电网的自动发电控制可以把这点误差消化掉。但分布式电源叠加在负荷曲线上之后净负荷曲线变成“负荷减新能源出力”波动幅度成倍放大预测误差也跟着成倍放大。这三个前提被打破之后单阶段确定性模型的软肋就暴露了它把所有预测值当成真值算出一个在“预测完全正确的平行世界”里最优的方案。一旦真实运行数据偏离预测这个方案就处处掣肘。1.2 一组数字看懂单阶段模型的“账面最优”用一组数据模拟会更直观。假设10号母线接了1MW光伏24号母线接了0.5MW风机32号母线接了0.6MW/1.2MWh的储能。某天上午10点的日前预测是光伏出力0.9MW负荷1.2MW。单阶段模型为了让购电成本最低把联络线交换功率压到0.5MW储能同时在充电。到了实际运行那一刻天上云层一厚光伏出力掉到0.45MW负荷又涨到1.35MW。这时候节点净负荷是0.9MW日前安排的0.5MW购电根本不够用。差额0.4MW只能到实时市场高价购买实时电价常常是日前电价的1.5到2倍要是联络线容量有限还可能直接触发限负荷。更麻烦的是光伏出力骤降的同时局部电压会在很短时间内往下掉单阶段模型没有为这种场景预留调节裕度实际运行电压经常越限。这就是“账面最优、实际亏损”的直观来源。解决问题的关键是不要把一天的调度计划一次性锁死而是拆成“日前计划日内修正”两个阶段来做。日前阶段基于预测做出不可轻易更改的粗计划日内阶段利用更新的预测信息对粗计划做局部修正既保证经济性又保留对不确定性的应对能力。2. 两阶段模型的数学骨架从日前计划到日内修正2.1 第一阶段锁定哪些“不反悔”的决策第一阶段在日前执行决策时间尺度是小时级通常覆盖24个时段。这一阶段要定下来的是那些“事后难反悔”的变量主要包括联络线交换功率。这是和上级电网之间的购电协议一旦签了日前合同日内大幅调整就要付出较高的实时市场差价。燃气轮机的启停状态和出力计划。机组启停是0-1整数变量日内频繁启停不仅成本高对设备寿命也不利。储能充放电的基准曲线。虽然日内还可以微调储能但日内调整幅度受限于SOC和功率上限所以日前就要给一个大体合理的充放电框架。无功补偿装置投切计划如果模型需要细化无功调压也会放在这一阶段。第一阶段的目标函数一般写成min Σ_t [ c_buy(t)·P_buy(t) c_gas·P_mt(t) Σ c_om·P_dg(t) c_ess·(P_ch(t) P_dis(t)) c_loss·P_loss(t) ]其中c_buy是分时购电价P_buy是联络线购电功率c_gas是燃气轮机燃料成本c_om是分布式电源运行维护成本P_loss是网损项。需要注意第一阶段目标不能只盯着购电成本网损和DG运行成本都要放进去否则模型会倾向于把功率送到很远的末端节点产生不切实际的调度方案。2.2 第二阶段在最新预测驱动下“小步快跑”进入实际运行日每过15分钟或1小时就能拿到未来几个小时的超短期预测。超短期预测的误差远小于日前预测这时候把第一阶段给出的计划当成“基准”在它附近做修正就是第二阶段做的事情。第二阶段的决策变量基本都是连续量燃气轮机出力的微调量ΔP_mt储能充放电功率的修正量ΔP_ch和ΔP_dis可中断负荷的调节量光伏和风电的弃用功率ΔP_curt。第二阶段的目标函数以偏差惩罚为主min Σ_m [ α⁺·ΔP_buy⁺(m) α⁻·ΔP_buy⁻(m) β_curt·ΔP_curt(m) β_load·ΔP_loadcut(m) ]其中ΔP_buy⁺和ΔP_buy⁻是相对日前购电计划的向上/向下调整量α⁺和α⁻是相应的惩罚单价β_curt是弃风弃光惩罚β_load是切负荷惩罚。这个目标函数的核心含义是允许你偏离日前计划但偏离要付出代价允许你弃风弃光或切负荷但代价更高。这样优化器会优先动用储能和燃气轮机的调节能力而不是被动承受预测误差。第二阶段不是从头重新优化一遍而是在第一阶段计划的基础上做修正这一点非常关键。两阶段的本质区别不是变量多寡而是决策所依据的信息不同日前用的是粗预测日内用的是精预测。2.3 两阶段如何耦合基准点、偏差惩罚与不确定性两个阶段不是孤立跑的。第一阶段的结果会作为一组基准值传递到第二阶段第二阶段把实际运行变量限制在基准值附近偏离就产生惩罚成本。耦合关系可以写成P_buy_real(t) P_buy_ref(t) ΔP_buy⁺(t) − ΔP_buy⁻(t) P_mt_real(t) P_mt_ref(t) ΔP_mt(t) soc_real(t) soc_ref(t) Δsoc(t)这里的下角标ref表示第一阶段的日前计划值。第二阶段的决策变量就是各个Δ。这种“基准点偏差修正”的结构保证了日内运行不会完全推翻日前计划因为日前计划已经涉及购电协议、机组启停等硬约束同时它又给日内运行留了足够的调整空间。在不确定性建模上常见的有三条路线。随机规划给预测误差生成多个场景目标函数变成期望成本鲁棒优化考虑最坏场景目标是让最坏情况下的成本最小滚动时域/模型预测控制MPC每个时段滚动刷新预测并重新求解一个有限时域优化问题。Matlab工程实现中MPC思路最常用因为它不需要维护大量场景也不需要求解复杂的对偶问题只需要把两阶段模型封装成一个函数在各个时段滚动调用就行。前提是每个时段的求解要在几分钟内完成而线性化DistFlow加上LP/MILP的规模完全能够满足。2.4 潮流约束的选择线性DistFlow与SOCP的取舍配电网潮流约束是优化模型里最容易让人头疼的部分。完整交流潮流方程是非凸的放进优化模型会让问题变成NP-hard33节点小系统还勉强能算系统一大基本没法在可接受时间内求到最优解。工程上最常用的替代方案是DistFlow方程。DistFlow是为辐射状配电网专门设计的潮流简化形式原始版本还是包含P²Q²项的非线性方程。如果忽略网损项就得到线性化DistFlowP_ij(t) Σ P_jk(t) P_j,load(t) − P_j,DG(t)右边第一项是所有以j为首端的子支路有功功率之和。无功功率的方程结构完全一样。电压方程是V_j²(t) V_i²(t) − 2( r_ij·P_ij(t) x_ij·Q_ij(t) )这样处理之后整个模型变成LP或MILP求解速度快、数值稳定性好。代价是网损被忽略了。如果系统不算重载网损通常只占总负荷的百分之几线性化带来的误差完全在工程可接受范围内。如果项目要求精确计及网损可以升级到SOCP二阶锥松弛。具体做法是把DistFlow原始方程中的vV²和l(P²Q²)/V²用锥约束松弛为不等式拟合精度显著提升但模型复杂度和求解时间也会上升。我个人的习惯是先用线性DistFlow把模型跑通、把逻辑理顺再根据项目是否需要网损精度决定要不要升到SOCP。绝大多数研究场景和工程预研场景线性化已经够用。3. Matlab实现YALMIP建模、求解器与可复现代码3.1 环境准备与版本坑Matlab里做优化建模我不建议手写约束矩阵然后直接调linprog/intlinprog。约束一多、索引一复杂手写矩阵容易把人绕晕改一个参数就得重新梳理半天索引。推荐用YALMIP做建模层底层求解器用Gurobi。学术用户直接申请一个Gurobi license完全免费。环境配置上有几个版本坑值得提前说。YALMIP必须和当前Matlab版本兼容太老的YALMIP在新版Matlab上会出现莫名其妙的语法错误Gurobi版本太新而YALMIP太老优化时会提示找不到求解器如果装了Gurobi但Matlab提示Java相关错误多半是JVM路径或环境变量问题跟模型本身无关。验证环境是否正常的办法很简单在Matlab命令行跑一下yalmiptest能列出可用的求解器列表就说明环境OK。如果实在装不上Gurobi备选方案是Cplex或Mosek。小算例用Matlab自带的intlinprog也能跑只是速度慢一些对学习验证来说完全够用。3.2 IEEE 33节点算例与数据改造代码落地最常用的算例是IEEE 33节点系统。基准电压12.66kV基准容量通常取10MVA一共33个节点、32条支路1号节点是变电站根节点总负荷约3.715MWj2.3Mvar。数据网上很容易找到也可以用Matpower直接加载case33bw。DG接入位置按文献常见配置来设置10号节点接1MW光伏17号节点接0.8MW风机24号节点接0.6MW/1.2MWh储能8号节点接1MW微型燃气轮机。这样改造后在馈线中后段形成了分布式电源多点接入的格局最能体现两阶段模型的调节效果。负荷和新能源时序数据需要自己生成。常用做法是用一条典型日负荷率曲线乘以各节点额定峰值负荷得到24小时负荷光伏出力用光照强度曲线换算风电出力用风速曲线换算。为了模拟预测误差可以在日前预测值基础上叠加正态分布随机数标准差取5%到20%不等视你要研究的场景而定。3.3 核心代码从变量定义到求解下面这段代码是我常用的框架根据自己的算例改参数就能用。变量定义部分T 24; % 日前调度时段数 nb 33; % 节点数 nl 32; % 支路数 Pij sdpvar(nl, T, full); % 支路有功功率 Qij sdpvar(nl, T, full); % 支路无功功率 Vi2 sdpvar(nb, T, full); % 节点电压平方 Pbuy sdpvar(1, T); % 联络线购电功率 Pmt sdpvar(1, T); % 燃气轮机出力 soc sdpvar(1, T1, full); % 储能SOC Pch sdpvar(1, T, full); % 储能充电功率 Pdis sdpvar(1, T, full); % 储能放电功率约束构建部分。DistFlow约束用循环逐支路、逐时段构建Constraints []; % 根节点电压固定为1.0标幺 Constraints [Constraints, Vi2(1,:) 1.0]; % 逐时段构建DistFlow约束 for t 1:T for k 1:nl i BusFrom(k); % 支路首端节点 j BusTo(k); % 支路末端节点 child find(BusFrom j); % 以j为首端的子支路集合 Constraints [Constraints, ... Pij(k,t) sum(Pij(child,t)) Pd(j,t) - Ppv(j,t) - Pwt(j,t)]; Constraints [Constraints, ... Vi2(j,t) Vi2(i,t) - 2*(R(k)*Pij(k,t) X(k)*Qij(k,t))]; end % 电压上下限约束取平方后比较 Constraints [Constraints, 0.95^2 Vi2(:,t), Vi2(:,t) 1.05^2]; endPd是节点负荷矩阵Ppv和Pwt是节点光伏和风电注入矩阵这些数据要在进入优化前提前算好。燃气轮机和储能的功率项可以加到对应节点的功率平衡方程右侧我这里为了展示核心结构先省略。储能约束和SOC周期耦合soc_min 0.1*1.2; % 单位MWh soc_max 1.2; Pch_max 0.6; % 单位MW Pdis_max 0.6; eta_c 0.95; % 充电效率 eta_d 0.95; % 放电效率 for t 1:T Constraints [Constraints, ... soc(t1) soc(t) eta_c*Pch(t) - (1/eta_d)*Pdis(t)]; Constraints [Constraints, 0 Pch(t) Pch_max]; Constraints [Constraints, 0 Pdis(t) Pdis_max]; Constraints [Constraints, soc_min soc(t1) soc_max]; end % 一个运行日下来SOC回到初始值 Constraints [Constraints, soc(T1) soc(1)];目标函数与求解Objective sum(Pbuy .* c_buy) sum(Pmt .* c_gas) ... c_pv_curt * sum(Ppv_curt) ... % 弃光惩罚 c_wt_curt * sum(Pwt_curt); % 弃风惩罚 ops sdpsettings(solver,gurobi,verbose,2,mipgap,1e-4); sol optimize(Constraints, Objective, ops); if sol.problem ~ 0 disp(sol.info); end代码里的Ppv_curt和Pwt_curt是弃光弃风变量正式建模时需要把它们加入光伏和风机的功率平衡约束形如Ppv_real(t) Ppv_curt(t) Ppv_forecast(t)。这里为了控制代码篇幅没有完整展开但逻辑上是必须的。使用YALMIP时有个维度问题要特别小心1×24和24×1混用经常报维度不一致建议所有一维变量统一用1×T不要一会行向量一会列向量。3.4 模型求解时间与性能调优纯LP模型没有整数变量在33节点24时段规模下Gurobi通常几秒内解完。如果加入燃气轮机启停变量变成MILP求解时间可能升到几十秒甚至几分钟。几个性能调优经验约束尽量一次性用矩阵拼接构建。在循环里不断追加ConstraintsYALMIP每次append都要重建符号表达式循环次数一多会明显变慢。先跑LP松弛版本验证物理逻辑再加整数变量不要一开始就上MILP。mipgap设1e-4足够工程使用默认1e-9很多时候只是为了最后那一点精度多等几十分钟。如果求解时间实在压不下来把日前时间粒度从1小时改成2小时或把日内滚动窗口从8小时压到4小时速度会显著提升。4. 结果怎么读曲线、指标与方案对比4.1 必画的六类调度结果图两阶段模型跑完之后别直接看一眼目标值就收工。我建议把下面这些图全部画出来逐张检查。第一张是净负荷曲线。把负荷减去光伏和风电预测出力得到的净负荷画出来能一眼看出一天里哪个时段系统最紧张哪个时段需要储能充电。第二张是联络线交换功率对比图把第一阶段的日前计划和第二阶段修正后的实际值画在同一张图上偏差大的时段就是预测误差影响最明显的时段。第三张是DG出力图光伏、风电、燃机的出力逐时段画出来重点看大风和强光时段有没有明显弃风弃光。第四张是储能SOC曲线和充放电功率图用来验证储能是否真的在低谷充电、高峰放电有没有出现整天都在放电的不合理结果。第五张是电压包络图把24个时段每个节点最高和最低电压画成上下两条包络线一旦碰到0.95或1.05的限值线就说明那个节点或时段需要重点处理。第六张是系统总网损曲线如果模型包含网损项这条曲线能反映整体运行质量。4.2 单阶段与两阶段的量化差距我拿一个典型的33节点算例做过对比结果大致如下。具体数字和参数强相关主要看趋势指标单阶段确定性模型两阶段优化模型运行总成本含日内调整成本基准值低8%~12%24小时内电压越限次数5~8次0~1次弃风弃光率15%左右5%以下联络线功率波动幅度较大明显缓和单阶段模型不是成本一定更高而是它的成本里根本没有计入日内实时调整的偏差惩罚和越限风险。一旦把实际运行中的不平衡成本算进去两阶段模型的优势就很明显。另一个关键改善是电压质量因为第二阶段能根据最新预测提前调整储能和燃机出力电压越限次数大幅减少。4.3 预测误差对两阶段效果的影响把预测误差的标准差从5%逐步加到20%会看到一个规律误差越大两阶段模型相对单阶段的成本优势越明显。原因很直接误差越大日内需要修正的空间越大固定日前计划的代价越高。但这里有个边界。当误差大到储能和燃气轮机的可调容量都覆盖不了时两阶段模型同样会出现电压越限或切负荷只是比单阶段模型来得晚一些。所以做对比实验时建议用Monte Carlo方式跑多组误差场景取平均结果再比较不要单跑一两次就下结论。只跑一次很可能因为随机数运气好或差得出与真实规律相反的结论。5. 调试心得那些让求解器“翻车”的细节5.1 不可行解按这个顺序排查YALMIP返回“Infeasible problem”时最常见的原因有三个。维度错配排第一。很多新手把sdpvar定义成1×T却用T×1的数组初始化YALMIP会静默生成一个非预期维度的变量然后约束比较时直接报错。排查办法是在构建约束前打印size(变量)确认。储能约束自相矛盾排第二。比如soc_min设成0.1、充电功率上限又小、充电时间短同时要求首末SOC相等这种约束无论怎么优化都找不到可行域。联络线功率上限设得太死排第三。储能同时充电加负荷高峰的情况下购电功率很容易突破上级变压器容量上限。建议的调试顺序是先注释掉全部约束然后逐个打开第一次出现不可行的那条约束就是问题所在。这个过程虽然有点笨但对中小规模算例非常快比盯着报错信息猜高效得多。5.2 储能SOC周期耦合一个让结果“白嫖”电量的Bug如果第一阶段模型里少了soc(T1) soc(1)这条约束优化器通常会做一件事在一天结束前把SOC放空把存储的电能全部换算成电费收益因为成本模型里没有为“明天还要继续运行”买单。拿到结果的人会看到储能一天到晚都在放电SOC曲线一路跌到下限表面收益很漂亮实际根本没法滚动运行。解决的办法除了加首末耦合约束也可以把SOC终止值设为软约束让末时刻SOC尽量接近初始值偏离时加惩罚。这样做的好处是给日内运行留一点灵活性而不是把终止SOC死死钉在初始值上尤其适合日内预测偏差较大的场景。5.3 第二阶段不动作惩罚系数设置不当第二阶段模型跑出来所有ΔP变量全是0目标函数值和第一阶段完全一样。这不是代码bug而是偏差惩罚系数设得太小。优化器算了一笔账与其花燃气轮机的燃料成本和储能损耗去修正偏差不如什么都不做、交点惩罚费更便宜。解决办法是把α⁺和α⁻从实时电价差值往上调至少要高于“实时购电价格减日前购电价格”才会触发修正动作。更稳妥的做法是先人为制造一个大的预测偏差用测试场景验证第二阶段会不会动作确认修正逻辑正常后再放回真实数据。这一步测试虽然简单但我见过不少同行跳过它最后花很长时间才定位到是系数问题。5.4 弃风弃光惩罚系数的取值经验弃风弃光惩罚系数不能乱设。设得太低优化器会为了降低运行成本选择大量弃风弃光设得太高比如1e9优化器又会为了规避天价惩罚去压其他约束产生不合理的运行策略。我一般把弃光惩罚系数取在当前实时电价的0.5到2倍之间或者参考度电成本上限来确定然后根据结果微调。做敏感性分析时把惩罚系数从低到高扫描一遍看新能源利用率的变化曲线就能找到一个比较合理的取值区间。最后再分享一个我自己坚持到现在的习惯任何两阶段模型跑通之后别急着拿结果写报告先把过去某几天的预测数据和实际数据拿出来做一次回放测试。用日前预测跑第一阶段再用实际数据模拟第二阶段修正逐时段核对系统是否越限、成本是否和模型预测值接近。这个习惯帮我抓出过不少只有在数据回放时才会暴露的建模问题。如果你也是第一次搭这类模型强烈建议把回放测试当成验收环节而不是可选项。
