基于勢(shì)流理論的內(nèi)孤立波與立管作用數(shù)值研究1)

胡英杰 鄒 麗,?,2) 孫 哲 金國(guó)慶 馬鑫宇

* (大連理工大學(xué)船舶工程學(xué)院,遼寧大連 116024)

? (工業(yè)裝備結(jié)構(gòu)分析國(guó)家重點(diǎn)實(shí)驗(yàn)室,遼寧大連 116024)

引言

海洋內(nèi)波是一種發(fā)生在海面以下的波浪,具有波幅大、周期長(zhǎng)的特征,內(nèi)波在海洋中經(jīng)常以大波幅的內(nèi)孤立波及波列的形式傳播,在世界范圍的海域內(nèi)廣泛存在[1-3].在傳播過(guò)程中,內(nèi)孤立波的流場(chǎng)會(huì)在海洋表面產(chǎn)生明顯的特征,因而使用合成孔徑雷達(dá)(SAR)對(duì)海面特征的捕捉成為對(duì)其觀測(cè)的重要依據(jù)[4-6].同時(shí),使用聲學(xué)多普勒流速儀(ADCP)等儀器的海洋實(shí)測(cè)也在多個(gè)區(qū)域觀測(cè)到了內(nèi)孤立波的存在[7-8].通過(guò)遙感和實(shí)地觀測(cè)表明我國(guó)南海也是內(nèi)孤立波頻發(fā)的海域[9].在海洋資源開(kāi)發(fā)過(guò)程中,海洋石油和深海礦產(chǎn)輸送立管長(zhǎng)期在海上作業(yè),海洋立管在內(nèi)孤立波流場(chǎng)的作用下會(huì)產(chǎn)生劇烈而大幅度的動(dòng)力響應(yīng),對(duì)管道的構(gòu)形、結(jié)構(gòu)強(qiáng)度、作業(yè)效率與安全都產(chǎn)生嚴(yán)重影響[10].海洋內(nèi)孤立波不僅對(duì)水下潛器和結(jié)構(gòu)物影響巨大,還會(huì)對(duì)海洋平臺(tái)等水面結(jié)構(gòu)物產(chǎn)生影響.內(nèi)孤立波在遠(yuǎn)距離傳播過(guò)程中,會(huì)在海洋表面引導(dǎo)出明顯的海水流動(dòng),因而會(huì)造成海洋平臺(tái)的長(zhǎng)距離漂移,對(duì)錨泊結(jié)構(gòu)是一個(gè)重大考驗(yàn),影響正常生產(chǎn)作業(yè),造成嚴(yán)重的經(jīng)濟(jì)損失[11].在安德曼海作業(yè)的一個(gè)石油鉆井機(jī)曾被海洋內(nèi)孤立波的流場(chǎng)推移了30.48 m,嚴(yán)重影響作業(yè)安全[12].

