A transient stability constrained optimal power flow calculation method for wind power access to power system

Abstract

The transient stability constrained optimal power flow calculation method for wind power access to the power system, establishes a mathematical model of wind power output based on Weibull distribution; generates a large number of wind power scenarios based on Latin hypercube sampling; reduces the generated large number of wind power scenarios based on the synchronous back-substitution reduction method; constructs the transient stability constraints of the power system based on the input convex neural network; adds the constructed transient stability constraints of the power system to the optimal power flow model, and constructs the transient stability constrained optimal power flow model containing wind power access to the power system; uses the search space reduction algorithm to solve and calculate the constructed transient stability constrained optimal power flow model. The present invention provides a transient stability constrained optimal power flow calculation method for wind power access to the power system, which effectively guarantees the stability of the power system containing wind power access.

Description

Technical Field

[0001] The present invention relates to the technical field of transient stability calculation of power systems, and in particular to a transient stability constrained optimal power flow calculation method for wind power access to a power system. Background Art

[0002] As new energy enters a rapid development stage, the impact of new energy access on grid stability must be considered. Wind power output is characterized by uncertainty, and its large-scale access poses a higher challenge to the stable operation of the power system. Therefore, the transient stability constrained optimal power flow calculation method taking into account wind power uncertainty is of great significance.

[0003] In the technical solution disclosed in the Chinese patent "A trusted data-driven solution method for optimal power flow with strong generalization ability" (application number: 202310570200.X), the fixed relationship model of the optimal power flow output and input is optimized to obtain the gradient relationship model of the optimal power flow output and input, and the input data is decoded to obtain the optimal power flow input characteristics. In the technical solution disclosed in the Chinese patent "A data-driven optimal power flow calculation method considering topological feature learning" (application number: 202211690444.3), this method constructs a graph convolutional neural network based on the optimal power flow model, and trains the graph convolutional neural network to obtain a neural network model for optimal power flow calculation. In the technical solution disclosed in the Chinese patent "Optimal power flow calculation method and device for flexible DC traction power supply system" (application number: 202211384111.8), this method can be compatible with both centralized and decentralized power supply systems, taking into account more complete system operation feature constraints and typical fault sets, and proposes multiple optimization objectives to meet different needs.

[0004] However, in the above-disclosed prior art, there is a lack of a simple and practical optimal power flow calculation method that takes transient stability constraints into consideration for a power system with wind power access. Summary of the invention

[0005] To solve the above technical problems, the present invention considers the influence of uncertainty of wind power access in the transient stability constrained optimal power flow model, and proposes a transient stability constrained optimal power flow calculation method for wind power access power system, which effectively ensures the stability of the power system containing wind power access.

[0006] The technical solution adopted by the present invention is:

[0007] The transient stability constrained optimal power flow calculation method for wind power access to the power system includes the following steps:

[0008] Step 1: Establish a mathematical model of wind power output based on Weibull distribution;

[0009] Step 2: Generate a large number of wind power scenarios based on Latin hypercube sampling;

[0010] Step 3: Reduce the large number of wind power scenarios generated in step 2 based on the synchronous back generation reduction method;

[0011] Step 4: Construct power system transient stability constraints based on the input convex neural network;

[0012] Step 5: Add the transient stability constraint of the power system constructed in step 4 to the optimal power flow model to construct the transient stability constraint optimal power flow model with wind power access to the power system;

[0013] Step 6: Use the search space reduction algorithm to solve the transient stability constrained optimal power flow model constructed in step 5.

[0014] In step 1, the Weibull wind power modeling method is a method for analyzing the relationship between wind speed and power in a power system including wind power access. It describes the wind speed probability density function and probability distribution function based on the Weibull distribution. By analyzing the historical wind speed data, it is possible to obtain parameters that conform to the Weibull distribution, and use these parameters to establish the probability density function of the wind speed, as follows:

[0015]

[0016] In formula (1), v is the wind speed; h and z are the shape parameter and scale parameter of the Weibull distribution, respectively.

[0017] The output model of wind power generation using wind speed is as follows:

[0018]

[0019] In formula (2): P ω is the output of the wind turbine; P E is the rated power of wind power; v, v β 、v R and v ο They are actual wind speed, cut-in wind speed, rated wind speed and cut-out wind speed respectively.

