Method and system for determining inertia margin of high-proportion power electronic system

Abstract

The invention discloses a high-proportion power electronic system inertia margin determination method and system, and the method comprises the steps: firstly constructing an optimization model with the minimum difference value between the minimum inertia margin of a system and the equivalent inertia of the system as a target, and then solving the optimization model based on a chaos snow ablation optimization algorithm, and obtaining the optimal value of the minimum inertia margin of the system. According to the method, the stability and reliability when the power system is disturbed can be improved, the accuracy of the minimum inertia margin value of the system is improved, and the optimization speed and precision are effectively improved.

Description

Technical Field

[0001] The present invention belongs to the field of power system inertia assessment, and in particular relates to a method and system for determining inertia margin of a high-proportion power electronic system. Background Art

[0002] A reasonable and accurate inertia margin assessment system and method helps the power system maintain a certain inertia level, which plays an important role in the safety and stability of the power system frequency. With the transformation of the power grid to a new power system, a high proportion of power electronic systems is an inevitable trend in the development of future power systems.

[0003] However, due to the large-scale access of new energy and DC transmission systems to the power system, new challenges have emerged in the frequency stability of the power system. Due to the low rotational inertia characteristics of new energy generators and energy storage devices, a series of frequency stability problems caused by the decline in system inertia levels have become increasingly prominent. Inertia, as an important characteristic parameter that determines the dynamic behavior of the frequency after the power grid is disturbed, has complex and diverse manifestations in the new energy power system. Therefore, a power system inertia assessment method has been introduced to guide power system planning by assessing the inertia margin of the power system. This method has been applied in power system stable operation planning and other occasions.

[0004] Although this method can evaluate the inertia margin of the power system to guide power system planning, it still has the following defects:

[0005] 1. The existing inertia margin assessment method does not take into account the characteristics of high proportion of power electronics in future power systems. In future power systems, there will be a large number of power electronic equipment in new energy stations and energy storage devices. New energy stations and energy storage devices considering virtual inertia control can also provide certain inertia support for the power system.

[0006] 2. The existing technology considers incomplete constraints, generally considering the maximum frequency change rate constraint or the maximum frequency deviation constraint, that is, an inertia margin evaluation system has not been formed.

[0007] 3. The existing technology uses solvers or intelligent optimization algorithms to solve the inertia evaluation problem, which has the problems of slow solution speed and inaccurate solution. Summary of the invention

[0008] The object of the present invention is to provide a method and system for determining inertia margin of a high-ratio power electronic system in view of the above-mentioned problems existing in the prior art.

[0009] To achieve the above objectives, the technical solution of the present invention is as follows:

[0010] In a first aspect, the present invention provides a method for determining an inertia margin of a high-ratio power electronic system, comprising:

[0011] S1. Construct an optimization model with the goal of minimizing the difference between the system minimum inertia margin and the system equivalent inertia;

[0012] S2. Solve the optimization model to obtain the optimal value of the system's minimum inertia margin.

[0013] In S1, the objective function of the optimization model includes:

[0014] min|H minc -H sys |

[0015] H sys =H G (K NE (K NE-ES )+1)

[0016] In the above formula, H minc is the minimum inertia margin of the system, H sys is the system equivalent inertia, H G is the inertia time constant of the running traditional synchronous unit, K NE is the penetration rate of new energy, K NE-ES It is the virtual inertia penetration rate of new energy unit-energy storage combination;

[0017] The constraints of the optimization model include the maximum frequency change rate constraint, the maximum frequency deviation constraint, the equivalent inertia constraint, and the primary frequency modulation recovery frequency difference constraint.

[0018] The K NE , K NE-ES Calculated using the following formula:

[0019]

[0020] In the above formula, P NE,j , H NE,j are the rated capacity and inertia time constant of the running new energy unit j, P G,i , H G,i are the power and inertia time constant of the running traditional synchronous unit i, P ES,k , H ES,k are the power and inertia time constant of the system energy storage device k, J ES,k is the virtual inertia of the system energy storage device k, w s is the angular velocity of the synchronous unit, J G is the moment of inertia of the synchronous unit, k ES is the ratio of the battery state of charge change rate to the generator speed change rate, u B_ES is the rated voltage of the system energy storage device, P N_ES is the rated power of the battery, γ soc_initial is the initial state of charge, EG is the kinetic energy of the synchronous unit.

[0021] The maximum frequency change rate constraint includes:

[0022] |RoCoF max| ≤|RoCoF lim |

[0023]

[0024] In the above formula, RoCoF max 、RoCoF lim are the system maximum frequency change rate and its limit, ΔP d is the system power disturbance, f0 is the system initial frequency;

