一種非完全覆蓋多鍍層微圓板諧振器熱彈性阻尼通用計(jì)算框架

楊龍飛， 李 普， 葉一舟， 方昱斌

(1.南京理工大學(xué) 智能制造學(xué)院, 南京 210094； 2. 東南大學(xué) 機(jī)械工程學(xué)院, 南京 211189； 3. 重慶大學(xué) 光電工程學(xué)院, 重慶 400094)

隨著5 G物聯(lián)網(wǎng)時(shí)代飛速發(fā)展，體積小、精度高、耗能小的微機(jī)電系統(tǒng)(micro-electro-mechanical system, MEMS)呈現(xiàn)爆炸式增長(zhǎng)[1]，微機(jī)械諧振器是MEMS傳感器和執(zhí)行器的核心器件[2]。品質(zhì)因子(quality-factor, Q-factor)是衡量微諧振器精度和分辨率的關(guān)鍵參數(shù)之一，其倒數(shù)可理解為系統(tǒng)阻尼[3]。微機(jī)械諧振器運(yùn)行時(shí)主要存在內(nèi)部阻尼和外部阻尼兩種[4]。外部阻尼主要包括空氣阻尼和支撐阻尼，空氣阻尼可通過(guò)超低真空封裝消除；支撐阻尼可通過(guò)優(yōu)化錨區(qū)的結(jié)構(gòu)阻斷應(yīng)力波。但內(nèi)部阻尼只能降低，不能被完全消除，主要表現(xiàn)為熱彈性阻尼，并已被試驗(yàn)證實(shí)[5]。TED本質(zhì)是振動(dòng)過(guò)程中不均勻應(yīng)力場(chǎng)產(chǎn)生的溫度梯度，在弛豫過(guò)程中熵增造成的能量損失。通過(guò)研究微諧振器中TED的機(jī)理研究可為結(jié)構(gòu)優(yōu)化設(shè)計(jì)提供理論支撐，有助于提高其Q-factor值。因此，近幾年來(lái)研究人員和工程技術(shù)人員對(duì)TED的關(guān)注度一直很高。

Zener[6]于1937年首次提出了薄梁的TED解析模型，并給出了簡(jiǎn)單模型形式，被科研和工程設(shè)計(jì)領(lǐng)域廣為使用。2000年，Lifshitz等[7]采用新思路推導(dǎo)出基于復(fù)頻率法的熱彈性阻尼解析模型(以下簡(jiǎn)稱L-R模型)。Zener模型和L-R模型為后續(xù)研究各種微諧振器件TED奠定了基礎(chǔ)。2009年，Sun等[8]提出了軸對(duì)稱單層微圓板諧振器TED解析模型。2012年，Li等[9]提出了考慮z向應(yīng)變的微矩形板和圓板的TED模型，研究發(fā)現(xiàn)矩形板和圓板的TED表達(dá)式完全相同。2021年，Zhou等[10]提出了考慮非局部效應(yīng)和雙遲滯效應(yīng)微圓板TED模型。

隨著鍍膜工藝日趨完善，MEMS采用多層結(jié)構(gòu)滿足更多的功能性和多樣性需求，極大擴(kuò)展了其應(yīng)用場(chǎng)景。比如金屬膜經(jīng)常用于電極、質(zhì)量檢測(cè)器、光學(xué)反射、磁性單元和熱導(dǎo)體等。2018年，Liu等[11]提出了考慮雙層微圓板諧振器軸向和厚度方向熱傳導(dǎo)的TED借些模型。2020年，Li等[12]推導(dǎo)出功能梯度微圓板諧振器TED理論模型。

以上雙層微圓板TED模型均以鍍層與基底完全覆蓋為出發(fā)點(diǎn)進(jìn)行研究。實(shí)際應(yīng)用中，機(jī)械夾緊和電絕緣導(dǎo)致鍍層幾乎不可能完全與基底重合。Sandberg等[13]通過(guò)試驗(yàn)證實(shí)，即使很薄的金屬鍍層也會(huì)導(dǎo)致Q-factor大幅度下降，并建議在基底進(jìn)行部分區(qū)域鍍膜，以提高Q-factor。本文提出了一種通用理論框架，可計(jì)算具有多個(gè)非完全覆蓋鍍層微圓板諧振器TED。該模型可退化為完全覆蓋雙層微圓板TED模型[14]，并利用ANSYS驗(yàn)證了當(dāng)前模型的有效性。此外，本文提出了可用工程領(lǐng)域快速計(jì)算的簡(jiǎn)單模型，并分析簡(jiǎn)單模型的適用范圍。最后研究了具有2個(gè)非完全鍍層微圓板TED頻率譜特性。

