不同雷諾數(shù)下傾斜圓柱繞流三維數(shù)值模擬研究

高洋洋，張演明，劉 彩，王 濱

(1. 浙江大學(xué) 海洋學(xué)院，浙江 舟山 316021; 2. 中國(guó)電建集團(tuán)華東勘測(cè)設(shè)計(jì)研究院有限公司，浙江 杭州 311122)

在海洋工程中，鋼懸鏈線立管因其結(jié)構(gòu)簡(jiǎn)單、經(jīng)濟(jì)可靠及對(duì)上部平臺(tái)運(yùn)動(dòng)更好的適應(yīng)性越來越多的被應(yīng)用于深水油氣資源開發(fā)中。鋼懸鏈線立管近觸底區(qū)由于傾斜角度變化大，在流體繞流傾斜立管的過程中，尾渦相互干擾機(jī)制更加復(fù)雜，更易導(dǎo)致立管發(fā)生疲勞破壞。然而目前基于豎直圓柱的尾渦特性研究因忽略軸向流影響，其尾流形態(tài)特征、受力特性及旋渦脫落頻率與豎直圓柱相比存在明顯差異。近年來，針對(duì)傾斜圓柱繞流的研究受到海洋工程界和學(xué)術(shù)界的廣泛關(guān)注。

Shirakasi等[1]在傾斜圓柱繞流的風(fēng)洞試驗(yàn)中發(fā)現(xiàn)了沿圓柱軸向的渦運(yùn)動(dòng)，且隨傾角的增大旋渦脫落頻率逐漸減小。李壽英和顧明[2]通過對(duì)升力系數(shù)進(jìn)行功率譜密度分析發(fā)現(xiàn)豎直圓柱與傾斜圓柱繞流存在顯著差異，在傾斜圓柱尾流中存在低頻諧波。Hogan & Hall[3]研究表明隨著傾斜角度的增大，傾斜圓柱繞流流場(chǎng)的三維效應(yīng)愈發(fā)明顯，不存在規(guī)則脫落的旋渦。Franzini等[4]通過物理模型試驗(yàn)進(jìn)一步探討了圓柱傾斜角度對(duì)尾流渦街的影響，研究發(fā)現(xiàn)當(dāng)傾斜角度達(dá)到α=45°時(shí)，圓柱后方的軸向流會(huì)對(duì)圓柱表面自由剪切層之間的相互作用造成一定影響，使得尾流的旋渦脫落變得不規(guī)則。Mstsumoto等[5-6]指出，圓柱后方軸向流的作用類似于分流板，即軸向流會(huì)影響自由剪切層的相互作用，降低旋渦脫落的強(qiáng)度。為了簡(jiǎn)化計(jì)算，通常采用獨(dú)立性原則(independence principle，簡(jiǎn)稱IP)進(jìn)行傾斜圓柱繞流特性研究，該原則認(rèn)為，在不同傾斜角度下，通過垂直于圓柱軸線的速度分量計(jì)算得到的升阻力系數(shù)及旋渦脫落頻率趨于一致[7]。Vakil & Green[8]開展了低雷諾數(shù)下不同長(zhǎng)細(xì)比及傾斜角度條件下圓柱繞流的數(shù)值模擬研究，研究發(fā)現(xiàn)當(dāng)傾斜角度大于45°時(shí)，獨(dú)立性原則不再適用。Zhao等[9]研究發(fā)現(xiàn)當(dāng)雷諾數(shù)達(dá)到Re=1 000時(shí)，在小傾角情況下，阻力系數(shù)均值與無量綱渦脫頻率St數(shù)適用于獨(dú)立性原則，相反，升力系數(shù)均方根值不符合獨(dú)立性原則。Thapa等[10]開展了近壁面傾斜圓柱繞流的數(shù)值模擬研究，研究結(jié)果表明由于邊界效應(yīng)的影響，隨著傾斜角度的不斷增大，流場(chǎng)的三維特性逐漸減弱。梁輝等[11]研究了傾斜圓柱繞流的展向特性，歸納了不同展向截面處的不同渦脫形態(tài)。Hu等[12]采用高速粒子測(cè)速技術(shù)(PIV)開展了方柱繞流的試驗(yàn)研究，發(fā)現(xiàn)方柱后側(cè)的軸向流會(huì)對(duì)尾流流場(chǎng)特征產(chǎn)生顯著的影響，隨著傾斜角度的增大，順流向方柱出現(xiàn)自底端流向頂端的軸向流動(dòng)，而逆流向方柱出現(xiàn)自頂端流向底端的軸向流動(dòng)。

