Power grid multi-region frequency deviation spatio-temporal differentiation machine analysis and collaborative control method

Abstract

The application relates to the technical field of estimation of maximum frequency deviation of power grid nodes and is used for solving the problem that the existing estimation method of maximum frequency deviation of power grid nodes cannot reflect the time-space difference characteristics of frequency response of an interconnected power system. First, a two-region frequency response model TAFR considering the time-space difference characteristics of frequency of power grid nodes is established, a time-domain analytical expression describing the node frequency difference mechanism is obtained by adopting the mode analysis and average system frequency model aggregation idea, and a maximum frequency deviation estimation method reflecting the time-space difference is proposed on the basis of the analytical expression. The established TAFR model considers the distribution of tie-line parameters, system inertia and damping, and equivalently processes the new energy frequency modulation parameters, compared with the existing frequency dynamic response model, the model has improved accuracy and applicability while ensuring the time-space difference.

Description

TECHNICAL FIELD

[0001] The present application relates to the technical field of estimation of maximum frequency deviation of power grid nodes, in particular to a method for analyzing and cooperatively controlling spatial and temporal differences of power grid multi-region frequency deviation. BACKGROUND

[0002] With more and more renewable energy and power electronic devices accessing the power grid, the operation safety of the power grid is increasingly prominent. As one of the key indicators in the operation of the power grid, the stability of the frequency directly affects the operation safety of the power grid. Due to the transformation of the power supply structure and the form of the power grid, the inertia distribution of the system becomes uneven, and the frequency time-space distribution characteristics become more prominent. It is more important to more accurately and quickly predict the frequency changes in different inertia level regions. Among them, the minimum point of the frequency of the nodes of the power grid as the key parameter of the frequency monitoring and the decision basis of the frequency control, the maximum frequency deviation must be accurately predicted. At present, the estimation methods of the maximum frequency deviation of the nodes of the power grid mainly include time domain simulation, model analysis and artificial intelligence algorithm.

[0003] The time-domain simulation method has high precision and detailed information, but the construction of the model is complex, the model parameters are numerous, and the calculation amount is large. When the system size increases, the simulation time also increases. The artificial intelligence algorithm is a sample-based learning estimation method, which has fast training speed and does not require specific model structure and parameters, but the accuracy of the estimation result will decrease when the network structure changes. The model analysis method is to establish a simplified model for the system, and the time-domain analytical expression of the frequency response is obtained by analyzing the model. The most commonly used models at present are the ASF model and the SFR model. The ASF model ignores the influence of the system network structure and retains the model of each group of governors. The traditional ASF model uses the detailed analytical expression of the governor model, which has high precision but is difficult to obtain and has high order. Moreover, the order of the ASF model will continuously increase with the increase of the number of groups, making the calculation difficult and increasing the calculation cost. An open-loop linear model based on the ASF model for estimating the frequency minimum point is proposed in the prior art, which equivalent each generator speed regulation system to a first-order model, reducing the system model order and improving the calculation efficiency, but the accuracy of the first-order model of the speed regulation system is insufficient. The frequency minimum point prediction FNP model based on the polynomial model of the speed regulation system is proposed on this basis, which increases the model accuracy, but due to the insufficient consideration of the response characteristics of the speed regulation system by the polynomial fitting method, the universality is relatively poor, i.e. when modeling different speed regulation systems, additional steps are needed to determine the polynomial order, otherwise the accuracy may be affected. The SFR model equivalent all generator frequency modulation controls to a reheating link speed regulation system, which converts the power grid into a single machine model, has low model order and fast calculation. However, the SFR model requires that all generator speed regulation systems must be the same, otherwise the accuracy will decrease significantly, which has poor practicability. The improved model of SFR considers the equivalent governor dynamic characteristics, but does not consider the new energy station. The universal SFR model generalizes the speed regulation system model and is suitable for frequency dynamic analysis of power grids combining thermal power generation, hydroelectric power generation and renewable energy power generation, but the unique equivalent speed regulation system of the SFR model limits the accuracy and flexibility of the model, resulting in large error in frequency dynamic analysis under large disturbance. The influence of wind power integration on system frequency is analyzed based on the SFR model, but the generators in the model are all the same steam turbine, which cannot guarantee the analysis accuracy when the system structure changes. The generalized inertia response system and analysis method based on SFR are not verified on specific models.

[0004] The model analysis method applied most in the above method is basically limited to the idea of using an overall frequency response to represent the frequency responses at all places in the whole network. However, with the continuous development of the UHV AC-DC hybrid power grid, the regional power grid is continuously expanding and gradually connected into a large interconnected power system through long-distance grid connection lines. The dynamic frequency response after the system disturbance has obvious time-space differentiation characteristics. The time-space differentiation characteristics of the frequency response refer to the frequency oscillation appearing in the dynamic frequency response. This oscillation phenomenon is usually related to the swing of a group of generators in one region of the system or a group of generators in one or more other regions of the system. The traditional SFR model and the ASF model cannot reflect the time-space differentiation characteristics of the frequency response of the interconnected power system, and the time-domain simulation is too complex. SUMMARY

[0005] The purpose of the present application is to provide a power grid multi-region frequency deviation time-space differentiation machine analysis and collaborative control method, which is used to solve the problem that the existing power grid node maximum frequency deviation estimation method cannot reflect the frequency response time-space differentiation characteristics of the interconnected power system.

[0006] The technical scheme adopted by the present application to solve its technical problems is: a power grid multi-region frequency deviation time-space differentiation machine analysis and collaborative control method, comprising the following steps:

[0007] S1, the parameters of each region are equivalent aggregated.

[0008] The inertia time constant H, the equivalent damping D and the equivalent droop coefficient 1 / R are aggregated according to formula (1), Wherein, S i represents the total rated capacity of the generators in each region, S B is the selected reference capacity, H i is the inertia time constant of each region, D i is the damping coefficient, and the value of i is 1 or 2.

[0009] The original machine high-pressure cylinder work proportion F H and the original machine reheating time constant T R are solved according to formula (2): Wherein, F Hi is the original machine high-pressure cylinder work proportion, T Ri is the original machine reheating time constant of each region, and λ i is the conversion coefficient of each region during aggregation, R i is the droop coefficient of the governor.

