一种天然气管道多状态故障模型及电-气耦合连锁故障模拟方法

1.本发明涉及综合能源系统安全性及可靠性领域，尤其涉及系统故障在电、气两系统间的传播以及初始故障导致的一系列连锁故障。


背景技术：

2.现如今，新能源发展迅速，世界各国都在大力发展新能源以期转变能源发展方式、减轻环境污染以及优化传统能源的结构，未来新能源将发挥越来越多的作用。社会发展对清洁能源的需求大大增加了天然气系统与电力系统的耦合程度，这就会使天然气系统的故障影响到电力系统，这对综合能源系统的安全运行产生了巨大威胁。
3.y ao r等人(2017)主要研究了电力系统的连锁故障，而忽略了天然气系统带给耦合系统的影响。目前，对于电力系统的故障分析已经较为成熟。在新能源发展背景下，其它系统对于电力系统的影响亟待研究。王伟亮等人(2016)对典型电-气-热综合能源系统划分了层次，建立了系统的稳态模型，分析了小干扰事件对综合能源系统安全性及稳定性的影响。然而，耦合系统的运行是一个复杂的动态过程，使用稳态模型会使结果偏离实际。并且这种方法只考虑了小扰动而忽略了故障对系统的影响，对于系统的安全性分析而言是远不够的。jiang zhai等人(2021)建立了电-气一体化能源系统中耦合装置及天然气系统动态过程模型，用特征法求解双曲偏微分方程，分析了耦合装置调整引起的天然气系统暂态过程。然而这种方法忽略了天然气系统暂态过程对电力系统的影响。陈娟伟等人(2018)建立了双向耦合的电-气综合能源系统，通过解析法分析了天然气供应不足对供电可靠性的影响。这种方法能够揭示当燃气轮机故障供气不足时，对电力系统发电产生的一系列影响。但是，其研究的故障比较单一，且只探究了故障对于供电的影响。没有涉及到天然气系统关键元件故障时对两个系统整体的影响。同样的，包铭磊等人(2019)虽然分析了天然气系统随机故障造成的两个系统的连锁故障过程，但其仅考虑了气源随机故障，故障种类较少。并且随着天然气系统的发展，气源的供应已经成熟稳定，可靠性较高。所以气源的随机故障代表性不强。laube u等人(2016)用ford-fulkerson法量化了天然气系统风险，还分析了天然气系统气源、管道等关键元件随机发生故障对于电力系统的影响。这种方法虽然考虑到了天然气系统关键元件的故障，但是其对于故障的模拟过于简单，选择固定的故障来分析两系统的运行状态会使实验结果与现实存在偏差，不能充分探究故障传播机理及其造成的影响。目前，大部分的研究都是选择的固定的系统故障或是单一故障来分析耦合系统的连锁反应。然而，系统故障是一个十分复杂的动态过程，故障与故障之间会相互影响、相互传播，并且故障状态是会发生转变的。另外，系统元件的故障率不是一成不变的，会随所处环境及外部条件的变化而变化。以上方法都忽略了这些因素。


技术实现要素：

