張樹海, 李 虎, 王益民

(中國空氣動(dòng)力研究與發(fā)展中心 空氣動(dòng)力學(xué)國家重點(diǎn)實(shí)驗(yàn)室, 四川 綿陽 621000)

0 引 言

激波和旋渦是可壓縮流動(dòng)的基本結(jié)構(gòu)。在許多流動(dòng)中，比如超聲速混合層、超聲速射流和燃燒等，激波與旋渦共存，它們之間會(huì)發(fā)生相互作用，如激波與激波之間、激波與旋渦之間以及旋渦與旋渦之間的相互作用，這種作用，不僅會(huì)引起激波和旋渦的變形，還會(huì)產(chǎn)生噪聲。因此，吸引了很多學(xué)者開展理論、實(shí)驗(yàn)和數(shù)值模擬研究。

1955年，Hollingsworth 和Richards[1]首先研究了激波與旋渦的相互作用，他們利用激波管產(chǎn)生平面激波，該激波掃過一個(gè)具有一定攻角的機(jī)翼后，產(chǎn)生一個(gè)旋渦，激波碰到激波管端部的固壁后反射回來，與旋渦發(fā)生干擾。該實(shí)驗(yàn)發(fā)現(xiàn)，激波與旋渦干擾后，伴隨著激波與旋渦變形會(huì)產(chǎn)生聲波。Dosanjh 和Week[2]對此問題進(jìn)行了較細(xì)致的實(shí)驗(yàn)，測量了第一道聲波的強(qiáng)度，成為后來實(shí)驗(yàn)研究和數(shù)值模擬研究的一個(gè)標(biāo)準(zhǔn)。Barbosa 和Skews[3]設(shè)計(jì)了一個(gè)分叉型激波管，激波管所產(chǎn)生的激波被分成兩支，一支通過激波管內(nèi)的拐角產(chǎn)生旋渦，另一支直接與旋渦發(fā)生干擾，研究發(fā)現(xiàn)，激波與旋渦干擾后，在旋渦中心附近一個(gè)很小的區(qū)域，出現(xiàn)壓力脈沖區(qū)，其脈沖壓力是周圍壓力的兩倍以上。

為了解釋聲波的產(chǎn)生機(jī)理，早在20世紀(jì)60年代，Ribner[4]就發(fā)展了線性分析理論。他利用傅立葉變換將旋渦分解成平面正弦剪切波，這種單正弦波與激波相互作用后產(chǎn)生聲波，經(jīng)過合成后得到的激波與旋渦干擾產(chǎn)生的聲波在垂直于激波的方向反對稱分

布，其理論結(jié)果與Dosanjh和Week[2]的測試結(jié)果的反對稱部分一致。后來，Week和Dosanjh[5]發(fā)展了Lightill[6]的聲比擬理論，將周向壓力分布展開成四極、兩極和單極聲源之和。單極的作用使得聲波偏離了單純反對稱態(tài)，他們所得到的聲波與其實(shí)驗(yàn)一致。

到20世紀(jì)90年代，伴隨著高階精度捕捉激波格式的成熟，采用數(shù)值模擬研究激波和旋渦干擾并捕捉聲波已成為研究聲波產(chǎn)生機(jī)理的重要手段。Ellzey等[7-9]采用四階差分格式，通過求解非定?？蓧嚎sEuler方程，系統(tǒng)地模擬了激波與旋渦的干擾，得到的第一道聲波強(qiáng)度與實(shí)驗(yàn)測得的聲波強(qiáng)度基本一致。Inoue等[10-12]利用六階緊致格式直接模擬了不同強(qiáng)度激波與旋渦相互作用及聲波的產(chǎn)生過程，給出了弱激波與旋渦相互作用聲波的產(chǎn)生及傳播機(jī)理。他們通過精細(xì)的模擬發(fā)現(xiàn)，當(dāng)激波碰到旋渦的瞬間，旋渦邊緣發(fā)生局部壓縮和膨脹，形成雙極的第一道聲波，隨著激波與旋渦干擾的繼續(xù)，出現(xiàn)第二組壓縮和膨脹區(qū)，此壓縮和膨脹區(qū)與前面出現(xiàn)的壓縮和膨脹區(qū)域交替出現(xiàn)，形成第一道完整的聲波，第一道聲波演變成四極波。激波通過旋渦后，在旋渦兩側(cè)，形成兩道反射激波，反射激波后伴隨第二道和第三道聲波。Grasso和Pirrozli[13-16]也系統(tǒng)地模擬了較大范圍激波與旋渦的干擾問題，并根據(jù)激波變形的形狀，將激波與旋渦的干擾分成三類：弱干擾，具有規(guī)則反射的強(qiáng)干擾,具有馬赫反射的強(qiáng)干擾。

除了激波與旋渦之間的相互作用外，在包含激波和旋渦的復(fù)雜流場中，還存在旋渦與旋渦之間的相互作用，兩個(gè)或多個(gè)旋渦之間的相互作用也是噪聲產(chǎn)生的重要機(jī)制。比如，兩個(gè)同向旋轉(zhuǎn)的旋渦的合并過程，是剪切層噪聲的重要部分。

這種包含激波、旋渦的復(fù)雜流場的非定常運(yùn)動(dòng)產(chǎn)生的噪聲是激波噪聲。目前，雖然氣動(dòng)聲學(xué)理論已有相當(dāng)?shù)陌l(fā)展，但是關(guān)于激波噪聲的產(chǎn)生機(jī)制認(rèn)識(shí)還很不夠。伴隨高階精度數(shù)值格式和大型并行計(jì)算機(jī)的快速發(fā)展，數(shù)值模擬已成為研究激波噪聲產(chǎn)生機(jī)制的重要手段，但是激波噪聲的計(jì)算仍然是計(jì)算氣動(dòng)聲學(xué)所面臨的相當(dāng)嚴(yán)峻的挑戰(zhàn)，具體表現(xiàn)在：

