計算氣動聲學(xué)CAA若干學(xué)習(xí)經(jīng)驗分享，強烈推薦！

2016-10-26 by:CAE仿真在線來源:互聯(lián)網(wǎng)

在論壇上看到越來越多的人也在做氣動聲學(xué)相關(guān)的東西,頗有得遇同道中人的喜悅。本人在碩士階段就開始接觸一些氣動聲學(xué)相關(guān)的東西,工作后主要的研究內(nèi)容就更專一了:航空聲學(xué)。工作一年后,通過各種亂七八糟的學(xué)習(xí)過程,對計算氣動聲學(xué)有了更多的理解。受版主水若無痕的影響(他是我的同學(xué)),因此打算在此寫個與計算氣動聲學(xué)(CAA)相關(guān)的東西,和大家交流交流。

對氣動聲學(xué)的關(guān)注始于上世紀的50年代,原因就是當時渦噴式航空發(fā)動機的噴流噪聲實在是太嚇人了。于是,牛逼的萊特希爾(Lighthill)坐在火車上,在一個信封上一頓寫,就把N-S方程給改寫成了波動方程的形式。方程的左邊是一個經(jīng)典聲學(xué)的波動方程,而右邊則是一個主要與湍流相關(guān)的源項,被后人稱為萊特希爾應(yīng)力張量。這就是所謂的萊特希爾方程了,氣動聲學(xué)的開山之作。萊爾希爾方程的聲源為四極子聲源,也就是湍流噪聲源,主要適用于高速、湍流為主要噪聲源的情況,如高速噴流。方程的聲源項未知,需要采用CFD或者試驗來獲取。

再后來,柯爾(Curler)同志對萊特希爾方程進一步發(fā)展,得出了考慮了固壁影響的柯爾方程。柯爾方程主要適用于低速情況下的固壁繞流噪聲計算,如低速的圓柱繞流、機翼繞流等。此時,氣動噪聲源主要為偶極子聲源,聲源的強度為聲源表面對流體的作用力。這種作用力不單是壓力,還包括表面動量流量。當然,對于固壁來說,法向速度為零,也就沒有動量流量了,因此采用固壁表面作為聲源面時,只需要壁面的壓力脈動即可。而在采用通流面作為積分面時,則需要考慮動量流量了,這在后面會有介紹。

福茨威廉斯與霍金斯(Ffcows Williams&Hawkings)兩位在萊特希爾方程的基礎(chǔ)上,發(fā)展出FW-H方程。FW-H方程的發(fā)展主要是針對運動壁面的發(fā)聲情況。這里說的運動壁面指的是在來流中的運動,也就是說壁面具有加速度,如螺旋槳。FW-H方程包含了所有的噪聲源,單極子、偶極子和四極子。這三種聲源的發(fā)聲效率遞減,指向性差異很大。一般來說,FW-H方程能夠描述所有的氣動噪聲問題,只不過你需要根據(jù)你計算問題的具體情況,來確定哪種噪聲源為主,哪種噪聲源可以忽略?，F(xiàn)在主流的氣動聲學(xué)計算軟件基本上都用的是FW-H方程。

上面大概介紹了一下氣動聲學(xué)理論方面你的東西。具體的方程形式復(fù)雜,推導(dǎo)困難,我是不會的。不過隨便找本相關(guān)的書都有這方面的介紹,大家可以好好看看。這三個方程有個一致的假設(shè),就是聲場與流場不存在相互影響。這三個方程的主要作用有兩個:一是告訴了我們聲源的發(fā)聲機理,以及怎么由流場參數(shù)去求聲源參數(shù);二是方程的積分解可以用來解決一些簡單的氣動聲學(xué)問題,后面會提及。

有了這些方程后,我們就應(yīng)該想著去計算氣動噪聲了。一個完整的氣動噪聲計算應(yīng)該包括以下三個部分:聲源計算、聲傳播計算和聲輻射計算。如下面這張圖片所示。

1 聲源計算

