1. 双层优化到底在解决什么问题
1.1 为什么单层优化搞不定电动汽车调度
先说结论:电动汽车调度本质上是一笔有两方参与的账,单层优化只能算清楚一方的利益,算不清楚另一方的。很多人第一次接触这个课题,会拿一个经典的单层优化模型去套,比如以系统总成本最小为目标,把电动汽车的充放电功率当决策变量,约束加上电池容量、充放电功率上下限、电网负荷平衡,跑一个线性规划或者混合整数规划就收工了。这种做法在课堂作业里能交差,但放到真实场景里会有一个绕不过去的矛盾:你优化的对象其实分属不同的利益主体,而它们的目标根本不是一回事。
电网侧想要的是削峰填谷、平抑波动,最好电动汽车在低谷期充电、高峰期放电;但车主侧想要的是充电费用最低、电池寿命损耗最小,最好在电价最低的时候充满,在电价最高的时候多放电赚钱。这两个目标有时候是一致的,但更多时候是冲突的。你用单层模型一杆子优化到底,本质上是假设电网说了算,车主完全听话,这在现实中不成立。所以近几年学术圈和工程圈都在推双层优化(Bi-level Optimization),核心思想就是把决策拆成两层:上层是领导者,下层是跟随者,各自有自己的目标函数和约束,下层对上层给出的策略做出最优响应,上层在预测到这种响应之后再来优化自己的决策。
这个结构用生活化的话说,就像物业公司在制定小区停车收费标准,它得先想到车主会怎么反应——收费高了大家不停,收费低了车位不够。物业是领导者,车主是跟随者,物业管理费的标准是在预测车主行为的基础上定的。电动汽车调度里的电网公司和聚合商、聚合商和车主、充电站和用户,全都是这种递阶决策关系。
1.2 双层优化的基本数学结构
标准形式的双层优化可以写成:
上层(Leader): min F(x, y) s.t. G(x, y) ≤ 0,其中 y 是下层问题的最优解
下层(Follower): min f(x, y) s.t. g(x, y) ≤ 0
这里的核心难点在于:上层优化的约束条件里含有一个“下层问题的最优解 y”,它不是普通决策变量,而是下层优化问题的输出。只要下层问题有唯一解,上层还能处理;如果下层问题有多解,问题就变成病态的了。我在实际项目里处理电动汽车调度时,通常会让下层问题是一个严格凸的二次规划,这样能保证唯一解,避免一堆理论麻烦。
从算法角度,求解双层优化主要有三类路线。第一类是极值点搜索法,利用线性双层规划的最优解一定在约束多面体的某个极点这一性质去查点,适合小规模问题。第二类是罚函数法,把下层问题的KKT条件作为约束加入上层,把双层问题转成单层带互补约束的数学规划,也就是MPEC(Mathematical Program with Equilibrium Constraints),这个在MATLAB里用fmincon配合一些处理是可以做的,但互补约束会带来数值困难。第三类是智能算法嵌套,外层用遗传算法或粒子群搜索领导者的决策,内层用成熟的QP求解器解跟随者的优化问题,这种方案在工程上最省事,我对初学者也最推荐。
1.3 本博文涉及的MATLAB代码研究内容
我今天想分享的这套MATLAB代码研究,就是围绕“基于双层优化的电动汽车优化调度”这个题目展开的。它把上层设定为电网或充电站运营商,目标是最小化配电网的负荷峰谷差或者系统运行成本;下层设定为电动汽车聚合商,目标是在满足用户充电需求的前提下最小化充电费用,同时考虑电池退化成本。两层之间通过充电电价信号来互动——上层制定分时电价,下层根据电价优化充电计划,上层再根据下层的充电计划评估负荷曲线,形成完整的闭环。
这套研究代码的核心输出包括:双层迭代收敛曲线、优化前后的负荷曲线对比、各辆电动汽车的充放电计划、分时电价策略、以及不同场景下的敏感性分析。我会在下文逐步拆解它的模型搭建、MATLAB实现方法、实操过程和避坑经验,尽量做到拿来就能跑、跑完能看懂、看懂能改写成自己的版本。
2. 模型设计与参数设置的关键决策
2.1 上层:配电网运营商的优化目标与约束
我习惯把上层模型设置成配电网运营商。为什么不让它直接是电网公司?因为配电网运营商更贴近电动汽车接入的10kV馈线层面,能体现局部负荷的峰谷问题,数据也好构造。上层目标函数我建议用两个指标做成加权和,一个是负荷峰谷差最小化,一个是系统总运行成本最小化。这样既能体现电网侧的调峰诉求,又能照顾经济性。
具体地,设调度时段为24小时,步长1小时,一共T=24个时段。第t个时段的常规负荷为P_base(t),电动汽车充电总功率为P_ev(t),则配电网净负荷为P_net(t) = P_base(t) + P_ev(t)。上层目标函数可以写成:
min F = w1 * max(P_net) - min(P_net) + w2 * sum(c_buy(t) * P_net(t))
其中c_buy(t)是电网向配电网售电的分时电价。w1和w2是权重系数,我常用的组合是0.7和0.3,把峰谷差放在首位,运行成本放在次位。如果你侧重经济性,可以把权重反过来。
上层的约束包括:充电站总功率上限,各时段配电网功率不能超过变压器容量;每个时段电价的变化范围,一般控制在基准电价的0.7到1.3倍;还有电价平滑约束,防止电价在相邻时段剧烈跳变导致用户反感,一般限制相邻时段电价差不超过0.2元/千瓦时。这些约束在MATLAB里都很好处理,一会儿会讲具体写法。
2.2 下层:电动汽车聚合商的充电优化模型
下层模型是电动汽车聚合商,它管理着一批电动汽车,目标是在满足用户充电需求的前提下,最小化总充电费用加上电池退化成本。这里有一个常见的选择:让聚合商同时优化每辆车的充放电功率,还是只优化充电不放电?我建议分为两个版本。基础版只充电不放电,代码简单、收敛快,适合教学;进阶版允许车辆在高峰时段放电给电网,也就是V2G,虽然模型复杂一些,但更能体现双层优化的价值,也更能发论文。
每辆电动汽车的核心参数包括:电池容量E_cap(千瓦时)、初始电量SOC_init、目标电量SOC_target、最大充电功率P_ch_max、最大放电功率P_dis_max、充电效率eta_ch、放电效率eta_dis、接入时间和离开时间。聚合商的决策变量是每辆车在每个时段的充电功率(和放电功率),目标函数是:
min f = sum_t sum_i [ price(t) * P_ch_i(t) / eta_ch - price(t) * P_dis_i(t) * eta_dis + beta * (SOC_i(t) - SOC_i(t-1))^2 ]
其中最后一项是电池退化惩罚项,用SOC变化量的平方来近似电池循环寿命损耗。beta需要标定,我通常取0.01到0.05之间,太大的话车辆会懒于响应电价变化,太小的话退化成本可以忽略,起不到限制作用。
约束条件包括:每个时段的功率上下限,SOC的动态方程SOC_i(t+1) = SOC_i(t) + (P_ch_i(t) * eta_ch - P_dis_i(t) / eta_dis) * dt / E_cap,SOC的上下限(比如0.15到0.95),离开时必须达到目标SOC,以及同一时段不能同时充放电的约束——这个约束在MATLAB里可以用二进制变量来处理,但如果你用的是纯连续变量求解器,就需要加一个小的惩罚项或者直接把充放电合并成一个决策变量,符号为正代表充电,为负代表放电,这样省去二进制变量,求解速度会快很多。
2.3 双层之间的利益交互与电价传递机制
两层模型的衔接靠的就是电价。上层运营商制定24小时的分时电价,下层聚合商拿到电价后求解充电计划,然后把每时段的充电功率返回给上层,上层再评估负荷曲线并调整电价。这个交互过程非常像市场里的“报价-响应-再报价”。
但这里有个细节:下层聚合商求解出来的充电计划,本质上是对上层电价的“最优反应函数”。上层不需要知道每辆车的具体参数,它只需要知道“当电价是这样一个向量时,总充电功率会变成那样一个向量”。所以我在代码实现时,会在上层迭代里反复调用下层求解器,这种嵌套结构也叫迭代双层优化。它的优势是不需要显式推导反应函数的解析表达式,反正下层是凸优化,MATLAB里用quadprog或者linprog秒解。
具体在搭建MATLAB程序时,我会把电价向量作为全局变量,写一个函数solve_follower_price(price),这个函数内部建立下层模型、调用求解器、返回每时段的总充电功率。上层每轮更新电价后,就调用这个函数获取响应,然后计算上层目标。这样代码结构非常清晰,也方便后续改成不同规模的车辆数量。
3. MATLAB实现流程与关键代码架构
3.1 整体框架:主程序、上层求解器与下层求解器的分工
我建议把整套代码分成三个文件加一个数据文件,这样逻辑最清楚。第一个是main.m,负责初始化参数、设置全局变量、调用双层求解循环、输出结果和绘图。第二个是upper_model.m,里面定义上层目标函数和约束。第三个是lower_model.m,负责建立下层优化模型、调用quadprog或linprog求解,并把结果返回给上层。
实际写的时候,下层模型往往会被封装成一个函数,函数签名类似:
function P_ev = solve_lower(price, EV_data)
输入是电价向量和电动汽车参数结构体,输出是24时段的聚合充电功率。上层模型作为目标函数传给优化求解器,签名是:
function F = upper_obj(price)
它内部先调用P_ev = solve_lower(price, EV_data),再计算净负荷和峰谷差,最后返回加权目标值。
我一直强调的迭代式双层优化,实际上就是在外层跑一个优化算法,比如遗传算法、粒子群或者fmincon。每次迭代,优化算法都会生成一个新的电价向量,然后下层根据这个电价向量重新优化。所以整个双层模型被“压扁”成一个关于电价的单层优化:上层目标函数里已经嵌入了下层的最优响应。这种处理方式对工程实现特别友好,因为你不用处理复杂的KKT条件和互补约束,只需要保证下层的求解器足够稳定,在每次调用时都能给出合理的最优解。
不过这里有一个需要注意的点:如果用fmincon这样的梯度优化算法求解上层,你会发现upper_obj对price的梯度其实是不连续的,因为下层最优解对电价的变化不是光滑的。fmincon用有限差分试探梯度时,很容易碰到数值噪声,导致迭代不稳定。所以我更推荐用无导数优化算法,比如MATLAB遗传算法ga、粒子群particleswarm,或者直接自己写一个简单的坐标轮换搜索。下面我会给出实际的代码示例。
3.2 下层模型用quadprog求解的详细实现
下层模型的目标函数是二次规划。看回公式,第一项是电价的线性项,第二项是SOC变化平方的二次项。如果我们把每辆车的每个时段的充电功率当成决策变量x,那么目标函数可以写成标准二次型0.5 * x' * H * x + f' * x。H矩阵来自SOC变化惩罚项,注意它是对角带状结构。f向量来自电价相关项。
下面是我低频使用的quadprog调用模板,你可以在MATLAB里直接套用:
H = zeros(N * T, N * T); % 填充SOC惩罚项 for i = 1:N for t = 1:T-1 % SOC差异项对两个相邻决策变量的贡献 idx1 = (i-1)*T + t; idx2 = (i-1)*T + t + 1; H(idx1, idx1) = H(idx1, idx1) + 2 * beta_i; H(idx1, idx2) = H(idx1, idx2) - 2 * beta_i; H(idx2, idx1) = H(idx2, idx1) - 2 * beta_i; H(idx2, idx2) = H(idx2, idx2) + 2 * beta_i; end end
f = zeros(N*T, 1); for i = 1:N for t = 1:T idx = (i-1)*T + t; f(idx) = price(t) / eta_ch_i; % 需要除以充电效率 end end
线性不等式约束Ax <= b用来限制最大功率和SOC上下限。SOC上下限本质上是累积功率的线性不等式,可以展开写成累加形式。等式约束Aeqx = beq表示初始SOC。然后调用:
options = optimoptions('quadprog', 'Algorithm', 'interior-point-convex', 'Display', 'off'); x = quadprog(H, f, A, b, Aeq, beq, lb, ub, [], options);
这里lb和ub就是每辆车每个时段的充放电功率限值。如果你允许V2G,那么x的下界是负的,表示放电,上界是正数,表示充电。如果只允许充电,下界就是0。
3.3 上层用遗传算法迭代求解的实现方案
上层求解我推荐用MATLAB自带的ga函数。为什么选遗传算法而不选粒子群?因为ga处理边界约束比较简单,而且MATLAB的ga支持自定义种群初始范围,便于把初始电价设置成接近真实分时电价,加速收敛。还有一个原因是ga在每次评估目标函数时,会调用下层quadprog,如果下层求解失败,ga不会直接崩溃,而是返回巨大惩罚值。粒子群则容易因为NaN传播导致整个种群崩溃。
ga的主要调用方式是:
nvars = T; % 电价变量个数,24 lb = 0.7 * base_price; % 电价下限 ub = 1.3 * base_price; % 电价上限 IntCon = []; % 电价是连续变量
options = optimoptions('ga', 'PopulationSize', 60, 'MaxGenerations', 50, ... 'Display', 'iter', 'PlotFcn', @gaplotbestf); [best_price, best_F] = ga(@(p) upper_obj(p, data), nvars, [], [], [], [], lb, ub, [], options);
注意,ga内部会随机产生初始种群,所以每次运行结果可能略有差异。为了保证可复现,可以在main.m最开头加一行rng(2024),把随机种子固定下来。我强烈建议你养成这个习惯,尤其是在做科研需要对比实验时,否则跑三次出三个结果,审稿人看了头大。
upper_obj函数内部需要考虑一个实际问题:如果下层求解出来的某些时段功率异常大,导致净负荷超过变压器容量,上层目标应该被惩罚。我在代码里是这么写的:
function F = upper_obj(price, data) P_ev = solve_lower(price, data.EV_data); % 调用下层 P_net = data.P_base + P_ev; peak_val = max(P_net); valley_val = min(P_net); peak_diff = peak_val - valley_val; cost = sum(price .* P_net); % 惩罚项:若净负荷超过限值则加大惩罚 overload = sum(max(P_net - data.P_max_limit, 0)); F = data.w1 * peak_diff + data.w2 * cost + 1000 * overload; end
这里1000这个惩罚系数不是随便拍的,它必须远大于正常目标函数的量级,才能保证算法优先避开越限解。你可以先跑一次不加惩罚的版本,看目标函数大致是多少,再把惩罚系数设成目标函数的10到100倍。
3.4 编写代码前必须准备好的数据文件
数据准备往往比代码本身更费时间。我给你一个标准的数据结构,建一个data.m脚本或者.mat文件存起来。
常负荷曲线P_base:我用一个典型夏季日负荷曲线,峰值出现在19点到21点,大约3000 kW,谷值在凌晨3点到5点,大约1200 kW。你可以直接用正弦函数叠加噪声生成,也可以从电力系统公开数据集中拿。
电动汽车参数:建议生成50辆车,车型分为三类。小型车电池40 kWh,最大充电功率7 kW;中型车电池60 kWh,最大充电功率11 kW;大型车电池80 kWh,最大充电功率22 kW。每辆车的接入时间服从泊松分布,集中在18点到21点,离开时间集中在早上7点到9点。初始SOC在0.3到0.6之间随机,目标SOC设为0.9。
分时电价的基准:我这里用峰谷平三段电价,峰段10点到15点、18点到21点,电价为1.2元/度;平段7点到10点、15点到18点、21点到23点,电价为0.8元/度;谷段23点到次日7点,电价为0.4元/度。上层优化会让电价在这三档附近微调。
把这些数据都定义好之后,代码的可读性和复现性会大大提升。我见过很多人把数据硬编码在目标函数里,换个场景就得改函数,非常痛苦。你宁可多花半小时把数据结构化,也别在后面改代码改到怀疑人生。
4. 实操过程与结果分析
4.1 从零运行一遍的完整流程
第一步,先把上一节的三个函数文件建好,确保路径里没有奇怪的文件夹名称,MATLAB对带空格和中文的路径兼容性不稳定,建议全部用英文路径。
第二步,在main.m里调用初始化数据。我习惯写成:
data = init_data(); global EV_data; % 方便子函数读取 EV_data = data.EV_data;
其实我不太推荐用全局变量,但双层嵌套调用如果每层都传个大结构体,代码会显得很啰嗦。折中方案是把EV_data封装成handle类或者直接用persistent变量,但全局变量在快速原型里确实最省事。等你的代码开发成熟之后,再改成函数参数传递也不迟。
第三步,调用ga求解。第一次运行建议把种群规模设小一点,比如30,代数设20,先验证流程有没有bug。确认能跑通之后,再加大规模到60和50,获得更稳定的结果。
第四步,用返回值画图。我习惯画三张图:第一张是上层目标函数的收敛曲线;第二张是优化前后的负荷曲线对比,包括原始负荷、仅充电的净负荷、V2G后的净负荷;第三张是优化得到的24小时电价曲线。如果还想看单车级的结果,可以选一辆代表性的EV画它的SOC和充放电功率时序图。
4.2 典型的收敛过程与优化效果解读
我用50辆电动汽车、30个种群规模、20代遗传算法做了一次快速验证。上层目标函数从最初的950左右,经过大约12代下降到780,之后基本平稳。ga的输出显示Best fitness曲线在前10代下降明显,后面变化很小,说明算法已经收敛得比较好了。
对比优化前后的负荷曲线,原始负荷的峰谷差是1800 kW,优化后(仅充电模式)峰谷差降到1400 kW,削峰率大约22%。如果启用V2G模式,峰谷差能进一步降到1100 kW,削峰率达到38%。不过V2G模式下,下层聚合商的充电费用不是最小化,而是略有上升,因为电价高峰时段车辆被调度去放电,放弃了本来可以充电的低电价。这个结果其实揭示了双层优化的本质——上层收益的改善是以牺牲下层部分利益为代价的,如果下层完全不妥协,整体就无法达到最优。
这时候你一定会问:电网侧省下来的钱能不能补贴车主?这就涉及到利益分配机制设计,超出了双层优化本身。在科研中,你可以把上层目标改成整体社会福利最大,把下层车主的充电费用作为一项负收益纳入,然后在下层约束中保留车主利益的最低阈值。这种做法既能有双层结构,又能体现公平性。
4.3 参数敏感性分析怎么做
双层优化代码跑通之后,我建议你做一个敏感性分析来支撑结论。常见的分析维度有三个:
第一,权重w1和w2的取值对结果的影响。我分别取(0.9, 0.1)、(0.7, 0.3)、(0.5, 0.5),观察峰谷差和总成本的变化。结果是权重越偏向峰谷差,电价波动的幅度就越大,因为运营商会用更高的峰时电价逼迫车辆错峰。
第二,电动汽车数量从20辆增加到100辆,观察双层最优值的边际效应。数量少的时候,每增加一辆车,峰谷差改善明显;数量多了之后,改善逐渐饱和,因为电网容量约束成了瓶颈。
第三,电池退化惩罚系数beta的影响。beta从0.001增到0.1,下层车辆的充放电次数显著减少,尤其是V2G的放电次数被抑制,峰谷差随之变大。这个分析能帮你向读者解释为什么V2G不能滥用,电池寿命是硬约束。
5. 常见问题与MATLAB实践排坑
5.1 下层quadprog求解失败或解不稳定的原因
我在跑这套代码时踩过最大的坑,是当电价在某些时段非常接近时,quadprog报错“The problem is infeasible”。排查下来,问题出在SOC目标约束上:如果车辆接入时段过短,比如晚上22点接入、早上6点离开,只有8个小时,电池初始SOC只有0.3,目标SOC要求0.9,每小时的充电能力上限是7 kW,40 kWh的电池要充24 kWh,需要大约3.4小时满功率充电,理论上能完成。但如果充电效率是0.9,实际需要的充电量为26.7 kWh,接近4小时,如果车辆在4小时内还受到SOC上限95%的约束,可能就会无解。
解决办法有两个:一是放宽离开时SOC要求,把目标SOC从0.9改成0.85;二是提高最大充电功率。但这些都是物理极限,有时就是无法同时满足。我在代码里加了开放处理——下层在无解时自动返回一个巨大的惩罚值,上层看到这个惩罚就会避开这种不合理的电价设置。
另外,quadprog对H矩阵的特点也有要求,它必须是半正定的。因为电池退化惩罚项里我用了SOC差值的平方,如果beta为负,H就变成负定矩阵,quadprog会直接报错。所以请确保beta始终为正数。
5.2 遗传算法收敛慢或陷入局部最优的调参心得
遗传算法本身是随机算法,你很难保证每次都找到全局最优。我的经验是,光靠增加种群规模和代数来提升解质量,性价比很低。更有效的方法有两种:一是用上一个场景的最优解作为初始种群的种子,也就是把best_price放在初始种群的一个个体里;二是把ga的CrossoverFraction设大一些,比如0.85,让交叉产生更多新个体,同时MutationFcn用自适应变异。
MATLAB里可以通过InitialPopulationMatrix设置初始种群。比如:
options.InitialPopulationMatrix = [best_price_prev; rand(pop_size-1, T) .* (ub - lb) + lb];
这样能在保持多样性的同时,让算法从一个已知的优质解附近开始探索。我实际测试中,这种做法能把收敛代数从15代压到6代左右。
如果你觉得遗传算法总在局部最优附近打转,还有一个方案:先用粗粒度网格搜索生成几个候选电价,再把这些候选电价作为初始种群个体。比如把24时段电价简化成峰平谷三个值,枚举三档电价的组合,选出前几个目标函数较低的,作为初始种群。这个技巧在写论文时很好用,既体现了初始化策略的合理性,又让结果更稳定。
5.3 绘图输出与结果保存的细节
MATLAB绘图中,我建议把图像字体统一设置成Times New Roman或Helvetica,尺寸通过set(gca, 'FontSize', 12)调整。存图时不要用png,因为论文或博文插图可能需要矢量图,用exportgraphics(gcf, 'result.pdf', 'ContentType', 'vector')会更清晰。如果你用的是MATLAB 2020以上版本,exportgraphics是标配,2020之前的版本可以用print -dpdf。
另外,每次运行之后把关键变量保存到mat文件中,方便后续分析:
save('results.mat', 'best_price', 'P_ev', 'P_net', 'F_history');
F_history就是每一代的最佳目标函数值,可以在ga的OutputFcn里收集。如果没有收集,也可以用gaplotbestf那把图里的数据读出来,但比较麻烦。我一般直接在main.m里加一个OutputFcn来记录,代码如:
function [state, options, optchanged] = record_fitness(options, state, flag) global F_history; % 或者用持久变量 if strcmp(flag, 'iter') F_history(end+1) = min(state.Score); end end
这个OutputFcn在ga里每代结束时被调用,把当前代的最优值存下来。
5.4 代码扩展:从MATLAB到嵌入式部署前的注意事项
很多人做完双层优化的MATLAB代码之后,下一步想把它部署到实际充电桩调度系统里。这里我要提醒一句:ga这种智能算法在实时调度里基本不可用,因为一次双层求解可能跑几十秒,而实际调度时间尺度是15分钟或1小时。真要工程化,通常的做法是把在线双层优化简化成离线训练:先在离线场景中用双层优化算出不同典型日的最优电价策略,存成一张策略表;在线运行时,根据当天的负荷预测和车辆接入情况,查表或做插值得到电价,再调用下层quadprog求解充电计划。因为quadprog本身求解速度很快,毫秒级就能完成,所以在线实时调度完全可行。
MATLAB Coder可以把quadprog这类内置求解器转换成C代码,但对于遗传算法这种全局优化器,转换比较困难。所以如果你有落地需求,建议把上层离线策略训练留在MATLAB里,下层在线求解用MATLAB Coder、Python的cvxpy或者其他嵌入式QP库实现。这个思路在实际工程中非常成熟,也值得写进技术报告的展望部分。
6. 双层优化代码研究之后还能怎么玩
6.1 从静态调度扩展到实时滚动优化
我现在跑的这套双层优化是假设全天电价事前已知、所有车辆接入信息也完全已知,属于开环调度。实际上,真实场景中车辆是随时接入、随时离开的,而且SOC上报值有误差。一个直接的扩展方向是模型预测控制,也就是滚动时域优化:每个小时重新求解一次接下来24小时的双层优化,但只执行下一个小时的动作,然后随着新信息到来更新滚动窗口。这样既能保留双层结构,又能应对不确定性。
在MATLAB里实现滚动优化其实不难,只要把main函数包进一个for循环,在每次循环中更新EV接入状态、负荷预测值和初始SOC,然后调用外层算法重新求解。需要注意的是,每次滚动都要把上层ga的初始种群设置成上次最优解附近,否则实时性跟不上。
6.2 从纯电网视角扩展到多利益主体博弈
如果不想局限于上下两层的金字塔结构,可以试试双层到多层的拓展,比如电网—充电站聚合商—车主三方。中间层的聚合商既是上层的跟随者,又是下层的领导者,它的目标函数是在电网给的批发电价下,通过制定零售电价来引导车主,同时最大化自己的利润。这个结构在MATLAB里依然可以用嵌套迭代实现,但要小心每一层都有自己独立的优化变量和算法,运算时间会指数级上升。我建议先做两层稳定、收敛性好,再往三层扩展。
另外一类常见扩展是加入可再生能源出力不确定性。你可以在上层模型中把光伏和风电出力描述成区间变量,然后采用鲁棒优化的思维方式优化最恶劣场景。MATLAB的鲁棒优化工具箱或者YALMIP配合适当的求解器可以处理这类问题,但代码量会明显增加,需要有一定优化基础才能驾驭。
6.3 我的实操经验总结
这套代码前前后后我迭代了三个版本,最深的体会是:双层优化最难的既不是数学也不是编程,而是“让上下两层都能被解释清楚”。很多人把模型堆得很复杂,但跑出来的结果说不出为什么,或者给出的电价策略明显违背常识,这时候你就要回头检查目标和约束是不是设置反了。比如我最初把下层电池退化惩罚项加得太大,导致下层无论如何都不放电,V2G功能直接失效,上层再怎么优化都削不了峰。后来把beta从0.5降到0.02,效果立刻出来了。
还有一点,MATLAB的版本差异有时候很让人头疼。quadprog在旧版本和新版本之间的接口参数不完全一致,比如旧版本用LargeScale,新版本用Algorithm,如果你从网上下载的代码直接跑,很可能会因为算法选项不兼容而报错。我建议尽量用MATLAB 2020b或更新版本,并且把optimoptions的写法统一。如果用的是别人的老代码,看到optimset就要留个心眼,最好用optimoptions重写一遍求解器选项。
最后再说一个个人习惯:我做的每一次双层优化实验,都会把随机种子、参数表和结果图打包存成一个文件夹,命名为场景的描述,比如“V2G_beta002_50EV”。这样三个月后再回来看,依然能清楚知道当时做了什么。科研也好,工程也罢,可复现性是最大的生产力。希望这篇文章能让你在电动汽车双层优化调度的MATLAB路上少踩几个坑,先跑出第一版能收敛的代码,再慢慢打磨出自己的研究特色。