4.本发明主要是克服现有方法中的不足之处，本发明的目的是提供一种考虑天然气
系统管道故障过程的电-气耦合系统连锁反应分析方法。
5.为达到以上技术目的，本发明采用以下技术方案：
6.一种天然气管道多状态故障模型及电-气耦合连锁故障模拟方法，其特征在于，所述分析方法包括以下步骤：
7.s1：在天然气系统实际的可靠性数据统计工作中，通常采用年平均故障率，不能准确反应实际故障情况，这里根据不同故障因素拟合实时故障率。
8.统计某地区近几年来导致天然气管道故障的因素：天然气管道泄漏的因素有老化、极端事件、低温断裂、腐蚀以及其他原因(第三方破坏、违章操作、施工不当等)；天然气管道堵塞的原因主要有积水堵塞以及冰霜堵塞；因为管道脆性断裂发生泄漏和管道冰堵多发生在温度低的春冬季，具有明显的季节性，所以这里采用月故障率。
9.进一步的，计算故障因素权重。
10.对于管道泄漏，有前述5种因素，在泄漏中所占权重为ωj(1≤j≤5)。
11.n＝[n
1 n
2 n
3 n
4 n5]
[0012]
nj为第j种因素引起的管道泄漏事件数量，同理可用kj表示第j种因素引起的管道堵塞时间数量,不同因素在管道堵塞的权重为σj(1≤j≤2)。
[0013]
进一步的，归一化故障因素，得到不同故障因素对管道故障的权重。
[0014][0015][0016]
s2：计算故障因素修正系数：天然气管道的故障率随使用年限的变化符合浴盆曲线，利用改进的weibull分布来拟合。
[0017][0018]
式中：t为管道使用时间，t1为管道使用初期磨合中故障偶发时段，t1～t2表示管道使用中期，管道故障率趋于稳定，t2为管道使用后期故障多发。β可以用极大似然估计法求解，t1、t2由统计数据给出。
[0019]
极端事件修正模型：因为极端事件突发性强，只能根据某地连续n年的给出的极端事件数据，比较本月所处情况与历史情况来进行拟合。修正值如式所示
[0020][0021]
式中，n表示待评估月发生极端事件的次数，ni表示统计年中第m个月发生极端事件的次数。同样的，低温断裂的修正模型以同样的方法得到。
[0022]
由于管道腐蚀以及其他原因导致的管道故障难以统计且难寻其规律，故只考虑其权重。
[0023]
管道水堵修正模型：统计连续n年管道水堵次数，拟合水堵故障率修正值
[0024][0025]
s为待评估月管道水堵次数，sm为当年第m月发生管道水堵的次数。
[0026]
通过连续几年的统计数据发现，管道冰堵通常发生在寒冷的冬季，且管道的故障率基本相同，所以取管道冰堵修正系数e2＝1。
[0027]
s3：从历史统计数据中得到管道泄漏、堵塞连续n年的月平均故障率，管道泄漏主要因素有5种，管道堵塞主要因素有两种，基于s1、s2中计算得到的权重系数以及修正系数得到管道泄漏动态故障率管道堵塞动态故障率φ。
[0028]
管道泄漏动态故障率
[0029]
管道堵塞动态故障率
[0030]
式中：λj为统计数据得到的月平均管道泄漏率，γi为月平均管道堵塞率，1≤i≤12。
[0031]
s4:分别建立管道泄漏4态马尔可夫模型和管道堵塞单独的双状态模型；管道泄漏4态马尔可夫模型中，状态1为管道正常运行态，状态2为管道破裂态，状态3为管道小孔泄漏态，状态4为管道大孔泄漏态，其中每个状态发生的概率为pi(i＝1、2、3、4)，从状态i转移到状态j的状态转移率为(j＝1、2、3、4，且i≠j)，状态i经修复后进入状态j的修复率为μ
ij
；管道堵塞双状态模型中，状态1为管道正常运行态，状态2为管道堵塞态，其中管道堵塞故障概率为f，转移率为φ修复率为ν。
[0032]
元件停运概率可由下式求得：
[0033][0034]
f＝φ/(φ ν)
[0035]
外界环境变化对修复率的影响很小，故认为修复率不变，可从统计数据中得出。
[0036]
s5：根据管道泄漏与堵塞的状态空间图建立下式
[0037][0038][0039]
解得
[0040]
p1 p2 p3 p4＝1
[0041][0042]
[0043][0044]
其中mttr为故障平均修复时间，单位为小时。除了未知外，其他数据均可以从统计数据中计算得出。
[0045]
s6：根据频率持续时间法，列出下式
[0046]
p2μ
21
＝f2[0047][0048][0049]
其中fi为状态i的失效频率。
[0050]
将s5中的式子代入上式，可以得到式子中的值。
[0051]
s7：采用双层蒙特卡罗方法对天然气管道故障过程进行模拟。
[0052]
在第一层蒙特卡罗方法中，确定现在的状态是管道泄漏还是管道堵塞；
[0053]
如果选中状态为管道泄漏，则继续采样两个独立的随机数来确定泄漏孔的位置和大小。设管道长度为l，泄漏点与管道起点的距离为x，统一转化为圆孔计算，孔径为d，管道的直径为d，假设泄漏前管道初始点的压力为p1，末端压力为p4，泄漏点压力为p3，对应位置管道中线上的压力为p2，每个点的温度为ti，泄漏量为q，原来的管道流量为q0，全部转化为圆孔来进行计算。
[0054]
对0到1之间均匀分布的随机数u、u’进行采样
[0055]
x＝lu
[0056]
d＝du’[0057]
如果d/d≤0.2，此时为小孔模型，这种模型下泄漏对管道内的压力与流量的影响极小，故认为泄漏前后参数不变。
[0058]
如果0.2＜d/d＜0.5，此时为大孔模型，采用下式计算泄漏点后的管道流量
[0059][0060][0061]
管道上游流量q
so
＝q 0.5q，管道下游流量q
or
＝q-0.5q
[0062]
a为泄漏孔面积，cd为流量系数，取0.9～0.98，k为绝热指数，一般取1.3，pa为外界环境压力，是已知量；m是天然气的摩尔质量(kg/mol)，g为流过单位面积管道的质量流量(kg/(m2.s))，通用气体常数是r(j/(mol.k))，水力摩擦系数ξ；
[0063]
如果d/d≥0.5，认为管道断裂，管道内流量为0。
[0064]
进一步的，在第二层蒙特卡罗算法中采样当前状态和过渡到下一个状态的持续时间。根据从当前状态到下一个可能状态的转换速率，对所有可能的当前状态持续时间进行采样。选择持续时间最小的状态确定下一个过渡状态，并且将此时间作为当前状态的持续
时间。
[0065][0066]
进一步的，如果在第一层采样中选中的状态为管道堵塞，继续对管道堵塞后实际内径d1随机采样
[0067]
d1＝dln u
[0068]
采用下式计算堵塞后管道内流量
[0069][0070][0071][0072]
a为常数，ρa为正常工况下天然气密度，ta为环境温度，qo为堵塞前管道流量，qd为堵塞后管道内流量。
[0073]
第二层采样过程与管道堵塞情况相同。
[0074]
s8：通过蒙特卡罗模拟天然气管道故障，计算故障后的能流，若是能流越限，天然气系统将会调整气源出力以及切负荷，切负荷目标函数需使系统运行成本及切负荷成本最小。
[0075][0076]
mg作为天然气系统总的节点数量；c
i,g
与c
i,g,cut
分别为天然气系统在节点i处的供气成本以及削负荷成本；q
i,g
与δq
i,g,cut
分别为天然气系统在节点i处的供气量以及天然气削减量。
[0077]
天然气系统故障对电力系统的影响主要是对耦合节点的供气量减少。进行最优能流计算确定δq
i,g,cut
，正常运行时天然气系统供给燃气发电机节点的供气量为q
i,g,l
，故障后该节点的发电功率为
[0078]
p
i,g
＝(q
i,g,l-δq
i,g,cut
)
·
ηg[0079]
p
i,g
为燃气发电机组的有功功率，ηg为天然气转换效率。
[0080]
s9：天然气系统故障将是燃气机组的供气量减少，由此其发电功率相应减少，电力系统发用不平衡，电力系统需调节发电功率及电负荷来实现功率平衡，如果调整发用已经不能平衡功率，则根据电力系统备用容量ri来确定发电量和切负荷量。
[0081][0082]
式中：为电力系统节点m处传统发电机组f的最大功率，表示初始负荷，ns表示电力系统节点个数。
[0083]
若ri＞0，说明孤岛i内的备用容量能够满足负荷需求，这时只需调整发电功率与负荷达到发用平衡。若ri＜0，说明孤岛i内的备用容量无法满足负荷需求，此时就需要进行切负荷处理。负荷削减量通过故障后的供电缺额求得：
[0084]kic
＝min{δpi/∑l
dr,1
}
[0085]kic
为某时刻孤岛i的规定削减比例；δpi为孤岛i的供电缺额；∑l
dr,1
为可削减负荷点的负荷总量，若故障后的供电缺额比可削减负荷之和大，所有的可中断负荷就要以比例1完全削减，其他重要负荷的供电才能得到保证。
[0086]
s10：；用系统负荷损失率αg评估天然气系统故障对两系统的影响。
[0087][0088]
pi为正常运行时系统在节点i处的供电负荷/供气量；p
ks，i
为孤岛s第k次故障后系统在节点i处的供电负荷/供气量；nc为系统负荷节点的数量；ns为故障后系统的孤岛数量；n
sk
是孤岛s在发生k次故障后的节点数量。
[0089]
本发明所提供的一种天然气管道多状态故障模型及电-气耦合连锁故障模拟方法，以历史统计故障数据为基础，将造成管道故障的因素进行分类，分别计算不同故障因素在故障中的权重，并且根据实际情况修正故障系数，拟合了实时动态的故障率，克服了常用固定年平均故障率的不足，更加接近实际情况，使得模拟结果更加准确。结合马尔可夫模型以及蒙特卡罗模拟方法，建立了管道多状态故障模型，对于管道泄漏模型，随机选取泄漏点孔径、在管道中的位置，对于管道堵塞模型，随机选取管道堵塞程度，同时给出了不同故障的数学模型，克服了只采样单一故障的片面性，随机选取故障参数的方式更加接近实际情况，所得结果更加合理可靠。本发明除用于模拟天然气系统故障对电力系统的影响外，还可推广用于模拟天然气系统故障对热力系统的影响、电力系统连锁故障、电热系统连锁故障等等，具有广泛的应用价值。
[0090]
有益效果：
[0091]
与现有技术相比，本发明具有以下有益效果：
[0092]
以历史统计故障数据为基础，将造成管道故障的因素进行分类，分别计算不同故障因素在故障中的权重，并且根据实际情况修正故障系数，拟合了实时动态的故障率，克服了常用固定年平均故障率的不足，更加接近实际情况，使得模拟结果更加准确。结合马尔可夫模型以及蒙特卡罗模拟方法，建立了管道多状态故障模型，对于管道泄漏模型，随机选取泄漏点孔径、在管道中的位置，对于管道堵塞模型，随机选取管道堵塞程度，同时给出了不同故障的数学模型，克服了只采样单一故障的片面性，随机选取故障参数的方式更加接近实际情况，所得结果更加合理可靠。本发明除用于模拟天然气系统故障及电-气耦合系统连锁反应外，还可推广用于模拟天然气系统故障对热力系统的影响、电力系统连锁故障、电热系统连锁故障等等，具有广泛的应用价值。
附图说明
[0093]
图1为电-气耦合测试系统；
[0094]
图2为2016～2021年城市输配管线泄漏事故统计；
[0095]
图3为2016～2021年城市输配管线堵塞事故统计；
[0096]
图4为管道泄漏4态马尔可夫模型；
[0097]
图5为管道堵塞双状态模型；
[0098]
图6为天然气管道故障模拟及电-气耦合连锁故障模拟流程图；
[0099]
图7为不同场景故障下气转电负荷切除比例分布曲线；
[0100]
图8为不同场景故障下电力负荷负荷切除比例分布曲线；
具体实施方式
[0101]
为了使本发明的目的、技术方案及优点更加清楚明白，以下结合实施例，对本发明进行进一步详细说明。应当理解，此处所描述的具体实施例仅仅用以解释本发明，并不用于限定本发明。
[0102]
实施例：
[0103]
一种考虑天然气系统管道故障过程的电-气耦合系统连锁反应模拟方法，具体的方法包括以下步骤：
[0104]
s1：在天然气系统实际的可靠性数据统计工作中，通常采用年平均故障率，不能准确反应实际故障情况，这里根据不同故障因素拟合实时故障率。
[0105]
统计某地区近几年来导致天然气管道故障的因素：天然气管道泄漏的因素有老化、极端事件、低温断裂、腐蚀以及其他原因(第三方破坏、违章操作、施工不当等)；天然气管道堵塞的原因主要有积水堵塞以及冰霜堵塞；因为管道脆性断裂发生泄漏和管道冰堵多发生在温度低的春冬季，具有明显的季节性，所以这里采用月故障率。
[0106]
进一步的，计算故障因素权重。
[0107]
对于管道泄漏，有前述5种因素，在泄漏中所占权重为ωj(1≤j≤5)。
[0108]
n＝[n
1 n
2 n
3 n
4 n5]
[0109]
nj为第j种因素引起的管道泄漏事件数量，同理可用kj表示第j种因素引起的管道堵塞时间数量,不同因素在管道堵塞的权重为σj(1≤j≤2)。
[0110]
进一步的，归一化故障因素，得到不同故障因素对管道故障的权重。
[0111][0112][0113]
s2：计算故障因素修正系数：天然气管道的故障率随使用年限的变化符合浴盆曲线，利用改进的weibull分布来拟合。
[0114][0115]
式中：t为管道使用时间，t1为管道使用初期磨合中故障偶发时段，t1～t2表示管道使用中期，管道故障率趋于稳定，t2为管道使用后期故障多发。β可以用极大似然估计法求解，t1、t2由统计数据给出。
[0116]
极端事件修正模型：因为极端事件突发性强，只能根据某地连续n年的给出的极端事件数据，比较本月所处情况与历史情况来进行拟合。修正值如式所示
[0117][0118]
式中，n表示待评估月发生极端事件的次数，ni表示统计年中第m个月发生极端事件的次数。同样的，低温断裂的修正模型以同样的方法得到。
[0119]
由于管道腐蚀以及其他原因导致的管道故障难以统计且难寻其规律，故只考虑其权重。
[0120]
管道水堵修正模型：统计连续n年管道水堵次数，拟合水堵故障率修正值
[0121][0122]
s为待评估月管道水堵次数，sm为当年第m月发生管道水堵的次数。
[0123]
通过连续几年的统计数据发现，管道冰堵通常发生在寒冷的冬季，且管道的故障率基本相同，所以取管道冰堵修正系数e2＝1。
[0124]
s3：从历史统计数据中得到管道泄漏、堵塞连续n年的月平均故障率，管道泄漏主要因素有5种，管道堵塞主要因素有两种，基于s1、s2中计算得到的权重系数以及修正系数得到管道泄漏动态故障率管道堵塞动态故障率φ。
[0125]
管道泄漏动态故障率
[0126]
管道堵塞动态故障率
[0127]
式中：λj为统计数据得到的月平均管道泄漏率，υi为月平均管道堵塞率，1≤i≤12。
[0128]
s4:分别建立管道泄漏4态马尔可夫模型和管道堵塞单独的双状态模型；管道泄漏4态马尔可夫模型中，状态1为管道正常运行态，状态2为管道破裂态，状态3为管道小孔泄漏态，状态4为管道大孔泄漏态，其中每个状态发生的概率为pi(i＝1、2、3、4)，从状态i转移到状态j的状态转移率为(j＝1、2、3、4，且i≠j)，状态i经修复后进入状态j的修复率为μ
ij
；管道堵塞双状态模型中，状态1为管道正常运行态，状态2为管道堵塞态，其中管道堵塞故障概率为f，转移率为φ修复率为ν。
[0129]
元件停运概率可由下式求得：
[0130][0131]
f＝φ/(φ ν)
[0132]
外界环境变化对修复率的影响很小，故认为修复率不变，可从统计数据中得出。
[0133]
s5：根据管道泄漏与堵塞的状态空间图建立下式
[0134][0135][0136]
解得
[0137]
p1 p2 p3 p4＝1
[0138][0139][0140][0141]
其中mttr为故障平均修复时间，单位为小时。除了未知外，其他数据均可以从统计数据中计算得出。
[0142]
s6：根据频率持续时间法，列出下式
[0143]
p2μ
21
＝f2[0144][0145][0146]
其中fi为状态i的失效频率。
[0147]
将s5中的式子代入上式，可以得到式子中的值。
[0148]
s7：采用双层蒙特卡罗方法对天然气管道故障过程进行模拟。
[0149]
在第一层蒙特卡罗方法中，确定现在的状态是管道泄漏还是管道堵塞；
[0150]
如果选中状态为管道泄漏，则继续采样两个独立的随机数来确定泄漏孔的位置和大小。设管道长度为l，泄漏点与管道起点的距离为x，统一转化为圆孔计算，孔径为d，管道的直径为d，假设泄漏前管道初始点的压力为p1，末端压力为p4，泄漏点压力为p3，对应位置管道中线上的压力为p2，每个点的温度为ti，泄漏量为q，原来的管道流量为q0，全部转化为圆孔来进行计算。
[0151]
对0到1之间均匀分布的随机数u、u’进行采样
[0152]
x＝lu
[0153]
d＝du’[0154]
如果d/d≤0.2，此时为小孔模型，这种模型下泄漏对管道内的压力与流量的影响极小，故认为泄漏前后参数不变。
[0155]
如果0.2＜d/d＜0.5，此时为大孔模型，采用下式计算泄漏点后的管道流量
[0156][0157][0158]
管道上游流量q
so
＝q 0.5q，管道下游流量q
or
＝q-0.5q
[0159]
a为泄漏孔面积，cd为流量系数，取0.9～0.98，k为绝热指数，一般取1.3，pa为外界环境压力，是已知量；m是天然气的摩尔质量(kg/mol)，g为流过单位面积管道的质量流量
(kg/(m2.s))，通用气体常数是r(j/(mol.k))，水力摩擦系数ξ；
[0160]
如果d/d≥0.5，认为管道断裂，管道内流量为0。
[0161]
进一步的，在第二层蒙特卡罗算法中采样当前状态和过渡到下一个状态的持续时间。根据从当前状态到下一个可能状态的转换速率，对所有可能的当前状态持续时间进行采样。选择持续时间最小的状态确定下一个过渡状态，并且将此时间作为当前状态的持续时间。
[0162][0163]
进一步的，如果在第一层采样中选中的状态为管道堵塞，继续对管道堵塞后实际内径d1随机采样
[0164]
d1＝dln u
[0165]
采用下式计算堵塞后管道内流量
[0166][0167][0168][0169]
a为常数，ρa为正常工况下天然气密度，ta为环境温度，qo为堵塞前管道流量，qd为堵塞后管道内流量。
[0170]
第二层采样过程与管道堵塞情况相同。
[0171]
s8：通过蒙特卡罗模拟天然气管道故障，计算故障后的能流，若是能流越限，天然气系统将会调整气源出力以及切负荷，切负荷目标函数需使系统运行成本及切负荷成本最小。
[0172][0173]
mg作为天然气系统总的节点数量；c
i,g
与c
i,g,cut
分别为天然气系统在节点i处的供气成本以及削负荷成本；q
i,g
与δq
i,g,cut
分别为天然气系统在节点i处的供气量以及天然气削减量。
[0174]
天然气系统故障对电力系统的影响主要是对耦合节点的供气量减少。进行最优能流计算确定δq
i,g,cut
，正常运行时天然气系统供给燃气发电机节点的供气量为q
i,g,l
，故障后该节点的发电功率为
[0175]
p
i,g
＝(q
i,g,l-δq
i,g,cut
)
·
ηg[0176]
p
i,g
为燃气发电机组的有功功率，ηg为天然气转换效率。
[0177]
s9：天然气系统故障将是燃气机组的供气量减少，由此其发电功率相应减少，电力系统发用不平衡，电力系统需调节发电功率及电负荷来实现功率平衡，如果调整发用已经不能平衡功率，则根据电力系统备用容量ri来确定发电量和切负荷量。
[0178]
[0179]
式中：为电力系统节点m处传统发电机组f的最大功率，表示初始负荷，ns表示电力系统节点个数。
[0180]
若ri＞0，说明孤岛i内的备用容量能够满足负荷需求，这时只需调整发电功率与负荷达到发用平衡。若ri＜0，说明孤岛i内的备用容量无法满足负荷需求，此时就需要进行切负荷处理。负荷削减量通过故障后的供电缺额求得：
[0181]kic
＝min{δpi/∑l
dr,1
}
[0182]kic
为某时刻孤岛i的规定削减比例；δpi为孤岛i的供电缺额；∑l
dr,1
为可削减负荷点的负荷总量，若故障后的供电缺额比可削减负荷之和大，所有的可中断负荷就要以比例1完全削减，其他重要负荷的供电才能得到保证。
[0183]
s10：；用系统负荷损失率αg评估天然气系统故障对两系统的影响。
[0184][0185]
pi为正常运行时系统在节点i处的供电负荷/供气量；p
ks，i
为孤岛s第k次故障后系统在节点i处的供电负荷/供气量；nc为系统负荷节点的数量；ns为故障后系统的孤岛数量；n
sk
是孤岛s在发生k次故障后的节点数量。
[0186]
实施例1：
[0187]
本实施例所用数据来自公开文献(曹明涛等人，2015)，电-气耦合测试系统中天然气系统选择比利时20节点天然气系统，电力系统选择ieee30节点电力系统，如图1所示。该测试系统有管道19条，气源3个，压缩机3台，g1、g2、g6是燃气发电机组，传统发电机组是其余3个。
[0188]
这里以某市天然气系统在2021年12月发生故障的概率为例，统计该地区2016-2021年输配管线故障发生的次数，如图2、3所示，代入s1中式子得到权重情况如表1、表2所示
[0189]
表1管道泄漏不同故障因素所占权重
[0190][0191][0192]
表2管道堵塞不同故障因素所占权重
[0193]
管道堵塞故障因素权重积水堵塞0.4
冰霜堵塞0.6
[0194]
该地区天然气管道在12月份的月故障率由统计数据给出，管道泄漏故障率为0.276
×
10-4
/km管道堵塞故障率为0.041
×
10-4
/km。分别计算管道泄漏与管道堵塞因素的修正系数。将管道泄漏故障的历史数据代入weibull模型拟合，求得s2中β的值,结果如表3所示。
[0195]
表3weibull分布拟合数据
[0196]
βt1t2t0.670.465
[0197]
将表3的数据代入s2中求c1，可得c1的值为1。
[0198]
通过统计发现该地区12月发生5.0级以上地震1次，该地区连续5年该月总共发生4次自然灾害，代入s2中式子得到极端事件因素修正值c2＝1.25。该年12月管道因为低温发生6次脆性断裂，连续5年本月发生了共26次，计算得c3＝1.15。通过统计数据发现该地区12月没有出现管道水堵的情况，所以e1＝0。
[0199]
将上述计算得到的值代入s3中式子，得到修正的管道泄漏动态故障率管道堵塞动态故障率管道堵塞动态故障率基于修正值的设备实时故障率与年故障率的对比如表4所示。
[0200]
表4年故障率与修正后实时故障率对比
[0201][0202]
可以看出，得到的故障率相差很大，这可能是由于当月低温造成的管道脆性断裂和管道冰堵影响了管道的实时故障率，如果用年平均故障率这种差异性就得不到体现。
[0203]
分别建立管道泄漏故障模型以及管道堵塞故障模型，如图4、图5所示。
[0204]
管道四态马尔可夫模型的输入数据见表5，管道堵塞故障转移率即为上面所求的故障率，修复率为0.573
×
10-4
/月。
[0205]
表5管道四态马尔可夫模型输入数据
[0206][0207]
模型中的三个未知转移率和可以用公式计算。和过渡率分别为1.103
×
10-4
、3.310
×
10-4
和1.204
×
10-4
/月。
[0208]
采用蒙特卡罗法模拟天然气系统随机故障，流程如图6所示，一共设立了以下3种场景，每种场景均仿真10000次。天然气网及电力系统参数如表6、表7、表8所示。
[0209]
表6天然气网参数
[0210]
管道起点终点直径/mm长度/km1128904223890633489026455590.143566590.129647590.11974789055881489059910890201091089025111111890421212128904013131389051414148901015151589025161116395.510.5171717315.526181818315.598191910315.56
[0211]
表7天然气系统初始运行状态
[0212][0213]
表8电力系统初始运行状态
[0214][0215]
场景1：仅考虑管道断裂双状态故障影响，并采用年平均故障率；场景2：采用固定年平均故障率，管道故障不涉及状态转移；场景3：采用实时故障率，并考虑管道多状态故障影响。三种场景下的系统负荷损失率αg如图7、图8所示。
[0216]
可见，当仅考虑管道断裂故障情况时，气转电负荷切除比例最高，同样的电力负荷切除比例也是最高的。管道故障有多种状态，不同故障状态对天然气系统的影响是不同的。如果仅仅把管道断裂作为燃气管道的唯一故障，那么管道故障对天然气系统和电力系统的影响将被高估，使得系统的运维成本增大。当我们采用年平均故障率来对多模式故障进行模拟时，从图中可以发现电、气系统的负荷切除比例较小，使得天然气系统故障对电-气耦合系统的影响被低估。原因可能是12月的气温较低，比起温度较高的夏季，这个月天然气管道发生脆性断裂、冰堵的概率大大增加，而采用年平均故障率无法体现这样的变化。这里我
们拟合了实时故障率，并模拟了天然气管道故障的复杂过程，比起现有的其他方法，结果更加接近实际，更可靠，对电-气系统可靠性评估具有重要意义。
[0217]
本发明以历史统计故障数据为基础，将造成管道故障的因素进行分类，分别计算不同故障因素在故障中的权重，并且根据实际情况修正故障系数，拟合了实时动态的故障率，克服了常用固定年平均故障率的不足，更加接近实际情况，使得模拟结果更加准确。结合马尔可夫模型以及蒙特卡罗模拟方法，建立了管道多状态故障模型，对于管道泄漏模型，随机选取泄漏点孔径、在管道中的位置，对于管道堵塞模型，随机选取管道堵塞程度，同时给出了不同故障的数学模型，克服了只采样单一故障的片面性，随机选取故障参数的方式更加接近实际情况，所得结果更加合理可靠。本发明除用于模拟天然气系统故障及电-气耦合系统连锁反应外，还可推广用于模拟天然气系统故障对热力系统的影响、电力系统连锁故障、电热系统连锁故障等等，具有广泛的应用价值。
[0218]
以上所述，并非对本发明作任何形式上的限制，虽然本发明已通过上述实施例揭示，然而并非用以限定本发明，任何熟悉本专业的技术人员，在不脱离本发明技术方案范围内，当可利用上述揭示的技术内容作出些变动或修饰为等同变化的等效实施例，但凡是未脱离本发明技术方案的内容，依据本发明的技术实质对以上实施例所作的任何简单修改、等同变化与修饰，均仍属于本发明技术方案的范围内。