CN111059477A - 一种基于双层架构反向寻源化工管道泄漏检测与定位方法 - Google Patents

Abstract

本发明公开了一种基于双层架构反向寻源化工管道泄漏检测与定位方法，包括：采集化工管道的首末端和中间位置的流量值和压力值；在第一层架构的泄漏检测与定位计算模块中计算得到泄漏发生的大致位置x_layer1和泄漏量m_layer1；根据泄漏发生的大致位置x_layer1和泄漏量m_layer1计算确定泄漏孔径d_L的搜索范围；根据泄漏位置率k的搜索范围、泄漏量m的搜索范围以及泄漏孔径d_L的搜索范围在第二层架构的泄漏检测与定位计算模块中对化工管道模型进行反向寻源计算，最终获得精确的泄漏定位和泄漏量预测。本发明的方法比常规单层结构的反向寻源方法，具有收敛速度快、定位精度高的优点。

Description

一种基于双层架构反向寻源化工管道泄漏检测与定位方法

技术领域

本发明涉及化工安全生产管道泄漏检测与定位技术领域，具体涉及一种基于双层架构反向寻源化工管道泄漏检测与定位方法。

背景技术

随着化工行业的发展和化工装置的进步，化学工业发展十分迅速，其中管道泄漏检测与定位系统(LDS)是化学工业安全运行的重要保障，目前已有多种管道泄漏检测与定位方法，这些方法大致可分为硬件检测法和软件检测法两大类。其中，基于硬件的方法主要是指对于泄漏物进行直接检测，如采用光纤传感器、声学传感器、漏磁技术等；基于软件的方法又称为间接检测法，主要借助于计算机系统，通过检测因泄漏造成的影响，如流量、压力、温度等物理参数的变化来判断泄漏是否发生及对泄漏进行定位。早期的检测方法主要以负压波法(NPW)为主，其基本原理为：当某带压设备中某个位置发生泄漏，在泄漏口处会瞬间喷射处大量流体，根据能量守衡的原则，泄漏口周围的流体就会向泄漏口补充流体介质，导致泄漏口处产生一个瞬变的水锤减压波(即负压波)，这种能量会向管道的上下游传播，使管路沿线在一定程度上压力降低，通过采集首末端接收的压力变化进行泄漏判断，并通过接收时间差判断漏点位置。NPW方法检测精度高、响应快，并能进行长距离检测，但不适用于短距离的漏点检测，并且对缓慢的、微弱的泄漏无法进行有效的检验。针对这个问题，有人将统计方法序贯概率比(SPRT)应用于故障检测检测判断中，该方法基于从管道观测数据对泄漏或无泄漏的估计概率进行比较，当压力值较长时间偏离正常，则泄漏的概率就会增加，当概率远大于未发生泄漏的概率时，就会发生泄漏报警，该方法能够对小泄漏进行有效地检测，并能及时获得泄漏报警。SPRT对于瞬态过程也具有良好的检测性能。基于瞬态的泄漏检测方法得到了很多发展，该方法利用瞬态流动的水力学特性来检测管道中的泄漏。具体方法是对管道系统引入适当带宽和振幅的压力波，然后通过分析压力测量点的数据分析系统响应。这种基于瞬态的方法因其快速，易在线实施和大范围操作范围的优势而变得流行。目前，有四种典型的基于瞬态的方法，包括：(1)基于瞬态波反射的方法；(2) 基于瞬态波阻尼的方法；(3)基于瞬态频率响应的方法；(4)基于逆瞬态分析的方法。其中，基于ITA的方法是这四种方法中最鲁棒和最综合的方法，但其准确性和效率高度依赖于历史数据和质量，以及所采用的逆分析算法。 Liggett首次将逆瞬态分析方法(ITA)应用于管道泄漏检测和定位。 Sarkamaryan等引入高斯函数来模拟每个候选泄漏，并将高斯函数参数作为ITA决策变量，降低了传统ITA的维度问题。

尽管基于瞬态流的方法在一定程度上提高了缓慢、微小泄漏的检测能力和定位精度，但是仍然存在很多问题，基于瞬态的方法对于模型的精度也有较高的要求，然而化工过程管道由于具有直段短、分之多、异型组件复杂的特点，会对所建模型的精度有很大的影响；另外，基于瞬态方法的应用需要外部给定的瞬态激励，如，人为地调节阀门开度。这就带来人为的外部干扰会对管道的正常运行带来一定的影响，导致该方法在化工过程中的实际应用中受到了一定的限制。为了解决这个问题，Zhang等结合ITA和改进的PSO 算法，从准确性、稳定性、鲁棒性、误报率等四个指标来验证方法有效性。该方法有效避免了外部瞬态激励，具有计算精度高、稳定性好、误报率低、抗噪声能力强等优点。但该方法计算机较大，且容易在反向寻源中陷入局部最优。

发明内容

本发明针对化工过程管道提出了一种基于双层架构反向寻源化工管道泄漏检测与定位方法，能够有效地解决化工过程故障检测领域中，化工管道存在较多弯头和其他局部阻力附件的情况下，传统的检测与定位方法不能有效进行准确定位和高效检测的问题，通过本发明公开的一种基于双层架构反向寻源化工管道泄漏检测与定位方法，进一步提高了基于模型法在化工过程管道故障检测的可靠性。

