亚洲免费av电影一区二区三区,日韩爱爱视频,51精品视频一区二区三区,91视频爱爱,日韩欧美在线播放视频,中文字幕少妇AV,亚洲电影中文字幕,久久久久亚洲av成人网址,久久综合视频网站,国产在线不卡免费播放

        ?

        共晶和共混對(duì)CL-20/DNB感度和熱解機(jī)理影響的MD模擬

        2017-05-11 11:05:43付一政康志鵬郭志婧苗瑞珍孟瑞鴻楊潞霞劉亞青
        含能材料 2017年2期
        關(guān)鍵詞:結(jié)合能感度鍵長(zhǎng)

        付一政, 康志鵬, 郭志婧, 苗瑞珍, 孟瑞鴻, 楊潞霞, 劉亞青

        (1. 中北大學(xué)材料科學(xué)與工程學(xué)院, 山西 太原 030051; 2. 中北大學(xué)化工與環(huán)境學(xué)院, 山西 太原 030051;3. 山西大學(xué)商務(wù)學(xué)院信息學(xué)院, 山西 太原 030031)

        1 引 言

        六硝基六氮雜異伍茲烷[1](CL-20)是目前能量最高的單質(zhì)炸藥之一,具有廣闊的應(yīng)用前景,但高感度及高成本嚴(yán)重制約了它的發(fā)展應(yīng)用。發(fā)展低感度、低成本和高能量炸藥一直以來都是含能材料領(lǐng)域研究者們追求的目標(biāo),但單質(zhì)炸藥很難在能量、安全和價(jià)格上得以平衡[2]。共晶和共混均可以作為含能材料的改性方法,但與共混相比,共晶是一種新的改性技術(shù),其通過將兩種或兩種以上的炸藥在分子層面上通過分子間作用力加以結(jié)合,可克服普通改性方法不能改變炸藥內(nèi)部組成和晶體結(jié)構(gòu)的局限,形成具有獨(dú)特結(jié)構(gòu)的共晶炸藥,賦予CL-20炸藥新的性能[3-5]。含能材料主要通過在極端條件下發(fā)生復(fù)雜的化學(xué)反應(yīng)釋放能量,反應(yīng)速度非???給實(shí)驗(yàn)研究帶來了時(shí)間和空間上的挑戰(zhàn),當(dāng)前的實(shí)驗(yàn)條件很難從分子或原子尺度上提供含能材料的微觀信息。分子動(dòng)力學(xué)(MD)模擬可以描述原子水平和飛秒尺度上的詳細(xì)信息,有助于在原子水平上認(rèn)識(shí)含能材料在極端條件下的性質(zhì)和反應(yīng)機(jī)理,得到其能量釋放規(guī)律[6-8],對(duì)新材料的設(shè)計(jì)、開發(fā)、合成具有重要意義。CL-20/1,3-二硝基苯(DNB)共晶由于兼具CL-20高能和DNB鈍感廉價(jià)的優(yōu)勢(shì),具有優(yōu)異的使用價(jià)值,近年來得到極大的關(guān)注[9]。本課題組已經(jīng)采用反應(yīng)分子動(dòng)力學(xué)(Reactive Molecular Dynamics, RMD)模擬方法對(duì)CL-20/DNB共晶的熱解機(jī)理進(jìn)行了初步研究[10]。為了進(jìn)一步探討共晶和共混對(duì)CL-20/DNB結(jié)構(gòu)和熱解機(jī)理的影響,本研究在已有工作的基礎(chǔ)上首先基于COMPASS[11]力場(chǎng)采用MD模擬方法對(duì)CL-20、DNB、CL-20與DNB物理混合組成的共混物及CL-20/DNB共晶結(jié)構(gòu)進(jìn)行了研究,分析了共晶和共混對(duì)含能材料感度和力學(xué)性能等的影響,然后采用ReaxFF/lg力場(chǎng)[12],通過RMD模擬對(duì)其熱解過程進(jìn)行了研究,分析了共晶和共混對(duì)含能材料熱解機(jī)理的影響規(guī)律。

        2 模型構(gòu)建與模擬細(xì)節(jié)

        2.1 模型構(gòu)建

        CL-20/DNB共晶結(jié)構(gòu)(圖1a)的構(gòu)建見文獻(xiàn)[10],為了比較,同時(shí)構(gòu)建了CL-20的2×2×2(圖1b)、DNB的2×1×4晶體結(jié)構(gòu)(圖1c)和與共晶結(jié)構(gòu)配比相同的CL-20與DNB物理混合組成的共混物(圖1d)具體建模參數(shù)如表1所示。

        a. CL-20/DNB cocrystal

        c. DNB

        b. CL-20

        d. CL-20/DNB mixture

        圖1CL-20/DNB共晶、CL-20、DNB和CL-20/DNB共混物分子結(jié)構(gòu)

        Fig.1The molecular structures of CL-20/DNB cocrystal, CL-20,DNB and CL-20/DNB mixture

        表1模擬系統(tǒng)的晶胞參數(shù)、密度和原子數(shù)

        Table1Crystal parameters, density and number of atoms of the simulation systems

        systemspacegrouplatticeparametersofsupercell/?abcρ/g·cm-3numberofatomsnumberofCL?20numberofDNBCL?20P21/n26.426016.322029.79601.918115232-DNBPbn2126.514014.048015.22401.575512-32CL?20/DNBcocrystalPbca9.470313.458933.62001.88016643232CL?20?DNBmixture-20.493820.493844.44901.72616643232

        利用Materials Studio 6.0軟件[13]在COMPASS力場(chǎng)條件下對(duì)所構(gòu)建的分子模型進(jìn)行結(jié)構(gòu)優(yōu)化,將優(yōu)化后的模型在300 K下先進(jìn)行500 ps的等溫等容(NVT)MD模擬,再進(jìn)行500 ps溫度為300 K,壓力為0 GP的等溫等壓(NPT)MD模擬,目的是為了對(duì)超晶胞內(nèi)部的壓力進(jìn)行弛豫,獲得常溫常壓條件下系統(tǒng)的初始狀態(tài),使其更加接近實(shí)際情況。最后50 ps體系已經(jīng)平衡,用于分析鍵長(zhǎng)、結(jié)合能和力學(xué)性能等。模擬過程中采用Andersen控壓方法[14],Berendsen控溫方法[15],各分子起始速度按Maxwell分布取樣,velocity Verlet算法[16]進(jìn)行求解,范德華(vdW)和靜電(coulomb)作用分別用atom based[17]和Ewald[18]方法計(jì)算,非鍵截取半徑(cutoff distance)取0.95 nm,樣條寬度(spline width)取0.1 nm,緩沖寬度(buffer width)取0.05 nm,時(shí)間步長(zhǎng)1 fs。

        最后運(yùn)用LAMMPS(Large-scale Atomic/Molecular Massively Parallel Simulator)軟件[18]在ReaxFF/lg力場(chǎng)條件下進(jìn)行2000 K下50 ps的NVT分子動(dòng)力學(xué)計(jì)算,采用Berendsen方法對(duì)溫度和壓力進(jìn)行控制,耦合參數(shù)為100 fs,使其在設(shè)定的范圍內(nèi)波動(dòng),采用周期性邊界條件,步長(zhǎng)為0.1 fs,每隔50 fs記錄一次原子軌跡以及鍵級(jí)信息,鍵級(jí)為0.3。

        3 結(jié)果與討論

        3.1 共晶和共混對(duì)引發(fā)鍵鍵長(zhǎng)的影響

        研究表明引發(fā)鍵的最大鍵長(zhǎng)可以關(guān)聯(lián)感度[20-21],可作為感度相對(duì)大小的理論判據(jù),雖然最大鍵長(zhǎng)在MD平衡軌跡文件中出現(xiàn)的概率很小,但是它一旦被引發(fā),最容易引發(fā)初始的熱解和爆轟反應(yīng)。表2給出不同模擬體系中引發(fā)鍵的平均鍵長(zhǎng)(Lave)和最大鍵長(zhǎng)(Lmax)的MD模擬結(jié)果,可見與CL-20和DNB相比,無論是共晶還是共混都會(huì)對(duì)引發(fā)鍵的最大鍵長(zhǎng)產(chǎn)生影響,使得CL-20引發(fā)鍵N—NO2的最大鍵長(zhǎng)減小,感度降低,DNB引發(fā)鍵C—NO2的最大鍵長(zhǎng)增大,感度增加,而且共晶的影響效果比共混的大。

        表2CL-20 引發(fā)鍵(N—NO2)和DNB引發(fā)鍵(C—NO2)的最大鍵長(zhǎng)(Lmax)和平均鍵長(zhǎng)(Lave)

        Table2Maximum length(Lmax)and average length(Lave)of the trigger bond (N—NO2) of CL-20 and (C—NO2)of DNB for different systems

        triggerbondLmax/?Lave/?N—NO2ofCL?201.5431.395N—NO2ofCL?20ofmixture1.5251.395N—NO2ofCL?20ofcocrystal1.5171.393C—NO2ofDNB1.6001.454C—NO2ofDNBofmixture1.6051.455C—NO2ofDNBofcocrystal1.6181.453

        3.2 共晶和共混對(duì)體系結(jié)合能的影響

        結(jié)合能是組分間相互作用強(qiáng)弱的標(biāo)志,結(jié)合能越大,組分間的相互作用越強(qiáng),形成的體系越穩(wěn)定。表3給出共晶和共混模型得到的CL-20與DNB之間的結(jié)合能,及各能量項(xiàng)對(duì)結(jié)合能的貢獻(xiàn)情況,可見共晶結(jié)構(gòu)中CL-20與DNB的結(jié)合能遠(yuǎn)大于共混體系,說明共晶的結(jié)構(gòu)較共混更加穩(wěn)定,同時(shí)還發(fā)現(xiàn)對(duì)結(jié)合能的貢獻(xiàn)主要是由非鍵相互作用引起,其中靜電力所占比例較大,范德華能量項(xiàng)所占比例較小??梢姽簿ㄟ^改變炸藥的內(nèi)部組成和晶體結(jié)構(gòu),使組分間形成強(qiáng)烈的分子間作用力,增加了共晶炸藥分子體系的穩(wěn)定性,提高了共晶分子對(duì)機(jī)械外力的抗振性,從而提高了其安全性。

        表3 體系的結(jié)合能(Ebind)及庫侖相互作用能(Ecoulomb)和范德華能量(EvdW)

        3.3 共晶和共混對(duì)體系力學(xué)性能的影響

        力學(xué)性能是含能材料的最重要性能之一,其關(guān)系到材料的制備、加工和使用,表4給出基于MD模擬的原子運(yùn)動(dòng)軌跡分析所得體系的力學(xué)性能。由表4可見,與CL-20和DNB相比,CL-20/DNB共晶和共混體系的大部分彈性系數(shù)(Cij)、模量均比CL-20的低,比DNB的高,說明共晶和共混均可以使CL-20體系的剛性減小,柔性增強(qiáng),即體系變“軟”,當(dāng)體系受到外力作用時(shí),可以有效緩沖和分散外力作用,減小炸藥顆粒之間的摩擦,使其內(nèi)部應(yīng)力分布更均勻,從而減少“熱點(diǎn)”的形成,減低材料的感度,增加體系的安全性,但共混會(huì)使體系的力學(xué)性能下降過多。

        表4不同體系的力學(xué)性能

        Table4The mechanical properties of different systems

        parameterCL?20DNBcocrystalmixtureC1112.14.38.98.0C2214.74.47.05.1C3316.75.512.96.0C442.51.12.81.2C552.41.19.12.0C666.01.22.71.9C125.12.63.93.4C135.12.99.62.7C236.52.94.42.9(C12?C44)2.61.51.12.2tensilemodulus(E)9.22.52.93.3Poissonratio(v)0.30.370.40.4bulkmodulus(K)8.13.36.33.5shearmodulus(G)3.31.02.51.3K/G2.53.32.52.7

        Note:Cijis an element of the elastic coefficient matrix of 6×6. The unit ofCij,E,KandGis GPa.

        3.4 共晶和共混對(duì)體系熱解勢(shì)能和總物種數(shù)的影響

        圖2a給出了不同體系熱解過程中勢(shì)能隨模擬時(shí)間的變化情況,為了便于比較,縱坐標(biāo)設(shè)為勢(shì)能的變化量。從圖2a中可以看出由于CL-20感度較高,受熱很快發(fā)生熱解,放出大量的熱,勢(shì)能發(fā)生明顯下降,但DNB由于感度較低,幾乎在整個(gè)模擬時(shí)間內(nèi)都沒有發(fā)生明顯變化; 而共晶體系中CL-20與DNB分子間隔排列,由于DNB的阻隔作用使得相鄰CL-20分解的初級(jí)和次級(jí)產(chǎn)物很難再發(fā)生接觸,限制了CL-20的進(jìn)一步快速分解; 而共混體系中CL-20與DNB僅在兩者界面處有少量接觸,DNB的加入并不會(huì)明顯影響CL-20的反應(yīng),所以反應(yīng)初期共混體系中CL-20迅速分解,放出大量熱量,導(dǎo)致較共晶體系勢(shì)能下降的快; 但在反應(yīng)后期(25 ps以后),由于共晶中CL-20與DNB接觸的機(jī)會(huì)較共混的多,更容易促進(jìn)DNB的熱解,共晶體系中的DNB開始大量分解,所以勢(shì)能下降趨勢(shì)較共混體系的更加明顯。

        a. potential energy

        b. total species

        圖2體系勢(shì)能及反應(yīng)過程中總物種數(shù)量隨模擬時(shí)間的變化

        Fig.2Change in system potential energy and total species number in the reaction process with the simulation time

        圖2b是模擬反應(yīng)過程中總物種數(shù)量隨時(shí)間的變化情況,由于CL-20和DNB體系中初始原子個(gè)數(shù)較少(分別為1152和512個(gè)原子),所以最終產(chǎn)生的物種數(shù)比共晶和共混體系的少; DNB的感度較低,在前12.5 ps內(nèi)都沒有發(fā)生熱解反應(yīng),最終產(chǎn)生的物種數(shù)也很少,而CL-20的感度較高,反應(yīng)初期(5 ps內(nèi))就有大量的產(chǎn)物生成,隨后在17 ps左右發(fā)生波動(dòng)。將兩者通過機(jī)械共混后,由于低感度DNB的加入使得體系初期的分解速度和產(chǎn)生的物種數(shù)較純CL-20的低,但由于體系的總原子個(gè)數(shù)較純CL-20的多,后期產(chǎn)生的總物種數(shù)快速增加。而對(duì)于共晶體系,反應(yīng)開始分解產(chǎn)物的物種數(shù)就比CL-20的多,2.15 ps后可能要引發(fā)DNB的分解,使得反應(yīng)速率下降,產(chǎn)生的物種數(shù)較CL-20的少,但20 ps后伴隨DNB的大量分解,產(chǎn)物的總物種數(shù)迅速增加,而且高于共混體系。

        3.5 共晶和共混對(duì)體系主要產(chǎn)物分布的影響

        a. CL-20/DNB cocrystal

        c. CL-20

        d. DNB

        圖3CL-20/DNB共晶, CL-20/DNB共混, CL-20和DNB主要產(chǎn)物分布曲線隨模擬時(shí)間的變化

        Fig.3Change in distribution curves of main products of CL-20/DNB cocrystal, CL-20/DNB mixture, CL-20 and DNB with the simulation time

        圖3為不同體系熱解的主要產(chǎn)物,可以看到CL-20的分解速度很快(圖3a,圖3b和圖3c),幾乎在0.05 ps內(nèi)就分解完畢,但DNB的分解速度卻非常慢(圖3d),在20 ps內(nèi)都未見明顯變化,即使到50 ps仍有大量的殘余。4個(gè)體系的熱解產(chǎn)物主要有NO2、N2、NO、H2O、HONO、HON以及CO2等,其中最早生成的產(chǎn)物均為NO2,這主要是由于CL-20的熱解反應(yīng)中首先發(fā)生的是N—NO2鍵的斷裂,這是環(huán)狀硝胺類化合物熱解的特征引發(fā)反應(yīng),而DNB熱解的引發(fā)鍵為C—NO2鍵,它是硝基芳香化合物中最弱的鍵,這兩個(gè)鍵斷裂后的產(chǎn)物均為NO2,所以降解初期主要產(chǎn)物為NO2,由于C—NO2鍵斷裂所需的能量比N—NO2的多,所以DNB的分解速率比CL-20慢。隨后次級(jí)反應(yīng)會(huì)消耗掉大量的NO2,使其數(shù)量急劇減少,生成NO、HONO和HNO3等物質(zhì)。NO也是模擬體系中一種重要的產(chǎn)物,其變化情況與NO2相似,前期均快速生成,隨后發(fā)生次級(jí)反應(yīng),最終生成大量的N2。

        4 結(jié) 論

        在COMPASS力場(chǎng)和ReaxFF/lg力場(chǎng)條件下采用分子動(dòng)力學(xué)(MD)模擬方法對(duì)CL-20、DNB、CL-20與DNB組成的共晶和共混物的結(jié)構(gòu)進(jìn)行了研究,分析了材料的感度、力學(xué)性能和熱解過程,得到如下結(jié)論:

        (1) 共晶和共混會(huì)使CL-20引發(fā)鍵N—NO2最大鍵長(zhǎng)減小,感度降低,DNB引發(fā)鍵C—NO2最大鍵長(zhǎng)增大,感度增加,但共晶的效果更明顯。共晶結(jié)構(gòu)中CL-20與DNB的結(jié)合能遠(yuǎn)大于共混體系,說明共晶較共混結(jié)構(gòu)更加穩(wěn)定。共晶和共混均可以使得CL-20體系的剛性減小,柔性增強(qiáng),增加體系的安全性,但共混會(huì)使體系的力學(xué)性能下降過多。

        (2) 熱解過程初期由于共晶中DNB對(duì)CL-20次級(jí)反應(yīng)的阻隔作用使得共混體系勢(shì)能下降的趨勢(shì)較共混體系的慢,但后期共晶會(huì)促進(jìn)DNB的分解,使得體系勢(shì)能下降較共混體系明顯。共晶和共混體系的主要熱解產(chǎn)物主要有NO2、N2、NO、H2O、HONO、HON以及CO2等。綜上所述,相同配比條件下共晶對(duì)CL-20的改性效果比共混的明顯。

        參考文獻(xiàn):

        [1] Nielsen A T, Chafin A P, Christian S L, et al. Synthesis of polyazapolycyclic caged polynitramines[J].Tetrahedron, 1998, 54(39): 11793-11812.

        [2] Sikder A K, Sikder N. A review of advanced high performance, insensitive and thermally stable energetic materials emerging for military and space applications[J].JournalofHazardousMaterials, 2004, 112(1): 1-15.

        [3] Bolton O, Matzger A J. Improved stability and smart-material functionality realized in an energetic cocrystal[J].AngewandteChemieInternationalEdition, 2011, 50(38): 8960-8963.

        [4] Bolton O, Simke L R, Pagoria P F, et al. High power explosive with good sensitivity: A 2:1 cocrystal of CL-20: HMX[J].CrystalGrowth&Design, 2012, 12(9): 4311-4314.

        [5] Guo C, Zhang H, Wang X, et al. Crystal structure and explosive performance of a new CL-20/caprolactam cocrystal[J].JournalofMolecularStructure, 2013, 1048: 267-273.

        [6] Van Duin A C T, Dasgupta S, Lorant F, et al. ReaxFF: a reactive force field for hydrocarbons[J].TheJournalofPhysicalChemistryA, 2001, 105(41): 9396-9409.

        [7] Rom N, Zybin S V, van Duin A C T, et al. Density-dependent liquid nitromethane decomposition: molecular dynamics simulations based on ReaxFF[J].TheJournalofPhysicalChemistryA, 2011, 115(36): 10181-10202.

        [8] Han S, van Duin A C T, Goddard III W A, et al. Thermal decomposition of condensed-phase nitromethane from molecular dynamics from ReaxFF reactive dynamics[J].TheJournalofPhysicalChemistryB, 2011, 115(20): 6534-6540.

        [9] Wang Y, Yang Z, Li H, et al. A novel cocrystal explosive of HNIW with good comprehensive properties[J].Propellants,Explosives,Pyrotechnics, 2014, 39(4): 590-596.

        [10] 苗瑞珍, 劉偉帥, 王建, 等. CL-20/DNB共晶高溫?zé)峤獾腞eaxFF反應(yīng)分子動(dòng)力學(xué)模擬[J]. 含能材料, 2016, 24(2): 111-117.

        MIAO Rui-zhen, LIU Wei-shuai, WANG Jian, et al. ReaxFF reactive molecular dynamics simulation of thermal decomposition under high temperature for CL-20/DNB cocrystal[J].ChneseJounralofEnergeticMaterials(HannengCailiao), 2016, 24(2): 111-117.

        [11] Sun H. COMPASS: an ab initio force-field optimized for condensed-phase applications overview with details on alkane and benzene compounds[J].TheJournalofPhysicalChemistryB, 1998, 102(38): 7338-7364.

        [12] Liu L, Liu Y, Zybin S V, et al. ReaxFF-lg: Correction of the ReaxFF reactive force field for London dispersion,with applications to the equations of state for energetic materials[J].TheJournalofPhysicalChemistryA, 2011, 115(40): 11016-11022.

        [13] Materials Studio 6.0[CP]. Accelrys Software Inc, USA: 2013.

        [14] Andersen H C. Molecular dynamics simulations at constant pressure and/or temperature[J].TheJournalofChemicalPhysics, 1980, 72(4): 2384-2393.

        [15] Wiberg K B, Wasserman D J, Martin E. Enthalpies of hydration of alkenes. 2. The n-heptenes and n-pentenes[J].TheJournalofPhysicalChemistry, 1984, 88(16): 3684-3688.

        [16] Swope W C, Andersen H C, Berens P H, et al. A computer simulation method for the calculation of equilibrium constants for the formation of physical clusters of molecules: Application to small water clusters[J].TheJournalofChemicalPhysics, 1982, 76(1): 637-649.

        [17] Tosi M P. Cohesion of ionic solids in the Born model[J].SolidStatePhysics, 1964, 16: 1-120.

        [18] Ewald P P. Die Berechnung optischer and elektrostatischer Gitterpotentiale[J].AnnalenderPhysik, 1921, 369(3): 253-287.

        [19] Plimpton S J, Large-scale atomic/molecular massively parallel simulator[CP]. Sandia National Laboratories: 2003.

        [20] 肖繼軍, 朱衛(wèi)華, 朱偉, 等. 高能材料分子動(dòng)力學(xué)[M]. 北京: 科學(xué)出版社, 2012: 10-34.

        XIAO Ji-jun, ZHU Wei-hua, ZHU Wei, et al. Molecular Dynamics of Energetic Materials[M]. Beijing: Science Press, 2012: 10-34.

        [21] 朱偉, 肖繼軍, 鄭劍, 等. 高能混合物的感度理論判別——不同配比和不同溫度AP/HMX的MD研究[J]. 化學(xué)學(xué)報(bào), 2008, 66(23): 2592-2596.

        ZHU Wei, XIAO Ji-jun, ZHENG Jan, et al. A theoretical criterion for sensitivity of energetic composites—molecular dynamics studies on AP/HMX systems at various concentrations and temperatures[J].ActaChimSinica, 2008, 66(23): 2592-2596.

        猜你喜歡
        結(jié)合能感度鍵長(zhǎng)
        基于機(jī)器學(xué)習(xí)的RDX-CMDB推進(jìn)劑安全性能預(yù)測(cè)
        陰離子調(diào)控錳基鈣鈦礦中Mn─O的鍵長(zhǎng)和磁性
        晶體結(jié)合能對(duì)晶格動(dòng)力學(xué)性質(zhì)的影響
        借鑒躍遷能級(jí)圖示助力比結(jié)合能理解*
        密度泛函理論研究鎘的二鹵化合物分子的結(jié)構(gòu)和振動(dòng)頻率
        淺議鍵能與鍵長(zhǎng)的關(guān)系
        ε-CL-20/F2311 PBXs力學(xué)性能和結(jié)合能的分子動(dòng)力學(xué)模擬
        對(duì)“結(jié)合能、比結(jié)合能、質(zhì)能方程、質(zhì)量虧損”的正確認(rèn)識(shí)
        高感度活性稀釋劑丙烯酰嗎啉的合成研究
        FOX-7晶體形貌對(duì)感度的影響
        性色av免费网站| 午夜黄色一区二区不卡| 亚洲精品中文字幕91| 日韩 无码 偷拍 中文字幕| 桃花色综合影院| 中文字幕精品久久天堂一区| 久久久精品2019免费观看| 91华人在线| 一道本加勒比在线观看| 精品亚洲成a人在线观看| 日本阿v网站在线观看中文| 人人妻人人玩人人澡人人爽| 亚洲精品国产二区三区在线| 亚洲a级视频在线播放| 中文字幕在线日亚州9| 竹菊影视欧美日韩一区二区三区四区五区 | 国产极品美女高潮无套| 日本亚洲色大成网站www久久| 毛片无遮挡高清免费久久| 亚洲av手机在线一区| 久久aaaa片一区二区| 台湾佬娱乐中文22vvvv| 国产91在线|亚洲| 午夜视频一区二区三区播放 | 国产成人精品日本亚洲| 亚洲欧洲日产国码无码| 视频女同久久久一区二区| 疯狂添女人下部视频免费| 91av精品视频| 18禁成人免费av大片一区| 欧美性白人极品1819hd| 黑人巨大白妞出浆| 亚洲午夜无码久久久久软件| 在线观看国产激情视频| 亚洲av无码一区二区三区不卡| 在线人妻无码一区二区| 按摩偷拍一区二区三区| 男女交射视频免费观看网站| 欧美黑人性暴力猛交喷水黑人巨大| 亚洲欧美日韩国产一区二区精品| 亚洲一区二区日韩精品|