王 旭，尤云祥，黃文昊

(上海交通大學(xué)海洋工程國家重點實驗室，上海 200240)

內(nèi)孤立波是最大振幅發(fā)生在密度分層海水內(nèi)部的一種常見海洋動力學(xué)現(xiàn)象，不僅是海洋能量級串中的一個重要環(huán)節(jié)，也是對海洋工程結(jié)構(gòu)物安全性產(chǎn)生重要影響的一類海洋環(huán)境因素。在南海流花油田的早期延長測試期間，就曾發(fā)生過因內(nèi)孤立波產(chǎn)生的突發(fā)性強流而導(dǎo)致纜繩拉斷、船體碰撞，甚至拉斷和擠破漂浮軟管等事故［1］。在南海陸豐油田的早期延長測試期間，也曾發(fā)生過內(nèi)孤立波產(chǎn)生的突發(fā)性強流使半潛式鉆井船與錨定油輪在連接輸油管道時發(fā)生困難等問題［2］。

非線性和頻散效應(yīng)是密度分層海洋中影響內(nèi)孤立波傳播演化特性的兩個基本因素，非線性效應(yīng)趨向于使波形變陡，而頻散效應(yīng)則趨向于使波形變得平坦，內(nèi)孤立波的穩(wěn)定傳播是非線性與頻散效應(yīng)在一定尺度上達到平衡的結(jié)果，可以用KdV、eKdV和MCC等理論模型來描述［3］。設(shè)ε=a/h為非線性參數(shù)，μ=( h /λ)2為色散參數(shù)(其中a 和λ分別為內(nèi)孤立波振幅和特征寬度，h為垂向特征尺度)，KdV理論是建立在弱非線性、弱色散且兩者平衡(即ε=O(μ)?1)條件下的一類內(nèi)孤立波理論［4］。研究表明，在小振幅情況下KdV理論與實驗結(jié)果相符，但在大振幅情況KdV理論明顯與實驗結(jié)果不符［5-8］。

隨著內(nèi)孤立波振幅的增大，高階非線性項的影響會逐漸增大，而且非線性和色散性不再保持平衡［9］，這正是KdV理論不能適用于大振幅內(nèi)孤立波的主要原因之一。這個缺陷可以通過在KdV方程中加入立方非線性項等途徑解決，所得修正方程稱為eKdV方程［10］。但eKdV理論依然是建立在弱非線性和弱色散條件下的，為克服這兩類理論中需要弱非線性這個限制性條件的缺陷，從完全非線性歐拉方程出發(fā)，Choi和Camassa［11］建立了一個強非線性和弱色散的兩層Green-Naghdi模型，在定態(tài)條件下該理論可以簡化為一個常微分方程形式的非線性方程，稱為MCC理論。這是建立在強非線性和弱色散條件下的一類定態(tài)內(nèi)孤立波理論，研究表明在大振幅內(nèi)孤立波的情況，MCC理論與實驗結(jié)果相符［12］。

上述KdV、eKdV和MCC理論均是在弱色散條件下建立的，而對KdV和eKdV理論還需要弱非線性這個限制性條件，但關(guān)于弱非線性和弱色散這兩個限制性條件迄今仍是定性的。最近，黃文昊等［13］利用大型密度分層水槽，采用雙推板內(nèi)孤立波造波方法，研究了弱非線性和弱色散這兩個條件的定量表征問題，結(jié)果表明在給定密度分層及內(nèi)孤立波振幅條件下，KdV、eKdV和MCC理論的適用性與非線性和色散參數(shù)的某種組合有關(guān)。

各種海洋環(huán)境條件下深海浮式平臺載荷特性及其運動響應(yīng)的合理確定，是其設(shè)計與應(yīng)用中的一項關(guān)鍵性指標(biāo)。程友良等［14］和蔡樹群等［15-16］將Morison公式與KdV理論結(jié)合，而XIE等［17］則將Morison公式與MCC理論結(jié)合，研究了內(nèi)孤立波作用下直立貫底圓柱體的載荷特性問題。最近，尤云祥等［18-19］將Morison公式與eKdV理論結(jié)合，研究了內(nèi)孤立波作用下張力腿和半潛式平臺的載荷與動力響應(yīng)問題，而宋志軍等［20］將Morison公式與KdV理論結(jié)合，研究了內(nèi)孤立波作用下Spar平臺的載荷與動力響應(yīng)問題。在這些研究中，Morison公式中的慣性力和拖曳力系數(shù)是參照表面波的方法選取的，這種選取方法尚缺乏理論和實驗依據(jù)。

