基于ABAQUS的土层地震反应分析

李煜东，李 平，孙强强，乔 峰，王 亮

（1.防灾科技学院，河北三河 065201；2.中国地震局工程力学研究所，黑龙江哈尔滨 150080）



基于ABAQUS的土层地震反应分析

李煜东1，李 平1，孙强强1，乔 峰2，王 亮2

（1.防灾科技学院，河北三河 065201；2.中国地震局工程力学研究所，黑龙江哈尔滨 150080）

摘 要：基于一维波动理论和Spring/Dashpots单元实现了ABAQUS中粘性边界条件的施加，进而对弹性半空间的波源问题进行波动数值模拟，验证了ABAQUS软件中粘性边界条件的有效性。分别采用ABAQUS软件和RSLEIBM程序计算响嘡台阵3号测井的土层地震反应，并将计算结果与实际记录进行对比分析。结果表明，在地震动不大的情况下，ABAQUS计算结果能够较好地模拟地表地震反应，较RSLEIBM程序计算结果具有更高的准确性。

关键词：ABAQUS软件；粘性边界；RSLEIBM程序；响嘡台阵；土层地震反应分析

0 引言

震害表明，场地条件对地震动和结构震害的影响十分显著。开展土层地震反应分析工作，研究场地条件对地震动和结构震害的影响，为工程抗震提供合理的设计地震动参数具有重要意义［1］。土层地震反应分析是工程抗震的重要部分，目前我国常用的土层地震反应分析方法为等效线性化方法，其概念清晰，易于被广大工程师接受，在我国的工程建设中得到广泛应用。

然而等效线性化方法及其相应程序在大量的工程应用中也暴露出一些不足之处，如输入地震动较大、场地较软时计算结果误差较大；含软弱夹层的剖面计算结果不合理；存在死循环现象；高频部分计算误差大等［2-3］。随着计算机技术的发展，基于有限元的数值模拟逐渐成为研究场地效应的重要手段。ABAQUS作为一种功能强大的通用有限元数值模拟软件，拥有多种岩土本构模型，不仅能够处理复杂的边界条件问题，还能较好地处理非线性问题，是进行土层地震反应分析的有力工具。本文基于一维波动理论，实现了ABAQUS软件中粘性边界条件的施加，进而利用ABAQUS软件及RSLEIBM等效线性化程序分析响嘡场地的土层地震反应，将计算结果与地表实际记录进行对比分析，验证ABAQUS软件在土层地震反应分析中的可靠性。

1 边界处理

进行土层地震反应分析时需要从半无限介质中切取有限的计算域进行计算，若基于刚性基底假设在基底直接输入地震动加速度时程，会使得波动能量在计算域内往复反射，这会影响数值分析结果的可靠性。为保证散射波向无限地基传播，笔者基于一维波动理论，在ABAQUS软件中实现粘性边界条件的施加，以模拟无限地基的辐射阻尼效应。

1.1 粘性边界

粘性边界是一种局部人工边界，最早由lysmer和Kuhlemeyer基于一维波动理论推导提出［4］，其概念清晰，便于实现，是目前土层地震反应分析中应用最广泛的人工边界条件之一［5］。

粘性边界条件的原理可以由半无限弹性杆内的波动方程推导得出。如图1所示，对于在半无限弹性杆内传播的平面波，其一维标准波动方

根据达朗贝尔原理，上述一维标准波动方程（1）的通解为

式中，u为总波场，ui为内行波，uo为外行波。

将半无限弹性杆截断，获得有限的计算域。以S波为例，由边界上的应力条件可得：

式中，τ为剪应力，γ为剪应变，G为剪切模量，ρ为密度，cs为剪切波速

对于P波，同样由边界上的应力条件可得：

式中，cp为压缩波速。

由式（5）和式（6）可知，只要在人为截断的边界处加上合适的阻尼（如图2），便能使边界处的应力状态和截断前一样，从而模拟实际散射波透过边界时的应力状态。

对于波源问题，没有入射的内行波，总波场只有外行波，因此只要在边界处加上阻尼系数为ρc的粘滞阻尼器，便可以使边界处的应力条件与截断前一样，从而满足外行波向无限远处传播的应力条件；对于散射问题，除了外行波还有入射的内行波，这时除了在边界处加上阻尼系数为ρc的粘滞阻尼器外，还需要将内行波转化为的等效荷载进行输入，以便模拟截断前边界上的应力条件（其中对于P波，c为压缩波速；对于S波，c为剪切波速；可由输入地震动加速度时程积分得到）。

