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

    RDX及其衍生物高溫?zé)峤獾姆磻?yīng)分子動力學(xué)模擬

    2017-05-12 06:58:02彭莉娟王靜波李澤榮李象遠(yuǎn)
    物理化學(xué)學(xué)報 2017年4期
    關(guān)鍵詞:含碳單胞衍生物

    彭莉娟 姚 倩 王靜波 李澤榮 朱 權(quán),* 李象遠(yuǎn),*

    RDX及其衍生物高溫?zé)峤獾姆磻?yīng)分子動力學(xué)模擬

    彭莉娟1姚 倩2王靜波3李澤榮2朱 權(quán)3,*李象遠(yuǎn)3,*

    (1四川大學(xué)空天科學(xué)與工程學(xué)院,成都610065;2四川大學(xué)化學(xué)學(xué)院,成都610064;3四川大學(xué)化學(xué)工程學(xué)院,成都610065)

    本文采用ReaxFF反應(yīng)分子動力學(xué)方法,研究了RDX及其衍生物晶體在高溫條件下(2000、2500和3000 K)的熱分解機(jī)理以及主要產(chǎn)物隨時間的變化情況。結(jié)果表明:RDX及其衍生物晶體熱解的第一步反應(yīng)均為N―NO2鍵斷裂生成NO2分子,隨后反應(yīng)由于六元環(huán)上和側(cè)鏈基團(tuán)的不同導(dǎo)致側(cè)鏈脫除或開環(huán)反應(yīng)的順序不同。在這四種單胞體系熱解中,NO2和NO分子為共同的中間產(chǎn)物,形成之后迅速發(fā)生次級反應(yīng)并最終生成N2。各體系熱解終產(chǎn)物一致,均為N2、H2O和CO2,其中N2分子數(shù)最多,大于20個。由于原始分子結(jié)構(gòu)和組成中C/N比、H/O比的不同,各體系熱解后H2O和CO2分子數(shù)目相差較大。不同溫度下,各單胞體系熱分解生成的最大含碳團(tuán)簇中C原子數(shù)均較小。在進(jìn)一步超胞體系的模擬中,RDX和RDX-D2體系生成的含碳團(tuán)簇中C原子數(shù)分別達(dá)到約30和16個,遠(yuǎn)高于單胞模擬,且受溫度影響較大;而RDX-D1和RDX-D3單胞或超胞模擬結(jié)果相近,均未生成含碳團(tuán)簇,僅存在小分子含碳碎片。因此,初始分子的結(jié)構(gòu)和元素比對含碳團(tuán)簇的生成有明顯影響。

    RDX及其衍生物;高溫?zé)岱纸猓籖eaxFF;含碳團(tuán)簇;分子動力學(xué)

    1 引言

    含能材料是一類儲存高能量的特殊材料,在一定外界條件刺激下,自身能在短時間內(nèi)發(fā)生氧化還原反應(yīng),釋放大量能量,在軍用和民用領(lǐng)域廣泛應(yīng)用1。高溫?zé)岱纸膺^程和反應(yīng)機(jī)理的研究,不僅可提高對含能材料反應(yīng)特性的認(rèn)識,而且還可在原子分子水平上進(jìn)行改進(jìn)和設(shè)計,從而研制出更好性能的炸藥和推進(jìn)劑配方。

    黑索金(RDX)是一種綜合性能優(yōu)良的單質(zhì)炸藥和推進(jìn)劑。通常RDX存在α和β兩種晶型,常溫常壓下α-RDX性質(zhì)穩(wěn)定而β-RDX極不穩(wěn)定,因此α-RDX具有更高的實用價值2,本文研究均采用α晶型。文獻(xiàn)中,針對RDX的熱分解過程開展了大量的理論和實驗研究。Goddard等3采用從頭算方法研究了單分子RDX熱分解過程,提出了三種主要初始分解途徑:協(xié)同分解反應(yīng),N―NO2的斷裂反應(yīng)和連續(xù)HONO消去反應(yīng),其中連續(xù)HONO消去反應(yīng)被認(rèn)為是能量最低的反應(yīng)途徑。趙繼軍等4研究了溫度為3000 K下RDX的熱分解行為,表明熱解后的主要產(chǎn)物為N2、H2O和CO2。Schweigert5采用從頭算分子動力學(xué)研究氣相單分子RDX在高溫2000 K下的詳細(xì)裂解反應(yīng)機(jī)理。Irikura6采用密度泛函理論計算比較了RDX單分子和雙分子熱解反應(yīng)途徑,為凝聚相反應(yīng)機(jī)理的研究提供了理論依據(jù)。Behrens和Bulusu7從實驗上指出RDX的主要分解路徑為RDX分子去除HNO和兩個HONO分子生成1,3,5-三嗪(TAZ)。雖然,氣相單分子RDX熱解機(jī)理研究已經(jīng)比較成熟,但凝聚相RDX由于高溫?zé)峤庵蟹肿娱g復(fù)雜的相互作用,反應(yīng)機(jī)理研究仍不夠深入。

    ReaxFF反應(yīng)力場8-10通過鍵長、鍵級和鍵能之間關(guān)系,以及原子間電荷的動態(tài)傳遞來模擬原子間化學(xué)鍵的斷裂和生成,為研究較大體系的物理化學(xué)性質(zhì)提供了有效途徑。ReaxFF反應(yīng)力場已經(jīng)廣泛應(yīng)用于熱分解11-14、燃燒15,16和爆轟17,18等快速反應(yīng)過程。采用這一方法,Strachan等19研究了不同溫度和密度下RDX的熱分解過程,表明產(chǎn)物中CO和CO2分子數(shù)目很大程度上與系統(tǒng)密度有關(guān),而N2和H2O分子的數(shù)目幾乎不受密度的影響。Zhang等20研究了不同溫度和密度下奧克托金(HMX)和三氨基三硝基苯(TATB)晶體的熱分解過程,表明HMX的初始裂解為N―NO2鍵的斷裂反應(yīng)和HONO消去反應(yīng);TATB在低溫下幾乎沒有發(fā)現(xiàn)C―NO2鍵的斷裂,高溫下主要為C―NO2鍵斷裂反應(yīng)和HONO消去反應(yīng)。但在RDX的熱解模擬中,關(guān)于反應(yīng)物分子結(jié)構(gòu)和元素組成對裂解產(chǎn)物和反應(yīng)途徑的影響,一直缺乏系統(tǒng)的對比分析,難以對含能材料的結(jié)構(gòu)設(shè)計提供指導(dǎo)。

    本工作根據(jù)RDX分子,改變六元環(huán)上和側(cè)鏈的取代基團(tuán)設(shè)計衍生物分子,從而改變含能材料中碳、氮含量比,研究其對反應(yīng)過程中能量釋放和含碳團(tuán)簇形成過程的影響。通過RDX取代基團(tuán)的逐步替代設(shè)計了三種化合物,它們的結(jié)構(gòu)分別如圖1所示,這三種衍生物的分子式分別為:C2H5O6N7(RDX-D1)、C3H7O6N7(RDX-D2)、C2H4O8N8(RDXD3)。相對于RDX,RDX-D1采用NH取代RDX主環(huán)上的CH2,改變了主環(huán)結(jié)構(gòu),提高了N元素含量;相對于RDX-D1,RDX-D2用―CH3取代了主環(huán)上NH基團(tuán)中的氫原子,改變了C元素含量和側(cè)鏈結(jié)構(gòu);相對于RDX-D1,RDX-D3用―NO2取代了主環(huán)上NH基團(tuán)中的氫原子,進(jìn)一步提高了N元素含量和改變了側(cè)鏈結(jié)構(gòu)。本工作采用ReaxFFMD方法,在單胞水平上模擬研究了高溫?zé)峤獾某跏挤磻?yīng)過程和反應(yīng)產(chǎn)物隨時間變化,在超胞水平上深入分析了反應(yīng)過程中含碳團(tuán)簇的演化情況,從而探討分子元素組成和分子結(jié)構(gòu)對RDX類化合物高溫?zé)峤庑阅艿挠绊憽?/p>

    圖1 RDX及其衍生物的分子結(jié)構(gòu)Fig.1 Structures of the RDX molecule and its derivatives

    2 計算模型與方法

    2.1 初始構(gòu)型的構(gòu)建

    RDX模擬采用的初始體系為一個單胞,依據(jù)RDX的X射線衍射單晶數(shù)據(jù)21,運用Materials Studio(MS)程序包Visualizer模塊,構(gòu)建單胞體系,如圖2(a)所示,該單胞結(jié)構(gòu)包含8個RDX分子,共計168個原子,初始密度為1.86 g·cm-3。由于RDX-D1與RDX分子具有相似的分子結(jié)構(gòu),故采用和RDX相同的單晶參數(shù),采用MS程序包Visualizer模塊構(gòu)建單胞體系,密度為1.87 g· cm-3。依據(jù)肖鶴鳴等22采用密度泛函理論計算結(jié)果,RDX-D3密度選取為1.92 g·cm-3。RDX-D2與RDX-D3具有相似的分子結(jié)構(gòu),故選用相同的參數(shù)。各體系單胞結(jié)構(gòu)如圖2所示,其中灰色、白色、藍(lán)色和紅色原子分別表示碳原子、氫原子、氮原子和氧原子,模擬均采用周期性邊界條件。

    針對含碳團(tuán)簇生成過程模擬,對RDX及其衍生物分別搭建了(2×2×1)超胞結(jié)構(gòu)。RDX超胞結(jié)構(gòu)共含有32個分子,共計672個原子;衍生物RDX-D1,RDX-D2和RDX-D3超胞結(jié)構(gòu)所含原子個數(shù)分別為640個、736個和704個。

    2.2 計算方法

    分子動力學(xué)模擬采用ReaxFF力場20,選取Velocity Verlet算法,采用ReaxFF-MD方法分別模擬各晶胞體系的高溫?zé)崃呀膺^程。在ReaxFF反應(yīng)力場中,各初始晶胞體系分別在低溫10 K下進(jìn)行幾何優(yōu)化,獲得最小能量結(jié)構(gòu),然后在相同溫度采用等溫等容(NVT)系綜進(jìn)行弛豫,使得模擬體系達(dá)到平衡狀態(tài),弛豫時間為5 ps,時間步長選取為0.1 fs。其后,各體系得到的平衡結(jié)構(gòu)分別在2000、2500、3000 K三個溫度下采取NVT系綜進(jìn)行高溫分解模擬,模擬時間步長選取0.1 fs,總模擬時間為100 ps。由于升溫過程可能會導(dǎo)致化學(xué)鍵的斷裂,因此,本文不考慮升溫過程的模擬。在所有的NVT-MD模擬中,各系統(tǒng)溫度均用Berendsen方法23進(jìn)行調(diào)節(jié),溫度調(diào)整參數(shù)選取100 fs,使溫度在設(shè)定值附近波動。當(dāng)任意原子對的鍵級大于0.3時,判定新化學(xué)鍵形成。

    3 結(jié)果與討論

    3.1 初始反應(yīng)途徑

    在3000 K條件下對四個單胞體系的熱分解進(jìn)行模擬,通過分析原子成鍵情況及各體系在不同時刻下的裂解產(chǎn)物,最終得到高溫?zé)峤獾某跏挤磻?yīng)途徑,如圖3所示。RDX及其三種衍生物熱解過程,第一步反應(yīng)均為N―NO2鍵發(fā)生斷裂生成NO2分子,表明主環(huán)和側(cè)鏈基團(tuán)的變化并不影響側(cè)鏈硝基斷鍵的初始反應(yīng)步驟。RDX分子第一步首先發(fā)生N―NO2鍵的斷裂,從而生成大量的NO2,與Botcher等在實驗中觀察到RDX快速熱解初期生成大量NO2一致24,25。在隨后反應(yīng)中,RDX主要發(fā)生剩余N―NO2鍵斷裂反應(yīng)和母體環(huán)斷裂反應(yīng),這一結(jié)果與Goddard等3采用密度泛函理論計算得到的RDX中N―N鍵斷裂(解離能約83.68 kJ·mol-1)優(yōu)先于協(xié)同開環(huán)反應(yīng)的結(jié)論一致。

    圖2 RDX(a),RDX-D1(b),RDX-D2(c)和RDX-D3(d)單胞結(jié)構(gòu)示意圖Fig.2 Unit cell structure of RDX(a),RDX-D1(b),RDX-D2(c)and RDX-D3(d) color online

    圖3 RDX及其衍生物的初始反應(yīng)途徑Fig.3 Initial reaction path for RDX and its derivatives

    RDX由于分子本身的對稱性,側(cè)鏈第一和第二個NO2基的消去無選擇性,但在RDX-D1中,由于六元環(huán)上的CH2基團(tuán)被NH基團(tuán)取代,在熱解過程中優(yōu)先脫去NH基團(tuán)鄰位氮原子上的兩個側(cè)鏈NO2基,其后,不同于RDX的開環(huán)反應(yīng),RDX-D1發(fā)生HONO脫去反應(yīng),之后再開環(huán)。在RDX-D1基礎(chǔ)上將環(huán)上NH中的氫原子用甲基取代,RDXD2中兩個NO2基脫除也優(yōu)先發(fā)生于N―CH3鄰位氮原子上,其后發(fā)生的是甲基脫除反應(yīng),然后發(fā)生開環(huán)反應(yīng)。因此,與RDX-D1相比,CH3基的引入對鄰位N―NO2的斷鍵影響不大,而對對位NO2的脫除影響較大,使得后續(xù)的反應(yīng)途徑完全不同。將RDX-D2中的甲基用硝基取代,RDX-D3中相鄰的兩個NO2基優(yōu)先脫除,六元環(huán)上形成氮氮雙鍵,其后雙鍵間位上的NO2基發(fā)生脫除,然后發(fā)生開環(huán)反應(yīng)。

    在RDX及其衍生物晶胞熱解過程中主要觀察到三種初始分解機(jī)理:N―NO2鍵的斷裂反應(yīng),HONO消去反應(yīng)和分子主環(huán)的開環(huán)反應(yīng),與Goddard等采用密度泛函理論計算的結(jié)果一致3,22。對比這四種單胞,第一步反應(yīng)均為N―NO2鍵斷裂,表明這一通道是能量最低的反應(yīng)通道,與Imamura等研究的其它硝胺化合物的熱分解機(jī)理一致26-29。然后才出現(xiàn)HONO消去、甲基取代基脫除和開環(huán)反應(yīng),進(jìn)一步生成小分子碎片。在四種體系的高溫?zé)峤庵?,?cè)鏈取代基團(tuán)對N―NO2第一步斷鍵無影響,主要改變了隨后的側(cè)鏈脫除或開環(huán)反應(yīng)的順序。

    3.2 溫度對勢能和總分子碎片數(shù)的影響

    在不同溫度條件下,我們模擬了各單胞體系(經(jīng)過5 ps初始弛豫)勢能隨時間的變化,如圖4所示。在初始熱解時,系統(tǒng)迅速吸收熱量發(fā)生初級熱解反應(yīng),使得系統(tǒng)內(nèi)勢能升高,隨后的次級反應(yīng)生成一系列中間及最終產(chǎn)物,并放出大量的熱,因而,系統(tǒng)勢能發(fā)生明顯下降。溫度越高,系統(tǒng)達(dá)到化學(xué)反應(yīng)動態(tài)平衡所需時間越短,熱解反應(yīng)越完全,生成的產(chǎn)物越多,勢能下降也越多。因而,各體系勢能的改變與其熱分解反應(yīng)過程中的化學(xué)鍵斷裂和生成緊密相關(guān)30。

    從圖4(a)可見,系統(tǒng)溫度越低,初始反應(yīng)越慢,吸熱階段越長;隨著系統(tǒng)溫度的升高,各體系分子初始反應(yīng)加劇,從而導(dǎo)致次級反應(yīng)越早開始,系統(tǒng)達(dá)到化學(xué)反應(yīng)動態(tài)平衡的時間縮短。當(dāng)溫度為3000 K時,各晶胞體系勢能在20 ps后相對穩(wěn)定;而當(dāng)系統(tǒng)處于2000 K的較低溫度時,各體系勢能在100 ps內(nèi)依然緩慢減少,需要更長時間才能達(dá)到穩(wěn)定,故在低溫條件達(dá)到化學(xué)反應(yīng)平衡的時間更長。

    圖5為各單胞體系在不同溫度下熱解過程中總分子碎片數(shù)隨時間的變化情況。從圖5(a)可見,在高溫初始裂解階段,溫度越高RDX熱解釋放能量的速率越快,物種數(shù)增加越快、產(chǎn)生的總分子碎片數(shù)越多,在10 ps左右達(dá)到峰值后,這些快速熱解產(chǎn)生的分子碎片會進(jìn)一步轉(zhuǎn)化,但整體數(shù)量維持在50左右。其他衍生物具有類似的變化規(guī)律,各晶胞體系的總分子碎片數(shù)在前10 ps內(nèi)快速增加,隨后變化相對比較緩慢,且溫度越高,總分子碎片數(shù)增加速度越大。在50 ps后,各晶胞體系的總分子碎片數(shù)在3000 K溫度條件下基本維持恒定,而在較低的溫度2000 K條件下仍緩慢增加。綜上,隨熱解溫度升高,各晶胞體系產(chǎn)生的總分子碎片數(shù)變多,總分子碎片數(shù)增加到最大值的時間縮短,體系釋放的能量也更多。

    3.3 主產(chǎn)物及中間產(chǎn)物分析

    圖4 四個體系的勢能隨時間變化的關(guān)系Fig.4 Time evolution of potential energy for the four materials

    圖5 四個體系的熱解總分子碎片數(shù)隨時間變化關(guān)系Fig.5 Time evolution of total fragments for the four materials

    RDX及其衍生物在3000 K時的主要分解產(chǎn)物N2、H2O和CO2及中間產(chǎn)物NO和NO2分子數(shù)目隨時間的變化如圖6所示。在各體系主產(chǎn)物中,N2是生成數(shù)目最多的產(chǎn)物,主要來自次級反應(yīng)中NO2和NO的轉(zhuǎn)化,其主要機(jī)理為NO2或NO與含碳分子碎片反應(yīng)生成N2。在前20 ps內(nèi),NO2和NO均迅速增大然后減少,此過程中N2數(shù)目一直迅速增加到約20,說明消耗掉的大部分NO2和NO生成了N2。在這四種不同體系中,主產(chǎn)物分子數(shù)目相對較多的還有H2O和CO2,這三種主要產(chǎn)物數(shù)目的大小順序為:N2>H2O>CO2。

    在四種不同的晶胞體系中,主要分解產(chǎn)物H2O和CO2分子數(shù)目相差較大。表1列出了各體系分子的元素組成比以及主產(chǎn)物CO2和H2O分子的最終數(shù)目。從表1可以清晰看出RDX和RDX-D2分子中C/N比大于等于0.5,高于RDX-D1和RDX-D3分子中的C/N比;RDX和RDX-D2體系中碳原子含量相多較高,熱解穩(wěn)定后生成的CO2分子數(shù)均小于10,這說明體系中的碳原子的最終轉(zhuǎn)化途徑并不以生成穩(wěn)定產(chǎn)物CO2為主,而是轉(zhuǎn)變成了其它的含碳分子碎片。與此對應(yīng),RDX和RDX-D2由于其分子構(gòu)成的H/O元素比(≥1)高于另外兩種衍生物,熱解后生成H2O分子數(shù)大于20,H2O的生成消耗了大量的氧,同時也抑制了CO2的生成。說明在初始的分子設(shè)計中C/N比和H/O比,對熱解反應(yīng)路徑和最終產(chǎn)物分布起決定作用。

    仔細(xì)觀察圖6各體系中N2和CO2分子數(shù)目,發(fā)現(xiàn)它們幾乎在同一時間達(dá)到最大,這與N2和CO2的生成途徑有關(guān);反應(yīng)早期的主要步驟是N―NO2的斷裂生成NO2,所以剛開始體系中NO2數(shù)量最多,隨后NO2進(jìn)一步反應(yīng)生成N2,而NO2中的O則與中間物CO反應(yīng)生成CO2,因而N2和CO2達(dá)到最大值的時間一致。同時還表明,高溫下N―N兩原子主要以N2的形式存在,與Strachan等10,19的研究結(jié)果一致。此外,在反應(yīng)過程中還觀察到一些數(shù)量很少基團(tuán)出現(xiàn),如:H、HO、CO、HNO3、N2O和HONO等基團(tuán),這些小分子基團(tuán)相互聚合或反應(yīng)生成其它產(chǎn)物;因其數(shù)量較少,故在本文沒有對其進(jìn)行詳細(xì)討論。

    由于圖6很難清晰展示中間體NO2和NO的變化,而中間體的變化影響著整個熱分解反應(yīng)過程,故研究了前20 ps中間體NO2和NO隨時間的變化,如圖7所示。RDX-D2中的NO2峰值最高,RDX-D3中的NO2峰值最低,RDX和RDX-D2中的NO2峰值(>13)略高于RDX-D1和RDX-D3的NO2峰值(〈10)。從圖7(a)中可以看出,NO2的數(shù)量在前2.5 ps內(nèi)迅速達(dá)到峰值,表明四種物質(zhì)首先發(fā)生斷裂的鍵是N―NO2鍵,這與圖3中研究表明高溫下N―NO2鍵的斷裂是各體系最主要的初始反應(yīng)途徑相一致。隨后NO2的數(shù)量開始減少,這是由于中間產(chǎn)物NO2參與次級反應(yīng)不斷被消耗,生成NO,NO3和HONO等中間產(chǎn)物。

    圖6 四個體系的主要產(chǎn)物分子(N2,H2O,CO2,NO2,NO)數(shù)目隨時間變化關(guān)系Fig.6 Time evolution of key produced molecules(N2,H2O,CO2,NO2,NO)for the four materials

    同時,NO也是很活躍的中間產(chǎn)物,從圖7(b)可以看出,NO和NO2的分布曲線具有相似的演化特征,均是先增大后減小的過程,但各體系NO達(dá)到峰值前的增長速率以及隨后的衰減速率相比NO2較小。在整個反應(yīng)過程中,NO數(shù)量整體上明顯低于NO2,且NO分子數(shù)達(dá)到最大值的時間相對NO2滯后一些,這是因為后續(xù)反應(yīng)消耗掉一部分NO2生成NO。NO的產(chǎn)生主要有三種途徑,(1)N―NO2→N―ONO,然后NO―NO單鍵斷裂生成NO;(2)氫原子轉(zhuǎn)移反應(yīng)生成的HONO,隨后分解為HO和NO分子;(3)由NO2和NO3等相互反應(yīng)生成NO。在20 ps后,NO2和NO的數(shù)量很少,幾乎被完全消耗。

    表1 各體系的C/N和H/O比值及主產(chǎn)物CO2和H2O的生成量Table 1 C/N and H/O ratios of each system and the amount of the main products of CO2and H2O

    3.4 大分子含碳團(tuán)簇的形成

    實驗研究表明,富碳高能材料在爆轟過程中能形成大分子含碳團(tuán)簇30-34。大分子含碳團(tuán)簇的出現(xiàn)是深入了解含能材料性質(zhì)的重要方面,目前采用分子動力學(xué)模擬,已經(jīng)開展含能材料熱解過程含碳團(tuán)簇形成動力學(xué)機(jī)制的研究35。為了深入了解含碳團(tuán)簇形成過程,首先模擬了各單胞體系高溫?zé)峤膺^程最大含碳團(tuán)簇C原子數(shù)隨時間的演變,如圖8(a,b,c)所示。本文中我們將C原子數(shù)大于8定義為含碳團(tuán)簇,隨著溫度升高,RDX和RDX-D2單胞體系熱解過程中最大含碳團(tuán)簇C原子數(shù)有所增加,在3000 K時含碳團(tuán)簇最大碳原子數(shù)可達(dá)19。而在RDX-D1和RDX-D3單胞體系熱解過程中,均未生成含碳團(tuán)簇,只存在小分子含碳碎片(C原子數(shù)小于等于3),且在不同溫度下幾乎保持不變。由于大分子含碳團(tuán)簇的生成受初始模型的影響較大,因此,為了正確考慮其形成過程,我們進(jìn)一步采用超胞體系進(jìn)行模擬,如圖8(d)所示。

    圖7 四個體系熱解過程中NO2(a)和NO(b)的分子數(shù)量變化關(guān)系Fig.7 Quantities of NO2(a)and NO(b)for the thermal pyrolysis of four systems

    圖8 不同溫度下各體系單胞(a-c)和超胞(d)熱解過程中最大含碳團(tuán)簇所含C原子數(shù)隨時間變化關(guān)系Fig.8 Time evolution of C atomic number in the maximum molecule formed in the thermal pyrolysis of unit cell(a-c)and super cell(d)under different temperatures

    對比圖8(c,d)單胞和多胞體系可以看出,RDX和RDX-D2超胞體系中生成的最大含碳團(tuán)簇C原子數(shù)增加迅速,最終C原子數(shù)目分別達(dá)到約30和16個;RDX-D1和RDX-D3單胞和超胞體系的模擬結(jié)果相近,均未生成大分子含碳團(tuán)簇,生成的小分子含碳碎片中C原子數(shù)均小于等于3??紤]初始分子結(jié)構(gòu)中,RDX和RDX-D2分子中原子比C/N≥0.5,而RDX-D1和RDX-D3分子中C/N〈0.5,說明:當(dāng)分子中C/N≥0.5時,生成的最大含碳團(tuán)簇C原子數(shù)與模型大小有關(guān);而對C/N〈0.5的體系,含碳團(tuán)簇的形成過程基本不受模型大小的影響。

    因而,初始設(shè)計的分子中碳、氮含量比和反應(yīng)物分子結(jié)構(gòu)對RDX及其衍生物中最大含碳團(tuán)簇的形成過程影響較大,對最大含碳團(tuán)簇中C原子數(shù)目多少也影響顯著,且C/N比越大越容易形成大分子含碳團(tuán)簇。根據(jù)上述模擬結(jié)果,在熱分解的初級反應(yīng)中,環(huán)上脫落NO2、NO和CH3等自由基后,剩余部分發(fā)生聚合形成大分子含碳團(tuán)簇,這與實驗和理論上研究硝胺類化合物的高溫?zé)峤膺^程是一致的20,36。

    4 結(jié)論

    采用ReaxFF反應(yīng)力場分子動力學(xué)方法,模擬研究了RDX及其衍生物單胞和超胞在不同溫度下的熱分解行為。研究表明,這四種晶體高溫?zé)峤獾牡谝徊骄鶠镹―NO2鍵的斷裂生成NO2,然后發(fā)生取代基的脫去和環(huán)的斷裂反應(yīng)。RDX在斷裂兩個硝基后就發(fā)生開環(huán)反應(yīng),而RDX的衍生物在斷裂兩個硝基后并沒有發(fā)生開環(huán)反應(yīng),而是先發(fā)生取代基的脫去反應(yīng)再發(fā)生開環(huán)反應(yīng),這是由于各體系分子六元環(huán)上基團(tuán)和側(cè)鏈基團(tuán)的改變使得初始反應(yīng)途徑不同。

    通過對各體系熱分解過程主要產(chǎn)物和中間產(chǎn)物進(jìn)行詳細(xì)分析,發(fā)現(xiàn)NO2和NO是裂解過程中十分重要的中間產(chǎn)物。NO2主要來自N―NO2中的N―N鍵的斷裂產(chǎn)生;NO主要來自N―NO2→N―ONO,然后NO―NO單鍵斷裂生成NO,中間產(chǎn)物HONO分解為HO和NO以及次級反應(yīng)生成NO分子。熱分解后最穩(wěn)定的氣相產(chǎn)物為N2、H2O和CO2,其中分子數(shù)最多的是N2(>20),主要來源于NO2和NO參與的次級反應(yīng)中。各分子體系中的C/ N比和H/O比不同,使得主要分解產(chǎn)物H2O和CO2分子數(shù)目相差較大。C/N和H/O比越大的體系,熱解后生成的H2O分子數(shù)越多,H2O的生成消耗了大量的氧,故生成穩(wěn)定產(chǎn)物CO2分子數(shù)越少,生成含碳團(tuán)簇分子量越大。模型分子設(shè)計中C/N比和H/O比,對于最終熱解過程和產(chǎn)物分布起決定作用。

    單胞體系模擬分析表明:溫度升高有利于RDX和RDX-D2單胞體系熱解中最大含碳團(tuán)簇的形成,在3000 K時初步形成大分子含碳團(tuán)簇(C原子數(shù)約為10)。進(jìn)一步超胞體系模擬結(jié)果表明:相對于單胞模擬,RDX和RDX-D2超胞熱解過程中最大含碳團(tuán)簇C原子數(shù)顯著增加(約為30和16),初始模型大小和反應(yīng)溫度對含碳團(tuán)簇的生成影響明顯。而RDX和RDX-D2體系中的C/N比明顯高于RDX-D1和RDX-D3體系,對于后兩者的模擬,單胞和超胞模擬結(jié)果相近,均未生成含碳團(tuán)簇,僅存在小分子含碳碎片(碳原子數(shù)小于等于3)。含碳團(tuán)簇的大量生成,不利于含能材料能量的釋放,故分子設(shè)計中減小C/N比能很好地改善含能材料的熱解性能。

    (1) Shaw,R.W.;Brill,T.B.;Thompson,D.L.Overviews of Recent Research on Energetic Materials;World Scientific:NJ,2005, Vol.16.

    (2) Teipel,U.Energetic Materials:Particle Processing and Characterization;Wiley-VCH:Weinheim,Germany,2005.

    (3) Chakraborty,D.;Muller,R.P.;Dasgupta,S.;Goddard,W.A.J. Phys.Chem.A 2000,104,2261.doi:10.1021/jp9936953

    (5) Schweigert,I.V.J.Phys.Chem.A 2015,119,2747. doi:10.1021/jp510034p

    (6) Irikura,K.K.J.Phys.Chem.A 2013,117,2233.doi:10.1021/ jp310247z

    (7) Behrens,R.,Jr.;Bulusu,S.J.Phys.Chem.1992,96,8877. doi:10.1063/1.462245

    (8) Van Duin,A.C.;Dasgupta,S.;Lorant,F.;Goddard,W.A.J. Phys.Chem.A 2001,105,9396.doi:10.1021/jp004368u

    (9) Cohen,R.;Zeiri,Y.;Wurzberg,E.;Kosloff,R.J.Phys.Chem.A 2007,111,11074.doi:10.1021/jp072121s

    (10) Strachan,A.;van Duin,A.C.;Chakraborty,D.;Dasgupta,S.; Goddard,W.A.,III Phys.Rev.Lett.2003,91,098301. doi:10.1103/PhysRevLett.91.098301

    (11) Zhou,T.T.;Shi,Y.D.;Huang,F.L.Acta Phys.-Chim.Sin. 2012,28(11),2605.[周婷婷,石一丁,黃風(fēng)雷.物理化學(xué)學(xué)報, 2012,28(11),2605.].doi:10.3866/PKU.WHXB201208031

    (12) Zhang,L.;Chen,L.;Wang,C.;Wu,J.Y.Acta Phys.-Chim.Sin. 2013,29(6),1145.[張 力,陳 朗,王 晨,伍俊英.物理化學(xué)學(xué)報,2013,29(6),1145.].doi:10.3866/PKU. WHXB201303221

    (13) Rom,N.;Zybin,S.V.;van Duin,A.C.;Goddard,W.A.,III.; Zeiri,Y.;Katz,G.;Kosloff,R.J.Phys.Chem.A 2011,115, 10181.doi:10.1021/jp202059v

    (14) Liu,H.;Dong,X.;He,Y.H.Acta Phys.-Chim.Sin.2014,30, 232.[劉 海,董 曉,何遠(yuǎn)航.物理化學(xué)學(xué)報,2014,30,232.] doi:10.3866/PKU.WHXB201312101

    (15) Qian,H.J.;van Duin,A.C.;Morokuma,K.;Irle,S.J.Chem. Theory Comput.2011,7,2040.doi:10.1021/ct200197v

    (16) Liu,L.;Bai,C.;Sun,H.;Goddard,W.A.,III.J.Phys.Chem.A 2011,115,4941.doi:10.1021/jp110435p

    (17)Ge,N.N.;Wei,Y.K.;Ji,G.F.;Chen,X.R.;Zhao,F.;Wei,D. Q.J.Phys.Chem.B 2012,116,13696.doi:10.1021/jp309120t

    (18) Zhang,L.Z.;Zybin,S.V.;Van Duin,A.C.T.;Goddard,W.A., III.J.Energ.Mater.2010,28,92.doi:10.1080/ 07370652.2010.504682

    (19) Strachan,A.;Kober,E.M.;van Duin,A.C.;Oxgaard,J.; Goddard,W.A.,III.J.Chem.Phys.2005,122,054502. doi:10.1063/1.1831277

    (20) Zhang,L.;Zybin,S.V.;van Duin,A.C.;Dasgupta,S.; Goddard,W.A.,III;Kober,E.M.J.Phys.Chem.A 2009,113, 10619.doi:10.1021/jp901353a

    (21) Hakey,P.;Ouellette,W.;Zubieta,J.;Korter,T.Acta Crystallographica Section E:Structure Reports Online 2008,64 (8),o1428.doi:10.1107/S1600536808019727

    (22) Qiu,L.;Gong,X.;Ju,X.;Xiao,H.Sci.China Ser.B:Chem. 2008,51,1231.doi:10.1007/s11426-008-0141-1

    (23) Berendsen,H.J.;Postma,J.V.;van Gunsteren,W.F.;DiNola, A.R.H.J.;Haak,J.R.J.Chem.Phys.1984,81,3684. doi:10.1063/1.448118

    (24) Botcher,T.R.;Wight,C.A.J.Phys.Chem.1994,98,5441. doi:10.1021/j100072a009

    (25) Oyumi,Y.;Brill,T.B.Combust.Flame 1985,62,213. doi:10.1016/0010-2180(85)90147-6

    (26)Kohno,Y.;Ueda,K.;Imamura,A.J.Phys.Chem.1996,100, 4701.doi:10.1021/jp9503223

    (27)Prabhakaran,K.V.;Bhide,N.M.;Kurian,E.M.Thermochim. Acta 1995,249,249.doi:10.1016/0040-6031(95)90704-1

    (28) Oxley,J.C.;Kooh,A.B.;Szekeres,R.;Zheng,W.J.Phys. Chem.1994,98,7004.doi:10.1021/j100079a019

    (29) Brill,T.B.;Oyumi,Y.J.Phys.Chem.1986,90,6848. doi:10.1021/j100284a028

    (30) Mironov,E.V.;Petrov,E.A.;Korets,A.Y.Combust.Explos. Shock Waves.2004,40,473.doi:10.1023/B: CESW.0000033571.82326.6a

    (31)Anisichkin,V.F.;Dolgushin,D.S.;Petrov,E.A.Combust. Explos.Shock Waves.1995,31,106.doi:10.1007/BF00755966

    (32) Titov,V.M.;Anisichkin,V.F.;Mal′kov,I.Y.Combust.Explos. Shock Waves.1989,25,372.doi:10.1007/BF00788819

    (33) Baidakova,M.J.Phys.D:Appl.Phys.2007,40,6300. doi:10.1088/0022-3727/40/20/S14

    (34) Mal'Kov,I.Y.;Filatov,L.I.;Titov,V.M.;Litvinov,B.V.; Chuvilin,A.L.;Teslenko,T.S.Combust.Explos.Shock Waves. 1993,29,542.doi:10.1007/BF00782983

    (35) Chevrot,G.;Sollier,A.;Pineau,N.J.Chem.Phys.2012,136, 084506.doi:10.1063/1.3686750

    (36) Hobbs,M.L.;Kaneshige,M.J.;Gilbert,D.W.;Marley,S.K.; Todd,S.N.J.Phys.Chem.A 2009,113,10474.doi:10.1021/ jp906134f

    Pyrolysis of RDX and Its Derivatives via Reactive Molecular Dynamics Simulations

    PENG Li-Juan1YAO Qian2WANG Jing-Bo3LI Ze-Rong2ZHU Quan3,*LI Xiang-Yuan3,*
    (1School of Aeronautics&Astronautics,Sichuan University,Chengdu 610065,P.R.China; 2College of Chemistry,Sichuan University,Chengdu 610064,P.R.China; 3School of Chemical Engineering,Sichuan University,Chengdu 610065,P.R.China)

    The thermal decomposition mechanisms of cyclotrimethylene trinitramine(RDX)and its derivatives are investigated using ReaxFF reactive molecular dynamics simulation at high temperatures(2000,2500,and 3000 K).It is shown that the first pyrolysis step of RDX and its derivatives is the N―NO2homolytic cleavage to produce NO2;however,the subsequent reaction mechanism is completely different due to the different sixmembered rings and side chain groups.In these four model systems,NO2and NO molecules are common intermediates in the thermal decomposition process,which eventually transform to N2in the following reactions. The most stable products are N2,H2O,and CO2in the thermal process,among which N2has the maximum molecular number(more than 20).The numbers of cracked H2O and CO2molecules are very different in the four model systems due to the different C/N and H/O ratios.At different temperatures,the maximum numbers of carbon atoms in the carbon clusters are all small in the unit cell simulations of the four systems.In the further super cell simulation for RDX and RDX-D2 the numbers of carbon atoms reach about 30 and 16,respectively, which are higher than those from the unit cell simulation.In both the unit and super cell simulations for RDX-D1 and RDX-D3,the carbon cluster cannot be formed and there exist only small carbon molecular fragments. Therefore,the initial molecular structure and elemental ratio have a large effect on the production of carbon clusters.

    RDX and its derivatives;Thermal decomposition;ReaxFF;Carbon clusters;Molecular dynamics

    O643

    Xu,J.;Zhao,J.;Sun,L.Mol.Simul.2008,34,961.

    10.1080/ 08927020802162892

    doi:10.3866/PKU.WHXB201701161

    Received:October 22,2016;Revised:January 16,2017;Published online:January 16,2017.

    *Corresponding authors.ZHU Quan,Email:qzhu@scu.edu.cn;Tel:+86-28-85402951.LI Xiang-Yuan,Email:xyli@scu.edu.cn.

    The project was supported by the National Natural Science Foundation of China(91441132)and Program for New Century Excellent Talents in University of Ministry of Education of China(NCET-13-0398).國家自然科學(xué)基金(91441132)和教育部新世紀(jì)優(yōu)秀人才支持計劃(NCET-13-0398)資助項目

    猜你喜歡
    含碳單胞衍生物
    基于NURBS的點陣材料參數(shù)化建模方法
    復(fù)合材料周期結(jié)構(gòu)數(shù)學(xué)均勻化方法的一種新型單胞邊界條件
    考慮界面層影響的三維機(jī)織復(fù)合材料單胞模型研究
    黑龍江省造林樹種含碳率與土壤性質(zhì)研究
    森林工程(2018年4期)2018-08-04 03:23:20
    烴的含氧衍生物知識鏈接
    新型螺雙二氫茚二酚衍生物的合成
    鐵礦含碳球團(tuán)還原后性能的實驗研究
    Xanomeline新型衍生物SBG-PK-014促進(jìn)APPsw的α-剪切
    基于B-H鍵的活化對含B-C、B-Cl、B-P鍵的碳硼烷硼端衍生物的合成與表征
    內(nèi)蒙古某含碳低品位硫化鉛鋅礦石選礦試驗
    金屬礦山(2013年6期)2013-03-11 16:53:57
    午夜免费观看性视频| 高清日韩中文字幕在线| 国产又色又爽无遮挡免| 国产亚洲一区二区精品| 午夜精品在线福利| av卡一久久| 嫩草影院精品99| 丰满人妻一区二区三区视频av| 身体一侧抽搐| 成人欧美大片| 免费黄网站久久成人精品| 午夜激情福利司机影院| 视频中文字幕在线观看| 日日撸夜夜添| 日韩欧美精品v在线| 国产伦在线观看视频一区| 国产在视频线在精品| 建设人人有责人人尽责人人享有的 | 国产男女超爽视频在线观看| 精品99又大又爽又粗少妇毛片| 在现免费观看毛片| 久久午夜福利片| 成人午夜高清在线视频| 国产视频首页在线观看| 亚洲内射少妇av| 菩萨蛮人人尽说江南好唐韦庄| 精品久久久久久久久亚洲| 成人毛片a级毛片在线播放| 亚洲国产高清在线一区二区三| 国产一区二区三区综合在线观看 | 青春草视频在线免费观看| 一级av片app| 好男人在线观看高清免费视频| 欧美成人一区二区免费高清观看| 深爱激情五月婷婷| 国产精品三级大全| 日韩伦理黄色片| 国内少妇人妻偷人精品xxx网站| 免费大片黄手机在线观看| 欧美高清成人免费视频www| 久久这里有精品视频免费| 免费黄频网站在线观看国产| 韩国av在线不卡| 好男人在线观看高清免费视频| 美女内射精品一级片tv| 久久韩国三级中文字幕| 国产高清有码在线观看视频| 国产精品日韩av在线免费观看| 最近最新中文字幕免费大全7| 日本午夜av视频| 国产精品国产三级专区第一集| 国产永久视频网站| 22中文网久久字幕| 国产黄色免费在线视频| 美女脱内裤让男人舔精品视频| 国产熟女欧美一区二区| 欧美区成人在线视频| 大香蕉久久网| www.av在线官网国产| 九草在线视频观看| 国产探花极品一区二区| 亚洲在线自拍视频| 国产一区有黄有色的免费视频 | 免费av观看视频| 久久久久久久久中文| 国产午夜精品一二区理论片| 日韩视频在线欧美| 免费人成在线观看视频色| 日韩一区二区视频免费看| 午夜精品一区二区三区免费看| 熟妇人妻不卡中文字幕| 亚洲不卡免费看| 熟妇人妻不卡中文字幕| 精品99又大又爽又粗少妇毛片| 精品久久久久久久久久久久久| 国产一区有黄有色的免费视频 | 99久久九九国产精品国产免费| 18+在线观看网站| 亚洲,欧美,日韩| 日韩国内少妇激情av| 欧美人与善性xxx| 人人妻人人看人人澡| 日韩视频在线欧美| 久久精品熟女亚洲av麻豆精品 | 亚洲精品aⅴ在线观看| 身体一侧抽搐| 婷婷色麻豆天堂久久| 亚洲欧美精品专区久久| 2021天堂中文幕一二区在线观| 国产黄频视频在线观看| 3wmmmm亚洲av在线观看| 97人妻精品一区二区三区麻豆| 91aial.com中文字幕在线观看| 亚洲精品成人久久久久久| 看非洲黑人一级黄片| 美女黄网站色视频| 亚州av有码| 精品午夜福利在线看| 午夜激情久久久久久久| 国产成人福利小说| 一级毛片久久久久久久久女| 国产av国产精品国产| 日韩av在线大香蕉| 国产免费视频播放在线视频 | 联通29元200g的流量卡| 日韩电影二区| 免费观看精品视频网站| 久久久精品免费免费高清| 午夜免费男女啪啪视频观看| 欧美人与善性xxx| 18禁在线无遮挡免费观看视频| 亚洲av成人av| 久久99热这里只有精品18| 免费av毛片视频| 男女边吃奶边做爰视频| 国产av国产精品国产| 一个人观看的视频www高清免费观看| 国产精品久久久久久av不卡| 一个人观看的视频www高清免费观看| 日本免费a在线| 欧美日韩在线观看h| 国产永久视频网站| 波野结衣二区三区在线| 免费av不卡在线播放| 亚洲成人av在线免费| 久久久久久久久久久免费av| 免费无遮挡裸体视频| 国产熟女欧美一区二区| 夫妻性生交免费视频一级片| 国产乱人视频| 亚洲精品成人久久久久久| 亚洲av成人精品一二三区| 国产男人的电影天堂91| xxx大片免费视频| 色哟哟·www| 乱人视频在线观看| 久久6这里有精品| 大香蕉97超碰在线| 免费观看的影片在线观看| 午夜激情福利司机影院| 欧美3d第一页| 国产单亲对白刺激| 国产久久久一区二区三区| 又爽又黄a免费视频| 啦啦啦韩国在线观看视频| 免费看光身美女| 伦理电影大哥的女人| 午夜福利在线观看免费完整高清在| 久久久久久九九精品二区国产| 日韩av不卡免费在线播放| 精品人妻视频免费看| 亚洲欧洲国产日韩| 久久久久久久久大av| 国产成人免费观看mmmm| 亚洲自偷自拍三级| 舔av片在线| 亚洲aⅴ乱码一区二区在线播放| 国产成人a∨麻豆精品| 99久国产av精品国产电影| 麻豆av噜噜一区二区三区| 欧美成人精品欧美一级黄| 亚洲精品乱码久久久v下载方式| 日韩电影二区| 麻豆成人av视频| 久久久久久久久久久免费av| 国内精品美女久久久久久| av国产免费在线观看| 丰满人妻一区二区三区视频av| 麻豆国产97在线/欧美| 久久久精品欧美日韩精品| 亚洲aⅴ乱码一区二区在线播放| 春色校园在线视频观看| 成人鲁丝片一二三区免费| 亚洲不卡免费看| 欧美一级a爱片免费观看看| 亚洲性久久影院| 精品国产露脸久久av麻豆 | 2021天堂中文幕一二区在线观| 美女国产视频在线观看| 国产国拍精品亚洲av在线观看| 久久久精品94久久精品| 综合色av麻豆| 一个人观看的视频www高清免费观看| 精品一区二区三区视频在线| 内地一区二区视频在线| 免费人成在线观看视频色| 国产亚洲精品av在线| 插逼视频在线观看| 成人毛片a级毛片在线播放| 亚洲精品日韩av片在线观看| 97超视频在线观看视频| 免费大片18禁| 一级毛片我不卡| 建设人人有责人人尽责人人享有的 | 免费看美女性在线毛片视频| 亚洲国产最新在线播放| 秋霞在线观看毛片| av线在线观看网站| 在线播放无遮挡| 亚洲成人av在线免费| 亚洲精品亚洲一区二区| 国产成人免费观看mmmm| 久久精品综合一区二区三区| 国产高清有码在线观看视频| 国产人妻一区二区三区在| 日本黄色片子视频| 晚上一个人看的免费电影| 美女被艹到高潮喷水动态| 欧美最新免费一区二区三区| 22中文网久久字幕| 亚洲乱码一区二区免费版| 国产男人的电影天堂91| 久久午夜福利片| 97精品久久久久久久久久精品| kizo精华| 亚洲欧美成人精品一区二区| h日本视频在线播放| 免费观看在线日韩| 亚洲成人精品中文字幕电影| 久久久精品94久久精品| 亚洲精品影视一区二区三区av| 午夜精品在线福利| av线在线观看网站| 精品少妇黑人巨大在线播放| av.在线天堂| 日日啪夜夜撸| 九色成人免费人妻av| 最后的刺客免费高清国语| 国产乱来视频区| 国产综合懂色| 国产探花极品一区二区| 日韩av不卡免费在线播放| 一本—道久久a久久精品蜜桃钙片 精品乱码久久久久久99久播 | 一个人看视频在线观看www免费| 亚洲av日韩在线播放| 免费观看a级毛片全部| 色哟哟·www| 97热精品久久久久久| 免费看av在线观看网站| 可以在线观看毛片的网站| 伦理电影大哥的女人| 亚洲真实伦在线观看| 国产成人福利小说| 精品久久久久久久末码| 亚洲精品成人av观看孕妇| 一级爰片在线观看| 色综合色国产| 2022亚洲国产成人精品| 亚州av有码| 天天躁日日操中文字幕| 日本黄色片子视频| 日韩视频在线欧美| 欧美最新免费一区二区三区| 国产精品无大码| 最近视频中文字幕2019在线8| 亚洲欧美中文字幕日韩二区| videos熟女内射| 国产午夜精品论理片| 中文字幕av成人在线电影| 亚洲精品色激情综合| 毛片一级片免费看久久久久| 国产亚洲一区二区精品| kizo精华| 永久免费av网站大全| 美女黄网站色视频| 欧美成人a在线观看| 蜜桃亚洲精品一区二区三区| 亚洲av中文字字幕乱码综合| 人人妻人人澡人人爽人人夜夜 | 色哟哟·www| 日韩亚洲欧美综合| 又黄又爽又刺激的免费视频.| 欧美最新免费一区二区三区| 欧美bdsm另类| av在线播放精品| 日本黄色片子视频| 久久亚洲国产成人精品v| 国产黄色免费在线视频| 久久久精品94久久精品| 欧美xxxx性猛交bbbb| 久久这里有精品视频免费| 天天躁夜夜躁狠狠久久av| 国产精品人妻久久久影院| 免费观看的影片在线观看| 乱码一卡2卡4卡精品| av福利片在线观看| 五月伊人婷婷丁香| 国产黄色小视频在线观看| 又爽又黄无遮挡网站| 亚洲精华国产精华液的使用体验| 男女边摸边吃奶| 久久久久久国产a免费观看| 亚洲欧美一区二区三区国产| 成人午夜精彩视频在线观看| 国内少妇人妻偷人精品xxx网站| 2018国产大陆天天弄谢| 国产淫片久久久久久久久| 免费少妇av软件| 亚洲电影在线观看av| 又粗又硬又长又爽又黄的视频| 国内精品宾馆在线| 国产综合懂色| 亚洲欧美精品自产自拍| eeuss影院久久| 亚洲aⅴ乱码一区二区在线播放| 亚洲国产成人一精品久久久| 不卡视频在线观看欧美| 国产不卡一卡二| 好男人视频免费观看在线| 国产 一区 欧美 日韩| 熟女电影av网| 国产 亚洲一区二区三区 | 午夜亚洲福利在线播放| 久久草成人影院| 日韩,欧美,国产一区二区三区| 久久久久精品久久久久真实原创| 美女内射精品一级片tv| 国产在线一区二区三区精| 精品不卡国产一区二区三区| 美女黄网站色视频| 婷婷色综合www| 亚洲av成人精品一区久久| ponron亚洲| 成人特级av手机在线观看| 午夜激情欧美在线| 亚洲av中文字字幕乱码综合| 少妇高潮的动态图| 91精品一卡2卡3卡4卡| 在现免费观看毛片| 99久久精品热视频| 黄色日韩在线| 精品久久久噜噜| 亚洲高清免费不卡视频| 建设人人有责人人尽责人人享有的 | 一区二区三区乱码不卡18| 老女人水多毛片| 日韩一本色道免费dvd| 久久这里只有精品中国| av在线播放精品| 国产精品蜜桃在线观看| 欧美最新免费一区二区三区| 免费不卡的大黄色大毛片视频在线观看 | 欧美3d第一页| 成年版毛片免费区| 国产成人免费观看mmmm| 一级毛片我不卡| 国产探花在线观看一区二区| 国产亚洲5aaaaa淫片| 偷拍熟女少妇极品色| 特大巨黑吊av在线直播| 免费av不卡在线播放| 亚洲av男天堂| 一本一本综合久久| 天堂av国产一区二区熟女人妻| 国产黄色免费在线视频| 少妇裸体淫交视频免费看高清| 蜜桃亚洲精品一区二区三区| 成年版毛片免费区| 国产精品久久久久久精品电影| 五月天丁香电影| 国产午夜精品一二区理论片| 婷婷色综合大香蕉| 亚洲va在线va天堂va国产| 国产欧美日韩精品一区二区| 日韩av在线免费看完整版不卡| 日本av手机在线免费观看| 国产久久久一区二区三区| 欧美xxⅹ黑人| 狂野欧美激情性xxxx在线观看| 一级av片app| 成人亚洲欧美一区二区av| 五月玫瑰六月丁香| 一级黄片播放器| 黄色配什么色好看| 亚洲一级一片aⅴ在线观看| 日韩亚洲欧美综合| 亚洲性久久影院| 国产精品久久久久久精品电影| 禁无遮挡网站| 国产午夜精品一二区理论片| 亚洲国产精品sss在线观看| 欧美成人a在线观看| 精品一区二区三卡| 婷婷色麻豆天堂久久| 一二三四中文在线观看免费高清| 老女人水多毛片| 国产亚洲5aaaaa淫片| 亚洲aⅴ乱码一区二区在线播放| 免费黄频网站在线观看国产| 午夜福利高清视频| 青春草亚洲视频在线观看| 欧美极品一区二区三区四区| 精品国产露脸久久av麻豆 | 好男人在线观看高清免费视频| 午夜视频国产福利| 蜜臀久久99精品久久宅男| 大片免费播放器 马上看| 99热这里只有是精品50| 神马国产精品三级电影在线观看| 日日摸夜夜添夜夜爱| 青春草视频在线免费观看| 亚洲真实伦在线观看| av国产久精品久网站免费入址| 联通29元200g的流量卡| 国产一区二区三区av在线| 美女内射精品一级片tv| 国产探花在线观看一区二区| 亚洲乱码一区二区免费版| 日韩欧美精品免费久久| 小蜜桃在线观看免费完整版高清| 3wmmmm亚洲av在线观看| 久久久久精品性色| www.av在线官网国产| 亚洲18禁久久av| 蜜桃亚洲精品一区二区三区| 日本黄大片高清| 在线观看美女被高潮喷水网站| 男女那种视频在线观看| 国产黄频视频在线观看| 一本久久精品| 高清视频免费观看一区二区 | 亚洲精品456在线播放app| 精品不卡国产一区二区三区| 国产精品爽爽va在线观看网站| av线在线观看网站| 婷婷色麻豆天堂久久| 大香蕉97超碰在线| 黄色配什么色好看| 五月玫瑰六月丁香| 国产v大片淫在线免费观看| 三级国产精品欧美在线观看| 亚洲精品久久午夜乱码| 精品一区二区三卡| 国产伦精品一区二区三区视频9| 最近2019中文字幕mv第一页| 亚洲av二区三区四区| 国产亚洲午夜精品一区二区久久 | 久久久久久久国产电影| 久久热精品热| 午夜福利在线观看吧| 伦精品一区二区三区| 亚洲自偷自拍三级| 99久国产av精品| 国产爱豆传媒在线观看| 69av精品久久久久久| 波多野结衣巨乳人妻| 高清日韩中文字幕在线| 一级av片app| 男的添女的下面高潮视频| 两个人视频免费观看高清| 欧美潮喷喷水| 在线观看一区二区三区| 狂野欧美白嫩少妇大欣赏| 久久久久久久久久人人人人人人| 97在线视频观看| 国产综合懂色| av国产免费在线观看| 男人爽女人下面视频在线观看| 蜜桃久久精品国产亚洲av| av女优亚洲男人天堂| 免费大片黄手机在线观看| 亚洲精品国产成人久久av| 日本一本二区三区精品| 日本色播在线视频| 日韩 亚洲 欧美在线| 欧美区成人在线视频| 最近中文字幕高清免费大全6| 秋霞伦理黄片| 美女黄网站色视频| 国产精品伦人一区二区| 全区人妻精品视频| 搡老妇女老女人老熟妇| 一级毛片aaaaaa免费看小| 国产伦精品一区二区三区视频9| 国产乱来视频区| 边亲边吃奶的免费视频| 亚洲欧美清纯卡通| 99re6热这里在线精品视频| 久久久久久久午夜电影| 天堂影院成人在线观看| 国产极品天堂在线| 国内精品宾馆在线| .国产精品久久| 69av精品久久久久久| 亚洲国产日韩欧美精品在线观看| 国产国拍精品亚洲av在线观看| 国产久久久一区二区三区| 精品人妻视频免费看| 人妻少妇偷人精品九色| 亚洲成人av在线免费| 91在线精品国自产拍蜜月| 中文字幕av在线有码专区| 亚州av有码| 亚洲不卡免费看| 国产真实伦视频高清在线观看| 久久久色成人| 国产白丝娇喘喷水9色精品| 观看美女的网站| 欧美变态另类bdsm刘玥| 97超视频在线观看视频| 熟妇人妻不卡中文字幕| 婷婷色麻豆天堂久久| 久久久色成人| 大陆偷拍与自拍| 国产免费视频播放在线视频 | 国产精品国产三级国产av玫瑰| 亚洲精品乱码久久久v下载方式| 亚洲成人精品中文字幕电影| 狂野欧美白嫩少妇大欣赏| 女人久久www免费人成看片| 搡老乐熟女国产| 联通29元200g的流量卡| 久久久久久久久大av| 中文字幕亚洲精品专区| 亚洲精品久久久久久婷婷小说| 国产精品美女特级片免费视频播放器| 日本免费在线观看一区| 少妇猛男粗大的猛烈进出视频 | 国产亚洲精品av在线| 日韩av免费高清视频| 夜夜爽夜夜爽视频| 丰满人妻一区二区三区视频av| 99热这里只有精品一区| 中文欧美无线码| 久久人人爽人人片av| 日韩av在线免费看完整版不卡| 在线免费十八禁| 寂寞人妻少妇视频99o| 国产精品久久视频播放| 国产白丝娇喘喷水9色精品| 一个人免费在线观看电影| 国产伦理片在线播放av一区| 久久久久久久午夜电影| 一级毛片aaaaaa免费看小| 亚洲av成人av| 国产一区有黄有色的免费视频 | 天美传媒精品一区二区| 中国美白少妇内射xxxbb| 亚洲精品国产成人久久av| 日韩一区二区视频免费看| 久久久国产一区二区| 亚洲三级黄色毛片| 亚洲自偷自拍三级| 国产高清三级在线| 菩萨蛮人人尽说江南好唐韦庄| 国产精品蜜桃在线观看| 亚洲精品,欧美精品| 韩国高清视频一区二区三区| 精品人妻偷拍中文字幕| 51国产日韩欧美| 日韩精品青青久久久久久| 国产麻豆成人av免费视频| 国产永久视频网站| 国产黄色视频一区二区在线观看| 永久免费av网站大全| 国产淫片久久久久久久久| 高清午夜精品一区二区三区| 亚洲国产精品成人综合色| 大香蕉久久网| 久久韩国三级中文字幕| 国内精品一区二区在线观看| 免费观看无遮挡的男女| 成人高潮视频无遮挡免费网站| 嫩草影院新地址| 黄色一级大片看看| 久久6这里有精品| 久久鲁丝午夜福利片| 国精品久久久久久国模美| 岛国毛片在线播放| 女人十人毛片免费观看3o分钟| 精品少妇黑人巨大在线播放| 免费看美女性在线毛片视频| 久久久久久久久久成人| 噜噜噜噜噜久久久久久91| 中文资源天堂在线| 青春草亚洲视频在线观看| 中文资源天堂在线| 欧美一区二区亚洲| 国产伦理片在线播放av一区| 嫩草影院精品99| 91精品一卡2卡3卡4卡| 舔av片在线| 麻豆成人av视频| 一级毛片久久久久久久久女| 欧美+日韩+精品| 久久6这里有精品| 黄色日韩在线| 最后的刺客免费高清国语| 夫妻性生交免费视频一级片| 午夜老司机福利剧场| 插逼视频在线观看| 麻豆成人午夜福利视频| 日韩av不卡免费在线播放| 乱系列少妇在线播放| 亚洲精品久久久久久婷婷小说| 99热全是精品| 国产一区二区三区av在线| 中文字幕制服av| 搡女人真爽免费视频火全软件| 亚洲天堂国产精品一区在线| 日韩大片免费观看网站| 夜夜爽夜夜爽视频| 亚洲国产欧美人成| 少妇的逼水好多| 免费大片18禁| 九草在线视频观看| 国产成人a∨麻豆精品| 成人亚洲欧美一区二区av| 国产亚洲av嫩草精品影院| 永久免费av网站大全| 久久久久性生活片| 日本色播在线视频| 女的被弄到高潮叫床怎么办| eeuss影院久久| 99久久人妻综合|