非均質(zhì)結(jié)構(gòu)非平穩(wěn)隨機(jī)響應(yīng)的快速算法

非均質(zhì)結(jié)構(gòu)非平穩(wěn)隨機(jī)響應(yīng)的快速算法

陳玉震, 張盛, 陳飆松, 張洪武

(大連理工大學(xué)工業(yè)裝備結(jié)構(gòu)分析國家重點(diǎn)實(shí)驗(yàn)室運(yùn)載工程與力學(xué)學(xué)部工程力學(xué)系,遼寧大連116024)

摘要：將擴(kuò)展多尺度有限元法應(yīng)用于非平穩(wěn)隨機(jī)振動(dòng)時(shí)域顯式法中，實(shí)現(xiàn)了對(duì)非均質(zhì)結(jié)構(gòu)非平穩(wěn)隨機(jī)響應(yīng)的快速精確計(jì)算。首先，論文闡述了擴(kuò)展多尺度有限元法基本原理。其次，探討了時(shí)域顯式法在非平穩(wěn)隨機(jī)響應(yīng)分析中的優(yōu)勢。時(shí)域顯式法基于動(dòng)力響應(yīng)顯式表達(dá)式直接在時(shí)域進(jìn)行隨機(jī)振動(dòng)分析，較傳統(tǒng)反應(yīng)譜法，具有良好的計(jì)算精度和計(jì)算效率。最后，根據(jù)兩算法的特性和優(yōu)勢提出統(tǒng)一的多尺度算法框架，使其適用于非均質(zhì)結(jié)構(gòu)非平穩(wěn)隨機(jī)響應(yīng)分析的快速求解。數(shù)值算例驗(yàn)證了該算法的高效性和精確性。

關(guān)鍵詞：非均質(zhì)結(jié)構(gòu); 非平穩(wěn)隨機(jī)振動(dòng)； 擴(kuò)展多尺度有限元法; 基函數(shù)； 時(shí)域顯式法

中圖分類號(hào):TB122

文獻(xiàn)標(biāo)志碼：A

DOI:10.13465/j.cnki.jvs.2015.19.004

Abstract:A fast computation was implemented for non-stationary random responses of heterogeneous material structures by applying the extended multiscale finite element method (EMsFEM) into the time-domain explicit method(TDEM). Firstly, the fundamental principle of EMsFEM was presented. Then, the advantages of TDEM in non-stationary random response analysis were explored. Based on the explicit expressions of dynamic responses, the random response analysis was performed directly in time domain, it was shown that TDEM has good accuracy and efficiency compared with the response spectrum method. Finally, according to the characteristics and advantages of the two algorithms, a unified multiscale framework was proposed, it was applicable for non-stationary random response analysis of heterogeneous material structures. Numerical examples showed that the proposed method has high efficiency and accuracy.

基金項(xiàng)目：國家自然科學(xué)基金資助項(xiàng)目(51275540)；中央高?；究蒲袠I(yè)務(wù)費(fèi)資助項(xiàng)目(CDJZR13110001)

收稿日期：2014-08-07修改稿收到日期：2014-10-17

A fast algorithm for non-stationary random responses of heterogeneous material structures

CHENYu-zhen,ZHANGSheng,CHENBiao-song,ZHANGHong-wu(State Key Laboratory of Structural Analysis for Industrial Equipment, Department of Engineering Mechanics, Faculty of Vehicle Engineering and Mechanics, Dalian University of Technology, Dalian 116024, China)

Key words:heterogeneous material structures; non-stationary random response; extended multiscale finite element method (EMsFEM); base function; time-domain explicit method (TDEM)

自然存在以及人工形成的大部分材料都具有非均質(zhì)特性，例如地下巖土以及航空航天工業(yè)中廣泛使用的復(fù)合材料等，由這些非均質(zhì)材料構(gòu)成的結(jié)構(gòu)，在工程實(shí)踐中每時(shí)每刻都要受到各種形式的動(dòng)態(tài)載荷，其中，非平穩(wěn)隨機(jī)激勵(lì)是最為普遍的一種載荷形式，例如地震作用、陣風(fēng)載荷等，為了確保這些結(jié)構(gòu)的可靠性與安全性，研究非均質(zhì)結(jié)構(gòu)在非平穩(wěn)隨機(jī)激勵(lì)作用下的動(dòng)力響應(yīng)具有非常重要的意義。

