灯泡贯流式水电站厂房三维静动力分析(四)

摘要：国内水利水电工程建设目前正处于前所未有的蓬勃发展时期，许多低水头径流式水电站建设逐步在我国的江河上兴建，其中灯泡贯流式水电站由于流道平坦，机组过流量大、单位转速高、效率高、尺寸小、重量轻、能量及经济指标好等优.点成为目前比较普遍的一种开发型式。然而，由于灯泡贯流式水电站厂房独特的布置型式，致使应力分布有不同于常规水电站厂房的特点，特别是在高地震烈度区修建的灯泡贯流式水电站。因此，本项目的研究分析具有十分重要的现实意义。

关键词：灯泡贯流式水电站 静 动力计算分析

2.4 结构动力问题的有限元法

动力学问题在国民经济和科学技术的发展中有着广泛的应用领域。最经常遇到的是结构动力学问题，它有两类研究对象。一类是在运动状态下工作的结构，另一类是承受动力荷载作用的工程结构。结构受载荷处于平衡状态时，是静止不动的；结构有变形，而位移是不随时间而改变的，载荷和内部应力也不随时间而变化，这是静力问题。结构受载荷没达到平衡状态，或由于结构的弹性和惯性而围绕平衡位置振动时，其位移、应力等都是时间的函数，各点有位移还有速度和加速度，这是一种动力问题。有限元方法可以用来分析连续结构的动力问题[70]。

2.4.1结构动力学方程[71]

对于动态结构而言，所受的外力（包括体力、面力、集中力、惯性力和阻尼力）和产生的位移都是时间的函数。应用达伦贝尔原理，把结构的惯性力加入平衡方程中，就可以将弹性的结构的动力问题转化为静力平衡问题来处理。

用有限元法求解弹性结构的动力问题，也是把结构离散成有限个单元的集合体，并取出任意单元，此时单元上任意点的位移都是时间的函数，以表示单元上的节点位移向量，再利用单元的位移插值公式，写出单元的上任意点的位移函数：

(2-11)

其中，为形函数，是位移的插值函数，与时间无关。

则速度和加速度函数为：

(2-12)

(2-13)

其中，、为单元节点的速度和加速度列阵。

将单元内惯性力与阻力作为体积分布载荷分配到单元各节点上，分别记为、，有

将式（2-11）、（2-13）代入上式，有

令 (2-14)

称为单元质量矩阵；

令 (2-15)

称为单元阻尼矩阵。

按达伦贝尔原理，将惯性力、阻力作为载荷，单元叠加得到弹性结构的动力平衡方程：

(2-16)

令 、

则方程（2-16）改写为：

(2-17)

弹性结构的振动本身是连续体的振动，位移是连续的，具有无限多个自由度。经有限元离散化后，单元内的位移按假定的位移形式来变动，可用节点位移插值表示。这样，连续系统的运动就离散化为有限个自由度系统的运动。尽管如此，结构动力有限元计算量比静力的大得多。为保证计算的方便、快捷并满足一定计算精度的要求，可以采用合理的计算方法和计算程序；宜可从力学角度简化动力方程，如通过集中质量矩阵、静力缩聚、主副自由度、模态综合等方法已达到降阶和简化方程的目的。

一般的连续结构都可以用有限元方法化为有限自由度系统问题，并列出相应的动力方程。在给定的节点载荷作用下，求解动力方程，可归纳为两种方法。一是通过求解大型的矩阵特征值问题确定结构的动力特性，经模态矩阵变换，化为互不耦合的N个单自由度问题，逐个求解并迭加，称振型迭加法。这需要算出系统的各阶振型，而且也仅适用于线性系统和简单的阻尼情况。二是用数值计算直接积分多自由度系统的微分方程，写成矩阵形式用计算机逐步求解，这可用于一般阻尼的情况，并且可按增量法，用逐段线性化的方法求解非线性系统问题。

（1）振型迭加法

对于多个自由度系统，结构的动力反应可以用各个振型动力反应的线性组合来表示，即

(2-18)

式中，为位移向量；为广义的坐标向量；矩阵为振型矩阵，振型矩阵中第列向量即为系统的第个振型向量。将（2-18）式代入系统的动力方程式（2-17），并左乘振型向量后，可得

(2-19)

利用振型关于质量和刚度矩阵的正交性，并假定阻尼矩阵也满足正交性条件，可以得到：

(2-20)

式中、分别为振型质量和振型刚度，为振型阻尼，根据假定也满足正交性条件，即，当采用瑞利阻尼时，很明显，，这个条件是自然满足的；称为振型节点荷载。

逐个求解（2-20）式，即可得到个广义坐标，代入式（2-11），即将得到了结构系统的反应。用振型分解法求得的节点位移是时间的函数，由它插值的单元内部位移、应力、应变的计算与静力计算一样，不同的是这些量都是时间的函数。

用振型分解法求解结构系统的动力反应时有两个明显的优点：一是个相互耦连的方程利用振型正交性解耦后相互独立，变成了个自由度方程，使计算过程大大简化。二是只需按要求求解少数几个振型的方程，就可以得到满意的解答，因为在大多数情况下，结构的动力反应主要是前面几个低阶振型起控制作用。

（2）直接积分法

在结构动力计算中，常用的直接积分法有中心差分法、线性加速度法、法和法等等。

1）、中心差分法

中心差分法的基本思路是将动力方程式中的速度向量用位移的某种组合来表示，将微分方程组的求解问题转化为代数方程组的求解问题，并在时间历程内求出每个微小时段的递推公式，进而逐步求的整个时程的反应。

