改進的區間參數結構頻響函數迭代解法

范芷若， 姜　東， 董萼良， 費慶國

(1.東南大學 土木工程學院 工程力學系, 南京　210096;2. 江蘇省工程力學分析重點實驗室, 南京　210096; 3. 南京林業大學 機械電子工程學院, 南京　210096)

改進的區間參數結構頻響函數迭代解法

范芷若1, 2， 姜東3， 董萼良1, 2， 費慶國1, 2

(1.東南大學 土木工程學院 工程力學系, 南京210096;2. 江蘇省工程力學分析重點實驗室, 南京210096; 3. 南京林業大學 機械電子工程學院, 南京210096)

針對區間不確定性系統頻響函數求解提出了一種改進的迭代方法。采用區間分析方法描述不確定性結構參數，將頻響函數求解轉化為區間線性方程組形式，使用定點迭代方法，構造迭代格式；為了解決迭代方法中出現的收斂性問題，結合子區間概念提出改進的迭代方法，計算結構頻響函數的區間。以彈簧-質量系統為仿真算例，研究了子區間數目對計算結果的影響。最后，將方法應用于文獻中的工程結構，對比研究表明：該方法適用于全頻段頻響函數分析、具有較高的精度和效率。

不確定性；區間分析；區間線性方程；子區間；頻響函數

結構參數不確定性是工程中不可避免的問題。參數不確定性產生的原因很多[1]，結構物理量存在離散性，如楊氏模量、泊松比、密度等；在設計初始階段因掌握資料不足產生的不可預測性，如邊界條件和連接方式難以準確參數化。不確定性的存在對結構動態特性和結構動響應的影響不容忽視。

針對不確定性問題的概率方法研究較為廣泛[2]，這類方法在掌握不確定參數的概率密度函數或均值、方差等統計量后開展研究，其中最有效的方法是蒙特卡洛模擬法[3]，但這種方法要求大量確定性計算，在實際工程結構的分析應用中存在難度。此外，工程中難以得到不確定參數的概率分布或統計量，概率方法并不適用于描述某些特殊的不確定性。為了克服概率方法的缺點，出現了一些非概率的方法，其中區間分析與其他方法在理論上有密切聯系。Moore等[4]最早提出了區間分析的概念，用一個帶上下邊界的區間來定義參數的不確定性。區間有限元方法[5]是在考慮所有參數的不確定性后，分析結構響應的取值范圍。這類方法因其所需信息少、計算初始條件簡潔等特點，在不確定性結構響應分析中具有廣闊的應用前景。Ben-Haim等[6]提出了不確定性結構力學分析的凸模型方法。Chen等[7]采用區間量描述系統矩陣的不確定性，并提出區間攝動法，以減少區間參數的相關性在區間運算中導致的區間擴張問題。邱志平等[8]將區間攝動法應用于不確定性結構振動問題特征值求解中。在攝動法基礎上，李金平等[9]提出基于Neumann展開與攝動理論相結合的線性方程組的解法，提高了攝動法的準確性，但未解決攝動法在不確定性水平較高時的精度不足問題。蘇靜波等[10]將子區間方法與攝動法聯合應用，提高不確定性水平較高時的區間運算精度，并通過計算得出精度較高時子區間的數目。Muhanna等[11]提出一種結構靜力響應的區間線性方程解法，該算法是對線性方程形式的一種特定迭代解法。

頻響函數求解在結構動力學響應計算中至關重要。針對于帶區間參數結構頻響函數的求解，Moens等[12]提出基于模態疊加法原理的模態矩形法，用于計算帶區間參數結構頻響函數，其結果精度高，但是不適用于全頻段分析，且算法較復雜；Dessombz等[13]提出使用區間算法計算結構頻響函數區間。

本文提出一種改進的不確定性結構動態特性頻響函數的求解方法。將區間線性方程方法應用于帶區間參數結構頻響函數包絡值的求解，結合子區間概念提出改進解法以提高計算精度，采用彈簧-質量系統為仿真算例驗證方法的有效性，并將方法應用于工程結構分析。

1　區間頻響函數

1.1區間參數

區間參數是實數域R內的一個封閉集合，對于給定的區間參數xI∈ IR，可以將它表示為

(1)

(2)

令

(3)

得到區間參數的另一種常見表達形式

xI=x0(1+αI)

(4)

xI○yI=[min(xi○yi),max(xi○yi)],

○∈{+,-,×,/}

(5)

1.2不確定線性系統頻響函數