在隨機(jī)振動(dòng)分析中，功率譜法是求解隨機(jī)響應(yīng)的主要方法[1-2]，對(duì)于多點(diǎn)(m點(diǎn))非平穩(wěn)隨機(jī)響應(yīng)問題，為了求解各離散頻點(diǎn)處的響應(yīng)功率譜，功率譜法在每個(gè)頻點(diǎn)都需要進(jìn)行約m次時(shí)程分析，對(duì)于大規(guī)模結(jié)構(gòu)，計(jì)算量巨大，工程中難以接受。時(shí)域顯式法[3-4](TDEM)從時(shí)域出發(fā)，根據(jù)線性動(dòng)力系統(tǒng)輸入與輸出之間的線性關(guān)系，建立各時(shí)刻結(jié)構(gòu)響應(yīng)關(guān)于隨機(jī)激勵(lì)的顯式表達(dá)式，從而可以直接利用一階矩和二階矩的運(yùn)算規(guī)律計(jì)算各時(shí)刻結(jié)構(gòu)響應(yīng)的均值和方差，對(duì)于在m個(gè)非平穩(wěn)隨機(jī)激勵(lì)作用下的結(jié)構(gòu)，顯式表達(dá)式的建立只需要進(jìn)行2m次時(shí)程分析，計(jì)算量小，工程中容易接受。時(shí)域顯式法中結(jié)構(gòu)響應(yīng)已經(jīng)完全解耦，只需要對(duì)所關(guān)心的響應(yīng)求解即可，實(shí)現(xiàn)了空間尺度降維。當(dāng)隨機(jī)響應(yīng)的均值和方差為時(shí)間的慢變函數(shù)時(shí)，時(shí)域顯式法可以在較大時(shí)間間隔上計(jì)算響應(yīng)的均值和方差，實(shí)現(xiàn)了時(shí)間尺度降維。時(shí)域顯式法僅需要非平穩(wěn)隨機(jī)激勵(lì)時(shí)域相關(guān)函數(shù)信息，所以可以摒棄各種功率譜模型，適用范圍更廣。時(shí)域顯式法是一種高效、精確的非平穩(wěn)隨機(jī)振動(dòng)分析方法。

時(shí)域顯式法在求解非均質(zhì)結(jié)構(gòu)的非平穩(wěn)隨機(jī)響應(yīng)時(shí)，為了捕獲結(jié)構(gòu)的微觀非均質(zhì)特性，結(jié)構(gòu)網(wǎng)格需要剖分的非常細(xì)致，這使得結(jié)構(gòu)規(guī)模巨大，甚至無法計(jì)算。多尺度方法[5]是求解非均質(zhì)結(jié)構(gòu)力學(xué)性能的一種有效途徑，本文選用擴(kuò)展多尺度有限元法[6](EMsFEM)，且為了能夠求解結(jié)構(gòu)動(dòng)力響應(yīng)，與傳統(tǒng)擴(kuò)展多尺度有限元法不同，該方法[7]是以各時(shí)間步內(nèi)等效靜力平衡方程為基礎(chǔ)，數(shù)值構(gòu)造出滿足局部特性微分算子并能反映結(jié)構(gòu)動(dòng)力性能的多尺度基函數(shù)，從而在宏觀尺度上對(duì)原問題進(jìn)行求解，很大程度地減少計(jì)算量。時(shí)域顯式法的主要計(jì)算來自2m次時(shí)程分析，將擴(kuò)展多尺度有限元法引入時(shí)域顯式法中，能夠很大程度的減少每次時(shí)程分析的計(jì)算量，從而使得原問題能夠快速求解。數(shù)值算例表明，擴(kuò)展多尺度有限法與時(shí)域顯式法的結(jié)合應(yīng)用，在保證精度的同時(shí)大大提高了計(jì)算效率，為非均質(zhì)結(jié)構(gòu)的非平穩(wěn)隨機(jī)響應(yīng)分析提供了一種快速精確的計(jì)算方法。

1擴(kuò)展多尺度有限元法

擴(kuò)展多尺度有限元法的基本思想是通過局部求解宏觀單元域內(nèi)的平衡方程，數(shù)值構(gòu)造出宏觀單元的多尺度基函數(shù)(形函數(shù))，這些基函數(shù)能反映宏觀單元內(nèi)部的微觀非均質(zhì)特征，從而只需在宏觀尺度上對(duì)原問題進(jìn)行求解即可以獲得滿意的結(jié)果，大大減少了計(jì)算量。擴(kuò)展多尺度有限元法的計(jì)算流程主要可分為三個(gè)部分[8-9]，即多尺度基函數(shù)構(gòu)造、宏觀計(jì)算、降尺度計(jì)算，本節(jié)以二維問題為例簡要介紹一下這一流程。

