Method for determining the ultimate performance response of long flexible blades of large wind turbines

Abstract

The present invention discloses a method for determining the extreme performance response of long flexible blades of large wind turbines. First, a dynamic structural response model of long flexible blades is established based on geometric nonlinear beam theory, and a discretization solution is performed through a one-dimensional equivalent beam unit to accurately calculate the performance responses such as blade tip deformation, stress and strain. Next, the method comprehensively considers the joint uncertainty of the 10-minute average wind speed and turbulence intensity, performs a large number of 10-minute long flexible blade structural deformation simulations at discrete average wind speed and turbulence intensity sample points, extracts performance response extreme value samples, and fits their conditional cumulative distribution functions based on the maximum likelihood estimation method. According to the joint probability of the average wind speed and turbulence intensity discrete sample points, the performance response extreme value conditional cumulative distribution functions at each sample point are weighted summed, the exceedance probability of the performance response extreme value samples is calculated, and the 50-year limit performance response value is calculated using the statistical extrapolation method.

Description

Technical Field

[0001] The invention relates to a method for determining the limit performance response of a long flexible blade of a large wind turbine generator, and in particular to a method for determining the limit performance response of a long flexible blade of a large wind turbine generator taking into account the joint distribution of average wind speed and turbulence intensity. Background Art

[0002] As a key component of large wind turbines, accurate estimation of the 50-year extreme performance response of blades has a crucial impact on the long life cycle design results of the unit considering the extreme state. The current method for calculating the extreme performance response of large wind turbines is to select multiple discrete wind speed sample points, run the 10-minute random structure simulation model multiple times at each wind speed sample point, fit its performance response extreme value distribution, and then perform weighted summation based on the long-term probability of wind speed to obtain the 10-minute extreme performance response long-term distribution model, and achieve statistical extrapolation of the extreme performance response by calculating the 50-year recurrence probability. However, the existing method for estimating the extreme performance response of long flexible blades of large wind turbines has the following defects:

[0003] (1) With the trend of large-scale wind turbines, the slender elastic structure of the blades and the application of composite materials make the blades subject to deformation coupling in different directions, and their deformation nonlinear characteristics are significant, which leads to deviations in the calculation of the existing Euler-Bernoulli beam blade structure model, thus affecting the accuracy of the prediction of the ultimate performance response;

[0004] (2) The existing extreme performance response calculation method only considers the long-term distribution of wind speed, but ignores the joint distribution of other environmental factors such as turbulence intensity, which may lead to certain deviations in the prediction results. Summary of the invention

[0005] In order to solve the above-mentioned technical problems existing in the prior art, the purpose of the present invention is to propose a method for determining the ultimate performance response of long and flexible blades of large wind turbines, aiming at the problems of inaccurate blade structure simulation model and incomplete consideration of environmental factor uncertainty in the calculation process of the ultimate performance response of long and flexible blades of large wind turbines. The method improves the prediction accuracy of the ultimate performance response of blades by integrating the nonlinear structural analysis model of the blades and comprehensively considering the joint distribution of environmental uncertainties such as wind speed and turbulence intensity.

[0006] The technical solution for implementing the present invention is as follows:

[0007] A method for determining the limit performance response of a long flexible blade of a large wind turbine generator comprises the following steps:

[0008] S1: According to the design parameters of the wind turbine and the measured environmental data of the actual wind field, determine the unit operating condition X = [V 10 ,I 10] of discrete sample points X1,X2,…,X M , where V 10 Represents the average wind speed of the unit within 10 minutes, I 10 represents the turbulence intensity of the unit within 10 minutes, and M is the total number of discrete sample points;

[0009] S2: At each V 10 ,I 10 At discrete sample points, N 10-minute unit blade structure simulations are performed; the blade performance response time series data is extracted to obtain N sample points of performance response extreme values ​​L, and statistical distribution fitting is performed using the maximum likelihood estimation method to obtain the cumulative distribution function F of the performance response extreme value under this condition. L (l|X i );

