基于正交試驗的光聲池氣流性能模擬及參數(shù)優(yōu)化

程 剛, 曹亞南, 田 興, 曹 淵, 劉 錕

1. 安徽理工大學深部煤礦采動響應與災害防控國家重點實驗室, 安徽 淮南 232001 2. 中國科學院安徽光學精密機械研究所, 安徽 合肥 230031

引 言

利用光聲光譜原理實現(xiàn)對氣體濃度的檢測是光聲技術中最為典型的應用。 與其他光譜類氣體檢測方法相比, 光聲光譜氣體檢測技術主要有結構簡單、 探測器不受波長的限制、 零背景噪聲、 成本低等優(yōu)點, 在氣體檢測領域已得到了廣泛的認同與應用。 光聲池作為光聲光譜氣體檢測系統(tǒng)中最為核心的部件, 其性能的優(yōu)劣直接對系統(tǒng)檢測結果產(chǎn)生重要影響, 因而對光聲池開展優(yōu)化成為了該領域研究的熱點。 在光聲池性能優(yōu)化方面, 國內(nèi)外學者進行了大量的研究, 取得了許多有益的結論和結果。 Kost等通過數(shù)值模擬獲得了一種由多個軸對稱錐體連續(xù)連接形狀的光聲池結構[1]; Haouari等采用聲學三維拓撲數(shù)值計算與平滑處理技術獲得了一種形如“土豆”結構形狀的光聲池[2]; Mikael等從圓柱腔圓柱面開設腔, 并采用計算流體動力學分析了光聲池的流致噪聲[3]; Cheng等提出一種聲強得到增益的“喇叭口”形狀及階梯復合形光聲池[4-5]; Zhang等設計出一種高靈敏度的橢球形共振光聲池[6]; Li等設計了一種母線為雙曲線的曲體束腰型光聲池結構方案[7]; Liu設計了一種“T”型光聲池[8]; Shi等、 Zhao等及Jiao等對球形光聲池的性能開展相關研究[9-11]; Zhao等設計了一種拱形光聲池, 并確定了其最佳結構尺寸[12]; Yu等設計并試驗了一款基于耦合積分球光聲池的在線氣體檢測裝置[13]; Yin等設計了具有兩個相同諧振器的差分式光聲池[14]; Ma等提出一種利用棱鏡作用實現(xiàn)多通道反射增強型光聲池[15]。 通過優(yōu)化光聲池的結構, 上述研究取得了良好的預期效果, 為光聲光譜氣體檢測性能的提升以及光聲池的結構創(chuàng)新設計提供了大量的技術支持和參考。 然而, 上述大多數(shù)研究主要是在系統(tǒng)的靜態(tài)條件下進行的, 光聲光譜在線檢測系統(tǒng)中光聲池中氣體流動性能和動態(tài)時間響應等方面研究報道較少, 這無疑構成了進一步提高光聲光譜檢測性能的障礙。 因此, 在前人研究的基礎上, 探索動態(tài)條件下的光聲池參數(shù)影響規(guī)律與優(yōu)化組合進而提高光聲池的檢測性能具有重要意義。

基于三維流場數(shù)值模擬方法, 建立了光聲池流場的穩(wěn)態(tài)和瞬態(tài)模擬模型, 重點研究了光聲池腔內(nèi)氣體流動過程引起的動態(tài)壓力效應和時間響應, 利用正交試驗設計和熵權法分析了光聲池參數(shù)的影響規(guī)律, 通過對指標的權重平衡, 提出了光聲池的最佳參數(shù)組合, 相關結論可為光聲池的后續(xù)擴展和優(yōu)化研究提供借鑒。

1 計算模型與三維流場模擬

1.1 計算模型

