基于改進T-B法的寬帶高斯過程疲勞損傷分析

袁奎霖，孫卓成

（大連理工大學工業裝備結構分析國家重點實驗室船舶工程學院，遼寧 大連 116024）

0 引 言

船舶與海洋工程結構物服役于惡劣的海洋環境中，長期遭受如風載荷、波浪載荷、海流載荷等多種隨機交變載荷作用，極易產生疲勞損傷。在初期設計階段，通常采用基于功率譜的頻域方法計算隨機應力響應造成的疲勞損傷。當結構應力是一個零均值窄帶高斯隨機過程時，其應力范圍服從Ray?leigh 分布，Bendat 等[1]根據這種特性推導了在頻域內的窄帶疲勞損傷解析解。然而，實際情況下船海結構物在多種外部載荷和自身彈性變形影響下，結構應力往往是一個寬帶高斯隨機過程，其應力范圍的概率分布難以確定。工程上往往仍采用基于窄帶假設的Rayleigh 分布模型計算寬帶疲勞損傷，但當帶寬系數變大時將會高估疲勞損傷。因此，研究寬帶高斯隨機過程下疲勞損傷的頻域評估方法具有一定的理論意義和實際應用價值。

目前，針對寬帶高斯應力造成的疲勞損傷問題，學者們提出了一系列近似方法。Wirsching 和Light[2]假設寬帶應力過程的應力范圍仍服從Rayleigh 分布，在此基礎上引入一個可以考慮帶寬影響的修正系數對窄帶近似疲勞損傷解進行修正。W-L 修正系數已經在船舶與海洋工程結構物疲勞規范[3-4]中得到廣泛應用。由于寬帶高斯應力的雨流幅值概率密度函數難以從應力響應功率譜直接推導，也有學者試圖借助時域模擬仿真建立雨流幅值的經驗概率密度函數，以此建立更為可靠的寬帶疲勞損傷預報模型。Dirlik[5]采用不同參數形式的功率譜密度函數，通過蒙特卡羅時域模擬仿真對雨流計數后的應力范圍進行統計分析，提出了一個雨流幅值的經驗概率密度函數。Dirlik模型由一個指數分布和兩個Rayleigh 分布組成描述雨流幅值分布，由此推導出的疲勞損傷解析解的計算精度較W-L修正系數法有明顯改善。Zhao 和Baker[6]提出了由一個Weibull 分布和一個Rayleigh 分布線性組合而成的雨流幅值概率密度的近似模型，但該模型對于帶寬系數較大的功率譜而言，計算效果仍有待改善。Tovo和Benasciutti[7-8]提出了T-B 法即將雨流計數損傷表示為由一個權重系數bTB控制的范圍計數損傷和窄帶假設損傷的線性組合，并通過大量數值實驗得到了權重系數bTB的近似公式。已有研究[9-10]表明，T-B法具有適用范圍廣泛且計算精度較高的特點，對于不同譜型都具有較好的魯棒性。

然而，T-B 法中參數bTB的近似公式僅與譜矩有關，并沒有考慮S-N曲線斜率參數m的影響，因而隨著斜率參數m的增大，T-B 法與雨流法的差異亦逐漸增大。針對該問題，本文在原有T-B 法基礎上，在權重系數中引入S-N曲線斜率參數m加以修正，并通過大量數值實驗結果擬合得到一個新的權重系數bMTB的非線性函數模型。在此基礎上提出一種基于改進T-B法的寬帶高斯疲勞損傷分析方法，并對該方法的有效性進行驗證。

1 寬帶高斯隨機過程疲勞損傷分析方法

1.1 譜矩和帶寬系數

對于單邊功率譜密度函數為SXX(ω)的平穩高斯隨機過程X(t) 而言，其譜矩定義為

式中，ω為角頻率，單位為rad/s。對于高斯過程，其平均跨零率v0和平均峰值率vp可由譜矩表示：

此外，譜矩可以對隨機過程的帶寬進行表征，即帶寬系數可表示為

其中，最常用的兩個帶寬系數α1和α2定義如下：

此外，工程上常用的與α1和α2有關的另外兩種帶寬系數定義如下：