[0025] The maximum frequency deviation constraint includes:

[0026] |Δf t |≤Δf lim

[0027] Δf t =Δf t-1 +R t Δt

[0028]

[0029] In the above formula, Δf t is the system frequency deviation at time t, Δf lim is the system frequency deviation limit, R t is the rate of change of system frequency at time t, Δt is the time interval, ΔP m,t-1 is the primary frequency modulation power adjustment at time t-1, D sys is the system damping coefficient;

[0030] The equivalent inertia constraints include:

[0031] H min ≤H sys ≤H max

[0032]

[0033] In the above formula, H max , H min are the upper and lower limits of the system equivalent inertia, α is the frequency space distribution index, S sys is the system capacity, H k is the equivalent inertia at node k, R max,k is the maximum frequency change rate at node k;

[0034] The primary frequency modulation recovery frequency difference constraint includes:

[0035]

[0036] In the above formula, Δf pfd K is the frequency difference restored by a frequency modulation. sys is the primary frequency modulation coefficient of the system, f db K is the primary frequency modulation dead zone of the system, G It is the equivalent primary frequency regulation coefficient of the synchronous unit.

[0037] The S2 uses a chaotic snow melting optimization algorithm to solve the optimization model, including:

[0038] S21, randomly generate an initial solution in the feasible region, that is, the initial value of the minimum inertia margin of the system;

[0039] S22, using Logistic chaotic mapping to process the initial solution;

[0040] S23, by simulating the evaporation of liquid water and the Brownian motion of water vapor to traverse the feasible domain, generate a population P a , and generate population P by simulating the snow melting process b ;

[0041] S24. Calculate population P a and P b The fitness of each individual in the optimization model, that is, the objective function value of the optimization model;

[0042] S25, taking the population individual with the best fitness as the updated population;

[0043] S26, determine whether the iteration termination condition is met, if so, output the current population as the optimal result; if not, return to S23 for the next iteration.

[0044] In S22, the initial solution after the Logistic chaotic mapping process is determined based on the following formula:

[0045]

[0046] s1=μ×s×(1-s)

[0047] n=n 0,min +s1(n 0,max -n 0,min )

[0048] In the above formula, s and s1 are the intermediate quantities generated by the Logistic chaotic map, n0 is the randomly generated initial solution, and n 0,max 、n 0,minare the maximum and minimum values ​​of n0 respectively, μ is the chaotic mapping parameter, and n is the initial solution after Logistic mapping processing.

[0049] In a second aspect, the present invention provides a system for determining inertia margin of a high-proportion power electronic system, comprising a model building module and a model solving module;

[0050] The model building module is used to build an optimization model with the goal of minimizing the difference between the system minimum inertia margin and the system equivalent inertia;

[0051] The model solving module is used to solve the optimization model to obtain the optimal value of the minimum inertia margin of the system.

[0052] The optimization model objective function includes:

[0053] min|H minc -H sys |

[0054] H sys =H G (K NE (K NE-ES )+1)

[0055]

[0056]

[0057] In the above formula, H minc is the minimum inertia margin of the system, H sys is the system equivalent inertia, H G is the inertia time constant of the running traditional synchronous unit, K NE is the penetration rate of new energy, K NE-ES is the virtual inertia penetration rate of new energy unit-energy storage, P NE,j , H NE,j are the rated capacity and inertia time constant of the running new energy unit j, P G,i , H G,i are the power and inertia time constant of the running traditional synchronous unit i, P ES,k , H ES,k are the power and inertia time constant of the system energy storage device k, J ES,k is the virtual inertia of the system energy storage device k, w s is the angular velocity of the synchronous unit, J G is the moment of inertia of the synchronous unit, k ES is the ratio of the battery state of charge change rate to the generator speed change rate, u B_ES is the rated voltage of the system energy storage device, P N_ES is the rated power of the battery, γsoc_initial is the initial charge state, E G is the kinetic energy of the synchronous unit;

[0058] The constraints of the optimization model include:

[0059] Maximum frequency change rate constraint

[0060] |RoCoF max |≤|RoCoF lim |

[0061]

[0062] In the above formula, RoCoF max 、RoCoF lim are the system maximum frequency change rate and its limit, ΔP d is the system power disturbance, f0 is the system initial frequency;

[0063] Maximum frequency deviation constraint

[0064] |Δf t |≤Δf lim

[0065] Δf t =Δf t-1 +R t Δt

[0066]

[0067] In the above formula, Δf t is the system frequency deviation at time t, Δf lim is the system frequency deviation limit, R t is the rate of change of system frequency at time t, Δt is the time interval, ΔP m,t-1 is the primary frequency modulation power adjustment at time t-1, D sys is the system damping coefficient;