(1) 流場的強(qiáng)非定常特性。激波噪聲是激波與其它流場結(jié)構(gòu)相互作用產(chǎn)生的，其中存在激波變形和產(chǎn)生過程，與激波相互作用的流場結(jié)構(gòu)也隨之發(fā)生變化。比如，在激波與旋渦的相互作用過程中[17-20]，入射激波會(huì)變形并產(chǎn)生反射激波，在旋渦中心附近還產(chǎn)生小激波串，同時(shí)旋渦被壓扁，而且出現(xiàn)扭轉(zhuǎn)運(yùn)動(dòng)，當(dāng)激波強(qiáng)度達(dá)到一定程度之后，會(huì)引起旋渦的破碎。在激波與剪切層相互作用過程中，存在有激波的局部震蕩[21]。因此，傳統(tǒng)的定常流體動(dòng)力學(xué)問題計(jì)算方法或是大尺度問題的非定常方法已不再適用。盡管已發(fā)展了多種高階精度Runge-Kutta方法[22]，但大多Runge-Kutta方法是顯式格式，由于穩(wěn)定性條件的限制，允許的時(shí)間步長非常小，對包含邊界層的流動(dòng)計(jì)算效率非常低，因此仍需要發(fā)展時(shí)間離散精度和效率更高的計(jì)算方法。

(2) 流場的多尺度特性。激波噪聲是流場的一部分，但是在很多情況下，它的擾動(dòng)尺度比流體動(dòng)力學(xué)尺度小很多[23]。比如在馬赫數(shù)為1.5的噴流中，距噴流直徑40倍處測得的噪聲強(qiáng)度是124 dB，聲場的脈動(dòng)速度與噴流的速度之比是1.5×10-4。特別地，當(dāng)流場中存在強(qiáng)激波時(shí)，聲波引起的物理量脈動(dòng)比激波引起的變化尺度差距更大。這就要求不僅能精確計(jì)算宏觀尺度的流體動(dòng)力學(xué)量，能夠捕捉激波，還能精確計(jì)算微弱尺度的聲波，對計(jì)算方法和計(jì)算網(wǎng)格的要求遠(yuǎn)比動(dòng)力學(xué)計(jì)算高。

(3) 模擬激波噪聲所需的計(jì)算區(qū)域遠(yuǎn)大于模擬氣動(dòng)力所需的計(jì)算區(qū)域。研究激波噪聲既關(guān)心近場流體動(dòng)力學(xué)的變化，挖掘噪聲產(chǎn)生的機(jī)制，又關(guān)心遠(yuǎn)場激波噪聲的傳播特性，計(jì)算區(qū)域邊界與聲源的距離非常遠(yuǎn)，這是與傳統(tǒng)計(jì)算流體動(dòng)力學(xué)的鮮明差異。

(4) 激波與聲波的計(jì)算對格式需求存在矛盾的因素。捕捉激波需要格式具有一定的耗散，而聲波是近乎無色散、無耗散的等熵運(yùn)動(dòng)，聲波的模擬需要高分辨率和無耗散格式。因此，現(xiàn)在還不存在一種較為理想的高階精度、高分辨率、低耗散的激波捕捉格式。目前，聲波的計(jì)算普遍采用線性緊致格式，精度和分辨率非常高，耗散性小，適合于聲波的計(jì)算。但是線性格式不能捕捉流場中的強(qiáng)激波，而捕捉激波普遍采用WENO格式[24]。盡管WENO格式的魯棒性非常好，可以模擬很強(qiáng)的激波，但同所有的捕捉激波格式一樣，在激波下游會(huì)產(chǎn)生微弱的非物理波動(dòng)。雖然這種波動(dòng)不影響流場的整體特性，但其強(qiáng)度常常會(huì)比流場中的噪聲高。另外，雖然可以設(shè)計(jì)很高精度的WENO格式[25]，但格式的耗散性比較大，對聲波的分辨率不滿意。盡管已有很多種簡單的混合格式[26]、優(yōu)化WENO格式[27]、非線性加權(quán)緊致格式[28-29]以及DRP格式[30]，試圖提高捕捉激波類格式的分辨率，但實(shí)際應(yīng)用比WENO格式并沒有明顯的優(yōu)勢，因此目前還不存在一種既能捕捉激波又適用于聲波計(jì)算的較理想數(shù)值方法。

針對這些問題，我們近年開展了較為系統(tǒng)性的研究，解決了計(jì)算激波噪聲的一些方法上存在的問題，對包含有激波、旋渦和聲波的一系列問題進(jìn)行了直接數(shù)值模擬，揭示了激波噪聲產(chǎn)生機(jī)制。本文對這些研究做簡要回顧。

1 針對激波噪聲計(jì)算的數(shù)值方法的研究

1.1 WENO格式的不收斂問題[31-32]

WENO格式是一種高階精度捕捉激波的格式[24],被廣泛用于工程實(shí)際和科學(xué)計(jì)算。但是像所有其它的激波捕捉格式一樣，用WENO格式計(jì)算定常激波時(shí)，在激波下游會(huì)出現(xiàn)微弱的非物理波動(dòng)，導(dǎo)致殘差不能收斂到滿意程度。雖然這種波動(dòng)不會(huì)影響整體的氣動(dòng)特性，但是至少會(huì)造成兩方面的影響，第一，殘差是結(jié)束計(jì)算的重要判據(jù)，殘差無法收斂到很小的值會(huì)使計(jì)算無法判斷何時(shí)停止，只能靠經(jīng)驗(yàn)結(jié)束進(jìn)程；第二，雖然激波后的非物理波動(dòng)是微弱的，不影響氣動(dòng)力等宏觀流動(dòng)特性，但是對于激波噪聲，這種微弱的非物理波動(dòng)則非常大，有時(shí)比激波噪聲的幅度高幾個(gè)數(shù)量級(jí)，因此計(jì)算氣動(dòng)噪聲問題，應(yīng)該消除這種微弱的非物理波動(dòng)。經(jīng)過系統(tǒng)的數(shù)學(xué)物理分析和大量的數(shù)值實(shí)驗(yàn)，我們提出了兩種方法消除這種微弱的非物理波動(dòng)[31-32]。

