基于COMSOL 多物理場耦合仿真建模方法研究?

李淑君，王惠泉，趙文玉，孟文俊，文 豪

（1.山西三明重工機(jī)械制造有限公司，山西 晉中 030600；2.北京起重運(yùn)輸機(jī)械設(shè)計(jì)研究院，北京 100007；3.太重集團(tuán) 礦山設(shè)備分公司，山西 太原 030024；4.太原科技大學(xué) 機(jī)械工程學(xué)院，山西 太原 030024）

0 引言

目前，在工程實(shí)際中，很多多物理場耦合作用下的實(shí)驗(yàn)還不具備開展條件，且無具體的理論指導(dǎo)設(shè)計(jì)，必須采用數(shù)值仿真［1，2］的方法來研究和測評(píng)。因此，有必要建立一套仿真建模方法來進(jìn)行多物理場耦合分析。COMSOL Multiphysics具有高效的計(jì)算性能能力和獨(dú)特的多場全耦合分析能力，因此可以保證數(shù)值仿真的高度精確［3，4］。本文在相關(guān)理論研究的基礎(chǔ)上，提出基于COMSOL Multiphysics的聲振耦合、熱力耦合、熱聲耦合等的建模方法，并對(duì)軟件在耦合仿真中的前后處理著重進(jìn)行分析。

1 聲振耦合

1.1 聲振耦合仿真建模方法

在聲振耦合界面必須滿足兩個(gè)連續(xù)條件：聲壓連續(xù)、法向加速度連續(xù)［5］。在進(jìn)行聲振耦合仿真時(shí)必須考慮此邊界條件。

聲振耦合仿真建模方法如下：首先在COMSOL Multiphysics軟件中調(diào)用壓力聲學(xué)模塊和結(jié)構(gòu)力學(xué)模塊；然后在壓力聲學(xué)模塊中加載邊界載荷法向加速度，調(diào)用結(jié)構(gòu)力學(xué)中的參數(shù)，在結(jié)構(gòu)力學(xué)模塊中加載邊界條件壓力，調(diào)用壓力聲學(xué)中的壓力參數(shù)。完成設(shè)置后，軟件即可根據(jù)其內(nèi)部算法在代表兩個(gè)物理場的模塊間互相迭代，直至求解結(jié)束。

1.2 應(yīng)用實(shí)例

為驗(yàn)證此方法的正確性，選取一個(gè)簡單的算例應(yīng)用以上方法進(jìn)行聲振耦合仿真分析。將一個(gè)圓柱形的鋁塊置于水中，有一個(gè)平面波入射進(jìn)來，仿真計(jì)算平面波經(jīng)水入射到鋁柱外表面，然后鋁柱外表面反射，最終入射波與反射波疊加的過程。鋁柱浸入水中的模型如圖1 所示，此模型由兩個(gè)域構(gòu)成，一個(gè)為水聲域，一個(gè)為固體域。圓柱外表面即為耦合界面，被聲壓壓緊，圓柱高20mm，半徑為5mm。外圍的水域由一個(gè)球體來模擬無限大水域。分別計(jì)算圓柱做為硬聲場邊界和鋁材料兩種情況下，在60kHz入射波時(shí)的頻響。仿真結(jié)果如圖2 、圖3 所示。

圖1 鋁柱浸入水中的模型

圖2 聲壓、位移和加速度的分布云圖

由圖2 中可以清楚地看到波的傳播路徑，同時(shí)可以看出鋁柱的變形是非常小的，但是它的加速度足夠大，可以影響波的重新分布。

由圖3 中可以看到，在鋁柱表面聲波傳播的上游側(cè)，硬表面時(shí)比鋁表面時(shí)聲壓級(jí)要高，但是在下游側(cè)，硬表面時(shí)比鋁表面時(shí)的振動(dòng)幅值要低。由此可見硬聲場反射的波更多，穿過的波更少，可以得出一個(gè)結(jié)論，材料的吸聲特性對(duì)波形分布有很大的影響。

經(jīng)驗(yàn)證仿真分析結(jié)果正確無誤，驗(yàn)證了聲振耦合建模方法的正確性。

圖3 聲壓隨傳播方向的變化

2 熱力耦合

2.1 熱力耦合仿真建模方法

熱力耦合可以看作兩部分：熱膨脹和振動(dòng)發(fā)熱［6］。相對(duì)應(yīng)于這兩個(gè)條件，在COMSOL Multiphysics軟件中有對(duì)應(yīng)的邊界條件設(shè)置。結(jié)構(gòu)力學(xué)模塊中有熱膨脹邊界條件，設(shè)置好熱膨脹系數(shù)及參考溫度即可；傳熱模塊中有相應(yīng)的體熱源來對(duì)應(yīng)于振動(dòng)發(fā)熱。完成設(shè)置后軟件就可以仿真計(jì)算出相應(yīng)的熱膨脹和振動(dòng)生熱。此為熱力耦合仿真建模方法。

2.2 應(yīng)用實(shí)例

