一種適用于陡波條件下研究半潛平臺周圍波浪非線性爬升效應的計算方法

單鐵兵，楊建民，李 欣，肖龍飛

（上海交通大學 海洋工程國家重點實驗室，上海200240）

1 引 言

隨著深海油氣田的不斷開發，淺海開發技術已難以滿足要求，傳統的重力式平臺和導管架平臺因其自重和造價隨水深增加而大幅度增加，已經不適應深海油氣發展的要求，取而代之的是浮式生產系統。其中半潛平臺具有性能優良，運動響應較小，工作水深范圍廣，抗風浪能力強，甲板面積大及裝載能力強等優點，在海洋石油勘探和開采中得到了廣泛的應用[1]。隨著石油勘探開采朝深海以及超深海發展，復雜的深海環境條件以及近年來超強臺風等極具破壞力的極端海況的頻繁出現，對平臺的安全性提出了嚴峻挑戰。

諸多的安全性評估指標中，氣隙預報（海洋平臺下層甲板底部至波面間垂直距離的預報）一直是半潛平臺設計過程中一個關鍵問題[2]。在波浪與平臺的相互作用過程中，有一種現象特別引起關注，即波浪遇平臺立柱后將沿其表面往上攀升，有時其高度明顯大于其它區域的波浪，甚至能爬升至平臺下甲板，以致嚴重影響平臺的氣隙性能，該非線性現象稱為波浪爬升。波浪的爬升過程從理論上可解釋為當入射波浪沖擊浮體時，水波的動能轉化為勢能同時波高發生放大的過程。

導管架平臺由于相應管架直徑小，對入射波擾動小，故波浪作用機理簡單，采用線性理論即可對其周圍的波面分布進行預報。對于大型的浮式結構物如半潛平臺而言，因體型大，且下浮體與甲板之間通過若干立柱支撐，當波陡較大時，波浪非線性效應較強，柱體周圍波浪爬升現象尤為明顯，波高易發生放大，并表現出較強的非線性特征，甚至超過平臺下甲板發生越浪，此過程中產生的砰擊載荷將對平臺的安全性構成威脅。陡波條件下的波浪爬升問題近年來已成為水動力學研究的熱點問題之一，是平臺氣隙預報的一個重要方面。

很早以來就有學者對波浪的繞射和爬升效應進行了試探性的研究。McCamy等[3]1954年給出了單根圓柱周圍繞射波浪的經典線性解。Niedzwecki等[4]1992年針對該圓柱上的波浪爬升問題開展了模型試驗，試驗結果表明，波陡（H/λ，其中H為波高，λ為波長）較大時，采用線性散射理論大大低估了波浪的爬升高度。Huston等[5]對一座四立柱的張力腿平臺（TLP）在單色規則波下進行了一系列的試驗測量，旨在研究波陡、立柱間的間隙等對立柱上的波浪爬坡和甲板下部波浪放大的影響。Lee和Newman等[6]1994年采用基于二階邊界元的方法對平臺周圍的波浪效應進行了模擬。盡管二階理論能夠提高波浪爬升預報的精確度，但仍不能滿足平臺詳細設計的要求，且無法捕捉到波浪沖擊立柱時的高階非線性現象[7-8]。

計算機的發展使平臺立柱周圍的波浪散射、海浪沿立柱爬升、甚至砰擊平臺甲板后波浪發生變形和翻滾等強非線性現象的直接數值模擬成為可能。與二階理論不同的是，CFD理論基于N-S方程直接求解時域內的流場。近年來有學者針對一些自由液面大變形的強非線性問題提出了較為精確的追蹤重構自由面的數值方法。其中常見的有Volume-of-Fluid（VOF）、Smooth-Particle-Hydrodynamics（SPH）、Level-Set（LS）以及CIP方法等。Donald[9]采用兩種方法求解重力式平臺周圍的波面分布，第一種為二階散射程序（WIMIT），基于勢流和小波陡假設,采用截斷的攝動展開并結合邊界元理論求解波動場，考慮二階非線性特征。第二種為基于VOF理論的CFD方法，該方法直接求解不可壓縮N-S方程，采用的是全非線性自由面邊界條件。通過對比發現，線性波條件下，二階散射理論可較為準確地模擬平臺上的波浪爬升，但當波陡較大時，二階攝動理論限制了其預報的有效性。相比而言，VOF理論即使在入射波浪波陡較大時也能夠捕捉到波物相互作用時的高階非線性效應。Roberto等[10]基于雷諾平均的N-S方程，使用僅求解水相的LS方法模擬自由液面，研究粘性效應對立柱周圍波面分布的影響。

