2022, Vol. 54引用本文 |
2. 北京工业大学 建筑工程学院，北京 100124;
3. 中国电力工程顾问集团东北电力设计研究院有限公司，吉林 长春 130021;
4. 上海海事大学 海洋科学与工程学院，上海 201306;
5. 河海大学 岩土力学与堤坝工程教育部重点实验室，江苏 南京 210098
2. College of Architecture and Civil Eng., Beijing Univ. of Technol., Beijing 100124, China;
3. Northeast Electric Power Design Institute Co., Ltd. of China Power Eng. Consulting, Changchun 130021, China;
4. College of Ocean and Eng., Shanghai Maritime Univ., Shanghai 201306, China;
5. Key Lab. of Geomechanics and Embankment Eng. (Hohai Univ.), Ministry of Edu., Nanjing 210098, China
螺旋锚基础具有结构简单、施工方便、环境影响小等优点，在建筑、输电线路、海洋等工程中得到广泛的应用，其抗拔性能与平板锚相似，抗拔承载力预测是设计的关键。Ghaly[1]和Hao[2]等缩比尺模型试验和Hao等[3]离心机试验结果表明，螺旋锚锚盘型式对承载力的影响很小，预埋螺旋锚与平板锚的上拔承载力基本相同。锚板上拔承载力计算主要有极限平衡法[3-4]、极限分析法[5-6]、有限元法[7-8]、孔扩张分析法[9-10]等，其中极限平衡法的应用最为广泛，该分析方法将滑裂面内土体视为刚体，无需考虑土体应力–变形关系，同时对滑裂面形式和滑裂面上的应力分布及大小进行了假设。
一些学者通过室内模型试验、离心机试验、现场试验及数值模拟的方法对锚板上拔破坏滑裂面提出不同的简化，包含圆柱面[11]，倒锥台[12-14]及曲面[15-17]。对于沿滑裂面上的应力分布及大小，有研究者提出应力沿锚埋深呈线性分布，大小为侧压力系数Ku与自重应力的乘积，且Ku为常值，如Meyerhof等[13]基于Adams等[18]的试验结果，采用Caquot[19]曲面滑裂面时的土压力系数Ku，经过形状修正系数s计算圆锚的侧压力系数为sKu；Mitsch等[14]所用Ku值与Meryerhof等[13]相似，但考虑安装扰动的影响，Ku值减少了30%～40%；Murray等[16]假设Ku=K0=1–sinφ，φ为内摩擦角；Ghaly等[1]采用被动土压力系数
| $ {K_{\text{n}}} = 1 - \dfrac{{\sin \;{\phi _{{\text{cr}}}}\left( {1 + \cos 2\psi } \right)}}{2} $ | (1) |
式中，φcr为临界状态内摩擦角。
大多数极限平衡法包含的经验参数由缩比尺模型试验结果确定，在预测实际较大尺寸锚板上拔承载力时通常不准确，存在高估风险。因此，Cerfontaine等[8]基于2维有限元软件PLAXIS，采用小应变硬化土本构模型(HS small)模拟大直径锚板上拔过程，探讨浅埋圆锚或螺旋锚抗拔机理及半解析解的可靠性。但对于具有应变软化的砂土，HS模型不能捕捉大应变时峰值后出现的软化现象，数值模拟会高估锚板的上拔承载力，尤其对于较大埋深比情况。为考虑砂土应变软化特性，陈榕等[21]基于Abaqus平台采用改进的Mohr–Coulomb弹塑性本构模型，模拟密砂中圆锚上拔行为，揭示砂土渐进破坏过程，并分析尺寸效应产生的原因，但未对滑裂面上应力及上拔承载力计算理论等进行分析。
为更真实地模拟砂土压硬性、剪胀性及应变软化特性对锚板上拔破坏机理的影响，本文以Abaqus有限元软件为平台，基于陈榕等[21]改进的Mohr–Coulomb本构模型及强度参数对应力的依赖关系，采用欧拉拉格朗日耦合法(CEL法)模拟不同密实度砂土中不同埋深比圆或螺旋形锚板上拔过程，探讨浅埋锚上拔破坏机理，分析破坏滑裂面形式及滑裂面上法向应力及土体强度参数的分布规律，提出了更准确地预测圆锚或螺旋锚上拔承载力的方法，同时对比分析已有的理论计算公式，并建议各公式的适用条件。
1 有限元模型及验证 1.1 有限元模型为分析不同密实度砂土中圆锚或螺旋锚上拔承载机理，建立松、中密和极密砂中埋深比H:D=1.0～7.5的圆锚拉拔试验有限元模型。
1.1.1 计算模型圆锚上拔承载问题为轴对称问题，故在Abaqus/explicit分析模块中建立四分之一锚–土耦合体系的三维模型。将土体设置为欧拉体，锚设置为拉格朗日体。为与离心机试验结果[3]对比，锚盘直径D取为400 mm，螺距为0.25D，盘厚为0.05D，锚杆轴径为0.235D。考虑锚体刚度远大于土体，将锚设为刚体，忽略其变形，离散化为拉格朗日单元C3D8R，在锚顶部设置参考点，约束参考点的水平向位移和轴向转动。根据文献[22]，计算域从锚杆轴线径向取10D；锚盘底部向下取锚埋深H；为允许上拔过程中土体隆起和流入空的欧拉单元，地面以上竖向取5D，设置为零强度和刚度的空单元。计算域底部边界约束各方向位移；侧面边界约束速度和转角。计算域内土体网格规则划分，采用EC3D8R单元，锚板上下最大5D、水平向4D范围内网格密化，最小网格宽度ΔB在板边缘处，计算几何模型如图1所示。锚–土接触面采用自动识别，法向作用设为硬接触，切向作用采用罚函数法，摩擦系数设为0.3。在锚杆顶部参考点施加上拔位移1D，由于采用动态分析模块进行准静态分析，设置的加载时间应尽可能接近准静态；同时，为获得更稳定的结果，动力显示分析步中时间比例因子取0.6[23]。
 |