1) 針對五階WENO格式，我們提出了一種新的光滑測試因子如下：

(1)

具體形式為：

(2)

采用這種光滑測試因子后，消除了激波下游的微弱非物理波動(dòng)，對典型問題計(jì)算可以收斂到機(jī)器零。

2) 前面的方法只適合于五階精度的格式，對更高階精度格式不適用，而且會(huì)引起格式在高階極值點(diǎn)附近降階。經(jīng)過大量的數(shù)值實(shí)驗(yàn)，我們發(fā)現(xiàn)WENO格式在特征投影過程中，采用Roe平均方法計(jì)算半點(diǎn)上的物理量確定特征值和特征矩陣，是影響WENO格式不收斂的重要原因。為此，我們提出采用迎風(fēng)插值方法代替原來的Roe平均方法。即：

具體的迎風(fēng)插值公式如下：

一階插值：U(1)=UiU(2)=Ui+1

(3)

(4)

最優(yōu)加權(quán)插值：

(5)

(6)

最優(yōu)插值精度與WENO重構(gòu)過程精度一致，如五階精度的WENO格式采用五階插值，七階精度的WENO格式采用七階插值，九階精度的WENO格式采用九階插值。圖1和圖2是采用兩種方法對馬赫數(shù)2的一維定常激波計(jì)算得到的密度分布和收斂歷程，圖1中ResA為殘差的平均值。可以看到，改進(jìn)的WENO格式可以收斂到機(jī)器零。值得指出的是，采用最優(yōu)加權(quán)插值方法求解半點(diǎn)上的值，使WENO格式具有自相似特性，并一致高階精度，詳見文獻(xiàn)[31-32]。

1.2 一種類譜分辨率的高階精度緊致格式[34-35]

WENO格式是一種高階精度捕捉激波格式，雖然格式的精度可以很高，但是對聲波的分辨率并不理想，而且耗散性也比較大，采用WENO格式模擬湍流結(jié)構(gòu)、激波噪聲等復(fù)雜多尺度流動(dòng)并不理想。為了提高格式的分辨率，降低耗散，并保持捕捉激波的魯棒性，我們將WENO格式的設(shè)計(jì)思想應(yīng)用到緊致格式中，設(shè)計(jì)了一種類譜分辨率的高階精度緊致格式。這一格式是基于Lele的半點(diǎn)型線性緊致格式[33]設(shè)計(jì)的,Lele的半點(diǎn)型線性緊致格式形式是：

(7)

該格式左端是網(wǎng)格點(diǎn)上的函數(shù)導(dǎo)數(shù)值，右端是半點(diǎn)的函數(shù)值，它最大的優(yōu)點(diǎn)是分辨率非常高，接近譜方法分辨率。但是右端半點(diǎn)的函數(shù)值是未知的，Lele[33]設(shè)計(jì)了一種線性緊致插值的方法求得半點(diǎn)的值。圖3是Lele的半點(diǎn)型線性緊致格式的設(shè)計(jì)模板，可以看出，該格式有兩個(gè)缺點(diǎn)：(1)設(shè)計(jì)模板包括了整點(diǎn)和半點(diǎn)的函數(shù)值，但Lele僅用了部分信息，格式?jīng)]有達(dá)到該模板上的最優(yōu)精度；(2)半點(diǎn)上的值需要插值方法求得，任何插值方法都存在傳遞誤差，導(dǎo)致格式對高波數(shù)的分辨率大幅度下降。

為克服Lele格式的缺點(diǎn)，我們設(shè)計(jì)了一種類譜分辨率的高階精度緊致格式。在網(wǎng)格節(jié)點(diǎn)上，我們采用如下格式計(jì)算其一階導(dǎo)數(shù)：

(8)

采用同樣的格式計(jì)算半點(diǎn)上的一階導(dǎo)數(shù)值：

(9)

求得網(wǎng)格點(diǎn)和半點(diǎn)處的一階導(dǎo)數(shù)之后，通過求解相應(yīng)的控制方程，即可求得下一時(shí)刻的網(wǎng)格點(diǎn)和半點(diǎn)處函數(shù)值。

與Lele的格式相比，式(8)、式(9)給出的數(shù)值格式有兩個(gè)優(yōu)點(diǎn)：1) 格式的精度有一定的提高，Lele格式最高達(dá)到十階精度，我們設(shè)計(jì)的格式可以達(dá)到十四階；2) 網(wǎng)格點(diǎn)和半點(diǎn)上的函數(shù)值都是通過同一格式計(jì)算得到，消除了原格式緊致插值引起的傳遞誤差，保證了計(jì)算過程的精度和分辨率。圖4給出了六階和八階格式修正波數(shù)及其與 Lele 格式的比較，相比之下我們改進(jìn)的格式分辨率與精確值相差無幾，詳見文獻(xiàn)[34]。

線性格式不能捕捉流場中的強(qiáng)激波，為了使格式具有捕捉強(qiáng)激波的功能，我們將WENO重構(gòu)思想應(yīng)用到插值技術(shù)，發(fā)展了高階WENO插值方法，并將該方法應(yīng)用于式(8)、式(9)右端，初步發(fā)展了一種非線性類譜分辨率的緊致格式。具體的WENO插值方法是：在格式的設(shè)計(jì)模板S2r-1={xi-r+1,...,xi+r-1}上，采用插值函數(shù):

(10)

可以得到半點(diǎn)上的通量值：

(11)

與WENO格式重構(gòu)類似，將模板S2r-1分成三個(gè)子模板，在每個(gè)子模板上，存在r階插值：

(12)

采用加權(quán)思路，可以得到半點(diǎn)上的加權(quán)插值：

(13)

在式(8)的右端，采用如下混合插值方法替代原來的精確值fi+1/2，即：

(14)

而在式(9)中，采用類似的混合插值方法。得到的非線性格式，能夠保證格式的精度，具有較好的捕捉激波能力，同時(shí)格式的分辨率比現(xiàn)有的WENO明顯提高，耗散性也大幅度降低，詳見文獻(xiàn)[35]。

