毛虎平，高鵬飛

（中北大學(xué) 能源動(dòng)力工程學(xué)院,太原 030051）

復(fù)雜結(jié)構(gòu)動(dòng)力學(xué)分析需要對(duì)巨大自由度的模型進(jìn)行計(jì)算，當(dāng)高頻激勵(lì)力作用時(shí)，要求計(jì)算步長(zhǎng)非常小，這將造成計(jì)算耗時(shí)指數(shù)級(jí)增加。為了提高計(jì)算效率，可在保證一定精度的情況下，用少量自由度模型代替大量自由度模型，即模型縮減。所謂模型縮減是通過一定的變換，將對(duì)總體結(jié)構(gòu)動(dòng)力學(xué)分析影響較小的次自由度用對(duì)總體結(jié)構(gòu)動(dòng)力學(xué)分析影響較大的少量自由度表示，達(dá)到減少自由度的目的，其中少量自由度就是主自由度。然而如何從龐大的自由度中選擇主自由度，目前在結(jié)構(gòu)動(dòng)力學(xué)領(lǐng)域仍屬極具挑戰(zhàn)性的問題。不過學(xué)術(shù)界目前提出一些選擇主自由度的原則，最具代表性的有：

(1)將結(jié)構(gòu)振動(dòng)方向定為主自由度；

(2)在質(zhì)量或轉(zhuǎn)動(dòng)慣量相對(duì)較大而剛度又相對(duì)較小的位置選擇主自由度；

(3)在施加力或非零位移的位置選擇主自由度。

這些原則在具體選擇主自由度時(shí)，僅僅是指導(dǎo)思想，隨意性較大。主自由度的位置和數(shù)目直接影響模態(tài)分析縮減質(zhì)量矩陣的精度。針對(duì)無阻尼系統(tǒng)的各種有效合理PDOFS選擇方法不能直接推廣到阻尼系統(tǒng)中的問題，提出基于自由度能量分布比值的阻尼系統(tǒng)PDOFs選擇方法，并通過算例說明了該方法的有效性和可靠性[1]。由于結(jié)構(gòu)的動(dòng)力行為受低階模態(tài)控制，不需要計(jì)算系統(tǒng)的所有特征值，提出了一種單元級(jí)能量估計(jì)方法，建立了簡(jiǎn)化的有限元模型，該方法有效地節(jié)省了計(jì)算時(shí)間，并能準(zhǔn)確預(yù)測(cè)全局系統(tǒng)的特征值[2]。羅虹等[3]提出2種主自由度選擇方法，以單層懸臂梁為對(duì)象分析了這2種方法的特點(diǎn)和適用范圍。劉孝保等[4]對(duì)多種主自由度選擇方案進(jìn)行了分析研究，表明結(jié)構(gòu)靜態(tài)縮減模態(tài)分析中誤差的最大影響因素是主自由度數(shù)量和分布，特別是主自由度分布。包學(xué)海等[5]以轉(zhuǎn)向架為分析對(duì)象，提出了選取主自由度的部分準(zhǔn)則。

