An injection-type harmonic impedance measurement method based on non-injection-type priori results

Abstract

The present invention relates to the technical field of harmonic impedance calculation, and discloses an injection-type harmonic impedance measurement method based on non-injection priori results. First, relevant data on the PCC point are collected and processed, and a relevant impedance priori model is obtained and background noise is modeled by a non-injection impedance measurement method. Thereafter, based on the non-injection impedance measurement priori results, a multi-sinusoidal disturbance signal for injection impedance measurement is designed and then injected and fitted to an impedance curve. Then, the multi-sinusoidal signal is redesigned and iteratively fitted according to the estimated error, and finally an impedance curve with required accuracy is obtained. The present invention combines the advantages of non-injection measurement method and injection measurement method. Through multiple rounds of measurement process, a non-injection method is first used to obtain a preliminary estimate of the harmonic impedance on the system side, and then an injection method is used to further measure and correct specific harmonics with high attention, effectively reducing measurement errors and improving the reliability of measurement results and the rapidity of measurement.

Description

Technical Field

[0001] The invention relates to the technical field of harmonic impedance calculation, and in particular to an injection-type harmonic impedance measurement method based on non-injection-type priori results. Background Art

[0002] In recent years, with the continuous deepening of my country's energy system reform, the access of a high proportion of power electronic devices and nonlinear loads has caused increasingly severe harmonic pollution problems. In urban distribution networks, there are multiple users with high quality requirements for power quality; at the same time, the high cable rate of urban distribution networks potentially exacerbates power quality problems. Therefore, how to make power quality meet high quality requirements is one of the core issues of common concern to both power suppliers and users.

[0003] At present, a large number of power quality monitoring devices are installed in the distribution network to monitor harmonic data with high precision. According to the monitored harmonic data, the existing calculation method usually divides the entire system from a common connection point into the system side and the user side to solve the harmonic impedance. In the problem of harmonic impedance measurement, most of the solution methods are based on the calculation of the voltage and current of the harmonics at a common connection point PCC. At present, the implementation methods of impedance measurement are divided into two categories: injection measurement method and non-injection measurement method.

[0004] The non-injection measurement method analyzes the harmonic voltage measured at the PCC point and the harmonic current measured on the public line under normal operation of the power grid. Its biggest feature is that it does not need to inject external disturbance sources, but uses the inherent harmonics or noise of its own port for analysis. This method has no effect on the normal operation of the power grid, is low-cost, can be used repeatedly, and can even perform real-time evaluation when the algorithm calculation speed is high. However, this method can only use existing harmonic and noise signals, the accuracy is generally low, and it is not possible to independently design the spectrum distribution of the disturbance signal. In addition, the existing algorithm has high requirements for the algorithm when only the harmonic voltage measured at the PCC point and the harmonic current measured on the public line are given; and it is usually only effective in specific scenarios. When the on-site working conditions change, large measurement errors may occur. Therefore, how to improve the accuracy, rapidity and universality of the algorithm is a hot issue that scholars at home and abroad continue to study.

[0005] The injection measurement method artificially injects disturbance signals into the system, making the fluctuation of the harmonic source on one side significantly higher than that on the other side, thereby obtaining the harmonic impedance. Although the injection measurement method has good measurement accuracy and is highly applicable in various scenarios, the injected disturbance signal must have a certain intensity to reflect its effect, which may affect the safe and stable operation of the system and cause the system's operating point to shift. In addition, the quantitative results obtained by this method can only reflect the harmonics on both sides at the time. Due to the above reasons and the economic issues of injecting disturbance signals, it is not realistic to repeatedly use the active injection measurement method under different working conditions at different times. In summary, the existing methods still have certain limitations in terms of harmonic impedance measurement and cannot well reflect the harmonic impedance characteristics.

[0006] In summary, the existing technology has the following problems:

[0007] 1) The non-injection impedance measurement method can only use the existing harmonic and noise signals, the accuracy is generally low, and it is impossible to independently design the spectrum distribution of the disturbance signal. In addition, the existing algorithm has high requirements when only the harmonic voltage measured at the PCC point and the harmonic current measured on the public line are given; and it usually works better only in specific scenarios. When the on-site working conditions change, large measurement errors may occur.

[0008] 2) The disturbance signal injected by the injection impedance measurement must have a certain intensity to reflect its effect, which may affect the safe and stable operation of the system and cause the system's operating point to shift. In addition, the quantitative results obtained by this method can only reflect the harmonics on both sides at that time. Due to the above reasons and the economic issues of injecting disturbance signals, it is not realistic to repeatedly use the active injection measurement method under different working conditions at different times, and it cannot well reflect the harmonic impedance characteristics. Summary of the invention

[0009] In view of the above problems, the purpose of the present invention is to provide an injection-type harmonic impedance measurement method based on non-injection-type priori results, which can effectively reduce measurement errors and improve the reliability of measurement results and the rapidity of measurement. The technical solution is as follows:

[0010] An injection-type harmonic impedance measurement method based on non-injection-type priori results comprises the following steps:

[0011] Step 1: Non-injection impedance measurement based on complex independent component method: collect and process the voltage and current data at the PCC point under the hth harmonic, obtain the impedance prior model through the non-injection impedance measurement method, and obtain the non-injection impedance measurement prior result;

[0012] Step 2: Based on the prior results of non-injection impedance measurement, a multi-sinusoidal disturbance signal model for injection impedance measurement is designed, and a multi-sinusoidal disturbance signal for injection impedance measurement is generated;

[0013] Step 3: Allocate the number of measurement points, inject the multi-sinusoidal disturbance signal generated in step 2 into the system to be measured, collect the voltage and current data at the PCC point, calculate the injection impedance measurement result, combine the non-injection impedance measurement prior results obtained in step 1, calculate the estimation error, determine the impedance of the measurement point that meets the requirements, and obtain the impedance curve through vector fitting;

[0014] Step 4: According to the estimated error, the number of measurement points is reallocated in the next round of injection impedance measurement, and then the multi-sinusoidal signal is redesigned and iteratively fitted to finally obtain an impedance curve that meets the accuracy requirements.

[0015] Further, the step 1 specifically includes:

[0016] Step 1.1: Establish complex independent component analysis mathematical model:

[0017] X=AS (1)

[0018] Where X=[X1,X2,…,X mo ] T is the observation signal matrix composed of mo observation signals; S = [S1, S2, ..., S no ] T is a source signal matrix composed of no source signals; A is a mixing matrix, representing the proportion of the source signal matrix S in the observation signal matrix X;

[0019] Step 1.2: Based on the collected voltage and current data at the PCC point under the hth harmonic, let the observation signal matrix Source signal matrix The mixing matrix A is composed of the system side harmonic impedance Z S and user-side harmonic impedance Z C The unknown matrix is ​​formed; the impedance prior model is represented by the blind source separation model, as shown in the following formula:

[0020]

[0021] Among them, sc1 and sc2 represent the scale of scaling; and They are the fast-changing components of voltage and current measured at the PCC point respectively; and They are the fast-changing components of the system-measured and user-side currents, respectively;

[0022] Step 1.3: Define the mixing coefficient and calculate the system-side harmonic impedance Z S:

[0023]

[0024] Among them, kc1 and kc2 are mixing coefficients; A 11 , A 21 , A 12 and A 22 is the element in the mixing matrix A;

[0025] From the expression of the mixing matrix A in formula (2), we can know

[0026]

[0027] Step 1.4: Obtain the impedance under odd harmonics through non-injection harmonic impedance measurement, and use the vector fitting method to perform curve fitting to obtain the prior results of non-injection impedance measurement, including impedance curve and impedance measurement set.

[0028] Furthermore, the step 2 specifically includes:

[0029] Step 2.1: Amplitude selection

[0030] The amplitudes of the harmonics of the system in the background noise at the PCC are:

[0031]

[0032] Among them, A ni (f) is the amplitude of the noise current at frequency f, A nv (f) is the amplitude of the noise voltage at frequency f; P(f) is the noise power at frequency f;

[0033] In the i-th injection impedance measurement, at frequency f k The minimum amplitude A that meets the voltage signal-to-noise ratio requirement SNRv (f k ) (i) for:

[0034] A SNRv (f k ) (i) =|SNR r A nv (f k )| (6)

[0035] Among them, A nv (f k ) is the frequency f k The amplitude of the noise voltage; SNR r is the minimum signal-to-noise ratio; the superscript (i) indicates that it is the i-th injection impedance measurement;

[0036] Similarly, at frequency f k The minimum amplitude A that meets the current signal-to-noise ratio requirement SNRi (f k ) (i) for:

[0037] A SNRi (f k ) (i) =|SNR r A ni (f k )Z(f k ) (i-1) | (7)

[0038] Among them, A ni (f k ) is the frequency f k The amplitude of the noise current, Z(f k ) (i-1) is the frequency f in the i-1th injection harmonic impedance measurement k Impedance;

[0039] At the same time, the frequency f k The disturbance signal amplitude required for the voltage and current signal-to-noise ratio is:

[0040] A(f k ) (i) =max[A SNRv (f k ) (i) ,A SNRi (f k ) (i) ] (8)

[0041] Therefore, in the discrete time domain, the time domain expression of the sinusoidal signal is:

[0042]

[0043] Where K is the number of sinusoidal components of the multi-sinusoidal signal; the amplitude, frequency and phase of the kth sinusoidal component are A(f k ) (i) 、f k and T max Indicates the signal duration, N p Indicates the number of sampling points;

[0044] Step 2.2: Frequency Selection

[0045] When the injection frequency is selected as f k The signal frequency is f ban The signals of frequency f are not injected simultaneously; ban The expression is:

[0046]

[0047] Among them, f base is the fundamental frequency;

[0048] And let the forbidden frequency set f ban for:

[0049] f ban ={f|f=f ban} (11)

[0050] The spectral resolution of the acquired discrete-time signal {x(n)} is:

[0051]

[0052] Among them, f s is the sampling rate;

[0053] The frequency f that can be selected for harmonic impedance measurement is m , that is, the discrete frequency points that can be distinguished in the spectrum are:

[0054] f m =m×Δf m=1,2,...,N p / 2 (13)