本文以連續性方程和動量方程為控制方程，按照與深水試驗池相同的搖板造波方式，采用VOF方法捕捉自由液面，建立了三維數值波浪水槽，為了防止水池末端的波浪發生反射而對工作區域造成影響，在動量方程中添加源項進行消波。采用六面體結構化網格對數值模型進行離散，保證了網格的質量，提高了計算的精度。文中還對多種網格劃分形式進行了詳細研究，同時與試驗測量值進行對比，確定了最優的網格劃分方案，保證了數值方法的可行性。實際的模型試驗中，因浪高儀測量設備數量有限，模型試驗前確定平臺哪些位置安裝浪高儀對研究平臺周圍的波浪分布更具價值就顯得十分重要。此外，入射波陡較大時，波浪的非線性特征較強，與平臺之間的水動力作用加劇，將導致平臺立柱周圍的波浪爬升效應變得更為復雜，甚至產生砰擊、越浪等對平臺結構和安全性構成一定威脅的強非線性現象，相比平緩的波浪，研究波陡較大的波浪更有實際意義。介此，在該波浪水槽中，對陡波條件下半潛平臺周圍的波浪分布進行數值模擬。計算時，在平臺周圍設定多個波浪監測點，用以詳細預報不同位置處波浪的繞射和爬升效應，確定波浪非線性較強的區域，為下一步模型試驗中確定浪高儀的安裝位置提供重要依據。此外，文中還從物理現象上模擬并分析了波浪沿立柱的繞射和爬升過程，為今后進一步研究波物之間的非線性相互作用，如波浪參數對立柱周圍波浪爬升的影響、極限波浪下半潛平臺的強非線性砰擊和越浪、半潛平臺的氣隙響應等提供參考。

2 數值波浪

2.1 控制方程

半潛平臺波浪繞射和爬升的數值模擬是在一個類似于海洋深水試驗池的數值水池中進行的。數值模擬包括搖板造波的動態區域，過渡區域，工作區域以及消波區域。

設定流體不可壓，則流場的連續性方程為：

相應的動量方程：

式中：x1，x2，x3分別為x，y，z方向的坐標；u1，u2，u3分別為x，y，z方向的速度；f1=0，f2=0，f3=-g；S1，S2，S3分別為x，y，z方向上的附加源項；μ為動力學粘性系數；ρ為流體的密度。

自由液面采用VOF（Volume of Fluid）方法進行求解。該方法的基本原理是通過網格單元中的各流體和網格體積比函數aq來求解自由面的位置。aq=1為該網格單元全部由指定相流體所占據；aq=0說明該網格單元內無指定相流體；0＜aq＜1說明該網格單元內有部分指定相流體。氣液兩相流的自由面波動的求解方程可寫為：

其中：q=1表示空氣相，q=2表示水相。

2.2 搖板造波

數值造波和消波的技術是數值波浪水池模擬成功的關鍵，一直是國內外學者關注的問題。造波的方式有很多種，Boo等[11]通過在入口邊界給定速度或速度勢，然后基于時間步進方式造出二階和三階stokes波。Brorsen和Larsen等[12]提出了一種適合積分方程法的域內源造波方法，即在計算區域內沿水深方向設定一造波源，在源項兩邊同時產生方向相反的兩列波。周勤俊和王本龍等[13]從忽略粘性的歐拉方程出發，以N-S方程為控制方程，將源項添加到動量方程進行造波和消波。董志等[14]采用仿物理的推板造波方式，通過設定造波板邊界按簡諧運動實現規則波的模擬，該方法造波原理簡單，物理造波機理論能較容易地應用到數值造波中。