在選擇主自由度后，面臨的問題就是模型縮減。最早的模型縮減法是Guyan[6]提出的，稱為Guyan縮減法，也稱為靜態(tài)縮減法。該方法忽略了自由度相關(guān)的慣性項(xiàng)和阻尼項(xiàng)，質(zhì)量矩陣、剛度矩陣、狀態(tài)向量和載荷向量均被分為主自由度和次自由度兩個(gè)部分，然后經(jīng)過矩陣變換，用包含主自由度的部分表示包含次自由度的部分?？紤]慣性項(xiàng)，就產(chǎn)生了改進(jìn)的縮減方法，通過Guyan縮減法獲得慣性項(xiàng)，其結(jié)果與整體結(jié)構(gòu)模態(tài)更加接近[7-8]。通過能量估計(jì)選擇候選元素來構(gòu)建簡(jiǎn)化系統(tǒng)，在所選元素關(guān)聯(lián)的候選自由度中通過依次刪除來選擇主自由度，該方法有效地節(jié)省了計(jì)算成本，并能從最低頻率到截尾頻率范圍內(nèi)恢復(fù)整個(gè)系統(tǒng)的高精度特征值[9]。在模態(tài)疊加法中，所求的需要疊加的模態(tài)與載荷完全沒有關(guān)系，實(shí)質(zhì)上有些模態(tài)可能貢獻(xiàn)很小，因此可以考慮用里茲向量疊加，因?yàn)槔锲澫蛄扛荏w現(xiàn)結(jié)構(gòu)動(dòng)態(tài)特性[10]。將簡(jiǎn)化系統(tǒng)與子域方案相結(jié)合，對(duì)每個(gè)子域進(jìn)行尺寸與形狀優(yōu)化，由于在每個(gè)子域中都采用了約簡(jiǎn)方案，因此該方案對(duì)大規(guī)模問題的設(shè)計(jì)優(yōu)化問題非常有效[11]。針對(duì)Guyan縮減法的不足，提出了改進(jìn)的思路并推導(dǎo)了相關(guān)公式[12]。將里茲向量疊加與靜態(tài)子結(jié)構(gòu)法結(jié)合起來，形成了動(dòng)力分析的子結(jié)構(gòu)方法[13]。在工程結(jié)構(gòu)損傷識(shí)別中應(yīng)用逐級(jí)近似模型縮聚法，采用改進(jìn)Guyan遞推縮聚法的一級(jí)縮聚模型獲得識(shí)別精度最高的結(jié)果[14]。提出了一種求解特征問題迭代凝聚的加速方法，采用序列消除法或能量法選擇主自由度，并尋求了不同縮聚過程中方程的系統(tǒng)推導(dǎo)和比較，矩陣更新不僅包含了逆迭代，還隱含了子空間變換[15]。用結(jié)構(gòu)各階模態(tài)的DC增益作為其價(jià)值判斷的準(zhǔn)則，進(jìn)行模態(tài)截?cái)嗖?shí)現(xiàn)模型縮減[16]。

通過文獻(xiàn)分析可以看出，人們對(duì)結(jié)構(gòu)動(dòng)態(tài)縮減法以及對(duì)其影響最大的主自由度的選擇方法均進(jìn)行了大量研究，然而主自由度的選擇一直以來沒有一種確切的精度高的方法，鑒于此，本文提出一種基于節(jié)點(diǎn)里茲勢(shì)能主自由度的結(jié)構(gòu)動(dòng)態(tài)縮減方法，利用里茲向量與結(jié)構(gòu)自身動(dòng)態(tài)特性和結(jié)構(gòu)所承受的載荷分布形態(tài)相關(guān)聯(lián)的特點(diǎn)，定義了節(jié)點(diǎn)里茲勢(shì)能的概念，并在此基礎(chǔ)上給出加權(quán)系數(shù)的公式，將兩者點(diǎn)乘獲得了節(jié)點(diǎn)里茲勢(shì)能向量，以其為依據(jù)選擇主自由度，最后用改進(jìn)的動(dòng)態(tài)縮減方法獲得小規(guī)模的結(jié)構(gòu)動(dòng)力方程并求解。

1 節(jié)點(diǎn)里茲勢(shì)能計(jì)算及主自由度選擇

模態(tài)疊加法是計(jì)算結(jié)構(gòu)動(dòng)態(tài)響應(yīng)的一種有效方法，其利用振型的正交性將動(dòng)力學(xué)方程解耦，分別求解每一個(gè)方程再疊加起來。然而在參與計(jì)算的模態(tài)中，有些對(duì)響應(yīng)影響較小，這是因?yàn)檎裥团c結(jié)構(gòu)所受載荷分布方式無任何關(guān)系。里茲向量是一組正交的并與載荷空間分布有關(guān)的向量，其是通過將質(zhì)量矩陣歸一化并正交化獲得的，能夠反映慣性力的影響。里茲向量構(gòu)造過程如下：