[0010] S2, solving the analytical solution of the two-region frequency response model TAFR

[0011] The aggregated parameters H, D and 1 / R and the parameters F H , TR Substitute the TAFR model analytical solution expression of formula (3) into In the formula, Δf1(t) is the system frequency change of region one, Δf2(t) is the system frequency change of region two, Δf COI (t) is the inertial center frequency response expression, Δf OSC1 (t) is the oscillation frequency response expression of region one, Δf OSC2 (t) is the oscillation frequency response expression of region two.

[0012] S3, differentiate each part of the analytical solution to obtain formula (4), and obtain the inertial center frequency maximum deviation time and the oscillation frequency extreme point expression;

[0013] In formula (4), ΔP d is the disturbance power, α is the amplitude coefficient of Δf COI (t), ω r is the oscillation frequency of Δf COI (t), ω r ′ is the oscillation frequency of Δf OSC1 (t), φ is the initial phase angle of Δf COI (t), φ1 is the initial phase angle of the differential expression of the inertial center frequency response Δf′ COI (t), ω n is the natural frequency of Δf COI (t), ω n ′ is the oscillation frequency of Δf OSC1 (t), is the damping ratio of Δf COI (t), is the damping ratio of Δf OSC1 (t), φ is the initial phase angle of the differential expression of the oscillation frequency response Δf′ OSC1 (t), β is the amplitude coefficient of Δf OSC1 (t).

[0014] Let the two expressions in formula (4) be zero to obtain the inertial center frequency minimum point time t COI and the oscillation frequency extreme point time t OSC expressions: k is a non-negative integer;

[0015] S4, let t COI =t OSC， , k * is the integer part, denoted as [k * ], to obtain three extreme point times of the oscillation term, k *To meet t COI = t OSC The non-negative constant of the condition.

[0016] S5, maximum frequency offset calculation

[0017] When β < 0 and [k * ] is odd, or β > 0 and [k * ] is even, the maximum frequency offset of the two regions is respectively:

[0018] When β < 0 and [k * ] is even, or β > 0 and [k * ] is odd, the maximum frequency offset of the two regions is respectively:

[0019] In the formula, Δf 1max is the maximum frequency offset of region one, Δf 2max is the maximum frequency offset of region two, Δf1(t1) is the frequency change amount of region one at t1 time, Δf1(t2) is the frequency change amount of region one at t2 time, Δf1(t3) is the frequency change amount of region one at t3 time, Δf2(t1) is the frequency change amount of region two at t1 time, Δf2(t2) is the frequency change amount of region two at t2 time, and Δf2(t3) is the frequency change amount of region two at t3 time.

[0020] Further, the frequency response model TAFR of the two regions is:

[0021] (2H1+2H2)y″1+(D′1+L′2)y′1=ω0(ΔP F1 +ΔP F2 +ΔP d1 +ΔP d2 ) (12);

[0022] In the formula, ΔP di is the disturbance power of each region, ω0 is the rated angular velocity, ΔP Fi is the first-order inertia part of the unbalanced acceleration power of each region, k is the node admittance coefficient of the system, y represents the frequency response expression and the result after differentiation or integration, y1′ is the inertia center frequency response expression, y1″ is the first derivative of y1′, y2 represents the result after integration of y2′, y2′ is the oscillation frequency response expression, and y2″ is the first derivative of y2′.

[0023] Further, the calculation formula of Δf COI (t) is:

[0024] Further, Δf OSC1 The calculation formula of (t) is: In formula (15), ω' is the synthetic coefficient of node admittance and rated angular velocity, D1' is the equivalent damping coefficient of region one, and D2' is the equivalent damping coefficient of region two.

[0025] Further, Δf OSC2 The calculation formula of (t) is: In formula (17)

[0026] γ1 and γ2 are the amplitude coefficients of Δf OSC2 (t), is the initial phase angle of Δf OSC2 (t). Further, after adding control in the fan, the change amount of output power is (16), wherein ΔP w is the active power change amount of the generator, Δf is the system frequency change amount; K1 is the droop coefficient; K2 is the inertia coefficient; when the system has active power disturbance and frequency fluctuation, there is In formula (16), ΔP d is the disturbance power, ΔP G is the active power output of the synchronous machine, ΔP M is the unbalanced acceleration power in the generator frequency modulation process; formula (16) is substituted into formula (17) to obtain In formula (19), H W is the equivalent inertia time constant of the system when the fan participates in frequency modulation, D W is the equivalent damping coefficient when the fan participates in frequency modulation; for the interconnected regional system with wind power access, the rotor motion equation of the interconnected regional SFR model is represented as In formula (20), ΔP di is the disturbance power of each region, ΔP Mi is the unbalanced acceleration power in the generator frequency modulation process of each region, ΔP Gi is the active power output of the synchronous machine of each region, Δf i is the system frequency change amount of each region, Δδ i is the rotor angle change amount of each generator; according to formula (20), the two-region frequency response model TAFR is obtained. Further, the transfer function of TAFR is In formula (21), Δδ is a column vector containing the rotor angle change amounts of each generator, s is a complex variable in Laplace transform; H is a diagonal matrix of inertia constants, D is a diagonal matrix of damping coefficients, K is a node admittance coefficient matrix, ΔP M is an unbalanced acceleration power column vector; Δfi Substituting the center frequency Δf, formula (21) is transformed into a two-degree-of-freedom vibration equation: In the formula, ΔP F is a column vector of the unbalanced acceleration power first-order inertia part; formula (24) is solved by using the mode superposition method of solving two-degree-of-freedom vibration equations in structural dynamics, and My''+Cy'+Zy=F(t)(25) is obtained, wherein, Assuming that the system is uniformly damped, D'1 / H1=D'2 / H2 is allowed, formula (25) is decoupled to obtain formulas (12) and (13); in formula (26), M is a diagonal matrix containing inertia constants of each region, C is a matrix of damping coefficients and combined coefficients of inertia constants, Z is a diagonal matrix of combined coefficients of inertia constants and node admittance, and F(t) is a column vector of unbalanced acceleration power and disturbance power.