綜上所述，當(dāng)前針對(duì)圓柱繞流的研究主要集中于豎直圓柱繞流及單一雷諾數(shù)條件下同一傾斜方向圓柱繞流特性研究，不同傾斜方向及雷諾數(shù)對(duì)傾斜圓柱繞流的流場(chǎng)形態(tài)及受力特性的影響規(guī)律尚未闡明。因此，采用計(jì)算流體力學(xué)開源代碼OpenFOAM開展了不同雷諾數(shù)、傾斜方向及傾斜角度條件下傾斜圓柱繞流的三維數(shù)值模擬研究，主要分析了瞬時(shí)及時(shí)均三維流場(chǎng)形態(tài)、流線拓?fù)浼拜S向流分布特征，揭示了傾斜圓柱尾流特征、受力特性及旋渦脫落頻率隨傾斜角度、雷諾數(shù)及傾斜方向變化的影響規(guī)律，闡明了不同雷諾數(shù)條件下獨(dú)立性原則(IP)的適用范圍。

1 數(shù)學(xué)模型

1.1 控制方程

在直角坐標(biāo)系中，不可壓縮黏性流體在層流狀態(tài)下的連續(xù)方程和動(dòng)量守恒方程N(yùn)avier-Stokes(N-S)方程為：

·U=0

(1)

(2)

其中，U，p，ρ和υ分別表示速度矢量、壓力、流體密度和流體動(dòng)力黏度。

大渦模擬控制方程可由湍流狀態(tài)下將N-S方程進(jìn)行濾波處理得到[13-17]，如式(3)-(4)所示：

(3)

(4)

利用有限體積法進(jìn)行離散。在層流工況下采用icoFoam算法進(jìn)行速度與壓力耦合，壓力與拉普拉斯項(xiàng)采用二階高斯線性正交格式進(jìn)行離散，時(shí)間項(xiàng)采用一階隱式進(jìn)行離散，梯度與散度項(xiàng)采用高斯線性進(jìn)行離散。在湍流工況下采用pisoFoam算法進(jìn)行速度與壓力耦合，壓力項(xiàng)和拉普拉斯項(xiàng)采用二階高斯線性正交格式進(jìn)行離散，時(shí)間項(xiàng)采用二階隱式進(jìn)行離散，梯度與散度項(xiàng)采用高斯線性進(jìn)行離散。

1.2 計(jì)算域及邊界條件設(shè)置

