于秀霞 李 欣 凌賢長 劉俊偉

(①青島理工大學， 青島 266000， 中國)

(②青建集團股份公司， 青島 266000， 中國)

0 引 言

隨著經(jīng)濟的高速發(fā)展，對于能源的需求逐漸增多，人類面臨的問題也日益增多，如環(huán)境污染、資源短缺、人口膨脹等。海洋環(huán)境中擁有大量的可利用資源，對海洋能源的開采是目前海洋學者研究的重要方向，學者們對海洋資源的開發(fā)展開了不同的研究(陳林平等， 2020； 鄧檢良等， 2020； 李守定等， 2020； 劉曉磊等， 2020)，目前超過六成的海上風電基礎(chǔ)采用單樁基礎(chǔ)的型式(王國粹等， 2011)。海上風力機單樁基礎(chǔ)可能會受到強非線性甚至破浪的破壞，波浪荷載作用下，單樁基礎(chǔ)附近的土體可能發(fā)生液化，進而加劇海上風力機的振動。因此，了解這些機理并準確預(yù)測它們對單樁基礎(chǔ)的影響在工程設(shè)計中尤為重要。

國內(nèi)外的大量學者針對波浪-海床-結(jié)構(gòu)物相互作用的課題進行了大量的研究，取得了眾多的成果(Sumer et al.，2012； Zhang et al.，2012， 2013； Sumer， 2014； 段倫良等， 2017， 2018； Jeng， 2018)。以往研究大多針對坐落式結(jié)構(gòu)物，對單樁周圍海床動力響應(yīng)問題研究較少。波浪運動到大直徑單樁基礎(chǔ)附近時會產(chǎn)生反射和繞射現(xiàn)象，但Li et al. (2011)采用有限元軟件ABAQUS建立單樁模型得到單樁周圍海床中孔隙水壓力的變化規(guī)律沒有考慮波浪與單樁的相互作用。Chang et al.(2014)建立三維數(shù)值模型模擬東海風機群樁基礎(chǔ)(8個樁基基礎(chǔ))，考慮了樁與波浪的相互作用，該數(shù)值模型能夠有效地反映樁前波浪的反射和繞射，對群樁基礎(chǔ)周圍海床動力響應(yīng)進行研究。Sui et al. (2016)采用全動力模型建立三維數(shù)值模型，研究單樁周圍海床的動力響應(yīng)與土壤位移，研究表明，波浪的反射和繞射對海床的動力響應(yīng)具有不可忽略的影響。Zhang et al. (2015， 2016)利用三維綜合數(shù)值模型，研究了單樁基礎(chǔ)對其周圍孔隙水壓力、土體應(yīng)力和液化分布的影響，繼而在海床控制方程中引入非均勻分布梯度來模擬三維非均勻土壤特性下波浪引起的海床響應(yīng)，結(jié)果表明，粗粒土的存在會降低鄰近細粒土的液化深度。以往研究中的假設(shè)單樁基礎(chǔ)為實心柱體，本文采用開口管樁模型真實反映實際工況下的海床液化規(guī)律。

Jeng(2001)、Rahman(1991)基于莫爾-庫侖破壞準則對海床的穩(wěn)定性問題進行了研究，研究表明，在波浪荷載作用下，當海床中某點所受的剪應(yīng)力超過該點的剪切強度時，海床土體會發(fā)生剪切破壞。Rahman(1994)研究表明海床內(nèi)的超孔隙水壓力是造成海床液化的主要原因。Okusa(1985)、Tsai(1995)、Zen et al.(1990)、Jeng(1997)分別對孔隙水壓力和有效應(yīng)力判斷標準提出了相應(yīng)的液化判定準則，因海床中的孔隙水壓力難以直接測量，基于孔隙水壓力的液化判斷標準應(yīng)用更為廣泛。

本文采用的液化判定準則為，當海床上存在結(jié)構(gòu)物時，海床在土體和結(jié)構(gòu)物自重作用下會經(jīng)歷一段時間達到固結(jié)狀態(tài)，會造成海床中孔隙水壓力衰減和土體有效應(yīng)力的增高等現(xiàn)象。根據(jù)Zhao et al. (2014)平均法向有效正應(yīng)力修正液化準則：