[0027] The application has the beneficial effects that: the application aims at the minimum point estimation of the grid node frequency under the high proportion of new energy, considers the time-space difference characteristics of the frequency response in the power system, establishes a two-region frequency response model containing new energy and capable of characterizing the time-space differentiation, obtains a full response dynamic analytical solution of the model, and proposes a grid node frequency deviation estimation method under large power loss. The established TAFR model considers the distribution of tie line parameters, system inertia and damping, and equivalently processes the frequency modulation parameters of new energy, compared with the existing frequency dynamic response model, the model has improved accuracy and applicability while reflecting the time-space difference. The obtained frequency dynamic analytical expression has high accuracy, and the form of superposition of multiple frequency responses can accurately describe the frequency oscillation phenomenon between regions caused by time-space characteristics, thereby providing a train of thought for the analytical solution of the multi-region frequency response model. The proposed grid node frequency deviation estimation method has high accuracy and simple form, the time range of the maximum frequency deviation can be determined by once differentiating each part of the analytical solution, without considering the comparison of multiple solutions, thereby reducing the calculation. BRIEF DESCRIPTION OF DRAWINGS

[0028] Figure 1 It is a SFR model diagram of the single-region frequency response in the regional interconnected system with wind power access.

[0029] Figure 2 It is a two-region frequency response model TAFR diagram for characterizing the time-space differentiation of the interconnected system.

[0030] Figure 3 It is a model aggregation diagram of the application.

[0031] Figure 4 It is a diagram of the oscillation term Δf OSC1 (t) and Δf OSC2 (t).

[0032] Figure 5 is a three-machine two-area model diagram.

[0033] Figure 6 is a frequency response curve diagram of low inertia area two.

[0034] Figure 7 is a TAFR model and TAFRAS frequency response curve comparison diagram.

[0035] Figure 8 is a low frequency load shedding action state determination result comparison diagram.

[0036] Figure 9 is a multi-area system schematic diagram.

[0037] Figure 10 is a three-area system model diagram.

[0038] Figure 11 is a three-area frequency response comparison diagram before optimization.

[0039] Figure 12 is a three-area frequency response comparison diagram after optimization. DETAILED DESCRIPTION

[0040] The single-machine ASF model and SFR model cannot describe the mechanism of the time and space difference of the frequency response of the interconnected system, and the estimation of the frequency deviation of the grid node will have a large error. In the research of large interconnected power systems, for any control area, all other areas can be regarded as an equivalent area, thereby forming a two-area power system. If the inertia center of each area is used to describe the local frequency response of the area, each area can be equivalent to a generator, so the frequency dynamic change between two equivalent generators can be used to reflect the time and space difference characteristics of the frequency response in the system, and a two-area frequency response model based on the SFR model is used to solve this problem. In view of the above problem, the two-area frequency response model TAFR considering the time and space difference characteristics of the frequency of the grid node is first established, the modal analysis and the aggregation of the average system frequency model are adopted, and the time domain analytical expression describing the mechanism of the node frequency difference is obtained. Further, a maximum frequency deviation estimation method representing the time and space difference is proposed based on the time domain analytical expression. Finally, a two-area simulation model is established based on the MATLAB / Simulink simulation platform to verify the effectiveness of the proposed model and method. The following describes the time and space difference mechanism analysis and collaborative control method of the grid multi-area frequency deviation of the present application:

[0041] S1. Construct a two-area frequency response model TAFR representing the time and space difference mechanism.

[0042] S1.1 Two-area frequency response model TAFR and its simplification and decoupling.

[0043] In the conventional synchronous generator, there is a coupling relationship between the rotor speed and the system frequency, the moment of inertia. In order to make the fan simulate the frequency regulation characteristics of the synchronous generator, the frequency variation and its rate of change of the system are taken as inputs to change the output power of the fan. After adding control in the fan, the change of the output power of the fan is: Where, ΔP w is the active power variation of the generator, Δf is the system frequency variation; K1 is the droop coefficient; K2 is the inertia coefficient.

[0044] When the system has an active disturbance and the frequency fluctuates, there is: Where, H is the inertia time constant, ΔP d is the disturbance power, ΔP G is the output active power of the synchronous machine, D is the load damping coefficient, and ΔP M is the unbalanced acceleration power of the generator during frequency regulation.

[0045] Substituting equation (1) into equation (2) can obtain: In equation (4), H W is the equivalent inertia time constant of the system when the fan participates in frequency regulation, D W is the equivalent damping coefficient when the fan participates in frequency regulation.

[0046] Therefore, for a single region in the regional interconnected system with wind power access, the SFR model representing the frequency response of the region is as shown in Figure 1 Where, R is the governor droop coefficient, F H is the high-pressure cylinder work proportion of the prime mover, T R is the reheating time constant of the prime mover. The rotor motion equation of the SFR model of the interconnected region is:

[0047] According to equation (5), the two-region frequency response model TAFR representing the spatial and temporal differentiation of the interconnected system can be obtained, as shown in Figure 2 Where, H i is the inertia time constant of each region, D i is the damping coefficient of each region, R i is the droop coefficient of each region, F Hi is the high-pressure cylinder work proportion of the prime mover of each region, T Ri is the reheating time constant of the prime mover of each region, ΔP di is the disturbance power of each region, ΔP Mi is the unbalanced acceleration power of the generator during frequency regulation of each region, ΔP Gi is the output active power of the synchronous machine of each region, k ij is the node admittance matrix parameter, Δf iω0is the rated angular speed, Δδ i is the rotor angle variation of each generator, and i is 1 or 2. By Figure 2 The transfer function of TAFR is further obtained as follows: wherein Δδis a column vector containing the rotor angle variation of each generator, s is a complex variable in Laplace transform, H is a diagonal matrix of inertia constants, D is a diagonal matrix of damping coefficients, K is a matrix of nodal admittance coefficients, ΔP M is a column vector of unbalanced accelerating power. ΔP M can be converted into the sum of a proportional link and a first-order inertia, and the proportional link part can be combined with the diagonal matrix D of damping coefficients into an equivalent damping matrix D'. Δf i is replaced by the center frequency Δf, and the equation set is converted into a two-degree-of-freedom vibration equation: wherein ΔP F is a column vector of the first-order inertia part of unbalanced accelerating power.