隨著計算流體力學(xué)(CFD)技術(shù)的發(fā)展，采用CFD方法來研究內(nèi)孤立波的生成傳播及其與海上結(jié)構(gòu)物的相互作用問題，已成為目前國際上的研究熱點之一。由于采用CFD方法可以直接獲得內(nèi)孤立波與海洋結(jié)構(gòu)物相互作用過程中的速度場及壓力場的變化特性，因此可以直接獲得內(nèi)孤立波作用下海上結(jié)構(gòu)物的載荷等水動力特性。在這種方法中，通常在入口邊界條件中以定態(tài)內(nèi)孤立波的解析解(包括KdV、eKdV和MCC等)來驅(qū)動完全非線性模型，進而模擬內(nèi)孤立波的生成傳播及其與海上結(jié)構(gòu)物的作用特性。

陳鈺等［21］、關(guān)輝等［22］和李水娟等［23］利用 KdV 理論解作為入口邊界條件，而高原雪等［24］利用 MCC 理論解作為入口邊界條件，采用CFD方法研究了內(nèi)孤立波的生成傳播問題。付東明等［25］利用KdV理論解作為入口邊界條件，而陳杰等［26］和關(guān)輝等［27］利用eKdV理論解作為入口邊界條件，采用CFD方法研究了內(nèi)孤立波作用下水下潛體的載荷特性問題。劉碧濤等［28］利用eKdV理論解作為入口邊界條件，采用CFD方法研究了內(nèi)孤立波作用下深海立管的載荷特性問題。

根據(jù)文獻［13］的實驗結(jié)果，在給定密度分層及內(nèi)孤立波振幅條件下，由于KdV、eKdV和MCC理論的適用性與非線性和色散參數(shù)的某種組合有關(guān)，因此如果選擇不合適的內(nèi)孤立波理論解作為CFD數(shù)值模擬的入口邊界條件，則有可能會導(dǎo)致振幅及波形都不可控的數(shù)值模擬結(jié)果。因此，在給定密度分層及內(nèi)孤立波振幅條件下，如何選擇合適的內(nèi)孤立波理論解作為CFD數(shù)值模擬的入口邊界條件是利用CFD方法來研究內(nèi)孤立波的生成傳播及其與海上結(jié)構(gòu)物相互作用問題時需要研究和解決的關(guān)鍵之一。

有鑒于此，以有限深兩層流體中內(nèi)孤立波誘導(dǎo)上下層深度平均水平速度為入口邊界條件，以理想流體完全非線性歐拉方程為控制方程，結(jié)合VOF動態(tài)界面跟蹤方法，對內(nèi)孤立波的生成特性進行完全非線性數(shù)值模擬。在此基礎(chǔ)上，研究在給定密度分層及內(nèi)孤立波振幅條件下如何選擇合適的內(nèi)孤立波理論解作為完全非線性數(shù)值模擬入口邊界條件等問題。

1 數(shù)值模擬方法

設(shè)兩層流體均為理想不可壓縮流體，當(dāng)流體處于靜平衡狀態(tài)時，上層流體深度與密度分別為h1和ρ1，而下層流體深度與密度分別為h2和ρ2，水深為h=h1+h2。建立直角坐標(biāo)系oxz，其中ox軸位于流體靜止時兩層流體的界面上，oz軸垂直向上為正。對水表面作剛蓋假定，假定水底為不可滲透平坦的剛性固壁。設(shè)ζ為兩層流體界面位移，ui、wi為流場中的速度矢量，其中i=1、2分別表示上下層流體。

在上下層流體中，流場控制方程如下:

其中，g為重力加速度，各物理量的下標(biāo)t、x和z分別表示關(guān)于時間和空間坐標(biāo)的導(dǎo)數(shù)，pi為壓力。

在兩層流體界面上，要求滿足如下連續(xù)性條件:

在水面和水底，要求滿足如下固壁條件:

流場計算的控制區(qū)域如圖1所示，包括內(nèi)孤立波生成傳播區(qū)和消波區(qū)兩個區(qū)域。采用速度入口方法生成內(nèi)孤立波，當(dāng)形成穩(wěn)定的內(nèi)孤立波后，對所生成內(nèi)孤立波的傳播特性進行監(jiān)測分析。

圖1 流場計算控制區(qū)域Fig.1 Governing domain for flow computation

速度入口方法為:在上層流體入口處，入口速度取為u1;在下層流體入口處，入口速度取為。內(nèi)孤立波在生成與傳播過程中，兩層流體的界面會發(fā)生變化，采用VOF(volume of fluid)方法追蹤兩層流體界面的變化［29］。在VOF方法中，用體積分數(shù)aq表示單元內(nèi)第q相流體所占體積與該單元的體積之比，其中q=1表示上層流體，q=2表示下層流體。當(dāng)aq=1時，表示單元內(nèi)全部為第q相流體，而當(dāng)aq=0時，表示單元不存在第q相流體，對含界面的單元則有0＜aq＜1。在每一個流體單元中各相的體積分數(shù)之和為1，即=1，而且對于任一相的體積分數(shù)aq，要求滿足如下連續(xù)方程:

采用海綿層方法對水槽末端的內(nèi)孤立波進行消波處理，該方法通過添加人工粘性項實現(xiàn)阻尼消波之目的。為此，對消波區(qū)，在動量方程的右端添加源項，將其修改為

其中，海綿層衰減系數(shù)μ依據(jù)文獻［30］中的方法確定。

采用KdV、eKdV和MCC理論計算內(nèi)孤立波誘導(dǎo)上下層流體中的層深度平均水平速度。為此，設(shè)a和λ分別為內(nèi)孤立波振幅與特征寬度，定義非線性參數(shù)ε=a/h和色散參數(shù)μ=(h/λ)2，那么在弱色散和弱非線性(即ε?1而且μ?1)條件下，界面位移ζ滿足如下的eKdV方程［3］

其中

當(dāng)內(nèi)孤立波為定態(tài)平面前進波時，eKdV方程(10)有如下形式的解［3］:

其中

其中，cekdv為內(nèi)孤立波eKdV理論解的相速度。

當(dāng)ε=O(μ)?1，即當(dāng)內(nèi)孤立波是弱非線性、弱色散且兩者平衡時，在eKdV方程中的立方非線性項可以忽略，即可取c3=0，這時eKdV理論解(13)即退化為如下的KdV理論解［3］

其中，ckdv為內(nèi)孤立波KdV理論解的相速度。

在上述理論中，均要求內(nèi)孤立波是弱非線性和弱色散的，對KdV理論還要求兩者是平衡的。因此，這兩類內(nèi)孤立波理論一般只適用于振幅較小的情況。為克服這個缺陷，Choi和Camassa［11］建立了一個強非線性和弱色散的內(nèi)孤立波理論模型，稱為MCC理論模型，其定態(tài)解可表示為

其中，cMCC為內(nèi)孤立波MCC理論解的相速度，

a-和a+(a-＜a+)分別為下面方程的兩個根

式中:

對于KdV理論解來說，只要h1/h2≠ ( h1/h2)c，則不存在極限振幅，即無論振幅有多大，該理論均能給出對應(yīng)的內(nèi)孤立波解。Camassa等［12］指出，KdV理論無極限振幅的這個缺陷，正是其不能適用于大振幅內(nèi)孤立波的原因之一。與KdV理論的一個重要不同之處是，其它兩類內(nèi)孤立波理論均存在極限振幅，其中eKdV和MCC理論解的極限振幅分別為［3］

2 結(jié)果與分析

文中將結(jié)合該文獻中的相關(guān)實驗結(jié)果進行數(shù)值模擬與分析。為此，數(shù)值水槽主尺度、上下層流體密度及其深度比均與該文一致。在數(shù)值模擬中，內(nèi)孤立波生成傳播區(qū)的長度為24 m，而消波區(qū)的長度為6 m。生成傳播區(qū)內(nèi)網(wǎng)格的縱向間距為0.03 m，垂向網(wǎng)格劃分以底部向上0.4 m位置為界，上部區(qū)域取網(wǎng)格垂向間距為0.005 m，下部區(qū)域首層垂向間距取為0.005 m，為保證網(wǎng)格平滑過渡，后續(xù)網(wǎng)格間距按照1.02的比例逐漸加寬;消波區(qū)內(nèi)垂向網(wǎng)格間距同生成傳播區(qū)的情況一致，而縱向首層網(wǎng)格間距取0.03 m，后續(xù)網(wǎng)格間距按1.04的比例逐漸加寬;整個計算域的網(wǎng)格單元共計60 800個。

利用商業(yè)軟件Fluent進行數(shù)值模擬，初始時水槽中流場物理量均取為零，時間步長為0.01 s。計算區(qū)域采用結(jié)構(gòu)化網(wǎng)格進行離散，在兩層流體界面附近，網(wǎng)格需進行加密處理?？刂品匠滩捎糜邢摅w積法(FVM)進行離散，對流項采用二階迎風(fēng)格式離散，時間步推進選用二階隱式時間匹配格式，兩相界面的構(gòu)造方法選用幾何重構(gòu)法，應(yīng)用非定常壓力的隱式算子分割法(PISO)進行壓力-速度的耦合計算。

文獻［13］的實驗結(jié)果表明，存在一個臨界色散參數(shù)μ0=0.1，以此為基準(zhǔn)對弱色散條件μ?1進行定量化，可將其分為弱色散(μ＜μ0)和強色散(μ≥μ0)兩種情況;在弱色散(μ＜μ0)條件下，對非線性條件也進行定量化，可將其分為弱非線性(ε≤μ)、中等非線性(μ＜ε≤)和強非線性(ε＞)三種情況;在系列實驗基礎(chǔ)上，獲得了KdV、eKdV和MCC理論的如下適用性條件:KdV理論適用于弱非線性和弱色散情況，eKdV理論適用于中等非線性和弱色散情況，而MCC理論適用于強非線性和弱色散或強色散情況。

