环境参数对风电齿轮箱传动系统疲劳损伤的影响

伍源，朱才朝，谭建军，宋朝省，张会阳; WU Yuan; ZHU Caichao; TAN Jianjun; SONG Chaosheng; ZHANG Huiyang

网刊加载中。。。

使用Chrome浏览器效果最佳，继续浏览，你可能不会看到最佳的展示效果，

确定继续浏览么?

复制成功，请在其他浏览器进行阅读

环境参数对风电齿轮箱传动系统疲劳损伤的影响 PDF

- ORCID：
伍源 ¹
✉
- ORCID：
朱才朝 ¹
✉
- ORCID：
谭建军 ¹
- ORCID：
宋朝省 ¹
- ORCID：
张会阳 ²

1. 重庆大学机械传动国家重点实验室，重庆 400044； 2. 中国船舶重工集团海装风电股份有限公司，重庆 401122

中图分类号： TH113

最近更新：2024-04-01

DOI：10.11835/j.issn.1000-582X.2022.115

目录contents

摘要

在风电机组全寿命周期内，长期风速概率分布会使风电齿轮箱传动系统动载荷出现随机特性，影响其疲劳损伤预估精度。笔者提出了一种考虑长期风速概率分布特征的风电齿轮箱传动系统疲劳损伤预估方法，通过建立大功率海上风电机组OpenFAST-SIMPACK联合仿真模型，计算不同平均风速与湍流强度组合工况下的风电齿轮箱传动系统齿轮短期疲劳损伤，进而采用代理模型技术重构“平均风速、湍流强度-短期疲劳损伤”映射关系，预测齿轮长期疲劳损伤。研究结果表明：风电齿轮箱传动系统低速级太阳轮容易发生接触疲劳失效；在额定风速以下，低速级太阳轮短期疲劳损伤与平均风速呈正相关，在额定风速附近，平均风速与湍流强度的随机特性均会增大其长期疲劳损伤不确定性，增大其疲劳失效风险。

关键词

风电齿轮箱; 传动系统; 疲劳损伤; 随机特性

大功率海上风电机组是目前最有效地开发海上风能资源的重大海洋工程装备之一。风电齿轮箱是风电机组中传递力与运动的关键传动装置。据统计，海上风电齿轮箱故障率比陆上高5%以上^[

1]。这是由于海面粗糙度比陆地小，海表风速大、季节性变化明显，因而由风剪切效应造成的作用在风电齿轮箱的气动载荷波动幅值大、随机性强，造成海上风电齿轮箱疲劳失效间隔时间短^{[参考文献 2

百度学术}2]。同时，海上风电装备运维成本比陆上高10%以上^{[参考文献 1

百度学术}1]，对海上风电齿轮箱可靠性设计提出了更高要求。因此，开展考虑长期风速概率分布特征的风电齿轮箱传动系统疲劳损伤随机特性分析对指导海上风电齿轮箱设计、提高其可靠性具有重要意义。

近年来，国内外学者对风电齿轮箱传动系统疲劳损伤开展了大量研究。Dong等^[

3]分别建立风电机组整机动力学模型和齿轮箱动力学模型，利用解耦方法计算了不同风速作用下风电机组全局气动载荷，分析了不同气动载荷作用下风电齿轮箱传动系统齿轮接触应力概率分布特征；向东等^{[参考文献 4

百度学术}4]建立随机风载作用下风电齿轮箱传动系统短期疲劳损伤计算模型，分析了随机风速下各齿轮的短期弯曲疲劳与接触疲劳损伤；Nejad等^{[参考文献 5-7}5-7]针对风电齿轮箱传动系统疲劳损伤分析开展了大量研究，包括考虑长期风速概率分布的10 MW与5 MW风电齿轮箱传动系统齿轮长期疲劳损伤计算模型、考虑随机波浪作用下4种浮式平台对风电齿轮箱传动系统长期疲劳损伤影响等；熊中杰^{[参考文献 8

百度学术}8]建立考虑时变啮合刚度和阻尼的齿轮箱动力学模型，分析不同风速、不同载荷类型下的风电齿轮箱疲劳损伤规律。文献[3-8]报道了有效预测风电齿轮箱传动系统长期疲劳损伤的方法，通过等间距离散平均风速，计算各风速工况下的齿轮短期疲劳损伤后加权求和得到长期疲劳损伤，但离散间距取值主要依靠设计者经验，离散后所需仿真的环境参数组合工况数量多，计算耗时。

为了提高传动系统长期疲劳损伤的计算效率，部分学者将多项式响应面(polynomial response surface，PRS)、克里金(kriging，KRG)、人工神经网络(artificial neural network，ANN)和径向基函数(radial basis function，RBF)等代理模型方法引入到复杂系统动力学分析中^[

9-12]。Juan等^{[参考文献 13

百度学术}13]建立了环境参数与风电机组等效疲劳载荷和发电量之间的PRS代理模型，并开展环境参数全局灵敏度分析；Li等^{[参考文献 14

百度学术}14]考虑风浪等环境工况相关性，利用KRG代理模型预测漂浮式风电机组系泊缆、塔架底座和顶部横截面的长期疲劳损伤；Zhang等^{[参考文献 15

百度学术}15]利用ANN和KRG等代理模型计算动态环境载荷下浮式风电机组系泊线的疲劳失效概率；Wilkie等^{[参考文献 16

百度学术}16]利用高斯过程回归代理模型预测了海上风电机组疲劳损伤、疲劳可靠性等性能，分析风和波浪等环境参数对性能指标的影响。文献[13-16]报道了高效的风电机组疲劳损伤分析的计算方法，但忽略了平均风速与湍流强度概率分布对风电齿轮箱传动系统疲劳损伤随机特性的影响。

