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

雙[1,2-二(三氟甲基)乙烯-1,2-二硫基]鎳與丁二烯反應的溶劑效應

2010-11-30 10:49:12孫麗麗趙月紅韓清珍
物理化學學報 2010年12期

孫麗麗 趙月紅 韓清珍 溫 浩,*

(1中國科學院過程工程研究所多相復雜系統國家重點實驗室,北京 100190;2中國科學院研究生院,北京 100049)

雙[1,2-二(三氟甲基)乙烯-1,2-二硫基]鎳與丁二烯反應的溶劑效應

孫麗麗1,2趙月紅1韓清珍1溫 浩1,*

(1中國科學院過程工程研究所多相復雜系統國家重點實驗室,北京 100190;2中國科學院研究生院,北京 100049)

采用密度泛函理論(DFT)在B3LYP/6-31G(d)水平上研究了雙[1,2-二(三氟甲基)乙烯-1,2-二硫基]鎳(Ni[S2C2(CF3)2]2)與丁二烯的反應機理.采用極化連續介質模型(PCM),考察了溶劑對各反應駐點的電荷分布、偶極矩、溶劑化自由能的影響.計算結果表明:Ni[S2C2(CF3)2]2與丁二烯的反應為前線軌道對稱性匹配的協同反應,溶劑介電常數的增大有利于穩定各反應駐點.同時在同種溶劑中,過渡態和產物穩定的程度大于反應物,從而反應更加容易進行.

密度泛函理論;雙[1,2-二(三氟甲基)乙烯-1,2-二硫基]鎳;丁二烯;溶劑效應

近年的研究指出,雙[1,2-二(三氟甲基)乙烯-1,2-二硫基]鎳(Ni[S2C2(CF3)2]2)能夠選擇性并可逆地與乙烯、丙烯等小分子單烯烴或共軛二烯烴發生絡合反應,生成1:1加合產物[1-4].這類反應的活性中心是Ni[S2C2(CF3)2]2分子中的S原子,與金屬鹽溶液相比, Ni[S2C2(CF3)2]2具有不易被炔烴、CO、H2S等雜質毒化的特點[1].因此,人們希望將其用于烷烴與烯烴的絡合分離[5].由于Ni[S2C2(CF3)2]2與烯烴的反應是在溶液中進行的,溶劑將對參與反應的分子電荷分布和相互作用產生影響,進而影響整個反應過程[6].已有的研究顯示,Ni[S2C2(CF3)2]2與小分子單烯烴的反應是分兩步完成的環加成反應,極性溶劑有利于降低反應的活化能和提高產物的穩定性[7-8].

關于Ni[S2C2(CF3)2]2與共軛二烯烴反應的實驗研究指出,Ni[S2C2(CF3)2]2與共軛二烯烴的反應是熱力學可逆的,反應速率明顯高于Ni[S2C2(CF3)2]2與小分子單烯烴的反應.其中,Ni[S2C2(CF3)2]2與丁二烯的反應速率最高,反應產物最易分解[9-12].由于Ni[S2C2(CF3)2]2與丁二烯的反應是在溶液中進行的,研究該反應的溶劑效應對溶劑選擇將有重要的參考價值.

本文采用DFT方法討論Ni[S2C2(CF3)2]2與丁二烯的絡合反應歷程,采用PCM模型[13-14]考察該反應體系的溶劑效應,通過分析不同極性的溶劑中分子的電荷分布、偶極矩和溶劑化自由能,研究溶劑對Ni[S2C2(CF3)2]2與丁二烯反應的影響.

1 計算方法

本文討論的反應體系是在溫度T=298.15 K,壓強p=101.3 kPa條件下,雙[1,2-二(三氟甲基)乙烯-1, 2-二硫基]鎳與丁二烯在氣相和7種常用溶劑中的反應.雙[1,2-二(三氟甲基)乙烯-1,2-二硫基]鎳的分子結構如圖1所示.為簡化表達,分別用Ni[S2C2(CF3)2]2和C4H6·Ni[S2C2(CF3)2]2表示雙[1,2-二(三氟甲基)乙烯-1,2-二硫基]鎳及其與丁二烯的反應產物.

本文采用的7種常用溶劑及其介電常數ε分別為:庚烷(ε=1.920)、苯(ε=2.247)、甲苯(ε=2.379)、二乙基醚(ε=4.335)、氯苯(ε=5.621)、苯胺(ε=6.890)和四氫呋喃(ε=7.580).