| 图1 数值计算模型 Fig. 1 Numerical model |
1.1.2 土体模型及参数
土体本构模型采用改进的Mohr–Coulomb弹塑性模型[21]，模型参数包括变形模量E、泊松比v、初始侧压力系数K0、等效塑性应变阀值
| 表1 土性参数 Tab. 1 Soil parameters |
 |
1.2 模型及计算结果验证 1.2.1 加载速率及网格密度的影响
考虑计算的稳定性及时效性，在批量计算之前，进行加载速率和网格密度影响的变参数计算。额外选择极密砂中更深埋（H:D=9）的情况进行对比计算，其中，土性参数φp=42.5°，ψp=23°，E=53.3 MPa。为分析加载速率的影响，取加载速率V=0.025D/s、0.050D/s、0.100D/s，以最小网格尺寸比ΔB:D=0.1进行计算，不同加载速率的位移–荷载半对数坐标曲线如图2（a）所示。由图2（a）可见，加载速率对上拔承载力的影响较小。当计算至0.5D上拔位移时，速率V1和V2计算时间约为V3的12和2倍，从耗时和稳定性考虑选择V2加载速率进行后续计算。为分析网格密度的影响，选择3种密度网格进行计算，最小网格尺寸比ΔB:D分别为0.05，0.10和0.20，对应划分的网格数量分别为613 130、280 864、125 856，不同网格密度的计算结果如图2（b）所示。由图2（b）可见，网格越疏，计算上拔承载力越大，ΔB:D=0.05和0.10两种情况的计算结果接近，但0.05网格密度比的计算时间为0.10网格密度比的6.6倍，因此，后续计算中采用最小网格尺寸比ΔB:D=0.10的网格划分方法。
 |
| 图2 加载速率和网格密度影响 Fig. 2 Effects of loading rate and mesh density |
1.2.2 计算结果
图3为不同密实度砂土中不同埋深比锚板的上拔位移–荷载曲线，以半对数坐标形式绘制。埋深较浅时，位移–荷载曲线有明显峰值出现，峰值之后，承载力下降；随着埋深比增加，曲线特征发生变化，可分为初始直线段、弯曲段和后期平稳段（近似直线）3个阶段。锚板极限上拔承载力取曲线的峰值或曲线平稳段起点，图3中以空心圆表示，同时将极限上拔荷载时对应的破坏位移up、净极限上拔承载力Qu（扣除锚自重）及上拔力系数Nγ=Qu/γAH列于表2，其中，γ为土体重度，A为锚盘面积，A=πD2/4。
 |
| 图3 上拔位移–荷载曲线 Fig. 3 Curves of uplift resistance and displacement |
| 表2 数值计算结果 Tab. 2 Numerical results for different conditions |
 |
1.2.3 结果验证
为验证模型的准确性，将数值计算（FEM）结果与Hao等[3]密砂中螺旋锚离心机试验结果进行比较，如图4所示。
 |
| 图4 数值结果与离心机试验对比 Fig. 4 Comparison between numerical and centrifugal tests results |
由图4可见，随着埋深比增加，数值计算结果获得的Nγ与离心机试验结果趋势相同；离心机试验的砂土密实度在85.8%～96.5%之间，其试验结果介于中密与极密砂的数值计算结果之间，验证了计算模型的可靠性。但从严格角度，数值模拟结果较离心机试验结果偏高，可能是由于按照Bolton关系式[25]确定的剪胀角较实际情况偏大引起。
2 破坏模式 2.1 埋置深度及密实度对破坏模式的影响图5为不同密实度土中不同埋深比锚板上拔破坏时对应的土体等效塑性应变云图，其中，应变超过2%区域以黑色显示，红线为郝冬雪[2-3]及Giampa[4]等建议的浅埋时理论破坏面，即与竖轴夹角为剪胀角的倒锥台面。
 |
