A microgrid inertia constant estimation method based on African vulture algorithm

Abstract

The present invention discloses a microgrid inertia constant estimation method based on the African vulture algorithm. After Laplace transform, discretization and Z transform of the generator rotor motion equation, a transfer function expression containing the inertia constant H and damping coefficient D to be determined is obtained. The power and frequency change signals at the outlet of each unit during system operation are obtained through the PMU measurement device. Based on the measured data, the system inertia constant containing virtual inertia is finally obtained through algorithm optimization. This is used to solve the problem that existing inertia constant evaluation methods cannot be used for systems containing virtual inertia. First, a unified inertia constant analytical model of the equivalent virtual inertia of the new energy VSG and the rotational inertia of the synchronous generator is established. The problem of solving the inertia constant is converted into an algorithm optimization problem. The power-frequency time series data obtained by the PMU is used to solve the inertia constant in the unified analytical model using the African vulture algorithm. The present invention realizes the rapid identification of the system inertia constant.

Description

Technical Field

[0001] The present invention relates to the technical field of power systems, and in particular to a microgrid inertia constant estimation method based on an African vulture algorithm. Background Art

[0002] The annual increase in installed wind turbine capacity has led to a decrease in the inertia level of the power system, severely weakening the system's inertia support and frequency regulation capabilities under active power disturbances. To address this issue, applying virtual inertia control to grid-connected wind turbines can enable them to provide virtual inertia to the system, helping to improve the inertia level of new power systems. In this context, research on estimation methods that can accurately calculate the system inertia constant can quantitatively understand the contribution of grid-connected wind turbines to the system inertia level, provide a reference for subsequent comprehensive evaluation of the equivalent inertia level of new power systems, and provide reference data for dispatchers.

[0003] Most existing methods for evaluating inertia constants rely on constructing a dynamic error model and then using system identification to determine the inertia constant. System identification involves identifying a model equivalent to the measured system from a given set of model classes, based on a set of input and output data. Although traditional system identification methods have matured, they still suffer from the drawback of not being able to obtain ideal identification results for nonlinear systems. While some modern system identification methods have addressed these issues to some extent, they are generally tailored to specific model structures. Like traditional system identification methods, they are also unable to effectively address identification problems in complex system models where both the system model and parameters are unknown. Summary of the Invention

[0004] To solve the above problems, the present invention proposes to use the African vulture algorithm to evaluate the inertia constant of the new power system, and compares the results with the traditional system identification and evaluation method to solve the problem that the traditional system identification method cannot obtain ideal identification results when the model structure and parameters are unknown.

[0005] The present invention adopts the following technical solutions:

[0006] A microgrid inertia constant evaluation method based on the African vulture algorithm includes the following steps:

[0007] Step 1: Compare the new energy power source to a synchronous generator, equate it to a virtual synchronous generator, and propose the concept of virtual inertia.

[0008] Step 2: Construct a unified inertia constant analytical model including the VSG virtual inertia and the synchronous generator rotational inertia.

[0009] Step 3: The rotor motion equation of the unified inertia constant analytical model constructed in step 2 is subjected to Laplace transform, inverse Laplace transform, discretization, and Z transform in sequence, and finally a transfer function expression containing the unknown parameters H and D is obtained, where H is the inertia constant and D is the damping coefficient.

[0010] Step 4: The system operating power change data ΔP measured by the PMU o Substituting this into the calculated transfer function expression, we can obtain the corresponding frequency change Δf.

[0011] Step 5: Construct the fitness function y = ∑(Δf-Δf o ).

[0012] Step 6: Use the African vulture algorithm to optimize the parameters H and D to be determined. When the fitness function value y is within the required range, H and D calculated at that moment are taken as the parameters to be determined.

[0013] Furthermore, the unified inertia constant analytical model containing the VSG virtual inertia and the synchronous generator rotational inertia in step 2 is derived as follows:

[0014] According to the conservation of the overall inertia of the power system, the equivalent inertia constant of the power system can be calculated by the following formula:

[0015]

[0016] Among them: H gen,i represents the inertia constant of the i-th synchronous generator, S gen,i represents the rated capacity of the i-th synchronous generator, H vir,j represents the inertia constant of the jth VSG, S vir,j represents the rated capacity of the jth VSG, S sys Indicates the rated capacity of the entire system.

[0017] Furthermore, in step 3, to obtain a transfer function expression containing the parameters H and D to be determined, the following steps are included:

[0018] (1) The transient electromechanical characteristics of the system are expressed using the equivalent rotor motion equation:

[0019]