1 理論框架

1.1 動(dòng)力學(xué)方程求解

以軸對(duì)稱周邊固支微圓板諧振器為研究對(duì)象，在圓板徑向布置有m(m≥1)個(gè)同心鍍層，幾何模型如圖1所示。假設(shè)基底層和鍍層是均質(zhì)、各向同性和線性熱彈性體。

圖1 多鍍層微圓板諧振器結(jié)構(gòu)示意圖Fig.1 Schematic diagram of a circular microplate resonator partially covered by n coatings

鍍層改變了基底抗彎剛度，覆蓋區(qū)域變?yōu)殡p層結(jié)構(gòu)，未覆蓋區(qū)域仍為單層結(jié)構(gòu)。將結(jié)構(gòu)整體劃分為2m+1個(gè)子區(qū)域R{p}，上標(biāo){p}為各子區(qū)域索引，p=1，3，5，…，2m+1代表未覆蓋區(qū)，p= 2，4，6，…，2m表示覆蓋區(qū)。覆蓋區(qū)中基底層厚度為hs，半徑為rs。圓柱坐標(biāo)系(r，φ，z)固定在圓板圓心O處。下標(biāo)k表示第k個(gè)鍍層，k=1，2，3，…，m。第k個(gè)鍍層內(nèi)、外半徑分別為r{p=2k-1}和r{p=2k}，厚度為hk。子區(qū)域r{2k-2}≤r

子區(qū)域R{p}簡(jiǎn)諧線性軸對(duì)稱振動(dòng)的振型函數(shù)設(shè)為

(1)

Kirchhoff-Love板的轉(zhuǎn)動(dòng)慣量和剪切力對(duì)控制方程的影響可忽略，自由振動(dòng)控制方程可設(shè)為[15]

(2)

式中：ω=2πf為角頻率；D{p}為子區(qū)域等效抗彎剛度。

(3)

(4)

w{p}(r)=A{p}J0(γ{p}r)+B{p}Y0(γ{p}r)+
C{p}I0(γ{p}r)+M{p}K0(γ{p}r)

(5)

式中：J0和Y0為零階第一類(lèi)和第二類(lèi)Bessel函數(shù)；I0和K0為零階第一類(lèi)和第二類(lèi)修正Bessel函數(shù)。

對(duì)于未覆蓋區(qū)R{p=1}，為確保圓心位置位移函數(shù)是有限的，系數(shù)B{1}和M{1}應(yīng)為0，其解簡(jiǎn)化為

w{1}(r)=A{1}J0(γ{1}r)+C{1}I0(γ{1}r)

(6)

單層區(qū)與雙層區(qū)連接位置的結(jié)構(gòu)連續(xù)性需要滿足以下邊界條件[17]，主要包括：

(1) 位移連續(xù)性

(7)

(2) 坡度連續(xù)性

(8)

(3) 彎矩連續(xù)性

(9)

(4) 剪力連續(xù)性

(10)

此外，基底層外邊緣被夾緊，表示為

(11)

將邊界條件式(7)～式(11)分別代入振型函數(shù)式(5)和式(6)中，可得8m+2個(gè)齊次方程，與振型函數(shù)式(5)中包含得8m+2個(gè)未知參數(shù)A{p}，B{p}，C{p}和M{p}一一對(duì)應(yīng)。將該8m+2個(gè)方程組寫(xiě)為矩陣形式

M(γ{p})·u=0

(12)

為了獲得非平凡解，系數(shù)矩陣M(γ{p})的行列式必須為0

|M(γ{p})|=0

(13)

該式為高度非線性超越方程。

固有頻率ω和γ{p}之間的關(guān)系

(14)

經(jīng)分析，固有頻率ω是矩陣方程式(13)中唯一未知量，可通過(guò)迭代尋根法求解[18]。將得到的ω代入方程(12)求解得到模態(tài)振型參數(shù)u。

軸對(duì)稱振動(dòng)模式圓盤(pán)的位移分量為[19]

(15)

線性彈性理論下極坐標(biāo)系中應(yīng)變分量為

(16)

根據(jù)胡克定律，未覆蓋區(qū)應(yīng)力分量為

(17)

覆蓋區(qū)應(yīng)力分量表示為

(18)

1.2 溫度場(chǎng)求解

熱彈性耦合[20]導(dǎo)致圓板在簡(jiǎn)諧振動(dòng)時(shí)產(chǎn)生簡(jiǎn)諧溫度場(chǎng)，圓板中性面不受壓縮或拉伸，保持在平衡溫度T0?；讓雍湾儗拥臏囟葓?chǎng)增量可寫(xiě)為

