Wind power plant dynamic inertia partition evaluation method based on improved clustering

Abstract

The invention provides a wind power plant dynamic inertia partition evaluation method based on improved clustering, and relates to the technical field of wind power plants. The method comprises the following steps: taking a Canopy clustering algorithm as a preprocessing step through improved Canopy-Kmeans; the mean value point of the collected data set is calculated, and the measuring point closest to the mean value point is selected as the initial clustering center, so that the result is more reasonable and stable. The grouping precision is improved by considering the wake effect and the time delay effect and measuring multi-dimensional wind turbine generator state variables, in order to avoid a large amount of data dimensionality disasters, PCA principal component analysis is adopted for data dimensionality reduction, components with the variance contribution rate being 85% or above are selected as dominant variables, and the complexity of data collection is reduced. According to the method, the improved Canopy algorithm and the k-means algorithm are combined, the optimal classification number k value and the clustering center can be recognized more accurately, and then large-area wind turbine generator partitioning is achieved quickly and accurately.

Description

Technical Field

[0001] The present invention relates to the technical field of wind farms, and in particular to a wind farm dynamic inertia partitioning evaluation method based on improved clustering. Background Art

[0002] As the construction of new power systems based on renewable energy continues to advance, the penetration rate of wind power in regional power grids across the country will increase rapidly. However, large-scale wind farms cannot provide rotational inertia support, and even after the introduction of virtual inertia control, there is still a lack of partitioned and accurate virtual inertia assessment methods. Currently, many researchers only consider the impact of wind speed on the available inertia of wind farms, without considering the influence of different wind directions and wake effects at actual wind farms, which will also affect the available inertia of wind farms. Among clustering methods, traditional K-means clustering requires the number of clusters and the initial cluster centers to be determined in advance. The pre-setting of the K value will directly affect the accuracy of the clustering results. Summary of the Invention

[0003] The technical problem to be solved by the present invention is to address the deficiencies of the above-mentioned existing technologies and provide a wind farm dynamic inertia partition evaluation method based on improved clustering. By combining the improved Canopy-Kmeans clustering algorithm and the Copula function, the accuracy and adaptability of the dynamic equivalent modeling of the wind farm are improved, and the inertia partition calculation of a large-scale wind farm is optimized.

[0004] In order to solve the above technical problems, the technical solution adopted by the present invention is:

[0005] A wind farm dynamic inertia partition assessment method based on improved clustering includes the following steps:

[0006] S1: Data acquisition and preprocessing: Obtain wind turbine operating data and environmental data from the wind farm's control system and monitoring system, including 13 initial state variables and environmental data, and preprocess the acquired data;

[0007] S2: dimensionality reduction by principal component analysis;

[0008] Construct a data matrix, converting the preprocessed 13 initial state variable data into a data matrix X(l,n), where n is the number of samples and l is the number of features, i.e., 13. Calculate the covariance matrix to measure the linear relationship between the initial state variables and the distribution characteristics of the data. Solve for eigenvalues ​​and eigenvectors, and select the y eigenvectors corresponding to the first h components as the dominant variables, with the sum of the cumulative variance contributions of the first h components being greater than 85%.

[0009] S3: Determine the initial cluster center and threshold using the improved Canopy-Kmeans algorithm, select the point closest to the mean point as the cluster center, and establish the Canopy dataset;

[0010] S4: Based on the clustering results, the available inertia is evaluated. By analyzing the spatiotemporal distribution characteristics of wind speed, a wind speed probability model of the wind farm is constructed using a hybrid function. The available inertia of each wind farm group is superimposed to obtain the available inertia distribution of the wind farm.

[0011] Furthermore, the specific method of step S1 is:

[0012] S1.1: Acquire data from the wind farm's control and monitoring systems, including wind turbine operating and environmental data. Specifically, this includes 13 initial state variables: pitch angle, rotational speed, electromagnetic torque, mechanical torque, wind turbine torque, stator current dq-axis components, rotor current dq-axis components, rotor voltage dq-axis components, and potential dq-axis components. Environmental data includes wind speed, wind direction, temperature, and air pressure. Perform a preliminary analysis of the collected data to identify any abnormalities.

[0013] S1.2: For obviously erroneous data, delete or replace them directly; when replacing, interpolate the data from the previous and next moments, or use statistical methods to estimate;

[0014] S1.3: For missing data, interpolate the average value using normal data from the previous and next moments, or use time series analysis methods for prediction; if there are many missing data, use data mining methods to fill them;

[0015] S1.4: For each variable, calculate its maximum and minimum values ​​in the collected data for normalization calculation; use the following formula to normalize the data:

[0016]

[0017] Among them, x is the original data, x min and x max are the minimum and maximum values ​​of the variable, respectively, x norm The data are normalized.

[0018] Furthermore, the specific method of step 2 is:

[0019] S2.1: Construct the data matrix X(n,l) from the 13 preprocessed initial state variables as follows:

[0020]

[0021] Each row represents a sample, i.e. a wind turbine, and each column represents an initial state variable. Each column corresponds to 13 features, i.e. pitch angle β, speed s, electromagnetic torque T M , mechanical torque T T , wind turbine torque T e , stator current d-axis and q-axis components I ds with I qs , rotor current d-axis and q-axis components I dr with I qr , rotor voltage d-axis and q-axis components U dr with U qr , the d-axis and q-axis components of the electric potential E d With E q ;

[0022] S2.2: Calculate the covariance matrix C of the data matrix X to measure the linear relationship between the initial state variables and the distribution characteristics of the data; the covariance matrix C is calculated as follows:

[0023]

