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

基于新UNIFAC基團的尼龍66鹽溶解度的計算方法

2016-06-24 06:48:43顧雪萍田璐璐馮連芳張才亮浙江大學化學工程聯合國家重點實驗室浙江杭州310027
化工學報 2016年2期

顧雪萍,田璐璐,馮連芳,張才亮(浙江大學化學工程聯合國家重點實驗室,浙江 杭州 310027)

基于新UNIFAC基團的尼龍66鹽溶解度的計算方法

顧雪萍,田璐璐,馮連芳,張才亮
(浙江大學化學工程聯合國家重點實驗室,浙江 杭州 310027)

摘要:尼龍66(PA66)鹽在水中溶解度的準確計算是尼龍66聚合過程模型化基礎。選擇UNIFAC活度系數方法,針對PA66聚合體系對組分進行基團劃分,定義的新基團為羧基和氨基形成的締合基團(—CH2COO?·+H3NCH2—)、羧基與其相鄰的亞甲基(—CH2COOH)。利用PA66鹽-水的固液平衡實驗數據以及PA66鹽的熔點和熔化焓,回歸得到UNIFAC模型的基團交互作用參數。基于構建的UNIFAC物性模型預測了0~100℃范圍內PA66鹽在水中的溶解度,與實驗數據的平均相對誤差為2.1%;采用該模型進一步計算工業操作條件(120℃)下的鹽濃度,與Aspen自帶參數的UNIFAC模型比較,其誤差從20%左右降到5%以下。

關鍵詞:尼龍66鹽;水溶液;溶解度;UNIFAC模型;UNIFAC基團

2015-08-03收到初稿,2015-12-09收到修改稿。

聯系人:馮連芳。第一作者:顧雪萍(1972—),女,博士,副研究員。

Received date: 2015-08-03.

引 言

尼龍66(PA66)以其優越的性能而被廣泛應用于電子、汽車、紡織等領域,世界PA66的產能集中在歐美國家,國內產能最大的河南平頂山神馬集團僅為20萬~30萬噸,單線產能低,遠遠不能滿足需求[1]。對PA66聚合過程建立模型,消化吸收和改造所引進的工藝和裝置是提升產能的有效辦法。

根據PA66生產工藝流程,PA66鹽溶液在濃縮槽中升溫濃縮,若溫度過低,溶解鹽以固體鹽形式析出,若溫度過高,鹽溶液開始聚合;在高溫高壓聚合過程中,隨著溫度的升高和反應時間的延長,系統中的水和PA66鹽的含量均下降,PA66鹽在水中的溶解度決定了體系中水的脫除速率。反應體系中存在溶解鹽、二聚體、己二酸、產物PA66、低聚物、小分子水、CO2和NH3、易揮發物己二胺等組分。因此針對PA66聚合過程確定合適的溶解度計算模型及準確的模型參數尤為重要。

對組分的基團進行劃分,定義新基團并回歸參數,得到準確的尼龍66鹽在水中的溶解度,為實現尼龍66生產過程的工藝優化提供指導。

1 物性計算方法選擇和參數確定

熱力學性質計算模型有兩類,一是活度系數模型,二是狀態方程模型。前者不考慮液體的可壓縮性,適用于極性較強的體系;后者適用于非極性分子,可壓縮流體[2]。而逐步聚合體系的單體大都是極性類物質,故對于PA66聚合體系選擇活度系數模型比較合適。

活度系數模型采用活度系數修正了假設理想溶液狀態下的Raoult定律,Schaffer等[3]采用Flory-Huggins(FH)模型僅計算了尼龍熔融縮聚過程中水的活度系數,目前不能準確知道PA66體系聚合物-溶劑交互作用參數,該模型不適用。Chen[4]將NRTL模型[5]和FH模型組合起來應用于聚合物體系,提出了Poly-NRTL模型,模型涉及多個二元交互作用參數;Seavey等[6]采用此模型提出了適用于水/己內酰胺/PA6混合物的相平衡模型,該模型適用的溫度和組成范圍較寬,可利用原有的NRTL模型數據庫,但是其二元交互作用參數通常由平衡數據回歸或者預測模型得到,其應用有一定的局限性。UNIFAC模型由Fredenslund等[7]提出,對于PA66聚合體系,由于文獻相平衡數據的缺乏,且系統中存在多個組分,實驗測定存在困難,采用數據回歸關聯系數的模型受到很大限制,而基于基團貢獻法的UNIFAC物性模型具有預測功能,處理復雜體系更具優勢。UNIFAC模型應用廣泛[8-11],選擇UNIFAC物性模型計算物性和組分相平衡。