實際的深水試驗池大多采用搖板式造波方法生成波浪，即通過機械驅動使搖板繞固定軸擺動，強迫水質點運動，達到模擬相應波動的目的。波浪波高通過搖板的擺動幅度控制，其擺動頻率則決定波浪的周期大小。實驗室對規則波的模擬中，搖板驅動信號通過線性波浪理論計算得出。根據波浪理論，若搖板按固定角度做簡諧運動，在離開搖板一段距離后波浪將呈現二維規則波特性，即波浪周期與搖板簡諧運動相同，波高與擺幅有關。本文基于CFD理論，采用與上海交通大學深水試驗池相同的搖板造波方式生成波浪。

對于規則波浪，造波板自由液面處的運動方程可寫為：

其中：S0為造波板自由液面處的運動幅值，k為波數。

通過勢流理論和造波板的運動方程，可推導出造波板運動幅值和波高H之間的轉換關系[15-16]，如（5）式。已知規則波浪的相關參數，為了確定造波板的運動規律，需要知道搖板造波方式下的水力傳遞函數T（ω,d），該值同樣可由勢流理論得出，如（6）式所示。

其中：d為水深，k為波數，ω為波浪圓頻率。而k可由下面的色散關系求得。

造波板的運動方程可寫為：

根據線性造波機理論，在水深為d的水池中，波面方程為：

本文利用動網格模塊中的鋪層（Laying）功能中的網格更新技術，指定造波板為動邊界，并繞某一固定點轉動來實現數值造波。按照深水試驗池中所采用的線性造波機原理，推導出搖板擺動的角度信號，如（10）式所示。其中R為初始時刻造波板的水下長度。

2.3 源項消波

對于數值波浪水池，消波功能的實現有很多。如主動式消波、阻尼層消波、Sommerfeld輻射邊界消波、動量源項消波和多孔介質消波等。主動式消波[17]是指在造波的同時，通過實時的修正造波板的運動，在生成目標波浪的同時產生額外的波動抵消反射波浪。而被動式消波方法如動量源項消波，通過在動量方程中添加源項以消除反射波浪。海綿層阻尼消波[17]主要是通過在特定的計算區域設定1-2倍波長的人工阻尼區，對波浪垂向速度進行強迫衰減，同時開邊界處滿足Sommerfeld輻射條件。多孔介質消波方法[19]是一種仿物理消波法，通過在N-S方程中增加動量衰減的源項達到消波的目的。

模擬計算過程中，采用王本龍等[13]提出的消波技術，在消波區域內添加源項來消除水池尾部池壁的波浪反射，確保工作區域波浪不受干擾。消波段內波動場可寫為：

i=1方向，已加源項的動量方程按照忽略粘性的歐拉方程出發可離散為：

未添加源項的動量方程可離散成：

將（11）式代入方程聯立求解可得源項S1的表達式：

同理i=3方向的源項S3的表達式：

其中下標C代表計算值；上標N+1、N分別代表N+1、N時刻的值；λ為與空間位置有關的光滑過渡加權函數，避免速度的劇烈變化影響流場的求解精度。λ=λ（x）在數值水池消波段的表達式為：

式中：xb為消波區的起始位置，xe為消波區的末端，因此可知：

2.4 幾何模型的建立

本文采用ICEMCFD軟件建立了波浪模擬所需的幾何模型，其示意圖如圖1所示，相應的三維立體模型見圖8和圖9。設定坐標原點位于平臺重心在初始自由液面的投影處。x軸正向與波浪傳播方向一致，y軸正方向與平臺橫向方向一致，z軸沿平臺的吃水方向。該數值波浪水池長24 m，寬5 m，水深8 m，靜水面以上的高度為1 m，造波板位于X軸-10 m位置處，各區域長分別為搖板及前端過渡區5 m，工作區10 m，尾部過渡區5 m，尾部消波區4 m，在水池試驗區域的中央位置設置一虛擬浪高儀來獲取波浪的時歷。

圖1 數值波浪水池的示意圖Fig.1 Schematic show of numerical wave tank

圖2 電容式浪高儀測量系統Fig.2 Measurement system of capacitive wave probes

2.5 波浪時歷的試驗測量

