張憶州， 廖晨聰,b， 陳錦劍,b

(上海交通大學(xué) a. 土木工程系； b. 海洋工程國家重點實驗室， 上海 200240)

隨著海洋工程技術(shù)的不斷發(fā)展，波浪荷載作用下海洋結(jié)構(gòu)物與海床土體的相互作用問題已成為海洋工程與巖土工程學(xué)科交叉領(lǐng)域的重要研究課題，因其關(guān)系到工程結(jié)構(gòu)安全而引起了工程界與學(xué)術(shù)界的廣泛關(guān)注.單樁基礎(chǔ)是風(fēng)電機等海洋結(jié)構(gòu)物的重要基礎(chǔ)形式，近20年來相關(guān)研究者對波浪荷載作用下的單樁響應(yīng)以及樁-土相互作用問題進行了研究[1-2]，但其中主要考慮的是深水條件下的線性波荷載.

在實際工程中，很多風(fēng)機樁基礎(chǔ)位于近海岸淺水區(qū)域，不宜采用線性波模型模擬其波浪荷載特性.對于相對水深(水深與波浪波長之比d/L)小于1/8的淺水波浪，應(yīng)采用由Korteway等提出的橢圓余弦波理論.相關(guān)研究包括：Isobe[3]研究了一階橢圓余弦波理論的計算公式及其應(yīng)用；肖波等[4]推導(dǎo)了推板運動軌跡的計算公式，并在造波機上制造了橢圓余弦波；Xu等[5]研究了有限深度海床中由橢圓余弦波所導(dǎo)致的滲流問題；Zhou等[6]研究了在橢圓余弦波作用下海床的液化和穩(wěn)定性問題，但有關(guān)橢圓余弦波作用下海洋樁基的動力特性研究還很少見.

隨著理論方法與有限元方法的發(fā)展，許多學(xué)者對海洋結(jié)構(gòu)物周圍的海床響應(yīng)進行了研究[7-8].張衛(wèi)平等[9]研究了樁-土耦合作用下的樁柱響應(yīng)，發(fā)現(xiàn)樁-土耦合對不規(guī)則波作用下樁柱響應(yīng)的影響較大；Zhang等[10]研究了埋置單樁對周圍海床響應(yīng)的影響，發(fā)現(xiàn)單樁的埋置將會改變樁周的孔隙水壓及液化深度；文獻[11-12]中采用數(shù)值方法研究了波浪荷載作用下不透水防波堤周圍的海床土體孔隙水壓分布及瞬時液化響應(yīng).但是，這些研究大多將樁基礎(chǔ)等結(jié)構(gòu)物視為固定剛體，而實際工程中的單樁基礎(chǔ)經(jīng)常會受到海洋中波浪和風(fēng)等周期性荷載的作用而產(chǎn)生樁身振動，進而對樁周的海床響應(yīng)造成影響.近年來，一些學(xué)者針對風(fēng)和波浪等周期性荷載作用而產(chǎn)生的樁身振動開展了研究.Jeng等[13]研究了二維海床的樁身振動問題，并對不同參數(shù)下的孔隙水壓響應(yīng)進行分析；Cuéllar[14]通過模型實驗研究了風(fēng)機樁基礎(chǔ)分別受到長期和短期循環(huán)荷載作用下的樁-土響應(yīng)問題；Hansen[15]通過實驗和數(shù)值模擬方法研究了周期荷載作用下風(fēng)機單樁的樁-土相互作用問題，但其中并未充分考慮波浪荷載的作用.因此，本文針對位于海洋中淺水區(qū)域的單樁基礎(chǔ)，采用橢圓余弦波荷載作用形式并考慮樁身的振動作用，對樁-土相互作用問題進行了動態(tài)分析，建立了包含波浪和海床的三維數(shù)值模型，并對三維模型進行了合理驗證，以期為海洋工程中受到波浪及外荷載作用的單樁基礎(chǔ)設(shè)計提供參考.

1 數(shù)值模型與計算方法

圖1 波浪-海床-單樁數(shù)值模型示意圖Fig.1 Sketch of the wave-seabed-pile numerical model

對于波浪模型，在Flow3D軟件中建立了1個尺寸為196 m×30 m×8 m的含樁數(shù)值水槽，其中包含100 m的波浪區(qū)域以及96 m的波浪吸收層，單樁放置在距水槽左側(cè)入波邊界50 m處.此外，利用COMSOL有限元軟件建立海床-單樁模型，采用偏微分方程模塊，對多孔海床以及單樁設(shè)置了2組本構(gòu)關(guān)系方程進行模擬.