[0020] Where: Δw represents the angular frequency deviation at the unit outlet, ΔP m∑ , ΔP e∑ are the total mechanical power increment and load power increment of the system respectively, and D is the equivalent damping coefficient of the system.

[0021] (2) Performing Laplace transform on the equivalent rotor motion equation, the transfer function of the corresponding system can be obtained as:

[0022]

[0023] Where: Δw(s) is the angular frequency deviation in the frequency domain, S is the Laplace operator, and ΔP(s) is ΔP m∑ With ΔP e∑ After the difference is made, it is obtained by Laplace transform.

[0024] (3) Perform inverse Laplace transform on the transfer function, and its time domain form is:

[0025]

[0026] Obviously, g(t) decays exponentially, and its amplitude is 2H sys The negative reciprocal of .

[0027] (4) Discretize the transfer function in the time domain, replace t with NT, and discretize g(t) into a time series function. The transfer function becomes:

[0028]

[0029] (5) Perform Z-transformation on the discrete time series function g(NT) to obtain the following formula:

[0030]

[0031] Where: Z represents a discrete variable, and T is a unit time period.

[0032] Because the output signal sequence, when solved in the time domain, is equal to the convolution sum of the input signal sequence and the unit sample response sequence of the system being used, tedious convolution sum calculations are inevitably required if the input signal sequence is to be the output signal sequence after being processed by a certain system. However, the convolution sum characteristics of the Z transform can greatly simplify this process. Simply calculate the Z transform of the input signal sequence and the unit sample response sequence of the system, and then calculate the inverse Z transform of the product of the two to obtain the output signal sequence.

[0033] Furthermore, in step 6, to obtain the parameters H and D to be determined, the African vulture algorithm AVOA includes the following steps:

[0034] (1) Initialize the African vulture population RP, set the total number of African vultures to N, the maximum number of iterations to maxiterations, and randomly initialize the positions of the African vulture groups. After the initial population is formed, calculate the fitness of all solutions, select the best solution as the best vulture of the first group, select the second best solution as the best vulture of the second group, and move other solutions to the best solutions of the first and second groups using the following formula.

[0035]

[0036] Where: L1 and L2 are parameters given before the search operation, their values are between 0 and 1, and the sum of the two parameters is 1. Use the following formula to obtain the probability of selecting the best solution, and select each best solution for each group.

[0037]

[0038] (2) Use the following formula to mathematically model the hunger rate of vultures.

[0039]

[0040]

[0041] Among them: F represents the vulture's satiety rate, Iteration i represents the current number of iterations, maxiterations represents the maximum number of iterations, z is a random number between -1 and 1 that changes with each iteration, h is a random number between -2 and 2, and rand1 is a random number between 0 and 1. When the z value drops below 0, it means the vulture is hungry, and if the z value increases to 0, it means the vulture is full.

[0042] The proportion of vultures to the total number of vultures is decreasing, and the decrease is larger with each iteration. When the value of F is greater than 1, the vultures search for food in different areas and AVOA enters the exploration phase; if the value of F is less than 1, AVOA enters the exploitation phase and the vultures search for food near the optimal solution.

[0043] (3) Vultures use the following formula to explore:

[0044]

[0045] Where: rand P1 is a random number between [0, 1]. P1 is a preset exploration parameter used to control the exploration strategy. P(i+1) is the vulture position vector in the next iteration. F is the vulture satiety rate obtained in the current iteration. R(i) is the best vulture. rand2 and rand3 are both random numbers between [0, 1]. lb and ub are the upper and lower bounds of the optimization, respectively.

[0046] (4) When the value of F is between 0.5 and 1, AVOA enters the first stage of the development phase. In the first stage, two different foraging strategies are executed, namely, rotation flight and siege strategy. The choice of strategy is based on the value of P2. The specific process is shown as follows:

[0047]

[0048] Where: rand P2 , rand4, rand5 are all random numbers between [0,1], F is the satiety rate of the vulture obtained in the iteration, and R(i) is one of the best vultures.

[0049] (5) When the value of F is less than 0.5, the second stage of the algorithm is executed. The actions of the two vultures gather several types of vultures at the food source and start an aggressive struggle to besiege and compete for food. Different strategies are selected based on P3. The specific process is as follows:

[0050]

[0051] Where: BestVulture1(i) is the best vulture in the first group of the current iteration, BestVulture2(i) is the best vulture in the second group of the current iteration, F is the vulture's satiety rate, P(i) is the vulture's current position vector, Levy(d) is the Levy flight mechanism used to improve the effectiveness of the algorithm, and d(t) represents the distance between the vulture and the best vulture in the two groups.