如果你能看到上述的三個方程,那么你會發(fā)現(xiàn)方程右邊的源項主要是由流場中的參數(shù)構(gòu)成的,比如雷諾應(yīng)力、粘性應(yīng)力、壓力、動量流量、質(zhì)量流量等等。顯然,這些參數(shù)的獲取則需要靠CFD計算了,這也就是所謂的聲源計算了。因此,準確的CFD流場計算是CAA計算的基礎(chǔ)。一般來說,聲源的計算需要采用高級的湍流模型,如LES、DES,進行非穩(wěn)態(tài)計算來實現(xiàn)。當然,URANS模型也可,不過計算出來的結(jié)果往往只有主要特征的影響,得不到細節(jié)方面的東西。需要補充的是,穩(wěn)態(tài)計算,也能在一定的程度上得到聲源分布,這種牛逼的算法叫隨機噪聲產(chǎn)生與傳播(Stochastic Noise Generation and Radiation,SNGR),它是根據(jù)湍動能的概率分布來轉(zhuǎn)換噪聲源的。德國人在這方面做的比較多,FLUENT里面所謂的寬頻噪聲源法采用的也是這種理論。

2 聲傳播計算

問題的復(fù)雜性不同,聲傳播計算的重要性也就不同。對于一些簡單問題,聲傳播計算并不是必要的。比如說低速圓柱繞流噪聲。不需要考慮聲傳播計算的主要原因是結(jié)構(gòu)和流場對聲傳播的影響可以忽略。因此,對于這種簡單的情況,直接利用CFD的計算結(jié)果,對FW-H方程求積分解,就可以得到圓柱繞流的輻射噪聲了。FLUENT里面的聲學(xué)模塊就是這么干的。

這里說一下FW-H方程積分法。積分法的思想是比較簡單的,它實際上是一種解析方法,當然求起來就不一定簡單了。以圓柱繞流噪聲計算為例,我們把圓柱表面劃分成一個個的微元,每個微元當做一個偶極子聲源,聲源的強度就是脈動壓力大小,然后對每一偶極子求解波動方程,得到其在聲接受點處的輻射聲壓,再把每個偶極子所產(chǎn)生的聲壓相加,就得到了總聲壓,問題就解決了。不過積分法里面一個很重要的參數(shù)是延遲時間。所謂延遲時間指的是,聲源產(chǎn)生的聲波傳到聲接收點需要一定的時間。聲源位置不同,延遲時間也就不同,因此在編程計算中是一個比較麻煩的東西。

而對于復(fù)雜的問題,那么聲傳播的計算則是必須的了。舉幾個例子,一是管內(nèi)流噪聲問題。聲波產(chǎn)生后,在管道內(nèi)的傳播受到管道結(jié)構(gòu)的影響,必須以某種聲模態(tài)來傳播。這時候,你想用FLUENT來計算噪聲那是行不通的。第二個例子,航空發(fā)動機的風(fēng)扇噪聲的計算。航空發(fā)動機短艙在進氣段一般都安裝了消聲用的聲襯,也就是說邊界條件直接對聲傳播有影響,所以也必須進行聲傳播計算。第三個例子,航空發(fā)動機噴流噪聲計算。發(fā)動機尾流中存在強烈的剪切層,有著很多大的速度梯度和溫度梯度,將影響到聲的傳播,因此必須進行聲傳播計算?？偟脕碚f,當結(jié)構(gòu)和流場對聲的傳播有這不可忽略的影響,你就得做聲傳播計算了。

3 聲輻射計算

聲輻射計算的目的在于求取遠場的聲壓。當然,如果你的計算機很牛逼,你的聲傳播計算算法也很牛逼,那么你可以直接建立龐大的計算域,直接采用傳播計算來求取遠場的聲壓。否則的話,你就得用上面說過的FW-H方程積分法來求了。目前來說,這個難度不大,對算法沒什么要求,方程也比較好求解,花點時間基本上也能搗騰出來。或者直接用FLUENT現(xiàn)有的也可以,因此不多說了。