（1）求初始向量{x1}。其可通過對(duì)[M] 歸一化處理獲得，即

其中：

（2）構(gòu)造迭代式。

里茲向量與模態(tài)向量是相對(duì)應(yīng)的，在前k階里茲向量里分別取p1,p2,…,pk個(gè)最大里茲向量分量對(duì)應(yīng)的自由度，將其組合并刪除重合項(xiàng)，獲得最終的主自由度，這稱為基于里茲向量的主自由度選擇方法。

節(jié)點(diǎn)里茲勢(shì)能是指將模態(tài)空間轉(zhuǎn)換到里茲向量空間，用里茲向量與節(jié)點(diǎn)自由度質(zhì)量向量點(diǎn)乘而得到的勢(shì)能向量。節(jié)點(diǎn)質(zhì)量是將質(zhì)量矩陣每一行元素求和作為結(jié)構(gòu)有限元節(jié)點(diǎn)自由度的質(zhì)量。選取節(jié)點(diǎn)里茲勢(shì)能向量較大的分量對(duì)應(yīng)的自由度作為主自由度。然而這樣做會(huì)過分強(qiáng)調(diào)低階頻率，并忽略高階頻率，因此可以通過定義加權(quán)系數(shù)來減低計(jì)算結(jié)果過分集中于低階頻率的問題。節(jié)點(diǎn)里茲勢(shì)能計(jì)算公式如下

2 構(gòu)造縮減系統(tǒng)

其中：[K] 為n階剛度矩陣，[M] 為n階質(zhì)量矩陣，ω

其中：φp為主自由度對(duì)應(yīng)的特征向量；φs為次自由度對(duì)應(yīng)的特征向量；下標(biāo)p表示主自由度對(duì)應(yīng)的量，下標(biāo)s表示次自由度對(duì)應(yīng)的量。

結(jié)構(gòu)穩(wěn)態(tài)簡(jiǎn)諧響應(yīng)可表示為

其中：Ω為簡(jiǎn)諧激勵(lì)力頻率，f為激勵(lì)力。

如果fs=0，根據(jù)式（12）可以獲得精確的縮減關(guān)系

應(yīng)用二項(xiàng)式定理將式(13)展開并省略2階Ω以上的項(xiàng)可得

當(dāng)Ω=0時(shí)

那么，靜態(tài)縮減轉(zhuǎn)換矩陣為

通過以下近似來消除Ω

3 基于節(jié)點(diǎn)里茲能量主自由度的結(jié)構(gòu)動(dòng)態(tài)縮減實(shí)施

圖1說明了節(jié)點(diǎn)里茲能量主自由度的結(jié)構(gòu)動(dòng)態(tài)縮減流程，其實(shí)施步驟如下：

(1)構(gòu)造里茲向量。通過式（1）至式（5）構(gòu)造里茲向量，第1個(gè)里茲向量可以通過質(zhì)量矩陣對(duì)角元素與剛度矩陣逆矩陣相乘，并與質(zhì)量矩陣歸一化獲得。

(2)計(jì)算節(jié)點(diǎn)里茲勢(shì)能。節(jié)點(diǎn)里茲勢(shì)能由節(jié)點(diǎn)質(zhì)量與里茲向量對(duì)應(yīng)的分量乘積獲得，其表示結(jié)構(gòu)動(dòng)態(tài)特性的貢獻(xiàn)率，較大者說明貢獻(xiàn)大，較小者貢獻(xiàn)小，因此將此作為選擇主自由度的依據(jù)。然而，這樣做過分強(qiáng)調(diào)低階頻率，通過定義加權(quán)系數(shù)可以提高高頻的精度。

(3)主自由度選擇。通過式（6）計(jì)算獲得節(jié)點(diǎn)里茲勢(shì)能向量，選擇其分量較大者作為主自由度。