[0048] The mode superposition method for solving the two-degree-of-freedom vibration equation in structural dynamics is used to solve equation (10), and after arrangement, My" + Cy' + Zy = F(t) (12) is obtained. In equation (12), wherein M is a diagonal matrix containing inertia constants of each region, C is a matrix of damping coefficients and combined coefficients of inertia constants, Z is a diagonal matrix of combined coefficients of inertia constants and nodal admittance, and F(t) is a column vector of unbalanced accelerating power and disturbance power.

[0049] Assuming that the system is uniformly damped, i.e., D'1 / H1 = D'2 / H2, equation (12) is decoupled to obtain the frequency response model TAFR of two regions, i.e., equations (14) and (15), and y1 and y2 are subsequently solved. F1 + ΔP F2 + ΔP d1 + ΔP d2 )(14).

[0050] wherein ΔP Fi is the first-order inertia part of unbalanced accelerating power, k is the nodal admittance coefficient of the system, y1' is the inertia center frequency response expression, y1 represents the inertia center frequency response expression and the result after differentiation or integration thereof, y1" is the first-order derivative of y1', y2' is the oscillation frequency response expression, y2 represents the result after integration of y2', and y2" is the first-order derivative of y2'.

[0051] S1.2 Analytical solution of TAFR model

[0052] The decoupled TAFR model is analyzed to obtain its complete time-domain expression. It is known from equation (14) that the calculation of y'1 is the same in form as that of the ASF model, and the idea of aggregating two-machine SFR model can be used to solve it. In a two-area system, the frequency response curves of each area are not the same, and a synchronous generator that can reflect the average frequency response characteristics of the two areas is defined, which is called the equivalent inertia center in the ASF model. In this way, the two-area model is aggregated into a low-order simplified SFR model, and it is assumed that the power disturbance is a step disturbance with an amplitude of ΔP di The time-domain expression is expressed as: Equation (16) is the analytical solution SFRAS of the SFR model, where, α is the amplitude coefficient of Δf COI (t), ω r is the oscillation frequency of Δf COI (t), ω r ′ is the oscillation frequency of Δf OSC1 (t), φ is the initial phase angle of Δf COI (t), φ1 is the initial phase angle of Δf ′ COI (t), ω n is the natural frequency of Δf COI (t), ω n ′ is the oscillation frequency of Δf OSC1 (t), is the damping ratio of Δf COI (t), is the damping ratio of Δf OSC1 (t), φ is the initial phase angle of Δf ′ OSC1 (t), β is the amplitude coefficient of Δf OSC1 (t).

[0053] The superposition theorem is used to solve y'2, ΔP F and ΔP d are solved separately, and then added to obtain y'2. First, ΔP d is solved by Laplace transform and inverse transform when it acts alone: In equation (18), where ω' is the combined coefficient of the node admittance and the rated angular velocity, D1' is the equivalent damping coefficient of area one, and D2' is the equivalent damping coefficient of area two. Then, ΔP F is inverse Laplace transformed when it acts alone: where, There are 6 eigenvalues, 4 modes in equation (21), which is difficult to solve directly. The equation has a summation part and the expression of the center frequency of the ASF model, so the idea of ASF model aggregation is adopted, as shown in equation (22), which can simplify equation (21) to 2 pairs of conjugate complex roots, 2 modes, greatly simplifying the calculation. Figure 3

[0054] Figure 3 , where S i represents the total rated capacity of the generator in each area, equation (21) is converted to equation (22): In equation (23),

[0055] The Laplace inverse transform of equation (23) is obtained to get the time domain expression: In equation (25) (26). In the equation, γ1, γ2 are Δf OSC2 (t) amplitude coefficient, Δf OSC2 (t) initial phase angle. Finally, according to Δδ = Φy, the coordinate inverse transformation can obtain the analytical expression of the frequency response of the two-area interconnected system: In the equation, Δf1(t) is the system frequency variation of area 1, Δf2(t) is the system frequency variation of area 2, Δf COI (t) is the center frequency response expression of inertia, Δf OSC1 (t) is the oscillation frequency response expression of area 1, Δf OSC2 (t) is the oscillation frequency response expression of area 2. Equation (27) is the complete analytical solution of the two-area frequency response TAFRAS. Analysis of it shows that the frequency response of the two areas has a symmetrical form, and their frequency responses are superimposed on the basis of the system center frequency response with certain oscillation frequencies. Although the oscillation frequency is relatively complex in form, it is the superposition of multiple sinusoidal functions, but the oscillation frequencies of the two areas have the same period, the oscillation directions are opposite, and the oscillation amplitudes are inversely proportional to the equivalent moment of inertia of each area.

[0056] S2 considers the estimation of the frequency deviation of the power grid node considering the time and space differences.

[0057] ​Compared to traditional frequency response models, frequency response models that can characterize spatiotemporal differences consider the influence of inter-regional power oscillations. Analytically, this is reflected as an inertial center frequency component superimposed on the frequency oscillation component, together forming the multi-regional frequency response of the interconnected system. However, due to the presence of the oscillation component, the traditional method of determining the maximum frequency offset point—that is, determining the extreme point of the expression through a single derivative—cannot be used. This is because multiple extreme points exist after a single derivative. To obtain the maximum offset point, all possible solutions within the range must be compared, resulting in excessive computational complexity.

[0058] By analyzing the parameters of each part of the total response expression, it is noted that the inertial center frequency response expression Δf COI (t) The trigonometric function part has a large period, and this formula mainly exhibits exponential characteristics. The trigonometric function part of the oscillation part has a small period. Compared with the parameters of the exponential part, it mainly exhibits a set of decaying sine curves and a set of exponential curves with decaying oscillations. Therefore, it can be considered that the maximum offset point of the synthesized response curve is determined by the oscillation part within one period before and after the maximum offset point of the inertial center frequency. That is, the time of the extreme point of the oscillation frequency near the time of the lowest point of the center frequency is used to approximate the time of the maximum frequency offset of the whole response.

[0059] First, for Δf COI (t) and Δf OSC1 (t) Taking the differential, we get equation (28): In the formula: Δf′ COI (t) is the differential expression for the inertial center frequency response, Δf′ OSC1 (t) is the differential expression for the oscillation frequency response.

