Power distribution network dynamic harmonic suppression method based on model predictive control

Abstract

The invention discloses a power distribution network dynamic harmonic suppression method based on model predictive control, relates to the technical field of harmonic suppression, and solves the problems that only fixed subharmonic elimination or generated harmonic passive offset is performed, the dynamic change of the resonant frequency of a power distribution network is not considered, a windowed frequency band regulation and control mechanism is lacked, and the harmonic suppression efficiency is low. And the treatment effect of the harmonic waves of the power grid is influenced. A harmonic frequency band injected by the harmonic source avoids a system resonant frequency by adjusting window parameters; and combining a power distribution network equivalent impedance model to determine the system resonant frequency in real time. And a nonlinear equation set based on the SHE-PWM algorithm eliminates specific sub-harmonics, so that the treatment efficiency is improved. A PWM period is divided into a plurality of windows, a switching angle combination is independently optimized in each window, the influence of the switching angle combinations of different windows on harmonic output is visually reflected by a window correlation topological graph, a state space model is established to describe harmonic propagation, a rolling optimization problem is set, and a control sequence is generated according to a result of solving the rolling optimization problem.

Description

Technical Field

[0001] This invention belongs to the field of harmonic control, specifically a dynamic harmonic control method for distribution networks based on model predictive control. Background Technology

[0002] Harmonics in distribution networks refer to the phenomenon where the waveforms of voltage or current in a distribution network deviate from an ideal sine wave. Essentially, they are frequency components in non-sinusoidal periodic signals whose frequencies are integer multiples of the fundamental frequency, such as the 3rd, 5th, 7th…50th harmonics, corresponding to 150Hz, 250Hz, 350Hz…2500Hz respectively. They may also include non-integer "interharmonics." The main sources of harmonics are power electronic equipment and nonlinear loads, such as frequency converters, charging piles, new energy inverters, as well as electric arc furnaces and fluorescent lamps. The switching actions or nonlinear characteristics of these devices distort the current / voltage waveforms, decomposing harmonic components. Harmonics can lead to increased transformer and line losses, shortened equipment lifespan, relay protection malfunctions, and energy metering errors. In severe cases, they can resonate with the distribution network impedance, causing overvoltage / overcurrent and threatening the safe operation of the system. Therefore, a novel dynamic harmonic mitigation method for distribution networks is needed.

[0003] Traditional dynamic harmonic mitigation schemes for distribution networks only target the elimination of fixed harmonics or the passive cancellation of existing harmonics. They do not consider the dynamic changes in the resonant frequency of the distribution network, lack a windowed frequency band control mechanism, and are difficult to generate control sequences based on the results of window characteristics, thus affecting the effectiveness of harmonic mitigation in the power grid. Summary of the Invention

[0004] This invention aims to solve at least one of the technical problems existing in the prior art. To this end, this invention proposes a dynamic harmonic control method for distribution networks based on model predictive control, which is used to solve the technical problem that only fixed harmonics are eliminated or existing harmonics are passively canceled, without considering the dynamic changes of the distribution network resonant frequency and lacking a windowed frequency band control mechanism, thus affecting the effectiveness of harmonic control in the power grid.

[0005] To address the above problems, a first aspect of the present invention provides a method for dynamic harmonic mitigation in distribution networks based on model predictive control, comprising the following steps: Based on the line distribution parameters, transformer leakage reactance, and load impedance, an equivalent impedance model of the distribution network is established, and the system resonant frequency is determined. The nonlinear equations of the SHE-PWM algorithm are used to eliminate specific subharmonics, and the PWM period is divided into several windows, with a window function defined. By adjusting the window parameters, the harmonic frequency band injected by the harmonic source avoids the system resonant frequency. The harmonic energy density function is defined, and the optimization objective is set to minimize the harmonic energy within the resonant frequency band. The switching angle combination is independently optimized within each window by using the solution of the optimization objective. Based on the optimized window switching angle combination, the correlation between windows is analyzed, a window correlation topology graph is established, the degree centrality and betweenness centrality of windows are calculated based on the topology graph, and the window data influence coefficient is analyzed. Based on the influence coefficient of window data, a state-space model is established to describe harmonic propagation. A rolling optimization problem is set up, and a control sequence is generated based on the result of solving the rolling optimization problem. The equivalent impedance model of the distribution network is used to analyze the equivalent impedance of the power grid in real time to generate compensation current. The recursive least squares method is used to update the harmonic impedance online, and the compensation frequency band is adjusted in combination with the generated control sequence.

