高拱坝流固耦合数值分析方法研究

毕 业 设 计(论 文)高拱坝流固耦合数值分析方法研究专业年级 学 号 姓 名 指导老师 评 阅 人 II二 xx 年六月中国 南京I摘 要在水工抗震计算中，地震作用下的动水压力是水工建筑物必须考虑的重要荷载之一。研究地震作用下水工建筑物与水体的相互作用是水工抗震设计的一个重要课题。本文研究了流固耦合分析中的数值仿真方法，并基于 ANSYS 有限元分析软件，结合天花板水电站拱坝工程建立三维有限元分析模型，进行高拱坝的有限元抗震分析。主要内容如下：1、综述了拱坝的结构特点和发展概况、流固耦合分析以及抗震分析的研究现状。介绍了 ANSYS 软件及其 APDL 语言的二次开发技术、水工结构动力分析方法以及反应谱法和时程分析法在 ANSYS 软件中的实现过程。2、详细论述了流固耦合分析的基本理论及求解方法。重点研究了在ANSYS 软件中考虑动水压力的方法，并对其实现步骤作了说明。以一个二维重力坝水体地基有限元分析模型为例，分析重力坝在空库、死水位、满库三种情况下的自振特性，并且分析坝库系统承受水平位移激励作用时的水体动水压力分布规律。3、结合天花板水电站拱坝工程，建立了拱坝-地基的三维有限元模型，分析了正常蓄水位和死水位下的两种基本静力工况的位移及应力的分布情况；研究了拱坝的动力特性，并采用反应谱分析法和时程分析法分析了正常蓄水期拱坝在地震作用下的动位移以及动应力分布情况，最后对天花板拱坝抗震性能作出安全评价，成果合理，可以为工程设计提供参考。关键词:流固耦合；ANSYS ；拱坝；抗震分析IIAbstractIn the calculation of hydraulic seismic, hydrodynamic pressure under the earthquake is one of important loads must be considered.The research on the interaction between hydraulic structure and water is an important topic in the aseismic design. This paper studied numerical simulation method of fluid-solid coupling,based on ANSYS finite element analysis software, established a project 3-d finite element analysis model for high arch dams combining the ceiling hydropower dam.Main contents are as follows:1.Summarized the structural characteristics and development situation of arch dam survey, fluid-solid coupling analysis and the status quo of seismic analysis . Introduced the ANSYS software and its secondary development technique of language APDL, also introduced the realization process of hydraulic structure dynamic analysis method and response spectrum method and the time-history analysis in the ANSYS software.2.Discusses the basic theory and solving method of fluid-solid coupling analysis in details.Mainly studied considering the hydrodynamic pressure in the ANSYS software and explained its implementation steps.Take a two-dimensional gravity dam water-foundation model for example,analyzing the gravity dams free vibration under the three situation as followed: the empty reservoir ,dead water level and the full reservoir.Analyzed distribution law of hydrodynamic pressure when the dam-library system under horizontal displacement excitation.3.Established a arch dam-foundation three-dimensional finite element model combining with the ceiling hydropower station, and analysed displacement and distribution of stress under the normal storage level and the dead water level.Study the dynamic characteristics of arch dam, and used the response spectrum analysis and time-history analysis to analyze displacement and distribution of stress under the earthquake during the normal storage level period.Finally, made a safety evaluation to ceiling seismic performance.The achievement is reasonable, can be provided reference for the engineering design.IIIKeyword： fluid-solid coupling,ANSYS,arch dam,seismic analysis目录第一章 绪论 - 1 -1.1 选题的背景及意义 .- 1 -1.2 拱坝的类型、特点及其发展 .- 1 -1.2.1 拱坝的类型 - 2 -1.2.2 拱坝的主要特点 - 2 - 1.2.3 拱坝的发展 - 3 -1.3 流固耦合理论研究现状 .- 3 -1.3.1 流固耦合的定义及其分类 - 3 -1.3.2 库水流固耦合理论研究发展过程 - 5 -1.4 结构地震反应分析理论现状 .- 6 -1.5 本文的主要研究内容 .- 7 -第二章 动力分析理论 - 8 -IV2.1 ANSYS 软件简介 - 8 -2.1.1 ANSYS 软件概述 .- 8 -2.1.2 APDL 参数化语言概述 - 8 -2.2 动力分析有限元基本理论 .- 9 -2.2.1 结构动力平衡方程 - 9 -2.2.2 结构自振特性的计算 - 10 -2.3 反应谱分析理论与方法 - 12 -2.4 时程分析法理论与方法 - 16 -2.5 地震波输入与选取 - 19 -2.6 本章小结 - 21 -第三章 流固耦合分析的基本理论 - 22 -3.1 无粘小扰动流动的基本方程和表达形式 .- 22 -3.1.1 无粘小扰动流动的基本方程 .- 22 -3.1.2 以压力 p 为场变量的表达形式 .- 23 -3.2 流固耦合系统有限元分析的格式 .- 23 -3.2.1 流固耦合系统的动力学模型和基本方程及边界条件 - 23 -3.2.2 流固耦合的有限元方程 .- 25 -3.3 两种求解动水压力的方法 - 27 -3.3.1 基于 ANSYS 的附加质量法 - 27 -3.3.2 物理场耦合分析方法 - 31 -3.4 考虑流固耦合的重力坝振动特性研究 - 32 -3.4.1 模型的选取及计算参数的确定 .- 32 -3.4.2 重力坝流固耦合的附加质量求解 .- 33 -3.4.3 基于 ANSYS 的顺序弱耦合分析 - 34 -3.5 本章小结 - 36 -V第四章 天花板拱坝抗震分析 - 37 -4.1 工程简介 - 37 -4.1.1 工程概况 .- 37 -4.1.2 工程地质 - 37 -4.2 计算模型 .- 38 -4.2.1 计算模型与网格划分 .- 38 -4.2.2 物理力学参数 - 39 -4.3 计算荷载及工况组合 - 40 -4.3.1 计算荷载 .- 40 -4.3.2 工况组合 - 41 -4.4 拱坝静力分析 - 42 -4.4.1 基本工况一 - 42 -4.4.2 基本工况二 - 45 -4.4.3 基本工况分析小结 - 48 -4.5 拱坝地震响应反应谱法分析 .- 49 -4.5.1 模态分析 - 49 -4.5.2 拱坝在地震作用下的动力响应 - 53 -4.5.3 静动叠加分析 - 56 -4.6 拱坝地震响应时程分析 - 63 -4.6.1 拱坝动位移与应力响应 - 63 -4.6.2 静动叠加分析 .- 65 -4.7 拱坝抗震安全评价 .- 68 -4.8 本章小结 .- 69 -第五章 结论与展望 - 70 -5.1 本文工作总结 .- 70 -VI5.2 展望 .- 70 -参考文献 - 72 -致谢 - 74 - 1 -第一章 绪论1.1 选题的背景及意义随着我国经济的高速发展，我国水电建设迎来高峰时期，高坝大库不断涌现，而小湾、溪洛渡、拉西瓦等规模巨大的 300 米级高拱坝多数在水能资源最丰富的西部地区进行建设。但是，我国是一个地震多发的国家，地震范围广而且强度高，而西部地区地震带更加活跃，新构造运动强烈，地质环境极不稳定，属于高地震烈度区。在这一地震高烈度区，抗震安全多为工程建设中的关键之一。高拱坝作为大型水利工程，其安全运行必须得到保证，高拱坝抗震性能评价关系到下游广大地区工农业生产和人民生命财产的安全，所以具有特别重要的意义。地震作用必将引起坝体运动和库水的运动，导致水体荷载的分布和大小发生改变，而且受坝体变形的影响。而拱坝由于体型单薄，单位质量所承受的动水压力很大，地震时坝体对动水压力的反应比较敏感，动水压力的大小直接关系到大坝的抗震安全。因此，水库对上游坝面作用的动水压力是危害坝体安全的一个不可忽视的重要因素 1。因此，在水工结构抗震计算中，常将动水压力作为必须考虑的重要荷载之一。尤其对高拱坝而言，搞清楚动水压力的规律，从而研究减小动水压力的工程措施，对提高拱坝的抗震安全性具有重要的意义。1933 年韦斯特伽特(H.M.Westergaard)发表其著名论文 地震时作用于坝面上的动水压力 2，文中假定坝体和地基为刚性，首先研究了刚性垂直坝面动水压力。从这时起，研究地震作用下高拱坝与库水的相互作用开始成为水工抗震设计中的重要课题。在随后的 70 多年期间，国内外学者基于不同的分析模型对坝库系统流固耦合的理论和计算方法开展了广泛研究，取得了一些成果。但是由于流固耦合问题相当复杂，目前无论在理论分析还是数值计算方面，关于动水压力对坝体地震反应影响的问题仍然没有得到合理解决，常用的方法是将库水作为附加质量引入。本文期望利用 ANSYS 软件对流固耦合问题进行数值模拟，以得到更接近于现实的结果。1.2 拱坝的类型、特点及其发展 3拱坝是一种建筑在峡谷中的拦水坝,是平面上凹向上游三向固定的空间高次超静定结构，它可以看成由竖立于地基上的一根根悬臂梁和横跨两岸的一层层水平拱组成，它能把上游坝面水压力、风浪压力等荷载相当大的部分通过水平拱的作用传递给两岸岩体，而将另一部分荷载通过悬臂梁的作用传至坝底基岩。它不像重力坝那样全靠自重维持稳定，而是利用筑坝强度承担以轴向压力为主的拱内- 2 -力，并由两岸拱端岩体来支承拱端推力。我国自新中国成立以来修建了许多拱坝。据不完全统计，至 1985 年底，全国已建坝高 15m 以上的各种拱坝总数达 800 多座，约占全世界已建拱坝总数的1/4，我国已经成为世界上修建拱坝最多的国家。进入 21 世纪，我国在建的高拱坝有小湾(292m) 、溪洛渡 (273m)、拉西瓦(254m)、锦屏一级(305m)及正在规划设计的构皮滩(225m) 等工程，显示了中国 21 世纪高拱坝的勘测、设计、施工和科研方面已经达到一个新的水平。而毋庸置疑的是，拱坝这一坝型将在我国的水利事业和国民经济发展中发挥至关重要而无可替代的作用。1.2.1 拱坝的类型拱坝的体形可以按照拱弧半径和拱中心角，分为单曲拱坝和双曲拱坝。(1) 单曲拱坝：又称为定外半径定中心角拱坝，只在水平截面上呈拱形，而悬臂梁断面不弯曲或者曲率很小。对 U 型或矩形断面的河谷，其宽度上下相差不大，各高程中心角比较接近，外半径可保持不变，仅需下游半径变化以适应坝厚变化的要求。(2) 双曲拱坝：又称穹形拱坝，在水平和铅直截面内都呈拱形。在 V 形河谷或其他上宽下窄的河谷中，宜将水平拱圈的半径从上到下逐渐减小，以使上下各层拱圈的中心角基本相等，并在铅直向设计成一定曲率，形成变半径等中心角双曲拱坝。对于大多数河谷，很难做到上下层拱圈中心角相等，为此广泛采用变半径变中心角的双曲拱坝，各层拱圈内、外的圆心连线均为光滑的曲线。1.2.2 拱坝的主要特点拱坝与其他坝型比较，具有如下一些特点：(1) 利用拱结构特点，充分发挥利用材料强度。拱作用越大，坝体材料的抗压强度越能充分发挥，坝体的厚度也可以变得越薄。所以同一坝址相同坝高的拱坝与重力坝相比，拱坝体积可节省 1/32/3。(2) 利用两岸岩体维持稳定。与重力坝靠自重在坝基上产生摩阻力来维持稳定的方式不同，拱坝将外荷载的大部分通过拱作用传至两岸岩体，主要依靠两岸坝肩岩体维持稳定，坝体自重对拱坝稳定性影响不占主导作用。(3) 荷载特点。拱坝不设永久性伸缩缝，其周边通常固接于基岩上，属于超静定结构。因而温度变化、地基变形等对坝体应力有显著影响。此外，坝体自重- 3 -和扬压力对拱坝应力的影响较小。坝体越薄，上述特点更加明显。(4) 超载能力强，安全度高。可视为拱梁系统的拱坝结构，当外荷载增大或某一部位因拉应力过大而发生局部开裂时，能调整拱梁系统的荷载分配，改变应力分布状态，不致使坝全部丧失承载能力。所以，拱坝允许坝面局部开裂，从而拥有很强的超载能力。(5) 抗震性能好。由于拱坝是整体性空间结构，厚度薄，富有弹性，因而其抗震能力较好。综上所述，拱坝是能够抵御洪水和强烈地震作用的高性能坝工结构，具有较大的经济性和较高的综合安全性。1.2.3 拱坝的发展早在两千年以前，人们就意识到拱结构具有比板其他结构更强的拦蓄流水的能力，并开始修建十余米高的圆筒形拱坝。它最早起源于古罗马时代的欧洲，在现今的法国境内圣里米省南部建造了第一座鲍姆拱坝 4。到了十三世纪末，伊朗修建了一座高约 60 米高的砌石拱坝。到二十世纪初，较高的拱坝在美国修建，如 1910 年完建的巴菲比尔拱坝，高 99 米。到第一次世界大战以前，国外坝高15 米以上的拱坝大约有 40 座。第一次世界大战以后，拱坝逐渐在欧洲流行。1920 年建成的瑞士蒙特萨文斯拱坝，采用了坝体竖向弯曲，水平拱为悬链曲线的变截面体形。1939 年意大利建成坝高 76.8 米，设置垫座及边缘缝的奥西尔埃塔双曲薄拱坝。19231935年美国提出考虑切向变位与水平扭转变位的新概念，即为以后试载法的基础，并于 1936 年建成第一座坝高超过 200 米的鲍尔德重力拱坝。第二次世界大战以后，拱坝又有了新的发展，修建数量也显著增加。这一时期修建的拱坝中，有不少为双曲薄拱坝。如意大利 1960 年所建瓦扬特拱坝，为一座高 265.5 米的双曲薄拱坝。五十年代中期，西欧各国的拱坝分析方法，多数已经采用纯拱法或拱冠梁法，并主要通过结构模型试验来验证。现今拱坝已经成为当今大坝设计中的三大优选坝型之一。根据国际大坝委员会世界大坝登记(1998) 一书统计，截止 1986 年底，全世界己建成的坝高大于 15m 的拱坝 1608 座，占世界大坝总数 36235 座的 4.4%5。- 4 -1.3 流固耦合理论研究现状1.3.1 流固耦合的定义及其分类在船舶、航空、水利、化工和核动力等工程领域中，都会遇到流体与结构的相互作用问题，简称流固耦合。流固耦合是研究变形固体在流场作用下的各种行为及固体变形对流场的影响的学科，涉及到流体力学、固体力学、振动力学等学科知识，是典型的多学科综合问题。其显著的特征是流固两相介质的相互作用：变形固体在流体荷载作用下会产生变形或运动，而这种变形或运动反过来又对流场产生影响，改变流场荷载的大小及分布。反应到其耦合方程中表现为：未知变量含有描述流体的变量和含有描述固体的变量，方程的定义域同时有流体域与固体域，一般具有以下两个特征：(1) 流体域和固体域均不可能单独求解；(2) 无法显式地消去描述流体运动的独立变量或描述固体运动的独立变量。目前对于流固耦合问题，一般从其耦合机理和其控制方程求解方法上分类。(1) 按流固耦合的耦合机理可分为两大类。第一大类问题的特征是：固体域和流体域部分或者全部重叠在一起，难以明显地分开，使描述物理现象的方程，特别是本构方程需要针对具体的物理现象来建立，其耦合效应通过描述问题的微分方程来体现。比较典型的是土木工程中的渗流问题；第二大类问题的特征好似：耦合作用仅仅发生在固体域和流体域的交界面上，在方程上的耦合是由两相耦合面的力学平衡及运动协调关系引入的。根据文献6，第二大类问题可分为三种问题：a. 固体与流体之间有大的相对运动，比较典型的是机翼颤振或悬桥振荡中发生的气固耦合，被人们习惯上称为气动弹性力学问题。b. 流体有限位移的短期问题，这类问题由引起位形变化的流体中的爆炸或冲击引起，其相互作用在瞬间完成，总位移有限，但是流体的压缩性十分重要。c. 流体有限位移的长期问题，这类问题比较常见，如充液容器的固液耦合振动、近海结构对波或地震的响应、大坝库水的耦合振动等。一般将仅为界面耦合的流固耦合分为以下四个方面的研究分支 7：a. 土木工程中的坝水耦合。主要研究在动水压力作用下的坝体结构的自振特性及动力响应。- 5 -b. 流体引起管道振动的研究。主要研究输液管道的稳定性及动力特性，包括直管和曲管在多种支撑条件下的振动特性。c. 船水耦合动力学。该分支的一个非常活跃的研究领域是水下结构的声振问题，其研究方式为将水下结构(如潜艇) 简化成为一水下球壳、圆柱壳、非圆柱壳或加肋壳体，然后应用解析法和半解析法研究其声振特性。d. 储液容器的振荡问题。研究储液容器中的罐-液耦合系统的动力特性及地震响应、壳-液相互作用以及在地震作用下储液罐的强度和稳定性问题等。(2) 按流固耦合问题控制方程求解方法同样分为两大类 89：弱耦合和强耦合。弱耦合法是对流体域和固体域分别建立和求解运动方程，通过某种差值函数进行方程间数据( 流体压力、结构位移等)的传递，从而实现结构和流体的耦合。弱耦合法是将流体和结构用各自的求解器在时域内积分，交错时间推进。其典型过程是：a.在流体力作用下对结构响应进行积分，推进至下一时刻；b.把结构的边界位移和运动传递给流体系统；c.更新流体动态网格；d.流体积分，计算新的流体压力和应力场；e. 把流体压力和应力转换成结构荷载，传递给结构。这种方法只需利用计算流体动力学和计算结构动力学的方法和程序，将其作少量修改就可实现程序的模块化，但是也存在缺陷，由于积分是交错进行的，流体和结构的时间推进积分总是存在时间滞后，耦合界面上的能量不能保持守恒，故无法满足动平衡。强耦合法是对流体和结构分别建立运动方程，采用数值方法将两个方程耦合在一起形成一个耦合方程，在求解耦合方程时，将同时求解耦合界面上的压力分布及结构压力响应，因此避免了弱耦合法中出现的时间滞后问题。显然强耦合法能够准确地描述流体的运动，更真实的反映流体和结构的相互作用。对流固耦合方程的求解，不论是强耦合还是弱耦合，都可以分为解析法、半解析法和数值法。对于简单的问题可以考虑用解析法和半解析法求解，但是对于边界条件和几何构造复杂的问题，则需要借助计算机采用数值方法求解。目前数值计算方法常用的有：有限元法、边界元法、有限差分法及这几种方法的混合方法等。这几种方法中，使用最多的还是有限元法。1.3.2 库水流固耦合理论研究发展过程1933 年韦斯特伽特(H.M.Westergaard)首先对坝面动水压力的分布进行了研究 19。由于 Westergaard 解中假定了刚性坝面并给定了地面运动规律，因而其解本质上并未计及固液两相耦合作用，只是求解了一个给定边界条件的流体动力问题，但是该公式也反映了动水压力的一些本质特征，对随后动水压力的研究和坝体设计产生了重要的影响。因此，目前许多国家的水工结构抗震计算以及大坝的- 6 -动力分析研究中，有关动水压力的计算均采用 Westergaard 附加质量公式。1952 年 Zanger 应用电模拟试验研究了在水平地震作用时，不可压缩库水作用下有倾度的折线型上游坝面的动水压力分布。1967 年 Chopra 的研究 10表明Westergaard 给出的坝面动水压力的解答仅适用于简谐地面水平运动频率小于库水基频的情况。考虑坝面倾斜时，很难求得可压缩库水作用于坝面的动水压力的解析解，因此在解析分析时，一般不考虑水体的压缩性。1978 年 Chwang 用二维势能流理论得到了常倾角倾斜坝面动水压力分布的精确积分公式。1982 年Clough 教授对 Westergaard 公式进行了推广 11。同年 Hanna 和 Humar 考虑库水可压缩性，建立了分析由于坝面水平运动引起的库水动水压力的边界元法，分析了不同坝面倾角、库底倾角和库水上游边界时坝面动水压力。六十年代以前，流固耦合的研究解答都是采用解析法进行的。七十年代开始，数值计算方法和计算机技术的发展，使流固耦合可以采用有限元、边界元及其混合法等数值方法进行计算。Chopra 等人 12采用子结构法研究了重力坝的地震作用，将坝、水、地基当作相互独立的子结构，进行了重力坝-库水-地基的二维有限元分析。其后仍用子结构法将二维发展到三维，进行了拱坝的三维流固耦合分析。傅作新在挡水坝面动水压力方面作了较多的研究 1315，采用流体边界元法和结构有限元法分别讨论了刚性坝及弹性坝的动水压力，并考虑了地基与库水的相互作用及库底泥沙吸收的影响。1.4 结构地震反应分析理论现状 1619结构地震响应取决于地震动和结构特性，尤其是结构的动力特性。结构地震响应分析的水平也是随着人们对这两方面认识的逐步深入而提高的。近几十年来，人们对地震动的谱成分和各类结构的动力特性有了深入的认识。因此，结构的分析也随之有了相应的发展。结构的动力理论发展按照不同的计算理论大体可以分为三种：拟静力法、反应谱法和时程分析法。(1) 拟静力法日本的大森吉房教授在 20 世纪初期提出用拟静力法分析地震动对结构的影响，假定地震时结构振动加速度与地面运动加速度完全一样，地震引起的惯性力就可以当作静力作用在结构上，然后按静力学方法计算结构的响应。因为这种方法把结构当成绝对刚性，没有考虑结构的动力特性，所以通常称它为静力法。事实上，只有结构的固有周期远远小于地面运动周期时，结构才可以被视为刚体而不计自身变形。用拟静力法计算由地震引起的惯性力的方法对于刚性结构- 7 -是适用的，但是对于像拱坝这样高而薄的柔性结构，就会产生较大的误差。因此，我国水工建筑物抗震设计规范规定：对于高度超过 150m 的大坝应进行时程分析。(2) 反应谱法现行结构工程抗震设计规范普遍采用的是反应谱法。它是根据已取得的多个实测地震波，将其分别代入单自由度动力反应方程，计算出各自的加速度、速度、位移等最大弹性地震反应，这样即可得到结构最大地震反应与其自振周期的关系曲线，即反应谱。将结构的最大地震作用通过反应谱，转化为作用于结构的等效荷载，然后根据这一荷载用静力分析方法求得结构的地震反应和变形，从而使得十分复杂的多自由度体系地震反应求解变得简单。但是，地震作用是一个时间过程，反应谱法不能反映结构在地震过程中的经历，同时目前应用的加速度反应谱属于弹性分析范畴，当结构在强烈地震作用下进入塑性阶段时，反应谱分析将不能反映结构真实的地震响应。(3) 时程分析法随着电子计算机的大量普及而兴起的结构动力响应数值分析以及强震观测记录和震害经验的积累，人们逐步认识到像反应谱那样的等效静力法，其计算结果并不足以保证结构的抗震安全，考虑地震动全过程进行真正的结构响应动力分析，是非常必要的。时程分析法又称逐步积分法，其实质是用数值积分法来求解运动微分方程。该方法是由初始状态开始积分直至地震终止，求出结构在地震作用下从静止到振动，直至振动终止整个过程的地震反应 (位移、速度和加速度 )。 目前常用的数值积分方法主要有三种：线性加速度法、 Wilson- 法和Newmark- 法。 但其计算工作量大，计算步骤繁琐，只能依靠计算机编程计算，随着计算机和计算技术的飞速发展时程分析法已逐步在工程开始得到应用。1.5 本文的主要研究内容1、阅读国内外相关文献，综述了拱坝的结构特点和发展现状， 动力分析方法以及流固耦合理论研究现状。2、介绍了 ANSYS 软件及其 APDL 语言的二次开发技术，水工结构动力分析的理论与方法(包括结构动力平衡方程、结构自振特性计算方法)以及基于 ANSYS- 8 -软件的反应谱法和时程分析法的实现。3、详细论述了流固耦合分析的基本理论及求解方法。重点研究了在 ANSYS软件中考虑动水压力的方法，并对其实现步骤作出说明。以一二维重力坝水体地基模型为例，分析重力坝在空库，死水位，满库三种情况下的自振特性，并且分析坝库系统承受水平位移激励作用时的水体动水压力分布规律。4、结合天花板水电站拱坝工程，建立了拱坝-地基的三维有限元模型，分析了正常蓄水位和死水位下的两种基本静力工况的位移及应力的分布情况；研究了拱坝的动力特性，并采用反应谱分析法和时程分析法分析了正常蓄水期拱坝在地震作用下的动位移以及动应力分布情况，最后对天花板拱坝抗震性能作出安全评价。第二章 动力分析理论2.1 ANSYS 软件简介2.1.1 ANSYS 软件概述ANSYS2022软件是第一个通过 ISO9001 质量认证的大型分析设计类软件，是美国机械工程师协会(ASME)、美国核安全局(NQA) 及近二十种专业技术协会认证的标准分析软件，是融合结构、热、流体、电场、磁场、声场分析于一体的大型通用有限元分析软件，可广泛应用于土木工程、地质矿产、水利、铁道、汽车交通、国防军工、航天航空、船舶、机械制造、核工业、石油化工、轻工、电子、日用家电和生物医学等一般工业及科学研究。通过它可以对多数结构在各种外荷载条件下的受力、变形、稳定性及各种动力特性做出全面分析，在力学计算、- 9 -组合分析等方面提出全面的解决方案，从而为土木工程师提供功能强大且方便易用的分析手段。ANSYS 分析过程包含 3 个主要的步骤：前处理、分析计算、后处理 23。(1) 前处理。指创建实体模型以及有限元模型。它包括创建实体模型，定义单元属性，划分有限元网格，修正模型等几项内容。ANSYS 在实体建模时可以采用自底向上或自顶向下方式。自底向上方式首先定义关键点，然后利用这些关键点定义较高级的线、面、体等实体图元；自顶向上方式则是一开始就从较高级的实体图元构造模型，程序会自动生成所有从属于该体素的较低级图元。(2) 分析计算。在前处理阶段完成建模后，可以定义分析类型、分析选项、载荷数据和载荷步选项，然后进行有限元求解计算。ANSYS 软件提供的分析类型包括结构分析、热分析、电磁分析、流体分析、耦合场分析等。在 ANSYS 中可以定义的载荷包括边界条件(约束、支承)和其他外部或内部作用载荷，可分为六大类： DOF 约束；表面分布载荷；体载荷；惯性载荷；耦合场载荷；力。(3) 后处理。指检查分析的结果，这可能是分析中最重要的一环。后处理可使用两个后处理器：通用后处理器 POST1 和时间历程后处理器POST26。POST1 允许检查整个模型在某一载荷步和子步(或对某一特定时间点或频率) 的结果。POST26 可以检查模型的指定点的特定结果相对于与时间、频率或其他结果项的变化。2.1.2 APDL 参数化语言概述 ANSYS 软件提供了一种以命令流方式进行分析的功能，即 APDL 语言。它的全称是 ANSYS Parametric Design Language，是一种类似 FORTRAN 的解释语言，它提供一般程序语言的功能，如参数、宏、标量、向量及矩阵运算、分支、循环、重复以及访问 ANSYS 有限元数据库等，另外还提供简单界面定制功能，实现参数交互输入、消息机制、界面驱动和运行应用程序等。用户可以利用APDL 语言将 ANSYS 命令组织起来，编写出参数化的用户程序，从而实现有限元分析的全过程，即建立参数化的 CAD 模型、参数化的网格划分与控制、参数化的材料定义、参数化的载荷和边界条件定义、参数化的分析控制和求解以及参数化的后处理。因此 APDL 语言大大方便了建模、计算及数据处理，同时使有些凭手工无法想象的工作成为可能 24。APDL 参数化语言具有下列优点 25：- 10 -(1) 可以减少大量的重复工作，特别适用于经少许修改后需要多次重复计算的场合，可为设计人员节省大量的时间；(2) 便于保存和携带，一个 APDL 的 ASCII 文件一般只有几十 KB，最多也只有几百千字节，其文件的容量仅为 GUI 数据文件的千分之一；(3) 不受 ANSYS 软件的系统操作平台的限制，即用户使用 APDL 文件既可以在 Windows 平台进行交流运行，也可以在 UNIX 或其他的操作平台上运行。(4) 不受 ANSYS 软件版本限制。一般情况下，ANSYS 软件以 GUI 方式生成的数据文件只能向上兼容一个版本，而 APDL 文件则不存在这个限制，仅有个别命令会受影响；(5) 在进行优化设计和自适应网格分析时，则必须使用 APDL 文件系统；(6) 利用 APDL 方式，用户很容易建立参数化得零件库，以利于其快速生成有限元分析模型；(7) 利用 APDL 方式，可以编写一些常用命令流的集合即宏命令，或者是制作快捷键，并将其放在工具栏上；(8) 可以利用 APDL 对 ANSYS 进行二次开发。2.2 动力分析有限元基本理论2.2.1 结构动力平衡方程 26有限单元法应用于结构的动力学问题时，在基本概念上与静力有限元法完全相似。但是，这时除了考虑结构受到体力、面力等外荷载作用之外，还要考虑由于结构运动而引起的惯性力和阻尼力。采用有限元法分析结构的地震响应时，根据最小势能原理可以导出整个结构的动力平衡方程为(2-1)gKCM式中， 、 、 分别为结构的结点位移、结点速度和结点加速度列阵；、 分别为系统的整体劲度矩阵和整体阻尼矩阵； 是结构的集中质量KC 矩阵，可由下列单元劲度矩阵 、质量矩阵 、阻尼矩阵 集合而成。ekedmeceTvBDv- 11 -evTddvNmevc式中， 为应变矩阵； 为弹性矩阵； 为形函数矩阵； 为材料的密度；BD为材料的阻尼系数。 为地震时的地面运动加速度列阵。gTgggg wvuvu式中， 、 、 分别表示地面加速度沿着 3 个坐标轴方向的分量。根据基本guvw方程(2-1)，可以解决结构自振特性和结构在地震时的动力反应等问题。2.2.2 结构自振特性的计算 27结构的自振特性(频率和振型)是结构的动力计算中的主要内容之一。在动力方程(2-1) 中，令 为零，便得到自由振动方程。在实际工程问题中，阻尼对g结构的自振特性影响极小，所以求频率和振型是可以不考虑阻尼的影响，所以在动力方程(2-1)中， 为零，可得到结构的无阻尼自由振动方程：C(2-2) 0KM设在自由振动时，各质点作简谐振动，各结点的位移可以表示成为：(2-3)0cost式中， 为结点振幅列阵，即振型； 为与该振型对应的频率。0将式(2-3)代入方程(2-2)，可得广义特征方程：(2-4)20()KM结构在自由振动时，各结点的振幅不全为零，所以式(2-4)的系数行列式必须为零，由此可得求解结构自振频率的方程：(2-5)20方程(2-5)称为特征方程。因为结构的劲度矩阵 和质量矩阵 都是 阶KMn的方阵，其中 等于自由度的数目，所以式(2-5) 是关于 的 次代数方程，由该n 2 n方程解出的 个 的值可按升序排列为：2- 12 -(2-6)2210n第 个特征值 的算术平方根 称为结构的第 阶固有频率。i2i i i对于每个自振频率，由式(2-5)可确定一组各结点的振幅值。它们互相之间应保持固定的比值，但绝对值可任意变01020,Tiiini化，它们构成一个向量，称为特征向量，在工程上通常称为结构的振型。在大型复杂结构的振动分析中，系统的自由度可能多达几百甚至几千，但是需要用到的自振频率和振型往往是最低的 10 阶或 20 阶。因此，选择合适的分析方法将有助于提高动力的精度以及保证结构分析的准确性。本文采用子空间迭代法来分析有限元模型的动力特性。子空间迭代法是解决大型结构问题的有效方法之一，它可以按照任意要求的精度求得体系的自振周期与振型，因而在实践中对子空间迭代法甚为推荐。运用子空间迭代法可以很顺利地求解结构体系的最低的几阶频率及振型。子空间迭代法是对一组初始向量反复地应用幂法和李兹法。应用幂法师使其低阶振型成分增加，而高阶振型成分缩减。经过幂法迭代得到的一般不能满足正交性条件，如果不经过正交化处理，连续多次迭代后则都将收敛得到第一振型。应用李兹法，使 转化为 ，不仅满足了正交性条件，而且也得到了前 s个频率和振型的最佳近似值。经多次迭代后，将收敛于前 s 个自振频率和振型。子空间迭代法通常用于结构频率范围难以估计，同时无法选择主自由度的情况。对大型对称特征值求解，对于初始迭代向量的选择要求不高，计算特征对的个数不受限制。所以，综合对 ANSYS 中各模态提取方法所采用的求解器、计算精度、计算速度和适用范围进行比较后，选择子空间迭代法进行动力求解。2.3 反应谱分析理论与方法 28振型分解反应谱法是在振型叠加法的基础上推导出的一种近似方法。这个方法需要事先求出结构的若干个振型和频率，但是，可以直接利用标准的设计反应谱，求各振型的最大动力反应最大绝对加速度、最大相对速度和最大相对位移。振型分解反应谱法无须对动力方程(2-1)做数值积分，计算存贮量和运算时间最省。振型分解反应谱法的优点是可以采用由统计方法得到的标准反应谱，避免了选择地震加速度记录的困难，它的缺点是不能用于非线性振动情况。根据反应谱理论，结构各阶振型的最大地震反应与具有相同周期的单自由度- 13 -体系的最大反应成正比，即 maxmax0ii iAaaii iVmaxax0ii iU式中， 、 、 是体系第 阶振型的最大绝对加速度、最大相对maxiaxiai速度和最大相对位移； 、 、 分别是周期相同的单自由度弹性体系mAxVmax对同一地震波的响应； 、 分别是第 阶振型的振型向量和对应的振型参与i0ii系数 2728。由图 2-1 可见，反映结构振动加速度的动力系数 的大小与结构的自振周期有关，当结构较柔( 自振周期 较大)时，相应的 值将减小。TT同时，从反应谱也可看出，当结构的自振周期 和场地的特征周期 (即相TgT当于地震波的主震周期)相接近时，将发生共振，这时，结构的 值将明显增大。在反应谱中，设计加速度反应谱最大值 对阻尼比为 0.05 的水工钢筋混max凝土结构，如拱坝可取 。max2.50当用反应谱法求结构的地震反应时，应先求出结构的若干个低阶振型和周期，根据求得的周期并利用设计反应谱图求得各振型对应的设计反应谱值 ，便可i进而求出结构各振型的最大加速度向量、最大荷载向量、最大位移向量、和最大图 2-1 水工建筑物抗震设计反应谱- 14 -应力向量(2-7)max0iiikg(2-8)aii iFM(2-9)2max0/i ii iUkg(2-10)ax0/iiii式中， 为地震系数，在 7 度地震时， ； 是重力加速度； 为第 阶k 1.kgi0振型向量； 为与第 阶振型向量对应的振型应力向量； 为结构第 阶振型的i0i i参与系数，其分量形式为， ，2,()kixmxyz2,()kiymyxz2,()kizmzxy式中， 为 结点的质量； 、 、 为 结点的振型分量。k kkz地震作用效应按平方和求平方根法进行组合，即(2-11)21mEjjS式中， 为组合后的地震作用总效应； 为第 阶振型的地震作用效应； 为ESj m计算采用的振型数目。当两个振型的频率差的绝对值与其中一个较小的频率之比小于 0.1 时，地震作用效应采用完全二次型方根法进行组合。即(2-12)1mEijjijSS(2-13)322 28()(1)44()ijijij ij ij 式中， 为第 阶振型的地震作用效应； 为第 阶和第 阶的振型相关系数；iSij、 分别为第 阶、第 阶振型的阻尼比； 为圆频率比， ； 、ijij iji- 15 -分别为第 阶、第 阶振型的频率。jij当各振型的阻尼比相同时，式(2-13)可简化为(2-14)23/28(1)()4()ij 基于 ANSYS 软件的 APDL 语言，开发了用于反应谱分析的命令流如下：(1) 模态分析求解自振频率命令流/soluantype,2modopt,lanb,10 !用子空间迭代法求前 10 阶振型mxpand,10,yessolvefinish(2) 反应谱分析命令流/solu*dim,freq,18,1,1*dim,sv,18,1,1Tg=0.50 !定义特征周期betamax=2.3 !设计加速度反应谱最大值*do,i,1,16,1freq(i)=1/(3-(3-Tg)*(i-1)/15)sv(i)=0.143*9.81*betamax*(Tg/(3-(3-Tg)*(i-1)/15)*0.9*enddofreq(17)=10 freq(18)=1000 sv(17)=0.143*9.81*betamax sv(18)=0.143*9.81/soluantype,spectr !求振型和频率 - 16 -spopt,sprs,10,yes !单点谱分析 svtype,2 !加速度谱 sed,1,1,0 !定义地震激励方向freq,freq(1),freq(2),freq(3),freq(4),freq(5),freq(6),freq(7),freq(8),freq(9), sv,sv(1),sv(2),sv(3),sv(4),sv(5),sv(6),sv(7),sv(8),sv(9), freq,freq(10),freq(11),freq(12),freq(13),freq(14),freq(15),freq(16),freq(17),freq(18)sv,sv(10),sv(11),sv(12),sv(13),sv(14),sv(15),sv(16),sv(17),sv(18)allselsolvefinish/soluantype,modal !模态扩展expass,onmxpand,10,yes,0,0001/soluantype,spectr srss,0.15,disp !合并模态solvefinish2.4 时程分析法理论与方法 28结构地震反应分析的反应谱方法是将结构所受的最大地震作用通过反应谱，转化成作用于结构的等效侧向荷载，然后根据这一荷载用静力分析方法求得结构的地震内力和变形。因其计算简便，所以广泛为各国的规范所采纳。但地震作用是一个时间过程，反应谱法不能反映结构在地震动过程中的经历，同时，目前应用的加速度反应谱属于弹性分析范畴，当结构在强烈地震下进入塑性阶段时，用此法进行计算将不能得到真正的结构地震反应。对于长周期结构，地震动态作用下的地面运动速度和位移可能对结构的破坏具有更大影响，但是振型分解反应谱法对此无法作出估计。所谓时程分析法，是根据选定的地震波和结构恢复力特性曲线，采用逐步积- 17 -分的方法对动力方程进行直接积分，从而求得结构在地震过程中每一瞬时的位移、速度和加速度反应，以便观察结构在强震作用下从弹性到非弹性阶段的内力变化以及构件开裂、损坏直至结构倒塌的破坏全过程。这类方法是指不通过坐标变换，直接求解数值积分动力平衡方程。其实质是基于以下两种思想：第一，将本来在任何连续时刻都应满足动力平衡方程的位移 ，代之以仅在有限个离散时刻 ，t0t， ，可以满足这一方程的位移 ，从而获得有限个时刻上的近似动力平1t2 t衡方程;第二，在时间间隔 内，以假设的位移、速度和加速度的变化iit1规律代替实际未知的情况，所以真实解与近似解之间总有某种程度差异，误差决定于积分每一步所产生的截断误差和舍入误差以及这些误差在以后各步计算中的传播情况。其中前者决定了解的收敛性，后者则与算法本身的数值稳定性有关。实用中一般取等时间间隔，从初始时刻 到某一指定时刻 ，逐步积分0t Tt0求得动力平衡方程的解。把区间0，T进行 n 等分后有 ，相应的第 i 个n/离散时刻为 (i=1， 2，3) 。it此法费时较多，且确定计算参数尚有许多困难，因此目前仅在重要的、特殊的、复杂的以及高层建筑结构的抗震设计中应用，在拱坝结构的抗震分析中应用较少。此外，时程法亦用于结构在地震作用下破坏机理和改进抗震设计方法的研究。目前有限元法进行大型结构动力计算的解法有线性加速度法、Wilson法和 Newmark 法。ANSYS 程序中结构动力分析使用的是 Newmark 时间积分法，在离散的时间点上求解这些方程，在 时刻，系统动力方程式为tn(2-15)1111nnnnMCKf式中， 、 、 分别表示加速度向量、速度向量及位移向量；1n表示作用力矢量。f在 Newmark 方法中，应满足下列方程(2-16)1111nnnnMCKf(2-17)2t- 18 -(2-18)111nnnt其中，(2-19)21 1/2nnnntt (2-20)nnt式中， 、 、 是 时刻的位移向量、速度向量和加速度向量； 、nnt 为 Newmark 参数，控制方法的精度和稳定性，通常取 ， ； 25.0.、 是预估值； 、 是校正值。1n1n1n1n当给定初始位移 和初始速度 后，可以从下式求出初始加速度00 0(2-21)0000MfCK然后将动力问题变成“静力等效问题”求解(2-17)(2-22)式，并结合Newton-Raphson 方法求解，具体步骤如下(1) 预估值( 为迭代计算变量)i211 1/2innnnntt1t t211/0iinnnt(2) 计算残余力 1111i i i innnnfMCK (3) 变步时，用下式计算等效劲度矩阵 *2/Ktt(4) 求解位移增量 *ii- 19 -(5) 修正位移、速度、加速度 11iiinn21 /ii t11ii innnt(6) 如果 和 或两者之一不满足收敛条件，则令 转到第二ii 1i步，否则继续下去。(7) 在下一时间步内，令11inn1i1inn同时令 ，形成新的劲度矩阵 并开始下一时间步运算。K在 ANSYS 软件中利用 APDL 语言编制的时程分析命令流如下：Finish/config,nres,3000 !结果文件中允许的荷载步数或子步数的最大值*dim,acex,array,1001,1*creat,a !创建宏文件*vread,acex,acel,txt,jik,1,1001 !读取 acel 文件数据并产生数组参数矢量或矩阵(e10.3)*dim,acey,array,1001,1*vread,acey,acel,txt,jik,1,1001 !读取 acel 文件数据并产生数组参数矢量或矩阵(e10.3)*end- 20 -/input,a !为随后的命令打开文件/soluantype,trans*do,t,1,1001,1time,t*0.01alphad,betad,acel,acex(t,1),acey(t,1) !定义 X 向、Y 向地震加速度allsel,allsolve*enddooutres,all,allfinish2.5 地震波输入与选取 2930对重要的建筑物和其它特殊结构采用多节点多自由度的结构动力有限元分析时，往往是把地震强迫振动的激振地震加速度时程直接输入，对结构进行时程反应分析。这种动态时程的分析方法可以考虑各种不同因素，使结构抗震计算分析的结果更加符合实际震害现象，也使结构工程师更清楚结构地震动力破坏的机理和正确提高结构抗震能力的途径。但结构动力时程分析的可靠性不仅取决于结构计算模型的合理和计算方法的精确，而且与所使用的地震输入即地震波的选用有直接关系。因此，对于时程分析法，合理选择适宜的地震波是极为重要的。地震作为一种自然现象，无论其规模、产生概率以及地震的波形特征都具有其随机性。目前国内外已经积累了一定数量的强震记录可供时程分析选用。这些记录由于能真实反映地震动的特点，已被广泛应用于重要工程的抗震设计。(1) 地震波选取输入地震波选择的原则是使输入地震波的特性和建筑场地的条件相符合。主- 21 -要参数有：地震烈度、地震强度参数、场地的土壤类别、卓越周期和反应谱等。选择地震波时应选其主要周期与建筑场地卓越周期接近的地震波。此外还要满足地震活动三要素的要求：即频谱特性(可用地震影响系数曲线表征，依据所处的建筑场地类别和设计地震分组确定)、幅值(一般按规范所列地震加速度最大值采用)和地震加速度时程曲线持续时间(一般为结构基本周期的510 倍) ，这三者均要符合