基于N-R 法的某水闸桩基负载及变形仿真分析

程天竞

（淮安市水利工程建设管理服务中心，江苏 淮安 223300）

1 引言

桩土的相互作用十分复杂[1]，求解亦为非线性求解[2]，模型较大且庞杂的情况下，求解过程可能会严重迟缓[3]，收敛也会极其缓慢[4]，甚至会出现计算不收敛的情况[5]。造成非线性行为的原因可归结为3 种情况[6，7]：①几何非线性：在水力驱动下，构建模型出现较大的变形而引发结构的非线性响应。②材料非线性：材料会由于加载历史、周围环境以及加载的时间跨度不同等因素，对材料的应力应变性质产生影响。③状态非线性：桩土接触状态会引发结构非线性，结构受力会随着状态的变化而变化。桩土模型存在的问题属于第③种情况时，建模过程中若能够全面地去考虑问题、分析原因并及时调整一些设置，选择合理的求解器，可使计算更迅速、更易收敛。概而言之，择取合适的桩土非线性求解方法具有重要意义，而牛顿-拉普森迭代法(Newton-Raphson， N-R)能够克服传统方法的不足，提升计算性能。

某水闸工程为软土闸基。为控制建筑物水平位移和沉降，满足工程施工工期及导流要求。根据工程施工和闸址工程地质特点，确定了采用预制管桩基础处理方案。基于以上研究背景，依托该工程开展了基于N-R 法的预制管桩承载数值仿真的研究。

2 基于N-R 法桩土模型非线性求解方法

求解非线性结构时，传统的增量近似会出现误差积累，导致最终结果失衡[8]，传统增量近似与真实响应对比见图1。

牛顿-拉普森迭代法克服了传统方法的不足，其可以使解在每个载荷增量末尾均能达到收敛平衡。计算过程中，若收敛准则未达到要求，需进一步估计非平衡静载，并修改刚度矩阵，可再一次获得新的解，循环往复这种迭代直到收敛为止。牛顿-拉普森选代近似与真实响应对比如图2 所示。

图2 牛顿-拉普森迭代近似与真实响应结果对比

由图1、图2 可知，非线性静态分析解决一些不稳定的系统时，牛顿-拉普森迭代法适用性较强。

3 闸基桩模型建立

在借助有限元软件对桩基负载性能开展分析时，在确保模型的真实性和准确性的前提下，对模型进行了简化处理，主要为：①土体在试验前，由于自重作用下已经固结完毕，加载之前需进行地应力平衡，即通过施加自重应力产生初始位移，保留自重的重力场。②预制管桩为各向同性的弹性材料，且桩四周及底部范围内的土是均匀的。

3.1 模型尺寸

模型尺寸参考现场试验桩的尺寸：①竖向承载桩取M2 桩长30.00 m；②水平承载桩取N3 桩长27.00 m。两根单桩的外桩径均为0.80 m，内桩径均为0.58 m。将桩周土范围取为桩外直径的10 倍，为8.00 m；桩端土范围取到钻孔达到的位置，为42.00 m。

3.2 桩土参数

桩体参数弹性模量E=3.8×104MPa，泊松比ν=0.2，密度ρ=2 500 kg/m3，土体参数如表1 所示。

表1 土体参数

3.3 模型边界条件

模型尺寸横向取桩径10 倍，纵向尺寸也较长，因此模型边界条件对桩的影响很小。模型周边设横向及纵向约束，竖向不施加约束，底部设置全约束。

3.4 模型网格的剖分

采用自上向下建模方式直接生成实体。根据闸基实际情况先设定桩体与土体的每个线的段数，而后对构建的桩土模型进行单元剖分。M2 桩取1/4模型，实体单元为6 020 个，刚性接触单元为401 个，柔性接触单元为420 个，整个模型经过划分后单元个数为6 573 个；N3 桩取完整模型，实体单元为308 864 个，刚性接触单元为4 128 个，柔性接触单元为4 269 个，整个模型经过划分后单元个数为38 294 个。

3.5 模型接触单元

计算模型接触选用面-面接触，摩擦系数值取为0.20。桩体、土体分别选用刚、柔性接触面，其单元分别设为刚、柔性接触单元，即分别选取TARGE170、CONTA174 进行仿真。

3.6 模型加载方式

