基于混沌理論的地震誘導近海風力機動力響應研究

牛凱倫， 閆陽天， 李 春,2， 薛世成， 李志昊， 孫 康

(1.上海理工大學 能源與動力工程學院，上海 200093；2.上海市動力工程多相流動與傳熱重點實驗室，上海 200093)

《新時代的中國能源發展》白皮書指出：風電將在實現碳達峰和碳中和過程中發揮重要作用[1]。截至2020年底，我國風電裝機容量達2.8億kW，居世界首位[2]。在氣候雄心峰會上我國宣布到2030年，中國風電、太陽能發電總裝機容量將達到12億kW以上[3]。相較陸上風電，海上風電具有風速穩定、噪聲小及風電場可建設空間廣等優點，因此風電逐漸將發展重心轉向蘊藏豐富風能資源的近海地區[4]。

海上風力機運行過程中，面臨湍流風、波浪等復雜環境載荷影響，同時我國風能資源相對豐富的沿海地區多位于環太平洋地震帶，高強度地震載荷易對風力機的結構安全性和運行穩定性造成潛在危害[5]。因此，對于海上風力機結構設計研究，除環境載荷作用外，還須考慮地震載荷對其動力特性的影響。

Zuo等[6]通過分析5 MW風力機在湍流風與波浪載荷作用下的動力特性，發現氣動阻尼對風力機塔頂位移響應有較大影響。Wang等[7]研究發現，10 MW風力機在極端非線性波浪載荷及湍流風作用下，其動力特性由波浪載荷主導。Sajeer等[8]分析發現，當湍流風氣動力與波浪載荷水動力作用頻率及方向一致時，對風力機動力特性影響最大。但以上研究僅限于環境載荷對風力機動力特性的影響，未考慮地震作用效果。

針對地震誘導近海風力機動力響應特性研究，國內外學者開展了一系列研究工作。閆陽天等[9]發現10 MW風力機在不同強度地震作用下的動力學響應差異明顯。楊陽等[10]通過分析5 MW風力機在湍流風及地震載荷聯合作用下的動力學響應，表明塔基的振動響應主要由地震載荷決定。Meng等[11]基于5 MW風力機縮比實驗得出湍流風和地震的相互作用會使風力機結構整體響應降低。

以上研究均采用傳統數理方法對風力機動力響應進行分析，但作用于風力機葉片的非定常氣動載荷與結構相互耦合引起的幾何非線性效應和地震、波浪載荷作用下風力機基部動態響應之間的強耦合作用使整機運動響應呈顯著非線性特征，若僅采用常規方法對其在極端環境載荷作用下的動力響應時間序列數據進行分析，無法反映其內在動力特征。為此，部分學者提出將混沌理論應用于風力機設計中[12-13]，進而從多角度對多維度強非線性結構動力特性展開研究，以期實現對結構問題進行更透徹的剖析，為結構優化及控制系統建模提供更有效的理論基礎。綜上，筆者以DTU 10 MW單樁式近海風力機為研究對象，基于p-y曲線法及Winkler土-構相互作用模型建立土-構耦合模型，研究其在地震及環境載荷作用下的動力學響應，同時采用相空間重構定性分析和最大Lyapunov指數定量分析方法，揭示海上風力機動力響應混沌特性，為海上風力機結構抗震設計提供參考。

1 研究對象

1.1 大型單樁式風力機

采用Vestas公司與丹麥科技大學合作開發的DTU 10 MW單樁式近海風力機為研究對象。其主要結構和設計參數如圖1和表1[14]所示。

圖1 DTU 10 MW單樁式近海風力機模型

表1 風力機主要參數

風力機塔架壁厚隨高度的增加呈線性減小趨勢，塔基外徑為7.665 m，壁厚為34 mm；機艙與塔架連接處塔架直徑為5.5 m，壁厚為26 mm。模型材料選用密度為7 850 kg/m3、彈性模量為210 GPa、泊松比為0.3、材料屈服極限為380 MPa的A709結構鋼，由于構件油漆、焊接及法蘭等因素的影響，材料密度修正為8 500 kg/m3[15]。

1.2 土-構耦合

