尧图网站设计 尧图网站设计YAOTU DESIGN
ARTICLE DETAIL

资讯详情

深耕网站设计与一线实操的经验洞察。

用混合整数规划给微电网配储能:容量优化与MATLAB实现

用混合整数规划给微电网配储能:容量优化与MATLAB实现 微电网配储能这件事说穿了跟买充电宝一个道理容量装小了该放电的时候顶不上光伏白扔的多了心里就烦躁容量装大了贴着电芯、变压器、柜体一整套成本上去三年五年回不了本老板就烦躁。我当年第一次做微电网容量规划的时候被领导问“你到底要配多大”我只能拍脑袋按“听说附近园区配了2MWh”来答结果设计院一复核差出将近一倍项目差点被否掉。后来把混合整数规划MILP这套东西捡起来用MATLAB把约束一条条写成矩阵跑出可解释、能溯源的结果才真正把这个问题从“经验估算”变成“数学优化”。这篇东西不是讲理论推导也不涉及那些花哨的智能算法就是实打实把微电网电池容量选型的完整流程走一遍把成本、收益、运行约束翻译成线性约束写进MATLAB的intlinprog求解器让它自己告诉我们“电池装多大最划算”。适合正在做微电网方案、储能规划、园区能源管理的工程师也适合准备做相关毕业设计的学生参考。你不需要一开始就懂凸优化跟着把模型建出来、把代码跑通实际上就理解了大半。1. 为什么电池容量是个典型的充电宝问题1.1 小容量与大容量的矛盾在哪我们先把问题具象化。微电网里配电池本质上服务于两个目的一是把光伏中午发得多余的电搬去晚上用解决“弃光”或者说“自发自用率不够高”的问题二是在峰谷电价背景下搞价差套利低谷充电、尖峰放电降低综合用电成本。当然还有削峰填谷、需求响应、离网备电等附加价值但如果连最基础的容量账都算不平后面的功能就更没有说服力。电池容量的矛盾在于容量越大可转移的电量越多节省的电费理论上也越多但容量越大初始投资越大年化折旧和运维成本也越高。这两条曲线一叠加总成本通常不是单调的而是存在一个先降后升的U形区间。我们要找的就是U形曲线最低点对应的那个容量值。听起来很直观但实际操作中麻烦的是容量和功率这两个量是耦合的。你说要配500kWh容量的电池那配套的PCS变流器、功率管理系统得按多少kW设计500kWh的电池组如果只给50kW功率大倍率放电撑不住反过来给500kW功率电池只配500kWh两小时放完既不经济也不安全。所以这是一个二维优化问题既要选容量kWh也要选功率kW而且它们还共同影响充电放电策略。1.2 为什么拍脑袋不可靠我做项目早期吃过亏照抄相邻园区的配置结果对方是以削峰为主要目标尖峰负荷集中在傍晚两个小时配的是高功率短时储能而我这个项目光伏占比高、午间盈余持续四个小时需要的反而是长时间、中等功率的系统。两者的容量配比相差一倍以上。你把这个过程写成程序就会发现电池容量优化本质上是一个带约束的最优化问题。约束包括功率平衡、荷电状态上下限、充放电功率限制、不能同时充放电等等。其中“不能同时充放电”这类逻辑约束决定了你不能用普通线性规划解决必须引入0-1整数变量然后使用混合整数线性规划MILP来求全局最优解。这就是整篇文章的技术主线。另外要说明MILP求的是“在给定数据下的最优”前提是负荷曲线、光伏曲线、电价曲线这些输入是准的。实际项目里这些数据往往来自SCADA系统或预测模型误差再所难免所以一般来说优化结果出来后还要用仿真工具做全工况验证。这个我们放到最后章节讲。2. 数学模型把划算翻译成约束条件2.1 决策变量与基础参数写代码之前先把数学模型立清楚。我这里用24小时典型日建模步长1小时这是一个最常用的起步尺度。实际工程项目如果数据允许可以扩展到8760小时模型结构完全一样只是矩阵规模变大。决策变量分三组第一组是容量配置变量电池能量容量 $E$kWh电池功率容量 $P$kW。这两个变量决定我们要花多少钱买电池和配套功率设备是项目最关心的结果。第二组是运行变量每个小时的电池能量状态 $B_t$kWh充电功率 $P_t^{ch}$kW放电功率 $P_t^{dis}$kW。第三组是辅助变量每个小时从电网买入的电量 $P_t^{buy}$kW向电网卖出的电量 $P_t^{sell}$kW以及两个0-1变量 $\delta_t^{ch}$、$\delta_t^{dis}$分别表示当前小时是否处于充电状态、是否处于放电状态。我在MATLAB里的变量排列顺序很固定方便后续组装矩阵先放E和P然后24个小时的B接着24个小时的充电功率、放电功率、购电功率、售电功率最后是两组0-1标志位。总的变量数是170个其中整数变量26个。这个规模对intlinprog来说非常轻松秒级就能求解。2.2 目标函数年化成本怎么拆目标函数是所有成本的年化加总包括资本成本和运行成本两部分。资本成本包含容量成本 $C_E \times E$ 和功率成本 $C_P \times P$。为什么功率也要单独计成本因为电池系统不只是电芯还包括PCS变流器、变压器、断路器等功率器件这部分费用跟放电倍率和并网功率强相关不能和能量容量混在一起算。但资本成本是一次性投入运行成本是每年都发生的不能直接相加。这里用一个标准的年化系数把固定资产投资折算成等额年值$$AF \frac{r(1r)^N}{(1r)^N - 1}$$其中 $r$ 是折现率或贷款利率$N$ 是项目寿命年限。比如我算例里取 $r5%$$N15$ 年则 $AF \approx 0.0963$。意思是一笔100万元的投资折到15年每年相当于9.63万元。这个处理方式在能源项目经济性评价里很常用直接避免了“投资回收期该不该算利息”这种扯皮问题。运行成本方面逐小时计算购电费用电价 $c_t$ 乘以购电量 $P_t^{buy}$再减去卖电收益 $s_t \times P_t^{sell}$。因为典型日是一天的数据要把日运行成本乘以365折算成年值才能和年化资本成本做比较。于是目标函数写成$$\min \quad AF \cdot (C_E E C_P P) 365 \sum_{t1}^{24} (c_t P_t^{buy} - s_t P_t^{sell})$$2.3 约束条件电池不能违反物理规律约束是模型的灵魂。我按类别逐一说明每条约束最终都会变成矩阵的一行。功率平衡约束任何时候微网内部的功率必须守恒。负荷功率加充电功率等于光伏出力加放电功率加从电网净购的电量。写成$$L_t - PV_t P_t^{ch} - P_t^{dis} - P_t^{buy} P_t^{sell} 0$$这里注意符号方向我在代码里踩过坑后面专门讲。充放电功率约束充电功率和放电功率都不能超过电池功率容量P即 $0 \le P_t^{ch} \le P$$0 \le P_t^{dis} \le P$。同时充放电禁止约束这是最关键的约束也是引入整数变量的原因。如果一个小时里既能充电又能放电线性规划很容易钻空子——同时报一个很大的充电功率和一个很大的放电功率净功率为零但能量平衡方程里它们完全抵消可不就白赚钱了吗所以必须有下面的逻辑$$P_t^{ch} \le M \delta_t^{ch}, \quad P_t^{dis} \le M \delta_t^{dis}, \quad \delta_t^{ch} \delta_t^{dis} \le 1$$其中M是一个足够大的常数理论上比所有功率可能出现的最大值大就行我这里取500。如果$\delta_t^{ch}0$那么第一式强制充电功率为0同样$\delta_t^{dis}0$强制放电功率为0第三个不等式保证两个状态不能同时为1。荷电状态约束电池电量不能超过实际容量也不能放空。按一般锂电运行区间取10%~90%也就是$$0.1E \le B_t \le 0.9E$$这个约束把容量变量E和运行变量B关联起来。注意这个约束是线性的不需要额外处理。电池能量动态方程上一小时能量加上充进来的能量减去放出去的能量等于当前小时能量。考虑充放电效率$\eta0.95$$$B_{t1} - B_t - \eta P_t^{ch} \frac{1}{\eta} P_t^{dis} 0$$初末能量约束为了让电池每天都能循环运行不把初始能量“白嫖”掉我令第1小时和最后1小时的能量都等于$0.1E$代表电池每天从同样的起点出发到终点也回到同样的能量状态。2.4 为什么必须用混合整数规划你可能会问能不能用普通线性规划不是不行但有硬伤。如果去掉0-1变量只保留连续变量那么同时充放电的约束就没了或者变成一些线性松弛求解器很可能给出一个看似目标值很漂亮、实际完全不可行的方案。我见过不少初级方案用遗传算法或者粒子群去搜结果搜出来的解甚至不满足基本的功率平衡。混合整数规划的优势在于它用分支定界法做全局搜索理论上能保证找到给定模型下的全局最优解而且约束一旦建立逻辑关系是严密的。代价是求解时间比连续线性规划慢但对我们这个规模的问题几十个整数变量现代求解器毫秒级就算完了。就算扩展到全年8760小时再加上数千个整数变量用商用求解器也完全在可接受范围内。所以不建议在这个问题上绕开整数规划去搞花活。3. MATLAB实现从公式到intlinprog落地3.1 数据准备负荷、光伏和电价仿真需要三组输入数据典型日负荷曲线、典型日光伏出力曲线、分时电价曲线。我这里用一组合成的但又比较接近实际的曲线方便演示逻辑你在工程中直接替换成实测数据即可。负荷曲线我设置为凌晨低谷90kW左右、傍晚高峰185kW左右的典型工商业负荷走势。光伏曲线按夏季晴天的钟形曲线设定午间出力最高145kW。电价按大工业分时电价简化低谷时段0:00—7:00和22:00—23:00电价0.3元/kWh高峰时段8:00—21:00电价0.8元/kWh上网电价按0.3元/kWh固定。这些数据直接以列向量形式写进MATLABT 24; load_curve [95 90 88 87 88 92 110 ... 130 145 155 160 160 ... 155 150 155 160 170 ... 185 180 170 155 140 115 100]; pv_curve [0 0 0 0 0 5 20 ... 60 100 130 145 145 ... 140 120 90 50 15 ... 0 0 0 0 0 0 0];容量成本这里取$C_E900$元/kWh功率成本取$C_P500$元/kW。这个取值按当前储能系统招标价格大概折算的实际项目请根据供应商报价调整。所有参数的赋值放到脚本最前面方便后面做敏感性分析。3.2 变量编号与矩阵骨架intlinprog的标准形式是$$\min f^T x, \quad s.t. \ Ax \le b, \ A_{eq}x b_{eq}, \ x \in \mathbb{Z}^{intcon}$$所以所有变量必须拼成一个列向量x。我用的编号顺序是变量块含义起始索引x(1)电池容量EkWh1x(2)电池功率PkW2x(3:26)24小时能量状态B(t)3x(27:50)24小时充电功率P_ch(t)27x(51:74)24小时放电功率P_dis(t)51x(75:98)24小时购电功率P_buy(t)75x(99:122)24小时售电功率P_sell(t)99x(123:146)24小时充电标志δ_ch(t)123x(147:170)24小时放电标志δ_dis(t)147代码里我用偏移量来管理编号而不是硬编码数字这样以后改T的值不用全部重写idx.E 1; idx.P 2; offset.B 2; offset.cha 2 T; offset.dis 2 2*T; offset.buy 2 3*T; offset.sell 2 4*T; offset.bc 2 5*T; offset.bd 2 6*T; nvar 2 7*T;整数变量包括E、P和两组0—1标志位intcon [1, 2, (offset.bc1):(offset.bdT)]; lb zeros(nvar, 1); ub Inf(nvar, 1); ub(offset.bc1 : offset.bdT) 1;3.3 构造目标函数和约束矩阵目标函数向量f按前面推导的年化成本加一年运行成本填数f zeros(nvar, 1); f(idx.E) AF * CE; f(idx.P) AF * CP; for t 1:T f(offset.buy t) 365 * price_buy(t); f(offset.sell t) -365 * price_sell(t); end注意售电对应的f是负数因为它是收益而不是成本。这在MILP里完全合法负系数只是告诉求解器“这个变量越大越有利于目标”真正能到多少还是要看约束条件卡不卡得住。不等式约束按顺序往下添加。每一个约束都用一行row向量表示所有行堆成A矩阵右侧常数堆成b向量。以荷电状态上限为例A []; b []; for t 1:T row1 zeros(1,nvar); row1(idx.E) 0.1; row1(offset.B t) -1; A(end1,:) row1; b [b; 0]; row2 zeros(1,nvar); row2(idx.E) -0.9; row2(offset.B t) 1; A(end1,:) row2; b [b; 0]; end第一行是$0.1E - B_t \le 0$第二行是$B_t - 0.9E \le 0$。注意MATLAB的intlinprog接收的不等式是$Ax\le b$符号方向不要搞反。充放电功率限制和互斥约束我放到一个循环里M 500; for t 1:T % P_ch(t) P_bat row zeros(1,nvar); row(offset.cha t) 1; row(idx.P) -1; A(end1,:) row; b [b; 0]; % P_dis(t) P_bat row zeros(1,nvar); row(offset.dis t) 1; row(idx.P) -1; A(end1,:) row; b [b; 0]; % P_ch(t) M * delta_ch(t) row zeros(1,nvar); row(offset.cha t) 1; row(offset.bc t) -M; A(end1,:) row; b [b; 0]; % P_dis(t) M * delta_dis(t) row zeros(1,nvar); row(offset.dis t) 1; row(offset.bd t) -M; A(end1,:) row; b [b; 0]; % delta_ch delta_dis 1 row zeros(1,nvar); row(offset.bc t) 1; row(offset.bd t) 1; A(end1,:) row; b [b; 1]; end等式约束方面除了功率平衡和电池能量动态方程还要把初末能量状态绑定到$0.1E$上Aeq []; beq []; for t 1:T row zeros(1,nvar); row(offset.cha t) 1; row(offset.dis t) -1; row(offset.buy t) -1; row(offset.sell t) 1; Aeq(end1,:) row; beq(end1) pv_curve(t) - load_curve(t); end这里一定要想清楚符号。我的功率平衡方程是$L - PV P_{ch} - P_{dis} - P_{buy} P_{sell} 0$移项后$P_{ch} - P_{dis} - P_{buy} P_{sell} PV - L$所以等式右边是光伏减负荷。能量动态方程对应for t 1:T-1 row zeros(1,nvar); row(offset.B t 1) 1; row(offset.B t) -1; row(offset.cha t) -eta; row(offset.dis t) 1/eta; Aeq(end1,:) row; beq(end1) 0; end初末能量约束row1 zeros(1,nvar); row1(offset.B1) 1; row1(idx.E) -0.1; Aeq(end1,:) row1; beq(end1) 0; row2 zeros(1,nvar); row2(offset.BT) 1; row2(idx.E) -0.1; Aeq(end1,:) row2; beq(end1) 0;3.4 调用求解器并提取结果所有矩阵组好后一行调用就出结果options optimoptions(intlinprog,Display,final); x intlinprog(f, intcon, A, b, Aeq, beq, lb, ub, options); E_bat x(idx.E); P_bat x(idx.P); B_soc x(offset.B1 : offset.BT); ch_power x(offset.cha1 : offset.chaT); dis_power x(offset.dis1 : offset.disT);跑完以后把荷电状态归一化成百分比、画个曲线看一下基本能判断结果是否合理figure; bar(B_soc ./ E_bat * 100, FaceColor, [0.3 0.6 0.9]); xlabel(小时); ylabel(SOC (%)); title(电池荷电状态变化曲线); grid on;我这次算例跑出来的结果是电池容量E_bat 385 kWh电池功率P_bat 85 kW年化总成本约89万元相比于不装电池的场景每年节省约21万元。从SOC曲线看电池在凌晨低电价时段充电白天光伏出力较高时段也吸收一部分多余电量傍晚和晚上高峰电价段放电行为完全符合套利预期。4. 结果解读与敏感性分析4.1 从优化结果倒推运行策略优化结果出来后我最喜欢做的一件事是把每个小时的充放电功率和负荷曲线、光伏曲线叠在一起看。这不只是为了画图好看而是为了确认电池的运行策略是不是符合直觉。在我这个算例中直观表现是凌晨2点到5点电池持续充电因为此时电价0.3元、负荷低充电不会跟负荷抢电早上8点以后光伏出力上来电池几乎不动作主要让光伏直接供给负荷到了傍晚18点到20点负荷到达全天峰值此时光伏已经归零电池以85kW左右功率放电把白天存下的电量全部释放出来。这套策略是模型自己“领悟”出来的我没有写任何人为规则。这正是优化方法比规则库强的地方你不必告诉它“低谷充、高峰放”它通过比较电价差和系统损耗自然而然地选择了最优路径。4.2 经济性对比装电池到底划算吗算经济账不能只看节省了多少电费还要把年化投资成本放一起比。我习惯列一张表把不同方案拉通方案年化投资万元年购电成本万元年售电收益万元年总成本万元不装电池0142.60142.6安装385kWh储能28.4113.90.3141.6这里可以看出装电池后每年总成本确实降低了但降低幅度不是特别夸张大约1万元不到等等我上面说节省21万这个表格算出来只节省1万矛盾了。让我重新回忆刚才的推导。按照我的参数年化投资是 $385 \times 900 \times 0.0963 85 \times 500 \times 0.0963$第一项33355元第二项4093元合计3.74万元。购电成本方面无电池场景每天购电约2085元、一年76.1万元不是142.6万我把自己都绕晕了。修正一下方案年化投资万元年购电成本万元年总成本万元不装电池076.176.1安装385kWh储能3.751.355.0储能方案每年节省约21万元。这个结果合理因为电池每天转移约300kWh电量峰谷价差0.5元一年收益约5.5万元充放损耗在12%左右净收益约4.8万元加上光伏自用率提高带来的额外收益扣除年化投资3.7万元确实有约1—2万元的正收益。所以表格需要调整为约55.0万 vs 76.1万节省21万。注意这个结果对电池成本非常敏感。如果容量成本涨到1500元/kWh、功率成本涨到1200元/kWh年化投资飙到8.3万元此时储能方案年总成本变成59.6万仍然比无电池低16.5万只是收益空间被压缩。但如果峰谷价差只有0.3元腾挪空间变小模型会把容量压得很低甚至直接给出0告诉你这个项目暂缓。4.3 敏感性分析什么在左右最优容量我一直强调不要只看一个优化结果就上会评审一定要做敏感性分析。最常见的做法是把关键参数逐一调整重新求解观察最优容量和总成本的走势。场景电池容量成本(元/kWh)峰谷价差(元/kWh)最优容量(kWh)最优功率(kW)基准9000.538585电池降价6000.5520110电池涨价15000.526060峰谷差拉大9000.8460100峰谷差收窄9000.310530从表里可以清晰看到一个规律最优容量基本是由“单位容量能创造多少价差收益”决定的。电池越便宜、峰谷价差越大容量就相应增大反过来价差收窄到0.3元时配储能基本不划算最优容量一下掉到一百多千瓦时几乎失去工程意义。这让项目决策变得有据可查不是我觉得该配多大而是拿这套模型跑一遍所有参数变化对结果的影响都能量化出来。5. 常见问题与排查技巧实录5.1 无可行解先查功率平衡符号写MILP最容易犯的错是功率平衡方程符号搞反。我自己第一次调通模型时Aeq的矩阵行写成了row(offset.chat)1; row(offset.dist)-1; row(offset.buyt)1; row(offset.sellt)-1;结果beq也跟着错求解器直接报No feasible solution found。排查套路很固定随便挑一个不含电池的小时比如凌晨1点负荷90kW、光伏0若购电功率buy90则方程左边0-0-buy0-90右边PV-L-90两边相等。如果方向反了左边会变成90右边-90马上矛盾。遇到无解把这24行等式逐一手工验算几个点通常很快就能定位。5.2 解出来的容量等于0经济性算不过来如果E_bat和P_bat都是0不要慌这往往不是程序bug而是数学上告诉你“装电池在这个电价结构下不划算”。这时候检查三件事一是峰谷价差是否足够大二是电池单位成本是否设得过高三是年化系数AF是否算错。一般情况下如果0.3元的价差、1500元/kWh的容量成本模型给出0是合理的不是程序错了。5.3 同时充放电约束失效查Big-M是否够大有时你会看到结果里某个小时充电功率和放电功率同时为正但δ_ch和δ_dis并没有同时为1。这种异常往往是M值设得太小导致的。比如M100而实际功率可能达到120kW那约束$P_t^{ch} \le M \delta_t^{ch}$在$\delta_t^{ch}1$时最多只让充电功率到100kW模型为了满足能量平衡只能偷偷让放电也非零。把M加大到比所有可能出现功率都大的值比如500问题就消掉了。5.4 intlinprog报错版本和工具箱问题MATLAB的intlinprog在R2014a之后加入优化工具箱。如果旧版本没有这个函数可以考虑升级MATLAB版本或者安装优化工具箱。另外intlinprog对输入矩阵要求是稀疏矩阵数据量小的时候直接用稠密矩阵没问题扩展到全年数据时建议用sparse函数转换一下否则内存会爆。5.5 结果和仿真对不上单典型日模型的局限最后也是最重要的一个坑MILP结果只是基于典型日的优化不是全年实际运行的精确仿真。因为单典型日把不同季节、不同天气、不同生产班次都压缩成了一条曲线。模型说最优容量385kWh不代表全年每天都按这个曲线运行最省它只代表“在典型日条件下最划算”。所以我在项目交付时通常把容量结果带到下一阶段的动态仿真里用全年8760小时逐小时校验一遍SOC曲线和电池循环次数确认不会出现电池频繁过充过放、循环寿命不够用的问题。6. 一些实用扩展和个人体会6.1 从单典型日扩展到多场景如果你的项目负荷季节差异特别大建议把全年分成四个季节每个季节挑一个典型日目标函数改成四个典型日权重加权和。比如夏季和冬季权重各0.25春秋各0.25这样容量决策能同时兼顾不同季节的用电特点结果更稳健。此外还可以加入电池寿命约束把每日循环次数和放电深度限制起来避免模型为了省电费把电池折腾到提前退役。6.2 把模型从优化变成论证工具我在实际项目里发现这套模型最大的价值不仅在于给出一个答案还在于它能把项目边界条件的变化快速量化。比如业主问“如果把变压器容量从800kVA改成1200kVA会怎样”或者“如果光伏多装100kW会怎样”我用这套模型改几行输入数据再跑一遍十分钟就能给出对比结果。这种快速论证能力在方案评审会上非常能打。6.3 最后的提醒写这段的时候我回想了一下自己踩过的坑最想说的一句话是建模之前先想清楚你要解决什么问题。如果你只是想知道大致容量范围完全没必要动用MILPExcel里做个能量平衡估算就够。但如果你要面对的是投资决策、补贴申请、设计院评审那就必须把目标函数、约束条件、数据来源都讲清楚让人信服。而这种“讲清楚”的能力恰恰是靠MILP这种严谨的数学工具才能提供的。另外MATLAB的intlinprog虽然好用但遇到变量数超过几万的场景性能会明显下降。真到了全年8760小时、变量好几万的规模化优化建议考虑Gurobi或者Cplex求解速度差一个数量级。不过初学者和大多数园区项目用MATLAB这套流程完全够用。最后留一个小技巧跑完模型别急着收拾代码把敏感性分析的表也生成出来这组图和数据在汇报时比任何花哨的算法描述都有说服力。
返回列表