(4)構(gòu)造縮減系統(tǒng)。通過IRS方法，在靜態(tài)縮減法的基礎(chǔ)上，考慮結(jié)構(gòu)慣性力，該慣性項(xiàng)能使模態(tài)結(jié)果更加逼近完整模型的模態(tài)，它采用靜態(tài)縮減法獲得。式（20）是最重要的轉(zhuǎn)換矩陣，式（8）至式（19）為其推導(dǎo)過程。

(5)采用廣義Schur分解法[17]求解該矩陣。對(duì)于任意n階矩陣A，存在一個(gè)酉矩陣U，使得U′AU成為上三角矩陣，且該上三角矩陣的對(duì)角線元素為A矩陣的特征值，利用該性質(zhì)進(jìn)行縮減系統(tǒng)求解。

4 實(shí)例分析

4.1 圓柱曲板

圓柱形曲板半徑為100 mm，高為100 mm，兩側(cè)端固定，彈性模量為0.3 MPa，泊松比為0.3，密度為0.01 kg/m3，采用SHELL63單元，共有216單元，247個(gè)節(jié)點(diǎn)，1326個(gè)自由度，如圖2所示。

根據(jù)里茲向量和節(jié)點(diǎn)里茲勢(shì)能選擇主自由度，并將與主自由度相連的單元顯示出來（其中顏色較深的為主自由度對(duì)應(yīng)的單元），圖3為選擇90個(gè)主自由度時(shí)2種方法的不同效果。

圖1 基于節(jié)點(diǎn)里茲能量主自由度的結(jié)構(gòu)動(dòng)態(tài)縮減流程

從圖3可以看出，基于里茲向量法選擇的主自由度過于向結(jié)構(gòu)中心集中，這是過分強(qiáng)調(diào)低頻的一個(gè)重要表現(xiàn)；而基于節(jié)點(diǎn)里茲勢(shì)能法選擇的主自由度大部分也集中在結(jié)構(gòu)中心，但有一部分主自由度向兩側(cè)擴(kuò)散，這可能是提高高階固有頻率的一個(gè)重要現(xiàn)象。如圖4所示。

圖2 圓柱曲板幾何參數(shù)

圖3 圓柱曲板上選擇主自由度對(duì)應(yīng)的單元

圖4 根據(jù)不同主自由度選擇方法所得圓柱曲板計(jì)算結(jié)果相對(duì)誤差比較

從圖4可以看出，在選擇相同主自由度數(shù)的情況下，任意選擇主自由度誤差太大，而基于里茲向量法選擇主自由度誤差非常小，在前30階模態(tài)中，最大相對(duì)誤差不超過10%，然而基于節(jié)點(diǎn)里茲勢(shì)能法選擇主自由度誤差最小，在前30階模態(tài)中，最大誤差3%。

表1是對(duì)應(yīng)的固有頻率數(shù)值，從中可以看出本文方法的優(yōu)勢(shì)。如圖5所示。

從圖5可以看出，本文方法選擇不同主自由度數(shù)時(shí)，結(jié)構(gòu)模態(tài)的相對(duì)誤差變化隨著自由度數(shù)量的增加越來越小，當(dāng)主自由度數(shù)為450時(shí)，前30階模態(tài)相對(duì)誤差最大不超過2%；另一方面，自由度數(shù)量增加相同，而結(jié)構(gòu)模態(tài)相對(duì)誤差減小得越來越慢，450個(gè)自由度大約是全部自由度數(shù)的1/3，如果繼續(xù)增加主自由度數(shù)，雖然可以繼續(xù)減小誤差，然而對(duì)于工程應(yīng)用來說意義不大，因此在使用本方法時(shí)只要誤差滿足精度要求即可。

表1 根據(jù)不同算法所得圓柱曲板模態(tài)比較（400個(gè)主自由度）

4.2 曲軸

曲軸軸徑為6 mm，曲柄直徑為4 mm，曲柄長(zhǎng)4 mm，曲軸總長(zhǎng)34 mm，兩端固定，彈性模量為210 GPa，泊松比為0.3，密度為7850 kg/m3，采用十節(jié)點(diǎn)四面體結(jié)構(gòu)實(shí)體SOLID92單元，共有1022個(gè)單元，2127個(gè)節(jié)點(diǎn)，6231個(gè)自由度，如圖6所示。