由于粘性边界是基于一维的基础推导得出，如果简单地将其应用到多维情况会导致较大误差［6］。而实际上多数情况下土层为水平成层，根据斯奈尔定理，地震波到达近地表时近似于垂直入射，因此大多数情况下可以将土层地震反应分析作为一维问题来处理。在采用有限元软件ABAQUS进行水平成层土层场地地震反应分析时，通常将地震动加速度时程转化为力的时程施加到底边界，对应于在各节点施加粘滞阻尼，可实现粘性边界的施加。其实现方法主要有两种：一种方法是在ABAQUS的 CAE模块中，通过Spring/Dashpots单元在人工截断边界的各节点处施加一系列粘滞阻尼器来实现；另一种方法通过修改inp文件，在inp文件中添加*Element，type =Spring/Dashpots的语句，来定义边界处的弹簧-阻尼器单元，实现粘性边界条件的施加［7-8］。

1.2 数值算例

为检验ABAQUS中施加粘性边界的正确性，分别采用粘性边界和远置边界模型对弹性半空间的波源问题进行数值模拟。以远置边界模型的模拟结果作为验证条件［9］，计算模型如图3所示，对比分析模型底面水平向加速度时程。由于本文进行一维土层地震反应分析，因此边界验证采用一维剪切梁模型，取单位面积的竖向土柱进行研究，土柱高20m，宽1m，密度1700kg/m3，弹性模量1.7×106Pa，泊松比0.25，上部垂直输入剪切脉冲波，周期为0.5s，幅值为1g，采用修改inp文件的方法，在土柱底部定义弹簧-阻尼器单元，实现粘性边界条件的施加。

图4为粘性边界与远置边界计算结果对比图。可以看出，粘性边界模型和远置边界模型底部的计算加速度时程曲线吻合程度很好，说明粘性边界很好地吸收了散射源产生的散射波，未造成波动能量在计算域的往复反射，较好地模拟了实际情况下散射波向半无限地基传播。

2 土层地震反应分析

从上述波源问题的数值模拟结果可以看出，ABAQUS中施加的粘性边界条件能够很好地模拟散射波向弹性半无限空间传播的情况，与实际波动传播规律基本一致，证明了ABAQUS中基于一维波动理论施加的粘性边界条件的可靠性。在此基础上，收集唐山响嘡三维台阵土层剖面资料，建立响嘡台阵3号测井的土层计算模型，以实测基岩地震动作为输入，分别采用ABAQUS软件及RSLEIBM程序计算响嘡台阵土层计算模型的地震反应，计算得到地表加速度时程、加速度傅里叶谱及加速度反应谱，并将计算结果与地表实测记录进行对比分析。

2.1 模型及地震动选取

唐山响嘡台阵始建于1992年，目前共设置有3个测井，是我国首个研究场地效应的三维观测台阵，目的在于利用其详细的土层剖面资料以及相应的实测地震动加速度记录来研究场地条件对地震动的影响，进而检验现行土层地震反应分析的理论、计算模型和方法的正确性［10］。由于响嘡台阵建台初期未对1、2号测井内土样进行土动力学试验，因此本文选取拥有详细土层剖面资料的3号测井建立土层计算模型。3号测井土层剖面由地表往下依次为杂填土、砂土、细砂、粉砂、粘性土、粉土、细砂及粉粘土等，下卧岩层为花岗岩混合岩［11］，其详细土层剖面资料与土的非线性参数列于表1和表2［12］。

表1 响嘡台阵3号测井场地土层参数表Tab.1 Parameters of soil layers at Xiangtang array No.3 logging

响嘡台阵3号测井在地表及井下47m深处各布置了一个三分量测点，能够记录地表及井下47m处的东西向、南北向以及竖向地震动加速度时程，由于计算使用基岩地震波取自井下47m处，因此本文选择井下47m处实测记录点作为输入基岩面，分别选取2001年及2002年滦县地震所记录的GZ005.EW和GZ019.EW基岩加速度时程作为输入地震动，选用的地震动记录具体资料如表3所示，其加速度时程及傅里叶幅值谱如图5及图6所示。其中，GZ005.EW基岩地震动峰值为31.73gal，持时为2.17s；GZ019.EW基岩地震动峰值为22.99gal，持时为2.36s。

表2 响嘡台阵3号测井土样动剪切模量比和阻尼比Tab.2 Dynamic shear modulus ratio and damping ratio of soil layers at Xiangtang array No.3 logging

表3 地震动资料表Tab.3 Data of seismic ground motion

2.2 计算结果对比分析

分别采用ABAQUS软件及RSLEIBM程序计算得到地表加速度时程，对比分析其峰值、频谱、持时及反应谱与地表实际记录的差异，分析如下。

