基于Craig-Bampton方法的海上風電機組基礎結構超單元計算精度研究

鐘 耀，王文華，周 潤，李 昕

(1.中國電建集團中南勘測設計研究院有限公司，湖南 長沙 410014；2.大連理工大學海岸和近海工程國家重點實驗室，遼寧 大連 116024)

1 國內外研究進展綜述

與陸上風能資源開發不同，海上風電場建設面臨著復雜的環境條件，比如我國東部沿海地區的臺風、環渤海海域的地震及海冰等均給海上風電機組(海上風機)及基礎設計帶來了技術挑戰。上述復雜環境條件下如何建立海上數值仿真計算模型，并準確得到其結構響應已成為現階段海上風電結構領域的研究重點。

目前，海上風電結構分析模型普遍采用整體耦合模型和半整體模型，整體耦合模型即直接建立包含轉子—機艙—控制—塔筒—基礎結構的整體結構數值仿真模型，并開展耦合時域反應分析，已廣泛應用于浮式和固定式海上風機耦合機理及運行控制策略研究。Ma[1]、Li[2]、萬德成[3]和劉穎蓮[4]等均基于整體耦合模型揭示了Spar、TLP等不同類型浮式風機整體結構耦合反應機理。王洪哲等[5]采用整體耦合分析方法研究了極端風浪作用下海上風電結構響應規律，并提出了相應的運行控制策略。方通通[6]、Rezaei[7]和Koukoura[8]等采用整體方法開展了隨機風浪作用下海上風電疲勞分析，系統研究了環境荷載耦合效應、風浪聯合分布及控制策略等對于海上風電結構疲勞累積的影響。

由上述研究可得，基于整體耦合模型能夠得到更為合理的海上風電結構反應，尤其對于結構疲勞累積；不過，方通通[6]和Kvittem[9]等指出采用整體耦合疲勞計算模型將顯著增加數值仿真計算成本。與整體耦合模型不同，半整體模型通過基礎結構等效剛度、質量矩陣代替相應有限元模型，建立轉子—機艙—塔筒—等效特征矩陣的氣彈數值仿真計算模型；將氣彈分析所得風機荷載作用于塔筒與基礎結構交接點，并施加其他環境荷載(如波浪、海冰、地震等)于基礎結構有限元模型，計算得到海上風機結構響應。張戈[10]、夏露[11]、陳法波[12]和莫仁杰[13]等均采用該模型開展了風浪和地震等不同類型荷載作用下固定式海上風電結構動力反應計算，Minguez等[14]采用該方法系統對比了單樁、三樁和導管架基礎結構的極限承載力。王鵬等[15]以基礎結構用鋼量為優化目標，采用半整體方法對海上風機三樁基礎結構進行了優化設計。

本文對半整體模型關鍵技術基礎結構超單元矩陣計算展開研究，探討將Craig-Bampton方法應用于海上風機基礎結構剛度、質量超單元矩陣計算時，其計算參數選取對于所得超單元矩陣計算精度的影響。進一步，分別建立固定式海上風電整體耦合模型和半整體模型，對基于Craig-Bampton方法所得超單元矩陣計算精度進行驗證。

2 半整體模型運動方程

空氣動力荷載作用下海上風電轉子—機艙—塔筒—超單元的運動方程為

(1)

進一步，可得風機荷載作用下海上風機基礎結構運動方程為

(2)

基礎結構超單元計算普遍采用Craig-Bampton(C-B)方法[16]。Craig-Bampton方法將海上風機基礎結構節點劃分為邊界節點和內部節點，劃分后結構剛度和質量矩陣如下

(3)

基于該方法可得海上風機基礎結構超單元計算公式，即

(4)

(5)

其中，

(6)

(7)

(8)

(9)

3 樣本風機有限元模型

選用我國近海風電場軟弱地基所應用的五樁基礎結構，如圖 1所示，并結合NREL 5MW基準風機[17]得到本次研究的樣本海上風電機組，其中，NREL 5MW基準風機采用變速變槳式控制方式，輪子、輪轂直徑分別為126、3 m，切入、額定、切出風速分別為3、11.4、25 m/s、切入、額定轉速分別為6.9、12.1 rad/min。基于ANSYS[18]，選用PIPE單元建立基礎結構有限元模型。同時，選用COMBINE單元，基于p-y、t-z和Q-z承載力曲線建立非線性樁土相互作用模型。由此，所建立的基礎結構有限元模型如圖2所示。

圖1 五樁基礎結構基本參數(單位：mm)

圖2 五樁基礎結構有限元模型

4 基礎結構超單元計算分析

4.1 凝聚節點影響研究

基于圖 2所示有限元模型主要桿件節點編號，分別選取基礎結構主筒體模型節點1 000、2 000、3 000和4 000作為凝聚節點，開展圖 3所示風機荷載Fx作用下，所選取凝聚節點對于C-B方法所得基礎結構超單元計算精度影響研究。風機荷載Fx作用于基礎結構頂部位置。

圖3 風機荷載Fx

風機荷載作用下基礎結構頂部位移統計值相對誤差隨所選取凝聚節點變化規律如圖4所示，相對誤差計算公式為

(10)

式中，S1為荷載作用下基于圖2a所示有限元模型所得結構響應統計值；S2為考慮不同凝聚節點影響，基于基礎結構超單元矩陣所得結構響應統計值。

圖4 基礎頂部位移統計值相對誤差