(1)

本文將建立的波浪-管樁-海床相互作用三維數(shù)值模型應(yīng)用到實際工況中，考慮了單樁自重影響的海床液化問題，分析現(xiàn)實工況作用下管樁周圍波浪場的非線性變化以及海床的動力響應(yīng)，最終得到海床的液化范圍，對數(shù)值模型在實際工程中的海床響應(yīng)進行分析評估。

1 模型建立

1.1 波浪模型

假設(shè)波浪模型為不可壓縮性和連續(xù)性流體，波浪模型由Pablo(2015)改進的olaFoam求解器求解不可壓縮連續(xù)性方程和RANS方程，olaFoam求解器是由IHFOAM演變來的，可以邊界上主動造波和消波，重力作用于Z軸方向負方向，波浪在XY平面?zhèn)鞑ィ瑢λ涂諝鈨上嗖豢蓧嚎s流進行模擬(包括質(zhì)量守恒方程和動量守恒方程)。

基于Kissling et al. (2010)VOF方程的應(yīng)用，采用體積平均VOF方程描述不同流體的運動，用于追蹤流體的自由表面。

式中：α1是VOF的指示函數(shù)，是每個單元的單位體積的含水量； 1-α1為空氣的體積分數(shù)。

(3)

使用體積分數(shù)α1可以表示任何流體性質(zhì)的空間變化。

φ=α1φw+(1-α1)φa

(4)

式中：φw和φa分別是水和空氣的任何一種性質(zhì)。

1.2 海床模型

一般來說，Biot理論在飽和多孔介質(zhì)中固體-孔隙-流體相互作用的耦合模型中應(yīng)用最為廣泛。本文所用的Biot理論忽略了固體和流體的慣性項，將Biot固結(jié)理論簡化為準靜態(tài)方程，適用于小滲透或低頻荷載作用下的固結(jié)過程。

基于質(zhì)量守恒，可壓縮多孔介質(zhì)中可壓縮孔隙流體流動的三維方程可以表示為：

(5)

(6)

式中：us、vs、ws分別是x、y、z方向上的土壤位移。

(7)

與有效應(yīng)力和孔隙壓力有關(guān)的多孔彈性介質(zhì)的總平衡方程為：

(8)

(9)

(10)

基于廣義胡克定律得到土中力平衡方程：

(11)

(12)

(13)

1.3 邊界條件

圖1為波浪模型和海床模型的邊界條件示意圖，從圖中可以看出，管樁外側(cè)的海床表面，垂直有效法向應(yīng)力和剪應(yīng)力為零，海床表面的孔隙水壓力等于波浪模型求得的動波壓力； 管樁內(nèi)側(cè)的海床表面海床孔隙水壓力為零； 在海床底部，采用不透水剛性邊界條件，其中海床位移為零，不發(fā)生垂直流動； 在海床的4個橫向邊界處，假定材料是不透水和剛性的； 在管樁與海床邊界處，管樁結(jié)構(gòu)設(shè)置為剛性不透水，管樁內(nèi)外海床孔隙壓力的法向梯度為零。此外，土壤相對于結(jié)構(gòu)沒有相對位移，即土骨架與結(jié)構(gòu)物是同步移動的。

圖1 波浪模型和海床模型邊界條件示意圖

2 模型驗證

本文模擬Wang et al. (2019)在西南交通大學進行的波浪水槽實驗，并對比不同深度處單樁周圍海床的孔隙水壓力時程曲線。水槽尺寸為60m(長)×1.8m(寬)×2.0m(高)，單樁插在尺寸為7.0m(長)×2.0m(寬)×1.0m(深)的沉淀池中，由上部結(jié)構(gòu)固定。表1為水槽實驗波浪和海床的參數(shù)設(shè)置。

表1 水槽實驗波浪和海床參數(shù)

