王春凌 王寅

摘要：內(nèi)波是分層水體中普遍存在的水動(dòng)力現(xiàn)象，破碎產(chǎn)生的湍流是打破水體層化的重要機(jī)制之一。本文采用數(shù)值模擬技術(shù)模擬研究了內(nèi)孤立波與典型大陸坡的相互作用過(guò)程。

關(guān)鍵詞：破碎；內(nèi)孤立波；能量耗散；地形

基金項(xiàng)目：國(guó)家自然科學(xué)基金青年基金項(xiàng)目

項(xiàng)目編號(hào)：52009087、52109090

基金項(xiàng)目：浙江省自然科學(xué)基金公益項(xiàng)目

項(xiàng)目編號(hào)：LGF19E090005

基金項(xiàng)目：中國(guó)博士后科學(xué)基金面上項(xiàng)目

項(xiàng)目編號(hào)：2022M721426

一、研究背景

內(nèi)波是密度層化水體中內(nèi)界面波動(dòng)現(xiàn)象，在遭遇地形或障礙物時(shí)，內(nèi)波會(huì)發(fā)生不同程度的破碎，激發(fā)的紊流有助于局部水體混合，改善局部水質(zhì)，無(wú)規(guī)律的水質(zhì)點(diǎn)運(yùn)動(dòng)作用于結(jié)構(gòu)，使結(jié)構(gòu)體發(fā)生“渦振”，破壞結(jié)構(gòu)穩(wěn)定性。內(nèi)波破碎和誘發(fā)機(jī)制復(fù)雜，Saffarinia和Kao等人用數(shù)值方法模擬了兩層流體中內(nèi)孤立波在平坦地形上傳播和破碎，發(fā)現(xiàn)當(dāng)波面指點(diǎn)流速超過(guò)內(nèi)波相速度時(shí)，下層水體被卷躍至上層，內(nèi)波發(fā)生破碎且此時(shí)波幅為上層水深的2.65倍，指出波的破碎是對(duì)流不穩(wěn)定導(dǎo)致的。Vlasenko和Hutter對(duì)內(nèi)孤立波與大陸架作用過(guò)程模擬發(fā)現(xiàn)，淺水效應(yīng)使內(nèi)波波面不斷陡峭直至翻轉(zhuǎn)、破碎，此時(shí)的破碎由對(duì)流不穩(wěn)定引起，隨后Afansyev和Peltier、Legg和Adcroft、Klymak和Pinkel等也得到類(lèi)似結(jié)論。Barad 和Fringer等人的研究則指出，剪切不穩(wěn)定同樣可導(dǎo)致內(nèi)波發(fā)生破碎，在內(nèi)重力波和加速剪切流的傳播過(guò)程中，剪切不穩(wěn)定的發(fā)生取決于剪切流的加速度，即當(dāng)加速度較小時(shí)，內(nèi)波會(huì)因?yàn)榧羟胁环€(wěn)定發(fā)生破碎，且內(nèi)波誘發(fā)紊流產(chǎn)生的水體混合效率與不穩(wěn)定剪切相關(guān)。然而，自然條件下內(nèi)波破碎通常是對(duì)流不穩(wěn)定和剪切不穩(wěn)定共同作用的結(jié)果，不穩(wěn)定特征與波幅和傳播方向相關(guān)。劉國(guó)濤等人對(duì)內(nèi)波在白噪聲干擾下引起的不穩(wěn)定現(xiàn)象進(jìn)行分析，結(jié)果顯示，靜力不穩(wěn)定最先出現(xiàn)，持續(xù)時(shí)間較長(zhǎng)，隨后剪切不穩(wěn)定和對(duì)流不穩(wěn)定相繼出現(xiàn)。

