Method for determining key protection areas of wind turbine blades

Abstract

This invention provides a method for determining key protection areas of wind turbine blades, comprising the following steps: S1. Obtaining a wind turbine blade model under actual operating conditions and importing the model into Ansys Fluent software, where the blade model is divided into several grids; S2. Obtaining environmental parameters of the environment where the wind turbine is located under actual operating conditions, including ambient temperature, wind speed, humidity, and rainfall, and importing the environmental parameters into Ansys Fluent software to simulate wind turbine operation; S3. Determining the icing thickness of each grid of the wind turbine blade during simulated operation; S4. Determining the chord direction of the wind turbine blade and determining the average icing thickness of multiple regions of the blade along the chord direction; S5. Determining the lift coefficient and drag coefficient of each region of the wind turbine blade based on the average icing thickness; S6. Identifying the regions corresponding to the maximum lift coefficient and the minimum drag coefficient as key protection areas.

Description

TECHNICAL FIELD

[0001] The application relates to a fan blade protection treatment method, in particular to a fan blade key protection area determination method. BACKGROUND

[0002] Wind power generation is widely used as a renewable clean energy, and with the popularization of wind power generation, the icing problem of wind turbine generators is increasingly prominent.

[0003] The icing of the blades of the wind power generator changes the original aerodynamic shape, affects the aerodynamic performance of the blades, and reduces the power generation efficiency of the unit. At the same time, due to the different effects of icing in different areas of the wind turbine on the aerodynamic performance of the blades, it is difficult to determine the key protection area. Therefore, in the deicing protection of the wind power generator, a comprehensive protection is generally used, which is a waste of resources. Another way is to determine the key protection area according to experience, which leads to the fact that the blades of the wind power generator are not deiced and protected in a targeted manner, so that the wind power generator is still affected by icing.

[0004] Therefore, in order to solve the above technical problems, a new technical means is needed. SUMMARY

[0005] Therefore, the purpose of the present application is to provide a fan blade key protection area determination method, which simulates the operation of the actual wind turbine blade airfoil and the actual working condition environment of the wind turbine. In the simulation operation, the icing mass and the icing thickness are used to determine the key protection area, so as to effectively improve the accuracy of the determination of the key protection area, provide accurate reference for the deicing protection of the wind turbine, facilitate the proposal of targeted deicing scheme, reduce or even avoid resource waste, and have strong adaptability.

[0006] The application provides a fan blade key protection area determination method, which comprises the following steps:

[0007] S1. Obtain the blade model of the wind turbine in the actual working condition, and import the model into the Ansys Fluent software. In the Ansys Fluent software, the blade model is divided into a plurality of grids;

[0008] S2. Obtain the environment parameters of the environment where the wind turbine is located in the actual working condition, including the environment temperature, wind speed, humidity and rainfall, and import the environment parameters into the Ansys Fluent software to simulate the operation of the wind turbine;

[0009] S3. Determine the icing thickness of each grid of the blade of the wind turbine in the simulation operation;

[0010] S4. Determine the chord direction of the wind turbine blade, and determine the icing thickness average of the blade in the chord direction of the blade in multiple regions;

[0011] S5. Determine the lift coefficient and drag coefficient of each region of the wind turbine blade based on the icing thickness average;

[0012] S6. Find the region corresponding to the maximum value of the lift coefficient and the minimum value of the drag coefficient as the key protection region.

[0013] Further, the lift coefficient and the drag coefficient are determined by the following method:

[0014] C l = cosα(C dmax sinα+A2cotα);

[0015] C d =C dmax -cosα(C dmax cosα-B2);

[0016] Where: α represents the attack angle, C l represents the lift coefficient, C d represents the drag coefficient, C dmax is the maximum drag coefficient under the current wind turbine blade airfoil; A2 and B2 are both coefficients;

[0017]

[0018] Where: α s is the stall angle of attack under the current blade airfoil, C ls is the lift coefficient under the stall angle of attack α s , C ds is the drag coefficient under the stall angle of attack α s ;

[0019] AR is the aspect ratio under the current grid, C ice is the icing thickness average in the chord direction of the current grid, d i is the icing thickness of the chord direction i grid, n is the number of grids, and c is the chord length of the current wind turbine blade airfoil.

[0020] Further, the icing thickness d i is determined by the following method:

[0021]

[0022] Where: Δt is the icing time step, ρ ice represents the icing density, m ice is the icing mass, and A is the grid length.

[0023] Further, the ice density p ice is determined by the following formula:

[0024]

[0025] wherein: V ∞ represents the combined velocity of the incoming flow of the wind turbine; MVD is the median diameter of the water droplets;

