# Faulty Line Selection Based on Modified CEEMDAN Optimal Denoising Smooth Model and Duffing Oscillator for Un-Effectively Grounded System

Authors: Sizu Hou; Wei Guo
Journal: Mathematical Problems in Engineering
 (2020)

---

## Abstract

As the un-effectively grounded system fails, the zero-sequence current contains strong noise and nonstationary features. This paper proposes a novel faulty line selection method based on modified complete ensemble empirical mode decomposition with adaptive noise (MCEEMDAN) and Duffing oscillator. Here, based on multiscale permutation entropy, fuzzy c-means clustering, and general regression neural network for abnormal signal detection, the MCEEMDAN is proposed. The endpoint mirror method is used to suppress the endpoint effect problem in the decomposition stage. The proposed algorithm is able to decompose the original signal into a series of intrinsic mode functions, which can complete the first filtering. The research shows that it can efficiently suppress the mode confusing phenomenon of empirical mode decomposition (EMD) and is also more complete and orthogonal than ensemble empirical mode decomposition (EEMD) and complementary ensemble empirical mode decomposition (CEEMD). The optimal denoising smooth model is established for choosing optimal intrinsic mode functions to complete the second filtering. It can ensure that the reconstructed filtered signal has better smoothness and similarity. The optimal denoising smooth model of MCEEMDAN can not only keep useful details of the original signal but also reduce the noise and smooth signal. The bifurcation characteristic of the chaotic oscillator is applied in weak signal detection. The zero-sequence current’s denoising result is extracted as the input signal of the Duffing system. The faulty line could be selected by observing the phase diagram of the system. The research results verify the usability and effectiveness of the proposed method.


---

## Body

## 1. Introduction

The un-effectively grounded system is used in the distribution network. Data show that 80% of the fault occurring in the un-effectively grounded system is single-phase grounding fault. When a single-phase grounding fault occurs, the fault current is pretty small, and the zero-sequence current is nonstationary and nonlinear [1, 2]. If the faulty line remains unsolved, it could cause a serious threat to the insulation of the distribution network. It is important for stable operation and safe of the distribution network to select the faulty line accurately [3, 4]. But the adaptability and accuracy of traditional faulty line selection methods are not enough due to the weak fault signal and the complex working condition.Most existing faulty line selection methods use the injection signal method, steady-state information method, or transient information method. The injection signal method needs additional signal injection equipment, and its price is high, and its engineering operation is complex. The characteristic of a fault signal is weak in the steady-state information method, and it is greatly affected by the arc suppression coil, which makes the faulty line selection unreliable [5]. The characteristic of the fault signal in the transient information method is larger than the steady characteristic, and it is not affected by the arc suppression coil [6] and it does not require additional equipment. Liu et al. [7] used wavelet analysis to extract fault information of the zero-sequence current for faulty line selection, but the wavelet analysis is easily affected by noise, and the fault feature band which is extracted by wavelet analysis is inefficient, and the method needs to select appropriate wavelet basis function. Shu et al. [8] used time-frequency characteristics to the analysis of transient zero-sequence current for the faulty line selection, but it is invalid while a high-resistance grounding fault. Fault characteristics as an initial abrupt wave, real power, and transient are employed in the proposed fusion method [9], but its calculation amount is large. The paper [10] shows that the 5th harmonic fault signal at the exit of the faulty line is larger than any nonfaulty line 5th harmonic signal and towards the opposite, but the 5th harmonic fault signal is very weak, and the fault background signal interferes with the 5th harmonic fault signal. Zhang et al. [11] used EMD to decompose the zero-sequence current, and the first intrinsic mode function is extracted as the characteristic parameter of the zero-sequence current. But when the signal contains pulse signal, intermittent signal, and noisy signal, the mode mixing of EMD would arise and cause wrong judgment [12, 13]. Ensemble empirical mode decomposition (EEMD) can restrain modal aliasing of EMD at a certain level; however, the calculation amount grows, and the completeness loses due to residual noise [14, 15]. Complementary ensemble empirical mode decomposition (CEEMD) can decrease the reconstruction error by adding white noise to the target signal and using the ensemble mean to extract intrinsic mode functions [16, 17]. However, its calculation speed is slow, and the efficiency is low. Some effective characteristics will be lost when the high-frequency signal is filtered, which will affect the denoising effect. Complete ensemble empirical mode decomposition with adaptive noise (CEEMDAN) adds adaptive white noise to the original signal and calculates unique residual signal to obtain every intrinsic mode function (IMF) [18, 19]. The CEEMDAN algorithm can improve computational efficiency and decrease the reconstruction error. But CEEMDAN still has noise existing in the IMF region, and some IMFs may appear in hysteretic IMFs.In order to overcome the above deficiencies, modified complete ensemble empirical mode decomposition with adaptive noise (MCEEMDAN) is proposed. The modified method based on combing CEEMDAN, general regression neural network, fuzzy c-means clustering, and multiscale permutation entropy is established. The decomposition stage uses an endpoint mirror method to suppress end effects. MCEEMDAN, which restrains the mode mixing problem, has excellent adaptability. The approach can also sift out the pseudocomponents generated by the signal decomposition. It is complete and orthogonal than EEMD and CEEMD. The optimal denoising smooth mathematical model, which weights the contradiction between similarity and smoothness of the filtered signal, is established for completing signal reconstruction and effectively extracting useful information from the original signal. The zero-sequence current in each line, which is processed by the MCEEMDAN optimal denoising smooth mathematical model, is extracted as the periodic external dynamic of the Duffing system. The trisection symmetry phase estimation method is applied for searching the critical phase, and each input signal is moved according to the critical phase. The faulty line could be selected by observing the phase diagram of the Duffing system.


## 2. Multiscale Permutation Entropy

The permutation entropy (PE) can describe the random and dynamic behavior of the time series [20]. It is based on a comparison of neighborhood values that PE measures signal complexity, and it does not consider the specific value of the time series. PE includes many advantages, such as simple algorithms, high computational efficiency, good robustness, strong anti-interference ability, and so forth. The PE algorithm is described below. Given a time series {x(i) i = 1, 2, …, N} with the length N, the time series which is constructed by using the following equation is defined as [21](1)X1=x1,x1+τ,…,x1+m−1τXi=xi,xi+τ,…,xi+m−1τ⋮XN−m−1τ=xN−m−1τ,xN−m−2τ,…,xN,where τ represents the time delay and m represents the embedding dimension. Thus, X(i) is arranged in ascending sorting, and the permutation is defined as(2)xi+j1−1τ≤xj+j2−1τ≤⋯≤xi+i+jm−1τ.Forxi+ji1−1τ=xi+ji2−1τ, the permutation is ordered as xi+ji1−1τ≤xi+ji2−1τ for ji1≤ji2. Any vector X(i) has a unique symbol sequence as follows:(3)Sg=j1,j2,…,jm,where g = 1∼k and k ≤ m!. The ascending sequence is one of the m! permutations of m-tuple vector. The occurrence frequency of each permutation is defined as Pi(i = 1, 2, …, k). The PE is defined by the Shannon entropy, which is calculated as follows:(4)Hpm=−∑g=1kPglnPg,where 0≤Hp≤lnm!, and the normalized PE can be defined as(5)PE=Hpmlnm!,where 0 ≤ PE ≤ 1. PE is an effective method for reacting a slight change in the time series. In order to solve the single-scale problem of PE, multiscale permutation entropy (MPE) is proposed by the paper [22]. The most important step of calculating MPE is acquiring coarse-grained time series by using the following equation:(6)ysj=1s∑i=j−1s+1jsxi,1≤j≤Ns,where s is the scale factor and ys(j) is a coarse-grained time series. PE of each coarse-grained time series is calculated as MPE.Before computing MPE, three important parameters, including the embedding dimensionm, the time delay τ, and the scale factor s are needed to set. If the embedding dimension m is too small, it is difficult to distinguish the randomness of the signal. The method will be meaningless. If the m is too large, the algorithm could take a long time to calculate, and it is difficult to observe the subtle change of the signal. Zheng et al. [23] suggested that the embedding dimension m is set to 6 in this paper. The time delay τ is set to 1 in this paper because it has little effect on the permutation entropy. If the maximal scale factor smax is >10, the MPE can reflect the important information of the signal [24]. The scale factor s is set to 12 in this paper.


