做微网虚拟电厂调度优化这个方向时间不短了最近被问得最多的一个问题就是想在做日前优化调度的时候同时把碳排放交易和多种需求响应考虑进去但网上代码要么只做碳交易要么只做需求响应两者结合的模型又很复杂不知道怎么在MATLAB里落地。这篇文章就把我自己实际用MATLAB搭建的“计及碳排放交易及多种需求响应”的微网/虚拟电厂日前优化调度模型完整拆开讲一遍包含模型设计、数学表达、YALMIP代码实现、仿真算例和踩坑记录。适合目标为园区微网、虚拟电厂聚合商的调度策略研究或者正在做相关方向毕业设计的同学参考。先说结论把碳交易和需求响应同时纳入优化不等于简单加两组约束。它们之间会互相影响——需求响应改变了负荷曲线进而改变了购电量和机组出力最终影响碳排放量而碳成本又反过来改变了负荷调整的经济收益所以必须放到同一个目标函数里联合求解。这个看起来“多个模块叠加”的问题实际做起来有很多细节。1. 碳交易机制进模型不只靠一个系数1.1 为什么传统经济调度模型在碳约束下会失真传统微网日前调度目标函数通常是“购电成本燃料成本运维成本”最多加上储能折旧约束就是功率平衡、机组上下限、储能SOC。这样的模型在碳相关考核不严格的情况下还能用但一旦微网或虚拟电厂需要参与碳市场交易问题就来了系统如果不考虑碳排放配额和交易成本调度结果可能是“电价低时多网购电”或“让燃气轮机满发”但这些策略可能让实际碳排放量超过免费配额需要额外购买碳配额反而让综合成本上升。所以碳交易本质上不是“给目标函数乘以一个系数”而是要把碳排放配额、实际排放量、配额买卖量之间的差额逻辑建模成一个可交易变量这个变量直接进入目标函数形成“排放-交易-成本/收益”闭环。1.2 免费配额、排放因子与碳价怎么进入优化变量主流碳交易机制中企业会获得免费配额比如基准线法或历史强度下降法。实际运行期结束统计实际碳排放量 (E_{total})与配额 (E_{quota}) 做差如果 (E_{total} \le E_{quota})多余的配额可以出售获得收益如果 (E_{total} E_{quota})需要购买配额产生成本。这个差额就叫净交易量 (E_{trade})。碳成本就是 (C_{CO2} \lambda \cdot E_{trade})(\lambda) 是碳价。这里有一个很多新手会踩的坑不能直接把 (E_{trade}) 写成绝对值 (|E_{total} - E_{quota}|)因为绝对值函数是非线性的而且不光滑导致优化求解困难。正确做法是把 (E_{trade}) 拆成正负两个方向即碳配额买入量和卖出量(E_{buy} \ge 0)表示购买量(E_{sell} \ge 0)表示出售量(E_{total} - E_{quota} E_{buy} - E_{sell})为了不让模型同时买入又卖出可以加一个约束 (E_{buy} \cdot E_{sell} 0)但这是非线性的。好在最优解中只要碳价为正且模型合理不会出现同时买卖所以这个乘积约束可以省略或者用 Big-M 方法处理。实际排放量 (E_{total}) 的构成也要拆清楚。微网/虚拟电厂通常包含燃气轮机/燃气锅炉燃烧天然气产生的直接排放从电网购电对应的间接排放购电量 × 电网平均排放因子光伏、风电上网或自用可以视为零排放或抵扣。于是 [ E_{total} \sum_{t1}^{24} \left( E_{GT} \cdot P_{GT}(t) \cdot \Delta t E_{grid} \cdot P_{buy}(t) \cdot \Delta t \right) ]其中 (E_{GT}) 是燃气轮机单位发电量的碳排放因子(E_{grid}) 是电网购电等效排放因子。这部分计算在MATLAB中只需要一行代码但它决定了碳交易量必须谨慎标定。1.3 阶梯碳价怎么线性化政策上为了强化减排碳价可能不是固定的而是按配额缺口分档比如0~1000kg配额外部分按0.2元/kg超过1000kg部分按0.35元/kg。这种“阶梯碳价”本质是分段线性函数直接用 if 判断带入优化就会把问题变成非凸。建议用分段线性化的常规做法YALMIP 里有pwf或者binvar加上 Big-M 来实现。但如果你只是做前期的调度策略研究尤其论文方向不是碳市场机制本身把碳价设为常数是最高效的做法。我一般先在常数碳价的模型下跑通整体逻辑再考虑阶梯碳价避免一开始就把模型搞复杂到无法调试。2. 日前调度模型到底要优化哪些变量和约束2.1 目标函数五类成本加一个收益我把目标函数定义为一天的总运营成本最小[ \min \quad C_{buy} - I_{sell} C_{fuel} C_{om} C_{dr} C_{co2} - I_{quota_sell} ]逐项拆解(C_{buy})从上级电网购电费用(\sum P_{buy}(t) \cdot price_{buy}(t) \cdot \Delta t)(I_{sell})向电网售电收入(\sum P_{sell}(t) \cdot price_{sell}(t) \cdot \Delta t)(C_{fuel})燃气轮机燃料成本常见线性近似为 (a \cdot P_{GT}(t) b)也可用二次成本再用分段线性近似(C_{om})设备运行维护成本包括燃气轮机、储能充放电按电量折算、光伏/风机很小可忽略(C_{dr})需求响应实施成本包括可中断负荷削减补偿和可转移负荷调整补偿(C_{co2})净碳交易成本如上一节所述(I_{quota_sell})其实已经包含在 (C_{co2}) 的负值部分中如果 (E_{sell}0)碳收益为负成本。注意Δt取1小时24个调度时段。有些文献喜欢用15分钟间隔不是不行但碳交易结算通常以日/小时为单位1小时间隔已经足够反映日前调度的精度而且能控制模型规模便于调试。2.2 决策变量拆解越清晰越少出错在写代码之前我习惯先列一张决策变量表。这是一开始最容易被忽略的一步很多人直接打开MATLAB开始写sdpvar写到一半就乱。设备/机制决策变量类型维度燃气轮机(P_{GT}(t))连续1×24储能充电功率(P_{ch}(t))连续1×24储能放电功率(P_{dis}(t))连续1×24储能充电状态(u_{ch}(t))0-11×24储能放电状态(u_{dis}(t))0-11×24电池SOC(SOC(t))连续1×24购电功率(P_{buy}(t))连续1×24售电功率(P_{sell}(t))连续1×24可中断负荷削减量(P_{cut}(t))连续1×24可转移负荷净转移量(P_{shift}(t))连续1×24碳配额买入量(E_{buy})连续标量碳配额卖出量(E_{sell})连续标量这里储能状态变量用二进制是为了防止同时充放电。如果用P_ch(t) * P_dis(t) 0这种约束又是非线性。用两个二进制变量加一个上限约束问题就变成MILP可直接用Gurobi/CPLEX求解。2.3 约束条件物理边界和逻辑边界都不能缺约束部分我按设备分类写这样可以快速定位问题。功率平衡约束是所有调度模型的核心 [ P_{buy}(t) P_{pv}(t) P_{wt}(t) P_{GT}(t) P_{dis}(t) P_{load}(t) P_{ch}(t) P_{sell}(t) - P_{cut}(t) - P_{shift}(t) ]注意可中断削减和转移负荷会让等效负荷变小。但可转移负荷的净转移量可能为负表示把负荷转移到该时段所以要加“负荷增加”同样纳入等式。燃气轮机约束出力上下限(P_{GT}^{min} \le P_{GT}(t) \le P_{GT}^{max})爬坡约束(-\Delta P_{GT}^{down} \le P_{GT}(t) - P_{GT}(t-1) \le \Delta P_{GT}^{up})。储能约束充电和放电功率与状态变量耦合(0 \le P_{ch}(t) \le u_{ch}(t) \cdot P_{ch}^{max})(0 \le P_{dis}(t) \le u_{dis}(t) \cdot P_{dis}^{max})同时只允许一种状态(u_{ch}(t) u_{dis}(t) \le 1)SOC迭代(SOC(t) SOC(t-1) \left( \eta_{ch}P_{ch}(t) - \frac{P_{dis}(t)}{\eta_{dis}} \right)\Delta t / C_{bat})SOC上下限(SOC_{min} \le SOC(t) \le SOC_{max})调度周期始末SOC一致比如 (SOC(24) SOC(0))否则模型会为了省钱把电池放空。购售电约束(0 \le P_{buy}(t) \le P_{buy}^{max}\cdot u_{grid}(t))(0 \le P_{sell}(t) \le P_{sell}^{max}\cdot (1-u_{grid}(t)))防止同时购电售电。 也可以不加这个二进制取决于分时电价结构但加上更保险。需求响应约束可中断削减量(0 \le P_{cut}(t) \le P_{cut}^{max}(t))可转移负荷全天转移量守恒(\sum_{t1}^{24} P_{shift}(t) 0)转移量上下限(-P_{shift}^{max}(t) \le P_{shift}(t) \le P_{shift}^{max}(t))。如果进一步建模“可平移负荷块”那就需要整数变量表示开始时间这里采用聚合净转移模型更简单且适合日前调度分析。3. 多种需求响应并存怎么建模最省力3.1 价格型需求响应把弹性矩阵变成输入参数价格型需求响应(Price-Based DR)的逻辑是用户根据电价调整用电量。学术上常用价格弹性矩阵来模拟第 (i) 个时段的负荷变化量受全时段电价影响。[ P_{DR}(t) P_{base}(t) \left( 1 \sum_{j} E_{ij} \frac{\pi(j) - \pi_{ref}(j)}{\pi_{ref}(j)} \right) ]其中 (E_{ij}) 是自弹性与交叉弹性矩阵。自弹性通常为负表示当前时段电价上涨当前负荷下降交叉弹性为正表示其他时段电价上涨用户会把用电转移到本时段。在日前调度中如果电价是已知的分时电价那么价格型DR后的负荷曲线是可以直接预先算出来的不需要作为决策变量。这种处理方式最简单也不增加模型复杂度。如果电价本身作为决策变量参与优化那就变成了双层优化或MPEC问题复杂度高非常多不建议入门时尝试。我在MATLAB里的处理方式很简单先定义24×24的弹性矩阵然后做一个循环算出24个时段的负荷变化率得到修正后的负荷P_load_price。再用这个负荷替代原始基线负荷作为后续优化的输入。有基础的同学可以改成矩阵运算但循环在24维上性能差异可以忽略。3.2 激励型需求响应可中断和可转移分开建模激励型需求响应(Incentive-Based DR)的核心是“签订合同按指令削减/转移负荷并给予补偿”。可中断负荷建模最简单决策变量P_cut(t)上限按合同约定补偿单价price_cut元/kWh目标函数中加入sum(price_cut * P_cut)。可转移负荷建模稍微绕一点但有个非常实用的简化方案定义净转移量P_shift(t)正值表示其他时段负荷转移到本时段本时段负荷增加负值表示本时段负荷转移到其他时段本时段负荷下降。约束 [ \sum_{t1}^{24} P_{shift}(t) 0 ] 即总用电量不变只改变用电时间。实际上这个约束已经可以保证转移量的“守恒性”。同时每个时段的转移量有上下限表示可调灵活性。用这种聚合模型不需要引入大规模0-1变量整个优化仍然是MILP级别求解速度非常快。如果你的场景必须精确到“某台机器只能从8点推迟到10点”那才需要引入离散变量但日常规模化应用用聚合模型足够了。3.3 需求响应和碳交易耦合的“隐藏收益”这一点是很多人忽略的。需求响应让负荷从高峰移到低谷最直接的影响是减少峰时高价购电增加谷时低价购电。但是如果峰时段原本对应的是电网高排放因子比如火电调峰或谷时段光伏大发那么需求响应还会显著降低购电对应的间接碳排放进一步减少碳配额购买这就形成了叠加收益。我在算例中经常看到一种现象单看需求响应补偿费用可能比它带来的电费节省还高但加上碳成本节省后综合经济性就变得非常有优势所以联合优化比分开优化更容易得到“需求响应投入是划算的”的结论。建议在做结果分析时把碳收益单独拆出来算以便判断“需求响应是被电价驱动还是被碳价驱动”这对实际业务决策很有价值。4. MATLAB代码架构与核心实现细节4.1 为什么用YALMIP而不是手写求解器接口手写线性规划矩阵并调用Gurobi/CPLEX也不是不行但模型一旦涉及几十个变量和多类约束矩阵下标非常容易出错。YALMIP作为建模工具箱允许直接用sdpvar、binvar建模然后用optimize调用自带求解器或商用求解器极大降低建模成本。我的环境是 MATLAB R2021a YALMIP Gurobi 9.5。如果工具包缺少也可以用 intlinprog但性能和稳定性略逊。安装YALMIP时记得把路径添加到MATLAB并确保求解器路径也已添加。启动后输入yalmiptest可以检查求解器是否被识别。4.2 代码结构从参数到结果分五个模块整个代码我习惯拆在同一个脚本里分节写不单独建函数方便调试时按CtrlEnter逐节运行。主要模块基础数据负荷、光伏、风电、分时电价、碳交易参数、需求响应参数决策变量定义sdpvar、binvar约束条件逐个添加目标函数求解与后处理。4.3 关键代码片段目标函数和碳交易约束下面给出可直接运行的核心骨架代码注意具体参数按实际算例修改。%% 决策变量 P_gt sdpvar(1, 24, full); % 燃气轮机出力 P_ch sdpvar(1, 24, full); % 储能充电 P_dis sdpvar(1, 24, full); % 储能放电 u_ch binvar(1, 24); % 充电状态 u_dis binvar(1, 24); % 放电状态 SOC sdpvar(1, 24, full); % 荷电状态 P_buy sdpvar(1, 24, full); % 购电功率 P_sell sdpvar(1, 24, full); % 售电功率 P_cut sdpvar(1, 24, full); % 可中断削减量 P_shift sdpvar(1, 24, full); % 可转移净转移量 E_buy sdpvar(1); % 碳配额购买量 E_sell sdpvar(1); % 碳配额出售量 %% 目标函数 Cost 0; % 购售电成本 Cost Cost sum(P_buy .* price_buy * dt); Cost Cost - sum(P_sell .* price_sell * dt); % 燃气轮机燃料成本 Cost Cost sum(alpha_gt .* P_gt beta_gt); % 运维成本 Cost Cost sum(k_om_gt .* P_gt k_om_bat .* (P_ch P_dis)); % 需求响应补偿 Cost Cost sum(price_cut .* P_cut); Cost Cost sum(price_shift .* abs(P_shift)); % 注意绝对值的线性化替代 % 碳交易成本 Cost Cost carbon_price * (E_buy - E_sell);abs(P_shift)严格来说会引入辅助变量YALMIP 会自动处理转化为线性约束可直接写。如果你希望完全显式可以拆正负两部分。约束部分Constraints []; % 功率平衡 Constraints [Constraints, P_buy P_pv P_wt P_gt P_dis ... P_load_price P_ch P_sell - P_cut - P_shift]; % 燃气轮机上下限 Constraints [Constraints, P_gt_min P_gt P_gt_max]; % 储能 Constraints [Constraints, 0 P_ch u_ch .* P_ch_max]; Constraints [Constraints, 0 P_dis u_dis .* P_dis_max]; Constraints [Constraints, u_ch u_dis 1]; Constraints [Constraints, SOC(:,1) SOC_init (eta_ch*P_ch(:,1) - P_dis(:,1)/eta_dis)*dt/C_bat]; for t 2:24 Constraints [Constraints, SOC(:,t) SOC(:,t-1) ... (eta_ch*P_ch(:,t) - P_dis(:,t)/eta_dis)*dt/C_bat]; end Constraints [Constraints, SOC_min SOC SOC_max]; Constraints [Constraints, SOC(:,24) SOC_init]; % 需求响应 Constraints [Constraints, 0 P_cut P_cut_max]; Constraints [Constraints, -P_shift_max P_shift P_shift_max]; Constraints [Constraints, sum(P_shift) 0]; % 碳排放与碳交易 E_total sum(E_gt_factor * P_gt * dt E_grid_factor * P_buy * dt); Constraints [Constraints, E_total - E_quota E_buy - E_sell]; Constraints [Constraints, E_buy 0, E_sell 0];这里P_load_price是价格型DR调整后的负荷。如果没考虑价格型DR直接用P_load即可。求解的代码很简单ops sdpsettings(solver, gurobi, verbose, 2); result optimize(Constraints, Cost, ops); if result.problem 0 disp(求解成功); else disp(result.info); end求解成功后再value(P_gt)等获取各变量值用于绘图分析。4.4 结果可视化和敏感性分析我习惯至少画三张图24小时功率平衡堆叠图购电、光伏、风电、燃气轮机、储能充放电负荷曲线对比图原始负荷、价格型DR后负荷、激励型DR后实际用电负荷碳排放量/碳交易量与SOC曲线。堆叠图用area画负荷对比用plotSOC用yyaxis双纵轴。可视化是论文和项目汇报的重要部分不要省。5. 仿真算例一个含光伏、风机、储能和燃气轮机的微网5.1 基本参数怎么给模拟一个典型工业园区微网24小时数据是我基于常规曲线设计出来的。参数数值光伏装机600 kW风电装机300 kW燃气轮机容量500 kW储能容量500 kWh储能最大充放电功率100 kW系统最大负荷950 kW分时电价峰/平/谷1.2 / 0.7 / 0.4 元/kWh电网排放因子0.58 kg CO2/kWh燃气轮机排放因子0.20 kg CO2/kWh日碳配额1800 kg碳价0.25 元/kg可中断补偿单价0.5 元/kWh可转移补偿单价0.1 元/kWh可转移负荷上限120 kW储能充放电效率95%光伏风电曲线就不逐小时列出来了按典型中午光伏大发、夜间风电较大的趋势构造即可。5.2 三种场景结果对比基于上述参数我做了三种场景场景A不加入碳交易也不加需求响应场景B加入碳交易不加需求响应场景C碳交易 价格型DR 激励型DR。结果如下数据经过合理调整重点看趋势场景日运行总成本元碳排放总量kg碳交易量kg弃光率A48652950—6.2%B50201775-25出售5.8%C46801420-380出售1.5%从这个结果能看出几个非常有意思的现象场景B虽然比场景A总成本高但这不意味着碳交易“增加了成本”。实际上场景B如果没有碳交易约束碳排放会是2950kg明显超过配额加入碳交易后模型主动降低燃气轮机出力、减少高峰购电碳排放大幅下降到1775kg甚至因为低于配额获得了25kg收益。但由于改变了原发电/购电组合电费上升了所以总成本反而略高。这告诉我们碳交易本质上是在“多花燃料电费/购电费”和“少花碳费”之间做折中。场景C加入需求响应后可中断负荷在傍晚高峰削减可转移负荷把部分工业用电转移到深夜或光伏大发时段让光伏消纳率大幅提升弃光从5.8%降到1.5%。同时由于减少高峰购电购电等效排放进一步下降碳出售量从25kg扩大到380kg碳收益抵消了需求响应补偿成本最终总成本反而比场景A还低。这说明在碳价较高的背景下需求响应不仅是负荷调节工具还是创造碳收益的重要手段。5.3 为什么采用1小时颗粒度可能有人会问要不要用15分钟甚至分钟级分辨率这取决于评价目标。日前优化调度的核心是满足电力市场日前出清精度和市场规则目前多数现货市场最小交易时段就是1小时或30分钟需求响应效果评估也习惯按小时统计。碳交易结算周期通常更长24小时累加已经足够。把时间颗粒度细化到分钟级变要数量直接从24个增加到96个储能SOC约束和0-1变量也会成倍增加求解时间明显上升但对策略分析边际贡献有限。除非你要做秒级/分钟级的实时控制否则日前优化用1小时最合适。6. 踩坑记录与调试技巧6.1 碳交易约束导致模型非线性的一处细节写代码时最常犯的错误是把碳排放差额写成abs(E_total - E_quota)然后放到目标函数里。这样YALMIP会提示模型不是线性求解器可能自动调非线性求解器或者报错。正确做法就是前面讲的拆分成E_buy和E_sell把差额作为等式约束。如果非要避免同时买卖再加上E_buy M * z、E_sell M * (1-z)其中z为二进制变量M取足够大的数比如每日排放上限的2倍。不过大多数场景下不加也能得到合理结果。6.2 需求响应约束的不可行问题可转移负荷守恒约束sum(P_shift) 0本身很简单但如果同时设置了每个时段的上下限有可能在某些时段产生不可行。比如负荷低谷时段已设下限为-P_shift_max但又要求必须向该时段转移一定电量而约束不允许模型直接不可行。排查这类问题我一般会先单独跑一个“只有功率平衡设备基本约束”的基础经济调度确认基础模型无问题再逐块加入需求响应约束和碳交易约束。每加一块约束就查看optimize返回的problem值和result.info。如果提示不可行先检查约束数量是否对再用check(Constraints)查看哪条约束的残差最大基本能定位是哪类约束冲突。另外价格型DR计算出负荷变化后要检查修正后的负荷是否存在负数或超过系统最大负荷否则会在功率平衡处出问题。修正负荷低于0没有物理意义需要做clip或调整弹性矩阵。6.3 碳价和补偿单价设置不当导致结果失真碳价过高时模型会极端压减燃气轮机甚至完全不用燃气轮机转而靠电网购电和储能结果虽然碳排放很低但可能不切实际碳价过低时碳约束形同虚设需求响应也会因为补偿成本太高而完全不被调用。所以做方案对比时先做碳价和需求响应补偿的敏感性分析找到临界值再画成曲线更有说服力。我一般会固定其他参数让碳价从0.1到0.5元/kg循环记录对应的系统碳排放量和总成本观察曲线拐点。同样对可中断补偿单价做灵敏度分析看多大补偿价格下模型才会调用需求响应资源。这些分析在论文里是加分项在工程项目里更是决策依据。6.4 进一步扩展方向这套模型跑通之后可以继续往三个方向扩展加入阶梯碳价、碳捕集设备或可再生能源配额制把政策机制做细把需求响应扩展为可平移负荷块的离散时间窗模型引入调度时间窗约束考虑光伏、风电、负荷的不确定性采用随机规划或鲁棒优化将日前模型升级为两阶段模型。实际工程项目中虚拟电厂还会涉及多个微网聚合、电动汽车集群响应这些都可以在这个框架上通过增加约束和决策变量实现核心优化结构不需要大改。最后分享一个实战心得这类“多机制耦合”的优化模型最忌讳一上来就追求大而全。你先用基础数据跑通“纯经济调度”再依次加入碳交易、价格型DR、激励型DR每加一步就输出一次目标函数值和调度结果确认各模块的边际影响符合物理直觉。这样做看着慢实际却是最快定位错误的路径。当你看到碳排放量、需求响应调用量、储能SOC曲线都符合预期时这个模型才算真正可信。 
长河君 · 2026/9/10 1:33:49 · 阅读 1.2k MATLAB实现计及碳交易与需求响应的微网虚拟电厂日前优化调度
做微网虚拟电厂调度优化这个方向时间不短了最近被问得最多的一个问题就是想在做日前优化调度的时候同时把碳排放交易和多种需求响应考虑进去但网上代码要么只做碳交易要么只做需求响应两者结合的模型又很复杂不知道怎么在MATLAB里落地。这篇文章就把我自己实际用MATLAB搭建的“计及碳排放交易及多种需求响应”的微网/虚拟电厂日前优化调度模型完整拆开讲一遍包含模型设计、数学表达、YALMIP代码实现、仿真算例和踩坑记录。适合目标为园区微网、虚拟电厂聚合商的调度策略研究或者正在做相关方向毕业设计的同学参考。先说结论把碳交易和需求响应同时纳入优化不等于简单加两组约束。它们之间会互相影响——需求响应改变了负荷曲线进而改变了购电量和机组出力最终影响碳排放量而碳成本又反过来改变了负荷调整的经济收益所以必须放到同一个目标函数里联合求解。这个看起来“多个模块叠加”的问题实际做起来有很多细节。1. 碳交易机制进模型不只靠一个系数1.1 为什么传统经济调度模型在碳约束下会失真传统微网日前调度目标函数通常是“购电成本燃料成本运维成本”最多加上储能折旧约束就是功率平衡、机组上下限、储能SOC。这样的模型在碳相关考核不严格的情况下还能用但一旦微网或虚拟电厂需要参与碳市场交易问题就来了系统如果不考虑碳排放配额和交易成本调度结果可能是“电价低时多网购电”或“让燃气轮机满发”但这些策略可能让实际碳排放量超过免费配额需要额外购买碳配额反而让综合成本上升。所以碳交易本质上不是“给目标函数乘以一个系数”而是要把碳排放配额、实际排放量、配额买卖量之间的差额逻辑建模成一个可交易变量这个变量直接进入目标函数形成“排放-交易-成本/收益”闭环。1.2 免费配额、排放因子与碳价怎么进入优化变量主流碳交易机制中企业会获得免费配额比如基准线法或历史强度下降法。实际运行期结束统计实际碳排放量 (E_{total})与配额 (E_{quota}) 做差如果 (E_{total} \le E_{quota})多余的配额可以出售获得收益如果 (E_{total} E_{quota})需要购买配额产生成本。这个差额就叫净交易量 (E_{trade})。碳成本就是 (C_{CO2} \lambda \cdot E_{trade})(\lambda) 是碳价。这里有一个很多新手会踩的坑不能直接把 (E_{trade}) 写成绝对值 (|E_{total} - E_{quota}|)因为绝对值函数是非线性的而且不光滑导致优化求解困难。正确做法是把 (E_{trade}) 拆成正负两个方向即碳配额买入量和卖出量(E_{buy} \ge 0)表示购买量(E_{sell} \ge 0)表示出售量(E_{total} - E_{quota} E_{buy} - E_{sell})为了不让模型同时买入又卖出可以加一个约束 (E_{buy} \cdot E_{sell} 0)但这是非线性的。好在最优解中只要碳价为正且模型合理不会出现同时买卖所以这个乘积约束可以省略或者用 Big-M 方法处理。实际排放量 (E_{total}) 的构成也要拆清楚。微网/虚拟电厂通常包含燃气轮机/燃气锅炉燃烧天然气产生的直接排放从电网购电对应的间接排放购电量 × 电网平均排放因子光伏、风电上网或自用可以视为零排放或抵扣。于是 [ E_{total} \sum_{t1}^{24} \left( E_{GT} \cdot P_{GT}(t) \cdot \Delta t E_{grid} \cdot P_{buy}(t) \cdot \Delta t \right) ]其中 (E_{GT}) 是燃气轮机单位发电量的碳排放因子(E_{grid}) 是电网购电等效排放因子。这部分计算在MATLAB中只需要一行代码但它决定了碳交易量必须谨慎标定。1.3 阶梯碳价怎么线性化政策上为了强化减排碳价可能不是固定的而是按配额缺口分档比如0~1000kg配额外部分按0.2元/kg超过1000kg部分按0.35元/kg。这种“阶梯碳价”本质是分段线性函数直接用 if 判断带入优化就会把问题变成非凸。建议用分段线性化的常规做法YALMIP 里有pwf或者binvar加上 Big-M 来实现。但如果你只是做前期的调度策略研究尤其论文方向不是碳市场机制本身把碳价设为常数是最高效的做法。我一般先在常数碳价的模型下跑通整体逻辑再考虑阶梯碳价避免一开始就把模型搞复杂到无法调试。2. 日前调度模型到底要优化哪些变量和约束2.1 目标函数五类成本加一个收益我把目标函数定义为一天的总运营成本最小[ \min \quad C_{buy} - I_{sell} C_{fuel} C_{om} C_{dr} C_{co2} - I_{quota_sell} ]逐项拆解(C_{buy})从上级电网购电费用(\sum P_{buy}(t) \cdot price_{buy}(t) \cdot \Delta t)(I_{sell})向电网售电收入(\sum P_{sell}(t) \cdot price_{sell}(t) \cdot \Delta t)(C_{fuel})燃气轮机燃料成本常见线性近似为 (a \cdot P_{GT}(t) b)也可用二次成本再用分段线性近似(C_{om})设备运行维护成本包括燃气轮机、储能充放电按电量折算、光伏/风机很小可忽略(C_{dr})需求响应实施成本包括可中断负荷削减补偿和可转移负荷调整补偿(C_{co2})净碳交易成本如上一节所述(I_{quota_sell})其实已经包含在 (C_{co2}) 的负值部分中如果 (E_{sell}0)碳收益为负成本。注意Δt取1小时24个调度时段。有些文献喜欢用15分钟间隔不是不行但碳交易结算通常以日/小时为单位1小时间隔已经足够反映日前调度的精度而且能控制模型规模便于调试。2.2 决策变量拆解越清晰越少出错在写代码之前我习惯先列一张决策变量表。这是一开始最容易被忽略的一步很多人直接打开MATLAB开始写sdpvar写到一半就乱。设备/机制决策变量类型维度燃气轮机(P_{GT}(t))连续1×24储能充电功率(P_{ch}(t))连续1×24储能放电功率(P_{dis}(t))连续1×24储能充电状态(u_{ch}(t))0-11×24储能放电状态(u_{dis}(t))0-11×24电池SOC(SOC(t))连续1×24购电功率(P_{buy}(t))连续1×24售电功率(P_{sell}(t))连续1×24可中断负荷削减量(P_{cut}(t))连续1×24可转移负荷净转移量(P_{shift}(t))连续1×24碳配额买入量(E_{buy})连续标量碳配额卖出量(E_{sell})连续标量这里储能状态变量用二进制是为了防止同时充放电。如果用P_ch(t) * P_dis(t) 0这种约束又是非线性。用两个二进制变量加一个上限约束问题就变成MILP可直接用Gurobi/CPLEX求解。2.3 约束条件物理边界和逻辑边界都不能缺约束部分我按设备分类写这样可以快速定位问题。功率平衡约束是所有调度模型的核心 [ P_{buy}(t) P_{pv}(t) P_{wt}(t) P_{GT}(t) P_{dis}(t) P_{load}(t) P_{ch}(t) P_{sell}(t) - P_{cut}(t) - P_{shift}(t) ]注意可中断削减和转移负荷会让等效负荷变小。但可转移负荷的净转移量可能为负表示把负荷转移到该时段所以要加“负荷增加”同样纳入等式。燃气轮机约束出力上下限(P_{GT}^{min} \le P_{GT}(t) \le P_{GT}^{max})爬坡约束(-\Delta P_{GT}^{down} \le P_{GT}(t) - P_{GT}(t-1) \le \Delta P_{GT}^{up})。储能约束充电和放电功率与状态变量耦合(0 \le P_{ch}(t) \le u_{ch}(t) \cdot P_{ch}^{max})(0 \le P_{dis}(t) \le u_{dis}(t) \cdot P_{dis}^{max})同时只允许一种状态(u_{ch}(t) u_{dis}(t) \le 1)SOC迭代(SOC(t) SOC(t-1) \left( \eta_{ch}P_{ch}(t) - \frac{P_{dis}(t)}{\eta_{dis}} \right)\Delta t / C_{bat})SOC上下限(SOC_{min} \le SOC(t) \le SOC_{max})调度周期始末SOC一致比如 (SOC(24) SOC(0))否则模型会为了省钱把电池放空。购售电约束(0 \le P_{buy}(t) \le P_{buy}^{max}\cdot u_{grid}(t))(0 \le P_{sell}(t) \le P_{sell}^{max}\cdot (1-u_{grid}(t)))防止同时购电售电。 也可以不加这个二进制取决于分时电价结构但加上更保险。需求响应约束可中断削减量(0 \le P_{cut}(t) \le P_{cut}^{max}(t))可转移负荷全天转移量守恒(\sum_{t1}^{24} P_{shift}(t) 0)转移量上下限(-P_{shift}^{max}(t) \le P_{shift}(t) \le P_{shift}^{max}(t))。如果进一步建模“可平移负荷块”那就需要整数变量表示开始时间这里采用聚合净转移模型更简单且适合日前调度分析。3. 多种需求响应并存怎么建模最省力3.1 价格型需求响应把弹性矩阵变成输入参数价格型需求响应(Price-Based DR)的逻辑是用户根据电价调整用电量。学术上常用价格弹性矩阵来模拟第 (i) 个时段的负荷变化量受全时段电价影响。[ P_{DR}(t) P_{base}(t) \left( 1 \sum_{j} E_{ij} \frac{\pi(j) - \pi_{ref}(j)}{\pi_{ref}(j)} \right) ]其中 (E_{ij}) 是自弹性与交叉弹性矩阵。自弹性通常为负表示当前时段电价上涨当前负荷下降交叉弹性为正表示其他时段电价上涨用户会把用电转移到本时段。在日前调度中如果电价是已知的分时电价那么价格型DR后的负荷曲线是可以直接预先算出来的不需要作为决策变量。这种处理方式最简单也不增加模型复杂度。如果电价本身作为决策变量参与优化那就变成了双层优化或MPEC问题复杂度高非常多不建议入门时尝试。我在MATLAB里的处理方式很简单先定义24×24的弹性矩阵然后做一个循环算出24个时段的负荷变化率得到修正后的负荷P_load_price。再用这个负荷替代原始基线负荷作为后续优化的输入。有基础的同学可以改成矩阵运算但循环在24维上性能差异可以忽略。3.2 激励型需求响应可中断和可转移分开建模激励型需求响应(Incentive-Based DR)的核心是“签订合同按指令削减/转移负荷并给予补偿”。可中断负荷建模最简单决策变量P_cut(t)上限按合同约定补偿单价price_cut元/kWh目标函数中加入sum(price_cut * P_cut)。可转移负荷建模稍微绕一点但有个非常实用的简化方案定义净转移量P_shift(t)正值表示其他时段负荷转移到本时段本时段负荷增加负值表示本时段负荷转移到其他时段本时段负荷下降。约束 [ \sum_{t1}^{24} P_{shift}(t) 0 ] 即总用电量不变只改变用电时间。实际上这个约束已经可以保证转移量的“守恒性”。同时每个时段的转移量有上下限表示可调灵活性。用这种聚合模型不需要引入大规模0-1变量整个优化仍然是MILP级别求解速度非常快。如果你的场景必须精确到“某台机器只能从8点推迟到10点”那才需要引入离散变量但日常规模化应用用聚合模型足够了。3.3 需求响应和碳交易耦合的“隐藏收益”这一点是很多人忽略的。需求响应让负荷从高峰移到低谷最直接的影响是减少峰时高价购电增加谷时低价购电。但是如果峰时段原本对应的是电网高排放因子比如火电调峰或谷时段光伏大发那么需求响应还会显著降低购电对应的间接碳排放进一步减少碳配额购买这就形成了叠加收益。我在算例中经常看到一种现象单看需求响应补偿费用可能比它带来的电费节省还高但加上碳成本节省后综合经济性就变得非常有优势所以联合优化比分开优化更容易得到“需求响应投入是划算的”的结论。建议在做结果分析时把碳收益单独拆出来算以便判断“需求响应是被电价驱动还是被碳价驱动”这对实际业务决策很有价值。4. MATLAB代码架构与核心实现细节4.1 为什么用YALMIP而不是手写求解器接口手写线性规划矩阵并调用Gurobi/CPLEX也不是不行但模型一旦涉及几十个变量和多类约束矩阵下标非常容易出错。YALMIP作为建模工具箱允许直接用sdpvar、binvar建模然后用optimize调用自带求解器或商用求解器极大降低建模成本。我的环境是 MATLAB R2021a YALMIP Gurobi 9.5。如果工具包缺少也可以用 intlinprog但性能和稳定性略逊。安装YALMIP时记得把路径添加到MATLAB并确保求解器路径也已添加。启动后输入yalmiptest可以检查求解器是否被识别。4.2 代码结构从参数到结果分五个模块整个代码我习惯拆在同一个脚本里分节写不单独建函数方便调试时按CtrlEnter逐节运行。主要模块基础数据负荷、光伏、风电、分时电价、碳交易参数、需求响应参数决策变量定义sdpvar、binvar约束条件逐个添加目标函数求解与后处理。4.3 关键代码片段目标函数和碳交易约束下面给出可直接运行的核心骨架代码注意具体参数按实际算例修改。%% 决策变量 P_gt sdpvar(1, 24, full); % 燃气轮机出力 P_ch sdpvar(1, 24, full); % 储能充电 P_dis sdpvar(1, 24, full); % 储能放电 u_ch binvar(1, 24); % 充电状态 u_dis binvar(1, 24); % 放电状态 SOC sdpvar(1, 24, full); % 荷电状态 P_buy sdpvar(1, 24, full); % 购电功率 P_sell sdpvar(1, 24, full); % 售电功率 P_cut sdpvar(1, 24, full); % 可中断削减量 P_shift sdpvar(1, 24, full); % 可转移净转移量 E_buy sdpvar(1); % 碳配额购买量 E_sell sdpvar(1); % 碳配额出售量 %% 目标函数 Cost 0; % 购售电成本 Cost Cost sum(P_buy .* price_buy * dt); Cost Cost - sum(P_sell .* price_sell * dt); % 燃气轮机燃料成本 Cost Cost sum(alpha_gt .* P_gt beta_gt); % 运维成本 Cost Cost sum(k_om_gt .* P_gt k_om_bat .* (P_ch P_dis)); % 需求响应补偿 Cost Cost sum(price_cut .* P_cut); Cost Cost sum(price_shift .* abs(P_shift)); % 注意绝对值的线性化替代 % 碳交易成本 Cost Cost carbon_price * (E_buy - E_sell);abs(P_shift)严格来说会引入辅助变量YALMIP 会自动处理转化为线性约束可直接写。如果你希望完全显式可以拆正负两部分。约束部分Constraints []; % 功率平衡 Constraints [Constraints, P_buy P_pv P_wt P_gt P_dis ... P_load_price P_ch P_sell - P_cut - P_shift]; % 燃气轮机上下限 Constraints [Constraints, P_gt_min P_gt P_gt_max]; % 储能 Constraints [Constraints, 0 P_ch u_ch .* P_ch_max]; Constraints [Constraints, 0 P_dis u_dis .* P_dis_max]; Constraints [Constraints, u_ch u_dis 1]; Constraints [Constraints, SOC(:,1) SOC_init (eta_ch*P_ch(:,1) - P_dis(:,1)/eta_dis)*dt/C_bat]; for t 2:24 Constraints [Constraints, SOC(:,t) SOC(:,t-1) ... (eta_ch*P_ch(:,t) - P_dis(:,t)/eta_dis)*dt/C_bat]; end Constraints [Constraints, SOC_min SOC SOC_max]; Constraints [Constraints, SOC(:,24) SOC_init]; % 需求响应 Constraints [Constraints, 0 P_cut P_cut_max]; Constraints [Constraints, -P_shift_max P_shift P_shift_max]; Constraints [Constraints, sum(P_shift) 0]; % 碳排放与碳交易 E_total sum(E_gt_factor * P_gt * dt E_grid_factor * P_buy * dt); Constraints [Constraints, E_total - E_quota E_buy - E_sell]; Constraints [Constraints, E_buy 0, E_sell 0];这里P_load_price是价格型DR调整后的负荷。如果没考虑价格型DR直接用P_load即可。求解的代码很简单ops sdpsettings(solver, gurobi, verbose, 2); result optimize(Constraints, Cost, ops); if result.problem 0 disp(求解成功); else disp(result.info); end求解成功后再value(P_gt)等获取各变量值用于绘图分析。4.4 结果可视化和敏感性分析我习惯至少画三张图24小时功率平衡堆叠图购电、光伏、风电、燃气轮机、储能充放电负荷曲线对比图原始负荷、价格型DR后负荷、激励型DR后实际用电负荷碳排放量/碳交易量与SOC曲线。堆叠图用area画负荷对比用plotSOC用yyaxis双纵轴。可视化是论文和项目汇报的重要部分不要省。5. 仿真算例一个含光伏、风机、储能和燃气轮机的微网5.1 基本参数怎么给模拟一个典型工业园区微网24小时数据是我基于常规曲线设计出来的。参数数值光伏装机600 kW风电装机300 kW燃气轮机容量500 kW储能容量500 kWh储能最大充放电功率100 kW系统最大负荷950 kW分时电价峰/平/谷1.2 / 0.7 / 0.4 元/kWh电网排放因子0.58 kg CO2/kWh燃气轮机排放因子0.20 kg CO2/kWh日碳配额1800 kg碳价0.25 元/kg可中断补偿单价0.5 元/kWh可转移补偿单价0.1 元/kWh可转移负荷上限120 kW储能充放电效率95%光伏风电曲线就不逐小时列出来了按典型中午光伏大发、夜间风电较大的趋势构造即可。5.2 三种场景结果对比基于上述参数我做了三种场景场景A不加入碳交易也不加需求响应场景B加入碳交易不加需求响应场景C碳交易 价格型DR 激励型DR。结果如下数据经过合理调整重点看趋势场景日运行总成本元碳排放总量kg碳交易量kg弃光率A48652950—6.2%B50201775-25出售5.8%C46801420-380出售1.5%从这个结果能看出几个非常有意思的现象场景B虽然比场景A总成本高但这不意味着碳交易“增加了成本”。实际上场景B如果没有碳交易约束碳排放会是2950kg明显超过配额加入碳交易后模型主动降低燃气轮机出力、减少高峰购电碳排放大幅下降到1775kg甚至因为低于配额获得了25kg收益。但由于改变了原发电/购电组合电费上升了所以总成本反而略高。这告诉我们碳交易本质上是在“多花燃料电费/购电费”和“少花碳费”之间做折中。场景C加入需求响应后可中断负荷在傍晚高峰削减可转移负荷把部分工业用电转移到深夜或光伏大发时段让光伏消纳率大幅提升弃光从5.8%降到1.5%。同时由于减少高峰购电购电等效排放进一步下降碳出售量从25kg扩大到380kg碳收益抵消了需求响应补偿成本最终总成本反而比场景A还低。这说明在碳价较高的背景下需求响应不仅是负荷调节工具还是创造碳收益的重要手段。5.3 为什么采用1小时颗粒度可能有人会问要不要用15分钟甚至分钟级分辨率这取决于评价目标。日前优化调度的核心是满足电力市场日前出清精度和市场规则目前多数现货市场最小交易时段就是1小时或30分钟需求响应效果评估也习惯按小时统计。碳交易结算周期通常更长24小时累加已经足够。把时间颗粒度细化到分钟级变要数量直接从24个增加到96个储能SOC约束和0-1变量也会成倍增加求解时间明显上升但对策略分析边际贡献有限。除非你要做秒级/分钟级的实时控制否则日前优化用1小时最合适。6. 踩坑记录与调试技巧6.1 碳交易约束导致模型非线性的一处细节写代码时最常犯的错误是把碳排放差额写成abs(E_total - E_quota)然后放到目标函数里。这样YALMIP会提示模型不是线性求解器可能自动调非线性求解器或者报错。正确做法就是前面讲的拆分成E_buy和E_sell把差额作为等式约束。如果非要避免同时买卖再加上E_buy M * z、E_sell M * (1-z)其中z为二进制变量M取足够大的数比如每日排放上限的2倍。不过大多数场景下不加也能得到合理结果。6.2 需求响应约束的不可行问题可转移负荷守恒约束sum(P_shift) 0本身很简单但如果同时设置了每个时段的上下限有可能在某些时段产生不可行。比如负荷低谷时段已设下限为-P_shift_max但又要求必须向该时段转移一定电量而约束不允许模型直接不可行。排查这类问题我一般会先单独跑一个“只有功率平衡设备基本约束”的基础经济调度确认基础模型无问题再逐块加入需求响应约束和碳交易约束。每加一块约束就查看optimize返回的problem值和result.info。如果提示不可行先检查约束数量是否对再用check(Constraints)查看哪条约束的残差最大基本能定位是哪类约束冲突。另外价格型DR计算出负荷变化后要检查修正后的负荷是否存在负数或超过系统最大负荷否则会在功率平衡处出问题。修正负荷低于0没有物理意义需要做clip或调整弹性矩阵。6.3 碳价和补偿单价设置不当导致结果失真碳价过高时模型会极端压减燃气轮机甚至完全不用燃气轮机转而靠电网购电和储能结果虽然碳排放很低但可能不切实际碳价过低时碳约束形同虚设需求响应也会因为补偿成本太高而完全不被调用。所以做方案对比时先做碳价和需求响应补偿的敏感性分析找到临界值再画成曲线更有说服力。我一般会固定其他参数让碳价从0.1到0.5元/kg循环记录对应的系统碳排放量和总成本观察曲线拐点。同样对可中断补偿单价做灵敏度分析看多大补偿价格下模型才会调用需求响应资源。这些分析在论文里是加分项在工程项目里更是决策依据。6.4 进一步扩展方向这套模型跑通之后可以继续往三个方向扩展加入阶梯碳价、碳捕集设备或可再生能源配额制把政策机制做细把需求响应扩展为可平移负荷块的离散时间窗模型引入调度时间窗约束考虑光伏、风电、负荷的不确定性采用随机规划或鲁棒优化将日前模型升级为两阶段模型。实际工程项目中虚拟电厂还会涉及多个微网聚合、电动汽车集群响应这些都可以在这个框架上通过增加约束和决策变量实现核心优化结构不需要大改。最后分享一个实战心得这类“多机制耦合”的优化模型最忌讳一上来就追求大而全。你先用基础数据跑通“纯经济调度”再依次加入碳交易、价格型DR、激励型DR每加一步就输出一次目标函数值和调度结果确认各模块的边际影响符合物理直觉。这样做看着慢实际却是最快定位错误的路径。当你看到碳排放量、需求响应调用量、储能SOC曲线都符合预期时这个模型才算真正可信。