巨建民 劉新明? 張書娜

(1.大連交通大學(xué)土木與安全工程學(xué)院, 大連 116028)(2.邢臺職業(yè)技術(shù)學(xué)院 建筑工程系,邢臺 054035)

一般阻尼動(dòng)力系統(tǒng)非奇異攝動(dòng)法

巨建民1劉新明1?張書娜2

(1.大連交通大學(xué)土木與安全工程學(xué)院, 大連 116028)(2.邢臺職業(yè)技術(shù)學(xué)院 建筑工程系,邢臺 054035)

借助矩陣攝動(dòng)理論,將模態(tài)疊加法運(yùn)用于一般阻尼矩陣的動(dòng)力學(xué)方程求解結(jié)構(gòu)的動(dòng)響應(yīng)是一種較為理想的方法.但當(dāng)系統(tǒng)的外荷載激振頻率接近于系統(tǒng)的固有頻率時(shí),直接將阻尼矩陣作為攝動(dòng)矩陣,會(huì)使解產(chǎn)生奇異,并導(dǎo)致求解失敗或誤差過大,這是因?yàn)槟B(tài)坐標(biāo)下的動(dòng)力學(xué)方程是無阻尼方程.為了解決這一問題,本文考慮在模態(tài)坐標(biāo)的動(dòng)力學(xué)方程中保留一定的阻尼.即將阻尼做分解,代入振動(dòng)方程,得到不同階次攝動(dòng)方程,再將攝動(dòng)方程變換到模態(tài)坐標(biāo),即采用非奇異攝動(dòng)方法.最后通過數(shù)值算例,得到一階、二階攝動(dòng),將其與精確解進(jìn)行比較.精度明顯得到改善,基本趨于精確解.從而驗(yàn)證了本方法的精確性和有效性.

動(dòng)力學(xué), 矩陣攝動(dòng), 阻尼矩陣, 非奇異, 模態(tài)變換, 攝動(dòng)量

引言

在工程問題中經(jīng)常遇到參數(shù)小變化的情形,而結(jié)構(gòu)參數(shù)小變化引起結(jié)構(gòu)振動(dòng)特性變化的問題,對工程結(jié)構(gòu)的優(yōu)化設(shè)計(jì)有重要意義.矩陣攝動(dòng)法就是用來研究由于結(jié)構(gòu)參數(shù)小變化而引起的結(jié)構(gòu)振動(dòng)的固有頻率和模態(tài)向量的變化的方法.文獻(xiàn)[1-2]分別介紹了結(jié)構(gòu)振動(dòng)特征值問題的矩陣攝動(dòng)法和退化系統(tǒng)振動(dòng)分析的矩陣攝動(dòng)法.之后Chen等[3]和Qiu等[4]分別基于一階矩陣攝動(dòng)法與一階參數(shù)攝動(dòng)法提出了計(jì)算區(qū)間參數(shù)結(jié)構(gòu)動(dòng)力響應(yīng)的方法.雖然采用一階區(qū)間分析方法可以對區(qū)間參數(shù)結(jié)構(gòu)動(dòng)力響應(yīng)進(jìn)行分析,但是當(dāng)參數(shù)的不確定量稍大時(shí),采用一階區(qū)間分析方法對結(jié)構(gòu)動(dòng)力響應(yīng)范圍進(jìn)行估計(jì)可能會(huì)失效. Qiu等[5]基于二階泰勒展開與區(qū)間運(yùn)算,對具有區(qū)間參數(shù)非線性結(jié)構(gòu)動(dòng)力響應(yīng)范圍問題進(jìn)行了討論. 然而由于此方法涉及到區(qū)間運(yùn)算,得到的結(jié)果可能對動(dòng)力響應(yīng)范圍過分高估,使得這一方法失去實(shí)際的工程意義.文獻(xiàn)[6]基于二階攝動(dòng)展開法與DC算法,提出一種估計(jì)區(qū)間參數(shù)結(jié)構(gòu)動(dòng)力響應(yīng)范圍的方法,避免了區(qū)間運(yùn)算帶來的對結(jié)構(gòu)動(dòng)力響應(yīng)的過高估計(jì).