## 3. FCM-GRNN Identifies Abnormal Signal

### 3.1. Fuzzy C-Means Clustering

Many areas, such as data mining and pattern recognition, use fuzzy c-means clustering (FCM) to solve their problems. FCM, which is an unsupervised learning, uses a membership degree to determine the value of each sample objective, which belongs to a cluster. FCM continuously modifies the cluster centers and optimizes the cost function by calculating the Euclidean distance between the data points and the cluster center until the termination criterion is satisfied. In FCM, membership values range between 1 and 0. The label of each sample objective is based on the maximum coefficient of the membership. Similar sample objects are assigned to the same cluster.X = {x1, x2, …, xn} is a data set that is composed of n samples of dimension. It is divided into C = {c1, c2, …, cc} fuzzy clusters. The cost function is defined as [25](7)JU,c1,...,cc=∑i=1c∑j=1nuijmdij2s.t.∑i=1cuij=1,∀j=1,2,...,n,where m ∈[1, ∞) is the fuzzy weighted exponent, uij is the membership matrix of the jth sample point to the ith cluster, and dij is the Euclidean distance between xj and ci. We use Lagrange multipliers λ = [λ1, λ2, …, λn] with the length n to form a new cost function for equation (7) to reach a minimum as follows:(8)J¯U,c1,…,cc,λ1,…,λn=JU,c1,…,cc+∑j=1nλj∑i=1cuij−1=∑i=1c∑jnuijmdij2+∑j=1nλj∑i=1cuij−1,where J¯ is a new cost function, and the necessary conditions for equation (7) to reach its minimum are(9)ci=∑j=1nuijmxj∑j=1nuijm,uij=1∑k=1cdij/dkj2/m−1.We repeatedly update equation (9). If either J’s value is below a certain threshold value or J’s improvement value minus its previous iteration value is below a certain threshold, the constraint condition is satisfied, and the cyclic iteration is stopped. The cluster to which the individual belongs is determined by the membership degree of individual relative to each cluster center, and each cluster is labeled for representing the rule degree of the fuzzy cluster signal.


### 3.2. General Regression Neural Network

The generalized regression neural network (GRNN) is a critical branch of radial basis function neural network. GRNN includes many advantages, such as the strong nonlinear mapping ability, good global convergence, high fault tolerance, fast training speed, and so forth. The schematic diagram of the GRNN is shown in Figure1.Figure 1
The schematic diagram of GRNN.P, input vector; Q, number of input vectors; R, number of sample elements per group; KW1,1, weight matrix of input layer; LW2,1, weight matrix of specific linear layer; ||dist||, distance function; b1, threshold of hidden layer; n1, n2, each layer output vector; a1, a2, linear transfer function; nprod, weight function of output layer.GRNN consisted of three layers, including the input layer, radial basis layer, and specific linear layer [26]. The input layer is composed of neurons and does not participate in the operation. It only transfers the input vectors to the radial basis layer, and the number of layers is equal to the dimension of the input vector. The number of radial basis layers is equal to the number of the learning samples, and the neuron transfer function is Gaussian function pi as follows:(10)pi=expX−XiTX−Xi2δ2,i=1,2,…,n,where n is the number of radial basic neurons, X is the input variable, Xi is the sample observation corresponding to the ith neuron, and δ is the smoothness factor. The transfer function of the output layer is a linear function purelin. The final network output value is defined as(11)y=purelinLW2,1×a1suma1.To sum up, the structure of GRNN is simple, and the smoothness factorδ will affect the execution effect of the GRNN. If the smoothness factor δ is selected randomly by hand, it could have greater subjectivity and affect the output result of GRNN. In this paper, the fruit fly optimization algorithm (FOA) [27] is used to select the optimal smoothness factor δ adaptively.


### 3.3. FCM-GRNN Identifies Abnormal Signal Model

We collect 1600 signals. These signals include 492 regular signals and 1108 random signals. There are white noise signals, Gaussian white noise signals, intermittent signals, and other signals in random signals. The MPEs of the signals are clustered by the FCM algorithm. Because these signals have two states, the clustering center is set to 2. The weighted index is set to 2, and the iterative termination tolerance is set to 10e−3. The input is MPE of each signal. The clustering center of each cluster and fuzzy membership matrix are obtained. We calculate the within-class means of two categories and obtain the distance matrix from all samples in the corresponding category to within-class mean. The 100 samples with the smallest distance are, respectively, selected from the distance matrix of the two categories, and the corresponding network outputs are set to 1 and 2. The MPE of 12 × 200 groups of signals is taken as the input of GRNN, and the output is the corresponding signal category. The FOA algorithm optimizes the smoothness factor of GRNN. The iteration number is set to 100, and the population size is set to 50. The optimization process diagram of FOA-GRNN is shown in Figure 2. As shown in Figure 2, when the smoothness factor is 0.1212, the minimum RMSE is 0.0145.Figure 2
The optimization process diagram of FOA-GRNN.The smoothness factor 0.1212 is brought into the GRNN, and the above MPE of 12 × 200 sample signals is trained to obtain the prediction model. The prediction model is used to analyze the MPE. The MPE as the eigenvector cannot be visualized in a high dimension. In order to show the effect of the model, the principal component analysis (PCA) algorithm is used to reduce the dimensionality of the 12-dimensional MPE and the clustering center vector. Figure3 shows the 2-dimensional spatial distribution and 2-dimensional countermap of the random signals and regular signals after PCA dimension reduction. The MPE eigenvectors of signals are distributed around the two clustering centers according to regularity and randomness. The distinguishing effect of different types of signals is ideal, and the characteristics of the same type of signals are closely related. This shows that the signal’s MPE as the characteristic quantity combining with FCM unsupervised clustering and GRNN tutor learning can complete the random detection of the signal, and the detection rate is maintained at 100%. This method can be used to identify abnormal signals.Figure 3
The clustering contour diagram.


## 3.1. Fuzzy C-Means Clustering

Many areas, such as data mining and pattern recognition, use fuzzy c-means clustering (FCM) to solve their problems. FCM, which is an unsupervised learning, uses a membership degree to determine the value of each sample objective, which belongs to a cluster. FCM continuously modifies the cluster centers and optimizes the cost function by calculating the Euclidean distance between the data points and the cluster center until the termination criterion is satisfied. In FCM, membership values range between 1 and 0. The label of each sample objective is based on the maximum coefficient of the membership. Similar sample objects are assigned to the same cluster.X = {x1, x2, …, xn} is a data set that is composed of n samples of dimension. It is divided into C = {c1, c2, …, cc} fuzzy clusters. The cost function is defined as [25](7)JU,c1,...,cc=∑i=1c∑j=1nuijmdij2s.t.∑i=1cuij=1,∀j=1,2,...,n,where m ∈[1, ∞) is the fuzzy weighted exponent, uij is the membership matrix of the jth sample point to the ith cluster, and dij is the Euclidean distance between xj and ci. We use Lagrange multipliers λ = [λ1, λ2, …, λn] with the length n to form a new cost function for equation (7) to reach a minimum as follows:(8)J¯U,c1,…,cc,λ1,…,λn=JU,c1,…,cc+∑j=1nλj∑i=1cuij−1=∑i=1c∑jnuijmdij2+∑j=1nλj∑i=1cuij−1,where J¯ is a new cost function, and the necessary conditions for equation (7) to reach its minimum are(9)ci=∑j=1nuijmxj∑j=1nuijm,uij=1∑k=1cdij/dkj2/m−1.We repeatedly update equation (9). If either J’s value is below a certain threshold value or J’s improvement value minus its previous iteration value is below a certain threshold, the constraint condition is satisfied, and the cyclic iteration is stopped. The cluster to which the individual belongs is determined by the membership degree of individual relative to each cluster center, and each cluster is labeled for representing the rule degree of the fuzzy cluster signal.


