本发明涉及天然气系统稳态能流计算技术领域，具体是一种考虑多重不确定性的天然气概率-模糊能流分析方法。
背景技术：
能流计算作为求解天然气状态分布的基本计算，是天然气运行和规划的基础。而在实际系统中，往往存在诸多不确定性，如负荷波动、温度变化等。它们的变化通常会对系统运行造成不可忽视的影响。因此，有必要展开天然气系统不确定能流的相关研究。
现有关于不确定性能流研究主要为仅考虑以负荷为代表的系统注入能量不确定性，并评估其对于系统能流分布的影响，未考虑系统的参数不确定性。事实上，管道综合参数与环境温度、气质等多方面因素有关，当其大小随外界环境变化时，会改变天然气的能流分布，甚至影响系统的安全运行，因此也应被视为变量进行分析。
在对上述因素建模时，应充分考虑其分布特性，建立相应的不确定模型。对于负荷这类不确定因素，由于其统计数据多且随机分布明确，具有概率特性，可采用正态分布刻画其不确定性；而对于管道综合参数，由于缺乏充分的历史数据，且其分布范围和特性并不明显，具有显著的模糊特性，因而更适合采用隶属度函数表示。因此，对于负荷和管道综合参数，应分别建立其相应的概率模型和模糊模型。
综上所述，为充分研究天然气系统各类不确定因素的影响，对于负荷和管道综合参数两类不确定模型，应结合其分布特点，建立相应的概率模型和模糊模型，而现有研究没有同时考虑两类不确定因素对天然气系统能流的影响。
技术实现要素：
本发明的目的是解决现有技术中存在的问题。
为实现本发明目的而采用的技术方案是这样的，一种考虑多重不确定性的天然气概率-模糊能流分析方法，主要包括以下步骤：
1)建立天然气系统模型。
进一步，所述天然气系统模型中主要包括一个平衡节点、n个非平衡节点和m个管道。φn表示所述n个非平衡节点所在的集合。φm表示所述m个管道所在的集合。
所述天然气系统模型的输入为网络结构参数。所述网络结构参数主要包括气源参数、输气管道参数、气负荷、平衡节点承受的压力初值和节点负荷相关系数矩阵。
所述天然气系统模型的输出变量主要包括所有非平衡节点承受的压力和管道天然气流量。
2)建立所述非平衡节点负荷的概率模型。所述非平衡节点的负荷服从正态分布。所述非平衡节点负荷的概率函数如下所示：
式中，x为非平衡节点负荷的数值。μx为非平衡节点负荷的均值。σx为非平衡节点负荷的标准差。
3)利用拉丁超立方法对所述非平衡节点负荷的概率模型抽样，从而得到n×nl维负荷概率样本矩阵x'。
式中，矩阵x'的行下标为所述非平衡节点负荷对应的节点编号。矩阵x'的列下标为所述非平衡节点负荷概率的样本序列。
进一步，得到所述负荷概率样本矩阵x'的主要步骤如下：
3.1)定义n个非平衡节点的负荷为n维变量，即：x1，x2，···，xi，···，xn。非平衡节点负荷的相关系数矩阵cx表示如下：
式中，ρxij为第i个非平衡节点的负荷和第j个非平衡节点的负荷之间的相关系数。i、j∈φn。当i≠j时，0≤ρxij≤1。当i＝j时，ρxij＝1。
3.2)定义具有相关性的n维标准正态随机变量：y1，y2，···，yiyn。n维标准正态随机变量对应的变量样本矩阵为矩阵y。矩阵y的维度为n×nl。标准正态随机变量yi的累积概率分布函数表示如下：
式中，yi为任意标准正态随机变量。
矩阵y的相关系数矩阵为矩阵cy。矩阵cy表示如下：
式中，ρyij为第i个标准正态随机变量和第j个标准正态随机变量之间的相关系数。
3.3)利用nataf变换理论，得到相关系数ρyij。
式中，
利用二分法求解公式7，根据相关系数ρxij得到相关系数ρyij，从而得到矩阵cy中每个元素的值。
3.4)利用公式8对矩阵cy进行cholesky分解，得到下三角矩阵l。
cy＝llt。(8)
3.5)定义各变量彼此独立且服从标准正态分布的n维变量：z1,z2,···,zn。任意变量zi的累积概率分布φ(zi)的取值范围为[0,1]。
利用拉丁超立方抽样法，将φ(zi)等分为nl个子区间，即：[0,1/nl]，[1/nl,2/nl]，···，[(nl-1)/nl,1]。
第i个变量的第s个子区间中，将点值(2s-1)/nl作为采样点的累积概率值。第i个变量的第s个子区间中，采样值zis由标准正态分布的累积概率分布函数反函数得出。采样值zis表示如下：
式中，φ-1(·)表示标准正态分布的累积分布函数的反函数。s代表第i个变量的第s个子区间。s＝0,1,···,nl。nl为子区间总数。
3.6)根据步骤3.5得到n×nl维独立的标准正态分布随机样本矩阵z。矩阵z中的元素为zis。
3.7)根据矩阵l和矩阵z，得到矩阵y中所有元素的值。
y＝lz。(10)
式中，l为矩阵l。z为矩阵z。
3.8)将矩阵y中的元素和矩阵x'中的元素按照等概率原则转换，从而得到矩阵x'中元素xis的值。xis表示如下：
式中，fi-1(·)表示变量xi的累积分布函数反函数。ф(·)表示标准正态分布的累积分布函数。yis为矩阵y中的元素值。
4)建立管道综合参数的模糊模型。
设所述天然气系统模型中的管道m的综合参数为km，m∈φm。m为所述天然气系统模型中的管道总数。
管道m的综合参数的基准值km0表示如下：
式中，dm为管道m的直径。lm为管道m的长度。fm为摩擦系数。sg为天然气的相对密度。
5)利用α-截集法对所述管道综合参数的模糊模型进行抽样。
进一步，对所述管道综合参数的模糊模型进行抽样主要步骤如下：
5.1)取隶属度α＝0时，综合管道参数km的取值范围为
5.2)在综合管道参数km的取值范围内均匀抽样nα次，从而得到综合管道参数km的nα个样本。
5.3)重复步骤5.2，对所有管道综合参数进行抽样，从而获得隶属度α＝0时的管道综合参数模糊样本矩阵k。矩阵k表示如下：
式中，矩阵k的行数表示管道编号。矩阵k的列数表示抽样次数。
5.4)改变隶属度α，使得α＝1。重复步骤5.2和步骤5.3，获取隶属度α＝1时，所有管道综合参数的样本矩阵k'。
6)基于所述非平衡节点负荷的概率模型和所述管道综合参数的模糊模型，利用牛顿法计算所述天然气系统中的管道能流和管道非平衡节点承受的压力。
进一步，计算所述天然气系统模型中的管道能流和管道压力的主要步骤如下：
6.1)根据管道综合参数km和输入的平衡节点承受的压力初值，建立相应的管道天然气流量方程。
设管道两端节点分别为d和g(d、g∈φn)。φn表示所述n个非平衡节点所在的集合。将管道综合参数km改写为kdg，则管道天然气流量fm,dg计算公式如下：
式中，pd为管道非平衡节点d承受的压力。pg表示管道非平衡节点g承受的压力。kdg为管道综合参数。sm,dg为管道非平衡节点d和管道非平衡节点g之间的气流方向。若pd≥pg，sm,dg＝1。若pd<pg,则sm，dg＝-1。
所有管道天然气流量fm,dg由列向量f表示。
在所述天然气系统模型中，管道非平衡节点d的天然气流量平衡量表示如下：
式中，
6.2)计算节点天然气流量的不平衡量。
将管道非平衡节点d的天然气流量平衡量fd带入式14，从而获得节点d的天然气流量不平衡量δfd。不平衡量δfd表示如下：
式中，
不平衡量δfd由矩阵δf统一表示。矩阵δf表示如下：
δf＝[δf1,δf2,…,δfn]t。(16)
式中，δf1,δf2,…,δfn分别为第1个非平衡节点、第2个非平衡节点、…第n个非平衡节点的天然气流量不平衡量。n为管道非平衡节点总数。
6.3)建立雅克比矩阵j。
雅克比矩阵j的计算公式为：
式中，j为n×n的方阵。f1、…、fn分别为第1个节点、第2个节点、…第n个节点的天然气流量平衡量。p1、…、pn分别为第1个节点、第2个节点、…第n个节点承受的压力。
6.4)计算管道非平衡节点承受压力的修正量δп。
通过管道非平衡节点天然气流量不平衡量矩阵δf和雅克比矩阵j，计算所有节点承受压力的修正量δп：
δп＝-j-1δf。(18)
式中，j为雅克比矩阵。δf为管道天然气流量不平衡量矩阵。
δп表示如下：
δп＝[δp1,δp2,…,δpn]t。(19)
式中，δp1,δp2,…,δpn分别为第1个非平衡节点、第2个非平衡节点、…第n个非平衡节点承受压力的修正量。n为非平衡节点总数。
迭代更新节点承受的压力，即：
π(k+1)＝π(k)+δπ(k)。(20)
式中，k表示迭代次数，当k＝1时，δп由公式15得到。k≠1时，δп由式18得到。δп为节点承受的压力修正量。п为第k次迭代时，节点承受的压力。
6.5)重复步骤3和步骤4，利用牛顿法迭代计算雅克比矩阵j和管道非平衡节点承受压力的修正量。迭代次数为tmax。
6.6)收敛性判断。
当管道非平衡节点承受压力的修正量δп满足max(|δп|)<ε时，迭代结束。
当max(|δп|)<ε，且迭代次数k≥tmax时，停止迭代，输出“能流不收敛”，结束计算。
当max(|δп|)≥ε，且迭代次数k＜tmax时。令k＝k+1。返回步骤5继续迭代计算。
6.7)计算天然气系统能流和管道非平衡节点承受的压力。主要步骤如下：
6.7.1)统计隶属度为0下的输出变量范围。
依次将管道综合参数模糊样本矩阵k的第一列、第二列、…、第nα列作为管道综合参数样本，并代入管道天然气流量计算公式13中，记录计算结果。
比较所有管道综合参数样本下的输出结果，并保存所有管道非平衡节点承受压力的最大值п0,max、管道非平衡节点承受压力最小值п0,min、管道天然气流量的最大值f0,max和管道天然气流量最小值f0,min。
6.7.2)统计隶属度为1下的输出变量范围。
依次将管道综合参数模糊样本矩阵k'的第一列、第二列、…、第nα列作为管道综合参数样本，并代入管道天然气流量计算公式13中，记录计算结果。
比较所有管道综合参数样本下的输出结果，并保存所有管道非平衡节点承受压力的最大值п1,max、管道非平衡节点承受压力最小值п1,min、管道天然气流量的最大值f1,max和管道天然气流量最小值f1,min。
6.7.3)统计所述负荷概率样本下输出变量的可能性分布范围。
分别选取所述负荷概率样本矩阵x'的第1列、…、第nl列作为负荷概率样本，并分别计算在nl种负荷概率样本下的管道非平衡节点承受压力可能性分布
7)计算得到输出变量的可能性测度pos和必然性测度nec，从而得到输出变量的信任度函数bel和似真度函数pl。
进一步，计算得到输出变量的可能性测度pos、必然性测度nec、输出变量的信任度函数bel和似真度函数pl的主要步骤如下：
7.1)当管道非平衡节点承受的压力为pd时，得到nl个可能性分布
区间长度为dpd为：
dpd＝pd,max-pd,min。(21)
式中，pi,max为管道非平衡节点承受压力pd的最大值。pd,min为管道非平衡节点承受压力pd的最小值。
将分布区间[pd,min,pd,max]平分为c份，并建立c个子集，则第c个子集ec为：
式中，pd为管道非平衡节点承受压力。pd,max为管道非平衡节点承受压力pd的最大值。pd,min为管道非平衡节点承受压力pd的最小值。c为任意子集。c为子集总数。
则子集ec在pd的第l个可能性分布
式中，pd为管道非平衡节点承受压力。
则子集ec在pd的第l个可能性分布
式中，inf(·)表示取c个集合中的最小值。pd为管道非平衡节点承受压力。
7.2)子集ec的似真度pl表示如下：
式中，ec为任意子集。posl,pd为可能性测度。l既可表示可能性分布编号，也代表了对应的负荷概率样本编号。nl为负荷概率样本总数。dl为负荷第l个概率样本的采样概率。dl＝1/nl。ec为管道非平衡节点承受压力的可能性分布的子集。
子集ec的信任度bel表示如下：
式中，ec为任意子集。necl,pd为必然性测度。dl为负荷第l个概率样本的采样概率。dl＝1/nl。ec为管道非平衡节点承受压力的可能性分布的子集。
本发明的技术效果是毋庸置疑的。本发明结合概率和模糊两类不确定因素的影响，可以有效、准确的计算出天然气系统输出变量的变化范围。即，在天然气系统的不确定能流分析中，本发明同时考虑了概率不确定变量和模糊不确定变量，可以有效计及两类不确定性对能流的影响。
附图说明
图1为管道综合参数k的隶属度函数；
图2为能流方法的流程框图；
图3为基于nataf变换和拉丁超立方抽样建立负荷概率样本流程图；
图4为基于牛顿法的天然气系统能力计算流程框图；
图5为网络结构图；
图6为节点1压力仿真结果图；
图7为管道7流量仿真结果图。
具体实施方式
下面结合实施例对本发明作进一步说明，但不应该理解为本发明上述主题范围仅限于下述实施例。在不脱离本发明上述技术思想的情况下，根据本领域普通技术知识和惯用手段，做出各种替换和变更，均应包括在本发明的保护范围内。
实施例1：
参见图1至图4，一种考虑多重不确定性的天然气概率-模糊能流分析方法，主要包括以下步骤：
1)建立天然气系统模型。
进一步，所述天然气系统模型中主要包括一个平衡节点、n个非平衡节点和m个管道。φn表示所述n个非平衡节点所在的集合。φm表示所述m个管道所在的集合。
所述天然气系统模型的输入为网络结构参数。所述网络结构参数主要包括气源参数、输气管道参数、气负荷、管道内部平衡节点承受的压力初值和节点负荷相关系数矩阵。
所述天然气系统模型的输出变量主要包括所有管道内部非平衡节点承受的压力和管道天然气流量。
2)建立所述非平衡节点负荷的概率模型。所述非平衡节点的负荷服从正态分布。所述非平衡节点负荷的概率函数如下所示：
式中，x为非平衡节点负荷的数值。μx为非平衡节点负荷的均
值。σx为非平衡节点负荷的标准差。
3)利用拉丁超立方法对所述非平衡节点负荷的概率模型抽样，从而得到n×nl维负荷概率样本矩阵x'。
式中，矩阵x'的行下标为所述非平衡节点负荷对应的节点编号。矩阵x'的列下标为所述非平衡节点负荷概率的样本序列。
进一步，得到所述负荷概率样本矩阵x'的主要步骤如下：
3.1)定义n个非平衡节点的负荷为n维变量，即：x1，x2，···，xi，···，xn。非平衡节点负荷的相关系数矩阵cx表示如下：
式中，ρxij为第i个非平衡节点的负荷和第j个非平衡节点的负荷之间的相关系数。i、j∈φn。当i≠j时，0≤ρxij≤1。当i＝j时，ρxij＝1。
3.2)定义具有相关性的n维标准正态随机变量：y1，y2，···，yiyn。n维标准正态随机变量对应的变量样本矩阵为矩阵y。矩阵y的维度为n×nl。标准正态随机变量yi的累积概率分布函数表示如下：
式中，yi为任意标准正态随机变量。
矩阵y的相关系数矩阵为矩阵cy。矩阵cy表示如下：
式中，ρyij为第i个标准正态随机变量和第j个标准正态随机变量之间的相关系数。
3.3)利用nataf变换理论，得到相关系数ρyij。
式中，
利用二分法求解公式7，根据相关系数ρxij得到相关系数ρyij，从而得到矩阵cy中每个元素的值。
3.4)利用公式8对矩阵cy进行cholesky分解，得到下三角矩阵l。
cy＝llt。(8)
进一步，cholesky分解是把一个对称正定的矩阵表示成一个下三角矩阵l和其转置的乘积的分解。它要求矩阵的所有特征值必须大于零，故分解的下三角的对角元也是大于零的。cholesky分解法又称平方根法，是当a为实对称正定矩阵时，lu三角分解法的变形。
3.5)定义各变量彼此独立且服从标准正态分布的n维变量：z1,z2,···,zn。任意变量zi的累积概率分布φ(zi)的取值范围为[0,1]。
利用拉丁超立方抽样法，将φ(zi)等分为nl个子区间，即：[0,1/nl]，[1/nl,2/nl]，···，[(nl-1)/nl,1]。
进一步，设置拉丁超立方抽样次数nl＝2000。
拉丁超立方体抽样的步骤是：将每一维分成互不重迭的m个区间，使得每个区间有相同的概率(通常考虑一个均匀分布，这样区间的长度相同)。在每一维里的每一个区间中随机的抽取一个点。再从每一维里随机抽出(2)中选取的点，将它们组成向量。
第i个变量的第s个子区间中，将点值(2s-1)/nl作为采样点的累积概率值。第i个变量的第s个子区间中，采样值zis由标准正态分布的累积概率分布函数反函数得出。采样值zis表示如下：
式中，φ-1(·)表示标准正态分布的累积分布函数的反函数。s代表第i个变量的第s个子区间。s＝0,1,···,nl。nl为子区间总数。
3.6)根据步骤3.5得到n×nl维独立的标准正态分布随机样本矩阵z。矩阵z中的元素为zis。
3.7)根据矩阵l和矩阵z，得到矩阵y中所有元素的值。
y＝lz。(10)
式中，l为矩阵l。z为矩阵z。
3.8)将矩阵y中的元素和矩阵x'中的元素按照等概率原则转换，从而得到矩阵x'中元素xis的值。xis表示如下：
xis＝fi-1(φ(yis))。(11)
式中，fi-1(·)表示变量xi的累积分布函数反函数。ф(·)表示标准正态分布的累积分布函数。yis为矩阵y中的元素值。
4)建立管道综合参数的模糊模型。
设所述天然气系统模型中的管道m的综合参数为km，m∈φm。m为所述天然气系统模型中的管道总数。
管道m的综合参数的基准值
式中，dm为管道m的直径。lm为管道m的长度。fm为摩擦系数。sg为天然气的相对密度。
天然气是一种无毒无色无味的气体，其主要成份是甲烷，在能流计算中通常取为0.6094。
5)利用α-截集法对所述管道综合参数的模糊模型进行抽样。
进一步，对所述管道综合参数的模糊模型进行抽样主要步骤如下：
5.1)取隶属度α＝0时，综合管道参数km的取值范围为
5.2)在综合管道参数km的取值范围内均匀抽样nα次，从而得到综合管道参数km的nα个样本。
5.3)重复步骤5.2，对所有管道综合参数进行抽样，从而获得隶属度α＝0时的管道综合参数模糊样本矩阵k。矩阵k表示如下：
式中，矩阵k的行数表示管道编号。矩阵k的列数表示抽样次数。
5.4)改变隶属度α，使得α＝1。重复步骤5.2和步骤5.3，获取隶属度α＝1时，所有管道综合参数的样本矩阵k'。
6)基于所述非平衡节点负荷的概率模型和所述管道综合参数的模糊模型，利用牛顿法计算所述天然气系统中的管道能流和管道非平衡节点承受压力。牛顿法收敛精度ε＝10-6。牛顿法最大迭代次数tmax＝50。可能性分布区间数为c＝1000。管道能流即管道中天然气流量。
进一步，计算所述天然气系统模型中的管道能流和管道压力的主要步骤如下：
6.1)根据管道综合参数km和输入的平衡节点承受压力初值，建立相应的管道天然气流量方程。
设管道两端节点分别为d和g(d、g∈φn)。φn表示所述n个非平衡节点所在的集合。将管道综合参数km改写为kdg，则管道天然气流量fm,dg计算公式如下：
式中，pd为管道非平衡节点d承受的压力。pg表示管道非平衡节点g承受的压力。kdg为管道综合参数。sm,dg为管道非平衡节点d和管道非平衡节点g之间的气流方向。若pd≥pg，sm,dg＝1。若pd<pg,则sm，dg＝-1。
所有管道天然气流量fm,dg由列向量f表示。
在所述天然气系统模型中，管道非平衡节点d的天然气流量平衡量表示如下：
式中，
6.2)计算节点天然气流量的不平衡量。
将管道非平衡节点d的天然气流量平衡量fd带入式14，从而获得节点d的天然气流量不平衡量δfd。不平衡量δfd表示如下：
式中，
不平衡量δfd由矩阵δf统一表示。矩阵δf表示如下：
δf＝[δf1,δf2,…,δfn]t。(16)
式中，δf1,δf2,…,δfn分别为第1个非平衡节点、第2个非平衡节点、…第n个非平衡节点的天然气流量不平衡量。n为管道非平衡节点总数。
6.3)建立雅克比矩阵j。
雅克比矩阵j的计算公式为：
式中，j为n×n的方阵。f1、…、fn分别为第1个节点、第2个节点、…第n个节点的天然气流量平衡量。p1、…、pn分别为第1个节点、第2个节点、…第n个节点承受的压力。
在向量微积分中，雅可比矩阵是一阶偏导数以一定方式排列成的矩阵，其行列式称为雅可比行列式。雅可比矩阵的重要性在于它体现了一个可微方程与给出点的最优线性逼近。
6.4)计算管道非平衡节点承受压力的修正量δп。
通过管道非平衡节点天然气流量不平衡量矩阵δf和雅克比矩阵j，计算所有节点承受压力的修正量δп：
δп＝-j-1δf。(18)
式中，j为雅克比矩阵。δf为管道天然气流量不平衡量矩阵。
δп表示如下：
δп＝[δp1,δp2,…,δpn]t。(19)
式中，δp1,δp2,…,δpn分别为第1个非平衡节点、第2个非平衡节点、…第n个非平衡节点承受压力的修正量。n为非平衡节点总数。
迭代更新节点承受的压力，即：
π(k+1)＝π(k)+δπ(k)。(20)
式中，k表示迭代次数，当k＝1时，δп由公式15得到。k≠1时，δп由式18得到。δп为节点承受的压力修正量。п为第k次迭代时，节点承受的压力。
6.5)重复步骤3和步骤4，利用牛顿法迭代计算雅克比矩阵j和管道非平衡节点承受压力的修正量。迭代次数为tmax。
6.6)收敛性判断。
当管道非平衡节点承受压力的修正量δп满足max(|δп|)<ε时，迭代结束。
当max(|δп|)<ε，且迭代次数k≥tmax时，停止迭代，输出“能流不收敛”，结束计算。
当max(|δп|)≥ε，且迭代次数k＜tmax时。令k＝k+1。返回步骤5继续迭代计算。
6.7)计算天然气系统能流和管道非平衡节点承受压力。主要步骤如下：
6.7.1)统计隶属度为0下的输出变量范围。
依次将管道综合参数模糊样本矩阵k的第一列、第二列、…、第nα列作为管道综合参数样本，并代入管道天然气流量计算公式13中，记录计算结果。
比较所有管道综合参数样本下的输出结果，并保存所有管道非平衡节点承受压力的最大值п0,max、管道非平衡节点承受压力最小值п0,min、管道天然气流量的最大值f0,max和管道天然气流量最小值f0,min。
6.7.2)统计隶属度为1下的输出变量范围。
依次将管道综合参数模糊样本矩阵k'的第一列、第二列、…、第nα列作为管道综合参数样本，并代入管道天然气流量计算公式13中，记录计算结果。
比较所有管道综合参数样本下的输出结果，并保存所有管道非平衡节点承受压力的最大值п1,max、管道非平衡节点承受压力最小值п1,min、管道天然气流量的最大值f1,max和管道天然气流量最小值f1,min。
6.7.3)统计所述负荷概率样本下输出变量的可能性分布范围。
分别选取所述负荷概率样本矩阵x'的第1列、…、第nl列作为负荷概率样本，并分别计算在nl种负荷概率样本下的管道非平衡节点承受压力可能性分布
7)计算得到输出变量的可能性测度pos和必然性测度nec，从而得到输出变量的信任度函数bel和似真度函数pl。
进一步，计算得到输出变量的可能性测度pos、必然性测度nec、输出变量的信任度函数bel和似真度函数pl的主要步骤如下：
7.1)当管道非平衡节点承受的压力为pd时，得到nl个可能性分布
区间长度为dpd为：
dpd＝pd,max-pd,min。(21)
式中，pi,max为管道非平衡节点承受的压力pd的最大值。pd,min为管道非平衡节点承受压力pd的最小值。
将分布区间[pd,min,pd,max]平分为c份，并建立c个子集，则第c个子集ec为：
式中，pd为管道非平衡节点承受压力。pd,max为管道非平衡节点承受压力pd的最大值。pd,min为管道非平衡节点承受压力pd的最小值。c为任意子集。c为子集总数。
则子集ec在pd的第l个可能性分布πlpd下的可能性测度posl,pd表示如下：
式中，pd为管道非平衡节点承受压力。
则子集ec在pd的第l个可能性分布πlpd下的必然性测度necl,pd表示如下：
式中，inf(·)表示取c个集合中的最小值。pd为管道非平衡节点承受压力。
7.2)子集ec的似真度pl表示如下：
式中，ec为任意子集。posl,pd为可能性测度。l既可表示可能性分布编号，也代表了对应的负荷概率样本编号。nl为负荷概率样本总数。dl为负荷第l个概率样本的采样概率。dl＝1/nl。ec为管道非平衡节点承受压力的可能性分布的子集。
子集ec的信任度bel表示如下：
式中，ec为任意子集。necl,pd为必然性测度。dl为负荷第l个概率样本的采样概率。dl＝1/nl。ec为管道非平衡节点承受压力的可能性分布的子集。
实施例2：
参见图5至图7，利用一种考虑多重不确定性的天然气概率-模糊能流分析方法分析具有11个节点的天然气系统的能流的实验，主要包括以下步骤：
1)建立天然气系统模型。所述天然气系统具有1个平衡节点和10个非平衡节点。天热气系统管道数m＝14。
设置拉丁超立方抽样次数nl＝2000，α截集法每个截集抽样次数nα＝2000，牛顿法收敛精度ε＝10-6，牛顿法最大迭代次数tmax＝50，可能性分布区间数为c＝1000。
网络数据采用2016年第164卷《appliedenergy》中“steadystateanalysisofgasnetworkswithdistributedinjectionofalternativegas”一文中的11节点天然气系统的数据。修改节点编号，平衡节点编号为0，原来节点2～10编号改为1～10，管道编号保持不变，仍为1～14，但两端节点编号发生变化。
2)建立所述非平衡节点负荷的概率模型。所述非平衡节点的负荷服从正态分布。对于每个非平衡节点，其负荷服从正态分布，以节点2所接负荷为例，其负荷基准值为219m3/h，标准差取其大小的5％，则其累积概率分布函数为：
式中，x为非平衡节点负荷的数值。μx为非平衡节点负荷的均值。σx为非平衡节点负荷的标准差。
3)利用拉丁超立方法对所述非平衡节点负荷的概率模型抽样，从而得到n×nl维负荷概率样本矩阵x‘。
式中，矩阵x'的行下标为所述非平衡节点负荷对应的节点编号。矩阵x'的列下标为所述非平衡节点负荷概率的样本序列。
3.1)令节点负荷变量为x1,x2,···,x10，且变量相关系数矩阵cx已知。矩阵cx如下所示：
3.2)定义有相关性的n维标准正态随机变量样本矩阵为y，其相关系数矩阵为cy，由于xi服从正态分布，得到矩阵cy如下所示：
3.3)对矩阵cy进行cholesky分解，得到下三角矩阵l：
3.4)定义有n维变量z1,z2,···,zn，各变量独立且服从标准正态分布，应用拉丁超立方抽样原理，产生n×nl样本矩阵z。通过等式y＝lz可求取矩阵y。将矩阵y与矩阵x'中的各元素按照等概率原则转换即可得到矩阵x'中元素值，以第一列样本为例，给出x1'中各节点负荷样本值：
4)建立管道综合参数的模糊模型。
设所述天然气系统模型中的管道m的综合参数为km，m∈φm。m为所述天然气系统模型中的管道总数。
建立如图1所示的隶属度函数，其中基准值由已给数据计算而来。以管道1(连接节点1和节点2)为例，其直径d1＝160mm，长度l1＝50m，摩擦系数f1＝0.0056，sg＝0.6096kg/m3，则根据式33，计算出管道1的综合管道参数基准值为：
式中，dm为管道m的直径。lm为管道m的长度。fm为摩擦系数。sg为天然气的相对密度。
5)利用α-截集法对所述管道综合参数的模糊模型进行抽样。
对图1所示的隶属度函数，取隶属度α＝0，以管道1为例，则k1的取值范围为
6)基于所述非平衡节点负荷的概率模型和所述管道综合参数的模糊模型，利用牛顿法计算所述天然气系统中的管道能流和管道非平衡节点承受压力。
根据步骤(3)和(5)中抽取的负荷和管道综合参数的样本，将其作为已知量代入下一步的能流计算。
6.1)根据得到的管道综合参数及已有管道非平衡节点承受压力初值，建立相应的管道天然气流量方程，以管道1为例，算得管道天然气流量f12为：
6.2)计算不平衡量
根据节点承受压力、管道综合参数的实际取值，获得节点流量不平衡量：
6.3)雅克比矩阵
雅克比矩阵j为：
6.4)计算修正量
根据上述步骤获得的节点注入流量不平衡量△f和雅克比矩阵j，以第一次迭代后的结果为例，计算出各节点承受压力的修正量：
根据所得结果更新节点承受压力，则有：
6.5)收敛性判断
当不平衡量△п满足max(|△п|)＜ε，则计算结束。
当max(|△п|)<ε并且迭代次数k≥tmax，则停止迭代，输出“能流不收敛”，结束计算。
当max(|△п|)≥ε并且迭代次数k＜tmax，执行k＝k+1，并则返回步骤6.3继续迭代计算。
根据前面步骤，计算可知当迭代8次后，满足收敛条件，此时max(|△п|)＝3.3429e-8＜ε。
6.6)计算天然气系统能流分布
根据已知节点承受压力，获得系统能流分布。记入平衡节点承受压力p0＝75mbar，采用列向量п'表示节点0～10压力为：
管道天然气流量f为：
6.7)统计输出变量范围
6.7.1)统计隶属度为0下的输出变量范围。
依次将管道综合参数模糊样本矩阵k的第一列、第二列、…、第nα列作为管道综合参数样本，并代入管道天然气流量计算公式13中，记录计算结果。
比较所有管道综合参数样本下的输出结果，并保存所有管道非平衡节点承受压力的最大值п0,max、管道非平衡节点承受压力最小值п0,min、管道天然气流量的最大值f0,max和管道天然气流量最小值f0,min。
比较所有样本下的输出结果并保存各节点承受压力和管道天然气流量的最大值和最小值如下：
6.7.2)统计隶属度为1下的输出变量范围。
依次将管道综合参数模糊样本矩阵k'的第一列、第二列、…、第nα列作为管道综合参数样本，并代入管道天然气流量计算公式13中，记录计算结果。
比较所有管道综合参数样本下的输出结果，并保存所有管道非平衡节点承受压力的最大值п1,max、管道非平衡节点承受压力最小值п1,min、管道天然气流量的最大值f1,max和管道天然气流量最小值f1,min。
绘制出各输出变量的可能性分布，以节点1和管道2为例，其可能性分布分别如表1和表2所示。
表1节点1压力可能性分布
表2管道2流量可能性分布
6.7.3)统计所述负荷概率样本下输出变量的可能性分布范围。
分别选取所述负荷概率样本矩阵x'的第1列、…、第nl列作为负荷概率样本，并分别计算在nl种负荷概率样本下的管道非平衡节点承受压力可能性分布
7)计算得到输出变量的可能性测度pos和必然性测度nec，从而得到输出变量的信任度函数bel和似真度函数pl。
以节点1压力为例，取其在所有可能性分布中的最值，形成分布区间为[58.7372,69.9069]，其区间长度为dp1为：
dp1＝58.7428。(46)
当得到其nl个可能性分布
将分布区间平分为c份，并建立c个子集，则c个子集ec为：
根据式23和式24求得各子集ec在l个可能性分布
计算出各子集ec的似真度pl(ec)和信任度bel(ec)：
对所有节点承受压力和管道天然气流量进行上述操作，即可获得待求输出变量的似真度函数和信任度函数，并绘制相应似真度和信任度曲线。