在三種上下層流體深度比情況下，文獻［13］中系列實驗所得內(nèi)孤立波振幅如表1所示。在本文數(shù)值模擬中，初始內(nèi)孤立波的設(shè)計振幅ad也取為表中所列相應(yīng)的值。由式(21)可知，當(dāng)h1∶h2=3∶7和1∶9時，h=-0.182 6和 -0.211 8，而/h=-0.196 7和 -0.396 7;當(dāng) h1∶h2=2∶8時，/h =-0.234 1，而/h=-0.296 7。由此可知，當(dāng)h1∶h2=3∶7和1∶9時，各設(shè)計振幅均在eKdV和MCC理論解的極限振幅內(nèi);而當(dāng)h1∶h2=2∶8時，各設(shè)計振幅雖均在MCC理論解的極限振幅內(nèi)，但最后三個振幅已超出eKdV理論解的極限振幅。

表1 數(shù)值模擬中的無因次設(shè)計振幅ad/hTab.1 Non-dimensional design amplitudes for numerical simulations ad/h

在數(shù)值模擬中，對給定的設(shè)計振幅ad，利用KdV、eKdV和MCC理論解分別計算其誘導(dǎo)的上下層流體中的層深度平均水平速度(見式(6))，將其作為入口邊界條件來驅(qū)動完全非線性數(shù)值模型(見式(1)-(5))。在上述初始內(nèi)孤立波設(shè)計振幅ad下，這里重點研究如何選擇合適的定態(tài)內(nèi)孤立波理論解來驅(qū)動完全非線性數(shù)值模型的問題。

為此，設(shè)εm和μm分別為數(shù)值模擬所得波的非線性和色散參數(shù)，而εd和μd分別為在設(shè)計振幅ad下利用內(nèi)孤立波理論解計算得到的非線性和色散參數(shù)。在圖2中，給出了當(dāng)h1/h2=3∶7時，(εm，μm)和(εd，μd)的結(jié)果比較。其中:“▲”表示(εm，μm)，“2”表示(εd，μd)。在給定設(shè)計振幅ad下，點▲和2之間的橫向距離表示μm和μd之間的絕對誤差，而點▲和2之間的垂向距離表示εm和εd之間的絕對誤差。由圖可知，在給定設(shè)計振幅ad下，利用eKdV理論解作為入口邊界條件，所得μm和μd以及εm和εd之間的絕對誤差最小。因此，在表1所示各設(shè)計振幅下，利用eKdV理論解作為入口邊界條件，數(shù)值模擬所得內(nèi)孤立波振幅及設(shè)計振幅誤差最小，而且兩者波形最為吻合。

圖2 當(dāng) h1∶h2=3∶7 時 (εm，μm)和 (εd，μd)結(jié)果比較Fig.2 Comparisons between(εm，μm)and(εd，μd)when h1∶h2=3∶7

由圖2(a)可知，在各設(shè)計振幅ad下，KdV理論解和數(shù)值模擬所得波的非線性和色散參數(shù)均屬于弱色散和中等非線性區(qū)域，而文獻［13］的實驗表明，KdV理論適用于弱色散和弱非線性的情況，因此采用KdV理論解作為入口邊界條件本身就是不合適的。由圖2(b)可知，在各設(shè)計振幅ad下，eKdV理論解和數(shù)值模擬所得波的非線性和色散參數(shù)也均屬于弱色散和中等非線性區(qū)域，而文獻［13］的實驗表明，eKdV理論適用于弱色散和中等非線性的情況，兩者正好是一致的，即采用eKdV理論解作為入口邊界條件是合適的。由圖2(c)可知，MCC理論解和數(shù)值模擬所得波的非線性和色散參數(shù)均屬于弱色散和中等非線性區(qū)域，而文獻［13］的實驗表明，MCC理論適用于弱色散和強非線性或強色散的情況，因此采用MCC理論解作為入口邊界條件本身就是不合適的。