阻尼一直是結(jié)構(gòu)動(dòng)力學(xué)中很難處理的一個(gè)問題.例如,模態(tài)迭加法是求解結(jié)構(gòu)動(dòng)態(tài)響應(yīng)較為有效的方法之一,但在考慮結(jié)構(gòu)阻尼時(shí),模態(tài)迭加法往往不能直接求解[7],然而,對大部分結(jié)構(gòu)而言,忽略結(jié)構(gòu)阻尼勢必造成很大誤差[8].過去,在大部分情況下,對于連續(xù)體結(jié)構(gòu),阻尼往往都被忽略掉了,這樣的結(jié)論對于結(jié)構(gòu)的設(shè)計(jì)及強(qiáng)度分析也往往不是很有說服力[9-10].考慮到在大部分情況下,阻尼比較小,阻尼對結(jié)構(gòu)動(dòng)力性能的影響是有限的,因此,動(dòng)力學(xué)方程中的阻尼項(xiàng)可作為數(shù)學(xué)模型中的攝動(dòng)項(xiàng).對于大型結(jié)構(gòu)問題,可使用矩陣攝動(dòng)方法對此類問題進(jìn)行求解.對于弱阻尼問題結(jié)構(gòu)頻率的求解已有文獻(xiàn)論述[11].關(guān)于阻尼動(dòng)力系統(tǒng)的模態(tài)分析矩陣攝動(dòng)法和求解動(dòng)力響應(yīng)問題的矩陣攝動(dòng)模態(tài)疊加法,可參考文獻(xiàn)[12-13].

借助矩陣攝動(dòng)理論,將模態(tài)疊加法運(yùn)用于一般阻尼矩陣的動(dòng)力學(xué)方程求解結(jié)構(gòu)的動(dòng)響應(yīng)不失為一種較為理想的方法.在大部分情況下直接將阻尼矩陣作為攝動(dòng)矩陣,都不會(huì)出現(xiàn)問題.但是,當(dāng)系統(tǒng)的外荷載激振頻率接近于系統(tǒng)的固有頻率時(shí),這種方法會(huì)使解產(chǎn)生奇異,并導(dǎo)致求解失敗或誤差過大.這是因?yàn)槟B(tài)坐標(biāo)下的動(dòng)力學(xué)方程是無阻尼方程.因此,為了解決這一問題,需要考慮在模態(tài)坐標(biāo)的動(dòng)力學(xué)方程中保留一定的阻尼.即采用非奇異攝動(dòng)方法.

1 非奇異攝動(dòng)理論

設(shè)一般多自由度系統(tǒng)的振動(dòng)方程為:

(1)

對于常見的工程問題,阻尼陣[C]的范數(shù)小于質(zhì)量矩陣和剛度矩陣.文獻(xiàn)[13]直接將[C]作為攝動(dòng)矩陣,有時(shí)會(huì)出現(xiàn)問題.為了在模態(tài)坐標(biāo)的動(dòng)力學(xué)方程中保留阻尼,將阻尼矩陣做如下分解.即:

[C]=[Cd]+[Cr]

(2)

其中,矩陣[Cd]和[Cr]由以下方法確定:

設(shè):[Cdm]=Diag([Φ]T[C][Φ])

(3)

即[Cdm]為阻尼矩陣在使用模態(tài)矩陣[Φ]經(jīng)過模態(tài)變換后矩陣主對角元素矩陣.而[Crm]被定義為：

[Crm]=[Φ]T[C][Φ]-[Cdm]

(4)

也就是說[Crm]為[Φ]T[C][Φ]的非對角線元素矩陣.將矩陣[Cdm]和[Crm]變換到物理坐標(biāo)后為:

[C]=[M][Φ][Cdm][Φ]T[M]+ [M][Φ][Crm][Φ]T[M]

(5)

比較式(2)和(5),即可得矩陣[Cd]和[Cr]的表達(dá)式:

[Cd]=[M][Φ][Cdm][Φ]T[M]

(6)

[Cr]=[M][Φ][Crm][Φ]T[M]

(7)