為了驗證數值波浪水池的精確度，在上海交通大學海洋深水池內對波浪時歷進行了實際測量。用于該試驗研究的測量儀器主要為電容式精密浪高儀，該儀器由探頭、弓形支架、電纜以及信號輸出系統等組成，如圖2所示，其通常安裝在水池拖車下方，通過測量因水位引起的電容變化來獲取波浪時歷。試驗前，首先需對浪高儀進行標定，以確定電壓轉換系數。待水體充分平靜后，通過數據采集系統對浪高儀進行采零處理，正式采集得到的波浪時歷均減去該采零值。隨后，在控制程序中輸入波浪相關參數并轉化為搖板信號來驅動造波板實現物理造波。待波形穩定之后進行數據采集，該次試驗的采樣頻率選定為40 Hz，采樣時間3分鐘，采樣點數大約為7 200點。

2.6 網格尺度的測試

CFD計算中，網格的質量、尺度和布置方式等將直接影響數值計算結果的精度。網格質量較差勢必影響數值結果的收斂性，甚至可能造成解的發散。六面體結構網格具有存儲空間小、數據結構簡單、迭代運算速度快、網格質量高的優點，且對于VOF模型來說，排列整齊，形狀規則的六面體網格不僅可以提高求解精度，同時可以捕捉到更為光順的自由液面。因此，本文全部采用H-H型六面體結構化網格對數值波浪水槽進行離散，以確保網格質量良好且方便布置網格節點。此外，過密的網格將大大增加迭代計算所需的時間，過于稀疏的網格不僅無法捕捉到流場的相關細節，而且將引起嚴重的數值耗散，以致計算結果失真。一般而言，在進行網格的生成時需遵循的原則是：在不影響計算精度的情況下，盡量減少網格數量，提高計算的速度。因此本文在網格劃分時對一些流場梯度變化較大的區域進行了加密，而對物理量不敏感的區域采用了粗網格，如自由液面附近采用較細的網格，遠離這些區域網格間距按指數分布逐漸拉大，以起到疏密有致的效果。但即便如此，計算精度也不能完全得到保證，因此必須同試驗測量值進行對比，才能確定最終的網格劃分方式。

結合試驗室的造波能力，本文選取了陡度較大的規則波浪作為研究的對象，即波陡（H/λ，其中λ為入射波波長）為1/15，周期T為1.0 s，波高H為0.104 m，采用三種密度由疏到密的網格劃分方式對該波浪進行了數值模擬。考慮到AB段網格對數值造波是否能夠成功有重要的影響，為了減少波浪傳播過程中的數值耗散，對其等距離布置較細密的網格，并按一個波長內100個、50個和25個節點三種方式來檢驗網格的不同間距對計算結果帶來的影響。此外，考慮到自由液面附近，相應的物理量梯度變化大，在EF段同樣采用等距離的網格布置方式。且三種方案中，單個波高內的網格數分別為40、20和10個。BC段為模型區域，該段的網格參數無論對計算迭代的收斂時間和精度，還是對模型周圍的波動場細節的模擬等都起著至關重要的作用，因此該處的網格應更為細密。其余段的網格劃分均按指數分布的形式，尺度由梯度變化大的區域往梯度小的區域逐漸拉大，具體如表1所示。

表1 網格尺度的測試Tab.1 Grid size test

從圖4可以看出，網格C的計算結果與試驗值偏差嚴重，且隨著模擬時間的增加，波形因數值的耗散逐漸衰減，因此不符合高品質波浪的要求。而網格A和網格B所得到的結果均與試驗值較為接近，說明網格B的劃分方案可以使模擬結果達到要求。相比網格B，劃分最為精細的網格A并沒有帶來更好的精度，且由于網格數量巨大，將大大地消耗數值運算的時間和資源。因此，從計算的精確度和時間消耗兩方面權衡考慮，最終我們選取網格B作為數值波浪水池網格劃分的基礎。

圖5為網格B方案下，數值波浪水池中央位置處未經擾動的波形時歷曲線，造波板周圍的波形分布如圖7。通過與線性Airy波理論值對比發現，計算得到的波峰高而尖，波谷較為平坦，波浪出現了非線性的特征，已經不再滿足線性波理論。

圖3 數值波浪水池網格劃分方案的示意圖Fig.3 Schematic show of mesh strategy of numerical wave tank

圖4 不同網格劃分方案的波浪時歷分布Fig.4 Time series of wave elevation model in different mesh scheme

圖5 試驗值與線性波理論值對比Fig.5 Comparison of wave elevation between test and linear wave theory

圖6 計算值與線性波理論值的波浪譜密度對比Fig.6 Comparison of wave spectrum between calculation and linear theory

