斜入射非線性電離層Langmuir 擾動的電磁波傳播特性*

楊利霞 劉超 李清亮 閆玉波

1) (安徽大學(xué)電子信息工程學(xué)院,合肥 230031)

2) (安徽大學(xué),信息材料與智能感知安徽省實驗室,合肥 230031)

3) (中國電波傳播研究所,電波環(huán)境特性及?；夹g(shù)重點實驗室,青島 266107)

1 引言

電離層作為一種時變的非均勻磁化等離子體,會與入射到其中的電磁波發(fā)生復(fù)雜的相互作用[1].一般地,利用高功率尋常波(O 波)對電離層進(jìn)行加熱實驗,由于參量衰變不穩(wěn)定性,泵波可在反射高度附近激發(fā)出電子Langmuir 波和離子聲波[2],同時波粒相互作用導(dǎo)致O 波轉(zhuǎn)換為Z 波并向電離層更高區(qū)域傳播,進(jìn)一步形成Langmuir 擾動[3],Langmuir 擾動的形成可引起非相干散射雷達(dá)背向散射功率增強(qiáng),即離子線與等離子體線的增強(qiáng)現(xiàn)象.

電離層加熱是近年來電波傳播領(lǐng)域的熱點問題,電離層局部加熱性導(dǎo)致等離子體的非線性變化[4],可以極大地改變局部地區(qū)的電磁波傳播特性,具有重要的研究價值和軍事研究意義.而Langmuir擾動是電離層加熱效應(yīng)的重要研究方向之一[5,6],國外研究人員做了大量的Zakharov 模擬實驗來研究電離層中強(qiáng)擾動和弱擾動的特性,并將得到的能量譜與觀察到的受激電磁輻射以及雷達(dá)觀測結(jié)果進(jìn)行比較[7,8].Zakharov 模型對時間尺度進(jìn)行了劃分,得到高頻電磁場的復(fù)數(shù)包絡(luò)和等離子體低頻響應(yīng)的方程[9].Mj?lhus 等[10]建立了一個兩級模型來表示Langmuir 擾動和電磁波之間的耦合,對Langmuir 波使用Zakharov 模型,對電磁波使用格林函數(shù)法.Eliasson 等[11,12]進(jìn)行了一系列的研究,對電離層Langmuir 擾動建立了一維(1D)廣義Zakharov 數(shù)值模型.目前國內(nèi)對Langmuir 擾動本身的研究尚處于起步階段,如武漢大學(xué)的劉默然等[13]、何昉等[14]以及江蘇大學(xué)的朱婷[15]都對電離層垂直加熱下Langmuir 的擾動進(jìn)行了數(shù)值模擬研究.

對于斜入射[16]電離層加熱下Langmuir 的擾動目前國內(nèi)研究較少,為此,本文基于廣義Zakharov模型,通過將二維(2D)麥克斯韋方程等價地轉(zhuǎn)換為一維麥克斯韋方程[17,18],避免了用傳統(tǒng)的二維時域有限差分方法分析該散射問題,這樣將節(jié)省大量內(nèi)存,極大地提高計算效率,有效降低編程的復(fù)雜程度.同時,深入細(xì)致地研究計算背景為等離子體的斜入射問題,從而進(jìn)一步研究基于電離層非線性等離子體擾動的電磁波傳播特性的全尺度時域計算方法,并基于上述方法分析電離層擾動的電磁場分布和傳播特性.

2 斜入射非線性等離子體下Langmuir擾動分析

2.1 斜入射非磁化等離子體的迭代式推導(dǎo)

TEz斜入射非磁化等離子體的示意圖如圖1所示,其中介質(zhì)區(qū)間1 和3 表示真空,介質(zhì)區(qū)間2表示電離層等離子體.這里選取電磁波的波數(shù)k與z軸夾角為θ,T Ez波在區(qū)間1 中以θ角度斜入射到等離子體區(qū)域中.