用于光聲光譜氣體檢測的圓柱形光聲池結構示意圖與實物如圖1所示。 光聲池與待測氣體接觸的內(nèi)腔部分主要包括光聲池Z軸向腔內(nèi)正中間的圓柱形諧振腔和其兩端的圓柱形緩沖腔, 以及殼體的表面對稱開設有進、 出氣孔形成的孔腔。 當光聲池處于動態(tài)檢測過程中時, 被測氣體的流動路徑是從進氣孔流向靠近進口孔的緩沖腔, 然后沿著諧振腔流向靠近出氣孔的緩沖腔, 最后從出氣孔流出。 現(xiàn)有的光聲池結構原始尺寸為: 光聲池殼體總長L=220 mm, 其Z軸端面為正方形, 尺寸a×a=100 mm×100 mm; 諧振腔半徑Rre=5 mm, 長度Lre=100 mm; 緩沖腔半徑Rbuff=35 mm, 長度Lbuff=60 mm; 進、 出氣孔半徑Rin=Rin=5 mm, 進、 出氣孔中心距離Ljc=170 mm。 光聲池三維流場計算物理模型是光聲室腔內(nèi)形成的封閉氣體區(qū)域。

圖1 光聲池結構模型圖(a): 模型圖; (b): 實物圖Fig.1 Structure diagram of photoacoustic cell(a) Model diagram; (b): Physical photos

1.2 控制方程

光聲池腔內(nèi)氣體流場模擬采用計算流體動力學(computational fluid dynamics, CFD)方法, 腔內(nèi)氣體流動需要遵循流場基本控制方程的約束。 因不涉及熱量交換, 所以控制方程主要由質(zhì)量守恒方程、 動量守恒方程組成, CFD計算均以上述方程為基礎。 在笛卡爾坐標系中可寫出對應的氣體控制方程通用形式如式(1)[16]

(1)

式(1)中,φ為通用變量,Γ為廣義擴散系數(shù),S為廣義源項。

其質(zhì)量守恒方程為

(2)

式(2)中,ρ為流體密度,u為速度矢量,t為時間。 式(2)是瞬態(tài)三維可壓流體的質(zhì)量守恒方程, 若流體視為不可壓, 式(2)可變?yōu)?/p>

(3)

其動量守恒方程為

(4)

式(4)中,μ為流體的動力粘度,Su、Sv、Sw為動量守恒方程中的廣義源項。

1.3 基本假設與邊界條件

抽取光聲室腔內(nèi)形成的封閉氣體區(qū)域為流體求解域。 采用穩(wěn)態(tài)和瞬態(tài)流場模擬開展計算, 基本假設如下: ①簡化光聲池的細節(jié)結構, 忽略腔內(nèi)的小凸臺、 小凹槽等幾何模型細節(jié)特征, 不考慮因機械加工帶來的結構誤差; ②腔內(nèi)氣體視為理想性氣體, 為不可壓縮流體介質(zhì), 且不考慮濕度影響; ③視腔內(nèi)表面為理想光滑表面, 不考慮腔內(nèi)粗糙度對氣流的影響。

邊界條件: ①腔內(nèi)氣壓初始狀態(tài)為0.1MPa, 介質(zhì)為空氣, 體積分數(shù)為100%, 通入氣體為純甲烷氣體, 體積分數(shù)為100%; ②進氣孔處邊界設為氣流速度入口,uin=0.06～0.10 m·s-1, 并設置為“充分發(fā)展流動”, 出氣孔處邊界設為壓力出口, 為0.1 MPa, 考慮氣體產(chǎn)生的重力影響; ③腔內(nèi)溫度恒定為293 K, 腔內(nèi)壁面為絕熱無滑移邊界; ④腔內(nèi)氣流動態(tài)壓力模擬設置為CFD穩(wěn)態(tài)計算, 腔內(nèi)氣體濃度調(diào)節(jié)時間模擬設置為CFD瞬態(tài)計算, 設定模擬時間t=300 s, 計算步長為3 s。

1.4 網(wǎng)格獨立性驗證