UNIFAC模型提出活度系數作為多組分組成和溫度的函數方程,分為兩部分表示。

式中,Z為配位數。對于液體,Z = 10;純組分的參數ri和qi由基團的體積參數Rk和表面積參數Qk求得,而Rk和Qk由Bondi[12]給出的基團體積Vwk和表面積Awk計算得到。

式中,anm和amn為基團交互作用參數,anm≠amn,必須由相平衡的實驗數據求得,單位為K;T為體系的溫度,單位為K。

基團體積參數Rk和表面積參數Qk是UNIFAC物性模型的重要參數,與基團本身的大小和尺寸有關,其值取決于基團的劃分和參數值計算方法的可靠性。基團劃分時基團應盡可能小,但是劃分的基團又必須能夠影響分子的性質。Lohmann等[11]、孫曉波[13]和曾思[14]針對研究體系組分定義新的基團,修正后UNIFAC模型預測數據的相對誤差較小,預測結果較好。與以上經驗式的基團劃分不同,Wu 等[15]嘗試采用量子力學分子軌道計算軟件對分子中基團和原子的電荷進行檢測,根據所劃分的基團呈電中性和每個基團是最小單元的原則,從理論上給出了基團劃分的準則。Gong等[8]和Lü等[10]認為新劃分基團的Rk和Qk結構參數具有加和性,并且相似基團的結構參數相同,根據此假設可以由Bondi法[12]計算得到新定義基團的Rk和Qk;趙亞軍等[16]提出一種基于元素和化學鍵的快速計算分子體積和表面積的方法,將構成分子的元素和化學鍵作為貢獻單元,所計算分子體積和表面積參數的平均相對誤差分別為2.02%、2.84%,該方法簡單快速,應用范圍更廣。基團的劃分和準確合理的計算方法是獲得可靠的基團體積參數和表面積參數值的關鍵。

anm為UNIFAC物性模型的基團交互作用參數,由相平衡的實驗數據得到,其值取決于實驗數據是否可靠[17]。根據相平衡原理,忽略壓力和熱容差的影響,并用熔點和熔化焓代替三相點溫度及其焓變后,可得到溶質的溶解度方程[7]。

由式(12)可知,根據相平衡數據,只要有組分i的摩爾分數xi,結合純組分i的熔點及熔化焓,就可以計算活度系數γi,利用UNIFAC活度系數方程式(1)就可以回歸得到基團交互作用參數。

2 基礎數據的測定

2.1PA66鹽的溶解度測定

2.1.1實驗儀器及原料溶解釜,自制50 ml;CH1015型超級恒溫水浴槽;85-1型恒溫磁力攪拌器;精密溫度計,量程為0~50、50~100℃,精度為±0.05℃;AL204型電子天平,精度為±0.1 mg。實驗所用的原料PA66鹽為BASF公司提供產品,室溫下存放于干燥器中;己二酸和己二胺皆為市售,純度為分析純以上;水為實驗室自制去離子二次蒸餾水。

2.1.2實驗裝置及實驗步驟動態法測定溶解度的實驗裝置如圖1所示,溶解釜由玻璃制成,夾套中通循環水,由超級恒溫水浴控制其溫度,體系溫度由插入物料中的精密溫度計測定。溶解釜中的耐腐蝕磁力攪拌轉子連續攪拌以促進溶解并使固液兩相快速達到平衡,并且外接冷凝管以防止溶劑的揮發。

圖1 PA66鹽溶解度測定實驗裝置Fig.1 Solubility measurement apparatus of PA66-salt 1—thermostatic bath; 2—dissolution vessel; 3—thermometer; 4—condenser tube; 5—magnetic stirrer sub; 6—magnetic stirrer

用電子天平準確稱量一定量的溶質和一定量的去離子水加入到溶解釜中;依次開啟冷卻水、恒溫磁力攪拌器、超級恒溫水浴開關,充分攪拌固液兩相體系;控制超級恒溫水浴的升溫速率,用精密溫度計測定溶解釜中的溫度;當體系由渾濁變為澄清時,即為該溶質的溶解溫度。