為驗(yàn)證此方法的正確性同樣選取了簡單的模型來進(jìn)行熱力耦合仿真分析。取一個(gè)懸臂梁，材料為鋁合金，其幾何尺寸為10mm×10mm×100mm。固定其左端，在此梁上施加一個(gè)30℃的溫度差載荷，加載時(shí)間為10s。根據(jù)熱膨脹仿真方法進(jìn)行建模，仿真計(jì)算可得其熱變形，仿真結(jié)果如圖4 、圖5 所示。

圖4 10s時(shí)梁的熱變形

圖5 梁的一條邊（x＝9，y＝10）上的位移分布

由圖4 可知，梁的位移主要集中于x方向，屬于線膨脹。其余兩個(gè)方向由于幾何尺寸的限制，位移不大。由圖5 可以看到，在0mm～10mm范圍內(nèi)，位移變化不是線性的，說明，此處包含有體膨脹，同時(shí)也受到了固定端約束的影響。

鋁合金的熱膨脹系數(shù)約為1.8×10－5℃－1～2.4×10－5℃－1，按理論公式計(jì)算可得其理論位移約為0.054mm～0.072mm。計(jì)算結(jié)果與理論相符，在一定程度上說明了該仿真方法的正確性。

3 熱聲耦合

3.1 熱聲耦合仿真建模方法

熱聲耦合效應(yīng)［7］實(shí)際上就是熱與聲的相互轉(zhuǎn)化。熱量分布會(huì)引起傳聲介質(zhì)的密度變化，進(jìn)而影響聲場的分布，同時(shí)由于熱場中各處聲壓不同，熱場分布也會(huì)因此而產(chǎn)生變化。

熱聲耦合仿真建模方法如下：首先在COMSOL Multiphysics軟件中調(diào)用壓力聲學(xué)模塊和傳熱模塊，在壓力聲學(xué)模塊中調(diào)用傳熱學(xué)中的溫度分布參數(shù)，在傳熱模塊中添加聲壓邊界條件；接下來軟件會(huì)在代表熱場和聲場的兩個(gè)模塊之間來回迭代，每次運(yùn)算都要調(diào)用前一次的結(jié)果，進(jìn)而仿真出熱和聲之間的相互影響。

3.2 應(yīng)用實(shí)例

為驗(yàn)證此方法的正確性同樣選取了簡單的模型來進(jìn)行熱聲耦合仿真分析。建立一個(gè)正方形的空氣域模型，分兩種情況進(jìn)行了模擬，第一種情況下溫度場分布均勻，第二種情況下左側(cè)溫度比右側(cè)溫度高60℃，達(dá)到穩(wěn)態(tài)后溫度沿x軸為線性分布。兩種情況下，在左側(cè)加一個(gè)入射平面波，其頻率為500Hz，大小為2Pa。圖6 為熱場分布對(duì)聲場分布的影響。由圖6 中可以看出，有溫度場分布情況下聲場分布更密集一些。

根據(jù)理論我們可以知道，溫度高的地方氣體的密度會(huì)下降，在此處的聲速就會(huì)下降，在頻率不變的情況下，其波長就會(huì)變短，其聲場分布就會(huì)變得密集。

仿真結(jié)果與理論推導(dǎo)一致，說明該仿真建模方法的正確性。

圖6 熱場分布對(duì)聲場分布的影響

4 結(jié)論

根據(jù)由多場耦合理論得出的多場耦合機(jī)理，結(jié)合COMSOL軟件的特征，開發(fā)了多場耦合建模方法，包括聲振耦合建模方法、熱力耦合建模方法、熱聲耦合建模方法，并通過幾個(gè)簡單算例驗(yàn)證了文中建模方法的正確性，應(yīng)用此方法即可解決多場耦合下的仿真問題。

［1］楊鼎寧，鄒經(jīng)湘，蓋登宇.計(jì)算機(jī)輔助工程（CAE）及其發(fā)展［J］.力學(xué)與實(shí)踐，2005，27（3）：7-16.

［2］付德彬.數(shù)值仿真及其在航天發(fā)射技術(shù)中的應(yīng)用［M］.北京：國防工業(yè)出版社，2011.

［3］王堯.油膜軸承巴氏合金與鋼體的結(jié)合強(qiáng)度理論與試驗(yàn)研究［D］.太原：太原科技大學(xué)，2014：1-15.

［4］Teich L，Hütten A，Schr?der C.Utilization of COMSOL Multiphysics'JAVA API for the implementation of a micromagnetic modeling and simulation package with a customized user interface ［C］／／COMSOL Conference Europe.Mailand：［s.n.］，2012.

［5］白長青，周進(jìn)雄，閆桂榮.聲振耦合對(duì)薄壁圓柱結(jié)構(gòu)動(dòng)力特性的影響［J］.機(jī)械工程學(xué)報(bào)，2011，47（5）：78-84.

［6］張胤.材料傳熱學(xué)［M］.北京：化學(xué)工業(yè)出版社，2013.

［7］黃鑫，胡軍，李青，等.Rijke型熱聲自激振蕩機(jī)理研究進(jìn)展［J］.低溫工程，2010（1）：5-10.