Multi-mode ultrasonic oil-gas-water three-phase flow pattern identification method

Abstract

The invention provides a horizontal pipeline oil-gas-water three-phase flow pattern identification method based on multi-mode ultrasound, and the method comprises the steps: combining a single-crystal pulse wave ultrasonic sensor with a bicrystal continuous wave ultrasonic sensor, and obtaining echo intensity data reflecting phase interface distribution characteristics and Doppler frequency shift data reflecting flow velocity fluctuation characteristics; based on response characteristics of test data under different flow patterns, data fluctuation characteristics are extracted from multiple angles of a time domain, a frequency domain and a space domain, characteristic dimension reduction is carried out by using an equidistant characteristic mapping method, and multi-domain joint characterization and quantification of flow characteristics are completed; a classifier is established based on a nonlinear support vector machine, and accurate identification of the flow pattern of the oil-gas three-phase flow of the horizontal pipeline is realized.

Description

Technical Field

[0001] The present invention belongs to the technical field of multiphase flow measurement, and relates to a method for non-invasively identifying flow patterns of oil-gas-water three-phase flow in a horizontal pipeline by using an ultrasonic sensor with two different working modes, namely pulse wave and continuous wave, through effective feature extraction and dimensionality reduction. Background Art

[0002] Oil, gas and water three-phase flow is a common form of output and raw material in the oil and gas industry. It often appears in related pipelines in the energy and chemical industries, involving production processes such as mining, transportation, storage, and refining. The mixing of three-phase media with different dynamic properties such as density and viscosity causes the phase interface in the oil, gas and water three-phase flow to present a complex non-uniform spatiotemporal distribution structure, which is called the flow pattern. The flow pattern has an important influence on the pressure drop characteristics and mass transfer / heat transfer characteristics of the fluid. It is not only related to the efficiency improvement and safety assurance in the industrial production process, but also has important guiding significance for the accurate measurement of multiphase flow process parameters and the establishment of related mechanism models. Therefore, the accurate identification of the flow pattern of oil, gas and water three-phase flow is a basic issue that cannot be ignored in the field of multiphase flow.

[0003] Accurate identification of flow patterns often relies on high-precision sensing technology, effective feature extraction methods, and appropriate classification algorithms. Among them, high-precision sensing technology is used to obtain process fluctuation signals of the fluid, which are usually closely related to the distribution or fluctuation characteristics of the flow pattern, such as flow rate, content, pressure, etc. At present, many researchers at home and abroad have developed a variety of fluid measurement sensors based on the specific characteristics of the measured object and different sensitive principles, including electrical sensors (probes, wire mesh, planar electrodes, cross-sectional distributed electrodes, etc.), optical sensors (laser Doppler velocimeters, particle image velocimeters, etc.), radiation sensors (X-rays, gamma rays, etc.), microwave sensors (antennas, resonant cavities, etc.), ultrasonic sensors, etc. Different sensing technologies have their own characteristics and application limitations. For example, electrical probes and wire mesh sensors are invasive measurements that interfere with the flow field and bring potential threats of pressure loss and blockage. Optical sensors have strict requirements on the light transmittance of fluids and pipelines, and the test system is complex and costly. Radiation sensors are radioactive and require strict radiation isolation protection when used to ensure the safety of personnel and the environment. Microwave sensors are easily affected by fluid mineralization and flow patterns, and the hardware electronic circuits and sensor design are difficult. Ultrasonic sensors are easily affected by temperature fluctuations. In practical applications, measurement sensors should be reasonably selected according to specific usage scenarios. An effective feature extraction method is to extract objective feature parameters from the test signal based on the sensing technology, according to the response characteristics of the test signal to the flow pattern, and use certain information processing algorithms to effectively characterize the flow state. At present, time domain statistical features (mean, variance, maximum value, minimum value, etc.), frequency domain features (maximum frequency, average frequency, etc.), time-frequency domain features (band correlation coefficient, energy proportion, etc.), and nonlinear features (entropy, fractal dimension, etc.) have been used to characterize flow characteristics. Due to the complexity of the flow process, it is sometimes necessary to combine multiple features to achieve a comprehensive description of the flow characteristics, but the redundancy between features is an issue that needs to be considered. Finally, the extracted flow features are fed into the classifier to identify the flow pattern. Compared with traditional classification methods, machine learning-based methods, such as support vector machines, decision trees, random forests, artificial neural networks, etc., have shown obvious advantages in solving nonlinear and random problems such as multiphase flow pattern identification, and have been successfully applied in many cases.

[0004] In the past few decades, many domestic and foreign researchers have proposed a variety of flow pattern identification methods, but most of them focus on relatively simple gas-water or oil-water two-phase flow. Compared with two-phase flow, oil-gas-water three-phase flow has a more complex phase distribution structure and phase interaction. The accurate identification of its flow pattern needs to rely on richer test information, more effective feature expression and better classification algorithm. At present, multi-sensor fusion strategy has become an effective solution to improve the accuracy of flow pattern identification. However, in practical engineering applications, the increase of sensing modes will greatly increase the complexity, installation difficulty and application cost of the test system, and the most ideal flow pattern identification solution should be accurate and simple. Among many sensing technologies, the ultrasonic method can use different acoustic phenomena such as transmission, reflection, scattering, diffraction, Doppler, etc. emitted when sound waves propagate in multiphase flow to form different measurement modes, which are used to obtain information on various flow parameters such as content, interface position, flow velocity, etc. It has the advantages of rich information acquisition, simple sensor structure, no interference with flow state, convenient installation and maintenance, safe and radiation-free, low cost, and can measure non-transparent fluids in non-transparent pipes. Therefore, the present invention proposes a method for identifying oil-gas-water three-phase flow patterns based on multi-mode ultrasound, which involves an acoustic testing method combining continuous waves and pulse waves, a multi-dimensional feature extraction and dimensionality reduction method, and a support vector machine classification algorithm. Summary of the invention