本文的計算均采用密度泛函理論在Gaussian 03程序包[15]中完成.

對于氣相中Ni[S2C2(CF3)2]2與丁二烯反應體系的計算和分子幾何構型優化,Ni原子采用相對論有效核勢近似(ECP)[16-17]和修正的LANL2DZ基組,即將兩個最外層p軌道的指數和系數用優化的Ni(41)分裂價基[18]代替,并添加f極化函數[19];對反應中心的S原子添加極化函數[20],對反應中心的C原子增加彌散函數,采用6-31+G(d)基組提高計算精度;對遠離反應中心的C、F和H原子采用6-31G基組[21-22].

圖1 雙[1,2-二(三氟甲基)乙烯-1,2-二硫基]鎳的分子結構Fig.1 Molecular structure of bis[1,2-di(trifluoromethyl) ethylene-1,2-dithiolato]nickel

本文采用過渡態(TS)方法搜索反應的過渡態,并對過渡態進行振動頻率分析確認.進一步,對過渡態進行內稟反應坐標(IRC)計算,以驗證過渡態的正確性.

本文采用PCM模型,在溫度T=298.15 K,壓強p=101.3 kPa條件下,考察Ni[S2C2(CF3)2]2與丁二烯反應的溶劑效應.由于對于此反應體系,采用PCM模型時溶劑場中做幾何優化難以收斂,尤其是對于過渡態的優化.本文采用在氣相中優化的分子構型,在溶劑場中做單點能計算和頻率分析.

由反應物和產物在氣相和溶液中的Gibbs自由能可計算反應物和產物的溶劑化自由能ΔGsolv、反應的活化自由能Δ≠G和反應的自由能變化ΔG.

式中,(Gi)gasphas和(Gi)solution分別為氣相和溶液中物質i的Gibbs自由能,下角中的i為Ni[S2C2(CF3)2]2、丁二烯和C4H6·Ni[S2C2(CF3)2]2.

2 結果與討論

2.1 Ni[S2C2(CF3)2]2與丁二烯的反應歷程

圖2 Ni[S2C2(CF3)2]2和丁二烯分子的LUMOs和HOMOsFig.2 LUMOs and HOMOs of Ni[S2C2(CF3)2]2and butadiene moleculesLUMO:the lowest unoccupied molecular orbital; HOMO:the highest occupied molecular orbital

圖3 Ni[S2C2(CF3)2]2與丁二烯的反應歷程Fig.3 Reaction course of Ni[S2C2(CF3)2]2and butadiene

已有文獻[2,11]指出,Ni[S2C2(CF3)2]2與共軛二烯烴的反應屬前線軌道對稱性允許的反應.如圖2所示, Ni[S2C2(CF3)2]2分子的最低空軌道(LUMO)和最高占據軌道(HOMO)都主要由S原子的p軌道和C原子的p軌道組合成的垂直分子平面的π分子軌道構成, Ni[S2C2(CF3)2]2分子的LUMO的4個S原子的位相不完全相同,而HOMO的4個S原子的位相相同.丁二烯分子的LUMO和HOMO主要由C原子的p軌道組合而成,丁二烯分子的LUMO包含一個節面, HOMO則包含了兩個節面.Ni[S2C2(CF3)2]2分子的LUMO與丁二烯分子的HOMO的軌道為對稱性匹配,軌道能級差為ΔELH=0.571 eV.Ni[S2C2(CF3)2]2分子的HOMO與丁二烯分子的LUMO的軌道對稱性也匹配,但軌道能級差較大,ΔELH=6.667 eV.因此, Ni[S2C2(CF3)2]2與丁二烯的反應是電子從丁二烯的HOMO躍遷到Ni[S2C2(CF3)2]2的LUMO而成鍵的反應.Ni[S2C2(CF3)2]2的LUMO與丁二烯的HOMO沿著垂直于分子平面的方向接近,可以使二者有效重疊,形成兩個新的S—C鍵.

分子的幾何結構參數(鍵長、鍵角和二面角)在反應過程中的變化,能夠在一定程度上體現出Ni[S2C2(CF3)2]2分子與丁二烯分子的成鍵過程.圖3和表1給出Ni[S2C2(CF3)2]2與丁二烯的反應歷程和分子結構參數的變化.

