做微网/虚拟电厂的日前优化调度，把“计及碳排放交易及多种需求响应”一起塞进模型里，再用Matlab把这套优化跑通，听起来像是一个标准到不能再标准的课题组合。实际上这也是当前调度方向里很多人正在复现的内容：单纯追求运行成本最低已经不够，碳成本正成为一个真实的价格信号；负荷也不再适合当作固定边界，多种需求响应会让负荷从“刚性”变成“可调”。如果你刚入手这类调度题目，或者你已经有基础模型，想往里面加碳交易、加需求响应，这篇应该能帮你省不少时间。下面我会把问题拆开讲，从建模思路、数学约束、Matlab代码结构到常见坑位，一次性梳理清楚。
1. 项目整体设计：碳成本与负荷弹性如何进入调度
1.1 微网与虚拟电厂调度的共同骨架
很多人一开始会被“微网”和“虚拟电厂”这两个词搞混，其实在日前优化调度这一层，不需要把它们的差异想得太复杂。微网更侧重区域内的供能自洽，虚拟电厂更偏重把分散的分布式电源、储能、可控负荷聚合起来，像一个“虚拟机组”参与电网调度。两者落到优化模型上都差不多：时间尺度是未来24小时，步长通常取1小时或15分钟，给定光伏、风电和负荷预测曲线，求解各机组出力、储能充放电、与上级电网购售电计划，以及需求响应执行量。
这个模型的核心逻辑可以概括成一句话：在所有技术约束条件下，让系统运行总成本尽量低。但成本怎么定义，会直接决定调度结果。传统做法里，成本主要是燃料成本、购电成本、设备启停维护成本。加入碳排放交易后，成本列表里多了碳配额买卖成本；加入需求响应后，负荷侧不再是一个固定的曲线，可中断、可转移负荷都变成了决策变量，同时它们的调用会产生补偿成本。这些成本彼此耦合，才是这个项目的真正难点。
1.2 碳排放交易机制不是简单惩罚，而是一种价格信号
碳排放交易进入优化模型后，最常见也最好理解的处理方式是“基准线配额法”。调度中心根据历史电量、负荷峰值或者发电量，给系统分配一个日度免费碳排放配额。实际运行中，机组燃烧燃料会产生碳排放，从上级电网购电也可以按电网平均排放强度折算成间接排放。如果实际碳排放量小于免费配额，则可以在碳市场出售多余配额，获得收益；如果实际排放量超过免费配额，则需要额外购买碳排放配额，产生成本。
把这个机制写成数学表达式其实就是：
实际碳排放量 = ∑ 机组出力 × 机组碳排放强度 × Δt + 上级电网购电量 × 电网平均排放强度 × Δt 碳交易成本 = 碳价 × (实际碳排放量 - 免费配额)从中能看出，碳价和免费配额会共同影响调度决策。免费配额给得高，系统可能并不急于压碳；配额给得紧，模型就会自然减少高排放机组的出力，转向更清洁的光伏、风电、储能，或者增加从电网买入低碳电力。注意这里“实际排放量”不是运行后统计出来的，而是在优化前由决策变量表达出来的，所以碳交易成本才能进目标函数。如果要做更复杂的阶梯碳价，也就是排放量超过某个区间后碳价上涨，那还需要引入分段的碳价格函数，后面我再专门说坑。
1.3 多种需求响应：可中断、可转移、价格型三类先分清
需求响应种类很多，但在日前调度模型里，最常出现的是三类。
第一类是可中断负荷，也叫可削减负荷。用户在高峰时段接受调度指令降低一部分用电，系统根据削减量给予补偿。这类负荷进入模型时，通常是一个带上下限的连续变量，同时还可以限制每天最大中断次数或总中断时长。
第二类是可转移负荷，比如洗衣机、电动车充电、蓄热式电锅炉这类可在时间上移动的用电。它的特点是总用电量不变，但用电时段可以从高峰挪到低谷。模型里需要满足一天内“转入某时段的总电量等于转出某时段的总电量”这个守恒关系。
第三类是价格型需求响应。这类响应不依赖调度指令，而是用户自行对分时电价做出回应，高电价少用，低电价多用。工程上一般用需求价格弹性矩阵描述，在模型外面先做一次负荷预测修正。
“多种需求响应”这四个字，并不只是把变量都加进去就行，关键的是分清哪些在优化内被调度，哪些在优化前做预测修正。我的做法是：基础负荷先经过价格型需求响应修正，得到一条新的基准负荷曲线；在优化模型内部，只把可中断和可转移负荷作为决策变量。这样既满足多种需求响应，又不至于把模型搞到不可解。
| 需求响应类型 | 常见建模方式 | 补偿/成本特征 | 典型使用时段 |
|---|---|---|---|
| 可中断/可削减 | 0到上限的连续削减变量，可加中断次数约束 | 按削减电量给补偿 | 高峰时段削峰 |
| 可转移 | 转入、转出两个连续变量，要求总电量守恒 | 按转移量或调用次数给补偿 | 峰谷时段电量平移 |
| 价格型 | 需求价格弹性矩阵修正基础负荷 | 不直接在优化中计费，反映在电价套利中 | 全天负荷预测 |
2. 数学模型：把碳成本与需求响应写成可求解的约束
2.1 目标函数拆解：五个成本项一次说清
我建议把目标函数拆成五个部分，这样代码和调试都会清晰很多。
第一部分是机组发电成本。微型燃气轮机、柴油机这类可控机组通常用二次函数描述，形如：
C_g = a_g × P_g^2 + b_g × P_g + c_g其中P_g是该机组在时段t的出力，系数a_g、b_g、c_g由机组热效率和燃料价格拟合得到。开停机还有额外成本，所以也要在目标函数中反映，一般用启动成本乘以启动状态变量。
第二部分是与上级电网交互的购售电成本。假设分时电价已知，购电成本是购电量乘以购电价，售电收益是售电量乘以售电价。这里需要说明，如果只考虑经济性，系统会尽量低谷买电、高峰卖电来套利，所以通常要加购电和售电互斥约束，不能让同一时段又买又卖。
第三部分是碳交易成本。碳排放量由机组发电排放和电网购电间接排放两部分组成，免费配额作为常数向量放进模型中，碳价乘以排放缺口或者排放富余量，得到碳交易成本。碳价取值可以是固定常数，也可以是阶梯函数。
第四部分是需求响应补偿成本。可中断负荷削减了用户的用电，需要按削减电量支付补偿；可转移负荷虽然总用电量不变，但调整用电时段会影响用户生产，也需要给一定补偿。这部分成本一定要进目标函数，否则模型会“免费”滥用需求响应，得出的方案没有实际意义。
第五部分是储能和其他运维成本。储能运维成本可以用充放电电量乘以单位运维成本近似，风电、光伏运维成本通常比较低，如果项目不要求特别精确，可以简化成常数项，甚至直接忽略也不影响主结果。
最终目标函数可以统一写成：
min F = 机组发电成本 + 启停成本 + 购售电成本 + 碳交易成本 + 需求响应补偿成本 + 储能运维成本2.2 关键约束“四件套”：功率平衡、机组、储能、需求响应
约束条件才是这类项目最花时间的地方。我习惯分四大块写。
第一块是功率平衡约束。微网或虚拟电厂内部必须时刻满足有功平衡：
∑ P_G(t) + P_PV(t) + P_WT(t) + P_dis(t) - P_ch(t) + P_buy(t) - P_sell(t) = P_load0(t) - P_cut(t) + P_shift_in(t) - P_shift_out(t)等号左边是电源、储能放电、上级电网购电，等号右边是负荷需求。需求响应体现在右边：可中断负荷P_cut降低了该时段需求，可转移负荷转入P_shift_in增加该时段需求，转出P_shift_out减少该时段需求。左右两边同时调整才是真正的功率平衡，如果把削减量直接写到左边，虽然数学上仍可能凑出正确解，但物理含义容易混乱，排查问题也麻烦。
第二块是可控机组约束。包括出力上下限、爬坡约束和最小启停时间。出力上下限一般写为：
0 ≤ P_G(t) ≤ P_Gmax × on(t) P_Gmin × on(t) ≤ P_G(t)这里on(t)是0/1开停机状态。爬坡约束需要跨时段写，例如：
P_G(t) - P_G(t-1) ≤ ramp_up P_G(t-1) - P_G(t) ≤ ramp_down如果忽视了爬坡约束，结果模型也许能求出整数解，但执行层面根本做不到，所以这一步不能省。
第三块是储能约束。储能是日前调度里最重要的灵活调节资源之一，建模时重点关注电量平衡：
E_bat(t+1) = E_bat(t) + η_ch × P_ch(t) × Δt - P_dis(t) × Δt / η_dis其中E_bat是储能当前电量，单位MWh，P_ch、P_dis是充电功率和放电功率，η_ch、η_dis是充放电效率。此外还要限制电量上下限、充放电功率上下限，并让调度周期首末电量相等，这是为了保留储能在次日继续调节的能力。
第四块是需求响应约束。可中断负荷要限制削减量上限和每日最大可中断次数：
0 ≤ P_cut(t) ≤ P_cut_max(t) ∑ on_cut(t) ≤ N_cut_max可转移负荷要保证全天转移前后总用电量不变：
∑ P_shift_in(t) = ∑ P_shift_out(t) 0 ≤ P_shift_in(t) ≤ P_shift_in_max(t) 0 ≤ P_shift_out(t) ≤ P_shift_out_max(t)如果转移补偿单价为正，模型不会傻到同一时段既转进又转出，因为那只会增加成本，不会改善平衡，所以不需要额外加互斥二进制变量。要是你以后去掉转移补偿，只考虑电价套利，就必须补上互斥约束，否则会出现无意义的空转。
2.3 模型属性与求解难度：这本质是一个MILP/MIQP问题
把上面的目标函数和约束放在一起看，二进制变量主要来自机组启停、购售电互斥、可中断负荷计数，某些情况下还包括阶梯碳价的区间选择。连续变量包括机组出力、储能充放电功率、购售电量以及需求响应量。如果机组成本采用二次函数，那么目标函数会带二次项，模型属于混合整数二次规划；如果对二次成本做分段线性化，模型就退化为混合整数线性规划。
论文复现里常见的规模并不大：假设1台燃气轮机、1台柴油机、1套储能、24个时段，决策变量大概几百个，二进制变量几十到上百个，用Gurobi这种商用求解器很快就能出结果。但如果把碳交易和需求响应全部精细化到每个设备、每个分区，变量会快速膨胀，尤其是引入了阶梯碳价的多个0/1变量之后，求解时间会成倍增加。所以我建议初版模型先跑线性碳价、可中断和可转移需求响应，跑通后逐步往里面加复杂度，这样比较容易定位问题。
3. Matlab代码实现：从框架到可运行的调度程序
3.1 工具选择：YALMIP做建模，求解器负责算
Matlab里实现日前调度，最常用的是YALMIP加Gurobi或CPLEX的组合。YALMIP可以理解成一个中间建模层，专门处理sdpvar、binvar、optimize这类建模语法，让代码读起来非常接近数学表达式，极大减少手写线性规划标准矩阵的负担。
打开Matlab后，先确认YALMIP和求解器被正确识别：
% 检查YALMIP是否可用 if exist('yalmip','file') == 0 error('未找到YALMIP，请先下载并设置路径'); end % 查看YALMIP识别到的求解器 yalmiptest环境本身不难配，但要注意YALMIP只是建模语言，不是求解器。没有Gurobi或CPLEX时，很多带二进制的调度模型会跑得很吃力，甚至报错找不到求解器。如果手头没有商用求解器，也可以把目标函数和约束全部线性化，然后直接使用Matlab自带的intlinprog。不过代码可读性和后续扩展性都会差不少，碰到大规模案例时会很痛苦。
工程实现里还有一个容易忽略的点：统一单位。建议全部采用MW、MWh、元和小时四个基础单位。时段步长Δt取1小时时，功率直接乘以时间就是电量，温度、成本和排放约束写起来非常顺。如果单位一会儿用kW、一会儿用MW、一会儿用p.u.，求解器数值尺度过大会出现各种收敛困难的怪问题。
3.2 变量定义与目标函数代码
假设系统里有Ng台可控机组，一个典型YALMIP代码如下：
T = 24; Ng = 2; % 可控机组数量 P_G = sdpvar(Ng, T, 'full'); % 机组出力 on = binvar(Ng, T); % 开停机状态 up = binvar(Ng, T); % 启动标志 dn = binvar(Ng, T); % 停机标志 P_buy = sdpvar(1, T, 'full'); % 上级电网购电功率 P_sell = sdpvar(1, T, 'full'); % 向上级电网售电功率 x_buy = binvar(1, T); % 购电状态 x_sell = binvar(1, T); % 售电状态 P_cut = sdpvar(1, T, 'full'); % 可中断负荷削减量 P_sh_in = sdpvar(1, T, 'full');% 可转移负荷转入 P_sh_out = sdpvar(1, T, 'full');% 可转移负荷转出 E_bat = sdpvar(1, T + 1, 'full'); % 储能电量 P_ch = sdpvar(1, T, 'full'); % 充电功率 P_dis = sdpvar(1, T, 'full'); % 放电功率目标函数可以写得很直接：
Cost_gen = sum(sum(a_g .* P_G.^2 ... + b_g .* P_G ... + c_g .* on)); Cost_startup = sum(sum(start_cost .* up)); Cost_grid = sum(price_buy .* P_buy - price_sell .* P_sell); E_actual = eg_co2 * sum(P_G, 1) ... + egrid_co2 * sum(P_buy); Cost_carbon = co2_price * (E_actual - Allowance); Cost_dr = sum(price_cut .* P_cut ... + price_sh_in .* P_sh_in ... + price_sh_out .* P_sh_out); Objective = Cost_gen + Cost_startup ... + Cost_grid + Cost_carbon + Cost_dr;代码里的a_g、b_g、c_g是机组二次成本系数，start_cost是单次启动成本，price_buy/price_sell是分时购售电价，eg_co2是机组排放强度，egrid_co2是电网平均排放强度，Allowance是当日免费配额，co2_price是碳价。这些都属于基础数据，应当放在调度程序前部一次性定义好，不要把数字硬编码在目标函数里，否则后面做灵敏度分析会很麻烦。
3.3 约束搭建代码
功率平衡约束写在循环里会更直观。假设基础负荷P_load0已经经过价格型需求响应修正：
Constraints = []; for t = 1:T Constraints = [Constraints, ... sum(P_G(:, t)) + P_PV(t) + P_WT(t) ... + P_dis(t) - P_ch(t) ... + P_buy(t) - P_sell(t) ... == P_load0(t) - P_cut(t) ... + P_sh_in(t) - P_sh_out(t)]; end机组状态和启动停机之间需要用约束关联，YALMIP里可以写成：
% 第一时段，不强制启停 Constraints = [Constraints, ... P_G >= 0, ... P_G <= P_gmax .* on, ... P_G >= P_gmin .* on]; % 启动停机和开停机状态的关系 Constraints = [Constraints, ... on(:, t) - on(:, t - 1) == up(:, t) - dn(:, t)]; % 启动成本只需要关心发生了启动的时间点 Constraints = [Constraints, ... up(:, t) >= on(:, t) - on(:, t - 1)];储能电量平衡与充放电上下限：
Constraints = [Constraints, E_bat(1) == E_bat_init]; Constraints = [Constraints, E_bat(T + 1) == E_bat_init]; for t = 1:T Constraints = [Constraints, ... E_bat(t + 1) == E_bat(t) ... + ch_eff * P_ch(t) ... - P_dis(t) / dis_eff]; Constraints = [Constraints, ... E_bat_min <= E_bat(t + 1) <= E_bat_max, ... 0 <= P_ch(t) <= P_ch_max, ... 0 <= P_dis(t) <= P_dis_max]; end需求响应约束里，可转移负荷守恒最容易漏：
Constraints = [Constraints, sum(P_sh_in) == sum(P_sh_out)]; Constraints = [Constraints, 0 <= P_sh_in <= P_sh_in_max]; Constraints = [Constraints, 0 <= P_sh_out <= P_sh_out_max]; Constraints = [Constraints, 0 <= P_cut <= P_cut_max];最后加上购售电互斥约束，调用求解器：
Constraints = [Constraints, ... P_buy <= P_grid_max .* x_buy, ... P_sell <= P_grid_max .* x_sell, ... x_buy + x_sell <= 1]; options = sdpsettings('solver', 'gurobi', ... 'verbose', 2, ... 'gurobi.MIPGap', 0.01); sol = optimize(Constraints, Objective, options); if sol.problem == 0 P_G_opt = value(P_G); P_buy_opt = value(P_buy); ... else disp('求解失败'); disp(sol.info); end跑完后再用value函数把不同变量的结果取