一种基于双层架构反向寻源化工管道泄漏检测与定位方法，其特征在于，包括以下步骤：

步骤1.1：采集化工管道的首末端和中间位置的流量值和压力值；

步骤1.2：在第一层架构的泄漏检测与定位计算模块中计算得到泄漏发生的大致位置x_{layer 1}和泄漏量m_{layer 1}；

步骤1.3：根据步骤1.2获得的泄漏发生的大致位置x_{layer 1}和泄漏量 m_{layer 1}计算确定泄漏孔径d_L的搜索范围；

步骤1.4：根据步骤1.3中的泄漏位置率k的搜索范围、泄漏量m的搜索范围以及泄漏孔径d_L的搜索范围在第二层架构的泄漏检测与定位计算模块中对化工管道模型进行反向寻源计算，最终获得精确的泄漏定位和泄漏量预测。

以下作为本发明的优选技术方案：

步骤1.2中，在第一层架构的泄漏检测与定位计算模块中计算得到泄漏发生的大致位置x_{layer 1}和泄漏量m_{layer 1}，具体包括：

将化工管道的首末端的压力值和流量值输入到第一层架构的泄漏检测与定位计算模块中进行计算，其中，化工管道的首末端的压力值通过多序贯概率比检验法实时获取与正常序列不同的似然比来初步判断管道泄漏是否发生，若泄漏发生，则利用负压波定位(NPW)方法初步判断泄漏发生的大致位置x_{layer 1}；

同时，化工管道的首末端的流量值通过第一层架构的泄漏检测与定位计算模块中质量流量平衡式初步获得泄漏量m_{layer 1}；

步骤1.3中，根据步骤1.2获得的泄漏发生的大致位置x_{layer 1}和泄漏量 m_{layer 1}计算确定泄漏孔径d_L的搜索范围，具体包括：

首先，定义泄漏位置率k＝x_L/L，其中，L为化工管道的总长，单位m，x_L为泄漏点距离首端的距离，单位m，将步骤1.2第一层架构得到泄漏发生的大致位置x_{layer 1}，计算得到泄漏位置率k_{layer 1}＝x_{layer 1}/L，并将k的搜索范围的上限设置为泄漏位置率k_{layer 1}的0.5倍，下限为k_{layer 1}的2倍，即 k∈[0.5k_layer1，2k_layer1]若2k_layer1＞1，则2k_layer1＝1，另外，将泄漏量m的搜索范围按照步骤1.2中得到的泄漏量m_{layer 1}的0.5倍作为区间的左边界，将泄漏量 m_{layer 1}的2倍作为区间的右边界，即m∈[0.5m_layer1，2m_layer1]，若

则

最后将泄漏孔径d_L的搜索范围设置为d_L∈(0，D]，其中D为化工管道的管径。

步骤1.4中，根据步骤1.3中的泄漏位置率k的搜索范围、泄漏量m的搜索范围以及泄漏孔径d_L的搜索范围在第二层架构的泄漏检测与定位计算模块中对化工管道模型进行反向寻源计算，最终获得精确的泄漏定位和泄漏量预测，具体包括：

(1).确定管路系统基本参数：化工管道的总长L，化工管道的管径D，绝对粗糙度，流体密度dens，压力波速a和首末端流量、压力测量数据；

(2).划分计算时间网格和计算空间网格，计算时间网格由信号采集频率决定，计算空间网格根据化工管道的总长L划分为长度为Δx的n段，根据化工管道模型采用四阶龙格库塔法计算管路沿程的初始流量和压力分布；

(3).设定泄漏位置率k初始值，将化工管道的总长L分为L1和L2两段，其中，L1＝k*L、L2＝L-L1两段，由计算空间网格获得两段L1和L2的计算节点数分别为N1和N2，其中，两段管路L1和L2分别包含两个边界条件，L1段的边界条件为上游压力、流量实时数据和左泄漏边界条件，L2 的边界条件为右泄漏边界条件和下游压力、流量实时数据；

(4).根据式(4d)和(4e)分别计算管道L1段和L2段内的压头H和流速v分布，其中，当边界条件为H时，L1段左边界由式(4d)求得L1 段流速的值v_u，t+Δt，右边界由顺特征线经过变形得到式(4e)，计算出L1段泄漏点左端处的压头H_d，t+Δt，同时，计算得到L2段压头的值H_Lk+和L2段泄漏点右端处的流速v_Lk+；

v_Lk＝v_Lk-+v_Lk+ (4a)

H_Lk＝H_Lk-＝H_Lk+ (4b)

其中式4a中，v_Lk表示泄漏点出口处的流速，v_Lk-表示泄漏点左端管道的流速，v_Lk+表示泄漏点右端管道的流速，单位m/s；式4b中，H_Lk表示泄漏点出口处的压头，H_Lk-表示泄漏点左端管道口处的压头，H_Lk+表示泄漏点右端管道口处的压头，单位m；式4c中，Q_L表示泄漏量，单位m³/S，C_d表示小孔系数，g表示重力加速度，Z_L表示泄漏点处的高程，单位m；式4d中，v_u，t+Δt表示t+Δt时刻通过计算得到的管道上游的流速，v_1，t表示t时刻管道节点1处的流速，H_1，t表示t时刻管道节点1处的压头，a表示压力波速，

