振動體聲學靈敏度分析的Burton-Miller奇異邊界法及其MATLAB工具箱開發

張汝毅， 王發杰， 程隋福， 劉建政

(青島大學 機電工程學院， 山東 青島 266071)

0 引 言

在產品設計和開發過程中，產品的NVH性能是衡量產品競爭力的一項重要指標。聲學靈敏度分析是產品NVH性能分析的重要組成部分，可以為產品的優化設計提供方向和依據，降低開發成本。因此，發展聲學靈敏度分析[1]的準確高效數值算法和開發相關計算軟件平臺具有重要的理論價值和實際意義。

本文基于Burton-Miller奇異邊界法(BM-SBM)[2-3]，采用MATLAB軟件實現聲學靈敏度分析的相關代碼生成，并對聲學靈敏度相關的幾個標準算例進行分析測試。BM-SBM引入Burton-Miller公式，采用物理力學問題控制方程的奇異基本解作為插值基函數，可有效解決聲學問題計算中出現的虛假特征頻率問題，是一種半解析邊界型無網格方法。該方法可簡化傳統有限元法[4]中的網格劃分等前處理過程，也無須傳統邊界元法[5]中的奇異計算，具有數學理論簡單、計算精度高、易于數值實現等特點。BM-SBM通過源點強度因子避免基本解的奇異性[6-7]，因此如何有效計算源點強度因子尤為重要。目前，許多學者針對該問題提出不同的解決方案，主要包括反插值技術[8]、加減去奇異技術[9]、經驗公式[2,10]等。其中，經驗公式是最簡單直接的一種方法，具有節省計算成本、提高計算效率的優勢。因此，本文采用經驗公式計算源點強度因子。

MATLAB具有強大的計算能力，是功能種類豐富的商業軟件，在數據分析、算法開發、信息可視化、人機交互等方面應用廣泛。MATLAB還內置各類擴展工具箱，如Stateflow、RF Toolbox等，可以使軟件應用于不同領域。[11]本文利用MATLAB軟件實現BM-SBM聲學靈敏度分析代碼生成，開發相應的圖形用戶界面(graphical user interface, GUI)。用戶無須進行繁瑣的算法理論研究，可直接根據需要輸入參數，對各種二維和三維結構進行聲學靈敏度分析。該用戶界面具有可操作性強、運行可靠和界面友好等特點。

1 Burton-Miller奇異邊界法BM-SBM

奇異邊界法即邊界型無網格配點法，通過關于不同源點的基本解的線性組合逼近問題的解。該方法不需要設置虛擬邊界，而是將源點和場點重合布置，利用源點強度因子OIFs代替邊界奇異點。當邊界源點總數為N時，BM-SBM的基本形式為

αiuBM,ii,xi∈Γu,sj∈Γ

(1)

αiqBM,ii,xi∈Γq,sj∈Γ

(2)

(3)

式中：α為待求的未知系數；λ=i/(k+1)，i為虛數單位；uBM,ii和qBM,ii為源點強度因子，可以通過經驗公式獲得；Γu和Γq為邊界Γ上已知邊界條件的邊界部分；ns和nx分別為源點sj和邊界配點xi處的單位外法向量。

OIFs是BM-SBM算法中的關鍵，因此準確求解OIFs十分重要。在已有的OIFs有效算法中，經驗公式最為簡單實用。經驗公式可分為二維問題和三維問題2種情況，

(4)

(5)

(6)

(7)

式中：Lj和Aj分別為源點在二維和三維問題中的影響范圍(見圖1)；γ為歐拉常數；S為域的表面積；二維問題δji=Lj/Li，三維問題δji=Aj/Ai，Li和Ai分別為二維和三維問題場點影響區域的范圍；G0(xi,sj)為拉普拉斯方程的基本解，

(8)

邊界條件分為Dirichlet邊界和Neumann邊界2類，

(9)

(10)

將已知的邊界條件代入式(1)和(2)，得

Aα=b

(11)

式中：A為系數矩陣；α為源點處待求系數向量；b為已知向量。

(a)二維問題

通過求解式(11)可以得到未知系數向量α。將求解得到的未知系數α代入式(12)和(13)，可求得任一點的值

(12)

(13)

2 聲學靈敏度分析

聲音在均勻的各向同性介質中的傳播，可以通過Helmholtz方程描述，

?2p(x)+k2p(x)=0,x∈Ω

(14)

式中：?2為拉普拉斯算子；k為波數，k=ω/c，ω為角頻率，c為空氣中的聲速。

在Dirichlet邊界和Neumann邊界條件下，

(15)

(16)

本文的聲學靈敏度分析涉及輻射聲場，因此有必要為外部聲輻射問題引入無窮遠條件，即無窮遠點上的聲壓應滿足Sommerfeld輻射條件

(17)

式中：r為點x到聲場中心的距離；d為所研究問題的維數，d=2,3。

將已知的邊界條件代入式(1)和(2)，求解得到未知系數α，將α代入式(12)和(13)，并對設計變量進行直接求導，即可求得任一點聲壓關于設計變量的靈敏度值

(18)

(19)

綜上所述，利用BM-SBM求解聲學靈敏度的分析流程見圖2。

圖2 BM-SBM求解聲學靈敏度的分析流程

3 MATLAB程序實現

3.1 GUI界面的創建