[0005] The purpose of the present invention is to provide an accurate, non-invasive, low-cost method for identifying the flow pattern of oil, gas and water three-phase flow in a horizontal pipeline by using ultrasonic sensors with two different working modes: pulse wave and continuous wave. The technical solution is as follows:

[0006] A method for identifying flow patterns of oil-gas-water three-phase flow in a horizontal pipeline based on multi-mode ultrasound, which combines a single-crystal pulse wave ultrasonic sensor with a dual-crystal continuous wave ultrasonic sensor, and is used to obtain echo intensity data reflecting the distribution characteristics of the phase interface and Doppler frequency shift data reflecting the flow velocity fluctuation characteristics respectively; based on the response characteristics of the test data under different flow patterns, the data fluctuation characteristics are extracted from multiple angles of the time domain, frequency domain and space domain, and the feature dimension reduction is performed using the equidistant feature mapping method to complete the multi-domain joint characterization and quantification of the flow characteristics; a classifier is constructed based on a nonlinear support vector machine to realize the accurate identification of the flow pattern of the oil-gas three-phase flow in the horizontal pipeline; the method comprises the following steps:

[0007] (1) Using a pulse wave ultrasonic sensor and a continuous wave ultrasonic sensor to synchronously collect flow process information of oil and gas three-phase flow in a horizontal pipeline under different working conditions: the pulse wave ultrasonic sensor includes a piezoelectric ceramic chip installed vertically to the pipe wall, which is used to obtain echo intensity profile data reflecting the distribution characteristics of the fluid phase interface; the continuous wave ultrasonic sensor includes two obliquely installed piezoelectric ceramic chips, which are used to obtain Doppler frequency shift data reflecting the fluid velocity fluctuation characteristics;

[0008] (2) In the echo intensity profile data, the radial position of the pipeline where the maximum echo amplitude value appears is located in the profile obtained within each pulse repetition period, and its probability distribution is statistically analyzed; based on the statistical results, the kurtosis, skewness, weighted average value and the position corresponding to the maximum probability are extracted from the probability distribution diagram, which are used to quantify the flatness, symmetry, fluctuation distribution and extreme value of the probability distribution, respectively, as the spatial domain features of the flow pattern characterization;

[0009] (3) Perform short-time Fourier transform on the Doppler frequency shift data to obtain a time-frequency spectrum, and calculate the average frequency f(τ) of the spectrum at time τ and the corresponding average flow velocity u(τ); based on the calculation results, extract the mean, variance, and average zero-crossing rate from the average flow velocity fluctuation time series, which are used to quantify the central tendency, fluctuation degree, and morphological change intensity of the data, respectively, as time domain features for flow type characterization;

[0010] (4) The Doppler frequency shift data is decomposed into the sum of several intrinsic mode functions and residual terms using the ensemble empirical mode decomposition method; based on the decomposition results, the correlation coefficients between the intrinsic mode functions at each level and the original Doppler frequency shift data are calculated, and the first several intrinsic mode function components with larger correlation coefficients are selected as the principal components representing the significant characteristics of the Doppler frequency shift data, and the energy proportion corresponding to the principal components is used as the multi-scale domain feature of the flow type characterization;

[0011] (5) forming a multidimensional feature vector based on pulse wave and continuous wave ultrasonic test data;

[0012] (6) Taking the multidimensional feature vector as the input data set, reducing its dimension by using the isometric feature mapping method to obtain the feature vector after dimension reduction;

[0013] (7) adding the oil-gas-water three-phase flow pattern label to the feature vector after dimension reduction, and forming a data set with the feature vectors under different working conditions;

[0014] (8) Dividing the data set into a training data set and a test data set according to a certain ratio, the training data set is used to train the flow type classifier model, and the test data set is used to evaluate the classification error of the trained classifier model;

[0015] (9) Obtain a nonlinear support vector machine flow type classifier for outputting the flow type of unknown samples.

[0016] Furthermore, in step (1), the pulse wave ultrasonic sensor intermittently emits ultrasonic pulses into the fluid at a certain repetition frequency, and receives echoes within the pulse interval, and obtains echo intensity profile data A(y,t PRF ), where t PRF ∈(0,nT PRF) is the measurement time, T PRF is the pulse repetition period, n is the number of pulse excitations within the sampling time, y∈(0,D) is the distance from the bottom of the pipe, and D is the pipe diameter.

[0017] Furthermore, in step (1), the continuous wave ultrasonic sensor comprises two piezoelectric ceramic chips installed at an angle of θ (45°≤θ≤60°) with the horizontal direction, which are responsible for continuously transmitting and receiving ultrasonic waves into the fluid during operation; and the Doppler frequency shift data f reflecting the velocity fluctuation characteristics of the fluid are obtained by frequency demodulating the received ultrasonic waves. d (t), where t∈(0,T) is time and T is the sampling time.

[0018] Furthermore, the method of step (2) is:

[0019] In the echo intensity profile data A(y,t PRF ),t PRF ∈(0,nT PRF ), locate the position where the maximum echo amplitude occurs in the profile obtained within each pulse repetition period, and count its probability distribution:

[0020]

[0021] Where P(y) represents the maximum echo amplitude A within the sampling time. max The probability of appearing in the position interval [y,y+Δy], Num(y≤A max ≤y+Δy) represents the maximum value of the echo amplitude within the sampling time A max The number of occurrences in the position interval [y,y+Δy];

[0022] Based on the statistical results, the kurtosis P is extracted from the probability distribution diagram. kur , skewness P ske , weighted average P wei The position P corresponding to the maximum probability lo , which are used to quantify the flatness, symmetry, fluctuation distribution and extreme value of probability distribution as the spatial domain characteristics of flow pattern characterization:

[0023]

[0024] Among them, E(*) represents the average operation, P max is the maximum probability value in the probability distribution, P mean and P var are the mean and variance of the probability distribution statistics, respectively.

[0025] Furthermore, the method of step (3) is:

[0026] For the Doppler frequency shift data f d (t),t∈(0,T) is short-time Fourier transformed to obtain the time-frequency spectrum S d (τ,f), and for S d Each column of (τ,f) uses the following formula to calculate the average frequency of the spectrum at time τ And the corresponding average flow velocity of the fluid

[0027]

[0028] Where c is the speed of sound in the sound wedge, f0 is the frequency of the emitted sound wave, and θ is the angle between the central axis of the piezoelectric ceramic chip and the horizontal direction;

[0029] Based on the calculation results, in the average flow velocity fluctuation time series Extract the mean variance Average zero crossing rate They are used to quantify the central tendency, fluctuation degree and drastic degree of morphological change of the data as the time domain characteristics of flow pattern characterization:

[0030]

[0031] Among them, sgn(*) is the sign function, is the data after subtracting its own mean from the average flow velocity fluctuation time series, and M is the total number of data points in the average flow velocity fluctuation time series.

[0032] Furthermore, the method of step (4) is:

[0033] The Doppler frequency shift data f is transformed into d (t),t∈(0,T) is decomposed into several intrinsic mode functions IMF j (t)(j=1,2,...,N) and the residual term r N (t) and:

[0034]

[0035] Where N is the number of intrinsic mode functions obtained by decomposition;

[0036] Based on the decomposition results, the correlation coefficients between the intrinsic mode functions at each level and the original Doppler frequency shift data are calculated, and the first b intrinsic mode function components IMF with the largest correlation coefficient are selected. i ,i∈{j1,j2,j3,...,j b} as the principal component that characterizes the significant characteristics of the Doppler frequency shift data, and the energy proportion R corresponding to the principal component i Multi-scale domain features as flow pattern representation:

[0037]

[0038] Furthermore, in step (5), the multidimensional feature vector S is:

[0039]

[0040] Furthermore, in step (6), the multidimensional feature vector S is used as the input data set X = {x i |x1,x2,···,x b+7}, where x i Represent each feature parameter in the multidimensional feature vector S, use the isometric feature mapping method to reduce its dimension to m-dimensional space (m < b + 7), and obtain the feature vector after dimension reduction, which is achieved by the following sub-steps:

[0041] (6.1) Use the K nearest neighbor search algorithm to determine each data point x in the data set i The k nearest neighbors of i ={x i1 ,x i2 ,···,x ik},i∈(1,b+7);

[0042] (6.2) Construct a neighborhood connection graph G = [V, E], where V represents the nodes in the neighborhood connection graph, which is all the data points in the data set X, and E represents the edges in the neighborhood connection graph, whose length is the distance d between each data point in the data set X. G Measure: When data x i With x j When they belong to the same neighborhood space, d G Let the Euclidean distance between two points be d G =d2(x i ,x j )=||x i -x j ||2, when data x i With x j When they do not belong to the same neighborhood space, d G =∞, indicating no connection;

[0043] (6.3) Use the shortest path algorithm to determine the shortest distance between any two points in the neighborhood connection graph G and obtain the geodesic matrix: D G ={dist(x i ,x j )},x i ,x j ∈X, where

[0044]

[0045] (6.4) Construct low-dimensional embedding: Use the multidimensional scaling (MDS) algorithm to construct the kernel matrix Where S G =[dist(x i ,x j )] 2 , H = IJ / (b+7), I is the unit diagonal matrix, J is the matrix of all 1s, and τ(D G ) performs singular value decomposition and takes the first m eigenvalues ​​λ1≤λ2≤...≤λ m and the corresponding eigenvectors v1,v2,...,v m , then the new coordinates of the original data set X formed by the multidimensional feature vector S in the low-dimensional m-dimensional space are expressed as:

[0046]

[0047] Among them, Y is the eigenvector after dimensionality reduction, Λ m =diag(λ1,λ2,...,λ m ), V m =[v1,v2,…,v m ]

[0048] Furthermore, the shortest path algorithm is Dijkstra's or Floyds' algorithm.

[0049] The essential features of the present invention are: combining two ultrasonic sensors with different working modes, pulse wave and continuous wave, to realize non-invasive synchronous measurement of oil-gas-water three-phase flow, respectively obtaining echo intensity profile data reflecting the distribution characteristics of fluid phase interface and Doppler frequency shift data reflecting the fluid velocity fluctuation characteristics; deeply mining ultrasonic test data, extracting characteristic parameters from multiple angles in spatial domain, time domain and multi-scale domain to realize objective and comprehensive characterization of flow process; realizing comprehensive dimensionality reduction of nonlinear global characteristic relationship in high-dimensional feature space through equidistant feature mapping; constructing flow type classifier by nonlinear support vector machine method suitable for small sample size, realizing accurate identification of flow type of oil-gas-water three-phase flow. The beneficial effects and advantages of the present invention are as follows:

[0050] 1. This method adopts a multi-mode ultrasonic combined testing method, which has the advantages of non-invasiveness, simple structure and low cost;

[0051] 2. This method extracts data features from the spatial domain, temporal domain, and multi-scale domain, and can achieve a comprehensive and objective characterization of the flow process;

[0052] 3. This method uses the isometric feature mapping method to achieve feature dimensionality reduction, which helps to reduce the computational complexity, improve the generalization of the model, and reduce overfitting in the process of building a flow type classifier. BRIEF DESCRIPTION OF THE DRAWINGS

[0053] The following drawings and tables describe selected embodiments of the present invention, which are exemplary and not exhaustive or limiting, wherein:

[0054] Figure 1 Schematic diagram of pulse wave and continuous wave ultrasonic testing system in the measurement method of the present invention;