[0060] like Figure 4 As shown, the frequency response time period of interest in this invention is Δf OSC2 (t) relative to Δf COI (t) and Δf OSC1 The influence of (t) is relatively small and can be ignored when solving for the maximum frequency offset time.

[0061] Within the range of normal parameter values, there are always Setting both equations in (28) to zero yields the time t at the point of lowest inertial center frequency. COI and the extreme point of oscillation frequency and time t OSC expression: k is a non-negative integer. Let t COI =t OSC We can obtain: k * After rounding down, we have [k] * At this point, we have: In formula (32), t1, t2, t3 are three extreme point times of oscillation term, k * In order to satisfy t COI < t1 < t2 < t3 OSC The non-negative constant of the condition.

[0062] In the case of β < 0, when [k * ] is odd, t1 is a maximum point, and the minimum point of region one is t1, and the minimum point of region two can be t2 or t3; when [k * ] is even, t1 is a minimum point, and the minimum point of region two is t1, and the minimum point of region one can be t2 or t3.

[0063] In the case of β > 0, when [k * ] is odd, t1 is a minimum point, and the minimum point of region two is t1, and the minimum point of region one can be t2 or t3; when [k * ] is even, t1 is a maximum point, and the minimum point of region one is t1, and the minimum point of region two can be t2 or t3.

[0064] When β < 0 and [k * ] is odd or β > 0 and [k * ] is even, the maximum frequency offset of the two regions is respectively: In the formula, Δf 1max 1 is the maximum frequency offset of region one, Δf 2max 2 is the maximum frequency offset of region two, Δf1(t1) is the frequency change amount of region one at t1 time, and Δf2(t2) is the frequency change amount of region two at t2 time. Similarly, when β < 0 and [k * ] is even or β > 0 and [k * ] is odd, the maximum frequency offset of the two regions is respectively: In the formula, Δf1(t2) is the frequency change amount of region one at t2 time, Δf1(t3) is the frequency change amount of region one at t3 time, and Δf2(t1) is the frequency change amount of region two at t1 time.

[0065] The above is the estimated value of the maximum frequency offset of the two regions in various cases.

[0066] As described above, the maximum frequency offset estimation step of the application is: first, the parameters of each region are aggregated, the inertia time constant H, the equivalent damping D, and the equivalent adjustment coefficient 1 / R are aggregated according to formula (35), S i represents the total rated capacity of the generator of each region, and S Bfor the selected reference capacity.

[0067] For the high-pressure cylinder of the prime mover H , the prime mover reheating time constant T R Solved from equation (36): In the formula, λ i is the conversion factor of each area during polymerization.

[0068] Substitute the post-polymerization parameters and area parameters into the analytical equation (27), and determine the disturbance occurrence area and size, first, the maximum offset of the Δf COI (t) part is obtained by partial differentiation, then the extreme value point time expression of the Δf OSC1 (t) part is obtained, the k* value in the extreme value point time expression (31) is determined from equation (30) by the maximum offset time of the center of inertia frequency, and three possible maximum offset times t1, t2, t3 are obtained from equation (32). According to the state of β and [k * ], the time is substituted into the estimation equation (33) or equation (34) to obtain the maximum frequency offset estimation value of each area.

[0069] So far, the frequency response analytical equation of the two-area interconnected system and the estimation method for determining the maximum frequency offset have been obtained. Through analysis of the analytical equation and the oscillation frequency curve, it can be seen that the reason for the difference in frequency response between areas is the size of the oscillation part amplitude. The size of the oscillation part amplitude is affected by two factors, one is the amplitude of the oscillation frequency itself, which is mainly related to the electrical distance between areas and the disturbance occurrence area, and the other is the inertia level of each area. If in the case of multiple areas, the inertia level difference between areas is very large, and the oscillation of each area is different, the maximum frequency offset difference between areas will increase, which will further increase the risk of frequency limit of low inertia area.

[0070] In a large power grid system with a high proportion of new energy, due to uneven inertia distribution, it is often necessary to subdivide the system according to the inertia level to more accurately describe the frequency deviation of nodes in the area. When divided into a multi-area system, each area contains synchronous generators and new energy generators, as shown in Figure 9 According to the mechanism analysis of the time and space differentiation of the frequency response of the interconnected system, a large difference in inertia level between areas will lead to a significant time and space differentiation of the frequency deviation of the grid nodes, which will seriously challenge the frequency safety of the low inertia area of the system. Therefore, by optimizing the coordination of the frequency modulation capability of new energy generators between areas to improve the frequency deviation differentiation of the grid nodes, a multi-area oriented node frequency deviation optimization model is proposed. The optimization objective of the model is to reduce the frequency response differentiation of each area and narrow the maximum frequency deviation difference of the inertia center frequency response: In the formula, Δf imaxis the maximum frequency deviation of each region coimax is the maximum frequency deviation of the inertia center.

[0071] The decision variable is the frequency modulation capability of the new energy unit in each region, and the difference between regions is reduced by coordinating the frequency modulation capability of the new energy unit in each region. The application considers the case of wind farm access, at this time the frequency modulation parameter of the wind turbine in each region is taken as the decision variable, and the value range of the frequency modulation parameter is determined according to the wind farm parameter and the running wind condition: In the formula: K pimin is the lower limit of the virtual inertia coefficient of the wind farm in each region, K pimax is the upper limit of the virtual inertia coefficient of the wind farm in each region; K dimin is the lower limit value of the frequency droop coefficient, K dimax is the upper limit value of the frequency droop coefficient.

[0072] In order to meet the requirements of the convergence of the frequency response of each region and the maximum deviation of the overall system frequency not increasing, the maximum deviation of the inertia center frequency before and after optimization needs to be constrained: Δf coiinitmax ≤Δf coioptmax (39) In the formula: Δf coiinitmax is the maximum deviation of the inertia center frequency before optimization, Δf coioptmax is the maximum deviation of the inertia center frequency after optimization, respectively.

[0073] As described above, the optimization model is established with the objective function (38) and the constraint condition (39). The model takes the minimum difference between the maximum deviation of the frequency response of each region after power disturbance and the maximum deviation of the inertia center frequency response as the optimization target, considers the overall system frequency safety and the spatial and temporal difference between regions, and can improve the problem of excessive spatial and temporal difference of the node frequency deviation of the system power grid.