圖7 網格B方案下搖板附近的瞬時波形圖Fig.7 Instantaneous wave shape around flap wavemaker in MESH-B

采用FFT理論可將一段時歷的波列轉化為基于波頻的譜密度分布形式，從而得到波浪的各階諧頻成份以及相應的貢獻[22]。從圖6中可以看出，airy波為線性規則波，其僅包含一階諧頻成份，由于入射波周期T為1.0 s，故譜峰值出現在頻率ω約為6.28 rad/s處。此外，我們對網格B劃分方案計算得到的時歷曲線進行頻譜分析發現，譜密度曲線中出現兩個峰值：第一個峰值出現在頻率ω約6.28 rad/s位置，即主頻位置，表示一階諧頻成份（First-harmonic component）；第二個譜峰對應的頻率約為12.56 rad/s，其值的大小恰好為主頻的2倍，對應的譜峰值稱為二階諧頻成份（Second-harmonic component）。故認為該入射波具有較明顯的二階非線性特征。盡管波浪的一階線性成份所占的比例較大，其高階諧頻成份卻是影響波浪非線性爬升效應的重要因素，在以后的研究中將對此重點考察。

3 波浪繞射和爬升的數值計算

3.1 模型建立

本文在上述理論建立的數值波浪水池內，以中海油3 000 m深水半潛式鉆井平臺為研究對象，對其周圍的波浪繞射和波浪爬升現象進行了數值模擬。文章重點研究波浪散射效應（包括繞射、反射等效應）對平臺周圍波浪爬升的影響，因此設定平臺為剛性固定，以消除運動引起的輻射效應的干擾。

表2 中海油半潛平臺主要參數Tab.2 Main pariticulars of CNOOC’semisubmersible platform

首先建立了計算所需的幾何模型，如圖8和9所示。該半潛平臺包含水下浮箱、四根直立柱、橫撐和平臺甲板等，模型縮尺比為1:60，平臺實型和模型的具體尺度見表2。半潛平臺置于數值波浪水池工作區域的中央位置，平臺吃水為34.5 cm，初始氣隙（靜水面至平臺下甲板的垂直距離）為15.6 cm。

為了測量平臺周圍波面分布情況，在平臺前立柱周圍設置5排浪高儀，各排浪高儀間的角度為45°，后立柱周圍設有1排浪高儀。為了揭示波浪沿立柱爬升的過程，每排均設置四根浪高儀，其總長（R）為立柱的等效半徑。考慮到波浪越靠近立柱，波浪爬升效應越復雜，故設置浪高儀沿徑向的布置間距（r）由外到內逐漸變小，分別為0.4r/R，0.4r/R，0.193 75r/R和0.006 25r/R，且單排內浪高儀的編號依次為iA，iB，iC和iD（i為浪高儀的排號），如圖10所示。

3.2 邊界條件

圖8 半潛平臺三維數值模型Fig.8 3D numerical model of semi-submersible

圖9 計算區域的幾何模型Fig.9 Geometrical model of calculated domains

如圖9所示，水池頂部設置為壓力入口邊界條件，半潛平臺和造波板為無滑移壁面。因數值計算中，所選的波浪為迎浪狀態（與平臺艏部方向夾角為180°），且平臺沿中縱剖面對稱，因此水池內流場亦關于中縱剖面對稱，為了節省網格的數量，設定水池中縱剖面為對稱面，其余均默認為壁面。波動場采用層流模式（Laminar）求解，通過VOF方法追蹤自由面；壓力—速度耦合項采用SIMPLEC算法進行迭代計算，動量方程中的對流項和擴散項均采用二階迎風格式；VOF方程采用幾何重構法[20]（Geo-Reconstruct）求解。在數值模擬時，設定左邊界上部的搖板為運動邊界，使搖板繞其末端轉動，采用鋪層技術控制動態造波區域的網格變形、潰滅和再生來產生波浪。

圖10 虛擬浪高儀的布置位置Fig.10 Layout of virtual wave probes

3.3 網格劃分

圖11 計算區域網格分布Fig.11 Meshes of calculated domains

圖12 造波板和自由液面附近網格劃分Fig.12 Meshes around wavemaker and free surface

圖13 平臺周圍的網格劃分Fig.13 Meshes around semi-submersible