近年來(lái),內(nèi)孤立波對(duì)海洋立管的作用引起了更多學(xué)者的關(guān)注.文獻(xiàn)[13-15]針對(duì)內(nèi)孤立波對(duì)海洋立管的作用問(wèn)題開(kāi)展了系統(tǒng)深入的研究,通過(guò)KdV 和mKdV 方程計(jì)算了內(nèi)孤立波的流場(chǎng)特征,采用動(dòng)力有限元模型得到了頂張力立管的動(dòng)力響應(yīng),并討論了內(nèi)波參數(shù)對(duì)管道響應(yīng)的影響.文獻(xiàn)[16-18]采用試驗(yàn)和數(shù)值模擬方法研究了內(nèi)孤立波與海流共同作用下管道與海洋平臺(tái)樁柱的動(dòng)力響應(yīng)問(wèn)題,獲得了不同入射角度及不同波浪參數(shù)條件下的管道變形特征.Wang等[19]采用非線(xiàn)性有限元模型計(jì)算了弱約束深海礦產(chǎn)輸送管道在內(nèi)孤立波作用下的變形與運(yùn)動(dòng)特性,討論了內(nèi)孤立波不同入射角度及水面平臺(tái)運(yùn)動(dòng)對(duì)管道的影響.Wang等[20]采用試驗(yàn)方法研究了內(nèi)孤立波對(duì)小尺度樁柱的作用力分布特征,討論了不同水深條件下樁柱的受力情況,分析了水動(dòng)力系數(shù)選取與KC數(shù)的關(guān)系.Fan等[21]基于歐拉-伯努利梁理論,考慮了海流、表面波和內(nèi)孤立波的共同作用,發(fā)現(xiàn)內(nèi)孤立波的存在明顯增大了立管的變形包絡(luò),并且分層流體的密度差也會(huì)對(duì)立管變形產(chǎn)生影響.李朝瑋等[22]采用準(zhǔn)靜態(tài)分析法研究了非均勻海流、表面波和內(nèi)孤立波作用下隔水管的變形,發(fā)現(xiàn)內(nèi)孤立波對(duì)立管的橫向位移和應(yīng)力都有不可忽視的影響且表面波和內(nèi)波相位角分別在12°和0°時(shí)對(duì)立管的影響最大.Zhang等[23]采用數(shù)值模擬和試驗(yàn)方法研究了平臺(tái)運(yùn)動(dòng)和內(nèi)孤立波共同作用下頂張力立管的渦激振動(dòng)特性,發(fā)現(xiàn)內(nèi)孤立波的波幅越大、平臺(tái)升沉運(yùn)動(dòng)的頻率越高,立管的渦激振動(dòng)越明顯.Wang等[24]在分層流水槽中開(kāi)展模型試驗(yàn)分析了斜向入射的內(nèi)孤立波的波流結(jié)構(gòu)特征及其對(duì)水下細(xì)長(zhǎng)結(jié)構(gòu)的載荷作用特性.文獻(xiàn)[25-28]采用莫里森經(jīng)驗(yàn)公式、模態(tài)分解和回歸分析等方法對(duì)內(nèi)孤立波對(duì)柱體及管道的動(dòng)力載荷進(jìn)行了深入分析.

內(nèi)孤立波流場(chǎng)速度和加速度信息的準(zhǔn)確獲取是立管載荷確定的關(guān)鍵,進(jìn)而影響整個(gè)計(jì)算模型的準(zhǔn)確性.對(duì)于大尺度的海洋內(nèi)孤立波,無(wú)法開(kāi)展實(shí)驗(yàn)獲得流場(chǎng)信息,同時(shí)求解NS 方程的CFD 方法計(jì)算量較大.目前對(duì)大尺度內(nèi)孤立波與立管作用的研究中,內(nèi)孤立波流場(chǎng)信息多是基于解析模型近似給定的.本文工作基于勢(shì)流理論采用多域邊界單元法建立數(shù)值模型求解了內(nèi)孤立波的流場(chǎng)特征,較為高效又能實(shí)時(shí)準(zhǔn)確地計(jì)算流場(chǎng)信息.依據(jù)流場(chǎng)信息采用莫里森方程給定立管載荷,采用動(dòng)力學(xué)有限元方法計(jì)算了頂張力立管和弱約束立管的動(dòng)力響應(yīng),以期為海洋立管設(shè)計(jì)在考慮內(nèi)孤立波作用時(shí)提供參考.

1 計(jì)算模型

1.1 分層流體運(yùn)動(dòng)勢(shì)流計(jì)算模型

內(nèi)孤立波與立管作用如圖1 所示,兩層流體的深度分別為h1,h2;密度分別為ρ1,ρ2.下凹型內(nèi)孤立波位于密度界面分層處.

圖1 內(nèi)孤立波與立管作用示意圖Fig.1 Schematic diagram of the interaction between internal solitary wave and the riser

兩層流體都假定為無(wú)旋,無(wú)黏性,不可壓縮的,因此兩層流體中存在流動(dòng)勢(shì)函數(shù)φk滿(mǎn)足拉普拉斯方程

在分層流體界面處的運(yùn)動(dòng)學(xué)邊界條件為上下層流體沿界面法線(xiàn)方向的速度連續(xù),因此有

其中η2為內(nèi)孤立波波面,n為計(jì)算區(qū)域Ω1,Ω2在界面處的外法線(xiàn)方向矢量.

在分層流體界面處的動(dòng)力學(xué)邊界條件為上下層流體區(qū)域在界面處的壓力連續(xù),因此有

其中g(shù) 為重力加速度.海底滿(mǎn)足不可穿透條件

