分析三维土_结构动力相互作用体系等效输入的时域显式有限元法

1、分析三维土-结构动力相互作用体系等效输入的时域显式有限元法杨笑梅1关慧敏2张学胜2杨柏坡2（1.广东工业大学，广东广州 510006；2. 中国地震局工程力学研究所，黑龙江哈尔滨 150080）摘要：基于直接法中分析三维土-结构相互作用问题的两步法思想，提出关键性第一步等效输入计算的时域简化算法。本方法在采用集中质量显式有限单元法结合显式人工边界的同时，借鉴子结构法中刚性基础假设及动力 相互作用的概念，简化三维土-结构相互作用体系动力计算过程。和其他同类方法相比，本方法完全采用时域显 式递推的过程，计算效率较高，计算量少；等效输入的获得，可使土-结构复杂体系的动力分析等效为常规的上部 结构动力

2、反应计算，时域算法较常用的频域法在土体非线性方面存在更大潜力。本文还通过典型的三维土-结构 模型的建立，与相应的参考解相比，验证本文提出的计算方法的有效性和计算精度，并初步探讨土-结构动力相 互作用体系的一些基本规律。关键词：土-结构动力相互作用；等效输入；时域显式有限元；子结构法中图分类号：TU435文献标识码：A文章编号：1000-131X（2010）增-0535-06The explicit finite element method in time domain for analysis of equivalentinput of three-dimensional dynamic S

3、SI systemYang Xiaomei1Guan Huimin2Zhang Xuesheng2Yang Baipo2（1. Guangdong University of Technology, Guangzhou 510643, China；2. Institute of Engineering Mechanics, China Seismological Bureau, Harbin 150080, China）Abstract: For decreasing the computation cost, a new explicit numerical simplified proce

4、dure in time domainisproposed for the rigid boundary equivalent input computation for three-dimension dynamic SSI system in this paper.Based on the lumped-mass explicit finite element procedure with a local transmitting artificial boundary in time domain, the rigid foundation assumption and the dyna

5、mic interaction relationship usually used in substructure method are adopted to simplify the computation procedure. This method has advantages in computational feasibility and the applicability in engineering when the explicit finite element and step-by-step integration scheme is adopted. The lump-

6、beam model for upper structure and the rigid foundation assumption can make the cost of computation equal to that of free field computation. Compared to the computation procedure in frequency domain, the computation in time can be advanced in deal with the nonlinear soil body. In addition, the equiv

7、alent input of upper structure obtained by the simplified procedure make the complex SSI system analysis equal to the common upper structure dynamic analysis. To illustrate the efficacy and reliability of the procedure, several three-dimension SSI numerical examples are presented. The computational

8、predictions based on the assumption of presented procedures are in good agreement to the reference solutions. At the same time, numerical examples are also used to study for the influence on structural seismic response by the different seismic input cases usually used in design.Keywords: soil-struct

9、ure interaction; equivalent input; explicit FEM; substructure methodE-mail: xmyang_86yahoo 引言基金项目：广东省基金资助（06300263）和地震工程与工程振动开放实验室基金（2009A09） 作者简介：杨笑梅，博士，副教授 收稿日期：2010-08-24为正确估算结构在地震作用下的反应时，应建立可考虑土-结构动力相互作用的合理计算模型。最 直接的考虑是将土体与结构看成一个整体，建立统536土 木 工程 学 报2010 年一的动力学方程并求解，这就是通常所说的一步直接法或直接法1-2。为简化计算，直接法可

10、分两步完 成，即第1步将上部结构简化为集中质量-杆系模型， 模型应能反应总体结构的动力特性，计算后可以得 到结构与基础相连处的反应，称为复杂结构的“ 等 效输入”；第2 步将等效输入作为复杂模型（此模型 仅为一步法的上部结构模型）刚性基底处的地震输 入，计算结构的地震反应。这一过程简称两步法3。 两步法的优点在于计算量较之常规子结构法增加并 不多，且具有直接法计算过程简单的特点，因此成 为现时分析三维土-结构性互作用体系最有效的方 法之一。但目前两步法中等效输入的计算主要依赖 专用程序FLUSH4及SASSI5等完成，计算采用频域 内隐式的计算过程，因此不能直接考虑土的非线 性，只能用等效线性

11、化方法间接考虑土的非线性， 另外隐式的方法需要求解线性方程组，且计算时间 相对较长，很大程度上影响了它们的计算精度和解 题规模，因此模拟复杂的实际工程问题时受到很大 局限。为克服两步法中关键性第1 步频域隐式算法的 局限，本文提出了改进等效输入计算的时域简化计 算方法，借鉴常规子结构的刚性基础假设6及动力 的相互作用的概念3，6，在刚性基础处建立土体、结 构及基础的时域动力平衡关系，将土-结构体系的分 析统一起来，整个分析完全采用高效的显式有限元 数值过程。本文方法的优势体现在几个方面： 不 需要求取阻抗函数及自由场计算；计算直接在时 域完成，时域方法将在非线性的处理上存在巨大潜 力，克服采用

