组合受荷桩承载力传递矩阵通解方法与流程

本发明涉及一种组合受荷桩承载力传递矩阵通解方法。

背景技术：

桩基础一开始以承担上部结构竖向承载力为主要设计目的，随着工程需求的变化，出现了以承担水平荷载作用为主的水平受荷桩。然而，随着对基桩承载特性的研究的深入和理论的逐步完善，发现实际工程中的桩基通常不但承担桩顶竖向荷载作用，还同时受到水平荷载作用，即组合荷载的作用，如桥梁桩基础，海上风机单桩基础等。

组合荷载下水平力使桩身产生较大内力和位移，竖向分力也因桩身挠曲变形而产生附加弯矩，尤其是当地基土质较差、地面以上桩自由长度较大时，附加的桩身挠曲变形和弯矩更不可忽略，即所谓的P-Δ效应，因此组合荷载作用下基桩承载特性越来越受到重视。大量的学者开展了众多关于组合受荷桩承载特性的离心模型实验研究、室内模型试验研究以及二维、三维有限元数值模拟研究并得出诸多有益的结论。然而，这些研究成果只能反映特定情况下的桩身响应特性，且一般情况下这些研究成果反映的是承载特性定性结论，并不能为实际的工程设计人员提供很好的设计依据。因此，组合受荷桩桩身响应的解析模型的建立以及四阶微分控制方程的如何求解是得到桩身响应的关键。诸多学者相继开展了以有限差分法、幂级数法、变分法、传递矩阵法等方法为基础的组合受荷桩承载力计算方法研究，然而，有限差分法是一种纯数值方法，对编程能力要求非常高，不适合工程设计人员使用；幂级数法和变分法的推导过程和解的形式非常繁琐，极易出差且编程复杂，因而限制了其应用；现有的传递矩阵法其传递矩阵系数一般采用幂级数的方式表达且只能解决特定问题，因而限制了其应用范围。

技术实现要素：

本发明的目的在于克服现有技术中的不足，提供一种组合受荷桩承载力传递矩阵通解方法，解决现有技术中组合受荷桩承载力计算方法步骤繁琐、易出差、应用范围受限的技术问题。

为解决上述技术问题，本发明所采用的技术方案是：组合受荷桩承载力传递矩阵通解方法，包括如下步骤：

步骤一：建立组合受荷桩通解受力分析模型并提出便于求通解时的桩土相互作用关系的统一表达式；

步骤二：根据传递矩阵的一维线性特征将土抗力系数、桩身轴力分布、桩身水平荷载分布进行离散化、均匀化以及常数化处理；

步骤三：采用Laplace正逆变换推导得出了桩身自由段、弹性段以及塑性段的传递矩阵系数的解析通解；

步骤四：在传递矩阵法原理基础上代入桩顶桩端边界条件，从而得出组合受荷桩承载力通解。

步骤一的具体步骤如下：

将线弹性、非线性弹性、线弹性-塑性以及非线性弹性-塑性桩土相互作用模型进行归类，提出采用统一线弹性-塑性p-y曲线表达式模拟组合受荷桩桩土相互作用，如下式所述：

式中：p为土抗力；k为土抗力模量；y为桩身变形；p_u为极限土抗力。

步骤二的具体步骤如下：

根据土层分层以及桩身截面尺寸的变化，对地基土体重新进行n数量的分层，其中对第i层地基土体中的桩身再进行m_i数量的等分，则第j小段桩身平均轴向荷载为

其中：为第i层地基土中经过m_i等分后第j小段桩的桩身平均轴力荷载作用；

V_(i,j-1)，V_(i,j)分别为第i部分桩中第j小段桩的桩顶和桩端位置处的轴力荷载作用；

对于地面下桩土相互作用，根据其弹性、塑性状态分为以下两种：

①当桩土相互作用为弹性阶段时，第i层地基土中经过m_i数量的等分后第j小段桩的桩侧平均土抗力模量为：

