999精品在线视频,手机成人午夜在线视频,久久不卡国产精品无码,中日无码在线观看,成人av手机在线观看,日韩精品亚洲一区中文字幕,亚洲av无码人妻,四虎国产在线观看 ?

基于團簇結構模型的Na2O-A l2O3-SiO2三元體系熔體吉布斯混合摩爾自由能計算

2016-09-13 03:10:16李永康尤靜林王建王敏馬楠魏廣超
物理化學學報 2016年3期
關鍵詞:體系結構模型

李永康 尤靜林王建王 敏馬楠 魏廣超

(省部共建高品質特殊鋼冶金與制備國家重點實驗室,上海市鋼鐵冶金新技術開發應用重點實驗室,上海大學材料科學與工程學院,上海200072)

基于團簇結構模型的Na2O-A l2O3-SiO2三元體系熔體吉布斯混合摩爾自由能計算

李永康尤靜林*王建王敏馬楠魏廣超

(省部共建高品質特殊鋼冶金與制備國家重點實驗室,上海市鋼鐵冶金新技術開發應用重點實驗室,上海大學材料科學與工程學院,上海200072)

以三元硅酸鹽熔體團簇結構模型為基礎,選取了Na2O-Al2O3-SiO2體系不同成分的團簇結構,采用半經驗量子化學方法MNDO/d分別計算該三元體系熔體中不同結構的團簇基元在1473、1873、2000K溫度下的熵、焓、熱容和自由能等熱力學數據,計算得出不同團簇結構基元的混合自由能,并根據統計熱力學波爾茲曼分布定律,推導計算得出Na2O-Al2O3-SiO2三元體系各成分下的混合摩爾自由能。三元硅酸鹽熔體的熱力學性質與該熔體的微觀結構密切相關。

Na2O-Al2O3-SiO2;團簇結構模型;熔體;混合摩爾自由能

從Neuville10,Yang11和劉欽12等的研究中都可以得知鋁硅酸玻璃中鋁的配位數是四、五和六,四配位鋁占多數。長期以來人們一度認為玻璃繼承了熔體的結構,因而很多學者憑借對玻璃的研究來說明熔體結構。可是,玻璃的結構不能等同于熔體的結構。劉欽等12通過高溫拉曼光譜對鋁硅酸鹽的研究認為鋁在硅酸鹽熔體中可以發現四配位和六配位的結構而且主要的配位數是四,但并沒有證據證明有五配位鋁的結構存在于熔體結構中。說明雖然很多研究指出鋁硅酸鹽玻璃中鋁的配位數是四、五和六,但是鋁硅酸鹽熔體中只有含四配位鋁和六配位鋁的結構。類似的例子可以從鍺酸鹽和硼酸鹽的研究中發現。鍺酸鹽玻璃和晶體中鍺的配位數為四和六,但是熔體中只有四配位鍺13;硼酸鹽玻璃和晶體中硼的配位數為三和四,但是熔體中只有三配位硼14,15。

You等16,17利用高溫拉曼光譜和量子化學從頭計算對硅酸鹽晶體、玻璃和熔體結構進行了研究,并預測了硅酸鹽熔體中可能存在的團簇結構基元。除此之外,王威和尤靜林18基于團簇結構模型采用半經驗方法對堿金屬二元硅酸鹽的熱力學性質進行了計算,得出100-2000K時的熵、焓、吉布斯自由能及比定壓熱容,并由此計算出該二元體系的混合自由能,并與實驗數據相吻合。從硅酸鹽微觀精細結構出發計算得到某成分下硅酸鹽宏觀熱力學性質是新的研究方法,具有指導意義。因此本文基于Na2O-A l2O3-SiO2三元體系熔體團簇結構,采用量子化學中的半經驗方法,進行三元體系的高溫熱力學性質計算,給出了該體系的熱力學性質隨溫度的變化規律。

2 實驗部分

2.1鋁硅酸鹽熔體團簇結構模型

團簇的選取既要考慮相關晶體的結構,同時也包含了熔體和玻璃體系特有的結構19,20。尤靜林等12,17,21采用高溫拉曼光譜技術探測了硅酸鹽和鋁硅酸鹽類熔體,并將熔體里的局域結構與拉曼光譜位移建立了對應的關系,并用量子化學從頭計算進行相關有效的驗證,從而作為熔體團簇結構假設的必要基礎條件。在此基礎上得到如圖1所示的Na2O-A l2O3-SiO2三元體系熔體團簇結構。

2.2計算方法

