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

氣固兩相湍流場納米顆粒演變特性綜述

2022-01-10 07:55:00石瑞芳林建忠
航空學報 2021年12期

石瑞芳,林建忠

浙江大學 航空航天學院 流體工程研究所,杭州 310027

存在于流動氣體介質中的納米顆粒以氣溶膠的形式出現,多數情況下流動呈湍流狀態,如空氣凈化系統、噴射器、呼吸道流[1]、汽車尾氣噴流、地效飛機懸停流場[2]和燃燒排氣系統等。納米顆粒兩相湍流在自然界和實際應用中很普遍,了解湍流場內納米顆粒的尺度分布和演化機理,有助于人們在實際應用中對其控制從而達到預期的目標。

要有效地預測和控制納米顆粒在湍流場中的分布,就要了解顆粒在湍流場中的生成、對流、擴散、凝并、破碎等機理,為此人們進行了大量的研究[3-8]。對單顆粒而言,由于Stokes數(St)遠小于1,所以當顆粒的體積分數較小時,可以忽略顆粒對流場的影響而采取單向耦合的方法求解方程。但如果顆粒呈致密的類分形團塊,因顆粒具有獨特的傳熱傳質特性,故需考慮顆粒對流場的影響[9]。納米顆粒兩相湍流的湍動特性體現在流場速度脈動、顆粒數密度脈動、化學組分濃度脈動、溫度和飽和度脈動等[10]。與層流相比,湍流場對顆粒擴散的影響更大并顯著增強顆粒在流場中的混合。受湍流場影響,顆粒的分布高度不均勻[11-13],湍流會在流動界面卷起形成兼具大尺度和小尺度運動的擬序結構,該結構的存在對增強顆粒凝并和彌散起著重要的作用[5]。由于納米顆粒兩相流的顆粒數一般都非常多,所以對顆粒群的研究通常采用統計的方法。常見的研究納米顆粒動力學的方法有分子動力學、分區法、矩方法、蒙特卡洛方法等[14]。

以下從納米顆粒氣固兩相湍流場中最常見的顆粒生成、凝并、破碎和沉降4個方面敘述相關的研究狀況和進展。

1 顆粒生成

不同尺寸和結構的顆粒,呈現的物理性質往往有所差別。實際工業應用中,通過不同方式制備特定尺度范圍的功能性納米顆粒。圖1給出了納米顆粒的生成方式,由于本文涉及的是氣固兩相流,所以以下主要介紹與此相關的流場中顆粒的成核和表面反應生長(冷凝)過程。

圖1 納米顆粒的產生方式

1.1 納米顆粒在不同環境中的生成機理

納米顆??梢栽诖髿猸h境[15-16]、湍流噴霧火焰[17]和湍流反應流[18]等不同環境下生成。顆粒在氣體環境中的生成是氣相化學反應產生的可冷凝蒸汽物質因表面冷卻、絕熱膨脹或混合、湍流混合[19]或化學過程產生的過飽和所導致;而在液相中的生成與顆粒的粒徑、化學成分、表面和電荷性質等有關[20]。不同來源與生成機理的顆粒呈現出不同的顆粒尺度、形態和結構特性及散射性、吸附性和生長性等物理特性。

液相體系中得到納米顆粒如水/有機溶劑在聚合物做保護劑的前提下被還原劑還原生成金屬納米微粒,其形狀主要取決于保護劑的種類、劑量、聚合度及配料的速率和順序等,產物有納米線、棒、無規則、球、橢球或三角形等多種形式。太陽能系統加入不同形狀的銀納米顆粒可不同幅度提高溶液的吸光度和光熱轉換性[21]。航天復合材料,加入納米顆粒可降低航天器材料的質量損失和剝蝕率[22],激光加工制備的納米尺度顆粒和微孔陣列形成微觀機械咬合等可提高接頭的強度和黏附性[23]。

