考慮多場耦合的碰撞載荷下半潛式風機動力響應

王禹凱，宋 娜，王加夏，劉 昆

(江蘇科技大學 船舶與海洋工程學院，江蘇 鎮(zhèn)江 212003)

0 引 言

隨著經(jīng)濟的快速發(fā)展，清潔的可再生能源逐漸成為各個國家開發(fā)的重點[1]，其中，風能以其無污染且儲量大的特點成為近年來發(fā)展較快的新能源之一[2-3]。開發(fā)風能的方式包含陸上風機和海上風機。海上風機可分為固定式風機和海上浮式風機兩種類型。海上浮式風機通?？拷降溃趷毫犹鞖庀拢芭c其發(fā)生碰撞事故的概率會大幅提高。此外，風電場需要運維船定期維護，且風電場所在的海域常有較大的風浪，這也會增加碰撞風險。因此，研究半潛式風機在撞擊載荷下的運動規(guī)律及結(jié)構(gòu)響應情況對其耐撞性的優(yōu)化及使用壽命的延長有重要的理論價值，對保證海上風電場的安全可靠性具有至關(guān)重要的意義。

半潛式風機的動態(tài)響應包含風、浪、流和系泊系統(tǒng)等復雜的多場載荷耦合問題。WAYMAN[4]采用線性頻域水動力法開發(fā)一套在線性頻域下計算水動力、氣動力和結(jié)構(gòu)耦合響應的分析程序。WANG等[5]提出一種基于AQWA水動力軟件的浮式風機耦合響應的簡單頻域法，研究風、浪、流和系泊系統(tǒng)阻尼對浮式基礎(chǔ)響應的影響，但頻域法的局限性在于其不能捕捉非線性動態(tài)特性和瞬態(tài)效應。JONKMAN等[6]開發(fā)一套能全耦合考慮空氣動力學、水動力學、控制系統(tǒng)、結(jié)構(gòu)分析的程序以分析Barge型浮式風機的結(jié)構(gòu)動力響應。MA等[7]基于此程序分析Hywind Spar浮式風機系統(tǒng)在不同海況下的運動和系統(tǒng)載荷響應特性，討論風、浪對風機系統(tǒng)的影響，但忽略系泊系統(tǒng)的影響。為了更精確地模擬正常作業(yè)條件下浮式風機系統(tǒng)的性能，對浮式風機的風-浪-系泊進行全耦合數(shù)值模擬是必要的。

本文以DeepCwind浮式風機模型為研究對象，提出一種風-浪-系泊耦合場景下船舶-浮式風機碰撞動態(tài)響應數(shù)值分析方法。研究風-浪-系泊耦合場景下半潛式風機在船舶撞擊載荷作用下的動力響應，有效地模擬半潛式風機在沖擊載荷作用下的動力過程。

1 船舶-浮式風機碰撞耦合算法

1.1 三維勢流理論

對于大尺寸構(gòu)件(D/λ>0.2，其中：D為橫截面直徑；λ為波長)，入射波經(jīng)過該結(jié)構(gòu)時會受到嚴重的輻射和繞射，需要選取基于輻射和繞射理論的勢流理論預報浮體結(jié)構(gòu)的水動力系數(shù)和波浪力。當入射波與浮體結(jié)構(gòu)互相作用時，波浪的速度勢可分解為3個部分：入射波速度勢ΦI，繞射波速度勢ΦD和輻射波速度勢ΦR。其中：繞射波速度勢指的是由于浮體結(jié)構(gòu)(海上風機平臺)的存在，入射波流過浮體結(jié)構(gòu)時產(chǎn)生繞射波；輻射波速度勢是由海上風機平臺自身的運動而產(chǎn)生的波浪，這是有關(guān)于置于初始靜止流場產(chǎn)生微幅振蕩的結(jié)構(gòu)與這個流場之間的互相作用問題。

單色波入射情況下的一階速度勢為

Φ=Κ(ΦI+ΦD+ΦR)=Re{[ΦI,1(x,y,z)+ΦD,1(x,y,z)+ΦR,1(x,y,z)]eiω t}

(1)