[Cd]可以被稱為主阻尼矩陣,[Cr]稱為剩余阻尼矩陣.一般來說,剩余阻尼矩陣的范數(shù)更小一些,可將其作為攝動(dòng)矩陣,因此定義:

[Cr]=ε[Cp]

(8)

(9)

依據(jù)攝動(dòng)理論,當(dāng)參數(shù)發(fā)生微小變化時(shí),其解也必發(fā)生相應(yīng)變化.設(shè)響應(yīng)[X]的解在攝動(dòng)后為:

[X]=[X0]+ε[X1]+ε2[X2]+……

(10)

在保留同量級誤差的前提下,可得不同階次的動(dòng)力學(xué)方程:

(11)

(12)

(13)

其余各階攝動(dòng)方程依次類推為：

(14)

從求解過程可以看出,其計(jì)算過程實(shí)為一迭代運(yùn)算.對于常見的工程問題,一般來說經(jīng)過幾次攝動(dòng)計(jì)算,基本都可以達(dá)到一個(gè)滿意的精度.

2 模態(tài)變換解耦和求解過程

模態(tài)疊加法是在模態(tài)坐標(biāo)下進(jìn)行動(dòng)力學(xué)方程的求解.因此,在實(shí)際運(yùn)行中需將攝動(dòng)方程(11～14)變換到模態(tài)坐標(biāo).設(shè)[q(t)]是方程模態(tài)坐標(biāo)下的響應(yīng)解.則:

[X(t)]=[Φ][q(t)]

(15)

模態(tài)坐標(biāo)下的各階攝動(dòng)的動(dòng)力學(xué)方程可寫為:

(16)

(17)

(18)

實(shí)際的運(yùn)行過程可按以下步驟進(jìn)行：

(1)不考慮阻尼項(xiàng),對無阻尼動(dòng)力系統(tǒng)進(jìn)行模態(tài)分析,確定模態(tài)矩陣[Φ]和頻率矩陣[Λ].

(2)按照式(3)和式(4)確定主阻尼矩陣[Cdm]和剩余阻尼矩陣[Crm],并計(jì)算攝動(dòng)阻尼矩陣[Crmp].

(3)確定模態(tài)坐標(biāo)下的初始條件.

(19)

(20)

(4)由式(16)～(18)計(jì)算模態(tài)坐標(biāo)下響應(yīng)的各階攝動(dòng)值[q0],[q1],[q2]等等.

(5)疊加求解總響應(yīng)值:

[q]=[q0]+ε[q1]+ε2[q2]+…

(21)

(6)將響應(yīng)值返回到原始坐標(biāo):

[X]=[Φ][q]

(22)

3 數(shù)值算例

阻尼系統(tǒng)動(dòng)力學(xué)問題的矩陣攝動(dòng)方法是為較大規(guī)模動(dòng)力學(xué)問題設(shè)計(jì)的.但為了對其計(jì)算效果和精度進(jìn)行分析,這里先以一個(gè)簡單的二自由度動(dòng)力方程為例.

某二自由度問題動(dòng)力方程如下：

其中ω=0.85rad/s

可以看出,阻尼矩陣較質(zhì)量矩陣和剛度矩陣的范數(shù)相對較小.所以,可以把阻尼矩陣作為攝動(dòng)矩陣,經(jīng)過模態(tài)分析后,可計(jì)算出如表1.

表1 非奇異攝動(dòng)法各階攝動(dòng)響應(yīng)解比較Table 1 Comparison of perturbation response solutions for nonsingular perturbation

從表1可以看出,不考慮阻尼的解與精確解都有一定差距.而通過非奇異攝動(dòng)法改進(jìn)后,可以看出,零階攝動(dòng)與精確解相比,精度明顯得到改善.一、二階攝動(dòng)后,基本趨于精確解.

在本例中,由計(jì)算得到系統(tǒng)的第一自振頻率為0.85,此值與外荷載激振頻率0.85相等,但所得到的解并未產(chǎn)生奇異,誤差也很小,所以說,一般阻尼動(dòng)力系統(tǒng)的非奇異攝動(dòng)法可以很好地解決當(dāng)外荷載激振頻率接近固有頻率時(shí),直接將阻尼矩陣作為攝動(dòng)矩陣所造成的解的奇異性和誤差過大問題.