(19)

根據(jù)薄板熱彈性耦合理論，沿z軸的應(yīng)變?yōu)?/p>

(20)

由于體積應(yīng)變是產(chǎn)生溫度梯度的根源，因此圓板的總體積應(yīng)變?yōu)?/p>

(21)

熱彈性圓板中的熱傳導(dǎo)方程為[21]

(22)

(23)

微諧振器采用真空超低真空封裝可減少氣體阻尼，其上下表面處于絕熱態(tài)，熱邊界條件表示為

(24)

此外，磁控濺射鍍膜工藝可保證基底層與鍍層原子級(jí)別的結(jié)合，可假設(shè)界面處熱傳導(dǎo)是理想態(tài)，即溫度和熱通量連續(xù)

(25)

熱傳導(dǎo)方程式(22)本質(zhì)上屬于非齊次偏微分方程。為求解此方程，首先應(yīng)求解其齊次形式

(26)

式中，α=κ/CV為材料的熱擴(kuò)散系數(shù)。

對(duì)于雙層區(qū)域R{2}，R{4}…，R{2m}，采用分離變量方法求解，式(26)解可設(shè)為[22]

(27)

將式(27)代入式(26)得

(28)

根據(jù)熱邊界條件式(24)，可設(shè)

(29)

熱邊界條件式(24)和式(25)代入式(29)可得

(30)

(31)

為求解式(30)，其行列式應(yīng)為0

(32)

式(32)化簡(jiǎn)后可得超越方程

(33)

(34)

對(duì)于未覆蓋區(qū)域R{1}，R{3}…，R{2m+1}，其本質(zhì)是均質(zhì)單層板。熱傳導(dǎo)控制方程式(22)退化為

(35)

根據(jù)熱邊界條件式(24), 式(35)的解設(shè)為

(36)

(37)

時(shí)間變量函數(shù)Γ[j](t)的解為

(38)

(39)

至此，熱傳導(dǎo)的特征函數(shù)和特征值全部求解完畢。

(40)

式中，格林函數(shù)為

(41)

(42)

未覆蓋區(qū)域溫度場(chǎng)函數(shù)式(19)的解為

(43)

式中，格林函數(shù)為

(44)

(45)

將格林函數(shù)式(41)和式(44)和熱源項(xiàng)式(23)代入式(40)和式(44)，溫度場(chǎng)函數(shù)化簡(jiǎn)后可得

(46)

(47)

(48)

溫度場(chǎng)函數(shù)式(46)的虛部為

(49)

1.3 熱彈性阻尼模型

根據(jù)Bishop等[23]之前提出的計(jì)算多層板熱彈性阻尼的框架，針對(duì)不規(guī)則復(fù)雜結(jié)構(gòu)，本文提出將將整體結(jié)構(gòu)劃分為若干規(guī)則子區(qū)域的新思路。依據(jù)能量法，熱彈性阻尼可由式(50)求得

(50)

式中：ΔW{p}為一個(gè)周期的子區(qū)域能量損失；Wmax為存儲(chǔ)的最大應(yīng)變能。子區(qū)域耗散項(xiàng)ΔW{p}可通過(guò)式(51)計(jì)算

(51)

將應(yīng)力本構(gòu)式(17)、式(18)和溫度場(chǎng)虛部式(49)代入式(51)，化簡(jiǎn)后得到

(52)

式中，

(53)

(54)

(55)

一個(gè)振動(dòng)周期內(nèi)總能量損失為

(56)

在一個(gè)振動(dòng)周期中，圓板中存儲(chǔ)的最大勢(shì)能為

(57)

式中，

(58)

本文所有積分形式的結(jié)果均是顯式的。根據(jù)熱彈性阻尼的定義，利用式(56)和式(57)得到熱彈性阻尼解析解為

(59)

式中，參數(shù)ψ{p}

(60)

以及熱弛豫時(shí)間常數(shù)為

(61)

只保留第一項(xiàng)(n=1)可得熱彈性阻尼簡(jiǎn)單模型

(62)

式中，熱弛豫時(shí)間常數(shù)為

(63)

表1列出了非完全覆蓋多鍍層微圓板諧振器TED的完整模型和簡(jiǎn)單模型表達(dá)式。

表1 非完全覆蓋多鍍層微圓板諧振器熱彈性阻尼完整模型和簡(jiǎn)單模型表達(dá)式Tab.1 Comparison of the full model, the simple model of TED

2 模型分析與討論

為便于分析，將結(jié)構(gòu)參數(shù)進(jìn)行無(wú)參化設(shè)置：