式中：ΦI,1為一階線性入射波速度勢；ΦD,1為一階線性繞射波速度勢；ΦR,1為一階線性輻射波速度勢；ω為波浪頻率；t為時間；Κ為波浪斜率且定義為

(2)

式中：A為波幅;L為波長；k為波數(shù)。

式(1)的一階線性入射波速度勢可由線性波浪理論計算得出：

(3)

式中：g為重力加速度，取9.81 m/s2；z為波浪高度；d為水深；k=kexcosθ+keysinθ，其中θ為波浪入射角，ex和ey分別為沿x軸和y軸的單位向量；x=xex+xey。

(1)繞射速度勢需滿足如下邊界條件：

(4)

式中：n為物體表面的單位法向向量，方向朝外；r為物體表面的位置。

(2)輻射速度勢需滿足如下邊界條件：

(5)

式中：ξ1為一階平移運動；α1為一階轉(zhuǎn)向運動。

通常利用格林函數(shù)計算包含邊值條件的拉普拉斯方程，然后建立速度勢與格林函數(shù)之間的積分方程。需要利用邊界單元法計算該積分方程，將物體離散為一定數(shù)目的單元，之后假設(shè)每個單元內(nèi)的速度勢都可利用單元節(jié)點勢的多項式表達，基于伽遼金法或配點法建立節(jié)點勢的線性方程組，得出各節(jié)點處的速度勢。在得到繞射勢和輻射勢后，一階線性水動壓強和自由表面高程[8]為

(6)

(7)

式(6)和式(7)中：p1為一階線性水動壓強；ρ為海水密度，取1 025 kg/m3；Φ1為一階線性速度勢；η1為一階自由表面高程，是時間和空間的函數(shù)。

沿物體的濕表面進行積分就能得到總波浪力(F1、F2、F3)和力矩(M1、M2、M3)為

(8)

(9)

式(8)和式(9)中：S為物體的濕表面。另外，一階波浪力F1可定義為3個部分的和：

F1=FHS,1+FR,1+FEX,1

(10)

式中：FHS,1為一階剛性靜水回復力；FR,1為一階輻射力；FEX,1為一階波浪力(由繞射速度勢和入射波速度勢產(chǎn)生)。

1.2 多場耦合分析方法

圖1為DeepCwind浮式風機的頻域分析流程。在有限元軟件ANSYS Workbench和Hypermesh中分別對DeepCwind浮式風機進行建模與網(wǎng)格劃分。利用AQWA的Hydrodynamic Diffraction模塊對浮式風機進行靜水分析，得到附加質(zhì)量慣性矩和靜水回復剛度矩陣。由于水動力軟件AQWA采用勢流理論，無法計及黏性阻尼的影響，與實際情況會有較大差異，因此通過黏性阻尼經(jīng)驗公式計算臨界阻尼和黏性阻尼并將其添加至軟件中以達到阻尼修正的目的。關(guān)于半潛式風機基礎(chǔ)黏性阻尼的計算，目前還缺乏成熟的方法。垂蕩、橫搖和縱搖方向上等效黏性阻尼可取臨界阻尼的10%。臨界阻尼計算式[9]為

圖1 DeepCwind浮式風機頻域分析流程

(11)

式中：M為質(zhì)量；Kz為對應自由度的靜水回復剛度。式(11)變?yōu)?/p>

(12)

式中：I為質(zhì)量慣性矩；ΔI為附加質(zhì)量慣性矩。

之后利用頻域分析研究浮式結(jié)構(gòu)在水中的運動響應。應用AQWA對DeepCwind浮式風機進行水動力性能數(shù)值計算，采用經(jīng)典勢流理論的計算方法，利用頻域計算的快速性和線性化特點求解單位規(guī)則波下的半潛式風機平臺的水動力參數(shù)，主要包括：附加質(zhì)量系數(shù)、輻射阻尼系數(shù)以及不同頻率、不同浪向下的一階波浪力等。