[0055] Figure 2 Schematic diagram of the structure of pulse wave and continuous wave ultrasonic sensors in the measurement method of the present invention;

[0056] Figure 3 A flow chart of the short-time Fourier transform method in the measurement method of the present invention;

[0057] Figure 4 A flow chart of the ensemble empirical mode decomposition method in the measurement method of the present invention;

[0058] Figure 5 Algorithm flow chart of nonlinear support vector machine in the measurement method of the present invention;

[0059] Figure 6 This is a flow chart for identifying the flow pattern of the oil-gas-water three-phase flow in the measurement method of the present invention. DETAILED DESCRIPTION

[0060] The embodiments of the present invention are described in detail below with reference to the accompanying drawings and embodiments.

[0061] Figure 1 The figure is a schematic diagram of the pulse wave and continuous wave ultrasonic test system in the measurement method of the present invention. The test system is composed of three modules: an excitation unit I, a sensing unit II and an acquisition unit III. The sensing unit II includes a single-crystal pulse wave ultrasonic sensor 0 and a dual-crystal continuous wave ultrasonic sensor 1, which are installed at a certain interval at the bottom of the pipeline and work simultaneously, and the installation order is not limited; the excitation unit I includes a signal generator 2 and a power amplifier 3, which are used to generate continuous wave and pulse wave excitation signals required for the operation of the ultrasonic sensor; the acquisition unit III includes a data processing card 4, a data acquisition card 5 and a computer 6, wherein the data processing card 4 is used to process the sound wave signal modulated by the fluid received by the ultrasonic sensor and output the echo intensity data related to the distribution characteristics of the phase interface and the Doppler frequency shift data related to the flow velocity fluctuation characteristics, and the data acquisition card 5 synchronously collects the output data of the data processing card 4 and uploads it to the computer 6 to complete the real-time storage of the data.

[0062] Figure 2The figure is a schematic diagram of the structure of the pulse wave and continuous wave ultrasonic sensors in the measurement method of the present invention. The pulse wave ultrasonic sensor 0 and the continuous wave ultrasonic sensor 1 used in the present invention are installed at the bottom of the pipe 7 at a certain interval to ensure that the measurement spaces 8 of the two sensors do not interfere with each other. The pulse wave ultrasonic sensor 0 contains a piezoelectric ceramic chip 901 with a center frequency of 1MHz attached to the acoustic wedge 101. When the piezoelectric ceramic chip is working, it will intermittently emit ultrasonic pulses into the fluid at a certain repetition frequency and receive the echo within the pulse interval. The amplitude information of the received sound wave is related to the distribution characteristics of the phase interface of the measured fluid. The continuous wave ultrasonic sensor 1 contains two piezoelectric ceramic chips 902 and 903 with a center frequency of 1 MHz, attached to the acoustic wedges 102 and 103. The two piezoelectric ceramic chips are responsible for continuously transmitting and receiving ultrasonic waves into the fluid when working. The frequency offset between the transmitted and received sound waves is related to the speed of the measured fluid. The acoustic wedges 102 and 103 are cut into a fixed shape to ensure that the angle between the normal of the piezoelectric ceramic chip and the horizontal direction is θ (45°≤θ≤60°), and the beam path direction is opposite to the flow direction of the fluid 12. In addition, the back of the piezoelectric ceramic chip in the pulse wave ultrasonic sensor 0 and the continuous wave ultrasonic sensor 1 is filled with sound absorbing material 11 to effectively reduce the reflection and scattering of sound waves inside the sensor and improve the measurement signal-to-noise ratio and resolution.

[0063] Figure 3 The figure is a flow chart of the short-time Fourier transform method in the measurement method of the present invention. Short-time Fourier transform (STFT) is an important time-frequency analysis method. Its idea is to use a shorter window function to slide along the time axis and calculate the frequency characteristics of the signal at different times in different windows, which can effectively capture the transient changes and time-varying characteristics of the signal. Mathematically, the STFT of the time domain signal s(t) is defined as: Where h(*) is a window function, which is used to limit the time range of analysis. Commonly used window functions include rectangular window, Hamming window, Hanning window, Blackman window, etc. The implementation of STFT for a given time domain signal s(t) mainly includes the following steps: 1) Setting the window function: Select the window function and its length according to the application requirements, and set the frame shift. The frame shift determines the sliding step of the window function on the time axis. The common frame shift is 50% or 75% of the window function length; 2) Framing: Use the window function to divide the signal into multiple overlapping short time frames, and each frame is multiplied by the window function to obtain a windowed signal segment. The length of each signal segment is the length of the window function; 3) Perform Fourier transform on each signal segment to obtain its frequency domain representation; 4) Arrange the Fourier transform results of all signal segments in chronological order to form a time-frequency matrix S(τ,f).

[0064] Figure 4 The flowchart of the ensemble empirical mode decomposition method in the measurement method of the present invention. Ensemble Empirical Mode Decomposition (EEMD) is a data-adaptive multi-resolution technology that can decompose the fluctuation components of different frequencies contained in the signal into intrinsic mode functions of different scales. EEMD is proposed on the basis of the traditional empirical mode decomposition (EMD) method. It utilizes the frequency uniform distribution characteristics and zero mean characteristics of white noise. By adding multiple groups of white noise to the signal, it can overcome the modal aliasing problem that occurs in the EMD method and improve the accuracy of decomposition. Therefore, the EEMD method is also called a white noise auxiliary analysis method. For a given signal x(t), EEMD decomposition is performed, which mainly includes the following steps: 1) Adding Gaussian white noise ω with a mean of 0 and a standard deviation of 0.1 to 0.4 times the standard deviation of x(t) to x(t); i (t), construct a new time series X i (t) = x(t) + ω i (t), i = 1, 2…g; 2) for X i (t) is decomposed by EMD to obtain IMF ij is the jth intrinsic mode function of the ith time series, r in is the decomposition residual term of the i-th time series; 3) the g-time decomposition results are averaged to remove the influence of white noise, and the intrinsic mode functions IMF of each level that can characterize the fluctuation characteristics of different frequencies in the signal are obtained. j and the residual term r that characterizes the trend characteristics of the signal N : Finally, x(t) can be expressed as:

[0065] Figure 5 The flowchart is an algorithm of nonlinear support vector machine in the measurement method of the present invention. Support vector machine (SVM) is a semi-supervised machine learning method based on statistical learning theory and the principle of structural risk minimization. It is widely used to solve small sample, nonlinear and high-dimensional pattern recognition and regression problems. For nonlinear separable classification problems, SVM is essentially a binary classifier. Its basic idea is to use kernel function to map data from low-dimensional linear inseparable feature space to high-dimensional Hilbert space, and find a separating hyperplane in high-dimensional space that can separate samples of different categories in the data set and has the largest geometric interval. It can be achieved through the optimization algorithm of convex quadratic programming. The specific steps are as follows: 1) Input training data set X = {(x1, y1), (x2, y2),…, (x n,y n )}, where x i ∈R D , i = 1, 2, ..., n is the n samples in the D-dimensional feature space, y i ∈{-1,+1} is the sample category; 2) Select an appropriate kernel function L(*) and penalty factor C, construct and solve the convex quadratic programming problem: Get the optimal solution Select α * A positive component of Calculation parameters Then the separating hyperplane can be expressed as Use the separating hyperplane to construct the classification decision function: Where sgn(*) represents the sign function. For samples of unknown category, the output of the classification decision function can be used to determine its category.

[0066] Figure 6 This is a flow chart for identifying the flow pattern of the oil-gas-water three-phase flow in the measurement method of the present invention.

[0067] The oil-gas-water three-phase flow pattern identification method proposed in the present invention involves an acoustic test method combining continuous waves and pulse waves, a multidimensional feature extraction and dimensionality reduction method, and a support vector machine classification algorithm. The specific implementation steps are as follows:

[0068] Step 1: Exploitation Figure 1 , Figure 2 The ultrasonic testing system and the sensor shown synchronously obtain pulse wave and continuous wave ultrasonic measurement data under different flow patterns of oil, gas and water three-phase flow in a horizontal pipeline.

[0069] 1) The pulse wave ultrasonic sensor works in the self-transmitting and self-receiving mode, that is, the piezoelectric ceramic chip is measured at a repetition frequency f PRF Ultrasonic pulses are intermittently emitted into the fluid, and the echoes are received by the same piezoelectric ceramic chip during the pulse interval. By demodulating the amplitude of the received echoes, the echo intensity profile data A(y,t PRF ), where t PRF ∈(0,nT PRF ) is the measurement time, T PRF =1 / f PRF is the pulse repetition period, n is the number of pulse excitations within the sampling time, y∈(0,D) is the distance from the bottom of the pipe, and D is the pipe diameter.

[0070] 2) The continuous wave ultrasonic sensor works in a transmit-receive mode, that is, during measurement, one piezoelectric ceramic chip continuously transmits ultrasonic waves into the fluid, while the other piezoelectric ceramic chip continuously receives the echo. According to the Doppler effect, by frequency demodulating the received echo, the Doppler frequency shift data reflecting the velocity fluctuation characteristics of the fluid can be obtained. d (t), where t∈(0,T) is time and T is the sampling time.

[0071] Step 2: Process the echo intensity data and Doppler frequency shift data, and perform feature extraction from multiple angles to achieve flow type characterization.

[0072] 1) In the echo intensity profile data A(y,t PRF ),t PRF ∈(0,nT PRF ), locate the position where the maximum echo amplitude occurs in each pulse repetition period, and count its probability distribution:

[0073]

[0074] Among them, P(y) represents the probability that the maximum echo amplitude value within the sampling time appears in the position interval [y, y+Δy], Num(y≤A max ≤y+Δy) indicates the number of times the maximum echo amplitude appears in the position interval [y,y+Δy] within the sampling time.

[0075] The shape of the probability distribution graph is modulated by the spatial distribution characteristics of the fluid phase interface at the radial position of the pipeline. Therefore, the kurtosis P is extracted from the probability distribution graph. kur , skewness P ske , weighted average P wei The position P corresponding to the maximum probability lo , which are used to quantify the flatness, symmetry, fluctuation distribution and extreme value of probability distribution as the spatial domain characteristics of flow pattern characterization:

[0076]

[0077] Among them, E(*) represents the average operation, P mean and P var The mean and variance of the probability distribution statistics respectively.

[0078] 2) Utilize Figure 3 The short-time Fourier transform algorithm process shown in the figure processes the Doppler frequency shift data f d (t),t∈(0,T), obtain the time-frequency spectrum S d (τ,f). S d(τ,f) is a two-dimensional matrix, and each column vector represents the spectrum obtained by Fourier transform after windowing the signal at different time τ. d Each column of (τ,f) calculates the average frequency of the spectrum using the following formula:

[0079]

[0080] Furthermore, according to the Doppler effect, the average flow velocity u of the fluid at time τ is calculated using the following formula:

[0081]

[0082] Among them, c0 is the speed of sound in the fluid, f0 is the frequency of the emitted sound wave, and θ0 is the angle between the sound wave in the fluid and the horizontal direction.

[0083] Based on the calculation results, extract the average flow velocity fluctuation time series The mean variance Average zero crossing rate They are used to quantify the central tendency, fluctuation degree and drastic degree of morphological change of the data as the time domain characteristics of flow pattern characterization:

[0084]

[0085] Among them, sgn(*) is the sign function, is the data after subtracting its own mean from the average flow velocity fluctuation time series, and M is the total number of data points in the average flow velocity fluctuation time series.