圖1　四節(jié)點(diǎn)宏觀單元 Fig.1 Four-node coarse element

1.1多尺度基函數(shù)構(gòu)造

見圖1的某平面四節(jié)點(diǎn)宏觀單元，其內(nèi)部微觀節(jié)點(diǎn)的位移可以表示為

(1)

式中:N為宏觀單元的多尺度基函數(shù)矩陣；u單胞內(nèi)所有微觀節(jié)點(diǎn)的位移向量；u′E為宏觀單元四個(gè)節(jié)點(diǎn)的位移向量。它們可以表達(dá)成如下公式

u=[u1v1u2v2…unvn]T

N=[RT1RT2…RTn]T

(2)

其中

(3)

式中,n為子網(wǎng)格的微觀節(jié)點(diǎn)總數(shù)，Nixy表示在宏觀節(jié)點(diǎn)i發(fā)生x方向單位位移時(shí)，宏觀單元的子網(wǎng)格上所有微觀節(jié)點(diǎn)產(chǎn)生的y向位移值。

1.2宏觀計(jì)算

當(dāng)多尺度基函數(shù)構(gòu)造完成后，即可推導(dǎo)得到相應(yīng)宏觀單元的等效剛度陣。見圖1，考慮子網(wǎng)格內(nèi)任意一個(gè)細(xì)網(wǎng)格單元e，其節(jié)點(diǎn)編號(hào)為(e1，e2，e3，e4)，單元e的彈性應(yīng)變能Πe可表達(dá)為

(4)

式中:Ke是單元e的常規(guī)剛度矩陣；ue是單元e的節(jié)點(diǎn)位移向量。

由式(1)和式(2)可得細(xì)網(wǎng)格單元e的節(jié)點(diǎn)位移向量與相應(yīng)的宏觀單元節(jié)點(diǎn)位移向量的關(guān)系為

ue=Geu′E,Ge=[RTe1RTe2RTe3RTe4]T

(5)

式中:Ge為細(xì)網(wǎng)格單元e的轉(zhuǎn)換矩陣，其表征細(xì)網(wǎng)格單元e和對(duì)應(yīng)的宏觀單元之間的節(jié)點(diǎn)位移映射關(guān)系。

將宏觀單元內(nèi)部的所有細(xì)網(wǎng)格單元應(yīng)變能力相加，可得到該宏觀單元的等效剛度陣

(6)

式中:p為宏觀單元內(nèi)部細(xì)網(wǎng)格單元的總數(shù)。

宏觀單元的等效剛度陣獲取后，就可以組裝總體剛度陣，建立宏觀層次的求解方程，進(jìn)而進(jìn)行宏觀尺度的求解。

1.3降尺度計(jì)算

宏觀計(jì)算后，即可利用式(1)得到宏觀單元內(nèi)部任一細(xì)網(wǎng)格單元e的節(jié)點(diǎn)位移ue，進(jìn)而可以得到該單元的應(yīng)力應(yīng)變

εe=Beue=Teu′E,Te=BeGe

(7)

和

σe=Deεe=Seu′E,Se=DeBeGe

(8)

依次在各宏觀單元內(nèi)計(jì)算，即可得到整個(gè)結(jié)構(gòu)的微觀位移Yk+1、微觀應(yīng)力σk+1和微觀應(yīng)變?chǔ)舓+1。

2非平穩(wěn)隨機(jī)振動(dòng)時(shí)域顯式法

2.1動(dòng)力響應(yīng)顯式表達(dá)式

對(duì)于n個(gè)自由度的線性體系，在m個(gè)激勵(lì)作用下結(jié)構(gòu)體系的運(yùn)動(dòng)方程可寫為

(9)

文獻(xiàn)[3-4]從精細(xì)積分法出發(fā)推導(dǎo)動(dòng)力響應(yīng)顯式表達(dá)式，本文從更一般的時(shí)程積分法Newmark法出發(fā)推導(dǎo)動(dòng)力響應(yīng)顯式表達(dá)式。記時(shí)間步長為Δt，時(shí)間步數(shù)為l，外載F(t)在時(shí)刻t0,t1,…,tl處的值分別記做F0,F1,…,Fl。結(jié)構(gòu)第k+1步的響應(yīng)可由第k步的響應(yīng)遞推得到，遞推公式為

Vk+1=TVk+WFk+1

(10)

式中,

(11)

其中

(12)

(13)

(15)

記做

V0=QF0

(16)

則由遞推公式可寫出時(shí)刻ti=iΔt時(shí)刻，結(jié)構(gòu)動(dòng)力響應(yīng)Vi的表達(dá)式

