Introduction

Oil–water two-phase flow is a major two-phase flow (Gao et al. 2018), which is widely present in petroleum transportation pipelines. Its flow pattern is complex and changeable, and its dynamic evolution mechanism is extremely complicated. Therefore, it is very difficult to recognize the flow patterns of oil–water two-phase flow accurately. A large number of scientists have been concerned to identify the flow patterns of oil–water two-phase flow quickly and accurately.

Many scholars have been focusing on the identification of two-phase flow patterns. In the early stage, scholars mainly studied the flow patterns of oil–water two-phase flow by observation (Trallero et al. 1997; Nadler and Mewes 1997; Angeli and Hewitt 2000). The main methods included high-speed camera observation (Tan et al. 2018), micro-probe detection (Flores et al. 1997) and PIV technology (Wang et al. 2011). In recent years, many researchers have focused on the identification of flow patterns of two-phase flows by using indirect methods. These methods can not only extract dynamic characteristics of various flow patterns from experimental fluctuation signals, but also identify different flow patterns. The fluctuation signal is commonly used to identify various flow patterns that can reflect the pressure or conductance fluctuation of the mixed fluid. The time–frequency characteristics (Gao et al. 2016a, b) are used to describe the motion behavior of the flow pattern. In order to reflect the multi-scale dynamic characteristics of the flow pattern, Hilbert-Huang transform (Ding et al. 2007) and wavelet transform method (Gao et al. 2017) are performed to extract the characteristics of the flow pattern. Recently, complex network features (Gao et al. 2016a, b, 2018) have been proven to provide effective solutions for flow pattern identification. Many researchers have applied artificial neural networks to the process of flow pattern recognition and have achieved great results. A method for identifying gas–solid two-phase flow patterns in horizontal pneumatic conveying pipelines based on electrostatic sensor array (ESA) and artificial neural network (ANN) was proposed by Fu and Li (2018). The experimental results show that their proposed method can identify fully suspended flow, stratified flow, dune flow, and slug flow effectively. Electrical capacitance tomography (ECT) and neural network are proposed to identify two-phase flow patterns (Roman et al. 2016). This method can achieve the accuracy of liquid–vapor flow pattern recognition up to 98.1%. Fuzzy logic with principal component analysis (PCA) and support vector machine (SVM) are applied to improve the classification accuracy of gas–liquid two-phase flow regimes (Shanthi and Pappa 2017). During the experiment, they proved SVM with features reduced using PCA gives the better classification accuracy and computationally less intensive. A method combining artificial neural networks with related dimensions to construct novel gas–liquid flow pattern diagrams to distinguish between the bubble, bubble/plug transitional, plug, slug, and annular flows (Huang et al. 2017). A method based on combination of multi-beam gamma ray attenuation and dual-modal density measurement technology used radial basis function (RBF) neural network for identifying the flow pattern and determining the void fraction in gas–liquid two-phase flows independent of the liquid phase changes of gas–liquid two-phase flow (Roshani et al. 2017). Neural networks have also been successfully applied to electroencephalogram (EEG) signal analysis (Shrestha et al. 2019; Michielli et al. 2019; Tang et al. 2020), financial time series analysis (Araujo et al. 2019), and stock price prediction (Qiu et al. 2020).

So far in the field of multiphase flow, we have made some achievements in the identification of mixed fluid flow patterns. However, in terms of accuracy and speed of flow pattern recognition, we need to conduct more research deeply. Therefore, this paper proposes to perform MEMD on the four channel conductance signal and use the IMF normalized energy as the characteristic values and utilize the ELM for training to achieve accurate identification of oil–water two-phase flow patterns. In this paper, firstly, we verify the MEMD is more suitable for simultaneous decomposition of multi-channel signals than the EMD and the EEMD. Then, we perform the MEMD on multivariate conductance signal and divide the high- and low-frequency components of the IMF. Finally, we adopt high-frequency components eigenvalues for training with ELM. The experimentally results verify that the proposed method can identify the flow patterns of oil–water two-phase flow quickly and accurately.

Methodology

Multivariate empirical mode decomposition (MEMD) theory