圖14 立柱周圍的網格劃分Fig.14 Meshes around columns

CFD計算中，網格離散化的質量將對數值計算結果的精度產生重要影響，本文采用分塊耦合網格劃分技術，將整個流場劃分成若干區域，各流域均采用六面體結構化網格進行離散。多塊結構化網格是單塊結構網格和非結構網格的一個折中方案，它簡單而又能處理較復雜的形狀。按照上面闡述的方案生成流場網格，部分區域相應的網格分布如圖11-14所示，整個模型網格單元在360萬左右。

3.4 計算結果與分析

圖15為Mavrakos等[21]的模型試驗中，規則波浪條件下平臺立柱周圍的波浪爬升情形。本文也得到了類似的波面分布，如圖16所示，且波浪沿立柱往上爬升過程中形成的“水舌”也極為相似，基于VOF方法求解得到的自由液面分布也較光順，因此文中采用的數值計算方法能夠較好地模擬這一非線性特征。

圖15 模型試驗中波浪爬升情形Fig.15 Wave run-up in model tests

圖16 數值模擬中立柱周圍的波浪爬升Fig.16 Wave run-up in numerical simulation

圖17為前端立柱、后端立柱迎浪面中央位置，波浪沿柱體爬升時最大波面分布情況。從圖中可以看出，由于正面迎浪，兩處的波浪爬升都較大，且隨著波浪向前傳播，最大波面升高逐漸增大，在靠近立柱壁面時波浪迅速爬升至最高位置，這與Morris-Thomas[22]以及Nielsen等[23]研究結論相一致。

圖18給出的是前后立柱壁面附近的波浪爬升時歷分布。從圖中可以看出，相比平臺前立柱壁面附近的3D位置，6D處（后立柱的壁面附近）波浪時歷曲線顯示出更強的非線性特征，且波浪爬升幅值更大，部分波浪已經爬升至下層甲板，此時上端水體仍具有一定能量，將對甲板局部產生載荷沖擊。同時結合圖24（b）可知，文中采用的計算方法可以捕捉到波浪砰擊后立柱表面時所產生的水體分離和破碎等強非線性現象。由于該類砰擊載荷可能造成平臺甲板的結構破壞，因此已成為國際上研究的熱點問題之一，今后將對該強非線性現象展開詳細分析。

圖17 最大波面升高隨r/R的變化曲線Fig.17 Maximum wave elevation along wave probe row

圖18 立柱壁面附近的波浪爬升時歷Fig.18 Time history of wave elevation near column wall

圖19顯示的是虛擬浪高儀組2（虛擬浪高儀組3順時針方向45°位置）和4處（虛擬浪高儀組3逆時針方向45°位置）的最大波面分布情況。盡管波浪爬升幅度比柱體的迎浪面中央位置小，但仍較為明顯，并且表現出較強的非線性特性。虛擬浪高儀組2處最大的波面升高規律與虛擬浪高儀組3與6處基本一致，越靠近立柱，波浪爬坡現象越明顯。然而，虛擬浪高儀組4處，最大的波浪爬升并未發生在立柱壁面附近，而是離立柱較遠的4B位置，造成這種現象的原因可能是由于波浪遇前立柱后發生繞射和反射，部分反射的水體與入射波在遠離柱體位置發生疊加所致，該現象同時由Mavrakos等[21]在模型試驗中得以證實。

此外，通過波浪時歷對比發現，盡管4B與2D兩位置處最大波面升高相接近，但在4B處，由于立柱內側區域受到的干擾因素多，使得波谷位置單個周期內出現了明顯的二次波峰，具有比2D處更強的高階諧頻非線性特征。

當波浪繞過虛擬浪高儀組2和4處繼續往后傳播時，其爬升幅值迅速減小。且一個立柱半徑R范圍內，最大波面升高基本無變化，但其非線性特征卻越明顯,如圖21、22所示。

圖19 最大波面升高隨r/R的變化曲線Fig.19 Maximum wave elevation along wave probe row

圖20 立柱壁面附近的波浪爬升時歷Fig.20 Time history of wave elevation near column wall

圖21 最大波面升高隨r/R的變化曲線Fig.21 Maximum wave elevation along wave probe row

圖22 立柱壁面附近的波浪爬升時歷Fig.22 Time history of wave elevation near column wall