2.1.3實驗方法和裝置的可靠性驗證為驗證實驗方法和裝置測定數據的可靠性,如圖2所示,測定己二酸在水中的溶解度,并與文獻值[18]比較。由圖2可見,在實驗溫度范圍內,測定值與文獻值吻合,該實驗方法和裝置測得的溶解度數據可靠。

圖2 己二酸在水中的溶解度實驗值與文獻值的比較Fig.2 Experimental solubility of adipic acid in water compared with data in literature

2.1.4PA66鹽在水中的溶解度數據測定的PA66鹽在水中的溶解度數據見表1,其中x為溶質的溶解度,以摩爾分數表示。

表1 PA66鹽在水中的溶解度實驗測定值Table 1 Solubility of PA66-salt in water

2.2PA66鹽的熔點和熔化焓測定

文獻[19-20]給出PA66鹽的熔點為195℃±5℃,而其熔化焓數據并沒有相關文獻報道。采用DSC測定PA66鹽的熔點和熔化焓,為UNIFAC模型參數回歸提供基礎數據。

取一定量的PA66鹽置于真空烘箱中干燥12 h冷卻后,用精度為0.01 mg的分析天平準確稱取2~3 mg的樣品制樣,在PE-DSC7型差示掃描量熱儀中以5℃·min?1的升溫速率從100℃升溫到350℃,其結果如圖3所示。

圖3 PA66鹽的熱流曲線Fig.3 Heat flux curve of PA66-salt sample

當溫度升至207℃時,PA66鹽開始聚合,不能得到完整的PA66鹽DSC熔化熱曲線。將PA66 鹽DSC熔化熱曲線近似看成是對稱的,采用如圖3所示方法,先積分195.8~206.1℃溫度區間,得熔化熱H1=157.89 J·g?1;再積分186.9~206.1℃溫度區間,得熔化熱H2=204.58 J·g?1;最后積分100.9~186.9℃溫度區間,得熔化熱H3=4.34 J·g?1。采用以上假設,PA66鹽的熔化焓H=H1+ 2(H3+H2?H1)= 259.97 J·g?1,摩爾質量MPA66-salt= 262.35 g·mol?1,轉化為摩爾熔化焓ΔHMfus= HMPA66-salt=68200.0 J·mol?1。PA66鹽的熔點TM=202.5℃=475.6 K。

3 UNIFAC模型及參數確定

3.1基團的劃分和定義

反應體系中組分基團的劃分是否合理是UNIFAC模型估算物性準確性的關鍵,需根據聚合體系中組分的結構進行基團的合理劃分。涉及的化合物和官能團及其結構式列于表2。

表2 PA66聚合過程中所涉及的化合物及其官能團Table 2 Components and function groups involved in PA66 production process

現有的基團不能準確定義PA66鹽的結構。將羧基和氨基形成的締合基團(—CH2COO?·+H3NCH2-)、羧基與其相鄰的亞甲基(—CH2COOH)分別定義為新的UNIFAC基團。

UNIFAC模型中已存在羧基(—COOH)和酰胺鍵基團(—CONHCH2—),考慮到羧基對相鄰亞甲基的影響,需要重新定義新的基團—CH2COOH;而新定義的羧基和氨基形成的締合基團(—CH2COO?·+H3NCH2—)則為SOL-SALT 和DIS-SALT的特征基團。

3.2表面積參數和體積參數的計算

體系中涉及的常規組分含有的基團數目,及所劃分基團的體積參數和表面積參數見表3。

新定義基團的體積參數Rk和面積參數Qk值可由結構參數具有加和性的性質近似計算。

以上涉及的基團、元素和化學鍵的體積參數和表面積參數值采用Bondi法[12]及元素和化學鍵法[16]中的文獻值計算。

3.3基團交互作用參數的回歸

UNIFAC模型參數回歸可以在Aspen軟件的Regression模塊中完成,輸入的PA66鹽溶解平衡方程式如下

參數回歸過程中需要定義體系中所需要的基團結構和組分基團數目,以及PA66鹽-水固液相平衡數據和PA66鹽的熔點和熔化焓數據,回歸得到基團交互作用參數。

基團參數Qk和Rk來自Aspen plus數據庫和自定義計算值,以及組分基團的劃分見表3;二元體系PA66鹽-水固液相平衡數據采用本實驗所測數據,見表1,其溫度范圍為305.45~363.05 K。在PA66鹽物性估算的基礎上回歸基團交互作用參數aij值,見表4。

表3 組分基團的劃分及其體積參數和表面積參數Table 3 Group identifications of components and group volume and area parameters