相對而言,氣相體系下制備的納米顆粒不發生化學反應,僅通過熱源使可凝性物質在高溫下蒸發,并在惰性氣體氛圍下冷凝從而形成納米顆粒。顆粒大小取決于惰性氣體壓力和溫度,可制備精細表面清潔的顆粒,但效率低難控制產物形狀。產物形狀更多是顆?;蜴湢钸B接顆粒乃至納米線。一維的納米線,應用于傳感器原件,表現出極高的柔韌性、響應速度和靈敏度[24]。另一種氣體原料在氣相介質內發生化學反應得到基本粒子則會經成核、生長階段,更多以團聚體或點、角相接的團簇或小顆粒在大顆粒上附著的低維附聚體形式存在。源自氣相燃燒合成的顆粒,其顆粒除離散單體形式外,更多的是顆粒聚集形態,如原級粒子以面相接,難再分散的凝聚體。懸浮納米顆粒在界面堆積可改變沸點提高液體燃料著火率[25]。

顆粒在生成階段影響其尺度和形狀的主要控制參數是成核和生長速率。成核需要克服能壘形成的臨界核才能實現,成核過程中所需的成核能量決定了成核速率,隨著能壘的增加,成核速率呈指數下降[26]。顆粒成核和表面生長理論如表1[27-32]所示,顆粒成核方程為

表1 顆粒成核和生長理論

(1)

式中:n(v,t)表示以顆粒體積v和時間t為變量的顆粒尺度分布函數;v*為生成的穩定單體體積;δ為狄拉克函數;G(v)和J表示冷凝導致的顆粒生長核函數與成核速率。

1.2 顆粒生成的數值模擬

發動機尾氣排放是典型的顆粒生成過程,人們對與尾氣排放相關的各類射流以及與擴散火焰反應器相關的顆粒生成機理進行了一系列研究。Yin等用大渦模擬方法結合顆粒成核和演化模型先后研究了運動排氣雙噴流[33]和撞擊雙射流[34]中納米顆粒的形成。受行駛物尾部流場結構的影響,顆粒的生成取決于環境風速與排氣速度之比,較高的環境風速會降低納米顆粒的生成速率并增大顆粒的尺度。在撞擊射流中,顆粒的成核主要發生在射流界面及撞擊平面附近的區域,其中撞擊平面附近顆粒成核數最多;噴嘴到撞擊平面的距離越大越不利于顆粒成核。在平行雙射流中,顆粒成核區與撞擊射流的情形類似,硫含量、相對濕度、射流雷諾數能促進顆粒成核[35]。Yu等[36]數值模擬了擴散火焰反應器四異丙醇鈦分解均相成核和顆粒的形成過程。在此基礎上,Yu等[37]進一步給出了前驅物量值對擴散火焰反應器合成非球形二氧化鈦納米顆粒的影響。此后,Yu等[38]將大渦模擬方法與新提出的泰勒級數展開矩方法結合,對發動機尾氣中的顆粒生成進行了數值模擬,發現大渦主導著顆粒的演化,硫化物-水二元均質成核主要出現在排氣與周圍冷氣體的交界面;燃料含硫量和相對濕度的增加或環境溫度的降低均導致顆粒生成率以及平均粒徑的增加;隨著湍流度的增大,顆粒的分布變寬、平均直徑變小(見圖2)[38]。

圖2 湍流強度I對顆粒數密度和平均粒徑的影響[38]

Yu和Lin[39]研究了有和無背景顆粒情況下水-硫化物的混合物在大氣環境中二元均相成核的過程,發現顆粒成核及其后續的運動強烈依賴于硫化物的濃度。Lin和Liu[40]對混合層中顆粒成核的研究表明,溫度較低的區域顆粒更容易成核。Chan等[41]對飛行器近尾跡區顆粒的凝并、成核和湍流彌散過程的研究發現,顆粒凝并的時間尺度較大,而成核的時間尺度較小。Chan等[42]將雙峰泰勒級數展開矩方法與大渦模擬方法結合,證實了湍流的擬序結構增強了顆粒的擴散。Garrick[18]研究了二氧化鈦顆粒生成過程中冷凝和凝并所起的作用,如圖3[18]所示,在近射流區,成核和冷凝起主導作用,一旦射流核心區坍塌,凝并起主導作用。圖中1 nm顆粒(綠色)遍布整個流場,而凝并后的3 nm顆粒僅存在于下游區域。