在頻域分析中沒有考慮系泊系統(tǒng)的回復力和風載荷的影響，因此計算結(jié)果不能完全準確地反映半潛式風機的運動響應，而時域分析可用于船舶與海洋結(jié)構(gòu)物在特定的海洋環(huán)境條件下的一定運動范圍內(nèi)的運動情況仿真，能夠較真實地模擬船舶與海洋結(jié)構(gòu)物的運動響應，并且可以考慮在各時間步內(nèi)半潛式風機系統(tǒng)各部分之間的耦合作用，充分考慮前一時刻的運動響應對后一時刻的影響[10]。對于一些特定的問題，如結(jié)構(gòu)物之間的撞擊、海底擱淺都能清楚地反映其具體發(fā)生的時刻和全過程的時間歷程。圖2為DeepCwind浮式風機時域分析流程。將AQWA計算出來的水動力參數(shù)導入美國可再生能源實驗室(National Renewable Energy Laboratory，NREL)開發(fā)的FAST程序[11]中，使用FAST程序考慮風、浪以及系泊系統(tǒng)與半潛式風機之間的耦合作用，計算半潛式風機的運動響應、氣動推力、水動力和錨鏈力。

圖2 DeepCwind浮式風機時域分析流程

在LS-DYNA中無法完成與結(jié)構(gòu)動力響應相關(guān)的加載模式，即載荷無法隨著結(jié)構(gòu)響應變化而變化，無法達到載荷與結(jié)構(gòu)耦合的效果，因此通過子程序Loadud LS-DYNA加載。

圖3為載荷與結(jié)構(gòu)響應相互耦合的實現(xiàn)流程，先將FAST程序中計算出的氣動推力載荷和水動力載荷導入LS-DYNA進行結(jié)構(gòu)動力學計算。由于在FAST程序中模型的動力學參考點為水線面中心點，水動力載荷也作用在水線面中心點，因此在LS-DYNA中，將水動力載荷通過子程序Subroutine Loadud作用在水線面中心點，即原點處。氣動推力則進行適當簡化，時域載荷作為一個集中力直接施加在風機結(jié)構(gòu)上，具體作用點如圖4所示。將錨鏈力進行適當簡化后作用于半潛式風機平臺系泊點處。當T

圖3 耦合算法說明

根據(jù)半潛式風機的實際情況，將系泊系統(tǒng)作為邊界條件。在LS-DYNA中，如果準確地建立錨鏈模型并考慮其與水的相互作用，將會額外帶來較大的計算負擔。因此，為提高計算效率并保證一定的計算精度，對半潛式風機的系泊系統(tǒng)進行適當簡化。當半潛式風機小角度運動且位移較小時，可將系泊纜繩簡化成線性彈性系統(tǒng)。在FAST中，利用中心差分方法對偏導數(shù)進行計算，得到系泊纜繩的彈性剛度矩陣為

(13)

在LS-DYNA中，在半潛式風機的3個系泊點處均采用3個轉(zhuǎn)動和3個平移彈簧單元對系泊系統(tǒng)進行建模，如圖4和圖5所示。

圖4 載荷添加位置示例

圖5 彈簧單元有限元模型

2 數(shù)值計算模型

2.1 半潛式風機模型描述

以NREL提出的DeepCwind浮式風機模型為研究對象，該風機發(fā)電功率為5 MW，額定風速為11.4 m/s，基本參數(shù)[12]如表1所示。

表1 DeepCwind浮式風機模型主要參數(shù)

采用有限元軟件ANSYS Workbench對Deep-Cwind浮式風機進行建模，風機基礎(chǔ)主要由3個浮筒組成，并且有橫支撐桿和斜支撐桿。此外，對風機模型進行簡化，機艙及葉片采用集中質(zhì)量代替。半潛式風機會產(chǎn)生6個自由度的平移和旋轉(zhuǎn)運動(見圖6)，平移運動包括橫蕩、縱蕩和垂蕩，旋轉(zhuǎn)運動包括縱搖、橫搖和艏搖?？v坐標原點是在水線處由3根圓柱組成的等邊三角形的重心，z軸垂直向上，xy平面與海平面重合。

圖6 半潛式風機運動示例

2.2 撞擊船模型描述

考慮到深海風電場海域內(nèi)實際的過往船舶狀況，選用風電運維船作為撞擊船。撞擊船主要參數(shù)如表2所示。在計算中主要考慮DeepCwind浮式風機在碰撞中的結(jié)構(gòu)動力響應，為縮短計算時長，撞擊船采用剛性體且建模時對船體內(nèi)部結(jié)構(gòu)進行簡化，并且省略其上層建筑。船舶有限元模型及網(wǎng)格如圖7所示。