本文使用Materials Studio 5.5軟件22搭建不同成分的Na2O-Al2O3-SiO2三元系團簇模型,運用半經驗的VAMP(Vienna Ab initio Molecular-dynam ics Package)模塊對其進行結構優化和熱力學計算。計算參數設置為:哈密頓函數選擇MNDO/d23,24,收斂判據設為16.747 kJ?mol-1?nm-1,算法為Hartree-Fock(RHF)自洽分子軌道理論方法,自洽場(SCF)誤差選擇1.0×10-5,溫度設置為100-2000K,步長為25K,計算結果為熵、焓、熱容在不同溫度下的值。計算結果將用于計算不同溫度下該三元體系的混合自由能。

3 計算結果與討論

3.1Na2O-Al2O3-SiO2三元系團簇基元混合自由能

Na2O-A l2O3-SiO2三元系團簇基元的熵、熱容和焓值可由計算結果文件中直接給出。團簇基元的自由能G(T)由下式來計算:

根據VAMP模塊數據處理的需要,若要求得各團簇基元的自由能G(T),需要對其進行以下修正:

其中:G0(298 K)為計算時輸出文件中的生成熱(heatof formation)的值。

圖1  Na2O-A l2O3-SiO2三元系團簇基元的構型圖Fig.1 Diagram sofm odel clustersof the Na2O-A l2O3-SiO2ternary system Element labelsare all the sameasin the diagram of AlSi2O8Na5(a).

由半經驗方法的計算結果得到不同溫度下同一團簇基元的熵、焓和熱容,再通過公式(1)和(2)得出自由能。A lSi2O8Na5團簇基元的熵、焓、熱容、自由能與溫度的關系圖見圖2。結構不同的硅酸鹽團簇,它們的熵、焓、熱容和自由能隨溫度的變化規律是相同的,其中熵、焓和熱容隨著溫度升高而有不同程度的增加,自由能隨著溫度的升高而減小。

圖2  AlSi2O8Na5團簇基元(圖1(a))的熵、焓、熱容、自由能與溫度的關系圖Fig.2 Tem peratue dependententropy,enthalpy,heat capacity,free energy of A lSi2O8Na5cluster in Fig.1(a)

表1 歸一化后的Li2O、Na2O和SiO2團簇的體系自由能Tab le1 Normalized freeenergy of Li2O,Na2O,and SiO2clusters

為計算Na2O-A l2O3-SiO2三元系團簇基元在高溫下某一溫度的混合自由能,從100到2000K中選取了硅酸鹽熔體中常用的可用于比較的三個溫度點,分別是1475、1875和2000K。為了確定純物質Na2O、Al2O3和SiO2在1473、1873和2000K的自由能,從Na2O、SiO2和A l2O3的超晶胞中提取不同原子數的Na2O、SiO2和A l2O3團簇基元,結構優化后用半經驗方法MNDO/d計算并將團簇基元的自由能歸一化為單位摩爾分數的Na2O、SiO2和A l2O3,歸一化的自由能分別用、和表示。選擇Na2O、A l2O3和SiO2團簇基元自由能的相對低點分別記為Na2O、A l2O3和SiO2的穩定能級。并假設團簇基元的能級分布滿足準連續近似,根據統計熱力學波爾茲曼分布定律,通過式(3)將團簇基元的自由能轉化為體系的平均自由能。

上式中ΔGmin和ΔGsys分別為團簇基元和體系的摩爾自由能,計算結果見表1。從表1中可以發現,純物質Na2O、A l2O3和SiO2體系混合自由能的計算值要比實驗值25普遍偏高,而且同一純物質不同溫度的計算值與實驗值之間的差值幾乎一樣。說明目前采用VAMP模塊計算物質的絕對能量存在系統誤差。熔體的混合自由能是一個相對能量值,公式(4)中的各個項都是用同樣的計算方法計算得出,系統誤差在計算過程中會被抵消,對最后結果影響有限。

若以m Na2O?n A l2O3?p SiO2統一表示Na2OA l2O3-SiO2三元系團簇基元的組成(其中m、n和p分別為Na2O、A l2O3和SiO2的摩爾分數),則Na2OA l2O3-SiO2團簇基元的混合摩爾自由能可通過公式(4)歸一化求得。Na2O-Al2O3-SiO2團簇基元的混合摩爾自由能計算如圖3所示。