[0024] S2.3: Perform eigenvalue decomposition on the covariance matrix C to obtain eigenvalues ​​λ1, λ2, ..., λ 13 , and their corresponding eigenvectors v1, v2, ..., v 13 ; Sort the eigenvalues ​​in descending order, and select the first h components as the principal components based on the principle that the cumulative variance contribution rate reaches the preset threshold, where h < 13;

[0025] The variance contribution rate βi of each principal component is:

[0026] Add the variance contribution rates of h principal components to get the cumulative variance contribution rate, that is:

[0027] Among them, λ i represents the eigenvalue corresponding to the i-th principal component;

[0028] When the cumulative variance contribution rate is greater than 85%, the dimension can be effectively reduced. The contribution of the remaining principal components after h is ignored, and the number of eigenvectors contained in the first h principal components is y;

[0029] S2.4: The selected y eigenvectors are formed into a principal component matrix P with a dimension of 13×y; the original data matrix X is linearly transformed using the principal component matrix P to obtain the reduced-dimensional data matrix Y, i.e., Y=X*P; each row of the reduced-dimensional data matrix Y represents the projection of a sample in the principal component space, and each column corresponds to a principal component, thus realizing the dimensionality reduction from 13 initial state variables to h principal components; based on the size of the eigenvalues, the y eigenvectors corresponding to the first h principal components whose cumulative variance contribution rate reaches the preset threshold are selected as the dominant variables, and the clustering of the doubly fed wind turbines is based on the first y principal components as clustering indicators.

[0030] Furthermore, the specific method of step 3 is:

[0031] S3.1: Use the Canopy algorithm to cluster the data set and generate multiple Canopy sets. Each Canopy set contains a certain number of measurement points that have certain similarities in inertia characteristics. The specific method is as follows:

[0032] S3.1.1: In selecting the initial threshold, first calculate the mean of the first h data sets, and then calculate the Euclidean distance D from all data points to the mean. i =|mean-c i |, where c i Represents the eigenvector of the i-th sample point in the data set, selects the data point closest to the mean point as the initial cluster center, calculates the distance between other data points and the initial cluster center, and counts the maximum distance as L1 and the minimum distance as L2; ​​at the same time, according to the distribution characteristics of the data set, determines the initial thresholds T1 and T2, T1 = L1-L2,

[0033] S3.1.2: Generation of Canopy Sets;

[0034] Determine the distance between the initial cluster center and all remaining data points, and classify the data points with a distance less than T1 into a Canopy set; remove the data points with a distance less than T2 from the dataset;

[0035] S3.1.3: Repeat step S3.1.2 until the first h data sets are completely cleared, that is, each data point has been assigned to a corresponding canopy set; each canopy set contains a certain number of data points that have a certain degree of similarity in inertia characteristics;

[0036] S3.2: Improved Canopy-Kmeans clustering algorithm;

[0037] Directly use the center point of the Canopy cluster and use the K-means algorithm to further cluster the data. Through iterative optimization, the inertia values ​​of the data points in the same cluster are made closer, and the inertia values ​​between different clusters are more different. Finally, the inertia clustering results of each partition are determined. The specific method is as follows:

[0038] S3.2.1: First, directly use the center point of the Canopy set in step S3.1 as the initial cluster center of k-means to perform preprocessing operations on the initial wind turbine dataset;

[0039] S3.2.2: Calculate the distance f from each data point to the initial cluster center, and classify each point into the corresponding category based on the principle of the closest distance.

[0040] S3.2.3: Calculate the average value of the state variables for all measurement points within each category and use this average value as the new cluster center for that category. Based on the new cluster center location, recalculate the distance between each data point and the new cluster center, and then reclassify each data point.

[0041] S3.2.4: Repeat the above steps S3.2.2 and S3.2.3 until the partitioning results and cluster centers remain stable and no longer change. At this time, the final wind turbine group partitioning results are obtained.

[0042] Furthermore, the specific method of step 4 is:

[0043] S4.1: Wind speed modeling;

[0044] By analyzing the spatiotemporal distribution characteristics of wind speed of each wind turbine in the wind farm, a probability distribution model of wind speed in the wind farm is constructed using a mixed Copula function.

[0045] Wake effect analysis: Considering the wake effect between wind turbines in a wind farm, the Jensen wake model is used to calculate the wind speed relationship between wind turbines:

[0046]

[0047] Among them, V WT1 、V WT2 is the wind speed at the terminals of wind turbines WT1 and WT2; d is the wind speed drop coefficient; CT is the thrust coefficient; R is the blade radius; D is the distance between the two wind turbines; tanα is the wake constant, and α is the cone apex coefficient. However, in actual operation, due to the geographical location of the wind turbines and the change of wind direction, most wind turbines are located in areas where the wake effect is partially blocked. Therefore, the revised formula is:

[0048]

[0049] in, is the actual wind speed of the downstream wind turbine WT2 under partial wake shielding; is the projected shading area of ​​wind turbine WT1 at WT2 in the case of partial shading; S wt2 is the rotating area of ​​the blade of wind turbine WT2; h wt1 、h wt2 are the blade shaft altitudes of wind turbine WT1 at WT2; δ is the wind speed height coefficient; for the wake effect of the upstream wind turbine n', the wind speed V of the downstream wind turbine j is j Affected by the superposition of the upstream wind turbine wake effect, as shown in the following formula:

[0050]

[0051] Among them, V j0 is the wind speed of wind turbine j without wake effect;

[0052] Time delay effect analysis: Considering the time delay of wind speed transmission from upstream wind turbines to downstream wind turbines, the wind speed relationship between wind turbines is calculated using the time delay model. The wind speed v(t) at wind turbine j under the time delay effect is expressed as:

[0053]

[0054] Among them, A k is the shielding area of ​​wind turbine k on wind turbine j, τ is the time required for wind to be transmitted from upstream to downstream, v(t) is the actual wind speed of wind turbine j at time t, and v0(t-τ) is the original wind speed of the upstream wind turbine at time t-τ;

