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

轉盤表面黏性薄液膜穩態流動特性數值模擬

2017-06-05 01:21:40王東祥凌祥彭浩楊新俊崔政偉
化工學報 2017年6期

王東祥,凌祥,彭浩,楊新俊,崔政偉

?

轉盤表面黏性薄液膜穩態流動特性數值模擬

王東祥1,2,凌祥2,彭浩2,楊新俊1,2,崔政偉1

(1江南大學機械工程學院,江蘇省食品先進制造裝備技術重點實驗室,江蘇無錫214122;2南京工業大學江蘇省過程強化與新能源裝備技術重點實驗室,江蘇南京211800)

轉盤表面薄膜流廣泛存在于轉盤反應器、造粒和分子蒸餾等化工領域,其流動特性對造粒、反應以及熱質傳遞具有重要影響。建立了穩態薄膜流動特性理論模型,并對比實驗數據與數值模擬結果,研究了等效Froude數、量綱1特征厚度與長度比值和澆注尺寸r對薄膜厚度分布的影響,導出了水躍和同步區半徑模型并通過實驗驗證。結果表明:等效Froude數不會對薄膜厚度分布產生明顯影響;量綱1和澆注尺寸r是水躍發生與否的決定條件,增大或縮小r有助于水躍出現,水躍半徑穩定于=0.85;平均徑向速度呈典型三分區特征,縮小澆注尺寸r將引起澆注區范圍加大而加速區范圍顯著縮小,導致薄膜無明顯加速現象而直接進入同步區,但同步區半徑維持在=1.53。研究結果為轉盤反應器以及離心粒化器等的設計與優化提供了可借鑒的理論與應用基礎。

流體力學;模型;實驗驗證;流域;轉盤;穩態流

引 言

液體在轉盤表面因離心作用以薄層形式向外流動形成薄膜流,具有高熱質傳遞系數等獨特優點,廣泛存在于分子蒸餾[1]、轉盤反應器[2-8]、離心造粒[9-11]以及高溫熔融渣料余熱回收[12-14]等領域。薄膜流體自由液面形狀或者液膜厚度對造粒和熱質傳遞具有重要影響[15-18]。如可調控液膜厚度以避免傳熱表面發生干斑或干區而導致熱敏性物料變味變質、非熱敏性物料結焦,以至覆蓋傳熱界面引起傳熱能力或者反應速率的顯著下降[19-20];在轉盤粒化法回收熔融渣料余熱工藝中,液膜厚度決定了流體在轉盤邊緣的破碎形態[21-23],進而作用于熔融渣顆粒的凝固換熱和物相轉變[24],是影響余熱回收效率和渣粒資源化再利用的重要因素。

Emslie等[25]最早探索轉盤表面的薄膜流動行為,研究發現若忽略薄膜與轉盤間的速度滑移(即忽略科式力),薄膜流理論上由黏性力與離心力共同作用,較大Ekman數(>1)條件下的液膜厚度與轉盤徑向距離間滿足∝-2/3,提出了離心模型。在一些特殊條件下,某一轉盤區域的液膜厚度會顯著上升,形成水躍[26-29],但實現水躍半徑的理論預測還鮮有涉及。Wang等[30]在離心模型基礎上考慮了速度滑移的影響,研究發現離心模型只適用于水躍之后的薄膜流,科式力會顯著提高水躍區域的薄膜厚度。Woods[31]分析了強慣性力(<1)薄膜流,提出了Pigford模型。Burns等[32]采用電阻法對薄膜流的平均徑向速度進行測量,發現徑向速度分布可分為澆注區、加速區以及同步區,并提出黏性因子,對Pigford模型進行了修正。Sisoev等[33]研究發現空氣誘導剪切力引起薄膜表面可能伴隨波動,但平均液膜厚度不會顯著改變,轉盤軸心和邊緣區域的液膜近似呈穩態流動[34]。

