韓長志,陳淑玲,楊 帥,許思乾
(江蘇科技大學,江蘇 鎮江 212003)
從“十二五”到“十四五”,國家風電發展政策逐漸向海上發電傾斜。2020 年,我國海上風電累計裝機約900 萬千瓦,其中新增裝機達到306 萬千瓦。隨著我國風電技術的日漸成熟,風電的深水化、大功率化已成為必然的發展趨勢。為保證風機的運行及降低建造成本,基礎的形式也由傳統的固定式向漂浮式發展。在水深大于80~100 m 海域,相對于固定式基礎,浮式基礎無論經濟效益還是安全性能,都處于領先地位。但是浮式風機鋪設海域的海況更加復雜,更易受到大波浪的影響,從而其穩定性、可靠性和生存能力相對于淺海固定式風機更具挑戰性。
近年來,依靠計算流體力學(Computational Fluid Dynamics,CFD)方法,越來越多的學者對浮式風機系統進行耦合數值研究。萬德成等[1]基于開源軟件OpenFOAM 開發了用于風機氣動模擬的nao-FOAM-SJTU 求解器,對浮式風機系統的氣—水動力全耦合模型進行了數值分析。TRAN T T 等[2]對DeepCWind 半潛浮式風機平臺分別用CFD 方法和勢流方法進行研究,分析了不同湍流模型對自由衰減運動的影響,并對時間步長和網格數量進行了收斂性驗證。RODDIER D 等[3]采用數值分析和模型試驗相結合的方法對風機進行了水動力分析,其中水動力分析采用TIMEFLOAT,伺服彈性計算軟件FAST 與TIMEFLOAT 耦合計算平臺運動和風力機載荷。為了驗證FOWT 動力學仿真程序的準確性,開展了OC3 和OC4 項目,對不同海上結構的模擬響應進行代碼比對[4-6]。KOO B J 等[7]開發了自己的代碼MLTSIMFAST,以分析OC4 半潛式風力發電機,并將其數值結果與DeepCwind 模型試驗結果進行比較。KIM H C 等[8]利用FAST-CHARM3D 時域全耦合動力分析程序,對5 MW OC4 FOWT 在有或無定常、動力風的隨機波浪中進行了數值模擬。REN N等[9]使用FLUENT 分析了波浪—風耦合條件下的5 MW 張力平臺型風機,并將模擬結果與試驗數據進行了驗證。ZHANG Y 等[10]使用商業計算流體力學軟件STAR-CCM+,對DeepCWind 半潛式浮動平臺與美國國家可再生能源實驗室(National Renewable Energy Laboratory,NREL) 5MW 風力發電機模型在風浪聯合激勵環境條件下進行了完全耦合的動態分析。TRAN T T 等[11]使用動態流體相互作用(Dynamic Fluid Body Interaction,DFBI) 方法和風浪激勵條件下的重疊網格技術對一個OC4 半潛式FOWT進行建模,發現CFD 結果和FAST 數據之間有很好的整體一致性,兩個方法都使用了準靜態方法對系泊線進行建模。COULLING A J 等[12]使用DeepCWind半潛浮式平臺結合5 MW 風機進行了模型試驗,同時將試驗結果與FAST 計算結果對比,結果表明FAST 在耦合的浮式風力機動力學上能夠獲取到許多相關的物理現象。此外,提出了FAST 和試驗過程的差距和改進方案,以確保浮動風機系統的精確數值建模。LIANG L 等[13]使用MATLAB 編寫了用于海上浮式風力發電機動態模擬的空氣—水動力耦合計算代碼,其中水動力是基于勢流理論模擬,且通過試驗方法進行驗證。BENITZ M A 等[14]進一步對半潛浮式平臺的六自由度體運動進行模擬,并與DeepCwind 項目的試驗結果進行全面比較。YU Z Y等[15]使用qaleFoam 混合模型[16-17]計算分析了風機葉片氣動載荷、支撐平臺水動力特性和尾跡區流場信息的耦合效應。結果表明,該混合模型可以很好地用于FOWT 耦合建模計算。YAN X K 等[18]對帶有NREL 5 MW 風機的半潛式DeepCwind 平臺設計了相應的系泊系統,然后在AQWA-OrcaFlex-FAST 中建立了完全耦合的空氣—水—伺服彈性模型,分析和比較了均勻風和動態風對平臺運動、系泊系統的有效張力、葉片變形和塔架受力的影響。此外,為了驗證具有冗余設計的系泊系統的可行性和安全性,對系泊系統的影響進行了分析。
本文首先利用計算流體學軟件STAR-CCM+對10 MW 風機和OC4-DeepCWind 5 MW 半潛浮式平臺分別進行氣動力和水動力性能的數值驗證。根據計算結果,進一步對OC4-DeepCWind 10 MW 半潛浮式風機系統耦合動力性能開展數值研究,計算了不同葉輪俯仰角對浮式風機系統發電效率及平臺運動響應的影響。
本文基于CFD 技術研究風浪聯合作用下海上浮式風機系統的流體動力性能。在笛卡爾坐標系下,對于連續、不可壓縮、非定常流體的連續性方程和N-S 方程如下。