風力機入土樁基與地基之間彈性模量明顯不同，當地震載荷通過地基輸入樁基時，兩者之間所產生的變形量也存在較大差異，進而產生力的相互作用，即土-構耦合作用[16]。筆者根據美國石油協會(American Petroleum Institute, API)推薦的p-y曲線法[17]和Winkler土-構相互作用模型[18]建立土-構耦合模型，如圖2所示。其中,γ為有效重度，kN/m3；φ為砂礫內摩擦角，(°)。

當樁基發生橫向位移時，土壤作用力與側向位移呈非線性關系，且根據針對砂性土不同深度水平承載力推導出的經驗公式，砂性土的p-y曲線可表達為：

P=APutanh[kHy/(APu)]

(1)

Pus=(C1H+C2D)γH

(2)

Pud=C3DγH

(3)

Pu=min{Pus,Pud}

(4)

式中：A為經驗調整系數，A=(3-0.8H/D)，其中H為樁深，D為樁基外徑；P為土壤反力；y為側向位移；Pu為樁深至H時的極限承載力；Pus為樁深至H時的周向極限承載力；Pud為樁深至H時的徑向極限承載力；k為地基反力系數，其值取244.4 kN/m3；C1、C2、C3為系數，由砂礫內摩擦角φ確定，如圖3所示。

圖3 C1、C2與C3系數曲線

計算所得不同樁深處p-y曲線結果如圖4所示。

圖4 不同樁深處p-y曲線

2 地震及環境載荷

筆者構建了考慮湍流風、波浪及地震載荷作用下的極端環境多物理場非全耦合模型，氣動載荷部分通過OpenFAST/AeroDyn求解,波浪載荷在AQWA中生成,地震載荷基于實測近地位移數據產生。

2.1 地震載荷

為研究極端環境下不同地震強度與結構響應之間的關系，選用由國家強震動臺網中心提供的中國四川2017年九寨溝7.0級、2008年汶川6.0級及2014年越西5.2級地震數據。因我國海上風電場多位于東南沿海，故根據抗震設計規范及東南沿海地質條件[19]，確定以抗震設防烈度為7度，場地類型為Ⅳ，罕遇地震的反應譜對所選地震進行目標譜匹配，匹配后的位移時域曲線如圖5所示，7.0級、6.0級、5.2級地震對應的峰值地面位移(PGD)分別為36.20 cm、11.37 cm和7.99 cm。

2.2 湍流風載荷

參照IEC 61400-3設計準則[20]，基于IEC Kaimal湍流風譜模型，風在各方向上的功率譜密度為：

(5)

(6)

其中，湍流尺度參數ΛU為：

(7)

標準差σK等于湍流強度與平均風速之比，不同分量之間的關系定義如下：

(a) 7.0級地震

(b) 6.0級地震

(c) 5.2級地震

σv=0.8σu

(8)

σw=0.5σu

(9)

由于空間相干模型的存在，u分量的標準差會出現小幅波動。本文風場模型為IECKAI模型，其速度譜(以及標準差)穿過網格時保持恒定值。

通過TurbSim[21]生成平均風速為25 m/s的風場數據，輪轂高度處風速時域變化曲線如圖6所示。

圖6 湍流風場

2.3 波浪載荷

利用Morison方程計算作用在風力機下部的波浪與海流載荷，Morison方程可以有效地計算出振蕩流中作用在圓柱體上的力。

(10)

3 混沌及其判別方法

混沌是指確定性動力系統所表現出類似隨機的、不可預測的運動。目前，通常采用定性或定量方法對混沌現象展開分析，常見的定性研究方法有相空間重構法、龐加萊截面法及分頻采樣法等，常見的定量研究方法有最大Lyapunov指數法、自功率譜密度法和分形維數法等[22]。選取相對成熟簡便的相空間重構法和最大Lyapunov指數法對風力機運動響應分別進行定性和定量混沌分析。

3.1 相空間重構法