在不確定系統中，由于參數具有不確定性，將不確定性采用區間參數描述，則結構剛度矩陣KI∈IRn×n、質量矩陣MI∈IRn×n及阻尼矩陣CI∈IRn×n，IRn×n表示n×n維區間矩陣，代入多自由度線性系統的動力學方程

(-ω2MI+iωCI+KI)HI=I

(6)

式中：ω為圓頻率，HI為區間頻響函數，I為單位陣。

由此可將不確定線性系統頻響函數求解轉化為求解形如Ax=b的區間線性方程的解集

Σ(A,b)={x∈Rn|?A∈A,?b∈b:Ax=b}

(7)

其中

(8)

2　定點迭代解法

在區間線性方程的求解方面目前已有一些較成熟的方法[11]，例如區間高斯消元、區間Gauss-Seidel迭代、Krawczyk迭代和定點迭代等。本文采用定點迭代方法構造迭代格式求解區間線性方程。

2.1迭代格式

對方程進行如下處理，

R=inv(mid(A))=A0-1

(9)

解的估值為

x0=R×mid(b)

(10)

定義迭代初始條件

g=R×(b-A×x0)=-αRA1x0

(11)

G=[I]-R×A=-αRA1

(12)

迭代格式如下

(13)

設ε為迭代系數，對迭代格式做處理

(14)

將迭代格式進行展開，在進行了n次迭代后，解的格式可以表示為

(15)

(16)

由式(15)，(16)可看出，迭代最終結果中包含區間向量的乘方運算。區間運算規則決定了區間運算后會出現誤差累積，因此在迭代次數過多時，該方法會出現過度估計的結果。需要盡量減少迭代次數，以提高算法精度。這就要求迭代矩陣G具有強收斂性。

2.2改進的迭代解法

上述解決區間線性方程的迭代算法，是基于定點迭代方法的推廣，因此迭代矩陣G必須滿足收斂性條件ρ(|G|)<1，才能得到算法的保守解。

由式(12)可知，不確定參數區間寬度較小時，迭代矩陣G收斂性較好，線性方程迭代解法可以獲得較高精度，而對于參數區間寬度過大的系統，其迭代形式的收斂性就會降低，甚至出現不收斂的情況，迭代結果與真實結果相比誤差相應增大。

為解決這一問題，結合子區間概念對迭代格式進行改進，提出一種改進的迭代解法，可以有效的提高迭代精度，解決迭代收斂問題。

i=1,2,…,L

(17)

則區間參數矩陣A相應分為L個

(18)

對于每個Ai，區間寬度由α變為αi，A0的值變為A=A0+mid(αi)A1，A1不變。相應對于每個子區間線性方程，其迭代格式變為

Ri=inv(mid(Ai))

(19)

x0=Ri×mid(b)

(20)

gi=Ri×(b-Ai×x0)=-αiRiA1x0

(21)

Gi=[I]-Ri×Ai=-αiRiA1

(22)

(23)

利用結合子區間參數矩陣修改后的迭代格式進行求解，待求解的區間線性方程改寫為

Aixi=bx=∪xii=1,2,…,L

(24)

對于帶多個區間參數的系統，可將迭代中各參數做相應推廣(帶N個區間參數)

(25)

3　算例研究

3.1彈簧質量系統

考慮如圖1所示的帶有阻尼的6自由度彈簧質量系統，結構阻尼系數η=0.06。激勵點與激勵方向如圖1所示，確定性參數如表1所示。

表1　6自由度系統確定性參數

圖1　6自由度彈簧質量系統Fig.1 Spring mass system of 6 DOFs

結構系統的部分剛度k(N/m),質量m(kg)具有誤差或有界不確定性。設這些區間參數分別為

k1=[6 400, 9 600], m2=[0.009 6, 0.014 4]

對該系統計算質點1振動頻響函數。圖2、圖3分別給出了頻率從0～250 Hz時，區間迭代方法及其改進算法計算得到的結構頻響函數幅值、實部和虛部的上下界，以及蒙特卡洛方法計算得到的頻響函數上下界，其中改進算法中將每個區間參數分別劃分為16個子區間，蒙特卡洛抽樣樣本為20 000，蒙特卡洛的收斂解可被視為準確解[15]。

對比圖2(a)、圖3(a)系統頻響函數幅值圖，在頻響函數非峰值處，三種計算方法得到頻響函數幅值區間結果相近；在頻響函數峰值處，與蒙特卡洛方法對比，改進的迭代方法精度明顯高于區間線性方程法。