表4 UNIFAC模型的基團交互作用參數Table 4 Binary interaction parameters for UNIFAC model

4 模型驗證

4.1PA66鹽在水中的溶解度計算

為驗證模型及參數的準確性,對0~100℃范圍內PA66鹽在水中的溶解度進行模擬計算,并與實驗值進行比較,結果見表5,其平均相對誤差AADP 為2.1%。其中,AADP表示相對誤差百分數。

表5 PA66鹽-水溶解度模擬值與實驗值比較Table 5 Comparison of simulation data and experiment data of solubility of PA66-salt and H2O

4.2濃縮槽內鹽濃度的模擬

由PA66的聚合工藝過程知,PA66聚合過程涉及汽液平衡。為了驗證由PA66鹽-水固液平衡數據回歸參數后的UNIFAC模型能否應用于PA66聚合體系,考慮到鹽濃縮過程涉及PA66鹽-水的汽液平衡,模擬不同工業操作條件下的鹽濃度,并與工廠數據對比,其結果見表6。

表6 工業操作條件下鹽濃度計算值與測量值的比較Table 6 Comparison between calculation data and measurement data of salt concentration under industrial conditions

表6中,第2列為工廠濃縮槽出口鹽濃度數據;第3列為將PA66鹽重新劃分基團,由固液平衡數據回歸參數后UNIFAC物性模型的模擬結果;第5列為未考慮PA66鹽的離子絡合結構,采用數據庫中的基團和參數的模擬結果。采用Aspen中默認基團的UNIFAC模型計算結果誤差為20%左右,而采用劃分新基團回歸參數后的UNIFAC模型,計算結果相對誤差大大降低,其值小于5%。

5 結 論

PA66聚合體系為極性體系,選擇活度系數模型;由于文獻相平衡數據的缺乏,且系統中存在多個組分,實驗測定困難,基于基團貢獻法的UNIFAC物性模型具有預測功能,其更具吸引力,選擇UNIFAC活度系數模型計算物性和尼龍66鹽的溶解度。

UNIFAC模型中包含基團體積參數、表面積參數和基團交互作用參數3個參數。基團體積參數和表面積參數值與分子大小和形狀有關,其值取決于基團的劃分;基團交互作用參數由相平衡數據回歸得到。

實驗測定回歸參數所需的PA66鹽-水的固液平衡數據和PA66鹽的熔點及熔化焓;對體系組分進行基團劃分,定義的新基團為羧基和氨基形成的締合基團(—CH2COO?·+H3NCH2—)、羧基與其相鄰的亞甲基(—CH2COOH),計算新定義基團的體積參數和表面積參數;回歸UNIFAC模型的基團交互作用參數。

用回歸參數后的UNIFAC物性模型模擬0~100℃范圍內PA66鹽-水的固液平衡溶解度,與實驗值的平均相對誤差為2.1%;采用該模型進一步計算了工業操作條件下的鹽濃度,其相對誤差為5%以下。

符號說明

AADP ——相對誤差百分數

a ——基團交互作用參數,K

H ——熔化焓或熔化熱,J·g?1

ΔH ——摩爾熔化焓,J·mol?1

M ——摩爾質量,g·mol?1

Qk——基團的面積參數

q ——純組分的面積參數

R ——氣體常數,8.314 J·mol?1·K?1

Rk——基團的體積參數

r ——純組分的體積參數

T ——熱力學溫度,K

X——液相中基團的摩爾分數

x——液相中組分的摩爾分數

Z——配位數,對液體取為10

α——液相中組分的活度

Г——基團的剩余活度系數

γ——液相中組分的活度系數

θ——組分的面積分數

ν——基團數目

φ——組分的體積參數

上角標

c——結合項

exp——實驗

fus——熔化

(i)——組分i

r——剩余項

sim——模擬

下角標

i, j——組分

k, m, n——組分的基團

M——熔點

1,2,3——分別為195.8~206.1℃、186.9~

206.1℃、100.9~186.9℃區間

References

[1] 華陽, 劉振明, 劉權毅, 等. 尼龍66國內外生產現狀及發展建議[J]. 彈性體, 2010, 20(6): 78-82. HUA Y , LIU Z M, LIU Q Y, et al. Present situation and development suggestion of nylon-66 [J]. China Elastomerics, 2010, 20(6): 78-82.

