基于Chebyshev展開的區(qū)間穿孔板超材料分析1)

劉 堅(jiān) 雷濟(jì)榮 夏百戰(zhàn)

(湖南大學(xué)汽車車身先進(jìn)設(shè)計制造國家重點(diǎn)實(shí)驗(yàn)室，長沙410082)

基于Chebyshev展開的區(qū)間穿孔板超材料分析1)

劉 堅(jiān) 雷濟(jì)榮 夏百戰(zhàn)2)

(湖南大學(xué)汽車車身先進(jìn)設(shè)計制造國家重點(diǎn)實(shí)驗(yàn)室，長沙410082)

目前對于聲學(xué)超材料的傳輸特性分析和優(yōu)化大多是基于確定的數(shù)值和確定的模型，然而在實(shí)際工程和結(jié)構(gòu)設(shè)計中存在大量材料自身特性和幾何物理參數(shù)的不確定性.如果忽略這些不確定變量對聲學(xué)超材料傳輸特性分析和優(yōu)化過程的影響，得到的結(jié)果可能不正確.針對這一現(xiàn)狀，擬將切比雪夫區(qū)間模型引入多層穿孔板超材料，提出多層穿孔板超材料聲學(xué)透射率的區(qū)間切比雪夫展開--蒙特卡洛模擬法(interval Chebyshev expansion-Monte Carlo simulation method，ICE-MCSM).該方法采用截斷切比雪夫多項(xiàng)式近似擬合多層穿孔板超材料的聲學(xué)透射率響應(yīng)曲線，構(gòu)造聲學(xué)透射率響應(yīng)曲線的切比雪夫代理模型；然后采用蒙特卡洛模擬法(Monte Carlo simulation method，MCSM)隨機(jī)生成一定數(shù)量的不確定區(qū)間變量的樣本數(shù)據(jù)點(diǎn)，并將生成的不確定區(qū)間變量樣本數(shù)據(jù)點(diǎn)代入切比雪夫代理模型，預(yù)測單個不確定區(qū)間變量和多個不確定區(qū)間變量條件下的多層穿孔板超材料聲學(xué)透射率區(qū)間的上界和下界.數(shù)值分析結(jié)果表明,ICE-MCSM預(yù)測的聲學(xué)透射率變化區(qū)間的上界和下界與直接蒙特卡洛法(direct Monte Carlo simulation method，DMCSM)預(yù)測的聲學(xué)透射率上界和下界的結(jié)果非常接近.與DMCSM相比，ICE-MCSM具有更高的計算效率.因此，ICE-MCSM可有效且高效地分析不確定區(qū)間變量條件下多層穿孔板超材料聲學(xué)透射率傳輸特性，具有良好的工程應(yīng)用前景.

超材料,聲學(xué)透射率,切比雪夫展開,區(qū)間分析

引言

作為一個新興研究領(lǐng)域，超材料引起了極大的關(guān)注.超材料最開始起源于電磁復(fù)合材料.1968年，前蘇聯(lián)物理學(xué)家Veselago在理論上提出一種介電常數(shù)ε和磁導(dǎo)率μ同時為負(fù)的電磁材料[1].最近十多年，隨著實(shí)驗(yàn)條件的改善，研究者在電磁超材料的理論和實(shí)驗(yàn)上做出了一系列成果顯著的研究工作[2-3].根據(jù)聲波方程和麥克斯韋方程之間的相似性，研究人員通過類比電磁超材料，設(shè)計和制造出了聲學(xué)超材料[4-8].當(dāng)聲波或者彈性波在聲學(xué)超材料中傳播時，能夠表現(xiàn)出一些特殊的性質(zhì)，如負(fù)折射[9-10]聲聚焦[11]和超級聲通道[12]等.

聲學(xué)超材料的特性與其帶隙的產(chǎn)生有著密切聯(lián)系.帶隙的產(chǎn)生主要源于布拉格散射[13]和局域共振[14].因?yàn)椴祭裆⑸涞膸额l率一般出現(xiàn)于波長與晶格常數(shù)在同一數(shù)量級的頻率區(qū)域,所以其難點(diǎn)在于在較小的周期尺寸條件下得到頻率較低的帶隙.而基于局域共振機(jī)理的結(jié)構(gòu)[15]，其帶隙所對應(yīng)的波長可以遠(yuǎn)大于晶格常數(shù)，突破了布拉格散射機(jī)理的限制.最近Akozbek等[16]研究了多層穿孔板超材料的聲學(xué)特性，其特點(diǎn)在于保留了基于法布里--珀羅(Fabry Perot，F(xiàn)P)共振的亞波長穿孔板的超透射(extraordinary acoustic transmission)[17]特性，同時還產(chǎn)生了周期性材料所具有的帶隙特性[18].這種結(jié)構(gòu)可以廣泛應(yīng)用于帶通濾波器、醫(yī)學(xué)檢測、能量收集和噪聲控制.此外，多層穿孔板最大的優(yōu)點(diǎn)在于可以通過調(diào)節(jié)穿孔直徑的大小，改變穿孔板的材料參數(shù)來調(diào)節(jié)帶隙頻率和超透射頻率位置.

上述討論的聲學(xué)超材料的聲學(xué)傳輸特性都是基于確定的模型和確定的數(shù)值分析得到的.然而在實(shí)際工程設(shè)計中存在大量材料自身特性、幾何物理參數(shù)的不確定性，如果忽略這些不確定因素，對聲學(xué)超材料的聲學(xué)傳輸特性進(jìn)行分析，得到的結(jié)果可能不正確[19-21].