[0006] Optionally, in one example of the above aspects, an equivalent impedance model of the distribution network is established based on line distributed parameters, transformer leakage reactance, and load impedance, and the system resonant frequency is determined, including the following steps: Based on the line distribution parameters, transformer leakage reactance, and load impedance, an equivalent impedance model of the distribution network is established. Frequency-related terms are extracted from the imaginary part of the equivalent impedance Zsys(f) of the distribution network. The equivalent inductance value in the k-th mode is obtained by Taylor expansion or frequency domain fitting. The capacitive component is separated from the admittance part of the equivalent impedance Zsys(f) of the distribution network. The equivalent capacitance value in the k-th mode is determined by fitting the frequency domain response curve. Based on the obtained equivalent inductance Leq,k and equivalent capacitance Ceq,k in the k-th mode, the system resonant frequency fres,k in the k-th mode is determined by frequency domain scanning. .

[0007] Alternatively, in one example of the above aspects, specific subharmonics are eliminated by the nonlinear equations of the SHE-PWM algorithm; By solving the nonlinear equations of the SHE-PWM algorithm, the switching angle required to eliminate specific harmonics can be obtained.

[0008] Optionally, in one example of the above aspects, the PWM period is divided into several windows, and a window function is defined, including the following steps: The PWM cycle is divided into M windows, and the switching angle combination is independently optimized within each window. A window function is defined as follows: αi(m)(t) = αi(m) + Δαm(t) Where: αi(m)(t) is the switching angle of the m-th window at time t, αi(m) is the initial switching angle of the m-th window, and Δαm(t) is the dynamic adjustment term of the m-th window; The dynamic adjustment term for the m-th window is calculated using resonant frequency feedback:

[0009] Where Δαm(t) is the dynamic adjustment term of the m-th window at time t, Kpe(t) is the proportional gain parameter of the PI controller of the m-th window at time t, Ki is the integral gain parameter of the PI controller of the m-th window, e(t)=fres,k(t)- ftar(t), e(t) is the resonant frequency error at time t, fres,k(t) is the actual resonant frequency of the m-th window under the k-th mode at time t, and ftar(t) is the target frequency of the m-th window at time t.

[0010] Optionally, in one example of the above aspects, by adjusting the window parameters to ensure that the harmonic frequency band injected by the harmonic source avoids the system resonant frequency, the harmonic energy density function is defined, including the following steps: By adjusting the window parameters to ensure that the harmonic frequency band injected by the harmonic source avoids the system resonant frequency fres,k in the k-th mode, the harmonic energy density function is established as follows:

[0011] Where Eh(f) is the harmonic energy density of the m-th window, M is the total number of windows into which the PWM period is divided, tm is the initial time of the m-th window, e is the base of the natural logarithm, j is the imaginary unit, f is the frequency, vpwm(t) is the PWM voltage of the m-th window at time t, and Tv is the PWM period.

[0012] Optionally, in one example of the above aspects, the optimization objective is set to minimize the harmonic energy within the resonant frequency band. The switching angle combination is independently optimized within each window using the solution to the optimization objective, including the following steps: The optimization objective is set as follows: by optimizing the switching angle αi(m)(t) of the m-th window of the control parameter at time t, minimize the total harmonic energy of the n-th harmonic source in the frequency band near the k-th resonant frequency fres,k of the system, expressed as:

[0013] Where f is the frequency, Δf is the safety offset bandwidth of the resonant frequency, and Eh(f) is the harmonic energy density of the m-th window; By statistically analyzing the switching angle αi(m)(t) at time t corresponding to the minimum total harmonic energy of the nth harmonic source in the frequency band near the kth resonant frequency fres,k of the system, the switching angle combination is independently optimized within each window.

[0014] Optionally, in one example of the above aspects, based on the optimized combination of window opening and closing angles, the correlation between windows is analyzed, a window correlation topology graph is established, the degree centrality and betweenness centrality of the windows are calculated based on the topology graph, and the window data influence coefficient is analyzed, including the following steps: The switching angle combination of each window is normalized to obtain the normalized Xi; The difference in the average switching angle values ​​between windows is calculated to obtain the average difference value ΔXi. The maximum and minimum values ​​of the difference in switching angle values ​​between windows, MaX and MiX, are also calculated. Finally, the degree of correlation between windows is calculated.

[0015] Where ρ is the resolution coefficient, ρ=0.5; Treat windows as nodes, establish connection edges between windows with a connection degree greater than a threshold based on the connection degree between windows, and add weight coefficients to the connection edges, with the weight coefficients set to the connection degree between the corresponding windows; The degree centrality and betweenness centrality of the window are calculated based on the topological graph, and the influence coefficient of the window data is obtained by weighted averaging of the degree centrality and betweenness centrality.

[0016] Optionally, in one example of the above aspects, a state-space model describing harmonic propagation is established based on the window data influence coefficient, a rolling optimization problem is set up, and a control sequence is generated based on the result of solving the rolling optimization problem, including the following steps: Within the control period Ts, the influence coefficients of the window data within which the time falls are arranged in chronological order to form a vector of window data influence coefficients, and a state-space model is established to predict harmonic propagation:

[0017] Where x(Ts) is the state vector in control period Ts, x(Ts+1) is the state vector in control period Ts+1, A is the state transition matrix, u(Ts) is the input control vector in control period Ts, u(Ts)=[WSHE-PWM switching angle control value, APF compensation current control value]T; B is the input matrix, w(Ts) is the process noise vector, y(Ts) is the output vector in control period Ts, C is the output matrix, q(Ts) is the window data influence coefficient vector in control period Ts, and v(Ts) is the measurement noise vector in control period Ts; The "output tracking error" and "control quantity amplitude" are weighted and summed, and rolling optimization is performed to select the control sequence u(Ts)´ that minimizes the weighted summation result and apply it to the control cycle after Ts; Among them, "output tracking error" is the difference between the output vector and the output reference value, and "control magnitude" is the norm of the control sequence.

[0018] Optionally, in one example of the above aspects, the equivalent impedance of the power grid is analyzed in real time using the equivalent impedance model of the distribution network to generate a compensation current, the harmonic impedance is updated online using the recursive least squares method, and the compensation frequency band is adjusted in combination with the generated control sequence, including the following steps: The equivalent impedance model of the distribution network is used to analyze the equivalent impedance Zsys(f) of the power grid in real time, and the compensation current Iapf(f) = -Ginv(f) * Vh(f) is generated based on Zsys(f). Where Vh(f) is the actual observed voltage, Ginv(f) is the inverse transfer function, and Ginv(f) = 1 / Zsys(f); Based on the control sequence u(Ts)´, the APF compensation current control value is extracted, and the APF compensation current control value is weighted and averaged with the generated compensation current Iapf(f) to obtain the final compensation current; The equivalent impedance of the power grid is updated online using the RLS recursive least squares method; Based on the updated equivalent impedance of the power grid, the system resonant frequency fres,k and the final compensation current are updated. The final compensation current, combined with the control sequence u(Ts)´ that minimizes the weighted summation result, is applied to the control cycle after Ts.

[0019] Compared with the prior art, the beneficial effects of the present invention are: This invention adjusts window parameters to ensure that the harmonic frequency band injected by the harmonic source avoids the system's resonant frequency, reducing the likelihood of resonance triggering at the source. Compared to passive compensation-based mitigation, it fundamentally suppresses the generation of resonant overvoltage. The system resonant frequency is determined in real time using the distribution network's equivalent impedance model. The nonlinear equations based on the SHE-PWM algorithm eliminate specific harmonics, achieving precise suppression of characteristic harmonics without affecting the fundamental component, thus improving mitigation efficiency.

[0020] This invention divides the PWM cycle into several windows, independently optimizing the switching angle combination within each window. This allows for differentiated management of harmonics in different frequency bands, improving the flexibility and targeting of harmonic control. By actively adjusting the output frequency band of the harmonic source, the resonant coupling between the harmonic source and the distribution network impedance is reduced, lowering the risk of equipment failure due to resonance and improving the overall stability of the distribution network.

[0021] The window correlation topology diagram of this invention intuitively reflects the influence of different window switching angle combinations on each other's harmonic output. For example, a sudden change in the switching angle of a window with a high window data influence coefficient may trigger new harmonic distortion. Based on the window data influence coefficient, a state-space model is established to describe harmonic propagation, a rolling optimization problem is set, and a control sequence is generated based on the result of solving the rolling optimization problem to improve the quality of the PWM output waveform. Attached Figure Description

[0022] To more clearly illustrate the technical solutions in the embodiments of the present invention or the prior art, the drawings used in the description of the embodiments or the prior art will be briefly introduced below. Obviously, the drawings described below are only some embodiments of the present invention. For those skilled in the art, other drawings can be obtained based on these drawings without creative effort.

[0023] Figure 1 This is a schematic diagram of the system framework of the present invention. Detailed Implementation

[0024] The technical solution of the present invention will be clearly and completely described below with reference to the embodiments. Obviously, the described embodiments are only some embodiments of the present invention, and not all embodiments. Based on the embodiments of the present invention, all other embodiments obtained by those skilled in the art without creative effort are within the scope of protection of the present invention.

