時變計數法模型及其驗證

范文亮，魏剛毅，李正良，2

(1.重慶大學 土木工程學院，重慶 400045；2.重慶大學 山地城鎮建設與新技術教育部重點實驗室，重慶 400045)

疲勞失效是航空工程、機械工程、船舶工程、交通運輸工程和土木工程等諸多領域的普遍物理現象，由此造成的災害性事故屢見不鮮。如飛機、核壓力容器和橋梁等都存在嚴重的疲勞破壞風險[1]。眾所周知，降低結構或構件疲勞失效危險的有效措施就是對其進行疲勞損傷分析，進而實現抗疲勞設計。

常規疲勞分析中，通常需要對所關注部位的應力時間歷程進行統計，將此統計結果與疲勞累積損傷相結合即可獲得結構或構件的損傷狀態估計，由此可以進一步進行結構剩余壽命評估。其中，對應力時間歷程進行統計進行統計的過程就稱為計數，與之對應的方法即為計數法，是疲勞分析中非常關鍵的環節之一。

計數法種類繁多，主要包括水平跨越計數法(Level-Crossing Counting)、峰值計數法(Peak Counting)、單變程計數法(Simple-Range Counting)和雨流類計數方法(Rainflow Counting and Related Methods)等[2]。自1968年Matsuiski和Endo提出雨流計數法[3]以來，雨流類計數方法獲得了廣泛的關注，提出了三峰谷值法[4-6]和四峰谷值法[7-8]兩大類方法，并在此基礎上發展了很多改進措施[9-14]，目前已成為了應用最為廣泛的計數法。仔細考察上述各類雨流計數法不難發現，盡管早期雨流計數法的各類缺陷逐漸被新近方法所改進，但是它們均存在一個共同的特點：欲獲得某一時間歷程的計數結果，必須從該歷程的起始點開始運算。

隨著可持續發展理念和精細化設計思路的重要性日益凸顯，工程結構的全壽命周期設計近十年來引起了廣泛的興趣。顯然，結構性能的時變變化規律和特性是全壽命周期設計的基礎。對應于遭遇疲勞破壞風險的結構，結構的時變疲勞損傷分析就至關重要。換言之，需要了解結構疲勞損傷隨著時間歷程的變化。若采用傳統的計數方法，任一時刻疲勞損傷的計算都需要從疲勞荷載剛作用于結構時開始計數，從而導致了很多不必要的重復運算，降低了計算效率。當需要由構件的疲勞損傷評估整體結構的疲勞損傷時，效率的下降就更為明顯，對于復雜工程結構系統尤甚。

為了滿足全壽命周期分析和設計的需要，借鑒動力分析中的增量方法思路，本文力圖發展一種可以充分利用已有歷程計數結果，再結合增量時程計數，從而實現任意時刻計數的時變計數法，并對其有效性和高效性進行驗證。

1 時變計數法模型及實現

顯而易見，對于一般的時間歷程，若將其分為兩段，且采用雨流計數法對整體和兩分段分別計數，那么整體計數結果與兩分段計數結果之和并不一定相等。以圖1所示應力歷程簡單說明之。將該歷程在圖中虛線處分段，表1分別給出了圖1(a)的計數結果和圖1(b)中兩段時程計數結果之和，兩者存在顯著差異。

表1 不分段與分段累加計數結果的比較

圖1 歷程不分段和分段兩種情況

1.1 時變計數法模型

為解決上述問題，文中給出了時變計數法模型。詳細過程如下：

(1) 可根據分析要求對時程進行分段。

(4) 將第一段殘波Rt1的后三個點接入第二分段ΔSt2的首部，形成新的第二分段ΔSt2’。若Rt1的點數不足三個，就Rt1將整體接入ΔSt2。

1.2 時變計數法的實現

根據時變計數法模型的闡述可以發現，除獲得t1時刻的計數結果稍有差異外，其他時刻的計數結果的計算過程是完全相同的。因此，時變計數法的實現亦可分為兩部分：首先是獲得t1時刻的計數結果及其中間結果，然后逐步讀入增量歷程并采用建議方法獲得其他時刻的計數結果。時變計數法的詳細流程圖見圖2。

圖2 時變計數法的實現流程圖

2 時變計數法的計算性能分析

為了避免計算機運行狀態和程序算法結構對各計數法運算效率的影響，文中將對時變計數法與Amzallag法所包含的各類操作數量進行簡單的理論分析，并通過對比說明建議方法的高效性。