隨著反應的進行,丁二烯中共軛C=C雙鍵(C1=C2鍵和C3=C4鍵)的鍵長LC1—C2和LC3—C4將由0.134 nm逐漸伸長至產物C4H6·Ni[S2C2(CF3)2]2中的0.150 nm,形成新的C—C單鍵.丁二烯中的共軛C—C單鍵(C2—C3鍵)的鍵長LC2—C3則由0.147 nm逐漸縮短至0.134 nm,成為新的C=C雙鍵.在此環加成反應中,丁二烯分子的π電子將向Ni[S2C2(CF3)2]2分子轉移而成鍵.Ni[S2C2(CF3)2]2分子中,作為反應中心的S原子與丁二烯的C1和C4原子逐漸接近并成鍵,C1—S1鍵和C4—S2鍵的鍵長LC1—S1和LC4—S2將從過渡態的0.260 nm逐漸縮短至產物C4H6· Ni[S2C2(CF3)2]2中的0.190 nm.Ni[S2C2(CF3)2]2分子中的成鍵S原子和丁二烯分子中的成鍵C原子逐漸由sp2軌道雜化轉變為sp3軌道雜化.由過渡態到產物C4H6·Ni[S2C2(CF3)2]2,鍵角∠Ni1—S1—C1、∠Ni1—S2—C4、∠S1—C1—C2和∠S2—C4—C3均表現出一定程度的增大.由于Ni[S2C2(CF3)2]2分子與丁二烯分子的前線軌道對稱性匹配,二者相互平行地靠近即可實現軌道的有效重疊而成鍵.因此,在反應過程中,Ni[S2C2(CF3)2]2分子中的Ni、S和C原子基本處于同一平面,二面角∠S1—S2—S3—S4和∠S1—Ni1—S4—S3均保持在0°或180°附近.丁二烯以鄰位交叉或順式構象參與反應,在反應過程中丁二烯的4個C原子基本處于同一平面,二面角∠C1—C2—C3—C4在反應過程中保持在0°附近.

關于C4H6·Ni[S2C2(CF3)2]2的分子結構目前尚未見有關文獻報道,但可參考文獻[12]中關于Ni[S2C2(CF3)2]2與2,3-二甲基丁二烯的反應產物C6H10·Ni[S2C2(CF3)2]2的分子結構.其中,NiS4單元及其絡合的共軛二烯烴均為平面構型,成鍵的S原子和C原子均為四面體雜化.這些結構特征與本文得到的C4H6·Ni[S2C2(CF3)2]2的幾何構型特征是完全一致的.

表1 各反應駐點的幾何結構參數Table 1 Structural parameters of the stationary points

圖4 Ni[S2C2(CF3)2]2和丁二烯反應的能量曲線Fig.4 Energy curve of the reaction between Ni[S2C2(CF3)2]2 and butadiene

我們對Ni[S2C2(CF3)2]2與丁二烯反應的過渡態做IRC計算,驗證過渡態的可靠性,能量曲線如圖4所示,其中正的反應坐標方向趨近于反應物,負的反應坐標方向趨近于產物.

2.2 溶劑對原子的電荷分布的影響

表2給出氣相和溶劑中過渡態和反應產物中部分原子(Ni1、S1、C5、C1、C2等)和結構片段(—Ni[S2C2(CF3)2]2—、—C4H6—)的靜電勢(ESP)電荷分布.從中可以看出,在過渡態和反應產物中,結構片段—Ni[S2C2(CF3)2]2—均帶有部分負電荷,—C4H6—則相應地帶有等量正電荷.顯然,在Ni[S2C2(CF3)2]2與丁二烯的反應中,出現從結構片段—C4H6—向—Ni[S2C2(CF3)2]2—的電子轉移.同時,溶劑極性增強,將使結構片段—Ni[S2C2(CF3)2]2—所帶的負電荷數增大;相應地,結構片段—C4H6—所帶的正電荷數也出現等量的增加.這說明隨著溶劑極性的增強,電子在結構片段—Ni[S2C2(CF3)2]2—上出現的概率更大,使—Ni [S2C2(CF3)2]2—與—C4H6—之間的靜電相互作用增強,有利于過渡態和反應產物的穩定.

表3 各反應駐點的偶極矩(μ)Table 3 The dipole moment(μ)of the stationary points

