大規模單點基礎激勵隨機振動分析及并行計算

范宣華,于晨陽，王柯穎，方 葉，肖世富

(中國工程物理研究院 總體工程研究所， 四川 綿陽 621900)

大規模單點基礎激勵隨機振動分析及并行計算

范宣華,于晨陽，王柯穎，方 葉，肖世富

(中國工程物理研究院 總體工程研究所， 四川 綿陽 621900)

基礎激勵作用下的隨機振動分析為結構動力學有限元計算中的一種主要分析類型，在實際結構分析中得到了廣泛應用?；谀B疊加法，對基礎運動激勵下的隨機振動核心算法理論進行了系統推導?；贘AUMIN框架并行數據結構，在PANDA平臺下完成了大規模算法設計和并行程序研發。結合數值算例對隨機振動響應分析模塊進行了對比驗證和大規模并行可擴展性測試。研究結果表明：所研發的隨機振動分析軟件與商業軟件的計算結果完全吻合，在并行求解能力上可達到上億自由度、上萬核，遠超出商業有限元軟件的計算水平。

單點基礎激勵；并行計算；隨機振動；PANDA平臺；并行可擴展性；模態疊加法

隨機振動是自然界和工程結構分析中常見的一類振動。諸如大型建筑物因地震、車輛因路面不平等引起的振動都屬于隨機振動分析的范疇。隨機振動不同于其他確定性振動，通常沒有固定周期，無法用確定性函數進行描述，但具有一定的統計規律，一般可以用功率譜密度等統計量進行描述[1]。

大型復雜裝備或建筑物經常會受到來自約束部位的一致基礎激勵作用，有限元隨機振動分析多以平穩隨機振動為前提，主要采用功率譜密度對結構的響應進行描述[2]。在目前的商業有限元軟件中，串行分析求解能力基本維持在百萬至數百萬自由度量級，即便是國內部分開放的商業有限元并行軟件，受技術禁運(開放核數最多在數百個以內，而且價格異常昂貴)和多年來串行編程機制下形成的內核影響，僅計算規模略有增加，但并行性能往往難以提升，極大制約了復雜結構的精細數值模擬水平[3]。

并行計算是解決這類大規模有限元問題的一個重要途徑。作為近20年來隨著計算機硬件快速發展而形成的一個熱門研究方向，并行計算對提升計算分析的規模和精度、縮短數值模擬時間具有重要意義。本文以單點基礎激勵隨機振動分析為研究對象，對基于自主并行軟件平臺PANDA[4]的大規模有限元并行計算研究進行介紹，以期為大規模精細有限元分析提供借鑒。

1 基礎激勵隨機振動理論

隨機振動分析的輸入一般為作用于結構約束部位的載荷功率譜密度曲線，載荷可以是位移、速度、加速度等基礎激勵。本文從最為常用的基礎加速度激勵入手，對隨機振動的核心理論進行推導。

隨機振動算法理論的核心就是推導基礎激勵自功率譜密度與結構關注點自功率譜密度響應之間的傳遞關系。對于基礎運動情形，以結構和地基的相對位移xr作為變量，可建立如下多自由度系統的運動方程：

(1)

式(1)中xr為各節點自由度對應的相對位移向量，與絕對位移x以及基礎運動位移u之間滿足如下關系：

x=xr+du

(2)

對于運動方程(1)，采用模態疊加法[5]進行方程解耦。首先進行模態分析，與式(1)對應的廣義特征值方程為

Mφ=ω2Kφ

(3)

i=1,2,…,m

(4)

對于式(4)，在頻域內的穩態解可以寫為

(5)

其中Hi(ω)為第i階模態對應的頻響函數，即

(6)

對于平穩隨機振動，根據功率譜密度定義，相應的模態功率譜密度可以表示為：

(7)

式(7)給出了基礎運動加速度激勵自譜密度和模態功率譜密度之間的遞推關系。根據相對位移和模態位移的關系，可以得到有限元模型中第k個自由度相對位移功率譜密度Sxrkxrk和基礎激勵加速功率譜密度之間的關系：

φikφjkSηiηj=

(8)

其中φik為第i階振型在第k個自由度處的取值，上標*表示共軛。結合虛擬激勵法[6]思想，可對式(8)做進一步簡化，得到

(9)

再根據相對位移和絕對位移的關系表達式(2)，得到第k個自由度的絕對位移自功率譜密度：

(10)

式中dk為方向向量d在第k個自由度處的取值。式(9)(10)分別給出了第k個自由度相對位移和絕對位移自功率譜密度與基礎加速度激勵自功率譜密度之間的遞推表達式。根據求解得到的位移自譜密度，可以直接得到速度和加速度自功率譜密度，三者之間滿足如下轉換關系：