圖2給出了波浪荷載作用下不同深度處單樁周圍海床孔隙水壓力時程曲線和波浪水槽實驗結(jié)果對比結(jié)果。從圖中可以看出，在海床表面以下0.05m、0.2m和0.35m處海床孔隙水壓力隨時間變化曲線與實驗結(jié)果吻合較好。而距離海床表面0.7m處孔隙水壓力時程曲線在波谷處和實驗結(jié)果存在些許誤差。分析原因可能為海床深度較深，海床的滲透特性發(fā)生改變，而海床模型中海床性質(zhì)是均勻不變的。圖3給出了海床最大孔隙水壓力隨深度變化的數(shù)值結(jié)果和實驗結(jié)果對比，從圖中可看出，單樁存在的條件下本文數(shù)值結(jié)果與實驗結(jié)果吻合較好。在有結(jié)構(gòu)物存在的情況下，本文建立的波浪-單樁-海床模型得出的結(jié)論與實驗結(jié)果吻合度較高，能夠準確模擬波浪作用下結(jié)構(gòu)物周圍海床的動力響應(yīng)。

圖2 不同深度處海床孔隙水壓力隨時間變化數(shù)值結(jié)果與實驗結(jié)果對比

圖3 海床最大孔隙水壓力沿深度變化數(shù)值結(jié)果

3 模型設(shè)置

3.1 數(shù)值波浪水槽設(shè)置

本耦合WSSI模型被進一步應(yīng)用于研究由于波浪作用在中等水深下管樁基礎(chǔ)附近的海床動力響應(yīng)。采用來自丹麥“波浪載荷”項目(Paulsen et al.，2014)的波浪參數(shù)，波浪場與直徑為6m的管樁相互作用。波浪模型選取的波浪參數(shù)如表2所示。

表2 波浪、海床和管樁參數(shù)

根據(jù)丹麥“波浪載荷”項目(Paulsen et al.，2014)的波浪特性，本文波浪水槽數(shù)值模型采用笛卡爾坐標系，坐標零點位于波浪水槽左側(cè)中點，波浪產(chǎn)生于水槽左側(cè)并向右側(cè)傳播，右側(cè)設(shè)置消波邊界層，降低波浪反射的影響。圖4為本文建立波浪-管樁-海床相互作用三維數(shù)值模型(WSSI模型)三維示意圖。波浪水槽是沿坐標軸方向的長度為Lx、Ly和Lz，尺寸為Lx=4lw、Ly=8D、Lz=28m，管樁位于波浪水槽中心(2lw， 4D)，Chakrabarti(1994)表明當結(jié)構(gòu)物寬度與波浪水槽寬度之比小于0.2時，可忽略水槽寬度的邊界效應(yīng)。

圖4 波浪-管樁-海床相互作用三維數(shù)值模型示意圖

3.2 海床模型設(shè)置

如圖5所示，本文海床模型尺寸與波浪模型相同，Ye et al. (2012)提出海床模型的長度大于兩倍波長即可消除固定橫向邊界的影響。管樁仍位于海床模型的中心，海床模型的側(cè)向邊界并不會影響管樁基礎(chǔ)周圍的模擬結(jié)果。在x-y平面上采用和海床模型相同的加密方式，即在管樁周圍5倍半徑范圍內(nèi)進行網(wǎng)格加密(Δx=0.59m，Δy=0.15m)，管樁內(nèi)部海床采用中間方形網(wǎng)格以達到更好的網(wǎng)格質(zhì)量(Δx=Δy=0.025m)。海床的水平方向網(wǎng)格密度對于海床孔隙水壓力的影響較小，豎向網(wǎng)格密度對海床孔隙水壓力有較大影響，在x-z方向上，海床模型分為兩部分，分別采用不同的網(wǎng)格密度，考慮到孔隙水壓力和垂向有效應(yīng)力在海床表面的變化速度較快，所以在管樁基礎(chǔ)存在的上半部分海床采用上密下疏的網(wǎng)格形式，即在管樁長度范圍內(nèi)網(wǎng)格尺寸由0.03m漸變至0.16m，樁底以下的海床尺寸為0.3m。在管樁基礎(chǔ)周圍設(shè)置6個測點，研究在波浪荷載作用下管樁基礎(chǔ)附近及管樁內(nèi)部孔隙水壓力和垂向有效應(yīng)力的分布規(guī)律，其中測點A、E距離管樁基礎(chǔ)5倍半徑處，測點B、C、D位于管樁基礎(chǔ)外緣，測點F位于管樁內(nèi)部中心處，測點A、B、D、E、F位于一條直線上。

圖5 觀測點示意圖