同時在溶劑場中各原子的電荷分布也發生了變化,無論是過渡態還是產物,隨著溶劑極性的增強, Ni1原子的正電荷數增加,S1原子和S4原子所帶的正電荷數減少(過渡態中S1原子的負電荷數增加),說明由Ni原子向S原子轉移電子,同時由于溶劑介電常數增大,結構片段—C4H6—向—Ni[S2C2(CF3)2]2—的S1和S2原子轉移的電子數有顯著增加.然而對于過渡態和產物,S1原子的電荷數改變很小,說明在溶劑場中S1原子向C5原子轉移部分電子.對于C5原子和C8原子,C8原子所帶的負電荷數大于C5原子的負電荷數,這說明在溶劑場中,對于過渡態和產物出現了由靠近Ni[S2C2(CF3)2]2反應中心的Ni1、S1和C5原子向距離反應中心較遠的C8和S4原子的電子轉移.

表2 過渡態和反應產物中部分原子和結構片段的靜電勢(ESP)電荷分布(e)Table 2 Electrostatic potential(ESP)charge denstity(e)of some atoms and substructures in TS and product

表4 各駐點的溶劑化自由能(ΔGsolv)、反應的活化自由能(Δ≠G)和反應的自由能變(ΔG)Table 4 Solvation free energy(ΔGsolv)of the stationary points,activation free energy of the reaction(Δ≠G)and free energy change of the reaction(ΔG)

表3給出了溶劑對反應物、過渡態和產物分子偶極矩μ的影響.其中,由于Ni[S2C2(CF3)2]2分子具有D2h對稱性,分子的偶極矩為0.000,因此Ni[S2C2(CF3)2]2分子在溶劑中出現的電荷分布變化并不會影響分子偶極矩.與Ni[S2C2(CF3)2]2分子相比,丁二烯、過渡態和產物分子均表現出偶極矩隨溶劑極性的增強而增大的趨勢,這是由于丁二烯、過渡態和產物分子都具有C1對稱性,正負電荷中心不重合,它們將在溶劑中產生誘導偶極,從而使丁二烯、過渡態和產物分子的偶極矩隨溶劑極性的變化而發生改變.然而在同種溶劑中,反應物、過渡態和產物分子產生的誘導偶極矩不同.過渡態和產物分子的誘導偶極矩大于反應物分子的誘導偶極矩.這種誘導偶極與溶劑間的相互作用更有利于過渡態和產物分子的穩定.

2.3 溶劑化效應對自由能的影響

表4中給出了氣相和溶劑中各反應駐點的溶劑化自由能ΔGsolv、反應的活化自由能Δ≠G和反應的自由能ΔG.表4中數據顯示,雖然各反應駐點的溶劑化過程均為自發過程(ΔGsolv<0),但不同反應駐點的ΔGsolv有明顯差別.在同種溶劑中,過渡態和產物的溶劑化自由能均小于反應物的溶劑化自由能之和.同時,各反應駐點的ΔGsolv將隨溶劑極性的增強呈現出持續下降的趨勢.特別是當溶劑的介電常數1.000<ε≤4.335時,影響較為顯著.表4說明,選擇具有一定極性的溶劑將有利于反應的進行和產物的穩定.各反應駐點ΔGsolv的差別及其隨溶劑極性的變化,應可歸因于溶劑分子與各反應駐點間相互作用的差異.實際上,溶劑化自由能ΔGsolv包含靜電自由能ΔGel和非靜電自由能ΔGnon的貢獻[13]:

ΔGsolv=ΔGel+ΔGnon=ΔGel+ΔGcav+ΔGdis-rep(4)其中,靜電自由能ΔGel是溶質分子與連續介質(溶劑)因靜電作用導致的體系能量變化;非靜電自由能ΔGnon體現為空穴自由能ΔGcav和分散-排斥自由能ΔGdis-rep的共同影響,ΔGcav是連續介質中因騰出空穴以容納溶質分子而導致的體系能量變化,ΔGdis-rep是連續介質與溶質之間范德華相互作用而導致的體系能量變化.溶劑介電常數雖可在一定程度上體現溶質與溶劑的靜電相互作用,但并不能表現非靜電相互作用[13].

表4同時列出了在氣相和溶劑中,Ni[S2C2(CF3)2]2與丁二烯反應的活化自由能Δ≠G和反應自由能ΔG.表4中數據顯示,與氣相相比,所研究的溶劑均能夠降低反應的活化自由能Δ≠G和反應自由能ΔG,使Ni[S2C2(CF3)2]2與丁二烯的反應更易于進行.

3 結論