[0068] Equivalent inertia constraint:

[0069] H min ≤H sys ≤H max

[0070]

[0071] In the above formula, H max , H min are the upper and lower limits of the system equivalent inertia, α is the frequency space distribution index, S sys is the system capacity, H k is the equivalent inertia at node k, R max,kis the maximum frequency change rate at node k;

[0072] Frequency difference constraint after one frequency modulation:

[0073]

[0074] In the above formula, Δf pfd K is the frequency difference restored by a frequency modulation. sys is the primary frequency modulation coefficient of the system, f db K is the primary frequency modulation dead zone of the system, G It is the equivalent primary frequency regulation coefficient of the synchronous unit.

[0075] The model solving module is used to solve the optimization model using the chaotic snow melting optimization algorithm. The specific process includes:

[0076] A. Randomly generate an initial solution in the feasible region, that is, the initial value of the system's minimum inertia margin;

[0077] B. Use Logistic chaotic mapping to process the initial solution;

[0078] C. Generate population P by traversing the feasible domain by simulating the evaporation of liquid water and the Brownian motion of water vapor a , and generate population P by simulating the snow melting process b ;

[0079] D. Calculate the population P a and P b The fitness of each individual in the optimization model, that is, the objective function value of the optimization model;

[0080] E. The individuals in the population with the best fitness are used as the updated population;

[0081] F. Determine whether the iteration termination condition is met. If so, output the current population as the optimal result; if not, return to C for the next iteration.

[0082] The initial solution after the Logistic chaotic mapping process is determined based on the following formula:

[0083]

[0084] s1=μ×s×(1-s)

[0085] n=n 0,min +s1(n 0,max -n 0,min )

[0086] In the above formula, s and s1 are the intermediate quantities generated by the Logistic chaotic map, n0 is the randomly generated initial solution, and n 0,max 、n 0,minare the maximum and minimum values ​​of n0 respectively, μ is the chaotic mapping parameter, and n is the initial solution after Logistic mapping processing.

[0087] Compared with the prior art, the present invention has the following beneficial effects:

[0088] 1. The optimization model constructed by the method for determining the inertia margin of a high-proportion power electronic system proposed in the present invention takes the minimum difference between the system minimum inertia margin and the system equivalent inertia as the optimization goal, which can further guide the inertia level design during power system design and improve the stability and reliability of the power system when disturbed.

[0089] 2. The optimization model constraints constructed by the method for determining the inertia margin of a high-proportion power electronic system proposed in the present invention include the maximum frequency change rate constraint, the maximum frequency deviation constraint, the equivalent inertia constraint and the primary frequency modulation recovery frequency difference constraint. The model simultaneously considers the inertial response of the high-proportion power electronic system after disturbance and its primary frequency modulation response. By establishing a complete inertia margin evaluation system for a high-proportion power electronic system, the accuracy of the minimum inertia margin value of the system is effectively improved.

[0090] 3. The method for determining the inertia margin of a high-proportion power electronic system proposed in the present invention takes into account the influence of traditional synchronous units and new energy units based on virtual inertia control, the influence of electrical distance between nodes on node inertia changes, the influence of load fluctuations on system inertia, and the influence of energy storage devices based on virtual inertia control on system equivalent inertia when calculating the system equivalent inertia. That is, the influence of multiple factors such as "source-grid-load-storage" is considered at the same time, thereby ensuring the accuracy and reliability of the system equivalent inertia, and further improving the accuracy of the system minimum inertia margin value.

[0091] 4. A method for determining inertia margin of high-proportion power electronic system proposed in the present invention adopts a chaotic snow melting optimization algorithm to solve the optimization model. On the one hand, the algorithm introduces Logistic chaotic mapping to process the initial solution to avoid the initial solution falling into the local optimum; on the other hand, the algorithm adopts a dual population mechanism and uses the population P a The function of traversing the feasible domain is realized by simulating the evaporation of liquid water and the Brownian motion of water vapor. b By simulating the snow melting process, the function of exploring the global optimal solution near the local optimal solution is realized, which not only ensures that the population traverses the feasible domain, but also ensures that the population can explore the global optimal solution near the local optimal solution, thereby improving the optimization speed and accuracy. BRIEF DESCRIPTION OF THE DRAWINGS

[0092] Figure 1 This is a topological diagram of the IEEE39 node system used in Example 1.

[0093] Figure 2 This is an overall flow chart of the method described in Example 1.

[0094] Figure 3 This is a flow chart of the chaotic snow melting optimization algorithm in Example 1.

[0095] Figure 4 This is a convergence curve diagram of the chaotic snow melting optimization algorithm in Example 1.

