移动简谐荷载作用下层状道路结构的安定下限分析

林缘祥; 郑俊杰; 后如意; 方昊

doi:10.11779/CJGE202211008

移动简谐荷载作用下层状道路结构的安定下限分析 English Version

华中科技大学岩土与地下工程研究所, 湖北武汉 430074

基金项目:

国家自然科学基金项目 51878313

国家自然科学基金项目 52078236

详细信息

作者简介:
林缘祥(1996—)，男，博士研究生，主要从事交通岩土工程方面的研究。E-mail：ce_linyx@hust.edu.cn

通讯作者:
郑俊杰，E-mail: zhengjj@hust.edu.cn

中图分类号: TU43
计量
- 文章访问数: 132
- HTML全文浏览量: 15
- PDF下载量: 20
出版历程
- 收稿日期: 2021-09-06
- 网络出版日期: 2022-12-08
- 刊出日期: 2022-10-31

Lower shakedown limits of layered road structures under moving harmonic loads

Institute of geotechnical and underground engineering, Huazhong University of Science and Technology, Wuhan 430074, China

摘要

摘要: 为了研究移动简谐荷载作用下层状道路结构的安定性问题，先通过傅里叶变换方法和数值积分求得三维层状道路结构在时间–空间域内的动力响应，然后考虑饱和土层的有效应力场而不是总应力场，对现有的静力安定理论和安定极限求解方法进行改进，提出了有效安定极限的概念，并与考虑总应力的安定极限求解方法进行了对比分析。此外，针对饱和土层选取不同的有效内摩擦角，分别研究了荷载移动速度、荷载频率以及道路面层刚度对层状道路结构的有效安定极限和有效临界深度的影响规律。结果表明：有效安定极限与考虑总应力的安定极限求解方法得到的安定极限有明显差异；而且在荷载移动速度较大时，有效临界深度比考虑总应力的安定求解方法得到的临界深度更深。有效安定极限的求解方法更适用于包含饱和土层的层状道路结构设计和安全评估。
- 安定性分析 /
- 层状道路结构 /
- 饱和土 /
- 静力安定定理 /
- 有效应力
Abstract: The shakedown limits of layered road structures subjected to a moving harmonic load are studied. The inverse Fourier transform and the numerical integration are used to obtain the dynamic responses of a three-dimensional layered road structure in the time and space domain. Considering the effective stress field of saturated subsoil instead of the total stress field, the existing static shakedown theorem and the solving method for the shakedown limits are improved, and the concept of the effective shakedown limit is proposed and compared with the solving method for the shakedown limits considering the total stress. In addition, different effective internal friction angles are selected for the saturated soil layer, and the influences of load-moving speed, load frequency and pavement stiffness on the effective shakedown limit and effective critical depth of the layered road structure are studied respectively. The results show that there is a significant difference between the effective shakedown limits and the shakedown limits obtained by the solving method for the shakedown limits considering the total stress. Moreover, the effective critical depth is deeper than that obtained by the solving method for the shakedown limits considering the total stress when the load-moving speed is relatively high. The proposed method for solving the effective shakedown limits is more suitable for the design and safety assessment of layered road structures containing saturated soil layers.
- shakedown analysis /
- layered road structure /
- saturated soil /
- static shakedown theorem /
- effective stress

HTML全文

0. 引言

中国地处环太平洋地震带以及欧亚地震带之间，其地震活动呈现出高频性、高强度、浅震源等特点，当地震发生在人口稠密的地区时，会带来惨重的经济和社会损失。对工程进行抗震研究，减轻地震灾害的影响将变得十分重要。

在对工程进行地震响应分析过程中，考虑到现实情况的复杂性，常采用数值方法进行计算^[1]。动力时程分析法可以较好地揭示结构在地震全时段内的响应规律^[2]，随着计算机硬件设备的不断发展，该方法得到了一定的认可。应用有限元方法进行动力时程分析过程中，为降低计算成本，提高计算效率，往往要将地层取出规定大小的有限区域来进行计算，并在此区域外施加人工边界来模拟无限地基的情况^[3]。其中人工边界包括无限元边界、黏性边界和黏弹性边界等形式，而黏弹性人工边界因具有好的低频和高频稳定性，应用方便，被广泛使用^[4]。

在进行工程抗震分析时，根据波动理论，考虑到地壳密度随地层深度的增加而增加，普遍认为地震波在地表附近的传播方向接近竖向^[5]，故而常常探讨地震波在垂直入射下情况。但现实情况下，根据地层条件的不同，地震波往往存在一定的入射角度。Jin等^[6]和Takahiro等^[7]通过对各地震动记录进行反演分析，得出了不同的地震波入射角度。因此，若仅考虑垂直入射的情况对结构进行地震响应分析，或将带来不能忽视的误差^[8]。因此，考虑地震波斜入射也具有一定的现实意义。而在此方面，许多学者^[9-11]同样进行了很多相关的研究。

而对于水平层状场地地震波斜入射问题，刘晶波等^[12-13]提出了一维化时域有限元计算方法；赵密等^[14-15]和尹侯权等^[16]利用弹性介质的应力–位移本构关系建立了人工边界处应力与速度的阻抗边界条件，替代黏性边界条件，改进P-SV波斜入射时成层半空间自由场的时域算法；王笃国等^[17-18]应用频域传递矩阵结合等效线性化方法，求解了地震波斜入射下水平层状场地的地震响应。而直接运用黏弹性边界理论，因时间延迟计算较为繁琐，缺乏深入研究。