(11)

2 算法設計

根據本文的討論和理論推導，可以建立如下單點基礎加速度激勵作用下的隨機振動分析算法：

算法1 基礎加速度激勵單點隨機振動分析算法輸入:結構有限元模型,模態分析階數m,模態阻尼比ξi,加速度激勵功率譜密度曲線S¨u¨u(ω),激勵方向;輸出:指定節點自由度的位移、速度或加速度功率譜密度。1.根據結構有限元模型,在JAUMIN框架下進行離散,并行生成分布式質量矩陣M和剛度矩陣K;2.根據質量矩陣和剛度矩陣構造廣義特征值問題,進行模態分析,獲取的前m階特征對(ωi,?i);3.根據基礎運動作用方向在所有非約束自由度上構造單位方向向量d;4.結合整體質量矩陣、振型和方向向量求解基礎激勵對應的各階模態參與系數;5.以指定計算的節點自由度作為第一重循環,以加速度輸入自譜曲線離散后的頻率作為第二重循環: a.獲取每個計算自由度在各階模態振型中的取值和對應方向向量中的值; b.計算模態坐標系下的各階頻響函數Hi(ω)(式(6)); c.按照式(9)和式(10)進行模態疊加,分別計算相對位移和絕對位移功率譜密度曲線; d.根據需要按照式(11)分別計算速度和加速度的自譜密度曲線。6.功率譜密度曲線計算結果輸出。

在以上算法設計中，基礎激勵輸入可以是位移或速度的功率譜密度曲線，三者輸入之間同樣滿足式(11)的換算關系。此外，總體質量矩陣和剛度矩陣的生成、模態分析過程、響應求解等多個環節均可借助JAUMIN框架和PANDA平臺進行并行求解。

3 并行實現

3.1 JAUNMIN框架和PANDA平臺

以上單點基礎激勵隨機振動的并行實現主要借助中國工程物理研究院自主研發的JAUMIN并行計算框架[7]和PANDA平臺實現。JAUMIN是根據超大機群硬件結構特點面向大規模非結構網格計算而研發的并行計算框架，提供基本的底層并行數據結構和矩陣向量操作等數學運算，負責不同數據塊之間的并行通信和數據管理等操作，同時為各類應用軟件提供各類接口。有關JAUMIN的詳細介紹參見文獻[7]。

PANDA是筆者團隊基于JAUMIN框架并行數據結構研發的結構力學有限元并行分析軟件平臺。PANDA平臺基本架構以及與JAUMIN框架之間的關系如圖1所示。

圖1 PANDA平臺的基本結構

目前整個PANDA平臺包括靜力學、模態和振動、沖擊動力學、多物理場耦合分析等多個有限元分析模塊，主要采用 C++ 和MPI編寫語言，包含代碼10萬余行。對于復雜工程結構，JAUMIN框架結合前后處理軟件完成結構的有限元建模和網格區域分解，PANDA平臺根據結構材料屬性和外載荷條件等完成矩陣并行組裝和求解。單點基礎激勵隨機振動分析只是模態和振動分析軟件中的一個基本分析類型，此外還包含模態分析、地震響應譜分析、諧響應分析以及多點隨機振動分析等多個動力學分析類型。

3.2 并行實現概述

首先，根據工程結構特點建立結構有限元模型。在基于JAUMIN框架的并行計算中，多采用自主前處理軟件Supermesh進行建模，此外也支持商業有限元軟件建模方式。對于建立的有限元模型，JAUMIN框架采用圖剖分功能，將有限元網格進行區域分解，分成多個子區域。在區域分解時綜合應用負載平衡技術，將有限元模型網格信息均勻分配到各個CPU計算節點，在各個CPU內，PANDA平臺將結合模型本構和材料、邊界等物理參數，并行生成分布式質量矩陣和剛度矩陣。

其次，利用生成的質量矩陣和剛度矩陣，PANDA平臺調用模態分析模塊進行模態分析并行計算，獲取模態疊加所需的固有頻率和模態振型等。模態分析在整個動力學分析過程中是最耗時間和資源的環節，其計算能力也基本決定了后續振動分析的計算能力。由于已經在PANDA平臺下進行了相應并行求解的集成實現，故目前支持Krylov-Schur算法[8]和Jacobi-Davidson算法[9]等開展大規模模態分析并行求解[4,10]。

再次，根據模態分析計算結果，結合JAUMIN框架提供的矩陣向量等并行操作運算，獲取各階模態參與系數。

最后，將模態參與系數、計算自由度對應的各階振型值、模態固有頻率和模態阻尼比等作為響應計算輸入，設計相應的模態疊加C++類，完成從基礎激勵到結構計算自由度直接頻響函數的求解，并按照式(9)(10)設計隨機振動分析C++類，實現相對位移和絕對位移自功率譜密度曲線的求解和計算結果輸出。