1.1 控制方程

1.1.1波浪控制方程 波浪模型采用Navier-Stokes方程描述不可壓縮黏性流體運動，并選取κ-ε湍流模型模擬波浪與單樁的相互作用.由于Flow3D軟件中運用FAVOR(Fractional Area Volume Obstacle Representation)方法進行網(wǎng)格劃分而產(chǎn)生了流體面積和體積分數(shù)項，采用張量形式表達，其連續(xù)性方程和動量方程可分別寫為

式中：Ai(i,j=x,y,z)為3個方向的流體面積分數(shù)；qi(qj)表示3個方向的流體速度；VF為流體體積分數(shù)；ρ為流體密度；p為孔隙水壓；g為重力加速度；fi為3個方向的黏性加速度.

雙方程的κ-ε湍流模型表達式為

式中：κT為湍流動能；εT為湍流動能耗散率；PT為湍流動能產(chǎn)生項；GT為浮力產(chǎn)生項；DT為湍流動能擴散項；Dε為湍流動能耗散率的擴散項；C1、C2、C3均為無量綱計算系數(shù).

1.1.2橢圓余弦波 當波浪進入淺水區(qū)后，海底邊界的摩擦阻力對它的影響變得顯著，使得波浪的波形和波高等發(fā)生變化.此時，應(yīng)采用橢圓余弦波來描述波浪運動.1階橢圓余弦波的波面方程為

(5)

式中：H為波高；cn(·)為橢圓余弦函數(shù)；T為波浪周期；m為模數(shù)；K和E分別為第1和第2類橢圓積分，即

1.1.3海床控制方程 研究人員在波浪荷載作用下海床響應(yīng)的研究中發(fā)展出了基于Biot理論[16]的多種解析方法并得以應(yīng)用[17].Biot理論中的準靜態(tài)方程適用于低頻波浪荷載以及滲透系數(shù)較小的海床土體，其與動力方程計算結(jié)果的差別很小[18]，能夠保證計算精度.因此，本文采用Biot的準靜態(tài)方程研究波浪引起的海床響應(yīng)，假定海床為各向同性的且其各個方向的滲透系數(shù)相同，孔隙流體滿足達西定律.孔隙流體流動的質(zhì)量守恒方程可表示為

(6)

海床土體的變形采用多孔彈性本構(gòu)模型來描述，忽略體力和慣性項的作用，土體中力的平衡方程可寫為

(7)

式中：usi表示土體x、y和z方向的位移(即us、vs和ws)；εs為土體應(yīng)變；Gs為土體剪切模量；μs為土體泊松比.

1.1.4單樁控制方程 基于彈性理論，采用與海床控制方程類似的偏微分方程描述單樁的運動，但忽略孔隙水壓的作用，單樁的控制方程可寫為

(8)

式中：upi表示單樁x、y和z方向的位移(即up、vp和wp)；εp為單樁應(yīng)變；Gp為單樁剪切模量；μp為單樁泊松比.

1.2 邊界條件

1.2.1波浪邊界條件 如圖1所示，波浪模型的左側(cè)為入波邊界，模型的右側(cè)為出波邊界.為了保證計算穩(wěn)定性，在出波邊界前設(shè)置了波浪吸收層.上部為水與空氣的交界面，可認為其壓力等于標準大氣壓；水-土交界面采用墻面邊界，邊界上的流體法向速度為0.前后兩側(cè)面采用對稱邊界，即可認為邊界兩側(cè)存在流體，且流體的通量和剪應(yīng)力均為0.

1.2.2海床邊界條件 求解海床控制方程需設(shè)置合適的邊界條件.如圖1所示，在海床表面，孔隙水壓p等于由波浪模型所得作用于海床表面的動態(tài)波浪壓力pw，且豎向有效正應(yīng)力和剪應(yīng)力為0，即在z=0處，

(9)

海床底部是不透水的剛性邊界，其海床位移為0且無豎向流動，即在z=-h處，

(10)

海床的4個側(cè)向邊界也是不透水的，對應(yīng)方向的位移為0且無法向流動，即在x=0,W時，

(11)

在y=0,b(海床模型的寬度)時，

(12)

1.2.3單樁邊界條件 通過在樁頂產(chǎn)生周期性的水平位移來實現(xiàn)樁身振動.在模型中，單樁頂部的位移為

up=up,maxsinωt

(13)

式中：up,max為樁頂水平振動位移的峰值；ω為振動的角頻率；t為時間.