在圖3中，給出了當(dāng)h1/h2=3∶7時，數(shù)值模擬結(jié)果與初始內(nèi)孤立波及相關(guān)實驗結(jié)果的比較。由圖中第一列可知，當(dāng)利用KdV理論解作為入口邊界條件時，數(shù)值模擬所得波的振幅及其波形與初始輸入波的振幅及其波形差異都很大，而且隨著設(shè)計振幅的增大，數(shù)值模擬結(jié)果還會出現(xiàn)顯著的尾波現(xiàn)象。由圖中第二列可知，當(dāng)利用eKdV理論解作為入口邊界條件時，數(shù)值模擬所得波的振幅及其波形與初始輸入波和文獻［13］中實驗所得波振幅及其波形之間的差異都很小，而且數(shù)值模擬結(jié)果沒有明顯的尾波現(xiàn)象。由圖中第三列可知，與利用KdV理論解作為入口邊界條件的結(jié)果相比，當(dāng)利用MCC理論解作為入口邊界條件時，數(shù)值模擬所得波的振幅及其波形與初始輸入波的振幅及其波形之間的差異明顯減小，但這種差異要比利用eKdV理論解作為入口邊界條件的結(jié)果大，因此采用eKdV理論解作為入口邊界條件更為合適。

由圖2可知，當(dāng)h1/h2=3∶7時，在各設(shè)計振幅ad下，利用KdV理論解作為入口邊界條件，數(shù)值模擬所得波均屬于弱色散和中等非線性的。因此，根據(jù)文獻［13］中的實驗結(jié)果，如果該波為內(nèi)孤立波，則應(yīng)該屬于eKdV型的。為此，在圖4中，給出了當(dāng)ad/h=0.034 0、0.105 5和0.158 0時，利用KdV理論解作為入口邊界條件所得數(shù)值模擬結(jié)果與等振幅eKdV理論解的比較。由圖4(a)可知，在較小無因次振幅的情況，數(shù)值模擬結(jié)果的波形與等振幅eKdV理論解波形一致，這表明數(shù)值模擬所得波為eKdV型內(nèi)孤立波。由圖4(b)和4(c)可知，在較大無因次振幅的情況，數(shù)值模擬所得波的波形與等振幅eKdV理論解的波形差異很大，這表明在這種情況下利用KdV理論解作為入口邊界條件是不合適的。

圖3 當(dāng)h1∶h2=3∶7時數(shù)值模擬結(jié)果與初始內(nèi)孤立波及相關(guān)實驗結(jié)果比較Fig.3 Comparisons of numerical results with theoretical and experimental ones when h1∶h2=3∶7

圖4 當(dāng)h1∶h2=3∶7時利用KdV理論解作為入口邊界條件所得數(shù)值模擬結(jié)果與eKdV理論解比較Fig.4 Comparisons of numerical results with eKdV theoretical solutions based on the KdV input conditions when h1∶h2=3∶7

圖5給出了當(dāng)h1/h2=2∶8時，(εm，μm)和(εd，μd)的結(jié)果比較。由圖可知，對表1中前3個設(shè)計振幅ad，利用KdV理論解作為入口邊界條件，所得μm和μd以及εm和εd之間的誤差最小，而且利用KdV理論解和數(shù)值模擬所得波的非線性和色散參數(shù)均屬于弱色散和弱非線性區(qū)，與文獻［13］實驗結(jié)果一致，即這時采用KdV理論解作為入口邊界條件最為合適。對表1中第4到9個設(shè)計振幅ad，利用eKdV理論解作為入口邊界條件，所得μm和μd以及εm和εd之間的誤差最小，而且利用eKdV理論解和數(shù)值模擬所得波的非線性和色散參數(shù)均屬于弱色散和中等非線性區(qū)，與文獻［13］實驗結(jié)果一致，即這時采用eKdV理論解作為入口邊界條件最為合適。對表1中最后4個設(shè)計振幅ad，利用MCC理論解作為入口邊界條件，所得μm和μd以及εm和εd之間的誤差最小，而且利用MCC理論解和數(shù)值模擬所得波的非線性和色散參數(shù)均屬于弱色散和強非線性區(qū)域，與文獻［13］實驗結(jié)果一致，即這時采用MCC理論解作為入口邊界條件最為合適。

圖5 當(dāng) h1∶h2=2∶8 時 (εm，μm)和 (εd，μd)結(jié)果比較Fig.5 Comparisons between(εm，μm)and(εd，μd)when h1∶h2=2∶8

在圖6中，給出了當(dāng)h1/h2=2∶8時，數(shù)值模擬結(jié)果與初始內(nèi)孤立波及相關(guān)實驗結(jié)果的比較。由圖可知，當(dāng)ad/h=0.062 1時，與利用eKdV和MCC理論解作為入口邊界條件的結(jié)果相比，利用KdV理論解作為入口邊界條件所得數(shù)值模擬波的振幅及其波形與初始輸入波最為吻合，而且與文獻［13］中的實驗一致。當(dāng)ad/h=0.130 4時，與利用KdV和MCC理論解作為入口邊界條件的結(jié)果相比，利用eKdV理論解作為入口邊界條件所得數(shù)值模擬波的振幅及其波形與初始輸入波最為吻合，而且與文獻［13］中的實驗一致。當(dāng)ad/h=0.243 3時，與利用KdV和eKdV理論解作為入口邊界條件的結(jié)果相比，利用MCC理論解作為入口邊界條件所得數(shù)值模擬波的振幅及其波形與初始輸入波最為吻合，而且與文獻［13］中的實驗一致。