## 3.2. General Regression Neural Network

The generalized regression neural network (GRNN) is a critical branch of radial basis function neural network. GRNN includes many advantages, such as the strong nonlinear mapping ability, good global convergence, high fault tolerance, fast training speed, and so forth. The schematic diagram of the GRNN is shown in Figure1.Figure 1
The schematic diagram of GRNN.P, input vector; Q, number of input vectors; R, number of sample elements per group; KW1,1, weight matrix of input layer; LW2,1, weight matrix of specific linear layer; ||dist||, distance function; b1, threshold of hidden layer; n1, n2, each layer output vector; a1, a2, linear transfer function; nprod, weight function of output layer.GRNN consisted of three layers, including the input layer, radial basis layer, and specific linear layer [26]. The input layer is composed of neurons and does not participate in the operation. It only transfers the input vectors to the radial basis layer, and the number of layers is equal to the dimension of the input vector. The number of radial basis layers is equal to the number of the learning samples, and the neuron transfer function is Gaussian function pi as follows:(10)pi=expX−XiTX−Xi2δ2,i=1,2,…,n,where n is the number of radial basic neurons, X is the input variable, Xi is the sample observation corresponding to the ith neuron, and δ is the smoothness factor. The transfer function of the output layer is a linear function purelin. The final network output value is defined as(11)y=purelinLW2,1×a1suma1.To sum up, the structure of GRNN is simple, and the smoothness factorδ will affect the execution effect of the GRNN. If the smoothness factor δ is selected randomly by hand, it could have greater subjectivity and affect the output result of GRNN. In this paper, the fruit fly optimization algorithm (FOA) [27] is used to select the optimal smoothness factor δ adaptively.


## 3.3. FCM-GRNN Identifies Abnormal Signal Model

We collect 1600 signals. These signals include 492 regular signals and 1108 random signals. There are white noise signals, Gaussian white noise signals, intermittent signals, and other signals in random signals. The MPEs of the signals are clustered by the FCM algorithm. Because these signals have two states, the clustering center is set to 2. The weighted index is set to 2, and the iterative termination tolerance is set to 10e−3. The input is MPE of each signal. The clustering center of each cluster and fuzzy membership matrix are obtained. We calculate the within-class means of two categories and obtain the distance matrix from all samples in the corresponding category to within-class mean. The 100 samples with the smallest distance are, respectively, selected from the distance matrix of the two categories, and the corresponding network outputs are set to 1 and 2. The MPE of 12 × 200 groups of signals is taken as the input of GRNN, and the output is the corresponding signal category. The FOA algorithm optimizes the smoothness factor of GRNN. The iteration number is set to 100, and the population size is set to 50. The optimization process diagram of FOA-GRNN is shown in Figure 2. As shown in Figure 2, when the smoothness factor is 0.1212, the minimum RMSE is 0.0145.Figure 2
The optimization process diagram of FOA-GRNN.The smoothness factor 0.1212 is brought into the GRNN, and the above MPE of 12 × 200 sample signals is trained to obtain the prediction model. The prediction model is used to analyze the MPE. The MPE as the eigenvector cannot be visualized in a high dimension. In order to show the effect of the model, the principal component analysis (PCA) algorithm is used to reduce the dimensionality of the 12-dimensional MPE and the clustering center vector. Figure3 shows the 2-dimensional spatial distribution and 2-dimensional countermap of the random signals and regular signals after PCA dimension reduction. The MPE eigenvectors of signals are distributed around the two clustering centers according to regularity and randomness. The distinguishing effect of different types of signals is ideal, and the characteristics of the same type of signals are closely related. This shows that the signal’s MPE as the characteristic quantity combining with FCM unsupervised clustering and GRNN tutor learning can complete the random detection of the signal, and the detection rate is maintained at 100%. This method can be used to identify abnormal signals.Figure 3
The clustering contour diagram.


## 4. Modified CEEMDAN Optimal Denoising

### 4.1. Modified CEEMDAN Algorithm