The MEMD decomposes the signal into several IMF, and each IMF component contains the local characteristic signal of the original signal at different scales. Each component of IMF represents each frequency component in the original signal, and is arranged according to the high-frequency components and low-frequency components. MEMD can decompose multiple signals and obtain the same mode of different channels. We can use the MEMD to process the multi-channel signal of oil–water two-phase flow. This method solves the mode calibration problem of multi-channel signals.

Given an n-dimensional data vector \(\left\{ {a\left( t \right)\}_{t = 1}^{\text{T}} = \{ a_{1} (t),a_{2} (t), \ldots a_{n} (t)} \right\}\), where T is the length of the signal, \(\alpha^{k} = \left\{ {\alpha_{1}^{k} ,\alpha_{2}^{k} \ldots \alpha_{n - 1}^{k} } \right\}\) is an angle of \(n - 1\) dimensional sphere, \(x^{\alpha k} = \left[ {x_{1}^{k} ,x_{2}^{k} \ldots x_{n}^{k} } \right]\) is the set of direction vectors for \(\alpha^{k}\). We assume that K direction vectors are established on a sphere, \(k = 1,2,3, \ldots K\). The specific algorithm steps (Rehman and Mandic 2010) of the MEMD are as follows:

  1. (1)

    Choose a sample pointset for direction vectors on the \(n - 1\) sphere.

  2. (2)

    Compute the projection of the original input signal \(\left\{ {a\left( t \right)} \right\}_{t = 1}^{\text{T}}\) on the kth direction vector \(x^{\alpha k}\), recorded as \(\left\{ {p^{\alpha k} \left( t \right)} \right\}_{t = 1}^{\text{T}}\).

  3. (3)

    Find the extreme value corresponding to the instantaneous moment of the projection signal of direction vector, recorded as \(\left\{ {t_{i}^{\alpha k} } \right\}_{k = 1}^{K}\), \(i\) is the position of the extreme point, \(i \in \left[ {1,\;T} \right]\).

  4. (4)

    Interpolate \(\left[ {t_{i}^{\alpha k} ,\;a\left( {t_{i}^{\alpha k} } \right)} \right]\) to obtain multivariate envelope curves \(\left\{ {e^{\alpha k} \left( t \right)} \right\}_{k = 1}^{K}\).

  5. (5)

    The mean of the envelope curves of the K direction vectors can be obtained

    $$m\left( t \right) = \frac{1}{K}\sum\limits_{k = 1}^{K} {e^{\alpha k} \left( t \right)}$$
    (1)
  6. (6)

    Obtain the IMF by \(d\left( t \right) = a\left( t \right) - m\left( t \right)\), if \(d\left( t \right)\) meets the criterion of multivariate-modal function (Rilling et al. 2003), then substitute \(a\left( t \right) - m\left( t \right)\) to step (2) as the original input signal. We obtain new multivariate IMF components from step (2) to step (6) and repeatedly perform MEMD decomposition. The original multivariate signal \(\left\{ {a\left( t \right)\}_{t = 1}^{\text{T}} = \{ a_{1} (t),a_{2} (t) \ldots a_{n} (t)} \right\}\) is decomposed into multiple IMF components \(\left\{ {d\left( t \right)} \right\}_{i = 1}^{q}\) and the remainder \(r\left( t \right)\). As shown in Eq. (2).

    $$a\left( t \right) = \sum\limits_{i = 1}^{q} {d_{i} \left( t \right)} + r\left( t \right)$$
    (2)

    where q indicates the number of IMF.

    \(d_{i} (t) = \{ d_{i,1} (t),d_{i,2} (t),d_{i,n} (t)\}\) and \(r_{i} (t) = \{ r_{i,1} (t),r_{i,2} (t),r_{i,n} (t)\}\) correspond to n sets of IMF and n residual components of the n dimension original signal, respectively.

Extreme learning machine (ELM) theory

ELM is a typical fast and efficient algorithm in the neural network family (Huang et al. 2004). This algorithm has the unique advantages of saving time and high accuracy. The calculation structure of the ELM is consisted of m arbitrary samples \(\left( {x_{j} ,\;t_{j} } \right)\). Its structure can be expressed as an equation as:

$$\sum\limits_{i = 1}^{L} {\beta_{i} g\left( {W_{i} \cdot x_{j} + b_{i} } \right)} = o_{j} ,\quad j = 1, \ldots m$$
(3)

where \(W_{i} = \left[ {W_{i1} ,W_{i2} \cdots ,W_{im} } \right]^{\text{T}}\) represents the weight vector connecting the input nodes to the ith hidden node. \(\beta_{i}\) is the weight vector connecting the output nodes to the ith hidden node. \(b_{i}\) is the thresholds of the ith hidden node. During the training process, we can use an activation function with L hidden nodes. In this process, we employ the ELM to approximate the error between the output values of these m samples and the expected value to zero. We can get the following equation:

$$\sum\limits_{j = 1}^{m} {\left\| {o_{j} } \right.} - \left. {t_{j} } \right\| = 0$$
(4)

where \(o_{j}\) represents the actual value of the ELM output. We hope to get \(b_{i} ,\;W_{i}\) and \(\beta_{i}\) so that

$$\sum\limits_{i = 1}^{L} {\beta_{i} h\left( {W_{i} \cdot x_{j} + b_{i} } \right)} = t_{j} ,\quad j = 1, \ldots m,$$
(5)

Equation (5) can be abbreviated as

$$H\beta = T$$
(6)

where \(\beta\) is the output weight, \(T\) is the expected value of the output, and \(H\) is the output value of the hidden layer nodes. The detailed equation is as

$$\begin{aligned} & H\left( {W_{1} , \cdots W_{L} ,b_{1} , \cdots b_{L} ,x_{1} , \cdots x_{L} } \right) \\ & \quad = \left( {\begin{array}{*{20}c} {h\left( {W_{1} \cdot x_{1} + b_{1} } \right)} & \cdots & {h\left( {W_{L} \cdot x_{1} + b_{L} } \right)} \\ \vdots & \ddots & \vdots \\ {h\left( {W_{1} \cdot x_{N} + b_{1} } \right)} & \cdots & {h\left( {W_{L} \cdot x_{N} + b_{L} } \right)} \\ \end{array} } \right)_{N \times L} \\ \end{aligned}$$
(7)

During the ELM training process, we hope to get \(W_{i} ,\;b_{i}\) and \(\beta_{i}\) so that

$$\left\| {H\left( {\hat{W}_{i} ,\hat{b}_{i} } \right)\hat{\beta }_{i} - T} \right\| = \mathop {\hbox{min} }\limits_{W,b,\beta } \left\| {H\left( {W_{i} ,b_{i} } \right)\beta_{i} - T} \right\|$$
(8)

Equation (8) is equivalent to minimizing the loss function:

$$e = \sum\limits_{j = 1}^{N} {\left( {\sum\limits_{i = 1}^{L} {\beta_{i} h\left( {W_{i} \cdot x_{j} + b_{i} } \right) - t_{j} } } \right)}$$
(9)

For the training of a single hidden layer ELM neural network, it can be transformed into a linear system problem. Therefore, the value of the output weight can be determined

$$\hat{\beta } = H^{ + } T.$$
(10)

where \(H^{ + }\) is the generalized inverse of the matrix, and the resulting norm solution is unique and minimal. During the ELM algorithm learning process, once the bias and input weight of the hidden layer are determined, the output matrix that uniquely determines the hidden layer can be obtained. The main feature of the ELM is that the number of hidden layer nodes can be selected randomly or artificially. From this perspective, the learning process of the ELM only needs to calculate the value of the output weight. Therefore, the biggest advantages of the ELM are fast learning speed and good generalization performance.

IMF normalized energy

The time scale of the energy distribution is an important parameter in the signal analysis process. When water and oil phases with different flow rates are passed through the 4-electrode distributed conductivity sensor pipeline, the frequency component signals will also change accordingly. This means that the energy of each frequency component signal contains a large amount of information of different flow patterns. Theoretically, we obtain each group of signals of the MEMD decomposition which is arranged according to the high-frequency components and low-frequency components. Therefore, the IMF energy can be used as the characteristic vector for the conductance fluctuation signals under different flow patterns.