如果用實(shí)際復(fù)雜的模型來解決不確定性問題計算成本巨大，因?yàn)樵谇蠼獠淮_定性問題時涉及大量的計算.比如運(yùn)用有限元方法來計算汽車懸架系統(tǒng)的不確定性問題，每一次計算時間會達(dá)到幾個小時，而整個不確定性分析需要計算的次數(shù)達(dá)到上百次甚至上千次，整個不確定性分析過程非常耗時.因此，代理模型被廣泛用來代替復(fù)雜模型.近些年，一系列代理模型得到了快速的發(fā)展，比如多項(xiàng)式回歸模型、Kriging模型、徑向基函數(shù)(RBF)和多變量自適應(yīng)回歸樣條函數(shù)(MARS).Jin等[22]基于精確性、穩(wěn)健性、效率和操作簡便性這4個標(biāo)準(zhǔn)對比研究了這4種方法.結(jié)果顯示多項(xiàng)式回歸模型在效率和簡便性上具有較大的優(yōu)勢.Kriging,RBF和MARS對于非線性較強(qiáng)的問題具有較好的適應(yīng)性，但是其穩(wěn)定性不如多項(xiàng)式回歸模型.對于非線性較強(qiáng)的問題，傳統(tǒng)二階多項(xiàng)式回歸模型擬合精度不高，但是可以通過提高多項(xiàng)式的階數(shù)來提高精度.然而一般多項(xiàng)式擬合階數(shù)越高，容易產(chǎn)生龍格現(xiàn)象，但當(dāng)采用第一類切比雪夫多項(xiàng)式作為代理模型來近似擬合就能避免龍格現(xiàn)象.所以本文將采用第一類切比雪夫多項(xiàng)式來近似模擬聲學(xué)超材料的復(fù)雜模型有限元計算過程,從而獲得一個有效的代理模型.

目前處理不確定性問題的最常用辦法包括概率法和非概率分析法.由于概率法需要大量樣本數(shù)據(jù)點(diǎn)來構(gòu)造不確定變量的概率密度函數(shù)，因此當(dāng)樣本數(shù)據(jù)點(diǎn)較少時，該方法具有一定局限性.非概率分析法能通過較少的信息來獲得不確定性邊界，因此可以較好地處理樣本數(shù)據(jù)點(diǎn)較少的工程問題.區(qū)間模型是非概率可靠性模型中的一種.目前常用的區(qū)間模型分析方法包括：高斯消去法(Gaussian elimination method，GEM)、區(qū)間攝動法[23](interval perturbation method，IPM)、蒙特卡洛模擬抽樣法(MonteCarlosim-ulation method，MCSM)[24]和區(qū)間切比雪夫展開分析方法(interval Chebyshev expansion method，ICEM)[25].沈祖和[26]研究發(fā)現(xiàn)僅在矩陣對角元素占優(yōu)時，高斯消去法才具有較大的精度，且高斯消去法的精度受限于消去步驟，隨著消去步驟的增加，精度會逐漸降低.區(qū)間攝動法是一種高效的區(qū)間方法.將仿射法(affinearithmeticmethod，AAM)引入?yún)^(qū)間攝動法[27]，能在一定程度上緩解區(qū)間擴(kuò)張現(xiàn)象.但區(qū)間攝動法較大的近似誤差限制了它的應(yīng)用.MCSM是最簡單最穩(wěn)健的區(qū)間分析方法，根據(jù)MCSM概率收斂性，隨著樣本數(shù)據(jù)的增加，其計算精度會逐步提升.當(dāng)樣本數(shù)據(jù)趨于無窮大時，MCSM的計算結(jié)果逐漸收斂于真實(shí)解.近年來區(qū)間切比雪夫展開分析方法在不確定區(qū)間問題上已有較成熟的發(fā)展.Xia等[28]將切比雪夫區(qū)間方法應(yīng)用于不確定時域系統(tǒng)動態(tài)響應(yīng)分析中.Wu等[29]將區(qū)間切比雪夫展開分析方法應(yīng)用于不確定條件下汽車懸架的優(yōu)化.Yin等[30]將區(qū)間切比雪夫展開分析方法推廣到了中頻聲場響應(yīng)分析.

綜上所述，切比雪夫區(qū)間分析方法已經(jīng)成熟地運(yùn)用到了不確定時域系統(tǒng)、汽車懸架等的不確定性優(yōu)化以及中頻聲場響應(yīng)分析中.但是在聲學(xué)超材料領(lǐng)域內(nèi)，區(qū)間切比雪夫展開分析方法還沒有得到相應(yīng)的研究和應(yīng)用.本文將采取區(qū)間切比雪夫展開方法對聲學(xué)超材料進(jìn)行區(qū)間不確定分析.首先介紹了多層穿孔板超材料聲學(xué)透射率的有限元 (finitelement，F(xiàn)E)計算方法，并采用直接蒙特卡洛模擬法(direct Monte-Carlo simulation method，DMCSM)來預(yù)測多層穿孔板超材料聲學(xué)透射率的區(qū)間邊界.DMCSM是通過蒙特卡洛模擬法生成大量的樣本數(shù)據(jù)點(diǎn)，并將這些樣本數(shù)據(jù)點(diǎn)代入有限元方程預(yù)測多層穿孔板超材料聲學(xué)透射率的區(qū)間邊界.由于將所有的樣本數(shù)據(jù)點(diǎn)代入有限元方程預(yù)測聲學(xué)透射率區(qū)間邊界是一個計算成本較高的過程.因此，為了解決DMCSM計算成本高的問題，提出了區(qū)間切比雪夫展開--蒙特卡洛模擬法(ICE-MCSM).該方法采用截斷切比雪夫展開多項(xiàng)式近似擬合多層穿孔超材料聲學(xué)透射率響應(yīng)曲線，構(gòu)建聲學(xué)透射率響應(yīng)曲線的代理模型，代替有限元計算聲學(xué)透射率的方程；然后通過蒙特卡洛模擬方法生成大量的樣本數(shù)據(jù)點(diǎn)，并將這些樣本點(diǎn)代入代理模型來計算單個不確定區(qū)間變量和多個區(qū)間不確定變量聲學(xué)透射率的區(qū)間邊界.最后比較了DMCSM和ICE-MCSM兩種方法對聲學(xué)透射率響應(yīng)曲線進(jìn)行不確定區(qū)間分析的效率和精度.