[0086] 3) Utilize Figure 4 The algorithm flow of the ensemble empirical mode decomposition method shown in FIG. 1 is to convert the Doppler frequency shift data f d (t),t∈(0,T) is decomposed into several intrinsic mode functions IMF j (t)(j=1,2,...,N) and the residual term r N (t) and:

[0087]

[0088] Where N is the number of intrinsic mode functions obtained by decomposition. The intrinsic mode functions of each level obtained by ensemble empirical mode decomposition reflect the different frequency band components in the original signal, and the residual term represents the trend information in the original signal.

[0089] Based on the decomposition results, calculate the energy proportion R of each level of eigenmode function j :

[0090]

[0091] Furthermore, the correlation coefficient CC between each level of intrinsic mode function and the original Doppler frequency shift data is calculated using the following formula: IMFj :

[0092]

[0093] Where E[*] represents the averaging operation. Generally, the larger the correlation coefficient, the more the corresponding intrinsic mode function can reflect the most significant characteristics of the original signal.

[0094] Sort the eigenmode function components according to the size of the correlation coefficient The first 10 intrinsic mode function components with larger correlation coefficients are selected as the main components to characterize the significant characteristics of Doppler frequency shift data, and their corresponding energy proportions are used as multi-scale domain features for flow type characterization:

[0095] R j ,j={j1,j2,j3,...,j b} (14)

[0096] 4) Form a multi-dimensional feature vector S based on pulse wave and continuous wave ultrasonic test data:

[0097]

[0098] The dimension of the multidimensional feature vector S is 17 dimensions.

[0099] Step 3: Use the isometric feature mapping method to reduce the dimension of the multi-dimensional feature vector S.

[0100] Isometric Feature Mapping (Isomap) is a popular learning algorithm that can map high-dimensional data points to low-dimensional space to achieve data dimensionality reduction by constructing an adjacency graph and calculating the shortest path to approximate the geodesic distance matrix on the manifold. It has the advantages of high computational efficiency, global optimization, and good asymptotic convergence, and is widely used in comprehensive dimensionality reduction of nonlinear global feature relationships in high-dimensional space. The multidimensional feature vector S is used as the input data set X = {x i |x1,x2,···,x 17}, where x i Represent each feature parameter in the multidimensional feature vector and use the Isomap method to reduce its dimension to m-dimensional space (m < 17), which mainly includes the following sub-steps:

[0101] 1) Use the K nearest neighbor search algorithm to determine the nearest neighbor of each data point x in the data set. i The k nearest neighbors of i ={x i1 ,x i2 ,···,x ik},i∈(1,17);

[0102] 2) Construct a neighborhood connection graph G = [V, E], where V represents the nodes in the neighborhood connection graph, which is all the data points in the data set X, and E represents the edges in the neighborhood connection graph, whose length is the distance d between each data point in the data set X. G Measure: When data x i With x j When they belong to the same neighborhood space, d G Let the Euclidean distance between two points be d G =d2(x i ,x j )=||x i -x j ||2, when data x i With x j When they do not belong to the same neighborhood space, d G =∞, indicating no connection;

[0103] 3) Use the shortest path algorithm (Dijkstra or Floyds algorithm) to determine the shortest distance between any two points in the neighborhood connection graph G and obtain the geodesic matrix: D G ={dist(x i ,x j )},x i ,x j ∈X, where

[0104] 4) Construct low-dimensional embedding: Use the multidimensional scaling (MDS) algorithm to construct the kernel matrix in H = IJ / 17, I is the unit diagonal matrix, J is the matrix of all 1s, and τ(D G ) performs singular value decomposition and takes the first m eigenvalues ​​λ1≤λ2≤...≤λ m and the corresponding eigenvectors v1,v2,...,v m , then the new coordinates of the original data set X formed by the multidimensional feature vector S in the low-dimensional m-dimensional space can be expressed as:

[0105]

[0106] Among them, Y is the eigenvector after dimensionality reduction, Λ m =diag(λ1,λ2,…,λ m ), V m =[v1,v2,...,v m ].

[0107] Step 4: Using the nonlinear support vector machine method, the flow pattern classifier is constructed using the feature vector after dimension reduction to realize the flow pattern identification of oil-gas-water three-phase flow, which includes the following sub-steps:

[0108] 1) Add labels to the data set: According to the actual observation results, flow type labels are added to the feature vector data set after dimension reduction to form a data set W. In the embodiment of the present invention, a total of 610 ultrasonic measurement data of oil-gas-water three-phase flow under working conditions are collected, including a total of 8 flow type labels, namely water-based-oil-water separation-laminar flow, water-based-oil-water dispersion-laminar flow, water-based-oil-water separation-wavy flow, water-based-oil-water dispersion-wavy flow, water-based-oil-water dispersion-plug flow, water-based-oil-water dispersion-elastic flow, water-based-oil-water dispersion-through-elastic flow and water-based-oil-water dispersion-annular flow.

[0109] 2) Use the holdout method to divide the data set: divide the data set into two subsets Q and Z in a ratio of 7:3, satisfying Q∪Z=W and The subset Q will be used as a training data set for training the flow type classifier model, and the subset Z will be used as a test data set to evaluate the classification error of the trained classifier model.

[0110] 3) Constructing and testing a flow pattern classifier based on a nonlinear support vector machine method: The classic SVM algorithm is designed to solve binary classification problems. When solving multi-classification problems, a "one-to-one" or "one-to-many" method can be used to construct a multi-classifier. In the embodiment of the present invention, a "one-to-one" strategy is used to construct a multi-classifier for oil-gas-water three-phase flow pattern recognition. Specifically, during training, the data sets of any two types of flow pattern samples are trained using Figure 5 The algorithm flow shown in the figure trains a decision function, and a total of (8×7) / 2=28 decision functions can be obtained under 8 flow types, corresponding to 28 SVM binary classifiers, each of which only classifies between two types of flow types. When predicting the flow type of the test sample, the prediction results of the 28 SVM binary classifiers are voted and counted, and the result with the highest number of votes is the flow type prediction result of the test sample.