[0020] The step 2 comprises the following steps:

[0021] Step 2.1: Cumulative probability distribution function G of wind power generation N As shown in formula (3):

[0022]

[0023] In formula (3): G N is the cumulative probability distribution function of wind power; v R 、v β and vo are rated wind speed, cut-in wind speed and cut-out wind speed respectively; p ω is the output of the wind turbine; p r is the maximum output of the wind turbine; h and z are the shape parameter and scale parameter of the Weibull distribution respectively.

[0024] First, determine the sampling scale M and set G N The distribution range of is divided into M intervals, that is,

[0025] A sample is drawn in each sampling interval, and the cumulative probability of the sample in the lth interval is δ l The calculation formula is:

[0026] δ l =(1 / M)×m i +(l-1) / M (4);

[0027] In formula (4): l is the cumulative probability of the sample in the lth interval; m i is a random number between the interval [0,1]; M is the total number of intervals; l is the number of intervals.

[0028] Step 2.2: The calculated cumulative probability δ l Substitute G N The inverse function G N -1 In the above equation, we get the corresponding sampling value ω l ,Right now:

[0029] ω l =G N -1 (δ l ) (5);

[0030] In formula (5): l is the sampling value; δ l is the cumulative probability of the sample in the lth interval; is the inverse cumulative probability distribution function.

[0031] Step 2.3: Repeat the above steps M times to obtain M sampling values, and randomly sort the sampling values ​​to obtain M scenes.

[0032] The step 3 comprises the following steps:

[0033] Step 3.1: Based on the synchronous back-substitution method of Kantorovich distance, select the scene with the smallest Kantorovich distance from the original scene set Ω. The Kantorovich distance is defined as follows:

[0034]

[0035] In formula (6): C k (s i ,s j ) is scene s i and j Kantorovich distance; Ω is the reduced scene set; s i and j is a scene in the scene set Ω; d(s i ,s j ) is a non-negative, continuous, and symmetric distance function; ω(s i ,s j ) is scene s i and j The probability product of i and ρ j Respectively for scenes i and j The probability in Ω; T is the number of time periods in a day; Respectively for scenes i and j The power during period t.

[0036] Step 3.2: Eliminate scenes that meet the following conditions w :

[0037]

[0038]

[0039] In formula (7) and (8): C k (s w ,s i ) is scene s w and i Kantorovich distance; s w is the scenario with the closest probability distance to other scenarios and the worst representativeness; i is other scenes; M is the number of scenes; ρ i and ρ w Respectively for scenes i and w The probability of d(s w ,s i ) is the distance between two scenes; T is the number of time periods in a day; and Respectively for scenes w and i The power during period t.

[0040] Step 3.3: Change the number of original scenes and remove the distance scene sw Recent scenes w′ :

[0041]

[0042]

[0043] In formula (9) and (10): is the number of scenes after the change; M is the number of original scenes; s w For the scenes to be eliminated; w′ is the distance scene s w Recent scenes; d(s w ,s w′ ) is scene s w and w′ The distance between i For the scene set Ω, excluding scene s w Other scenes.

[0044] Step 3.4: Change and cull scenes w The closest scene w′ The probability of , so that the sum of the probabilities of all remaining scenarios is 1:

[0045]

[0046] In formula (11): For the updated scene w′ The probability of w′ To update the scene before w′ The probability of w For scenes w probability.

[0047] Step 3.5: Go to step 1, loop and iterate the calculation until the number of remaining scenarios meets the set target number of scenarios, and output the scenario results. The schematic diagram of the synchronous back-substitution method based on Kantorovich distance is as follows: Figure 2 shown.

[0048] In step 4, the mathematical expression of the input convex neural network ICNN is:

[0049]

[0050] In formula (12), a is the input layer vector; b is the hidden layer vector of a; ξ and β are the weights and biases of different layers respectively; h is the activation function; b i+1 is the i+1th hidden layer vector of a; is the weight of the i-th hidden layer; b i is the i-th hidden layer vector of a; is the weight of the input layer a in the i-th hidden layer; β i is the bias in the i-th hidden layer.

[0051] Based on the proposed input convex neural network ICNN model to predict the transient stability of the power system:

[0052] Input convex neural network ICNN input a = [P G ,Q G ,P L ,Q L ], where P G and Q G are the active power and reactive power of the generator respectively; P L and Q L They are active power demand and reactive power demand respectively;

