最近在梳理新能源电网优化方向的课题把一套“电气互联系统有功-无功协同优化模型”从头到尾实现了一遍代码用的是Matlab优化框架走的是YalmipCplex/Gurobi这条路。做完之后我最大的感受是这个项目不是单纯地把电力潮流计算和无功优化程序拼在一起而是要在两个异质网络的耦合边界上做文章同时还要把“碳中和”这个约束落到目标函数里。这篇内容我不打算念教科书直接讲思路、讲建模、讲代码实现把踩过的坑一并整理出来。整个模型做的事情可以用一句话概括在电力网络和天然气网络共同构成的多能流系统里同时优化有功出力和无功控制变量在满足潮流、气压、碳排放等约束的前提下让系统总运行成本含碳排放成本最低。听起来复杂拆开看就是“物理建模 优化建模 编程实现”三段式。我会把这几个环节逐一展开文末再给一版通用的问题排查清单想复现这个方向的同学可以直接对照使用。1. 项目思路拆解这个模型到底在优化什么1.1 碳中和背景下为什么不能只盯着调度侧的有功传统电力系统经济调度关注的核心是“有功怎么分配最省钱”机组爬坡、启停、燃料成本这些都围绕有功展开。但到了高比例新能源接入的场景情况明显变了。风电、光伏本身出力波动大而且大多通过电力电子变流器并网不像传统同步机组那样天然具备无功支撑能力。结果就是系统电压越限的频率明显上升线路损耗增大甚至影响设备安全。无功功率不直接对外做功但它决定了电压水平而电压水平又直接决定了输电网络的传输效率和稳定性。碳中和目标下新能源占比继续提升送端和受端电网都面临更大的电压控制压力。在这种背景下无功优化从过去的“锦上添花”变成了日常调度必须面对的基本问题——这不是某个局部电网的特殊需求而是系统性的普遍现象。1.2 把天然气网络拉进来的核心原因新能源出力波动大系统需要快速响应的灵活调节电源。燃气轮机启动快、调节范围宽、碳排放强度相对燃煤更低成了非常常见的配套电源。但燃气机组接入之后电力系统和天然气系统就不再是“各自为政”的独立系统了。燃气机组发电需要消耗天然气天然气管网的气压、管道容量、气源供气能力都会限制燃气机组的可用出力。我在做这个项目之前看过不少电气独立优化的案例最典型的失败场景就是电网调度给出的燃气机组出力计划在天然气侧校验时发现管道流量越限、节点气压跌破下限导致计划根本执行不了。所以要想得到可落地的运行方案必须把电力系统、天然气系统和耦合元件放到一个模型里统一求解。这就是“电气互联系统”这个表述的含义——这里的“电气”指的是电力-天然气互联而不是通常说的电气工程。进一步看还有P2G电转气技术的引入。富余风电可以转化为天然气注入管网相当于给系统增加了一条能量时移通道。P2G让电-气系统从“单向供气”变成“双向耦合”优化问题的可行域和运行方式都发生了本质变化。这种系统中的无功优化自然不能再沿用过去“固定有功出力、只调无功”的老思路。1.3 有功-无功协同优化与传统无功优化的本质区别传统无功优化的标准做法是先通过有功经济调度确定发电机有功出力然后在有功计划固定的前提下以发电机端电压、无功补偿装置容量、变压器分接头位置为控制变量以网损最小或电压偏差最小为目标做二次优化。这种解耦策略在运行工况稳定、系统耦合关系简单的场景下问题不大。但电气互联系统里有功和无功的耦合关系变得非常紧密。燃气机组的耗气量是有功出力的函数这是电-气耦合的核心而无功出力会影响节点电压、进而影响线路有功损耗线路损耗反过来又会影响系统的有功平衡和总购气量。换句话说有功决定无功的“边界条件”无功又通过网损影响有功的“结果”。把两者拆开优化很容易陷入局部最优。协同优化做得更彻底把有功出力、无功出力、天然气气源产气量、管道流量、P2G电功率等统一作为决策变量在同一个优化问题里求解。这样得到的不再是“先定有功再调无功”的次优方案而是整个电-气系统在当前约束下的全局最优运行点。代价也很明显——模型规模变大、非线性约束增多、求解难度成倍上升。怎么在精度和可解性之间做平衡是后面建模和代码阶段的核心问题。2. 核心数学模型电气互联系统怎么用方程描述2.1 电力系统潮流约束与无功控制变量建模电力系统的核心约束是潮流方程。我用的是极坐标形式的交流潮流模型节点i的有功和无功注入方程为$$P_i V_i \sum_{j \in N_i} V_j (G_{ij} \cos\theta_{ij} B_{ij} \sin\theta_{ij})$$$$Q_i V_i \sum_{j \in N_i} V_j (G_{ij} \sin\theta_{ij} - B_{ij} \cos\theta_{ij})$$其中 (P_i) 和 (Q_i) 是节点注入有功和无功(V_i) 是节点电压幅值(\theta_{ij}) 是节点i和j的相角差(G_{ij})、(B_{ij}) 是导纳矩阵的实部和虚部。模型中需要显式表达的无功控制手段主要有三类发电机端电压幅值 (V_{G})、无功补偿装置出力 (Q_C)、变压器分接头档位 (T_k)。它们的取值范围作为约束加入模型$$V_{G,min} \leq V_G \leq V_{G,max}$$$$0 \leq Q_C \leq Q_{C,max}$$$$T_{k,min} \leq T_k \leq T_{k,max}$$变压器分接头通常被建模为整数变量这会引入混合整数约束。如果只是做常规场景测试先把分接头档位离散化步长设置得大一些能显著降低求解压力。我看很多新手拿到项目就喜欢把分接头步长设置成0.5%甚至更小结果求解器跑到天荒地老也出不来。实际工程中1.25%或2.5%的步长已经很够用了。2.2 天然气系统网络方程与管道约束处理天然气系统的核心变量是节点气压 (p_i) 和管道流量 (f_{ij})。节点流量平衡方程可以写成$$\sum_{s \in S_i} F_s - \sum_{l \in L_i} F_l - \sum_{g \in G_i} F_{g} \sum_{p2g \in P2G_i} F_{p2g} L_i^{gas}$$其中 (F_s) 是气源供气量(F_l) 是管道流出流量(F_g) 是燃气机组耗气量(F_{p2g}) 是P2G注入流量(L_i^{gas}) 是天然气负荷。管道流量与两端气压的关系标准用Weymouth方程描述$$f_{ij} K_{ij} \cdot \text{sign}(p_i^2 - p_j^2) \cdot \sqrt{|p_i^2 - p_j^2|}$$这个方程天然是非凸的直接丢给求解器基本无解。实际做法有两种一是引入辅助变量 (u p^2) 把方程转化为关于 (u) 的线性形式后再处理符号函数二是在目标函数中引入惩罚项或对Weymouth方程做分段线性逼近。我在代码里使用的是把气压平方作为变量、然后对符号项做引入大M法线性化处理的方式这样整个模型可以直接交给MIP或SOCP求解器。压缩机的建模也比较关键。压缩机用来弥补管道压降它的运行需要消耗天然气或者电能两种方式在耦合模型里都可以考虑。为了控制规模我在基础模型中把压缩机简化成带压比约束的升压节点$$p_{out} \leq R_{c,max} \cdot p_{in}$$这样既能保证气压水平在可接受范围又不至于引入过多非线性项。2.3 燃气轮机和P2G等耦合元件的运行特征电-气系统的耦合元件是燃气轮机和P2G设备。燃气机组的耗气量通常用二次函数近似$$F_g \alpha_g P_g^2 \beta_g P_g \gamma_g$$为了降低非线性我在实际实现中使用了分段线性化逼近效果很稳定。分段数取3到5段已经能保证很高的精度如果再往上加精度提升的幅度非常有限求解时间却会成倍变长。P2G设备的耦合关系相对简单$$F_{p2g} \eta_{p2g} \cdot P_{p2g} / LHV_{gas}$$其中 (\eta_{p2g}) 是电转天然气效率(LHV_{gas}) 是天然气低位热值(P_{p2g}) 是P2G消耗的电功率。P2G可以制氢也可以进一步甲烷化注入天然气管网。这里按注入天然气管网来处理耦合关系更直接。2.4 约束汇总一张表看清约束体系分类约束名称数学描述处理方式电力约束节点有功平衡(P_i V_i \sum V_j (G \cos\theta B \sin\theta))交流潮流或SOCP松弛电力约束节点无功平衡(Q_i V_i \sum V_j (G \sin\theta - B \cos\theta))同上电力约束机组有功/无功上下限(P_{min} \leq P_g \leq P_{max})线性电力约束电压幅值范围(V_{min} \leq V \leq V_{max})线性电力约束无功补偿容量(Q_{C,min} \leq Q_C \leq Q_{C,max})线性天然气约束节点流量平衡(\sum F_s - \sum F_l - \sum F_g \sum F_{p2g} L)线性天然气约束管道流量特性(f K \cdot \text{sign}(\Delta p^2) \sqrt{|\Delta p^2|})平方变量线性化天然气约束节点气压上下限(p_{min} \leq p \leq p_{max})线性天然气约束气源供气上限(F_s \leq F_{s,max})线性耦合约束燃气机组耗气量(F_g \alpha P_g^2 \beta P_g \gamma)分段线性化耦合约束P2G产气量(F_{p2g} \eta P_{p2g} / LHV)线性碳排放约束总排放量限制(\sum E_{gen} \sum E_{gas} \leq E_{max})线性3. 目标函数设计碳排放成本怎么量化和嵌入3.1 碳排放的来源与量化方式在电气互联系统里碳排放来源主要有三块火电机组直接排放、燃气机组燃烧排放、气源侧开采和输送产生的间接排放。前两者与发电出力直接相关后者与天然气流量相关。火电机组的碳排放量和有功出力的关系近似线性$$E_{coal,i} e_{coal} \cdot P_i \cdot \Delta t$$燃气机组则根据耗气量折算$$E_{gas,i} e_{gas} \cdot F_{g,i} \cdot \Delta t$$气源端排放可以按供气量的固定比例折算。如果系统有外购电力还要考虑外购电力的隐含排放——这个在跨区互联场景里经常被忽略但实际占比不低。3.2 三种目标函数构建思路对比目标函数怎么选直接决定模型的性质和求解复杂度。我试验过三种构建方式各有适用场景方式目标函数形式特点适用场景经济成本碳成本(\min C_{fuel} C_{gas} C_{carbon})单目标方便求解碳成本是线性项模型保持凸性工程应用最常用碳排放最小化(\min E_{total})单目标但经济性完全被忽略实际难以落地低碳政策研究多目标Pareto(\min (C_{total}, E_{total}))需要NSGA-II等算法计算量大结果需要做决策偏好分析学术科研我个人最推荐第一种。碳价格合理设置后碳排放会通过成本传导机制被自动“惩罚”不需要额外引入复杂的多目标处理流程。而且线性碳成本不破坏模型的凸性求解器跑起来非常稳。代码调试阶段我用的是第一种后面做灵敏度分析时才把目标函数改成碳排放最小化做对比。3.3 碳配额、阶梯碳价和低碳激励的处理碳中和目标下比较现实的碳管控手段是碳交易和碳配额。模型里可以这样建模$$E_{total} - E_{quota} \leq \Delta E_{buy}$$目标函数中的碳成本项为$$C_{carbon} c_{co2} \cdot (E_{total} - E_{quota})$$如果 (E_{total} - E_{quota}) 为负说明系统低碳运行有余量可以把配额出售目标函数里体现为收益项。我在项目里还测试过阶梯碳价即超过配额越多、边际碳价越高相当于分档惩罚。这种设置下优化结果会更倾向于压低高碳机组的出力效果立竿见影。还有一个非常容易忽略的细节碳排放约束与有功-无功耦合的关系。无功出力增加会增大线路电流进而增大有功损耗这部分损耗也需要由发电机出力来弥补所以网损增加会带来额外的间接碳排放。这就是无功优化在低碳模型里不能被忽略的根本原因。把碳排放项放进目标函数之后模型会自动倾向于通过无功优化降低网损、减少碳排放——这是协同优化的核心价值。4. Matlab代码实现从参数表到求解器调通的完整路径4.1 求解策略选型SOCP松弛与模型凸化电气互联系统的原始模型是非凸的直接求解非常困难。我采用的策略是对电力系统部分做二阶锥松弛SOCP将潮流约束转换为锥约束对天然气系统部分把气压平方作为变量将Weymouth方程线性化。经过这两步处理整个模型转化成一个混合整数二阶锥规划MISOCP可以用成熟的商业求解器高效求解。这里必须提醒一点松弛不是免费的午餐。SOCP松弛在某些情况下会引入松弛间隙松弛后的最优解可能不满足原始潮流方程。实践经验是收敛之后务必把电压和潮流回代到原始方程里验证一遍算一下对偶间隙或松弛误差。如果误差超过1%通常需要给目标函数增加合适的惩罚项来压缩松弛间隙。我在项目里就用到了这个技巧把松弛变量的和加上一个很小的惩罚系数效果很好。4.2 代码结构设计与核心函数拆分整个Matlab工程我拆成了五个核心文件分工非常清楚|-- main.m % 主程序数据加载、模型构建、求解调用、结果展示 |-- data_case.m % 算例数据电网参数、气网参数、耦合元件参数、碳参数 |-- build_power_model.m % 电力系统约束构建 |-- build_gas_model.m % 天然气系统约束构建 |-- build_coupling_model.m % 耦合元件约束构建这样拆的好处是改数据不用动模型改电力部分不影响天然气部分。项目后期我加了P2G设备只需要在coupling模型里增加几行约束主程序完全不受影响。对于做课题研究、需要反复改动模型的人来说这种模块化结构能省大量时间。4.3 关键代码段逐段拆解下面直接上关键代码我按“变量定义 - 目标函数 - 约束构建 - 求解配置”的顺序来展示。第一步定义决策变量%% 电力侧变量 P_g sdpvar(n_gen, 1); % 发电机有功出力 Q_g sdpvar(n_gen, 1); % 发电机无功出力 V sdpvar(n_bus, 1); % 节点电压幅值 theta sdpvar(n_bus, 1); % 节点相角 Q_c sdpvar(n_cap, 1); % 无功补偿出力 %% 天然气侧变量 F_s sdpvar(n_source, 1); % 气源供气量 F_pipe sdpvar(n_pipe, 1); % 管道流量 p_sq sdpvar(n_gas_bus, 1); % 节点气压平方 F_g sdpvar(n_gu, 1); % 燃气机组耗气量 %% 耦合变量 P_p2g sdpvar(n_p2g, 1); % P2G消耗电功率 F_p2g_out sdpvar(n_p2g, 1); % P2G输出天然气流量第二步目标函数构建%% 目标函数发电成本 购气成本 碳排放成本 网损惩罚 obj_fuel coal_cost_coeff(:,1) .* P_g(filter_coal).^2 ... coal_cost_coeff(:,2) .* P_g(filter_coal) ... coal_cost_coeff(:,3); obj_gas gas_price * F_s; obj_carbon c_co2 * (sum(coal_emission_coeff .* P_g(filter_coal)) ... sum(gas_emission_coeff .* F_g) - E_quota) * delta_t; obj sum(obj_fuel) sum(obj_gas) obj_carbon;第三步电力系统约束构建Constraints []; %% 有功、无功平衡约束我这里按直流潮流简化和交流SOCP松弛两种方式都跑过 % 交流SOCP松弛方式 for i 1:n_bus % 节点有功平衡 Constraints [Constraints, P_g_i(i) P_renew_i(i) - P_load_i(i) ... V(i)^2*G(i,i) sum(V(i).*V(j).*(G(i,j)*cos(theta(i)-theta(j)) B(i,j)*sin(theta(i)-theta(j))))]; % 节点无功平衡 Constraints [Constraints, Q_g_i(i) Q_c_i(i) - Q_load_i(i) ... -V(i)^2*B(i,i) sum(V(i).*V(j).*(G(i,j)*sin(theta(i)-theta(j)) - B(i,j)*cos(theta(i)-theta(j))))]; end %% 机组出力约束 Constraints [Constraints, P_g_min P_g P_g_max]; Constraints [Constraints, Q_g_min Q_g Q_g_max]; %% 节点电压约束 Constraints [Constraints, V_min V V_max]; %% 无功补偿约束 Constraints [Constraints, 0 Q_c Q_c_max];第四步天然气系统约束构建%% 节点流量平衡 for i 1:n_gas_bus inflow sum(F_s(find_source(i))) sum(F_pipe_in(find_in_pipe(i))); outflow sum(F_pipe_out(find_out_pipe(i))) sum(F_g(find_gas_load(i))) gas_demand(i); Constraints [Constraints, inflow - outflow 0]; end %% 管道Weymouth方程线性化 % 使用气压平方作为变量向量化构建 for k 1:n_pipe from pipe_from(k); to pipe_to(k); Constraints [Constraints, F_pipe(k) pipe_C(k) * ... (p_sq(from) - p_sq(to)) / max(abs(p_sq(from) - p_sq(to)), 1e-6)]; end %% 节点气压上下限约束 Constraints [Constraints, p_min.^2 p_sq p_max.^2]; %% 气源供气约束 Constraints [Constraints, 0 F_s F_s_max];4.4 求解器配置与收敛性设置Matlab下用Yalmip建模求解器我切换到Gurobi或Cplex。Gurobi在求解MISOCP时表现稳定Cplex同样可以。关键是求解器选项要设置到位否则默认参数经常出问题。我常用的配置如下options sdpsettings(solver, gurobi, ... verbose, 2, ... gurobi.MIPGap, 0.01, ... gurobi.TimeLimit, 300, ... gurobi.NumericFocus, 1);MIPGap控制在1%以内既能让结果可信又能控制求解时间。TimeLimit设成300秒避免某些场景下求解器无限迭代。NumericFocus设置为1可以在模型尺度差异较大时提升数值稳定性——这点在做电气互联系统时特别重要因为电力标幺值在0.9到1.1之间而天然气气压平方可能在几百到上千的绝对量级数值范围差距非常大。5. 仿真结果分析优化效果与参数灵敏度5.1 算例设置与基础数据准备我用的是修改版IEEE 30节点系统耦合一个6节点天然气网络。电力系统部分挂了5台火电机组、2台燃气机组、2台风电场、1个无功补偿站天然气系统部分有2个气源、6条管道、1个P2G设备。碳配额设置为系统历史排放水平的90%碳价按照50元/吨起步。这个配置参考了当前多数综合能源系统文献的标准做法既能反映系统耦合关系又不至于因为网络规模太大导致调试困难。做方向验证的阶段节点数不是越多越好——我记得自己刚上手时就吃过亏直接拿118节点系统测试结果每次改动模型光等求解器就要十几分钟非常影响排查问题的效率。5.2 优化前后的系统运行对比指标独立优化电网单独、气网单独协同优化变化幅度系统总运行成本/万元234.5218.7-6.7%碳排放总量/tCO21124983-12.5%网损/MWh28.622.3-22.0%电压最低点/p.u.0.9410.9682.9%燃气机组出力/MW12815621.9%协同优化的效果非常直观总成本下降碳排放大幅降低网损也明显改善。原因不复杂——独立优化时电网侧只关心发电成本容易把出力压力加在燃煤机组上而协同优化可以统筹天然气价格和碳排放惩罚让燃气机组在更经济、更低碳的区间运行。同时无功优化的加入降低了网损系统有功需求减少各机组的出力都有下调空间。5.3 碳价和新能源渗透率对优化结果的影响碳价从20元/吨一路涨到200元/吨系统运行方式的变化趋势很清楚碳价越高燃气机组出力占比越大煤电出力被逐步替代碳排放总量持续下降。但当碳价超过150元/吨后碳排放下降速度明显趋缓原因是剩余排放基本来自系统必须开机的边际机组可调节空间已经不大。这说明碳价信号的边际效益有一个递减区间做政策仿真时可以重点关注这个拐点。新能源渗透率从15%提高到40%系统的特点是无功优化压力明显增大因为风电场的无功支撑能力有限电压越限风险增加需要投入更多无功补偿资源。同时在新能源出力较高的时段P2G设备的运行功率明显上升富余电力被转化为天然气存储起来极大改善了系统的消纳水平。协同优化模型能很好地捕捉到这些跨系统的转移关系这是只做电力系统优化看不到的结论。6. 调试实录这些坑我替你们先踩了6.1 常见报错与解决办法速查表报错现象根本原因解决方法Yalmip报“No suitable solver”未安装Cplex/Gurobi或未配置路径正确安装求解器并运行yalmip clear后重新检测Warning: Solver not applicable模型包含非凸二次约束转用SOCP松弛或分段线性化处理求解结果出现NaN初始值设置不合理或数值尺度失衡检查单位和标幺值统一量纲Gurobi提示Infeasible model约束之间矛盾用export模型到lp文件人工检查冲突约束求解速度极慢整数变量过多或Big-M值过大合理设置分接头步长缩小Big-M电压结果越界但不报错潮流松弛间隙过大对比SOCP松弛解与原始潮流方程的误差增加惩罚项6.2 模型不收敛的排查思路模型不收敛是这类多能流优化项目最常见的问题。我的排查顺序是固定的先看可行域是否为空再看凸化是否有效最后看数值条件是否良好。可行域为空的概率最高通常是某个约束参数填错了。我之前调试时遇到一次死活找不到可行解查了两天才发现是天然气系统的气压下限写成了绝对压力5MPa而系统正常运行时该节点压力只有0.8MPa——参数量纲搞混了。这提醒大家在建模时统一单位的习惯必须养好电力部分用标幺值天然气部分用MPa和kcf/h耦合环节再统一用热值折算。凸化失效的典型表现是模型报告求解成功但把结果代回去原始潮流方程误差超过5%。解决方法是检查SOCP松弛是否满足精确性条件必要时在目标函数里加松弛惩罚项。数值条件问题则相对容易处理把决策变量做归一化或者调整求解器的NumericFocus参数通常就能解决。6.3 让代码更健壮的几个小习惯第一自定义数据结构而不是散装数组。我在data_case.m里用的是struct数组比如bus_data.voltage_min、gas_data.pressure_min这样虽然代码行数略多但出错的概率大幅下降以后改参数也不用在代码里到处搜索替换。第二构建矩阵用矢量化而非逐点循环。上面的代码里为了可读性用了for循环但在实际工程版本中所有约束都是用矩阵拼接和矢量化运算完成的。Matlab的for循环处理大规模稀疏矩阵时性能差异非常大尤其IEEE 118节点这类中等规模系统矢量化能提速几十倍。第三求解完成之后一定要做可行性验证。不能只看solver返回的求解状态还要把最优解代回原始约束逐条检查。我习惯写一个check_results.m函数输出每条约束的最大违背量如果超过1e-4就认为是数值问题需要处理。这个习惯帮我在后期省了无数时间强烈建议保留。一点心得这个项目做下来最大的体会是多能流协同优化和单能流优化有本质区别。单能流系统里能量只走一条路约束结构相对规整一旦把电和气放在一个模型里面对的不仅是设备模型数量的叠加更是不同物理量纲、不同网络特性、不同时间常数之间的碰撞。代码写的再漂亮如果物理关系没理顺求解结果依然没有工程意义。最后分享一个实用技巧做这类项目时不要一开始就上完整模型先从最简单的“电力系统单台燃气机组”开始跑通了再加天然气网络再上P2G设备最后再引入碳排放约束和碳交易机制。每一步都做结果验证确认无误后再叠加下一步。这样即使中间出了问题也能迅速定位到是哪一层引入的矛盾。我的经验是超过一半的“模型不收敛”问题其实是前一步的模型还没完全调稳就急着加复杂度导致的。