傾斜圓柱繞流流場(chǎng)計(jì)算域如圖1所示，計(jì)算域長(zhǎng)度和寬度分別取L=30D和W=20D，其中D表示圓柱直徑，取D=0.01 m。圓柱底部橫截面的中心與計(jì)算域入口的距離為10D，同時(shí)根據(jù)Behr等[18]的研究結(jié)果，圓柱距離兩側(cè)壁的距離設(shè)置為10D。計(jì)算域展向高度設(shè)為Z=πD，與文獻(xiàn)[19-20]保持一致，從而保證流場(chǎng)在展向得到充分發(fā)展，準(zhǔn)確捕捉傾斜圓柱的三維渦脫。圓柱傾角α表示來流速度與圓柱法線的夾角，如圖1(b)所示。整個(gè)計(jì)算域的初始速度與壓力均設(shè)置為0，入口設(shè)置為自由來流邊界，入口速度設(shè)置分別為u=u0=0.01 m/s、0.15 m/s和0.39 m/s，對(duì)應(yīng)不同的雷諾數(shù)Re=100、1 500和3 900，而橫流向與展向速度設(shè)置為v=w=0；出口設(shè)置為自由出流邊界條件；兩側(cè)壁面設(shè)置為對(duì)稱邊界，即法向速度與壓力的梯度為0；上下兩壁面設(shè)置為零梯度邊界，即所有變量沿法向梯度為0；圓柱表面設(shè)置為無滑移壁面，即u=v=w=0。具體邊界條件設(shè)置詳見表1。在層流工況下，基于前人關(guān)于圓柱繞流研究[21-23]，雷諾數(shù)采用Re=100，在該雷諾數(shù)下圓柱繞流流場(chǎng)呈現(xiàn)明顯的二維特性。在湍流工況下，前人關(guān)于圓柱繞流流場(chǎng)的實(shí)驗(yàn)及數(shù)值模擬研究大多采用雷諾數(shù)Re=1 500[24-26]和3 900[4,13-14]，當(dāng)Re=1 500和3 900時(shí)，湍流流場(chǎng)均得到充分發(fā)展。因此，在本研究中采用直接模擬方法和Smargorinsky-Lilly亞格子模式的大渦模擬(LES)湍流模型開展Re=100、1 500和3 900下傾斜圓柱繞流研究，分析不同雷諾數(shù)對(duì)傾斜圓柱繞流流場(chǎng)的影響。傾斜角度范圍為-60° ≤α≤ 60°，每間隔15°設(shè)為一組，共9組工況，其中正值和負(fù)值分別表示圓柱沿順流向傾斜和逆流向傾斜。

圖1 計(jì)算域和傾斜角度αFig. 1 Computational domain and definition of the inclined angle α

表1 邊界條件設(shè)置
Tab. 1 Boundary conditions

邊界條件設(shè)定入口自由來流u=u0=0.01 m/s、0.15 m/s和0.39 m/s,v=w=0出口自由出流u/x=0, P=0兩側(cè)邊界對(duì)稱邊界?u/?y=0,?w/?y=0,v=0;?P/?y=0上下邊界零梯度邊界u/z=0,P/z=0圓柱表面無滑移邊界u=v=w=0

1.3 網(wǎng)格及時(shí)間步驗(yàn)證

(5)

(6)

St=fsD/u

(7)

(8)

圖2 傾斜圓柱繞流流場(chǎng)計(jì)算網(wǎng)格(α=30°)Fig. 2 Schematic view of the computational mesh(α=30°)

表2 豎直圓柱繞流模型網(wǎng)格信息
Tab. 2 The details of three different meshes for flow past a single vertical cylinder

網(wǎng)格密度粗網(wǎng)格中網(wǎng)格細(xì)網(wǎng)格網(wǎng)格總單元數(shù)666 7541 174 0081 259 218圓柱節(jié)點(diǎn)數(shù)100100120展向結(jié)點(diǎn)數(shù)255050

表3 不同網(wǎng)格結(jié)果對(duì)比(α=0°, Re=100、1 500和3 900)Tab. 3 Comparisons between the present results and previous literature (α=0°, Re=100, 1 500 and 3 900)

表4 不同時(shí)間步長(zhǎng)結(jié)果對(duì)比(α=0°,Re=100、1 500和3 900)

Tab. 4Comparison of the results for a single vertical circular cylinder with different time steps (α=0°,Re=100, 1 500 and 3 900)

雷諾數(shù)Δt/sCdC′lStRe=1000.021.4310.2290.1600.011.391 (2.88%)0.217 (5.53%)0.163 (1.88%)0.0051.385 (0.43%)0.215 (0.93%)0.163 (0.00%)Re=1 5000.0031.3630.5860.2070.0021.354 (0.66%)0.555 (5.59%)0.200 (3.50%)0.0011.361 (0.52%)0.573 (3.24%)0.198 (1.01%)Re=3 9000.0021.1390.2140.1900.0011.132 (0.62%)0.209 (2.39%)0.189 (0.53%)0.000 51.134 (0.18%)0.210 (0.48%)0.189 (0.00%)

