呂春靜 韓一平

1) (西安電子科技大學(xué)通信工程學(xué)院,西安 710071)

2) (西安電子科技大學(xué)物理與光電工程學(xué)院,西安 710071)

為了研究高斯光束在湍流等離子體鞘套中的傳輸特性,根據(jù)廣義惠更斯-菲涅耳原理,采用基于快速傅里葉變換的功率譜反演法,用多隨機相位屏來模擬湍流帶來的影響.根據(jù)超聲速飛行器繞流等離子體流場厚度在厘米級別的特點,光束在兩個相位屏之間的傳輸過程中采用菲涅耳衍射積分的兩次快速傅里葉變換算法(double fast Fourier transform algorithm),利用多隨機相位屏模擬等離子體鞘套湍流對光束傳輸產(chǎn)生的影響,解決了多隨機相位屏模擬湍流研究中的超短距離傳輸問題.當(dāng)飛行高度為45 km,飛行速度為18馬赫時,通過對超聲速飛行器繞流等離子體流場的統(tǒng)計分析,發(fā)現(xiàn)在此飛行條件下折射率起伏方差的強度范圍10-11—10-14.對高斯光束在湍流等離子體流場中的傳輸特性進行了數(shù)值仿真.結(jié)果表明：在等離子體鞘套湍流中折射率起伏強度、波長、傳輸距離等都是影響高斯光束質(zhì)量的重要因素.折射率方差越大,傳輸距離越長,光斑彌散越嚴(yán)重,光強起伏越大,光強減弱也越明顯.光束的波長越長,高斯光束抑制湍流的能力越強,光斑彌散程度越小,光強起伏也越小.

1 引 言

當(dāng)飛行器以高超聲速在大氣中飛行時,其動能會迅速轉(zhuǎn)化為大量的熱能導(dǎo)致飛行器周圍的空氣產(chǎn)生離解和電離,最終形成等離子體包覆的流場[1-4].等離子體的存在對于電磁波信號有很大的衰減作用,甚至導(dǎo)致信號中斷形成“黑障效應(yīng)”[5,6].目前已經(jīng)提出幾種有希望解決這一難題的幾種方法[7-11],例如通過提高電磁波頻率來減小信號衰減,或者在鞘套層中加入親電材料,這兩種方法的目的都是為了提高電磁波頻率與等離子體頻率的比值.同時有研究表明等離子體鞘套層是不穩(wěn)定的[12,13],所以湍流也是影響“黑障效應(yīng)”的重要因素.近些年,國內(nèi)外也早對湍流等離子體流場進行了較多的研究,例如Shi等[14]為了研究等離子體鞘套在無線信道上的動態(tài)效應(yīng),分析飛行條件變化(例如攻角、小尺度湍流)對介質(zhì)傳輸性能的影響,提出了一種新的自適應(yīng)馬爾可夫信道建模方法,利用該方法建立的信道的一階統(tǒng)計特性與無線電波的可變衰減具有良好的一致性,而且該方法具有良好的自適應(yīng)能力.He等[15]從通信角度采用透射矩陣法研究了再入等離子體鞘套對電磁波傳播的動態(tài)影響,所建立的模型中包含空間湍流和等離子體電子密度的變化,并對由等離子體的非線性特性引起的幅度和相位變化進行了數(shù)值仿真,發(fā)現(xiàn)幅度變化服從對數(shù)正態(tài)分布,相位變化服從瑞利分布.Yao等[16]針對等離子體流場的不穩(wěn)定性,首先在時-空-頻域揭示了電子密度的波動的規(guī)律,并基于這些規(guī)律提出了時變電子密度分層模型,利用Monte Carlo準(zhǔn)靜態(tài)電磁波數(shù)值方法計算了動態(tài)等離子鞘套層對Ka波段電磁波的影響.Lin和Sproul[17]重點研究了電子密度波動對電磁波傳播的影響,利用Navier-Stokes模型、三維拋物化Navier-Stokes模型、耦合邊界層和非黏性流動代碼估算了截頭錐體的電子密度的均值和方差.Li等[18]基于在高超聲速湍流實驗中測量的分維數(shù),建立了分形模型和等離子鞘套湍流功率譜,為研究等離子體鞘套湍流中的光束傳輸特性提供了重要參考.