現有研究工作主要針對轉盤軸心或邊緣區域的薄膜流動行為。本文針對整個轉盤表面薄膜流體的穩態流動特性,采用量綱1分析與積分平均法簡化Navier-Stokes(N-S)方程,建立穩態薄膜流動特性理論模型,研究薄膜自由液面形狀和液膜厚度分布,同時探討水躍和同步區演化的一般規律,為轉盤反應器以及離心粒化器等的設計與優化提供可借鑒的理論與應用基礎。

1 穩態薄膜流動模型

1.1 控制方程量綱1化

工質以恒定流量沿轉盤軸心澆注,液膜呈旋轉周期性流動。考慮一無限大轉盤,水平放置,以恒定角速度沿軸旋轉。工質澆注口半徑為in,流量為,牛頓流體,不可壓縮,完全潤濕轉盤表面。液膜厚度分布如圖1所示,r為水躍半徑。

參考系固定于轉盤,微團所受體積力為

式中,為工質密度,為轉盤轉速,為液膜與轉盤之間相對速度,為重力加速度。

將式(1)代入Navier-Stokes(N-S)方程并進行量綱1化處理化簡方程,引入薄膜特征厚度0以及特征長度0[29,34]

式中,為工質運動黏度。基于相似定律,速度尺度和壓力尺度分別為

0=0=l0,0=0,0=220(3)

結合式(1)~式(3)處理連續性方程和N-S方程可得式(4)~式(7)。、、、分別為量綱1速度和壓力,=0/0,=20/。

(5)

(6)

1.2 薄膜近似

澆注口尺寸一般滿足in?0,薄膜特征厚度0相比其特征長度0為微量,可采用薄膜近似[30],即=0/0?1,式(4)~式(7)中的2和4項均可忽略;采用積分平均法將式(4)~式(7)中的物理量按照處理。變換后,式(4)~式(7)進一步簡化為式(8)。為便于理解,式(8)中已略去所有物理量的上劃線。

式中,=s+,為等效Froude數,s為大氣壓強。

1.3 邊界條件

薄膜流表面需滿足運動學以及應力平衡邊界條件。Kim等[29]認為r處的切向速度可采用式(9)計算,但并未給出r的計算方法。而工質澆注條件(澆注尺寸或澆注速度)均有可能影響薄膜流特性,可定義r=in/0。若忽略工質沖擊轉盤引起的微弱能量損失,動能由軸向完全轉變為徑向,r處平均徑向速度可以近似為式(10),=r處的液膜厚度應滿足連續性方程式(11)。綜上,式(8)需滿足邊界條件式(9)~式(11)。

(10)

(11)

1.4 模型驗證

薄膜流厚度一般在數百微米,測量方法一般為探針法[28]和電阻法[26,32]。圖2、圖3為式(8)預測值分別與Leshev等[28]的探針法和Miyasaka[26]的電阻法測量值對比,工質均為甘油水溶液。Miyasaka[26]的實驗條件為=1200 kg·m-3,=(3~8)×10-3Pa·s,進口in為200,=0.1~0.3。可以看出,薄膜的厚度分布與實驗測量吻合得很好。Miyasaka在實驗中觀察到了水躍現象,從圖3可以看出,式(8)準確地預測到了水躍的發生半徑,并且薄膜厚度與測量值也完全吻合。

Burns等[32]采用電阻法對薄膜平均徑向速度進行了測量,發現徑向速度分布可分為澆注區、加速區以及同步區。由圖4可以看出,式(8)對薄膜流徑向速度的預測值與測量值吻合得很好,完全可以體現速度的三分區特征。

2 結果與討論

2.1 液膜厚度分布

