海上風電樁基水平承載力循環衰減二維簡化分析

夏華盛，張陳蓉，俞 劍，黃茂松

（1.同濟大學 巖土及地下工程教育部重點實驗室，上海 200092；2. 同濟大學 地下工程與建筑系，上海 200092）

1 引 言

海上風電項目近幾年在我國得到了快速的發展，目前海上風電選址主要還是在近海區域，從可行性和經濟性方面考慮，樁基礎仍是海上風電的首選基礎形式，海上風電樁基礎長期承受著河流、風機工作荷載等水平循環荷載的作用。針對水平循環受荷樁的問題，學術界開展了大量的工作。較早的有Matlock[1]和Reese等[2]分別對飽和軟黏土中鋼管樁和砂土中打入樁進行了水平靜荷和循環荷載試驗，提出了相應的經驗p-y曲線公式。近年來高明等[3]開展了飽和砂土中單樁模型試驗，在與現場試驗對比后，基于土的應力-應變關系給出了水平受荷樁的p-y曲線建議公式。田平等[4]根據試樁資料，提出了黏土中橫向周期荷載作用下樁基礎的p- y曲線統一法。Rajashree等[5]通過引入與循環次數相關的衰減系數，給出了循環加載后的水平極限承載力公式。Gerber等[6]進行了位移控制的現場水平循環樁試驗，發現循環荷載下樁側的最大土抗力小于對應的靜力加載，其認為循環p- y曲線中，樁側土抗力急劇上升點與前次循環加載幅值相關。

以往的研究大多集中在小數目循環荷載作用下的荷載-位移曲線，但對于海上風電樁基礎而言，設計單位更關心的是大數目長期荷載作用下樁基承載力的衰減。本文在商業有限元軟件基礎上，借鑒Zhou等[7]所建立的軟黏土不排水強度弱化模型概念，建立了水平循環受荷樁的有限元數值模擬方法。在與有限元計算結果對比驗證基礎上，提出了長期循環荷載作用下水平受荷樁基礎承載力衰減的簡化理論分析方法，從而解決了有限元法在大數目循環荷載計算中的局限性問題，為海上風電樁基礎長期安全預測提供了一定的理論依據。

2 不排水彈塑性有限元分析

針對本文中的軟黏土不排水條件，土體模型借鑒Zhou等[7]提出的基于修正Tresca屈服準則的理想彈塑性模型概念，通過累積塑性應變考慮循環荷載作用下樁周土體不排水抗剪強度的弱化。在實際的有限元數值模擬過程中，分開計算靜力加載與循環加載，靜力加載過程中，不考慮土體強度的弱化；循環加載過程中，在加載階段不考慮土體力學參數的弱化，僅累積其在加載階段產生的塑性應變，進入卸載階段，即反向加載的第一個增量步時，根據累積塑性應變對土體的彈性模量以及不排水抗剪強度進行弱化修正。隨著循環荷載的持續施加，重復上述工作，不斷累積土體產生的塑性應變和進一步修正土體參數，以此類推。本文采用的循環荷載作用下土體弱化衰減模型為

式中：Su為弱化后的土體不排水抗剪強度；Su0為土體初始不排水抗剪強度；Es為弱化后的土體彈性模量；Es0為土體初始彈性模量；δrem為完全重塑土的抗剪強度與初始土體的抗剪強度比值，即土體靈敏度的倒數；ε為循環加載時土體中累積的塑性應變；ε95為應變弱化參數，即土體發生 95%重塑時所對應的累積塑性應變。