關(guān)于內(nèi)波破碎判據(jù)，Helfrich將離岸最遠(yuǎn)的初始湍流斑點(diǎn)出現(xiàn)位置定義為內(nèi)波破碎位置，并得出當(dāng)a/d≥0.4±0.1（a是斜坡起點(diǎn)處波幅，d是破碎位置的下層水深）時(shí)，內(nèi)波發(fā)生破碎。Vlasenko和Hutter認(rèn)為：當(dāng)無(wú)量綱振幅a*=am /（Hb－Hm）滿足a*=0.8/γ+0.4時(shí)，如果Hs＜Hb，內(nèi)波的背風(fēng)面在斜坡上破碎，反之會(huì)以頻散尾波的形式繼續(xù)傳播。Fringer和Street認(rèn)為內(nèi)波破碎發(fā)生的不穩(wěn)定主要與界面厚度相關(guān)，并給出kδ＜0.56時(shí)（k為水平波數(shù)，δ為界面厚度），內(nèi)波波峰和波谷出現(xiàn)Kelvin－Helmholtz（K－H）不穩(wěn)定；當(dāng)0.56≤kδ＜2.33時(shí)，初始不穩(wěn)定為剪切不穩(wěn)定；當(dāng)kδ≥2.33時(shí)，由于界面厚度較大無(wú)法形成K－H不穩(wěn)定，為對(duì)流不穩(wěn)定。Troy的研究結(jié)論認(rèn)為當(dāng)0＜kδ＜0.66時(shí)，擾動(dòng)長(zhǎng)度僅隨界面厚度變化，具有最大增長(zhǎng)率的擾動(dòng)得以充分發(fā)展，形成K－H滾動(dòng)，內(nèi)波破碎主要是剪切不穩(wěn)定引起；當(dāng)0.66≤kδ＜0.3時(shí)，擾動(dòng)增長(zhǎng)率減小，剪切不穩(wěn)定不能充分發(fā)展；當(dāng)kδ≥1.3時(shí)，剪切不穩(wěn)定不可能發(fā)生，內(nèi)波不發(fā)生破碎。

Barad和Fringer分別給出不穩(wěn)定性的判據(jù)條件：Richardson數(shù)Ri小于0.1；σT＞5， T是Ri數(shù)＜0.25的時(shí)間長(zhǎng)度，σ是擾動(dòng)在此時(shí)間內(nèi)的平均增長(zhǎng)率；Lw/L＞0.86， Lw是Ri＜1/4區(qū)域的長(zhǎng)度的一半，L是內(nèi)波的半寬度。近年來(lái)，有研究者提出內(nèi)波破碎與地形坡度和內(nèi)波坡度的比值相關(guān)，比值越小，破碎摻混越劇烈。相對(duì)于斜坡，孤島地形導(dǎo)致的內(nèi)波破碎通常用阻塞系數(shù)ζ（ζ=（a+H1）/（H1+H2+hs），a是波幅，H1、H2分別為上下層水深，Hs是障礙物高度描述將內(nèi)波與障礙物作用程度分為弱作用（ζ＜0.55），適度作用（0.55＜ζ＜0.7）和內(nèi)波破碎（ζ＞0.7），各研究者根據(jù)阻塞系數(shù)的取法不一致，結(jié)論有所不同。對(duì)于三層分層流體中的第二模態(tài)內(nèi)波，當(dāng)阻塞系數(shù)0.2＜B＜2（B=am/?h， am為波幅，?h為內(nèi)界面到地形頂面的高度）時(shí)，內(nèi)波將完全破碎分解。

針對(duì)內(nèi)波能量計(jì)算，多數(shù)研究者的論文中采用下式估算內(nèi)波勢(shì)能：

； （1）

此法簡(jiǎn)單易行，只需勾勒出波面形態(tài)即可進(jìn)行估算，多用于物理水槽試驗(yàn)中。然而，在內(nèi)波破碎劇烈情況下，用此方法進(jìn)行勢(shì)能估算誤差較大，一方面，混合區(qū)界面無(wú)法準(zhǔn)確捕捉；另一方面，作用過(guò)程形成的小尺度和微尺度內(nèi)波難以被精確刻畫(huà)。

綜上所述，內(nèi)波破碎仍有兩個(gè)問(wèn)題有待進(jìn)一步探討：內(nèi)波破碎過(guò)程特征尚未確定；傳統(tǒng)內(nèi)波勢(shì)能計(jì)算公式存在一定局限，真實(shí)能量的損耗、轉(zhuǎn)化需進(jìn)一步分析。

