大规模电动汽车随机充放电优化：局部求解策略与MATLAB实现

最近帮一个园区做充电桩配套调度方案时最头疼的问题就是晚间六点到九点一两百辆车同时扎进电网。车主到达时间不确定、剩余电量不确定、第二天出发时间也不确定——这其实就是一个典型的大规模电动汽车随机充放电优化问题。起初我的想法比较天真优化目标就是充电费用最低约束就是电池SOC、充电功率上限、变压器容量扔给求解器跑就行。但等我把规模加到几千辆车之后才发现集中式全局优化的求解时间完全没法接受。后来我转向了“局部优化”的思路把车按集群拆开、把调度时段按滚动窗口切短、把随机场景用代表性样本替代。说白了就是不让一个超级大的模型一次性吃下所有变量而是通过分解、协调和反复修正来逼近全局最优。这套思路在工程里非常实用计算时间能从“按小时算”降到“按秒算”优化效果跟全局最优解的差距通常也就几个百分点。这篇文章就把建模思路、局部优化的三种切法、MATLAB里可以直接用的代码骨架以及我实测500辆和5000辆车时踩过的坑全部梳理一遍。1. 为什么“大规模随机局部优化”三个词会凑到一起1.1 先别急着堆模型这个问题的计算量是怎么涨上去的电动汽车充放电优化的基本形态并不复杂给每一辆车在每一个时段安排一个充电功率或者放电功率让它满足电池约束同时不冲击配电网。但一旦把规模做上去变量数量就非常夸张。假设有 N 辆车调度周期是24小时按15分钟一个时段就是 T96 个时段。再叠加 S 个随机场景做期望值评估优化变量规模就是 N×T×S。举个例子1000辆车、96个时段、50个随机场景那就是480万个连续变量。如果还要求充放电互斥可能还要加480万个0/1整数变量。这种规模的混合整数二次规划用Gurobi这种商业求解器都不一定能在一个小时内解出来更别提实际调度要求的是几分钟甚至几十秒内给出方案。所以当标题里出现“大规模”这个词时本质上就是在提醒你集中式全局优化大概率跑不动必须做某种形式的分解或近似。1.2 随机性才是这类问题真正难啃的地方很多人做电动汽车充放电优化时先假设所有车辆的接入时间、初始SOC、离开时间都是已知的把问题简化成确定性优化。但这种做法在真实场景里基本没用——车主的行为本身就是高度随机的。随机性主要体现在三个层面接入时间的随机性晚高峰集中在18点到21点之间但每辆车的到达时刻都有波动极端情况下前后差好几个小时。初始SOC的随机性有人到家还剩60%电有人只剩10%这个分布通常比较分散。离开时间和目标SOC的随机性有人第二天早上6点就走有人中午才走目标SOC也各不相同。如果忽略这些随机性你算出来的“最优充放电策略”很可能在实际执行中完全不成立。因为真实到达车辆和模型预设的不同功率分配自然需要推倒重来。但如果你严格把随机性全部嵌入到两阶段随机规划模型里计算复杂度又会再翻几倍。所以说随机性不是锦上添花的扩展而是这种调度问题的核心矛盾。你要么用场景法把不确定性纳入优化要么用滚动窗口边执行边修正反正不能假装所有信息都已知。1.3 局部优化的本质分解、求解、再协调局部优化的思路说起来并不玄乎既然一个全局大模型解不动那就把它切成若干个个可以快速求解的小问题每个小问题只负责一个局部区域、一段局部时间或者一个局部场景然后在上一层用协调机制把结果拼起来。电动汽车充放电问题之所以适合这样做是因为车辆与车辆之间的耦合其实很稀疏。车与车的可行域基本互不干扰它们的“交集”主要集中在一个点上同一台变压器或同一条馈线下的总功率不能超限。只要把这个全局耦合拆开剩下的单车问题几乎是一个带电池动态的线性规划或者二次规划规模很小一秒钟能算几百辆。协调的办法也很经典给每个集群施加一个“虚拟功率信号”或“虚拟电价”让子问题在优化时自动避开超容量时段。上层迭代调整这个虚拟信号直到所有子问题的加总功率都不越限。这个方法在数学上对应拉格朗日分解在工程上就像一个小型价格调节系统我后面会给出代码骨架。2. 数学底子先打好单台车、集群约束与随机变量2.1 单台车的SOC递推与行为约束单台车是整个优化问题的最小单元。先定义一个默认参数电池容量 E60kWh最大充电功率 P_max7kW常见交流慢充最大放电功率 D_max5kW充放电效率都取95%。SOC的递推关系很简单SOC(t1) SOC(t) (ηc × Pc(t) - Pd(t)/ηd) × Δt / E其中 Pc(t) 和 Pd(t) 分别是 t 时段的充、放电功率ηc 和 ηd 分别是充、放电效率Δt 是时段长度这里是0.25小时。注意第一项做除法还是乘法要看你的效率定义我在项目里统一按“充电按效率折算进电池放电按效率折算出发动机”来写免得来回改符号改出错。单台车的主要约束包括SOC上下限0.2 ≤ SOC(t) ≤ 0.9留一点安全余量也避免电池深度充放。充放电功率约束Pc(t) ≤ P_maxPd(t) ≤ D_max。充放电互斥同一时段要么充电要么放电数学上是 Pc(t)·Pd(t)0在求解时经常用0/1变量或者SOS1约束实现。离开时的目标SOC约束SOC(T_depart) ≥ SOC_target。这是硬约束代表车主第二天出行的基本需求。如果允许V2G放电还需要在目标函数里加入电池退化成本否则优化器会在所有高电价时段疯狂放电表面上看费用很低实际电池寿命损失远大于收益。2.2 集群容量约束唯一的“硬耦合”单个车辆的约束再完善也只是孤立问题。真正把所有车关联起来的是集群层面的功率上限。比如一个台区变压器容量是500kW接入100辆车每辆车都可能满功率7kW充电100辆车同时满充就是700kW变压器立刻过载。所以集群约束一般写成∑_{i∈G_k} Pc_i(t) P_base(t) ≤ P_cap(t)这里的 P_base(t) 是除电动汽车以外的常规负荷P_cap(t) 是该集群的可用容量上限G_k 是集群内的车辆集合。这个约束是跨车耦合的它决定了所有车的充电计划必须“错峰”。这个约束也是局部优化能成立的支点。因为它只有一个而且形式非常线性非常适合做拉格朗日分解把容量约束松弛掉加上一个虚拟价格 λ(t)每个子问题独立优化上层迭代调整 λ(t) 逼近真实容量约束。2.3 随机变量怎么建模更贴近现实随机变量建模如果做得太粗糙后面不管用什么算法都白搭。我在项目里的做法是先做基础调研统计再用截断正态分布或分段经验分布拟合。以接入时间为例晚高峰的到达时刻大致服从钟形分布我取均值18:30、标准差2小时然后截断在16:00到22:00之间。初始SOC我习惯用0.2到0.6之间的均匀分布加一点偏置因为大多数人不会把电跑光才回家。离开时间上午6:00到9:00标准差取1小时。生成随机场景之后常规的做法是跑蒙特卡洛抽样每个场景都是一种“可能发生的一天”。但如果场景数太多优化模型就解不动。所以后续一定要做场景缩减这个我在第三章详细说。2.4 目标函数充电费用、电池损耗与削峰填谷目标函数的设计直接影响最终策略。如果只有充电费用最小那么所有车都会挤在凌晨电价低谷时段充电虽然对用户便宜但对电网来说可能形成新的负荷尖峰。如果只做削峰填谷可能让用户承担更高充电费用没有经济性驱动力。我常用的目标函数是三部分加权min ∑_t price(t) × P_net(t) β ∑_i∑_t Pd_i(t) γ ∑_t (P_load(t) - P_avg)²第一项是购电费用price(t) 是分时电价第二项是V2G放电的电池退化惩罚β按照每kWh放电折算成电池损耗成本第三项是削峰填谷项用集群总负荷与平均负荷的方差来衡量。γ 的取值决定你更倾向经济性还是更倾向电网友好性我一般从小到大调几轮先看费用变化再定。3. 局部优化的三种“切法”时间、空间、场景3.1 时间切法滚动窗口加前馈修正最简单的局部优化思路就是别把24小时一次算完而是只优化未来一小段时间。这种方式和模型预测控制MPC的逻辑一致当前时刻基于最新观测信息对未来的 N 个时段比如4小时16个时段做一次优化只执行第一个时段的决策到下一个时段再重新求解。时间切法对随机性的包容性很强。因为车主实际到达时间和你预测的不一致你可以等车真正接入后再把它的信息纳入模型不需要提前假设所有车都准时出现。滚动窗口也天然具备反馈修正能力某辆车晚到了之前分配的功率已经没用重新求解时系统会自动把它从队列里拿掉。在MATLAB里实现滚动窗口其实不难核心是一个while循环外层推进时间内层每次只优化一个短时域子问题。3.2 空间切法按接入点分组用虚拟价格协调空间切法是处理“大规模”最直接的方式。一个几百辆车的园区通常不会只有一个接入点而是分散在若干台区、若干栋楼、若干根馈线下。这正好提供了天然的分组边界。每一组的内部可以独立优化但组与组之间共用一个上级变压器容量约束。我用的协调方式是虚拟电价迭代初始化 λ(t)0表示不施加任何额外约束。每个集群按照 price(t)λ(t) 作为自己的“综合电价”来求解子问题。汇总各集群的功率曲线检查是否超出上级总容量。如果某个时段超了就提高该时段的 λ 值让该时段的充电变得更“贵”各子问题自然会减少这个时段的功率。反复迭代直到总功率曲线收敛到容量限值以内。我通常迭20到40次就差不多了速度很快。这种方法的优点是不需要改造求解器子问题还可以并行计算。缺点是 λ 的调整参数需要试几次调整比例太大容易震荡太小收敛太慢。3.3 场景切法蒙特卡洛抽样加代表性场景聚类随机场景如果全部纳入优化计算量承受不住。我一般先用蒙特卡洛抽样生成1000个场景然后用K-means聚成10到20个代表性场景每个场景带一个概率权重。做法是把每个场景看成一条多维向量包含到达时间、初始SOC、离开时间等在标准化之后做K-means聚类。聚类完的每个簇中心就是一个“代表性场景”优化时只针对这10到20个场景做期望值计算。这样做的损失很有限因为大多数随机场景的行为模式是相似的最后体现出来的最优决策也相差不大。我实测把1000个场景缩减到15个之后优化结果的变化不超过2%但计算时间从几十分钟降到了几分钟。3.4 三种切法怎么组合实际项目里这三种切法不是互斥的而是层层配合。我常用的组合是外层按空间切分先按接入点分成若干集群每个集群一个子问题。中层按场景切分每个集群内部跑代表性随机场景按概率加权得到期望目标。内层按时间切分子问题时域不是96个时段全上而是用滚动窗口只算未来16个时段。三层切完原来的千万变量大模型就变成了一堆几千变量的小模型而且小模型和小模型之间还能并行。把三层组合起来之后我在MATLAB里用相对普通的硬件就能做到5000辆车秒级出方案。4. MATLAB代码骨架从数据模拟到协调循环4.1 随机接入数据生成下面的代码用于生成一晚上的车辆接入数据。重点是把随机性参数模型化方便后面反复替换数据集。% 基础参数 dt 15/60; % 15分钟单位小时 T 96; % 一天96个调度时段 n 1000; % 车辆数 E 60; % 电池容量 kWh P_max 7; % 最大充电功率 kW D_max 5; % 最大放电功率 kW eta 0.95; % 充放电效率 % 生成分时电价元/kWh示例曲线 price [zeros(1,32), 0.62*ones(1,24), ... 1.08*ones(1,28), 0.62*ones(1,12)]; % 随机变量生成 arrival_hour 18 2.5*randn(n,1); % 到达时刻 arrival_hour max(16, min(22, arrival_hour)); duration_h 6 3*randn(n,1); % 接入时长 duration_h max(2, min(14, duration_h)); soc_init 0.25 0.3*rand(n,1); % 初始SOC soc_target 0.85 0.1*rand(n,1); % 离开时目标SOC % 把到达时刻折算成时段编号 arrival_idx round((arrival_hour - 16) / dt) 1; depart_idx arrival_idx round(duration_h / dt) - 1;这里有两个容易踩的细节一是时段编号和小时之间的换算稍不注意就会差一个时段后面容量约束全乱二是截断分布一定要做不截断的话会出现凌晨3点到家甚至负数的奇葩数据。4.2 单车充放电子问题的优化函数单车子问题可以直接调用 fmincon也可以用 linprog 线性化。如果要求严谨的充放电互斥我建议用 intlinprog 加0/1变量。下面给一个适合改造成linprog形式的目标函数骨架function [Pc, Pd, soc] solve_single_ev(soc0, soc_target, arrival_idx, ... depart_idx, price, lambda, E, P_max, D_max, eta, dt, T) % 单车优化给定虚拟电价 lambda(t)求充电/放电计划 % 之所以把 lambda 直接加到电价上是因为上层协调就是靠这个信号 % 引导单车避开高峰时段 % 变量定义Pc(1:T), Pd(1:T)。实际优化窗口取 [arrival_idx, depart_idx] % 目标函数min sum( (pricelambda) .* (Pc - Pd) ) beta * sum(Pd) % 核心约束 % 1) SOC递推soc(t1) soc(t) (eta*Pc(t) - Pd(t)/eta) * dt / E % 2) SOC(1)soc0, SOC(depart_idx)soc_target % 3) 0PcP_max, 0PdD_max % 4) 到达前和离开后不调度arrival前与depart后功率为0 % 5) 可选互斥约束用二进制变量 % 该函数返回的Pc、Pd就是该车在协调循环中上报给集群的功率曲线 end注意我在函数里让 lambda 直接参与电价计算这样整个协调循环的代码非常干净子问题根本不知道全局容量约束的存在它只知道“现在这个时段的综合电价涨了我少充一点更划算”。4.3 集群功率分配的协调主循环主循环负责把多个单车子问题的结果汇总并不断调整 λ。比例积分形式的调节比纯比例调节更稳我实际用下来效果不错。% 集群参数 n_group 10; % 集群数 P_sub_limit 800 * ones(1,T); % 集群总功率上限 lambda zeros(1,T); % 虚拟电价增量 lambda_sum zeros(1,T); % 积分项 Kp 0.05; Ki 0.01; for iter 1:40 demand zeros(n_group, T); for k 1:n_group g groups{k}; % 第k个集群的车辆编号 for i 1:length(g) ev evs(g(i)); [Pc, Pd, ~] solve_single_ev(ev.soc0, ev.soc_target, ... ev.arrival_idx, ev.depart_idx, price, lambda, ... E, P_max, D_max, eta, dt, T); demand(k, :) demand(k, :) Pc - Pd; end end total_demand sum(demand, 1); overload total_demand - P_sub_limit; % 比例积分修正超限时段抬高虚拟电价 lambda_sum lambda_sum overload; lambda max(0, lambda Kp * overload Ki * lambda_sum); % 如果超限量已经很小可以提前退出 if max(abs(overload)) 1e-3 break; end end这个循环每迭代一轮所有车都要重新求解一次所以单车子问题必须写得足够轻。好在每辆车的优化窗口都很短问题规模小40轮迭代跑几千辆车在我的机器上也就是十几秒。4.4 求解器选型对比不同规模、不同约束形式适合不同的求解器。我做了个简单的对比表方便你按自己的场景选求解器/工具适用情况优点容易踩的坑fmincon连续变量目标可导MATLAB自带接口简单互斥约束要额外处理非线性约束多了容易陷入局部解linprog线性目标、线性约束速度快代码简单充放电互斥无法严格表达需用惩罚近似intlinprog需要0/1互斥变量能精确表达互斥整数变量一多求解时间暴涨YALMIPGurobi中小规模MIQP建模速度快求解器强需要单独安装License处理麻烦自写协调子问题大规模分解计算可控可并行调试λ参数需要经验代码量大我个人的建议是验证算法逻辑时用YALMIPGurobi写个小规模MIQP做基准真正跑大规模对比时用“协调循环单车intlinprog”这种分解结构。两条线互相印证既能保证正确性又能保证性能。5. 实测对比与踩坑实录5.1 50辆、500辆、5000辆的实测数据我在自己的算例环境中分别跑了三个规模环境是MATLAB R2023aCPU是普通的8核桌面级处理器。每个规模都跑同一套随机数据集和同一套参数只是车辆数不同。车辆数集中式MIQP求解时间局部优化求解时间与全局最优解的差距相比无序充电的费用下降5012.4s0.4s0.8%18.2%500约1180s3.6s2.1%15.6%5000内存溢出/超过3600s19.8s无法直接对比14.9%从数据能看出一个明显趋势车辆数越多集中式求解的时间呈超线性增长而局部优化时间的增长基本接近线性。50辆车时全局优化还能忍500辆时已经是按分钟起步5000辆时基本不可用。而局部优化的费用表现跟全局最优的差距一直在几个百分点以内这个代价换回几十倍的性能提升工程上完全划算。另外局部优化在削峰填谷上的效果也可以量化50辆规模下集群最大负荷从无序充电的约480kW降到约360kW降幅约25%5000辆规模下最大负荷的降幅大概在30%左右。这种效果主要来自两个机制一是把充电负荷往电价低谷时段搬二是通过虚拟电价把不同集群的充电时间错开。5.2 踩坑一效率不对称导致能量的“凭空出现”如果充放电效率都用同一个值比如0.95那还好。但如果充电效率是0.95、放电效率也是0.95而你的SOC递推公式写成 SOC(t1)SOC(t)(Pc-Pd)*dt/E那就会出现一个bug同一度电从电池放出来再充回去能量居然守恒了没有损耗。这在目标函数上表现为“钻空子”——优化器会让电池在一个低电价时段充电再在高电价时段放电加一点点损耗循环套利。正确的写法是充电和放电分开折算公式里我前面已经给了。实际执行时还需要注意效率放在哪一侧我习惯把充电效率乘在充电功率上把放电效率除在放电功率上这样SOC从电池角度始终是“净变化”。调试时可以在单车模拟里加一行断言SOC增加值和充电电量之比应该约等于充电效率避免低级错误。5.3 踩坑二Big-M和0/1互斥变量带来的数值抖动用整数变量处理充放电互斥时很多人喜欢写成 Pc(t) ≤ M × z(t)Pd(t) ≤ M × (1-z(t))。Big-M如果取得太大比如直接取10000会让求解器在数值上出现病态整数变量一直抖求解时间暴涨。我的经验是 Big-M 只需要比单时段最大可能能量值大一点点就行。比如功率上限7kW时段长度15分钟能量就是1.75kWhBig-M取10已经是绰绰有余。甚至可以不用Big-M改用intlinprog对每时段设置同一组功率变量加一个互斥约束代码更干净。另一种更工程化的处理是先忽略互斥约束求解连续问题再看结果里有没有同时充放电的小时段——一般出现这种情况是因为价格差和损耗惩罚之间的边际收益太小实际影响不大可以强行清零小功率值。这个近似在快速原型阶段很高效但在正式论文或报告里要说明清楚。5.4 踩坑三场景数量不削减随机性直接把计算打崩我第一次跑随机版本时直接用1000个场景叠加结果本地跑了一个多小时还没结束最后内存直接不够用。后来才意识到随机性场景不能“一锅端”进去必须先做场景削减。用K-means做场景削减时有个细节需要注意聚类前一定要对变量做标准化。因为到达时间、初始SOC、离开时间这三者的量纲差异很大不标准化会导致聚类结果完全被数值大的变量主导。标准化之后再聚类一般选15到20个簇就够用了。还有一个容易忽略的问题代表场景的初始SOC要直接用簇中心值而不是把簇内所有场景的SOC平均后再舍入。平均后一般会产生小数本身没问题但一定要保证场景内部的时序约束一致否则后续优化会反反复复不收敛。5.5 关于随机种子和结果复现这个坑我觉得值得单独拎出来说。两套算法做对比时如果随机种子不一致两批车辆数据本身就差了一截最终费用差距自然包含了数据差异算法差异反而看不出来。我现在的习惯是全项目统一用 rng(2024) 这样的固定种子生成所有随机场景跑不同算法时用完全相同的EV接入数据。这样算出来的费用差距才能干净地反映算法本身的差异。这个习惯也让我在排查优化器是否发疯时省了很多时间——分配方案突然异常先检查数据是不是被覆盖了再去怀疑算法。6. 后续演进从离线优化到在线MPC再到分布式协调这套局部优化框架做好之后我很快就想到了下一步。当前代码是把“未来一天的预测信息”当作已知来做计划但真实场景中预测不断更新所以最自然的扩展是把滚动窗口和集群协调结合起来做成在线MPC每分钟或每15分钟重新读取当前接入车辆状态重算未来4小时计划只执行下一步。这样司机临时不来、临时早走、充电桩故障等情况都能被系统自动吸收。另一个可以接的方向是分布式ADMM。如果园区规模大到单台机器都扛不住5000辆车或者担心单点计算瓶颈可以把每个集群的优化任务下放到边缘设备上集群之间只交换“边界功率”或者“虚拟电价”这类少量信息。ADMM的收敛速度比纯对偶分解快不少代价是代码复杂度又上了一个台阶。最后还想说一句关于数据质量的体会。优化算法做到后来决定结果上限的往往不是求解器而是输入数据的准确性。车主到达时间分布、SOC初始值的估计、常规负荷曲线这些如果偏差太大任何高级优化都救不回来。所以我在每个项目里都会先花时间做数据分析和分布拟合确认随机变量模型靠谱之后才开始调优化器参数。我的个人体会是先把车的数据分布估计准比先调优一个复杂的优化器更值得投入精力。