1 多層穿孔板超材料有限元方法

1.1 多層穿孔板超材料基本模型

圖1是穿孔鋁板為與空氣層交替排列組成的多層穿孔板超材料周期性結(jié)構(gòu)，其鋁板單元晶胞的尺寸為dx和dy，面積為S2=dxdy，亞波長孔半徑為r，面積為S1=πr2.鋁板厚度為la，空氣層厚度為lb.cair和ρa(bǔ)ir分別是空氣聲速和空氣密度，θ是平面波與Z軸的夾角,即入射角.具體的多層穿孔板超材料模型結(jié)構(gòu)和材料參數(shù)見表1.

圖1 多層穿孔板超材料模型Fig.1 Model of multi-layer metmaterials

表1 多層穿孔板超材料模型結(jié)構(gòu)和材料參數(shù)Table 1 Parameters of multi-layer metmaterials

1.2 多層穿孔板超材料聲學(xué)透射率有限元計算方法

角頻率為 ω的諧振彈性波在線彈性、各向同性、小變形且無源均勻介質(zhì)中的傳播行為可用以下波動方程描述[31]

其中，r=(x,y,z)為位置矢量；u(r)為位移矢量，為矢量微分算符；為拉普拉斯算子；λ和 μ為材料的Lame常數(shù)，與縱波波速cl和橫波波速ct的關(guān)系為

方程(1)也可以寫成分量形式

整理式(3)可得

其中i,j=x,y,z，ui為位移矢量u(r)的分量.

平面波的傳播介質(zhì)是空氣，即在傳播過程中不存在橫波，只存在縱波.故式(1)的波動方程可簡化為

其中p為壓力，ρ為傳播介質(zhì)的密度.將式(5)寫成分量形式

其中i=x,y,z.

當(dāng)用有限元方法求解式(5)時，可將其寫成離散形式的廣義特征值方程

其中,K為整體結(jié)構(gòu)的剛度矩陣，可寫成K=B為應(yīng)變矩陣；Ve表示單胞的整個區(qū)域；M為質(zhì)量矩陣，可寫成為形函數(shù)矩陣；為單胞的位移矩陣.

對于有限周期結(jié)構(gòu)，通常用頻率響應(yīng)函數(shù)來描述其傳輸特性.如果輸入輸出端同時對其聲壓的平方進(jìn)行面積積分，則輸出輸入端的比值為一個無量綱的比例系數(shù)，即聲學(xué)透射率.

接下來可以利用有限元軟件直接求解式(7)的特征值方程.將特征值代入式(5)求出任意位置的聲壓值.多層穿孔板超材料聲學(xué)透射率t可表述為

其中，p0和p1分別表示輸入和輸出端的聲壓值；s0和s1分別表示輸入和輸出端的表面積.故可以利用式(8)求出任意頻率下的聲學(xué)透射率t.

當(dāng)n分別等于1,2,3和4時(n表示穿孔鋁板與空氣層組成的周期數(shù))，聲學(xué)透射率與頻率之間的響應(yīng)曲線如圖2所示.從圖2可以看出，當(dāng)周期數(shù)n大于等于2時，透射率中出現(xiàn)了帶隙.該帶隙對應(yīng)的聲學(xué)透射率系數(shù)會隨層數(shù)的增加而變小.另外，當(dāng)周期數(shù)n=2時，透射率在15kHz得到明顯的加強(qiáng)，當(dāng)周期數(shù)n增加到4時，由于FP共振，超透射波峰會隨之增加到3個.

圖2 周期數(shù)n分別等于1,2,3和4時聲學(xué)透射率響應(yīng)曲線Fig.2 Response curve of transmittance for periodicityn=1,2,3 and 4

為了更好地說明聲波在多層穿孔板超材料中的傳播效果，本文給出了在周期數(shù)n=2，入射平面聲波大小為1Pa時兩種不同頻率下(11kHz,15kHz)的聲波分布圖，如圖3.從圖3(a)可以看出其出射聲壓接近0.25Pa，聲學(xué)透射率為0.2，這與圖2所示的結(jié)果吻合.從圖3(b)可以看出,其出射聲壓為-1Pa，把其入射波和出射波聲壓值代入式(8)可以計算其聲學(xué)透射率為1，結(jié)果與圖2所示的有限元結(jié)果吻合.

圖3 入射波頻率分別為11kHz(a)和15kHz(b)時的聲壓分布圖(周期數(shù)n=2)Fig.3 Sound pressure fiel for plane wave incidence under the frequency of 11kHz(a)and 15kHz(b)(n=2)

1.3 不確定條件下多層穿孔板超材料聲學(xué)透射率有限元區(qū)間分析方法

本小節(jié)在聲學(xué)透射率有限元計算方法的基礎(chǔ)上，采用直接蒙特卡洛模擬法對多層穿孔板超材料的聲學(xué)透射率進(jìn)行區(qū)間分析.

假設(shè)不確定區(qū)間變量的個數(shù)為k，區(qū)間變量可表達(dá)為

