鄧露 ，劉建祥

（1.工程結(jié)構(gòu)損傷診斷湖南省重點實驗室，湖南長沙410082；2.湖南大學(xué)土木工程學(xué)院，湖南長沙410082）

日益增長的能源需求促進了海上風(fēng)電行業(yè)的發(fā)展，半潛式浮式風(fēng)機平臺適用水深范圍廣，是一種開發(fā)具有儲量豐富、風(fēng)速大等優(yōu)點的深水風(fēng)能資源的重要海洋平臺[1].長期以來，柱體繞流一直是海洋工程領(lǐng)域的重點問題[2].流體流經(jīng)浮式風(fēng)機平臺的柱體結(jié)構(gòu)時，柱體后方周期性脫落的漩渦會導(dǎo)致平臺大幅往復(fù)運動，這將加劇平臺系泊結(jié)構(gòu)的疲勞損傷，甚至導(dǎo)致系泊結(jié)構(gòu)疲勞破壞[3].因此開展半潛式浮式風(fēng)機平臺繞流研究有著重要的工程意義.

在單柱繞流方面，國內(nèi)外學(xué)者針對單柱繞流的力學(xué)特性、流場特性及其影響因素等已進行了大量的研究，并取得了豐碩的成果[4-7].在雙柱和多柱式結(jié)構(gòu)繞流方面，賈曉荷和劉樺[8]采用大渦模擬方法分析了并列或串列排布的雙圓柱的繞流特性，結(jié)果表明串列排布時，上游圓柱對下游圓柱受力有影響，并列排布時，兩柱的阻力系數(shù)均值基本相等.鄭宇華和顧杰[9]的研究表明，間距比（雙柱中心的間距與柱寬的比值）影響并列雙柱開縫下方的漩渦個數(shù)、形態(tài)等.Kitagawa T 和Ohtab H[10]采用大渦方法模擬了串列雙圓柱的繞流特性，結(jié)果表明當(dāng)間距比小于某一臨界值時，只有下游圓柱有渦脫現(xiàn)象.李聰洲等[11]研究了串列雙柱繞流的力學(xué)特性，結(jié)果表明存在一個臨界間距比使得圓柱的升力和阻力發(fā)生跳躍性變化.劉為民等[12]的研究表明，0°和 45°來流時，陣列四立柱柱群中立柱的渦脫頻率大于單柱的渦脫頻率，渦脫模式表現(xiàn)為2S 模式.劉正浩等[13]的研究表明，來流角對陣列四柱式海洋平臺的立柱的升力和阻力系數(shù)影響較大，不同來流角下平臺底部的流場不同.

雙柱和多柱式結(jié)構(gòu)繞流時，柱體尾流之間的相互干擾會導(dǎo)致立柱的力學(xué)特性和流場結(jié)構(gòu)更復(fù)雜[11].然而，以上針對雙柱和多柱式結(jié)構(gòu)繞流的研究，主要分析了流速、來流角和間距比等因素對結(jié)構(gòu)的受力、流場形態(tài)和渦脫模式等的影響，而深入分析流場結(jié)構(gòu)以及尾流之間干擾機理的研究極少.因此，本文采用DDES 方法對不同來流角和流速下多柱式半潛式浮式風(fēng)機平臺的繞流進行了三維數(shù)值模擬，分析了平臺立柱的阻力系數(shù)平均值、升力系數(shù)頻率等力學(xué)指標(biāo)，從流場的相干結(jié)構(gòu)層面研究了尾流的干擾機理，并分析了尾流的空間相關(guān)性.研究結(jié)果闡釋了多柱式結(jié)構(gòu)繞流的尾流相互干擾機理，為深入理解多柱式結(jié)構(gòu)繞流現(xiàn)象提供了理論依據(jù).

1 數(shù)值建模

1.1 DDES 方法

Spalart 等[14]提出的 DES 方法（Detached Eddy Simulation）是一種雷諾時均法（Reynolds-averaged Navier-Stokes，RANS）和大渦模擬（Large Eddy Simulation，LES）相結(jié)合的方法，其基本處理原則是：在近壁區(qū)采用雷諾時均法進行模擬，在遠離壁面的區(qū)域采用LES 來模擬分離流動.DES 方法通過長度尺度來實現(xiàn)RANS 到LES 的轉(zhuǎn)變，其定義為：

式中：d 為離最近壁面的距離，Cdes取 0.65，Δmax為局部網(wǎng)格的最大尺寸.