[0053] The output is based on the transient stability index (TSI), which reflects the maximum value of the rotor angle difference of each generator during the transient period. The TSI formula is as follows:

[0054]

[0055] In formula (13): max is the maximum rotor angle difference between any two generators.

[0056] Based on ICNN, the transient stability constraint can be expressed as:

[0057] TSI>TSI B (14);

[0058] In formula (14): TSI B is the transient stability boundary TSI of the power system obtained by the ICNN model, which is used to distinguish transient stability from instability.

[0059] In step 5, the transient stability constraint is added to the optimal power flow model to construct a transient stability constrained optimal power flow (TSCOPF) model for a power system with wind power access. The TSCOPF model is a tool for coordinating operating costs and stability requirements through rescheduling. The transient stability constrained optimal power flow model for a power system with wind power access consists of four parts: an objective function that minimizes the total fuel cost of system power generation, an equality constraint, an inequality constraint, and a transient stability constraint.

[0060] The transient stability constrained optimal power flow model is as follows:

[0061] 1) Objective function:

[0062]

[0063] In formula (15): best is the objective function, i.e. the minimum total cost of power generation; m i 、n i 、h i The generator G i Cost coefficient of For the generator G i Active output; S G Represents all generators.

[0064] 2) Equality constraints:

[0065]

[0066] In formula (16): and The generator G i Active output and reactive output; and are the active power demand and reactive power demand of load i respectively; is the voltage amplitude of node i; S A is the total number of nodes; θ ij is the phase angle difference between nodes i and j; G ij and B ij are the conductance and susceptance between node i and node j respectively.

[0067] 3) Inequality constraints:

[0068]

[0069] In formula (17): The generator G i Active output and reactive output; are the lower and upper limits of the active output of generator i respectively; are the lower and upper limits of reactive output of generator i respectively; V i , are the lower and upper limits of the voltage at node i respectively; V i is the voltage of node i; S li and are the apparent power and apparent power upper limit of node i respectively; n G 、n m and n l They are respectively the generator set, node set and line set of the power system.

[0070] 4) Transient stability constraints:

[0071] TSI>TSIB (18).

[0072] In step 6, the TSCOPF model is solved and calculated based on the Search Space Reduction Algorithm (SSRA); SSRA is a population-based metaheuristic algorithm used to solve optimization problems. The main idea of ​​the algorithm is to effectively explore and utilize the search space, and then reduce the search space as the iteration proceeds. The TSCOPF model is solved and calculated based on the Search Space Reduction Algorithm, including the following steps:

[0073] Step 6.1: Based on the Gaussian penalty function method, an unconstrained objective function is formed;

[0074]

[0075]

[0076] In formula (19) and (20): is the penalty function; λ is the penalty factor; n is the number of variables; x i is the ith discrete variable; b i is the center of the ith neighborhood; c i is the i-th neighborhood factor; best is the objective function; b ' est is an unconstrained objective function.

[0077] Step 6.2: Define the general parameters of the optimization problem, such as defining the dimension as D; the population of search agents as pop; and the maximum number of iterations as n max ; The upper and lower boundaries of the search agent are XB and YB;

[0078] Step 6.3: Define algorithm-specific parameters, such as defining C as the number of optimal positions for evaluating midpoints; R f The reduction factor represents the speed at which the search space range is reduced in each iteration;

[0079] Step 6.4: Initialize the search range within the upper and lower boundaries of the first iteration;

[0080] Step 6.5: For the initial iteration, randomly generate a search population within the entire search range, and the search population position is V = [v1, v2…v pop ] T Among them, v i represents the position of the i-th search agent in the D-dimensional search space; v1,v2…v pop They respectively represent the positions of the 1st, 2nd, …, pop search agents in the D-dimensional search space; T is the transposition symbol.

[0081] Step 6.6: Evaluate each search agent position v i χ b ' est , χ b ' est , is the unconstrained objective function given in step 6.1;

[0082] Step 6.7: Sort the objective function value of each search agent in descending order Among them, χ best1 ,x best2 … bestpop They represent the minimum value of the objective function, the second smallest value of the objective function, …, and the maximum value of the objective function, respectively, and determine the C global optimal solutions V obtained in this iteration gbestC ;

