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

矩形微通道內液相黏度對氣泡界面的作用機制

2023-08-16 09:06:44陳蔚陽宋欣殷亞然張先明朱春英付濤濤馬友光
化工進展 2023年7期
關鍵詞:平尾界面

陳蔚陽,宋欣,殷亞然,張先明,朱春英,付濤濤,馬友光

(1 浙江理工大學材料科學與工程學院,紡織纖維材料與加工技術國家地方聯合工程研究中心,浙江 杭州 310018;2 化學工程聯合國家重點實驗室,天津大學化工學院,天津 300072)

作為化工過程強化設備,微通道反應器憑借其提供的高濃度/熱梯度和巨大比表面積,已在生物化學、有機合成、分離、鹵化等涉及傳質受限的氣-液、液-液體系中展現潛在的應用前景[1-5]。闡釋微通道內流體動力學行為是實現微尺度受限空間內傳質強化和反應性能有效調控的重要基礎。

彈狀流流型因分布均勻、易于調控和重復性好等優點在諸多體系中受到青睞,而通道幾何尺寸、兩相流速和流體物性等因素卻影響彈狀流中氣泡/液滴的形狀與流動,因而探究氣泡/液滴流動與形變機理對于實現氣泡/液滴行為和熱/質傳遞強化的有效控制具有重要意義[6-10]。目前,已有學者針對彈狀液滴界面形變規律、液膜厚度和流場結構等內容展開研究[11-13]。研究發現,液滴形狀與液滴尺寸、毛細管數(Ca)和黏度比等緊密相關[14]。Kovalev等[15]利用Ca數、雷諾數(Re)和韋伯數(We)繪制液滴形狀分布圖,并揭示了不同液滴形狀間的轉變機制,表示連續相的高剪切速率和流體結構變化導致液滴由彈狀轉變為啞鈴狀。相較于液-液體系,氣-液體系中高的界面張力和液/氣黏度比將導致氣泡界面形變與流動規律顯著不同[16],而如何闡釋液相黏度對氣泡界面的作用機理成為相關研究的重點。Yao等[17]發現在Ca數小于0.07范圍內,隨著Ca的增大,氣泡在液相剪切作用下頭尾曲率增大,而增加氣液相We數,液彈壓力升高,導致氣泡尾部更加扁平甚至出現凹陷[18]。Fei等[19]在矩形微通道中發現彈狀、啞鈴狀和手榴彈狀三種氣泡形狀,并發現增加液相黏度氣泡變形的臨界流速降低。此外,相界面的表面活性劑在連續相強剪切下將遷移至氣泡/液滴尾部,產生馬蘭戈尼效應,引起尾部界面張力降低甚至破裂[20]。盡管以上研究指出了氣泡不規則界面及尾部破裂現象[20-21],而關于液相黏度對氣泡界面作用機制的認識尚不充分。Wu等[20]表明微通道內氣泡/液滴界面形變和破裂是一個時間過程,但對界面演變過程卻未做深入分析。

本文針對微通道中不同黏度體系下氣泡界面及演變規律進行實驗研究。利用高速攝像機追蹤并獲取氣泡形狀和長度等特征尺寸,重點闡釋氣泡界面形狀轉變機理。利用氣液流率比和Ca數繪制氣泡形狀流型分布圖,并建立流型轉變判別模型。根據氣泡形狀特征長度研究了氣泡不同形狀間的演變過程,并分析了操作條件對流動過程氣泡由棒槌狀向平尾狀和棒槌狀向尖尾狀轉變的影響。最后,結合氣液流率比和Ca數分析了尖尾狀氣泡的破裂機理與臨界條件,并提出尖尾破裂條件的預測模型。

1 材料和方法

1.1 實驗裝置