表示t+Δt时刻测量得到的管道上游压头，α表示管道倾斜角，λ表示管道的流体阻力系数，D表示化工管道的管径；其中式4e中，H_d，t+Δt表示 t+Δt时刻通过计算得到的管道下游的压头，H_n，t表示t时刻管道第n个节点处的压头，v_n，t表示t时刻管道第n个节点处的流速，

表示t+Δt时刻通过测量得到的管道下游的流速；

(5).泄漏点处的压头和流速关系满足式(4a)(4b)(4c)关系，其中，式(4c)为泄漏位置压头与泄漏量的关系，由式(6)计算得到t+Δt时刻的泄漏量Q_Lc，t＝(v_L-，t-v_L+，t)×A_dL，A_dL为泄漏孔截面积，根据收敛条件定义式(5)可以得到t时刻的相对误差之和Erro_VH；

其中，式(5)中，H_0，t+Δt表示t+Δt时刻首端计算的压头，

表示 t+Δt时刻首端实际测量的压头，V_n+1，t+Δt表示t+Δt时刻末端计算的流速，

表示t+Δt时刻末端实际测量的流速，Q_Lc，t+Δt表示t+Δt时刻的泄漏量Q_Lc，t+Δt，η_VH表示阈值；

式(6)中，H_L表示泄漏点处的高程的压头，Z_L表示泄漏点处的高程， C_d表示孔流系数；

(6).判断式(5)中Erro_VH是否满足收敛条件，即小于阈值η_VH，若满足，则t+Δt时刻反向寻源结束，得到泄漏信息：泄漏位置率k、泄漏孔径的大小d_L和泄漏量Q_Lc，t+Δt；若不满足，则使用混沌粒子群进行优化寻优，混沌粒子群将在k∈[0.5k_layer1，2k_layer1]、m∈[0.5m_layer1，2m_layer1]、d_L∈(0，D]，(D为管径)该搜索区域内寻优直到得到满足收敛条件的泄漏位置率k、泄漏孔径大小d_L和泄漏量Q_Lc，t+Δt值。

与现有技术相比，本发明具有如下优点：

本发明的一种基于双层架构反向寻源化工管道泄漏检测与定位方法比常规单层结构的反向寻源方法，具有收敛速度快、定位精度高的优点。首先利用多序贯概率比(M-SPRT)结合负压波定位法和质量流量平衡式能够初步快速地对管路进行泄漏点定位和泄漏量判断，然后根据本发明提出的搜索空间设计方法，能够有效提高第二层中利用实时瞬态模型(RTTM)进行反向寻源的效率，同时利用混沌粒子群(CPSO)智能优化算法进行快速寻优，减少了误报率和漏报率的发生，能够在第一层泄漏检测和定位的基础上，进一步提高了定位的精度和泄漏量与预测的准确性。

附图说明

图1为双层架构反向寻源化工管道泄漏检测与定位流程图；

图2a为1.2mm泄漏孔径下，泄漏率在不同流量下的变化情况图；

图2b为3mm泄漏孔径下，泄漏率在不同流量下的变化情况图；

图2c为6mm泄漏孔径下，泄漏率在不同流量下的变化情况图；

图3为上下游压力差在三种孔径下的分布图；

图4为管道首末端测量的流量和压力变化趋势图图；

图5为M-SPRT在ln A＝4.6对首末端压力的检测效果对比图；

图6为首末端流量和压头的测量值和计算值随时间的分布图，其中，图 6(a)为上游流量的测量值和计算值随时间的分布图，图6(b)为下游压头的测量值和计算值随时间的分布图；

图7a为位置率k的寻源结果图；

图7b为泄漏孔径dl的寻源结果图；

图7c为泄漏量ql的泄漏结果图。

具体实施方式

如图1所示，一种基于双层架构反向寻源化工管道泄漏检测与定位方法，其特征在于，包括以下步骤：

步骤1.1：采集化工管道的首末端和中间位置的流量值和压力值；

步骤1.2：将化工管道的首末端的压力值和流量值输入到第一层架构的泄漏检测与定位计算模块中进行计算，其中，化工管道的首末端的压力值通过多序贯概率比检验法实时获取与正常序列不同的似然比来初步判断管道泄漏是否发生，若泄漏发生，则利用负压波定位(NPW)方法初步判断泄漏发生的大致位置；同时，化工管道的首末端的流量值通过质量流量平衡式初步获得泄漏量；得到泄漏发生的大致位置x_{layer 1}和泄漏量m_{layer 1}；

其中，多序惯概率比的备择假设Δμ的设计主要通过实验分析备择假设与泄漏前后首末端压力的关系，备则假设的设计要最大程度满足各种规模的泄漏，通过实验平台分析不同规模(微型泄漏、正常规模泄漏和大规模泄漏) 下，管道首末端的压力差情况，图2a、2b、2c分别表示在1.2mm、3mm和 6mm三种泄漏孔径下，泄漏率随着管道流量增大(3.1m3/h～4.6m3/h)的关系。

