CDP模型是ABAQUS中专门为混凝土设计的非线性材料模型，适用于模拟混凝土结构的渐进式破坏行为，包括拉压异性的损伤演化机制、非关联塑性流动法则。

拉压异性的损伤演化机制通过独立定义拉伸损伤d_t和压缩损伤d_c来精确描述混凝土桩在受拉和受压的不同手里状态下刚度退化与破坏，其中受拉行为以微裂纹为主导，损伤表现为刚度逐渐丧失（软化），最终完全开裂（d_t=1）；压缩行为则以塑性变形和局部压溃为主，损伤体现为承载力下降（硬化后软化），但残余强度仍可能存在（d_c<1）。

采用非关联流动法则（塑性势函数≠屈服函数），避免了传统关联法则导致的不合理体积膨胀，更符合混凝土的真实变形特征。其中通过膨胀角（ψ）控制塑性体积膨胀程度，模拟混凝土受压时的横向膨胀（如箍筋约束效应），偏心率（ϵ）调整塑性势形状，避免数值奇异性，影响收敛。

为确定塑性损伤模型（CDP）参数，本文根据中国《混凝土结构设计规范》（GB 50010—2010）中计算公式计算确定^[12]。

以C40混凝土为混凝土桩材料，通过规范所给计算公式绘制Excel表格计算最终得出CDP参数取值（表1）。

表1 CDP参数Table 1 CDP parameters

其中：E_c为弹性模量；ν为泊松比；ψ为膨胀角；ϵ为偏心率；f_b0/f_c为双轴抗压强度与单轴抗压强度的比值；K为屈服面的形状参数；μ为粘性参数。

受压行为下应力和非弹性应变以及损伤曲线图及受拉行为下应力与开裂应变以及损伤曲线图见图1（a）和（b）：

图1 CDP损伤参数曲线图Fig.1 CDP damage parameter curves

本文并未对土层材料进行变量研究，采用砂土作为土体研究，本构模型为摩尔库伦塑性模型^[13]。

表2 土体参数Table 2 Soil parameters

其中：ρ为密度；E_c为弹性模量；ν为泊松比；φ为摩擦角。

为减少有限元分析计算量，模型采用四分之一对称模型，在土体底部采用固定约束，桩体、抱夹器以及土体的X与Y方向采用对称约束。

选取混凝土桩长为12 m，桩径d桩=0.5 m，管壁厚为0.1 m。网格尺寸为0.05 m，网格单元个数为384个，采用C3D8R单元。网格划分如图2。

图2 桩体网格划分Fig.2 Pile meshing

选取抱夹器外径D=1.4 m，内径d=0.5 m，高度为0.3 m。网格尺寸为0.05 m，网格单元个数为54个，采用C3D8R单元。网格划分如图3。

图3 抱夹器网格划分Fig.3 Clamping device meshing

为模拟土体特征性，土体单元采用欧拉单元，总尺寸为半径5 m，土厚28 m，其中上部预留3 m空气域，在空气域中赋予空气材料属性，使其参与能量交换，可作为缓冲层，缓解土体表面因突变荷载引起的数值振荡进而保持网格的几何规律性，提高收敛性，能够真实反应土体和空气界面的动态变化。

土体网格尺寸为0.2 m，网格单元个数119 140，网格单元类型为欧拉单元EC3D8R。网格划分如图4。

图4 土体网格划分Fig.4 Soil meshing

由于沉桩过程是一个冲击过程，所以使用显式分析。在分析步长选择上，过长的分析步长度会导致惯性效应失真以及能量误差积累；过短的分析步会导致大量增量步，计算时间激增。根据前人计算公式和经验，在求解冲击问题中，即使外荷载已经移除，模型仍会因弹性回弹或惯性效应产生数值振荡，其原因是显式算法并无隐式算法的自动阻尼机制，需要人工设置或者自然衰减，因此需要延长分析步时间。同时为了研究多次锤击，本模型设置了4个分析步，每个分析步长度为0.6 s^[14]。

前人已经证明通过轴对称和平面对称问题的数值算例中，通过对称约束以及荷载作用下，对比四分之一模型与完整模型的精度差异通常在误差1％内，但为有限元计算量的减少起到较大的帮助^[15]。本文主要研究抱夹器材料属性以及抱夹器高度对能量传递效率的影响，所以具体变量参数如下表（其中①号为顶部锤击）。

表3 对照组设计Table 3 Finite element control model design