其中：为第i层地基土中经过m_i数量的等分后第j小段桩的桩侧平均土抗力模量；

k_(i,j-1)，k_(i,j)分别为第i部分桩中第j小段桩的桩顶和桩端位置处的土抗力模量；

②当桩土相互作用为塑性阶段时，第i层地基土中桩经过m_i数量的等分后任意j小段桩侧极限土抗力为

其中：为第i层地基土中桩任意j小段桩侧平均极限土抗力值；

p_u(i,j-1)，p_u(i,j)分别为第i部分桩中第j小段桩的桩顶和桩端位置处的极限土抗力。

步骤三的具体步骤如下：

采用Laplace正逆变换对桩身自由段、弹性段以及塑性段的微分控制方程进行求解，得出以下传递矩阵系数解析解：

①对于考虑桩顶竖向荷载作用时的地面上自由段桩，其自由段第w小段桩的桩身传递矩阵系数表达式如下：

其中：l_a为自由段桩进行m_a数量等分后的每小段长度；为桩身自由段经过m_a等分之后任意w小段桩身平均轴向荷载作用；为桩身自由段经过m_a等分之后任意w小段桩侧所受的平均水平荷载作用；ψ_a(w)＝α_a(w)l_a，EI_a为桩身抗弯刚度；

②对于不考虑桩顶竖向荷载作用时的地面上自由段桩，其自由段第w小段桩的桩身传递矩阵系数表达式如下：

③对于考虑桩顶竖向荷载作用时的地面下弹性段桩，其弹性阶段的第i部分桩中第j小段桩的传递矩阵系数表达式如下：

式中：l_e(i,j)为弹性段进行m_i数量的等分后的第j小段桩的长度；ω_x(i,j)＝±(γ_i,j±ξ_i,ji)，x＝1、2、3和4，EI_i为第i层地基土位置处桩身抗弯刚度；

④对于不考虑桩顶竖向荷载作用时的地面下弹性段桩，其弹性阶段的第i部分桩中第j小段桩的传递矩阵系数表达式如下：

式中：t_e(i,j)＝β_(i,j)l_e(i,j)；

⑤对于考虑桩顶竖向荷载作用时的地面下塑性段桩，其塑性阶段的第i部分桩中第j小段桩的传递矩阵系数表达式如下：

式中：t_p(i,j)＝α_p(i,j)l_p(i,j)，l_p(i,j)为塑性段等分后的第j小段桩的长度，EI_i为第i层地基土位置处桩身抗弯刚度；

⑥对于不考虑桩顶竖向荷载作用时的地面下塑性段桩，其塑性阶段的第i部分桩中第j小段桩的传递矩阵系数表达式如下：

根据传递矩阵原理，结合所推导的传递矩阵系数可得整个桩传递矩阵方程为：

式中：和S₀＝[y₀ θ₀ M₀ Q₀ 1]^T分别为整个桩的桩端和桩顶的位移、转角、弯矩和剪力状态量；

U_(i,j)为自由段矩阵传递系数或弹性段传递矩阵系数或塑性段传递矩阵系数，U为整个桩身的总传递矩阵系数，是从整个桩的桩端到桩顶的每小段桩的传递矩阵系数的连续乘积；

根据沿深度的桩土相互作用状态以及是否存在自由段桩，总传递矩阵系数U求解过程如下：

①如果地面以上存在自由段桩时，则自由段桩身总传递矩阵系数U_a如下式：

则总传递矩阵系数U根据沿深度的桩土相互作用状态按照下式进行计算：

②如果地面以上不存在自由段桩时，则总传递矩阵系数U根据沿深度的桩土相互作用状态按照下式进行计算：

联合桩顶、桩端的边界条件和总传递矩阵系数U，即可解得桩顶变形、转角、弯矩和剪力的状态量S₀，进而按照下式计算任意位置处的桩身响应：

步骤四的具体步骤如下：