[0055] S4.2: Optimal Copula function selection;

[0056] The Euclidean distance between the empirical Copula function and the theoretical Copula function is used to select the optimal Copula function to describe the wind speed correlation between each partition.

[0057] S4.2.1: Draw a probability histogram of the joint distribution of wind speed sequences between the central unit of each cluster and the wind tower, analyze its tail characteristics and symmetry characteristics, select the function with the smallest Euclidean distance to the empirical Copula function as the optimal Copula function, and obtain a probability histogram of the joint distribution of wind speed sequences between the central unit of each cluster and the wind tower;

[0058] S4.2.2: Construct an empirical Copula function. Suppose the wind speed data at the central cluster and the wind tower constitute a two-dimensional population (Z w1 ,Z w2 ) is a sample of (Z w1,i ,Z w2,i ), where i = 1, 2, ..., q; Z w1 and Zw2 The empirical distribution functions are F q (Z w1 ) and H q (Z w2 ), then the empirical Copula function is defined as:

[0059]

[0060] Where u and v are the cumulative probability values ​​of the marginal distribution; u = F(Z w1 ), indicating Z w1 The probability that the observed value is less than or equal to a certain threshold; v = H (Z w2 ), indicating Z w2 The probability that the observed value is less than or equal to a certain threshold; u,v∈[0,1], is the indicator function, when F q (Z w1,i )≤u, Otherwise, it is 0; where F q (Z w1,i ) is the variable Z in the i-th sample w1 The empirical marginal distribution function value is the upstream wind speed, H q (Z w2,j ) is the variable Z in the i-th sample w2 The empirical marginal distribution function value of is the downstream wind speed; when H q (Z w2,j )≤v, Otherwise, 0;

[0061] The maximum likelihood estimation method is used to calculate the relevant unknown parameters of each theoretical Copula function; the square of the Euclidean distance between the theoretical Copula function and the empirical Copula function is shown as follows:

[0062]

[0063] Where p is the theoretical Copula function type; Reflects the square of the Euclidean distance between the theoretical Copula function and the empirical Copula function, The Copula function corresponding to the minimum value is the optimal Copula;

[0064] After comparing and selecting the optimal Copula function using the above-mentioned Euclidean distance method, the binary joint probability distribution function between wind turbine i and wind tower m is:

[0065] F i (v i ,v m )=C i,m (F i(v i ),F m (v m ));

[0066] Among them, F i (v i ,v m ) represents the joint cumulative distribution function of wind speed at wind turbine i and wind tower m; C i,m (F i (v i ),F m (v m )) is the Copula function connecting wind turbine i and wind tower m, which is used to describe the dependency structure between the two; F i (v i ) and F m (v m ) are the marginal cumulative distribution functions of wind turbine i and wind tower m, respectively;

[0067] When estimating inertia, assuming that the wind speed data of the wind tower is known, the optimal Copula function is selected based on the binary joint probability distribution function in the above formula to describe the wind speed correlation between the wind tower and each group of turbines. The instantaneous wind speed probability density function of each group is obtained as follows:

[0068]

[0069] Among them, f(v i ∣v m ) is the given reference wind speed v m Under the condition of i The probability density function of i ∣v m ) is the given reference wind speed v m Under the condition, the wind speed of fan i is less than or equal to v i The cumulative probability of

[0070] S4.2.3: Based on the optimal Copula function, construct an inertia distribution model for the wind farm and evaluate the total available inertia of the wind farm and its distribution;

[0071] f(v i ∣v m )Integrate to obtain the probability distribution function of available inertia The binary search-numerical integration method is used to obtain the inertia values ​​corresponding to the wind turbine cluster at the probability values ​​of α / 2 and 1-α / 2 of the distribution function. The upper and lower limits of the available inertia interval of the wind turbine cluster are recorded respectively, where α represents the significance level and is used to define the confidence interval of the available inertia estimate. The available inertia of each cluster in the wind farm is superimposed to obtain the confidence interval of the available inertia of the entire wind farm. At the same time, the wake effect is used to analyze the wind speed of each unit in the cluster to obtain the available inertia distribution of each unit. The total available inertia and inertia distribution diagram are output.

[0072] The beneficial effect of adopting the above technical solution is that: the wind farm dynamic inertia partition evaluation method based on improved clustering provided by the present invention uses the Canopy clustering algorithm as a preprocessing step through the improved Canopy-Kmeans to optimize the blindness of the k-means algorithm in selecting the k value. The Canopy algorithm can obtain a rough clustering range through calculation without specifying the exact number of classifications, thereby reducing the uncertainty of the initial value selection in the k-means algorithm. The present invention calculates the mean point of the collected data set and selects the measurement point closest to the mean point as the initial cluster center, making the result more reasonable and stable. By combining the improved Canopy algorithm with the k-means algorithm, the present invention can more accurately identify the optimal classification number k value and cluster center, and then quickly and accurately realize the partitioning of large-area wind turbines. The improved Canopy-Kmeans clustering algorithm not only improves the partitioning efficiency, but also improves the reliability of the results. The present invention improves the accuracy of clustering by taking into account the wake effect and time lag effect and measuring the multi-dimensional wind turbine state variables. In order to avoid the dimensionality disaster of a large amount of data, PCA principal component analysis is adopted to reduce the data dimension, and the components with a variance contribution rate of more than 85% are selected as the dominant variables, which reduces the complexity of data collection. BRIEF DESCRIPTION OF THE DRAWINGS

[0073] Figure 1 A flow chart of a wind farm dynamic inertia partitioning assessment method based on improved clustering provided by an embodiment of the present invention;

[0074] Figure 2 A flowchart of the improved Canopy-Kmeans algorithm provided in an embodiment of the present invention;