網(wǎng)格劃分對計算收斂值具有重要影響, 其精細度也決定著計算機運算時間。 驗證網(wǎng)格獨立性可以在計算精度和時間之間達成折衷。 設置計算模型網(wǎng)格劃分為三維自動劃分, 最小縫隙尺寸為2 mm, 每次運算均提高初始網(wǎng)格的精度級別。 模型網(wǎng)格計算圖如圖2所示, 諧振腔、 進氣口、 出氣口以及過渡處均進行網(wǎng)格自動加密, 整個模型采用了結構化網(wǎng)格進行處理。 設定uin=0.06 m·s-1, 計算目標為諧振腔軸線中點處動壓值Pj0, 其隨網(wǎng)格數(shù)量變化的影響規(guī)律如圖3所示, 當網(wǎng)格數(shù)量在分別在42萬和14萬左右時, 諧振腔軸線中點處動壓值Pj0變化率僅為3.7%, 趨于穩(wěn)定狀態(tài), 繼續(xù)增加網(wǎng)格數(shù)量對結果影響不明顯, 因而綜合確定光聲池流場模擬計算時的網(wǎng)格精度初始級別為6級, 網(wǎng)格數(shù)量為14萬左右。

圖2 計算網(wǎng)格圖(a): 整體網(wǎng)格圖; (b): 進氣孔處網(wǎng)格圖; (c): 諧振腔與過渡處網(wǎng)格圖Fig.2 Calculation grid diagram(a): Overall grid diagram; (b): Grid diagram at air inlet hole; (c): Grid diagram of resonator and transition

圖3 網(wǎng)格獨立性驗證Fig.3 Grid independence verification

1.5 模擬結果

1.5.1 穩(wěn)態(tài)模擬

模擬求解穩(wěn)態(tài)條件下光聲池縱向?qū)ΨQ截面的流場分布如圖4所示, 可知光聲池腔內(nèi)進、 出氣孔以及諧振腔處的動壓值過大。 過大的動壓會引起較強的氣流擾動, 從而增加信號檢測噪聲, 特別是諧振腔軸線中部的動壓波動會直接對聲學傳感器的檢測造成負面影響; 此外, 從圖4中還可以看出, 光聲池兩端緩沖腔中存在氣流渦旋現(xiàn)象, 這也會增強流體噪聲影響。

圖4 腔內(nèi)動壓與流線分布切面圖Fig.4 Section diagram of dynamic pressure and streamline distribution in the cavity

模擬計算進氣速度uin分別為0.06、 0.08和0.10 m·s-1時, 諧振腔軸線上動壓分布規(guī)律, 如圖5所示, 可知諧振腔軸線上動壓大小隨著進氣速度uin的增大而增大, 由文獻[3, 17]可知, 流速越大腔內(nèi)渦旋回流將會越大, 導致的氣體噪聲越強, 因而動壓與流速可以等同地表征氣體的擾動影響。 氣體從緩沖腔流進和流出諧振腔時動壓波動較為明顯, 且氣流進入諧振腔時的動壓變化梯度強于氣體流出諧振腔時的動壓變化梯度, 但在諧振腔中動壓表現(xiàn)的相對平穩(wěn), 上述是由于腔體截面直徑的突擴和突縮流體局部阻力變化引起的。 因此可以預見, 減少腔內(nèi)氣流流速及優(yōu)化光聲池中的過渡結構將會改善氣流引發(fā)的動壓波動程度。

圖5 不同進氣速度下諧振腔軸線動壓值分布Fig.5 Dynamic pressure distribution at the axis of the resonator at different inlet speeds

1.5.2 瞬態(tài)模擬