在本节研究中，通过最终进入土体能量的曲线图分析每次锤击下进入土体能量变化规律。根据Spencer在《Continuum Mechanics》中从连续介质的基本能量守恒方程出发，并引入适用于打桩问题的简化条件：绝热过程，冲击时间短，热量来不及散失；准静态近似，当桩停止运动后，动能趋于零。简化后方程：

式（5）中， ρ为密度； e为单位质量内能； ta为时间；V为所研究土体的空间区域；De/Dta表示随材料质点运动的材料导数；b为单位质量体积力向量；t为作用于边界上的表面牵引力向量；v为质点速度向量。

随后通过应力功率（Stress Power）将外力做功与内能变化联系起来：

式（6）中，σ为柯西应力张量； 为应变率张量。

最后Spencer基于土体弹塑性本构模型，将应变率分解为弹性与塑性部分：

D=D^e+D^p（7）

式（7）中，D^e为弹性应变率张量；D^p为塑性应变率张量。

对应的应力功率为：

σ:D=σ:D^e+σ:D^p（8）

其中弹性部分是可逆的，存储为弹性应变能；而塑性部分为不可逆的，耗散为内能。

总结打桩过程的具体分析为冲击阶段：锤击力通过桩传递到土体表面，土体产生高应变率（D较大），应力功率σ:D显著；能量转化过程：土体颗粒发生塑性滑梯、破碎及孔隙压缩，塑性耗散（σ:D^p）直接转化为内能；平衡阶段：桩停止运动，动能趋近于0，外部功完全转化为土体内能，即∆e=∫σ:D^pdt^[16]。

综上所述，最终传入土体能量可近似看作土体最终的内能，即在Abaqus中Visualization模块中提取土体内能ALLIE，并绘制曲线。如图5.

图5 土体内能整体曲线Fig.5 Overall curve of soil internal energy

对图5曲线进行分析，每组曲线都呈现先上升→轻微下降→最终平衡的变化趋势。再上升阶段是由于在打桩初期，土体经历高应变率加载，材料迅速进入塑性流动状态，塑性耗散能快速积累，直接贡献了内能的增长，同时应力波传播需要时间覆盖整个模型，初期能量积累表现为单调上升， Belytschko指出显式分析中塑性耗散与时间积分同步增长，初期内能由不可逆变形主导16；当桩停止运动后，土体中未发生塑性屈服的区域释放弹性应变能产生弹性回弹并激发局部动态响应，部分内能转化为动能，导致内能出现短暂下降；最后土体塑性变形趋于稳定，不产生新的塑性功，塑性耗散饱和，系统的阻尼将残余动能完全耗散为热力学平衡状态内能，达到了系统阻尼平衡，热力学平衡时内能变化率满足 ，内能趋稳^[17]。

在能量效率的计算中，本文选取每一击结束时土体能量的增量与单次输入能量相比作为评价指标，而不是直接采用能量曲线在峰值处的读数。这是因为能量曲线在冲击峰值阶段往往会受到动能与弹性应变能的相互转换影响，曲线在峰后会出现轻微下降，直到动能逐渐耗散、塑性耗散趋于稳定后才达到真实的“净吸收”水平。因此，采用每一击结束时的能量增量，能够更准确地反映该次冲击对土体的实际能量贡献。

因此为方便进行比较分析，提取每段曲线中每个周期土体内能稳定时数值作为每个周期最终的能量值，简化曲线图，如图6（a）和（b）。

图6 简化土体能量曲线Fig.6 Simplified soil energy curves

根据图6分析可知，无论是抱夹器材料属性为Q345还是Q235传入土体能量整体随着高度的增加而减小；在（a）中，h=1/4处比h=1/2处传入土体能量多33.2%，h=1/2处比h=3/4处传入土体能量多22.6%，h=1/2处比h=95%处（近似顶部）传入能量多20.8%；在（b）中，h=1/4处比h=1/2处传入土体能量多34.6%，h=1/2处比h=3/4处传入土体能量多21.1%，h=1/2处比h=95%处（近似顶部）传入能量多21.8%，可知虽然位置都相差1/4的高度，但是当抱夹器处于桩体高度的下半部分时，输入能量大小变化对高度更为敏感，而在上半部分抱夹器每下降一个单位传入土体能量值就会较上一个单位提高21％上下，同时无论抱夹器材料属性为Q345还是Q235都不会改变曲线整体的变化趋势以及不同高度之间的相对大小关系。

从时间维度来看，随着周期数增加，除h=1/4处由于传入能量显著过大导致土体压制后剪切强度降低，后续加载更易变形，能量传入需求减少，其余位置随着周期数增加，每次传入能量值没有明显变化，比较稳定。