在采用隨機相位屏研究湍流問題的文獻中,兩個隨機相位屏之間的傳輸過程大多采用單次快速傅里葉變換(single fast Fourier transform,S-FFT)算法[19-21],但是這種菲涅耳衍射積分的S-FFT算法適用于長距離傳輸,例如大氣湍流[22-24]和海洋湍流[25-27].飛行器周圍等離子體包覆的流場的厚度很小,在厘米或毫米級別[28,29].對于厘米和毫米級別的短距離傳輸問題,S-FFT算法會出現(xiàn)欠采樣的問題,最終將會導(dǎo)致信號失真.針對等離子體鞘套厚度較小這一特點,采用菲涅耳衍射積分的兩次快速傅里葉變換(double fast Fourier transform,D-FFT)算法研究光束在兩個相位屏之間的傳輸過程.本文第2部分首先在高度為45 km、飛行速度為18馬赫和0攻角的飛行條件,數(shù)值模擬了鈍頭體飛行目標(biāo)的高超音速等離子體鞘套湍流流場.然后利用Gladstone-Dale公式根據(jù)等離子體密度計算折射率,并對折射率的方差進行了估算;第3部分對比研究了菲涅耳衍射積分的D-FFT算法和S-FFT算法,說明D-FFT算法在厘米級別的短距離傳輸問題中的顯著優(yōu)勢,并利用多隨機相位屏模擬了高斯光束在等離子體鞘套湍流中的傳輸特性.

2 高超聲速鈍頭體等離子體鞘套湍流流場計算

2.1 鈍頭體飛行目標(biāo)繞流流場數(shù)值模擬

利用三維Navier-Stokes (N-S)方程和湍流兩方程k-ε模型數(shù)值模擬湍流流場.在笛卡爾坐標(biāo)系中,考慮1989年提出的雙溫?zé)崃W(xué)模型[30]和7組分空氣化學(xué)反應(yīng)模型(七組分化學(xué)反應(yīng)模型如表1所列),并對三維N-S方程進行適當(dāng)封閉,相應(yīng)模型如下[31-33]：

其中,Q和S分別是守恒變量和源項矢量;以x方向為例,E和Ev分別是對流通量矢量和黏性通量矢量;κ和κt分別是熱傳導(dǎo)系數(shù)和湍流熱傳導(dǎo)系數(shù);u,v和w分別是3個方向的瞬時速度分量;k是湍動能;是振動源項,是每個組分所占的質(zhì)量分數(shù);τij(ij=x,y,z) 為各個方向的剪應(yīng)力張量,需要利用湍流模型進行封閉[34];ev,i,hi和ρi分別是每一種組分的振動能量、顯焓和密度,下標(biāo)i代表每一種組分 (i=1,2,3,···,7);ρ,p和Ua分別為混合物的密度、壓強和平均速度;Eint,Et和Tint分別是混合物的內(nèi)部能量、湍流內(nèi)部能量和內(nèi)部溫度.

基于有限體積法,利用AUSMPW+(the advection upstream splitting method by pressurebased weight function)方法[35]求解N-S方程,其中時間推進采用LU-SGS (lower-upper symmetric Gauss-Seidel relaxation)方法[36].值得注意的是,文中采用的湍流k-ε模型是利用不同的方程來控制湍流能量傳輸并適合高雷諾數(shù)環(huán)境的兩方程模型,其中控制湍動能k和湍動能耗散率ε的兩方程分別為

表1 化學(xué)反應(yīng)模型中的方程Table 1.Reactions considered in chemistry model.

其中,μ,μt,μi分別是等離子體的動力學(xué)黏性系數(shù)、湍流黏性系數(shù)和質(zhì)量平均速度;常數(shù)Cε1=1.44,Cε2=1.92,σk=1.0,σε=1.3;P被定義為ui,uj和um為張量表示法中時均速度,δij是克羅內(nèi)克爾符(當(dāng)i=j時,δij=1,否則δij=0 ).

2.2 對流場數(shù)據(jù)做統(tǒng)計分析,估算等離子體折射率起伏方差

基于Gladstone-Dale公式[37],根據(jù)2.1節(jié)中計算的等離子體密度求解其折射率,Gladstone-Dale公式為