實驗所用微通道芯片由兩片聚甲基丙烯酸甲酯(PMMA)板制成,下板刻有微通道結構,上板為光滑平板,上下板通過螺栓進行密封。如圖1(a)所示,微通道為對流T形結構,進口分支和主通道的長度分別為5mm 和25mm,均為寬(w)0.6mm、深(h)0.4mm 的矩形截面。微通道的制造精度為±6μm,占微通道深度1.5%,根據已有研究[22-23],此粗糙度對流體流動幾乎不產生影響。實驗所用裝置及流程如圖1(b)所示,由兩臺微流注射泵(Harvard PHD 2000,USA,0.35%穩定精度)以設定的流量分別將氣液兩相送至T形微通道的分支進口。微流注射泵能夠精確控制氣液流速,并在T形結構處生產穩定的氣泡[氣泡長度(LB)誤差小于6%,如圖2 所示]。液相流率(QL)和氣相流率(QG)范圍分別設定為40~100mL/h 和20~500mL/h。氣液兩相在T形交叉口處交匯并產生氣泡,隨后通過主通道流入收集瓶中。待流型和壓力穩定后,在冷光燈照明下,采用高速攝像儀(FastCam SA1.1,Photron,Japan)以每秒2000 幀的拍攝頻率(精度為0.5μs)記錄微通道內氣液兩相的流動狀態。采用壓差傳感器(Honeywell ST3000,USA)測量主通道的進口壓力,出口壓力為大氣壓。所有實驗在(293.15±2)K和室內壓力(101kPa)下進行。

圖1 微通道結構和實驗裝置流程示意圖

圖2 氣泡初始長度的穩定性分析

1.2 材料與物理性質

實驗所用氮氣(N2,純度>99.99%),杭州今工特種氣體有限公司;甘油(甘油,純度>99%)和十二烷基硫酸鈉(SDS,純度>99%),上海麥克林生化科技有限公司,藥品使用前均未經任何處理。氣液兩相分別為氮氣和不同甘油含量(c,0、20%、40%、60%、80%、90%,質量分數)的水溶液。水溶液中加入質量分數0.4%的表面活性劑SDS以提高通道壁的潤濕性,產生穩定流型。實驗采用加熱方法,在甘油水溶液中溶解0.4%SDS粉末至溶液澄清,之后溶液放置在室內冷卻至20℃,并在使用中依舊澄清,因而SDS質量質量分數確定為0.4%。經接觸角測定儀(JY-82B,承德鼎盛試驗機測試設備有限公司)測量,不同甘油質量分數的水溶液與PMMA 通道芯片表面的接觸角均小于46°(圖3),表明溶液對材料表面具有很好潤濕性。水溶液的密度(ρL)、界面張力(σ)和黏度(μL)分別 由 密 度 計 (DMA-4500-M,Anton Paar,Austria)、表 面 張 力 儀(OCAH200,Dataphysics instruments GmbH,Germany)和 黏 度 儀(iVisc,LAUDA,Germany)測定。所有物理特性至少重復測量5 次,實驗數據采用平均值,誤差小于5%。相應物性數據見表1。

表1 甘油水溶液的物理性質

圖3 甘油水溶液與PMMA通道芯片表面的接觸角

2 結果與討論

2.1 氣泡形狀

實驗觀察表明,氣泡在微通道內運動過程中呈現四種氣泡形狀,分別為子彈狀[圖4(a)]、棒槌狀[圖4(b)]、平尾狀[圖4(c)]和尖尾狀[圖4(d)]。理論上,氣泡形狀是由液相壓力和黏性剪切力克服界面張力形成的。為闡明氣泡界面變形機理,對氣泡頭尾部壓差進行分析。根據Young-Laplace 方程,氣泡尾部和頭部的拉普拉斯壓降分別為式(1)、式(2)。

圖4 微通道中典型的氣泡形狀

式中,P1為氣泡上游液相的壓力;P2為氣泡尾部內部壓力;P3為氣泡頭部內部壓力;P4為氣泡下游液相壓力;σ為界面張力;Rrear和rrear分別為氣泡尾部在y軸和z軸方向上的界面曲率半徑;Rhead和rhead分別為氣泡頭部在y軸和z軸方向上的界面曲率半徑。研究表明,液滴內部壓力在不同流型下幾乎保持定值[24],而氣相黏度比液相黏度低2~4個數量級,因此氣泡內部壓力也可認為是均勻的,即P2=P3。氣泡相鄰兩個液彈的壓差可由式(3)進行表示。

