1. 项目概述与两阶段调度思路1.1 为什么配电网调度要分“两阶段”先聊点背景。以前搞配电网调度大家面对的基本是单向潮流的被动网络电源侧就是变电站、馈线出口负荷曲线虽然预测不准但整体波动可控。但现在分布式电源大量接入风电、光伏的出力谁都知道难伺候——早上的光伏爬坡、傍晚的骤降、台风天前后风电的忽大忽小这还只是天气层面的不确定。再加上负荷预测本身就有误差如果还按传统的“一条预测曲线定全天计划”的做法要么计划过于保守导致经济性差要么计划过于激进导致实际运行时不安全。这就是“两阶段”调度模型出现的核心原因把决策拆成“日前阶段”和“日内/实时阶段”两层来处理。第一阶段在日前做决策使用预测数据提前决定机组启停、储能充放电计划、联络线交换功率等“慢变量”第二阶段在日内根据实际/更准确的风光出力修正在已有日前计划的基础上做“再调度”处理“快变量”和偏差量。一句话概括第一阶段定骨架第二阶段补肌肉这样既保证了经济性又留出了应对不确定性的调节空间。1.2 这个模型到底解决什么问题在含分布式电源的配电网里光伏和风电接入会带来三个典型问题潮流反向中午光伏大发时馈线末端电压可能越上限甚至向上一级电网倒送功率。电压波动分布式电源出力的随机性直接反映在节点电压上尤其在弱电网末端电压闪变问题突出。经济性退化如果调度方式跟不上分布式电源的消纳率低弃光弃风严重投资回报周期变长。两阶段优化调度模型正是为了在经济性购电成本、网损、储能损耗与安全性电压约束、支路潮流约束之间找最优平衡点。具体来说日前阶段以预测场景下的总运行成本最小为目标通过优化常规机组出力、储能充放电策略、联络线功率等先给出一个基准方案日内阶段以实际场景或更精细的预测场景为准最小化调整成本在基准方案基础上做修正确保安全约束始终满足。这套逻辑从思路上讲并不晦涩真正难的是建模细节和Matlab实现时的各种坑。2. 分布式电源与不确定性建模2.1 风光出力的数学建模分布式电源最典型的就是光伏PV和风电WT。做日前调度时第一步就是把它们的出力预测曲线转换成可计算的数学模型。光伏出力的简化模型通常基于光照强度[ P_{PV}(t) P_{STC} \cdot \frac{G(t)}{G_{STC}} \cdot [1 k(T(t) - T_{STC})] ]其中 (P_{STC}) 是标准测试条件下的额定功率(G(t)) 是实际光照强度(G_{STC}) 是标准光照强度通常取1000 W/m²(k) 是温度系数一般取 -0.0045 左右(T(t)) 是光伏板表面温度。实际工程中如果你没有光照和温度数据也可以用典型日出力曲线加随机扰动来近似——但扰动必须符合一定的概率分布否则后续场景分析就失真了。风电出力模型则基于风速与功率的转换关系[ P_{WT}(v) \begin{cases} 0, v v_{in} \text{ 或 } v v_{out} \ \frac{v - v_{in}}{v_{rated} - v_{in}} P_{rated}, v_{in} \leq v \leq v_{rated} \ P_{rated}, v_{rated} v \leq v_{out} \end{cases} ]风速数据通常用Weibull分布描述(v_{in}) 是切入风速、(v_{out}) 是切出风速、(v_{rated}) 是额定风速。这个分段函数看着简单但在Matlab里向量化实现时要注意索引边界尤其是风速恰好落在临界点附近时容易出现数组越界或除零问题。我的建议是单独写一个wind_power_curve.m函数用逻辑索引一次性处理不要用循环逐点判断速度会快很多倍。2.2 场景生成与削减从“一条曲线”到“一组场景”既然预测不准干脆别假装准。最常用的做法是用蒙特卡洛采样生成大量可能的风光出力场景然后用场景削减技术挑出概率最大、最具代表性的少数几个场景参与优化——这比直接对预测值做鲁棒优化计算量小得多也比单场景优化结果更稳健。场景生成的主要步骤对历史预测误差做统计分析得到误差的概率分布通常假设为正态分布或Beta分布。对每个时刻的预测值叠加随机误差生成N个完整的日出力场景。用同步回代消除法scenario reduction或K-means聚类把N个场景削减成K个典型场景K一般取5~10个。同步回代法的核心逻辑是每次迭代找到距离最近的一对场景将其中概率较小的一个删除并将其概率累加到较大概率场景上直到剩下指定数量的场景。Matlab里可以实现如下function [scen_reduced, prob_reduced] scenario_reduction(scen, prob, K) % scen: 原始场景矩阵, 每列一个场景 % prob: 每个场景的概率 % K: 目标场景数量 n_scen size(scen, 2); dist zeros(n_scen); for i 1:n_scen for j i1:n_scen dist(i,j) norm(scen(:,i) - scen(:,j)); dist(j,i) dist(i,j); end end while n_scen K % 找距离最小的场景对 min_dist inf; for i 1:n_scen for j i1:n_scen if dist(i,j) min_dist min_dist dist(i,j); idx_pair [i, j]; end end end % 删除概率较小的那个概率加到另一个上 [~, idx_del] min(prob(idx_pair)); idx_keep idx_pair(3 - idx_del); prob(idx_keep) prob(idx_keep) prob(idx_del); % 删除场景 scen(:, idx_del) []; prob(idx_del) []; dist(:, idx_del) []; dist(idx_del, :) []; n_scen n_scen - 1; end scen_reduced scen; prob_reduced prob / sum(prob); end这段代码在场景数量大比如1000个时计算较慢主要是两层循环算距离矩阵太耗时。实践中可以先随机抽样50个初筛场景再做削减效果差别不大但速度提升明显。3. 两阶段优化调度模型构建3.1 第一阶段日前调度主问题第一阶段的决策变量是常规机组启停状态、各时段出力、储能充放电功率与状态、联络线交换功率、分布式电源的日前计划出力。目标函数是让全天总成本最小[ \min ; C \sum_{t1}^{T} \left[ \sum_{g \in G} (a_g P_{g,t}^2 b_g P_{g,t} c_g) \lambda_t^{grid} P_{t}^{grid} C_{ESS,t} \right] ]其中 (T24)第一项是常规机组的燃料成本通常简化为二次函数第二项是从上级电网购电的成本(\lambda_t^{grid}) 是分时电价第三项是储能充放电损耗成本。这个目标函数没有考虑弃风弃光惩罚如果项目要求优先消纳分布式电源需要在目标里增加一项[ C_{curtail} \sum_{t1}^{T} \sum_{d \in DG} c_{curtail} \cdot (P_{d,t}^{forecast} - P_{d,t}^{schedule}) ]惩罚系数 (c_{curtail}) 应大于购电电价否则优化器会为了省钱而主动弃风弃光违背了分布式电源优先消纳的初衷这是我做项目时踩过的坑调参时需要格外注意。第一阶段的约束包括潮流约束DistFlow 线性化方程或交流潮流方程节点电压上下限支路电流/功率上限常规机组出力上下限与爬坡约束储能SOC约束与充放电功率约束联络线交换功率约束3.2 第二阶段日内再调度子问题第二阶段的思路是第一阶段的日前计划已经确定但当实际场景或更精确的预测场景出来后某些约束可能不满足了比如光伏实际出力比预测低导致部分负荷需要切掉或电压越限。此时需要在“最小化调整量”的目标下重新分配出力[ \min ; \sum_{t1}^{T} \left( \sum_{g \in G} c_g^ \Delta P_{g,t}^ c_g^- \Delta P_{g,t}^- \sum_{d \in DG} c_d^{curtail} \Delta P_{d,t}^{curtail} \right) ]式中 (\Delta P^) 和 (\Delta P^-) 分别表示常规机组的向上/向下调整量(\Delta P^{curtail}) 是分布式电源的削减量。注意 (\Delta P^) 和 (\Delta P^-) 不能同时非零这需要引入辅助变量或约束处理。我用YALMIP建模时直接把它写成两个非负变量再加一个互斥约束 (x_ x_- \leq 1)其中 (x_) 和 (x_-) 是0-1变量。不过这样会增加求解难度如果对精度要求不高也可以用大M法把互斥约束松弛掉——因为目标函数是正成本优化器自然倾向于不同时调整。第二阶段的目标函数除了调整成本也可以加权一个安全约束越限惩罚项比如电压越限的软约束成本这样可以在极端场景下避免无解。这个技巧在工程中非常实用因为完全满足所有硬约束的可行解不一定存在软约束可以保证求解器总能返回一个可用的次优解。3.3 配电网潮流约束的线性化处理配电网潮流计算比输电网复杂的地方在于网络拓扑通常是辐射状R/X比值较大传统的PQ分解法几乎不收敛必须用前推回代法或者牛拉法。但在优化模型里直接嵌入非线性潮流方程会导致模型变成MINLP求解极其困难。工程上最主流的做法是采用 DistFlow 分支潮流方程并对其做二阶锥松弛SOCP或线性化近似。DistFlow 方程如下考虑有功、无功和电压降[ P_{ij} \sum_{k \in N(j)} P_{jk} r_{ij} \cdot \frac{P_{ij}^2 Q_{ij}^2}{V_i^2} P_{j}^{load} - P_{j}^{DG} ][ Q_{ij} \sum_{k \in N(j)} Q_{jk} x_{ij} \cdot \frac{P_{ij}^2 Q_{ij}^2}{V_i^2} Q_{j}^{load} - Q_{j}^{DG} ][ V_j^2 V_i^2 - 2(r_{ij}P_{ij} x_{ij}Q_{ij}) (r_{ij}^2 x_{ij}^2) \cdot \frac{P_{ij}^2 Q_{ij}^2}{V_i^2} ]这里面的非线性项 ( \frac{P_{ij}^2 Q_{ij}^2}{V_i^2} ) 是难点。二阶锥松弛的做法是引入辅助变量 (l_{ij} \frac{P_{ij}^2 Q_{ij}^2}{V_i^2})将约束松弛为[ | \begin{bmatrix} 2P_{ij} \ 2Q_{ij} \ l_{ij} - V_i^2 \end{bmatrix} |2 \leq l{ij} V_i^2 ]这样就把非线性潮流转化成了一个凸约束整个模型变成了混合整数二阶锥规划MISOCPYALMIP配合Gurobi或Mosek都能高效求解。如果不追求精确潮流结果只做日前调度方案对比也可以直接用线性化潮流[ P_{ij} \approx \sum_{k \in N(j)} P_{jk} P_{j}^{load} - P_{j}^{DG} ][ V_j \approx V_i - r_{ij}P_{ij} - x_{ij}Q_{ij} ]这种方法忽略了网损项在最重负载场景下误差可能达到百分之几但在优化方案预筛选阶段完全够用。我在仿真里对比过SOCP方案和线性化方案在最终调度结果上的差别通常不到2%但SOCP的求解时间可能是线性化的几十倍。如果只是做演示或快速验证线性化足够如果要做精密分析或发表论文建议用SOCP。4. Matlab代码实现与求解配置4.1 模型初始化与数据准备Matlab里实现这类优化调度模型我的标准套路如下基础数据节点、线路、负荷存成结构体数组单独写一个脚本load_case_data.m分布式电源参数和预测曲线存成.mat文件用load命令读入优化模型用 YALMIP 建模求解器用 Gurobi如果只有Mosek也可以但Mosek对二阶锥支持不如Gurobi顺手数据准备阶段比较容易忽略的是电价的时段划分。很多论文里的分时电价用的是峰、平、谷三段但实际仿真时如果用的是真实电网数据电价可能是半小时一个点甚至15分钟一个点。这会给日前调度带来一个问题模型的时间尺度到底取多少常见做法是日前阶段取1小时分辨率日内阶段取15分钟分辨率两阶段之间通过滑动窗口对接。如果模型只是演示性质统一用1小时分辨率也可以但要在论文里说明。4.2 YALMIP建模核心代码以第一阶段为例YALMIP建模关键代码如下% 定义变量 P_g sdpvar(n_g, T, full); % 常规机组出力 u_g binvar(n_g, T, full); % 机组启停状态 P_ch sdpvar(n_ess, T, full); % 储能充电功率 P_dis sdpvar(n_ess, T, full); % 储能放电功率 SOC sdpvar(n_ess, T, full); % 储能的荷电状态 P_grid sdpvar(1, T, full); % 联络线交换功率 % 目标函数 objective 0; for t 1:T objective objective sum(a_g * P_g(:,t).^2 b_g * P_g(:,t) c_g * u_g(:,t)); objective objective grid_price(t) * P_grid(t); objective objective sum(C_ess * (P_ch(:,t) P_dis(:,t))); end % 约束条件集合 constraints []; % 机组出力上下限 constraints [constraints, P_min .* u_g P_g P_max .* u_g]; % 储能SOC递推 for t 2:T constraints [constraints, SOC(:,t) SOC(:,t-1) eta_ch * P_ch(:,t) - P_dis(:,t)/eta_dis]; end % 功率平衡 for t 1:T constraints [constraints, sum(P_g(:,t)) sum(P_dis(:,t)) - sum(P_ch(:,t)) P_grid(t) sum(P_wt(:,t)) sum(P_pv(:,t)) total_load(t)]; end % 求解 ops sdpsettings(solver, gurobi, verbose, 2); optimize(constraints, objective, ops);有几个细节要提醒新手sum(P_wt(:,t)) sum(P_pv(:,t))用的是预测值不是决策变量。如果要做两阶段这里应该把分布式电源出力也设为决策变量并受到预测值上限约束。储能SOC公式里eta_ch和eta_dis分开设置因为充放电效率不同。有些简化模型直接用一个效率但实际电池在充和放时损耗确实不一样尤其在SOC较高时充电效率下降明显。目标函数里的P_g(:,t).^2是非线性项如果Gurobi处理不了二次目标MIQP需要先把目标线性化。一个常见做法是用分段线性近似替代二次成本函数YALMIP里可以用pwf函数处理不过直接用Gurobi的MIQP求解器通常也没问题。4.3 第二阶段的模型与衔接处理第二阶段的基本思路是把第一阶段确定的机组状态 (u_g^*)、储能基准功率 (P_{ess}^{base})、分布式电源基准出力 (P_{dg}^{base}) 当作已知参数然后在新场景下求解再调度模型。% 第二阶段场景s下做再调度 % 固定第一阶段决策中的整数变量 constraints_2nd [constraints_2nd, u_g u_g_star]; % 整数变量保持不变 % 再调度调整量 delta_Pg_up sdpvar(n_g, T, full); delta_Pg_down sdpvar(n_g, T, full); delta_Pdg sdpvar(n_dg, T, full); % 分布式电源削减量 % 功率平衡场景s下的实际值 for t 1:T constraints_2nd [constraints_2nd, ... sum(P_g_star delta_Pg_up(:,t) - delta_Pg_down(:,t)) ... sum(P_dis_star(:,t) - P_ch_star(:,t)) ... % 这里简化处理储能 P_grid_2nd(t) sum(P_wt_scen(:,t)) - sum(delta_Pdg_wt(:,t)) ... sum(P_pv_scen(:,t)) - sum(delta_Pdg_pv(:,t)) total_load(t)]; end这里的关键是“固定整数变量”。在Matlab中直接把u_g替换成u_g_star常数值即可不要让求解器再优化启停。否则第二阶段变成一个完整的MIP问题那就失去了两阶段“快速再调度”的意义。第二阶段求解完成后把各场景下的调整成本加权求和加上第一阶段的基准成本就是这个日前调度方案的期望总成本。至此两阶段模型的计算闭环完成。4.4 求解器选型与参数调优Matlab里求解这类问题可选的求解器不少我的经验是Gurobi对于MISOCP和MIQP速度和稳定性都没得说学术免费强烈推荐。Mosek对锥规划的数值稳定性比Gurobi略好但整数变量的MIP能力稍弱。CPLEX老牌求解器目前新版本对配电网优化支持也不错但在学生群体里用得少了。Matlab内置的求解器如intlinprog只能解MILP处理不了SOCP。如果只是做线性化潮流可以用一旦涉及非线性潮流约束就束手无策了。调参方面最常碰到的坑是sdpsettings(solver, gurobi)之后YALMIP有时会把变量自动转换为其他格式导致Gurobi报错 “Model is infeasible”。这种情况下我一般是先检查约束里是不是存在矛盾条件比如某个节点既要求电压不超过1.0又要求注入功率满足某个下限两者可能互斥。把约束逐个注释掉跑一次模型找到引起不可行的那一组约束通常很快就能定位问题。5. 常见问题与排查技巧实录5.1 模型求解速度慢怎么优化两阶段模型如果直接用全网所有节点建立约束负荷节点多了后约束数量会爆炸式增长。比如一个33节点系统24小时每个时段的潮流约束就有几十条加上机组、储能、配网约束变量总数轻松上万。求解MISOCP可能要几分钟甚至更久。优化手段用场景削减缩小第二阶段场景数量10个场景以内是比较合理的。电压约束用软约束只要在目标函数中惩罚越限就可以改用单纯形快速求解。如果只是研究调度策略而不关注潮流细节可以把配电网等效为一个大节点只保留馈线出口的功率平衡约束求解速度能提升几个数量级。我这里有一个实际案例某个33节点的配电网模型初始建模用了完整DistFlow加SOCP松弛加上8个场景Gurobi求解耗时约200秒。后来把分布式电源接入点压缩到3个关键节点潮流网络做了等效简化求解时间降到8秒而调度成本只多了不到3%。对前期方案比选来说这个效率提升非常值得。5.2 出现 “Infeasible problem” 怎么排查这个问题我在教学中被问得最多。排查思路很重要一定要按顺序来第一检查功率平衡约束。很多新手在写功率平衡时忘了把线路损耗、储能自损耗算进去导致每个时段的总发电大于总负荷或小于无解。可以先用线性潮流把网损忽略看看模型是否可解如果可解再逐步加入网损。第二检查储能SOC递推公式。SOC的初值如果不合理比如要求首末SOC相等但电池容量又不够大就会导致无解。我一般设置首末SOC相同并给出一个较小的充放电功率上限避免约束过紧。第三检查爬坡约束。机组爬坡能力设置过小而负荷波动又大可能会导致特定时段无法满足功率平衡。把爬坡约束系数放大10倍如果模型变可解了说明问题出在这儿。第四检查电压约束。分布式电源接入容量越大电压越限越容易发生。如果约束是 0.95 ≤ V ≤ 1.05可以先放宽到 0.9~1.1 试试如果能解再逐渐收紧。5.3 分布式电源渗透率过高时调度结果不合理有时候模型求解正常但结果里会看到在分布式电源出力高、负荷轻的时段联络线功率为0甚至为负数向上一级电网倒送但常规机组仍然以最小出力运行。这看起来反直觉实际上是目标函数里没有给常规机组设置启停惩罚导致的。解决方案是增加机组的启停成本或者给常规机组设置最小运行时间约束。否则优化器会在某个时段让机组停机以降低成本但下一时段又需要它启动了启停太频繁工程上完全不可接受。这个问题在IEEE 33节点系统里非常典型我建议直接加入最小启停时间约束即使求解时间增加一点也值得。5.4 常见错误速查表下面这张表是我做Matlab配电网调度模型时总结出来的高频问题分享给大家现象可能原因处理方法求解器提示“Infeasible”功率平衡约束不平衡检查发电/负荷/网损是否闭合储能SOC越界递推公式效率参数方向错误确认充电/放电效率是否对应正确方向优化结果里机组频繁启停缺少启停成本或最小运行时间约束在目标中加入启停惩罚项电压越限但优化器不处理电压约束写成了软约束但惩罚系数太小增大惩罚权重或改为硬约束第二阶段无解第一阶段基准方案太激进放宽部分约束或加入松弛变量和惩罚项求解时间过长场景数过多或SOCP约束过多削减场景简化潮流为线性近似MATLAB报错“Undefined function”工作路径没有包含函数文件检查当前文件夹及路径设置5.5 代码调试的独家技巧最后分享几个实战小技巧这些是常规教程里不会细讲的第一个技巧是写一个check_constraints.m函数在优化求解前把所有约束的残差打印出来这样能非常直观地判断是哪个约束导致无解。我见过不少人用YALMIP时只盯着优化器返回状态结果状态码又不明确白白浪费了好几个小时。其实只需在optimize之后加一行check(constraints)YALMIP 就会逐一列出所有约束的 OK 状态和残差排查效率瞬间翻倍。第二个技巧是在模型测试阶段先固定所有整数变量只用测试数据跑连续优化验证连续模型的可行性再去处理整数变量的MIP问题。这样可以快速区分是数学模型本身的问题还是整数变量的组合爆炸问题。第三个技巧是注意YALMIP和Gurobi版本兼容性。我踩过一次坑YALMIP版本太老传给Gurobi的二阶锥约束在求解器里被识别成了非凸二次约束结果求解时间增加了10倍换了新版本YALMIP后立刻恢复。如果你的模型突然变慢优先怀疑版本匹配问题。从我个人体会来说做这类含分布式电源的配电网两阶段优化调度模型真正难的不是算法原理而是把原理落地的过程中那些零零碎碎的细节——场景怎么生成、约束怎么松弛、求解器怎么调、无解时怎么定位。这些经验不是看几篇论文就能获得的需要在反复的调试中慢慢积累。希望这篇文章能帮你少走一些弯路把时间花在真正有价值的问题上。 
长河君 · 2026/9/19 5:28:37 · 阅读 1.2k 配电网分布式电源两阶段优化调度模型及Matlab实现详解
1. 项目概述与两阶段调度思路1.1 为什么配电网调度要分“两阶段”先聊点背景。以前搞配电网调度大家面对的基本是单向潮流的被动网络电源侧就是变电站、馈线出口负荷曲线虽然预测不准但整体波动可控。但现在分布式电源大量接入风电、光伏的出力谁都知道难伺候——早上的光伏爬坡、傍晚的骤降、台风天前后风电的忽大忽小这还只是天气层面的不确定。再加上负荷预测本身就有误差如果还按传统的“一条预测曲线定全天计划”的做法要么计划过于保守导致经济性差要么计划过于激进导致实际运行时不安全。这就是“两阶段”调度模型出现的核心原因把决策拆成“日前阶段”和“日内/实时阶段”两层来处理。第一阶段在日前做决策使用预测数据提前决定机组启停、储能充放电计划、联络线交换功率等“慢变量”第二阶段在日内根据实际/更准确的风光出力修正在已有日前计划的基础上做“再调度”处理“快变量”和偏差量。一句话概括第一阶段定骨架第二阶段补肌肉这样既保证了经济性又留出了应对不确定性的调节空间。1.2 这个模型到底解决什么问题在含分布式电源的配电网里光伏和风电接入会带来三个典型问题潮流反向中午光伏大发时馈线末端电压可能越上限甚至向上一级电网倒送功率。电压波动分布式电源出力的随机性直接反映在节点电压上尤其在弱电网末端电压闪变问题突出。经济性退化如果调度方式跟不上分布式电源的消纳率低弃光弃风严重投资回报周期变长。两阶段优化调度模型正是为了在经济性购电成本、网损、储能损耗与安全性电压约束、支路潮流约束之间找最优平衡点。具体来说日前阶段以预测场景下的总运行成本最小为目标通过优化常规机组出力、储能充放电策略、联络线功率等先给出一个基准方案日内阶段以实际场景或更精细的预测场景为准最小化调整成本在基准方案基础上做修正确保安全约束始终满足。这套逻辑从思路上讲并不晦涩真正难的是建模细节和Matlab实现时的各种坑。2. 分布式电源与不确定性建模2.1 风光出力的数学建模分布式电源最典型的就是光伏PV和风电WT。做日前调度时第一步就是把它们的出力预测曲线转换成可计算的数学模型。光伏出力的简化模型通常基于光照强度[ P_{PV}(t) P_{STC} \cdot \frac{G(t)}{G_{STC}} \cdot [1 k(T(t) - T_{STC})] ]其中 (P_{STC}) 是标准测试条件下的额定功率(G(t)) 是实际光照强度(G_{STC}) 是标准光照强度通常取1000 W/m²(k) 是温度系数一般取 -0.0045 左右(T(t)) 是光伏板表面温度。实际工程中如果你没有光照和温度数据也可以用典型日出力曲线加随机扰动来近似——但扰动必须符合一定的概率分布否则后续场景分析就失真了。风电出力模型则基于风速与功率的转换关系[ P_{WT}(v) \begin{cases} 0, v v_{in} \text{ 或 } v v_{out} \ \frac{v - v_{in}}{v_{rated} - v_{in}} P_{rated}, v_{in} \leq v \leq v_{rated} \ P_{rated}, v_{rated} v \leq v_{out} \end{cases} ]风速数据通常用Weibull分布描述(v_{in}) 是切入风速、(v_{out}) 是切出风速、(v_{rated}) 是额定风速。这个分段函数看着简单但在Matlab里向量化实现时要注意索引边界尤其是风速恰好落在临界点附近时容易出现数组越界或除零问题。我的建议是单独写一个wind_power_curve.m函数用逻辑索引一次性处理不要用循环逐点判断速度会快很多倍。2.2 场景生成与削减从“一条曲线”到“一组场景”既然预测不准干脆别假装准。最常用的做法是用蒙特卡洛采样生成大量可能的风光出力场景然后用场景削减技术挑出概率最大、最具代表性的少数几个场景参与优化——这比直接对预测值做鲁棒优化计算量小得多也比单场景优化结果更稳健。场景生成的主要步骤对历史预测误差做统计分析得到误差的概率分布通常假设为正态分布或Beta分布。对每个时刻的预测值叠加随机误差生成N个完整的日出力场景。用同步回代消除法scenario reduction或K-means聚类把N个场景削减成K个典型场景K一般取5~10个。同步回代法的核心逻辑是每次迭代找到距离最近的一对场景将其中概率较小的一个删除并将其概率累加到较大概率场景上直到剩下指定数量的场景。Matlab里可以实现如下function [scen_reduced, prob_reduced] scenario_reduction(scen, prob, K) % scen: 原始场景矩阵, 每列一个场景 % prob: 每个场景的概率 % K: 目标场景数量 n_scen size(scen, 2); dist zeros(n_scen); for i 1:n_scen for j i1:n_scen dist(i,j) norm(scen(:,i) - scen(:,j)); dist(j,i) dist(i,j); end end while n_scen K % 找距离最小的场景对 min_dist inf; for i 1:n_scen for j i1:n_scen if dist(i,j) min_dist min_dist dist(i,j); idx_pair [i, j]; end end end % 删除概率较小的那个概率加到另一个上 [~, idx_del] min(prob(idx_pair)); idx_keep idx_pair(3 - idx_del); prob(idx_keep) prob(idx_keep) prob(idx_del); % 删除场景 scen(:, idx_del) []; prob(idx_del) []; dist(:, idx_del) []; dist(idx_del, :) []; n_scen n_scen - 1; end scen_reduced scen; prob_reduced prob / sum(prob); end这段代码在场景数量大比如1000个时计算较慢主要是两层循环算距离矩阵太耗时。实践中可以先随机抽样50个初筛场景再做削减效果差别不大但速度提升明显。3. 两阶段优化调度模型构建3.1 第一阶段日前调度主问题第一阶段的决策变量是常规机组启停状态、各时段出力、储能充放电功率与状态、联络线交换功率、分布式电源的日前计划出力。目标函数是让全天总成本最小[ \min ; C \sum_{t1}^{T} \left[ \sum_{g \in G} (a_g P_{g,t}^2 b_g P_{g,t} c_g) \lambda_t^{grid} P_{t}^{grid} C_{ESS,t} \right] ]其中 (T24)第一项是常规机组的燃料成本通常简化为二次函数第二项是从上级电网购电的成本(\lambda_t^{grid}) 是分时电价第三项是储能充放电损耗成本。这个目标函数没有考虑弃风弃光惩罚如果项目要求优先消纳分布式电源需要在目标里增加一项[ C_{curtail} \sum_{t1}^{T} \sum_{d \in DG} c_{curtail} \cdot (P_{d,t}^{forecast} - P_{d,t}^{schedule}) ]惩罚系数 (c_{curtail}) 应大于购电电价否则优化器会为了省钱而主动弃风弃光违背了分布式电源优先消纳的初衷这是我做项目时踩过的坑调参时需要格外注意。第一阶段的约束包括潮流约束DistFlow 线性化方程或交流潮流方程节点电压上下限支路电流/功率上限常规机组出力上下限与爬坡约束储能SOC约束与充放电功率约束联络线交换功率约束3.2 第二阶段日内再调度子问题第二阶段的思路是第一阶段的日前计划已经确定但当实际场景或更精确的预测场景出来后某些约束可能不满足了比如光伏实际出力比预测低导致部分负荷需要切掉或电压越限。此时需要在“最小化调整量”的目标下重新分配出力[ \min ; \sum_{t1}^{T} \left( \sum_{g \in G} c_g^ \Delta P_{g,t}^ c_g^- \Delta P_{g,t}^- \sum_{d \in DG} c_d^{curtail} \Delta P_{d,t}^{curtail} \right) ]式中 (\Delta P^) 和 (\Delta P^-) 分别表示常规机组的向上/向下调整量(\Delta P^{curtail}) 是分布式电源的削减量。注意 (\Delta P^) 和 (\Delta P^-) 不能同时非零这需要引入辅助变量或约束处理。我用YALMIP建模时直接把它写成两个非负变量再加一个互斥约束 (x_ x_- \leq 1)其中 (x_) 和 (x_-) 是0-1变量。不过这样会增加求解难度如果对精度要求不高也可以用大M法把互斥约束松弛掉——因为目标函数是正成本优化器自然倾向于不同时调整。第二阶段的目标函数除了调整成本也可以加权一个安全约束越限惩罚项比如电压越限的软约束成本这样可以在极端场景下避免无解。这个技巧在工程中非常实用因为完全满足所有硬约束的可行解不一定存在软约束可以保证求解器总能返回一个可用的次优解。3.3 配电网潮流约束的线性化处理配电网潮流计算比输电网复杂的地方在于网络拓扑通常是辐射状R/X比值较大传统的PQ分解法几乎不收敛必须用前推回代法或者牛拉法。但在优化模型里直接嵌入非线性潮流方程会导致模型变成MINLP求解极其困难。工程上最主流的做法是采用 DistFlow 分支潮流方程并对其做二阶锥松弛SOCP或线性化近似。DistFlow 方程如下考虑有功、无功和电压降[ P_{ij} \sum_{k \in N(j)} P_{jk} r_{ij} \cdot \frac{P_{ij}^2 Q_{ij}^2}{V_i^2} P_{j}^{load} - P_{j}^{DG} ][ Q_{ij} \sum_{k \in N(j)} Q_{jk} x_{ij} \cdot \frac{P_{ij}^2 Q_{ij}^2}{V_i^2} Q_{j}^{load} - Q_{j}^{DG} ][ V_j^2 V_i^2 - 2(r_{ij}P_{ij} x_{ij}Q_{ij}) (r_{ij}^2 x_{ij}^2) \cdot \frac{P_{ij}^2 Q_{ij}^2}{V_i^2} ]这里面的非线性项 ( \frac{P_{ij}^2 Q_{ij}^2}{V_i^2} ) 是难点。二阶锥松弛的做法是引入辅助变量 (l_{ij} \frac{P_{ij}^2 Q_{ij}^2}{V_i^2})将约束松弛为[ | \begin{bmatrix} 2P_{ij} \ 2Q_{ij} \ l_{ij} - V_i^2 \end{bmatrix} |2 \leq l{ij} V_i^2 ]这样就把非线性潮流转化成了一个凸约束整个模型变成了混合整数二阶锥规划MISOCPYALMIP配合Gurobi或Mosek都能高效求解。如果不追求精确潮流结果只做日前调度方案对比也可以直接用线性化潮流[ P_{ij} \approx \sum_{k \in N(j)} P_{jk} P_{j}^{load} - P_{j}^{DG} ][ V_j \approx V_i - r_{ij}P_{ij} - x_{ij}Q_{ij} ]这种方法忽略了网损项在最重负载场景下误差可能达到百分之几但在优化方案预筛选阶段完全够用。我在仿真里对比过SOCP方案和线性化方案在最终调度结果上的差别通常不到2%但SOCP的求解时间可能是线性化的几十倍。如果只是做演示或快速验证线性化足够如果要做精密分析或发表论文建议用SOCP。4. Matlab代码实现与求解配置4.1 模型初始化与数据准备Matlab里实现这类优化调度模型我的标准套路如下基础数据节点、线路、负荷存成结构体数组单独写一个脚本load_case_data.m分布式电源参数和预测曲线存成.mat文件用load命令读入优化模型用 YALMIP 建模求解器用 Gurobi如果只有Mosek也可以但Mosek对二阶锥支持不如Gurobi顺手数据准备阶段比较容易忽略的是电价的时段划分。很多论文里的分时电价用的是峰、平、谷三段但实际仿真时如果用的是真实电网数据电价可能是半小时一个点甚至15分钟一个点。这会给日前调度带来一个问题模型的时间尺度到底取多少常见做法是日前阶段取1小时分辨率日内阶段取15分钟分辨率两阶段之间通过滑动窗口对接。如果模型只是演示性质统一用1小时分辨率也可以但要在论文里说明。4.2 YALMIP建模核心代码以第一阶段为例YALMIP建模关键代码如下% 定义变量 P_g sdpvar(n_g, T, full); % 常规机组出力 u_g binvar(n_g, T, full); % 机组启停状态 P_ch sdpvar(n_ess, T, full); % 储能充电功率 P_dis sdpvar(n_ess, T, full); % 储能放电功率 SOC sdpvar(n_ess, T, full); % 储能的荷电状态 P_grid sdpvar(1, T, full); % 联络线交换功率 % 目标函数 objective 0; for t 1:T objective objective sum(a_g * P_g(:,t).^2 b_g * P_g(:,t) c_g * u_g(:,t)); objective objective grid_price(t) * P_grid(t); objective objective sum(C_ess * (P_ch(:,t) P_dis(:,t))); end % 约束条件集合 constraints []; % 机组出力上下限 constraints [constraints, P_min .* u_g P_g P_max .* u_g]; % 储能SOC递推 for t 2:T constraints [constraints, SOC(:,t) SOC(:,t-1) eta_ch * P_ch(:,t) - P_dis(:,t)/eta_dis]; end % 功率平衡 for t 1:T constraints [constraints, sum(P_g(:,t)) sum(P_dis(:,t)) - sum(P_ch(:,t)) P_grid(t) sum(P_wt(:,t)) sum(P_pv(:,t)) total_load(t)]; end % 求解 ops sdpsettings(solver, gurobi, verbose, 2); optimize(constraints, objective, ops);有几个细节要提醒新手sum(P_wt(:,t)) sum(P_pv(:,t))用的是预测值不是决策变量。如果要做两阶段这里应该把分布式电源出力也设为决策变量并受到预测值上限约束。储能SOC公式里eta_ch和eta_dis分开设置因为充放电效率不同。有些简化模型直接用一个效率但实际电池在充和放时损耗确实不一样尤其在SOC较高时充电效率下降明显。目标函数里的P_g(:,t).^2是非线性项如果Gurobi处理不了二次目标MIQP需要先把目标线性化。一个常见做法是用分段线性近似替代二次成本函数YALMIP里可以用pwf函数处理不过直接用Gurobi的MIQP求解器通常也没问题。4.3 第二阶段的模型与衔接处理第二阶段的基本思路是把第一阶段确定的机组状态 (u_g^*)、储能基准功率 (P_{ess}^{base})、分布式电源基准出力 (P_{dg}^{base}) 当作已知参数然后在新场景下求解再调度模型。% 第二阶段场景s下做再调度 % 固定第一阶段决策中的整数变量 constraints_2nd [constraints_2nd, u_g u_g_star]; % 整数变量保持不变 % 再调度调整量 delta_Pg_up sdpvar(n_g, T, full); delta_Pg_down sdpvar(n_g, T, full); delta_Pdg sdpvar(n_dg, T, full); % 分布式电源削减量 % 功率平衡场景s下的实际值 for t 1:T constraints_2nd [constraints_2nd, ... sum(P_g_star delta_Pg_up(:,t) - delta_Pg_down(:,t)) ... sum(P_dis_star(:,t) - P_ch_star(:,t)) ... % 这里简化处理储能 P_grid_2nd(t) sum(P_wt_scen(:,t)) - sum(delta_Pdg_wt(:,t)) ... sum(P_pv_scen(:,t)) - sum(delta_Pdg_pv(:,t)) total_load(t)]; end这里的关键是“固定整数变量”。在Matlab中直接把u_g替换成u_g_star常数值即可不要让求解器再优化启停。否则第二阶段变成一个完整的MIP问题那就失去了两阶段“快速再调度”的意义。第二阶段求解完成后把各场景下的调整成本加权求和加上第一阶段的基准成本就是这个日前调度方案的期望总成本。至此两阶段模型的计算闭环完成。4.4 求解器选型与参数调优Matlab里求解这类问题可选的求解器不少我的经验是Gurobi对于MISOCP和MIQP速度和稳定性都没得说学术免费强烈推荐。Mosek对锥规划的数值稳定性比Gurobi略好但整数变量的MIP能力稍弱。CPLEX老牌求解器目前新版本对配电网优化支持也不错但在学生群体里用得少了。Matlab内置的求解器如intlinprog只能解MILP处理不了SOCP。如果只是做线性化潮流可以用一旦涉及非线性潮流约束就束手无策了。调参方面最常碰到的坑是sdpsettings(solver, gurobi)之后YALMIP有时会把变量自动转换为其他格式导致Gurobi报错 “Model is infeasible”。这种情况下我一般是先检查约束里是不是存在矛盾条件比如某个节点既要求电压不超过1.0又要求注入功率满足某个下限两者可能互斥。把约束逐个注释掉跑一次模型找到引起不可行的那一组约束通常很快就能定位问题。5. 常见问题与排查技巧实录5.1 模型求解速度慢怎么优化两阶段模型如果直接用全网所有节点建立约束负荷节点多了后约束数量会爆炸式增长。比如一个33节点系统24小时每个时段的潮流约束就有几十条加上机组、储能、配网约束变量总数轻松上万。求解MISOCP可能要几分钟甚至更久。优化手段用场景削减缩小第二阶段场景数量10个场景以内是比较合理的。电压约束用软约束只要在目标函数中惩罚越限就可以改用单纯形快速求解。如果只是研究调度策略而不关注潮流细节可以把配电网等效为一个大节点只保留馈线出口的功率平衡约束求解速度能提升几个数量级。我这里有一个实际案例某个33节点的配电网模型初始建模用了完整DistFlow加SOCP松弛加上8个场景Gurobi求解耗时约200秒。后来把分布式电源接入点压缩到3个关键节点潮流网络做了等效简化求解时间降到8秒而调度成本只多了不到3%。对前期方案比选来说这个效率提升非常值得。5.2 出现 “Infeasible problem” 怎么排查这个问题我在教学中被问得最多。排查思路很重要一定要按顺序来第一检查功率平衡约束。很多新手在写功率平衡时忘了把线路损耗、储能自损耗算进去导致每个时段的总发电大于总负荷或小于无解。可以先用线性潮流把网损忽略看看模型是否可解如果可解再逐步加入网损。第二检查储能SOC递推公式。SOC的初值如果不合理比如要求首末SOC相等但电池容量又不够大就会导致无解。我一般设置首末SOC相同并给出一个较小的充放电功率上限避免约束过紧。第三检查爬坡约束。机组爬坡能力设置过小而负荷波动又大可能会导致特定时段无法满足功率平衡。把爬坡约束系数放大10倍如果模型变可解了说明问题出在这儿。第四检查电压约束。分布式电源接入容量越大电压越限越容易发生。如果约束是 0.95 ≤ V ≤ 1.05可以先放宽到 0.9~1.1 试试如果能解再逐渐收紧。5.3 分布式电源渗透率过高时调度结果不合理有时候模型求解正常但结果里会看到在分布式电源出力高、负荷轻的时段联络线功率为0甚至为负数向上一级电网倒送但常规机组仍然以最小出力运行。这看起来反直觉实际上是目标函数里没有给常规机组设置启停惩罚导致的。解决方案是增加机组的启停成本或者给常规机组设置最小运行时间约束。否则优化器会在某个时段让机组停机以降低成本但下一时段又需要它启动了启停太频繁工程上完全不可接受。这个问题在IEEE 33节点系统里非常典型我建议直接加入最小启停时间约束即使求解时间增加一点也值得。5.4 常见错误速查表下面这张表是我做Matlab配电网调度模型时总结出来的高频问题分享给大家现象可能原因处理方法求解器提示“Infeasible”功率平衡约束不平衡检查发电/负荷/网损是否闭合储能SOC越界递推公式效率参数方向错误确认充电/放电效率是否对应正确方向优化结果里机组频繁启停缺少启停成本或最小运行时间约束在目标中加入启停惩罚项电压越限但优化器不处理电压约束写成了软约束但惩罚系数太小增大惩罚权重或改为硬约束第二阶段无解第一阶段基准方案太激进放宽部分约束或加入松弛变量和惩罚项求解时间过长场景数过多或SOCP约束过多削减场景简化潮流为线性近似MATLAB报错“Undefined function”工作路径没有包含函数文件检查当前文件夹及路径设置5.5 代码调试的独家技巧最后分享几个实战小技巧这些是常规教程里不会细讲的第一个技巧是写一个check_constraints.m函数在优化求解前把所有约束的残差打印出来这样能非常直观地判断是哪个约束导致无解。我见过不少人用YALMIP时只盯着优化器返回状态结果状态码又不明确白白浪费了好几个小时。其实只需在optimize之后加一行check(constraints)YALMIP 就会逐一列出所有约束的 OK 状态和残差排查效率瞬间翻倍。第二个技巧是在模型测试阶段先固定所有整数变量只用测试数据跑连续优化验证连续模型的可行性再去处理整数变量的MIP问题。这样可以快速区分是数学模型本身的问题还是整数变量的组合爆炸问题。第三个技巧是注意YALMIP和Gurobi版本兼容性。我踩过一次坑YALMIP版本太老传给Gurobi的二阶锥约束在求解器里被识别成了非凸二次约束结果求解时间增加了10倍换了新版本YALMIP后立刻恢复。如果你的模型突然变慢优先怀疑版本匹配问题。从我个人体会来说做这类含分布式电源的配电网两阶段优化调度模型真正难的不是算法原理而是把原理落地的过程中那些零零碎碎的细节——场景怎么生成、约束怎么松弛、求解器怎么调、无解时怎么定位。这些经验不是看几篇论文就能获得的需要在反复的调试中慢慢积累。希望这篇文章能帮你少走一些弯路把时间花在真正有价值的问题上。