2 包含激波、旋渦和聲波的復(fù)雜多尺度流動(dòng)數(shù)值模擬

2.1 兩個(gè)同向旋轉(zhuǎn)的Gaussian渦合并過程[20]

Gaussian渦是一種典型的旋渦，在一定的條件下，兩個(gè)同向旋轉(zhuǎn)的Gaussian渦會(huì)合并[36]，對于合并條件，已有很多研究[37]。1995年，Mitchell等[36]采用直接數(shù)值模擬方法，研究了兩個(gè)同向旋轉(zhuǎn)的

(a)r=λ/2

(b)r=2λ

圖7兩個(gè)同向旋轉(zhuǎn)的Gaussian渦合并過程中
產(chǎn)生的聲波及其比較
Fig.7Far-fieldpressuretracesatr=λ/2andr=2λandthecomparisonbetweenourdirectnumericalsimulation(DNS)andthatbytheM?hring’sequation

值得指出的是，我們的結(jié)果與Mitchell的結(jié)果并不一致，其原因可能是由于初始參數(shù)的微弱不一致性，導(dǎo)致兩個(gè)旋渦合并時(shí)間的差異。為了驗(yàn)證我們結(jié)果的正確性，我們進(jìn)行了網(wǎng)格收斂性研究，初始條件的影響以及采用不同方法對此問題進(jìn)行了模擬，結(jié)果一致，詳見文獻(xiàn)[20]。

2.2 兩個(gè)Taylor旋渦的相互作用[20]

Taylor渦也是一種典型的旋渦，與Gaussian渦結(jié)構(gòu)不同，Taylor渦具有雙層結(jié)構(gòu)。為了深入了解渦聲的產(chǎn)生機(jī)理，我們系統(tǒng)研究了兩個(gè)Taylor旋渦的相互作用，并根據(jù)旋渦特性，將兩個(gè)Taylor旋渦的相互作用分成四類：(1) 兩個(gè)反向旋轉(zhuǎn)的Taylor渦的相互作用；(2) 兩個(gè)同向旋轉(zhuǎn)Taylor渦的相互作用；(3) 兩個(gè)同向旋轉(zhuǎn)Taylor渦的合并；(4) 兩個(gè)尺度或強(qiáng)度相差較大的旋渦的相互作用。

2.2.1 兩個(gè)反向旋轉(zhuǎn)的旋渦的相互作用

Taylor渦具有雙層結(jié)構(gòu)，內(nèi)層外層渦量符號(hào)相反。如果兩個(gè)反向旋轉(zhuǎn)的旋渦距離足夠近，它們之間會(huì)發(fā)生耦合，結(jié)果導(dǎo)致內(nèi)層外層分離，兩個(gè)渦核逐漸靠攏，形成一個(gè)渦偶極子，兩個(gè)Taylor渦的外層逐漸靠攏，形成一個(gè)反向運(yùn)動(dòng)的渦偶極子，如圖8所示,是強(qiáng)度為Mv=0.5、初始距離為4倍旋渦半徑的兩個(gè)反向渦相互作用的渦量演化歷程。詳見文獻(xiàn)[20]。

2.2.2 兩個(gè)同向旋轉(zhuǎn)Taylor旋渦的相互作用

圖9是兩個(gè)強(qiáng)度為Mv=0.5、初始旋渦間距是4倍旋渦半徑的兩個(gè)同向旋轉(zhuǎn)的Taylor旋渦相互作用過程的渦量演化歷程。這種相互作用存在兩種特性，1) 兩個(gè)同向旋轉(zhuǎn)的旋渦將演化成兩個(gè)不對稱的旋渦區(qū)，每一個(gè)渦區(qū)存在三個(gè)渦核結(jié)構(gòu)；2) 與兩個(gè)反向旋轉(zhuǎn)的旋渦相互作用具有本質(zhì)區(qū)別的是，每一個(gè)渦區(qū)的兩個(gè)渦核來自同一個(gè)初始旋渦，強(qiáng)的渦核來自初始旋渦的內(nèi)層，而初始旋渦的外層逐漸從初始旋渦中分離形成圍繞渦核的兩個(gè)獨(dú)立的渦核。

2.2.3 兩個(gè)同向旋轉(zhuǎn)Taylor旋渦的合并過程

與兩個(gè)同向旋轉(zhuǎn)的Gaussian旋渦相似，在一定條件下兩個(gè)同向旋轉(zhuǎn)的Taylor旋渦也會(huì)合并。如圖10所示，是兩個(gè)強(qiáng)度為Mv=0.5的同向旋轉(zhuǎn)的Taylor旋渦的合并過程，兩個(gè)旋渦的初始距離是旋渦半徑的2倍。兩個(gè)渦核合并形成主渦核，而在渦核的外部，形成兩個(gè)旋轉(zhuǎn)臂結(jié)構(gòu)，兩個(gè)渦的外層逐漸集中于旋轉(zhuǎn)臂附近，形成一個(gè)三核渦結(jié)構(gòu)。

2.2.4 兩個(gè)強(qiáng)度或尺度相差較大的旋渦的相互作用

兩個(gè)強(qiáng)度相差很大的旋渦的耦合過程是最令人感興趣的。1) 兩個(gè)強(qiáng)度或尺度相差較大的旋渦的相互作用會(huì)產(chǎn)生多核渦結(jié)構(gòu)；2) 強(qiáng)度或尺度的較大差異，導(dǎo)致它們之間相互作用所產(chǎn)生的多核渦結(jié)構(gòu)，具有非對稱特性。這種相互作用包括有兩種典型狀態(tài)，一是強(qiáng)度相差較大的兩個(gè)Taylor渦的相互作用，另一是尺度相差較大的Taylor渦的相互作用。

