基于希尔伯特－黄变换的水平管气液两相流流型识别方法

1、一种基于希尔伯特-黄变换的水平管气液两相流流型识别方法，其特征在 于：它是在一个由文丘利管、差压变送器与计算机构成的流型信号处理系统中， 通过采集文丘利管上的差压信号，引入希尔伯特-黄变换方法计算差压信号不同 固有模态所占的能量比得到识别流型的流态图，再通过与残差均值的结合得到 流型识别规则。
2、根据权利要求1所述的一种基于希尔伯特-黄变换的水平管气液两相流流 型识别方法，它的识别步骤为：
1).选取与不同流型相对应的水平管道气水、气油两相差压流数据作为样本 集，对差压信号进行固有模态分解，得到描述信号的不同分量，然后通过定义 的能量比公式对不同的模态进行能量比与平均瞬时频率计算，得到不同频率带下 能量比与流型的关系，通过对3Hz与30Hz频率范围的能量比的分析得到气液两 相流流型识别的流态图；
2).根据样本数据计算得到的流态图识别出泡状流或单相水、弹状流、混状 流，通过残差均值的计算把单相水流从小的泡状流中区别出来，利用3Hz与30Hz 附近的能量比与残差均值的结合总结出流型识别规则，识别出单相水、泡状流、 弹状流、混状流；
3).上面得到的流态图与规则作为知识存储于计算机中，对现场采集的数据 通过计算得到能量比与残差得均值，然后利用流态图或规则进行判断，得到两 相流的流型。
3、根据权利要求2所述的一种基于希尔伯特-黄变换的水平管气液两相流流型 识别方法，其特征在于：所说通过定义的能量比公式对不同的模态进行能量比 与平均瞬时频率计算，得到不同频率带下能量比与流型的关系：是通过计算得到 在3 Hz附近混状流所占的能量比为75％～98％，弹状流为40％～90％，泡状流为 45％以下，在30Hz附近混状流占的能量比为10％以下，弹状流为8％～55％，泡 状流为40％～90％，通过计算这两个频率范围内的能量比，然后判断能量比的范 围，识别出混状流、弹状流、泡状流或单相水。
4、根据权利要求2所述的一种基于希尔伯特-黄变换的水平管气液两相流流 型识别方法，其特征在于：所说通过对3Hz与30Hz频率范围的能量比的分析得 到气液两相流流型识别的流态图：是把50mm与40mm管径气液两相流在3Hz附 近的能量比E1与30Hz附近的能量比E3画成流态识别图，通过流态识别图识别 出不同管径下泡状流、弹状流与混状流，阴影部分分别为泡-弹与弹-混状流型 的过渡区间。
5、根据权利要求2所述的一种基于希尔伯特-黄变换的水平管气液两相流流 型识别方法，其特征在于：所说通过残差均值的计算把单相水流从小的泡状流 中区别出来，利用3Hz与30Hz附近的能量比与残差均值的结合总结出流型识别 规则，识别出单相水、泡状流、弹状流、混状流：是通过对残差均值的计算， 得到单相流时残差的均值在0.94以上，泡状流在0.65～0.91之间，弹状流在 0.40～0.65，而混状流分布在0.4以下的区域，利用残差均值与流态识别图相结合 识别出单相水、泡状流、弹状流、混状流流型，另外如果同时满足两个规则， 则可以判别为该流型为过渡流型。