2 結(jié)果與討論

2.1 瞬時(shí)流場(chǎng)

采用Hunt等[27]推薦的Q準(zhǔn)則繪制了不同雷諾數(shù)及不同傾角情況下圓柱繞流的三維無量綱瞬時(shí)渦量等勢(shì)面圖，如圖3～圖5所示。在低雷諾數(shù)Re=100和高雷諾數(shù)Re=1 500和3 900時(shí)Q取值分別為0.1和100。在低雷諾數(shù)Re=100情況下，在-45° ≤α≤ 45°范圍內(nèi)，傾斜圓柱三維尾渦結(jié)構(gòu)均較為規(guī)則，尾流渦管與圓柱近乎平行，隨著圓柱傾角的變化而變化，幾乎無展向變形，表現(xiàn)出極強(qiáng)的二維特性。然而，在大傾角α=± 60°條件下，旋渦脫落完全被抑制，圓柱后方已不存在渦管，無旋渦脫落發(fā)生。當(dāng)雷諾數(shù)Re增大至1 500和3 900時(shí)，對(duì)于豎直圓柱(α=0°)而言，圓柱尾流并未出現(xiàn)明顯的展向摻混，但尾渦結(jié)構(gòu)在展向上存在一定的相位差異，圓柱尾流渦街較為規(guī)則；當(dāng)傾角增大至|α|=15°時(shí)，流場(chǎng)的展向摻混現(xiàn)象明顯，尾流渦街仍較為規(guī)則；而當(dāng)傾角繼續(xù)增大至|α| ≥ 30°時(shí)，圓柱尾流渦街的形態(tài)會(huì)發(fā)生明顯變化，沿圓柱軸向方向會(huì)產(chǎn)生較強(qiáng)的軸向流。在大傾角 |α|=45°和60°條件下，圓柱尾流渦脫會(huì)受到較大程度的抑制，從瞬時(shí)渦量圖中可以看出，圓柱后方尾流僅存在沿展向分布的渦管。在雷諾數(shù)Re=3 900情況下，流場(chǎng)的展向摻混更為強(qiáng)烈，三維特性更加明顯。

綜上所述，隨著圓柱傾角的增大，尾流旋渦的脫落與渦街的形成會(huì)逐漸受到抑制，細(xì)碎旋渦逐漸消失，尾流寬度也隨之減小。隨著雷諾數(shù)的增加，圓柱尾流渦管會(huì)發(fā)生明顯的變形，展向摻混隨著雷諾數(shù)的增大而愈發(fā)明顯，導(dǎo)致出現(xiàn)大量細(xì)碎旋渦，呈現(xiàn)出明顯的三維特性。在相同傾斜角度下，順流向傾斜圓柱與逆流向傾斜圓柱后方的尾渦結(jié)構(gòu)之間的差異非常小，不存在明顯的區(qū)別,說明沿順流向傾斜圓柱自底端流向頂端的軸向流動(dòng)與沿逆流向傾斜圓柱自頂端流向底端的軸向流動(dòng)對(duì)流場(chǎng)的影響較為一致。

圖3 Re=100時(shí)不同傾角瞬時(shí)無量綱渦量Q =0.1等勢(shì)面圖Fig. 3 Instantaneous iso-surface of vorticity contours at different inclined angles for Re=100 (Q =0.1)

圖4 Re=1 500時(shí)不同傾角瞬時(shí)無量綱渦量Q =100等勢(shì)面圖Fig. 4 Instantaneous iso-surface of vorticity contours at different inclined angles for Re=1 500 (Q =100)

圖5 Re=3 900時(shí)不同傾角瞬時(shí)無量綱渦量Q =100等勢(shì)面Fig. 5 Instantaneous iso-surface of vorticity contours at different inclined angles for Re=3 900 (Q =100)

