0 引言

能源综合利用能够有效提高能源利用效率[1]，在能源转型的迫切需求下，能源互联网[2]和综合能源系统[3]等是实现能源综合利用的有效方式。燃气轮机（gas turbine, GT）和电转气（power-to-gas, P2G）技术的快速发展使得电力系统和天然气系统的耦合程度不断加深[4]，近年来电-气综合能源系统（integrated electricity and natural gas systems, IEGS）关键技术[5-8]受到了高校和科研机构的广泛关注。

潮流计算是IEGS运行安全分析的基础环节，目前已有不少文献进行了相关研究。借鉴电力系统稳态潮流计算的思路，文献[7]对天然气系统和能源集站分别建模，并基于牛顿拉夫逊法实现了IEGS稳态潮流的分立求解，文献[9]基于牛顿拉夫逊法进行了IEGS的统一潮流稳态分析。目前，用于IEGS潮流计算的核心算法仍为传统的牛顿拉夫逊法，其用于电力系统潮流计算时收敛性良好，但是用于天然气系统潮流计算时存在初值敏感问题。

在IEGS中，电力系统和天然气系统联合运行时，两个系统中的扰动会通过耦合元件相互影响，甚至会威胁整个耦合系统的安全运行[10-11]。文献[10]基于IEGS稳态潮流计算分析了两个子系统的交互作用。文献[12]指出天然气系统中的扰动会传播到电力系统并影响电力系统的安全运行。文献[13]指出电力系统中燃气轮机天然气需求的频繁变化会影响天然气系统的安全。因此，IEGS的相互作用机理分析对于保障耦合系统的安全运行具有重要意义。上述研究侧重从能量流动的角度分析系统的相互作用，没有定量分析系统中的不确定因素对子系统的影响，也没有对系统薄弱环节进行辨识。文献[10]基于统一潮流模型提出了IEGS的静态灵敏度分析方法，通过节点气压-节点注入功率灵敏度矩阵辨识了系统的薄弱环节。文献[14]采用灵敏度分析法定义了两个用于电、气系统交互作用机理分析的灵敏度矩阵，辨识了系统薄弱环节并给出了提升系统安全的措施。上述研究的灵敏度矩阵只能描述两个变量间的交互关系，无法评估耦合系统中的多种扰动对系统运行安全的影响。

针对上述问题，本文首先建立了电力系统、天然气系统和耦合元件的稳态模型，从而建立电-气综合能源系统的统一潮流模型，并基于改进的牛顿拉夫逊法实现潮流计算。基于统一潮流模型和灵敏度分析法，推导出关于电压/气压灵敏度矩阵的通用形式，并定义节点电压/气压-注入功率综合灵敏度指标，从而定量分析了电-气耦合系统的相互作用机理。最后，通过仿真算例分析了GT和P2G耦合作用下的电-气耦合系统的相互影响，从而为辨识系统薄弱环节提供理论指导。

1 考虑双向能量流动的IEGS稳态建模

本文研究的电-气综合能源系统考虑电、气子系统间的双向能量流动，如图1所示，由电力系统输电网、天然气输气系统、燃气轮机以及P2G装置组成。

图1 电-气综合能源系统结构示意图
Fig.1 Schematic diagram of IEGS

1.1 IEGS变量分类

对电-气综合能源系统进行建模时，其变量主要可以分为控制变量和状态变量两类，具体变量参见表1。电-气综合能源系统的状态变量描述了系统的节点运行信息，是安全分析和控制的主要变量。

表1 电-气综合能源系统变量分类
Table 1 Variable classification in IEGS

系统 节点类型 控制变量 状态变量电力平衡节点，PV节点，PQ节点平衡节点的电压幅值和相角，PV节点的电压幅值和有功功率，PQ节点的有功和无功功率平衡节点的有功和无功功率，PV节点的电压相角和无功功率，PQ节点的电压幅值和相角天然气 定气压节点，定流量节点定气压节点的气压，定流量节点的流量定气压节点的流量，定流量节点的气压

1.2 耦合元件模型

本文考虑的耦合元件为GT和P2G装置，两类耦合元件的能量转换关系[15]为

式中：PiGT 和PiP2G 分别表示电节点i 的GT产生和P2G消耗的电功率；fpGT 和fp P2G 分别表示气节点p 的GT消耗和P2G产生的天然气量；ηGp T和ηiP2G 分别为GT和P2G的转换效率；GHV 为天然气热值。