除矩陣組裝和模態分析階段的并行環節以外，在隨機振動分析實現中的并行環節主要有2個：① 求解模態參與系數; ② 計算各自由度的自譜密度曲線。在求解模態參與系數時，主要借助模態分析的質量矩陣以及模態振型在各個進程的分布式數據進行并行求解，而在計算各個節點自由度的自功率譜密度時，其計算過程是一個天然并行模式，每個CPU進程只負責本進程內節點自由度的自功率譜密度計算，相鄰進程之間不需要任何數據通信。

4 數值算例

為驗證研發模塊的正確性和并行可擴展性，本節以光機靶球隨機振動分析作為算例開展并行計算研究。光機靶球的有限元模型如圖2所示。采用四面體單元進行劃分，初始模型自由度數為340萬，模態分析提取前100階模態，采用Jacobi-Davidson算法求解，各階模態阻尼比取0.01，在支架4個腳底施加沿水平方向的基礎加速度激勵，激勵曲線為0～50 Hz范圍內的白直譜(加速度功率譜密度為1(m/s2)2/Hz的直線譜)。計算靶球上方桿件最上端的水平方向的位移功率譜密度。

圖2 光機靶球有限元模型

對于340萬自由度初始模型，采用商業有限元軟件ANSYS進行對比驗證，ANSYS和PANDA計算得到的前100階固有頻率保持小數點后3位有效數字一致。圖3給出了靶球頂端關注點的絕對位移功率譜密度對比情況，可以看出，ANSYS和PANDA計算得到的功率譜密度曲線幾乎完全重合，驗證了PANDA單點基礎激勵程序模塊的正確性。

為進一步驗證PANDA單點激勵隨機振動分析的并行可擴展性，在百萬億次大型機群上對靶球340萬自由度模型進行了自適應網格加密后的隨機振動計算，經過1次網格自適應加密后達到 1 700 萬自由度，經過2次加密后達到1.3億自由度。對于如此上千萬乃至上億自由度的有限元模型，已超出國內通用商業有限元軟件的計算能力，目前只能借助并行計算軟件完成。

圖3 絕對位移自功率譜密度曲線對比

對于靶球模型加密一次后得到的1 700萬自由度規模，CPU核數從64個一直測試到 8 192個，得到不同并行CPU核數上的并行計算時間和并行效率，如表1所示?？梢钥闯觯簩τ谠撚嬎阋幠?，在8 192核內，隨著核數的增加，計算時間持續下降，并未出現計算拐點。8 192核相對于64核的并行效率為24%，說明PANDA平臺下的單點隨機振動程序具有優異的并行可擴展性。

表1 靶球1 700萬自由度模型隨機振動分析計算時間和并行效率

對于加密2次得到的1.3億自由度規模，分別進行1 024核、2 048核、4 096核、8 192核以及10 000核的并行可擴展性測試，均成功算出。相應計算時間情況如表2所示?？梢钥闯觯涸?8 192核內隨著核數的增加，計算時間近乎線性減少，8 192核較1 024核的并行效率提高達84%，表現出非常優異的并行可擴展性；超過8 192核后，由于大型機群系統自身網絡分組配置原因導致跨域并行通信，并行效率有所下降，出現計算拐點，但萬核級相對于千核級并行效率仍高達60%以上。為了更直觀地表示分析模塊的并行可擴展性，繪制1 700萬自由度和1.3億自由度的并行加速比曲線，如圖4所示。

表2 靶球1.3億自由度模型隨機振動分析計算時間和并行效率

圖4 兩種規模的計算加速比曲線

從表1～ 2以及圖4可以看出：

1) 計算規模越大，PANDA隨機振動分析功能模塊的并行可擴展性越優異。這主要是因為隨著模型的增大，每個節點計算區域內的計算自由度數較邊界通信部分的自由度數大大增加，通信帶來的計算量相對變小。同理，對于同一計算規模，隨CPU核數的增加，各區域帶來的并行通信量增大，并行效率也逐漸下降。

2) 對于光機靶球之類復雜的算例，千萬自由度計算時間從64核的26 h減少到8 192核的1 h左右。上億自由度規模從1 024核的40余h減少到上萬核的6 h左右。這相對于串行計算而言是幾乎無法實現的，并行計算則不但大幅提升了計算規模，還大大縮短了數值模擬的計算時間，凸顯了并行計算的優越性；

3) 本文的數值算例表明隨機振動分析模塊具備“上億自由度、上萬核”的并行可擴展能力。

5 結束語

