能源系统方向的硕士论文近几年翻来覆去就是那几个热门题“双碳”目标下的电力系统优化运行、分布式能源接入、需求响应、电动汽车有序充电……而“可再生能源发电与电动汽车的协同调度策略”几乎是把几个热点方向缝在了一块儿。我去年完整复现过几篇这个题目的论文从Matlab代码到论文里的所有图都重新跑了一遍。这个题目听起来玄乎其实就是两件事处理风电光伏出力的不确定性以及把电动汽车充电负荷变成可调度的弹性资源然后把它们放进同一个优化框架里求解。这篇博文不是贴代码就完事。我会把复现时整个思路拆开讲先搞清楚论文模型里到底优化什么再讲Matlab里怎么建模型、怎么调用求解器、怎么把论文里的图复现出来最后把调试时遇到的那些坑都列出来。正在选这个题的同学、准备做仿真的研究生可以参考一下能少走不少弯路。1. 为什么选择这个题目协同调度的需求背景与核心思路1.1 风电光伏的“出力看天吃饭”特性风电和光伏接入电网后最让人头疼的问题就是出力不确定性。风电出力取决于风速具有明显的季节性和时段性风速忽高忽低风机出力跟着剧烈波动光伏更“看天吃饭”白天有光才有出力晚上直接归零遇到阴天出力可能只有晴天的20%。这跟传统水火电机组“想发多少就发多少”的控制方式完全不一样。电力系统的基本要求是发电和负荷实时平衡。新能源大规模接入以后系统净负荷负荷减去新能源出力的波动被放大了传统的“源随荷动”模式越来越吃力。结果就是两个典型问题一是弃风弃光尤其是夜间风电大发而负荷低谷时只能限制风电出力二是需要火电机组频繁调整出力运行成本和碳排放都上去了。复现这类论文时我建议第一步就把新能源出力模型搞明白。硕士论文里常见做法是用Weibull分布模拟风速再用风机功率曲线计算风电出力光伏则用Beta分布描述光照强度或者直接用某地典型日辐照度数据。搞清楚这一层后面写Matlab代码才不会乱。1.2 电动汽车是“时间尺度上的移动储能”电动汽车的充电负荷跟传统负荷最大的区别在于时间弹性。车主的出行规律决定了充电需求的总量但充电起始时间、充电功率、充电时长是可以调整的。一台家用电动汽车一天24小时里实际行驶时间通常只有1到2小时其余时间都在停着。这意味着充电负荷在时间轴上有很大的平移空间。如果所有电动汽车都在傍晚下班后同时开始充电配电网的峰值负荷会被瞬间拉高变压器过载风险剧增。反过来如果让一部分车错峰充电或者利用电价信号引导车主把充电时间挪到夜间低谷期就能起到削峰填谷的作用。更进一步V2G技术还能让电动汽车在高峰时段反向放电相当于一个个小型储能电站。当然V2G在硕士论文里更多是作为展望方向实际建模通常只考虑有序充电。这个特性就是协同调度的核心机会点。要理解的是电动汽车本身不是“额外负担”而是可以被调度的灵活性资源。论文题目里的“协同”本质上就是把这种灵活性利用起来。1.3 协同调度到底在“协同”什么把前面两点结合起来看协同调度的本质是在一个统一的优化框架下同时考虑新能源出力的不确定性和电动汽车充电负荷的可调性求出一组最优的调度策略。这里的“协同”体现在三个层面。第一个层面是时间协同新能源出力高峰时引导电动汽车多充电新能源出力不足时控制充电负荷甚至让车反向送电。第二个层面是源网协同发电机组的出力计划和电动汽车充电计划要同时优化而不是分开做。第三个层面是多目标协同调度策略往往要同时兼顾运行成本最低、新能源消纳最大化、负荷波动最小等多个目标彼此之间可能存在冲突需要用加权或者罚函数的方式折中。读论文的时候最好先把作者的目标函数拆出来看。有的论文是成本最小化有的论文是负荷方差最小化有的把弃风弃光作为惩罚项塞进了成本函数。目标函数定了约束条件定了剩下的就是求解方法的选择——这就是后面Matlab实现的主要内容。2. 复现前的技术准备先读懂模型再打开Matlab2.1 拿到论文先做的不是找代码而是列公式很多人拿到源码就开始跑跑通了就觉得复现成功。实际上这是本末倒置。硕士论文复现的核心是理解模型代码只是模型的载体。我自己的习惯是先把论文的第三章通常是“调度模型”或“优化模型”逐句读一遍把每个公式抄在纸上标注清楚哪些是决策变量、哪些是已知参数、哪些是目标函数项、哪些是约束条件。一篇典型的协同调度论文模型结构大概是这样决策变量各时段火电机组出力、电动汽车充电功率可能还有储能充放电功率、弃风弃光量目标函数运行成本最小化包括火电燃料成本、弃风弃光惩罚成本、电动汽车充电成本等约束条件功率平衡约束、机组出力上下限、爬坡约束、EV充电功率约束、SOC约束、传输线容量约束等把这张表列出来之后再去对照代码就能看出作者是怎么把一个数学模型翻译成Matlab语言的。这一步做扎实了后面改参数、改约束、扩展到自己的场景都不难。2.2 Matlab环境配置与工具箱选型Matlab环境这块步骤其实不复杂但有几个容易栽跟头的地方。先说工具箱。做调度优化至少需要以下工具箱Optimization Toolbox提供linprog、fmincon、intlinprog等标准求解函数Statistics and Machine Learning Toolbox生成随机数、拟合分布需要它比如wblrnd、normrndGlobal Optimization Toolbox如果论文用了遗传算法、粒子群算法就会用到ga、particleswarm版本上R2020a以上基本就够用不追求最新版。我复现时用的R2022b整体比较稳定。再说求解器。很多论文用的是YALMIP框架加CPLEX或Gurobi求解器。YALMIP是一个建模工具箱可以把你用数学公式写出来的优化模型自动转换成求解器需要的格式比直接调用linprog要直观得多。CPLEX和Gurobi是商业求解器学生可以申请免费学术许可求解混合整数线性规划MILP的效率远高于Matlab内置求解器。安装YALMIP很简单去GitHub或者官网下压缩包解压后把整个文件夹添加到Matlab路径即可。CPLEX则需要注册学术账号下载对应版本安装后在Matlab里配置路径。2.3 代码文件结构不要堆成一个巨型脚本复现论文不是只写一个脚本就完事。一个大几百行的脚本你当时能看懂三天后回来看就完全不知道哪里是哪里了。我建议按这个结构组织代码project/ ├── data/ # 输入数据负荷数据、风速数据、EV参数 ├── code/ │ ├── generate_wind.m # 风电出力生成 │ ├── generate_pv.m # 光伏出力生成 │ ├── generate_ev_load.m # EV充电负荷生成 │ ├── build_model.m # 构建优化模型 │ ├── solve_model.m # 求解 │ └── plotting.m # 出图 └── results/ # 运行结果这个做法的好处是边界清晰。改参数不会动到模型逻辑改模型不会动到数据生成逻辑。后面做多场景对比时特别方便——你只需要改data里的参数其他代码不用碰。3. 核心代码实现从风光出力到优化求解的完整拆解3.1 风电与光伏出力建模的常用套路风电出力建模的常见流程是先用Weibull分布抽样得到风速序列再按风机功率曲线计算功率。Weibull分布的两个参数是形状参数k和尺度参数c典型取值是k2、c8m/s这个组合在论文里非常常见。T 24; % 调度时段数按小时 k 2; c 8; % Weibull 参数示意 v wblrnd(c, k, [1, T]); % 生成风速序列 % 风机功率曲线参数 v_cutin 3; % 切入风速 m/s v_rated 12; % 额定风速 m/s v_cutout 25; % 切出风速 m/s P_w_rated 200; % 风电场额定容量单位略 P_w zeros(1, T); for t 1:T if v(t) v_cutin || v(t) v_cutout P_w(t) 0; elseif v(t) v_cutin v(t) v_rated P_w(t) P_w_rated * (v(t)^3 - v_cutin^3) / (v_rated^3 - v_cutin^3); else P_w(t) P_w_rated; end end这段代码里有个物理常识必须写对风速低于切入风速或者高于切出风速时风机出力为0在切入风速和额定风速之间风功率与风速的三次方成正比。很多初版复现会漏掉切出风速的判断导致出现“大风出力为0、中等风速出力反而高”的荒谬曲线。光伏出力建模相对直接一点用辐照度数据乘上光电转换效率和面积就行。硕士论文里通常用Beta分布抽样辐照度或者直接给一条典型日曲线。前者更学术后者更省事。需要注意辐照度为零的夜间时段光伏出力一定为零这个判断逻辑要写在代码里否则后面优化模型里容易出现负的光伏出力或者莫名其妙的极端值。3.2 电动汽车充电负荷的蒙特卡洛模拟EV负荷这一块主流做法是蒙特卡洛模拟。每台车的数据包括到家/接入充电桩时间、第二天离家时间、起始SOC、充电功率、电池容量。论文里一般会假定三个随机变量的分布规律到家时间近似正态分布均值18点左右标准差2到4小时离家时间近似正态分布均值8点左右标准差1到2小时日行驶里程影响起始SOC常用对数正态分布或直接给统计值num_ev 200; % 电动汽车数量 P_ev zeros(num_ev, T); % 每台车各时段充电功率 SOC_min 0.2; SOC_max 0.9; battery_capacity 60; % kWh charge_power 7; % kW for i 1:num_ev % 采样接入时间、离开时间 t_arrive round(normrnd(18, 3)); t_depart round(normrnd(8, 2)) 24; % 次日上午 % 采样起始SOC s0 max(0.2, normrnd(0.4, 0.15)); s0 min(s0, SOC_max); need_energy (SOC_max - s0) * battery_capacity; % 需要充入的电量 need_time ceil(need_energy / charge_power); % 需要的充电时长 % 无序充电接入即充 for j t_arrive : min(t_arrive need_time - 1, t_depart - 1) if j T P_ev(i, j) charge_power; elseif j T % 跨天处理映射到第二天凌晨 P_ev(i, j - T) charge_power; end end end P_ev_load sum(P_ev, 1); % 汇总为系统EV总负荷这段代码里有个细节值得注意跨天时段。车辆如果在23点之后接入充电桩充电行为会延续到第二天凌晨直接循环到第二天凌晨的时段。很多初版复现代码没有处理这个边界导致凌晨时段的EV负荷少算了。蒙特卡洛的收敛性是另一个坑。如果车辆数太少比如50辆单次随机样本的波动会很大导致结果不稳定。一般建议至少跑1000辆车或者多次重复采样取平均。但考虑到调度模型本身计算量论文里常见的是直接给定EV负荷曲线或者用少量场景代表这也可以在代码里通过设置随机种子保证结果可复现。3.3 优化模型用YALMIP建思路更清晰优化模型的实现是复现工作的核心。我推荐用YALMIP因为它的语法和课本上的数学公式几乎一一对应。先声明决策变量再写目标函数再写约束最后调用求解器。以一个相对简单的模型为例目标函数为系统运行成本最小化包括火电燃料成本、弃风弃光惩罚成本、EV充电成本T 24; P_g sdpvar(1, T); % 火电机组出力 P_w_use sdpvar(1, T); % 实际使用的风电 P_pv_use sdpvar(1, T); % 实际使用的光伏 P_ev sdpvar(1, T); % EV充电总功率 % 火电成本系数 a 0.001; b 20; c 100; % 弃风弃光惩罚系数 penalty 500; Objective sum(a * P_g.^2 b * P_g c) ... penalty * sum(P_w_available - P_w_use) ... penalty * sum(P_pv_available - P_pv_use) ... sum(price_ev .* P_ev); % 功率平衡约束 Constraints [sum(P_g) sum(P_w_use) sum(P_pv_use) sum(P_load) sum(P_ev)]; % 出力上下限 Constraints [Constraints, P_g_min P_g P_g_max]; Constraints [Constraints, P_w_use 0, P_w_use P_w_available]; Constraints [Constraints, P_pv_use 0, P_pv_use P_pv_available]; % EV充电功率约束 Constraints [Constraints, 0 P_ev P_ev_max]; % 调用求解器 ops sdpsettings(solver, cplex, verbose, 0); optimize(Constraints, Objective, ops);留意几个容易出问题的地方。一是目标函数里的二次项。火电成本函数一般写成ax²bxc这是个凸二次函数YALMIP能直接处理但对应的求解器必须是能解二次规划的比如CPLEX、Gurobi的QP模块或者Matlab自带的quadprog。二是功率平衡约束到底按时段写还是按总能量写。调度模型通常要求每个时段都满足功率平衡而不是一天的总发电等于总用电。上面代码里为了简洁写了总功率平衡但实际论文里几乎都是逐时段约束。如果你按总能量约束会得到一个整体可行但每个时段严重失配的结果。三是罚函数系数。弃风弃光惩罚系数如果设得太大求解器会优先保证新能源全额消纳设得太小则可能出现“弃风比干活划算”的荒谬结果。这个系数需要根据论文的定价逻辑来调我一般先按论文原文取跑不通再微调。3.4 结果可视化哪些图必须复现出来出图这一步决定你复现出来的效果符不符合发表水平。标准的几张图包括系统各时段负荷与新能源出力曲线对比图优化前后EV充电负荷的柱状图或曲线图火电出力与新能源出力构成堆积图弃风弃光量的柱状图Matlab绘图的常用套路figure; plot(1:T, P_load, k-o, LineWidth, 1.2); hold on; plot(1:T, P_load P_ev_optimized, r--s, LineWidth, 1.2); plot(1:T, P_load P_ev_uncontrolled, b--^, LineWidth, 1.2); legend(基准负荷, 有序充电后总负荷, 无序充电后总负荷); xlabel(时段/h); ylabel(功率/MW); grid on;把这张图画出来基本上就能清楚看到协同调度和普通无序充电的区别无序充电在晚间出现明显的负荷尖峰有序充电则把充电负荷转移到了新能源出力充足的时段。这个对比图也是论文里最能说明问题的核心图很多审稿人一眼就能看出来你做没做出来。4. 实操过程与结果复现从参数设置到算例分析4.1 算例参数设置的四个关键表复现任何一个调度模型之前先复制原论文的算例数据表。我见过很多失败的复现问题不是出在代码而是出在参数跟论文对不上。常见的待设置参数包括系统基准负荷通常引用某地区典型日负荷数据调整到论文的量纲风电场、光伏电站容量按论文给出的装机容量设置火电机组数量与容量单机模型或多机模型多机要分别给爬坡约束EV数量论文里一般是几百到几千辆决定了可调度容量电池容量家用轿车的60kWh或商用的更大容量充电功率慢充7kW、快充50kW把论文里的表2、表3一个个对应到代码里。如果原文给的是标幺值要先把基准功率算出来再反推有名值。很多同学在这里栽跟头数据量纲错了结果曲线的量级完全对不上。4.2 调试时最折磨人的几个坎调试时最容易卡壳的几个地方我按出现频率排个序第一个是维度匹配错误。YALMIP里sdpvar(1, T)声明的是1行T列的变量跟T×1的常量矩阵做运算时Matlab会隐式扩展但结果可能完全不是你想要的。最好在每个约束加进Constraints之前用size()检查维度。第二个是求解器返回不可行。这说明约束之间本来就矛盾了。排查方法很简单先把EV充电约束放宽把功率平衡约束改成不等式留出余量逐步收紧看哪个约束导致不可行。第三个是收敛精度问题。CPLEX和Gurobi默认精度都比较高但如果模型里存在大数值比如成本系数达到10^6求解器会报数值问题。解决办法是重新标度把功率、费用都除以一个基准值让数值落在10^-2到10^3之间。第四个是随机数的问题。蒙特卡洛模拟每次跑的结果都不一样。复现论文时一定要固定随机种子否则你永远无法得到一个稳定可对比的基线结果。4.3 结果要对比、分析、记录三项一起做模型跑通之后输出要做四件事对比、分析、记录、入图。对比是指把有序充电和无序充电两种场景的结果放在一起看系统性指标的变化系统峰值负荷是否下降了负荷方差是否减小了新能源消纳率是否提高了运行成本是否降低了很多论文在摘要里给出的一组百分比数字就是从这三个指标里来的。分析是要解释“为什么是这个结果”。以削峰填谷为例有序充电的本质是把高峰期充电需求转移到谷时段。看最终曲线时要能指出EV负荷的转移方向是否与新能源出力曲线一致。如果光伏大发时段在中午而你的充电负荷大部分被安排到了凌晨2点那说明目标函数里没有把“促进消纳”考虑进去或者约束条件设置有问题。记录是指把每次实验的配置和结果存下来。我的习惯是运行脚本最后自动把关键指标写成表格存成Excel或CSV方便后面写论文时直接引用数据而不是每次临时截图。场景峰值负荷/MW负荷方差新能源消纳率运行成本/元无序充电15233560091.2%125400协同调度13872180096.8%113200指标改善-8.9%-38.8%6.1%-9.7%5. 常见问题与排查技巧实录5.1 YALMIP与求解器相关的问题YALMIP最常见的坑是路径添加后仍然报“未定义函数或变量”。多数情况下是解压位置的问题把整个文件夹放在当前工作目录就可以但代码里出现相对路径变化后就找不到工具包了。建议把工具包放在固定的绝对路径下在启动脚本里统一addpath。求解器找不到的另一个原因是许可证。CPLEX和Gurobi学术许可证申请的是站点许可证安装后需要配置环境变量。YALMIP调用时如果报“No solver available”用yalmiptest命令检查哪些求解器可用。如果yalmiptest里CPLEX显示为不可用多半是版本匹配问题建议把YALMIP更新到最新版再检查CPLEX路径。5.2 求解速度慢的几个优化思路多时段模型加上整数变量比如EV充电状态的0/1变量之后求解时间可能从几秒暴涨到几分钟。常用处理思路有几个线性化是最有效的手段把二次目标函数分段线性化把非线性约束用大M法转成线性约束把整数规划变成纯LP会快一个数量级。简化时间粒度也能缓解把1小时一个时段改成15分钟一个时段或者反过来在精度和速度之间做取舍。减少场景数同样见效如果不确定性用了多个场景先试着把场景数减到1看能不能跑通再逐步加。最后是调整求解器参数比如设置gap容差0.1%到0.5%的精度损失对工程问题来说完全可接受。5.3 曲线结果不合理时从哪里查起跑出来的曲线颠三倒四先检查数据再检查模型最后检查求解器。数据层面看新能源出力序列有没有负值负荷曲线有没有瞬时跳变EV到达时间有没有超出1到24的范围。模型层面逐条约束检查有没有漏掉“每个时段”的限定。最容易漏的是爬坡约束火电机组在两个相邻时段之间的出力变化不能超过爬坡速率漏了这个约束偶尔会出现相邻时段出力跳变的锯齿形曲线。求解层面看求解器返回的求解状态是Infeasible、Unbounded还是Zero。状态为Zero的常见原因是所有目标函数系数和约束右侧都是一类很小的数数值求解器直接把它当成了平凡问题。5.4 报错速查表报错信息原因解决办法Index exceeds number of array elements数组越界检查循环变量的上限尤其是跨天充电的时段索引Conversion to double from sdpvar is not possible把YALMIP变量当普通数值处理确认对所有变量都先赋值optimize再取value()Solver not found求解器未安装或未配置路径运行yalmiptest查看可用求解器列表NaN in output数据包含NaN或除零检查原始数据确保没有空值避免除零Infeasible problem约束矛盾逐个放宽约束列定位矛盾来源Undefined function wblrnd缺少统计工具箱安装Statistics and Machine Learning Toolbox6. 写在后面几点个人体会到这里整个复现流程就讲完了。最后分享几点我在实际项目里体会最深的经验。其一代码能不能跑通是表层问题模型是不是合理才是核心。复现论文的时候一定要先有自己的判断——这个目标函数合理吗这个约束符合物理规律吗把问题理解透了代码改起来才有方向。其二数据管理和实验记录比代码本身更重要。我因为前期没保存好参数配置后期补算一组对比实验时跑了两天才发现初始场景参数弄错了。建议一开始就建一个参数记录表每个实验跑完顺手更新。其三这个题目还有很多可以扩展的方向加入储能系统、考虑分时电价下的用户响应行为、引入多目标优化算法、用深度强化学习方法替代传统优化……想水一篇出来不难想把这个框架吃透做深够你研究一年半载。就我个人经验把这篇复现真正搞明白的人后面做任何能源系统优化方向的课题都能少踩很多坑。