[0083] Step 6.8: Use the C optimal solutions obtained to evaluate the midpoint, that is, calculate the C global optimal solutions V obtained in steps 6.4 to 6.7 gbestC The average value of

[0084] Step 6.9: Update the search range z, the formula is as follows:

[0085]

[0086] In formula (21): R f is the reduction factor, which determines the speed of search range reduction during the entire iteration process; V is the search population position; n max is the maximum number of iterations; z(n+1) is the search range after the n+1th iteration; z(n) is the search range after the nth iteration.

[0087] The search population is randomly regenerated in the updated search area based on the following relationship:

[0088]

[0089] In formula (22), i = 1, 2…pop, pop is the population of search agents; ψ1 and ψ2 represent random numbers in the interval [0, 1]; v i (n+1) is the search population after the n+1th iteration; mid(n) is the C global optimal solution V of the nth iteration gbestC The average value of; z(n) is the search range after the nth iteration.

[0090] Step 6.10: Verify that the regenerated search population is within the range of XB and YB:

[0091] v i (n) = YB + (XB - YB) * ψ, vi <YB or v i >XB (23);

[0092] In Equation (23): v i (n) is the search population after the nth iteration; XB and YB are the upper and lower bounds of the search agent respectively; ψ is a random number within the interval [0, 1]; v i is the search population.

[0093] Step 6.11: Use the updated search agent to go to Step 6.6 and loop again until a predetermined termination condition is reached. The present invention uses the number of iterations as the termination condition.

[0094] The transient stability constrained optimal power flow calculation method for a wind power integrated power system of the present invention has the following technical effects:

[0095] 1) TSCOPF is a non-convex optimization problem with a large computational amount. In order to find the global optimal solution and improve the convergence, the present invention introduces ICNN into the transient stability constrained optimal power flow model of the power system, which can accurately and quickly calculate the transient stability index of the power system boundary.

[0096] 2) The TSCOPF scheme based on SSRA of the present invention coordinates the requirements of system operation economy and stability by changing the active power of each generator set in the power system, and at the same time satisfies the constraint conditions. Compared with other optimization algorithms, SSRA has superior convergence characteristics, which can make the calculation speed faster and the accuracy higher when solving the TSCOPF model.

[0097] 3) The advantage of Step 1 of the present invention is that wind speed is one of the main factors affecting the output of wind power generation. The Weibull distribution can well describe the variation characteristics of wind speed, so as to more accurately simulate the wind power output, with easy-to-estimate parameters and strong reliability.

[0098] 4) In Step 2 of the present invention, the advantage of Latin hypercube sampling is that it can ensure that the distribution of samples in each dimension is uniform, thus avoiding the concentration and skewness of samples. In addition, Latin hypercube sampling usually requires fewer samples in a multi-dimensional space to achieve a similar coverage range.

[0099] 5) In Step 3 of the present invention, the synchronous back substitution method based on the Kantorovich distance is a powerful metric method, which can accurately quantify the differences between scenarios, can efficiently process a large number of wind power scenarios, and reduces the computational cost and time consumption.

[0100] 6) In step 4 of the present invention, the transient stability of the power system is evaluated based on the input convex neural network. The objective function of the input convex neural network is a convex function, so there is a global optimal solution, and the training process is more stable and not easy to fall into a local optimal solution.

[0101] 7) In step 5 of the present invention, compared with the traditional transient stability constraint based on time domain simulation, the data-driven transient stability constraint based on input convex neural network proposed by the present invention has the advantages of less computing resource consumption and short operation time.

[0102] 8) In step 6 of the present invention, SSRA is used to solve the TSCOPF model. SSRA has excellent convergence characteristics, which can speed up the calculation when solving the TSCOPF model. Compared with other algorithms, SSRA has higher accuracy. BRIEF DESCRIPTION OF THE DRAWINGS

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

[0104] Figure 1 It is a Latin hypercube diagram;

[0105] Figure 2 It is a schematic diagram of synchronous back-substitution method;

[0106] Figure 3 is a flow chart of the method of the present invention;

[0107] Figure 4 It is a flow chart of the solution algorithm of the present invention;

[0108] Figure 5 is a schematic diagram of the IEEE 39-node system of the present invention;