本文基于模態疊加法，利用中國工程物理研究院自主研發的JAUMIN框架和PANDA平臺，對單點基礎激勵隨機振動分析進行了算法設計和并行實現，結合光機靶球數值算例驗證了隨機振動分析模塊的正確性和并行可擴展性，取得了商業有限元分析軟件無法達到的大規模并行計算能力。

本文研發的單點基礎激勵隨機振動分析模塊僅是PANDA平臺下模態和振動分析軟件中的一個基本分析模塊，旨在說明基于JAUMIN框架研發的PANDA平臺具有超強的并行可擴展性。在單點激勵基礎上，筆者所在團隊近期已將隨機振動分析從單點分析擴展到了多點激勵分析，并將分析規模進一步提升到了10億自由度以上，限于篇幅，在此不做深入介紹。

[1] CHRISTIAN L..Random Vibration:Mechanical Vibration and Shock Analysis[C]//3rdedition.Published in Great Britain and the United States by ISTE Ltd and John Wiley & Sons,2014.

[2] CHOPRA A K.Dynamics of Structures-Theory and Applications to Earthquake Engineering[M]//4th edition.New Jersey:Prentice Hall,2012.

[3] FAN X,WANG K,XIAO S,et al.Some progress on parallel modal and vibration analysis using the JAUMIN framework[J].Mathematical Problems in Engineering,2015(2):1-8.

[4] 范宣華.基于Panda框架的大規模有限元模態分析并行計算及應用[D].北京：北京大學,2013.

[5] ITOH T.Damped vibration mode superposition method for dynamic response analysis[J].Earthquake Engineering & Structure dynamics,1973,2(1):47-57.

[6] 林家浩，張亞輝.隨機振動的虛擬激勵法[M].北京：科學出版社,2006.

[7] LIU Q K,ZHAO W B,CHENG J,et al.A programming framework for large scale numerical simulations on unstructured mesh[C]//IEEE 2nd International Conference on High Performance and Smart Computing.New York:[s.n.]， 2016:310-315.

[8] STEWART G W.A Krylov-Schur algorithm for large eigenproblems[J].SIAM Journal on Matrix Analysis and Applications,2001,23(3):601-614.

[9] SLEIJPEN G L G，VORST H A V.A Jacobi-Davidson iteration method for linear eigenvalue problems[J].SIAM Journal on Matrix Analysis and Applications,1996,17(2):401- 425.

[10] FAN X H,CHEN P,WU R,et al.Parallel computing study for the large-scale generalized eigenvalue problems in modal analysis[J].Science China Physics,Mechanics and Astronomy,2014,57(3):477- 489.

(責任編輯楊黎麗)

ParallelComputationofLarge-ScaleRandomVibrationAnalysisUnderSingle-PointMotion-BasedExcitation

FAN Xuanhua， YU Chenyang, WANG Keying, FANG Ye, XIAO Shifu

(Institute of Systems Engineering, China Academy of Engineering Physics, Mianyang 621900, China)

The random vibration analysis under the basic excitation is a kind of main analysis type in the finite element calculation of structural dynamics, and it has been widely used in practical structural analysis. Based on the modal superposition method, the random vibration algorithm theories under the basic motion excitation were systematically deduces. Based on the parallel data structure of JAUMIN framework, large-scale algorithm design and parallel program development were carried out under the PANDA platform. Numerical examples are given to verify the rightness as well as the large-scale parallel scalabilities of the random vibration response analysis. The results show that the random vibration analysis software developed by us is in good agreement with the calculated results of commercial software, and can reach hundreds of millions of degrees of freedom, tens of thousands of CPU cores in the parallel solution scalability, going far beyond the calculation level of commercial finite element software.

single-point base excitation; parallel computation; random vibration; PANDA platform; parallel scalability; mode superposition method

2017-06-21

國家自然科學基金面上資助項目(11472256)；科技部“高性能計算”重大專項課題(2016YFB0201005)；國防基礎科研計劃項目(C1520110002)；中國工程物理研究院院長基金、院發展基金和雙百人才基金資助項目(YZ2015011，2014B0202025，ZX04003)

范宣華(1981—),男,山東日照人,博士,研究員,主要從事結構動力學并行計算研究,E-mail：fanxh@caep.cn。

范宣華,于晨陽，王柯穎，等.大規模單點基礎激勵隨機振動分析及并行計算[J].重慶理工大學學報(自然科學)，2017(10):56-61,89.

formatFAN Xuanhua，YU Chenyang, WANG Keying, et al.Parallel Computation of Large-Scale Random Vibration Analysis Under Single-Point Motion-Based Excitation[J].Journal of Chongqing University of Technology(Natural Science)，2017(10):56-61,89.

10.3969/j.issn.1674-8425(z).2017.10.009

TB132

A

1674-8425(2017)10-0056-06