技术领域\n本发明涉及气液两相流测量技术，尤其涉及一种基于希尔伯特-黄变换 (Hirbert Huang Trasform)的水平管气液两相流流型识别方法。\n背景技术\n两相流动广泛存在于自然界和许多工业过程中，如电力、能源、石油、核 工业等。由于两相流各相间存在着界面反应和相对速度，相界面在时间和空间 上都是随机变化的，致使其流动结构极其复杂。因此，两相流参数检测难度较 大。描述两相流的参数主要有：流型、空隙率、分相流量等。流型的研究是两 相流的机理及其运动规律的重要组成部分，从工程角度看，流型直接影响到两 相流系统的管程压降，决定了动力工具的设计依据，流型又影响两相流系统的 稳定性与安全性。\n由于应用流态图识别流型的局限性，目前人们开始研究利用仪器设备直接 识别管道内气液两相流的流型的方法。主要要两类，一类是直接根据两相流的 流动形式确定流型，如观测法、高速摄影法等，该方法具有一定的主观性，对 于不同的观测者可能得到不同的结果；另一类是通过反映两相流动特性的波动 信号进行处理，提取特征量，进行流型识别。该方法通过提取信号的特征量， 实现对流型的定量识别。\n采用信号处理方法实现流型的识别是一个研究的热点，目前的流型识别大 多数人采用阻抗、电容测量的空隙率作为分析信号，然后提取特征值通过模糊 算法与神经网络技术进行流型识别。这些方法的最大缺点就是成本高，传感器 设计复杂。\n另外，也有很多研究者采用管道的压力或者差压信号，通过计算功率谱密 度或概率密度函数的方法进行流型识别，功率谱密度与概率分布函数要求信号 是平稳的，而二相流动信号具有强烈的非线性特性，难以满足平稳性要求。因 此，难以得到通用的流型识别方法。\n为了突破传统信号分析方法的限制，Huang等1998年提出了一个新的数据 处理方法：希尔伯特-黄变换方法，该方法将原始信号用固有模态展开，展开后 的每个模态有实际的物理意义，而且展开后使瞬时频率具有了实际意义。该方 法可以处理非线性或非平稳的信号，是一种真正意义上的时频分析方法。该方 法具有完整性，近乎正交性，更重要的是具有自适应性。目前希尔伯特-黄变换 方法已经在海洋、机械等领域获得成功的应用，在解决非线性问题时体现出很 好的性能。\n为了准确提取两相流流型信息，利用信号处理方法提取特征量，信号的选 取也十分重要的。管道直接取差压的缺点是信号的变化幅度小，容易丢失有用 的信息，而文丘利管的差压信号具有提取容易、信号幅度大、对流型影响小适 于用作流型识别的信号。\n两相流动是一个复杂的多变量随机过程，随着随机过程理论和信号处理技术 的不断完善和发展，应用数理统计、参数估计、模式识别等理论和技术，进行 多相流参数的软测量是一个很有前途的研究方向。\n采用多传感器组合，进行多参数组合测量确定流量、流型与干度也获得了很 多的应用。采用文丘利管与测量空隙率的仪表的结合实现两相分相流量的测量 也取得了很好的结果。另外，利用文丘利管上的差压信号进行空隙率(干度) 的估计也取得了一些进展。采用文丘利管计算两相流流量的经典公式与利用文 丘利管的差压信号识别流型的结合也是进一步研究的方向。\n发明内容\n本发明的目的在于提出一种基于希尔伯特-黄变换的水平管气液两相流流 型识别方法。\n它是在一个由文丘利管、差压变送器与计算机构成的流型信号处理系统中， 通过采集文丘利管上的差压信号，引入希尔伯特-黄变换方法计算差压信号不同 固有模态所占的能量比得到识别流型的流态图，再通过与残差均值的结合得到 流型识别规则。\n本发明针对目前两相流研究的现状与存在的问题，提出了利用传统的文丘利 管，通过提取差压信号，使用希尔伯特-黄变换方法实现气液两相流的流型识别 方法。该方法克服了传统的流型识别方法通用性差、受实验条件应响大、难于 实际应用的不足。\n附图说明\n图1是水平管气液两相流流型识别系统示意图；\n图2是利用希尔伯特-黄变换实现流型识别的过程示意图；\n图3是单相水信号及其固有模态分解图；\n图4是泡状流信号及其固有模态分解图；\n图5是弹状流信号及其固有模态分解图；\n图6是混状流信号及其固有模态分解图；\n图7是50mm管道的瞬时频率和能量比与流型的关系图；\n图8是40mm管道的瞬时频率和能量比与流型的关系图；\n图9是50mm管道的残差均值和液相流量与流型的关系图；\n图10是40mm管道的残差均值和液相流量与流型的关系图；\n图11是水平管气液两相流不同能量比识别流型的流态图。\n具体实现方式\n基于希尔伯特-黄变换的水平管气液两相流流型识别方法的步骤为：\n1).选取与不同流型相对应的水平管道气液(气水、气油)两相差压流数据作为 样本集，对差压信号进行固有模态分解，得到描述信号的不同分量，然后通过 定义的能量比公式对不同的模态进行能量比与平均瞬时频率计算，得到不同频率 带下能量比与流型的关系，通过对3Hz与30Hz频率范围的能量比的分析得到气 液两相流流型识别的流态图；\n2).根据样本数据计算得到的流态图识别出泡状流(单相水)、弹状流、混状 流，通过残差均值的计算把单相水流从小的泡状流中区别出来，利用3Hz与30Hz 附近的能量比与残差均值的结合总结出流型识别规则，识别出单相水、泡状流、 弹状流、混状流；\n图1为气液两相流流型识别系统示意图，图2为具体实现的过程。通过测量 反映流型信息的差压信号，然后通过计算机对信号进行希尔伯特-黄变换，得到 不同模态的能量比与残差的均值，利用得到的流态图与规则进行不同流型的识 别。希尔伯特-黄变换原理与具体算法如下所述：\n希尔伯特-黄变换假设任一信号都是由许多固有模态函数组成，如果固有模 态之间相互重叠，便形成复合信号。在希尔伯特-黄变换中描述信号的基本量是 瞬时频率，这是与其它信号处理方法的本质区别，而瞬时频率对每一个固有模 态函数都有实际意义。瞬时频率可以通过Hilbert变换求得，即先将对信号s(t) 作Hilbert变换得到信号的解析形式：\nz(t)＝s(t)+jH(s(t))＝a(t)ejΦ(t)    (1)\n其中幅值函数为\n$a\left(t\right)=\sqrt{s^2\left(t\right)+H^2\left(s\left(t\right)\right)}---\left(2\right)$\n相位函数为\n$\Phi \left(t\right)=\arctan \frac{H\left[s\right(t\left)\right]}{s\left(t\right)}---\left(3\right)$\n对相位函数求导数即得到瞬时频率为\n$f\left(t\right)=\frac{1}{2\pi }\frac{\mathrm{d\Phi }\left(t\right)}{\mathrm{dt}}---\left(4\right)$\n其具体过程如下：\n①求出信号s(t)的局部极大值和极小值\n②用三次样条插值分别构成信号的上下包络线\n③求上下包络线的平均包络，记为m1(t)\n④用原信号减去m1(t)，记为h1(t)\nh1(t)＝s(t)-m1(t)    (5)\n⑤.判断h1(t)是否满足IMF条件，如满足记\nc1(t)＝h1(t)\n为第一个IMF\n若不满足，把h1(t)看作s(t)，求\nh11(t)＝h1(t)-m11(t)    (6)\nh1k(t)=h1(k-1)(t)-m1k(t)    (7)\n直到h1k(t)满足IMF条件\n⑥设c1(t)=h1k(t)为第一个IMF\n⑦记    r1(t)=s(t)-c1(t)    (8)\n将r1(t)作为信号s(t)重复①到⑥求出c2(t)，直到rn(t)一个单调函数或直流分 量。经过EMD后，信号可以表示为\n$s\left(t\right)=\stackrel{n}{\Sigma }{c}_i\left(t\right)+r_n\left(t\right)---\left(9\right)$\n对式(9)中的每个IMF作Hilbert变换可得\n$s\left(t\right)=\mathrm{Re}\underset{i=1}{\overset{n}{\Sigma }}{a}_i{\left(t\right)e}^{j{\Phi }_i\left(t\right)}$\n$=\mathrm{Re}\underset{i=1}{\overset{n}{\Sigma }}{a}_i\left(t\right)e^{{j\int \omega }_j\left(t\right)\mathrm{dt}}---\left(10\right)$\n这里省略了残差函数r，Re表示取实部。式(10)称为Hilbert幅值谱，简 称Hilbert-Huang谱，记作\n$H(\omega ,t)=\mathrm{Re\Sigma }{a}_i{\left(t\right)e}^{j\int {\omega }_l\left(t\right)\mathrm{dt}}---\left(11\right)$\n每个模态的瞬时频率fi(t)可以通过下式求出：\n$f_i\left(t\right)=\frac{1}{2\pi }{\omega }_i\left(t\right)---\left(12\right)$\n为了利用希尔伯特-黄变换方法处理两相流动信号，现定义分解后各层平均 能量为\n$E_i=\underset{i=1}{\overset{N}{\Sigma }}{\left|a_i\left(t\right)\right|}^2/N---\left(13\right)$\n其中N为每个模态的数据点数。\n分解后各层的总平均能量为\n$E=\underset{i=1}{\overset{n}{\Sigma }}{E}_i---\left(14\right)$\n定义各层所占的能量比为\n$M_i=\frac{E_i}{E}---\left(15\right)$\n图3-6为原始的不同流型信号及其固有模态分解结果，第一个子图为原始差 压信号，imf1-imf8为分解后的固有模态，最后是残差函数。表1为计算得到四 种不同的典型流型的分解后各模态的平均瞬时频率与能量比，其它情况的计算 相同。从分解结果可以看出，不同的固有模态基本上是反映信号的波动成分， 而残差则体现了信号的平稳部分。从图3-6可以看出，由单相流到混状流变化时， 信号的均值与残差逐渐减小。从信号的时域图可以看出，用文丘利管的差压信 号可以看出不同的流型信息，但是用常规的信号处理方法却很难提取出定量的 描述。\n表1.不同模态的频均瞬时频率与能量比\n  Imf1   Imf2   Imf3   Imf4   Imf5   Imf6   Imf7 单相水Ei(％) 47.1761 26.2428 13.7704 5.8723 2.6686 2.3676 1.9022 fi(Hz) 28.6273 13.9663 6.4761 3.6056 0.8736 0.8584 0.3662 泡状流Ei(％) 75.6701 8.8507 3.7417 3.4433 4.4808 1.8425 1.9709 fi(Hz) 31.5066 15.7983 8.1838 4.1270 0.9305 0.8963 1.4684 弹状流Ei(％) 26.2864 10.5941 7.2009 12.1268 22.3949 18.1222 3.2748 fi(Hz) 30.8179 16.4172 8.6372 4.1701 0.9837 0.9089 0.8047 混状流Ei(％) 2.0634 2.2615 16.3537 27.4682 29.3461 9.7548 12.7522 fi(Hz) 29.9238 14.9524 6.6244 5.5428 1.9994 1.3501 0.6859\n利用上述过程对不同流型的文丘利管的差压信号进行分解与计算，在50mm 与40mm气水两相流的实验结果分别如图7、图8所示，可以看出在3Hz附近， 能量比从混状流到弹状流、泡状流依次减小，而在30Hz附近则刚好相反。在15 Hz附近虽然与30Hz附近有相似的结果，但不是很明显。由此可以把整个分解 结果当成三个特征向量来处理，即3Hz附近所占的总的能量比E1，15Hz附近 所占的能量比E2与30Hz附近所占的能量比E3。通过计算得到在3Hz附近混状 流所占的能量比为75％～98％，弹状流为40％～90％，泡状流为45％以下，在30Hz 附近混状流占的能量比为10％以下，弹状流为8％～55％，泡状流为40％～90％。 通过计算这两个频率范围内的能量比，然后判断能量比的范围就可以识别出混 状流、弹状流、泡状流(单相水)，由于单相水的能量比与气泡很小时泡状流的 能量比在相同范围内，不能识别出来。\n为了把单相水流从泡状流中识别出来，我们分解后的残差进行分析，得到 了残差均值与流型的关系，结果如图9、图10所示，可以看出，单相流时残差 的均值在0.94以上，泡状流在0.65～0.91之间，弹状流在0.40～0.65，而混状流 时分布在0.4以下的区域。残差是把原始信号的波动成份经过EMD后得到的反 映信号变化趋势的成分，通过对50mm和40mm的差压信号的分析看出二者的 变化范围几乎是一致的，通过残差的均值可以十分容易地把单相水与泡状流区 别开。该结果可以补充利用不同的能量比不能把小泡流与单相流区别开的不足。\n利用相同的实验过程，我们在40mm和50mm的管道上进行油水两相流实 验，实验用油为0#柴油。油流量范围为4.2m3/h～7.0m3/h。其它条件与气水实验 条件一样。\n为了更加明确地反映出能量比与流型的关系，我们把50mm与40mm管径气 液两相流在3Hz附近的能量比E1与30Hz附近的能量比E3画在图11中，以便 找到与流型的关系的直观描述。图中主要有4个区域，每个区域包含了一种流 型。其中，单相流可以看作是泡状流的一种特殊形式。也就是通过该图可以很 好的识别出不同管径下泡状流、弹状流与混状流。再通过与计算的残差得均值 的结合，把单相水从泡状流中分离出来，从而实现单相流、泡状流、弹状流与 混状流的有效识别。其中，阴影部分分别为泡-弹与弹-混状流型的过渡区间。\n通过以上的分析，总结出以下的规则来实现对水平管气液两相流流型的识 别，如表2所示，通过这些规则就可以十分容易的识别出不同的流型，另外如 果同时满足两个规则，则可以判别为该流型为过渡流型，如同时满足2和3则 可判别该流型为泡-弹过渡流型。\n表2、利用能量比与残差均值识别气液两相流流型的规则\n  规则     E1(％)     E3(％)     残差均值(V)     流型   1     10～35     45～70     ＞0.94     单相水   2     3～45     40～90     ＜0.92     泡状流   3     40～90     8～55     弹状流   4     75～98     ＜15     混状流\n结果表明，利用文丘利管上的差压信号采用希尔伯特-黄变换方法进行流型识 别，可以针对不同的管道直径、不同的液相流量下得到的单相流、泡状流、弹 状流、混状流流型，利用相同的规则进行识别。该方法识别流型受实验条件影 响小，即使采用不同量程的差压变送器得到结果也是一致的。\n参考文献\n1.Huang，N.E.，Steven，Z.S.，Long，R.，et al，1998.The Empirical mode decomposition and Hilbert spectrum for non-stationary time series analysis. Proceedings of the Royal Society of London，Series A454，903-995.\n2.Leisk，G.G.Application of the Hilbert-Huang Transform to machine tool condition/health monitoring.AIP Conference Proceedings.615，1711-1718.\n3.Montesinosa，M.E.，2003.Hilbert-Huang analysis of BWR neutron detector signals： application to DR calculation and to corrupted signal analysis.Annals of Nuclear Energy 30，715-727.\n4.Veltcheva，A.D.，2002.Wave and group transformation by a Hilbert spectrum. Coastal Engineering Journal 44，283-300.\n5.熊学军，郭炳火等，2002.EMD方法和Hilbert谱分析法的应用与探讨。黄 渤海海洋.20，12-21.