[0075] Figure 3 A flow chart for evaluating the available inertia of a wind farm provided by an embodiment of the present invention;

[0076] Figure 4 A schematic diagram of the wake effect provided by an embodiment of the present invention;

[0077] Figure 5 This is a joint distribution probability histogram of wind speed sequences between the central unit of a cluster and a wind tower provided by an embodiment of the present invention. DETAILED DESCRIPTION

[0078] The following embodiments of the present invention are described in further detail with reference to the accompanying drawings and examples. The following examples are used to illustrate the present invention but are not intended to limit the scope of the present invention.

[0079] like Figure 1 As shown, the method of this embodiment is as follows.

[0080] S1: Data acquisition and preprocessing: Obtain the wind turbine operating data and environmental data from the wind farm control system and monitoring system, including 13 initial state variables and environmental data, and preprocess the acquired data. The specific method is:

[0081] S1.1: Acquire data from the wind farm's control and monitoring systems, including wind turbine operating and environmental data. Specifically, this includes 13 initial state variables: pitch angle, rotational speed, electromagnetic torque, mechanical torque, wind turbine torque, stator current dq-axis components, rotor current dq-axis components, rotor voltage dq-axis components, and potential dq-axis components. Environmental data includes wind speed, wind direction, temperature, and air pressure. Perform a preliminary analysis of the collected data to identify any abnormalities.

[0082] S1.2: For obviously erroneous data, delete or replace them directly; when replacing, interpolate the data from the previous and next moments, or use statistical methods to estimate;

[0083] S1.3: For missing data, interpolate the average value using normal data from the previous and next moments, or use time series analysis methods for prediction; if there are many missing data, use data mining methods to fill them;

[0084] S1.4: For each variable, calculate its maximum and minimum values ​​in the collected data for normalization calculation; use the following formula to normalize the data:

[0085]

[0086] Among them, x is the original data, x min and x max are the minimum and maximum values ​​of the variable, respectively, x norm The data are normalized.

[0087] S2: Dimensionality reduction by principal component analysis.

[0088] Construct a data matrix. The 13 pre-processed initial state variable data are constructed into a data matrix X(n,l), where n is the number of samples and l is the number of features, i.e., 13. Calculate the covariance matrix to measure the linear relationship between the initial state variables and the distribution characteristics of the data. Solve for eigenvalues ​​and eigenvectors, and select the y eigenvectors corresponding to the first h components as the dominant variables. The sum of the cumulative variance contributions of the first h components is greater than 85%. The specific method is:

[0089] S2.1: Construct the data matrix X(n,l) from the 13 preprocessed initial state variables as follows:

[0090]

[0091] Each row represents a sample, i.e. a wind turbine, and each column represents an initial state variable. Each column corresponds to 13 features, i.e. pitch angle β, speed s, electromagnetic torque T M , mechanical torque T T , wind turbine torque T e , stator current d-axis and q-axis components I ds with I qs , rotor current d-axis and q-axis components I dr with I qr , rotor voltage d-axis and q-axis components U dr with U qr , the d-axis and q-axis components of the electric potential E d With E q ;

[0092] S2.2: Calculate the covariance matrix C of the data matrix X to measure the linear relationship between the initial state variables and the distribution characteristics of the data; the covariance matrix C is calculated as follows:

[0093]

[0094] S2.3: Perform eigenvalue decomposition on the covariance matrix C to obtain eigenvalues ​​λ1, λ2, ..., λ 13 , and their corresponding eigenvectors v1, v2, ..., v 13 ; Sort the eigenvalues ​​in descending order, and select the first h components as the principal components based on the principle that the cumulative variance contribution rate reaches the preset threshold, where h < 13;

[0095] The variance contribution rate βi of each principal component is:

[0096] Add the variance contribution rates of h principal components to get the cumulative variance contribution rate, that is:

[0097] Among them, λi represents the eigenvalue corresponding to the i-th principal component;

[0098] When the cumulative variance contribution rate is greater than 85%, the dimension can be effectively reduced. The contribution of the remaining principal components after h is ignored, and the number of eigenvectors contained in the first h principal components is y;

[0099] S2.4: The selected y eigenvectors are formed into a principal component matrix P with a dimension of 13×y; the original data matrix X is linearly transformed using the principal component matrix P to obtain the reduced-dimensional data matrix Y, i.e., Y=X*P; each row of the reduced-dimensional data matrix Y represents the projection of a sample in the principal component space, and each column corresponds to a principal component, thereby achieving a dimensionality reduction process from 13 initial state variables to h principal components; based on the size of the eigenvalues, the y eigenvectors corresponding to the first h principal components whose cumulative variance contribution rate reaches a preset threshold are selected as the dominant variables, and the clustering of the doubly-fed wind turbine is based on the first y principal components as clustering indicators. In this embodiment, based on the size of the eigenvalues, the first four principal components whose cumulative contribution rate reaches a preset threshold are selected as the dominant variables, so the clustering of the doubly-fed wind turbine can be based on the first four principal components as clustering indicators.

[0100] S3: Determine the initial cluster center and threshold through the improved Canopy-Kmeans algorithm, select the point closest to the mean point as the cluster center, and establish the Canopy data set.

[0101] In the traditional Canopy algorithm, T1 and T2 are manually selected, which requires experience and trial and error and is inefficient. In view of the shortcomings of the traditional algorithm, this embodiment makes the following improvements in the selection of the initial threshold, such as Figure 2 As shown, the specific method is:

[0102] S3.1: Use the Canopy algorithm to cluster the data set and generate multiple Canopy sets. Each Canopy set contains a certain number of measurement points that have certain similarities in inertia characteristics. The specific method is:

[0103] S3.1.1: In selecting the initial threshold, first calculate the mean of the first h data sets, and then calculate the Euclidean distance D from all data points to the mean. i =|mean-c i |, where c i Represents the eigenvector of the i-th sample point in the data set, selects the data point closest to the mean point as the initial cluster center, calculates the distance between other data points and the initial cluster center, and counts the maximum distance as L1 and the minimum distance as L2; ​​at the same time, according to the distribution characteristics of the data set, determines the initial thresholds T1 and T2, T1 = L1-L2,

[0104] S3.1.2: Generation of Canopy Sets;

[0105] Determine the distance between the initial cluster center and all remaining data points, and classify the data points with a distance less than T1 into a Canopy set; remove the data points with a distance less than T2 from the dataset;

[0106] S3.1.3: Repeat step S3.2 until the first h data sets are completely cleared, that is, each data point has been classified into a corresponding Canopy set; each Canopy set contains a certain number of data points, and these data points have certain similarities in inertia characteristics.

[0107] S3.2: Improved Canopy-Kmeans clustering algorithm.

[0108] Initial cluster centers are selected from the Canopy set, and the K-means algorithm is used to further cluster the data. Through iterative optimization, the inertia values ​​of data points within the same cluster are made closer, and the inertia values ​​between different clusters are more different. Finally, the inertia clustering results of each partition are determined. The specific method is:

[0109] S3.2.1: First, directly use the center point of the Canopy set in step S3.1 as the initial cluster center of k-means to perform preprocessing operations on the initial wind turbine dataset;

[0110] S3.2.2: Calculate the distance f from each data point to the initial cluster center, and classify each point into the corresponding category based on the principle of the closest distance.

[0111] S3.2.3: Calculate the average value of the state variables for all measurement points within each category and use this average value as the new cluster center for that category. Based on the new cluster center location, recalculate the distance between each data point and the new cluster center, and then reclassify each data point.

[0112] S3.2.4: Repeat the above steps S3.2.2 and S3.2.3 until the partitioning results and cluster centers remain stable and no longer change. At this time, the final wind turbine group partitioning results are obtained.

[0113] S4: Based on the clustering results, the available inertia is evaluated. By analyzing the spatiotemporal distribution characteristics of wind speed, a wind speed probability model of the wind farm is constructed using a mixed function. The available inertia of each wind farm group is superimposed to obtain the available inertia distribution of the wind farm. The improved Canopy-Kmeans clustering algorithm is introduced, and the Copula function is used to evaluate the available inertia of the wind farm under different wind conditions. Figure 3 As shown, the specific method is:

[0114] S4.1: Wind speed modeling.

[0115] By analyzing the spatiotemporal distribution characteristics of wind speed of each wind turbine in the wind farm, a probability distribution model of wind speed in the wind farm is constructed using a mixed Copula function.

[0116] Wake effect analysis: Consider the wake effect between wind turbines in the wind farm. The wake effect principle is as follows: Figure 4 As shown, the Jensen wake model is used to calculate the wind speed relationship between wind turbines:

[0117]

[0118] Among them, V WT1 、V WT2 is the wind speed at the terminals of wind turbines WT1 and WT2; d is the wind speed drop coefficient; CT is the thrust coefficient; R is the blade radius; D is the distance between the two wind turbines; tanα is the wake constant, and α is the cone apex coefficient. However, in actual operation, due to the geographical location of the wind turbines and the change of wind direction, most wind turbines are located in areas where the wake effect is partially blocked. Therefore, the revised formula is:

[0119]

[0120] in, is the actual wind speed of the downstream wind turbine WT2 under partial wake shielding; is the projected shading area of ​​wind turbine WT1 at WT2 in the case of partial shading; S wt2 is the rotating area of ​​the blade of wind turbine WT2; h wt1 、h wt2 are the blade shaft altitudes of wind turbine WT1 at WT2; δ is the wind speed height coefficient; for the wake effect of the upstream wind turbine n', the wind speed V of the downstream wind turbine j is j Affected by the superposition of the upstream wind turbine wake effect, as shown in the following formula:

[0121]

[0122] Among them, V j0 is the wind speed of wind turbine j without wake effect;

[0123] Time delay effect analysis: Considering the time delay of wind speed transmission from upstream wind turbines to downstream wind turbines, the wind speed relationship between wind turbines is calculated using the time delay model. The wind speed v(t) at wind turbine j under the time delay effect is expressed as:

[0124]

[0125] Among them, A kis the shielding area of ​​wind turbine k on wind turbine j, τ is the time required for wind to be transmitted from upstream to downstream, v(t) is the actual wind speed of wind turbine j at time t, and v0(t-τ) is the original wind speed of the upstream wind turbine at time t-τ.

[0126] For the spatiotemporal modeling of wind speed, it is necessary to construct a mixed copula function and select the optimal copula function. Copula function is a method for describing the dependence of multiple variables. Sklar's theorem states that there exists a copula probability distribution function C(·) such that F(x1,x2,…,x n )=C(F1(x1),F2(x2),…,F n (x n ))f(x1,x2,…,x n ) holds true, and the correlation structure of the mixed Copula function is more flexible than that of a single Copula function, and can reflect the correlation pattern between different random variables. By weightedly combining multiple functions for any mixed Copula, the spatiotemporal correlation of wind speed can be more comprehensively reflected.

[0127] S4.2: Optimal Copula function selection;

[0128] The Euclidean distance between the empirical Copula function and the theoretical Copula function is used to select the optimal Copula function to describe the wind speed correlation between each partition. The specific method is as follows:

[0129] S4.2.1: Draw a joint distribution probability histogram of the wind speed sequence between the central unit of each cluster and the wind tower, analyze its tail characteristics and symmetry characteristics, select the function with the smallest Euclidean distance to the empirical Copula function as the optimal Copula function, and obtain the joint distribution probability histogram of the wind speed sequence between the central unit of each cluster and the wind tower; taking cluster 1 (assuming an average wind speed of 0.5) and wind tower 1 (assuming an average wind speed of 0.5) as an example, the joint distribution probability histogram of the wind speed sequence between the central unit of each cluster and the wind tower is obtained as follows: Figure 5 As shown;

[0130] S4.2.2: Construct an empirical Copula function, assuming that the wind speed data at the central cluster and the wind tower constitute a two-dimensional population (Z w1 ,Z w2 ) is a sample of (Z w1,i ,Z w2,i ), where i = 1, 2, ..., q; Z w1 and Z w2 The empirical distribution functions are F q (Z w1 ) and H q (Z w2), then the empirical Copula function is defined as:

[0131]

[0132] Where u and v are the cumulative probability values ​​of the marginal distribution; u = F(Z w1 ), indicating Z w1 The probability that the observed value is less than or equal to a certain threshold; v = H (Z w2 ), indicating Z w2 The probability that the observed value is less than or equal to a certain threshold; u,v∈[0,1], is the indicator function, when F q (z w1,i )≤u, Otherwise, it is 0; where F q (z w1,i ) is the variable z in the i-th sample w1 The empirical marginal distribution function value is the upstream wind speed, H q (z w2,j ) is the variable z in the i-th sample w2 The empirical marginal distribution function value of is the downstream wind speed; when H q (Z w2,j )≤v, Otherwise, 0;

[0133] The maximum likelihood estimation method is used to calculate the relevant unknown parameters of each theoretical Copula function; the square of the Euclidean distance between the theoretical Copula function and the empirical Copula function is shown as follows:

[0134]

[0135] Where p is the theoretical Copula function type; Reflects the square of the Euclidean distance between the theoretical Copula function and the empirical Copula function, The Copula function corresponding to the minimum value is the optimal Copula;

[0136] In order to more accurately select the Copula function, the above-mentioned Euclidean distance method is used for comparison. After selecting the optimal Copula function, the binary joint probability distribution function between wind turbine i and wind tower m is:

[0137] F i (v i ,v m )=C i,m (F i (v i ),F m (v m ));

[0138] Among them, F i (v i ,v m ) represents the joint cumulative distribution function of wind speed of wind turbine i and wind tower m; C i,m (F i (v i ),F m (v m )) is the Copula function connecting wind turbine i and wind tower m, which is used to describe the dependency structure between the two; F i (v i ) and F m (v m ) are the marginal cumulative distribution functions of wind turbine i and wind tower m, respectively;

[0139] When estimating inertia, assuming that the wind speed data of the wind tower is known, the optimal Copula function is selected based on the binary joint probability distribution function in the above formula to describe the wind speed correlation between the wind tower and each group of turbines. The instantaneous wind speed probability density function of each group is obtained as follows:

[0140]

[0141] Among them, f(v i ∣v m ) is the given reference wind speed v m Under the condition of i The probability density function of i ∣v m ) is the given reference wind speed v m Under the condition, the wind speed of fan i is less than or equal to v i The cumulative probability of

[0142] S4.2.3: Based on the optimal Copula function, construct an inertia distribution model for the wind farm and evaluate the total available inertia of the wind farm and its distribution;

[0143] f(v i ∣v m )Integrate to obtain the probability distribution function of available inertia The binary search-numerical integration method is used to obtain the inertia values ​​corresponding to the wind turbine cluster at the probability values ​​of α / 2 and 1-α / 2 of the distribution function. The upper and lower limits of the available inertia interval of the wind turbine cluster are recorded respectively, where α represents the significance level and is used to define the confidence interval of the available inertia estimate. The available inertia of each cluster in the wind farm is superimposed to obtain the confidence interval of the available inertia of the entire wind farm. At the same time, the wake effect is used to analyze the wind speed of each unit in the cluster to obtain the available inertia distribution of each unit. The total available inertia and inertia distribution diagram are output.

[0144] Finally, it should be noted that the above embodiments are only used to illustrate the technical solutions of the present invention, rather than to limit it. Although the present invention has been described in detail with reference to the above embodiments, those skilled in the art should understand that they can still modify the technical solutions described in the above embodiments, or make equivalent replacements for some or all of the technical features therein. However, these modifications or replacements do not cause the essence of the corresponding technical solutions to deviate from the scope defined by the claims of the present invention.

Claims

1. A wind farm dynamic inertia zoning assessment method based on improved clustering, characterized by: The method comprises the following steps: S1: Data acquisition and preprocessing: Obtain wind turbine operating data and environmental data from the wind farm's control system and monitoring system, including 13 initial state variables and environmental data, and preprocess the acquired data; S2: dimensionality reduction by principal component analysis; Construct a data matrix, converting the preprocessed 13 initial state variable data into a data matrix X(l,n), where n is the number of samples and l is the number of features, i.e., 13. Calculate the covariance matrix to measure the linear relationship between the initial state variables and the distribution characteristics of the data. Solve for eigenvalues ​​and eigenvectors, and select the y eigenvectors corresponding to the first h components as the dominant variables, with the sum of the cumulative variance contributions of the first h components being greater than 85%. S3: Determine the initial cluster center and threshold using the improved Canopy-Kmeans algorithm, select the point closest to the mean point as the cluster center, and establish the Canopy dataset; S4: Based on the clustering results, the available inertia is evaluated. By analyzing the spatiotemporal distribution characteristics of wind speed, a wind speed probability model of the wind farm is constructed using a hybrid function. The available inertia of each wind farm group is superimposed to obtain the available inertia distribution of the wind farm.

