A Harmonic Impedance Estimation Method Based on the Cauchy Mixed Model

In this paper, a new method without any tradition assumption to estimate the utility harmonic impedance of a point of common coupling (PCC) is proposed. But, the existing estimation methods usually are built on some assumptions, such as, the background harmonic is stable and small, the harmonic impedance of the customer side is much larger than that of utility side, and the harmonic sources of both sides are independent. However these assumptions are unpractical to modern power grid, which causes very wrong estimation.)e proposed method first uses a CauchyMixedModel (CMM) to express the Norton equivalent circuit of the PCC because we find that the CMM can right fit the statistical distribution of the measured harmonic data for any PCC, by testing and verifying massive measured harmonic data. Also, the parameters of the CMM are determined by the expectation maximization algorithm (EM), and then the utility harmonic impedance is estimated by means of the CMM’s parameters. Compared to the existingmethods, the main advantages of our method are as follows: it can obtain the accurate estimation results, but it is no longer dependent of any assumption and is only related to the measured data distribution. Finally, the effectiveness of the proposed method is verified by simulation and field cases.


Introduction
With the continual increase of the capacity of renewable energy and power electronic equipments in power grid, the utility harmonic problems are more and more serious, and the public requirements for power quality are more and more raised, and thus, the harmonic control becomes more significant [1,2]. e harmonic suppression requires accurately acquiring the harmonic impedance of the utility side [3][4][5][6][7].
At present, the mainstream methods for estimating the utility harmonic impedance are noninvasive, which mainly include the linear regression method [8][9][10][11][12][13], fluctuation quantity (FQ) method [14,15], independent random vector (IRV) method [16], and independent component analysis (ICA) method [17][18][19][20][21][22]. e linear regression method uses linear regression to calculate Norton equivalent circuit equation coefficient and obtain the harmonic impedance. e FQ method uses the sign characteristics of the ratio of harmonic voltage fluctuation to harmonic current fluctuation at a point of common coupling (PCC) to estimate the utility harmonic impedance. e linear regression method and FQ method can effectively estimate the utility harmonic impedance under the condition that the background harmonic is stable and small, and the customer harmonic impedance is much larger than that of utility side. e IRV method assumes the harmonic current of the PCC is weakcorrelation with the utility harmonic voltage, i.e., assumes implicitly that the harmonic impedance of customer side is much larger than that of utility side, which means that the harmonic current of the PCC is dominated by the customer. us, the covariance of the harmonic current of the PCC and the utility harmonic voltage can be approximated to zero, and then by solving the covariance equation, the utility harmonic impedance can be obtained. e main characteristic of ICA is that it can simultaneously figure out the harmonic impedances of both sides. However, ICA requires that the harmonic sources of both sides are independent.
Generally, the existing methods need one or one more assumptions, which includes that the background harmonic is stable and small, the harmonic impedance of the customer side is much larger than that of utility side, and the harmonic sources of both sides are independent. However, for modern power grid, these assumptions are hard to be held up. For instance, the harmonic sources caused by renewable energy and power electronic equipment in power grid widely distribute, the harmonic sources among power electronic equipments interact, and harmonics amplify for transmission lines and cables. us, it is obviously unreasonable that the background harmonic is assumed to be very small and stable [17]. Besides, filters or reactive compensation capacitors are usually installed on the customer side, which may cause the harmonic impedances of both sides to be approximately equivalent. As for, the independence hypothesis of the harmonic sources of both sides is much idealized, and many researchers have pointed out that the assumption may not meet the practice [4,5,17,20,22].
In this paper, a novel method without any assumption is proposed, which employs a Cauchy Mixed Model (CMM) to fit the distribution characteristic of the measured data at the PCC because we have found that the statistical distribution of the harmonic measurements of the PCC can be depicted with CMM through the analysis of the field harmonic data. So, a fitting CMM will be chosen to match the distribution relationship between the harmonic voltage, the harmonic current, and the utility harmonic impedance of the PCC. Also, then the utility harmonic impedance can be obtained by the expectation maximization algorithm (EM).
is paper mainly studies from the following four aspects.
(1) e principle of CMM is introduced, and the Norton equivalent circuit of the PCC is expressed as a CMM (2) e feasibility of CMM is verified by field measured data (3) e CMM estimation processes are derived in detail (4) e effectiveness and robustness of the proposed method are verified by comparing the estimation results of simulation and field data with existing algorithms is paper is organized as follows. Section 2 presents the analysis for the validity of CMM. Section 3 derives the process of estimating the utility impedance. Sections 4 and 5 verify the validity of the proposed method by simulation and field case. e paper is concluded in Section 6.