若分別考慮系統頻響函數的實部與虛部區間，如圖2(b)、2(c)、3(b)、3(c)，在對應頻響函數峰值的頻率點處，區間線性方程法與蒙特卡洛方法計算結果差距極大，而改進的迭代方法結果為蒙特卡洛方法吻合度較高。

圖2　區間線性方程法與蒙特卡洛方法結果對比Fig.2 Results comparison between Liner interval equation method and Monte-carlo method

圖3　改進的迭代方法(16子區間)與蒙特卡洛方法對結果比Fig.3 Results comparison between modified iteration method and Monte-carlo method

使用子區間攝動法[10]對質點1振動頻響函數區間進行求解作為參照。該將頻響函數峰值點處四種方法數值區間提取得表2。

在f=20 Hz頻率點處，隨著子區間劃分數目變化，改進的迭代方法結果與蒙特卡洛方法對比的誤差如圖4所示。

表2　四種方法在特定頻率點處幅值區間

由圖4可知，在子區間劃分數目超過8以后，改進的區間線性方程方法計算所得結果與蒙特卡洛方法計算結果誤差降低到5%以內，且隨著子區間數目增多，計算結果收斂。

表3　計算時間比較

圖4　20 Hz處誤差隨子區間劃分數目變化趨勢Fig.4 Change of deviation against the number of sub-intervals at 20 Hz

計算所需時間對比如表3，由時間對比可得改進的區間線性方程迭代解法運算時間遠少于蒙特卡洛方法。

綜上算例結果可得，區間線性方程解法的改進解法可以得到足夠的精度，計算時間相較于蒙特卡洛方法大幅減少，可以有效地提高計算效率。

3.2工程結構

將區間線性方程方法應用于文獻[12]中的工程結構，進行頻響函數區間包絡值求解。圖5為相應的確定性結構，為一個27自由度的二維梁模型。模型的各參數的屬性在表3中列出。其中節點3的質量以及單元6的截面積為區間不確定參數，取值區間為

m3=[144, 216]kg,A6=[3.0, 5.0]×10-4m2

表4　單元參數

表5　節點參數

對該模型節點3豎向自由度的頻響函數進行求解。圖6、圖7分別將區間線性方程法及其改進方法求出的頻響函數包絡值與蒙特卡洛方法結果進行對比。由于此為無阻尼系統，因此對頻響函數最終結果進行處理，去掉結果在(0, 10-11)以及(10-5,5)數值段的數據，保留其他數據進行對比分析。

圖6、圖7進一步證明了改進的區間線性方程方法的精度。與文獻[12]方法相比，改進的迭代方法在保證精度的情況下，不會出現文獻[12]中結構部分頻段頻響函數不能求解的情況，因此在區間參數結構頻響函數求解上具有更大的優越性。

圖5　二維梁模型Fig.5 Two dimension beam model

圖6　區間線性方程法與蒙特卡洛方法結果對比Fig.6 Results comparison between Liner interval equation method and Monte-Carlo method

圖7　改進的迭代方法(16子區間)與蒙特卡洛方法對結果比Fig.7 Results comparison between modified iteration method and Monte-Carlo method

4　結　論

本文對區間不確定工程結構的動力特性頻響函數的求解，給出了一種改進的迭代方法。將求解頻響函數的問題轉化為區間線性方程的形式，并按照區間的運算規則，運用定點迭代的方法，得出結構頻響函數上、下界的求解方法。運用子區間的思想，對迭代格式進行改進，保證迭代的精度和收斂性。通過對算例計算結果的對比和分析表明，文中的改進方法能求出較精確的解、計算效率高且能得到全頻段的響應區間。

[1] 蘇靜波, 邵國建. 基于區間分析的工程結構不確定性研究現狀與展望[J]. 力學進展, 2005, 35(3):338-344.

SU Jing-bo, SHAO Guo-jian. Current research and prospects on interval analysis in engineering structure uncertainty analysis [J]. Advances in Mechanics, 2005, 35(3): 338-344.

[2] 邱志平, 王靖. 不確定參數結構特征值問題的概率統計方法和區間分析方法的比較[J]. 航空學報, 2007, 28(3):590-592.

QIU Zhi-ping, WANG Jing. Comparison of probabilistic Statistical method and interval analysis method for Eigenvalue problem of structure system with uncertain parameters [J]. Acta Aeronautica Et Astronautica Sinica-Series,2007, 28(3):590-592.

[3] Metropolis N, Ulam S. The monte carlo method [J]. Journal of the American statistical association,1949,44(247):335-341.