12、频域法4-5在非线性处理上的麻烦； 利用简化的上部结构及子结构法中刚性基础假设及 力的相互作用关系，使计算可以在土体、上部结构 分别进行，再通过刚性基础的平衡条件，统一实现 土-结构体系的动力反应，这在很大程度上克服了三 维土-结构相互作用分析中计算量过大的局限；显 式有限元计算过程不需要求解方程组，提高了三维 土-结构相互作用分析中计算规模及计算效率，局部 透射人工边界在保证无限域问题简化为有限域问 题过程的计算精度。本文还以典型土-结构模型为例 初步验证了提出方法基所编程序的可靠性，并分析 典型体系的地震输入差异。计算，将分析区域分为土体加人工边界（图1（b）、刚性基础（图1（d）及上部结

13、构（图1（c）三部分，分 别用上标g，b，s表示。具体推导过程介绍如下。1.1土体运动如图1（b）所示，考虑阻尼时，采用时域集中质 量显式有限单元法，土体内点i解耦的运动方程可以 写成 1+ Ksp s p- 1 t2ggil ul- t Kil ulup+1pp- 1= 2ui - ui -iMgtill（1）其中：Mi 为集中于内点i的质量；Kil 为反应内点i 与相 邻节点l 之间弹性恢复力大小的刚度系数；u 为节点时域显式有限元计算结构基础处等效输入1的位移向量；g 为与u 成正比的黏性阻尼系数；t 为l为确定结构等效输入，计算如图1模型基础处的动力反应即为结构刚性边界处的等效输入。为

14、简化时间步长；p为非负整数，满足时间t=pt。人工边界上的节点无法直接计算，采用文献3的方法计算。第 43 卷 增刊杨笑梅等分析三维土-结构动力相互作用体系等效输入的时域显式有限元法537土体与基础相连的结点，由基础运动确定。1.2结构运动如图1（c）所示，考虑阻尼时，将上部结构单元结 点划分为两类、即与基础相连的节点及内点。结构内点仍采用集中质量显式有限单元法，考虑瑞利阻 尼系数分别为s及s时，其解耦的运动方程可以写成：部结构各点的惯性力向基础形心简化，即：MsTz =AiFsi （4）TTTF s = FsTxFsTy FsTz MsTx MsTyn其中FTs 为上部结构作用力向刚性基础形

15、心简化后的合力，这样的处理是因为Fs i可由式（5）直接确定Fsi = - 1+ Ksu p + Msu p-ssil b li s ittib 1+t Ksil ubp- 1 -M i usi tt2sKil ul +up+1pp- 1s s= 2ui - ui -s p（5）isp- 1MslitisibMi ui - M u ss Ks up- 1+ s p- 1s p（2）及tili ittl 10100000010000- ZiYi100Zi0- Xi010- YXi0001i式中各量如前述。1.3基础运动土-结构相互作用的连接部分是基础，基础传递 土-结构相互作用力，因此，基础在土

16、-结构相互作用分析中占着非常重要的位置。由于核电厂或高层建筑基础底板较厚，在计算时通常假设基础是刚性 的，基础的运动只有6 个分量，即3 个水平分量，3 个 转动分量。土和结构作用在基础上的合力使基础产 生刚体运动，继而使土和结构与基础相连点的运 动，从而影响着土和结构的其他部分的运动。因此 研究基础的运动，必须先分析土和结构对基础作用 的合力，求出基础的动力反应。文献3及文献7 采用这一思想，分别用对上部结构用模态分析的方 法计及通过商业软件中的时域计算过程计算结构 作用力。由于本模型上部结构简单，且自由度少，本 文采用与土体统一的时域显式有限元过程计算上部 结构的作用力，可直接实现，且计算

17、量较小，不必像 文献7一样通过商业软件计算复杂结构体系的反 应，再求合力及返回计算土体，这样的计算过程和 一步法相比，自由度数目并不减少且更为复杂。因 此本文提出如下的过程简化计算。根据刚性基础假设，土-结构相互作用的受力分 析力学模型见图1（d）。为求基础及上部结构对基础 的作用力，将基础部分的动平衡方程表示为：00Ai =（6）000 其中：Xi、Yi 及 Zi 分别为节点i 对基础形心的相对坐标。考虑土体对刚性基础产生的作用时，若土体采 用黏性阻尼Cbg= g Kbg，直接可得土体在其与基 础相连的结点k处作用的弹性恢复力：Fgk=FgxFgy FgzkT=- km m + 1+ ggK