[0109] Figure 6 It is a power angle diagram of the IEEE 39-node system when a fault occurs in the present invention;

[0110] Figure 7 is a power angle diagram of the IEEE 39-node system of the present invention after being solved based on the optimization method of the present invention;

[0111] Figure 8 It is the fuel cost graph of SSRA and Grey Wolf optimization algorithm in IEEE 39-bus system. DETAILED DESCRIPTION

[0112] The present invention will be further described below in conjunction with the accompanying drawings and embodiments.

[0113] The present invention provides a transient stability constrained optimal power flow calculation method for wind power access to a power system, the method comprising the following steps: Figure 3 As shown:

[0114] Step 1: Modeling the uncertainty of wind power generation based on Weibull distribution:

[0115] The present invention mainly considers the uncertainty of wind power generation, and the wind speed probability density function is:

[0116]

[0117] Where: v is the wind speed; h and z are the shape parameter and scale parameter of the Weibull distribution respectively.

[0118] The output model of wind power generation using wind speed is as follows:

[0119]

[0120] Where: P E is the rated power of wind power; v, v β 、v R and v ο They are actual wind speed, cut-in wind speed, rated wind speed and cut-out wind speed respectively.

[0121] It can be seen from equations (1) and (2) that the estimation accuracy of the shape parameter h and scale parameter z of the Weibull distribution is very important in wind power uncertainty modeling.

[0122] Step 2: Generate a large number of wind power scenarios based on the Latin Hypercube Sampling (LHS) method. The specific steps are as follows:

[0123] Step 2-1: Based on the idea of ​​stratified sampling, the value range of the probability distribution function is evenly divided into n sub-areas;

[0124] Step 2-2: Conduct a random sampling in each sub-area;

[0125] Step 2-3: Obtain the wind power output situation in reverse based on the sampling values ​​and generate a large number of scenarios.

[0126] Step 3: Reduce the large number of generated scenes based on the synchronous back generation reduction method. The specific steps are as follows:

[0127] Step 3-1: Calculate the probability distance between each scene and eliminate the scene with the smallest probability distance;

[0128] Step 3-2: Change the total number of scenarios and the corresponding probabilities of the scenarios;

[0129] Step 3-3: Repeat the above steps until the number of retained scenes is reduced to the specified number.

[0130] Step 4: Construct power system transient stability constraints based on the Input Convex Neural Networks (ICNN) method.

[0131] ICNN can parameterize the relationship between output and input features in the form of nonlinear convex functions. ICNN can effectively avoid the model convergence difficulties or local optimal problems caused by training neural networks with non-convex boundary constraints. Since ICNN has the hyperparametric boundary convex property, objects modeled based on ICNN can obtain the global optimal solution through a fast gradient decision algorithm.

[0132] The mathematical expression of ICNN is:

[0133]

[0134] f(a;θ)=b;

[0135] Where a is the input layer vector; b is the hidden layer vector of a; ξ and β are the weights and biases of different layers respectively; h is the activation function; f is the ICNN model; θ is all the parameters in ICNN.

[0136] The transient stability of the power system is predicted based on the proposed ICNN model. The input of ICNN is a = [P G ,Q G ,P L ,Q L ], where P G and Q G are the active power and reactive power of the generator respectively; P L and Q L The output is based on the transient stability index (TSI), which reflects the maximum value of the rotor angle difference of each generator during the transient period. The formula is as follows:

[0137]

[0138] Where: max is the maximum rotor angle difference between any two generators.

[0139] Based on ICNN, the transient stability constraint can be expressed as:

[0140] TSI>TSI B (14);

[0141] Where: TSI B is the transient stability boundary TSI of the power system obtained by the ICNN method, which is used to distinguish between transient stability and instability.

[0142] Step 5: Add the transient stability constraint of the power system to the optimal power flow model and construct a transient stability constrained optimal power flow (TSCOPF) model with wind power access to the power system. The specific steps are as follows:

[0143] 1) Objective function:

[0144]

[0145] Where: m i 、n i 、h i The generator G i Cost coefficient of For the generator G i Active output; S G Represents all generators.

[0146] 2) Equality constraints:

[0147]

[0148] Where: and The generator G i Active output and reactive output; and are the active power demand and reactive power demand of load i respectively; is the voltage amplitude of node i; S A is the total number of nodes; θ ij is the phase angle difference between nodes i and j; G ij and B ij are the conductance and susceptance between node i and node j respectively.