圖11是兩個(gè)強(qiáng)度相差較大的Taylor渦的相互作用過程，其中強(qiáng)旋渦的強(qiáng)度為Mvu=-0.8，順時(shí)針旋轉(zhuǎn)，弱旋渦的強(qiáng)度為Mvd=0.25，逆時(shí)針旋轉(zhuǎn)，它們之間的初始距離為半徑的4倍?？梢钥吹?，與兩個(gè)等強(qiáng)度的反向旋渦相互作用類似，兩個(gè)強(qiáng)度相差較大的Taylor渦相互作用會(huì)產(chǎn)生兩個(gè)渦偶極子。兩個(gè)旋渦的渦核逐漸靠攏，形成一個(gè)強(qiáng)的渦偶極子,兩個(gè)旋渦的外層逐漸從原旋渦中分離，形成一個(gè)弱的渦偶極子。由于其強(qiáng)度相差較大，兩個(gè)渦偶極子都呈現(xiàn)非對稱特性，導(dǎo)致弱的渦核圍繞強(qiáng)渦核旋轉(zhuǎn)，較弱的渦偶極子圍繞強(qiáng)渦偶極子旋轉(zhuǎn)。

圖12是兩個(gè)尺度相差較大的反向旋渦相互作用的過程，其中兩個(gè)旋渦的強(qiáng)度都是0.5，一個(gè)逆時(shí)針旋轉(zhuǎn)，一個(gè)順時(shí)針旋轉(zhuǎn)。兩個(gè)旋渦的尺寸相差較大，上面的旋渦半徑為1， 而下面的旋渦半徑為0.2， 兩個(gè)旋渦的初始距離為2.2。可以看到，小旋渦的渦核穿過大旋渦的外層，進(jìn)入到大旋渦的渦核外圍，并圍繞大旋渦的渦核不斷旋轉(zhuǎn)，最終形成一個(gè)橢圓形的渦核，而旋渦的外層形成兩個(gè)渦核結(jié)構(gòu)，圍繞大的渦核旋轉(zhuǎn)，形成一個(gè)三核渦結(jié)構(gòu)，并持續(xù)做蛙跳運(yùn)動(dòng)。

2.2.5 兩個(gè)Taylor渦相互作用過程中聲波產(chǎn)生機(jī)理分析

兩個(gè)Taylor渦的相互作用包括豐富的動(dòng)力學(xué)特性，這種過程會(huì)產(chǎn)生噪聲。圖13和圖14 是兩個(gè)典型狀態(tài)兩個(gè)Taylor渦相互作用的聲壓Δp=(p-p0)/p0等值線和其徑向分布。旋渦的強(qiáng)度均是0.5，兩個(gè)旋渦初始距離都是半徑的4倍。其中圖13是兩個(gè)反向Taylor渦相互作用所產(chǎn)生的噪聲，圖14是兩個(gè)同向Taylor渦相互作用所產(chǎn)生的噪聲。可以看到，兩種旋渦相互作用產(chǎn)生的噪聲具有很大區(qū)別，兩個(gè)反向旋渦相互作用僅產(chǎn)生幾個(gè)聲脈沖，而兩個(gè)同向旋轉(zhuǎn)的Taylor渦相互作用，與兩個(gè)同向旋轉(zhuǎn)的Gaussian渦相互作用很類似，會(huì)產(chǎn)生持續(xù)的強(qiáng)噪聲。

以兩個(gè)尺度相差較大的旋渦相互作用為例，分析噪聲產(chǎn)生機(jī)理。我們計(jì)算了Lamb量·(ρω×u)，它代表著噪聲產(chǎn)生的源。如圖15所示，是Lamb量的等值線演化過程?？梢钥吹?，旋渦相互作用具有三個(gè)明顯的運(yùn)動(dòng):1) 初始圓形旋渦被扭曲，弱旋渦被撕扯成兩個(gè)渦帶，他們與強(qiáng)旋渦繼續(xù)相互作用；2) 盡管強(qiáng)旋渦的渦心位置幾乎不變，但整個(gè)渦核被壓成橢圓形；3) 兩個(gè)渦帶與強(qiáng)旋渦的渦核之間相互作用產(chǎn)生一種蛙跳式運(yùn)動(dòng)。圖16是監(jiān)測點(diǎn)(x,y)=(100,0)和(0,100)的聲壓隨時(shí)間變化曲線，壓力峰值P1、P2、P3、P4、P5、P6和P7位于時(shí)間軸103、121、139、163、190、223和262?？紤]到聲從渦心附近產(chǎn)生聲波到達(dá)監(jiān)測點(diǎn)所需要的時(shí)間，它們對應(yīng)于長軸在y軸的方向，而聲波谷M1、M2、M3、M4、M5、M6和M7位于時(shí)間軸114、130、151、176、206、241和285，它們對應(yīng)于長軸位于x的方向。這意味著，噪聲主要由旋渦的蛙跳運(yùn)動(dòng)產(chǎn)生。

2.3 激波與旋渦的相互作用[17-19]

激波與旋渦相互作用是激波噪聲的一個(gè)簡化模型，通過研究激波與旋渦相互作用，可以深入了解激波噪聲產(chǎn)生機(jī)制。通過對激波與旋渦相互作用系統(tǒng)的模擬，我們發(fā)現(xiàn)激波與強(qiáng)旋渦相互作用存在多級(jí)干擾。第一級(jí)干擾是入射激波和初始旋渦的相互作用，第二級(jí)干擾是反射激波和變形旋渦的相互作用，第三級(jí)干擾是渦心附近產(chǎn)生的小激波串和變形旋渦的干擾,如圖17所示。與此同時(shí)，旋渦的變形也存在多級(jí)特征。在第一級(jí)干擾中，初始圓形旋渦被壓成橢圓形；在第二級(jí)干擾中，橢圓形旋渦被壓回成圓形旋渦；在第三級(jí)干擾中，圓形旋渦又被壓縮成橢圓形。此外，這種多級(jí)干擾會(huì)引起旋渦在水平線附近做振蕩運(yùn)動(dòng)。根據(jù)激波與強(qiáng)旋渦相互作用過程中的多級(jí)干擾過程，當(dāng)時(shí)我們預(yù)測，這種多級(jí)干擾過程會(huì)產(chǎn)生更為復(fù)雜的聲波。這一預(yù)測被Chatterjee和Vijayaraj[39]的研究所證實(shí)，如圖18所示。此外，我們還系統(tǒng)研究了激波與渦列相互作用過程流場結(jié)構(gòu)和聲波產(chǎn)生機(jī)制[17-19]。