式中,Gk表示由速度梯度引起的湍動能生成項;Gω表示ω 方程生成項;Γk和Γω表示k 和ω 的有效擴散項;Yk和Yω表示由于湍流運動引起的和的耗散項。
自由液面的捕捉采用流體體積法(Volume of Fluid,VOF),該方法引入流體體積分數α,α 表示一個單元格內某種流體所占比例,本次計算中的取值如下[21]。

本文使用滑移網格技術和重疊網格技術分別模擬風機旋轉和浮式風機系統運動。
重疊網格也稱作“嵌套”網格,用于以任意方式相互重疊的多個不同網格的計算區域,因此可以良好處理相互獨立網格之間產生的無約束位移運動[21]。任何涉及重疊網格區域的研究都具有封閉整個求解域的背景區域,以及包含域中的體的一個或更多較小區域,如圖1 所示。

圖1 重疊區域示意圖
在重疊網格中,網格單元分組為活動、不活動或受體網格單元。在活動網格單元中,將對離散控制方程進行求解。在不活動網格單元中,不會對任何方程進行求解,但是如果重疊區域在移動,這些網格單元可以變為活動狀態。
滑移網格技術具體來說,整個計算區域分為兩個子區域。葉片區域周圍的網格隨葉片移動,并與葉片保持相對靜止,而其他區域保持靜止狀態。它只需要處理局部網格的整體運動和界面的插值。它不涉及網格的變形和重疊網格,從而節省了計算時間和計算機內存[22]。
使用DFBI 模塊模擬剛性物體在流體域中壓力和剪切力作用下的運動,并考慮系泊線的恢復力。STAR-CCM+通過計算由于各種影響作用在物體上的合力和力矩,進而求解剛體運動的控制方程,以確定剛體的位置變化。圖2 顯示了說明DFBI 方法的工作流程圖。

圖2 DFBI 工作流程
本節使用STAR-CCM+軟件建立數值風洞和數值水池模型,分別對實際尺寸風機和半潛浮式風機系統開展氣動性能和水動力性能分析,完成對STAR-CCM+軟件數值建模的可靠性驗證。
2.1.1 風機參數
本風機葉片源自丹麥科技大學(Technical University of Denmark)設計的DTU 10 MW 風機[23]發電機葉片模型,如圖3 所示,風輪直徑178 m,葉片數量3 片,符合我們數值耦合研究的尺寸和要求。具體風機參數如表1 所示。

圖3 DTU 10 MW 風機模型

表1 風機模型參數
2.1.2 半潛浮式平臺參數
本文釆用美國國家可再生能源實驗室設計的5 MW OC4-DeepCwind 半潛浮式平臺對數值水池進行驗證工作[12],采用10 MW OC4-DeepCwind 半潛浮式平臺進行半潛浮式風機系統的耦合計算工作。半潛浮式平臺各部分結構是通過桁架結構進行剛性連接,風機由中部立柱支撐。按照表2 模型參數建立圖4 所示OC4-DeepCwind 半潛浮式平臺模型。

表2 半潛浮式平臺模型參數

圖4 OC4-DeepCwind 半潛浮式平臺
2.2.1 風機氣動性能驗證
根據網格收斂性驗證結果,兼顧到計算效率,最終確定網格總體數量5×106。在葉片表面設定邊界層為5 層、增長率為1.3,表面最小網格尺寸為0.008 m,網格如圖3 所示。網格運動區域采用滑移網格技術,使用SST k-ω 湍流模型,在速度入口定義恒定的均勻風速。翼型網格劃分如圖5 所示,風場區域和邊界劃分如圖6 所示。

圖5 葉片截面網格

圖6 風場區域劃分示意圖
選取了8 m/s、10 m/s、11 m/s、12 m/s、16 m/s、25 m/s 六組作業風況對風機的推力和功率進行數值驗證。在數值模擬時來流風速和葉輪轉速要保持匹配關系。表3 所示DTU 10 MW 風機在不同風速下對應的額定轉速和仰角。