然而，DES 方法具有局限性，在某些情況下，當(dāng)RANS 轉(zhuǎn)變到LES 時，會出現(xiàn)模型應(yīng)力損耗現(xiàn)象，并導(dǎo)致“網(wǎng)格誘導(dǎo)分離”問題.因此，Spalart 等[15]對DES方法作了改進，提出了含有延遲函數(shù)的DDES 方法，解決了“網(wǎng)格誘導(dǎo)分離”問題.DDES 方法通過引入延遲函數(shù)對長度尺度的定義進行修正，修正后的定義為：

式中：fd為延遲函數(shù)，其表達式為：

式中：rd為延遲因子，其表達式為：

式中：κ 為卡門常數(shù)，取 0.41；ui，j為速度梯度；vt為動力黏性系數(shù)；v 為分子黏性系數(shù).

DDES 方法是一種改進的RANS 與LES 相結(jié)合的方法，其在模擬流體分離流動時十分有效[16]，已得到廣泛的應(yīng)用.該方法具有多個優(yōu)點，如不需要精細的網(wǎng)格也能得到較好的結(jié)果，適合復(fù)雜的結(jié)構(gòu)等.因此，本文采用該方法進行研究.

1.2 平臺模型

選取美國國家可再生能源實驗室設(shè)計的半潛式浮式風(fēng)機平臺作為分析對象，平臺的立柱（包括邊柱和中柱）、連桿和浮箱底座的結(jié)構(gòu)布置如圖1 所示（圖中括號內(nèi)為構(gòu)件編號）.平臺的結(jié)構(gòu)尺寸，吃水深度等設(shè)計參數(shù)詳見文獻[17].

圖1 半潛式浮式風(fēng)機平臺模型圖Fig.1 Model of semi-submersible floating wind turbine platform

1.3 計算域與網(wǎng)格劃分

流場計算域及模型平面布置如圖2 所示.θ 為來流角，即中柱與邊柱2 的連線與x 軸負向的夾角.流域入口邊界距離邊柱2 中心為10D，出口邊界距離邊柱2 中心為25D，左右邊界距離中柱中心均為10D.D 為邊柱直徑.

圖2 計算域與平面布置圖Fig.2 Computational domain and plane configuration

邊界條件設(shè)定為：入口為速度入口，出口為出流邊界，頂部為自由邊界，其余邊界為對稱邊界，平臺表面采用無滑移固定壁面.

網(wǎng)格劃分示意圖如圖3.由于平臺幾何結(jié)構(gòu)復(fù)雜，因此，將流場分成兩種類型的網(wǎng)格進行劃分，在平臺附近區(qū)域采用四面體網(wǎng)格以適應(yīng)復(fù)雜的幾何結(jié)構(gòu)，在其他區(qū)域采用六面體網(wǎng)格以加快計算速度.為控制網(wǎng)格密度，將流域分成九個區(qū)域，在平臺附近和尾流區(qū)劃分密集網(wǎng)格，其他區(qū)域則劃分較稀疏的網(wǎng)格.根據(jù)張新曙等[18]的建議，經(jīng)過驗算，并綜合考慮計算耗時成本，確定網(wǎng)格方案為：加密區(qū)（平臺四周1.4D 的區(qū)域）的網(wǎng)格大小為0.12D；尾流區(qū)的初始網(wǎng)格大小為0.25D，且隨著尾流遠離平臺，網(wǎng)格逐漸稀疏.

圖3 模型網(wǎng)格劃分示意圖Fig.3 Mesh generation of model

1.4 計算設(shè)置

流域采用Ansys Fluent 進行求解.時間離散采用二階隱式格式，空間離散采用二階迎風(fēng)格式，壓力與速度耦合方程采用PISO 方法（Pressure-Implicit with Splitting of Operators）.時間步長取0.05 s，每個時間步長內(nèi)最大迭代次數(shù)為20，模擬時長為6 000 ～10 000 s，以保證獲得穩(wěn)定的結(jié)果.數(shù)據(jù)處理時，選取最后約0.6 倍模擬時長內(nèi)較穩(wěn)定的計算結(jié)果進行處理.

1.5 計算工況

本文選取 3 種典型來流角：θ=0°、30°、60°.根據(jù)中國南海的實際情況[19]，選取12 種來流速度：u=0.32 m/s、0.42 m/s、0.53 m/s、0.63 m/s、0.74 m/s、0.85 m/s、0.95 m/s、1.06 m/s、1.16 m/s、1.27 m/s、1.37 m/s、1.48 m/s，來流方向沿x 軸正向.

2 計算結(jié)果及分析

2.1 力學(xué)特性