[2] 顧凱, 黃繼紅. 聚合過程模擬與優化——基于Polymer Plus[M].北京: 化學工業出版社, 2010: 125. GU K, HUANG J H . Simulation and Optimization of Polymerization Process-Based on Polymer Plus [M]. Beijing: Chemical Industry Press, 2010: 125.

[3] SCHAFFER M A, MCAULEY K B, CUNNINGHAM M F, et al. Experimental study and modeling of nylon polycondensation in the melt phase [J]. Industrial & Engineering Chemistry Research, 2003, 42(13): 2946-2959.

[4] CHEN C C. A segment-based local composition model for the Gibbs energy of polymer solutions [J]. Fluid Phase Equilibria, 1993, 83: 301-312.

[5] RENON H, PRAUSNITZ J M. Local compositions in thermodynamic excess functions for liquid mixtures [J]. AIChE Journal, 1968, 14(1): 135-144.

[6] SEAVEY K C, KHARE N P, LIU Y, et al. A new phase-equilibrium model for simulating industrial nylon-6 production trains [J]. Industrial & Engineering Chemistry Research, 2003, 42(17): 3900-3913.

[7] FREDENSLUND A, JONES R L, PRAUSNITZ J M. Group-contribution estimation of activity coefficients in nonideal liquid mixtures [J]. AIChE Journal, 1975, 21(6): 1086-1099.

[8] GONG X, Lü Y, ZHANG Y, et al. Liquid-liquid equilibria of the quaternary system water + caprolactam + 1-octanol + ammoniumsulfate [J]. Journal of Chemical & Engineering Data, 2007, 52(3): 851-855.

[9] GMEHLING J G, ANDERSON T F, PRAUSNITZ J M. Solid-liquid equilibria using UNIFAC [J]. Industrial & Engineering Chemistry Fundamentals, 1978, 17(4): 269-273.

[10] Lü Y C, LIN Q, LUO G S, et al. Solubility of berberine chloride in various solvents [J]. Journal of Chemical & Engineering Data, 2006, 51(2): 642-644.

[11] LOHMANN J, R?PKE T, GMEHLING J. Solid-liquid equilibria of several binary systems with organic compounds [J]. Journal of Chemical & Engineering Data, 1998, 43(5): 856-860.

[12] BONDI A. van der Waals volumes and radii [J]. The Journal of Physical Chemistry, 1964, 68(3): 441-451.

[13] 孫曉波. 混合二元酸綜合利用工程基礎研究[D]. 鄭州: 鄭州大學, 2007. SUN X B. Study on engineering foundation of comprehensive utilization of mixed dibisic acids [D]. Zhengzhou: Zhengzhou University, 2007.

[14] 曾思. 工業聚酯裝置酯化過程的建模與流程模擬[D]. 杭州: 浙江大學, 2010. ZENG S. Modeling and simulation of esterification section in industrial poly(ethylene terephthalate) process [D]. Hangzhou: Zhejiang University, 2010.

[15] WU H S, SANDLER S I. Use of ab initio quantum mechanics calculations in group contribution methods (Ⅰ): Theory and the basis for group identifications [J]. Industrial & Engineering Chemistry Research, 1991, 30(5): 881-889.

[16] 趙亞軍, 司繼林, 夏力, 等. 基于元素和化學鍵的快速計算分子體積和表面積的方法 [J]. 計算機與應用化學, 2013, 30(7): 739-742. ZHAO Y J, SI J L , XIA L, et al. A new method based on elements and chemical bonds for fast calculation of moleculer volume and surface [J]. Computers and Applied Chemistry, 2013, 30(7): 739-742.

[17] 陳明鳴. 對苯二甲酸生產中相關相平衡的測定與研究[D]. 天津:天津大學, 2003. CHEN M M. Determination and study on the corresponding phase equilibria in terephthalic acid manufacture [D]. Tianjin: Tianjin University, 2003.

[18] STEPHEN H, STEPHEN T. Solubilities of Inorganic and Organic Compounds [M]. Oxford: Pergamon Press, 1963.

[19] 張朝暉. 國內外尼龍66鹽市場分析及需求預測[J]. 化工技術經濟, 2002, 20(4): 26-28. ZHANG Z H. Market analysis and forecast of nylon66 salt [J]. Chemical Techno-Economics, 2002, 20(4): 26-28.

[20] 己二酸、己二胺、尼龍66鹽[J]. 合成纖維工業, 1983, (3): 60-63. Adipic acid, hexamethylene diamine and PA66-salt [J]. China Syntehtic Fiber Industry, 1983, (3): 60-63.