圖6～圖8分別表示雷諾數(shù)Re=100、1 500和3 900情況下不同傾角條件下瞬時(shí)三維流線圖。當(dāng)雷諾數(shù)Re=100時(shí)，除|α|=60°外，圓柱尾流流場(chǎng)幾乎不存在展向摻混，流線十分規(guī)則且近乎平行，而當(dāng)傾角|α|=15°時(shí)，在圓柱后方的軸向流動(dòng)呈現(xiàn)較為明顯的螺旋形態(tài)。當(dāng)雷諾數(shù)增大至Re=1 500和3 900時(shí)，圓柱后方的三維流線形態(tài)較為近似且三維特性明顯：對(duì)于豎直圓柱(α=0°)而言，在圓柱后方存在回流區(qū)，隨著流動(dòng)逐漸往下游發(fā)展，旋渦消失，流線變得較為規(guī)則；當(dāng)傾角增大至|α|=15°，在圓柱后方存在回流區(qū)，此外，在圓柱后方觀測(cè)到較為明顯的呈螺旋運(yùn)動(dòng)形態(tài)的軸向流動(dòng)；當(dāng)傾角進(jìn)一步增大至 |α| ≥ 30°時(shí)，圓柱后方的旋渦逐漸消散，軸向流更加明顯，軸向流到達(dá)圓柱下端后繼續(xù)向下游發(fā)展。在相同傾斜角度下，沿順流向傾斜圓柱與逆流向傾斜圓柱的軸向流強(qiáng)度幾乎相同，沿順流向傾斜圓柱出現(xiàn)自底端流向頂端的軸向流動(dòng)，而沿逆流向傾斜圓柱出現(xiàn)自頂端流向底端的軸向流動(dòng)。

圖6 Re=100時(shí)不同傾角瞬時(shí)三維流線Fig. 6 Instantaneous three-dimensional streamline at different inclined angles for Re=100

圖7 Re=1 500時(shí)不同傾角瞬時(shí)三維流線Fig. 7 Instantaneous three-dimensional streamline at different inclined angles for Re=1 500

圖8 Re=3 900時(shí)不同傾角瞬時(shí)三維流線Fig. 8 Instantaneous three-dimensional streamline at different inclined angles for Re=3 900

2.2 時(shí)均流場(chǎng)

綜上所述，隨著雷諾數(shù)的增大，傾斜圓柱繞流的尾流流線愈加混亂，流場(chǎng)的三維特性愈加明顯。隨著傾角的變化，尾流回流區(qū)的存在范圍也會(huì)發(fā)生改變，當(dāng)雷諾數(shù)Re=100和1 500時(shí)，在小傾角|α| ≤ 15°條件下，傾斜圓柱尾流中存在一定范圍的回流區(qū)，在大傾角 |α| ≥ 30°條件下，回流區(qū)完全消失；而當(dāng)雷諾數(shù)Re=3 900時(shí)，在|α| ≤ 30°條件下存在回流區(qū)。在相同傾斜角度下，順流向傾斜圓柱與逆流向傾斜圓柱后方的回流區(qū)寬度一致，時(shí)均流線形態(tài)幾乎對(duì)稱，沿順流向傾斜圓柱出現(xiàn)自底端流向頂端的軸向流，而沿逆流向圓柱出現(xiàn)自頂端流向底端的軸向流。

圖9 Re=100時(shí)在y/D=0斷面處不同傾角的時(shí)均流線圖Fig. 9 Time-averaged streamline topologies at different inclined angles for Re=100

圖10 Re=1 500時(shí)在y/D=0斷面處不同傾角的時(shí)均流線圖Fig. 10 Time-averaged streamline topologies at different inclined angles for Re=1 500

圖11 Re=3 900時(shí)在y/D=0斷面處不同傾角的時(shí)均流線圖Fig. 11 Time-averaged streamline topologies at different inclined angles for Re=3 900