老板娘見我滿臉疑惑的樣子，便有點惋惜地說：“你不知道，原來我們柳江鎮(zhèn)還沒屁股那么大，前幾年上面說要擴建，要發(fā)展什么旅游業(yè)，就把鎮(zhèn)子建成了現(xiàn)在這個樣子。當時要擴建，就把原來挨著老鎮(zhèn)周圍的耕地和住房全占用了。被占用土地的農(nóng)民，凡兩畝以上的可以補貼點錢在鎮(zhèn)子上分一套六十平米的住房?！?/p>

單樁與海床的交界面是無滑移且不透水的，樁-土交界面上的土體位移與樁身位移相等且無法向流動，即

2 數(shù)值模型驗證

Chang等[19]研究了橢圓余弦波通過水下矩形障礙物的響應(yīng)特性，并設(shè)計了模型實驗以測量其波形.本文將上述波浪模型的計算結(jié)果與Chang等的實驗數(shù)據(jù)進行對比，以驗證所提出模型的準確性.Chang等的模型實驗中，橢圓余弦波的波浪參數(shù)分別為T=2.0 s，相對波高H/d=0.15，d=0.24 m，對比結(jié)果如圖2所示.由圖2可見，波浪模型的計算結(jié)果與實測值較吻合，從而驗證了波浪模型的準確性.

圖2 橢圓余弦波模型的計算值與實測波形數(shù)據(jù)對比Fig.2 Comparison of wave surface between present cnoidal wave model and experimental results

圖3 海床模型計算結(jié)果與海床響應(yīng)解析解的對比Fig.3 Comparsion of seabed response between present seabed model and previous analytical solution

3 數(shù)值模擬與結(jié)果分析

本文基于三維數(shù)值模型的結(jié)果對樁身水平位移、樁周海床響應(yīng)進行分析，并對波浪和海床的主要參數(shù)及樁周液化特性等進行分析.選取的標準計算參數(shù)：① 橢圓余弦波的T=6.0 s，H=1.5 m，d=5.0 m，L=43 m；② 多孔海床的h=20 m，浮重度γ′=10.0 kN/m3，k=1.0×10-4m/s，Gs=15 MPa，Sr=0.99，e=0.46，μs=0.35，靜止土壓力系數(shù)K0=0.5；③ 單樁基礎(chǔ)的D=2 m，總長度S=25 m，l=15 m，Ep=30 GPa，Gp=12.5 GPa，μp=0.2；④ 單樁樁頂?shù)膗p,max=5 mm，Tp=6 s，ω=π/3 rad/s.

3.1 樁身水平位移

采用數(shù)值模型計算了7個波浪周期內(nèi)的海床響應(yīng).圖4所示為在t處于 34.5～40.0 s的1個周期內(nèi)單樁頂部水平位移隨時間變化的示意圖.其中，位置1和位置7分別對應(yīng)于正、負位移峰值，位置3和位置5分別為正、負1/2位移峰值處，位置4為單樁位移的原點，點A表示樁前側(cè)距樁表面 0.1 m的位置，以下分析所用計算響應(yīng)均采用點A處的結(jié)果.

圖5所示為不同時刻樁身水平位移的豎向分布.由于位移的變化呈對稱分布，圖中只給出了1/2周期內(nèi)的變化情況.計算中，取up,max值為5 mm，而水-土交界面的最大水平位移約為1 mm.由圖5可見，樁身水平位移在水-土交界面以上的變化幅度較大，而在水-土交界面以下埋置部分的變化幅度很小，這是受到海床土體約束作用的緣故.

圖4 單樁頂部水平位移隨時間變化的示意圖Fig.4 Illustration of horizontal displacement on the mono-pile head in time series

圖5 不同時刻樁身水平位移的豎向分布Fig.5 Vertical distributions of horizontal displacement of mono-pile at different time

3.2 樁周的海床響應(yīng)

圖6所示為1個周期(34.5～40.0 s)內(nèi)樁周孔隙水壓的豎向分布.其中，t=34.5 s時波谷正好位于單樁頂部.由圖可見，從整體上看，樁周孔隙水壓值隨豎向深度增大而逐漸減小，且其在樁頂附近的變化速率較快，而在z=-6 m以下的變化趨勢明顯減緩.當波峰和波谷位于樁頂時，孔隙水壓的峰值分別為 4.65 和 -5.57 kPa.在樁身振動的作用下，樁周孔隙水壓值在部分時刻出現(xiàn)了局部變化，其變化范圍處于z=-0.5 m至海床表面的區(qū)域.這是由于樁身振動而導(dǎo)致樁周海床土體發(fā)生了位移，不同的位移方向使得土體受到擠壓或吸力的作用，最終使得孔隙水壓發(fā)生了局部變化.以t=40.0 s為例，海床頂部的波浪壓力為 -5.25 kPa，而在z=-0.4 m處的孔隙水壓值局部增大至 -5.57 kPa.此時，單樁由位置2向位置1運動，樁前土體產(chǎn)生了x正方向的位移，使得樁前海床土體受到吸力作用，從而導(dǎo)致局部孔隙水壓值上升.