(1)子彈狀氣泡[圖4(a)] 出現在低的氣液流率比和液相黏度條件,其形狀特征為具有穩定的弧形頭尾部,且尾部曲率明顯小于頭部,表明氣泡相鄰液彈內壓力P4<P1,氣泡和與通道壁間的液膜將從氣泡尾部流向頭部。由于壁面摩擦作用,通道中心氣泡比通道壁附近流體流動快,液膜內產生速度梯度,而氣泡將受到液相流體反向剪切,從而導致氣泡頭部曲率大于尾部[25]。由圖5(a)可知,隨著液相黏度的增加,頭部曲率半徑變小而尾部弧形趨于平緩。這是由于較大黏度增強矩形截面角區液膜的流動,進而增強對氣泡頭部側面的反向剪切作用[17]。因此,在受限空間內,氣泡通過改變頭尾曲率以提高頭尾壓差(P1-P4)來平衡增大的剪切應力[25]。

圖5 液相黏度對微通道中氣泡形狀的影響

(2)棒槌狀氣泡[圖4(b)] 隨著氣液流率比升高,氣泡頭部曲率幾乎保持不變,而氣泡尾部(a區域)出現凹陷[圖4(b)][26],尾部y軸方向上曲率逐漸減大,但始終保持Rrear>Rhead,即P4<P1[圖5(b)],因而凹陷處液膜可能來自氣泡尾部[27]。由于矩形通道深度小于寬度(h<w),氣泡受上下壁摩擦作用更為嚴重,因而凹陷部分液膜更傾向來自左右壁的角區。在較大氣液流率比下氣泡凹陷界面處在液體剪切作用被進一步拉長,從而形成棒槌狀[15]。隨著液相黏度的增加,液相黏性剪切力和液彈壓力增大,液體更易流向左右兩端而拉長氣泡,導致氣泡尾部附近液膜變厚、尾部曲率變大[圖5(b)][28]。

(3)平尾狀氣泡[圖4(c)] 進一步增大氣液流率比,氣泡頭部曲率仍維持不變,但尾部曲率較棒槌狀氣泡而減小,形狀為扁平或凹面狀,氣泡相鄰液彈內壓力仍滿足P4<P1。隨著氣/液流率比的增加,氣泡上游液彈壓力等于甚至高于氣泡內壓力(P1≥P2),在上游液彈擠壓作用下氣泡尾部曲面最終由凸變平甚至變凹[29]。如圖6(a)所示,平尾狀氣泡是由棒槌狀氣泡在流動過程演變而來。增加液相黏度,棒槌狀氣泡將提前向平尾狀轉變[圖5(c)]。氣泡由棒槌狀向平尾狀的演變過程將在第2.3 節展開討論。

圖6 微通道中氣泡形狀的演變過程

(4)尖尾狀氣泡[圖4(d)] 高的液相黏度和氣液流率比下,氣泡在被拉長的同時尾部曲率半徑銳減,其特征為具有穩定的尖尾。隨著液相黏度和氣液流率比的增加,y軸方向氣泡尾部曲率半徑減小,而z軸方向上氣泡尾部形狀由凸變凹。研究發現,方形截面通道中液滴或氣泡在靠近通道壁處形成2 處或4 處尖尾[20,29],而在矩形截面通道中則發現1處尖尾[15]。矩形微通道內,液彈在通道中心流速最大[30],且受壁面摩擦效應顯著,因而尖尾極有可能為2 個,形成于上下通道壁面中心[20,29][圖4(d)]。