圖3 顆粒在流場中的瞬時分布[18]

2 顆粒凝并

顆粒凝并是顆粒相互碰撞并黏附形成大顆粒的過程,如表2[43]所示,導致顆粒碰撞的原因包括布朗運動、湍流場剪切作用、速度梯度、差異沉降等。顆粒的凝并取決于顆粒的尺度和流場的特性[44-45]。

表2 顆粒凝并機制[43]

2.1 布朗凝并

布朗凝并通常發生在小于1 μm的顆粒,不同角度、不同直徑和不同碰撞類型的顆粒發生碰撞時,其碰撞凝并有不同的表達形式[46-50]。Smoluchowski最早提出凝并模型,將顆粒尺度譜的演化轉化為顆粒碰撞頻率函數的演化,后來Müller[51]發展了連續顆粒尺度分布的數密度方程:

(2)

式中:β(v,v1)表示體積分別為v1和v的兩顆粒凝并的顆粒凝并核函數,取決于顆粒間的碰撞和凝并率。顆粒布朗凝并取決于顆粒布朗運動引起的碰撞,當顆粒粒徑遠大于氣體分子平均自由程時,顆粒的碰撞受擴散制約;而當粒徑遠小于分子平均自由程時,顆粒的碰撞可由分子運動理論確定。在近連續區和自由分子區,顆粒的布朗凝并核函數可分別表示為[52]

(3)

(4)

式中:λ為氣體分子平均自由程;K和T分別為Boltzmann常數與溫度;μ和ρ為流體動力黏度和密度;vp和rp為初始顆粒的體積和半徑;Df為顆粒分形維數,表征顆粒的形態。

2.2 湍流凝并

湍流凝并通常發生在粒徑為1 μm以上的顆粒或者是高湍流情況下的亞微米顆粒。Camp和Stein[53]用湍動能耗散與流體黏度比值的平方根代替湍流中的速度梯度,將用于層流剪切凝并的表達式推廣到湍流剪切凝并。粒徑為1 μm以下的顆粒通常位于黏性區,Saffman和Turner[54]導出了黏性區內粒徑遠小于Kolmogorov尺度的顆粒凝聚核,其中由剪切引起的球形顆粒凝并核為

(5)

式中:ν和ε分別為流體黏度和湍動能耗散率;r1和r2為顆粒半徑。Flesch等[55]將式(5)拓展到分形結構體:

(6)

總的顆粒凝并核可以由布朗凝并βB和湍流凝并βT構成:

(7)

在實際應用中,不同湍流場對顆粒凝并有不同的作用。馮鵬等[56]比較了不同類型的繞流場中顆粒的凝并,發現Y型繞流場由于有更強的渦度場,所以凝并最為明顯。Yu等[57]研究了湍旋流中的顆粒凝并,分析了渦流發生器類型、結構參數和運行參數對亞微米顆粒凝并的影響,給出了當螺旋葉片內外徑比為0.3時顆粒濃度和速度對100 nm顆粒凝并的影響(如圖4[57]所示,N為顆粒濃度,U為流場速度,dp為顆粒直徑,r為顆粒體積分數),說明顆粒濃度呈單峰形式分布且峰值隨著流速的增加而迅速增加。進一步的研究表明,如果流場的速度增加,則小顆粒的數量增加,而大顆粒的數量減少,其原因一是較高的流速會影響顆粒在濃度梯度方向上的運動,顆粒擴散受到一定程度的抑制,使得顆粒凝并削弱;原因二是較高流速使得顆粒停留在流場中的時間更短,顆粒碰撞的概率降低。若保持流速不變,分布曲線向右移動,表明顆粒凝并增強,因為較高的顆粒濃度導致顆粒之間更多的相互作用和碰撞。