对于动力方程（2-17）各阶微分可以用中心差分表示为

(2-21)

(2-22)

式中为均匀的时间步长，、和分别为时刻及其前、后时刻的节点位移向量。将式b、c代入a式后可得到一个递推公式如下：

(2-23)

上式即为中心差分法的计算公式，在求得结构的和后，就可以根据t时刻及t-Δt时刻的结点位移，按（2-23）式推算出t+Δt时刻的结点位移；并可逐步推出t+2Δt，t+3Δt，…，tend各时刻的结点位移。 式（2-23）对于t=0的时刻并不适用，因为一般运动的初始条件给出的是初始位移和初始速度，而难以给出前一个Δt时刻的位移，无法直接按式（2-23）进行第一步的计算，因此，这时就要利用其他条件建立中心差分的计算公式，

= (2-24)

(2-25)

再利用t=0时刻的动力方程：

(2-26)

由（2-24）、（2-25）、（2-26）三式，可以求得、和。求解的方程式如下：

(2-27)

这个方程式中的、和都是已知的，因此可以解出。而后就可以按式（2-24）解出和，…。这是一种将时间段划分为若干个相同的时段后的逐步求解方法，求解出的量均是每个时刻结点的位移，因此，很适合于像有限元方法这样以结点位移来计算单元内部位移、应力和应变的各种数值求解问题。

2） 线性加速度法

这个方法的基本思路是把整个振动时程分成很多个时间间隔，并假定在范围内加速度按直线规律变化，在此基础上计算出时刻内的增量位移、增量速度和增量加速度，一步一步地求得整个时程的反应。

将动力方程式写成增量形式的方程：

(2-28)

用时刻的和表示和，代入（2-28）并整理后得

在求出后，及可按下式求出：

(2-30)

这样，t时刻的位移、速度和加速度可按下式求出：

(2-31)

重复上述步骤，可根据体系的初始条件，一步一步地求得各时刻（1，2，…，n）时系统的动力位移、速度和加速度反应。

数值计算方法的一个基本要求是算法的收敛性好，上一节介绍的线性加速度法当体系自振周期较短而计算步长较大时，有可能出现计算过程发散的情况，即计算的反应数值越来越大，直至溢出（overflow），对于多自由度系统，其最小的自振周期可能很小，此时，计算步长Δt必须取得很小才能保证计算不发散。对于结构抗震分析来说，Δt需要选得比地面运动中高频分量的周期以及结构的自振周期小很多（例如10倍以上），才能保证必要的精确度。因此，线性加速度法是一种条件收敛的算法。

Wilson-θ法是在线性加速度法基础上改进得到的一种无条件收敛的数值方法，它的基本假定仍然是加速度按线性变化但其范围延伸到时间步长为θΔt的区段，只要参数θ取得合适（θ≥1.37），就可以取得收敛的计算结果。当然，Δt取得较大时，计算误差也将较大。

在时刻t+θΔt，多自由度系统的运动方程式为

[M]{(t+Δt)}+[C]{ (t+Δt)}+[K]{ (t+Δt)}={P(t+Δt)}

(2-32)

根据Wilson-θ法的基本假定，加速度反应在[t，t+θΔt]上线性变化，即在此区段上运用线性加速度法得到的公式，并将时间步长改为θΔt，即可求得时刻t+θΔt时的加速度反应为

{(t+Δt)}=

(2-33)

在[t，t+θΔt]时段内采用内插法，可以求得t+Δt时刻的加速度为

{(t+Δt)}={(t)}+

={(t+Δt)}+

= (2-34)

根据线性加速度法的基本关系式，利用{(t+Δt)}可得

(2-35)

{} (2-36)

式（2-35）、（2-36）即为用Wilson-θ法计算结构动力反应的公式。

4）Newmark-β法

Newmark-β法的基本假定是：

{δ(t+Δt)}={δ(t)}+ (2-37)

其中，γ和β是按积分的精度和稳定性要求而调整的参数。研究表明，当γ>=0.5，β>=0.25(0.5+γ)2时，Newmark-β法是无条件稳定的。

由式（2-37），可利用{：

{(t+t)}=

(2-38)

{}

(2-39)

考虑到t+Δt时刻的动力方程，有：

[M]{（t+Δt）}＋[C] {(t+Δt)}+[K]{}={P(t+t)} (2-40)

将式（2-39）代入上式，可得：

（2-41）

式中

求解方程（2-41），可得{δ（t+Δt）}，然后由式（2-39）可解出{}和{}。以此类推，可求出各时刻的位移、速度和加速度。

2.4.3结构体系自振周期、振型计算

结构的自由振动问题可以归纳为求解广义特征值问题[66，76]，广义特征值为1/ω2，广义特征向量为结构的固有振型。

忽略结构的阻尼影响，结构的自由振动方程为：

(2-42)

假设位移向量，由上式得：

(2-43)

式中：[K]、[M]分别为结构的整体刚度矩阵、质量矩阵；

、分别为结构各质点的位移、加速度；

ω为结构自由振动的圆频率。

一般地振型向量≠0，由齐次线性方程组解的理论得：

（2-44）

由式（2-44）得到n个不同的圆频率ω1、ω2、ω3、…、ωn，将圆频率代入方程（2-43）可得到固有振型{A}1、{A}2、…、{A}n。

由于非对称框架结构隔震系统的质量矩阵、刚度矩阵为非对角矩阵，程序中求解自振频率及振型采用广义雅可比法。