[0025] Please see Figure 1 The first aspect of this invention provides a method for dynamic harmonic mitigation in distribution networks based on model predictive control, comprising the following steps: Based on the line distribution parameters, transformer leakage reactance, and load impedance, an equivalent impedance model of the distribution network is established, and the system resonant frequency is determined. The nonlinear equations of the SHE-PWM algorithm are used to eliminate specific subharmonics, and the PWM period is divided into several windows, with a window function defined. By adjusting the window parameters, the harmonic frequency band injected by the harmonic source avoids the system resonant frequency. The harmonic energy density function is defined, and the optimization objective is set to minimize the harmonic energy within the resonant frequency band. The switching angle combination is independently optimized within each window by using the solution of the optimization objective. Based on the optimized window switching angle combination, the correlation between windows is analyzed, a window correlation topology graph is established, the degree centrality and betweenness centrality of windows are calculated based on the topology graph, and the window data influence coefficient is analyzed. Based on the influence coefficient of window data, a state-space model is established to describe harmonic propagation. A rolling optimization problem is set up, and a control sequence is generated based on the result of solving the rolling optimization problem. The equivalent impedance model of the distribution network is used to analyze the equivalent impedance of the power grid in real time to generate compensation current. The recursive least squares method is used to update the harmonic impedance online, and the compensation frequency band is adjusted in combination with the generated control sequence.

[0026] Specifically, by adjusting window parameters, the harmonic frequency band injected by the harmonic source avoids the system's resonant frequency, reducing the possibility of resonance triggering at the source. Compared to passive compensation methods, such as traditional APF which only cancels existing harmonics, this approach fundamentally suppresses the generation of resonant overvoltages. By combining the distribution network's equivalent impedance model to determine the system's resonant frequency in real time, the window parameters can adaptively adjust with changes in the resonant frequency, ensuring effective avoidance of the resonant frequency band regardless of changes in the distribution network's operating state, such as load switching or new energy integration. The nonlinear equations based on the SHE-PWM algorithm eliminate specific harmonics, achieving precise suppression of characteristic harmonics from the 2nd to the 50th order without affecting the fundamental component, resulting in high control accuracy. Dividing the PWM period into several windows and independently optimizing the switching angle combination within each window allows for differentiated control of harmonics in different frequency bands, improving the flexibility and targeting of harmonic control. By actively adjusting the output frequency band of the harmonic source, the resonant coupling between the harmonic source and the distribution network impedance is reduced, lowering the risk of equipment failures due to resonance, such as transformer overheating and capacitor damage, thus improving the overall stability of the distribution network.

[0027] The window correlation topology diagram intuitively reflects the influence of different window switching angle combinations on each other's harmonic output. For example, a sudden change in the switching angle of a window with a high window data influence coefficient may trigger new harmonic distortion.

[0028] Based on the influence coefficient of window data, a state-space model is established to describe harmonic propagation. A rolling optimization problem is set up, and a control sequence is generated based on the result of solving the rolling optimization problem to improve the quality of PWM output waveform.

[0029] The state-space model incorporates key variables such as harmonic amplitude, phase, resonant frequency, equivalent impedance of the distribution network, and window data influence coefficients into the state vector, comprehensively describing the propagation law of harmonics in the distribution network. Compared with static models, it can more accurately predict the harmonic state at future moments, providing a reliable mathematical basis for optimization decisions.

[0030] Within each control cycle, the rolling optimization problem is solved again. Based on the harmonic state at the current moment, the window correlation analysis results, and the prediction information for the next N steps, the control sequence is dynamically generated. This "prediction-optimization-execution" rolling mechanism can quickly respond to dynamic changes in the distribution network (such as resonant frequency shifts and abrupt changes in harmonic sources), avoiding the lag of traditional fixed control strategies.

[0031] In one embodiment of the present invention, an equivalent impedance model of the distribution network is established based on line distribution parameters, transformer leakage reactance, and load impedance, and the system resonant frequency is determined, including the following steps: Based on the line distribution parameters, transformer leakage reactance, and load impedance, an equivalent impedance model for the distribution network is established:

[0032] Where Zsys(f) is the equivalent impedance of the power grid, j is the imaginary unit, f is the harmonic frequency, Cshunt is the line-to-ground distributed capacitance, Rline(f) is the equivalent resistance of the distribution network transmission line, Xline(f) is the equivalent reactance of the distribution network transmission line, Cshunt is the parallel capacitance, Rload is the load equivalent resistance, Xload is the load equivalent reactance, and Lload is the load equivalent inductance. Frequency-dependent terms are extracted from the imaginary part of Zsys(f), and the equivalent inductance value in the k-th mode is obtained by Taylor expansion or frequency domain fitting. The capacitive component is separated from the admittance part of Zsys(f), and the equivalent capacitance value in the k-th mode is determined by fitting the frequency domain response curve. Based on the obtained equivalent inductance Leq,k and equivalent capacitance Ceq,k in the k-th mode, the system resonant frequency fres,k in the k-th mode is determined by frequency domain scanning. .

[0033] k represents the order of the resonant mode, corresponding to the resonant frequency of different modes in the system, such as the fundamental frequency k=1, the second harmonic k=2, etc.

[0034] In one embodiment of the present invention, specific subharmonics are eliminated by the nonlinear equations of the SHE-PWM algorithm, using the following formula:

[0035] Where n is the harmonic order, which is an integer multiple of the discrete value h, n=(h,2h,…), h is the fundamental frequency or minimum frequency interval, αi is the switching angle of the i-th switch in half a cycle, and N is the total number of switches in half a cycle. By solving this set of equations, the switching angle required to eliminate a specific harmonic can be obtained.

[0036] In one embodiment of the present invention, the PWM period is divided into several windows, and a window function is defined, including the following steps: The PWM cycle is divided into M windows, and the switching angle combination is independently optimized within each window. A window function is defined as follows: αi(m)(t) = αi(m) + Δαm(t) Where: αi(m)(t) is the switching angle of the m-th window at time t, αi(m) is the initial switching angle of the m-th window, and Δαm(t) is the dynamic adjustment term of the m-th window; The dynamic adjustment term for the m-th window is calculated using resonant frequency feedback:

[0037] Where Δαm(t) is the dynamic adjustment term of the m-th window at time t, Kpe(t) is the proportional gain parameter of the PI controller of the m-th window at time t, Ki is the integral gain parameter of the PI controller of the m-th window, e(t)=fres,k(t)- ftar(t), e(t) is the resonant frequency error at time t, fres,k(t) is the actual resonant frequency of the m-th window under the k-th mode at time t, and ftar(t) is the target frequency of the m-th window at time t.

[0038] In one embodiment of the present invention, by adjusting the window parameters to ensure that the harmonic frequency band injected by the harmonic source avoids the system resonant frequency, a harmonic energy density function is defined, including the following steps: By adjusting the window parameters to ensure that the harmonic frequency band injected by the harmonic source avoids the system resonant frequency fres,k in the k-th mode, the harmonic energy density function is established as follows:

[0039] Where Eh(f) is the harmonic energy density of the m-th window, M is the total number of windows into which the PWM period is divided, tm is the initial time of the m-th window, e is the base of the natural logarithm, j is the imaginary unit, f is the frequency, vpwm(t) is the PWM voltage of the m-th window at time t, and Tv is the PWM period.

[0040] It is a twitch factor in complex exponential form, which can be used to represent the complex form of a sinusoidal signal with frequency f. In Fourier transform, it is used to decompose the frequency components of the signal. It can also describe the phase change of the signal over time, with a phase of -2πft, which reflects the periodic phase rotation of the signal.

[0041] In one embodiment of the present invention, the optimization objective is set as minimizing the harmonic energy within the resonant frequency band. The switching angle combination is independently optimized within each window using the solution to the optimization objective, including the following steps: The optimization objective is set as follows: by optimizing the switching angle αi(m)(t) of the m-th window of the control parameter at time t, minimize the total harmonic energy of the n-th harmonic source in the frequency band near the k-th resonant frequency fres,k of the system, expressed as:

[0042] Where f is the frequency, Δf is the safety offset bandwidth of the resonant frequency, and Eh(f) is the harmonic energy density of the m-th window; By statistically analyzing the switching angle αi(m)(t) at time t corresponding to the minimum total harmonic energy of the nth harmonic source in the frequency band near the kth resonant frequency fres,k of the system, the switching angle combination is independently optimized within each window. The switching angle combination is thus the switching angle sequence.

[0043] f∈[fres,k-Δfk,fres,k+Δfk]: Defines a frequency window, indicating that the frequency f falls within the interval centered at the system's k-th resonant frequency fres,k, with a bandwidth of 2Δfk. Δfk is the safe offset bandwidth of the resonant frequency, which is usually taken as 0.2fres,k to ensure that resonance is avoided.

[0044] In one embodiment of the present invention, based on the optimized combination of window opening and closing angles, the correlation between windows is analyzed, a window correlation topology graph is established, the degree centrality and betweenness centrality of windows are calculated based on the topology graph, and the window data influence coefficient is analyzed, including the following steps: The switching angle combination of each window is normalized to obtain the normalized Xi; The difference in the average switching angle values ​​between windows is calculated to obtain the average difference value ΔXi. The maximum and minimum values ​​of the difference in switching angle values ​​between windows, MaX and MiX, are also calculated. Finally, the degree of correlation between windows is calculated.

[0045] Where ρ is the resolution coefficient, ρ=0.5; Treat windows as nodes, establish connection edges between windows with a connection degree greater than a threshold based on the connection degree between windows, and add weight coefficients to the connection edges, with the weight coefficients set to the connection degree between the corresponding windows; The degree centrality and betweenness centrality of the window are calculated based on the topological graph, and the influence coefficient of the window data is obtained by weighted averaging of the degree centrality and betweenness centrality.

[0046] In this embodiment, the degree centrality between windows is determined by the ratio of the degree of the corresponding node of the window to the total number of nodes; Betweenness centrality between windows is determined by calculating the number of shortest paths between each pair of nodes in the network, and for a specific node, calculating how many pairs of nodes have shortest paths passing through that specific node. The betweenness centrality between windows is obtained by the ratio of the calculated number of shortest paths to the number of shortest paths passing through that specific node.