笔者考虑随机风速与风电机组拓扑结构，建立大功率海上风电机组OpenFAST-SIMPACK联合仿真模型；考虑长期风速概率分布特征，利用Copula函数构建平均风速-湍流强度联合概率分布，建立平均风速-湍流强度关键工况集合，计算对应的风电齿轮箱传动系统齿轮短期疲劳损伤；采用代理模型技术重构“平均风速、湍流强度-疲劳损伤”映射关系，预测齿轮长期疲劳损伤。

1 海上风电机组OpenFAST-SIMPACK联合仿真建模

1.1 风电机组结构及运行原理

图1所示为某型风电机组齿轮箱传动系统结构及原理。海上风电机组主要由风轮、主轴、齿轮箱、发电机、塔筒及塔架等组成。在随机风速作用下，叶片将风能转变为机械能，并通过主轴驱动齿轮箱转动。齿轮箱主要包括低速级、中间级和高速级三级传动结构，其中低速级和中间级为行星斜齿轮传动，高速级为平行轴斜齿轮传动。发电机通过联轴器与高速轴相连。海上风电机组基本设计参数如表1所示，风电齿轮箱主要设计参数如表2所示。

图1 某型风电机组齿轮箱传动结构及原理

Fig. 1 Transmission structure and principle of a wind turbine gearbox

表1 某型海上风电机组主要参数

Table 1 Main parameters of an offshore wind turbine

项目	参数	项目	参数
额定功率/MW	5.8	轮毂高度/m	105.50
类型	上风向、3叶片	风轮扫掠直径/m	171.44
切入、额定、切出风速/(m⋅s^-1)	3.50、9.58、25.00	齿轮箱总传动比	1:120.7
风轮切入、额定转速/(r⋅min^-1)	5.4、10.1	发电机额定转速/(r⋅min^-1)	1 219.07

表2 某型风电齿轮箱设计参数

Table 2 Design parameters of a wind turbine gearbox

速度等级	齿轮箱构件	齿数	模数/mm	螺旋角/( $°$ )	压力角/( $°$ )
低速级	内齿圈	93	24	5	20
	行星轮	29	24	5	20
	太阳轮	32	24	5	20
中间级	内齿圈	118	17	8.5	20
	行星轮	47	17	8.5	20
	太阳轮	23	17	8.5	20
高速级	大齿轮	121	12	9	20
高速级	小齿轮	24	12	9	20

1.2 OpenFAST-SIMPACK联合仿真模型

图2所示为大功率海上风电机组OpenFAST-SIMPACK联合仿真模型，主要包括2类不同层级的子模型，即风电机组整机全局耦合模型与风电齿轮箱传动系统动力学模型。首先，利用OpenFAST^[

17]建立大功率海上风电机组整机全局耦合模型，控制策略采用变速-变桨控制^{[参考文献 18

百度学术}18]；为了提高计算效率，风电齿轮箱传动系统简化为传动比。然后，利用SIMPACK^{[参考文献 19

百度学术}19]建立风电齿轮箱传动系统动力学模型，其中齿轮副采用切片理论建模、传动轴与行星架简化为刚体，轴承利用6×6刚度矩阵模拟；发电机和联轴器采用集中质量单元建模^{[参考文献 20

百度学术}20]。

图2 大功率海上风电机组OpenFAST-SIMPACK联合仿真模型

Fig. 2 OpenFAST-SIMPACK co-simulation model for high-power offshore wind turbines

为了实现OpenFAST-SIMPACK联合仿真，利用Matlab动态修改OpenFAST风电机组整机全局耦合模型工况.inp文件，计算任意给定平均风速-湍流强度环境参数下的风电机组轮毂处6自由度气动载荷；利用Matlab将气动载荷文件格式转为SIMPACK风电齿轮箱传动系统动力学模型载荷.afs文件，并同时调用其宏命令.sjs文件进行仿真，计算风电齿轮箱传动系统齿轮副动态啮合力。

2 考虑长期风速分布的风电齿轮箱传动系统疲劳损伤预估

2.1 长期风速概率分布

根据IEC61400-1标准^[

21]，将10 min时序风速的平均值作为平均风速，其概率密度函数通常服从威布尔分布

f_{v} (u) = \frac{θ}{λ} {(\frac{u}{λ})}^{θ - 1} e^{- {(\frac{u}{λ})}^{θ}}

，

(1)

式中： $u$ 为平均风速； $θ$ 为威布尔分布的形状参数； $λ$ 为威布尔分布的尺度参数。

将10 min时序风速的变化情况作为湍流强度，其概率密度函数通常服从伽马分布^[

22]

f_{I} (ϑ) = \frac{β^{α}}{Γ (α)} ϑ^{α - 1} e^{(- β ϑ)}

，

(2)

式中： $ϑ$ 为湍流强度； $α$ 为伽马分布的尺度参数； $β$ 为伽马分布的形状参数； $Γ (α)$ 为伽马函数； $Γ (α) = \int_{0}^{\infty} x^{α - 1} e^{- x} d x$ 。式(1)和(2)中各项特征参数均采用最小二乘法估计^[

23]。