| 图5 破坏时等效塑性应变云图 Fig. 5 Contours of equivalent plastic strain at failures |
浅埋时，塑性区从锚盘边缘开始，向两侧发展，直至地表，形成贯通的剪切破坏带，而锚板上方有较大的弹性区，此弹性区会随锚埋深增加逐渐变小；除H:D=1外，上拔位移至一定值后，地表处杆周出现塑性区，并沿杆向下发展，锚埋深越大，杆周塑性区向下发展得越深，表明在一定范围内杆侧摩阻得到发挥，但由于杆轴径较小，侧摩阻力的影响很小。根据文献[26-27]，由最大剪应变值近似确定破坏面，则浅埋时，数值计算获得的滑裂面可视为直线，这种破坏模式为浅埋破坏模式；对于松砂和中密砂，数值滑裂破坏面与郝冬雪[2-3]及Giampa[4]等理论破坏面（红线）基本一致，而极密砂土中数值破坏面（绿虚线）与理论破坏面偏差稍大，H:D=2～6时，理论破坏面较数值破坏面外扩约5°。
随着埋深增加，土体塑性区向锚盘外侧发展受限制，开始向锚盘正上方发展，其上方弹性区逐渐消失，数值破坏面与理论直线破坏面偏差增大，破坏面由斜面形式向朝杆轴弯曲的曲面发展，但依然会持续开展至地表，此破坏模式为过渡破坏模式，如松砂H:D=4～5，中密砂H:D=6和极密砂H:D=7.5。
当埋深较大时，由锚盘边缘开始的塑性区不会连续发展至地表，而是局限在地表以下，呈梨状封闭面，出现了深埋破坏模式，如松砂H:D=6.0、7.5，塑性区开展高度为5D～6D。锚盘以上一定高度范围内应变等值线呈气泡状（深色区域），超过此高度，应变等值线近似呈直线向轴线收缩。
在不同密实度砂土中，随着锚埋深的增加，其破坏模式均会出现上述的浅、过渡及深破坏模式，但其出现过渡的埋深比会受土体密实度影响，土体越密实，过渡破坏模式出现的埋深越大。根据图5，松砂、中密砂和极密砂3种情况过渡破坏模式出现的埋深比分别为4.0、6.0和7.5。
2.2 浅埋滑裂面假设数值计算结果表明，砂土中浅埋锚的破坏滑裂面基本为倒锥台型，倾角可近似为剪胀角。过渡埋深破坏面虽然在近地表时出现收缩，但按浅埋破坏模式开展至地表时，破坏面与地表的交点仍处于塑性范围，为简化后续计算，将过渡破坏模式按照浅埋破坏模式处理。因此，在理论计算时假设浅及过渡埋深时的滑裂面均为倾斜角度等于砂土剪胀角的倒锥台面。
3 滑裂面上内摩擦角及应力变化规律本文采用的模型较理想弹塑性本构模型更能够真实地反映锚上拔破坏过程[21]，破坏滑裂面上土体处于不同的应变状态，故沿滑裂面的应力比及土体强度不同，即滑动面上土体侧压力系数并非常值，内摩擦角亦非峰值内摩擦角。
3.1 滑裂面上内摩擦角变化图6为自锚盘边缘开始沿破坏滑裂面上各单元的内摩擦角，其中z:D为各单元中点至锚盘的垂直距离与盘径之比。
 |
| 图6 破坏时滑裂面上内摩擦角分布 Fig. 6 Distribution of friction angles along failure surface |
由图6可见：从锚盘边缘开始内摩擦角逐渐增大至峰值内摩擦角，之后又开始减小，即滑裂面上土体强度发挥的程度不同；埋深越大，峰值内摩擦角的位置离锚盘越远。内摩擦角的发展规律与等效塑性应变相关，图5中盘边缘两侧剪切带一定高度范围的塑性应变超过2%，表明土体的强度经历了初始值至峰值的过程，进入强度降低的某个阶段，该区域内摩擦角小于峰值内摩擦角；对于较大埋深的锚，破坏时盘边缘等效塑性应变超过20%，内摩擦角进入临界状态内摩擦角，如松砂H:D≥5.0，中密砂和极密砂H:D=7.5。
松砂在H:D≥5.0，中密砂H:D=7.5和极密砂H:D≥2时，在接近地表处滑裂面上的内摩擦出现异常增大，这是因为它们的理论滑裂面与数值滑裂面在接近地表处的偏差所致；偏差越大，内摩擦角异常增大越明显，如松砂H:D≥5.0，中密砂H:D≥6.0及极密砂H:D=7.5的过渡和深埋破坏情况。但对于过渡破坏的埋深，仅在靠近地表的很小范围内出现内摩擦角增大，且增量不超过3°，故对于过渡破坏模式，可采用浅破坏模式假设。
3.2 滑裂面上应力分布规律提取沿滑裂面上各单元积分点的大小主应力，并根据各单元的内摩擦角值按Mohr–Coulomb平衡条件计算滑裂面上的法向应力σn及土压力系数Ku，Ku=σn/γ(H–z)，z为计算点至锚盘距离。将Ku以Hao[3]等提出的土压力系数Kn进行标准化，获得了埋深范围内的标准化土压力系数Ku:Kn，如图7所示。自锚盘边缘起，侧压力系数Ku先近似呈线性增加至峰值Ku,peak，而后非线性快速衰减至接近Kn值，这与Cerfontaine等[8]基于硬土模型的数值模拟结果及Schiavon等[28]利用光弹法测量透明土中上拔螺旋桩的应力规律相同，但峰值及峰值点的位置有所不同。
 |