好了,通過上面的介紹,我們知道,一個問題的難度主要就在于是否要考慮聲傳播計算了。為什么不能用傳統(tǒng)的CFD算法直接來計算聲傳播呢?主要的原因有兩個:一是CFD的控制方程,N-S方程也好,Euler方程也好,都是非線性的,而聲壓相對于流場壓力是一個非常小的量,在計算過程總很快就會被忽略掉;二是聲壓信號是一個非常小的量,傳統(tǒng)的低階格式對聲壓來說具有很大的耗散性,將導(dǎo)致聲壓信號迅速被耗散掉。因此,對CAA來說,真正有技術(shù)含量的工作都集中在聲傳播計算這一塊。

接下來將針對簡單問題和復(fù)雜問題分別介紹一些我做的計算。

簡單問題

由上面的介紹我們知道,簡單問題只需要采用FW-H方程的積分解,就可以解決。這對一些簡單的工程問題已經(jīng)足夠了。對于簡單問題,計算的重點有兩個:一是CFD流場的計算,你必須算對了算準了才行;二是聲源面的選取了。前一個問題在于你對CFD的理解和掌握程度了,這里主要說一下第二個問題。

聲源面的選取有兩種方法,一是直接選取固體壁面作為聲源面,二則是選取虛擬的通流面作為聲源面,即所謂的積分面。積分面就是采用通流面,對應(yīng)著FLUENT里面的interior邊界條件,將聲源區(qū)域給包起來。采用積分面的原因是為了簡化計算當中的體積分問題。我們知道,對于固體邊界不是主要噪聲源的情況,如噴流噪聲,噪聲主要由湍流產(chǎn)生,而湍流的分布則是空間分布了。如果我們采用積分法來求解湍流噪聲源的聲輻射,那么就必須采用體積分了,這將導(dǎo)致計算量的劇增。而采用積分面,則可將體積分轉(zhuǎn)換為面積分,從而大大減小計算量。

積分面的選取方法有多種,用的比較多的主要由基爾霍夫(KirchhoffMethod)法和可通流面FW-H方程(PorousSurface FW-H Equation)。FLUENT當中的積分面法采用的應(yīng)當是后一種(不確定,里面介紹的不多)。積分面的選取比較麻煩。一個比較重要的注意事項是:應(yīng)當避免尾流或者渦直接穿過積分面,導(dǎo)致聲信號的污染。此外,積分面的位置不同,對計算結(jié)果有一定的影響。

下面的一個例子是利用FLUENT,采用固壁面和積分面兩種方法來計算二維圓柱繞流噪聲。圓柱直徑50mm,來流馬赫數(shù)0.12,采用的湍流模型為SST。網(wǎng)格如下圖所示,其中綠色線框就是兩個積分面了,記為IS100和IS200。而紅色的面就是積分面的封閉面,由于有尾流和渦穿過,其聲信號是被污染的,計算時不能取用,但你可以用它來對比一下,研究它的影響。

計算的過程是比較簡單的。畫好網(wǎng)格之后,首先做穩(wěn)態(tài)計算,然后在穩(wěn)態(tài)計算的基礎(chǔ)上開始非穩(wěn)態(tài)計算。當非穩(wěn)態(tài)計算達到一個動態(tài)穩(wěn)定的狀態(tài)后,開啟氣動聲學(xué)模塊,設(shè)定聲源面,輸出聲源文件。這里有一個比較重要的問題。當你采用interior作為聲源面的時候,必須將流域設(shè)成多個計算域,其中積分面內(nèi)設(shè)成一個域,積分面外設(shè)成另一個域,如此你才能選取兩個域之間的interior面作為聲源面,并選取相應(yīng)的聲源計算域。聲源數(shù)據(jù)輸出后,指定聲接收點,然后使用積分法計算聲信號,計算就完成了。我的計算結(jié)果如下圖所示。