Vi=Ai,0F0+Ai,1F1+…+Ai,iFi

(i=1,2,…,l)

(17)

式中,Ai,0,Ai,1,…Ai,i的值為

(18)

各時(shí)刻響應(yīng)遞推公式展開后形式如下

V1=A1,0F0+A1,1F1

V2=A2,0F0+A2,1F1+A2,2F2

?

Vl=Al,0F0+Al,1F1+…+

Al,l-1Fl-1+Al,lFl

(19)

矩陣形式為

V=AR

(20)

式中:V=[VT1VT2…VTl]T，R=[FT0FT1…FTl]T。

第i時(shí)刻的響應(yīng)，可進(jìn)一步表示為

Vi=AiRi,(i=1,2,…,l)

(21)

其中

(22)

式(17)或式(21)即為結(jié)構(gòu)動(dòng)力響應(yīng)顯式表達(dá)式。

2.2響應(yīng)均值和方差

結(jié)構(gòu)動(dòng)力響應(yīng)顯式表達(dá)式只與結(jié)構(gòu)相關(guān)，與外載性質(zhì)無關(guān)，當(dāng)外載F(t)為隨機(jī)激勵(lì)時(shí)，動(dòng)力響應(yīng)顯式表達(dá)式依然成立。此時(shí)，隨機(jī)激勵(lì)F(t)在各時(shí)刻的離散量F0,F0,…,Fl，分別為一隨機(jī)激勵(lì)。根據(jù)一階矩與二階矩的運(yùn)算規(guī)律，由式(21)可得時(shí)刻ti處響應(yīng)Vi的均值向量和協(xié)方差矩陣分別為

E[Vi]=AiE[Ri](i=1,2,…,l)

(23)

cov[Vi,Vi]=Aicov[Ri,Ri]ATi

(i=1,2,…,l)

(24)

式中,E[Ri]和cov[Ri,Ri]分別為隨機(jī)激勵(lì)向量的均值向量和協(xié)方差矩陣，可分別由F(t)的均值函數(shù)向量和互相關(guān)函數(shù)矩陣求出。

式(21)中，結(jié)構(gòu)響應(yīng)已完全解耦，可只計(jì)算所關(guān)心的響應(yīng)量，即

vi=aiRi(i=1,2,…,l)

(25)

在結(jié)構(gòu)隨機(jī)振動(dòng)分析中，一般只對(duì)結(jié)構(gòu)某些關(guān)鍵部位的響應(yīng)感興趣，并不需要求解所有響應(yīng)的均值和方差。由式(25)即可實(shí)現(xiàn)對(duì)關(guān)鍵部位響應(yīng)均值和方差的單獨(dú)計(jì)算，如下式

μvi=E(vi)=aiE(Ri)(i=1,2,…,l)

(26)

(i=1,2,…,l)

(27)

這樣實(shí)現(xiàn)了對(duì)結(jié)構(gòu)的空間尺度降維，大大節(jié)省了計(jì)算量，提高了計(jì)算效率。

為了滿足積分精度和激勵(lì)描述的需要，時(shí)間步長為Δt不能太大，因此動(dòng)力響應(yīng)顯式表達(dá)式的推導(dǎo)是在一系列間隔較小的離散時(shí)間點(diǎn)上進(jìn)行的。但當(dāng)隨機(jī)響應(yīng)的均值和方差為時(shí)間t的慢變函數(shù)時(shí)，并不需要在上述每一個(gè)時(shí)刻計(jì)算響應(yīng)的均值和方差，即可在較大時(shí)間間隔上應(yīng)用式(26)和(27)，這樣實(shí)現(xiàn)了對(duì)結(jié)構(gòu)的時(shí)間尺度降維，進(jìn)一步提高了計(jì)算效率。

由上述過程可以看出，時(shí)域顯式法求解響應(yīng)均值和方差時(shí)，直接采用隨機(jī)激勵(lì)的時(shí)域相關(guān)函數(shù)信息，無需用到在時(shí)頻域中表征隨機(jī)激勵(lì)的時(shí)變功率譜密度函數(shù)信息，可以摒棄非平穩(wěn)隨機(jī)激勵(lì)的各種功率譜模型，適用范圍更廣。

2.3時(shí)域顯式法系數(shù)矩陣的快速計(jì)算