EEMD and CEMMD make the distribution of the extreme points of the original signal more uniform by adding white noise, and the added noise covers the abnormal signals, such as high frequency, intermittent, or noise in the original signal. Therefore, the problem of modal confusion can be reduced. However, the algorithms need to limit the number of iterations, and the decomposed components do not necessarily satisfy the definition of the IMF. From the perspective of the application, it is meaningless to sacrifice the accuracy of the component and the physical meaning of the instantaneous frequency of the component for the adaptive decomposition.In fact, after the abnormal signal is decomposed, it is not necessary to use EMD to decompose the noises completely, and as long as the decomposition completeness of the added noise signal can be ensured. In order to solve these problems, MCEEMDAN is proposed in this paper. The MCEEMDAN can avoid unnecessary ensemble mean in EEMD and CEEMD. It can make the decomposition results satisfy the definition of inherent mode function. The proposed method can reduce the reconstruction error caused by adding white noise and ensure the completeness of decomposition. The key step of MCEEMDAN is the detection of abnormal signals. Given the advantages of the MPE-FCM-GRNN, MCEEMDAN can be designed. Here, based on MPE-FCM-GRNN for abnormal signal detection, the adaptability of MCEEMDAN also can be improved.Supposing the operatorEk(·) is the kth intrinsic mode function, which is decomposed by EMD. IMFq′ is the qth intrinsic mode function, which is decomposed by CEEMDAN. I(t) is the original signal. ni(t) stands for the ith white noise with the mean value of zero and the variance of one (i = 1, 2, …, NR, and NR is the times of performing EMD); βi−1 stands for adaptive coefficient (i = 1, 2, …, NR). EMD adopts the end point mirror image method during the decomposition stage. The steps of MCEEMDAN are described as follows:(1)
The first mode of CEEMDAN is obtained by EMD NR times to decompose the signalI(t) + β0ni(t) as follows:(12)IMF1′=1NR∑i=1NRE1It+β0nit.(2)
The first residue is as follows:(13)r1t=It−IMF1′.(3)
The second mode of CEEMDAN is obtained by EMD NR times to decompose the signalr1(t) + E1(ni(t)) as follows:(14)IMF2′=1NR∑i=1NRE1r1t+β1E1nit.(4)
Whenk = 2, 3, …, N, the kth residue and the (k + 1)th mode are as follows:(15)rkt=rk−1t−IMFk′,IMFk+1′=1NR∑i=1NRE1rkt+βkE1nit.(5)
Step (4) is repeated until either the mean value of the upper and lower extreme envelope of the residue is zero or the difference between the number of the zero points and the number of the extreme points is less than two.(6)
The MPE of each intrinsic mode function of CEEMDAN decomposition is analyzed by using trained GRNN. When the output classification result is the same as the random signal, it is considered to be an abnormal signal and is removed from the original signal. Supposing thepth mode is an abnormal signal, then the residue signal R(t) is defined as follows:(16)Rt=It−∑i=1pIMFi′.(7)
TheR(t) is decomposed by EMD to obtain the MCEEMDAN’s intrinsic mode functions, which are arranged in the order from high frequency to low frequency.The signalx is composed of sinusoidal signal x1 = 8sin (2π40t + π/2), amplitude modulation signal x2 = (t + 1) sin (2π10t + π/3), and intermittent random signal x3. The sampling frequency is 1000 Hz, and the total time is 2 s. Then, the waveforms of simulation signal and its components are shown in Figure 4. Because the results of CEEMD and EEMD are the same, only the result of EEMD decomposition is given. The x is, respectively, decomposed by EMD, EEMD, CEEMDAN, and MCEEMDAN, and the results are shown in Figure 5. It can be seen from waveform in Figure 5 that there is an obvious problem of modal aliasing in the EMD algorithm. EEMD and CEEMD can restrain modal aliasing of EMD at a certain level, but there are still a lot of pseudocomponents. The decomposition effect of the MCEEMDAN algorithm is better. The noise signal with strong randomness is filtered out. Both the residual component and the pseudocomponents are close to zero. The end effect problem can be suppressed. The decomposition results are regular and stable.Figure 4
The waveforms of simulation signal and its components.Figure 5
The decomposition results of four methods. (a) Decomposition results of EMD. (b) Decomposition results of EEMD. (c) Decomposition results of CEEMDAN. (d) Decomposition results of MCEEMDAN.
(a)
(b)
(c)
(d)Figure6 shows the Hilbert–Huang spectrum of four decomposition methods. As shown in Figure 6, EMD, EEMD, and CEEMDAN have the problem of instantaneous frequency instability, but the MCEEMDAN decomposition method can effectively extract the frequency information of the signal, and the frequency is more stable. The quantitative analysis of the four decomposition methods is based on the correlation between the components of the original signal and the intrinsic mode functions of four decomposition methods as follows:(17)corr=C1,2C1,1∗C2,2,C=covIMFj,xi,where j is the numbering of intrinsic modal function and i = 1, 2, 3 denotes signal x1, x2, x3, and C is a covariance matrix. As the correlation degree |corr| is close to 1, the correlation is high. The IMF correlation analysis of four decomposition methods can be known from Tables 1–4. EMD cannot decompose the effective components of the original signal. Although x1, x2 can be extracted by EEMD, both IMF3-IMF4 and x2 have a high correlation, indicating that there is modal aliasing in EEMD. The CEEMDAN has a high correlation in both IMF4-IMF5 and x1, which indicates that CEEMDAN has the problem of model mixing. MCEEMDAN can effectively extract x1 and x2 from x, and the only IMF1 is highly correlated with x1, and the only IMF2 is highly correlated with x2. As described above, MCEEMDAN can restrain modal aliasing. The orthogonality index can evaluate the endpoint effect of the algorithm. The smaller the orthogonal index is, the better the orthogonality of each intrinsic modal function is, and the lower the influence of the endpoint effect on the algorithm is. Table 5 shows orthogonality indexes for different methods. As shown in Table 5, the endpoint effect of the MCEEMDAN decomposition method is effectively suppressed.Figure 6
Hilbert–Huang spectrums of four methods. (a) EMD. (b) EEMD. (c) CEEMDAN. (d) MCEEMDAN.
(a)
(b)
(c)
(d)Table 1
EMD decomposition IMF correlation.
x
1
x
2
x
3
IMF1
0.7449
0.0051
0.0282
IMF2
0.6273
0.1157
0.0416
IMF3
0.0013
0.7748
0.0105
IMF4
0.0003
0.0034
−0.0055
IMF5
−0.0006
0.0112
0.0356
IMF6
0.0003
−0.0124
0.0261
RS7
0.0010
−0.0143
0.0087Table 2
EEMD decomposition IMF correlation.
x
1
x
2
x
3
IMF1
0.0463
0.0084
0.7371
IMF2
0.9994
−0.0006
−0.0042
IMF3
0.2515
0.9442
−0.0401
IMF4
0.0110
0.9545
−0.0323
IMF5
0.0023
0.5837
−0.0413
IMF6
0.0008
0.0070
0.0364
IMF7
0.0001
0.0027
−0.0009
IMF8
0.0006
−0.0109
−0.0375
IMF9
0.0001
−0.0092
−0.0232
RS10
−0.0003
−0.0020
−0.0180Table 3
CEEMDAN decomposition IMF correlation.
x
1
x
2
x
3
IMF1
0.0329
0.0086
0.5025
IMF2
0.0566
0.0310
0.1209
IMF3
0.1060
0.0081
0.0323
IMF4
0.9886
−0.0018
0.0341
IMF5
0.9913
0.0085
0.0234
IMF6
0.0347
0.6953
0.0216
IMF7
0.0010
0.9761
0.0411
IMF8
0.0005
0.1301
0.0049
IMF9
−0.0006
−0.0022
0.0055
IMF10
−0.0009
0.0093
0.0258
IMF11
−0.0011
0.0153
0.0415
RS12
−0.0009
0.0192
−0.0112Table 4
MCEEMDAN decomposition IMF correlation.
x
1
x
2
x
3
IMF1
0.9990
0.0043
0.0301
IMF2
0.0002
0.9781
0.0066
IMF3
−0.0001
0.0354
−0.0061
IMF4
−0.0012
0.0211
−0.0091
IMF5
−0.0010
0.0140
0.0352
RS6
−0.0009
0.0190
−0.0194Table 5
Orthogonality indexes for different methods.
Method type
Orthogonality index
EMD
0.0276
EEMD
0.0582
CEEMDAN
0.2459
MCEEMDAN
9.5775e−4The simulation results show that the proposed algorithm has better decomposition results than the existing algorithms. MCEEMDAN can restrain model mixing and overcome the problem of the endpoint effect. Its frequency of the Hilbert–Huang spectrum is stable. The orthogonality of the intrinsic mode function is good. These advantages provide a good basis for signal reconstruction.


### 4.2. Optimal Denoising Smooth Model

Different combinations of intrinsic mode functions can construct different reconstructed signals. The optimal denoising smooth model is constructed by the similarity degree and smoothness index [28–30]. The standard deviation of the difference between the reconstructed signal HP(t) and the original signal x(t) is calculated as the similarity degree as follows:(18)Xt=xt−HPt,Amse=stdXt=∑i=1NXi−X¯2N−1,where X(t) is the difference between x(t) and HP(t), N is the number of sampling points, and Amse represents the standard deviation. The smaller the Amse is, the closer the reconstructed signal is to the original signal. When point x1 is on the curve f(x), the left and right curvatures at the point x1 are as follows:(19)K−=f−x1″1+f−x1′23/2,K+=f+x1″1+f+x1′23/2.Whenf(x) is smooth and differentiable at point x1, the left and right curvature at point x1 are equal as follows:(20)K−=K+⟶f−x1″=f+x1″.The second derivative is expanded according to the discrete formula as follows:(21)f−x1″≈fx1−2h−2fx1−h+fx1h2,f+x1″≈fx1+2h−2fx1+h+fx1h2,where h is the step length. Then, the smoothness SNx=x1 of the curve f(x) at point x1 is defined as(22)SNx=x1=fx1+2h−fx1−2h−2fx1+h−fx1−h.Supposing smoothness index Asmse is the standard deviation of SN of all points excluding the left and right endpoints in the sample. The smaller the Asmse is, the smoother the signal is.Considering the smoothness of the filtering algorithm and the error standard deviation between the reconstruction result and the original signal, a weight coefficientμ(0 ≤ μ ≤ 1) is set, and then the objective function is defined as(23)Aminf=μAmse+1−μAsmse.When Aminf is the minimum value, the objective function obtains the optimal solution. The algorithm model corresponding to the optimal solution is the optimal denoising smooth model.


## 4.1. Modified CEEMDAN Algorithm