對比壁面和積分面的計算結(jié)果,低頻區(qū)基本一致,而在高頻區(qū),由于積分面法考慮了尾流中湍流的影響,因此相差較大。但是由于主要噪聲源是偶極子噪聲源,即低頻段,因此兩種方法在總聲壓級相差并不大。而噪聲的指向性應(yīng)該如第二張圖中紫色曲線所示,像一個豬腰子,這就是偶極子聲源的指向特性,說明低速圓柱繞流主要的噪聲源是偶極子噪聲源。

好了,對于一個簡單的問題,上面的方法基本能夠滿足工程上的需要了。總結(jié),對于簡單問題,重點在于你的CFD計算結(jié)果,以及聲源面的選取了。接下來介紹復(fù)雜問題。

復(fù)雜問題

CAA真正有難度的在于復(fù)雜問題的聲傳播這一塊。聲傳播的計算有兩種方法,一種是頻域法,另一種是時域法。頻域法是計算聲學(xué)最常用的方法,如sysnoise、ACTRAN等它們采用聲學(xué)有限元的方法來求解頻域上的波動方程,來實現(xiàn)聲傳播的計算。而時域法則是真正有技術(shù)含量的工作了,它需要好的算法和高階格式。這兩種方法理論上應(yīng)當是等效的,但是在實際情況中還是應(yīng)該有所偏重的。

當對聲傳播產(chǎn)生影響的主要因素是結(jié)構(gòu)邊界時,這時候采用頻域法是比較適宜的。這主要是因為頻域法對邊界條件的處理比較方便。比如說,對與聲襯或者有消聲材料的情況,采用時域法時該邊界條件需要一個復(fù)雜的方程來描述,而對于頻域法來說,只需要給定阻抗即可。

當對聲傳播產(chǎn)生影響的主要因素是流場時,這時候采用時域法是比較合適的,頻域法往往都無能為力。比如說發(fā)動機尾流,存在著強烈的剪切層時。下面的介紹則針對頻域法與時域法分別展開

頻域法

對于比較喜歡現(xiàn)有軟件的同志們來說,有兩個軟件可以作為頻域法求解的工具: ACTRAN和LMS VirtualLab。這兩個軟件采用的都是頻域法,過程大致相同:首先由CFD計算結(jié)果轉(zhuǎn)換出聲源參數(shù),得到聲源的時域分布,然后對聲源參數(shù)做DFT變換,得到聲源的頻域分布,最后采用計算聲學(xué)的方法來求取聲源作用下的聲壓。兩個軟件有著很大的差異。ACTRAN的聲源取的是體聲源,VL取的是面聲源;ACTRAN采用聲學(xué)有限元計算近場聲壓分布,在采用聲學(xué)無限元對遠場插值,求取輻射噪聲,而VL則是直接采用邊界元法來求取聲壓分布?？陀^的說,ACTRAN在氣動聲學(xué)方面的計算能力要優(yōu)于VL。

下面給出的是我采用ACTRAN計算的一個分支管氣動噪聲結(jié)果。

首先,做非穩(wěn)態(tài)CFD計算,結(jié)果如下所示:

接下來,根據(jù)Lighthill方程,做聲源轉(zhuǎn)換,得到聲源的時域分布

對聲源做DFT變換,得到聲源的頻域分布:

最后采用聲學(xué)有限元法,求管內(nèi)聲壓分布

由管內(nèi)聲壓分布可以看出,這種方法能夠得出管內(nèi)聲傳播的高階模態(tài)。指定聲接受點,取各個頻率下的聲壓級,則可得到聲壓的頻譜曲線了。

這就是采用ACTRAN做管內(nèi)氣動噪聲的完全過程了。答案對不對不好說,但是過程肯定是完整的,呵呵。當然,ACTRAN同樣可以用來做自由場的聲學(xué)問題,下圖是我用ACTRAN做的一個噴流噪聲計算結(jié)果,用的是ACTRAN自帶的CFD結(jié)果和網(wǎng)格文件。

時域法