根据传递矩阵原理，代入桩顶与桩端的边界条件进行建立桩身响应的传递矩阵表达式，并进行桩身变形与内力的求解，计算收敛关键步骤如下：

①根据计算模型确定相应的基本计算参数并进行离散化与常数化处理；对于桩土相互作用模型为非线性模型时需假定一个初始土抗力模量k；在计算初始阶段需先假定桩土相互作用均为弹性状态；

②根据上一次迭代的计算结果，如果存在非线性的桩土相互作用p-y曲线，则需要判断经过r次迭代后割线刚度是否满足公式的要求，如果满足割线刚度迭代精度要求，则判断经过s次迭代后整体精度是否满足公式的要求；

其中，表示经过第r次迭代后的第i部分桩中第j小段桩的桩土相互作用p-y曲线割线刚度；ε表示迭代收敛标准推荐值；表示经过s次迭代后从地表往桩身方向发展的塑性区深度；

如果满足要求则进行桩身响应最终状态求解，如果不满足上述2种迭代要求中的任意一种，则需要进行再次迭代直到两种状态的迭代精度均满足要求。

自由段以及地面下桩身等分数量m_a和m_i确定原则为其等分后每小段的桩身长度要小于等于相应位置处桩身直径的二十分之一。

对于桩土相互作用模型为非线性模型时，取该p-y曲线模型中y＝0.01mm所对应的割线刚度作为初始土抗力模量k。

迭代收敛标准推荐值为ε＝0.005。

与现有技术相比，本发明所达到的有益效果是：

本发明从组合受荷桩解的通用性角度出发，结合传递矩阵法的一维线性特征，提出离散化与均匀化的处理方法将每一微段单元的所有计算参数状态转换成常数量，并采用Laplace正逆变换求解变系数四阶微分方程进而得出桩身自由段、弹性段以及塑性段的传递矩阵系数解析解表达式；本发明所阐述的计算方法能够考虑任意的桩顶竖向荷载、水平荷载、弯矩荷载以及桩身自由段水平分布荷载的组合作用、桩身轴力荷载作用、桩土相互作用p-y曲线模型、桩身截面尺寸和长度以及任意多土层数量，同时传递矩阵系数解析解表达式非常简洁，便于编程计算，具有非常好的适应性和通用性，能够很好地解决岩土工程中遇到的各种工况下的组合受荷桩承载特性问题。

附图说明

图1为组合受荷桩受力分析通解模型；

图2为桩土相互作用关系模型，其中图2(a)为线弹性模型，图2(b)为非线性弹性模型，图2(c)为线弹性-塑性模型，图2(d)为非线性弹性-塑性模型；

图3为地面上自由段桩身微分单元受力模型；

图4为地面下弹性段桩身微分单元受力模型；

图5为地面下塑性段桩身微分单元受力模型；

图6为沿桩身深度的桩土相互作用状态示意图，其中图6(a)为塑性-弹性模型，图6(b)为塑性-弹性-塑性-弹性模型，图6(c)为塑性-弹性-塑性模型，图6(d)为塑性-塑性模型；

图7为基于Laplace正逆变换的组合受荷桩承载力传递矩阵通解算法计算流程图；

图8为线弹性桩土相互作用模型对比验证；

图9为非线弹性桩土相互作用模型对比验证；

图10为线弹性-塑性桩土相互作用模型对比验证；

图11为非线弹性-塑性桩土相互作用模型对比验证，其中图11(a)为桩身变形对比，图11(b)为桩身弯矩对比；

图12为组合受荷桩桩身弯矩对比验证，其中图12(a)为第一组桩，图12(b)为第二组桩，图12(c)为第三组桩，图12(d)为第四组桩。

具体实施方式

下面结合附图对本发明作进一步描述。以下实施例仅用于更加清楚地说明本发明的技术方案，而不能以此来限制本发明的保护范围。

