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

    磁約束聚變等離子體的大規(guī)?;匦齽?dòng)理學(xué)模擬

    2016-06-17 03:23:10
    關(guān)鍵詞:核聚變等離子體

    肖 湧

    (浙江大學(xué) 聚變理論與模擬中心,浙江 杭州 310027)

    ?

    磁約束聚變等離子體的大規(guī)?;匦齽?dòng)理學(xué)模擬

    肖湧

    (浙江大學(xué) 聚變理論與模擬中心,浙江 杭州 310027)

    【摘要】文章介紹了磁約束核聚變能源科學(xué)及其密切相關(guān)的聚變等離子體物理,以及研究磁約束聚變等離子體非線性物理所必需的大規(guī)?;匦齽?dòng)理學(xué)模擬的方法及其前沿進(jìn)展.

    【關(guān)鍵詞】等離子體;核聚變;磁約束;大規(guī)模數(shù)值模擬;回旋動(dòng)理學(xué)

    1聚變能源科學(xué)

    1.1聚變能源

    核聚變能源是資源取之不盡、清潔的理想能源,將徹底解決人類社會(huì)能源問題.如圖1,太陽通過氫原子聚變生成氦原子釋放能量,在太陽核心部位,每秒鐘約有6.2億噸的氫原子發(fā)生聚變反應(yīng). 核聚變能源具有以下幾個(gè)特點(diǎn).

    1.1.1能量巨大

    核能是能量巨大的能源,一次核反應(yīng)釋放的能量是MeV量級,而平時(shí)我們廣泛使用的煤、石油、天然氣等化石能源涉及的是化學(xué)反應(yīng),釋放的能量是eV量級,只有前者的千分之一.核能有重核裂變和輕核聚變兩種方式,其中核聚變比核裂變釋放出更多的能量.例如,鈾235的裂變反應(yīng),將千分之一的物質(zhì)變成了能量;而氘的聚變反應(yīng),將近千分之四的物質(zhì)變成了能量.

    1.1.2資源豐富

    重核裂變使用的主要原料是鈾,目前探明的儲(chǔ)量僅夠使用約1 000年;而輕核聚變使用的燃料是海水中的氘,1升海水能提取30毫克氘,在聚變反應(yīng)中能產(chǎn)生約等于300升汽油的能量,即“1升海水約等于300升汽油”,地球上海水中就有45萬億噸氘,足夠人類使用百億年.聚變所需的另外一種燃料是氚,在地球上儲(chǔ)量并不多,但是可以通過鋰轉(zhuǎn)變獲得,而地球上尤其是月球上都有大量的鋰礦.一旦氘氚聚變實(shí)現(xiàn)商用,我們也可以采用其他聚變過程.總之,聚變能基本可以算是取之不盡的能源.

    1.1.3成本低廉

    1千克氘的價(jià)格只為1千克濃縮鈾的1/40.

    1.1.4安全、無污染

    核聚變不會(huì)產(chǎn)生長半衰期的放射性污染物,萬一發(fā)生事故,反應(yīng)堆會(huì)自動(dòng)冷卻而停止反應(yīng),不會(huì)發(fā)生爆炸.

    圖1 太陽的聚變反應(yīng)示意圖Figure 1 Fusion reaction schematic of the sun

    1.2磁約束聚變

    1.2.1約束方式

    聚變原料發(fā)生聚變反應(yīng)時(shí),最佳溫度高達(dá)一千萬度.在這種高溫下,電子被完全剝離出原子,形成物質(zhì)第四態(tài)等離子體.當(dāng)沒有約束時(shí),帶電粒子的熱運(yùn)動(dòng)會(huì)使等離子體迅速擴(kuò)散,無法發(fā)生核聚變反應(yīng),所以需要外界約束使高溫等離子體聚在一起,保持一定的密度,從而發(fā)生可觀的核聚變反應(yīng).等離子體目前有兩種約束方式,磁場約束和慣性約束.

    1) 磁約束

    磁場可以約束粒子,使其沿磁力線運(yùn)動(dòng),如圖2.等離子體作為帶電粒子,會(huì)受到電磁場的影響,所以可以通過外加磁場約束等離子體.磁場是目前已知最好的絕緣體,其隔熱性能比航天飛機(jī)的瓦片還要強(qiáng)1 000倍.

    圖2 磁場對粒子的約束Figure 2 Constraint of magnetic field on the particles

    2) 磁約束中存在的問題

    磁場只能約束粒子在垂直磁場方向的運(yùn)動(dòng),但并不能約束粒子在平行磁場方向的運(yùn)動(dòng),如圖3,帶電粒子在垂直磁場的方向做回旋運(yùn)動(dòng),磁場不均勻性會(huì)產(chǎn)生梯度漂移.為此可以將磁場彎曲成環(huán)形,首尾相接,如圖4,這樣就可以把粒子約束在環(huán)中運(yùn)動(dòng).

    圖3 帶電粒子在磁場中的運(yùn)動(dòng)圖Figure 3 Movement of charged particles in magnetic field

    圖4 環(huán)形磁場中螺旋形狀磁場線的方向Figure 4 Direction of magnetic field lines with circular spiral shape in magnetic field

    但是,當(dāng)把磁場彎曲成環(huán)形時(shí),磁場中存在垂直于磁場方向的梯度,粒子會(huì)發(fā)生垂直于磁場方向的梯度漂移.同時(shí)粒子沿著環(huán)形磁力線上的運(yùn)動(dòng),會(huì)受到沿磁場曲率方向的離心力作用,從而產(chǎn)生一個(gè)垂直于離心力和外加磁場為方向的曲率漂移.因此需要再增加一個(gè)沿極向的磁場,使得粒子在不同位置處的漂移作用相互抵消,從而保證約束的穩(wěn)定.

    目前國際上比較主流的磁約束聚變裝置有托卡馬克和仿星器等,如圖5(a).其中托卡馬克由于其磁場設(shè)計(jì)的相對簡單和約束性能的良好成為各國進(jìn)行磁約束聚變研究的首選.仿星器有在穩(wěn)態(tài)運(yùn)行方面有其獨(dú)特的優(yōu)點(diǎn),但是其復(fù)雜磁場的設(shè)計(jì)和建造也極有挑戰(zhàn)性.

    3) 慣性約束

    除了磁約束的方法之外,還可以使用慣性來約束等離子體.宇宙中的許多星體都是依靠自身巨大質(zhì)量產(chǎn)生的引力作用來約束粒子的,如圖5(b)所示,上面為太陽依靠其巨大的慣性引力約束等離子體,下面為星云依靠其巨大的慣性引力約束等離子體.在實(shí)驗(yàn)室里可以通過強(qiáng)激光壓縮氘氚靶丸實(shí)現(xiàn)慣性約束,如美國的國家點(diǎn)火裝置(NIF)和我國的神光系列裝置都是采用這一原理的.

    圖5 等離子體的約束Figure 5 Plasma confinement

    1.2.2磁約束裝置

    1) 托卡馬克(Tokamak)磁約束聚變實(shí)驗(yàn)裝置

    托卡馬克是前蘇聯(lián)科學(xué)家發(fā)明的一種利用磁約束來實(shí)現(xiàn)磁約束聚變的環(huán)性容器,依靠環(huán)形線圈和等離子體電流產(chǎn)生的環(huán)狀螺旋形強(qiáng)磁場,將極高溫等離子體狀態(tài)的聚變物質(zhì)約束在環(huán)形容器里,以此實(shí)現(xiàn)聚變反應(yīng).

    2) 國際磁約束核聚變實(shí)驗(yàn)裝置

    自從20世紀(jì)50年代以來,人類開始研究聚變能源及其密切相關(guān)的等離子體物理,測試了很多不同的實(shí)驗(yàn)裝置,發(fā)現(xiàn)到目前為止托卡馬克還是表現(xiàn)最好的.評判聚變裝置約束效果的重要參數(shù)包括聚變等離子體的溫度、密度以及約束時(shí)間,理論分析表明只有當(dāng)這三者的乘積超過某一臨界值時(shí),燃料才能燃燒發(fā)生聚變,這就是著名的勞遜判據(jù).自20世紀(jì)90年代以來,歐盟的JET、美國的TFTR、DIII-D和日本的JT-60、JT-60SA這幾個(gè)大型托卡馬克裝置在磁約束核聚變研究中獲得許多重要成果.圖6為不同托卡馬克試驗(yàn)裝置大半徑對比圖.我們現(xiàn)在的實(shí)驗(yàn)裝置越來越靠近勞遜判據(jù)的臨界值,如圖7和圖8.另外,一個(gè)經(jīng)驗(yàn)判定是,托卡馬克尺寸越大越穩(wěn)定,約束時(shí)間越長,越接近勞遜判據(jù)的臨界值,這就是為什么現(xiàn)在托卡馬克越做越大的原因.

    圖6 不同托卡馬克試驗(yàn)裝置大半徑對比Figure 6 Major radius comparison among different Tokamak experiment devices

    圖7 不同托卡馬克實(shí)驗(yàn)裝置約束時(shí)間對比Figure 7 Confinement time comparison among different Tokamak experiment devices

    圖8 不同托卡馬克實(shí)驗(yàn)裝置約束效果對比Figure 8 Confinement effect comparison among different Tokamak experiment devices

    隨著研究進(jìn)展的深入,各國之間的學(xué)術(shù)交流愈加頻繁,國際熱核聚變實(shí)驗(yàn)反應(yīng)堆(International Thermonuclear Experimental Reactor,ITER)將成為世界上最大的磁約束等離子體物理實(shí)驗(yàn)裝置,如圖9. 這是目前正在建設(shè)的世界上最大的實(shí)驗(yàn)性托卡馬克核聚變反應(yīng)堆,由七個(gè)成員實(shí)體資助和運(yùn)行,包括歐盟、印度、日本、中國、俄羅斯、韓國和美國.它建立在TFTR、歐洲聯(lián)合環(huán)形加速器(JET)、JT-60和T-15等裝置的研究之上,并將顯著的超越所有前者.ITER實(shí)驗(yàn)成功后將是人類歷史上第一個(gè)長時(shí)間實(shí)現(xiàn)核聚變反應(yīng)的實(shí)驗(yàn)裝置,對整個(gè)聚變能源科學(xué)的影響是巨大的.

    國際熱核反應(yīng)試驗(yàn)堆計(jì)劃(ITER)是商用化核聚變能源的關(guān)鍵一步.作為ITER成員,中國在國家中長期科技規(guī)劃中把ITER計(jì)劃作為重點(diǎn)高技術(shù),投入近百億人民幣進(jìn)行研究.我國關(guān)于ITER采購包的研制在我國科研和工程技術(shù)人員的協(xié)同努力下進(jìn)展順利,在七方成員中排在前列,體現(xiàn)了我國的科學(xué)技術(shù)實(shí)力和綜合國力.

    圖9 國際熱核反應(yīng)試驗(yàn)堆(ITER)結(jié)構(gòu)示意圖Figure 9 Structure of ITER

    ITER計(jì)劃中的等離子體的參數(shù):

    ·等離子體大半徑:6.2 m

    ·等離子體小半徑:2.0 m

    ·等離子體容量:840 m3

    ·等離子體電流:15.0 MA

    ·軸向環(huán)形磁場強(qiáng)度:5.3T

    ·聚變功率:500 MW

    ·等離子體維持時(shí)間:>400s

    ·能量倍增因數(shù)(Q值):>10

    (3) 國內(nèi)磁約束核聚變實(shí)驗(yàn)裝置

    目前,國內(nèi)主要的托卡馬克實(shí)驗(yàn)裝置有中科院等離子體所的EAST(圖10)和核工業(yè)部西南物理研究院的HL-2A(圖11).

    圖10 中科院等離子體所EAST超導(dǎo)托卡馬克Figure 10 EAST superconducting tokamak in Institute of Plasma Physics, Chinese Academy of Sciences

    圖11 核工業(yè)部西南物理研究院HL-2A托卡馬克Figure 11 HL-2A tokamak in Southwestern Institute of Physics, Department of Nuclear Industry

    先進(jìn)實(shí)驗(yàn)超導(dǎo)托卡馬克實(shí)驗(yàn)裝置(Experimental Advanced Superconducting Tokamak,EAST)是一個(gè)具有非圓小截面的全超導(dǎo)托卡馬克.它可以實(shí)現(xiàn)高參數(shù)和長脈沖的實(shí)驗(yàn)運(yùn)行;科研人員在以上條件下可以開發(fā)先進(jìn)運(yùn)行模式,并進(jìn)行熱流平衡和粒子流平衡控制的實(shí)驗(yàn)研究.

    環(huán)流器二號A(HL-2A)是我國第一個(gè)具有先進(jìn)偏濾器位形的非圓截面的托卡馬克核聚變實(shí)驗(yàn)研究裝置,可以進(jìn)行高參數(shù)等離子體條件下的改善約束實(shí)驗(yàn),并利用其獨(dú)特的大體積封閉偏濾器結(jié)構(gòu),開展核聚變領(lǐng)域許多前沿物理課題以及相關(guān)工程技術(shù)的研究.

    1.3磁約束聚變的目標(biāo)和難點(diǎn)

    磁約束聚變的目標(biāo)是使粒子在擴(kuò)散之前發(fā)生聚變,并讓聚變產(chǎn)生的alpha粒子加熱燃料,從而實(shí)現(xiàn)自持燃燒.實(shí)現(xiàn)這些目標(biāo)需要對磁約束聚變發(fā)生的等離子體物理有深刻的理解,主要涉及等離子體的約束、輸運(yùn)、加熱和電流控制等復(fù)雜物理過程.同時(shí)這些聚變等離子體物理需要研究波-粒子、波-波以及粒子-粒子等多時(shí)空尺度下復(fù)雜的非線性行為,需要設(shè)計(jì)復(fù)雜的計(jì)算機(jī)模擬來輔助裝置設(shè)計(jì),物理實(shí)驗(yàn)設(shè)計(jì)以及實(shí)驗(yàn)解釋和預(yù)測.從20世紀(jì)90年代開始,磁約束聚變等離子體模擬就開始使用日益興起的超級計(jì)算技術(shù)進(jìn)行大規(guī)模并行模擬,開發(fā)了針對聚變等離子體宏觀微觀不穩(wěn)定性、湍流輸運(yùn)和加熱物理研究的大規(guī)模數(shù)值模擬程序,并興起了以回旋動(dòng)理學(xué)模擬為主導(dǎo)的模擬科學(xué).

    2大規(guī)模并行數(shù)值模擬

    2.1大規(guī)模并行數(shù)值模擬的起源

    大規(guī)模并行數(shù)值模擬是指將一個(gè)復(fù)雜的計(jì)算問題分成若干部分,通過以網(wǎng)絡(luò)(光纖)相互聯(lián)接的大量的計(jì)算機(jī)核(幾百到20萬核)的超級計(jì)算機(jī),同時(shí)計(jì)算同一個(gè)任務(wù)的不同部分.其目的是提高問題求解速度,同時(shí)擴(kuò)大求解問題規(guī)模.

    大規(guī)模并行數(shù)值模結(jié)合了最尖端的計(jì)算機(jī)科學(xué)、應(yīng)用數(shù)學(xué)和具體科技應(yīng)用領(lǐng)域的跨學(xué)科全新研究方法,適應(yīng)于那些對科學(xué)研究和技術(shù)創(chuàng)新有重大革新性意義的課題,適合于團(tuán)隊(duì)作業(yè).

    2.2大規(guī)模并行數(shù)值模擬的發(fā)展

    數(shù)值模擬已經(jīng)成為除理論、實(shí)驗(yàn)外,支持科學(xué)探索和技術(shù)革新的第三個(gè)主要研究手段.大規(guī)模數(shù)值模擬能夠幫助理解和預(yù)測那些科學(xué)技術(shù)上極其重要的系統(tǒng)的行為,從而加速科技創(chuàng)新的速度.超級計(jì)算能力已成為衡量一個(gè)國家(地區(qū))科技實(shí)力和創(chuàng)新能力的重要標(biāo)志.

    2.3超級計(jì)算機(jī)的世界分布

    大規(guī)模數(shù)值模擬從20世紀(jì)90年代中葉開始興起,最早的超級計(jì)算機(jī)計(jì)算能力約為1G flops,其計(jì)算能力以每年2~3倍的速度增加.現(xiàn)在世界上最強(qiáng)的超級算計(jì)機(jī)的速度在1015flop量級,正在向1018flop進(jìn)軍.目前世界上最強(qiáng)的超級算計(jì)主要集中在美國、中國、歐洲和日本.其中美國的超級計(jì)算資源占到全世界一半左右.

    1) 國際超級計(jì)算機(jī)

    近幾年,世界上排名前500的超級計(jì)算機(jī)中,較著名的有美國橡樹嶺國家實(shí)驗(yàn)室的“泰坦”,美國勞倫斯利福摩爾國家實(shí)驗(yàn)室的“紅杉”,日本理化學(xué)研究所的“京”,美國阿貢國家實(shí)驗(yàn)室的“米拉”等.

    2) 國內(nèi)超級計(jì)算機(jī)

    中國硬件發(fā)展上已經(jīng)成為僅次于美國的超級計(jì)算大國,深圳星云超級計(jì)算機(jī),曾排名世界第二;天津天河超級計(jì)算機(jī),2011排名世界第一;廣州天河-2超級計(jì)算機(jī),2014排名世界第一.

    現(xiàn)階段我國超級計(jì)算機(jī)發(fā)展中的主要問題是軟件發(fā)展和硬件發(fā)展不相稱.超級計(jì)算平臺(tái)輔助軟件和大規(guī)模數(shù)值模擬應(yīng)用程序的發(fā)展嚴(yán)重滯后于超級計(jì)算機(jī)的建造.

    2.4大規(guī)模數(shù)值模擬的應(yīng)用

    大規(guī)模數(shù)值模擬廣泛應(yīng)用于新能源、新材料、物理、化工、生物、氣候環(huán)境、物流、金融等領(lǐng)域,對國民經(jīng)濟(jì)和科技創(chuàng)新產(chǎn)生廣泛而深刻的影響,如圖12.

    圖12 大規(guī)模數(shù)值模擬應(yīng)用舉例Figure 12 Application example of large-scale numerical simulation

    3回旋動(dòng)理學(xué)模擬

    基于磁約束聚變等離子物理的高度復(fù)雜性和多時(shí)空尺度特性,回旋動(dòng)理學(xué)模擬已成為研究聚變中不穩(wěn)定性和湍流輸運(yùn)的主要工具.近幾十年來的研究表明,托卡馬克中普遍存在的電磁場微擾是造成托卡馬克反常輸運(yùn)的主要原因,一般也稱之為湍流輸運(yùn).這些電磁場微擾是由各種宏觀和微觀的電磁不穩(wěn)定性激發(fā)的.聚變等離子體約束的好壞直接決定于其輸運(yùn)過程,其產(chǎn)生和控制的物理機(jī)制是目前研究的前沿課題.

    3.1等離子體的輸運(yùn)過程

    3.1.1輸運(yùn)過程的種類

    在等離子體中,輸運(yùn)過程可以分為碰撞輸運(yùn)與湍流輸運(yùn).

    1) 碰撞輸運(yùn)

    等離子體中粒子之間會(huì)存在相互作用,一般稱為碰撞,可以用碰撞算符來表示.碰撞可以使粒子回旋運(yùn)動(dòng)的導(dǎo)心或者沿著磁力線運(yùn)動(dòng)的香蕉軌道中心產(chǎn)生偏移,從而導(dǎo)致碰撞輸運(yùn),如圖13.碰撞輸運(yùn)即使在沒有電磁場擾動(dòng)的情況下仍然存在,是托卡馬克裝置最低水平的輸運(yùn)過程.

    圖13 帶有相反電荷的粒子由于導(dǎo)心漂移發(fā)生碰撞Figure 13 collision of particles with opposite charge due to guide heart drift

    2) 湍流輸運(yùn)

    湍流是托卡馬克中粒子和熱輸運(yùn)的主要原因.湍流輸運(yùn)可以分為兩種,靜電性和電磁性.根據(jù)其不穩(wěn)定性來源又可細(xì)分為各種不同模,比如離子梯度模(ITG)、電子梯度模(ETG)、捕獲電子模(TEM)、動(dòng)理學(xué)氣球模(KBM),等等.每個(gè)模線性不穩(wěn)定性原因不一樣,同時(shí)飽和機(jī)制和輸運(yùn)機(jī)制也往往不一用,從而導(dǎo)致輸運(yùn)系數(shù)不一樣.一種最常見湍流輸運(yùn)原因就是帶電粒子的在靜電場微擾下的EXB漂移,這種垂直于磁場的漂移由于電場的退相關(guān)過程而導(dǎo)致隨機(jī)運(yùn)動(dòng),從而導(dǎo)致徑向的湍流輸運(yùn),如圖14.

    圖14 導(dǎo)心的EXB漂移導(dǎo)致湍流輸運(yùn)Figure 14 Turbulent transport because of EXB drift of guide heart

    3.2回旋動(dòng)理學(xué)簡介

    3.2.1Vlasov方程與Maxwell方程

    等離子體有動(dòng)理學(xué)和流體力學(xué)兩種描述方式,其中流體力學(xué)描述可通過對動(dòng)理學(xué)方程求矩獲得,但是它忽略了朗道阻尼等重要?jiǎng)永韺W(xué)內(nèi)容,因而不能準(zhǔn)確描述聚變等離子體中一些重要的動(dòng)理學(xué)不穩(wěn)定性及其湍流輸運(yùn).最一般的動(dòng)理學(xué)方程可以從非平衡態(tài)統(tǒng)計(jì)力學(xué)中的BBGKY遞推關(guān)系獲得如下的帶電磁場的玻爾茲曼方程

    (1)

    -▽2φ=4π(qi∫d3vfi-e∫d3vfe).

    (2)

    其中,f表示單粒子的相空間分布函數(shù).第一個(gè)方程在等離子體物理被稱之為含碰撞的Vlasov方程,方程右邊為代表兩粒子相關(guān)函數(shù)的碰撞算符.為簡單起見,我們用第二個(gè)方程即泊松方程代替求解電磁場的Maxwell(麥克斯韋)方程組.這組方程可以視為靜電極限下的Vlasov-Maxwell方程組.上面這組方程可以描述時(shí)間尺度從等離子體震蕩到輸運(yùn),空間尺度從德拜長度到系統(tǒng)尺寸等大范圍多尺度物理.依靠現(xiàn)在的超級計(jì)算機(jī)的計(jì)算能力是不可能模擬這么大范圍的多尺度物理.

    但在托卡馬克中引起湍流輸運(yùn)最危險(xiǎn)的模是低頻漂移波湍流,比如ITG、TEM、KBM和ETG等,在燃燒等離子體中可能會(huì)激發(fā)比較高頻的Alfven波,但是仍然低于離子回旋頻率,如圖15.

    圖15 托卡馬克中不穩(wěn)定性的時(shí)空頻率特點(diǎn)Figure 15 Space and time frequency characteristics of instability in Tokamak

    (3)

    在這種條件下,我們可以將式(3)的Vlasov方程進(jìn)行坐標(biāo)變換:

    (x,v)→(R,μ,v||),

    (4)

    將方程從粒子坐標(biāo)移到導(dǎo)心坐標(biāo),然后再回旋平均得到如下的回旋動(dòng)理學(xué)(gyrokinetics)-麥克斯韋方程組[1-6]:

    =C[f].

    (5)

    相應(yīng)的回旋動(dòng)理學(xué)泊松方程整理如下:

    (6)

    其中

    代表從導(dǎo)心到粒子位置的反演[7].強(qiáng)磁場下的導(dǎo)心運(yùn)動(dòng)如圖16.

    圖16 強(qiáng)磁場下的導(dǎo)心運(yùn)動(dòng)Figure 16 Movement of guide heart in strong magnetic field

    3.3回旋動(dòng)理學(xué)模擬

    3.3.1模擬方法和算法特點(diǎn)

    回旋動(dòng)理學(xué)模擬是研究磁約束等離子體微觀不穩(wěn)定性中非線性物理的重要工具.它將六維相空間的動(dòng)理學(xué)方程變成五維相空間,大大提高計(jì)算效率.相對流體模擬而言,它又保留了磁約束聚變等離子體物理中如朗道阻尼重要的動(dòng)理學(xué)效應(yīng),可以比較真實(shí)地描述等離子體物理.等離子體物理模擬一般通過將模擬結(jié)果與理論分析進(jìn)行對比來驗(yàn)證程序真實(shí)反映了理論模型,再用模擬與實(shí)驗(yàn)結(jié)果進(jìn)行比較來驗(yàn)證物理模型的可信度,最終達(dá)到模擬的可預(yù)測性.

    回旋動(dòng)理學(xué)湍流模擬的理論基礎(chǔ)起源于20世紀(jì)70年代,從90年代以來隨著超級計(jì)算機(jī)和大規(guī)模數(shù)值模擬的發(fā)展而獲得蓬勃發(fā)展.回旋動(dòng)理學(xué)-麥克斯韋方程體系有兩種主要數(shù)值方法可求解其動(dòng)力學(xué)行為,一種叫粒子模擬方法,即用蒙特卡洛的方法對相空間取樣;另外一種是用解多維偏微分方程的網(wǎng)格法,也叫Vlasov方法或者連續(xù)法.在算法上,為了去除模擬中的粒子數(shù)或網(wǎng)格數(shù)與物理要求相差很遠(yuǎn)而導(dǎo)致的噪聲,我們一般采用所謂部分分布函數(shù)法,關(guān)注分布函數(shù)中偏離平衡分布的部分,從而提高信噪比.

    3.3.2模擬程序

    從20世紀(jì)90年代以來,涌現(xiàn)了很多回旋動(dòng)理學(xué)程序.這些程序的發(fā)展代表了磁約束聚變模擬這二十年來最富代表性的發(fā)展歷程.代表性的程序美國有GTC、GYRO、GS2、GEM、XGC等,歐洲有GENE、ORB5等,日本有GT3D等.

    其中,GTC(Gyrokinetic Toroidal Code)程序[8]是美國加州大學(xué)的林志宏教授90年代從普林斯頓讀研究生時(shí)開始發(fā)展的三維全域性粒子回旋動(dòng)理學(xué)程序.歷經(jīng)二十年的發(fā)展,已成長為一個(gè)包含了微觀湍流、動(dòng)理學(xué)磁流體、射頻加熱等多功能的模擬程序或者模擬平臺(tái).近年來浙江大學(xué)、中科院物理所和北京大學(xué)也積極參與了其程序開發(fā).

    GTC程序的創(chuàng)新之一體現(xiàn)在其采用磁坐標(biāo),構(gòu)造了沿著磁力線的三維電磁場網(wǎng)格,如圖17,圖18和圖19.在某一個(gè)極面上GTC采用非結(jié)構(gòu)化網(wǎng)格,保證徑向和極向計(jì)算精度大致相近;在同一個(gè)磁面上,格點(diǎn)沿著磁力線成為一條直線.沿著磁力線的格點(diǎn)設(shè)置契合磁約束聚變等離子體各項(xiàng)異性的特點(diǎn).回旋平均和導(dǎo)心到粒子位置的反演都采用了所謂4點(diǎn)平均的辦法[9].這些計(jì)算算法的創(chuàng)新使得GTC能夠進(jìn)行真正三維的實(shí)空間的全域模擬,同時(shí)大大提高了計(jì)算效率.

    圖17 GTC磁面上的格點(diǎn)設(shè)置Figure 17 Lattice point set on magnetic surface of GTC

    圖18 GTC極面上的格點(diǎn)設(shè)置Figure 18 Lattice point set on polar surface of GTC

    圖19 GTC的3維格點(diǎn)設(shè)置和湍流特點(diǎn)Figure 19 3D point set and turbulence characteristics of GTC

    3.3.3粒子推動(dòng)方程

    回旋動(dòng)理學(xué)粒子模擬方法主要由兩部分組成,導(dǎo)心(或者回旋中心)運(yùn)動(dòng)積分器和電磁場解法器.導(dǎo)心的運(yùn)動(dòng)方程可以描述如下:

    (7)

    (8)

    假設(shè)磁場擾動(dòng)可以寫成如下形式:

    δB=▽×αB,

    (9)

    ρc=ρ||+α,

    (10)

    D=ρc(gI′-Ig′)+I+qg,

    (11)

    磁坐標(biāo)下的運(yùn)動(dòng)方程可以表述為[10]:

    (12)

    (13)

    (14)

    (15)

    式(12)~(15)這組導(dǎo)心運(yùn)動(dòng)方程可以用來描述導(dǎo)心在托卡馬克中運(yùn)動(dòng),一般用2階龍格庫塔方法來積分導(dǎo)心軌道,因?yàn)樗且环N顯式算法,可以和描述碰撞的蒙特卡洛方法相匹配.

    3.3.4泊松方程

    (16)

    式(16)的關(guān)鍵在于計(jì)算分子中的拉普拉斯算符(Laplacian).

    3.3.5一般幾何位型下的拉普拉斯算符

    磁坐標(biāo)下的拉普拉斯算符(Laplacian)展開可得

    (17)

    其中,協(xié)變張量gxy可以通過B樣條方法計(jì)算得出.逆變張量gxy可以通過協(xié)變張量gxy計(jì)算得出.

    3.4基于回旋動(dòng)理學(xué)模擬的發(fā)現(xiàn)

    回旋動(dòng)理學(xué)模擬能夠幫助解釋實(shí)驗(yàn)現(xiàn)象和促進(jìn)理論的發(fā)展,并為實(shí)驗(yàn)建模提供必需的參數(shù).

    3.4.1帶狀流對湍流的抑制

    回旋動(dòng)理學(xué)模擬最有影響的結(jié)果之一是模擬發(fā)現(xiàn)湍流會(huì)自發(fā)激發(fā)帶狀流(zonal flow),從而反過來抑制湍流[7],如圖20.這一模擬上的發(fā)現(xiàn)極大促進(jìn)人們對于湍流理論和實(shí)驗(yàn)的理解,比如促進(jìn)了漂移波參量不穩(wěn)定性理論的發(fā)展[11].我國的實(shí)驗(yàn)也多次發(fā)現(xiàn)了帶狀流對湍流的調(diào)制作用[12-13].

    3.4.2湍流輸運(yùn)系數(shù)

    關(guān)于磁約束聚變等離子體輸運(yùn)過程的集成建模迫切需要在給定等離子體徑向剖面的條件下湍流輸運(yùn)的強(qiáng)度,基于第一性原理的回旋動(dòng)理學(xué)模擬能夠提供湍流輸運(yùn)的系數(shù)[15],如圖21.而對于湍流強(qiáng)度的徑向分布目前很難用現(xiàn)有理論預(yù)測.

    圖20 GTC模擬發(fā)現(xiàn)捕獲電子模湍流中帶狀流對湍流有抑制作用[14]Figure 20 GTC simulation found that capture banded in electronic mode turbulence flow has inhibitory effect on turbulence

    圖21 回旋動(dòng)理學(xué)模擬提供的熱傳導(dǎo)系數(shù)的徑向分布Figure 21 Radial distribution of heat transport coefficient by gyrokinetic simulation

    4結(jié)語

    本文介紹了磁約束核聚變中的基于回旋動(dòng)理學(xué)的大規(guī)模數(shù)值模擬.最近20年來這個(gè)領(lǐng)域已經(jīng)取得很大進(jìn)步,但是還有很多關(guān)鍵性問題沒有得到解決,比如聚變等離子體邊界電磁模擬和長時(shí)間尺度模擬,仍需要研究人員的持續(xù)努力.同時(shí),這個(gè)領(lǐng)域隨著超級計(jì)算技術(shù)的進(jìn)步,將會(huì)進(jìn)一步促進(jìn)可控磁約束聚變的實(shí)現(xiàn).

    (致謝:本文作者感謝袁泉準(zhǔn)備部分初稿并仔細(xì)校對文稿,同時(shí)感謝P.J. Catto和Z. Lin博士提供部分插圖.)

    【參考文獻(xiàn)】

    [1]CATTO P J. Linear gyrokinetics[J]. Plasma Physics,1978,20:719.

    [2]ANTONSEN T, LANE B. Kinetic equations for low frequency instabilites in inhomogeneous plasmas[J].Physics of Fluids,1980,23(6):061205.

    [3]FRIEMAN E A, CHEN L. Nonlinear gyrokinetic equations for low-frequency electromagnetic waves in general plasma equilibria[J]. Physics of Fluidsvol,1982,25(3):502-508.

    [4]HAHM T S. Nonlinear gyrokinetic equations for tokamak microturbulence[J]. Physics of Fluids,1988,31(9):2670-2673.

    [5]QIN H W, TANG M, REWOLDT G. Gyrokinetic theory for arbitrary wavelength electromagnetic modes in tokamaks[J]. Physics of Plasmas,1998,5(4):1035-1049.

    [6]LEE W W. Nonlinear gyrokinetic equations[J]. Physics of Fluids,1983,26(12):3524.

    [7]LIN Zhihong, HAHM T S, LEE W W, et al.Turbulent transport reduction by zonal flows: Massively parallel simulations[J]. Science,1998,281(5384):1835-1837.

    [8]LEE W W. Gyrokinetic particle simulation model[J]. Journal of Computational Physics,1987,72(1):243-269.

    [9]LIN Z, LEE W W. Method for solving the gyrokinetic Poisson equation in general geometry[J]. Phys Rev E,1995,52(5):5646-5652.

    [10]WHITE R B, CHANCE M S. Hamiltonian guiding center drift orbit calculation for plasmas of arbitrary cross-section[J]. Physics of Fluids,1984,27(10):2455-2467.

    [11]CHEN L, LIN Z, WHITE R B. Excitation of zonal flow by drift waves in toroidal plasmas[J].Physics of Plasmas,2000,7(8):3129.

    [12]ZHAO K J, LAN T, DONG J Q, et al. Toroidal symmetry of the geodesic acoustic mode zonal flow in a tokamak plasma[J]. Physical Review Letters,2006,96(25):255004.

    [13]XU G, WAN B N, SONG M, et al. Direct measurement of poloidal long-wavelength E× B flows in the HT-7 tokamak[J]. Physical Review Letters,2003,91(12):125001.

    [14]XIAO Y, LIN Z. Turbulent transport of trapped-electron modes in collisionless plasmas[J]. Physical Review Letters, 2009,103(8):085004.

    [15]XIAO Y, HOLOD I, WANG Z X, et al. Gyrokinetic particle simulation of microturbulence for general magnetic geometry and experimental profiles[J]. Physics of Plasmas,2015(22):022516.

    Massively parallel gyrokinetic simulation for magnetic fusion plasmas

    XIAO Yong

    (Institute for Fusion Theory and Simulation, Zhejiang University, Hangzhou 310027, China)

    Abstract:This article introduces the magnetic fusion energy science and its associated fusion plasma physics. Massively parallel gyrokinetic simulation has emerged as a major tool to study the nonlinear physics. Its research frontier is also discussed.

    Key words:plasma; fusion; magnetic confinement; large scale data simulation; gyrokinetics

    【文章編號】1004-1540(2016)01-0001-10

    DOI:10.3969/j.issn.1004-1540.2016.01.001

    【收稿日期】2016-01-18《中國計(jì)量學(xué)院學(xué)報(bào)》網(wǎng)址:zgjl.cbpt.cnki.net

    【作者簡介】肖湧(1977- ),男,湘南省湘潭人,教授,主要研究方向?yàn)榇偶s束核聚變等離子體中的不穩(wěn)定性,輸運(yùn)和加熱以及大規(guī)模并行回旋動(dòng)理學(xué)模擬等. E-mail:yxiao@zju.edu.cn

    【中圖分類號】O532+.11

    【文獻(xiàn)標(biāo)志碼】A

    猜你喜歡
    核聚變等離子體
    人造太陽,追夢無限能源
    趨勢
    核聚變能應(yīng)用前景現(xiàn)一絲曙光
    連續(xù)磁活動(dòng)對等離子體層演化的影響
    基于低溫等離子體修飾的PET/PVC浮選分離
    氣道內(nèi)氬等離子體凝固術(shù)治療應(yīng)用進(jìn)展
    等離子體對燃?xì)庠谘a(bǔ)燃室中燃燒特性的影響
    等離子體種子處理技術(shù)介紹
    百科探秘·航空航天(2016年9期)2016-12-01 03:04:17
    德國仿星器:聚變裝置造出氫等離子體
    新天地(2016年3期)2016-05-30 10:48:04
    国产精品久久视频播放| 亚洲av第一区精品v没综合| 国产高清激情床上av| 免费看a级黄色片| 亚洲 国产 在线| 亚洲成人久久性| 每晚都被弄得嗷嗷叫到高潮| 国产真人三级小视频在线观看| 中文字幕人妻丝袜一区二区| 18禁裸乳无遮挡免费网站照片 | 可以在线观看毛片的网站| 在线十欧美十亚洲十日本专区| 999久久久精品免费观看国产| 在线观看一区二区三区| 90打野战视频偷拍视频| 久99久视频精品免费| 国产精品 欧美亚洲| 国产亚洲精品一区二区www| 亚洲精品在线美女| 18禁黄网站禁片免费观看直播| 哪里可以看免费的av片| 精品欧美国产一区二区三| 一卡2卡三卡四卡精品乱码亚洲| 一级毛片高清免费大全| 亚洲国产欧美日韩在线播放| 在线免费观看的www视频| 日本a在线网址| 国产极品粉嫩免费观看在线| 脱女人内裤的视频| 妹子高潮喷水视频| 1024香蕉在线观看| 啪啪无遮挡十八禁网站| 欧美色欧美亚洲另类二区| 欧美久久黑人一区二区| 国产精品98久久久久久宅男小说| 一级a爱片免费观看的视频| 国产午夜精品久久久久久| 视频在线观看一区二区三区| 色av中文字幕| 激情在线观看视频在线高清| 巨乳人妻的诱惑在线观看| 亚洲无线在线观看| 99精品久久久久人妻精品| 日本在线视频免费播放| 婷婷亚洲欧美| 18禁黄网站禁片免费观看直播| av有码第一页| 精品熟女少妇八av免费久了| 亚洲 欧美一区二区三区| 婷婷精品国产亚洲av| 伊人久久大香线蕉亚洲五| 国产欧美日韩一区二区三| 久久精品国产亚洲av高清一级| 免费在线观看成人毛片| 国产精品野战在线观看| 悠悠久久av| 亚洲 国产 在线| 伊人久久大香线蕉亚洲五| 亚洲 国产 在线| 99国产极品粉嫩在线观看| 此物有八面人人有两片| 精品国产一区二区三区四区第35| 88av欧美| 国产精品自产拍在线观看55亚洲| 国产精品av久久久久免费| 亚洲人成77777在线视频| videosex国产| 国产成人精品久久二区二区91| 欧美成人性av电影在线观看| 免费看十八禁软件| 久久这里只有精品19| 无人区码免费观看不卡| 久久精品国产99精品国产亚洲性色| 久久国产亚洲av麻豆专区| 熟女电影av网| 女同久久另类99精品国产91| 久久狼人影院| 香蕉丝袜av| 1024手机看黄色片| 日韩一卡2卡3卡4卡2021年| 国产日本99.免费观看| 国产高清有码在线观看视频 | 女性被躁到高潮视频| 一级黄色大片毛片| 亚洲成人久久性| 黄色成人免费大全| 熟妇人妻久久中文字幕3abv| 亚洲中文日韩欧美视频| 99热6这里只有精品| 黄片大片在线免费观看| 日韩欧美在线二视频| 日日摸夜夜添夜夜添小说| 成人av一区二区三区在线看| 99在线人妻在线中文字幕| 99riav亚洲国产免费| 国产成人啪精品午夜网站| 国产伦在线观看视频一区| 欧美激情久久久久久爽电影| 无遮挡黄片免费观看| 女性生殖器流出的白浆| 免费搜索国产男女视频| 国产aⅴ精品一区二区三区波| 天天一区二区日本电影三级| 两个人免费观看高清视频| 俄罗斯特黄特色一大片| 欧美日韩亚洲国产一区二区在线观看| 色婷婷久久久亚洲欧美| 波多野结衣av一区二区av| 免费看美女性在线毛片视频| 国产极品粉嫩免费观看在线| 黑人操中国人逼视频| 亚洲人成网站在线播放欧美日韩| 久久久水蜜桃国产精品网| 这个男人来自地球电影免费观看| 精品久久久久久成人av| 欧美绝顶高潮抽搐喷水| 久久国产精品男人的天堂亚洲| 99国产综合亚洲精品| 免费在线观看黄色视频的| 可以在线观看毛片的网站| 欧美日韩黄片免| 丰满的人妻完整版| 一边摸一边抽搐一进一小说| 欧美久久黑人一区二区| 亚洲国产欧美日韩在线播放| 成人免费观看视频高清| 精品一区二区三区四区五区乱码| 美国免费a级毛片| 国产极品粉嫩免费观看在线| 淫妇啪啪啪对白视频| 亚洲欧美精品综合一区二区三区| 久久国产乱子伦精品免费另类| 99在线人妻在线中文字幕| 中国美女看黄片| 久久人妻av系列| 国产精品98久久久久久宅男小说| 国产不卡一卡二| 国产黄色小视频在线观看| x7x7x7水蜜桃| 在线观看舔阴道视频| 日韩免费av在线播放| 国产精品一区二区三区四区久久 | 狂野欧美激情性xxxx| 看免费av毛片| 亚洲 国产 在线| 美女高潮喷水抽搐中文字幕| 国产av又大| 制服人妻中文乱码| 亚洲精品中文字幕在线视频| 日韩 欧美 亚洲 中文字幕| 亚洲七黄色美女视频| aaaaa片日本免费| 看片在线看免费视频| 悠悠久久av| 哪里可以看免费的av片| aaaaa片日本免费| 9191精品国产免费久久| 人人妻,人人澡人人爽秒播| 啦啦啦 在线观看视频| 久久亚洲真实| 高清在线国产一区| 欧美亚洲日本最大视频资源| 国产一区二区三区在线臀色熟女| 欧美中文日本在线观看视频| 国产精品亚洲一级av第二区| 中文字幕人成人乱码亚洲影| 日本 欧美在线| 午夜久久久久精精品| 一二三四在线观看免费中文在| av欧美777| 久久久久久国产a免费观看| 欧美色欧美亚洲另类二区| 99久久99久久久精品蜜桃| 成人国产一区最新在线观看| 黄色视频,在线免费观看| 亚洲欧美激情综合另类| 国产精品国产高清国产av| 精品国产美女av久久久久小说| 99久久无色码亚洲精品果冻| 久久久久久久久免费视频了| 国产亚洲精品av在线| 亚洲中文字幕一区二区三区有码在线看 | 欧美性猛交╳xxx乱大交人| 亚洲人成伊人成综合网2020| 韩国av一区二区三区四区| 成人三级黄色视频| 国产成人精品久久二区二区91| 美女扒开内裤让男人捅视频| 亚洲一区中文字幕在线| 18禁国产床啪视频网站| 又大又爽又粗| 一级毛片精品| 久久久国产欧美日韩av| 亚洲av第一区精品v没综合| 久久久久精品国产欧美久久久| 亚洲五月色婷婷综合| 婷婷六月久久综合丁香| 黄色片一级片一级黄色片| 国产精品野战在线观看| 国产在线观看jvid| 午夜福利视频1000在线观看| 亚洲国产日韩欧美精品在线观看 | cao死你这个sao货| 亚洲第一电影网av| 国产极品粉嫩免费观看在线| 天天一区二区日本电影三级| 变态另类丝袜制服| 国产精品亚洲av一区麻豆| 国产乱人伦免费视频| 精品国产美女av久久久久小说| 日本免费a在线| 成人国语在线视频| 成年版毛片免费区| 久久精品夜夜夜夜夜久久蜜豆 | 色精品久久人妻99蜜桃| 男女那种视频在线观看| 婷婷精品国产亚洲av在线| 俄罗斯特黄特色一大片| 亚洲第一电影网av| 人人妻,人人澡人人爽秒播| 亚洲av电影在线进入| 99久久国产精品久久久| 色老头精品视频在线观看| 午夜福利欧美成人| 国产成人影院久久av| 国产成人精品久久二区二区免费| 亚洲第一av免费看| 99精品久久久久人妻精品| 亚洲 欧美一区二区三区| 国产精品久久久久久亚洲av鲁大| a级毛片在线看网站| 高清在线国产一区| 亚洲专区字幕在线| 欧美黄色片欧美黄色片| 亚洲成人久久爱视频| 日韩精品中文字幕看吧| 欧美黄色片欧美黄色片| 无人区码免费观看不卡| 一边摸一边做爽爽视频免费| 国产精品自产拍在线观看55亚洲| 精品日产1卡2卡| 在线观看免费视频日本深夜| 天天添夜夜摸| 国产国语露脸激情在线看| 观看免费一级毛片| 久久国产亚洲av麻豆专区| 亚洲国产高清在线一区二区三 | 欧美人与性动交α欧美精品济南到| 在线看三级毛片| 50天的宝宝边吃奶边哭怎么回事| 免费在线观看视频国产中文字幕亚洲| 亚洲自偷自拍图片 自拍| www日本在线高清视频| 欧美精品亚洲一区二区| 手机成人av网站| 悠悠久久av| 国产不卡一卡二| 淫秽高清视频在线观看| 男女之事视频高清在线观看| 日韩av在线大香蕉| 午夜免费激情av| 母亲3免费完整高清在线观看| 夜夜爽天天搞| 国产精品久久视频播放| 日本一区二区免费在线视频| 午夜久久久久精精品| 天天一区二区日本电影三级| 最好的美女福利视频网| 变态另类丝袜制服| 禁无遮挡网站| 欧美国产日韩亚洲一区| 欧美绝顶高潮抽搐喷水| av中文乱码字幕在线| 亚洲熟妇熟女久久| 高清毛片免费观看视频网站| 久久热在线av| 黄色成人免费大全| 亚洲欧美日韩无卡精品| 免费高清在线观看日韩| 亚洲天堂国产精品一区在线| 午夜影院日韩av| a级毛片在线看网站| 午夜免费鲁丝| 中文字幕人成人乱码亚洲影| 18禁黄网站禁片午夜丰满| 曰老女人黄片| 听说在线观看完整版免费高清| 亚洲欧美精品综合久久99| 女同久久另类99精品国产91| 免费观看人在逋| 91大片在线观看| 97碰自拍视频| 国产亚洲av嫩草精品影院| 久久久久免费精品人妻一区二区 | 麻豆国产av国片精品| 中文字幕av电影在线播放| 国产亚洲精品第一综合不卡| 后天国语完整版免费观看| 人妻丰满熟妇av一区二区三区| 久久国产精品男人的天堂亚洲| 国产欧美日韩一区二区精品| 在线观看免费视频日本深夜| 国产黄色小视频在线观看| 91老司机精品| 男人舔女人下体高潮全视频| 亚洲成人久久爱视频| 欧美成人午夜精品| 91av网站免费观看| 首页视频小说图片口味搜索| 中文字幕精品免费在线观看视频| 精品久久久久久,| 两个人视频免费观看高清| 久久久久久亚洲精品国产蜜桃av| 亚洲一区二区三区不卡视频| 亚洲国产精品999在线| 欧美日韩亚洲综合一区二区三区_| 可以在线观看的亚洲视频| 久久国产亚洲av麻豆专区| 亚洲中文字幕一区二区三区有码在线看 | 丁香欧美五月| 大香蕉久久成人网| 脱女人内裤的视频| 国产亚洲精品久久久久久毛片| 国产一卡二卡三卡精品| 成年人黄色毛片网站| 一卡2卡三卡四卡精品乱码亚洲| 亚洲一区二区三区不卡视频| 亚洲国产精品sss在线观看| 国产欧美日韩精品亚洲av| 在线观看日韩欧美| 麻豆一二三区av精品| 神马国产精品三级电影在线观看 | 一个人观看的视频www高清免费观看 | 久久久国产欧美日韩av| 一本一本综合久久| 日韩精品中文字幕看吧| 亚洲av电影在线进入| 精品国产乱子伦一区二区三区| 一级a爱片免费观看的视频| 成人手机av| 国产精品1区2区在线观看.| 变态另类丝袜制服| 美女扒开内裤让男人捅视频| 97人妻精品一区二区三区麻豆 | 黄色毛片三级朝国网站| 18禁国产床啪视频网站| 亚洲第一欧美日韩一区二区三区| 91字幕亚洲| cao死你这个sao货| 好看av亚洲va欧美ⅴa在| 国产亚洲欧美在线一区二区| 美女扒开内裤让男人捅视频| 免费看a级黄色片| 日韩欧美免费精品| 一级毛片精品| 婷婷六月久久综合丁香| xxxwww97欧美| 丝袜在线中文字幕| 免费看十八禁软件| 国产精品亚洲美女久久久| 视频区欧美日本亚洲| 国产精品 欧美亚洲| 人人妻人人澡欧美一区二区| 青草久久国产| 欧美日本视频| 两个人免费观看高清视频| 精品日产1卡2卡| 亚洲 国产 在线| 在线观看午夜福利视频| 日韩精品中文字幕看吧| 免费人成视频x8x8入口观看| 黑人操中国人逼视频| 黄片小视频在线播放| 精品欧美一区二区三区在线| 国产精品免费视频内射| 人成视频在线观看免费观看| 一级毛片高清免费大全| 亚洲精品美女久久久久99蜜臀| 亚洲avbb在线观看| 亚洲全国av大片| 视频在线观看一区二区三区| 色播在线永久视频| 欧美性猛交黑人性爽| 国产在线精品亚洲第一网站| 国产免费av片在线观看野外av| 女警被强在线播放| 身体一侧抽搐| 禁无遮挡网站| 成人av一区二区三区在线看| 亚洲国产欧美一区二区综合| 国产精品99久久99久久久不卡| 制服丝袜大香蕉在线| 国产激情偷乱视频一区二区| 欧美黑人欧美精品刺激| 午夜免费鲁丝| 国产私拍福利视频在线观看| 日本三级黄在线观看| avwww免费| 欧美黄色淫秽网站| www.自偷自拍.com| 免费在线观看成人毛片| 欧美+亚洲+日韩+国产| 国产极品粉嫩免费观看在线| 亚洲精品中文字幕在线视频| 18禁裸乳无遮挡免费网站照片 | АⅤ资源中文在线天堂| 亚洲天堂国产精品一区在线| 精华霜和精华液先用哪个| 成年女人毛片免费观看观看9| 男人舔女人的私密视频| 这个男人来自地球电影免费观看| 成年版毛片免费区| 久久久久久久精品吃奶| 激情在线观看视频在线高清| 18禁美女被吸乳视频| 高潮久久久久久久久久久不卡| 国产av一区二区精品久久| 久久久久国产一级毛片高清牌| 2021天堂中文幕一二区在线观 | 757午夜福利合集在线观看| 91av网站免费观看| 成人永久免费在线观看视频| 91av网站免费观看| 久久精品亚洲精品国产色婷小说| 免费高清在线观看日韩| 欧美乱色亚洲激情| 老汉色av国产亚洲站长工具| 首页视频小说图片口味搜索| 国产野战对白在线观看| 欧美色欧美亚洲另类二区| 大型黄色视频在线免费观看| 久久香蕉激情| 日本免费一区二区三区高清不卡| 在线观看免费午夜福利视频| 国产一级毛片七仙女欲春2 | 夜夜看夜夜爽夜夜摸| 午夜免费激情av| 一级a爱视频在线免费观看| 久久香蕉激情| 一级毛片精品| 午夜福利一区二区在线看| 精品久久蜜臀av无| 久久国产乱子伦精品免费另类| 香蕉国产在线看| 国产一级毛片七仙女欲春2 | 亚洲一区中文字幕在线| 久久精品国产亚洲av高清一级| 一本综合久久免费| 男女之事视频高清在线观看| 日本撒尿小便嘘嘘汇集6| 天天躁夜夜躁狠狠躁躁| 欧美成人一区二区免费高清观看 | 亚洲性夜色夜夜综合| 亚洲欧美精品综合久久99| 精品福利观看| 变态另类丝袜制服| 法律面前人人平等表现在哪些方面| 欧美一区二区精品小视频在线| 香蕉丝袜av| 三级毛片av免费| www日本在线高清视频| 亚洲精品美女久久久久99蜜臀| 亚洲九九香蕉| 国内精品久久久久久久电影| 精品福利观看| 国产三级黄色录像| 成人亚洲精品av一区二区| 不卡一级毛片| 日韩欧美免费精品| 在线av久久热| 熟女电影av网| 国产精品av久久久久免费| 亚洲天堂国产精品一区在线| 欧美又色又爽又黄视频| 最近最新免费中文字幕在线| 看片在线看免费视频| 黄色女人牲交| √禁漫天堂资源中文www| 一a级毛片在线观看| 久久草成人影院| 99国产极品粉嫩在线观看| 91麻豆av在线| 欧美在线黄色| svipshipincom国产片| 午夜影院日韩av| 日本撒尿小便嘘嘘汇集6| 国产男靠女视频免费网站| 男人的好看免费观看在线视频 | 久久中文字幕人妻熟女| 91九色精品人成在线观看| 99国产综合亚洲精品| 国产精品国产高清国产av| 精品久久久久久久毛片微露脸| 国产亚洲精品av在线| 极品教师在线免费播放| 夜夜躁狠狠躁天天躁| 天天添夜夜摸| 波多野结衣高清无吗| 香蕉丝袜av| 男人舔奶头视频| 国产伦人伦偷精品视频| 很黄的视频免费| 精品少妇一区二区三区视频日本电影| 国产伦人伦偷精品视频| 国产激情欧美一区二区| 国产一区在线观看成人免费| 一级片免费观看大全| 中文资源天堂在线| 欧美激情极品国产一区二区三区| 精品久久久久久成人av| 999久久久国产精品视频| 亚洲五月天丁香| 俺也久久电影网| 两个人看的免费小视频| 免费看a级黄色片| 国产私拍福利视频在线观看| 国产亚洲欧美98| 悠悠久久av| 久久中文字幕一级| 欧美最黄视频在线播放免费| 男女做爰动态图高潮gif福利片| 99re在线观看精品视频| 国产成人精品久久二区二区91| 9191精品国产免费久久| 91麻豆精品激情在线观看国产| 久久久国产成人精品二区| 麻豆久久精品国产亚洲av| 俺也久久电影网| 中文字幕高清在线视频| 欧美激情极品国产一区二区三区| e午夜精品久久久久久久| 妹子高潮喷水视频| 又大又爽又粗| 国产成+人综合+亚洲专区| 啦啦啦 在线观看视频| 日韩大码丰满熟妇| 精品国产国语对白av| 久久久久国产精品人妻aⅴ院| 1024香蕉在线观看| 中出人妻视频一区二区| 老司机靠b影院| 一个人观看的视频www高清免费观看 | 曰老女人黄片| 国产精品,欧美在线| 人妻久久中文字幕网| 久久精品国产亚洲av香蕉五月| 一卡2卡三卡四卡精品乱码亚洲| 亚洲男人天堂网一区| 亚洲人成网站高清观看| 国产三级黄色录像| 国产伦在线观看视频一区| 神马国产精品三级电影在线观看 | cao死你这个sao货| 精品少妇一区二区三区视频日本电影| 亚洲av美国av| 好看av亚洲va欧美ⅴa在| 亚洲aⅴ乱码一区二区在线播放 | 久久精品国产综合久久久| 18禁美女被吸乳视频| 久久国产精品影院| 母亲3免费完整高清在线观看| ponron亚洲| 国产区一区二久久| 久久亚洲精品不卡| 日日夜夜操网爽| 成人亚洲精品一区在线观看| 成人18禁高潮啪啪吃奶动态图| 欧美性猛交黑人性爽| 国产精品亚洲av一区麻豆| 精品国产美女av久久久久小说| 免费观看精品视频网站| 国产精品自产拍在线观看55亚洲| 精品福利观看| 一级片免费观看大全| 久久婷婷成人综合色麻豆| 别揉我奶头~嗯~啊~动态视频| 亚洲一区二区三区色噜噜| 热99re8久久精品国产| 国内揄拍国产精品人妻在线 | 亚洲男人的天堂狠狠| 欧美日韩精品网址| 亚洲国产精品成人综合色| 亚洲电影在线观看av| 精品熟女少妇八av免费久了| 精品国产一区二区三区四区第35| 亚洲天堂国产精品一区在线| 90打野战视频偷拍视频| 日韩国内少妇激情av| 久久久久九九精品影院| 日韩成人在线观看一区二区三区| 在线观看午夜福利视频| 国产亚洲精品第一综合不卡| 91老司机精品| 天天躁狠狠躁夜夜躁狠狠躁| 国产人伦9x9x在线观看| 99热只有精品国产| 精品久久久久久久人妻蜜臀av| 午夜福利在线观看吧| 久久精品人妻少妇| 亚洲第一欧美日韩一区二区三区| 国产精品免费一区二区三区在线| 韩国av一区二区三区四区| 精品久久久久久久末码| 国产男靠女视频免费网站| 在线播放国产精品三级| 丝袜美腿诱惑在线| 日日摸夜夜添夜夜添小说| 日韩免费av在线播放| 久久久久国产一级毛片高清牌| 香蕉av资源在线| 国产99久久九九免费精品| 免费女性裸体啪啪无遮挡网站| 久久 成人 亚洲| 久久狼人影院| 日韩高清综合在线|