時(shí)域顯式法在隨機(jī)響應(yīng)分析時(shí)主要包含兩部分計(jì)算：一是動(dòng)力響應(yīng)顯式表達(dá)式中系數(shù)矩陣的計(jì)算；二是基于動(dòng)力響應(yīng)顯式表達(dá)式的隨機(jī)響應(yīng)計(jì)算。當(dāng)只求解結(jié)構(gòu)某些關(guān)鍵部位的響應(yīng)時(shí)，時(shí)域顯式法的計(jì)算量主要來自前者。而由式(18)可以發(fā)現(xiàn)，整個(gè)系數(shù)矩陣除第一列Ai,0與第二列Ai,1(i=1,2,…,l)需要單獨(dú)計(jì)算外，其他各列均可由前兩列遞推得到。

當(dāng)結(jié)構(gòu)受單點(diǎn)激勵(lì)作用(例如隨機(jī)地震作用)時(shí)，Ai,0和Ai,1為3n維向量。由式(17)及式(19)可知，Ai,0的值為在t=0時(shí)刻施加單位脈沖載荷F0=1，其他時(shí)刻載荷為0時(shí)，結(jié)構(gòu)在各時(shí)刻的響應(yīng)Vi。而Ai,1的值為在t=Δt時(shí)刻施加單位脈沖載荷F1=1，其他時(shí)刻載荷為0時(shí)，結(jié)構(gòu)在各時(shí)刻的響應(yīng)Vi。此時(shí)，僅需要2次時(shí)程分析即可得到系數(shù)矩陣前兩列，進(jìn)而遞推得到整個(gè)系數(shù)矩陣。

當(dāng)結(jié)構(gòu)受m(m≥2)點(diǎn)激勵(lì)作用時(shí)，Ai,0和Ai,1為3n×m維矩陣。此時(shí)Ai,0每列的值分別為t=0時(shí)刻在各激勵(lì)點(diǎn)源依次施加單位脈沖載荷F0,j=1(j表示第j個(gè)激勵(lì)點(diǎn)源，j=1,2,…,m)，其他時(shí)刻載荷為0時(shí)，結(jié)構(gòu)在各時(shí)刻的響應(yīng)Vi,j。而Ai,1每列的值分別為t=Δt時(shí)刻在各激勵(lì)點(diǎn)源依次施加單位脈沖載荷F1,j=1(j=1,2,…,m)，其他時(shí)刻載荷為0時(shí)，結(jié)構(gòu)在各時(shí)刻的響應(yīng)Vi,j。此時(shí)，僅需要2m次時(shí)程分析即可得到系數(shù)矩陣前兩列，進(jìn)而遞推得到整個(gè)系數(shù)矩陣。

3時(shí)域顯式法的多尺度框架

由第2節(jié)內(nèi)容可知，采用時(shí)域顯式法計(jì)算非平穩(wěn)隨機(jī)振動(dòng)的主要計(jì)算量來自2m次特定脈沖激勵(lì)下的時(shí)程分析，提高單次時(shí)程分析的計(jì)算效率成為提高整個(gè)算法效率的關(guān)鍵。對(duì)于非均質(zhì)結(jié)構(gòu)，采用多尺度方法是減小問題規(guī)模，提高計(jì)算效率的有效途徑，為了能夠求解非均質(zhì)結(jié)構(gòu)的動(dòng)力響應(yīng)，本文所采用的擴(kuò)展多尺度方法以各時(shí)間步內(nèi)等效靜力平衡方程為基礎(chǔ)來構(gòu)造多尺度基函數(shù)，其余過程與第1節(jié)相同。為了提高計(jì)算精度，本文采用多節(jié)點(diǎn)宏觀單元策略[10]。本節(jié)將詳細(xì)介紹一下將多尺度方法引入時(shí)域顯式法，構(gòu)建時(shí)域顯式法多尺度框架的基本流程?？蚣軋D見圖2，基本流程如下:

(1)讀取粗細(xì)網(wǎng)格信息，基于各時(shí)間步等效靜力平衡方程計(jì)算各宏觀單元?jiǎng)恿瘮?shù)Nd。

(4)構(gòu)造2m個(gè)脈沖激勵(lì)，1:m個(gè)激勵(lì)在t=0時(shí)刻施加，m+1:2m個(gè)激勵(lì)在t=Δt時(shí)刻施加，其幅值均為：F=JE，J為一n×m階常數(shù)矩陣，用于定位激勵(lì)，E為m×m階單位矩陣。

(5)初始狀態(tài)計(jì)算，根據(jù)構(gòu)造激勵(lì)計(jì)算宏觀初始狀態(tài)向量，并降尺度得到微觀初始狀態(tài)向量。

(6)各時(shí)間步內(nèi)等效載荷為

(28)