1.3 电力系统模型

本文采用经典的交流潮流模型[16]对电力系统进行描述，其节点注入功率方程和节点功率平衡方程如下：

式中：Pi 和Qi 分别为节点i 的注入有功和无功功率；Vi 为节点i 的电压幅值；θij 为节点i 和j 的相角差；Gij 和Bij 分别表示支路ij 的电导和电纳；ΔPi 和ΔQi 分别表示节点i 的有功和无功不平衡量；PiEG 和QiEG 分别表示节点i 非燃气发电机的有功和无功输出功率；PiES 和QiES 分别表示节点i 储电装置释放的有功和无功功率；QiGT 为节点i 燃气轮机的无功输出功率；PiED 和QiED 分别为节点i 负荷消耗的有功功率和无功功率。

1.4 天然气系统模型

天然气系统常用的潮流计算方法有节点法和回路法，本文采用节点法对天然气系统进行分析。

压缩机作为输气系统的重要元件，用于补偿长距离输送中的压力损失。对于不含压缩机的管道，管道流量与管道两端节点压力的关系满足一定的规律，输气系统一般为高压系统，因此本文采用Weymouth方程[17]对不含压缩机的管道流量和气压关系进行建模，具体表达式见式（5）和（6）。对于含压缩机的管道，其管道流量与压缩机马力有关[17]，其关系表达式见式（7）。

式中：f Ppq为管道pq的气流量；Cpq为管道pq的管道系数，其值由管道长度、直径、效率系数等参数共同决定；πp为节点p的气压；f Cpq为压缩机管道pq的气流量；Hp、Zp、κp 和Bp分别为入节点为p的压缩机的马力、气体可压缩系数、比热容比和压缩机常数。

压缩机通常有三种运行模式，定入口压力、定出口压力和定压缩比，本文选择定压缩比的运行模式。当压缩比为Kpq时，压缩机出、入口压力关系为

压缩机通常有燃气驱动型和电驱动型，本文考虑燃气驱动型压缩机，压缩机的马力与其消耗天然气量的关系为[18]

式中：表示管道pq的压缩机消耗的气流量；αp、βp和γp分别表示入节点为p的压缩机的天然气消耗系数。

考虑GT和P2G装置，节点流量平衡方程为

式中：Δfp 为节点p 气流量的不平衡量；fpGS 和fpGD 分别表示节点p 的气源和气负荷；fpGST 表示节点p 储气装置释放的气流量。

2 基于牛顿下山法的IEGS潮流计算

2.1 统一潮流建模

本文将电力系统的潮流建模方法拓展到IEGS，将电、气系统的潮流方程统一建模，联立求解。

在潮流计算中，电力系统的发电机节点分为平衡节点和PV节点，而天然气系统的气源节点也分为定气压节点和定流量节点，其中定气压节点类似于电力系统的平衡节点，可以平衡气流量。

考虑GT和P2G两类耦合元件，当GT接于电力系统平衡节点时，它的有功功率是变化的，因此GT所需要的天然气量是变化的，从而影响天然气系统的运行状态；当P2G接于天然气系统的定气压气源节点时，该P2G的电功率需求随着天然气系统运行状态的变化而变化。而GT不接于电力系统平衡节点，P2G不接于天然气系统定气压气源节点时，多能流模型的式（1）和式（2）可以消去，即GT和P2G的出力为定值。

综上，考虑GT/P2G接于电力系统/天然气系统的平衡节点，IEGS的统一潮流模型为

式中：X表示系统状态变量向量，包括电压相角θ、电压幅值V、节点气压π和压缩机马力H；ΔF表示系统的不平衡量向量，其值与系统状态变量X有关；ΔP、ΔQ和Δf 分别表示有功、无功和气流量的不平衡量向量。

电力系统有1个平衡节点，EA个PV节点和EB个PQ节点，则式（11）中有EA+EB个有功功率不平衡方程，EB个无功功率不平衡方程。天然气系统中有G个节点（含GA个定气压气源节点），GB个压缩机，则式（11）中有G-GA+GB个气流量不平衡方程。

2.2 潮流求解

式（11）给出的统一潮流模型是一组非线性方程组，本文采用牛顿拉夫逊法[16]进行求解。IEGS统一潮流模型的修正方程为