[0149] 3) Inequality constraints:

[0150]

[0151] Where: The generator G i Active output and reactive output; are the lower and upper limits of the active output of generator i respectively; are the lower and upper limits of reactive output of generator i respectively; V i , are the lower and upper limits of the voltage at node i respectively; S li and are the apparent power and apparent power upper limit of node i respectively; n G 、n m and nl They are respectively the generator set, node set and line set of the power system.

[0152] 4) Transient stability constraints:

[0153] TSI>TSI B (18);

[0154] In step 6: the TSCOPF model is solved and calculated based on the Search Space Reduction Algorithm (SSRA):

[0155] SSRA is a population-based metaheuristic algorithm for solving optimization problems. The main idea of ​​the algorithm is to effectively explore and utilize the search space and narrow the search space as it iterates. In each iteration, SSRA is different from previous population-based metaheuristic algorithms in that it gradually narrows the search space and randomly creates search agents within the narrowed search range. This SSRA idea helps to get rid of the dilemma of local optimality. SSRA has the advantages of simple implementation and low computational complexity. The TSCOPF algorithm based on SSRA includes the following steps, such as Figure 4 As shown:

[0156] Step 6-1: Based on the penalty function method, change the form of the TSCOPF model;

[0157] Step 6-2: Define the general parameters of the optimization problem, such as defining the dimension as D; the population of search agents as pop; and the maximum number of iterations as n max ; The upper and lower boundaries of the search agent are XB and YB;

[0158] Step 6-3: Define algorithm-specific parameters, such as defining C as the number of optimal positions for evaluating midpoints; R f The reduction factor represents the speed at which the search space range is reduced in each iteration;

[0159] Step 6-4: Initialize the search range within the upper and lower boundaries of the first iteration. The initial search range is as follows:

[0160]

[0161] Where: XB and YB are the upper and lower bounds of the search agent respectively;

[0162] Step 6-5: For the initial iteration, randomly generate a search population within the entire search range. The search population position is V = [v1, v2…v pop ] T Among them, v i represents the position of the i-th search agent in the D-dimensional search space;

[0163] Step 6-6: Evaluate each search agent position v i The cost function of

[0164] Step 6-7: Sort the objective function values ​​of each search agent in descending order to determine the C global optimal solutions V obtained in this iteration gbestC ;

[0165] Step 6-8: Use the C optimal solutions obtained to evaluate the midpoint, that is, calculate the C global optimal solutions V obtained in steps 4-7 gbestC The average value of

[0166] Step 6-9: Update the search range z, the formula is as follows:

[0167]

[0168] Where: R f is the reduction factor, which determines the speed of search range reduction during the entire iteration process; V is the search population position; n max is the maximum number of iterations;

[0169] The search population is randomly regenerated in the updated search area based on the following relationship:

[0170]

[0171] Where: i = 1, 2…pop; ψ1 and ψ2 represent random numbers in the interval [0, 1];

[0172] Step 6-10: Verify that the regenerated search population is within the range of XB and YB;

[0173] Step 6-11: Go to step 6-6 using the updated search agent and re-loop until a predetermined termination condition is reached.

[0174] Verification example:

[0175] To verify the effectiveness of the present invention, a test was conducted on an IEEE 39-node system. The schematic diagram of the IEEE 39-node system is shown in FIG. Figure 5 In the IEEE 39-node system, the synchronous generators at nodes 34 and 38 are replaced by wind turbines W6 and W10, and the power generation of the wind turbines is the same as that of the original generators.

[0176] In a typical wind power scenario, a three-phase short circuit occurs on line 10-13 in the IEEE 39-bus system, and the fault duration is 0.25s. Figure 6 As shown. Based on SSRA, the system is processed, and the generator power angle after the system is stable is as follows Figure 7 As shown. Figure 7 It can be seen that after the system operation mode is optimized based on the TSCOPF model proposed in the present invention, although the generator power angle undergoes certain disturbances, it eventually tends to be flat, and the maximum power angle does not exceed the threshold of 360 degrees, indicating that the system maintains a transient stable state.

