CN117744409A - Method and system for predicting blade deformation and impeller hub load of offshore floating fan - Google Patents
Method and system for predicting blade deformation and impeller hub load of offshore floating fan Download PDFInfo
- Publication number
- CN117744409A CN117744409A CN202410182336.8A CN202410182336A CN117744409A CN 117744409 A CN117744409 A CN 117744409A CN 202410182336 A CN202410182336 A CN 202410182336A CN 117744409 A CN117744409 A CN 117744409A
- Authority
- CN
- China
- Prior art keywords
- blade
- impeller
- model
- wind speed
- aerodynamic
- Prior art date
- Legal status (The legal status is an assumption and is not a legal conclusion. Google has not performed a legal analysis and makes no representation as to the accuracy of the status listed.)
- Pending
Links
- 238000000034 method Methods 0.000 title claims abstract description 96
- 238000007667 floating Methods 0.000 title claims abstract description 29
- 238000004458 analytical method Methods 0.000 claims abstract description 45
- 238000010168 coupling process Methods 0.000 claims abstract description 22
- 230000008878 coupling Effects 0.000 claims abstract description 19
- 238000005859 coupling reaction Methods 0.000 claims abstract description 19
- 238000012545 processing Methods 0.000 claims abstract description 17
- 241001130943 Phyllanthus <Aves> Species 0.000 claims abstract description 5
- 239000002184 metal Substances 0.000 claims abstract description 4
- 239000011159 matrix material Substances 0.000 claims description 81
- 230000000694 effects Effects 0.000 claims description 56
- 230000033001 locomotion Effects 0.000 claims description 42
- 230000008569 process Effects 0.000 claims description 36
- 239000013598 vector Substances 0.000 claims description 30
- 230000008859 change Effects 0.000 claims description 24
- 238000009826 distribution Methods 0.000 claims description 21
- 238000006073 displacement reaction Methods 0.000 claims description 19
- 230000001133 acceleration Effects 0.000 claims description 15
- 230000005484 gravity Effects 0.000 claims description 13
- 238000010008 shearing Methods 0.000 claims description 12
- 230000010354 integration Effects 0.000 claims description 9
- 238000004364 calculation method Methods 0.000 description 24
- 238000013016 damping Methods 0.000 description 12
- 238000010586 diagram Methods 0.000 description 7
- 230000006698 induction Effects 0.000 description 7
- 230000009471 action Effects 0.000 description 6
- 238000012937 correction Methods 0.000 description 5
- 238000013461 design Methods 0.000 description 5
- KFLQGJQSLUYUBF-WOJBJXKFSA-N Phyllanthin Chemical compound C([C@H](COC)[C@@H](COC)CC=1C=C(OC)C(OC)=CC=1)C1=CC=C(OC)C(OC)=C1 KFLQGJQSLUYUBF-WOJBJXKFSA-N 0.000 description 4
- KFLQGJQSLUYUBF-PMACEKPBSA-N Phyllanthin Natural products C([C@@H](COC)[C@H](COC)CC=1C=C(OC)C(OC)=CC=1)C1=CC=C(OC)C(OC)=C1 KFLQGJQSLUYUBF-PMACEKPBSA-N 0.000 description 4
- 230000005540 biological transmission Effects 0.000 description 4
- 238000012905 input function Methods 0.000 description 4
- 230000000737 periodic effect Effects 0.000 description 4
- 238000005452 bending Methods 0.000 description 3
- 238000006243 chemical reaction Methods 0.000 description 3
- 238000013500 data storage Methods 0.000 description 3
- 238000011161 development Methods 0.000 description 3
- 230000008676 import Effects 0.000 description 3
- 230000004048 modification Effects 0.000 description 3
- 238000012986 modification Methods 0.000 description 3
- 235000017166 Bambusa arundinacea Nutrition 0.000 description 2
- 235000017491 Bambusa tulda Nutrition 0.000 description 2
- 241001330002 Bambuseae Species 0.000 description 2
- 235000015334 Phyllostachys viridis Nutrition 0.000 description 2
- 230000003044 adaptive effect Effects 0.000 description 2
- 239000011425 bamboo Substances 0.000 description 2
- 230000007547 defect Effects 0.000 description 2
- 230000001934 delay Effects 0.000 description 2
- 239000000463 material Substances 0.000 description 2
- 238000011160 research Methods 0.000 description 2
- 238000005070 sampling Methods 0.000 description 2
- 238000004088 simulation Methods 0.000 description 2
- 238000009825 accumulation Methods 0.000 description 1
- 230000009286 beneficial effect Effects 0.000 description 1
- 238000012512 characterization method Methods 0.000 description 1
- 230000001351 cycling effect Effects 0.000 description 1
- 230000005284 excitation Effects 0.000 description 1
- 238000010438 heat treatment Methods 0.000 description 1
- 244000144972 livestock Species 0.000 description 1
- 230000014759 maintenance of location Effects 0.000 description 1
- 230000007246 mechanism Effects 0.000 description 1
- 238000002620 method output Methods 0.000 description 1
- 239000000203 mixture Substances 0.000 description 1
- 238000005457 optimization Methods 0.000 description 1
- 238000005381 potential energy Methods 0.000 description 1
- 238000010248 power generation Methods 0.000 description 1
- 238000005096 rolling process Methods 0.000 description 1
- 239000007787 solid Substances 0.000 description 1
- 238000001228 spectrum Methods 0.000 description 1
- 230000004083 survival effect Effects 0.000 description 1
- 239000000725 suspension Substances 0.000 description 1
- 238000012546 transfer Methods 0.000 description 1
- 230000009466 transformation Effects 0.000 description 1
Landscapes
- Wind Motors (AREA)
Abstract
The invention discloses a method and a system for predicting blade deformation and impeller hub load of an offshore floating fan, wherein the predicting method comprises the following steps: s1, acquiring initial parameters; s2, processing the blade parameters based on a geometric precise beam theory considering geometric large deformation, and constructing an impeller model; s3, processing the wind speed parameter based on an unsteady phyllanthus metal momentum theory, and establishing a aerodynamic model, wherein the aerodynamic model is used for calculating aerodynamic load; s4, constructing a pneumatic-elastic coupling model by combining the impeller model and the aerodynamic model, applying the pneumatic load to the impeller, analyzing pneumatic-elastic characteristics, and outputting analysis results, wherein the analysis results comprise blade load, blade deformation and impeller load. The invention has the characteristics of high precision and high efficiency.
Description
Technical Field
The invention relates to the technical field of wind driven generators, in particular to a method and a system for predicting blade deformation and impeller hub load of an offshore floating fan.
Background
The energy source is an important material guarantee for human survival and civilization development, and is an important mark of the motive power and scientific technical level of national economic development. Wind energy is caused by uneven heating of the atmosphere, so that the wind energy is clean and environment-friendly, and has huge reserves. Wind power generators have also gained substantial development as a primary device for converting wind energy into electrical energy. With the enlargement of the offshore floating fan, the diameter of the impeller is increased, the flexibility of the blades is increased, the rotational inertia and the mass are increased, and in addition, the working environment of the fan is bad (unsteady pneumatic load and platform movement), so that the deformation prediction difficulty of the blades of the fan is increased, and the situation is particularly good under the marine environment. Meanwhile, under the action of non-uniform wind speed, three blades of the fan impeller generate larger impeller hub non-torque load, and a larger challenge is brought to the operation of a fan transmission system. The deformation of the large-scale offshore floating fan blade is accurately predicted, so that a solid foundation is laid for the load analysis and design of the blade; the method can accurately predict the hub load of the impeller and effectively ensure the internal structural layout and design analysis of the fan transmission system. The existing fan blade deformation analysis method is mostly based on finite element or rigid blade analysis, wherein the finite element method is low in calculation efficiency. In addition, the rigidity blade ignores the flexibility of the blade, the prediction precision of the large-scale blade is low, the existing model is often analyzed for a single blade, the combined action among three blades is ignored, and the hub load of the impeller cannot be predicted accurately.
Disclosure of Invention
The invention aims to: the invention aims to provide a method and a system for predicting the deformation of a blade and the hub load of an impeller of an offshore floating fan, which have high efficiency and high precision.
The technical scheme is as follows: the invention relates to a method for predicting the deformation of a blade and the hub load of an impeller of a marine floating fan, which comprises the following steps:
s1, acquiring initial parameters, wherein the initial parameters comprise blade parameters and wind speed parameters;
s2, processing the blade parameters based on a geometric precise beam theory considering geometric large deformation, and constructing an impeller model, wherein the impeller model comprises a plurality of blade models; the impeller model corresponds to one impeller, the blade model corresponds to one blade, and each impeller is provided with a plurality of blades;
s3, processing the wind speed parameter based on an unsteady phyllanthus metal momentum theory, and establishing a aerodynamic model, wherein the aerodynamic model is used for calculating aerodynamic load;
s4, constructing a pneumatic-elastic coupling model by combining the impeller model and the aerodynamic model, applying the pneumatic load to the impeller, analyzing pneumatic-elastic characteristics, and outputting analysis results, wherein the analysis results comprise blade load, blade deformation and impeller load.
Further, the geometric precise beam model comprises a plurality of geometric precise beam units and beam nodes, each blade is scattered into a plurality of geometric precise beam units, and the plurality of blades are connected by the beam nodes; the geometric precise beam unit is constructed by adopting three nodes, and the formula is as follows:
;
In the method, in the process of the invention,as a unit-shaped function;/>Is a third-order unit vector; />Is a node 1 shape function; />Is a node 2-shaped function; />Is a node 3-shaped function.
Further, a geometrically accurate beam model is established to simulate the deformation of the blade, the geometrically accurate beam model formula is as follows:
in the method, in the process of the invention,is a beam unit mass matrix;the acceleration matrix is a beam node acceleration matrix;the gyro matrix is a beam unit gyro matrix;is a beam node velocity matrix;is a beam unit stiffness matrix;the displacement matrix is a beam node displacement matrix;the gravity matrix is born by the beam unit;liang Shan of a shape of Liang ShanA matrix of external forces;is a beam cell internal force matrix.
Further, in the unsteady phyllanthin momentum theory, wind speed distribution on the blades is uneven due to unsteady factors.
Further, the wind speed parameters comprise a wind speed size parameter and a wind speed space distribution parameter, and the wind speed is divided into a basic stable wind speed, a gust wind speed, a gradual change wind speed and a turbulent wind speed according to the wind speed size parameter; dividing the unsteady factors into wind shearing effect, tower shadow effect, yaw effect and basic movement of the wind driven generator according to the wind speed space distribution parameters.
Further, in step S3, calculating the aerodynamic load includes the following formula:
;
In the method, in the process of the invention,representing air density; />Representing the cross-sectional chord length of the blade; />A drag coefficient representing a blade section; />A lift coefficient representing a blade section; />A twist coefficient representing a cross section of the blade; />Is corresponding axial aerodynamic force; />For corresponding cuttingAerodynamic force to the direction; />Torque caused by the offset of the aerodynamic center relative to the shear center; />Lift per unit blade length; />Representing the resistance per unit blade length; />Is a relative velocity vector; />The final real wind speed and the included angle of the rotating plane are obtained.
Further, a loose coupling method is adopted to establish a pneumatic-elastic coupling model, and data between the impeller model and the aerodynamic model are exchanged at each time step.
Further, in step S4, the reaction is based on the generalized-And the time integration method is used for analyzing the aeroelastic characteristics and solving by adopting a self-adaptive time step strategy.
Further, when the aeroelastic characteristic analysis is carried out, the gravity factors of the blades and the influence of a control system are considered.
The technical scheme is as follows: the invention relates to a marine floating fan blade deformation and impeller hub load prediction system, which comprises:
the data acquisition module is used for acquiring initial parameters, wherein the initial parameters comprise blade parameters and wind speed parameters;
The impeller model building module is used for processing the blade parameters based on a geometric precise beam theory considering geometric large deformation to build an impeller model, and the impeller model comprises a plurality of blade models; the impeller model corresponds to one impeller, the blade model corresponds to one blade, and each impeller is provided with a plurality of blades;
the aerodynamic model building module is used for processing the wind speed parameter based on the unsteady leaf element momentum theory and building an aerodynamic model which is used for calculating aerodynamic load;
and the aeroelastic analysis module is used for combining the impeller model and the aerodynamic model to construct an aeroelastic coupling model, applying the aerodynamic load to the impeller, carrying out aeroelastic characteristic analysis, and outputting analysis results, wherein the analysis results comprise blade load, blade deformation and impeller load.
The beneficial effects are that: the invention has the following remarkable effects: 1. the precision is high: on one hand, the constructed large-scale fan blade model takes the influence of unsteady inflow (wind shearing, tower shadow effect, yaw inflow and the like) into consideration, so that the typical running wind environment of the existing fan is comprehensively covered; on the other hand, the constructed fan blade aeroelastic coupling dynamic model considers the influence of the motion of the offshore floating fan matrix, and meanwhile, the constructed fan blade coupling analysis model has excellent data processing flow, high precision and strong usability; 2. the efficiency is high: on one hand, when modeling the blades, a complete three-blade impeller model is directly built, so that the defects of the existing single-blade modeling are overcome, the impeller hub load is directly output, and input conditions are provided for the design of a transmission system in the fan; on the other hand, the provided prediction method can operate only by providing basic blade parameters and wind speed parameters, is not limited by the technical foundation of operators, and is provided with corresponding operation interfaces in a prediction system, so that the calculation efficiency is greatly improved on the premise of ensuring the accuracy determination.
Drawings
FIG. 1 is a general flow chart of a prediction method provided by the invention;
FIG. 2 is a schematic diagram of a geometrically precise beam state;
FIG. 3 is a schematic illustration of a three-bladed impeller model;
FIG. 4 is a schematic diagram of aerodynamic flow rate and aerodynamic force under a fixed airfoil;
FIG. 5 is a schematic illustration of tower shadow effect on wind speed;
FIG. 6 is a schematic diagram of wake deflection under yaw effect;
FIG. 7 is a six degree of freedom motion schematic of a blower;
FIG. 8 is a schematic diagram of a pneumatic-elastic coupling model structure;
FIG. 9 is a schematic diagram of software functionality presented based on a predictive system;
FIG. 10 is a schematic diagram of an NREL 5-MW fan blade airfoil distribution;
FIG. 11 is a welcome function interface schematic;
FIG. 12 is a schematic view of a blade structural attribute interface;
FIG. 13 is a schematic view of a blade aerodynamic interface;
FIG. 14 is a schematic view of an airfoil aerodynamic performance interface;
FIG. 15 is a steady wind parameter setting interface;
FIG. 16 is an gust parameter setting interface;
FIG. 17 is a gradual wind parameter setting interface;
FIG. 18 is a turbulent wind parameter setting interface;
FIG. 19 is a wind shear parameter setting interface;
FIG. 20 is a tower shadow effect parameter setup interface;
FIG. 21 is a yaw effect parameter set interface;
FIG. 22 is a matrix motion parameter setting interface;
FIG. 23 is a calculation parameters function interface;
FIG. 24 is a data retention interface;
FIG. 25 is a pneumatic properties interface in a results display;
FIG. 26 is a blade load attribute interface in a results display;
FIG. 27 is a blade deformation interface in a results display;
fig. 28 is an impeller load interface in a results display.
Detailed Description
The invention is further elucidated below in connection with the drawings and the detailed description.
Referring to FIG. 1, the invention discloses a method and a system for predicting blade deformation and impeller hub load of an offshore floating fan. The method for predicting the deformation of the blade and the hub load of the impeller of the offshore floating fan comprises the following steps:
s1, acquiring initial parameters, wherein the initial parameters comprise blade parameters and wind speed parameters;
s2, processing the blade parameters based on a geometric precise beam theory considering geometric large deformation, and constructing an impeller model, wherein the impeller model comprises a plurality of blade models; the impeller model corresponds to one impeller, the blade model corresponds to one blade, and each impeller is provided with a plurality of blades;
s3, processing the wind speed parameter based on an unsteady phyllanthus metal momentum theory, and establishing a aerodynamic model, wherein the aerodynamic model is used for calculating aerodynamic load;
S4, constructing a pneumatic-elastic coupling model by combining the impeller model and the aerodynamic model, applying the pneumatic load to the impeller, analyzing pneumatic-elastic characteristics, and outputting analysis results, wherein the analysis results comprise blade load, blade deformation, impeller load and pneumatic attribute.
An offshore floating wind turbine blade deformation and impeller hub load prediction system comprising:
the data acquisition module is used for acquiring initial parameters, wherein the initial parameters comprise blade parameters and wind speed parameters;
the impeller model building module is used for processing the blade parameters based on a geometric precise beam theory considering geometric large deformation to build an impeller model, and the impeller model comprises a plurality of blade models; the impeller model corresponds to one impeller, the blade model corresponds to one blade, and each impeller is provided with a plurality of blades;
the aerodynamic model building module is used for processing the wind speed parameter based on the unsteady leaf element momentum theory and building an aerodynamic model which is used for calculating aerodynamic load;
and the aeroelastic analysis module is used for combining the impeller model and the aerodynamic model to construct an aeroelastic coupling model, applying the aerodynamic load to the impeller, carrying out aeroelastic characteristic analysis, and outputting analysis results, wherein the analysis results comprise blade load, blade deformation, impeller load and aerodynamic properties.
In this embodiment, the blade is large-sized and flexible.
In the step S1, in the data acquisition module, the blade parameters comprise blade structure attribute parameters, blade aerodynamic attribute parameters and airfoil aerodynamic performance parameters. The blade structure attribute parameters mainly define the structure of the blade and are used for constructing a blade model; the blade aerodynamic property parameters are mainly used for defining the wing profiles of the blades at the lengths of all axes of the blades, the wing profile aerodynamic property parameters are mainly used for defining aerodynamic variables of all wing profiles, and the blade aerodynamic property parameters and the wing profile aerodynamic property parameters lay a foundation for calculation of aerodynamic forces in aerodynamic force models. The wind speed parameters include a wind speed magnitude parameter and a wind speed spatial distribution parameter. The wind speed is divided into a basic stable wind speed, a gust wind speed, a gradual change wind speed and a turbulent wind speed according to the wind speed size parameter; according to the wind speed space distribution parameters, the unsteady factors (concepts in the unsteady leaf element momentum theory) are divided into wind shearing effect, tower shadow effect, yaw effect and basic movement of the wind driven generator.
S2, in the impeller model building module, the geometric precise beam model comprises a plurality of geometric precise beam units and beam nodes, each blade is scattered into a plurality of geometric precise beam units, and the plurality of blades are connected through the beam nodes; the geometrically accurate beam unit is constructed by adopting three nodes, and the formula is as follows:
In the method, in the process of the invention,is a unit-shaped function;is a third-order unit vector;is a node 1 shape function;is a node 2-shaped function;is a node 3-shaped function.
Establishing a geometrically accurate beam model to simulate the deformation of the blade, wherein the geometrically accurate beam model has the following formula:
in the method, in the process of the invention,is a beam unit mass matrix;the acceleration matrix is a beam node acceleration matrix;the gyro matrix is a beam unit gyro matrix;is a beam node velocity matrix;is a beam unit stiffness matrix;the displacement matrix is a beam node displacement matrix;the gravity matrix is born by the beam unit;an external force matrix of the beam unit;is a beam cell internal force matrix.
The process for establishing the geometrically accurate beam model is as follows:
step 1), referring to fig. 2, the original state and the deformed state of the geometrically precise beam are first determined. First, an inertial coordinate system o is established 1 i 1 i 2 i 3 . Then build up in the initial state of the beamInitial coordinate system O 1 b 1 b 2 b 3 And defines a local coordinate b i To describe the original state of the beam. Wherein b 1 Is the unit tangent vector of the reference curve. Then establishing a deformation coordinate system O in a deformation state 2 B 1 B 2 B 3 And defining the deformed beam from the original local coordinates b i Changing to a new local coordinate B i New local coordinates B i From displacement vectors u (x 1 ) And the cross-sectional rotation tensor R. Notably, in the deformation coordinate system O 2 B 1 B 2 B 3 In O 2 B 1 The direction of (a) is not necessarily exactly the tangent vector of the deformed reference curve, the rotation tensor R 0 And R, all in inertial coordinates (oi 1 i 2 i 3 ) Is established for reference. According to the geometrically accurate beam, determining the position vector of any point P' in the deformed beam is as follows:
(1)
wherein x is 2 And x 3 Is a local coordinate system O 2 B 1 B 2 B 3 The coordinate values below; w (w) i For a small warping displacement function, calculating according to the geometric characteristic and the rigidity characteristic of the cross section, namely determining by solving a two-dimensional finite element problem on the cross section; r is the reference point O 1 Is used for the initial position vector of (a); r is R 0 Is O 1 The rotation tensor of the beam section at the position is not an identity matrix because the beam model generally has initial bending and torsion; x is x 1 To reference the line arc length along the beam with initial bending and torsion;i is 1,2 and 3 which are original local coordinates and represent three directions of coordinates; b (B) i I is 1,2 and 3, representing three directions of coordinates; u (x) 1 ) Is a displacement vector; and R is a section rotation tensor.
Step 2), calculating the section strain and curvature of the beam model according to the variable amount of the formula (1), and further solving the strain energy variation contained in the beam unit; calculating the inertial speed of the beam section by the variable of the formula (1), and further obtaining the functional variation of the beam unit; and calculating the virtual work of the beam section. And processing strain energy variation, kinetic energy variation and virtual work based on the Hamiltonian principle, and obtaining a motion control equation of the geometrically accurate beam through integration. The process is specifically as follows.
Step 21): the section strain, curvature and strain energy variation formulas of the beam model are calculated as follows:
(2)
(3)
wherein, κ is the curvature of the section;for r to x 1 Seeking a derivative;for u to x 1 Seeking a derivative;strain energy per unit length;is of unit length;An antisymmetric tensor of r;is thatIs an antisymmetric tensor of (2);an antisymmetric tensor that is a virtual rotation vector;ϵis the section strain; 1 is the direction vector of the integral inertial coordinate system 1; the symbol ()' represents the pair x 1 Seeking a derivative; the symbol axial (∙) is denoted as an antisymmetric tensorWhen vector a= (a) 1 , a 2 , a 3 ,) T The corresponding antisymmetric tensor is;Is a virtual rotation vector, andthe method comprises the steps of carrying out a first treatment on the surface of the F is Liang Jie face force matrix under the whole coordinate system; m is a beam section mass matrix under the whole coordinate system. Wherein the tensor is antisymmetricCan be expressed as:
(4)
(5)
wherein F is s And M s Respectively the force and moment of the beam section; c (C) * The beam section rigidity matrix is obtained by adopting two-dimensional finite element analysis for the beam section rigidity matrix under a local material coordinate system;is a beam section stiffness matrix under an integral inertial coordinate system.
Step 22): calculating the inertial speed of the beam section, and further obtaining the kinetic energy variation of the beam unit as follows:
(6)
wherein h and g are respectively the linear momentum and the angular momentum of the lower beam section of the inertial coordinate system; Kinetic energy per unit length;u is the time derivative per unit length;an antisymmetric tensor of the time derivative of u;is the time derivative of the virtual rotation vector. () means deriving time. Wherein h and g are represented as:
(7)
wherein M is a beam section mass matrix under an integral coordinate system;is an angular velocity vector. Wherein, the beam section mass matrix M is as follows:
(8)
wherein m is the mass per unit length of the beam unit;is thatIs represented by a tensor of (c),,η * to the centroid of the cross-section relative to the point O in the local coordinate system 2 Is a position vector of (2);,Γ * the moment of inertia of the lower section of the local coordinate system;is a third order unit vector.
Step 23): and calculating the virtual work of the beam section, wherein the formula is as follows:
(9)
wherein f ext And m ext Respectively externally distributed forces and moments on the beam section;is the imaginary work of the beam section.
Step 24): and processing the strain energy variation, the kinetic energy variation and the virtual work generated in the steps 21) to 23) based on the Hamiltonian principle, and obtaining a motion control equation of the geometrically accurate beam through integration:
(10)
the hamiltonian principle is formulated as follows:
(11)
in the method, in the process of the invention,is the end time;for a start time;is the unit length;is the kinetic energy of the beam unit;is the potential energy of the beam unit.
Step 25): for ease of application, a matrix form of the motion control equations of geometrically accurate beams is proposed. In this embodiment, a generalized- α numerical integration method is used to solve, and linearize the inertial force, elastic force, and external force in a matrix form. Wherein the matrix form of the motion control equation is as in equation 12:
(12)
Wherein F is I Is an inertial force; f (F) C And F D Is elastic force; f (F) ext Is an external force (including gravity);is F C For x 1 And (5) deriving.
In the linearization process, firstly, incremental displacement Δq, incremental speed Δv and incremental acceleration Δa are set as follows:
(13)
in the method, in the process of the invention,is incremental translational displacement;is incremental rotational displacement;is the incremental translational velocity;is the incremental rotational speed;is the incremental translational acceleration;is the incremental rotational acceleration.
Substituting equations (7) and (8) into equation (14) for linearization, and then performing matrix transformation to rewrite the linearized inertial force into a form containing generalized coordinates (equation 15):
(14)
(15)
wherein K is I Is a rigidity matrix; g I Is a gyro matrix; m is M I Is a quality matrix;to increase the force. Wherein, the rigidity matrix, the gyro matrix and the mass matrix are as follows:
(16)
in the linearization process, the elastic force can be linearized by substituting the formulas (2) and (6) into the formula (14) (formula 17):
(17)
in the method, in the process of the invention,is the increment of damping force;the rigidity matrix is a beam section rigidity matrix under an integral inertial coordinate system;for incremental displacementFor x 1 Seeking a derivative;is one of the coefficient matrices;is energy-increasing displacement;is one of the coefficient matrices;is one of the coefficient matrices. Wherein each coefficient may be expressed as:
(18)
Wherein E is 1 =r' + u';C 11 、C 21 Is part of matrix C; o (O) 12 As well as a part of matrix O. Wherein, matrix C is as follows:
(19)
step 3), in many cases, will generally note the effect of beam unit quality. Meanwhile, in order to avoid convergence of the high-frequency vibration acceleration model, damping force is often introduced. Wherein, when considering beam cell mass, a geometrically accurate Liang Zhiliang matrix is proposed:
(20)
wherein g is a gravitational acceleration vector;is the unit gravity;is the mass per unit length of the beam unit.
When a damping force is introduced, there are:
(21)
wherein mu is a damping coefficient;is a damping force;is translational damping force;is a bending damping force. Referring to the form of equation (14), the form of writing the damping force as a matrix is:
(22)
in the method, in the process of the invention,is a damping force matrix;for an increased damping force matrix.
The damping force increment is:
(23)
wherein the variables are as follows:
(24)
after introducing the damping force and gravity of equation (13), the motion control equation of the geometrically accurate beam becomes:
(25)
discretizing the control equation (25) using a finite element method:
(26)
wherein l is the length of the beam unit; n is a unitary function matrix. The function (26) is a nonlinear function and is linearized as follows:
(27)
the displacement, velocity and acceleration delta within a cell are represented by values at nodes:
(28)
In the method, in the process of the invention,is the displacement at the node;for unit displacement vs x 1 Seeking a derivative;is the node speed;the acceleration is node acceleration;as a unit form function pair x 1 And (5) deriving.
Step 4), substituting formulas (15), (17), (23) and (28) into formula (27), and finally (27) may be finished as:
(29)
in the method, in the process of the invention,is a quality matrix;the node acceleration matrix is formed;is a gyro matrix;is a node speed matrix;is a rigidity matrix;is a node displacement matrix;is a gravity matrix;is an external force matrix;is an internal force matrix. Wherein the mass matrix, gyro matrix, stiffness matrix and force matrix of the beam unit can be expressed as:
(30)
in the formula, the variable is shown in the formula (24).
Thus, a complete geometrically accurate beam model is built. The geometric precise beam model comprises a plurality of geometric precise beam units and beam nodes, each blade is scattered into the plurality of geometric precise beam units, and the plurality of blades are connected by the beam nodes. In the model, the geometric precise beam unit is constructed by adopting three nodes, and the formula is as follows:
(31)
in the method, in the process of the invention,is a unit-shaped function;is a third-order unit vector;is a node 1 shape function;is a node 2-shaped function;is a node 3-shaped function. Wherein the shape function h 1 、h 2 And h 3 The method comprises the following steps:
(32)
where s is a Gaussian point.
Wherein, according to the Euler rotation theorem: in three dimensions, this displacement is equivalent to a rotation about a fixed axis containing a fixed point, provided that the inside of a rigid body is at least somewhat stationary while it is rotating. The Wiener-Milenkovic parameter is used in this embodiment to define the finite rotation c of the beam node:
(33)
in the method, in the process of the invention,is the angle of rotation;is the unit rotation axis vector.
Rotation matrix of beam nodesCan be expressed as:
(34)
wherein c= [ c ] 1 , c 2 , c 3 ] T ;c 0 = 2-1/8×c T c;c 1 A first column vector that is a finite rotation c; c 2 A second column vector that is a finite rotation c; c 3 The third column vector, which is a finite rotation c.
Referring to fig. 3, since the blade is typically of an elongated structure, in order to balance the computational efficiency and the computational accuracy, a blade model is constructed using the geometrically accurate beam model described above. Each blade is discretized into 12 three-node geometrically precise beam units, while the unit lengths of the geometrically precise beam units at the blade root and tip are small relative to the middle portion of the blade. This is mainly because the root structural properties of the blade change faster, and smaller units are used in order to be able to describe the root structural properties more accurately; similarly, the deformation of the blade tip is most severe, the aeroelastic effect is more obvious, and the blade tip adopts smaller unit length for more accurately describing the deformation of the blade. In this embodiment, she Pianxuan uses three. In the impeller model, in order to connect three blades as a whole, the hub of the impeller is regarded as a part of the blades, and a simulation is performed with a high rigidity (rigid unit). And simultaneously, the mass of the hub is equally divided into three parts and added to the root of the blade. The support of the rotor is modeled using linear springs in-plane. The blade model is directly connected at the intermediate node o, so that the combined load (six opposite loads including three translational loads and three rotational moments) of the three blades under the action of non-uniform aerodynamic forces, namely the reaction force at the hub of the point o, can be directly obtained. The method can effectively overcome the defect that the influence of deformation of the blades is ignored by direct integral accumulation of each blade of the existing model.
In the step S3, in the aerodynamic model building module, an unsteady phyllanthin momentum theory is adopted to build an aerodynamic model, the phyllanthin momentum theory can be used for very rapidly calculating the aerodynamic load suffered by the blade, and the aerodynamic model has satisfactory calculation accuracy, and the blade is generally influenced by unsteady inflow, so that the unsteady phyllanthin momentum theory is adopted to calculate the aerodynamic load of the blade. In the unsteady phyllin momentum theory, wind speed distribution on the blades is uneven due to unsteady factors. Calculating the aerodynamic load comprises the following process of steps 1) to 5):
step 1), defining the blade section and the state of aerodynamic inflow: (a) The speed of each section of all the blades is calculated, and the calculation includes the speed of the whole rotation of the blades and the speed of the local vibration of the blades. (b) And calculating the local free flow speed of each section according to the position of the section of the blade and the spatial distribution of the wind speed at the moment. (c) Initializing or reading the induction speed of each section from the previous time step, and correcting the induction speed of the blade by taking the tip loss correction, the hub loss correction and the Ge Laowo special correction into consideration. Wherein, the formula of the induction speed of the correction blade is as follows:
(35)
(36)
in the method, in the process of the invention, The number of blades of the wind driven generator;is the radius of rotation of the vane element;negative for rotor plane edgeNormal vector of the shaft;is Ge Laowo specially modified and is an empirical expression;is lost by the blade tipAnd hub lossLoss factor coefficient of common composition:is the rotation radius of the blade tip; superscript、Respectively the tangential and axial direction of the blade element;is the induction speed;is air density;the final real wind speed and the included angle of the rotating plane are obtained;is the wind speed;is the radius of the hub; l represents the lift per unit blade length.
Step 2), referring to fig. 4, the relative speed of each section of the blade is calculated according to the calculated speed of each section of the blade in (a):
(37)
in the superscript、Respectively the tangential and axial direction of the blade element;is a relative velocity vector; the angle between the relative velocity vector and the rotor plane is defined as the inflow angle phi;is a local velocity vector of the blade, including the local vibration velocity of the blade;is the blade unit speed in the global coordinate system, including the blade rotation and vibration speeds.
Step 3), calculating attack angles alpha of all sections of the blade e :
(38)
In the method, in the process of the invention,is the total torsion angle of one section;is the initial torsion angle of the section; Is the blade pitch angle;is the blade deformation angle; the superscript ax represents the axial direction and the superscript tan represents the tangential direction. And the total torsion angle of the section is formed by combining the initial torsion angle of the section, the pitch angle of the blade and the deformation angle of the blade at the section.
Step 4), calculating the aerodynamic load of each section of the blade:
(39)
in the method, in the process of the invention,representing air density;representing the cross-sectional chord length of the blade;a drag coefficient representing a blade section;a lift coefficient representing a blade section;a twist coefficient representing a cross section of the blade;for the corresponding axial aerodynamic force (superscript ax stands for axial);aerodynamic forces for the corresponding tangential directions (superscript tan stands for tangential direction);torque caused by the offset of the aerodynamic center relative to the shear center;lift per unit blade length;representing the resistance per unit blade length;is the relative velocity vector.
The drag coefficient, lift coefficient and torsion coefficient of the section are determined by the wing profile corresponding to the drag coefficient, lift coefficient and torsion coefficient, and the drag coefficient, lift coefficient and torsion coefficient of the section show different characteristics under different attack angles of the section. In this embodiment, a table look-up method is used to determine the three coefficients according to the attack angle of the cross section. Three coefficient data corresponding to different attack angles of each airfoil of the NREL 5-MW fan blade are disclosed, and can be conveniently obtained.
Step 5), it can be known from the calculation correction formula (35) of the induced velocity and the calculation formula (37) of the relative velocity that the two equations are involved in each other and influence each other, so the equations must be solved iteratively. However, the induction speed varies relatively slowly in time, so when the calculation time step is sufficiently small relative to the aerodynamic variation time scale, the iteration is replaced by a time evolution. Thus, instead of iteratively cycling the calculated induction rate at each time step, time evolution is employed. And finally, storing and recording the pneumatic load proposed in the step 4), so that the later pneumatic-elastic coupling application is facilitated.
In addition, in the unsteady phyllanthus momentum theory, wind speed distribution on the blades is uneven due to unsteady factors. In fact, the spatial distribution of the air flow of the whole impeller is related to the atmosphere and also to the structural characteristics of the blower itself, mainly expressed by the wind shearing effect, the tower shadow effect, the yaw effect and the like. These unsteady flows can cause maldistribution of gas flow rates across the wheel space, resulting in the hub being subjected to periodic non-torque loads. The characteristic of flow velocity maldistribution caused by each unsteady factor is described in detail below.
The wind shearing effect reflects the change of wind speed along with the vertical height, so that the wind speed on the fan blade is unevenly distributed, and especially the difference between the upper wind speed and the lower wind speed of the impeller is maximum. As the rotor diameter increases, the wind shear effect becomes more pronounced. A simple and reliable wind shear model is shown below:
(40)
In the method, in the process of the invention,andwind speed at hub height and hub height, respectively;the wind shearing experience index is related to the environment of the fan, represents the resistance degree of the ground environment to the wind speed, is relatively flat to the sea level,taking 0.1, but for mountain areas or more tree areas,taking 0.24, since the analysis is of an offshore floating fan, in the subsequent analysisTaking 0.1;is the position of the blade section in the vertical direction in the inertial coordinate system.
The tower shadow effect is: for upwind wind generators, the incoming wind speed will change due to the presence of the tower. It introduces periodic unsteady flows on the rotor plane, which has an important impact on the aerodynamic load. The invention adopts a typical tower shadow model for analysis. Please refer to fig. 5, which is a schematic diagram illustrating the influence of the tower shadow effect on the wind speed. Considering blade deformation, the axial and tangential inflow velocities of the blade section Q under the effect of the tower shadow can be expressed as:
(41)
in the method, in the process of the invention,velocity for wind in axial direction;is the velocity of the wind in the horizontal direction;the speed of wind along the radial direction of the tower drum is given;the speed of wind along the tangential direction of the tower drum;the azimuth angle is the rotation angle of the blade; r is (r) t Is a tower having a height equal to the position of the blade element Cylinder radius. Wherein r is t The calculation formula is as follows:
(42)
in the method, in the process of the invention,the radii of the top and bottom of the tower, respectively, it is noted that in general the tower shadow effect only acts on the lower half of the impeller rotation plane, i.e. 90 degrees<θ wing <270 degree area.
The yaw effect is: due to the frequent changes in wind speed direction, even if the wind turbine has a yaw control system, the wind turbine rotor often operates out of alignment with the wind speed due to delays in the wind turbine yaw system. Referring to fig. 5, the yaw effect causes the wake to tilt. In addition, the induction speed is also caused to change at different azimuth angles. In the present invention, a yaw model describing an induced speed skew distribution proposed by Ge Laowo is employed:
(43)
wherein χ is wake deflection; θ 0 Is the reference position (theta) where the blade enters the wake deepest wing >At 0 is pi/2, theta wing <3 pi/2 at 0);is yaw-induced speed;the total length of the blade;the azimuth angle is the blade rotation.
In addition, the axial and tangential velocities also vary due to yaw inflow, which can be expressed as:
(44)
in the method, in the process of the invention,is the yaw angle.
The basic movement of the wind driven generator is as follows: the motion of a basic platform of a fan (wind driven generator) causes a fan impeller to move along with the fan foundation, so that the wind speed is linearly changed on the impeller, and the relative wind speeds of different parts on the impeller are changed. Referring to fig. 7, the offshore floating wind turbine foundation motion is considered to be a rigid motion with six degrees of freedom (heave, roll, heave, pitch). OXYZ is the global coordinates of the fan and OXYZ is the local rotor plane coordinates. Angle θ between y-axis and blade wing Is the azimuth angle of the blade; v (V) 0 For the inflow wind speed. Since the direction of movement of the fan wheel caused by the pitching and pitching movement of the fan foundation is opposite to the direction of wind speed, the aerodynamic force and thrust of the impeller can be obviously influenced by the two movements. The pitching motion of the fan foundation causes the fan impeller to deflect while rolling, and the influence on the aerodynamic characteristics of the fan is more complex. Thus, the present invention only analyzes the pitching motion of the fan foundation. In order to highlight the effect of the fan foundation movement, the influence of wind shear, tower shadow effect, wind turbulence and other unsteady effects are ignored here. The invention assumes the basic motion of the fan as a sine function:
(45)
in θ pp And t is the pitch angle and time of the fan foundation platform respectively; a is that pitch And f pitch Is the amplitude and frequency of the pitching of the platform. The initial phase angle phi in the present embodiment 0 Is set to 0, and for a real offshore floating fan, the amplitude and the frequency are respectively A in pitching pitch = 0° ~ 4° 、 f pitch =0 to 0.2 Hz, in this embodiment, f pitch Take 0.1 Hz.
In the basic movement of wind power generation, for any section Q on a blade, the relative speed caused by the basic movement of a fan basic platformExpressed as:
(46)
wherein V is p An additional wind speed caused by pitching movement of the base platform; h is the distance from the height of the hub to the pitching center of the foundation platform; omega pitch Angular velocity for pitch of the platform; y is Q Is the y-axis coordinate in the oxyz coordinate system of blade section Q, where, for a rigid blade, y Q Expressed as:
(47)
wherein f rotor The rotation frequency of the impeller is set; phi (phi) wing The initial phase angle of the blade is here set to zero, i.e. initially the blade 1 is just at the very top of the impeller.
However, for a real large fan blade, the flexibility of the blade is often larger, and errors are caused when the formula is adopted for calculation. Therefore, the method outputs the coordinates of each section of the blade after deformation directly through real-time calculation based on the proposed blade model. The frequency of the basic movement is very small and the frequency of the rotation of the gear box is far different, and meanwhile, the amplitude of the basic movement is smaller than 4 degrees, which is equivalent to the fan transmission system working in a quasi-static environment. Therefore, when analyzing, the basic excitation of the gear box by the basic movement of the fan is ignored, and only aerodynamic force changes caused by the basic movement of the fan are considered.
And S4, in a pneumatic-elastic coupling model building module, a pneumatic-elastic coupling model is built by adopting a loose coupling method, and data between the impeller model and the aerodynamic model are exchanged in each time step. Based on broad-meaning-And the time integration method is used for analyzing the aeroelastic characteristics and solving by adopting a self-adaptive time step strategy. The aeroelastic characteristic analysis is performed by selectively considering Blade gravity factors and control system effects. The process of building the aeroelastic coupling model is specifically explained below.
Referring to fig. 8, at each time step, the position, velocity and torsional deformation of each blade section are transferred from the impeller model (structural module) to the aerodynamic model (aerodynamic module). These data are the basis for calculating the relative wind speed and angle of attack for each blade section. In contrast, aerodynamic loads are exerted on the discrete blade beam units as externally distributed forces. In each calculation time step, keeping the aerodynamic load unchanged until the residual error of the impeller model solution in the time step converges, and adopting a generalized-A time integration method. It should be noted in particular that the present invention designs and employs an adaptive variable time step strategy to increase computational efficiency during the computation process. The adaptive time-step strategy is: when the iteration frequency limit is reached in a certain time step, the residual error requirement is not met, and the solving time step is reduced; when the residual convergence requirement is met and the adopted solving time step does not reach the specified maximum time step, the solving time step is further increased. The variable step length solving strategy can well improve the calculation efficiency on the premise of ensuring the solving precision.
Referring to fig. 9, a software embodiment of a system for predicting blade deformation and wheel hub load of an offshore floating wind turbine according to the present invention is specifically described below. The software mainly comprises five functional interfaces of welcome use, blade parameter input, wind speed parameter input, calculation parameter input and result output. Wherein, welcome the function interface to be used for appointing the working path, giving out the operation step, controlling the operation of the whole procedure; the blade parameter input functional interface is used for defining the blade structure attribute, the blade aerodynamic attribute and the airfoil aerodynamic performance; the wind speed parameter input function interface is used for defining the wind speed spatial distribution when the integral wind speed change (stable wind, gust, gradual change wind and turbulent wind) and the unstable influencing factors (wind shearing, tower shadow effect, yaw effect and matrix movement) act; meter with a meter bodyThe parameter calculation functional interface is mainly used for generalized-The time integration method is used for controlling, and meanwhile, whether the influence of a control system and the gravity factor of the blade are considered or not is selected; in the result output functional interface, the calculated data result is mainly stored for further joint analysis, and the result is directly displayed for individual blade aeroelastic characteristic analysis. In this example, the U.S. energy laboratory NREL 5-MW fan is used as the data source, which is a typical large three-blade windward megawatt fan with an impeller radius up to 63m and a tower height up to 90m. Therefore, the structural parameters of the wind turbine blade can effectively represent the large-scale and flexible characteristics of the modern large-scale offshore floating fan blade. Referring to fig. 10, the airfoil of the blade is of an elongated configuration. In order to achieve optimal blade aerodynamic distribution, the blade airfoils at different radial positions of the impeller are different. The main performance parameters and operating conditions of the NREL 5-MW fan are shown in table 1.
TABLE 1 Main Performance parameters and operating conditions of NREL 5-MW wind generators
Parameters (parameters) | Numerical value |
Rated power | 5 MW |
Fan direction, number of blades (N) b ) | Upwind 3 blades |
Impeller radius (R), hub radius (R hub ) Tower height (H) | 63 m, 1.5 m, 90 m |
Cut-in wind speed, rated wind speed and cut-out wind speed | 3 m/s, 11.4 m/s, 25 m/s |
Impeller cutting-in rotation speed and rated rotation speed of impeller | 6.9 r/min, 12.1r/min |
Hub cantilever length(s) 0 ) Blade pretilt angle (θ) cone ) | 5 m, 2.5º |
Tip and root tower radius (r top , r bottom ) | 1.935 m, 3 m |
Fan power transfer efficiency | 94.4% |
Hub mass, hub polar moment of inertia | 56780 kg, 115926 kg·m2 |
Referring to fig. 11, a welcome function interface is specifically illustrated below. After the software is opened, a welcome function interface is entered first. The execution button of the present interface, which is used to clear the last running result, is used to clear the previous running result and prepare for the running of the software. The execution button of "reset all parameters" resets all parameters in the present software to the initial parameters, so this button needs to be used carefully, otherwise the settings of the previous user will be all cleared, and the default will be restored. The "working path" is used to create a working path for the software, which is used to store temporary files generated during the running process of the software. The operation steps of the software are given at the left bottom of the page and are used for guiding the primary user to use the software, and mainly comprise: (1) inputting/importing blade performance parameters; (2) inputting wind speed parameters; (3) setting calculation parameters; (4) result display/output. The right side of the interface gives a schematic view of the aerodynamic coupling of a single blade, in which the blade generates lift, thrust and torsional moments under the action of wind, resulting in a buckling coupling deformation of the blade.
Referring to fig. 12 to 14, a specific explanation of the blade parameter input function interface is provided below, where the function interface includes a blade structure attribute, a blade aerodynamic attribute, and an airfoil aerodynamic performance interface. Referring to fig. 12, the blade structure attribute interface mainly defines blade structure parameters critical for modeling such as blade unit length, blade mass, flapping, shimmy direction, section stiffness, section torsional stiffness, axial tensile stiffness, etc. (other blade structure attributes may be displayed by a cross bar below the interface). These blade structure parameters are dynamic characterization parameters of the blade structure, and are mainly used for constructing geometrically accurate beam models. The software supports the direct import of the blade structure parameters from EXCLE, TXT and other types of files, and is convenient and quick. Meanwhile, the software is provided with buttons of adding, deleting and saving, supports manual input/modification of each blade structural parameter, can save the input/modified blade structural parameter, and is convenient for direct calling next time. Referring to FIG. 13, the blade aerodynamic parameters of the blade aerodynamic interface mainly include the twist angle of the blade, the chord length of the blade, and the airfoil of the blade. In the blade aerodynamic attribute interface, the left side is an input table of blade aerodynamic parameters, and the right image shows airfoil data to be input in the interface. The function key at the lower right corner has the same function as the key of the blade structure attribute function page, supports manual input/modification of pneumatic parameters of each blade, can store the input/modified pneumatic parameter data of the blade, and also supports data import. Referring to fig. 14, the aerodynamic parameters of the airfoil in the aerodynamic performance interface of the airfoil mainly include lift coefficient, drag coefficient and twist coefficient, which are determined by the characteristics of the airfoil itself and the angle of attack. The left side is a parameter main input area, the right side is a related parameter description, and a function key at the lower right corner supports manual input/modification of each airfoil aerodynamic parameter, can save the input/modified airfoil aerodynamic parameter data, and also supports data import.
Referring to fig. 15 to 22, the following specifically illustrates a wind speed parameter input function interface. Wind speed parameters are used to describe the magnitude of wind speed and its spatial distribution across the impeller, which affects the aerodynamic forces on the individual blades of the impeller. The software respectively carries out parameter setting on the whole wind speed and the space unstable wind speed distribution, and tries to cover the wind speed possibly suffered by the impeller in the actual running process. Setting a wind speed parameter setting interface and a wind speed space distribution setting interface.
Typically, wind speed magnitude parameters include substantially steady wind speed, gust wind speed, wind ramp speed, and turbulent wind speed. The software can respectively provide the four wind speeds for analyzing the aeroelastic characteristics of the fan blade. The setting interface of the stable wind speed is shown in fig. 15, and the corresponding stable wind speed value can be directly input, and whether the set wind speed change condition is correct can be observed through the display button. The setting interface of gust wind speed is shown in fig. 16, the gust wind speed variation form is defined by setting the basic wind speed, the wind speed increment, the starting time and the duration time, and a display window is provided on the right side of the interface and is used for checking parameter setting conditions. The setting interface of the gust wind speed of the software is shown in fig. 17, and the gradual change wind speed change form is defined by setting the basic wind speed, the maximum increment, the starting time, the peak time and the ending time, and a display window is provided on the right side of the interface and is used for checking parameter setting conditions. The interface for setting the wind speed of the gust of the software is shown in fig. 18, turbulent wind is generated by setting the average wind speed, the harmonic sampling number, the harmonic amplitude and the harmonic frequency, and a display window is provided on the right side of the interface and is used for checking parameter setting conditions. In the software, four wind speeds, namely, a stable wind speed, a gust wind speed, a gradual change wind speed and a turbulent wind speed, have mutual exclusivity. In the using process, only one wind speed type is required to be set, and meanwhile, data in other wind speeds are zeroed, so that software misjudgment is prevented. All wind speed parameters can be set to zero without additional setting by welcome the reset parameters of the function interface in the early stage of software operation.
The stable wind speed is the most basic wind speed, the wind speed does not change along with the time, corresponding stable wind speed numerical values can be directly input in the software, and whether the set wind speed change condition is correct can be observed through the display button. Gusts refer to sudden increases in wind speed, which are the result of air being disturbed, and can be described by the following formula:
(48)
in the method, in the process of the invention,andrespectively representing the time of beginning gust and the period of gust;is the wind speed increment caused by gust, is described by sine function, and in order to improve the calculation efficiencyAssuming a constant;is gust wind speed; t is time; t is the conversion period.
The gradual change wind refers to a phenomenon that the wind speed gradually changes with time, and can be described by the following equation:
(49)
in the method, in the process of the invention,、 、 the starting time, the peak value arrival time and the gradual change ending time of the gradual change of the wind speed are respectively;is the maximum wind value of the gradual change wind;is a gradual change wind speed;is the base wind speed.
Turbulent wind is the most common wind speed variation, in which the wind speed varies in a random continuous fashion, described by the sum of harmonics:
(50)
wherein V is 0 、N、A i And omega i The average wind speed, the harmonic sampling number, the harmonic amplitude and the harmonic frequency are respectively. The software calculates the harmonic amplitude by adopting a Dryden frequency spectrum, and the value range of the harmonic frequency is 0.1-10 Hz. After calculation, a small frequency corresponds to a larger amplitude, and a larger frequency corresponds to a smaller fluctuation amplitude. Therefore, the wind speed simulated in this way has a high degree of reliability.
The air flow space distribution of the whole impeller is related to the atmosphere and the structural characteristics of the fan, and mainly shows wind shearing effect, tower shadow effect, yaw effect, basic movement of the fan and the like. These unsteady flows can cause maldistribution of gas flow rates across the wheel space, resulting in the hub being subjected to periodic non-torque loads. The safety and reliability of the internal structure of the fan are seriously affected by the non-torque load, the design service life of equipment is shortened, and the failure rate of the equipment is increased. The wind speed spatial distribution parameters are thus derived from wind shear effects, tower shadow effects, yaw effects and wind turbine base motions.
The wind shearing effect reflects the change of wind speed along with the vertical height, so that the wind speed on the fan blade is unevenly distributed, and especially the difference between the upper wind speed and the lower wind speed of the impeller is maximum. As the rotor diameter increases, the wind shear effect becomes more pronounced. Although the sea friction is small relative to land, the shearing effect is reduced. On the other hand, the impeller of the offshore wind turbine tends to be large, so that the shearing effect of the whole impeller is large. The setting interface of the software wind shear is shown in fig. 19, and control parameters of a switching key and a wind shear factor of whether the wind shear influence is considered are given on the left side of the interface: wind speed at hub height, hub height and wind shear empirical index.
The tower shadow effect is: for upwind wind generators, the incoming wind speed will change due to the presence of the tower. It introduces periodic unsteady flows on the rotor plane, which has an important impact on the aerodynamic load. The tower shadow effect cannot be avoided, and the influence of the tower shadow effect can be reduced as much as possible only through optimization of parameters. The software tower shadow effect parameter setting interface is shown in fig. 20, and on-off key of whether the tower shadow effect is considered and control parameters of the tower shadow effect are given on the left side of the interface: the tower top radius, the tower bottom radius, the hub height and the suspension distance.
The yaw effect is: due to frequent changes in wind speed direction, even if the wind turbine has a yaw control system, the wind turbine is often operated out of alignment with the wind speed due to delays in the yaw system of the wind turbine. Yaw effects can cause the wake to tilt and cause the induced speed to change at different azimuth angles. The yaw effect parameter setting interface of the software is shown in fig. 21, and on-off key of whether the yaw effect is considered and control parameters of the yaw effect are given on the left side of the interface: yaw angle. The yaw angle may be negative, indicating that the direction of yaw is the same.
The basic movement of the wind driven generator is as follows: the floating fan on sea is for the fan on land, and its basis is in continuous motion, and the basis motion of platform can produce the influence to the relative inflow speed of fan blade, along with the fan is big, and the height of tower section of thick bamboo is in continuous extension too for matching impeller diameter, has finally led to the basis to swing the motion and has enlarged the influence to the impeller is whole through tower section of thick bamboo length, especially to the pitching motion of platform. The software basic parameter setting interface is shown at 22, and control parameters of whether the on-off key and basic movement of the livestock are considered are given on the left side of the interface. The software considers only two basic movement forms (pitch and heave) that have the greatest influence on the relative wind speed, and the control parameters of the basic movement include the wind speed at hub height, hub height and pitch/heave control parameters. Only pitch control parameters (amplitude and frequency) are shown in fig. 22, which are similar to pitch as the amplitude and frequency of the pitch motion.
Referring to fig. 23, the calculation parameter input function interface will be specifically explained. The calculation parameters mainly control a generalized-alpha time integration method in the calculation process, wherein the control parameters are total simulation time, maximum solving step length, iteration precision and maximum iteration times. The right side of the interface is provided with internal calculation logic of the aeroelastic coupling, and simultaneously, the interface is provided with two switch buttons for selecting whether to consider the gravity factors of the blades and the influence of a control system. Neglecting the blade gravity can better observe the action mechanism of different influencing factors on the blade aeroelastic characteristics, so the user can conveniently analyze the blade aeroelastic characteristics by setting a switch key.
Referring to fig. 24 to 28, the result output function interface is specifically explained below. The result output functional page comprises a data storage interface and a result display interface, and is mainly used for storing the calculated data result for further deep analysis and directly displaying the result for analyzing the aeroelastic characteristics of the independent blade.
As shown in fig. 24, the software classifies and stores data commonly used for blade aeroelastic characteristic analysis, so as to meet the personalized requirements of users. The stored data mainly comprises: (1) Aerodynamic attribute results-tangential relative wind speed, axial relative wind speed, and angle of attack. Wherein the tangential relative wind speed and the axial relative wind speed are v rel The tangential and axial components directly affect the final angle of attack. The software outputs the change rule of the section of the blade along with time, can accurately reveal the wind speed distribution rule of the impeller, and provides the most original data for the aeroelastic characteristic analysis of the blade. (2) Blade load results-tangential load, axial load and twistAnd (5) turning the load. The tangential load, the axial load and the torsional load are pneumatic loads borne by the blades, are key input parameters of load analysis of single blades, and facilitate later stress deformation analysis of the blades. (3) Blade deformation results-shimmy deformation, flapping deformation and torsional deformation. The deformation of the blade under the action of pneumatic load is conveniently evaluated as a result of the shimmy deformation, the flapping deformation and the torsional deformation, and whether the deformation of the blade meets the requirements or not is judged, and whether the deformation of the blade meets the distance of the blade deviating from the tower barrel or not is judged. (4) Impeller load results-thrust, horizontal force, vertical force, overturning moment, yaw moment and torque. The invention adopts a data storage form to facilitate analysis and research rules: all data used for analysis can be derived in the data storage mode, and research designers can conduct comparison analysis on different data according to own requirements. In other embodiments, the relevant graphic comparison module can be formulated for comparison analysis, which is convenient for use.
Referring to fig. 25 to 28, the results show that the interface is mainly used for users who independently use the software or who perform basic debugging analysis, and can directly check aerodynamic properties, blade loads, blade deformations and impeller loads, and primarily understand basic aeroelastic properties of the analyzed blades. The aerodynamic properties, blade loading, blade deformation, and impeller loading results are shown in fig. 25, 26, 27, 28, respectively. The beginning and end of the result presentation time segment may be set at each result presentation interface. As the two-dimensional image display result is adopted, only one axial position of the blade can be selected for display at a time. The results of the first level under aerodynamic properties, blade loading, blade deformation and impeller loading can be selected autonomously by the hooking function, and different types of results can be displayed simultaneously. In other embodiments, corresponding three-dimensional drawing modules and animation modules are developed, so that designers can be helped to observe the calculation results and analysis rules more clearly and intuitively.
Claims (10)
1. The method for predicting the deformation of the blade and the hub load of the impeller of the offshore floating fan is characterized by comprising the following steps of:
s1, acquiring initial parameters, wherein the initial parameters comprise blade parameters and wind speed parameters;
S2, processing the blade parameters based on a geometric precise beam theory considering geometric large deformation, and constructing an impeller model, wherein the impeller model comprises a plurality of blade models; the impeller model corresponds to one impeller, the blade model corresponds to one blade, and each impeller is provided with a plurality of blades;
s3, processing the wind speed parameter based on an unsteady phyllanthus metal momentum theory, and establishing a aerodynamic model, wherein the aerodynamic model is used for calculating aerodynamic load;
s4, constructing a pneumatic-elastic coupling model by combining the impeller model and the aerodynamic model, applying the pneumatic load to the impeller, analyzing pneumatic-elastic characteristics, and outputting analysis results, wherein the analysis results comprise blade load, blade deformation and impeller load.
2. The method for predicting the deformation and impeller hub load of an offshore floating wind turbine blade of claim 1, wherein the geometrically precise beam model comprises a plurality of geometrically precise beam units and beam nodes, each blade is discretized into a plurality of geometrically precise beam units, and the plurality of blades are connected by the beam nodes; the geometric precise beam unit is constructed by adopting three nodes, and the formula is as follows:
;
in the method, in the process of the invention, Is a unit-shaped function; />Is a third-order unit vector; />Is a node 1 shape function; />Is a node 2-shaped function; />Is a node 3-shaped function.
3. The method for predicting the deformation of a blade and the hub load of an impeller of a floating offshore wind turbine according to claim 2, wherein in step S2, a geometrically accurate beam model is established to simulate the deformation of the blade, and the geometrically accurate beam model is formulated as follows:
;
in the method, in the process of the invention,is a beam unit mass matrix; />The acceleration matrix is a beam node acceleration matrix; />The gyro matrix is a beam unit gyro matrix; />Is a beam node velocity matrix; />Is a beam unit stiffness matrix; />The displacement matrix is a beam node displacement matrix; />The gravity matrix is born by the beam unit;an external force matrix of the beam unit; />Is a beam cell internal force matrix.
4. The method for predicting the deformation and the hub load of an impeller of a floating offshore wind turbine according to claim 1, wherein in the unsteady principle of momentum, the wind speed on the blade is unevenly distributed due to unsteady factors.
5. The method for predicting the deformation and impeller hub loading of the offshore floating wind turbine blade according to claim 4, wherein the wind speed parameters comprise a wind speed size parameter and a wind speed space distribution parameter, and the wind speed is divided into a basic stable wind speed, a gust wind speed, a gradual change wind speed and a turbulent wind speed according to the wind speed size parameter; dividing the unsteady factors into wind shearing effect, tower shadow effect, yaw effect and basic movement of the wind driven generator according to the wind speed space distribution parameters.
6. The method of predicting the blade deformation and the wheel hub load of an offshore floating wind turbine of claim 1, wherein in step S3, calculating the aerodynamic load comprises the following formula:
;
in the method, in the process of the invention,representing air density; />Representing the cross-sectional chord length of the blade; />A drag coefficient representing a blade section; />A lift coefficient representing a blade section; />A twist coefficient representing a cross section of the blade; />Is corresponding axial aerodynamic force; />Is a corresponding tangential aerodynamic force; />Torque caused by the offset of the aerodynamic center relative to the shear center; />Lift per unit blade length; />Representing the resistance per unit blade length; />Is a relative velocity vector; />The final real wind speed and the included angle of the rotating plane are obtained.
7. The method of predicting the deformation of a floating offshore wind turbine blade and the hub load of an impeller of claim 1, wherein a loose coupling method is used to build a pneumatic-elastic coupling model, and data between the impeller model and the aerodynamic model are exchanged at each time step.
8. The method for predicting the blade deformation and the hub load of an offshore floating wind turbine of claim 1, wherein in step S4,based on broad-meaning-And the time integration method is used for analyzing the aeroelastic characteristics and solving by adopting a self-adaptive time step strategy.
9. The method for predicting blade deformation and impeller hub load of an offshore floating wind turbine of claim 1, wherein the blade gravity factor and control system effects are considered when performing the aeroelastic characteristic analysis.
10. An offshore floating wind turbine blade deformation and impeller hub load prediction system, comprising:
the data acquisition module is used for acquiring initial parameters, wherein the initial parameters comprise blade parameters and wind speed parameters;
the impeller model building module is used for processing the blade parameters based on a geometric precise beam theory considering geometric large deformation to build an impeller model, and the impeller model comprises a plurality of blade models; the impeller model corresponds to one impeller, the blade model corresponds to one blade, and each impeller is provided with a plurality of blades;
the aerodynamic model building module is used for processing the wind speed parameter based on the unsteady leaf element momentum theory and building an aerodynamic model which is used for calculating aerodynamic load;
and the aeroelastic analysis module is used for combining the impeller model and the aerodynamic model to construct an aeroelastic coupling model, applying the aerodynamic load to the impeller, carrying out aeroelastic characteristic analysis, and outputting analysis results, wherein the analysis results comprise blade load, blade deformation and impeller load.
Priority Applications (1)
Application Number | Priority Date | Filing Date | Title |
---|---|---|---|
CN202410182336.8A CN117744409A (en) | 2024-02-19 | 2024-02-19 | Method and system for predicting blade deformation and impeller hub load of offshore floating fan |
Applications Claiming Priority (1)
Application Number | Priority Date | Filing Date | Title |
---|---|---|---|
CN202410182336.8A CN117744409A (en) | 2024-02-19 | 2024-02-19 | Method and system for predicting blade deformation and impeller hub load of offshore floating fan |
Publications (1)
Publication Number | Publication Date |
---|---|
CN117744409A true CN117744409A (en) | 2024-03-22 |
Family
ID=90277657
Family Applications (1)
Application Number | Title | Priority Date | Filing Date |
---|---|---|---|
CN202410182336.8A Pending CN117744409A (en) | 2024-02-19 | 2024-02-19 | Method and system for predicting blade deformation and impeller hub load of offshore floating fan |
Country Status (1)
Country | Link |
---|---|
CN (1) | CN117744409A (en) |
Cited By (1)
Publication number | Priority date | Publication date | Assignee | Title |
---|---|---|---|---|
CN117973157A (en) * | 2024-04-01 | 2024-05-03 | 上海勘测设计研究院有限公司 | Simulation method for action of composite cylinder foundation and soil body, storage medium and electronic equipment |
Citations (5)
Publication number | Priority date | Publication date | Assignee | Title |
---|---|---|---|---|
CN101882177A (en) * | 2010-06-18 | 2010-11-10 | 北京航空航天大学 | Aeroelastic stability fluid-structure interaction prediction method of turbo-machine changed interblade phase angles |
CN102708266A (en) * | 2012-06-12 | 2012-10-03 | 中国科学院工程热物理研究所 | Method for predicting and calculating limit load of horizontal-axis wind turbine blade |
CN106768917A (en) * | 2016-11-23 | 2017-05-31 | 中国科学院工程热物理研究所 | A kind of pneumatic equipment bladess scene load test and appraisal procedure |
CN115544883A (en) * | 2022-10-08 | 2022-12-30 | 浙江大学 | Online measurement method and system for load and platform deformation of floating type offshore wind turbine generator |
CN117436322A (en) * | 2023-12-21 | 2024-01-23 | 浙江远算科技有限公司 | Wind turbine blade aeroelastic simulation method and medium based on phyllin theory |
-
2024
- 2024-02-19 CN CN202410182336.8A patent/CN117744409A/en active Pending
Patent Citations (5)
Publication number | Priority date | Publication date | Assignee | Title |
---|---|---|---|---|
CN101882177A (en) * | 2010-06-18 | 2010-11-10 | 北京航空航天大学 | Aeroelastic stability fluid-structure interaction prediction method of turbo-machine changed interblade phase angles |
CN102708266A (en) * | 2012-06-12 | 2012-10-03 | 中国科学院工程热物理研究所 | Method for predicting and calculating limit load of horizontal-axis wind turbine blade |
CN106768917A (en) * | 2016-11-23 | 2017-05-31 | 中国科学院工程热物理研究所 | A kind of pneumatic equipment bladess scene load test and appraisal procedure |
CN115544883A (en) * | 2022-10-08 | 2022-12-30 | 浙江大学 | Online measurement method and system for load and platform deformation of floating type offshore wind turbine generator |
CN117436322A (en) * | 2023-12-21 | 2024-01-23 | 浙江远算科技有限公司 | Wind turbine blade aeroelastic simulation method and medium based on phyllin theory |
Non-Patent Citations (8)
Title |
---|
ZHANWEI LI 等: "Aerodynamic and aeroelastic characteristics of flexible wind turbine blades under periodic unsteady inflows", JOURNAL OF WIND ENGINEERING & INDUSTRIAL AERODYNAMICS, vol. 197, 16 December 2019 (2019-12-16), pages 104057 * |
ZHANWEI LI 等: "Dynamic Performance Evaluation of Variable‑Speed Wind Turbine Transmission System in Non‑Stationary Conditions", JOURNAL OF VIBRATION ENGINEERING & TECHNOLOGIES, 26 June 2023 (2023-06-26), pages 1 - 21 * |
中国海洋工程学会: "《第十九届中国海洋(岸)工程学术讨论会会论文集 上》", 30 September 2019, 海洋出版社, pages: 249 * |
刘利琴 等: "浮式垂直轴风机的动力学建模、仿真与实验研究", 力学学报, vol. 49, no. 02, 15 December 2016 (2016-12-15), pages 299 - 307 * |
刘格梁;胡志强;段斐;: "海上浮式风机在支撑平台运动影响下的气动特性研究", 海洋工程, no. 01, 30 January 2017 (2017-01-30), pages 42 - 50 * |
吴有生 等: "《第二十七届全国水动力学研讨会文集 下》", 31 October 2015, 海洋出版社, pages: 952 * |
姜劲 等: "半潜型垂直轴风机纵荡下气动载荷特性", 哈尔滨工程大学学报, vol. 40, no. 06, 18 December 2018 (2018-12-18), pages 1084 - 1089 * |
徐浩然;杨华;朱卫军;刘超;: "偏航工况下水平轴风力机叶片气动变形数值模拟研究", 可再生能源, no. 05, 20 May 2013 (2013-05-20), pages 51 - 55 * |
Cited By (2)
Publication number | Priority date | Publication date | Assignee | Title |
---|---|---|---|---|
CN117973157A (en) * | 2024-04-01 | 2024-05-03 | 上海勘测设计研究院有限公司 | Simulation method for action of composite cylinder foundation and soil body, storage medium and electronic equipment |
CN117973157B (en) * | 2024-04-01 | 2024-06-04 | 上海勘测设计研究院有限公司 | Simulation method for action of composite cylinder foundation and soil body, storage medium and electronic equipment |
Similar Documents
Publication | Publication Date | Title |
---|---|---|
Ponta et al. | Effects of rotor deformation in wind-turbine performance: the dynamic rotor deformation blade element momentum model (DRD–BEM) | |
Chen et al. | Numerical analysis of unsteady aerodynamic performance of floating offshore wind turbine under platform surge and pitch motions | |
Liu et al. | Aeroelastic analysis of a floating offshore wind turbine in platform‐induced surge motion using a fully coupled CFD‐MBD method | |
Dai et al. | Analysis of wind turbine blades aeroelastic performance under yaw conditions | |
Gebhardt et al. | Non-linear aeroelasticity: an approach to compute the response of three-blade large-scale horizontal-axis wind turbines | |
Rezaei et al. | Development of a reduced order model for nonlinear analysis of the wind turbine blade dynamics | |
CN117744409A (en) | Method and system for predicting blade deformation and impeller hub load of offshore floating fan | |
Tang et al. | Study on the aeroelastic responses of a wind turbine using a coupled multibody-FVW method | |
Jokar et al. | Dynamic modeling and free vibration analysis of horizontal axis wind turbine blades in the flap-wise direction | |
Lupton | Frequency-domain modelling of floating wind turbines | |
Chen et al. | Study on composite bend-twist coupled wind turbine blade for passive load mitigation | |
Matha et al. | Efficient critical design load case identification for floating offshore wind turbines with a reduced nonlinear model | |
Liu et al. | Development of a fully coupled aero-hydro-mooring-elastic tool for floating offshore wind turbines | |
Ma et al. | Analysis of the coupled aeroelastic wake behavior of wind turbine | |
Andersen et al. | Load alleviation on wind turbine blades using variable airfoil geometry | |
CN106777499B (en) | Complete machine dynamic modeling method for double-fed asynchronous wind generating set | |
CN115525988B (en) | Wind turbine autonomous load simulation calculation and correction system | |
Sessarego et al. | Development of an aeroelastic code based on three‐dimensional viscous–inviscid method for wind turbine computations | |
Vitsas et al. | Multiscale aeroelastic simulations of large wind farms in the atmospheric boundary layer | |
Wang et al. | A study of nonlinear aeroelastic response of a long flexible blade for the horizontal axis wind turbine | |
Qi et al. | Data-driven control of the coupled loads for floating wind turbine systems based on deformable trailing edge flaps | |
Jeong et al. | Numerical investigation of optimal yaw misalignment and collective pitch angle for load imbalance reduction of rigid and flexible HAWT blades under sheared inflow | |
El-Baklish et al. | Nonlinear model predictive pitch control of aero-elastic wind turbine blades | |
Santoni et al. | Coupling turbulent flow with blade aeroelastics and control modules in large-eddy simulation of utility-scale wind turbines | |
Shen et al. | Fast online reinforcement learning control of small lift-driven vertical axis wind turbines with an active programmable four bar linkage mechanism |
Legal Events
Date | Code | Title | Description |
---|---|---|---|
PB01 | Publication | ||
PB01 | Publication | ||
SE01 | Entry into force of request for substantive examination | ||
SE01 | Entry into force of request for substantive examination |