从实验数据可以看出，在不同流量下，同一种泄漏孔径下，不同流量对泄漏率的影响不大，随着流量的增大，泄漏率依旧在一个区间内变化，没有明显的上升或者下降的趋势；而对于不同的泄漏孔径，泄漏率却有很大的差别，1.2mm泄漏孔径的泄漏率平均在2％，3mm泄漏孔径的泄漏率平均为6％， 6mm泄漏孔径的泄漏率平均为20％。因此，可以得出泄漏规模单向取决于泄漏孔径大小，和管道主流量并无决定性影响，因此选择合适的备则假设只需要考虑泄漏孔径的大小造成的泄漏变化，由于本节主要研究不同泄漏引起的压力变化，因此本文针对以上三种泄漏孔径，得出了不同流量下，上游压力和下游压力变化量的分布图，如图3所示。

从图中可以看出在1.2mm的泄漏孔径，泄漏率约为2％的管道中，上、下游压力变化量在0.5m范围内；在3mm的泄漏孔径，泄漏率约为6％的管道中，上、下游压力变化量在1.5m范围内；在6mm的泄漏孔径，泄漏率约为20％的管道中，上、下游压力变化量在4m范围内。三种泄漏孔径分别代表微型泄漏、常规规模泄漏和大规模泄漏的水平，根据实验数据，可以对M-SPRT进行备则假设的设计，分别设计如下7种SPRT的备择假设： H₁∶Δμ₁＝0.5m、H₂∶Δμ₂＝1m、H₃∶Δμ₃＝2m、H₄∶Δμ₄＝3m、H₅∶Δμ₅＝4m、 H₆∶Δμ₆＝5m和H₇∶Δμ₇＝6m。

另外，步骤1.2中，根据负压波定位(NPW)方法初步获得泄漏的大致位置，结合质量流量平衡式初步获得泄漏量，即当泄漏发生时候，泄漏产生的压力波传播到入口处所用的时间为

判断泄漏发生的时间为t_in；同样压力波传播到出口处所用的时间为

判断泄漏发生的时间为t_out，则时间差为

最终可以根据NPW方法，即式(3)初步获得管道发生泄漏的大概位置；

对于质量流量平衡的计算，主要通过采集的首末端流量，建立质量流量平衡式

当通过M-SPRT检测出首末端压力序列可能有泄漏发生的信号，可以根据质量流量平衡式判断

(

为测量误差)是否成立，若成立则判断泄漏；进一步地，当

稳定后，可以初步确定该值为泄漏量

步骤1.3：首先，定义泄漏位置率k＝x_L/L，其中，L为化工管道的总长，单位m，x_L为泄漏点距离首端的距离，单位m。将步骤1.2第一层架构得到泄漏发生的大致位置x_{layer 1}，计算得到泄漏位置率k_{layer 1}＝x_{layer 1}/L，并将k 的搜索范围的上限设置为泄漏位置率k_{layer 1}的0.5倍，下限为k_{layer 1}的2倍，即k∈[0.5k_layer1，2k_layer1]若2k_layer1＞1，则2k_layer1＝1，另外，将泄漏量m的搜索范围按照步骤1.2中得到的泄漏量m_{layer 1}的0.5倍作为区间的左边界，将泄漏量 m_{layer 1}的2倍作为区间的右边界，即m∈[0.5m_layer1，2m_layer1]，若

则

最后将泄漏孔径d_L的搜索范围设置为 d_L∈(0，D]，(D为管径)；

步骤1.4：根据步骤1.3中的泄漏位置率k的范围、泄漏量m的范围以及泄漏孔径d_L的搜索范围在第二层架构的泄漏检测与定位计算模块中对化工管道模型，如式(4a，4b，4c，4d，4e)进行反向寻源计算，最终获得精确的泄漏定位和泄漏量预测。

所述的第一层架构，即步骤1.2，1.3，泄漏检测与定位方法中使用多序贯概率比检验法判断泄漏时刻、负压波方法(NPW)定位，质量流量平衡式计算泄漏量包括以下步骤：

step1.将实时采集的管道首端和末端压头序列

经过卡尔曼滤波之后输入到多序贯概率比(简称M-SPRT)检验式(1)，可以得到第j 个SPRT的似然比；

其中P_0，t和P_j，t表示通过首末端压头序列计算的第j个序贯概率比检验式的概率密度，其计算如式(2)所示

式中σ表示H_u，1(或H_d，1)至H_u，m(或H_d，m)压头序列的方差，μ₀表示H_u，1(或 H_d，1)至H_u，m(或H_d，m)的正常序列的均值，Δμ_j表示第j个SPRT的备择假设，即原假设：μ＝μ₀表示管道正运行，无泄漏，当μ＝μ₀+Δμ表示第i个备则假设，发生泄漏。其中Δμ＞0表示压力上升，Δμ＜0表示压力下降。当λ_j满足不等式

那么接受零假设μ＝μ₀，判断管道运行正常；如果任意一个备择假设的似然比满足不等式

那么拒绝零假设，则判断泄漏发生；如果这个比值满足不等式

不能作出接受还是拒绝零假设的判断，必须继续观察直到这个比值满足前两种情况为止，规定每个SPRT的α和β分别为犯第一类和第二类错误的概率，即允许的误报率和漏报率；

step2.当多序贯概率比检验法分别由step1判断出泄漏发生的时刻后，根据负压波定位法(NPW)初步获得泄漏的大致位置如式(3)，其中L表示管道总长，a为负压波波速，Δt表示泄漏产生的压力波传播到入口处所用的时间为

判断泄漏发生的时间为t_in，同样压力波传播到出口处所用的时间为