相空間重構法是研究非線性系統結構特性的成熟方法。筆者基于相空間重構法的基本思想，根據Takens理論和DTU 10 MW風力機的動力學響應時間序列信號，選取合適的延遲時間τ和嵌入維數m，構建出m維相空間，通過相空間中的軌線分布即可反映系統運動規律[23]：(1)若相軌跡趨于一點，則系統處于穩定狀態；(2)若相軌跡為一封閉曲線，則系統作周期運動；(3)若相軌跡雜亂無章分布在相空間，則系統作隨機運動；(4)若相軌跡具有一定規則，表明系統運動具有混沌特性。

在Takens理論中，對于理想的無限長和無噪聲的一維時間序列，嵌入維數m和延遲時間τ可以取任意值，但在實際應用中，因時間序列通常是含有噪聲的有限長序列，故其并不能任意取值，否則會影響相空間重構質量，進而影響后續數據處理結果[24]。筆者采用可反映數據點一般相關性的互信息函數法求取延遲時間τ，采用所需數據量相對較小的假近鄰法(FNN)計算嵌入維數m。

3.2 最大Lyapunov指數法

針對非線性信號，可通過引入Lyapunov指數λ來判斷其是否具有混沌特性。當λ>0時，系統具有混沌特性；當λ<0時，系統不具有混沌特性。λ表達式為：

(11)

4 結果與分析

4.1 時域特征分析

4.1.1 塔架位移

為研究不同地震強度及其與環境載荷聯合作用時風力機結構動力學響應特征，計算各條件下塔架位移最大值隨塔架高度的變化趨勢及塔頂和塔底的位移最大值，如圖7和表2所示。

(a) 僅地震載荷

(b) 環境載荷+地震載荷

表2 塔頂和塔底位移最大值

由圖7(a)可知，僅有地震載荷作用時，塔架各高度處位移最大值隨地震強度的增加而增加，這是由于風力機塔頂質量集中且塔架為細長柔性結構，隨地震強度的增加，在慣性及地震載荷共同作用下，塔架擺動程度加劇。

如圖7(b)所示，與其他條件下相比，僅有環境載荷作用時塔頂位移最大值最大，為1.519 7 m。當地震及環境載荷聯合作用時，隨著地震強度的增加，塔頂位移逐漸減小，說明塔頂位移主要為氣動彈性效應所導致的柔性變形量，當有地震載荷作用時，可有效緩解環境載荷所造成的塔頂振動。與其他條件下相比,僅有環境載荷作用時塔底位移最大值最小，為0.143 4 m。當地震及環境載荷聯合作用時，隨著地震強度的增加,塔底位移逐漸增加，說明地震載荷主要影響塔架中部及底部的動力學響應，其誘導塔頂振動幅值較小。

4.1.2 等效應力

在極端環境下，風力機在受到地震及環境載荷等外力作用時會產生較大形變，形變位置會產生應力集聚現象，進而使風力機結構遭受嚴重損壞。因此，計算塔架在地震及環境載荷聯合作用下等效應力達到峰值時其表面等效應力分布情況，結果見圖8。

圖8 塔架等效應力分布

如圖8所示，無地震載荷作用時，塔底等效應力最大，其隨塔架高度的增加而逐漸減小。等效應力在流向相較其他方向有更明顯的集聚現象,這是由于在環境載荷作用下，塔架變形主要為前后向彎曲變形，因此等效應力在流向較大。有地震載荷作用時，地面震動會緩解一部分因波浪載荷所造成的應力集聚，因此塔架底部等效應力明顯減小。隨著地震強度增加，塔架變形量逐漸增大，因此塔架中部的等效應力峰值也隨之增大。

4.1.3 塔頂加速度

風力機塔頂加速度為反映塔架振動的重要參數，隨著風力機朝大型化發展，塔架結構柔性增加，使其在極端環境中易產生較大形變，為檢測復雜環境下風力機塔架振動響應，計算不同地震及環境載荷下風力機塔頂流向及側向的加速度時域曲線，結果如圖9所示。

(a) 流向加速度

(b) 側向加速度

由圖9可知，地震載荷對塔頂側向加速度的影響比對塔頂流向加速度大。塔頂流向加速度主要決定因素為環境載荷，當有地震發生時，塔頂流向加速度沿僅有環境載荷時的加速度曲線上下波動；側向加速度則以0 m/s2為基線上下波動。因此塔頂流向加速度變化對地震的敏感性相對較小，側向加速度對地震的敏感性較大。