[0177] In order to understand the advantages of SSRA over existing algorithms, the TSCOPF model is calculated based on SSRA and Grey Wolf Optimization Algorithm (GWOA). Figure 8 As shown, compared with GWOA, the fuel cost and iteration number of the power system transient stability constrained optimal power flow model obtained based on SSRA are less, which shows that the SSRA proposed in the present invention has the characteristics of fast convergence.

Claims

1. The transient stability constrained optimal power flow calculation method for wind power access to power system is characterized by The following steps are involved: Step 1: Establish a mathematical model of wind power output based on Weibull distribution; Step 2: Generate a large number of wind power scenarios based on Latin hypercube sampling; Step 3: Reduce the large number of wind power scenarios generated in step 2 based on the synchronous back generation reduction method; Step 4: Construct power system transient stability constraints based on the input convex neural network; Step 5: Add the transient stability constraint of the power system constructed in step 4 to the optimal power flow model to construct the transient stability constraint optimal power flow model with wind power access to the power system; Step 6: Use the search space reduction algorithm to solve the transient stability constrained optimal power flow model constructed in step 5; In step 4, the mathematical expression of the input convex neural network ICNN is: （12）； In formula (12): is the input layer vector; yes The hidden layer vector of ; and are the weights and biases of different layers respectively; is the activation function; yes No. layer hidden layer vector; For the The weights of the hidden layers; yes No. layer hidden layer vector; For the Layer Hidden Layer Input Layer The weight of For the Biases in the hidden layers; Input to the Convex Neural Network ICNN = ,in, and are the active power and reactive power of the generator respectively; and They are active power demand and reactive power demand respectively; Output based on transient stability index , which reflects the maximum value of the rotor angle difference of each generator during the transient period, The formula is as follows: （13）； In formula (13): is the maximum rotor angle difference between any two generators; Based on ICNN, the transient stability constraint can be expressed as: （14）； In formula (14): is the transient stability boundary of the power system obtained by the ICNN model , used to distinguish transient stability from instability; In step 5, the transient stability constraint is added to the optimal power flow model to construct a transient stability constraint optimal power flow model with wind power access to the power system. The transient stability constraint optimal power flow model is specifically as follows: 1) Objective function: （15）； In formula (15): The objective function is to minimize the total cost of power generation; , , Generator Cost coefficient of For generator Active output; represents all generators; 2) Equality constraints: （16）； In formula (16): and Generator Active output and reactive output; and Load Active power demand and reactive power demand; For Node The voltage amplitude of is the total number of nodes; For Node and The phase angle difference between and Node and nodes Conductance and susceptance between 3) Inequality constraints: （17）； In formula (17): , Generator Active output and reactive output; , Generator The lower and upper limits of active output; , Generator The lower and upper limits of reactive output; , Node The lower and upper voltage limits; For Node Voltage; and Node The apparent power and the upper limit of apparent power; , and They are the power system generator set, node set and line set respectively; 4) Transient stability constraints: （18）； In step 6, solving and calculating the TSCOPF model based on the search space reduction algorithm includes the following steps: Step 6.1: Based on the Gaussian penalty function method, an unconstrained objective function is formed; （19）； （20）； In formula (19) and (20): is the penalty function; is the penalty factor; is the number of variables; For the discrete variables; For the Neighborhood centers; For the Neighborhood factors; is the objective function; is the unconstrained objective function; Step 6.2: Define the general parameters of the optimization problem and define the dimension as ; The population of search agents is ; The maximum number of iterations is ; The upper and lower bounds of the search agent are and ; Step 6.3: Define algorithm-specific parameters, define is the number of optimal positions for evaluating the midpoint; The reduction factor represents the speed at which the search space range is reduced in each iteration; Step 6.4: Initialize the search range within the upper and lower boundaries of the first iteration; Step 6.5: For the initial iteration, randomly generate a search population within the entire search range, and use the search population position ;in, Respectively represent the 1st, 2nd, ..., Search agents in The position in the dimensional search space; is the transpose symbol; Step 6.6: Evaluate each search agent location of , ; Step 6.7: Sort the objective function value of each search agent in descending order ,in, They represent the minimum value of the objective function, the second smallest value of the objective function, …, and the maximum value of the objective function, respectively, to determine the value obtained in this iteration. Global optimal solution ; Step 6.8: Use the obtained The optimal solution evaluates the midpoint, that is, calculates the C global optimal solutions obtained in steps 6.4 to 6.7 The average value of Step 6.9: Update the search scope , the formula is as follows: （21）； In formula (21): is the reduction factor, which determines the speed at which the search range is reduced during the entire iteration; To search for population locations; is the maximum number of iterations; For the The search range after iterations; For the The search range after iterations; The search population is randomly regenerated in the updated search area based on the following relationship: （22）； In formula (22): , is the population size of the search agent; and Represents a random number in the interval [0,1]; For the Search population after iterations; For the Iterate C global optimal solutions The average value of For the The search range after iterations; Step 6.10: Verify that the regenerated search population is and Within the scope of: （23）； In formula (23): For the Search population after iterations; and are the upper and lower bounds of the search agent, respectively; For interval Random numbers within To search for population; Step 6.11: Go to step 6.6 with the updated search agent and re-loop until a predetermined termination condition is reached.

