劉逸凡,劉紅軍,2*,劉 燦,潘光來

(1.中國海洋大學(xué)環(huán)境科學(xué)與工程學(xué)院，山東 青島 266100；2.山東省海洋環(huán)境地質(zhì)工程重點實驗室，山東 青島 266100)

引 言

基于Biot固結(jié)理論，Madsen[1]與Yamamoto[2]最早提出線性波浪作用下多孔彈性無限厚度海床土體的靜力解析解，Hsu和Jeng[3]則進一步推出三維解析解。楊少麗等[4]計算得到了有關(guān)粉砂海床土孔壓響應(yīng)的有限差分解及極端風(fēng)浪下的液化深度。王棟等[5]考慮了土骨架和孔隙流體加速度項，提出了適用于非均質(zhì)土層及復(fù)雜邊界條件下的有限元數(shù)值分析方法。漆文剛等[6]通過解析解計算了無粘性砂質(zhì)海床土瞬時液化深度，并與Mory等[7]針對海洋結(jié)構(gòu)物周圍的瞬時孔壓原位觀測試驗相驗證。Okusa[8]、 ZEN[9]、Tsai[10]分別推導(dǎo)了一維、二維至三維等不同條件下的波致海床土瞬時液化判別準(zhǔn)則，何影等[11]、葉建紅等[12]對上述判別準(zhǔn)則的液化區(qū)分布特征作出對比分析，并分析了各判別準(zhǔn)則在海洋工程中的適用性。

有關(guān)單樁與海床相互作用方面的研究，李小均等[13]通過ABAQUS分析了樁體周圍海床土體孔隙水壓力響應(yīng)，但模型中并無實質(zhì)單樁，僅設(shè)置海床的樁土交界面處為不排水邊界。該研究的不足之處在于將樁體的作用僅歸結(jié)為一種接觸屬性的改變。胡翔等[14]在ABAQUS中進一步細化了單樁的屬性參數(shù)，但作者并未將波浪荷載通過子程序以完整、連續(xù)的函數(shù)形式施加在海床表面，而是將動荷載分為數(shù)個分析步，以靜力荷載的形式依次施加，該種加載方式在一定程度上會影響結(jié)果的準(zhǔn)確性。朱斌等[15]綜合考量了單樁及海床所受到的波浪力，對比了修正劍橋、運動硬化兩種本構(gòu)模型下樁周海床土的超孔壓累積增長過程及液化分析，但文中樁體為小尺寸(1 m)實心方樁基礎(chǔ)，彈性模量設(shè)置為30 GPa，類似混凝土性質(zhì)，同時，其對單樁施加的波浪荷載為簡化的集中荷載形式。

雖然以往學(xué)者對于殘余孔壓的累積增長過程有了一定研究，但對樁周土體瞬時孔壓響應(yīng)的研究仍相對較少。本文從極端海況條件出發(fā)，細化波浪與大直徑單樁的相互作用過程，使單樁和海床分別同步受到來自波浪的水平和豎向荷載作用，分析在一個波浪周期內(nèi)，單樁周圍土體的瞬時孔壓響應(yīng)規(guī)律及瞬時液化的可能性。

1 數(shù)值模型與驗證

1.1 數(shù)值模型參數(shù)

為分析波浪-土體-單樁三者之間的相互作用，在有限元軟件ABAQUS中建立三維土體和單樁模型，其次將波浪以荷載的形式施加在單樁和土體表面。

參照海上風(fēng)電大直徑單樁實際尺寸建立單樁模型，其中長60 m，外徑5 m，壁厚5 cm。參照鋼管樁設(shè)置其物理、力學(xué)屬性，其中密度為7.85 g/cm3，彈性模量為210 GPa，泊松比為0.2。按照C3D8R單元類型對單樁模型進行網(wǎng)格劃分。