圖4 濃度和速度對凝并的影響[57]

2.3 數值求解方法

Yu等[58-61]提出了對顆粒數密度方程進行求解的泰勒級數展開矩方法并用于多種流場的計算,在保證精度的同時提高了計算效率。在此基礎上,Lin和Chen[62]基于分區法的漸近解,改進了泰勒級數展開矩方法,使得當顆粒尺度分布接近自相似或凝并時間足夠長時,能更精確地預測零階和二階矩的演化。Chen等[63]在泰勒級數展開矩方法的基礎上,提出了直接展開矩方法并用于求解顆粒數密度方程,給出了較大顆粒尺度范圍內顆粒布朗凝并的解[64]。Yu和Lin[65]提出了一種混合矩法來求解顆粒數密度方程,該方法采用三階泰勒級數展開來逼近具有附加性形式的凝聚核,用拉格朗日插值來逼近隱式矩,從而避免了針對特定核的推導,其有效性在應用中得到了證實。

2.4 湍流場的影響

Garrick[66]在平面混合層中研究了湍流對顆粒凝并的影響,發現在有顆粒的區域,流場大尺度與亞網格尺度的相互作用能促進顆粒的凝并。Das和Garrick[67]對平面射流中二氧化鈦顆粒凝并的研究表明,小尺度的波動對顆粒的凝并既有促進作用,也有抑制作用,但后者更為明顯。納米顆粒在充分發展的湍流管道內會向管道中心移動,單分散顆粒會因布朗運動的存在而凝并成不同大小的顆粒簇,在管道中心則產生最大的顆粒團[68]。Lin和Chen[69]發現,在充分發展的湍流邊界層中,顆粒的凝并與擴散使得顆粒的粒徑從壁面向外側增加,顆粒在向外側遷移中因較高的局部濃度而促進了凝并,導致較大的顆粒在外側積聚。Lin和Liu[40]的研究表明,混合層中的納米顆粒凝并受初始顆粒分布及湍流擴散控制,圖5[40]給出了顆粒在不同流向位置的無量綱平均直徑dave分布,可見在向下游發展過程中,顆粒的直徑分布逐漸變寬。

圖5 顆粒在不同流向位置的平均直徑分布[40]

湍流場對顆粒凝并的影響除了湍流強度的因素外,還體現在由湍流脈動所引發的顆粒數密度脈動。Lin等[70]以湍動能與平均動能的比值來表征顆粒數密度脈動對顆粒凝并的影響,并通過圓管湍流場的計算驗證了其有效性,數值模擬結果表明,在小Schmidt 數和高Reynolds數情況下,顆粒呈現出更大的平均直徑和更強的尺度分散性。然而,當流場存在回流區或湍動能遠大于平均動能時,以上方法不適用,對此Yang等[71]在以上方法基礎上提出了另一種形式,即以湍動能與總動能的比值來表征顆粒數密度脈動對顆粒凝并的影響。Shi等[72]采用這一方法對通風室內非球形顆粒的凝并進行了數值模擬,得到了與實驗[73]比較吻合的結果。在有顆粒凝并的流場中通常伴隨著顆粒破碎,Gan等[74]數值模擬了平面湍射流中存在顆粒凝并和破碎情況下的顆粒數密度方程,發現存在凝并與破碎的動態平衡。Barthelmes等[75]從理論上研究了非球形顆粒的動力學行為,建立了剪切誘導的顆粒凝并和破碎的群體平衡模型。

3 顆粒破碎

顆粒凝并后形成尺度較大的團聚體容易在流場剪切和其他因素作用下發生解體破碎,與顆粒破碎相關的方程為

(8)

式中:a(v)為破碎核函數,表示體積為v的顆粒的破碎頻率,有指數和冪律等形式;b(v|v1)為碎片分布函數,表3給出了常見碎片的分布形式。Marchisio等[76]在考慮不同破碎核與子分布形式下,用積分矩方法模擬了顆粒的凝并和破碎,并給出了存在凝并和破碎時的漸近解[77]。

