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

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

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

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

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

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

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

1 理論框架

1.1 動力學方程求解

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

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

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

子區域R{p}簡諧線性軸對稱振動的振型函數設為

(1)

Kirchhoff-Love板的轉動慣量和剪切力對控制方程的影響可忽略，自由振動控制方程可設為[15]

(2)

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

(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為零階第一類和第二類Bessel函數；I0和K0為零階第一類和第二類修正Bessel函數。

對于未覆蓋區R{p=1}，為確保圓心位置位移函數是有限的，系數B{1}和M{1}應為0，其解簡化為

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

(6)

單層區與雙層區連接位置的結構連續性需要滿足以下邊界條件[17]，主要包括：

(1) 位移連續性

(7)

(2) 坡度連續性

(8)

(3) 彎矩連續性

(9)

(4) 剪力連續性

(10)

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

(11)

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

M(γ{p})·u=0

(12)

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

|M(γ{p})|=0

(13)

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

固有頻率ω和γ{p}之間的關系

(14)

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

軸對稱振動模式圓盤的位移分量為[19]

(15)

線性彈性理論下極坐標系中應變分量為

(16)

根據胡克定律，未覆蓋區應力分量為

(17)

覆蓋區應力分量表示為

(18)

1.2 溫度場求解

熱彈性耦合[20]導致圓板在簡諧振動時產生簡諧溫度場，圓板中性面不受壓縮或拉伸，保持在平衡溫度T0。基底層和鍍層的溫度場增量可寫為

(19)

根據薄板熱彈性耦合理論，沿z軸的應變為

(20)

由于體積應變是產生溫度梯度的根源，因此圓板的總體積應變為

(21)

熱彈性圓板中的熱傳導方程為[21]

(22)

(23)

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

(24)

此外，磁控濺射鍍膜工藝可保證基底層與鍍層原子級別的結合，可假設界面處熱傳導是理想態，即溫度和熱通量連續

(25)

熱傳導方程式(22)本質上屬于非齊次偏微分方程。為求解此方程，首先應求解其齊次形式

(26)

式中，α=κ/CV為材料的熱擴散系數。

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

(27)

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

(28)

根據熱邊界條件式(24)，可設

(29)

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

(30)

(31)

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

(32)

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

(33)

(34)

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

(35)

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

(36)

(37)

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

(38)

(39)

至此，熱傳導的特征函數和特征值全部求解完畢。

(40)

式中，格林函數為

(41)

(42)

未覆蓋區域溫度場函數式(19)的解為

(43)

式中，格林函數為

(44)

(45)

將格林函數式(41)和式(44)和熱源項式(23)代入式(40)和式(44)，溫度場函數化簡后可得

(46)

(47)

(48)

溫度場函數式(46)的虛部為

(49)

1.3 熱彈性阻尼模型

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

(50)

式中：ΔW{p}為一個周期的子區域能量損失；Wmax為存儲的最大應變能。子區域耗散項ΔW{p}可通過式(51)計算

(51)

將應力本構式(17)、式(18)和溫度場虛部式(49)代入式(51)，化簡后得到

(52)

式中，

(53)

(54)

(55)

一個振動周期內總能量損失為

(56)

在一個振動周期中，圓板中存儲的最大勢能為

(57)

式中，

(58)

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

(59)

式中，參數ψ{p}

(60)

以及熱弛豫時間常數為

(61)

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

(62)

式中，熱弛豫時間常數為

(63)

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

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

2 模型分析與討論

為便于分析，將結構參數進行無參化設置：

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

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

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

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

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

表2 材料物理參數Tab.2 Material physical coefficients

2.1 模型驗證與收斂速度分析

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

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

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

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

圖3 位置比fp1=0.2、半徑比fr1=0.6、厚度比fh1=0.1時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時SiO2/Si微圓板中TED頻率譜的收斂性。fr1=1表示基底層Si圓板被鍍層SiO2完全覆蓋。從圖4中可看出，頻率譜中可觀察到兩個顯著的Debye峰。若使用簡單模型(n=1)計算則高頻率處Debye峰將被忽略，存在嚴重誤差。

圖4 半徑比fr1 = 1、厚度比fh1=1時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 簡單模型和完整模型的對比分析

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

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

圖6分析了不同厚度比fh1=h1/hs的SiO2/Si圓板TED頻率譜特性。可注意到，只有fh1=1時出現了雙峰現象，隨著fh1減小，低頻峰逐漸消失，但高頻峰的峰值逐漸增大。可知，SiO2層可降低高頻峰值，但在低頻處會增大TED值。

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

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

圖7 不同半徑比fr1=(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時，低頻峰消失。同圖7一樣，此時用簡單模型計算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個局部鍍層的微圓板TED頻率譜分析

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

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

圖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 結 論

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

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

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

(3) 在Si/DLC/SiO2微圓板的TED頻率譜中發現3個顯著的Debye峰。此時，只能使用完整模型計算，使用簡單模型會造成嚴重誤差。