We can obtain a set of IMF normalized energy as eigenvalues to represent different flow patterns. The calculation method is as follows

$${\text{IER}}\left( i \right) = \frac{{{\text{IE}}_{{{\text{IMF}}_{i} \left( j \right)}} }}{{{\text{IE}}_{x\left( j \right)} }} = \sqrt {\frac{{{\raise0.7ex\hbox{$1$} \!\mathord{\left/ {\vphantom {1 {N_{\text{s}} }}}\right.\kern-0pt} \!\lower0.7ex\hbox{${N_{\text{s}} }$}}\sum\nolimits_{j = 1}^{{N_{\text{s}} }} {\left| {{\text{IMF}}_{i} \left( j \right)} \right|^{2} } }}{{{\raise0.7ex\hbox{$1$} \!\mathord{\left/ {\vphantom {1 N}}\right.\kern-0pt} \!\lower0.7ex\hbox{$N$}}\sum\nolimits_{j = 1}^{N} {\left| {x\left( j \right)} \right|^{2} } }}}$$
(11)

where \({\text{IE}}_{{{\text{IMF}}_{i} \left( j \right)}}\) represents the energy of the ith IMF component of the jth conductance signal. \({\text{IE}}_{x\left( j \right)}\) indicates the energy of the jth conductance signal, \(x\left( j \right)\) is the jth conductance signal, \({\text{IER}}\left( i \right)\) represents the ith normalized energy, and \(N,\;N_{\text{s}}\) are the length of the conductance signal and the length of the IMF data, respectively.

In general, the length of the IMF data and the length of the conductance signal are both equal to N. Therefore, Eq. (11) can be expressed as

$${\text{IER}}\left( i \right) = \frac{{{\text{IE}}_{{{\text{IMF}}_{i} \left( j \right)}} }}{{{\text{IE}}_{x\left( j \right)} }} = \sqrt {\frac{{\sum\nolimits_{j = 1}^{N} {\left| {{\text{IMF}}_{i} \left( j \right)} \right|^{2} } }}{{\sum\nolimits_{j = 1}^{N} {\left| {x\left( j \right)} \right|^{2} } }}}$$
(12)

It is obvious that the normalized energy can extract the value between the specific root value IMF component and the signal from the root mean square. Therefore, the IMF component is highly correlated with various flow pattern conductance signals. That is, the greater the normalized energy, the greater the energy of the IMF, and the more information it contains.

MEMD applied to Lorenz systems

In order to verify the MEMD is more suitable for simultaneous decomposition of multi-channel signals than EMD (Huang and Wu 2008) and EEMD (Wu and Huang 2009), we perform the EMD, EEMD, and MEMD on the Lorenz system. The expression of Lorzen system is as follows