模擬求解瞬態(tài)條件下光聲池腔內(nèi)待測氣體變化規(guī)律如圖6所示, 調(diào)節(jié)時間t設定為當光聲池中甲烷氣體體積占比到達并保持在穩(wěn)態(tài)值的c=95%所需的最短時間, 探測點設置為諧振腔軸線的中點。 由圖6可知, 在不同的進氣速度下, 腔內(nèi)氣體的體積比例隨著時間的增加而增加, 并逐漸達到穩(wěn)定的狀態(tài), 在早期, 腔體中待測氣體的體積比例增長速度較快, 在后期, 其增長速度逐漸減慢。uin為0.06、 0.08、 0.10 m·s-1時, 調(diào)節(jié)時間t模擬結果分別為171、 132和108 s, 這是由于氣體流速增大時, 待測氣體在腔內(nèi)的擴散速度將隨之增大, 因此可以預見, 提高腔內(nèi)的氣體流動速度可以縮短調(diào)節(jié)時間, 這與實際也非常吻合, 進而提高系統(tǒng)檢測的動態(tài)響應速度。 結合穩(wěn)定模擬的動壓結果可知, 腔內(nèi)氣流流速將是光聲池內(nèi)動壓波動與動態(tài)響應的重要影響參量。

圖6 光聲池腔內(nèi)氣體濃度調(diào)節(jié)時間Fig.6 Adjustment time vs gas concentration in photoacoustic cell cavity

2 結構方案設計

通過對光聲池腔內(nèi)流場的模擬分析可知, 腔內(nèi)的動壓、 濃度調(diào)節(jié)時間與流速及結構突變等因素密切相關, 因此改善光聲池的內(nèi)腔結構具有重要意義。 為充分利用現(xiàn)有的光聲池, 最大限度減少其加工成本, 在不改變原光聲池主體結構的情況下, 采用模塊化設計策略, 如圖7所示, 具體如下: ①將光聲池殼體沿緩沖腔的軸線并按其直徑尺寸將其加工成通孔; ②將原諧振腔結構改為諧振腔組件, 諧振腔組件為獨立圓柱體結構, 其外徑大小與緩沖腔相同, 諧振腔組件可沿Z軸進行移動, 處于中間位置時, 光聲池Z軸兩端可以自然形成緩沖腔; ③根據(jù)流場模擬結果及文獻[17]可知, 對腔體截面突變處進行適當過渡處理可以進一步提高氣流局部阻力變化; ④為進一步改善諧振腔中的氣體流速的平穩(wěn)性, 在諧振腔組件上沿Z軸嘗試開設均勻分布的輔助孔。

圖7 光聲池結構改進方案示意圖Fig.7 Schematic diagram of structural improvement scheme of photoacoustic cell

3 正交試驗分析

3.1 試驗方案與計算

正交試驗是利用數(shù)理統(tǒng)計學和正交性原理開展合理試驗的一種科學方法, 特別適合研究和處理多因素優(yōu)化分析。 采用正交試驗的方法對光聲池的結構改進參數(shù)開展研究。 將緩沖腔與諧振腔過渡處圓角Ra(A)(2.0～5.0 mm)、 輔助孔數(shù)量為Nf(B) (2～8個)、 輔助孔半徑Rf(C)(2.5～4.0 mm)、 輔助孔中心圓半徑Rz(D)(20.0～27.5 mm)以及進氣速度uin(E)(0.06～0.08 m·s-1)5個參數(shù)選為考察因素, 每個因素取4個水平, 其他參數(shù)保持不變, 不考慮因素之間的交互影響, 以諧振腔軸線中點處動壓值Pj0、 腔內(nèi)氣體濃度調(diào)節(jié)時間t為試驗評價指標, 因素-水平表如表1所示。

表1 因素-水平表Table 1 Factor level table

采用正交表L16(45)安排試驗, 需要開展16次試驗, 試驗方案及計算如表2所示。

表2 試驗方案及計算結果Table 2 Test scheme and calculation results

3.2 極差分析