圖1(a)和圖1(b)分別為飛行目標(biāo)的區(qū)域劃分圖和五個區(qū)域的折射率統(tǒng)計方差圖.圖2(a)和圖2(b)分別為飛行目標(biāo)頭部和側(cè)部選取的一條傳輸路徑,頭部路徑三維坐標(biāo)范圍分別為：x=[0.0047 m,0.0170 m],y= [—0.0581 m,—0.0475 m],z= [0.0627 m,0.0512 m],側(cè)部路徑三維坐標(biāo)范圍分別為：x= [0.0674 m,0.0788 m],y= [—0.1555 m,—0.1225 m],z= [0.0628 m,0.0527 m].對兩條傳輸路徑上的折射率進行方差統(tǒng)計,統(tǒng)計方法為：以傳輸路徑中的每個三維流場坐標(biāo)點為中心,以超聲速湍流中典型的相關(guān)長度[38]—L/10為半徑(L為等離子體鞘套的厚度)建立球模型,對球內(nèi)所有數(shù)據(jù)進行統(tǒng)計,計算每個網(wǎng)格點處的折射率方差,流場邊界值采用外推法求解.通過對比圖2(a)和圖2(b)發(fā)現(xiàn),在接近飛行目標(biāo)壁面時,頭部折射率方差較大,比兩側(cè)的折射率方差要高出一個數(shù)量級,這是因為等離子鞘套中的湍流狀態(tài)還會受到激波、高溫、高壓等因素的影響,在靠近頭部的區(qū)域內(nèi)溫度、壓強較高且起伏較大.由圖1和圖2的計算結(jié)果顯示,在高度為45 km、飛行速度為18馬赫和0攻角的飛行條件下,飛行目標(biāo)的繞流流場中折射率起伏方差的強度的為10-11—10-14.

圖1 (a)五個區(qū)域;(b)五個區(qū)域的折射率方差Fig.1.(a) Five regions;(b) refractive index variance in five regions.

圖2 折射率方差 (a)頭部區(qū)域1;(b)側(cè)部區(qū)域2Fig.2.Refractive index variance in different regions：(a) Head area 1;(b) side area 2.

3 湍流等離子體鞘套中的隨機相位屏模擬

3.1 相鄰兩屏之間傳輸算法討論及仿真實驗驗證

S-FFT算法和D-FFT算法都是從菲涅耳衍射公式出發(fā),但對于衍射距離d的要求不同.SFFT算法主要適用于衍射距離較大的情況,而DFFT算法可適用于衍射距離較小的問題,下面對兩種算法進行理論分析和數(shù)值仿真,說明選用DFFT算法的原因.

3.1.1 菲涅耳衍射積分的S-FFT算法

在傍軸近似情況下,菲涅耳衍射積分可以表示為

其中,λ為光波波長,kλ=2π/λ為波數(shù),d是觀察面到衍射屏的距離,(x0,y0)和(x,y) 分別是光場在衍射面上和觀察面上某一點的復(fù)振幅分布.對(10)式進行適當(dāng)變換,并根據(jù)傅里葉變換法則,上式可改寫為

由于這種算法只需要進行一次傅里葉變換,所以稱為S-FFT算法.由(11)式可知其空間頻率為,即觀察屏上不同的點(x,y) 其實對應(yīng)不同的空間頻率(方向)的光.又因為光傳播方向與空間頻率的關(guān)系為fμ=1/X=cosα/λ,所以fμ=sin(π/2-α)/λ≈1/λ·(x/d) .

為了確定觀察面的尺寸,以x方向為例進行討論說明,y方向的情況與之類似.設(shè)衍射面在x方向上的寬度為Lx0,共有Nx個抽樣點,則按照抽樣定理,在用FFT變換計算時得到的空間頻率最高為fμmax=Nx/2Lx0,空間頻率為fμmax的這些衍射光對應(yīng)的傳播方向為fμmax=1/X=cosαmax/λ.這些衍射光出現(xiàn)在衍射觀察屏上坐標(biāo)的位置為fμmax=sin(π/2-α)/λ≈1/λ·(Lxmax/2/d) ,即衍射屏上的最大坐標(biāo)為Lxmax=Nxλd/Lx0.由此可見在使用S-FFT算法時,衍射觀察屏的大小不僅是波長的函數(shù),而且還是再現(xiàn)距離和抽樣數(shù)的函數(shù),當(dāng)d很小,如果保持取樣數(shù)不變,則觀察屏上的再現(xiàn)結(jié)果只是鄰近光軸的很小區(qū)域,所以該算法主要適用于衍射距離較大的情況.

