亚洲免费av电影一区二区三区,日韩爱爱视频,51精品视频一区二区三区,91视频爱爱,日韩欧美在线播放视频,中文字幕少妇AV,亚洲电影中文字幕,久久久久亚洲av成人网址,久久综合视频网站,国产在线不卡免费播放

        ?

        城際列車氣動性能分析與評估

        2016-09-08 03:20:12李文化尚克明楊明智
        鐵道科學與工程學報 2016年7期
        關鍵詞:湍流氣動阻力

        李文化,尚克明,楊明智

        (1.中車青島四方機車車輛股份有限公司,山東青島 266111;2.中南大學 交通運輸工程學院,湖南 長沙 410075;3.軌道交通安全教育部重點實驗室,湖南 長沙 410075)

        ?

        城際列車氣動性能分析與評估

        李文化1,尚克明1,楊明智2,3

        (1.中車青島四方機車車輛股份有限公司,山東青島 266111;2.中南大學 交通運輸工程學院,湖南 長沙 410075;3.軌道交通安全教育部重點實驗室,湖南 長沙 410075)

        城際列車的運行速度一般處于100~200km/h之間,對于這類列車氣動性能的研究基本處于空白狀態(tài),而隨著城際列車的不斷開通,其安全問題日益加深,因此有必要對城際列車進行氣動性能分析。采用數(shù)值模擬計算的方法,對溫州S1線列車的氣動性能進行研究。主要分析列車明線運行時阻力和噪聲、列車過隧道車體表面壓力變化和車廂內部壓力變化幾個方面。通過對結果的分析對比,得出其變化規(guī)律:隨著速度的增大,阻力增大,車體表面壓力增大,車內壓力增大,噪聲增大;列車過隧道時,40m2圓形隧道要比33m2方形隧道氣動性能好;對于單車過隧道,動態(tài)密封指數(shù)要大于10s;列車以160km/h速度運行時,線邊25m遠、沿列車軸線方向的最大噪聲級約80dB,在環(huán)境噪聲允許范圍內。

        數(shù)值模擬;城際列車;氣動性能

        隨著我國高速列車的蓬勃發(fā)展,高速列車的研制日益加快,列車的運行速度不斷攀上高峰。國內外很多研究人員針對高速列車的氣動性能進行了研究。梁習鋒等[1-8]研究了明線運行時高速列車的氣動性能;SUZUKI等[9-11]研究了高速列車過隧道時的氣動性能;劉加利等[12-18]研究了高速列車的氣動噪聲特性[15-20]。然而,對于中速度級列車氣動性能的研究明顯不足。城際列車的運行速度一般處于100~200km/h之間,對于城際列車的研究基本處于空白狀態(tài),而隨著城市間城際列車的不斷開通,城際列車的安全運行問題越來越受重視。因此應工程實際的要求,需要評估城際列車的整體性能,其中列車氣動性能是評價列車整體性能的主要參數(shù)。

        對列車氣動性能的研究方法主要有實車試驗、列車動模型試驗和數(shù)值模擬計算三大類。本文采用數(shù)值模擬計算的方法,對溫州S1線列車的氣動性能進行了研究和評估。計算分析了列車空氣阻力、車體表面壓力、車內壓力和噪聲四個主要參數(shù)。通過對以120,140和160km/h3種不同速度運行的列車數(shù)值模擬分析,得出了列車氣動性能變化規(guī)律,分析了溫州S1線城際列車行的氣動特性。

        1 基本理論

        1.1湍流模型

        與列車相關的流動現(xiàn)象絕大多數(shù)是湍流流動。因此,我們在討論列車周圍空氣流場流動的數(shù)值模擬時,自然也離不開如何模擬湍流現(xiàn)象的問題[19-20]。

        湍流的數(shù)值模擬方法主要有直接數(shù)值模擬湍流(DirectNumericalSimulationTurbulence-DNST)、大渦模擬(LargeEddySimulation-LES)和湍流模型(包括渦粘性模型、雷諾應力模型(ReynoldsStressModel)。

        目前,國際上所作的直接模擬僅限于較低的雷諾數(shù)和有簡單外形的問題,而復雜的列車周圍空氣流場流動數(shù)值模擬還不適應;大渦模擬雖然已經開始應用于列車的湍流和氣動噪聲的計算,但由于其要求網(wǎng)格尺度非常小、網(wǎng)格規(guī)模非常大,對計算機的速度和內存仍然有較高的要求,還不能完全滿足工程設計和應用的需要;雷諾應力模型和渦黏性模型中的k-ε兩方程模型,特別是k-ε兩方程模型,在列車周圍空氣流場流動的湍流數(shù)值模擬中應用最為廣泛。本文湍流的模擬采用k-ε湍流模型。

        (1)

        式中:ui為渦黏性系數(shù);k為湍流動能;ε湍流耗散率;Cu為湍流常數(shù),一般情況下取Cu=0.09。

        1.2氣動聲學方程

        Lighthill根據(jù)N-S方程和連續(xù)性方程推導了經典的氣動聲學傳播方程:

        (2)

        式中:Tij為Lighthill張量;eij為黏性應力張量;δij為單位張量;ρ′為流體密度的擾動量;ρ′為流場中壓力的脈動量;c為聲速。

        2 氣動阻力數(shù)值模擬及結果分析

        2.1幾何模型

        計算模型為溫州S1線列車,為4車連掛,及頭車+兩節(jié)中間車+尾車,中間以風擋連接,見圖1。列車模型比較復雜,包括轉向架、風擋、車鉤、防撞裝置、車下設備等復雜結構。

        圖1 計算區(qū)域Fig.1 Calculation zone

        氣動阻力計算來流區(qū)域長度取100m,尾流區(qū)域長度取200m,左右寬為60m,高60m,計算區(qū)域見圖2(a)。

        列車過隧道時,流場計算由于列車與隧道之間存在相對運動,為此,計算區(qū)域采用分區(qū)對接網(wǎng)格技術。為真實模擬列車通過隧道全過程,選擇列車頭部距隧道口50m作為初始運動點,尾流區(qū)域長度取200m,左右寬為60m,高60m,計算區(qū)域見圖2(b)。隧道長度由列車最不利隧道長度計算公式確定,隧道長度見表1。根據(jù)實際情況要求,計算時采用兩種不同隧道截面,一種是33m2方形隧道,一種是40m2圓形隧道,隧道截面形狀見圖3。

        單位:m(a)阻力計算區(qū)域;(b)過隧道計算區(qū)域圖2 計算區(qū)域尺寸圖Fig.2 Calculation zone size

        速度/(km·h-1)120140160隧道長度/m274320441588

        (a)40 m2圓形隧道;(b)33 m2方形隧道圖3 隧道截面形狀Fig.3 Tunnel cross section shape

        2.2邊界條件

        列車明線運行時,在入口ABFE給定第一類邊界條件:入口來流取理論上無窮遠處的均勻來流,速度為列車運行速度,方向與列車運行方向相反;在出口DCGH給定第二類邊界條件:壓力固定為參考壓力,靜壓為零。外邊界AEHD、EFGH、BFGC按光滑壁面處理,給定對稱面;列車表面為無滑移的壁面邊界條件;地面ABCD按光滑壁面處理,為消除地面效應,給出的是滑移邊界條件,大小和方向與入口來流一致[5-6]。

        列車過隧道時,對于列車表面,給定運動邊界條件,即X方向速度分量等于給定的列車運行速度V,Y、Z向速度分量等于0;流域底面和隧道壁面給定無滑移邊界條件;流域兩側面、頂面給定對稱面條件;出口給定壓力出口邊界條件,相對壓強為0;進口給定壓力進口邊界條件,相對壓強為0。

        2.3計算網(wǎng)格

        模型的復雜性給網(wǎng)格劃分帶來了極大的挑戰(zhàn)。網(wǎng)格劃分時,在列車頭部以及周圍流動量變化都較大的位置均采用了較細的網(wǎng)格,在列車轉向架位置和受電弓位置也采用了比較細的網(wǎng)格,以免因為網(wǎng)格過大原因模糊掉了復雜位置幾何結構。為了捕捉列車頭部來流和尾部尾流流場特性,還在來流區(qū)和尾流區(qū)進行了加密。圖4為計算網(wǎng)格。

        列車明線運行時,采用OpenFoam軟件對模型進行離散,物面為四邊形網(wǎng)格,空間為六面體網(wǎng)格??紤]邊界層的影響,對車體附近區(qū)域進行了局部網(wǎng)格加密,空間體網(wǎng)格約3 000萬。第一層網(wǎng)格的厚度為0.001 25m。圖4(a)為列車頭部網(wǎng)格示意圖。

        列車過隧道時,隧道和地面用六面體結構化網(wǎng)格離散;列車由于其外形為三維自由曲面,用四面體非結構化網(wǎng)格離散??臻g體網(wǎng)格約為1 000萬,最小網(wǎng)格尺寸為0.007m。圖4(b)為列車過隧道網(wǎng)格示意圖。

        2.4計算結果分析

        2.4.1氣動阻力數(shù)值模擬計算結果分析

        利用商用流體軟件Fluent,對溫州S1線列車進行了數(shù)值模擬計算。通過計算數(shù)值求解不可壓縮連續(xù)性方程、動量方程及湍流模型方程得到不同速度列車氣動阻力特性。表2為氣動阻力系數(shù)計算結果,圖5為列車壓力分布云圖,圖6為列車阻力分布圖,圖8為列車整車氣動阻力隨速度變化曲線。

        (a)列車頭部車體表面網(wǎng)格示意圖(b)列車過隧道表面網(wǎng)格示意圖圖4 計算網(wǎng)格Fig.4 Calculation grid

        圖5 列車壓力分布云圖Fig.5 Train pressure distribution diagram

        從圖5中可以看出,列車的頭車鼻錐端為駐點所在位置,此處流速為零,正壓壓力最大,經過該點后,氣流速度逐漸加快,壓力下降。其沿縱剖面的流動情況,到達頭部與車頂過渡處負壓降至最大值,而尾車開始時為較小的負壓,氣流到達尾部與車頂過渡處負壓降至最大值,這均是由于過渡弧面變化很大,空氣繞流速度加快,從而使這一區(qū)域的壓力急劇降低所致,至車體頂面壓力又再次回升成為平穩(wěn)的較小的負壓。

        表2不同速度列車氣動阻力系數(shù)

        Table2Aerodynamicdragcoefficientofdifferentspeedtrain

        速度/(km·h-1)120140160頭車0.27690.27290.2737中間車10.14320.14280.1425中間車20.14110.14130.1415尾車0.35110.34650.3428總阻力0.91230.90350.9005

        圖6 列車氣動阻力分布圖Fig.6 Aerodynamic drag distribution of train

        圖7 列車整車氣動阻力隨速度變化曲線Fig.7 Curve of aerodynamic drag along with the speed

        從上述圖表中可以看出,不同車速各節(jié)車氣動阻力系數(shù)基本相同;不同車速列車的氣動阻力變化規(guī)律一致,即尾車阻力最大,約占總阻力的30%,兩節(jié)中間車阻力相近最小,約占總阻力的15.5%,頭車阻力小于尾車大于中間車,約占總阻力的39%;所有車速的阻力均為正值;列車氣動阻力與車速的關系近似成平方,符合規(guī)律。

        以120km/h速度為例進一步分析。120km/h速度的列車頭車和尾車阻力比重最大的部位在流線型車頭位置,各占總阻力的60.8%和73.4%,這主要是因為溫州S1線列車為鈍頭車,使的列車阻力過大。而頂部空調位置加大了迎風面積,亦使阻力變大。為了減小列車阻力,需要對列車外形進行優(yōu)化,使其頭部流線位置變長,減小空調迎風面積,提高列車氣動性能。

        2.4.2列車過隧道數(shù)值模擬計算結果分析

        由于隧道壁面的限制,列車通過隧道的過程必須考慮空氣的壓縮性,同時,列車通過隧道的過程,壓力和速度等流場參數(shù)均隨時間在變化,因此,需要求解可壓縮和非定常的連續(xù)性方程、動量方程及湍流模型方程,以得到列車通過隧道的氣動特性。

        由于列車運行速度的增大,隧道通過性能隨之降低,列車通過單線隧道時,車體表面壓力幅值隨著速度的增大而增大;列車前部車鉤中心位置壓力最大,壓力變化幅值也最大,具體數(shù)值見表3。圖8為列車以120km/h速度運行時,該位置壓力變化曲線。因為隧道橫截面不同,圓形隧道的阻塞比要比方形隧道的小,所以列車通過圓形隧道時的氣動性能比通過方形隧道時好。

        表3 列車壓力最大值及幅值/Pa

        (b)40m2圓形隧道

        圖8 列車最大壓力變化曲線Fig.8 Maximum pressure curves of the train

        2.4.3列車過隧道車內壓力計算結果分析

        列車在隧道內運行時,車外的壓力變化會引起車廂內壓力的變化,從而會帶來乘客舒適性問題。

        一般而言,我們通過數(shù)值計算可以得到車體壁面上任意一點的壓力變化,然后根據(jù)車體表面壓力變化求出車廂內部的壓力變化,車外的壓力向車廂內部傳播主要取決于車體的剛度和車輛的氣密性。列車的剛度比較大,可以忽略其影響。因此,本文主要考慮車輛的氣密性的影響。

        假定車內部壓力變化率與內外壓差成正比:

        (3)

        式中:pi為車內部壓力;pε為外部壓力;τ為密封指數(shù)。由氣密性試驗測量得到:

        (4)

        根據(jù)以上公式,由前文計算所得到的列車車體表面的壓力變化,計算得到車廂內部的壓力變化,進而可以得到車廂內部的壓力變化率。

        表4為列車以120,140和160km/h的速度在隧道內單車運行時,車廂內部的壓力3s變化率,其中密封指數(shù)為τ=3,5和10s。

        表4列車車廂內部壓力變化率/(Pa/3s)

        Table4Internalpressurechangerateofthetrain

        速度/(km·h-1)120140160τ=0s方形隧道153921312816圓形隧道110615162004τ=3s方形隧道107114461825圓形隧道85511861556τ=5s方形隧道79610611328圓形隧道6428681105τ=10s方形隧道475628788圓形隧道387515653

        計算結果表明,車廂內部壓力變化最大來自頭車,當動態(tài)密封指數(shù)超過10s時,司機室內部的壓力變化率已經小于800/3Pa/s。因此,對于單車過隧道,動態(tài)密封指數(shù)大于10s即可。因為隧道橫截面不同,圓形隧道的阻塞比要比方形隧道的小,列車在過圓形隧道時要比過方形隧道乘客感覺更加舒適。

        3 噪聲水平分析研究

        噪聲會對人體的不良影響和危害,鐵路噪聲對沿線居民的影響最大。歐盟和日本等均出臺了高速列車噪聲聲壓級限制標準,隨著高速列車運行速度的提升,該噪聲標準的規(guī)定也越來越低[13]。因此有必要對列車的噪聲水平進行分析研究。本文結合大渦模擬與FW-H聲學模型分析了溫州S1線列車遠場氣動噪聲的分布規(guī)律,對關鍵位置的噪聲進行了評估。

        遠場噪聲測點在距軌道中心線25m和距軌面高度3.5m處,沿列車縱向平行等間距選取,在升弓和降弓位置也選取了兩個測點,測點編號從車頭前20mZ1開始依次到車尾后20mZ17。圖9和圖10分別給出了沿列車軸線方向聲壓級分布和最大值點聲壓頻譜分布。

        圖9 沿列車軸線方向聲壓級分布Fig.9 Sound pressure distribution in the axial direction of train

        圖10 最大值點聲壓頻譜分布Fig.10 Sound pressure spectral distribution of the Maximum point

        從圖中可以看出,溫州S1線列車以160km/h速度運行時,所有測點遠場聲壓級幅值分布比較均勻,集中分布在76-81dB之間,整體水平較低,噪聲水平最高位置處于升弓處。

        4 結論

        1)隨著速度的增大,列車的氣動阻力增大,車體表面壓力增大,車內壓力增大,噪聲增大。

        2)列車過隧道時,40m2圓形隧道要比33m2方形隧道氣動性能好;動態(tài)密封指數(shù)要大于10s。

        3)列車以160km/h速度運行時,線邊25m遠、沿列車軸線方向的最大噪聲級約80dB,在環(huán)境噪聲允許范圍內。

        [1] 袁磊,李人憲.高速列車氣動噪聲及影響巨[J].機械工程與自動化,2013,180(5):31-35.

        YUANLei,LIRenxian.Aerodynamicnoiseofhighspeedtrainanditsimpact[J].MechanicalEngineeringAutomaytion,2013,180(5):31-35.

        [2] 劉加利,張繼業(yè),張衛(wèi)華.考慮地面效應的高速列車遠場氣動噪聲計算方法研究[J].計算力學學報,2013,30(1); 94-100.

        LlUJiali,GHANUJiye,GHANUWeihua.Studyofcomputationmethodoffar-fieldaerodynamicnoiseofhigh-speedtrainconsideringgroundeffect[J].ChineseJournalofComputationalMechanics, 2013,30(1):94-100.

        [3] 田紅旗.列車空氣動力學[M].北京:中國鐵道出版社,2007:26-28,160-161.

        TIANHongqi.TrainAerodynamics[M].Beijing:ChinaRailwayPublishingHouse, 2007:26-28,160-161.

        [4]JOHND.Anderson.計算流體力學基礎及其應用[M].北京:機械工業(yè)出版社,2012:3.

        JOHND.Anderson.ComputationalFluidDynamic[M] .Beijing:ChinaMachinePress, 2012,3.

        [5] 梁習鋒,田紅旗.列車氣動性能評估參數(shù)研究[J].中國鐵道科學,2003,2(1):38-42.

        LIANGXifeng,TIANHongqi.Researchonevaluatingparametersoftrainaerodynamic[J].ChinaRailwayScience,2003,24(1):38-42.

        [6]TIANHongqi.Formationmechanismofaerodynamicdragofhigh-speedtrainandsomereductionmeasure[J].JournalofCentralSouthUniversity, 2009, 16(1):166-171.

        [7]JOSEPHAS.Aerodynamicofhigh-speedtrains[J].AnnualReviewofFluidMechanics, 2001, 33(1): 371-414.

        [8]Raghuathans,KimHD,SetoguchiT.Aerodynamicsofhigh-speedrailwaytrain[J].ProgressinAerospaceSciences,2002,38(1):469-514.

        [9] 肖云華,李行,熊小慧.明線運行時動車組空調裝置表面壓力分布數(shù)值分析[J].鐵道科學與工程學報,2013,10(3):108-111.

        XIAOYunhua,LIHang,XIONGXiaohui.NumericalsimulationofsurfacepressuredistributionofairconditioningequipmentwhenEMUrunninginopenair[J].JournalofRailwayScienceandEngineering. 2013,10(3):108-111.

        [10] 蔡軍爽,張潔,劉堂紅.基于DES的高速列車氣動阻力分布特性研究[J].鐵道科學與工程學報,2014,11(5):64-69.

        CAIJunshuang,ZHANGJie,LIUTanghong.Studyondistributioncharacteristicsofaerodynamicdragofhigh-speedtrainbasedonDES[J].JournalofRailwayScienceandEngineering,2014,11(5):64-69.

        [11] 姚拴寶,郭迪龍,楊國偉,等.高速列車氣動阻力分布特性研究[J].鐵道學報,2012,3(7):18-23.

        YAOShuanbao,GUODilong,YANGGuowei,etal.Distributionofhigh-speedtrainaerodynamicdrag[J].JournaloftheChinaRailwaySociety, 2012, 34(7):18-23.

        [12] 陳燕榮,肖友剛.高速列車空氣動力學性能計算[J].鐵道車輛,2009,47(1):14-16.

        CHENYanrong,XIAOYougang.Calculationinaerodynamicperformanceofhighspeedtrain[J].RollingStock,2009,47(1):14-16.

        [13]SuzukiM.Unsteadyaerodynamicforceactingonhighspeedtrainsintunnel[J].QRofRTRI, 2001,42(2): 89-93.

        [14] 梅元貴,趙海恒,劉應清.鐵路隧道壓力波數(shù)值分析[J].西南交通大學學報,1995,30 ( 6 ):667 -672.

        MEIYuangui,ZHAOHaiheng,LIUYingqing.Thenumericalpredictionofthepressuretransientgeneratedbyahighspeedtrainpassingatunnel[J].JournalofSouthwestJiaotongUniversity. 1995,30 (6):667-672.

        [15] 趙晶,李人憲.高速列車進入隧道的氣動作用數(shù)值模擬[J].西南交通大學學報,2009,44(1):96-100.

        ZHAOJing,LIRenxian.Numericalanalysisofaerodynamicsofhighspeedtrainsrunningintotunnels[J].JournalofSouthwestJiaotongUniversity. 2009,44 (1):96-100.

        [16] 馬大猷. 現(xiàn)代聲學理論基礎[M]. 北京:科學出版社, 2004: 296-306.

        MADayou.Modernacousticstheory[M].Beijing:SciencePress,2004:296-306.

        [17] 劉加利,張繼業(yè),張衛(wèi)華. 高速列車車頭的氣動噪聲數(shù)值分析[J]. 鐵道學報. 2011,33(09): 19-26.

        LIUJiali,ZHANGJiye,ZHANGWeihua.Numericalanalysisonaerodynamicnoiseofthehigh-speedtrainhead[J].JournaloftheCHINARailwaySociety. 2011,33(9):19-26.

        [18]BrancatiA,AliabadiMH.Boundaryelementsimulationsforlocalactivenoisecontrolusinganextendedvolume[J].EngineeringAnalysiswithBoundaryElements, 2012, 36(2): 190-202.

        [19]PapadopoulosChristosI.Developmentofanoptimized,standard-compliantproceduretocalculatesoundtransmissionloss:Numericalmeasurements[J].AppliedAcoustics, 2003, 64(11): 1069-1085.

        [20] 肖友剛,康志成.高速列車車頭曲面氣動噪聲的數(shù)值預測[J].中南大學學報(自然科學版), 2008, 39(6): 1267-1272.

        XIAOYougang,KANUZhicheng.Numericalpredi-ctionofaerodynamicnoiseradiatedfromhighspeedtrainheadsurface[J].JournalofCenterSouthUniversity(ScienceandTechnology),2008,39(6):1267-1272.

        Analysis and evaluation on characteristics of aerodynamics of intercity trains

        LIWenhua1,SHANGKeming1,YANGMingzhi2,3

        (1. CSR Qingdao Sifang Co. Ltd. , Qingdao 266111, China;2.SchoolofTransportationEngineering,CentralSouthUniversity,Changsha410075,China;3.KeyLaboratoryofTrafficSafetyonTrackofMinistryofEducation,Changsha410075,China)

        Thespeedofintercitytrainislimitedbetween100km/hand200km/h,buttherearefewaerodynamicstudiesforthisvelocitymagnitude.Withtheoperationoftheintercitytrain,thesafetyproblemisbecomingmoreandmoreimportant.So,itisnecessarytostudythisproblem.ThepaperstudiestheaerodynamiccharacteristicsofWenzhouS1-linetrainbynumericalsimulationmethod.Itmainlyanalyzedthedragandnoisewhenthetrainrunsinopenair.Thesurfacepressureandinsidepressurewhenthetraintravelsthroughthetunnelwerealsoinvestigated.Throughtheanalysisandcomparison,theregulationwasgiven:thedrag,surfacepressure,insidepressureandnoisewerelargerwiththeincreaseoftrainspeed;Theaerodynamicperformanceof40m2-circulartunnelisbetterthan33m2-squaretunnel;Thedynamicsealingindexshouldbegreaterthan10s;Trainsrunat160km/hspeed,lineedge25maway,alongtheaxialdirectionofthetrainthemaximumnoiselevelofabout80dB,thenoiseintheenvironmentwithintheallowablerange.

        numericalsimulationmethod;intercitytrain;aerodynamiccharacteristics

        2015-12-15

        國家自然科學基金資助項目(11372360)

        楊明智(1972-),男,湖南望城人,副教授,博士,從事列車空氣動力學研究;E-mail:yqyymz@126.com

        U270

        A

        1672-7029(2016)07-1407-07

        猜你喜歡
        湍流氣動阻力
        中寰氣動執(zhí)行機構
        鼻阻力測定在兒童OSA診療中的臨床作用
        基于NACA0030的波紋狀翼型氣動特性探索
        零阻力
        英語文摘(2020年10期)2020-11-26 08:12:12
        重氣瞬時泄漏擴散的湍流模型驗證
        基于反饋線性化的RLV氣動控制一體化設計
        測控技術(2018年9期)2018-11-25 07:44:24
        別讓摩擦成為學習的阻力
        阻力不小 推進當循序漸進
        “青春期”湍流中的智慧引渡(三)
        “青春期”湍流中的智慧引渡(二)
        一区二区三区国产天堂| 51国偷自产一区二区三区| 无码少妇a片一区二区三区| 一本一道波多野结衣av中文| 亚洲熟妇av日韩熟妇av| 国产在线视频一区二区三区| 在线a人片免费观看高清| 亚洲国内精品一区二区在线| 亚洲最新精品一区二区| 日韩有码中文字幕在线观看| 久久久噜噜噜久久| 日本高清aⅴ毛片免费| 亚洲精品永久在线观看| 国产爆乳无码一区二区在线| 亚洲日韩区在线电影| 亚洲国产成人资源在线桃色| 国产成人精品视频网站| 日本国产精品高清在线| 伊人久久这里只有精品| 三年片大全在线观看免费观看大全| 日本japanese丰满多毛| 人妻少妇精品专区性色av| 国产主播一区二区三区在线观看 | 麻豆婷婷狠狠色18禁久久| 九色九九九老阿姨| 日本做受高潮好舒服视频| 天码av无码一区二区三区四区| 国产成人亚洲综合小说区| 国产一区二区三区在线影院| 在线精品国产亚洲av麻豆| 人妻无码一区二区三区| 内射爽无广熟女亚洲| 久久精品国产一区二区电影| 女同性恋精品一区二区三区| 亚洲天堂一区二区三区视频| 国产av丝袜熟女丰满一区二区| 白白在线视频免费观看嘛| 性高朝大尺度少妇大屁股| 亚洲国产精品国自产电影| 久久久高清免费视频| 亚洲高清国产拍精品熟女|