本文通過(guò)數(shù)值模擬方法，系統(tǒng)研究典型地形與內(nèi)波作用過(guò)程，分析水體剪切不穩(wěn)定和對(duì)流不穩(wěn)定發(fā)生機(jī)制，確定內(nèi)波與大陸架地形破碎的基本特征。同時(shí)，為探求更為精確的內(nèi)波能量損失值，對(duì)用式（1）的估算值進(jìn)行分析和修正，為理論模型用于實(shí)際工程提供可靠的參考及修正依據(jù)。

二、數(shù)學(xué)模型

（一）控制方程

基于連續(xù)介質(zhì)假設(shè)，描述不可壓縮粘性流體的連續(xù)方程和動(dòng)量方程為：

； （2）

；? （3）

式中，ui 和xi分別是i方向的速度和坐標(biāo)；ρ是壓強(qiáng)；是流體密度，本研究中清水密度取值1000kg/m3，鹽水密度取值為1030kg/m3；t是時(shí)間；μ是流體動(dòng)力黏度；gi是i方向加速度，豎直方向取值為－9.8m/s2，其他方向取值為0；si是動(dòng)量源項(xiàng)。

通常將內(nèi)波概化為兩層分層流體的多相流問(wèn)題，數(shù)值計(jì)算中通過(guò)求解體積分?jǐn)?shù)對(duì)當(dāng)?shù)匚锢砹窟M(jìn)行更新：

且? ?； （4）

式中，λ為更新的物性值，例如密度，Φm和λm分別為第m相的體積分?jǐn)?shù)和相屬性值。對(duì)于不可壓縮兩相流，只需計(jì)算某一相體積分?jǐn)?shù)，余下即為另外一相。

本文采用大渦模型求解控制方程，濾波后的控制方程和體積分?jǐn)?shù)輸運(yùn)方程為：

（5）

τij為亞格子應(yīng)力項(xiàng)，是可解尺度運(yùn)動(dòng)與亞格子尺度運(yùn)動(dòng)之間的動(dòng)量輸運(yùn)，需進(jìn)行模化處理，采用擬序渦結(jié)構(gòu)模型封閉方程。

（二）求解方法和邊界條件

采用結(jié)構(gòu)化網(wǎng)格對(duì)計(jì)算區(qū)域進(jìn)行離散，通過(guò)分步法（Fractional Method）分別求解對(duì)流項(xiàng)、擴(kuò)散項(xiàng)和壓力項(xiàng)，方程如下：

； （6）

； （7）

； （8）

式中，“*”和“**”分別代表對(duì)流步和擴(kuò)散步求解的中間物理量，“n”代表當(dāng)前時(shí)刻，“n+1”代表新時(shí)刻。在得到流場(chǎng)后，求解體積分?jǐn)?shù)（Volume Of Fraction，簡(jiǎn)稱(chēng)“VOF”）輸運(yùn)方程：

； （9）

同時(shí)，由式（3）更新本地密度。

邊界是否設(shè)置合理直接關(guān)系計(jì)算收斂和數(shù)值準(zhǔn)確性。本文中水槽自由液面采用垂向速度為零，其他方向速度梯度為零；底部采用固體壁面邊界；入流邊界采用滑移邊界；出流采用Sommerfeld邊界，避免入射波反射。

應(yīng)用兩層流體內(nèi)波關(guān)系式：（h1+a） U1 = （h2-a） U2 = c0η， 得到上下層水體的理論流速值。其中，η為波面高; h1和h2分別為上下層水深；U1和U2分別為上下層理論平均流速?；诓婧瘮?shù)和流速理論解，采用邊界造波方法進(jìn)行數(shù)值造波，與傳統(tǒng)方法相比，該方法保證了完整的波形，便于波能計(jì)算。

三、模型實(shí)驗(yàn)及結(jié)果

（一）模擬試驗(yàn)

