• <tr id="yyy80"></tr>
  • <sup id="yyy80"></sup>
  • <tfoot id="yyy80"><noscript id="yyy80"></noscript></tfoot>
  • 99热精品在线国产_美女午夜性视频免费_国产精品国产高清国产av_av欧美777_自拍偷自拍亚洲精品老妇_亚洲熟女精品中文字幕_www日本黄色视频网_国产精品野战在线观看 ?

    IMF 北向時太陽風粒子向磁層輸運的試驗粒子模擬研究*

    2023-01-14 12:49:10趙明現(xiàn)
    空間科學學報 2022年6期
    關(guān)鍵詞:磁層太陽風步長

    趙明現(xiàn)

    (中國氣象局空間天氣重點開放實驗室/國家衛(wèi)星氣象中心(國家空間天氣監(jiān)測預警中心) 北京 100081)

    (許健民氣象衛(wèi)星創(chuàng)新中心 北京 100081)

    0 引言

    太陽風和行星際磁場(IMF)與地球內(nèi)稟磁場相互作用,形成了向陽面被壓縮、背陽面被拉伸的地球磁層。磁層中的等離子體主要有內(nèi)外兩個來源,內(nèi)源為電離層,外源為太陽風[1]。一般來說,太陽風等離子體對磁層中的遠地區(qū)域起支配作用,而來自電離層的粒子對磁層中的近地區(qū)域起主導作用[2,3]。

    太陽風等離子體向磁層輸運是太陽風-磁層相互作用的一個重要過程,對這個過程物理機制的探討一直是一個重要的研究方向。Dungey[4]提出的磁場重聯(lián)機制是目前普遍接受的一種理論,在向陽面磁層頂日下點發(fā)生磁場重聯(lián),從而剝離向陽面磁層中的磁力線,并將之拖向背陽面磁尾。日下點的磁場重聯(lián)受IMF 控制,IMF 南向時的重聯(lián)率要顯著的高于IMF北向時[5],此時太陽風粒子主要通過晨側(cè)磁層頂直接進入磁層中的等離子體片和環(huán)電流區(qū)域[2]。在IMF北向期間,磁重聯(lián)發(fā)生在高緯的南北極尖區(qū)背陽面一側(cè)[6,7],這個過程會在向陽面磁層中形成新的閉合磁力線,閉合磁力線攜帶太陽風粒子沉入磁層內(nèi)部,并沿方位角方向往晨昏兩側(cè)形成低緯邊界層(Low Latitude Boundary Layer,LLBL)。此 外,IMFBx和By在太陽風-磁層能量耦合過程中也起著重要的作用,同樣條件下IMF 徑向時進入磁層的能量比IMF北向時要大[8]。

    Axford 和Hines[9]提出,太陽風與磁層頂側(cè)面等離子體的黏性相互作用會在磁層中形成對流,從而使太陽風動量傳入磁層內(nèi)。這種機制不受IMF 方向的影響,不過對太陽風向磁層的能量傳輸貢獻并不大,大概在10%~20%[10]。K-H 不穩(wěn)定(Kelvin-Helmholtz Instability,KHI)也是一種可能的傳輸機制。觀測表明,在低緯的晨昏兩側(cè),太陽風等離子體相對于磁層頂內(nèi)的等離子體具有很高的速度,K-H 波在邊界層會形成局地的小尺度渦旋[11,12]。雖然具體的物理過程尚無定論,但是從衛(wèi)星觀測數(shù)據(jù)分析,在IMF 北向時,磁尾晨昏兩側(cè)的KHI 會導致磁鞘中的太陽風向磁層內(nèi)的傳輸[13]。另外,磁鞘中的動力阿爾芬波也能讓磁鞘等離子體進入磁層[14,15],如果磁鞘中存在動量足夠大的不均勻結(jié)構(gòu)體,也可以直接穿過磁層頂進入磁層[16,17]。不同的傳輸機制之間并不互斥,在實際的物理過程中,可能同時起作用[18]。

    試驗粒子方法已經(jīng)在空間物理和天文學中得到廣泛的應用,常用于日地空間中粒子的起源、傳播和損耗等過程的研究。試驗粒子法需要給定背景場,包括電場和磁場等,通過跟蹤單個或多個粒子的動力學過程,來研究復雜系統(tǒng)中的物理過程。由于不考慮試驗粒子對背景場的反饋以及各試驗粒子之間的相互作用,所以這個方法是不自洽的。但是在研究大尺度的復雜問題時,自洽模型會帶來巨大的運算量,而作為系統(tǒng)一階近似的試驗粒子方法足夠研究很多粒子的動力學問題,是大尺度磁流體模擬的有益補充[19]。也可以在計算試驗粒子軌道時,使用時間回推方法,研究粒子起源[20]。

    試驗粒子法又可以分為粒子軌道法、蒙特卡羅法等[19],粒子軌道法通過單個粒子或一群粒子的運動軌道,來研究系統(tǒng)中的粒子輸運過程,本文研究使用的是試驗粒子軌道法。Speiser[21]使用試驗粒子軌道法基于Dungey 磁層模型研究了太陽風粒子向磁層電流片中的輸運和極光粒子的加速,但其使用的磁層背景場設定太過簡化,與真實物理過程相去甚遠。Delcourt[22]通過單粒子的軌道變化研究了磁層中氧離子在地磁亞暴期間的加速過程。Li等[23,24]研究亞暴期間能量粒子注入磁層的過程,Tang等[25]對磁暴期間能量粒子進入磁層做了數(shù)值模擬研究,這些研究著重于能量粒子,而不是能量較低的太陽風粒子。Moore等[2]使用四階龍格庫塔方法求解全洛倫茲方程,以LFM 磁層模型的磁流體(MHD)模擬結(jié)果為背景場,研究了太陽風和極區(qū)電離層向磁層中的粒子輸運特征,不過其所使用的MHD 背景場是選擇了IMF 南向和北向兩個時刻,是固定的背景場,此外該MHD模擬并未考慮地磁軸傾角的影響。Ding[3,26]使用試驗粒子法模擬計算了太陽風粒子向磁層等離子體片的輸運,背景場使用BATS-R-US 的MHD 模擬結(jié)果,背景場隨時間變化且考慮了地磁傾角,不過MHD 模擬使用的是人工設定的理想太陽風條件,并非實際觀測得到的太陽風數(shù)據(jù),另外粒子軌道的計算采用了引導中心近似,而引導中心近似方法得到的粒子運動軌跡與全洛倫茲方程相比會存在一定的誤差[3]。Guo等[27]基于自主開發(fā)的磁層MHD 模型和試驗粒子法,研究了磁尾等離子體片的粒子注入機制,其中磁層背景是給定行星際條件下計算得到的穩(wěn)態(tài)結(jié)果,通過對19600 個粒子的追蹤,討論了一些典型的粒子運動軌跡。Cao等[28]研究了非磁暴期間IMF 南向時太陽風粒子向磁層的注入,不過對IMF 北向的情景,只討論了單個粒子的運動軌道特征,沒有做大規(guī)模粒子運動的模擬研究,而且使用的BATS-R-US 磁層MHD 模擬結(jié)果是理想太陽風狀況下的穩(wěn)態(tài)結(jié)果,不隨時間變化。

    本文以ACE 衛(wèi)星實時觀測數(shù)據(jù)驅(qū)動的BATSR-US 全球MHD 模擬為背景場,使用試驗粒子法對北向IMF 時太陽風粒子向磁層中的輸運過程進行模擬,研究北向IMF 時太陽風粒子注入磁層的過程以及粒子在磁層中的分布特征。

    1 數(shù)值模擬方法

    雖然使用固定的(或人工設定的隨時間變化的)行星際參數(shù)驅(qū)動磁層MHD 模型可以得到IMF 北向時的磁層狀態(tài)[3,26,28],但是使用衛(wèi)星實時觀測的太陽風和IMF 參數(shù)作為驅(qū)動,并將磁傾角等參數(shù)設定為實際數(shù)值,得到的MHD 模擬結(jié)果會更加接近實際情況。在2003 年10 月22-24 日,有一次持續(xù)時間較長的北向IMF 事件,ACE 衛(wèi)星觀測的太陽風和IMF數(shù)據(jù)如圖1 所示。選用2003 年10 月22 日00:00 UT作為太陽風-磁層MHD 模擬的起始時間,每4 min保存一次MHD 模擬結(jié)果。通過對空間進行三維插值,可以求得模擬區(qū)域內(nèi)任意位置的MHD 背景場參數(shù),再通過線性時間插值,就可以得到一個隨時間變化的太陽風-磁層相互作用的背景場。太陽風-磁層MHD 數(shù)值模擬通過CCMC(Community Coordinated Modeling Center)*https://ccmc.gsfc.nasa.gov/實現(xiàn),磁層模型選用BATS-R-US。

    圖1 2003 年10 月22-25 日ACE 衛(wèi)星觀測到的太陽風和IMF 數(shù)據(jù)Fig.1 ACE satellite observation data on 22-25 October 2003

    在研究磁層粒子分布時常以離子為對象,這是因為電子具有更復雜多樣的起源機制,并且很難區(qū)分電子的來源[3]。太陽風中的離子主要是質(zhì)子,所以本文以質(zhì)子為追蹤粒子,研究太陽風粒子向磁層輸運的過程。帶電粒子在行星際和磁層中的運動可以用以下洛倫茲方程來描述:

    其中,q和m為 粒子電荷和質(zhì)量,v為粒子速度,p為粒子動量,γ為洛倫茲因子,c為真空光速,E和B分別為粒子所在位置的電場強度和磁感應強度,其之間滿足E=-vMHD×B,這里vMHD為磁流體速度。

    在滿足絕熱近似的條件時,可以用引導中心近似的方法描述粒子的運動軌跡,這樣可以將求解6 個未知量的洛倫茲方程簡化為求解4 個未知量的引導中心近似方程,同時引導中心近似方法也使得求解時可以使用更大的時間步長,從而大大地降低運算量,但是與全洛倫茲方程相比,引導中心近似方法得到的粒子運動軌跡會存在一定的誤差[3],而且由于使用了絕熱條件,所以并不能反映粒子的加熱過程。本文選用4 階龍格庫塔方法求解全洛倫茲方程以得到粒子運動軌跡,并使用變時間步長的方法,根據(jù)粒子的回旋周期,自動調(diào)整時間步長。時間步長的選擇也會對粒子的運動軌跡造成影響。Ding[3]使用引導中心近似求解中,認為應該將時間步長設置為約1/10 個回旋周期,對于全洛倫茲方法,則需要更小的步長。Peroomian 和EI-Alaoui[29]認為求解高能粒子的全洛倫茲方程時應該將步長設置為離子回旋周期的0.002 倍。Birdsall 和 Langdon[30]認為時間步長應該設定為εP/F,其中P為粒子動量,F(xiàn)為洛倫茲力,ε為一個預設很小的值[31]。Pulkkinen等[32]認為時間步長可以設定為粒子回旋周期的約1/180。

    為了選擇合適的時間步長,對單個粒子設定不同時間步長,比較其運動軌跡的差異。選擇一個穩(wěn)態(tài)的BATS-R-US 模擬結(jié)果作為背景場,粒子初始位置位于磁層內(nèi)部,坐標為(8.80,0.04,1.99)Re(Re為地球半徑)。粒子初始速度接近所在位置的磁流體速度(-228.4,70.5,206.8)km·s-1。分別計算時間步長 Δt為(1/20)T,(1/50)T,(1/100)T,(1/200)T(T為質(zhì)子所在位置的回旋周期)時質(zhì)子的運動軌跡,質(zhì)子位置坐標的存儲時間間隔為1 s,質(zhì)子運行軌跡如圖2 所示。Δt=(1/20)T時,質(zhì)子的軌跡與其他情況相比差異明顯,其位置超出模擬區(qū)域的時間(逃逸時間)為30917 s。Δt=(1/50)T時,運動軌跡與(1/100)T相比在磁尾區(qū)域稍有差異,其逃逸時間為43898 s。Δt=(1/100)T和 Δt=(1/200)T兩種情況的質(zhì)子運動軌跡非常接近,其逃逸時間分別為44548 s 和44352 s。綜合考慮計算精度和計算量,在本文的大規(guī)模粒子模擬中,時間步長設定為(1/100)T。

    圖2 使用不同時間步長計算得到的粒子運動軌跡Fig.2 Trajectories of a single particle using different settings of time step

    2 數(shù)值模擬結(jié)果

    CCMC 得到的太陽風-磁層MHD 背景場使用的是非均勻網(wǎng)格,為了提高插值運算速度,先將所有時刻的背景場插值轉(zhuǎn)換為間隔為1Re的均勻網(wǎng)格。背景場采用GSM 坐標系,在x方向的范圍為(-80~20)Re,在y和z方向的范圍均為(-40~40)Re。試驗粒子數(shù)量為107個,入射位置隨機分布于x=19Re的平面上,粒子入射初始速度由觀測數(shù)據(jù)插值計算得到,溫度設定為 10 eV,熱速度滿足麥克斯韋分布。初始時每秒入射104個,后續(xù)逃逸出模擬區(qū)域的粒子放在入射位置重新注入模擬區(qū)域。由于每次入射的粒子數(shù)量會存在差異,所以對每次的入射粒子均會賦予一個權(quán)重系數(shù),從而使得折算后的粒子密度接近于實際的太陽風粒子密度。權(quán)重系數(shù)的關(guān)系式為

    其中,W(t)為t時刻入射粒子的權(quán)重系數(shù),n(t)為t時刻太陽風粒子數(shù)密度,A為入射位置的截面面積,N(t)為t時刻的入射粒子個數(shù)。

    每隔100 s,輸出一次試驗粒子分布的計算結(jié)果。將模擬區(qū)域內(nèi)的每個粒子,按其權(quán)重系數(shù)和粒子所在網(wǎng)格中的位置與網(wǎng)格頂點的距離,線性分配給網(wǎng)格的8 個頂點。每個格點周圍的粒子對其貢獻之和,即為該點處的粒子數(shù)密度。

    2.1 與觀測數(shù)據(jù)的對比

    為驗證模擬結(jié)果中粒子分布情況,將模擬結(jié)果與Geotail 衛(wèi)星觀測結(jié)果進行對比。由于Geotail 衛(wèi)星在10 月22-23 日處于本文設定的模擬區(qū)域之外(GSM 坐標x> 20Re),所以只能比較24-25 日期間衛(wèi)星進入模擬區(qū)域的結(jié)果。IMF 長時間北向是在23 日,與Geotail 觀測數(shù)據(jù)的比較無法獲取IMF 北向期間的信息,但模擬結(jié)果與觀測結(jié)果的對比依然有助于確定模擬結(jié)果的可信度。

    如圖3 所示,圖3(a)為Geotail 衛(wèi)星觀測的低能粒子(Low Energy Particle)數(shù)密度,圖3(b)為粒子模擬結(jié)果中Getotail 衛(wèi)星所處位置的粒子數(shù)密度,圖3(c)為Geotail 衛(wèi)星在GSM 坐標系中的x坐標位置。

    圖3 Geotail 衛(wèi)星觀測粒子密度與模擬結(jié)果對比Fig.3 Comparison between the particle density observed by Geotail and that from our simulation

    太陽風粒子從模擬區(qū)域的邊界(x=20Re處)進入模擬區(qū)域,首先進入的是弓激波上游的太陽風區(qū)域,模擬結(jié)果中粒子密度比較穩(wěn)定。太陽風進入磁鞘后,粒子密度有明顯的增強和波動。最后在磁層內(nèi)部,粒子密度相比磁鞘區(qū)域顯著降低。觀測數(shù)據(jù)的粒子密度比模擬結(jié)果要大數(shù)倍,這主要是因為本文的模擬粒子數(shù)量有限所致,另外,模擬中只計算了質(zhì)子,未考慮其他粒子的影響。Geotail 的觀測數(shù)據(jù)存在一些缺失,但是其變化趨勢與模擬結(jié)果基本吻合,表明本文的模擬結(jié)果能夠反映太陽風-磁層相互作用過程中粒子的分布情況。

    2.2 太陽風粒子在磁層中的三維分布

    在IMF 北向期間,太陽風粒子可以通過多種機制和途徑進入內(nèi)磁層和磁尾區(qū)域[18],由于IMF 北向時的太陽風進入近地磁尾并未經(jīng)過顯著的加熱加速過程,所以與IMF 南向時相比,在磁尾中的等離子體片溫度偏低,密度偏高,從而形成冷而密的等離子體片(Cold Dense Plasma Sheet,CDPS)[33]。在模擬時間t=66600 s,對應物理時間為2003 年10 月22 日18:30 UT 左右,此時是長時間IMF 北向的起始階段,進入磁層的太陽風粒子密度分布如圖4 所示(圖中磁層頂位置是通過流線方法識別得到的[34],弓激波位置是利用磁流體速度突變?yōu)榕卸藴首R別得到的)。此時進入磁層中的粒子數(shù)量較少,進入磁層的粒子主要分布在南北極尖區(qū)和環(huán)電流區(qū)域。在磁層頂附近,越靠近磁尾,越過磁層頂進入磁層的太陽風粒子越多,不過這些粒子主要分布在磁層頂內(nèi)側(cè)附近,很難深入磁層內(nèi)部。

    圖4 太陽風粒子在磁層中的三維密度分布(黑色實線為磁層頂位置,白色實線為弓激波位置)Fig.4 Three dimensional density distribution of solar wind particles in magnetosphere (Black solid lines mark the location of magnetopause,and the white solid lines mark the location of bow shock)

    再選取模擬時間為170000 s 的模擬結(jié)果,對應物理時間為2003 年10 月23 日23:13 UT 左右,此時IMF 長時間的北向過程已經(jīng)快要結(jié)束,太陽風粒子在磁層中的三維密度分布如圖5 所示。與IMF 北向起始階段不同,此時已經(jīng)有大量的太陽風粒子分布在磁層內(nèi)部,南北極尖區(qū)和環(huán)電流區(qū)域的粒子數(shù)顯著增多,在磁尾區(qū)域,出現(xiàn)了粒子聚集區(qū)。為了更清晰分析磁層內(nèi)部粒子密度分布,取xy,xz,yz三個截面進行分析。

    圖6 為圖5 所示的赤道截面,由于考慮了磁傾角,所以地球赤道面與GSM 坐標的xy平面并不重合,為了便于觀察近地磁層中的粒子分布,選取的截面接近赤道面。在晨昏兩側(cè)的磁層頂內(nèi),均有太陽風粒子的注入。等離子體片晨昏兩翼中低能粒子密度的增加由磁鞘粒子注入引起[35,36],粒子可以穿越側(cè)翼的磁層頂進入,也可以通過極尖區(qū)進入磁層。Sorathia等[37]認為通過側(cè)翼磁層頂?shù)腒H 渦旋進入磁層的粒子在晨昏兩側(cè)的分布是對稱的,而通過極尖區(qū)進入磁層的粒子表現(xiàn)出晨昏不對稱性,進入晨側(cè)等離子體片的粒子要多于昏側(cè)。圖6 中磁尾等離子體片靠近磁層頂位置的晨昏兩翼存在晨昏不對稱,晨側(cè)粒子密度要大于昏側(cè),這與Wing等[36]的統(tǒng)計結(jié)果相符。近地空間的環(huán)電流區(qū)域捕獲了大量的太陽風粒子,在地球赤道面形成完整的環(huán)狀分布。環(huán)電流區(qū)域的粒子分布也表現(xiàn)出晨昏不對稱性,晨側(cè)粒子分布更寬,這可能是由于晨向電場與向陽面磁層內(nèi)離子的西向漂移共同作用產(chǎn)生的[37],Yao等[38]統(tǒng)計結(jié)果表明晨側(cè)磁層頂附近的動理學阿爾芬波能量要高于昏側(cè),這也可以導致晨側(cè)穿越磁層頂進入磁層的粒子比昏側(cè)多。

    圖5 太陽風粒子在磁層中的三維密度分布Fig.5 Three dimensional density distribution of solar wind particles in magnetosphere

    圖6 太陽風粒子在xy 平面的密度分布Fig.6 Density distribution of solar wind particles in xy plane

    子午面(xz)中粒子密度分布如圖7 所示,圖7 中可以清晰分辨出弓激波和磁層頂?shù)奈恢?,尤其是北?cè)磁層頂邊界清晰,說明從北側(cè)磁層頂進入磁層的粒子極少。南側(cè)磁層頂內(nèi)側(cè)有少量粒子分布,特別是在遠磁尾區(qū)域,說明有少量粒子可以通過南側(cè)磁層頂進入磁層,這種南北兩側(cè)磁層頂?shù)牟町惪赡苁鞘艿卮排紭O傾角的影響。觀察南北極尖區(qū)中的粒子密度,南半球的極尖區(qū)中粒子密度更高,說明太陽風粒子更容易進入南半球極尖區(qū),這也與地磁偶極傾角有關(guān)[3]。

    圖7 太陽風粒子在xz 平面的密度分布Fig.7 Density distribution of solar wind particles in xz plane

    在磁尾-15Re<z< 5Re,x> -50Re的區(qū)域,存在粒子聚集區(qū),也就是長時間IMF 北向時,在磁尾形成的冷而密等離子體片。Bai等[33]統(tǒng)計CDPS 事件主要分布在x< -10Re的區(qū)域,Li等[39]認為大部分的CDPS 觀測位于x> -30Re的近地區(qū)域,Terasawa等[40]認為北向IMF 時磁尾粒子密度增加發(fā)生在-50Re<x< -15Re,的區(qū)域,這與圖7 中的分布范圍基本符合。同時也注意到,在x≈ -75Re附近,存在另一個粒子密度較大的區(qū)域。

    CDPS 事件的基本特征為粒子密度較大(>1 cm-3),10 月22 日IMF 北向持續(xù)3 h 后,Cluster 衛(wèi)星觀測到粒子密度達到1.5 cm-3,在這個水平維持了約30 h,DMSP 觀測到此期間的粒子密度平均為 2~3 cm-3[41]。圖7 中磁尾等離子體片區(qū)域密度較大處的lgN略大于6,也即粒子密度略大于 1 cm-3。雖然密度最大的區(qū)域滿足CDPS 對密度的要求,但是與觀測數(shù)據(jù)相比還是略小,這是因為模擬中的粒子數(shù)量較少所致,雖然采取了給粒子賦予權(quán)重的方法,但與實際觀測相比,還存在粒子數(shù)量偏少的問題。

    圖8 所示為磁尾x=-20Re處yz平面內(nèi)的粒子密度分布。磁尾南北尾瓣粒子分布極少,位于赤道附近的低緯邊界層增厚,這與Sckopke等[42]的結(jié)果相一致。LLBL 區(qū)域的粒子密度也比較高,lgN可達到6~6.5 cm-3左右,也即粒子密度在 1~3 cm-3左右。Ding[3]計算得到IMF 北向期間LLBL 密約為3.2 cm-3,而觀測結(jié)果一般在0.1~10 cm-3,均與本文的計算結(jié)果吻合。除了赤道附近的LLBL,還存在一個從左上(北-晨側(cè))到右下(南-昏側(cè))的粒子分布結(jié)構(gòu)。Guo等[43]分析認為,在IMF 北向時,磁尾橫斷面會出現(xiàn)顯著的由北-晨側(cè)和南-昏側(cè)指向磁層內(nèi)部的粒子流,可以用來解釋IMF 北向時磁鞘粒子對等離子體片更高的注入效率。

    圖8 太陽風粒子在 x=-20 Re 處的yz 平面的密度分布Fig.8 Density distribution of solar wind particles in yz plane at x=-20 Re

    2.3 磁尾粒子數(shù)量的時間演化

    太陽風粒子進入磁層的過程并不是一個穩(wěn)態(tài)的過程,粒子進入磁層的速率也不是一成不變的。在長時間的IMF 北向期間,磁尾電流片中的粒子總數(shù)也應該隨時間有所變化,為了研究磁尾粒子數(shù)量的變化情況,簡單劃定一個長方體形狀的區(qū)域(-40Re<x<-10Re,-15Re<y< 15Re,-10Re<z< 10Re),分析其中粒子平均密度和最大密度值隨時間的演化特征。如圖9 所示,圖9(a)為磁尾劃定區(qū)域內(nèi)的粒子平均密度,圖9(b)為區(qū)域內(nèi)粒子密度最大值。起始時間為IMF長時間北向的開始(t=66600 s,22 日18:30 UT),結(jié)束時間為IMF長時間北向結(jié)束(t=181800 s,24 日02:30 UT)的時刻,豎虛線為前文分析三維密度分布所使用的時刻(t=170000 s)。

    圖9 磁尾選定區(qū)域(-40 Re < x < -10 Re,-15 Re < y < 15 Re,-10 Re < z < 10 Re )內(nèi)的粒子平均密度和密度最大值在IMF 北向期間隨時間的變化Fig.9 Variation of total particle average density and particle max density maximum during northward IMF in the selected magnetotail region (-40 Re < x < -10 Re,-15 Re < y < 15 Re,-10 Re < z < 10 Re)

    粒子平均密度變化曲線中,可以看到在大部分時間都存在一個幅度較小的準周期性變化,周期大約為20 min,在平均密度曲線沒有劇烈變化時比較明顯。最大密度值曲線變化趨勢與平均密度曲線相仿,不過由于密度最大值存在較大的偶然性,所以曲線抖動比較明顯。Sorathia等[37]分析了約1 h 的粒子注入速率,發(fā)現(xiàn)晨昏側(cè)翼的粒子注入存在約幾分鐘的短時變化,極尖區(qū)注入的粒子存在約8 min 的準周期變化。粒子平均密度存在的約20 min 的短時變化特征,可能跟磁層頂晨昏兩翼以及南北極尖區(qū)等的注入機制存在短周期變化有關(guān),是多種機制共同作用的結(jié)果。

    從兩個曲線都可以看出,在IMF 長時間北向期間,磁尾區(qū)域粒子數(shù)量存在持續(xù)增加的趨勢,在這個趨勢上,平均密度曲線還存在一個幅度較大的準周期變化,周期約為 5~6 h。這表明,在IMF 長時間北向期間,磁尾等離子體片中的粒子并不是處于穩(wěn)定狀態(tài)的,而是可能存在一種從聚集到釋放,然后再聚集的循環(huán)過程。

    3 結(jié)論

    以ACE 衛(wèi)星實時觀測數(shù)據(jù)驅(qū)動的BATS-R-US全球MHD 模擬作為背景場,選取2003 年10 月22-24 日的長時間IMF 北向事件,使用試驗粒子法對太陽風粒子向磁層輸運的過程進行模擬研究,分析討論北向IMF 時太陽風粒子注入磁層的過程以及粒子在磁層中的分布和時間演化特征。

    使用全洛倫茲方程對粒子軌道進行求解,并選用變時間步長的方法,根據(jù)粒子的回旋周期,自動調(diào)整時間步長。通過對比分析幾種不同的時間步長取值對粒子運動軌跡的影響,綜合考慮求解精度和運算量,最終選用時間步長為回旋周期的1/100。將一千萬個粒子從磁層外太陽風上游隨機選擇位置入射進入模擬區(qū)域,跑出模擬區(qū)域的粒子重新回到入射位置,使得模擬區(qū)域內(nèi)一直保持一千萬粒子的規(guī)模。模型中粒子密度與Geotail 衛(wèi)星觀測數(shù)據(jù)具有一致的變化趨勢,表明模擬結(jié)果可以反映太陽風-磁層相互作用中的粒子分布情況。

    在IMF 長時間北向的初期,進入磁層的粒子數(shù)量較少,主要分布在環(huán)電流區(qū)域和南北極尖區(qū)。IMF 長時間持續(xù)北向后,不僅南北極尖區(qū)和環(huán)電流區(qū)域的粒子數(shù)量大量增加,在磁尾也聚集了大量的粒子。磁尾等離子體片粒子密度分布存在晨昏不對稱,晨側(cè)粒子密度高于昏側(cè);環(huán)電流區(qū)域的粒子分布也存在晨昏不對稱,在晨側(cè)粒子分布的范圍更廣。磁尾南北側(cè)磁層頂邊界清晰,南北尾瓣粒子密度極小,而低緯邊界層區(qū)域粒子密度大,表明磁鞘中的太陽風粒子基本不會通過南北磁層頂進入磁層,低緯邊界層是磁鞘粒子進入磁層的一個主要通道。

    長時間的IMF 北向使得太陽風粒子進入磁層并積累在磁尾,形成了冷而密的等離子體片。磁尾等離子體片的粒子聚集區(qū)在x> -50Re,-15Re<z< 5Re的區(qū)域,粒子密度略大于1 cm-3,達到了CDPS 事件的標準。在磁尾存在從北-晨側(cè)向南-昏側(cè)的粒子密度分布結(jié)構(gòu),這與IMF 北向時,磁鞘粒子向磁層的注入增強有關(guān)。

    本文分析了磁尾一個指定區(qū)域內(nèi)的粒子數(shù)隨時間的演化特征。在IMF 長時間北向期間,磁尾粒子數(shù)量呈現(xiàn)出明顯的增長趨勢,并存在約20 min 的小幅度準周期變化和約 5~6 h 的較大幅度的準周期變化。小幅度的短時變化可能是由粒子注入機制存在短周期變化引起,較大幅度的長周期變化則表明磁尾等離子體片中的粒子并不是穩(wěn)態(tài)的,而是應該存在一種從聚集到釋放,然后再聚集的循環(huán)過程。

    致謝ACE數(shù)據(jù)由ACE SWEPAM儀器團隊和ACE 科學中心提供。模擬結(jié)果由位于戈達德太空飛行中心的CCMC(http://ccmc.gsfc.nasa.gov)提供。BATS-RUS 模型由密歇根大學空間環(huán)境建模中心(CSEM)Tamas Gombosi 博士等開發(fā)。Geotail 等離子體數(shù)據(jù)由Saito Y 通過日本宇宙科學研究所的數(shù)據(jù)檔案和傳輸系統(tǒng)(DARTS)提供。

    猜你喜歡
    磁層太陽風步長
    基于Armijo搜索步長的BFGS與DFP擬牛頓法的比較研究
    基于對抗神經(jīng)網(wǎng)絡有限角度的磁層邊界CT重構(gòu)技術(shù)
    磁鞘快速流引起的磁層頂凹陷事件
    等離子體鞘層效應對磁層探測電場儀設計的影響
    多種觀測數(shù)據(jù)驅(qū)動的三維行星際太陽風MHD模擬
    地球磁層亞暴統(tǒng)計分析
    基于ACE飛船觀測的銀河宇宙線與太陽風變化的統(tǒng)計研究
    在太空聽聲音
    戈爾的夢想飛船啟航:美國發(fā)射首顆位于L1點的太陽風監(jiān)測衛(wèi)星
    太空探索(2015年5期)2015-07-12 12:52:26
    基于逐維改進的自適應步長布谷鳥搜索算法
    成人午夜精彩视频在线观看| 久久国产亚洲av麻豆专区| 高潮久久久久久久久久久不卡| 人人妻,人人澡人人爽秒播 | 这个男人来自地球电影免费观看| 人人妻人人澡人人爽人人夜夜| 国产在线免费精品| 国产三级黄色录像| 男女免费视频国产| 欧美精品一区二区免费开放| 国产成人91sexporn| 一本色道久久久久久精品综合| 久久久国产一区二区| 亚洲av电影在线观看一区二区三区| 日日爽夜夜爽网站| 婷婷色综合www| 丰满少妇做爰视频| 黄网站色视频无遮挡免费观看| 精品福利观看| 成人国产一区最新在线观看 | 亚洲中文av在线| 成人亚洲精品一区在线观看| 欧美亚洲 丝袜 人妻 在线| av福利片在线| 2018国产大陆天天弄谢| 高潮久久久久久久久久久不卡| 99九九在线精品视频| 国产精品一区二区精品视频观看| 亚洲精品一区蜜桃| 精品人妻在线不人妻| 啦啦啦 在线观看视频| 国产精品av久久久久免费| 电影成人av| 欧美日韩福利视频一区二区| 日日夜夜操网爽| 天堂俺去俺来也www色官网| 日韩精品免费视频一区二区三区| 欧美在线一区亚洲| av又黄又爽大尺度在线免费看| 亚洲中文日韩欧美视频| 99久久人妻综合| 99国产精品免费福利视频| 久久精品成人免费网站| 成人18禁高潮啪啪吃奶动态图| 亚洲欧美中文字幕日韩二区| avwww免费| 视频区欧美日本亚洲| 久久免费观看电影| 国产成人免费观看mmmm| 免费一级毛片在线播放高清视频 | 免费少妇av软件| 91麻豆精品激情在线观看国产 | 精品国产一区二区久久| 国产高清国产精品国产三级| 国产片特级美女逼逼视频| 国产高清视频在线播放一区 | 亚洲国产成人一精品久久久| 午夜福利一区二区在线看| 国产三级黄色录像| 人妻一区二区av| 国产日韩欧美视频二区| 免费在线观看日本一区| 精品久久久久久电影网| 叶爱在线成人免费视频播放| 欧美变态另类bdsm刘玥| 18禁观看日本| 丝袜脚勾引网站| 国产成人精品久久二区二区免费| 国产视频首页在线观看| 七月丁香在线播放| 国产主播在线观看一区二区 | 91国产中文字幕| 久久久久视频综合| av国产久精品久网站免费入址| 亚洲五月色婷婷综合| 好男人电影高清在线观看| 国产xxxxx性猛交| 国产亚洲欧美在线一区二区| 香蕉国产在线看| 日本vs欧美在线观看视频| av天堂久久9| 美国免费a级毛片| 欧美性长视频在线观看| 久久久久久久大尺度免费视频| 日本五十路高清| 亚洲成人国产一区在线观看 | 国产极品粉嫩免费观看在线| 91麻豆av在线| 高清欧美精品videossex| 999精品在线视频| 国产真人三级小视频在线观看| 国产欧美日韩精品亚洲av| 亚洲欧美一区二区三区黑人| 菩萨蛮人人尽说江南好唐韦庄| 国产高清videossex| 国产1区2区3区精品| 黄色视频在线播放观看不卡| 人人妻,人人澡人人爽秒播 | 亚洲av日韩精品久久久久久密 | 美女午夜性视频免费| 91麻豆精品激情在线观看国产 | 99久久99久久久精品蜜桃| 99热全是精品| 无遮挡黄片免费观看| 99re6热这里在线精品视频| 51午夜福利影视在线观看| 欧美日韩视频精品一区| 各种免费的搞黄视频| 天堂8中文在线网| 男男h啪啪无遮挡| 午夜激情av网站| 99久久综合免费| 无限看片的www在线观看| 亚洲 欧美一区二区三区| 曰老女人黄片| 天天躁夜夜躁狠狠躁躁| 欧美黑人精品巨大| 欧美日韩亚洲综合一区二区三区_| 国产麻豆69| 伦理电影免费视频| 亚洲欧美激情在线| 美女高潮到喷水免费观看| 午夜老司机福利片| 欧美精品啪啪一区二区三区 | 欧美97在线视频| 欧美精品人与动牲交sv欧美| 我的亚洲天堂| 亚洲av在线观看美女高潮| 精品国产一区二区三区久久久樱花| 操美女的视频在线观看| a 毛片基地| 久久久欧美国产精品| 最新的欧美精品一区二区| 精品一区二区三区av网在线观看 | 中国国产av一级| 亚洲成人免费电影在线观看 | 精品久久久精品久久久| 国产精品久久久人人做人人爽| 久久99精品国语久久久| 色播在线永久视频| 亚洲第一av免费看| 只有这里有精品99| av国产久精品久网站免费入址| 国产一区有黄有色的免费视频| 高清黄色对白视频在线免费看| 午夜久久久在线观看| 亚洲欧美日韩另类电影网站| 大香蕉久久网| xxx大片免费视频| 九色亚洲精品在线播放| 国产精品九九99| 国产成人一区二区三区免费视频网站 | 亚洲精品国产色婷婷电影| 在线观看www视频免费| 啦啦啦 在线观看视频| 亚洲av男天堂| 九色亚洲精品在线播放| 日韩熟女老妇一区二区性免费视频| 一级毛片 在线播放| 男人爽女人下面视频在线观看| a级片在线免费高清观看视频| 成人黄色视频免费在线看| 免费一级毛片在线播放高清视频 | 亚洲精品成人av观看孕妇| 亚洲成av片中文字幕在线观看| 又大又黄又爽视频免费| 老司机在亚洲福利影院| 欧美激情极品国产一区二区三区| 2018国产大陆天天弄谢| 999久久久国产精品视频| 少妇裸体淫交视频免费看高清 | 久久久久国产一级毛片高清牌| 日韩 亚洲 欧美在线| www.999成人在线观看| 亚洲国产av新网站| 亚洲欧美日韩另类电影网站| 最近中文字幕2019免费版| 亚洲国产欧美在线一区| 老司机靠b影院| 欧美亚洲 丝袜 人妻 在线| 中文字幕高清在线视频| 午夜免费男女啪啪视频观看| 欧美黄色淫秽网站| 蜜桃国产av成人99| 人成视频在线观看免费观看| 国产成人免费无遮挡视频| 天天躁狠狠躁夜夜躁狠狠躁| 国产精品一国产av| www日本在线高清视频| 亚洲精品一二三| 亚洲一卡2卡3卡4卡5卡精品中文| 丝袜脚勾引网站| 人妻一区二区av| 一级毛片 在线播放| 国产亚洲av高清不卡| 性高湖久久久久久久久免费观看| 青春草亚洲视频在线观看| 成年av动漫网址| 国产成人精品久久久久久| 精品一区二区三区av网在线观看 | 欧美日韩国产mv在线观看视频| 人妻 亚洲 视频| 婷婷色综合大香蕉| 一二三四在线观看免费中文在| 女人久久www免费人成看片| 亚洲综合色网址| 亚洲精品在线美女| 欧美精品一区二区大全| 亚洲国产欧美网| 午夜福利乱码中文字幕| 国产成人免费观看mmmm| 久久久久久久大尺度免费视频| 亚洲国产精品一区二区三区在线| 中文字幕av电影在线播放| 日韩 欧美 亚洲 中文字幕| 亚洲色图 男人天堂 中文字幕| 国产精品久久久久成人av| 亚洲av美国av| 又粗又硬又长又爽又黄的视频| 欧美 日韩 精品 国产| 亚洲av日韩在线播放| 50天的宝宝边吃奶边哭怎么回事| 青春草视频在线免费观看| av福利片在线| 丰满饥渴人妻一区二区三| 国产av精品麻豆| 国产精品一区二区在线不卡| 欧美在线黄色| 精品久久久久久电影网| 欧美成人精品欧美一级黄| 亚洲一码二码三码区别大吗| 国产av国产精品国产| 久久狼人影院| 在线观看免费午夜福利视频| 久久久精品区二区三区| 日韩av不卡免费在线播放| 91精品国产国语对白视频| avwww免费| 又大又黄又爽视频免费| 侵犯人妻中文字幕一二三四区| 亚洲av日韩精品久久久久久密 | 精品人妻一区二区三区麻豆| 99精品久久久久人妻精品| a 毛片基地| 亚洲人成77777在线视频| 伦理电影免费视频| avwww免费| 日韩制服骚丝袜av| 亚洲精品中文字幕在线视频| 国产女主播在线喷水免费视频网站| 亚洲精品久久午夜乱码| 热99久久久久精品小说推荐| 国产精品亚洲av一区麻豆| 久久久精品区二区三区| 电影成人av| 脱女人内裤的视频| 亚洲国产看品久久| 亚洲精品乱久久久久久| 乱人伦中国视频| 超色免费av| 2021少妇久久久久久久久久久| 亚洲国产日韩一区二区| 国产男女内射视频| 人人妻,人人澡人人爽秒播 | 久久久久国产精品人妻一区二区| 国产成人精品久久二区二区91| 国产成人免费无遮挡视频| 国产不卡av网站在线观看| 久久久久久亚洲精品国产蜜桃av| 婷婷成人精品国产| 欧美日韩成人在线一区二区| 777久久人妻少妇嫩草av网站| 国产男女超爽视频在线观看| 在线观看免费高清a一片| 黄片小视频在线播放| 精品欧美一区二区三区在线| 18禁观看日本| 国产精品成人在线| 亚洲成人免费电影在线观看 | 欧美日韩黄片免| 亚洲黑人精品在线| 在线 av 中文字幕| 一级毛片黄色毛片免费观看视频| 午夜福利视频精品| 国产三级黄色录像| 亚洲av电影在线观看一区二区三区| 巨乳人妻的诱惑在线观看| 又大又黄又爽视频免费| 中文字幕制服av| 欧美日韩成人在线一区二区| 亚洲五月色婷婷综合| 国产成人精品久久久久久| 午夜日韩欧美国产| 黄色一级大片看看| 亚洲欧美日韩高清在线视频 | 日本五十路高清| 一区在线观看完整版| 亚洲欧美一区二区三区国产| 老鸭窝网址在线观看| 天天躁日日躁夜夜躁夜夜| 国产在线视频一区二区| 校园人妻丝袜中文字幕| 91成人精品电影| 天天添夜夜摸| 蜜桃国产av成人99| 热99久久久久精品小说推荐| 日韩伦理黄色片| 99久久99久久久精品蜜桃| 国产午夜精品一二区理论片| 首页视频小说图片口味搜索 | 男女下面插进去视频免费观看| 国产日韩一区二区三区精品不卡| 99久久人妻综合| 丝袜人妻中文字幕| 久久青草综合色| 中文字幕制服av| 国产男女超爽视频在线观看| 国产亚洲欧美精品永久| 黄色视频在线播放观看不卡| 亚洲国产精品999| 国产精品秋霞免费鲁丝片| 热99国产精品久久久久久7| 国产精品国产av在线观看| 男人操女人黄网站| 丝袜美足系列| 男女无遮挡免费网站观看| 香蕉国产在线看| 欧美乱码精品一区二区三区| 成年动漫av网址| av福利片在线| 欧美性长视频在线观看| 国产欧美日韩一区二区三 | 成人三级做爰电影| 亚洲欧美色中文字幕在线| 操美女的视频在线观看| 精品亚洲成国产av| 色综合欧美亚洲国产小说| 亚洲国产成人一精品久久久| 欧美日韩视频精品一区| 悠悠久久av| 国产成人影院久久av| 黄色视频在线播放观看不卡| 亚洲一码二码三码区别大吗| 9热在线视频观看99| 天堂中文最新版在线下载| av电影中文网址| 国产男女超爽视频在线观看| 一级毛片我不卡| 夫妻性生交免费视频一级片| 大片免费播放器 马上看| av有码第一页| 啦啦啦视频在线资源免费观看| 免费少妇av软件| 午夜福利免费观看在线| 日韩电影二区| 亚洲精品久久久久久婷婷小说| 精品一区二区三区av网在线观看 | 波野结衣二区三区在线| 国产一区亚洲一区在线观看| 水蜜桃什么品种好| 久久久国产精品麻豆| 国产精品三级大全| 精品人妻1区二区| 色精品久久人妻99蜜桃| 丰满少妇做爰视频| 国产精品二区激情视频| 啦啦啦在线观看免费高清www| 一级黄色大片毛片| 色94色欧美一区二区| 青春草亚洲视频在线观看| 日韩电影二区| 精品国产一区二区久久| 亚洲欧美激情在线| 男女下面插进去视频免费观看| 美国免费a级毛片| 国产一区二区三区综合在线观看| 男女边吃奶边做爰视频| 国产成人系列免费观看| 日本一区二区免费在线视频| 亚洲熟女毛片儿| 欧美精品人与动牲交sv欧美| 亚洲成人手机| 国产精品国产av在线观看| 99久久精品国产亚洲精品| 三上悠亚av全集在线观看| 汤姆久久久久久久影院中文字幕| 超碰成人久久| 19禁男女啪啪无遮挡网站| 国产精品欧美亚洲77777| 国产日韩欧美亚洲二区| 视频区图区小说| 亚洲欧美激情在线| 国产亚洲精品久久久久5区| 精品少妇内射三级| 精品少妇黑人巨大在线播放| 日本一区二区免费在线视频| 精品福利观看| 国产伦人伦偷精品视频| 在线观看国产h片| 中文字幕亚洲精品专区| 汤姆久久久久久久影院中文字幕| 我要看黄色一级片免费的| 久久久久精品人妻al黑| 亚洲精品乱久久久久久| 欧美成人精品欧美一级黄| 这个男人来自地球电影免费观看| 又黄又粗又硬又大视频| 亚洲国产日韩一区二区| 亚洲av电影在线观看一区二区三区| 久久久国产欧美日韩av| 欧美日韩福利视频一区二区| 亚洲中文字幕日韩| 日韩av不卡免费在线播放| 天天躁日日躁夜夜躁夜夜| 中文字幕人妻丝袜一区二区| xxx大片免费视频| 狠狠精品人妻久久久久久综合| 两性夫妻黄色片| 777久久人妻少妇嫩草av网站| 叶爱在线成人免费视频播放| 黄片播放在线免费| 一区二区日韩欧美中文字幕| 中文字幕色久视频| 欧美日韩黄片免| 亚洲精品久久久久久婷婷小说| 久久精品亚洲熟妇少妇任你| 99国产综合亚洲精品| 久久99热这里只频精品6学生| 80岁老熟妇乱子伦牲交| 一本色道久久久久久精品综合| 欧美 日韩 精品 国产| 精品福利永久在线观看| 午夜免费鲁丝| 男人添女人高潮全过程视频| 最新在线观看一区二区三区 | 日韩电影二区| 久久久精品国产亚洲av高清涩受| 啦啦啦啦在线视频资源| 国产有黄有色有爽视频| 国产成人av激情在线播放| 一个人免费看片子| 婷婷色综合大香蕉| 亚洲国产日韩一区二区| 日韩人妻精品一区2区三区| 99久久99久久久精品蜜桃| 午夜福利免费观看在线| 亚洲,一卡二卡三卡| 欧美人与善性xxx| 亚洲av电影在线观看一区二区三区| 久久国产精品男人的天堂亚洲| 日韩伦理黄色片| 亚洲中文字幕日韩| 伊人亚洲综合成人网| 国产精品一区二区精品视频观看| 国产主播在线观看一区二区 | 欧美+亚洲+日韩+国产| 黄色片一级片一级黄色片| 女人被躁到高潮嗷嗷叫费观| 亚洲精品日韩在线中文字幕| 又黄又粗又硬又大视频| 国产亚洲欧美在线一区二区| 丰满迷人的少妇在线观看| av视频免费观看在线观看| 午夜激情久久久久久久| 午夜免费成人在线视频| 极品人妻少妇av视频| 91字幕亚洲| 国产国语露脸激情在线看| 欧美+亚洲+日韩+国产| 超色免费av| 国产日韩欧美视频二区| 日韩中文字幕欧美一区二区 | 久久国产亚洲av麻豆专区| kizo精华| 青春草亚洲视频在线观看| 91麻豆精品激情在线观看国产 | 久久国产精品大桥未久av| 两个人看的免费小视频| 日本91视频免费播放| 青草久久国产| 精品福利观看| 精品国产乱码久久久久久男人| 高清不卡的av网站| 99re6热这里在线精品视频| 中国国产av一级| 两个人免费观看高清视频| 欧美日韩亚洲国产一区二区在线观看 | 我的亚洲天堂| 国产欧美日韩一区二区三 | 日本午夜av视频| 欧美亚洲日本最大视频资源| 国产黄频视频在线观看| 亚洲欧美精品综合一区二区三区| 大话2 男鬼变身卡| 国产精品久久久久成人av| 中文字幕另类日韩欧美亚洲嫩草| 国产成人a∨麻豆精品| 欧美黑人欧美精品刺激| 精品国产国语对白av| 天天操日日干夜夜撸| 蜜桃国产av成人99| 狂野欧美激情性xxxx| 久久天堂一区二区三区四区| 色播在线永久视频| 天天躁夜夜躁狠狠久久av| 欧美 日韩 精品 国产| 少妇的丰满在线观看| 宅男免费午夜| av国产久精品久网站免费入址| 母亲3免费完整高清在线观看| 每晚都被弄得嗷嗷叫到高潮| 91麻豆av在线| 久久久久精品人妻al黑| 99国产精品一区二区三区| 一区二区三区精品91| 99久久人妻综合| svipshipincom国产片| 欧美日本中文国产一区发布| 我要看黄色一级片免费的| 欧美日韩成人在线一区二区| av在线app专区| 久久久精品94久久精品| 午夜福利免费观看在线| 亚洲精品一卡2卡三卡4卡5卡 | 国产97色在线日韩免费| 国产精品av久久久久免费| 一本大道久久a久久精品| 黄色 视频免费看| 日韩av在线免费看完整版不卡| 日日摸夜夜添夜夜爱| 精品人妻熟女毛片av久久网站| 真人做人爱边吃奶动态| 久久精品人人爽人人爽视色| 好男人视频免费观看在线| 19禁男女啪啪无遮挡网站| 91精品伊人久久大香线蕉| 天天添夜夜摸| 亚洲国产最新在线播放| 天天添夜夜摸| 女人爽到高潮嗷嗷叫在线视频| 色网站视频免费| 又黄又粗又硬又大视频| 肉色欧美久久久久久久蜜桃| 9191精品国产免费久久| 美女脱内裤让男人舔精品视频| 大片免费播放器 马上看| 美女脱内裤让男人舔精品视频| 亚洲精品国产av成人精品| 天堂中文最新版在线下载| 无限看片的www在线观看| 亚洲欧美精品自产自拍| 亚洲国产精品成人久久小说| 久久久亚洲精品成人影院| 99热全是精品| 欧美精品一区二区免费开放| 丝袜在线中文字幕| 午夜老司机福利片| 在线观看人妻少妇| 自拍欧美九色日韩亚洲蝌蚪91| 国产av一区二区精品久久| 青春草亚洲视频在线观看| 亚洲成av片中文字幕在线观看| 亚洲情色 制服丝袜| 亚洲精品美女久久久久99蜜臀 | 男女午夜视频在线观看| 精品欧美一区二区三区在线| 久久久精品免费免费高清| 国产人伦9x9x在线观看| 亚洲成av片中文字幕在线观看| 国产成人a∨麻豆精品| 90打野战视频偷拍视频| 午夜91福利影院| 亚洲欧美一区二区三区黑人| 狂野欧美激情性bbbbbb| 成人手机av| 国产亚洲一区二区精品| 久9热在线精品视频| 丝袜在线中文字幕| 国产熟女午夜一区二区三区| 18禁观看日本| 国产成人欧美在线观看 | 免费少妇av软件| 老司机影院毛片| 中文字幕精品免费在线观看视频| 日本a在线网址| 黄色视频在线播放观看不卡| 1024视频免费在线观看| 丝袜在线中文字幕| 色婷婷av一区二区三区视频| 男人添女人高潮全过程视频| 一级毛片电影观看| 一二三四社区在线视频社区8| 9191精品国产免费久久| 亚洲七黄色美女视频| 水蜜桃什么品种好| 日本av免费视频播放| 日韩制服丝袜自拍偷拍| 可以免费在线观看a视频的电影网站| 免费少妇av软件| 欧美大码av| 国产亚洲欧美在线一区二区| 国产精品久久久久成人av| 老鸭窝网址在线观看| 免费看av在线观看网站| 亚洲伊人久久精品综合| 国产免费又黄又爽又色| 国产精品一区二区在线不卡| 黄色一级大片看看| 国产精品久久久久成人av| 男的添女的下面高潮视频| 欧美日韩黄片免| 精品国产乱码久久久久久小说| 国产野战对白在线观看| 丰满人妻熟妇乱又伦精品不卡| 亚洲久久久国产精品|