[0074] In order to verify the correctness and accuracy of the time and space difference frequency deviation estimation method proposed in the application, the three-machine two-region system based on Simulink is used to compare the two-region model.

[0075] I. Accuracy verification of TAFR model representing spatial and temporal differences

[0076] As Figure 5 shown, region 1 is an equivalent synchronous generator, region 2 is an equivalent synchronous generator and an equivalent wind turbine, the rated power of the equivalent synchronous generator G1 is 500MW, the rated power of the equivalent synchronous generator G2 is 500MW, the rated power of an equivalent double-fed wind turbine G W is 200MW, the rated wind speed of the wind turbine is 12m / s, and the load shedding level of the wind turbine is 20%. The parameters of the region aggregation are shown in Table 1.

[0077] Table 1 - two regional model parameter aggregation results

[0078]

[0079] As Figure 6 shown, the frequency response curve of the low inertia region two, Table 2 records the values of the key indicators in the curve, including the lowest point of the system frequency, the time of the lowest point of the system frequency and the quasi-steady frequency. As can be seen from the curve of the detailed model, the power shortage makes the frequency response of region two produce obvious oscillation, which is particularly obvious within 5 seconds after the disturbance occurs, and the maximum frequency deviation and its time of the two regions have obvious differences. The SFR model can only describe the overall trend of the central frequency response of the system, and cannot describe the time and space difference characteristics of the system, and does not reflect the oscillation phenomenon of the frequency response, while the establishment of the two-region model can well describe the time and space difference characteristics of the regional power grid.

[0080] Table 2 - key parameters of three models

[0081]

[0082] As can be seen from Table 2, the maximum frequency deviation of the SFR model is 0.25Hz, the maximum deviation of the detailed system is 0.29Hz, and the error is 13.8%. The maximum frequency deviation of the TAFR model is 0.28Hz, which is very close to the actual value, and the error is only 3.5%. In terms of time estimation, the TAFR model is also more accurate, and the time of the maximum frequency deviation is also closer to the actual value, and the error is 3.3%, while the error of the SFR model is 5.4%. Compared with the two, it is obvious that the TAFR model has higher precision and advantages in off-line setting of unit protection action parameters and on-line frequency safety analysis.

[0083] On-line frequency safety analysis is based on the power grid operation data provided by the energy management system, and according to the current statistical synchronous machine inertia of the system, combined with the frequency modulation model parameter data of the system, the frequency safety analysis under the node power disturbance of the expected fault of the current operation mode can help the dispatching operation personnel to timely master the operation state of the power grid and take effective control measures. The model of the application needs to respectively count the inertia of each subarea, and for large interconnected power grids with obvious time and space distribution characteristics of system frequency, the system of the near area of the fault, but has a larger frequency change amplitude, such as assuming that the system uses a unified inertia center frequency calculation result is optimistic, and the model of the application can obtain more accurate and reasonable results. As Figure 6 shown, the frequency minimum point estimated by the SFR model is 49.75Hz, which is more optimistic than the actual value, while the TAFR model considering the time and space difference estimates the result as 49.72Hz, which is close to the actual value of 49.71Hz, and can better assist the frequency safety control decision of the dispatching level.

[0084] II. Verification of the accuracy of TAFRAS and the frequency deviation estimation method

[0085] The simulation results of the TAFR model are compared to verify the accuracy of the TAFRAS and the frequency deviation estimation method. The specific parameters of the two-area and SFR model are shown in Table 3, and all the parameters in Table 3 are in per unit.

[0086] Table 3 - Aggregated results of model parameters

[0087]

[0088] According to the method of step S2, the parameters of TAFRAS are calculated as follows:

[0089] Assuming that the load in the low-inertia area suddenly increases by 0.2 p.u., the frequency response curves calculated by the TAFR model and TAFRAS are shown in Figure 7 The results obtained by the estimation method are shown in Table 4.

[0090] As can be seen from Figure 7 , within 2-7 seconds after the disturbance occurs, the two curves are almost identical, and the analytical expression is very accurate in the early stage of frequency response. When the system enters quasi-steady state, the steady-state values of the two are slightly different, but they will also converge to the same after a period of time. From the above analysis, it can be seen that the accuracy of the analytical solution is very high.

[0091] Table 4 gives the five key indicators estimated by the analytical expression, and compares them with the simulation results of the TAFR model. It can be seen that the maximum frequency deviation and its time estimated by the analytical expression have very small errors compared with the results obtained by the TAFR model, so this method can completely replace SFRAS to estimate the maximum frequency deviation after the disturbance occurs.

[0092] Table 4 - Error analysis of TAFR model and TAFRAS

[0093]

[0094] As an important defense line in power systems, low-frequency load shedding should be able to meet the serious situation of 0-40% power shortage. In the process of low-frequency load shedding, in order to fully utilize the rotating reserve capacity of the system, the frequency starting value should not exceed 49.25 Hz, the frequency setting interval should be 0.2 Hz, and the basic wheel number can be 3-8 wheels.

[0095] ​The traditional low-frequency load shedding scheme generally determines the load shedding scheme by simulating the maximum power shortage that can occur in the actual system in advance to ensure the recovery requirements of the system frequency after various accidents. Due to the time and space distribution characteristics of the system frequency, the estimated results of the frequency deviation of the two models are quite different, and there is a large difference between the estimated results of the load shedding action based on the model and the actual load shedding action, which can cause misjudgment of the decision-making at the dispatching level. Figure 8 The comparison between the actual action and the estimated action when the system has a 0.5 p.u. power disturbance is shown. When the first round starting value of the low-frequency load shedding is 49.0 Hz, it can be seen that the action results are different under the estimated results of the two methods, and in the actual situation, the first round load shedding has started to act at 0.84 s, while the SFRAS estimated result shows that the frequency value is 49.16 Hz at this time, and the load shedding action is not triggered, and falls to the trigger value at 1.11 s. The TAFRAS estimated result shows that the frequency has fallen to 48.94 Hz at 0.84 s, the frequency has exceeded the limit, and the first round load shedding has acted, which is consistent with the actual action.