[4] Moore R E, Kearfott R B, Cloud M J. Introduction to interval analysis [M]. Society for Industrial and Applied Mathematics, Philadelphia, 2009.

[5] Moens D, Hanss M. Non-probabilistic finite element analysis for parametric uncertainty treatment in applied mechanics: Recent advances [J]. Finite Elements in Analysis and Design, 2011, 47(1): 4-16.

[6] Ben-Haim Y, Elishakoff I. Convex models of uncertainty in applied mechanics [M]. Elsevier, 2013.

[7] Chen S H, Qiu Z P, Song D. A new method for computing the upper and lower bounds on frequencies of structures with interval parameters [J]. Mechanics Research Communications, 1994, 21(6): 583-592.

[8] 邱志平, 陳塑寰, 劉中生. 區間參數結構振動問題的矩陣攝動法[J]. 應用數學和力學, 1994(6): 519-527.

QIU Zhi-ping, CHEN Su-huan, LIU Zhong-sheng. Matrix perturbation method of structure vibration problems with interval parameters [J]. Applied Mathematics and Mechanics, 1994(6): 519-527.

[9] 李金平, 陳建軍, 朱增青,等. 結構區間有限元方程組的一種解法[J]. 工程力學, 2010, 27(4): 79-83.

LI Jin-ping, CHEN Jian-jun, ZHU Zeng-qing,et al. A method for solving the structural interval finite element equations [J]. Engineering Mechanics, 2010, 27(4): 79-83.

[10] 蘇靜波, 邵國建, 劉寧. 基于單元的子區間攝動有限元方法研究[J]. 計算力學學報, 2007, 24(4): 524-528.

SU Jing-bo, SHAO Guo-jian, LIU Ning. Static sub-interval perturbed finite element method based on the elements [J]. Chinese Journal of Computational Mechanics, 2007, 24(4): 524-528.

[11] Muhanna R L, Zhang H, Mullen R L. Interval finite elements as a basis for generalized models of uncertainty in engineering mechanics [J]. Reliable Computing, 2007, 13(2): 173-194.

[12] Moens D, Vandepitte D. An interval finite element approach for the calculation of envelope frequency response functions [J]. International Journal for Numerical Methods in Engineering, 2004, 61(14): 2480-2507.

[13] Dessombz O, Thouverez F, Lané J P, et al. Analysis of mechanical systems using interval computations applied to finite element methods [J]. Journal of Sound and Vibration, 2001, 239(5): 949-968.

[14] Rump S M. On the solution of interval linear systems [J]. Computing, 1992, 47(3/4): 337-353.

[15] 梁震濤, 陳建軍, 王小兵. 不確定性結構區間分析的改進MonteCarlo方法[J]. 系統仿真學報,2007,19(6):1220-1223.

LIANG Zhen-tao, CHEN Jian-jun, WANG Xiao-bing. Improved Monte Carlo method for interval analysis of uncertain structures [J]. Journal of System Simulation, 2007, 19(6): 1220-1223.

A modified iteration method to solve frequency response function of structures with interval parameters

FAN Zhi-ruo1,2, JIANG Dong3, DONG E-liang1,2, FEI Qing-guo1,2

(1. Department of Engineering Mechanics, Southeast University, Nanjing 210096, China; 2. Jiangsu Key Laboratory of Engineering Mechanics, Nanjing 210096, China;3. College of Mechanical and Electronic Engineering, Nanjing Forestry University, Nanjing 210096, China)

A modified iteration method for solving frequency response function of a structure with uncertain parameters was presented. Interval parameters were adopted to represent structural uncertain parameters. The equation to solve a structural frequency response function was converted to a group of interval linear equations. The iteration matrix was obtained by using the fixed point iteration theorem. Due to a convergence problem in iteration, sub-intervals were introduced and a modified iteration method was proposed to get the upper and lower bounds of the frequency response function. The accuracy of the solution was improved and the convergence problem was solved after modification. A spring-mass system was adopted in simulation to explore the influence of sub-interval’s quantity. At last, a comparative study was conducted using a structure in reference. The results showed that the proposed method is suitable for a full frequency band frequency response analysis with better accuracy and computation efficiency.

uncertainty; interval analysis; interval linear equation; sub-interval; frequency response function

10.13465/j.cnki.jvs.2016.13.004

國家自然科學基金(10902024);教育部新世紀優秀人才支持計劃(NCET-11-0086)

2015-04-27修改稿收到日期：2015-06-26

范芷若 女，碩士生，1991年生

董萼良 男，副教授，1966年生

TB122

A