圖1 TEz 波斜入射電離層等離子體的傳播模型Fig.1.Propagation model of TEz wave obliquely incident ionospheric plasma.

由圖1 可知,在斜入射情況下,電磁波在電離層中的傳播過程是一個二維問題.從宏觀角度來看,將麥克斯韋旋度方程在直角坐標(biāo)系下展開,可得 T Ez波在區(qū)間1 的時域麥克斯韋方程:

根據(jù)文獻(xiàn)[19]中的推導(dǎo)方法,最終(1)式可化為時域里等價的一維形式:

式中k1為區(qū)間1 (真空)中的電磁波波數(shù),k1z為k1在z軸方向上的投影,即k1z=k1cosθ.

同理,在區(qū)間2 (等離子體)中,一維斜入射情況下電磁波同樣滿足(2)式形式的一維解,并轉(zhuǎn)換到頻域:

式中k2為區(qū)間2 (電離層等離子體)中的電磁波波數(shù),k2z為k2在z軸方向上的投影,即k2z=k2cos?.

為了解決區(qū)間1,3 與區(qū)間2 的電磁波傳播過程中遇到的相位及邊界問題,由邊界條件和相位匹配原理知將其代入(3)式得

對方程(4)中的第2 式進(jìn)行整理,可得

式中磁化等離子體的相對介電常數(shù)為矩陣形式,可寫為

其中σ表示電導(dǎo)率,將(6)式代入(5)式并通過逆傅里葉變換到時域,可得

在等離子體介質(zhì)中電流密度為

其中,下標(biāo)“ i ”表示離子,下標(biāo)“e”表示電子.

由于離子的質(zhì)量遠(yuǎn)大于電子的質(zhì)量,離子成份的流體速度遠(yuǎn)小于電子成份,離子可以被視為是靜止的.忽略離子部分后將(8)式代入(7)式,最后聯(lián)立(4)式的第一個式子并進(jìn)行逆傅里葉變換,可得

同理推導(dǎo)Ey1D,Hx1D具有相似的形式:

電離層等離子體可被看成是由電子流體和離子流體組成的導(dǎo)電的連續(xù)介質(zhì),因此描述等離子體運(yùn)動采用流體力學(xué)方程.在雙流體力學(xué)描述中,帶電粒子的數(shù)密度和流體速度分別滿足連續(xù)性方程和動量方程:

式中nα為粒子數(shù)密度,υα為粒子的流體速度,qα表示粒子電荷;磁感應(yīng)強(qiáng)度B=B0+B1,其中B0是地磁場,B1是波磁場;mα是粒子質(zhì)量;Pα是壓強(qiáng);γα是與粒子間碰撞等因素有關(guān)的阻尼算子.

將(11)式中的電子數(shù)密度連續(xù)性方程,經(jīng)過高斯定理推導(dǎo)出如下形式:

對于高頻波,(11)式中的電子動量方程變?yōu)?/p>

(9)式、(10)式、(12)式和(13)式構(gòu)成了高頻電磁波斜入射電離層等離子體的非線性耦合方程組,也是此研究的理論基礎(chǔ).

時域有限差分方法可以直接求出Maxwell 方程組顯式下的離散形式,(9)式和(10)式進(jìn)行離散后得到磁場強(qiáng)度的迭代方程為

從迭代式(14)—(16)可以看出,電場的垂直分量E⊥與磁場的垂直分量以及等離子體密度、流體速度密切相關(guān),而電場的平行分量Ez僅與其本身、電子的數(shù)密度和z方向的流體速度有關(guān).

對(12)式關(guān)于空間的一階偏導(dǎo)數(shù)選取后向差分近似,得到電子密度的離散迭代方程為

對(13)式等號右側(cè)的第一項作近似處理,再通過矩陣求解得到電子流體速度的離散迭代方程為

電磁波在電離層中斜入射傳播過程的時域數(shù)值完整迭代步驟如圖2 所示.

圖2 FDTD 電磁場時域的逐步迭代流程圖Fig.2.Step by step iterative flow chart of FDTD electromagnetic field in time domain.