無論是建議的時變計數法和Amzallag法，計算過程均包含如下操作：數據檢測和壓縮、對給定時長的時程數據進行計數、殘波復制和移位處理。假設上述各操作與對象的時間跨度呈正比，且單位時長的數據檢測和壓縮耗時為a，單位時長時程計數耗時為b，單位時長殘波復制和移位耗時為c，且假設單位時長數據計數所剩殘波的比例為d。

若存在一時長為[0,tn]的應力時程，需要獲得在時刻t1，t2，…，tn的結構損傷。為此，必須先得到各時刻的計數結果?；谏鲜龊唵渭僭O，考察分別采用時變計數法和Amzallag法所需的理論計算時間。

不難發現，采用Amzallag法對單位時長歷程計數耗時為a+b+c+bd，于是所有n個時刻計數結果的耗時為(t1+t2+…+tn)×(a+b+c+bd)。

若采用時變計數法，則獲得t1時刻的計數結果需耗時t1×(a+b+c+bd)；獲得ti時刻的計數結果需耗時為(ti-ti-1)×(a+b)+Ai(a+b)+dAi(a+b+c)，其中Ai=di-1t2+di-2(t3-t2) +…+d(ti-ti-1)，總耗時即為各時刻耗時之和。

為便于比較，不妨再假設各時刻均勻分布，且間隔時長為τ，那么采用Amzallag法所需總耗時為n(n+1)τ/2×[a+(1+d)×b+c]；而采用時變計數法所需耗時為 [nτ+(1+d)A0]×(a+b)+(dA0+τ)c+dτb，其中A0=(n-1)dτ/(1-d)+[(1-dn-1)(1-2d)dτ]/(1-d)2。通常，d遠小于1，隨著n的增大，dn-1逐漸趨近于0，于是A0≈(n-1)dτ/(1-d)+(1-2d)dτ/(1-d)2。顯然，Amzallag法的總耗時隨n呈非線性增長，而建議方法耗時基本與n呈正比，且d越小，比例系數越小。因此，隨著n的增大時變計數法的高效性逐漸顯著。以n=100，d=0.5為例，Amzallag法的總耗時為5 050τ(a+1.5b+c)，時變計數法的總耗時為137.125τa+137.625τb+13.375τc，不足前者的1/37。若n=10 000，d=0.5，則建議方法耗時不足Amzallag法的1/3 600，高效性較n=100情形更為明顯。

3 時變計數法模型的算法驗證

文中針對兩類簡單算例，將建議方法與常用的雨流類計數法展開了詳細的對比分析，以驗證時變計數法的準確性。此處，常用雨流類方法選取兩種，分別為Bannantine等提出的三峰谷值法(簡稱為Bannantine法)與Amzallag法。

3.1 隨機應力時程

為簡便，文中構造了一類虛擬的隨機應力歷程，其構造過程如下：首先生成m個服從N(0,1)分布的隨機數據點，并將其視為有序數列，然后將該數列中的各元素依次交替乘以“-1”和“1”得到新的有序數列，最后以該有序數列為縱坐標(即應力值)、以從0開始的等間隔時刻為橫坐標得到隨機應力歷程。

對于構造的變幅應力歷程，可將其分為n段，由時變計數法給出n個時刻對應的計數結果，同時亦可由Bannantine法和Amzallag法給出最終時刻的計數結果。僅從驗證時變計數法準確性的角度出發，只需給出最終時刻的計數結果并將其與常用計數法比較即可。

為了驗證建議方法的普適性，文中生成了3條隨機應力時程。采用上述各計數法對此3條歷程進行計數的結果分別示于圖3～圖5。圖中，x軸為滯回環的幅值，y軸為其均值，z軸為各幅值-均值對出現的次數，下同。

圖3 隨機應力時程1的計數結果(m=25 000)

圖4 隨機應力時程2的計數結果(m=50 000)

圖5 隨機應力時程3的計數結果(m=75 000)

不難發現，對于任意構造的3個隨機應力時程，時變計數法的結果均與常用雨流類計數法的結果吻合，有效地驗證了建議模型的準確性和合理性。此外，由上述3個算例亦可發現Bannantine法和Amzallag法的計數結果也是一致的，進一步驗證了三峰谷值法與四峰谷值法的等效性，與McInnes等[15]的理論分析結果相互佐證。

