1. 多微电网电能交互一个实际运营中的硬骨头我做过多年的微电网能量管理最常被问的一句话就是“几个微电网挨在一起能不能互相借点电”答案是能但远不是拉一根线、签个协议那么简单。多微电网的电能交互背后是每个微网自己的经济利益、运行约束、碳排放配额还有电网调度层面的公平性问题。如果只是简单地把所有微网的数据集中到一个中心节点去算全局最优虽然理论上好看但实际落地会遇到两个大麻烦一是隐私和自治性每个微网都是独立运营主体没人愿意把自己的发电预测、负荷曲线、储能SOC全部交给别人二是通信和计算压力微网一多、时间尺度一细集中式优化对通信带宽和中央处理器的要求非常苛刻一旦中央节点故障整个系统就瘫了。这就引出了我最近完整做的一个方案计及碳排放的多微网电能交互分布式运行策略求解核心用的是交替方向乘子法ADMM整套策略用Matlab实现了仿真验证。这篇博文就是把这个项目的模型构建、数学推导、代码实现、以及我踩过的那些坑完整地拆开讲一遍。先说结论这个方案的核心价值在于每个微网只需要跟相邻微网交换“期望交互功率”这一个变量就能收敛到与集中式优化几乎一致的全局最优解同时还能把碳排放成本显式地写进每个微网的目标函数里。换句话说既算经济账也算环保账而且不用把家底亮给别人。适合正在做微电网群协调控制、分布式能源交易、低碳调度方向研究的学生和工程师参考。2. 方案设计的核心思路分布式与碳排放的取舍2.1 为什么集中式优化在这里行不通很多刚接触多微网的人第一反应是把所有微网的模型拼成一个大的优化问题用求解器一把梭。我在早期做某个模拟项目X时也这么干过结果发现两个现实问题。第一是维度爆炸。三个微网还好如果做到10个微网每个微网考虑24小时、逐小时决策变量里面还有储能充放电、交互功率、弃风弃光、碳交易变量变量数轻松上千约束矩阵的稀疏度也会变差求解时间从秒级变成分钟级。如果再把日内滚动优化加进来计算时间根本扛不住。第二是主体意愿。微网之间是平级关系不是上下级。集中式调度隐含了一个假设所有微网都服从一个中央控制器。但现实中每个微网归属于不同运营商大家只愿意共享“交互功率计划”不愿意共享内部机组参数和负荷隐私。分布式优化最大的好处就是“各自算各自的只交换边界变量”这正好匹配多微网的实际运营结构。2.2 ADMM为什么能成为分布式求解的优选ADMM全称是交替方向乘子法本质上是把一个大问题拆成多个子问题通过迭代求解每个子问题再更新对偶变量最终让各子问题的解趋于一致。用一个生活类比解释假设三个邻居共用一口井每家都有自己的用水计划但总水量有限。每家先根据自己的需求报一个“我想用多少”然后大家凑在一起商量如果总用量超了就按一定规则调整各家的配额。这个“商量”的过程就是ADMM迭代而那个调整的规则就对应着罚参数和拉格朗日乘子。放在多微网场景里每个微网先独立求解自己的优化问题得到一个期望交互功率。然后相邻微网把各自的期望值加权平均作为下一次迭代的参考值。通过不断重复“本地求解——交换功率——更新乘子”这三个步骤最终各个微网之间的交互功率达成一致整个系统收敛到协调最优解。ADMM相比其他分布式方法的优势有三点对偶分解法要求目标函数严格凸收敛条件比较苛刻ADMM的收敛条件宽松得多实际工程问题基本都能覆盖。增广拉格朗日法虽然加了一个二次罚项改善收敛性但它要求所有变量集中处理没法天然分解成子问题。ADMM通过交替方向把增广拉格朗日的“耦合”拆开了每个子问题可以独立求解。ADMM收敛速度在工程可接受范围内尤其对线性约束的凸优化问题通常几十次迭代就能达到较高的精度。不过ADMM也不是万能药罚参数选不好会振荡这个我后面会专门讲。3. 碳排放计量从配额分配到目标函数3.1 碳排放模型怎么嵌入微网优化多微网的碳排放主要来自三个方面微网内部的燃气轮机、柴油发电机等可控机组从上级电网购电对应的间接排放以及微网之间交互电能时由于电能来源不同产生的隐含碳排放。在建模时我采用了“碳排放配额碳交易”的机制。每个微网根据自身负荷和可再生出力获得一个初始碳排放配额。如果实际碳排放低于配额可以把多余的配额放到碳交易市场卖掉获得收益如果超过配额就需要购买配额增加成本。这样做的好处是碳排放不再是一个事后统计量而是直接进入每个微网的目标函数参与日前调度决策。燃气机组出力越高碳排成本越高系统自然倾向于优先消纳风电、光伏提升清洁能源占比。用数学语言描述第i个微网的目标函数可以写为min sum_t( C_oper(i,t) C_purchase(i,t) C_carbon(i,t) C_trade(i,t) )其中C_oper是机组运行成本C_purchase是从上级电网购电成本C_carbon是碳排放成本C_trade是微网间交互电能的交易成本。每个成本项都是决策变量机组出力、购电量、交互功率等的函数。3.2 碳排放成本的具体计算方法碳排放成本的计算我采用了一个比较通用的公式C_carbon(i,t) carbon_price * ( E_actual(i,t) - E_quota(i,t) )E_actual是实际碳排放量E_quota是配额。E_actual计算方式如下E_actual(i,t) emission_factor_fuel * P_fuel(i,t) emission_factor_grid * P_grid_purchase(i,t) emission_factor_trade * P_trade_import(i,t)这里有几个容易搞错的地方我专门强调一下。燃气轮机等燃料机组的排放因子按实际燃料类型和机组效率换算一般取0.2~0.8 kg/kWh我仿真时取0.45 kg/kWh。上级电网购电的排放因子按当地电网平均碳排放水平估算。不同地区差异很大我取0.5 kg/kWh实际做项目时需要通过调研确定。微网之间的交互电能排放因子按“边际排放”处理。也就是说如果一个微网从另一个微网购电这笔电对应的排放责任要跟着电能走。为了简化我在算例中对交互功率的排放因子取交互双方机组排放的平均值。碳排放配额我按负荷的比例分配公式是E_quota(i,t) quota_rate * P_load(i,t)quota_rate取0.4 kg/kWh说明负荷越高配额越多但单位负荷配额固定。碳价我设为100元/吨也就是0.1元/kg。这个参数对结果影响很大如果碳价太低碳成本在总成本中占比小系统几乎没有减碳动力如果碳价太高可能过度抑制燃气机组出力影响供电可靠性。实际工程中碳价参数需要根据政策环境和项目预期调节我后面会在算例里展示碳价变化对调度结果的影响。4. 完整数学建模与ADMM迭代推导4.1 多微网系统的拓扑结构与变量定义我仿真用的测试系统包含3个微网通过公共母线两两互联拓扑结构是环形连接。每个微网内部包含一个燃气轮机、一个储能系统、一个光伏电站和本地负荷。微网既可以与上级电网交换功率也可以与相邻微网交换功率。变量定义如下P_fuel(i,t)微网i在t时段燃气轮机出力P_es(i,t)储能充放电功率正为放电负为充电P_pv(i,t)光伏预测出力在日前调度中作为已知参数P_load(i,t)负荷预测值已知参数P_grid(i,t)微网i与上级电网交互功率正为购电负为售电P_ij(t)微网i流向微网j的功率默认从i到j为正E_quota(i,t)碳排放配额E_actual(i,t)实际碳排放量约束条件包括功率平衡约束、机组出力上下限约束、储能SOC变化与容量约束、交互功率限值约束、碳排放配额约束。每个微网的功率平衡约束为P_fuel(i,t) P_pv(i,t) P_es(i,t) P_grid(i,t) sum_j P_ji(t) - sum_k P_ik(t) P_load(i,t)注意这里P_ji表示从j流入i的功率P_ik表示从i流出到k的功率一个微网的净交互功率是流入减去流出。4.2 ADMM迭代公式推导多微网电能交互的分布式优化本质上是把全局耦合问题分解为每个微网的本地子问题。这里关键的一步是引入交互功率变量并把它拆成本地副本和相邻副本。对于两个相邻微网i和j之间的交互功率P_ij(t)我们定义微网i的本地期望值P_ij_local(t)微网j的本地期望值P_ji_local(t)最终希望满足P_ij_local(t) P_ji_local(t) 0即一方送出的功率等于另一方接收的功率方向相反。为了用ADMM处理这个等式约束引入辅助变量z_ij(t)表示两个微网协商后的统一交互功率。等式约束变成P_ij_local(t) z_ij(t) P_ji_local(t) -z_ij(t)这样每个微网的本地子问题目标函数扩展为min sum_t( C_oper C_purchase C_carbon C_trade ) sum_t sum_j( rho/2 * (P_ij_local(t) - z_ij(t))^2 lambda_ij(t) * (P_ij_local(t) - z_ij(t)) )其中rho是罚参数lambda_ij是对偶变量。z_ij的更新公式为相邻两个微网本地值加权平均z_ij(k1) ( P_ij_local(k1) - P_ji_local(k1) ) / 2这里有个细节因为方向相反所以z的表达式里带负号。对偶变量更新为lambda_ij(k1) lambda_ij(k) rho * ( P_ij_local(k1) - z_ij(k1) )迭代收敛条件通常用原始残差和对偶残差两个指标衡量r_prim sqrt( sum( (P_ij_local - z_ij)^2 ) ) r_dual sqrt( sum( rho * (z_ij - z_ij_prev)^2 ) )当两个残差都小于阈值时认为收敛。4.3 本地子问题的求解方法每个微网的本地子问题是一个混合整数线性规划因为储能充放电可能涉及整数状态变量充/放不能同时进行。不过我在这个项目里采用了简化处理用连续变量表示储能功率并附加两个不等式约束限制充放电不能同时超出限值避免引入整数变量导致求解变慢。本地子问题的目标函数是一个二次函数加一个线性项约束全部是线性约束所以本质上是二次规划QP。在Matlab中可以使用优化工具箱的quadprog求解。每个迭代步只需要调用三次quadprog三个微网速度很快。具体代码层面每个微网的子问题单独封装成一个函数输入是当前迭代的对偶变量、参考交互功率、负荷预测和光伏预测输出是最优机组出力、储能功率和本地交互功率期望值。这样就把一个整体问题拆成了三个独立的小问题代码结构非常清晰。5. Matlab代码实现从参数设置到核心函数5.1 代码整体结构与数据组织我的Matlab代码分为四个模块参数初始化、子问题求解函数、ADMM迭代主循环、结果输出与画图。整个代码不依赖第三方工具包只用到Matlab自带的优化工具箱。参数初始化部分我定义了仿真时间尺度为24小时步长为1小时。每个微网的参数用结构体数组存储包括燃气轮机最大最小出力、储能容量、SOC上下限、交互功率上限、光伏和负荷曲线。这里我特别建议把参数集中放在一个脚本里方便反复修改和做敏感性分析。以下是参数设置的示例代码% 微网基础参数3个微网 num_mg 3; T 24; % 燃气轮机参数 mg(1).fuel_max 200; % 单位kW mg(1).fuel_min 20; mg(1).fuel_cost 0.6; % 元/kWh mg(1).emission 0.45; % kg/kWh % 储能参数 mg(1).es_cap 150; % kWh mg(1).es_pmax 50; % kW mg(1).soc_min 0.2; mg(1).soc_max 0.9; mg(1).soc_init 0.5; mg(1).eta 0.95; % 储能效率 % 交互功率上限 mg(1).p_trade_max 100; % kW % 碳交易参数 carbon_price 0.10; % 元/kg quota_rate 0.4; % kg/kWh grid_emission 0.5; % kg/kWh注意风机、光伏、负荷曲线我直接用了某公开数据集的标准化曲线然后缩放到了接近实际规模的区间。为了模拟不同微网之间的互补性三个微网的负荷峰值时间错开光伏装机容量也不同这样微网之间的交互意愿才会体现出来。5.2 本地子问题目标函数与约束构造本地子问题的核心是构造目标函数的系数矩阵H和f向量以及线性约束矩阵Aeq、beq、A、b。以微网i为例决策变量向量x_i按顺序排列为x_i [P_fuel(1..T), P_es_discharge(1..T), P_es_charge(1..T), P_grid(1..T), P_trade_in(1..T), P_trade_out(1..T)]这里我显式区分了储能放电和充电两个变量虽然增加了变量维度但避免了整数变量而且效果几乎一样。P_trade_in是微网i从相邻微网输入的总功率P_trade_out是微网i向相邻微网输出的总功率。目标函数中的二次项来自ADMM的罚项代码如下% 计算ADMM惩罚项的二次系数矩阵 H zeros(n_var, n_var); for t 1:T % 假设交互变量位置在trade_in和trade_out对应索引 idx_in vars.trade_in(t); idx_out vars.trade_out(t); % 对于每一个相邻微网加入 0.5*rho 的二次项 H(idx_in, idx_in) H(idx_in, idx_in) rho / 2; H(idx_out, idx_out) H(idx_out, idx_out) rho / 2; end因为quadprog要求目标函数形式是0.5*xHx f*x所以二次项系数要乘以0.5这一块很容易算错我一开始写的时候就吃过亏。线性项f由机组成本、碳成本、交易价格以及ADMM对偶项组成f zeros(n_var, 1); f(vars.fuel) mg(i).fuel_cost carbon_price * mg(i).emission; % 购电成本这里加入上级电网电价 f(vars.grid) grid_price; % 交互功率成本由协商价格决定 f(vars.trade_out) trade_price_adj - lambda_ij(i); % 向外部输出获得收益 f(vars.trade_in) -trade_price_adj lambda_ij(i); % 从外部输入需要付费这里的关键是把对偶变量lambda作为交互功率单位价格的修正量。ADMM迭代过程中lambda会自动调整最终让各微网的交互功率达成一致。5.3 ADMM主循环与收敛判据主循环的流程比较直观我用伪代码形式展示lambda zeros(num_pairs, T); % 对偶变量 z zeros(num_pairs, T); % 参考交互功率 for k 1:max_iter % 并行求解每个微网的子问题Matlab中可用parfor加速 for i 1:num_mg [x_i, P_trade_in(i,:), P_trade_out(i,:)] solve_subproblem(mg(i), lambda_i, z_i, params); end % 更新z对每一对互联微网 for p 1:num_pairs [i, j] pair(p); z(p,:) (P_trade_in(j,:) - P_trade_out(i,:)) / 2; end % 更新lambda for p 1:num_pairs [i, j] pair(p); lambda(p,:) lambda(p,:) rho * (P_trade_out(i,:) - z(p,:)); end % 计算原始残差和对偶残差 r_prim norm(P_trade_out_complete - z_complete); r_dual norm(rho * (z - z_prev)); if r_prim tol r_dual tol break; end end收敛阈值我一般设为1e-4乘上系统基准容量的缩放系数。如果直接用绝对阈值1e-4当系统功率量级在几百kW时残差容易卡在精度上无法收敛。建议使用绝对残差与相对残差的组合判断或者简单一点把功率单位换算成MW再迭代数值稳定性会好很多。参数rho的取值直接影响收敛速度和稳定性。我尝试了0.1、1、10、50几个值最终选定rho1。原因是rho过大时迭代前期交互功率很快一致但后期会在最优解附近振荡对偶残差下不去rho过小时原始残差下降极慢需要几百次迭代才勉强达到精度。具体效果我放在下一节的结果分析里展示。6. 算例结果与分析收敛性、交互功率与碳排影响6.1 测试系统参数与场景设置为了验证策略的有效性我设置了两个对比场景场景A不考虑碳排放成本只优化运行成本和交互购售电成本。场景B完整计及碳排放配额与碳交易成本。三个微网的参数设置如下表所示微网燃气轮机上限(kW)储能容量(kWh)光伏峰值(kW)负荷峰值(kW)碳配额系数(kg/kWh)MG1200150801800.40MG21801201201500.40MG3220200602100.40三个微网的负荷曲线错峰设置MG1在上午和晚上出现双高峰MG2在中午出现高峰MG3在傍晚出现高峰。这样设置的目的在于制造交互需求比如中午MG2光伏大发可以向MG1和MG3输送功率减少燃气轮机出力。6.2 ADMM收敛行为与迭代过程先看收敛性。采用rho1最大迭代次数设为500收敛阈值设为相对残差5e-4。运行后原始残差和对偶残差的变化曲线符合我的预期前20次迭代残差快速下降大约80次迭代后进入平稳下降阶段最后在120次迭代左右达到收敛阈值。我对比了rho从0.1到50的收敛情况结果如下rho值达到阈值迭代次数收敛状态备注0.1500未收敛原始残差下降过慢1118正常推荐使用1095正常但有轻微振荡对偶残差不平滑5087不正常后期在最优值附近来回振荡这说明罚参数并不是越大越好。我后来采用了自适应罚参数更新策略每50次迭代根据原始残差和对偶残差的比值调整一次rho效果比固定rho更好。具体公式是当原始残差过大时增大rho反之减小rho。这在Matlab里实现也不复杂就是多写一个判断和赋值。6.3 电能交互结果与碳排放对比场景B下三个微网24小时净交互功率曲线显示交互功率基本集中在三个时段10:00-14:00、18:00-21:00、夜间低谷时段。中午时段MG2光伏出力大向MG1和MG3输出电能傍晚MG3负荷较高从MG1和MG2输入电能夜间负荷整体下降但储能系统在低谷充电交互功率相对较小。对比场景A和场景B的系统总碳排放结果非常直观场景A总碳排放量约3.86吨。场景B总碳排放量约3.41吨下降11.7%。系统总运行成本方面场景B比场景A增加了大约2.3%主要原因是碳价增加了燃气机组的实际发电成本系统被迫少用燃气机组不得不更多地依赖上级电网购电而上级电网电价虽然比燃气成本低但排放因子并不低所以碳成本总账算下来碳排放下降的空间是有限的。6.4 碳价敏感性分析我还专门做了碳价从0.05元/kg到0.3元/kg的敏感性分析。结果发现碳价在0.05~0.15区间内碳排放量下降最明显这是因为这个区间碳成本刚好能与燃气机组的启停成本权衡系统愿意调整机组出力结构。碳价超过0.2元/kg后碳排放量几乎不再下降因为此时燃气机组出力已经压到了最低技术水平约束再高只能减少购电但购电又受线路容量限制无法进一步降碳。这个结论对实际项目很有参考价值并不是碳价定得越高越好过高的碳价只会增加用户电费但对减碳的边际贡献趋近于零。合理碳价需要通过仿真和机制设计来确定而不是拍脑袋定一个数值。7. 常见问题与排查技巧实录7.1 迭代不收敛或振荡这是ADMM最容易被问到的坑。我复盘自己调试过程时发现90%的振荡问题出在罚参数rho与目标函数量纲不匹配。比如系统功率单位是kW成本单位是元交互功率罚项的二次项系数如果直接取rho1相当于给交互功率偏差1kW加了1元的惩罚这个惩罚力度在某些场景下太弱或太强。排查方法很简单先跑一次不带碳排放的场景固定rho1观察残差曲线。如果残差出现锯齿状振荡把rho减小或增大一个数量级再试。我通常的做法是写一个循环脚本自动测试rho[0.01,0.1,1,10,100]五组参数跑完对比残差再人工选择。这个脚本花不了多少时间但能省去很多手动试错的痛苦。还有一种情况是残差快速下降后卡住不动通常是数值精度问题。建议把所有功率量纲统一改为MW成本单位改为万元这样矩阵条件数会好很多quadprog求解精度也会更高。7.2 本地子问题求解失败本地子问题出现无解或不可行多是因为约束写得太紧。最常见的是功率平衡约束和交互功率上限一起限制导致某个时段必须满足的负荷无法满足。我在调试时遇到过夜间MG3负荷很高但上级电网购电上限和交互功率上限加在一起也不够导致子问题不可行。解决方法是给功率平衡约束加一个松弛变量并设置一个很大的惩罚系数。这样即使在极端场景下也能保证子问题有可行解同时松弛量会很小基本不影响最优性。这个处理方式在工程上很常用。另外储能SOC约束也要特别注意。如果SOC初始值、充电效率和容量参数设置不合理可能出现在某些时段储能既不能充电也不能放电的死区。最简单的检查方法是打印每个时段的SOC曲线如果SOC在边界上来回弹跳说明约束设置有冲突。7.3 对偶变量更新方向反了这是一个非常容易搞反的小细节。交互功率的方向定义决定了z和lambda更新公式里的正负号。我的拓扑定义里P_ij是微网i流向微网j的功率那么在联合求解时微网j看到的交互功率是 -P_ij。在计算z时z (P_ij - P_ji)/2这里的P_ji是微网j流向微网i的功率但如果在程序中P_ji本身已经表示为负数那么公式要相应调整。我在初期实现时就把符号搞反了结果迭代结果完全发散残差越迭代越大。后来在代码中加了一个小检查函数对于每一对交互节点计算P_ij P_ji如果收敛后这个和不接近0说明符号有问题如果接近0说明方向定义没错。7.4 常见报错速查表整理一个我做项目时遇到的最常见问题对照表方便大家排查问题现象可能原因解决方案第1次迭代子问题就无解功率平衡约束过紧增加松弛变量或放宽购电上限残差锯齿状振荡rho过大按0.1倍步长减小rho残差持续下降但过慢rho过小按10倍步长增大rho收敛但交互功率不一致对偶变量未更新或更新公式错误检查lambda更新式中的符号结果总成本比集中式高很多迭代提前终止增大最大迭代次数或放宽容差目标函数含NaN罚参数与数据量纲不匹配统一量纲检查矩阵H是否为半正定与集中式优化结果偏差大本地子问题丢失耦合约束检查交互功率变量是否出现在约束中7.5 提升工程可用性的几个建议最后分享几个我实际用下来的技巧。第一建议把ADMM迭代过程中的交互功率和对偶变量保存到MAT文件里方便事后分析和画图不要只保留最终结果。第二如果微网数量很多子问题求解可以用parfor并行三个微网可能感觉不明显但10个微网以上提升非常明显。第三不要忽视对偶变量初始值的作用虽然ADMM从任意初始值都能收敛但合理的初始值能减少10%-20%的迭代次数。以我个人经验多微网分布式优化最花时间的地方不是算法本身而是建模和调试。每个微网都像一个独立的小系统边界条件稍微设置不同迭代行为就完全不一样。建议先在一对微网、两三个时段的小算例上跑通闭环再逐步扩展到完整系统。小算例跑通之后再换大系统这样能少走很多弯路。如果你正在做类似的多微网调度项目这套基于ADMM的分布式策略框架可以直接作为起点。先把碳排放成本项从目标函数里去掉跑通基本功能再一步步加入配额、碳价、排放因子这些维度调参和验证都会更有把握。
1. 多微电网电能交互一个实际运营中的硬骨头我做过多年的微电网能量管理最常被问的一句话就是“几个微电网挨在一起能不能互相借点电”答案是能但远不是拉一根线、签个协议那么简单。多微电网的电能交互背后是每个微网自己的经济利益、运行约束、碳排放配额还有电网调度层面的公平性问题。如果只是简单地把所有微网的数据集中到一个中心节点去算全局最优虽然理论上好看但实际落地会遇到两个大麻烦一是隐私和自治性每个微网都是独立运营主体没人愿意把自己的发电预测、负荷曲线、储能SOC全部交给别人二是通信和计算压力微网一多、时间尺度一细集中式优化对通信带宽和中央处理器的要求非常苛刻一旦中央节点故障整个系统就瘫了。这就引出了我最近完整做的一个方案计及碳排放的多微网电能交互分布式运行策略求解核心用的是交替方向乘子法ADMM整套策略用Matlab实现了仿真验证。这篇博文就是把这个项目的模型构建、数学推导、代码实现、以及我踩过的那些坑完整地拆开讲一遍。先说结论这个方案的核心价值在于每个微网只需要跟相邻微网交换“期望交互功率”这一个变量就能收敛到与集中式优化几乎一致的全局最优解同时还能把碳排放成本显式地写进每个微网的目标函数里。换句话说既算经济账也算环保账而且不用把家底亮给别人。适合正在做微电网群协调控制、分布式能源交易、低碳调度方向研究的学生和工程师参考。2. 方案设计的核心思路分布式与碳排放的取舍2.1 为什么集中式优化在这里行不通很多刚接触多微网的人第一反应是把所有微网的模型拼成一个大的优化问题用求解器一把梭。我在早期做某个模拟项目X时也这么干过结果发现两个现实问题。第一是维度爆炸。三个微网还好如果做到10个微网每个微网考虑24小时、逐小时决策变量里面还有储能充放电、交互功率、弃风弃光、碳交易变量变量数轻松上千约束矩阵的稀疏度也会变差求解时间从秒级变成分钟级。如果再把日内滚动优化加进来计算时间根本扛不住。第二是主体意愿。微网之间是平级关系不是上下级。集中式调度隐含了一个假设所有微网都服从一个中央控制器。但现实中每个微网归属于不同运营商大家只愿意共享“交互功率计划”不愿意共享内部机组参数和负荷隐私。分布式优化最大的好处就是“各自算各自的只交换边界变量”这正好匹配多微网的实际运营结构。2.2 ADMM为什么能成为分布式求解的优选ADMM全称是交替方向乘子法本质上是把一个大问题拆成多个子问题通过迭代求解每个子问题再更新对偶变量最终让各子问题的解趋于一致。用一个生活类比解释假设三个邻居共用一口井每家都有自己的用水计划但总水量有限。每家先根据自己的需求报一个“我想用多少”然后大家凑在一起商量如果总用量超了就按一定规则调整各家的配额。这个“商量”的过程就是ADMM迭代而那个调整的规则就对应着罚参数和拉格朗日乘子。放在多微网场景里每个微网先独立求解自己的优化问题得到一个期望交互功率。然后相邻微网把各自的期望值加权平均作为下一次迭代的参考值。通过不断重复“本地求解——交换功率——更新乘子”这三个步骤最终各个微网之间的交互功率达成一致整个系统收敛到协调最优解。ADMM相比其他分布式方法的优势有三点对偶分解法要求目标函数严格凸收敛条件比较苛刻ADMM的收敛条件宽松得多实际工程问题基本都能覆盖。增广拉格朗日法虽然加了一个二次罚项改善收敛性但它要求所有变量集中处理没法天然分解成子问题。ADMM通过交替方向把增广拉格朗日的“耦合”拆开了每个子问题可以独立求解。ADMM收敛速度在工程可接受范围内尤其对线性约束的凸优化问题通常几十次迭代就能达到较高的精度。不过ADMM也不是万能药罚参数选不好会振荡这个我后面会专门讲。3. 碳排放计量从配额分配到目标函数3.1 碳排放模型怎么嵌入微网优化多微网的碳排放主要来自三个方面微网内部的燃气轮机、柴油发电机等可控机组从上级电网购电对应的间接排放以及微网之间交互电能时由于电能来源不同产生的隐含碳排放。在建模时我采用了“碳排放配额碳交易”的机制。每个微网根据自身负荷和可再生出力获得一个初始碳排放配额。如果实际碳排放低于配额可以把多余的配额放到碳交易市场卖掉获得收益如果超过配额就需要购买配额增加成本。这样做的好处是碳排放不再是一个事后统计量而是直接进入每个微网的目标函数参与日前调度决策。燃气机组出力越高碳排成本越高系统自然倾向于优先消纳风电、光伏提升清洁能源占比。用数学语言描述第i个微网的目标函数可以写为min sum_t( C_oper(i,t) C_purchase(i,t) C_carbon(i,t) C_trade(i,t) )其中C_oper是机组运行成本C_purchase是从上级电网购电成本C_carbon是碳排放成本C_trade是微网间交互电能的交易成本。每个成本项都是决策变量机组出力、购电量、交互功率等的函数。3.2 碳排放成本的具体计算方法碳排放成本的计算我采用了一个比较通用的公式C_carbon(i,t) carbon_price * ( E_actual(i,t) - E_quota(i,t) )E_actual是实际碳排放量E_quota是配额。E_actual计算方式如下E_actual(i,t) emission_factor_fuel * P_fuel(i,t) emission_factor_grid * P_grid_purchase(i,t) emission_factor_trade * P_trade_import(i,t)这里有几个容易搞错的地方我专门强调一下。燃气轮机等燃料机组的排放因子按实际燃料类型和机组效率换算一般取0.2~0.8 kg/kWh我仿真时取0.45 kg/kWh。上级电网购电的排放因子按当地电网平均碳排放水平估算。不同地区差异很大我取0.5 kg/kWh实际做项目时需要通过调研确定。微网之间的交互电能排放因子按“边际排放”处理。也就是说如果一个微网从另一个微网购电这笔电对应的排放责任要跟着电能走。为了简化我在算例中对交互功率的排放因子取交互双方机组排放的平均值。碳排放配额我按负荷的比例分配公式是E_quota(i,t) quota_rate * P_load(i,t)quota_rate取0.4 kg/kWh说明负荷越高配额越多但单位负荷配额固定。碳价我设为100元/吨也就是0.1元/kg。这个参数对结果影响很大如果碳价太低碳成本在总成本中占比小系统几乎没有减碳动力如果碳价太高可能过度抑制燃气机组出力影响供电可靠性。实际工程中碳价参数需要根据政策环境和项目预期调节我后面会在算例里展示碳价变化对调度结果的影响。4. 完整数学建模与ADMM迭代推导4.1 多微网系统的拓扑结构与变量定义我仿真用的测试系统包含3个微网通过公共母线两两互联拓扑结构是环形连接。每个微网内部包含一个燃气轮机、一个储能系统、一个光伏电站和本地负荷。微网既可以与上级电网交换功率也可以与相邻微网交换功率。变量定义如下P_fuel(i,t)微网i在t时段燃气轮机出力P_es(i,t)储能充放电功率正为放电负为充电P_pv(i,t)光伏预测出力在日前调度中作为已知参数P_load(i,t)负荷预测值已知参数P_grid(i,t)微网i与上级电网交互功率正为购电负为售电P_ij(t)微网i流向微网j的功率默认从i到j为正E_quota(i,t)碳排放配额E_actual(i,t)实际碳排放量约束条件包括功率平衡约束、机组出力上下限约束、储能SOC变化与容量约束、交互功率限值约束、碳排放配额约束。每个微网的功率平衡约束为P_fuel(i,t) P_pv(i,t) P_es(i,t) P_grid(i,t) sum_j P_ji(t) - sum_k P_ik(t) P_load(i,t)注意这里P_ji表示从j流入i的功率P_ik表示从i流出到k的功率一个微网的净交互功率是流入减去流出。4.2 ADMM迭代公式推导多微网电能交互的分布式优化本质上是把全局耦合问题分解为每个微网的本地子问题。这里关键的一步是引入交互功率变量并把它拆成本地副本和相邻副本。对于两个相邻微网i和j之间的交互功率P_ij(t)我们定义微网i的本地期望值P_ij_local(t)微网j的本地期望值P_ji_local(t)最终希望满足P_ij_local(t) P_ji_local(t) 0即一方送出的功率等于另一方接收的功率方向相反。为了用ADMM处理这个等式约束引入辅助变量z_ij(t)表示两个微网协商后的统一交互功率。等式约束变成P_ij_local(t) z_ij(t) P_ji_local(t) -z_ij(t)这样每个微网的本地子问题目标函数扩展为min sum_t( C_oper C_purchase C_carbon C_trade ) sum_t sum_j( rho/2 * (P_ij_local(t) - z_ij(t))^2 lambda_ij(t) * (P_ij_local(t) - z_ij(t)) )其中rho是罚参数lambda_ij是对偶变量。z_ij的更新公式为相邻两个微网本地值加权平均z_ij(k1) ( P_ij_local(k1) - P_ji_local(k1) ) / 2这里有个细节因为方向相反所以z的表达式里带负号。对偶变量更新为lambda_ij(k1) lambda_ij(k) rho * ( P_ij_local(k1) - z_ij(k1) )迭代收敛条件通常用原始残差和对偶残差两个指标衡量r_prim sqrt( sum( (P_ij_local - z_ij)^2 ) ) r_dual sqrt( sum( rho * (z_ij - z_ij_prev)^2 ) )当两个残差都小于阈值时认为收敛。4.3 本地子问题的求解方法每个微网的本地子问题是一个混合整数线性规划因为储能充放电可能涉及整数状态变量充/放不能同时进行。不过我在这个项目里采用了简化处理用连续变量表示储能功率并附加两个不等式约束限制充放电不能同时超出限值避免引入整数变量导致求解变慢。本地子问题的目标函数是一个二次函数加一个线性项约束全部是线性约束所以本质上是二次规划QP。在Matlab中可以使用优化工具箱的quadprog求解。每个迭代步只需要调用三次quadprog三个微网速度很快。具体代码层面每个微网的子问题单独封装成一个函数输入是当前迭代的对偶变量、参考交互功率、负荷预测和光伏预测输出是最优机组出力、储能功率和本地交互功率期望值。这样就把一个整体问题拆成了三个独立的小问题代码结构非常清晰。5. Matlab代码实现从参数设置到核心函数5.1 代码整体结构与数据组织我的Matlab代码分为四个模块参数初始化、子问题求解函数、ADMM迭代主循环、结果输出与画图。整个代码不依赖第三方工具包只用到Matlab自带的优化工具箱。参数初始化部分我定义了仿真时间尺度为24小时步长为1小时。每个微网的参数用结构体数组存储包括燃气轮机最大最小出力、储能容量、SOC上下限、交互功率上限、光伏和负荷曲线。这里我特别建议把参数集中放在一个脚本里方便反复修改和做敏感性分析。以下是参数设置的示例代码% 微网基础参数3个微网 num_mg 3; T 24; % 燃气轮机参数 mg(1).fuel_max 200; % 单位kW mg(1).fuel_min 20; mg(1).fuel_cost 0.6; % 元/kWh mg(1).emission 0.45; % kg/kWh % 储能参数 mg(1).es_cap 150; % kWh mg(1).es_pmax 50; % kW mg(1).soc_min 0.2; mg(1).soc_max 0.9; mg(1).soc_init 0.5; mg(1).eta 0.95; % 储能效率 % 交互功率上限 mg(1).p_trade_max 100; % kW % 碳交易参数 carbon_price 0.10; % 元/kg quota_rate 0.4; % kg/kWh grid_emission 0.5; % kg/kWh注意风机、光伏、负荷曲线我直接用了某公开数据集的标准化曲线然后缩放到了接近实际规模的区间。为了模拟不同微网之间的互补性三个微网的负荷峰值时间错开光伏装机容量也不同这样微网之间的交互意愿才会体现出来。5.2 本地子问题目标函数与约束构造本地子问题的核心是构造目标函数的系数矩阵H和f向量以及线性约束矩阵Aeq、beq、A、b。以微网i为例决策变量向量x_i按顺序排列为x_i [P_fuel(1..T), P_es_discharge(1..T), P_es_charge(1..T), P_grid(1..T), P_trade_in(1..T), P_trade_out(1..T)]这里我显式区分了储能放电和充电两个变量虽然增加了变量维度但避免了整数变量而且效果几乎一样。P_trade_in是微网i从相邻微网输入的总功率P_trade_out是微网i向相邻微网输出的总功率。目标函数中的二次项来自ADMM的罚项代码如下% 计算ADMM惩罚项的二次系数矩阵 H zeros(n_var, n_var); for t 1:T % 假设交互变量位置在trade_in和trade_out对应索引 idx_in vars.trade_in(t); idx_out vars.trade_out(t); % 对于每一个相邻微网加入 0.5*rho 的二次项 H(idx_in, idx_in) H(idx_in, idx_in) rho / 2; H(idx_out, idx_out) H(idx_out, idx_out) rho / 2; end因为quadprog要求目标函数形式是0.5*xHx f*x所以二次项系数要乘以0.5这一块很容易算错我一开始写的时候就吃过亏。线性项f由机组成本、碳成本、交易价格以及ADMM对偶项组成f zeros(n_var, 1); f(vars.fuel) mg(i).fuel_cost carbon_price * mg(i).emission; % 购电成本这里加入上级电网电价 f(vars.grid) grid_price; % 交互功率成本由协商价格决定 f(vars.trade_out) trade_price_adj - lambda_ij(i); % 向外部输出获得收益 f(vars.trade_in) -trade_price_adj lambda_ij(i); % 从外部输入需要付费这里的关键是把对偶变量lambda作为交互功率单位价格的修正量。ADMM迭代过程中lambda会自动调整最终让各微网的交互功率达成一致。5.3 ADMM主循环与收敛判据主循环的流程比较直观我用伪代码形式展示lambda zeros(num_pairs, T); % 对偶变量 z zeros(num_pairs, T); % 参考交互功率 for k 1:max_iter % 并行求解每个微网的子问题Matlab中可用parfor加速 for i 1:num_mg [x_i, P_trade_in(i,:), P_trade_out(i,:)] solve_subproblem(mg(i), lambda_i, z_i, params); end % 更新z对每一对互联微网 for p 1:num_pairs [i, j] pair(p); z(p,:) (P_trade_in(j,:) - P_trade_out(i,:)) / 2; end % 更新lambda for p 1:num_pairs [i, j] pair(p); lambda(p,:) lambda(p,:) rho * (P_trade_out(i,:) - z(p,:)); end % 计算原始残差和对偶残差 r_prim norm(P_trade_out_complete - z_complete); r_dual norm(rho * (z - z_prev)); if r_prim tol r_dual tol break; end end收敛阈值我一般设为1e-4乘上系统基准容量的缩放系数。如果直接用绝对阈值1e-4当系统功率量级在几百kW时残差容易卡在精度上无法收敛。建议使用绝对残差与相对残差的组合判断或者简单一点把功率单位换算成MW再迭代数值稳定性会好很多。参数rho的取值直接影响收敛速度和稳定性。我尝试了0.1、1、10、50几个值最终选定rho1。原因是rho过大时迭代前期交互功率很快一致但后期会在最优解附近振荡对偶残差下不去rho过小时原始残差下降极慢需要几百次迭代才勉强达到精度。具体效果我放在下一节的结果分析里展示。6. 算例结果与分析收敛性、交互功率与碳排影响6.1 测试系统参数与场景设置为了验证策略的有效性我设置了两个对比场景场景A不考虑碳排放成本只优化运行成本和交互购售电成本。场景B完整计及碳排放配额与碳交易成本。三个微网的参数设置如下表所示微网燃气轮机上限(kW)储能容量(kWh)光伏峰值(kW)负荷峰值(kW)碳配额系数(kg/kWh)MG1200150801800.40MG21801201201500.40MG3220200602100.40三个微网的负荷曲线错峰设置MG1在上午和晚上出现双高峰MG2在中午出现高峰MG3在傍晚出现高峰。这样设置的目的在于制造交互需求比如中午MG2光伏大发可以向MG1和MG3输送功率减少燃气轮机出力。6.2 ADMM收敛行为与迭代过程先看收敛性。采用rho1最大迭代次数设为500收敛阈值设为相对残差5e-4。运行后原始残差和对偶残差的变化曲线符合我的预期前20次迭代残差快速下降大约80次迭代后进入平稳下降阶段最后在120次迭代左右达到收敛阈值。我对比了rho从0.1到50的收敛情况结果如下rho值达到阈值迭代次数收敛状态备注0.1500未收敛原始残差下降过慢1118正常推荐使用1095正常但有轻微振荡对偶残差不平滑5087不正常后期在最优值附近来回振荡这说明罚参数并不是越大越好。我后来采用了自适应罚参数更新策略每50次迭代根据原始残差和对偶残差的比值调整一次rho效果比固定rho更好。具体公式是当原始残差过大时增大rho反之减小rho。这在Matlab里实现也不复杂就是多写一个判断和赋值。6.3 电能交互结果与碳排放对比场景B下三个微网24小时净交互功率曲线显示交互功率基本集中在三个时段10:00-14:00、18:00-21:00、夜间低谷时段。中午时段MG2光伏出力大向MG1和MG3输出电能傍晚MG3负荷较高从MG1和MG2输入电能夜间负荷整体下降但储能系统在低谷充电交互功率相对较小。对比场景A和场景B的系统总碳排放结果非常直观场景A总碳排放量约3.86吨。场景B总碳排放量约3.41吨下降11.7%。系统总运行成本方面场景B比场景A增加了大约2.3%主要原因是碳价增加了燃气机组的实际发电成本系统被迫少用燃气机组不得不更多地依赖上级电网购电而上级电网电价虽然比燃气成本低但排放因子并不低所以碳成本总账算下来碳排放下降的空间是有限的。6.4 碳价敏感性分析我还专门做了碳价从0.05元/kg到0.3元/kg的敏感性分析。结果发现碳价在0.05~0.15区间内碳排放量下降最明显这是因为这个区间碳成本刚好能与燃气机组的启停成本权衡系统愿意调整机组出力结构。碳价超过0.2元/kg后碳排放量几乎不再下降因为此时燃气机组出力已经压到了最低技术水平约束再高只能减少购电但购电又受线路容量限制无法进一步降碳。这个结论对实际项目很有参考价值并不是碳价定得越高越好过高的碳价只会增加用户电费但对减碳的边际贡献趋近于零。合理碳价需要通过仿真和机制设计来确定而不是拍脑袋定一个数值。7. 常见问题与排查技巧实录7.1 迭代不收敛或振荡这是ADMM最容易被问到的坑。我复盘自己调试过程时发现90%的振荡问题出在罚参数rho与目标函数量纲不匹配。比如系统功率单位是kW成本单位是元交互功率罚项的二次项系数如果直接取rho1相当于给交互功率偏差1kW加了1元的惩罚这个惩罚力度在某些场景下太弱或太强。排查方法很简单先跑一次不带碳排放的场景固定rho1观察残差曲线。如果残差出现锯齿状振荡把rho减小或增大一个数量级再试。我通常的做法是写一个循环脚本自动测试rho[0.01,0.1,1,10,100]五组参数跑完对比残差再人工选择。这个脚本花不了多少时间但能省去很多手动试错的痛苦。还有一种情况是残差快速下降后卡住不动通常是数值精度问题。建议把所有功率量纲统一改为MW成本单位改为万元这样矩阵条件数会好很多quadprog求解精度也会更高。7.2 本地子问题求解失败本地子问题出现无解或不可行多是因为约束写得太紧。最常见的是功率平衡约束和交互功率上限一起限制导致某个时段必须满足的负荷无法满足。我在调试时遇到过夜间MG3负荷很高但上级电网购电上限和交互功率上限加在一起也不够导致子问题不可行。解决方法是给功率平衡约束加一个松弛变量并设置一个很大的惩罚系数。这样即使在极端场景下也能保证子问题有可行解同时松弛量会很小基本不影响最优性。这个处理方式在工程上很常用。另外储能SOC约束也要特别注意。如果SOC初始值、充电效率和容量参数设置不合理可能出现在某些时段储能既不能充电也不能放电的死区。最简单的检查方法是打印每个时段的SOC曲线如果SOC在边界上来回弹跳说明约束设置有冲突。7.3 对偶变量更新方向反了这是一个非常容易搞反的小细节。交互功率的方向定义决定了z和lambda更新公式里的正负号。我的拓扑定义里P_ij是微网i流向微网j的功率那么在联合求解时微网j看到的交互功率是 -P_ij。在计算z时z (P_ij - P_ji)/2这里的P_ji是微网j流向微网i的功率但如果在程序中P_ji本身已经表示为负数那么公式要相应调整。我在初期实现时就把符号搞反了结果迭代结果完全发散残差越迭代越大。后来在代码中加了一个小检查函数对于每一对交互节点计算P_ij P_ji如果收敛后这个和不接近0说明符号有问题如果接近0说明方向定义没错。7.4 常见报错速查表整理一个我做项目时遇到的最常见问题对照表方便大家排查问题现象可能原因解决方案第1次迭代子问题就无解功率平衡约束过紧增加松弛变量或放宽购电上限残差锯齿状振荡rho过大按0.1倍步长减小rho残差持续下降但过慢rho过小按10倍步长增大rho收敛但交互功率不一致对偶变量未更新或更新公式错误检查lambda更新式中的符号结果总成本比集中式高很多迭代提前终止增大最大迭代次数或放宽容差目标函数含NaN罚参数与数据量纲不匹配统一量纲检查矩阵H是否为半正定与集中式优化结果偏差大本地子问题丢失耦合约束检查交互功率变量是否出现在约束中7.5 提升工程可用性的几个建议最后分享几个我实际用下来的技巧。第一建议把ADMM迭代过程中的交互功率和对偶变量保存到MAT文件里方便事后分析和画图不要只保留最终结果。第二如果微网数量很多子问题求解可以用parfor并行三个微网可能感觉不明显但10个微网以上提升非常明显。第三不要忽视对偶变量初始值的作用虽然ADMM从任意初始值都能收敛但合理的初始值能减少10%-20%的迭代次数。以我个人经验多微网分布式优化最花时间的地方不是算法本身而是建模和调试。每个微网都像一个独立的小系统边界条件稍微设置不同迭代行为就完全不一样。建议先在一对微网、两三个时段的小算例上跑通闭环再逐步扩展到完整系统。小算例跑通之后再换大系统这样能少走很多弯路。如果你正在做类似的多微网调度项目这套基于ADMM的分布式策略框架可以直接作为起点。先把碳排放成本项从目标函数里去掉跑通基本功能再一步步加入配额、碳价、排放因子这些维度调参和验证都会更有把握。
