多物質(zhì)ALE方法在海嘯波沖擊方柱問題中的研究

翟 秋 ，張書劍，王華坤，張稼昊，嚴(yán)思寒

(河海大學(xué) 港口海岸與近海工程學(xué)院，南京210098)

海嘯是一種破壞性極強(qiáng)的自然災(zāi)害。一些海嘯災(zāi)害調(diào)查報(bào)告[1-2]顯示，海嘯襲擊后受災(zāi)地區(qū)建筑物的損害程度差異較大，保存完好的建筑物往往位于大型建筑物的正陸側(cè)，損壞嚴(yán)重的建筑物通常位于街道的盡頭。以上兩種現(xiàn)象被稱為屏蔽效應(yīng)和聚焦效應(yīng)，它們?cè)跒I海城市規(guī)劃、海嘯防災(zāi)減災(zāi)等方面有重要啟示作用。

目前大部分物理模型試驗(yàn)側(cè)重于海嘯波與單個(gè)樁柱或墻體構(gòu)件相互作用研究，只有少量工作涉及樁柱群和復(fù)雜建筑物的海嘯波水力特性分析。Hayatdavoodi等[3]對(duì)比研究了平板與 T梁在海嘯作用下的受力，探討了水深、波高、淹沒深度和抬高高度變化對(duì)波浪力的影響，發(fā)現(xiàn)平板上的受力一般是線性分布；楊萬理等[4]研究了門窗及屋面板洞口對(duì)低矮房屋海嘯作用力的影響，發(fā)現(xiàn)開洞率越大海嘯水平力越??；荀東亮等[5]研究了整流板、帶切縫的整流板和帶空氣孔的整流板等3種抗海嘯措施，結(jié)果表明整流板和帶切縫的整流板能減小海嘯水平力，帶空氣孔的整流板能有效減小豎向力；陳橙[6]通過模型實(shí)驗(yàn)考察了潰壩波對(duì)由平板-樁-斜坡組成的簡(jiǎn)易高樁碼頭的沖擊過程，分析了面板壓強(qiáng)隨海嘯波高和斜坡角度等因素的變化規(guī)律。

由于模型實(shí)驗(yàn)成本較高，學(xué)者們更多用數(shù)值模擬手段來開展研究工作。景旭斌等[7]采用 ALE方法，分析了漂浮物在海嘯作用下對(duì)陸上建筑的作用力，發(fā)現(xiàn)漂浮物的質(zhì)量、撞擊速度對(duì)作用力起決定性作用，而被撞建筑的剛度對(duì)作用力的影響并不大；楊志瑩等[8]用數(shù)值模擬的方法，分析了海嘯和颶風(fēng)作用下波高及淹沒系數(shù)對(duì)橋梁中T梁、箱梁受力的影響，發(fā)現(xiàn)相比于颶風(fēng)波，海嘯對(duì)橋梁主梁的作用力更大，對(duì)橋梁安全性威脅更大。Pringgana[9]采用SPH方法研究了方柱自身旋轉(zhuǎn)角度對(duì)屏蔽效應(yīng)及聚焦效應(yīng)的影響，發(fā)現(xiàn)方柱正對(duì)海嘯波時(shí)受到的波浪力最大；Yang等[10]利用ANSYS Fluent軟件，研究了潰壩波沖擊下前屋對(duì)后屋的防護(hù)作用，發(fā)現(xiàn)如果前后屋之間的間隙足夠大，則前屋對(duì)后屋幾乎沒有防護(hù)作用；Wei等[11]用GPUSPH的方法，研究了海嘯沖擊下副公路橋及防波堤對(duì)主橋的沖擊減緩作用，發(fā)現(xiàn)防波堤與主橋之間的最佳減緩距離約為當(dāng)?shù)厮畹?倍或來襲海嘯高度的13倍。