值得指出的是，盡管各方法得到的幅值-均值對循環的出現頻次完全相同，但各方法獲得的幅值-均值對序列的次序不同。

3.2 一般應力時程

為了驗證建議方法對更一般時間歷程的適用性，文中截取了某輸電塔關鍵構件在風致振動下的應力時程為計數對象。該時程共包含3 004個峰谷點，時間間隔0.25 s，總時長為751 s。采用時變計數法時，將其分為3段，即1～1 000、1 001～2 000和2 001～3 004這三段。由于Bannantine法和Amzallag法的等效性，此處僅給出Amzallag法和建議方法得到的最終計算結果，如圖6所示。經比較，建議方法的準確性得到了進一步驗證。

圖6 輸電塔應力時程的最終計數結果

圖7 輸電塔應力時程的時變計數結果

此外，在采用時變計數法獲得最終計數結果的過程中可以很方便地給出前兩段增量歷程終點時刻的計數結果，如圖7所示。

4 結 論

為了提高時變疲勞損傷分析的效率，文中借鑒動力分析中增量法的思路，提出了一種可以逐步讀入增量時程進行計數的時變計數法。與傳統計數法相比，該方法不僅可以獲得最終時刻的計數結果，而且可以給出各增量時程終點時刻的計算效果。由于不必每次均從零時刻開始計算，時變計數法顯著地改善了傳統計數法在時變分析中的效率，且所關注的時點越多，效率提高就越為明顯。最終，文中通過一個虛擬算例和一個實際工程算例驗證了建議計數法的準確性。

此外，文中通過算例的形式驗證了三峰谷值法與四峰谷值法的等效性，與McInnes等的理論分析結果是相吻合的。

[1] 張行，趙軍.金屬構件應用疲勞損傷力學[M].北京: 國防工業出版社，1998，1-2.

[2] ASTM E1049-85.Standard practices for cycle counting in fatigue analysis [S].West Conshohocken，Pennsylvania: ASTM International，2005.

[3] Matsuishi M，Endo T.Fatigue of metals subjected to varying stress [C].In Proceeding of the Kyushu Branch of Japan Society of Mechanics Engineering，Fukuoka，Japan (in Japanese)，1968，37-40.

[4] Downing S，Socie D F.Simple rainflow counting algorithms [J].International Journal of Fatigue，1982，4(1):31-40.

[5] Glinka G，Kam J C P.Rainflow counting algorithm for very long stress histories [J].International Journal of Fatigue，1987，9(3):223-228.

[6] Bannantine J A，et al.Fundamentals of metal fatigue analysis [M].New Jersey: Prentice-Hall，1990，193-196.

[7] Amzallag C，et al.Standardization of the rainflow counting method for fatigue analysis [J].International Journal of Fatigue，1994，16:287-293.

[8] 蔣薈.實時雨流計數法的“三變程”計數原則[J].航空計算技術，2008,38 (5):5-7.

JJIANG Hui.Three-range rule of real-time rain-flow counting method [J].Aeronautical Computing Technique，2008，38(5):5-7.

[9] Rychlik I.A new definition of the rainflow cycle counting method [J].International Journal of Fatigue，1987，9(2):119-121.

[10] Hong N.A modified rainflow counting method [J].International Journal of Fatigue，1991，13(6):465-469.

[11] Anthes R J.Modified rainflow counting keeping the load sequence [J].International Journal of Fatigue，1997,19(7):529-537.

[12] 朱亦鋼.一種隨機變化載荷的疲勞損傷累計方法[J].機械強度，2004，26(S):32-35.

ZHU Yi-gang.Damage counting algorithm for random load fatigue analysis [J].Journal of Mechanical Strength，2004，26(S):32-35.

[13] 田軍，李強.改進的雨流法實時計數模型[J].北京交通大學學報，2009,33(1):28-31.

TIAN Jun，LI Qiang.Improved model of rain-flow real-time counting method [J].Journal of Beijing Jiao Tong University，2009，33(1):28-31.

[14] 蔣東方.雨流計數的遞歸算法[J].航空學報，2009，30(1):99-103.

JIANG Dong-fang.Recursive algorithm for rain-flow counting [J].Acta Aeronautica et Astronautica Sinica，2009，30(1):99-103.

[15] McInnes C H，et al，Equivalence of four-point and three-point rainflow cycle counting algorithms [J].International Journal of Fatigue，2008，30:547-559.