3.3 時間步設(shè)置

一般說來，波浪模型收斂所用的時間間隔和網(wǎng)格密度遠小于海床模型所需的時間間隔和網(wǎng)格密度。波浪模型的最小時間間隔遠小于海床模型，但考慮到單樁的存在，很難捕捉到波浪的運動，為了確保計算精度，本文采用了時間匹配的方案，波浪模型和海床模型均采用時間步長為0.001s，每50個時間步輸出為一個時間文件夾。

4 數(shù)值結(jié)果及分析

4.1 大直徑單樁周圍波浪場變化

管樁基礎(chǔ)的截面為圓弧形，當波浪運動到管樁附近時很容易產(chǎn)生圓柱繞流現(xiàn)象，波浪流速加大，會產(chǎn)生波面升高的現(xiàn)象。管樁周圍波浪場發(fā)生劇烈變化。海床的動力響應(yīng)與海床表面波浪壓力的變化具有重大聯(lián)系，波浪場的差異造成海床表面孔隙水壓力的差異，首先對管樁周圍波浪場進行分析研究。圖6為波峰和波谷運動到管樁基礎(chǔ)附近的波浪形態(tài)，由圖6可以發(fā)現(xiàn)，由于單樁的存在，波浪運動到管樁前側(cè)時會發(fā)生反射和繞射，波面形態(tài)發(fā)生明顯擾動，當波峰運動到樁前時，這種現(xiàn)象更加明顯。

圖6 波高和波谷運動到管樁基礎(chǔ)附近的波浪形態(tài)

圖7給出了管樁周圍各測點在有無結(jié)構(gòu)物波浪自由表面時程曲線對比。如圖所示，在管樁基礎(chǔ)迎浪側(cè)，越靠近管樁，波浪形態(tài)發(fā)生的變化越大，具體表現(xiàn)為：波峰變尖，波谷加深，在波谷處發(fā)生“跳躍”現(xiàn)象； 管樁周圍垂直于波浪傳播方向處波峰減小，波谷加深； 管樁背浪側(cè)波谷變緩，波峰變陡，管樁基礎(chǔ)附近波浪發(fā)生明顯的非線性變化。

圖7 管樁周圍各測點處波浪自由表面時程曲線對比

4.2 大直徑單樁周圍海床動力響應(yīng)

圖8和圖9分別為管樁周圍各測點孔隙水壓力最大值和垂向有效應(yīng)力最大值隨深度變化曲線，由圖8可知，管樁周圍各測點和管樁內(nèi)部的孔隙水壓力隨深度的變化規(guī)律與上一節(jié)規(guī)律相似，管樁外側(cè)孔隙水壓力在淺層海床中(z/h<0.1)衰減速度較快，當z/h>0.1時，孔隙水壓力衰減速度明顯減小，幅值約降為海床表面波浪壓力的1/9。管樁附近樁后D點和樁側(cè)處C點的孔隙水壓力衰減速度較快，而樁前處B點的孔壓衰減速度較慢。由于單樁的存在，B點和D點處的孔隙水壓力衰減速度大于A點和E點處。根據(jù)液化判斷準則可知，海床液化發(fā)生的必要條件是波浪作用下產(chǎn)生向上的壓力梯度，壓力梯度的產(chǎn)生與孔隙水壓力的衰減速度成正相關(guān)，孔隙水壓力隨深度的衰減變化速度越快，產(chǎn)生向上的壓力梯度越大，液化的可能性也越大。由圖8可知，當波谷運動到管樁附近時，樁后和樁側(cè)發(fā)生液化的可能性大。由圖9可知，由于單樁的存在，靠近管樁樁壁的位置，在海床表面存在較大的有效應(yīng)力，在管樁底部附近存在應(yīng)力集中現(xiàn)象，這種現(xiàn)象在距離樁底較遠的B點和E點并不明顯，此現(xiàn)象很可能是因為樁底處海床土體和樁底的彈性模量變化較大造成的。

圖8 各測點處孔隙水壓力最大值隨深度變化規(guī)律

圖9 各測點處有效應(yīng)力最大值隨深度變化規(guī)律