[0096] Figure 5 This is a graph showing the relationship between the system's minimum inertia margin and the new energy penetration rate.

[0097] Figure 6 It is a system topology diagram including synchronous units, wind-storage combined system and load.

[0098] Figure 7 This is the result diagram of minimum inertia margin under different fluctuation scenarios.

[0099] Figure 8 This is a graph showing the relationship between the minimum inertia margin and the virtual inertia permeability.

[0100] Fig. 9 This is a structural diagram of the system described in Example 2. DETAILED DESCRIPTION

[0101] The present invention is further described in detail below in conjunction with specific implementations and drawings.

[0102] Embodiment 1:

[0103] This embodiment is for the modified IEEE39 node system with new energy access built in MATLAB (its topology diagram is shown in Figure 1 ), the system includes 39 nodes, 45 transmission corridors, and 10 power sources.

[0104] A method for determining inertia margin of high-proportion power electronic system, such as Figure 2 As shown, follow the steps below to implement it:

[0105] 1. Build an optimization model

[0106] The model aims to minimize the difference between the system's minimum inertia margin and the system's equivalent inertia, and its objective function is:

[0107] min|H minc -H sys |

[0108] H sys =H G (K NE (K NE-ES )+1)

[0109]

[0110] In the above formula, H minc is the minimum inertia margin of the system, H sys is the system equivalent inertia, H G is the inertia time constant of the running traditional synchronous unit, K NE is the penetration rate of new energy, K NE-ES is the virtual inertia penetration rate of new energy unit-energy storage, P NE,j , H NE,j are the rated capacity and inertia time constant of the running new energy unit j, P G,i , H G,i are the power and inertia time constant of the running traditional synchronous unit i, P ES,k , H ES,k are the power and inertia time constant of the system energy storage device k, J ES,k is the virtual inertia of the system energy storage device k, w s is the angular velocity of the synchronous unit, J G is the moment of inertia of the synchronous unit, k ES is the ratio of the battery state of charge change rate to the generator speed change rate, u B_ES is the rated voltage of the system energy storage device, P N_ES is the rated power of the battery, γ soc_initial is the initial state of charge, E G is the kinetic energy of the synchronous unit;

[0111] The constraints include:

[0112] Maximum frequency change rate constraint

[0113] |RoCoF max |≤|RoCoF lim |

[0114]

[0115] In the above formula, RoCoF max 、RoCoF lim They are the maximum frequency change rate of the system and its limit. It is usually believed that the frequency change rate of the system is the largest when it is disturbed. d is the system power disturbance, f0 is the system initial frequency;

[0116] Maximum frequency deviation constraint

[0117] |Δf t |≤Δf lim

[0118] Δf t =Δf t-1 +Rt Δt

[0119]

[0120] In the above formula, Δf t is the system frequency deviation at time t, Δf lim is the system frequency deviation limit, R t is the rate of change of system frequency at time t, Δt is the time interval, ΔP m,t-1 is the primary frequency modulation power adjustment at time t-1, D sys is the system damping coefficient;

[0121] Equivalent inertia constraint:

[0122] H min ≤H sys ≤H max

[0123]

[0124]

[0125] In the above formula, H min is the lower limit of the system equivalent inertia, taking the system equivalent inertia without virtual inertia support, H max is the upper limit of the system equivalent inertia, and the system equivalent inertia when the virtual inertia is the maximum is taken. In order to consider the influence of electrical distance on node inertia, the frequency space distribution index α, which represents the electrical distance between nodes, is introduced. sys is the system capacity, H k is the equivalent inertia at node k, R max,k is the maximum frequency change rate at node k;

[0126] Frequency difference constraint after one frequency modulation:

[0127]

[0128] In the above formula, Δf pfd K is the frequency difference restored by a frequency modulation. sys is the primary frequency modulation coefficient of the system, f db K is the primary frequency modulation dead zone of the system, G It is the equivalent primary frequency regulation coefficient of the synchronous unit.

[0129] In this embodiment, the parameter values ​​of the above model are shown in Table 1:

[0130] Table 1 Parameters of the minimum inertia margin evaluation optimization model

[0131] content Numeric Load fluctuation 36.6MW System equivalent inertia 7s New energy penetration rate 0.3 Virtual inertia permeability 0 System damping coefficient 0.1

[0132] 2. The chaotic snow melting optimization algorithm is used to solve the optimization model and obtain the optimal value of the system's minimum inertia margin.

[0133] The snow melting algorithm randomly generates the initial solution in the feasible domain, which may cause the initial solution to fall into the local optimum directly, causing premature maturity. Therefore, this step introduces the Logistic chaotic map and uses the highly random characteristics of the chaotic map to improve the snow melting algorithm. The specific process is as follows: Figure 3 As shown, including:

[0134] 2.1. Initialization parameters, including randomly generating the initial solution in the feasible domain, i.e., the initial value of the minimum inertia margin of the system, setting the population size to 20 and the maximum number of iterations to 20.

[0135] 2.2. Use Logistic chaotic mapping to process the initial solution:

[0136]

[0137] s1=μ×s×(1-s)

[0138] n=n 0,min +s1(n 0,max -n 0,min )

[0139] In the above formula, s and s1 are the intermediate quantities generated by the Logistic chaotic map, n0 is the randomly generated initial solution, and n 0,max 、n 0,min are the maximum and minimum values ​​of n0 respectively, μ is the chaotic mapping parameter, and n is the initial solution after Logistic mapping processing.

[0140] 2.3. Generate population P by traversing the feasible domain by simulating the evaporation of liquid water and the Brownian motion of water vapor a , and generate population P by simulating the snow melting process b ;

[0141] Generate population P a In the process, when the molecule exists in the form of water vapor, it traverses the solution space according to Brownian motion. The position update expression of Brownian motion is as follows:

[0142]

[0143] In the above formula, Z i (t+1) represents the position of the i-th molecule at the t+1th time; Elite(t) is a random variable, which is obtained from G(t), Z second (t), Z third (t) and Z c (t) Randomly select from four quantities; BM i(t) is a Gaussian distribution based on Brownian motion, θ1 is a random number taken from the interval [0,1], and G(t) is the optimal position of the molecule at time t; is the average position of all molecules, Z second (t) is the suboptimal position of the molecule at time t; Z third (t) is the third position of the molecule at time t; Z c (t) is the average position of the molecule at time t; N is the population of all molecules, and N1 is the population of molecules with the top 50% fitness;

[0144] Generate population P b In the process, when the molecules exist in the form of snow, the process of snow conversion into liquid water is simulated according to the snow melting process, and its expression is as follows:

[0145] M=DDF×(T-T1)

[0146]

[0147] In the above formula, DDF is the snow ablation coefficient in the snow ablation model, which ranges from [0.35, 0.6], T is the daily average temperature, and T1 is the basic temperature;

[0148] The expression for simulating the snow melting process to update the population position is as follows:

[0149]

[0150] In the above formula, M is the snow melting rate, and θ2 is a random number ranging from [-1, 1].

[0151] 2.4. Calculate the population P a and P b The fitness of each individual in is the objective function value of the optimization model.

[0152] 2.5. The individuals in the population with the best fitness are taken as the updated population.

[0153] 2.6. Determine whether the iteration termination condition is met, that is, whether the maximum number of iterations is reached. If so, the population at this time is output as the optimal result; if not, return to 2.3 for the next iteration.

[0154] The algorithm convergence curve of this embodiment is as follows Figure 4 shown.

[0155] Effect test analysis:

[0156] 1. For comparative analysis, the following three cases are set (load disturbance occurs at 1.5s in all three cases):

[0157] (1) Case 1, the load at node 15 is reduced by 36.6MW. The system inertia is set according to the optimal value of the system minimum inertia margin obtained by the method described in this embodiment, and the inertia change of node 16 is observed;

[0158] (2) Case 2, the load at node 15 is reduced by 36.6MW. The system inertia is set according to the optimal value of the system minimum inertia margin obtained by the method described in this embodiment, and the inertia change of node 23 is observed;

[0159] (3) Case 3, the load at node 15 is reduced by 36.6 MW. The system inertia is set according to the optimal value of the system minimum inertia margin obtained by the method described in this embodiment, and the inertia change of node 8 is observed.

[0160] The minimum inertia margin values ​​for the above three cases are shown in Table 3:

[0161] Table 3 System minimum inertia margin values ​​under different cases

[0162] Case node Minimum inertia margin Case 1 16 4.8 Case 2 23 4.5 Case 3 8 3.5

[0163] As shown in Table 3, the time for Case 1 system to recover balance is 4.8s, the time for Case 2 system to recover balance is 4.5s, and the time for Case 3 system to recover balance is 3.5s. Within the error range, they are all consistent with the simulation results of the IEEE 39-node system example. The electrical distances of the three cases from the node 15 where the load disturbance occurs are from near to far. By comparing Case 1, Case 2 and Case 3, it can be seen that the farther the electrical distance from the load disturbance node, the smaller the minimum inertia margin value of the system, and the less affected by the load disturbance.

[0164] 2. In order to explore the impact of new energy penetration rate on the system minimum inertia margin, for Case 1, the system minimum inertia margin under different new energy penetration rates was investigated. The results are as follows: Figure 5 shown.