[0010] S3: According to each V 10 ,I 10 The probability P of discrete sample points X (X i ), the cumulative distribution function of the extreme value of the performance response F L (l|X i ) is weighted summed as shown below:

[0011]

[0012] S4: For all performance response extreme value sample points l j , the exceedance probability is calculated according to the following formula

[0013] P F (L>l j )=1-F L (l j ) (2)

[0014] After taking the logarithm of all sample points, numerical fitting and extrapolation are performed to obtain the exceedance probability P F =3.8×10 -7 The corresponding limit performance response l 50 The calculation result is the ultimate performance response of the wind turbine blades with a return period of 50 years.

[0015] Furthermore, the 10-minute unit blade structure simulation in S2 adopts a nonlinear geometric deformation model. The model construction and solution steps are as follows:

[0016] S2.1: Discretize the blade into several one-dimensional beam units, and obtain the 6×6 mass matrix M and stiffness matrix S at the cross section of each beam unit by referring to the unit design documents or composite lamination theory;

[0017] S2.2: Based on the blade element theory, the distributed aerodynamic performance response of each beam unit of the blade is calculated. In addition, the gravity and centrifugal force of each beam unit of the blade are calculated. According to the blade pitch angle, the aerodynamic and structural torsion angle at the previous time step, and the deformation data of the previous time step, the coordinate system of each performance response at the beam unit is transformed to obtain the externally applied distributed force f ext and distributed moment m ext ;

[0018] S2.3: Based on the Hamiltonian principle, the constitutive equation of the blade one-dimensional beam model is constructed:

[0019]

[0020] Among them, h and g are the linear momentum and angular momentum of the beam element respectively, u is the linear displacement of the beam element, F and T represent the response force and response moment generated by the cross section of the beam element respectively, and r represents the initial position of the beam element. represents the time derivative operator, (·)′ represents the derivative of the coordinates of the center axis of the one-dimensional beam, Represents the transformation of a vector into its corresponding antisymmetric tensor.

[0021] Furthermore, let u = [u1,u2,u3] T ,but

[0022]

[0023] h, g, F, T are obtained through the continuity equation:

[0024]

[0025] Among them, ω is the angular velocity of the cross section, γ and κ are the strain and bending vector of the beam respectively. The Newton-Raphson method is used to solve equation (3) to obtain the dynamic response of the blade.

[0026] Furthermore, in S2, the conditional cumulative distribution function F of the blade performance response extreme value is L (l|X i )Statistical fitting, select the distribution type such as logarithmic distribution, Weibull distribution, or generalized extreme value distribution.

[0027] Furthermore, each V in S3 10 ,I 10 Discrete sample point probability P X (X i ) calculation includes the following steps:

[0028] S3.1: Based on the wind farm measured data or historical database, use the maximum likelihood estimation method to fit the 10-minute wind speed mean V 10 With standard deviation ∑10 The parameters C, k and a, b of the marginal probability density function are:

[0029]

[0030]

[0031] S3.2: Based on the fitting results in S3.1, calculate V 10 and∑ 10 The marginal cumulative distribution functions u and v of , and calculate V 10 and∑ 10 The Gumbel Copula density function c v∑ (v 10 ,σ 10 ; τ), as shown in the following formula:

[0032]

[0033] Where τ is the Kendall parameter, w = (-lnu) 1 / (1-τ) +(-lnv) 1 / (1-τ) ;

[0034] S3.3: Calculate the 10-minute average wind speed V of the group according to the following formula: 10 and turbulence intensity I 10 Joint probability density function:

[0035]

[0036] S3.4: According to the determined discrete working condition measurement point resolution ΔV 10 ,ΔI 10 The probability P corresponding to each working condition measurement point X (X i )for:

[0037] P X (X i ) = f 10 (v 10 ,i 10 ; C,k,a,b,τ)ΔV 10 ΔI 10 (3).