本工作是基于團簇模型作為熔體結構的基本假設。從Na2O-A l2O3-SiO2三元體系相圖26(見圖3)可知,編號為(g)、(h)、(j)、(l)和(m)的團簇模型所在成分點的材料在1473 K下并未熔解,因此會表現較大差異。但這五個團簇結構的成分點材料在1873和2000K高溫下皆處于熔體狀態,因此它們適用于1873和2000K下的模擬計算。

圖3  Na2O-A l2O3-SiO2三元相圖26Fig.3 Phasc diagram of Na2O-A l2O3-SiO226

團簇基元混合自由能的計算值如圖4、圖5和圖6所示呈散點分布。通過將各成分的最低點篩選擬合出一個相對光滑且無奇點的光滑曲面(如圖4中實線網格曲面),并稱該曲面為團簇基元最低能量前沿(CULEF),這意味著該體系中不會存在一個團簇基元的混合自由能低于該曲面。處在CULEF曲面附近的團簇基元混合自由能較低,這些團簇基元在熔體中將會較穩定存在并且其分布權重也大。又由計算過程可知,由于團簇基元的混合自由能是通過差減得到的相對能量,有效地降低了因計算方法而帶來的系統誤差,因而更具有實際參考價值。

圖4  Na2O-A l2O3-SiO2團簇基元在1473 K的混合自由能和其最低能量前沿Fig.4Mixing freeenergy of Na2O-A l2O3-SiO2cluster and cluster unit low est energy frontier CULEF at1473 K

從以上內容可知,團簇基元(f)、(g)、(h)、 (j)、(k)、(l)、(m)、(r)是最接近最低能量前沿CULEF曲面的團簇模型,因此它們是在高溫下各自的成分點上的分布權重較大的團簇結構。而且我們得到了穩定團簇結構的同時,也得到了該三元體系的團簇結構在不同溫度和成分下的最低混合自由能值。

圖5  Na2O-A l2O3-SiO2團簇基元在1873K的混合自由能和其最低能量前沿Fig.5Mixing freeenergy of Na2O-A l2O3-SiO2cluster and cluster unit lowestenergy frontier CULEF at1873K

圖6  Na2O-A l2O3-SiO2團簇基元在2000K的混合自由能和其最低能量前沿Fig.6Mixing freeenergy of Na2O-Al2O3-SiO2cluster and cluster unit lowest energy fron tier CULEF at 2000K

3.2Na2O-Al2O3-SiO2三元體系平均混合自由能

在較高的溫度下,團簇基元的鍵長和鍵角會在一定范圍內變化,這導致團簇基元的混合自由能發生連續變化。由此可知,團簇基元混合自由能的能級分布依然滿足準連續近似,同樣通過(3)式,可由最低能量前沿的混合自由能求得Na2OA l2O3-SiO2三元體系的平均混合自由能(如圖7中的曲面)。Na2O-A l2O3-SiO2三元體系的平均混合自由能與HSC chem istry熱學數據27進行對比,結果見圖7。

圖7  Na2O-Al2O3-SiO2三元系在1473、1873和2000K的平均混合自由能和HSC chem istry熱學數據的比較Fig.7 Com parison of averagem ixing freeenergy of Na2O-A l2O3-SiO2system with HSC chem istry data calculated at 1473,1873,and 2000K

從圖7可以看出平均混合自由能和HSC chem istry熱學數據有一定的誤差,誤差小于10%。1473 K下NaA lSiO4通過本工作計算得出的平均混合摩爾自由能與HSC計算值有較大的誤差。原因是NaA l-SiO4在1473 K下仍未熔解(見圖3),晶體狀態下的NaAlSiO4的平均混合自由能顯然比以計算熔體的方法模擬計算得到的平均混合自由能低。相反,在1473 K處于熔體狀態的其他三個成分點計算得出的結果與HSC計算值吻合較好,并且在1873和2000K下處于熔體狀態的NaA lSiO4的HSC計算值依然與本工作模擬計算得出的平均混合摩爾自由能吻合較好。這種差異恰恰驗證了本工作模擬計算方法的合理性。

4 結論