式中：ΔX表示状态变量的变化量；J为统一雅可比矩阵，对角块表示子系统自身潮流与状态变量的偏微分关系，非对角块表示不同子系统间的耦合关系，下标代表变量关系，如JPθ为有功功率P与电压相角θ的偏微分关系。

对于电力系统，传统牛顿拉夫逊法具有良好的收敛性，其状态变量的迭代公式为

式中：k 表示迭代次数。

然而天然气系统存在对初值敏感的问题，用传统牛顿拉夫逊法进行求解时容易不收敛。因此，本文采用改进的牛顿拉夫逊法，即牛顿下山法计算天然气系统的潮流，其状态变量的迭代公式为[19]

式中：λ（k） 为第k 次迭代的步长修正因子，其具体计算过程见表2。

表2 步长修正因子的计算
Table 2 Calculation of the Down-hill factor

步长修正因子λ（k）的计算输入：λ（k） =1；ΔX == F (X（k） )2；1：计算R= X X（k）1；gas F1 gas gas（）（）（）k k k −λ Δ ；2：WhileF R gas gas gas 1( ) 2 ＞F ；（）（）k k λ λ=/2；（）（）（）k k k R= X X−λ Δ ；gas gas End输出：λ（k）。

3 IEGS灵敏度分析与综合灵敏度指标构建

3.1 基于统一多能流模型的灵敏度分析

灵敏度分析法可以量化分析在当前运行点的基础上，系统受到扰动或输入变化时其运行状态的变化，灵敏度指标可以用来描述影响程度。

在IEGS中，电、气子系统的相互影响是通过耦合元件传播的，一个子系统受到的扰动可能会传递到另一个子系统。因此，有必要评估IEGS中扰动在不同系统间的传播机理，定义统一的标准形式来量化耦合系统的灵敏度关系。

基于统一多能流模型，获取潮流收敛后的统一雅可比矩阵J，基于灵敏度分析法，保留状态方程的Taylor展开一次项可得：

式中：Δf2表示压缩机出节点气流量的偏差量；Δf1表示除去定气压气源节点和压缩机出节点的其他节点的气流量的偏差量，雅可比矩阵块也进行相应的分块；Δθ、ΔV、Δπ和ΔH分别表示电压相角、电压幅值、节点气压和压缩机马力的变化量。

在IEGS中，用于进行系统安全分析的状态变量主要是节点电压和节点气压。式（16）中，如果存在GT接于电力系统的平衡节点，则JPπ和JPH不为0，如果存在P2G接于天然气系统的平衡节点，则Jf1θ和Jf1V不为0。

3.2 节点电压的综合灵敏度指标

以节点电压-注入功率灵敏度为表征量，定义节点电压的综合灵敏度指标，采用高斯消去法消去式（16）中的Δθ、Δπ 和ΔH，可得到：

式中：M1至M6为中间变量。

除了综合灵敏度指标外，还可以获取关于节点电压的灵敏度矩阵：

式中：SVP为节点电压-注入有功功率的灵敏度矩阵；SVQ为节点电压-注入无功功率的灵敏度矩阵；SVf1、SVf2和SVf为节点电压-注入气流量的灵敏度矩阵；灵敏度矩阵的具体计算对照式（17）。

3.3 节点气压的综合灵敏度指标

以节点气压-注入功率灵敏度为表征量，定义节点气压的综合灵敏度指标，采用高斯消去法消去式（16）中的Δθ、ΔV 和ΔH，可得到：

式中：R1~R5 为中间变量。

除了综合灵敏度指标外，还可以获取关于节点气压的灵敏度矩阵：

式中：Sπ P 为节点气压-注入有功功率的灵敏度矩阵；SπQ 为节点气压-注入无功功率的灵敏度矩阵；Sπ f1、Sπ f2 和Sπ f 为节点气压-注入气流量的灵敏度矩阵；灵敏度矩阵的具体计算对照式（20）。

3.4 分析流程

在IEGS系统当前状态已知时，通过求解灵敏度矩阵和综合灵敏度指标，可以定量分析系统中的不确定因素对系统的影响。本文IEGS灵敏度分析的流程见图2。

图2 IEGS灵敏度分析流程图
Fig.2 Flowchart of sensitivity analysis for IEGS

综合灵敏度指标的值越大，说明该点的状态变量受系统扰动的影响越大，节点电压/气压受到的影响越大。其中，综合灵敏度指标值较大的节点即为系统薄弱环节。同理，灵敏度矩阵的元素值也体现了状态变量和负荷扰动间的影响程度。