由于內(nèi)孤立波傳播時(shí)在海洋表面產(chǎn)生的波動(dòng)較小,為了提高計(jì)算效率,本文對(duì)上表面采用剛蓋假定,因此上表面的邊界條件為

通過(guò)格林第二定理,可將上下兩層流體區(qū)域Ω1,Ω2的控制方程拉普拉斯方程轉(zhuǎn)化為兩個(gè)邊界上的曲線(xiàn)積分方程

其中r為源點(diǎn),q為場(chǎng)點(diǎn),G(r,q)=-1/(2π)ln|r-q|為格林函數(shù),c(r)=1-α(r)/(2π),α(r)為場(chǎng)點(diǎn)r處的單元夾角,Γk為上下層區(qū)域邊界.

由于內(nèi)孤立波位于密度分層的流場(chǎng)中,所以計(jì)算過(guò)程需要涉及兩個(gè)流體區(qū)域,進(jìn)而要將單獨(dú)區(qū)域的邊界單元法擴(kuò)展到多個(gè)區(qū)域[29-30].對(duì)于計(jì)算區(qū)域Ω1,Ω2的邊界分別采用n個(gè)線(xiàn)單元進(jìn)行離散,從而積分方程 (6) 可以寫(xiě)成如下的形式

根據(jù)兩個(gè)區(qū)域的關(guān)系,可以將邊界點(diǎn)分為公共邊界點(diǎn)和自身邊界點(diǎn),從而將積分方程進(jìn)一步處理為

其中,上標(biāo)c和s分別指公共邊界點(diǎn)和自身邊界點(diǎn).

波面位置和速度勢(shì)的更新格式分別根據(jù)運(yùn)動(dòng)學(xué)和動(dòng)力學(xué)邊界條件給定

內(nèi)孤立波初始波形的給定依據(jù)強(qiáng)非線(xiàn)性的MCC 模型,實(shí)現(xiàn)對(duì)大波幅內(nèi)孤立波的準(zhǔn)確模擬.

1.2 管道動(dòng)力響應(yīng)模型

在內(nèi)孤立波由遠(yuǎn)方向立管傳播過(guò)程中,立管處的流場(chǎng)速度及加速度是實(shí)時(shí)變化的,因此不能用簡(jiǎn)單的均勻流或剪切流靜力假定求解,需要采用動(dòng)力學(xué)進(jìn)行求解.

海洋立管的動(dòng)力學(xué)控制方程

其中,M,C,K分別為質(zhì)量陣、剛度陣和阻尼陣,x(t),分別為管道的變形位移、速度和加速度.Q(t)為載荷向量,根據(jù)內(nèi)孤立波流場(chǎng)的速度、加速度以由莫里森方程給定.

在計(jì)算中需要考慮管道軸向張力的影響,因此管道的剛度陣由兩個(gè)部分組成,分別是彈性剛度陣KE和幾何剛度陣KG,即

其中L為管道單元的長(zhǎng)度.

阻尼陣采用瑞利阻尼,質(zhì)量陣與剛度陣的線(xiàn)性組合形式

其中 α 和 β 為比例系數(shù),根據(jù)模態(tài)頻率給定,分別為

ω1和ω2分別為管道的一階和二階頻率,ζ為管道的阻尼比.

采用Newmark-β法求解動(dòng)力學(xué)控制方程 (11),獲得管道的動(dòng)力響應(yīng).

2 內(nèi)孤立波載荷

在通過(guò)對(duì)邊界積分方程的求解獲得區(qū)域邊界上的速度勢(shì)及其法向?qū)?shù)值后,流場(chǎng)內(nèi)部速度勢(shì)可以用方程 (15) 表達(dá).對(duì)積分方程 (15) 進(jìn)行微分,即可將速度勢(shì)的積分方程轉(zhuǎn)化為關(guān)于流場(chǎng)速度的積分方程 (16),對(duì)該方程進(jìn)行離散求解即可獲得內(nèi)孤立波流場(chǎng)內(nèi)部的速度分布特征,即