[0047] In one embodiment of the present invention, a state-space model is established to describe harmonic propagation based on the window data influence coefficient, a rolling optimization problem is set, and a control sequence is generated based on the result of solving the rolling optimization problem, including the following steps: Within the control period Ts, the influence coefficients of the window data within which the time falls are arranged in chronological order to form a vector of window data influence coefficients, and a state-space model is established to predict harmonic propagation:

[0048] Here, x(Ts) is the state vector in control period Ts, a column vector containing key variables describing the internal characteristics of the system, such as harmonic amplitude, phase, and frequency. x(Ts+1) is the state vector in control period Ts+1, i.e., the system state in the next time step. A is the state transition matrix, describing the autonomous evolution of the system state between time steps, i.e., how the state changes from control period Ts to control period Ts+1 when there is no input. u(Ts) is the input control vector in control period Ts, u(Ts) = [WSHE-PWM switching angle control value, APF compensation current control value]T; B is the input matrix, describing how the input u(Ts) affects the state change. w(Ts) is the process noise vector, a random vector representing the influence of system modeling errors, external unknown disturbances, etc., on state evolution. y(Ts) is the output vector in control period Ts, which are observable quantities of the system, such as harmonic voltage and current measurements at the point of common coupling. C is the output matrix, describing how the system state x(Ts) is mapped to the output y(Ts). q(Ts) is the window data influence coefficient vector within the control period Ts. The window data influence coefficients within the control period Ts are arranged in chronological order to form the window data influence coefficient vector. v(Ts) is the measurement noise vector within the control period Ts, which is a random vector representing the influence of sensor measurement errors, etc., on the output.

[0049] Specifically, in this embodiment, the state vector x(Ts) = [harmonic current, harmonic voltage, resonant frequency]. T ; In the control period Ts, the input control vector u(Ts) = [WSHE - PWM switching angle control value, APF compensation current control value]. T ; The "output tracking error" and "control quantity amplitude" are weighted and summed, and rolling optimization is performed to select the control sequence u(Ts)´ that minimizes the weighted summation result and apply it to the control cycle after Ts; Among them, "output tracking error" is the difference between the output vector and the output reference value, and "control magnitude" is the norm of the control sequence.

[0050] In this embodiment, when the "output tracking error" and "control amplitude" are weighted and summed, the weight of the "output tracking error" is set to 0.5, and the weight of the "control amplitude" is set to 0.5.

[0051] In one embodiment of the present invention, a compensation current is generated by analyzing the equivalent impedance of the power grid in real time using the equivalent impedance model of the distribution network, and the harmonic impedance is updated online using the recursive least squares method. The compensation frequency band is then adjusted in conjunction with the generated control sequence. The process includes the following steps: The equivalent impedance model of the distribution network is used to analyze the equivalent impedance Zsys(f) of the power grid in real time, and the compensation current Iapf(f) = -Ginv(f) * Vh(f) is generated based on Zsys(f). Where Vh(f) is the actual observed voltage, Ginv(f) is the inverse transfer function, and Ginv(f) = 1 / Zsys(f); Based on the control sequence u(Ts)´, the APF compensation current control value is extracted, and the APF compensation current control value is weighted and averaged with the generated compensation current Iapf(f) to obtain the final compensation current; The equivalent impedance of the power grid is updated online using the RLS recursive least squares method:

[0052] in, This is the updated estimate of the equivalent impedance of the power grid at the (g+1)th iteration. Let P(g) be the updated estimate of the equivalent impedance of the power grid at the g-th iteration, and let P(g) be the covariance matrix used to measure the regression vector. The correlation of each component and the statistical characteristics of the estimation error. For the regression vector, =[Ih(g),Vh(g)] T ; For the regression vector The transpose vector, Vh(g), is the actual observed voltage of the power grid at the g-th iteration, and Ih(g) is the actual observed current of the power grid at the g-th iteration; Based on the updated equivalent impedance of the power grid, the system resonant frequency fres,k and the final compensation current are updated. The final compensation current, combined with the control sequence u(Ts)´ that minimizes the weighted summation result, is applied to the control cycle after Ts.

[0053] The above embodiments are only used to illustrate the technical methods of the present invention and are not intended to limit it. Although the present invention has been described in detail with reference to preferred embodiments, those skilled in the art should understand that modifications or equivalent substitutions can be made to the technical methods of the present invention without departing from the spirit and scope of the technical methods of the present invention.

Claims