判断泄漏发生的时间为t_out，则时间差为

根据序惯概率比泄漏发生的负压波传递到首末端位置的时间差；

step3.将实时采集的管道首端和末端流速序列

输入到质量流量平衡式

中，判断

(

为测量误差)是否成立，若成立则判断泄漏；进一步地，当

稳定后，可以初步确定该值为泄漏量

step4.通过上述步骤初步得到泄漏位置和泄漏量，为了第二层寻源的计算，需要确定搜索空间。首先定义泄漏位置率k_{layer 1}＝x_{layer 1}/L，并将k的搜索范围的上限设计k_{layer 1}的0.5倍，下限为k_{layer 1}的2倍，即k∈[0.5k_layer1，2k_layer1]， (若2k_layer1＞1，则2k_layer1＝1)，另外，将泄漏量m的搜索范围按照步骤1.2中初步得到的结果的0.5倍作为区间的左边界，将该值的2倍作为区间的右边界，即m∈[0.5m_layer1，2m_layer1]，(若

则

)，最后将泄漏孔径的搜索范围设置为d_L∈(0，D]，(D为管径)。

步骤1.4中，对化工管道模型进行反向寻源计算具体包括：

Step1.确定管路系统基本参数：化工管道的总长L，管径D，绝对粗糙度，流体密度dens，压力波速a和首末端流量、压力测量数据；

Step2.划分计算时间网格和计算空间网格，计算时间网格由信号采集频率决定，计算空间网格根据化工管道的总长L划分为长度为Δx的n段，根据化工管道模型采用四阶龙格库塔法计算管路沿程的初始流量和压力分布；

Step3.设定泄漏位置率k初始值，将化工管道的总长L分为L1和L2两段，其中，L1＝k*L、L2＝L-L1两段，由计算空间网格获得两段L1和L2的计算节点数分别为N1和N2，其中，两段管路L1和L2分别包含两个边界条件，L1段的边界条件为上游压力、流量实时数据和左泄漏边界条件，L2 的边界条件为右泄漏边界条件和下游压力、流量实时数据；

Step4.根据式(4d)和(4e)分别计算管道L1段和L2段内的压头H和流速v分布，其中，当边界条件为H时，L1段左边界由式(4d)求得L1 段流速的值v_u，t+Δt，右边界由顺特征线经过变形得到式(4e)，计算出L1段泄漏点左端处的压头H_d，t+Δt，同时，计算得到L2段压头的值H_Lk+和L2段泄漏点右端处的流速v_Lk+；

v_Lk＝v_Lk-+v_Lk+ (4a)

H_Lk＝H_Lk-＝H_Lk+ (4b)

其中式4a中，v_Lk表示泄漏点出口处的流速，v_Lk-表示泄漏点左端管道的流速，v_Lk+表示泄漏点右端管道的流速，单位m/s；式4b中，H_Lk表示泄漏点出口处的压头，H_Lk-表示泄漏点左端管道口处的压头，H_Lk+表示泄漏点右端管道口处的压头，单位m；式4c中，Q_L表示泄漏量，单位m³/S，C_d表示小孔系数，Z_L表示泄漏点处的高程，单位m；式4d中，v_u，t+Δt表示t+Δt时刻通过计算得到的管道上游的流速，v_1，t表示t时刻管道节点1处的流速，H_1，t表示t时刻管道节点1处的压头，a表示压力波速，

表示t+Δt时刻测量得到的管道上游压头，α表示管道倾斜角，λ表示管道的流体阻力系数，D 表示管道直径；其中式4e中，H_d，t+Δt表示t+Δt时刻通过计算得到的管道下游的压头，H_n，t表示t时刻管道第n个节点处的压头，v_n，t表示t时刻管道第n 个节点处的流速，

表示t+Δt时刻通过测量得到的管道下游的流速。

Step5.泄漏点处的压头和流速关系满足式(4a)(4b)(4c)关系，其中，式(4c)为泄漏位置压头与泄漏量的关系，由式(6)计算得到t+Δt时刻的泄漏量Q_Lc，t＝(v_L-，t-v_L+，t)×A_dL，A_dL为泄漏孔截面积。根据收敛条件定义式(5)可以得到t时刻的相对误差之和Erro_VH；

VH边界条件：

其中，式(5)中，H_0，t+Δt表示t+Δt时刻首端计算的压头，

表示 t+Δt时刻首端实际测量的压头，V_n+1，t+Δt表示t+Δt时刻末端计算的流速，

表示t+Δt时刻末端实际测量的流速，Q_Lc，t+Δt表示t+Δt时刻的泄漏量Q_Lc，t+Δt，η_VH表示阈值；

式(6)中，H_L表示泄漏点处的高程的压头，Z_L表示泄漏点处的高程， C_d表示孔流系数；

Step6.判断式(5)中Erro_HV是否满足收敛条件，即小于阈值η_VH，若满足，则t+Δt时刻反向寻源结束，得到泄漏信息：泄漏位置率k、泄漏孔径的大小d_L和泄漏量Q_Lc，t+Δt；若不满足，则使用混沌粒子群进行优化寻优，混沌粒子群将在k∈[0.5k_layer1，2k_layer1]、m∈[0.5m_layer1，2m_layer1]、d_L∈(0，D]，(D为管径)该搜索区域内寻优直到得到满足收敛条件的泄漏位置率k、泄漏孔径大小d_L和泄漏量Q_Lc，t+Δt值。