通過分子模擬半經驗量化計算,利用Materials Studio 5.5計算軟件中半經驗的MNDO/d計算方法計算得出Na2O-A l2O3-SiO2三元系的不同團簇基的熵、焓、熱容和自由能隨溫度變化的能量曲線。進而計算得出Na2O-A l2O3-SiO2三元系不同的團簇基元的混合自由能,從而推導出該成分點在體系中的1473、1873和2000K的平均混合自由能。VAMP模塊的計算方法仍存在系統誤差,但是熔體的混合自由能是一個相對值,公式中的各個項都是用同樣的計算方法計算得出,系統誤差在這樣的計算過程中會被消除掉,對混合自由能值沒有太大的影響。Na2O-A l2O3-SiO2三元體系的平均混合自由能計算結果與HSC chem istry數據庫的熱學數據吻合較好,表明了Na2O-A l2O3-SiO2三元體系的熱力學性能不僅與硅氧四面體的種類有關,而且與硅氧四面體的近鄰連接情況和鋁的配位數有關,即與硅酸鹽的精細結構密切相關。搭建硅酸鹽團簇模型需要預先了解硅酸鹽熔融后的結構形態,掌握其高溫特性。本研究工作表明團簇模型可以描述熔體中的微觀結構及獲取相關的熱力學性能。

Refe ren ces

(1)Ma,X.C.;Yu,C.M.;Cheng,G.G.Journalof Iron and Steel Research 2011,5,3.[馬小春,于春梅,成國光.鋼鐵研究學報,2011,5,3.]doi:10.3969/j.issn.1001-7208.2011.02.010

(2)Zhang,J.JournalofBeijing University of Iron and Steel Techno logy 1988,10,412.[張鑒.北京鋼鐵學院學報, 1988,10,412.]doi:10.13374/j.issn1001-053x.1988.04.023

(3)Lü,N.N.;Yu,J.K.;Su,C.;Wang,H.Z.Journalof Northeastern University Natural Science 2013,34(12),1743. [呂寧寧,于景坤,蘇暢,王洪章.東北大學學報(自然科學版),2013,34(12),1743.]

(4)Mysen,B.O.;Lucier,A.;Cody,G.D.AmericanMineralogist 2003,88(11-12),1668.doi:10.2138/am-2003-11-1206

(5)Lin,P.;Pelton,A.Metallurgical TransactionsB 1979,10(4), 667.doi:10.1007/BF02662569

(6)Masson,C.Journalofthe American Ceramic Society1968,51 (3),134.doi:10.1111/jace.1968.51.issue-3

(7)Ottonello,G.JournalofNon-Crystalline So lids 2001,282(1), 72.doi:10.1016/S0022-3093(01)00330-1

(8)Zaitsev,A.I.;A rutyunyan,N.A.;Shaposhnikov,N.G.; Zaitseva,N.E.;Burtsev,V.T.Russ.J.Phys.Chem.2006,80(3),335.doi:10.1134/S0036024406030058

(9)Glibin,V.P.;King,P.L.Calphad-Comput.Coupling Ph. Diagrams Thermochem.2015,49,19.doi:10.1016/j. calphad.2015.02.001

(10)Neuville,D.R.;Corm ier,L.;Massiot,D.ChemicalGeology 2006,229(1),173.doi:10.1016/j.chemgeo.2006.01.019

(11)Yang,X.Y.;Vantomme,A.;Xiao,F.S.;Su,B.L.Catalysis Today 2007,128(3),123.doi:10.1016/j.cattod.2007.07.017

(12)Liu,Q.;You,J.L.;Wang,Y.Y.;Wang,C.Y.;Wang,J.;Liu,X. W.Spectroscopy and SpectralAnalysis 2013,10,31.[劉欽,尤靜林,王媛媛,王晨陽,王靜,劉曉偉.光譜學與光譜分析,2013,10,31.]doi:10.3964/j.issn.1000-0593(2013)10-2705-06

(13)Wang,Y.Y.;You,J.L.;Jiang,G.C.Chinese Journalof Inorganic Chemistry 2008,24(5),765.[王媛媛,尤靜林,蔣國昌.無機化學學報,2008,24(5),765.]

(14)Wang,J.;You,J.L.;Wang,Y.Y.Journal ofSynthetic Crystals 2013,3,4.[王靜,尤靜林,王媛媛.人工晶體學報,2013,3, 4.]

(15)Wu,J.;Stebbins,J.F.JournalofNon-Crystalline Solids2010, 356(41),2097.doi:10.1016/j.jnoncrysol.2010.08.015

(16)You,J.;Jiang,G.;Hou,H.Isomorphic Representations of Hyperfine Structureof Binary Silicatesby Interior Stress, VibrationalWavenumber and Special Fractional Dimension.In JournalofPhysics:Conference Series;IOPPublishing Ltd.: Singapore2006;p 25.

(17)You,J.L.;Jiang,G.C.;Hou,H.Y.Journal ofRaman Spectroscopy 2005,36(3),237.doi:10.1002/jrs.1287