流體流經(jīng)半潛式浮式風(fēng)機平臺時，在水平面內(nèi)，平臺立柱受到兩個力的作用，一個是順流向的阻力Fx，又叫拖曳力，另一個是橫流向的升力Fy.定義阻力系數(shù)為，升力系數(shù)為其中，u 為來流速度，ρ 為海水密度，取為1 025 kg/m3，A 為柱體迎流投影面積.

2.1.1 阻力系數(shù)平均值

圖4 ～圖7 分別為中柱、邊柱1、邊柱2、邊柱3的阻力系數(shù)平均值隨來流角和流速變化的關(guān)系圖.從這四幅圖中可知，在0°來流下，中柱的阻力系數(shù)平均值隨流速的變化而小幅波動，其余立柱的阻力系數(shù)平均值隨流速的改變而基本不變.在30°來流下，邊柱2 和邊柱3 的阻力系數(shù)平均值隨流速的改變而基本不變，但中柱和邊柱1 的阻力系數(shù)平均值隨流速的改變而變化較大，離散性較強.這是因為流體流經(jīng)上游立柱后發(fā)展成為湍流，湍流本身具有較強的隨機波動性，所以，受上游尾流的影響，下游中柱和邊柱1 的阻力系數(shù)平均值隨流速的變化而波動較大[12].在60°來流下，中柱、邊柱2 和邊柱3 的阻力系數(shù)平均值隨流速的改變而基本不變，而邊柱1 的阻力系數(shù)平均值出現(xiàn)了跳躍波動的現(xiàn)象.這也是因為邊柱1 處于下游位置，其阻力系數(shù)平均值受上游尾流的影響較大.

圖4 中柱的阻力系數(shù)平均值Fig.4 Average drag force coefficient of middle column

圖5 邊柱1 的阻力系數(shù)平均值Fig.5 Average drag force coefficient of side column 1

圖6 邊柱2 的阻力系數(shù)平均值Fig.6 Average drag force coefficient of side column 2

圖7 邊柱3 的阻力系數(shù)平均值Fig.7 Average drag force coefficient of side column 3

由圖4 可知，中柱的阻力系數(shù)平均值隨來流角的增大而減小.阻力是漩渦脫落引起的立柱前后表面的壓力差和立柱表面的摩擦力共同作用所形成的.隨著來流角的改變，立柱和立柱、立柱和底座在垂直于來流方向的相對距離發(fā)生了改變.由圖2 的平面布置圖可知，在0°來流時，中柱后方尾流距離立柱和底座較遠，下游柱體的阻擋作用小，上游中柱的尾流漩渦能夠充分地發(fā)展，所以阻力系數(shù)最大.在30°來流時，邊柱1 和底座1 靠近中柱尾流，下游柱體的阻擋作用增大，下游中柱后方的漩渦脫落不能充分發(fā)展，故阻力系數(shù)減小.在60°來流時，邊柱1 和底座1位于中柱的正后方，完全阻擋了中柱尾流漩渦的發(fā)展，故阻力系數(shù)最小.由圖5 可知，邊柱1 的阻力系數(shù)平均值隨來流角的增大呈現(xiàn)先減小后增大的趨勢.由圖6 和圖7 可知，邊柱2 和邊柱3 的阻力系數(shù)平均值均隨來流角的增大而增大.這些也均是因為隨著來流角的改變，柱體之間在垂直于來流方向的相對距離發(fā)生了改變，柱體與柱體之間的相互影響所導(dǎo)致的.

2.1.2 升力系數(shù)時程曲線及頻譜分析

限于篇幅，僅給出了 θ=0°、u=1.06 m/s 時各個立柱的升力系數(shù)時程曲線，如圖8 所示.從圖8 中可以看出，邊柱2、邊柱3 和中柱的升力系數(shù)時程變化表現(xiàn)出較好的諧波特性，而邊柱1 的升力系數(shù)變化波動較大.邊柱2 和邊柱3 的升力系數(shù)變化相反，中柱的升力系數(shù)變化滯后于邊柱2.當(dāng)邊柱2 的升力系數(shù)達到波峰時，中柱的升力系數(shù)稍有增大，且曲線變化趨勢發(fā)生了改變.以圖8 中t1-t2時間段為例，此時間段內(nèi)中柱的升力系數(shù)未按圖中粗劃虛線趨勢變化，而是按實線趨勢變化.這可能是因為中柱處于邊柱2 的尾流中，當(dāng)上游邊柱2 的升力系數(shù)達到波峰時，其后方的漩渦脫落，隨后撞擊到中柱上，導(dǎo)致中柱的升力系數(shù)稍微增大而改變了曲線的變化趨勢.