海嘯波與結(jié)構(gòu)物相互作用數(shù)值模擬涉及到流固耦合和自由液面問題，應(yīng)用較多的有SPH方法[9，11]和VOF方法[12-13]等。然而，SPH方法有著計(jì)算量大、耗費(fèi)機(jī)時(shí)、自由表面計(jì)算精度較低、張力不穩(wěn)定等問題[14]；VOF方法根據(jù)體積比函數(shù)F來構(gòu)造和追蹤自由面，在處理F的變化時(shí)稍顯繁瑣[15]。相較而言，多物質(zhì)ALE法在流固耦合及自由液面問題上具有計(jì)算效率高、穩(wěn)定性好等優(yōu)勢(shì)。景旭斌[16]曾采用多物質(zhì)ALE法對(duì)孤立波夾雜漂浮物撞擊防波堤的過程進(jìn)行了數(shù)值模擬，體現(xiàn)出該方法在流體沖擊問題中的適用性，然而其在沿海結(jié)構(gòu)物海嘯防災(zāi)減災(zāi)方面的應(yīng)用較少。

本文基于多物質(zhì)ALE方法建立了海嘯波與結(jié)構(gòu)物相互作用數(shù)值模型，研究了海嘯波對(duì)不同結(jié)構(gòu)布置形式下方柱的沖擊問題，以探明海嘯災(zāi)害中的聚焦效應(yīng)及屏蔽效應(yīng)。

1 數(shù)值模型的構(gòu)建

流體運(yùn)動(dòng)及自由液面的處理采用多物質(zhì)ALE方法，而在固體域離散時(shí)使用拉格朗日描述，固體運(yùn)動(dòng)由線彈性結(jié)構(gòu)動(dòng)力方程或剛體運(yùn)動(dòng)方程來控制，流體與固體之間的耦合采用罰函數(shù)法，基于上述算法構(gòu)建海嘯波與結(jié)構(gòu)物相互作用數(shù)值模型。

1.1 ALE描述流體運(yùn)動(dòng)控制方程

流體運(yùn)動(dòng)Navier-Stokes方程和連續(xù)方程的ALE形式可表述如下

(1)

(2)

式中：ξ為ALE坐標(biāo)，xi和xj為空間坐標(biāo)，t為時(shí)間，ρ為流體密度，vi為流體速度，ci是流體質(zhì)點(diǎn)相對(duì)于網(wǎng)格點(diǎn)的對(duì)流速度，σij為流體應(yīng)力張量，可用應(yīng)力偏量和壓力來表示

(3)

1.2 ALE多物質(zhì)材料模型

在多物質(zhì)ALE算法中，不同的材料能在計(jì)算區(qū)域的網(wǎng)格中自由流動(dòng)，一個(gè)網(wǎng)格可以存在兩種及以上的材料，在計(jì)算過程中跟蹤每種材料的邊界，并在相應(yīng)的單元中進(jìn)行物質(zhì)交換及輸送，在處理多相流與結(jié)構(gòu)之間相互作用的問題上表現(xiàn)較好。本文涉及水體及空氣兩種流體材料，水的狀態(tài)方程采用 Gruneisen方程

(4)

式中：p為壓力；E為單位體積內(nèi)能；C為us-up曲線截距；μ為體積變化率；S1，S2和S3為us-up曲線斜率系數(shù)；γ0為Gruneisen系數(shù)；α為對(duì)γ0的一階體積修正系數(shù)。

空氣的狀態(tài)方程則采用線性多項(xiàng)式狀態(tài)方程

p=C0+C1μ+C2μ2+C3μ3+(C4+C5μ+C6μ2)E0

(5)

式中：E0為初始比內(nèi)能，C0～C6為自定義常數(shù)。相關(guān)參數(shù)均可參考景旭斌[7]。

1.3 流固耦合處理方法

采用罰函數(shù)法處理流體與固體之間的耦合接觸。流體與結(jié)構(gòu)體相互接觸時(shí)，設(shè)流體為主面，結(jié)構(gòu)體為從面。每一個(gè)時(shí)步先檢查各從節(jié)點(diǎn)是否穿透主表面，沒有穿透則不對(duì)該從節(jié)點(diǎn)做任何處理；如果穿透，在該從節(jié)點(diǎn)與主表面間、主節(jié)點(diǎn)與從表面間引入大小與穿透深度和接觸剛度成正比的截面接觸力，其物理意義相當(dāng)于在其中放置一系列限制穿透的法向彈簧。