試驗條件下的均值與極差計算結果如表3所示, 從試驗結果中可直觀看出, 第7組試驗組合A2B3C4D1E2中動壓值Pj0最小, 此時Pj0=9.1×10-4Pa, 調(diào)節(jié)時間t=171 s; 第5組試驗組合A2B1C2D3E4中調(diào)節(jié)時間t最小, 此時Pj0=1.1×10-2Pa, 調(diào)節(jié)時間t=117 s。 極差R表示各個因素對指標結果的影響程度,R=max(k1,k2,k3,k4)-min(k1,k2,k3,k4),kj為第j列同一水平試驗結果的平均值, 極差R越大, 表明該因素對試驗結果影響越大, 其重要性越強。 從表3極差結果可知, 光聲池的各參數(shù)對動壓值Pj0影響的主次順序為: 輔助孔半徑C>輔助孔數(shù)B>進氣速度E>過渡圓角A>輔助孔中心圓半徑D, 且各因素的影響程度區(qū)分度不大; 各參數(shù)對調(diào)節(jié)時間t影響的主次順序依次為: 進氣速度E>輔助孔半徑C>輔助孔數(shù)B=輔助孔中心圓半徑D>過渡圓角A, 且進氣速度遠高于其他因素對調(diào)節(jié)時間的影響。

表3 均值與極差Table 3 Mean and range

動壓值Pj0和調(diào)節(jié)時間t均期望越小越好, 因而選取因素水平應該取使指標小的水平。 在所做試驗的范圍內(nèi), 僅考慮動壓值Pj0時, 光聲池的結構參數(shù)最優(yōu)組合為A1B4C4D1E1&A1B4C4D4E1, 上述參數(shù)組合在正交試驗中未出現(xiàn), 通過仿真計算獲得此時Pj0=4.6×10-4Pa &Pj0=4.6×10-4Pa, 調(diào)節(jié)時間t=216 s &225 s; 僅考慮調(diào)節(jié)時間t時, 光聲池的結構參數(shù)最優(yōu)組合為A2B1C3D2E4, 此時Pj0=8.0×10-3Pa, 調(diào)節(jié)時間t=114 s。 上述結果表明, 通過正交試驗極差的分析, 在參數(shù)試驗范圍內(nèi)可實現(xiàn)光聲池最佳的動壓值和調(diào)節(jié)時間的參數(shù)組合, 顯示了所采用正交試驗的優(yōu)越性和可行性。

3.3 熵權法及最優(yōu)參數(shù)組合

光聲池的最優(yōu)參數(shù)組合需要同時考慮動壓值和調(diào)節(jié)時間2個評價指標, 屬于多目標參數(shù)優(yōu)化, 同時兩者為望小特征。 從上述正交試驗中可知, 動壓值和調(diào)節(jié)時間兩者是互相沖突的, 即不能同時取到兩者的最小值。 為兼顧兩個指標的影響, 找到最佳指標下的參數(shù)組合, 需將多目標參數(shù)優(yōu)化問題轉(zhuǎn)化成單目標優(yōu)化問題。 常見的方法采用單目標權重法來進行處理, 為客觀地給出動壓值和調(diào)節(jié)時間的權重, 采用熵權法建立評價模型, 利用信息熵計算兩者的權重。 數(shù)據(jù)分析對象為正交試驗中的16組試驗計算結果。 同時, 將正交試驗中僅考慮動壓值或調(diào)節(jié)時間的計算結果作為補充數(shù)據(jù), 共有19組Pj0&t指標數(shù)據(jù), 補充數(shù)據(jù)如表4所示。

表4 補充數(shù)據(jù)Table 4 Supplementary data

由于動壓值和調(diào)節(jié)時間兩者的量綱不同, 需對19組Pj0&t數(shù)據(jù)進行標準化處理, 因兩者均為負向指標, 標準化公式如式(7)[18]

(7)

(8)

Pj0&t的信息熵eij計算公式如式(9)

(9)

Pj0&t的信息熵ωij計算公式如式(10)

(10)

數(shù)據(jù)集中每個方案的綜合評價得分si計算方式如式(11), 綜合評價得分的值越大越好。

(11)

動壓值和調(diào)節(jié)時間權重計算結果如表5所示。

表5 熵權法計算結果Table 5 Calculation results of entropy weight method

