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 = i+1: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 = i+1: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_{t=1}^{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] ]
其中 (T=24),第一项是常规机组的燃料成本(通常简化为二次函数),第二项是从上级电网购电的成本((\lambda_t^{grid}) 是分时电价),第三项是储能充放电损耗成本。这个目标函数没有考虑弃风弃光惩罚,如果项目要求优先消纳分布式电源,需要在目标里增加一项:
[ C_{curtail} = \sum_{t=1}^{T} \sum_{d \in DG} c_{curtail} \cdot (P_{d,t}^{forecast} - P_{d,t}^{schedule}) ]
惩罚系数 (c_{curtail}) 应大于购电电价,否则优化器会为了省钱而主动弃风弃光,违背了分布式电源优先消纳的初衷,这是我做项目时踩过的坑,调参时需要格外注意。
第一阶段的约束包括:
- 潮流约束(DistFlow 线性化方程或交流潮流方程)
- 节点电压上下限
- 支路电流/功率上限
- 常规机组出力上下限与爬坡约束
- 储能SOC约束与充放电功率约束
- 联络线交换功率约束
3.2 第二阶段:日内再调度子问题
第二阶段的思路是:第一阶段的日前计划已经确定,但当实际场景(或更精确的预测场景)出来后,某些约束可能不满足了,比如光伏实际出力比预测低,导致部分负荷需要切掉或电压越限。此时需要在“最小化调整量”的目标下重新分配出力:
[ \min ; \sum_{t=1}^{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 ]
这样就把非线性潮流转化成了一个凸约束,整个模型变成了混合整数二阶锥规划(MISOCP),YALMIP配合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后立刻恢复。如果你的模型突然变慢,优先怀疑版本匹配问题。
从我个人体会来说,做这类含分布式电源的配电网两阶段优化调度模型,真正难的不是算法原理,而是把原理落地的过程中那些零零碎碎的细节——场景怎么生成、约束怎么松弛、求解器怎么调、无解时怎么定位。这些经验不是看几篇论文就能获得的,需要在反复的调试中慢慢积累。希望这篇文章能帮你少走一些弯路,把时间花在真正有价值的问题上。