1. A method for dynamic harmonic mitigation in distribution networks based on model predictive control, characterized in that, Includes the following steps: Based on the line distribution parameters, transformer leakage reactance, and load impedance, an equivalent impedance model of the distribution network is established, and the system resonant frequency is determined. The nonlinear equations of the SHE-PWM algorithm are used to eliminate specific subharmonics, and the PWM period is divided into several windows, with a window function defined. By adjusting the window parameters, the harmonic frequency band injected by the harmonic source avoids the system resonant frequency. The harmonic energy density function is defined, and the optimization objective is set to minimize the harmonic energy within the resonant frequency band. The switching angle combination is independently optimized within each window by using the solution of the optimization objective. Based on the optimized window switching angle combination, the correlation between windows is analyzed, a window correlation topology graph is established, the degree centrality and betweenness centrality of windows are calculated based on the topology graph, and the window data influence coefficient is analyzed. Based on the influence coefficient of window data, a state-space model is established to describe harmonic propagation. A rolling optimization problem is set up, and a control sequence is generated based on the result of solving the rolling optimization problem. The equivalent impedance model of the distribution network is used to analyze the equivalent impedance of the power grid in real time to generate compensation current. The recursive least squares method is used to update the harmonic impedance online, and the compensation frequency band is adjusted in combination with the generated control sequence.

2. The method for dynamic harmonic control of distribution networks based on model predictive control according to claim 1, characterized in that, Based on the line distribution parameters, transformer leakage reactance, and load impedance, an equivalent impedance model of the distribution network is established, and the system resonant frequency is determined, including the following steps: Based on the line distribution parameters, transformer leakage reactance, and load impedance, an equivalent impedance model for the distribution network is established: Where Zsys(f) is the equivalent impedance of the power grid, j is the imaginary unit, f is the harmonic frequency, Cshunt is the line-to-ground distributed capacitance, Rline(f) is the equivalent resistance of the distribution network transmission line, Xline(f) is the equivalent reactance of the distribution network transmission line, Cshunt is the parallel capacitance, Rload is the load equivalent resistance, Xload is the load equivalent reactance, and Lload is the load equivalent inductance. Frequency-dependent terms are extracted from the imaginary part of Zsys(f), and the equivalent inductance value in the k-th mode is obtained by Taylor expansion or frequency domain fitting. The capacitive component is separated from the admittance part of Zsys(f), and the equivalent capacitance value in the k-th mode is determined by fitting the frequency domain response curve. Based on the obtained equivalent inductance Leq,k and equivalent capacitance Ceq,k in the k-th mode, the system resonant frequency fres,k in the k-th mode is determined by frequency domain scanning. .

3. The method for dynamic harmonic mitigation of distribution networks based on model predictive control according to claim 1, characterized in that, Eliminating specific subharmonics using the nonlinear equations of the SHE-PWM algorithm is achieved through the following formula: Where n is the harmonic order, which is an integer multiple of the discrete value h, n=(h,2h,…), h is the fundamental frequency or minimum frequency interval, αi is the switching angle of the i-th switch in half a cycle, and N is the total number of switches in half a cycle. By solving this set of equations, the switching angle required to eliminate a specific harmonic can be obtained.

4. The method for dynamic harmonic control of distribution networks based on model predictive control according to claim 1, characterized in that, Divide the PWM period into several windows and define a window function, including the following steps: The PWM cycle is divided into M windows, and the switching angle combination is independently optimized within each window. A window function is defined as follows: αi(m)(t) = αi(m) + Δαm(t) Where: αi(m)(t) is the switching angle of the m-th window at time t, αi(m) is the initial switching angle of the m-th window, and Δαm(t) is the dynamic adjustment term of the m-th window; The dynamic adjustment term for the m-th window is calculated using resonant frequency feedback: Where Δαm(t) is the dynamic adjustment term of the m-th window at time t, Kpe(t) is the proportional gain parameter of the PI controller of the m-th window at time t, Ki is the integral gain parameter of the PI controller of the m-th window, e(t)=fres,k(t)- ftar(t), e(t) is the resonant frequency error at time t, fres,k(t) is the actual resonant frequency of the m-th window under the k-th mode at time t, and ftar(t) is the target frequency of the m-th window at time t.

5. The method for dynamic harmonic control of distribution networks based on model predictive control according to claim 1, characterized in that, By adjusting the window parameters to ensure that the harmonic frequency band injected by the harmonic source avoids the system's resonant frequency, the harmonic energy density function is defined, including the following steps: By adjusting the window parameters to ensure that the harmonic frequency band injected by the harmonic source avoids the system resonant frequency fres,k in the k-th mode, the harmonic energy density function is established as follows: Where Eh(f) is the harmonic energy density of the m-th window, M is the total number of windows into which the PWM period is divided, tm is the initial time of the m-th window, e is the base of the natural logarithm, j is the imaginary unit, f is the frequency, vpwm(t) is the PWM voltage of the m-th window at time t, and Tv is the PWM period.

