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

    激波/湍流邊界層干擾壓力脈動特性數(shù)值研究1)

    2021-11-09 06:26:08童福林段俊亦周桂宇李新亮
    力學(xué)學(xué)報 2021年7期
    關(guān)鍵詞:物面外層邊界層

    童福林 段俊亦 周桂宇 李新亮

    * (中國空氣動力研究與發(fā)展中心空氣動力學(xué)國家重點實驗室,四川綿陽 621000)

    ? (中國科學(xué)院力學(xué)研究所高溫氣體動力學(xué)國家重點實驗室,北京 100190)

    ** (中國空氣動力研究與發(fā)展中心計算空氣動力研究所,四川綿陽 621000)

    ?? (中國科學(xué)院大學(xué)工程科學(xué)學(xué)院,北京 100049)

    引言

    激波與湍流邊界層的相互作用問題廣泛存在于各類高速飛行器中.已有研究表明[1-2],在強激波作用下,干擾區(qū)內(nèi)壓力脈動峰值急劇升高,分離激波的非定常運動導(dǎo)致物面脈動壓力低頻能量的急劇增強,這將使得飛行器結(jié)構(gòu)出現(xiàn)振動疲勞問題,進而嚴重影響飛行安全.因此,深入研究激波/湍流邊界層干擾區(qū)內(nèi)壓力脈動特性對認識和理解激波與湍流邊界層的相互作用機制非常重要.

    自20 世紀70 年代以來,國內(nèi)外學(xué)者一直致力于對該問題的研究.早期研究主要以風(fēng)洞試驗為主,受限于實驗測量技術(shù),風(fēng)洞試驗研究多局限于物面壓力脈動信號,并以此獲得干擾區(qū)內(nèi)典型流動結(jié)構(gòu)的特征信息,如分離激波非定常運動特性、湍流結(jié)構(gòu)、速度場等.Settle 等[3]研究了激波強度對壓縮拐角物面壓力脈動強度的影響規(guī)律.Ringuette 等[4]研究表明,低雷諾數(shù)下,干擾區(qū)內(nèi)物面壓力脈動強度降低而間歇性增強.在分離激波低頻振蕩方面,Erengil等[5]采用條件取樣和變窗口系綜平均技術(shù)分析了物面壓力分布與分離激波非定常運動的相關(guān)性.結(jié)果表明,分離泡內(nèi)的低頻壓力脈動與激波位置密切相關(guān):隨著激波往下游移動,分離泡內(nèi)壓力升高,反之則降低.Dolling 等[6]比較了不同構(gòu)型誘導(dǎo)的激波/湍流邊界層干擾區(qū)物面壓力脈動特性,研究發(fā)現(xiàn),分離激波的低頻振蕩機制主要來源于下游流場的局部或全局脈動,其過零頻率與上游湍流邊界層的特征頻率無關(guān).此外,Selig 等[7]基于物面壓力與邊界層質(zhì)量通量的時空關(guān)聯(lián)分析,給出了干擾區(qū)內(nèi)湍流結(jié)構(gòu)的演化規(guī)律.Piponniau 等[8]采用本征正交分解和線性隨機估計方法,通過實驗測得的物面壓力脈動信息對速度場進行了重構(gòu).

    由于風(fēng)洞試驗?zāi)軌虻玫降膲毫π畔O為有限,為此國內(nèi)外學(xué)者針對干擾區(qū)內(nèi)壓力脈動特性也開展了大量的高精度數(shù)值模擬研究.特別是近些年,直接數(shù)值模擬方法(direct numerical simulation,DNS)已經(jīng)在激波/湍流邊界層干擾問題方面取得了長足進展,其中Pirozzoli 等[9-11]的工作最具代表性.Pirozzoli等[9]分析比較了干擾區(qū)內(nèi)物面壓力脈動與剪切層壓力脈動的差異,并基于聲共振機制成功解釋了分離激波的非定常運動現(xiàn)象.Bernardini 等[10]著重探討了干擾區(qū)內(nèi)物面壓力脈動特性的演化規(guī)律,如脈動強度、時空關(guān)聯(lián)、對流速度和頻譜等,發(fā)現(xiàn)物面壓力脈動對流速度明顯降低而其展向尺度則顯著增強.隨后,Volpiani 等[11]進一步考察了壁面溫度對物面壓力脈動特性的影響規(guī)律.與此同時,國內(nèi)童福林等針對物面壓力脈動特性也開展了較為詳細的直接數(shù)值模擬研究,歸納了激波強度[12]、壁面溫度[13]、膨脹效應(yīng)[14-15]等因素的影響規(guī)律.

    盡管在激波/湍流邊界層干擾區(qū)壓力脈動特性研究方面,國內(nèi)外已開展了大量的風(fēng)洞試驗和數(shù)值研究工作,但是這些研究大多關(guān)注物面壓力脈動特性及其與分離激波非定常運動的相互關(guān)系方面,而關(guān)于邊界層內(nèi)壓力脈動在干擾區(qū)內(nèi)演化機制的認識仍然是很有限的.實際上,Duan 等[16]的研究結(jié)果表明,邊界層內(nèi)壓力脈動結(jié)構(gòu)特征與物面壓力脈動完全不同,其特征頻率沿法向?qū)⒔档投卣鞒叨葘⒓眲≡龃?然而,目前學(xué)術(shù)界對邊界層內(nèi)外層大尺度壓力脈動結(jié)構(gòu)與分離激波的相互作用機制仍然缺乏充分的認識,亟待開展相關(guān)研究工作.

    本文采用直接數(shù)值模擬方法研究M∞=2.25,Reθ=3567 來流下激波角β=33.2°的入射激波/平板湍流邊界層干擾區(qū)內(nèi)壓力脈動特性.為了便于比較和驗證結(jié)果,計算參數(shù)的選取與Dupont 等[17]的實驗和Fang 等[18]的DNS 相近.計算結(jié)果與前人的DNS 結(jié)果及風(fēng)洞實驗數(shù)據(jù)進行了仔細驗證,著重探討干擾區(qū)外層壓力脈動與物面壓力脈動演化機制的差異,特別是激波干擾對脈動強度、功率譜、特征尺度及時空關(guān)聯(lián)特性等的影響規(guī)律.

    1 DNS 計算

    控制方程為曲線坐標系(τ,ξ,η,ζ)下的三維可壓縮無量綱Navier-Stokes 方程組

    其中J?1為直角坐標系(x,y,z)變換為曲線坐標系(ξ,η,ζ) 的Jacobian 矩陣,U為守恒變量,F,G和H分別為(x,y,z)方向上對應(yīng)的無黏通量,Fv,Gv和Hv分別為(x,y,z)方向上對應(yīng)的黏性通量,具體表達式參見文獻[19].本文采用的解算器是高精度有限差分軟件OPENCFD-SC,該求解器目前已在壓縮拐角[19-20]、膨脹角入射激波/湍流邊界層干擾[14-15,21]等復(fù)雜流動問題中得到廣泛應(yīng)用和驗證確認,可以保證計算結(jié)果的可靠性.

    對流項計算時,采用Martin 等[22]優(yōu)化構(gòu)造的WENO_SYMBO_LMT 格式以及Steger-Warming流通量分裂方法,這樣能夠在極大抑制激波間斷區(qū)的強數(shù)值振蕩同時也保證了對邊界層內(nèi)不同尺度湍流脈動結(jié)構(gòu)的高分辨率捕捉.黏性項采用八階中心差分格式進行離散,時間方向上推進采用三階Runge-Kutta 方法.

    圖1 為計算模型示意圖,來流方向為圖中從左往右,來流參數(shù)如下:馬赫數(shù)M∞=2.25,基于單位長度來流雷諾數(shù)Re∞=2.5×104mm?1,靜溫T∞=169.44 K.模型長度Lx×Ly×Lz=137.6 mm×12.7 mm×4.4 mm,這里L(fēng)x,Ly和Lz分別為流向、法向和展向長度,坐標系原點取為計算域入口處.采用與文獻[21]類似的方法生成入射激波,在計算域上邊界xin=82.2 mm左右兩側(cè)分別設(shè)置為波前和波后參數(shù),即:在xin左側(cè)流動參數(shù)取為自由來流參數(shù),而在xin右側(cè)則按照激波關(guān)系式Rankine-Hugoniot 給出.這里激波角β取為33.2°,圖中入射激波在壁面上的名義入射點為xs=101.6 mm.

    圖1 計算模型Fig.1 Computational model

    計算網(wǎng)格點數(shù)為Nx×Ny×Nz=3700×300×250 (流向×法向×展向).圖2 為計算網(wǎng)格示意圖,為了便于展示,流向間隔了10 個點,法向間隔了5 個點.網(wǎng)格生成采用代數(shù)解析方法,這里流向網(wǎng)格點采用三段分布:上游轉(zhuǎn)捩區(qū)0 mm

    圖2 計算網(wǎng)格示意圖Fig.2 Sketch of the computational grid

    邊界條件設(shè)置如下:計算域入口處取為相同來流下的層流解(laminar inflow),具體分布如圖3 所示.壁面取為無滑移和等溫壁,壁溫Tw=321.9 K,上邊界采用簡單無反射邊界條件.出口處采用超聲速出口無反射條件及緩沖區(qū),以便消除擾動波反射帶來的影響.計算域展向兩側(cè)采用周期性邊界條件.

    圖3 入口處層流剖面Fig.3 Laminar profile at inlet

    計算時,通過在上游壁面添加多頻吹吸擾動從而轉(zhuǎn)捩生成充分發(fā)展湍流邊界層,進而與下游入射激波產(chǎn)生相互作用.這里,法向擾動速度Vbs的具體表達式如下[18]

    DNS 計算時取lmax=10,mmax=5,擾動幅值A(chǔ)=0.2,擾動頻率?=0.628U∞/δ,這里δ表示參考點處邊界層厚度,下文類似.擾動帶起始點分別為xa=7.62 mm 和xb=20.32 mm.式中相位差φl和φm取為0 到1 的隨機數(shù).

    表1 給出了參考點xref處湍流邊界層參數(shù)與Fang 等[18]DNS 結(jié)果比較情況,表中Reδ,Reδ*,Reθ分別為基于邊界層厚度、位移厚度和動量厚度的雷諾數(shù),Cf為摩阻系數(shù).兩者較為接近,這也表明本文計算結(jié)果的準確性.為了進一步驗證本文DNS 結(jié)果的可靠性,還與以往數(shù)值模擬結(jié)果和風(fēng)洞實驗數(shù)據(jù)進行了比較,包括湍流邊界層的速度剖面、脈動強度,物面壓力頻譜以及平均壓力分布.

    表1 參考點xref 處湍流邊界層參數(shù)Table 1 Turbulent boundary layer parameters at xref

    首先,圖4~ 圖6 分別給出了上游湍流邊界層在參考點xref處的統(tǒng)計特性.圖4 為平均流向速度剖面,定義為

    圖4 平均流向速度剖面Fig.4 Mean streamwise velocity profile at xref

    式中ρ和u分別為密度和流向速度,上標 ? 表示平均量.需要特別說明的是,本文中平均指的是時間和展向平均.速度分布在黏性底層y+<10 符合線性律,在30

    圖5 湍流脈動強度Fig.5 Turbulence intensities at xref

    圖6 物面壓力功率譜Fig.6 Power spectral density for wall pressure at xref

    此外,圖7 還給出了干擾區(qū)內(nèi)物面壓力分布情況,其中物面壓力P*和流向坐標x*分別定義為

    圖7 物面壓力分布Fig.7 Distribution of wall pressure

    式中,P1是反射激波的波后壓力值,xr*為轉(zhuǎn)換后的反射激波腳對應(yīng)的平均流向位置,這里采用與Fang等[18]類似的方法,即(P1+P∞)/2 對應(yīng)的流向位置.圖中紅色虛線代表壓力分布的無黏解.可見,計算得到的壓力分布與Fang 等[18]的DNS 結(jié)果及Dupont 等[17]的實驗結(jié)果相符較好.在黏性干擾作用下,干擾區(qū)內(nèi)壓力梯度變化較為緩和且在下游逐步趨近于無黏解.

    2 脈動強度

    圖8 給出了上游湍流邊界層(turbulent boundary layer,TBL)壓力脈動強度沿法向分布情況.可見,峰值位置出現(xiàn)在近壁區(qū),約為y/δ=0.04,其量值沿物面法向急劇降低,在邊界層外緣,約為近壁峰值的28%.采用當(dāng)?shù)匚锩婕羟袘?yīng)力τw進行無量綱化后,計算得到的脈動分布與Bernardini 等[28]和Duan 等[16]的DNS 結(jié)果符合較好.

    圖8 上游湍流邊界層壓力脈動強度分布Fig.8 Distribution of pressure fluctuation intensity for the incoming TBL

    從圖9 可以看出,與上游湍流邊界層相比,干擾區(qū)內(nèi)壓力脈動強度呈現(xiàn)截然不同的分布規(guī)律.從量值上來看,脈動峰值呈現(xiàn)先急劇升高,這主要是激波的增強作用,隨后在下游再附區(qū)伴隨著擾動邊界層的恢復(fù)則呈現(xiàn)緩慢衰減的流向演化趨勢.

    圖9 干擾區(qū)壓力脈動強度云圖Fig.9 Contour of pressure fluctuation intensity in the interaction region

    為了定量描述這一現(xiàn)象,圖10 給出了物面壓力脈動強度沿流向的分布情況.計算得到脈動強度與Fang 等[20]的DNS 結(jié)果符合良好,兩者在整體往下偏移0.02 的情況下與Dupont 等[17]的實驗數(shù)據(jù)也較為吻合.這主要是由于壓力傳感器截止頻率太小導(dǎo)致實驗數(shù)據(jù)明顯偏低[17].圖中符號S和R分別表示分離點和再附點位置,通過平均摩阻的過零點確定.可見,分離區(qū)內(nèi)物面壓力脈動急劇增強,其峰值約為上游的4.3 倍,隨后在下游再附區(qū)脈動峰值緩慢降低,在x*?xr*=12 時,仍約為上游的2.2 倍.另外,從峰值位置來看,其在分離區(qū)呈現(xiàn)逐漸遠離壁面的態(tài)勢,這與激波誘導(dǎo)生成分離區(qū)剪切層相關(guān),如圖9中平均音速線所示(粉色點劃線).強脈動主要出現(xiàn)在干擾區(qū)內(nèi)壓縮波系、剪切層附近等強剪切區(qū)域.此后,在下游再附邊界層,其峰值位置仍集中在外層區(qū)域.

    圖10 物面壓力脈動強度分布Fig.10 Distribution of wall pressure fluctuation intensity

    3 功率譜

    沿外層y/δ=0.8 和物面分別選取5 個流向站位進行對比分析,其中站位1 位于上游湍流邊界層內(nèi),站位2 位于分離區(qū)內(nèi),站位3~ 5 位于下游再附邊界層內(nèi),各站位具體位置見圖9.如無特別說明,下文中外層均為y/δ=0.8 處.

    圖11 比較了上游湍流邊界層站位1 壓力脈動的預(yù)乘功率譜ωφ(ω),其中ω和φ(ω)分別為角頻率和功率譜密度.功率譜密度的計算采用Pwelch 方法和Hamming 窗函數(shù),同時對分段數(shù)據(jù)按50%重疊處理.為了便于定量比較,圖中還采用當(dāng)?shù)孛}動強度的平方進行了歸一化處理.顯而易見,外層壓力脈動的時間尺度明顯要大于物面處,前者峰值頻率約為ωδ/U∞=3.0,后者約為ωδ/U∞=12.3,這與Duan 等[13]的研究結(jié)論是一致的.這是由于外層壓力脈動主要表征為大尺度脈動結(jié)構(gòu).

    圖11 站位1 壓力脈動預(yù)乘功率譜Fig.11 Pre-multiplied power spectral density for pressure fluctuations at station 1

    圖12 分別給出了干擾區(qū)外層和物面壓力脈動預(yù)乘功率譜的變化情況.從兩者的比較來看,以下幾個方面值得特別關(guān)注.首先,相較于上游站位1,分離區(qū)站位2 處兩者都出現(xiàn)了峰值頻率向更低頻區(qū)偏移的現(xiàn)象,這與激波的低頻振蕩運動密切相關(guān).

    圖12 干擾區(qū)內(nèi)不同站位壓力脈動預(yù)乘功率譜Fig.12 Pre-multiplied power spectral density for pressure fluctuations at various streamwise locations

    如圖12(a)和圖12(b)所示,外層和物面壓力峰值頻率分別約為ωδ/U∞=0.7 和ωδ/U∞=4.4.以往研究表明[29],在分離激波低頻振蕩的作用下,分離區(qū)內(nèi)物面壓力脈動的時間尺度通常比上游湍流邊界層低1~ 2 個數(shù)量級.本文計算結(jié)果進一步驗證了該結(jié)論,同時結(jié)果還表明,外層壓力脈動也同樣符合這一規(guī)律.在功率譜的定性分布方面,顯然,外層壓力脈動的變化較物面壓力更為劇烈,這很可能是由于外層站位非常靠近分離激波,受其非定常運動的影響更為直接,而物面壓力處于分離泡內(nèi)(見圖9),影響則要弱得多.

    從再附區(qū)站位3~ 5 的結(jié)果來看,外層壓力脈動呈現(xiàn)快速恢復(fù)的態(tài)勢,且峰值頻率向高頻區(qū)偏移,逐步趨近于上游站位1,這與物面壓力脈動的演化過程則完全不同.從圖12(a)中可以清楚看到,再附邊界層內(nèi)物面壓力脈動的低頻能量占比仍相對較高,峰值頻率依然維持在較低頻區(qū),并沒有出現(xiàn)往高頻區(qū)恢復(fù)的現(xiàn)象.造成兩者該差異的原因可以認為主要有兩方面:一方面,上游激波誘導(dǎo)分離形成了剪切層,分離泡上方剪切層內(nèi)不穩(wěn)定渦結(jié)構(gòu)沿流向急劇發(fā)展,并進入到下游再附邊界層內(nèi),這些不穩(wěn)定結(jié)構(gòu)極大地抑制了物面壓力脈動結(jié)構(gòu)的恢復(fù);另一方面,外層壓力脈動法向位置相對較高,盡管在站位2 處激波干擾對其影響顯著,但是過激波后外層壓力脈動受近壁區(qū)流動結(jié)構(gòu)的影響則要小得多,主要表征為外層大尺度脈動結(jié)構(gòu)的破碎過程,更為詳細的定量分析可見下一節(jié)的兩點相關(guān)統(tǒng)計結(jié)果.

    4 兩點相關(guān)

    通過兩點相關(guān)性考察干擾區(qū)外層壓力脈動擬序結(jié)構(gòu)長度尺度的變化規(guī)律,且定量分析比較其與物面壓力脈動結(jié)構(gòu)尺度的差異.兩點相關(guān)系數(shù)Cpp(Δx,Δz)定義為

    式中p′為雷諾平均后得到的壓力脈動,Δx和Δz分別代表流向間距和展向間距,x0為特征點流向位置,yref分別取物面和外層,上橫線代表時間和展向平均.

    圖13 分別給出了上游湍流邊界層外層和物面脈動壓力兩點相關(guān)系數(shù)的分布情況.以往的零壓力梯度平板結(jié)果表明[28],小間距時脈動壓力相關(guān)性曲線可近似為圓形,這是小尺度壓力脈動結(jié)構(gòu)的各向同性所決定的;大間距時,大尺度脈動結(jié)構(gòu)的各向異性導(dǎo)致相關(guān)性曲線在展向出現(xiàn)了急劇拉伸,近似為橢圓形分布,表征了展向拉長的大尺度結(jié)構(gòu).顯然,圖13(a)中的結(jié)果也完全符合該規(guī)律.然而,對于外層壓力脈動,如圖13(b)所示,整體來看,外層脈動結(jié)構(gòu)的流向尺度和展向尺度都明顯大于物面壓力脈動.這里以相關(guān)性系數(shù)0.3 為參考確定其特征長度,外層結(jié)構(gòu)的流向和展向尺度分別約為0.7δ和1.0δ,而物面壓力脈動則分別約為0.14δ和0.2δ.此外,盡管外層相關(guān)性等值線小間距時仍趨近于圓形,但隨著間距增大,等值線在流向和展向都被急劇拉伸,展向明顯強于流向,這表明外層壓力脈動擬序結(jié)構(gòu)仍以展向大尺度結(jié)構(gòu)為主,但其特征尺度要大得多.

    圖13 站位1 壓力脈動兩點相關(guān)系數(shù)Fig.13 Two-point correlation coefficient of pressure fluctuations at station 1

    圖14(a)和圖14(b)分別給出了站位2 物面和外層壓力脈動的兩點相關(guān)等值線,其中等值線從外到內(nèi)依次為0.3,0.5,0.7 和0.9.為了便于比較說明,圖中也給出了站位1 的結(jié)果.從分布規(guī)律上來看,等值線形態(tài)基本類似于上游湍流邊界層,仍然以展向拉伸的橢圓形分布為主,但整體呈現(xiàn)向外擴張的趨勢,這表明激波干擾作用增大了外層和物面壓力脈動結(jié)構(gòu)尺度,特別是在展向.需要特別關(guān)注的是,從圖14(b)中可以清楚看到,同樣以等值線0.3 為參考,站位1 處的展向尺度約為1.0δ,而站位2 處增大到約2.0δ,而流向尺度兩者基本維持在0.7δ.這很可能是由于站位2 外層位于入射激波附近(見圖9),激波的展向結(jié)構(gòu)增強了壓力脈動的展向相關(guān)性.更為細致的定量分析還有待深入研究.

    圖14 站位2 兩點相關(guān)系數(shù)(虛線:站位1)Fig.14 Two-point correlation coefficient of pressure fluctuations at station 2(dashed lines:station 1)

    如圖15 所示,再附區(qū)站位3~ 5 外層和物面壓力脈動的兩點相關(guān)等值線與上游湍流邊界層也基本類似.小間距時,曲線接近于圓形,這表明小尺度結(jié)構(gòu)并沒有發(fā)生實質(zhì)改變,仍以各向同性為主.從大尺度結(jié)構(gòu)的演化歷程來看,兩者差異顯著.研究發(fā)現(xiàn),大間距時,圖15(a)~ 圖15(c)中等值線以向外擴張為主,這表明再附區(qū)內(nèi)物面壓力大尺度脈動結(jié)構(gòu)的特征尺度仍在增大,而圖15(d)~ 圖15(f)等值線則以向內(nèi)收縮為主,對應(yīng)的是外層壓力大尺度脈動結(jié)構(gòu)特征尺度的減小.

    圖15 站位3~ 5 兩點相關(guān)系數(shù)(上:物面;下:y/δ=0.8;虛線:站位1)Fig.15 Two-point correlation coefficient at station 3~ 5 (top:wall;down:y/δ=0.8;dashed lines:station 1)

    為了定量描述該現(xiàn)象,圖16 進一步給出了干擾區(qū)內(nèi)物面及外層壓力脈動積分尺度的分布.積分尺度可以用來表征脈動場最大尺度結(jié)構(gòu)的特征尺度,這里流向和展向積分尺度分別定義為[10]

    可見,外層壓力脈動流向和展向積分尺度遠大于物面壓力脈動,這與Duan 等[16]的研究結(jié)果是一致的.從整體變化規(guī)律上來看,外層壓力Λx和Λz在站位2 處急劇增大,隨后在站位3~ 5 則又緩慢減小,而物面壓力的Λx和Λz整體上呈現(xiàn)逐步增大的趨勢,這與圖15 中相關(guān)性分析結(jié)果是相符的.從量值變化上來看,以圖16(b)中的展向積分尺度Λz為例,對于外層壓力,在激波的增強作用下,站位1 處Λz=0.73δ增大到站位2 處Λz=1.17δ,隨后在再附區(qū)又緩慢衰減到站位5 的Λz=1.02δ;而物面壓力,從站位1 處Λz=0.26δ逐步增大到站位5 處的Λz=0.82δ,后者約為前者的3.15 倍.可以看到,在站位5 處,此時再附區(qū)內(nèi)物面壓力脈動結(jié)構(gòu)的展向和流向積分尺度已非常接近于外層的積分尺度.綜上所述,可以認為這是由于外層壓力脈動大尺度結(jié)構(gòu)在過激波后呈現(xiàn)快速破碎的過程,而物面上大尺度結(jié)構(gòu)仍在緩慢增長,這也解釋了圖12 中物面壓力脈動低頻能量的相對升高而外層壓力脈動峰值頻率往高頻區(qū)的偏移現(xiàn)象.

    圖16 壓力脈動積分尺度分布Fig.16 Distribution of integral scale for pressure fluctuations

    5時空關(guān)聯(lián)

    進一步分析干擾區(qū)外層和物面壓力脈動的時空關(guān)聯(lián)特性,從而得到向下游傳播的壓力波的對流速度.

    圖17 分別給出了站位1 處物面和外層壓力的時空關(guān)聯(lián)云圖,圖中時空關(guān)聯(lián)系數(shù)Rpp(Δx,Δt)定義為

    式中Δt為延遲時間.以往研究結(jié)果表明[28,30],由于脈動結(jié)構(gòu)往下游傳播,物面壓力時空關(guān)聯(lián)云圖呈現(xiàn)傾斜的類橢圓形分布,如圖17(a)所示,本文計算結(jié)果也完全符合該結(jié)論.相較于物面壓力,從圖17(b)中可以看到,外層壓力時空關(guān)聯(lián)云圖的橢圓率和傾斜角都明顯增大.

    圖17 站位1 壓力脈動時空關(guān)聯(lián)系數(shù)Fig.17 Space-time correlation coefficient for pressure fluctuations at station 1

    圖18 和圖19 分別給出了分離區(qū)站位2 和再附區(qū)站位3~ 5 壓力脈動時空關(guān)聯(lián)系數(shù)與站位1 (黑色虛線)的比較情況.如圖18(a)所示,可以發(fā)現(xiàn),相較于上游湍流邊界層,分離區(qū)物面脈動壓力時空關(guān)聯(lián)曲線更為聚集,同時傾斜角也急劇減小,但整體形態(tài)仍滿足橢圓形(只是等值線0.1 較為雜亂).值得注意的是,對于站位2 外層脈動壓力,顯然,圖18(b)中等值線形態(tài)變化非常劇烈,已完全不符合橢圓形分布,特別是等值線0.1 和0.3,在時間方向上可近似為直線分布.造成這個現(xiàn)象的原因很可能是由于站位2 外層處于強激波區(qū)(見圖9),在分離激波的強干擾下,此時泰勒凍結(jié)假設(shè)完全失效,脈動壓力在時空上表征為強非線性.此外,從圖19 中還可以看到,無論是再附區(qū)外層還是物面,其脈動壓力時空關(guān)聯(lián)等值線整體形態(tài)與上游站位1 都較為類似.從站位3~站位5,強相關(guān)等值線被急劇拉長壓縮,且傾斜角增大,呈現(xiàn)逐步逼近上游站位1 的演化過程.從定性比較的角度來看,外層脈動壓力時空關(guān)聯(lián)曲線的恢復(fù)速度明顯要快于物面.

    圖18 站位2 壓力脈動時空關(guān)聯(lián)系數(shù)(虛線:站位1)Fig.18 Space-time correlation coefficient for pressure fluctuations at station 2 (dashed lines:station 1)

    圖19 站位3~ 5 壓力脈動時空關(guān)聯(lián)系數(shù)(上:物面;下:y/δ=0.8;虛線:站位1)Fig.19 Space-time correlation coefficient at station 3~ 5 (top:wall;down:y/δ=0.8;dashed lines:station 1)

    壓力脈動對流速度的計算方法如下[16,31]:首先給定延遲時間Δt,隨后依據(jù)時空關(guān)聯(lián)曲線通過

    獲得流向間距rx,進而得到對流速度Uc=rx/Δt.采用上述方法計算得到了上游湍流邊界層物面和外層脈動壓力對流速度隨延遲時間的分布情況.如圖20 和圖21 所示.總體來看,外層脈動壓力對流速度要遠大于物面脈動壓力,這與Duan 等[16]的研究結(jié)果也較為符合.另外,對于外層脈動壓力,可見,其對流速度對延遲時間的依賴性較小,基本維持在Uc=0.94U∞,與當(dāng)?shù)仄骄飨蛩俣冉咏?然而,物面脈動壓力對流速度隨延遲時間則變化劇烈,在小延遲時間(小尺度脈動)時,對流速度約為Uc=0.56U∞,而在大延遲時間(大尺度脈動) 時,其增大到Uc=0.73U∞.已有研究表明[13],這是由于物面大尺度脈動壓力與外層大尺度結(jié)構(gòu)相關(guān),而小尺度脈動壓力則主要受近壁區(qū)結(jié)構(gòu)影響,故前者對流速度明顯大于后者.

    從圖20 和圖21 中還可以看到,分離區(qū)物面脈動壓力對流速度降低,約為Uc=0.2U∞,隨后下游再附區(qū)則又逐漸增大,但仍小于上游湍流邊界層.Na等[32]在強逆壓力梯度平板的DNS 結(jié)果中也發(fā)現(xiàn)了類似規(guī)律,本文的計算結(jié)果進一步證實了該現(xiàn)象在入射激波/湍流邊界層干擾區(qū)內(nèi)的適用性.從圖20(b)中可以看到,采用當(dāng)?shù)赝鈱铀俣萓e對Uc進行歸一化后,分離區(qū)和再附區(qū)各站位處物面脈動壓力的對流速度分布均吻合較好,約為Uc=0.25Ue,與Bernardini 等[10]在跨聲速激波/湍流邊界層干擾下的研究結(jié)論也較為接近,但后者為初始分離(incipient separation)工況.如圖21 所示,干擾區(qū)外層不同站位處脈動壓力對流速度分布與上游結(jié)果較為類似,其對延遲時間的依賴性依然較弱,但其量值減小到約為Uc=0.4Ue.結(jié)合之前的功率譜和兩點相關(guān)分析,可以推測外層大尺度脈動結(jié)構(gòu)特征尺度的減小可能是造成該現(xiàn)象的重要因素.

    圖20 干擾區(qū)內(nèi)各站位物面壓力脈動對流速度Fig.20 Convection velocity for wall pressure fluctuations at various stations in the interaction region

    圖21 干擾區(qū)內(nèi)各站位外層壓力脈動對流速度Fig.21 Convection velocity for pressure fluctuations in the outer layer

    6 結(jié)論

    本文采用直接數(shù)值模擬方法研究了來流馬赫數(shù)2.25 下激波角33.2°的入射激波對平板湍流邊界層壓力脈動特性的影響規(guī)律,詳細地分析了外層和物面壓力脈動場的典型統(tǒng)計特征,如脈動強度、功率譜密度、兩點相關(guān)和時空關(guān)聯(lián)特性等,得到以下結(jié)論:

    (1) 數(shù)值模擬準確捕捉到了干擾區(qū)內(nèi)脈動壓力.研究發(fā)現(xiàn),分離區(qū)外層脈動壓力以低頻能量為主,高頻能量占比較小.與物面脈動壓力不同的是,再附區(qū)外層脈動壓力峰值頻率逐漸往高頻區(qū)偏移,而物面壓力脈動能量仍維持在相對較低頻率.

    (2) 激波干擾對壓力脈動結(jié)構(gòu)有顯著的增強作用,干擾區(qū)內(nèi)展向積分尺度明顯大于流向積分尺度.積分尺度在外層和物面的演化歷程差異明顯,前者經(jīng)歷了先急劇增長后緩慢衰減的過程,而后者整體上呈現(xiàn)增長趨勢.

    (3) 激波干擾極大地降低了壓力脈動場的對流速度.時空關(guān)聯(lián)結(jié)果表明,干擾區(qū)下游外層壓力波的對流速度仍明顯高于物面,前者約為0.4Ue,而后者約為0.25Ue.

    猜你喜歡
    物面外層邊界層
    一種溶液探測傳感器
    傳感器世界(2022年4期)2022-11-24 21:23:50
    基于HIFiRE-2超燃發(fā)動機內(nèi)流道的激波邊界層干擾分析
    讓吸盤掛鉤更牢固
    一種購物袋
    科技資訊(2016年6期)2016-05-14 13:09:55
    新型單面陣自由曲面光學(xué)測量方法成像特性仿真
    一類具有邊界層性質(zhì)的二次奇攝動邊值問題
    專題Ⅱ 物質(zhì)構(gòu)成的奧秘
    “人”字變身
    賽車空動優(yōu)化的秘密
    汽車之友(2014年21期)2014-11-03 17:33:45
    非特征邊界的MHD方程的邊界層
    中文字幕免费在线视频6| 国产 一区精品| 97超碰精品成人国产| 亚洲欧美一区二区三区国产| 日韩三级伦理在线观看| 日韩av在线免费看完整版不卡| 蜜桃在线观看..| 大香蕉97超碰在线| 大码成人一级视频| 九九在线视频观看精品| 久久精品aⅴ一区二区三区四区 | 日本猛色少妇xxxxx猛交久久| 99香蕉大伊视频| 亚洲少妇的诱惑av| 日韩视频在线欧美| 在线观看人妻少妇| 亚洲人成网站在线观看播放| 人妻一区二区av| 性高湖久久久久久久久免费观看| av.在线天堂| av在线观看视频网站免费| 日韩伦理黄色片| 精品一区二区免费观看| 国产成人免费观看mmmm| 一本—道久久a久久精品蜜桃钙片| 十分钟在线观看高清视频www| 精品人妻偷拍中文字幕| a级毛色黄片| 午夜免费鲁丝| 久久午夜综合久久蜜桃| av又黄又爽大尺度在线免费看| 国产熟女午夜一区二区三区| 91午夜精品亚洲一区二区三区| 蜜臀久久99精品久久宅男| 亚洲一区二区三区欧美精品| 亚洲色图 男人天堂 中文字幕 | 观看av在线不卡| 亚洲av日韩在线播放| 亚洲一码二码三码区别大吗| 国产极品天堂在线| 日韩制服骚丝袜av| 精品国产一区二区久久| 一本—道久久a久久精品蜜桃钙片| 一本—道久久a久久精品蜜桃钙片| 美女国产视频在线观看| 自拍欧美九色日韩亚洲蝌蚪91| 成人黄色视频免费在线看| 交换朋友夫妻互换小说| www.av在线官网国产| 久热这里只有精品99| 香蕉精品网在线| 在线免费观看不下载黄p国产| 亚洲,一卡二卡三卡| www.色视频.com| 两性夫妻黄色片 | 天天操日日干夜夜撸| 不卡视频在线观看欧美| 亚洲三级黄色毛片| 在线观看www视频免费| 女性被躁到高潮视频| 看免费av毛片| 国产欧美另类精品又又久久亚洲欧美| av播播在线观看一区| 高清黄色对白视频在线免费看| 国产精品一区二区在线观看99| 免费看不卡的av| 久久久久久人人人人人| 成年人午夜在线观看视频| 亚洲欧洲日产国产| 最新中文字幕久久久久| 日韩av在线免费看完整版不卡| 日韩熟女老妇一区二区性免费视频| 两个人免费观看高清视频| 欧美成人午夜精品| 人妻人人澡人人爽人人| 2021少妇久久久久久久久久久| 99热6这里只有精品| 中文字幕亚洲精品专区| 免费黄频网站在线观看国产| 国产精品久久久av美女十八| 老女人水多毛片| 狠狠精品人妻久久久久久综合| 一个人免费看片子| 久久亚洲国产成人精品v| 亚洲精品一二三| 看免费成人av毛片| 午夜福利在线观看免费完整高清在| 大话2 男鬼变身卡| 国产熟女欧美一区二区| 最新的欧美精品一区二区| 成人黄色视频免费在线看| 欧美人与性动交α欧美精品济南到 | 亚洲美女视频黄频| 寂寞人妻少妇视频99o| 免费女性裸体啪啪无遮挡网站| 日韩av不卡免费在线播放| 免费观看a级毛片全部| 少妇人妻久久综合中文| 国产精品人妻久久久影院| 爱豆传媒免费全集在线观看| 少妇的丰满在线观看| 一级片'在线观看视频| 一边亲一边摸免费视频| 久久久久视频综合| 国产男人的电影天堂91| 97在线视频观看| 一级毛片电影观看| 久久久久久久精品精品| 亚洲内射少妇av| 99热网站在线观看| 日韩中字成人| 七月丁香在线播放| 精品久久国产蜜桃| 色婷婷av一区二区三区视频| 国产精品一国产av| 69精品国产乱码久久久| 少妇人妻久久综合中文| 少妇人妻久久综合中文| 国产一区亚洲一区在线观看| 久久人人爽人人爽人人片va| 夫妻午夜视频| 久久午夜综合久久蜜桃| 色视频在线一区二区三区| 国产成人a∨麻豆精品| 王馨瑶露胸无遮挡在线观看| 亚洲成人手机| 夜夜骑夜夜射夜夜干| 在线 av 中文字幕| 亚洲国产精品成人久久小说| 建设人人有责人人尽责人人享有的| 秋霞在线观看毛片| 亚洲成人手机| 色哟哟·www| 黑人猛操日本美女一级片| 久久综合国产亚洲精品| av在线播放精品| 国产av国产精品国产| 一区二区三区乱码不卡18| 亚洲欧美清纯卡通| a级毛片黄视频| 国产亚洲一区二区精品| 久久久久久人人人人人| 人妻 亚洲 视频| 深夜精品福利| 国产女主播在线喷水免费视频网站| 日韩伦理黄色片| 十分钟在线观看高清视频www| 777米奇影视久久| 国产高清国产精品国产三级| 啦啦啦中文免费视频观看日本| 欧美人与性动交α欧美软件 | 另类亚洲欧美激情| 中文字幕免费在线视频6| 新久久久久国产一级毛片| 国内精品宾馆在线| 精品卡一卡二卡四卡免费| 午夜激情久久久久久久| av免费在线看不卡| 精品一区二区免费观看| 国产高清不卡午夜福利| 国产欧美亚洲国产| 久久久久久久久久成人| 亚洲国产精品999| a级片在线免费高清观看视频| 一区二区日韩欧美中文字幕 | 乱人伦中国视频| 国产日韩欧美在线精品| 日韩一区二区视频免费看| 美女脱内裤让男人舔精品视频| videosex国产| 亚洲人成77777在线视频| 久久国产精品大桥未久av| 大陆偷拍与自拍| 国产熟女欧美一区二区| 国产成人精品在线电影| 女人久久www免费人成看片| 黄色配什么色好看| 欧美激情 高清一区二区三区| 考比视频在线观看| 青春草视频在线免费观看| 久久久久精品人妻al黑| 精品熟女少妇av免费看| 亚洲色图 男人天堂 中文字幕 | 国产免费视频播放在线视频| 国产亚洲精品久久久com| 丝袜脚勾引网站| 久久精品久久久久久久性| 亚洲av在线观看美女高潮| 久久婷婷青草| 中文字幕精品免费在线观看视频 | 免费人成在线观看视频色| 免费av不卡在线播放| 久久精品亚洲av国产电影网| 欧洲精品卡2卡3卡4卡5卡区| 亚洲专区字幕在线| 成人国产一区最新在线观看| 男女床上黄色一级片免费看| 99riav亚洲国产免费| 多毛熟女@视频| 欧美不卡视频在线免费观看 | 国产免费男女视频| 99国产极品粉嫩在线观看| 日韩三级视频一区二区三区| 国产精品 国内视频| 久久久久视频综合| 国产精品乱码一区二三区的特点 | 手机成人av网站| 日本五十路高清| 精品国产乱码久久久久久男人| 成熟少妇高潮喷水视频| 丝瓜视频免费看黄片| 黄频高清免费视频| 婷婷丁香在线五月| 女性被躁到高潮视频| av在线播放免费不卡| 久久国产精品大桥未久av| 国产精品免费大片| av天堂在线播放| 母亲3免费完整高清在线观看| 99精品在免费线老司机午夜| 日韩欧美一区二区三区在线观看 | 国产精品亚洲一级av第二区| 少妇裸体淫交视频免费看高清 | 欧美中文综合在线视频| 老熟女久久久| a级毛片在线看网站| 真人做人爱边吃奶动态| 国产精品国产av在线观看| 美女 人体艺术 gogo| 国产精品99久久99久久久不卡| 十八禁网站免费在线| 91麻豆av在线| 乱人伦中国视频| 免费看十八禁软件| 99国产精品一区二区三区| 99精品欧美一区二区三区四区| 国产精品亚洲av一区麻豆| 国产亚洲av高清不卡| 午夜福利在线免费观看网站| 人妻丰满熟妇av一区二区三区 | 亚洲一卡2卡3卡4卡5卡精品中文| 国产精品98久久久久久宅男小说| 亚洲性夜色夜夜综合| 国产人伦9x9x在线观看| 无遮挡黄片免费观看| 欧美日韩成人在线一区二区| 不卡av一区二区三区| 亚洲精品美女久久av网站| 老司机福利观看| 婷婷成人精品国产| 一边摸一边抽搐一进一小说 | 精品福利永久在线观看| 亚洲一区二区三区不卡视频| 黄片小视频在线播放| 18禁裸乳无遮挡动漫免费视频| 精品一区二区三区视频在线观看免费 | 精品高清国产在线一区| 欧美黄色片欧美黄色片| 国产精品亚洲av一区麻豆| 超碰97精品在线观看| 满18在线观看网站| 国产亚洲精品一区二区www | 欧美日韩瑟瑟在线播放| 精品久久久久久久毛片微露脸| av一本久久久久| 十八禁高潮呻吟视频| 黄片大片在线免费观看| 精品福利观看| 少妇裸体淫交视频免费看高清 | 看黄色毛片网站| 日韩欧美一区视频在线观看| 老司机午夜福利在线观看视频| 12—13女人毛片做爰片一| 亚洲国产欧美一区二区综合| 99国产精品一区二区三区| 黄色视频不卡| 亚洲在线自拍视频| 精品国产一区二区三区四区第35| aaaaa片日本免费| 热re99久久精品国产66热6| 男女之事视频高清在线观看| 美女福利国产在线| 十八禁网站免费在线| www.精华液| 午夜免费成人在线视频| 欧美精品一区二区免费开放| 一二三四在线观看免费中文在| av有码第一页| 久久人妻av系列| 亚洲国产看品久久| 中文字幕色久视频| 大型黄色视频在线免费观看| 国产熟女午夜一区二区三区| 精品国产亚洲在线| 18禁裸乳无遮挡免费网站照片 | 91精品三级在线观看| 国产激情久久老熟女| 国产男靠女视频免费网站| 国产色视频综合| 又黄又爽又免费观看的视频| 久久人妻熟女aⅴ| 丁香六月欧美| 国产一区二区三区综合在线观看| 9色porny在线观看| 国产蜜桃级精品一区二区三区 | 久久久精品免费免费高清| 亚洲欧美激情综合另类| 久久久久国产精品人妻aⅴ院 | 成年女人毛片免费观看观看9 | 91麻豆精品激情在线观看国产 | 日韩一卡2卡3卡4卡2021年| 曰老女人黄片| 老汉色av国产亚洲站长工具| 亚洲免费av在线视频| 午夜两性在线视频| 女性生殖器流出的白浆| 大片电影免费在线观看免费| 亚洲欧美日韩另类电影网站| 久久久久久久午夜电影 | 久久香蕉激情| 午夜福利影视在线免费观看| 久久中文看片网| 国产麻豆69| 女同久久另类99精品国产91| 欧美日韩中文字幕国产精品一区二区三区 | 国产高清国产精品国产三级| 少妇粗大呻吟视频| 日韩视频一区二区在线观看| 一级片'在线观看视频| 日本黄色视频三级网站网址 | 怎么达到女性高潮| videos熟女内射| 夜夜爽天天搞| 91老司机精品| 免费黄频网站在线观看国产| 1024视频免费在线观看| 我的亚洲天堂| 18禁观看日本| 国产成人精品无人区| 久久午夜综合久久蜜桃| 18禁观看日本| 男女下面插进去视频免费观看| 黑人欧美特级aaaaaa片| 久久精品国产a三级三级三级| 51午夜福利影视在线观看| 国产人伦9x9x在线观看| 黑人巨大精品欧美一区二区蜜桃| 亚洲av日韩在线播放| 看免费av毛片| 欧美在线黄色| 午夜精品久久久久久毛片777| 亚洲一区中文字幕在线| 日韩三级视频一区二区三区| 老司机午夜福利在线观看视频| 欧美日韩黄片免| 天堂中文最新版在线下载| 精品国产乱子伦一区二区三区| 久久国产精品男人的天堂亚洲| 亚洲熟女精品中文字幕| 两性午夜刺激爽爽歪歪视频在线观看 | 欧美国产精品一级二级三级| 如日韩欧美国产精品一区二区三区| 性色av乱码一区二区三区2| 久久国产精品人妻蜜桃| 精品久久久久久,| 两性夫妻黄色片| 欧美国产精品va在线观看不卡| 91精品三级在线观看| 99久久国产精品久久久| 国产精品免费视频内射| 久久久久久久午夜电影 | 亚洲精品国产色婷婷电影| 精品福利观看| 99精品欧美一区二区三区四区| 老司机午夜十八禁免费视频| 亚洲午夜精品一区,二区,三区| 在线免费观看的www视频| 日本wwww免费看| 亚洲欧美日韩另类电影网站| 欧美日韩中文字幕国产精品一区二区三区 | 最新在线观看一区二区三区| 国产精品1区2区在线观看. | 国产精品自产拍在线观看55亚洲 | 久久久久国产精品人妻aⅴ院 | 国产精品成人在线| 久久久久久人人人人人| 天堂动漫精品| 国产色视频综合| 免费在线观看日本一区| 婷婷成人精品国产| 欧美黄色片欧美黄色片| а√天堂www在线а√下载 | x7x7x7水蜜桃| 丰满迷人的少妇在线观看| 老司机在亚洲福利影院| 精品亚洲成a人片在线观看| 深夜精品福利| 99riav亚洲国产免费| 91大片在线观看| 国产精品免费视频内射| 精品人妻在线不人妻| 色精品久久人妻99蜜桃| 亚洲国产中文字幕在线视频| 成年人午夜在线观看视频| 欧美日本中文国产一区发布| 熟女少妇亚洲综合色aaa.| 乱人伦中国视频| 欧美日韩成人在线一区二区| 成年人午夜在线观看视频| 露出奶头的视频| 亚洲在线自拍视频| 精品视频人人做人人爽| 午夜老司机福利片| 精品一区二区三区视频在线观看免费 | av天堂在线播放| 亚洲一区高清亚洲精品| av欧美777| 国产高清videossex| 亚洲,欧美精品.| av一本久久久久| 男女床上黄色一级片免费看| 一级毛片高清免费大全| 中文字幕另类日韩欧美亚洲嫩草| 久久久久国产精品人妻aⅴ院 | 天天躁夜夜躁狠狠躁躁| 亚洲免费av在线视频| 国产真人三级小视频在线观看| 国产精品欧美亚洲77777| 1024香蕉在线观看| 中文字幕av电影在线播放| 久久国产精品大桥未久av| 欧美 亚洲 国产 日韩一| 国产成人精品无人区| 在线观看66精品国产| 亚洲五月婷婷丁香| 少妇粗大呻吟视频| 国产97色在线日韩免费| 精品人妻在线不人妻| 欧美久久黑人一区二区| 午夜福利免费观看在线| 757午夜福利合集在线观看| a级毛片在线看网站| 激情视频va一区二区三区| 在线观看免费午夜福利视频| 在线视频色国产色| 成人18禁在线播放| 国产欧美亚洲国产| 午夜福利影视在线免费观看| 中亚洲国语对白在线视频| 精品一区二区三区视频在线观看免费 | 久久国产乱子伦精品免费另类| 99热国产这里只有精品6| 天天躁日日躁夜夜躁夜夜| 91成年电影在线观看| 高清欧美精品videossex| 啪啪无遮挡十八禁网站| 热re99久久精品国产66热6| 一进一出抽搐动态| 美女午夜性视频免费| 在线播放国产精品三级| 老汉色av国产亚洲站长工具| 搡老熟女国产l中国老女人| 久久精品国产亚洲av高清一级| 国产精品99久久99久久久不卡| 成年人黄色毛片网站| av线在线观看网站| www.精华液| 人人妻人人爽人人添夜夜欢视频| 9191精品国产免费久久| 大型av网站在线播放| 一区二区三区国产精品乱码| 19禁男女啪啪无遮挡网站| 日日摸夜夜添夜夜添小说| 成在线人永久免费视频| 亚洲国产精品合色在线| 国产色视频综合| 妹子高潮喷水视频| 久久九九热精品免费| 在线十欧美十亚洲十日本专区| 亚洲九九香蕉| 国产亚洲精品一区二区www | 国产在线精品亚洲第一网站| 大香蕉久久网| 欧美激情久久久久久爽电影 | 看片在线看免费视频| 女性被躁到高潮视频| 婷婷成人精品国产| 超碰成人久久| 999久久久国产精品视频| 90打野战视频偷拍视频| 国产熟女午夜一区二区三区| 中出人妻视频一区二区| 一级黄色大片毛片| 高清黄色对白视频在线免费看| 一二三四在线观看免费中文在| 高清黄色对白视频在线免费看| 女同久久另类99精品国产91| 久久九九热精品免费| 欧美成人午夜精品| 免费在线观看黄色视频的| 国产乱人伦免费视频| 在线观看免费日韩欧美大片| 亚洲精品中文字幕在线视频| 免费观看a级毛片全部| 99国产极品粉嫩在线观看| 欧美黄色淫秽网站| 999久久久国产精品视频| 久久久水蜜桃国产精品网| 亚洲精品久久午夜乱码| 露出奶头的视频| 久久久久国产精品人妻aⅴ院 | 久久精品亚洲精品国产色婷小说| 久久亚洲真实| 欧美精品亚洲一区二区| 国产熟女午夜一区二区三区| 妹子高潮喷水视频| 久久午夜亚洲精品久久| 精品亚洲成a人片在线观看| 国产区一区二久久| 一a级毛片在线观看| 丰满迷人的少妇在线观看| 欧美激情 高清一区二区三区| av片东京热男人的天堂| 亚洲中文日韩欧美视频| 一级黄色大片毛片| 久久中文看片网| 午夜亚洲福利在线播放| 黄片大片在线免费观看| 国产精品 欧美亚洲| 九色亚洲精品在线播放| 久久久久久久国产电影| 19禁男女啪啪无遮挡网站| 午夜老司机福利片| 一级毛片高清免费大全| 国产高清视频在线播放一区| 成人免费观看视频高清| 别揉我奶头~嗯~啊~动态视频| 亚洲精品av麻豆狂野| 村上凉子中文字幕在线| 午夜福利在线免费观看网站| 国产麻豆69| 老司机深夜福利视频在线观看| 三上悠亚av全集在线观看| 国产精品.久久久| 色综合欧美亚洲国产小说| 丰满饥渴人妻一区二区三| 亚洲av熟女| 亚洲精品av麻豆狂野| 宅男免费午夜| 成人手机av| 女人爽到高潮嗷嗷叫在线视频| 1024香蕉在线观看| 国产成人av教育| 亚洲精品国产区一区二| 亚洲人成77777在线视频| 中出人妻视频一区二区| 在线观看免费视频日本深夜| 一进一出好大好爽视频| 在线av久久热| 亚洲欧美精品综合一区二区三区| 十八禁网站免费在线| a级毛片在线看网站| av片东京热男人的天堂| 好看av亚洲va欧美ⅴa在| 女人高潮潮喷娇喘18禁视频| 国产成人精品久久二区二区91| 久久久久精品国产欧美久久久| 女人被躁到高潮嗷嗷叫费观| 国产精品九九99| 一区二区三区精品91| 亚洲 欧美一区二区三区| 色在线成人网| 亚洲熟女精品中文字幕| 一本综合久久免费| 欧美日本中文国产一区发布| 国产97色在线日韩免费| 12—13女人毛片做爰片一| 亚洲精品美女久久av网站| 1024视频免费在线观看| 亚洲 欧美一区二区三区| 久久人人97超碰香蕉20202| 久久中文字幕人妻熟女| 免费一级毛片在线播放高清视频 | 欧美在线黄色| 桃红色精品国产亚洲av| 国产在视频线精品| 精品视频人人做人人爽| 中文欧美无线码| 亚洲成人手机| 亚洲av美国av| 中文字幕人妻熟女乱码| 欧美精品啪啪一区二区三区| av欧美777| 操美女的视频在线观看| 久久 成人 亚洲| 涩涩av久久男人的天堂| 国产三级黄色录像| 大码成人一级视频| 亚洲精品国产区一区二| 99国产精品一区二区蜜桃av | 欧美国产精品va在线观看不卡| 国产淫语在线视频| 国产97色在线日韩免费| 国产亚洲欧美精品永久| 中文字幕制服av| 18在线观看网站| 亚洲av成人一区二区三| 国产91精品成人一区二区三区| bbb黄色大片| 韩国精品一区二区三区| 欧美黄色淫秽网站| av中文乱码字幕在线| 久久久国产一区二区| 日韩欧美一区二区三区在线观看 | 在线观看www视频免费|