2.4 激波與剪切層相互作用[40]

激波與剪切層相互作用是一個(gè)典型的問題，也可以看做是研究超聲速噴流激波噪聲的一個(gè)簡化模型。為了了解超聲速噴流激波噪聲的產(chǎn)生機(jī)制，我們設(shè)計(jì)了一個(gè)激波與剪切層相互作用的模型。如圖19所示，一道斜激波打到一個(gè)剪切層上，下邊界是一個(gè)反射固壁，透射激波打到固壁上后，形成反射激波，反射激波再次與剪切層相互作用發(fā)生反射，形成類似于超聲速噴流的激波串結(jié)構(gòu)。通過較為系統(tǒng)的數(shù)值模擬，我們揭示了兩種噪聲產(chǎn)生機(jī)制，一是激波與旋渦相互作用，發(fā)生在圖20所示的區(qū)域1， 另一是激波泄露機(jī)制，發(fā)生在圖20所示的區(qū)域2。

圖21是第一個(gè)區(qū)域入射激波與剪切層中旋渦相互作用過程的一個(gè)演化周期中幾個(gè)典型時(shí)刻的脹量圖，其中紅色虛線標(biāo)記的是剪切層下層流體中被追蹤的兩道聲波。我們看到由于流體以超聲速向下游傳播，聲波也隨著流體向下游方向?qū)α?，并且聲波的輻射半徑在逐漸增大。當(dāng)t=91.69時(shí)，產(chǎn)生了第二道聲波，可以明確地看到聲波的產(chǎn)生源在渦與小激波相互作用的位置，這和張樹海等[17-19]的激波與旋渦相互作用的研究結(jié)果相同。剪切層的渦列穿過激波，每個(gè)渦與激波發(fā)生干擾時(shí)，根據(jù)這種激波噪聲產(chǎn)生和傳播的機(jī)制，在噪聲源處向外輻射出一道道的聲波，同時(shí)我們可以看到當(dāng)聲波傳播到下壁面時(shí)聲波反射回內(nèi)場，從而聲波再次與激波和剪切層相互干擾。

圖22是第二個(gè)區(qū)域反射激波與剪切層相互作用過程中一個(gè)演化周期中幾個(gè)典型時(shí)刻的脹量圖，其中紅色虛線標(biāo)記的是剪切層下層流體中的一道聲波。我們看到由于流體以超聲速向下游傳播，聲波也隨著流體向下游方向?qū)α?，并且聲波的輻射半徑在逐漸增大。當(dāng)t=91.69時(shí)，噪聲源處產(chǎn)生了第二道聲波，可以明確地看到聲波的產(chǎn)生源在旋渦之間的鞍點(diǎn)處，即辮子區(qū)的位置，這和Manning等[41-43]的激波泄漏機(jī)制完全相同。在剪切層的渦列穿過激波時(shí)，渦對與激波發(fā)生干擾，在渦對之間的鞍點(diǎn)處激波透過剪切層，同時(shí)在鞍點(diǎn)處激波被辮子區(qū)的旋臂所阻礙，激波與辮子區(qū)的相互干擾使得激波強(qiáng)度變?nèi)?，部分激波穿過剪切層，還有一部分激波沒有穿過剪切層，而是在鞍點(diǎn)處以聲波的形式發(fā)生泄漏。這說明了強(qiáng)激波與剪切層相互作用中也存在激波泄漏機(jī)制，在鞍點(diǎn)處激波泄漏并向外輻射聲波。根據(jù)這種激波噪聲產(chǎn)生和傳播的機(jī)制，當(dāng)剪切層穿過激波時(shí)，在噪聲源處向外輻射出一道道的聲波，同時(shí)我們可以看到當(dāng)聲波傳播到下壁面時(shí)聲波反射回內(nèi)場，從而聲波再次與激波和剪切層相互干擾。

2.5 軸對稱超聲速噴流[44]

超聲速噴流噪聲是一種典型的激波噪聲，與很多工程密切相關(guān)。采用五階WENO格式直接求解軸對稱Navier-Stokes方程，數(shù)值模擬了射流馬赫數(shù)(完全膨脹)為Mj=1.19(嘯聲呈現(xiàn)軸對稱模態(tài)，三維效應(yīng)不明顯)的軸對稱欠膨脹射流，其對應(yīng)的射流聲速馬赫數(shù)為Ma=1.05，雷諾數(shù)為Re=6.216×105，噴管直徑為D=25.4 mm，噴管厚度為0.4D，環(huán)境溫度為288.15 K。

圖23給出了噴管唇口壁面[0.0,0.642D]和噴管出口平面[0.0,2.0D]上壓力監(jiān)測點(diǎn)處噪聲信號(hào)的頻譜信息，包括頻率和聲壓級(jí)。結(jié)果顯示，在大于5000 Hz的頻率區(qū)間內(nèi)，有4個(gè)尖銳的聲壓級(jí)突起：前兩個(gè)聲壓級(jí)突起的頻率為6567 Hz(128 dB)和8637 Hz(123 dB)，分別是嘯聲A1模態(tài)和嘯聲A2模態(tài)；另外兩個(gè)聲壓級(jí)突起的頻率為12087 Hz(125 dB)和14310 Hz(123 dB)，分別是嘯聲B模態(tài)的諧頻和嘯聲A0模態(tài)。

(a) 噴管唇口壁面

(b) 噴管出口平面

圖23監(jiān)測點(diǎn)處壓力信號(hào)的譜分析
Fig.23Soundpressurelevel

綜合流場結(jié)構(gòu)和聲場信息可以得到的嘯聲模態(tài)的圖像及其產(chǎn)生位置，其中流場結(jié)構(gòu)用數(shù)值紋影(密度梯度的模)表示，聲場由脹量場表示。根據(jù)嘯聲頻率可以得到嘯聲的波長，利用波長可以在脹量場中清晰地分辨出嘯聲A1模態(tài)和A2模態(tài)，相應(yīng)的波長分別為λ=2.04D和λ=1.51D。從圖24中可以清晰地分辨出向上游傳播的嘯聲，并且嘯聲產(chǎn)生位置都是在第三個(gè)激波柵格和第五個(gè)激波柵格之間的區(qū)域。