6. The method for dynamic harmonic control of distribution networks based on model predictive control according to claim 5, characterized in that, The optimization objective is set as minimizing harmonic energy within the resonant frequency band. The switching angle combination is independently optimized within each window using the solution to the objective, including the following steps: The optimization objective is set as follows: by optimizing the switching angle αi(m)(t) of the m-th window of the control parameter at time t, minimize the total harmonic energy of the n-th harmonic source in the frequency band near the k-th resonant frequency fres,k of the system, expressed as: Where f is the frequency, Δf is the safety offset bandwidth of the resonant frequency, and Eh(f) is the harmonic energy density of the m-th window; By statistically analyzing the switching angle αi(m)(t) at time t corresponding to the minimum total harmonic energy of the nth harmonic source in the frequency band near the kth resonant frequency fres,k of the system, the switching angle combination is independently optimized within each window.

7. The method for dynamic harmonic mitigation of distribution networks based on model predictive control according to claim 1, characterized in that, Based on the optimized window switching angle combinations, the correlation between windows is analyzed, a window correlation topology graph is established, the degree centrality and betweenness centrality of windows are calculated based on the topology graph, and the window data influence coefficient is analyzed, including the following steps: The switching angle combination of each window is normalized to obtain the normalized Xi; The difference in the average switching angle values ​​between windows is calculated to obtain the average difference value ΔXi. The maximum and minimum values ​​of the difference in switching angle values ​​between windows, MaX and MiX, are also calculated. Finally, the degree of correlation between windows is calculated. Where ρ is the resolution coefficient, ρ=0.5; Treat windows as nodes, establish connection edges between windows with a connection degree greater than a threshold based on the connection degree between windows, and add weight coefficients to the connection edges, with the weight coefficients set to the connection degree between the corresponding windows; The degree centrality and betweenness centrality of the window are calculated based on the topological graph, and the influence coefficient of the window data is obtained by weighted averaging of the degree centrality and betweenness centrality.

8. The method for dynamic harmonic mitigation of distribution networks based on model predictive control according to claim 1, characterized in that, Based on the influence coefficients of the window data, a state-space model is established to describe harmonic propagation. A rolling optimization problem is set up, and a control sequence is generated based on the results of solving the rolling optimization problem, including the following steps: Within the control period Ts, the influence coefficients of the window data within which the time falls are arranged in chronological order to form a vector of window data influence coefficients, and a state-space model is established to predict harmonic propagation: Where x(Ts) is the state vector in control period Ts, x(Ts+1) is the state vector in control period Ts+1, A is the state transition matrix, u(Ts) is the input control vector in control period Ts, u(Ts)=[WSHE-PWM switching angle control value, APF compensation current control value]T; B is the input matrix, w(Ts) is the process noise vector, y(Ts) is the output vector in control period Ts, C is the output matrix, q(Ts) is the window data influence coefficient vector in control period Ts, and v(Ts) is the measurement noise vector in control period Ts; The "output tracking error" and "control quantity amplitude" are weighted and summed, and rolling optimization is performed to select the control sequence u(Ts)´ that minimizes the weighted summation result and apply it to the control cycle after Ts; Wherein, "output tracking error" is the difference between the output vector and the output reference value, and "control magnitude" is the norm of the control sequence.

9. The method for dynamic harmonic control of distribution networks based on model predictive control according to claim 8, characterized in that, The compensation current is generated by analyzing the equivalent impedance model of the distribution network in real time, and the harmonic impedance is updated online using the recursive least squares method. The compensation frequency band is then adjusted in conjunction with the generated control sequence. The process includes the following steps: The equivalent impedance model of the distribution network is used to analyze the equivalent impedance Zsys(f) of the power grid in real time, and the compensation current Iapf(f) = -Ginv(f) * Vh(f) is generated based on Zsys(f). Where Vh(f) is the actual observed voltage, Ginv(f) is the inverse transfer function, and Ginv(f) = 1 / Zsys(f); Based on the control sequence u(Ts)´, the APF compensation current control value is extracted, and the APF compensation current control value is weighted and averaged with the generated compensation current Iapf(f) to obtain the final compensation current; The equivalent impedance of the power grid is updated online using the RLS recursive least squares method: in, This is the updated estimate of the equivalent impedance of the power grid at the (g+1)th iteration. Let P(g) be the updated estimate of the equivalent impedance of the power grid at the g-th iteration, and let P(g) be the covariance matrix used to measure the regression vector. The correlation of each component and the statistical characteristics of the estimation error. For the regression vector, =[Ih(g),Vh(g)] T ; For the regression vector The transpose vector, Vh(g), is the actual observed voltage of the power grid at the g-th iteration, and Ih(g) is the actual observed current of the power grid at the g-th iteration; Based on the updated equivalent impedance of the power grid, the system resonant frequency fres,k and the final compensation current are updated. The final compensation current, combined with the control sequence u(Ts)´ that minimizes the weighted summation result, is applied to the control cycle after Ts.

Images