此時(shí)，載荷直接施加在微觀網(wǎng)格上，這將產(chǎn)生不平衡節(jié)點(diǎn)力，本文采用位移分解技術(shù)[10]來處理，真實(shí)位移u被分解成宏觀整體響應(yīng)u′和局部擾動(dòng)響應(yīng)u″，u′由宏觀等效力驅(qū)動(dòng)，u″則由微觀擾動(dòng)力驅(qū)動(dòng)，本步驟即計(jì)算第k步宏觀整體等效載荷及各單胞域內(nèi)微觀等效載荷。

(9)判斷是否達(dá)到最大時(shí)間步，完成整個(gè)時(shí)程分析計(jì)算。

(10)提取2m個(gè)時(shí)程響應(yīng)結(jié)果，組裝得到系數(shù)矩陣前兩列Ai,0與Ai,1(i=1,2,…,l)。

(11)由系數(shù)矩陣前兩列，遞推得到整個(gè)系數(shù)矩陣A。

(12)通過式(26)或(27)計(jì)算響應(yīng)均值和方差。

圖2　時(shí)域顯式法多尺度框架 Fig.2 The multiscale framework for TDEM

4數(shù)值算例

通過2個(gè)數(shù)值算例來驗(yàn)證本文算法的有效性。為了便于比較，兩算例均以虛擬激勵(lì)法的計(jì)算結(jié)果作為參考解。

兩算例均采用下列形式的均勻調(diào)制非平穩(wěn)激勵(lì)

F(t)=g(t)f(t)

(29)

式中,g(t)為均勻調(diào)制函數(shù)，有兩種形式

g1(t)=1

(30)

(31)

f(t)的自譜密度均取過濾白噪聲模型

(32)

式中,tb=8.5s,tc=20.0s,c1=0.1572,S0=142.75,ωg=19.07,ζg=0.544。

兩算例中，結(jié)構(gòu)尺寸及材料屬性等數(shù)據(jù)均為無量綱量。計(jì)算所取時(shí)間和頻率分別為t∈[0,40],Δt=0.5和ω∈[0,50],Δω=0.5。均采用瑞利阻尼，系數(shù)為α=0.015和β=0.02。

4.1算例1周期性非均質(zhì)結(jié)構(gòu)受非平穩(wěn)水平地面加速度響應(yīng)分析

見圖3(a)T型結(jié)構(gòu)，無量綱尺寸圖中已給出，細(xì)網(wǎng)格單元尺寸為10×10，粗網(wǎng)格單元尺寸為100×100，粗網(wǎng)格剖分及單胞見圖3(b)。結(jié)構(gòu)底部全約束，整體受一非平穩(wěn)水平地面加速度作用，均勻調(diào)制函數(shù)由式(31)給定，激勵(lì)自譜密度由式(32)給定。結(jié)構(gòu)由5種無量綱材料周期性構(gòu)成，各材料屬性見表1，計(jì)算右上頂點(diǎn)處的位移均方值。

圖3　T型結(jié)構(gòu)模型圖 Fig.3 The T-shape structure model

材料密度彈性模量泊松比E150003.0e110.4E230001.0e110.3E338002.0e110.35E432002.1e110.32E548002.5e110.38

粗網(wǎng)格單元邊界宏觀節(jié)點(diǎn)數(shù)依次取2、4、6和8，四種情況下虛擬激勵(lì)法(PEM)、時(shí)域顯式法(TDEM)以及擴(kuò)展多尺度時(shí)域顯式法(EMsTDEMm)的計(jì)算結(jié)果見圖4，其中m表示粗網(wǎng)格單元邊界宏觀節(jié)點(diǎn)數(shù)目。以虛擬激勵(lì)法(PEM)作為數(shù)值參考解，其他方法相對(duì)數(shù)值誤差通過L2范數(shù)定義，其公式如下：

(33)

式中,va表示PEM參考解，v表示TDEM及EMsTDEMm計(jì)算結(jié)果；n為時(shí)間步數(shù)。

圖4　端點(diǎn)位移均方值 Fig.4 Mean square of the endpoint displacement

由圖4及表2中的L2范數(shù)可以看出，TDEM與PEM的計(jì)算精度一致。隨著邊界宏觀節(jié)點(diǎn)數(shù)目的增加，EMsTDEMm的計(jì)算精度逐漸趨近PEM，當(dāng)m=8時(shí)，基本一致。

由表2中的計(jì)算時(shí)間可得，TDEM計(jì)算用時(shí)遠(yuǎn)遠(yuǎn)小于PEM。隨著邊界宏觀節(jié)點(diǎn)數(shù)目的增加，EMsTDEMm的計(jì)算用時(shí)逐漸增加，但效率較TDEM依然有極大的提高。

