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

    入流激勵(lì)下可壓縮剪切層中Kelvin-Helmholtz渦的響應(yīng)特性*

    2020-02-18 03:17:56張冬冬譚建國(guó)姚霄
    物理學(xué)報(bào) 2020年2期
    關(guān)鍵詞:大渦渦量流向

    張冬冬 譚建國(guó) 姚霄

    (國(guó)防科技大學(xué),高超聲速?zèng)_壓發(fā)動(dòng)機(jī)技術(shù)重點(diǎn)實(shí)驗(yàn)室,長(zhǎng)沙 410073)

    通過直接數(shù)值求解Navier-Stokes方程,研究了入流激勵(lì)下可壓縮剪切層中Kelvin-Helmholtz (KH)渦結(jié)構(gòu)的響應(yīng)特性,結(jié)果清晰地展示了KH渦的獨(dú)特演化方式.基于流動(dòng)可視化數(shù)據(jù),采用兩點(diǎn)相關(guān)性分析獲得了流場(chǎng)擬序結(jié)構(gòu)的空間尺寸和結(jié)構(gòu)角分布.通過分析不同激勵(lì)頻率下渦結(jié)構(gòu)的動(dòng)態(tài)特性,揭示了入流激勵(lì)下可壓縮剪切層中KH渦結(jié)構(gòu)的獨(dú)特演化機(jī)理.研究結(jié)果表明,低頻入流激勵(lì)(f=5 kHz)下KH渦尺寸在遠(yuǎn)場(chǎng)區(qū)域達(dá)到飽和后呈現(xiàn)鎖頻狀態(tài),KH渦量厚度穩(wěn)定在12—14 mm之間;與自由剪切層渦結(jié)構(gòu)通過配對(duì)合并的方式實(shí)現(xiàn)生長(zhǎng)的機(jī)理不同,低頻入流激勵(lì)下剪切層的發(fā)展是通過中間渦核順時(shí)針吞噬KH不穩(wěn)定波誘導(dǎo)的一串外圍小渦結(jié)構(gòu)來實(shí)現(xiàn)生長(zhǎng).此外,針對(duì)高頻激勵(lì) (f=20 kHz)下的剪切層流動(dòng),研究了渦結(jié)構(gòu)特性和入流激勵(lì)參數(shù)之間的定量關(guān)系,發(fā)現(xiàn)均勻分布渦結(jié)構(gòu)的尺寸近似等于對(duì)流速度與入流激勵(lì)頻率之比.

    1 引 言

    近年來隨著超燃沖壓發(fā)動(dòng)機(jī)技術(shù)以及高速導(dǎo)彈成像制導(dǎo)技術(shù)的發(fā)展,可壓縮剪切層流動(dòng)得到了廣泛的關(guān)注和持續(xù)不斷的研究[1,2].超燃沖壓發(fā)動(dòng)機(jī)有限的燃燒室尺寸內(nèi),超聲速空氣與燃料剪切混合的時(shí)間尺度為毫秒量級(jí),如何在較短的時(shí)間內(nèi)實(shí)現(xiàn)空氣與燃料的充分混合成為其中最為核心的一項(xiàng)關(guān)鍵技術(shù)[3,4].此外,對(duì)于紅外成像制導(dǎo)技術(shù)而言,導(dǎo)彈以高馬赫數(shù)飛行時(shí),其光學(xué)側(cè)窗周圍存在嚴(yán)重的氣動(dòng)光學(xué)效應(yīng),使得目標(biāo)圖像出現(xiàn)偏移、抖動(dòng)以及模糊.這其中為消除氣動(dòng)熱效應(yīng)而采用的制冷劑噴流與外部來流形成的可壓縮剪切層流場(chǎng)結(jié)構(gòu)是產(chǎn)生氣動(dòng)光學(xué)效應(yīng)的一個(gè)重要來源[5].一系列工程應(yīng)用上的需求使得可壓縮剪切層流場(chǎng)結(jié)構(gòu)的調(diào)控變得尤為迫切.在科學(xué)價(jià)值方面,剪切層流動(dòng)存在速度拐點(diǎn),其誘導(dǎo)的Kelvin-Helmholtz (KH)不穩(wěn)定波促成了流場(chǎng)中典型大尺度渦—KH渦的卷起,KH渦的發(fā)展和演化是剪切層實(shí)現(xiàn)增長(zhǎng)的重要方式.此外,剪切層不受壁面干擾,結(jié)構(gòu)相對(duì)簡(jiǎn)單,深入研究其流場(chǎng)中KH渦結(jié)構(gòu)的相互作用及渦結(jié)構(gòu)的生長(zhǎng)特性具有重要的學(xué)術(shù)意義和理論價(jià)值.

    Brown和Roshko[6]于1974年首次在剪切層中發(fā)現(xiàn)KH不穩(wěn)定波誘導(dǎo)的擬序結(jié)構(gòu),且這些結(jié)構(gòu)主導(dǎo)著剪切層的能量交換及動(dòng)力學(xué)行為.隨后Ortwerth和Shine[7]通過實(shí)驗(yàn)發(fā)現(xiàn)在可壓縮狀態(tài)下剪切層中同樣存在著大尺度KH擬序渦結(jié)構(gòu),但是可壓縮效應(yīng)的存在使得KH渦結(jié)構(gòu)的增長(zhǎng)受到顯著抑制.為了衡量壓縮性對(duì)剪切層的影響,Papamoschou和Roshko[8]基于大尺度KH渦結(jié)構(gòu)首次提出了對(duì)流馬赫數(shù)(convective Mach number,Mc)的概念,這一無量綱量后來成為國(guó)內(nèi)外學(xué)者開展可壓縮剪切層穩(wěn)定性分析,混合增強(qiáng)以及流動(dòng)控制研究采用的重要參數(shù).

    針對(duì)剪切層主動(dòng)流動(dòng)控制的研究,最早由Ho和Huang[9]在低速剪切層中采用人工擾動(dòng)的方式分析了流場(chǎng)中KH渦結(jié)構(gòu)的演化及剪切層的增長(zhǎng)情況,其研究表明在流場(chǎng)中加入人工擾動(dòng)可以實(shí)現(xiàn)對(duì)流場(chǎng)結(jié)構(gòu)的控制.隨后McLaughlin等[10]在可壓縮剪切層中引入電火花激勵(lì),通過電火花在流動(dòng)入口處產(chǎn)生高頻溫度擾動(dòng),激發(fā)出流場(chǎng)中的大尺度KH渦結(jié)構(gòu),從而實(shí)現(xiàn)了可壓縮條件下剪切層的流動(dòng)控制,但是受到實(shí)驗(yàn)條件的限制,其研究結(jié)果并沒有觀察到人工激勵(lì)下流場(chǎng)的精細(xì)結(jié)構(gòu)以及KH渦結(jié)構(gòu)的生長(zhǎng)方式.Guo等[11]基于氣動(dòng)光學(xué)研究背景,采用大渦模擬方法研究了脈沖激勵(lì)下可壓縮剪切層的演化機(jī)理,其研究基于渦核位置提取方法,對(duì)流場(chǎng)中渦結(jié)構(gòu)的空間尺寸和瞬時(shí)對(duì)流速度的動(dòng)態(tài)特性進(jìn)行了定量計(jì)算.Freeman和Catrakis[12]采用等離子體激勵(lì)的方式對(duì)剪切層渦結(jié)構(gòu)演化特性進(jìn)行了研究,結(jié)果表明選擇恰當(dāng)?shù)牡入x子體激勵(lì)頻率,可以有效抑制渦結(jié)構(gòu)的生長(zhǎng).在周期性激勵(lì)控制流動(dòng)混合方面,Yu等[13]的研究表明入口人工擾動(dòng)下可壓縮剪切層的增長(zhǎng)率可以提高2—3倍;然而Feng[14]采用大渦模擬方法研究表明,入口人工激勵(lì)對(duì)剪切層演化的影響呈現(xiàn)出先促進(jìn)增長(zhǎng),后抑制增長(zhǎng)的特性.

    由此可見,雖然國(guó)內(nèi)外學(xué)者對(duì)剪切層進(jìn)行控制時(shí)采用的人工激勵(lì)方式有差別,但其本質(zhì)上都是通過人工激勵(lì)來調(diào)控流場(chǎng)中KH擬序渦結(jié)構(gòu)的發(fā)展.在周期性的人工激勵(lì)下,渦結(jié)構(gòu)必然會(huì)呈現(xiàn)出某些特有的規(guī)律,進(jìn)而影響流場(chǎng)的發(fā)展和演化方式.然而過去的研究更多的是對(duì)人工激勵(lì)下可壓縮剪切層中KH渦結(jié)構(gòu)進(jìn)行可視化分析[12],相關(guān)的定量分析如結(jié)構(gòu)空間尺寸,結(jié)構(gòu)角分布等還缺乏深入的研究,進(jìn)而使得對(duì)人工激勵(lì)下KH渦結(jié)構(gòu)特有的響應(yīng)特性和生長(zhǎng)機(jī)理也缺乏足夠的認(rèn)識(shí).此外,對(duì)比Yu等[13]和Feng[14]的研究可以發(fā)現(xiàn),在入口激勵(lì)影響剪切層摻混方面,前人的研究對(duì)于剪切層的增長(zhǎng)特性也尚未達(dá)成共識(shí).

    基于此,本文研究采用入流周期性激勵(lì)的方式對(duì)可壓縮剪切層添加人工擾動(dòng),通過直接數(shù)值求解可壓縮Navier-Stokes方程,借助于KH渦量厚度分析,空間相關(guān)性分析等研究手段對(duì)入流激勵(lì)下剪切層的增長(zhǎng)特性和KH擬序渦結(jié)構(gòu)的空間尺寸,結(jié)構(gòu)角分布進(jìn)行定量研究.通過分析不同激勵(lì)頻率下流場(chǎng)中渦結(jié)構(gòu)的動(dòng)態(tài)演化過程,獲得了KH渦結(jié)構(gòu)演化特性和入流激勵(lì)參數(shù)之間的定量關(guān)系,揭示了入流激勵(lì)下可壓縮剪切層KH渦結(jié)構(gòu)的演化機(jī)理和獨(dú)特生長(zhǎng)特性.

    2 計(jì)算方法和模型

    2.1 數(shù)值計(jì)算方法

    通過直接數(shù)值求解可壓縮Navier-Stokes方程獲得流動(dòng)的瞬時(shí)流場(chǎng),即各種尺度的隨機(jī)運(yùn)動(dòng),從而獲得流場(chǎng)的全部時(shí)空信息,結(jié)合可視化和相關(guān)后處理方法揭示流場(chǎng)結(jié)構(gòu)的發(fā)展演化特性.計(jì)算采用的高精度數(shù)值計(jì)算方法由本課題組自行開發(fā),具有高精度,高分辨率和強(qiáng)穩(wěn)定性的特點(diǎn),已經(jīng)在可壓縮流動(dòng)的研究中得到了一系列應(yīng)用[15,16].限于篇幅,這里僅對(duì)格式的構(gòu)造作簡(jiǎn)要介紹.對(duì)于可壓縮Navier-Stokes方程,對(duì)流項(xiàng)采用具有5階精度的WENO (weighted essentially non-oscillatory) 格式進(jìn)行求解,WENO格式能夠在克服ENO(essentially non-oscillatory) 格式帶來的舍入誤差大及經(jīng)濟(jì)性差的缺陷的同時(shí)保持較好的計(jì)算剛性和精度.黏性項(xiàng)采用6階精度中心差分格式離散,時(shí)間推進(jìn)上采用三階精度具有TVD保持性質(zhì)的Runge-Kutta方法,庫(kù)朗特?cái)?shù)(CFL)選為0.1.相關(guān)格式的具體構(gòu)造方式可以參見[17,18].

    圖1 計(jì)算模型示意圖Fig.1.Schematic of computational model.

    由于本研究并未關(guān)注KH渦轉(zhuǎn)捩形成湍流的過程,因此考慮到計(jì)算成本,采用了二維模型進(jìn)行計(jì)算,如圖1所示.前人的研究亦表明,如果僅僅關(guān)注剪切層中KH不穩(wěn)定及KH渦的演化,采用二維的計(jì)算模型是可行的[19].計(jì)算中取混合區(qū)的流向長(zhǎng)度為L(zhǎng)=600 mm,橫向長(zhǎng)度2H=200 mm,流向長(zhǎng)度和橫向長(zhǎng)度的比值為3.在文獻(xiàn)[19]中,通過設(shè)置流向長(zhǎng)度為橫向長(zhǎng)度的4倍來保證橫向區(qū)域足夠大,從而有效減小上下邊界對(duì)主流的影響.本文取這一比值為3,從而將影響進(jìn)一步降低.

    本次計(jì)算的來流參數(shù)如表1所列.上下來流的速度分別為473 m/s和284 m/s,當(dāng)?shù)芈曀倏赏ㄟ^關(guān)系式獲得

    其中γ為比熱比,γ=1.4.R為氣體常數(shù),取R=287 J·kg—1K—1.由關(guān)系式 Ma=U/a 求得上下兩層超聲速來流的馬赫數(shù)分別為2.15和1.12.

    表1 數(shù)值計(jì)算來流參數(shù)Table 1.Inflow parameters of numerical simulation.

    對(duì)于邊界條件的處理,由于上下兩股來流均為超聲速流動(dòng),因此入口直接給定超聲速入口邊界條件,即只提物理邊界條件,給定上下兩層的來流參數(shù).出口只提數(shù)值邊界條件.上下邊界均給無反射邊界條件,以保證邊界的信息不會(huì)對(duì)剪切層主流產(chǎn)生影響.對(duì)于5階WENO格式,3個(gè)子模板共涉及5個(gè)節(jié)點(diǎn)的信息.對(duì)于邊界部分的節(jié)點(diǎn),需要人為補(bǔ)充3個(gè)鏡像點(diǎn)(ghost cells)來提供足夠的網(wǎng)格點(diǎn)構(gòu)造子模板.

    2.2 算例與網(wǎng)格無關(guān)性驗(yàn)證

    采用Goebel-Dutton[20]經(jīng)典實(shí)驗(yàn)算例對(duì)本計(jì)算采用的數(shù)值格式進(jìn)行驗(yàn)證.為了將計(jì)算結(jié)果和Goebel-Dutton的實(shí)驗(yàn)結(jié)果相對(duì)比,所有條件均與實(shí)驗(yàn)條件相對(duì)應(yīng),Goebel-Dutton實(shí)驗(yàn)來流條件如表2所列.

    表2 Goebel-Dutton超聲速混合層實(shí)驗(yàn)來流參數(shù)Table 2.Inflow parameters of supersonic mixing layer experiments conducted by Goebel and Dutton.

    圖2(a)對(duì)比了湍流發(fā)展到自相似區(qū)域時(shí)流場(chǎng)的時(shí)均速度剖面與Goebel-Dutton的實(shí)驗(yàn)結(jié)果,其中ΔU=U1— U2,U1和U2分別代表上下兩層流體的速度.y*為縱向無量綱坐標(biāo),y*=(y-y0)/b,y0為流場(chǎng)y方向中心坐標(biāo),b為y方向?qū)?yīng)速度值U=U1— 0.1ΔU和 U=U2+0.1ΔU之間的距離.可以發(fā)現(xiàn)時(shí)均速度剖面與實(shí)驗(yàn)結(jié)果吻合的較好.圖2(b)對(duì)計(jì)算獲得的流向脈動(dòng)強(qiáng)度與實(shí)驗(yàn)結(jié)果進(jìn)行了對(duì)比,發(fā)現(xiàn)其結(jié)果同樣具有較好的一致性.通過對(duì)比流場(chǎng)的一階量和二階量驗(yàn)證了本次數(shù)值計(jì)算采用的DNS方法在模擬超聲速混合層流場(chǎng)方面的有效性.

    圖2 數(shù)值與實(shí)驗(yàn)對(duì)比 (a) 時(shí)均速度;(b) 流向湍流強(qiáng)度Fig.2.Comparison between numerical and experimental results:(a) Mean velocity;(b) turbulent intensity.

    此外為了進(jìn)行網(wǎng)格無關(guān)性的驗(yàn)證,本文采用了三套網(wǎng)格進(jìn)行數(shù)值計(jì)算.網(wǎng)格數(shù)量分別為960×240,1440×240 以及 1920×240.網(wǎng)格分布為流向均勻分布,橫向網(wǎng)格采用正切函數(shù)方式進(jìn)行加密.圖3對(duì)比了不同網(wǎng)格條件下混合層的渦量厚度沿流向的變化情況.渦量厚度 δω定義為δω=ΔU/|du/dy|max.可以發(fā)現(xiàn)在 1440×240 以及 1920×240 網(wǎng)格分布下,混合層的渦量厚度沿流向的變化大體一致,而在差網(wǎng)格條件下(960×240),渦量厚度的變化與中等網(wǎng)格及精細(xì)網(wǎng)格條件下相差較大.

    此外,對(duì)于一個(gè)可信的結(jié)果而言,其最小網(wǎng)格尺度不能超過Kolmogorov尺度的幾倍,即最小網(wǎng)格尺度與Kolmogorov尺度至少在一個(gè)量級(jí)上[19].可壓縮剪切層剪切運(yùn)動(dòng)最劇烈的部分發(fā)生在剪切層的核心區(qū)域,對(duì)于本文研究,基于獲得的剪切層中心處的最小Kolmogorov尺度約為0.004 mm.這里 ν和ε 分別為動(dòng)力黏性系數(shù)和湍動(dòng)能耗散系數(shù).本文采用的網(wǎng)格在剪切層核心處的最小尺度為0.03 mm,約為Kolmogorov尺度的7.5倍.因此,綜合考慮計(jì)算精度和計(jì)算耗費(fèi)資源,本次計(jì)算采用1440×240的網(wǎng)格分布進(jìn)行計(jì)算.

    圖3 不同網(wǎng)格條件下混合層的渦量厚度沿流向的變化Fig.3.Vorticity thickness variations versus stream wise direction for different mesh distributions.

    2.3 激勵(lì)方式

    采用連續(xù)激勵(lì)控制的方式來模擬對(duì)可壓縮剪切層流場(chǎng)的周期性擾動(dòng).采用的激勵(lì)控制函數(shù)如下:

    式中,A為入口激勵(lì)幅值,f為激勵(lì)頻率,φ 為位于[— π,π]之間的隨機(jī)相位.G (y)為y方向服從Y~N (0.1,0.0052)的正態(tài)分布函數(shù),作用是來控制入流激勵(lì)的影響區(qū)域,C為常數(shù),其取值依賴于G(y)的參數(shù).此處G(y)其表達(dá)式為

    圖4 連續(xù)控制信號(hào)分布Fig.4.Distribution of input continuous signal.

    3 結(jié)果與分析

    3.1 流場(chǎng)可視化分析

    為了和入流激勵(lì)下剪切層的生長(zhǎng)特性進(jìn)行對(duì)比分析,首先研究了可壓縮自由剪切層流場(chǎng)結(jié)構(gòu)的演化特性.圖5定性展示了自由剪切層中KH渦結(jié)構(gòu)的演化過程.流場(chǎng)區(qū)域采用橫向半高H進(jìn)行無量綱化,渦量值采用來流對(duì)流速度Uc和入口動(dòng)量厚度δ(0)進(jìn)行無量綱化.由于剪切層渦量具有極大值,速度分布具有拐點(diǎn),所以剪切層流動(dòng)對(duì)于KH不穩(wěn)定波的擾動(dòng)是無黏不穩(wěn)定的[20].當(dāng)流場(chǎng)中最不穩(wěn)定模態(tài)幅值達(dá)到最大值時(shí),即在流場(chǎng)下游X/H=0.3處,KH渦結(jié)構(gòu)完成卷起.同時(shí),KH渦結(jié)構(gòu)的生長(zhǎng)是通過配對(duì)合并的方式來進(jìn)行.在時(shí)間間隔為200 μs內(nèi),流場(chǎng)中的渦結(jié)構(gòu)P和Q運(yùn)動(dòng)到P'和Q' 的位置,并且完成了配對(duì)與合并,這種對(duì)并的方式是自由剪切層中KH渦結(jié)構(gòu)的生長(zhǎng)機(jī)理,主導(dǎo)了剪切層的發(fā)展,學(xué)者Olsen和Dutton[22]以及Zhang等[23,24]也證實(shí)了這一生長(zhǎng)機(jī)理的存在.此外,圖5也定性地反映了自由剪切層是通過近似線性增長(zhǎng)的方式來實(shí)現(xiàn)上下兩層流體的摻混.

    圖5 自由剪切層流場(chǎng)渦結(jié)構(gòu)分布Fig.5.Distribution of vortex structures of free shear layers.

    圖6為高頻入流激勵(lì)(f=20 kHz)條件下剪切層渦結(jié)構(gòu)的分布,上下兩幅圖的時(shí)間間隔為20 μs.在高頻激勵(lì)下,剪切層在下游X/H=0.4處完成渦結(jié)構(gòu)的卷起,并且渦結(jié)構(gòu)的尺寸很快達(dá)到飽和,在流場(chǎng)向下游演化過程中,渦結(jié)構(gòu)始終保持飽和狀態(tài),這種流動(dòng)現(xiàn)象的出現(xiàn)說明流場(chǎng)中渦結(jié)構(gòu)與入口擾動(dòng)波發(fā)生了共振,這一現(xiàn)象最早由Ho和Huang[9]在不可壓低速剪切層流動(dòng)中發(fā)現(xiàn),此處在可壓縮剪切層中也發(fā)生了這一共振效應(yīng).此外,在時(shí)間間隔為20 μs范圍內(nèi),渦結(jié)構(gòu)并沒有出現(xiàn)自由剪切層中通過“對(duì)并”來生長(zhǎng)的現(xiàn)象,而是自第一個(gè)渦結(jié)構(gòu)卷起后,即維持其飽和尺寸狀態(tài)往下游發(fā)展.

    圖6 高頻入流激勵(lì)下流場(chǎng)渦結(jié)構(gòu)的分布 (f=20 kHz)Fig.6.Vortex structures distribution under inlet forcing with high frequency (f=20 kHz).

    由于渦結(jié)構(gòu)的尺寸在整個(gè)流場(chǎng)中保持一致,因此有必要對(duì)其進(jìn)行分析,探究渦尺度與入流激勵(lì)之間是否存在定量關(guān)系.由于渦核處壓力存在極小值[14],因此只要獲得流場(chǎng)中的壓力分布即可得到流場(chǎng)的渦結(jié)構(gòu)尺寸.這里定義渦尺寸為兩個(gè)相鄰渦核之間的距離λ.圖7為n=10個(gè)渦核區(qū)域范圍內(nèi)的壓力分布,采樣點(diǎn)的位置均位于橫向Y/H=0.6處.渦尺寸可以用下式獲得:

    其中Xp為10個(gè)渦核之間的無量綱距離,Xp=1.706,通過計(jì)算可得基于流場(chǎng)可視化得到的渦尺寸 λ=18.96 mm.此外,Ghoniem和 Ng[25]研究表明,剪切層流場(chǎng)內(nèi)的特征渦尺度λv可以表達(dá)為

    其中

    式中,Uc為剪切層的對(duì)流速度,a1和a2分別為上下兩層的當(dāng)?shù)芈曀?其值可通過 (1)式的聲速關(guān)系式求得.fv為流動(dòng)的某種特征頻率,這里取fv為入流激勵(lì)頻率,fv=20 kHz.通過計(jì)算可得特征渦尺度λv為19.25 mm,可以發(fā)現(xiàn)基于入流激勵(lì)頻率計(jì)算得到的渦尺度λv與λ相差在2%范圍之內(nèi),因此研究結(jié)果表明均勻分布的渦結(jié)構(gòu)的尺寸近似等于對(duì)流速度與入流激勵(lì)頻率之比.

    圖7 流場(chǎng)渦核之間的壓力分布Fig.7.Pressure distribution between the vortex core in the flow field.

    為了進(jìn)一步驗(yàn)證此結(jié)論的正確性,采用快速傅里葉變換(fast Fourier transform,FFT)提取了剪切層流場(chǎng)下游X/H=1和4處的頻譜分布,采樣點(diǎn)的位置均位于橫向Y/H=0.6處,結(jié)果如圖8所示.與自由剪切層中存在各種諧波擾動(dòng)不同,高頻激勵(lì)下剪切層的流場(chǎng)特征頻率鎖定在20 kHz,這正是入流激勵(lì)施加的頻率.在X/H=1處,流場(chǎng)中仍然存在著少許呈規(guī)律分布的高頻擾動(dòng)波,如圖8(a)中黑色虛線框中所示.隨著流場(chǎng)向下游發(fā)展,在X/H=4處,這種高頻成分已經(jīng)完全消失.其原因在于,入流高頻激勵(lì)在和剪切層固有的KH不穩(wěn)定波相互作用過程中占據(jù)主導(dǎo)地位,使得KH不穩(wěn)定誘導(dǎo)的各種不穩(wěn)定模態(tài)波在來流初始剪切處就開始被抑制.隨著流動(dòng)的發(fā)展,高頻激勵(lì)完全控制著渦結(jié)構(gòu)的發(fā)展并實(shí)現(xiàn)了鎖頻的作用.

    為了進(jìn)一步探究完善渦結(jié)構(gòu)生長(zhǎng)特性與不同入流激勵(lì)條件之間的關(guān)系,有必要分析低頻激勵(lì)下渦結(jié)構(gòu)的運(yùn)動(dòng)特性.圖9展示了f=5 kHz時(shí)剪切層流場(chǎng)的可視化結(jié)果,圖9(a)和圖9(b)時(shí)間間隔為 100 μs.在 100 μs的時(shí)間間隔內(nèi),并沒有出現(xiàn)自由剪切層中通過渦結(jié)構(gòu)配對(duì)合并來實(shí)現(xiàn)剪切層增長(zhǎng)的現(xiàn)象,而是出現(xiàn)了一種新的渦結(jié)構(gòu)生長(zhǎng)方式:在流場(chǎng)下游 X/H<2范圍內(nèi),出現(xiàn)了一串由KH不穩(wěn)定波誘導(dǎo)的渦結(jié)構(gòu)剪切帶,在X/H=2處,入流激勵(lì)誘導(dǎo)的渦核開始出現(xiàn),外圍剪切帶圍繞著渦核旋轉(zhuǎn)并逐漸被渦核吞噬 (engulfment),實(shí)現(xiàn)剪切層渦結(jié)構(gòu)的生長(zhǎng).這種吞噬的生長(zhǎng)機(jī)制使得渦結(jié)構(gòu)在剪切層近場(chǎng)處就能實(shí)現(xiàn)其尺寸的快速增長(zhǎng).吞噬過程完成后,渦結(jié)構(gòu)達(dá)到飽和,在向下游演化過程中其尺度不再變化,剪切層的厚度近似呈現(xiàn)出先線性增長(zhǎng),后保持不變的特性.

    圖8 高頻入流激勵(lì)下頻譜分布 (f=20 kHz) (a) X/H=1;(b) X/H=4Fig.8.Frequency spectrum distribution under inlet forcing with high frequency (f=20 kHz):(a) X/H=1;(b) X/H=4.

    圖9 低頻入流擾動(dòng)下剪切層流場(chǎng)渦結(jié)構(gòu)分布 (f=5 kHz)Fig.9.Vortex structures distribution under inlet forcing with low frequency (f=5 kHz).

    圖10提取了低頻激勵(lì)下流向X/H=1,4,橫向均為Y/H=0.6兩處的流場(chǎng)速度振蕩頻譜分布.在X/H=1處,入流激勵(lì)頻率及其倍頻成分的幅值都很大,基頻5 kHz并沒有完全突顯出來,在流場(chǎng)下游X/H=4處,入流激勵(lì)實(shí)現(xiàn)了對(duì)流場(chǎng)渦結(jié)構(gòu)的控制,渦結(jié)構(gòu)運(yùn)動(dòng)的頻率鎖定在5 kHz.此外,f=20 kHz的高頻激勵(lì)下流動(dòng)近場(chǎng)處 (X/H=1)的高頻擾動(dòng)波成分較少,而f=5 kHz時(shí)流動(dòng)近場(chǎng)(X/H=1)處的高頻擾動(dòng)波成分顯著增多,表明相比于低頻激勵(lì),入流高頻激勵(lì)能夠更快地實(shí)現(xiàn)對(duì)流場(chǎng)結(jié)構(gòu)演化的主導(dǎo)和控制.

    圖10 低頻入流激勵(lì)下頻譜分布 (f=5 kHz) (a) X/H=1;(b) X/H=4Fig.10.Frequency spectrum distribution under inlet forcing with low frequency (f=5 kHz):(a) X/H=1;(b) X/H=4.

    為了更好地揭示低頻激勵(lì)下中間渦核對(duì)于外圍剪切帶的吞噬作用,圖11給出了流場(chǎng)中渦結(jié)構(gòu)的動(dòng)態(tài)演化過程.t=0至t=Tω對(duì)應(yīng)于一個(gè)完整的吞噬過程.在中間高速旋轉(zhuǎn)的大渦M的作用下,t=0時(shí)刻從剪切帶脫落的KH渦N進(jìn)入大渦M的作用區(qū)域,直至t=Tω被完全吞噬,從而實(shí)現(xiàn)了渦結(jié)構(gòu)尺寸的增長(zhǎng).與自由剪切層中KH渦結(jié)構(gòu)通過配對(duì)合并實(shí)現(xiàn)增長(zhǎng)不同,在整個(gè)吞噬過程中,渦N的形狀逐漸變得狹長(zhǎng),而中間渦M的形狀則保持的較好.這表明中間大渦結(jié)構(gòu)對(duì)于整個(gè)吞噬過程起主導(dǎo)作用.此外,另一個(gè)重要現(xiàn)象是中間大渦可以同時(shí)吞噬外圍剪切帶的多個(gè)小渦結(jié)構(gòu).在t=0.2Tω時(shí)刻,渦結(jié)構(gòu)R進(jìn)入了中間大渦的影響區(qū)域,高速旋轉(zhuǎn)的大渦結(jié)構(gòu)使得渦結(jié)構(gòu)N和R在這個(gè)周期內(nèi)均完成了吞噬過程.

    圖11 渦結(jié)構(gòu)的動(dòng)態(tài)吞噬過程Fig.11.Dynamic engulfment process of vortex structures.

    3.2 KH渦增長(zhǎng)特性

    入流激勵(lì)下KH渦結(jié)構(gòu)展現(xiàn)出了獨(dú)特的生長(zhǎng)特性,而與流場(chǎng)中擬序結(jié)構(gòu)密切相關(guān)的剪切層必然也會(huì)呈現(xiàn)出不同的增長(zhǎng)特性,本文采用KH渦量厚度指標(biāo)研究了渦結(jié)構(gòu)生長(zhǎng)對(duì)剪切層厚度的影響.圖12為自由剪切層和入流激勵(lì)下剪切層渦量厚度隨流向距離的分布.在自由剪切層中,渦量厚度呈現(xiàn)出近似線性增長(zhǎng)的特點(diǎn),這與圖5體現(xiàn)出的趨勢(shì)相一致.高頻激勵(lì)下(f=20 kHz),剪切層在經(jīng)歷入流的剪切過程后,其渦量厚度值在渦結(jié)構(gòu)達(dá)到飽和后,穩(wěn)定在3—4 mm之間,表明高頻擾動(dòng)下,無論是剪切層的近場(chǎng)區(qū)域還是遠(yuǎn)場(chǎng)區(qū)域,其混合效果都沒有得到改善.低頻激勵(lì)下(f=5 kHz),在剪切層的近場(chǎng)區(qū)域,剪切層的厚度快速增長(zhǎng),渦核的吞噬作用明顯促進(jìn)了來流的混合.在圖12中藍(lán)色虛線框中的分界處,剪切層渦量厚度達(dá)到飽和.在剪切層的遠(yuǎn)場(chǎng)處,剪切層的渦量厚度穩(wěn)定在12—14 mm之間,入流周期性低頻激勵(lì)下(f=5 kHz)剪切層的演化呈現(xiàn)出先加速增長(zhǎng),后平穩(wěn)過渡的特性.

    圖12 剪切層渦量厚度隨流向距離的分布Fig.12.Vorticity thickness distribution versus streamwise direction.

    3.3 空間相關(guān)性分析

    為了對(duì)可壓縮剪切層的渦結(jié)構(gòu)大小、形狀以及結(jié)構(gòu)角進(jìn)行定量分析,進(jìn)一步明晰入流激勵(lì)下渦結(jié)構(gòu)的動(dòng)力學(xué)特性,本文基于計(jì)算結(jié)果,采用空間兩點(diǎn)相關(guān)性分析的方法,研究了流向速度脈動(dòng)在X-Y平面的空間相關(guān)系數(shù)分布,采用的計(jì)算公式為

    式中,(x0,y0)為相關(guān)性計(jì)算的參考點(diǎn),u'代表流向速度的脈動(dòng)均方根值.Bourdon和Dutton[26]對(duì)兩點(diǎn)相關(guān)性分析參數(shù)的研究表明,當(dāng)樣本數(shù)為500左右時(shí)能夠確保統(tǒng)計(jì)結(jié)果的穩(wěn)定性.

    圖13為自由剪切層取不同流向和豎向參考點(diǎn)(X0,Y0)位置時(shí)的相關(guān)性分布,圖13(a)—(d)的參考點(diǎn)無量綱坐標(biāo)分別為(1.5,0.6),(2.4,0.6),(3.5,0.6),(5.1,0.6).定義中心點(diǎn)的值為1,等值線從中心點(diǎn)向外圍以0.05的間隔遞減,最外圍等值線的數(shù)值為0.8.

    在自由剪切層相關(guān)系數(shù)的等值線分布圖中,等值線的輪廓呈現(xiàn)出飽滿的橢圓形.同時(shí)隨著流向距離的增加,流場(chǎng)相關(guān)區(qū)域的面積增加.Bourdon和Dutton[26]通過實(shí)驗(yàn)發(fā)現(xiàn),等值線代表擬序結(jié)構(gòu)的平均結(jié)構(gòu),因此隨著流場(chǎng)往下游發(fā)展,自由剪切層渦流場(chǎng)結(jié)構(gòu)尺寸近似呈現(xiàn)出線性增長(zhǎng)的趨勢(shì),這與圖12中自由剪切層的渦量厚度沿流向的分布相一致.

    為了進(jìn)行對(duì)比分析,圖14給出了入流低頻激勵(lì)下(f=5 kHz)不同參考點(diǎn)處的空間相關(guān)系數(shù)分布,圖14(a)—(d)的參考點(diǎn)無量綱坐標(biāo)分別為(2.0,0.6),(2.7,0.6),(3.5,0.6),(5.1,0.6).與自由剪切層不同的是,在(2.0,0.6)和(5.1,0.6)處,流場(chǎng)結(jié)構(gòu)呈現(xiàn)出飽滿的橢圓形,而在(2.7,0.6)和(3.5,0.6)處,流場(chǎng)結(jié)構(gòu)則趨于狹長(zhǎng)的橢圓形.對(duì)比圖9可以發(fā)現(xiàn),在流向區(qū)域X/H=2—4之間,中間大渦結(jié)構(gòu)通過吞噬一串外圍剪切帶來實(shí)現(xiàn)渦結(jié)構(gòu)的生長(zhǎng).與自由剪切層相比,在這種新的渦結(jié)構(gòu)生長(zhǎng)機(jī)制下,相關(guān)性系數(shù)分布的橢圓結(jié)構(gòu)特性和結(jié)構(gòu)角有著本質(zhì)的變化.對(duì)比圖14(a)和圖14(d),發(fā)現(xiàn)經(jīng)過中間大渦的吞噬作用后,渦結(jié)構(gòu)的尺寸實(shí)現(xiàn)了增長(zhǎng),剪切層的厚度實(shí)現(xiàn)了增加.

    圖13 自由剪切層中不同流向位置處空間相關(guān)性分布 (a) (1.5,0.6);(b) (2.4,0.6);(c) (3.5,0.6);(d) (5.1,0.6)Fig.13.Spatial correlation distributions of free shear layers in different streamwise positions:(a) (1.5,0.6);(b) (2.4,0.6);(c) (3.5,0.6);(d) (5.1,0.6).

    圖14 低頻激勵(lì)下不同流向位置處空間相關(guān)性分布 (a) (2.0,0.6);(b) (2.7,0.6);(c) (3.5,0.6);(d) (5.1,0.6)Fig.14.Spatial correlation distributions with low frequency forcing in different streamwise positions:(a) (2.0,0.6);(b) (2.7,0.6);(c) (3.5,0.6);(d) (5.1,0.6).

    圖15 擬序結(jié)構(gòu)的大小和結(jié)構(gòu)角示意圖Fig.15.Schematic of vortex size and structure angle.

    Wu等[27]的研究指出,可以采用最小二乘擬合的方法對(duì)相關(guān)性分布的外圍等值線進(jìn)行橢圓擬合,從而獲得渦結(jié)構(gòu)的大小和結(jié)構(gòu)角分布.本文對(duì)外圍0.8等值線進(jìn)行橢圓擬合,如圖15所示.得到的橢圓長(zhǎng)軸2d為相關(guān)性計(jì)算得到的擬序結(jié)構(gòu)大小,α為結(jié)構(gòu)角,α定義為擬合后橢圓長(zhǎng)軸與流向的夾角.

    圖16為通過橢圓擬合得到的自由剪切層和入流低頻激勵(lì)下 (f=5 kHz)剪切層在不同參考點(diǎn)處無量綱結(jié)構(gòu)大小δnon的分布,以自由剪切層在(1.5,0.6)處的結(jié)構(gòu)尺寸為1.可以發(fā)現(xiàn)自由剪切層中,由相關(guān)性分析得到的結(jié)構(gòu)大小沿流向近似呈線性增長(zhǎng)的趨勢(shì),這與圖12得到的自由剪切層渦量厚度的變化趨勢(shì)相一致.而低頻激勵(lì)下,流向下游無量綱位置X/H=2和2.7處結(jié)構(gòu)尺寸的變化是先從1.31快速增長(zhǎng)到3.17左右,然后在X/H>2.7時(shí),流場(chǎng)結(jié)構(gòu)的大小逐步減小,并在遠(yuǎn)場(chǎng)處維持在2.0左右.結(jié)合上文分析的低頻激勵(lì)下渦結(jié)構(gòu)的生長(zhǎng)機(jī)理不難發(fā)現(xiàn),X/H=2.7位置附近正是中間大渦吞噬外圍剪切帶的區(qū)域,吞噬作用是導(dǎo)致渦結(jié)構(gòu)尺寸快速增長(zhǎng)的重要原因.

    圖16 無量綱結(jié)構(gòu)大小沿流向的分布Fig.16.Non-dimensional structure size distribution versus streamwise direction.

    圖17 結(jié)構(gòu)角沿流向的分布Fig.17.Structure angle distribution versus streamwise direction.

    圖17給出了不同流向位置處自由剪切層和低頻激勵(lì) (f=5 kHz)剪切層的結(jié)構(gòu)角分布.自由剪切層中結(jié)構(gòu)角沿流向的變化起伏較小,這是由于無外加擾動(dòng)的自由剪切層完全依靠渦結(jié)構(gòu)的配對(duì)與合并來實(shí)現(xiàn)剪切層的增長(zhǎng),這種近似線性的增長(zhǎng)方式使得布置在剪切層核心發(fā)展區(qū)的參考點(diǎn)結(jié)構(gòu)的旋轉(zhuǎn)特性較為微弱.而低頻激勵(lì)下結(jié)構(gòu)角呈現(xiàn)出不同的變化趨勢(shì).在近場(chǎng)X/H=2處,結(jié)構(gòu)角較小;在下游X/H=2.7以及X/H=3.5處,流場(chǎng)結(jié)構(gòu)表現(xiàn)出強(qiáng)烈的傾斜和旋轉(zhuǎn)特性,結(jié)構(gòu)角分別達(dá)到157°和24°.結(jié)合圖14(b)和圖14(c)中相關(guān)性系數(shù)在(2.7,0.6)和(3.5,0.6)處的分布可以發(fā)現(xiàn),在流場(chǎng)吞噬作用發(fā)生的區(qū)域,外圍剪切帶被中間大渦吞噬使得流場(chǎng)結(jié)構(gòu)的傾斜特征變得明顯.同時(shí),結(jié)構(gòu)角的大小在X/H=2.7,X/H=3.5以及X/H=5.1的變化趨勢(shì)表明,中間大渦對(duì)外圍剪切帶的順時(shí)針吞噬作用是入流低頻激勵(lì)下剪切層實(shí)現(xiàn)增長(zhǎng)的重要機(jī)理.

    4 結(jié) 論

    通過直接數(shù)值求解Navier-Stokes方程,研究了入流激勵(lì)下可壓縮剪切層中KH渦結(jié)構(gòu)的響應(yīng)特性.計(jì)算結(jié)果清晰地展示了自由剪切層和入流激勵(lì)下剪切層的流場(chǎng)擬序結(jié)構(gòu)及其動(dòng)態(tài)演化過程.采用KH渦量厚度指標(biāo)定量分析了剪切層的增長(zhǎng)特性.基于流動(dòng)可視化結(jié)果,結(jié)合空間相關(guān)性分析探究了流場(chǎng)結(jié)構(gòu)尺寸和結(jié)構(gòu)角分布,揭示了入流激勵(lì)下剪切層的混合特性及其渦結(jié)構(gòu)的獨(dú)特生長(zhǎng)機(jī)理.

    自由剪切層中,流場(chǎng)KH渦結(jié)構(gòu)通過配對(duì)與合并的方式實(shí)現(xiàn)渦結(jié)構(gòu)的生長(zhǎng),上下兩層流體是通過近似線性增長(zhǎng)的方式來實(shí)現(xiàn)摻混.在入流高頻激勵(lì)下(f=20 kHz),渦結(jié)構(gòu)在下游X/H=0.3處卷起后很快達(dá)到飽和狀態(tài),并且維持其飽和尺度往下游發(fā)展,剪切層的厚度值穩(wěn)定在3—4 mm之間,這種由渦結(jié)構(gòu)與入流擾動(dòng)波發(fā)生共振效應(yīng)導(dǎo)致的現(xiàn)象最早在Ho和Huang[9]研究的不可壓低速剪切層中發(fā)現(xiàn).此外高頻激勵(lì)下,均勻分布飽和渦結(jié)構(gòu)的尺寸近似等于對(duì)流速度與入流激勵(lì)頻率之比.同時(shí)在流場(chǎng)下游,入流高頻擾動(dòng)占據(jù)主導(dǎo)地位,剪切層的流場(chǎng)特征頻率鎖定在20 kHz.

    入流低頻激勵(lì)下(f=5 kHz),流場(chǎng)可視化結(jié)果表明擬序結(jié)構(gòu)的生長(zhǎng)是通過中間大渦的吞噬作用來實(shí)現(xiàn),這種吞噬的生長(zhǎng)機(jī)制使得在剪切層近場(chǎng)處就能實(shí)現(xiàn)渦結(jié)構(gòu)尺寸的快速增長(zhǎng).空間相關(guān)性分析表明,與自由剪切層中等值線輪廓呈現(xiàn)出飽滿的橢圓形不同,在吞噬作用發(fā)生的區(qū)域 (2.7,0.6)和(3.5,0.6),等值線輪廓趨于狹長(zhǎng)的橢圓形.同時(shí)結(jié)構(gòu)角分析表明,在下游X/H=2.7以及X/H=3.5處,低頻激勵(lì)下流場(chǎng)結(jié)構(gòu)表現(xiàn)出強(qiáng)烈的傾斜和旋轉(zhuǎn)特性,結(jié)構(gòu)角分別達(dá)到157°和24°.結(jié)構(gòu)角在流向的變化趨勢(shì)表明,中間大渦對(duì)于外圍剪切帶的順時(shí)針吞噬作用是入流低頻激勵(lì)下剪切層實(shí)現(xiàn)增長(zhǎng)的重要機(jī)理.

    猜你喜歡
    大渦渦量流向
    含沙空化對(duì)軸流泵內(nèi)渦量分布的影響
    小溪?。×飨蜻h(yuǎn)方
    井岡教育(2020年6期)2020-12-14 03:04:42
    基于壁面射流的下?lián)舯┝鞣欠€(wěn)態(tài)風(fēng)場(chǎng)大渦模擬
    自由表面渦流動(dòng)現(xiàn)象的數(shù)值模擬
    軸流風(fēng)機(jī)葉尖泄漏流動(dòng)的大渦模擬
    十大漲幅、換手、振副、資金流向
    基于大渦模擬的旋風(fēng)分離器錐體結(jié)構(gòu)影響研究
    流向逆轉(zhuǎn)的啟示
    航態(tài)對(duì)大型船舶甲板氣流場(chǎng)的影響
    The application of numerical simulation of delta wing with blunt leading edge using RANS/LES hybrid method
    一区二区三区激情视频| 亚洲精华国产精华液的使用体验| 亚洲精品aⅴ在线观看| 日韩精品有码人妻一区| 亚洲精品久久午夜乱码| 啦啦啦啦在线视频资源| 黄片小视频在线播放| 边亲边吃奶的免费视频| 97在线视频观看| 免费日韩欧美在线观看| 天天影视国产精品| 交换朋友夫妻互换小说| 国产一区亚洲一区在线观看| a级毛片黄视频| av在线老鸭窝| 99精国产麻豆久久婷婷| 七月丁香在线播放| 美国免费a级毛片| 国产成人欧美| 国产精品.久久久| 最近最新中文字幕大全免费视频 | 国产精品久久久久久精品电影小说| 国产精品久久久久久久久免| 女性生殖器流出的白浆| 少妇人妻久久综合中文| 最近手机中文字幕大全| 亚洲三级黄色毛片| 母亲3免费完整高清在线观看 | 久久久久久人妻| 国产伦理片在线播放av一区| 制服诱惑二区| 精品少妇内射三级| 精品亚洲成国产av| 肉色欧美久久久久久久蜜桃| 久久久国产一区二区| 2021少妇久久久久久久久久久| 高清在线视频一区二区三区| 永久网站在线| 大片电影免费在线观看免费| av女优亚洲男人天堂| 欧美精品高潮呻吟av久久| 美女脱内裤让男人舔精品视频| 搡女人真爽免费视频火全软件| 天美传媒精品一区二区| 久久久久久久大尺度免费视频| 欧美+日韩+精品| 免费人妻精品一区二区三区视频| 久久青草综合色| 亚洲国产欧美日韩在线播放| 欧美人与性动交α欧美软件| 一级爰片在线观看| 久久久久久久久久久久大奶| 国产一区二区 视频在线| 建设人人有责人人尽责人人享有的| 国产又爽黄色视频| 中文字幕另类日韩欧美亚洲嫩草| 亚洲一级一片aⅴ在线观看| 人妻人人澡人人爽人人| 国产成人精品一,二区| 免费女性裸体啪啪无遮挡网站| 桃花免费在线播放| 亚洲av免费高清在线观看| 国产日韩欧美在线精品| 亚洲精品自拍成人| 精品第一国产精品| 久久久久久人人人人人| 成人亚洲欧美一区二区av| 日韩制服骚丝袜av| 精品国产国语对白av| 黑人巨大精品欧美一区二区蜜桃| 欧美日韩一区二区视频在线观看视频在线| www.精华液| 黄色视频在线播放观看不卡| 人妻一区二区av| 两性夫妻黄色片| 人人澡人人妻人| 蜜桃国产av成人99| 国产高清不卡午夜福利| 可以免费在线观看a视频的电影网站 | 国产亚洲午夜精品一区二区久久| 亚洲欧美一区二区三区国产| 99re6热这里在线精品视频| 亚洲内射少妇av| 多毛熟女@视频| 只有这里有精品99| 午夜福利在线免费观看网站| 午夜福利网站1000一区二区三区| 欧美bdsm另类| 国产人伦9x9x在线观看 | 久久鲁丝午夜福利片| 精品福利永久在线观看| 性少妇av在线| 国产精品欧美亚洲77777| 精品少妇内射三级| 天堂俺去俺来也www色官网| 亚洲av在线观看美女高潮| 午夜福利乱码中文字幕| 国产成人精品福利久久| 成人毛片60女人毛片免费| 美女主播在线视频| 免费久久久久久久精品成人欧美视频| 看免费av毛片| 亚洲欧美一区二区三区国产| 国产日韩欧美亚洲二区| 亚洲在久久综合| 男女免费视频国产| 精品国产乱码久久久久久小说| 欧美国产精品一级二级三级| 卡戴珊不雅视频在线播放| 99热全是精品| 男人添女人高潮全过程视频| 人成视频在线观看免费观看| 欧美日韩视频高清一区二区三区二| 国产精品久久久久久久久免| 免费高清在线观看日韩| 精品午夜福利在线看| 亚洲国产色片| 99re6热这里在线精品视频| 免费播放大片免费观看视频在线观看| 在线 av 中文字幕| 99国产精品免费福利视频| 国产极品粉嫩免费观看在线| 亚洲欧美中文字幕日韩二区| 亚洲精品美女久久久久99蜜臀 | 欧美亚洲 丝袜 人妻 在线| 青青草视频在线视频观看| 天天躁夜夜躁狠狠躁躁| 香蕉丝袜av| 有码 亚洲区| 国产精品一二三区在线看| www.av在线官网国产| 99热网站在线观看| 亚洲欧美中文字幕日韩二区| 欧美激情极品国产一区二区三区| 欧美日韩亚洲国产一区二区在线观看 | 大香蕉久久成人网| 看免费av毛片| 国产无遮挡羞羞视频在线观看| 国产色婷婷99| 老汉色av国产亚洲站长工具| 男女免费视频国产| 最新中文字幕久久久久| 精品一区在线观看国产| 一级毛片我不卡| 亚洲国产欧美日韩在线播放| 欧美老熟妇乱子伦牲交| 精品福利永久在线观看| 免费不卡的大黄色大毛片视频在线观看| 国产av国产精品国产| 青春草国产在线视频| 一区二区日韩欧美中文字幕| 免费观看av网站的网址| 搡女人真爽免费视频火全软件| 免费观看性生交大片5| 国产精品一二三区在线看| 亚洲视频免费观看视频| 日韩av不卡免费在线播放| 亚洲内射少妇av| 欧美日韩av久久| 中文字幕制服av| 国精品久久久久久国模美| 日韩免费高清中文字幕av| 中文字幕亚洲精品专区| 亚洲国产看品久久| 欧美日韩综合久久久久久| 免费播放大片免费观看视频在线观看| 在线观看美女被高潮喷水网站| 99热全是精品| 久久毛片免费看一区二区三区| av网站免费在线观看视频| 久久久久久久国产电影| 18+在线观看网站| 电影成人av| 午夜久久久在线观看| www日本在线高清视频| 国产日韩欧美亚洲二区| 久久精品人人爽人人爽视色| 黄片无遮挡物在线观看| 激情五月婷婷亚洲| 97在线人人人人妻| 两性夫妻黄色片| 天天躁日日躁夜夜躁夜夜| 日韩欧美一区视频在线观看| 精品久久久久久电影网| 99热全是精品| 亚洲av日韩在线播放| 9色porny在线观看| 欧美成人午夜免费资源| 超色免费av| 丰满少妇做爰视频| 中文字幕人妻丝袜一区二区 | 国产高清国产精品国产三级| 国产精品嫩草影院av在线观看| 91久久精品国产一区二区三区| 日韩av免费高清视频| 高清欧美精品videossex| 欧美成人午夜免费资源| 亚洲久久久国产精品| 天天操日日干夜夜撸| 国产 一区精品| 久久人人爽av亚洲精品天堂| 少妇熟女欧美另类| 日韩在线高清观看一区二区三区| 久久久亚洲精品成人影院| 国产精品三级大全| 97人妻天天添夜夜摸| 国产午夜精品一二区理论片| 亚洲精品aⅴ在线观看| 亚洲精品国产一区二区精华液| 成年人午夜在线观看视频| 国产xxxxx性猛交| 午夜福利网站1000一区二区三区| 国产激情久久老熟女| 晚上一个人看的免费电影| 男女下面插进去视频免费观看| 啦啦啦中文免费视频观看日本| 久久国产精品男人的天堂亚洲| 国产精品偷伦视频观看了| 春色校园在线视频观看| 性高湖久久久久久久久免费观看| 日本猛色少妇xxxxx猛交久久| 女的被弄到高潮叫床怎么办| 亚洲av国产av综合av卡| 99久国产av精品国产电影| 亚洲精品美女久久久久99蜜臀 | av网站免费在线观看视频| 国产成人免费无遮挡视频| 午夜福利乱码中文字幕| 少妇的丰满在线观看| 亚洲精品在线美女| 精品人妻偷拍中文字幕| 9191精品国产免费久久| 欧美激情极品国产一区二区三区| 亚洲国产精品一区三区| 80岁老熟妇乱子伦牲交| 亚洲国产欧美在线一区| 国产亚洲精品第一综合不卡| 看非洲黑人一级黄片| 九色亚洲精品在线播放| 国产免费福利视频在线观看| 午夜老司机福利剧场| 尾随美女入室| 国产片内射在线| 国产日韩欧美亚洲二区| 欧美日韩精品成人综合77777| 人妻少妇偷人精品九色| 国产成人精品久久久久久| 欧美日韩视频高清一区二区三区二| 精品国产一区二区三区四区第35| 亚洲欧美成人综合另类久久久| 亚洲欧美日韩另类电影网站| 亚洲av欧美aⅴ国产| 嫩草影院入口| 亚洲国产日韩一区二区| 成人漫画全彩无遮挡| 欧美av亚洲av综合av国产av | 人人妻人人添人人爽欧美一区卜| 国产不卡av网站在线观看| 成人国语在线视频| 丝瓜视频免费看黄片| 丝袜人妻中文字幕| 狂野欧美激情性bbbbbb| 欧美xxⅹ黑人| 好男人视频免费观看在线| 久久亚洲国产成人精品v| 免费不卡的大黄色大毛片视频在线观看| 999精品在线视频| 精品人妻偷拍中文字幕| 亚洲国产精品成人久久小说| 中文字幕人妻丝袜制服| 精品国产国语对白av| 99国产精品免费福利视频| 男女边吃奶边做爰视频| 日本欧美视频一区| 一本—道久久a久久精品蜜桃钙片| 晚上一个人看的免费电影| 亚洲激情五月婷婷啪啪| 97人妻天天添夜夜摸| 国产精品av久久久久免费| 午夜免费鲁丝| 欧美+日韩+精品| 国产精品三级大全| 国产精品香港三级国产av潘金莲 | 免费观看a级毛片全部| 一二三四在线观看免费中文在| 在线天堂最新版资源| av国产精品久久久久影院| 久久久欧美国产精品| 国产男人的电影天堂91| 在线观看一区二区三区激情| 成人国产av品久久久| videos熟女内射| 亚洲成色77777| 爱豆传媒免费全集在线观看| 亚洲欧美日韩另类电影网站| 国产精品二区激情视频| 国产成人免费观看mmmm| 久久久久精品性色| 国产在线视频一区二区| 777米奇影视久久| 国产乱来视频区| 国产av码专区亚洲av| 国产日韩一区二区三区精品不卡| 午夜精品国产一区二区电影| 高清黄色对白视频在线免费看| 搡老乐熟女国产| 久久人人97超碰香蕉20202| 欧美日韩综合久久久久久| xxx大片免费视频| 麻豆精品久久久久久蜜桃| 美女高潮到喷水免费观看| 十分钟在线观看高清视频www| 亚洲欧美精品综合一区二区三区 | 国产精品国产三级专区第一集| 水蜜桃什么品种好| 国产精品二区激情视频| av片东京热男人的天堂| 一二三四在线观看免费中文在| 美女高潮到喷水免费观看| 精品午夜福利在线看| 韩国精品一区二区三区| 国产精品秋霞免费鲁丝片| 亚洲欧美精品综合一区二区三区 | 丰满少妇做爰视频| 人体艺术视频欧美日本| 午夜老司机福利剧场| 大香蕉久久网| 91精品国产国语对白视频| 91精品三级在线观看| 亚洲欧洲精品一区二区精品久久久 | 老汉色∧v一级毛片| 国产黄色免费在线视频| 少妇熟女欧美另类| 免费高清在线观看视频在线观看| 日韩一区二区三区影片| 国产极品天堂在线| 最近的中文字幕免费完整| 国产福利在线免费观看视频| 亚洲精品日本国产第一区| 欧美激情高清一区二区三区 | 极品少妇高潮喷水抽搐| 久久久久久伊人网av| 狠狠精品人妻久久久久久综合| 制服诱惑二区| 晚上一个人看的免费电影| 中文欧美无线码| 日韩一区二区三区影片| 91国产中文字幕| 秋霞在线观看毛片| 18在线观看网站| 久久久久国产精品人妻一区二区| 男女边摸边吃奶| 下体分泌物呈黄色| 成年美女黄网站色视频大全免费| 91久久精品国产一区二区三区| 久热久热在线精品观看| 黄频高清免费视频| 97在线视频观看| 国产一级毛片在线| 亚洲av欧美aⅴ国产| 免费观看性生交大片5| 新久久久久国产一级毛片| 1024视频免费在线观看| 国产免费一区二区三区四区乱码| 久久这里有精品视频免费| 中文字幕人妻熟女乱码| a 毛片基地| 伦理电影大哥的女人| 免费av中文字幕在线| 久久精品国产亚洲av天美| 久久久亚洲精品成人影院| 亚洲欧美一区二区三区久久| 免费看不卡的av| 女人精品久久久久毛片| 国产精品久久久久久精品古装| 寂寞人妻少妇视频99o| av国产精品久久久久影院| 欧美变态另类bdsm刘玥| 亚洲av电影在线观看一区二区三区| a级片在线免费高清观看视频| 91久久精品国产一区二区三区| 热99久久久久精品小说推荐| 日韩制服骚丝袜av| 一本大道久久a久久精品| 肉色欧美久久久久久久蜜桃| 女性被躁到高潮视频| 满18在线观看网站| 一二三四在线观看免费中文在| 亚洲精品视频女| 少妇 在线观看| 18在线观看网站| 亚洲三区欧美一区| 亚洲天堂av无毛| 午夜精品国产一区二区电影| av片东京热男人的天堂| 欧美在线黄色| 国产人伦9x9x在线观看 | 欧美日韩视频高清一区二区三区二| 高清视频免费观看一区二区| 亚洲精品av麻豆狂野| 亚洲中文av在线| 亚洲国产欧美网| 国产av精品麻豆| 一本大道久久a久久精品| 亚洲欧美精品自产自拍| 最近2019中文字幕mv第一页| 街头女战士在线观看网站| 亚洲国产欧美日韩在线播放| 极品人妻少妇av视频| 久久久a久久爽久久v久久| 国产男人的电影天堂91| 国产又色又爽无遮挡免| 久久久久久久精品精品| 国产高清不卡午夜福利| 亚洲精品日本国产第一区| 久久久亚洲精品成人影院| 一边摸一边做爽爽视频免费| 最黄视频免费看| 亚洲精品一区蜜桃| 午夜激情久久久久久久| 熟女av电影| 在线观看美女被高潮喷水网站| www.av在线官网国产| 国产免费一区二区三区四区乱码| 男女免费视频国产| 韩国av在线不卡| 交换朋友夫妻互换小说| 性色avwww在线观看| 亚洲av日韩在线播放| 777久久人妻少妇嫩草av网站| 最新的欧美精品一区二区| 啦啦啦啦在线视频资源| 国产一区二区在线观看av| 国产精品麻豆人妻色哟哟久久| 自拍欧美九色日韩亚洲蝌蚪91| 精品一品国产午夜福利视频| 久久久久久久久免费视频了| 亚洲经典国产精华液单| 欧美国产精品va在线观看不卡| 亚洲精品aⅴ在线观看| 亚洲五月色婷婷综合| 如日韩欧美国产精品一区二区三区| 18禁观看日本| 大码成人一级视频| 亚洲欧美成人精品一区二区| 国产免费视频播放在线视频| 一级毛片黄色毛片免费观看视频| 狂野欧美激情性bbbbbb| 美女主播在线视频| tube8黄色片| 欧美少妇被猛烈插入视频| 高清视频免费观看一区二区| 99九九在线精品视频| 国产成人精品无人区| 麻豆精品久久久久久蜜桃| 久久精品国产自在天天线| 91成人精品电影| www.精华液| 欧美人与善性xxx| 午夜日韩欧美国产| 亚洲av欧美aⅴ国产| 久久久精品区二区三区| 日韩制服骚丝袜av| 在线观看免费视频网站a站| 久久综合国产亚洲精品| 久久国产精品大桥未久av| 久久久久网色| 欧美日韩视频精品一区| 女人高潮潮喷娇喘18禁视频| 国产亚洲午夜精品一区二区久久| 久久久久精品人妻al黑| 高清欧美精品videossex| 91午夜精品亚洲一区二区三区| 国产精品久久久久久精品古装| 精品一区二区三区四区五区乱码 | 亚洲国产精品成人久久小说| 国产精品一国产av| 国产亚洲av片在线观看秒播厂| 高清av免费在线| 亚洲,欧美精品.| 久久韩国三级中文字幕| 波野结衣二区三区在线| 水蜜桃什么品种好| 久久精品国产自在天天线| 少妇精品久久久久久久| 美女xxoo啪啪120秒动态图| 五月开心婷婷网| 五月天丁香电影| 一区二区三区精品91| 99久久人妻综合| 日本猛色少妇xxxxx猛交久久| 91精品三级在线观看| 国语对白做爰xxxⅹ性视频网站| 天天躁夜夜躁狠狠久久av| 老汉色∧v一级毛片| 十分钟在线观看高清视频www| 少妇的逼水好多| 久久97久久精品| 午夜福利在线观看免费完整高清在| 少妇 在线观看| 日韩精品有码人妻一区| 久久婷婷青草| 日日摸夜夜添夜夜爱| 久久精品国产亚洲av高清一级| 久久久精品区二区三区| 中国三级夫妇交换| 亚洲国产精品999| videossex国产| 男女边摸边吃奶| 伊人久久大香线蕉亚洲五| 日韩av免费高清视频| 国产在线视频一区二区| 亚洲精品日本国产第一区| 如日韩欧美国产精品一区二区三区| 国产野战对白在线观看| 成人漫画全彩无遮挡| 色网站视频免费| 亚洲欧美日韩另类电影网站| 极品人妻少妇av视频| 80岁老熟妇乱子伦牲交| 少妇猛男粗大的猛烈进出视频| 成人毛片60女人毛片免费| 国产精品久久久久久久久免| 国产成人精品婷婷| 在线天堂中文资源库| 在线 av 中文字幕| 中文字幕人妻熟女乱码| 日韩精品免费视频一区二区三区| 啦啦啦视频在线资源免费观看| 在线看a的网站| 黄色 视频免费看| 精品午夜福利在线看| 久久精品国产亚洲av涩爱| 婷婷色综合www| 街头女战士在线观看网站| 免费黄频网站在线观看国产| 欧美日韩视频高清一区二区三区二| 国产有黄有色有爽视频| 免费看av在线观看网站| 自拍欧美九色日韩亚洲蝌蚪91| 国产在线视频一区二区| 亚洲国产色片| 亚洲欧美成人综合另类久久久| 久久久精品免费免费高清| 亚洲成人手机| 婷婷色麻豆天堂久久| 老熟女久久久| 亚洲av电影在线观看一区二区三区| 叶爱在线成人免费视频播放| 热99国产精品久久久久久7| 亚洲国产精品一区二区三区在线| 亚洲,一卡二卡三卡| 亚洲成人手机| 国产高清不卡午夜福利| 国产在视频线精品| 高清黄色对白视频在线免费看| 97人妻天天添夜夜摸| 久久免费观看电影| 国产黄色视频一区二区在线观看| 日韩人妻精品一区2区三区| 国产精品久久久久久精品古装| 在线观看三级黄色| 岛国毛片在线播放| 久久精品夜色国产| 国产97色在线日韩免费| 一本一本久久a久久精品综合妖精 国产伦在线观看视频一区 | 国产在线一区二区三区精| 青春草亚洲视频在线观看| 黄频高清免费视频| 18禁动态无遮挡网站| 黑人猛操日本美女一级片| 99精国产麻豆久久婷婷| 婷婷色av中文字幕| 亚洲内射少妇av| 午夜福利在线免费观看网站| 韩国高清视频一区二区三区| 国产男人的电影天堂91| 日韩熟女老妇一区二区性免费视频| 免费黄色在线免费观看| 一本色道久久久久久精品综合| 捣出白浆h1v1| 97人妻天天添夜夜摸| 菩萨蛮人人尽说江南好唐韦庄| 亚洲,欧美精品.| 国产精品嫩草影院av在线观看| 久久久久精品久久久久真实原创| 亚洲人成网站在线观看播放| 精品久久久精品久久久| 一区在线观看完整版| 久热这里只有精品99| 一区二区三区乱码不卡18| 亚洲熟女精品中文字幕| 久久久精品免费免费高清| av电影中文网址| 99香蕉大伊视频| 国产成人午夜福利电影在线观看| 一本大道久久a久久精品| 少妇被粗大猛烈的视频| kizo精华| 亚洲精品日本国产第一区| 亚洲国产av新网站| 午夜福利,免费看| 亚洲欧美中文字幕日韩二区| 亚洲综合色惰| 国产在线视频一区二区| 侵犯人妻中文字幕一二三四区| av免费在线看不卡| 免费大片黄手机在线观看| 99热全是精品| 亚洲综合色网址| 夫妻性生交免费视频一级片| 婷婷成人精品国产| www.熟女人妻精品国产| 国产av码专区亚洲av| 国产欧美日韩一区二区三区在线|