在与顶锤式作比较下，1/4位置组、1/2位置组以及3/4位置组都比顶锤式方法有较好的能量传递结果，其原因减少了能量在桩内传递的距离，使得降低能量传播中在桩体材料的内部耗散以及应力波的几何扩散与衰落；而95%位置组比顶锤式方法传递效果略低，尽管减少了能量在桩内传播的距离从而减少消耗，但同时能量经过抱夹器传入桩时会发生应力波反射导致部分能量的损失高过前者减少的消耗，以至于无法有比顶锤式方法更好的沉桩效果。综上可抱夹式锤击方法的沉桩效果与顶锤式方法沉桩效果等效高度在3/4处~顶部之间。

为了进一步揭示抱夹器抱夹高度对能量传递规律的影响，本文在土体能量增量统计的基础上，将每一击结束时的土体能量增量与单击外部输入能量进行对比，得到分击能量传递效率随击次的变化情况。这一处理方式不仅能够反映不同安装高度下能量传递的差异，还能揭示多次锤击过程中能量利用的边际变化规律上述单从能量绝对值方面比较，现分析各模型的能量传递效率，公式为：

其中初始地应力内能为596 kJ，初始动能为1.67×10⁴ kJ。取平衡后土体内能计算后绘制曲线图7。

如图7，逐一分析可知，1/4处组的效率随着锤击次数增加，效率有着比较明显的下降趋势，近似从72%下降至46%；1/2处组的效率曲线波动较大，呈现先减小再增大最后降低的规律；而3/4处组与顶锤式方法的效率曲线相对比较平稳，没有明显的波动；95%处组在第一次锤击后效率与顶锤式方法一致，随后出现降低并平稳。

图7 土体能量效率曲线Fig.7 Soil energy efficiency curves

为进一步分析抱夹器材料强度对能量传递效率的影响，本文选取Q235与Q345两种钢材作为代表性材料进行对比。两者弹性模量相近（E≈210 GPa），但屈服强度差异显著，Q345较Q235具有更高的承载力与抗变形能力，可在冲击加载过程中减少塑性耗能与结构变形，从而影响能量传递特性。不同抱夹高度下，Q235与Q345两种材料对应的传入土体能量如图8所示。

图8 同一位置不同抱夹器材料属性土体能量对比图Fig.8 Comparison of soil energy under different clamp material properties at the same position

图8所示为1/4位置处、1/2位置处、3/4位置处以及95％位置处下，抱夹器金属材料属性分别为Q345与Q235下土体内能简化曲线图，由四幅曲线图分析得，每一位置的两条曲线十分接近甚至出现部分重合情况。在抱夹器材料属性（Q345与Q235）差异下土体内能图未发生明显变化的根本原因在于，土体自身的能量耗散机制主导了能量传递过程，而抱夹器作为能量传递媒介的作用为因材料强度差异发生显著变化。

以下从土体响应主导性与能量传递路径两方面展开分析：

1）土体在动态加载中能量吸收与耗散主要以下途径：①土体颗粒的不可逆重排/压实以及剪切的形成消耗了大部分能量（占输入能量的60%~80%）。土体塑性变形对输入能量幅值和频率高度敏感，但对传递介质的材料强度不敏感；②土体的率相关特性导致能量随应变率升高而耗散。此过程取决于土体本构模型而非抱夹器材料属性；③土体非均质性引发应力波散射，能量随传播距离按几何衰减和材料阻尼衰减规律耗散。此类衰减由土体属性决定，与抱夹器材料无关。综上三点对土体能量的影响远高于抱夹器材料属性差异带来的影响；

2）若抱夹器在动态加载中未达到屈服，两者刚度相同，应力波传播速度 相近，导致输入土体的波形和能量密度无显著差异；且土体对低频能量的吸收效率高于高频能量。抱夹器材料差异未改变激荡频率特性，土体内能分布模式将保持一致。

综上，土体内能曲线图中未因抱夹器材料有明显变化而改变的核心原因是：土体自身的耗散机制主导了能量吸收过程，而抱夹器材料差异未显著改变输入土体的能量特征。

尽管曲线上十分接近，但无论哪个高度位置下，在效率方面，Q345抱夹器会比Q235高3‰~5‰。其主要原因是屈服强度差异导致Q345在相同动荷载下不易发生局部塑性变形，从而减少能量在抱夹器自身塑性耗散的损失。平均能量传递效率见表4。

表4 平均能量传递效率Table 4 Average energy transfer efficiency