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

    高超聲速進氣道強制轉(zhuǎn)捩流動的大渦模擬

    2016-11-14 00:41:17趙曉慧鄧小兵毛枚良楊偉趙慧勇
    航空學報 2016年8期
    關(guān)鍵詞:大渦進氣道邊界層

    趙曉慧, 鄧小兵,*, 毛枚良, 楊偉, 趙慧勇

    1.中國空氣動力研究與發(fā)展中心 空氣動力學國家重點實驗室, 綿陽 621000 2.中國空氣動力研究與發(fā)展中心 計算空氣動力研究所, 綿陽 621000 3.中國空氣動力研究與發(fā)展中心 高超聲速沖壓發(fā)動機重點實驗室, 綿陽 621000

    ?

    高超聲速進氣道強制轉(zhuǎn)捩流動的大渦模擬

    趙曉慧1, 鄧小兵1,*, 毛枚良2, 楊偉1, 趙慧勇3

    1.中國空氣動力研究與發(fā)展中心 空氣動力學國家重點實驗室, 綿陽621000 2.中國空氣動力研究與發(fā)展中心 計算空氣動力研究所, 綿陽621000 3.中國空氣動力研究與發(fā)展中心 高超聲速沖壓發(fā)動機重點實驗室, 綿陽621000

    吸氣式高超聲速飛行器常在進氣道邊界層內(nèi)布置粗糙顆?;驕u流發(fā)生器強制流動轉(zhuǎn)捩為湍流以確保發(fā)動機正常啟動。為了清晰認識強制轉(zhuǎn)捩過程,采用隱式大渦模擬方法,對強制轉(zhuǎn)捩問題開展了數(shù)值模擬研究。針對鉆石形和斜坡形渦流發(fā)生器,計算得到渦流發(fā)生器誘導(dǎo)的流動結(jié)構(gòu),顯示出強制轉(zhuǎn)捩流動由渦流發(fā)生器產(chǎn)生的反向旋轉(zhuǎn)流向渦對主導(dǎo)。擾動沿流向增長和發(fā)展,導(dǎo)致流向渦對以偶模式或奇模式失穩(wěn),偶模式失穩(wěn)產(chǎn)生對稱形式的渦對破碎,而奇模式失穩(wěn)則導(dǎo)致非對稱(彎曲)形式的渦對破碎。流向渦對破碎后產(chǎn)生一系列發(fā)卡渦并最終促使邊界層轉(zhuǎn)捩為湍流。最后就計算網(wǎng)格和數(shù)值耗散對隱式大渦模擬結(jié)果的影響以及計算的收斂性進行了討論。

    高超聲速; 進氣道; 強制轉(zhuǎn)捩; 渦流發(fā)生器; 大渦模擬

    吸氣式高超聲速飛行器通常巡航在高度比較高,而密度則相對較低的大氣環(huán)境中,此時飛行雷諾數(shù)相對較小,自然轉(zhuǎn)捩發(fā)生的位置通常超出了進氣道長度。為了滿足發(fā)動機設(shè)計對邊界層湍流流態(tài)的要求,需要在進氣道采用強制轉(zhuǎn)捩措施。2004年,Hyper-X (X-43A)采用渦流發(fā)生器在飛行試驗中成功實現(xiàn)了強制轉(zhuǎn)捩[1],測試結(jié)果顯示,對于馬赫數(shù)7和10的飛行條件,飛行器前體需要采用粗糙帶來保證進氣道上游的全湍流狀態(tài)。

    渦流發(fā)生器強制轉(zhuǎn)捩的作用機理圖像尚未完全清晰,最早的一些研究者試圖用Tollmien-Schlichting (T-S)波解釋粗糙顆粒強制轉(zhuǎn)捩機理[2],Reshotko和Tumin[3]則認為瞬時增長在粗糙顆粒強制轉(zhuǎn)捩中起了重要作用。數(shù)值計算結(jié)果顯示強制轉(zhuǎn)捩過程與渦流發(fā)生器產(chǎn)生的流向渦對的發(fā)展有關(guān),這與G?rtler渦的二次失穩(wěn)相似,在文獻[4]中這種粗糙顆粒引起的流向反向旋轉(zhuǎn)的渦對也被稱作G?rtler渦。Li等[5-6]在研究低速和高速流動G?rtler渦時指出其二次失穩(wěn)主要有奇模式(Odd mode)和偶模式(Even mode)兩種,分別導(dǎo)致了彎曲(Sinuous)和對稱(Varicose)兩種不同形式的流向渦對破碎,并且其不穩(wěn)定模式與G?rtler渦波長有關(guān)。Tullio等[7]利用線性穩(wěn)定性理論以及求解Navier-Stokes方程研究超聲速平板利用顆粒強制轉(zhuǎn)捩時也指出,對稱和彎曲是兩種最容易出現(xiàn)的不穩(wěn)定形態(tài)。

    高超聲速邊界層的渦流發(fā)生器強制轉(zhuǎn)捩過程受渦流發(fā)生器的形狀、尺寸、分布和環(huán)境噪聲等諸多因素的影響。Berry等[8-10]對五種形狀的渦流發(fā)生器作了試驗比較,發(fā)現(xiàn)鉆石形和斜坡形渦流發(fā)生器具有較好的轉(zhuǎn)捩觸發(fā)效果。趙慧勇等[11]在中國空氣動力研究與發(fā)展中心的FL-31(直徑0.5 m)超聲速風洞中,開展了鉆石形和斜坡形渦流發(fā)生器強制轉(zhuǎn)捩研究,給出了鉆石形渦流發(fā)生器觸發(fā)轉(zhuǎn)捩的有效高度。Danehy等[12]采用平面激光誘導(dǎo)熒光(Plane Laser-Induced Fluorescence, PLIF) 技術(shù)研究了平板表面半球凸起物對邊界層的影響,岡敦殿等[13]采用基于納米示蹤的平面激光散射技術(shù)(Nano-tracer Planar Laser Scattering, NPLS)進行了超聲速平板上直立圓臺凸起物繞流流場研究,張慶虎[14]采用NPLS技術(shù)對三角楔渦流發(fā)生器產(chǎn)生的流動結(jié)構(gòu)給過較為細致的展示,這些研究都給出邊界層中突起物尾跡某一平面上的精細結(jié)構(gòu),對研究突起物尾跡的發(fā)展有很大幫助。Borg和Schneider[15]在靜聲風洞中開展了X-51前體強制轉(zhuǎn)捩研究,指出來流噪聲會使轉(zhuǎn)捩位置提前;鄧小兵[16]采用隱式大渦模擬(Implicit Large Eddy Simulation, ILES)方法,通過在計算來流中加噪聲信號來模擬試驗條件,發(fā)現(xiàn)相同的現(xiàn)象。這表明,要通過試驗準確預(yù)測轉(zhuǎn)捩位置需采用靜聲風洞以消除風洞噪聲的影響。

    相比于試驗研究,數(shù)值計算可以獲得更為全面的流場信息,而且不受風洞噪聲的影響。研究強制轉(zhuǎn)捩問題采用的數(shù)值模擬方法主要有直接數(shù)值模擬(Direct Numerical Simulation, DNS),大渦模擬(Large Eddy Simulation, LES)和雷諾平均Navier-Stokes(RANS)方程方法。其中,DNS方法最準確,Marxen[17]、Iyer[18]和Duan[19-20]等利用DNS對突起物強制轉(zhuǎn)捩流動進行研究,獲得了清晰的流動結(jié)構(gòu)和轉(zhuǎn)捩圖像。DNS方法的缺點是所需的計算資源過多,在當前的計算機硬件條件下還難以推廣應(yīng)用。涂國華等[21]曾經(jīng)利用RANS方法對強制轉(zhuǎn)捩流動進行模擬,結(jié)果表明RANS方法難以預(yù)測轉(zhuǎn)捩位置,且對轉(zhuǎn)捩后壓縮面上熱流分布的預(yù)測也不理想。周玲等[22-23]改進了k-ω-γ轉(zhuǎn)捩模式,并對高超聲速飛行器進氣道前體邊界層強制轉(zhuǎn)捩進行數(shù)值計算,指出改進模型對轉(zhuǎn)捩位置有較好的預(yù)測能力。然而采用RANS方法難以準確捕捉強制轉(zhuǎn)捩的流動結(jié)構(gòu)。LES方法能夠以遠低于DNS的計算開銷給出較為準確的非定常流動結(jié)構(gòu)。Sayadi和Moin[24]利用大渦模擬方法對馬赫數(shù)0.2可壓縮平板邊界層的H-type和K-type轉(zhuǎn)捩進行了模擬,指出渦黏性在層流和轉(zhuǎn)捩區(qū)的存在會阻礙擾動的發(fā)展,不利于轉(zhuǎn)捩預(yù)測。ILES不加入亞格子模型,在層流區(qū)不會引入亞格子耗散,能夠更準確地模擬轉(zhuǎn)捩區(qū)的擾動增長。Orlik和Fedioun[25]采用基于5階WENO格式的隱式大渦模擬方法計算了射流控制的高超聲速邊界層強制轉(zhuǎn)捩流動,結(jié)果顯示ILES能夠捕捉到豐富的流動細節(jié)。本文重點關(guān)注渦流發(fā)生器誘導(dǎo)強制轉(zhuǎn)捩過程的模擬,因而采用ILES是合適的。

    通常認為,大渦模擬方法需要采用高階精度格式,鄧小兵和金玲[26]通過對湍流數(shù)值模擬中誤差傳播和發(fā)展的動力學行為的分析,認為二階精度為基礎(chǔ)的數(shù)值方法可以應(yīng)用于湍流大渦模擬。本文采用二階有限體積隱式大渦模擬方法,對鉆石形和斜坡形渦流發(fā)生器激發(fā)的高超聲速進氣道邊界層轉(zhuǎn)捩現(xiàn)象開展數(shù)值模擬研究,計算得到了流向渦失穩(wěn)的兩種基本模式:偶模式與奇模式,同時給出了脈動量沿流向的增長規(guī)律,最后對數(shù)值耗散和網(wǎng)格分辨率共同作用下ILES方法的收斂性作了探討。

    1 研究方法

    1.1隱式大渦模擬方法

    采用基于二階有限體積的ILES方法研究強制轉(zhuǎn)捩問題。由于ILES采用數(shù)值黏性模擬亞格子耗散,數(shù)值耗散的控制就十分重要。對流項空間離散格式的數(shù)值耗散對大渦模擬結(jié)果有明顯的影響,過大的數(shù)值耗散往往導(dǎo)致能譜過度耗散,抑制不穩(wěn)定波的初始增長,所以數(shù)值耗散的控制對轉(zhuǎn)捩位置的準確預(yù)測至關(guān)重要。采用基于MUSCL( Monotone Upstream-centered Scheme for Conservation Laws)重構(gòu)的ROE格式,在光滑區(qū),該重構(gòu)表達式退化為

    (1)

    式中:qL,i+1/2和qR,i+1/2分別為i和i+1單元界面變量左側(cè)和右側(cè)重構(gòu)值;當κ=1/3時為常用的三階格式,此時格式對ILES是過耗散的,κ=1時格式變?yōu)槎A中心格式,是無耗散的,通過調(diào)整κ可以將數(shù)值耗散控制在合理水平。高超聲速流動中存在強激波,計算中對差商Δ-、Δ+采用min-mod限制器保證計算的魯棒性[27]。黏性項離散采用二階中心格式。非定常時間推進采用雙時間步方法[28],其中物理時間推進采用二階歐拉后差格式,隱式時間推進采用LU-SGS方法[29]。

    1.2算例驗證

    Duan等[30]開展了來流馬赫數(shù)Ma=2.99的超聲速湍流平板邊界層DNS研究,這里用此算例驗證所采用軟件和方法的有效性。

    算例來流密度ρ∞=0.089 1 kg/m3,來流溫度T∞=218.2 K。計算域尺寸為:流向69 mm,展向35 mm,法向110 mm。計算網(wǎng)格單元量約900萬,其中流向dx+=12,展向dy+=4.5,法向第一層網(wǎng)格dz+=0.3。壁面邊界條件為無滑移等溫壁,壁面溫度Tw=567.3 K。計算入口邊界層厚度為7.2 mm,入口邊界采用時間平均剖面疊加瞬時脈動的方法,其中時間平均剖面采用RANS計算的結(jié)果,瞬時脈動采用RRM(Rescaling-Recycling Method)[31-32]將下游x=52 mm位置處的瞬時脈動按照邊界層相似關(guān)系重新標定后疊加到入口邊界。

    物理模型采用基于改進延遲脫體渦模擬(Improved Delayed Detached-Eddy Simulation, IDDES)[33]方法的壁面?;鬁u模擬(Wall-Modeled LES, WMLES),其中湍流雷諾應(yīng)力模型采用一方程Spalart-Allmaras(S-A)模型[34],MUSCL重構(gòu)耗散參數(shù)κ=0.8。

    圖1為計算得到的邊界層平均速度剖面(時間平均加展向平均)。模擬結(jié)果與理論曲線以及文獻[30]DNS模擬結(jié)果符合較好,這表明本文所用的軟件和方法是可靠的。

    圖1 超聲速平板邊界層平均速度剖面計算結(jié)果Fig.1 Computational results of averaged velocity profile in supersonic flat boundary layer

    2 模型與網(wǎng)格說明

    計算采用的進氣道模型是在趙慧勇等[11]試驗?zāi)P突A(chǔ)上簡化得到的,如圖2所示。進氣道為分段壓縮結(jié)構(gòu),三個壓縮面均與水平方向有一定夾角。水平方向從前緣指向下游定義為x正方向,豎直方向定義為y向,渦流發(fā)生器(Vortex Generators, VG)凸起方向大致為y的負方向,展向定義為z方向。z向取了試驗?zāi)P椭芯€附近的一小段,x向長度約300 mm。渦流發(fā)生器采用了鉆石形和斜坡形兩種,后緣大致放在x=88 mm的位置(當?shù)剡吔鐚雍穸燃s為1.2 mm)。計算域展向尺寸以及模擬渦流發(fā)生器的個數(shù)和高度見表1。

    圖2 進氣道與渦流發(fā)生器模型及其網(wǎng)格(G1)Fig.2 Models and grids of inlet and vortex generators (G1)

    表1 計算模型

    Table 1 Computational models

    CaseVGshapeComputationalwidth/mmVGnumberVGheight/mm1Diamond6.021.02Ramp5.730.5

    考慮計算量,本文對鉆石形的渦流發(fā)生器只模擬了兩顆,對斜坡形的模擬了三顆。展向采用周期邊界條件。為了更好地分辨關(guān)鍵區(qū)域的流動結(jié)構(gòu),加密了渦流發(fā)生器附近和尾跡等區(qū)域的網(wǎng)格,并采用G1和G2粗細兩套網(wǎng)格計算考察了網(wǎng)格分辨率對計算結(jié)果的影響,這兩套網(wǎng)格在包含主要流動結(jié)構(gòu)的區(qū)域配置見表2。其中dx+為基于壁面單位的流向網(wǎng)格尺度,dz+為展向,dy+為法向,dy2+為法向第一層網(wǎng)格尺度。按照文獻[35]給出的LES網(wǎng)格分辨率準則,表2給出的粗網(wǎng)格G1已經(jīng)滿足了LES模擬要求。

    表2 計算網(wǎng)格

    3 數(shù)值計算和分析

    3.1計算條件

    根據(jù)趙慧勇等的試驗[11],選用馬赫數(shù)6試驗中的兩組試驗參數(shù)作為計算工況。來流馬赫數(shù)Ma=5.96,迎角α=1°,來流壓力p∞,來流溫度T∞以及第一、第二、第三壓縮面的物面溫度Tw1、Tw2、Tw3見表3。

    表3 計算條件

    3.2渦流發(fā)生器產(chǎn)生的渦結(jié)構(gòu)

    利用大渦模擬,對兩個工況進行模擬,獲取了非定常流場結(jié)構(gòu)。渦流發(fā)生器周圍的渦結(jié)構(gòu)如圖3 所示。其中Q值定義為

    (2)

    式中:Ω為基于特征速度為來流聲速和特征長度為1 mm的無量綱渦強;S為變形率張量。

    圖3 渦流發(fā)生器附近的流線和流場Q等值面Fig.3 Streamlines and Q iso-surfaces of flow around vortex generators

    根據(jù)工況1的瞬時流場可以看出,在鉆石形渦流發(fā)生器前緣根部產(chǎn)生了一個“馬蹄”渦,同時在后緣產(chǎn)生一對反向旋轉(zhuǎn)的流向渦。后緣產(chǎn)生的這對渦,其旋轉(zhuǎn)方向與外側(cè)的“馬蹄”渦反向,而位置比外側(cè)的“馬蹄”渦要遠離壁面一些,其強度也高于“馬蹄”渦。工況2的瞬時流場顯示,斜坡形的渦流發(fā)生器產(chǎn)生的渦結(jié)構(gòu)相對簡單,其主要結(jié)構(gòu)是渦流發(fā)生器側(cè)緣產(chǎn)生的一對反向旋轉(zhuǎn)的流向渦,這對渦的結(jié)構(gòu)、轉(zhuǎn)向與工況1中鉆石形渦流發(fā)生器后緣產(chǎn)生的流向渦相似。由于渦流發(fā)生器形狀、高度和分布密度不同,這里兩個工況中,鉆石形渦流發(fā)生器產(chǎn)生的流向渦更強,斜坡形渦流發(fā)生器產(chǎn)生的渦對更密集。

    3.3強制轉(zhuǎn)捩過程分析

    圖4顯示出鉆石形渦流發(fā)生器下游的主導(dǎo)結(jié)構(gòu)是尾緣產(chǎn)生的流向渦對。這對流向渦相互抬升,并在第二壓縮面上破碎并脫落出一系列的“發(fā)卡”渦,“發(fā)卡”渦很大程度上加劇了流動內(nèi)層和外層的流體交換,破壞了邊界層的穩(wěn)定性,促成轉(zhuǎn)捩。圖5顯示出斜坡形渦流發(fā)生器主要流向渦對結(jié)構(gòu)的發(fā)展,與鉆石形渦流發(fā)生器強制轉(zhuǎn)捩相似,斜坡形渦流發(fā)生器產(chǎn)生的流向渦對在第二壓縮面破碎脫落產(chǎn)生“發(fā)卡”渦,促成轉(zhuǎn)捩,明顯有別于工況1的是,工況2中流動結(jié)構(gòu)有展向的搖擺和彎曲,搖擺的頻率與“發(fā)卡”渦脫落的頻率是相關(guān)的。計算結(jié)果顯示,渦流發(fā)生器促使流動轉(zhuǎn)捩過程是其產(chǎn)生的流向渦對主導(dǎo)的,流向渦對結(jié)構(gòu)失穩(wěn)、破碎產(chǎn)生一系列“發(fā)卡”渦,使得邊界層流動轉(zhuǎn)捩為湍流。

    圖4 鉆石形渦流發(fā)生器強制轉(zhuǎn)捩瞬時流場Q等值面圖(G1, κ=0.9)Fig.4 Q iso-surfaces of forced transition instantaneous flow by diamond shaped vortex generators (G1, κ=0.9)

    圖5 斜坡形渦流發(fā)生器強制轉(zhuǎn)捩瞬時流場Q等值面圖(G2, κ=0.9)Fig.5 Q iso-surfaces of forced transition instantaneous flow by ramp shaped vortex generators (G2, κ=0.9)

    圖6 沿流向渦物理量脈動的均方根Fig.6 Root-mean-square of variable fluctuations along a stream-wise vortex

    兩種渦流發(fā)生器強制轉(zhuǎn)捩過程都與流向渦對結(jié)構(gòu)的失穩(wěn)有關(guān)。圖7為瞬時流場中壓力脈動p′沿流向渦對結(jié)構(gòu)兩個渦心的分布。圖7(a)顯示工況1瞬時流場中壓力脈動相位相同,表明工況1流向渦對失穩(wěn)由偶模式主導(dǎo),這種模式導(dǎo)致流向渦對對稱破碎;圖7(b)顯示工況2瞬時流場中壓力脈動相位相反,表明工況2流向渦對失穩(wěn)由奇模式主導(dǎo),這種模式導(dǎo)致流向渦對彎曲破碎。

    圖7 瞬時流場沿流向渦對兩個渦心的壓力脈動分布Fig.7 Pressure fluctuation of instantaneous field distribution along two vortex centers of counter-rotating vortices

    渦流發(fā)生器的幾何參數(shù)和分布密度會影響流向渦對結(jié)構(gòu)不穩(wěn)定增長的主導(dǎo)模式、頻率、破碎位置等,這些參數(shù)的具體影響有待深入研究。Li等[5-6]在研究G?rtler渦二次失穩(wěn)時指出,G?rtler渦波長(展向)較大時二次失穩(wěn)由偶模式主導(dǎo),波長較小時更容易出現(xiàn)奇模式。文中工況1渦對之間的距離較遠,對應(yīng)的波長為3 mm,工況2中渦對之間的距離較近,波長為1.9 mm,文中兩個工況在失穩(wěn)模式和波長的對應(yīng)關(guān)系上與Li等對于G?rtler渦二次失穩(wěn)的結(jié)論相符合。

    3.4計算收斂性討論

    對于本文研究的渦流發(fā)生器強制轉(zhuǎn)捩流動,流向渦對結(jié)構(gòu)破碎的位置基本決定了轉(zhuǎn)捩位置,準確模擬渦破碎位置取決于對擾動增長過程的準確模擬,而準確模擬擾動增長的關(guān)鍵則是避免數(shù)值耗散對其產(chǎn)生明顯的抑制。因此,網(wǎng)格密度和數(shù)值方法所引入耗散的大小,是準確計算轉(zhuǎn)捩位置的關(guān)鍵。本節(jié)通過對數(shù)值耗散和網(wǎng)格密度的調(diào)節(jié),考察了本文所采用的數(shù)值方法的收斂性。

    圖8和圖9給出了在G1和G2兩套網(wǎng)格上計算結(jié)果隨網(wǎng)格和耗散參數(shù)κ的變化,其中σp為壓力均方根。本文以流向渦對破碎位置作為計算收斂的判斷依據(jù)。渦破碎位置不僅受網(wǎng)格的影響,同樣受到數(shù)值格式耗散的影響。大體趨勢是網(wǎng)格越密、數(shù)值格式耗散越低,渦破碎位置越靠前,這是由于數(shù)值耗散能夠抑制擾動的增長。隨著格式耗散的降低和網(wǎng)格的加密,渦對破碎位置前移并趨于不變。如前所述,本文采用的粗網(wǎng)格G1滿足LES模擬的需求[35],如果一定的格式耗散(取決于κ值)下,渦破碎位置不再隨網(wǎng)格加密變化,則可以認為數(shù)值方法是合適的。

    從圖8和圖9中流向渦對破碎位置以及流場中脈動壓力的增長看,在格式耗散較小(κ較大)的情況下,計算更容易達到網(wǎng)格收斂;κ=0.9時,G1和G2兩套網(wǎng)格下計算得到的渦破碎位置幾乎一致,可認為實現(xiàn)了網(wǎng)格收斂;對于文中兩個工況,在G1網(wǎng)格下,采用κ=0.9計算基本得到了收斂的流向渦對破碎位置。

    圖8 第二壓縮面上流場Q等值面圖 Fig.8 Q iso-surfaces of flow on the second compress surface

    圖9 沿流向渦壓力脈動的均方根Fig.9 Root-mean-square of pressure fluctuations along a stream-wise vortex

    4 結(jié) 論

    1) 渦流發(fā)生器強制轉(zhuǎn)捩過程的物理圖像可以描述為:渦流發(fā)生器在其下游產(chǎn)生一對流向渦,這對流向渦失穩(wěn)產(chǎn)生一系列“發(fā)卡”渦并逐漸演化直至破碎,最終使邊界層轉(zhuǎn)捩為湍流。

    2) 渦流發(fā)生器下游流向渦的失穩(wěn)包括兩種模式:偶模式和奇模式,前者產(chǎn)生對稱形式的渦破碎,后者產(chǎn)生彎曲形式的渦對破碎。本文開展的數(shù)值模擬研究同時捕捉到了這兩種失穩(wěn)模式。

    3) 采用隱式大渦模擬方法計算強制轉(zhuǎn)捩問題時,數(shù)值耗散對計算結(jié)果有明顯影響,文中對網(wǎng)格和數(shù)值耗散影響作了測試:在較小的格式耗散情況下,計算更容易達到網(wǎng)格收斂;隨著網(wǎng)格加密、耗散變小,計算得到的流向渦對破碎前移并趨于不變。

    [1]BERRY S A, DARYABEIGI K, WURSTER K. Boundary-layer transition on X-43A[J]. Journal of Spacecraft and Rockets, 2010, 47(6): 922-934.

    [2]KLEBANOFF P S, TIDSTROM K D. Mechanism by which a two-dimensional roughness element induces boundary-layer transition[J]. Physics of Fluids, 1972, 15(7): 1173-1188.

    [3]RESHOTKO E, TUMIN A. Role of transient growth in roughness-induced transition[J]. AIAA Journal, 2004, 42(4): 766-770.

    [4]SESCU A, PENDYALAY R K, THOMPSON D. On the growth of G?rtler vortices excited by distributed roughness elements: AIAA-2014-2885[R]. Reston: AIAA, 2014.

    [5]LI F, MALIK M R. Fundamental and subharmonic secondary instabilities of G?rtler vortices[J]. Journal of Fluid Mechanics, 1995, 297: 77-100.

    [6]LI F, CHOUDHARI M, CHANG C L, et al. Development and breakdown of G?rtler vortices in high speed boundary layers: AIAA-2010-705[R]. Reston: AIAA, 2010.

    [7]TULLIO N D, PAREDES P, SANDHAM N D, et al. Laminar-turbulent transition induced by a discrete roughness element in a supersonic boundary layer[J]. Journal of Fluid Mechanics, 2013, 735: 613-646.

    [8]BERRY S A, HORVATH T J. Discrete-roughness transition for hypersonic flight vehicles[J]. Journal of Spacecraft and Rockets, 2008, 45(2): 216-227.

    [9]BERRY S A, AUSLENDER A H, DILLEY A D, et al. Hypersonic boundary layer trip development for Hyper-X[J]. Journal of Spacecraft and Rockets, 2001, 38(6): 853-864.

    [10]BERRY S A, NOWAK R J, HORVATH T J. Boundary layer control for hypersonic airbreathing vehicles: AIAA-2004-2246[R]. Reston: AIAA, 2004.

    [11]趙慧勇, 周瑜, 倪鴻禮, 等. 高超聲速進氣道邊界層強迫轉(zhuǎn)捩試驗[J]. 實驗流體力學, 2012, 26(1): 1-6.

    ZHAO H Y, ZHOU Y, NI H L, et al. Test of forced boundary-layer transition on hypersonic inlet[J]. Journal of Experiments in Fluid Mechanics, 2012, 26(1): 1-6 (in Chinese).

    [12]DANEHY P M, BATHEL B, IVEY C, et al. NO PLIF study of hypersonic transition over a discrete hemispherical roughness element: AIAA-2009-0394[R]. Reston: AIAA, 2009.

    [13]岡敦殿, 易仕和, 陳植, 等. 超聲速平板上直立圓臺突起物繞流流場研究[C]//中國力學大會-2013. 北京: 中國力學學會, 2013.

    GANG D D, YI S H, CHEN Z, et al. Flow field investigations of supersonic flow over a circular protuberance mounted on a flat plate[C]//CCTAM-2013. Beijing: The Chinese Society of Theoretical and Applied Mechanics, 2013 (in Chinese).

    [14]張慶虎. 超聲速流動分離及其控制的試驗研究[D]. 長沙: 國防科學技術(shù)大學, 2013.

    ZHANG Q H. Experimental investigation of supersonic flow separation and its micro-ramp control[D]. Changsha: National University of Defense Technology, 2013 (in Chinese).

    [15]BORG M P, SCHNEIDER S P. Effect of freestream noise on roughness-induced transition for the X-51A forebody[J]. Journal of Spacecraft and Rockets, 2008, 45(6): 1106-1116.

    [16]鄧小兵. 來流噪聲對高超聲速進氣道強制轉(zhuǎn)捩的影響[C]//中國力學大會-2013. 北京: 中國力學學會, 2013.

    DENG X B. Freestream noise impact on forced transition of hypersonic boundary layer[C]//CCTAM-2013. Beijing: The Chinese Society of Theoretical and Applied Mechanics, 2013 (in Chinese).

    [17]MARXEN O, IACCARINO G, SHAQFEH E S G. Numerical simulations of hypersonic boundary-layer instability with localized roughness: AIAA-2011-0567[R]. Reston: AIAA, 2011.

    [18]IYER P S, MUPPIDI S, MAHESH K. Roughness-induced transition in high speed flows: AIAA-2011-0566[R]. Reston: AIAA, 2011.

    [19]DUAN Z W, XIAO Z X, FU S. Direct numerical simulation of hypersonic transition induced by ramp roughness elements: AIAA-2014-2037[R]. Reston: AIAA, 2014.

    [20]DUAN Z W, XIAO Z X. Direct numerical simulation of geometrical parameter effects on the hypersonic ramp-induced transition: AIAA-2014-2495[R]. Reston: AIAA, 2014.

    [21]涂國華, 燕振國, 趙曉慧, 等. SA和SST湍流模型對高超聲速邊界層強制轉(zhuǎn)捩的適應(yīng)性[J]. 航空學報, 2015, 36(5): 1471-1479.

    TU G H, YAN Z G, ZHAO X H, et al. SA and SST turbulence models for hypersonic forced boundary layer transition[J]. Acta Aeronautica et Astronautia Sinica, 2015, 36(5): 1471-1479 (in Chinese).

    [22]周玲, 閆超, 孔維萱. 高超聲速飛行器前體邊界層強制轉(zhuǎn)捩數(shù)值模擬[J]. 航空學報, 2014, 35(6): 1487-1495.

    ZHOU L, YAN C, KONG W X. Numerical simulation of forced boundary layer transition on hypersonic vehicle forebody[J]. Acta Aeronautica et Astronautica Sinica, 2014, 35(6): 1487-1495 (in Chinese).

    [23]周玲, 閆超, 郝子輝, 等. 轉(zhuǎn)捩模式與轉(zhuǎn)捩準則預(yù)測高超聲速邊界層流動[J]. 航空學報, 2016, 37(4): 1092-1102.

    ZHOU L, YAN C, HAO Z H, et al. Transition model and transition criteria for hypersonic boundary layer flow[J]. Acta Aeronautica et Astronautica Sinica, 2016, 37(4): 1092-1102 (in Chinese).

    [24]SAYADI T, MOIN P. Large eddy simulation of controlled transition to turbulence[J]. Physics of Fluids, 2012, 24(11): 923-938.

    [25]ORLIK E, FEDIOUN I. Hypersonic boundary-layer transition forced by wall injection: A numerical study[J]. Journal of Spacecraft and Rockets, 2014, 51(4): 1306-1318.

    [26]鄧小兵, 金玲. 讓誤差飛一會兒——湍流大渦模擬中的動態(tài)自適應(yīng)迎風方法[C]//中國力學大會-2013. 北京: 中國力學學會, 2013.

    DENG X B, JIN L. Let the error fly for a while—A dynamic adaptive upwind method for large eddy simulations of turbulent flow[C]//CCTAM-2013. Beijing: The Chinese Society of Theoretical and Applied Mechanics, 2013 (in Chinese).

    [27]朱自強, 吳子牛, 李津, 等. 應(yīng)用計算流體力學[M]. 北京: 北京航空航天大學出版社, 1998: 54-88.

    ZHU Z Q, WU Z N, LI J, et al. Application of computational fluid dynamics[M]. Beijing: Beihang University Press, 1998: 54-88 (in Chinese).

    [28]JAMESON A. Time dependent calculations using multigrid with application to unsteady flows past airfoils and wings: AIAA-1991-1596[R]. Reston: AIAA, 1991.

    [29]YOON S, JAMESON A. Lower-upper symmetric-Gauss-Seidel method for the Euler and Navier-Stokes equations[J]. AIAA Journal, 1988, 26(9): 1025-1026.

    [30]DUAN L, BEEKMAN I, MARTN M P. Direct numerical simulation of hypersonic turbulent boundary layers with varying freestream Mach number: AIAA-2010-0353[R]. Reston: AIAA, 2010.

    [31]LUND T S, WU X, SQUIRES K D. Generation of turbulent inflow data for spatially-developing boundary layer simulations[J]. Journal of Computational Physics, 1998, 140(2): 233-258.

    [32]EDWARDS J R, CHOI J, BOLES J A. Hybrid LES/RANS simulation of a mach 5 compression-corner interaction: AIAA-2008-0718[R]. Reston: AIAA, 2008.

    [33]SHUR M L, SPALART P R, STRELETS M K, et al. A hybrid RANS-LES approach with delayed-DES and wall-modelled LES capabilities[J]. International Journal of Heat and Fluid Flow, 2008, 29(6): 1638-1649

    [34]SPALART P R, ALLMARAS S R. A one-equation turbulence model for aerodynamic flows: AIAA-1992-0439[R]. Reston: AIAA, 1992.

    [35]GEORGIADIS N J, RIZZETTA D P, FUREBY C. Large-eddy simulation: current capabilities, recommended practices, and future research[J]. AIAA Journal, 2010, 48(8): 1772-1784.

    趙曉慧男, 博士研究生, 助理研究員。主要研究方向: 湍流數(shù)值模擬, 高超聲速空氣動力學。

    Tel.: 0816-2463219

    E-mail: xhzhao@skla.cardc.cn

    鄧小兵男, 博士, 副研究員。主要研究方向: 計算流體力學, 湍流數(shù)值模擬。

    Tel.: 0816-2463219

    E-mail: xbdeng@skla.cardc.cn

    毛枚良男, 博士, 研究員, 博士生導(dǎo)師。主要研究方向: 計算流體力學, 高階精度格式, 高超聲速空氣動力學。

    Tel.: 0816-2463296

    E-mail: mm1219@163.com

    Large eddy simulation for forced transition flow at hypersonic inlet

    ZHAO Xiaohui1, DENG Xiaobing1,*, MAO Meiliang2, YANG Wei1, ZHAO Huiyong3

    1. State Key Laboratory of Aerodynamics, China Aerodynamics Research and Development Center,Mianyang621000, China 2. Computational Aerodynamics Institute, China Aerodynamics Research and Development Center,Mianyang621000, China 3. Science and Technology on Scramjet Laboratory, China Aerodynamics Research and Development Center,Mianyang 621000, China

    Forced transition is commonly used for a hypersonic air breathing vehicle to ensure its scramjet’s normal start by placing roughness elements or vortex generators in the inlet boundary layer. To get a clear image of forced transition process, an implicit large eddy simulation method is used for laminar-turbulent forced transition flows in the boundary layer of a hypersonic inlet configuration. Main structures are obtained for the forced transition flows induced by the vortex generators shaped of diamonds and ramps, which show that the flows are dominated by stream-wise counter-rotating vortices. Fluctuations grow in the stream-wise direction and cause instabilities. Two fundamental modes of the instabilities are found in the simulations, the even mode and the odd one. The even mode leads to a varicose way of breaking down of the counter-rotating vortices, while the odd one leads to a sinuous way. Strings of hairpin vortices are generated after the breaking down and finally cause the transition. A discussion is carried out on the computation convergence together with the influence of the grids and numerical dissipation.

    hypersonic; inlet; forced transition; vortex generator; large eddy simulation

    2016-02-16; Revised: 2016-02-17; Accepted: 2016-06-21; Published online: 2016-06-2715:34

    National Natural Science Foundation of China (11402289)

    . Tel.: 0816-2463219E-mail: xbdeng@skla.cardc.cn

    2016-02-16; 退修日期: 2016-02-17; 錄用日期: 2016-06-21;

    時間: 2016-06-2715:34

    www.cnki.net/kcms/detail/11.1929.V.20160627.1534.008.html

    國家自然科學基金 (11402289)

    .Tel.: 0816-2463219E-mail: xbdeng@skla.cardc.cn

    10.7527/S1000-6893.2016.0200

    V211.3; O355

    A

    1000-6893(2016)08-2445-09

    引用格式: 趙曉慧, 鄧小兵, 毛枚良, 等. 高超聲速進氣道強制轉(zhuǎn)捩流動的大渦模擬[J]. 航空學報, 2016, 37(8): 2445-2453. ZHAO X H, DENG X B, MAO M L, et al. Large eddy simulation for forced transition flow at hypersonic inlet[J]. Acta Aeronautica et Astronautica Sinica, 2016, 37(8): 2445-2453.

    http://hkxb.buaa.edu.cnhkxb@buaa.edu.cn

    URL: www.cnki.net/kcms/detail/11.1929.V.20160627.1534.008.html

    猜你喜歡
    大渦進氣道邊界層
    基于AVL-Fire的某1.5L發(fā)動機進氣道優(yōu)化設(shè)計
    基于輔助進氣門的進氣道/發(fā)動機一體化控制
    基于HIFiRE-2超燃發(fā)動機內(nèi)流道的激波邊界層干擾分析
    基于壁面射流的下?lián)舯┝鞣欠€(wěn)態(tài)風場大渦模擬
    軸流風機葉尖泄漏流動的大渦模擬
    基于大渦模擬的旋風分離器錐體結(jié)構(gòu)影響研究
    一類具有邊界層性質(zhì)的二次奇攝動邊值問題
    The coupling characteristics of supersonic dual inlets for missile①
    非特征邊界的MHD方程的邊界層
    某柴油機進氣道數(shù)值模擬及試驗研究
    汽車零部件(2014年2期)2014-03-11 17:46:30
    精品国产超薄肉色丝袜足j| 老熟妇仑乱视频hdxx| 特大巨黑吊av在线直播 | 国产亚洲欧美98| 亚洲精品久久成人aⅴ小说| ponron亚洲| 久久婷婷人人爽人人干人人爱| 免费女性裸体啪啪无遮挡网站| 日本一区二区免费在线视频| www.熟女人妻精品国产| 亚洲午夜精品一区,二区,三区| 欧美日本亚洲视频在线播放| 久久精品夜夜夜夜夜久久蜜豆 | 亚洲av成人不卡在线观看播放网| 精品久久久久久,| 国产真人三级小视频在线观看| 美女免费视频网站| 久久国产精品男人的天堂亚洲| 成在线人永久免费视频| 成年免费大片在线观看| 亚洲欧美日韩高清在线视频| 变态另类成人亚洲欧美熟女| 男女做爰动态图高潮gif福利片| 少妇裸体淫交视频免费看高清 | 在线观看免费午夜福利视频| 淫妇啪啪啪对白视频| 日韩精品中文字幕看吧| av在线播放免费不卡| 美女高潮到喷水免费观看| 欧美日韩乱码在线| 日日干狠狠操夜夜爽| 亚洲黑人精品在线| 中亚洲国语对白在线视频| 国产伦人伦偷精品视频| 午夜两性在线视频| 欧美日韩黄片免| 一进一出抽搐gif免费好疼| 国产精品98久久久久久宅男小说| 国产精品1区2区在线观看.| 亚洲欧美日韩无卡精品| а√天堂www在线а√下载| 老司机午夜十八禁免费视频| 国产精品野战在线观看| 亚洲精品中文字幕在线视频| 99在线人妻在线中文字幕| 嫁个100分男人电影在线观看| 一卡2卡三卡四卡精品乱码亚洲| 好看av亚洲va欧美ⅴa在| 999精品在线视频| 丰满的人妻完整版| 免费看十八禁软件| 两个人看的免费小视频| 丝袜在线中文字幕| 老汉色av国产亚洲站长工具| 久久久久国内视频| 欧美另类亚洲清纯唯美| 国产精品国产高清国产av| 国产成人系列免费观看| 黄频高清免费视频| 亚洲男人的天堂狠狠| 又黄又爽又免费观看的视频| 精品国产一区二区三区四区第35| 这个男人来自地球电影免费观看| 怎么达到女性高潮| 不卡一级毛片| 男女做爰动态图高潮gif福利片| 天天添夜夜摸| 一夜夜www| 两个人视频免费观看高清| 免费av毛片视频| 一区二区日韩欧美中文字幕| 91麻豆精品激情在线观看国产| 亚洲欧美一区二区三区黑人| 欧美中文综合在线视频| 精品少妇一区二区三区视频日本电影| 日韩欧美一区二区三区在线观看| 中亚洲国语对白在线视频| 在线播放国产精品三级| 久热爱精品视频在线9| 在线永久观看黄色视频| 男人操女人黄网站| 一区二区三区国产精品乱码| 日本五十路高清| 90打野战视频偷拍视频| 中出人妻视频一区二区| 婷婷亚洲欧美| 一区二区日韩欧美中文字幕| 久久性视频一级片| 国产精品野战在线观看| 欧美又色又爽又黄视频| 99国产综合亚洲精品| 女人被狂操c到高潮| 搞女人的毛片| 宅男免费午夜| 麻豆成人午夜福利视频| 在线观看舔阴道视频| 欧美日本亚洲视频在线播放| 母亲3免费完整高清在线观看| 伦理电影免费视频| 精品国产超薄肉色丝袜足j| 亚洲精品在线观看二区| www.999成人在线观看| 一级片免费观看大全| 日本一区二区免费在线视频| 岛国在线观看网站| 精品人妻1区二区| 香蕉丝袜av| 色老头精品视频在线观看| 少妇裸体淫交视频免费看高清 | 啦啦啦 在线观看视频| 亚洲五月天丁香| 正在播放国产对白刺激| av免费在线观看网站| 97碰自拍视频| 欧美激情极品国产一区二区三区| 精品无人区乱码1区二区| 欧美乱色亚洲激情| 亚洲国产毛片av蜜桃av| 桃红色精品国产亚洲av| 欧美黄色淫秽网站| 成人特级黄色片久久久久久久| 97人妻精品一区二区三区麻豆 | 午夜精品在线福利| 国产v大片淫在线免费观看| 国产精品国产高清国产av| 国产午夜福利久久久久久| 欧美av亚洲av综合av国产av| 国产精品98久久久久久宅男小说| 国产精品日韩av在线免费观看| 视频区欧美日本亚洲| 亚洲人成电影免费在线| 久久热在线av| 国产亚洲精品第一综合不卡| 亚洲精品美女久久久久99蜜臀| 中文字幕人妻熟女乱码| 天天一区二区日本电影三级| 18禁黄网站禁片免费观看直播| 色精品久久人妻99蜜桃| 中文字幕人妻熟女乱码| 亚洲色图av天堂| 亚洲真实伦在线观看| av在线天堂中文字幕| 国产伦人伦偷精品视频| 香蕉丝袜av| 欧美性猛交╳xxx乱大交人| 一级黄色大片毛片| 不卡av一区二区三区| 亚洲色图av天堂| 免费在线观看黄色视频的| 黑丝袜美女国产一区| 久久中文字幕人妻熟女| a级毛片a级免费在线| 免费在线观看亚洲国产| 亚洲第一av免费看| 欧美日本亚洲视频在线播放| 2021天堂中文幕一二区在线观 | 这个男人来自地球电影免费观看| 日本三级黄在线观看| 窝窝影院91人妻| 国产野战对白在线观看| 大香蕉久久成人网| 一本大道久久a久久精品| 国产精品影院久久| 99riav亚洲国产免费| 国产人伦9x9x在线观看| 亚洲真实伦在线观看| 欧美性猛交黑人性爽| 亚洲国产精品成人综合色| 最近最新中文字幕大全电影3 | 国产伦在线观看视频一区| 夜夜看夜夜爽夜夜摸| 亚洲美女黄片视频| 午夜免费激情av| 欧美国产日韩亚洲一区| 男人舔奶头视频| 国产一区二区在线av高清观看| 久久草成人影院| 午夜日韩欧美国产| 国产高清激情床上av| 男女之事视频高清在线观看| 桃色一区二区三区在线观看| 午夜激情av网站| 欧美+亚洲+日韩+国产| 日韩欧美一区二区三区在线观看| 欧美日本视频| 国产极品粉嫩免费观看在线| 午夜福利成人在线免费观看| 国产激情偷乱视频一区二区| 亚洲av中文字字幕乱码综合 | 99国产精品99久久久久| 免费在线观看视频国产中文字幕亚洲| 天天躁狠狠躁夜夜躁狠狠躁| 日韩av在线大香蕉| 久久 成人 亚洲| 欧美亚洲日本最大视频资源| 俄罗斯特黄特色一大片| 校园春色视频在线观看| 精品电影一区二区在线| а√天堂www在线а√下载| 精品久久久久久久久久免费视频| 成人av一区二区三区在线看| 欧美黑人巨大hd| 国产激情久久老熟女| 亚洲av五月六月丁香网| 久久亚洲精品不卡| 久久精品91无色码中文字幕| 亚洲av片天天在线观看| 日韩成人在线观看一区二区三区| 欧美zozozo另类| 日日爽夜夜爽网站| 亚洲真实伦在线观看| www.www免费av| 在线播放国产精品三级| 亚洲美女黄片视频| 欧美黑人巨大hd| 男男h啪啪无遮挡| 欧美+亚洲+日韩+国产| 午夜福利在线观看吧| 日本在线视频免费播放| 国产午夜福利久久久久久| 亚洲国产欧洲综合997久久, | 999久久久国产精品视频| 18禁美女被吸乳视频| 亚洲avbb在线观看| 久久久久亚洲av毛片大全| 久久中文字幕人妻熟女| 国产97色在线日韩免费| 国产精品一区二区精品视频观看| 国产精品自产拍在线观看55亚洲| 亚洲avbb在线观看| 又紧又爽又黄一区二区| 日本a在线网址| 精品日产1卡2卡| 老熟妇仑乱视频hdxx| 淫秽高清视频在线观看| 我的亚洲天堂| 亚洲美女黄片视频| 久9热在线精品视频| 精品久久久久久久人妻蜜臀av| av欧美777| 一级片免费观看大全| netflix在线观看网站| 亚洲精品一区av在线观看| 国产av不卡久久| 人妻久久中文字幕网| 熟女电影av网| 亚洲aⅴ乱码一区二区在线播放 | 天堂影院成人在线观看| 精品不卡国产一区二区三区| 黄色视频,在线免费观看| 精品一区二区三区四区五区乱码| 亚洲一区二区三区色噜噜| 亚洲av第一区精品v没综合| 精品久久久久久成人av| 一边摸一边做爽爽视频免费| 欧美zozozo另类| 国产片内射在线| 最好的美女福利视频网| 欧美黑人巨大hd| 人人妻人人看人人澡| 久久久国产欧美日韩av| 国产麻豆成人av免费视频| 美国免费a级毛片| 欧美黑人欧美精品刺激| 日本撒尿小便嘘嘘汇集6| 日本黄色视频三级网站网址| 欧美性猛交╳xxx乱大交人| 夜夜躁狠狠躁天天躁| ponron亚洲| 亚洲狠狠婷婷综合久久图片| 欧美色欧美亚洲另类二区| 色哟哟哟哟哟哟| 亚洲黑人精品在线| 又黄又爽又免费观看的视频| 久久中文字幕一级| 久久久久久久久中文| 国产野战对白在线观看| 在线看三级毛片| 久久天堂一区二区三区四区| 日本三级黄在线观看| 一边摸一边做爽爽视频免费| 麻豆av在线久日| 天堂影院成人在线观看| 少妇被粗大的猛进出69影院| av福利片在线| 国产一级毛片七仙女欲春2 | 黑人操中国人逼视频| 亚洲国产精品sss在线观看| 男女之事视频高清在线观看| 亚洲天堂国产精品一区在线| 韩国精品一区二区三区| 麻豆成人午夜福利视频| 国产精品亚洲美女久久久| 免费av毛片视频| 日韩欧美免费精品| 99久久99久久久精品蜜桃| 男女床上黄色一级片免费看| www国产在线视频色| 18禁观看日本| 精品免费久久久久久久清纯| 欧美久久黑人一区二区| 搡老熟女国产l中国老女人| 国内精品久久久久精免费| 亚洲国产精品成人综合色| 日韩大码丰满熟妇| 动漫黄色视频在线观看| 成人特级黄色片久久久久久久| 亚洲av成人av| 欧美中文日本在线观看视频| 女人高潮潮喷娇喘18禁视频| 亚洲五月婷婷丁香| 欧美精品亚洲一区二区| 男女床上黄色一级片免费看| 一卡2卡三卡四卡精品乱码亚洲| 动漫黄色视频在线观看| 精品免费久久久久久久清纯| 日韩大尺度精品在线看网址| 亚洲中文字幕日韩| 超碰成人久久| 成人18禁高潮啪啪吃奶动态图| 99精品在免费线老司机午夜| 久久久久国产精品人妻aⅴ院| 嫩草影视91久久| 精品国产美女av久久久久小说| 亚洲国产欧洲综合997久久, | 好男人在线观看高清免费视频 | 91九色精品人成在线观看| 国产乱人伦免费视频| 黄网站色视频无遮挡免费观看| 国产又爽黄色视频| 中国美女看黄片| 狂野欧美激情性xxxx| 亚洲精品av麻豆狂野| 亚洲中文字幕日韩| 精品国产一区二区三区四区第35| 久久久国产精品麻豆| 在线观看午夜福利视频| 成人国语在线视频| 怎么达到女性高潮| 亚洲最大成人中文| 亚洲精品一区av在线观看| 哪里可以看免费的av片| 香蕉国产在线看| 天堂影院成人在线观看| 亚洲一区中文字幕在线| 精品人妻1区二区| 欧美黑人精品巨大| 免费高清在线观看日韩| 天天添夜夜摸| 亚洲av电影不卡..在线观看| 人人澡人人妻人| 亚洲av成人不卡在线观看播放网| 88av欧美| 黄片小视频在线播放| 色播亚洲综合网| 久9热在线精品视频| 欧美色欧美亚洲另类二区| 亚洲欧洲精品一区二区精品久久久| 亚洲av片天天在线观看| 一区二区三区高清视频在线| 国产高清激情床上av| 免费人成视频x8x8入口观看| 午夜久久久在线观看| 亚洲精品中文字幕在线视频| 十八禁网站免费在线| 国产亚洲欧美98| 久久精品亚洲精品国产色婷小说| 亚洲一区二区三区色噜噜| 一级a爱片免费观看的视频| 丰满的人妻完整版| 日本 欧美在线| 亚洲欧美精品综合久久99| 久久青草综合色| 超碰成人久久| 无限看片的www在线观看| 国语自产精品视频在线第100页| 亚洲精品国产区一区二| 黄色女人牲交| 亚洲欧美一区二区三区黑人| 特大巨黑吊av在线直播 | 精品国产超薄肉色丝袜足j| 人人妻人人澡欧美一区二区| 变态另类丝袜制服| 老司机午夜福利在线观看视频| 美女国产高潮福利片在线看| 在线免费观看的www视频| 老司机午夜十八禁免费视频| 伦理电影免费视频| 嫩草影视91久久| 久久久国产精品麻豆| 国产亚洲av嫩草精品影院| 性欧美人与动物交配| 免费一级毛片在线播放高清视频| 亚洲中文日韩欧美视频| 手机成人av网站| 国产爱豆传媒在线观看 | 后天国语完整版免费观看| 国产成人一区二区三区免费视频网站| 男男h啪啪无遮挡| 欧美不卡视频在线免费观看 | cao死你这个sao货| 特大巨黑吊av在线直播 | 亚洲国产精品合色在线| 国产成人欧美| 亚洲自拍偷在线| 久久久久久久午夜电影| 天堂影院成人在线观看| 一区二区三区高清视频在线| 国产精品av久久久久免费| 欧美国产日韩亚洲一区| 此物有八面人人有两片| 久久国产精品男人的天堂亚洲| 精品久久久久久久末码| 亚洲国产中文字幕在线视频| 一级黄色大片毛片| 一进一出抽搐gif免费好疼| 日本免费a在线| 一区福利在线观看| 免费在线观看完整版高清| 中文字幕人妻丝袜一区二区| 亚洲av成人不卡在线观看播放网| 久久久久久久久久黄片| 黄片大片在线免费观看| 久久午夜综合久久蜜桃| 欧美日韩精品网址| 亚洲专区中文字幕在线| 黄色视频,在线免费观看| 三级毛片av免费| 午夜激情av网站| 国产不卡一卡二| cao死你这个sao货| aaaaa片日本免费| 成人欧美大片| 亚洲一区二区三区色噜噜| 欧美性猛交黑人性爽| 成人手机av| 中文资源天堂在线| 97碰自拍视频| 天天添夜夜摸| 观看免费一级毛片| 色哟哟哟哟哟哟| 亚洲精品国产精品久久久不卡| 国产熟女午夜一区二区三区| 精品久久久久久久久久免费视频| 变态另类成人亚洲欧美熟女| 国产av一区二区精品久久| 老鸭窝网址在线观看| 日本三级黄在线观看| 在线国产一区二区在线| 精品国产超薄肉色丝袜足j| 国产不卡一卡二| 久久久久久久久久黄片| 18美女黄网站色大片免费观看| 黄色视频,在线免费观看| 香蕉丝袜av| 国产区一区二久久| 国产精品野战在线观看| 亚洲国产高清在线一区二区三 | 久久草成人影院| 亚洲av电影不卡..在线观看| 成人国产综合亚洲| 亚洲激情在线av| 可以免费在线观看a视频的电影网站| 麻豆一二三区av精品| 久久精品91蜜桃| 俄罗斯特黄特色一大片| 青草久久国产| 国产片内射在线| 亚洲av五月六月丁香网| av欧美777| 久久亚洲真实| 成人三级做爰电影| 两性午夜刺激爽爽歪歪视频在线观看 | 欧美成人午夜精品| 国产精品爽爽va在线观看网站 | www.www免费av| 91老司机精品| 久久婷婷成人综合色麻豆| 久久久久久久午夜电影| 国产单亲对白刺激| 国产精品99久久99久久久不卡| 亚洲,欧美精品.| 在线观看日韩欧美| 欧美乱码精品一区二区三区| 久久精品国产亚洲av香蕉五月| 高清毛片免费观看视频网站| 人人妻人人看人人澡| 一区二区日韩欧美中文字幕| 国产熟女午夜一区二区三区| 女人高潮潮喷娇喘18禁视频| 免费在线观看成人毛片| 日本 av在线| 久久人妻av系列| 久久精品夜夜夜夜夜久久蜜豆 | 欧美成狂野欧美在线观看| 精品国产国语对白av| 视频区欧美日本亚洲| 99国产精品99久久久久| 他把我摸到了高潮在线观看| 国产成人精品久久二区二区免费| 欧美zozozo另类| 一本综合久久免费| 桃红色精品国产亚洲av| av电影中文网址| 中国美女看黄片| 欧美一区二区精品小视频在线| 1024手机看黄色片| 美女免费视频网站| 婷婷丁香在线五月| 黄片小视频在线播放| 在线播放国产精品三级| 久久久久久大精品| 亚洲av成人一区二区三| 1024视频免费在线观看| 黑人巨大精品欧美一区二区mp4| bbb黄色大片| 18禁观看日本| 人人妻人人澡欧美一区二区| 日韩一卡2卡3卡4卡2021年| 在线观看一区二区三区| 欧美成人免费av一区二区三区| 中文字幕人妻丝袜一区二区| 精品第一国产精品| 日韩欧美免费精品| 久久久久久免费高清国产稀缺| 久久精品国产亚洲av高清一级| 不卡一级毛片| 欧美日韩一级在线毛片| 久久欧美精品欧美久久欧美| www.精华液| 欧美黑人精品巨大| 亚洲中文日韩欧美视频| АⅤ资源中文在线天堂| 美女国产高潮福利片在线看| 好男人在线观看高清免费视频 | 亚洲精品美女久久久久99蜜臀| 好看av亚洲va欧美ⅴa在| 叶爱在线成人免费视频播放| 日本黄色视频三级网站网址| 国产av又大| 12—13女人毛片做爰片一| 色av中文字幕| 哪里可以看免费的av片| 手机成人av网站| 午夜精品在线福利| 成人三级黄色视频| 人人澡人人妻人| 黄频高清免费视频| 丝袜在线中文字幕| 欧美另类亚洲清纯唯美| 一进一出抽搐gif免费好疼| 一卡2卡三卡四卡精品乱码亚洲| 国产精品久久久av美女十八| 亚洲男人的天堂狠狠| 亚洲国产高清在线一区二区三 | 欧美精品啪啪一区二区三区| 欧美激情 高清一区二区三区| 在线av久久热| 99国产极品粉嫩在线观看| 成人手机av| 999久久久精品免费观看国产| 精品免费久久久久久久清纯| 中文字幕av电影在线播放| 色综合站精品国产| 淫妇啪啪啪对白视频| 日韩精品免费视频一区二区三区| 精品欧美国产一区二区三| 精品电影一区二区在线| www.www免费av| 波多野结衣巨乳人妻| а√天堂www在线а√下载| 无限看片的www在线观看| 欧美色视频一区免费| 久久精品夜夜夜夜夜久久蜜豆 | 国产日本99.免费观看| 两个人免费观看高清视频| 久久久久国产一级毛片高清牌| 一二三四社区在线视频社区8| 人妻丰满熟妇av一区二区三区| 日日干狠狠操夜夜爽| 欧美精品亚洲一区二区| 可以在线观看的亚洲视频| 免费av毛片视频| 我的亚洲天堂| 一个人观看的视频www高清免费观看 | 不卡一级毛片| 18禁黄网站禁片午夜丰满| 久久精品国产清高在天天线| 美国免费a级毛片| 黑丝袜美女国产一区| 亚洲精品在线美女| 精品欧美一区二区三区在线| 色在线成人网| 中文字幕av电影在线播放| 久久久久国内视频| 久久久久久免费高清国产稀缺| 欧美另类亚洲清纯唯美| 国产欧美日韩一区二区三| 亚洲av片天天在线观看| 一区二区三区精品91| 亚洲电影在线观看av| 欧美另类亚洲清纯唯美| 宅男免费午夜| 又黄又爽又免费观看的视频| 国产国语露脸激情在线看| 一级a爱片免费观看的视频| 99国产精品99久久久久| 亚洲狠狠婷婷综合久久图片| 999精品在线视频| 亚洲国产中文字幕在线视频| 亚洲精品美女久久av网站| 国产乱人伦免费视频| 亚洲精品久久国产高清桃花| 亚洲第一av免费看| 亚洲av成人一区二区三|