[0165] Depend on Figure 5 It can be seen that the higher the penetration rate of new energy, the higher the minimum inertia margin of the system.

[0166] 3. In order to test the effectiveness and applicability of the method proposed in the present invention in a power system with a high proportion of power electronics, a system including synchronous units, wind-storage combined system and loads is constructed in MATLAB (its topology is shown in the figure below). Figure 6 As shown in Table 4, multiple wind-storage combined system fluctuations and different load fluctuation scenarios are set, and the minimum inertia margin under each scenario is determined by the method described in the present invention. The influence of the change of virtual inertia permeability on the minimum inertia margin is investigated. The results are shown in Table 4. Figure 7-8 As shown:

[0167] Table 4 System equivalent inertia and minimum inertia margin results under different scenarios

[0168]

[0169] Comparing Scenario 1 and Scenario 2, when the output of the wind-storage combined system fluctuates, the greater the output fluctuation, the smaller the system minimum inertia margin; comparing Scenario 3 and Scenario 4, the greater the load fluctuation, the greater the system minimum inertia margin; comparing Scenario 5 and Scenario 4, when the system damping coefficient is larger, the system minimum inertia margin is larger; comparing Scenario 6 and Scenario 4, when the system equivalent inertia increases, the required inertia of the system is smaller.

[0170] Figure 8 The relationship between the minimum inertia margin and the virtual inertia penetration rate is shown in Figure 2. It can be seen that the greater the virtual inertia penetration rate, the smaller the minimum inertia margin. This shows that the wind-storage joint system based on virtual inertia control has a certain effect on improving the inertia level of high-proportion power electronic systems.

[0171] Embodiment 2:

[0172] A system for determining inertia margin of a high-proportion power electronic system, such as Fig. 9 As shown, it includes a model building module and a model solving module.

[0173] The model building module is used to build the following optimization model:

[0174] Objective function:

[0175] min|H minc -H sys |

[0176] H sys =H G (K NE (K NE-ES )+1)

[0177]

[0178] In the above formula, H minc is the minimum inertia margin of the system, H sys is the system equivalent inertia, H G is the inertia time constant of the running traditional synchronous unit, K NE is the penetration rate of new energy, K NE-ES is the virtual inertia penetration rate of new energy unit-energy storage, P NE,j , H NE,j are the rated capacity and inertia time constant of the running new energy unit j, P G,i , H G,i are the power and inertia time constant of the running traditional synchronous unit i, PES,k , H ES,k are the power and inertia time constant of the system energy storage device k, J ES,k is the virtual inertia of the system energy storage device k, w s is the angular velocity of the synchronous unit, J G is the moment of inertia of the synchronous unit, k ES is the ratio of the battery state of charge change rate to the generator speed change rate, u B_ES is the rated voltage of the system energy storage device, P N_ES is the rated power of the battery, γ soc_initial is the initial charge state, E G is the kinetic energy of the synchronous unit;

[0179] Constraints:

[0180] Maximum frequency change rate constraint

[0181] |RoCoF max |≤|RoCoF lim |

[0182]

[0183] In the above formula, RoCoF max 、RoCoF lim are the system maximum frequency change rate and its limit, ΔP d is the system power disturbance, f0 is the system initial frequency;

[0184] Maximum frequency deviation constraint

[0185] |Δf t |≤Δf lim

[0186] Δf t =Δf t-1 +R t Δt

[0187]

[0188] In the above formula, Δf t is the system frequency deviation at time t, Δf lim is the system frequency deviation limit, R t is the rate of change of system frequency at time t, Δt is the time interval, ΔP m,t-1 is the primary frequency modulation power adjustment at time t-1, D sys is the system damping coefficient;

[0189] Equivalent inertia constraint:

[0190] H min ≤Hsys ≤H max

[0191]

[0192] In the above formula, H max , H min are the upper and lower limits of the system equivalent inertia, α is the frequency space distribution index, S sys is the system capacity, H k is the equivalent inertia at node k, R max,k is the maximum frequency change rate at node k;

[0193] Frequency difference constraint after one frequency modulation:

[0194]

[0195] In the above formula, Δf pfd K is the frequency difference restored by a frequency modulation. sys is the primary frequency modulation coefficient of the system, f db K is the primary frequency modulation dead zone of the system, G It is the equivalent primary frequency regulation coefficient of the synchronous unit.

[0196] The model solving module is used to solve the optimization model using the chaotic snow melting optimization algorithm to obtain the optimal value of the minimum inertia margin of the system. The specific process includes:

[0197] A. Randomly generate an initial solution in the feasible region, that is, the initial value of the system's minimum inertia margin;

[0198] B. Use Logistic chaotic mapping to process the initial solution:

[0199]

[0200] s1=μ×s×(1-s)