EEMD and CEMMD make the distribution of the extreme points of the original signal more uniform by adding white noise, and the added noise covers the abnormal signals, such as high frequency, intermittent, or noise in the original signal. Therefore, the problem of modal confusion can be reduced. However, the algorithms need to limit the number of iterations, and the decomposed components do not necessarily satisfy the definition of the IMF. From the perspective of the application, it is meaningless to sacrifice the accuracy of the component and the physical meaning of the instantaneous frequency of the component for the adaptive decomposition.In fact, after the abnormal signal is decomposed, it is not necessary to use EMD to decompose the noises completely, and as long as the decomposition completeness of the added noise signal can be ensured. In order to solve these problems, MCEEMDAN is proposed in this paper. The MCEEMDAN can avoid unnecessary ensemble mean in EEMD and CEEMD. It can make the decomposition results satisfy the definition of inherent mode function. The proposed method can reduce the reconstruction error caused by adding white noise and ensure the completeness of decomposition. The key step of MCEEMDAN is the detection of abnormal signals. Given the advantages of the MPE-FCM-GRNN, MCEEMDAN can be designed. Here, based on MPE-FCM-GRNN for abnormal signal detection, the adaptability of MCEEMDAN also can be improved.Supposing the operatorEk(·) is the kth intrinsic mode function, which is decomposed by EMD. IMFq′ is the qth intrinsic mode function, which is decomposed by CEEMDAN. I(t) is the original signal. ni(t) stands for the ith white noise with the mean value of zero and the variance of one (i = 1, 2, …, NR, and NR is the times of performing EMD); βi−1 stands for adaptive coefficient (i = 1, 2, …, NR). EMD adopts the end point mirror image method during the decomposition stage. The steps of MCEEMDAN are described as follows:(1)
The first mode of CEEMDAN is obtained by EMD NR times to decompose the signalI(t) + β0ni(t) as follows:(12)IMF1′=1NR∑i=1NRE1It+β0nit.(2)
The first residue is as follows:(13)r1t=It−IMF1′.(3)
The second mode of CEEMDAN is obtained by EMD NR times to decompose the signalr1(t) + E1(ni(t)) as follows:(14)IMF2′=1NR∑i=1NRE1r1t+β1E1nit.(4)
Whenk = 2, 3, …, N, the kth residue and the (k + 1)th mode are as follows:(15)rkt=rk−1t−IMFk′,IMFk+1′=1NR∑i=1NRE1rkt+βkE1nit.(5)
Step (4) is repeated until either the mean value of the upper and lower extreme envelope of the residue is zero or the difference between the number of the zero points and the number of the extreme points is less than two.(6)
The MPE of each intrinsic mode function of CEEMDAN decomposition is analyzed by using trained GRNN. When the output classification result is the same as the random signal, it is considered to be an abnormal signal and is removed from the original signal. Supposing thepth mode is an abnormal signal, then the residue signal R(t) is defined as follows:(16)Rt=It−∑i=1pIMFi′.(7)
TheR(t) is decomposed by EMD to obtain the MCEEMDAN’s intrinsic mode functions, which are arranged in the order from high frequency to low frequency.The signalx is composed of sinusoidal signal x1 = 8sin (2π40t + π/2), amplitude modulation signal x2 = (t + 1) sin (2π10t + π/3), and intermittent random signal x3. The sampling frequency is 1000 Hz, and the total time is 2 s. Then, the waveforms of simulation signal and its components are shown in Figure 4. Because the results of CEEMD and EEMD are the same, only the result of EEMD decomposition is given. The x is, respectively, decomposed by EMD, EEMD, CEEMDAN, and MCEEMDAN, and the results are shown in Figure 5. It can be seen from waveform in Figure 5 that there is an obvious problem of modal aliasing in the EMD algorithm. EEMD and CEEMD can restrain modal aliasing of EMD at a certain level, but there are still a lot of pseudocomponents. The decomposition effect of the MCEEMDAN algorithm is better. The noise signal with strong randomness is filtered out. Both the residual component and the pseudocomponents are close to zero. The end effect problem can be suppressed. The decomposition results are regular and stable.Figure 4
The waveforms of simulation signal and its components.Figure 5
The decomposition results of four methods. (a) Decomposition results of EMD. (b) Decomposition results of EEMD. (c) Decomposition results of CEEMDAN. (d) Decomposition results of MCEEMDAN.
(a)
(b)
(c)
(d)Figure6 shows the Hilbert–Huang spectrum of four decomposition methods. As shown in Figure 6, EMD, EEMD, and CEEMDAN have the problem of instantaneous frequency instability, but the MCEEMDAN decomposition method can effectively extract the frequency information of the signal, and the frequency is more stable. The quantitative analysis of the four decomposition methods is based on the correlation between the components of the original signal and the intrinsic mode functions of four decomposition methods as follows:(17)corr=C1,2C1,1∗C2,2,C=covIMFj,xi,where j is the numbering of intrinsic modal function and i = 1, 2, 3 denotes signal x1, x2, x3, and C is a covariance matrix. As the correlation degree |corr| is close to 1, the correlation is high. The IMF correlation analysis of four decomposition methods can be known from Tables 1–4. EMD cannot decompose the effective components of the original signal. Although x1, x2 can be extracted by EEMD, both IMF3-IMF4 and x2 have a high correlation, indicating that there is modal aliasing in EEMD. The CEEMDAN has a high correlation in both IMF4-IMF5 and x1, which indicates that CEEMDAN has the problem of model mixing. MCEEMDAN can effectively extract x1 and x2 from x, and the only IMF1 is highly correlated with x1, and the only IMF2 is highly correlated with x2. As described above, MCEEMDAN can restrain modal aliasing. The orthogonality index can evaluate the endpoint effect of the algorithm. The smaller the orthogonal index is, the better the orthogonality of each intrinsic modal function is, and the lower the influence of the endpoint effect on the algorithm is. Table 5 shows orthogonality indexes for different methods. As shown in Table 5, the endpoint effect of the MCEEMDAN decomposition method is effectively suppressed.Figure 6
Hilbert–Huang spectrums of four methods. (a) EMD. (b) EEMD. (c) CEEMDAN. (d) MCEEMDAN.
(a)
(b)
(c)
(d)Table 1
EMD decomposition IMF correlation.
x
1
x
2
x
3
IMF1
0.7449
0.0051
0.0282
IMF2
0.6273
0.1157
0.0416
IMF3
0.0013
0.7748
0.0105
IMF4
0.0003
0.0034
−0.0055
IMF5
−0.0006
0.0112
0.0356
IMF6
0.0003
−0.0124
0.0261
RS7
0.0010
−0.0143
0.0087Table 2
EEMD decomposition IMF correlation.
x
1
x
2
x
3
IMF1
0.0463
0.0084
0.7371
IMF2
0.9994
−0.0006
−0.0042
IMF3
0.2515
0.9442
−0.0401
IMF4
0.0110
0.9545
−0.0323
IMF5
0.0023
0.5837
−0.0413
IMF6
0.0008
0.0070
0.0364
IMF7
0.0001
0.0027
−0.0009
IMF8
0.0006
−0.0109
−0.0375
IMF9
0.0001
−0.0092
−0.0232
RS10
−0.0003
−0.0020
−0.0180Table 3
CEEMDAN decomposition IMF correlation.
x
1
x
2
x
3
IMF1
0.0329
0.0086
0.5025
IMF2
0.0566
0.0310
0.1209
IMF3
0.1060
0.0081
0.0323
IMF4
0.9886
−0.0018
0.0341
IMF5
0.9913
0.0085
0.0234
IMF6
0.0347
0.6953
0.0216
IMF7
0.0010
0.9761
0.0411
IMF8
0.0005
0.1301
0.0049
IMF9
−0.0006
−0.0022
0.0055
IMF10
−0.0009
0.0093
0.0258
IMF11
−0.0011
0.0153
0.0415
RS12
−0.0009
0.0192
−0.0112Table 4
MCEEMDAN decomposition IMF correlation.
x
1
x
2
x
3
IMF1
0.9990
0.0043
0.0301
IMF2
0.0002
0.9781
0.0066
IMF3
−0.0001
0.0354
−0.0061
IMF4
−0.0012
0.0211
−0.0091
IMF5
−0.0010
0.0140
0.0352
RS6
−0.0009
0.0190
−0.0194Table 5
Orthogonality indexes for different methods.
Method type
Orthogonality index
EMD
0.0276
EEMD
0.0582
CEEMDAN
0.2459
MCEEMDAN
9.5775e−4The simulation results show that the proposed algorithm has better decomposition results than the existing algorithms. MCEEMDAN can restrain model mixing and overcome the problem of the endpoint effect. Its frequency of the Hilbert–Huang spectrum is stable. The orthogonality of the intrinsic mode function is good. These advantages provide a good basis for signal reconstruction.


## 4.2. Optimal Denoising Smooth Model

