-
在修正剑桥模型(MCC)的基础上,姚仰平等[36-42]引入了当前屈服面和参考屈服面,超固结参数R,潜在强度应力比Mf,统一硬化参数H等概念,开发了统一硬化(UH)本构模型。该模型使用与修正剑桥模型完全相同的6个参数,分别为:极限状态应力比M(三轴压缩路径),有效应力泊松比ν,等向压缩线斜率λ,等向回弹线斜率κ,e-lnp空间各向同性正常压缩线(NCL)截距N(p=1 kPa时的孔隙比e)以及初始孔隙比e0。以上参数都可以通过常规的土工试验获得,故该模型具有较强的实际应用潜力。因为UH模型是一个有效应力本构模型,文章中若无特殊声明所有应力均为有效应力。
-
UH本构模型共设计了2个屈服面,如图1所示。其中参考屈服面即修正剑桥模型屈服面,该屈服面以塑性体积应变$ \varepsilon _v^p $作为硬化参数;当前屈服面与参考屈服面几何相似,但以统一硬化参数H作为硬化参数。图1中点A称为当前应力点,与之对应的点$\bar A$称为参考应力点,即当前应力点所对应的应力比延长线与参考屈服面的交点。
该模型遵从相关联流动法则,即屈服函数$\bar F$与塑性势函数$\bar G$相同。因此,参考屈服面可以表达为:
$$ \bar F = \bar G = \ln \frac{{\bar p}}{{{{\bar p}_x}}} + \ln \left( {1 + \frac{{{{\bar q}^2}}}{{{M^2}{{\bar p}^2}}}} \right) = 0 $$ (1) 式中:
$\bar p$ ——参考应力点$\bar A$的平均主应力(kPa);
${\bar p_x}$ ——参考屈服面与平均主应力轴(p轴)交点的平均主应力(kPa);
$\bar q$ ——参考应力点$\bar A$的广义剪应力(kPa);
$M$ ——极限状态应力比。
引入硬化参数$ \varepsilon _{\mathrm{v}}^p $,参考屈服函数可以改写为:
$$ \bar F = \bar G = \ln \frac{{\bar p}}{{{{\bar p}_0}}} + \ln \left( {1 + \frac{{{{\bar q}^2}}}{{{M^2}{{\bar p}^2}}}} \right) - \frac{1}{{{c_p}}}\varepsilon _{\mathrm{v}}^p = 0 $$ (2) 式中:
${\bar p_0}$ ——先期固结平均主应力(kPa);
${c_p}$ ——调整系数,可以表示为${c_p} = \left( {\lambda - \kappa } \right)/\left( {1 + {e_0}} \right)$,其中${e_0}$为初始应力比;
$ \varepsilon _{\mathrm{v}}^p $ ——塑性体应变,即参考屈服面的硬化参数。
而当前屈服面与参考屈服面几何相似,其屈服函数F可以表达为:
$$ F = G = \ln \frac{p}{{{p_x}}} + \ln \left( {1 + \frac{{{q^2}}}{{{M^2}{p^2}}}} \right) = 0 $$ (3) 式中:
p ——当前应力点A的平均主应力(kPa);
px ——当前屈服面与平均主应力轴(p轴)交点的平均主应力(kPa);
q ——当前应力点A的广义剪应力(kPa)。
同样地,引入硬化参数H,当前屈服面的屈服函数可以改写为:
$$ F = G = \ln \frac{p}{{{p_0}}} + \ln \left( {1 + \frac{{{q^2}}}{{{M^2}{p^2}}}} \right) - \frac{1}{{{c_p}}}H = 0 $$ (4) 式中:
p0——初始平均主应力(kPa);
H ——统一硬化参数,将在1.4节具体介绍。
-
超固结参数R的作用与常用的超固结比OCR相似,用于表示土体的超固结程度。该参数由图1中的当前应力点A和参考应力点$\bar A$的平均主应力之比表示:
$$ R = \frac{p}{{\bar p}} = \frac{q}{{\bar q}} $$ (5) 由式(2)可知$\bar A$的平均主应力$\bar p$为:
$$ \bar p = {\bar p_0}\left( {\frac{{{M^2}}}{{{M^2} + {\eta ^2}}}} \right){{\mathrm{e}}^{\frac{{\varepsilon _{\mathrm{v}}^p}}{{{c_p}}}}} $$ (6) 式中:
$\eta $——调动应力比,可以表示为$\eta = {q \mathord{\left/ {\vphantom {q p}} \right. } p}$。
将$\bar p$的表达式代入式(5)可得R的表达式为:
$$ R = \frac{p}{{{p_0}}}\left( {1 + \frac{{{\eta ^2}}}{{{M^2}}}} \right){{\mathrm{e}}^{ - \frac{{\varepsilon _{\mathrm{v}}^p}}{{{c_p}}}}} $$ (7) -
姚仰平等[36, 40]将传统的两个直线强度边界线:零拉应力线和Hvorslev线使用一条连续的抛物线统一表达。这条抛物线即UH本构模型的剪切强度包络线,称为改进的Hvorslev线,如图2所示。
确定UH本构模型的剪切强度包络线后,潜在强度应力比Mf可以由超固结参数R表示:
$$ {M_{\mathrm{f}}} = 6\left[ {\sqrt {\frac{k}{R}\left( {1 + \frac{k}{R}} \right)} - \frac{k}{R}} \right] $$ (8) 式中:
k ——极限状态应力比M的函数,$k = {M^2}/ \left[ {2\left( {3 - M} \right)} \right]$。
-
在本构模型的理论体系中,硬化参数是一个重要的概念。硬化参数的巧妙设计可以将屈服面的流动(即应力发展)与土体变形情况(即应变发展)之间的关系准确的表达出来。姚仰平等[37]在搭建UH本构模型时,设计了由潜在强度应力比Mf,当前调动应力比$\eta $和极限状态应力比M 3个参数共同参与的硬化参数。
其中潜在强度应力比Mf表达的是土体潜在的抵抗最大剪切破坏的能力,该参数决定了峰值应力比的大小,进而决定了应变硬化和软化的分界线;调动应力比$\eta $则刻画了当前应力状态的具体位置;极限状态应力比M决定了极限状态出现的时机。若调动应力比$\eta $高于极限状态应力比M则应力状态落在屈服面顶点的左侧,土处于剪胀状态;若调动应力比$\eta $低于极限状态应力比M则应力状态落在屈服面顶点的右侧,土处于剪缩状态。
结合了以上3个应力比参数的统一硬化参数H可以表示为:
$$ H = \int {\mathrm{d}}{H = \int {\frac{{M_{\mathrm{f}}^4 - {\eta ^4}}}{{{M^4} - {\eta ^4}}}{\mathrm{d}}\varepsilon _{\mathrm{v}}^p} } $$ (9) 为方便引用,令${\varOmega} = \left( {{M^4} - {\eta ^4}} \right)/\left( {{M_{\mathrm{f}}^4} - {\eta ^4}} \right)$。
当η<M<Mf时,土处于剪缩状态,此时式(9)的分数项为正数,且$ {\mathrm{d}}\varepsilon _{\mathrm{v}}^p $为正,故dH为正,当前屈服面和参考屈服面都处于硬化状态,将持续膨胀,当前屈服面追赶参考屈服面;当M<η<Mf时,土进入剪胀阶段,此时式(9)中分数项为负数,且$ {\mathrm{d}}\varepsilon _{\mathrm{v}}^p $为负,故dH为正,当前屈服面处于硬化状态,继续发生膨胀,而参考屈服面处于软化状态,发生收缩,此时当前屈服面继续追赶参考屈服面;当M<Mf<η时,式(9)中分数项为正数,且$ {\mathrm{d}}\varepsilon _{\mathrm{v}}^p $为负,故dH为负,当前屈服面和参考屈服面都处于软化阶段,都发生收缩。
-
UH本构模型的强度各向异性特征,是由其选用的空间滑动面强度准则(SMP准则)[34]决定的。这一准则在π平面上的破坏包络线如图3实线所示,呈一圆角三角形。
矢径$ \overrightarrow {OA} $的大小会跟随应力洛德角$\theta $发生变化,而$\overrightarrow {OA} $的大小即代表了A点所对应的应力路径的破坏强度,由此SMP准则可以很好的对材料的强度各向异性做出描述。但是将这样一个各向异性强度准则直接应用于岩土材料本构模型是较困难的,需要通过应力变换将真实应力空间具有各向异性的SMP强度准则转为变换应力空间中各向同性的广义Mises强度准则。变换应力可以表示为:
$$ {\sigma }_{变换}=p+\frac{{q}_{{\mathrm{s}}}}{q}\left(\sigma -p\right) $$ (10) $$ {q_{\mathrm{s}}} = \frac{{{I_1}{I_2} - 9{I_3} + 3\sqrt {\left( {{I_1}{I_2} - {I_3}} \right)\left( {{I_1}{I_2} - 9{I_3}} \right)} }}{{4{I_2}}} $$ (11) 式中:
σ变换 ——变换应力,遵从变换应力空间中各向同性的广义Mises强度准则(kPa);
σ ——真实应力,遵从真实应力空间中各向异性的SMP强度准则(kPa);
qs ——真实应力空间中用SMP准则表示的三轴压缩路径($\theta = {0^ \circ }$)的偏应力(kPa);
$ {I_i}{\text{ }}(i = 1,2,3) $ ——真实应力空间的应力不变量。
-
在数值模拟过程中,首先需要为其6个本构参数赋值,这些本构参数具体物理意义和表示符号见表1。
表 1 UH本构参数物理意义
Table 1. Physical meanings of UH constitutive parameters
表示符号 物理意义 M 极限状态应力比(三轴压缩应力路径) ν 泊松比 κ 等向回弹线斜率 λ 等向压缩线斜率 e0 初始孔隙比 N e-lnp空间正常压缩线(NCL)截距 -
如上文所述,UH本构模型继承了修正剑桥模型的基本特点,并且借助当前屈服面和参考屈服面之间的动态关系以及潜在强度Mf的变化来描述超固结土的硬化、软化、剪胀、剪缩特性和不同应力路径下的应力-应变过程。为更加便捷地将该本构模型应用于岩土工程边值问题,本研究对其进行了下述面向实际工程需求的修正。
-
模型参数“初始孔隙比e0”需用户根据土体的先期固结应力和当前应力状态进行计算,同一土层初始孔隙比设置为常数。而实际工程中土体的孔隙比随深度变化。为模拟孔隙比随深度变化,需将地层划分为若干子层,对每个子层赋恒定孔隙比。使用该方法非常麻烦且不精确。为此,我们将参数“初始孔隙比e0”替换为“超固结比OCR”。超固结比可以方便地从一维固结试验中测得,也可以从静力触探解析中得到。任一深度处土体的初始孔隙比e0的换算由用户子程序(UMAT)按照式(12)自动进行。
$$ \begin{gathered} e = N - \lambda \ln ({p_{{\mathrm{OCR}}}}) - (\lambda - \kappa )\ln (1 + \frac{{{\eta ^2}}}{{{M^2}}})+ \\ \kappa \ln (\frac{{1 + 2{K_0}}}{{1 + 2{{({K_0})}_{{\mathrm{OCR}}}}}} \cdot {\mathrm{OCR}}) \\ \end{gathered} $$ (12) 式中:
N ——e-lnp空间正常压缩线(NCL)截距,即p=1 kPa时的孔隙比e;
${p_{{\mathrm{OCR}}}}$ ——超固结土先期固结应力状态下的平均主应力,${p_{{\mathrm{OCR}}}} = \left[ {{{\left( {1 + 2{K_0}} \right)} \mathord{\left/ {\vphantom {{\left( {1 + 2{K_0}} \right)} {\left( {1 + 2{{\left( {{K_0}} \right)}_{{\mathrm{OCR}}}}} \right)}}} \right. } {\left( {1 + 2{{\left( {{K_0}} \right)}_{{\mathrm{OCR}}}}} \right)}}} \right] \cdot {\mathrm{OCR}} \cdot p$,其中p为超固结土当前固结应力状态下的平均主应力;
$\eta $ ——正常固结状态应力比,$\eta = \left[ {3\left( {1 - {K_0}} \right)} \right]/ (1 + 2{K_0}) $对正常固结粘土可假设${K_0} = 0.6$,$\eta = 0.545$;
$({K_0})_{\mathrm{OCR}}$——超固结状态静止土压力系数,$({K_0})_{\mathrm{OCR}} = {K_0} \cdot {{\mathrm{OCR}}^{0.4}}$ [43]。
式(12)即一般应力状态下初始孔隙比计算方法。对于各项同性正常固结状态,式(12)会退化为$e = N - \lambda \ln p$;对于各向异性正常固结土,即OCR=1的K0固结状态,需要在各向同性固结状态的基础上进行K0固结修正,即式(12)中的第二项“$ - (\lambda - \kappa )\ln (1 + {\eta ^2}/{M^2}) $”,孔隙比$ e = N - \lambda \ln (p) - (\lambda - \kappa ) \ln (1 + {\eta ^2}/{M^2}) $;对于各向异性固结超固结状态,则需要在各向异性正常固结的状态基础上再进行超固结度修正,即式(12)中的第三项“$ + \kappa \ln \left[ {\left( {1 + 2{K_0}} \right)}/ \left( {1 + 2{({K_0})}_\mathrm{OCR}} \right) \cdot {\mathrm{OCR}} \right] $”,这样就得到了一般应力状态下的孔隙比计算方法。
-
如式(13)、式(14)所示,与修正剑桥模型一致,UH模型的弹性体变模量K和剪切模量G是参数ν和κ的函数:
$$ K = \frac{{(1 + e) \cdot p}}{\kappa } $$ (13) $$ G = \frac{{3(1 - 2\nu )}}{{2(1 + \nu )}}K $$ (14) 由式(13)、式(14)可得到应力水平相关的弹性体变模量和剪切模量。然而,实践中发现由式(14)计算的剪切模量与土体真实的初始剪切模量存在较大差距。以参数κ为例,虽其具有明确的物理意义,但是通过一维压缩试验确定时常会有较大的误差,在模型参数标定过程中通常将其当作一个可以调整的参数。此外,有效应力泊松比ν的取值也会显著影响剪切模量的计算结果。因此,若要将UH模型应用于对初始剪切模量G0非常敏感的工程问题,需格外小心。
在本研究中,参照Andersen[44]基于大量土工试验数据得到的初始剪切刚度经验公式对UH模型进行修正,将剪切模量G与本构参数ν和κ解耦,修正后的模型不需要再输入本构参数ν。弹性体变模量依然用式(13)计算,而弹性剪切模量由式(15)、式(16)所示:
$$ G = (30 + \frac{{75}}{{\frac{{{I_{\mathrm{p}}}}}{{100}} + 0.03}}) \cdot {\mathrm{OC}}{{\mathrm{R}}^{0.5}} \cdot {\sigma _{{\mathrm{ref}}}} $$ (15) $$ {\sigma _{{\mathrm{ref}}}} = 100 \cdot {(\frac{{{\sigma _{{\mathrm{vc}}}}}}{{100}})^{0.9}} $$ (16) 式中:
${I_{\mathrm{p}}}$ ——塑性指数(%);
$ {\sigma _{{\mathrm{ref}}}} $ ——参考竖向应力;
$ {\sigma _{{\mathrm{vc}}}} $ ——竖向固结应力。
采用上述修正后,用户无需输入有效应力泊松比ν,但此时需要输入塑性指数Ip。
为验证修正结果,本节将进行土单元模拟试验并对比修正剪切模量的UH模型和原始UH模型的预测结果。数值模拟过程中UH模型取姚仰平等[41]由Weald Clay试验数据标定得到的本构参数(表2中前6项),而UH_G0模型添加本构参数Ip(表2中第7项)。所有土单元试验中,有效竖向固结应力${\sigma _{{\mathrm{{\mathrm{vc}}}}}} = 100\;{\text{kPa}}$,有效水平固结应力${\sigma _{{\mathrm{rc}}}} = {K_0}{\sigma _{{\mathrm{vc}}}}$。本研究中取正常固结静止土压力系数${K_{{\mathrm{0NC}}}} = 0.6$,不同超固结度下静止土压力系数${K_0} = {K_{{\mathrm{0NC}}}} \cdot {\mathrm{OC}}{{\mathrm{R}}^{0.4}}$。
参数 取值 极限状态应力比M 0.87 泊松比ν 0.2 等向回弹线斜率κ 0.035 等向压缩线斜率λ 0.093 e-lnp空间正常压缩线截距N 1.06 超固结比OCR 见图中标注 塑性指数Ip 25 修正剪切模量后的UH模型与原始UH模型在各超固结度水平下(OCR=1, 2, 4, 8, 12, 16, 20, 50)的不排水三轴压缩(TC)和三轴拉伸(TE)的模拟结果对比如图4所示,为方便区分,将修正剪切模量后的UH模型记为UH_G0模型。图4(c)中p-q坐标经归一化处理,将p和q分别除以等效应力pe[42]。
图 4 UH模型与UH_G0模型三轴压缩(拉伸)模拟对比
Figure 4. Comparison between UH model and UH_G0 model for triaxial compression (tensile) simulations
从图4(a),图4(b)对比结果可以看出,UH模型和UH_G0模型在不同的超固结度水平下强度和孔压模拟结果差异很小;从图4(c)可以发现修正后模型对应力路径的预测也没有受到影响。
图4的模拟结果对比表明对剪切模量G的修正没有改变模型对强度和孔压的预测结果,仅影响了应力-应变曲线的刚度。图5展示了较小应变范围内应力-应变关系的对比,图5中各图横坐标为剪切应变,纵坐标为使用各自剪切模式下不排水剪切强度(Su)的归一化剪切应力。
图5中灰色虚线(图例G0)表示按照式(15)、式(16),取竖向固结应力${\sigma _{{\mathrm{vc}}}} = 100{\text{ kPa}}$计算得到的初始剪切刚度理论值。从图5的不同应力路径刚度对比结果可发现,在不同的应力路径上UH_G0模型对土体的初始刚度预测都得到了提升,并且可以发现UH_G0模型对应力-应变曲线初始刚度的影响在超固结比水平较高时更为显著。
整合了Andersen初始剪切模型公式[44]的UH模型提高了小应变阶段刚度模拟的准确性,对变形(刚度)敏感的边值问题具有重要意义。
-
为说明UH_G0模型对实际工程问题的适用性,开展了水平加载大直径单桩桩-土相互作用数值模拟试验。同时,为便于比较,使用原始UH模型和Abaqus内嵌的Tresca模型开展平行对比数值试验。
-
有限元分析依托商用有限元软件ABAQUS,有限元模型如图6所示。模型由2部分组成,分别为大直径单桩,桩径D=10 m,桩长L=3D=30 m,共划分为5760个单元,单元类型C3D8;土体,土域直径为6D=60 m,深度为4D=40 m,共划分为21120个单元,单元类型C3D8P。桩身与桩周土体完全绑定。为方便研究取模型的1/2(沿直径所在平面划分)为研究对象。土体的上表面设为自由排水边界,加载过程为完全不排水。桩头荷载包括沿X轴正方向的水平荷载H=20 MN和绕Y轴的顺时针方向弯矩M=600 MN·m,荷载施加在桩顶中心位置,如图6所示。
模拟中为简化模型,将壁厚$t = 100\;{\text{mm}}$,直径$D = 10 \;{\text{m}}$的钢管桩等效为实心桩,依据抗弯刚度(EI)等效计算实心桩的等效弹性模量,计算方法如下式:
$$ \left. \begin{gathered} {E_{{\mathrm{eq}}}} = \frac{{E \cdot ({D^4} - {d^4})}}{{{D^4}}} \\ d = D - 2t \\ \end{gathered} \right\} $$ (17) 式中:
${E_{{\mathrm{eq}}}}$ ——等效弹性模量(MPa);
$E$ ——钢材弹性模量,取$E = 201\;{\text{GPa}}$;
$d$ ——桩的内径(m)。
按式(19)计算可得实心桩弹性模量${E_{{\mathrm{eq}}}} = 16\;{\text{GPa}}$。
土体设置为两层,上层土(层厚1.5D=15 m)为正常固结土,有效重度$\gamma ' = 8\;{\text{kN/}}{{\text{m}}^{\text{3}}}$,渗透系数$k = 5 \times {10^{ - 7}}{\text{m/s}}$;下层土(层厚2.5D=25 m)为超固结土(OCR=4),有效重度$\gamma ' = 10\;{\text{kN/}}{{\text{m}}^{\text{3}}}$,渗透系数与第一层相同。UH模型和UH_G0模型参数见表2。确定模型参数后进行土单元试验计算在不同应力路径下的抗剪强度随深度的变化,模拟结果如图7所示。
从图7可以看出,UH_G0模型预测得到的强度剖面中,单剪模式下的不排水抗剪强度处于中间水平。在Tresca本构的对比性单桩模拟试验中,由于Tresca模型为各向同性本构,不排水抗剪强度剖面取UH_G0模型预测得到的单剪模式抗剪强度$S_{\mathrm{u}}^{{\mathrm{DSS}}}$,且应力-应变曲线为UH_G0预测的单剪模式应力-应变曲线。弹性刚度随深度的变化与UH_G0模型的初始弹性刚度保持一致,Tresca模型具体参数见表3。
表 3 Tresca模型本构参数
Table 3. Constitutive parameters of Tresca model
参数 取值 不排水强度Su 与UH_G0模型的单剪路径不排水强度SuDSS保持一致 弹性刚度E 与UH_G0模型的初始弹性刚度保持一致 泊松比ν 0.495 摩擦角φ/(°) 0 剪涨角δ/(°) 0 -
图8为小荷载水平下桩顶水平位移随桩头荷载的变化,从图中可以发现:UH_G0模型与Tresca模型预测的桩-土相互作用结果吻合良好,较原始UH模型刚度提升较大。当桩头水平荷载达到4 MN时,UH_G0模型和Tresca模型预测桩顶水平位移分别为0.0070 m和0.0069 m,而原始UH模型预测桩顶水平位移达到0.055 m,为UH_G0模型预测数值的7.8倍。
图 8 小荷载水平下桩顶水平位移随桩头荷载的变化
Figure 8. Variation of small-load horizontal displacement with load at pile head
图9为本模拟最大荷载(H=20 MN,M=600 MN·m)加载过程中桩头水平位移随桩头荷载的变化,从图中可以发现:随着荷载水平增大,UH_G0预测的桩-土相互作用刚度较Tresca模型更小,但始终高于原始UH模型预测的桩-土相互作用刚度。当桩头水平荷载达到20 MN时,Tresca模型预测桩头位移为0.51 m,UH_G0模型预测桩头位移为1.15 m,原始UH模型预测桩头位移为1.53 m。
图10为桩顶水平荷载H=20 MN,弯矩M=600 MN·m时的桩身位移曲线以及对应的剪力和弯矩曲线。对比曲线结果可以发现:(1)Tresca模型预测桩身整体转角最小,其次为UH_G0模型,原始UH模型预测桩身整体转动角度最大;此时Tresca模型和UH_G0模型预测的桩头水平位移分别为原始UH模型预测结果的33%和75%;(2)受同样的外荷载时,UH模型和UH_G0模型预测的沿桩长的弯矩与剪力剖面几乎重合,表明两个本构模型所预测的土反力随深度的分布相同,但是所需要的桩身位移差距较大,这说明对UH模型进行的剪切模量修正并未影响其强度预测;(3)Tresca模型截面剪力衰减更快,表明Tresca模型调动的水平土抗力更大。这是由于Tresca模型采用各向同性的抗剪强度(等同于相同深度位置UH_G0模型的单剪模式抗剪强度$S_{\mathrm{u}}^{{\mathrm{DSS}}}$),且在第一层土中,UH模型预测的单剪强度与三轴压缩强度基本一致,但高于三轴拉伸强度(如图7所示),而UH_G0模型考虑土体不排水抗剪强度的各向异性,故UH_G0模型预测得到的土反力应小于同深度处Tresca模型的模拟结果。
-
3.2小节从宏观角度分析了UH_G0本构模型在大直径单桩水平桩-土相互作用问题模拟中得到的桩头和桩身模拟结果,以及这些结果产生的原因。为进一步分析UH_G0本构模型所得结果与原始UH模型模拟结果之间的差异,验证模型修正的效果,以下将追踪大直径单桩有限元模型中不同特征单元在单桩水平加载过程中所调动的应力-应变关系,孔隙水压力与有效应力路径。特征单元位置如图11所示:单元A-D处于正常固结土层(OCR=1),单元E-G处于超固结土层(OCR=4)。表4总结了特征单元在加载结束时的应变情况(采用UH_G0本构):单元A近似处于三轴拉伸状态;单元B近似处于三轴压缩状态;单元C近似处于平面拉伸状态;单元D近似处于平面压缩状态;单元E, F, G近似处于单剪状态。
表 4 特征单元应变(UH_G0模型结果)
Table 4. Strain of characteristic elements (from UH_G0 model results)
单元 深度/m εx εy εz γxy γxz γyz A −4.75 −0.1075 0.0517 0.0474 −0.0266 −0.0081 −0.0117 B −4.75 0.1153 −0.0566 −0.0530 −0.0219 0.0022 −0.0566 C −14.75 −0.0372 0.0116 0.0201 −0.0087 0.0249 −0.0021 D −14.75 0.0498 −0.0149 −0.0547 −0.0087 0.0335 0.0030 E −31.25 0.0391 0.0013 −0.0321 0.0057 −0.1613 0.0015 F −31.25 −0.0348 −0.0015 0.0554 0.0047 −0.1899 −0.0008 G −31.25 −0.0067 0.0000 0.0029 0.0055 −0.1595 0.0001 图12为边值问题模拟中各特征单元的应力-应变图、孔压-应变图和有效应力路径图。图中平均主应力$p = \left( {{\sigma _x} + {\sigma _y} + {\sigma _z}} \right)/3$,剪应力$q = {\sigma _z} - {\sigma _x}$。
图 12 单桩水平加载过程中特征单元调动的应力-应变、孔压-应变响应及应力路径
Figure 12. Stress-strain, pore pressure-strain response and stress path of characteristic elements during horizontal loading of monopile
对比原始UH模型和UH_G0模型的模拟结果可以发现:对于不同超固结度,不同应力路径的特征单元,两个模型给出的应力路径模拟契合度很高。孔压和抗剪强度的终值很接近,区别是达到终值时的应变不同以及模拟过程中的刚度不同,这与单桩桩-土相互作用宏观上的响应是一致的。
Application of Unified Hardening Model in Soil-Pile Analysis of Offshore Wind Power Large-Diameter Monopiles
-
摘要:
目的 文章旨在推动海上风电基础的优化设计,将理论更加完善,模拟结果更加符合实际的高级岩土本构模型推广至实际工程应用。 方法 首先介绍了一种高级本构模型:UH模型(Unified Hardening model,统一硬化本构模型)的理论框架,随后对其进行了面向实际工程应用需求的修正,最后将修正后的模型应用在大直径单桩桩-土分析中,并对结果进行了宏观和单元层面的分析以验证其有效性与实用性。 结果 经过修正和验证,得到了适用于实际工程问题的UH_G0模型。该模型将原始UH模型中的参数“初始孔隙比”替换为“超固结比”;并引入Andersen经验公式,将剪切模量与本构参数κ解耦,显著提升了原始UH模型的初始剪切刚度,使得UH模型对变形(刚度)敏感的边值问题具有适用性。 结论 UH模型作为理论框架清晰,参数物理意义明确的高级本构模型具有较好的实际工程应用潜力。研究进行的面向工程设计需求的两项修正是必要的,二者分别提高了UH模型实际应用于边值问题的便捷性和初始刚度计算的准确性。修正后的UH_G0模型在有限元模拟过程中表现出对实际工程问题良好的有效性与实用性。 Abstract:Introduction In order to facilitate the design optimisation of foundation for offshore wind turbines, it is necessary to promote application of advanced geotechnical constitutive models that are theoretically sound and capable of modelling realistic soil behaviours in real engineering scenarios. Method Firstly, the theoretical framework of the Unified Hardening model (UH model), an advanced constitutive model, was introduced. Then, further modifications oriented to requirements in real engineering applications were presented. Finally, the modified model was applied to the pile-soil analysis of large-diameter monopoles, and the results were analyzed at the macroscopic and element level to verify its effectiveness and practicability. Result The UH_G0 model obtained through modification and verification is applicable to addressing issues in real engineering practices. The parameter "initial void ratio" in the original UH model is replaced by "overconsolidation ratio". By incorporating the Andersen empirical formula, the shear modulus is decoupled from the constitutive parameter κ, significantly improving the initial shear stiffness of the original UH model, and thus allowing the UH model to be applicable to the deformation (stiffness) sensitive boundary value problems. Conclusion As an advanced constitutive model with a well-defined theoretical framework and physically meaningful parameters, the UH model demonstrates high potential for practical engineering applications. The two modifications tailored to engineering design requirements are necessary, which improve the convenience in application of the UH model to address boundary value problems and the accuracy in calculating initial stiffness, respectively. The modified UH_G0 model is proven to be effective and practical for solving practical engineering problems in the process of finite element simulation. -
表 1 UH本构参数物理意义
Tab. 1. Physical meanings of UH constitutive parameters
表示符号 物理意义 M 极限状态应力比(三轴压缩应力路径) ν 泊松比 κ 等向回弹线斜率 λ 等向压缩线斜率 e0 初始孔隙比 N e-lnp空间正常压缩线(NCL)截距 参数 取值 极限状态应力比M 0.87 泊松比ν 0.2 等向回弹线斜率κ 0.035 等向压缩线斜率λ 0.093 e-lnp空间正常压缩线截距N 1.06 超固结比OCR 见图中标注 塑性指数Ip 25 表 3 Tresca模型本构参数
Tab. 3. Constitutive parameters of Tresca model
参数 取值 不排水强度Su 与UH_G0模型的单剪路径不排水强度SuDSS保持一致 弹性刚度E 与UH_G0模型的初始弹性刚度保持一致 泊松比ν 0.495 摩擦角φ/(°) 0 剪涨角δ/(°) 0 表 4 特征单元应变(UH_G0模型结果)
Tab. 4. Strain of characteristic elements (from UH_G0 model results)
单元 深度/m εx εy εz γxy γxz γyz A −4.75 −0.1075 0.0517 0.0474 −0.0266 −0.0081 −0.0117 B −4.75 0.1153 −0.0566 −0.0530 −0.0219 0.0022 −0.0566 C −14.75 −0.0372 0.0116 0.0201 −0.0087 0.0249 −0.0021 D −14.75 0.0498 −0.0149 −0.0547 −0.0087 0.0335 0.0030 E −31.25 0.0391 0.0013 −0.0321 0.0057 −0.1613 0.0015 F −31.25 −0.0348 −0.0015 0.0554 0.0047 −0.1899 −0.0008 G −31.25 −0.0067 0.0000 0.0029 0.0055 −0.1595 0.0001 -
[1] 国家发展改革委, 国家能源局. “十四五”现代能源体系规划 [EB/OL]. (2022-03-22) [2023-01-26]. https://www.ndrc.gov.cn/xwdt/tzgg/202203/t20220322_1320017.html. National Development and Reform Commission, National Energy Administration. "Fourteenth Five-Year Plan" modern energy system planning [EB/OL]. (2022-03-22) [2023-01-26]. https://www.ndrc.gov.cn/xwdt/tzgg/202203/t20220322_1320017.html. [2] 元国凯, 方辉, 马兆荣,等. 船撞海上风电单桩基础结构剩余强度研究 [J]. 内蒙古电力技术, 2022, 40(5): 7-15. DOI: 10.19929/j.cnki.nmgdljs.2022.0074. YUAN G K, FANG H, MA Z R, et al. Study on residual strength of offshore wind turbine monopile foundations after ship collision [J]. Inner Mongolia electric power, 2022, 40(5): 7-15. DOI: 10.19929/j.cnki.nmgdljs.2022.0074. [3] 曾雨欣, 施伟, 张礼贤, 等. 10 MW大型单桩式海上风机桩土作用研究 [J]. 南方能源建设, 2023, 10(1): 1-12. DOI: 10.16516/j.gedi.issn2095-8676.2023.01.001. ZENG Y X, SHI W, ZHANG L X, et al. Research on pile-soil interaction of 10 MW large monopile offshore wind turbine [J]. Southern energy construction, 2023, 10(1): 1-12. DOI: 10.16516/j.gedi.issn2095-8676.2023.01.001. [4] 刘晋超. 海上大直径单桩基础沉桩施工关键技术研究 [J]. 南方能源建设, 2022, 9(1): 47-51. DOI: 10.16516/j.gedi.issn2095-8676.2022.01.007. LIU J C. Research on key technologies of pile driving construction for monopile [J]. Southern energy construction, 2022, 9(1): 47-51. DOI: 10.16516/j.gedi.issn2095-8676.2022.01.007. [5] ZDRAVKOVIĆ L, TABORDA D M G, POTTS D M, et al. Numerical modelling of large diameter piles under lateral loading for offshore wind applications [C]//3rd International Symposium on Frontiers in Offshore Geotechnics, Oslo, Norway, June 10-12, 2015. Boca Raton: CRC Press, 2015: 759-764. [6] BURD H J, TABORDA D M G, ZDRAVKOVIĆ L, et al. PISA design model for monopiles for offshore wind turbines: application to a marine sand [J]. Géotechnique, 2020, 70(11): 1048-1066. DOI: 10.1680/jgeot.18.p.277. [7] BYRNE B W, HOULSBY G T, BURD H J, et al. PISA design model for monopiles for offshore wind turbines: application to a stiff glacial clay till [J]. Géotechnique, 2020, 70(11): 1030-1047. DOI: 10.1680/jgeot.18.P.255. [8] TERZAGHI K. Theoretical soil mechanics [M]. New York: John Wiley & Sons, Inc. , 1943. [9] DRUCKER D C, PRAGER W. Soil mechanics and plastic analysis or limit design [J]. Quarterly of applied mathematics, 1952, 10(2): 157-165. DOI: 10.1090/qam/48291. [10] ROSCOE K H, SCHOFIELD A N, THURAIRAJAH A. Yielding of clays in states wetter than critical [J]. Géotechnique, 1963, 13(3): 211-240. DOI: 10.1680/geot.1963.13.3.211. [11] ROSCOE K H, BURLAND J B. On the generalized stress-strain behavior of ‘wet’ clays [M]. Cambridge: Cambridge University Press, 1968. [12] LADE P V, DUNCAN J M. Elastoplastic stress-strain theory for cohesionless soil [J]. Journal of the geotechnical engineering division, 1975, 101(10): 1037-1053. DOI: 10.1061/AJGEB6.0000204. [13] DAFALIAS Y F, POPOV E P. Plastic internal variables formalism of cyclic plasticity [J]. Journal of applied mechanics, 1976, 43(4): 645-651. DOI: 10.1115/1.3423948. [14] DAFALIAS Y F. Bounding surface plasticity. I: mathematical foundation and hypoplasticity [J]. Journal of engineering mechanics, 1986, 112(9): 966-987. DOI: 10.1061/(ASCE)0733-9399(1986)112:9(966). [15] DAFALIAS Y F. Plastic spin: necessity or redundancy? [J]. International journal of plasticity, 1998, 14(9): 909-931. DOI: 10.1016/S0749-6419(98)00036-9. [16] DAFALIAS Y F, MANZARI M T, PAPADIMITRIOU A G. SANICLAY: simple anisotropic clay plasticity model [J]. International journal for numerical and analytical methods in geomechanics, 2006, 30(12): 1231-1257. DOI: 10.1002/nag.524. [17] DAFALIAS Y F, POPOV E P. A model of nonlinearly hardening materials for complex loading [J]. Acta mechanica, 1975, 21(3): 173-192. DOI: 10.1007/BF01181053. [18] DAFALIAS Y F, MANZARI M T, AKAISHI M. A simple anisotropic clay plasticity model [J]. Mechanics research communications, 2002, 29(4): 241-245. DOI: 10.1016/S0093-6413(02)00252-5. [19] HASHIGUCHI K, CHEN Z P. Elastoplastic constitutive equation of soils with the subloading surface and the rotational hardening [J]. International journal for numerical and analysis methods in geomechanics, 1998, 22(3): 197-227. DOI: 10.1002/(SICI)1096-9853(199803)22:33.0.CO;2-T. [20] HASHIGUCHI K. Generalized plastic flow rule [J]. International journal of plasticity, 2005, 21(2): 321-351. DOI: 10.1016/j.ijplas.2003.12.003. [21] HASHIGUCHI K, OZAKI S, OKAYASU T. Unconventional friction theory based on the subloading surface concept [J]. International journal of solids and structures, 2005, 42(5/6): 1705-1727. DOI: 10.1016/j.ijsolstr.2004.08.006. [22] ASAOKA A, NODA T, FERNANDO G S K. Consolidation deformation behavior of lightly and heavily overconsolidated clay foundations [J]. Soils and foundations, 1998, 38(2): 75-91. DOI: 10.3208/sandf.38.2_75. [23] ASAOKA A, NAKANO M, NODA T. Superloading yield surface concept for highly structured soil behavior [J]. Soils and foundations, 2000, 40(2): 99-110. DOI: 10.3208/sandf.40.2_99. [24] ASAOKA A, NODA T, YAMADA E, et al. An Elasto-plastic description of two distinct volume change mechanisms of soils [J]. Soils and foundations, 2002, 42(5): 47-57. DOI: 10.3208/sandf.42.5_47. [25] 黄文熙. 土的弹塑性应力-应变模型理论 [J]. 岩土力学, 1979(1): 1-20. DOI: 10.16285/j.rsm.1979.01.002. HUANG W X. Theory of elastoplastic stress-strain model for soil [J]. Rock and soil mechanics, 1979(1): 1-20. DOI: 10.16285/j.rsm.1979.01.002. [26] 黄文熙. 硬化规律对土的弹塑性应力-应变模型影响的研究 [J]. 岩土工程学报, 1980, 2(1): 1-11. DOI: 10.3321/j.issn:1000-4548.1980.01.001. HUANG W X. The influence of the hardening law on the formulation of the Elasto-plastic model of soil [J]. Chinese journal of geotechnical engineering, 1980, 2(1): 1-11. DOI: 10.3321/j.issn:1000-4548.1980.01.001. [27] 黄文熙, 濮家骝, 陈愈炯. 土的硬化规律和屈服函数 [J]. 岩土工程学报, 1981, 3(3): 19-26. DOI: 10.3321/j.issn:1000-4548.1981.03.003. HUANG W X, PU J L, CHEN Y J. Hardening rule and yield function for soils [J]. Chinese journal of geotechnical engineering, 1981, 3(3): 19-26. DOI: 10.3321/j.issn:1000-4548.1981.03.003. [28] 沈珠江. 土的弹塑性应力应变关系的合理形式 [J]. 岩土工程学报, 1980, 2(2): 11-19. DOI: 10.3321/j.issn:1000-4548.1980.02.002. SHEN Z J. The rational form of stress-strain relationship of soils based on Elasto-plasticity theory [J]. Chinese journal of geotechnical engineering, 1980, 2(2): 11-19. DOI: 10.3321/j.issn:1000-4548.1980.02.002. [29] 沈珠江. 土的三重屈服面应力应变模式 [J]. 固体力学学报, 1984(2): 163-174. DOI: 10.19636/j.cnki.cjsm42-1250/o3.1984.02.002. SHEN Z J. A stress-strain model for soils with three yield surfaces [J]. Acta mechanica solida sinica, 1984(2): 163-174. DOI: 10.19636/j.cnki.cjsm42-1250/o3.1984.02.002. [30] 沈珠江. 几种屈服函数的比较 [J]. 岩土力学, 1993, 14(1): 41-50. DOI: 10.16285/j.rsm.1993.01.005. SHEN Z J. Comparison of several yield function [J]. Rock and soil mechanics, 1993, 14(1): 41-50. DOI: 10.16285/j.rsm.1993.01.005. [31] 沈珠江. 三种硬化理论的比较 [J]. 岩土力学, 1994, 15(2): 13-19. DOI: 10.16285/j.rsm.1994.02.002. SHEN Z J. Comparison among three kinds of hardening theories [J]. Rock and soil mechanics, 1994, 15(2): 13-19. DOI: 10.16285/j.rsm.1994.02.002. [32] 沈珠江. 粘土的双硬化模型 [J]. 岩土力学, 1995, 16(1): 1-8. DOI: 10.16285/j.rsm.1995.01.001. SHEN Z J. A double hardening model for clays [J]. Rock and soil mechanics, 1995, 16(1): 1-8. DOI: 10.16285/j.rsm.1995.01.001. [33] 殷宗泽, 邓肯J M. 剪胀土与非剪胀土的应力应变关系 [J]. 岩土工程学报, 1984, 6(4): 24-40. DOI: 10.3321/j.issn:1000-4548.1984.04.003. YIN Z Z, DUNCAN J M. A stress-strain relationship for dilative and non-dilative soils [J]. Chinese journal of geotechnical engineering, 1984, 6(4): 24-40. DOI: 10.3321/j.issn:1000-4548.1984.04.003. [34] 殷宗泽. 一个土体的双屈服面应力-应变模型 [J]. 岩土工程学报, 1988, 10(4): 64-71. DOI: 10.3321/j.issn:1000-4548.1988.04.007. YIN Z Z. A double yield surface stress-strain model of soil [J]. Chinese journal of geotechnical engineering, 1988, 10(4): 64-71. DOI: 10.3321/j.issn:1000-4548.1988.04.007. [35] 殷宗泽. 土体本构模型剖析 [J]. 岩土工程学报, 1996, 18(4): 95-97. DOI: 10.3321/j.issn:1000-4548.1996.04.016. YIN Z Z. Analysis of soil constitutive model [J]. Chinese journal of geotechnical engineering, 1996, 18(4): 95-97. DOI: 10.3321/j.issn:1000-4548.1996.04.016. [36] YAO Y P, LU D C, ZHOU A N, et al. Generalized non-linear strength theory and transformed stress space [J]. Science in China series E:technological sciences, 2004, 47(6): 691-709. DOI: 10.1360/04ye0199. [37] YAO Y P, ZHOU A N, LU D C. Extended transformed stress space for geomaterials and its application [J]. Journal of engineering mechanics, 2007, 133(10): 1115-1123. DOI: 10.1061/ASCE0733-93992007133:101115. [38] 侯伟. 超固结土的本构模型 [D]. 北京: 北京航空航天大学, 2007. DOI: 10.26937/d.cnki.gbhtu.2007.000001. HOU W. A consititutive model for overconsolidated clays [D]. Beijing: Beihang University, 2007. DOI: 10.26937/d.cnki.gbhtu.2007.000001. [39] YAO Y P, SUN D A, MATSUOKA H. A unified constitutive model for both clay and sand with hardening parameter independent on stress path [J]. Computers and geotechnics, 2008, 35(2): 210-222. DOI: 10.1016/j.compgeo.2007.04.003. [40] YAO Y P, HOU W, ZHOU A N. UH model: three-dimensional unified hardening model for overconsolidated clays [J]. Géotechnique, 2009, 59(5): 451-469. DOI: 10.1680/geot.2007.00029. [41] YAO Y P, GAO Z W, ZHAO J D, et al. Modified UH model: constitutive modeling of overconsolidated clays based on a parabolic Hvorslev envelope [J]. Journal of geotechnical and geoenvironmental engineering, 2012, 138(7): 860-868. DOI: 10.1061/(ASCE)GT.1943-5606.0000649. [42] 姚仰平, 罗汀, 侯伟. 土的本构关系 [M]. 北京: 人民交通出版社, 1991. YAO Y P, LUO T, HOU W. Soil constitutive models [M]. Beijing: China Communications Press Co. , Ltd. , 1991. [43] MAYNE P W, KULHAWY F H. K-OCR relationships in soil [J]. International journal of rock mechanics and mining sciences & geomechanics abstracts, 1983, 20(1): A2. DOI: 10.1016/0148-9062(83)91623-6. [44] ANDERSEN K H. Cyclic soil parameters for offshore foundation design [C]//The 3rd McClelland Lecture-Frontiers in Offshore Geotechnics Ⅲ, London, England, May 15, 2015. London: Taylor & Francis Group, 2015: 5-82.