2. The wind farm dynamic inertia zoning assessment method based on improved clustering according to claim 1, characterized in that: The specific method of step S1 is: S1.1: Acquire data from the wind farm's control and monitoring systems, including wind turbine operating and environmental data. Specifically, this includes 13 initial state variables: pitch angle, rotational speed, electromagnetic torque, mechanical torque, wind turbine torque, stator current dq-axis components, rotor current dq-axis components, rotor voltage dq-axis components, and potential dq-axis components. Environmental data includes wind speed, wind direction, temperature, and air pressure. Perform a preliminary analysis of the collected data to identify any abnormalities. S1.2: For obviously erroneous data, delete or replace them directly; when replacing, interpolate the data from the previous and next moments, or use statistical methods to estimate; S1.3: For missing data, use the normal data of the previous and next moments to interpolate the average value, or use time series analysis methods to make predictions; If there are many missing data, data mining methods are used to fill them; S1.4: For each variable, calculate its maximum and minimum values ​​in the collected data for normalization calculation; The data were normalized using the following formula: Among them, x is the original data, x min and x max are the minimum and maximum values ​​of the variable, respectively, x norm The data are normalized.

3. The wind farm dynamic inertia zoning assessment method based on improved clustering according to claim 2, characterized in that: The specific method of step 2 is: S2.1: Construct the data matrix X(n,l) from the 13 preprocessed initial state variables as follows: Each row represents a sample, i.e. a wind turbine, and each column represents an initial state variable. Each column corresponds to 13 features, i.e. pitch angle β, speed s, electromagnetic torque T M , mechanical torque T T , wind turbine torque T e , stator current d-axis and q-axis components I ds with I qs , rotor current d-axis and q-axis components I dr with I qr , rotor voltage d-axis and q-axis components U dr with U qr , the d-axis and q-axis components of the electric potential E d With E q ; S2.2: Calculate the covariance matrix C of the data matrix X to measure the linear relationship between the initial state variables and the distribution characteristics of the data; the covariance matrix C is calculated as follows: S2.3: Perform eigenvalue decomposition on the covariance matrix C to obtain eigenvalues ​​λ1, λ2, ..., λ 13 , and their corresponding eigenvectors v1, v2, …, v 13 ; Sort the eigenvalues ​​in descending order, and select the first h components as the principal components based on the principle that the cumulative variance contribution rate reaches the preset threshold, where h < 13; The variance contribution rate βi of each principal component is: Add the variance contribution rates of h principal components to get the cumulative variance contribution rate, that is: Among them, λ i represents the eigenvalue corresponding to the i-th principal component; When the cumulative variance contribution rate is greater than 85%, the dimension can be effectively reduced. The contribution of the remaining principal components after h is ignored, and the number of eigenvectors contained in the first h principal components is y; S2.4: The selected y eigenvectors are formed into a principal component matrix P with a dimension of 13×y; the original data matrix X is linearly transformed using the principal component matrix P to obtain the reduced-dimensional data matrix Y, i.e., Y=X*P; each row of the reduced-dimensional data matrix Y represents the projection of a sample in the principal component space, and each column corresponds to a principal component, thus realizing the dimensionality reduction from 13 initial state variables to h principal components; based on the size of the eigenvalues, the y eigenvectors corresponding to the first h principal components whose cumulative variance contribution rate reaches the preset threshold are selected as the dominant variables, and the clustering of the doubly fed wind turbines is based on the first y principal components as clustering indicators.

4. The wind farm dynamic inertia zoning assessment method based on improved clustering according to claim 3 is characterized by: The specific method of step 3 is: S3.1: Use the Canopy algorithm to cluster the data set and generate multiple Canopy sets. Each Canopy set contains a certain number of measurement points that have certain similarities in inertia characteristics. The specific method is as follows: S3.1.1: In the selection of the initial threshold, first calculate the mean of the first h data sets, and then calculate the Euclidean distance D from all data points to the mean. i =|mean-c i |, where c i Represents the eigenvector of the i-th sample point in the data set, selects the data point closest to the mean point as the initial cluster center, calculates the distance between other data points and the initial cluster center, and counts the maximum distance as L1 and the minimum distance as L2; ​​at the same time, according to the distribution characteristics of the data set, determines the initial thresholds T1 and T2, T1 = L1-L2, S3.1.2: Generation of Canopy Sets; Determine the distance between the initial cluster center and all remaining data points, and classify the data points with a distance less than T1 into a Canopy set; remove the data points with a distance less than T2 from the dataset; S3.1.3: Repeat step S3.1.2 until the first h data sets are completely cleared, that is, each data point has been assigned to a corresponding canopy set; each canopy set contains a certain number of data points that have a certain degree of similarity in inertia characteristics; S3.2: Improved Canopy-Kmeans clustering algorithm; Directly use the center point of the Canopy cluster and use the K-means algorithm to further cluster the data. Through iterative optimization, the inertia values ​​of the data points in the same cluster are made closer, and the inertia values ​​between different clusters are more different. Finally, the inertia clustering results of each partition are determined. The specific method is as follows: S3.2.1: First, directly use the center point of the Canopy set in step S3.1 as the initial cluster center of k-means to perform preprocessing operations on the initial wind turbine dataset; S3.2.2: Calculate the distance f from each data point to the initial cluster center, and classify each point into the corresponding category based on the principle of the closest distance. S3.2.3: Calculate the average value of the state variables for all measurement points within each category and use this average value as the new cluster center for that category. Based on the new cluster center location, recalculate the distance between each data point and the new cluster center, and then reclassify each data point. S3.2.4: Repeat the above steps S3.2.2 and S3.2.3 until the partitioning results and cluster centers remain stable and no longer change. At this time, the final wind turbine group partitioning results are obtained.