本发明从组合受荷桩解的通用性角度出发，结合传递矩阵法的一维线性特征，提出离散化与均匀化的处理方法将每一微段单元的所有计算参数状态转换成常数量，并采用Laplace正逆变换求解变系数四阶微分方程进而推导得出桩身自由段、弹性段以及塑性段的传递矩阵系数解析解表达式，最后通过公知的传递矩阵原理进行求解任意位置处的桩身响应。具体计算步骤如下：

(1)建立组合受荷桩通用受力分析模型，如图1所示，明确相应的各计算参数状态。图中L_a、L_b分别为桩位于地面上和地面下的长度；H_i、d_i、m_i、EI_i分别为地面下每层地基土中桩身的长度、直径、等分数量以及抗弯刚度；d_a、m_a、EI_a分别为地面上桩身的直径、等分数量以及抗弯刚度；自由段以及地面下桩身等分数量m_a和m_i确定原则为其等分后每小段的桩身长度要小于等于相应位置处桩身直径的二十分之一。

(2)如图2所示，将目前常用的线弹性、非线性弹性、线弹性-塑性以及非线性弹性-塑性桩土相互作用模型进行归类，并提出模拟组合受荷桩桩土相互作用的统一p-y曲线表达式，

其中：p为土抗力；k为土抗力模量；y为桩身变形；p_u为极限土抗力。因此，对于弹性模型而言(图2(a)和2(b))，采用p＝ky进行计算；当考虑塑性影响时(图2(c)和2(d))，其弹性阶段采用p＝ky计算，塑性阶段采用p＝p_u计算。

(3)结合传递矩阵法的一维线性特征，提出离散化与均匀化的处理方法将每一微段单元的所有计算参数状态转换成常数量。因此对于地面上自由段桩而言(见图3)，任意第w小段桩侧所受的水平分布荷载按下式计算：

其中为任意w小段桩侧所受的平均水平荷载作用，q_a(w-1)和q_a(w)为任意w小段桩顶和桩端位置处的水平荷载作用。

任意第w小段桩身所受平均轴力为：

其中为任意w小段桩身平均轴向荷载作用，V_a(w-1)和V_a(w)分别为任意w小段桩顶和桩端位置处的轴向荷载作用。

对于地面下弹性段桩而言(如图4)，第i层地基土中经过m_i等分后第j小段桩的桩侧平均土抗力模量为：

其中，为第i层地基土中经过m_i等分后第j小段桩的桩侧平均土抗力模量；k_(i,j-1)，k_(i,j)分别为第i部分桩中第j小段桩的桩顶和桩端位置处的土抗力模量。

对于地面下塑性段桩而言(如图5)，第i层地基土中桩的任意j小段桩侧平均极限土抗力采用下式计算：

其中为第i层地基土中桩任意j小段桩侧平均极限土抗力值；p_u(i,j-1)，p_u(i,j)分别为第i部分桩中第j小段桩的桩顶和桩端位置处的极限土抗力。

对于地面下桩而言(图4和图5)，对第i层地基土体中的桩身任意第j小段桩身平均轴向荷载采用下式计算：

其中，为第i层地基土中经过m_i等分后第j小段桩的桩身平均轴力荷载作用；V_(i,j-1)，V_(i,j)分别为第i部分桩中第j小段桩的桩顶和桩端位置处的轴力荷载作用。

则可得桩身响应变系数四阶微分控制方程为：

地面上自由段：

地面下弹性段：

地面下塑性段：

(4)为采用Laplace正逆变换求解自由桩微分方程，将公式(7)采用矩阵方式表述：

式中，S_a(w)、A_a(w)和f_a(w)为矩阵变量，其中S_a(w)为自由段桩中第w小段桩的桩端状态量，各矩阵变量如下所示：

S_a(w)＝[y_a(w) θ_a(w) M_a(w) Q_a(w)]^T (11)

式中：y_a(w)、z_a(w)、θ_a(w)、M_a(w)和Q_a(w)分别为地面上自由段桩经过等分m_a段后任意w小段桩的桩端位移、深度、转角、弯矩和剪力。

