1. 为什么选这个课题P2G在电-气综合能源系统中的关键位置最近帮一位师弟做硕士论文复现题目就是“计及P2G厂站的电-气综合能源系统规划研究Matlab代码实现”。说实话我刚拿到这个题目的时候脑子里冒出来的第一个问题是P2G到底解决了什么问题值得单独撑起一篇硕士论文P2G全称是Power to Gas中文叫电转气。简单地说就是把电力网络里多余的电能通过电解水变成氢气再进一步甲烷化生成天然气然后注入天然气管网。这个环节最直接的价值是给两个能源系统之间装了一座“双向桥”电网电多了可以往气网里存气网不够用了又能把气转回电通过燃气轮机。很多论文里把P2G当成综合能源系统的核心耦合环节并不是因为它多高科技而是因为它让原本割裂的电力平衡和天然气平衡变成了一个可以互相调节的整体。这类课题适合谁我觉得主要是这几类人一是正在选题的电气工程、能源系统工程方向的硕士生二是做新能源消纳、多能互补研究的工程师三是想验证自己优化建模能力、想看看Matlab和求解器怎么配合使用的初学者。不同基础的人读这个题目的收获完全不一样。这个课题真正吸引我的地方在于它不是在做一个孤立的设备模型而是把电力网络、天然气网络、P2G厂站、风电出力、负荷曲线全部放进一个优化框架里。写代码的时候要用到Matlab的矩阵操作要把物理约束写成数学约束还要选择一个能算得动优化问题的求解器。整套流程走下来基本等于把“电力系统分析、燃气系统基础、最优化理论、Matlab编程”这四门课串起来做了一次实战。而且我注意到现在很多人复现论文最大的障碍不是看不懂公式而是不知道怎么把一个优化模型从论文里的数学表达式变成能在Matlab里跑出结果的代码。学科交叉的地方尤其如此电力专业的同学不熟悉天然气网络燃气专业的同学又对潮流方程和fmincon陌生。这篇博文我就按我实际复现的时候的步骤把模型怎么建、代码怎么组织、坑在哪里一五一十讲清楚。复现这个课题的产出也挺直观你会得到一套完整的Matlab程序输入电网拓扑、气网拓扑、风电和负荷数据输出的是最优规划方案包括P2G厂站的容量、要不要新增线路、典型的运行出力安排以及总成本。这套东西无论用来发论文、写大作业还是把它改造成自己课题的基础模块都是很有价值的。2. 模型搭建从目标函数到约束条件的完整链条2.1 系统框架与决策变量设定先别急着写代码建模之前必须把系统是什么样子画清楚不然后面的约束条件很容易越写越乱。这个课题里的“电-气综合能源系统”一般包含三块电力系统、天然气系统、耦合环节。电力系统里有常规机组、风电场、电力负荷、输电线路如果做规划还要考虑新增线路和变电站扩容。天然气系统里有气源节点、天然气负荷、输气管道、压缩机如果需要还有储气库。耦合环节就是P2G厂站它既是电力系统的负荷耗电也是天然气系统的注入源产气。我复现时采用的是比较常见的“规划-运行联合优化”框架外层决定投资建设方案也就是建多大容量的P2G、要不要新增某条线路内层在每个典型日里做运行调度。这类问题在数学上通常被建模成一个混合整数优化问题整型变量管“建还是不建”连续变量管“功率和流量怎么分配”。决策变量要分两类来说。投资类变量是P2G建设容量、线路增容决策、气网新建管道决策这些变量直接决定了系统的物理结构。运行类变量是每个时刻的常规机组出力、风电消纳量、P2G耗电功率和产气功率、气源购气量、节点压力等这些变量是在给定结构下让系统以最低成本运行的结果。从硕士学位论文的角度这个框架已经足够撑起完整的研究逻辑了。如果是我自己搭建我会在论文里专门用一个表格列出所有变量及其单位、类型、维度避免后面建模的时候混单位。这里有一点必须提醒单位是复现过程中最容易翻车的地方后面我会单独说。2.2 目标函数与成本建模绝大多数这类论文的目标函数都是总成本最小化。总成本通常由三块组成投资建设成本、系统运行成本、惩罚项。投资成本这一块核心是P2G厂站单位容量投资成本乘以规划容量然后是电力网络扩容、天然气管道扩容的一次性投入。学过工程经济的都知道一次性投资不能直接加到年运行成本里要把它们折算成年值用到一个年值折算系数C_annual C_total * r * (1r)^N / ((1r)^N - 1)其中r是折现率N是设备寿命。公式的原理很简单相当于把未来N年要分摊的等额成本折算到每年。很多初学的同学直接拿总投资和年运行成本相加这样在数学上能求解但是规划结果会偏向少建设备因为投资被高估了。我建议还是按年值折算来做这是论文里能体现专业度的地方。运行成本主要是购电成本和购气成本。购电成本按照电网节点电价或上网电价乘以购电量购气成本是气源价格乘以购气量。如果系统里已经明确了火电机组的煤耗也可以把煤耗成本写进去。风电和光伏在运行阶段的边际成本可以近似为0但在目标函数里可以加一个对弃风弃光的惩罚项。目标函数的写法看起来不复杂但我遇到过一个问题如果直接把投资成本的量纲和运行成本加在一起数值会差个数量级。每千瓦P2G投资可能几千元每兆瓦时电价只有几百元算出来的目标函数值完全被投资项覆盖求解器对运行部分的优化就非常不敏感。我的处理办法有两种一是把投资成本除以统计年限后按小时折算到每个时段再统一相加二是通过权重系数调节让投资项和运行项处于相近的数值范围。这两种方法都能跑通但建议在论文里明确说明折算方式否则盲审容易被打回。2.3 约束条件与P2G环节的建模细节这个课题最核心的约束条件就是电网潮流平衡、气网流量平衡、P2G转换关系、设备运行上下限。电力网络部分如果做直流潮流每个节点的有功平衡约束可以写成该节点注入的有功机组出力、风电场出力、从上级电网购入功率减去该节点流出的部分负荷、P2G耗电功率、注入其他节点的线路潮流等于零。线路潮流本身受传输极限约束。直流潮流的本质是忽略了无功和电压只考虑有功对大多数综合能源规划研究完全够用。天然气网络部分要复杂一些。每个节点的流量平衡包括气源注入、P2G注入、气负荷以及管道流入流出。管道流量和两端节点压力之间通常用Weymouth方程描述这个方程是非线性的会让问题变成非线性规划。处理方式无非两种一是直接用非线性规划求解器比如Ipopt硬解二是把管道流量方程在某个工作点附近线性化转成线性约束。论文往往采用第二种用分段线性逼近来处理管道气压关系。P2G环节的建模是整个模型的灵魂。最简单的建模方式是两个功率流关系和一个容量约束电解槽耗电功率 P_p2g(t) 和产氢速率成正比氢气再进入甲烷反应器生成天然气P2G产出气体的热值功率 G_p2g(t) η_p2g * P_p2g(t)P2G的耗电功率上限由装机容量决定P_p2g(t) ≤ S_p2g。η_p2g 是整个电转气链路的综合效率一般取0.5到0.7具体要看论文中假设的是从电力到氢、还是从电力到甲烷。这里有一个极其容易踩坑的单位问题电力功率单位是kW或MW天然气流量的单位如果是m³/h或者kg/h两者无法直接相加。正确做法是把天然气流量按热值换算成功率单位。比如1标准立方米天然气低位热值大约是35.6MJ那么1m³/h的天然气流量大约对应9.89kW的热功率。这个问题不搞清楚气网平衡方程写出来就是错的再怎么调求解器也得不到合理结果。除开这些还要考虑风电出力上限、常规机组的爬坡约束、气源的供气能力上限、节点压力上下限。如果研究里考虑了储能还要再加一个储气罐的容量约束和充放平衡约束。说实话每一条约束单独拿出来都不难但二三十条约束堆在一起符号标错一个整个模型就解不出来了。3. Matlab实现求解器选择、代码架构与关键函数3.1 为什么用Yalmip而不是裸写fminconMatlab里实现优化问题有几条路最直接的是直接用fmincon或linprog。但我的经验是这种复杂模型不要裸写fmincon。为什么因为fmincon要求你把目标函数和约束条件都写成标准的函数形式尤其是不等式约束要写成C(x) 0你得手工把一堆矩阵和向量拼起来后期改一个变量维度都要回去翻代码效率非常低。我推荐的方式是使用Yalmip工具箱。Yalmip本质上是一个建模语言它的作用是让你用接近数学表达式的语法定义变量和目标函数然后由它自动转换成求解器需要的格式。你只需要负责把模型写对具体矩阵生成和调用求解器的脏活累活都交给Yalmip。用Yalmip的另一个好处是求解器可以灵活切换。同一个模型今天你可以用Gurobi解线性规划明天想试试Ipopt解非线性规划只需要改一行 sdpsettings 的solver参数不需要重写模型。这对复现论文特别友好因为你经常需要在不同求解设置之间做对比确认解不是求解器本身的伪最优。当然还有纯Matlab的优化工具箱比如intlinprog如果模型是MILP也能用。但我个人在实际操作中还是更习惯Yalmip为外层框架再根据问题类型指定底层求解器。下面是一段我最常用的主程序框架基本复用了好几次这个课题%% main.m 主程序框架 clear; clc; close all; % 读取数据 mpc load_case(case9); % 电网数据 gas load_gas_data(gas6); % 气网数据 profile load_load_profile(); % 负荷和风电时序数据 param set_param(); % 经济参数 % 构建优化模型 [model, cons, obj] build_ies_model(mpc, gas, profile, param); % 配置求解器 ops sdpsettings(solver, gurobi, verbose, 2, savesolveroutput, 1); % 求解 result optimize(cons, obj, ops);这里每个函数各司其职数据读取、模型构建、求解、结果显示分开调试的时候非常方便。可以这么理解你不需要把两千行代码堆在一个脚本里要容忍自己在建模过程中反复调整所以代码结构清晰比什么都重要。3.2 构建模型的核心Matlab代码模型构建函数是整套代码最需要仔细的地方。我用Yalmip定义变量和约束的方式如下%% build_ies_model.m 简化版核心代码 % 时间尺度T个时段 T 24; % 投资决策变量 S_p2g sdpvar(1, 1, full); % P2G建设容量(MW) x_line binvar(1, n_line, full); % 线路增容0-1变量 % 运行决策变量 P_import sdpvar(1, T, full); % 从上级电网购入功率 P_wind sdpvar(1, T, full); % 风电消纳功率 P_p2g sdpvar(1, T, full); % P2G耗电功率 G_p2g sdpvar(1, T, full); % P2G产气功率(热值功率,MW) G_import sdpvar(1, T, full); % 气网购入气量(热值功率,MW) % 目标函数投资年值 运行成本 C_p2g_inv c_p2g_unit * S_p2g * crf(r, N); C_line_inv sum(c_line * d_line .* x_line) * crf(r, N); C_op sum(C_gas_price * G_import) sum(C_elec_price * P_import); objective C_p2g_inv C_line_inv C_op; % 约束条件 cons []; % 电功率平衡 for t 1:T cons [cons, P_import(t) P_wind(t) P_thermal(t) ... P_load(t) P_p2g(t)]; % P2G转换关系 cons [cons, G_p2g(t) eta_p2g * P_p2g(t)]; % P2G容量限制 cons [cons, P_p2g(t) S_p2g]; % 风力上限 cons [cons, 0 P_wind(t) P_wind_avail(t)]; end % 气网节点平衡热值单位 for t 1:T cons [cons, G_import(t) G_p2g(t) G_load(t)]; end % 求解 ops sdpsettings(solver, gurobi, verbose, 2); sol optimize(cons, objective, ops);小提示上面这个是非常简化的网络模型没有考虑多节点潮流和天然气管道压力。你真做论文复现的时候肯定要加上节点矩阵形式用矩阵运算代替for循环否则48个时段、几十个节点的问题会跑得非常慢。我的经验是用向量化写法建立一个平衡节点映射矩阵用矩阵乘的方式来代表每一个节点的功率和流量平衡。Yalmip里sdpvar是可以直接用矩阵参与的效率比for循环高太多。我第一次复现这个课题的时候用的是for循环7节点系统还能忍换到9节点气网和24时段就明显吃力改写成矩阵形式之后计算速度快了不止两倍。3.3 结果可视化和多方案横向对比优化求解完成之后最大的问题是怎么把结果讲清楚。评审老师最关心三个东西系统配置变化、运行特性曲线、经济性提升。系统配置变化直接看决策变量取值比如P2G容量建了多少、哪条线路需要增容。这个可以直接输出成表格我一般会加一列对照参考把“无P2G”和“含P2G”两个方案并列。运行特性曲线建议画两张图。第一张是24小时风电出力、风电消纳、P2G耗电功率的关系时序曲线这张图能直观看出P2G主要在哪个时间段吃掉了原本要弃掉的风电。第二张是气网侧的气源供气、P2G产气、气负荷平衡堆栈图这张图能说明P2G对气网调峰的作用。这两张图画完整个方案的运行逻辑基本就一目了然了。经济性对比则用表格整理总投资、年运行成本、弃风率、系统总成本。如果你额外做了敏感性分析比如气价发生变化对P2G规划容量的影响可以画一条P2G最优容量随气价变化的折线图这张图是我认为全篇论文里最有说服力的东西。实际上我在复现里做过三方案对比一种是不装P2G、一种是装P2G但不扩容电网、一种是装P2G同时优化电网扩容。结果并不意外不装P2G的时候弃风率最高装了P2G之后弃风率明显下降但只有同时考虑电网扩容P2G的利用率才真正提升。这三行结论虽然简单却足以体现规划研究的价值不是单独看某一项技术而是看整个系统怎么配合。4. 复现过程中的问题排查与心得4.1 最容易踩的坑这个课题我前后完整跑通不下三遍踩过的坑很多挑最典型的几个分享。第一个坑是单位不统一。电网功率用的是MW天然气流量用的是m³/h两者直接在一个平衡约束里相加求解器没报错但结果乱七八糟。我后来把整个模型全部统一成热值功率气流量乘以热值系数再参与计算这才算对上。这个错误的特点是模型能求解约束也满足但结果与论文数值完全对不上。第二个坑是Yalmip的变量维度问题。S_p2g定义为1×1变量但在循环里给约束的时候如果写的是 cons [cons, P_p2g(t) S_p2g(1,t)]而S_p2g是1×1就会出现维度不匹配。我建议从一开始就统一所有变量的维度含义尽量定义成 1×T 的sdpvar然后用括号索引单个时刻这样最不容易错。第三个坑是求解器对非线性约束的处理。如果保留Weymouth方程模型成了MINLP用Gurobi或CPLEX直接报错因为这两个求解器不支持非线性。我最初不知道这一点折腾了很久才发现要嘛换Ipopt要嘛对管道方程做线性化。在硕士论文场景下我倾向于做简化线性化用分段线性函数近似然后用Gurobi直接解MILP稳定性高很多。第四个坑是二进制变量与连续变量的初始化问题。intlinprog或Gurobi在解混合整数问题时对初始解不敏感但如果你想让它在合理时间内收敛最好自己加一个启发式初值。我给x_line赋初值的时候是读取现有无P2G方案的结果把它当成初始可行解实测能减少不少求解时间。下面是问题速查表直接抄作业常见现象可能原因解决办法模型求解后目标值离谱单位不一致统一转换为热值功率MW求解器报错“No suitable solver”模型含非线性约束但指定了线性求解器换用ipopt或线性化长时间不收敛二进制变量多、初值差提供可行初值、减少典型日数量约束全满足但弃风率反而上升P2G容量和线路扩容耦合没有建模检查容量约束与运行约束的索引对应关系从论文复现但结果对不上参数取值不同核对折现率、设备寿命、气价、效率取值4.2 运行时间和收敛性优化四五十行约束的模型如果只有连续变量Gurobi算起来是秒级。但一旦加上二进制变量尤其是有几十条线路增容决策、多个P2G站址选择问题规模会迅速膨胀求解时间可能从几秒涨到几十分钟。我的优化手段主要有三个。第一是削减典型日不要用365天调度而是通过K-means聚类选出几个典型日每个典型日赋予权重代表整年的运行情况。我试过用4个典型日就能稳住整体趋势计算速度提升非常明显。第二是去掉冗余约束比如有些时段风电出力上限已经确定了P2G最大耗电没必要保留部分宽松的线路容量约束。第三是给目标函数设置一个合理的可行性容差sdpsettings里可以调整gurobi的MIPGap默认1e-4有时过严调成1e-2之后求解速度快非常多结果差异基本可以忽略。如果你像我一样是在普通笔记本上跑复现建议用MIPGap0.01作为默认配置先得到一套近似最优把全部流程跑通之后再把精度调高做最终计算。这个节奏能避免你花了一个小时等求解器收敛结果发现模型里还有个小错误。4.3 我的实操心得这个课题复现下来我个人最大的感受是不要把建模和代码当成两件事。建模的时候就要想着代码怎么写比如变量维度和单位在做公式推导的时候就直接确定好。反过来写代码的时候也要不断回到公式上确认每个约束是不是真的代表物理过程。单纯对着论文抄公式或者单纯照着网上代码改都会让你迷失。还有一点就是这个课题特别适合扩展。你只需要改动很少的地方就能把它延伸成计及碳交易的电-气综合能源系统低碳规划、考虑源荷不确定性的两阶段随机规划、或者多目标帕累托优化。我后来就试过在目标函数里加入碳排放惩罚项整个模型的求解逻辑基本不变但研究视角立刻不一样了。如果你准备拿相关方向继续做研究这个复现题目是一个很好的起点。最后再分享一个小技巧在跑完结果之后一定要把节点电价、气源出力和P2G出力这三条曲线放到同一张图里看。如果P2G耗电高峰恰好落在气价低谷或者弃风高峰时段说明模型传出来的信号是合理的再往下展开分析就有底气了。否则就算求解器说收敛了你也得怀疑是不是某个时序数据错位了。