[0052] Beneficial effects of the present invention:

[0053] Compared to building a dynamic error model and then using system identification methods to estimate the system's inertia constant, this paper proposes an African vulture algorithm for inertia constant estimation, providing a new approach for evaluating the inertia constant of new power systems. This method can estimate the inertia time constant over a specific time period based on data points measured by the PMU. It can quickly track local trends in noise-like signals and better track changes in the overall system inertia. BRIEF DESCRIPTION OF THE DRAWINGS

[0054] Figure 1 is a flow chart of the present invention;

[0055] Figure 2 This is the topological structure diagram of the four-machine two-zone system;

[0056] Figure 3 Optimize the iteration curve for the African vulture algorithm. DETAILED DESCRIPTION

[0057] The technical solutions in the embodiments of the present invention will be described clearly and completely below with reference to the accompanying drawings in the embodiments of the present invention.

[0058] like Figure 1 、 2 As shown in FIG3 , an implementation example of the present invention discloses a microgrid inertia constant evaluation method based on the African vulture algorithm, comprising the following steps:

[0059] Step 1: Compare the new energy power source to a synchronous generator, equate it to a virtual synchronous generator, and propose the concept of virtual inertia.

[0060] Step 2: Construct a unified inertia constant analytical model including the VSG virtual inertia and the synchronous generator rotational inertia.

[0061] Step 3: The rotor motion equation of the unified inertia constant analytical model constructed in step 2 is subjected to Laplace transform, inverse Laplace transform, discretization, and Z transform in sequence, and finally a transfer function expression containing the unknown parameters H and D is obtained, where H is the inertia constant and D is the damping coefficient.

[0062] Step 4: The system operating power change data ΔP measured by the PMU o Substituting this into the calculated transfer function expression, we can obtain the corresponding frequency change Δf.

[0063] Step 5: Construct the fitness function y = ∑(Δf-Δf o ).

[0064] Step 6: Use the African vulture algorithm to optimize the parameters H and D to be determined. When the fitness function value y is within the required range, H and D calculated at that moment are taken as the parameters to be determined.

[0065] Furthermore, the unified inertia constant analytical model containing the VSG virtual inertia and the synchronous generator rotational inertia in step 2 is derived as follows:

[0066] According to the conservation of the overall inertia of the power system, the equivalent inertia constant of the power system can be calculated by the following formula:

[0067]

[0068] Among them: H gen,i represents the inertia constant of the i-th synchronous generator, S gen,i represents the rated capacity of the i-th synchronous generator, H vir,j represents the inertia constant of the jth VSG, S vir,j represents the rated capacity of the jth VSG, S sys Indicates the rated capacity of the entire system.

[0069] Furthermore, in step 3, to obtain a transfer function expression containing the parameters H and D to be determined, the following steps are included:

[0070] Step 3.1, use the equivalent rotor motion equation to express the transient electromechanical characteristics of the system:

[0071]

[0072] Where: Δw represents the angular frequency deviation at the unit outlet, ΔP m∑ , ΔP e∑ are the total mechanical power increment and load power increment of the system respectively, and D is the equivalent damping coefficient of the system.

[0073] Step 3.2, perform Laplace transform on the equivalent rotor motion equation, and the transfer function of the corresponding system can be obtained as:

[0074]

[0075] Where: Δw(s) is the angular frequency deviation in the frequency domain, S is the Laplace operator, and ΔP(s) is ΔP m∑ With ΔP e∑ After the difference is made, it is obtained by Laplace transform.

[0076] Step 3.3, perform inverse Laplace transform on the transfer function, and its time domain form is:

[0077]

[0078] Obviously, g(t) decays exponentially, and its amplitude is 2H sys The negative reciprocal of .

[0079] In step 3.4, the transfer function in the time domain is discretized. Replace t with NT and discretize g(t) into a time series function. The transfer function becomes:

[0080]

[0081] Step 3.5: Perform Z-transformation on the discrete time series function g(NT) to obtain the following formula:

[0082]

[0083] Where: Z represents a discrete variable, and T is a unit time period.

[0084] Because the output signal sequence, when solved in the time domain, is equal to the convolution sum of the input signal sequence and the unit sample response sequence of the system being used, tedious convolution sum calculations are inevitably required if the input signal sequence is to be the output signal sequence after being processed by a certain system. However, the convolution sum characteristics of the Z transform can greatly simplify this process. Simply calculate the Z transform of the input signal sequence and the unit sample response sequence of the system, and then calculate the inverse Z transform of the product of the two to obtain the output signal sequence.