土體模型為半徑50 m，厚50 m的圓柱體。本文研究對象為黃河三角洲區(qū)域粉質(zhì)海床土，參照馮秀麗等[16]，賈永剛等[17]所研究的埕北淺灘區(qū)域地層物理力學(xué)性質(zhì)參數(shù)，設(shè)置土體模型屬性，其中土體濕密度為1.9 g/cm3，干密度1.5 g/cm3，泊松比為0.3，彈性模量為3.5 MPa，內(nèi)摩擦角為30°，粘聚力為20 kPa，采用摩爾庫倫本構(gòu)模型表示。設(shè)置土體孔隙比為0.8，孔隙水容重為10 kN/m3,滲透系數(shù)為10-6m/s，土體模型整體屬性保持一致。土體網(wǎng)格單元類型為C3D8P(八節(jié)點六面體位移孔壓耦合單元)。單樁入土深度為30 m，樁土接觸屬性采用法向硬接觸，切向摩擦系數(shù)設(shè)置為0.3。

單樁與土體數(shù)值模型尺寸及網(wǎng)格劃分如圖1所示。

1.2 荷載邊界條件

首先為單樁和土體模型設(shè)置初始應(yīng)力狀態(tài)：對單樁和土體模型施加重力荷載，同時為土體設(shè)置初始應(yīng)力場以平衡地應(yīng)力；設(shè)置土體初始孔壓場為零，限制土體側(cè)面和底面位移。

海上風(fēng)電單樁基礎(chǔ)及周圍土體所受到的荷載包括：

荷載一：上部風(fēng)機結(jié)構(gòu)自重對單樁基礎(chǔ)產(chǎn)生的豎向荷載。本文參照5 MW風(fēng)機結(jié)構(gòu)重量為500 t，設(shè)置豎向荷載為5 MN，以集中荷載的形式施加在單樁頂部表面。

荷載二：海床表面受到的波浪豎向荷載作用。其隨波面起伏不斷變化，已知波浪表面方程如式：

η=αcos(kx-ωt)

(1)

式中，η為波面位置；α為波幅；x、t為波浪的水平位置及時間；k為波數(shù)；ω為角速度。

結(jié)合隨水深的波浪衰減系數(shù)，得到作用在海床表面的動波壓力如式(2)，該式為凈波壓力，不含靜水壓力。

(2)

式中，P0為波壓力；ρ為水體密度；h為水深；H為波高。式(1)(2)以水面為基準(zhǔn)面。

與陸地上土體表面受到的壓力不同，波浪對海床表面不僅有波壓力作用，同時賦予海床表面相應(yīng)的孔隙水壓力邊界，該孔壓邊界值與波壓力值一致，同樣采用式(2)設(shè)置。根據(jù)土體內(nèi)不同位置的孔壓分布，可以判斷孔隙水在土體中的滲流方向。三維海床泥面上所受到的波壓力按式(2),作為X軸方向上位置和時間的函數(shù)施加在單樁基礎(chǔ)周圍，因表達式(2)并未表現(xiàn)單樁的影響，因此土體表面所受波壓力與波浪遇到單樁產(chǎn)生的反射和繞流現(xiàn)象無關(guān)，更多地體現(xiàn)的是樁土接觸面帶來的影響。上述豎向波浪荷載及孔壓邊界條件以子程序的形式在ABAQUS中添加。

荷載三：單樁受到的波浪水平荷載作用。本文以線性波形式，采用莫里森方程計算波浪力。以往的研究中針對該水平荷載多以集中力的形式[18]施加在樁頂剛片上，表達式如式(3)所示。但該式為積分方式計算得到的集中力荷載，同時對式(4)中的參數(shù)進行了近似和簡化，從而忽略了波面起伏、水深等因素對波浪力的影響。

(3)

(4)

式中，F(xiàn)H為集中力荷載；D為單樁直徑；CD和CM為拖曳力系數(shù)和慣性力質(zhì)量系數(shù)，分別取1.2和2.0[19]；θ表示kx-ωt。