为使用Laplace正逆变换，令：

L[S_a(w)(z_a(w))]＝F_a(w)(s_a(w))，L[f_a(w)(z_a(w))]＝g_a(w)(s_a(w)) (14)

其中，符号“L”表示Laplace正变换；F_a(w)和g_a(w)均表示Laplace方程，s_a(w)表示Laplace方程的变量。根据Laplace正变换原理，对公式(10)进行Laplace正变换得：

F_a(w)(z_a(w))＝(I×s_a(w)-A_a(w))^-1S_a(w-1)+(I×s_a(w)-A_a(w))^-1g_a(w)(s_a(w))

对公式(15)采用Laplace逆变换进行变换，可得：

S_a(w)＝L^-1[(I×s_a(w)-A_a(w))^-1]S_a(w-1)+L^-1[(I×s_a(w)-A_a(w))^-1g_a(w)(s_a(w))] (16)

式中：符号“L^-1”表示Laplace逆变换。S_a(w-1)为自由段桩中第w小段桩顶状态量：

S_a(w-1)＝[y_a(w-1) θ_a(w-1) M_a(w-1) Q_a(w-1)]^T (17)

令：可得：

则公式(18)的Laplace逆变换为：

其中，ψ_a(w)＝α_a(w)l_a，l_a为自由段进行m_a数量等分后的每小段长度。

将公式(20)和(21)代入公式(16)，并根据传递矩阵法的原理可得：

S_a(w)＝U_a(w)S_a(w-1) (22)

式中，U_a(w)为自由段第w小段桩的传递矩阵系数，根据公式(20)和(21)可得：

当不考虑桩身轴力荷载作用时，对公式(23)求极限可得：

(5)为采用Laplace正逆变换求解弹性桩微分方程，将公式(8)采用矩阵方式表述：

式中，S_e(i,j)、A_e(i,j)和f_e(i,j)为矩阵变量，其中S_e(i,j)为位于弹性阶段的第i部分桩中第j小段桩的桩端状态量，各矩阵变量如下所示：

S_e(i,j)＝[y_e(i,j) θ_e(i,j) M_e(i,j) Q_e(i,j)]^T (26)

f_e(i,j)＝[0 0 0 0]^T (27)

对公式(25)进行Laplace逆变换可得：

S_e(i,j)＝L^-1[(I×s_e(i,j)-A_e(i,j))^-1]S_e(i,j-1)+L^-1[(I×s_e(i,j)-A_e(i,j))^-1g_e(i,j)(s_e(i,j))] (29)

式中，符号“L”表示Laplace正变换；F_e(i,j)和g_e(i,j)均表示Laplace方程，s_e(i,j)表示Laplace方程的变量。S_e(i,j-1)为位于弹性阶段的第i部分桩中第j小段桩的桩顶状态量，也是第(j-1)小段桩的桩端状态量：

S_e(i,j-1)＝[y_e(i,j-1) θ_e(i,j-1) M_e(i,j-1) Q_e(i,j-1) 1]^T (30)

式中：y_e(i,j-1)、θ_e(i,j-1)、M_e(i,j-1)和Q_e(i,j-1)为第j小段桩的桩顶的位移、转角、弯矩和剪力。

其中，

(I×s_e(i,j)-A_e(i,j))^-1g_e(i,j)(s_e(i,j))＝[0 0 0 0]^T (32)

在采用Laplace逆变换求解过程中，需要解得如下方程的根：

实际上大部分桩的桩身轴力则可得方程(33)根为：

ω_x(i,j)＝±(γ_i,j±ξ_i,ji)(x＝1,2,3,4) (34)

式中：符号“i”(正体)表示虚数，且令l_e(i,j)为弹性段进行m_i数量的等分后的第j小段桩的长度，则可得：

L^-1[(I×s_e(i,j)-A_e(i,j))^-1g_e(i,j)(s_e(i,j))]＝[0 0 0 0]^T (36)