式中，δ稱為Vanmarcke 帶寬系數[11]，取值范圍為0 ?δ?1；ε稱為Wirsching 帶寬系數[2]，取值范圍為0 ?δ?1。當α1和α2越趨近于0，隨之δ和ε越趨近于1 時，該隨機過程越趨近于理想寬帶隨機過程，反之則越趨近于窄帶隨機過程。工程上一般認為δ大于0.1時，可將一個隨機過程視作寬帶隨機過程。

1.2 高斯隨機疲勞損傷分析

目前工程上的疲勞損傷分析主要是基于S-N曲線方法進行的。材料的S-N曲線表示了該材料在恒幅循環載荷作用下應力范圍與疲勞失效循環次數之間的關系，其一般表達式如下：

式中，S代表應力范圍，N代表在特定的應力水平作用下疲勞失效時的循環次數，m和K分別是材料參數。

當結構遭受多級應力載荷時，工程上常采用Miner線性累積損傷理論計算總疲勞損傷：

式中，ni表示在應力范圍Si作用下的循環次數，Ni表示在相同應力范圍Si作用下疲勞失效時的循環次數，n代表載荷作用時長T內的總循環次數為Sm的數學期望。

當結構應力是一個隨時間變化的隨機過程時，其雨流循環的概率密度函數為fs(S)，則在隨機響應作用時間T下的疲勞損傷可表示為

當隨機過程X(t)是一個嚴格的零均值窄帶高斯過程時，其峰值和谷值在每一個循環中對稱出現，由此可以認為其雨流幅值分布與其峰值分布均服從Rayleigh分布：

其應力范圍S是應力幅值R的2倍，也服從Rayleigh 分布，并且平均峰值率vp等于平均跨零率v0。由式（8）可得到作用時間T范圍內的窄帶疲勞損傷解析解為

式中，Γ(?)表示gamma函數。目前,工程中常用Rayleigh分布模型來近似代替寬帶過程，將一個零階譜矩和平均跨零率與原始隨機過程相等的窄帶過程來近似代替實際寬帶過程，即窄帶近似方法。

1.3 寬帶高斯疲勞損傷分析方法

對于寬帶高斯隨機過程，學者們根據不同功率譜生成大量隨機時間歷程，力爭頻域結果與雨流計數時域結果的誤差最小化，提出了多種寬帶高斯疲勞損傷的計算公式。

1.3.1 Wirsching-Light方法（W-L）

Wirsching 和Light[2]基于等效近似窄帶過程的概念，提出了一種估算寬帶高斯疲勞損傷的經驗模型：

式中，DWL為W-L法估算的寬帶高斯疲勞損傷，DNB為式（10）計算的窄帶高斯疲勞損傷，ρWL是W-L法的修正因子，它是關于帶寬系數ε和S-N曲線斜率參數m的函數，即

式中：

由上式可知，W-L方法認為寬帶高斯過程的疲勞損傷僅與λ0、λ2和λ4三個譜矩參數和單斜率S-N曲線斜率參數m有關。

1.3.2 Dirlik方法（DK）

Dirlik[5]提出了雨流幅值分布的半經驗公式，該分布由一個指數型分布和兩個Rayleigh 型分布組成：

式中：

進而可以得到寬帶高斯過程的疲勞損傷為

由上式可知，DK 方法的雨流幅值分布僅與兩個帶寬參數α1和α2有關，而α1和α2又與λ0、λ1、λ2和λ4四個譜矩參數有關。

1.3.3 Zhao-Baker方法（Z-B）

Zhao和Baker[6]提出了一個由Weibull分布和Rayleigh分布線性組合的雨流幅值分布，公式如下：

式中：

由式（8）和式（17）得到疲勞損傷公式如下：

1.3.4 Tovo-Benasciutti方法（T-B）

Tovo 和Benasciutti[7-8]通過理論分析發現，對于一個穩態高斯隨機過程，其雨流計數疲勞損傷DRFC總是介于基于范圍計數的疲勞損傷DRC和窄帶假設疲勞損傷DNB之間，即

式中，DNB可根據公式（10）計算，DRC的近似解如下：

基于這種思想，Tovo 和Benasciutti 提出雨流計數損傷可表示為由一個權重系數bTB控制的范圍計數損傷和窄帶假設損傷的線性組合，即T-B模型：

為了正確評估DRFC需要確定權重系數bTB，Tovo 和Benasciutti 經過大量數值實驗得到bTB的近似公式如下：