圖8 0°來流角下各柱升力系數(shù)時程曲線Fig.8 Time history curves of lift force coefficients of four columns under 0°heading

圖9 ～圖11 為u=1.06 m/s 時，三種來流角下平臺各立柱的升力系數(shù)頻譜圖.立柱后方漩渦周期性地脫落是導(dǎo)致升力系數(shù)周期性變化的原因.由圖9可知，0°來流角下，中柱的頻譜除主峰外，還存在多個高頻的小峰，小峰的幅度遠小于主峰的幅度，說明中柱后方漩渦脫落以主峰頻率為主，并伴隨多個高頻成分.邊柱1 的頻譜圖中存在兩個較大的峰值，低頻對應(yīng)的幅度顯著大于高頻對應(yīng)的幅度，說明邊柱1 后方低頻的渦脫引起了較大的升力.邊柱2 和邊柱3 的頻譜圖中均只有一個峰值，說明這兩立柱后方的的漩渦為單頻脫落.此外，4 個立柱的頻譜圖中均有一峰值對應(yīng)的頻率為f=0.018 5 Hz.綜上可知，平臺立柱后方漩渦脫落的主頻率為f=0.018 5 Hz，部分立柱還伴有高頻或低頻成分.

由圖10 可知，30°來流角下，中柱、邊柱 2 和邊柱3 的升力系數(shù)頻譜雖然有多個峰值，但對應(yīng)的頻率極低.受上游尾流的干擾，下游邊柱1 的升力系數(shù)頻譜沒有明顯的峰值.可見在此來流角下，平臺立柱后方漩渦脫落不規(guī)則.

圖9 0°來流角下各立柱升力系數(shù)頻譜Fig.9 Lift coefficient spectrum of each column under 0°heading

圖10 30°來流角下各立柱升力系數(shù)頻譜Fig.10 Lift coefficient spectrum of each column under 30°heading

圖11 60°來流角下各立柱升力系數(shù)頻譜Fig.11 Lift coefficient spectrum of each column under 60°heading

由圖11 可知，60°來流角下，中柱頻譜沒有峰值，這是因為流體繞過中柱后，撞擊在下游的底座1和邊柱1 上，沒有形成周期性的渦脫.其余立柱的頻譜均有多個峰值，但沒有相同頻率的成分.由此可見，平臺立柱后方的漩渦呈多頻脫落模式.其他流速下，各柱的升力系數(shù)頻譜也具有類似的性質(zhì)，限于篇幅，不再贅述.

2.2 流場分析

2.2.1 流場相干結(jié)構(gòu)與尾流干擾

湍流場中有組織性的漩渦結(jié)構(gòu)稱之為相干結(jié)構(gòu)，其在湍流的發(fā)展演化中起著重要的作用[20-21].u=1.06 m/s 時，0°、30°來流角下的渦量等值面圖（按 x方向速度 u 著色）分別如圖12（a）、（b），渦量等值面采用 Hunt 等[22]建議的 Q 判據(jù)（Q-criterion）表示.從圖12 中可以看出，流體流經(jīng)平臺后，形成了多種形式和尺度不一的漩渦.圖12（a）和（b）中的弧狀渦包成月牙形分布在浮箱底座的底面和頂面邊緣，反映出流體流經(jīng)浮箱底座時，在底面和頂面邊緣發(fā)生流動分離，并在隨后的區(qū)域形成了馬蹄渦.圖12（a）中，邊柱1、邊柱2 和邊柱3 后方有成對的、錯位排布的管狀渦包，說明這些立柱后方形成了交替脫落、錯位排布的尾渦，尾渦流經(jīng)中柱后，逐漸拉長而發(fā)展成為流向渦.從圖12（b）中可以看出，邊柱3 后方形成了斜向下的錐形渦包，說明該處漩渦結(jié)構(gòu)的半徑隨著水深的增加而逐漸縮小.浮箱底座3 后方呈螺旋狀上升的螺旋渦包，反映出其后方形成了發(fā)夾渦[21]，該漩渦結(jié)構(gòu)由渦頭和兩條渦腿組成，渦頭受到渦腿旋轉(zhuǎn)向上的抬升力而上揚.靠近浮箱底座的發(fā)夾渦的尺度明顯小于遠離底座的發(fā)夾渦的尺度，說明隨著湍流的發(fā)展，發(fā)夾渦逐漸發(fā)展壯大.

流體流經(jīng)平臺時，平臺的立柱、連桿和浮箱底座的尾流之間會相互干擾.如圖12（a）中所示，邊柱2后方尺度較大的漩渦結(jié)構(gòu)與下游斜桿2 后方尺度較小的漩渦結(jié)構(gòu)相互融合，最終匯合成為一個整體，相