| 图7 滑裂面上标准化侧压力系数分布 Fig. 7 Distribution of normalized coefficients of lateral pressure along failure surface |
图8为不同密实度土中不同埋深比圆锚破坏滑裂面上的峰值土压力系数Ku,peak，图9为不同埋深比时峰值点至锚盘距离与盘径比zp:D及其与Cerfontaine等[8]结果的对比。由图8可见，埋深比相对于土体密实度对Ku,peak有更重要的影响；在埋深比较小时，Ku,peak基本保持不变，之后，随着埋深比增加基本呈现线性增加；埋深比相同时，砂土密实度越小，土压力系数峰值Ku,peak越大，这是由于土体越松，压缩性越高，拉拔过程中挤密程度越大，从而获得的土压力系数越大。与Cerfontaine等[8]结果的对比发现，密实度相近的砂土，基于硬土模型的有限元分析获得的Ku,peak基本上比本文采用软化模型获得的计算结果高，并且采用硬土模型的结果表现出土体密实度对Ku,peak的影响比本文获得的密实度影响更大。
 |
| 图8 侧压力系数峰值Ku,peak对比 Fig. 8 Comparison of the peak lateral coefficients Ku,peak |
 |
|
图9 峰值侧压力系数对应的位置
|
由图9可见：Cerfontaine等[8]获得的规律为，zp:D明显受埋深比和相对密实度的影响，密实度越小，zp:D越小；在较小埋深时，随着埋深比H:D增加，zp:D线性增加至某一限值，而后保持不变，但zp:D达到限值对应的埋深比与临界埋深比无关。而本文获得的zp:D与埋深比H:D及密实度Dr的规律为，锚浅埋及过渡埋深时，峰值土压力系数位置zp:D基本不受H:D和Dr的影响，即极密砂中H:D≤7.5，中密砂H:D≤6.0，松砂H:D≤5.0，zp:D≈0.2；对于松砂深埋H:D=6.0、7.5及中密砂深埋H:D=7.5情况，zp:D≈0.33。
4 理论计算 4.1 基于数值结果的极限上拔承载力极限上拔承载力Qu为滑裂面上法向应力σn及剪切强度
| $ \begin{aligned}[b] &{Q_{\tau 1}} = \int_0^H {\left( {{\tau _{\rm{f}}}\cos \;{\psi _{\rm{p}}} - {\sigma _{\rm{n}}}\sin \;{\psi _{\rm{p}}}} \right) \cdot \text{π}\left( {D + 2\textit{z}\tan \;{\phi _{\rm{p}}}} \right)\frac{{{\rm{d}}\textit{z}}}{{\cos \;{\psi _{\rm{p}}}}}} {\rm{ = }}\\ &\;\;\;\;\;\;\;\;\;\; \displaystyle\sum\limits_{i = 1}^n \{ \text{π} \gamma '\;{K_{{\rm{u}}i}}\left( {\tan\; {\phi _i} - \tan\; {\psi _{\rm{p}}}} \right)\left( {H - {\textit{z}_i}} \right)\left( {{\textit{z}_{i{\rm{ + }}1}} - {\textit{z}_i}} \right)\cdot \\&\;\;\;\;\;\;\;\;\;\; \left[ {D + \left( {{\textit{z}_i} + {\textit{z}_{i{\rm{ + }}1}}} \right)\tan \;{\psi _{\rm{p}}}} \right] \}\\[-10pt] \end{aligned} $ | (2) |
| $ W=\text{π}\gamma 'H\left( {\frac{{{D^2}}}{4} + \frac{{HD\tan\; {\psi _{\text{p}}}}}{2} + \frac{{{H^2}{{\tan }^2}\;{\psi _{\text{p}}}}}{3}} \right) $ | (3) |
式中，γ′为土体有效重度，H为锚埋深，D为锚盘直径，ψp为峰值剪胀角，Kui为单元i上土侧压力系数，zi+1和zi分别为i单元上、下边界与锚盘的垂直距离，φi为i单元的内摩擦角。
将不同密实度砂土中埋深比H:D≤7.5锚滑裂面上抗力Qτ1及滑裂面内土重W计算结果列于表3。
| 表3 基于FEM结果的极限上拔力理论值和FEM结果对比 Tab. 3 Comparison between theoretical ultimate uplift resistance based on FEM and numerical results |
 |
由表3可见：对于浅埋情况，松及中密砂，采用由假设滑裂面上提取的应力和内摩擦角数值结果计算的极限上拔承载力Qu1与数值模拟获得的QuFEM误差在–10%～7%之间，表明浅埋假设的破坏滑裂面与实际情况基本一致；对于Dr=100%的情况（H:D≤7.5），Qu1与QuFEM误差在–0.3%～15%之间，除H:D=1外，Qu1均大于QuFEM，因为密砂中假设破坏滑裂面普遍偏大，尤其对于H:D=7.5过渡情况，见图5（c）。
对于深埋情况，由于采用浅埋破坏模式假设进行的计算，计算结果高于有限元结果，并且随着埋深增加，偏差越大，如Dr=30%，H:D=6.0和7.0情况，分别高估19%和35%，再次印证松砂H:D>5.0，中密砂H:D>6.0，已表现为深埋破坏模式。
4.2 浅埋锚极限上拔承载力理论公式及比较 4.2.1 等效内摩擦角沿着破坏滑裂面土体塑性应变不同，土体强度被激发的程度则不同，部分土体处于峰值状态，部分土体强度经历峰值后开始下降或仍未达到峰值强度，沿滑裂面各单元内摩擦角并非常值φp。Davis等[29-31]考虑土体非关联流动性或应变软化，采用式（4）计算极限荷载时破坏面（间断面）上的残余强度参数φ*。
| $ \tan \;{\phi ^*} = \dfrac{{\sin \;{\phi _{\text{p}}}\cos \;{\psi _{\text{p}}}}}{{1 - \sin \;{\phi _{\text{p}}}\sin \;{\psi _{\text{p}}}}} $ | (4) |
φ*实际上为在利用极限分析法确定极限荷载时，考虑土体非关联流动性而在间断面上对土体内摩擦角进行的折减[29-30]或考虑应变软化材料极限荷载时破坏面上的残余应力比[31]。Drescher等[30]指出对于承载力问题，采用该强度参数进行极限分析获得的理论值会低于基于非关联理想弹塑性本构模型的有限元分析结果，即意味着以理想弹塑性有限元计算获得的滑裂面上的应力比σt:σn应比由式（4）计算的结果高。
下面考察采用改进Mohr–Coulomb模型获得的浅埋锚滑裂面上的等效内摩擦角
| 表4 等效内摩擦角及对比 Tab. 4 Equivalent friction angles and comparisons |
 |
| $\begin{aligned}[b]& \tan \;{\phi _1^*} - \tan \;\psi =\\& \frac{{{Q_{\tau 1}}}}{{\displaystyle\sum\limits_{i = 1}^n {\left\{ {{\text{π}}\gamma '{K_{{\text{u}}i}}\left( {H - {\textit{z}_i}} \right)\left( {\textit{z}_{i+1}}- {{\textit{z}_i}} \right)\left[ {D + \left( {{\textit{z}_i} + {\textit{z}_{i{\text{ + }}1}}} \right)\tan \;{\psi _{\text{p}}}} \right]} \right\}} }} \end{aligned}$ | (5) |
由表4可见：松砂条件，采用Davis的折减内摩擦角
根据第3.2节分析的滑裂面上侧压力分布规律构造Ku函数，用于极限承载力理论分析。峰值前采用线性拟合，从某一值开始增加至Ku,peak；峰值后，Ku非线性快速下降直至接近Kn，采用指数函数进行拟合。为简化，z=0的初始点Ku可按Kn取值，Kn如式（1）所示；对于浅埋情况，峰值点的位置zp可统一取为0.2D。滑裂面上侧压力系数Ku的函数可表达成式（6），至此，确定变化函数的关键在于确定Ku,peak值及指数衰减函数的参数A，A表示衰减速度，A越大，衰减越快。
| $ {K}_{\text{u}}=\left\{\begin{array}{l}{K}_{\text{n}}\text+k\textit{z},\textit{z}\le \textit{z}_{\text{p}};\hfill \\ {K}_{\text{n}}\text+\left({K}_{\text{u,peak}}-{K}_{\text{n}}\right){\text{e}}^{-A\left(\textit{z}-\textit{z}_{\text{p}}\right)},\textit{z} > \textit{z}_{\text{p}}\hfill \end{array} \right.$ | (6) |
式中，k为直线段斜率，k=(Ku,peak–Kn)/zp。
图10为3种密实度砂土不同埋深比时峰值侧压力系数Ku,peak/Kn。图11为3种密实度砂土中浅埋锚滑裂面上Ku的变化及拟合结果（实线），对于Dr=30%、60%、100%时，A值分别取为1.0、1.5、3.0。
 |
| 图10 不同密实度土中不同埋深比时滑裂面上侧压力系数比 Fig. 10 Normalized peak earth lateral pressure coefficient for different conditions |
 |
| 图11 土侧压力系数Ku拟合结果 Fig. 11 Fitting results of lateral earth pressure coefficient Ku |
4.2.3 理论公式计算结果
将滑裂面上等效内摩擦角表达式（4）及滑裂面上土应力分布函数（6）代入式（2），则滑裂面上土抗力Qτ可表达成式（7）～（9）：
| $\begin{aligned}[b] &{Q}_{\tau }={\displaystyle {\int }_{0}^{H}\left({\tau }_{\text{f}}\mathrm{cos}\;\psi -{\sigma }_{\text{n}}\mathrm{sin}\;\psi \right)\cdot \text{π}\left(D+2\textit{z}\mathrm{tan}\;\psi \right)\dfrac{\text{d}\textit{z}}{\mathrm{cos}\;\psi }}\text{=}\\&\;\;\;\;\;\;\;\;\;{Q}_{\tau 直线段}\text+{Q}_{\tau 指数段}\\[-10pt] \end{aligned}$ | (7) |
| $ {Q}_{\tau 直线段}\text{=π}\gamma \left(\mathrm{tan}\;{\phi }^{\ast }-\mathrm{tan}\;{\psi }_{\text{p}}\right)\cdot {\left({B}_{1}\textit{z}+{B}_{2}{\textit{z}}^{2}+{B}_{3}{\textit{z}}^{3}+{B}_{4}{\textit{z}}^{4}\right)\bigr|}_{0}^{{\textit{z}}_{\text{p}}} $ | (8) |
式（7）～（8）中，
| $ \begin{array}{l}{Q}_{\tau 指数段}=\text{π}\gamma \left(\mathrm{tan}\;{\phi }^{\ast }-\mathrm{tan}\;{\psi }_{\text{p}}\right)\cdot \left({C}_{0}+{C}_{1}\textit{z}+{C}_{2}{\textit{z}}^{2}+{C}_{3}{\textit{z}}^{3}\right)\bigr|_{\textit{z}_{\rm{p}}}^H \end{array} $ | (9) |
式中，
对于过渡埋深的情况，由于假设的滑裂面比实际大会引起上拔力高估，松砂H:D=5.0，中密砂H:D=6.0，理论公式分别较数值计算结果高7.92%和7.17%。在这种情况，可以考虑采用下一等级密实度的侧压力衰减系数A进行计算，如对于松及中密砂过渡埋深情况，可取中密砂和密砂的A=1.5和3.0进行计算。
表5为针对浅埋锚采用理论公式计算的极限上拔力Qτ及与有限元计算结果的比较，其中，松砂及中密砂过渡埋深计算时采用了上述A的取值方法。结果表明，对于不同密实度的砂土，通过直线–指数型函数表达的侧压力系数Ku计算的上拔承载力与数值模拟结果吻合较好。
| 表5 理论公式计算结果与有限元计算结果对比 Tab. 5 Comparisons between theoretical and numerical results |
 |
采用Hao等[3]和Giampa等[4]的理论公式及Cerfontaine等[8]基于FEM结果的理论方法估算本文有限元算例。前两个理论公式的差别在于采用了不同的侧向土压力系数公式及是否考虑非关联法则的影响。对任意密实度砂土中的浅埋锚极限上拔承载力，Giampa等[4]结果比Hao等[3]结果高约20%～35%，埋深越大，两者结果相差越大。图12为浅埋锚极限上拔承载力不同计算公式的误差对比。其中，由于文献[8]中仅给出相对密实度Dr=50%～90%时侧压力系数图，因此，对于松砂的算例，未采用Cerfontaine等理论方法计算，对于Dr=60%的算例，侧压力系数插值确定，对于Dr=100%的算例，采用Dr=90%的侧压力系数计算。
 |
| 图12 浅埋锚极限上拔承载力不同计算公式误差对比 Fig. 12 The errors of different formulas for ultimate uplift capacity of shallow circular anchors |
由图12可见：松砂时Giampa等[4]和Hao等[3]理论公式均低估锚的上拔承载力，低估量范围分别在–4%～–23%、–51%～–25%，基于本文FEM结果的理论公式估算相对较好；对于中密砂，Giampa等[4]和本文基于FEM的理论公式估算结果较好，误差在–11%～6%之间，而Hao等结果仍低估锚上拔承载力–13%～–34%，Cerfontaine等[8]结果高估承载力16%～60%；对于密砂，Hao等[3]及本文基于FEM的理论公式计算结果误差在–5～7%以内，而Giampa等结果高估了16%～32%，Cerfontaine等结果则高估了19%～105%，且埋深越大误差越大，主要由于Giampa等[4]和Cerfontaine等[8]高估了密砂中滑裂面附近土体的侧压力，且未考虑应变软化带来的土体强度降低。
综上，从安全和经济的角度，砂土中浅锚上拔承载力计算可采用本文基于FEM结果的理论方法计算不同密实度的情况，但计算相对繁琐；亦可采用Giampa等[4]和Hao等[3]公式，建议当为松砂和中密砂时，采用Giampa等[4]公式计算，内摩擦角不折减；当为密砂时，采用Hao等[3]公式计算，并考虑内摩擦角的折减。
5 结 论本文利用Abaqus中3维欧拉大变形有限元分析方法，同时采用能够考虑应变软化的改进Mohr–Coulomb模型对圆锚拉拔过程进行数值模拟。通过与离心机试验结果进行对比，验证了数值模型的有效性，进而探讨了不同密实度砂土中锚的拉拔破坏机理，并基于数值分析结果提出了浅埋圆锚或螺旋锚上拔承载力理论公式，同时对比分析已有的理论计算公式，建议各公式的适用条件，为工程应用提供参考。具体结论如下：
1）随着埋深增加，锚周围土体的破坏模式不同，从浅埋时的倾斜角度接近于剪胀角的倒锥台型滑裂面向底部臌胀顶部渐缩的曲面型过渡模式发展。对于松砂、中密砂和极密砂，埋深比H:D分别小于等于4.0、5.0、6.0时发生浅埋破坏模式；松砂和中密砂，当H:D分别大于6.0和7.5时会发生深埋破坏模式；当H:D=7.5时极密砂仍为过渡破坏模式。
2）对于浅埋及过渡埋深比时，破坏滑裂面均可按直线表示，沿着滑裂面土体侧压力系数Ku并非常值，自锚盘边缘到地表，呈先线性增加至峰值Ku,peak，而后近似以指数形式衰减至稳定值；埋深比相对于土体密实度对Ku,peak有更重要的影响；Ku,peak位置基本不受H:D和Dr的影响，Ku衰减速率明显受密实度影响，密实度越高，衰减越快，Ku最后的稳定值接近Hao等[3]提出的土侧压力系数。
3）上拔破坏时，沿着滑裂面上，自锚盘边缘开始内摩擦角逐渐增大至峰值内摩擦角，之后又开始减小。当进行极限平衡分析时，采用等效内摩擦角φ*将滑裂面上变化的内摩擦角以常值表示，对中密砂土，φ*满足Davis等[29]提出的内摩擦角折减公式；对于松砂，Davis等公式略低估φ*；而对于极密砂，该公式高估φ*约3°。
4）通过各理论计算方法的比较，表明，对于任意密实度砂土中浅埋圆锚或螺旋锚，采用本文基于FEM的理论公式计算极限上拔承载力，既安全又经济，但计算相对繁琐；亦可采用Giampa等[4]和Hao等[3]公式，建议当为松砂和中密砂条件时，采用Giampa等[4]公式计算，内摩擦角不折减；当为密砂时，采用Hao等[3]计算公式，并考虑内摩擦角的折减。
致谢：感谢中国海洋大学王栋教授和郑敬宾副教授为本研究提供相关程序。
| [1] |
Ghaly A M,Hanna A M,Hanna M S. Uplift behavior of screw anchors in sand.I:Dry sand[J]. Journal of Geotechnical Engineering,ASCE, 1991, 117(5): 773-793. DOI:10.1061/(ASCE)0733-9410(1991)117:5(773) |
| [2] |
Hao Dongxue,Fu Shengnan,Chen Rong,et al. Experimental investigation of uplift behavior of anchors and estimation of uplift capacity in sands[J]. Chinese Journal of Geotechnical Engineering, 2015, 37(11): 2101-2106. [郝冬雪,符胜男,陈榕,等. 砂土中锚板拉拔模型试验及其抗拔力计算[J]. 岩土工程学报, 2015, 37(11): 2101-2106. DOI:10.11779/CJGE201511023] |
| [3] |
Hao Dongxue,Wang Dong,O’Loughlin C D,et al. Tensile monotonic capacity of helical anchors in sand:interaction between helices[J]. Canadian Geotechnical Journal, 2019, 56(10): 1534-1543. DOI:10.1139/cgj-2018-0202 |
| [4] |
Giampa J R,Bradshaw A S,Schneider J A. Influence of dilation angle on drained shallow circular anchor uplift capacity[J]. International Journal of Geomechanics, 2017, 17(2): 04016056. DOI:10.1061/(ASCE)GM.1943-5622.0000725 |
| [5] |
Merifield R S,Sloan S W. The ultimate pullout capacity of anchors in frictional soils[J]. Canadian Geotechnical Journal, 2006, 43(8): 852-868. DOI:10.1139/T06-052 |
| [6] |
Bhattacharya P,Kumar J. Uplift capacity of anchors in layered sand using finite-element limit analysis:formulation and results[J]. International Journal of Geomechanics, 2016, 16(3): 04015078. DOI:10.1061/(ASCE)GM.1943-5622.0000560 |
| [7] |
Tagaya K,Tanaka A,Aboshi H. Application of finite element method to pullout resistance of buried anchor[J]. Soils and Foundations, 1983, 23(3): 91-104. DOI:10.3208/sandf1972.23.3_91 |
| [8] |
Cerfontaine B,Knappett J A,Brown M J,et al. Effect of soil deformability on the failure mechanism of shallow plate or screw anchors in sand[J]. Computers and Geotechnics, 2019, 109: 34-45. DOI:10.1016/j.compgeo.2019.01.007 |
| [9] |
Vesić A S.Expansion of Cavities in Infinite Soil Mass[J].Journal of the Soil Mechanics and Foundations Division,1972,98(3):265-290.
|
| [10] |
Zhuang P Z,Yu H S. Uplift resistance of horizontal strip anchors in sand:A cavity expansion approach[J]. Géotechnique Letters, 2018, 8(4): 284-289. DOI:10.1680/jgele.18.00083 |
| [11] |
Majer J. Zur berechnung von zugfundamenten[J]. Osterreichische Bauzeitgschrift, 1955, 10(5): 85-90. |
| [12] |
Bobbitt D E,Clemence S P.Helical anchors:application and design criteria[C].Proceedings of the 9th Southeast Asian Geotechechnical Conference,Bangkok:Southeast Asian Geotechnical Society,1987:105–120.
|
| [13] |
Meyerhof G G,Adams J I. The ultimate uplift capacity of foundations[J]. Canadian Geotechnical Journal, 1968, 5(4): 225-244. DOI:10.1139/t68-024 |
| [14] |
Mitsch M P,Clemence S P.The uplift capacity of helix anchors in sand[C].Proceedings of Uplift Behavior of Anchor Foundations in Soil.New York:American Society of Civil Engineers(ASCE),1985:26-47.
|
| [15] |
Baker W H,Konder R L. Pullout load capacity of a circular earth anchor buried in sand[J]. Highway Research Record, 1966, 108(1): 1-10. |
| [16] |
Murray E J,Geddes J D.Uplift anchor of plates in sand[J].Journal of Geotechnical Engineering,1987,113(3):202–215.
|
| [17] |
Saeedy H S. Stability of circular vertical earth anchors[J]. Canadian Geotechnical Journal, 1987, 24(3): 452-456. DOI:10.1139/t87-056 |
| [18] |
Adams J I,Klym T W. A study of anchorages for transmission tower foundations[J]. Canadian Geotechnical Journal, 1972, 9(1): 89-104. DOI:10.1139/t72-007 |
| [19] |
Caquot A,Kérisel J,Bec M A.Tables for the calculation of passive pressure,active pressure and bearing capacity of foundations[M].Paris:Gauthier-Villars,1948.
|
| [20] |
White D J,Cheuk C Y,Bolton M D. The uplift resistance of pipes and plate anchors buried in sand[J]. Géotechnique, 2008, 58(10): 771-779. DOI:10.1680/geot.2008.3692 |
| [21] |
Chen Rong,Fu Shengnan,Hao Dongxue,et al. Scale effects of uplift capacity of circular anchors in dense sand[J]. Chinese Journal of Geotechnical Engineering, 2019, 41(1): 78-85. [陈榕,符胜男,郝冬雪,等. 密砂中圆形锚上拔承载力尺寸效应分析[J]. 岩土工程学报, 2019, 41(1): 78-85.] |
| [22] |
Chen Zongrui,Tho K K,Leung C F,et al. Influence of overburden pressure and soil rigidity on uplift behavior of square plate anchor in uniform clay[J]. Computers and Geotechnics, 2013, 52: 71-81. DOI:10.1016/j.compgeo.2013.04.002 |
| [23] |
Qiu Gang,Henke Sascha,Grabe Jürgen. Application of a Coupled Eulerian–Lagrangian approach on geomechanical problems involving large deformations[J]. Computers and Geotechnics, 2011, 38(1): 30-39. DOI:10.1016/j.compgeo.2010.09.002 |
| [24] |
Hsu S T,Liao H J. Uplift behaviour of cylindrical anchors in sand[J]. Canadian Geotechnical Journal, 1998, 35(1): 70-80. DOI:10.1139/cgj-35-1-70 |
| [25] |
Bolton M D. The strength and dilatancy of sands[J]. Géotechnique, 1986, 36(1): 65-78. DOI:10.1680/geot.1986.36.1.65 |
| [26] |
Liu Jinyuan,Liu Mingliang,Zhu Zhende. Sand deformation around an uplift plate anchor[J]. Journal of Geotechnical and Geoenvironmental Engineering, 2012, 138(6): 728-737. DOI:10.1061/(ASCE)GT.1943-5606.0000633 |
| [27] |
Shi Danda,Mao Yiyao,Yang Yong,et al. Experimental study on the deformation characteristics of soils around uplift circular plate anchors using digital image correlation technology[J]. Rock and Soil Mechanics, 2020, 41(10): 3201-3213. [史旦达,毛逸瑶,杨勇,等. 基于DIC技术的砂土中圆形锚板上拔土体变形特性试验研究[J]. 岩土力学, 2020, 41(10): 3201-3213.] |
| [28] |
Schiavon J,Tsuha C,Esquivel E.Experimental stress analysis in helical pile foundations by the photoelastic method[M].ICPMG 2014–Physical Modelling in Geotechnics.Boca Raton:CRC Press, 2013:757-762.
|
| [29] |
Davis E H.Theories of plasticity and the failure of soil masses[M].London:Butterworth,1968:341-380.
|
| [30] |
Drescher A,Detournay E. Limit load in translational failure mechanisms for associative and non-associative materials[J]. Géotechnique, 1993, 43(3): 443-456. DOI:10.1680/geot.1993.43.3.443 |
| [31] |
Vermeer P A. The orientation of shear bands in biaxial tests[J]. Géotechnique, 1990, 40(2): 223-236. DOI:10.1680/geot.1990.40.2.223 |