[0038] The beneficial effects of the present invention are: (1) the present invention integrates the nonlinear structural simulation model of the blade, and the 10-minute structural simulation of the blade is more in line with the actual operation of the wind turbine; (2) the joint distribution of environmental uncertainties such as wind speed and turbulence intensity is comprehensively considered to achieve accurate evaluation of the extreme performance response. BRIEF DESCRIPTION OF THE DRAWINGS

[0039] The present invention is further described below in conjunction with the accompanying drawings and embodiments;

[0040] Figure 1 Schematic diagram of the inventive method flow chart;

[0041] Figure 2 Schematic diagram of the long flexible blade and its coordinate system in the invention;

[0042] Figure 3 The invention shows a curve diagram of the tip deformation of the long flexible blade in a static wind field as the wind speed changes;

[0043] Figure 4 The 10-minute average wind speed and turbulence intensity joint probability density distribution diagram in the embodiment of the invention;

[0044] Figure 5 Yes V 10 =19m / s,I 10 =0.14 Next is a typical 10-minute long flexible blade tip deformation time series data graph;

[0045] Figure 6 Yes V 10 =19m / s,I 10 = 0.14, the statistical distribution result diagram of the extreme value samples of blade tip deformation and their corresponding fitting;

[0046] Figure 7 It is the result diagram of blade tip deformation exceedance probability predicted by this method; DETAILED DESCRIPTION

[0047] The technical solutions in the embodiments of the present application will be clearly and completely described below in conjunction with the drawings in the embodiments of the present application. Obviously, the described embodiments are only a part of the embodiments of the present application, rather than all of the embodiments. Based on the embodiments in the present application, all other embodiments obtained by ordinary technicians in the field without creative work are within the scope of protection of the present application. In order to make the above-mentioned purposes, features and advantages of the present application more obvious and easy to understand, the present application will be further described in detail below in conjunction with the drawings and specific implementation methods.

[0048] like Figure 1 As shown, the method for determining the limit performance response of a long flexible blade of a large wind turbine of the present invention comprises the following steps:

[0049] S1: Within the normal wind speed and turbulence intensity range of the fan, uniformly determine the unit operating condition X = [V 10 ,I 10 ] of discrete measurement points X1,X2,…,X M , where V 10 Represents the average wind speed of the unit within 10 minutes, I 10It represents the turbulence intensity of the unit within 10 minutes, and M is the total number of operating measurement points.

[0050] S2: At each V 10 ,I 10 At discrete sample points, perform N 10-minute unit blade structure simulations. Extract the blade performance response time series data to obtain N sample points with the maximum performance response value L. According to the distribution types such as lognormal distribution, Weibull distribution or generalized extreme value distribution, perform statistical distribution fitting by the maximum likelihood estimation method to obtain the cumulative distribution function F of the performance response under this condition. L (l|X i ). Among them, constructing the blade unit structure simulation model specifically includes the following steps:

[0051] S2.1: Discretize the blade into several one-dimensional beam units, and obtain the 6×6 mass matrix M and stiffness matrix S at the cross section of each beam unit by referring to the unit design documents or the Classical Laminate Theory (CLT) and other methods.

[0052] S2.2: Based on the blade element theory, the distributed aerodynamic performance response of each beam unit of the blade is calculated. In addition, the gravity and centrifugal force of each beam unit of the blade are calculated. According to the blade pitch angle, the aerodynamic and structural torsion angle at the previous time step, and the deformation data of the previous time step, the coordinate system of each performance response at the beam unit is transformed to obtain the externally applied distributed force f ext and distributed moment m ext .

[0053] S2.3: Based on the Hamiltonian principle, the constitutive equation of the blade one-dimensional beam model is constructed:

[0054]