當(dāng)i=j 時，表示該點的脈動量自相關(guān)，當(dāng)i ≠j時，表示兩點的脈動量互相關(guān).相關(guān)系數(shù)的大小表征了脈動量在流場不同位置的關(guān)聯(lián)程度，R=±1 表示完全相關(guān)，R=0 表示不相關(guān).

圖13 為0°來流角下y=0 平面結(jié)構(gòu)的剖面圖.為考察上下游立柱尾流間的空間相關(guān)性，選擇x 方向速度u 作為脈動量，計算了y=0 平面，z=1.17D高度處邊柱2（bz2）尾流點與中柱（zz）尾流點的相關(guān)系數(shù)，計算結(jié)果如圖14 所示，其中，圖14（b）為圖14（a）的水平投影圖.當(dāng)于大漩渦“吞并”了小漩渦.如圖12（b）中所示，浮箱底座3 后方形成的發(fā)夾渦在上揚的過程中，由于渦腿的旋轉(zhuǎn)效應(yīng)，誘導(dǎo)了邊柱3 的尾渦向斜下方運動，最終與浮箱底座3 的尾渦匯合在一起，并導(dǎo)致了次生流向渦的形成.此外，由于該發(fā)夾渦渦腿的旋轉(zhuǎn)效應(yīng)，也誘導(dǎo)了中柱后方的流向渦發(fā)生偏轉(zhuǎn)彎折，從而改變了中柱尾流的流動方向.由以上分析可知，尾流干擾的機理正是相干結(jié)構(gòu)之間的相互作用.

2.2.2 空間相關(guān)性分析

流場內(nèi) i、j 兩點的脈動量 Xi、Xj在空間上的相關(guān)性可用相關(guān)系數(shù)R 來表示，其定義[23]表示為：

圖13 y=0 平面結(jié)構(gòu)剖面圖Fig.13 Stucture section at y=0

由圖14 可知，邊柱2 尾流與中柱尾流的x 方向速度在一定區(qū)間內(nèi)具有較強的相關(guān)性.經(jīng)計算，x 坐標(biāo)在 [0.85D，1.1D]區(qū)段的邊柱2 尾流與x 坐標(biāo)在[3.1D，3.4D]區(qū)段的中柱尾流的x 方向速度的相關(guān)系數(shù)為0.4 ＜R ＜0.55，表現(xiàn)出較強的正相關(guān)性，而x 坐標(biāo)在[1.7D，2D]區(qū)段的邊柱2 尾流與x 坐標(biāo)在[3.1D，3.4D]區(qū)段的中柱尾流的x 方向速度的相關(guān)系數(shù)為-0.55 ＜R ＜-0.4，表現(xiàn)出較強的負相關(guān)性.這說明雖然邊柱2 后方交替脫落的尾渦流經(jīng)中柱后，變成了流向渦，但兩立柱的尾流之間仍有著較強的空間相關(guān)性.

圖14 x 方向速度u 的相關(guān)系數(shù)分布圖Fig.14 Correlation coefficient distribution of velocity u in x direction

3 結(jié) 論

本文通過 DDES 方法研究了 0°、30°和 60°來流下半潛式浮式風(fēng)機平臺繞流的力學(xué)特性，并對平臺的尾流進行了深入分析，得出以下結(jié)論：

1）由于尾流的干擾，下游立柱的阻力系數(shù)平均值隨流速變化而波動，離散性大，而上游立柱的阻力系數(shù)平均值變化不大.在0°來流時，平臺立柱后方的漩渦脫落以某一主頻率為主，并伴有其他高頻或低頻成分；60°來流時，平臺立柱后方的漩渦脫落呈多頻模式.

2）流場相干結(jié)構(gòu)具有多尺度性，且類型多樣，包括發(fā)夾渦、流向渦、次生流向渦、馬蹄渦等.相干結(jié)構(gòu)之間的相互作用是尾流干擾的內(nèi)在原因.尾流干擾存在兩種模式：處于上下游位置的立柱與連桿之間的尾渦相互交融，導(dǎo)致它們的尾流成為一個整體；浮箱后的發(fā)夾渦誘導(dǎo)立柱漩渦改變運動方向，導(dǎo)致立柱尾流匯入浮箱尾流或者運動發(fā)生偏轉(zhuǎn).

3）流體流經(jīng)成串列排布的立柱時，即使流場的相干結(jié)構(gòu)發(fā)生了轉(zhuǎn)變，上下游立柱尾流之間仍有著較強的空間相關(guān)性.