所研究極端海況取自黃河三角洲海域50 a一遇[20]的波浪條件，其具體波浪參數(shù)為波高6.7 m，波長77.4 m，周期8.6 s，單樁所處位置水深假設(shè)為10 m。

對于極端海況下的波浪，不僅波浪周期長，同時具有大波高，使得波峰與波谷時刻水面位置相差較大，且持續(xù)時間較長。為更好地體現(xiàn)一個周期過程中，波浪對單樁產(chǎn)生的水平荷載變化，本文根據(jù)莫里森方程[21]得到水中某一深度處單樁受到的波浪力原始表達式，如式(5)。

(5)

(6)

式中，fH為分布式荷載；ux為水質(zhì)點速度；z為計算水質(zhì)點至海床表面的距離。

對于該表達式，首先是有關(guān)波浪的函數(shù)，受波浪運動到單樁時的相位時刻、單樁所處的水平位置影響，其次同樣受到單樁的影響，根據(jù)波面起伏位置確定單樁受力范圍。因而，本文通過ABAQUS荷載自定義子程序，基于有限元原理，將該表達式所述荷載對應(yīng)施加在單樁表面所有深度位置處；最后設(shè)置條件函數(shù)，保證單樁只有位于水中的部分受到波浪水平荷載作用。如此，則單樁表面受到的水平波浪力可以隨著波面的起伏而不斷變換作用范圍與作用大小，同時波浪相位與海床所受波浪相位始終保持一致。其中部分子程序如圖2所示。

圖2 單樁所受波浪力的子程序部分

根據(jù)荷載條件共設(shè)置三組計算方案，分別是：

(1) 對照方案，在不含單樁的土體模型上施加波浪荷載二；

(2)方案一，對單樁-海床體系施加荷載一和荷載三，觀察單樁前后運動對樁周土的影響；

(3)方案二，在方案一基礎(chǔ)上，再對體系施加荷載二，觀察波浪的水平與豎向荷載共同作用下樁周土孔壓分布。其中方案一與方案二的荷載示意圖如圖3所示。

1.3 本文數(shù)值模型的驗證

針對波浪作用下海床土體內(nèi)孔隙水壓力響應(yīng)，參照以往研究中Hsu和Jeng[3]提出的解析解、Liu和Jeng[22]得到的試驗結(jié)果，在相同的波浪、土體條件參數(shù)下計算本模型的數(shù)值結(jié)果。將所得的三組數(shù)據(jù)繪制曲線如圖4所示，通過分析波致孔壓(相對)最值隨相對入土深度的變化，可以發(fā)現(xiàn)本文所使用的數(shù)值模型計算結(jié)果與以往學(xué)者得到的解析解數(shù)據(jù)和試驗數(shù)據(jù)結(jié)果具有較好的吻合性，驗證了本文數(shù)值模型計算方法的可靠性。

2 計算結(jié)果分析

2.1 未受單樁影響處海床孔壓響應(yīng)

為觀察單樁對海床土孔壓響應(yīng)帶來的影響，首先僅對波浪作用下的海床土體(不含單樁)孔壓響應(yīng)進行分析，以更好地與后期(包含單樁結(jié)果)進行對比。如圖5所示，波浪運動到該位置時為波谷狀態(tài)，因而首先出現(xiàn)負(fù)孔隙水壓力，之后隨時間呈周期性變化。如圖6所示，隨著深度的增加，孔壓峰值在有規(guī)律地遞減，以波峰時為例，由海床表面處的26.7 kPa降低至30 m處為3.2 kPa。

圖5 不同深度處海床孔壓隨時間變化曲線

圖6 孔壓峰值隨深度變化曲線

2.2 僅受波浪水平荷載影響下樁周孔壓響應(yīng)