圖10和圖11分別為波峰和波谷時管樁周圍孔隙水壓力場和有效應(yīng)力場分布，由圖10和圖11可知，當波峰和波谷作用在單樁附近時，海床中的孔隙水壓力場和有效應(yīng)力場產(chǎn)生較大的變化趨勢，在海床深度方向上產(chǎn)生較大的變化梯度； 孔隙水壓力場變化明顯的區(qū)域主要集中在淺層海床； 在深層海床，孔隙水壓力場變化并不明顯。隨著海床深度的增加，孔隙水壓力場和有效應(yīng)力場變化趨勢逐漸減小，但同一深度處差異依然明顯。管樁底部區(qū)域出現(xiàn)明顯的應(yīng)力集中現(xiàn)象，分析原因為樁底附近的楊氏模量變化較大。另外，管樁內(nèi)部海床孔隙水壓力基本不變，有效應(yīng)力主要表現(xiàn)為橫向的差異，垂向的幾乎沒有變化，分析原因為模型沒有完全耦合。

圖10 波峰作用下樁周孔隙水壓力場和有效應(yīng)力場分布

圖11 波谷作用下樁周孔隙水壓力場和有效應(yīng)力場分布

4.3 液化深度分析

圖12給出了管樁周圍海床液化范圍的變化過程，圖12a和圖12b分別表示波谷即將傳播到管樁位置處管樁附近海床的液化范圍，圖12c表示波谷恰好運動到管樁位置處管樁附近海床的液化范圍，圖12d表示波谷離開管樁位置處管樁附近海床的液化范圍。由圖12a可知當波谷即將傳播至管樁位置處時，此時管樁周圍海床的液化范圍由管樁左右分開，即樁前為液化區(qū)域，樁后為非液化區(qū)域。由圖12c可知當波谷正好運動到管樁位置時，管樁附近大部分范圍內(nèi)海床均發(fā)生液化，此結(jié)果與波谷處容易發(fā)生液化規(guī)律吻合。隨著波浪的傳播，受波浪壓力的變化，管樁附近的海床土體液化深度先增大后減小。另外發(fā)現(xiàn)，樁基礎(chǔ)周圍液化帶的分布是對稱的，管樁樁壁附近一定范圍內(nèi)存在非液化區(qū)域，管樁附近海床的最大液化深度達到10m，樁底下方海床不會發(fā)生液化，最大液化深度發(fā)生在垂直于波浪運動方向。

圖12 管樁周圍海床液化范圍發(fā)展過程

5 結(jié) 論

本文將建立的波浪-海床-管樁相互作用數(shù)值模型應(yīng)用到工程實際中，研究海上風電管樁基礎(chǔ)周圍海床的振蕩響應(yīng)； 從樁基礎(chǔ)附近及管樁內(nèi)部提取6個特征點，對管樁附近的波浪形態(tài)、孔隙水壓力和有效應(yīng)力方面進行研究，最終探索管樁周圍海床的液化深度變化規(guī)律，得到以下結(jié)論：

(1)波浪運動到管樁附近時，波浪形態(tài)會發(fā)生巨大變化，在管樁基礎(chǔ)迎浪側(cè)的波浪形態(tài)變化最大，具體表現(xiàn)為：波峰變尖，波谷加深，在波谷處發(fā)生“跳躍”現(xiàn)象。

(2)分析各個測點的孔隙水壓力沿深度的垂向分布發(fā)現(xiàn)，在波谷運動到管樁前側(cè)時，管樁基礎(chǔ)前側(cè)的孔隙水壓力衰減速度最快，產(chǎn)生的孔隙水壓力最大，發(fā)生液化的可能性最大。

(3)各個測點的有效應(yīng)力隨深度變化曲線研究發(fā)現(xiàn)，管樁周圍海床有效應(yīng)力隨著深度的增加逐漸增大，這剛好與孔隙水壓力的變化趨勢相反?？紫端畨毫龊陀行?yīng)力場在水平方向隨時間變化不大，深度方向隨時間具有較大變化，特別是在波峰和波谷處變化明顯。

(4)管樁基礎(chǔ)周圍海床的液化規(guī)律為：波谷運動到管樁附近時，海床發(fā)生液化，最大液化深度達到10m，最大液化深度發(fā)生在垂直于波浪運動方向，管樁附近海床液化深度先增大后減小。