試驗(yàn)在一數(shù)值水槽中進(jìn)行，水槽長(zhǎng)10m、寬0.5m、高0.4m，圖1為數(shù)值水槽示意圖。在入流邊界定義流速進(jìn)行造波，速度如何取值見(jiàn)上文。hs為大陸架高度，h1和h2分別為清水層和鹽水層厚度，清水密度為1000kg/m3，鹽水密度為1030kg/m3，水槽右側(cè)為自由開(kāi)邊界，防止反射波影響試驗(yàn)觀測(cè)區(qū)。

（二）模擬結(jié)果

模擬實(shí)驗(yàn)揭示了下凹型內(nèi)孤立波與大陸架地形相互作用過(guò)程，該過(guò)程可分為三種情況：當(dāng)?shù)匦胃叨扰c內(nèi)分界面距離較大時(shí)，內(nèi)孤立波不受影響繼續(xù)保持原有波形前行（圖2a）。隨著地形高度增加，在內(nèi)孤立波接近地形過(guò)程中，由于局部水流增大波面發(fā)生不同程度的扭曲變形甚至破碎產(chǎn)生紊流。因此，內(nèi)孤立波也被稱(chēng)為水體的天然攪拌器。局部的混合不僅利于上下層水體打破分層將富含氧的水體輸送到底層，其形成的第二模態(tài)內(nèi)波（圖2b）還可實(shí)現(xiàn)界面處的物質(zhì)遠(yuǎn)距離輸運(yùn)。除此之外，模擬研究發(fā)現(xiàn)與地形作用過(guò)程中存在內(nèi)孤立波轉(zhuǎn)變極性的情況（圖2c），這主要取決于地形上方的上下層水深比。根據(jù)KdV理論，h1/h2的值偏離1，其非線性越明顯，表現(xiàn)在波形上是波幅越大；反之，表現(xiàn)為彌散作用占主導(dǎo)，即波形越變得扁平。

當(dāng)?shù)匦胃叨鹊扔诨蚋哂诜謱咏缑鏁r(shí)，內(nèi)孤立波幾乎完全反射（圖2d），但反射系數(shù)與大陸架坡度有一定相關(guān)性，這一點(diǎn)在Chen等人（2007）的研究中已得到證實(shí)。圖3中為直角大陸架和非直角大陸架上內(nèi)孤立波誘發(fā)的流場(chǎng)，圖中內(nèi)孤立波向右傳播，引起的下層水體流向與波的傳播方向相反，很明顯在直角大陸架的迎波面有逆時(shí)針的漩渦形成，原因是直角處的水體受剪切；當(dāng)坡面較緩更貼合水流特性，坡底水體受剪切不易形成漩渦。漩渦的存在大大增加了對(duì)底部的剪切，造成地形不同程度的侵蝕，消耗的能量也相對(duì)較多。圖4為直角大陸架和非直角大陸架對(duì)內(nèi)波能量消耗的對(duì)比圖，從圖中可以看出，直角大陸架消耗的內(nèi)波能量高于帶坡度大陸架地形約1.5%，整體能量損失隨地形高度先增后減且符合線性規(guī)律，峰值位置出現(xiàn)在km≈0.3處，大小約為初始能量的20%～25%。

四、結(jié)語(yǔ)

本文采用數(shù)值試驗(yàn)的方法對(duì)內(nèi)波與典型大陸架地形的相互作用過(guò)程進(jìn)行了模擬，通過(guò)分析流場(chǎng)和能量計(jì)算得到如下結(jié)論：

內(nèi)波與大陸架作用的演變過(guò)程與地形高度相關(guān)，當(dāng)?shù)匦胃叨刃∮谙聦铀?/6時(shí)，內(nèi)波幾乎不受地形影響保持原有形態(tài)前行；當(dāng)高度超過(guò)此數(shù)值，內(nèi)波的演變還與地形上的上下層水深比相關(guān)，根據(jù)該比值偏離1的程度發(fā)生裂變和轉(zhuǎn)變極性。同時(shí)，伴隨內(nèi)波第二模態(tài)出現(xiàn)，地形上方水深比大于1時(shí)，第二模態(tài)波跟隨反射波；反之，模態(tài)二波跟隨透射波。