[0055] Among them, h and g are the linear momentum and angular momentum of the beam element respectively, u is the linear displacement of the beam element, F and T represent the response force and response moment generated by the cross section of the beam element respectively, and r represents the initial position of the beam element. represents the time derivative operator, (·)′ represents the derivative of the coordinates of the center axis of the one-dimensional beam, Indicates converting a vector into its corresponding antisymmetric tensor. For example, if u = [u1,u2,u3] T ,but

[0056]

[0057] h, g, F, T can be obtained through the continuity equation:

[0058]

[0059] Among them, ω is the angular velocity of the cross section, γ and κ are the strain and bending vector of the beam respectively. The Newton-Raphson method is used to solve equation (3) to obtain the dynamic response of the blade.

[0060] S3: According to each V 10 ,I 10 The probability P of discrete sample points X (X i ), the cumulative distribution function of the extreme value of the performance response F L (l|X i ) is weighted summed as shown below:

[0061]

[0062] Each V 10 ,I 10 The probability calculation of discrete sample points follows the following steps:

[0063] S3.1: Based on the wind farm measured data or historical database, use the maximum likelihood estimation method to fit the 10-minute wind speed mean V 10 With standard deviation ∑ 10 The parameters C, k and a, b of the marginal probability density function are:

[0064]

[0065]

[0066] S3.2: Based on the fitting results in S3.1, calculate V 10 and∑ 10 The marginal cumulative distribution functions u and v of , and calculate V 10 and∑ 10 The Gumbel Copula density function c v∑ (v 10 ,σ 10 ; τ), as shown in the following formula:

[0067]

[0068] Where τ is the Kendall parameter, w = (-lnu) 1 / (1-τ) +(-lnv) 1 / (1-τ) .

[0069] S3.3: Calculate the 10-minute average wind speed V of the group according to the following formula: 10 and turbulence intensity I 10 Joint probability density function:

[0070]

[0071] S3.4: According to the determined discrete working condition measurement point resolution ΔV 10 ,ΔI 10 The probability P corresponding to each working condition measurement point X (X i )for:

[0072] P X (X i ) = f 10 (v 10 ,i 10 ; C,k,a,b,τ)ΔV 10 ΔI 10 (9)

[0073] S4: For all performance response extreme value sample points l j , the exceedance probability is calculated according to the following formula

[0074] P F (L>l j )=1-F L (l j ) (2)

[0075] After taking the logarithm of these sample points and performing numerical fitting and extrapolation, the exceedance probability is obtained as P F =3.8×10 -7 The corresponding limit performance response l 50 As the prediction result of the ultimate performance response of wind turbine blades with a return period of 50 years.

[0076] Example

[0077] The process of the present invention is applied to the 5MW reference model proposed by the American Renewable Energy Laboratory, thereby proving the effectiveness of the method of the present invention. The 5MW model has a blade length of 61.5m, a cut-in wind speed of 3m / s, a cut-out wind speed of 25m / s, a rated wind speed of 11.4m / s, a rated speed of 12.1rpm, a hub height of 90m, and an integrated variable speed and pitch control system. Under normal power generation conditions of the unit, the method of the present invention predicts the tip swing deformation value of the unit once in 50 years.

[0078] First, in order to demonstrate the accuracy of the blade one-dimensional equivalent beam model established by S2 in the simulation of nonlinear large deformation, within the wind speed range of 11 to 25 m / s, wind speed sample points were determined every 2 m / s to establish a static wind field, and the unit yaw angle was fixed to 0°, the pitch angle was 2°, and the wind turbine speed was 12.1 rpm. The blade tip swing deformation was simulated using the one-dimensional equivalent beam model of the blade in this method and the Euler-Bernoulli beam, and the results were as follows: Figure 3As shown. The results show that the simulation model of the blade structure of the present invention is significantly smaller than that of the Euler-Bernoulli beam, and the deformation of the blade tip is reduced by 14.3% at a wind speed of 25m / s. This is because the establishment of the Euler-Bernoulli beam theory needs to obey the small deformation assumption, and in the application scenario of long and flexible blades of large wind turbines, the blades often undergo large nonlinear deformations, resulting in inaccurate prediction results of the Euler-Bernoulli beam theory, which further illustrates that the method of the present invention can effectively simulate the nonlinear deformation of the blade structure and lay the foundation for the accurate simulation of the extreme performance response.