[0111] Table 1 Flow pattern identification results of oil-gas-water three-phase flow

[0112]

[0113]

[0114] Table 1 shows the recognition results of the oil-gas-water three-phase flow patterns in the measurement method of the present invention. It can be seen from the table that the recognition rates of eight flow patterns, namely, water-based-oil-water separation-laminar flow, water-based-oil-water dispersion-laminar flow, water-based-oil-water separation-wavy flow, water-based-oil-water dispersion-wavy flow, water-based-oil-water dispersion-plug flow, water-based-oil-water dispersion-elastic flow, water-based-oil-water dispersion-through-elastic flow and water-based-oil-water dispersion-annular flow, are 93.33%, 100%, 88.89%, 94.11%, 90.00%, 93.75%, 93.33% and 100% respectively, and the average recognition rate of all flow patterns is 93.44%, indicating that the proposed method is effective.

[0115] Finally, it should be noted that the above embodiments are only used to illustrate the technical solutions of the present invention rather than to limit it. Although the present invention has been described in detail with reference to the above embodiments, ordinary technicians in the relevant field should understand that the specific implementation methods of the present invention can still be modified or replaced by equivalents. Any modification or equivalent replacement that does not depart from the spirit and scope of the present invention should be covered within the scope of protection of the claims of the present invention.

Claims

1. A method for identifying flow patterns of oil, gas and water three-phase flow in a horizontal pipeline based on multi-mode ultrasound, which combines a single-crystal pulse wave ultrasonic sensor with a dual-crystal continuous wave ultrasonic sensor to obtain echo intensity data reflecting the distribution characteristics of the phase interface and Doppler frequency shift data reflecting the flow velocity fluctuation characteristics respectively; based on the response characteristics of the test data under different flow patterns, the data fluctuation characteristics are extracted from multiple angles of the time domain, frequency domain and space domain, and the feature dimension reduction is performed using the equidistant feature mapping method to complete the multi-domain joint characterization and quantification of the flow characteristics; a classifier is constructed based on a nonlinear support vector machine to realize the accurate identification of the flow pattern of oil and gas three-phase flow in the horizontal pipeline; the method comprises the following steps: (1) Using a pulse wave ultrasonic sensor and a continuous wave ultrasonic sensor to synchronously collect flow process information of oil and gas three-phase flow in a horizontal pipeline under different working conditions: the pulse wave ultrasonic sensor includes a piezoelectric ceramic chip installed vertically to the pipe wall, which is used to obtain echo intensity profile data reflecting the distribution characteristics of the fluid phase interface; the continuous wave ultrasonic sensor includes two obliquely installed piezoelectric ceramic chips, which are used to obtain Doppler frequency shift data reflecting the fluid velocity fluctuation characteristics; (2) In the echo intensity profile data, the radial position of the pipeline where the maximum echo amplitude value appears is located in the profile obtained within each pulse repetition period, and its probability distribution is statistically analyzed; based on the statistical results, the kurtosis, skewness, weighted average value and the position corresponding to the maximum probability are extracted from the probability distribution diagram, which are used to quantify the flatness, symmetry, fluctuation distribution and extreme value of the probability distribution, respectively, as the spatial domain features of the flow pattern characterization; (3) Perform short-time Fourier transform on the Doppler frequency shift data to obtain a time-frequency spectrum, and calculate the average frequency f(τ) of the spectrum at time τ and the corresponding average flow velocity u(τ); based on the calculation results, extract the mean, variance, and average zero-crossing rate from the average flow velocity fluctuation time series, which are used to quantify the central tendency, fluctuation degree, and morphological change intensity of the data, respectively, as time domain features for flow type characterization; (4) The Doppler frequency shift data is decomposed into the sum of several intrinsic mode functions and residual terms using the ensemble empirical mode decomposition method; based on the decomposition results, the correlation coefficients between the intrinsic mode functions at each level and the original Doppler frequency shift data are calculated, and the first several intrinsic mode function components with larger correlation coefficients are selected as the principal components representing the significant characteristics of the Doppler frequency shift data, and the energy proportion corresponding to the principal components is used as the multi-scale domain feature of the flow type characterization; (5) forming a multidimensional feature vector based on pulse wave and continuous wave ultrasonic test data; (6) Taking the multidimensional feature vector as the input data set, reducing its dimension by using the isometric feature mapping method to obtain the feature vector after dimension reduction; (7) adding the oil-gas-water three-phase flow pattern label to the feature vector after dimension reduction, and forming a data set with the feature vectors under different working conditions; (8) Dividing the data set into a training data set and a test data set according to a certain ratio, the training data set is used to train the flow type classifier model, and the test data set is used to evaluate the classification error of the trained classifier model; (9) Obtain a nonlinear support vector machine flow type classifier for outputting the flow type of unknown samples.

2. The method for identifying flow patterns of oil, gas and water three-phase flow in a horizontal pipeline according to claim 1 is characterized in that: In step (1), the pulse wave ultrasonic sensor intermittently emits ultrasonic pulses into the fluid at a certain repetition frequency, and receives echoes within the pulse interval, and obtains echo intensity profile data A(y,t PRF ), where t PRF ∈(0,nT PRF ) is the measurement time, T PRF is the pulse repetition period, n is the number of pulse excitations within the sampling time, y∈(0,D) is the distance from the bottom of the pipe, and D is the pipe diameter.