協(xié)議書附圖是描述邊界線地理位置的法律圖件，是協(xié)議書的附圖，由鄉(xiāng)級行政區(qū)域界線信息與修測后的邊界地形圖疊加一起制作形成，與界線協(xié)議書具有同等的法律作用。

3.1.2 菲涅耳衍射積分D-FFT算法

同S-FFT算法相類似,在直角坐標(biāo)系中,傍軸近似下,當(dāng)衍射距離滿足菲涅耳衍射條件時,根據(jù)衍射積分表達式(10)可以得到其脈沖響應(yīng)：

顯然,菲涅耳衍射下的脈沖響應(yīng)是不變的,所以可以用卷積運算來完成,即

由于在空域完成卷積運算需要花費很長的機器時間,所以選擇到頻域完成.假設(shè)菲涅耳衍射的傳遞函數(shù)為(fμ,fv) ,則菲涅耳積分可以分為兩次傅里葉變換完成,所以稱其為D-FFT算法.可以通過傅里葉變換的性質(zhì)直接得到菲涅耳衍射傳遞函數(shù)(fμ,fv) 的解析表達式：

綜上所述,菲涅耳衍射公式可以這樣計算：

D-FFT算法中根據(jù)已經(jīng)確定的傳遞函數(shù)(fμ,fv)的解析表達式,只需要將頻域中x方向和y方向的頻率fμ和fv代入解析式,就可以得到傳遞函數(shù)的數(shù)值.假設(shè)實際問題是在寬度為Lx0×Ly0的空間采樣,并且有Nx×Ny個抽樣點,最高空間頻率分別為.所以,fμ的取值范圍從到取Nx個點,同樣,fv的取值范圍從到取Ny個點.

3.2 S-FFT算法和D-FFT算法的數(shù)值模擬結(jié)果的比較

在衍射計算過程中,應(yīng)該特別注意欠采樣的問題,在保證抽樣數(shù)和光束波長不變的情況下,改變衍射距離d,比較S-FFT算法和D-FFT算法抑制欠采樣發(fā)生的能力,說明D-FFT算法在處理短距離衍射積分問題的優(yōu)勢.

圖3 在不同衍射距離下S-FFT和D-FFT兩種算法的光強分布 (a) S-FFT,d = 100 mm;(b) D-FFT,d = 100 mm;(c) S-FFT,d = 10 mm;(d) D-FFT,d = 10 mmFig.3.Intensity distribution of S-FFT and D-FFT algorithm at different diffraction distances：(a) S-FFT,d = 100 mm;(b) D-FFT,d = 100 mm;(c) S-FFT,d = 10 mm;(d) D-FFT,d = 10 mm.

圖3(a)和圖3(b)分別為衍射距離為100 mm時利用S-FFT算法和D-FFT算法得到的衍射光強分布,圖3(c)和圖3(d)分別為衍射距離為10 mm時利用S-FFT算法和D-FFT算法得到的衍射光強分布.通過對比圖3(a)和圖3(b)發(fā)現(xiàn)：當(dāng)衍射距離為100 mm時,利用S-FFT和D-FFT兩種算法得到的衍射光強分布均未出現(xiàn)欠采樣現(xiàn)象.通過對比圖3(c)和圖3(d)可知：隨著衍射距離的減小,采用S-FFT算法得到的衍射光強出現(xiàn)了明顯的欠采樣現(xiàn)象,但是D-FFT算法未出現(xiàn),可見D-FFT算法的抗欠采樣能力要明顯高于S-FFT算法.結(jié)果表明,D-FFT算法不懼怕計算衍射距離近的問題,即使衍射距離d= 1 mm時,用D-FFT算法依然可以得到很好的結(jié)果.超高音速飛行器周圍形成的等離子體鞘套的厚度屬于厘米數(shù)量級.由此可見,用FFT算法數(shù)值模擬光束在相鄰兩屏之間的傳輸過程時,選用D-FFT算法更合適.

利用譜反演法生成等離子體鞘套湍流隨機相位屏,首先需要在頻域生成一個均值為1,方差為0的N×N復(fù)隨機矩陣H(kx,ky) ,利用等離子體鞘套湍流功率譜[18]對其進行濾波,然后再進行傅里葉逆變換便可以得到等離子體鞘套湍流的隨機相位屏Φ(kx,ky) .等離子體鞘套湍流功率譜為