[0201] n=n 0,min +s1(n 0,max -n 0,min )

[0202] In the above formula, s and s1 are the intermediate quantities generated by the Logistic chaotic map, n0 is the randomly generated initial solution, and n 0,max 、n 0,min are the maximum and minimum values ​​of n0, μ is the chaotic mapping parameter, and n is the initial solution after Logistic mapping processing;

[0203] C. Generate population P by traversing the feasible domain by simulating the evaporation of liquid water and the Brownian motion of water vapor a , and generate population P by simulating the snow melting process b;

[0204] D. Calculate the population P a and P b The fitness of each individual in the optimization model, that is, the objective function value of the optimization model;

[0205] E. The individuals in the population with the best fitness are used as the updated population;

[0206] F. Determine whether the iteration termination condition is met. If so, output the current population as the optimal result; if not, return to C for the next iteration.

Claims

1. A method for determining inertia margin of a high-proportion power electronic system, characterized in that: The method comprises: S1. Construct an optimization model with the goal of minimizing the difference between the system minimum inertia margin and the system equivalent inertia; S2. Solve the optimization model to obtain the optimal value of the system's minimum inertia margin.

2. The method for determining inertia margin of a high-proportion power electronic system according to claim 1, characterized in that: In S1, the objective function of the optimization model includes: min|H minc -H sys | H sys =H G (K NE (K NE-ES )+1) In the above formula, H minc is the minimum inertia margin of the system, H sys is the system equivalent inertia, H G is the inertia time constant of the running traditional synchronous unit, K NE is the penetration rate of new energy, K NE-ES It is the virtual inertia penetration rate of new energy unit-energy storage combination; The constraints of the optimization model include the maximum frequency change rate constraint, the maximum frequency deviation constraint, the equivalent inertia constraint, and the primary frequency modulation recovery frequency difference constraint.

3. The method for determining inertia margin of a high-proportion power electronic system according to claim 2, characterized in that: The K NE , K NE-ES Calculated using the following formula: In the above formula, P NE,j , H NE,j are the rated capacity and inertia time constant of the running new energy unit j, P G,i , H G,i are the power and inertia time constant of the running traditional synchronous unit i, P ES,k , H ES,k are the power and inertia time constant of the system energy storage device k, J ES,k is the virtual inertia of the system energy storage device k, w s is the angular velocity of the synchronous unit, J G is the moment of inertia of the synchronous unit, k ES is the ratio of the battery state of charge change rate to the generator speed change rate, u B_ES is the rated voltage of the system energy storage device, P N_ES is the rated power of the battery, γ soc_initial is the initial state of charge, E G is the kinetic energy of the synchronous unit.

4. The method for determining inertia margin of a high-proportion power electronic system according to claim 2, characterized in that: The maximum frequency change rate constraint includes: |RoCoF max |≤|RoCoF lim | In the above formula, RoCoF max 、RoCoF lim are the system maximum frequency change rate and its limit, ΔP d is the system power disturbance, f0 is the system initial frequency; The maximum frequency deviation constraint includes: |Δf t |≤Δf lim Δf t =Δf t-1 +R t Δt In the above formula, Δf t is the system frequency deviation at time t, Δf lim is the system frequency deviation limit, R t is the rate of change of system frequency at time t, Δt is the time interval, ΔP m,t-1 is the primary frequency modulation power adjustment at time t-1, D sys is the system damping coefficient; The equivalent inertia constraints include: H min ≤H sys ≤H max In the above formula, H max , H min are the upper and lower limits of the system equivalent inertia, α is the frequency space distribution index, S sys is the system capacity, H k is the equivalent inertia at node k, R max,k is the maximum frequency change rate at node k; The primary frequency modulation recovery frequency difference constraint includes: In the above formula, Δf pfd K is the frequency difference restored by a frequency modulation. sys is the primary frequency modulation coefficient of the system, f db K is the primary frequency modulation dead zone of the system, G It is the equivalent primary frequency regulation coefficient of the synchronous unit.

5. The method for determining inertia margin of a high-proportion power electronic system according to claim 1, characterized in that: The S2 uses a chaotic snow melting optimization algorithm to solve the optimization model, including: S21, randomly generate an initial solution in the feasible region, that is, the initial value of the minimum inertia margin of the system; S22, using Logistic chaotic mapping to process the initial solution; S23, by simulating the evaporation of liquid water and the Brownian motion of water vapor to traverse the feasible domain, generate a population P a , and generate population P by simulating the snow melting process b ; S24. Calculate population P a and P b The fitness of each individual in the optimization model, that is, the objective function value of the optimization model; S25, taking the population individual with the best fitness as the updated population; S26, determine whether the iteration termination condition is met, if so, output the current population as the optimal result; if not, return to S23 for the next iteration.

