基于最小二乘數據擬合的微小壓力傳感器仿真設計?

宮凱勛雷 程?董志超趙珠杰熊繼軍齊 蕾

(1.中北大學儀器科學與動態測試教育部重點實驗室，山西 太原030051；2.北方自動控制技術研究所，山西 太原030051)

在MEMS器件的設計加工中，通常具有周期長、費用高等特點，可靠的仿真計算在器件的設計中十分重要[1]。硅壓阻式壓力傳感器因其靈敏度高、精度高、體積小、動態頻響高等優點一直作為MEMS壓力傳感器領域的熱門應用與研究對象[2]。壓阻式壓力傳感器敏感膜上的力學敏感器件受壓后阻值的變化決定著傳感器的輸出[3－5]，而阻值的變化與應力有關。普遍的對于仿真計算得到的壓敏電阻上應力的處理是取其中心一點的應力值進行估算或者通過定義路徑計算路徑上壓敏電阻的應力積分，再除以壓敏電阻的長度求得壓敏電阻所在位置的平均應力。

上述算法都有其局限性，忽視了壓敏電阻整體的應力分布不均勻的性質。在對壓阻式壓力傳感器的仿真設計中通常將承壓膜看作薄板結構[6]，施加垂直于板面的應力時，其應力的分布均以中面為對稱面，且在中面處降至最低。壓力傳感器特別是微壓傳感器因其膜厚與壓敏電阻的高度比值較小，不可忽視壓敏電阻延高度方向上應力的變化。為解決上述問題，本文采用基于最小二乘法的多項式擬合方法對壓敏電阻橫截面方向上的應力進行函數擬合，求得在固定z軸數值的平均應力，通過對不同z軸坐標下壓敏電阻平均應力的計算，擬合出z軸對應力的函數關系，再計算壓敏電阻的整體平均應力，為壓阻式壓力傳感器提供更科學的設計計算方法。

1 算法原理

1.1 彈性薄板理論

壓阻式壓力傳感器上敏感膜作為彈性感壓部件，一般情況下，敏感膜在理論分析時被看作是在厚板與薄膜之間的薄板結構。根據彈性力學理論，其在分析時滿足小撓度理論的模型，彈性薄板小撓度理論的基本假設:板彎曲時，板中面內各點只存在垂直方向上的位移沒有平行中面方向上的位移；彎曲變形前垂直于薄板中面的直線，變形后仍垂直于彈性曲面且長度不變；薄板各層的纖維在變形前后均互不擠壓[7]。

基于這些假設，對于薄板問題的求解最重要的是對于薄板撓曲函數ω(x，y)的求解，圖一所示為在方形薄板施加面載荷。

在薄板上建立直角坐標系，其撓度變化可表示為:

式中:D為彎曲剛度，可用式(2)表示；p為薄板所受的壓力載荷。

式中:E為材料的彈性模量；v為泊松比；t為薄板厚度。

薄板的邊界條件為x＝0，x＝a，y＝0，y＝a，即:

代入式(1)可解得薄板的撓度為:

薄板上任一點(x，y)延x方向應力δx與延y方向的應力δy為:

由式(5)可知，主要的應力分量延z軸線性分布，故可先對(x，y)坐標與應力值擬合，再對不同Z坐標下的平均應力值做擬合積分得到更為準確的仿真壓敏電阻平均應力計算值。

1.2 最小二乘法擬合

數據擬合的目的是利用最小二乘方法，將基于模型函數所計算出的數值與觀測量之間的殘差平方和最小化[8－11]。

以本文采用的三元二次多項式方程擬合為例，對于方程:

有N組具有相同精度的數據(xi，yi)，其中i＝1，…，N；xi，yi都是精確的，對與之相對應的zi值進行預測。用最小二乘法估計參數時，要求觀測值zi的偏差加權平方和最小。用泛函誤差δ表示:

a0，a1，a2，a3，a4，a5的值影響著δ的大小，求取δ的最小值，即對a0，a1，a2，a3，a4，a5分別進行求偏導，轉化為一個泛函求極值問題。

2 有限元仿真

COMSOL Multiphysics是一款功能強大的多物理場仿真軟件，用于仿真模擬工程、制造和科研等各個領域的設計、設備及過程。其界面簡潔，模型建立簡單且仿真計算后的數據導入與導出都極為便捷。本文采用comsol軟件對壓敏電阻應力值進行仿真計算。