2.6 聲波穿過激波的畸變過程

聲波產(chǎn)生在近流場區(qū)域，流體介質(zhì)的非定常運(yùn)動(dòng)產(chǎn)生的噪聲經(jīng)過復(fù)雜流場結(jié)構(gòu)傳播達(dá)到遠(yuǎn)場。當(dāng)噪聲穿過復(fù)雜流場結(jié)構(gòu)過程中，會(huì)與流場相互作用，發(fā)生畸變。其中，聲波穿過激波是一個(gè)典型問題，我們采用直接數(shù)值模擬方法，研究了聲波穿過激波的畸變特性。

初始條件為：

(15)

(16)

(17)

計(jì)算區(qū)域?yàn)閤∈[-10,10]，馬赫數(shù)M=2.0，聲波幅值ε=1.0d-5，頻率ω=0.6π，比熱比取γ=1.4，計(jì)算結(jié)果如圖25所示，由圖可以看出，聲波穿過馬赫數(shù)M=2.0的激波，振幅放大了約3.64倍多。

2.7 聲波與旋渦相互作用

當(dāng)聲波穿過旋渦時(shí)，聲波的幅值和頻率都會(huì)發(fā)生變化。采用我們發(fā)展的線性緊致格式，通過求解Navier-Stokes方程，我們研究了平面聲波穿過等熵渦的散射特性。

圖26是計(jì)算模型，平面聲波與一個(gè)二維Taylor渦相互作用。計(jì)算區(qū)域的左端入射聲波為：

(18)

二維Taylor渦置于計(jì)算區(qū)域的中心。圖27是入射聲波頻率ω=0.6π與強(qiáng)度為Mv=0.25的Taylor渦相互作用的瞬態(tài)脈動(dòng)壓力場云圖。可以看到，旋渦對聲波的影響非常明顯。在旋渦下游，形成兩條強(qiáng)噪聲干擾條帶，在強(qiáng)條帶外側(cè)，有二次條帶，由于旋渦以逆時(shí)針方向旋轉(zhuǎn)，干擾條帶關(guān)于y=0不對稱。圖28是散射聲壓的均方根(RMS)等值線圖，圖29是離渦心r=10的圓周上散射聲壓的指向性分布。

3 結(jié)束語

近年來，針對激波噪聲計(jì)算的需求，我們開展了數(shù)值方法研究，通過對典型問題的直接數(shù)值模擬，揭示了激波噪聲產(chǎn)生機(jī)理。

(1) 在數(shù)值方法方面，針對高階精度捕捉激波的WENO格式不收斂問題，提出了一種新的光滑測試因子，發(fā)展了一種迎風(fēng)插值技術(shù)計(jì)算半點(diǎn)上物理量，用以確定WENO重構(gòu)過程中特征投影所需物理量。所改進(jìn)的WENO格式，消除了激波下游的微弱非物理波動(dòng)，對典型問題計(jì)算可以收斂到機(jī)器零。發(fā)展了一種類譜方法的緊致格式，這一格式包括線性和非線性兩部分，其中線性格式具有高階精度、高分辨率和低耗散特性，是模擬低速氣動(dòng)聲學(xué)問題的理想算法。非線性格式捕捉激波能力與WENO格式相當(dāng)，分辨率明顯提高，耗散性大幅度降低。

(2) 采用直接數(shù)值模擬方法，系統(tǒng)研究了兩個(gè)旋渦相互作用、激波與旋渦/渦列相互作用、激波與剪切層相互作用以及超聲速噴流激波噪聲。發(fā)現(xiàn)激波與強(qiáng)旋渦相互作用具有多級(jí)干擾特性，提出每級(jí)干擾都有自己的聲波產(chǎn)生。

當(dāng)然，對激波噪聲的認(rèn)識(shí)還相當(dāng)有限，有許多問題有待進(jìn)一步解決。

參 考 文 獻(xiàn)：

[1]Hollingsworth M A, Richards E J. On the sound generated by the interaction of a vortex and a shock wave[R]. British Aeronautical Research Council, Fluid Motion Sub-Committee, 18257, FM2371, 1956.

[2]Dosanjh D S, Weeks T M. Interaction of a starting vortex as well as a vortex street with a travelling shock wave[J]. AIAA J, 1965, (3): 216-223.

[3]Barbosa F J, Skews B W. Shock wave interaction with a spiral vortex[J]. Phys Fluids, 2001, 13: 3049-3060.

[4]Ribner H S. The sound generated by interaction of a single vortex with a shock wave[R]. UTIA Report No. 61, 1959.

[5]Weeks T M, Dosanjh D S. Sound generation by shock vortex interaction[J]. AIAA J, 1967, (5): 660-669.

[6]Lighthill M J. On sound generated aerodynamically I: general theory[J]. Proceedings of the Royal Society of London, 1952, 211: 564-587.

[7]Ellzey J L, Henneke M R, Picone J M, et al. The interaction of a shock with a vortex: shock distortion and the production of acoustic waves[J]. Phys Fluids, 1995, 7: 172-184.

[8]Ellzey J L, Henneke M R. The shock-vortex interaction: the origins of the acoustic wave[J]. Fluid Dyn Res, 1997, 21: 171-184.

[9]Ellzey J L, Henneke M R. The acoustic wave from a shock-vortex interaction: comparison between theory and computation[J]. Fluid Dyn Res, 2000, 27: 53-64.

[10]Inoue O, Hattori Y. Sound generation by shock-vortex interaction[J]. J Fluid Mech, 1999, 380: 81-116.

[11]Inoue O. Propagation of sound generated by weak-shock interaction[J]. Phys Fluids, 2000, 12: 1258-1261.

[12]Inoue O, Takahashi T, Hatakeyama N. Separation of reflected shock waves due to secondary interaction with vortices: Another mechanism of sound generation[J]. Phys Fluids, 2002, 14: 3733-3744.