6. A method for determining inertia margin of a high-proportion power electronic system according to claim 5, characterized in that: In S22, the initial solution after the Logistic chaotic mapping process is determined based on the following formula: s1=μ×s×(1-s) n=n 0,min +s1(n 0,max -n 0,min ) In the above formula, s and s1 are the intermediate quantities generated by the Logistic chaotic map, n0 is the randomly generated initial solution, and n 0,max 、n 0,min are the maximum and minimum values ​​of n0 respectively, μ is the chaotic mapping parameter, and n is the initial solution after Logistic mapping processing.

7. A high-proportion power electronic system inertia margin determination system, characterized in that: The system includes a model building module and a model solving module; The model building module is used to build an optimization model with the goal of minimizing the difference between the system minimum inertia margin and the system equivalent inertia; The model solving module is used to solve the optimization model to obtain the optimal value of the minimum inertia margin of the system.

8. The high-ratio power electronic system inertia margin determination system according to claim 7, characterized in that: The optimization model objective function includes: min|H minc -H sys | H sys =H G (K NE (K NE-ES )+1) In the above formula, H minc is the minimum inertia margin of the system, H sys is the system equivalent inertia, H G is the inertia time constant of the running traditional synchronous unit, K NE is the penetration rate of new energy, K NE-ES is the virtual inertia penetration rate of new energy unit-energy storage, P NE,j , H NE,j are the rated capacity and inertia time constant of the running new energy unit j, P G,i , H G,i are the power and inertia time constant of the running traditional synchronous unit i, P ES,k , H ES,k are the power and inertia time constant of the system energy storage device k, J ES,k is the virtual inertia of the system energy storage device k, w s is the angular velocity of the synchronous unit, J G is the moment of inertia of the synchronous unit, k ES is the ratio of the battery state of charge change rate to the generator speed change rate, u B_ES is the rated voltage of the system energy storage device, P N_ES is the rated power of the battery, γ soc_initial is the initial state of charge, E G is the kinetic energy of the synchronous unit; The constraints of the optimization model include: Maximum frequency change rate constraint |RoCoF max |≤|RoCoF lim | In the above formula, RoCoF max 、RoCoF lim are the system maximum frequency change rate and its limit, ΔP d is the system power disturbance, f0 is the system initial frequency; Maximum frequency deviation constraint |Δf t |≤Δf lim Δf t =Δf t-1 +R t Δt In the above formula, Δf t is the system frequency deviation at time t, Δf lim is the system frequency deviation limit, R t is the rate of change of system frequency at time t, Δt is the time interval, ΔP m,t-1 is the primary frequency modulation power adjustment at time t-1, D sys is the system damping coefficient; Equivalent inertia constraint: H min ≤H sys ≤H max In the above formula, H max , H min are the upper and lower limits of the system equivalent inertia, α is the frequency space distribution index, S sys is the system capacity, H k is the equivalent inertia at node k, R max,k is the maximum frequency change rate at node k; Frequency difference constraint after one frequency modulation: In the above formula, Δf pfd K is the frequency difference restored by a frequency modulation. sys is the primary frequency modulation coefficient of the system, f db K is the primary frequency modulation dead zone of the system, G It is the equivalent primary frequency regulation coefficient of the synchronous unit.

9. The high-ratio power electronic system inertia margin determination system according to claim 7, characterized in that: The model solving module is used to solve the optimization model using the chaotic snow melting optimization algorithm. The specific process includes: A. Randomly generate an initial solution in the feasible region, that is, the initial value of the system's minimum inertia margin; B. Use Logistic chaotic mapping to process the initial solution; C. Generate population P by traversing the feasible domain by simulating the evaporation of liquid water and the Brownian motion of water vapor a , and generate population P by simulating the snow melting process b ; D. Calculate the population P a and P b The fitness of each individual in the optimization model, that is, the objective function value of the optimization model; E. The individuals in the population with the best fitness are used as the updated population; F. Determine whether the iteration termination condition is met. If so, output the current population as the optimal result; if not, return to C for the next iteration.

10. The high-ratio power electronic system inertia margin determination system according to claim 9, characterized in that: The initial solution after the Logistic chaotic mapping process is determined based on the following formula: s1=μ×s×(1-s) n=n 0,min +s1(n 0,max -n 0,min ) In the above formula, s and s1 are the intermediate quantities generated by the Logistic chaotic map, n0 is the randomly generated initial solution, and n 0,max 、n 0,min are the maximum and minimum values ​​of n0 respectively, μ is the chaotic mapping parameter, and n is the initial solution after Logistic mapping processing.

Images