圖5 根據(jù)不同主自由度數(shù)所得圓柱曲板計(jì)算結(jié)果相對(duì)誤差比較

圖6 曲軸幾何參數(shù)

圖7是根據(jù)2種方法選擇的500個(gè)主自由度對(duì)應(yīng)單元（顏色較深的單元），從中可以看出，基于里茲向量法選擇的主自由度集中在結(jié)構(gòu)中部，而基于節(jié)點(diǎn)里茲勢(shì)能法選擇的主自由度從中部向兩端擴(kuò)散，而且非常明顯。

圖7 所選擇曲軸主自由度對(duì)應(yīng)的單元

圖8是采用隨機(jī)選擇、基于里茲向量法和基于節(jié)點(diǎn)里茲勢(shì)能法分別選擇1800個(gè)主自由度并用IRS法構(gòu)造的縮減系統(tǒng)，求解獲得的曲軸模態(tài)相對(duì)誤差，可以看出，隨機(jī)選擇不可行，基于節(jié)點(diǎn)里茲勢(shì)能法選擇主自由度結(jié)果相對(duì)誤差最小，前30階模態(tài)中，最大相對(duì)誤差不超過10%，而基于里茲向量法選擇主自由度雖然優(yōu)于隨機(jī)選擇方法，但前30階的相對(duì)誤差中，最大誤差超過20%。

圖9為基于節(jié)點(diǎn)里茲勢(shì)能法選擇不同數(shù)量的主自由度的計(jì)算結(jié)果，進(jìn)一步說明隨著主自由度數(shù)量的不斷增加，結(jié)構(gòu)模態(tài)相對(duì)誤差不斷減小，但是減小的幅度越來越小，而主自由度數(shù)為2000時(shí)，在前30階中，結(jié)構(gòu)模態(tài)相對(duì)誤差最大不超過5%。

表2提供圖8對(duì)應(yīng)的具體數(shù)據(jù)，更能從數(shù)值上看到本文方法的優(yōu)勢(shì)。

表2 模態(tài)計(jì)算比較（2000個(gè)主自由度）

圖8 根據(jù)不同選擇主自由度方法所得曲軸計(jì)算結(jié)果相對(duì)誤差比較

圖9 根據(jù)不同主自由度數(shù)所得曲軸計(jì)算結(jié)果相對(duì)誤差比較

5 結(jié)語

本文通過里茲向量法和節(jié)點(diǎn)里茲勢(shì)能法選擇主自由度，并采用改進(jìn)縮減法構(gòu)造縮減系統(tǒng)，最后應(yīng)用廣義Schur分解法求解，通過圓柱曲板和曲軸分析說明了基于節(jié)點(diǎn)里茲勢(shì)能法選擇主自由度方法的有效性，并獲得以下結(jié)論：

(1)將模態(tài)空間轉(zhuǎn)換到里茲向量空間，可以用很少的幾個(gè)里茲向量捕捉到非常精確的動(dòng)態(tài)特性。由于其過分強(qiáng)調(diào)低階頻率，需進(jìn)一步通過定義加權(quán)系數(shù)來提高高階頻率的精度。

(2)節(jié)點(diǎn)里茲勢(shì)能法比里茲向量法能捕捉到更加合適的主自由度，從而能更好地反映結(jié)構(gòu)動(dòng)態(tài)特性，在同等條件下，能獲得更高的精度，為下一步高效的動(dòng)態(tài)響應(yīng)優(yōu)化奠定基礎(chǔ)。

(3)在結(jié)構(gòu)縮減中，用大約1/3的主自由度數(shù)量比較合適，如圓柱曲面板總自由度數(shù)為1326，而取450個(gè)主自由度獲得了前30階模態(tài)相對(duì)誤差不超過2%的精度；曲軸總自由度為6231，取2000個(gè)主自由度獲得前30階模態(tài)相對(duì)誤差不超過5%的精度。