2.1 建模與網格劃分

壓阻式壓力傳感器的敏感膜一般為方型膜與圓型膜，方形薄板的最大應力比圓形薄板的最大應力大，有利于獲得更高靈敏度的壓力傳感器，所設計的微壓傳感器采用方形膜作為感壓膜面。

在comosl軟件中，一般的通過布爾操作對幾何體進行結構建模的方法對于計算后網格數據提取與處理較為困難。為解決這一問題，采用各結構獨立建模，最后通過構成聯合體的方式形成完整的器件結構，為網格劃分與數據的后處理提供了便利。

敏感膜模型尺寸設計為460μm×460μm×7μm，壓敏電阻尺寸為65μm×13μm×1.5μm。壓敏電阻位置在膜邊緣10μm處。單晶硅材料屬性包括楊氏模量170 GPa、泊松比0.3、密度為2 329 kg/m3。傳感器采用惠斯通半橋結構，建模時對垂直方向與平行方向的一個壓敏電阻進行網格加密剖分。

各部分單獨建模的好處是在劃分網格的過程中可較好地控制局部網格的密度，網格劃分采用用戶控制網格，將壓敏電阻網格延Z軸方向劃分10層，X、Y軸方向自由剖分，得到延Z軸分層結構的四面體網格，結構剩余部分采用自由四面體網格自適應劃分，圖2為網格刨分結果圖。應力從敏感膜面到與支撐體相連邊界處后迅速衰減，將實際模型簡化為下圖結構以減少不必要的計算。邊界條件:敏感膜上表面施加40 kPa的壓力，結構底部采用固定約束方式。

圖2 網格刨分圖

2.2 有限元分析結果

從圖3可以看出，垂直膜面位置壓敏電阻應力分布結果與彈性力學分析結果一致，應力延中面對稱分布，且中面位置應力降到最低。在距離膜面1.5μm位置上的壓敏電阻應力下降為膜表面處應力的一半，應力延z軸方向上的變化已不能忽略。

圖3 沿膜厚方向von Mises應力分布

從圖4垂直膜邊壓敏電阻表面應力等值線圖可以得出，壓敏電阻表面不同位置應力大小也不同，故在壓力傳感器設計中必須對壓敏電阻整體結構的應力分布進行考慮。

圖4 壓敏電阻縱向應力橫向應力差值分布圖

3 多項式數據擬合處理

常用到的數據擬合主要有線性回歸擬合與多項式擬合，線性回歸擬合要求自變量與因變量之間具有線性關系，多項式擬合較線性回歸擬合總體上更為靈活，可以模擬一些相當復雜的關系[12]。在數據范圍內，多項式擬合方法可以提供很好的擬合效果，但在數據范圍外可能會大范圍發散。壓敏電阻不存在對數據范圍以外進行預測，其上各點位置的應力值都可以通過仿真計算得到，使用MATLAB軟件通過基于最小二乘法的多項式擬合方式來對壓敏電阻平均應力進行計算。

在COMSOL軟件的后處理部分，通過域選擇導出壓敏電阻部分的計算數據集，并對不同Z坐標的數據進行分類，將處理好的數據通過importdata函數導入到MATLAB中。

對垂直于膜邊壓敏電阻應力進行處理，在Z坐標為0即壓敏電阻頂面位置時，圖5為采用二次多項式擬合方式得到的結果，其中離散點為COMSOL計算得到的原始數據，曲面為擬合得到的結果。

圖5 x，y坐標值與應力值擬合

擬合得到的多項式函數為:

擬合確定系數R-square:0.999 8，接近于1，表明該方程可以很好的對散點數據進行表達。

將原始數據坐標代入到擬合得到的函數進行計算，通過MATLAB里plot3函數作圖與原始數據進行對比，從圖6來看，原始數據與擬合數據在同一曲面上，可以忽略由離擬合得到的函數的計算值與仿真得到的數據誤差。

圖6 擬合計算值與仿真計算值對比

在MATLAB語言中，求符號函數的定積分是使用int函數，其調用格式如下:int(F，x，a，b)。a表示定積分的下限；b表示定積分的上限；使用int函數對表達式進行積分，再除以電阻條截面積可得到固定z坐標下壓敏電阻截面上的平均應力，表1為計算結果。