18、g up u （7）Kgp- 1km mttmnmn式中：mn 为与点k（土体与基础相连点）相连的土体内点总数，这里Fgk为3阶向量。若将Fgk 向基础形 心简化，得：MgTz =AkFgk （8）F Tg = FgTxFgTy FgTz MgTx MgTyTTm其中：F Tg 为土体作用力向刚性基础形心简化后的合力；m表示土体与基础相连的节点总数；Ak为式（6）的前3行子矩阵。考虑刚体性质，将与基础相连的每个结构或土体节点的惯性力（Msibbus ，Mg ug ）向刚性基础ii ibbbb形心简化，则刚性基础动平衡方程（3）简化表达为：Mbu +Cbu +Kbu =F +F （3）bbbgs

19、Mub=F T+F T（9）其中：质量阵Mb由基础质量、与基础相连的结构bgs及土体点的质量组成；加速度向量u 由刚性基础其中：ub为基础质心加速度向量；F T及F T分别由bbsg质心加速度、以及与基础相连的结构及土体各点的加速度组成；Fs及Fg分别表示为上部结构及土体 对基础的作用的合力，但通常并不作用于基础形 心。此时刚性基础上受力分析如图1（d），若忽略刚 性基础的阻尼，对于结构部分产生的弹性恢复力Fs，可根据力学定律向基础形心简化，其等效于上式（4）及式（8）计算；M为Mbb的简化矩阵，包含刚性基础质量、转动惯量及与基础相连的土体及结 构点相对刚性基础形心的简化后的集中质量及转动 惯

20、量。需要指出，文献7-8中在此步处理上没有考 虑M中与基础相连的结构点及土体点的对刚性体系质量矩阵的贡献。由于土体的密度与刚性基础的 质量密度差异不会超过一个量级，这样的简化可造538土 木 工程学报2010 年成不小的误差。需要说明的是本文在此假设与基础相接触处的土体为均匀介质，并忽略了厚度对基础 刚体平面运动时加速度的影响。仍采用中心差分格式得基础反应的递推公式：ubp+1 = 2ubp - ubp- 1 - t2 M-（1 F T+F T） （10）bbbsg由刚体运动决定了土体和结构与基础相连点t 的运动，其算法用下式表示：=Atub（11）ubp+1b p+1t式中At取式（6）。数

21、值计算2鉴于本文并不着重于上部结构的动力特性分析，因此算例的选取主要针对于本文两个重要简化 即基于刚性基础假设及等效输入的概念，用以说明 本文方法的合理性、有效性并验证相关精度。另外， 还对现行工程中常用的结构基底地震动输入法作了 初步的比较，指出各类方法存在差异。所有算例在 人工边界处输入垂直入射水平单位脉冲加速度时 程，宽度取0.1。算例1 均匀半空间有质量的刚性表面矩形基础 的时程反应本例选取弹性半空间上作用一有质量刚性基础 模型，用以验证本文刚性基础假设在时域显式计算 时的合理性及精度。具体模型为: 基础尺寸L W H=78m77m6.5m，容重2.4t/m3，半空间土体参数 为剪切波

22、速为Vs=200m/s，容重=1.7t/m3。根据显式有限元的网格精度计算要求及稳定性条件选取网格尺寸及时间步距，计算所得形心处的水平u 及基础边 缘处的竖直w 加速度时程反应曲线，如图2所示。图 中同时给出了采用一步法模型进行验证的计算结 果，该模型应用已被广泛验证3，9的三维显式有限元 结合透射边界方法，对于本算例，基础的剪切波速Vs 分别取为：3000m/s，5000m/s，8000m/s。从两种方法结 果的对比，可以看出两者水平分量结果吻合较好， 竖直分量在刚性基础Vs达8000m/s时，时程反应曲线 果也吻合较好，说明本文方法中的刚性基础假设是 合理及可靠的。算例2性验证土-结构相互

23、作用体系等效输入的可靠通过上节算例，可认为土与刚性基础的反应是合理的，在此基础上，本例选取图1（a）所示完整模 型，以验证用本文方法计算上部结构反应的可靠 性，一并验证等效输入概念的合理性。如图1（a）上 部结构简化为一集中质量杆系模型，由两根梁单元 和两个质点组成，梁单元截面性质及质点的质量和 转动惯量见表1，结构的杨氏模量E=3.0 1010MPa，剪切模量G=1.2 10 MPa，基础为无质量及厚度的10刚性表面基础，下面为一覆盖土层，厚50m，覆盖土层下面为均质基岩，土和基岩的物理参数见表2。表 1 梁单元截面性及质点信息Table 1 Properties of the beams