[13]Grasso F, Pirozzoli S. Shock-wave thermal inhomogeneity week shock-vortex interaction in a mixing zone[J]. AIAA J, 1999, 33: 1797-1802.

[14]Grasso F, Pirozzoli S. Shock-wave-vortex interactions: shock and vortex deformations, and sound production[J]. Theroet Comput Fluid Dyn, 2000, 13: 421-456.

[15]Grasso F, Pirozzoli S. Simulations and analysis of the coupling process of compressible vortex pairs: Free evolution and shock induced coupling[J]. Phys Fluids, 2001, 13: 1343-1366.

[16]Pirozzoli S, Grasso F, D’Andrea A. Interaction of a shock wave with two counter-rotating vortices: Shock dynamics and sound production[J]. Phys Fluids, 2001, 13: 3460-3481.

[17]Zhang S, Zhang Y T, Shu C W. Multistage interaction of a shock wave and a strong vortex[J]. Phys Fluids, 2005, 17: 116101.

[18]Zhang S, Zhang Y T, Shu C W. Interaction of an oblique shock wave with a pair of parallelvortices: Shock dynamics and mechanism of sound generation[J]. Phys Fluids, 2006, 18(12): 126101.

[19]Zhang S, Jiang S, Yong-Tao Zhang, et al. The mechanism of sound generation in the interaction between a shock wave and two counter-rotating vortices[J]. Phys Fluids, 2009, 21: 076101.

[20]Zhang S, Li H, Liu X, et al. Classification and sound generation of two-dimensional interaction of two Taylor vortices[J]. Physics of Fluids, 2013, 25: 056103.

[21]Lui C V. C. M. A numerical investigation of shock-associated noise[D]. Stanford University, 2003.

[22]Gottlib S, Shu C W. Total variation diminishing Runge-Kuttaschemes[J]. Mathematics of Computation, 1998, 67: 73-85.

[23]Tim Colonius, SanjivaK Lele. Computational aeroacoustics: progress on nonlinear problems of sound generation[J]. Progress in Aerospace Sciences, 2004, 40: 345-416.

[24]Jiang G S, Shu C W. Efficient implementation of weighted ENO schemes[J]. Journal of Computational Physics, 1996, 126: 202-228.

[25]Balsara D S, Shu C W. Monotonicity preserving weighted essentially non-oscillatory schemes with increasingly high order accuracy[J]. Journal of Computational Physics, 2000, 160: 405-452.

[26]Pirozzoli S. Conservative hybrid compact-WENO schemes for shock-turbulence interaction[J]. J Comput Phys, 2002, 178: 81-117.

[27]Wang Z J, Chen R F. Optimized weighted essentially non-oscillatory schemes for linear waves with discontinuity[J]. J Comput Phys, 2001, 174: 381-404.

[28]Deng X, Zhang H. Developing high-order weighted compact nonlinear schemes[J]. J Comput Phys, 2000, 165: 22-44.

[29]Zhang S, Jiang S, Shu C W. Development of nonlinear weighted compact schemes with increasingly high order accuracy[J]. J Comput Phys, 2008, 227: 7294-7321.

[30]Tam C K W, Webb J C. Dispersion-Relation-Preserving finite difference schemes for computational acoustics[J]. J Comput Phys, 1993, 107: 262-281.

[31]Zhang S, Shu C W. A new smoothness indicator for the WENO schemes and its effect on the convergence to steady state solutions[J]. Journal of Scientific Computing, 2007, 31: 273-305.

[32]Zhang S, Jiang S, Shu C W. Improvement of convergence to steady state solutions of Euler equations with the WENO schemes[J]. Journal of Scientific Computing, 2011, 47: 216-238.

[33]Lele S K. Compact finite difference schemes with spectral-like resolution[J]. Journal of Computational Physics, 1992, 103: 16-42.

[34]Liu X, Zhang S, Zhang H, et al. A new class of central compact schemes with spectral-like resolution I: Linear schemes[J]. Journal of Computational Phys, 2013, 248: 235-256.

[35]Liu X, Zhang S, Zhang H, et al. A new class of central compact schemes with spectral-like resolution II: Hybrid weighted nonlinear schemes[J]. Journal of Computational Phys, 2015, 284: 133-154.

[36]Mitchell B E, Lele S K, Moin P. Direct computation of the sound from a compressible co-rotating vortex pair[J]. Journal of Fluid Mechanics, 1995, 285: 181-202.

[37]Cerretelli and Williamson C H K. The physical mechanism for vortex merging[J]. J Fluid Mech, 2003, 475: 41-77.

[38]Mohring W. On vortex sound at low Mach number[J]. J Fluid Mech, 1978, 85: 685-691.

[39]Chatterjee A, Vijayaraj S. Multiple sound generation in interaction of shock wave with strong vortex[J]. AIAA Journal, 2008, 46: 2558-2567.

[40]Liu X, Zhang S. Direct numerical simulation of the interaction of 2D shock wave and shear layer[J]. Chinese Journal of Theoretical and Applied Mechanics, 2013, 45: 61-75. (in Chinese)劉旭亮, 張樹海. 二維激波與剪切層相互作用的直接數(shù)值模擬研究[J]. 力學(xué)學(xué)報(bào), 2013, 45: 61-75.

[41]Manning T A. A numerical investigation of sound generation in supersonic jet screech[D]. Stanford University, 1999.

[42]Suzuki T, Lele S K. Shock leakage through an unsteady vortex-laden mixing layer: application to jet screech[J]. Journal of Fluid Mechanics, 2003, 490: 139-167[43]Lui C C M. A numerical investigation of shock-associated noise[D]. Stanford University, 2003.

[44]Li H. Numerical study of shock-associated noise in axisymmetric supersonic jet[D]. China Aerodynamics Research and Development Center, 2016. (in Chinese)李虎. 軸對稱超聲速射流激波噪聲數(shù)值模擬研究[D]. 中國空氣動(dòng)力研究與發(fā)展中心, 2016.