其中為折射率方差,L0為湍流外尺度,kx和ky分別為不同方向的頻率.由折射率頻譜可得到相位頻譜為Φ(kx,ky)=2πk2Φn(kx,ky) ,則相位頻譜的方差為σ2(kx,ky)=(2πd/NΔx)2Φ(kx,ky) ,其中d是每兩個隨機相位屏的間距,N為采樣數(shù),Δx是柵格間距.通過傅里葉變化便可得到空域中的隨機相位屏[39]：Φ(x,y)=F[H(kx,ky)σ(kx,ky)] .

利用菲涅耳衍射積分的D-FFT算法求解高斯光束在相鄰兩屏之間的傳輸過程,利用功率譜反演法生成的多隨機相位屏模擬等離子鞘套湍流帶來的影響,多隨機相位屏示意圖如圖4所示.

圖4 多隨機相位屏示意圖Fig.4.Schematic diagram of multi-random phase screen.

圖5 不同傳輸條件時的高斯光束光強數(shù)值仿真 (a)λ = 1.55 × 10-6 m,ω0 = 40 mm,d = 13 mm;(b)ω0 = 40 mm,d =13 mm, = 10-12;(c)λ = 1.55 × 10-6 m,ω0 = 40 mm, = 10-12Fig.5.Numerical simulation of intensity of Gaussian beam in different propagation conditions：(a)λ = 1.55 × 10-6 m,ω0 = 40 mm,d = 13 mm;(b)ω0 = 40 mm,d = 13 mm, = 10-12;(c)λ = 1.55 × 10-6 m,ω0 = 40 mm, = 10-12.

3.3 等離子體鞘套湍流中影響高斯光束質(zhì)量的因子分析

圖5的傳輸參數(shù)為：束腰半徑ω0= 40 mm,相位屏間隔d= 0.013 m,相位屏采樣數(shù)為 1024×1024.分析折射率方差、波長和傳輸距離對高斯光束質(zhì)量的影響.

圖5(a)是高斯光束在折射率起伏方差分別為=10-13,10-12,10-11時通過10個隨機相位屏后系統(tǒng)歸一化的三維(3D)和二維(2D)光強數(shù)值仿真結(jié)果.由圖5(a)可知,在湍流折射率方差不斷增大的過程中,光斑分裂成形狀不規(guī)則的若干小光斑,光斑彌散程度加劇,而且光強起伏逐漸增大,強度逐漸減小.

圖5(b)是高斯光束波長分別為λ=0.6328×10-6,1.05×10-6,1.55×10-6m時通過10個隨機相位屏后系統(tǒng)歸一化的3D和2D光強數(shù)值仿真結(jié)果.通過對比分析可知,隨著波長的不斷增加,高斯光束抵抗湍流的能力逐漸增強,具體表現(xiàn)在光斑彌散程度減小、光強起伏減弱上,所以在等離子體鞘套通信系統(tǒng)中選用波長較長的激光光束作為信息載體比較好.

圖5(c)是高斯光束沿z軸傳輸時距離分別為6.5,13,19.5 cm,即分別通過5,10,15個隨機相位屏后系統(tǒng)歸一化的3D和2D光強數(shù)值仿真結(jié)果.結(jié)果反映出：傳輸距離越大,光斑彌散越嚴(yán)重,光強起伏也越明顯.

4 結(jié) 論

首先利用菲涅耳衍射積分的D-FFT算法解決了湍流問題中的短距離傳輸問題,然后利用多隨機相位屏數(shù)值模擬了高斯光束在等離子體鞘套湍流中的傳輸特性.研究表明折射率方差越大,高斯光束的光斑彌散程度越大,光強越弱,光強起伏越明顯;隨著波長的不斷增加,高斯光束抵抗湍流的能力逐漸增強,所以在等離子體鞘套通信系統(tǒng)中選用波長較長的激光光束作為信息載體比較好;傳輸距離越大,光斑彌散越嚴(yán)重,光強起伏也越明顯.綜上所述,高斯光束在等離子體鞘套湍流中的傳輸質(zhì)量與外界因素折射率起伏方差、傳輸距離及自身因素波長、束腰半徑等密切相關(guān).本文對研究激光在湍流等離子體鞘套中的傳輸問題具有一定的意義.