建筑物结构动力分析用自振频率和振型的计算方法

摘要

本发明涉及建筑物结构动力分析用自振频率和振型的计算方法，该方法包括有如下步骤：建立多自由度的建筑结构体系自振频率方程，并将其转化为特征方程问题；从自振频率方程中求解建筑结构体系的自振频率，先采用带“原点位移”的Schmidt正交化QR分解算法来求解自振频率方程的全部特征值，然后计算自振频率；根据自振频率方程解出的特征值来计算建筑结构体系的振型，通过1次求解非齐次线性方程组的方法直接求解特征向量，特征向量标准化后就是所求的振型。本发明的计算方法计算速度快，计算过程稳定，不会出现特征值遗漏的情况，该算法不仅可以在桥梁等工程中应用，还可以求解有阻尼情况下的振动频率及对应振型。

说明书

技术领域

本发明涉及土木建筑工程领域，特别涉及到多自由度体系下建筑物及高层建筑进行结构动力分析时自振频率和振型的计算方法。

背景技术

自振频率和振型是建筑物体系本身的固有性质，它们是设计和研究结构体系的抗风、抗震的基础。在多自由度下，自振频率和振型不止一个。自振频率与本身的刚度系数及其质量分布有关，而与外部荷载无关。而振型可以理解为结构体系振动的特定形式，通过标准化方法可以唯一确定结构体系的相对振幅。振动频率和振型是建筑物进行结构动力分析的基本设计参数，它们是其他动力学分析的起点，例如瞬态动力学分析、谱响应分析和谱分析。瞬态动力学分析亦称时间历程分析，是用于确定承受任意随时间变化载荷结构的动力学响应的一种方法。谱响应分析和谱分析在建筑尤其是高层建筑设计中，常用来计算等效地震荷载。

按多自由度体系来计算建筑物的自由振动问题，如多层厂房的侧向振动、不等高排架、高层建筑的振动等，常采用在建筑（包括高层建筑）中应用十分广泛的“刚度法”或“柔度法”求解，前者是通过建立力的平衡方程求解，后者通过建立位移协调方程求解。

现有技术中对多自由度体系自振频率和振型的求解，一般采用近似的或迭代的方法，虽然总是能产生最低的几个特征值及自振频率，但这一点并没有得到数学上的严格证明，有时还会出现不收敛问题。所以，还需要通过序列检验法来检验是否有遗漏的特征值及自振频率。同时，有了特征值及自振频率之后，仍是通过迭代法来求解，一般需要多次求解线性方程组，造成计算过程非常繁琐，计算结果不稳定。

发明内容

本发明的目的在于克服上述现有技术中存在的不足，提供一种新的建筑物特别是高层建筑进行结构动力分析时自振频率和振型向量的计算方法。本发明的计算方法要能够做到计算速度快，计算结果稳定，不会出现特征值遗漏。

为了达到上述发明目的，本发明专利提供的技术方案如下：

一种建筑物结构动力分析用自振频率和振型的计算方法，其特征在于，该计算方法包括有如下步骤：

第一步，对于选定的建筑物，多自由度的建筑结构体系在无阻尼自主振动时，以柔度法建立建筑结构体系的自振频率方程，并将其转化为特征方程，假定建筑结构体系只有平移自由度，其刚度矩阵[K]表示为[K]=k[K₀]，k为常系数，质量矩阵[M]表示为[M]=m[M₀]，m为常系数；柔度矩阵[δ]为刚度矩阵[K]的逆矩阵，[δ]=[K]^-1=δ[δ₀]，δ=1/k；

则自振频率方程为：

（1）

即

（2）

该自振频率方程是求解自振频率和振型的方程，令

，，上述自振频率方程就转化为代数特征值问题：

（3）

其中A为n阶实矩阵，为振型；

第二步，从自振频率方程中求解全部特征值 ，采用带“原点位移”的Schmidt正交化QR分解算法求解，对于非奇异的n阶实矩阵A，可以通过Schmidt正交化QR分解方法分解成正交矩阵Q和上三角矩阵R的乘积：

>

“带原点位移”的Schmidt正交化QR分解法具体算法如下：

① QR分解：将进行QR分解，即

② RQ构造：

其中为单位矩阵，这里按如下取值：

（5）

其中是不为零的数，非奇异矩阵A经过QR算法后的矩阵序列，本质收敛于上三角矩阵或块上三角矩阵，且对角块为1阶子块和2阶子块，1阶子块就是A的实特征值，2阶子块含有A的一对共轭复根或实重根；计算出全部特征值后，可直接计算建筑结构体系的自振频率：

第三步，根据自振频率方程解出的特征值来计算建筑结构体系的振型：通过1次求解非齐次线性方程组的方法，直接求解特征向量，特征向量标准化后就是所求的振型，该振型为质点的相对振幅。

在本发明建筑物结构动力分析用自振频率和振型的计算方法中，将已经计算出的一个实和特征向量代入式（3），计算结果一般存在误差项：

>