(18)Wang,W.;You,J.L.Thermodynam ic Calculation for LiO2-SiO2and Na2O-SiO2Binary MeltBased on Cluster Structure ModelofMolten Silicates.ShanghaiUniversity,Shanghai, 2009.[王威,尤靜林.基于團簇結構模型的Li2O-SiO2、Na2O-SiO2二元系熔體的熱力學性質計算[D].上海:上海大學,2009.]

(19)Park,S.Y.;Lee,S.K.Geochimica EtCosmochimica Acta 2014,147,26.doi:10.1016/j.gca.2014.10.019

(20)Neuville,D.;Cormier,L.;Ligny,D.American Mineralogist 2008,93(1),228.doi:10.2138/am.2008.2646

(21)You,J.L.;Wang,Y.Y.;Jiang,G.C.;Xie,J.;Yu,X.H.;Mo,X. X.Quantum Chemistry ab initio Calculation Study of Structure and Raman Spectraof A lum inosilicates.In Xxii InternationalConference on Raman Spectroscopy, Proceedingsof the AIPConference,Boston,America,Aug 8-13,2010;Champion,P.M.,Ziegler,L.D.Eds.;Elsevier: USA,2010.

(22)Segall,M.;Linda,P.;Probert,M.Materials Studio,Version 5.5;Accelrys:American,2010.

(23)Thiel,W.;Voityuk,A.A.Theoretica Chimica Acta 1992,81 (6),391.doi:10.1007/BF01134863

(24)Thiel,W.;Voityuk,A.A.The JournalofPhysicalChemistry 1996,100(2),616.doi:10.1021/jp952148o

(25)Farnan,I.;Stebbins,J.Journalofthe American Chemical Society 1990,112(1),32.doi:10.1021/ja00157a008

(26)Rui,L.PhysicalChemistry of Inorganic Materials;Chinese Construction Industry Press:Beijing,1986;p 194.[瑞倫.無機材料物理化學.北京:中國建筑工業出版社,1986:194.]

(27)Roine,A.HSCChemistry,Ver6.1;OutokumpuOy:Finland, 2002.

Molar Gibbs Mixing Free Energy Calcu lation for a Na2O-A l2O3-SiO2Ternary Melt Based on the Cluster Model

LIYong-Kang YOU Jing-Lin*WANG Jian WANG Min MA Nan WEIGuang-Chao
(State Key Laboratory ofAdvanced Special Steel,ShanghaiKey Laboratory ofAdvanced Ferrometallurgy, SchoolofMaterials Science and Engineering,ShanghaiUniversity,Shanghai200072,P.R.China)

Based on the clustermode lofmolten ternary silicates,the thermodynam ic properties,including entropy,enthalpy,and hea t capacity o f Na2O-A l2O3-SiO2a t 1473,1873,and 2000K w ere ca lculated by the modified neglecto fdifferentia loverlap(MNDO/d)sem i-em piricalmethod based on the prim itive assum ption of clusters in themelt.Themixing free energies of the Na2O-Al2O3-SiO2ternary system were derived.Them ixing free energy of the Na2O-Al2O3-SiO2ternary system is the sum ofa ll the clusterunits thatexistaccording to the Boltzmann distribution law.The thermodynam ic properties of this ternary silicate melt depend on its m icrostruc ture.

Na2O-A l2O3-SiO2;Cluste rmode l;Me lt;Mola rm ixing free ene rgy

1 引言

硅酸鹽及含鋁硅酸鹽熔體的熱力學性質一直是冶金、地學、玻璃和陶瓷等領域重要的研究課題之一。而對于硅酸鹽的熔體熱力學性質的研究,已經有文獻分別從分子模型、分子離子共存模型、聚合模型等方面計算和討論硅酸鹽活度和組分之間的關系,進而利用活度求出熱力學數據1-5。分子離子共存模型和聚合模型是應用比較廣泛的兩個熱力學模型。分子離子共存理論模型最初是由前蘇聯丘意柯教授基于實驗事實提出,張鑒2用此模型對二元硅酸鹽熔體、三元硅酸鹽熔體和熔鹽做了驗證。很多考慮到硅酸鹽中的聚合反應的模型都應用到硅酸鹽的熱力學性質預測中并得到后續的發展,包括Lin和Pelton5,Masson6,Ottonello7,Zaitsev等8的研究。聚合模型是假設在給定的溫度壓強下各成分下的硅酸鹽熔體由金屬陽離子、離子氧化物和單元均衡分布聚合而成。關于硅酸鹽聚合理論的研究已經擴展到三元體系當中。Glibin和King9基于當前對硅酸鹽熔體結構的認識,以Guggenheim的準化學模型預測了二元體系的硅酸鹽和鋁硅酸鹽的混合熱力學函數,他們利用氧化物晶體、含氧分子和離子中已有的結構和熱力學數據計算模型中所需的參數,計算結果與實驗數據吻合較好。