此外,表面活性劑存在情況下,Marangoni 效應也可促進尖尾形成[31]。由于彈狀流中氣泡速度高于周圍流體平均速度,氣泡界面上的表面活性劑在壁面剪切效應下將從氣泡頭部剪切至尾部,導致氣泡界面產生張力梯度,而尾部界面張力最小[32],因而在較大壁面剪切和上游液彈中心擠壓下,氣泡尾部變成近壁面尖而中間凹的尖尾形狀[圖5(d)][33]。如圖6(b)所示,尖尾狀氣泡通常也由棒槌狀氣泡演變而來。這是因為氣泡在生成過程中,氣泡界面處表面活性劑轉移和吸附飽和時間要長于氣泡生成時間(亞秒級),剛生成氣泡界面的表面活性劑濃度梯度小,其尾部不容易發生變形,從而維持棒槌狀[34]。而隨著氣泡向下游流動,氣泡界面表面活性劑在剪切作用下由頭部向尾部轉移,導致氣泡界面活性劑濃度梯度增大,尾部界面張力最小,從而逐漸形成尖尾狀[29]。隨著液相黏度和兩相流率的增大,氣泡界面處增大的黏性剪切力將使表面活性劑的轉移和吸附速度加快[35-36],氣泡界面容易形成較大張力梯度,導致尖尾狀提前出現[圖5(d)],甚至出現尖尾破裂[圖5(e)]。

微尺度限域內氣泡與其周圍液相流場相互耦合,因而氣泡形狀將密切影響周圍液相流場和傳質特性。彈狀流型下,相鄰彈狀氣泡間的液彈內循環已被證實能夠增強氣液傳質[21];棒槌狀氣泡的凹部表明,與之相鄰的液膜中存在較高壓力和壓力梯度,這有助于加速壁面流體的流動,進而提高氣泡主體的界面更新與傳質[15,26];平尾狀氣泡兩端較大的壓降有利于相鄰液彈間的交匯,從而增強氣泡主體與頭尾部的傳質速率;在尖尾狀氣泡流型下,液膜厚且不易飽和,因而液膜的傳質貢獻不容忽略[15]。此外,較大的界面剪切和Marangoni 效應也有助于傳質增強[37]。目前,關于彈狀氣泡流型下的氣液傳質特性已有大量報道,而其他形狀氣泡周圍的液相流場和傳質特性尚需進一步研究。

2.2 氣泡形狀分布和預測

為進一步探究氣泡形狀分布及轉變機理,系統分析黏度體系和兩相流速對氣泡形狀的影響[38]。利用兩相Ca數(Ca=uμL/σ,表示剪切力與界面張力的比值,其中u為兩相表觀速度)和QG/QL繪制氣泡形狀分布圖(圖7)。子彈狀氣泡出現在低Ca(0.002~0.203)區域內,且隨QG/QL的增加而變窄。因為在Ca>0.0058,由擠壓和剪切共同控制下的氣泡生成尺寸較小,抗衡此狀態下氣泡界面張力需要更大的剪切力[17]。增加Ca,氣泡形狀由子彈狀轉變為棒槌狀,而棒槌狀向平尾狀和尖尾狀的轉變需要更高Ca。增加液相黏度,流體流動更為穩定,導致棒槌狀向尖尾狀轉變需要更大剪切力,從而造成棒槌狀-尖尾狀的轉變線相比于棒槌狀-平尾狀轉變線向更高的Ca偏移。氣泡形狀區域的轉變線可根據QG/QL和Ca數進行擬合預測,如式(4)所示。

圖7 氣泡形狀分布圖與轉變線

式中,參數a、b和c為擬合值,可根據最小二乘法計算得到。子彈狀與棒槌狀之間的轉變線,見式(5),棒槌狀與平尾狀之間的轉變線見式(6),棒槌狀與尖尾狀之間的轉變線見式(7)。

如圖7所示,轉變線的計算值與氣泡形狀過渡的實驗值之間具有良好的一致性。該轉變線適用于矩形截面微通道(600μm×400μm)、0.002<Ca<2.588 和0.5<QG/QL<5 的實驗條件。由以上分析可知,液相黏度和QG/QL的增加有利于形成尖尾狀氣泡。