2 數(shù)值模型的驗(yàn)證

2.1 潰壩波沖擊彈性板算例

根據(jù)Rafiee和Thiagarajan[17]提供的算例建立數(shù)值水槽，如圖1所示。整個(gè)水槽由蓄水體、空氣域、壁面三部分構(gòu)成，在空氣域中放置一彈性板，其底端固定在床面。彈性板采用拉格朗日單元，材料為線彈性模型；蓄水體及空氣域采用多物質(zhì)ALE單元；水槽壁面按滑移邊界處理；對(duì)整個(gè)系統(tǒng)施加重力。

具體參數(shù)如下：蓄水體長(zhǎng)L取14.6 cm，高2L，整個(gè)水槽長(zhǎng)4L；彈性板與水體相距L，厚b取1.2 cm，高20/3b，密度ρ為2 500 kg/m3，楊氏模量E為106N/m2，重力加速度g取9.8 m/s2。當(dāng)水體右側(cè)閘門被迅速抽離之后，生成潰壩波沖擊彈性板。

圖1 潰壩波沖擊彈性板數(shù)值模型 圖2 潰壩波沖擊彈性板的SPH與多物質(zhì)ALE法結(jié)果對(duì)比Fig.1 Numerical model of elastic plate impacted by dam break Fig.2 Comparison of SPH and multi-material ALE simulation for breaking dam on elastic plate

圖3 彈性板自由端位移歷時(shí)曲線比較Fig.3 Comparison between time histories of displacement of the free end of the elastic plate

圖2給出了潰壩波自由液面的歷時(shí)變化。圖左為文獻(xiàn)[17]SPH結(jié)果，右為本文結(jié)果，可見兩者的潰壩波液面形態(tài)吻合很好。將彈性板自由端的位移歷時(shí)曲線與文獻(xiàn)[17-20]結(jié)果進(jìn)行比較，見圖3?？梢钥闯?，不同方法得到的位移在0.25 s達(dá)到最大值，之后振動(dòng)相位互有偏差，但大體趨勢(shì)基本一致，即在振動(dòng)幅值逐漸衰減波動(dòng)。與各文獻(xiàn)數(shù)據(jù)相比，本文結(jié)果處于合理的范圍內(nèi)，證明了本文數(shù)值模型的可靠性。

2.2 潰壩波沖擊立柱算例

利用Gomez-Gesteira等[21]和Crespo等[22]提供的模型實(shí)驗(yàn)數(shù)據(jù)進(jìn)一步驗(yàn)證本文算法，建立的水槽模型如圖4所示。數(shù)值水槽構(gòu)建方法同上一小節(jié)，閘門左側(cè)蓄水體高30 cm，右側(cè)床面以上設(shè)置高1 cm的淺水層，在空氣域中放置一剛性方柱，高75 cm，底端全約束。

潰壩波沖擊方柱的自由液面演變?nèi)鐖D5所示。圖6比較了柱體所受的潰壩波沖擊力?？梢钥闯?，本文數(shù)值結(jié)果與實(shí)驗(yàn)數(shù)據(jù)接近，柱體所受沖擊力計(jì)算值為34 N，實(shí)驗(yàn)結(jié)果為32.5 N。但在數(shù)值模擬中潰壩波沖擊柱體時(shí)刻比實(shí)驗(yàn)中提前約0.07 s。這是因?yàn)楸緮?shù)值模型將水槽底面作滑移邊界處理，而實(shí)驗(yàn)中存在底床摩擦，雖在閘門右側(cè)區(qū)域設(shè)置了淺薄水層以減輕水槽摩阻的影響，但無法將其完全消除，實(shí)驗(yàn)中摩擦的存在消耗了部分潰壩水體動(dòng)能，從而導(dǎo)致潰壩波沖擊立柱時(shí)刻稍晚于數(shù)值模擬，且峰值沖擊力也略有降低?？傮w上來看，基于多物質(zhì)ALE方法構(gòu)建的數(shù)值模型能較為準(zhǔn)確地模擬結(jié)構(gòu)物受到的潰壩波沖擊力。