3. The method for identifying flow patterns of oil, gas and water three-phase flow in a horizontal pipeline according to claim 2 is characterized in that: In step (1), the continuous wave ultrasonic sensor comprises two piezoelectric ceramic chips installed at an angle of θ (45°≤θ≤60°) to the horizontal direction, which are responsible for continuously transmitting and receiving ultrasonic waves into the fluid during operation; and the Doppler frequency shift data f reflecting the velocity fluctuation characteristics of the fluid are obtained by frequency demodulating the received waves. d (t), where t∈(0,T) is time and T is the sampling time.

4. The method for identifying flow patterns of oil, gas and water three-phase flow in a horizontal pipeline according to claim 3 is characterized in that: The method of step (2) is: In the echo intensity profile data A(y,t PRF ),t PRF ∈(0,nT PRF ), locate the position where the maximum echo amplitude occurs in the profile obtained within each pulse repetition period, and count its probability distribution: Where P(y) represents the maximum echo amplitude A within the sampling time. max The probability of appearing in the position interval [y,y+Δy], Num(y≤A max ≤y+Δy) represents the maximum value of the echo amplitude within the sampling time A max The number of occurrences in the position interval [y,y+Δy]; Based on the statistical results, the kurtosis P is extracted from the probability distribution diagram. kur , skewness P ske , weighted average P wei The position P corresponding to the maximum probability lo , which are used to quantify the flatness, symmetry, fluctuation distribution and extreme value of probability distribution as the spatial domain characteristics of flow pattern characterization: P lo =y(when P=P max ) Among them, E(*) represents the average operation, P max is the maximum probability value in the probability distribution, P mean and P var are the mean and variance of the probability distribution statistics, respectively.

5. The method for identifying flow patterns of oil, gas and water three-phase flow in a horizontal pipeline according to claim 4, characterized in that: The method of step (3) is: For the Doppler frequency shift data f d (t),t∈(0,T) is short-time Fourier transformed to obtain the time-frequency spectrum S d (τ,f), and for S d Each column of (τ,f) uses the following formula to calculate the average frequency of the spectrum at time τ And the corresponding average flow velocity of the fluid Where c is the speed of sound in the sound wedge, f0 is the frequency of the emitted sound wave, and θ is the angle between the central axis of the piezoelectric ceramic chip and the horizontal direction; Based on the calculation results, in the average flow velocity fluctuation time series Extract the mean variance Average zero crossing rate They are used to quantify the central tendency, fluctuation degree and drastic degree of morphological change of the data as the time domain characteristics of flow pattern characterization: Among them, sgn(*) is the sign function, is the data after subtracting its own mean from the average flow velocity fluctuation time series, and M is the total number of data points in the average flow velocity fluctuation time series.

6. The method for identifying flow patterns of oil, gas and water three-phase flow in a horizontal pipeline according to claim 5, characterized in that: The method of step (4) is: The Doppler frequency shift data f is transformed into d (t),t∈(0,T) is decomposed into several intrinsic mode functions IMF j (t)(j=1,2,...,N) and the residual term r N (t) and: Where N is the number of intrinsic mode functions obtained by decomposition; Based on the decomposition results, the correlation coefficients between the intrinsic mode functions at each level and the original Doppler frequency shift data are calculated, and the first b intrinsic mode function components IMF with the largest correlation coefficient are selected. i ,i∈{j1,j2,j3,…,j b } as the principal component that characterizes the significant characteristics of the Doppler frequency shift data, and the energy proportion R corresponding to the principal component i Multi-scale domain features as flow pattern representation:

7. The method for identifying flow patterns of oil, gas and water three-phase flow in a horizontal pipeline according to claim 6, characterized in that: In step (5), the multidimensional feature vector S is:

8. The method for identifying flow patterns of oil, gas and water three-phase flow in a horizontal pipeline according to claim 7, characterized in that: In step (6), the multidimensional feature vector S is used as the input data set X = {x i |x1,x2,···,x b+7 }, where x i Represent each feature parameter in the multidimensional feature vector S, use the isometric feature mapping method to reduce its dimension to m-dimensional space (m < b + 7), and obtain the feature vector after dimension reduction, which is achieved by the following sub-steps: (6.1) Use the K nearest neighbor search algorithm to determine each data point x in the data set i The k nearest neighbors of i ={x i1 ,x i2 ,···,x ik },i∈(1,b+7); (6.2) Construct a neighborhood connection graph G = [V, E], where V represents the nodes in the neighborhood connection graph, which is all the data points in the data set X, and E represents the edges in the neighborhood connection graph, whose length is the distance d between each data point in the data set X. G Measure: When data x i With x j When they belong to the same neighborhood space, d G Let the Euclidean distance between two points be d G =d2(x i ,x j )=||x i -x j ||2, when data x i With x j When they do not belong to the same neighborhood space, d G =∞, indicating no connection; (6.3) Use the shortest path algorithm to determine the shortest distance between any two points in the neighborhood connection graph G and obtain the geodesic matrix: D G ={dist(x i ,x j )},x i ,x j ∈X, where (6.4) Construct low-dimensional embedding: Use the multidimensional scaling (MDS) algorithm to construct the kernel matrix Where S G =[dist(x i ,x j )] 2 , H = IJ / (b+7), I is the unit diagonal matrix, J is the matrix of all 1s, and τ(D G ) performs singular value decomposition and takes the first m eigenvalues ​​λ1≤λ2≤…≤λ m and the corresponding eigenvectors v1,v2,...,v m , then the new coordinates of the original data set X formed by the multidimensional feature vector S in the low-dimensional m-dimensional space are expressed as: Among them, Y is the eigenvector after dimensionality reduction, Λ m =diag(λ1,λ2,...,λ m ), V m =[v1,v2,...,v m ].

9. The method for identifying flow patterns of oil, gas and water three-phase flow in a horizontal pipeline according to claim 8, characterized in that: The shortest path algorithm is Dijkstra or Floyds algorithm.

Images