4 結(jié)論

(1)大部分工程問題,其材料均為連續(xù)介質(zhì)體(如金屬材料、混凝土等),阻尼為結(jié)構(gòu)阻尼,具有弱阻尼特性.運(yùn)用矩陣攝動(dòng)理論,提出了依據(jù)攝動(dòng)解計(jì)算有阻尼問題動(dòng)力響應(yīng)的方法.該方法有效地克服了按常規(guī)無阻尼計(jì)算時(shí)造成的計(jì)算誤差.

(2)攝動(dòng)法采用了經(jīng)典的多自由度動(dòng)態(tài)模態(tài)分析理論,簡單而實(shí)用.與無阻尼問題模態(tài)迭加法相比,因?yàn)槭褂昧艘延械哪B(tài)和固有頻率計(jì)算,增加的計(jì)算量很小.所以,可以適用大型結(jié)構(gòu)的動(dòng)力計(jì)算.

(3)由算例可以看出,采用非奇異攝動(dòng)法,即使激振頻率等于自振頻率,也并未產(chǎn)生奇異,并且精度明顯得到提高,二階攝動(dòng)后基本趨于精確解.

(4)由非奇異攝動(dòng)法,將阻尼分解,使模態(tài)坐標(biāo)下的動(dòng)力學(xué)方程轉(zhuǎn)化為有阻尼方程,有效地解決了由于激振頻率等于系統(tǒng)固有自振頻率時(shí),直接將阻尼陣作為攝動(dòng)矩陣所造成的解的奇異性以及誤差過大問題.

1 陳塑寰. 退化系統(tǒng)振動(dòng)分析的矩陣攝動(dòng)法. 吉林工業(yè)大學(xué)學(xué)報(bào), 1981,11(4):11～18 (Chen S H. Matrix perturbation method for vibration analysis of degenerate system.JournalofJilinUniversityofTechnology, 1981,11(4):11～18 (in Chinese))

2 劉應(yīng)力,陳塑寰. 結(jié)構(gòu)振動(dòng)特征值問題的二階矩陣攝動(dòng)法. 吉林工業(yè)大學(xué)學(xué)報(bào), 1983,29(3):1～11 (Liu Y L, Chen S H. The second order matrix perturbation method for structural vibration eigenvalue problem.JournalofJilinUniversityofTechnology, 1983,29(3):1～11 (in Chinese))

3 Chen S H, Lian H D, Yang X W. Dynamic response analysis for structures with interval parameters.StructuralEngineeringandMechanics, 2002,13(3):299～312.

4 Qiu Z P, Wang X J. Parameter perturbation method for dynamic responses of structures with uncertain-but-bounded parameters based on interval analysis.InternationalJournalofSolidsandStructures, 2005,42(18-19):4958～4970.

5 Qiu Z P, Ni Z. Interval modal superposition method for impulsive response of structures with uncertain-but-bounded external loads.AppliedMathematicalModelling, 2011,35:1538～1550.

6 李琦,邱志平,張旭東. 基于二階攝動(dòng)法求解區(qū)間參數(shù)結(jié)構(gòu)動(dòng)力響應(yīng). 力學(xué)學(xué)報(bào), 2015,47(1):147～153 (Li Q, Qiu Z P,Zhang X D. Second-order parameter perturbation method for dynamic structures with interval parameters.ChineseJournalofTheoreticalandAppliedMechanics, 2005,47(1):147～153 (in Chinese))

7 陳國平,韋勇. 有阻尼結(jié)構(gòu)現(xiàn)行振動(dòng)系統(tǒng)的模態(tài)綜合. 振動(dòng)工程學(xué)報(bào), 2003,16(4):442～445 (Chen G P,Wei Y. Modal synthesis of current vibration system with damping structure.JournalofVibrationEngineering, 2003,16(4):442～445 (in Chinese))