根据式(1)和(2),通过积分可得平均风速 $f_{v} (u)$ 与湍流强度 $f_{I} (ϑ)$ 的累积分布函数分别为 $F_{v} (u)$ 和 $F_{I} (ϑ)$ 。根据Sklar定理^[

24]，利用二元Copula函数构建平均风速-湍流强度联合概率密度函数

f_{v I} (u, ϑ)

：

f_{v I} (u, ϑ) = c (F_{v} (u), F_{I} (ϑ)) f_{v} (u) f_{I} (ϑ)

。

(3)

式中， $c (F_{v} (u), F_{I} (ϑ)) = \frac{\partial C (F_{v} (u), F_{I} (ϑ))}{\partial F_{v} (u) \partial F_{I} (ϑ)}$ 。 $C (F_{v} (u), F_{I} (ϑ))$ 函数通常包括Gumbel Copula、Clayton Copula、Frank Copula等^[

25]，其中，Gumbel Copula表达式为

C (F_{v} (u), F_{I} (ϑ)) = e x p \{- {[{(- l n F_{v} (u))}^{τ} + {(- l n F_{I} (ϑ))}^{τ}]}^{\frac{1}{τ}}\} ，

(4)

式中，采用核密度对 $τ$ 进行参数估计^[

26]。

2.2 关键环境工况

为了从 $f_{v I} (u, ϑ)$ 中高效地抽选权重占比较高的平均风速-湍流强度环境参数组合工况，结合拉丁超立方抽样法^[

27]与最大差异法^{[参考文献 28

百度学术}28]两者优点，建立图3所示关键环境工况分析流程，以保证在抽样样本数量较少的条件下，同时保持各抽样维度均匀性与差异最大性，具体步骤如下。

图3 关键环境工况分析流程

Fig. 3 Analysis flow of key environmental conditions

步骤1：利用Nataf变换^[

29]对平均风速和湍流强度进行独立化与标准正态化预处理。

步骤2：采用拉丁超立方抽样方法^[

27]对步骤1预处理后的平均风速和湍流强度2个维度进行随机抽样，生成

M

组拉丁超立方样本

X_{i} ， i = 1,2, \dots, M

。

步骤3：利用Nataf变换的逆变换过程^[

30]对步骤2生成的抽样样本数据逆处理，得到