其中为误差，为确保方程求解的稳定性，可按奇偶列的不同按如下方法取值（其中是一个很小的数）：

（7）

利用公式（6）求解非齐次线性方程组直接求解得出特征向量，并同样方式求解建筑物中建筑结构体系的所有振型。

在本发明建筑物结构动力分析用自振频率和振型的计算方法中， 所述的刚度矩阵和质量矩阵中k为常系数，m也为常系数，其由建筑物中建筑结构体系的平面布置、质量分布、材料实际性能、施工质量和地基基础确定。

基于上述技术方案，本发明专利与现有技术相比具有如下技术优点：

1．与现有技术的迭代法（如瑞利-里兹法、反幂法等）相比，在求解自振频率时，本发明方法的第二步具有算法稳定性，不会出现遗漏的特征值，但计算精度是依赖于特定的计算机系统的一种算法，其计算精度依赖于计算机的位数和编译系统的有效数字位数。

2. 在求解振型时，本发明的计算方法为直接法，除了算法稳定外，计算速度也非常快。迭代法常采用的反幂法，假设迭代为s次，那么至少需要求解2s次线性方程组；而采用第三步的方法，仅需要求解1次线性方程组，所用时间是反幂法的1/2s，当模型的节点很大时，本发明的计算方法能节约大量的计算时间。

附图说明

图1是本发明建筑物结构动力分析用自振频率和振型的计算方法的流程示意图。

图2是本发明建筑物结构动力分析用自振频率和振型的计算方法中实施例1所针对的建筑物刚度矩阵和质量矩阵示意图。

图3是本发明专利的实施例1中建筑物的第一层楼板产生单位位移（其他各质点的位移保持为零）时结构的刚度系数示意图。

图4是本发明专利的实施例1中建筑物的第二层楼板产生单位位移（其他各质点的位移保持为零）时结构的刚度系数示意图。

图5是本发明专利的实施例1中建筑物的第三层楼板产生单位位移（其他各质点的位移保持为零）时结构的刚度系数示意图。

图6是本发明专利的实施例1计算出的建筑物结构第一主振型示意图。

图7是本发明专利的实施例1计算出的建筑物结构第二主振型示意图。

图8是本发明专利的实施例1计算出的建筑物结构第三主振型示意图。

具体实施方式

下面我们结合附图和具体的实施例来对本发明专利建筑物结构动力分析用自振频率和振型的计算方法做进一步的详细阐述，以求更为清楚明了地理解其计算原理和步骤，但不能以此来限制本发明专利的保护范围。

本发明是一种建筑物结构动力分析用全新的自振频率和振型的计算方法。如图1所示，该计算方法包括有如下步骤：

第一步，建立多自由度的建筑结构体系自振频率方程，并将其转化为特征方程问题。

第二步，从自振频率方程中求解建筑结构体系的自振频率，先采用带“原点位移”的Schmidt正交化QR分解算法来求解自振频率方程的全部特征值，然后计算自振频率。

第三步，根据自振频率方程解出的特征值来计算建筑结构体系的振型，通过1次求解非齐次线性方程组的方法直接求解特征向量，特征向量标准化后就是所求的振型。

在上述计算方法的第一步中，通过“柔度法”建立建筑结构体系自振频率方程，并将其转化为特征方程问题。按多自由度体系来计算建筑物的自由振动问题（如多层厂房的侧向振动、不等高排架、高层建筑的振动等），常在建筑中（包括高层建筑）应用十分广泛的采用“刚度法”或“柔度法”求解，前者是通过建立力的平衡方程求解，后者通过位移协调方程求解。本发明方法采用的是“柔度法”，而“柔度法”的自振频率方程可由“刚度法”的方程间接地导出，之后将自振频率的求解转化为特征方程问题的求解。

多自由度的建筑结构体系在无阻尼自由振动时，假设建筑结构体系只有平移自由度，其刚度矩阵[K]和质量矩阵[M]就可表示为：[K]=k[K₀]、[M]=m[M₀],k和m分别为常系数。也就是说，在假设建筑结构体系的平面布置、质量分布、材料实际性能、施工质量、地基基础等情况均能准确确定，这样一旦建筑物确定，k和m均为恒定常数。[K₀]和[M₀]分别为质点间的相对刚度和相对质量，而柔度矩阵[δ]是刚度矩阵的逆矩阵（[δ]=[K]^-1=δ[δ₀]，δ=1/k）。[K]中的元素k_ij为结构的刚度系数，是使点j产生单位位移（其他各点的位移保持为零）时在点i所需施加的力。结构的刚度系数，可以根据达朗伯原理，通过结构所受的力和结构的位移之间满足的刚度方程求解。本发明采用的“柔度法”所用刚度矩阵，与求内力、位移等建筑结构体系静力计算所用的刚度矩阵是一致的。

有了[K]、[M]和[δ]后，就可以建立建筑结构体系的“柔度法”自振频率方程（是求解自振频率和振型的方程）：

（1）

即

（2）

令 ，，则自振频率方程就转化为了求解代数特征值问题：