(1)設(shè)不確定區(qū)間變量的個數(shù)為k，并確定不確定區(qū)間變量的上、下界，初始化循環(huán)次數(shù)i=1；

(2)假設(shè)k個不確定區(qū)間變量在給定區(qū)間內(nèi)服從均勻分布，并對k個不確定區(qū)間變量分別隨機(jī)生成n組樣本數(shù)據(jù)點(diǎn)

(3)通過有限元方法(式(8))計算不確定區(qū)間變量每一組數(shù)據(jù)樣本點(diǎn)的聲學(xué)透射率ti(f,xi)的值.循環(huán)次數(shù)i=i+1；

(4)判斷循環(huán)次數(shù)i.如果i≤n，返回步驟3繼續(xù)計算下一組樣本數(shù)據(jù)的值.如果i＞n，比較所有聲學(xué)透射率的值t1(f,x1),t2(f,x2),··,tn(f,xn)，得到所有頻率下最大聲學(xué)透射率的值tmax(f,xmax)和最小聲學(xué)透射率的值tmin(f,xmin)，即為聲學(xué)透射率區(qū)間邊界值.

圖4 DMCSM聲學(xué)透射率區(qū)間邊界分析流程Fig.4 Flow chart of DMCSM for bounds analysis of transmittance

MCSM具有概率收斂的特性.如果產(chǎn)生的樣本數(shù)據(jù)足夠大，DMCSM所求得的聲學(xué)透射率上、下界將收斂于真實(shí)區(qū)間的上、下界.因此DMCSM可以作為聲學(xué)透射率區(qū)間分析的參考解.DMCSM是直接將樣本數(shù)據(jù)點(diǎn)代入有限元方程來預(yù)測聲學(xué)透射率的區(qū)間邊界，這是一個高成本的計算過程.為了提高計算聲學(xué)透射率區(qū)間邊界的效率，下文將采用截斷切比雪夫級數(shù)來近似擬合有限元計算聲學(xué)透射率響應(yīng)曲線，再通過蒙特卡洛模擬方法隨機(jī)生成樣本數(shù)據(jù)點(diǎn)，然后將這些隨機(jī)樣本點(diǎn)代入聲學(xué)透射率切比雪夫代理模型來預(yù)測聲學(xué)透射率的區(qū)間邊界.相比DMCSM復(fù)雜的有限元計算方程，ICE-MCSM構(gòu)建的聲學(xué)透射率切比雪夫代理模型是關(guān)于區(qū)間變量的簡單函數(shù)，聲學(xué)透射率區(qū)間邊界的計算效率會得到顯著的提高.

2 不確定條件下多層穿孔板超材料聲學(xué)透射率ICE-MCSM

2.1 聲學(xué)透射率切比雪夫代理模型的構(gòu)建

當(dāng)不確定區(qū)間變量個數(shù)為k且其不確定區(qū)間范圍為xi∈[-1,1],(i=1,2,··,k)，k階切比雪夫多項(xiàng)式可以表示為

Michalewicz函數(shù)三維圖形如圖5所示.

圖5 Michalewicz函數(shù)Fig.5 Michalewicz function

本文主要通過平均誤差指數(shù)emean來描述切比雪夫代理模型的精度

式中，y和分別表示原函數(shù)聲學(xué)透射率的真值和用切比雪夫多項(xiàng)式計算出來的聲學(xué)透射率值.為了保持測試樣本點(diǎn)的一致性，本文采用Hamersley序列測試點(diǎn)產(chǎn)生1000個測試點(diǎn).

從圖 6可以看出隨著切比雪夫截斷階數(shù)的增加，切比雪夫代理模型的誤差值會隨著減小.誤差值越低說明其切比雪夫代理模型的精度越高.當(dāng)切比雪夫截斷階數(shù)為2階時，其平均誤差emean為0.8，當(dāng)截斷階數(shù)增加到6階時，其平均誤差emean為0.079.這說明當(dāng)截斷階數(shù)增加到一定時，切比雪夫代理完全可以滿足較高的精度要求.

圖6 切比雪夫代理模型平均誤差Fig.6 Average error of Chebyshev approximation

2.2 聲學(xué)透射率ICE-MCSM計算流程

當(dāng)用截斷切比雪夫多項(xiàng)式近似計算聲學(xué)透射率時，式(11)就被認(rèn)為是聲學(xué)透射率切比雪夫代理模型.因?yàn)閤是不確定區(qū)間變量，所以代理模型式(11)是關(guān)于區(qū)間的函數(shù).在求得式(12)中的切比雪夫代理模型的常系數(shù)向量后，聲學(xué)透射率切比雪夫代理模型就是一個關(guān)于區(qū)間變量的簡單區(qū)間函數(shù).相比有限元方法，其計算效率會有較大的提高.在構(gòu)建切比雪夫代理模型后，再把蒙特卡洛模擬法引入到聲學(xué)透射率切比雪夫代理模型中，就能非常高效地計算得到聲學(xué)透射率區(qū)間邊界.ICE-MCSM的主要流程如圖7所示.具體步驟如下：

(1)設(shè)不確定區(qū)間變量的數(shù)量為k，并確定不確定區(qū)間變量的上、下界，初始化循環(huán)次數(shù)i=1;

(2)假設(shè)k個不確定區(qū)間變量在給定區(qū)間內(nèi)服從均勻分布，并對k個不確定區(qū)間變量分別隨機(jī)生成n組樣本數(shù)據(jù)點(diǎn)

(3)根據(jù)式(12)求得聲學(xué)透射率截斷切比雪夫多項(xiàng)式的常系數(shù)向量，并構(gòu)建式(11)所示的聲學(xué)透射率切比雪夫代理模型；

