粉土地震液化动力反应的数值分析

罗 强

(攀枝花学院土木与建筑工程学院，四川攀枝花617000)

土体地震动力反应的数值计算分析是求解土体在动荷载作用下，任一时刻土中各点反应值的过程，是定量判别场地液化和估计地震引起变形的有效方法。通过动力分析，有助于对土体振动和波动进行分析，也有助于正确评价土体的抗液化性能。随着计算技术和量测技术的发展，动力分析方法也从最初的线性总应力法，发展为基于非线性有限元单元法基础上的有效应力动力分析方法和采用复杂弹塑性模型并考虑水土耦合作用的动力分析方法［1］。本次使用的有限元程序采用等效线性化模型来考虑土的静力及动力非线性特性。该程序可以进行初始应力分析、动力反应分析，而且也可以估算土体中的孔隙水压力。

1 土的动力本构关系

在等效线性化模型中，土的动力性能用Gmax、G/Gmax～γ曲线和λ～γ曲线表示［2］。Gmax用现场波速试验或室内试验确定，可表示为：

式中：ρ为土的质量密度，Vs为土的剪切波速，或：

式中：K1、n1为试验参数，和土的类型有关。Pa为大气压力，σ0为静平均正应力。

G/Gmax～γ关系较多采用Hardin-Drnevich模型：

式中：γr为参考剪应变。其中：

这里，K2、n2为试验参数。

阻尼比表达式改写成下面形式：

M为试验参数，与土的类型有关。

2 静力有限元分析

进行静力有限元分析，目的是了解土体内的静应力分布情况，确定土体内各单元的初始静围压，为动力分析提供初始条件。其有限元求解方程为：

式中：［K］为求解体系的总刚度矩阵，由单元刚度叠加而成，单元刚度矩阵取决于单元的尺寸和材料的模量；{U}为在静荷载向量作用下结点发生的静位移向量，是待求量，求得位移向量后，即可求出土体内各计算单元的应变和应力分量；{R}为结点的静荷载向量，在形成静荷载向量时考虑材料的自重、渗透力和其它外荷载的作用。

考虑到土体的本构关系是非线性的，因此在程序计算求解时采用增量迭代法，一般3～5次迭代即可。

3 粉土的动力反应分析

3.1 动力有限元分析

选定所要分析的二维场地区域后，进行有限元离散，可得计算体系的动力有限元求解方程：

式中：［M］为体系的质量矩阵，采用集中质量法合成；［K］为体系的总刚度矩阵，由单元刚度矩阵合成，取决于土体单元剖分情况和土壤模量；{E}x、{E}y分别为水平向地震荷载和竖向地震荷载的质量列阵，{E}x的奇数元素为相应结点的质量，偶数元素为0，{E}y的奇数元素为0，偶数元素为相应的结点质量；{U}、{}和{} 分别为结点相对基岩运动的位移、速度和加速度向量，为待求量；［C］为体系的阻尼矩阵，由单元阻尼矩阵合成，单元阻尼矩阵采用瑞利阻尼形式：

其中：

式中：ω为体系的自振圆频率，通常采用主阵圆频率ω1，ξ为材料的阻尼比。

动力有限元方程的求解采用wilson-θ法，θ取值为1.4。由式(7)求出每一时刻位移后，就可求出相应时刻土单元的应力分量，进而可以计算粉土层的孔隙水压力。

3.2 边界的处理

对于底边界采用静止底边界，认为基岩输入地震动为边界不动而在计算域各点上作用着惯性力的问题。对于侧边界认为进行静力计算时横向固定、竖向自由；而在动力计算时为竖向固定、横向自由。

3.3 孔压和应变计算模型

粉土孔隙水压力的增长模型可采用如下函数表示：

式中：αu为孔压比；λ为时间参数。αu0、λ0分别为第一点的孔压比和时间参数，a、b为拟合所需的模型参数，这4个参数是最大剪切作用面上往返剪应力比的函数［3］。

3.4 液化判别方法

在此计算程序中，地震液化采用孔隙水压力比来判别。如果孔隙水压力比达到1，则认为发生液化。

4 数值计算

4.1 有限元离散

本次计算区域选取的二维计算区域的深度和宽度分别为1.5 m和3.8 m。假设土层分布均匀，用矩形单元对计算区域进行有限元离散，整个计算区域共包含570个单元，624个节点。

4.2 地震动输入

本次计算工作输入地震动是通过将取土场地50年超越概率10%的地震加速度反应时程按一定相似律压缩放大后得到的。HZ1、HZ2的地震动时程如图1。

图1 计算输入地震加速度时程曲线

4.3 数值计算结果与分析

分别输入地震动HZ1、HZ2进行计算，不同深度土层的计算结果见图2～图5。图2为输入地震波计算得到的孔压比，图3为输入HZ1波计算得到的加速度和水平剪应力时程，图4为输入HZ2波计算得到的加速度和水平剪应力时程，图5为输入HZ1波作用时产生的竖向加速度时程。

图2 输入地震波计算得到的孔压比

图3 输入HZ1波计算得到的水平加速度和水平剪应力时程

5 结论

(1)通过以上计算工作可以发现，所输入的地震动都将造成粉土不同程度液化的发生。

(2)从孔隙水压力增长时程曲线可见，三个地震动输入情况下上部土层都发生了液化，下部的粉土只是产生接近液化的孔压或者不发生液化。这表明上部粉土层更容易液化，抵抗液化的能力比埋深大的土弱。

图4 输入HZ2波计算得到的水平加速度和水平剪应力时程

图5 输入HZ1波作用时产生的竖向加速度

(3)地震动作用时，液化后的土层将对地震动产生一定的滤波和放大作用。从输出点的加速度时程可见，深度不同土层的加速度时程曲线幅值由下向上逐渐增大，并存在一定的相位差，但其形状几乎完全相同；剪应力时程曲线具有相反的规律，由下向上逐渐减小。

(4)尽管只是输入水平向地震加速度，但在计算过程中发现输出了一定的垂向加速度，且垂向加速度在土层中也是由下向上逐渐增大。