（3）

在上述公式中，A为n阶实矩阵，为振型。

在上述计算方法的第二步中，从式（2）的自振频率方程中求解全部特征值，是采用带“原点位移”的Schmidt正交化QR分解算法来求解的。对于非奇异的n阶实矩阵A，可以通过Schmidt正交化方法分解成正交矩阵Q和上三角矩阵R的乘积：

>

“带原点位移”的Schmidt正交化QR分解法具体算法如下：

① QR分解：将进行QR分解，即

② RQ构造：

其中为单位矩阵，这里按如下取值：

>

其中是不为零的数；

非奇异矩阵A经过QR算法后的矩阵序列，本质收敛于上三角矩阵或块上三角矩阵，且对角块为1阶子块和2阶子块，1阶子块就是A的实特征值，2阶子块含有A的一对共轭复根或实重根。所谓矩阵非奇异，就是指A的行列式不为零。

有了全部特征值后，就可以直接计算建筑结构体系的自振频率：

在上述计算方法的第三步中，根据自振频率方程解出的特征值来计算建筑结构体系的振型：通过1次求解非齐次线性方程组的方法，直接求解特征向量，特征向量标准化后就是所求的振型，这里的振型就是质点的相对振幅；

若将已经求得的一个实和特征向量代入式（3），计算结果一般存在误差项：

>

其中为误差。

在基于本发明计算方法进行的计算机程序设计中，为确保方程求解的稳定性，可按奇偶列的不同按如下方法取值（其中是一个很小的数）：

（7）

式（6）可采用求解非齐次线性方程组的方法直接求解得出，同理可求得其他全部特征值对应的特征向量。

实施例1

本实施例是利用本发明的计算方法来求图2所示多自由度刚架的自振频率和主振型。假设横梁的变形略去不计，第一、第二和第三层的层间刚度系数分别为和。刚架的质量都集中在楼板上，第一、第二和第三层楼板处的质量分别为2m、m和m。本实施例是采用“柔度法”来进行计算。

利用本发明的计算方法求解的过程如下：

第一步：刚架的刚度系数如图3、图4和图5所示，其中图3、图4和图5对应的(a)、(b)和(c)分别表示当第一层、第二层和第三层楼板产生单位位移（其他各层楼板的位移保持为零）时在各层楼板所需施加的力，即结构的刚度系数。刚度系数的计算过程从略,刚度矩阵和质量矩阵分别为：

[K]中的元素k_ij为结构的刚度系数，是使点j产生单位位移（其他各点的位移保持为零）时在点i所需施加的力。而刚架的柔度矩阵是刚度矩阵的逆矩阵：

其中

于是有：

所求自振频率方程为：

即

令

自振频率方程就转化为代数特征值问题：

第二步：从自振频率方程中求解自振频率

采用带“原点位移”的Schmidt正交化QR分解算法求解自振频率方程的全部特征值如下：

从自振频率方程中求解出全部特征值后，就可以根据下式直接计算自振频率：

自振频率为：

第三步：根据自振频率方程的特征值计算主振型

，求解对应特征向量，这里。这里没用按式（7）对方程进行预处理。

①列非齐次线性方程组

列非齐次线性方程组如下：

待求解的非齐次线性方程组如下：

②求解特征向量并标准化后得到振型

采用Schmidt正交化QR分解法求解非齐次线性方程组，求解的未单位化的向量及其范数如下：

单位化后有：

若按结构力学教材的标准化方法，只要将第三个元素归一化，第一主振型的计算结果与结构力学教材一致，如图6所示：

同理，可求得、对应的第二和第三主振型。

对于第二主阵型，待求解的非齐次线性方程组如下：

求解的未单位化的向量及其范数如下：

单位化后有：

若按文献的标准化方法，只要将第三个元素归一化，第二主振型的计算结果与结构力学教材一致，如图7所示：

对于第三主阵型，待求解的非齐次线性方程组如下：

求解的未单位化的向量及其范数如下：

单位化后有：

若按文献的标准化方法，只要将第三个元素归一化，第三主振型的计算结果与结构力学教材一致，如图8所示：

采用本发明的计算方法，其算法的稳定性体现在第二步，计算速度的提高体现在第三步，即在求解振型上，因为只需要求解1次非齐次线性方程组，所以它所用的时间是采用反幂法计算所用时间（假设迭代s次）的1/2s。若采用现有技术的反幂法计算n个主振型需要迭代s=1000次，每次求解线性方程组需要0.01秒，那么反幂法计算振型需要运行t=2*1000*0.01=20秒，而采用本发明的方法仅需要0.01秒（注：两种方法的线性方程组的求解时间是相同的），是反幂法计算振型所用时间的。

本发明的计算方法计算速度快，计算过程稳定，不会出现特征值遗漏的情况，但计算精度依赖于特定的计算机系统，其精度依赖于计算机的位数和编译系统的有效数字位数。同时，修改本发明的算法，不仅可以在桥梁等工程中应用，还可以求解有阻尼情况下的振动频率及对应振型。