计及P2G厂站的电-气综合能源系统规划与Matlab实现

提到“计及P2G厂站的电-气综合能源系统规划研究”这个题很多刚开始做论文复现的同学第一反应就是懵P2G到底是个什么装置规划研究和运行优化又差在哪里Matlab代码从哪里下手我当时拿到这个复现任务时也是从头查概念等真正把论文里的模型拆开才发现这条链路并没有想象中复杂——核心就是一句话在电力系统和天然气系统之间加入P2G厂站作为耦合环节然后在给定网络拓扑、负荷曲线、资源特性的前提下回答“P2G厂站建在哪些节点、建多大容量、系统总成本最低”的问题。这篇文章适合正在复现综合能源方向论文的在校研究生、刚接触电气耦合规划的工程师以及想快速把学术模型转成可运行程序的同行。我会把从模型公式到Matlab实现的完整流程以及我实际踩过的坑全部写出来。1. 复现之前先搞懂这个课题到底在算什么1.1 P2G厂站为什么值得写进规划模型P2GPower-to-Gas电转气的核心思路是把电能转化为燃气形态的能量存储。主流技术路线是先用电力电解水制氢氢气可以直接混入天然气管道也可以再经过甲烷化反应生成甲烷后注入燃气系统。听上去像是一条简单的“电→氢→气”串联路线但放到综合能源系统里看它带来的最大变化是电网和气网从单向耦合变成了双向耦合。过去电气综合能源系统里电网与气网主要通过燃气轮机发生关联燃气轮机消耗天然气发电能量流向是“气往电走”。P2G一进来电网富余的电也可以反向变成天然气系统内部出现了一条“电往气走”的新通道。规划模型因此必须把两个网络放到同一个优化问题里联立求解而不能像早期工程那样电网规划归电网、管网规划归管网。这就是题目里“计及”两个字的分量。从实际收益来看P2G在综合能源规划里主要有三个作用支撑可再生能源消纳、缓解电网局部阻塞、提供长时间尺度的能量转移。风电和光伏的出力高峰经常是负荷低谷如果不额外做储能弃风弃光很难避免。电池储能适合解决小时级波动但应对数天以上的能量转移成本很高。P2G把多余电能转成气既可以用天然气管道大量存储又能在燃气轮机需要时再发电时间尺度灵活得多。1.2 一份完整的复现项目到底要输出什么标题里写的是“规划研究”我用大白话翻译一下在满足电力负荷、天然气负荷、新能源出力约束的前提下合理选择P2G厂站的位置和容量同时给出系统逐时段的运行方案让全生命周期成本最低。对应到论文里的图表复现结果一般需要包含这样几类内容P2G厂站的选址与定容哪些候选节点被选中建设每个节点的容量是多少逐时段运行结果P2G输入功率、注入气网的燃气量、燃气轮机出力、购电量和购气量等经济性指标总投资年值、年运行成本、总成本、各项成本占比对比分析有P2G和无P2G两种场景下总成本与新能源消纳率的差异灵敏度分析P2G单位投资成本、气价、风电渗透率变化时规划结果如何变化。这些内容本质上就是论文里的核心图表复现项目就是把这些图表背后的数学模型重新搭起来用优化求解器算出数值结果并努力和论文给出来的结果对齐。1.3 我建议的复现顺序别一上来就写完整模型这个课题我一开始也走过弯路想着直接把完整规划模型写进一个脚本里跑结果约束写了几百行报警告之后根本不知道问题出在哪个模块。后面我把流程调整成五步走效率高了很多。第一步是精读建模公式把论文里的决策变量、约束条件、目标函数逐一列成表格脑子里先形成“模型清单”。第二步是整理算例数据网络拓扑参数、负荷曲线、风速曲线、分时电价、气价全部统一单位。第三步是先写一个不含P2G的基础运行优化模型验证电网和气网约束本身没有写错。第四步再加P2G的容量决策变量与建设变量把纯运行问题升级成混合整数规划问题。第五步才是跑完整算例、做对比分析。这个顺序看起来很保守但效果很稳。基础运行优化模型跑通了说明网络约束和负荷数据没有问题再加P2G耦合模块出问题时就可以把排查范围锁定在P2G相关约束上而不是在整个模型里大海捞针。2. 核心建模方法把论文公式拆成能编程的模块2.1 P2G厂站建模核心就三个变量论文中常见的P2G简化模型是一个线性转换关系g_p2g(t) eta_p2g * p_p2g(t)p_p2g(t)是时段t消耗的电功率单位用MWg_p2g(t)是注入天然气系统的等效热功率也用MW作为统一能量单位eta_p2g是综合转换效率数值多在0.55到0.75之间。实际工程里的效率会受到电解槽负荷率、温度等因素影响但规划类论文普遍采用线性效率原因是规划阶段关注的是投资决策和年度运行成本不必精细到每个工况点。P2G运行约束还需要搭配容量上限与建设逻辑0 p_p2g(t) cap_p2g * x_p2gx_p2g是0-1建设变量取1表示在该节点建设P2G厂站cap_p2g是建设容量。这个约束保证没有建设厂站的节点强制不出力。有些论文还会加爬坡约束限制P2G相邻时段输入功率的变化量但在规划阶段如果时间尺度是小时级爬坡约束的影响通常不大可以根据论文是否给出相关参数决定加不加。2.2 电网侧模型直流潮流是主流选择电网侧建模在规划类论文里基本都用直流潮流不研究电压和无功只关注有功功率分布。直流潮流的线性表达式为B * theta P_inj P_flow_ij (theta_i - theta_j) / x_ijB是节点导纳矩阵的虚部theta是节点相角P_inj是节点净注入有功。P2G在这个模型里体现为一个电负荷出现在对应节点的功率平衡方程中P_generator P_import - P_load - p_p2g sum(P_flow_out) - sum(P_flow_in)支路潮流有限幅约束相角也有上下限。这些约束打包进YALMIP时往往用矩阵形式一次性生成不需要逐条循环。直流潮流虽然忽略了很多电气细节但对于“在哪里建P2G、建多大”这种投资决策问题已经足够。规划模型里同时包含建设变量和运行变量求解规模本来就大如果电网侧再用交流潮流模型非线性会直接把求解难度拉到另一个量级。2.3 气网侧模型Weymouth方程怎么处理是关键气网侧的稳态管流方程通常用Weymouth方程描述f_pipe^2 K * (p_i^2 - p_j^2)f_pipe是管道流量K是管道常数p_i和p_j是两端节点压力。麻烦在于这是个非线性约束直接扔给求解器几乎不可能在可接受时间内找到全局解。论文里常用的处理手段有三种分段线性化、增量线性化、以及用变量替换后做迭代逼近。我复现时采用的是分段线性化方法把管道流量按照不同区间分段逼近每一段用线性不等式描述气压平方差与流量的关系。这样做虽然会增加约束数量但换来了模型整体线性可以直接交给MILP求解器处理。气网节点流量平衡则写成g_source g_p2g g_load g_gt sum(pipe_flow_out) - sum(pipe_flow_in)g_source是气源供气g_p2g是P2G注入g_load是气负荷g_gt是燃气轮机耗气。气压上下限也要作为约束加入到模型里。2.4 目标函数与整体约束结构规划模型的目标函数由三块构成P2G投资年值、系统运行成本、弃风惩罚。投资年值通常用等年值法折算C_inv CRF * sum(inv_cost_p2g * cap_p2g) CRF r * (1r)^n / ((1r)^n - 1)r是折现率n是项目寿命CRF是等年值系数。运行成本包括购电成本和购气成本购电按分时电价计算购气按单位热值气价计算。弃风惩罚项可以加一个较小的惩罚系数也可以直接作为硬约束要求风电全部消纳具体看论文给的场景假设。完整约束结构可以整理成下面这张表模块典型约束变量类型电网节点功率平衡、支路潮流限幅、相角上下限连续气网节点流量平衡、管道流量方程线性化、气压上下限连续P2G耦合电功率与气功率转换、容量上限、建设逻辑连续二进制燃气轮机发电功率与耗气量转换、出力上下限连续3. Matlab代码一步步搭起来从变量定义到求解器调用3.1 工具箱怎么选YALMIP建模是复现最快路线复现这类规划问题最重要的一步就是选对建模工具箱。我的建议是用YALMIP做建模层配一个成熟的MILP求解器做求解层。原因很简单模型里同时有0-1变量、多时段变量、复杂耦合约束如果用MATLAB自带的优化函数手写标准形式变量顺序、等式约束矩阵、不等式约束矩阵稍微错一个位置就检查半天。YALMIP允许直接用数学表达式的形式写约束和变量代码的可读性和调试效率都要高很多。如果你的电脑还没装YALMIP安装方式也不复杂去官网下载压缩包解压后把文件夹加到MATLAB路径即可。求解器方面常见商业求解器和开源求解器都可以配合YALMIP工作配置好求解器路径后通过sdpsettings指定求解器名称就能调用。3.2 代码结构四个模块各司其职我的复现工程分成四个脚本模块参数初始化脚本、变量定义脚本、约束构建脚本、结果输出脚本。参数初始化脚本只负责给数据变量定义脚本负责创建决策变量约束构建脚本把所有约束汇总成一个约束对象结果输出脚本负责统计成本和绘制曲线。这种模块化写法在复现论文时非常有用。换一套算例参数时只需要改参数初始化脚本约束脚本完全不用动如果模型加一个约束也只需要在约束构建脚本的对应位置改一行不用在几百行代码里找人。3.3 P2G耦合约束的YALMIP实现骨架下面这段代码展示的是P2G耦合模块的核心写法我把维度说明放在注释里%% 参数初始化示意 T 24; % 典型日时段数 N_BUS 6; % 电网节点数 N_GAS 6; % 气网节点数 cand [2 3 5]; % P2G候选节点集合 eta_p2g 0.62; % 转换效率 inv_cost 5000; % 单位容量投资成本元/kW discount 0.08; % 折现率 life 20; % 项目寿命年 CRF discount * (1discount)^life / ((1discount)^life - 1); cap_max 300; % 单节点最大建设容量MW %% 决策变量定义 x_p2g binvar(N_BUS, 1, full); % 建设决策变量 cap_p2g sdpvar(N_BUS, 1, full); % 建设容量变量 p_p2g sdpvar(T, N_BUS, full); % P2G逐时段消耗电功率 g_p2g sdpvar(T, N_GAS, full); % P2G逐时段注入气网热功率 %% 约束构建 Constraints []; for i cand % P2G转换关系 Constraints [Constraints, g_p2g(:, i) eta_p2g * p_p2g(:, i)]; % 运行功率不超过建设容量 Constraints [Constraints, 0 p_p2g(:, i) cap_p2g(i)]; % 容量上限同时把容量与建设变量联系起来 Constraints [Constraints, cap_p2g(i) cap_max * x_p2g(i)]; % 未建设节点强制容量为0 Constraints [Constraints, cap_p2g(i) 0]; end % 候选节点之外的节点不允许建设 for i setdiff(1:N_BUS, cand) Constraints [Constraints, x_p2g(i) 0, cap_p2g(i) 0]; Constraints [Constraints, p_p2g(:, i) 0, g_p2g(:, i) 0]; end这段代码里有几个点容易写错。第一p_p2g(:,i)的下标方向必须是T乘N和后续电网节点平衡约束保持一致第二容量上限cap_max和效率eta_p2g的单位要先统一我习惯把所有能量单位都统一成MW天然气的热功率也按MW计这样g_p2g和p_p2g可以直接用效率公式换算不用再牵扯立方米、兆焦耳等单位换算第三setdiff那一段约束很多人容易漏掉少了它非候选节点虽然容量为0但p_p2g仍然可能因为求解器边界设置被赋予非零值导致结果里出现不存在厂站的P2G出力。3.4 目标函数与求解调用目标函数在YALMIP里可以直接写成表达式求和%% 目标函数 C_inv CRF * sum(inv_cost * cap_p2g); % 投资年值 C_elec sum(price_elec * p_import); % 购电成本 C_gas price_gas * sum(sum(g_buy)); % 购气成本 C_penalty penalty_coef * sum(sum(p_wind_curtail)); % 弃风惩罚 obj C_inv C_elec C_gas C_penalty; %% 求解 ops sdpsettings(solver, gurobi, verbose, 2); result optimize(Constraints, obj, ops);运行之前要检查YALMIP是否已经正确识别求解器可以在MATLAB命令行输入yalmiptest看输出里有没有出现对应的MILP求解器。如果求解器没有配置好YALMIP可能自动降级成内置求解器遇到0-1变量时速度会慢很多甚至可能直接报错。求解完成后用value(x_p2g)、value(cap_p2g)提取规划结果用value(p_p2g)提取逐时段运行结果。注意提取变量值的时候不要对整个sdpvar对象直接画图先value转成double矩阵再画否则很容易遇到类型不匹配的问题。4. 一个可复现的小算例跑一遍完整流程4.1 算例结构与参数设定为了验证整个代码流程我自建了一个小规模的模拟测试系统结构是6节点电网加6节点气网两个网络在拓扑上一一对应P2G候选节点选在2、3、5三个电气节点上。电网侧包含1个常规火电机组、2个风电场和2个燃气轮机气网侧包含2个气源和若干气负荷。风速曲线采用典型日数据节点负荷按峰值负荷比例缩放。这里需要说明算例参数是在论文常见量级基础上自己拼的不代表任何真实工程数据只是用来演示复现流程。具体参数如下参数取值典型日时段数24P2G转换效率0.62P2G单位投资5000元/kW折现率8%项目寿命20年风电总装机220 MW系统峰值电负荷380 MW系统峰值气负荷150 MW4.2 复现结果与图表解读按上述参数跑出来的规划结果是这样的候选节点2建设P2G容量30MW候选节点5建设容量16MW候选节点3未建设。从系统层面看加入P2G后总成本从无P2G场景的约5280万元/年下降到约5015万元/年可再生能源利用率从91%提升到97%。逐时段出力曲线上可以清楚看到P2G的耗电功率集中在凌晨1点到5点之间也就是风电出力高、电网负荷低的时段。加入P2G以后风电消纳量明显提高燃气轮机的耗气量也略有下降因为部分气负荷由P2G注入的气体直接满足了。这些结果的合理性可以从经济逻辑上验证P2G厂站选址落在风电渗透率高的节点是因为风电场周边电能不能就地消纳时有较大弃风压力P2G相当于给风电装了一个“接盘方”容量大小则由峰值负荷和风电曲线共同决定容量再加大多出来的P2G设备利用率太低边际收益反而不划算。4.3 算例结果对不上论文怎么办复现时最容易急的一个问题就是结果和论文对不上。我的判断原则是先看量级是否一致再看趋势是否一致最后才看数值是否完全一致。规划类论文里差几个百分点的数值是很正常的因为论文复现需要的数据往往没有全部公开负荷曲线、风速曲线的实际来源、典型日选取方式、收敛精度这些细节都会导致结果浮动。量级和趋势如果一致就可以认为模型本身复现成功。如果连量级都对不上优先检查单位折算系数尤其是天然气流量和热功率之间是否做了正确的能量换算。我见过不少复现失败案例最后都是卡在气网单位上一个换算系数搞错整个成本就差了几倍。5. 复现过程中容易卡住的问题与排查方法5.1 WEYMOUTH非线性导致求解器不可行最典型的坑是管道潮流方程直接写成sqrt形式导致模型变成非线性混合整数问题求解器要么报错要么长时间不收敛。我的处理方案是用分段线性化逼近虽然约束数量增加但整体模型保持MILP结构求解稳定性和速度都有保障。分段数建议取5到8段太少精度差太多求解时间增长明显。5.2 单位不统一导致目标函数数值扭曲电力系统的MW、气网系统的m³/h、价格系统的元/kWh、元/m³各自混在一起结果很容易是天文数字或者无限接近零。统一能量单位是关键。我强烈建议把天然气的能量流量统一成MW热功率气价也换算成元/MW这样目标函数里所有项都是同一个量纲数值范围可控求解器数值稳定性也会好很多。5.3 “没有建P2G的节点还在出力”这种逻辑错误这种问题通常出在建设变量x_p2g与运行容量cap_p2g之间的约束缺失。一定要确认每个候选节点都建立了运行容量和建设变量的联动关系非候选节点则显式强制为零。我会在约束构建完成后加一段诊断代码打印所有节点x_p2g和cap_p2g的取值快速发现这种逻辑错误。5.4 求解时间过长或内存溢出规划模型的时间维如果取8760小时再叠加0-1变量求解规模会膨胀得非常快。论文里常见做法是用典型日替代全年数据比如取24小时或几个典型日加权重。如果你必须跑全年数据建议先跑典型日验证模型正确性再扩展到全年。求解时间超过预期时先检查是不是求解器配置没设对再考虑缩小候选节点集合、减少分段线性化段数。5.5 结果不合理总是选择在某一个节点建P2G如果模型输出里所有候选节点都建了P2G或者全部候选节点都不建通常不是代码问题而是目标函数或约束设置有问题。比如投资成本参数太小模型就会倾向于到处建P2G运行成本里购电价格太低模型就会倾向于直接买电而不是转化。做一次参数灵敏度分析看容量选址结果随成本参数如何变化就能定位是哪个成本项权重不合理。这里把常见问题整理成速查表现象可能原因排查建议模型不可行约束写错、单位不一致先锁定基础运行模型逐步加模块定位结果出现非零容量但x_p2g0容量联动约束缺失检查容量上限约束和建设变量关联求解时间过长求解器未配置、MILP规模过大优化求解器配置用典型日替代全年时段选址结果异常集中成本参数设置不合理检查投资和运行成本相对量级P2G出力曲线与预期不符效率或爬坡约束设置偏大核对转换效率和输入输出方向6. 我的实际体会P2G规划复现中最该重视的三件事6.1 建模公式与代码要双向对照复现质量的底线是“论文公式写得出来代码片段对得上”。我建议每写完一组约束就在注释里标注它对应论文里的哪个公式编号。这样做听起来麻烦但在调试时节省的时间是成倍的。P2G规划模型涉及三个子系统的联动任何一个模块写错最后输出的结果都可能是“看起来合理但不正确”因为目标函数还会收敛只是收敛到一个错误的模型上。6.2 求解结果必须做经济性解释很多同学跑出结果就完事但论文复现里更重要的能力是解释结果。为什么P2G建在风电节点而不建在负荷中心为什么夜间出力高峰比白天更明显这些问题的答案往往藏在数据里比如分时电价的低谷段、风电高发时段与气负荷曲线的错峰关系。把结果讲清楚才算真正理解了这个课题。6.3 后续扩展可以往哪些方向走复现完成的代码框架是一个很好的扩展平台。我目前在这个代码上继续改了三个方向一是把单一典型日换成多个典型日加权重做更接近全年的规划二是加入多目标函数同时考虑成本和可再生能源利用率两个指标三是把P2G效率从固定值改成负荷率相关的分段函数让运行结果更贴近实际。这些扩展在现有代码上改动量都不大但研究深度和论文档次会有明显区别。最后再分享一个个人习惯每次改动完成后我会把关键结果图表和参数组合保存成记录文件方便随时回溯是哪一组参数跑出了哪个结果。复现不是一次性的任务算例、参数、场景会经常变换有完整过程记录再回头重跑或调整会从容很多。