圖4 潰壩波沖擊剛體柱數(shù)值模型 圖5 潰壩波沖擊剛體柱的各時(shí)刻狀態(tài)Fig. 4 Numerical model of rigid cylinder impacted by breaking dam Fig.5 Time evolution of dam break wave impacting a rigid cylinder

圖6 潰壩波沖擊力數(shù)值模擬與實(shí)驗(yàn)數(shù)據(jù)對(duì)比Fig.6 Comparison between numerical and experimental results of tsunami impact force

3 海嘯波沖擊不同排列下的多柱體研究

3.1 雙柱順排布置

采用多物質(zhì)ALE法對(duì)海嘯波沖擊雙柱順排布置下的剛性方柱進(jìn)行探究，數(shù)值模型如圖7所示。以潰壩波來模擬海嘯波，前柱(方柱1)距潰壩水體保持0.5 m不變。后柱(方柱2)與后水槽壁保持0.56 m的距離不變，后柱與前柱之間的距離L分別取0.2 m、0.5 m、0.8 m、1 m四種工況，即無量綱間距比L/d(d為方柱邊長(zhǎng))分別為1.67、4.17、6.67、8.33。工況1～4及工況5(單柱)如圖8所示。

圖7 雙柱順排布置數(shù)值模型Fig.7 Numerical model for two cylinders in tandem arrangement

方柱1受到的海嘯波沖擊力如圖8所示，五種工況下的沖擊力幾乎相同，最大值均為24 N左右。各工況下的海嘯波沖擊力歷時(shí)曲線差別不大，這表明方柱1受后柱影響很小。

方柱2和單柱受到的波浪沖擊力如圖9所示，四種工況下方柱2受到的波浪沖擊力最大值分別為10.5 N、21.8 N、24 N、24 N；單柱受到的最大沖擊力為25 N。工況1下的最大沖擊力比其它工況下的減小了50%以上，當(dāng)L/d由1.67增加到4.17時(shí)，最大沖擊力由單柱的42%上升到單柱的87%；當(dāng)L/d繼續(xù)增加到6.67時(shí)，最大沖擊力與單柱的最大沖擊力25 N幾乎相同，即屏蔽效應(yīng)在L/d為6.67時(shí)基本消失。

圖8 方柱1受到的沖擊力 圖9 方柱2受到的沖擊力Fig.8 Impact force on cylinder 1 Fig.9 Impact force on cylinder 2

不同工況下方柱2受力最大時(shí)刻水體自由表面形態(tài)如圖10所示。在工況1及工況2中，因受到方柱1的屏蔽作用，水體繞過方柱1向后傳播，前后柱之間有一明顯的低洼水域，沖擊在方柱2上的水體動(dòng)量比方柱1要小。在工況4中，沖擊在方柱2上的水體動(dòng)量與單柱相比幾乎相同，前后柱之間觀察不到明顯的低洼區(qū)域。

10-a 單柱10-b 工況110-c 工況210-d 工況4圖10 不同工況下方柱2受力最大時(shí)刻水體自由表面形態(tài)Fig.10 Free surface of water corresponding to the maximum force on cylinder 2