表3 來流速度和風機轉速
圖7 所示不同風況下DTU 10 MW 風機推力和功率與CFD 參考值[23]的對比結果。

圖7 不同工況下風機推力和功率
風機推力結果最大誤差不超過3%,功率結果最大誤差不超過7%,驗證結果滿足后續耦合計算要求,同時證明滑移網格技術對風機氣動性能研究的可靠性。
2.2.2 自由衰減運動驗證
圖8 是半潛式平臺的網格劃分,包括自由液面區域、背景區域和重疊網格區域。考慮到計算效率和準確性的要求,網格總體數量2.8×106,表面最小網格尺寸為0.006 25 m,邊界層5 層,增長率1.3[24]。圖8 所示的是平臺的網格劃分結果。

圖8 平臺截面網格

表4 錨鏈布置及參數
系泊載荷是依據懸鏈線方程計算的準靜態模型,其中包括由于懸鏈線質量和拉伸力引起的動態效應,這意味在整個計算過程中錨鏈和半潛浮式平臺完全耦合的。在STAR-CCM+中的Catenary Coupling 模塊中可以設置準靜態懸鏈線耦合,在不考慮錨鏈躺底部分的條件下,懸鏈線耦合模塊等效于浮式風機的系泊錨鏈的作用[22]。懸鏈線懸掛在兩個端點之間,承受重力場導致的自身重量。在局部笛卡爾坐標系中,懸鏈線的形狀由以下參數方程組給定。

參數值u1和u2表示懸鏈線端點在參數空間中的位置,表達式給定。


式中,g 為重力加速度;λ0和Leq分別是無力條件下懸鏈線的單位長度質量和松弛長度;D 為懸鏈線剛度;α 和β 為積分常數,取決于兩個端點的位置和懸鏈線的總質量。錨鏈與平臺連接示意圖如圖9 所示。

圖9 平臺系泊示意圖
浮式風機系統的自由衰減能夠獲得自由運動的固有周期。本節兼顧計算效率的影響,只使用半潛浮式平臺部分計算,但模型質量屬性設置均與浮式風機系統相同。給平臺一個初始的位移,然后釋放,使其從初始位置自由移動。在靜水條件下,使用SST k-ω 湍流模型分別對平臺縱蕩、垂蕩、縱搖3 個自由度進行數值驗證。圖10 是半潛浮式平臺自由衰減運動的時歷曲線,結合下列方程計算出平臺固有周期。

式中,t1時間點對應的運動幅值為A1;t2時間點對應的運動幅值為A2;t3時間點對應的運動幅值為A3。將圖10 中的數值結果,帶入上述方程便可以計算出對應周期T。

圖10 半潛浮式平臺運動衰減時歷曲線
根據表5 中半潛浮式平臺固有周期對比結果,可以看出縱蕩的固有周期與試驗結果[9]相差較大,這主要是由于準靜態的懸鏈線模型在計算過程中是分段迭代求解,縱蕩的位移遠大于縱搖和垂蕩,準靜態模型在迭代求解的過程中會積累較大誤差,但是縱蕩誤差沒超過14%;垂蕩和縱搖與試驗值、FAST 計算結果[9]均小于5%,在誤差允許的范圍內。

表5 半潛浮式平臺固有周期對比
本節著重考慮風浪聯合作用下,浮式風機系統的運動響應。計算過程中并未考慮塔柱對流場的影響,但考慮到了整體系統的質量屬性。選定額定計算工況,其風速為11 m/s,波高為3.66 m,波周期為9.7 s,風機葉片的額定轉速為8.837 rpm。根據相似理論[25-26],采用1 ∶45 的縮尺模型進行計算。考慮到實際工程應用中,風力機有一個范圍在0°至5°的預置仰角,所以本文通過改變風機葉輪俯仰角在0°和5°情況下,在風浪聯合條件下對風機系統在水動力和氣動力性能進行分析。
為了確保浮式系統能夠準確地在波浪模擬中運動,本節所需的規則波工況是由計算域邊界直接輸入速度和壓力完成的,STAR-CCM+的直接輸入造波方法是基于波浪理論的精準方法,在開始數值計算之前對數值造波的準確性進行了驗證。采用縮尺工況為周期1.45 s、波高0.081 3 m 的規則波,數值水池后端設有兩倍波長消波區域。通過對比數值解和解析解(圖11),計算發現其最大誤差率仍低于5%,說明數值造波的準確度較高,可以用于后續的數值計算。