2. The transient stability constrained optimal power flow calculation method for wind power access to the power system according to claim 1 is characterized in that In step 1, a probability density function of wind speed is established based on Weibull distribution, as follows: （1）； In formula (1): is the wind speed; and are the shape parameter and scale parameter of the Weibull distribution, respectively; The output model of wind power generation using wind speed is as follows: （2）； In formula (2): is the output of the wind turbine; is the rated power of wind power; , , and They are actual wind speed, cut-in wind speed, rated wind speed and cut-out wind speed respectively.

3. The transient stability constrained optimal power flow calculation method for wind power access to the power system according to claim 1 is characterized in that : The step 2 comprises the following steps: Step 2.1: Cumulative probability distribution function of wind power generation As shown in formula (3): （3）； In formula (3): is the cumulative probability distribution function of wind power; , and They are rated wind speed, cut-in wind speed and cut-out wind speed respectively; is the output of the wind turbine; is the maximum output of the wind turbine; and are the shape parameter and scale parameter of the Weibull distribution, respectively; First, determine the sampling size ,Will The distribution range is divided into interval, that is , , …, ; In each sampling interval, a sample is drawn. The cumulative probability of the interval samples The calculation formula is: （4）； In formula (4): For the The cumulative probability of the interval sample; is a random number between the interval [0,1]; is the total number of intervals; is the interval number; Step 2.2: The calculated cumulative probability Substitution The inverse function of The corresponding sampling value is obtained ,Right now: （5）； In formula (5): is the sampling value; For the The cumulative probability of the interval sample; is the inverse cumulative probability distribution function; Step 2.3: Repeat the above steps times, get sampling values, sort the sampling values ​​randomly, and get A scene.

4. The transient stability constrained optimal power flow calculation method for wind power access to the power system according to claim 1 is characterized in that : Step 3 comprises the following steps: Step 3.1: Based on the Kantorovich distance, the synchronous back-substitution method is used to start from the original scene set Select the scene with the smallest Kantorovich distance from the others. The Kantorovich distance is defined as follows: （6）； In formula (6): For the scene and Kantorovich distance; To reduce the front scene set; and For scene set Scenes in; is a non-negative, continuous, and symmetric distance function; For the scene and The probability product of and Respectively for scenes and exist The probability of The number of time periods into which the day is divided; , Respectively for scenes and middle Power of the time period; Step 3.2: Eliminate scenes that meet the following conditions : （7）； （8）； In formula (7) and (8): For the scene and Kantorovich distance; The scenario with the closest probability distance to other scenarios and the worst representativeness; For other scenes; is the number of scenes; and Respectively for scenes and The probability of is the distance between the two scenes; The number of time periods into which the day is divided; and Respectively for scenes and middle Power of the time period; Step 3.3: Change the number of original scenes and remove the distance scenes Recent scenes : （9）； （10）； In formula (9) and (10): is the number of scenes after the change; is the number of original scenes; For the scenes to be eliminated; For distance scene Recent scenes; For the scene and The distance between For scene set In addition to the scene Other scenes of Step 3.4: Changing and culling scenes The closest scene The probability of , so that the sum of the probabilities of all remaining scenarios is 1: （11）； In formula (11): For the updated scene The probability of To update the scene before The probability of For the scene The probability of Step 3.5: Go to step 1 and iterate the calculation until the number of remaining scenarios meets the set target number of scenarios, and output the scenario results.

Images