在進(jìn)行計(jì)算求解時(shí),當(dāng)內(nèi)部的待求點(diǎn)與邊界節(jié)點(diǎn)接近時(shí),格林函數(shù)會(huì)表現(xiàn)出近奇異特征,導(dǎo)致計(jì)算精度降低甚至使計(jì)算失效.本文的研究中采用了單元子分法[31-32]來(lái)解決這一問(wèn)題,得到了良好的計(jì)算結(jié)果.為了驗(yàn)證計(jì)算模型的準(zhǔn)確性,與文獻(xiàn)中的內(nèi)孤立波水平速度垂向分布試驗(yàn)結(jié)果[33]進(jìn)行了對(duì)比.試驗(yàn)在分層流體水槽中開(kāi)展,上層流體深度h1=15 cm,下層深度h2=62 cm,密度分別為ρ1=0.999 kg/m3,ρ2=1022 kg/m3,內(nèi)孤立波無(wú)量綱波幅為a/h1=0.22.計(jì)算結(jié)果與分層流體內(nèi)孤立波傳播試驗(yàn)流場(chǎng)的測(cè)量結(jié)果吻合良好如圖2 所示.計(jì)算得到的內(nèi)孤立波流場(chǎng)速度分布如圖3 所示,計(jì)算結(jié)果表明內(nèi)孤立波流場(chǎng)速度具有明顯的剪切特征,上層流體的質(zhì)點(diǎn)速度與內(nèi)孤立波的傳播方向相同,而下層流體則與之相反,且上層流體質(zhì)點(diǎn)速度的絕對(duì)值要明顯大于下層流體.圖4 為內(nèi)孤立波流場(chǎng)的水平速度在上層流體區(qū)域z/h1=0.2 處沿水平方向的分布情況.可以發(fā)現(xiàn)水平速度的最大值位于波谷位置(x/h1=25),且隨著與波谷距離的增大而迅速衰減.圖5為z/h1=0.2 處流體質(zhì)點(diǎn)加速度在波長(zhǎng)方向的分布情況,加速度的最大值位于波谷兩側(cè),波谷位置處的加速度為零.

圖2 水平速度垂向分布與文獻(xiàn)[33]結(jié)果對(duì)比Fig.2 Comparison between the calculation result about the vertical distribution of horizontal velocity with the result of Ref.[33]

圖3 內(nèi)孤立波流場(chǎng)速度分布Fig.3 Velocity distribution of internal solitary wave flow field

圖4 內(nèi)孤立波水平速度沿水平方向的分布Fig.4 Horizontal velocity distribution of internal solitary wave

圖5 內(nèi)孤立波水平加速度分布Fig.5 Horizontal acceleration distribution of internal solitary wave

由于深海礦產(chǎn)輸送管道具有大長(zhǎng)細(xì)比特征,因此采用經(jīng)典的莫里森公式作為橋梁將非線(xiàn)性?xún)?nèi)孤立波數(shù)值計(jì)算模型和管道響應(yīng)的計(jì)算模型連接起來(lái),從而實(shí)現(xiàn)內(nèi)孤立波對(duì)管道作用的準(zhǔn)確求解.

由于內(nèi)孤立波對(duì)管道的作用過(guò)程中管道變形較大,所以要考慮管道的位移以及管道傾角的影響,因此內(nèi)孤立波對(duì)管道作用力的莫里森公式寫(xiě)成如下形式

其中,Cm為慣性力系數(shù),本文計(jì)算取為2.0;Cd為拖曳力系數(shù),取為1.2;Un為垂直于軸向的內(nèi)孤立波流場(chǎng)速度,為對(duì)應(yīng)的加速度;ρ為海水密度,x˙ 和x¨ 分別為管道運(yùn)動(dòng)的速度和加速度.根據(jù)前面內(nèi)孤立波流場(chǎng)速度和加速度的計(jì)算結(jié)果及莫里森方程可知,當(dāng)管道的直徑相對(duì)較小時(shí),對(duì)管道響應(yīng)起主要作用的載荷為內(nèi)孤立波對(duì)管道的拖曳力,慣性力的影響相對(duì)較小.

3 管道動(dòng)力響應(yīng)

3.1 頂張力立管動(dòng)力響應(yīng)