由式(8)~式(11)可以看出,薄膜流的流動特性與等效Froude數、澆注口半徑r以及有關。圖5為=0.01與r=0.1時不同的薄膜厚度(或自由液面形狀)沿徑向的分布。可以看出,隨著由0增至0.04,薄膜厚度自澆注口后先呈微弱下降,而后在水躍半徑處匯聚,在水躍之后的薄膜流同樣先微弱下降而后完全重合。水躍半徑r則不受控制始終維持在0.85左右,另外對于>1.5的轉盤區域,薄膜流自由液面完全重合。因此只對水躍附近的薄膜流產生了微弱影響,而不會改變水躍半徑以及遠離轉盤軸心的薄膜流特性,忽略后在=0.48處引入的最大誤差小于3%。從式(8)可以看出,等效Froude數對薄膜流的影響主要通過薄膜流體壓強來體現,但由于薄膜流厚度一般在數百微米,其內部壓強近似等于自由液面的大氣壓強而沒有顯著改變,即使Froude數在0.25~∝范圍時也不會對薄膜流自由液面形狀產生顯著影響。

圖6為=0與r=0.1時不同值的薄膜厚度徑向分布。隨著由0.001增至0.2,薄膜流從無水躍現象逐漸過渡到明顯的發生水躍,水躍半徑最終穩定于=0.85的轉盤位置;當<1.5時,薄膜厚度呈顯著下降直至>0.05后趨于穩定,而>1.5轉盤區域的薄膜流則不受影響,完全重合。由于相同澆注條件下主要體現了轉速對薄膜流的影響,因此水躍現象是否發生與轉盤轉速直接相關。圖7所示為=0與=0.01時不同澆注口尺寸r的薄膜厚度分布。其變化規律與圖6類似,澆注口半徑r<0.1時水躍現象開始于=0.85處出現,當r降至0.02以下時,水躍穩定且自由液面開始趨于一致。Wang等[30]采用CFD方法對不同澆注口尺寸in時的薄膜流特性進行了數值分析,其中=45 ml·s-1,=2590 kg·m-3,=0.7 Pa·s,=0.478 N·m-1,=209.4 r·s-1。從圖8同樣可以看出,相同工況條件下澆注口尺寸顯著影響水躍是否發生和水躍前的薄膜厚度,而對水躍后區域則無影響。

結合圖6和圖7可以看出,相同工況下,工質進口條件(澆注口尺寸、澆注流量)以及轉盤轉速直接決定了水躍現象是否發生,并且影響<1.5轉盤區域的薄膜流特性。但水躍的發生半徑始終維持在=0.85左右,>1.5轉盤區域的薄膜流特性也不受影響。原因主要在于以=1.5為界,>1.5薄膜流徑向流動特性主要受離心力和黏性力控制,而<1.5的區域作用力主要為慣性力和黏性力。對于給定的澆注流量,澆注口尺寸縮小引起工質沖擊轉盤軸心速度的顯著上升,在黏性力作用下澆注區速度急劇縮小,工質由急流轉變為緩流,部分動能轉化為位能而導致水躍現象的發生。而增大轉盤轉速時,慣性力影響區域減小,離心作用區域加大而引起薄膜快速進入加速區而無水躍現象發生。

2.2 同步區半徑

Burns等[32]采用電阻法對薄膜流的平均徑向速度進行測量,發現徑向速度分布可分為澆注區、加速區以及同步區。以圖9中r=0.15為例,工質澆注至轉盤軸心后由于黏性力作用引起速度下降,在水躍附近速度降至最低,而后離心力引起薄膜切向速度逐漸與轉盤一致,徑向速度也同時上升直至=1.53處達到最大,最終進入同步區。但從圖9也可以看出,隨著澆注口尺寸的縮小或薄膜流慣性力的增大,澆注區范圍逐漸加大,工質速度因黏性減速效果也越明顯,而加速區范圍顯著縮小,薄膜無明顯加速現象,但不影響同步區尺寸。圖10所示為根據文獻測量數據得到的同步區半徑回歸方程。可以看出,量綱1同步區半徑=1.47,與本文的理論分析結果基本一致,這也驗證了本文的研究結論。

3 結 論

本文以轉盤表面薄膜流穩態流動特性為研究對象,采用量綱1分析與積分平均法化簡N-S方程,建立了量綱1薄膜流理論模型,研究了等效Froude數、量綱1和澆注口尺寸r對量綱1薄膜厚度分布的影響,探討了薄膜流水躍半徑以及同步區半徑與操作參數的一般規律。得到以下結論。