24、for the simplified structure转动惯量（kgm2）截面惯性矩（m4）梁号（质点号）梁高度（m）质点质量（kg）横截面积（m2）X轴Y轴Z轴X轴Y轴Z轴126.06.21704019930274.7274.712.1114.1712.1114.1710101679001679001679001679001000010000第 43 卷 增刊杨笑梅等分析三维土-结构动力相互作用体系等效输入的时域显式有限元法539表 2 地基物理参数 Table 2 Physical parameters of soil medium and bedrock 算例3 软弱地基，本文的等效输

25、入和现行常规方法计算土-结构相互作用体系地震输入结果的比 较计算模型和算例2相同，取截止频率f=25Hz，为比较土-结构模型与目前设计常用计算模型的差异。 受到垂直入射的剪切脉冲型波作用，考虑以下3种情 况，分别再在ABAQUS刚性基底模型底边界处输入， 计算上部结构最高点的时程反应。（1）均匀半空间基岩地表反应（基岩参数取自表2），此输入可视为直接在结构输入基岩地震记录工 况。（2）不考虑基础及上部结构的自由场反应加速 度时程作为输入，此情况可视为目前采用地震安全性评价规范情况获得的地震输入工况。（3）考虑土-结构体系相互作用的本文方法得到 的等效输入（算例2结果）。上述3 种情况的输入时程

26、和计算结果分别绘于图4和图5。可以看到软土介质时，考虑相互作用后， 结构基底输入出现摆动分量，上部结构的加速度时程反应也明显不同于直接输入地震波（工况1）或仅考虑土层反应情况的地震输入模型（工况2），说明 建立合理的土-结构计算模型的必要性。土类剪切波速V（s m/s）密度（t/m3）泊松比阻尼比土基岩5007001.71.80.420.410.070.05本文方法可直接得到刚性基础处和结构各质量点处的动力反应时程，并将最高点计算所得时 程绘于图3 。 由 于本算 例 很 难 直 接 利 用 通 用 软 件 ABAQUS10进行计算，且本例仅需检验本方法上部 结构反应的合理性，因此仅对上部结构

27、即2 节点集 中质量模型建立ABAQUS刚性基底边界模型，只要 在ABAQUS结构基底处输入由本方法计算得到的刚 性基础处反应加速度时程（水平及摆动分量），再由 ABAQUS 刚性基底边界模型计算其上部结构最高点 时程反应即可，其结果也绘于图3。从图3 中可以看 到，由本文直接计算结果及ABAQUS刚性基底边界 上部结构模型计算的反应吻合的非常好，即同样基 础动力输入下，两算法结构反应计算一致，进一步 说明计算是合理可靠的。此例由本方法计算所得刚 性基础处反应加速度时程（水平及摆动分量）即为 仅含上部结构模型的等效输入，即两步法的中第2 步的复杂体系地震动输入。本例还说明了，当上部 结构模型一

28、致时，两种算法结果一致，可验证等效 输入概念的可行性。540土 木 工程 学报2010 年参 考文 献1刘晶波，谷音，李彬. 结构地基开放系统动力相互作用问题的高效解法J. 土木工程学报，2006, 39（5）:112- 116 （Liu JingBo, Gu Yin, Li Bin. An efficient method for dynamic interaction of open structure-foundation systemsJ. China Civil Engineering Journal , 2006 , 39（5）: 112- 116（in Chinese）2 荣峰，

29、汪家春，何树延，等. CARR 堆厂房土壤结构相互作用与楼层反应谱分析J. 核动力学报, 2006, 27（5）: 19- 23（Rong Feng, Wang Jiachun, He Shuyan, et al. Analysis of soil-structure interaction and floor response spectrum of reactor building for china advanced research reactorJ. Nuclear Power Engineering, 2006, 27（5）:19- 23（in Chinese）34廖振鹏. 工程波

30、动理论M. 北京: 科学出版社，2002Berger E , Lysmer J , Seed H B . FLUSH a computer program for seismic response analysis of axisymmetric soil-structure systemR. Berkeley: University of California,1975Lysmer J, Symes B, TabataBaie-Raissi M, et al. SASSI/S analysis of soil-structure problemsR.Berkeley: University of California, 1981Wolf J P. Soil structure in the time domainM. NJ : Prentice-Hall, 1988需要说明的是图5 仅为特定场地条件及特定时程情况下的单一计算结果，并不能反应地震记录及 反应谱特征，对不同地基参数及结构参数在各类地 震记录下的土-结构地震反应研究，以及反应谱特性将是进一步需考虑的工作。563结论7 唐敢，陈少林，王法武，等. 空间结构-地基动力相互作用的三维时域数值分析方法J. 工程力学, 2009, 26（8）: 143- 149（Tang Gan, Chen Shaolin, Wang Faw