可見，在雙柱順排布置的形式下，兩柱的間距越小，前面方柱1的屏蔽作用越明顯，方柱2受到的最大沖擊力明顯減?。划?dāng)兩柱的間距超過一定范圍后，將不再產(chǎn)生屏蔽效應(yīng)。這種現(xiàn)象是因?yàn)楹[波沖擊方柱1后受到阻礙，在其后方一定區(qū)域內(nèi)形成遮蔽區(qū)域，此區(qū)域內(nèi)波速、波高都顯著減?。辉谶h(yuǎn)離方柱1一定距離之后波浪能量有所恢復(fù)，此時(shí)前后兩柱受到的最大沖擊力幾乎相同。

3.2 三柱錯(cuò)排布置

對(duì)海嘯波沖擊三柱錯(cuò)排布置下的方柱結(jié)構(gòu)物進(jìn)行探究，數(shù)值模型如圖11所示。工況1只設(shè)置1根方柱，工況2、工況3、工況4均設(shè)置3根方柱，方柱的底面中心點(diǎn)連接起來為一個(gè)正三角形，邊長(zhǎng)L分別為0.5 m、0.7 m、0.9 m，即L/d分別為2.5、3.5、4.5。對(duì)各工況下方柱1和方柱3受到的海嘯波沖擊力進(jìn)行分析，由于方柱1和方柱2對(duì)稱布置，此處不再給出方柱2的受力結(jié)果。

圖11 三柱錯(cuò)排布置數(shù)值模型Fig.11 Numerical model for three cylinders in staggered arrangement

工況1～4下方柱1受到的沖擊力如圖12所示。可見，各工況峰值沖擊力基本相同，均在106 N左右，各沖擊力歷時(shí)曲線也非常接近。工況2～4的結(jié)果顯示，隨著L/d的增加，方柱1的受力變化較小，這說明方柱1受到的最大沖擊力受周圍柱體干擾很小。

四種工況下方柱3受到的沖擊力如圖13所示，沖擊力最大值分別為106 N、112 N、120 N、106 N。工況2和工況3中方柱3的沖擊力最大值比單柱分別增大了5.6%和13.2%，工況4中方柱3的沖擊力最大值與單柱相近。基本規(guī)律是：隨著L/d的增加，方柱3的受力先增大再減小，最后和單柱接近。當(dāng)L/d在一定范圍內(nèi)時(shí)，方柱1、2的存在會(huì)使方柱3受到的峰值沖擊力增大。這是因?yàn)樗髟诜街?、2之間形成很強(qiáng)的間隙流，導(dǎo)致沖擊方柱3的水體動(dòng)量增加。當(dāng)L/d超過一定值時(shí)，間隙流逐漸消失，方柱3的受力接近單柱。

圖12 方柱1受到的沖擊力 圖13 方柱3受到的沖擊力Fig.12 Impact force on cylinder 1 Fig.13 Impact force on cylinder 3

不同工況下方柱3受力最大時(shí)刻潰壩波態(tài)如圖14所示。在工況2和工況3下，能觀察到方柱1、2之間較強(qiáng)的間隙流，但在工況4下方柱1、2間隙較大，水流的聚集效果已不明顯。

14-a 工況114-b 工況214-c 工況314-d 工況4圖14 不同工況下方柱3受力最大時(shí)刻水體自由表面形態(tài)Fig.14 Free surface of water corresponding to the maximum force on cylinder 3

4 結(jié)論

本文基于多物質(zhì)ALE法，綜合考慮了空氣、水體與結(jié)構(gòu)的耦合作用，對(duì)海嘯波沖擊彈性板及剛體柱兩個(gè)經(jīng)典算例進(jìn)行了驗(yàn)證，證明了該方法不僅能夠很好地模擬自由液面，而且能準(zhǔn)確地計(jì)算水流沖擊力和結(jié)構(gòu)動(dòng)力響應(yīng)。在此基礎(chǔ)上，應(yīng)用該方法探究了海嘯波對(duì)雙柱順排布置及三柱錯(cuò)排布置形式下方柱的沖擊力規(guī)律。結(jié)果表明，前柱受力與結(jié)構(gòu)布置形式無關(guān)，當(dāng)L/d較小時(shí)(L/d=1.67)，雙柱順排布置形式下后柱受到的波浪沖擊力相比單柱能降低50%以上，隨著L/d增大，后柱受力逐漸增大，且最終接近于單柱受到的沖擊力(L/d≥6.67)，此時(shí)意味著屏蔽效應(yīng)基本消失；三柱錯(cuò)排布置形式下后柱受到的波浪沖擊力相比單柱會(huì)增大，且隨著L/d的增加呈現(xiàn)先增大再減小的趨勢(shì)，其受力相比單柱最多可增長(zhǎng)10%以上(L/d=3.5)，此時(shí)聚焦效應(yīng)顯著。因此，濱海城市規(guī)劃及海工建筑物設(shè)計(jì)均應(yīng)考慮結(jié)構(gòu)物的布置形式，盡量利用順排布置下前方建筑物的遮擋效應(yīng)，避免錯(cuò)排布置下聚焦水流的沖擊，以使海嘯災(zāi)害減到最小。