對比沿平臺內側布置的浪高儀組3、4、5發現，3處的波浪爬升最為顯著，隨著波浪繞立柱繼續向前傳播，波高逐漸減小，波浪爬升效應明顯減弱，這與Nielsen等[23]在對半潛平臺的單立柱展開模型試驗時所得出的結論相同。值得注意的是，3處波浪爬升幅度雖相對較大，但波形較為規整，波谷的形狀較為對稱，隨著波浪繞立柱往前傳播至4處，波谷位置出現了二次波峰，在波浪傳至5處時，二次波峰幅值不斷的增大，同時主波峰（一個波浪周期內第一次出現的波峰）的形狀更尖瘦，該過程將使波浪非線性效應增強，與平臺間的相互作用變得更為復雜。

為了更為清晰地理解波物之間的相互作用，文中對沿平臺立柱的繞射和爬升現象進行了模擬。圖23為平臺前立柱周圍波浪的爬升過程。該情景可描述為：波浪遇立柱時，部分水流繞過立柱繼續朝前傳播，而另一部分水體因立柱的阻礙開始沿立柱爬升，此時波浪的能量很大。隨著波浪繼續往上爬升，其能量不斷耗散，興起的波峰不斷增高，波峰周圍的速度逐漸減小。待波峰增至一定位置后停止升高，在重力的作用下開始回落，但此時緊貼立柱表面的水體卻繼續往上爬升，如圖（c）所示，直至其速度接近0時，波浪開始下滑，此時可以觀察到，回落的水體沿立柱繞射至側面時將興起波浪，且其波峰周圍的速度顯著增大，如圖（d）所示。

圖24顯示的是后立柱迎浪面波浪的爬升現象。從圖（a）可以看出，后立柱周圍的波浪已上升到平臺下甲板并與之發生非線性砰擊。砰擊結束后，下甲板上的部分水流由于重力作用直接脫落，剩余部分沿立柱回落，同時有小部分水體與尾部水舌發生分離，如圖（b），這種現象可能是因水流下滑過程中，兩部分水體速度不同引起的。相比前立柱，該規則波浪下，后立柱周圍的波浪爬升情況更為嚴重，非線性效應更強，原因可能因后立柱位于前立柱的尾流中，其周圍的波面更易受擾動影響。

圖23 前立柱迎浪面波浪的爬升過程Fig.23 Wave run-up process on front side of forward columns

圖24 后立柱迎浪面波浪爬升現象Fig.24 Wave run-up process on front side of aft columns

4 結 論

文中基于N-S方程和VOF方法建立了三維數值波浪水槽，運用類似于深水試驗池的搖板造波方式產生波浪。采用H-H型六面體結構化網格對含平臺模型的數值波浪水槽進行離散，保證了網格的質量，提高了計算的精度。此外，文中還對多種網格劃分形式進行了詳細研究，同時與試驗測量值進行對比，確定了最優的網格劃分方案，并保證了數值方法的可行性。隨后，在該水池中對陡波條件下半潛平臺周圍的波浪繞射和爬升效應進行了數值模擬。結果表明，本文采用的計算方法能夠對平臺周圍的波面分布給出詳細的數值模擬，同時可預報出波浪沿平臺立柱的繞射以及從波浪開始爬升至回落的完整過程，為下一步有關平臺氣隙的模型試驗前確定浪高儀的安裝位置提供了重要參考。計算獲得的波浪爬升情形和模型試驗相似，且文中得到的重要結論與相關文獻一致，進一步證明該方法模擬波浪繞射和爬升等現象是可行的。通過計算和分析可以得到以下結論：

（1）立柱正面迎浪位置，波浪爬升較大，且越靠近立柱，波峰越尖瘦，波浪爬升越為嚴重。因此在下一步的模型試驗中，需在該位置附近安裝浪高儀進行進一步研究。

（2）一般而言，越靠近立柱波浪爬升效應越明顯，但也有例外，如偏離前立柱中縱剖面45°的內側位置附近（浪高儀組4附近），可能的原因是由于內側位置發生繞射和反射，部分反射的水體與入射波在遠離柱體位置發生疊加所致。

（3）平臺立柱的內側，波面受干擾的因素多，波浪爬升現象比相應的外側更為明顯，爬升幅度更大，成因復雜，應加以詳細研究。