對海上風電樁基礎中單個樁受水平循環荷載作用進行分析（見圖1）。本文有限元模型為水平方向的二維切片，對應的是水平受荷樁深層不考慮三維效應的力學狀態。考慮模型以及荷載的對稱性，建模選取一半模型。Bransby[8]建議，為消除邊界影響，水平受荷樁的土體外邊界尺寸取樁半徑的25～50倍。據此，本文水平向寬度取30倍樁徑，豎向取15倍樁徑，采用平面應變8節點四邊形單元，有限元模型包括411個單元，1300個節點，網格劃分如圖2所示，其中紅色標示部分為樁身。采用位移控制加載方式施加雙向水平循環荷載。樁土界面不考慮相對滑移，單樁截面半徑r為0.5 m，土體初始不排水抗剪強度為4.5 kPa，泊松比為0.49，靈敏度取5，即δrem取0.2。考慮剛度指數Es0/Su0、位移幅值yc/r以及應變弱化參數ε95的變化對水平受荷樁的影響，共進行了8組算例計算，循環加載次數均為15次，詳細參數見表1。

圖1 單樁結構示意圖Fig.1 Sketch of single pile

圖2 二維有限元網格Fig.2 2D FEM model

表1 各算例參數取值Table1 Parameters for 8 cases

圖3為各算例得到的歸一化p- y曲線滯回圈。圖中，pu為初始狀態靜力加載的單樁樁側水平極限抗力。由應變弱化系數ε95定義可知，ε95越大，則土體強度衰減同樣程度所需的累積塑性應變越大，土體強度的衰減速度越慢，從圖3(a)、3(g)、3(f)、3(h)可以看出，其他條件相同時，ε95越小則土體強度衰減越快，循環后樁側殘余抗力越小。由圖3(a)、3(b)可見，雖然土體剛度指數較大時p- y曲線上升速度較快，但對樁的水平抗力影響不大。

