CN113569360B - 一种风力机叶片抗颤振翼型簇设计方法 - Google Patents
一种风力机叶片抗颤振翼型簇设计方法 Download PDFInfo
- Publication number
- CN113569360B CN113569360B CN202110962098.9A CN202110962098A CN113569360B CN 113569360 B CN113569360 B CN 113569360B CN 202110962098 A CN202110962098 A CN 202110962098A CN 113569360 B CN113569360 B CN 113569360B
- Authority
- CN
- China
- Prior art keywords
- airfoil
- flutter
- inertia
- profile
- wind turbine
- 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.)
- Active
Links
- 238000000034 method Methods 0.000 title claims abstract description 39
- 230000002068 genetic effect Effects 0.000 claims abstract description 9
- 238000004422 calculation algorithm Methods 0.000 claims abstract description 8
- 241001264313 Semioptera wallacii Species 0.000 claims description 3
- 238000010206 sensitivity analysis Methods 0.000 claims description 3
- 238000000556 factor analysis Methods 0.000 claims description 2
- 238000004364 calculation method Methods 0.000 abstract description 25
- 238000013461 design Methods 0.000 abstract description 23
- 238000005457 optimization Methods 0.000 abstract description 12
- 239000012530 fluid Substances 0.000 abstract description 5
- 238000011160 research Methods 0.000 abstract description 2
- 239000010410 layer Substances 0.000 description 8
- 238000010586 diagram Methods 0.000 description 5
- 230000008569 process Effects 0.000 description 4
- 238000004458 analytical method Methods 0.000 description 3
- 238000012804 iterative process Methods 0.000 description 3
- 230000005484 gravity Effects 0.000 description 2
- 230000010354 integration Effects 0.000 description 2
- 238000012805 post-processing Methods 0.000 description 2
- 230000002829 reductive effect Effects 0.000 description 2
- 238000004088 simulation Methods 0.000 description 2
- 239000012790 adhesive layer Substances 0.000 description 1
- 230000009286 beneficial effect Effects 0.000 description 1
- 230000005540 biological transmission Effects 0.000 description 1
- 230000008859 change Effects 0.000 description 1
- 239000002131 composite material Substances 0.000 description 1
- 238000007405 data analysis Methods 0.000 description 1
- 230000007423 decrease Effects 0.000 description 1
- 230000003111 delayed effect Effects 0.000 description 1
- 238000012938 design process Methods 0.000 description 1
- 230000000694 effects Effects 0.000 description 1
- 238000005516 engineering process Methods 0.000 description 1
- 230000008676 import Effects 0.000 description 1
- 230000006872 improvement Effects 0.000 description 1
- 230000002401 inhibitory effect Effects 0.000 description 1
- 230000000670 limiting effect Effects 0.000 description 1
- 239000000463 material Substances 0.000 description 1
- 238000013178 mathematical model Methods 0.000 description 1
- 238000012986 modification Methods 0.000 description 1
- 230000004048 modification Effects 0.000 description 1
- 230000036961 partial effect Effects 0.000 description 1
- 238000010248 power generation Methods 0.000 description 1
- 238000000926 separation method Methods 0.000 description 1
- 238000010008 shearing Methods 0.000 description 1
- 238000012546 transfer Methods 0.000 description 1
- 230000007704 transition Effects 0.000 description 1
- 230000000007 visual effect Effects 0.000 description 1
Classifications
-
- G—PHYSICS
- G06—COMPUTING; CALCULATING OR COUNTING
- G06F—ELECTRIC DIGITAL DATA PROCESSING
- G06F30/00—Computer-aided design [CAD]
- G06F30/10—Geometric CAD
- G06F30/17—Mechanical parametric or variational design
-
- F—MECHANICAL ENGINEERING; LIGHTING; HEATING; WEAPONS; BLASTING
- F03—MACHINES OR ENGINES FOR LIQUIDS; WIND, SPRING, OR WEIGHT MOTORS; PRODUCING MECHANICAL POWER OR A REACTIVE PROPULSIVE THRUST, NOT OTHERWISE PROVIDED FOR
- F03D—WIND MOTORS
- F03D1/00—Wind motors with rotation axis substantially parallel to the air flow entering the rotor
- F03D1/06—Rotors
- F03D1/065—Rotors characterised by their construction elements
- F03D1/0675—Rotors characterised by their construction elements of the blades
-
- G—PHYSICS
- G06—COMPUTING; CALCULATING OR COUNTING
- G06F—ELECTRIC DIGITAL DATA PROCESSING
- G06F30/00—Computer-aided design [CAD]
- G06F30/20—Design optimisation, verification or simulation
- G06F30/27—Design optimisation, verification or simulation using machine learning, e.g. artificial intelligence, neural networks, support vector machines [SVM] or training a model
-
- G—PHYSICS
- G06—COMPUTING; CALCULATING OR COUNTING
- G06F—ELECTRIC DIGITAL DATA PROCESSING
- G06F30/00—Computer-aided design [CAD]
- G06F30/20—Design optimisation, verification or simulation
- G06F30/28—Design optimisation, verification or simulation using fluid dynamics, e.g. using Navier-Stokes equations or computational fluid dynamics [CFD]
-
- G—PHYSICS
- G06—COMPUTING; CALCULATING OR COUNTING
- G06N—COMPUTING ARRANGEMENTS BASED ON SPECIFIC COMPUTATIONAL MODELS
- G06N3/00—Computing arrangements based on biological models
- G06N3/12—Computing arrangements based on biological models using genetic models
- G06N3/126—Evolutionary algorithms, e.g. genetic algorithms or genetic programming
-
- G—PHYSICS
- G06—COMPUTING; CALCULATING OR COUNTING
- G06F—ELECTRIC DIGITAL DATA PROCESSING
- G06F2111/00—Details relating to CAD techniques
- G06F2111/04—Constraint-based CAD
-
- G—PHYSICS
- G06—COMPUTING; CALCULATING OR COUNTING
- G06F—ELECTRIC DIGITAL DATA PROCESSING
- G06F2111/00—Details relating to CAD techniques
- G06F2111/10—Numerical modelling
-
- G—PHYSICS
- G06—COMPUTING; CALCULATING OR COUNTING
- G06F—ELECTRIC DIGITAL DATA PROCESSING
- G06F2113/00—Details relating to the application field
- G06F2113/08—Fluids
-
- G—PHYSICS
- G06—COMPUTING; CALCULATING OR COUNTING
- G06F—ELECTRIC DIGITAL DATA PROCESSING
- G06F2119/00—Details relating to the type or aim of the analysis or the optimisation
- G06F2119/14—Force analysis or force optimisation, e.g. static or dynamic forces
-
- Y—GENERAL TAGGING OF NEW TECHNOLOGICAL DEVELOPMENTS; GENERAL TAGGING OF CROSS-SECTIONAL TECHNOLOGIES SPANNING OVER SEVERAL SECTIONS OF THE IPC; TECHNICAL SUBJECTS COVERED BY FORMER USPC CROSS-REFERENCE ART COLLECTIONS [XRACs] AND DIGESTS
- Y02—TECHNOLOGIES OR APPLICATIONS FOR MITIGATION OR ADAPTATION AGAINST CLIMATE CHANGE
- Y02E—REDUCTION OF GREENHOUSE GAS [GHG] EMISSIONS, RELATED TO ENERGY GENERATION, TRANSMISSION OR DISTRIBUTION
- Y02E10/00—Energy generation through renewable energy sources
- Y02E10/70—Wind energy
- Y02E10/72—Wind turbines with rotation axis in wind direction
Landscapes
- Engineering & Computer Science (AREA)
- Physics & Mathematics (AREA)
- Theoretical Computer Science (AREA)
- General Physics & Mathematics (AREA)
- Evolutionary Computation (AREA)
- Geometry (AREA)
- General Engineering & Computer Science (AREA)
- Life Sciences & Earth Sciences (AREA)
- Health & Medical Sciences (AREA)
- Biophysics (AREA)
- Computer Hardware Design (AREA)
- Bioinformatics & Cheminformatics (AREA)
- Pure & Applied Mathematics (AREA)
- Mathematical Optimization (AREA)
- Artificial Intelligence (AREA)
- Evolutionary Biology (AREA)
- Software Systems (AREA)
- Computing Systems (AREA)
- Bioinformatics & Computational Biology (AREA)
- Mathematical Physics (AREA)
- Mathematical Analysis (AREA)
- Algebra (AREA)
- Data Mining & Analysis (AREA)
- Fluid Mechanics (AREA)
- Medical Informatics (AREA)
- Computer Vision & Pattern Recognition (AREA)
- Physiology (AREA)
- Genetics & Genomics (AREA)
- Biomedical Technology (AREA)
- Computational Linguistics (AREA)
- Computational Mathematics (AREA)
- General Health & Medical Sciences (AREA)
- Molecular Biology (AREA)
- Sustainable Development (AREA)
- Sustainable Energy (AREA)
- Chemical & Material Sciences (AREA)
- Combustion & Propulsion (AREA)
- Mechanical Engineering (AREA)
- Wind Motors (AREA)
- Structures Of Non-Positive Displacement Pumps (AREA)
Abstract
本发明公开了一种风力机叶片抗颤振翼型簇设计方法，属于风力机叶片翼型设计技术领域，包括以下步骤：S1：翼型表达；S2：颤振临界速度因素分析；S3：选择约束条件；S4：确定翼型目标函数。本发明通过在原有的翼型基础上添加翼型扰动函数来修改翼型的几何外形，提高其抗颤振稳定性能并增强其气动性能(升阻比更大)；通过经典抗颤振理论的研究，得出翼型的极惯性矩对其抗颤振性能具有突出的影响，其值可通过高斯‑勒让德积分计算得出；并将最大升阻比和极惯性矩作为优化目标，在isight平台中融合翼型扰动函数、高斯‑勒让德积分计算、网格划分、流体力学计算、多目标遗传算法来获取优化后的新翼型簇,其抗颤振性能及气动性能均得到了显著的提高。
Description
技术领域
本发明涉及属于风力机叶片翼型设计技术领域，具体涉及一种风力机叶片抗颤振翼型簇设计方法。
背景技术
叶片主要是由复合材料制成的薄壁空心结构，随着叶片的长度增加，其刚度降低，柔性特征突出，易形成气动弹性不稳定性，即颤振。颤振会使叶片发生疲劳破坏，产生裂纹，严重时会导致叶片直接断裂。为了保证大型风力发电设备的正常工作，必须重视风力机叶片的抗颤振性能。
风力机叶片的气动外形作为捕获风能的关键因素，通过对风机叶片翼型的优化设计，提高其颤振临界速度，是改善抗颤振性能的有效方法；风能作为一种清洁、取之不尽的能源在世界中被广泛应用，伴随着风力发电技术不断升级，提高风力机叶片的抗颤振性能成为叶片设计的重点。为此，提出一种风力机叶片抗颤振翼型簇设计方法。
发明内容
本发明所要解决的技术问题在于：如何通过对风机叶片翼型的优化设计，提高其颤振临界速度，改善抗颤振性能，提供了一种风力机叶片抗颤振翼型簇设计方法。
本发明是通过以下技术方案解决上述技术问题的，本发明包括以下步骤：
S1：翼型表达
对原有翼型添加形状扰动函数来表达新的翼型；
S2：颤振临界速度影响因素分析
对风机叶片进行敏感性分析，得到颤振临界速度影响因素，提高颤振临界影响速度因素中的扭转刚度以提高抗颤振性能；
S3：选择约束条件
确定翼型的约束条件为厚度约束条件与升阻比约束条件；
S4：确定翼型目标函数
将翼型在攻角α0下的升阻比与翼型的极惯性矩I作为目标函数，通过采用多目标遗传算法对风力机叶片翼型进行优化设计。
更进一步地，在所述步骤S1中，所述原有翼型为荷兰DU翼型系列中的DU35。
更进一步地，在所述步骤S1中，在原有翼型基础上添加形状扰动函数，保证X轴坐标不变，得到新翼型的型值点坐标：
上翼面翼型的型值点坐标：
yu(i)＝you(i)+Δyu(i)
xu(i)＝xou(i)
下翼面翼型的型值点坐标：
yl(i)＝yol(i)+Δyl(i)
xl(i)＝xol(i)
其中，yu(i)代表新翼型上翼面第i个点的纵坐标，you(i)代表原翼型上翼面第i个点的纵坐标，Δyu(i)代表上翼面第i个点的形状扰动函数，yl(i)代表新翼型下翼面第i点的纵坐标，yol(i)代表原始翼型上翼面第i个点的横坐标，Δyl(i)代表下翼面第i个点的形状扰动函数；xu(i)代表新翼型上翼面第i个点的横坐标，xou(i)代表原始翼型上翼面第i个点的横坐标，xol(i)代表原始翼型下翼面第i个点的横坐标。
更进一步地，上、下翼面形状扰动函数如下：
上翼面形状扰动函数：
下翼面形状扰动函数：
其中，下标u和l分别代表上翼面和下翼面，i代表翼型型值点索引，k与f(k)分别是控制翼型厚度与弯度分布关键点变量的个数与系数，g为给点矢量，f(k)代表Hicks-Henne型函数，ξ和η为指数因子。
更进一步地，在所述步骤S2中，翼型的扭转刚度表达式如下：
其中，G为剪切模量，R为翼型所在界面距离叶根的半径，I为翼型极惯性矩。
更进一步地，所述极惯性矩I的计算过程如下：
S21：计算出翼型的面积和质心坐标点：
翼型的面积：
其中，S为翼型的面积；b为翼型弦长，单位标准翼型值弦长取1；u代表横坐标，v代表纵坐标；θ1(u)表示翼型下翼面曲线函数；θ2(u)为上翼面曲线函数；
翼型质心坐标：
其中，ρ(u，v)代表不同坐标点处的密度，假设叶片是均匀的，则ρ(u，v)为常数；
则上式可进一步推导：
S22：坐标原点位于翼型的前缘点，求出翼型对x，y方向的惯性矩：
S23：平移坐标系，使原点与翼型质心重合，通过平行轴定理，求出翼型相对质心轴系的惯性矩：
I1＝Ix-Ay2
I2＝Iy-Ax2
其中：I1代表Y轴相对于质心的惯性矩；
A代表翼型的面积；
y代表Y轴与质心轴的垂直距离；
I2代表X轴相对于质心的惯性矩；
x代表X轴与质心的垂直距离；
则翼型极惯性矩为：
I＝I1+I2。
更进一步地，在所述步骤S3中，厚度约束条件为t1≤t≤t2，式中t为翼型厚度，t1、t2为设定值。
更进一步地，在所述步骤S3中，优化后的翼型在不同攻角下的升阻比不低于原翼型的升阻比，选择[α1，α2，…，αi]攻角下的升阻比作为约束条件，升阻比约束条件为：
cl′/cd′≥cl/cd
式中cl'为优化后的升力系数，cd'为优化后的阻力系数，cl为原有翼型的升力系数，cd为原有翼型的阻力系数。
更进一步地，在所述步骤S4中，目标函数为：
maximize＝[f1(x)，f2(x)]
f1(x)＝I
f2(x)＝cl/cd。
更进一步地，在所述步骤S4中，多目标遗传算法为NSGA-II。
本发明相比现有技术具有以下优点：
1)通过采用三次样条拟合翼型的型值点，得到了上下翼面的曲线函数，将翼型分割成为50等分，得到51个点的横坐标值，带入到上下翼面曲线函数中，得到51个点横坐标对应的纵坐标，其中上下翼面各51个点，确保了重新拟合的翼型曲线与原曲线一致，采用分段三次样条来拟合这51个点，从而得到50个子域；利用此方法可以得到极惯性矩的具体数值，并且实现优化寻优过程中的自动计算。
2)基于isight优化设计平台实现了在迭代进程中自动更新翼型几何外形、网格模型、自动求解目标及约束函数值，将翼型几何外形更新程序、翼型流场ug几何生成模块、icem网格划分模块、cfx模块集成到一起，进行数据传递，从而实现了优化数学模型的迭代进程。
3)提出了风力机叶片抗颤振翼型多目标优化设计方法，建立了翼型截面极惯性矩与翼型几何外形的关系，通过多目标遗传算法得到了优化后的翼型簇集，其不仅提高了抗颤振稳定性能并且增大了其气动性能(升阻比更大)。
4)通过在原翼型上添加形状扰动函数，不断在原有翼型基础上更新，直到寻找到帕累托最优解集。
5)优化后的翼型集合与原翼型相比较，在大攻角条件下具有良好的抑制层流分离的能力，实现转捩控制，使得翼型失速延期。
附图说明
图1是本发明实施例一中的原始翼型图(DU35)；
图2是本发明实施例一中翼型表面压力及剪力分布示意图；
图3是本发明实施例一中优化设计的帕累托最优解集；
图4a是本发明实施例一中的优化后第一种翼型与原始翼型的比较图；
图4b是本发明实施例一中的优化后第二种翼型与原始翼型的比较图；
图4c是本发明实施例一中的优化后第三种翼型与原始翼型的比较图；
图5是本发明实施例一中优化后翼型与原始翼型的升阻比比较图；
图6是本发明实施例一中的isight平台流程图；
图7是本发明实施例二中优化设计中的翼型流场几何域；
图8是本发明实施例二中优化设计中的网络拓扑图；
图9a是本发明实施例二中优化设计的翼型全局计算网格；
图9b是本发明实施例二中优化设计的翼型局部计算网格；
图10a是本发明实施例二中优化设计中优化翼型绕流流线图；
图10b是本发明实施例二中初始翼型绕流流线图。
具体实施方式
下面对本发明的实施例作详细说明，本实施例在以本发明技术方案为前提下进行实施，给出了详细的实施方式和具体的操作过程，但本发明的保护范围不限于下述的实施例。
实施例一
本实施例以DU35为设计原始翼型，但该方法不限于此翼型，本发明从翼型的形状拟合、极惯性矩的计算方法、数值计算的方法选取等方面进行了创新设计。
翼型的拟合曲线
本实施例是选取函数线性叠加的方法来表达新的翼型，通过在原有翼型基础上添加形状扰动函数，以幅值参数和指数因子为设计变量，对原有翼型进行更新。
极惯性矩的计算
在对极惯性矩计算过程中，首先建立极惯性矩与几何外形的关系，然后计算出翼型的面积和质心坐标；平移坐标系，使坐标原点和质心重合，根据平行轴原理，可求出翼型相对于质心轴系的惯性矩，最终得到翼型的极惯性矩。
数值计算方法
前面虽然给出了极惯性矩计算的具体方法，但要得到具体数值，并实现优化过程中寻找最优解的自动计算，需要采用数值计算的方法，采用高斯-勒让德积分对求解域进行积分计算。
(1)借助ug生成计算域的确定几何模型，流场区域采用半椭圆型流场，这样做的目的是方便引入不同的入流攻角。利用ugupdate，此为ug和ansys的集成计算算程序，将ug的prt文件先储存成x_t文件格式，再导入ansys进行计算，既可以更新prt文件，又可以存储成x_t文件。
(2)利用CAD软件ug导出的parasolid格式的文件，并将其导入到icem中进行网格划分。应该注意网格划分不能过密或者过疏，否则会造成计算精度过差或者计算资源浪费。
(3)在cfxpre中应该合理的创建计算域和边界条件，有利于控制数值计算的精度。
本实施例的风力机叶片抗颤振翼型多目标优化设计方法，具体包括以下步骤：
步骤一：翼型的表达方式是进行优化的先决条件，本发明的翼型的几何形状是在原有翼型基础上添加形状扰动函数，本发明选取解析函数线性叠加法来表达新的翼型。既继承了原翼型的总体特征，又可以节省时间。原始翼型如图1所示。
在原有翼型基础上添加形状扰动函数，保证X轴坐标不变，得到新翼型的型值点坐标。
上翼面翼型的型值点坐标：
yu(i)＝you(i)+Δyu(i)
xu(i)＝xou(i)
下翼面翼型的型值点坐标：
yl(i)＝yol(i)+Δyl(i)
xl(i)＝xol(i)
其中，yu(i)代表新翼型上翼面第i个点的纵坐标，you(i)代表原翼型上翼面第i个点的纵坐标，Δyu(i)代表上翼面第i个点的形状扰动函数，yl(i)代表新翼型下翼面第i点的纵坐标，yol(i)代表原始翼型上翼面第i个点的横坐标，Δyl(i)代表下翼面第i个点的形状扰动函数；xu(i)代表新翼型上翼面第i个点的横坐标，xou(i)代表原始翼型上翼面第i个点的横坐标，xol(i)代表原始翼型下翼面第i个点的横坐标；
其中，上翼面形状扰动函数：
下翼面形状扰动函数：
式中下标u和l分别代表上翼面和下翼面，i代表翼型型值点索引，k代表形状模式索引，g为给点矢量，f(k)代表Hicks-Henne型函数，ξ和η指数因子。
取g(k)＝[0.1，0.2，0.3，0.5，1，3，5，8]；k与f(k)分别是控制翼型厚度与弯度分布关键点变量的个数与系数。采用Hicks-Henne型函数，对翼型上下表面各关键点形状的影响平滑均衡，既可以保证翼型关键点的形状调整，又可以保证关键点与翼型形状平滑。
所述步骤一的原翼型选择为荷兰DU翼型系列中的DU35。
步骤二：翼型极惯性矩计算：
对于风机叶片进行敏感性分析，得到扭转刚度、重心位置、空气密度及弹性轴位置这四个因素对于颤振临界速度影响较大。将重心位置，弹性轴位置向前缘移动或者将叶片放在低密度区域均可以提高抗颤振性能，但这些都是不易实现的。我们对翼型的扭转刚度进行分析，发现提高翼型的扭转刚度可以提高抗颤振性能，翼型的扭转刚度与三个因素有关；
其中，为扭转刚度；G为剪切模量，与材料的选择有关，剪切模量高的材料其扭转刚度好，则抗颤振性能越高；
R为翼型所在界面距离叶根的半径，随着半径的增大，其扭转刚度变小，则抗颤振性能减小；
I为翼型极惯性矩，与翼型外形密切相关，翼型抗颤振性能随极惯性矩增大而增加；
极惯性矩I取决于翼型几何外形，但极惯性矩I为几何外形的隐式表达，需要建立几何外形与极惯性矩之间的关系。首先，计算出翼型的面积和质心坐标点。
求翼型的面积：
其中，S为翼型的面积；b为翼型弦长，单位标准翼型值弦长取1；u代表横坐标，v代表纵坐标；θ1(u)表示翼型下翼面曲线函数；θ2(u)为上翼面曲线函数。
翼型质心坐标：
其中，ρ(u，v)代表不同坐标点处的密度，故假设叶片是均匀的，则ρ(u，v)为常数。
则上式可进一步推导：
坐标原点位于翼型的前缘点，可求出翼型对x，y方向的惯性矩：
平移坐标系，使原点与翼型质心重合，通过平行轴定理，求出翼型相对质心轴系的惯性矩
I1＝Ix-Ay2
I2＝Iy-Ax2
其中：I1代表Y轴相对于质心的惯性矩；
A代表翼型的面积；
y代表Y轴与质心轴的垂直距离；
I2代表X轴相对于质心的惯性矩；
x代表X轴与质心的垂直距离；
翼型极惯性矩为：
I＝I1+I2
步骤三：约束条件的选择，翼型沿着风力机叶片展向分布时，其在不同叶片处存在着厚度要求，所以要对翼型厚度t进行约束，同时要保证优化后的翼型在不同攻角下的升阻比不低于原翼型的升阻比，选择[α1，α2，…，αi]攻角下的升阻比作为约束条件；
厚度约束条件 t1≤t≤t2
升阻比约束条件 cl′/cd′≥cl/cd
式中cl'为优化后的升力系数，cd'为优化后的阻力系数，cl为原翼型的升力系数，cd为原翼型的阻力系数。
本实施例选择DU35翼型为案例，其厚度约束条件：
0.3≤t≤0.4
在追求极惯性矩的同时，肯定会影响翼型的气动性能。升阻比作为衡量翼型气动性能的重要指标，因为DU35翼型在8°攻角附近其升阻比最大(在步骤四中将8°攻角下的升阻比作为一个目标函数)，所以选择[0°，2°，4°，6°，10°，12°，14°]攻角下的升阻比作为约束条件。
升阻比约束条件：
cl′/cd′≥cl/cd
其中，cl′为优化后的升力系数，cd′为优化后的阻力系数，cl为原翼型的升力系数，cd为原翼型的阻力系数；
采用计算流体力学数值计算方法计算翼型上下翼面的压力及剪切力，从而得到翼型上下表面的升阻力系数，并进一步得到翼型的升阻比，优化后翼型升阻比不低于原翼型。翼型表面压力和剪力分布示意图如图2所示。
步骤四：
通过对步骤2进行分析，随着叶片增长，其半径R必然增加，则叶片扭转刚度减小，增加了叶片的柔性，使得叶片更容易发生颤振。翼型极惯性矩越大，其扭转刚度越高，则抗颤振性能越强，故将翼型的极惯性矩I作为一个设计目标；不同的翼型在不同的攻角下具有可变的升阻比，这里，先定义一个攻角α0，在此攻角下翼型具有最大的初始升阻比，将翼型在α0下的升阻比作为另一个目标函数。
在isight平台中通过采用多目标遗传算法NSGA-II对风力机叶片翼型进行优化设计。利用isight平台对matlab、ug、icem、cfx等软件进行集成运算，建立仿真分析流程，设定和修改设计变量以及约束条件和目标函数，自动进行多次分析循环。可以自动的进行翼型坐标点的更新，自动进行几何建模、网格划分以及cfx的自动化计算。
在本实施例中，将极惯性矩和8°攻角下的升阻比作为目标函数：
maximize＝[f1(x)，f2(x)]
f1(x)＝I
f2(x)＝cl/cd
在isight平台中可以对CAD/CAE软件进行集成，用户在isight平台中自行修改设计变量、约束条件和升阻比，采用多目标遗传因子算法NSGA-II对仿真优化流程自动进行循环。多目标优化的帕累托最优解集如图3所示。
在帕累托最优解集中选择三组翼型与原有翼型比较，如图4a，4b，4c所示；
选择的三组翼型与原有翼型的升阻比进行比较，如图5所示。
步骤五：通过步骤一到步骤四，实现了对翼型抗颤振多目标的设计。此优化流程设计过程中，重点在于如何实现迭代进程中自动更新翼型的几何模型，网格划分，自动求解目标以及求解约束函数。本案例中将翼型几何外形更新程序(matlab文件)、翼型流场几何生成模块(ug文件)、网格划分模块(icem)、数值计算模块(cfx)共同集成到isight平台中，实现数据传递并进行优化迭代进程。同时isight平台可以将设计变量，约束条件，目标函数值的变化历程以可视化的方式直观的显示出来，还可以提供多种数据分析和后处理工具。在isight平台的流程图如图6所示。
实施例二
本实施例以DU35为原始翼型，该翼型最大相对厚度为0.35，位于27.50％翼弦处，该翼型弧度为-1.96％，位于30％的翼弦处。
计算域：
借助ug生成确定计算域的几何模型，流场区域采用半椭圆几何形状，完整的椭圆长轴为100b短轴长度为50b，采用布尔运算从半椭圆流场中减去翼型几何曲线构成的几何体，便可得到翼型绕流流场的计算域的几何模型。如图7所示。
计算网格：
从ug导出parasolid格式文件，然后导入到网格划分软件ICEM中。在翼型边界层区域，考虑到近壁面粘滞作用，分别采用内O-block和外O-block拓扑生成边界层网格，首先，使用内O-block确定翼型几何拓扑，并将内部block删除，其次，使用外O-block确定翼周边界层拓扑区域，以生成边界层网格。近壁层第一层网格高度为0.2mm，确保粘性层第一层至少有一层网格。对网格数目进行划分，网格数为82852。网格拓扑图如图8所示，网格图如图9a、9b所示。
计算方法：
在cfx-pre中对流体分析进行前处理，设置参数。设置流体类型为理想气体，收敛精度为1e-6。cfx提供了丰富的后处理选项，通过后处理可以直观地显示出流场计算结果，部分流线图如图10a、10b所示。
综上所述，上述实施例的风力机叶片抗颤振翼型簇设计方法，通过在原有的翼型基础上添加翼型扰动函数来修改翼型的几何外形，提高其抗颤振稳定性能并增强其气动性能(升阻比更大)；通过经典抗颤振理论的研究，得出翼型的极惯性矩对其抗颤振性能具有突出的影响，其值可通过高斯-勒让德积分计算得出；并将最大升阻比和极惯性矩作为优化目标，在isight平台中融合翼型扰动函数、高斯-勒让德积分计算、网格划分、流体力学计算、多目标遗传算法来获取优化后的新翼型簇,其抗颤振性能及气动性能均得到了显著的提高，值得被推广使用。
尽管上面已经示出和描述了本发明的实施例，可以理解的是，上述实施例是示例性的，不能理解为对本发明的限制，本领域的普通技术人员在本发明的范围内可以对上述实施例进行变化、修改、替换和变型。
Claims (9)
1.一种风力机叶片抗颤振翼型簇设计方法，其特征在于，包括以下步骤：
S1：翼型表达
对原有翼型添加形状扰动函数来表达新的翼型；
在所述步骤S1中，在原有翼型基础上添加形状扰动函数，保证X轴坐标不变，得到新翼型的型值点坐标：
上翼面翼型的型值点坐标：
yu(i)＝you(i)+Δyu(i)
xu(i)＝xou(i)
下翼面翼型的型值点坐标：
yl(i)＝yol(i)+Δyl(i)
xl(i)＝xol(i)
其中，yu(i)代表新翼型上翼面第i个点的纵坐标，you(i)代表原有翼型上翼面第i个点的纵坐标，Δyu(i)代表上翼面第i个点的形状扰动函数，yl(i)代表新翼型下翼面第i点的纵坐标，yol(i)代表原有翼型上翼面第i个点的横坐标，Δyl(i)代表下翼面第i个点的形状扰动函数；xu(i)代表新翼型上翼面第i个点的横坐标，xou(i)代表原有翼型上翼面第i个点的横坐标，xol(i)代表原有翼型下翼面第i个点的横坐标；
S2：颤振临界速度影响因素分析
对风机叶片进行敏感性分析，得到颤振临界速度影响因素，提高颤振临界速度影响因素中的扭转刚度以提高抗颤振性能；
S3：选择约束条件
确定翼型的约束条件为厚度约束条件与升阻比约束条件；
S4：确定翼型目标函数
将翼型在攻角α0下的升阻比与翼型的极惯性矩I作为目标函数，通过采用多目标遗传算法对风力机叶片翼型进行优化设计。
2.根据权利要求1所述的一种风力机叶片抗颤振翼型簇设计方法，其特征在于：在所述步骤S1中，所述原有翼型为荷兰DU翼型系列中的DU35。
3.根据权利要求1所述的一种风力机叶片抗颤振翼型簇设计方法，其特征在于：上、下翼面形状扰动函数如下：
上翼面形状扰动函数：
下翼面形状扰动函数：
其中，下标u和l分别代表上翼面和下翼面，i代表翼型型值点索引，k与f(k)分别是控制翼型厚度与弯度分布关键点变量的个数与系数，g为给点矢量，f(k)代表Hicks-Henne型函数，ξ和η为指数因子。
4.根据权利要求1所述的一种风力机叶片抗颤振翼型簇设计方法，其特征在于：在所述步骤S2中，翼型的扭转刚度表达式如下：
其中，G为剪切模量，R为翼型所在界面距离叶根的半径，I为翼型极惯性矩。
5.根据权利要求4所述的一种风力机叶片抗颤振翼型簇设计方法，其特征在于：所述极惯性矩I的计算过程如下：
S21：计算出翼型的面积和质心坐标点：
翼型的面积：
其中，S为翼型的面积；b为翼型弦长，单位标准翼型值弦长取1；u代表横坐标，v代表纵坐标；θ1(u)表示翼型下翼面曲线函数；θ2(u)为上翼面曲线函数；
翼型质心坐标：
其中，ρ(u，v)代表不同坐标点处的密度，假设叶片是均匀的，则ρ(u，v)为常数；
则上式可进一步推导：
S22：坐标原点位于翼型的前缘点，求出翼型对x，y方向的惯性矩：
S23：平移坐标系，使原点与翼型质心重合，通过平行轴定理，求出翼型相对质心轴系的惯性矩：
I1＝Ix-Ay2
I2＝Iy-Ax2
其中：I1代表Y轴相对于质心的惯性矩；
A代表翼型的面积；
y代表Y轴与质心的垂直距离；
I2代表X轴相对于质心的惯性矩；
x代表X轴与质心的垂直距离；
则翼型极惯性矩为：
I＝I1+I2。
6.根据权利要求5所述的一种风力机叶片抗颤振翼型簇设计方法，其特征在于：在所述步骤S3中，厚度约束条件为t1≤t≤t2，式中t为翼型厚度，t1、t2为设定值。
7.根据权利要求6所述的一种风力机叶片抗颤振翼型簇设计方法，其特征在于：在所述步骤S3中，优化后的翼型在不同攻角下的升阻比不低于原翼型的升阻比，选择[α1，α2，…，αi]攻角下的升阻比作为约束条件，升阻比约束条件为：
cl′/cd′≥cl/cd
式中cl'为优化后的升力系数，cd'为优化后的阻力系数，cl为原有翼型的升力系数，cd为原有翼型的阻力系数。
8.根据权利要求7所述的一种风力机叶片抗颤振翼型簇设计方法，其特征在于：在所述步骤S4中，目标函数为：
maximize＝[f1(x)，f2(x)]
f1(x)＝I
f2(x)＝cl/cd。
9.根据权利要求1所述的一种风力机叶片抗颤振翼型簇设计方法，其特征在于：在所述步骤S4中，多目标遗传算法为NSGA-II。
Priority Applications (1)
| Application Number | Priority Date | Filing Date | Title |
|---|---|---|---|
| CN202110962098.9A CN113569360B (zh) | 2021-08-20 | 2021-08-20 | 一种风力机叶片抗颤振翼型簇设计方法 |
Applications Claiming Priority (1)
| Application Number | Priority Date | Filing Date | Title |
|---|---|---|---|
| CN202110962098.9A CN113569360B (zh) | 2021-08-20 | 2021-08-20 | 一种风力机叶片抗颤振翼型簇设计方法 |
Publications (2)
| Publication Number | Publication Date |
|---|---|
| CN113569360A CN113569360A (zh) | 2021-10-29 |
| CN113569360B true CN113569360B (zh) | 2024-03-22 |
Family
ID=78172308
Family Applications (1)
| Application Number | Title | Priority Date | Filing Date |
|---|---|---|---|
| CN202110962098.9A Active CN113569360B (zh) | 2021-08-20 | 2021-08-20 | 一种风力机叶片抗颤振翼型簇设计方法 |
Country Status (1)
| Country | Link |
|---|---|
| CN (1) | CN113569360B (zh) |
Families Citing this family (3)
| Publication number | Priority date | Publication date | Assignee | Title |
|---|---|---|---|---|
| CN114065423B (zh) * | 2021-11-12 | 2023-03-28 | 西北工业大学 | 快速评估航空发动机风扇叶片颤振的方法 |
| GB2623958A (en) * | 2022-10-31 | 2024-05-08 | Smallspark Space Systems | Design of an aerodynamic component |
| CN116776600B (zh) * | 2023-06-21 | 2024-04-12 | 安徽工业大学 | 基于自适应代理模型的风力机叶片优化设计方法及系统 |
Citations (2)
| Publication number | Priority date | Publication date | Assignee | Title |
|---|---|---|---|---|
| CN104612892A (zh) * | 2014-12-30 | 2015-05-13 | 中国科学院工程热物理研究所 | 一种风力机翼型的多学科优化设计方法 |
| CN112507471A (zh) * | 2020-10-23 | 2021-03-16 | 湖北工业大学 | 一种大攻角范围条件下垂直轴风力机翼型设计方法 |
-
2021
- 2021-08-20 CN CN202110962098.9A patent/CN113569360B/zh active Active
Patent Citations (2)
| Publication number | Priority date | Publication date | Assignee | Title |
|---|---|---|---|---|
| CN104612892A (zh) * | 2014-12-30 | 2015-05-13 | 中国科学院工程热物理研究所 | 一种风力机翼型的多学科优化设计方法 |
| CN112507471A (zh) * | 2020-10-23 | 2021-03-16 | 湖北工业大学 | 一种大攻角范围条件下垂直轴风力机翼型设计方法 |
Non-Patent Citations (2)
| Title |
|---|
| 基于EI加点准则与代理模型的风力机专用翼型气动与结构优化设计;张强等;中国电机工程学报;20210811;第1-11页 * |
| 风力机叶片三维模型颤振问题;王伟等;噪声与振动控制;20131031;第33卷(第5期);第90-96页 * |
Also Published As
| Publication number | Publication date |
|---|---|
| CN113569360A (zh) | 2021-10-29 |
Similar Documents
| Publication | Publication Date | Title |
|---|---|---|
| CN113569360B (zh) | 一种风力机叶片抗颤振翼型簇设计方法 | |
| WO2021077855A1 (zh) | 一种直升机旋翼翼型确定方法及系统 | |
| CN104612892B (zh) | 一种风力机翼型的多学科优化设计方法 | |
| CN111079228B (zh) | 一种基于流场预测的气动外形优化方法 | |
| WO2020134388A1 (zh) | 一种基于随机等几何分析的叶片高刚度设计方法 | |
| CN109190233B (zh) | 一种结构拓扑优化方法 | |
| CN106126791B (zh) | 一种考虑几何不确定性的高超声速机翼气动力/热分析方法 | |
| CN107391891A (zh) | 一种基于模型融合方法的大展弦比机翼优化设计方法 | |
| CN106407577B (zh) | 一种模拟风力机尾流的改进致动面模型建立方法 | |
| CN114528634A (zh) | 高隐身高机动布局飞机的弹性机翼气动隐身优化设计方法 | |
| CN108829970A (zh) | 基于涡动力学的轴流风扇叶片优化设计方法及优化系统 | |
| CN111400834A (zh) | 风力发电机组叶片翼型气动优化设计方法、模型及装置 | |
| CN114154434A (zh) | 一种无尾飞翼布局的多约束精细化气动优化设计方法 | |
| Song et al. | Airfoil optimization using a machine learning-based optimization algorithm | |
| Liu et al. | Optimization of nano-rotor blade airfoil using controlled elitist NSGA-II | |
| CN117272870B (zh) | 基于自适应重叠网格的动态编队飞行数值仿真方法 | |
| CN112231836B (zh) | 一种基于遗传算法和数值仿真的翼型优化方法 | |
| CN109992893A (zh) | 一种螺旋桨气动外形优化设计方法 | |
| CN112199777B (zh) | 一种适用于模化仿生前缘流场特征的方法 | |
| CN117744242A (zh) | 一种倾转旋翼桨叶气动外形确定方法、系统及电子设备 | |
| CN117494611A (zh) | 一种基于有粘当地流活塞理论的超声速颤振分析方法 | |
| CN108664742A (zh) | 短舱参数化设计的多目标优化方法 | |
| CN116306088B (zh) | 基于共形参数化的多尺度扫掠曲面结构设计方法及设备 | |
| CN114692455B (zh) | 一种基于分割等距重组曲线的复合材料层合板的建模方法 | |
| Pang et al. | Intelligence algorithm for optimization design of the low wind speed airfoil for wind turbine |
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 | ||
| GR01 | Patent grant | ||
| GR01 | Patent grant |