圖6 當(dāng)h1∶h2=2∶8時數(shù)值模擬結(jié)果與初始內(nèi)孤立波及相關(guān)實驗結(jié)果比較Fig.6 Comparisons of numerical results with theoretical and experimental ones when h1∶h2=2∶8

由圖5可知，當(dāng)h1/h2=2∶8時，對表1中最后10個設(shè)計振幅ad，利用KdV理論解作為入口邊界條件，數(shù)值模擬所得波均屬于弱色散和中等非線性的。因此，根據(jù)文獻［13］中實驗結(jié)果，如果該波為內(nèi)孤立波，則應(yīng)該屬于eKdV型的。為此，在圖7中，給出了當(dāng)ad/h=0.130 4和0.243 3時，利用KdV理論解作為入口邊界條件所得數(shù)值模擬結(jié)果與等振幅eKdV理論解的比較。由圖可知，數(shù)值模擬結(jié)果的波形與等振幅eKdV理論解波形一致，這表明數(shù)值模擬所得波為eKdV型內(nèi)孤立波。

圖7 當(dāng)h1∶h2=2∶8時，利用KdV理論解作為入口邊界條件所得數(shù)值模擬結(jié)果與eKdV理論解比較Fig.7 Comparisons of numerical results with eKdV theoretical solutions based on the KdV input conditions when h1∶h2=2∶8

圖8 當(dāng) h1∶h2=1∶9 時 (εm，μm)和 (εd，μd)結(jié)果比較Fig.8 Comparisons between(εm，μm)and(εd，μd)when h1∶h2=1∶9

在圖8中，給出了當(dāng)h1∶h2=1∶9時，(εm，μm)和 (εd，μd)的結(jié)果比較。由圖可知，在給定設(shè)計振幅 ad下，利用MCC理論解作為入口邊界條件，所得μm和μd以及εm和εd之間的誤差最小。因此，在表1所示各設(shè)計振幅下，利用MCC理論解作為入口邊界條件，數(shù)值模擬所得內(nèi)孤立波振幅與設(shè)計振幅誤差最小，而且兩者波形最為吻合。由圖8(a)和8(b)可知，在各設(shè)計振幅ad下，KdV和eKdV理論解和數(shù)值模擬所得波的非線性和色散參數(shù)均屬于強弱色散區(qū)域。文獻［13］的實驗表明，KdV理論適用于弱色散和弱非線性的情況，而eKdV理論適用于弱色散和中等非線性的情況，因此采用KdV和eKdV理論解作為入口邊界條件本身就是不合適的。由圖8(c)可知，在各設(shè)計振幅ad下，利用MCC理論解和數(shù)值模擬所得波的非線性和色散參數(shù)均屬于強弱色散區(qū)域，而文獻［13］的實驗表明，MCC理論適用于強弱色散的情況，兩者正好是一致的，即采用MCC理論解作為入口邊界條件是合適的。

在圖9中，給出了當(dāng)h1∶h2=1∶9時，數(shù)值模擬結(jié)果與初始內(nèi)孤立波及相關(guān)實驗結(jié)果的比較。由圖中第一列和第二列可知，當(dāng)利用KdV和eKdV理論解作為入口邊界條件時，隨著設(shè)計振幅的增大，數(shù)值模擬所得波的振幅及其波形與初始輸入波的振幅及其波形之間的差異也隨著增大。由圖中第三列可知，當(dāng)利用MCC理論解作為入口邊界條件時，數(shù)值模擬所得波的振幅及其波形與初始輸入波和文獻［13］中實驗所得波振幅及其波形之間的差異都很小。

圖9 當(dāng)h1∶h2=1∶9時數(shù)值模擬結(jié)果與初始內(nèi)孤立波及相關(guān)實驗結(jié)果比較Fig.9 Comparisons of numerical results with theoretical and experimental ones when h1∶h2=1∶9

圖10 當(dāng)h1∶h2=1∶9時利用KdV和eKdV理論解作為入口邊界條件所得數(shù)值模擬結(jié)果與MCC理論解比較Fig.10 Comparisons of numerical results with MCC theoretical solutions based on the KdV and eKdV input conditions when h1∶h2=1∶9