将公式(35)和(36)代入公式(29)，并根据传递矩阵法的原理可得：

S_e(i,j)＝U_e(i,j)S_e(i,j-1) (37)

式中，U_e(i,j)为位于弹性阶段的第i部分桩中第j小段桩的传递矩阵系数，根据公式(35)和(36)可得：

当不考虑桩身轴力荷载作用时，对公式(38)求极限可得：

式中：t_e(i,j)＝β_(i,j)l_e(i,j)。

(6)为采用Laplace正逆变换求解弹性桩微分方程，将公式(9)采用矩阵方式表述：

式中，S_p(i,j)、A_p(i,j)和f_p(i,j)为矩阵变量，其中S_p(i,j)为位于塑性阶段的第i部分桩中第j小段桩的桩端状态量，各矩阵变量如下所示：

S_p(i,j)＝[y_p(i,j) θ_p(i,j) M_p(i,j) Q_p(i,j)]^T (41)

公式(40)的Laplace逆变换为：

S_p(i,j)＝L^-1[(I×s_p(i,j)-A_p(i,j))^-1]S_p(i,j)+L^-1[(I×s_p(i,j)-A_p(i,j))^-1g_p(i,j)(s_p(i,j))] (44)

其中，符号“L”表示Laplace正变换；F_p(i,j)和g_p(i,j)均表示Laplace方程，s_p(i,j)表示Laplace方程的变量。

根据自由段、弹性段的基于Laplace正逆变换求解过程可得塑性段桩身传递关系为：

S_p(i,j)＝U_p(i,j)S_p(i,j-1) (45)

式中，U_p(i,j)为位于塑性阶段的第i部分桩中第j小段桩的传递矩阵系数，为：

式中：t_p(i,j)＝α_p(i,j)l_p(i,j)，l_p(i,j)为塑性段等分后的第j小段桩的长度，EI_i为第i层地基土位置处桩身抗弯刚度。

当不考虑桩身轴力荷载作用时，对公式(46)求极限可得：

(7)根据传递矩阵原理，结合所推导的传递矩阵系数可得整个桩传递矩阵方程为：

式中：和S₀＝[y₀ θ₀ M₀ Q₀ 1]^T分别为整个桩的桩端和桩顶的位移、转角、弯矩和剪力状态量。U_(i,j)为自由段矩阵传递系数或弹性段传递矩阵系数或塑性段传递矩阵系数，根据具体的桩土相互作用状态以及是否存在自由段桩进行决定。U为整个桩身的总传递矩阵系数，是从整个桩的桩端到桩顶的每小段桩的传递矩阵系数的连续乘积。根据沿深度的桩土相互作用状态(如图6)以及是否存在自由段桩，总传递矩阵系数U求解过程如下：

①如果地面以上存在自由段桩时，则自由段桩身总传递矩阵系数U_a如下式：

则总传递矩阵系数U根据沿深度的桩土相互作用状态(见图6)按照下式进行计算：

②如果地面以上不存在自由段桩时，则总传递矩阵系数U根据沿深度的桩土相互作用状态(见图6)按照下式进行计算：

联合桩顶、桩端的边界条件和总传递矩阵系数U，即可解得桩顶变形、转角、弯矩和剪力如下所示：

桩顶自由-桩端自由：

式中：U_ij表示总传递矩阵系数U(5×5阶)中的第i行第j列元素。

桩顶自由-桩端嵌固：

桩顶嵌固-桩端自由：

桩顶嵌固-桩端嵌固：

在桩顶状态量S₀确定之后，进而按照下式计算任意位置处的桩身响应：

(8)根据传递矩阵原理，代入桩顶与桩端的边界条件进行建立桩身响应的传递矩阵表达式，并编制相应的核心计算程序进行桩身变形与内力的求解，该程序计算收敛关键步骤如下(如图7)：