表2　L 2范數(shù)及計(jì)算用時(shí)

4.2算例2 隨機(jī)非均質(zhì)懸臂梁結(jié)構(gòu)受多點(diǎn)部分相干隨機(jī)激勵(lì)響應(yīng)分析

見圖5(a)T型結(jié)構(gòu)，無量綱尺寸圖中已給出，細(xì)網(wǎng)格單元尺寸為5×2.5，粗網(wǎng)格單元尺寸為50×25，粗網(wǎng)格剖分及單胞見圖5(b)。懸臂梁左端全約束，上邊界中點(diǎn)與右端點(diǎn)處施加兩部分相干隨機(jī)激勵(lì)。結(jié)構(gòu)由3種材料構(gòu)成(見表3)，且3種材料在結(jié)構(gòu)中完全隨機(jī)分布，計(jì)算右下端點(diǎn)處的位移均方值。

圖5　懸臂梁結(jié)構(gòu)模型圖 Fig.5 The cantilever structure model

材料密度彈性模量泊松比E130001.0e110.31E240002.0e110.35E350003.0e110.38

兩隨機(jī)激勵(lì)的均勻調(diào)制函數(shù)分別由式(31)和(30)給定，激勵(lì)自譜密度分別取式(32)及S0，兩激勵(lì)相干矩陣為

(34)

粗網(wǎng)格單元邊界宏觀節(jié)點(diǎn)數(shù)依次取2、5、8和11，四種情況下虛擬激勵(lì)法(PEM)、時(shí)域顯式法(TDEM)以及擴(kuò)展多尺度時(shí)域顯式法(EMsTDEMm)的計(jì)算結(jié)果見圖6，其中m表示粗網(wǎng)格單元邊界宏觀節(jié)點(diǎn)數(shù)目，L2范數(shù)及計(jì)算用時(shí)見表4。

圖6　結(jié)構(gòu)端點(diǎn)位移方差 Fig.6 The variance of the structure’s endpoint

方法自由度L2范數(shù)計(jì)算時(shí)間/sPEM8282—1325.00TDEM82820.000049.67EMsTDEM21100.14025.06EMsTDEM56740.05896.24EMsTDEM812380.01838.48EMsTDEM1118022.39e-712.80

由圖6及表4中的L2范數(shù)可以看出，TDEM與PEM的計(jì)算精度一致。隨著邊界宏觀節(jié)點(diǎn)數(shù)目的增加，EMsTDEMm的計(jì)算精度逐漸趨近PEM，當(dāng)m=11時(shí)，基本一致。

表4中的計(jì)算時(shí)間可得，TDEM計(jì)算用時(shí)遠(yuǎn)遠(yuǎn)小于PEM。隨著邊界宏觀節(jié)點(diǎn)數(shù)目的增加，EMsTDEMm的計(jì)算用時(shí)逐漸增加，但效率較TDEM依然有極大的提高。

此外，本文提出的擴(kuò)展多尺度時(shí)域顯式法(EMsTDEMm)在計(jì)算單點(diǎn)、多點(diǎn)完全相干或者激勵(lì)點(diǎn)源不多的部分相干非平穩(wěn)隨機(jī)振動(dòng)時(shí)，具有顯著優(yōu)勢，計(jì)算效率很高；但在處理激勵(lì)點(diǎn)源較多(m個(gè))的部分相干非平穩(wěn)隨機(jī)振動(dòng)時(shí)，需要計(jì)算2m個(gè)脈沖激勵(lì)下結(jié)構(gòu)的時(shí)程分析，計(jì)算量巨大，該方法會(huì)受到一定的限制，因此，需要進(jìn)一步從并行程序設(shè)計(jì)角度來解決此問題。

5結(jié)論