由圖8可知，當(dāng)h1/h2=1∶9時，在各設(shè)計振幅ad下，利用KdV和eKdV理論解作為入口邊界條件，數(shù)值模擬所得波均屬于強色散的。因此，根據(jù)文獻［13］中的實驗結(jié)果，如果該波為內(nèi)孤立波，則應(yīng)該屬于MCC型的。為此，在圖10中，給出了當(dāng)ad/h=0.041 8，0.137 3和0.151 3時，利用KdV和eKdV理論解作為入口邊界條件所得數(shù)值模擬結(jié)果與等振幅MCC理論解的比較。由圖可知，數(shù)值模擬結(jié)果的波形與等振幅MCC理論解波形一致，這表明數(shù)值模擬所得波為MCC型內(nèi)孤立波。

3 結(jié)語

這里為有限深水內(nèi)孤立波非線性數(shù)值模擬中如何選擇合適的理論解作為入口邊界條件提供了理論依據(jù)。所得入口邊界條件選擇方法僅限于有限深水的情況，針對無限深水條件下應(yīng)如何選擇合適的內(nèi)孤立波理論解作為入口邊界條件等問題，尚需做進一步研究。

［1］ 陳景輝.南海流花11-1深海油田開發(fā)工程［J］.中國海洋平臺，1996，11(1):43-45.(CHEN Jing-hui.Liuhua 11-1 deep sea oil-field development project in South China Sea［J］.China Offshore Platform，1996，11(1):43-45.(in Chinese))

［2］ Bole J B，Ebbesmeyer C C，Romea R D.Soliton currents in the South China Sea:measurements and theoretical modeling［C］//The 16thOffshore Technology Conference.1994:367-376.

［3］ Helfrich K R，Melville W K.Long nonlinear internal waves［J］.Ann.Rev.Fluid Mech.，2006，38:395-425.

［4］ Ostrovsky L A，Stepanyants Yu A.Do internal solitons exist in the ocean［J］.Rev.Geophys.，1989，27(3):293-310.

［5］ Michallet H，Barthelemy E.Experimental study of interfacial solitary waves［J］.J.Fluid Mech.，1998，366:159-177.

［6］ Grue J，Jensen A，Rusas P E，et al.Properties of large-amplitude internal waves［J］.J.Fluid Mech.，1999，380:257-278.

［7］ Sveen J K，Guo Y，Davies P A，et al.On the breaking of internal solitary waves at a ridge［J］.J.Fluid Mech.，2002，469:161-188.

［8］ Walker S A，Martin A J，Easson W J，et al.Comparison of laboratory and theoretical internal solitary wave kinematics［J］.Journal of Waterway，Port，Coastal and Ocean Engineering，2003，129(5):210-218.

［9］ Brandt P，Rubino A，Alpers W，et al.Internal waves in a strait of messina studied by a numerical model and synthetic aperture radar form the ERS1/2 satellites［J］.J.Phys.Oceanography.，1997，27(5):648-663.

［10］Djordjevic V D.The fission and disintegration of internal solitary waves moving over two-dimensional topography［J］.J.Phys.Oceanography，1978，8(6):1016-1024.

［11］Choi W，Camassa R.Fully nonlinear internal waves in a two-fluid system［J］.J.Fluid Mech.，1999，396:1-36.

［12］Camassa R，Choi W，Michallet H，et al.On the realm of validity of strongly nonlinear asymptotic approximations for internal waves［J］.J.Fluid Mech.，2006，549:1-23.

［13］黃文昊，尤云祥，王旭，等.有限深兩層流體中內(nèi)孤立波造波實驗及其理論模型［J］.物理學(xué)報，2013，62(8):1-14.(HUANG Wen-hao，YOU Yun-xiang，WANG Xu，et al.Wave-making experiments and theoretical models for internal solitary waves in a two-layer fluid of finite depth［J］.Acta Phys.Sin.，2013，62(8):1-14.(in Chinese))

［14］CHENG Y L，LI J C，LIU Y F.The induced flow field by internal solitary wave and its action on cylindrical piples in the stratified ocean［C］//The 4th International Conference on Fluid Mechanics.2004:296-299.

［15］CAI S Q，LONG X M，GAN Z J.A method to estimate the forces exerted by internal solitons on cylinder piles［J］.Ocean Engineering，2003，30(5):673-689.

［16］CAI S Q，LONG X M，WANG S G.Force and torques exerted by internal solitions in shear flows on cylinder piles［J］.Applied Ocean Research，2008，30(1):72-77.

［17］XIE J S，JIANG Y J，YANG L G.Strongly nonlinear internal solution load on a small vertical circular cylinder in two-layer fluids［J］.Applied Mathematical Modelling，2010，34(8):2089-2101.

［18］尤云祥，李巍，何景異，等.海洋內(nèi)孤立波中張力腿平臺的水動力特性［J］.上海交通大學(xué)學(xué)報，2010，44:12-17.(YOU Yun-xiang，LI Wei，HE Jing-yi，et al.Hydrodynamic characteristics of tension leg platforms in ocean internal solitary waves［J］.Journal of Shanghai Jiao Tong University，2010，44:12-17.(in Chinese))