在海洋石油開(kāi)采等領(lǐng)域,頂張力立管得到了廣泛的應(yīng)用,本文計(jì)算了內(nèi)孤立波對(duì)頂張力立管的作用特性.為了驗(yàn)證本文計(jì)算模型的準(zhǔn)確性,將計(jì)算結(jié)果與文獻(xiàn)[34]結(jié)果進(jìn)行了對(duì)比.上下層流體深度分分別為60 m 和412 m,密度分別為ρ1=1025 kg/m3,ρ2=1028 kg/m3,管道內(nèi)外徑分別為0.26 m 和0.21 m,內(nèi)孤立波波幅為75 m.對(duì)比結(jié)果如圖6 所示,可以發(fā)現(xiàn)兩者吻合較好.由于本文內(nèi)孤立波流場(chǎng)的給定方法與文獻(xiàn)有一定區(qū)別,所以動(dòng)力響應(yīng)結(jié)果稍有偏差,管道順流向最大位移處的偏差約為3.34%.

圖6 管道流向最大位移與文獻(xiàn)[34]計(jì)算結(jié)果對(duì)比Fig.6 Comparison of the maximum displacement of the pipe in flow direction with the calculation results in the Ref.[34]

本文計(jì)算的頂張力立管的幾何參數(shù)和物理參數(shù)見(jiàn)表1.為了討論頂張力大小對(duì)管道動(dòng)力響應(yīng)的影響,頂張力T分別取為400 kN,500 kN,600 kN,700 kN,800 kN.

表1 頂張力立管參數(shù)Table 1 Parameters of the top tension riser

計(jì)算中采用的分層流體參數(shù)見(jiàn)表2.為了研究?jī)?nèi)孤立波參數(shù)對(duì)管道動(dòng)力響應(yīng)的影響,內(nèi)孤立波的波幅分別給定為40 m,50 m,60 m,70 m,80 m.

表2 流體分層參數(shù)Table 2 Parameters of stratified fluids

不同內(nèi)孤立波波幅時(shí)管道的變形特征如圖7 所示,可以發(fā)現(xiàn)立管的動(dòng)力響應(yīng)對(duì)波高敏感,隨著波高的增加,管道的流向位移明顯增大.同時(shí),由于上層流體速度明顯大于下層,且在所研究問(wèn)題中拖曳力遠(yuǎn)大于慣性力,因此管道順流向的最大位移發(fā)生在上層區(qū)域.

圖7 不同內(nèi)孤立波入射波幅時(shí)頂張力立管的流向變形特征Fig.7 Characteristics deformation of top tension riser in flow direction with different incident amplitudes of internal solitary waves

當(dāng)內(nèi)孤立波波幅為60 m 時(shí)不同頂張力條件下的立管變形特征如圖8 所示,可以發(fā)現(xiàn)頂張力對(duì)管道的流向變形具有重要影響.在相同內(nèi)孤立波參數(shù)的條件下,施加的頂張力越大,管道的變形越小.圖9展示了不同張力時(shí)管道流向最大變形的變化趨勢(shì).這是由于隨著頂張力T的增加,作用在管道單元上的有效張力Te也逐漸增大,通過(guò)方程 (13) 可知,管道的幾何剛度陣KG值將隨之增大,進(jìn)而提高管道的整體剛度.因此,在相同參數(shù)的內(nèi)孤立波作用下,頂張力越大,管道的順流向變形越小.

圖8 不同頂張力時(shí)管道的變形特征Fig.8 Deformation characteristics of the riser under different top tensions

圖9 管道最大流向位移隨頂張力大小的變化關(guān)系Fig.9 The relationship between the maximum displacement of the riser in flow direction and the value of the top tension

圖10 展示了z=-40 m 位置處管道的應(yīng)力在不同內(nèi)孤立波波幅條件下的變化特征.管道應(yīng)力隨著波幅的增大而顯著增大,同時(shí)波幅較小的內(nèi)孤立波作用時(shí)管道的最大應(yīng)力時(shí)刻出現(xiàn)滯后.根據(jù)分層流體K d V 理論,內(nèi)孤立波的波速可以表示為,其中c0和c1是與分層流體厚度和密度有關(guān)的量,因此當(dāng)分層參數(shù)確定之后,內(nèi)孤立波的波速與波幅密切相關(guān),波幅越大的內(nèi)孤立波具有更大的波速.因此,當(dāng)波幅較小時(shí),內(nèi)孤立波傳播到管道位置處的時(shí)間顯著增加,管道的最大應(yīng)力出現(xiàn)的時(shí)刻明顯滯后.