4 算例分析

本文采用由IEEE 30节点电力系统[20]和15节点天然气系统[17]构成的IEGS进行算例分析，其中电节点1为平衡节点，气节点1和2为定气压气源节点，系统结构见图3。

图3 测试系统示意图
Fig.3 Schematic diagram of the test system

4.1 灵敏度分析

关于IEGS的灵敏度分析，设置以下两种扰动。

扰动1：仅考虑电力系统所有负荷的波动，波动为每个节点的负荷增加该节点负荷的2%，用于研究电力系统中的扰动对天然气系统的影响。

扰动2：仅考虑天然气系统所有负荷的波动，波动为每个节点的负荷增加该节点负荷的2%，用于研究天然气系统中的扰动对电力系统的影响。

分别从GT和P2G耦合作用下的相互影响出发，设置以下场景。

场景1：研究GT耦合作用下的影响机理，考虑三个GT的出力，不考虑P2G的接入，其中GT1接于电力系统平衡节点，GT总出力占所有机组出力的50%，计算有关气压的灵敏度矩阵和综合灵敏度指标。

场景2：在场景1的基础上，改变GT出力总占比，分别为30%、50%和70%，计算扰动1下的气压综合灵敏度指标。

场景3：在场景1的基础上，改变GT1的接入位置，分别为气负荷节点3、4和13，计算扰动1下的气压综合灵敏度指标。

场景4：研究P2G耦合作用下的影响机理，考虑P2G的接入，不考虑GT1的接入（另外两个GT正常接入），即电力系统平衡节点为常规发电机，设置GT出力总占比为30%，计算有关电压的灵敏度矩阵和综合灵敏度指标。

场景5：在场景4的基础上，在天然气系统薄弱负荷节点（气压综合灵敏度值高的节点）增加储气罐，容量为负荷的30%，计算扰动2下的电压综合灵敏度指标。

场景6：在场景1的基础上，在电力系统薄弱负荷节点（电压综合灵敏度值高的节点）增加储电装置，容量为负荷的30%，计算扰动1下的气压综合灵敏度指标。

场景1的计算结果如图4和图5所示。场景2和场景3的结果分别见图6和图7。场景4的计算结果如图8和图9所示。场景5和6的结果分别见图10和图11。

图4 GT作用下的气压灵敏度矩阵
Fig.4 Pressure sensitivity matrices under the coupling of GT

图5 GT作用下的气压综合灵敏度指标
Fig.5 Composite sensitivity index of pressure under the coupling of GT

从图4（a）可知，在气压-注入气流量灵敏度矩阵中，对角线区域的灵敏度值较大，其中末端负荷节点13、14和15的值最大，说明其气压受电负荷变化的影响较大。从图4（b）的气压-注入有功功率灵敏度矩阵和图5的气压综合灵敏度指标可以看出，天然气系统末端负荷节点13、14和15的气压受电力系统有功扰动的影响最大，为天然气系统的薄弱环节。

图6和图7分别为GT耦合作用下，不同GT容量占比和不同GT接入位置的气压综合灵敏度指标的计算结果。图6是GT容量占比为50%和70%相比GT容量占比为30%的气压灵敏度指标的变化量，可以看出随着GT容量占比的增加，扰动1下气压综合灵敏度指标的值不断增大，是因为GT所消耗的天然气增加，天然气系统的运行受到的影响更大。图7是GT1分别接入气节点4和13相比GT1接入气节点3的气压灵敏度指标的变化量，可以看出GT1接入节点4的指标值变化量较小，其中节点4和7有小幅度增加而其他节点有所减小。这是因为节点4和节点3在该天然气系统中的位置相似，GT1接入节点4时节点4的气压受到的影响会变大，因此节点4和同侧的节点7的指标值有所增加。从图6可以看出，当GT1接入节点13时，指标值明显增加，因为节点13是天然气系统的末端负荷节点，电力系统中的扰动通过节点13传递到天然气系统的影响更大。从图5和图6指标变化的具体数值可知，该测试系统中电力系统平衡机组GT1接入位置的变化比GT容量变化对天然气系统的影响更大。

图6 气压综合灵敏度指标在不同GT容量占比的变化量
Fig.6 Variation of composite sensitivity index of pressure at different GT capacities