2.3 傾斜圓柱受力特性及旋渦脫落頻率

通過對(duì)不同雷諾數(shù)及傾斜角度下圓柱所受升力系數(shù)進(jìn)行快速傅里葉變換(FFT)，開展頻譜分析得到無量綱旋渦脫落頻率Strouhal數(shù)，如圖13～圖15所示。在雷諾數(shù)Re=100、1 500和3 900情況下，豎直圓柱對(duì)應(yīng)的St數(shù)分別為0.16、0.20和0.19，與Zhou等[23]和Vu等[25]所得結(jié)果較為接近，對(duì)應(yīng)的頻譜峰值均單一且明顯。當(dāng)雷諾數(shù)Re=100時(shí)，隨著圓柱傾角的不斷增大，St數(shù)呈單調(diào)遞減的趨勢(shì)，當(dāng)|α|=60°時(shí)旋渦脫落已被完全抑制，因此不存在St數(shù)。當(dāng)雷諾數(shù)Re=1 500時(shí)，傾角|α| ≥ 15°情況下，雖然圓柱尾渦結(jié)構(gòu)會(huì)受到軸向流的影響，但對(duì)應(yīng)的頻譜峰值仍比較單一，傾角|α|=45°和60°時(shí)對(duì)應(yīng)的St數(shù)分別為0.13與0.09。當(dāng)雷諾數(shù)Re=3 900時(shí)，傾角|α|=15°情況下，由于軸向流的存在導(dǎo)致軸向渦脫和尾流不規(guī)則渦脫，在頻譜主峰值的附近存在幅值較低的諧波。當(dāng)傾角進(jìn)一步增大到|α|=30°時(shí)，St數(shù)在0

圖13 Re=100時(shí)不同傾角的St數(shù)Fig. 13 Strouhal number at different inclined angles for Re=100

圖14 Re=1 500時(shí)不同傾角的St數(shù)Fig. 14 Strouhal number at different inclined angles for Re=1 500

圖15 Re=3 900時(shí)不同傾角的St數(shù)Fig. 15 Strouhal number at different inclined angles for Re=3 900

圖16 無量綱渦脫頻率St數(shù)隨傾角的變化曲線Fig. 16 Strouhal number as a function of the inclined angle

2.4 獨(dú)立性原則驗(yàn)證

圖17 獨(dú)立性原則驗(yàn)證Fig. 17 Validation of IP

3 結(jié) 語

對(duì)不同雷諾數(shù)(Re=100、1 500和3 900)和不同傾斜角度(-60°≤α≤ 60°)工況下沿順流向與逆流向傾斜的圓柱繞流流場(chǎng)開展了三維數(shù)值模擬研究，主要揭示了三維瞬時(shí)及時(shí)均尾流流場(chǎng)、流線拓?fù)洹⑸枇ο禂?shù)及旋渦脫落頻率隨雷諾數(shù)及傾斜角度變化的規(guī)律，驗(yàn)證了獨(dú)立性原則(IP)的適用性，得到的主要結(jié)論如下：

1) 隨著圓柱傾斜角度的增大，傾斜圓柱尾流產(chǎn)生較為明顯的軸向流，尾流旋渦脫落受到明顯抑制，細(xì)碎旋渦逐漸消失，尾流寬度也隨之減??；隨著雷諾數(shù)的增大，圓柱尾流渦管發(fā)生明顯的變形，展向摻混使得大量細(xì)碎旋渦產(chǎn)生，呈現(xiàn)出明顯的三維特性。

2) 隨著傾斜角度的增大，傾斜圓柱尾流回流區(qū)逐漸消失，當(dāng)雷諾數(shù)Re=100和1 500時(shí)，在傾角 |α| ≥ 30°條件下，回流區(qū)完全消失；當(dāng)雷諾數(shù)Re=3 900時(shí)，在|α| ≥ 45°條件下，回流區(qū)消失。