2.2 數(shù)值算法的穩(wěn)定性條件

1)時間離散間隔的穩(wěn)定性要求

在使用FDTD 算法對電磁波在電離層中的傳播進(jìn)行數(shù)值仿真時,需要考慮等離子體與電磁波相互作用的影響.等離子體中空間和時間離散間隔之間應(yīng)當(dāng)滿足的充分條件[20]為

其中wpe為等離子體特征頻率.(19)式表明一維情況下,時間間隔必須等于或小于波以光速通過Yee元胞所需的時間.

2)數(shù)值色散對空間離散間隔的穩(wěn)定性要求

有限差分代替波動方程的二階導(dǎo)數(shù)過程中,離散處理將會導(dǎo)致波的色散,此時需要對空間離散間隔做出要求來降低色散,一維情況下要求

式中λ對應(yīng)介質(zhì)中的最小波長.

本次數(shù)值模擬選取λ=60 m ,Δz=1 m,并使Δt滿足Courant 條件,本文在50—400 km 全范圍內(nèi)以1 m 步進(jìn)計算.

2.3 數(shù)值模型

本文建立的電磁波斜入射電離層的數(shù)學(xué)模型如圖3.數(shù)值模型采取的實際參數(shù)來自EISCAT 在挪威的電離層加熱實驗[21,22],考慮到參量不穩(wěn)定性發(fā)生的反射區(qū)域高度范圍一般在200—340 km,假設(shè)沿z軸垂直向上為正方向,水平方向正北為y軸正方向,水平向西為x軸正方向.計算區(qū)域的范圍在垂直方向延伸400 km.為了降低計算量對計算機(jī)內(nèi)存的要求,本文的數(shù)值模擬將電離層剖面的尺寸壓縮為實際尺寸的1/20.

圖3 電磁波加熱電離層傳播模型Fig.3.Ionospheric propagation model heated by electromagnetic wave.

2.4 數(shù)值仿真參量

在z=50 km 處通過連接邊界,采用線極化的調(diào)制高斯脈沖作為入射源,Ex=1.5 sin(ω0t)exp[-(z-50×103)/108]V/m.電磁波在電離層的傳播過程中發(fā)生模式轉(zhuǎn)換進(jìn)而產(chǎn)生Langmuir 擾動現(xiàn)象,必須滿足兩個條件:1)入射波頻率應(yīng)當(dāng)略小于 F2層的最大臨界頻率[23,24];2)在臨界層內(nèi)存在等離子體密度或者電場強(qiáng)度的擾動[25].根據(jù)所給參數(shù),計算出 F2層最大截止頻率fmax≈6.33 MHz,這里選取加熱電磁波頻率f0=5 MHz 以滿足第1 個條件.第2 個條件是電磁波對電離層加熱引起了參量的不穩(wěn)定性,會使等離子體中存在靜電波擾動.

在模擬計算區(qū)域中,上下邊界的截斷處采用一維的Mur 吸收邊界.根據(jù)吸收邊界原理和連接邊界理論[26],在斜入射情況下空間步長改為 Δz在z方向上的投影,具體修正為如下形式:

左截斷邊界處

TEz波只需在原有的連接邊界形式上乘以角度因子 c osθ,具體修正為如下形式:

2.5 O 波、X 波和Z 波

對于高頻電磁波,考慮垂直傳播的波,假設(shè)電磁波傳播方向k沿E1×B1方向,因此電場可能有兩種基本方向(如圖4 所示):E1//B0和E1⊥B0.當(dāng)E1平行于磁場B0時,這種波被稱為尋常(O)波[27],得到O 波的色散關(guān)系為結(jié)果表明尋常波的傳播不受磁場影響.對于E1垂直于磁場B0的情況,電子響應(yīng)隨電場而運(yùn)動,但洛倫茲力的作用使電子不能沿固定方向運(yùn)動,而要產(chǎn)生另一方向的速度分量.為了區(qū)別于和B0無關(guān)的尋常波,把和B0有關(guān)的波動稱為非尋常(X)波,非尋常波是由部分橫波和部分縱波構(gòu)成的電磁波.非尋常波的色散方程為