表3 顆粒破碎后的碎片分布

在攪拌流場中,顆粒容易受攪拌的外力場作用而破碎[78],攪拌器的速度和流動速度是影響顆粒破碎的2個主要因素。Pandya和Spielman[79]根據顆粒破碎和侵蝕機制,提出了一個控制湍流中顆粒團聚體粒徑分布的平衡方程式。此后,Hill和Ng[80]用新的離散化方法,推導出了適用于分散體系顆粒群的破碎方程,解決了區間交互問題,明顯提高了對顆粒粒徑分布的預測精度。

剪切破碎是導致顆粒破碎的主要因素,由剪切導致的有效破碎系數取決于剪切率和顆粒的體積分數,體積分數越大,顆粒的碰撞次數增加,剪切率越高意味著碰撞能量越強,顆粒越容易破碎[81]。當剪切率一定時,顆粒的體積分數增加,穩態的顆粒尺度減小。

當流場進出口的差壓較大時也容易導致顆粒破碎。Yuan等[82]研究了顆粒從高壓容器的噴嘴中排出的情形,通過控制壓力使得顆粒破碎后的尺度滿足所需的要求。Badyga等[83]對轉-定子混合器和高壓噴嘴粉碎機流場的數值模擬結果表明,高壓噴嘴粉碎機對顆粒的破碎產生更明顯的效果。圖6[83]給出了壓差(Δp)對顆粒平均粒徑(L30)的影響,可見進出口的壓差越大,顆粒的破碎程度越高。

圖6 壓差對平均粒徑的影響[83]

4 顆粒沉降

在氣體凈化設備、氣溶膠測量儀器以及人體呼吸這樣的納米顆粒氣力輸送過程中,預測顆粒的沉降速率至關重要。顆粒的沉降取決于顆粒的尺度、形狀和流體性質等因素。導致顆粒沉降的因素有重力、擴散、慣性撞擊、電場和熱遷移等[19]。在湍流場中,顆粒的沉降與黏性底層和過渡層中的參數密切相關,在充分發展的湍流邊界層中,粒徑較小顆粒的沉降主要因布朗擴散和湍流擴散所致;粒徑較大的顆粒主要受其慣性和重力支配,當地流場的影響較小。

顆粒沉降到壁面上會導致在壁面上的沉積,沉積率或輸運顆粒損失率η可表示為

(9)

式中:Ci和Co為進出口顆粒的質量濃度或數密度;fp為顆粒流過管道的通過率,對充分發展的圓管湍流場有:

(10)

其中:顆粒的沉降速度vd定義為單位時間、單位面積上的顆粒通量Jd和截面平均濃度Cave的比值;l和D分別為管道長度和直徑;U為流場速度。用壁速度uw無量綱化,可得無量綱沉降速度:

(11)

Lai和Nazaroff[84]考慮了布朗擴散和湍流擴散及重力因素對顆粒沉降的影響,提出了新的半經驗沉積模型(也稱“三層模型”)用于預測沉降到光滑表面的顆粒:

(12)

式中:εp為顆粒的渦流擴散率;DB為顆粒在邊界層內的布朗擴散系數;vs為重力引起的沉降速度;i表示曲面方向(表面向上取1,向下取-1,垂直為0);C為顆粒的質量濃度或數密度。

Chen和Lai[85]基于簡化的三層模型,在考慮了布朗擴散和湍流擴散、重力、庫侖力以及虛假慣性力的基礎上,建立了一個修正的費克定律方程,并提出一種新的顆粒漂移通量模型用于數值模擬室內顆粒的分布與沉降[86]。付崢嶸[87]提出的解析法模型能很好地預測通風空調管道充分發展湍流段中顆粒的沉降規律。擴散慣性模型[88]能預測湍流中顆粒的擴散和沉降規律。兩層分區的湍流模型能很好地預測圓管內納米顆粒的沉降速度[89]。Guichard等[90]提出了用粒徑分布矩表示的顆粒輸運和沉降模型。