GUI是MATLAB軟件附帶的工具箱，GUI提供文本框、按鍵、滑動條、圖形窗口等多種控件供用戶選擇，用戶可通過光標拖動編輯各類控件的位置，并通過回調功能自定義其屬性和功能，從而設計出所需要的圖形界面。GUIDE在用戶進行頁面設置時自動生成可以運行的.m文件，可以大大簡化應用程序的創建過程，用戶可以直接在框架中進行代碼的編寫。

軟件的主要操作流程包括初始化、計算工作和結果輸出3個過程。使用前單擊“幫助”按鈕可獲取使用信息。在初始化過程中，用戶可以進行算例選擇和相關參數(如節點個數、設計變量個數、法向振速、聲速等)設置，單擊“繪制圖形”按鈕可以將模型在二維或三維坐標系中繪制出來。待所有數據都輸入后，單擊“計算”按鈕可完成計算工作。計算結果以圖像方式輸出，也可以單擊“輸出靈敏度數據和輸出誤差數據”按鈕輸出文本文件，以查看測試點的具體計算結果。計算完成后點擊“復位”按鈕可實現軟件重置，以重新選擇算例和參數。軟件界面見圖3。

圖3 聲學靈敏度分析界面

3.2 操作步驟

算例屬性模塊示意見圖4。在算例屬性模塊中，選擇需要計算的算例，并設置好相應的參數，包括算例尺寸、邊界節點數、法向振速、聲速、空氣密度等。若算例為聲散射模型，還需設置入射波方向。

圖4 算例屬性模塊示意

計算設置模塊示意見圖5。在計算設置模塊中，選擇設計變量并輸入設計變量步長，同時設置測試點坐標。若設計變量為波數k，則需要設置波數的取值范圍。

圖5 計算設置模塊示意

在計算結果模塊中，單擊“計算”按鈕，待計算完成后在模塊中顯示計算時間、實部和虛部最大絕對誤差以及實部和虛部全局誤差，同時可以根據需要選擇導出靈敏度和誤差數據。計算結束后進行復位操作，準備下一次運算。

4 數值算例

用2個經典的標準算例[12-13]展示本文開發的工具箱具有計算準確、結果直觀以及易于操作的特點。對一個復雜的二維車腔模型進行聲學靈敏度分析，展示該工具箱在實際工程中的應用。

4.1 二維無限長脈動圓柱聲輻射

取空氣中聲速c=343 m/s，空氣密度ρ=1.2 kg/m3，給定波數k的取值范圍為0.1≤k≤10.0，步長dk=0.1，測試點坐標取(2,2)。使用BM-SBM對測試點聲學靈敏度進行求解，本文開發的工具箱計算結果見圖6，部分k值下的相對誤差見表1。

(a)實部

表 1 部分k值下的相對誤差

由圖6和表1可以看出，該算法的計算精度較高，不同k值下得出的計算結果與精確解之間差值很小，相對誤差穩定在1%以下，因此認為BM-SBM可以準確求解出二維-無限長脈動圓柱(聲輻射)關于波長k的聲學靈敏度問題。與此同時，計算過程反映出該工具箱具有操作簡單、使用便捷、計算速度快、效率高的特點。

4.2 三維振動球聲輻射

(a)實部

圖8 三維振動球的聲學靈敏度誤差

由圖7和8可以看出，不同k值下得到的聲學靈敏度數值結果與相對應的精確解結果高度一致，可見該工具箱能夠精確模擬三維振動球聲輻射關于波長k的聲學靈敏度問題。

4.3 二維車腔聲散射

除對一些特定形狀的模型進行分析外，本文開發的工具箱還可以對任意形狀物體進行聲學靈敏度計算。某車腔模型示意見圖9。將二維車腔模型的邊界節點信息導入到主程序的根目錄下，輸入給定參數為初始聲壓p0=1 Pa，入射波為沿x軸方向pi=p0exp(ikrcosθ)，波數0.1≤k≤10.0，步長dk=0.1，計算外部測試點(4,1)關于波長k的聲學靈敏度，結果見圖10。

圖9 車腔模型示意

(a)實部

5 結束語

建立聲學靈敏度分析的Burton-Miller奇異邊界法BM-SBM數值離散模型，在不同維度下求解聲學靈敏度問題。與傳統方法相比，本文方法能夠避免復雜的網格劃分和奇異積分計算，在簡化計算過程的同時保證較高的精度。在MATLAB軟件的GUI環境下創建對應的計算工具箱，為聲學靈敏度分析提供一種簡單高效的數值方法計算平臺。該工具箱界面友好、操作簡單，不僅可以用于學術研究，還可用于求解實際工程問題。

奇異邊界法是一種邊界型無網格配點法，所建立的矩陣方程為稠密矩陣。隨著計算規模和結構復雜度的增加，本文方法將需要更多的存儲空間和計算時間。因此，未來需要對算法進行加速，建立大規模聲學靈敏度分析的快速奇異邊界法，以優化程序、提高效率。

該用戶界面的開發過程完全基于MATLAB軟件，其中復雜結構的邊界配點借助于三維建模軟件SolidWorks和有限元軟件HyperMesh實現。用戶想要順利使用，無須深入了解奇異邊界法的理論，但可能會用到專業的建模和網格劃分軟件。如果能夠實現建模、節點布置、數值仿真一體化功能，那么會更加方便用戶操作，也會提升軟件的實用價值。