[0079] In this embodiment, the values ​​of the parameters of the probability density function of the joint distribution of wind speed and turbulence intensity are as follows: C = 6.5856, k = 2.5178, a = 3.1570, b = 0.4123, τ = 0.6826, and the obtained V 10 ,I 10 The joint probability density function is Figure 4 shown.

[0080] In this embodiment, within the wind speed range of 3m / s to 25m / s, a wind speed measurement point is selected every 2m / s; the turbulence intensity is between 0.02 and 0.40, and the interval is 0.02 as the turbulence intensity measurement point, that is, ΔV 10 =2m / s, ΔI = 0.02. Thus, a total of 12×20 = 240 working condition measurement points were determined. At each working condition measurement point, 100 10-minute blade structure simulations were performed to extract the tip deformation extreme value samples. The tip deformation extreme value cumulative distribution function F at each working condition measurement point was fitted by log-normal distribution combined with maximum likelihood estimation method. L (l|V 10 ,I 10 ). Figure 5 It shows a typical time series data of the tip deformation of a long flexible blade. Figure 6 The statistical distribution results of the extreme value samples of blade tip deformation at a measuring point in an operating condition and their corresponding fitting are shown.

[0081] Figure 7 The sample points of blade tip deformation exceedance probability distribution are shown, as well as the exceedance probability results fitted and extrapolated based on these sample points, and P is obtained. F =3.8×10 -7 The corresponding blade tip deformation value once in 50 years is 9.944m.

[0082] In order to verify the effectiveness of the proposed method, the results obtained in this embodiment are compared with the existing statistical extrapolation method process that only considers wind speed uncertainty (the method process refers to Yang Q, Li Y, Li T, et al. Statistical extrapolation methods and empirical formulae for estimating extreme loads on operating wind turbine towers [J]. Engineering Structures, 2022, 267: 114667.), and the experiment is repeated 100 times to compare the predicted mean and variance of the extreme performance response of the two. At the same time, the results of the simulation using the Monte-Carlo method of Sandia National Laboratory in the United States (refer to Barone MF, Paquette JA, Resor BR, et al. Decades of Wind Turbine Load Simulation [R]. Sandia National Lab. (SNL-NM), Albuquerque, NM (United States), 2011.) are used as reference values ​​for comparison, and the results are shown in Table 1. The results show that the method for predicting the extreme performance response of large long flexible blades proposed in the present invention takes into account the uncertainty distribution of turbulence intensity. At wind speeds where high-performance response extremes may occur, performance response extreme sample points corresponding to greater turbulence intensity can be obtained. Therefore, the prediction result is greater than the extreme performance response statistical extrapolation method that only considers wind speed uncertainty and is closer to the reference value obtained by the Monte-Carlo method, which proves the effectiveness of the present invention.

[0083]

[0084] Table 1 Comparison of calculation results of real-time cases (100 test results).

[0085] Those skilled in the art can understand that the above are only preferred examples of the invention and are not intended to limit the invention. Although the invention is described in detail with reference to the above examples, those skilled in the art can still modify the technical solutions recorded in the above examples or replace some of the technical features therein with equivalents. Any modification, equivalent replacement, etc. made within the spirit and principle of the invention shall be included in the protection scope of the invention.

Claims