[0055] Wherein, m represents the mth sampling point;

[0056] The frequency interval [f min ,3f base ] is evenly divided into K subintervals according to the logarithm, so that f k Located in the kth subinterval, that is

[0057]

[0058] Among them, f min is the minimum value of the selected frequency;

[0059] Step 2.3: Phase and number of frequency components selection

[0060] The phase and number of each frequency component of the multi-sine signal are set to satisfy:

[0061] v peak =max(v(n))<10%V nom (15)

[0062] Among them, v peak is the maximum value of the multi-sinusoidal voltage signal v(n); V nom is the rated voltage;

[0063] The phase of each frequency component corresponding to the minimum amplitude of the disturbance signal is obtained by optimizing the calculation through the improved particle swarm algorithm. and the number of frequency components K;

[0064] The objective function is:

[0065]

[0066] Here, min|x(n) means minimizing the magnitude of x(n).

[0067] Furthermore, the optimization calculation performed by the improved particle swarm algorithm in step 2.3 is specifically as follows:

[0068] In the particle swarm algorithm, the particle swarm position update formula and velocity update formula are as follows:

[0069]

[0070] in, are the d-th dimension components of the flight velocity vector of particle p in the k-th and k+1-th iterations respectively; ω is the inertia weight; c1 and c2 are the individual and group learning factors respectively; r1 and r2 are random numbers between [0,1]; is the d-th dimension component of the optimal position found by particle p so far; is the d-th dimension component of the optimal position found for all particles so far; and are the d-th dimension components of the position vector of particle p at the k-th and k+1-th iterations respectively;

[0071] Each particle represents a possible phase solution, so each particle has K dimensions, that is, dimension d = K;

[0072] The fitness function is:

[0073] Fitness=-|x(n)| (19)

[0074] Let the inertia weight ω be:

[0075]

[0076] Among them, ω max and ω min are the maximum inertia weight and the minimum inertia weight respectively, t is the current iteration number; N is the total iteration number;

[0077] Let the learning factors c1 and c2 be:

[0078]

[0079]

[0080] Among them, c 1max and c 1min are the maximum and minimum values ​​of the learning factor c1 respectively; c 2max and c 2min are the maximum and minimum values ​​of the learning factor c2 respectively;

[0081] The phase combination represented by the particle with the highest fitness, that is, the smallest amplitude, is the optimal phase solution and frequency solution to the problem.

[0082] Furthermore, the step 3 specifically includes:

[0083] The multi-sinusoidal signal obtained in step 2 is injected into the system under test through the impedance measurement device, and the voltage v is collected at the PCC point m 、Current i m Waveform, and apply DFT to obtain the voltage and current at the frequency of the disturbance signal, that is,

[0084]

[0085] Where V m (f k ), I m (f k ) represents the frequency f k DFT results below;

[0086] By using the wave quantity method, the impedance Z corresponding to the frequency is obtained m (f k ),Right now

[0087]

[0088] The K impedance measurement results are converted into Z m = {Z m (f k )|k=1,2,...,K}, and divide it into two complementary sets, namely

[0089]

[0090] in, The frequency is f k Impedance measurement set; is the frequency divided by f k External impedance measurement set;

[0091] Will The impedance measurement set obtained in the non-injected impedance measurement Perform vector fitting on each point in the equation to obtain the fitted impedance measurement curve Z fit(s), then the estimated error E k for:

[0092]

[0093] Among them, Z fit (f k ) is the frequency f on the fitting curve k The corresponding impedance value;

[0094] The estimated error E k ≤E r The impedance of the measurement point is compared with the impedance measurement set obtained in the non-injected impedance measurement Accurate set of components Right now

[0095]

[0096] Among them, E r is the error limit; is the accurate measurement set of the last round of impedance measurement;

[0097] Accurate set Contains all measurement points that meet the error requirements in the jth round and previous rounds; finally, the accurate set The impedance curve is obtained by vector fitting of all impedances in the equation.

[0098] Furthermore, the step 4 specifically includes:

[0099] In step 3, if there is an estimation error E k >E r The next round of injection impedance measurement is at Assign more measurement points to the top; define the attention level P k And normalize it to get:

[0100]

[0101] Among them, P k ' is the normalized attention;

[0102] In the next round of allocation The number of measurement points N on k for:

[0103] N k =[P k 'K] (30)

[0104] Where K is the number of sinusoidal components of the multi-sinusoidal signal determined in step 2, and [·] represents a rounding function;

[0105] If ΣN appears k≠K, then for N k The largest N max Adjust until ΣN k =K;

[0106] After that, the measurement points of the next round meet the frequency selection conditions of step 2.2 and are logarithmically uniformly distributed in their respective intervals;

[0107] The amplitude determination method adopts the method for determining the amplitude of multiple sinusoidal signal components mentioned in step 2.1;

[0108] The phase is calculated using the particle swarm optimization algorithm mentioned in step 2.3;

[0109] Finally, the next round of multi-sinusoidal signal model is injected into the system to be tested, and the impedance curve is solved by the method of step 3; if the maximum estimation error of the fitting curve is greater than the accuracy requirement, step 4 is repeated; otherwise, the fitting is completed and the corresponding impedance curve is obtained.

[0110] The beneficial effects of the present invention are:

[0111] The present invention collects and processes relevant data on the PCC point, obtains a priori model of relevant impedance and models background noise through a non-injection impedance measurement method, and then designs a multi-sinusoidal disturbance signal for injection impedance measurement based on the prior results of non-injection impedance measurement, injects and fits the impedance curve, and then redesigns and iteratively fits the multi-sinusoidal signal according to the estimated error, and finally obtains the impedance curve with the required accuracy. Compared with simple injection and non-injection impedance measurement, it is more precise and accurate.

[0112] The present invention combines the advantages of non-injection measurement method and injection measurement method. Through multiple rounds of measurement process, the non-injection method is first used to obtain a preliminary estimate of the harmonic impedance on the system side, and then the injection method is used to further measure and correct specific harmonics with high concern. Through this injection impedance measurement method based on non-injection prior results, the measurement error can be effectively reduced, and the reliability of the measurement results and the rapidity of the measurement can be improved. BRIEF DESCRIPTION OF THE DRAWINGS

[0113] Figure 1 The flowchart of the injection-type harmonic impedance measurement method based on non-injection-type priori results of the present invention.

[0114] Figure 2 This is a flowchart of step 1 of non-injection impedance measurement based on the complex independent component method. DETAILED DESCRIPTION

[0115] The present invention is further described in detail below with reference to the accompanying drawings and specific embodiments.

[0116] The present invention proposes an injection impedance measurement based on non-injection priori results. The basic flow chart is as follows: Figure 1 As shown, it is divided into four steps S1-S4:

[0117] S1: Non-injected impedance measurement based on the complex independent component method.

[0118] Non-injection impedance measurement based on complex independent component method: collect and process the voltage and current data at the PCC point under the hth harmonic, obtain the impedance prior model through the non-injection impedance measurement method, and obtain the non-injection impedance measurement prior result.

[0119] The mathematical model of CICA (Complex Independent Component Analysis) is:

[0120] X=AS (1)

[0121] Where X=[X1,X2,…,X mo ] T is the observation signal matrix composed of mo observation signals; S = [S1, S2, ..., S no ] T is a source signal matrix composed of no source signals, and in this embodiment, mo=no; A is a mixing matrix, representing the proportion of the source signal matrix S in the observation signal matrix X.

[0122] After the observed signal is preprocessed by decentralization and whitening, it is iteratively separated. According to the central limit theorem, the Gaussianity of the linear combination of independent source signals is enhanced (that is, the non-Gaussianity is weakened). When the non-Gaussianity of the transformed signal is the strongest, the obtained signal is considered to be the source signal. Therefore, the iterative idea of ​​the CICA algorithm is to find an unmixing matrix W (w is a row vector in the unmixing matrix W) so that w T Maximize the non-Gaussian property of X. According to the Kunhn-Tucker condition, it can be transformed into an unconstrained optimization problem, and its objective function F(w) is

[0123] F(w)=E[G(w T X)]+β(||w|| 2 -1) (2)

[0124] Where β is an arbitrary constant; E[·] means finding the mean, ||·|| means finding the second norm of the vector, and G(·) represents the following nonlinear equation:

[0125] G(s)=-exp(-s 2 / 2)

[0126] At this time, Newton's method can be used to solve F(w) to obtain the unmixing matrix and act on the observation signal to obtain the source signal.

[0127] When using the CICA algorithm to estimate source signals, the most important premise is that the source signals are statistically independent of each other. and user side current As source signals, there is a certain correlation between the two, and they cannot be directly separated using the CICA algorithm. and The fast-changing components of are approximately independent, and their fast-changing components can be obtained by median filtering technology. Therefore, let The mixing matrix A is composed of the system side impedance Z S and user side Z C The unknown matrix formed can be obtained by the least squares method.

[0128] Another point to note is that the solution of blind source separation has two uncertainties: one is the uncertainty of scaling, and the other is the uncertainty of the corresponding relationship. Using the blind source separation model to represent the impedance prior model, we get:

[0129]

[0130] Among them, sc1 and sc2 represent the scale of scaling; and They are the fast-changing components of voltage and current measured at the PCC point respectively; and They are the fast-changing components of the system-measured and user-side currents respectively.

[0131] Define the mixing coefficient and calculate the system-side harmonic impedance Z S :

[0132]

[0133] From the expression of the mixing matrix A in formula (3), we can know that:

[0134]

[0135] In summary, the system side impedance Z S The calculation steps are as follows: Figure 2 shown.

[0136] After obtaining the impedance under 3-19 odd harmonics through non-injection harmonic impedance measurement, the vector fitting method is used to directly perform curve fitting to obtain the non-injection impedance measurement priori results, that is, the impedance curve and impedance measurement set

[0137] S2: Design an injected multi-sinusoidal signal model based on non-injected prior results.

[0138] Based on the prior results of non-injection impedance measurement, a multi-sinusoidal disturbance signal model for injection impedance measurement is designed, and a multi-sinusoidal disturbance signal for injection impedance measurement is generated.

[0139] The present invention uses a multi-sinusoidal signal as the injection signal. The multi-sinusoidal signal is composed of multiple single-sinusoidal signals. The frequency components, amplitude and phase of each component can be flexibly designed according to the measurement requirements, and it has excellent spectrum energy distribution design capabilities. The time domain expression of the multi-sinusoidal signal x(t) is as follows:

[0140]

[0141] Among them, a i , f i , is the amplitude, frequency and initial phase of the i-th frequency component of the multi-sine signal; is the number of sinusoidal signals contained in the multi-sine signal. When designing a multi-sine signal, these four parameters need to be optimized. The following introduces the parameter selection of the multi-sine signal.

[0142] (1) Amplitude selection

[0143] Before signal design, the background noise needs to be measured at the PCC. The noise is mainly composed of Gaussian white noise and system harmonics with random phases. The recorded background noise power spectral density PSD can be obtained by first obtaining the autocorrelation function of the collected background noise and then performing FFT (Fast Fourier Transform) on it. Therefore, the amplitude of each harmonic of the system harmonics in the background noise is:

[0144]

[0145] Taking into account the amplitude error and phase error caused by noise, the minimum signal-to-noise ratio SNR required is met r and the corresponding error limit E r (%) and E p The relationship between (°) is:

[0146]

[0147] in, is the ceiling function.

[0148] Therefore, in the i-th impedance measurement, at frequency f k The minimum amplitude A that meets the voltage signal-to-noise ratio requirement SNRv (f k ) (i) for:

[0149] A SNRv (f k ) (i) =|SNR r A nv (f k )| (9)

[0150] Among them, A nv (f k ) is the frequency f k The amplitude of the noise voltage.

[0151] Similarly, at frequency f k The minimum amplitude A that meets the current signal-to-noise ratio requirement SNRi (f k ) (i) for:

[0152] A SNRi (f k ) (i) =|SNR r A ni (f k )Z(f k ) (i-1) | (10)

[0153] Among them, A ni (f k ) is the frequency f k The amplitude of the noise current, Z(f k ) (i-1) is the frequency f in the i-1th injection harmonic impedance measurement k In the initial measurement, the corresponding value on the impedance curve obtained by non-injection impedance measurement can be used as Z(f k ) 0 The value of. At the same time, the frequency f k The disturbance signal amplitude required for the voltage and current signal-to-noise ratio is:

[0154] A(f k ) (i) =max[A SNRv (f k ) (i) ,A SNRi (f k ) (i) ] (11)

[0155] Therefore, in the discrete time domain, equation (6) can be rewritten as:

[0156]

[0157] Where K is the number of sinusoidal components of the multi-sinusoidal signal; the amplitude, frequency, and phase of the kth sinusoidal component are A(f k ) (i) 、f k , The superscript (i) indicates the i-th injection impedance measurement; T max Indicates the signal duration, represented by N p sampling points (sampling rate f s =N p / T max ). K sinusoidal components must have K phases

[0158] (2) Frequency selection

[0159] First, in the impedance measurement process, affected by the frequency coupling characteristics, when a pair of frequencies are equal and differ by 2f base About f base Symmetrical, or about 2f base When symmetrical disturbance components are injected simultaneously, signals of the same frequency will be generated, which may make impedance measurement impossible. Therefore, when the injection frequency is selected as f k The signal frequency is f ban The signals of frequency f should not be injected simultaneously. ban The expression is:

[0160]

[0161] And let f ban Forbidden frequency set f ban for:

[0162] f ban ={f|f=f ban} (14)

[0163] Therefore, when selecting the frequency f for impedance measurement k When other frequencies are selected, f ban The frequencies contained in .

[0164] Secondly, when using Discrete Fourier Transformation (DFT), if affected by the fence effect, it will cause spectrum leakage, resulting in a decrease in the accuracy of DFT calculation. The spectrum resolution of the collected discrete time signal {x(n)} is:

[0165]

[0166] The spectral resolution is the smallest frequency interval that can be resolved in the DFT, so the frequency f used for impedance measurement is selectedm Should be:

[0167] f m =m×Δf m=1,2,...,N p / 2 (16)

[0168] Finally, considering the uniformity of frequency distribution and the fact that there is no frequency point distribution in the range of 0-150 Hz on the logarithmic axis, the frequency interval [f min ,3f base ] is evenly divided into K subintervals according to the logarithm, so that f k Located in the kth subinterval (this is because non-injection impedance measurement can only obtain the harmonic impedance existing in the power grid, and there are many 3-19 odd-order harmonics in the power grid), that is:

[0169]

[0170] (3) Phase and number of components selection

[0171] Since the superposition of multiple sine waves may produce a large peak, the accuracy of impedance measurement is affected. At the same time, the phase of each component of the multi-sine signal will affect the amplitude of the multi-sine signal. Therefore, it is necessary to set the phase and number of each frequency component of the multi-sine signal to meet the following requirements:

[0172] v peak =max(v(n))<10%V nom (18)

[0173] Among them, v peak is the maximum value of the multi-sinusoidal voltage signal v(n); V nom Here, the improved particle swarm algorithm is used to optimize the calculation and obtain the phase of each frequency component corresponding to the minimum amplitude of the disturbance signal. And the number of components K.

[0174] The objective function is:

[0175]

[0176] Here, min|x(n) means minimizing the magnitude of x(n).

[0177] It can be seen that this is a dual-objective optimization problem, but it can be converted into a single-objective optimization problem by determining the number of components K in advance and gradually increasing it.

[0178] In the particle swarm algorithm, the particle swarm position update formula and velocity update formula are as follows:

[0179]

[0180] in, are the d-th dimension components of the flight velocity vector of particle p in the k-th and k+1-th iterations respectively; ω is the inertia weight; c1 and c2 are the individual and group learning factors respectively; r1 and r2 are random numbers between [0,1]; is the d-th dimension component of the optimal position found by particle p so far; is the d-th dimension component of the optimal position found for all particles so far; and are the d-th dimension components of the position vector of particle p at the k-th and k+1-th iterations, respectively. Each particle represents a possible phase solution, so each particle has K dimensions, and each dimension corresponds to a That is, dimension d = K. Since the goal is to minimize the amplitude of the synthesized multi-sine signal, the fitness function is:

[0181] Fitness=-|x(n)| (22)

[0182] In particle swarm optimization (PSO), the number of particles is a key setting that affects the diversity of the initial solution and the breadth of the solution space explored during each iteration. However, increasing the number of particles will increase the computational burden of each iteration and may cause the algorithm performance to degrade and become a simple parallel random search. Here, the number of particles is set to 5.

[0183] The inertia weight ω can significantly affect the algorithm optimization.

[0184]

[0185] Among them, ω max and ω min are the maximum inertia weight and the minimum inertia weight, which are 0.9 and 0.4 respectively; t is the current iteration number; N is the total iteration number, which is set to 30. As the number of iterations increases, the ω value will gradually decrease, thereby improving the algorithm's search ability and convergence speed for the global optimal solution, while enhancing the stability of the solution, overcoming the limitations that the algorithm may encounter during the search process, and achieving efficient and simple search.

[0186] The learning factor adjusts the position of the particle swarm in the search space and determines its optimization result. Therefore, let c1 and c2 be

[0187]

[0188] Among them, c 1max and c 1min are the maximum and minimum values ​​of the learning factor c1 respectively; c 2max and c 2minare the maximum and minimum values ​​of the learning factor c2, respectively; and in this embodiment, the maximum value of the learning factor is 1.4, and the minimum value of the learning factor is 0.6. Finally, the phase combination represented by the particle with the highest fitness (ie, the smallest amplitude) is the solution to the problem.

[0189] S3: Inject multi-sinusoidal signals and fit the curve

[0190] The multi-sinusoidal disturbance signal generated by S2 is injected into the system to be tested, and then the voltage and current data at the PCC point are collected to calculate the injection impedance measurement result. Combined with the non-injection impedance measurement prior results obtained in S1, the estimation error is calculated to determine the measurement point impedance required by the sign, and the impedance curve is obtained by vector fitting.

[0191] The multi-sinusoidal signal designed in S2 is injected into the system under test through the impedance measurement device, and the voltage v is collected at the PCC point m 、Current i m Waveform, and apply DFT to obtain the voltage and current at the frequency of the disturbance signal, that is,

[0192]

[0193] Where V m (f k ), I m (f k ) represents the frequency f k The DFT result is shown below. The impedance Z at the corresponding frequency is obtained by the wave quantity method. m (f k ),Right now

[0194]

[0195] The K impedance measurement results are converted into Z m = {Z m (f k )|k=1,2,...,K}, and divide it into two complementary sets, namely

[0196]

[0197] Will Impedance measurement results obtained in non-injected impedance measurement Perform vector fitting on each point in the equation to obtain the fitted impedance measurement curve Z fit (s), then the estimated error E k for

[0198]

[0199] Among them, Z fit (f k) is the frequency f on the fitting curve k The corresponding impedance value.

[0200] The estimated error E k ≤E r The impedance of the measurement point and the impedance obtained by non-injection in S1 form an accurate set Right now

[0201]

[0202] Contains all measurement points that meet the error requirements in the jth round and previous rounds. The impedance curve is obtained by vector fitting of all impedances in the equation.

[0203] S4: Redesign and injection fitting of multi-sinusoidal signals

[0204] According to the estimated error, the number of measurement points is reallocated in the next round of injection impedance measurement, and then the multi-sinusoidal signal is redesigned and iteratively fitted to finally obtain an impedance curve that meets the accuracy requirements.

[0205] In S3, if there is an estimation error E k >E r If the measurement point is too small, it is necessary to measure around it in the next round of injection impedance measurement. Assign more measurement points. Define the attention P k And normalize it to

[0206]

[0207] Therefore, in the next round of allocation The number of measurement points N on k for

[0208] N k =[P k 'K] (33)

[0209] Where K is the number of sinusoidal components of the multi-sinusoidal signal determined in S2, and [·] represents the rounding function. k ≠K, then for N k The largest N max Adjust until ∑N k =K.

[0210] After that, the measurement points of the next round are evenly distributed logarithmically in their respective intervals while satisfying the frequency selection conditions mentioned in S2. The amplitude determination method is the same as the method for determining the amplitude of the multi-sinusoidal signal components in S2. The phase is calculated by the optimization algorithm mentioned in S2. Finally, the determined multi-sinusoidal signal model for the next round is injected into the system, and the impedance curve solution method is the same as S3. If the maximum estimated error of the fitting curve is greater than the accuracy requirement, repeat step S4; otherwise, the fitting is completed and the corresponding impedance curve is obtained.

Claims

1. An injection-type harmonic impedance measurement method based on non-injection-type priori results, characterized in that: The following steps are involved: Step 1: Non-injection impedance measurement based on complex independent component method: collect and process the voltage and current data at the PCC point under the hth harmonic, obtain the impedance prior model through the non-injection impedance measurement method, and obtain the non-injection impedance measurement prior result; Step 2: Based on the prior results of non-injection impedance measurement, a multi-sinusoidal disturbance signal model for injection impedance measurement is designed, and a multi-sinusoidal disturbance signal for injection impedance measurement is generated; Step 3: Allocate the number of measurement points, inject the multi-sinusoidal disturbance signal generated in step 2 into the system to be measured, collect the voltage and current data at the PCC point, calculate the injection impedance measurement result, combine the non-injection impedance measurement prior results obtained in step 1, calculate the estimation error, determine the impedance of the measurement point that meets the requirements, and obtain the impedance curve through vector fitting; Step 4: According to the estimated error, the number of measurement points is reallocated in the next round of injection impedance measurement, and then the multi-sinusoidal signal is redesigned and iteratively fitted to finally obtain an impedance curve that meets the accuracy requirements.

2. The injection-type harmonic impedance measurement method based on non-injection-type priori results according to claim 1, characterized in that: The step 1 specifically includes: Step 1.1: Establish complex independent component analysis mathematical model: X=AS (1) Where X=[X1,X2,…,X mo ] T is the observation signal matrix composed of mo observation signals; S = [S1, S2, ..., S no ] T is a source signal matrix composed of no source signals; A is a mixing matrix, representing the proportion of the source signal matrix S in the observation signal matrix X; Step 1.2: Based on the collected voltage and current data at the PCC point under the hth harmonic, let the observation signal matrix Source signal matrix The mixing matrix A is composed of the system side harmonic impedance Z S and user-side harmonic impedance Z C The unknown matrix is ​​formed; the impedance prior model is represented by the blind source separation model, as shown in the following formula: Among them, sc1 and sc2 represent the scale of scaling; and They are the fast-changing components of voltage and current measured at the PCC point respectively; and They are the fast-changing components of the system-measured and user-side currents, respectively; Step 1.3: Define the mixing coefficient and calculate the system-side harmonic impedance Z S : Among them, kc1 and kc2 are mixing coefficients; A 11 , A 21 , A 12 and A 22 is the element in the mixing matrix A; From the expression of the mixing matrix A in formula (2), we can know Step 1.4: Obtain the impedance under odd harmonics through non-injection harmonic impedance measurement, and use the vector fitting method to perform curve fitting to obtain the prior results of non-injection impedance measurement, including impedance curve and impedance measurement set.

3. The injection-type harmonic impedance measurement method based on non-injection-type priori results according to claim 2, characterized in that: The step 2 specifically includes: Step 2.1: Amplitude selection The amplitudes of the harmonics of the system in the background noise at the PCC are: Among them, A ni (f) is the amplitude of the noise current at frequency f, A nv (f) is the amplitude of the noise voltage at frequency f; P(f) is the noise power at frequency f; In the i-th injection impedance measurement, at frequency f k The minimum amplitude A that meets the voltage signal-to-noise ratio requirement SNRv (f k ) (i) for: A SNRv (f k ) (i) =|SNR r A nv (f k )| (6) Among them, A nv (f k ) is the frequency f k The amplitude of the noise voltage; SNR r is the minimum signal-to-noise ratio; the superscript (i) indicates that it is the i-th injection impedance measurement; Similarly, at frequency f k The minimum amplitude A that meets the current signal-to-noise ratio requirement SNRi (f k ) (i) for: A SNRi (f k ) (i) =|SNR r A ni (f k )Z(f k ) (i-1) | (7) Among them, A ni (f k ) is the frequency f k The amplitude of the noise current, Z(f k ) (i-1) is the frequency f in the i-1th injection harmonic impedance measurement k Impedance; At the same time, the frequency f k The disturbance signal amplitude required for the voltage and current signal-to-noise ratio is: A(f k ) (i) =max[A SNRv (f k ) (i) ,A SNRi (f k ) (i) ] (8) Therefore, in the discrete time domain, the time domain expression of the sinusoidal signal is: Where K is the number of sinusoidal components of the multi-sinusoidal signal; A(f k ) (i) 、f k and are the amplitude, frequency and phase of the kth sinusoidal component respectively; T max Indicates the signal duration, N p Indicates the number of sampling points; Step 2.2: Frequency Selection When the injection frequency is selected as f k When the signal frequency is f ban The signals of frequency f are not injected simultaneously; ban The expression is: Among them, f base is the fundamental frequency; And let the forbidden frequency set f ban for: f ban ={f|f=f ban } (11) The spectral resolution of the acquired discrete-time signal {x(n)} is: Among them, f s is the sampling rate; The frequency f that can be selected for harmonic impedance measurement is m , that is, the discrete frequency points that can be distinguished in the spectrum are: f m =m×Δf m=1,2,...,N p / 2 (13) Wherein, m represents the mth sampling point; The frequency interval [f min ,3f base ] is evenly divided into K subintervals according to the logarithm, so that f k Located in the kth subinterval, that is Among them, f min is the minimum value of the selected frequency; Step 2.3: Phase and number of frequency components selection The phase and number of each frequency component of the multi-sine signal are set to satisfy: v peak =max(|v(n)|)<10%V nom (15) Among them, v peak is the maximum value of the multi-sinusoidal voltage signal v(n); V nom is the rated voltage; The phase of each frequency component corresponding to the minimum amplitude of the disturbance signal is obtained by optimizing the calculation through the improved particle swarm algorithm. and the number of frequency components K; The objective function is: Here, min|x(n)| means minimizing the magnitude of x(n).

4. The injection-type harmonic impedance measurement method based on non-injection-type priori results according to claim 3 is characterized in that: The optimization calculation performed by the improved particle swarm algorithm in step 2.3 is specifically as follows: In the particle swarm algorithm, the particle swarm position update formula and velocity update formula are as follows: in, are the d-th dimension components of the flight velocity vector of particle p in the k-th and k+1-th iterations respectively; ω is the inertia weight; c1 and c2 are the individual and group learning factors respectively; r1 and r2 are random numbers between [0,1]; is the d-th dimension component of the optimal position found by particle p so far; is the d-th dimension component of the optimal position found for all particles so far; and are the d-th dimension components of the position vector of particle p at the k-th and k+1-th iterations respectively; Each particle represents a possible phase solution, so each particle has K dimensions, that is, dimension d = K; The fitness function is: Fitness=-|x(n)| (19) Let the inertia weight ω be: Among them, ω max and ω min are the maximum inertia weight and the minimum inertia weight respectively, t is the current iteration number; N is the total iteration number; Let the learning factors c1 and c2 be: Among them, c 1max and c 1min are the maximum and minimum values ​​of the learning factor c1 respectively; c 2max and c 2min are the maximum and minimum values ​​of the learning factor c2 respectively; The phase combination represented by the particle with the highest fitness, that is, the smallest amplitude, is the optimal phase solution and frequency solution to the problem.

5. The injection-type harmonic impedance measurement method based on non-injection-type priori results according to claim 3, characterized in that: The step 3 specifically includes: The multi-sinusoidal signal obtained in step 2 is injected into the system under test through the impedance measurement device, and the voltage v is collected at the PCC point m 、Current i m Waveform, and apply DFT to obtain the voltage and current at the frequency of the disturbance signal, that is, Where V m (f k ), I m (f k ) represents the frequency f k DFT results below; By using the wave quantity method, the impedance Z corresponding to the frequency is obtained m (f k ),Right now The K impedance measurement results are converted into Z m = {Z m (f k )|k=1,2,...,K}, and divide it into two complementary sets, namely in, The frequency is f k Impedance measurement set; is the frequency divided by f k External impedance measurement set; Will The impedance measurement set obtained in the non-injected impedance measurement Perform vector fitting on each point in the equation to obtain the fitted impedance measurement curve Z fit (s), then the estimated error E k for: Among them, Z fit (f k ) is the frequency f on the fitting curve k The corresponding impedance value; The estimated error E k ≤E r The impedance of the measurement point is compared with the impedance measurement set obtained in the non-injected impedance measurement The exact set Z j acc ,Right now Among them, E r is the error limit; is the accurate measurement set of the last round of impedance measurement; Accurate set Contains all measurement points that meet the error requirements in the jth round and previous rounds; finally, the accurate set The impedance curve is obtained by vector fitting of all impedances in the equation.

6. The injection-type harmonic impedance measurement method based on non-injection-type priori results according to claim 5, characterized in that: The step 4 specifically includes: In step 3, if there is an estimation error E k >E r The next round of injection impedance measurement is at Assign more measurement points to the top; define the attention level P k And normalize it to get: Among them, P k ' is the normalized attention; In the next round of allocation The number of measurement points N on k for: N k =[P k 'K] (30) Where K is the number of sinusoidal components of the multi-sinusoidal signal determined in step 2, and [·] represents a rounding function; If ∑N appears k ≠K, then for N k The largest N max Adjust until ∑N k =K; After that, the measurement points of the next round meet the frequency selection conditions of step 2.2 and are logarithmically uniformly distributed in their respective intervals; The amplitude determination method adopts the method for determining the amplitude of multiple sinusoidal signal components mentioned in step 2.1; The phase is calculated using the particle swarm optimization algorithm mentioned in step 2.3; Finally, the next round of multi-sinusoidal disturbance signal model is injected into the system to be tested, and the impedance curve is solved by the method of step 3; if the maximum estimation error of the fitting curve is greater than the accuracy requirement, step 4 is repeated; otherwise, the fitting is completed and the corresponding impedance curve is obtained.

Images