对于地震波斜入射中时间延迟的计算，目前普遍采用立体几何的方式对地震波路径逐一进行求解，或采用坐标变换^[19]的方式对计算过程进行一定成度的简化。然而，当面对水平层状场地等特殊情况时，存在地震波多次反射和透射的情况，此方法不易迭代，计算成本较高。

对此，本文尝试建立一种新的时间延迟计算方法，以当前地震波属性表征其生成路径全过程上的时间延迟；在此基础上，推导黏弹性边界下水平层状场地地震波斜入射方法，进而实现合理有效的地震动数值计算。最后建立有限元计算模型，利用此方法对单脉冲作用下的数值解与理论解进行对比，验证方法的准确性和适用性。

1. 黏弹性边界地震动输入

黏弹性人工边界的思想就是将弹簧和黏滞阻尼器并联，形成一种弹簧-阻尼器物理元件，并将其施加在边界上。要使用这种边界，首先要将地震波转化成人工边界上各处所对应的等效荷载，其中包括产生内行场反应所需要抵抗的人工边界物理元件上的结点力，以及因产生内行场反应而需抵抗近场介质所需的结点力^[20-21]。当黏弹性边界可以将计算区域内向外传播的散射波完全吸收时，人工边界的节点上便会形成因地震作用而产生的自由场运动。

这样，地震动的输入问题便转化为在人工边界的节点上所作用的自由场运动问题，通等效节点力的方式施加在所求模型的边界上。

等效节点力的计算公式为

${F_{{\text{b}}}} = ({{{\boldsymbol{ K}}}_{{\text{b}}}}{{\boldsymbol{u }}_{{\text{b}}}} + {{{\boldsymbol{C }}}_{{\text{b}}}}{\dot {\boldsymbol{ u}}_{{\text{b}}}} + {\sigma _{{\text{b}}}}n){A_{{\text{l}}}} 。$

(1)

在弹簧刚度矩阵 ${{\boldsymbol{K }}_{\text{b}}}$ 和阻尼系数矩阵 ${{\boldsymbol{ C}}_{\text{b}}}$ 的取值方面，许多学者都进行了研究，同时给出了一定的建议，但不同方法对结果影响不大。对此，本文选用文献[22]的方法进行计算。

在水平层状场地的计算过程中，式(1)的自由场位移向量 ${{\boldsymbol{u }}_{{\text{b}}}}$ 和自由场速度向量 ${\dot {\boldsymbol{ u}}_{{\text{b}}}}$ 需要根据斯涅尔方程得到地震波在成层界面上的转换情况，并分析其与入射波之间的时空关系，进而利用波的叠加原理进行求解，这其中便涉及到时间延迟的计算问题。在此方面，当前的计算方法较为复杂，求解较为困难。下面针对水平层状场地时间延迟的计算进行分析。

2. 时间延迟计算方法

针对水平层状场地，本文假定各层地层均匀弹性，在水平方向上无限延伸。根据地震波波动理论可知，平面P波斜入射时，在地层界面处生成反射P波、反射SV波、透射P波和透射SV波；平面SV波斜入射时，如入射角度不超过临界角，将生成反射P波、反射SV波、透射P波和透射SV波。因此，为将时间延迟计算公式进行统一，在这里假设入射波为m，反射波n，透射波l，其中m，n，l为平面P波或平面SV波，入射波、反射波及透射波的方向角分别为 $\alpha$ ， $\beta$ ， $\gamma$ ， ${T_0}$ 时刻入射波波前与入射波传播方向垂直，仅区分反射和透射两种情况进行讨论。

对反射波，如所示， $\alpha$ 和 $\beta$ 满足：

$\frac{{{{\text{sin}}}\alpha }}{{{v_1}}} = \frac{{{{\text{sin}}}\beta }}{{{v_2}}}{\text{，}}$

(2)

图 1 反射波时间延迟计算示意图

Figure 1. Time delay calculation of reflected waves

下载: 全尺寸图片幻灯片

式中， ${v_1}$ ， ${v_2}$ 分别为入射波m和反射波n在对应地层内的传播速度。

假设 ${X_{{\text{f}}}}({x_0}, {y_0})$ 为所求模型边界上的一点，AB为零时刻的波前面，则反射波n在该点处的时间延迟为

$\Delta t = {t_1} + {t_2} {\text{，}}$

(3)

${t_1} = \frac{{OD}}{{{v_1}}}{\text{，}}$

(4)

${t_2} = \frac{{O{X_{{\text{f}}}}}}{{{v_2}}}{\text{，}}$

(5)

式中， ${t_1}$ ， ${t_2}$ 分别为计算所需入射波m和反射波n的时间延迟。

在此，过A点做ON垂线，与ON延长线交于点C。由几何关系得

$OD = AO\sin \alpha {\text{，}}$

(6)

$OC = AO\sin \beta {\text{，}}$

(7)

${t_1} = \frac{{OD}}{{{v_1}}} = \frac{{OC}}{{{v_2}}}{\text{，}}$

