張汝毅， 王發(fā)杰， 程隋福， 劉建政

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

0 引 言

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

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

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

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

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

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

(1)

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

(2)

(3)

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

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

(4)

(5)

(6)

(7)

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

(8)

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

(9)

(10)

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

Aα=b

(11)

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

(a)二維問題

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

(12)

(13)

2 聲學(xué)靈敏度分析

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

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

(14)

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

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

(15)

(16)

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

(17)

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

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

(18)

(19)

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

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

3 MATLAB程序?qū)崿F(xiàn)

3.1 GUI界面的創(chuàng)建

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

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

圖3 聲學(xué)靈敏度分析界面

3.2 操作步驟

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

圖4 算例屬性模塊示意

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

圖5 計算設(shè)置模塊示意

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

4 數(shù)值算例

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

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

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

(a)實部

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

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

4.2 三維振動球聲輻射

(a)實部

圖8 三維振動球的聲學(xué)靈敏度誤差

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

4.3 二維車腔聲散射

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

圖9 車腔模型示意

(a)實部

5 結(jié)束語

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

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

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