由此可知，T-B 法的權重參數bTB僅與帶寬系數α1和α2即λ0、λ1、λ2和λ4四個譜矩參數有關，與S-N曲線斜率參數m無關。需要說明的是，Tovo和Benasciutti[7-8]通過數值實驗確定bTB時將S-N曲線斜率參數固定為m=3。然而，已有研究[13-16]表明，隨著S-N曲線斜率參數m的增大，T-B法計算的疲勞損傷與雨流法結果的誤差明顯增大。如圖1所示，本文在前期研究中分析了m=3和m=5時，不同帶寬（α1和α2）條件下由式（23）反推的bsim值與bTB近似公式即式（24）的差別。由圖1 可知，當帶寬系數α2=0.1時，bsim值受m值影響較小且與bTB近似公式吻合較好；隨著α2的增大，m=3和m=5對應的bsim值之間的差別變得更加明顯，表明T-B法權重參數bTB的近似公式應該考慮m的影響。

圖1 權重系數的模擬結果bsim與擬合公式bTB的對比Fig.1 Comparison of the weighting factor bsim from numeri?cal simulation and approximate formula of bTB

2 基于改進T-B法的寬帶高斯疲勞損傷分析方法

2.1 參數化功率譜

針對T-B法的權重系數提出新的擬合公式，本文采用逆傅里葉變換技術，根據不同參數形式的功率譜模擬生成大量的時間序列。采用時域的雨流計數法計算疲勞損傷，考慮不同的S-N曲線斜率參數m的影響，代入式（23）中反推出新的權重系數。為了能夠得到更大范圍的帶寬參數，同時考慮到不同譜形對疲勞損傷的影響，本文采用了如圖2所示的六種具有不同形狀的功率譜。

圖2 不同形狀的功率譜Fig.2 Illustration of spectral shapes considered in this study

圖2為六種參數化功率譜的示意圖。以圖2（a）中的對稱二次型功率譜為例，其表達式如下：

圖2（b）～（f）分別為反對稱二次型、線性以及常數型功率譜，表達式可參照式（25）獲得。圖2（a）～（e）五種譜型的ω1和ω3為固定值，分別取2π/1000 rad/s和2π rad/s，而ω2介于ω1和ω3之間，通過改變ω2、h1和h2的值可以得到不同形狀的功率譜。為了使得所有譜型的零階譜矩λ0為常數，各參數有如下關系：

其中，圖2（a）～（c）的M=3，圖2（d）的M=2，圖2（e）的M=1。

此外，圖2（f）為分離式矩形雙模態譜，常被用于高斯雙模態疲勞損傷分析[12-13]的數值模擬中，其形狀參數定義如下:

式中，B為高頻模態與低頻模態功率譜下的面積（能量）之比，R為兩個模態的中心頻率之比。為了保證兩個矩形譜是分離的，R還需要滿足如下條件：

本文中ωa取5 rad/s，B=0.01～10，R=3～20，A1+A2=1。c1和c2決定了低頻和高頻模態的帶寬特性，本文中c1和c2都取為1.1。

2.2 改進的權重系數bMTB

如式（24）所示，原T-B 法的權重系數bTB的近似公式僅是與帶寬系數α1和α2相關的函數。本文通過數值模擬方法重新建立權重參數b與帶寬系數α1、α2以及S-N曲線斜率參數m的函數關系。船舶與海洋工程結構為鋼質焊接結構，其S-N曲線一般為以N=107為拐點的雙斜率曲線。但是，為了便于研究，本文將采用單斜率的S-N曲線進行疲勞損傷計算。具體方法為選取不同的α1、α2和m值，利用2.1節中的六種譜型生成大量隨機時間序列進行數值模擬試驗，利用雨流計數法計算疲勞損傷DsRimFC。由式（23）可得到通過數值模擬確定權重系數bsim的表達式為

式中，DNB和DRC可利用式（10）和式（22）直接計算。由此，生成了一系列與不同α1、α2和m對應的bsim值，通過非線性擬合技術確定新的權重系數函數模型bMTB：