采用密度泛函理論和極化連續介質模型考察溶劑對Ni[S2C2(CF3)2]2與丁二烯反應的影響.結果表明:Ni[S2C2(CF3)2]2與丁二烯的反應是前線軌道對稱性匹配的環加成反應.隨溶劑介電常數的增大,過渡態和產物中結構片段—C4H6—向結構片段—Ni[S2C2(CF3)2]2—轉移電子數增加,二者的靜電相互作用增強,有利于穩定過渡態和反應產物.同時電子更傾向從靠近—Ni[S2C2(CF3)2]2—反應中心的Ni、C和S原子向遠離反應中心的C原子和S原子轉移.在同種溶劑中,過渡態和反應產物在溶劑場中產生的誘導偶極矩大于反應物的誘導偶極矩,因此溶劑對過渡態和反應產物的穩定程度大于反應物.同時在所考察的每種溶劑中,過渡態和產物的溶劑化自由能小于反應物的溶劑化自由能之和,因此相比于氣相,此反應在溶劑中的活化自由能壘和反應的自由能均降低,反應更容易進行.特別是當溶劑的介電常數為1.000<ε≤4.335時,影響較為顯著.同時考慮到絡合劑Ni[S2C2(CF3)2]2為非極性分子,在極性弱的溶劑中的溶解性更好,此反應宜選用弱極性溶劑.

致謝:作者向得克薩斯州立大學的Yubo Fan博士給予的指導致以誠摯的謝意!

1 Wang,K.;Stiefel,E.I.Science,2001,291:106

2 Fan,Y.B.;Hall,M.B.J.Am.Chem.Soc.,2002,124:12076

3 Harrison,D.J.;Nguyen,N.;Lough,A.J.;Fekl,U.J.Am.Chem. Soc.,2006,128:11026

4 Szilagyi,R.K.;Lim,B.S.;Glaser,T.;Holm,R.H.;Hedman,B.; Hodgson,K.O.;Solomon,E.I.J.Am.Chem.Soc.,2003,125: 9158

5 Smith,R.S.;Herrera,P.S.;Henderson,J.F.;Spence,R.E.V.H. Selective chemical binding for olefins/paraffins separation.CA, 2415064[P].2004-06-23

6 Field,M.J.;Bash,P.A.;Karplus,M.J.Comput.Chem.,1990,11: 700

7 Han,Q.Z.;Geng,C.Y.;Zhao,Y.H.;Qi,C.S.;Wen,H.Acta Phys.Sin.,2008,57:96 [韓清珍,耿春宇,趙月紅,戚傳松,溫 浩.物理學報,2008,57:96]

8 Han,Q.Z.;Zhao,Y.H.;Wen,H.Mol.Simulat.,2008,34:631

9 Schrauze,G.;Ho,R.K.Y.;Murillo,R.P.J.Am.Chem.Soc., 1970,92:3508

10 Baker,J.R.;Hermann,A.;Wing,R.M.J.Am.Chem.Soc.,1971, 93:6486

11 Clark,G.R.;Waters,J.M.;Whittle,K.R.J.Chem.Soc.Dalton Trans.,1973:821

12 Herman,A.;Wing,R.M.J.Organomet.Chem.,1973,63:441

13 Cossi,M.;Barone,V.;Cammi,R.;Tomasi,J.Chem.Phys.Lett., 1996,255:327

14 Barone,V.;Cossi,M.J.Phys.Chem.A,1998,102:1995

15 Frisch,M.J.;Trucks,G.W.;Schlegel,H.B.;et al.Gaussian 03. RevisionA.01.Pittsburgh,PA:Gaussian Inc.,2003

16 Hay,P.J.;Wadt,W.R.J.Chem.Phys.,1985,82:270

17 Hay,P.J.;Wadt,W.R.J.Chem.Phys.,1985,82:299

18 Couty,M.;Hall,M.B.J.Comput.Chem.,1996,17:1359

19 Ehlers,A.W.;Bohme,M.;Dapprich,S.;Gobbi,A.;Hollwarth, A.;Jonas,V.;Kohler,K.F.;Stegmann,R.;Veldkamp,A.; Frenking,G.Chem.Phys.Lett.,1993,208:111

20 Check,C.E.;Faust,T.O.;Bailey,J.M.;Wright,B.J.;Gilbert,T. M.;Sunderlin,L.S.J.Phys.Chem.A,2001,105:8111

21 Petersson,G.A.;Allaham,M.A.J.Chem.Phys.,1991,94:6081

22 Petersson,G.A.;Bennett,A.;Tensfeldt,T.G.;Allaham,M.A.; Shirley,W.A.;Mantzaris,J.J.Chem.Phys.,1988,89:2193