(8)

$\Delta t = {t_1} + {t_2} = \frac{{OC + O{X_{{\text{f}}}}}}{{{v_2}}} = \frac{{C{X_{{\text{f}}}}}}{{{v_2}}}。$

(9)

设点 $U({x_1}, {y_1}) \in AC$ ，则式（9）可改写为

$\Delta t = \frac{1}{{{v_2}}}\left[ {({x_1} - {x_0})\cos \beta + ({y_1} - {y_0})\sin \beta } \right]。$

(10)

对透射波，如所示， $\alpha$ 和 $\gamma$ 满足：

$\frac{{{{\text{sin}}}\alpha }}{{{v_1}}} = \frac{{{{\text{sin}}}\gamma }}{{{v_3}}}{\text{，}}$

(11)

图 2 透射波时间延迟计算示意图

Figure 2. Time delay calculation of transmission waves

下载: 全尺寸图片幻灯片

式中， ${v_3}$ 为透射波l在对应地层内的传播速度。

假设 ${X_t}({x_2}, {y_2})$ 为所求模型边界上的一点，AB为零时刻的波前面，则透射波l在该点处的时间延迟为

$\Delta t = {t_1} + {t_3}{\text{，}}$

(12)

${t_1} = \frac{{OD}}{{{v_1}}}{\text{，}}$

(13)

${t_3} = \frac{{O{X_{{\text{t}}}}}}{{{v_3}}}{\text{，}}$

(14)

式中，t₃为计算所需透射波的时间延迟。

在此，过A点做ON垂线，与ON延长线交于点 $C'$ 。由几何关系得

$OD = AO\sin \alpha {\text{，}}$

(15)

$OC' = AO\sin \gamma {\text{，}}$

(16)