圖7 ICE-MCSM聲學(xué)透射率區(qū)間邊界分析流程Fig.7 Flow chart of ICE-MCSM for bounds analysis of transmittance

(4)通過式(11)計算不確定區(qū)間變量每一組數(shù)據(jù)樣本點(diǎn)的聲學(xué)透射率ti(f,xi)的值.循環(huán)次數(shù)i=i+1；

(5)判斷循環(huán)次數(shù)i.如果i≤n，返回步驟4繼續(xù)計算下一組樣本數(shù)據(jù)的值.如果i＞n，比較所有的聲學(xué)透射率的值t1(f,x1),t2(f,x2),··,tn(f,xn)，得到所有頻率下聲學(xué)透射率的最大值tmax(f,xmax)和聲學(xué)透射率的最小值tmin(f,xmin)，即為聲學(xué)透射率的區(qū)間邊界值.

3 多層穿孔板超材料聲學(xué)透射率的ICEMCSM區(qū)間分析

3.1 單個不確定區(qū)間變量的聲學(xué)透射率的ICEMCSM區(qū)間分析

當(dāng)單個不同的結(jié)構(gòu)參數(shù)la和lb以及cair和穿孔半徑r分別作為區(qū)間變量時，其不確定范圍如表2所示.其中case 1,case 2和case 3分別表示三種不同的區(qū)間范圍不確定度，且其不確定度逐漸增加.分析頻帶范圍為1～20kHz，選擇4階切比雪夫級數(shù)進(jìn)行分析.在case 1，case 2和case 3三種不同區(qū)間范圍條件下，對4個不確定區(qū)間變量的聲學(xué)透射率區(qū)間分析結(jié)果如圖8～圖10所示.

表2 不確定區(qū)間變量的變化范圍Table 2 Uncertainty range of interval variables

從圖8～圖10可以看出，在case 1情況下，ICE-MCSM對所有的區(qū)間變量均與參考解DMCSM接近，在case 2和case 3情況下，在EAT頻率附近ICE-MCSM計算得到的聲學(xué)透射率區(qū)間上下界略微偏離參考解DMCSM，這是因?yàn)镋AT附近的頻率為系統(tǒng)特征頻率，與不確定區(qū)間變量有很強(qiáng)的非線性關(guān)系，所以聲學(xué)透射率在特征頻率附近對不確定區(qū)間變量較為敏感.

圖8 在case 1區(qū)間條件下區(qū)間變量la,lb,cair,r的聲學(xué)透射率區(qū)間邊界Fig.8 Interval bounds of transmittance forla,lb,cair,rin case 1

圖9 在case 2區(qū)間條件下區(qū)間變量la,lb,cair,r的聲學(xué)透射率區(qū)間邊界Fig.9 Interval bounds of transmittance forla,lb,cair,rin case 2

圖10 在case 3區(qū)間條件下區(qū)間變量la,lb,cair,r的聲學(xué)透射率區(qū)間邊界Fig.10 Interval bounds of transmittance forla,lb,cair,rin case 3

圖10 在case 3區(qū)間條件下區(qū)間變量la,lb,cair,r的聲學(xué)透射率區(qū)間邊界(續(xù))Fig.10 Interval bounds of transmittance forla,lb,cair,rin case 3(continued)

為了進(jìn)一步說明ICE-MCSM計算聲學(xué)透射率區(qū)間邊界的精度.這里以板厚區(qū)間la在case 3情況下為例，截斷階數(shù)取為4，5和6階.ICE-MCSM方法在EAT頻率附近(f=15kHz，f=15.5kHz)計算得到的聲學(xué)透射率區(qū)間邊界值與DMCSM參考解產(chǎn)生的相對誤差如表3所示.從表3可以明顯的看出,當(dāng)截斷階數(shù)為4時，相對誤差在3%～7%之間，隨著切比雪夫截斷階數(shù)的增加，ICE-MCSM精度會有所上升.當(dāng)截斷階數(shù)達(dá)到6時，用切比雪夫代理模型計算聲學(xué)透射率響應(yīng)曲線區(qū)間邊界基本與參考解重合.原因是隨著切比雪夫截斷階數(shù)的提高，聲學(xué)透射率響應(yīng)曲線的切比雪夫代理模型擬合精度會隨著提高.

計算成本是評估ICE-MCSM性能的另一個指標(biāo). ICE-MCSM的計算成本主要由兩部分構(gòu)成，第一部分是聲學(xué)透射率切比雪夫代理模型的構(gòu)建的計算成本，第二部分是將蒙特卡洛模擬方法均勻生成的樣本點(diǎn)代入切比雪夫代理模型計算聲學(xué)透射率區(qū)間邊界的計算成本.DMCSM的計算成本主要集中在直接將樣本點(diǎn)代入有限元方程進(jìn)行聲學(xué)透射率區(qū)間分析的過程，故樣本點(diǎn)數(shù)量直接決定了DMCSM的計算成本.為了得到一個較為合理的樣本點(diǎn)數(shù)量，本文對DMCSM進(jìn)行收斂性分析.在單個不確定區(qū)間變量條件下，以板厚區(qū)間la在case 3(15kHz)情況為例.當(dāng)數(shù)據(jù)樣本點(diǎn)的個數(shù)分別為100，300，500，700和900時，單個不確定區(qū)間變量條件下多層穿孔板超材料的聲學(xué)透射率上、下界在15kHz與數(shù)據(jù)樣本點(diǎn)的收斂關(guān)系如圖11所示.從圖11中可以明顯看出，當(dāng)抽樣點(diǎn)大于500時，聲學(xué)透射率上、下界分別收斂于0.9998和0.5339.