DOI:10.11949/j.issn.0438-1157.20151235

中圖分類號:TQ 028.8

文獻標志碼:A

文章編號:0438—1157(2016)02—0435—07

Corresponding author:Prof. FENG Lianfang, fenglf@zju.edu.cn

Thermodynamic modeling with new UNIFAC groups for solubility of nylon66-salt in water system

GU Xueping, TIAN Lulu, FENG Lianfang, ZHANG Cailiang
(State Key Laboratory of Chemical Engineering, College of Chemical Engineering and Biological Engineering, Zhejiang University, Hangzhou 310027, Zhejiang, China)

Abstract:Modeling of PA66 polymerization process requires accurate solubility computation of nylon66-salt in water system. The UNIFAC activity coefficient method with new groups was applied to predict the solubility. The new groups were —CH2COO?·+H3NCH2— and —CH2COOH, which characterized the component’s special structure. The pure parameters of new UNIFAC groups and interaction parameters between new groups were obtained by regression with melting point and molar fusion enthalpy of PA66-salt as well as the solubility data of PA66-salt in water obtained by experiments. The UNIFAC model with new groups was used to simulate the solubility of PA66-salt in aqueous solution, giving a relative error of 2.1% from the experimental data. The model was further employed to simulate the concentration of PA66-salt under the industrial conditions (120℃) and gave a relative error as low as 5%, compared with 20% from the UNIFAC model with default groups in Aspen.

Key words:nylon66-salt; aqueous solution; solubility; UNIFAC model; UNIFAC groups

主站蜘蛛池模板: 午夜视频免费试看| 一本二本三本不卡无码| 国产乱子伦视频三区| 日韩专区第一页| 欧美人人干| 国产视频只有无码精品| 亚洲日韩高清在线亚洲专区| 亚洲无码精彩视频在线观看| 岛国精品一区免费视频在线观看| 国产三级精品三级在线观看| 欧美中文字幕一区| 亚洲综合中文字幕国产精品欧美| 综1合AV在线播放| 亚洲IV视频免费在线光看| 国产一区二区三区精品久久呦| 亚洲成人一区二区| 一级爆乳无码av| 中文字幕亚洲乱码熟女1区2区| 一级爆乳无码av| 欧美成一级| 91精品视频网站| 中文无码精品a∨在线观看| 日韩在线网址| 午夜国产理论| 色综合久久久久8天国| 亚洲毛片一级带毛片基地| 国产在线观看精品| 久草视频精品| 色老头综合网| 嫩草国产在线| 麻豆精选在线| 亚洲欧美日韩色图| 亚洲精品动漫| 91无码人妻精品一区| 免费看a级毛片| 天堂在线亚洲| 久久免费精品琪琪| 亚洲天堂视频网站| 在线播放国产99re| 日日拍夜夜嗷嗷叫国产| 五月婷婷丁香综合| 99伊人精品| 日韩人妻无码制服丝袜视频| 特级做a爰片毛片免费69| 亚洲av无码人妻| 亚洲国产综合精品一区| 色欲色欲久久综合网| 国产乱子伦无码精品小说| 午夜限制老子影院888| 亚洲男人的天堂久久香蕉网| 国产熟睡乱子伦视频网站| 2021国产精品自产拍在线观看| 狠狠综合久久久久综| 国内精品免费| 日韩欧美色综合| 国产无人区一区二区三区| 国产免费精彩视频| 亚洲女同一区二区| 日韩精品一区二区三区swag| 四虎永久在线精品国产免费| 色综合天天娱乐综合网| 国产极品粉嫩小泬免费看| 国产一级在线观看www色| 日韩高清中文字幕| 久久精品中文无码资源站| 欧洲极品无码一区二区三区| 四虎永久在线精品影院| 又大又硬又爽免费视频| 一级毛片中文字幕| 亚洲午夜福利精品无码不卡| 免费人成视网站在线不卡| 欧美日韩精品综合在线一区| 四虎永久在线| 免费人欧美成又黄又爽的视频| 内射人妻无码色AV天堂| 无码一区18禁| 亚洲侵犯无码网址在线观看| 精品伊人久久久大香线蕉欧美 | 91精品久久久久久无码人妻| 精品欧美日韩国产日漫一区不卡| 国产一区二区三区在线无码| 国产产在线精品亚洲aavv|