當顆粒Stokes數很小且無外力作用時,顆粒因布朗擴散和湍流擴散遷移到壁面[6, 91-92]。當存在溫度梯度時,熱泳力對顆粒沉降也起到重要作用。研究結果表明,熱泳力的存在能增強對湍流擴散敏感的顆粒的沉降[93]。

Lin等[94]研究了不同Reynolds數下圓球形納米顆粒通過彎管時的穿透率,給出了8~550 nm顆粒的穿透率與綜合參數的關系式,該參數包括了Schmidt數、Dean數和管道長度的作用。他們還進一步研究了圓柱形納米顆粒通過彎管時的穿透率[95],其中穿透率(PE)與Stokes數(St)的關系如圖7[95]所示,可見隨著Stokes數的增加,穿透率先增大然后減小;同時,穿透率隨著Dean數(De)、Reynolds數(Re)、顆粒長徑比的減小而增大?;跀抵的M的數據,他們建立了穿透率與Dean數、Reynolds數、顆粒長徑比、Stokes數之間的關系。

圖7 穿透率與St和Re的關系(長徑比為8,De=1 862)[95]

Ounis等[96]對槽道兩相湍流場的直接數值模擬表明,顆粒布朗運動效應和擬序渦流形成的往壁面的流動導致了亞微米顆粒在壁面的沉積,顆粒的沉積數量和沉積速度隨著顆粒直徑的減小而增加。Wu和Young[97]研究了小顆粒在壁面上的沉積,圖8[97]給出了顆粒密度的分布,可見直徑為0.01 μm的最小顆粒占據了整個流場,因為它們具有很好的跟隨性,這種尺度顆粒的沉積主要是由擴散引起的,顆粒沉積在葉片的吸力面和壓力面上。然而,即使這么小的顆粒,慣性也起著重要的作用,因為壓力面的顆粒密度高于吸力面。

圖8 顆粒密度等值線圖[97]

呼吸道內顆粒的沉積也是人們關注的熱點。Mina等[98]研究了納米顆粒在人體呼吸系統上部的沉積,其結果如圖9[98]所示,研究結果表明,顆 粒越小其擴散性越強、顆粒濃度衰減就越快即顆粒的沉積量就越大。

圖9 不同尺度顆粒的質量分數分布[98]

5 結論與展望

本文回顧了含納米顆粒氣固兩相湍流場中顆粒的生成、凝并、破碎、沉降的部分研究狀況,主要研究結果歸納如下。

1)顆粒在氣體環境中的生成是氣相化學反應產生的可冷凝蒸汽物質因表面冷卻、絕熱膨脹或混合、湍流混合或化學過程產生的過飽和所導致;燃料含硫量和相對濕度的增加或環境溫度的降低均導致顆粒生成率以及平均粒徑增加;隨著湍流度的增大,顆粒的分布變寬、平均直徑變小;溫度較低的區域更容易導致顆粒成核。

2)導致顆粒碰撞凝并的原因包括布朗運動、湍流剪切作用、速度梯度、差異沉降;顆粒的凝并取決于顆粒的尺度和流場的特性;布朗凝并通常發生在小于1 μm的顆粒,湍流凝并通常發生在粒徑為1 μm以上的顆粒或者是高湍流情況下的亞微米顆粒;小尺度的流場脈動對顆粒的凝并既有促進也有抑制作用,但后者更為明顯;顆粒凝并受初始顆粒分布及湍流擴散控制;湍流場對顆粒凝并的影響除了湍流強度的因素外,還體現在由湍流脈動所引發的顆粒數密度的脈動。

3)顆粒凝并后形成尺度較大的團聚體容易在流場剪切和其他因素作用下發生解體破碎;剪切破碎是導致顆粒破碎的主要因素;有效破碎系數取決于剪切率和顆粒的體積分數,體積分數越大、剪切率越高,顆粒越容易破碎;當剪切率一定時,顆粒的體積分數增加,穩態的顆粒尺度減?。贿M出口壓差越大,顆粒破碎程度越高。