①根据计算模型确定相应的基本计算参数并进行离散化与常数化处理，自由段以及地面下桩身等分数量m_a和m_i确定原则为其等分后每小段的桩身长度要小于等于相应位置处桩身直径的二十分之一；对于桩土相互作用模型为非线性模型时，取该p-y曲线模型中y＝0.01mm所对应的割线刚度作为初始土抗力模量k；在计算初始阶段需先假定桩土相互作用均为弹性状态；

②根据上一次迭代的计算结果，如果存在非线性的桩土相互作用p-y曲线，则需要判断经过r次迭代后割线刚度是否满足公式的要求，如果满足割线刚度迭代精度要求，则判断经过s次迭代后整体精度是否满足公式的要求；其中，表示经过第r次迭代后的第i部分桩中第j小段桩的桩土相互作用p-y曲线割线刚度；ε表示迭代收敛标准推荐值；表示经过s次迭代后从地表往桩身方向发展的塑性区深度；如果满足要求则进行桩身响应最终状态求解，如果不满足上述2种迭代要求中的任意一种，则需要进行再次迭代直到两种状态的迭代精度均满足要求。

以上所述，仅为本发明较佳的具体实施方式，但本发明的具体保护范围并不局限于此，任何熟悉本技术领域的技术人员在本发明揭露的技术范围内，可轻易得到变换和替换，都应涵盖在本发明的保护范围之内。因此，本发明的保护范围应该以权利要求的保护范围为准。

实施例

1、统一p-y曲线模型正确性验证

鉴于不考虑竖向荷载作用时的水平受荷桩承载力计算研究文献相比组合受荷桩要全面，因此，基于组合受荷桩不考虑轴向荷载作用的特例情况，对p-y曲线统一表达式正确性验证。

①线弹性模型验证

一根等直变截面桩嵌入砂层土和软岩中。砂层厚度H₁＝7.925m，岩层厚度H₂＝12.192m。第一层砂土的土抗力模量为k₁＝1.04×10⁴z kN/m²(z从地面开始计算)；第二层软岩的土抗力模量分别为k₂＝5.1386×10⁴+5.91×10³(z-7.925)kN/m²(根据经验公式反算得出)和k₂＝4.9573×10⁴+7.419×10⁴(z-7.925)kN/m²(根据p-y曲线反算得出)。位于砂层中基桩直径(d₁)和抗弯刚度(EI₁)分别为2.59m和9.9324×10⁷kN·m²；位于软岩层中基桩直径(d₂)和抗弯刚度(EI₂)分别为2.44m和7.8327×10⁷kN·m²。桩顶和桩端边界均为自由边界。国外某学者采用变分法提出了该问题的解析解。计算结果如图8所示，可见本文所提出的基于Laplace正逆变换的传递矩阵解与国外学者基于变分法解析解的结果几乎完全一致，对比结果证明了本发明方法推导的正确性，也验证了本发明提出的p-y曲线统一表达式适用于线弹性桩土相互作用模型。

②非线弹性模型验证

国外某学者开展了一根长45m、直径1.83m的桩垂直嵌入8层黏性地基土中的水平承载特性的研究。该桩实际嵌入长度为43m(2m位于地面以上)，土体的弹性模量E_s＝500S_u，桩的弹性模量E_p＝30GPa，该学者采用有限差分法对该水平受荷进行计算，分析过程中桩土相互作用的p-y曲线采用双曲线型p-y曲线，如下式所示：

式中，k_s＝0.943(z/z_ref)^0.016(d/d_ref)v_s^-0.078E_s^1.036E_p^-0.031，p_u取(3+γz/S_u+Jz/d)S_u或10S_u两者中的最小值。桩顶水平的水平荷载F_t＝3000kN位于地表0.6m位置处，桩顶和桩端边界条件均为自由边界，土层参数详见表1所示。

表1双曲线模型基本计算参数