2.3 氣泡形狀演變過程

由圖6可知,氣泡平尾狀或尖尾狀是由棒槌狀演變而來,并且演變過程與時間相關[20]。通過追蹤氣泡特性長度(LB)和形狀轉變位置對氣泡形狀的演變機理進行研究。各種氣泡形狀特征長度的丈量方式如圖4所示。

固定甘油質量分數為40%,探究低黏體系中氣液流率比對氣泡由棒槌狀向平尾狀演變的影響,如圖8 所示。其中,坐標x表示氣泡中心與氣液入口通道交叉口間的距離。氣泡在運動初期為棒槌狀,隨著流體逐漸向通道出口流動,液膜內泄漏流減小[39],氣泡尾部受液彈壓力影響逐漸顯著而演變為平尾狀。氣泡長度增加同樣不利于液膜泄漏,從而增強上游液彈擠壓作用,更加易于平尾狀氣泡形成,因而增大QG/QL,氣泡形狀轉變位置不斷向通道入口偏移。QL的增加也會使氣泡形狀轉變位置向通道入口偏移。這是因為,氣泡上游液彈壓力隨著QL的增大而增大,氣泡尾部壓力及早地滿足P1≥P2。此外,在QG/QL>4時,氣泡趨向由慣性力主導的環狀流過渡,氣泡長度隨氣液流率比的增大而略微增加[40],并導致氣泡形狀轉變位置略有不同。

圖8 微通道內氣泡由棒槌狀向平尾狀的演變過程(c=40%)

固定液相甘油質量分數80%,高黏液相體系下氣泡由棒槌狀向尖尾狀的演變規律如圖6(b)和圖9所示。這一演變過程伴隨初期階段LB較大幅度增長,直至氣泡尖尾發展完全,LB幾乎保持不變。QG/QL增加同樣有利于氣泡向尖尾狀過渡,導致轉變點向通道入口偏移。并且,較大QL下較強液相剪切將縮短表面活性劑在氣泡界面的轉移和吸附時間,從而增大尖尾狀氣泡流動范圍。研究表明,氣泡圓弧狀向尖尾狀轉變需大于臨界黏性剪切力[20,29]。經研究,本實驗氣泡由橢圓形狀向尖尾狀轉變的臨界Ca值為0.08。進一步增大Ca,剪切效應增強且氣泡尾部界面張力因表面活性劑聚集而進一步減小,導致氣泡尖尾處微氣泡剪切脫離[33]。氣泡尖尾發生破裂時的臨界LB不隨QG/QL而改變,而隨QL增大而減小,說明尖尾破裂現象與液相剪切密切相關。

圖9 微通道內氣泡由棒槌狀向尖尾狀的演變過程(c=80%)

2.4 氣泡形狀轉變規律

兩相流率和液相黏度對微通道內氣泡形狀轉變位置的影響如圖10 所示,轉變位置(x)與QG/QL滿足冪律關系(QG/QL~axb)。由以上分析可知,較長氣泡易受液相壓力和剪切力作用而發生尾部變形,因而QG/QL越大,各種形變情況下臨界轉變位置x越小,即越靠近通道入口。值得注意的是,甘油質量分數為80%時,氣泡由棒槌狀向尖尾狀的轉變距離并未隨QL的增加而單一遞減[圖10(c)]。當QL從40mL/h 增至60mL/h 時,高黏體系(0.13<Ca<0.68)下由剪切機理控制下的氣泡生成尺寸顯著減小[圖5(e)][41],氣泡形狀較為穩定而不易變形[30],因而形狀轉變位置向通道出口偏移。而當QL進一步增加,液相黏性剪切力足以抗衡界面張力,致使轉變線不斷靠近通道入口。

圖10 不同液相黏度下微通道內氣泡形狀轉變規律

2.5 尖尾狀氣泡破裂的臨界條件