$$\left\{ {\begin{array}{*{20}l} {\dot{x} = - ax + ay} \hfill \\ {\dot{y} = bx - xz - y} \hfill \\ {\dot{z} = xy - cz} \hfill \\ \end{array} } \right.$$
(13)

where \(a = 10,\;b = 28,\;c = \frac{8}{3}\), and initial value is \(\left( {1,\;1,\;1} \right)\).

The decomposition result of the EMD is depicted in Fig. 1a–c that X and Y can get equal number of the IMFs, while Z has 8 IMFs. From Fig. 2a–c, it can be seen that the EEMD can obtain the same number of IMF components, but there may be a problem of inconsistent frequency bands of IMF component signals of the same scale. Therefore, the EEMD cannot be used for multi-channel simultaneous decomposition. As shown in Fig. 3a–c, the MEMD can be used to analyze multi-channel signals simultaneously, and it can generate the same number of IMF components. We can draw a conclusion that the EMD and the EEMD are only suitable for processing single-channel data, cannot be used to analyze data generated by multi-channel signals at the same time, the MEMD can be used to analyze multi-channel signals.

Fig. 1
figure 1

EMD results of the Lorenz system

Fig. 2
figure 2

EEMD results of the Lorenz system

Fig. 3
figure 3

MEMD results of the Lorenz system

Experiments and results analysis

We observed 6 flow patterns during the experiment: stratified flow pattern (ST), stratified flow with mixing at an interface pattern (ST&MI), dispersion of oil in water and water flow pattern (DO/W&W), dispersion of water in oil and oil in water flow pattern (DW/O&D O/W), dispersion of oil in water flow pattern (DO/W), dispersion of water in oil flow pattern (DW/O). ST, ST&MI, DO/W&W, DO/W and DW/O&DO/W widely exist in the actual oil pipeline transportation process. Therefore, this paper focuses on this five flow patterns. Figure 4 shows the conductance signals of four electrodes of five flow patterns measured experimentally.

Fig. 4
figure 4

Four-electrode conductivity signals for five typical flow patterns

We employ the DO/W flow pattern as an example to perform the MEMD decomposition, and the results are shown in Fig. 5. It can be seen from that the conductance signals of four channels are decomposed into seven IMF components and a residual component \(r\). We obtain the IMF components by the MEMD and divide them into high-frequency components, low-frequency components, and trend components. The high-frequency components are IMF1–IMF3, and IMF4–IMF7 are low-frequency components, and r is classified as a trend component.

Fig. 5
figure 5

MEMD decomposition of DO/W

We decompose the conductance signals of five flow patterns through the MEMD to obtain the IMF components, and then calculate the high-frequency component energy and low-frequency component energy, we select high-frequency normalized IMF and low-frequency normalized IMF of different channels as the eigenvalues. Experimental results are shown in Fig. 6.

Fig. 6
figure 6

Recognition rate of 4 channels with different feature values

From Fig. 6, we can see that the conductance signals of the five typical flow patterns obtained from the A channel make it easier to distinguish the various flow patterns. The sensitivity order of the ELM for using the energy normalized by four electrodes as the feature vector is A, D, B, and C. The reason for this phenomenon is that the A channel contains the most amount of information to reflect the characteristics of various flow patterns, so it is the most sensitive to the identification of flow patterns. We can conclude that selecting high-frequency normalized energy is far better than the low-frequency normalized energy as the eigenvalue. The low-frequency components contain characteristic information of the flow patterns less than the high-frequency component by the analysis, the high-frequency components can better reflect the characteristics of the flow patterns.

Figure 7 shows that choosing different activation functions has a great impact on the recognition rate of the test results. It can be seen from Fig. 7 that the choice of sigmoid (sig) is better than the sine (sin) and Hardlim as activation functions for the recognition of the five flow patterns.

Fig. 7
figure 7

Recognition rate of test set with different activation functions

We select the normalized energy of IMF1–IMF3 after the MEMD decomposition of A channel as the feature vector, select sig as the activation function. Then, we train and recognize five typical flow patterns. Each flow pattern has 100,000 data, MEMD decomposition is performed every two thousand, and the first 80% of the normalized energy is used as the training set and the last 20% is used as the test set. The experimental results are shown in Table 1.

Table 1 Test results of test samples

As can be seen from Table 1 that the ELM is used to identify five typical flow patterns, and the overall recognition rate of the convection pattern reaches 94%. Especially flow patterns recognition accuracy of ST, DW/O&DO/W and DO/W three types reaches 100%. Therefore, performing ELM to identify the flow pattern of horizontal oil–water two-phase flow is not only fast, but also has high generalization performance, and the accuracy of flow pattern recognition is high.

Conclusions

Considering the complexity and changeable of oil–water two-phase flow, the accuracy and speed of flow pattern recognition, we propose a method combining the ELM and the MEMD to identify flow patterns of oil–water two-phase flow. According to the analysis, when select the normalized high-frequency component of A channel as the eigenvalues can reach great accuracy. The experimental results indicate that the method increase the speed and improve accuracy of flow pattern recognition of horizontal oil–water two-phase flow. This paper provides a quick strategy and a new vision on the flow patterns characteristics analysis of oil–water two-phase flow.