1. A method for determining the limit performance response of a long flexible blade of a large wind turbine, characterized in that The following steps are involved: S1: According to the design parameters of the wind turbine and the measured environmental data of the actual wind field, determine the unit operating condition X = [V 10 ,I 10 ] of discrete sample points X1,X2,…,X M , where V 10 Represents the average wind speed of the unit within 10 minutes, I 10 represents the turbulence intensity of the unit within 10 minutes, and M is the total number of discrete sample points; S2: At each V 10 ,I 10 At discrete sample points, N 10-minute unit blade structure simulations are performed; the blade performance response time series data is extracted to obtain N sample points of performance response extreme values ​​L, and statistical distribution fitting is performed using the maximum likelihood estimation method to obtain the cumulative distribution function F of the performance response extreme value under this condition. L (l|X i ); S3: According to each V 10 ,I 10 The probability P of discrete sample points X (X i ), the cumulative distribution function of the extreme value of the performance response F L (l|X i ) is weighted summed as shown below: S4: For all performance response extreme value sample points l j , the exceedance probability is calculated according to the following formula P F (L>l j )=1-F L (l j ) (2) After taking the logarithm of all sample points, numerical fitting and extrapolation are performed to obtain the exceedance probability P F =3.8×10 -7 The corresponding limit performance response l 50 As the calculation result of the ultimate performance response of the wind turbine blades under the 50-year return period; The 10-minute unit blade structure simulation in S2 uses a nonlinear geometric deformation model. The model construction and solution steps are as follows: S2.1: Discretize the blade into several one-dimensional beam units, and obtain the 6×6 mass matrix M and stiffness matrix S at the cross section of each beam unit by referring to the unit design documents or composite lamination theory; S2.2: Based on the blade element theory, the distributed aerodynamic performance response of each beam unit of the blade is calculated. In addition, the gravity and centrifugal force of each beam unit of the blade are calculated. According to the blade pitch angle, the aerodynamic and structural torsion angle at the previous time step, and the deformation data of the previous time step, the coordinate system of each performance response at the beam unit is transformed to obtain the externally applied distributed force f ext and distributed moment m ext ; S2.3: Based on the Hamiltonian principle, the constitutive equation of the blade one-dimensional beam model is constructed: Among them, h and g are the linear momentum and angular momentum of the beam element, u is the linear displacement of the beam element, F and T represent the response force and response moment generated by the beam element section, and r represents the initial position of the beam element. represents the time derivative operator, (·)′ represents the derivative of the coordinates of the center axis of the one-dimensional beam, Indicates converting a vector into its corresponding antisymmetric tensor; Let u=[u1,u2,u3] T ,but h, g, F, T are obtained through the continuity equation: Among them, ω is the angular velocity of the cross section, γ and κ are the strain and bending vector of the beam, respectively. The Newton-Raphson method is used to solve equation (3) to obtain the dynamic response of the blade; Each V in S3 10 ,I 10 Discrete sample point probability P X (X i ) calculation includes the following steps: S3.1: Based on the wind farm measured data or historical database, use the maximum likelihood estimation method to fit the 10-minute wind speed mean V 10 With standard deviation ∑ 10 The parameters C, k and a, b of the marginal probability density function are: S3.2: Based on the fitting results in S3.1, calculate V 10 and∑ 10 The marginal cumulative distribution functions u and v of , and calculate V 10 and∑ 10 The Gumbel Copula density function c v∑ (v 10 ,σ 10 ; τ), as shown in the following formula: Where τ is the Kendall parameter, w = (-lnu) 1 / (1-τ) +(-lnv) 1 / (1-τ) ; S3.3: Calculate the 10-minute average wind speed V of the group according to the following formula: 10 and turbulence intensity I 10 Joint probability density function: S3.4: According to the determined discrete working condition measurement point resolution ΔV 10 ,ΔI 10 The probability P corresponding to each working condition measurement point X (X i )for: P X (X i )=f 10 (v 10 ,i 10 ;C,k,a,b,τ)ΔV 10 ΔI 10 (9)。 2. The method for determining the limit performance response of a long flexible blade of a large wind turbine according to claim 1, characterized in that: In S2, the extreme value conditional cumulative distribution function F of the blade performance response L (l|X i )Statistical Fitting, select Logarithmic, Weibull, or Generalized Extreme Value Distribution type.

Images