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

    基于勢(shì)流理論的回轉(zhuǎn)體并聯(lián)入水雙空泡演化動(dòng)力學(xué)研究?

    2018-12-18 05:58:18宋武超魏英杰路麗睿王聰盧佳興
    物理學(xué)報(bào) 2018年22期
    關(guān)鍵詞:空泡壁面并聯(lián)

    宋武超 魏英杰 路麗睿 王聰 盧佳興

    (哈爾濱工業(yè)大學(xué)航天學(xué)院,哈爾濱 150001)

    (2018年7月17日收到;2018年9月11日收到修改稿)

    1 引 言

    運(yùn)動(dòng)體以一定的速度撞擊自由液面,將自身動(dòng)能傳遞給了附近水域,引起流體質(zhì)點(diǎn)向運(yùn)動(dòng)體外側(cè)運(yùn)動(dòng),致使空泡敞開及自由液面噴濺,形成入水空泡.作為典型跨界質(zhì)過程的入水過程,涉及運(yùn)動(dòng)體、空氣和水三者的耦合作用,具有高度非線性和非定常性,運(yùn)動(dòng)體的運(yùn)動(dòng)特性和流場特性都極為復(fù)雜.

    早在19世紀(jì)末,Worthinghton和Cole[1,2]利用閃光攝影技術(shù)對(duì)球體等垂直入水時(shí)產(chǎn)生的入水空泡進(jìn)行了大量的試驗(yàn)研究,并對(duì)入水空泡流動(dòng)現(xiàn)象做了開創(chuàng)性的描述.二戰(zhàn)以后,在美國海軍支持下,學(xué)者們針對(duì)入水過程開展了系統(tǒng)性的試驗(yàn)研究工作,其中代表性的有Waugh[3],Birkho ff和Caywood[4],Glibarg和Anderson[5]以及May和Woodhull[6]等.21世紀(jì)以來,隨著試驗(yàn)設(shè)備、測量技術(shù)的革新,從運(yùn)動(dòng)體入水過程的試驗(yàn)中獲得的信息更加全面,為入水空泡流體特性的研究提供了更加有效的方法和途徑[7].Holfeld等[8]和Tabuteau等[9]分別在微重力、不同黏度液體等特殊條件下開展了入水實(shí)驗(yàn)研究.對(duì)于入水空泡形態(tài)預(yù)測的理論方法,Logvinovich[10]首次提出了基于能量守恒定律的空泡獨(dú)立膨脹原理,為之后空泡發(fā)展理論研究提供了重要基礎(chǔ).Lee等[11,12]基于能量守恒定律,提出了描述高速入水空泡形態(tài)發(fā)展演化的理論方法,對(duì)入水空泡擴(kuò)展、收縮和閉合過程進(jìn)行了較為精確的描述.Truscott[13]將Logvinovich提出的獨(dú)立膨脹原理進(jìn)行了擴(kuò)展,準(zhǔn)確地預(yù)測了運(yùn)動(dòng)體傾斜入水過程中入水空泡形體的演化規(guī)律.Yan等[14]基于勢(shì)流理論建立了低速入水空泡動(dòng)力學(xué)模型.Gordillo等[15]采用勢(shì)流理論推導(dǎo)了高雷諾數(shù)入水條件下空泡閉合點(diǎn)流動(dòng)特征.Aristo ff和Bush[16]同樣基于勢(shì)流理論得到了入水空泡形態(tài)的理論解.國內(nèi)方面,葉取源[17]基于歐拉-拉格朗日混合邊界元法,求解獲得了圓錐垂直入水過程中入水空泡形成、發(fā)展以及閉合的變化規(guī)律.路中磊[18]基于Yan等[14]的分析方法,建立了適用于開放空腔入水空泡波動(dòng)預(yù)測的理論模型.張阿漫等[19]基于勢(shì)流理論分析了氣槍氣泡遠(yuǎn)場壓力子波特性.

    在現(xiàn)有關(guān)于入水空泡動(dòng)力的研究中,多以單個(gè)運(yùn)動(dòng)體入水所形成入水空泡的動(dòng)力學(xué)特性作為研究重點(diǎn),主要考慮了空泡的演化過程,并對(duì)其流動(dòng)方式進(jìn)行了簡化,形成了較為完善的單獨(dú)低速軸對(duì)稱入水空泡動(dòng)力學(xué)計(jì)算模型.然而,對(duì)于雙運(yùn)動(dòng)體入水過程所形成的入水空泡的動(dòng)力特性,尚無人涉及.雙回轉(zhuǎn)體入水過程中所形成的入水空泡,呈現(xiàn)出了較強(qiáng)的非對(duì)稱性;兩回轉(zhuǎn)體軸線內(nèi)側(cè)區(qū)域空泡演化過程中,流體形成了相對(duì)流動(dòng),從而使得流動(dòng)不再具有軸對(duì)稱特性,而發(fā)展成復(fù)雜的三維流動(dòng),且具有強(qiáng)非線性,使得理論模型的建立難度更大.

    本文首先簡要介紹單獨(dú)入水空泡與并聯(lián)入水空泡各自的演化特性;而后基于勢(shì)流理論,在現(xiàn)有的單獨(dú)入水二維軸對(duì)稱計(jì)算模型的基礎(chǔ)上,將并聯(lián)入水過程中空泡之間的擾動(dòng)簡化為擾動(dòng)速度勢(shì),給出擾動(dòng)最強(qiáng)區(qū)域的空泡演化擾動(dòng)模型;基于非線性假設(shè),引入三維影響函數(shù),通過影響函數(shù)建立空泡自由擴(kuò)張區(qū)與擾動(dòng)擴(kuò)張區(qū)的過渡關(guān)系,得到三維非對(duì)稱空泡演化計(jì)算模型,并基于該模型分析同步并聯(lián)入水過程空泡形態(tài)的三維演化特性.

    2 單獨(dú)入水空泡與并聯(lián)入空泡演化特性分析

    2.1 實(shí)驗(yàn)系統(tǒng)與模型參數(shù)

    試驗(yàn)系統(tǒng)由試驗(yàn)水槽、光學(xué)測試系統(tǒng)和釋放機(jī)構(gòu)三部分組成,如圖1所示.試驗(yàn)水槽采用尺寸為300 mm×300 mm×900 mm的小型水槽,水槽四壁為鋼化玻璃,壁厚10 mm,槽底墊有緩沖橡膠墊,防止射彈入水對(duì)缸體造成破壞.試驗(yàn)光源采用4盞發(fā)光二極管(LED)作為背景光,并用柔光屏對(duì)背景光進(jìn)行柔化,使背景光均勻撒布,以達(dá)到良好的拍照效果;采用2盞LED燈作為頂光,以增強(qiáng)空泡內(nèi)部光照條件,進(jìn)而增強(qiáng)圓柱體邊界清晰度,便于圖片數(shù)據(jù)的后期處理.圓柱體運(yùn)動(dòng)過程通過Photron FASTCAM SA-X型高速攝像機(jī)采集,并通過自編程序?qū)\(yùn)動(dòng)軌跡及運(yùn)動(dòng)姿態(tài)進(jìn)行提取,以獲得圓柱體運(yùn)動(dòng)特性.相機(jī)采集幀率為2000 fps,相機(jī)上安裝廣角鏡頭,以擴(kuò)大視野,捕捉更多運(yùn)動(dòng)信息.

    釋放機(jī)構(gòu)由推拉式電磁鐵和豎直導(dǎo)軌組成,如圖2所示.其中豎直導(dǎo)軌起導(dǎo)向作用,以確保試驗(yàn)?zāi)P腿胨畷r(shí)保持豎直姿態(tài);推拉式電磁鐵通過觸發(fā)開關(guān)控制,共用同一個(gè)電信號(hào),以確保觸發(fā)的同步性.當(dāng)電源接通時(shí),電磁鐵推拉桿同時(shí)收縮,釋放模型.為防止電磁鐵收縮引起的支架振動(dòng)對(duì)釋放產(chǎn)生影響,導(dǎo)軌和電磁鐵分別安裝在不同的支架上,導(dǎo)軌支架和電磁鐵支架相互獨(dú)立.

    試驗(yàn)?zāi)P腿鐖D3所示.模型直徑D=10 mm,長L=60 mm,模型采用6063鋁合金加工,密度為2.7 g/cm3.

    圖1 試驗(yàn)系統(tǒng)示意圖Fig.1.Experiment system of water-entry.

    圖2 釋放機(jī)構(gòu)示意圖Fig.2.The mechanism to release the model.

    圖3 模型尺寸Fig.3.Sketch of water-entry model.

    2.2 實(shí)驗(yàn)結(jié)果分析

    圖4(a)給出了回轉(zhuǎn)體單獨(dú)入水空泡形態(tài)時(shí)歷演化過程.回轉(zhuǎn)體在入水過程中經(jīng)歷了入水撞擊(t=0 ms)、開空泡(t=10—30 ms)、空泡閉合(t=10—50 ms)、空泡隨動(dòng)及空泡潰滅(t=60—80 ms)等五個(gè)階段.在入水撞擊階段,回轉(zhuǎn)體通過入水沖擊將其自身動(dòng)能傳遞給附近水流體,水流體形成向外排開運(yùn)動(dòng),入水空泡形成;在開空泡階段,回轉(zhuǎn)體通過阻力作用不斷將動(dòng)能傳遞給附近水流體,致使回轉(zhuǎn)體底面邊緣處形成持續(xù)的流動(dòng)分離,在宏觀上表現(xiàn)為空泡敞開(t=10—30 ms);在空泡閉合階段,氣水界面處流體質(zhì)點(diǎn)在運(yùn)動(dòng)過程中受到壓力與表面張力的作用,做減速運(yùn)動(dòng),當(dāng)速度減小為0時(shí)開始反向加速,使得空泡形成“頸縮”現(xiàn)象,空泡發(fā)生閉合(t=40—50 ms);空泡完成閉合后,回轉(zhuǎn)體底面仍不斷形成流動(dòng)分離,分離所形成的空腔需要空氣進(jìn)行填充,使得封閉在空泡內(nèi)部的氣體隨著分離所形成的空腔位置的變化而流動(dòng),空泡形成隨動(dòng)(t=60—80 ms);在隨動(dòng)空泡形成后,閉合點(diǎn)處空泡內(nèi)外存在壓力梯度,使得水流體沿回轉(zhuǎn)體壁面向泡內(nèi)流動(dòng),形成回射流,空泡界面受回射流撞擊而逐漸云化潰滅.

    圖4(b)給出了回轉(zhuǎn)體同步并聯(lián)入水空泡形態(tài).對(duì)比單獨(dú)入水空泡形態(tài)與并聯(lián)入水空泡形態(tài)可以看出,單個(gè)回轉(zhuǎn)體擾動(dòng)下,流場流動(dòng)具有軸對(duì)稱性,因而空泡形態(tài)表現(xiàn)出良好的軸對(duì)稱特征;在同步并聯(lián)入水過程中,兩空泡在回轉(zhuǎn)體軸線內(nèi)側(cè)區(qū)域形成擾動(dòng),單空泡形態(tài)演化呈現(xiàn)出明顯的不對(duì)稱特征.在開空泡階段,空泡向兩回轉(zhuǎn)體軸線中間區(qū)域的擴(kuò)張受到了限制,以致兩回轉(zhuǎn)體外側(cè)空泡半徑明顯大于內(nèi)側(cè),中間區(qū)域形成的噴濺高度亦較外側(cè)高;在空泡閉合階段,兩回轉(zhuǎn)體軸線外側(cè)區(qū)域空泡界面率先收縮,空泡閉合點(diǎn)偏向兩回轉(zhuǎn)體軸線內(nèi)側(cè),且閉合時(shí)間變晚.在空泡隨動(dòng)階段,隨動(dòng)空泡亦呈現(xiàn)出不對(duì)稱特征.總之,在整個(gè)入水過程中,單獨(dú)入水空泡始終保持很強(qiáng)的軸對(duì)稱特性,而并聯(lián)入水空泡則始終呈現(xiàn)較強(qiáng)的鏡面對(duì)稱特征.

    圖4 單獨(dú)入水與并聯(lián)入水空泡演化對(duì)比(V0=1.9 m/s) (a)單獨(dú)入水空泡;(b)同步并聯(lián)入水空泡Fig.4.Comparison of parallel water-entry cavities and single water-entry cavities:(a)Single water-entry cavities;(b)parallel water-entry cavities.

    結(jié)合圖5對(duì)空泡形態(tài)不對(duì)稱原因解釋如下:在開空泡階段,外側(cè)空泡自由擴(kuò)張,而內(nèi)側(cè)空泡的發(fā)展受到了相鄰圓柱體的影響.空泡在兩圓柱體中間區(qū)域擴(kuò)張時(shí),空泡邊界附近流體形成了對(duì)流,流體質(zhì)點(diǎn)間動(dòng)能傳遞受到限制而轉(zhuǎn)變成流場勢(shì)能,致使中間區(qū)域噴濺高度增加而空泡半徑減小,空泡形態(tài)呈現(xiàn)出不對(duì)稱特性.同時(shí),由于中間區(qū)域較少的流體獲得了來自兩個(gè)圓柱體的動(dòng)能,流動(dòng)速度較外側(cè)區(qū)域大,致使該區(qū)域流場壓力低于外側(cè)空泡界面附近的流場壓力,導(dǎo)致在空泡閉合階段,外側(cè)空泡界面在較大壓力作用下率先收縮,且收縮速率大于內(nèi)側(cè)空泡界面,最終使得閉合點(diǎn)偏向兩圓柱體軸線內(nèi)側(cè).

    圖5 并聯(lián)入水空泡流動(dòng)示意圖Fig.5.The flow schematic of water entry cavity.

    3 單獨(dú)入水過程動(dòng)力學(xué)方程

    3.1 單獨(dú)入水過程初邊值分析

    分析可知,在單獨(dú)入水過程中入水空泡形態(tài)始終保持軸對(duì)稱特性,因此空泡演化動(dòng)力學(xué)模型建立在如圖6所示的柱坐標(biāo)系下.柱坐標(biāo)原點(diǎn)取在回轉(zhuǎn)體軸線與自由液面的交點(diǎn)處,r軸取回轉(zhuǎn)體徑向,z軸取回轉(zhuǎn)體軸線,其正方向指向回轉(zhuǎn)體運(yùn)動(dòng)方向.

    圖6 回轉(zhuǎn)體單獨(dú)入水空泡示意圖Fig.6.Schematic of single water-entry cavity.

    對(duì)于單獨(dú)回轉(zhuǎn)體入水過程,可假設(shè)流體處于無旋狀態(tài).因此,流場的速度可以用速度勢(shì)函數(shù)?(r,z,t)進(jìn)行描述,速度勢(shì)?(r,z,t)在流域內(nèi)滿足Laplace方程:

    此時(shí),流域內(nèi)的速度分布可以用下式表示:

    在自由液面處SF(r,t)處,速度勢(shì)滿足歐拉方程:

    其中ps為自由液面處的壓力,ρ為液體密度.在空泡閉合前,該壓力恒等于大氣壓.在回轉(zhuǎn)體表面,速度勢(shì)滿足物面條件:

    式中,v為回轉(zhuǎn)體運(yùn)動(dòng)速度,n=(nz,nr).速度勢(shì)在遠(yuǎn)場滿足遠(yuǎn)場邊界條件:

    在初始時(shí)刻,運(yùn)動(dòng)體僅僅觸碰到了自由液面,此時(shí),自由液面尚處于靜止?fàn)顟B(tài),初始速度勢(shì)滿足如下初始條件:

    3.2 回轉(zhuǎn)體動(dòng)力學(xué)方程及空泡演化方程

    在回轉(zhuǎn)體過程中,其動(dòng)力學(xué)方程可表達(dá)如下:

    式中M為回轉(zhuǎn)體質(zhì)量,F(t)為回轉(zhuǎn)體受到的流體動(dòng)力.回轉(zhuǎn)體在入水過程中,僅頭部沾濕,因此流體動(dòng)力可視為回轉(zhuǎn)體頭部所受到的流體壓力.該流體壓力由動(dòng)壓和靜壓兩部分組成,其表達(dá)式為

    式中,pw為回轉(zhuǎn)體頭部沾濕面上受到的壓力,其滿足伯努利方程:

    將(8)和(9)式代入(7)式得到回轉(zhuǎn)體單獨(dú)入水過程的動(dòng)力學(xué)方程如下:

    式中,Ab為回轉(zhuǎn)體底部面積,Ab=πR2,R為底部圓形面半徑;lc為空泡長度,其滿足

    對(duì)空泡進(jìn)行細(xì)長體假設(shè)[20],即dr/lc=ε?o(1),其中dr為空泡直徑,此時(shí),在空泡壁面附近的流動(dòng)可以假設(shè)為二維流動(dòng):

    此時(shí),(1)式可簡化為

    對(duì)于空泡發(fā)展過程,在入水空泡初生時(shí)刻,其半徑等于回轉(zhuǎn)體半徑.因此入水空泡半徑演化的隱式方程:

    其中,rc為空泡直徑.

    3.3 速度勢(shì)分析

    基于建立的入水空泡半徑演化方程(13)可知,若給出速度勢(shì)?(r,z,t),則可得到入水空泡半徑的演化規(guī)律,因此本節(jié)基于漸進(jìn)理論對(duì)速度勢(shì)函數(shù)?(r,z,t)進(jìn)行分析.

    為了描述空泡界面的運(yùn)動(dòng),入水過程的速度勢(shì)可以看作由一個(gè)位于回轉(zhuǎn)體底面中心處的點(diǎn)源?b和一個(gè)位于空泡軸線處的線源?c組成:

    空泡形成擴(kuò)張的速度主要來源于回轉(zhuǎn)體底部的空泡分離點(diǎn).因此,空泡的擴(kuò)張運(yùn)動(dòng)主要受三維點(diǎn)源?b的作用[14]:

    上式第二項(xiàng)對(duì)空泡演化的作用為抑制作用,由于對(duì)空泡進(jìn)行了細(xì)長假設(shè)[20],則(16)式的第二項(xiàng)為一小量,因此可忽略.根據(jù)文獻(xiàn)[14],點(diǎn)源強(qiáng)度由下式給出:

    式中,C1為修正系數(shù),文獻(xiàn)中給出的C1=2.78,本文經(jīng)過與試驗(yàn)結(jié)果對(duì)比,將其修正為2.76,得出點(diǎn)源速度勢(shì)為

    空泡的收縮主要受空泡軸線上分布的線源影響,其在某一空泡界面上具有二維作用效果,在整個(gè)空泡域,又具有三維作用效果.因此,空泡收縮的線源速度勢(shì)可表示為

    文獻(xiàn)[14]中給出了線源的強(qiáng)度表達(dá)式,由此可以確定總速度勢(shì):

    3.4 計(jì)算模型求解

    將速度勢(shì)代入(10)式得到回轉(zhuǎn)體單獨(dú)入水過程的完整動(dòng)力學(xué)方程,并采用經(jīng)典四階龍格-庫塔法對(duì)其進(jìn)行數(shù)值求解.在求解過程中,時(shí)間步長取?t=0.0001 s,在z方向(入水方向,如圖2所示)的離散取?z=0.001 m.圖7給出了求解所得z方向位移隨時(shí)間的變化與試驗(yàn)結(jié)果的對(duì)比,從圖中可以看出,計(jì)算結(jié)果與試驗(yàn)結(jié)果符合良好.

    圖7 位移結(jié)果對(duì)比(V0=1.9 m/s)Fig.7.Comparison of theory model calculative results and experimental results.

    圖8給出了入水空泡形態(tài)計(jì)算結(jié)果與試驗(yàn)結(jié)果的對(duì)比.從圖中可以看出,整體上看計(jì)算結(jié)果與試驗(yàn)結(jié)果符合度較好,證明了該理論模型計(jì)算結(jié)果的正確性.在靠近自由液面處,計(jì)算模型結(jié)果較試驗(yàn)結(jié)果差別較大,這是由于計(jì)算模型中忽略了自由液面處的非線性噴濺流動(dòng),但并不影響入水空泡整體形態(tài)預(yù)測.

    圖8 入水空泡形態(tài)計(jì)算結(jié)果與試驗(yàn)結(jié)果對(duì)比(V0=1.9 m/s)Fig.8.Comparison of water-entry cavities of calculative and experimental results(V0=1.9 m/s).

    4 同步并聯(lián)入水過程雙空泡演化擾動(dòng)模型

    4.1 動(dòng)力學(xué)建模簡化方法

    在回轉(zhuǎn)體同步并聯(lián)入水過程中,回轉(zhuǎn)體空泡外側(cè)的發(fā)展呈自由擴(kuò)張狀態(tài),內(nèi)側(cè)發(fā)展受到來自相鄰空泡的干擾,而呈現(xiàn)出受迫演化狀態(tài).由第2節(jié)對(duì)并聯(lián)入水空泡演化試驗(yàn)結(jié)果的分析,發(fā)現(xiàn)空泡在空間上呈對(duì)稱特性,因此可將內(nèi)側(cè)區(qū)域空泡演化對(duì)稱面作為擾動(dòng)勢(shì)對(duì)內(nèi)側(cè)空泡的發(fā)展演化進(jìn)行分析.該擾動(dòng)勢(shì)的作用范圍為?,如圖9所示.

    圖9 對(duì)稱面擾動(dòng)勢(shì)作用示意圖Fig.9.Schematic of symmetry plane disturbance potential.

    引入對(duì)稱面處的擾動(dòng)勢(shì),將入水空泡的演化通過二維軸對(duì)稱的方法進(jìn)行分析,其在三維空間中相當(dāng)于引入了一個(gè)柱面約束,空泡在該柱面約束下進(jìn)行軸對(duì)稱受迫發(fā)展,如圖10所示,該約束柱面可稱為有勢(shì)壁面.

    圖10 二維軸對(duì)稱模型有勢(shì)壁面約束本質(zhì)示意圖Fig.10.Schematic of potential wall constraint of the planar axisymmetrical model.

    4.2 有勢(shì)壁面速度勢(shì)分析

    由第2節(jié)體并聯(lián)入水空泡演化特性分析可知,回轉(zhuǎn)體入水過程約束勢(shì)?s可認(rèn)為由兩部分組成,一部分是對(duì)擴(kuò)張的約束,另一部分是對(duì)收縮的約束.結(jié)合3.3節(jié)關(guān)于速度勢(shì)分析可知,空泡的擴(kuò)張主要受位于回轉(zhuǎn)體底部中心處的點(diǎn)源速度勢(shì)?b的影響,而空泡的收縮主要受位于空泡軸線處的線源速度勢(shì)?c作用.因此,對(duì)于空泡擴(kuò)張與收縮的抑制,我們同樣引入一個(gè)隨回轉(zhuǎn)體底部中心運(yùn)動(dòng)的點(diǎn)匯?b,s和位于空泡軸線處的線源?c,s.有勢(shì)壁面的速度勢(shì)則可以表示如下:

    對(duì)于?b,s的分析,根據(jù)勢(shì)流疊加原理,通過分析?b,s與?b的合速度勢(shì)?b,sum來給出?b,s,?b,sum可由下式表示:

    合速度勢(shì)?b,sum可用三維點(diǎn)匯表示:

    其點(diǎn)匯強(qiáng)度σsum(t)與擾動(dòng)勢(shì)無量綱作用范圍??(??=?/R)呈負(fù)相關(guān),且根據(jù)試驗(yàn)結(jié)果,隨著入水時(shí)間的增加,空泡最大直徑減小,內(nèi)外側(cè)空泡直徑差異減小,表明其對(duì)空泡擴(kuò)展的約束與時(shí)間呈負(fù)相關(guān).

    根據(jù)有勢(shì)壁面的不可穿越性,合速度勢(shì)需要滿足如下邊界條件:

    其中,tp為空泡的無量綱閉合時(shí)間,tp=Tv/R,T為空泡閉合時(shí)間,可按照文獻(xiàn)[21]給出的經(jīng)驗(yàn)公式求得:根據(jù)試驗(yàn)修正給出.(24)式給出了滿足上述條件的點(diǎn)匯強(qiáng)度:

    其中C1=2.76.

    將(24)式代入(26)式得到合速度勢(shì)為

    根據(jù)勢(shì)流疊加原理推得有勢(shì)壁面的點(diǎn)匯如下:

    線源?c,s可以通過在r=2?處引入一個(gè)相同強(qiáng)度的線源進(jìn)行表示,其源強(qiáng)度由(19)式給出,得到?c,s如下:

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

    速度勢(shì)代入回轉(zhuǎn)體動(dòng)力學(xué)方程,采用經(jīng)典四階龍格-庫塔法求解,得到在有勢(shì)壁面影響下的位移,并將其與單獨(dú)入水位移和試驗(yàn)結(jié)果進(jìn)行比較,如圖11所示.從圖中可以看出,有勢(shì)壁面對(duì)回轉(zhuǎn)體z方向的位移影響很小.

    圖11 位移結(jié)果對(duì)比(V0=1.9 m/s)Fig.11.Comparison of displacement results.

    圖12 計(jì)算結(jié)果與試驗(yàn)結(jié)果對(duì)比(V0=1.9 m/s)Fig.12.Comparison of water-entry cavities of calculative and experimental results.

    基于第3.4節(jié)中對(duì)回轉(zhuǎn)體位移的求解結(jié)果,開展在有勢(shì)壁面約束下的空泡半徑演化計(jì)算,在z方向的離散取?z=0.001 m,圖12給出了有勢(shì)壁面約束下,空泡半徑與試驗(yàn)工況下兩回轉(zhuǎn)體軸線內(nèi)側(cè)空泡半徑的對(duì)比.從圖中可以看出,理論結(jié)果與試驗(yàn)結(jié)果符合度較好,證明了該擾動(dòng)模型的正確性.

    5 三維空泡形態(tài)演化模型

    5.1 非線性影響函數(shù)建立

    由4.1節(jié)分析可知,在空泡的演化中,空泡沿回轉(zhuǎn)體的徑向擴(kuò)張,對(duì)稱面對(duì)空泡的影響區(qū)域如圖13所示,白色區(qū)域空泡均呈現(xiàn)自由擴(kuò)張,僅在灰色區(qū)域空泡會(huì)受到對(duì)稱面的影響.

    圖13 對(duì)稱面對(duì)空泡演化的影響區(qū)域Fig.13.Schematic of the influence region of symmetry plane on the cavity evolution.

    基于上述的影響區(qū)域,我們將柱坐標(biāo)系下的柱面有勢(shì)壁面展開為笛卡爾直角坐標(biāo)系下的有勢(shì)平面,其對(duì)空泡的影響通過影響函數(shù)β來確定,其原理示意如圖14所示.

    設(shè)空泡的實(shí)際半徑為rc,在有勢(shì)壁面完全影響下的空泡半徑為rc1,自由擴(kuò)張半徑為rc2,則其三者的關(guān)系可用下式表示:

    其中,β為影響函數(shù),其具體形式如下文所述.

    圖14 有勢(shì)壁面展開示意圖Fig.14.Schematic of the potential wall expansion.

    以有勢(shì)壁面完全影響下的空泡邊界與展開后的有勢(shì)壁面之間的距離作為影響函數(shù)權(quán)重,在如圖15所示的極坐標(biāo)系下進(jìn)行推導(dǎo),其中l(wèi)1為有勢(shì)壁面,l2為有勢(shì)壁面完全約束下的空泡半徑,r1為回轉(zhuǎn)體中心到圓柱有勢(shì)壁面的距離,r2為回轉(zhuǎn)體中心距離展開后的有勢(shì)壁面的徑向距離.

    圖15 極坐標(biāo)系Fig.15.Polar coordinate system.

    (29)和(30)式給出了l1和l2在極坐標(biāo)系下的曲線方程:

    兩條曲線間的距離如下:

    圖15中,i為有勢(shì)壁面的法向,對(duì)空泡壁面的影響其主要作用的分量為空泡邊界的法向分量in,因此,影響函數(shù)的非線性假設(shè)如下:

    其中,k,b和q均為未知參數(shù).

    上式需要滿足邊界條件:

    此外,在θ=π/2處,即未受影響區(qū)域和受影響區(qū)域交界處,空泡截面邊界應(yīng)滿足連續(xù)可導(dǎo)條件.在笛卡爾坐標(biāo)系下,取1/2空泡截面,未受影響區(qū)域,空泡截面為一1/4圓,其在x=0處?y/?x為0,因此受影響區(qū)域在該處也應(yīng)滿足:

    其中

    得到影響函數(shù)β的表達(dá)式如下:

    5.2 計(jì)算結(jié)果驗(yàn)證

    基于5.1節(jié)中對(duì)影響函數(shù)的推導(dǎo),對(duì)基于影響函數(shù)的三維空泡演化理論模型進(jìn)行驗(yàn)證.圖16給出了基于影響函數(shù)的空泡同步并聯(lián)入水三維空泡形態(tài)計(jì)算結(jié)果及其與試驗(yàn)結(jié)果的對(duì)比.從圖中可以看出,計(jì)算結(jié)果與試驗(yàn)結(jié)果符合較好,說明了該空泡演化模型的正確性.

    圖16 三維空泡理論結(jié)果與試驗(yàn)結(jié)果對(duì)比 (a)t=0.01 s;(b)t=0.02 s;(c)t=0.03 s;(d)t=0.04 s;(e)t=0.05 sFig.16.Comparison of three-dimensional water-entry cavities of calculative and experimental results:(a)t=0.01 s;(b)t=0.02 s;(c)t=0.03 s;(d)t=0.04 s;(e)t=0.05 s.

    圖17 并聯(lián)入水空泡試驗(yàn)與數(shù)值計(jì)算對(duì)比(v0=1.9 m/s) (a)試驗(yàn)結(jié)果;(b)數(shù)值計(jì)算結(jié)果Fig.17.The water-entry cavity in paraller of experimental results and numerical results.

    圖18 并聯(lián)入水位移特性數(shù)值計(jì)算結(jié)果與試驗(yàn)結(jié)果對(duì)比Fig.18.Comparison of the numerical results and experimental results.

    由于試驗(yàn)過程中難以捕捉空泡橫截面形狀,因此我們采用數(shù)值模擬結(jié)果對(duì)三維空泡演化模型的正確性進(jìn)行驗(yàn)證[22].為驗(yàn)證數(shù)值計(jì)算結(jié)果的有效性,圖17給出了回轉(zhuǎn)體以1.9 m/s速度并聯(lián)入水過程中,入水空泡形態(tài)實(shí)驗(yàn)與數(shù)值計(jì)算結(jié)果的對(duì)比.圖18分別給出了回轉(zhuǎn)體并聯(lián)入水過程空泡形態(tài)時(shí)歷演化與角位移、線位移時(shí)歷曲線的對(duì)比.從圖中可以看出,數(shù)值計(jì)算所得的空泡形態(tài)在兩彈軸線內(nèi)側(cè)區(qū)域的發(fā)展程度較試驗(yàn)結(jié)果略大,但空泡形態(tài)輪廓的整體演化與試驗(yàn)符合度較好.位移特性曲線與試驗(yàn)結(jié)果符合較好.綜上分析,可認(rèn)為本文數(shù)值方法適用于回轉(zhuǎn)體并聯(lián)入水過程的研究.

    圖19給出了t=0.05 s時(shí),四個(gè)不同深度處空泡界面理論值與數(shù)值仿真結(jié)果的對(duì)比.從圖中可以看出,計(jì)算結(jié)果與數(shù)值結(jié)果具有較高的符合度,表明引入影響函數(shù)后的三維空泡演化計(jì)算模型的正確性,基于該模型可較為準(zhǔn)確地預(yù)測空泡在三維空間的演化過程.

    圖19 空泡截面形狀理論結(jié)果與數(shù)值結(jié)果對(duì)比(t=0.05 s)Fig.19.Comparison of the cavity cross-section of calculative and experimental results.

    6 并聯(lián)入水空泡三維形態(tài)演化特性研究

    不同于單獨(dú)入水過程的軸對(duì)稱空泡,并聯(lián)入水空泡半徑沿回轉(zhuǎn)體軸線呈非均一分布,且在空泡的不同發(fā)展階段,空泡隨周向的變化亦呈現(xiàn)出不同的變化規(guī)律.圖20給出了極坐標(biāo)系下t=0.05 s時(shí),不同深度處空泡半徑隨極角(范圍取0—π)的變化規(guī)律,此時(shí)空泡已在回轉(zhuǎn)體表面形成閉合,圖20(a)給出了閉合點(diǎn)上方空泡半徑隨極角的變化,圖20(b)給出了閉合點(diǎn)下方空泡半徑隨極角的變化.從圖中可以看出,在閉合點(diǎn)上方,隨著入水深度的增加,空泡半徑隨極角的變化呈現(xiàn)出從正相關(guān)到負(fù)相關(guān)的變化;在閉合點(diǎn)下方,隨著深度的增加,空泡半徑呈現(xiàn)出從負(fù)相關(guān)到正相關(guān)的變化.該現(xiàn)象的產(chǎn)生主要是受空泡頸縮的影響,在閉合點(diǎn)上方,隨著入水深度增加,空泡截面逐漸靠近閉合點(diǎn),此時(shí),自由演化區(qū)空泡(θ=π/2—π)隨著入水深度增加,收縮量增大;抑制演化區(qū),有勢(shì)壁面對(duì)空泡的收縮亦存在抑制作用,空泡收縮速率低于自由演化區(qū),使得在閉合點(diǎn)附近,抑制演化區(qū)的空泡半徑大于自由演化區(qū),進(jìn)而形成了隨著入水深度的增加,空泡半徑隨極角的變化呈現(xiàn)出從正相關(guān)到負(fù)相關(guān)變化的現(xiàn)象.在閉合點(diǎn)下方,靠近閉合點(diǎn)處,即H=6.1D處,空泡截面處在收縮過程中,其半徑隨極角的變化呈負(fù)相關(guān)的原因同上;隨著深度的增大,空泡截面的運(yùn)動(dòng)從收縮過渡至擴(kuò)張,在該過程中,抑制演化區(qū)受有勢(shì)壁面影響,擴(kuò)張受限,使得自由演化區(qū)空泡半徑大于抑制演化區(qū),因而空泡半徑隨極角變化呈現(xiàn)出正相關(guān)特性.

    從圖20中還可以發(fā)現(xiàn),在閉合點(diǎn)上方,無論是自由演化區(qū)還是抑制演化區(qū),空泡半徑隨入水深度的增大呈單調(diào)減小趨勢(shì);在閉合點(diǎn)下方,抑制演化區(qū)空泡半徑隨入水深度的增大單調(diào)減小,而自由演化區(qū)空泡半徑隨深度增加呈先增大后減小的趨勢(shì).該現(xiàn)象與空泡截面的發(fā)展程度有關(guān),在閉合點(diǎn)上方,各空泡截面均已充分發(fā)展并進(jìn)入收縮階段,因此空泡半徑隨入水深度具有較好的單調(diào)性;而在閉合點(diǎn)下方,自由演化區(qū)尚有未完全發(fā)展的空泡截面,因而出現(xiàn)了空泡半徑先增大后減小的趨勢(shì).

    為分析有勢(shì)壁面對(duì)并聯(lián)入水空泡半徑影響的空間變化規(guī)律,定義并聯(lián)入水空泡加權(quán)半徑如下:

    其中,r為并聯(lián)入水空泡半徑,rc2為未受有勢(shì)壁面影響的單獨(dú)入水空泡半徑.

    圖21給出了空泡加權(quán)半徑隨極角的變化規(guī)律,同樣,圖21(a)曲線給出了閉合點(diǎn)上方不同深度處的半徑變化,圖21(b)曲線給出了閉合點(diǎn)下方不同深度處的半徑變化.從圖21(a)中可以看出,在閉合點(diǎn)上方,隨著入水深度的增加,抑制演化區(qū)空泡半徑與自由演化區(qū)空泡半徑的比值逐漸增大,出現(xiàn)該現(xiàn)象的主要原因在于抑制演化區(qū)的空泡收縮率遠(yuǎn)小于自由演化區(qū).而對(duì)于閉合點(diǎn)下方,圖21(b)中可以看出隨著入水深度的增大,抑制演化區(qū)空泡半徑與自由演化區(qū)空泡半徑的比值隨深度增加先減小后增大,且當(dāng)深度較大時(shí),其值非常接近.

    圖20 空泡半徑沿周向變化規(guī)律(t=0.05 s) (a)閉合點(diǎn)上方;(b)閉合點(diǎn)下方Fig.20.The cavity radius of circumferential direction.

    圖21 抑制演化區(qū)空泡加權(quán)半徑沿周向變化規(guī)律(t=0.05 s) (a)閉合點(diǎn)上方;(b)閉合點(diǎn)下方Fig.21.The of the suppressed region.

    圖22給出了t=0.05 s時(shí)不同深度處空泡截面形態(tài),其中圖22(a)為閉合點(diǎn)上方不同深度空泡截面形態(tài),圖22(b)為閉合點(diǎn)下方不同深度空泡截面形態(tài).在閉合點(diǎn)上方,空泡經(jīng)歷了較為充分的演化,使得空泡有了明顯的擴(kuò)張、收縮運(yùn)動(dòng),隨著深度增大,空泡截面逐漸縮小.在H=1.1D處,空泡截面形成了過度演化,即空泡截面邊界在x方向超越了空泡的理論x最大值——空泡對(duì)稱線處的空泡形態(tài)x值.分析其形成原因,在該深度處,收縮勢(shì)較弱,空泡擴(kuò)張程度較大,然而隨著極角的增大,有勢(shì)壁面對(duì)擴(kuò)張的抑制逐漸減弱,在該位置處,空泡的擴(kuò)張程度與抑制程度的比值超過了對(duì)稱線處空泡擴(kuò)張程度與抑制程度的比值,使得空泡形成過度演化.空泡的過度演化只形成于深度較淺的位置,在深度較深的位置,空泡擴(kuò)張程度總體較弱的地方難以形成過度演化,如圖22(b)所示.

    圖22 不同深度處空泡截面形態(tài)(t=0.05 s) (a)閉合點(diǎn)上方;(b)閉合點(diǎn)下方Fig.22.The cross section shape of various depth.

    圖23分別給出了H=1.3D,H=3.3D,H=5.3D三個(gè)深度處極坐標(biāo)系下不同時(shí)刻的空泡形態(tài).從圖中可以看出,隨著時(shí)間的增加,極坐標(biāo)系下的空泡半徑呈先增大后減小的趨勢(shì),且增大速率隨著極角的增大而增大.該過程反映了空泡界面的經(jīng)典演化特性.在空泡發(fā)生收縮時(shí),由于影響函數(shù)在極角越大處產(chǎn)生的收縮抑制越小,使得極角較大的位置率先發(fā)生收縮,而此時(shí)極角為0處的空泡尚未形成收縮.隨著深度增大,極角較大位置的空泡收縮速率越快.

    圖23 極坐標(biāo)系下空泡截面半徑時(shí)歷演化 (a)H=1.3D;(b)H=3.3D;(c)H=5.3DFig.23.Time evolution of cavity cross section shape in polar coordinate system:(a)H=1.3D;(b)H=3.3D;(c)H=5.3D.

    圖24給出了極坐標(biāo)系下不同時(shí)刻的空泡加權(quán)半徑隨極角的變化規(guī)律.從圖中可以看出,空泡加權(quán)半徑與空泡截面半徑隨時(shí)間的變化規(guī)律基本一致.此外,從圖中還可以看出,在空泡從擴(kuò)張向收縮過渡的過程中,空泡加權(quán)半徑逐漸靠近,在H=3.3D處,t=0.03 s與t=0.04 s的加權(quán)半徑已接近完全重合.

    圖25給出了笛卡爾坐標(biāo)系下的空泡形態(tài)時(shí)歷演化,從圖25(a)中可以發(fā)現(xiàn),在空泡形成之初,空泡尚未出現(xiàn)過度演化現(xiàn)象,隨著時(shí)間的推移,空泡過度演化逐漸產(chǎn)生,其程度也逐漸增大,在空泡的收縮過程中,過度演化逐漸得到緩和.而隨著深度的增大,在空泡演化的整個(gè)過程中,均沒有過度演化現(xiàn)象產(chǎn)生.

    圖24 極坐標(biāo)系下空泡加權(quán)半徑時(shí)歷演化 (a)H=1.3D;(b)H=3.3D;(c)H=5.3DFig.24.Time evolution of in polar coordinate system:(a)H=1.3D;(b)H=3.3D;(c)H=5.3D.

    圖25 笛卡爾坐標(biāo)系下空泡截面形態(tài)時(shí)歷演化 (a)H=1.3D;(b)H=3.3D;(c)H=5.3D.Fig.25.Time evolution of in Cartesian coordinate system:(a)H=1.3D;(b)H=3.3D;(c)H=5.3D.

    7 結(jié) 論

    本文基于勢(shì)流理論,通過二維軸對(duì)稱方法分析了回轉(zhuǎn)體單獨(dú)入水過程空泡動(dòng)力學(xué)模型.基于回轉(zhuǎn)體單獨(dú)入水空泡的動(dòng)力學(xué)模型,引入有勢(shì)壁面和影響函數(shù),建立有勢(shì)壁面影響下的空泡受迫演化模型,并基于該模型分析了并聯(lián)入水空泡三維演化特性,得到如下結(jié)論.

    1)回轉(zhuǎn)體入水過程流場速度勢(shì)可以看作由一個(gè)隨回轉(zhuǎn)體運(yùn)動(dòng)的點(diǎn)源和位于空泡軸線處的線源疊加產(chǎn)生,其中空泡的擴(kuò)張運(yùn)動(dòng)主要受點(diǎn)源的影響,而收縮運(yùn)動(dòng)主要受線源影響.

    2)在并聯(lián)入水過程中,雙空泡演化在空間呈鏡面對(duì)稱,空泡間的相互擾動(dòng)可以通過引入有勢(shì)壁面進(jìn)行分析.有勢(shì)壁面對(duì)空泡的演化呈抑制作用,該抑制作用包含兩個(gè)方面,分別為對(duì)擴(kuò)張的抑制和對(duì)收縮的抑制.

    3)并聯(lián)入水空泡半徑隨極角的變化與空泡截面所處深度有關(guān),在靠近閉合點(diǎn)附近的抑制演化區(qū)空泡截面半徑隨極角的增大而逐漸減小,遠(yuǎn)離閉合點(diǎn)處的抑制演化區(qū)空泡截面半徑隨極角的增大而增大,空泡與半徑規(guī)律相反.

    4)在深度較淺處,空泡擴(kuò)張過程中會(huì)形成過度演化,空泡收縮階段過度演化逐漸被削弱并消失.

    猜你喜歡
    空泡壁面并聯(lián)
    二維有限長度柔性壁面上T-S波演化的數(shù)值研究
    識(shí)別串、并聯(lián)電路的方法
    水下航行體雙空泡相互作用數(shù)值模擬研究
    審批由“串聯(lián)”改“并聯(lián)”好在哪里?
    壁面溫度對(duì)微型內(nèi)燃機(jī)燃燒特性的影響
    基于LPV的超空泡航行體H∞抗飽和控制
    并聯(lián)型APF中SVPWM的零矢量分配
    一種軟開關(guān)的交錯(cuò)并聯(lián)Buck/Boost雙向DC/DC變換器
    基于CFD的對(duì)轉(zhuǎn)槳無空泡噪聲的仿真預(yù)報(bào)
    船海工程(2015年4期)2016-01-05 15:53:28
    SPH在水下高速物體空泡發(fā)展模擬中的應(yīng)用
    国产一区二区三区在线臀色熟女| x7x7x7水蜜桃| 国产成人系列免费观看| 一级黄色大片毛片| 小蜜桃在线观看免费完整版高清| 色综合欧美亚洲国产小说| 在线观看日韩欧美| 视频区欧美日本亚洲| 日韩免费av在线播放| 亚洲av片天天在线观看| av黄色大香蕉| 成年女人毛片免费观看观看9| 男女视频在线观看网站免费| 国产成人系列免费观看| 亚洲精品在线美女| 久久久久久久午夜电影| 欧美日本亚洲视频在线播放| 精品一区二区三区视频在线观看免费| 国产精品乱码一区二三区的特点| 怎么达到女性高潮| 亚洲第一电影网av| 国内精品久久久久久久电影| a在线观看视频网站| 免费看美女性在线毛片视频| 国产av不卡久久| 男女床上黄色一级片免费看| 国产免费av片在线观看野外av| 99国产极品粉嫩在线观看| 亚洲中文字幕一区二区三区有码在线看 | 亚洲狠狠婷婷综合久久图片| 国产一区二区在线av高清观看| 久久精品夜夜夜夜夜久久蜜豆| 好男人电影高清在线观看| 亚洲熟女毛片儿| 午夜日韩欧美国产| 精品国产乱子伦一区二区三区| 亚洲avbb在线观看| 亚洲人与动物交配视频| 一本一本综合久久| 欧美成狂野欧美在线观看| 午夜免费激情av| 国产成人aa在线观看| 99久久无色码亚洲精品果冻| 国产一区二区三区在线臀色熟女| 成人特级av手机在线观看| 国产欧美日韩精品亚洲av| 亚洲国产欧洲综合997久久,| 97超视频在线观看视频| 久久精品aⅴ一区二区三区四区| av国产免费在线观看| 国产主播在线观看一区二区| 国产av在哪里看| 国产精品1区2区在线观看.| 一级毛片精品| 亚洲天堂国产精品一区在线| 99热这里只有精品一区 | 91老司机精品| 真人做人爱边吃奶动态| 一区二区三区国产精品乱码| 夜夜看夜夜爽夜夜摸| 久久精品91蜜桃| 又紧又爽又黄一区二区| 成人永久免费在线观看视频| 久久久精品大字幕| 亚洲男人的天堂狠狠| 熟女少妇亚洲综合色aaa.| 亚洲成a人片在线一区二区| 亚洲人成伊人成综合网2020| 久久久久国内视频| 搡老妇女老女人老熟妇| 女生性感内裤真人,穿戴方法视频| 丰满的人妻完整版| 国产一区二区三区视频了| 婷婷亚洲欧美| 国产伦一二天堂av在线观看| 黄色丝袜av网址大全| 亚洲国产色片| 日本一二三区视频观看| 黄频高清免费视频| 69av精品久久久久久| 一个人观看的视频www高清免费观看 | 色综合站精品国产| 亚洲欧美一区二区三区黑人| cao死你这个sao货| 母亲3免费完整高清在线观看| а√天堂www在线а√下载| 国产一区二区三区在线臀色熟女| 91在线精品国自产拍蜜月 | 91久久精品国产一区二区成人 | 丰满的人妻完整版| 亚洲国产欧洲综合997久久,| 国产精品女同一区二区软件 | 国产日本99.免费观看| 亚洲精品456在线播放app | 成人三级做爰电影| 精品午夜福利视频在线观看一区| 成人三级做爰电影| 性欧美人与动物交配| 99久国产av精品| 美女高潮的动态| 级片在线观看| 亚洲av五月六月丁香网| 欧美中文综合在线视频| www国产在线视频色| 国产精品美女特级片免费视频播放器 | 国产高清有码在线观看视频| 久久久久久国产a免费观看| 日本精品一区二区三区蜜桃| 久久亚洲精品不卡| 国产一区在线观看成人免费| 男女下面进入的视频免费午夜| 亚洲av免费在线观看| www.999成人在线观看| 一二三四社区在线视频社区8| 麻豆久久精品国产亚洲av| 国产毛片a区久久久久| 国产精品永久免费网站| 丰满人妻熟妇乱又伦精品不卡| 制服人妻中文乱码| 精品久久久久久久久久免费视频| 国产高清三级在线| 脱女人内裤的视频| 色综合婷婷激情| 99久久99久久久精品蜜桃| 网址你懂的国产日韩在线| 亚洲中文字幕一区二区三区有码在线看 | 天天躁狠狠躁夜夜躁狠狠躁| 日韩欧美精品v在线| 国产精品香港三级国产av潘金莲| 久久精品aⅴ一区二区三区四区| 深夜精品福利| 国产一区二区三区在线臀色熟女| 一个人免费在线观看电影 | 脱女人内裤的视频| 曰老女人黄片| 一进一出抽搐gif免费好疼| 国产爱豆传媒在线观看| 午夜影院日韩av| 免费看a级黄色片| 91av网一区二区| 变态另类丝袜制服| 成在线人永久免费视频| 免费看十八禁软件| 亚洲国产色片| 桃色一区二区三区在线观看| 美女高潮的动态| 久久久久精品国产欧美久久久| 欧洲精品卡2卡3卡4卡5卡区| 成人欧美大片| 又黄又爽又免费观看的视频| av天堂中文字幕网| 亚洲一区二区三区色噜噜| 老汉色∧v一级毛片| 国产主播在线观看一区二区| 在线免费观看的www视频| 亚洲18禁久久av| 国产精品日韩av在线免费观看| 一个人免费在线观看电影 | 中文字幕高清在线视频| 欧美日韩瑟瑟在线播放| 手机成人av网站| 日本黄色视频三级网站网址| 日韩欧美 国产精品| 久久久久亚洲av毛片大全| 老熟妇仑乱视频hdxx| 99国产极品粉嫩在线观看| 日本五十路高清| 最好的美女福利视频网| 88av欧美| 亚洲精品国产精品久久久不卡| 免费观看人在逋| 日韩高清综合在线| 无人区码免费观看不卡| 两人在一起打扑克的视频| 伦理电影免费视频| 黑人欧美特级aaaaaa片| 国产伦一二天堂av在线观看| 国产精品98久久久久久宅男小说| 在线观看免费视频日本深夜| 国产99白浆流出| 伊人久久大香线蕉亚洲五| 久久人人精品亚洲av| 少妇人妻一区二区三区视频| 午夜a级毛片| 色综合站精品国产| 免费大片18禁| 国产精品美女特级片免费视频播放器 | 中文字幕精品亚洲无线码一区| 国产单亲对白刺激| 国产又黄又爽又无遮挡在线| 小蜜桃在线观看免费完整版高清| 少妇熟女aⅴ在线视频| 成人av一区二区三区在线看| 热99在线观看视频| 国产精品自产拍在线观看55亚洲| 麻豆av在线久日| 老司机午夜十八禁免费视频| 国产高清三级在线| 网址你懂的国产日韩在线| 精品久久久久久,| 日韩欧美三级三区| 欧美激情在线99| x7x7x7水蜜桃| 亚洲av片天天在线观看| 美女被艹到高潮喷水动态| 小说图片视频综合网站| 村上凉子中文字幕在线| 亚洲熟妇中文字幕五十中出| 欧美成人性av电影在线观看| 男插女下体视频免费在线播放| 国产成人欧美在线观看| 99久久99久久久精品蜜桃| 此物有八面人人有两片| 国产精品自产拍在线观看55亚洲| 美女高潮的动态| 色老头精品视频在线观看| 天天躁狠狠躁夜夜躁狠狠躁| 欧美黄色片欧美黄色片| 91在线精品国自产拍蜜月 | 噜噜噜噜噜久久久久久91| 狠狠狠狠99中文字幕| 性色avwww在线观看| 国内精品久久久久久久电影| 国产久久久一区二区三区| 国产激情欧美一区二区| 国产免费男女视频| 黑人巨大精品欧美一区二区mp4| 999久久久国产精品视频| 国产真人三级小视频在线观看| 亚洲av美国av| 长腿黑丝高跟| 欧美成人性av电影在线观看| 1024手机看黄色片| 亚洲av电影不卡..在线观看| 中文字幕久久专区| 老司机深夜福利视频在线观看| 999久久久精品免费观看国产| 成人三级黄色视频| www.999成人在线观看| 麻豆av在线久日| 成人亚洲精品av一区二区| 国产伦精品一区二区三区四那| 国语自产精品视频在线第100页| 日韩中文字幕欧美一区二区| 美女大奶头视频| 国产午夜精品久久久久久| 最近最新免费中文字幕在线| 在线免费观看不下载黄p国产 | 亚洲第一电影网av| 国产伦在线观看视频一区| 啦啦啦免费观看视频1| 日韩av在线大香蕉| 午夜激情欧美在线| 亚洲第一电影网av| 亚洲av成人一区二区三| 夜夜看夜夜爽夜夜摸| 日韩欧美精品v在线| 午夜福利免费观看在线| 最近在线观看免费完整版| 亚洲欧美精品综合久久99| 看免费av毛片| 在线a可以看的网站| 2021天堂中文幕一二区在线观| 99精品欧美一区二区三区四区| 色综合站精品国产| tocl精华| 中出人妻视频一区二区| 亚洲 欧美一区二区三区| 小说图片视频综合网站| av国产免费在线观看| 九九热线精品视视频播放| www.www免费av| 久久久久国产一级毛片高清牌| 天堂av国产一区二区熟女人妻| 全区人妻精品视频| www日本黄色视频网| av在线天堂中文字幕| 亚洲国产色片| 最新美女视频免费是黄的| 特级一级黄色大片| 91字幕亚洲| 99在线人妻在线中文字幕| 精品久久久久久久末码| 成人三级做爰电影| 熟妇人妻久久中文字幕3abv| 国产精品影院久久| 国产一区二区三区在线臀色熟女| 熟女人妻精品中文字幕| 亚洲av成人精品一区久久| 九九在线视频观看精品| 99久久精品热视频| 亚洲aⅴ乱码一区二区在线播放| 99热6这里只有精品| 天堂网av新在线| 麻豆成人av在线观看| 日本 av在线| www.www免费av| 99热精品在线国产| av中文乱码字幕在线| 亚洲中文av在线| 久久久色成人| 两个人看的免费小视频| 九九在线视频观看精品| 亚洲美女视频黄频| 亚洲欧美日韩卡通动漫| 国产精品亚洲av一区麻豆| 亚洲成av人片在线播放无| 欧美一区二区精品小视频在线| 波多野结衣高清作品| 日韩高清综合在线| 最近最新中文字幕大全电影3| 精品福利观看| 巨乳人妻的诱惑在线观看| 岛国在线免费视频观看| 特大巨黑吊av在线直播| 国产精品国产高清国产av| 国产真人三级小视频在线观看| 国产伦人伦偷精品视频| x7x7x7水蜜桃| 嫁个100分男人电影在线观看| 日本在线视频免费播放| 手机成人av网站| 日韩 欧美 亚洲 中文字幕| a级毛片在线看网站| 高清在线国产一区| 在线观看免费视频日本深夜| 黄片小视频在线播放| 国产亚洲av嫩草精品影院| 国产综合懂色| 欧美日韩黄片免| 国产欧美日韩精品亚洲av| 一本一本综合久久| 欧美大码av| 在线a可以看的网站| 国产高清视频在线观看网站| 美女cb高潮喷水在线观看 | 日韩精品青青久久久久久| 久久久水蜜桃国产精品网| 香蕉国产在线看| 日本一本二区三区精品| 999久久久精品免费观看国产| 色吧在线观看| 后天国语完整版免费观看| 亚洲欧美日韩高清在线视频| 国产男靠女视频免费网站| 久久久久久久精品吃奶| 人妻久久中文字幕网| 欧美精品啪啪一区二区三区| 99久久精品热视频| 欧美色视频一区免费| 国产熟女xx| 男女下面进入的视频免费午夜| a级毛片在线看网站| 黑人欧美特级aaaaaa片| 精品电影一区二区在线| 在线观看美女被高潮喷水网站 | 搞女人的毛片| 国产伦精品一区二区三区视频9 | 亚洲激情在线av| 亚洲人成电影免费在线| 9191精品国产免费久久| www日本黄色视频网| 桃红色精品国产亚洲av| 99精品在免费线老司机午夜| 又粗又爽又猛毛片免费看| 最近视频中文字幕2019在线8| 免费看日本二区| 白带黄色成豆腐渣| 精品国产美女av久久久久小说| 给我免费播放毛片高清在线观看| 久久久精品大字幕| 搡老岳熟女国产| 欧美日韩乱码在线| 精品午夜福利视频在线观看一区| 国产精品98久久久久久宅男小说| 亚洲熟妇熟女久久| aaaaa片日本免费| 看免费av毛片| 最好的美女福利视频网| 欧美成人性av电影在线观看| 少妇的丰满在线观看| 欧美不卡视频在线免费观看| 免费电影在线观看免费观看| 麻豆成人午夜福利视频| 亚洲中文日韩欧美视频| 日韩精品青青久久久久久| 亚洲国产精品成人综合色| 亚洲人成网站高清观看| 久久久久亚洲av毛片大全| 亚洲欧美一区二区三区黑人| 国产精品永久免费网站| 国产爱豆传媒在线观看| 久久久久久久午夜电影| 女警被强在线播放| 亚洲一区二区三区不卡视频| 成人三级黄色视频| 最近最新中文字幕大全免费视频| 久久久久精品国产欧美久久久| 欧美一区二区精品小视频在线| 日本免费a在线| 99久久成人亚洲精品观看| 动漫黄色视频在线观看| 视频区欧美日本亚洲| 老司机福利观看| 中文亚洲av片在线观看爽| 床上黄色一级片| 国产久久久一区二区三区| ponron亚洲| 男女视频在线观看网站免费| 国产精品野战在线观看| 男女下面进入的视频免费午夜| 国产亚洲av嫩草精品影院| 搡老熟女国产l中国老女人| 久久亚洲精品不卡| 香蕉国产在线看| 欧美日韩国产亚洲二区| 色综合亚洲欧美另类图片| 欧美日韩一级在线毛片| 又黄又粗又硬又大视频| 久久精品国产清高在天天线| 色视频www国产| av福利片在线观看| 亚洲国产欧美人成| 亚洲欧美精品综合久久99| 99视频精品全部免费 在线 | 免费看日本二区| 99视频精品全部免费 在线 | 夜夜躁狠狠躁天天躁| 99国产极品粉嫩在线观看| 性色avwww在线观看| 真人做人爱边吃奶动态| 91麻豆精品激情在线观看国产| 国产视频内射| 成人av一区二区三区在线看| 久久精品91无色码中文字幕| 99久久精品一区二区三区| 美女午夜性视频免费| 操出白浆在线播放| 亚洲欧美精品综合一区二区三区| 美女黄网站色视频| e午夜精品久久久久久久| 精品一区二区三区视频在线观看免费| 午夜亚洲福利在线播放| 最新美女视频免费是黄的| 午夜免费成人在线视频| 十八禁人妻一区二区| 99国产精品99久久久久| av欧美777| 午夜激情福利司机影院| 国产精品99久久99久久久不卡| 变态另类丝袜制服| 99久久成人亚洲精品观看| 女人高潮潮喷娇喘18禁视频| 欧美乱妇无乱码| 成人国产综合亚洲| 色综合站精品国产| 久久中文字幕一级| 非洲黑人性xxxx精品又粗又长| 亚洲人成网站在线播放欧美日韩| 我要搜黄色片| 99热这里只有精品一区 | 男人舔女人下体高潮全视频| 成人av在线播放网站| av福利片在线观看| 一个人免费在线观看电影 | 宅男免费午夜| 人妻久久中文字幕网| 18禁观看日本| 51午夜福利影视在线观看| 夜夜夜夜夜久久久久| 在线观看免费视频日本深夜| АⅤ资源中文在线天堂| 最近在线观看免费完整版| 欧美激情久久久久久爽电影| 床上黄色一级片| 丰满人妻熟妇乱又伦精品不卡| 久久久久亚洲av毛片大全| 黄片大片在线免费观看| 老熟妇乱子伦视频在线观看| 男女之事视频高清在线观看| 国产伦精品一区二区三区四那| 国产精品,欧美在线| 亚洲avbb在线观看| 婷婷精品国产亚洲av在线| 国产亚洲av嫩草精品影院| 亚洲欧美日韩东京热| 国产精品九九99| 国产综合懂色| 国产蜜桃级精品一区二区三区| 校园春色视频在线观看| 97碰自拍视频| 给我免费播放毛片高清在线观看| 岛国在线观看网站| 淫秽高清视频在线观看| 观看免费一级毛片| 少妇人妻一区二区三区视频| 国产真实乱freesex| 国产精品1区2区在线观看.| 色综合站精品国产| 天堂影院成人在线观看| 亚洲乱码一区二区免费版| 久久久色成人| 激情在线观看视频在线高清| 欧美性猛交╳xxx乱大交人| 亚洲人成网站高清观看| 国产真人三级小视频在线观看| 18禁黄网站禁片午夜丰满| 国产视频一区二区在线看| 日本a在线网址| 日韩欧美在线乱码| 久久精品影院6| www.熟女人妻精品国产| 神马国产精品三级电影在线观看| 性色avwww在线观看| 日本熟妇午夜| 亚洲狠狠婷婷综合久久图片| 高潮久久久久久久久久久不卡| 久久精品aⅴ一区二区三区四区| 桃色一区二区三区在线观看| 国产高清视频在线播放一区| 丝袜人妻中文字幕| 亚洲最大成人中文| 琪琪午夜伦伦电影理论片6080| 最新美女视频免费是黄的| 99re在线观看精品视频| 精品一区二区三区四区五区乱码| 国产精品久久久久久久电影 | 午夜免费观看网址| 特大巨黑吊av在线直播| av天堂在线播放| 免费观看的影片在线观看| 亚洲av五月六月丁香网| 我要搜黄色片| 小蜜桃在线观看免费完整版高清| 99精品欧美一区二区三区四区| 熟女少妇亚洲综合色aaa.| 国产三级在线视频| 99精品久久久久人妻精品| 久久香蕉精品热| 夜夜爽天天搞| 在线十欧美十亚洲十日本专区| 亚洲狠狠婷婷综合久久图片| 看免费av毛片| 日韩欧美三级三区| 高清在线国产一区| 最新中文字幕久久久久 | 天堂av国产一区二区熟女人妻| 97超视频在线观看视频| 黄色丝袜av网址大全| 亚洲激情在线av| 国内久久婷婷六月综合欲色啪| 丁香欧美五月| 欧美色欧美亚洲另类二区| 不卡一级毛片| 午夜精品在线福利| 亚洲中文字幕一区二区三区有码在线看 | 欧美在线一区亚洲| 又黄又粗又硬又大视频| 黄频高清免费视频| 亚洲国产精品合色在线| 亚洲欧美一区二区三区黑人| 欧美精品啪啪一区二区三区| 亚洲av第一区精品v没综合| 婷婷亚洲欧美| 18禁裸乳无遮挡免费网站照片| 又紧又爽又黄一区二区| 啪啪无遮挡十八禁网站| 国产精品一区二区三区四区免费观看 | aaaaa片日本免费| 日本黄色视频三级网站网址| 精品久久久久久久末码| 黄色日韩在线| 一个人观看的视频www高清免费观看 | 中文资源天堂在线| 久久精品综合一区二区三区| 欧美最黄视频在线播放免费| 亚洲精品粉嫩美女一区| 免费一级毛片在线播放高清视频| 欧美日韩中文字幕国产精品一区二区三区| 香蕉久久夜色| 免费人成视频x8x8入口观看| 麻豆成人av在线观看| 欧美xxxx黑人xx丫x性爽| 国产精品一区二区精品视频观看| 久久精品91无色码中文字幕| 亚洲人成网站在线播放欧美日韩| 欧美成狂野欧美在线观看| 国产一级毛片七仙女欲春2| 毛片女人毛片| 老熟妇仑乱视频hdxx| 午夜激情福利司机影院| 男人舔奶头视频| 色视频www国产| 97人妻精品一区二区三区麻豆| 欧美色欧美亚洲另类二区| 97人妻精品一区二区三区麻豆| 精品99又大又爽又粗少妇毛片 | 亚洲成人中文字幕在线播放| 欧美成人免费av一区二区三区| 男女下面进入的视频免费午夜| 日本三级黄在线观看| 午夜免费观看网址| 国产成人av激情在线播放| 变态另类丝袜制服| 亚洲国产日韩欧美精品在线观看 | 国产成人精品无人区| 国产激情欧美一区二区| 天堂√8在线中文| 免费一级毛片在线播放高清视频| 亚洲中文av在线| 国产亚洲精品av在线| 亚洲精品中文字幕一二三四区| 欧美日韩综合久久久久久 | 国产乱人视频| 狂野欧美白嫩少妇大欣赏|