土体在未进行沉桩之前未受到大扰动，所有位移为零，而所受的应力并非为零。因此在施加竖向力与水平力计算之前，需先进行地应力平衡。具体步骤为：在完成建模后，对模型施加重力加速度，使构建的桩土模型产生自重应力及位移，得出结果并保存，然后保留应力，清除重力产生的位移。如图3和图4 所示。

图3 自重下土体位移

图4 平衡后土体位移

在施加竖向静载时，将竖直驱动力平均施加于桩体顶部各个节点；在施加水平静载时，将水平驱动力平均施加于桩顶部表面积受力侧各个节点。

4 计算结果分析

4.1 竖向静载下单桩有限元分析

管桩在竖向静载驱动力下，由于桩与桩周土体各方面条件均对称于桩的中心轴线。为了减少桩土模型总体计算量，计算时选用1/4 模型开展计算。计算完成后，提取桩体竖向位移及应力值。

由计算可知，桩身静载值沿竖向呈逐渐减小规律。这是因为桩顶受竖向静载的驱动，桩顶部压力沿桩身竖向进行传递，相对于桩周土体而言，则桩身出现了向下的位移；与之相反，桩身受土体向上的摩阻力驱动，使得桩身向下竖向静载与向上的摩擦力部分相互抵消，减小了竖向的桩身静载。此外，桩侧摩阻力与桩端阻力二者并非协同驱动发展，而是桩侧摩阻驱动力先沿着桩体从上而下逐步传递极限；进一步地，端阻驱动力逐级传递至极限。桩侧摩阻驱动力与桩端驱动阻力传递步调与桩土的相对位移量相关。

图5 竖向静载-位移

图6 桩身轴力

图5、图6 分别为M2 桩体在竖向力驱动下的位移图和轴力图。由图4、图5 可知，与试验测值相比，桩体位移沿桩身数值仿真值变化规律与其保持一致，呈逐渐变小规律，最大静载为720 kN 时，桩顶位移值约为13.50 mm。

4.2 水平静载下单桩有限元分析

由于水平静载不对称，管桩在水平静载的驱动下，不能只对局部模型进行计算，因此，水平静载下单桩有限元分析选用整体模型。计算完成后，提取桩体水平位移与弯矩值。

沿桩身长度方向提取各个节点的水平静载位移。需要说明的是，由于桩身弯矩不能直接进行提取，而是由单元节点力求和法求得。该方法通过对所求截面的节点和单元进行选取，将单元节点力对所选截面中的一点弯矩进行求和，可求出该截面的弯矩。该方法求解的截面内力精确度较低，但能满足一般工程所需精度要求，尤其适用于桩或者梁截面弯矩的求解。

由计算可知，水平静载下单桩位移向附近呈逐渐减小趋势，桩体前、后方土体出现的水平位移值均较大。这是可能是因为在水平静载驱动力下，桩体对土体挤压而形成了水平位移；反之，被挤压的土体反过来会对桩形成反作用驱动力，进而使另一侧土体产生了位移。土体随着桩体所受静载驱动力的逐渐增加，水平静载位移呈逐渐增加趋势。

图7 水平静载-位移

图8 桩身弯矩

图7、图8 分别为N3 桩体在水平驱动力下的仿真位移图和弯矩图。由图6、图7 可知，与试验测值相比，桩体水平静载位移沿桩身数值仿真值变化规律与其保持一致，呈逐渐减小趋势，在最大静载值为130 kN 时，除去最上面1.00 m 的桩头，其位移值约为17.40 mm，桩身最大弯矩值约为432.29 kN·m。

由于三维有限元分析中做出较多假设，例如：土体本构模型的限制、实际土体存在不均匀性以及所选取的参数等众多因素均会引起有限元软件计算结果和桩基载荷试验得到的结果存有差异，但总体变化规律趋于一致，该方法验证了有限元模拟的可行性。

5 结论

（1）在借助有限元软件仿真分析水力静载驱动下单桩的承载性能状况时，基于牛顿—拉普森迭代法非线性求解方法较传统求解方法仿真复杂桩土问题更迅速、更易收敛。

（2）仿真分析时众多假设因素均会引起有限元软件仿真结果与桩基载荷试验得到的结果存有微小差异，但总体变化规律趋于一致，该方法验证了有限元模拟的可行性。