(1)等效Froude數在0~0.04范圍內時,對薄膜流的影響主要通過薄膜流體壓強來體現,但由于薄膜流厚度一般在數百微米,其內部壓強近似等于自由液面的大氣壓強而沒有顯著改變,不會對量綱1薄膜厚度產生明顯影響,進而水躍半徑與無關。

(2)量綱1以及澆注口尺寸r直接決定了水躍現象是否發生。澆注口尺寸縮小引起工質沖擊轉盤軸心速度顯著上升,黏性力引起澆注區速度由急流急劇向緩流轉變,流體位能增大而有利于水躍發生;而慣性力影響區域隨轉盤轉速增大而減小,離心作用區域加大而引起薄膜快速進入加速區而無水躍現象發生。增大或者縮小r均有利于水躍的發生,量綱1水躍半徑始終維持在=0.85。

(3)>1.53區域的薄膜流徑向流動特性主要受離心力和黏性力控制,而<1.53的區域作用力主要為慣性力和黏性力,由此導致薄膜的平均徑向速度呈現典型三分區特征,縮小澆注口尺寸r將引起澆注區范圍加大而加速區范圍顯著縮小,導致薄膜無明顯加速現象而直接進入同步區,但同步區半徑始終維持在=1.53,與實驗測量結果=1.47基本一致。

符 號 說 明

F——力,N g——重力加速度,m·s-2 h——薄液膜厚度,m l0——薄膜特征長度,m p——量綱1薄膜流壓強 ps——大氣壓強,Pa Q——工質流量,m3·s-1 r——量綱1轉盤半徑 ri——量綱1進口半徑 rin,rj,rs——分別為工質進口半徑、水躍半徑與同步區半徑,m u,v,w——分別為薄膜徑向、切向與軸向速度,m·s-1 u0,v0,w0——分別為量綱1薄膜徑向、切向與軸向速度 d0——薄膜特征厚度,m e——薄膜特征厚度與長度比值 μ——工質動力黏度,Pa·s ρ——工質密度,kg·m-3 υ——工質運動黏度,m2·s-1 w——轉盤轉速,rad·s-1 下角標 i——澆注進口 j——水躍半徑 0——特征參數

References

[1] BATISTELLA C B, MORAES E B, MACIEL W M R. Comparing centrifugal and falling film molecular stills using reflux and cascade for fine chemical separations[J]. Computers & Chemical Engineering, 1999, 23: S767-S770.

[2] BOIARKINA I, NORRIS S, PATTERSON D A. Investigation into the effect of flow structure on the photocatalytic degradation of methylene blue and dehydroabietic acid in a spinning disc reactor[J]. Chemical Engineering Journal, 2013, 222(4): 159-171.

[3] 韓進, 朱彤, 今井剛, 等. 基于高速轉盤法的剩余污泥可溶化處理[J]. 化工學報, 2008, 59(2): 478-483. HAN J, ZHU T, JIN J G,. Solubilization of excess sludge by high speed rotary disk[J]. Journal of Chemical Industry and Engineering (China), 2008, 59(2): 478-483.

[4] 李丹, 張大皓, 譚天偉. 轉盤反應器固定根霉重復批次發酵生產脂肪酶[J]. 化工學報, 2007, 58(9): 2347-2351. LI D, ZHANG D H, TAN T W. Product ion of lipase by repeated batch fermentation with immobilizedusing rotating biological contactor[J]. Journal of Chemical Industry and Engineering (China), 2007, 58(9): 2347-2351.

[5] 任南琪, 周顯嬌, 郭婉茜, 等. 染料廢水處理技術研究進展[J]. 化工學報, 2013, 64(1): 84-94. REN N Q, ZHOU X J, GUO W Q,. A review on treatment methods of dye wastewater[J]. CIESC Journal, 2013, 64(1): 84-94.