由表5可知, 在試驗范圍內(nèi), 從客觀角度可知光聲池調(diào)節(jié)時間的權重略大于動壓值的權重。 結合式(11)可計算出19組Pj0&t指標數(shù)據(jù)的綜合評價得分前3組的依次為第8組、 第4組、 第14組,s8=0.849,s4=0.848,s14=0.774。 根據(jù)以上分析, 并結合光聲池制作成本及實際需求等, 可以合理選擇其參數(shù)組合, 選擇綜合評分最高的第8組參數(shù)組合為其最優(yōu)參數(shù)組合。

3.4 優(yōu)化方案對比

設定其進氣速度uin均為0.06 m·s-1, 對比光聲池優(yōu)化前后的結果。 圖8為第8組參數(shù)組合下的腔內(nèi)動壓與流線分布切面圖, 結合圖4可知, 優(yōu)化前后的光聲池腔內(nèi)最大流速基本保持不變, 優(yōu)化后的腔體中的流速更加均勻、 平滑, 腔體中的動壓大幅降低。 圓角的處理使得諧振腔兩端處的腔體截面直徑的突擴和突縮引起的動態(tài)變化變得更加緩和, 同時, 出氣孔端的緩沖腔內(nèi)的氣流渦旋現(xiàn)象也得到了改善。 圖9優(yōu)化后腔內(nèi)球狀粒子流線圖, 可知光聲池腔內(nèi)氣體均勻地從輔助孔通道中分散流動, 以致氣流在諧振腔中更加低速、 平滑地流動, 如此可以降低光聲池腔內(nèi)的氣流擾動噪聲, 從而最終改善光聲池的信號檢測噪聲。

圖8 優(yōu)化后腔內(nèi)動壓與流線分布切面圖Fig.8 Section diagram of dynamic pressure distribution and streamline distribution in optimized cavity

圖9 優(yōu)化后腔內(nèi)球狀粒子流線圖Fig.9 Streamline diagram of spherical particles in cavity after optimization

圖10為光聲池優(yōu)化前后振腔軸線中點處動壓值Pj0、 腔內(nèi)氣體濃度調(diào)節(jié)時間t兩個指標的結果對比圖, 可知優(yōu)化后的動壓值Pj0為9.4×10-4Pa, 調(diào)節(jié)時間t為141s, 相較于優(yōu)化前的指標, 動壓值相對降低了88.1%, 調(diào)節(jié)時間相對降低了17.5%, 兩項評價指標均得到優(yōu)化提升, 優(yōu)化效果較為理想。

圖10 優(yōu)化前后結果對比Fig.10 Comparison of results before and after optimization

4 結 論

(1)基于三維流場數(shù)值模擬方法, 建立了光聲池流場的穩(wěn)態(tài)和瞬態(tài)模擬模型, 計算并獲得了光聲池腔內(nèi)氣體流場分布及其氣體濃度平衡響應情況和規(guī)律。 結果表明, 降低腔內(nèi)氣流流速、 優(yōu)化光聲池中的過渡結構可改善氣流引發(fā)的動壓波動程度以及縮短調(diào)節(jié)時間;

(2)運用正交試驗設計與熵權法, 研究了光聲池結構參數(shù)對諧振腔軸線中點處動壓值和氣體濃度調(diào)節(jié)時間的影響規(guī)律。 在研究參數(shù)范圍內(nèi), 獲得了其最佳參數(shù)組合為: 緩沖腔與諧振腔過渡處圓角為3.0 mm、 輔助孔數(shù)量為8個、 輔助孔半徑為3.5 mm、 輔助孔中心圓半徑為22.5 mm、 進氣速度為0.06 m·s-1;

(3)優(yōu)化后的光聲池諧振腔軸線中點處動壓值為9.4×10-4Pa, 腔內(nèi)氣體濃度調(diào)節(jié)時間為141 s, 相較于優(yōu)化前的指標, 動壓值相對降低了88.1%, 調(diào)節(jié)時間相對降低了17.5%。