September28,2015;Revised:December28,2015;Published onWeb:December28,2015.

O642

10.3866/PKU.WHXB201512281

*Corresponding author.Email:jlyou@163.com;Tel:+86-13331992297.

Theprojectwas supported by the Key Projectof NationalNatural Science Foundation of China(50932005),NationalNaturalScience Foundation of China(20973107,40973046),ShanghaiCommitteeof Scienceand Technology,China(12520709200),and CSIROMinerals Down Under Flagship of Australia.

國家自然科學基金重點項目(50932005),國家自然科學基金(20973107,40973046),上海科學技術委員會(12520709200)和澳大利亞Flagship的科學與工業研究組織CSIRO資助項目?Editorialofficeof Acta Physico-Chimica Sinica

猜你喜歡
體系結構模型
一半模型
《形而上學》△卷的結構和位置
哲學評論(2021年2期)2021-08-22 01:53:34
構建體系,舉一反三
重要模型『一線三等角』
重尾非線性自回歸模型自加權M-估計的漸近分布
論結構
中華詩詞(2019年7期)2019-11-25 01:43:04
論《日出》的結構
3D打印中的模型分割與打包
創新治理結構促進中小企業持續成長
現代企業(2015年9期)2015-02-28 18:56:50
“曲線運動”知識體系和方法指導
主站蜘蛛池模板: 韩国自拍偷自拍亚洲精品| 免费一级成人毛片| 欧美成人国产| 国产精品原创不卡在线| 久久香蕉国产线看观看式| 国产欧美日韩视频怡春院| 玩两个丰满老熟女久久网| 无码视频国产精品一区二区| 秋霞午夜国产精品成人片| 国产亚洲欧美日韩在线一区二区三区| 国内熟女少妇一线天| 2021亚洲精品不卡a| 午夜无码一区二区三区| 无码内射在线| 国产黄在线免费观看| 真实国产乱子伦视频| 97se亚洲综合在线天天| 真人高潮娇喘嗯啊在线观看 | 婷婷六月在线| 国产一区二区色淫影院| 欧美国产在线精品17p| 欧美伦理一区| 久久青草精品一区二区三区| 97超爽成人免费视频在线播放| 国内精品视频在线| 免费国产小视频在线观看| 少妇高潮惨叫久久久久久| 亚洲精品午夜天堂网页| 国产精品3p视频| 国产18在线播放| 精品无码人妻一区二区| 亚洲一区二区三区国产精华液| 久草视频一区| 91精品aⅴ无码中文字字幕蜜桃| 亚洲国产精品不卡在线 | 国产色网站| 亚洲AV电影不卡在线观看| 97se亚洲综合在线韩国专区福利| 无码网站免费观看| 国产又粗又猛又爽视频| 真人高潮娇喘嗯啊在线观看| 国产中文在线亚洲精品官网| 色妞永久免费视频| 色综合五月| 成年免费在线观看| 天天干天天色综合网| 国产无码精品在线播放| 国产视频大全| 国产精品xxx| 丝袜美女被出水视频一区| 有专无码视频| 国产人妖视频一区在线观看| 中文字幕欧美日韩| 欧美国产在线精品17p| 亚洲中文无码av永久伊人| 国产日韩丝袜一二三区| 亚洲资源站av无码网址| 国产小视频a在线观看| 国产一级二级在线观看| 国内视频精品| 麻豆AV网站免费进入| 91丝袜美腿高跟国产极品老师| 久久永久视频| 无码专区在线观看| 欧美午夜视频| 中文字幕第1页在线播| 国产一区二区人大臿蕉香蕉| 国产导航在线| 国产精品99久久久| 亚洲二区视频| 91精品啪在线观看国产91九色| 午夜一级做a爰片久久毛片| 2021天堂在线亚洲精品专区| 992tv国产人成在线观看| 日韩在线欧美在线| 自拍偷拍一区| 青青草原偷拍视频| 色婷婷视频在线| aaa国产一级毛片| 日韩国产黄色网站| 久久免费观看视频| 成人综合久久综合|