图7、图8分别是GZ005.EW及GZ019.EW基岩地震动输入下实测和计算绘出的地表加速度时程曲线，由地表地震动加速度时程曲线形态可以看出，两种地震动输入下，ABAQUS计算结果和RSLEIBM计算结果与地表实际记录十分接近，均在3s左右达到峰值点，之后迅速衰减，10s之后地表地震反应基本结束。

使用意大利seismosoft公司的地震波专用处理软件SeismoSignal得到两种不同地震动输入下的地表加速度峰值，列于表4。

可以看出，ABAQUS软件及RSLEIBM程序均高估了地表地震反应，但RSLEIBM程序对地表地震动的放大效应较为严重，ABAQUS软件计算峰值则与实际记录更为接近。

其中，当输入地震动为GZ005.EW时，基岩记录峰值为31.73gal，对应地表实测加速度峰值为96.03gal，放大倍数为3.03倍。ABAQUS软件计算地表加速度峰值为107.27gal，放大3.38倍，峰值和放大倍数都与地表实际记录接近。RSLEIBM程序计算地表加速度峰值为141.39gal，放大4.46倍，其加速度峰值比地表实测记录增大47%，计算结果与实际记录相差较大。

表4 计算地表加速度峰值对比表Tab.4 Seismic data

GZ019.EW输入下地表实测加速度峰值为48.8gal，放大倍数为2.12倍。ABAQUS软件计算地表加速度峰值为66.02gal，放大2.87倍，峰值和放大倍数都与地表实际记录较为接近。RSLEIBM程序计算地表加速度峰值为77.8gal，放大3.38倍，结果与实际记录相差较大。

将两种基岩地震动输入下计算得到的加速度时程进行傅里叶变换，得到加速度傅里叶幅值谱，分别绘于图9、图10中。

从图9和图10中可以看出，两个地表记录的主要频率成分均分布于2～20Hz之间，高于20Hz的高频成分较少。在两个基岩地震动输入下，ABAQUS软件及 RSLEIBM程序计算结果均对2～6Hz频率范围内的地震波有所抑制，且未能反映出地表实际记录在8～10Hz内出现的傅里叶谱峰值。而RSLEIBM程序计算得到的地震动均出现了严重的高频放大作用，由图中可以明显看出，两个基岩地震动输入下RSLEIBM程序计算结果均在15Hz附近出现高峰值点，与实际记录差别较大。综上所述，在GZ005.EW和GZ019.EW输入下，ABAQUS软件及RSLEIBM程序计算结果在2～12Hz范围内基本相同，而RSLEIBM程序计算结果对高频成分存在明显放大作用，ABAQUS软件计算结果则在15～25Hz的频范围内比RSLEIBM程序计算结果更接近实际记录。

图11、图12分别为GZ005.EW及GZ019. EW基岩地震动输入下的计算地表加速度反应谱对比图。可以看出，两个地表实测加速度反应谱峰值主要分布在周期为0.05～0.2s的范围内。当输入地震动为GZ005.EW时，实际记录的反应谱在0.12s处出现峰值，而 ABAQUS软件及RSLEIBM程序计算的反应谱都呈双峰特征，且反应谱峰值都较实际记录偏小，但在小于0.09s的短周期范围内，ABAQUS软件计算结果更接近实际记录。当输入地震动为 GZ019.EW 时，ABAQUS软件及RSLEIBM程序计算的反应谱曲线形式及走势基本一致，都呈现单峰特征，且都整体向短周期段偏移。但在小于0.09s的短周期部分，RSLEIBM程序计算结果存在严重的放大现象，计算结果较 ABAQUS软件计算结果误差更大。

为对比分析地表地震动持时，选择加速度峰值的十分之一为阀值［13］，计算得到地表地震动括弧持时。由表5可以看出，两个地表实测记录的持时分别为4.42s、5.67s，对应的ABAQUS软件计算的地表地震动持时分别为3.33s、4.69s，比实测记录稍少，但与实际记录持时相近；而RSLEIBM程序计算持时均小于2s，与地表实测地震动持时相差较大。

表5 计算地表加速度持时对比表Tab.5 Seismic data

综上，ABAQUS软件计算结果能够较好地反映地表的地震反应，较RSLEIBM程序计算结果而言，与实际地震记录更接近。

3 结论

本文在推导粘性边界原理的基础上，实现了ABAQUS中粘性边界条件的施加。分别采用ABAQUS软件和RSLEIBM程序计算响嘡台阵3号测井土层地震反应，对比分析地表加速度峰值、傅里叶幅值谱、持时及反应谱，最后得出以下结论：