[6] 孫志成, 韓進, 張明楊, 等. 機械轉盤聯合超聲破解剩余污泥[J]. 化工學報, 2016, 67(12): 5229-5236. SUN Z C, HAN J, ZHANG M Y,. Disintegration of excess sludge by mechanical rotary disk combination with ultrasonication[J]. CIESC Journal, 2016, 67(12): 5229-5236.

[7] BOIARKINA I, NORRIS S, PATTERSON D A. Investigation into the effect of flow structure on the photocatalytic degradation of methylene blue and dehydroabietic acid in a spinning disc reactor[J]. Chemical Engineering Journal, 2013, 222 (4): 159-171.

[8] 魏冰, 張巧玲, 劉有智, 等. 旋轉盤反應器光催化降解苯酚廢水動力學[J]. 環境工程學報, 2016, 10(3): 1305-1309. WEI B, ZHANG Q L, LIU Y Z,. Kinetics of phenol by photocatalysis degradation in spinning disk reactor[J]. Chinese Journal of Environmental Engineering, 2016, 10(3): 1305-1309.

[9] WANG D X, LING X, PENG H,. Experimental investigation of ligament formation dynamics of thin viscous liquid film at spinning disk edge[J]. Industrial & Engineering Chemistry Research, 2016, 55(34): 9267-9275.

[10] 吳君軍, 王宏, 朱恂, 等. 轉盤離心粒化中絲狀成粒特性[J]. 化工學報, 2015, 66(7): 2474-2480. WU J J, WANG H, ZHU X,. Characteristic of ligament in centrifugal granulation by spinning disc[J]. CIESC Journal, 2015, 66(7): 2474-2480.

[11] AHMED M, YOUSSEF M S. Influence of spinning cup and disk atomizer configurations on droplet size and velocity characteristics[J]. Chemical Engineering Science, 2014, 107(14): 149-157.

[12] ZHANG H, WANG H, ZHU X,. A review of waste heat recovery technologies towards molten slag in steel industry[J]. Applied Energy, 2013, 112(4): 956-966.

[13] BARATI M, ESFAHANI S, UTIGARD T A. Energy recovery from high temperature slags[J]. Energy, 2011, 36(9): 5440-5449.

[14] WANG D X, LING X, PENG H. Simulation of ligament mode breakup of molten slag by spinning disk in the dry granulation process[J]. Applied Thermal Engineering, 2015, 84: 437-447.

[15] JACOBSEN N C, HINRICHSEN O. Micromixing efficiency of a spinning disk reactor[J]. Industrial & Engineering Chemistry Research, 2012, 51(36): 11643-11652.

[16] AHMED M, YOUSSEF M S. Characteristics of mean droplet size produced by spinning disk atomizers[J]. Journal of Fluid Engineering, 2012, 134(7): 1-9.

[17] BOIARKINA I, NORRIS S, PATTERSON D A. The case for the photocatalytic spinning disc reactor as a process intensification technology: comparison to an annular reactor for the degradation of methylene blue[J]. Chemical Engineering Journal, 2013, 225(3): 752-765.

[18] CAPRARIIS B D, RITA D M, STOLLER M,. Reaction-precipitation by a spinning disc reactor: influence of hydrodynamics on nanoparticles production[J]. Chemical Engineering Science, 2012, 76(28): 73-80.

[19] 閻維平, 葉學民, 李洪濤. 液體薄膜流的流動和傳熱特性[J]. 華北電力大學學報, 2005, 32(1): 59-65. YAN W P, YE X M, LI H T. Hydrodynamics and heat transfer of liquid thin films[J]. Journal of North China Electric Power University, 2005, 32(1): 59-65.

[20] CHIANG C, CHANG M, LIU H,. Process intensification in the production of photocatalysts for solar hydrogen generation[J]. Industrial & Engineering Chemistry Research, 2012, 51(14): 5207-5215.

[21] PENG H, WANG N, WANG D X,. Experimental study on critical characteristics of liquid atomization by spinning disk[J]. Industrial & Engineering Chemistry Research, 2016, 55(21): 6175-6185.