针对本发明的一种基于双层架构反向寻源化工管道泄漏检测与定位方法做详细描述。本发明的实施方法如下：

首先搭建实验平台，包括管道设备平台、支撑固定框架、信号检测采集模块和信号加工处理模块。管道设备平台主要包括水槽、多级泵、主管路和若干球阀等，管道主体由32个标准90°弯头和总长约为45m的直管段交叉组成，管材为管径0.025m、厚度为1.6mm的不锈钢管，管内输送介质为自来水；信号检测采集模块主要由数据采集卡和信号传感器组成，数据采集卡型号为ATR，信号传感器指管道首末两端安装的一个涡轮流量计 (0.6-6m3/h，0.5级)和一个压力传感器(0-0.6Mpa，0.2级)以及管道中间距首端位置20m地方的一个压力传感器(0-0.6Mpa，0.2级)。

本发明使用ATR数据采集卡对管道的流量和压力进行实时监控，数据采集时间为30s，首先将管道流量流量设定为3.89m3/h，首末端压力和流量传感器开始采集数据，15.5s时开始泄漏，泄漏孔径为3mm，泄漏位置距首端的13m处，当数据采集30s后结束采样，得到管道首末端的流量、压力变化趋势如图4所示。

本发明第一层架构将利用上文设计的备择假设对实验平台的管道泄漏进行分析。对首末端压力数据采取以上设计的7种SPRT的备择假设： H₁∶Δμ₁＝0.5m、H₂∶Δμ₂＝1m、H₃∶Δμ₃＝2m、H₄∶Δμ₄＝3m、H₅∶Δμ₅＝4m、H₆∶Δμ₆＝5m 和H₇：Δμ₇＝6m进行泄漏检测，误报率和漏报率均为1％，门限取值为ln A＝4.6，分析结果如图5所示

从图中可以得出，当备择假设设置的过小(SPRT1、SPRT2)时，所求最大似然比一开始就超过门限值(ln A＝4.6)，不会产生报警；当备则假设设置到合适范围(SPRT3、SPRT4、SPRT5)时，SPRT检测到首末端的压力序列约在T＝15.5s处达到门限值并开始报警；当备则假设设置得过大时(SPRT6、 SPRT7)时，所求最大似然比几乎为零，SPRT始终不会报警，表1列出了7 组SPRT检测获得的报警情况。

因此，设置不同程度的备则假设，有利于SPRT在对压力信号进行泄漏检测时，更能有效地对不同规模的未知的泄漏进行检测判断。当M-SPRT设计不同规模的备则假设时，当前几组备择假设由于设计了较小的备择假设，导致了一开始检测的最大似然比就超过规定门限值，而无法正常报警；另外，设计较大的备择假设会使得系统的最大似然比一直处于很低水平(几乎为零)，而无法达到门限，始终不会报警，而产生漏报警。因此，可以选择介于前面几组和漏报警中间的备择假设作为正常报警的泄漏判断，其检测得到的首末端时间可以初步判断为泄漏发生的时间，最终通过式(3)初步确定泄漏位置，再结合质量平衡式可以初步计算出泄漏量，如表1所示为本次实验通过7组SPRT检测首末端发生泄漏的时刻、泄漏位置和泄漏量。

表1 M-SPRT检测的泄漏时刻及泄漏位置

根据上一层利用M-SPRT结合NPW和质量流量平衡法初步获得泄漏位置约为22m，即k约为0.52、泄漏量约为0.214m³/h、泄漏时刻约为15.6s。当判断泄漏发生后，进入第二层基于混沌粒子群的反问题寻源进行泄漏定位和泄漏量预测，首先将泄漏位置率k的搜索范围的上限设计为第一层预测值 k的0.5倍，下限为k的2倍(若2k＞1，则k＝1)，即k∈[0.26，1]，另外，将泄漏量的搜索范围按照第一层泄漏量粗测结果的0.5倍作为区间的左边界，将该值的2倍作为区间的右边界，即Q_L∈(0.107，0.428)，最后将泄漏孔径的搜索范围设置为d_L∈(0，D]，D＝0.023m。确定三个源参数的搜索范围后就可以进行混沌粒子群的反向寻源定位与泄漏量预测。

已知泄漏发生前，管道内流量为3.889m3/h，当泄漏发生后，泄漏口上游的管内流量增加为3.965m3/h，泄漏口下游的管内流量降低为3.731m3/h，则泄漏量为0.234m3/h，泄漏位置为距管道首端约13m处(k＝0.3)，泄漏孔径的大小为0.003m。对管路利用混沌粒子群进行反向寻源时，规定粒子群算法的粒子数nPSO＝30，最大迭代次数为200，混沌算法的混沌变量数 nChaos＝20，以HV边界条件计算的管道首末端的流量、压头所满足的收敛条件η_HV＝0.1，采样时间为0.1s。当判断泄漏发生后，经过混沌粒子群优化后的首末端流量和压头计算值随时间分布如图6所示