圖4 O 波與X 波Fig.4.O wave and X wave.

其中n為折射率,等離子體對于擾動的響應(yīng)是一個頻率為的振蕩,該頻率稱為高混雜頻率.

Z 波模式是特殊(X)模式的低頻分支的空間物理符號,它是等離子體的一個內(nèi)?；蚍@模,其頻率限定在右旋截止頻率WL和高混雜頻率WUH之間[28].

3 模擬結(jié)果

3.1 截止和共振

為了分析波的傳播特性,引入兩個重要的概念:截止和共振[29].當(dāng)波的折射率變?yōu)榱慊虿〝?shù)趨于零時,在等離子體中出現(xiàn)波的截止.由O 波的色散關(guān)系式知,波數(shù)趨于0 時其截止頻率等于電子等離子體頻率wpe,當(dāng)w＞wpe時,O 波可以在等離子體中傳播.由(24)式得,X 波在n2=0 ,k=0 情況下截止,因為n2＜0 使得波傳播因子變?yōu)檎穹p因子,意味著波在介質(zhì)中傳播時很快衰減,最終被截 止.而在n2→∞,k→∞條件下w與k無關(guān),這樣相速度、群速度都為0,波不能傳播,出現(xiàn)共振.在截止點和共振點,波都不能傳播,一般來說波在截止點被反射,在共振點被吸收.

令n2=0 求解(24)式,w方程應(yīng)該有4 個根,其中只有兩個根是合理的(w＞ 0),它們的截止頻率分別稱為左旋、右旋截止頻率,如圖5 所示.

圖5 w-k 色散曲線Fig.5.w-k dispersion curve.

通過對O 波和X 波的色散分析可知,對于O 波,當(dāng) 0＜w ＜wpe時不能傳播,w＞wpe時才能傳播.對于X 波,在 0＜w ＜wL和wUH＜w ＜wR兩 個頻率區(qū)域不能傳播,不能傳播的頻區(qū)稱為截止帶;當(dāng)wL＜w ＜wUH和w＞wR時,X 波可以傳播,對應(yīng)的頻區(qū)稱為傳播帶.

3.2 斜入射下的電磁波傳播特性

通過數(shù)值仿真模擬斜入射角為7°時的電離層傳播情況,結(jié)果如圖6 所示.本文數(shù)值模型的計算是在激勵源(50 km,圖6 粉色圓點處)加熱一段時間后開始進(jìn)行,電場的初始值在下邊界是一個沿著x方向線性極化的波,在進(jìn)入電離層等離子體區(qū)域前先經(jīng)過一段中性大氣(50—85 km,紅色虛線區(qū)域),該傳播過程中包絡(luò)形狀不變.

圖6 背景電子密度及t=0.249 ms 時刻電場 E x 分量幅值分布Fig.6.Background electron density and amplitude distribution of electric field E x component at t=0.249 ms.

入射波在電離層等離子體中傳播一段距離后出現(xiàn)了y方向上的分量,如圖7 所示,這是由法拉第旋轉(zhuǎn)現(xiàn)象引起的,并且電場y分量的包絡(luò)形狀與入射波基本一致.

圖7 背景電 子密度 及t=0.821 ms 時刻電 場 E x 和 E y 分量幅值分布Fig.7.Background electron density and distribution of electric field E x and E y component amplitude at t=0.821 ms.

如圖8 所示,隨著等離子體密度的增加,電磁波的波形開始發(fā)生變化,X 波反射點大約在271 km處(紅點所示),電磁波逐步接近X 波反射點的過程中,會導(dǎo)致波包出現(xiàn)分界點.

圖8 背景電 子密度 及t=0.90 ms 時刻電 場 E x 和 E y 分量幅值分布Fig.8.Background electron density and distribution of electric field E x and E y component amplitude at t=0.90 ms.