8 蒯行成,李永紅,任恒山,等. 一類振動(dòng)問題中阻尼矩陣的等效方法及應(yīng)用. 湖南大學(xué)學(xué)報(bào), 2004,31(1):95～98 (Kuai X C, Li Y H, Ren H S. A damping matrix equivalent method and its application in a kind of vibration problems.JournalofHunanUniversity, 2004,31(1):95～98 (in Chinese))

9 楊育臣. 非比例阻尼結(jié)構(gòu)動(dòng)力響應(yīng)分析. 哈爾濱：哈爾濱工業(yè)大學(xué), 2009 (Yang Y C. Dynamic response analysis of non-proportional damping structures. Harbin：Harbin Institute of Technology, 2009 (in Chinese))

10劉明輝. 非比例阻尼結(jié)構(gòu)動(dòng)力學(xué)計(jì)算方法及在衛(wèi)星設(shè)計(jì)中的應(yīng)用. 哈爾濱：哈爾濱工業(yè)大學(xué), 2009 (Liu M H. Computational methods of non-proportionally damped structures with applications to satellite Design. Harbin：Harbin Institute of Technology, 2009 (in Chinese))

11周君求,彭躍社. 求解弱非比例阻尼系統(tǒng)實(shí)模態(tài)解的阻尼矩陣攝動(dòng)法. 動(dòng)力學(xué)與控制學(xué)報(bào), 2008,6(5):322～326 (Zhou J Q, Peng Y S. Damping matrix perturbation method for solving real mode solutions of weakly non-proportional damping system.JournalofDynamicsandControl, 2008,6(5):322～326 (in Chinese))

12郭永剛,曲乃泗,董毓新. 在非比例阻尼情況下結(jié)構(gòu)動(dòng)力響應(yīng)的攝動(dòng)分析方法. 地震工程與工程振動(dòng), 1995,15(4):48～54 (Guo Y G, Qu N S, Dong Y X. A perturbation analysis method for structural dynamic response in the case of non-proportional damping.EarthquakeEngineeringandEngineeringVibration, 1995,15(4):48～54 (in Chinese))

13巨建民,蒙秋男. 具有一般阻尼的結(jié)構(gòu)動(dòng)力學(xué)模態(tài)分析的矩陣攝動(dòng)法. 計(jì)算力學(xué)研究與進(jìn)展. 中國力學(xué)學(xué)會(huì)青年工作委員會(huì)第三屆學(xué)術(shù)年會(huì)論文集.萬國學(xué)術(shù)出版社, 2000:152～156 (Ju J M, Meng Q N. Matrix perturbation method for structural dynamic modal analysis with general damping. Computational Mechanics Research and Development. Proceedings of the 3rd Annual Conference of the Chinese Youth Mechanics Association. Wanguo Academic Publishing House, 2000:152～156 (in Chinese))

? Corresponding author E-mail: liuxinming361@163.com

6 November 2016,revised 30 November 2016.

NON-SINGULAR PERTURBATION METHOD FOR GENERAL DAMPED DYNAMICAL SYSTEM*

Ju Jianmin1Liu Xinming1?Zhang Shuna2

(1.SchoolofCivilEngineeringandSafetyEengineering,DalianJiaotongUniversity,Dalian116028)

Using the matrix perturbation theory, it is an ideal method to apply the modal superposition method to the dynamic equation of general damping matrix to solve dynamic response of structure. However, when the frequency of the system is close to the natural frequency of the system, the damping matrix as the perturbation matrix make the solution singular and lead to the solution failure or large error, because the dynamic equation under modal coordinates is the undamped equation. In order to solve this problem, th consider to retain some damping in the dynamic equation of the modal coordinateshe damping is decomposed and substituted into the vibration equationobtain the different order perturbation equation, and the perturbation equation is transformed to modal coordinatesthe nonsingular perturbation method. Finally, the first-order and second-order perturbation solutions are obtained by numerical examples, which are compared with the exact solutions. The accuracy is improved, basically tend to be the exact solution, which verifies the accuracy and effectiveness of the method.

dynamics, matrix perturbation, damping matrix, non-singularity, modal transformation, perturbation

10.6052/1672-6553-2016-057

2016-11-06收到第1稿,2016-11-30收到修改稿.

? 通訊作者 E-mail: liuxinming361@163.com