尖尾狀氣泡在微通道內的流動可分為不破裂和破裂兩類,氣泡尖尾破裂存在臨界條件[20,29]。Wu等[20]通過研究微通道中氣泡的尖尾破裂現象,利用毛細管數Ca與量綱為1 的長度L(L=LB/w)建立氣泡尖尾破裂臨界條件的預測模型(Lc=aCab,Lc為臨界無量綱長度)。Wang等[29]利用量綱為1的Ca、Re和奧內佐格數(Oh)預測液滴尖尾破裂臨界條件。相關模型對本實驗中氣泡尖尾破裂臨界條件的預測結果如圖11所示。Wu等[20]模型可很好模擬破裂臨界條件趨勢,但卻與實驗值存在偏差[圖11(a)],這可能由于通道截面幾何結構不同所致。尖尾狀氣泡尾部破裂得益于上游液彈壓力積累與壁面摩擦,而矩形截面通道內角區相較于方形截面更厚,液膜易于在角區流通,不利于液彈內壓力積累,導致氣泡尖尾破裂需要更大Ca數(即較大液相流速和黏度)。此外,模型中Lc為氣泡流動過程的一項因變量,方程使用具有不便。Wang 等[29]模型不能很好描述本實驗氣泡尖尾破裂的臨界條件,兩者存在很大偏差[圖11(b)],可能是由于液滴與氣泡的物理性質,以及所用表面活性劑的濃度和種類不同。兩者都會影響相界面的張力梯度,導致尖尾破裂的臨界剪切力存在差異[20]。

圖11 尖尾狀氣泡破裂分布圖及轉變線

利用操作條件QG/QL和Ca繪制氣泡尖尾形狀破裂與不破裂分布圖(圖12),發現破裂臨界條件遵循式(8)關系。

圖12 尖尾狀氣泡破裂分布圖及轉變線

當QG/QL<0.83Ca-0.61,氣泡尾部在微通道內趨于穩定,不發生破裂,而當QG/QL>0.83Ca-0.61,氣泡尾部發生尖尾破裂。此關系式適用于實驗條件為0.5<QG/QL<5、0.08<Ca<2.59。

3 結論

本文利用高速攝像對矩形微通道(600μm×400μm)內不同黏度體系下氣泡形狀及其演變現象進行研究,揭示了0.002<Ca<2.588 和0.5<QG/QL<5條件下矩形通道構型、液相壓力和黏性對運動氣泡界面及其分布、演變和破裂的影響機制,為微通道反應器內氣-液兩相流調控提供理論指導。主要結論如下:

(1)通道內觀察到子彈狀、棒槌狀、平尾狀和尖尾狀四種氣泡形狀,平尾狀和尖尾狀氣泡界面分別由液相壓力和黏性剪切力控制,利用Ca數和氣液流率比可很好地預測氣泡形狀。

(2)高黏體系(0.13<Ca<0.68)下由剪切機理控制的氣泡生成尺寸影響運動過程的界面穩定性。在相對較低剪切效應下小氣泡的高界面張力穩定界面形狀,而當剪切力升高至一定值,并足以能對抗小氣泡界面張力,氣泡運動過程界面發生變形。通過探究微通道內氣泡特征長度的演變過程,發現平尾狀和尖尾狀氣泡均由棒槌狀氣泡在運動過程演變而來,并且氣泡形狀轉變距離與Ca數和氣/液流率比呈負相關。

(3)提出了矩形截面微通道內尖尾狀氣泡破裂的臨界條件模型。當QG/QL<0.83Ca-0.61,氣泡尾部在微通道內趨于穩定,不發生破裂。相反,當QG/QL>0.83Ca-0.61,氣泡尾部發生尖端破裂。

符號說明

Ca——兩相毛細管數(Ca=μLu/σ)

c——質量濃度,%

h——深度,m

L——量綱為1氣泡長度(L=LB/w)

LB——氣泡長度,m

Oh——奧內佐格數[Oh=μ(ρwσ)-0.5]

P——壓力,Pa

Q——體積流率,mL/s

R——半徑,m