表1 施加40 kPa壓力不同高度下壓敏電阻截面的平均應力

Z軸坐標與垂直膜邊壓敏電阻計算出對應的平均應力進行擬合，圖7為擬合得出的曲線。

圖7 Z坐標值與應力值擬合

擬合得到的函數為:

擬合確定系數為R-square:1。

對擬合得到的函數積分可求得在壓力為40 kPa時，垂直膜邊壓敏電阻平均應力值為12 261 kPa。重復上述計算，得到平行膜邊壓敏電阻平均應力值為19 787 kPa。

4 傳感器輸出計算

基于上述方法，計算得到不同壓力載荷下各壓敏電阻的平均應力值，表2為計算結果。

表2 不同壓力載荷下各壓敏電阻的平均應力值 單位:kPa

對于采用半橋方式實現壓敏電阻內部連接的壓力傳感器，其電壓輸出公式為:

圖8所示為在輸入電壓為3.3 V時靈敏度的仿真計算結果。

圖8 分層截面積分法壓力與輸出電壓曲線

擬合得到的函數為:

確定系數R-square:1

得到在電壓為3.3 V時，所設計傳感器靈敏度為0.4491 mv/KPa。

若采用中心點算法，提取電阻條表面中心點位置應力值，在40 kPa下，σa＝25 543 kPa，σb＝14 671 kPa，根據式(14)計算得輸出電壓為22.5 mV，圖9中虛線為中心點法壓力與輸出電壓曲線，得到采用中心點法設計的傳感器靈敏度為0.562 8 mV/KPa。采用路徑積分法，提取電阻條表面中心線處應力做積分，計算得到在40 kPa下，σa＝25 298 kPa，σb＝15 775 kPa，根據式(14)計算得輸出電壓為23.0 mV。為與中心點法計算輸出曲線區分，將其輸出曲線零點變為10 mV，圖9中實線為路徑積分法壓力與輸出電壓曲線，得到采用路徑積分法設計的傳感器靈敏度為0.575 mV/KPa。

圖9 中心點法與路徑積分法壓力與輸出電壓曲線

5 芯片加工與測試

芯片選用頂層硅6.5μm、埋氧層厚度0.5μm、片厚300μm的SOI作為襯底；采用離子注入方式摻雜形成壓敏電阻；擴散爐內溫度1 000℃、時間20 min重摻雜構建歐姆接觸區；磁控濺射鈦、鉑、金并在溫度500℃、時間5 min高真空下退火形成歐姆接觸；通過刻蝕至SOI埋氧層來精準控制敏感膜厚度，采用陽極鍵合工藝溫度380℃、壓力1 200 N、電壓800 V下實現硅－玻璃鍵合形成表壓參考腔。圖10所示為芯片正面結構顯微鏡觀察。

圖1 方形板直角坐標系

圖10 芯片正面結構

壓力傳感器芯片封裝在特制的不銹鋼測試管座，安裝在常溫壓力測試平臺氣壓泵接口，使用高精度萬用表測量記錄實驗數據，圖11所示為實驗測試結果，靈敏度為0.455 2 mv/kPa，與設計計算結果偏差為1.34%，采用中心點法與路徑積分法的設計偏差分別為23.64%與26.32%。采用分層截面積分方法可極大提高仿真設計精度。

圖11 不同壓力下電壓輸出結果

6 總結

為減小壓阻式壓力傳感器特別是微壓傳感器設計結果與實際制作測試的誤差，本文提出一種分層截面積分的方法。通過對壓敏電阻與敏感膜面分別建模實現對局部網格的精確加密與局部數據的導出，基于最小二乘法的多項式擬合方法實現對壓敏電阻離散點應力值的擬合計算，將傳統的中心點算法(零維)、路徑積分法(一維)提高到了分層截面積分(三維)的程度，提高了仿真計算結果的可信度，所設計的傳感器靈敏度與實際測得靈敏度偏差為1.34%。該方法為壓阻式壓力傳感器特別是微小壓力傳感器的靈敏度計算制作提供了更為可靠計算方式，在壓阻式壓力傳感器設計中具有現實意義。