Different combinations of intrinsic mode functions can construct different reconstructed signals. The optimal denoising smooth model is constructed by the similarity degree and smoothness index [28–30]. The standard deviation of the difference between the reconstructed signal HP(t) and the original signal x(t) is calculated as the similarity degree as follows:(18)Xt=xt−HPt,Amse=stdXt=∑i=1NXi−X¯2N−1,where X(t) is the difference between x(t) and HP(t), N is the number of sampling points, and Amse represents the standard deviation. The smaller the Amse is, the closer the reconstructed signal is to the original signal. When point x1 is on the curve f(x), the left and right curvatures at the point x1 are as follows:(19)K−=f−x1″1+f−x1′23/2,K+=f+x1″1+f+x1′23/2.Whenf(x) is smooth and differentiable at point x1, the left and right curvature at point x1 are equal as follows:(20)K−=K+⟶f−x1″=f+x1″.The second derivative is expanded according to the discrete formula as follows:(21)f−x1″≈fx1−2h−2fx1−h+fx1h2,f+x1″≈fx1+2h−2fx1+h+fx1h2,where h is the step length. Then, the smoothness SNx=x1 of the curve f(x) at point x1 is defined as(22)SNx=x1=fx1+2h−fx1−2h−2fx1+h−fx1−h.Supposing smoothness index Asmse is the standard deviation of SN of all points excluding the left and right endpoints in the sample. The smaller the Asmse is, the smoother the signal is.Considering the smoothness of the filtering algorithm and the error standard deviation between the reconstruction result and the original signal, a weight coefficientμ(0 ≤ μ ≤ 1) is set, and then the objective function is defined as(23)Aminf=μAmse+1−μAsmse.When Aminf is the minimum value, the objective function obtains the optimal solution. The algorithm model corresponding to the optimal solution is the optimal denoising smooth model.


## 5. Detection Principle of Duffing Oscillator Model

### 5.1. Duffing Oscillator Model

Duffing equation, which is widely used in weak signal detection, is a typical mathematical model to describe the chaotic property of the system [31, 32]. The chaotic system is sensitive to the periodic driving force. When the total periodic driving force of the system changes, the system will show different trajectories. The Duffing equation is defined as(24)x″+kx′−x3+x5=rcosωt,where k is the damping ratio, r and ω are the amplitude and the angular frequency of the inner periodic driving force, respectively, and −x3 + x5 represents the nonlinear restoring force. Equation (26) is expressed as the differential equation:(25)x′=y,y′=x3−x5−ky+rcosωt.When the damping ratiok and the angular frequency ω of the system are fixed, the trajectory of the system changes with the change of driving force. When r exceeds a certain threshold rc, the trajectory of the system enters the chaotic state. If r continues increasing over a certain threshold rd, the trajectory of the system enters the large-scale periodic motion state [33, 34]. When the signal is detected by the Duffing oscillator, the amplitude of the inner periodic driving force is adjusted to rd. The system is in a critical state of chaos, and then the cofrequency signal to be tested is added to the Duffing system as follows:(26)x″+kx′−x3+x5=rdcosωt+acosωt+φ,where acos(ωt + φ) is a signal to be measured. Duffing’s total driving force F(t) is composed of the internal periodic driving force and external periodic driving force as follows:(27)Ft=rdcosωt+acosωt+φ=rcosωt+θ.Then, the following equations come into existence:(28)r=rd2+2rdacosφ+a2,θ=arctgasinφrd+acosφ.Multiply the signal to be measured by a detection factor, and the inequationa≪rd can be obtained. The value of θ will tend to 0, and its impact on the system will be negligible. At this time, the state of the system is mainly determined by the value of φ. When the value of φ is in the following range as(29)π−arccosa2rd<φ<π+arccosa2rd,the r < rd and the system does not transition to the large-scale periodic motion state. If the value of φ is not in this range, the r >rd and the system can transition to the large-scale periodic motion state.After the fault of the distribution power network, the phases diagram trajectories of the fifth harmonic zero-sequence current of the normal line and faulty line are opposite, so the faulty line can be selected by the phase diagram trajectory of the chaotic oscillator.The above-discussed case is that the internal power angular frequency isω = 1 rad/s. But the frequency of the signal to be detected is a high-frequency signal in the power grid, and it is necessary to transform the scale of the Duffing system as follows:(30)xt=xωτ,x′t=dxtdt=1wx′ωτ,x″t=d2xtdt2=1w2x″ωτ.The detection equation can be defined as(31)1ω2x″ωτ+kωx′ωτ−x3ωτ+x5ωτ=rcosωτ.Equation (31) is expressed as the differential equation:(32)x′=ωy,x″=−ωkx′+ω2x3−x5+rcosωτ,y′=ω−ky+x3−x5+rcosωτ.The above change process is an equal amount of mathematical transformation, but only the scale of the chaotic system is different, and the characteristic of the system equation is not changed. After the change of scale, the state and change rule of the chaotic oscillator are not changed. As long as the value ofω is adjusted, the detection of the high-frequency signal can be completed. Because the field operation is complex and subject to electromagnetic interference, the phase diagram trajectory of the chaotic system will become rough, and excessive noise will affect the detection of the chaotic system. The signal needs to be processed by the optimal denoising smooth model of MCEEMDAN before it is input to the Duffing system.


### 5.2. Trisection Symmetry Phase Estimation

The trisection symmetry phase estimation is used to determine the phase of the signal to be measured relative to the internal driving force [35]. The phase is divided into three equal parts within a 2π circle. The MCEEMAND optimal denoising smooth model of any line is selected as the signal to be measured. When the signal is input into the chaotic oscillator by phase shift at 0, 2π/3, and 4π/3, there must be a phase diagram state which is different from the other two phase diagram states. For example, the Duffing system is a periodic state after the phase shift is 0° and 4π/3, but it is a large-scale periodic motion state after the phase shift is 2π/3, and then the critical phase must be within the (0, 2π/3) or (2π/3, 4π/3). The dichotomy is used to reduce the search area of the critical phase within the (0, 2π/3) and (2π/3, 4π/3) until the critical phase is searched. The results of MCEEMDAN optimal denoising of other lines that are moved according to the critical phase are input into the chaotic system to complete the signal detection. When the output state diagram of one of the lines is different from that of other lines, it is considered to be the faulty line. When all output state diagrams are the same, the bus bar is deemed to be the fault.


## 5.1. Duffing Oscillator Model

Duffing equation, which is widely used in weak signal detection, is a typical mathematical model to describe the chaotic property of the system [31, 32]. The chaotic system is sensitive to the periodic driving force. When the total periodic driving force of the system changes, the system will show different trajectories. The Duffing equation is defined as(24)x″+kx′−x3+x5=rcosωt,where k is the damping ratio, r and ω are the amplitude and the angular frequency of the inner periodic driving force, respectively, and −x3 + x5 represents the nonlinear restoring force. Equation (26) is expressed as the differential equation:(25)x′=y,y′=x3−x5−ky+rcosωt.When the damping ratiok and the angular frequency ω of the system are fixed, the trajectory of the system changes with the change of driving force. When r exceeds a certain threshold rc, the trajectory of the system enters the chaotic state. If r continues increasing over a certain threshold rd, the trajectory of the system enters the large-scale periodic motion state [33, 34]. When the signal is detected by the Duffing oscillator, the amplitude of the inner periodic driving force is adjusted to rd. The system is in a critical state of chaos, and then the cofrequency signal to be tested is added to the Duffing system as follows:(26)x″+kx′−x3+x5=rdcosωt+acosωt+φ,where acos(ωt + φ) is a signal to be measured. Duffing’s total driving force F(t) is composed of the internal periodic driving force and external periodic driving force as follows:(27)Ft=rdcosωt+acosωt+φ=rcosωt+θ.Then, the following equations come into existence:(28)r=rd2+2rdacosφ+a2,θ=arctgasinφrd+acosφ.Multiply the signal to be measured by a detection factor, and the inequationa≪rd can be obtained. The value of θ will tend to 0, and its impact on the system will be negligible. At this time, the state of the system is mainly determined by the value of φ. When the value of φ is in the following range as(29)π−arccosa2rd<φ<π+arccosa2rd,the r < rd and the system does not transition to the large-scale periodic motion state. If the value of φ is not in this range, the r >rd and the system can transition to the large-scale periodic motion state.After the fault of the distribution power network, the phases diagram trajectories of the fifth harmonic zero-sequence current of the normal line and faulty line are opposite, so the faulty line can be selected by the phase diagram trajectory of the chaotic oscillator.The above-discussed case is that the internal power angular frequency isω = 1 rad/s. But the frequency of the signal to be detected is a high-frequency signal in the power grid, and it is necessary to transform the scale of the Duffing system as follows:(30)xt=xωτ,x′t=dxtdt=1wx′ωτ,x″t=d2xtdt2=1w2x″ωτ.The detection equation can be defined as(31)1ω2x″ωτ+kωx′ωτ−x3ωτ+x5ωτ=rcosωτ.Equation (31) is expressed as the differential equation:(32)x′=ωy,x″=−ωkx′+ω2x3−x5+rcosωτ,y′=ω−ky+x3−x5+rcosωτ.The above change process is an equal amount of mathematical transformation, but only the scale of the chaotic system is different, and the characteristic of the system equation is not changed. After the change of scale, the state and change rule of the chaotic oscillator are not changed. As long as the value ofω is adjusted, the detection of the high-frequency signal can be completed. Because the field operation is complex and subject to electromagnetic interference, the phase diagram trajectory of the chaotic system will become rough, and excessive noise will affect the detection of the chaotic system. The signal needs to be processed by the optimal denoising smooth model of MCEEMDAN before it is input to the Duffing system.