[0085] Furthermore, in step 6, to obtain the parameters H and D to be determined, the following steps are included:

[0086] Step 6.1, initialize the African vulture population RP, set the total number of African vultures to N, the maximum number of iterations to maxiterations, randomly initialize the positions of the African vulture group, and after the initial population is formed, calculate the fitness of all solutions, select the best solution as the best vulture of the first group, select the second best solution as the best vulture of the second group, and move other solutions to the best solutions of the first and second groups using the following formula.

[0087]

[0088] Where: L1 and L2 are parameters given before the search operation, their values are between 0 and 1, and the sum of the two parameters is 1. Use the following formula to obtain the probability of selecting the best solution, and select each best solution for each group.

[0089]

[0090] In step 6.2, the hunger rate of vultures is mathematically modeled using the following formula.

[0091]

[0092]

[0093] Among them: F represents the vulture's satiety rate, Iteration i represents the current number of iterations, maxiterations represents the maximum number of iterations, z is a random number between -1 and 1 that changes with each iteration, h is a random number between -2 and 2, and rand1 is a random number between 0 and 1. When the z value drops below 0, it means the vulture is hungry, and if the z value increases to 0, it means the vulture is full.

[0094] The proportion of vultures to the total number of vultures is decreasing, and the decrease is larger with each iteration. When the value of F is greater than 1, the vultures search for food in different areas and AVOA enters the exploration phase; if the value of F is less than 1, AVOA enters the exploitation phase and the vultures search for food near the optimal solution.

[0095] In step 6.3, the vulture performs exploration using the following formula:

[0096]

[0097] Where: rand P1 is a random number between [0, 1]. P1 is a preset exploration parameter used to control the exploration strategy. P(i+1) is the vulture position vector in the next iteration. F is the vulture satiety rate obtained in the current iteration. R(i) is the best vulture. rand2 and rand3 are both random numbers between [0, 1]. lb and ub are the upper and lower bounds of the optimization, respectively.

[0098] In step 6.4, when the value of F is between 0.5 and 1, AVOA enters the first stage of the development phase. In the first stage, two different foraging strategies are implemented: rotation flight and siege strategy. The choice of strategy is based on the value of P2. The specific process is as follows:

[0099]

[0100] Where: rand P2 , rand4, rand5 are all random numbers between [0,1], F is the satiety rate of the vulture obtained in the iteration, and R(i) is one of the best vultures.

[0101] In step 6.5, if the value of F is less than 0.5, the second phase of the algorithm is executed. The actions of the two vultures attract several types of vultures to the food source, and they engage in a siege and aggressive struggle for food. Different strategies are selected based on P3. The specific process is as follows:

[0102]

[0103] Where: BestVulture1(i) is the best vulture in the first group of the current iteration, BestVulture2(i) is the best vulture in the second group of the current iteration, F is the vulture's satiety rate, P(i) is the vulture's current position vector, Levy(d) is the Levy flight mechanism used to improve the effectiveness of the algorithm, and d(t) represents the distance between the vulture and the best vulture in the two groups.

[0104] The algorithm parameters in this implementation example are as follows:

[0105] AVOA algorithm: the vulture group size is N=30, the dimension Dim=10, and the maximum number of iterations Max_iter=100.

[0106] The African vulture algorithm (AVOA) proposed in the present invention is used to estimate the inertia constant of the microgrid. Compared with the traditional method of system parameter identification, the method adopted by the present invention is simpler to calculate and has higher optimization accuracy.

[0107] Finally, it should be noted that the above description only describes specific embodiments of the present invention in detail. However, the present invention is not limited to the specific embodiments described above. Equivalent modifications and substitutions made by those skilled in the art are also within the scope of the present invention. Therefore, equivalent changes and modifications made without departing from the spirit and scope of the present invention are encompassed within the scope of the present invention.

Claims