5. The method for evaluating wind farm dynamic inertia based on improved clustering according to claim 4, characterized in that: The specific method of step 4 is: S4.1: Wind speed modeling; By analyzing the spatiotemporal distribution characteristics of wind speed of each wind turbine in the wind farm, a probability distribution model of wind speed in the wind farm is constructed using a mixed Copula function. Wake effect analysis: Considering the wake effect between wind turbines in a wind farm, the Jensen wake model is used to calculate the wind speed relationship between wind turbines: Among them, V WT1 、V WT2 is the wind speed at the terminals of wind turbines WT1 and WT2; d is the wind speed drop coefficient; CT is the thrust coefficient; R is the blade radius; D is the distance between the two wind turbines; tanα is the wake constant, and α is the cone apex coefficient. However, in actual operation, due to the geographical location of the wind turbines and the change of wind direction, most wind turbines are located in areas where the wake effect is partially blocked. Therefore, the revised formula is: in, is the actual wind speed of the downstream wind turbine WT2 under partial wake shielding; is the projected shading area of ​​wind turbine WT1 at WT2 in the case of partial shading; S wt2 is the rotating area of ​​the blade of wind turbine WT2; h wt1 、h wt2 are the blade shaft altitudes of wind turbine WT1 at WT2; δ is the wind speed height coefficient; for the wake effect of the upstream wind turbine n', the wind speed V of the downstream wind turbine j is j Affected by the upstream wind turbine wake effect, the following formula is used: Among them, V j0 is the wind speed of wind turbine j without wake effect; Time delay effect analysis: Considering the time delay of wind speed transmission from upstream wind turbines to downstream wind turbines, the wind speed relationship between wind turbines is calculated using the time delay model. The wind speed v(t) at wind turbine j under the time delay effect is expressed as: Among them, A k is the shielding area of ​​wind turbine k on wind turbine j, τ is the time required for wind to be transmitted from upstream to downstream, v(t) is the actual wind speed of wind turbine j at time t, and v0(t-τ) is the original wind speed of the upstream wind turbine at time t-τ; S4.2: Optimal Copula function selection; The Euclidean distance between the empirical Copula function and the theoretical Copula function is used to select the optimal Copula function to describe the wind speed correlation between each partition. S4.2.1: Draw a probability histogram of the joint distribution of wind speed sequences between the central unit of each cluster and the wind tower, analyze its tail characteristics and symmetry characteristics, select the function with the smallest Euclidean distance to the empirical Copula function as the optimal Copula function, and obtain a probability histogram of the joint distribution of wind speed sequences between the central unit of each cluster and the wind tower; S4.2.2: Construct an empirical Copula function. Suppose the wind speed data at the central cluster and the wind tower constitute a two-dimensional population (Z w1 ,Z w2 ) is a sample of (Z w1,i ,Z w2,i ), where i = 1, 2, ..., q; Z w1 and Z w2 The empirical distribution functions are F q (Z w1 ) and H q (Z w2 ), then the empirical Copula function is defined as: Where u and v are the cumulative probability values ​​of the marginal distribution; u = F(Z w1 ), indicating Z w1 The probability that the observed value is less than or equal to a certain threshold; v = H (Z w2 ), indicating Z w2 The probability that the observed value is less than or equal to a certain threshold; u,v∈[0,1], is the indicator function, when F q (Z w1,i )≤u, Otherwise, it is 0; where F q (Z w1,i ) is the variable Z in the i-th sample w1 The empirical marginal distribution function value is the upstream wind speed, H q (Z w2,j ) is the variable Z in the i-th sample w2 The empirical marginal distribution function value of is the downstream wind speed; when H q (Z w2,j )≤v, Otherwise, 0; The maximum likelihood estimation method is used to calculate the relevant unknown parameters of each theoretical Copula function; the square of the Euclidean distance between the theoretical Copula function and the empirical Copula function is shown as follows: Where p is the theoretical Copula function type; Reflects the square of the Euclidean distance between the theoretical Copula function and the empirical Copula function, The Copula function corresponding to the minimum value is the optimal Copula; After comparing and selecting the optimal Copula function using the above-mentioned Euclidean distance method, the binary joint probability distribution function between wind turbine i and wind tower m is: F i (v i ,v m )=C i,m (F i (v i ),F m (v m )); Among them, F i (v i ,v m ) represents the joint cumulative distribution function of wind speed at wind turbine i and wind tower m; C i,m (F i (v i ),F m (v m )) is the Copula function connecting wind turbine i and wind tower m, which is used to describe the dependency structure between the two; F i (v i ) and F m (v m ) are the marginal cumulative distribution functions of wind turbine i and wind tower m, respectively; When estimating inertia, assuming that the wind speed data of the wind tower is known, the optimal Copula function is selected based on the binary joint probability distribution function in the above formula to describe the wind speed correlation between the wind tower and each group of turbines. The instantaneous wind speed probability density function of each group is obtained as follows: Among them, f(v i ∣v m ) is the given reference wind speed v m Under the condition of i The probability density function of i ∣v m ) is the given reference wind speed v m Under the condition, the wind speed of fan i is less than or equal to v i The cumulative probability of S4.2.3: Based on the optimal Copula function, construct an inertia distribution model for the wind farm and evaluate the total available inertia of the wind farm and its distribution; f(v i ∣v m )Integrate to obtain the probability distribution function of available inertia The binary search-numerical integration method is used to obtain the inertia values ​​corresponding to the wind turbine cluster at the probability values ​​of α / 2 and 1-α / 2 of the distribution function. The upper and lower limits of the available inertia interval of the wind turbine cluster are recorded respectively, where α represents the significance level and is used to define the confidence interval of the available inertia estimate. The available inertia of each cluster in the wind farm is superimposed to obtain the confidence interval of the available inertia of the entire wind farm. At the same time, the wake effect is used to analyze the wind speed of each unit in the cluster to obtain the available inertia distribution of each unit. The total available inertia and inertia distribution diagram are output.

Images