[22] 王東祥, 凌祥, 彭浩. 轉盤離心粒化熔渣液膜流動特性數值模擬研究[J]. 南京工業大學學報(自然科學版), 2015, 37(3): 67-73. WANG D X, LING X, PENG H. Numerical simulation of film flow characteristics of molten slag on spinning disk in centrifugal atomization process[J]. Journal of Nanjing Tech University (Natural Science Edition), 2015, 37(3): 67-73.

[23] WANG D X, LING X, PENG H,. Experimental investigation of ligament formation dynamics of thin viscous liquid film at spinning disk edge[J]. Industrial & Engineering Chemistry Research, 2016, 55(34): 9267-9275.

[24] 何先琰, 王宏, 朱恂, 等. 鉛錫合金熔融顆粒風冷相變換熱特性實驗研究[J]. 工程熱物理學報, 2015, 36(8): 1748-1751.HE X Y, WANG H, ZHU X,. Experiment study on air-cooled phase change heat transfer characteristics of Sn-Pb alloy droplets[J]. Journal of Engineering Thermophysics, 2015, 36(8): 1748-1751.

[25] EMSLIE A G, BONNER F T, PECK L G. Flow of a viscous liquid on a rotating disk[J]. Journal of Applied Physics, 1958, 29(5): 858-862.

[26] MIYASAKA Y. On the flow of a viscous free boundary jet on a rotating disk[J]. Bulletin of the JSME, 1974, 40(331): 1461-1468.

[27] THOMAS S, FAGHRI A, HANKEY W. Experimental analysis and flow visualization of a thin liquid film on a stationary and rotating disk[J]. Journal of Fluids Engineering, 1991, 113(1): 73-80.

[28] LESHEV I, PEEV G. Film flow on a horizontal rotating disk[J]. Chemical Engineering and Processing, 2003, 42 (11): 925-929.

[29] KIM T S, KIM M U. The flow and hydrodynamic stability of a liquid film on a rotating disc[J]. Fluid Dynamics Research, 2009, 41(3): 475-479.

[30] WANG D X, LING X, PENG H. Theoretical analysis of free-surface film flow on the rotary granulating disk in waste heat recovery process of molten slag[J]. Applied Thermal Engineering, 2014, 63(1): 387-395.

[31] WOODS W P. The hydrodynamics of thin liquid films flowing over a rotating disk[D]. United Kingdom: University of Newcastle, 1995.

[32] BURNS J R, RAMSHAW C, JACHUCK R J. Measurement of liquid film thickness and the determination of spin-up radius on a rotating disc using an electrical resistance technique[J]. Chemical Engineering Science, 2003, 58 (11): 2245-2253.

[33] SISOEV G M, MATER O K, LAWRENCE C J. Axisymmetric wave regimes in viscous liquid film flow over a spinning disk[J]. Journal of Fluid Mechanics, 2003, 495(495): 385-411.

[34] PRIELING D, STEINER H. Unsteady thin film flow on spinning disks at large ekman numbers using an integral boundary layer method[J]. International Journal of Heat and Mass Transfer, 2013, 65(7): 10-22.

Numerical simulation of stable flow dynamics of viscous film on spinning disk surface

WANG Dongxiang1,2, LING Xiang2, PENG Hao2, YANG Xinjun1,2, CUI Zhengwei1

(1Jiangsu Key Laboratory of Advanced Food Manufacturing Equipment & Technology, School of Mechanical Engineering, Jiangnan University, Wuxi 214122, Jiangsu, China;2Jiangsu Key Laboratory of Process Enhancement and New Energy Equipment Technology, Nanjing Tech University, Nanjing 211800, Jiangsu, China)