${t_1} = \frac{{OD}}{{{v_1}}} = \frac{{OC'}}{{{v_3}}}{\text{，}}$

(17)

$\Delta t = {t_1} + {t_3} = \frac{{OC' + O{X_{{\text{t}}}}}}{{{v_3}}} = \frac{{C{X_{{\text{t}}}}}}{{{v_3}}}。$

(18)

设点 $V({x_3}, {y_3}) \in AC'$ ，式（18）可改写为

$\Delta t = \frac{1}{{{v_3}}}(({x_3} - {x_2})\cos \gamma + ({y_3} - {y_2})\sin \gamma )。$

(19)

至此，反射波n和透射波l的时间延迟计算公式中已经不含有其对应入射波m的任何参数，所有参数均为当前地震波的参数，且形式较为统一，便于进行迭代计算。

由于点 $U({x_1}, {y_1}) \in AC$ ，点 $V{{\text{(}}}{x_3}, {y_3}{{\text{)}}} \in AC'$ ，不妨将其均取为波前与地层界面相交的A点，作为所求地震波的等效零时刻点，由此便可以得到时间延迟的计算公式为

$\Delta {t_i} = \frac{1}{{{v_i}}}{{\text{(}}}\Delta {x_i}\cos {\theta _i} + \Delta {y_i}\sin {\theta _i}{{\text{)}}}{\text{，}}$

(20)

式中， $\Delta {x_i}$ ， $\Delta {y_i}$ 为所求点与A点坐标差， ${\theta _i}$ 为所求地震波传播方向与竖向的夹角， ${v_i}$ 为所求地震波在对应地层中的传播速度。

将其引申至三维模型，如图 3所示，波前与地层界面相交于A₁A₂，故A点可为A₁A₂上的任意一点，为便于计算可将其设为特殊点。此时时间延迟计算公式：

$\Delta {t_i} = \frac{1}{{{v_i}}}\left[ {\begin{array}{*{20}{c}} {\Delta {x_i}}&{\Delta {y_i}}&{\Delta {z_i}} \end{array}} \right]\left[ {\begin{array}{*{20}{c}} {\cos {\alpha _i}} \\ {\cos {\beta _i}} \\ {\cos {\gamma _i}} \end{array}} \right]{\text{，}}$

(21)

图 3 三维模型时间延迟计算示意图

Figure 3. Time delay calculation of 3D model

下载: 全尺寸图片幻灯片

式中， ${\alpha _i}$ ， ${\beta _i}$ ， ${\gamma _i}$ 分别为所求地震波传播方向与3个坐标轴的夹角。

3. 水平层状场地地震波斜入射

对于水平层状场地而言，地震波斜入射边界条件的计算重点在于等效节点力的求解上，而地震波在成层界面的反射透射导致其计算过程较为困难。因此，通过将上述时间延迟计算思路应用于水平层状场地的情况，采用迭代方式实现水平层状场地地震动的计算问题。

利用文献[]中地震波在成层界面上的对应关系，在计算过程中，对任意地震波，可通过将其产生过程中每一次转换所得到的振幅比值进行连乘来记录其与原入射波之间的幅度比值，即 ${E_i} = \prod {A{}_i{B_i}{C_i}{D_i}}$ 。其中A，B，C，D分别表示同类反射波、转换反射波、同类透射波、转换透射波与入射波幅值的比值。同时，由于地震波遇到成层界面时将会发生反射和透射现象，故在实际计算过程中存在终止条件的选择问题。在此取计算精度 $\delta$ ，忽略所得地震波与原入射波的比值小于 $\delta$ 的分量（即 ${E_i} < \delta$ ）及其后的地震波。考虑到地震是扩散衰减的，应用此方法进行计算无需担心收敛问题。 $\delta$ 的大小将在一定程度上影响计算精度和计算时效，故应在对二者进行权衡，但与动力有限元计算相比，此处计算量占比较小，无需过多斟酌。对地震波求解完成后，即可进一步计算黏弹性边界上各个节点的节点力，并施加在模型上。计算流程如图 4所示。

图 4 计算流程图

Figure 4. Flow chart of calculation

下载: 全尺寸图片幻灯片

以第k层地层边界面上的一点 $N({x_n}, {y_n}, {z_n})$ 为例，求解黏弹性边界水平层状场地地震波斜入射下的等效节点力。

假设Ⅰ₁为该层地层中与Z轴正向夹角大于90°的平面P波，Ⅰ₂为该层地层中与Z轴正向夹角小于90°的平面P波，Ⅰ₃为该层地层中与Z轴正向夹角大于90°的平面SV波，Ⅰ₄为该层地层中与Z轴正向夹角小于90°的平面SV波。计算过程中，应注意地震波传播过程中正方向的连续性，以保证公式的准确性。

假设 ${u_0}$ 为基岩入射波的位移时程，则N点处自由场位移：

${u_{Nx}} = \mathop \sum \limits_{i \in {{{\text{I}}}_1}, {{{\text{I}}}_2}} {E_i}{u_0}(t - \Delta {t_i}){{\text{sin}}}{\alpha _i} + \mathop \sum \limits_{i \in {{{\text{I}}}_3}, {{{\text{I}}}_4}} {E_i}{u_0}(t - \Delta {t_i}){{\text{cos}}}{\alpha _i} {\text{，}}$

(22)

${u_{Ny}} = 0{\text{，}}$

(23)

${u_{Nz}} = - \mathop \sum \limits_{i \in {{{\text{I}}}_1}} {E_i}{u_0}(t - \Delta {t_i}){{\text{cos}}}{\alpha _i} + \mathop \sum \limits_{i \in {{{\text{I}}}_2}} {E_i}{u_0}(t - \Delta {t_i}){{\text{cos}}}{\alpha _i} +$

$\sum\limits_{i\in {{\text{I}}}_{3}}{E}_{i}{u}_{0}(t-\Delta {t}_{i}){\text{sin}}{\alpha }_{i}-\sum\limits_{i\in {{\text{I}}}_{4}}{E}_{i}{u}_{0}(t-\Delta {t}_{i}){\text{sin}}{\alpha }_{i}{\text{ }}。$

(24)

将式（22）~（24）对时间求导可以得到N点处自由场速度：

${\dot u_{Nx}} = \mathop \sum \limits_{i \in {{{\text{I}}}_1}, {{{\text{I}}}_2}} {E_i}{\dot u_0}(t - \Delta {t_i}){{\text{sin}}}{\alpha _i} + \mathop \sum \limits_{i \in {{{\text{I}}}_3}, {{{\text{I}}}_4}} {E_i}{\dot u_0}(t - \Delta {t_i}){{\text{cos}}}{\alpha _i}{\text{，}}$

(25)

${\dot u_{Ny}} = 0{\text{，}}$

(26)

${\dot u_{Nz}} = - \mathop \sum \limits_{i \in {{{\text{I}}}_1}} {E_i}{\dot u_0}(t - \Delta {t_i}){{\text{cos}}}{\alpha _i} + \mathop \sum \limits_{i \in {{{\text{I}}}_2}} {E_i}{\dot u_0}(t - \Delta {t_i}){{\text{cos}}}{\alpha _i} +$

$\sum \limits_{i\in {{\text{I}}}_{3}}{E}_{i}{\dot{u}}_{0}(t-\Delta {t}_{i}){\text{sin}}{\alpha }_{i}-\sum \limits_{i\in {{\text{I}}}_{4}}{E}_{i}{\dot{u}}_{0}(t-\Delta {t}_{i}){\text{sin}}{\alpha }_{i}。$

(27)

在计算自由场应力时，可以根据波的叠加效应分别对模型的每个侧面及底面进行分析。

当对模型任意一点 $(x, y, z)$ ，都有x_n≥x时，N点处自由场应力：

${\sigma _{Nx}} = - \mathop \sum \limits_{i \in {{{\text{I}}}_1}, {{{\text{I}}}_2}} \frac{{{\lambda _k} + 2{G_k}{{\text{si}}}{{{\text{n}}}^2}{\alpha _i}}}{{{c_{pk}}}}{E_i}{\dot u_0}(t - \Delta {t_i}) -$

$\sum \limits_{i\in {{\text{I}}}_{1}, {{\text{I}}}_{2}}\frac{{G}_{k}{\text{sin}}2{\alpha }_{i}}{{c}_{sk}}{E}_{i}{\dot{u}}_{0}(t-\Delta {t}_{i}){\text{ }}{\text{，}}$

(28)

${\sigma _{Ny}} = 0{\text{，}}$

(29)

${\sigma _{Nz}} = \mathop \sum \limits_{i \in {{{\text{I}}}_1}} \frac{{{G_k}{{\text{sin}}}2{\alpha _i}}}{{{c_{pk}}}}{E_i}{\dot u_0}(t - \Delta {t_i}) - \mathop \sum \limits_{i \in {{{\text{I}}}_2}} \frac{{{G_k}{{\text{sin}}}2{\alpha _i}}}{{{c_{pk}}}} \cdot {E_i}{\dot u_0}(t - \Delta {t_i}) + \mathop \sum \limits_{i \in {{{\text{I}}}_3}} \frac{{{G_k}{{\text{cos}}}2{\alpha _i}}}{{{c_{sk}}}}{E_i}{\ddot u_0}(t - \Delta {t_i}) -$

$\sum \limits_{i\in {{\text{I}}}_{4}}\frac{{G}_{k}{\text{cos}}2{\alpha }_{i}}{{c}_{sk}}{E}_{i}{\dot{u}}_{0}(t-\Delta {t}_{i}){\text{ }}。$

(30)

当对模型任意一点 $(x, y, z)$ ，都有x_n≤x时，N点处自由场应力：

${\sigma _{Nx}} = \mathop \sum \limits_{i \in {{{\text{I}}}_1}, {{{\text{I}}}_2}} \frac{{{\lambda _k} + 2{G_k}{{\text{si}}}{{{\text{n}}}^2}{\alpha _i}}}{{{c_{pk}}}}{E_i}{\dot u_0}(t - \Delta {t_i}) +$

$\sum \limits_{i\in {{\text{I}}}_{1}, {{\text{I}}}_{2}}\frac{{G}_{k}{\text{sin}}2{\alpha }_{i}}{{c}_{sk}}{E}_{i}{\dot{u}}_{0}(t-\Delta {t}_{i}){\text{ }}{\text{，}}$

(31)

${\sigma _{Ny}} = 0{\text{，}}$

(32)

${\sigma _{Nz}} = - \mathop \sum \limits_{i \in {{{\text{I}}}_1}} \frac{{{G_k}{{\text{sin}}}2{\alpha _i}}}{{{c_{pk}}}}{E_i}{\dot u_0}(t - \Delta {t_i}) + \mathop \sum \limits_{i \in {{{\text{I}}}_2}} \frac{{{G_k}{{\text{sin}}}2{\alpha _i}}}{{{c_{pk}}}} \cdot$

${E_i}{\dot u_0}(t - \Delta {t_i}) - \mathop \sum \limits_{i \in {{{\text{I}}}_3}} \frac{{{G_k}{{\text{cos}}}2{\alpha _i}}}{{{c_{sk}}}}{E_i}{\dot u_0}(t - \Delta {t_i}) +$

$\sum \limits_{i\in {{\text{I}}}_{4}}\frac{{G}_{k}{\text{cos}}2{\alpha }_{i}}{{c}_{sk}}{E}_{i}{\dot{u}}_{0}(t-\Delta {t}_{i}){\text{ }}。$

(33)

当对模型任意一点 $(x, y, z)$ ，都有y_n≥y时，N点处自由场应力：

${\sigma _{Nx}} = 0{\text{，}}$

(34)

${\sigma _{Ny}} = - \mathop \sum \limits_{i \in {{{\text{I}}}_1}, {{{\text{I}}}_2}} \frac{{{\lambda _k}}}{{{c_{pk}}}}{E_i}{\dot u_0}(t - \Delta {t_i}){\text{，}}$

(35)

${\sigma _{Nz}} = 0。$

(36)

当对模型任意一点 $(x, y, z)$ ，都有y_n≤y时，N点处自由场应力：

${\sigma _{Nx}} = 0{\text{，}}$

(37)

${\sigma _{Ny}} = \mathop \sum \limits_{i \in {{{\text{I}}}_1}, {{{\text{I}}}_2}} \frac{{{\lambda _k}}}{{{c_{pk}}}}{E_i}{\dot u_0}{{\text{(}}}t - \Delta {t_i}{{\text{)}}}{\text{，}}$

(38)

${\sigma _{Nz}} = 0。$

(39)

当对模型任意一点 $(x, y, z)$ ，都有z_n≤z时，N点处自由场应力：

${\sigma _{Nx}} = - \mathop {\sum }\limits_{i \in {{{\text{I}}}_1}} \frac{{{G_k}{{\text{sin}}}2{\alpha _i}}}{{{c_{pk}}}}{E_i}{\dot u_0}(t - \Delta {t_i}) + \mathop \sum \limits_{i \in {{{\text{I}}}_2}} \frac{{{G_k}{{\text{sin}}}2{\alpha _i}}}{{{c_{pk}}}} \cdot$

${E_i}{\dot u_0}(t - \Delta {t_i}) - \mathop \sum \limits_{i \in {{{\text{I}}}_3}} \frac{{{G_k}{{\text{cos}}}2{\alpha _i}}}{{{c_{sk}}}}{E_i}{\dot u_0}(t - \Delta {t_i}) + {{\text{ }}}$

$\sum \limits_{i\in {{\text{I}}}_{4}}\frac{{G}_{k}{\text{cos}}2{\alpha }_{i}}{{c}_{sk}}{E}_{i}{\dot{u}}_{0}(t-\Delta {t}_{i}){\text{ }}{\text{，}}$

(40)

${\sigma _{Ny}} = 0{\text{，}}$

(41)

${\sigma _{Nz}} = \mathop \sum \limits_{i \in {{{\text{I}}}_1}, {{{\text{I}}}_2}} \frac{{{\lambda _k} + 2{G_k}{{\text{co}}}{{{\text{s}}}^2}{\alpha _i}}}{{{c_{pk}}}}{E_i}{\dot u_0}(t - \Delta {t_i}) -$

$\sum \limits_{i\in {{\text{I}}}_{3}, {{\text{I}}}_{4}}\frac{{G}_{k}{\text{sin}}2{\alpha }_{i}}{{c}_{sk}}{E}_{i}{\dot{u}}_{0}(t-\Delta {t}_{i}){\text{ }}。$

(42)

式（22）~（42）中的 $\Delta {t_i}$ 可根据式（21）进行计算。将式（22）~（42）代入等效节点力的计算公式，即可求得黏弹性边界下分层地层地震波斜入射过程中施加在模型边界节点上的等效节点力。在此基础上，依托现有的黏弹性边界地震动输入理论，利用常用的有限元计算软件和计算方法，即可对水平层状场地进行地震动有限元分析。

4. 方法验证

为验证本文方法的有效性和准确性，使用相关有限元计算软件建立如所示的计算模型。考虑两层水平层状场地结构进行模拟，模型上层 ${H_1} = 100\;{{\text{ m}}}$ ，下层 ${H_2} = 150\;{{\text{ m}}}$ 。长度L和宽度D均取250 m。采用T=0.25 s的狄拉克脉冲函数作为入射波，其位移的时程关系曲线如所示，先考虑模型两层密度 ${\rho _1} = {\rho _2} =$ 2.0×10³ ${{\text{kg/}}}{{{\text{m}}}^{{\text{3}}}}$ ，弹性模量 ${E_1} = {E_2} = 1{{\text{ GPa}}}$ ，泊松比 ${\mu _1} =$ ${\mu _2} = 0.2$ 的特殊情况，取入射角 $\alpha$ 为10°，20°平面P波和 $\alpha$ 为5°，15°平面SV波作为入射波，输入的时间间隔取0.005 s，对比本文方法与文献[19]方法的异同。

图 5 计算模型示意图

Figure 5. Computational model

下载: 全尺寸图片幻灯片

图 6 输入位移时程曲线

Figure 6. Time histories of input displacement

下载: 全尺寸图片幻灯片

分别采用两种方法得到有限元计算输入文件，图 7给出了文件的MD5计算结果，两种方法得到的MD5完全相同，说明文件一致，因此本文方法在有限次终止时与文献[19]方法的结果保持一致。

图 7 MD5计算结果

Figure 7. Results of MD5

下载: 全尺寸图片幻灯片

考虑一般情况，取上层密度 ${\rho _1} =$ 2.0×10³ ${{\text{kg/}}}{{{\text{m}}}^{{\text{3}}}}$ ，弹性模量 ${E_1} = 1\;{{\text{ GPa}}}$ ，泊松比 ${\mu _1} = 0.2$ ；下层密度 ${\rho _2} =$ 2.0×10³ ${{\text{kg/}}}{{{\text{m}}}^{{\text{3}}}}$ ，泊松比 ${\mu _2} = 0.2$ ，弹性模量 ${E_2} = 1.5$ GPa，分别取入射角 $\alpha$ 为10°，20°平面P波和 $\alpha$ 为5°，15°平面SV波作为入射波。

取 $\delta = 0.01$ ，图 8给出了模型在不同角度及不同类型入射波下地表中点处的水平（X方向）和竖向（Z方向）的位移时程，并与理论解进行对比。理论解采用频域传递矩阵法^[24]结合快速傅里叶变换得到^[14]。从图中可以看出，模型解与理论解吻合较好，相对误差较小，说明此方法在求解水平层状场地地震波斜入射方面具有良好的精度和准确性。

图 8 地表中点处的位移时程

Figure 8. Time histories of displacement at midpoint of surface

下载: 全尺寸图片幻灯片

以入射角 $\alpha =$ 20°的平面P波为例给出了模型在不同时刻下的位移场云图，从图中可以看出，本文的方法所引起地层内部变形规律与边界一致，变化较为均匀合理，且云图中的位移角度与地震波的斜入射的理论值基本上是吻合的，验证了本方法在计算水平层状场地地震波斜入射过程中的可靠性。

图 9 不同时刻位移场云图

Figure 9. Cloud diagram of displacement field at different time

下载: 全尺寸图片幻灯片

5. 结论

（1）提出了新的时间延迟计算方法。使用解析方式推导时间延迟的计算公式，使公式更加简化。在进行时间延迟的求解过程中，可通过所求目标地震波的具体信息及其所在地层的具体条件直接求解，摒弃以往计算方式中对所求地震波的全部来源依次进行计算的方式，使时间延迟公式统一，增强时间延迟计算方法的适用性。

（2）推导了基于等效零时刻点迭代的水平层状场地地震波斜入射方法。对于黏弹性边界下水平层状场地地震波波动输入问题，推导了基于本文时间延迟计算方法的人工边界上等效节点力的计算公式，进而得到了黏弹性边界下水平层状场地地震动输入的计算公式，拓展了黏弹性边界理论在水平层状场地斜入射地震动上的应用。

（3）进行方法验证。使用该方法建立水平层状场地三维有限元模型，求解两层地层结构地震波斜入射下的响应。将结果与理论值进行对比，验证了该方法的准确性和实用性。

图 1 确定有效安定极限的流程图

Figure 1. Flow chart for determining effective shakedown limit

下载: 全尺寸图片幻灯片

图 2 道路结构示意图

Figure 2. Diagram of road structure

下载: 全尺寸图片幻灯片

图 3 层状饱和多孔介质计算模型

Figure 3. Computational model for layered saturated porous media

下载: 全尺寸图片幻灯片

图 4 层状半空间计算模型

Figure 4. Computational model for layered half space

下载: 全尺寸图片幻灯片

图 5 本文计算结果与文献[17]对比

Figure 5. Comparison between calculated results and Reference[17]

下载: 全尺寸图片幻灯片

图 6 荷载移动速度对安定极限的影响

Figure 6. The influence of load moving speed on shakedown limit

下载: 全尺寸图片幻灯片

图 7 荷载移动速度对有效安定极限的影响

Figure 7. Influences of load-moving speed on effective shakedown limit

下载: 全尺寸图片幻灯片

图 8 荷载移动速度对临界深度的影响

Figure 8. The influence of load moving speed on critical depth

下载: 全尺寸图片幻灯片

图 9 荷载移动速度对有效临界深度的影响

Figure 9. The influence of load moving speed on effective critical depth

下载: 全尺寸图片幻灯片

图 10 道路面层刚度对有效安定极限的影响

Figure 10. The influence of pavement stiffness on effective shakedown limits

下载: 全尺寸图片幻灯片

图 11 道路面层刚度对有效临界深度的影响

Figure 11. The influence of pavement stiffness on effective critical depth

下载: 全尺寸图片幻灯片

图 12 荷载频率对有效安定极限的影响

Figure 12. The influence of load frequency on effective shakedown limit

下载: 全尺寸图片幻灯片

图 13 荷载频率对有效临界深度的影响

Figure 13. The influence of load frequency on effective critical depth

下载: 全尺寸图片幻灯片

表 1 层状道路结构的物理力学参数

Table 1 Physical and mechanical parameters of layered road structures

道路结构	μ/(N·m^-2)	λ/(N·m^-2)	M/(N·m^-2)	α	ϕ	ρ_s/(kg·m^-3)	ρ_f/(kg·m^-3)	α_∞	b_p/(N·s²m^-3)	H/m	c/kPa	$c'$ /kPa	$\varphi$ /(°)	$\varphi$ ′/(°)
道路面层	50.0×10⁷	33.30×10⁷	1.0×10^-4	1.0×10^-4	1.0×10^-4	2.5×10³	1.0×10^-4	1.0×10^-4	1.0×10^-4	0.2	1000	—	30	—
道路基层	20.0×10⁷	30.00×10⁷	1.0×10^-4	1.0×10^-4		2.0×10³	1.0×10^-4	1.0×10^-4	1.0×10^-4	0.2	400	—	30	—
路基	1.0×10⁷	2.33×10⁷	2.4×10⁸	9.7×10^-1	4.0×10^-1	2.0×10³	1.0×10³	1.0	2.0×10⁸	5.0	20	16	10~40	12~48
路基	2.0×10⁷	4.67×10⁷	2.4×10⁸	9.7×10^-1	3.0×10^-1	2.0×10³	1.0×10³	1.0	2.0×10⁸	10.0	20	16	10~40	12~48

下载: 导出CSV

参考文献(26)

[1]	ZHUANG Y, WANG K Y, LI H X. Shakedown solutions for ballasted track structure under multiple uniform loads[J]. Transportation Geotechnics, 2020, 22: 100298.
[2]	王娟, 余海岁. 道路安定理论的进展及其应用[J]. 岩土力学, 2014, 35(5): 1255–1262, 1268. doi: 10.16285/j.rsm.2014.05.026 WANG Juan, YU Hai-sui. Development and its application of shakedown theory for road pavements[J]. Rock and Soil Mechanics, 2014, 35(5): 1255–1262, 1268. (in Chinese) doi: 10.16285/j.rsm.2014.05.026
[3]	MELAN E. Der spannungsgudstand eines Henky-Mises schen kontinuums bei verlandicher belastung[J]. Sitzungberichte der Ak Wissenschaften Wie (Ser. 2A), 1938, 147: 73.
[4]	KOITER W T. General theorems for elastic-plastic solids[M]// Progress in Solid Mechanics. Amsterdam: North-Holland Publishing Company, 1960.
[5]	王永刚, 钱建固. 移动荷载下三维半空间动力安定性下限分析[J]. 岩土力学, 2016, 37(增刊1): 570–576. doi: 10.16285/j.rsm.2016.S1.074 WANG Yong-gan, QIAN Jian-gu. Dynamic shakedown lower-bound analysis of three-dimensional half-space under moving load[J]. Rock and Soil Mechanics, 2016, 37(S1): 570–576. (in Chinese) doi: 10.16285/j.rsm.2016.S1.074
[6]	SHARP R W, BOOKER J R. Shakedown of pavements under moving surface loads[J]. Journal of Transportation Engineering, ASCE, 1984, 110(1): 1–14.
[7]	SHIAU S H. Numerical Methods for Shakedown Analysis of Pavements Under Moving Surface Loads[D]. Newcastle: University of Newcastle, 2001.
[8]	YU H S. Three-dimensional analytical solutions for shakedown of cohesive-frictional materials under moving surface loads[J]. Proceedings of the Royal Society A: Mathematical, Physical and Engineering Sciences, 2005, 461(2059): 1951–1964.
[9]	YU H S, WANG J. Three-dimensional shakedown solutions for cohesive-frictional materials under moving surface loads[J]. International Journal of Solids and Structures, 2012, 49(26): 3797–3807.
[10]	WANG J, YU H S. Shakedown analysis for design of flexible pavements under moving loads[J]. Road Materials and Pavement Design, 2013, 14(3): 703–722.
[11]	WANG J, YU H S. Three-dimensional shakedown solutions for anisotropic cohesive-frictional materials under moving surface loads[J]. International Journal for Numerical and Analytical Methods in Geomechanics, 2014, 38(4): 331–348.
[12]	LIU S, WANG J, YU H S, et al. Shakedown solutions for pavements with materials following associated and non-associated plastic flow rules[J]. Computers and Geotechnics, 2016, 78: 218–226.
[13]	QIAN J G, WANG Y G, LIN Z G, et al. Dynamic shakedown analysis of flexible pavement under traffic moving loading[J]. Procedia Engineering, 2016, 143: 1293–1300.
[14]	QIAN J G, WANG Y G, WANG J, et al. The influence of traffic moving speed on shakedown limits of flexible pavements[J]. International Journal of Pavement Engineering, 2019, 20(2): 233–244.
[15]	ZHUANG Y, WANG K Y. Three-dimensional shakedown analysis of ballasted railway structures under moving surface loads with different load distributions[J]. Soil Dynamics and Earthquake Engineering, 2017, 100: 296–300.
[16]	ZHUANG Y, WANG K Y, LI H X, et al. Application of three-dimensional shakedown solutions in railway structure under multiple Hertz loads[J]. Soil Dynamics and Earthquake Engineering, 2019, 117: 328–338.
[17]	LU Z, QIAN J G, ZHOU R Y. Shakedown analysis of flexible pavement on saturated subgrade under moving traffic loading[C]// Advances in Environmental Vibration and Transportation Geodynamics, 2020. Singapore.
[18]	BIOT M A. Mechanics of deformation and acoustic propagation in porous media[J]. Journal of Applied Physics, 1962, 33(4): 1482–1498.
[19]	WANG J, YU H S. Residual stresses and shakedown in cohesive-frictional half-space under moving surface loads[J]. Geomechanics and Geoengineering, 2013, 8(1): 1–14.
[20]	LU J F, JENG D S. A half-space saturated poro-elastic medium subjected to a moving point load[J]. International Journal of Solids and Structures, 2007, 44(2): 573–586.
[21]	HALLONBORG U. Super ellipse as tyre-ground contact area[J]. Journal of Terramechanics, 1996, 33(3): 125–132.
[22]	XU B, LU J F, WANG J H. Dynamic response of a layered water-saturated half space to a moving load[J]. Computers and Geotechnics, 2008, 35(1): 1–10.
[23]	李广信. 高等土力学[M]. 2版. 北京: 清华大学出版社, 2016. LI Guang-xin. Advanced Soil Mechanics[M]. 2nd ed. Beijing: Tsinghua University Press, 2016. (in Chinese)
[24]	ACHENBACH J D, THAU S A. Wave propagation in elastic solids[J]. Journal of Applied Mechanics, 1974, 41(2): 544.
[25]	周仁义, 钱建固, 黄茂松. 不平顺路面的车辆动载诱发饱和地基的动应力响应[J]. 振动与冲击, 2016, 35(11): 93-101, 122. https://www.cnki.com.cn/Article/CJFDTOTAL-ZDCJ201611016.htm ZHOU Ren-yi, QIAN Jian-gu, HUANG Mao-song. Influences of vehicle dynamic load on dynamic stress in saturated poro-elastic ground[J]. Journal of Vibration and Shock, 2016, 35(11): 93–101, 122. (in Chinese) https://www.cnki.com.cn/Article/CJFDTOTAL-ZDCJ201611016.htm
[26]	公路工程技术标准: JTG B01—2014[S]. 北京: 人民交通出版社, 2015. Technical Standard of Highway Engineering: JTG B01—2014[S]. Beijing: China Communications Press, 2015. (in Chinese)

施引文献(2)

期刊类型引用(2)

1.	杜闯，高启辉，宋帅. 地震波斜入射下港珠澳大桥沉管隧道地震响应分析. 河北工业大学学报. 2024(03): 76-86 . 百度学术
2.	邓捷程，吕龙，孙海忠，吕玺琳. 基于等效线性土体模型的地基与结构整体地震响应模拟. 结构工程师. 2023(03): 85-91 . 百度学术