［19］尤云祥，李巍，胡天群，等.內(nèi)孤立波中半潛平臺動力響應(yīng)特性［J］.海洋工程，2012，30(2):1-7.(YOU Yun-xiang，LI Wei，HU Tian-qun，et al.Dynamic responses of a semi-submersible platform in internal solitary waves［J］.The Ocean Engineering，2012，30(2):1-7.(in Chinese))

［20］宋志軍，勾瑩，滕斌，等.內(nèi)孤立波作用下Spar平臺的運動響應(yīng)［J］.海洋學(xué)報，2010，32:12-19.(SONG Zhi-jun，GOU Ying，TENG Bin，et al.The motion responses of a Spar platform under internal solitary wave［J］.Acta Oceanologica Sinica，2010，32(2):12-19.(in Chinese))

［21］陳鈺，朱良生.基于 FLUENT的海洋內(nèi)孤立波數(shù)值水槽模擬［J］.海洋技術(shù)，2009，28(4):72-75.(CHEN Yu，ZHU Liang-sheng.Numerical wave tank simulation of oceanic internal solitary waves based on FLUENT［J］.Ocean Technology，2009，28(4):72-75.(in Chinese))

［22］關(guān)暉，蘇曉冰，田俊杰.分層流體中內(nèi)孤立波數(shù)值造波方法及其應(yīng)用［J］.力學(xué)季刊，2011，32(2):218-223.(GUAN Hui，SU Xiao-bing，TIAN Jun-jie.Numerical method of solitary wave generation in a stratified flow and its applications［J］.Chinese Quarterly Of Mechanics，2011，32(2):218-223.(in Chinese))

［23］李水娟，余真真，羅龍洪，等.三維數(shù)值波浪水槽內(nèi)孤立波特性分析［J］.水電能源科學(xué)，2012，30(7):96-99.(LI Shui-juan，YU Zhen-zhen，LUO Long-hong，et al.Analysis of internal solitary waves characteristics of 3d numerical wave flume［J］.Water Resources and Power，2012，30(7):96-99.(in Chinese))

［24］高原雪，尤云祥，王旭，等.基于 MCC理論的內(nèi)孤立波數(shù)值模擬［J］.海洋工程，2012，30(4):29-36.(GAO Yuan-xue，YOU Yun-xiang，WANG Xu，et al.Numerical simulation for the internal solitary wave based on MCC theory［J］.The Ocean Engineering，2012，30(4):29-36.(in Chinese))

［25］付東明，尤云祥，李巍.兩層流體中內(nèi)孤立波與潛體相互作用數(shù)值模擬［J］.海洋工程，2009，27(3):38-44.(FU Dong-ming，YOU Yun-xiang，LI Wei.Numerical simulation of internal solitary waves with submerged body in a two-layer fluid［J］.The Ocean Engineering，2009，27(3):38-44.(in Chinese))

［26］陳杰，尤云祥，劉曉東，等.內(nèi)孤立波與有航速潛體相互作用數(shù)值模擬［J］.水動力學(xué)研究與進展，2010，25(3):344-350.(CHEN Jie，YOU Yun-xiang，LIU Xiao-dong，et al.Numerical simulation of interaction of internal solitary waves with a moving submarine［J］.Chinese Journal of Hydrodynamics，2010，25(3):344-350.(in Chinese))

［27］關(guān)暉，魏崗，杜輝.內(nèi)孤立波與潛艇相互作用的水動力學(xué)特性［J］.解放軍理工大學(xué)學(xué)報:自然科學(xué)版，2012，13(5):577-582.(GUAN Hui，WEI Gang，DU Hui，et al.Hydrodynamic properties of interactions of three-dimensional internal solitary waves with submarine［J］.Journal of PLA University of Science and Technology，Natural Science Edition，2012，13(5):577-582.(in Chinese))

［28］劉碧濤，李巍，尤云祥，等.內(nèi)孤立波與深海立管相互作用數(shù)值模擬［J］.海洋工程，2011，29(4):1-7.(LIU Bi-tao，LI Wei，YOU Yun-xiang，et al.Numerical simulation of interaction of internal solitary waves with deep-sea risers［J］.The Ocean Engineering，2011，29(4):1-7.(in Chinese))

［29］Hirt C W，Nichols B D.Volume of fluid(VOF)method for the dynamics of free boundaries［J］.J.Computational Physics，1981，39(1):201-225.

［30］韓鵬.基于VOF方法的不規(guī)則波阻尼消波研究［D］.大連:大連理工大學(xué)，2008.(HAN Peng.The study of damping absorber for irregular waves based on VOF method［D］.Dalian:Dalian University of Technology，2008.(in Chinese))