从图中可以看出，蓝色实线为从t＝0～30s内采集的首末端流量和压头的数值，红色实线为经过混沌粒子群反向寻源得到的在判断泄漏发生后的首端流量和末端压力值。当第一层M-SPRT和流量平衡法检测到泄漏发生后(t＝15.8s)，第二层的反向寻源开始运算，两幅图中从泄漏时刻t＝15.8s至 t＝19s为寻源的过渡时间，该段数值均远远超过目标值，混沌粒子群无法搜寻到最佳目标，直到满足超过最大迭代次数的条件才进入下一个采样时刻的搜索，因此大约t＝19s后数值才稳定在目标值附近。经过统计，混沌粒子群优化的首末端流量与压头的值相比于实际测量的值，在首端流量下最大的相对误差为0.61％，平均相对误差为3.89％；在末端压头下最大的相对误差为 9.5％，平均相对误差为0.75％。

最终，通过本发明公开的双层架构就能有效实现对化工过程管道泄漏的检测与定位，图7a，7b，7c分别表示第二层架构经过混沌粒子群反向寻源得到的三个源参数k、dl和ql的计算结果。

从图7中可以看出，当第一层检测机构判断泄漏后，第二层架构的反向寻源算法开始计算源参数，在t＝15.8s～19s之间的过渡阶段无法寻到匹配的数值，寻源算法迭代达到最大收敛次数，直到t＝19s之后反向寻源算法才获得较为稳定的结果。泄漏点距离位于距首端12.8m处，占实际管道总长的30％，即k＝0.3，通过反向寻源算法对源参数k计算的结果为0.3224，相对误差为 7.45％；泄漏孔径为0.003m，通过反向寻源算法对源参数dl搜索的结果为 0.0034m，相对误差为13.31％；当泄漏稳定后，泄漏量为0.214m³/h，通过反向寻源算法对源参数ql搜索的结果为0.2317m³/h，相对误差为8.33％。

本发明公开了一种基于双层架构反向寻源化工管道泄漏检测与定位方法，针对化工过程管道直段短、分之多、异型组件复杂的特点，首先有效利用了M-SPRT方法对微小和缓慢的泄漏有较好的检测能力，并确定了 M-SPRT备择假设的选择可以根据泄漏孔径大小(1.2mm、3mm和6mm)，设计7种SPRT的备择假设：H₁∶Δμ₁＝0.5m、H₂∶Δμ₂＝1m、H₃∶Δμ₃＝2m、H₄∶Δμ₄＝3m、H₅∶Δμ₅＝4m、H₆∶Δμ₆＝5m和H₇∶Δμ₇＝6m，同时结合负压波方法能初步快速获得泄漏定位，利用质量流量平衡法能初步快速确认泄漏流量，实现了这类管道结构的初步故障检测和定位，然后第二层架构结合 RTTM模型和CPSO组成一种反向寻源的泄漏检测方法，该层利用混沌粒子群算法的特点能在整个搜索域中，搜索全局最优值而不会陷入局部最优，提高了泄漏定位和泄漏量预测的精度，同时具备高效的检测能力，最后有效地实现了对化工过程管道在线检测与定位的功能。总体来说，本发明公开的一种双层架构反向寻源化工管道泄漏检测与定位方法能够有效实现化工过程管道的故障检测与定位。

Claims (4)

1.一种基于双层架构反向寻源化工管道泄漏检测与定位方法，其特征在于，包括以下步骤：

步骤1.1：采集化工管道的首末端和中间位置的流量值和压力值；

步骤1.2：在第一层架构的泄漏检测与定位计算模块中计算得到泄漏发生的大致位置x_layer1和泄漏量m_layer1；

步骤1.3：根据步骤1.2获得的泄漏发生的大致位置x_layer1和泄漏量m_layer1计算确定泄漏孔径d_L的搜索范围；

步骤1.4：根据步骤1.3中的泄漏位置率k的搜索范围、泄漏量m的搜索范围以及泄漏孔径d_L的搜索范围在第二层架构的泄漏检测与定位计算模块中对化工管道模型进行反向寻源计算，最终获得精确的泄漏定位和泄漏量预测。

2.根据权利要求1所述的一种基于双层架构反向寻源化工管道泄漏检测与定位方法，其特征在于，步骤1.2中，在第一层架构的泄漏检测与定位计算模块中计算得到泄漏发生的大致位置x_layer1和泄漏量m_layer1，具体包括：

将化工管道的首末端的压力值和流量值输入到第一层架构的泄漏检测与定位计算模块中进行计算，其中，化工管道的首末端的压力值通过多序贯概率比检验法实时获取与正常序列不同的似然比来初步判断管道泄漏是否发生，若泄漏发生，则利用负压波定位方法初步判断泄漏发生的大致位置x_layer1；

同时，化工管道的首末端的流量值通过第一层架构的泄漏检测与定位计算模块中质量流量平衡式初步获得泄漏量m_layer1。

3.根据权利要求1所述的一种基于双层架构反向寻源化工管道泄漏检测与定位方法，其特征在于，步骤1.3中，根据步骤1.2获得的泄漏发生的大致位置x_layer1和泄漏量m_layer1计算确定泄漏孔径d_L的搜索范围，具体包括：

