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

    非線性擾動方程研究平面混合層發(fā)聲機制

    2016-04-01 07:26:57張程遠方一紅
    空氣動力學學報 2016年5期
    關(guān)鍵詞:特征函數(shù)等值線馬赫數(shù)

    張程遠,方一紅

    (天津大學機械工程學院力學系,天津300072)

    非線性擾動方程研究平面混合層發(fā)聲機制

    張程遠,方一紅*

    (天津大學機械工程學院力學系,天津300072)

    在流動過程中由脈動引起的噪聲問題在自然界和工業(yè)界中廣泛存在,研究流動過程中不穩(wěn)定波的發(fā)聲機制對于理解、預(yù)測并最終控制氣動噪聲有著重要意義。本文以二維Blasius相似性解作為基本流場,以超聲速混合層作為研究對象,采用非線性擾動方程(NLDE)研究超聲速混合層不穩(wěn)定波的近場動力學特性與遠場聲輻射之間的內(nèi)在聯(lián)系。針對不同類型的對流馬赫數(shù),分別研究對流馬赫數(shù)Mc=0.5(<1)和Mc=1.2(≥1)兩種情況下的擾動發(fā)聲機制,結(jié)果表明:當對流馬赫數(shù)Mc<1時,流場中的發(fā)聲機制主要由大尺度結(jié)構(gòu)的渦條發(fā)聲,且基本波與亞諧波之間的非線性作用能增強輻射強度;當對流馬赫數(shù)Mc≥1時,根據(jù)擾動相速度是相對于上層還是下層自由流速度為超聲速,可以進一步分為快慢兩種模態(tài),分別對快慢兩種模態(tài)以及其相互作用模態(tài)進行了數(shù)值模擬研究,計算得到流場中的發(fā)聲機制是以馬赫波形式輻射向遠場,即馬赫波輻射,其輻角的計算結(jié)果與理論值相符。

    非線性擾動方程;超聲速混合層;對流馬赫數(shù);渦條發(fā)聲;快慢模態(tài);馬赫波輻射

    0 引言

    流動中由脈動導(dǎo)致的噪聲問題在航空和自然工業(yè)中廣泛存在,例如發(fā)動機射流噪聲、機翼和旋翼噪聲、燃氣噪聲等,伴隨著航空航天飛行器的發(fā)展,由非定常氣動力產(chǎn)生的氣動噪聲對環(huán)境造成了嚴重地污染,氣動噪聲問題開始得到大家的重視,在此期間,相關(guān)的科研人員開始著手氣動噪聲問題的研究,隨之形成了氣動聲學這門學科[1]。流動中的強剪切區(qū)域一般具有復(fù)雜的動力學特性,最終將導(dǎo)致大范圍內(nèi)的時間和空間尺度的湍流結(jié)構(gòu)。湍流噪聲是流動中的雷諾應(yīng)力輻射的噪聲,從穩(wěn)定性理論的角度來看,對于可壓縮流動來說,流體中的不穩(wěn)定波與聲音的產(chǎn)生存在著直接的聯(lián)系,因此研究流動過程中不穩(wěn)定波的動力學特性對于理解、預(yù)測并最終控制氣動噪聲是十分必要的。這也引起了學術(shù)界和工業(yè)界的日益重視[2]。

    流動中不穩(wěn)定波與聲輻射之間是如何聯(lián)系在一起的?近年來已經(jīng)有許多學者開始研究二者之間的內(nèi)在聯(lián)系。Colonius等[3]和Bogey等[4]采用直接數(shù)值模擬研究二維平面混合層聲輻射問題,Mitchell等[5]利用數(shù)值模擬的方法研究對稱渦結(jié)構(gòu)產(chǎn)生的聲場,李曉東和高軍輝[6]使用線性Euler方程研究二維平行剪切層聲波產(chǎn)生和輻射,邢超和羅紀生[7]使用直接數(shù)值模擬研究聲波與剪切層的相互作用,Tam等[8-9]等使用漸進展開的方法在無粘和局部平行流假設(shè)的前提下研究剪切層中不穩(wěn)定波與聲輻射之間的關(guān)系,Wu[10-11]使用漸進分析的方法結(jié)合三層結(jié)構(gòu)理論,將流場分為粘性亞層、主層和勢流層,研究了剪切層中的聲輻射問題,Bertolotti等[12]使用線性拋物化穩(wěn)定性方程(Parabolized Stability Equations,PSE)研究射流聲輻射問題,這些研究主要還是從單個擾動波的演化性質(zhì)著手,入口加入的擾動初始幅值較小,對不穩(wěn)定波的非線性作用考慮較少。

    相比于線性不穩(wěn)定波的特性,對于非線性不穩(wěn)定波的近場動力學特性和遠聲場輻射之間的內(nèi)在聯(lián)系方面的研究較少,Cheung等[13]使用線性和非線性PSE研究二維可壓縮剪切層聲輻射問題,在亞聲速剪切層中出現(xiàn)了渦條結(jié)構(gòu),但是并沒有很好地捕捉到遠聲場的特性,可能的原因是PSE在法向上采用的是齊次邊界條件,而且PSE方法假設(shè)沿流向是緩慢變化,但混合層的增長速度較快,擾動在向下游演化過程中會受到強非平行性影響,這些都會影響到流場和聲場信息的捕捉。

    為了進一步研究剪切層中不穩(wěn)定波與聲輻射之間的內(nèi)在關(guān)聯(lián),本文采用非線性擾動方程(Nonlinear Disturbance Equations,NLDE)研究二維超聲速可壓縮混合層的聲輻射問題。同時考慮了非線性相互作用和非平行性的影響,數(shù)值模擬了對流馬赫數(shù)Mc<1和Mc≥1兩種情況下不穩(wěn)定波的擾動演化,給出了不同情況下擾動演化的計算結(jié)果,并在此基礎(chǔ)上研究擾動波的發(fā)聲機制。

    1 控制方程

    對當前混合層的穩(wěn)定性分析,采用非線性擾動方程的分析方法,推導(dǎo)擾動方程的基本思想是:將二維可壓縮N-S方程中的守恒量分解為基本流與擾動量之和的形式并帶到N-S方程中,去掉滿足定常流N-S方程的基本流項,經(jīng)過化簡可得到擾動方程,方程形式如下:

    粘性擾動項:

    2 差分格式和邊界條件

    計算過程中,粘性項采用六階中心差分格式,具體形式如下:

    對流項采用通量分裂,對正通量E+和負通量E-在空間上采用五階弱迎風緊致格式:

    時間項采用非TVD特性的四階四步Runge-Kutta格式:

    對于每一個邊界,都有四個特征波穿過,其特征波的傳播方向由各自的特征值正負號決定。對于超聲速流,四個特征值都為正,其對應(yīng)的特征波都向右傳;對于亞聲速流,有一個特征值為負,其它特征值為正,即三個波向右,一個波向左。本文考慮的是超聲速流動。入口邊界條件由線性穩(wěn)定性理論給出,其擾動形式為出口邊界條件由單邊差分給出,同時,為了減弱出流反射波的強度,在出口處添加一段嵌邊區(qū)[14],其中預(yù)先指定的嵌邊區(qū)出流處的擾動值設(shè)定為0,即Uc=0,上下邊界采用完全無反射邊界條件[15]。

    3 數(shù)值算例

    3.1 對流馬赫數(shù)Mc=1.2(≥1)的情況

    基本流通過求解定常可壓縮二維邊界層方程[16]得到,入口取無量綱長度δω0=1,對應(yīng)的流向長度x0≈337.50,圖1給出了入口處的基本流剖面。

    圖1 入口處基本流剖面Fig.1 Base flow profile at the inflow

    由線性穩(wěn)定性理論得到,對于Mc≥1.0的超聲速混合層流動,同一擾動頻率ω0對應(yīng)兩個增長波-αi,根據(jù)擾動相速度相對于自由流速度Mr,j= |ω/Re{α}-Ui|/ai(i,j=1,2)是否大于1,可以把擾動形式分為快慢兩種模態(tài)[17],當擾動相速度相對于快速流為亞聲速,相對于慢速流為超聲速,即Mr,1<1.0,Mr,2>1.0,稱其為快模態(tài);當擾動相速度相對于快速流為超聲速,相對于慢速流為亞聲速,即Mr,1>1.0,Mr,2<1.0,稱其為慢模態(tài)。

    下面計算在入口處分別加入不同擾動模態(tài)下的擾動演化過程,選取的擾動頻率為慢模態(tài)情況下最不穩(wěn)定波頻率,其具體的計算參數(shù)如表1。

    表1 入口加入不穩(wěn)定波的具體參數(shù)Table 1 Parameters of inflow instability wave

    圖2、圖3分別是工況A和工況B條件下擾動特征函數(shù)圖像,分別對應(yīng)于擾動快慢模態(tài)兩種情況,同一擾動頻率ω0=0.2,分別對應(yīng)著兩種不穩(wěn)定擾動波,它們的增長率-αi分別為0.0114和0.0131,圖2給出的是超聲速擾動波快模態(tài)情形,擾動相速度相對于上層快速流為亞聲速Mr,1=0.904(<1.0),而相對于下層慢速流為超聲速Mr,2=1.496(>1.0),擾動的特征函數(shù)沿上層快速流方向呈指數(shù)形式衰減,而擾動沿下層慢速流方向呈正弦波的形式衰減;與圖2的情況相反,圖3給出的是超聲速擾動波慢模態(tài)情形,擾動相速度相對于上層快速流為超聲速Mr,1=1.52 (>1.0),而相對于下層慢速流為亞聲速Mr,2=0.88 (<1.0),擾動的特征函數(shù)沿上層快速流方向呈正弦波的形式衰減,而擾動沿下層慢速流方向呈指數(shù)形式衰減。

    為了驗證程序,當前用擾動方程計算了慢模態(tài)下最不穩(wěn)定頻率ω=0.2的擾動波的演化,初始幅值A(chǔ)0=0.0001,圖4給出了擾動方程與LST的計算結(jié)果,二者吻合的很好,表明擾動方程程序的正確性。

    圖2 快模態(tài)的特征函數(shù)圖像Fig.2 Eigenfunction for fast mode

    圖3 慢模態(tài)的特征函數(shù)圖像Fig.3 Eigenfunction for slow mode

    圖4擾動波幅值演化,ω=0.2,A0=0.0001Fig.4 Amp lification of the disturbance,ω=0.2,A0=0.0001

    圖5、圖6分別給出在入口處加入A、B工況下的不穩(wěn)定波,初始幅值都給A0=0.005,其擾動向下游演化過程中的壓力擾動等值線p1。

    比較A、B兩種工況的壓力等值線圖,可以看出對于擾動相速度相對于自由流速度為超聲速(Mr,j>1.0)的部分,壓力擾動以一定的角度向遠場輻射,形成馬赫波輻射,通過馬赫角的計算公式μ= arcsin(1/Mr,j)可以求出快模態(tài)擾動的理論馬赫角為μ≈41.9°,由圖5的壓力等值線求得輻角μ≈40.7°,與理論值相符。同理,由圖6的慢模態(tài)等值線可以求得馬赫角為μ≈40°,也與理論值μ≈41.0°相符。

    圖5 加入快模態(tài)擾動壓力等值線Fig.5 Pressure contours for fast mode

    圖6 加入慢模態(tài)擾動壓力等值線Fig.6 Pressure contours for slow mode

    圖7 加入快慢模態(tài)的擾動壓力等值線Fig.7 Disturbance pressure contours for fast and slow mode

    在入口處同時加入快慢模態(tài)的擾動波,得到圖7所示的近場展向渦量和遠場擾動壓力等值線,它兼有二者的特性,擾動壓力向自由流兩側(cè)方向以馬赫波的形式輻射,從近場展向渦量ωz等值線看,沒有出現(xiàn)大尺度渦結(jié)構(gòu),圖8對應(yīng)著基本波能量幅值E10與其倍頻波的能量幅值E20。

    圖8 基本波與其倍頻波的能量幅值曲線Fig.8 Energy for the fundamental and first harmonic mode

    因此,對于Mc≥1的超聲速混合層流動,其發(fā)聲機制主要是不穩(wěn)定波以馬赫波的形式向遠場輻射。

    3.2 對流馬赫數(shù)Mc=0.5(<1.0)的情況

    利用線性穩(wěn)定性理論,可以求出入口處所加的擾動的特征值和特征函數(shù),計算求出的具體參數(shù)見表2。

    表2 入口加入不穩(wěn)定波的具體參數(shù)Table 2 Parameters of inflow instability wave

    圖9給出了最不穩(wěn)定波ω0=0.578對應(yīng)的擾動特征函數(shù)幅值曲線,可以看出擾動特征函數(shù)幅值以指數(shù)形式向自由流兩側(cè)快速衰減。

    圖10、圖11分別是在入口處加入工況D和E情況下的擾動波,使用非線性擾動方程(NLDE)計算得到擾動向下游演化過程中的遠場壓力等值線p和近場展向渦量云圖ωz,圖12是其兩組工況下所對應(yīng)的能量幅值曲線E。

    圖9 特征函數(shù)幅值曲線Fig.9 Amp litudes curve of the eigenfunction

    從圖10的近場渦量等值線可以看到明顯的大尺度渦結(jié)構(gòu),伴隨著渦結(jié)構(gòu)向下游演化,其擾動壓力直接輻射向自由流遠場,渦結(jié)構(gòu)形成位置在x≈100處,從圖12的能量幅值曲線上看,在此處,能量幅值達到最大,且正好對應(yīng)于遠場擾動壓力開始輻射的位置。

    圖10 工況D遠場壓力和近場渦量云圖Fig.10 Contours of the far-field pressure and nearfield spanw ise vorticity for case D

    圖11是在入口處加入基本波與亞諧波的計算結(jié)果,與在入口處只加入基本波的情況相似,其近場區(qū)出現(xiàn)大尺度渦結(jié)構(gòu),渦結(jié)構(gòu)形成的位置與遠場壓力輻射位置相對應(yīng),不同的是,渦結(jié)構(gòu)在更下游的位置處發(fā)生了渦與渦之間的合并,其渦合并發(fā)生在x≈190處,渦合并現(xiàn)象是基本波與亞諧波的非線性相互作用導(dǎo)致,圖12中的能量幅值曲線顯示在渦合并處,基本波與亞諧波的能量都處于峰值,此時,它們之間的非線性作用是最強的,從遠場壓力等值線圖可以看出它們的壓力強度也同時達到最大。

    圖11 工況E遠場擾動壓力和近場渦量云圖Fig.11 Contours of the far-field pressure and near-field spanw ise vorticity for case E

    圖12 基本波與其倍頻波的能量幅值曲線Fig.12 Am p litudes of energy for the fundamental and subharmonic mode

    因此,對于對流馬赫數(shù)Mc<1.0的超聲速混合層,其發(fā)聲機制主要是由大尺度渦結(jié)構(gòu)及其合并引起,且基本波與亞諧波之間的非線性相互作用能增強聲輻射強度。

    4 結(jié)論

    本文采用非線性擾動方程(NLDE)研究了超聲速情況下擾動的動力學特性和遠場發(fā)聲機制間的關(guān)聯(lián),分析了對流馬赫數(shù)Mc=0.5(<1)和Mc=1.2 (≥1)兩種情況下不穩(wěn)定波的發(fā)聲機制。得到了以下幾點結(jié)論:

    1)對流馬赫數(shù)Mc≥1.0的情況,同一頻率ω0的擾動對應(yīng)著兩個不同的不穩(wěn)定波-αi,分別對應(yīng)著快慢兩種模態(tài),其擾動特征函數(shù)呈現(xiàn)兩種完全不同的形式。

    2)分別對快慢兩種擾動模態(tài)及其相互作用進行數(shù)值模擬研究,計算得到流場中不穩(wěn)定波的發(fā)聲機制:壓力擾動以馬赫波形式向遠場輻射,即馬赫波輻射。

    3)對于Mc<1.0的情況,非線性擾動方程(NLDE)能夠準確捕捉到大尺度渦結(jié)構(gòu),并在下游處捕捉到渦對合并現(xiàn)象。

    4)對于Mc<1.0的情況,不穩(wěn)波產(chǎn)生遠聲場的機理主要是大尺度渦運動及渦對的合并,而基本波與亞諧波之間的非線性相互作用能明顯增強聲輻射的強度。

    [1]Hu Guoqing.Numerical simulation of aeroacoustic field generated by compressible flow[D].Beijing:Institute of Mechanics,Chinese Academy of Sciences,2001.(in Chinese)胡國慶.可壓縮流氣動聲場的數(shù)值模擬[D].北京:中國科學院力學研究所,2001.

    [2]Ju Hongbin,Shen Mengyu.Problem about computational aeroacoustics,methods and progress[J].Mechanics in Engineering,1995,17(5):1-10.(in Chinese)居鴻賓,沈孟育.計算氣動聲學的問題,方法與進展[J].力學與實踐,1995,17(5):1-10.

    [3]Colonius T,Lele S K,Moin P.Sound generation in a mixing layer[J].Journal of Fluid Mechanics,1997,330:375-409.

    [4]Bogey C,Bailly C,Juvé D.Numerical simulation of sound generated by vortex pairing in a mixing layer[J].AIAA Journal,2000,38(12):2210-2218.

    [5]Mitchell B E,Lele S K,Moin P.Direct computation of the sound generated by vortex pairing in an axisymmetric jet[J].Journal of Fluid Mechanics,1999,383:113-142.

    [6]Li Xiaodong,Gao Junhui.Numerical simulation of sound generation and radiation from a 2D parallel shear layer[J].Journal of Aerospace,2003,24(1):6-9.(in Chinese)李曉東,高軍輝.二維平行剪切層聲波產(chǎn)生和輻射的數(shù)值模擬[J].航空學報,2003,24(1):6-9.

    [7]Xing Chao,Luo Jisheng.Numercial simulation of mixing layer effect on sound waves[J].Journal of Aerospace Power,2005,20(1):8-12.(in Chinese)邢超,羅紀生.聲波與剪切流作用的數(shù)值模擬[J].航空動力學報,2005,20(1):8-12.

    [8]Tam C K W,Morris P J.The radiation of sound by the instability waves of a compressible plane turbulent shear layer[J].Journal of Fluid Mechanics,1980,98(02):349-381.

    [9]Tam C K W,Burton D E.Sound generated by instability waves of supersonic flows.Part 1:two-dimensional mixing layers[J].Journal of Fluid Mechanics,1984,138:249-271.

    [10]Wu Xuesong.Generation of sound and instability waves due to unsteady suction and injection[J].Journal of Fluid Mechanics,2002,453:289-313.

    [11]Wu Xuesong.Mach wave radiation of nonlinearly evolving supersonic instability modes in shear layers[J].Journal of Fluid Mechanics,2005,523:121-159.

    [12]Bertolotti F,Colonius T.On the noise generated by convected structures in a Mach 0.9 hot,turbulent jet[J].AIAA Paper,2003,1062:2003.

    [13]Cheung L C,Lele S K.Linear and nonlinear processes in twodimensional mixing layer dynamics and sound radiation[J].Journal of Fluid Mechanics,2009,625:321-351.

    [14]Qian Lin.The numerical study of the enhancement effect in 3D supersonic mixing layer with swirl[D].Tianjin University,2007.(in Chinese).錢琳.三維超聲速旋轉(zhuǎn)混合層增混作用的數(shù)值研究[D].天津大學,2007.

    [15]Ju Hongbin.Non-reflective boundary conditions and sound source modeling in computational aeroacoustics[J].Shanghai Jiaotong University,1998,32(7):36-39.(in Chinese)居鴻賓.氣動聲學數(shù)值模擬中無反射邊界處理及聲源模型的建立[J].上海交通大學學報,1998,32(7):36-39.

    [16]Day M J,Mansour N N,Reynolds W C.Nonlinear stability and structure of compressible reacting mixing layers[J].Journal of Fluid Mechanics,2001,446:375-408.

    [17]Fu Dexun,Ma Yanwen,Li Xinliang,Direct numerical simulation of compressible turbulent flow[M].Beijing:Science Press,2010.(in Chinese)傅德薰,馬延文,李新亮.可壓縮湍流直接數(shù)值模擬[M].北京:科學出版社,2010.

    Studing sound generation mechanism s of m ixing layer using nonlinear disturbance equations

    Zhang Chengyuan,F(xiàn)ang Yihong*
    (Department of Mechanics,School of Mechanical Engineering,Tianjin University,Tianjin 300072,China)

    The noise problem caused by pulse exists in the aviation and natural industry,studying sound generation mechanism of the unstable wave is of important significance for understanding,predicting and controlling the aerodynamic noise.In this paper,Nonlinear Disturbance Equation(NLDE)is used to study the relationship between near field dynamic characteristic and the far field sound radiation for supersonic mixing layer.For different types of convective Mach number,disturbance pressure are selected and calculated respectively.It is gained two different sound generation mechanisms.When the convective Mach number of supersonic mixing layer is less than 1,the sound is generated by the largescale vortex structure in near-field.Otherwise,the generation of sound is depended on whether the phase velocity of disturbance is supersonic relative to the upper flow or not.It turned out that the supersonic part of the disturbance pressure radiate to far-field at a fix angle,that is so called Mach wave radiation.The calculated angle of Mach wave is consistent with the theoretical one.

    nonlinear disturbance equations;supersonic mixing layer;convective Mach number; noise generation by vortex pairing;fast and slow mode;Mach wave radiation

    V212.1

    A

    10.7638/kqdlxxb-2014.0117

    0258-1825(2016)05-0659-07

    2015-01-15;

    2015-04-15

    國家自然科學基金(11332007,10972156)

    張程遠(1990-),男,福建,碩士研究生,主要從事邊界層轉(zhuǎn)捩預(yù)測的研究.E-mail:zhangcy9798@163.com

    方一紅*,副教授,主要從事流動穩(wěn)動性和空氣動力學.E-mail:yhfang@tju.edu.cn

    張程遠,方一紅.非線性擾動方程研究平面混合層發(fā)聲機制[J].空氣動力學學報,2016,34(5):659-665.

    10.7638/kqdlxxb-2014.0117 Zhang C Y,F(xiàn)ang Y H.Studing sound generation mechanisms of mixing layer using nonlinear disturbance equations[J].Acta Aerodynamica Sinica,2016,34(5):659-665.

    猜你喜歡
    特征函數(shù)等值線馬赫數(shù)
    一維非等熵可壓縮微極流體的低馬赫數(shù)極限
    基于規(guī)則預(yù)計格網(wǎng)的開采沉陷等值線生成算法*
    礦山測量(2020年6期)2021-01-07 04:52:06
    載荷分布對可控擴散葉型性能的影響
    亞純函數(shù)的Borel方向與Tsuji特征函數(shù)
    隨機變量的特征函數(shù)在概率論中的應(yīng)用
    黑龍江科學(2020年5期)2020-04-13 09:14:04
    特征函數(shù)的性質(zhì)在實變函數(shù)中的應(yīng)用
    特征函數(shù)在伽瑪分布中一個恒等式的證明及推廣
    等值線“慣性”變化規(guī)律的提出及應(yīng)用
    地理教學(2016年19期)2016-11-21 05:01:49
    基于Kriging插值的等值線生成算法研究
    等值線分析系統(tǒng)實際應(yīng)用之等值線填充
    97人妻精品一区二区三区麻豆| 99久久精品国产亚洲精品| 别揉我奶头 嗯啊视频| 国产激情偷乱视频一区二区| 国内精品美女久久久久久| 好看av亚洲va欧美ⅴa在| 波野结衣二区三区在线| 国产午夜精品久久久久久一区二区三区 | 日韩欧美三级三区| 99国产精品一区二区蜜桃av| 男女之事视频高清在线观看| 老司机午夜十八禁免费视频| 91在线观看av| 一进一出抽搐gif免费好疼| 色吧在线观看| 亚洲成av人片免费观看| 九色国产91popny在线| 精品福利观看| 久久99热这里只有精品18| 岛国在线免费视频观看| 成人亚洲精品av一区二区| 国产伦一二天堂av在线观看| 深夜a级毛片| av在线天堂中文字幕| 一进一出抽搐动态| 人人妻人人看人人澡| 12—13女人毛片做爰片一| 日韩欧美三级三区| 精品欧美国产一区二区三| 老女人水多毛片| 黄色日韩在线| 国产精品,欧美在线| 嫩草影院入口| 国产aⅴ精品一区二区三区波| 国产白丝娇喘喷水9色精品| 精品人妻视频免费看| 极品教师在线免费播放| 国产精品伦人一区二区| 欧美不卡视频在线免费观看| 最新中文字幕久久久久| 日韩高清综合在线| 夜夜躁狠狠躁天天躁| 亚洲成人免费电影在线观看| 男女床上黄色一级片免费看| 日韩有码中文字幕| 午夜日韩欧美国产| 国产午夜精品久久久久久一区二区三区 | 特级一级黄色大片| 在线观看一区二区三区| 色精品久久人妻99蜜桃| 老女人水多毛片| 熟女电影av网| 精品久久久久久久久av| 老司机午夜十八禁免费视频| 欧美日韩国产亚洲二区| 婷婷色综合大香蕉| 色精品久久人妻99蜜桃| 九色国产91popny在线| 一区二区三区激情视频| 精品久久久久久久久av| 欧美黄色淫秽网站| h日本视频在线播放| 国产不卡一卡二| 搡老妇女老女人老熟妇| 在线播放国产精品三级| 精品国产亚洲在线| 精品久久久久久久久av| 女人被狂操c到高潮| 色5月婷婷丁香| 午夜福利在线在线| 最近中文字幕高清免费大全6 | 中文字幕人妻熟人妻熟丝袜美| 欧美最新免费一区二区三区 | 蜜桃亚洲精品一区二区三区| 亚洲精品一卡2卡三卡4卡5卡| 一区福利在线观看| 久久九九热精品免费| 又黄又爽又刺激的免费视频.| 999久久久精品免费观看国产| 日韩欧美三级三区| 好男人在线观看高清免费视频| 精品一区二区三区av网在线观看| 日韩免费av在线播放| 亚洲国产欧美人成| 午夜激情欧美在线| 国产精品不卡视频一区二区 | 日韩欧美国产在线观看| 免费在线观看影片大全网站| 美女xxoo啪啪120秒动态图 | 亚洲国产欧美人成| 99久国产av精品| 激情在线观看视频在线高清| 黄片小视频在线播放| 欧美成狂野欧美在线观看| 免费人成在线观看视频色| 成人av一区二区三区在线看| av在线天堂中文字幕| 老司机午夜福利在线观看视频| 757午夜福利合集在线观看| 久久精品国产自在天天线| 色播亚洲综合网| 久久精品国产亚洲av天美| 三级男女做爰猛烈吃奶摸视频| 国产欧美日韩一区二区三| 午夜老司机福利剧场| 精品久久久久久久久久免费视频| 婷婷精品国产亚洲av在线| 日韩欧美免费精品| 亚洲久久久久久中文字幕| 久久香蕉精品热| 哪里可以看免费的av片| 国产精品伦人一区二区| 亚洲电影在线观看av| 天堂动漫精品| 久久精品国产99精品国产亚洲性色| 欧美激情国产日韩精品一区| 免费电影在线观看免费观看| 无人区码免费观看不卡| 欧美不卡视频在线免费观看| 99在线人妻在线中文字幕| 非洲黑人性xxxx精品又粗又长| 91在线观看av| 高清日韩中文字幕在线| 少妇人妻精品综合一区二区 | 九九热线精品视视频播放| 夜夜夜夜夜久久久久| 伊人久久精品亚洲午夜| 在线免费观看的www视频| 久99久视频精品免费| 欧美高清成人免费视频www| 午夜福利18| 欧美zozozo另类| av黄色大香蕉| 人人妻,人人澡人人爽秒播| 少妇裸体淫交视频免费看高清| 国产精品久久久久久精品电影| 欧美激情在线99| 日韩欧美国产在线观看| 2021天堂中文幕一二区在线观| 岛国在线免费视频观看| 国产av麻豆久久久久久久| 色在线成人网| 久久6这里有精品| 亚洲第一电影网av| 免费在线观看影片大全网站| 丝袜美腿在线中文| 日日摸夜夜添夜夜添小说| 午夜影院日韩av| 男女视频在线观看网站免费| 小蜜桃在线观看免费完整版高清| www.色视频.com| 成人av在线播放网站| 好看av亚洲va欧美ⅴa在| 精品国产三级普通话版| 亚洲人与动物交配视频| 欧美最新免费一区二区三区 | 日韩欧美国产在线观看| 久久人人精品亚洲av| 免费人成视频x8x8入口观看| 国产伦精品一区二区三区四那| 深夜a级毛片| 国模一区二区三区四区视频| 色综合婷婷激情| av国产免费在线观看| 日韩高清综合在线| 国产探花极品一区二区| 亚洲av免费在线观看| 国产精品美女特级片免费视频播放器| 国产精品1区2区在线观看.| 国产精品一区二区三区四区久久| 国产亚洲欧美98| av福利片在线观看| av在线天堂中文字幕| 日本在线视频免费播放| 亚洲真实伦在线观看| 毛片一级片免费看久久久久 | 欧美绝顶高潮抽搐喷水| 波多野结衣高清作品| 99热这里只有是精品在线观看 | 99在线视频只有这里精品首页| 亚洲成av人片在线播放无| 欧美在线黄色| 在线免费观看不下载黄p国产 | 国产主播在线观看一区二区| 国产视频内射| 婷婷亚洲欧美| 成人三级黄色视频| 长腿黑丝高跟| 欧美日韩乱码在线| 国产亚洲av嫩草精品影院| 亚洲av电影不卡..在线观看| 亚洲午夜理论影院| 久久精品人妻少妇| 成人特级av手机在线观看| 青草久久国产| 精品99又大又爽又粗少妇毛片 | 欧美成人免费av一区二区三区| 窝窝影院91人妻| 午夜福利高清视频| 18禁裸乳无遮挡免费网站照片| av专区在线播放| 亚洲精品成人久久久久久| 嫁个100分男人电影在线观看| 亚洲成av人片在线播放无| www.色视频.com| 午夜激情欧美在线| 亚洲精品在线美女| 免费人成视频x8x8入口观看| 在线观看舔阴道视频| avwww免费| 可以在线观看的亚洲视频| 亚洲va日本ⅴa欧美va伊人久久| 亚洲av免费在线观看| 首页视频小说图片口味搜索| 日本黄大片高清| 超碰av人人做人人爽久久| av黄色大香蕉| 五月伊人婷婷丁香| 欧美xxxx性猛交bbbb| 亚洲av五月六月丁香网| av在线天堂中文字幕| 国产又黄又爽又无遮挡在线| 色在线成人网| 国产成人av教育| 久久久久国产精品人妻aⅴ院| 五月玫瑰六月丁香| 精品99又大又爽又粗少妇毛片 | 国产精品爽爽va在线观看网站| 色尼玛亚洲综合影院| 国产三级中文精品| 老鸭窝网址在线观看| 日韩欧美在线乱码| 国产欧美日韩一区二区三| 久久性视频一级片| 久久精品国产亚洲av天美| 99久久九九国产精品国产免费| 国产乱人视频| www.999成人在线观看| 免费观看的影片在线观看| 亚洲第一电影网av| 在线天堂最新版资源| eeuss影院久久| 亚洲性夜色夜夜综合| 欧美xxxx黑人xx丫x性爽| 色av中文字幕| 日韩欧美在线乱码| 国产午夜精品久久久久久一区二区三区 | 不卡一级毛片| 免费人成视频x8x8入口观看| 国产熟女xx| 欧美黄色片欧美黄色片| www.熟女人妻精品国产| 97人妻精品一区二区三区麻豆| 简卡轻食公司| 亚洲中文字幕一区二区三区有码在线看| 51午夜福利影视在线观看| 欧美最黄视频在线播放免费| 亚洲激情在线av| 美女高潮的动态| 每晚都被弄得嗷嗷叫到高潮| 男插女下体视频免费在线播放| 欧美激情国产日韩精品一区| 国产精品自产拍在线观看55亚洲| av专区在线播放| 久久久久久久久大av| 欧美中文日本在线观看视频| 亚洲精品亚洲一区二区| 能在线免费观看的黄片| 麻豆久久精品国产亚洲av| 免费搜索国产男女视频| 内射极品少妇av片p| 亚洲国产欧美人成| 国产综合懂色| 欧美午夜高清在线| 直男gayav资源| 村上凉子中文字幕在线| 日本撒尿小便嘘嘘汇集6| 精品国产亚洲在线| 高清毛片免费观看视频网站| 国产白丝娇喘喷水9色精品| 中文字幕熟女人妻在线| 日韩欧美 国产精品| 精品欧美国产一区二区三| 中文在线观看免费www的网站| 精品国产亚洲在线| 久久久久久久久久黄片| 在现免费观看毛片| 亚洲一区二区三区不卡视频| av欧美777| 亚洲欧美日韩卡通动漫| 天堂网av新在线| 三级国产精品欧美在线观看| 嫩草影院新地址| 淫妇啪啪啪对白视频| 日韩免费av在线播放| 国产高清有码在线观看视频| 日韩欧美国产在线观看| 久久精品人妻少妇| 成年女人永久免费观看视频| 欧美中文日本在线观看视频| 婷婷精品国产亚洲av| 国产一区二区在线av高清观看| 成熟少妇高潮喷水视频| av在线蜜桃| or卡值多少钱| 观看免费一级毛片| 日本五十路高清| 日韩欧美精品免费久久 | 男女下面进入的视频免费午夜| 一本久久中文字幕| 搞女人的毛片| 亚洲激情在线av| 可以在线观看的亚洲视频| 波野结衣二区三区在线| 日韩欧美一区二区三区在线观看| 亚洲国产色片| av天堂中文字幕网| www.999成人在线观看| 在线看三级毛片| 悠悠久久av| 青草久久国产| 身体一侧抽搐| 午夜精品在线福利| 永久网站在线| 1024手机看黄色片| 不卡一级毛片| 亚洲国产色片| 精品人妻偷拍中文字幕| 精品人妻熟女av久视频| 色哟哟哟哟哟哟| 亚洲精品久久国产高清桃花| 麻豆成人av在线观看| 男女视频在线观看网站免费| 波野结衣二区三区在线| 亚洲精品久久国产高清桃花| 日韩中文字幕欧美一区二区| 亚洲成av人片免费观看| 在线天堂最新版资源| 欧美成人一区二区免费高清观看| 成人高潮视频无遮挡免费网站| 18禁在线播放成人免费| 亚洲一区二区三区不卡视频| 舔av片在线| 真人做人爱边吃奶动态| 国产一区二区亚洲精品在线观看| 日日摸夜夜添夜夜添小说| 天堂√8在线中文| 真实男女啪啪啪动态图| 宅男免费午夜| 亚洲在线观看片| 国产私拍福利视频在线观看| 亚洲人成网站高清观看| 99在线人妻在线中文字幕| 成人美女网站在线观看视频| 99在线人妻在线中文字幕| 人妻丰满熟妇av一区二区三区| 亚洲av美国av| 免费观看精品视频网站| 国产黄a三级三级三级人| 久久久精品大字幕| 亚洲av成人av| 日日夜夜操网爽| 变态另类丝袜制服| 亚洲av中文字字幕乱码综合| 中文字幕熟女人妻在线| 国产精品亚洲美女久久久| 听说在线观看完整版免费高清| 真人做人爱边吃奶动态| 琪琪午夜伦伦电影理论片6080| 成人欧美大片| 综合色av麻豆| 中文资源天堂在线| 少妇裸体淫交视频免费看高清| 国产成+人综合+亚洲专区| 美女xxoo啪啪120秒动态图 | 久久性视频一级片| 国产高潮美女av| 亚洲国产色片| 国产精品久久久久久亚洲av鲁大| 非洲黑人性xxxx精品又粗又长| 亚洲精品亚洲一区二区| 美女高潮的动态| 少妇人妻一区二区三区视频| 人妻久久中文字幕网| 18禁黄网站禁片免费观看直播| 亚洲一区高清亚洲精品| 3wmmmm亚洲av在线观看| 丰满人妻熟妇乱又伦精品不卡| 亚洲综合色惰| 夜夜躁狠狠躁天天躁| 一二三四社区在线视频社区8| 免费在线观看日本一区| 99在线人妻在线中文字幕| 国产av麻豆久久久久久久| 欧美一区二区国产精品久久精品| 91狼人影院| 精品免费久久久久久久清纯| 免费无遮挡裸体视频| 嫩草影院入口| 国产一区二区三区视频了| 男人和女人高潮做爰伦理| 每晚都被弄得嗷嗷叫到高潮| 能在线免费观看的黄片| 国产三级在线视频| 欧美成人免费av一区二区三区| 国内揄拍国产精品人妻在线| 最近最新中文字幕大全电影3| 天堂动漫精品| 18+在线观看网站| 中文字幕熟女人妻在线| 欧美精品啪啪一区二区三区| 中亚洲国语对白在线视频| 三级男女做爰猛烈吃奶摸视频| 亚洲欧美清纯卡通| 国内少妇人妻偷人精品xxx网站| 亚洲人与动物交配视频| 亚洲国产精品sss在线观看| 一卡2卡三卡四卡精品乱码亚洲| 九九久久精品国产亚洲av麻豆| 在线播放国产精品三级| 天堂√8在线中文| 亚洲中文日韩欧美视频| а√天堂www在线а√下载| 国产乱人伦免费视频| 俺也久久电影网| 两个人的视频大全免费| 老司机午夜十八禁免费视频| 国产一区二区在线观看日韩| 又爽又黄a免费视频| 精品午夜福利在线看| 国产激情偷乱视频一区二区| 一个人免费在线观看电影| 国内精品一区二区在线观看| 俺也久久电影网| 一区二区三区高清视频在线| 亚洲欧美激情综合另类| avwww免费| 九色国产91popny在线| 国产成人福利小说| 中文资源天堂在线| 欧美黑人巨大hd| 国产精品久久久久久久电影| 国产真实乱freesex| 又粗又爽又猛毛片免费看| 欧美色欧美亚洲另类二区| 中文字幕人成人乱码亚洲影| 午夜久久久久精精品| 亚洲熟妇中文字幕五十中出| 国产精品人妻久久久久久| 免费观看人在逋| 婷婷精品国产亚洲av| 欧美乱妇无乱码| 麻豆久久精品国产亚洲av| 给我免费播放毛片高清在线观看| 亚洲avbb在线观看| 国产精品久久久久久久久免 | 男人和女人高潮做爰伦理| 亚洲人成伊人成综合网2020| 变态另类成人亚洲欧美熟女| 午夜福利在线观看吧| 国产乱人伦免费视频| 日韩av在线大香蕉| 国产白丝娇喘喷水9色精品| 久久亚洲真实| 亚洲va日本ⅴa欧美va伊人久久| 最近中文字幕高清免费大全6 | 日本三级黄在线观看| 一a级毛片在线观看| 久久久国产成人精品二区| 波多野结衣高清无吗| 日韩有码中文字幕| 国产高潮美女av| 真人一进一出gif抽搐免费| 动漫黄色视频在线观看| 一区二区三区免费毛片| 欧美激情久久久久久爽电影| 日韩欧美免费精品| 久久久久免费精品人妻一区二区| 露出奶头的视频| 中文字幕高清在线视频| 欧美一区二区国产精品久久精品| 成人欧美大片| 亚洲精品粉嫩美女一区| 欧美最黄视频在线播放免费| АⅤ资源中文在线天堂| 怎么达到女性高潮| 亚洲熟妇中文字幕五十中出| 夜夜躁狠狠躁天天躁| 日本精品一区二区三区蜜桃| 日韩欧美国产在线观看| 一进一出好大好爽视频| 又爽又黄a免费视频| 国产欧美日韩一区二区三| 久久久久精品国产欧美久久久| 亚洲综合色惰| 国产又黄又爽又无遮挡在线| 88av欧美| 国产av麻豆久久久久久久| 国产一区二区在线av高清观看| 99热精品在线国产| 波多野结衣高清作品| 脱女人内裤的视频| 变态另类丝袜制服| 国产伦精品一区二区三区四那| 国产一区二区激情短视频| 午夜视频国产福利| 久久久久久国产a免费观看| 麻豆一二三区av精品| 国产精品三级大全| 可以在线观看毛片的网站| 久久久久久国产a免费观看| 小蜜桃在线观看免费完整版高清| 九九热线精品视视频播放| 国产亚洲精品综合一区在线观看| 麻豆国产97在线/欧美| 看片在线看免费视频| 精品久久国产蜜桃| 看黄色毛片网站| 欧美性猛交╳xxx乱大交人| 国产亚洲精品av在线| 日韩欧美免费精品| 在线十欧美十亚洲十日本专区| 亚洲aⅴ乱码一区二区在线播放| 激情在线观看视频在线高清| 亚洲五月天丁香| 白带黄色成豆腐渣| 国模一区二区三区四区视频| 在线观看66精品国产| 精品久久久久久久末码| 全区人妻精品视频| 俺也久久电影网| 亚洲一区二区三区色噜噜| 麻豆成人午夜福利视频| 一区二区三区高清视频在线| 成人性生交大片免费视频hd| 国产爱豆传媒在线观看| 少妇人妻精品综合一区二区 | 欧美在线黄色| 欧美3d第一页| 亚洲国产精品999在线| 精品午夜福利视频在线观看一区| 在线观看一区二区三区| 国产久久久一区二区三区| 色综合亚洲欧美另类图片| 亚洲天堂国产精品一区在线| 我的女老师完整版在线观看| 级片在线观看| 高清日韩中文字幕在线| 国产色爽女视频免费观看| 精品午夜福利在线看| 精品久久久久久久人妻蜜臀av| 一个人观看的视频www高清免费观看| 日本免费一区二区三区高清不卡| 亚洲美女搞黄在线观看 | 国产av不卡久久| 亚洲欧美日韩东京热| 最近在线观看免费完整版| 18禁黄网站禁片免费观看直播| 啦啦啦韩国在线观看视频| 啪啪无遮挡十八禁网站| 好男人电影高清在线观看| 又爽又黄a免费视频| 欧美一区二区亚洲| 精品一区二区三区人妻视频| 性欧美人与动物交配| 久久精品国产亚洲av香蕉五月| 国产色爽女视频免费观看| 日本精品一区二区三区蜜桃| 高潮久久久久久久久久久不卡| 少妇的逼好多水| 美女被艹到高潮喷水动态| 黄色视频,在线免费观看| 在线十欧美十亚洲十日本专区| 色精品久久人妻99蜜桃| 精品一区二区三区视频在线观看免费| 在现免费观看毛片| 在线观看av片永久免费下载| 国产v大片淫在线免费观看| 床上黄色一级片| 免费av不卡在线播放| 夜夜夜夜夜久久久久| 一二三四社区在线视频社区8| 白带黄色成豆腐渣| 很黄的视频免费| 国产精品,欧美在线| 亚洲最大成人中文| 51国产日韩欧美| 亚洲av.av天堂| 成人永久免费在线观看视频| 亚洲国产色片| 国产精品久久久久久精品电影| 亚洲,欧美,日韩| 国产人妻一区二区三区在| www.色视频.com| 一个人看视频在线观看www免费| 成人性生交大片免费视频hd| 丰满的人妻完整版| 老司机深夜福利视频在线观看| 能在线免费观看的黄片| 亚洲av.av天堂| 在线十欧美十亚洲十日本专区| 狂野欧美白嫩少妇大欣赏| 宅男免费午夜| 又黄又爽又免费观看的视频| 亚洲人成伊人成综合网2020| 波多野结衣高清作品| 国产精品1区2区在线观看.| 欧美黄色淫秽网站| 美女黄网站色视频| 亚洲中文日韩欧美视频| 69人妻影院| 成年版毛片免费区| 真人一进一出gif抽搐免费| 欧美又色又爽又黄视频| 一级av片app|