Re——雷諾數(Re=ρQ/wμ)

r——半徑,m

u——兩相表觀流速,m/s

w——寬度,m

μ——黏度,Pa·s

ρ——密度,kg/m3

σ——界面張力,N/m

下角標

c——臨界值

f——液膜

G——氣相

head——氣泡頭部

L——液相

rear——氣泡尾部

猜你喜歡
平尾界面
基于流場顯示技術的冰污染平尾失速試飛技術
雙層平尾對旋翼/平尾干擾的抑制機理研究
直升機技術(2021年4期)2022-01-12 13:18:12
國企黨委前置研究的“四個界面”
當代陜西(2020年13期)2020-08-24 08:22:02
民用飛機平尾載荷的不確定性及全局靈敏度分析
平尾損傷計算中特征載荷的算法研究
直升機技術(2018年3期)2018-10-09 10:44:12
全動式水平尾翼
大飛機(2018年1期)2018-05-14 15:59:08
基于FANUC PICTURE的虛擬軸坐標顯示界面開發方法研究
空間界面
金秋(2017年4期)2017-06-07 08:22:16
電子顯微打開材料界面世界之門
人機交互界面發展趨勢研究
主站蜘蛛池模板: 精品无码人妻一区二区| 麻豆国产在线观看一区二区| 国产精品亚洲а∨天堂免下载| 国产精品久久久久久影院| 久久精品一卡日本电影| 欧美色图第一页| 日韩精品一区二区三区视频免费看| 日本一区二区不卡视频| 九九九精品成人免费视频7| 91精品福利自产拍在线观看| 九九线精品视频在线观看| 又爽又黄又无遮挡网站| 3344在线观看无码| 爽爽影院十八禁在线观看| 少妇被粗大的猛烈进出免费视频| 好紧太爽了视频免费无码| 亚洲狼网站狼狼鲁亚洲下载| 婷婷色一二三区波多野衣| 91欧美亚洲国产五月天| 国产成人在线无码免费视频| 欧洲亚洲欧美国产日本高清| 国产欧美专区在线观看| 免费无码AV片在线观看国产| 青青草国产一区二区三区| 全午夜免费一级毛片| 国产综合亚洲欧洲区精品无码| 日本一区中文字幕最新在线| 曰韩免费无码AV一区二区| 伊人久久影视| 国产免费精彩视频| 在线中文字幕网| 亚洲女同欧美在线| 92午夜福利影院一区二区三区| 亚洲国产精品一区二区第一页免 | 中文字幕亚洲第一| 中文字幕色在线| 伊人成人在线视频| 国产日本视频91| 国产99视频精品免费观看9e| 精品国产免费观看一区| 国内老司机精品视频在线播出| 欧美视频在线观看第一页| 毛片大全免费观看| 久草性视频| 久久无码高潮喷水| 青青草原偷拍视频| 亚洲国产成人自拍| www亚洲天堂| 亚洲精品欧美日本中文字幕| 亚洲av色吊丝无码| 亚洲欧美一级一级a| 国产真实二区一区在线亚洲| 亚洲成人高清无码| 最新国产午夜精品视频成人| 国产午夜不卡| 99国产精品国产| 久久99精品久久久久纯品| 欧美中文字幕无线码视频| 在线观看精品国产入口| 无码精品国产VA在线观看DVD| 福利一区在线| 国产精品福利社| 青青操国产| 国产女人在线视频| 中美日韩在线网免费毛片视频 | 亚洲欧美日韩动漫| 麻豆精品在线播放| 久久久久九九精品影院| 国产高清又黄又嫩的免费视频网站| 日韩亚洲综合在线| 久久国产热| 亚洲人成影院午夜网站| 好吊色妇女免费视频免费| 日本AⅤ精品一区二区三区日| 久久99久久无码毛片一区二区| 美女黄网十八禁免费看| 久久综合色视频| 伊人久久大线影院首页| 国产日本视频91| 久久久亚洲色| 日韩区欧美区| 日韩无码一二三区|