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

大氣中苯甲酸參與硫酸-氨的成核機制研究

2025-07-01 00:00:00陳沙沙衛書琴袁佳怡張夢潔王睿張田雷
關鍵詞:生長體系

中圖分類號:064 文獻標識碼:A 文章編號:2096-3998(2025)03-0034-09

作為大氣氣溶膠的重要來源,新粒子形成(NPF)不僅是全球云凝結核數量的主要貢獻者[1],而且對云的形成和氣候系統產生重要影響[2]。此外,NPF 是霧霾事件的起源,對人體健康會產生嚴重的損害。基于此,深人理解NPF機制對于準確評估其對氣候和人體健康的影響具有至關重要的意義。

硫酸(SA)作為大氣中最重要的成核前驅體,其在NPF 過程中扮演著重要的角色[3]。先前的研究[45]表明,SA驅動的NPF對大氣顆粒物濃度的貢獻是不可忽視的。氨(A)具有穩定SA團簇的能力[6],能夠形成穩定的 SA-A二元團簇。然而,傳統的二元成核理論并不足以解釋 NPF 事件的發生[7]。Zhao 等[8發現三元體系比二元體系具有更強的氫鍵和更多的質子轉移來穩定團簇。雖然目前已對 SA-A驅動的三元成核進行了大量的研究,但仍然無法解釋外場觀測到的高NPF速率。

因有機酸能夠與 SA-A形成氫鍵增強穩定性,其參與 SA-A 成核過程的研究備受關注[9]。作為大氣中主要的人為酸之一,苯甲酸( C6H5COOH ,BA)在城市、工業、農業及植被覆蓋地區大量存在,主要由汽車尾氣排放、工業生產、農藥噴灑、生物質燃燒、植物代謝和芳香化合物光化學反應等產生[12]。據報道,苯、甲苯和鄰苯二甲酸酯等芳香化合物能夠存在于對流層[13]。而在日本富士山和北海道、孟加拉灣的西海岸以及中國的北京、天津、新疆等地區的氣溶膠中均檢測出 BA[14-15] ,其在大氣中的濃度可達 1010 molecules ?cm-3 。Zhang等[16在BA參與SA和 H2O 團簇形成過程中,發現BA能夠顯著提高 SA-H2O 氣溶膠成核的形成速率。Wang 等[17-I8]采用密度泛函理論研究BA與SA/A/DMA(二甲胺)等成核前驅體在大氣中的相互作用,發現當BA濃度較高時, 團簇的濃度可達到與(SA) 2 團簇相同量級或更高,這在一定程度上揭示了BA具有促進NPF的重要作用。然而,BA促進SA-A成核機制目前仍不清楚。

本文采用團簇構型搜索、密度泛函理論結合大氣團簇動力學模型,對 z?3 )體系在不同溫度和成核前體濃度下的成核過程進行模擬。具體來說:(i)闡明不同溫度和各單體濃度條件下BA對SA-A體系形成速率的增強作用;(ii)追蹤( SA)x(A)y(BA)z 體系的主要生長路徑并探討BA的參與機制。本研究可為消除大氣二次氣溶膠污染、改善大氣環境等提供理論指導。

1 計算方法

1.1 構型搜索

采用ABCluster程序結合M06-2X方法對 (SA)x(A)y(BA)z ( 0?y?x+z?3 )體系的構型進行搜索和優化。本文選擇M06-2X方法是由于該方法在團簇優化中能夠準確描述分子間非共價相互作用,并且能夠準確預測團簇的熱力學性質。具體而言,首先,采用ABCluster程序[19]和CHARMM36力場中的人工蜂群算法得到 ( 0?y?x+z?3 )體系的初始構型 n×1000(0

1.2大氣團簇動力學模型

為了闡明BA在大氣顆粒形成中的作用,利用大氣團簇動力學模型(ACDC)[25]研究了不同溫度和單體濃度對 SA-A-BA 團簇的影響。利用Matlab- ?R2014a[26] 中提供的ode15s算法,通過數值求解生-滅方程來計算團簇濃度的時間變化,生-滅方程式的具體細節如下:

其中, ci 表示團簇 i 的濃度,團簇 i 與團簇 j 的碰撞系數用 βi,. 表示, γ(i+j)i 表示蒸發系數,描述了組合的團簇 i+j 分解成團簇 i 和團簇 j 的過程, Qi 包含團簇 i 的所有附加源項, Si 表示團簇 i 的另一個潛在項。此外,外部損失通過使用固定的冷凝匯(CS)值 2.0×10-2 進行估算[27]。

2 結果與討論

2.1 靜電勢分析

為了研究SA、A和BA分子之間的相互作用,對其單體的靜電勢(ESP)進行了計算,分析了其潛在的相互作用位點。

如圖1所示,BA分子中-OH基團周圍的紅色區域表明存在正的ESP值 (?+53.10kcal?mol-1 ,說明-OH基團為有效的氫鍵供體位點,可與SA分子中的末端O原子和A分子中的N原子形成氫鍵。與此同時,BA分子中末端O原子周圍的藍色區域表明存在負的ESP值 ),可以作為有效的氫鍵受體位點,能夠與SA分子中的H原子和A分子中的H原子通過氫鍵結合。在這里以(SA)1?(A)1?(BA)1 團簇為例,發現BA在與SA和A分子結合時,主要以氫鍵相互作用為主,而SA與A結合時不僅存在氫鍵相互作用還發生了質子轉移。此外,在 (SA)x(A)y(BA)z 團簇中也存在著類似的相互作用。

2.2 熱力學穩定性分析

為得到 (SA)x(A)y(BA)z 團簇最穩定的結構,采用ABCluster程序進行系統的結構采樣,并通過一系列量子化學計算對結構進行逐步優化。最終在 M06-2X/6-311++G(2df,2pd) 水平上,獲得了27個(SA)x(A)y(BA)z 團簇的穩定構型。圖2展示了16個具有代表性的含BA團簇的最穩定結構。由圖可發現,隨著BA的數目的增加, (SA)x(A)y(BA)z 團簇中氫鍵的數量和環的數量增多,使得團簇的穩定性增強。

圖1SA、A和BA單體分子的靜電勢圖

氫鍵的長度單位為A;顏色代碼:藍色 Σ=Σ 氮,黃色 Σ=Σ 硫,紅色 Σ=Σ 氧,灰色 Σ=Σ 碳,白色 σ=σ 氫圖2 M06-2X/6-311++G(2df,2pd) 水平上 (SA)x(A)y(BA)z L 0?y?x+z?3, 團簇的最穩定結構

為了評估 (SA)x(A)y(BA)z 團簇的熱力學穩定性,研究了不同溫度(218.15、238.15、258.15、278.15、298.15K) 下團簇的形成Gibbs自由能 (ΔG,kcal?mol-1 )。結果表明,隨著溫度降低團簇的 ΔG 也逐漸降低,表明在低溫下團簇的熱力學穩定性更好。基于此,在圖3中列出了在 218.15K 和1atm條件下團簇形成過程的吉布斯自由能圖,可以看出,除了在 (BA)2 團簇中加入BA分子的過程是吸收能量的 (0.96kcal?mol-1 ),其余團簇在逐步形成過程中均為自發過程。值得注意的是,在除了團簇( 之外的其他團簇中加人SA分子會比加入BA分子更穩定,這可能是由于SA的酸性比BA更強。與BA相比,團簇與SA的結合更穩定, (SA)3?(A)3 團簇的 ΔG 值為 -66.75kcal?mol-1 ,在相應的二元SA-A和含BA團簇中最低。因此,選擇 (SA)4?(A)3 和 (SA)4?(A)4 團簇作為ACDC動力學模擬的邊界條件。

圖3 218.15K 和1 atm F(SA)x(A)y(BA)z(0-1 )圖

2.3 溫度與前驅體濃度的影響分析

為了進一步確定BA對SA-A成核過程的影響,通過一系列ACDC 模擬了不同溫度和前驅體濃度條件下有、無BA參與SA-A團簇形成速率。在此通過BA對SA-A體系形成的促進作用定義了增強強度(202 (R) 。 R 的計算方法為 (SA)x(A)y(BA)z(0-3?s-1, 與(SA)x(A)y(0

其中, x 和 y 分別表示大氣中[SA]和[A],范圍分別為 106~108 molecules ·cm-3 和 107~1011 molecules·cm-3 z 表示[BA],范圍為 106~1010 molecules ?cm-3

圖4(a)、(b)分別研究了 [SA]=107 molecules·cm-3,[A]=109 molecules·cm-3 和不同溫度(218.15,238.15,258.15,278.15,298.15K) 條件下, J 和 R 隨[BA]的變化情況。從圖4(a)中可以看出,SA-A-BA體系的 J 始終大于或等于SA-A體系的 J ,并且 J 隨著溫度的降低和[BA]的增加而升高。具體而言,在 218.15K 時, [BA]=1010 molecules ?cm-3 時的 J 比 [BA]=106 molecules·cm-3 時的 J 增加了2個數量級,這證明了BA在低溫和高濃度時對SA-A團簇的形成起著顯著的增強作用。此外,如圖4(b)所示,當 [BA]lt;108 molecules ?cm-3 時,在不同溫度下 RBA 基本保持不變, [BA]gt;108 molecules ? (20cm-3 時, RBA 在低溫( 218.15K 和 238.15K )下開始增加,而對于相對高溫 (258.15,278.15,298.15K) 時,只有當 molecules ?cm-3 時 RBA 才開始有微弱的增加。因此,溫度和[BA]對SA-A-BA體系的 J 和 R 都有著顯著的影響。

除溫度和[BA]之外,[SA]和[A]也是影響團簇 J 和 R 不可忽視的因素。圖4(c)、(d)示出了最優溫度(218.15K)和[BA]( 1010 molecules ?cm-3 )條件下, J 和 R 隨[SA]或[A]的變化情況。由圖4(c)可知,在 218.15K 時由于 J 與前驅體的濃度成正比, J 隨著[SA]和[A]的增加而升高。在238.15、258.15、278.15、298.15K時 J 與[SA]和[A]的變化規律與 218.15K 類似。與 J 隨[SA]和[A]變化趨勢不同,如圖4(d)所示 RBA 隨著[SA]和[A]的降低而升高,呈現出負相關效應。這可能是由于SA分子具有較強的酸性進而更容易與堿性分子A結合。上述分析表明,在低溫、高[BA」、低[SA」和[A]條件下,BA對SA-A體系的增強作用最為明顯。

(a)不同溫度下形成速率隨[BA]的變化
(b)不同溫度下增強強度隨[BA]的變化
(d)增強強度隨不同[SA]和[A]的變化
(c)形成速率隨不同[SA]和[A]的變化圖4BA對SA-A體系形成速率和增強強度的影響

2.4 團簇生長路徑分析

為進一步確定BA在SA-A-BA體系中的作用,追蹤了不同溫度下團簇形成的主要生長路徑。

如圖5所示,團簇生長路徑包含兩種途徑,一種是黑色箭頭表示的純SA-A團簇生長路徑,其中(SA)2?A 團簇的形成是第一步,有兩種形成途徑, SA?SA?A?(SA)2?A 和 SA?(SA)2?(SA)2 :A。而隨后團簇主要是通過與SA、A單體或與 SA?A 團簇碰撞形成,例如, (SA)2?A?(SA)3?A? (SA)3?(A)2?(SA)3?(A)3 和 (SA)2?A?(SA)3?(A)2?(SA)3?(A)3 。另一種是黃色箭頭表示的BA參與的生長路徑,在成核的初始階段由 和 SA?(SA)2?(SA)2?BA 路徑形成 (SA)2?BA 團簇。隨后 (SA)2?BA 團簇會與A碰撞形成 (SA)2?A?BA 團簇,由于BA的蒸發系數較高,能夠從含BA團簇中蒸發,使得BA成為促進SA-A團簇生長的“催化劑”。 SA?A?SA?A?A? (BA) )2?SA?(BA)2?SA?BA?SA?A?BA?SA?A 和 (SA)2?A?(SA)3?A?BA?(SA)3?A 過程也證明 BA對SA-A團簇的形成具有“催化”作用。通過對團簇生長路徑的進一步分析發現,BA單體或含BA的團簇可以通過與已有的SA-A團簇碰撞,直接“參與”團簇的形成 ((SA)3?(A)3?(SA)4?(A)3 ·(BA)2- Fluxout和 (SA)3?(A)3?(SA)4?(A)3?BA?Fluxout) ,并最終從模擬系統中生長出來。此外,在不同溫度(238.15、258.15、278.15、298.15K)條件下團簇的生長路徑有所不同,在238.15K時,BA在團簇的形成過程中只起到“催化劑”的作用,而在258.15~298.15K時,團簇的生長路徑為純 SA-A團簇。

注:黑色箭頭和黃色箭頭分別代表純SA-A團簇和含BA團簇的生長路徑

2.5 貢獻率分析

為了確定BA對SA-A團簇形成的貢獻率,研究了不同溫度和前驅體濃度對團簇聚集過程中生長途徑貢獻的影響。如圖6(a)所示,在相對高溫時 (258.15,278.15,298.15K) ,團簇中只存在純SA-A的生長路徑,這可能是由于在高溫條件下BA的蒸發速率比低溫下更高,導致其無法直接參與成核。當 T= 238.15K時,BA能夠作為“催化劑”催化SA-A團簇的形成,并未直接參與臨界團簇的形成過程,使得團簇的生長路徑仍然是以純SA-A團簇的形成為主。當溫度降低至 218.15K 時,由于BA的蒸發速率降低,使得BA能夠直接參與到團簇的形成過程中,并且貢獻率可達 52% 。由此得出,含BA參與的生長路徑與溫度呈負相關,這與圖4(b)所示的 RBA 與溫度的變化規律一致。如圖6(b)所示,在 [BA]= 106~108 molecules ?cm-3 時,只存在純 SA-A團簇生長路徑,在 [BA]=109 molecules ?cm-3 時,BA作為“參與者”參與SA-A團簇的形成中且對團簇形成的貢獻率達到 10% 。隨著[BA]升高至 1010 molecules ? (204號cm-3 時,含BA分子的生長路徑對團簇形成的貢獻逐漸增加至 52% ,相反,純SA-A團簇的形成路徑的貢獻從 88% 降低至 39% 。值得注意的是,由于SA與BA之間的競爭關系,含BA團簇對生長路徑的貢獻率隨著[SA]的降低而升高,從 41% 增加至 63% ,如圖6(c)。由于A與 RBA 之間呈負相關,含BA團簇對生長路徑的貢獻率隨著[A]的降低而升高。值得注意的是,在重污染區域的上空,[SA]與[A]隨著高度的增加而逐漸降低[28-30],該趨勢與BA 促進 SA 和A 團簇形成的貢獻率的結果一致。因此,本研究預測,在如中國東部、日本、德國及朝鮮南部等工業發達和人口密集的城市中,BA促進SA-A成核可能在秋末和冬季占主導地位。

圖5體系中團簇生長的主要途徑
圖6純SA-A團簇路徑和含BA路徑對系統通量的相對貢獻

3結論

本文通過采用團簇構型搜索、密度泛函理論結合大氣團簇動力學模型研究了BA參與SA-A成核的NPF機制,得到如下結論。

1)BA能夠通過氫鍵與SA和A分子結合并且還能促進SA-A發生質子轉移來穩定團簇。

2)在低溫、高[BA]、低[SA]和[A]條件下BA對團簇形成速率的增強最強,其中SA-A-BA體系的形成速率比純SA-A體系升高2個數量級。

3)在團簇形成過程中,BA既能作為“催化劑”促進SA-A團簇的形成,又能作為“參與者”直接參與SA-A 團簇的形成。與此同時,BA對團簇形成的參與程度隨著溫度的降低、[BA]的增加、[SA]和[A]的減少而增加。

[參考文獻]

[1]CAI J,SULOJ,GUY,etal.Elucidating the mechanismsofatmospheric new particle formation inthehighlypolutedPo Valley,Italy[J]. Atmospheric Chemistry and Physics,2024,24(4):2423-2441.

[2]BLICHNER SM,YLI-JUUTIT,MIELONENT,etal.Processevaluationofforestaerosol-cloud-climate fedbackshows clear evidence from observations and large uncertainty in models[J].Nature Communications,2O24,15:969. ence,2023,382:1308-1314.

[4]DAIL,ZHAOY,ZHANG L,et al.Particle number size distributions andformationand growth ratesof dferent newparticle formation types of a megacity in China[J].Journal of Environmental Sciences,2O23,131:11-25.

[5]YINR,YANC,CAI R,et al.Acid-baseclusters during atmospheric new particle formation inurban Beijing[J].Environmental Science amp; Technology,2021,55:10994-11005.

[6]LIY,FUTM,YUJZ,etal.Dissctingthecontributions oforganicnitrogen aerosols to globalatmosphericnitrogendeposition and implications for ecosystems[J].National Science Review,2023,10(12):nwad244.

[7]LIUL,YUFQ,DUL,et al. Rapid sulfuric acid-dimethylamine nucleation enhanced by nitricacid in polluted regions[J]. Proceedings of the National Academy of Sciences of the United States of America,2021,118 :e2108384118.

[8]ZHAOY,LIUYR,JIANG S,et al. Volatileorganiccompounds enhancing sulfuric acid-based ternaryhomogeneous nucleation:the important role of synergistic effect[J].Atmospheric Environment,2020,233:117609-117619.

[9] FANG X,HUM,SHANG DJ,etal.Observational evidence for the involvementof dicarboxylicacids in particle nucleation[J].Environmental Science amp; Technology Lettrs,2020,7:388-394.

[10]TCHINDA N T,DU L,LIU L,et al. Pyruvic acid,an effcient catalyst in SO 3 hydrolysis and effective clustering agent in sulfuric-acid-based new particle formation[J]. Atmospheric Chemistry and Physics,2O22,22:1951-1963.

[11]YANG Y,LIU L,WANG H X,et al. Molecular-Scale Mechanism of Sequential Reaction of Oxalic Acid with SO 3 : Potential Participator inAtmospheric Aerosol Nucleation[J]. Journal of Physical ChemistryA,2021,125:4200-4208.

[12]TIANY,IU X,HUOR,etal.Organiccompoundsource profilesof PM2.5fromtraffc emisions,coalcombustnindustrial processes and dust[J]. Chemosphere,2021,278:130429.

[13]IYER S,KUMAR A,SAVOLAINEN A,et al. Molecular rearrngementof bicyclic peroxy radicals is a keyroute to aerosol fromaromatics[J].Nature Communications,2O23,14:4984.

[14]KUMARI VR,SARMA VVS S,KUMAR M D.Spatial variabilityinaerosol compositionand itsseawater acidification potential in coastal watersof the western coastal BayofBengalJ].Journal ofEarth System Science,2O22,131(4):1-14.

[15]HO K F,HUANG RJ,KAWAMURA K,et al. Dicarboxylic acids,ketocarboxylic acids, α -dicarbonyls,fatty acids and benzoic acid inPM2.5aerosol coleced during CAREBeijing-2007:an effctof traffc restrictiononairquality[J].Atmospheric Chemistry and Physics,2015,15(6) :3111-3123.

[16] ZHANG R,SUHI,ZHAO J,et al. Atmospheric new particle formation enhanced by organic acids[J].Science,2004, 304(5676) :1487-1490.

[17]YAZDANI KACHOEI M,JALALI ASADABADI S. Thermoelectric properties of heavy fermion CeRhIn5 using density functional theory combined with semiclassical Boltzmann theory[J].RSC Advances,2019,9(62):36182-36197.

[18] WANG H,WANG X,HUY,et al.The efect of water molecules and air humidityon the cluster formation from benzoic acid withsulfuricacid/ammonia/dimethylamine:Amolecular-scale investigation[J].Journalof Molecular Liquids,2023, 390:123001.

[19]ZHANG J,DOLG M.ABCluster:the artificial beecolony algorithm for cluster global optimization[J].Physical Chemistry Chemical Physics,2015,17(37) :24173-24181.

[20] ELM J,MYLLYS N,HYTTINEN N,et al. Computational studyof theclustering of a cyclohexene autoxidation product (2號 C6H8O7 with itself and sulfuric acid[J]. Journal of Physical Chemistry A,2015,119(30):8414-8421.

[21]UNALY,NASSIF W,OZAYDIN B C,etal.Scale factordatabase for the vibration frequencies calculated in MO6-2X,one of the DFT methods[J]. Vibrational Spectroscopy,2021,112:103189.

[22] THIRUVANGOTH S,THAYYIL M S,NEENU KRISHNA P U,et al. Unveiling the structural properties,spectral information,scavenging mechanism,andNLOresponses of2-(4-Hydroxyphenyl)Ethylisothiocyanate:A theoretical insight[J]. Next Research,2024,1:100019.

[23] FRISCH M,TRUCKS G,SCHLEGEL H,et al. Gaussian O9[M]. Wallingford CT:Gaussian Inc,2009.

[24] LU T,CHENF.Bond order analysis basedon the Laplacian of electron density in fuzzyoverlapspace[J]. Journal of Physical Chemistry A,2013,117 :3100-3108.

[25]MCGRATH MJ,OLENIUS T,ORTEGA I K,etal. Atmospheric cluster dynamics dode:a flexible method forsolutionof thebirth-death equations[J].Atmospheric Chemistry and Phvsics,2012,12(5):2345-2355.

[26] SHAMPINE LF,REICHELTMW.The Matlab obe Suite[J].SIAMJournal on Scientific Computing,1997,18(1):1-22.

[27] YAO L,GARMASH O,BIANCHI F,et al.Atmospheric new particle formation from sulfuric acid and amines in a Chinese megacity[J].Science,2018,361:278-281.

[28] ROLLINS AW,THORNBERRYTD,WATTSLA,etal.Theroleofsulfurdioxide instratosphericaerosolformation evaluated byusing insitu measurements inthetropical lower stratosphere[J].Geophysical Research Leters,2O17,44:4280- 4286.

[29]MURPHY DM,FROYDKD,SCHWARZ JP,etal.Observations of thechemical compositionof stratospheric aerosolparticles[J]. Quarterly Journal of the Royal Meteorological Society,2014,140:1269-1278.

[30] WU K,ZHUSP,LIUYM,etal.Modeling ammoniaand itsuptake bysecondaryorganicaerosolover china[J].Journal ofGeophysical Research:Atmospheres,2021,126:e2020JD034109.

[責任編輯:魏強]

Abstract:Organic acids are considered important components of atmospheric aerosols and play a significant role inthe new particle formation(NPF)process.As the simplest aromatic organic acid in the atmosphere,benzoic acid (BA)has been shown to promote NPF. However,its role and nucleation mechanism in sulfuric acid(SA) and ammonia(A)-driven NPF processs remain unclear.Herein,the nucleation mechanism of the SA-A-BA system was investigated using density functional theory combined with an atmospheric cluster dynamics code. The results show that BA can enhance the formation rate of SA-A clusters,particularly atlow temperatures(218.15K),and in regions with high [BA]( 1010 molecules·cm-3 ),low[SA] ( 106 molecules·cm-3 ),and low[A]( 107 molecules·cm-3 ). The NPF mechanism suggests that in cold regions,BA can directly participate in the formation of SA-A clusters as a“participant”. This finding contributes to the understanding of the role of BA in the SA-A nucleation process,particularly in industrialized and densely populated urban areas.

Key words: benzoic acid; formation rates;enhanced strength;nucleation mechanism

猜你喜歡
生長體系
碗蓮生長記
小讀者(2021年2期)2021-03-29 05:03:48
構建體系,舉一反三
共享出行不再“野蠻生長”
生長在哪里的啟示
華人時刊(2019年13期)2019-11-17 14:59:54
探索自由貿易賬戶體系創新應用
中國外匯(2019年17期)2019-11-16 09:31:14
野蠻生長
NBA特刊(2018年21期)2018-11-24 02:48:04
生長
文苑(2018年22期)2018-11-19 02:54:14
《生長在春天》
如何建立長期有效的培訓體系
現代企業(2015年1期)2015-02-28 18:43:18
“曲線運動”知識體系和方法指導
主站蜘蛛池模板: 97亚洲色综久久精品| 国产波多野结衣中文在线播放 | 亚洲人网站| 国产区在线看| 日韩欧美国产区| 久久综合色视频| 天天做天天爱夜夜爽毛片毛片| 亚洲第一区欧美国产综合 | 夜夜操国产| 久久精品国产精品青草app| 精品丝袜美腿国产一区| 亚洲人成人无码www| 天天爽免费视频| 毛片网站免费在线观看| 激情亚洲天堂| 国产国产人免费视频成18| 真实国产乱子伦高清| 综合色区亚洲熟妇在线| 中文字幕第4页| 亚洲一区网站| 国产视频入口| 亚洲成人一区在线| 日韩中文精品亚洲第三区| 国产欧美精品一区aⅴ影院| 欧美国产日韩在线| 狠狠干综合| 制服丝袜无码每日更新| 中文字幕av一区二区三区欲色| 日本免费新一区视频| 在线国产资源| 国产精品久线在线观看| 人妻一区二区三区无码精品一区 | 九九免费观看全部免费视频| 美女无遮挡免费网站| 天堂在线视频精品| 久操中文在线| 国产成人做受免费视频| 毛片久久久| 精品国产成人高清在线| 久久伊人操| 亚洲第一在线播放| 欧美激情视频二区| 高清免费毛片| 国产精品女熟高潮视频| 无码专区在线观看| 国产凹凸一区在线观看视频| 538国产在线| 一区二区理伦视频| 久草视频一区| 国产丝袜91| 日韩黄色在线| 亚洲欧美综合在线观看| 欧美亚洲综合免费精品高清在线观看 | 日韩精品一区二区三区大桥未久 | 精品无码一区二区三区电影| 国产乱人伦偷精品视频AAA| 精品成人一区二区| 国产精品成人AⅤ在线一二三四| 熟妇丰满人妻| 免费国产高清精品一区在线| 久久永久免费人妻精品| 老熟妇喷水一区二区三区| 精品一区二区三区波多野结衣| 99久久精品免费观看国产| 国产91久久久久久| 欧美啪啪网| 一级毛片在线直接观看| 欧美成人综合在线| 伊人查蕉在线观看国产精品| 久久情精品国产品免费| 免费jizz在线播放| 在线观看av永久| 97超级碰碰碰碰精品| 亚洲精品免费网站| 亚洲精品国产精品乱码不卞 | 精品伊人久久久香线蕉 | 国产成人综合日韩精品无码首页 | 国产精品妖精视频| 国产区精品高清在线观看| 97狠狠操| 69精品在线观看| 亚洲青涩在线|