[0026] ΔT s-f = T s -T f , T s is the surface temperature of the blade of the wind turbine, and T f is the freezing temperature of the water droplets.

[0027] Further, the ice mass m ice is determined by the following method:

[0028] The mass balance equation in the grid is constructed as follows:

[0029]

[0030] wherein: m im is the mass of the supercooled water droplets impacting the grid, m in and m out are the mass of the water flowing into the grid and the mass of the water flowing out of the grid, respectively, and m eva is the mass of the evaporated water;

[0031] m im = β · LWC · V ∞ · A (2) ;

[0032] m out = (1-f) ·∑m Add -m eva (3) ;

[0033]

[0034] wherein: f is the ice formation ratio of the water droplets in the grid; L ew is the Lewis number, P sat-Ts and P sat-T are the saturated water vapor pressures corresponding to the blade grid surface temperature T s and the ambient temperature T

[0035] The energy balance equation in the grid is constructed as follows:

[0036] ∑Q grid-in = Q im + Q in (5) ;

[0037] ∑Q Loss =Q eva +Q i +Q out +Q f +Q c (6);

[0038] ∑Q grid-in =∑Q Loss (7);

[0039] Qim represents the energy generated by the impact of supercooled water droplets on the mesh, expressed as:

[0040]

[0041] Q in The energy brought by the water flowing into the grid is represented as:

[0042] Q in =m in C w (T s(i-1) -T f (9);

[0043] Q eva The energy carried away by water evaporation is represented as:

[0044] Q eva =m eva [C w ΔT s-f +L e (10);

[0045] Q out The energy carried away by the water flowing out of the grid is represented as:

[0046] Q out =m out C w ΔT s-f (11);

[0047] Q f The heat lost due to friction is expressed as:

[0048]

[0049] Q c The heat carried away by natural convection is represented as:

[0050] Q c =-h c AΔT s (13);

[0051] T f T is the freezing temperature of water droplets. sWhere is the blade surface temperature, T is the ambient temperature, and ΔT s-f =T s -T f ,△T s =TT s C a C w C i Specific heats of air, water, and ice, respectively; L f L e These are the latent heat of fusion and the latent heat of vaporization, respectively; h c Let be the convective heat transfer coefficient. r is the restitution coefficient for kinetic heating, expressed as:

[0052] When Pr is the Prandtl number, n0 = 1 / 2 when the grid boundary layer is in the laminar region and 1 / 3 when it is in the turbulent region.

[0053] By combining formulas (1) and (13), the mass of ice accumulation m can be calculated. ice .

[0054] The beneficial effects of this invention are as follows: By simulating the operation of the actual wind turbine blades and the actual working environment of the wind turbine, the key protection areas can be determined by the icing quality and icing thickness during the simulation operation. This can effectively improve the accuracy of the determination of key protection areas, provide accurate reference for the de-icing protection of wind turbines, facilitate the development of targeted de-icing solutions, reduce or even avoid resource waste, and has strong adaptability. Attached Figure Description

[0055] The present invention will be further described below with reference to the accompanying drawings and embodiments:

[0056] Figure 1 This is a flowchart of the present invention.

[0057] Figure 2 This is a schematic diagram of the wind turbine blades of the present invention.

[0058] Figure 3 This is a diagram illustrating the distribution of cells. Detailed Implementation

[0059] The present invention will be further described in detail below:

[0060] The present invention provides a method for determining the key protection area of ​​wind turbine blades, comprising the following steps:

[0061] S1. Obtain the blade model of the wind turbine under actual working conditions and import the model into Ansys Fluent software. In Ansys Fluent software, divide the blade model into several meshes.

[0062] S2. Obtain environmental parameters of the wind turbine's environment under actual operating conditions, including ambient temperature, wind speed, humidity, and rainfall, and import these environmental parameters into Ansys Fluent software to simulate wind turbine operation;

[0063] S3. Determine the icing thickness of the wind turbine blades for each grid during simulated operation;

[0064] S4. Determine the chord length direction of the wind turbine blades and determine the average icing thickness of multiple regions of the blades along the chord length direction; (e.g.) Figure 2 As shown, the chord length direction represents the straight-line distance from the tip to the root of the wind turbine blade. However, since some blades are not regularly shaped, when dividing the cells, a group of cells along the blade length direction is actually distributed along a curve. However, during calculation, the chord length direction, represented by the straight line in the figure, is used to divide the region, resulting in the following... Figure 3 As shown: Figure 3 The image shows a partial schematic of the blade surface unfolded into a plane. Each square represents a grid. Assuming the left side represents the blade root and the right side represents the blade tip, the unfolded surface along the actual chord length resembles a curve. Therefore, the area corresponding to this chord length line is... Figure 3 When calculating the average ice thickness of the grid area through which the curve passes, the average ice thickness of the grid area through which the curve passes is used as the standard. Therefore, starting from the root of the leaf and moving to the top, there will be multiple such corresponding areas. The corresponding calculation is performed on each area to determine the final key area.

[0065] S5. Determine the lift coefficient and drag coefficient of each region of the wind turbine blade based on the average icing thickness;

[0066] S6. Identify the areas corresponding to the maximum lift coefficient and minimum drag coefficient as key protection areas. Using the above method, simulated operation is conducted based on the airfoil of the actual wind turbine blades and the actual operating environment of the wind turbine. During the simulation, the key protection areas are determined by the icing quality and thickness, effectively improving the accuracy of key protection area identification. This provides an accurate reference for wind turbine de-icing protection, facilitates the development of targeted de-icing solutions, reduces or even avoids resource waste, and demonstrates strong adaptability.

[0067] In this embodiment, the lift coefficient and drag coefficient are determined by the following method:

[0068] C l =cosα(C dmax sinα+A2cotα);

[0069] Cd =C dmax -cosα(C dmax cosα-B2);

[0070] Where: α represents the angle of attack, C l C represents the lift coefficient. d C represents the drag coefficient. dmax A1 represents the maximum drag coefficient under the current wind turbine blade airfoil; A2 and B2 are both coefficients.

[0071]

[0072] Where: α s For the current blade airfoil, C is the stall angle of attack. ls Stall angle of attack α s The lift coefficient, C ds Stall angle of attack α s The drag coefficient below;

[0073] AR is the aspect ratio under the current grid, and C is... ice This represents the average ice thickness along the chord length of the current grid. d i Let represent the icing thickness of the i grids along the chord direction, n be the number of grids, and c be the current airfoil chord length of the wind turbine blade.

[0074] In this embodiment, the ice thickness d is determined by the following method. i :

[0075]

[0076] Where: Δt is the icing time step, ρ ice The density of the ice layer is expressed in m. ice Let A be the icing mass and A be the grid length.

[0077] In this embodiment, the ice density ρ ice Determined by the following formula:

[0078]

[0079] Where: V ∞ The incoming flow velocity of the wind turbine is represented; MVD is the median diameter of the water droplet.

[0080] ΔT s-f =T s -T f T s T represents the surface temperature of the wind turbine blades. f This is the freezing temperature of the water droplet.

[0081] In this embodiment, the ice mass m is determined by the following method. ice :

[0082] Construct the mass balance equation within the mesh:

[0083]

[0084] Where: m im For the mass of the supercooled water droplet impacting the grid, m in and m out These represent the water mass flowing into the grid and the water mass flowing out to the next grid, respectively, m. eva The mass of water evaporated;

[0085] m im =β·LWC·V ∞ •A (2);

[0086] m out = (1-f)·∑m Add -m eva (3);

[0087]

[0088] Where: f is the rate of water droplet freezing within the grid; L ew Let P be a Lewis number. sat-Ts and P sat-T The surface temperature T of the blade mesh is respectively s The saturated water vapor pressure corresponding to the ambient temperature T;

[0089] Construct the energy balance equation within the grid:

[0090] ∑Q grid-in =Q im +Q in (5);

[0091] ∑Q Loss =Q eva +Q i +Q out +Q f +Q c (6);

[0092] ∑Q grid-in =∑Q Loss (7);

[0093] Qim represents the energy generated by the impact of supercooled water droplets on the mesh, expressed as:

[0094]

[0095] Q inThe energy brought by the water flowing into the grid is represented as:

[0096] Q in =m in C w (T s(i-1) -T f (9);

[0097] Q eva The energy carried away by water evaporation is represented as:

[0098] Q eva =m eva [C w ΔT s-f +L e (10);

[0099] Q out The energy carried away by the water flowing out of the grid is represented as:

[0100] Q out =m out C w ΔT s-f (11);

[0101] Q f The heat lost due to friction is expressed as:

[0102]

[0103] Q c The heat carried away by natural convection is represented as:

[0104] Q c =-h c AΔT s (13);

[0105] T f T is the freezing temperature of a water droplet. s Where is the blade surface temperature, T is the ambient temperature, and ΔT s-f =T s -T f ,△T s =TT s C a C w C i Specific heats of air, water, and ice, respectively; L f L e These are the latent heat of fusion and the latent heat of vaporization, respectively; h c Let be the convective heat transfer coefficient. r is the restitution coefficient for kinetic heating, expressed as:

[0106] When Pr is the Prandtl number, n0 = 1 / 2 when the grid boundary layer is in the laminar region and 1 / 3 when it is in the turbulent region.

[0107] By combining formulas (1) and (13), the mass of ice accumulation m can be calculated. ice .

[0108] Finally, it should be noted that the above embodiments are only used to illustrate the technical solutions of the present invention and are not intended to limit it. Although the present invention has been described in detail with reference to preferred embodiments, those skilled in the art should understand that modifications or equivalent substitutions can be made to the technical solutions of the present invention without departing from the spirit and scope of the technical solutions of the present invention, and all such modifications or substitutions should be covered within the scope of the claims of the present invention.

Claims

1. A method for determining the key protection area of ​​wind turbine blades, characterized in that: Includes the following steps: S1. Obtain the blade model of the wind turbine under actual working conditions and import the model into Ansys Fluent software. In Ansys Fluent software, divide the blade model into several meshes. S2. Obtain environmental parameters of the wind turbine's environment under actual operating conditions, including ambient temperature, wind speed, humidity, and rainfall, and import these environmental parameters into Ansys Fluent software to simulate wind turbine operation; S3. Determine the icing thickness of the wind turbine blades for each grid during simulated operation; S4. Determine the chord length direction of the wind turbine blades and determine the average ice thickness of multiple regions of the blades along the chord length direction; S5. Determine the lift coefficient and drag coefficient of each region of the wind turbine blade based on the average icing thickness; S6. Identify the areas corresponding to the maximum lift coefficient and the minimum drag coefficient as key protection areas; The lift coefficient and drag coefficient are determined by the following method: ； ； in: Indicates angle of attack. Indicates the lift coefficient. Indicates the drag coefficient. This represents the maximum drag coefficient under the current wind turbine blade airfoil. and All are coefficients; ； ； in: This represents the stall angle of attack for the current blade airfoil. For stall angle of attack The lift coefficient at the bottom, For stall angle of attack The drag coefficient below; AR represents the aspect ratio under the current grid. This represents the average ice thickness along the chord length of the current grid. , Let n be the ice thickness of the i grid cells along the chord length, and n be the number of grid cells. This refers to the chord length of the current wind turbine blade airfoil.

2. The method for determining the key protection area of ​​wind turbine blades according to claim 1, characterized in that: The ice thickness is determined using the following method. : ； in: For the icing time step, Indicates the density of the ice layer. For icing quality, This represents the grid length.

3. The method for determining the key protection area of ​​wind turbine blades according to claim 2, characterized in that: The ice density Determined by the following formula: ； in: The incoming flow velocity of the wind turbine is represented; MVD is the median diameter of the water droplet. , The surface temperature of the wind turbine blades. This is the freezing temperature of the water droplet.

4. The method for determining the key protection area of ​​wind turbine blades according to claim 2, characterized in that: The ice accumulation quality is determined using the following method. : Construct the mass balance equation within the mesh: （1）； Where: m im For the mass of the supercooled water droplet impacting the grid, m in and m out These represent the water mass flowing into the grid and the water mass flowing out to the next grid, respectively, m. eva The mass of water evaporated; （2）； （3）； （4）； in: The freezing ratio of water droplets within the grid; Let P be a Lewis number. sat-Ts and P sat-T The surface temperature T of the blade mesh is respectively s The saturated water vapor pressure corresponding to the ambient temperature T; Construct the energy balance equation within the grid: （5）； （6）； （7）； Qim represents the energy generated by the impact of supercooled water droplets on the mesh, expressed as: （8）； Q in The energy brought by the water flowing into the grid is represented as: （9）； Q eva The energy carried away by water evaporation is represented as: （10）； Q out The energy carried away by the water flowing out of the grid is represented as: （11）； Q f The heat lost due to friction is expressed as: （12）； Q c The heat carried away by natural convection is represented as: （13）； T f T is the freezing temperature of a water droplet. s Where is the blade surface temperature, T is the ambient temperature, and ΔT s-f =T s -T f , △T s =TT s C a C w C i Specific heats of air, water, and ice, respectively; L f L e These are the latent heat of fusion and the latent heat of vaporization, respectively; h c Let be the convective heat transfer coefficient, and r be the restitution coefficient for kinetic heating, expressed as: ; For Prandtl number, when the grid boundary layer is in the laminar region, n0 = 1 / 2; when it is in the turbulent region, n0 = 1 / 3. By combining formulas (1) and (13), the mass of ice accumulation can be calculated. .

Images