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

    地磁感應(yīng)電流研究綜述

    2022-11-08 08:24:52步宛燁張慧茜黃清華
    關(guān)鍵詞:方法模型

    步宛燁,張慧茜,黃清華,2*

    1北京大學(xué)地球與空間科學(xué)學(xué)院地球物理學(xué)系,北京100871

    2河北紅山地球物理國家野外科學(xué)觀測研究站,北京大學(xué),北京100871

    0 引言

    太陽活動所引起的地磁擾動會在地球表面引起感生地電場,從而在地面導(dǎo)體系統(tǒng)中產(chǎn)生低頻電流,稱為地磁感應(yīng)電流(geomagnetically induced currents,GICs).

    GICs作為造成空間天氣災(zāi)害的重要因素之一,會對導(dǎo)體設(shè)施產(chǎn)生威脅.早在一個多世紀(jì)前,人類就觀測到了GICs對電報通信系統(tǒng)的影響(Barlow,1849).1859年9月1日至2日觀測到的卡林頓超級磁暴(Carrington,1859),甚至曾導(dǎo)致電報線起火(Boteler,2006).除了干擾通信系統(tǒng),GICs還能腐蝕油氣管道,給油氣運輸帶來長期的潛在威脅(Allison and Huddleston,1952;Pulkkinen et al.,2001).GICs也會影響軌道的自動控制系統(tǒng),對鐵路系統(tǒng)產(chǎn)生危害.例如,1982年7月,瑞典鐵路的軌道電路受到GICs的影響而出現(xiàn)信號燈“閃紅”現(xiàn)象(Wik et al.,2009).1989年3月的磁暴使前蘇聯(lián)Gorky鐵路的信號控制系統(tǒng)出現(xiàn)了異常(Eroshenko et al.,2010).隨著電力系統(tǒng)的不斷發(fā)展,電網(wǎng)系統(tǒng)成為受GICs影響最大的導(dǎo)體系統(tǒng).早在1940年,美國和加拿大電網(wǎng)受GICs干擾的現(xiàn)象就 已 被 發(fā) 現(xiàn)(Davidson,1940; McNish,1940).1989年3月13日,在強GICs的干擾下,加拿大魁北克省的電力系統(tǒng)出現(xiàn)大范圍的電壓崩潰,約600萬人持續(xù)斷電9個小時. 美國新澤西州的一個發(fā)電機(jī)的升壓變壓器也在該事件中被燒毀(Allen et al., 1989;Boteler, 2001;Bolduc,2002).

    事實上,GICs甚至可以造成更大程度的災(zāi)害.1859年的卡林頓級別的超級磁暴,是1989年導(dǎo)致魁北克斷電的磁暴強度的1.5倍(Siscoe et al.,2006).如果卡林頓級別的磁暴發(fā)生在導(dǎo)體系統(tǒng)高度發(fā)達(dá)的今天,可能會引起難以預(yù)估的后果(Kelbert,2020).2012年7月23日發(fā)生過一次強太陽活動事件,如果該事件襲擊地球,引起的感應(yīng)地電場強度甚至可能會超過1989年的磁暴(Ngwira et al.,2013).基于概率估計的研究也表明卡林頓級別的地磁暴存在重復(fù)發(fā)生的可能,例如,Riley(2012)基于冪律分布得出結(jié)論:卡林頓磁暴事件(Dst

    在早期,研究者們普遍認(rèn)為只有高緯度地區(qū)會遭受GICs的災(zāi)害,比如加拿大(Boteler,2001)、瑞典(Wik et al.,2008,2009)、芬蘭(Pirjola et al.,2007)等,因此早期的GICs研究主要集中在高緯度區(qū)域.然而,隨著電力系統(tǒng)的發(fā)展,在中低緯度地區(qū)也出現(xiàn)了GICs影響電網(wǎng)的事件,比如中國(劉連光等,2008)、南非(Koen and Guant,2003)、西班牙(Torta et al.,2012)和新西蘭(Béland and Small,2004).因此GICs研究在低緯度地區(qū)也開始發(fā)展,比如中國(Liu et al.,2009;Zhang et al.,2015,2020)、日本(Watari et al.,2009;Watari, 2015)、新 西 蘭(Marshall et al.,2012;Divett et al.,2017,2018)、巴西(Trivedi et al.,2007;Barbosa et al.,2015)、南非(Matandirotya et al.,2015,2016;Lotz and Danskin,2017)、澳大利亞(Marshall et al.,2019,2020)和西班牙(Torta et al.,2012,2014,2017)等.

    本文將對GICs的研究進(jìn)行概述.首先簡述空間天氣對GICs的驅(qū)動原理,然后分別介紹計算GICs的三個步驟的相關(guān)原理和研究進(jìn)展,包括地球表面地磁場的重建,感應(yīng)地電場的計算,以及地面導(dǎo)體系統(tǒng)中GICs的計算.最后對當(dāng)前GICs的研究方向加以總結(jié),并對其未來研究提出展望.

    1 驅(qū)動GICs的空間天氣鏈

    空間天氣這一概念首次出現(xiàn)于1950年代,并在1990年代普遍使用(Cade and Chan-Park,2015).空間天氣是指由太陽內(nèi)部、太陽大氣、太陽風(fēng)以及地球磁層、電離層和熱層組成的復(fù)雜系統(tǒng)的物理條件(Abda et al.,2020).從空間天氣的角度,GICs實際上是太陽—太陽風(fēng)—磁層—電離層—地球這一空間天氣鏈在地表的最終體現(xiàn)(Pulkkinen et al.,2017),如圖1所示.

    圖1 地磁感應(yīng)電流的“太陽—太陽風(fēng)—磁層—電離層—地球”空間天氣鏈(修改自Abda et al.,2020)Fig.1 The"sun-solar wind-magnetosphere-ionosphere-earth"space weather chain of geomagnetic induced currents(modified from Abda et al.,2020)

    近幾十年來,人們對太陽、太陽風(fēng)、磁層、電離層以及固體地球,特別是它們之間的耦合關(guān)系的認(rèn)識取得了重大進(jìn)展,因而對驅(qū)動GICs的關(guān)鍵過程的了解也日趨成熟(Pulkkinen et al.,2017).劇烈的太陽活動,例如太陽黑子、耀斑爆發(fā)、日冕物質(zhì)拋射(CME)等,會導(dǎo)致大量等離子體噴出,并且形成攜帶著巨大能量的太陽風(fēng).當(dāng)太陽風(fēng)高速運動并傳播到地球附近后,會與地球磁場發(fā)生復(fù)雜的相互作用,從而影響磁層以及電離層內(nèi)的物理過程,并產(chǎn)生一系列電流體系.這些電流源驅(qū)動了地球表面的電磁感應(yīng)現(xiàn)象,從而導(dǎo)致地磁擾動,形成地磁暴. 時變的磁場會在地表感應(yīng)出地電場,進(jìn)而在導(dǎo)體系統(tǒng)中產(chǎn)生GICs(劉春明,2009;Abda et al.,2020).

    地磁感應(yīng)電流是由地球表面的地磁擾動直接驅(qū)動的,因而GICs的計算通常從地磁場出發(fā),首先根據(jù)地磁觀測或空間天氣參數(shù)重建地球表面的地磁場,然后根據(jù)大地電導(dǎo)率結(jié)構(gòu)計算地表的感應(yīng)地電場,最終基于導(dǎo)體系統(tǒng)的參數(shù)和結(jié)構(gòu)計算其中流動的GICs.下面將分別介紹GICs計算的三個步驟的相關(guān)原理和發(fā)展.

    2 地磁場的重建

    2.1 地磁場的地面觀測和插值

    計算GICs的第一步是重建地球表面的地磁場.要獲得地面總地磁場變化的信息,最直接的方法就是地面觀測.人類對地磁場的觀測已經(jīng)開展了一個多世紀(jì).目前,很多國家都布設(shè)了長期的地磁臺站,形成了全球覆蓋的地磁觀測網(wǎng).但是,地磁觀測臺站的數(shù)目在全球范圍內(nèi)依舊稀少,并且分布不均,有些臺站的記錄也有所缺失,需要通過地球物理的方法對其在空間上進(jìn)行插值.

    最常用的地磁記錄插值方法是Amm和Viljanen(1999)提出的SECS(spherical elementary current systems)方法,該方法利用電離層電流中的無旋分量(不影響地球表面磁場,因而可被忽略)這一性質(zhì),將無散分量作為等效電流,根據(jù)已知的地磁場觀測數(shù)據(jù),計算無散的等效電離層電流分布,如圖2所示.而基于畢奧—薩伐爾定律,如果能確定磁層和電離層的電流體系,就可以準(zhǔn)確地計算出地面上相應(yīng)的外部地磁場的變化(Rast?tter et al.,2014).相比于傳統(tǒng)的傅里葉變換、球諧分析法和球冠諧分析法,SECS方法不需要對頻譜進(jìn)行截斷,并且可以任意設(shè)置等效電流系統(tǒng)的位置以適應(yīng)地磁測點的分布和所要研究問題的電流體系(Amm and Viljanen,1999).目前大部分運用地磁觀測數(shù)據(jù)來計算GICs的研究都是采用SECS方法進(jìn)行插值(例如,Caraballo et al.,2020;Lucas et al.,2020;Wang et al.,2020;Alves Ribeiro et al.,2021).

    圖2 SECS示意圖.(a)無旋的電流體系(與相關(guān)的場向電流);(b)無散的電流體系(修改自Amm and Viljanen,1999)Fig.2 Sketch of spherical elementary current systems(SECS).(a)Curl-free elementary system(with associated FACs);(b)Divergence-free elementary system(modified from Amm and Viljanen,1999).

    2.2 基于空間天氣的地磁場模擬

    除了地面觀測的方式,地球表面的地磁場還可以根據(jù)第一性原理,采用全球三維磁流體力學(xué)(magnetohydrodynamics, MHD)模型,根據(jù)空間天氣參數(shù)來模擬(Pulkkinen et al.,2013;Glocer et al.,2016).基于太陽風(fēng)和行星際磁場等參數(shù),利用結(jié)合流體力學(xué)方程和麥克斯韋方程組而得到的MHD模型,可以直接預(yù)測磁層和電離層的電流體系,進(jìn)而計算地表磁場的擾動情況.比較常用的MHD模型有BATS-R-US(the block-adaptive-treesolarwind-roe-upwind-scheme)模型(Powell et al.,1999)、GUMICS(grand unified ionospheremagnetosphere coupling simuation)模型(Palmroth et al.,2001)、OpenGGCM(open geospace general circulation model)模型(Raeder et al.,2001)、LFM(the lyon–fedder–mobarry)模型(Lyon et al., 2004)、SWMF(space weather modeling framework)模型(Tóth et al.,2005)和PPMLR-MHD(lagrangian version of piecewise parabolic method-magnetohydrodynamic)模型(Hu et al.,2005)等.

    目前已有很多GICs的相關(guān)研究利用MHD模型來計算地磁場變化.Pulkkinen等(2013)選取了12個地磁臺站的數(shù)據(jù),分別采用兩種經(jīng)驗?zāi)P停╓eimer和Weigel經(jīng)驗?zāi)P停┮约叭NMHD模型(LFM-MIX,OpenGGCM和SWMF模型),對6次地磁暴事件的地磁擾動進(jìn)行預(yù)測,計算了地磁場的時間導(dǎo)數(shù)(dB/dt),并通過列聯(lián)表評估各個模型的性能.Ngwira等(2014)利用全球三維MHD模型(SWMF模型)計算了卡林頓磁暴時的磁層—電離層響應(yīng)以及地磁場和感應(yīng)地電場.Zhang等(2020)利用SWMF模型與一維電導(dǎo)率模型,模擬了磁暴期間我國低緯度變電站的地磁擾動以及GICs. 模擬結(jié)果表明,該模型對Bx分量的擬合效果好于By分量.

    盡管在這一領(lǐng)域已經(jīng)取得了很多進(jìn)展,但要準(zhǔn)確地描述全球磁層和電離層電流的變化仍然非常有挑戰(zhàn)性.驅(qū)動地電場的地磁擾動的時空差異性很大,外部電流體系的微小變化就可能導(dǎo)致對局部地磁場預(yù)測的較大誤差,從而影響地電場和GICs的模擬(Ngwira et al.,2015;Pulkkinen et al., 2015,2017).另外,某些物理現(xiàn)象,比如亞暴電流體系或赤道電射流等,無法用單一的理想MHD模型來描述(Ngwira et al.,2014).同時,MHD方法依賴于對空間天氣參數(shù)的觀測,包括行星際磁場的三分量值、太陽風(fēng)的速度、密度和溫度等,因此未來更高精度的空間天氣探測計劃可能會有助于提高M(jìn)HD模擬的精度.

    3 感應(yīng)地電場的計算

    地電場是評估災(zāi)害和研究極端天氣事件的重要參數(shù),也是GICs研究的關(guān)鍵.在重建了地面地磁場以后,計算GICs的第二步是基于大地電阻率模型,通過重建的地磁場變化來計算地電場的時空變化.

    3.1 基于平面波假設(shè)的地電場計算方法

    幾十年來,Cagniard(1953)提出的“平面波法”一直是計算感應(yīng)地電場的基礎(chǔ)方法. 該方法假設(shè)來自磁層—電離層源的磁場可以近似為平面波.根據(jù)法拉第電磁感應(yīng)定律和麥克斯韋方程,地電場E可以表示成是通過地磁場B與大地阻抗張量Z在時域的卷積.經(jīng)過傅里葉變換,可以轉(zhuǎn)換為在頻率域的乘積(Lucas et al.,2020):

    式中,μ是磁導(dǎo)率,E是水平地電場,B是水平地磁場,Z是大地阻抗張量.

    目前應(yīng)用最廣泛的方法是一維平面波法,認(rèn)為地球電導(dǎo)率分布可以近似為僅隨深度變化的一維結(jié)構(gòu),從而采用一維半空間層狀電導(dǎo)率結(jié)構(gòu)模型來計算地表阻抗Z.

    阻抗張量Z在一維情況下可以寫為:

    將Z1D代入公式(1)即可計算地電場的水平分量.

    本文我們通過一些模擬來探究一維層狀電導(dǎo)率結(jié)構(gòu)和地磁暴頻率對感應(yīng)地電場的影響.假設(shè)地下電導(dǎo)率結(jié)構(gòu)由兩個厚度均為10 km的均勻?qū)蛹暗撞康木鶆虬肟臻g組成(圖3).取半空間電導(dǎo)率σ3=0.01 S/m保持不變,頂部兩層電導(dǎo)率分別在0.000 1~1 S/m范圍內(nèi)變化,當(dāng)改變某一層的電導(dǎo)率時,另一層電導(dǎo)率固定在0.01 S/m.

    圖3 三層大地電導(dǎo)率模型Fig.3 Three-layer earth conductivity model

    我們?nèi)〈疟╊l率分別為0.01 Hz、0.001 Hz和0.000 1 Hz,并假設(shè)地磁場的By分量為500 nT,分別計算對應(yīng)的地面地電場Ex分量的絕對值,結(jié)果如圖4所示.

    從圖4可以看出,地電場會隨任一層中電導(dǎo)率的增大而降低,并且高頻磁場會感生出更大幅度的地電場. 因此,高頻短周期的磁暴會在低電導(dǎo)率區(qū)域產(chǎn)生更強的地電災(zāi)害.同時,對比同一頻率的計算結(jié)果,當(dāng)兩層電導(dǎo)率均在0.000 1~1 S/m范圍內(nèi)變化時,地電場值隨著淺層電導(dǎo)率的變化更明顯,因此淺層電導(dǎo)率對地電場的影響更大.隨著頻率的降低,兩層的電導(dǎo)率對地電場的影響差異逐漸減小.因此,高頻磁暴期間所造成的地電災(zāi)害主要受到淺層電導(dǎo)率的影響.

    圖4 地電場隨第一、二層電導(dǎo)率以及頻率的變化曲線.藍(lán)線、紅線和橘線的頻率分別為0.01 Hz、0.001 Hz和0.000 1 Hz.實線表示隨著第一層電導(dǎo)率變化時的地電場,虛線表示隨著第二層電導(dǎo)率變化時的地電場Fig.4 The variation curves of the geoelectric field value with the conductivity and frequency of the first layer and the second layer.The blue, red and orange lines have frequencies of 0.01 Hz,0.001 Hz and 0.000 1 Hz,respectively.The solid lines show the geoelectric field as the conductivity of the first layer changes,and the dashed lines show the geoelectric field as the conductivity of the second layer changes

    上述模擬中假定地表磁場是均勻的,但在實際情況下,基于麥克斯韋方程,大地電導(dǎo)率模型的區(qū)域性差異會導(dǎo)致感應(yīng)地磁場的空間變化,因此利用均勻的地磁場來模擬磁暴在空間變化的電導(dǎo)率模型下的電磁響應(yīng)會產(chǎn)生誤差.為了探究電阻率模型的影響,我們考慮感應(yīng)地磁場的影響,根據(jù)波的傳播和衰減的特性,將地磁場和地電場看成是由入射波(Ei和Bi)和反射波(Er和Br)組成,將地磁場和地電場寫成如下形式(李金銘,2005):

    基于上述公式,選用與圖3同樣的電阻率模型,取入射磁場By為500 nT,計算頻率分別為0.01 Hz、0.001 Hz和0.000 1 Hz的地電場,如圖5所示.

    從圖5可以看出,地電場依然保持著隨電導(dǎo)率增大而降低的規(guī)律,并且地電場的幅值受到淺層電阻率的影響更大.同時,對比不同頻率,高頻磁場會感生出更大的地電場.因此,我們依舊得到了與方法一(不考慮感應(yīng)場)相同的結(jié)論,即在一維層狀介質(zhì)和平面波假設(shè)下,高頻磁暴和低電導(dǎo)率結(jié)構(gòu)會導(dǎo)致更強的地電災(zāi)害,并且相對淺層的電導(dǎo)率結(jié)構(gòu)對電場值的影響更大.

    圖5 考慮感應(yīng)場時,地電場隨第一、二層電導(dǎo)率以及頻率的變化曲線.藍(lán)線、紅線和橘線的頻率分別為0.01 Hz、0.001 Hz和0.000 1 Hz.實線表示隨著第一層電導(dǎo)率變化時的地電場,虛線表示隨著第二層電導(dǎo)率變化時的地電場Fig.5 The variation curves of geoelectric field value with the conductivity and frequency of the first layer and the second layer considering the induction field.The blue,red and orange lines have frequencies of 0.01 Hz,0.001 Hz and 0.000 1 Hz,respectively.The solid lines show the geoelectric field as the conductivity of the first layer changes,and the dashed lines show the geoelectric field as the conductivity of the second layer changes

    我們進(jìn)一步對比兩種方法所計算的地電場差異,當(dāng)我們考慮感應(yīng)地磁場(方法二)時,地電場隨第一層電導(dǎo)率的變化幅度遠(yuǎn)遠(yuǎn)超過了不考慮感應(yīng)場(方法一)的情況;而方法二計算得到的地電場隨著第二層電導(dǎo)率的變化幅度低于方法一. 這一結(jié)果表明,當(dāng)考慮地下電導(dǎo)率模型對感應(yīng)地磁場的影響時,淺部的低導(dǎo)層甚至?xí)?dǎo)致更嚴(yán)重的地電災(zāi)害,而較深部的高導(dǎo)結(jié)構(gòu)對地電場的影響相對較小.

    除了最簡單的一維半空間模型,還可以采用一維分段電導(dǎo)率模型來計算地電場.比較常用的分段模型有歐洲的一維模型(ádám et al.,2012)和美國的一維Fernberg(2012)電阻率模型,如圖6.Viljanen等(2012)采用歐洲的一維電導(dǎo)率模型,并引入球面坐標(biāo),計算了歐洲多個國家的地電場和GICs以驗證方法的有效性.Wei等(2013)使用Fernberg模型估算了歷史上兩次地磁暴時美國和加拿大南部的地電場,并得到結(jié)論:高電阻率結(jié)構(gòu)的低緯度區(qū)域也面臨著GICs風(fēng)險.一維分段模型的優(yōu)點是可以較好地描述區(qū)域電阻率分布,提高地電場模擬在較大范圍內(nèi)的精度,但是一維分段方法本質(zhì)上與傳統(tǒng)一維方法的區(qū)別不大,并不能描述真實電導(dǎo)率結(jié)構(gòu)的橫向變化.

    圖6 一維分段電導(dǎo)率模型.(a)歐洲(ádám et al.,2012,http://real.mtak.hu/id/eprint/2957);(b)美國(修改自Fernberg,2012)Fig.6 One-dimensional piecewise conductivity models of (a)Europe(ádám et al.,2012,http://real.mtak.hu/id/eprint/2957);(b)The United States (modified from Fernberg,2012)

    盡管上述提到的一維平面波法已經(jīng)有很多應(yīng)用,但是由于真正的地面電導(dǎo)率結(jié)構(gòu)是三維橫向變化的,因此采用簡單的一維近似可能會在電性結(jié)構(gòu)突變處產(chǎn)生較大的誤差.董博等(2015)采用有限元法模擬出在大地電導(dǎo)率橫向突變處感應(yīng)地電場會出現(xiàn)鄰近效應(yīng).Honkonen等(2018)認(rèn)為三維模型與局部一維模型所計算的地電場在高緯度海岸分界處會產(chǎn)生巨大差異.所以,基于三維電導(dǎo)率模型的地電場模擬是GICs研究的重要方向(Püthe et al.,2014).

    三維結(jié)構(gòu)中的阻抗張量可以寫為:只要已知三維阻抗張量Z3D,就可以代入公式(2)計算地電場水平分量.

    獲得Z3D的方法通常有兩種.第一種方法是已知研究區(qū)域的三維電導(dǎo)率模型,通過大地電磁的正演模擬來計算阻抗張量(例如,Wang et al.,2016).然而,詳細(xì)的三維電導(dǎo)率模型在大部分地區(qū)是未知的.全三維模型的有效替代是三維薄板模型,即在一維層狀模型上添加橫向變化的電導(dǎo)率薄層(Vasseur and Weidelt,1977).為了探究非一維電導(dǎo)率分布的影響,薄板模型已被應(yīng)用于感應(yīng)地電場和GICs的研究中(Beggan et al.,2013;Püthe and Kuvshinov,2013;Beggan,2015;Divett et al.,2017;Ivannikova et al.,2018).薄板模型可以考慮地球電導(dǎo)率的橫向變化,并模擬電導(dǎo)率結(jié)構(gòu)突變處的地電場異常.

    第二種方法是依據(jù)目標(biāo)區(qū)域已有的大地電磁測量的局部電場E和磁場B,直接計算該測點處的阻抗張量(例如,Lucas et al.,2020).局部大地電磁測量是獲得良好的區(qū)域阻抗張量,并將地磁場的變化映射到地電場的最佳方法(Pulkkinen et al.,2017).目前,在全球范圍內(nèi)已經(jīng)進(jìn)行了一些密集的MT勘探工作,比如北美的USArray項目(Schultz,2010)、澳大利亞的AusLAMP項目(Duan and Kyi,2018)和中國的SinoProbe計劃(Dong et al.,2013).但許多地區(qū)并沒有足夠密集的MT觀測數(shù)據(jù),這也限制了三維阻抗在計算地電場中的應(yīng)用.

    上述介紹的平面波方法大多使用了頻率域的公式,但時域公式也已得到開發(fā)和應(yīng)用(Viljanen and Pirjola,1989;Marti et al.,2014;Kelbert et al.,2017).時域方法的開發(fā)有助于實現(xiàn)GICs的實時估計,為GICs預(yù)警提供理論基礎(chǔ).

    盡管目前大部分GICs研究都是基于平面波假設(shè),但是該方法忽略了電離層電流的分布及其強度差異,可能會導(dǎo)致在大磁暴期間對地電場估算的不準(zhǔn)確(Kelbert et al.,2017).此外,不同緯度的空間電流體系模式不同,目前普遍認(rèn)為在中低緯區(qū)域產(chǎn)生地磁擾動的重要因素是環(huán)電流,而在高緯度地區(qū)產(chǎn)生地磁擾動的主要因素是磁層亞暴導(dǎo)致的極光電集流(Zheng et al.,2013;Ngwira and Pulkkinen,2018).平面波假設(shè)可以很好地描述低緯度區(qū)域的電流源,但并不適合高緯度地區(qū).

    3.2 地電場的其他計算方法

    Püthe和Kuvshinov(2013)提出將全球電流源表示為在地球上方流動的橫向變化的球面薄層等效電流,采用球諧系數(shù)展開法進(jìn)行求解.由于頻率域內(nèi)的磁場和地電場可以表示為關(guān)于源的線性組合,因此可以將電流源表示為球諧系數(shù)的組合,并利用麥克斯韋方程計算出頻率域內(nèi)每個球諧系數(shù)對應(yīng)的“單位”磁場和“單位”電場的分布,之后再進(jìn)行線性疊加,從而計算出全球的地電場分布.他們在沿海地區(qū)的地磁臺站驗證了該方法,并討論了一維和三維薄板電導(dǎo)率模型所計算的地電場的區(qū)別.但是上述方法具有一定的局限性,他們假設(shè)的電流源只有大尺度結(jié)構(gòu),因此僅使用次數(shù)和階次最多為3的球諧函數(shù)來描述源,僅適合于中低緯度區(qū)域的電磁模擬,并不足以描述高緯度區(qū)域的電離層電流源.

    Honkonen等(2018)沒有采用地面地磁場觀測數(shù)據(jù),而是結(jié)合了MHD方法(SWMF模型),利用太陽風(fēng)參數(shù)模擬地球表面的外源磁場的時空分布,再采用上述Püthe和Kuvshinov(2013)提出的球諧系數(shù)展開法模擬出薄層等效電流源,并計算全球地電場的時空分布.他們得到結(jié)論:考慮地球三維電導(dǎo)率結(jié)構(gòu)會導(dǎo)致橫向電導(dǎo)率驟變區(qū)域的地電場顯著增強(比如沿海區(qū)域),尤其是對于緯度較高的地區(qū).

    Marshalko等(2021)針對高緯度地區(qū),對比了基于平面波法、SECS法和SWMF模型的三維電磁場的模擬結(jié)果.圖7展示了在ABK臺站處采用三種方法所計算的三維地磁場和地電場的模擬值以及觀測值.他們認(rèn)為相比于SWMF模型,使用SECS方法構(gòu)建的源可以更準(zhǔn)確地再現(xiàn)研究區(qū)域的磁場擾動.使用SECS方法與平面波法所計算的地電場之間存在差異,特別是在電導(dǎo)率橫向?qū)Ρ葟娏业牡貐^(qū)(例如海岸)和更高緯度的區(qū)域.

    圖7 采用三種方法計算的三維電磁場的模擬值以及觀測值.(a)~(c)Abisko(ABK)地磁觀測站于2017年9月7日20:00 UT至2017年9月8日03:59 UT地磁場的模擬結(jié)果和觀測值.(d)~(e)地電場的模擬結(jié)果.黑線表示地磁場的觀測值;綠線表示基于平面波法的三維地電場模擬值;紫線表示基于SECS法的三維電磁場模擬值;紅線表示基于SWMF模型的三維電磁場模擬值(修改自Marshalko et al.,2021)Fig.7 The time series of the modeled and observed 3-D electromagnetic field calculated by three methods.(a)~(c)The modeled and observed time series of geomagnetic field at Abisko(ABK)geomagnetic observatory from 20:00 UT,September 7,2017 to 03:59 UT,September 8,2017.(d)~(e)The modeled time series of geoelectric field.The black lines represent the observed value of the geomagnetic field,the green lines represent the modeled time series of the 3-D geomagnetic field based on the plane-wave method,the purple lines represent the modeled time series of the 3-D electromagnetic field based on the SECS method,and the red lines represent the modeled time series of the 3-D electromagnetic field based on the SWMF model(modified from Marshalko et al.,2021)

    4 導(dǎo)體系統(tǒng)中GICs的計算

    在完成了地面地磁場重建并計算了感應(yīng)地電場后,就可以利用地電場來計算導(dǎo)體系統(tǒng)中的GICs值.地面導(dǎo)體系統(tǒng)可以分為離散接地(例如電網(wǎng))和連續(xù)接地(例如油氣管道)兩種類型.前者在變壓器處有離散的接地點,后者則通過導(dǎo)電性很小的絕緣材料連續(xù)接地.鐵路系統(tǒng)則是兩者兼具:鐵軌中的電力設(shè)備是離散接地,而鐵軌是連續(xù)接地(Pirjola,2002).

    兩種接地類型需要采用不同的方法來計算GICs:離散接地是利用矩陣的形式來計算(Lehtinen and Pirjola,1985),而連續(xù)接地則是采用分布式源傳輸線理論(Distributed Source Transmission Line,DSTL)計算(Boteler, 1997).

    4.1 離散接地系統(tǒng)中的GICs計算

    對于離散接地系統(tǒng),比如電網(wǎng),地磁暴期間產(chǎn)生的GICs在電網(wǎng)中流通的等效模型如圖8所示.

    圖8 GICs在電網(wǎng)中流通的等效模型(修改自劉同同等,2011)Fig.8 Equivalent model of GICs in power grids(modified from Liu et al.,2011)

    Lehtinen和Pirjola(1985)提出流入變壓器中性點的電流可以寫為矩陣的形式:

    式中,1是單位矩陣,Y和Z分別是網(wǎng)絡(luò)導(dǎo)納矩陣和接地阻抗矩陣,J是各節(jié)點的入地GICs電流組成的列矩陣.導(dǎo)納矩陣Y利用網(wǎng)絡(luò)中的導(dǎo)體電阻率來計算:

    式中,Rij和Rik表示兩節(jié)點之間的等效電阻.

    Vki是節(jié)點k和i之間輸電線路上的等效電壓源,可以表示成隨時間變化的沿輸電線路的地電場積分:

    若定義Ji為從電網(wǎng)流入節(jié)點i的所有電流Iki之和.根據(jù)歐姆定律,Ji可以寫為:

    將公式(6)和(8)代入(5)中,即可以計算變壓器中性點的GICs.已有很多基于上述公式的GICs模擬(例如,Beggan,2015;Marshall,et al.,2017;Viljanen and Pirjola,2017).

    從以上公式中可以看出,GICs是地電場沿輸電線路積分的時間函數(shù).GICs的計算并不依賴于單點的電場,而是取決于輸電線路上電場的積分.正是因為這種積分效應(yīng),盡管感應(yīng)地電場的空間差異很大(Bedrosian and Love,2015),真實的三維地電場所產(chǎn)生的線電壓仍可以用一維平均電場的等效電壓來近似(Viljanen et al.,2006,2012).因此,一維近似電阻率結(jié)構(gòu)對于大部分GICs研究而言是足夠準(zhǔn)確的.然而,這也并不意味著一維與三維方法是完全等價的,一維的方法有其固有的局限性,在橫向電導(dǎo)率驟變的區(qū)域,比如海岸,由于海陸電阻率差異而導(dǎo)致的海岸效應(yīng)無法通過輸電線積分來描述(Honkonen et al.,2018; Ivannikova et al.,2018;Liu et al.,2018).因此,采用三維電阻率的方法可以有效地改進(jìn)特定地質(zhì)環(huán)境下GICs的模擬效果.

    在實際的工程學(xué)計算中,通常會采用近似方法來計算電網(wǎng)中的GICs.在平面波假設(shè)下,如果假定地下電阻率結(jié)構(gòu)是橫向均勻的,那么GICs的計算可以簡化為關(guān)于電場水平分量的線性疊加(Viljanen and Pirjola,1994):

    式中,a和b系數(shù)取決于電力系統(tǒng)的電阻和幾何形狀.

    系數(shù)a和b通??梢愿鶕?jù)電網(wǎng)模型來估算(Wik et al.,2008;Viljanen et al.,2012),如果電網(wǎng)模型參數(shù)未知,但是有變壓器的GICs觀測值,則可以采用最小二乘法來擬合a和b系數(shù)(如

    Pulkkinen et al.,2007;Liu et al.,2009;Zhang et al.,2015,2020).但是需要注意的是,最小二乘法擬合有一定的局限性,由于變壓器的a和b系數(shù)由其各自的幾何參數(shù)和網(wǎng)絡(luò)配置決定,因此僅僅適用于有GICs測量值的變壓器,不能將擬合得到的系數(shù)應(yīng)用到周圍的其他變壓器上(Kelbert,2020).

    4.2 連續(xù)接地系統(tǒng)中的GICs計算

    對于連續(xù)接地系統(tǒng),比如油氣管道,由于管道與大地之間隔著絕緣涂層,管道內(nèi)的GICs可用分布源傳輸線理論來計算(Boteler,1997).如圖9所示,該理論將長距離的管道分割成許多小段,假設(shè)每小段管道的感應(yīng)電場是均勻的.可用阻抗Z和并聯(lián)導(dǎo)納Y組成的傳輸線電路來等效每段管道內(nèi)的電路(劉連光等,2015).

    圖9 管道GICs的分布源傳輸線模型(劉連光等,2015)Fig.9 Distribution source transmission line model of pipelines GICs(Liu et al.,2015)

    基于麥克斯韋方程可以建立管道內(nèi)的傳輸線方程,用于描述每一段管道中由均勻感應(yīng)電場E產(chǎn)生的電壓U和電流I的分布:

    電壓U的解析式為:

    目前已有很多基于分布源傳輸線理論的油氣管道GICs研究.Viljanen等(2006)在芬蘭創(chuàng)建了天然氣管道的監(jiān)測和實時預(yù)報系統(tǒng).Yu等(2019)首次實現(xiàn)了地磁場對管道的電磁干擾效應(yīng)的監(jiān)測,包括GICs和管地電位(pip-soil potential,PSP),并且認(rèn)為在中低緯度地區(qū),中小型磁暴就會對管道造成腐蝕.Ma和Liu(2019)建立了含漏點的管網(wǎng)模型,討論了不同漏點半徑和位置對GICs和PSP分布的影響.Liu等(2019)就我國油氣管道受GICs和管地電位的影響進(jìn)行了觀測實驗,證實了潮汐作用對磁暴中的GICs和PSP效應(yīng)的影響.

    5 總結(jié):GICs研究現(xiàn)狀與未來

    本文著重介紹了關(guān)于GICs計算方法的開發(fā)和應(yīng)用,除此之外,極值災(zāi)害分析也是GICs研究的重要內(nèi)容.極值災(zāi)害分析是基于統(tǒng)計模型,利用長期地磁觀測數(shù)據(jù)來預(yù)測地磁場或地電場極值(比如百年極值),從而評估該區(qū)域的地電災(zāi)害的研究方向.全球很多區(qū)域已開展了相關(guān)工作,比如加拿大(Boteler,2001;Nikitina,2016)、日本(Watari,2015)、英國(Beggan et al.,2013)、澳大利亞(Marshall et al.,2020)、西班牙(Torta et al.,2014)、南非(Lotz and Danskin,2017)以及美國(Love et al.,2016,2018a,2018b,2019;Lucas et al.,2020)等.

    對歷史事件的GICs分析只能評估區(qū)域的地電災(zāi)害風(fēng)險,GICs研究的最終目的是監(jiān)測、預(yù)報并規(guī)避極端天氣對基礎(chǔ)設(shè)施造成的災(zāi)害,因此對GICs的實時估計以及預(yù)報是研究的關(guān)鍵目標(biāo).GICs的實時估計需要基于時間域的感應(yīng)地電場模擬方法(Viljanen and Pirjola,1989;Marti et al.,2014;Kelbert et al.,2017),并且對計算效率有較高要求.GICs預(yù)報實際上是地磁場預(yù)報和GICs實時算法的結(jié)合,其中地磁場預(yù)報涉及到空間天氣的觀測和相關(guān)研究. 在拉格朗日L1點附近的先進(jìn)成分探測器(advanced composition explorer,ACE)和深太空氣候天文臺(deep space climate observatory,DSCOVR)觀測到的太陽風(fēng)數(shù)據(jù)比地面遭受磁暴的時間提前20 min至1小時(Kelbert,2020),可以用于MHD模型的輸入,并進(jìn)行GICs的短期預(yù)報(Pulkkinen et al.,2013;Glocer et al.,2016).此外,機(jī)器學(xué)習(xí)方法也是很有前景的GICs預(yù)報方法(如Bobra and Ilonidis,2016;Camporeale et al.,2018;Gruet et al.,2018;Camporeale,2019).

    盡管目前GICs的研究已經(jīng)取得了很大進(jìn)展,但仍然存在著許多重大的挑戰(zhàn).目前,我們對太陽風(fēng)—磁層—電離層耦合關(guān)系的了解還不夠透徹.明確誘發(fā)極端GICs的空間天氣因素并確定磁層—電離層對其的響應(yīng),有助于更好地模擬空間電流源,從而更準(zhǔn)確地計算感應(yīng)地電場和GICs.并且,我們對地電場和GICs的預(yù)測能力仍然非常有限.雖然MHD方法具備預(yù)測的能力,但是其模擬效果不夠精確,仍需進(jìn)一步改進(jìn).同時,高質(zhì)量的地球空間測量對于GICs研究至關(guān)重要.獲取更多的地球空間觀測也將有助于提高我們建模能力的準(zhǔn)確性.另外,量化GICs的危害需要多學(xué)科的交叉和合作,從地球物理學(xué)到電氣工程學(xué),如何更好地融合不同學(xué)科的知識,也是GICs研究領(lǐng)域的重要任務(wù).

    致謝

    感謝兩位審稿人和責(zé)任編輯對本文提出的建設(shè)性意見.感謝北京大學(xué)課程“地球物理研究的論文寫作指導(dǎo)”(課程號01201110)為本綜述的完成提供了支持.

    附中文參考文獻(xiàn)

    董博,王澤忠,劉連光,等.2015.大地電導(dǎo)率橫向突變處磁暴感應(yīng)地電場的鄰近效應(yīng)[J].地球物理學(xué)報,58(1):238-246.

    李金銘.2005.地電場與電法勘探[M].北京:地質(zhì)出版社.

    劉春明.2009.中低緯電網(wǎng)地磁感應(yīng)電流及其評估方法研究[D].北京:華北電力大學(xué).

    劉連光,劉春明,張冰,等.2008.中國廣東電網(wǎng)的幾次強磁暴影響事件[J].地球物理學(xué)報,2008(4):976-981.

    劉連光,張鵬飛,王開讓,等.2015.基于大地電導(dǎo)率分層模型的油氣管網(wǎng)地磁暴干擾評估方法[J].電網(wǎng)技術(shù),39(6):1556-1561.

    劉同同,劉連光,鄒明.2011.基于分層大地電導(dǎo)率模型的電網(wǎng)GIC算法研究[J].電網(wǎng)與清潔能源,27(5):26-30.

    猜你喜歡
    方法模型
    一半模型
    重要模型『一線三等角』
    重尾非線性自回歸模型自加權(quán)M-估計的漸近分布
    學(xué)習(xí)方法
    可能是方法不對
    3D打印中的模型分割與打包
    用對方法才能瘦
    Coco薇(2016年2期)2016-03-22 02:42:52
    FLUKA幾何模型到CAD幾何模型轉(zhuǎn)換方法初步研究
    四大方法 教你不再“坐以待病”!
    Coco薇(2015年1期)2015-08-13 02:47:34
    賺錢方法
    亚洲精品在线美女| 我的亚洲天堂| 51午夜福利影视在线观看| 成人亚洲精品一区在线观看| 18禁国产床啪视频网站| 日韩欧美一区二区三区在线观看| 午夜免费成人在线视频| 日韩欧美 国产精品| 国产一级毛片七仙女欲春2 | 中文亚洲av片在线观看爽| 日韩国内少妇激情av| 中文资源天堂在线| 久久精品亚洲精品国产色婷小说| 中文字幕高清在线视频| 色播在线永久视频| 91九色精品人成在线观看| 精品卡一卡二卡四卡免费| 精品高清国产在线一区| netflix在线观看网站| 在线国产一区二区在线| 成人一区二区视频在线观看| 91大片在线观看| 成人一区二区视频在线观看| 亚洲av中文字字幕乱码综合 | 午夜免费观看网址| 久久99热这里只有精品18| 日本五十路高清| av欧美777| 男人舔女人下体高潮全视频| 欧美黑人精品巨大| 男女午夜视频在线观看| 91成人精品电影| 色在线成人网| 色在线成人网| 国产又爽黄色视频| 亚洲中文字幕日韩| 不卡一级毛片| 国产精品亚洲美女久久久| 搞女人的毛片| 又黄又爽又免费观看的视频| 又黄又爽又免费观看的视频| 又黄又爽又免费观看的视频| 欧美激情久久久久久爽电影| 岛国在线观看网站| 一级a爱视频在线免费观看| 亚洲欧美精品综合一区二区三区| 一级黄色大片毛片| 亚洲 欧美 日韩 在线 免费| 身体一侧抽搐| 亚洲中文字幕一区二区三区有码在线看 | 欧美成人免费av一区二区三区| 国产亚洲欧美精品永久| 两人在一起打扑克的视频| 亚洲成国产人片在线观看| 国产成人av教育| 精品一区二区三区四区五区乱码| 老司机午夜十八禁免费视频| 日本一区二区免费在线视频| aaaaa片日本免费| 亚洲免费av在线视频| 夜夜夜夜夜久久久久| 成熟少妇高潮喷水视频| 日本三级黄在线观看| 两个人看的免费小视频| 久久热在线av| 国产av一区在线观看免费| 国产熟女午夜一区二区三区| 欧美日韩乱码在线| 国产精品自产拍在线观看55亚洲| 亚洲av电影在线进入| 90打野战视频偷拍视频| 欧美日韩乱码在线| 老鸭窝网址在线观看| 久久天躁狠狠躁夜夜2o2o| a在线观看视频网站| 叶爱在线成人免费视频播放| 国产伦人伦偷精品视频| 啪啪无遮挡十八禁网站| 精品免费久久久久久久清纯| 老鸭窝网址在线观看| 一卡2卡三卡四卡精品乱码亚洲| 免费看日本二区| 国产精品永久免费网站| 精品第一国产精品| 可以在线观看毛片的网站| 看黄色毛片网站| 欧美乱色亚洲激情| 亚洲av五月六月丁香网| 亚洲国产高清在线一区二区三 | 中文字幕人成人乱码亚洲影| 人人妻人人澡欧美一区二区| 国产不卡一卡二| 久久婷婷人人爽人人干人人爱| 亚洲第一青青草原| 天天添夜夜摸| 18美女黄网站色大片免费观看| 亚洲一卡2卡3卡4卡5卡精品中文| 亚洲免费av在线视频| 成人三级做爰电影| 制服人妻中文乱码| 一区二区三区国产精品乱码| 亚洲最大成人中文| 国产亚洲精品久久久久久毛片| 久久香蕉精品热| 午夜福利高清视频| 在线视频色国产色| 香蕉av资源在线| 国产精品一区二区三区四区久久 | 国产成人精品无人区| 18禁裸乳无遮挡免费网站照片 | 国产av一区在线观看免费| 脱女人内裤的视频| 男男h啪啪无遮挡| 亚洲精品一卡2卡三卡4卡5卡| 国产一级毛片七仙女欲春2 | 老司机午夜十八禁免费视频| 亚洲成a人片在线一区二区| 丰满人妻熟妇乱又伦精品不卡| 国产1区2区3区精品| 无限看片的www在线观看| 757午夜福利合集在线观看| 亚洲熟女毛片儿| 人人妻人人澡欧美一区二区| 国产一卡二卡三卡精品| 国产亚洲欧美在线一区二区| 欧美zozozo另类| 免费在线观看黄色视频的| 欧美中文日本在线观看视频| 成人18禁在线播放| 淫妇啪啪啪对白视频| 亚洲专区字幕在线| 亚洲熟女毛片儿| 亚洲自拍偷在线| 波多野结衣高清作品| 99国产综合亚洲精品| 成人手机av| xxxwww97欧美| 午夜久久久久精精品| 成人精品一区二区免费| 国产精品亚洲一级av第二区| 天天躁夜夜躁狠狠躁躁| 亚洲av成人不卡在线观看播放网| 亚洲国产精品sss在线观看| 国产不卡一卡二| a级毛片a级免费在线| 免费搜索国产男女视频| 非洲黑人性xxxx精品又粗又长| 国产片内射在线| 好看av亚洲va欧美ⅴa在| 午夜视频精品福利| 亚洲,欧美精品.| 国产成人精品无人区| 亚洲人成伊人成综合网2020| 日本在线视频免费播放| 麻豆一二三区av精品| 老司机深夜福利视频在线观看| 亚洲 欧美 日韩 在线 免费| 天天躁狠狠躁夜夜躁狠狠躁| 亚洲成人国产一区在线观看| 国产精品av久久久久免费| 久久久久国内视频| 久久久久久亚洲精品国产蜜桃av| 一本精品99久久精品77| 一区二区三区国产精品乱码| 婷婷精品国产亚洲av| 亚洲专区中文字幕在线| 久久亚洲精品不卡| 黄频高清免费视频| 国产真实乱freesex| 美女扒开内裤让男人捅视频| 看黄色毛片网站| 精品国产亚洲在线| 丁香欧美五月| 哪里可以看免费的av片| 69av精品久久久久久| 三级毛片av免费| 黄色a级毛片大全视频| 亚洲av五月六月丁香网| 婷婷精品国产亚洲av在线| 欧美午夜高清在线| 日韩视频一区二区在线观看| 丝袜人妻中文字幕| 日韩一卡2卡3卡4卡2021年| 免费在线观看亚洲国产| 久久午夜综合久久蜜桃| 99国产精品一区二区三区| 伊人久久大香线蕉亚洲五| 中文字幕人妻丝袜一区二区| 久久久久久久精品吃奶| 99精品在免费线老司机午夜| 99久久99久久久精品蜜桃| 国产亚洲精品久久久久久毛片| 91成人精品电影| 亚洲精品美女久久av网站| 国产伦人伦偷精品视频| 操出白浆在线播放| 一本精品99久久精品77| 岛国在线观看网站| ponron亚洲| 91麻豆精品激情在线观看国产| 国产亚洲精品久久久久5区| 女警被强在线播放| www.熟女人妻精品国产| 中文资源天堂在线| 久久精品91蜜桃| 男人舔女人下体高潮全视频| 丝袜美腿诱惑在线| 一进一出抽搐动态| 日本三级黄在线观看| 中文字幕精品免费在线观看视频| 亚洲专区中文字幕在线| 男女床上黄色一级片免费看| 午夜福利在线观看吧| 久久久久久亚洲精品国产蜜桃av| 国产精品电影一区二区三区| 国产精品亚洲一级av第二区| 欧美黑人精品巨大| 国产成人精品久久二区二区91| 久久久精品国产亚洲av高清涩受| 精品久久久久久久久久免费视频| 久热爱精品视频在线9| 亚洲真实伦在线观看| 欧美黑人欧美精品刺激| 免费人成视频x8x8入口观看| 国产区一区二久久| 国产成人精品久久二区二区免费| 每晚都被弄得嗷嗷叫到高潮| 国产精品免费一区二区三区在线| 日本免费a在线| 成年人黄色毛片网站| 午夜成年电影在线免费观看| 美国免费a级毛片| 日韩免费av在线播放| 岛国视频午夜一区免费看| 欧美性猛交╳xxx乱大交人| 99精品久久久久人妻精品| 欧美久久黑人一区二区| 国产精品1区2区在线观看.| 亚洲国产毛片av蜜桃av| or卡值多少钱| 操出白浆在线播放| 我的亚洲天堂| 久久国产精品影院| 亚洲片人在线观看| 亚洲精品中文字幕一二三四区| 欧美成人免费av一区二区三区| 国产精品一区二区精品视频观看| 黄片小视频在线播放| 后天国语完整版免费观看| 成人一区二区视频在线观看| www.自偷自拍.com| 国产精品日韩av在线免费观看| 这个男人来自地球电影免费观看| 欧美黄色淫秽网站| 精品久久久久久久人妻蜜臀av| 少妇粗大呻吟视频| 国产精品久久久人人做人人爽| 久久久久久久精品吃奶| 特大巨黑吊av在线直播 | 成人免费观看视频高清| 一级作爱视频免费观看| 成人18禁高潮啪啪吃奶动态图| 777久久人妻少妇嫩草av网站| 亚洲国产看品久久| 香蕉丝袜av| 久久这里只有精品19| 亚洲av电影不卡..在线观看| 欧美国产精品va在线观看不卡| 午夜福利18| 日韩欧美免费精品| 女人被狂操c到高潮| 欧美激情久久久久久爽电影| 国内揄拍国产精品人妻在线 | 欧美另类亚洲清纯唯美| 亚洲最大成人中文| 男女下面进入的视频免费午夜 | 精品久久久久久久久久久久久 | 美女免费视频网站| 日韩欧美免费精品| 在线av久久热| 两性午夜刺激爽爽歪歪视频在线观看 | 观看免费一级毛片| 亚洲精品中文字幕在线视频| 亚洲精品美女久久久久99蜜臀| bbb黄色大片| 超碰成人久久| 99re在线观看精品视频| 嫩草影院精品99| 久久久久国产一级毛片高清牌| 中文字幕精品免费在线观看视频| 欧美在线黄色| 亚洲五月天丁香| 国产日本99.免费观看| 脱女人内裤的视频| 午夜福利在线在线| 嫩草影院精品99| 国产成人啪精品午夜网站| 日本成人三级电影网站| 亚洲男人天堂网一区| 欧美人与性动交α欧美精品济南到| 日日爽夜夜爽网站| 神马国产精品三级电影在线观看 | 丁香欧美五月| 亚洲av熟女| 一级片免费观看大全| 黄色丝袜av网址大全| 不卡av一区二区三区| 亚洲精品久久成人aⅴ小说| 十八禁人妻一区二区| 亚洲欧美日韩高清在线视频| 久久婷婷成人综合色麻豆| 亚洲全国av大片| 久久精品国产清高在天天线| 91国产中文字幕| 在线国产一区二区在线| 波多野结衣av一区二区av| 69av精品久久久久久| 日本一本二区三区精品| 91字幕亚洲| 91av网站免费观看| 欧美一级毛片孕妇| 亚洲自偷自拍图片 自拍| 午夜影院日韩av| 日韩av在线大香蕉| 国产成人啪精品午夜网站| 999精品在线视频| 成人免费观看视频高清| 91麻豆av在线| 侵犯人妻中文字幕一二三四区| 亚洲精品美女久久久久99蜜臀| 午夜精品久久久久久毛片777| 国产精品野战在线观看| 男人舔奶头视频| 日日干狠狠操夜夜爽| 午夜成年电影在线免费观看| 老汉色∧v一级毛片| 淫秽高清视频在线观看| 欧美另类亚洲清纯唯美| 一级黄色大片毛片| 天堂√8在线中文| 97超级碰碰碰精品色视频在线观看| 欧美 亚洲 国产 日韩一| 免费av毛片视频| 99久久国产精品久久久| 美国免费a级毛片| 草草在线视频免费看| 男人的好看免费观看在线视频 | 天堂影院成人在线观看| 一区二区三区高清视频在线| 国产一区二区在线av高清观看| 色综合欧美亚洲国产小说| 又紧又爽又黄一区二区| 啦啦啦 在线观看视频| 欧美一级毛片孕妇| 亚洲中文字幕日韩| 亚洲av电影在线进入| 在线播放国产精品三级| 婷婷精品国产亚洲av| 老汉色av国产亚洲站长工具| 一级作爱视频免费观看| 波多野结衣巨乳人妻| 久久狼人影院| 在线观看午夜福利视频| 免费观看精品视频网站| 国产精品久久久久久亚洲av鲁大| 欧美zozozo另类| 中文字幕久久专区| 日本五十路高清| 麻豆av在线久日| 国产精品野战在线观看| 18禁黄网站禁片免费观看直播| 久久精品成人免费网站| 国产激情欧美一区二区| 无遮挡黄片免费观看| 18美女黄网站色大片免费观看| 最新美女视频免费是黄的| 亚洲av成人不卡在线观看播放网| 狠狠狠狠99中文字幕| 久久午夜综合久久蜜桃| 欧美+亚洲+日韩+国产| 久久中文字幕一级| 一本精品99久久精品77| www.自偷自拍.com| 69av精品久久久久久| 啦啦啦观看免费观看视频高清| 精品人妻1区二区| 性色av乱码一区二区三区2| 亚洲国产中文字幕在线视频| 亚洲熟女毛片儿| 黑丝袜美女国产一区| 亚洲av电影不卡..在线观看| 身体一侧抽搐| 午夜亚洲福利在线播放| 黑人欧美特级aaaaaa片| 午夜影院日韩av| 波多野结衣高清作品| 国产精品影院久久| 午夜亚洲福利在线播放| 国产伦一二天堂av在线观看| 视频在线观看一区二区三区| 久久中文看片网| 国产三级黄色录像| 国产精品永久免费网站| 亚洲成人久久性| 亚洲aⅴ乱码一区二区在线播放 | 免费观看精品视频网站| 国产一区在线观看成人免费| 欧美日韩福利视频一区二区| 欧美性长视频在线观看| 中文字幕人妻丝袜一区二区| www国产在线视频色| 成人特级黄色片久久久久久久| 淫秽高清视频在线观看| 成人亚洲精品av一区二区| 国产高清视频在线播放一区| 日本 欧美在线| 欧美日本视频| 亚洲成av片中文字幕在线观看| 久久精品国产99精品国产亚洲性色| 精品少妇一区二区三区视频日本电影| 成人亚洲精品一区在线观看| 久久午夜综合久久蜜桃| 久久国产乱子伦精品免费另类| 一本一本综合久久| 岛国在线观看网站| 国产三级黄色录像| 国产精品 欧美亚洲| 亚洲av成人不卡在线观看播放网| 欧美激情高清一区二区三区| 亚洲欧洲精品一区二区精品久久久| 大型黄色视频在线免费观看| 久久精品成人免费网站| 久久久久久久久中文| 国产免费av片在线观看野外av| aaaaa片日本免费| 精品少妇一区二区三区视频日本电影| 久久中文字幕人妻熟女| 久久久久久久久免费视频了| 美女午夜性视频免费| 操出白浆在线播放| 精品不卡国产一区二区三区| 天天添夜夜摸| 看免费av毛片| 在线永久观看黄色视频| 精品久久蜜臀av无| 欧美不卡视频在线免费观看 | 国产亚洲欧美在线一区二区| 国产精品久久电影中文字幕| 久久久久久九九精品二区国产 | 99riav亚洲国产免费| 久久久久免费精品人妻一区二区 | 国产亚洲精品av在线| 亚洲久久久国产精品| 国产久久久一区二区三区| 2021天堂中文幕一二区在线观 | 欧美亚洲日本最大视频资源| 搡老妇女老女人老熟妇| 男人操女人黄网站| 曰老女人黄片| 亚洲性夜色夜夜综合| 国产在线精品亚洲第一网站| 国产野战对白在线观看| 亚洲人成电影免费在线| 国产精品免费视频内射| 制服丝袜大香蕉在线| 看片在线看免费视频| 精品乱码久久久久久99久播| 久久伊人香网站| 国产又爽黄色视频| 国产日本99.免费观看| 亚洲中文字幕一区二区三区有码在线看 | 国产精品久久电影中文字幕| 最好的美女福利视频网| 欧美日韩一级在线毛片| 亚洲三区欧美一区| 人成视频在线观看免费观看| 国产成人精品久久二区二区免费| 亚洲国产欧洲综合997久久, | 国产免费av片在线观看野外av| 亚洲第一欧美日韩一区二区三区| 不卡一级毛片| 99riav亚洲国产免费| 在线十欧美十亚洲十日本专区| 黄色视频,在线免费观看| 国产成人av激情在线播放| 亚洲激情在线av| 欧美成人免费av一区二区三区| 国产野战对白在线观看| 人人妻人人澡人人看| 妹子高潮喷水视频| 在线观看午夜福利视频| 国产乱人伦免费视频| 亚洲av第一区精品v没综合| 男男h啪啪无遮挡| 男女床上黄色一级片免费看| 中文字幕最新亚洲高清| www日本黄色视频网| 久久久国产成人精品二区| www.自偷自拍.com| 久9热在线精品视频| 侵犯人妻中文字幕一二三四区| 大型av网站在线播放| 好看av亚洲va欧美ⅴa在| a级毛片a级免费在线| 免费观看精品视频网站| 深夜精品福利| 观看免费一级毛片| 亚洲激情在线av| 成人18禁高潮啪啪吃奶动态图| 免费电影在线观看免费观看| 正在播放国产对白刺激| 国产一区在线观看成人免费| 国产视频内射| 女人高潮潮喷娇喘18禁视频| 精品一区二区三区四区五区乱码| 一进一出抽搐gif免费好疼| 欧美乱妇无乱码| 亚洲欧美精品综合一区二区三区| 欧美成人性av电影在线观看| 好看av亚洲va欧美ⅴa在| 嫩草影院精品99| 亚洲午夜理论影院| 国产三级黄色录像| 黄色a级毛片大全视频| 十八禁网站免费在线| 国产又黄又爽又无遮挡在线| 久久99热这里只有精品18| 正在播放国产对白刺激| 日本 av在线| 首页视频小说图片口味搜索| 波多野结衣高清无吗| 美女 人体艺术 gogo| 韩国精品一区二区三区| 黄色片一级片一级黄色片| ponron亚洲| 免费在线观看亚洲国产| 久久久久亚洲av毛片大全| aaaaa片日本免费| 国产亚洲精品第一综合不卡| 国产区一区二久久| 韩国精品一区二区三区| 免费在线观看成人毛片| 神马国产精品三级电影在线观看 | 日韩一卡2卡3卡4卡2021年| 日韩精品青青久久久久久| 午夜老司机福利片| 亚洲av成人不卡在线观看播放网| 欧美人与性动交α欧美精品济南到| 精品国产超薄肉色丝袜足j| 给我免费播放毛片高清在线观看| 久久久久久免费高清国产稀缺| 制服诱惑二区| 一区二区三区高清视频在线| 99久久国产精品久久久| 男人舔奶头视频| 精品一区二区三区四区五区乱码| 国产精品影院久久| 日韩一卡2卡3卡4卡2021年| 欧美又色又爽又黄视频| 亚洲成人久久爱视频| 日本免费一区二区三区高清不卡| 国产国语露脸激情在线看| 男女午夜视频在线观看| 非洲黑人性xxxx精品又粗又长| 波多野结衣高清无吗| 亚洲免费av在线视频| 久久国产精品男人的天堂亚洲| 法律面前人人平等表现在哪些方面| 国产高清激情床上av| 亚洲天堂国产精品一区在线| 欧美黄色片欧美黄色片| 两个人看的免费小视频| 中文字幕人妻丝袜一区二区| 一级毛片女人18水好多| 久久精品影院6| 亚洲国产欧美网| 午夜亚洲福利在线播放| 婷婷亚洲欧美| 日韩高清综合在线| 精品一区二区三区视频在线观看免费| tocl精华| 日韩精品青青久久久久久| 欧美中文日本在线观看视频| 男女下面进入的视频免费午夜 | 亚洲av五月六月丁香网| 麻豆av在线久日| 国产午夜福利久久久久久| bbb黄色大片| 国产亚洲精品第一综合不卡| 狂野欧美激情性xxxx| 97人妻精品一区二区三区麻豆 | 一级黄色大片毛片| 午夜精品在线福利| 日韩精品青青久久久久久| 久久精品国产99精品国产亚洲性色| 日韩 欧美 亚洲 中文字幕| 国产片内射在线| 久久中文看片网| svipshipincom国产片| 啦啦啦免费观看视频1| 日本三级黄在线观看| 国产免费av片在线观看野外av| 欧美一级毛片孕妇| 日韩成人在线观看一区二区三区| 精品日产1卡2卡| 1024香蕉在线观看| 精品午夜福利视频在线观看一区| 男女下面进入的视频免费午夜 | 手机成人av网站| 男人的好看免费观看在线视频 | 丝袜在线中文字幕| 免费在线观看黄色视频的| 欧美av亚洲av综合av国产av| 高清在线国产一区| 男人的好看免费观看在线视频 |