由于本文有限元計算循環荷載是位移控制加載方式，在其他條件相同時，位移幅值越大，所對應的累積塑性剪應變越大，土體強度的衰減也就越多。由圖3(a、3(e)、3(c)、3(g)可見，位移加載幅值較小時土體強度下降較慢，而位移加載幅值較大時，土體強度迅速下降并趨于穩定。

圖4為循環加載后進行靜力加載得到的p- y曲線與初始靜力加載p- y曲線對比圖。由圖可見，當加載位移幅值較大時，由于樁周大部分土體進入塑性區并繼續發展，累積的塑性剪應變較大，因此樁側水平抗力急劇下降，并且趨于0.2倍的初始極限抗力；當加載幅值較小時，由于樁周較多土體尚未進入塑性狀態，累積的塑性剪應變較少，水平抗力下降不明顯。此外，應變弱化系數對樁側循環后的水平抗力影響，加載幅值較大時影響相對更小。當加載幅值較大時，由于累積塑性應變較大，樁側水平抗力在本文循環次數下已開始趨于穩定。

圖3 水平循環受荷單樁滯回曲線（本圖縱坐標由pc/pu修改為p/pu）Fig.3 Hysteretic curves for a laterally cyclic loaded pile

圖4 初始靜力加載和循環后靜力加載p-y曲線對比Fig.4 Comparison between p-y curves on for monotonic load and post monotonic load

3 二維簡化分析方法

靜力加載時，Gerolymos等[9]建議采用p- y曲線形式為

式中：k為樁側切線模量；pu為樁側極限抗力；n為p- y曲線形態控制參數，取值范圍為 0～∞，當n＞10時p- y曲線近似為二折線。

樁側切線模量k的確定，本文借鑒 Baguelin等[10]提出的模型，將樁周土體區域分為內外兩個圓環，外層圓環半徑范圍為 r1→R，內層圓環半徑范圍為 r→r1（r為樁的半徑）。樁側外層圓環區域土體處于彈性狀態，力學特性不受樁的影響，彈性模量為Es0，內層圓環區域因樁受循環荷載作用，土體彈性模量弱化為Es。針對軟黏土不排水情況，土體泊松比v=0.5，根據Baguelin等[10]的理論推導得到單樁的彈性位移y和水平抗力p之間的關系為

式中：v為土體泊松比；β=Es/Es0；α=r1/r；α′=1/α。

由 Bransby[8]的建議，得到R/r取值在 25～50之間。式（4）得到R/r足夠大時對計算結果沒有影響，為與有限元計算建模一致，本文采用R/r=30。Klar[11]通過極限分析認為，水平受荷樁取α=2.85較為合理，并得到樁側塑性位移與擾動區域土體塑性應變關系為

樁側土體的塑性衰減仍采用與有限元一樣的表達式。在計算中實現樁側p- y曲線的弱化，需要將位移y區分為彈性位移和塑性位移。具體為：加載后對式（2）、（5）同時進行迭代運算，可得到滿足以上兩個式子的塑性位移和相對應的塑性應變。

水平循環受荷樁的簡化計算分為兩個階段：（1）加載時采用的p- y曲線不考慮土體彈性模量以及不排水抗剪強度的弱化，僅累積其在加載時產生的塑性應變ε。（2）當一次加載結束后根據前期累積的塑性應變ε對下一次加載時土體的彈性模量以及不排水抗剪強度進行弱化修正后，進行彈性卸載并用弱化后的參數進行再次加載計算。為與有限元數值模擬對應，簡化方法的計算參數取 Es0/Su0=500，ε95=30，δrem=0.2。

圖5為通過調整式（3）的形態參數n，得到與有限元較為接近的p- y曲線，確定了本文簡化方法p- y曲線的表達式。在此基礎上，分別采用簡化理論方法和有限元數值模擬進行了等荷載幅值條件下的水平循環受荷樁的計算，在循環荷載后對單樁施加靜力荷載，得到其循環后的樁側極限抗力，計算結果對比如圖6所示。圖中，pc為水平循環荷載幅值；pu為靜載條件下樁側的極限抗力；pus為循環加載結束后再次進行靜力加載得到的樁側極限抗力。有限元計算對應的循環荷載次數分別為5、10、50、100次，共進行了3組不同的循環荷載幅值條件下的分析，分別為pc/pu=0.5、0.6、0.7。

圖5 靜力加載p-y曲線Fig.5 p-y curves for monotonic load

圖6 循環次數和殘余樁側抗力的關系Fig.6 Relationships between residual lateral capacity and cycle numbers

由圖6可見，針對不同的循環荷載幅值以及循環荷載加載次數，對循環荷載后靜載條件下樁側極限抗力的弱化預測，有限元數值模擬以及簡化方法給出的結果較為符合。如 pc/pu=0.7的情況，有限元計算結果為 69次循環后樁出現了破壞，而簡化方法則為 85次，所謂破壞即指由于樁側土體的弱化，樁側極限抗力已小于循環荷載加載幅值，水平受荷樁已無法承受對應的循環荷載。通過對比3種循環荷載幅值條件下樁側極限抗力衰減與循環荷載次數的關系曲線發現，隨著水平循環荷載幅值的增加，樁所能承受的循環荷載次數在逐漸減少。這主要有兩個原因，（1）荷載幅值的增加使得土體累積的塑性應變相應增大，導致土體加速弱化；（2）循環荷載幅值的增加使得破壞時所對應的樁側極限抗力同比例增加。

鑒于有限元數值模擬時間成本大，難以進行長期大數目循環荷載次數下樁側水平極限抗力的衰減預測，本文采用提出的簡化理論分析方法計算，得到了一系列循環荷載幅值比較小的情況下樁側極限抗力衰減與循環荷載次數的關系曲線，如圖 7所示。由圖可見，當 pc/pu＞0.2時，隨著循環次數的增加，樁均出現了破壞；當 pc/pu≤0.2時，循環荷載下樁的累積位移會不斷增加，但從承載能力角度來看，樁卻趨于穩定，沒有發生破壞。這是由于本文計算中δrem取值為0.2，其定義為完全重塑土的抗剪強度與初始土體的抗剪強度比值，因此，當循環荷載幅值小于δrem時，樁側土體即使由于循環荷載作用衰減完全弱化，樁側最小極限抗力也會大于循環加載幅值，不同的只是隨著水平循環荷載幅值減小，樁側極限抗力的衰減速度會越慢，如pc/ pu=0.05的情況，循環次數為 105次時，極限承載力為靜力加載極限承載力的82%。

圖7 循環次數和殘余樁側抗力的關系Fig.7 Relationships between residual lateral capacity and cycle numbers

4 結 語

通過塑性應變的積累考慮土體不排水強度的弱化，模擬軟黏土中單樁受水平循環荷載時樁側極限抗力的衰減弱化特性。有限元分析表明，當采用等幅位移加載方式施加水平循環荷載時，土體剛度指數對于樁側水平抗力衰減影響不大，應變弱化系數對樁側水平抗力衰減的影響在位移加載幅值較小時影響更大。二維簡化理論分析發現，當采用等幅荷載加載方式施加水平循環荷載時，當荷載加載幅值小于等于0.2倍初始樁側抗力時，單樁在長期水平循環荷載作用下會趨于穩定，不會發生破壞，0.2倍系數對應的是本文算例采用的土體靈敏度的倒數。在長期水平循環荷載作用下，樁側抗力能否趨于穩定，以及所能承受的不致破壞的最大循環荷載次數與循環荷載的加載幅值以及樁側土體的力學特性均有聯系。

在海上風電設計中，必須對選址區域風電樁基所承受的長期循環荷載以及土體的力學特性進行關注和優化分析，以保證海上風電的長期安全可靠性。鑒于本文的分析僅僅是二維切片深層模型的探討，沒有進入水平受荷樁的三維分析，對樁長、樁徑等設計參數還未涉及，也是后續工作的方向。

[1]MATLOCK H. Correlation for design of laterally loaded pile in soft clay[C]//Proceedings of 2nd Offshore Technology Conference. Houston: [s. n.],1970: 577－594.

[2]REESE L C,COX W R,KOOP F D. Analysis of laterally loaded piles in sand[C]//Offshore Technology in Civil Engineering: Hall of Fame Papers from the Early Years.Houston: [s. n.],1974: 473－483.

[3]高明,陳錦珍,鄭國芳,等. 樁在側向靜、動、循環荷載下的性能研究及 p-y曲線建議公式[J]. 海洋工程,1988,6(3): 34－43.GAO Ming,CHEN Jin-zhen,ZHENG Guo-fang,et al.Studies on behavior of lateral loaded piles under static,dynamic and cyclic loadings and proposal of p-y curve formulation[J]. The Ocean Engineering,1988,6(3): 34－43.

[4]田平,王惠初. 黏土中橫向周期性荷載樁的 p-y曲線統一法[J]. 河海大學學報,1993,21(1): 9－14.TIAN Ping,WANG Hui-chu. A unified method of p-y curves for laterally and periodically loaded piles in clay[J].Journal of Hohai University,1993,21(1): 9－14.

[5]RAJASHREE S S,SUNDARAVADIVELU R. Degradation model for one-way cyclic lateral load on piles in soft clay[J]. Computers and Geotechnics,1996,19(4):289－300.

[6]GERBER T M,ROLLINS K M. Cyclic p-y curves for a pile in cohesive soil[C]//Proceedings of the Geotechnical Earthquake Engineering and Soil Dynamics－IV.California: Geotechnical Special Publication,2008: 1－10.

[7]ZHOU H,RANDOLPH M F. Computational techniques and shear band development for cylinder and spherical penetrometers in strain-softening clay[J]. InternationalJournal of Geomechanics,2007,7(4): 287－295.

[8]BRANSBY M F. Selection of p-y curves for the design of single laterally loaded piles[J]. International Journal for Numerical and Analytical Methods in Geomechanics,1999,23: 1909－1926.

[9]GEROLYMOS N,GAZETAS G. Constitutive model for 1-D cyclic soil behavior applied to seismic analysis of layered deposits[J]. Soils and Foundations,2005,45(3):147－159.

[10]BAGUELIN F,FRANK R,SAID Y H.Theoretical study of lateral reaction mechanism of piles[J]. Geotechnique,1977,27(3): 405－434.

[11]KLAR A. Upper bound for cylinder movement using“elastic” fields and its possible applicaton to pile deformation analysis[J]. International Journal of Geomechanics,2008,8(2): 162－167.