fhk=hk/hs——鍍層k與基層s厚度比；

fpk=r[2k-1]/rs——鍍層k內(nèi)半徑r[2k-1]與基層半徑rs位置比；

frk=(r[2k]-r[2k-1])/rs——鍍層k外半徑r[2k]與內(nèi)半徑r[2k-1]之差與基層半徑rs半徑的比；

fαk=αk/αs——fαk代表鍍層k熱擴(kuò)散率αk與基層s熱擴(kuò)散率αs之比。

本文設(shè)定基底圓板幾何參數(shù)保持不變，半徑rs=500 μm，厚度hs=10 μm。表2給出了基底層和鍍層物理參數(shù)。

表2 材料物理參數(shù)Tab.2 Material physical coefficients

2.1 模型驗(yàn)證與收斂速度分析

可從兩方面驗(yàn)證本模型有效性。首先，完全覆蓋雙層圓板是非完全覆蓋雙層圓板的特殊情況。以覆蓋有一個(gè)鍍膜的微圓板諧振器為例，結(jié)構(gòu)示意如圖2所示。假設(shè)r{1}=0和r{2}=rs，當(dāng)前TED模型可完全退化為完全覆蓋雙層微圓板諧振器模型。

其次，可利用商業(yè)有限元軟件ANSYS進(jìn)行驗(yàn)證。為減少計(jì)算量，采用8節(jié)點(diǎn)二維軸對(duì)稱單元plane223。進(jìn)行網(wǎng)格無(wú)關(guān)性驗(yàn)證后，將網(wǎng)格尺寸設(shè)置為小于0.5 μm。沿圓板厚度方向施加均布力F0sin(ωt)，F(xiàn)0= 0.1 MPa。圖2給出了微圓板諧振器的有限元模型、諧響應(yīng)分析后的振型與溫度場(chǎng)云圖。

圖2 有限元模型、諧響應(yīng)分析后的振型與溫度場(chǎng)云圖Fig.2 FEM,mode shape and temperature field by the harmonic response analysis of circular microplate resonators

圖3給出了位置比f(wàn)p1=0.2，半徑比f(wàn)r1=0.6，厚度比f(wàn)h1= 0.1時(shí)Ag/Si微圓板諧振器TED頻率譜。同樣給出了有限元結(jié)果作為對(duì)比，驗(yàn)證了本模型的有效性。此外，在低于Debye峰對(duì)應(yīng)的臨界頻率，即f≤1.1×106Hz，簡(jiǎn)單模型(n=1)與n=15時(shí)的TED值誤差小于5%。但隨著頻率變大，其誤差逐漸增大，簡(jiǎn)單模型不再適用。n=10與n=15之間的誤差小于5%，因此完整模型取n=15即可保證足夠的精度。

圖3 位置比f(wàn)p1=0.2、半徑比f(wàn)r1=0.6、厚度比f(wàn)h1=0.1時(shí)Ag/Si微圓板諧振器TED頻率譜Fig.3 TED frequency spectra of Ag/Si circular microplate resonators in the case of position ratio fp1=0.2, radius ratio fr1=0.6, thickness ratio fh1 = 0.1

圖4分析了fr1=1和fh1=1時(shí)SiO2/Si微圓板中TED頻率譜的收斂性。fr1=1表示基底層Si圓板被鍍層SiO2完全覆蓋。從圖4中可看出，頻率譜中可觀察到兩個(gè)顯著的Debye峰。若使用簡(jiǎn)單模型(n=1)計(jì)算則高頻率處Debye峰將被忽略，存在嚴(yán)重誤差。

圖4 半徑比f(wàn)r1 = 1、厚度比f(wàn)h1=1時(shí)SiO2/Si微圓板諧振器TED頻率譜Fig.4 TED frequency spectra of SiO2/Si circular microplate resonators in the case of radius ratio fr1=1, thickness ratio fh1=1

2.2 簡(jiǎn)單模型和完整模型的對(duì)比分析

本節(jié)從鍍層的厚度、半徑出發(fā)分析簡(jiǎn)單模型的適用范圍。圖5給出了不同厚度比f(wàn)h1=h1/hs的Ag/Si圓板TED頻率譜圖。從圖5可看出，隨著fh1減小，模型收斂性逐漸增大。當(dāng)fh1=0.01時(shí)，TED簡(jiǎn)單模型與完整模型的值相對(duì)誤差小于2%。此時(shí)，可用簡(jiǎn)單模型代替完整模型進(jìn)行快速計(jì)算，以降低計(jì)算量。

圖5 不同厚度比f(wàn)h1=h1/hs的Ag/Si圓板TED頻率譜Fig.5 TED frequency spectra of Ag/Si microplates with different thickness ratio fh1=h1/hs