Harmonic Equivalent Analysis
Model. Norton equivalent circuit of the utility side and customer side is shown in Figure 1.
In Figure 1, _ I s is the utility harmonic current, and Z s is the harmonic impedance on the utility side. _ U PCC and _ I PCC are the harmonic voltage and current measured at the PCC, respectively.
According to the Norton equivalent circuit, the equation can be written as We have where, n is the location of the sampling point, and N is the number of sampling points. Δ _ I s (n) is the fluctuation of the utility harmonic current source, Δ _ I PCC (n) is the harmonic current fluctuation at the PCC, and Δ _ U PCC (n) is the harmonic voltage fluctuation at the PCC. en, equation (1) can be expressed as where us, equation (3) can be rewritten as where ΔI s−x and ΔI s−y are the real and imaginary parts of the utility harmonic current, respectively. ΔU PCC−x and ΔU PCC−y are the real and imaginary parts of the harmonic voltage measured at the PCC, respectively. ΔI PCC−x and ΔI PCC−y are the real and imaginary parts of the harmonic current measured at the PCC, respectively. Z s−x and Z s−y are the real and imaginary of the harmonic impedance on the utility side, respectively. erefore, where Customer side Utility side  Mathematical Problems in Engineering

Feasibility Analysis of the Model.
is section introduces the fitting methods of data distribution and the test of fitting results by CMM. First, we introduce the Copula function and Sklar's theorem [23,24]: Definition 1. A copula is a function C of D variables on the unit D-cube [0,1] D with the following properties: e key theorem due to Sklar clarifies the relations of dependence and the copula of a distribution: Therefore, the C is unique, that is, where To verify whether the CMM can fit the probability distribution of the harmonic voltage _ U PCC and current _ I PCC data, we combine the two-dimensional t-copula function [25] with the Monte Carlo method to generate a set of random data conforming to CMM. en, the Cumulative Distribution Functions (CDF) of field data and generated data are calculated and expressed as curves, respectively. Finally, we use the Kolmogorov-Smirnov (KS) test [26] to detect the maximum distance between two CDF curves.
The steps of generating the CMM data are as follows: (1) K models are selected, and the CMM is adopted to PCC data. en, the corresponding weight α k and covariance matrix k are obtained. (2) Convert the obtained covariance matrix k into a correlation matrix ρ k . (3) e t-Copula function is used to generate random two-dimensional Cauchy number of K groups associated with ρ k . (4) Generate N random numbers n i from 0 to max α k . If n i < α k , select the random number set of the kth twodimensional Cauchy model.
According to the abovementioned method, the CMM is used to fit the distribution of two groups' field data. e first set of field data is conducted at the site where a 150 kV bus bar feeds a 100 MW Electric Arc Furnace (EAF). e total sample time is 10 hours, and the calculated 3rd harmonic sampling values are averaged per minute, then the number of the all samples is 600. CMM random numbers are generated corresponding to the EAF data. e CDF of the data generated by CMM and the EAF data are shown in Figure 2.
The results of the KS test are shown in Table 1.
In Table 1, H � 0 indicates acceptance of the hypothesis at the significance level of 5%, P is the probability of acceptance, and KSSTAT is the value of the test statistic. It can be seen from Figure 2 and Table 1 that the CMM can fit EAF data.
The second set of field data is taken from a city power grid with multiple DC terminals, and the measured data are collected from 600 data points (20 sampled data per minute) at a 500 kV bus of the DC terminal for the 11th harmonic voltage and current. CMM random numbers are generated corresponding to the DC terminal data. e CDF of the data generated by CMM and DC terminal data are shown in Figure 3.
The results of the KS test are shown in Table 2.
The meanings of the parameters in Table 2 are the same as those in Table 1. According to Figure 3 and Table 2, it is Mathematical Problems in Engineering obvious that the CMM can fit DC terminal data even though there is a slight deviation between the CDF curve of the CMM and that of DC terminal.

Utility Harmonic Impedance Estimation
In the previous section, the feasibility of the CMM method is verified. In this section, the impedance estimation method based on CMM is described.
For the observation data , the joint probability function of the CMM under utility im- . By taking the logarithm likelihood function for P(θ), we know   It is very difficult to directly find the value of θ from equation (14). Since the EM algorithm has stable estimation ability for the mixture model, the EM algorithm is applied here.
e EM algorithm consists of two steps, namely, E-step: calculating the expectation of the log likelihood function of the complete data and M-step: obtaining a new parameter value by maximizing the value of the intermediate amount.
E-step: according to Bayesian formula, we can find the posterior probability of the nth group of samples (ΔU PCC (n), ΔI PCC , (n)) generated by the kth Cauchy model.
For equation (15), c nk is a known probability value and satisfies equation (16).
en, we obtain the following function l(θ) based on the principle of Jensen's inequality:   By adjusting the lower bound function of equation (17), the maximum value of the lower bound function can be obtained, so that equation (17) can gradually approach the optimal solution or obtain the local optimal solution as follows:

M-
Step is to calculate a new round of the iterative parameter model. We assume that So, the partial derivative of equation (18) for μ k is Let equation (20) be equal 0; we have Since equation (21) is nonlinear, it is impossible to separate μ k . erefore, in this step, equation (21) can be iterated until convergence. e partial derivative of equation (18) for k is obtained as follows: Let (22) be equal to 0, and after simplification, the update equation of k is Similarly, since (23) is also nonlinear, we cannot separate k . us, in this step, (23) can be iterated until convergence. Since Z s (k) is a negative diagonal matrix, its elements cannot be directly updated as a whole, and each element needs to be updated separately. e following transformation can be performed: where (25) e partial derivative of equation (25) for R k is We can obtain For equations (16) and (29), e μ k , k , R k , and α k pairs of parameters are updated by iterating over equations (21), (23), (27), and (29) until convergence. Since the harmonic data corresponding to each model affects the utility harmonic impedance, we use a weighted value α k for R k as follows: e steps of this method are in Algorithm 1.

Simulation Analysis
Simulations are conducted to verify the effectiveness of the proposed method in accordance of Figure 1. ey are programmed by Matlab2017a. Both the customer side and the utility side adopt a parallel model of harmonic current source and harmonic impedance. e specific parameters of the simulation model are set as follows: the amplitude of the utility harmonic impedance Z s is 10 Ω, and the phase angle is 80°(radian value is 1.3963). e amplitude of the customer harmonic impedance Z c is 10m(m � |Z c /Z s |) Ω, | · | is the module, and the phase angle is 33°(radian value is 0.5760). _ I s and _ I c are _ I s (n) � 200q · 1 + 0.5 sin(n) + λ 1 r s (n) · exp −30 · 1 + 0.1c s (n) − 0.05 π 180 j, where, λ 1 and λ 2 are the disturbance factors belonging to [0,1], within the value range of λ 1 and λ 2 , we have q ≈ I s /I c , sin(·) is a sine function, r s (n) and r c (n) are random numbers that obey the standard normal distribution, and c s (n) and c c (n) are random numbers that follow a uniform distribution. Also, the unit is Ampere. According to the abovementioned parameters, 300 harmonic samples are generated randomly, and the utility harmonic impedance is estimated by the CMM, Fluctuation Quantity (FQ), Binary Regression (BR), Independent Random Vector (IRV), and Independent Component Analysis (ICA) methods, respectively.

Influence of Background Harmonic on Utility Impedance
Estimation. To analyze the influence of background harmonic on the utility impedance estimation, we set q from 0.5 to 1.2 in steps of 0.1, and λ 1 � λ 2 � 0.5. When m � 3, the impedance of both sides is approximately equal. Here, we define the correlation coefficient ρ of the current of both sides as follows: where E(·) is the expected value of the data and * is the conjugate plural of a complex data. In the simulation, the ρ is about 0.43, and the relative errors of the harmonic impedance amplitude and the phase relative errors are shown in Figure 4, respectively. It can be seen from Figure 4, when the impedance of both sides is close, the relative errors of amplitude and phase of the five methods gradually increase with the increase of background harmonics. However, compared with the other four methods, the CMM method is the most reliable estimate of the utility harmonic impedance. When q ≥ 0.7, the FQ, BR, IRV, and ICA methods fail to estimate the utility impedance.
When m � 10, the consumer impedance is much larger than the utility impedance, the value of impedance on each side is approximately equal. e amplitude relative errors and phase relative errors of the five methods are shown in Figure 5.
According to Figure 5 and related analysis in the previous section, the five methods' estimation errors increase with the increase of amplitude ratio q. In addition, with the increase of current ratio q, the estimation accuracy of the Mathematical Problems in Engineering 7 CMM method is higher than that of other four methods. erefore, the CMM method can effectively estimate the utility harmonic impedance when the background harmonic increases.  impedance on utility impedance estimation, set q ranges from 0.5 to 1.2 and the impedance ratio of both sides ranges from 1.5 to 15.0, and for every q (its step size is 0.05.) and m (its step size is 0.25). And set the disturbance factors λ 1 � λ 2 � 0.4, the correlation coefficient ρ of current of both sides is about 0.46. e relative errors of amplitude and phase of utility impedance are shown in Figure 6.