表2 撞擊船主尺度

圖7 碰撞船有限元模型及網(wǎng)格劃分

將碰撞船有限元模型與DeepCwind浮式風機裝配在一起，為模擬浮式風機在運動狀態(tài)下遭受船舶撞擊的效果，撞擊船與浮式風機之間在發(fā)生碰撞前有一定的距離，撞擊船與風機的撞擊點位于風機吃水線附近。

此外，選取主-從接觸面方法分別定義在半潛式風機浮筒和撞擊船上，選用LS-DYNA有限元軟件提供的自動面面接觸(*CONTACT_AUTOMATIC_SURFACE_TO_SUFACE)，動、靜摩擦因數(shù)均取0.2[13]。

在撞擊船作用下，被撞半潛式風機平臺的立柱和支撐構(gòu)件可能會發(fā)生塑性變形，對于被撞半潛式風機采用常用的塑性動態(tài)模型(*MAT_PLASTIC_KINEMATIC)，預先假定一個最大塑性失效應變，若所計算的單元達到最大塑性失效應變，則該單元即刻失效，所對應的構(gòu)件就認為已損壞而不再承受任何外力。為節(jié)約計算時間，撞擊船舶采用剛性體(*MAT_RIGID)，半潛式風機和撞擊船的材料分別為船用高強鋼(Q345)和低碳鋼(Q235)，其材料參數(shù)[14]如表3所示。

表3 鋼材性能參數(shù)

撞擊船在碰撞前的運動中，需要考慮船舶結(jié)構(gòu)與水之間的相互作用，一般通過附加質(zhì)量模型和流固耦合模型解決?？紤]到流固耦合模型會使計算時間過長，且主要研究對象為半潛式風機，而附加質(zhì)量模型方法已得到過許多學者[15]的使用與驗證，因此選用附加質(zhì)量模型處理船體結(jié)構(gòu)與水之間的相互作用。此外，附加質(zhì)量為撞擊船質(zhì)量的0.05倍。

2.3 碰撞場景

選取典型碰撞場景，在天氣惡劣的條件下，海上運維船失控后與半潛式風機發(fā)生碰撞。根據(jù)DNV GL-ST-0126[16]和DNV-ST-0119[17]規(guī)范，當海上風機與船舶發(fā)生碰撞時，船舶速度不低于0.5 m/s，因此碰撞船舶的撞擊速度均取1.0 m/s，撞擊船與半潛式風機的撞擊點位于風機吃水線附近。海洋環(huán)境則取半潛式風機額定風速11.4 m/s，波浪高度為1.5 m，風向和浪向均取y軸負方向，如圖8所示。時間步長取0.05 s，模擬總時長為150 s。

圖8 簡化的碰撞場景示例

3 多場耦合下碰撞仿真結(jié)果及分析

3.1 碰撞力

圖9為撞擊船首部撞擊DeepCwind浮式風機碰撞力時程曲線，該曲線總體呈明顯的非線性狀態(tài)，碰撞是在非常短的時間內(nèi)發(fā)生的動力沖擊問題，整個碰撞過程持續(xù)時間約0.4 s，在船首與DeepCwind浮式風機的浮筒接觸開始，從0.05～0.20 s碰撞力急劇增加，這是由于隨著風機浮筒參與的撞擊范圍不斷增加，撞擊船首部所受抵抗力不斷增加，并于0.25 s時達峰值3.723×106N。此后，在0.25～0.45 s，碰撞力出現(xiàn)卸載情況迅速減小，這是由于在撞擊船撞擊風機浮筒的過程中，船速不斷降低，動能下降，直至碰撞過程結(jié)束，碰撞力也減小為零。

圖9 碰撞力時程曲線

3.2 碰撞能量