在達(dá)到X 波反射點后還有O 波繼續(xù)向前傳播,O 波在接近其反射點的過程中逐漸從左旋圓極化變成沿地磁場方向的線極化,圖9 所示的O 波反射點附近激發(fā)出了明顯的z向局部振蕩的電場分量.從圖9 可以明顯看出,大約在278 km 附近有明顯的截止特征,該處為O 波的反射點,與根據(jù)Appleton-Hartree 公式計算所得的O 波反射高度相符.

圖9 背景電 子密度 及t=1.03 ms 時刻電 場 E x ,E y 和Ez分量幅值分布Fig.9.Background electron density and distribution of electric field E x ,E y and E z component amplitude at t=1.03 ms.

由圖10 可以看出,z向分量電場大約在1.2 ms之后到達(dá)O 波反射區(qū)域,并形成駐波.圖11 是在O 波轉(zhuǎn)換區(qū)域(277—280 km)以及更高區(qū)域(310—320 km)隨時間變化的靜電擾動及電子密度(ns=ne-ni0)圖像.

圖10 背景電子密度及t=1.27 ms 時刻電場 E x ,E y 和Ez分量幅值分布Fig.10.Background electron density and distribution of electric field E x ,E y and E z component amplitude at t=1.27 ms.

從圖11(a)可以看出,大約在1.2 ms,有參量不穩(wěn)定性引起的靜電波擾動導(dǎo)致O 波轉(zhuǎn)換為Z 波,并向電離層的更高區(qū)域傳播.圖9 中電場在O 波反射點附近的快速振蕩也是由于參量不穩(wěn)定性過程中的Langmuir 波造成的,即高頻Langmuir波的連續(xù)振蕩引起電子運(yùn)動.如圖11(b)所示,在相互作用區(qū)域周圍越來越多的電子參與波的相互作用.從圖11(c)可以明顯看到O 波以上,位于314—316 km 之間依然有波的存在,此處為Z 波的轉(zhuǎn)換區(qū)域[30].從圖11(c)可以看出在Z 波反射區(qū)域,其電場分量的高度非常窄,這里是Z 波經(jīng)過參量不穩(wěn)定性激發(fā)出的短波段Langmuir 波,對應(yīng)的電子密度如圖11(d)所示.

圖11 O 波轉(zhuǎn)換 區(qū)域的 靜電擾動及 ns 擾 動 (a) O 波轉(zhuǎn)換區(qū)域 | Ez| 變 化;(b) O 波轉(zhuǎn)換區(qū)域 | ns| 變 化;(c) Z 波轉(zhuǎn)換區(qū)域 | Ez| 變 化;(d) Z 波轉(zhuǎn)換區(qū)域 | ns| 變 化Fig.11.Electrostatic disturbance and ns disturbance in O-wave conversion region:(a) Variation of O-wave conversion region with respect to | Ez| ; (b) variation of O-wave conversion region with respect to | ns| ;(c) variation of Z-wave conversion region with respect to | Ez| ; (d) variation of Z-wave conversion region with respect to | ns| .

4 結(jié)論

通過二維麥克斯韋方程等價轉(zhuǎn)換成一維麥克斯韋方程組,與雙流體力學(xué)方程進(jìn)行聯(lián)立,建立了斜入射電離層等離子體的數(shù)值模型,模擬了大功率電磁波在較小傾斜角度加熱電離層情況下O 波的參量不穩(wěn)定性.結(jié)果表明:在毫秒量級的時間尺度內(nèi),大功率高頻電磁波在小角度的斜入射加熱下,O 波反射區(qū)域有參量不穩(wěn)定性引起的靜電波擾動,O 波發(fā)生明顯的衰減并激發(fā)了Langmuir 波的形成以及Langmuir 擾動所引起的O 波轉(zhuǎn)換為Z 波現(xiàn)象.其中,Z 波在1.2 ms 左右到達(dá)O 波反射區(qū)域并形成駐波.