X_{i}^{'} ， i = 1,2, \dots, M

。

步骤4：计算抽样样本数据中任意2组数据的欧式距离 $d i s t (X_{i}^{'}, X_{j}^{'}) ， i = 1,2, \dots, M ， j = 1,2, \dots, M ， i \neq j$ ，以欧式距离最大的抽样样本作为训练样本集 $Q$ ，剩余部分作为测试样本集 $D$ ；以 $Q$ 与 $D$ 中各样本欧式距离最小作为不相似度 $d (X_{i}^{'}, D) ， X_{i}^{'} \in Q$ ，将 $D$ 中 $d (X_{i}^{'}, D)$ 最大的样本加入 $Q$ ，直至训练样本集 $Q$ 的样本数量达到设定值。

通过步骤1~4实现从平均风速-湍流强度联合分布中获得指定数量且权重占比较大的关键环境工况，用于风电齿轮箱传动系统疲劳损伤计算。

2.3 齿轮长期疲劳损伤

根据大功率海上风电机组OpenFAST-SIMPACK联合仿真模型，计算各关键工况下风电齿轮箱传动系统齿轮副时序动态啮合力 $F_{i}$ 。假设时序动态啮合力在 $F_{1}$ 和 $F_{2}$ 之间有 $T_{i}$ 个时间段，利用式(5)计算 $F_{1}$ 循环次数^[

31]；得到啮合力块后，根据ISO 6336-6标准计算齿轮接触应力(

δ_{H}

)与弯曲应力(

δ_{F}

)^{[参考文献 31

百度学术}31]，得到齿轮接触与弯曲应力块^{[参考文献 32-33}32-33]，最终应力载荷谱，如图4所示。图中

δ

为应力，

δ_{1}

和

δ_{2}

分别第1和第2个应力块的幅值。

n_{i} = \sum_{j = 1}^{T_{i}} k \frac{t_{j} w_{j}}{2 π}

，

（5）

δ_{H} = Z_{H} Z_{E} Z_{ε} Z_{β} \sqrt[]{\frac{F_{t}}{b d_{1}} \frac{i_{g e a r} + 1}{i_{g e a r}} K_{A} K_{V} K_{H β} K_{H α}}

，

（6）

δ_{F} = (\frac{F_{t}}{b m_{n}} Y_{F} Y_{S} Y_{β} Y_{B} Y_{D T}) K_{A} K_{V} K_{F β} K_{F α} 。

（7）

式中：对于内齿圈和太阳轮， $k$ 等于行星轮数量； $n_{i}$ 为啮合力块 $i$ 循环次数; $T_{i}$ 为啮合力块 $i$ 的总时间段数； $t_{j}$ 为啮合力块 $i$ 的第 $j$ 个时间段； $w_{j}$ 为啮合力块 $i$ 在时间段 $j$ 内齿轮平均转速； $F_{t}$ 为端面分度圆上的名义切向力； $b$ 为工作齿宽； $d_{1}$ 为小齿轮分度圆直径； $i_{g e a r}$ 为齿轮传动比； $Z_{H}$ 为节点区域系数； $Z_{E}$ 为材料弹性系数； $Z_{ε}$ 为重合度系数； $Z_{β}$ 为接触应力计算的螺旋角系数； $K_{A}$ 为使用系数； $K_{V}$ 为动载系数； $K_{H β} 、 K_{F β}$ 分别为接触应力计算和弯曲应力计算的齿向载荷分布系数； $K_{H α}$ 、 $K_{F α}$ 分别为接触应力计算和弯曲应力计算的齿间载荷分配系数； $Y_{F}$ 为齿形系数； $Y_{S}$ 为应力修正系数； $Y_{β}$ 为弯曲应力计算的螺旋角系数； $Y_{B}$ 为轮毂厚度系数； $Y_{D T}$ 为重合度系数； $m_{n}$ 为法向模数。

图4 齿轮接触与弯曲应力块统计

Fig. 4 Gear contact and bending stress block statistics

基于式(5)~(7)，根据线性疲劳累积损伤理论^[

34]，得到各关键工况下齿轮短期疲劳损伤

D_{c} (u, ϑ)

：

D_{c} (u, ϑ) = \sum_{i} \frac{n_{i} (u, ϑ)}{C δ_{i}^{- m}} = \frac{N_{b}}{C} \sum_{i} δ_{i}^{m} \cdot f (δ) \cdot ∆ δ \approx \frac{N_{b}}{C} \int_{0}^{\infty} δ_{i}^{m} \cdot f (δ) d δ

。

(8)

式中： $N_{b}$ 为应力块数量； $δ_{i}$ 为第 $i$ 个应力块幅值； $C$ 和 $m$ 分别为齿轮S-N曲线参数，利用齿轮的国际标准ISO 6336-2^[

32]、ISO 6336-3^{[参考文献 33

百度学术}33]、ISO 6336-5^{[参考文献 35

百度学术}35]计算得到；

f (δ)

为应力块分布函数。当

f (δ)

服从正态分布时，

D_{c} (u, ϑ)

可改写为

D_{c} (u, ϑ) = \frac{N_{b}}{C} \int_{0}^{\infty} δ_{i}^{m} \cdot f (δ) d δ = \frac{N_{b}}{C} \cdot \int_{0}^{\infty} δ_{i}^{m} \cdot \frac{1}{\sqrt[]{2 π} σ} e^{(- \frac{(δ - μ)}{2 σ^{2}})} d δ 。

(9)

在式(9)基础上，考虑各关键工况出现的概率及持续时间，可得齿轮长期疲劳损伤 $D_{t o l}$ 为：

D_{t o l} = \int_{u_{1}}^{u_{2}} \int_{ϑ_{1}}^{ϑ_{2}} T \cdot f_{v I} (u, ϑ) \cdot D_{c} (u, ϑ) d u d ϑ = \sum_{u_{1}}^{u_{2}} \sum_{ϑ_{1}}^{ϑ_{2}} T \cdot f_{v I} (u, ϑ) \cdot D_{c} (u, ϑ) ∆ u ∆ ϑ

。

(10)

式中： $u_{i}$ 为第 $i$ 个点风速； $ϑ_{i}$ 为第 $i$ 个点湍流强度；T为风电机组设计寿命，T=25 a。

为了高效计算 $D_{t o l}$ ，采用代理模型方法^[

36]计算各关键工况下齿轮短期疲劳损伤

D_{c} (u, ϑ)

，建立“平均风速、湍流强度-短期疲劳损伤”映射函数。为了确定拟合精度较高的代理模型，分别采用克里金模型(KRG)^{[参考文献 37

百度学术}37]、多项式响应面模型(PRS)^{[参考文献 38

百度学术}38]和径向基函数模型(RBF)^{[参考文献 39

百度学术}39]构建齿轮短期疲劳损伤函数，分别为

D_{c}^{K R G} 、 D_{c}^{P R S}

和

D_{c}^{R B F}

，如式(11)所示。利用全局精度指标均方根误差（root mean squared error, RMSE）、误差平方(

R^{2}

)和局部精度指标最大绝对误差（maximum absolute error，MAE）评价3种代理模型拟合精度^{[参考文献 40

百度学术}40]，分别如式(12)~(14)所示。

\{\begin{array}{l} D_{c}^{K R G} (u, ϑ) = \sum_{i = 1}^{n} κ_{i} p_{i} (u, ϑ) + Z (u, ϑ) ， \\ D_{c}^{P R S} (u, ϑ) = β_{0} + \sum_{i = 1}^{n} β_{i} (u^{i} + ϑ^{i}) + \sum_{i = 1}^{n} \sum_{j = 1}^{n} β_{i j} u^{i} ϑ^{j} ， \\ D_{c}^{R B F} (u, ϑ) = \sum_{i = 1}^{N_{s}} w_{i} ϕ (‖(u, ϑ) - (u_{i}, ϑ_{i})‖) 。 \end{array}

(11)

式中： $\sum_{i = 1}^{n} κ_{i} p_{i} (u, ϑ)$ 为多项式全局近似模型； $Z (u, ϑ)$ 是均值为0、方差为 $σ^{2}$ 、协方差非零的局部偏差模型； $β_{0}, β_{i}, β_{i j}$ 分别为多项式函数待定系数； $‖(u, ϑ) - (u_{i}, ϑ_{i})‖$ 为训练样本工况之间的欧氏距离； $w_{i}$ 为权重系数； $ϕ (\cdot)$ 为径向基函数； $N_{s}$ 为训练样本个数。

R M S E = \sqrt[]{\frac{\sum_{i = 1}^{n_{t}} {(D_{c i} (u, ϑ) - \hat{D_{c i}^{s} (u, ϑ)})}^{2}}{n_{t}}} ，

（12）

M A E = m a x |D_{c i} (u, ϑ) - \hat{D_{c i}^{s} (u, ϑ)}|

，

（13）

R^{2} = 1 - \frac{\sum_{i = 1}^{n_{t}} {(D_{c i} (u, ϑ) - \hat{D_{c i}^{s} (u, ϑ)})}^{2}}{\sum_{i = 1}^{n_{t}} {(D_{c i} (u, ϑ) - \bar{D_{c} (u, ϑ)})}^{2}} 。

（14）

式中： $n_{t}$ 为测试样本数量； $D_{c i}^{s} (u, ϑ)$ 为第 $i$ 组测试样本下的齿轮短期疲劳损伤； $\hat{D_{c i}^{s} (u, ϑ)}$ 为代理模型预测第 $i$ 组测试样本下的齿轮短期疲劳损伤； $\bar{D_{c}^{s} (u, ϑ)}$ 为所有测试样本下的齿轮短期疲劳损伤均值。

图5所示为3种代理模型对低速级太阳轮短期疲劳损伤的预测结果，对应的预测精度对比如表3所示。可以看出，径向基代理模型对齿轮短期疲劳损伤的预测精度最好，齿轮接触与弯曲疲劳损伤 $R^{2}$ 分别为0.990和0.989，因此，选用径向基代理模型构建“平均风速、湍流强度-短期疲劳损伤”映射函数。

图5 3种代理模型对低速级太阳轮短期疲劳损伤的预测效果

Fig. 5 Prediction results of 3 surrogate models on short-term fatigue damage of the low-speed sun gear

表3 3种代理模型对低速级太阳轮短期疲劳损伤的预测精度对比

Table 3 Prediction accuracies of 3 surrogate models for short-term fatigue damage of the low-speed sun gear

代理模型类型	弯曲疲劳			接触疲劳
代理模型类型	RMSE	MAE	$R^{2}$	RMSE	MAE	$R^{2}$
克里金模型	0.511	1.070	0.900	0.745	1.500	0.890
多项式响应面模型	0.339	0.850	0.955	0.466	1.158	0.960
径向基模型	0.163	0.341	0.989	0.220	0.371	0.990

3 结果讨论与分析

根据中国某海上风电场历年风速统计数据，建立平均风速-湍流强度联合概率密度函数，并选取关键环境工况，分别计算各关键环境工况下风电齿轮箱传动系统齿轮副动态啮合力(各工况仿真时间100 s)，并计算齿轮短期和长期疲劳损伤；最后分析环境参数对齿轮短期和长期疲劳损伤的影响规律。

3.1 环境参数对风电齿轮箱传动系统齿轮短期疲劳损伤影响

图6与图7分别为平均风速概率密度函数和湍流强度概率密度函数，其中平均风速威布尔分布函数为 $f_{v} (u) = \frac{u^{1.64}}{92.04} e^{- {(\frac{u}{8.01})}^{2.64}}$ ，湍流强度伽马分布函数为 $f_{Ι} (ϑ) = \frac{33 . 33^{2.49}}{Γ (2.49)} ϑ^{1.49} e^{- 33.33 ϑ}$ 。

图6 平均风速概率密度函数

Fig. 6 Density function of average wind speed

图7 湍流强度概率密度函数

Fig. 7 Density function of turbulence intensity

图8所示为平均风速与湍流强度之间的相关性。从图中可以看出，平均风速与湍流强度彼此之间存在明显的“上尾相关性”，符合Gumbel Copula函数特征，因此，根据式(3)和(4)可得平均风速-湍流强度联合概率密度函数，如图9所示。

图8 平均风速-湍流强度相关性

Fig. 8 Correlation between mean wind speed

图9 平均风速-湍流强度联合概率密度函数

Fig. 9 Joint function of mean wind speed

图10所示为平均风速、湍流强度分别对低速级太阳轮接触与弯曲应力均值与方差的影响。从图10(a)中可以看出，当平均风速低于额定风速时，低速级太阳轮接触与弯曲应力均值随平均风速增加而增大；当平均风速高于额定风速时，接触与弯曲应力均值基本稳定；低速级太阳轮接触应力均值远大于弯曲应力均值。从图10(b)中可以看出，接触应力方差在10~40 MPa波动，弯曲应力方差在0~25 MPa波动，总体上，低速级太阳轮接触应力方差大于弯曲应力方差。

图10 环境参数对低速级太阳轮接触与弯曲应力的影响

Fig. 10 Effects of environmental parameters on contact and bending stress of the low speed sun gear

图11所示为在额定风速(平均风速9.58 m/s，湍流强度0.11)作用下风电齿轮箱传动系统各级齿轮短期接触与弯曲疲劳损伤。可以看出，低速级太阳轮的短期接触疲劳损伤最大，同时各级齿轮短期接触疲劳损伤均大于弯曲疲劳损伤。

图11 风电齿轮箱传动系统齿轮短期接触与弯曲疲劳损伤

Fig. 11 Short-term contact and bending fatigue damage of gears in the wind turbine gearbox transmission system

图12所示为平均风速、湍流强度对低速级太阳轮短期接触与弯曲疲劳损伤的影响。从图12(a)中可知，当平均风速低于额定风速时，低速级太阳轮短期接触与弯曲疲劳损伤随平均风速增加而增大；当平均风速高于额定风速时，平均风速对其影响较小。其主要原因是当风速低于额定风速时，风电机组以追踪最佳风能利用系数运行，载荷变化明显；当风速高于额定风速时，由于变桨系统限制了风电机组功率超发，载荷较为稳定。从图12(b)中可以看出，当平均风速低于8.0 m/s时，低速级太阳轮短期接触与弯曲疲劳损伤主要受平均风速影响，湍流强度影响较小；当平均风速在8.0~9.5 m/s时，其损伤会随着湍流强度的增加而增大；当平均风速高于9.5 m/s时，由于控制策略作用，湍流强度对其损伤影响较小。

图12 环境参数对低速级太阳轮短期接触与弯曲疲劳损伤的影响

Fig. 12 Effects of environmental parameters on short-term contact and bending fatigue damage of the low-speed sun gear

3.2 环境参数对风电齿轮箱传动系统齿轮长期疲劳损伤影响

为了验证本文方法计算结果的精度与效率，将本文方法与常规方法计算的齿轮长期疲劳损伤进行对比分析，如式(15)所示。在常规计算方法中，在平均风速3~25 m/s、湍流强度0~0.25的区间内，以平均风速间隔0.75 m/s、湍流强度间隔0.05选取115个环境参数组合工况，直接根据概率累加计算风电齿轮箱传动系统齿轮长期疲劳损伤。

e = |\frac{D_{t o l}^{o l d} - D_{t o l}^{n e w}}{D_{t o l}^{n e w}}| \times 100 % ，

(15)

式中： $D_{t o l}^{o l d}$ 为常规方法计算的齿轮长期疲劳损伤； $D_{t o l}^{n e w}$ 为本文方法计算的齿轮长期疲劳损伤； $e$ 为误差百分比。

图13所示为本文方法与常规方法计算的风电齿轮箱传动系统齿轮长期疲劳损伤计算结果与误差对比。从图中可以看出，本文方法与常规方法计算的结果误差小于8%。整个仿真分析在一台CPU型号为i7-4790K，主频为4.00 GHz的台式计算机上完成，其中常规方法计算耗时约72 h，而本文方法计算仅耗时约20 h，效率提高260%。

图13 风电齿轮箱传动系统齿轮长期疲劳损伤计算结果与误差对比

Fig. 13 Long-term fatigue damage calculation results and errors of gears in the wind turbine gearbox transmission system

图14所示为平均风速、湍流强度对低速级太阳轮长期接触与弯曲疲劳损伤随机特性的影响。其中，黄色阴影部分表示当平均风速确定，湍流强度概率分布如图14(a)时，或当湍流强度确定，平均风速概率分布如图14(b)时，95%置信区间内的齿轮长期疲劳损伤。从图14(a)中可以看出，低速级太阳轮长期接触与弯曲疲劳损伤随着平均风速的增加而先增后减，同时由于湍流强度的随机性，尤其当平均风速位于10~11 m/s时齿轮长期疲劳损伤出现明显的随机不确定性。从图14(b)中可以看出，低速级太阳轮长期接触与弯曲疲劳损伤随湍流强度的增加也呈现出先增后减的趋势；当湍流强度在0.14附近时齿轮长期疲劳损伤最大，同时其随机不确定性也越明显，疲劳失效风险增大。

图14 环境参数对低速级太阳轮长期疲劳损伤的影响

Fig. 14 Effects of environmental parameters on long-term fatigue damage of the low-speed sun gear

4 结论

笔者考虑长期风速概率分布特征，建立了大功率海上风电机组OpenFAST-SIMPACK联合仿真模型，通过关键环境工况分析，计算对应的风电齿轮箱传动系统齿轮短期疲劳损伤，基于代理模型预测齿轮长期疲劳损伤，得出结论如下：

1) 通过关键环境工况选取与基于代理模型的风电齿轮箱传动系统齿轮长期疲劳损伤计算，相对于常规方法，其计算效率可以提高260%，且误差小于8%。

2) 风电齿轮箱传动系统各级齿轮接触疲劳损伤均高于弯曲疲劳损伤，其中低速级太阳轮接触疲劳损伤最大，易发生疲劳失效。

3) 在额定风速以下时，低速级太阳轮短期疲劳损伤与平均风速呈正相关，而在额定风速附近时，其主要受湍流强度影响；平均风速与湍流强度的随机特性会增大低速级太阳轮长期疲劳损伤不确定性。

参考文献

李垚, 朱才朝, 陶友传, 等. 风电机组可靠性研究现状与发展趋势[J]. 中国机械工程, 2017, 28(9): 1125-1133. [百度学术]

Li Y, Zhu C C, Tao Y C, et al. Research status and development tendency of wind turbine reliability[J]. China Mechanical Engineering, 2017, 28(9): 1125-1133. (in Chinese) [百度学术]

王磊. 海上风电机组系统动力学建模及仿真分析研究[D]. 重庆: 重庆大学, 2011. [百度学术]

Wang L. Study on systematic dynamic model and simulation for offshore wind turbine[D]. Chongqing: Chongqing University, 2011. (in Chinese) [百度学术]

Dong W B, Xing Y H, Moan T, et al. Time domain-based gear contact fatigue analysis of a wind turbine drivetrain under dynamic conditions[J]. International Journal of Fatigue, 2013, 48: 133-146. [百度学术]

向东, 蒋李, 沈银华, 等. 风电齿轮箱在随机风载下的疲劳损伤计算模型[J]. 振动与冲击, 2018, 37(11): 115-123. [百度学术]

Xiang D, Jiang L, Shen Y H, et al. Fatigue damage calculation model for wind turbine gearboxes under random wind loads[J]. Journal of Vibration and Shock, 2018, 37(11): 115-123. (in Chinese) [百度学术]

Nejad A R, Gao Z, Moan T. On long-term fatigue damage and reliability analysis of gears under wind loads in offshore wind turbine drivetrains[J]. International Journal of Fatigue, 2014, 61: 116-128. [百度学术]

Nejad A R, Bachynski E E, Kvittem M I, et al. Stochastic dynamic load effect and fatigue damage analysis of drivetrains in land-based and TLP, spar and semi-submersible floating wind turbines[J]. Marine Structures, 2015, 42: 137-153. [百度学术]

Wang S S, Nejad A R, Moan T. On design, modelling, and analysis of a 10-MW medium-speed drivetrain for offshore wind turbines[J]. Wind Energy, 2020, 23(4): 1099-1117. [百度学术]

熊中杰. 随机风速下风力发电机组齿轮箱疲劳断裂寿命研究[D]. 南京: 南京理工大学, 2019. [百度学术]

Xiong Z J. Research on fatigue fracture life of wind turbine gearbox under random wind[D]. Nanjing: Nanjing University of Science and Technology, 2019. (in Chinese) [百度学术]

Heidebrecht A, MacManus D G. Surrogate model of complex non-linear data for preliminary nacelle design[J]. Aerospace Science and Technology, 2019, 84: 399-411. [百度学术]

Palmer K, Realff M. Metamodeling approach to optimization of steady-state flowsheet simulations: model generation[J]. Chemical Engineering Research and Design, 2002, 80(7): 760-772. [百度学术]

Jia Z Y, Davis E, Muzzio F J, et al. Predictive modeling for pharmaceutical processes using kriging and response surface[J]. Journal of Pharmaceutical Innovation, 2009, 4(4): 174-186. [百度学术]

Milovanović S, von Sydow L. A high order method for pricing of financial derivatives using radial basis function generated finite differences[J]. Mathematics and Computers in Simulation, 2020, 174: 205-217. [百度学术]

Murcia J P, Réthoré P E, Dimitrov N, et al. Uncertainty propagation through an aeroelastic wind turbine model using polynomial surrogates[J]. Renewable Energy, 2018, 119: 910-922. [百度学术]

Li X A, Zhang W. Probabilistic fatigue evaluation of floating wind turbine using combination of surrogate model and copula model[C]//Proceedings of the AIAA Scitech 2019 Forum, San Diego, California. Reston, Virginia: AIAA, 2019: AIAA2019-0247. [百度学术]

Zhao Y L, Dong S. Probabilistic fatigue surrogate model of bimodal tension process for a semi-submersible platform[J]. Ocean Engineering, 2021, 220: 108501. [百度学术]

Wilkie D, Galasso C. Impact of climate-change scenarios on offshore wind turbine structural performance[J]. Renewable and Sustainable Energy Reviews, 2020, 134: 110323. [百度学术]

Jonkman B, Jonkman J. FAST v8. 16.00 a-bjj[EB/OL]. 2016-07-16[2021-12-20]. https://openfast.readthedocs.io/en/v3.5.2/_downloads/5f2ddf006568adc9b88d8118dc3f1732/FAST8_README.pdf. [百度学术]

陈旭, 朱才朝, 宋朝省, 等. 紧急停机工况下风力发电机系统动态特性分析[J]. 机械工程学报, 2019, 55(5): 82-88. [百度学术]

Chen X, Zhu C C, Song C S, et al. Dynamic characteristics analysis of wind turbine under emergency shutdown events[J]. Journal of Mechanical Engineering, 2019, 55(5): 82-88. (in Chinese) [百度学术]

陈岩松, 朱才朝, 谭建军, 等. 多工况下兆瓦级海上风电齿轮箱均载性能优化设计[J]. 重庆大学学报, 2022, 45(09): 1-14. [百度学术]

Chen Y S, Zhu C C, Tan J J, et al. Optimal design of load sharing performance of megawatt level offshore wind turbine gearbox under multi-operating conditions[J]. Journal of Chongqing University, 2022, 45(09): 1-14.(in Chinese) [百度学术]

刘华朝, 朱才朝, 柏厚义. 轮齿修形对兆瓦级风电齿轮箱NVH性能的影响[J]. 振动与冲击, 2016, 35(24): 158-163, 188. [百度学术]

Liu H C, Zhu C C, Bai H Y. The effect of gear modification on the NVH characteristics of a megawatt level wind turbine gearbox[J]. Journal of Vibration and Shock, 2016, 35(24): 158-163, 188. (in Chinese) [百度学术]

International Electrotechnical Commission. Wind turbines - Part 1: design requirements : IEC 61400-1:2005 [S]. IEC, 2006. [百度学术]

Li H X, Cho H, Sugiyama H, et al. Reliability-based design optimization of wind turbine drivetrain with integrated multibody gear dynamics simulation considering wind load uncertainty[J]. Structural and Multidisciplinary Optimization, 2017, 56(1): 183-201. [百度学术]

龚伟俊, 李为相, 张广明. 基于威布尔分布的风速概率分布参数估计方法[J]. 可再生能源, 2011, 29(6): 20-23. [百度学术]

Gong W J, Li W X, Zhang G M. The estimation algorithm on the probabilistic distribution parameters of wind speed based on Weibull distribution[J]. Renewable Energy Resources, 2011, 29(6): 20-23. (in Chinese) [百度学术]

Sklar A. Fonctions de repartition an dimensions et leurs marges[J]. Publ. inst. statist. univ. Paris, 1959, 8: 229-231. [百度学术]

涂志斌, 黄铭枫, 楼文娟, 等. 基于Copula函数的风浪多方向极限状态曲线[J]. 振动与冲击, 2021, 40(14): 1-9, 46. [百度学术]

Tu Z B, Huang M F, Lou W J, et al. Dimensional environmental contour lines of the wind and wave based on Copula functions[J]. Journal of Vibration and Shock, 2021, 40(14): 1-9, 46. (in Chinese) [百度学术]

侯亚楠. Copula函数的估计及其应用[D]. 武汉: 华中科技大学, 2013. [百度学术]

Hou Y N. The estimation and application of copula function[D]. Wuhan: Huazhong University of Science and Technology, 2013. (in Chinese) [百度学术]

Khowaja K, Shcherbatyy M, Härdle W K. Surrogate models for optimization of dynamical systems[EB/OL]. 2021: arXiv: 2101.10189. https://arxiv.org/abs/2101.10189.pdf. [百度学术]

Martini M, Guanche R, Armesto J A, et al. Met-ocean conditions influence on floating offshore wind farms power production[J]. Wind Energy, 2016, 19(3): 399-420. [百度学术]

周金宇, 谢里阳, 韩文钦, 等. 基于Nataf变换的载荷相关系统风险预测方法[J]. 机械工程学报, 2009, 45(8): 137-141. [百度学术]

Zhou J Y, Xie L Y, Han W Q,et al. Method for syetem risk prediction with load dependency based on nataf transformation[J]. Journal of Mechanical Engineering, 2009, 45(8): 137-141. (in Chinese) [百度学术]

张立波, 程浩忠, 曾平良, 等. 基于Nataf逆变换的概率潮流三点估计法[J]. 电工技术学报, 2016, 31(6): 187-194. [百度学术]

Zhang L B, Cheng H Z, Zeng P L, et al. A three-point estimate method for solving probabilistic load flow based on inverse Nataf transformation[J]. Transactions of China Electrotechnical Society, 2016, 31(6): 187-194. (in Chinese) [百度学术]

International Organization for Standardization. Calculation of load capacity of spur and helical gears - Part 6: calculation of service life under variable load: ISO 6336-6:2006 [S]. Switzerland: International Organization for Standardization, 2006． [百度学术]

International Organization for Standardization. Calculation of load capacity of spur and helical gears - Part 2: calculation of surface durability (pitting): ISO 6336-2:2006 [S]. Switzerland: International Organization for Standardization, 2006． [百度学术]

International Organization for Standardization. Calculation of load capacity of spur and helical gears - Part 3: calculation of tooth bending strength: ISO 6336-3:2006 [S]. Switzerland: International Organization for Standardization, 2006. [百度学术]

Miner M A. Cumulative damage in fatigue[J]. Journal of Applied Mechanics, 1945, 12(3): A159-A164. [百度学术]

International Organization for Standardization. Calculation of load capacity of spur and helical gears - Part 5: strength and quality of materials: ISO 6336-5:2003 [S]. Switzerland: International Organization for Standardization, 2003． [百度学术]

龙腾, 刘建, Wang G G, 等. 基于计算试验设计与代理模型的飞行器近似优化策略探讨[J]. 机械工程学报, 2016, 52(14): 79-105. [百度学术]

Long T, Liu J, Wang G G, et al. Discuss on approximate optimization strategies using design of computer experiments and metamodels for flight vehicle design[J]. Journal of Mechanical Engineering, 2016, 52(14): 79-105. (in Chinese) [百度学术]

Okpokparoro S, Sriramula S. Uncertainty modeling in reliability analysis of floating wind turbine support structures[J]. Renewable Energy, 2021, 165: 88-108. [百度学术]

Rashki M, Azarkish H, Rostamian M, et al. Classification correction of polynomial response surface methods for accurate reliability estimation[J]. Structural Safety, 2019, 81: 101869. [百度学术]

Wang H, Li W B, Qian Z H, et al. Reconstruction of wind pressure fields on cooling towers by radial basis function and comparisons with other methods[J]. Journal of Wind Engineering and Industrial Aerodynamics, 2021, 208: 104450. [百度学术]

Bhosekar A, Ierapetritou M. Advances in surrogate based modeling, feasibility analysis, and optimization: a review[J]. Computers & Chemical Engineering, 2018, 108: 250-267. [百度学术]

Home