系統(tǒng)總體能量轉(zhuǎn)化如圖10所示。隨著撞擊船與DeepCwind浮式風機的浮筒接觸，撞擊船的動能主要轉(zhuǎn)化為浮筒碰撞部位的變形內(nèi)能和浮式風機整體動能兩部分，還會產(chǎn)生少量的沙漏能、阻尼能和摩擦能。阻尼能和摩擦能的數(shù)量級極低，在圖10中略去。由圖10可知，碰撞始于25.95 s，由于浮式風機靠系泊系統(tǒng)進行定位，浮式基礎(chǔ)有較大的靈活性，因此撞擊船的動能大部分轉(zhuǎn)換為浮式風機的動能，動能最大值達1.37×108J，少量轉(zhuǎn)化為風機的內(nèi)能和沙漏能，這與船-船、船-固定式風機等碰撞問題不同。這說明在船-浮式風機的碰撞運動中，撞擊對浮式風機的損傷變形影響較小，但對浮式風機的運動響應影響較大。

圖10 系統(tǒng)能量轉(zhuǎn)化曲線

3.3 結(jié)構(gòu)應力應變

圖11為25.95 s時船舶與DeepCwind浮式風機發(fā)生正面碰撞的整體結(jié)構(gòu)應力云圖。由圖11可知，在船舶碰撞浮式風機的過程中，所產(chǎn)生的應力主要分布于浮筒的撞擊區(qū)域，存在明顯的局部性。

圖11 DeepCwind浮式風機與船首碰撞整體結(jié)構(gòu)應力云圖

圖12為DeepCwind浮式風機遭受撞擊的浮筒在不同時刻的結(jié)構(gòu)等效應力云圖。由圖12可知，在25.95～26.05 s，風機浮筒的結(jié)構(gòu)應力峰值逐漸增大，并于26.05 s達最大值3.626×108Pa，隨后隨著碰撞力的卸載而減小，且結(jié)構(gòu)應力分布有明顯的局部性。

圖12 撞擊過程中不同時刻浮筒結(jié)構(gòu)應力云圖

圖13為DeepCwind浮式風機遭受撞擊的浮筒在不同時刻的結(jié)構(gòu)塑性應變云圖。由圖13可知，結(jié)構(gòu)塑性應變在25.95～26.05 s呈逐漸增大的趨勢，結(jié)構(gòu)塑性應變區(qū)域集中于撞擊部位且損傷面積較小，呈現(xiàn)明顯的局部性，說明撞擊載荷對DeepCwind浮式風機的結(jié)構(gòu)損傷較小，撞擊船的動能只有一小部分轉(zhuǎn)化為內(nèi)能。

圖13 撞擊過程中不同時刻浮筒結(jié)構(gòu)塑性應變云圖

3.4 運動響應

圖14為 DeepCwind浮式風機在有無撞擊載荷下的橫蕩、橫搖、艏搖和垂蕩運動響應對比圖。風、浪和船舶撞擊方向皆為y軸負方向，縱搖和縱蕩因量級較低而忽略。由圖14可知，各曲線都呈現(xiàn)明顯的非線性狀態(tài)，垂蕩自由度的運動響應與其他不同，在有撞擊載荷下垂蕩自由度運動幅度低于無撞擊載荷下的運動幅度，但是有撞擊載荷下的橫蕩、橫搖和艏搖運動幅度明顯大于無撞擊載荷下的運動幅度，說明撞擊載荷對DeepCwind浮式風機的運動響應影響較大，這是由于撞擊船所產(chǎn)的動能大部分轉(zhuǎn)化為風機的動能，這些動能使DeepCwind浮式風機的橫蕩、橫搖和艏搖運動幅度加大。

圖14 有無撞擊載荷下不同工況運動響應對比圖

4 結(jié) 論

通過有效地模擬半潛式風機遭受撞擊的動態(tài)過程，得到如下結(jié)論：

(1)在船舶撞擊半潛式風機時，由于浮式風機采用系泊系統(tǒng)，浮式平臺具有一定的靈活性，撞擊船的動能大部分轉(zhuǎn)化為浮式風機的動能，使風機的橫搖、橫蕩和艏搖運動幅度皆大于無撞擊載荷下的運動幅度。

(2)在船舶碰撞半潛式風機的過程中，風機浮筒因吸收了較小的內(nèi)能而未有較大的損傷，所產(chǎn)生的應力應變主要分布于浮筒的撞擊區(qū)域，存在明顯的局部性，且風機浮筒的結(jié)構(gòu)應力呈現(xiàn)先增大后減小的趨勢。