Influence of the Continuous Change of Background
According to Figure 6, when the correlation coefficient ρ of the current amplitudes on each side is about 0.63, the relative error of the amplitude and phase angle of the CMM method is less than 12% and 10%, respectively, even under the condition that the impedances of both sides are approximately equivalent and the background harmonic is large, which verifies the robustness and effectiveness of the CMM method.

Influence of Correlation of Harmonic Currents of Both
Sides on Utility Harmonic Impedance Estimation. To further verify the effectiveness of the CMM method when the currents of both sides are not independent or are strongly correlated, we conducted another simulation. Since changing the interference of both sides of the current can change the independence of the current, we first set the disturbance factors λ 1 � λ 2 � λ and then set λ from 0.1 to 1.0 with a step size of 0.1. With the change in value of λ, the correlation coefficient ρ of the current of both sides is as shown in Table 3 (the current ratio q is fixed at 0.5, and the impedance ratio m is fixed at 5).
As can be seen from Table 3, the correlation coefficient ρ decreases with the increase of disturbance parameter λ. Also, the amplitude relative errors and phase relative errors of the five methods are shown in Figure 7.
From Figure 7, when the correlation is strong between the currents of the consumer side and the utility side, the estimation error ratio of the ICA method on the utility harmonic impedance is the highest, followed by the FQ, BR, and IRV methods. In addition, it can be seen that the performance of the CMM method is significantly more stable, which does not depend on the correlation of currents on each side. erefore, it proves that the CMM method provides reliable utility harmonic impedance regardless of the strong or weak correlation between the consumer side current and the utility side current.

Practical Cases
e feasibility and effectiveness of the proposed method are verified by simulation in Section 4. In this section, the practicability of the proposed method is verified by two practical cases under different backgrounds.

Case 1.
To verify the practicability of the proposed method, the field data of the EAF are used for verification, and the harmonic voltage and current amplitudes of the EAF are shown in Figure 8.
e estimation values of the five methods with 600 sample points are shown in Table 4.
In Table 4, the estimation results of the CMM, BR, IRV, and ICA methods are equivalent, while the amplitude estimation result of the FQ method is quite different from those of other four methods. To further verify the reliability of the proposed method in estimating the utility harmonic impedance, the abovementioned 600 data are divided into 11 subintervals. e first 10 subintervals, each of which has 54 sample points, and the last subinterval have 60 sample  points.
e five methods are used to estimate the utility impedance in per hour. e effectiveness of impedance estimation of the five methods is shown in Figure 9.
As can be seen from Figure 9, compared with the other four methods, the CMM method produces the most stable estimation results of the amplitude and phase radian of the utility harmonic impedance in each time period. e FQ method is unstable. Compared with the CMM method, the BR method still shows great fluctuations on its performances. e performance of the IRV method is roughly equivalent to the one of the CMM method, but in the 9th time period, the fluctuation of the IRV method is larger than that of the CMM method. is result proves the accuracy and robustness of the CMM method.

Case 2.
e amplitudes of harmonic voltage and current of the DC terminal data are shown in Figure 10. In the absence of load, the utility side is directly connected into    Table 5.
And then, the 11th harmonic data are divided into 10 subintervals, and each subinterval has 60 sample points. Figure 11 shows the estimation values with per subinterval.
According to the short-circuit capacity, the approximate value of the utility harmonic impedance is 79.5 Ω. It can be seen from Table 5 that the impedance amplitude estimated by the CMM method is closest to the rough estimation of short-circuit capacity. As can be seen from Figure 11(a), the estimation results of the CMM method in each data segment are more stable than those of the other four methods. e BR, IRV, and FQ methods have the same amplitude estimation fluctuation, but the amplitude estimated by the FQ method is the smallest. e fluctuation of the ICA method is relatively large. Figure 11  robust.
e phase estimation of the FQ method is the smallest. e ICA method is unstable for phase estimation.

Conclusions
In this paper, a novel method based on CMM for estimating utility harmonic impedance is proposed. In theory, CMM can approach any probability model, and in practice, it is easy to verify that field harmonic data of PCC can effectively be matched with CMM. erefore, we convert the Norton equivalent circuit model at the PCC into a CMM, and through adjusting the parameters of the CMM with the EM algorithm, the optimal estimation of utility harmonic impedance is calculated. e distinguishing feature of the method can obtain accurate and stable estimation results without any assumption conditions. Simulations and field cases verify the effectiveness of the proposed method.

Data Availability
e data used to support the findings of this study are available from the corresponding author upon request.