July 7,2010;Revised:August 20,2010;Published on Web:November 10,2010.

Solvent Effect in the Reaction between Bis[1,2-di(trifluoromethyl) ethylene-1,2-dithiolato]Nickel and Butadiene

SUN Li-Li1,2ZHAO Yue-Hong1HAN Qing-Zhen1WEN Hao1,*
(1State Key Laboratory of Multiphase Complex System,Institute of Process Engineering,Chinese Academy of Sciences,Beijing 100190,P.R.China;2Graduate University of Chinese Academy of Sciences,Beijing 100049,P.R.China)

We studied the reaction mechanism for the reaction between bis[1,2-di(trifluoromethyl) ethylene-1,2-dithiolato]nickel(Ni[S2C2(CF3)2]2)and butadiene by density functional theory(DFT)at the B3LYP/6-31G(d)level.The solvent effect on the charge distribution,dipole moment,and solvation free energies of the stationary points were investigated using the polarizable continuum model(PCM).The calculation results showed that this reaction was orbital symmetry allowed and concerted.The reaction stationary points become more stable with an increase of solvent dielectric constant.Additionally,the degree of stabilization for the transition state and the product is larger than that of the reactants in the same solvent,which means that the reaction occurs more easily.

Density functional theory;Bis[1,2-di(trifluoromethyl)ethylene-1,2-dithiolato]nickel; Butadiene; Solvent effect

O641

?Corresponding author.Email:hwen@home.ipe.ac.cn;Tel:+86-10-62626704.

The project was supported by the National Natural Science Foundation of China(20703047,20821092).國家自然科學基金(20703047,20821092)資助項目

主站蜘蛛池模板: 在线综合亚洲欧美网站| 91久久夜色精品| 91无码人妻精品一区二区蜜桃| 视频二区国产精品职场同事| 久草热视频在线| 国产国拍精品视频免费看| 亚洲AV无码不卡无码| 青草娱乐极品免费视频| 不卡视频国产| 日本欧美中文字幕精品亚洲| 欧美人在线一区二区三区| 国产日韩欧美视频| 日本爱爱精品一区二区| а∨天堂一区中文字幕| 亚洲精品国产精品乱码不卞| 亚洲成aⅴ人片在线影院八| a毛片在线| 直接黄91麻豆网站| 日本午夜三级| 99青青青精品视频在线| 国产精品视频猛进猛出| 亚洲天堂高清| 欧美视频在线第一页| 综合网久久| 色欲国产一区二区日韩欧美| 精品国产欧美精品v| 色悠久久综合| h网站在线播放| 国产成人1024精品下载| 老司机aⅴ在线精品导航| 97青青青国产在线播放| 国产欧美中文字幕| 91美女视频在线| 亚洲最大综合网| 2021国产精品自产拍在线观看| 欧美一道本| 欧美一区二区三区国产精品| 色窝窝免费一区二区三区| 国产精品自在拍首页视频8| 日韩高清无码免费| 超碰精品无码一区二区| 亚洲精品在线影院| 小蝌蚪亚洲精品国产| 亚洲美女一区| 91黄色在线观看| 久久情精品国产品免费| 国产精品欧美激情| 97人妻精品专区久久久久| 久草视频一区| 亚洲精品中文字幕午夜| 成人一级黄色毛片| 亚洲香蕉久久| 中文字幕66页| 日韩精品免费一线在线观看 | 亚洲一级无毛片无码在线免费视频| 日日摸夜夜爽无码| 精品91在线| 亚洲色图欧美激情| 欧美视频在线播放观看免费福利资源 | 日韩小视频在线播放| 国产午夜看片| 久久亚洲AⅤ无码精品午夜麻豆| 曰AV在线无码| 亚洲综合在线网| 国产精品不卡永久免费| 亚洲va欧美va国产综合下载| 中日韩欧亚无码视频| 国产白浆一区二区三区视频在线| 操操操综合网| 久久精品亚洲热综合一区二区| 婷婷亚洲视频| 2021精品国产自在现线看| 很黄的网站在线观看| 婷婷色丁香综合激情| 免费jjzz在在线播放国产| aaa国产一级毛片| 欧美国产精品不卡在线观看| 99久久精品久久久久久婷婷| 亚洲开心婷婷中文字幕| 91九色最新地址| 丁香亚洲综合五月天婷婷| 成年看免费观看视频拍拍|