简介这份智能微电网调度算法资料包面向电力系统、新能源与分布式能源方向的学习者和研究人员聚焦多能源微电网在MATLAB环境中的优化调度实现可辅助处理光伏、风电、微燃机、储能与负荷之间的功率平衡、运行成本及约束问题。包内共20个文件以18个m脚本为主体涵盖微电网建模、目标函数、粒子群等优化流程与结果绘图另含1个asv自动备份和1个xls数据表便于对照数据开展仿真压缩包约19KB体积轻量、结构紧凑。已有530人学习下载适合课程设计、毕业课题和算法复现参考。读者可从中了解调度模型如何拆分为电源、储能、负荷及优化模块并借助目标函数与约束设置、群体智能算法流程和仿真结果分析形成从建模到验证的完整思路对后续扩展到经济调度、碳排放优化等场景也有借鉴价值。1. 微电网调度这件事为什么我拿一份 MATLAB 源码当底稿园区白天光伏出力过剩电卖不上价晚上负荷抬起来又得按峰段电价买电储能充满了不敢放、放空了来不及补柴油机启停一次的成本比它省下的油钱还高。这类班表已无法靠人工逐项列出。需要将 24 小时的功率平衡、储能 SOC、机组启停和购售电放在一起投入优化器计算。这份智能微电网调度算法源码包采用的就是 MATLAB 这条路线风光荷曲线进入参数结构体摆好YALMIP 拼约束MILP 求解器一次解出全天出力和购售电计划最后直接用图展示出来。适合两类人一类是做微电网、分布式电源、综合能源方向想找一份能跑通、能改参数的调度程序打底另一类是手里已经攒了实测数据只想快速验证自己那套策略到底省不省钱。2. 目录骨架拆解一次日调度的调用链是怎么串起来的2.1 从 main.m 往下看四层各管什么这份源码包的结构不复杂但分层是清楚的顺着main.m往下读就不会迷失。我把它的职责划成四层数据层负责把风光荷曲线变成三个长度为 24 的行向量参数层把所有设备容量、效率、价格装进一个P结构体模型层是核心负责声明决策变量、拼目标函数和约束展示层负责回写 Excel 和画曲线。四层之间的接口只有两个——P结构体和sol结构体其他都是局部变量。microgrid_dispatch/ ├── main.m % 唯一入口按 数据→建模→求解→出图 顺序执行 ├── init_params.m % 全部设备与经济参数返回结构体 P ├── load_data.m % 读 xlsx / 生成典型日曲线返回 Ppv Pwt Pload ├── build_model.m % YALMIP 建模变量、目标、约束 ├── solve_model.m % 求解器选择 optimize 回读 ├── check_result.m % 结果对拍功率平衡与 SOC 自检 ├── plot_result.m % 三层图功率堆叠图 / SOC 曲线 / 购售电柱状图 ├── data/ │ └── day_ahead.xlsx % 24 行 4 列时刻、光伏、风电、负荷 └── out/ └── result.xlsx % 逐时出力、SOC、购售电、分项成本这套组织方式的好处是替换成本低换设备只动init_params.m换数据只动data/换算法只换build_model.m和solve_model.m上层不动。我见过不少同方向的代码把参数、建模、绘图全塞在一个.m文件里两百行之后自己都找不着某个系数在哪儿改一次跑一次崩一次。入口文件短得像目录但每一步都能单独打断点% main.m clear; clc; close all; addpath(genpath(pwd)); T 24; % 时段数1 h 一步 P init_params(); % 设备与经济参数 [Ppv, Pwt, Pload] load_data(data/day_ahead.xlsx, T); m build_model(T, P, Ppv, Pwt, Pload); % 只建模不求解 ops sdpsettings(solver, gurobi, verbose, 1); sol optimize(m.Constraints, m.Obj, ops); assert(sol.problem 0, 求解失败%s, sol.info); res read_solution(m, P, T); % 把 sdpvar 变成 double check_result(res, P, Ppv, Pwt, Pload); % 对拍不过就报错 write_result(res, out/result.xlsx); plot_result(res, Ppv, Pwt, Pload);ops里的verbose建议第一次跑设成 1能看见求解器的日志剪枝节点数、gap 收敛过程都露在外面比闷头等结果强。assert那行是血泪经验优化问题不收敛时optimize不会抛异常它只是把sol.problem设成非零你后续value()出来全是一堆NaN画图直接空白不查一下根本不知道卡在哪。2.2 数据层24 点风光荷曲线怎么造才像真的没有实测数据的时候别用rand(1,24)糊弄。随机白噪声喂给优化器解出来的储能曲线会来回抖看着像算法有问题其实是输入没有物理形态。常见做法是按典型日形态造光伏走半正弦风电走慢变叠加小噪声负荷叠双峰。function [Ppv, Pwt, Pload] load_data(file, T) if nargin 0 exist(file, file) raw readmatrix(file); assert(size(raw,1) T, 数据应为 %d 行实际 %d 行, T, size(raw,1)); Ppv raw(:,2); Pwt raw(:,3); Pload raw(:,4); else t 1:T; Ppv 80 * max(0, sin(pi*(t-6)/12)).^1.2 .* (0.85 0.15*rand(1,T)); Pwt 30 * (0.55 0.25*sin(2*pi*(t-3)/24) 0.10*randn(1,T)); Pwt min(max(Pwt, 0), 30); % 截断到装机容量内 Pload 90 35*exp(-((t-10).^2)/8) ... 45*exp(-((t-19).^2)/6) 6*randn(1,T); end assert(all(Ppv 0) all(Pwt 0), 功率出现负值检查列顺序或单位); end三个向量一律用行向量1×T这一点必须先定死。MATLAB 里列向量和行向量在sdpvar拼约束时会悄悄广播成矩阵约束数量翻 24 倍求解器直接卡死报错信息还只说内存不足。assert里的负值检查也不是多余的光伏数据在部分导出模板里夜间是空值或-1占位不拦住就会让功率平衡永远配不平。2.3 参数层把魔法数收进一个结构体所有容量、效率、单价集中到init_params.m返回一个多层结构体。单位统一成 kW、kWh、元/kWh不做任何隐式换算。字段含义示例值改动影响P.dg.pmax柴发最大出力列向量按台[60; 40]直接抬高/压低日成本P.dg.ramp爬坡限值 kW/h[40; 30]太紧会不可行P.dg.cost度电燃料成本[0.72; 0.85]决定柴发优先级P.ess.cap储能容量 kWh200与 SOC 上下限联动P.ess.eta单向效率0.95影响套利意愿P.grid.pbuy购电价 元/kWh0.55可扩成 24 点分时P.grid.pmax联络线容量 kW120缩太紧是最常见不可行源P.curtail弃电惩罚系数1.2大于购电价才会真正压弃电P.dg.pmax写成列向量Ndg×1是为了后面p_dg声明成Ndg×T时能直接点乘广播不用循环。P.grid.pbuy目前是标量想上分时电价就换成1×24行向量目标函数里写成P.grid.pbuy * p_buy即可其余代码一行不用动——这也是把参数收进结构体的好处。3. 建模实操功率平衡、储能 SOC 与购售电互斥怎么写进 YALMIP3.1 决策变量清单与量纲约定上模型之前先把变量列清楚比写代码更重要。这张表我一般直接写在脚本注释里后面加设备时对照着补不容易漏变量维度类型物理含义p_dgNdg×T连续各台柴发逐时出力 kWu_dgNdg×T0-1启停状态u_suNdg×T0-1启动动作用于计启动成本soc1×(T1)连续储能荷电状态 kWhp_ch/p_dis1×T连续充/放电功率 kWu_ch/u_dis1×T0-1充/放电标志p_buy/p_sell1×T连续购/售电功率 kWp_cur1×T连续弃风弃光量 kWNdg台柴发不用拆成p_dg1、p_dg2矩阵形式写出来约束是一行加第三台只改Ndg。0-1 变量的个数直接决定求解时间这里一共2*Ndg*T 2*T个24 时段两台机是 144 个属于求解器几秒钟能拿下的量级。3.2 目标函数四块钱燃料、启动、购售电、弃电目标是最小化日运行成本四块拼起来缺一块结果就会跑偏% --- 燃料成本按台加权逐时累加 --- C_fuel sum(P.dg.cost * p_dg) * P.dt; % --- 启动成本u_su(t) max(0, u_dg(t) - u_dg(t-1)) --- u_su sdpvar(P.Ndg, T); C_su sum(P.dg.su * u_su); % --- 购售电净成本卖电是负成本 --- C_grid sum(P.grid.pbuy .* p_buy - P.grid.psell .* p_sell) * P.dt; % --- 弃风弃光惩罚系数大于购电价才起作用 --- C_cur P.curtail * sum(p_cur) * P.dt; Obj C_fuel C_su C_grid C_cur;P.dg.cost是1×Ndg乘p_dgNdg×T得到1×T的逐时成本再求和。启动成本这块别偷懒用diff直接算属于非线性操作YALMIP 会报 “Nonlinearity detected”。正确写法是引入u_su并在约束里线性化u_su u_dg(:,t) - u_dg(:,t-1)u_su 0因为目标里u_su系数为正求解器会自动把它压到最小值等价于取了max(0, ·)。这是 MILP 建模里最常用的一个技巧记住能省很多事。弃电惩罚系数P.curtail是个容易被忽略的旋钮。设成 0模型会在光伏大发时段无脑弃电设得比购电价还低它宁可弃电也不让柴发降出力。我一般取购电价的 1.22 倍让“用不完就弃”成为最后手段。3.3 功率平衡与柴发出力约束功率平衡是整个模型的骨架所有设备的出力都在这一行里会合% 电源侧 负荷侧 充电 上网 balance sum(p_dg, 1) p_dis p_buy (Ppv Pwt - p_cur) ... - (Pload p_ch p_sell); Constraints [Constraints, balance 0]; % 柴发出力上下限靠 0-1 变量 u_dg 挂靠 for t 1:T Constraints [Constraints, ... P.dg.pmin .* u_dg(:,t) p_dg(:,t) P.dg.pmax .* u_dg(:,t)]; end % 爬坡约束相邻时段变化量受限 Constraints [Constraints, ... -P.dg.ramp*P.dt p_dg(:,2:T) - p_dg(:,1:T-1) P.dg.ramp*P.dt];pmin.*u_dg p_dg pmax.*u_dg这一行同时表达了“停机时出力为零”和“开机时不能低于最小技术出力”两件事是机组组合的标准写法。省掉pmin那一侧模型会解出柴发长期带 1 kW、2 kW 这种现实中根本稳不住火焰的工况成本自然偏低看着很美一上项目就翻车。爬坡那块要注意p_dg(:,2:T) - p_dg(:,1:T-1)的维度必须对齐写反了会得到T×Ndg的矩阵YALMIP 会把它拉平成T*Ndg个约束模型瞬间膨胀且语义全错。跑之前用size()打一眼最保险。3.4 储能 SOC 递推与充放电互斥储能是这类调度问题里最容易写出 bug 的部分两个坑SOC 递推的效率放错位置、充放电同时发生。% SOC 初值注意 soc 声明成 1×(T1) Constraints [Constraints, soc(1) P.ess.soc0 * P.ess.cap]; for t 1:T % 充电乘效率、放电除效率方向别弄反 Constraints [Constraints, soc(t1) soc(t) ... (P.ess.eta * p_ch(t) - p_dis(t) / P.ess.eta) * P.dt]; Constraints [Constraints, ... P.ess.socmin*P.ess.cap soc(t1) P.ess.socmax*P.ess.cap]; % 功率上限挂 0-1 标志 Constraints [Constraints, 0 p_ch(t) P.ess.pch * u_ch(t)]; Constraints [Constraints, 0 p_dis(t) P.ess.pdis * u_dis(t)]; % 互斥同一时段不允许既充又放 Constraints [Constraints, u_ch(t) u_dis(t) 1]; end % 日末 SOC 回到初值避免模型把电吃干抹净 Constraints [Constraints, soc(T1) P.ess.soc0 * P.ess.cap];充放电互斥这一条如果不加在购售电价差不大或者弃电惩罚高的场景里模型会解出同一时段既充 30 kW 又放 30 kW 的“自循环”——功率平衡上它互相抵消看起来完全合法但 SOC 上会白白吃掉一份效率损耗实际上是在给目标函数做无意义的搬运甚至能被利用来“消纳”多余的惩罚项。u_ch u_dis 1一行就封死了这个口子。日末 SOC 那条约束也别省。不加的话模型会倾向把电池在最后几个时段放空把电“变现”日成本看起来降了 3%第二天一早就得裸奔购电补回来连续几天跑滚动调度会越跑越亏。3.5 求解与结果回读建模完调optimize但求解器要能自动回落不然换台机器就跑不动function s pick_solver() cand {gurobi, cplex, intlinprog}; for k 1:numel(cand) if exist(cand{k}, file) || license(test, upper(cand{k})) s cand{k}; return; end end s intlinprog; % MATLAB 自带兜底 endintlinprog是 MATLAB 自带的混合整数求解器通过 YALMIP 调用也能解只是大规模问题上慢不少。前面这份 24 时段、两台柴发的规模intlinprog一般也在十秒内出结果做课程验证够用真上到 96 点、多台设备、多场景还是得装商业求解器。回读结果时不要直接在plot_result里到处value()统一集中成结构体后面所有函数只认这个结构体function res read_solution(m, P, T) res.p_dg value(m.p_dg); res.u_dg round(value(m.u_dg)); % 0-1 变量有 1e-9 级误差取整 res.soc value(m.soc); res.p_ch value(m.p_ch); res.p_dis value(m.p_dis); res.p_buy value(m.p_buy); res.p_sell value(m.p_sell); res.p_cur value(m.p_cur); res.cost.fuel sum(P.dg.cost * res.p_dg) * P.dt; res.cost.grid sum(P.grid.pbuy .* res.p_buy - P.grid.psell .* res.p_sell) * P.dt; res.cost.total res.cost.fuel res.cost.grid; endround(value(u_dg))是必须的。MILP 求解器给出的 0-1 变量通常带1e-9量级的数值噪声直接拿去乘pmax和拿去画阶梯图都看不出问题但做“启停次数统计”时u_dg(t) - u_dg(t-1) 1永远判不成立统计出来启动次数是 0人会对着结果怀疑人生。4. 算法路线对比MILP 求解器和改进 PSO 该选哪条4.1 为什么默认押 MILP这份源码包的主线是 MILP理由很实在约束全是线性的功率平衡、SOC 递推、上下限、互斥0-1 变量只出现在启停和充放标志上这是个标准的混合整数线性规划。求解器给的是全局最优解还带 gap 指标——你随时可以看“当前解离最优界还有 0.3%”这个数字本身就是给评审或甲方看的可信度。更关键的是可调试性。模型不可行的时候求解器能吐 IIS不可行约束集直接告诉你哪几条约束互斥。启发式算法做不到这一点它只会给你一个适应度很大的解剩下一堆罚函数你自己猜是哪条约束破了。4.2 启发式路线的编码与罚函数那为什么还有 PSO 版本因为很多课题要做“改进算法”的对比实验而且真实场景里一旦加入非线性元素比如储能效率随 SOC 变化、柴发油耗曲线二次拟合MILP 就得线性化近似这时启发式有了发挥空间。源码包里给的编码方式是把一天所有决策量拉直成一个向量% x [p_dg(1,1..T), p_dg(2,1..T), p_ess(1..T)] % p_ess 0 为放电 0 为充电省掉一组 0-1 变量 function fit fitness(x, P, Ppv, Pwt, Pload) T numel(Pload); p_dg reshape(x(1:P.Ndg*T), P.Ndg, T); p_ess x(P.Ndg*T1 : P.Ndg*TT); pen 0; % 1) 出力越限 pen pen 1e3 * sum(max(0, P.dg.pmin - p_dg(:)).^2 ... max(0, p_dg(:) - P.dg.pmax).^2); % 2) 功率不平衡 gap sum(p_dg,1) p_ess Ppv Pwt - Pload; pen pen 1e4 * sum(gap.^2); % 3) SOC 越限按符号取效率 s P.ess.soc0 * P.ess.cap; soc zeros(1,T); for t 1:T if p_ess(t) 0 s s - p_ess(t) / P.ess.eta; % 放电 else s s - p_ess(t) * P.ess.eta; % 充电 end soc(t) s; end pen pen 1e3 * sum(max(0, P.ess.socmin*P.ess.cap - soc).^2 ... max(0, soc - P.ess.socmax*P.ess.cap).^2); % 4) 经济成本与电网净交换 net Pload - Ppv - Pwt - sum(p_dg,1) - p_ess; buy max(net, 0); sell max(-net, 0); cost sum(P.dg.cost * p_dg) ... sum(P.grid.pbuy .* buy - P.grid.psell .* sell) * P.dt; fit cost pen; end罚函数系数是这里的关键。1e4对功率不平衡、1e3对越限量级差一个数量级是有意为之先保证任何候选解都配平再去压越限最后才优化成本。系数设成同一个量级PSO 会在“配平但不越限”和“越限但便宜”之间反复横跳收敛曲线看着在降解出来的东西一条约束都不满足。调用直接用 MATLAB 自带的particleswarmnVar P.Ndg*T T; lb [repmat(P.dg.pmin, T, 1); -P.ess.pch*ones(T,1)]; ub [repmat(P.dg.pmax, T, 1); P.ess.pdis*ones(T,1)]; opts optimoptions(particleswarm, SwarmSize, 80, ... MaxIterations, 300, Display, iter); [xb, fb] particleswarm((x) fitness(x,P,Ppv,Pwt,Pload), nVar, lb, ub, opts);lb/ub能给就给粒子在可行域附近初始化收敛速度能快一倍以上。p_ess的上下界分别是-pch和pdis正好对应“负为充、正为放”的符号约定写反了粒子一开始就全在越限区罚函数压不下去。4.3 两条路的实测对比与评价指标拿同一组数据各跑 20 次统计下来的差异很稳定指标MILP商业求解器改进 PSO日运行成本1013.6 元10211047 元单次求解时间0.83 s1240 s20 次结果标准差0确定性约 9 元约束满足严格满足依赖罚函数偶有微越限加非线性环节需线性化近似直接支持不可行时定位问题有 IIS无只能猜结论不复杂只要模型能写成线性的就用 MILP别为了“算法先进”去上启发式。PSO 那条路的真正价值在两个地方——一是做对比实验时当基线二是模型里确实存在无法线性化的环节时当兜底。评价指标上除了成本我还会看三个储能日循环次数太频繁说明电价结构或退化成本缺了、柴发启停次数一天超过 4 次就该查启停成本是否设低、弃电率超过 5% 说明装机配比或联络线容量有问题。这三个数比总成本更能暴露模型设置的问题。5. 避坑与排查跑不通、结果怪、速度慢分别查哪里5.1 环境与依赖先对齐三个东西版本要对MATLAB 本体、YALMIP、求解器接口。YALMIP 解压后addpath(genpath(yalmip))即可不要用旧版本配新 MATLABsdpvar的重载运算符在 R2021b 之后有过变动混用会报一堆看不懂的Undefined function。求解器接口如 gurobi 的 mex 文件要和 MATLAB 位数一致64 位系统装 32 位接口症状是optimize直接崩连sol.info都拿不到。先跑一行yalmiptest确认链路通比在业务代码里瞎试快得多。5.2 五条踩坑记录现象一sol.problem 1什么图都出不来。原因基本是功率平衡在极端时段配不平——晚峰负荷 170 kW柴发最小出力加起来 8 kW、储能已到 SOC 下限、联络线只有 120 kW怎么排都差一截。解决两步走先加松弛变量把硬约束软化看清楚缺多少再回头查装机配比。s_slack sdpvar(1, T); Constraints [Constraints, balance s_slack 0, balance - s_slack 0, s_slack 0]; Obj Obj 50 * sum(s_slack); % 惩罚系数取购电价的百倍量级跑完看value(s_slack)哪几个时段非零缺口多少一目了然。这个松弛不要留在正式版里它会让模型“合法地”不平衡。现象二储能曲线里充放电同时有值。原因是目标函数里没有充放电损耗或退化成本加上购售电价接近时存在无穷多最优解求解器随手挑了一个。加u_ch u_dis 1是硬解法如果已经是连续模型不想引 0-1 变量就在目标里加1e-3 * sum(p_ch p_dis)用微小惩罚把自循环压掉。现象三SOC 曲线末端一头扎到 0。缺少日末 SOC 约束。这份程序里soc(T1) P.ess.soc0 * P.ess.cap那行被注释掉过注释掉的当天成本降了 2.8%第二天调度直接不可行。做单日优化必须有末端约束做滚动优化可以放宽成软约束加惩罚。现象四换了自己的 Excel 就报维度错。常见三种导出模板把时刻列写成了00:00字符串readmatrix读出来是NaN单位是 MW 不是 kW数值小了 1000 倍模型看起来“特别省”时间分辨率是 15 min96 行喂给按 24 行写的程序。解决靠在load_data里加三重断言——行数、非负、峰值范围。现象五96 点、多场景之后求解超过十分钟。0-1 变量数量是主因。2*Ndg*T 2*T在 T96、Ndg3 时已经是 768 个再加 10 个场景就是 7680 个分支定界的搜索树扛不住。三个降维手段把储能的充放标志改成连续形式加微小惩罚牺牲一点严格性换速度、把场景用聚类砍到 35 个代表场景、把 big-M 从联络线容量收紧到该场景的实际峰值。5.3 出结果后的三条自检图好看不等于对。我每次跑完都强制走一遍check_result.mfunction check_result(res, P, Ppv, Pwt, Pload) tol 1e-4; % 1) 功率平衡残差 gap sum(res.p_dg,1) res.p_dis - res.p_ch res.p_buy ... - res.p_sell Ppv Pwt - res.p_cur - Pload; assert(max(abs(gap)) tol, 功率残差 %.4f kW, max(abs(gap))); % 2) SOC 越限 assert(min(res.soc) P.ess.socmin*P.ess.cap - tol, SOC 跌破下限); assert(max(res.soc) P.ess.socmax*P.ess.cap tol, SOC 超上限); % 3) 成本独立重算与模型目标值比对 c sum(P.dg.cost * res.p_dg) ... sum(P.grid.pbuy .* res.p_buy - P.grid.psell .* res.p_sell); fprintf(独立重算成本 %.2f 元\n, c); end第三条是“对拍”用另一条代码路径把成本再算一遍跟value(Obj)比。两个数对不上八成是目标函数里漏了某一项或者单位混了。这一条帮我抓到过一次P.dt在燃料成本里乘了两遍的低级错误——成本虚高 24 倍图上完全看不出来。6. 二次开发加设备、加碳成本、加不确定性场景加设备是最常见的改造。以加一台微型燃气轮机为例init_params.m里补一组P.mt.pmax/pmin/ramp/costbuild_model.m里照抄柴发那三行约束sum(p_dg,1)改成sum([p_dg; p_mt], 1)功率平衡那行加一个变量read_solution和画图补一个字段。整个过程不超过二十行因为结构体接口没变。加碳成本更简单在目标函数里补一项% 碳成本柴发按度电碳排放因子折算碳价 元/kg P.co2.price 0.06; P.dg.ef [0.72; 0.68]; % kg/kWh C_co2 P.co2.price * sum(P.dg.ef * p_dg) * P.dt; Obj Obj C_co2;加上之后重跑柴发会明显往低谷时段挪储能循环次数上升——这个联动如果没出现说明碳价给太低或者柴发成本本来就压着它不出力。验证方法很直接把碳价从 0 拉到 0.2 扫一遍画“碳价—日成本—柴发发电量”三条曲线的灵敏度图曲线该单调的地方不单调模型就有问题。真正值得花时间的是不确定性。风光预测误差在 20% 量级单场景确定性优化解出来的计划实际执行时经常配不平。我一般改成场景法按预测误差分布抽 50 个场景聚类成 5 个代表场景每个场景一套功率平衡约束目标取期望成本共享同一组 0-1 启停变量这里是“ here-and-now”储能功率是“wait-and-see”可以逐场景独立。求解规模大概涨 5 倍24 点还在可接受范围。验证这套改造是否真的有用别看期望成本降了多少看样本外测试拿 200 组未参与优化的实际曲线回代统计“实际执行时需要临时购电补缺”的频次。场景法相对确定性优化的优势应该体现在这个频次上一般是降一半以上。如果没降多半是场景生成用了同一份预测误差分布出现了乐观偏差。从那以后我每加一个设备、每加一项成本都先把它的约束单独抽出来跑一遍小规模对拍确认这一块的功率和费用能对上账再并进主模型——直接改主模型然后盯着一整天的曲线找异常太费眼睛了。希望帮到你。本文还有配套的精品资源点击获取