圖6 1個周期內(nèi)樁周孔隙水壓的豎向分布Fig.6 Vertical distributions of pore pressure around pile at different time in one period

3.3 孔隙水壓分析

本文選取不同的橢圓余弦波波高、海床土體的滲透系數(shù)和飽和度來分析t=34.5 s時刻孔隙水壓的變化情況.圖7(a)所示為不同的橢圓余弦波波高下樁周孔隙水壓的豎向分布.可見：橢圓余弦波波高越大，海床表面孔隙水壓值越大，并使得樁周上部孔隙水壓沿深度方向減小的速率加快；3條曲線在z=-6 m處的孔隙水壓值基本相等，且隨著海床深度繼續(xù)增大，孔隙水壓的變化趨勢幾乎相同.

圖7(b)所示為不同的海床土體滲透系數(shù)下樁周孔隙水壓的豎向分布.由圖可見，k值的減小會使得樁頂局部孔隙水壓的變化更為明顯.這是由于k值的減小會使孔隙水壓的消散速度變慢，使其受到樁身振動產(chǎn)生的擠壓或吸力作用更加明顯.當k=1×10-3m/s時，孔隙水壓的消散速度較快，其受到樁身振動的影響較弱，因而不會出現(xiàn)孔隙水壓的局部變化.

圖7 不同的波浪和海床土體參數(shù)下樁周孔隙水壓的豎向分布Fig.7 Vertical distributions of pore pressure around pile under various waves and seabed characteristics

圖7(c)所示為不同的海床飽和度下樁周孔隙水壓的豎向分布.由圖可知：當Sr=1.00, 即海床土體完全飽和時，樁頂局部孔隙水壓的變化更加明顯，即受到樁身振動的影響更加顯著；當Sr=0.98 時，樁頂孔隙水壓幾乎不產(chǎn)生局部變化；同時，海床土體飽和度的降低會加快樁周孔隙水壓沿豎向的減小速率.

3.4 樁周海床液化分析

在波浪荷載作用下，海床土體中超孔隙水壓的增大會減小土體有效應(yīng)力，最終導(dǎo)致結(jié)構(gòu)物周圍的海床液化，使得上部結(jié)構(gòu)失穩(wěn).根據(jù)Jeng[21]提出的液化勢理論，海床液化準則可表示為

(14)

式中：γs為土體重度.

圖8所示為在樁前和樁后距單樁表面 0.1 m處，31.5～37.5 s的1個完整波浪周期內(nèi)樁周海床液化深度隨時間的變化曲線.由圖可見，樁前的最大液化深度出現(xiàn)在 34.5 s時，即波谷恰好位于單樁頂部，其液化深度的峰值位于海床表面以下 1.41 m.受樁身振動的影響，樁后的最大液化深度出現(xiàn)在 35.5 s時，其液化深度的峰值位于海床表面以下 0.98 m.

圖8 樁前和樁后的海床液化深度隨時間的變化曲線Fig.8 Time series of the liquefaction depth variation in seabed at the front and back side of mono-pile

4 結(jié)論

(1) 樁周海床的孔隙水壓在樁頂區(qū)域出現(xiàn)了局部變化，這主要是樁身振動導(dǎo)致樁周土體受到擠壓或吸力作用的緣故.同時，樁身振動還將使孔隙水壓沿豎向深度方向的減小趨勢變緩，導(dǎo)致樁周上部整體的孔隙水壓值偏大.

(2) 橢圓余弦波波高越大，樁周孔隙水壓值越大，并使樁周上部孔隙水壓減小的速率加快；海床土體滲透系數(shù)對樁周孔隙水壓分布的影響較大，土體滲透系數(shù)的減小將使樁頂附近孔隙水壓的局部變化更為明顯；當海床土體完全飽和時，受到樁身振動作用所產(chǎn)生的孔隙水壓局部變化更加明顯，土體飽和度的降低將加快樁周孔隙水壓沿深度方向的減小速率.

(3) 樁周液化深度的峰值一般出現(xiàn)在波浪波谷位于單樁頂部時.受到樁身振動的影響，樁前和樁后的最大液化深度出現(xiàn)在不同時刻，且樁前與樁后的最大深度有所不同.