（工程力学专业论文）基于ABAQUS软件的地下洞室群抗震分析研究.pdf

摘要地下结构的抗震研究是当今地震工程界的重要研究方向之一，其研究方法与地面结构相比有较大的不同。近年来，随着地下结构震害的频繁发生，地下结构抗震问题日益受到世界各国地震工作者的重视，并取得了一定成就。目前的抗震研究绝大部分仍集中在地面建筑，对地下结构的抗震的研究也限于地铁、生命线工程等线形结构，对大断面的大型地下结构尤其是大型地下厂房等地下工程的抗震问题探索不足。本文首先介绍目前地下结构抗震分析的主要方法及存在的问题，在借鉴地下结构拟静力方法思想的基础上，应用a b a q u s 大型有限元软件实现了一套以拟静力设计思想为基础且适应于大型地下厂房结构抗震设计使用的分析方法。其次，详细阐述了文中拟静力方法的实施步骤与特点：先利用a b a q u s 软件的动力时程分析结合无限元边界进行一维岩体的地震动响应分析，然后基于前人提出的水平等效惯性加速度的求解方法，进行地下结构在地震作用下的全过程模拟分析，给出地下结构在输入地震波作用下的峰值变形与应力。本文在a b a q u s 平台上实现的拟静力方法概念清晰、操作简单，可以避免对岩( 土) - 结构整体模型进行复杂的动力相互作用分析。最后，结合实际工程民治水电站大型地下洞室群与基于a b a q u s 软件的动力响应谱分析进行对比研究，验证本文在a b a q u s 平台上实现的地下结构拟静力分析方法的可靠性及较好的模拟精度，给实际工程的抗震分析和设计提供参考。【关键词】抗震分析；地下洞室；拟静力法；a b a q u s ；动力分析a b s t r a c ts e i s m i cr e s e a r c ho ft h eu n d e 唱r o u n ds t r u c t u r e si st h em a i nr e s e a r c ho r i e n t a t i o ni nt h es e i s m i ce n g i n e e r i n gc i r c l e t h er e s e a r c hm e t h o d sa r eq u i t ed i 仃- e r e n tf 而mt h em e t h o d si nr e s e a r c h i n gt h ea b o v e g m u n ds t r u c t u r e s r e c e n t l y ，s e i s m i cr e s p o n s ea n a l y s e so fu n d e 唱r o u n ds h u c t u r e sa r ef o c u s e do nw i t ht h e 仔e q u e n to c c u r r e n c e so fs e i s m i cd a m a g e sa n da l s oo b t a i nc e r t a i na c h i e v e m e n t s n o v 吒t 1 1 es e i s m i cr e s p o n s er e s e a r c h e sa r em a i n l yf o c u s e do nt h ea b o v e g r o u n ds t l l l c t u r e s t h er e s e a r c h e so fu n d e 唱，o u n ds t m c t u r e sa r er e s t r i c t e dt ot h el i n e a rs t n l c t u r e ss u c ha sm e t r o ，l i 佬l i n ee n g i n e e r i n ge t c b u tt h es e i s m i ca n a l y s i so fl a 唱es e c t i o n ，e s p e c i a l l yl a r g eu n d e i ? g r o u n dp o w e r h o u s e si sd e 氍c i e n c y i nt h i sp a p e r ，m a i np r o b l e m so ft h em e t h o d sc u r r e n t l yu s e df o rs e i s m i cr e s p o n s ea n a l y s i so fu n d e f 伊o u n ds t t l l c t u r e sa r ed i s c u s s e d b yu s i n gt h el a 唱e6 n i t ee l e m e n ts o 胁a r ea b a q u st or e a l i z et h ep s e u d o - s t a t i cm e t h o d sf o rs e i s m i ca n a l y s i so fu n d e 唱r o u n ds t m c t u r e s ，a n dt h i sm e t h o dc a nb eu s e dt oa n a l y s i sa n dd e s i g nt h eu n d e r g r o u n ds 咖c t u r e s t h e nt h ei m p l e m e n t a t i o np r o c e d u r ea n ds p e c i a lf e a t u r e so ft h i sm e t h o da r ei n t r o d u c e db ys o m ed e t a i l s t h em e t h o do fd y n a m i c a lt i m e - h i s t o ua n a l y s i si na b a q u sc o m b i n e dw i t hi n f i n i t ee l e m e n tb o u n ( 1 a 叫i su s e dt o 锄a l y s i st h es e i s m i cr e s p o n s eo fo n e - d i m e n s i o n a lr o c km a s s t h e nb a s e do nt h es o l u t i o nm e t h o do fe q u i v a l e n th o r i z o n t a l i n e r t i a la c c e l e r a t i o ni nr e f s【1 9 】c o m p l e t es e i s m i ca n a i y s i so fu n d e 唱r o u n ds t r u c t l l r e si sc a r r i e do u t f i n a l l y ，t h ep e a kd e f o n n a t i o na n di n t e m a if o r c eo fu n d e 唱r o u n ds t m c t u r e sc 卸b eo b t a i n e d t h ep s e u d o s t a t i cm e t h o d sb a s e do na b a q u si sac l e a rc o n c e p ta n de a s yf o ra p p l i c a t i o n u s i n gt h i sm e t h o dw ec a na v o i dt h ec o m p l i c a t e da n a l y s i so fd y n a m i cs o i l - s t m c t u r ei n t e r a c t i o n s f i n a l l y ；i nt h i sp a p e r b yc o m p a r i n gt h ec a l c u l a t i o no ft h ep r a c t i c a le n g i n e e r i n g - m i n g z h ih y d r o p o w e rs t a t i o nw i t ht h er e s u l t sb a s e do na b a q u s ，w ec a nv e r i 矽t h a tt h eu n d e 略r o u n ds t r u c t u r ep s e u d o s t a t i cm e t h o db a s e do na b a q u si sr e l i a b l ea n dt h eb e t t e rs i m u l a t i o np r e c i s i o nc a nb ea c h i e v e d s ot h i sm e t h o dc a no f r e rr e 诧r e n c ef o rs e i s m i ca n a l y s i sa n dd e s i g no fp r a c t i c a le n g l n e e n n g 【k e y w o r d s 】：1 1 1 l d e 靴u i l dc a v 哪s ；q u a k e p r o o f 柚a l y s i s ；p s e u d o s t a t i cm e t h o d ；a b a q u s ：d y n 锄i ca n a l y s i 3学位论文独创性声明：本人所呈交的学位论文是我个人在导师指导下进行的研究工作及取得的研究成果。尽我所知，除了文中特别加以标注和致谢的地方外，论文中不包含其他人已经发表或撰写过的研究成果。与我一同工作的同事对本研究所做的任何贡献均已在论文中作了明确的说明并表示了谢意。如不实，本人负全部责任。论文作者c 签名) 壅盔7 9 诸年占月序日学位论文使用授权说明河海大学、中国科学技术信息研究所、国家图书馆、中国学术期刊( 光盘版) 电子杂志社有权保留本人所送交学位论文的复印件或电子文档，可以采用影印、缩印或其他复制手段保存论文。本人电子文档的内容和纸质论文的内容相一致。除在保密期内的保密论文外，允许论文被查阅和借阅。论文全部或部分内容的公布( 包括刊登) 授权河海大学研究生院办理。论文作者( 签名) ：乏叁蠡懈年6月占日第一章绪论第一章绪论1 1 概述随着科学技术和人类文明的发展，地下人工洞室的发展和现在地下空间利用的形态已是千姿百态，远远超出了为个人生活服务的利用领域，而扩大到为了保持作为集团的居民的生活需要空间。由于世界范围内能源及其他原材料的短缺，尤其是基于环境保护的需要，近年来，世界许多国家也开始日益重视地下空间的开发和利用。国际上普遍认为：2 1 世纪将是人类开发和利用“地下空间”的世纪【1 1 。当前，地下设施已经成为现代社会下部基础构造的主要部分，被广泛应用到地铁、公路、原材料储藏库以及输水发电系统等方面，位于地震活动较频繁的区域中的地下结构设施必须能够经受地震荷载与静力荷载的双重考验【z j 。随着工农业生产的发展和城市化程度的不断提高，地下结构的重要性也日益明显。尤其对大型和超大型地下洞室的设计研究日益迫切，在9 0 年代后，我国动工的峡谷高坝基本上都采用地下厂房设计方案，如三峡、二滩、溪洛渡等。在建的地下厂房中，三峡水电站尺寸为3 5 2 6x3 1 6x8 8 6 2 m ，为世界之最。在水电站地下厂房系统中，除地下厂房外，还包括主变洞、交变洞、调压井、尾水洞等大型洞室。随着一些干流水利资源的相继开发，水电工程中地下厂房及相关洞室等在向大型化或超大型化的方向发展，特别是伴随着装机容量的增大，地下厂房规模越来越大，以大跨度、高边墙、多交叉以及结构复杂为特征。在很多地下厂房系统中，大大小小的隧洞和地下空间交织在一起，形成庞大而复杂的地下洞室群p 】。地下结构的抗震研究是当今地震工程界的重要研究方向。地下结构抗震的研究方法与地面结构相比有较大的不同。近年来，随着地下结构震害的频繁发生，地下结构抗震问题日益受到世界各国地震工作者的重视，并取得了一定成就。目前的抗震研究绝大部分仍集中在地面建筑，对地下结构的抗震的研究也限于地铁、生命线工程等线形结构，对大断面的大型地下结构尤其是大型地下厂房等地下工程的抗震问题探索不足阡】。1 2 地下空间结构的震害现象在国内，1 9 7 6 年7 月2 8 日，唐山发生了7 8 级地震，一些地下通道、煤矿巷道和l河海大学硕士学位论文人防工程产生了轻微破坏，如在9 度区的天津宁河县和汉沽区，地下人防工程尤其是通道部分普遍出现有规律的l 3c m 环向裂缝，在接头转角处多发生断裂并有所错动造成漏水，个别底部有喷砂冒水和局部倒塌现象，在8 度区的天津塘沽地区，地下人防工程基本完好，但不少通道也出现许多环向裂缝。工程主体的裂缝多发生在丁字接头、拐角、出入口、不同结构的交接入口和断面变化部位。1 9 9 9 年9 月2 1 日台湾发生7 3 级地震，台中地区5 7 座山岭隧道中有4 9 座遭受不同程度的损坏，如衬砌开裂、剥落、移位笙【5 ，6 】寸。在国外，1 9 0 6 年美国旧金山发生8 3 级大地震，严重损坏了横穿圣安德烈斯断层的两座隧道，损坏主要是木梁柱折断，轨道上凸，奈特一号隧洞水平错动1 3 7 m 。1 9 2 3年日本关东地区发生8 3 级大地震，城市地下管网损坏严重，震区铁路线上1 1 6 座隧道有8 2 座遭到破坏或发生变形、边墙塌落、衬砌开裂等。东京附近距震中较近的2 4 座隧道产生了严重的损坏，损害主要是由地震振动所引起，隧道未穿越断层，破坏情况主要为拱部和边墙坍落、衬砌裂缝和变形错动，以及洞门砖石墙碎裂。1 9 3 0 年日本伊豆地震，丹那断裂带再次活动造成丹那隧道水平错动2 3 9 m ，竖向错动0 6 m 。1 9 5 2 年，美国的克恩县地震，使南太平洋铁路上的4 座隧道遭到了严重破坏。1 9 8 5 年墨西哥地震中，建在软弱地基一卜的地铁侧墙与地表结构相交部分发生分离破坏现象。特别是1 9 9 5年1 月1 7 日在日本兵库县南部发生里氏7 2 级“阪神地震”，这次地震首次出现地铁主体结构破坏。神户地铁的大开车站震后出现9 5 m 长的塌陷，其上方最大沉陷3 m ，有约3 0 个钢筋混凝土柱和3 5 个支撑平台倒塌：在向长田车站方向5 3 0 m 隧道中，有2 2 0 根立柱损坏，在神户的三宫站，约2 0 个钢筋混凝土柱破坏。新干线中有两个隧道严重破坏，六甲隧道在新大阪和新神户站约1 6 5 k m 之间，裂缝长达8 k m ，六甲隧道东端1 8 0 0 个栈道支撑墩中有4 0 0 个破坏，其中1 0 0 个倒塌【7 一。尽管地下结构具有良好的抗震性，在地震作用下结构处于不同的地质条件，采用不同的设计、施工方法，地下结构的抗震性能也各不相同。以上这些地震都使不少地下结构遭受震害。地震波动使地下结构开裂、坍塌和破坏。我国西部强震地区将有多座巨型水电站即将修建，其地下厂房受地形和地质条件的限制，所处的地质环境非常复杂，这些工程一旦遭受震害，不但将造成工程本身的损失，而且会带来巨大的次生灾害，后果将十分严重。因此研究超大型地下洞室群的抗震安全性十分重要。2第一章绪论1 3 地下结构抗震设计的研究现状1 3 1 地下结构的动力特性及研究方法从以往的震害报道中可以看出，地下结构与地面结构的振动特性有很大的不剐9 】：1 、地下结构的振动变形受周围地基土壤的约束作用显著，结构的动力反应一般不明显表观出自振特性的影响：2 、地下结构的存在对周罔地基震动的影响一般很小( 指地下结构的尺寸相对于地震波长的比例较小的情况) ；3 、地下结构的振动形态受地震波入射方向的影响很大，地震波的入射方向发生不大的变化，地下结构各点的变形和应力可以发生较大的变化；4 、地下结构在振动中各点的相位差别十分明显，而地面结构各点在振动中的相位差不很明显；5 、地下结构在振动中的应变一般与地震加速度的大小联系不很明显；6 、地下结构的地震反应随埋深发生的变化不很明显；7 、对地下结构和地面结构来说，它们与地基的相互作用都对它们的动力反应产生重要影响，但影响的方式和影响的程度则是不相同的。地下结构和地面结构动力反应特点的不同，决定了它们抗震分析方法的不同。地下结构抗震问题的研究方法大致分为三种【l o 】：地震观测、实验研究和理论分析。地震观测就是通过实测地下结构在地震时的动力特性来了解地下结构的地震特点。实验研究分为人工震源实验和振动台实验、人工震源实验为实地研究结构动力特性或求得地基弹簧阻尼特性等，进行现场激振试验。由于人工起振力偏小，很难真实地反映出建筑物的非线性性质和地基断裂等因素对地下结构地震反应的影响，所以一般不多采用。振动台实验法能够较好地把握地下结构的地震反应特性以及地下结构与地基之间的相互作用特性等问题，因此更受重视。通过模型实验使人们能更好地了解和掌握地下结构的工作特性，为抗震理论的发展奠定了基础。地下结构抗震设计计算方法，从力学特性上可分为拟静力计算方法和动力反应分析方法两类。由于地震是动力作用，显然动力反应分析方法应是理想的抗震设计方法。地震动力反应分析( 时程分析法) 和动力模型试验，能详细研究地震对地下结构的振动作用，但需要详细掌握有关的参数，如地震加速度、地层的动杨氏模量和阻尼系数等。借助于动力有限元法可以比较埋想的解决地下结构抗震问题，应用此方法时，除了需要解决有限元理论中的基本问题外，还要解决模拟振动过程中岩( 土) 体剪切模量的变化规律，模拟岩( 土) 体阻尼、边界问题这三个基本问题1 1 1 1 。3河海大学硕士学位论文但是，动力反应分析方法消耗计算机资源大，计算时间长，更重要的是动力反应分析需要高度的多方面专业知识和技能，对技术人员素质要求很高，且操作繁杂，对计算结果的评价也不容易。因此，目前国际上通常使用拟静力计算方法。我国现行地铁设计规范中关于抗震设计部分，也是规定基本采用反应位移法或地震系数法的拟静力计算方法进行分析【1 2 1 。1 3 2 地下结构抗震设计中的拟静力方法根据文献【1 0 】中所述，最早形成拟静力计算方法雏形是在上世纪6 0 年代，6 0 年代初，前苏联学者在抗震研究中将弹性理论应用于地下结构，以此求解均匀介质中关于单连通和多连通域的应力应变状态，得出了地下结构地震力的精确解和近似解。6 0 年代末，美国在建设快速地铁运输系统( b a 蛐时，对地下结构抗震进行了较深入的研究，提出了地下结构并不抵御惯性力而是具有强加变形的延性，同时还不散失其承受静荷载力等新的设计思想。到了7 0 年代后期和8 0 年代，日本学者依据地震观测和模型试验资料，结合波动理论，提出了反应位移法、应变传递法等实用的拟静力计算方法，大大丰富了地下结构的抗震理论。随着地下结构抗震理论的不断发展和工程实践的不断进步，在地下结构抗震设计中最具代表性、最常用的拟静力计算方法是反应位移法和有限元反应加速度法。但是对于大型地下洞室而言，这些拟静力计算方法仍然应用极少。水电站大型地下洞室群多位于围岩介质中，这与以往的地下结构抗震多考虑土一结构相互作用的分析模型有很大区别。由于岩体是赋存于自然界中的十分复杂的介质，岩体呈现明显的非线性、不连续性、不均质性和各向异性等复杂特性，因此，如何改进和完善拟静力计算方法，使其既保持操作上的简便性，又能够在较复杂的岩石介质条件下亦具有足够理想的计算精度，是大型地下洞室抗震设计中需要研究的问题。一、s t j o h n 法【1 3 】该法以弹性地基梁模型来考虑土一结构的相互作用问题，但忽略了土体与结构之间的动力相互作用，是一种拟静力分析力方法。该法认为在地震荷载作用下，隧道截面产生与自由场的轴向、弯曲和剪切变形相对应的轴向、弯曲和剪切应变。l 、地震波传播方向对地下结构产生变形的影响由地震波引起的地下结构承受的地震荷载，按其受载方式可分为3 种。4第一章绪论( 1 ) 弯曲荷载：它是由地震波沿着与隧道轴线平行的方向传播而产生的，其结果导致隧道在纵向和横向平面内的变形。( 2 ) 横向荷载：它是由横向传播的剪切( s h ) 波、压缩( p ) 波、r a y l e i g h 波、l 0 v e波等在隧道外壁产生的动力剪应力和门应力引起的。( 3 ) 轴向荷载：它是由平行于隧道轴向的质点运动产生的。s t j o h n 法在是否考虑土与结构的相互作用时，引入了p e c k 教授提出的柔度比概念，即：如果柔度比f 2 0 ，则认为衬砌是完全柔性的，土体与结构不发生相互作用，地下结构屈从于周围介质一起运动；如果柔度比f ( 1 2 2 )式中，【k 】为刚度矩阵；【u 】为结点位移列阵；【用为等效到结点的外荷载列阵。进行动力分析时，由于阻尼力的作用非常复杂，一般用粘性阻尼力来表示体系中的阻尼力，这种阻尼力和运动速度成正比。单元的体积力可以写成： p = p ) 一刁 矗) 一p 甜)( 1 - 2 3 )式中， p ) 为原有的体积力；7 为阻尼系数；p 为岩体密度。根据有限元的插值关系，单元中任意一点位移 “) = 【( 五j ，z ) 】 u ( f ) ) 。；【】为形1 0第一章绪论函数，是坐标( x ，y ，z ) 的函数，单元中任意一点的速度和加速度可以表示为： 矗) = 【】 d ) 。代入( 1 2 3 ) 式，得：= 【】 矽) 。( 1 2 4 )( 1 2 5 ) p ) = p ) - 巩】 【7 ) 8 一烈】 矽) 。( 1 - 2 6 )将 p ) 按照有限元中体力等效到结点的公式，可得到 p ) 的等效结点荷载，对于一个单元而言 p ) 8 = 肌r p 抄+ q ) 。+ f ) 8( 1 - 2 7 )式中， 奶。为作用在单元表面的等效结点荷载； f ) 8 为集中力的等效结点荷载。将( 1 2 7 ) 式代入( 1 2 6 ) 式中，得 p ) 8 = r ) 。一【c 】。 u ) 。一【m 】。 u ) 。( 1 2 8 )式中， r ) 。= j m r p ) 西+ q ) 。+ f ) 8( 1 2 9 )匕式代表作用于单元体的外荷载等效结点荷载列阵，而阿= j j 【】7 叩【协【m 】8 = 肌】r 烈】咖( 1 3 0 )( 1 3 1 )分别为体系的质量矩阵和阻尼矩阵，将式( 卜2 3 ) 代入式( 卜2 9 ) ，得到单元的运动方程为：【j l 彳】8 u ) 。+ 【c 】8 u ) 。+ 【k 】。 u ) 。= 尺 。( 1 3 2 )将各单元的运动方程按照常规的方法组装在一起，就可以得到结构的整体运动方程，为：【m 】 矽) + 【c 】 d ) + 【k 】 u ) = r )( 1 3 3 )式中，【肘】、【c 】、陋】分别为结构的质量、阻尼和刚度矩阵。初始条件为，= 乇， u ) = 砜) ， u ) = u 。)( 1 3 4 )1 1河海大学硕士学位论文而在地震时，没有外力直接作用在各质点上，地震力使体系所固着的地而发生一定的运动，经过系统将地震力传给各质点，使各质点相对地面进行振动。在地震作用下地面相对于固定坐标系可写出系统的振动方程【m 】 掰) + 【c m ) + 【k m ) = - 【m 】【g 】 ( f ) ( 1 - 3 5 )方程式( 1 3 5 ) 就是多自由度体系在地震地面运动影响下的振动方程。其中g 1 是从地面运动的三维笛卡尔坐标空间到”个自由度体系的刀维空间的转换矩阵。 口。( ，) )是地基边界均匀输入的地震加速度时程。地震波输入是进行抗震分析的关键环节，直接影响计算结果的精度和可信度。现有的输入模型有单一输入和多点输入。单一输入适用于小型结构，对大型结构则精度较差。多点输入需要模拟多点地面运动，而模拟多点地面运动的方法又可分为行波法、随机法和内差法【2 2 1 。地震波输入不能用半无限域模型，必须从半无限的地球介质中切取计算区，用人工边界把计算区与非计算区分离来实现地震波的输入，如图1 1 。因此，需要引入各种人工边界才能得到正确的结构响应【2 3 1 。数人：i ：边弊农卜人l ：边界图1 - 1 地下结构计算模型用有限元法分析地下岩体与结构动力相互作用时，需从半无限的地球介质中切取出感兴趣的有限计算区。在切取的边界上需建立人工边界以模拟连续介质的辐射，从计算区内部穿过人工边界而不发生反射。建立人工边界的方法，可广义地分为两大类2 4 】：精确边界和局部边界。第一类方法使人工边界满足无限介质的场方程、物理边界条件和辐射条件。这类边界在有限元的意义上是精确的，并能设置在不规则构造物和周围介质的界面上。这类精确方法在许多场合是有效的，然而，除了其它局限性外，这类方法最大的缺点是使人工边界上所有结点运动耦联，这将导致对计算机存储量提出更高要求并第一章绪论耗费较长时间。而局部边界的显著特征是其良好的实用性，人工边界上任一结点的运动与其它结点( 除邻近结点外) 解耦，因而计算机存储量小，计算时间短。现在比较成熟的人工边界有以下几种：粘性边界，一致边界，叠加边界，旁轴边界，透射边界，动力映射无限元等。文献 2 5 ，2 6 曾对以上的人工边界进行了对比和分析，在这些边界中，粘性边界、旁轴边界、透射边界属时域局部人工边界，而透射边界则具有较高的精度。另外，还可以采用无限元与有限元耦合的方法来模拟无限域或半无限域的平面问题。无限单元最早是由r f u n g l e s s 和r b e t t e s s 提出的，这种单元在1 个或2 个方向上是无限的，主要用来模拟无限域或半无限域的平面问题。其形函数不是原来的简单分片多项式插值函数，而是沿无限方向延伸的半解析函数。在无限元坐标形函数确定之后，单元刚度矩阵的形成过程与有限元相刚z 7 】。无限元方法的突出优点是无需涉及解析解表达式，使无限元成为有限元的一个部分，可以给出统一的求解格式。实际问题的边界条件往往比较复杂，难以给出解析解或者基本解，所以有限元与无限元的结合方法成为岩土工程计算中最有效的方法之一【2 8 】。本文在进行动力分析时就采用了无限元与有限元联合方法模拟半无限的地下岩体。1 4 问题的提出在实际地下结构工程的抗震计算中，拟静力法得到广泛应用，而作为其核心步骤的一维岩( 土) 层地震反应分析，在以往文献中多采用等效线性化程序shake91 、ee ra 、rslnlm 等完成。在完成后，进行二维有限元反应分析计算则多使用s u p e r f l u s h 、t d a p - 3 或msc marc 等软件。这样在各个步骤之间计算数据的转换应用会增加相当繁复的工作量，十分不便。而a b a q u s 大型有限元软件由美国h k s 公司( h i b b i t t ，k a r l s s o n & s o r e n s e n ，i n c ) 研制开发，在国际上被公认为是功能最强的有限元软件之一，可以分析各种固体力学、结构力学系统，特别是能够处理非常复杂的问题和模拟高度非线性问题。近些年来，a b a q u s 优秀的分析能力和模拟复杂系统的可靠性使得a b a q u s 被各国的工业和研究机构广泛采用。自从1 9 9 7 年进入中国以来，越来越多的国内企业和高校也采用a b a q u s 作为产品研发和科学研究的工具，内容涉及机械、材料、土木、电子等专业。为方便用户开发自己需要或感兴趣的模型，a b a q u s 提供了方便灵活的二次开发1 3河海大学硕士学位论文平台，包括若干用户子程序( u s e rs u b r o u t i n e s ) 以及在编程时可以调用的实用程序( u t i l 时r o u t i n e s ) 。a b a q u s 中的用户材料子程序为用户提供了自定义材料本构模型的程序接口，用户可通过f o r t r a n 或者v c 编程来实现所要解决的工程科学问题中所需要的特定材料模型【2 9 1 。在此基础上笔者应用大型有限元分析软件一a b a q u s ，结合结构抗震设计的拟静力计算方法，尝试在统一的a b a q u s 平台下实现其各个步骤，并且应用于超大型地下洞室的抗震安全性研究。1 5 本文的主要工作本文的核心内容是在借鉴地下结构抗震设计的拟静力计算方法思想的基础上，应用a b a q u s 大型有限元分析软件并结合实际工程民治水电站地下厂房洞室群进行计算，在统一的a b a q u s 平台下实现拟静力分析的各个步骤。围绕这一点主要做了以下工作：1 、应用无限元与有限元耦合的方法建立岩体自由场模型，给定岩体材料参数，在a b a q u s 中实现输入地震波作用下的一维岩体地震反应分析。2 、在完成一维地震反应分析后，采用文献 1 9 中的加速度求解方法来得到自由场岩体沿深度分布的水平等效惯性加速度。3 、建立地下洞室体系有限元模型，首先模拟洞室开挖，得到开挖后的地下洞室区域应力及位移分布( 在此步骤中对于岩体弱化的处理没有采用a b a q u s 传统的场变量控制材料参数变化的方法，而是应用单元的添加和移除功能进行模拟，获得了良好的分析效果) ，接着按照土层所在位置逐步施加第2 步中得到的沿深度分布的水平惯性加速度，然后按照静力有限元方法进行求解。值得指出，本文实现了拟静力计算方法各个步骤之间的数据在统一平台下的链接，省去了数据在不同软件间的转换，从而大大简化其操作过程。4 、基于a b a q u s 分析软件，实现基于响应谱方法的大型地下洞室抗震动力分析，并且进行地下结构拟静力分析方法和基于响应谱动力分析方法结果的对比研究，以验证本文提出的将地下结构拟静力分析方法应用于大型地下洞室的抗震设计的合理性，并且针对分析结果，对深埋地下厂房洞室群的抗震安全以及加固措施进行研究。1 4第二章结构动力分析的基本理论第二章结构动力分析的基本理论无论是一维岩体动力分析还是对于整个地下洞室群的地震动分析，都要借助于结构的有限元动力分析，本章就针对有限元的动力分析理论并结合其在a b a q u s 中的具体应用做进一步的研究。2 1结构动力学的基本方程有限单元法可以分析结构振动问题以及动力响应问题，即在动荷载作用下结构的动力响应，如位移、应力等。动力学问题的有限元也同结构静力学问题一样，要把物体离散为有限个数的单元体。不过此时在考虑单元特性时，结构所受到的荷载还要考虑单元的惯性力和阻尼力等因素。采用有限元分析结构的动力响应时，首先要建立有限元运动方程。建立有限元运动方程式要用到动态问题的变分原理，有势能原理、余能原理、广义变分原理等。就结构动力分析而言，采用与位移有关的能量变分原理是合适的，下面主要介绍一下哈密顿( h a m 订t o n ) 原理。该原理叙述如下：在满足协调性条件、约束条件或运动边界条件以及在时间t ，与t ：的条件下所有可能的位移随时间变化的形式中，是真实解的那种变化形式使拉格朗日泛函数取其驻值3 伪。这一理论中的拉格朗日泛函数定义为m 1 ：工= 丁一u 一一形( 2 - 1 )式中，丁为物体的动能；u 为物体的应变能；为阻尼力势能；形为外力势能。丁= 睁p y p 是质量密度即单位体积质量，p ) 是位移列向量p ) = 缸，v w rp j 是速度列向量，表示位移对时间的一阶导数俗) = 弛，帚，矽1 s( 2 2 )( 2 3 )( 2 4 )河海大学硕士学位论文在动力问题中，位移甜、1 ，、w 均为时间f 的函数u = 胪 ，y ，w 谚y( 2 5 )其中，么是应变能密度，即彳= 三p y p ) = 三p ) r 【d 怡)( 2 - 6 )式中，p ) 为结构应变列阵，p ) = k ，q ，占：，比) r ；p ) 为结构应力列阵，p ) = 蠢，盯y ，吒，p ；p 】为结构弹性矩阵。【。】= 揣lj lj lo0l 一l 一1 上00l 一loo生o2 ( 1 一)生o2 ( 1 一)1 2 2 ( 1 一)( 2 - 7 )其中，e 是材料的弹性模量，是材料的泊松比。是阻尼力势能：吸2 睁侈yp 炒( 2 - 8 )式中，c 为粘性阻尼系数。形是外力势能，它包括体积力势能和表面力势能，体积力势能为：= 胎) r 依沙y( 2 - 9 )式中，识) 是体积力列向量，即 e ) = 岁，】，z ) ，其中，ky 、z 分别为物体体积矿域内沿坐标x 、y 、z 方向单位体积的体积力。表面力势能为：1 6第二章结构动力分析的基本理论形：= 胎厂以洒( 2 - 1 0 )式中，以) 是表面力列向量，识) = 讧，】，z y ，其中，瓜kz 分别为物体表面s上沿三个坐标轴方向单位面积的表面力分量。由哈密顿原理知，使拉格朗日泛函数为极小的位移才是真实的，所以有：万l 三衍= o( 2 1 1 )先建立一个单元的运动微分方程式，等参元分析中有：= 【 ；p ) = 【】侥 ；p ) = 陋) ；p ) = 【d 括)( 2 1 2 )式中，谴) 为结点位移列向量； 为形状函数矩阵； b 为形状函数导数矩阵( 又称应变矩阵) 。对于不同的单元形式，上述矩阵是不完全相同的，将( 2 1 2 ) 代入( 2 一1 ) ，即得到拉格朗日泛函：丢儿帖兢r 防r 【胞) 一概y 【b r 陋p ) 一c 侥厂【r 【 一2 概) 7 【】r 以瘳yf 陋) 7 【】r 把洒( 2 一1 3 )应用哈密顿原理，在时间区间 t 。，t 。 上对三积分，并使其变分等于零，考虑到矩阵 d 的对称性后，有j ：锄= 小蚶 ( 妒r 【d p p y p g 侄h 眇m f p y p拗妙m p 少7 1 r 陆均一p 鼢 ( 妒】7 舡泅肛。将第二项应用分部积分公式，有1 7河海大学硕士学位论文j ：2g 垃玎眇时【p y 侥妇：陋 ( 妙眦一r ( 妙2 。1 5 按照哈密尔顿原理，上式中第一项的值为零，因为万慨“) ) = 万 疋o ：) ) = o ，于是只剩下第二项。同样，对于第三项有：f6 镜y 1 肝【】7 【p yl 概冲= 一j ：2p 缱) 7 l 肚【r 【p 矿l 经协( 2 1 6 )yr引入下列符号：医。】= 胪】7 p p 抄阻。】= 肛【】r 【妙矿【c 。】= p 陋】r 【p y忸。) = p r 亿扫矿+ 肛】r 以扫y( 2 1 7 )矿s式中，医。】为单元刚度矩阵；阻。】为单元质量矩阵；k 】为单元阻尼矩阵；忸。) 为瞬变的结点力列向量。则，式( 2 1 4 ) 变成：j ：2p 谴) 7 姚。) + 阻。琏) + 【c 。】侥) - 他妇= o( 2 1 8 )由于单元结点位移的变分万缱y 是任意的，故得到动态中单元的运动方程式为：阻。】髋) + 虹】硅) + k 。) = 亿)( 2 1 9 )结构整体的运动方程式，可以由单元的运动方程式作为基础进行组装得到：阻。】经) + 匠】娩j + k 。) = f r 。)( 2 2 0 )式中，般) 为结构整体所有结点位移的列向量；棱 为结构整体所有结点位移对时间一阶导数的列向量；概 为结构整体所有结点位移对时间二阶导数的列向量： m 为1 8第二章结构动力分析的基本理论结构整体的质量矩阵； c 为结构整体的阻尼矩阵； 幻为结构整体的刚度矩阵：忸。)为结构所有结点的瞬变的结点力列向量。其中阻尼矩阵 c 通常用瑞利阻尼表示，如下式：【c 】= 口】+ k 】( 2 2 1 )上式中系数口及由以下两式联立确定：毒：旦+ 堕磊2i + 了缶：芸+ 譬q 以2 式中参，白是结构第j 、阶振型的阻尼比；峨、彩，是第j 、阶振型的角频率。应用哈密尔顿原理推出了结构的运动微分方程，它是一个二阶微分方程组，反映了动力问题的两个方面3 劲。( 1 ) 方程的右端俅。( f ) ) ，表示外力对结构的激励，是随时间变化的函数，对特定的某个结构的动力问题，它是一个已知函数。于是，可在初始条件下求出方程的解髋( f ) ) ，即求解结构的动力响应问题。( 2 ) 对方程的左端，若右端项扭，( f ) ) = o ，在非零的初始条件下，方程也有非零解，这时结构处于自由振动状态。由于没有外载作用，方程的解反映了结构本身固有的特性，这就是结构的固有频率与振型，在数学问题上称之为特征值与特征向量，于是就导出了结构的动力特征值问题。2 2 结构自振分析2 2 1 自振频率和振型求解结构的频率及振型在工程的实际应用以及在求解动力响应方面，都具有很重要的意义，因而它在结构动力学中占有很重要的地位。在实际工程中，由于阻尼对结构自振频率及固有振型影响不大，因而在讨论结构的固有特性时，常不计阻尼作用，于是方程变为3 3 】：1 9河海大学硕士学位论文舷) + k ) = o )( 2 2 3 )此时体系将做简谐运动，有：慨 = 移) s i n ( 刎+ 口)( 2 2 4 )戗j = 一缈2 劬) s j n ( 刎+ 秒)( 2 - 2 5 )式中，彩一圆频率；秒一初始相位角；移) _ 一与时间t 无关的位移向量。将以上两式代入( 2 2 3 ) ，得：僻卜国2 阻蜘) = o )( 2 - 2 6 )或写成k 桫) = a 阻协)( 2 - 2 7 )其中，彳= 缈2 。上式就是结构动力分析中的广义特征值问题。它的核心是求解满足上式的元和非零解向量移) 。很显然，由上式求出的兄和移) 值，只取决于结构本身的刚度矩阵k 】和质量矩阵瞰】，即它们是结构的固有值。万= 万就是结构自振的圆频率，名称为结构的特征值，与万相应的空间振动形态( 即振型或模态) 称为特征向量。在三维动力分析中，质量矩阵阻】为一对称正定矩阵，而刚度矩阵k 】亦为对称正定矩阵，因此，所有特征值都是正( 或非负) 的实数。通常，在广义特征值问题中，由方程( 2 2 6 ) 可得到n 个特征值 ，如，以，即可得到结构的n 个自振频率。这些自振频率从小到大排列可形成频率向量，其中最小的频率叫做基本频率或第一频率。其相应的特征向量溉) ，钕) ，娩) 即为结构的n 阶振型。结构的刀个自振圆频率由小到大排列起来，依次记为q 吐 鸭 缈。其中功。称为基本频率，相应的 办) 称为基本振型。对于每个自振频率，由式( 2 - 2 6 ) 可确定一组各结点的振幅值 谚) = 谚。，谚：，丸 t ，它们构成一个向量，在工程上通常称为结构的振型。在每个振型中，各结点的振幅是相对的，其绝对值可取任意数值，在实际工作中，常用下列方法来决定振型的具体数值：( 1 ) 规准化振型：取 谚) 的某一项，例如取第刀项为1 ，即丸= l ，于是 谚) = 谚。，谚：，1 】t( 2 ) 正则化振型：设已求得一振型 荔 = 【荔。，石：，无】t ，令第二章结构动力分析的基本理论办= 无c ，c = 【 万) 7 【m 】砺) 】l ，2 ，则得到的 谚 = 【谚。，谚：，吮】t 为正则化振型。2 2 2 特征值的求解方法在求解特征值问题时，大致有两种情况，一种是求解结构系统特征方程的全部特征值问题，即所有的特征值和对应的特征向量：另一种是求解部分特征值问题，即部分( 通常是最小的一些) 特征值和对应的特征向量。这是因为，在结构动力学中，往往矩阵的阶数都很高，有时不可能也没有必要求解全部特征值和特征向量。依据这个原则，在求解方法上也分为两大类：一是直接求解法( d i r e c tm e t h o d s ) ，另一类是向量迭代求解法( i t e r a t i v em e t h o d s ) 。直接求解法包括h o u s e h o l d e rq r 和g i v e n s 两种求解技术。这些方法可以用来研究需要求解所有特征值和特征向量的自由度较少的系统，计算工作量较大。计算过程中结构的特征矩阵可以被转换为一个对称的对角阵形式，可以很容易地求解出所有特征值，向量迭代方法包括反向迭代法( 1 n v e r s ei t e r a t i o n ) 、子空间迭代法( s u b s p a c ei t e r a t i o n ) 等。使用迭代法可以求解较大结构系统的少数特征值问题，计算时间依据结构自由度的大小和需要提取的结构特征值的个数而定。在a b a q u s 中提供了两种特征值的求解器：子空间迭代求解器和l a n c z o s 特征值求解器。下面分别对两种求解器的原理和方法做一下介绍：1 、子空间迭代法求解3 钔子空间迭代法实质是用一组假设的初始特征向量将一个高阶特征值问题转化为一个低阶特征值问题，反复进行迭代求解。n 个自由度系统的运动方程：【d ) = 五 矽)式中，【d 】= 【k 】- l 【m 】，五2 吉选取s ( s ，+ 出= t “，+ t 幻，f + ( 丢一口) t 历l + 口t 舀，+ ， ，2c 2 3 5 ，2 7河海大学硕士学位论文其中，口( o ，号) 和万( o ，l 】是由积分精度和稳定性要求来决定的参数，称为权重因子。不同的学者通过定义不同的权重因子口、万的取值，建立了一系列经典的动力学方程时间积分方法。当权重因子万= 1 2 时，n e w m a r k 时间积分方法保持二阶精度，方程( 2 3 4 ) 导出如下形式 西) ，+ ，= 西) ，+ 丢 订) ，+ 丢 订) ，+ ， ，c 2 3 6 ，值得指出的是：n e 硼a r k 时间积分方法中，当万l 2 时算法可以考虑二阶精度，但是当占 ，+ m + ( 1 + 历( 【，】，+ m 一【r 】，+ a r ) 一厦【，l 一【r 1 ) = o( 2 4 8 )将( 2 3 4 ) 、( 2 3 5 ) 代入( 2 4 8 ) ，将权重因子口，万表征为权重因子的函数，如下所示口= 扣俐协4 9 ，万此时任意给定一个权重因子，就可以得到一组相对应的权重因子口，万，然后通过( 2 3 4 ) 、( 2 3 5 ) 式可以导出f + f 时刻的位移和速度。上述积分方法即为h h t时0 刚积分算法，它实际卜是n e w m a r k 方法的发展，这种算法在a b a q u s s t a n d a r d 中得到很好的运用。2 4 小结本章在利用前人有限元的动力分析理论的基础上结合其在a b a q u s 中的具体应用做了详细说明，分别对a b a q u s 软件的结构自振分析和求解运动方程的理论做了阐述，并且对其中不同的处理方法做了比较，总结出每个方法的特点及适用范围。针对本章所作的工作，需要对地下岩体以及地下洞室进行动力分析，在自振分析部分，选用l a n c z o s 向量法进行结构的固有频率和振型的提取，而对于求解结构运动方程，则使用隐式求解方法即选择a b a q u s s t a n d a r d 求解模块进行动力分析计算。第三章水平成层场地地震反应分析方法第三章水平成层场地地震反应分析方法3 1 前言地震的破坏力是震撼人心的。为了减轻地震灾害，人们早在上个世纪三十年代就开始了地震区划的研究，以指导抗