圖11 波幅時歷曲線
結合前面章節中網格劃分的方案,整體計算域長14 m、寬9 m、高10 m,整個耦合系統中心距離速度入口5 m、速度出口9 m,距離頂部邊界3 m、底部邊界7 m。在風機葉片和浮式平臺周圍分別再創建網格加密區域來捕捉附近復雜的不穩定流動。整個網格區域共分成三大部分,即背景區域、風機葉片組成的滑移網格區域,以及浮式風機系統整體組成的重疊網格區域。風和波浪通過在速度入口施加速度和壓力值實現。同時在水池前端施加波浪力,尾端施加阻尼消波,降低計算過程中背景波浪的影響。風機系統由DFBI 驅動重疊網格實現自由運動,風機葉片采用滑移網格技術套嵌在重疊內部,通過給定轉速驅動葉片旋轉。圖12 所示耦合區域及邊界劃分,圖13 所示耦合模型的計算流場。

圖12 耦合區域及邊界

圖13 耦合計算的流場
本節預設風機葉輪仰角0°和5°兩種結構,對風機的推力和扭矩進行了檢測,根據扭矩數值結果進一步計算出風機功率,推力和功率結果如圖14 所示。

圖14 風機推力和功率時歷曲線
根據圖14,我們可以發現隨著角度的增大,其推力和功率也相應下降,5°俯仰角的推力均值為9.45 N,功率均值為3.4 W;0°俯仰角推力均值為10.2 N,功率均值為4.4 W。5°俯仰角相較于0°俯仰角推力減小7.35%,功率減小25%。考慮到風載荷會對風機造成一個初始的縱搖角度,所以風機工作時的俯仰角都需要在設定角度的基礎上,再加上縱搖均值角度。由于風機的受風面積隨著葉輪俯仰角度的增大而減小,所以風機受到的推力減小,自然導致了整體功率的下降。
本文在數值模擬過程中,對浮式風機系統3 個自由度運動(縱蕩、縱搖、垂蕩)及系泊張力進行了檢測。由于風載荷和波浪載荷帶來的瞬態效應,使平臺在計算前期產生了較大非線性運動,因此本文主要針對數值計算后期較為平穩階段進行研究分析。
根據圖15、圖16 和圖17 時歷曲線可以看出隨著角度增加,幅值和均值都減小了。主導因素是葉輪角度改變導致風載荷對葉盤作用面積減小。同時由于波浪載荷帶來的瞬態效應,使平臺在計算前期產生了較大周期性的振蕩,因此下面主要對計算后期較為平穩階段進行分析研究。相較于0°,俯仰角為5°時的縱蕩幅值減少了5.26%,幅值均值減少了3.56%;縱搖幅值減少了1.22%,均值差別不大;但是垂蕩方向上幅值變化不大,幾乎一致。

圖15 縱蕩運動時歷曲線

圖16 垂蕩運動時歷曲線

圖17 縱搖運動時歷曲線
根據圖18 錨鏈1 張力時歷曲線可以看出,5°俯仰角比0°俯仰角的錨鏈響應幅值和均值都減少了0.43%。根據圖19 錨鏈2 張力時歷曲線可以看出,5°俯仰角比0°俯仰角幅值減少了8.2%,均值增加了1.8%。錨鏈的運動響應和縱蕩運動關聯較大,仰角的增加勢必會減小錨鏈的響應波動,因此適當增加俯仰角度有利于提高錨鏈壽命錨鏈幅值的波動容易影響錨鏈的整體使用壽命,故適當增大俯仰角可以提高錨鏈的壽命。

圖18 錨鏈1 張力時歷曲線

圖19 錨鏈2 張力時歷曲線
本文基于浮式風機耦合數值研究的現狀,使用軟件STAR-CCM+驗證了DTU-10MW 發電風機的推力和功率,以及對DeepCwind 半潛式浮動平臺進行了固有周期的驗證。在此研究的基礎上,通過對比不同葉輪俯仰角對半潛浮式風機系統的水動力和氣動力性能進行分析,得到如下結論。
(1)基于SST k-ω 湍流模型、滑移網格技術、重疊網格技術及準靜態懸鏈線系泊模塊建立的數值模型,可以很好地滿足半潛浮式風機系統的動態耦合模擬。
(2)根據動態耦合模擬結果,隨著風機葉輪俯仰角度從0°到5°的改變,整個浮式風機系統在穩性上沒有明顯的變化;但在風機功率上,相較于推力的減小,風機功率減小比例更大。
(3)風機葉輪俯仰角度增加后,錨鏈1 和錨鏈2 的張力波動均減小,但錨鏈響應周期沒有明顯變化,所以在工程意義上適當調節俯仰角可以提高錨鏈的使用壽命。