## 5.2. Trisection Symmetry Phase Estimation

The trisection symmetry phase estimation is used to determine the phase of the signal to be measured relative to the internal driving force [35]. The phase is divided into three equal parts within a 2π circle. The MCEEMAND optimal denoising smooth model of any line is selected as the signal to be measured. When the signal is input into the chaotic oscillator by phase shift at 0, 2π/3, and 4π/3, there must be a phase diagram state which is different from the other two phase diagram states. For example, the Duffing system is a periodic state after the phase shift is 0° and 4π/3, but it is a large-scale periodic motion state after the phase shift is 2π/3, and then the critical phase must be within the (0, 2π/3) or (2π/3, 4π/3). The dichotomy is used to reduce the search area of the critical phase within the (0, 2π/3) and (2π/3, 4π/3) until the critical phase is searched. The results of MCEEMDAN optimal denoising of other lines that are moved according to the critical phase are input into the chaotic system to complete the signal detection. When the output state diagram of one of the lines is different from that of other lines, it is considered to be the faulty line. When all output state diagrams are the same, the bus bar is deemed to be the fault.


## 6. Case Study

MATLAB-SIMULINK is used in this paper to simulate the un-effectively grounded system. The simulation model consists of four lines. The length of a single line is generally 2 km to 20 km in the un-effectively grounded system. Linel1 is set as 16 km overhead line; line l2 is set as 13 km cable-overhead line including 5 km overhead line and 8 km cable line; line l3 is set as 20 km overhead line; line l4 is set as 10 km cable line. Table 6 shows unit length line parameters. The overcompensation method is adopted in the arc suppression coil. When the compensation degree is set to 8% in this paper, the inductive reactance of the arc suppression coil is 0.4765H. The resistance value of the arc suppression coil takes 2% of the inductive reactance and is calculated to be 2.9943 Ω.Table 6
Unit length line parameters.
Line type
Parameters
R (Ω/km)
L (mH/km)
C (μF/km)
Overhead line
Positive sequence
0.178
1.25
0.061
Zero sequence
0.25
4.56
0.038
Cable line
Positive sequence
0.27
0.256
0.337
Zero sequence
2.7
1.039
0.278Case 1.
(low-resistance grounding). We usel3, 0.02 s, and 150 Ω to indicate the fault occurrence in line number 3 when its faulty resistance is 150 Ω, and the fault occurs at 0.02 s. Fault of line l3 occurs at the point 6 km from the bus. In order to perfectly simulate the actual situation in the field, Gaussian white noise is artificially added to the zero-sequence current in each line, and the noise intensity is −20 dBW. Zero-sequence current with noise in each line while a low-resistance grounding fault is shown in Figure 7. The zero-sequence current in each line is analyzed by Fast Fourier Transform (FFT), and the result is shown in Figure 8. Figure 8 shows that due to the compensation effect of the arc suppression coil, the fundamental component of the zero-sequence current of the fault line is not the largest component, and the zero-sequence current fundamental amplitude method for faulty line selection fails. In addition, the fifth harmonic is little affected by the arc suppression coil, but its content is less and the amplitude is small, which makes it difficult to select the faulty line. In a word, the application of FFT is limited.
The zero-sequence current of four cycles after fault is decomposed by MCEEMDAN proposed in this paper. Set the number of adding noise as 100 times, the amplitude of adding noise as 0.3, and the maximum number of iterations as 100 times. The results of MCEEMDAN decomposition in each line are shown in Figure9. It can be known from Figure 9 that the decomposed intrinsic mode functions are stable, and the noise component is suppressed effectively.
According to the optimal denoising smooth model, the optimal filtering algorithm is selected. The results of optimal denoising smooth model in each line are shown in Figure10. It can be known from Figure 10 that the optimal denoising signal is smooth and consistent with the change of the original signal. In short, not only can the optimal denoising smooth model of MCEEMDAN filter out noise but it also retains the useful information of the original signal.
When the external driving force is not added to the Duffing system, the chaotic oscillator is adjusted to the critical state. Set the damping ratio of the chaotic oscillator as 0.3; set the amplitude of internal driving force as 0.464602; set the angular frequency of internal driving force as 500π rad/s, and the calculated step size is 5 × 10−5. When the amplitude of the external driving force signal is much smaller than that of the internal driving force of the system, the θ value of formula (29) can be ignored. The external driving force needs to be multiplied by a detection factor before it is input into the Duffing system. Set the detection factor as 0.001 by the experiment many times. The fourth-order Runge–Kutta method is used to solve the equation of state of a chaotic system. It can be seen from the critical state diagram in Figure 11 that the critical state is a chaotic state.
It is first input into the Duffing system that MCEEMDAN optimal denoising result of linel4 is multiplied by the detection factor 0.001. The critical phase is 0.208π, which is calculated by the trisection symmetry phase estimation. The results of MCEEMDAN optimal denoising of other lines are moved according to the critical phase, and then they are, respectively, input into the Duffing system. The phase diagram trajectory while a low-resistance grounding fault is shown in Figure 12.
It can be known from Figure12 that the phase diagram trajectory of line l3 changes from the chaotic state to the large-scale periodic motion state, moving around the focal point (±1, 0) and the saddle point (0, 0). The phase diagram trajectories of other lines change slightly relative to the critical state, but they are still chaotic state. Therefore, it can be judged that the faulty line is line l3, which is consistent with the experimental setting.Figure 7
Zero-sequence current with noise in each line while a low-resistance grounding fault.Figure 8
Zero-sequence current amplitude-frequency diagram for different lines.Figure 9
Results of MCEEMDAN decomposition in each line. (a) Linel1 decomposition. (b) Line l2 decomposition. (c) Line l3 decomposition. (d) Line l4 decomposition.
(a)
(b)
(c)
(d)Figure 10
Results of optimal denoising in each line. (a) Linel1 optimal denoising. (b) Line l2 optimal denoising. (c) Line l3 optimal denoising. (d) Line l4 optimal denoising.
(a)
(b)
(c)
(d)Figure 11
Critical state diagram of the Duffing system.Figure 12
Phase diagram trajectory while a low-resistance grounding fault. (a) Linel1 diagram phase trajectory. (b) Line l2 diagram phase trajectory. (c) Line l3 diagram phase trajectory. (d) Line l4 diagram phase trajectory.
(a)
(b)
(c)
(d)Case 2.
(high-resistance grounding fault). We usel2, 0.02 s, and 1000 Ω to indicate the fault occurrence in line number 2 when its faulty resistance is 1000 Ω, and the fault occurs at 0.02 s. Fault of line l2 occurs at the point 2 km from the bus. In order to perfectly simulate the actual situation in the field, Gaussian white noise is artificially added to the zero-sequence current in each line, and the noise intensity is −1 dBW. Zero-sequence current with noise in each line while a high-resistance grounding fault is shown in Figure 13. It can be known from Figure 13 that the zero-sequence current signal is submerged in noise, and it is difficult to distinguish the changing trend of the signal.
In order to reflect the performance of MCEEMDAN optimal denoising, it is compared with the wavelet threshold noise reduction method, which is widely used at present. The wavelet function which uses the soft threshold method adopts db2 wavelet, and the number of decomposition layers is three.
The waveforms of different denoising methods are shown in Figure14. It can be known from Figure 14 that there is still noise in wavelet threshold denoising, but the effect of the MCEEMDAN optimal denoising in this paper is better, and the noise content is less than the result of wavelet threshold denoising, and the waveform surface is smooth. It is calculated that the correlation coefficient between MCEEMDAN optimal denoising and zero-sequence current without noise is 0.9932. There are still some problems in the wavelet threshold denoising method, such as the selection of appropriate wavelet basis function and the number of decomposition layers. In a word, the method proposed in this paper has some advantages over the wavelet threshold noise reduction method. It is first input into the Duffing system that MCEEMDAN optimal denoising of line l4 is multiplied by the detection factor 0.001. The critical phase is π/3, which is calculated by the trisection symmetry phase estimation. The final phase diagram trajectory is shown in Figure 15. It can be known from Figure 15 that the faulty line is line l2, which is consistent with the experimental setting. In summary, the proposed method is also suitable for high-resistance grounding fault.Figure 13
Zero-sequence current with noise in each line while a high-resistance grounding fault.Figure 14
Waveforms of different denoising methods. (a) Linel2 zero-sequence current with noise. (b) Line l2 zero-sequence current without noise. (c) Wavelet threshold denoising method. (d) MCEEMDAN optimal denoising.Figure 15
Phase diagram trajectory while a high-resistance grounding fault. (a) Linel1 diagram phase trajectory. (b) Line l2 diagram phase trajectory. (c) Line l3 diagram phase trajectory. (d) Line l4 diagram phase trajectory.
(a)
(b)
(c)
(d)Case 3.
(intermittent arc grounding fault). Make the intermittent arc grounding fault occur in line number 3. Make fault of linel3 occur at a point 5 km away from the bus, and the arc resistance is 10 Ω. The arc burns at 0.02 s, 0.04 s, 0.06 s, and 0.08 s, and the arc goes out at 0.03 s, 0.05 s, 0.07 s, and 0.09 s. The noise intensity is −10 dBW. Zero-sequence current with noise in each line while an intermittent arc grounding fault is shown in Figure 16.
It is first input into the Duffing system that the MCEEMDAN optimal denoising result of linel3 is multiplied by the detection factor 0.001. The critical phase is π/6, which is calculated by the trisection symmetry phase estimation. The results of MCEEMDAN optimal denoising of other lines are moved according to the critical phase, and then they are, respectively, input into the Duffing system. The phase diagram trajectory while an intermittent arc grounding fault is shown in Figure 17. It can be known from Figure 17 that the phase diagram trajectory of line l3 is the chaotic state, but the phase diagram trajectories of other lines are of the large-scale periodic motion state. Therefore, it can be judged that the faulty line is line l3, which is consistent with the experimental setting. In summary, the proposed method is also suitable for intermittent arc grounding fault.Figure 16
Zero-sequence current with noise in each line while an intermittent arc grounding fault.Figure 17
Phase diagram trajectory while an intermittent arc grounding fault. (a) Linel1 diagram phase trajectory. (b) Line l2 diagram phase trajectory. (c) Line l3 diagram phase trajectory. (d) Line l4 diagram phase trajectory.
(a)
(b)
(c)
(d)Case 4.
(fault conditions in different situations). Simulation experiments are carried out for different fault distance, fault line, and fault resistance, and the critical phase is determined by the denoising result of linel3. The fault line selection results are shown in Table 7. As can be seen in Table 7, the method proposed in this paper can accurately select fault lines under different fault conditions.Table 7
Fault line selection results as each line has different faults.
Fault situation
Linel1
Linel2
Linel3
Linel4
Judgment result
(Linel1, 20 Ω, 3 km)
Periodical
Chaotic
Chaotic
Chaotic
Linel1 fault
(Linel1, 50 Ω, 5 km)
Periodical
Chaotic
Chaotic
Chaotic
Linel1 fault
(Linel1, 200 Ω, 10 km)
Periodical
Chaotic
Chaotic
Chaotic
Linel1 fault
(Linel1, 500 Ω, 12 km)
Periodical
Chaotic
Chaotic
Chaotic
Linel1 fault
(Linel2, 1 Ω, 1 km)
Chaotic
Periodical
Chaotic
Chaotic
Linel2 fault
(Linel2, 10 Ω, 6 km)
Chaotic
Periodical
Chaotic
Chaotic
Linel2 fault
(Linel2, 100 Ω, 9 km)
Chaotic
Periodical
Chaotic
Chaotic
Linel2 fault
(Linel2, 600 Ω, 3 km)
Chaotic
Periodical
Chaotic
Chaotic
Linel2 fault
(Linel3, 8 Ω, 5 km)
Periodical
Periodical
Chaotic
Periodical
Linel3 fault
(Linel3, 10 Ω, 12 km)
Periodical
Periodical
Chaotic
Periodical
Linel3 fault
(Linel3, 100 Ω, 8 km)
Periodical
Periodical
Chaotic
Periodical
Linel3 fault
(Linel3, 200 Ω, 16 km)
Periodical
Periodical
Chaotic
Periodical
Linel3 fault
(Linel4, 6 Ω, 5 km)
Chaotic
Chaotic
Chaotic
Periodical
Linel4 fault
(Linel4, 8 Ω, 3 km)
Chaotic
Chaotic
Chaotic
Periodical
Linel4 fault
(Linel4, 50 Ω, 8 km)
Chaotic
Chaotic
Chaotic
Periodical
Linel4 fault
(Linel4, 100 Ω, 2 km)
Chaotic
Chaotic
Chaotic
Periodical
Linel4 fault
(Bus bar, 10 Ω)
Chaotic
Chaotic
Chaotic
Chaotic
Bus bar fault
(Bus bar, 100 Ω)
Chaotic
Chaotic
Chaotic
Chaotic
Bus bar fault
(Bus bar, 200 Ω)
Chaotic
Chaotic
Chaotic
Chaotic
Bus bar fault


## 7. Conclusion

In this paper, MCEEMDAN, which is based on combining GRNN, FCM, and MPE, has strong adaptability. The MCEEMDAN is more complete and orthogonal than the commonly used EMD, EEMD, and CEEMDAN in the adaptive decomposition stage. The modified method can also suppress the endpoint effect problem and mode confusing phenomenon at a certain level. Different intrinsic mode functions are obtained by using the MCEEMDAN scheme, and a series of filtering algorithms are obtained, respectively. The superior filtering algorithm is selected by the optimal denoising smooth model. Not only can the optimal denoising smooth model of MCEEMDAN filter out noise, but it also retains the useful information of the original signal. In addition, its denoising ability is better than that of the wavelet threshold denoising method. Duffing system trisection symmetry phase estimation is used to detect MCEEMDAN optimal denoising results in each line. The faulty line can be selected by observing the trajectory of the phase diagram. Compared with the traditional faulty line selection method, the result of this paper is more intuitive. A large number of experimental studies show that the proposed method can accurately select the faulty line under different fault situations. The research results verify the usability and effectiveness of the proposed method. It has a good application prospect and engineering value in an un-effectively grounded system single-phase grounding protection that this paper proposes the faulty line selection method.



---

*Source: 5761642-2020-04-06.xml*