由圖4可知，凝聚節點對于C-B方法所得基礎結構超單元矩陣計算精度具有顯著影響，相對于其他節點，當指定基礎頂部節點(節點4 000)為凝聚節點時，所得超單元矩陣具有較高計算精度，此時所得基礎結構頂部位移統計值最大相對誤差僅為2.3%。

4.2 內部節點主模態數量影響研究

選取基礎結構頂部節點4 000作為凝聚節點，同時分別選取基礎結構內部節點前4、14、34、54、94階模態作為主模態，基于C-B方法開展基礎結構超單元計算。基于超單元計算結果，施加風機荷載Fx，所得基礎結構頂部位移時程統計值與有限元模型計算結果時頻域對比如圖5所示。

圖5 風機荷載作用下基礎頂部位移時頻域反應對比

(11)

(12)

由圖5a可知，當所保留內部節點主模態數量為4階時，此時C-B方法所得超單元矩陣具有較高的計算精度，如所得基礎結構頂部位移統計值最大相對誤差僅為2.5%。隨著擬保留內部節點主模態數量的增加，超單元矩陣計算精度逐步提高。當所保留主模態數量達到34階時，此時基于基礎結構超單元矩陣所得風機荷載Fx作用下結構響應與有限元模型計算結果基本一致。

由圖5b可知，對比包含不同階數主模態向量的基礎結構超單元矩陣，均能保證風機荷載Fx作用下基礎結構位移反應控制頻率與有限元模型計算結果保持一致，如所得結構反應主要控制頻率均為0.319 Hz，故所選取超單元矩陣均具有較高計算精度。

因此，對于所選取風機荷載工況，當選用C-B方法開展基礎結構超單元計算時，如若選取基礎結構頂部節點為凝聚節點，同時內部節點主模態數量≥34階，此時得到的超單元具有足夠的計算精度。

5 基于整體耦合模型的海上風機超單元矩陣驗證

依據樣本風電結構參數，基于FAST[19]建立樣本海上風電機組的整體耦合模型。同時，依據第4部分研究所得C-B方法計算參數選取原則，選取塔筒與基礎結構交界點(節點4 000)為擬凝聚節點，將基礎結構總剛度及質量矩陣凝聚至該節點，所得與該節點自由度相耦合超單元矩陣如式(11)和式(12)所示。

基于式(11)和式(12)所示超單元矩陣，依據半整體模型運動方程，建立樣本的轉子—機艙—塔筒—超單元矩陣半整體計算模型。

分別選取8、11.4 m/s和14 m/s平均風速，基于IEC Kaimal譜擬合得到隨機風速1、2和3如圖6所示，開展上述隨機風速時程作用下海上風電整體耦合與半整體模型對比研究。隨機風速1作用下，基于整體耦合模型和半整體模型所得樣本海上風電機組塔筒頂部推力荷載、加速度時程響應統計值對比如圖7所示。

圖6 隨機風速時程

圖7 隨機風速作用下海上風電機組運動響應統計值計算誤差

由圖7可知，基于整體耦合模型和半整體模型所得海上風電結構運動響應基本相同。由圖7a可知，對于塔筒頂部推力荷載，隨機風速2作用下，基于上述模型所得結構響應平均值計算誤差最為顯著。隨機風速2的平均風速為額定風速，即在該風速作用下風機將啟動變槳控制策略。由此，進一步對比可得，由于變槳控制策略影響，隨機風速3作用下推力時程均值相對誤差相對于隨機風速2作用工況明顯減小。由圖7b可知，與塔筒頂部推力荷載變化規律稍有不同，隨機風速1作用下，半整體模型所得塔筒頂部加速度時程極大值和極小值相對誤差最為明顯，此時計算誤差將分別達到11.6%和6.0%。需指出，在隨機風速2和3作用下整體與半整體模型所得加速度時程統計值基本一致，比如上述工況下半整體模型所得加速度統計值最大相對誤差僅為5.8%。

綜上所述，將塔筒與基礎結構連接節點指定為基礎結構總剛度和質量矩陣凝聚節點，基于C-B方法所得超單元矩陣建立海上風電半整體計算模型所得結構時頻域響應與整體耦合模型計算結果基本一致。此外，海上風機運行狀態及控制策略對于半整體模型計算精度具有顯著影響，在達到額定風速后，啟動相應變槳控制策略將顯著提高半整體模型塔筒頂部推力荷載計算精度。

6 結 論

本文通過與整體耦合模型的對比，驗證了基于C-B方法計算得到基礎結構超單元矩陣建立海上風電半整體模型的計算精度。通過本次研究可得以下結論：

(1)基于有限元模型研究可得，凝聚節點選取對于Craig-Bampton方法所得超單元矩陣計算精度具有顯著影響；對于所選取工況，當凝聚節點選取為基礎結構頂部節點時，所得超單元矩陣具有較高計算精度。

(2)內部結點主模態數量對于超單元矩陣計算精度的影響亦不可忽視；對于所列工況，當內部主模態數量達到34階時，此時所得超單元矩陣計算結果與有限元模型基本一致。

(3)選取塔筒與基礎結構交界點為凝聚節點，基于C-B方法所得基礎結構超單元矩陣建立半整體模型所得海上風電結構運動響應與整體耦合模型計算結果基本一致。

(4)半整體模型計算精度與海上風機運行狀態及控制策略密切相關，在額定風速以上啟動變槳控制策略后，將顯著提高半整體模型塔筒位置推力荷載計算精度。