表3 不同截斷階數(shù)條件下ICE-MCSM的相對誤差(f=15kHz，f=15.5kHz)Table 3 Relative error of ICE-MCSM for dif f erent truncations(f=15kHz，f=15.5kHz)

圖11 基于DMCSM聲學(xué)透射率區(qū)間(a)上界(b)下界的收斂圖(15kHz)Fig.11 The covergence plot of transmittance for interval upper bound(a)and lower bound(b)based on DMCSM(15kHz)

ICE-MCSM和DMCSM兩種方法的具體計算成本如表4所示.單次有限元方法計算聲學(xué)透射率的時間為56s(3.30GHz Xeon(R)CPU E3 1230 v3),故DMCSM求出聲學(xué)透射率上、下界的時間成本為30000s.從表4可以看出，ICE-MCSM的總耗時是遠(yuǎn)遠(yuǎn)低于DMCSM的，這是因?yàn)镮CE-MCSM的計算成本主要集中在式 (12)和式(13)切比雪夫--高斯積分插值點(diǎn)的計算和常系數(shù)向量的計算.由于ICE-MCSM在構(gòu)建代理模型需要的插值點(diǎn)數(shù)量遠(yuǎn)小于DMCSM，故ICE-MCSM的計算效率遠(yuǎn)高于DMCSM.另外從表中可以看出在得到截斷切比雪夫級數(shù)后，再通過MCSM計算聲學(xué)透射率的區(qū)間上、下界的耗時只占用ICE-MCSM總耗時非常小的部分，這是因?yàn)槁晫W(xué)透射率切比雪夫級數(shù)代理模型是關(guān)于區(qū)間變量的簡單函數(shù)，其計算成本較低.最后，隨著切比雪夫截斷階數(shù)的提高，ICE-MCSM的總耗時是逐漸增加.但是相比于DMCSM的總耗時，是完全可以接受的.

表4 單個區(qū)間變量下ICE-MCSM和DMCSM聲學(xué)透射率區(qū)間邊界計算成本對比Table 4 Comparison of computational cost between ICE-MCSM and DMCSM for interval bounds of transmittance in the case of single variables

3.2 多個不確定區(qū)間變量的聲學(xué)透射率的ICEMCSM區(qū)間分析

本小節(jié)討論多個不確定區(qū)間變量條件下，多層穿孔板超材料聲學(xué)透射率的區(qū)間分析.假設(shè)4個不確定區(qū)間變量分別為鋁板的厚度la，空氣層的厚度lb,空氣聲速cair以及穿孔半徑r，其區(qū)間范圍大小為case 1,case 2和case 3.通過ICE-MCSM計算得到多層穿孔板超材料透射率響應(yīng)曲線區(qū)間上、下界如圖12～圖14所示.從圖12～圖14可以看出，在case 1和case 2情況下，當(dāng)截斷切比雪夫采用5階時，其聲學(xué)透射率區(qū)間上、下界能與DMCSM很好的吻合.在case 3情況下，在EAT頻率附近ICE-MCSM的邊界與DMCSM的邊界有略微偏差.接下來為了提高case 3情況下ICE-MCSM在EAT頻率附近的透射率響應(yīng)曲線的精度，把切比雪夫多項(xiàng)式截斷階數(shù)提高到6階，得到透射率響應(yīng)曲線區(qū)間上、下界如圖15所示.可以明顯看出，采用ICE-MCSM分析得到的聲學(xué)透射率區(qū)間上、下界與DMCSM分析得到的區(qū)間上、下界基本吻合.

圖12 Case 1條件下聲學(xué)透射率區(qū)間邊界(5階)Fig.12 Interval bounds of transmittance in case 1(5 order)

圖13 Case 2條件下聲學(xué)透射率區(qū)間邊界(5階)Fig.13 Interval bounds of transmittance in case 2(5 order)

圖14 Case 3條件下聲學(xué)透射率區(qū)間邊界(5階)Fig.14 Interval bounds of transmittance in case 3(5 order)

圖15 Case 3條件下聲學(xué)透射率區(qū)間邊界(6階)Fig.15 Interval bounds of transmittance in case 3(6 order)

為了進(jìn)一步說明多區(qū)間不確定變量條件下ICEMCSM的計算效率，以case 3為例，對比ICE-MCSM和DMCSM在計算聲學(xué)透射率區(qū)間上、下界的耗時(采用 10000個樣本點(diǎn))，其計算成本列在表 5中.從總耗時可以看出，ICE-MCSM的計算耗時遠(yuǎn)小于DMCSM.其主要原因是相比DMCSM中有限元方法計算聲學(xué)透射率區(qū)間邊界，ICE-MCSM中采用切比雪夫代理模型來計算分析聲學(xué)透射率的區(qū)間邊界，其效率會得到很大的提升.

表5 多個區(qū)間變量下ICE-MCSM和DMCSM聲學(xué)透射率區(qū)間邊界計算成本對比Table 5 Comparison of computational cost between ICE-MCSM and DMCSM for interval bounds of transmittance in the case of multi-variables

4 結(jié)論