[0096] The above comparison results show that the model of the application and the traditional SFR have significant advantages in model structure and simulation estimation accuracy. In terms of modeling method, the SFR model uses the system inertia center frequency response instead of the frequency response of the grid node, and has poor applicability to systems with uneven inertia distribution; the TAFR model considers the time and space characteristics of inertia distribution, damping difference and tie line parameters, and can accurately represent the time and space difference of the frequency response of the grid node, and has good applicability to various types of systems. At the same time, the frequency response dynamic analytical solution based on the model of the application and the frequency deviation estimation method accurately highlight the frequency oscillation phenomenon caused by the time and space difference, can accurately estimate the key indicators of the frequency response of each region, and through further analysis of the obtained full response frequency dynamic analytical expression, a universal conclusion suitable for systems with significant time and space difference can be obtained to guide the decision-making of frequency control.

[0097] Three, optimization of the time and space difference of the frequency deviation of three regions.

[0098] Taking a three-region system as an example, the particle swarm algorithm is used to analyze the node frequency deviation optimization model, and the three-region model is as shown in Figure 10 Each region is provided with one synchronous generator and one wind turbine, and the aggregated parameters of each region are shown in Table 5.

[0099] Table 5-Three-region model parameter aggregation results

[0100]

[0101] Figure 10In the middle, region 1 is the left part of bus 7, the inertia level is the highest; region 2 is the right part of bus 8, the inertia level is the second; region 3 is the lower part of bus 9, the inertia level is the lowest.

[0102] Table 6 - Comparison of maximum frequency deviation difference results before and after optimization

[0103]

[0104] When the load of 80 MW is suddenly increased at bus 9, the results before and after optimization are shown in Table 6, where f1 is the difference between the maximum frequency response of region 1 and the maximum frequency response of the inertia center, f2 is the difference between the maximum frequency response of region 2 and the maximum frequency response of the inertia center, and f3 is the difference between the maximum frequency response of region 3 and the maximum frequency response of the inertia center. The initial frequency response is shown in Figure 11 It can be seen that the frequency response of the three regions is very different from the inertia center frequency response before 4s, and the frequency oscillation of region 3 is the most severe, with a maximum frequency deviation of 0.042 Hz from the inertia center frequency, and the inertia level of region 1 is high, with a smaller difference between the response trend and the inertia center frequency, with a maximum frequency deviation of 0.014 Hz from the inertia center frequency. The oscillation amplitude of region 2 is smaller than that of region 3 at the beginning of the response, but the difference between the minimum point and the inertia center frequency response is still 0.042 Hz, the overall frequency drop of the system is close to 0.6 Hz, and the maximum frequency deviation time is only about 2s.

[0105] The optimized frequency response is shown in Figure 12 The overall maximum frequency deviation of the system is improved by 0.1 Hz, and the maximum frequency deviation time is delayed to 2.8s. Due to the high inertia level of region 1, the optimized response curve is very close to the inertia center frequency response curve, with only slight oscillation at the maximum frequency deviation, and the difference between the optimized maximum frequency deviation and the inertia center frequency response is only 0.007 Hz, which is 50% smaller than before optimization. The frequency oscillation of region 3, which is a low inertia region, is significantly improved, and the difference between the optimized maximum frequency deviation and the inertia center frequency response is 0.016 Hz, which is 62% smaller than before optimization. The inertia level of region 2 is between that of region 1 and region 3, and the difference is reduced to 0.018 Hz, which is 57% smaller than before optimization. From the above data, it can be seen that the difference in inertia level between regions is large, and the time and space difference in frequency deviation of nodes is more obvious, and the frequency oscillation of low inertia region is more severe.

[0106] The above example shows that the frequency response dynamic analytical solution and the frequency deviation estimation method based on the model of the application accurately depict the frequency oscillation phenomenon caused by the time and space differences, can accurately estimate the key indicators of the frequency response of each region, and through further analysis of the obtained full response frequency dynamic analytical expression, a universal conclusion suitable for systems with significant time and space differences can be obtained, guiding the formulation of frequency control decisions, and the optimization model can better improve the problem of excessive time and space differences of multi-region frequency deviation caused by uneven inertia distribution.

[0107] The application is aimed at the problem of estimating the minimum point of the grid node frequency under high proportion of new energy, considers the time and space difference characteristics of the frequency response in the power system, establishes a two-region frequency response model of new energy containing time and space difference characteristics, obtains the full response dynamic analytical solution of the model, and proposes a grid node frequency deviation estimation method under large power loss. The established TAFR model considers the distribution of tie-line parameters, system inertia and damping, and equivalent new energy frequency modulation parameters, compared with the existing frequency dynamic response model, the model improves the accuracy and applicability while ensuring the time and space difference. The obtained frequency dynamic analytical expression has high accuracy, adopts the form of superposition of multiple frequency responses, and can accurately describe the frequency oscillation phenomenon between regions caused by time and space characteristics, providing a train of thought for the analytical solution of the multi-region frequency response model. The proposed grid node frequency deviation estimation method is accurate and simple in form, through the first-order differentiation of each part of the analytical solution, the time range of the maximum frequency deviation can be determined, without considering the comparison problem of multiple solutions, and the calculation is reduced.

Claims