圖6分析了不同厚度比f(wàn)h1=h1/hs的SiO2/Si圓板TED頻率譜特性?？勺⒁獾?，只有fh1=1時(shí)出現(xiàn)了雙峰現(xiàn)象，隨著fh1減小，低頻峰逐漸消失，但高頻峰的峰值逐漸增大?？芍?，SiO2層可降低高頻峰值，但在低頻處會(huì)增大TED值。

圖6 不同厚度比f(wàn)h1=h1/hs的SiO2/Si圓板TED頻率譜Fig.6 TED frequency spectra of SiO2/Si microplates with different thickness ratio fh1=h1/hs

圖7給出了不同半徑比f(wàn)r1=(r[2k]-r[2k-1])/rs的Ag/Si圓板中TED頻率譜值。從圖7可看出，隨著fr1減小，頻率譜峰值逐漸減小，收斂性也逐漸變快。當(dāng)fr1=0.2時(shí)，可用簡(jiǎn)單模型代替完整模型快速估算TED值。

圖7 不同半徑比f(wàn)r1=(r[2k]-r[2k-1])/rs的Ag/Si圓板中TED頻率譜Fig.7 TED frequency spectra of Ag/Si microplates with different radii ratio fr1=(r[2k]-r[2k-1])/rs

圖8研究了不同半徑比(r[2k]-r[2k-1])/rs的SiO2/Si圓板中TED頻率譜特性。從圖8可知，隨著fr1逐漸變小，低頻峰值逐漸減小，但高頻峰值逐漸增大。直到fr1=0.2時(shí)，低頻峰消失。同圖7一樣，此時(shí)用簡(jiǎn)單模型計(jì)算TED值即可保證較高精度。

圖8 不同半徑比(r[2k]-r[2k-1])/rs的SiO2/Si圓板中TED頻率譜Fig.8 TED frequency spectra of SiO2/Si microplates with different radii ratio (r[2k]-r[2k-1])/rs

2.3 具有2個(gè)局部鍍層的微圓板TED頻率譜分析

圖9分析了Si/DLC/SiO2圓板TED頻率譜特性。DLC (diamond-like carbon)是類(lèi)金剛石典型特征的一類(lèi)材料的統(tǒng)稱。鍍層1(Si)、鍍層2(DLC)與基底層(SiO2)的厚度比分別為fh1=1,fh2=0.4；位置比分別為fp1=0,fp2=0.8；半徑比為fr1=0.8,fr2=0.2。

如圖9所示，TED頻率譜中存在3個(gè)顯著的Debye峰。如果使用簡(jiǎn)單模型計(jì)算，將丟失第2和第3個(gè)峰值，造成嚴(yán)重誤差。因此，此時(shí)不可使用簡(jiǎn)單模型計(jì)算TED。此外，當(dāng)n=20計(jì)算時(shí)，可觀察到2個(gè)Debye峰，第3個(gè)Debye峰丟失。只有使用完整模型取n=100時(shí)，才可捕捉到3個(gè)Debye峰。n=100時(shí)TED值與n=400時(shí)的結(jié)果幾乎完全重合。

圖9 Si/DLC/SiO2圓板的TED頻率譜特性曲線。厚度比為fh1=1,fh2=0.4；位置比為fp1=0,fp2=0.8；半徑比為fr1=0.8,fr2=0.2Fig.9 TED frequency spectra behaviors of Si/DLC/SiO2 in circular microplates with the thickness ratio fh1=1, fh2=0.4, position ratio fp1=0, fp2=0.8 and radii ratio fr1=0.8, fr2=0.2

3 結(jié) 論

本文首次提出一種可計(jì)算具有多個(gè)非完全覆蓋鍍層微圓板諧振器TED的通用計(jì)算框架。根據(jù)分析結(jié)果，概括出以下幾個(gè)重要結(jié)論：

(1) 該解析模型僅用于微薄板諧振器，可退化為完全覆蓋雙層圓板模型，與有限元計(jì)算結(jié)果吻合很好。

(2) 該解析模型收斂性隨厚度比和半徑比減小而逐漸增大。在fh1≤0.01或fr1≤0.2時(shí)可用簡(jiǎn)單模型替代完整模型計(jì)算TED值，以減少計(jì)算量。

(3) 在Si/DLC/SiO2微圓板的TED頻率譜中發(fā)現(xiàn)3個(gè)顯著的Debye峰。此時(shí)，只能使用完整模型計(jì)算，使用簡(jiǎn)單模型會(huì)造成嚴(yán)重誤差。