基于切比雪夫理論基礎(chǔ)和多層穿孔超材料聲學(xué)透射率的有限元計算方法，本文提出了切比雪夫區(qū)間有限元方法(ICE-MCSM)來進(jìn)行不確定條件下超材料的區(qū)間分析.首先用截斷切比雪夫級數(shù)展開擬合聲學(xué)透射率響應(yīng)曲線，從而構(gòu)建聲學(xué)透射率切比雪夫代理模型；接著通過MCSM生成一定數(shù)量的隨機(jī)點(diǎn)，并將這些點(diǎn)代入聲學(xué)透射率切比雪夫代理模型,預(yù)測單個和多個不確定區(qū)間變量的聲學(xué)透射率的區(qū)間邊界.為了說明ICE-MCSM在不確定性條件下對于聲學(xué)透射率區(qū)間分析的有效性，本文把DMCSM方法作為參考解，對比分析了ICE-MCSM的計算精度和計算效率.對比結(jié)果表明隨著切比雪夫多項(xiàng)式截斷階數(shù)的提高，ICE-MCSM對聲學(xué)透射率區(qū)間分析的精度會提高.當(dāng)切比雪夫截斷級數(shù)提升到6階時，采用ICE-MCSM得到的聲學(xué)透射率區(qū)間上下界基本與采用DMCSM得到的上下界幾乎完全重合.這說明通過提高切比雪夫多項(xiàng)式截斷階數(shù)可以提高ICE-MCSM對于聲學(xué)透射率區(qū)間分析的精度.另外，在多個不確定區(qū)間變量條件下，對于區(qū)間范圍不確定度較大的情況，同樣通過提高切比雪夫多項(xiàng)式截斷階數(shù)可以提高ICE-MSCM對于聲學(xué)透射率的區(qū)間分析精度.雖然隨著切比雪夫多項(xiàng)式截斷階數(shù)的增加，ICE-MCSM的計算成本也會增加.但是與獲得的高精度相比，計算成本的增加可以接受.此外，ICEMCSM的計算成本遠(yuǎn)小于DMCSM方法，因此ICEMCSM可以推廣到多層穿孔板超材料的不確定區(qū)間分析問題中.

從不確定變量的個數(shù)分析上來看，隨著不確定變量個數(shù)的增加，ICE-MCSM的計算成本會呈指數(shù)增長.所以ICE-MCSM只適用于不確定區(qū)間變量適中的情況.對于多個不確定區(qū)間變量的問題，可以通過抽樣方法生成更加均勻的樣本點(diǎn)來減少構(gòu)建高階切比雪夫多項(xiàng)式的樣本點(diǎn),從而提高分析效率.

1 Veselago VG.The electrodynamics of substances with simultaneously negative values of ε andμ.Soviet Physics Uspekhi,1968, 10(4)：509

2 Schurig D,Mock JJ,Justice BJ,et al.Metamaterial electromagnetic cloak at microwave frequencies.Science,2006,314(5801)：977-980

3 Smith DR,Padilla WJ,Vier DC,et al.Composite medium with simultaneously negative permeability and permittivity.Physical Review Letters,2000,84(18)：4184

4 Sheng P,Mei J,Liu ZY,et al.Dynamic mass density and acoustic metamaterials.Physica B:Condensed Matter,2007,394(2)：256-261

5 Fang N,Xi DJ,Xu JY,et al.Ultrasonic metamaterials with negative modulus.Nature Materials,2006,5(6)：452-456

6 Graci′a-Salgado R,Garc′?a-Chocano VM,Torrent D,et al.Negative mass density and ρ-near-zero quasi-two-dimensional metamaterials：design and applications.Physical Review B,2013,88(22)：224305

7 Lee SH,Park CM,Seo YM,et al.Composite acoustic medium with simultaneously negative density and modulus.Physical Review Letters,2010,104(5)：054301

8 Christensen J,de Abajo FJG.Negative refraction and backward waves in layered acoustic metamaterials.Physical Review B,2012, 86(2)：024301

9 Zhu R,Liu XN,Hu GK,et al.Negative refraction of elastic waves at the deep-subwavelength scale in a single-phase metamaterial.Nature Communications,2014,24(5)：5510

10 Qiu CY,Zhang XD,Liu ZY.Far-fiel imaging of acoustic waves by a two-dimensional sonic crystal.Physical Review B,2005,71(5)：054302

11 劉宸，孫宏祥，袁壽其等.基于溫度梯度分布的寬頻帶聲聚焦效應(yīng).物理學(xué)報,2016,65(4)：044303(Liu Chen,Sun Hongxiang, Yuan Shouqi,et al.Broadband acoustic focusing ef f ect based on temperature gradient distribution.Acta Physica Sinica,2016,65(4)：044303(in Chinese))

12 Xia BZ,Dai HQ,Yu DJ.Symmetry-broken metamaterial for blocking,cloaking,and supertunneling of sound in a subwavelength scale.Applied Physics Letters,2016,108(25)：251902

13 Kushwaha MS,Halevi P,Dobrzynski L,et al.Acoustic band structure of periodic elastic composites.Physical Review Letters,1993,71(13)：2022

14 Liu Z,Zhang X,Mao Y,et al.Locally resonant sonic materials.Science,2000,289(5485)：1734-1736

15 Oudich M,Djafari-Rouhani B,Pennec Y,et al.Negative ef f ective mass density of acoustic metamaterial plate decorated with low frequency resonant pillars.Journal of Applied Physics,2014,116(18)：184504

17 Wang XL.Theory of resonant sound transmission through small apertures on periodically perforated slabs.Journal of Applied Physics,2010,108(6)：064903

18 溫激鴻，韓小云，王剛等.聲子晶體研究概述.功能材料,2003, 34(4)：364-367(Wen Jihong,Han Xiaoyun,Wang Gang,et al. Overview of phononic crystals.Journal of Functional Materials, 2003,34(4)：364-367(in Chinese))

19 Xia BZ,Chen N,Xie LX,et al.Temperature-controlled tunable acoustic metamaterial with active band gap and negative bulk modulus.Applied Acoustics,2016,112：1-9