如图9所示为本文基于p-y曲线统一表达式的传递矩阵半解析解与国外学者基于p-y双曲线的有限差分解对比。从图中可以看出本发明通解计算结果与国外学者的结果几乎一致，对比结果证明了本发明方法推导的正确性，也验证了本发明提出的p-y曲线统一表达式适用于非线性弹性桩土相互作用模型。

③线弹性-塑性模型验证

国外某学者报道了一个在砂土中桩长12m、桩径0.5m的水平受荷桩试验，桩身抗弯EI_p＝1.02×10⁵kN·m²，砂土重度为10kN/m³，摩擦角为30度。该学者采用变分法求解了水平受荷桩桩身响应，计算过程中采用线弹性-塑性模型，该p-y模型的弹性段的土抗力模量k＝10000kN/m²，极限土抗力采用p_u＝3γzK_pd计算，其中桩顶施加水平荷载为F_t＝200kN，桩顶和桩端边界条件均为自由边界。则桩身变形对比如图10所示，本发明方法解得的桩身变形与国外学者解几乎吻合，对比结果证明了本发明方法推导的正确性，也验证了本发明提出的p-y曲线统一表达式适用于线弹性-塑性桩土相互作用模型。

④非线弹性-塑性模型验证

国外某学者报道了一个在双层黏性地基土中桩长26.6m、桩径1.02m的钢管桩水平受荷试验。钢管桩壁厚为16mm，相应的桩身抗弯刚度为1.26×10⁶kN·m²。桩真实嵌入深度为25.6m(地面上长度为1m)。计算过程中，该学者采用经典的Matlock粘土p-y曲线模型模拟桩土相互作用，即API粘土p-y曲线模型，如下所示：

式中，y₅₀＝ρε₅₀d，ρ为相关系数，取2.5；ε₅₀为三轴仪试验中最大主应力差一半时的应变值，计算参数如表2所示。

表2API粘土p-y曲线模型基本参数

桩身变形和弯矩对比如图11所示，本发明方法所得桩身变形与弯矩结果与国外学者所得结果几乎吻合，这证明了本文方法推导的正确性，也验证了本文提出的p-y曲线统一表达式适用于线弹性-塑性桩土相互作用模型。

2、组合受荷桩承载力传递矩阵通解正确性验证

国内某学者开展了组合荷载作用下单桩承载特性试验研究。桩采用铝合金材质的钢管桩，桩自重γ＝27.1kN/m³，桩外径为d＝0.016m，钢管桩壁厚为t＝0.002m；桩总长度为H＝0.8m，桩身抗弯刚度为EI＝0.189kN·m²。试验中，土体采用均质的砂性土，砂层表面以上存在一定长度的自由段桩，剩余部分桩全部嵌入均质的砂土中。由于是钢管桩且壁厚较薄，因此假定桩身轴力为常数分布。采用如下的桩土相互作用表达方式：

其中，n_hs为水平土抗力系数，单位为kN/m³；z₀为等效深度，单位为m；b为深度的指数，单位为一；y_u(＝3d/80)为土体屈服位移，单位为m。

由于该学者开展试验时所提供的参数不足以确定公式(54)中所有参数，因此根据图12(a)所示的第一组试验结果进行反算，得出n_hs＝5.84×103kN/m³、z₀＝0.3m和b＝0.8，其中y_u＝3d/80＝0.6×10^-3m。根据反算所得参数计算剩余三组试验，所得桩身弯矩结果如图12所示。通过图12对比可以发现本发明所提出的组合受荷桩承载力传递矩阵通解结果与实测结果非常接近，且桩身弯矩最大值吻合非常好，这说明了本发明提成的一种基于Laplace正逆变换的组合受荷桩承载力传递矩阵通解算法的正确性以及计算的合理性。

以上所述仅是本发明的优选实施方式，应当指出，对于本技术领域的普通技术人员来说，在不脱离本发明技术原理的前提下，还可以做出若干改进和变形，这些改进和变形也应视为本发明的保护范围。