• <tr id="yyy80"></tr>
  • <sup id="yyy80"></sup>
  • <tfoot id="yyy80"><noscript id="yyy80"></noscript></tfoot>
  • 99热精品在线国产_美女午夜性视频免费_国产精品国产高清国产av_av欧美777_自拍偷自拍亚洲精品老妇_亚洲熟女精品中文字幕_www日本黄色视频网_国产精品野战在线观看 ?

    中空液滴碰撞水平壁面數(shù)值分析?

    2017-07-31 01:37:48鄭志偉李大樹(shù)仇性啟崔運(yùn)靜
    物理學(xué)報(bào) 2017年1期
    關(guān)鍵詞:實(shí)心中空浸潤(rùn)性

    鄭志偉 李大樹(shù) 仇性啟? 崔運(yùn)靜

    1)(中國(guó)石油大學(xué)(華東)化學(xué)工程學(xué)院,青島 266580)2)(中海油研究總院,北京 100028)3)(中國(guó)石油大學(xué)(華東)機(jī)電工程學(xué)院,青島 266580)(2016年6月4日收到;2016年9月26日收到修改稿)

    中空液滴碰撞水平壁面數(shù)值分析?

    鄭志偉1)李大樹(shù)2)?仇性啟1)?崔運(yùn)靜3)

    1)(中國(guó)石油大學(xué)(華東)化學(xué)工程學(xué)院,青島 266580)2)(中海油研究總院,北京 100028)3)(中國(guó)石油大學(xué)(華東)機(jī)電工程學(xué)院,青島 266580)(2016年6月4日收到;2016年9月26日收到修改稿)

    采用耦合水平集-體積分?jǐn)?shù)法并綜合考慮傳熱及接觸熱阻作用建立了中空液滴碰撞水平壁面數(shù)值模型,并驗(yàn)證了模型的可靠性.通過(guò)分析計(jì)算結(jié)果,獲得了中空液滴與實(shí)心液滴撞壁的動(dòng)力學(xué)特征差異,揭示了中空液滴撞壁流動(dòng)傳熱機(jī)理和中心射流形成機(jī)制,探索了碰撞速度和壁面浸潤(rùn)性對(duì)中空液滴撞壁動(dòng)力學(xué)和傳熱特性的影響.研究表明:中空液滴撞壁后中心射流特征明顯,并伴隨有射流收縮和液殼破碎等現(xiàn)象.中空液滴內(nèi)部壓力梯度是液滴鋪展、中心射流產(chǎn)生和發(fā)展的主要原因;撞壁過(guò)程中中心射流表面溫度分布較為均勻,破碎液殼表面溫度分布波動(dòng)較大.碰撞速度與中空液滴撞壁最大鋪展系數(shù)的相關(guān)性較小,但其對(duì)無(wú)量綱射流長(zhǎng)度和壁面平均熱流密度的影響較大;壁面浸潤(rùn)性與中空液滴撞壁后期鋪展系數(shù)的相關(guān)性較大,但其對(duì)無(wú)量綱射流長(zhǎng)度和壁面平均熱流密度的影響較小.

    中空液滴,液滴撞壁,中心射流

    1引 言

    液滴撞壁是自然界和工業(yè)技術(shù)領(lǐng)域中普遍存在的現(xiàn)象,如噴霧燃燒、噴霧干燥、噴墨打印和熱噴涂等[1,2].考慮到液滴撞壁過(guò)程流動(dòng)機(jī)理復(fù)雜,相界面拓?fù)浣Y(jié)構(gòu)變化較大,研究液滴撞壁過(guò)程不僅能夠豐富自由界面氣液兩相流動(dòng)理論,而且對(duì)噴霧技術(shù)的工程應(yīng)用有重要指導(dǎo)意義[3].

    近年來(lái),國(guó)內(nèi)外學(xué)者對(duì)實(shí)心液滴撞壁研究做出了許多貢獻(xiàn).Rioboo等[4]采用高速攝像儀對(duì)實(shí)心液滴碰撞水平干壁面現(xiàn)象進(jìn)行了可視化研究,發(fā)現(xiàn)實(shí)心液滴撞壁后主要呈現(xiàn)出黏附、鋪展、破碎和飛濺等特征.畢菲菲等[5]實(shí)驗(yàn)觀測(cè)了不同實(shí)心液滴撞擊固體表面的動(dòng)力學(xué)形態(tài)變化,并探討了撞擊參數(shù)對(duì)液滴形態(tài)的影響,研究表明液滴黏度和表面張力的共同作用決定了液滴的振蕩特性,液滴最大鋪展因子隨撞擊速度的增大而增大,且液滴達(dá)到最大鋪展因子所需時(shí)間與表面張力密切相關(guān).李大樹(shù)等[6]采用液體體積(VOF)法數(shù)值分析了柴油液滴碰撞水平干壁面的過(guò)程,揭示了液滴撞壁流動(dòng)傳熱及破碎機(jī)制,并建立了柴油液滴碰撞水平干壁面最大鋪展系數(shù)半經(jīng)驗(yàn)理論解析模型.Kown[7]采用Level set法對(duì)液滴碰撞水平壁面進(jìn)行了數(shù)值模擬,獲得了液滴撞壁特性隨韋伯?dāng)?shù)的變化規(guī)律及液滴飛濺的臨界條件.Yokoi[8]采用耦合水平集-體積分?jǐn)?shù)法(CLSVOF)研究了水滴碰撞水平干壁面的飛濺現(xiàn)象,研究發(fā)現(xiàn)壁面前進(jìn)接觸角是液滴撞壁飛濺的關(guān)鍵因素,壁面前進(jìn)角越大,液滴越容易飛濺.

    現(xiàn)有液滴撞壁研究主要集中在實(shí)心液滴,然而,在一些液滴撞壁控制技術(shù)中,如熱噴涂領(lǐng)域,中空液滴與實(shí)心液滴相比存在較為明顯的優(yōu)勢(shì),主要體現(xiàn)在[9,10]:1)相同初始直徑下,中空液滴質(zhì)量較小,容易獲得較大的初始動(dòng)能;2)中空液滴導(dǎo)熱率可以通過(guò)改變液殼厚度以及液滴內(nèi)部氣體物性加以有效控制.然而,目前有關(guān)中空液滴撞壁研究十分缺乏,且集中在國(guó)外,在國(guó)內(nèi)尚未見(jiàn)相關(guān)報(bào)道.Gulyaev等[11]采用高速攝像儀首次觀測(cè)到甘油中空液滴碰撞水平壁面的法向中心射流現(xiàn)象.隨后,進(jìn)一步對(duì)甘油中空液滴撞壁特性進(jìn)行了分析,獲得了鋪展系數(shù)和中心射流速度隨不同雷諾數(shù)的變化[12].Solonenko等[13]和Shinoda等[14]實(shí)驗(yàn)對(duì)比分析了氧化鋯中空液滴和實(shí)心液滴的撞壁沉積現(xiàn)象,研究表明與實(shí)心液滴相比,中空液滴撞壁后沉積液膜更為平滑.在數(shù)值研究方面,Kumar等[15]采用VOF法數(shù)值分析了金屬錫中空液滴撞壁過(guò)程,研究發(fā)現(xiàn)中空液滴內(nèi)部氣相分布與其撞壁動(dòng)力學(xué)形態(tài)存在較大相關(guān)性.

    目前,現(xiàn)有研究對(duì)中空液滴撞壁動(dòng)力學(xué)特征的分析還存在較大不足,液滴撞壁流動(dòng)、傳熱機(jī)制尚不明確,撞壁特性規(guī)律仍有待進(jìn)一步探索[16].同時(shí),現(xiàn)有液滴撞壁研究多忽略傳熱作用,但液滴熱物性會(huì)對(duì)其動(dòng)力學(xué)行為存在一定影響.此外,在涉及傳熱的撞壁過(guò)程中,接觸熱阻是描述液固耦合的重要參數(shù)[17],而目前采用接觸熱阻的液滴撞壁模型還極其缺乏.因此,本文采用CLSVOF法并綜合考慮傳熱及接觸熱阻作用建立液滴碰撞水平壁面數(shù)值模型,對(duì)比分析中空液滴和實(shí)心液滴撞壁動(dòng)力學(xué)特征的差異,揭示中空液滴撞壁流動(dòng)的傳熱機(jī)理和中心射流的形成機(jī)制,探索碰撞速度和壁面浸潤(rùn)性對(duì)中空液滴撞壁動(dòng)力學(xué)和傳熱特性的影響.

    2數(shù)值模型

    中空液滴碰撞水平壁面的幾何模型如圖1所示.定義液滴的參數(shù)如下:D0是液滴初始直徑,Ds是液滴鋪展直徑,U0是液滴碰撞速度,δp是液殼厚度,則液滴的鋪展系數(shù)f和無(wú)量綱液殼厚度δ?p分別可以表示成

    數(shù)值模型采用甘油中空液滴作為研究對(duì)象,液滴初始直徑D0=5.25mm,無(wú)量綱液殼厚度=0.082,初始溫度Td=325 K,以一定的碰撞速度U0垂直沖擊恒溫不銹鋼壁面,壁面與環(huán)境溫度均為300 K,環(huán)境壓力為1個(gè)大氣壓.針對(duì)液滴的鋪展和收縮過(guò)程,計(jì)算模型分別采用前進(jìn)和后退接觸角以求解液滴與壁面間的液-固潤(rùn)濕作用.此外,數(shù)值模型考慮Marangoni效應(yīng)以精準(zhǔn)求解相界面的熱應(yīng)力,采用熱物性參數(shù)替代常數(shù)參數(shù),并考慮壁面接觸熱阻對(duì)液滴撞壁過(guò)程的影響.利用計(jì)算模型對(duì)控制方程進(jìn)行求解,并假設(shè)撞壁過(guò)程中熱量傳遞為空氣對(duì)流換熱以及液滴與基板表面間的導(dǎo)熱和對(duì)流換熱.

    圖1 中空液滴碰撞水平壁面幾何模型 (a)撞壁前;(b)撞壁后Fig.1. Geometry model of hollow droplet impact on a fl at surface:(a)Characteristics before impact;(b)characteristics after impact.

    根據(jù)CLSVOF方法,連續(xù)性方程、動(dòng)量方程和能量方程分別可以表示成

    式中,v是速度矢量,F是表面張力源項(xiàng),p是壓力,g是重力矢量,?是Level set函數(shù),ρ(?)是計(jì)算單元內(nèi)密度,μ(?)是計(jì)算單元內(nèi)動(dòng)力黏度,cp是計(jì)算單元內(nèi)比熱容,T是計(jì)算單元內(nèi)溫度,λ是計(jì)算單元內(nèi)導(dǎo)熱系數(shù).動(dòng)量方程中表面張力采用Brackbill等[18]提出的連續(xù)表面力(CSF)模型來(lái)求解,并通過(guò)引入Heaviside函數(shù)進(jìn)行連續(xù)化處理.表面張力源項(xiàng)可以表示成

    式中,κ(?)是曲率; σ是表面張力系數(shù);H(?)是Heaviside函數(shù);a是相界面處過(guò)渡區(qū)域厚度,a=1.5h(h是最小網(wǎng)格尺寸)[19?21].連續(xù)化后的密度和黏度分別為

    式中,下標(biāo)g和l分別表示氣相和液相.

    由Marangoni效應(yīng)引起的應(yīng)力τ可以表示為

    接觸角模型以壁面邊界條件的形式施加到數(shù)值模型中,可以表示成

    式中,θdy是壁面接觸角,Ucl是接觸線速度,θe是靜態(tài)接觸角,θr和θb分別是前進(jìn)和后退接觸角.

    計(jì)算采用二維軸對(duì)稱模型,計(jì)算區(qū)域取15mm×30mm,為提高相界面的計(jì)算精度,對(duì)壁面附近區(qū)域網(wǎng)格進(jìn)行自適應(yīng)加密,經(jīng)網(wǎng)格無(wú)關(guān)性檢驗(yàn),選取網(wǎng)格尺寸為0.03mm進(jìn)行劃分,既可保證計(jì)算精度,又可節(jié)省計(jì)算時(shí)間和成本,網(wǎng)格數(shù)為349800.采用有限體積法對(duì)控制方程進(jìn)行離散,壓力速度耦合采用SIMPLE方法,壓力求解采用PRESTO!方法,計(jì)算單元液相體積分?jǐn)?shù)采用CICSAM[22]方法離散求解,控制方程采用QUICK格式進(jìn)行離散,以減少假擴(kuò)散提高精度,時(shí)間步長(zhǎng)Δt=10?6s.本文經(jīng)收斂性檢驗(yàn)選取Δt內(nèi)迭代次數(shù)為100,在該次數(shù)內(nèi)可滿足收斂精度要求,計(jì)算結(jié)果收斂.

    3計(jì)算模型驗(yàn)證

    為驗(yàn)證所建立中空液滴撞壁數(shù)值模型的準(zhǔn)確性,將甘油中空液滴撞壁數(shù)值結(jié)果與文獻(xiàn)[11]中實(shí)驗(yàn)觀測(cè)結(jié)果進(jìn)行對(duì)比,如圖2所示,其中,液滴初始直徑D0=5.25mm,碰撞速度U0=5.94m/s,無(wú)量綱液殼厚度=0.082.圖2(a)是不同時(shí)刻中空液滴撞壁運(yùn)動(dòng)形態(tài)實(shí)驗(yàn)和數(shù)值模擬結(jié)果,圖2(b)是實(shí)驗(yàn)測(cè)量和數(shù)值模擬中心射流長(zhǎng)度隨時(shí)間的變化.

    圖2 實(shí)驗(yàn)和模擬中空液滴撞壁動(dòng)力學(xué)特征對(duì)比 (a)運(yùn)動(dòng)形態(tài)對(duì)比;(b)中心射流長(zhǎng)度對(duì)比Fig.2.Comparison of experiment and simulation results for hollow droplet impact:(a)Morphology of droplet;(b)length of jet.

    從圖2(a)可以看出,在整個(gè)撞壁過(guò)程中,數(shù)值模擬與實(shí)驗(yàn)觀測(cè)中空液滴運(yùn)動(dòng)形態(tài)較為一致.液滴接觸壁面后迅速鋪展,液滴高度不斷下降,液滴動(dòng)能逐漸轉(zhuǎn)化為鋪展過(guò)程中表面能和黏性耗散能,當(dāng)液滴鋪展到一定程度時(shí),碰撞中心處產(chǎn)生垂直于壁面豎直向上的射流(0.67ms),中心射流高度隨時(shí)間逐漸增大,當(dāng)其到達(dá)液滴頂部時(shí),對(duì)液殼產(chǎn)生較大的沖擊作用,液殼出現(xiàn)破碎(1ms).隨后,在慣性力和表面張力共同作用下,中心射流外邊緣迅速收縮,射流與壁面接觸面積有所減小,液滴逐漸達(dá)到完全射流狀態(tài)(3ms).可見(jiàn),中空液滴撞壁后主要呈現(xiàn)出鋪展、中心射流和液殼破碎等現(xiàn)象.從圖2(b)可以看出,不同時(shí)刻實(shí)驗(yàn)和數(shù)值模擬中心射流長(zhǎng)度L較為接近,均隨時(shí)間近似線性遞增.此外,實(shí)驗(yàn)與模擬液膜鋪展結(jié)果有所差異,分析認(rèn)為,與現(xiàn)有數(shù)字圖像處理技術(shù)和相界面追蹤方法的局限性有關(guān).

    4結(jié)果分析

    4.1 中空液滴與實(shí)心液滴撞壁對(duì)比

    圖3是數(shù)值模擬中空液滴與實(shí)心液滴以4m/s碰撞速度沖擊水平壁面的動(dòng)力學(xué)形態(tài).其中,中空液滴初始直徑D0=5.25mm,無(wú)量綱液殼厚度δ?p=0.082,壁面前進(jìn)接觸角和后退接觸角分別為θr=130°和θb=150°,靜態(tài)接觸角θe=140°. 為獲得相同的初始沖擊動(dòng)能,實(shí)心液滴采用與中空液滴具有相同質(zhì)量的液滴,初始直徑為3.92mm.

    從圖3(a)可以看出,撞壁初始階段,慣性力驅(qū)使液滴高度不斷下降,液滴內(nèi)空氣腔被壓縮,液殼厚度有所減小,液滴與壁面接觸邊緣的流體分別向內(nèi)外兩側(cè)運(yùn)動(dòng).外側(cè)液體克服表面張力及黏性力作用向外鋪展,逐漸形成鋪展液膜.內(nèi)側(cè)液體流向碰撞中心逐漸匯聚形成法向中心射流(0.45—1.45ms),由于慣性力作用,中心射流高度逐漸增加,當(dāng)其到達(dá)液滴頂部時(shí),液殼產(chǎn)生破碎,且外側(cè)液膜與中心射流連接處發(fā)生斷裂(2.75ms).隨后,射流邊緣在慣性力和表面張力的共同作用下迅速向碰撞中心收縮,表面能轉(zhuǎn)化為動(dòng)能,促使中心射流繼續(xù)向上擴(kuò)展,中空液滴逐漸發(fā)展成完全射流狀態(tài)(5.75—8.75ms).當(dāng)中心射流慣性力足夠克服壁面黏附力和液滴重力作用時(shí),其底部完全剝離壁面(13.75—21.75ms).需要指出,中心射流反彈與實(shí)心液滴撞擊非浸潤(rùn)性壁面的反彈機(jī)制不同,實(shí)心液滴撞擊非浸潤(rùn)性壁面反彈現(xiàn)象是由于壁面的非浸潤(rùn)性阻礙了液滴鋪展,減小了鋪展過(guò)程中黏性耗散,從而保證液滴具有足夠動(dòng)能可使其沿壁面回縮并從壁面反彈[23,24].而中空液滴撞壁后中心射流反彈是由于撞壁點(diǎn)處液體向內(nèi)外兩側(cè)流動(dòng),促使液體在碰撞中心匯聚形成豎直向上的中心射流,并逐漸剝離壁面.

    圖3 中空液滴與實(shí)心液滴動(dòng)力學(xué)形態(tài)對(duì)比 (a)中空液滴;(b)實(shí)心液滴Fig.3.Comparison of morphology of hollow and dense droplet impact:(a)Morphology of hollow droplet;(b)morphology of dense droplet.

    從圖3(b)可以看出,實(shí)心液滴撞壁后沒(méi)有產(chǎn)生中心射流,可見(jiàn),中心射流是中空液滴與實(shí)心液滴撞壁的主要?jiǎng)恿W(xué)特征差異.撞壁初期,實(shí)心液滴鋪展特征較為明顯,其沖擊動(dòng)能逐漸轉(zhuǎn)化成鋪展動(dòng)能、表面能和黏性耗散能.2.46ms時(shí),實(shí)心液滴達(dá)到最大鋪展直徑,此時(shí),液滴鋪展邊緣三相接觸線速度減小到零,在表面張力的作用下,三相接觸線向碰撞中心回縮.但由于甘油黏度較大,回縮過(guò)程中黏性耗散較大,導(dǎo)致液滴的回縮動(dòng)能較小.

    圖4將實(shí)心液滴和中空液滴撞壁鋪展系數(shù)表示成無(wú)量綱時(shí)間t?(t?=tU0/D0)的函數(shù).從圖4可以看出,撞壁初始階段,實(shí)心液滴鋪展系數(shù)隨無(wú)量綱時(shí)間逐漸增加到最大值,隨后有所減小,說(shuō)明實(shí)心液滴撞壁過(guò)程中產(chǎn)生了回縮現(xiàn)象.中空液滴撞壁鋪展系數(shù)變化則較為復(fù)雜,液滴接觸壁面后,鋪展系數(shù)逐漸增加.當(dāng)t?=1.1時(shí),鋪展系數(shù)達(dá)到最大值,同時(shí),液殼產(chǎn)生破碎,外側(cè)液膜與中心射流連接處發(fā)生斷裂(如圖3(a)中2.75ms),由于慣性力和表面張力的作用射流邊緣液體向碰撞中心收縮,鋪展系數(shù)迅速減小.當(dāng)t?=6.5時(shí),中空液滴鋪展系數(shù)減小為0,說(shuō)明中心射流與壁面產(chǎn)生了剝離.此外,依據(jù)中空液滴的流動(dòng)特征,取特征長(zhǎng)度D?=D0?d0(d0是中空氣體直徑),得到中空液滴局部雷諾數(shù)(Reh=49.8)小于實(shí)心液滴局部雷諾數(shù)(Rec=303.98),由此可知,中空液滴撞壁過(guò)程中黏性作用相對(duì)較大,阻礙了液膜鋪展,導(dǎo)致鋪展階段中空液滴的鋪展系數(shù)逐漸小于實(shí)心液滴.

    圖4 中空液滴與實(shí)心液滴鋪展系數(shù)對(duì)比Fig.4.Comparison of spreading factors for dense and hollow droplet impact.

    4.2 中空液滴的流動(dòng)傳熱

    為了探索中空液滴撞壁過(guò)程中流動(dòng)傳熱機(jī)理和中心射流的形成機(jī)制,圖5對(duì)中空液滴碰撞水平壁面不同時(shí)刻的壓力和速度分布進(jìn)行分析,其中,D0=5.25mm,δ?p=0.082,U0=4m/s,靜態(tài)接觸角θe=140°,液滴壓力為靜壓.

    由圖5(a)可以看出,0.45ms時(shí),由于液滴對(duì)壁面的沖擊作用,撞擊能量在液滴與壁面接觸位置A處聚集,出現(xiàn)壓力極值點(diǎn),為9269 Pa.此時(shí),鋪展邊緣B處壓力值為1639 Pa,小于A處的壓力極值,液滴內(nèi)部沿鋪展方向具有較大的壓力梯度,因此液體克服表面張力和黏性力影響,迅速向外側(cè)鋪展.由速度分布可得,此時(shí)鋪展邊緣B處速度最大(6.7m/s),大于液滴的碰撞速度(4m/s).同時(shí),從圖5(a)還可以看出,碰撞中心C處壓力(2600 Pa)也小于A處,因此,在上述壓力梯度的作用下,液滴內(nèi)部存在由A向C的回流,液滴最大回流速度為1.09m/s,小于液滴最大鋪展速度(6.7m/s).

    由圖5(b)可以看出,2.75ms時(shí),法向中心射流特征明顯,碰撞中心C處壓力最大(995 Pa),射流頂部D處壓力較小(92 Pa),射流內(nèi)部沿垂直壁面法向存在較大壓力梯度,促使中心射流繼續(xù)向上伸展.同時(shí),在慣性力和表面張力的共同作用下,射流邊緣沿壁面向碰撞中心收縮,射流動(dòng)量逐漸由徑向轉(zhuǎn)化為軸向,進(jìn)一步促進(jìn)了中心射流發(fā)展.由速度分布可以看出,射流收縮邊緣E處出現(xiàn)收縮速度最大值(1.73m/s),大于其頂部D處的法向射流速度(1.15m/s).此外,從速度分布還可以看出,液殼與中心射流間存在較大的氣流漩渦,漩渦氣體對(duì)破碎液殼邊界層產(chǎn)生較大的剪切作用,因此,促進(jìn)了二次霧化液滴的形成.

    由圖5(c)可以看出,8.25ms時(shí),液滴處于完全射流狀態(tài),碰撞中心C處壓力仍較大,因此射流繼續(xù)維持一定的法向運(yùn)動(dòng)能力,此時(shí),最大射流速度為0.93m/s.

    圖6是中空液滴碰撞水平壁面不同時(shí)刻的溫度分布,其中,液滴溫度高于壁面溫度,液滴向壁面?zhèn)鳠?

    從圖6可以看出,0.45ms時(shí),液滴表面及其所覆蓋壁面區(qū)域內(nèi)溫度分布較為均勻,碰撞點(diǎn)C處的溫度梯度(631 K/mm)大于液滴表面T處(17.24 K/mm)和鋪展邊緣B處(625 K/mm).這是由于碰撞中心C和鋪展邊緣B處以熱傳導(dǎo)為主,液滴表面T處以空氣對(duì)流傳熱為主,而鋪展邊緣B處液體質(zhì)量較少,液膜表面積較大,相對(duì)散熱較多,因此,碰撞點(diǎn)C處溫度梯度最大,液滴表面T處溫度梯度最小,而鋪展邊緣B處溫度梯度在兩者之間.2.75ms時(shí),射流表面溫度分布較為均勻,但由于射流沿壁面向碰撞中心收縮,壁面溫度分布不再均勻,同時(shí),液殼破碎產(chǎn)生二次霧化液滴,破碎液滴溫度較高,環(huán)境空氣溫度較低,破碎液滴與環(huán)境空氣之間進(jìn)行傳熱,因此破碎液殼表面溫度分布存在較大波動(dòng).此時(shí),射流具有較大軸向速度,因此射流頂部D處溫度梯度較大(69 K/mm).8.25ms時(shí),射流表面溫度分布仍較為均勻,但壁面附近溫度存在劇烈波動(dòng),這是由于壁面附近存在大量二次破碎液滴所造成的.

    圖5 (網(wǎng)刊彩色)中空液滴撞壁過(guò)程不同時(shí)刻壓力和速度分布 (a)0.45ms;(b)2.75ms;(c)8.25msFig.5.(color online)Pressure and velocity distribution of hollow droplet at di ff erent time during impact:(a)0.45ms;(b)2.75ms;(c)8.25ms.

    圖6 (網(wǎng)刊彩色)中空液滴撞壁過(guò)程不同時(shí)刻溫度分布 (a)0.45ms;(b)2.75ms;(c)8.25msFig.6.(color online)Temperature distribution of hollow droplet at di ff erent time during impact:(a)0.45ms;(b)2.75ms;(c)8.25ms.

    4.3 中空液滴撞壁特性

    圖7對(duì)比分析了不同碰撞速度下中空液滴撞擊不同浸潤(rùn)性壁面鋪展系數(shù)f隨無(wú)量綱時(shí)間t?的變化.其中,壁面靜態(tài)接觸角分別為θe=30°和θe=140°.

    從圖7可以看出,中空液滴碰撞相同浸潤(rùn)性壁面時(shí),不同碰撞速度下最大鋪展系數(shù)較為接近.當(dāng)碰撞速度分別為4,6和10m/s時(shí),最大鋪展系數(shù)分別是3,3.15和3.4.可見(jiàn),中空液滴撞壁最大鋪展系數(shù)隨碰撞速度的增加略有增大,這與實(shí)心液滴撞壁最大鋪展系數(shù)的碰撞速度效應(yīng)完全不同:實(shí)心液滴撞壁最大鋪展系數(shù)隨碰撞速度的增加顯著增大[25,26].這是由于中空液滴撞壁鋪展液膜中的液體質(zhì)量在液滴總質(zhì)量中占有比例較小,因此,不同碰撞速度下鋪展液體初始動(dòng)能較為接近,導(dǎo)致中空液滴最大鋪展系數(shù)與碰撞速度的相關(guān)性較小.

    圖7 中空液滴碰撞不同浸潤(rùn)性壁面鋪展系數(shù)的變化Fig.7.Spreading factor of hollow droplet impact on surface with di ff erent wettabilities.

    從圖7還可以看出,碰撞速度相同時(shí),中空液滴撞壁初始階段不同浸潤(rùn)性壁面的鋪展系數(shù)較為接近,均在t?大約為1.1時(shí)達(dá)到最大值,同時(shí)產(chǎn)生液殼破碎現(xiàn)象,可見(jiàn),壁面浸潤(rùn)性對(duì)中空液滴撞壁初始階段的鋪展過(guò)程和液殼破碎所需無(wú)量綱時(shí)間的影響較小.液殼破碎后,中心射流沿不同浸潤(rùn)性壁面均產(chǎn)生了收縮,因此鋪展系數(shù)隨無(wú)量綱時(shí)間均逐漸減小.但中空液滴碰撞浸潤(rùn)壁面(θe=30°)的鋪展系數(shù)大于憎浸壁面(θe=140°),說(shuō)明中心射流沿憎浸壁面的收縮速度較快.當(dāng)t?大約為6.5時(shí),不同碰撞速度下憎浸壁面的鋪展系數(shù)均減小到0,而浸潤(rùn)壁面的鋪展系數(shù)逐漸穩(wěn)定,說(shuō)明,不同碰撞速度下憎浸壁面中心射流均與壁面產(chǎn)生了剝離現(xiàn)象,而浸潤(rùn)壁面中心射流仍具有一定的鋪展特征.考慮到中心射流沿壁面收縮是慣性力、黏性力和表面張力共同作用的結(jié)果,射流沿壁面收縮過(guò)程表面能變化量ΔE可以表示成

    其中,S1和S2是收縮過(guò)程中的潤(rùn)濕面積.可以看出,中心射流沿壁面收縮表面能變化與接觸角有關(guān),當(dāng)射流潤(rùn)濕面積為一定值時(shí),接觸角越大,射流收縮需要克服的表面能越小,因此,中心射流沿憎浸壁面更容易產(chǎn)生收縮現(xiàn)象.此外,從圖7還可以看出,在中空液滴撞壁后期,當(dāng)碰撞速度較大時(shí)(U0=10m/s),浸潤(rùn)性壁面鋪展系數(shù)逐漸接近憎浸壁面,說(shuō)明隨著碰撞速度增大,壁面浸潤(rùn)性對(duì)中空液滴撞壁鋪展系數(shù)的影響逐漸減小.

    圖8 不同碰撞速度下中空液滴無(wú)量綱射流長(zhǎng)度的變化Fig.8.Dimensionless jet length for hollow droplet impact at di ff erent impacting velocity.

    圖8將不同碰撞速度下中空液滴撞擊不同浸潤(rùn)性壁面無(wú)量綱射流長(zhǎng)度L?表示成無(wú)量綱時(shí)間t?的函數(shù),其中,L?=L/D0,L為中心射流長(zhǎng)度.可以看出,不同工況下,無(wú)量綱射流長(zhǎng)度L?均隨t?近似線性遞增.碰撞速度相同時(shí),不同浸潤(rùn)性壁面的L?較為接近,說(shuō)明壁面浸潤(rùn)性對(duì)無(wú)量綱射流長(zhǎng)度的影響較小.當(dāng)中空液滴碰撞相同浸潤(rùn)性壁面時(shí),碰撞速度越大,無(wú)量綱射流長(zhǎng)度L?越大.由此可見(jiàn),無(wú)量綱射流長(zhǎng)度的碰撞速度效應(yīng)明顯.但隨著碰撞速度增大,無(wú)量綱射流長(zhǎng)度增幅逐漸減小,說(shuō)明中空液滴無(wú)量綱射流長(zhǎng)度與碰撞速度的相關(guān)性逐漸減小.

    進(jìn)一步地研究中空液滴撞壁傳熱特性,圖9將不同碰撞速度下中空液滴撞擊不同浸潤(rùn)性壁面平均熱流密度qw表示成無(wú)量綱時(shí)間t?的函數(shù).可以看出,不同碰撞速度下qw均隨t?先增加后減小.碰撞速度越大,qw越大,這是由于具有較大碰撞速度的液滴對(duì)壁面的沖擊擾動(dòng)較大,液滴與壁面間的傳熱能力也相應(yīng)地增大.但不同碰撞速度下壁面達(dá)到最大qw所需的無(wú)量綱時(shí)間較為接近,t?大約為0.16.可見(jiàn),中空液滴碰撞水平壁面qw與碰撞速度的相關(guān)性較大,但其達(dá)到最大qw所需無(wú)量綱時(shí)間與碰撞速度的相關(guān)性較小.結(jié)合圖7可知,中空液滴撞壁達(dá)到最大鋪展系數(shù)所需無(wú)量綱時(shí)間約為1.1,大于壁面達(dá)到最大qw所需無(wú)量綱時(shí)間,可見(jiàn),中空液滴撞壁鋪展特征明顯滯后于傳熱特征.此外,從圖9還可以看出,碰撞速度相同時(shí),不同浸潤(rùn)性壁面的qw較為接近,說(shuō)明壁面浸潤(rùn)性對(duì)中空液滴撞壁壁面平均熱流密度的影響較小.

    圖9 不同碰撞速度下壁面平均熱流密度的變化Fig.9.Average wall heat fl ux for hollow droplet impact at di ff erent impacting velocity.

    4.4 最大鋪展系數(shù)理論分析

    從能量守恒的角度對(duì)中空液滴撞壁過(guò)程進(jìn)行分析:液滴撞壁前后能量主要包括撞壁前液滴的動(dòng)能E1、表面能Es1以及撞壁后液膜的動(dòng)能E2、表面能Es2和黏性耗散能W,則中空液滴撞壁能量守恒方程可以表示成

    E2,Es2和W采用Gulyaev[12]所建立的關(guān)聯(lián)式,分別可以表示成

    式中,R(t)為鋪展半徑,h(t)為液膜厚度.結(jié)合(11)—(14)式并加入邊界條件[12]可得中空液滴撞壁最大鋪展系數(shù)表達(dá)式

    其中,β=(1?2δ?p)3,理論解析模型的適用范圍為12<Re<1400,230<We<15000.

    圖10 中空液滴fmax理論解析值與數(shù)值結(jié)果對(duì)比Fig.10.Comparison of fmaxfor hollow droplet impact between theoretical and numerical model.

    圖10將理論模型解析值與數(shù)值結(jié)果進(jìn)行對(duì)比,可以看出,不同Re下,數(shù)值計(jì)算fmax與理論解析值較為接近,均隨雷諾數(shù)的增加逐漸增大,但當(dāng)Re大約為500時(shí),fmax的增幅不再明顯,隨雷諾數(shù)的增加逐漸趨于穩(wěn)定.從圖10還可以看出,與數(shù)值計(jì)算值相比,理論解析fmax偏低.分析認(rèn)為,產(chǎn)生上述偏差主要有以下原因:1)理論模型在求解鋪展液膜厚度時(shí),對(duì)液滴形態(tài)進(jìn)行了假設(shè);2)理論模型在求解黏性耗散時(shí),沒(méi)有考慮熱量傳遞對(duì)液體黏度和流動(dòng)邊界層厚度的影響;3)現(xiàn)有相界面追蹤方法對(duì)接觸線移動(dòng)的精確描述具有局限性.

    5結(jié) 論

    采用CLSVOF法對(duì)中空液滴碰撞水平壁面進(jìn)行了數(shù)值分析,研究了中空液滴撞壁動(dòng)力學(xué)形態(tài)演變規(guī)律,揭示了中空液滴撞壁流動(dòng)傳熱機(jī)理和中心射流形成機(jī)制,并獲得了碰撞速度和壁面浸潤(rùn)性對(duì)中空液滴撞壁動(dòng)力學(xué)和傳熱特性的影響規(guī)律.得到主要結(jié)論如下:

    1)中空液滴撞壁后表現(xiàn)出明顯的鋪展和中心射流特征,并伴隨有中心射流收縮和液殼破碎等現(xiàn)象,實(shí)心液滴撞壁后主要呈現(xiàn)出鋪展和回縮特征;

    2)中空液滴內(nèi)部壓力梯度是液滴鋪展、中心射流產(chǎn)生和發(fā)展的主要原因,撞壁過(guò)程中中心射流表面溫度分布較為均勻,破碎液殼表面溫度分布存在較大波動(dòng);

    3)碰撞速度與中空液滴最大鋪展系數(shù)的相關(guān)性較小,但其對(duì)無(wú)量綱射流長(zhǎng)度和壁面平均熱流密度的影響較大;壁面浸潤(rùn)性與中空液滴撞壁后期鋪展系數(shù)的相關(guān)性較大,但其對(duì)無(wú)量綱射流長(zhǎng)度和壁面平均熱流密度的影響較小;中空液滴撞壁達(dá)到最大鋪展系數(shù)和最大熱流密度所需無(wú)量綱時(shí)間與碰撞速度和壁面浸潤(rùn)性均幾乎不相關(guān),且中空液滴撞壁鋪展特征明顯滯后于傳熱特征;

    4)中空液滴最大鋪展系數(shù)數(shù)值結(jié)果與理論解析值較為接近,均隨雷諾數(shù)的增加逐漸增大,但當(dāng)Re大約為500時(shí),fmax的增幅不再明顯,隨雷諾數(shù)的增加逐漸趨于穩(wěn)定.

    [1]Moreira A L N,Moita A S,Panao M R 2010 Prog.Energ.Combust.36 554

    [2]Li D S,Qiu X Q,Yu L,Xu J,Duan X L,Zheng Z W 2014 Ind.Heating 43 1(in Chinese)[李大樹(shù),仇性啟,于磊,許京,段小龍,鄭志偉2014工業(yè)加熱43 1]

    [3]Liang G T,Guo Y L,Shen S Q 2013 Acta Phys.Sin.62 024705(in Chinese)[梁剛濤,郭亞麗,沈勝?gòu)?qiáng)2013物理學(xué)報(bào)62 024705]

    [4]Rioboo R,Tropea C,Marengo M 2001 Atomization Spray 11 155

    [5]Bi F F,Guo Y L,Shen S Q,Chen J X,Li Y Q 2012 Acta Phys.Sin.61 184702(in Chinese)[畢菲菲,郭亞麗,沈勝?gòu)?qiáng),陳覺(jué)先,李熠橋2012物理學(xué)報(bào)61 184702]

    [6]Li D S,Qiu X Q,Cui Y J,Zheng Z W,Ma P Y,Qi F L 2014 T.Chin.Soc.Agric.Mach.45 25(in Chinese)[李大樹(shù),仇性啟,崔運(yùn)靜,鄭志偉,馬培勇,祁風(fēng)雷2014農(nóng)業(yè)機(jī)械學(xué)報(bào)45 25]

    [7]Kwon T J 2003 Ph.D.Dissertation(West Lafayette:Purdue University)

    [8]Yokoi K 2011 Soft Matter 7 5120

    [9]Solonenko O P,Gulyaev I P,Smirnov A V 2008 Tech.Phys.Lett.34 1050

    [10]Shukla R K,Kumar A 2015 J.Therm.Spray Techn.24 1368

    [11]Gulyaev I P,Solonenko O P,Gulyaev P Y,Smirnov A V 2009 Tech.Phys.Lett.35 885

    [12]Gulyaev I P,Solonenko O P 2013 Exp.Fluids 54 1

    [13]Solonenko O P,Smirnov A V,Gulyaev I P 2008 The 5th International Workshop on Complex Systems Sendai,Japan,September 26–-28,2008 982 561

    [14]Shinoda K,Murakami H 2010 J.Therm.Spray Techn.19 602

    [15]Kumar A,Gu S 2012 Int.J.Heat Fluid Fl.37 189

    [16]Kumar A,Gu S,Tabbara H,Kamnis S 2013 Surf.Coat.Techn.220 164

    [17]Tabbara H,Gu S 2012 Int.J.Heat Mass Trans.55 2081

    [18]Brackbill J U,Kothe D B,Zemach C 1992 J.Comput.Phys.100 335

    [19]Yokoi K 2013 J.Comput.Phys.232 252

    [20]Guo Y,Wei L,Liang G,Shen S 2014 Int.Commun.Heat Mass 53 26

    [21]Yang B H 2013 Ph.D.Dissertation(Chongqing:Chongqing University)(in Chinese)[楊寶海 2013博士學(xué)位論文(重慶:重慶大學(xué))]

    [22]Ubbink O,Issa R I 1999 J.Comput.Phys.153 26

    [23]Liang C,Wang H,Zhu X,Chen R,Ding Y D,Liao Q 2013 CIESC J.64 2745(in Chinese)[梁超,王宏,朱恂,陳蓉,丁玉棟,廖強(qiáng)2013化工學(xué)報(bào)64 2745]

    [24]Liu D W,Ning Z,Lü M,Yan K,Sun C H 2016 Chin.J.Comput.Mech.33 3027(in Chinese)[劉冬薇,寧智,呂明,閻凱,孫春華2016計(jì)算力學(xué)學(xué)報(bào)33 3027]

    [25]Song Y C,Ning Z,Sun C H,Lü M,Yan K,Fu J 2013 J.Combust.Sci.Technol.19 549(in Chinese)[宋云超,寧智,孫春華,呂明,閻凱,付娟 2013燃燒科學(xué)與技術(shù) 19 549]

    [26] ?ikalo ?,Marengo M,Tropea C,Ganic E N 2002 Exp.Therm.Fluid Sci.25 503

    PACS:47.55.db,47.55.nd,47.55.drDOI:10.7498/aps.66.014704

    *Project supported by the National Science and Technology Major Project of the Ministry of Science and Technology of China(Grant No.GZH201200602-01)and the Fundamental Research Funds for the Central Universities of Ministry of Education of China(Grant No.15CX06052A).

    ?Corresponding author.E-mail:lstax01@163.com

    ?Corresponding author.E-mail:apvshi@upc.edu.cn

    Numerical analysis of hollow droplet impact on a fl at surface?

    Zheng Zhi-Wei1)Li Da-Shu2)?Qiu Xing-Qi1)?Cui Yun-Jing3)

    1)(College of Chemical Engineering,China University of Petroleum,Qingdao 266580,China)2)(CNOOC Research Institute,Beijing 100028,China)3)(College of Machinery and Electrical Engineering,China University of Petroleum,Qingdao 266580,China)(Received 4 June 2016;revised manuscript received 26 September 2016)

    Many researches of a dense droplet impacting on a fl at surface have been reported in the literature.However,the mechanism of a hollow droplet impacting on a fl at surface has not yet been well addressed.A mathematical model is developed in the present research to resolve this impacting process.The model couples level set and volume of fl uid method,and considers heat transfer and contact resistance between the droplet and surface.The validation of the model is carried out by comparing simulation results with experiment data.Di ff erent impact behaviors are observed in the impacting processes of both the dense droplet and the hollow droplet on a fl at surface,obtained from the simulation result.The hydrodynamics and heat transfer behaviors of the hollow droplet impacting on a fl at surface and the formation of central jetting are also explored.The e ff ects of impact velocity and surface wettability on the impacting behavior of the hollow droplet are also analyzed.The results show that in the impacting process,the hollow droplet presents a spread and central jetting pattern,accompanying liquid shell contraction and breakup,while only spread and liquid shell contraction are observed in the dense droplet impacting process.It is also observed that the central jetting of the hollow droplet peels o ffthe surface in the fi nal impacting stage.The dimensionless spread factor for the hollow droplet is less than that of the dense droplet with the same initial kinetic energy in spread stage.The pressure gradient inside the hollow droplet is the main factor resulting in the spread and central jetting.The temperature distribution in the liquid shell and the surface is more uniform than in the central jetting,which is caused by the secondary breakup of the liquid shell.The spread factor of the hollow droplet remains unchanged as the impact velocity increases but is closely related to the surface wettability.The spread factor of the hydrophilic surface is larger than that of the hydrophobic surface.The e ff ects of the surface wettability on the spread factor gradually reduce with the increase of the impact velocity.The e ff ects of the impact velocity on the dimensionless jet length and the average wall heat fl ux are signi fi cant,while the surface wettability plays a negligible role in them.Improving the impact velocity increases the dimensionless length of the central jetting and the average wall heat fl ux,but this in fl uence diminishes under a high impact velocity condition.Neither the dimensionless time spans of reaching the maximum spread factor nor the maximum average wall heat fl ux for the hollow droplet is in fl uenced by the impact velocity and surface wettability and the development of the spread falls behind the heat transfer.Furthermore,the maximum spread factor increases with Reynolds number,and when Reynolds number is higher than 500,the increase in the maximum spread factor is no longer signi fi cant.

    hollow droplet,droplet impact,central jetting

    10.7498/aps.66.014704

    ?國(guó)家科技重大專項(xiàng)(批準(zhǔn)號(hào):GZH201200602-01)和中央高校基本科研業(yè)務(wù)費(fèi)專項(xiàng)資金(批準(zhǔn)號(hào):15CX06052A)資助的課題.

    ?通信作者.E-mail:lstax01@163.com

    ?通信作者.E-mail:apvshi@upc.edu.cn

    猜你喜歡
    實(shí)心中空浸潤(rùn)性
    浸潤(rùn)性乳腺癌超聲及造影表現(xiàn)與P63及Calponin的相關(guān)性
    一種實(shí)心輪胎檢測(cè)裝置及其系統(tǒng)
    GFRP管實(shí)心混凝土組合柱的研究現(xiàn)狀與發(fā)展
    江西建材(2018年4期)2018-04-10 12:36:52
    Q22、Q25 mmCr- Ni-Mo、Cr-Ni-W系列正七邊形中空釬鋼的研發(fā)
    某圓端型實(shí)心橋墩動(dòng)力彈塑性時(shí)程分析
    乳腺浸潤(rùn)性微乳頭狀癌的研究進(jìn)展
    乳腺浸潤(rùn)性導(dǎo)管癌組織β-catenin、cyclinD1、CDK4蛋白表達(dá)及臨床意義
    球磨機(jī)中空軸裂縫處理的新方法
    新疆鋼鐵(2015年1期)2015-11-07 03:11:03
    實(shí)心輪胎均勻硫化工藝的研究
    中空碳化硅微球的制備及其在催化NaBH4制氫中的應(yīng)用
    寂寞人妻少妇视频99o| 国产亚洲一区二区精品| 麻豆成人午夜福利视频| 欧美成人a在线观看| 性插视频无遮挡在线免费观看| 欧美激情在线99| 男女视频在线观看网站免费| 国产 亚洲一区二区三区 | 中文字幕亚洲精品专区| 2021少妇久久久久久久久久久| av一本久久久久| 久久久久久国产a免费观看| 亚洲精品国产成人久久av| 嫩草影院新地址| 亚洲国产欧美人成| 免费看美女性在线毛片视频| av在线天堂中文字幕| 丝袜美腿在线中文| 国精品久久久久久国模美| 一区二区三区免费毛片| 国产精品久久久久久精品电影小说 | 日日摸夜夜添夜夜爱| 两个人视频免费观看高清| 美女内射精品一级片tv| 亚洲四区av| 亚洲欧美一区二区三区国产| av在线播放精品| 99视频精品全部免费 在线| 亚洲国产高清在线一区二区三| 国产黄频视频在线观看| 最新中文字幕久久久久| 亚洲国产欧美在线一区| 久久久精品欧美日韩精品| 综合色av麻豆| 18禁在线无遮挡免费观看视频| 禁无遮挡网站| 亚洲最大成人手机在线| 欧美成人一区二区免费高清观看| 美女xxoo啪啪120秒动态图| 国产亚洲午夜精品一区二区久久 | 一级二级三级毛片免费看| 国产精品精品国产色婷婷| 丰满少妇做爰视频| 在线免费十八禁| 99久国产av精品| 国产色爽女视频免费观看| 亚洲,欧美,日韩| 免费观看av网站的网址| av播播在线观看一区| 中文欧美无线码| 日韩三级伦理在线观看| 国产91av在线免费观看| 日本免费a在线| 日韩强制内射视频| 亚洲内射少妇av| av一本久久久久| 欧美性猛交╳xxx乱大交人| 亚洲精品日韩av片在线观看| 久久99精品国语久久久| 久久这里只有精品中国| 亚洲内射少妇av| 免费观看无遮挡的男女| 亚洲高清免费不卡视频| 麻豆成人午夜福利视频| 久热久热在线精品观看| 晚上一个人看的免费电影| 国产在线一区二区三区精| 国产亚洲一区二区精品| 女人十人毛片免费观看3o分钟| 六月丁香七月| 亚洲国产精品成人综合色| 男插女下体视频免费在线播放| 精品人妻熟女av久视频| 韩国av在线不卡| 国产精品蜜桃在线观看| 国产大屁股一区二区在线视频| 国产免费又黄又爽又色| 日日摸夜夜添夜夜添av毛片| 18+在线观看网站| 三级经典国产精品| 肉色欧美久久久久久久蜜桃 | 国产高潮美女av| 搡女人真爽免费视频火全软件| 成人漫画全彩无遮挡| 亚洲四区av| 三级毛片av免费| av线在线观看网站| 日韩成人伦理影院| 最近中文字幕2019免费版| 久久99精品国语久久久| 成人国产麻豆网| 日韩强制内射视频| 五月天丁香电影| 亚洲美女搞黄在线观看| 免费少妇av软件| 久久久久久久久久人人人人人人| 婷婷色综合大香蕉| 日本av手机在线免费观看| 亚洲精品国产av成人精品| 99热全是精品| 欧美潮喷喷水| 插逼视频在线观看| 99久久精品热视频| 嫩草影院入口| 高清毛片免费看| 国产精品一区二区在线观看99 | 国产精品女同一区二区软件| 国产人妻一区二区三区在| 欧美xxxx性猛交bbbb| 午夜精品一区二区三区免费看| 欧美最新免费一区二区三区| 亚洲成人精品中文字幕电影| 激情 狠狠 欧美| 狂野欧美白嫩少妇大欣赏| 国产 一区 欧美 日韩| 欧美日本视频| 日韩人妻高清精品专区| 精品久久久久久久久av| 性色avwww在线观看| 老司机影院成人| 如何舔出高潮| 国产高清三级在线| 人妻制服诱惑在线中文字幕| 亚洲欧美日韩无卡精品| 久久综合国产亚洲精品| 欧美xxxx性猛交bbbb| 国产精品人妻久久久久久| 国产精品国产三级专区第一集| 日本猛色少妇xxxxx猛交久久| 日日干狠狠操夜夜爽| 好男人在线观看高清免费视频| 晚上一个人看的免费电影| 日韩在线高清观看一区二区三区| 人妻夜夜爽99麻豆av| 美女国产视频在线观看| 亚洲欧美中文字幕日韩二区| 久久久欧美国产精品| 国产一区二区三区av在线| 91精品伊人久久大香线蕉| 一级爰片在线观看| 黄片无遮挡物在线观看| 亚洲精品国产av成人精品| 日韩视频在线欧美| 婷婷色综合大香蕉| 国产亚洲av嫩草精品影院| 国产综合懂色| 日日撸夜夜添| 舔av片在线| 51国产日韩欧美| 精品久久久噜噜| 日本免费在线观看一区| 国产亚洲5aaaaa淫片| 亚洲综合色惰| 久久99热这里只频精品6学生| 日本欧美国产在线视频| 国产91av在线免费观看| 国产伦在线观看视频一区| 欧美精品一区二区大全| 久久久久久久久久黄片| 欧美成人一区二区免费高清观看| 又爽又黄无遮挡网站| 国产午夜精品一二区理论片| 免费观看在线日韩| 亚洲精品乱久久久久久| 一级二级三级毛片免费看| 国产精品美女特级片免费视频播放器| 国产黄片美女视频| 亚洲av男天堂| 少妇熟女欧美另类| 熟女电影av网| 一个人免费在线观看电影| 别揉我奶头 嗯啊视频| 国产精品久久久久久久电影| 午夜亚洲福利在线播放| 中文在线观看免费www的网站| 久久久亚洲精品成人影院| 最近中文字幕高清免费大全6| 中文字幕制服av| 午夜福利在线观看免费完整高清在| 久久久精品94久久精品| 国产亚洲av片在线观看秒播厂 | 欧美性猛交╳xxx乱大交人| 亚洲av在线观看美女高潮| 成人毛片60女人毛片免费| 欧美性感艳星| 夜夜爽夜夜爽视频| 99久久精品热视频| 非洲黑人性xxxx精品又粗又长| 精品酒店卫生间| 97超视频在线观看视频| 亚洲精品国产av蜜桃| 精品熟女少妇av免费看| 欧美zozozo另类| 亚洲在线观看片| 免费高清在线观看视频在线观看| 女人十人毛片免费观看3o分钟| 99热这里只有是精品50| 久久久久精品久久久久真实原创| 亚洲精品成人av观看孕妇| 免费大片黄手机在线观看| 亚洲国产欧美在线一区| 狠狠精品人妻久久久久久综合| 国产成人精品久久久久久| 久久久精品欧美日韩精品| 亚洲精品国产av蜜桃| 亚洲精品aⅴ在线观看| 国产综合精华液| 国产一区有黄有色的免费视频 | 国产精品久久久久久精品电影| 乱系列少妇在线播放| 十八禁网站网址无遮挡 | 国产在线男女| 亚洲自偷自拍三级| 久久精品国产亚洲av天美| av在线蜜桃| 免费观看性生交大片5| 高清在线视频一区二区三区| 只有这里有精品99| 精品一区二区三区视频在线| 夜夜爽夜夜爽视频| 午夜福利高清视频| 亚洲内射少妇av| 亚洲国产最新在线播放| 国产精品久久久久久精品电影小说 | 午夜免费观看性视频| 日本黄大片高清| 亚洲一级一片aⅴ在线观看| 能在线免费看毛片的网站| 国产在视频线在精品| 国产老妇女一区| 床上黄色一级片| 亚洲欧美一区二区三区国产| 在线播放无遮挡| 国产一区二区三区av在线| 国产黄片视频在线免费观看| 91精品国产九色| 少妇熟女aⅴ在线视频| 最新中文字幕久久久久| 性插视频无遮挡在线免费观看| 日韩欧美国产在线观看| 在线天堂最新版资源| 久久久久久久久大av| 国产老妇伦熟女老妇高清| 欧美不卡视频在线免费观看| 天美传媒精品一区二区| 97超碰精品成人国产| 国产 亚洲一区二区三区 | 少妇熟女欧美另类| 精品熟女少妇av免费看| 亚洲欧美日韩卡通动漫| 亚洲无线观看免费| 淫秽高清视频在线观看| 午夜免费观看性视频| 国产成人免费观看mmmm| 亚洲久久久久久中文字幕| 亚洲精品久久久久久婷婷小说| 欧美日韩国产mv在线观看视频 | 秋霞伦理黄片| 日韩欧美三级三区| 亚洲av电影不卡..在线观看| 亚洲自拍偷在线| 久久久久精品性色| 国产免费福利视频在线观看| 一区二区三区乱码不卡18| 成人高潮视频无遮挡免费网站| 久久久国产一区二区| 国产在线一区二区三区精| 久久99蜜桃精品久久| 亚洲av免费高清在线观看| 一个人看的www免费观看视频| 舔av片在线| 成人综合一区亚洲| 免费电影在线观看免费观看| 熟妇人妻不卡中文字幕| 黄片wwwwww| 久久久久精品性色| 亚洲欧美一区二区三区黑人 | 男女啪啪激烈高潮av片| 国产亚洲精品久久久com| 免费大片黄手机在线观看| 嫩草影院精品99| 国产 一区 欧美 日韩| 在线观看一区二区三区| 建设人人有责人人尽责人人享有的 | 97精品久久久久久久久久精品| 免费观看av网站的网址| 最后的刺客免费高清国语| 欧美人与善性xxx| 97在线视频观看| 成人美女网站在线观看视频| 最近最新中文字幕免费大全7| 精品99又大又爽又粗少妇毛片| 国产爱豆传媒在线观看| 午夜老司机福利剧场| 日韩,欧美,国产一区二区三区| 国产高清国产精品国产三级 | 女的被弄到高潮叫床怎么办| 国产三级在线视频| 国产精品嫩草影院av在线观看| 精品少妇黑人巨大在线播放| 深爱激情五月婷婷| 久久久欧美国产精品| 亚洲自偷自拍三级| 亚洲在线观看片| 亚洲av.av天堂| 欧美日韩亚洲高清精品| 亚洲熟女精品中文字幕| 国产中年淑女户外野战色| 国内精品美女久久久久久| 国产精品久久视频播放| 欧美激情久久久久久爽电影| 久久6这里有精品| 亚洲国产日韩欧美精品在线观看| 中国国产av一级| 天堂av国产一区二区熟女人妻| 在线免费观看的www视频| 蜜桃久久精品国产亚洲av| 91av网一区二区| 人妻系列 视频| 久久精品人妻少妇| 非洲黑人性xxxx精品又粗又长| 中文天堂在线官网| 久久久久网色| 边亲边吃奶的免费视频| 亚洲成人久久爱视频| 亚洲av日韩在线播放| 欧美日韩国产mv在线观看视频 | 欧美激情国产日韩精品一区| 搡老乐熟女国产| 久久精品熟女亚洲av麻豆精品 | 国产久久久一区二区三区| 欧美日韩精品成人综合77777| 干丝袜人妻中文字幕| 国产av在哪里看| 国产成人aa在线观看| 两个人的视频大全免费| 精品一区二区三区视频在线| 午夜福利高清视频| 在线免费观看不下载黄p国产| 永久免费av网站大全| 一本一本综合久久| 日本免费在线观看一区| 日韩亚洲欧美综合| 成人毛片60女人毛片免费| 国产精品嫩草影院av在线观看| 欧美日韩一区二区视频在线观看视频在线 | 2018国产大陆天天弄谢| 成人欧美大片| 麻豆国产97在线/欧美| 精品国内亚洲2022精品成人| 久久精品综合一区二区三区| 欧美变态另类bdsm刘玥| 国产成人免费观看mmmm| 观看免费一级毛片| 国产欧美另类精品又又久久亚洲欧美| 欧美97在线视频| 亚洲av免费在线观看| 十八禁国产超污无遮挡网站| 国产成年人精品一区二区| 亚洲成人一二三区av| 两个人的视频大全免费| av在线播放精品| 亚洲精品一区蜜桃| 日韩人妻高清精品专区| 两个人视频免费观看高清| 成年av动漫网址| 成人毛片a级毛片在线播放| 搡老妇女老女人老熟妇| 青春草视频在线免费观看| 日本wwww免费看| 99热这里只有是精品在线观看| 免费av不卡在线播放| 国产亚洲91精品色在线| 99久久人妻综合| 亚洲av中文av极速乱| 插阴视频在线观看视频| 亚洲国产日韩欧美精品在线观看| 97超视频在线观看视频| 久久草成人影院| 91午夜精品亚洲一区二区三区| 欧美成人a在线观看| 久久久欧美国产精品| 我的女老师完整版在线观看| 国产一级毛片在线| 中文天堂在线官网| 亚洲婷婷狠狠爱综合网| 中文资源天堂在线| 亚洲综合色惰| 久久精品久久久久久久性| 亚洲第一区二区三区不卡| 嫩草影院新地址| 秋霞伦理黄片| 97超视频在线观看视频| 极品教师在线视频| 91久久精品电影网| 亚洲欧美清纯卡通| 色综合色国产| 听说在线观看完整版免费高清| 人妻一区二区av| 一夜夜www| 黄色一级大片看看| 日韩欧美精品免费久久| 两个人的视频大全免费| 久久精品国产亚洲网站| 深夜a级毛片| 狠狠精品人妻久久久久久综合| 婷婷六月久久综合丁香| 亚洲在线观看片| 精品久久国产蜜桃| 麻豆成人av视频| 欧美zozozo另类| 三级男女做爰猛烈吃奶摸视频| 视频中文字幕在线观看| 黄色欧美视频在线观看| videos熟女内射| 久久久久久九九精品二区国产| 午夜福利在线观看吧| 一个人看视频在线观看www免费| 天堂√8在线中文| 日日啪夜夜爽| 少妇人妻一区二区三区视频| 国产乱来视频区| 22中文网久久字幕| 一区二区三区乱码不卡18| 白带黄色成豆腐渣| 偷拍熟女少妇极品色| 亚洲在线观看片| 国产一区二区三区av在线| 丝袜美腿在线中文| 久久久久免费精品人妻一区二区| 日韩大片免费观看网站| 午夜爱爱视频在线播放| 亚洲av成人精品一区久久| 免费观看a级毛片全部| 少妇人妻精品综合一区二区| 熟女电影av网| 啦啦啦中文免费视频观看日本| 午夜精品一区二区三区免费看| 国产永久视频网站| 嫩草影院入口| 日韩国内少妇激情av| 国产成人91sexporn| 啦啦啦啦在线视频资源| 久久久久久九九精品二区国产| 69av精品久久久久久| 亚洲,欧美,日韩| 97在线视频观看| 最近最新中文字幕大全电影3| 少妇丰满av| 2021少妇久久久久久久久久久| 久久久久久久久久成人| 综合色av麻豆| 精品亚洲乱码少妇综合久久| 亚洲精品久久久久久婷婷小说| 午夜福利视频1000在线观看| 白带黄色成豆腐渣| 91在线精品国自产拍蜜月| a级一级毛片免费在线观看| 日韩三级伦理在线观看| 国产视频内射| 久久国产乱子免费精品| 精品一区二区三区视频在线| 亚洲欧美精品自产自拍| 欧美zozozo另类| 精品酒店卫生间| 69人妻影院| 亚洲在久久综合| 深爱激情五月婷婷| 熟女电影av网| 免费av毛片视频| 亚洲av电影不卡..在线观看| 久久热精品热| 婷婷色综合www| 伦理电影大哥的女人| 校园人妻丝袜中文字幕| 午夜免费男女啪啪视频观看| 少妇的逼好多水| 日日啪夜夜爽| 日韩精品有码人妻一区| 欧美 日韩 精品 国产| 亚洲精品日韩在线中文字幕| 一本一本综合久久| 亚洲国产高清在线一区二区三| 一本一本综合久久| 亚洲成人精品中文字幕电影| 国产综合懂色| 中国国产av一级| 人人妻人人澡人人爽人人夜夜 | 欧美极品一区二区三区四区| 亚洲av福利一区| 亚洲精品一二三| 卡戴珊不雅视频在线播放| 寂寞人妻少妇视频99o| 天堂俺去俺来也www色官网 | 欧美日韩亚洲高清精品| 777米奇影视久久| 国产成人一区二区在线| 国产在视频线精品| 国产极品天堂在线| 久久精品国产亚洲av天美| 九色成人免费人妻av| 国产探花极品一区二区| 国内精品宾馆在线| 深夜a级毛片| 观看免费一级毛片| 免费看美女性在线毛片视频| 一级毛片黄色毛片免费观看视频| 久久精品国产亚洲网站| 99热网站在线观看| 免费看不卡的av| 成人漫画全彩无遮挡| 国产伦一二天堂av在线观看| 午夜亚洲福利在线播放| 国内精品宾馆在线| 成人亚洲精品一区在线观看 | 夜夜爽夜夜爽视频| 国产精品人妻久久久影院| 国产探花在线观看一区二区| 日韩强制内射视频| 特级一级黄色大片| 91精品伊人久久大香线蕉| 久久久久久久久久人人人人人人| 最后的刺客免费高清国语| 国产午夜福利久久久久久| 在现免费观看毛片| 在线观看av片永久免费下载| 精品久久久精品久久久| 看十八女毛片水多多多| 亚洲高清免费不卡视频| 亚洲国产高清在线一区二区三| 日本欧美国产在线视频| 午夜福利网站1000一区二区三区| 中文乱码字字幕精品一区二区三区 | 18禁在线播放成人免费| 国产精品久久久久久久电影| 免费看av在线观看网站| 亚洲熟妇中文字幕五十中出| 九草在线视频观看| 乱系列少妇在线播放| 搡女人真爽免费视频火全软件| 成人午夜精彩视频在线观看| 国产精品久久视频播放| 国产亚洲5aaaaa淫片| 日韩欧美国产在线观看| 免费大片18禁| 亚洲欧美精品自产自拍| 国产精品蜜桃在线观看| 高清av免费在线| 成年版毛片免费区| 国产美女午夜福利| 国产成人91sexporn| 97超碰精品成人国产| 精品久久久精品久久久| 亚洲精品,欧美精品| 一级爰片在线观看| 欧美zozozo另类| 日韩av不卡免费在线播放| 精品亚洲乱码少妇综合久久| 深爱激情五月婷婷| 日本黄大片高清| 亚洲精品国产av蜜桃| 色综合站精品国产| 欧美最新免费一区二区三区| 国产精品一二三区在线看| 国产一区亚洲一区在线观看| 岛国毛片在线播放| 亚洲精品第二区| 美女国产视频在线观看| av在线亚洲专区| 男人和女人高潮做爰伦理| 少妇人妻一区二区三区视频| 可以在线观看毛片的网站| 18禁在线无遮挡免费观看视频| 亚洲真实伦在线观看| 国产 一区 欧美 日韩| 亚洲精品久久午夜乱码| 两个人的视频大全免费| 久久精品国产亚洲av天美| av在线天堂中文字幕| 日本黄大片高清| 九九爱精品视频在线观看| 国产精品爽爽va在线观看网站| 亚洲国产高清在线一区二区三| 99久国产av精品| 色综合站精品国产| 日本免费在线观看一区| 亚洲国产成人一精品久久久| 国产一区二区亚洲精品在线观看| 大陆偷拍与自拍| 久久草成人影院| 一级毛片电影观看| 肉色欧美久久久久久久蜜桃 | 亚洲av二区三区四区| 人妻系列 视频| 国产精品一区二区三区四区久久| 久久久久久久久久久免费av| 秋霞在线观看毛片| 精品午夜福利在线看| 三级经典国产精品| 美女高潮的动态| 国产成人福利小说| 欧美区成人在线视频| 极品教师在线视频| 淫秽高清视频在线观看| 亚洲欧美精品自产自拍| 久久精品国产亚洲av天美| 综合色丁香网| 久久久亚洲精品成人影院| 国产色爽女视频免费观看| 色播亚洲综合网| 亚洲乱码一区二区免费版| 日韩视频在线欧美| 欧美成人一区二区免费高清观看| 亚洲精品日韩av片在线观看| 小蜜桃在线观看免费完整版高清|