同時，塔頂側向加速度峰值大于流向加速度峰值，說明側向相比流向的振動更劇烈。其主要原因是風輪重心懸垂于塔架中心線，導致塔頂質量在俯仰方向的慣性矩更大，且流向所具有的氣動阻尼對塔架流向加速度及位移具有較強的抑制作用，從而在側向具有更大的加速度。

4.2 混沌特性分析

由4.1.3節可知，風力機塔頂振動伴有明顯的隨機性和結構不穩定性，因此采用混沌方法對風力機運動特性進行非線性分析。在不同地震及環境載荷聯合作用下，風力機塔頂流向及側向加速度的響應最佳延遲時間τ如表3所示。

混沌的產生是由于系統的整體穩定性和局部的不穩定性共同作用的結果，系統局部不穩定性使其具有對初值敏感的特性，而整體穩定性使其在相空間表現出一定的分形結構，這種結構就稱為混沌吸引子[26]。根據相空間重構法繪制風力機在不同地震及環境載荷下流向及側向加速度時域數據的三維混沌吸引子，如圖10和圖11所示,其中t為當前迭代時間。

表3 最佳延遲時間

圖10 流向加速度混沌相圖

由圖10和圖11可知，塔頂加速度在不同地震及環境載荷下混沌吸引子形態各異，其分別代表了不同狀態下塔頂加速度的動力特性。

如圖10(a)所示，在僅有環境載荷作用時，塔頂流向加速度的吸引子軌跡反復纏繞形成一種細長“環狀”結構。當有地震載荷作用時，隨地震強度的增加,吸引子環狀區域逐漸稠密，最終形成“團絮狀”結構。此現象與圖9(a)所示的當有地震發生時，塔頂流向加速度沿僅有環境載荷時的加速度曲線上下波動現象相對應，進一步說明塔頂流向加速度主要決定因素為環境載荷。且由吸引子結構變化規律可知，隨著地震強度的增加，塔頂流向加速度響應的混沌特性逐漸增強，塔頂穩定性減弱。

如圖11(a)所示，當僅有環境載荷作用時，塔頂側向加速度吸引子軌跡反復纏繞成“棒狀”結構。當有地震載荷作用時，其吸引子軌跡糾纏發散成盤狀，且隨著地震強度的增加，其波動軌跡逐漸擴散加密，說明隨地震強度增加，其塔頂側向加速度混沌特性逐漸增強，塔頂穩定性減弱。

圖11 側向加速度混沌相圖

由以上分析可知，塔頂流向及側向加速度響應既非完全周期性響應也非完全隨機響應，而是具有鮮明的混沌特性。為進一步研究這一特性，計算不同地震及環境載荷聯合作用下塔頂流向及側向加速度響應的最大Lyapunov指數，結果見表4。

表4 最大Lyapunov指數

由表4可知，各種情況下塔頂流向及側向加速度響應的最大Lyapunov指數均大于0，說明風力機塔頂運動響應具有混沌特性，且隨地震強度的增加，最大Lyapunov指數逐漸變大，混沌特性增強。

5 結 論

(1) 地震與環境載荷聯合作用下,塔底位移最大值較僅有環境載荷時明顯變大，說明地震載荷加劇了風力機塔底振動。地震載荷主要影響塔架中部及底部動力學響應，其誘導塔頂振動幅值較小，塔頂振動主要為氣動彈性效應導致的柔性變形量。

(2) 塔頂流向加速度主要決定因素為環境載荷，該加速度變化對地震的敏感性相對較小，塔頂側向加速度對地震的敏感性較大。

(3) 塔頂加速度響應既非完全周期性響應也非完全隨機響應，而是具有鮮明的混沌特性。在不同載荷條件下其混沌吸引子形態各異，分別代表了不同狀態下塔頂加速度的動力特性。隨地震強度的增加，塔頂流向及側向加速度響應的吸引子逐漸擴散，最大Lyapunov指數逐漸變大，混沌特性逐漸增強，塔頂穩定性逐漸減弱。