20 夏百戰(zhàn)，覃緣，于德介等.區(qū)間模型下聲學(xué)超材料的可靠性優(yōu)化.機(jī)械工程學(xué)報,2016,13：094-102(Xia Baizhan,Qin Yuan,Yu Dejie,et al.Reliability-based optimization of the acoustic metamaterial under the interval model.Journal of Mechanical Engineering, 2016,13：094-102(in Chinese))

21 夏百戰(zhàn).不確定聲固耦合系統(tǒng)的數(shù)值分析與優(yōu)化方法研究.[博士論文].長沙：湖南大學(xué),2015(Xia Baizhan.A research on the numerical analysis and optimization method for the nondeterministic acoustic-structural coupled system.[PhD Thesis].Changsha：Hunan University,2015(in Chinese))

22 Jin RC,Wei C,Simpson TW.Comparison studies of metamodeling techniques under multiple modelling criteria.Structural and Multi-disciplinary Optimization,2001,23(1)：1-13

23 Xia BZ,Yu DJ.Modifiesub-interval perturbation finit element method for 2D acoustic fielprediction with large uncertain-butbounded parameters.Journal of Sound and Vibration,2012, 331(16)：3774-3790

24 Edgecombe S,Linse P.Monte Carlo simulation of two interpenetrating polymer networks：structure,swelling,and mechanical properties.Polymer,2008,49(7)：1981-1992

25 Wu JL,Luo Z,Zheng J,et al.Incremental modeling of a new highorder polynomial surrogate model.Applied Mathematical Modelling,2016,40(7)：4681-4699

26 沈祖和.區(qū)間分析方法及其應(yīng)用.應(yīng)用數(shù)學(xué)與計算數(shù)學(xué)學(xué)報，1983，2(1)：1-29(ShenZuhe.Intervalanalysismethodanditsapplication.Journal of Applied Mathematics and Computational Mathematics,1983，2(1)：1-29(in Chinese))

27 Muscolino G,SofA.Stochastic analysis of structures with uncertain-but-bounded parameters via improved interval analysis.Probabilistic Engineering Mechanics,2012,28：152-163

28 Xia B,Qin Y,Yu D,et al.Dynamic response analysis of structure under time-variant interval process model.Journal of Sound and Vibration,2016,381：121-138

29 Wu JL,Luo Z,Zhang YP,et al.An interval uncertain optimization method for vehicle suspensions using Chebyshev metamodels.Applied Mathematical Modelling,2014,38(5)：3706-3723

30 Yin SW,Yu DJ,Huang Y,et al.Hybrid chebyshev interval finite element and statistical energy analysis method for midfrequency analysis of built-Up systems with interval uncertainties.Journal of Engineering Mechanics,2016,142(10)：04016071

31 王艷鋒.含共振單元聲子晶體的帶隙特性及設(shè)計.[博士論文].北京：北京交通大學(xué),2015(Wang Yanfeng.Bandgap properties and design of phononic crystals with resonant units.[PhD Thesis].Beijing：Beijing Jiaotong University,2015(in Chinese))

THE INTERVAL ANALYSIS OF MULTILAYER-METAMATERIALS WITH PERFORATED APERTURES BASED ON CHEBYSHEV EXPANSION1)

Liu Jian Lei Jirong Xia Baizhan2)

(State Key Laboratory of Advanced Design and Manufacturing for Vehicle Body，Hunan University，Changsha410082，China)

Up to now,the response analysis and optimization of acoustic metamaterials are mostly based on deterministic parameters and deterministic models.However,in the real engineering world and structure design field where there are many uncertainties,such as the uncertain of material properties and geometric parameters.If the ef f ects of uncertain variables on analysis and optimization of acoustic metamaterials are not taken into account,the analysis and optimization results may not true.Aiming at this problem and situation,in this paper where the interval model is introduced into multilayer-metamaterials with perforated apertures,and interval Chebyshev expansion-Monte Carlo simulation method(ICE-MCSM)of multilayer-metamaterials with perforated apertures for transmission characteristic is proposed. In ICE-MCSM,truncated Chebyshev polynomials is applied to surrogate the transmittance of multilayer-metamaterials with perforated apertures;therefore,the surrogate model of transmittance is constructed.The samples of interval vari-ables are produced by Monte Carlo method,then the values of these samples are substituted into the surrogate model of transmittance to predict the interval bounds of multilayer-metamaterials with perforated apertures under single-interval variable and multi-interval variables.The results of numerical analysis show that the upper interval bounds and lower interval bounds calculated by ICE-MCSM match the response yields by direct Monte Carlo simulation method(DMCSM). Compared to DMCSM,ICE-MCSM achieves a higher accuracy in interval bound analysis,so ICE-MCSM can ef f ectively and efficiently analyze the interval bounds of multilayer-metamaterials with perforated apertures under uncertain interval variables.Thus,this proposed method in this paper can have promising prospects in real world engineering applications.

metamaterials,transmittance,Chebyshev expansion,interval analysis

O241

A doi：10.6052/0459-1879-16-254

2016-09-09收稿，2016-11-20錄用,2016-11-24網(wǎng)絡(luò)版發(fā)表.

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

2)夏百戰(zhàn)，副教授，主要研究方向：汽車振動和噪聲的不確定數(shù)值分析與優(yōu)化方法.E-mail：xiabz2013@hnu.edu.cn

劉堅(jiān)，雷濟(jì)榮，夏百戰(zhàn).基于Chebyshev展開的區(qū)間穿孔板超材料分析.力學(xué)學(xué)報,2017,49(1)：137-148

Liu Jian,Lei Jirong,Xia Baizhan.The interval analysis of multilayer-metamaterials with perforated apertures based on Chebyshev expansion.Chinese Journal of Theoretical and Applied Mechanics,2017,49(1)：137-148