（1）粘性边界能吸收散射的波动能量，避免散射波的反射，模拟散射波向地基无穷远处传播的情况。采用ABAQUS进行波动的数值模拟时，利用Spring/Dashpots单元设置接地阻尼，同时将输入地震波转化为等效节点力时程进行输入，能够方便地实现粘性边界条件的施加，提高了ABAQUS进行土层地震反应分析的适用性。

（2）通过ABAQUS软件及RSLEIBM程序计算响嘡台阵3号测井的地震反应，对比分析计算结果可见，小幅值输入情况下，ABAQUS软件计算结果能够较好地模拟地表地震反应，较RSLEIBM程序计算结果具有更高的准确性。

（3）RSLEIBM等效线性化程序已在工程中得到广泛应用，其在大多数情况下计算结果较为可靠，然而等效线性化方法毕竟只是考虑土体非线性的一种简化的线性方法，不能真实反映土体在地震动下的非线性特性。ABAQUS则是一种功能强大的有限元软件，不仅能较好地处理非线性问题，还能够方便地应用于工程之中，在现行等效线性化方法及相应程序存在诸多不足的情况下，基于有限元的ABAQUS软件不失为一种更加有效地估计土层地震反应的途径。

参考文献

［1］ 薄景山，李秀领，李山有.场地条件对地震动影响研究的若干进展［J］.世界地震工程，2003，19 （2）：11-15.

［2］ 齐文浩，薄景山.土层地震反应等效线性化方法综述［J］.世界地震工程，2007，23（4）：221-226.

［3］ 刘德东，齐文浩，张宇东，等.现行土层地震反应分析存在的问题［J］.防灾科技学院学报，2009，11（3）：34-37.

［4］ Lysmer J，Kuhlemeyer R L.Finite Dynamic Model for Infinite Media［J］.Journal of the Engineering Mechanics Division，1969，95：859-877.

［5］ 刘晶波，谷音，杜义欣.一致粘弹性人工边界及粘弹性边界单元［J］.岩土工程学报，2006，28 （9）：1070-1075.

［6］ 白建方.场地地震反应分析导论［M］.北京：中国铁道出版社，2010.

［7］ 何建涛，马怀发，张伯艳，等.黏弹性人工边界地震动输入方法及实现［J］.水利学报，2010，41 （8）：960-969.

［8］ 蒋云峰.速度大脉冲作用下土—结构相互作用问题的初步研究［D］.哈尔滨：中国地震局工程地震研究所，2009.

［9］ 廖振鹏.工程波动理论导论［M］.北京：科学出版社，2002.

［10］ 谢礼立，章文波，李沙白.唐山响堂三维场地影响观测台阵［J］.地震工程与工程振动，1999，19 （2）：1-8.

［11］ 卢滔，周正华，霍敬妍.土层非线性地震反应一维时域分析［J］.岩土力学，2008，29（8）：2170-2176.

［12］ 齐文浩.土层地震反应分析方法的比较研究［D］.中国地震局工程力学研究所，2004.

［13］ 袁一凡，田启文.工程地震学［M］.北京：地震出版社，2012.

中图分类号：P315.9

文献标识码：A

文章编号：1673-8047（2016）02-0064-08

收稿日期：2016-02-29

基金项目：国家自然科学基金项目（51408559，51508096）；中央高校基本科研业务费研究生科技创新基金项目（ZY20150302）；中国地震局地震科技星火计划项目（XH15067Y）

作者简介：李煜东（1991—），男，硕士研究生，主要从事工程地震方向的研究。程为

Seismic response analysis of soil layer based on ABAQUS

LI Yudong1，LI Ping1，SUN QiangQiang1，QIAO Feng2，WANG Liang2
（1.Institute of Disaster Prevention，Sanhe 065201，China；2.Institute of Engineering Mechanics，China Earthquake Administration，Harbin 150080，China）

Abstract：Based on wave motion theories and the Spring/Dashpots elements of ABAQUS software，the development of viscous boundary condition is realized.Then the viscous boundary condition is used to analyze the wave response of source problem in elastic half space.The records of Xiangtang borehole array（3#）site are also used to validate the seismic response with the ABAQUS software and the RSLEIBM procedures，respectively.The results of this study suggest that，when the earthquake motion is not strong，the results of ABAQUS software，which can simulate seismic response of soil layer appropriately，are more accurate than the results of RSLEIBM procedures.

Keywords：ABAQUS software；viscous boundary condition；RSLEIBM procedures；Xiangtang array；soil layer seismic response analysis