在方案一中，輸出樁后不同深度處土體的孔壓變化曲線如圖7所示，水平位移隨時間變化曲線如圖8所示，作用時間為一個波浪周期。以20 m土深處為分界點，單樁繞該位置隨波浪發(fā)生前后運動。將圖7與圖8對應(yīng)分析可得，隨著單樁的運動，當(dāng)某一時刻土體受到樁體壓力作用，發(fā)生正向位移時，此時該部分土體呈現(xiàn)出正孔隙水壓力；若土體受到樁體拉力作用，出現(xiàn)負(fù)位移時，此時土體內(nèi)呈現(xiàn)負(fù)孔隙水壓力。

圖7 樁后不同深度土體孔壓隨時間變化曲線

圖8 樁后不同深度土體水平位移隨時間變化曲線

2.3 波浪水平與豎向荷載共同影響下樁周孔壓響應(yīng)

在方案一的基礎(chǔ)上施加海床表面波浪荷載，得到方案二。已知對于周圍土體而言，單樁的存在類似于一種邊界限制，一方面阻礙了波浪豎向荷載作用下周圍土體的水平位移幅度，使得樁周土的孔壓變化相較于圖5出現(xiàn)較大改變；另一方面則受單樁的水平運動影響，如2.2所示。兩種作用效果相疊加后，得到樁后土體的孔壓分布曲線如圖9所示。如圖9(a)所示，對于1～8 m深度范圍的土體，其孔壓峰值絕對值，均超過了表面海床孔壓峰值，維持在30～50 kPa之間，其中負(fù)孔壓峰值接近-50 kPa，出現(xiàn)在時間軸的10 s位置附近。隨著深度的增加，在2～8 s內(nèi)，觀察到孔壓曲線逐漸向負(fù)孔壓方向靠攏，并在一個波浪作用過程中觀察到在13～20 m深度范圍內(nèi)，土體內(nèi)均為負(fù)孔隙水壓力。對比圖9(a)與圖9(b)，可以發(fā)現(xiàn)隨著深度增加，孔壓幅值不僅大幅波動，相位同樣出現(xiàn)明顯的變化。

圖9 樁后不同深度土體孔壓隨時間變化曲線

同樣的方法輸出樁體前側(cè)的土體孔壓分布，其整體分布規(guī)律與樁后土體一致。如圖10(a)所示，樁前土體的負(fù)孔隙水壓力可達-83 kPa，遠小于樁后土體相同深度位置的負(fù)孔壓(-50 kPa)；如圖10(b)所示，在樁土交界面(30 m)深度位置，樁周土孔壓呈現(xiàn)明顯的波動增加，隨著遠離單樁底面，孔壓曲線則逐漸趨于平穩(wěn)。

圖10 樁前不同深度土體孔壓隨時間變化曲線

2.4 不同計算參量對樁周孔壓響應(yīng)的影響

如表1所示，以增加單樁的穩(wěn)定性為出發(fā)點，分別從土體的彈性模量(依據(jù)已有勘察數(shù)據(jù)，大幅提高土體彈性模量)、單樁入土深度、單樁頂部所受豎向荷載等角度出發(fā)，研究計算參量改變對單樁周圍土體內(nèi)孔壓分布的影響。

表1 不同變量條件下的參數(shù)設(shè)置

采用方案二進行計算，取單樁后側(cè)土體孔隙水壓力峰值為研究對象，繪制出不同條件下其沿深度的分布曲線，如圖11所示。

圖11 不同條件下樁后孔壓峰值隨深度變化曲線

以圖中原始條件曲線為對照組。如條件一所示，隨著單樁入土深度的增加，使得土體對單樁的約束作用加強，表層深度范圍樁周土體的水平位移幅值減小，從而使得該部分土體孔隙水壓力有所降低；如條件二所示，上部荷載的增加對整體深度范圍內(nèi)土體的孔壓均帶來一定降低效應(yīng)，且孔壓分布規(guī)律未受影響，可見在單樁基礎(chǔ)已有的自重條件下，不同功率參數(shù)的風(fēng)機帶來的上部結(jié)構(gòu)重量差異相比于基礎(chǔ)結(jié)構(gòu)的自重比例較小，因而在條件二下孔壓分布變化并未出現(xiàn)十分顯著的變化；如條件三所示，當(dāng)增大土體彈性模量后，土體孔壓的降低幅度最大，在10 m深度范圍內(nèi)降低了約25 kPa，可見彈性模量的增加在提高土體強度、減小周圍土體運動幅度上具有重要作用，進而對孔壓的累積表現(xiàn)出一定的抑制效果。