大型復(fù)雜結(jié)構(gòu)在非平穩(wěn)隨機(jī)激勵(lì)作用下的響應(yīng)分析仍是目前面臨的一個(gè)難點(diǎn)問題，當(dāng)結(jié)構(gòu)材料具有非均質(zhì)特性時(shí)，分析變得更加困難。針對(duì)這一問題，本文從非平穩(wěn)隨機(jī)響應(yīng)分析算法和非均質(zhì)結(jié)構(gòu)動(dòng)力響應(yīng)算法兩個(gè)方面著手，提出一種適用于非均質(zhì)結(jié)構(gòu)平穩(wěn)隨機(jī)響應(yīng)分析的快速算法。非平穩(wěn)隨機(jī)響應(yīng)分析算法方面采用時(shí)域顯式法，該方法從時(shí)域出發(fā)建立結(jié)構(gòu)動(dòng)力響應(yīng)顯式表達(dá)式，直接通過一階矩和二階矩計(jì)算結(jié)構(gòu)響應(yīng)的均值和方差。對(duì)于單點(diǎn)非平穩(wěn)隨機(jī)響應(yīng)分析，主要計(jì)算僅為2次時(shí)程分析，較功率譜法幾十至幾百次的時(shí)程分析，具有較高的計(jì)算效率。非均質(zhì)結(jié)構(gòu)動(dòng)力響應(yīng)算法方面，采用擴(kuò)展多尺度有限元法，該法以時(shí)程分析各時(shí)間步的等效靜力平衡方程為基礎(chǔ)，構(gòu)造能夠反映結(jié)構(gòu)動(dòng)力效應(yīng)的廣義多尺度基函數(shù)，從而在宏觀尺度上對(duì)原動(dòng)力問題進(jìn)行求解，較常規(guī)有限元法，大大減小了問題規(guī)模，提高了計(jì)算效率。本文通過采用兩種方法相互耦合計(jì)算，為非均質(zhì)復(fù)雜結(jié)構(gòu)的非平穩(wěn)隨機(jī)響應(yīng)分析提供了一條可行途徑。數(shù)值算例表明，所提出算法在保證計(jì)算精度的前提下，具有很高的計(jì)算效率。

參考文獻(xiàn)

[1]Lin J H, Zhang W S, Li J J. Structural responses to arbitrarily coherent stationary random excitations[J]. Computers & structures, 1994, 50: 629-633.

[2]Lin J H, Li J J, Zhang W S, et al. Random seismic responses of multi-support structures in evolutionary inhomogeneous random fields[J]. Earthquake Engineering & Structural Dynamics, 1997, 26: 135-145.

[3]蘇成, 徐瑞. 非平穩(wěn)激勵(lì)下結(jié)構(gòu)隨機(jī)振動(dòng)時(shí)域分析法[J]. 工程力學(xué), 2010, 12:77-83.

SU Cheng, XU Rui. Random vibration analysis of structures subjected to non-stationary excitations by time domain method[J]. Engineering Mechanics, 2010,12:77-83.

[4]蘇成, 徐瑞, 劉小璐,等. 大跨度空間結(jié)構(gòu)抗震分析的非平穩(wěn)隨機(jī)振動(dòng)時(shí)域顯式法[J]. 建筑結(jié)構(gòu)學(xué)報(bào), 2011, 11:169-176.

SU Cheng, XU Rui, LIU Xiao-lu,et al. Non-stationary seismic analysis of large-span spatial structures by time-domain explicit method[J]. Journal of Building Structures, 2011,11:169-176.

[5]Babuska I, Caloz G, Osborn E. Special finite element methods for a class of second order elliptic problems with rough coefficients[J]. SIAM J. Numer. Anal, 1994,31(4): 945-981.

[6]Zhang H W, Fu Z D, Wu J K. Coupling multiscale finite element method for consolidation analysis of heterogeneous saturated porous media[J]. Advances in Water Resources, 2009, 32(2): 268-279.

[7]Zhang H W, Lu M K, Zheng Y G, et al. General coupling extended multiscale FEM for elasto-plastic consolidation analysis of heterogeneous saturated porous media[J].International Journal for Numerical and Analytical Methods in Geomechanics,2015,39(1):63-95.

[8]張洪武, 吳敬凱, 劉輝,等. 擴(kuò)展的多尺度有限元法基本原理[J]. 計(jì)算機(jī)輔助工程, 2010, 19(2): 3-9.

ZHANG Hong-wu, WU Jing-kai, LIU Hui, et al. Basic theory of extended multiscale finite element method[J]. Computer Aided Engineering, 2010, 19(2): 3-9.

[9]張洪武, 吳敬凱, 付振東. 周期性點(diǎn)陣桁架材料力學(xué)性能分析的一種新的多尺度計(jì)算方法[J]. 固體力學(xué)學(xué)報(bào), 2011, 32(2): 109-118.

ZHANG Hong-wu, WU Jing-kai, FU Zhen-dong. A new multiscale computational method for mechanical analysis of periodic truss materials[J]. Chinese Journal of Solid Mechanics, 2011, 32(2): 109-118.

[10]Zhang S, Yang D S, Zhang H W, et al. Coupling extended multiscale finite element method for thermoelastic analysis of heterogeneous multiphase materials[J]. Computers and Structures, 2013, 121: 32-49.

第一作者楊洋女，碩士，1988年生

通信作者褚志剛男，博士，副教授，碩士生導(dǎo)師，1978年生