式中，系數ai為隨不同S-N曲線斜率參數m變化的待定系數。對于船舶與海洋工程領域常用的S-N曲線[4]，斜率參數m的范圍為3～6，相應的擬合系數ai的具體值見表1。圖3 以m=3 和m=5 為例，顯示了數值模擬得到的權重系數bsim和兩種擬合函數模型bTB和bMTB隨α1和α2的變化規律。可以看出，新的權重系數模型bMTB與原T-B 法的權重系數模型bTB相比，與數值模擬結果的擬合效果更好，且能反映不同斜率參數m的影響。

表1 擬合系數ai(i = 1,2,3,4)Tab.1 Fitted coefficients ai(i = 1,2,3,4)

圖3 權重系數b的數值模擬結果與擬合公式結果對比Fig.3 Comparison of the weighting factors of b from numerical simulation and empirical formulae

2.3 不同寬帶高斯疲勞損傷評估方法的比較

為了將改進T-B 法與現有的寬帶高斯疲勞損傷評估方法進行對比，本文分別采用決定系數R2和Tovo-Benasciutti[8]定義的誤差指標EI作為衡量指標來評估各方法的準確性和魯棒性。R2和EI的定義如下：

針對2.1 節中的六種譜型，S-N曲線斜率參數分別取m=3 和m=5 進行數值模擬試驗，將改進T-B法與W-L 法[2]、DK 法[5]、Z-B 法[6]和T-B 法[7-8]的頻域疲勞損傷結果與相應的時域雨流結果進行對比，各方法的決定系數R2和誤差指標EI如表2和表3所示。

表2 不同頻域方法的決定系數R2Tab.2 Determination coefficients R2 for different frequency-domain methods

表3 不同頻域方法的誤差指標EITab.3 Error indexes EI for different frequency-domain methods

由表2 可知，對于六種參數化功率譜的數值模擬結果而言，當m=3 時，改進T-B 法的決定系數R2的平均值為0.9996，與其他四種方法相比更接近于1；當m=5時，五種方法的決定系數R2均有不同程度的降低，而改進T-B 法的決定系數R2為0.9963，表明當斜率參數m增大時該方法仍能保持良好的計算精度。

由表3可知，當采用誤差指標EI作為五種方法計算精度衡量指標時，改進T-B法的計算結果最為準確，其次為T-B法與DK法，且兩者精度相當，而W-L法和Z-B法的誤差最大。

為了更加直觀地將改進T-B法與DK法和T-B法進行對比，采用式（36）定義的指標γ，分析三種方法的計算精度隨Vanmarcke帶寬系數（0.1 ?δ?0.95）的變化規律，如圖4所示。

圖4 DK法、T-B法與改進T-B法疲勞損傷評估結果對比Fig.4 Comparison of fatigue damage evaluated by DK,T-B and modified T-B methods

由圖4 可知，對于m=3，DK 法、T-B 法和改進T-B 法三種方法得到的疲勞損傷與雨流損傷的比值均在1.0 附近，T-B 法的相對誤差范圍為-5%到5%，改進T-B 法的相對誤差范圍為-4%到4%；DK 法的損傷結果通常小于真實雨流損傷，并且當Vanmarcke 帶寬系數大于0.7 時，最大誤差超過-10%。由于疲勞損傷與應力循環的非線性關系會隨著m的增大而增強，當m=5時，三種方法的誤差都隨著Van?marcke 帶寬系數的遞增而明顯增大，改進T-B 法的最大誤差為-21.19%，但仍然小于T-B 法的最大誤差-25.25%和DK 法的最大誤差-33.6%。結果表明，在不同參數化譜型、帶寬范圍以及S-N曲線斜率參數的條件下，本文所提出的新方法的計算精度較已有的寬帶高斯疲勞損傷方法有了明顯改善。

3 算例驗證

本章選取船舶與海洋工程結構常遭受的一些典型寬帶響應譜，對本文所提出的改進T-B 法的有效性進行進一步驗證。

3.1 雙模態功率譜

對于船海結構物而言，其結構應力響應在功率譜上常常呈現出兩個明顯的峰值[14]，即高斯雙模態特征，因此高斯雙模態隨機過程是高斯寬帶過程的一種特例。在此，討論一個貼近于實際海洋結構物響應的雙模態功率譜[15]，其表達式如下：