4)顆粒的沉降取決于顆粒尺度、形狀和流體性質等因素;導致顆粒沉降的因素有重力、擴散、慣性撞擊、電場和熱遷移等;湍流場中顆粒的沉降與黏性底層和過渡層中的參數密切相關;粒徑較小顆粒的沉降主要因布朗擴散和湍流擴散所致,粒徑較大的顆粒主要受其慣性和重力支配;當存在溫度梯度時,熱泳力對顆粒沉降起到重要作用。

如前所述,雖然含納米顆粒氣固兩相湍流場的研究已取得了許多有價值的成果,但仍有一些問題有待于進一步研究。

1)湍流脈動對于納米的擴散有明顯的影響,需要建立由湍流脈動導致的納米顆粒擴散項的表達式。

2)湍流脈動不僅對顆粒的速度有影響,也會導致顆粒數密度的變化,所以要建立由湍流脈動導致的納米顆粒數密度增減項的表達式。

3)納米顆粒兩相流場的動力學特性覆蓋較大的空間和時間尺度,需要建立具有不同特征尺度作用的納米顆粒兩相湍流場較完整的數理模型并付諸于實際應用。

主站蜘蛛池模板: aaa国产一级毛片| 国产精品伦视频观看免费| 国产福利不卡视频| 538精品在线观看| 午夜视频www| 精品无码一区二区三区在线视频| 天天综合网色中文字幕| 国产精品亚洲αv天堂无码| 宅男噜噜噜66国产在线观看| 特级aaaaaaaaa毛片免费视频| 一本无码在线观看| 亚洲无码免费黄色网址| 九色91在线视频| 亚洲综合第一区| 亚洲精品男人天堂| 国产精品9| 国产成人综合欧美精品久久| 日韩无码黄色网站| 在线看片免费人成视久网下载| 国产精品亚欧美一区二区| 夜夜操国产| 国产精品久久自在自2021| 日本免费a视频| 女人爽到高潮免费视频大全| 色亚洲成人| 亚瑟天堂久久一区二区影院| 久草性视频| 久久窝窝国产精品午夜看片| 亚洲av无码人妻| 91精品人妻一区二区| 亚洲日本中文字幕乱码中文| 麻豆精品视频在线原创| 91精品国产福利| 国产精品极品美女自在线网站| 欧美日韩一区二区在线免费观看| 国产欧美日韩专区发布| 三级视频中文字幕| 久久免费视频播放| 久久国产精品影院| 激情网址在线观看| 国产成人精品亚洲77美色| 国产欧美性爱网| 91丝袜美腿高跟国产极品老师| 国产区免费精品视频| 国产精品久久国产精麻豆99网站| 波多野结衣国产精品| 亚洲日韩在线满18点击进入| 毛片免费网址| 亚洲美女一区二区三区| 无码人妻热线精品视频| 综合天天色| 婷婷色在线视频| 国产美女在线观看| 四虎精品黑人视频| 久久青青草原亚洲av无码| 久久不卡国产精品无码| 国产精品刺激对白在线| 性视频久久| 久久99精品久久久久纯品| 91po国产在线精品免费观看| 国产靠逼视频| 国产精品亚洲片在线va| 国产办公室秘书无码精品| 成人毛片免费在线观看| 欧美精品在线看| 久久久亚洲国产美女国产盗摄| 华人在线亚洲欧美精品| 成人免费午间影院在线观看| 超碰精品无码一区二区| 午夜电影在线观看国产1区| 亚洲色中色| 91色老久久精品偷偷蜜臀| 看你懂的巨臀中文字幕一区二区| 免费jjzz在在线播放国产| 亚洲色图在线观看| 色亚洲成人| 米奇精品一区二区三区| 91福利免费| 免费在线一区| 黄色污网站在线观看| 日韩精品资源| 国产另类乱子伦精品免费女|