2.5 樁周土體瞬時液化分析

根據(jù)Zen等[9]所提出的液化判別準(zhǔn)則，當(dāng)波浪作用下土體內(nèi)部的孔隙水壓力與土體表面的孔隙水壓力之間的差值所形成的向上的“超靜孔隙水壓力”，大于初始狀態(tài)下的土體有效自重時，上部的土體便會發(fā)生瞬時液化。具體表達式如下：

(P-P0)≥γ′z

(7)

式中，P表示波浪在深度為z處的土體位置產(chǎn)生的孔隙水壓力；P0表示同樣水平位置處海床表面的波致孔隙水壓力；γ′表示土體的初始有效容重。

根據(jù)判別準(zhǔn)則，瞬時液化應(yīng)當(dāng)出現(xiàn)在波浪的波谷作用時，此時土體內(nèi)的滲流力方向向上。如圖12所示，分別輸出樁周土體的初始自重應(yīng)力以及超靜孔隙水壓力隨深度的變化曲線。根據(jù)式(7)，判斷樁體前側(cè)土體瞬時液化最大深度可達160 cm，樁后土體液化深度最大為90 cm。隨著與單樁距離的增加，土體的液化深度范圍明顯降低，以距離單樁外表面2 m水平距離為例，此時樁前土體液化深度范圍降低為85 cm，樁后土體則降低至50 cm。

圖12 樁周不同位置土體超靜孔隙水壓力隨深度變化曲線

3 結(jié)論

本文通過ABAQUS子程序的方式細化了海床與單樁所受到的波浪荷載作用，其中單樁受到的水平波浪荷載能夠隨波面的起伏而改變作用范圍，與海床表面受到的豎向荷載保持同相位，與實際情況更接近，從而分析了極端海況下單樁基礎(chǔ)周圍土體瞬時孔壓響應(yīng)。

單樁基礎(chǔ)的存在使得周圍土體的孔壓分布在幅值和相位方面均發(fā)生巨大變化。單樁的前后運動會影響周圍土體的水平位移，從而生成相應(yīng)的正負(fù)孔隙水壓力。土體彈性模量、單樁入土深度的增加均能有效降低表層深度范圍內(nèi)土體的孔壓峰值水平。樁周土體的瞬時孔壓分布具有不對稱性，可見在波浪的水平作用下，單樁的前后運動對海床土原有的孔壓響應(yīng)機制產(chǎn)生了較大影響：在波浪的水平作用與豎向作用保持相位一致的情況下，結(jié)合2.2中單樁前、后水平位移幅值可以發(fā)現(xiàn)，相比于向后運動，單樁向前運動幅值相對更大，對周圍海床土體孔壓響應(yīng)產(chǎn)生的影響也更明顯。因此，在波浪的水平與豎向共同作用下，樁前土體的負(fù)孔隙水壓力峰值明顯小于樁體后側(cè)，使得前者的超靜孔隙水壓力明顯大于后者，從而使得樁前土體的最大瞬時液化深度超過樁體后側(cè)土體。同時，瞬時液化深度范圍隨著遠離樁體而有所降低。

本文的研究尚有一定不足，如波浪以荷載形式的方法施加尚不能描述波浪沿單樁繞流的影響，因而本文并未對樁體側(cè)面的土體響應(yīng)展開分析。后續(xù)的研究會繼續(xù)深入波浪流體性質(zhì)方面的研究，進而更好體現(xiàn)波浪-樁-土三者的相互作用。