式中，A是用于調整譜方差為1 的比例因子，TW=2π/ωw表示海浪周期，ωN=2π/TN是結構物的一階固有頻率，ξ為阻尼比。圖5 為雙模態功率譜算例，其中TW取9 s，結構振動的一階固有周期TN為2 s，阻尼比取0.04，該功率譜的帶寬系數δ= 0.5496。根據圖5 中的功率譜，模擬生成高斯隨機過程時域信號，由雨流法得到的疲勞損傷作為參考值，對五種頻域方法的計算精度進行對比，如表4所示。表中的相對誤差定義為

表4 各方法計算結果相對于雨流法的相對誤差（%，雙模態功率譜時）Tab.4 Relative errors of different methods compared with RFC method(%)

圖5 雙模態功率譜數值模擬算例（TW=9 s，TN=2 s，x=0.04）Fig.5 Bimodal spectrum used in numerical simulation in case of TW=9 s,TN=2 s,x=0.04

從表4 中可以看出，在m=3～6 時，W-L 法、DK 法、Z-B 法和T-B 法的誤差隨著m的增大而變大；當m=6 時，W-L 法的誤差為-16.61%，DK 法為-7.01%，Z-B 法為-7.99%，T-B 法為-13.63%。本文提出的改進T-B法在m=3～6時都較為準確，所有的誤差都未超過4%。

3.2 三模態功率譜

結構響應功率譜上存在三個明顯峰值的三模態隨機過程在海洋工程領域也是一種典型的寬帶高斯隨機過程[16]。圖6 為一個Spar 型漂浮式風機塔柱的彎矩響應譜[17]，三個模態的特征頻率分別為0.18 rad/s、0.48 rad/s 和2.4 rad/s，分別對應于縱搖固有頻率、波浪頻率和塔柱振動的一階固有頻率，該功率譜的帶寬系數δ= 0.6994。根據圖6 中的功率譜，模擬生成高斯隨機過程時域信號，由雨流法得到的疲勞損傷作為參考值，對五種頻域方法的計算精度進行對比，如表5 所示。可以看出在該工況下，相對于雨流法結果，W-L法和Z-B 法的結果更加保守，且相對誤差隨著S-N曲線斜率參數m的增大而降低，最大誤差分別為18.96%和14.40%；相比之下，DK 法、T-B 法和改進T-B法低估了疲勞損傷，且誤差隨S-N曲線斜率參數m的增大而增大，最大誤差分別為-19.33%、-21.71%和-12.16%。

表5 各方法計算結果相對于雨流法的相對誤差（%，三模態功率譜時）Tab.5 Relative errors of different methods compared with RFC method(%)

圖6 漂浮式Spar型風機三模態彎矩響應譜[17]Fig.6 Tri-modal bending moment spectrum of a floating Spar wind turbine[17]

以上的驗證算例表明，相對于已有的寬帶高斯疲勞損傷評估方法，本文所提出的改進T-B法對于船海結構物領域常見的真實響應功率譜具有更好的計算精度。需要說明的是，采用W-L 方法和Z-B方法計算本節中的雙模態和三模態功率譜算例時，隨著Vanmarcke 帶寬系數的不同，會存在高估或低估疲勞損傷的情況，在以往研究[9]中也有類似結論。

4 結 論

針對船舶與海洋工程領域中常見的寬帶高斯隨機過程，本文在原有T-B 法的基礎上考慮了S-N曲線斜率參數m的影響，建立了新的權重系數函數模型，以此提出了一種基于改進T-B法的寬帶高斯疲勞損傷頻域分析方法，并利用不同的參數化功率譜和真實響應譜對所提出方法的有效性進行了驗證，得到了以下結論：

（1）相比于原T-B法，本文提出新的權重系數函數模型bMTB能夠更為準確地反映帶寬參數a1、a2以及S-N曲線斜率參數m對于雨流損傷的影響，從而可獲得更為準確的疲勞損傷預報結果。

（2）隨著斜率參數m的增加，現有的通用寬帶譜方法計算的疲勞損傷相對于時域雨流損傷的相對誤差逐漸增大。改進T-B 法通過考慮斜率參數m的影響，能夠更加準確地評估在不同m值條件下的疲勞損傷。

（3）值得說明的是，在各船級社規范中鋼材S-N曲線的斜率參數取值主要為m=3。經過大量數值驗證表明，在m=3時改進T-B法計算的疲勞損傷與時域雨流法之間的相對誤差在±4%的范圍內，表明該方法在實際工程應用中具有較好的應用價值。