1. A microgrid inertia constant estimation method based on the African vulture algorithm, characterized in that: The following steps are involved: Step 1: Compare the new energy power source to a synchronous generator and make it equivalent to a virtual synchronous generator; Step 2: Construct a unified inertia constant analytical model including the virtual inertia of the VSG and the rotational inertia of the synchronous generator; Step 3: Perform Laplace transform, inverse Laplace transform, discretization, and Z transform on the rotor motion equation of the unified inertia constant analytical model constructed in step 2, and finally obtain a transfer function expression containing the unknown parameters H and D, where H is the inertia constant and D is the damping coefficient; Step 4: The system operating power change data measured by the PMU Substituting into the calculated transfer function expression, the corresponding frequency change is obtained ; Step 5: Construct fitness function ; Step 6: Use the African vulture algorithm to optimize the parameters H and D to be determined. When the fitness function value y is within the required range, the H and D calculated at that moment are taken as the parameters to be determined. The unified inertia constant analytical model containing the VSG virtual inertia and the synchronous generator rotational inertia in step 2 is derived as follows: According to the conservation of the overall inertia of the power system, the equivalent inertia constant of the power system is calculated by the following formula: ； in: represents the inertia constant of the i-th synchronous generator, represents the rated capacity of the i-th synchronous generator, represents the inertia constant of the jth VSG, represents the rated capacity of the jth VSG, Indicates the rated capacity of the entire system; In step 3, to obtain the transfer function expression containing the parameters H and D to be determined, the following steps are included: Step 3.1, use the equivalent rotor motion equation to express the transient electromechanical characteristics of the system: ； Where: Indicates the angular frequency deviation at the unit outlet, 、 are the total mechanical power increment and load power increment of the system respectively, is the equivalent damping coefficient of the system; Step 3.2, perform Laplace transform on the equivalent rotor motion equation, and the transfer function of the corresponding system is: ； Where: is the angular frequency deviation in frequency domain, S is the Laplace operator, for and After the difference is made and Laplace transform is performed, we can get: Step 3.3, perform inverse Laplace transform on the transfer function, and its time domain form is: ； Obviously, It decays exponentially, and its amplitude is The negative reciprocal of Step 3.4, discretize the transfer function in the time domain and use replace ,Will Discretized into a time series function, the transfer function becomes: ； Step 3.5, the time series function obtained by discretization Perform Z transform and get the following formula: ； Where: Z represents a discrete variable, and T is a unit time period.

2. The microgrid inertia constant estimation method based on the African vulture algorithm according to claim 1 is characterized in that: In step 6, to obtain the parameters H and D to be determined, the African vulture algorithm AVOA includes the following steps: Step 6.1: Initialize the African vulture population RP, set the total number of African vultures to N, the maximum number of iterations to maxiterations, and randomly initialize the positions of the African vulture groups. After the initial population is formed, calculate the fitness of all solutions, select the best solution as the best vulture of the first group, select the second best solution as the best vulture of the second group, and move other solutions toward the best solution of the first and second groups using the following formula; ； Where: L1 and L2 are parameters given before the search operation, their values are between 0 and 1, and the sum of the two parameters is 1; the probability of selecting the best solution is obtained using the following formula, and each best solution is selected for each group; ； Step 6.2, mathematically model the vulture's hunger rate using the following formula; ； ； Among them: F represents the vulture's satiety rate, Indicates the current iteration number, represents the maximum number of iterations, z is a random number between -1 and 1 that changes with each iteration, h is a random number between -2 and 2, is a random number between 0 and 1; when the z value drops below 0, it means the vulture is hungry, and if the z value increases to 0, it means the vulture is full; The proportion of the total number of vultures is decreasing, and the decrease is greater with each iteration; when the value of F is greater than 1, the vultures search for food in different areas, and AVOA enters the exploration phase; if the value of F is less than 1, AVOA enters the exploitation phase, and the vultures search for food near the optimal solution; In step 6.3, the vulture performs exploration using the following formula: ； in: is a random number between [0,1], P1 is the preset exploration parameter used to control the exploration strategy; P(i+1) is the vulture position vector in the next iteration, F is the vulture satiety rate obtained in the current iteration, and R(i) is one of the best vultures; and are all random numbers between [0,1]; lb and ub are the upper and lower bounds of the optimization respectively; In step 6.4, when the value of F is between 0.5 and 1, AVOA enters the first stage of the development phase. In the first stage, two different foraging strategies are executed, namely, rotating flight and siege strategy. The choice of strategy is based on the value of P2. The specific process is shown as follows: ； in: , , are all random numbers between [0,1], F is the satiety rate of the vulture obtained in the iteration, and R(i) is one of the best vultures; In step 6.5, when the value of F is less than 0.5, the second stage of the algorithm is executed. The actions of the two vultures gather several types of vultures at the food source, and they engage in a siege and aggressive struggle for food. Different strategies are selected based on P3. The specific process is as follows: ； in: is the best vulture in the first group of the current iteration, is the best vulture in the second group of the current iteration, F is the vulture's satiety rate, P(i) is the vulture's current position vector, Levy(d) is the Levy flight mechanism used to improve the effectiveness of the algorithm, and d(t) represents the distance between the vulture and the best vulture in the two groups.

Images