時域法的難度可想而知,反正我是沒本事去試水了。所以這方面的東西我只能做些基礎(chǔ)的介紹了。各位同仁如果想玩玩這個的話,就得多花點力氣了。

前面說過,傳統(tǒng)的CFD方法不能應(yīng)用在聲傳播計算的主要原因是控制方程和離散格式不給力,因此,要想做聲傳播的時域計算,就必須需要新的控制方程和高階格式了。

在聲傳播計算中用的比較多的控制方程主要有線性歐拉方程(LinearlisedEuler Equatuon,LEE)和聲擾動方程(AcousticPerturbation Equation,APE)。

LEE方程是對歐拉方程進行變換得出的。它將流場參數(shù)分成平均量與脈動量兩部分,代入到歐拉方程中,略掉高階小量后,就得到了脈動量的控制方程。所謂線性,指的就是脈動量的偏微分項前的系數(shù)不包含脈動量本身。在網(wǎng)上隨便搜一下LEE,你就能得到LEE方程的推導(dǎo)過程。在LEE方程中,脈動量前面的系數(shù)都是平均量的函數(shù),因此,方程就包含了流場對脈動量傳播的影響,也就是說反映了流場對聲傳播的影響。

同樣的,APE方程則是對N-S方程進行線性化后得出的,也有很多文獻可參考,不再多說。

有了控制方程,接下來要考慮的就是怎么去解這些方程了。而方程的解法則很多了,如有限差分法、有限體積法、間斷迦遼金法等等。每種算法都有很多文獻可查,這里就不多說了,一則是確實不怎么明白,二則是那玩意講起來也不是一會就能完事的。

解決了控制方程的問題后,接下來的就需要解決高階離散格式(四階及以上)的問題了。而高階離散格式則更多了,如高階Rungr-Kutta格式,Adams-Bashforth格式,DRP格式,ENO格式,等等。我的理解,離散格式的基礎(chǔ)基本上都是泰勒級數(shù),只不過采用的級數(shù)的階數(shù)不同,同時計算時取的節(jié)點數(shù)不同而已。離散格式相關(guān)的文獻也是多如牛毛。搞算法的那些人隨便改改就能給弄出種格式來,至于實不實用則另當別論了。受水平限制,對高階格式這玩意我是講不出多少東西了,所以同樣不多說了,大家針對感興趣的格式可以多看文獻。

而在聲學(xué)計算當中另一個重要的事情就是計算的邊界條件了。通常來說,我們比較關(guān)心噪聲對環(huán)境的影響,也就是噪聲污染。因此,求遠場的輻射噪聲往往是主要目的之一。在需要考慮聲傳播計算的情況下,我們的計算域不可能取的很大,最終還是需要取聲傳播計算域的邊界作為聲源面,然后采用積分法來求遠場的聲輻射。此時,我們必須考慮計算域邊界引起的偽反射,壁面反射聲波對聲信號的污染。此時我們就需要引入無反射邊界條件。而能夠?qū)崿F(xiàn)無反射條件的邊界主要有特征邊界條件(Characteristic boundary conditions)、吸收層(absorbinglayers)、完全匹配層(Perfectly Matched Layers)等多種,具體的情概況則需要閱讀文獻了。

下面這張圖片給出了則是一個發(fā)動機尾流場中聲傳播的時域計算結(jié)果了。由圖中可以看出,由于剪切層的存在,聲的傳播方向發(fā)生了偏轉(zhuǎn)。

對于不想編程的同志來說,ACTRAN軟件的DGM模塊是你一個不錯的選擇。這個模塊是采用間斷迦遼金法DGM來求解LEE。下面一個圖片給出了采用ACTRAN的DGM模塊來求解發(fā)動機尾流中聲傳播的模型設(shè)定和計算結(jié)果。需要說明的是,ACTRAN的試用版軟件里面沒有DGM這個模塊,而破解版,貌似現(xiàn)在還沒有。要是哪位大哥手上有的話,記得多和我交流一下!!