1. A method for analyzing and coordinating control of power grid multi-region frequency deviation, characterized in that, The method comprises the following steps: S1, equivalent aggregation of parameters of each region The inertia time constant H, the equivalent damping D, and the equivalent regulation coefficient 1 / R are aggregated according to formula (1), where S i represents the total rated capacity of the generators in each area, S B is a selected reference capacity, H i is the inertia time constant of each area, D i is the damping coefficient, and i is 1 or 2; Solve the proportion of the work of the high-pressure cylinder of the prime mover F according to formula (2) H , the reheating time constant of the prime mover T R : Wherein, F Hi is the proportion of the work of the high-pressure cylinder of the prime mover in each area, T Ri is the reheating time constant of the prime mover in each area, R i is the regulating coefficient of the governor, λ i is the conversion coefficient of each area during polymerization; S2, solving an analytical solution of a two-region frequency response model TAFR The post-polymerization parameters H, D, 1 / R and the area parameters F H , T R are substituted into the TAFR model analytical solution expression of formula (3), In the formula, Δf1(t) is the system frequency variation of area one, Δf2(t) is the system frequency variation of area two, Δf COI (t) is the inertia center frequency response expression, Δf OSC1 (t) is the oscillation frequency response expression of area one, Δf OSC2 (t) is the oscillation frequency response expression of area two; S3, differentiating each part of the analytical solution to obtain formula (4), and solving time expressions of a maximum frequency deviation of the center of inertia and an extreme point of an oscillation frequency; In equation (4), ΔP d is the disturbance power, α is the amplitude coefficient of COI (t), ω r is the oscillation frequency of COI (t), ω r ' is the oscillation frequency of OSC1 (t), φ is the initial phase angle of COI (t), φ1 is the initial phase angle of COI (t), ω n is the natural frequency of COI (t), ω n ' is the oscillation frequency of OSC1 (t), is the damping ratio of COI (t), is the damping ratio of OSC1 (t), φ is the oscillation frequency of OSC1 (t), β is the amplitude coefficient of OSC1 (t). Let the two equations in formula (4) be zero to get the lowest point of the center of inertia frequency time t COI and the oscillation frequency extreme point time t OSC Expression: k is a non-negative integer; S4, let t COI = t OSC， get k * The floor function is denoted by [k * ], which gives the three extreme points in time of the oscillating term, k * To satisfy the condition t COI = t OSC a non-negative constant; S5, maximum frequency deviation calculation When β < 0 and [k * ] is odd, or β > 0 and [k * ] is even, the maximum frequency offset of the two regions is respectively: when β < 0 and [k * ] is even, or β > 0 and [k * ] is odd, the maximum frequency deviation of the two regions are respectively: where Δf 1max is the maximum frequency deviation of region one, Δf 2max is the maximum frequency deviation of region two, Δf1(t1) is the frequency variation of region one at time t1, Δf1(t2) is the frequency variation of region one at time t2, Δf1(t3) is the frequency variation of region one at time t3, Δf2(t1) is the frequency variation of region two at time t1, Δf2(t2) is the frequency variation of region two at time t2, and Δf2(t3) is the frequency variation of region two at time t3.

2. The power grid multi-region frequency offset spatiotemporal differentiation mechanism understanding and collaborative control method according to claim 1, characterized in that, The two-region frequency response model TAFR is: (2H1+ 2H2)y1"+ (D1' + D2')y1' = ω0(ΔP F1 + ΔP F2 + ΔP d1 + ΔP d2 )(12); In the formula, ΔP di is the disturbance power for each area, ω0is the rated angular velocity, ΔP Fi is the first order inertia part of the unbalanced accelerating power for each area, k is the nodal admittance coefficient of the system, y represents the frequency response expression and the results after differentiation or integration, y1' is the inertia center frequency response expression, y1" is the first order derivative of y1', y2 represents the results after integration of y2', y2' is the oscillation frequency response expression, and y2" is the first order derivative of y2'.

3. The power grid multi-region frequency deviation spatiotemporal differentiation mechanism analysis and collaborative control method according to claim 2, characterized in that, Δf COI The calculation formula of (t) is:

4. The power grid multi-region frequency deviation spatio-temporal differentiation mechanism analysis and cooperative control method according to claim 3, characterized in that, Δf OSC1 The calculation formula of (t) is: In formula (15), ω' is the combined coefficient of the nodal admittance and the rated angular velocity, D1' is the equivalent damping coefficient of region one, and D2' is the equivalent damping coefficient of region two.

5. The power grid multi-region frequency deviation spatiotemporal differentiation machine explanation and collaborative control method according to claim 4, characterized in that, Δf OSC2 The calculation formula of (t) is: In formula (17) γ1, γ2 are Δf OSC2 the amplitude coefficient of (t), is Δf OSC2 the initial phase angle of (t).

6. The power grid multi-region frequency deviation spatio-temporal differentiation mechanism analysis and cooperative control method according to claim 5, characterized in that, After adding control to the fan, the change in output power is: In the formula, ΔP w Let K1 be the change in active power of the generator, Δf be the change in system frequency, K1 be the droop coefficient, and K2 be the inertia coefficient. When the system experiences an active power disturbance and the frequency fluctuates, we have... In the formula, ΔP d For the disturbance power, ΔP G For the synchronous machine to output active power, ΔP M The unbalanced acceleration power during generator frequency regulation; substituting equation (19) into equation (20) yields In equation (22), H W D is the equivalent inertial time constant of the system when the wind turbine participates in frequency regulation mode. W The equivalent damping coefficient when the wind turbine participates in frequency regulation; for the regional interconnection system with wind power access, the rotor motion equation of the SFR model of the interconnection area is expressed as: In the formula, ΔP di For the disturbance power in each region, ΔP Mi For the unbalanced acceleration power during the frequency regulation process of generators in each region, ΔP Gi The active power output of the synchronous machines in each region, Δf i Δδ represents the frequency variation of the system in each region. i The rotor angle change of each generator is given; the two-region frequency response model TAFR is obtained according to formula (23).

7. The power grid multi-region frequency deviation spatiotemporal differentiation mechanism analysis and cooperative control method according to claim 6, characterized in that, The transfer function of TAFR is where, Δδ is a column vector containing the rotor angle variation of each generator, s is a complex variable in Laplace transform; H is a diagonal matrix of inertia constant, D is a diagonal matrix of damping coefficient, K is a matrix of nodal admittance coefficient, ΔP M is a column vector of unbalanced accelerating power; replace Δf i with the center frequency Δf, formula (24) is converted into a two-degree-of-freedom vibration equation: where, ΔP F is a column vector of unbalanced accelerating power first-order inertia part; formula (27) is solved by using the mode superposition method of solving two-degree-of-freedom vibration equation in structural dynamics, and My″+Cy′+Zy=F(t)(28) is obtained, in formula (28), Assume that the system is uniformly damped, let D′1 / H1=D′2 / H2, decouple formula (28) to obtain formulas (12) and (13); in formula (29), M is a diagonal matrix containing inertia constants of each region, C is a matrix of damping coefficients and combined coefficients of inertia constants, Z is a diagonal matrix of combined coefficients of inertia constants and nodal admittance, and F(t) is a column vector of unbalanced accelerating power and disturbance power.

Images