做综合能源的硕士生，看到“计及P2G厂站的电-气综合能源系统规划”这种题目，第一反应应该是：这又是个要啃硬骨头的复现。P2G（Power to Gas）的本质是电转气，把电网里富余的电能变成天然气，让电力系统和天然气系统真正意义上双向互动。这篇复现的核心，就是用Matlab搭建一套“电网+气网+P2G选址定容”的联合规划模型，输出可跑的代码和清晰的求解逻辑。我当时拿到这个题目，第一反应是找作者要源码，结果跟大多数情况一样——没要到。自己吭哧吭哧复现了一周多，踩了不少坑，才把模型跑透。这篇记录就把整个过程拆开讲讲，重点是代码怎么组织、非线性约束怎么处理、以及那些容易让模型“跑飞”的细节。适合正在做电-气耦合规划课题的硕士生、刚入门综合能源方向的研究生，也适合想了解P2G建模要点的工程师。
1. 先理解P2G为什么是电-气系统规划的关键
1.1 电-气耦合的现实背景：为什么规划要绑定在一起
传统规划思路里，电网和天然气网是各管各的：电网规划只关心电源和线路，气网规划只关心气源和管道。但现实早就不是这样了，燃气轮机大规模接入后，电网的出力上限很多时候不是由机组决定的，而是由气网能供多少气决定的。反过来，气网里加压站、调压站的驱动能源也可能来自电网。两个系统在运行层面已经深度耦合，如果规划阶段还分开做，很容易出现“电网规划落地后气网供不上气”或者“气网扩建后电网送不出电”这种尴尬局面。
P2G在中间扮演的角色很特殊。它是目前唯一能让能量从电力系统反向流入天然气系统的设备，相当于在两个系统之间装了一个双向阀门。电力系统富余时，把多余的电变成甲烷注入气网储存；气网缺气时，这个“人造气源”顶上去。对规划来说，P2G的存在意味着电网和气网可以互为备用，比如某个区域电网输电阻塞，与其扩建线路，不如就地建一个P2G厂站把电转成气送走，这往往更经济。所以规划阶段就必须把P2G作为一种可选方案放进优化框架里，和线路扩容、气源扩建一起比选，这才是“计及P2G厂站”的真正含义。
另一个现实驱动力是可再生能源消纳。风电光伏大发时段电网消化不了，弃风弃光成了常态。P2G恰好提供了一个消纳出口：把本来要丢弃的电能转成气存起来，既缓解了弃电压力，又给气网增加了气源。我在自己的算例里配了30%比例的风电，基准场景弃风率大概12%，引入P2G之后可以降到5%以下。这就是这个课题最核心的应用场景：高比例新能源下的多能互补。
1.2 P2G厂站的内部结构与效率边界
很多刚接触P2G的人以为它是一台“电变气”的单一设备，其实不是。P2G厂站内部是两段工艺的串联：第一段是电解水制氢，碱性电解槽或者PEM电解槽，靠电把水分子拆成氢气和氧气；第二段是甲烷化反应，氢气在催化剂作用下与二氧化碳反应生成甲烷和水。之所以非要费劲做成甲烷而不是直接用氢气，是因为现有天然气管道、储气库和燃气设备对氢含量有严格限制，氢气比例过高会腐蚀管材、影响调压器寿命，而甲烷可以直接以一定比例掺混甚至完全替代天然气。
两段工艺耦合下来，整体能量转换效率大约在0.55到0.65之间，我代码里取0.6。这个数字很关键，它几乎决定了P2G在经济上什么时候划算。做个简单换算：1 MWh电力按0.6效率转换，得到0.6 MWh热值的天然气。标准天然气热值约36 MJ/m3，0.6 MWh对应2160 MJ，折算下来就是大约60标准立方米的天然气。换句话说，P2G厂站每消耗1 MWh电，大概能产出60 Nm3的甲烷。如果电价正常水平0.5元/kWh，光电费成本就相当于每立方米天然气8元多，远高于常规气价，所以P2G只有在弃电价格极低或者有额外政策补偿时才具备经济性。这也是为什么几乎所有的P2G规划文献都把场景设在“高比例新能源+高弃风”背景下。
这个效率边界还会直接影响规划结果：P2G容量配置大了，弃风消纳多，但投资和运维成本也高；容量配置小了，又起不到明显作用。所以模型里效率、投资成本、弃电价格这些参数都要做敏感性分析，不能只跑一组数据就下结论。
2. 规划建模的完整拆解：从物理网络到数学表达
2.1 电网侧建模：为什么规划场景下一定要用直流潮流
电网模型在规划问题里没有太多悬念，几乎清一色用直流潮流。原因很简单：规划问题本身就带有0-1整数变量，比如某个候选节点建不建P2G，这已经让模型变成混合整数规划了。如果电网侧再用交流潮流，支路潮流是电压相角的非线性函数，整个问题就成了混合整数非线性规划（MINLP），求解难度和耗时直接上升几个数量级，而且很难保证收敛到全局最优。
直流潮流的数学形式是线性的：支路有功潮流近似等于两端相角差除以支路电抗，节点功率平衡也变成线性等式约束。整个电网侧约束组就变成了一组线性方程，可以顺利嵌入MILP框架。
直流潮流牺牲掉的精度主要在网络损耗和电压支撑方面，对规划阶段来说误差通常在5%到10%区间，完全够用。但要知道它的边界：算出来的线路潮流和相角是近似值，电压不是变量而是被忽略的假设。如果算例里含大量风电接入，电压支撑问题不能靠直流潮流评估。我自己的做法是：MILP先求出规划方案，再用MATPOWER跑一次交流潮流做事后校验，确认节点电压不越限。这个“规划粗算+运行精算”的双校验流程，写论文时也很加分。
2.2 气网侧建模：Weymouth方程是第一个真正的坑
气网的建模比电网粗暴得多，但也复杂得多。核心原因在于它的物理约束是强非线性的。
节点流量平衡是线性的：每个气节点的进气量加气源注入量等于负荷量。但管道潮流与节点压力的关系遵循Weymouth方程，这是输气管道在稳态下的经典经验公式：
F_ij = K_ij · sign(π_i - π_j) · sqrt(|π_i - π_j|)
这里π_i是节点i气压的平方值，K_ij是管道常数，由管径、管长、摩擦系数决定。式子里那个sign符号和根号组合，导致函数既非线性又不平滑。π_i本身又是变量，所以约束里有变量相乘、开方、符号判断，直接扔给求解器基本无解。
初学者最常犯的错误是忽略符号处理。当节点i压力低于节点j时，气压平方差为负，流量方向要反向，如果直接把根号里的数当成非负处理，模型逻辑就错了。稳妥的做法是引入辅助变量Δπ=|π_i-π_j|，再配一个0-1变量表示流向，但这样会增加整数变量数量。
处理Weymouth方程的标准做法是分段线性化。因为π的取值范围是已知的——气压有上下限，所以Δπ的运行区间也是已知的。把sqrt函数在这个区间内切成若干线段，用增量模型或SOS2约束把它近似成线性关系。我在代码里每根管道切6段，精度误差大约1%，对规划结果的影响很小。分段数太少精度差，太多则整数变量爆炸，6到8段是实践中不错的平衡点。
2.3 P2G厂站建模：核心就四个约束
P2G厂站在规划模型里的数学表达不复杂，但细节容易出错。核心约束可以归成四类。
第一，选址与建设约束。候选节点i要么建设P2G厂站，要么不建，用0-1变量x_i表示。建设容量P_cap不能超过上限，并且只有x_i=1时才能有正容量，写成约束就是0 ≤ P_cap ≤ x_i · P_cap_max。这里P_cap_max可以取一个合理上限，比如该节点最大可接入容量。
第二，运行约束。P2G厂站在每个时段消耗的电功率不能超过建设容量。注意这里有个很多文献会忽略的细节：P2G设备存在最小运行负荷比，比如电解槽一般不能低于20%负荷运行，否则效率急剧恶化甚至损坏设备。于是约束变成P_min · x_i ≤ P_p2g(t) ≤ P_cap，其中P_min是额定容量的某个比例。
第三，能量转换约束。P2G的产气量与耗电功率成正比，G_p2g(t) = η · P_p2g(t) · GAS_CONVERT。这里GAS_CONVERT是单位换算系数，把MWh电量转换成标准立方米天然气。按照前面的换算，η=0.6时，1 MWh电对应60 Nm3气，所以转换系数的数值要仔细核对。
第四，耦合位置约束。P2G厂站在电网侧是电负荷，在气网侧是气源。建模时要明确它接入哪个电网节点、接入哪个气网节点，两个接口的功率必须由同一个决策变量耦合起来。
这四类约束都是线性的。只要把气网Weymouth方程线性化处理好，整个模型就可以保持MILP结构。
2.4 目标函数与约束体系汇总
目标函数通常是总成本最小化，包含三块：P2G厂站投资成本、系统运行成本、其他扩展成本。
投资成本按年化处理。P2G厂站一次性投资C_inv_total = C_fix · x_i + C_unit · P_cap，C_fix是固定建设成本（场地、并联设备等），C_unit是单位容量成本。年化要把寿命期内总成本折成等年值，公式涉及贴现率。不同文献贴现率取6%到10%，我代码里取8%，整个目标函数里贴现率是影响投资决策的一个敏感参数，改起来要谨慎。
运行成本包括：向上级电网购电成本、向气网购气成本、P2G厂站运维成本，可能还有弃风惩罚成本。这里有个时间尺度统一的问题：规划模型用典型日代表全年，运行成本要乘以典型日代表的天数再年化。比如冬季典型日和夏季典型日分别代表90天和275天，那么对应运行成本要乘以对应天数，不能直接加总。
约束体系汇总起来就是：电网节点功率平衡、线路潮流上下限、发电机出力上下限、气网节点流量平衡、管道流量约束、节点气压上下限、气源出力范围、P2G厂站建设与运行约束，加上所有变量的取值范围。整个约束体系写出来大约有上千行，但借助Matlab的YALMIP工具箱，构建起来比想象中省力。
3. Matlab代码实现：从模型到可运行代码
3.1 非线性环节的处理思路对比
在写代码之前，先得把求解策略想清楚。气网Weymouth非线性的处理，主流有三种思路，各有适用场景。
第一种是MILP分段线性化，这也是我代码最终采用的方式。把Weymouth方程分段逼近之后，整个模型变成混合整数线性规划，能保证全局最优，求解器可以直接给出结果。缺点是需要引入额外整数变量，分段越细整数变量越多，求解时间相应上升。但对规划问题来说，这个代价是完全可以接受的。
第二种是顺序迭代法，也叫解耦迭代。先假设P2G规划方案固定，电网问题与气网问题交替求解，通过耦合变量迭代修正。这个方法的好处是每一侧都可以用很成熟的商业软件，比如电网侧用MATPOWER、气网侧用专门的气网潮流工具，但缺点是整个过程无法保证收敛到全局最优，甚至可能出现振荡不收敛。
第三种是锥松弛。把气压的平方作为新变量，对非凸项做二阶锥松弛，原问题变成二阶锥规划（SOCP）。精度比线性化高，但YALMIP和Gurobi对SOCP的支持虽然很好，遇到大规模算例时数值稳定性有时候不如纯MILP。而且锥松弛在边界条件比较紧的情况下可能松弛间隙较大，结果需要通过后验收紧。
三类方法对比如下：
| 方法 | 优点 | 缺点 | 适用场景 |
|---|---|---|---|
| MILP分段线性化 | 全局最优、实现简单 | 分段误差、整数变量多 | 规划问题主流方案 |
| 顺序迭代 | 可复用成熟软件 | 无法保证全局最优 | 运行模拟、在线调度 |
| 锥松弛 | 精度较高 | 数据要求高、稳定性一般 | 学术研究、特殊网络 |
刚复现这个题目时，不用纠结直接用哪个，先走MILP线性化把完整链路跑通，后面再对比其他方法。
3.2 代码目录结构与运行环境准备
我的代码组织很清晰，五个部分：
- main.m —— 主入口，定义算例并调用建模与求解流程；
- data/ —— 电网节点数据、气网节点数据、负荷曲线、风电曲线、P2G候选集；
- model_build.m —— 用YALMIP定义全部决策变量和约束；
- solve_and_report.m —— 求解并输出报告与图表；
- utils/ —— 通用子函数，比如数据读取、单位换算、Weymouth系数计算。
环境配置方面，Matlab版本建议R2020a以上，YALMIP必须装好。YALMIP是Matlab环境下的建模工具箱，能把优化问题用接近自然语言的方式写出来，翻译成求解器能识别的标准形式。求解器推荐Gurobi，MILP性能目前没有任何替代品能比；学校没有Gurobi许可证的话，也可以用Cplex学术版，或者先用Matlab自带的intlinprog跑小算例验证逻辑。
安装YALMIP之后，强烈建议先跑一遍yalmiptest，确认建模环境与求解器能正常通信。很多“模型没问题但求解报错”的怪事，最终查出来都是YALMIP路径配置或者求解器识别的问题。
3.3 核心建模代码实现细节
下面把我代码里P2G相关的核心建模部分贴出来，这段是整个模型最有代表性的地方：
% 决策变量定义 x_p2g = binvar(n_cand, 1); % 候选厂站是否建设（0-1） P2G_cap = sdpvar(n_cand, 1); % 建设容量，单位MW P2G_pin = sdpvar(n_cand, T); % 各时段输入电功率 P2G_gout = sdpvar(n_cand, T); % 各时段输出天然气流量（Nm3/h） % 选址-容量耦合约束 M_cap = 120; % 容量上限取120MW C = [C, 0 <= P2G_cap <= x_p2g * M_cap]; C = [C, 0 <= P2G_pin <= repmat(P2G_cap, 1, T)]; % 最小运行负荷约束（20%额定容量） C = [C, repmat(0.20 * P2G_cap, 1, T) <= P2G_pin]; % 能量转换约束（效率0.6，换算系数60 Nm3/MWh） GAS_CONVERT = 60; % 0.6效率下1MWh电对应60Nm3气 C = [C, P2G_gout == P2G_eff * P2G_pin * GAS_CONVERT];注意repmat(P2G_cap, 1, T)这种写法，目的是把容量向量扩展到每个时段。很多人会在这里遇到维度报错，其实是因为YALMIP对矩阵隐式扩展支持不如一些其他语言，老老实实用repmat最稳妥。
电网侧直流潮流约束：
% 节点功率平衡 C = [C, sum(Pgen, 2) - P_load - P_p2g_elec == B_theta * theta]; % 线路潮流约束 flow = B_line * theta; % 直流潮流线性表达式 C = [C, -Flow_max <= flow <= Flow_max];这里B_theta是节点导纳矩阵的虚部部分，B_line是支路-节点关联矩阵与电抗对角阵的组合。构建之前先单独验证这两个系数矩阵，可以用简单两节点系统对比手算结果。
气网管道流量分段线性化的核心代码：
% 管道流量-压力关系的分段线性化（增量模型） pieces = 6; bp = linspace(0, DeltaPi_max, pieces + 1); flow_at_bp = K_pipe * sqrt(bp); % 各分段端点对应流量 delta = sdpvar(pieces, 1, 'full'); % 分段增量 zbin = binvar(pieces - 1, 1); % 分段激活指示变量 DeltaPi = 0; Flow = 0; for k = 1:pieces DeltaPi = DeltaPi + delta(k) * (bp(k+1) - bp(k)); Flow = Flow + delta(k) * (flow_at_bp(k+1) - flow_at_bp(k)); end这段是示意结构，真正的完整实现还要加分段激活约束，防止增量变量“跳过”中间段直接填满后面的分段。缺失这条约束时，线性化的分段可能被错位填充，流量数值会出现物理上不合理的跳跃。我用implies配合zbin来保证分段顺序填充，这是整个MILP实现里最容易被忽略的地方。
3.4 求解器参数与收敛控制
模型构建完成之后，求解设置有几个关键参数直接影响成败。Gurobi的MIPGap默认是1e-4，若规划算例规模大，整数变量数量超过几千个，可以把MIPGap放宽到1e-2，求解时间能缩短一个数量级，而规划结果的实用性差异微乎其微。除了MIPGap，还要设置一个合理的时间上限，比如7200秒，防止求解器无限“钻”下去。
代码里求解部分的骨架：
ops = sdpsettings('solver', 'gurobi', 'gurobi.MIPGap', 0.01, ... 'gurobi.TimeLimit', 7200, 'verbose', 2); sol = optimize(C, Objective, ops); if sol.problem == 0 % 求解成功 recon_error = check(C); else % 分析不可行原因 disp(sol.info); endcheck(C)是对所有约束进行残差检查的关键工具。求解结束后如果有约束残差明显偏离0，说明数值稳定性有问题，需要回到M值或者参数精度上去排查。这套“求解+校验+输出”三连流程，能让绝大部分隐蔽问题现形。
4. 算例复现与结果分析：怎么判断代码对不对
4.1 测试系统设计与参数设置
复现论文时，算例选择很讲究。太大跑不动，太小说明不了问题。我默认算例用的是IEEE 14节点电网系统改进版，加了一个风电场；气网用6节点输气系统，两个系统通过两个候选P2G厂站位置互联。风电场装机占系统总负荷约30%，基准场景弃风率约12%。P2G效率0.6，单位投资成本1.2万元/kW，生命周期按20年，贴现率8%。
这个中等规模算例的优点在于：MILP求解大概30秒内一定能跑完，特别适合调试和复现。等逻辑验证通过后，我再把算例替换成大型系统，比如IEEE 118节点电网配合20节点气网，代码框架完全不用动，只需要替换数据文件。
关于候选节点选择，不是所有节点都能建P2G。我的经验是选三类位置：一是风电送出受阻的节点，二是电网负荷中心附近且靠近气网的节点，三是气网气压偏低的节点。这三类位置P2G的价值体现最明显，也最容易跑出有说服力的结果。
4.2 从三个层次解读规划结果
算例跑完，怎么判断结果是否合理？我习惯分三个层次看。
第一层看规划方案本身：P2G落在哪些候选节点，容量配置多大。合理的结果应该落在我前面说的三类节点附近。如果P2G出现在莫名其妙的位置，先检查候选集设置是否合理，再检查目标函数是否写错——比如投资成本忘记年化，会导致P2G容量虚高。
第二层看系统运行指标：弃风率、购电成本、购气成本、总成本的变化。引入P2G后，弃风率应当显著下降，但购气成本可能上升，因为P2G产出的气也要计入系统气源。总成本是否下降取决于弃电价格与气价的相对高低。
第三层看气网压力分布：P2G注入点附近节点压力应当升高，但绝对不能越上限。如果压力越限，说明P2G容量配置超过了气网承受能力，需要调整容量或增加压缩机。这一层校验能直接暴露“电侧OK但气侧崩了”的问题，很多人漏掉这一步导致结果被审稿人质疑。
为了验证P2G的规划效果，我通常跑三个场景做对比：不建设P2G的基准场景、允许P2G规划的场景、P2G强制配置高容量的场景。三个场景总成本和弃风率放一起，规律一目了然。
4.3 敏感性分析与合理性校验
复现论文最怕结果和文献对不上。我的建议是先跑一个无P2G的基准场景，把电网购电成本、气网购气成本分别与单独规划结果核对，这一步能确认网络参数没偏。然后再开P2G，逐步增加候选容量上限，观察目标函数变化是否合理——总成本应该随容量上限先下降后平稳或略微上升，如果出现先上升后下降的“凹”形，多半是约束写错了。
敏感性分析重点扫两个参数：气价和P2G效率。气价从2元/m3扫到5元/m3时，P2G经济性逐渐变差，规划容量应该单调下降；P2G效率从0.5提到0.7时，规划容量应该增加。这种单调性检验虽然基础，但能过滤掉绝大部分建模错误。我甚至建议把敏感性分析结果画成曲线放论文里，比单一算例结果更有说服力。
还有一个容易被质疑的细节：P2G产气与气网负荷的时间匹配。P2G在风电大发时段大量产气，此时气网负荷可能处于低谷，多出来的气如果无法存储，实际消纳效果就要打折扣。所以算例里最好加入储气罐或者允许P2G产气替代其他气源的比例约束，这样结果更贴近工程实际。
5. 复现一路踩过的坑：问题排查速查
5.1 模型不可行的排查路径
复现过程中最高频的问题就是“模型无解”。第一次遇到时，直接懵了，怀疑是不是求解器坏了。后来发现90%的情况是约束维度写错或者耦合变量符号反了，另外10%才是数值问题。
我的排查流程很固定：先用sol.info看求解器返回的状态，如果是“Infeasible”，就用YALMIP的check(C)定位残差最大的约束，然后把Constraints列表逐个注释，用二分法定位哪一组约束导致不可行。一般十分钟内能找到问题。最容易出问题的就三处：电网节点功率平衡里的P2G耗电项符号、气网节点平衡里的P2G产气项符号、容量约束里矩阵维度不匹配。
5.2 大M值的选择：教训很深刻
P2G容量约束里需要一个M值，很多初学者的习惯是越大越好，干脆写1e6。我在某次调试时就吃过这个亏——模型求解特别慢，偶尔报数值警告，后来把M值从1e6缩小到容量上限的1.2倍，模型几秒钟就解完了，数值稳定性也恢复正常。
原因是求解器处理线性规划时，M值过大会导致约束的数值尺度差异太大，预求解器的缩放机制容易失效，进而引发舍入误差。规范做法是：对所有含M的约束，M取该变量实际业务边界的1.2到2倍。比如容量上限120MW，M取150以内足够。这条经验同样适用于其他0-1变量相关的“big-M”约束。
5.3 线性化精度校验这步不能省
分段线性化完成之后，如果你直接拿结果写论文，风险很大。我第一版代码就没做校验，后来手工核对发现：某条管道分段线性化后的流量与真实Weymouth方程计算值偏差达到6%。虽然选址结论没变，但敏感性分析曲线的形状完全错了。
从那以后，我固定会在代码里加一个校验子函数：把规划得到的P2G产气和气负荷回代到真实Weymouth方程里，对比线性化流量与真实流量的相对偏差，要求控制在2%以内。如果超过，就增加分段数或调整分段区间。这步校验的成本很低，但能让结果的可信度上一个台阶。写论文时还可以把这部分作为“模型精度验证”小节，直接加分。
5.4 环境配置与软件版本相关坑
最后说几个环境层面的坑。第一，YALMIP和Gurobi的版本匹配问题。新版Gurobi发布后，旧版YALMIP有时无法正确识别，求解器调用会出幺蛾子。建议先跑一个小的随机线性规划，确认求解器链路通畅，再跑大模型。第二，MATLAB版本太旧也可能导致YALMIP某些函数行为异常。我遇到过一次莫名其妙的check报错，升级MATLAB版本后自动消失。
第三，如果学校没有Gurobi，用Matlab自带的intlinprog也不是不行。但要注意它的性能极限：本文中等算例大约600多个0-1变量，用intlinprog能跑但明显慢，如果算例规模更大，建议还是装个商业求解器。学术许可证对学生免费，申请流程也不复杂，别在这上面耽误时间。
关于常见问题的汇总，整理成速查表：
| 现象 | 可能原因 | 处理办法 |
|---|---|---|
| 模型无解 | 约束维度错、耦合符号反 | check(C)定位残差，二分排查约束 |
| 求解特别慢 | 大M值过大、分段太密 | M取边界1.2倍，分段控制在6-8段 |
| 结果物理上不合理 | 候选集设置问题、P2G符号反 | 校验候选节点与平衡方程符号 |
| 线性化误差偏大 | 分段数不足或区间设置不当 | 增加分段数并回代校验Weymouth |
| 求解器报数值警告 | 大M值过大、参数尺度差异大 | 缩小M值、检查变量量纲 |
写在最后，复现论文这件事，10%的精力花在研究模型公式上，90%都在“伺候”数据和代码。模型跑不通别慌，把问题拆小、逐项验证，比到处找别人的代码模板有用得多。尤其是气网的Weymouth线性化、P2G的耦合符号、M值选择这三个点，这个模型里数值层面的坑基本都集中在这几个位置。把这几关过了，整个代码框架就变成你自己的了。后面你想加储能、加氢能、加碳交易，都能在这个底座上快速扩展，这也是复现类工作最大的价值——不是拿到一份能跑的代码，而是真正理解这个系统是怎么被数学表达出来的。