Film flow on spinning disk surface exists extensively in chemical engineering operations, such as centrifugal graining, molecular distillation and spinning disk reactors. Flow dynamics of the film flow has a major impact on graining, reaction, heat and mass transfer rate. Theoretical model of stable film flow dynamics was established and verified by comparison of experimental results and numerical simulations. The effect of equivalent Froude number, dimensionless ratio of characteristic thickness over length,and casting sizeron dimensionless film thickness distribution was studied, which a model of dimensionless hydraulic jump and synchronized zone radius was derived and verified by experiments. The results show that equivalent Froude number has little effect on film thickness distribution and occurrence of hydraulic jump phenomenon mainly depends uponandr. Increasingor decreasingrleads to appearance of hydraulic jump with a dimensionless hydraulic jump radius always at=0.85. Mean radial velocity exhibited features of typical three zone distribution. Reducing casting size would extend injection zone and shrink acceleration zone, such that film flow goes directly to synchronized zone at dimensionless synchronized radius of 1.53 without evident acceleration. The results will provide theoretical reference for design and optimization of spinning disk reactors and centrifugal pelletizers.

fluid mechanics; model; experimental validation; flow regimes; spinning disk; steady flow

10.11949/j.issn.0438-1157.20170040

TQ 021.1

A

0438—1157(2017)06—2321—07

凌祥。

王東祥(1985—),男,博士,講師。

國家自然科學基金項目(51606086,51406078)。

2017-01-10收到初稿,2017-02-16收到修改稿。

2017-01-10.

Prof. LING Xiang, xling@njtech.edu.cn

supported by the National Natural Science Foundation of China (51606086, 51406078).

主站蜘蛛池模板: 波多野结衣一区二区三区AV| 国产女人水多毛片18| 18禁黄无遮挡网站| 亚洲精品国产精品乱码不卞| 成人国产精品一级毛片天堂 | 婷婷99视频精品全部在线观看| 亚洲高清中文字幕在线看不卡| 国产一级裸网站| 国产亚洲欧美在线专区| 日本日韩欧美| 欧美在线视频不卡| 毛片基地美国正在播放亚洲 | 亚洲国产av无码综合原创国产| 美女被躁出白浆视频播放| 国产精品亚洲综合久久小说| 久久亚洲精少妇毛片午夜无码| 永久毛片在线播| 国产黄在线免费观看| 日韩欧美国产三级| 综1合AV在线播放| 国产在线自在拍91精品黑人| 色婷婷电影网| 中文字幕无码中文字幕有码在线| 91po国产在线精品免费观看| 一级成人a毛片免费播放| 日韩在线网址| 一级一毛片a级毛片| 成人国产精品一级毛片天堂| 伊人国产无码高清视频| 日本国产在线| 成人在线观看不卡| 国产乱人伦精品一区二区| 67194成是人免费无码| 日韩欧美在线观看| 午夜视频日本| 国产成人午夜福利免费无码r| 国产成人高精品免费视频| 粗大猛烈进出高潮视频无码| 91精品啪在线观看国产60岁 | 中文字幕亚洲精品2页| 亚洲国产中文综合专区在| 欧洲亚洲一区| 日韩在线第三页| 69国产精品视频免费| 性喷潮久久久久久久久| 亚洲av日韩av制服丝袜| 久久精品国产999大香线焦| 精品亚洲国产成人AV| 成人噜噜噜视频在线观看| 欧美日韩精品在线播放| 国模粉嫩小泬视频在线观看| 亚洲丝袜第一页| 国产麻豆va精品视频| 国产成人一级| 欧美精品另类| 国产男人天堂| 久久96热在精品国产高清| 精品视频在线一区| 精品免费在线视频| 五月婷婷精品| www.亚洲一区| 亚洲AV无码不卡无码 | 最近最新中文字幕免费的一页| 国产人人射| 欧美成人二区| 无码日韩视频| 国内精品伊人久久久久7777人| 国产福利观看| 国产成人综合久久| 91探花在线观看国产最新| 国产精品视频猛进猛出| 成人福利在线看| 在线亚洲小视频| 色精品视频| 国产精品自在在线午夜| 激情综合激情| 91久久精品日日躁夜夜躁欧美| 亚洲第一成网站| 好紧太爽了视频免费无码| 国产精品男人的天堂| 成人av专区精品无码国产| 综合色区亚洲熟妇在线|