1. 项目概述
搞电力系统优化调度的同行,最近大概率都躲不开“高比例新能源”和“电动汽车大规模接入”这两个关键词。单看任何一个,优化模型尚能应付;一旦把它们揉进同一个并网型微电网里,还要面对风电、光伏、负荷、充电行为的多重不确定性,传统确定性调度模型就明显力不从心了。这篇博文要拆解的项目,正是用Matlab实现的一套考虑不确定性的含集群电动汽车并网型微电网随机优化调度代码,核心解决“源—荷—车”三方不确定条件下,如何经济、安全地安排微电网内各机组出力与电动汽车充放电策略的问题。
先交代一下这个项目能干什么。它面向的是含风、光、储能、微型燃气轮机以及集群电动汽车的并网型微电网,通过场景法描述不确定性,构建两阶段或单阶段随机优化模型,决策变量涵盖机组出力、储能充放电功率、电动汽车充放电状态以及与主网的交换功率。最终输出的是各时段的最优调度方案、运行成本曲线以及电动汽车聚合体的充放电行为结果。适合的人群非常明确:正在做微电网优化调度方向毕业设计或论文复现的电气工程研究生、刚入门综合能源系统优化且需要一套可运行参考代码的工程师,以及想在Matlab环境里快速验证随机优化思路的科研人员。
为什么这类项目这几年特别“香”?因为并网型微电网本身就有“大电网兜底”的优势,可以灵活买卖电;但电动汽车集群又给系统增加了一个全新的可平移负荷维度——它既是负荷,又是潜在的储能资源。如果调度策略不考虑充电行为的随机性,高峰时段可能面临变压器过载,低谷时段又可能闲置容量。而随机优化调度,本质上就是“为多种可能发生的未来场景各准备一套预案,然后选一个综合下来最不亏的方案”。这套代码的价值,正是把这种思路从理论公式落地成可运行的Matlab程序。
整体来说,这套项目并不算“大而全”的商用电网调度系统,而更像一个聚焦于算法与模型验证的科研型工具。它的建模粒度集中在日前调度层面,时间尺度通常为24小时、以1小时为间隔,重点回答以下几个问题:风电和光伏出力有哪些典型场景?电动汽车集群的充放电功率边界怎么聚合?随机优化相比确定性优化到底省了多少成本?每一问在代码里都有对应的模块实现。下面我从设计思路到核心代码逻辑,再到运行调试的避坑经验,完整把这个项目拆开讲透。
2. 调度模型整体设计与思路拆解
2.1 为什么必须引入随机优化而不是继续用确定性模型
很多初接触这个方向的人,第一反应是把风电、光伏和负荷全部取预测期望值,然后套一个混合整数线性规划(MILP)或者粒子群算法去求解。这种做法在只有单一不确定性源的小系统里或许问题不大,一旦加入电动汽车集群,漏洞就暴露出来了:EV充电负荷的时间分布和功率大小受用户出行习惯影响极大,早高峰和晚高峰的充电需求可能相距悬殊,如果只按一个固定预测曲线安排调度,下午时段可能会因为实际充电功率超过预期而被迫从主网购高价电,甚至触发切负荷。
随机优化解决的核心矛盾是“决策必须做出,但不确定性尚未实现”。以风电为例,确定性模型默认风速预测准,据此安排燃气轮机出力;随机模型则生成多个风电场景——风大的场景、风小的场景、平风的场景——每个场景对应不同的出力调整需求。目标函数变成各场景成本的概率加权期望。这种“以场景求期望”的思路,非常类似投资组合里“不把鸡蛋放在一个篮子”的理念:宁可今天多花一点备用成本,也不要在恶劣场景下付出高昂的惩罚代价。
具体的数学表达,简单说就是:外层决策(如机组启停、EV充放电状态)先定下来,内层再根据每个随机场景调整出力与交换功率。如果采用两阶段随机规划,第一阶段变量通常被称为“here and now”,必须在不确定性实现前决定;第二阶段变量则是“wait and see”,可以随着场景变化。这个项目代码里,我见过的最常见实现是随机场景约束下的单阶段期望成本最小化,把场景集直接展开到约束中,等价于一个大规模线性规划;也有做得更细的,加入机组启停0-1变量,变成一个大规模MILP。Matlab环境下这类问题通常用Yalmip工具箱建模,求解器选Cplex或Gurobi,如果是线性问题单纯用linprog也能跑。
2.2 集群电动汽车建模:从单体到聚合的关键处理
电动汽车不能简单当成“一堆充电桩负荷”。因为每辆车的电池容量、初始SOC、到达离开时间、充电功率都不同,如果逐车建模,决策变量会爆炸增长,场景数稍微一多,求解器直接内存溢出。因此代码里普遍采用集群聚合模型,即将电动汽车集群看作一个“可调度的电池储能系统”,只是在充放电功率边界和能量边界上,需要根据车辆统计信息动态折算。
具体有两种聚合思路。第一种是基于蒙特卡洛模拟的聚合:先抽样生成每辆车的起始SOC、到达时间和驻留时长,统计得到集群在各时段的可用充放电功率上限和可调度容量范围。第二种是基于等效电池模型:把整个集群的聚合容量定义为所有车辆电池容量之和乘以一个同时率系数,充放电功率上限由配电网或停车场容量限制决定。实际项目中两种方法各有优劣,前者统计信息丰富但计算量大,后者简洁但精度稍逊。这套代码我建议采用第一种,因为论文里往往需要画出“集群充放电可行域”的图,而蒙特卡洛聚合正好能给出上下边界曲线,图表丰富度对写作十分有利。
再强调一个细节:电动汽车集群在调度模型里通常是“双状态”的,即既能充电又能放电(V2G),但实际运行中电池寿命损耗是个经济性问题。因此优化模型里的EV充放电成本不能按零处理,至少要在目标函数中加上电池退化成本项,或者设置充放电服务费系数。否则求解器会疯狂让EV放电来赚取峰谷价差,得到看似成本极低、实则完全不可行的调度方案。这一条很多人会漏掉,后面我在代码实操部分还会再提。
2.3 并网型微电网架构:各单元的角色与配合逻辑
一个典型的含集群EV并网型微电网,物理架构由这么几块组成:风力发电机(WT)、光伏阵列(PV)、微型燃气轮机(MT)、储能系统(BESS)、集群电动汽车(EV),以及并网点(PCC)连接上级配电网。调度目标通常是最小化总运行成本,总成本里包含MT燃料成本、储能折旧成本、EV充放电损耗成本、向主网购电成本;如果有弃风弃光惩罚或切负荷惩罚,也要加进去。
各单元的角色定位值得展开说。风电和光伏是“看天吃饭”的不可控电源,在随机优化里作为场景输入;微型燃气轮机是可控主力,响应速度快但燃料成本高,适合承担基荷和爬坡调节;储能是“缓冲器”,白天存过剩新能源、晚高峰放出,同时也能平抑EV充电带来的短时功率冲击;EV集群是“柔性负荷兼移动储能”,调度得好可以在电价低谷充电、高峰放电,但必须尊重用户出行需求,比如离开时SOC不能低于用户期望值。
这几者的配合逻辑,用一个生活场景类比就是:新能源发电像“看心情上班的兼职员工”,微型燃气轮机是“随时待命的正式员工”,储能是“公司仓库”,EV是“既能进货也能出货的运输车队”。随机优化的目标不是让某个单元最省钱,而是让整个团队在任何“天气剧本”下都能把活干完,且总工资支出期望最小。
3. 不确定性建模与场景生成的核心细节
3.1 风电、光伏出力场景怎么生成才靠谱
不确定性建模的第一大坑,是直接把历史数据拿来加噪声当场景。风电出力的概率分布是偏态的,且具有明显的时序相关性,简单用正态分布扰动会生成大量不合理场景——比如正午光伏出力为负、夜晚风电满发的极端样本。这套项目代码里,我常用的做法是拉丁超立方采样(LHS)+ Cholesky分解处理相关性,或者直接用历史典型日的出力曲线叠加预测误差分布。
具体步骤拆开是这样的:先获得风电、光伏预测出力的期望曲线,设定预测误差服从正态分布或Beta分布(光伏用Beta分布更贴切),然后通过LHS在概率空间分层采样,得到一组初始误差样本。由于风电和光伏之间存在负相关性(有风时往往云多,光伏出力低),还要用协方差矩阵结合Cholesky分解对样本进行相关性矫正。最后把期望曲线加上误差样本,截断到[0,额定功率]区间,就得到一组初始场景。
场景数量直接决定优化模型规模。如果初始生成500个场景,直接带入约束求解,变量数和约束数都会膨胀到难以收敛。因此必须做场景削减。最经典的是基于概率距离的同步回代削减法,核心思想是:计算各场景之间的 Kantorovich 距离,反复合并距离最近的两个场景,直到剩余场景数满足设定值。Matlab实现中,可以调用SceneReduction这类函数包,或者自己写一个30行的迭代程序。削减后典型场景一般取10~20个,既能反映不确定性分布,又不至于让求解器卡死。
3.2 电动汽车充电负荷的不确定性与可调度性边界
EV集群的不确定性体现在两个层面:一是充电负荷出现的时间和位置不确定,二是每辆车可提供的调节能力不确定。在日前调度框架下,通常用蒙特卡洛模拟生成EV群体的到达时间、起始SOC、目标SOC和驻留时间,统计得到集群在各时段的平均充电功率与可放电功率边界。
这里的关键参数是EV可调度容量时变曲线。举个例子,假设一个住宅区停车场有300辆车,下午5点下班陆续接入,晚8点达到接入峰值,次日早上7点陆续离开。那么晚上8点到次日6点这个时间段,集群的整体电池容量是可用的;白天的可调度容量就非常低。代码里应该输出两条曲线:一条是各时段集群可充电的最大功率P_ev_cha_max(t),一条是可放电的最大功率P_ev_dis_max(t),这两条曲线由停车场容量、接入车辆数和平均充电功率共同决定。
不少初学者会犯一个错误:把EV集群当成普通储能,忽略“SOC必须在用户离开前达到目标值”这一约束。这会导致调度方案在凌晨让EV大量放电,早上用户开车时电量不足,完全脱离实际。正确做法是在聚合模型中增加一个日总能量平衡约束:集群在调度周期结束时的总SOC不得低于初始值,或者不得低于用户期望阈值。同时,每个时段的荷电状态上下限要跟随接入车辆数动态调整,不能简单用0~1。对于这种细节,你可以在代码里用矩阵方式一次性生成24小时的可调度上下限,避免循环里反复判断。
3.3 场景削减的参数选择与Matlab实现建议
场景削减虽然是个“配角”环节,却经常决定随机优化的成败。削减后的场景数太少先,优化结果对不确定性描述不足;太多则计算时间成倍增加。我在实际测试中,初始场景500个、削减到20个时,优化结果与500个场景的全模型相比,成本偏差通常能控制在1%以内,但求解时间从半小时降到几十秒,性价比极高。
Matlab实现场景削减,有两种路线。一是用现成的fastForecast或SyncTrap等第三方函数;二是自己写一个简单版本。自己写的思路是:每个场景有概率权重,初始均等;循环里找到欧氏距离最小的两个场景,把概率较小的场景删掉,它的概率叠加到概率较大的那个;反复执行直到场景数达标。这个过程本质上和K-means聚类中的最近邻合并类似,但保留的是原始场景点而不是聚类中心,因此物理意义更明确。
对风电、光伏联合场景削减时,有一个容易踩的坑:距离度量不能只考虑同一时刻的出力差,还要考虑时序形状。如果两个场景在24个时段上的出力整体偏高,但走势相反,用欧氏距离会误判为“相似”。更好的办法是给距离矩阵加上一阶差分项,或者使用动态时间弯曲(DTW)距离来度量。当然,这会让计算量上升;对于论文级项目,用带加权系数的欧氏距离通常也够用,因为调度约束会过滤掉部分不合理结果。
4. 随机优化调度模型构建与Matlab实操
4.1 目标函数与约束条件的完整数学形式
在代码里,目标函数通常写作所有场景下运行成本期望的和。为了便于线性化,成本项会被拆成几块分别计算。以单阶段随机优化为例,目标函数可以表示为:
总成本 = 燃气轮机燃料成本 + 储能折旧成本 + EV充放电损耗成本 + 购电成本 - 售电收益 + 弃风弃光惩罚 + 失负荷惩罚
每项都要乘以对应场景的概率,然后对场景求和。这里燃气轮机的燃料成本通常用二次函数拟合,但在MILP里需要分段线性化。如果直接用非线性求解器(如fmincon),收敛性和全局最优性都会受到挑战,所以我更推荐在Matlab里用Yalmip建模、Gurobi或Cplex求解,这两款求解器对分段线性化的支持都很完善。
约束条件分为几大类:功率平衡约束、机组出力上下限约束、爬坡约束、储能SOC递推及容量约束、EV集群功率与能量约束、与主网交换功率约束。功率平衡约束是核心,要求任一时刻所有电源出力加购电功率,等于所有负荷加EV充电功率加售电功率,注意EV放电时要作为电源侧处理。由于是随机模型,每个场景下功率平衡都要满足,这意味着约束矩阵的行数要乘以场景数,建模时需特别注意索引对应关系。
4.2 Yalmip建模的关键代码片段解析
直接上一个简化的代码骨架,方便大家对照自己的项目修改。注意这段代码偏向示意,实际运行时需要根据你的具体参数调整维度。
% 假设已生成场景数据: P_wind(s,t), P_pv(s,t), P_load(s,t) % 决策变量 P_mt = sdpvar(1, 24, 'full'); % 燃气轮机出力 P_bess = sdpvar(1, 24, 'full'); % 储能充放电功率(正充负放) P_ev = sdpvar(1, 24, 'full'); % EV集群充放电功率(正充负放) P_grid = sdpvar(1, 24, 'full'); % 与主网交换功率(正购负售) SOC_bess = sdpvar(1, 25, 'full'); % 储能SOC,时段0~24 SOC_ev = sdpvar(1, 25, 'full'); % EV集群SOC % 场景相关的决策变量(第二阶段) P_mt_s = sdpvar(n_scene, 24, 'full'); P_grid_s = sdpvar(n_scene, 24, 'full'); P_ev_s = sdpvar(n_scene, 24, 'full'); Constraints = []; for t = 1:24 for s = 1:n_scene % 功率平衡 Constraints = [Constraints, ... P_mt_s(s,t) + P_bess(t) + P_ev_s(s,t) + P_grid_s(s,t) ... == P_load(s,t) + P_wind(s,t) + P_pv(s,t)]; % 注意符号定义 end end % 再补充机组上下限、爬坡、SOC递推约束... % 目标函数为场景概率加权成本期望 Objective = sum(prob .* sum(C_fuel + C_grid + C_bess + C_ev_deg));上面这个骨架里,我把储能出力设为与场景无关(第一阶段决策),EV出力设为场景相关,这是一种折中处理。实际中储能也可以在不确定性实现后再调整,需根据你的模型设定决定。有一点必须指出:代码里的功率平衡约束我写了等号左边是“发电侧”还是“负荷侧”?这个符号问题是最容易出错的地方。建议你在草稿纸上先画清楚节点功率流向,再写代码;否则正负号错了,求解器会给出极其离谱的调度结果,而且这种错误非常隐蔽,排查时令人崩溃。
4.3 求解器配置与运行调参经验
无论你用Yalmip还是直接用MATLAB优化工具箱,求解器配置都是影响成败的关键。如果模型是线性的,Gurobi通常只需几秒到几十秒;如果含有0-1变量(比如机组启停),则属于MILP,Gurobi的MIP gap设置尤其重要。我在跑类似项目时,一般会设置optimize(Constraints, Objective, sdpsettings('solver','gurobi','gurobi.MIPGap',0.01)),这里的1%相对间隙足够应付论文级精度,又能避免求解器在最优解附近死磕浪费时间。
另一个容易被忽视的参数是求解时间限制。遇到大规模场景时,MILP求解时间可能指数增长,跑半小时不出结果也是常事。建议设置gurobi.TimeLimit为300秒,并开启gurobi.NumericFocus为2或3,这能有效缓解约束矩阵尺度差异带来的数值问题。什么是数值问题?比如功率平衡约束的量级是兆瓦,SOC约束的量级是0到1,购电价格是几百元每兆瓦时,这些数值相差悬殊,求解器在计算内点时会遇到病态矩阵。解决办法是对模型进行归一化,或者设置合适的变量缩放。实践里我会把功率变量统一用兆瓦为单位,成本用千元为单位,SOC保留0~1,这样能大幅减少数值警告。
4.4 结果可视化的几个必画曲线
随机优化做完,光得到一个总成本数字不算完事。我写这类论文和报告时,有四张图是必须画的,这既是分析需要,也是向导师或审稿人展示建模充分性的直接证据。
第一张是场景削减前后的概率分布图,用折线或区域图展示风电/光伏典型场景集,旁边标注每个场景的概率。这张图能直观说明你的场景生成和削减不是拍脑袋。第二张是最优调度结果堆叠图,以时间为横轴,把燃气轮机出力、风电光伏消纳量、储能充放电、EV充放电、购售电功率堆叠起来,一看就能分析高峰时段谁在扛、低谷时段谁在吸收。第三张是EV集群SOC与充放电功率曲线,重点检查EV有没有在用户需求时段之外被强行放电,SOC是否始终处于合理区间。第四张是不同场景下的总成本对比柱状图,以及随机优化与确定性优化的成本差距对比,这是论证“考虑不确定性价值”的核心图表。
画图用Matlab自带plot、area、bar就够用。颜色建议统一风格,不要花哨;坐标轴标签一定写清楚单位,否则自己过两天回看都会懵。还有一个小窍门:在堆叠图里,用不同透明度区分预测值和场景值,论文里可以节省大量解释成本。
5. 常见问题与调试实录
5.1 求解结果荒谬:成本为负或者EV疯狂放电
这个问题我见过太多次,基本可以确定是目标函数中EV放电收益系数设置不当,或者缺少EV日总能量约束。如果EV放电没有惩罚成本,优化器会发现“低价充电、高价放电”是无本套利,于是让EV在低谷充满、高峰放空,甚至一天循环多次,最后成本低得离谱。解决办法有两种:要么在目标函数中加入EV充放电的电池损耗成本系数(比如0.2元/kWh),要么硬性约束EV集群SOC在调度周期末等于初始值,防止“能量凭空创造”。如果两种都做了还出现负成本,那要检查购售电价曲线是否出现倒挂,并网型微电网如果售电价高于购电价,模型会疯狂套利,这种情况需要给售电功率设上限。
5.2 场景削减后调度结果对场景数量异常敏感
有些用户发现,场景数从10个增加到30个,优化成本跳变非常剧烈,曲线不收敛。这通常是场景削减方法的问题,比如合并准则里概率更新算错,或者削减后没有重新归一化概率。更隐蔽的原因是削减后的场景集里出现了“极端场景被抹平”的现象——本来该保留的高风电低负荷的尖峰场景被合并掉了,导致优化方案对突发情况的准备不足。
我的建议是,削减前先查看各场景的概率分布直方图,确认极端场景概率非零;削减后再对比削减前后的期望出力曲线,如果偏差超过5%,就需要调整削减策略。另外,可以用多个随机种子生成场景集,观察优化结果的波动范围。如果波动太大,说明你的不确定性描述本身有问题,而不是调度模型的错。
5.3 求解时间过长,无法在合理时间内完成
当场景数达到100以上、且含0-1变量时,Matlab里单纯靠Yalmip+Gurobi也容易卡到怀疑人生。我踩过几次坑之后总结出四招:第一,尽量把连续变量决策保持线性,避免使用abs或min等非线性函数,改用引入辅助变量的方式线性化;第二,把非关键约束的“大M”值设得尽量紧,避免求解器在无效分支上浪费大量节点;第三,使用模型切割,比如先固定机组启停方案,求解连续优化得到下界,再在分支定界框架内加速收敛;第四,考虑把24小时模型压缩成典型时段,或者使用滚动时域方法,牺牲一点全局最优性换取计算速度。
Matlab的prob系列语法(optimproblem)在R2021b以后对MILP的建模也很友好,调试更为顺滑,只是表达复杂约束时没有Yalmip方便。个人习惯是:标准流程研究用Yalmip,快速原型验证用optimproblem都行。关键是求解器选项里打开Display迭代过程,这样你能实时看到gap下降情况,判断是继续等还是赶紧改模型。
5.4 数据输入格式导致矩阵维度不匹配
随机优化代码里的维度地狱,是新手最容易崩溃的地方。风电场景数组是n_scene×24,EV聚合功率边界是1×24,燃气轮机爬坡约束需要P_mt(t) - P_mt(t-1),而场景相关变量则是n_scene×24。在Matlab中,sdpvar的维度与普通数组的维度必须严格一致,任何ones(1,24)和ones(n_scene,1)的误用都会直接报错。
遇到这种问题,我建议先花15分钟用size()检查所有变量的维度,写注释标明每个变量的行含义和列含义。一种很实用的编码习惯是:所有场景都放在第一维,所有时间都放在第二维,这样转置操作会少很多。另外,在设计约束循环时,尽量用矩阵操作代替双层循环,否则代码运行时间会指数上升。Matlab的矢量化技巧在这里用得上:对场景维度的求和、对时间维度的递推,可以用repmat、sum(,1)、cumsum等函数组合完成。
5.5 不确定场景的稳定性:怎么验证调度方案的鲁棒性
很多人在写完随机优化后不知道如何验证结果的可靠程度。至少应该做这样一件事:用未参与优化的原始场景集(500个)重新评估优化得到的调度方案,计算这个固定决策在全部原始场景下的期望成本和最大成本。如果最大成本远超期望成本,说明方案虽然平均表现好,但存在“爆表”风险,可能需要增加备用约束或缩短调度时段。
更严谨的做法是进行样本外测试:把历史数据分成训练集和测试集,用训练集生成场景参与优化,再用测试集数据模拟实际运行,记录成本分布。这项测试在论文里属于加分的“鲁棒性分析”,在代码里实现也不复杂,多写一个评估脚本即可。我自己在项目交付时,除了调度代码,还会附带一个evaluate.m脚本,专门做这种后评估,给使用者的体验会好很多。
6. 扩展思考与个人体会
聊到这儿,这套随机优化调度项目的主体已经完整呈现出来了。不过我还想多啰嗦几句从实际项目中得到的体会,虽然不一定直接写在代码注释里,但对真正想把这套方法用好的人来说,这些比单纯的代码片段更有价值。
关于不确定性建模的深度,要有分寸感。随机优化不是场景越多越好,也不是分布越复杂越好。很多论文把风速分布改成韦布尔三参数、把负荷相关性用Copula刻画,看似高端,但落到实际调度中,收益往往被计算代价抵消。初学者我建议先用正态误差+LHS+同步回代削减,把主流程跑通,再逐步增加复杂度。这个项目的代码实现也是如此,先保证模型能解出来、结果符合物理直觉,再去谈高阶概率工具。
关于电动汽车集群的角色定位,我觉得接下来几年会越来越重要。现在的模型大多把EV当作价格响应型负荷或者V2G储能,但真实世界中用户充电意愿、电池健康状态、充电桩调度协议都会影响集群的灵活性。也就是说,随机优化里的“不确定性”远不只风电光伏,EV本身的行为不确定性同样值得单独建模。如果你打算在这个方向上继续深耕,可以考虑把EV充电决策建模成马尔可夫决策过程,或者引入博弈论描述充电运营商与微电网之间的互动,这些都是有意思的延伸方向。
最后分享一个调试小技巧:在正式跑大规模场景前,先把场景数设为3个,其中一个设成极端高价场景,一个设成极端低价场景,一个设成正常场景。这样既能快速验证模型约束是否正确,也方便你手算核对部分结果。用这种“人工小场景”把逻辑理顺后,再切换到完整场景集,往往能省下一整天的排查时间。我个人几乎每个优化项目都会这么做,屡试不爽。希望这篇拆解能帮你把随机优化调度模型从“看得懂公式”推进到“跑得通代码”,真到了那一步,你回头看这套项目的价值感会完全不一样。