與地形作用過(guò)程損失的能量與地形高度相關(guān)，隨地形高度的增加出現(xiàn)先增后減的趨勢(shì)，且在無(wú)量綱參數(shù)km=0.3時(shí)，耗散能量達(dá)到峰值，約占初始能量的20%～25%。

與非直角大陸架相比，直角大陸架消耗的能量偏多，主要原因是與地形作用過(guò)程中在坡腳產(chǎn)生旋渦結(jié)構(gòu)，該結(jié)構(gòu)加大了水體剪切，從而加速能量損耗。

參考文獻(xiàn)：

[1]Aghsaee，P.and L.Boegman.Experimental Investigation of Sediment Resuspension Beneath Internal Solitary Waves of Depression.Journal of Geophysical Research：Oceans，2015，

120（05）： 3301-3314.

[2]Zhu，H.，L.Wang，and H.Tang. Large-eddy Simulation of the Generation and Propagation of Internal Solitary Waves.Science China Physics，Mechanics & Astronomy，2014，57（06）：1128-1136.

[3]Saffarinia，K.and T.W.Kao，A Numerical Study of the Breaking of an Internal Soliton and Its Interaction with a Slope.Dynamics of Atmospheres and Oceans，1996，

23：379-391.

[4]Vlasenko，V.and K.Hutter.Numerical Experiments on the Breaking of Solitary Internal Waves. Journal of Physical Oceanography，2002，32：1779-1793.

[5]Vlasenko，V.I.and K.Hutter.Generation of Second Mode Solitary Waves by the Interaction of a First Mode Soliton with a Sill.Nonlinear Processes in Geophysi-cs，2001，08（4/5）：223-239.

[6]Barad，M.F.and O.B.Fringer.Simulations of Shear Instabilities in Interfacial Gravity Waves. Journal of Fluid Mechanics，2010.644：61-95.

[7]劉國(guó)濤，尚曉東.海洋內(nèi)波能量耗散的數(shù)值研究[J].海洋石油，2019，39（01）：82-87.

[8]Helfrich，K.R.Internal Solitary Wave Breaking and Run-up on a Uniform Slope. Journal of Fluid Mechanics，1992，243：133-154.

[9]Arthur，R.S.and O.B.Fringer.The Dynamics of Breaking Internal Solitary Waves on Slopes. Journal of Fluid Mechanics，2014，761：360-398.

[10]Maderich，V，et al.Numerical Simulations of the Nonhydrostatic Transformation of Basin-Scale Internal Gravity Waves and Wave-Enhanced Meromixis in Lakes， in Nonlinear Internal Waves in Lakes，2012：193-276.

[11]Talipova，T，et al.Internal Solitary Wave Transformation Over a Bottom Step： Loss of Energy. Physics of Fluids，2013，25（03）：032110.

[12]Liu，A.K，et al.Generation and Evolution of Mode-two Internal Waves in the South China Sea.Continental Shelf Research，2013（59）：18-27.

[13]Chen，C.Y.Amplitude Decay and Energy Dissipation Due to the Interaction of Internal Solitary Waves with a Triangular Obstacle in a Two-layer Fluid System：the Blockage Parameter.Journal of Marine Science and Technology，2009，14（04）：499-512.

[14]Wang，C.，etal.Numerical Modeling of Energy Dissipation of Internal Solitary Waves Encountering Step Topography.Ocean Engineering，2022，259：111853.

[15]Wang，C.L.，L.L.Wang，andH.W.Tang，Depression Internal Solitary Wave Interaction with Submerged Obstacle.Journal of Environmental Biology，2016，37（05）：1217-1223.

[16]Chen，Y.C.，et al.Wave Propagation at the Interface of a Two-layer Fluid System in the Laboratary.Journal of Marine Science and Technology，2007，15（01）：8-16.

作者簡(jiǎn)介：王春凌（1985），男，安徽省安慶市人，中級(jí)講師，博士學(xué)歷，主要研究方向?yàn)樗W(xué)。