（4）由于波浪沿前立柱的繞射、遇后立柱的反射以及波浪疊加等效應的相互作用，平臺前立柱的背面附近以及前后立柱之間區域的波面分布將變得較為復雜，波浪的非線性效應增強，應給予足夠重視。

[1]王世圣,謝 彬,曾恒一,馮 瑋,李曉平,張海濱.3 000米深水半潛式鉆井平臺運動性能研究[J].中國海上油氣,2007,19(4):277-280.

[2]Stansberg C T,Baarholm R,Kristiansen T.Extreme wave amplification and impact loads on offshore structures[C].Offshore Technology Conference,OTC 17487,2004.

[3]McCamy R,Fuchs R.Wave forces on piles:A diffraction theory[R].Tech.Memo No.69,Beach Erosion Board,US Army Corps of Engineers,Washington DC,USA,1954.

[4]Niedzwecki J N,Duggal A S.Wave run-up and forces on cylinders in regular and random waves[J].J Waterways,Port,Coastal,Ocean Eng,1992,118:615-634.

[5]Niedzwecki J M,Huston J R.Wave interaction with tension leg platforms[J].Ocean Engineering,1992,19(1):21-37.

[6]Lee C H,Newman J N.Second-order wave effects on offshore structures[C].Proceedings of the Seventh International Conference on Behaviour of Offshore structures,1994,2.

[7]Donald G Danmeier,Robert K M Seah,Timothy Finnigan,Dominique Roddier.Validation of wave run-up calculation methods for a gravity based structure[C].Proceedings of the 27th OMAE,2008.

[8]Teigen P,Trulsen K.Numerical investigation of non-linear wave effects around multiple cylinders[C].11th International Offshore and Polar Engineering Conference,2001.

[9]Danmeier G D,Robert K M S,Finnigan T,et al.Timothy Finnigan,Dominique Roddier.Validation of wave run-up calculation methods for a gravity based structure[C].Proceedings of the 27th OMAE,2008.

[10]Roberto M,Andrea D M,Riccardo B.Numerical simulation of the flow around an array of free-surface piercing cylinders in waves[C].Proceedings of the 26th OMAE,2007.

[11]Boo S Y,Kim C H.Simultation of fully nonlinear irregular waves in a 3-D numerical wave tank[C].Proceedings of the 17th ISOPE,1994.

[12]Brorsen M,Larsen.Source generation of nonlinear gravity waves with the boundary integral equation method[J].Coastal Engineering,1987,11(2):93-113.

[13]周勤俊,王本龍,蘭雅梅,劉 樺.海堤越浪的數值模擬[J].力學季刊,2005,26(4):629-633.

[14]董 志,詹杰民.基于VOF方法的數值波浪水槽以及造波、消波方法研究[J].水動力研究與進展A輯,2009,24(1):15-21.

[15]孫昭晨.推搖混合式造波機理論曲線[J].港口工程,1988(4):30-32.

[16]陶建華.水波的數值模擬[M].天津:天津大學出版社,2004.

[17]李雪臨,任 冰,王永學.VOF方法中主動吸收式無反射數值造波研究[C].第二十屆全國水動力學研討會文集,2006.

[18]韓 朋,任 冰,李雪臨,王永學.基于VOF方法的不規則波數值波浪水槽的阻尼消波研究[J].水道港口,2009,30(1).

[19]詹杰民,董 志.黏性數值波浪水槽的多孔介質消波方法[C].第十四屆中國海洋（岸）工程學術討論會論文集（上冊）,2009.

[20]Youngs D L.Time-dependent multi-material flow with large fluid distortion[C]//Proceedings of Numerical Methods for Fluid Dynamics.New York,1982.

[21]Mavrakos S A,Chatjigeorgiou I K,Grigoropoulos G,Maron A.Scale experiment of motions and wave run-up on a TLP model,subjected to monochromatic waves[C].Proceedings of the 23rd OMAE,2004.

[22]Morris-Thomas M T,Thiagarajan K P.The run-up on a cylinder in progressive surface gravity waves:Harmonic components[J].Applied Ocean Research,2004,26(3-4):98-113.

[23]Nielsen F G.Comparative study on airgap under floating platforms and run-up along platform columns[J].Marine Structures,2003,16(2):97-134.