圖10 不同波幅內(nèi)孤立波作用下管道的應(yīng)力變化Fig.10 Stress changes of the riser under the action of solitary waves with different amplitudes

3.2 弱約束立管動(dòng)力響應(yīng)

本文討論的弱約束立管是指上端與水面設(shè)備剛固連接,下端為自由端的海洋立管,一般在深海采礦等領(lǐng)域得到應(yīng)用.計(jì)算所采用的管道幾何和物理參數(shù)如表3 所示.

表3 弱約束立管參數(shù)Table 3 Parameters of the weakly constrained riser

為了討論內(nèi)孤立波參數(shù)對(duì)弱約束立管的作用規(guī)律,計(jì)算了不同波幅內(nèi)孤立波作用下的立管動(dòng)力響應(yīng)如圖11 所示.從圖中可以發(fā)現(xiàn)內(nèi)孤立波波幅的變化對(duì)立管變形作用明顯,立管的流向位移隨初始內(nèi)孤立波波幅的增大而顯著增大.

圖11 不同波幅內(nèi)孤立波作用下弱約束立管的變形特征Fig.11 Characteristics deformation of weakly constrained riser in flow direction with different incident amplitudes of internal solitary waves

對(duì)于深海礦產(chǎn)輸送管道,輸送礦物濃度的變化會(huì)導(dǎo)致內(nèi)部流體密度的變化,因此本文討論了內(nèi)部流體不同密度對(duì)管道響應(yīng)特征的影響.不同內(nèi)部流體濃度條件下管道變形特征如圖12 所示.從圖中發(fā)現(xiàn),管道內(nèi)部流體密度對(duì)弱約束管道的流向變形影響較小.

圖12 內(nèi)部流體密度對(duì)弱約束管道位移的影響Fig.12 Influence of inner fluid densities on displacement of weakly constrained riser

4 結(jié)論

本文基于勢(shì)流理論采用多域邊界單元法建立了內(nèi)孤立波流場(chǎng)的計(jì)算模型,內(nèi)孤立波波面采用全非線(xiàn)性動(dòng)力學(xué)和運(yùn)動(dòng)學(xué)邊界條件,準(zhǔn)確獲得了內(nèi)孤立波流場(chǎng)的速度分布特征.上下層的流體質(zhì)點(diǎn)速度方向相反,對(duì)于下凹型內(nèi)孤立波,上層流體速度與波的傳播方向相同,下層速度與之相反,且上層流體質(zhì)點(diǎn)速度明顯大于下層.

根據(jù)計(jì)算獲得的內(nèi)孤立波流場(chǎng)信息,采用莫里森方程確定了內(nèi)孤立波對(duì)海洋立管的載荷輸入,計(jì)算了頂張力立管和弱約束立管在內(nèi)孤立波作用下的動(dòng)力響應(yīng)特征.管道的順流向位移和管道應(yīng)力隨內(nèi)孤立波波幅的增大而顯著,同時(shí)頂張力通過(guò)改變幾何剛度陣而對(duì)管道的動(dòng)力響應(yīng)產(chǎn)生明顯影響.由于上層流體速度明顯大于下層,且在所研究問(wèn)題中拖曳力遠(yuǎn)大于慣性力,因此管道順流向的最大位移發(fā)生在上層區(qū)域.同時(shí)計(jì)算了不同波幅和不同內(nèi)部流體密度條件下弱約束管道的動(dòng)力響應(yīng),發(fā)現(xiàn)內(nèi)孤立波波幅對(duì)弱約束管道的動(dòng)力響應(yīng)影響顯著,而管道內(nèi)部流體密度的影響則相對(duì)較小.

本文中內(nèi)孤立波對(duì)立管的作用載荷根據(jù)流場(chǎng)信息采用莫里森方程給定,沒(méi)有考慮管道的渦激振動(dòng)響應(yīng)特征.為了更加全面分析管道在內(nèi)孤立波作用下的動(dòng)力響應(yīng)規(guī)律,需要在今后的研究中進(jìn)一步改進(jìn),為實(shí)際工程中的海洋立管設(shè)計(jì)和安全性評(píng)估提供更加充分的依據(jù).