首先，定义泄漏位置率k＝x_L/L，其中，L为化工管道的总长，单位m，x_L为泄漏点距离首端的距离，单位m，将步骤1.2第一层架构得到泄漏发生的大致位置x_layer1，计算得到泄漏位置率k_layer1＝x_layer1/L，并将k的搜索范围的上限设置为泄漏位置率k_layer1的0.5倍，下限为k_layer1的2倍，即k∈[0.5k_layer1,2k_layer1]，若2k_layer1>1,则2k_layer1＝1，另外，将泄漏量m的搜索范围按照步骤1.2中得到的泄漏量m_layer1的0.5倍作为区间的左边界，将泄漏量m_layer1的2倍作为区间的右边界，即m∈[0.5m_layer1,2m_layer1]，若

则

最后将泄漏孔径d_L的搜索范围设置为d_L∈(0,D]，其中D为化工管道的管径。

4.根据权利要求1所述的一种基于双层架构反向寻源化工管道泄漏检测与定位方法，其特征在于，步骤1.4中，根据步骤1.3中的泄漏位置率k的搜索范围、泄漏量m的搜索范围以及泄漏孔径d_L的搜索范围在第二层架构的泄漏检测与定位计算模块中对化工管道模型进行反向寻源计算，最终获得精确的泄漏定位和泄漏量预测，具体包括：

(1).确定管路系统基本参数：化工管道的总长L，化工管道的管径D，绝对粗糙度，流体密度dens，压力波速a和首末端流量、压力测量数据；

(2).划分计算时间网格和计算空间网格，计算时间网格由信号采集频率决定，计算空间网格根据化工管道的总长L划分为长度为Δx的n段，根据化工管道模型采用四阶龙格库塔法计算管路沿程的初始流量和压力分布；

(3).设定泄漏位置率k初始值，将化工管道的总长L分为L1和L2两段，其中，L1＝k*L、L2＝L-L1两段，由计算空间网格获得两段L1和L2的计算节点数分别为N1和N2，其中，两段管路L1和L2分别包含两个边界条件，L1段的边界条件为上游压力、流量实时数据和左泄漏边界条件，L2的边界条件为右泄漏边界条件和下游压力、流量实时数据；

(4).根据式(4d)和(4e)分别计算管道L1段和L2段内的压头H和流速v分布，其中，当边界条件为H时，L1段左边界由式(4d)求得L1段流速的值v_u,t+Δt，右边界由顺特征线经过变形得到式(4e)，计算出L1段泄漏点左端处的压头H_d,t+Δt，同时，计算得到L2段压头的值H_Lk+和L2段泄漏点右端处的流速v_Lk+；

v_Lk＝v_Lk-+v_Lk+ (4a)

H_Lk＝H_Lk-＝H_Lk+ (4b)

其中式4a中，v_Lk表示泄漏点出口处的流速，v_Lk-表示泄漏点左端管道的流速，v_Lk+表示泄漏点右端管道的流速，单位m/s；式4b中，H_Lk表示泄漏点出口处的压头，H_Lk-表示泄漏点左端管道口处的压头，H_Lk+表示泄漏点右端管道口处的压头，单位m；式4c中，Q_L表示泄漏量，单位m³/s，C_d表示小孔系数，g表示重力加速度，Z_L表示泄漏点处的高程，单位m；式4d中，v_u,t+Δt表示t+Δt时刻通过计算得到的管道上游的流速，v_1,t表示t时刻管道节点1处的流速，H_1,t表示t时刻管道节点1处的压头，a表示压力波速，

表示t+Δt时刻测量得到的管道上游压头，α表示管道倾斜角，λ表示管道的流体阻力系数，D表示化工管道的管径；其中式4e中，H_d,t+Δt表示t+Δt时刻通过计算得到的管道下游的压头，H_n,t表示t时刻管道第n个节点处的压头，v_n,t表示t时刻管道第n个节点处的流速，

表示t+Δt时刻通过测量得到的管道下游的流速；

(5).泄漏点处的压头和流速关系满足式(4a)(4b)(4c)关系，其中，式(4c)为泄漏位置压头与泄漏量的关系，由式(6)计算得到t+Δt时刻的泄漏量Q_Lc,t＝(v_L-,t-v_L+,t)×A_dL，A_dL为泄漏孔截面积，根据收敛条件定义式(5)可以得到t时刻的相对误差之和Erro_VH；

VH边界条件：

其中，式(5)中，H_0,t+Δt表示t+Δt时刻首端计算的压头，

表示t+Δt时刻首端实际测量的压头，V_n+1,t+Δt表示t+Δt时刻末端计算的流速，

表示t+Δt时刻末端实际测量的流速，Q_Lc,t+Δt表示t+Δt时刻的泄漏量Q_Lc,t+Δt，η_VH为阈值；

式(6)中，H_L表示泄漏点处的高程的压头，Z_L表示泄漏点处的高程，C_d表示孔流系数；

(6).判断式(5)中Erro_VH是否满足收敛条件，即小于阈值η_VH，若满足，则t+Δt时刻反向寻源结束，得到泄漏信息：泄漏位置率k、泄漏孔径的大小d_L和泄漏量Q_Lc,t+Δt；若不满足，则使用混沌粒子群进行优化寻优，混沌粒子群将在k∈[0.5k_layer1,2k_layer1]、m∈[0.5m_layer1,2m_layer1]、d_L∈(0,D],(D为管径)该搜索区域内寻优直到得到满足收敛条件的泄漏位置率、泄漏孔径大小d_L和泄漏量Q_Lc,t+Δt值。

Images