图7 气压综合灵敏度指标在不同GT1接入位置的变化量
Fig.7 Variation of composite sensitivity index of pressure at different GT1 access positions

由图8 （a）可知，在电压-注入有功功率灵敏度矩阵中，对角线位置的值较大，说明节点电压受节点自身扰动的影响较大。从图8 （b）可知，节点3、4、6、7、8和28的灵敏度值较大，其中节点7和8的负荷值较大，节点2的负荷值也较大，但节点2为发电机节点，所以其相邻的节点3和4有较大的灵敏度值。由图9可以看出，节点3、4、6、7、8和28的电压综合灵敏度指标值较大，说明这些节点的电压受天然气系统中扰动的综合影响较大，为电力系统的薄弱环节。

图8 P2G作用下的电压灵敏度矩阵
Fig.8 Voltage sensitivity matrices under the coupling of P2G

图9 P2G作用下的电压综合灵敏度指标
Fig.9 Composite sensitivity index of voltage under the coupling of P2G

图10为气节点13、14和15增加容量为节点负荷容量30%的储气罐后的电压综合灵敏度指标的变化值，增加储气装置后综合灵敏度指标有所降低，说明电力系统的运行状态有所改善。

图10 增加储气罐后的电压综合灵敏度指标的变化量
Fig.10 Variation of composite sensitivity index of voltage after adding gas storage

图11为在电节点7和8增加容量为节点负荷容量30%的储电装置后的气压综合灵敏度指标的变化值，增加储电装置后每个节点的气压综合灵敏度指标明显降低，因为储电装置可以对P2G电能消耗的变化进行缓冲，从而提升天然气系统的运行安全性。

图11 增加储电装置后的气压综合灵敏度指标的变化量
Fig.11 Variation of composite sensitivity index of pressure after adding electric power storage

4.2 牛顿下山法有效性分析

在4.1节的基础上，研究传统牛顿拉夫逊法和牛顿下山法求解天然气系统潮流的收敛特性。对4.1节的场景1采用两种潮流计算方法计算，结果见图12（a）和（b）；调整天然气系统的气压初始值，使其更接近潮流收敛后的气压值，然后再用两种潮流计算方法进行求解，结果见图12（c）和（d）。

图12 两种方法的不平衡量变化曲线
Fig.12 Evolutions of the mismatch value of the two calculation methods

从图12可知，采用牛顿拉夫逊法求解场景1的原始算例时，不平衡量随着迭代次数的增加而变大，即潮流不收敛，而调整初值后其结果是收敛的。图12（a）的迭代次数为9次，时间约为0.072 s；图12 （c）和 （d）的迭代次数分别为7次和18次，时间分别为0.052 s和0.095 s，显然牛顿下山法的计算时间更短。从图12（d）可以看出，调整初值后牛顿拉夫逊法是震荡收敛的，与牛顿下山法相比，迭代次数更多，计算时间更长。这是因为牛顿下山法要求迭代过程中不平衡量的值比前一次小，所以不平衡量逐渐变小直至收敛。而牛顿拉夫逊法没有该要求，只是根据计算结果更新状态变量，若初值远离收敛值，则不平衡量可能发散。综上，牛顿下山法求解天然气系统潮流的收敛性能更好，可以改善传统牛顿拉夫逊法的初值敏感问题。

5 结论

电-气耦合系统是一种重要的能源系统，本文基于统一潮流模型和灵敏度分析法研究了IEGS中电力系统和天然气系统的相互影响。首先，建立IEGS的统一潮流模型并基于牛顿拉夫逊法实现潮流求解，其中天然气系统的状态变量采用牛顿下山法更新。随后，定义了电压/气压-节点注入功率的综合灵敏度指标，分别用于分析系统在P2G和GT耦合作用下的影响关系。仿真结果表明，所定义的综合灵敏度指标能够同时定量评估系统多种扰动对状态变量的影响，从而辨识系统薄弱环节并给出安全提升措施建议。通过通用的算例对比说明了牛顿拉夫逊法求解天然气系统潮流因初值敏感性可能会不收敛，而牛顿下山法具有良好的收敛性能。

本文重点考虑的是电-气耦合系统的稳态建模，耦合元件建立的是线性模型。而在实际系统运行中耦合元件具有非线性特性，气网的动态特性也是影响系统运行的重要因素之一。未来可重点关注系统中非线性特性和动态特性的精确建模，进而分析耦合系统的相互作用机理。