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

        ?

        飛行過(guò)載下燃燒室凝相粒子沉積特征數(shù)值研究*

        2017-06-19 19:09:41劉長(zhǎng)猛余貞勇趙金萍
        固體火箭技術(shù) 2017年3期
        關(guān)鍵詞:發(fā)動(dòng)機(jī)

        劉長(zhǎng)猛, 余貞勇, 李 侃, 趙金萍

        (1.中國(guó)航天科技集團(tuán)公司四院四十一所,西安 710025; 2.中國(guó)航天科技集團(tuán)公司第四研究院,西安 710025)

        飛行過(guò)載下燃燒室凝相粒子沉積特征數(shù)值研究*

        劉長(zhǎng)猛1, 余貞勇2, 李 侃1, 趙金萍1

        (1.中國(guó)航天科技集團(tuán)公司四院四十一所,西安 710025; 2.中國(guó)航天科技集團(tuán)公司第四研究院,西安 710025)

        固體發(fā)動(dòng)機(jī)燃燒室部分凝相粒子在飛行過(guò)載作用下產(chǎn)生沉積,嚴(yán)重影響發(fā)動(dòng)機(jī)工作性能。通過(guò)確定燃燒室粒子參數(shù)和建立燃燒室內(nèi)兩相流場(chǎng)數(shù)值方法,獲得了發(fā)動(dòng)機(jī)不同軸向過(guò)載下粒子運(yùn)動(dòng)及沉積規(guī)律。與試驗(yàn)數(shù)據(jù)進(jìn)行對(duì)比分析,驗(yàn)證了計(jì)算方法的準(zhǔn)確性。數(shù)值結(jié)果表明,隨著軸向過(guò)載增大,后封頭及噴管潛入段粒子沉積質(zhì)量逐漸增大。沉積粒子粒徑大于凝相粒子平均值,即粒子粒徑越大,越容易沉積。軸向過(guò)載增大,減小了粒子在發(fā)動(dòng)機(jī)內(nèi)部的駐留時(shí)間,凝相粒子平均駐留時(shí)間均大于0.13 s。

        固體火箭發(fā)動(dòng)機(jī);飛行過(guò)載;粒子沉積

        0 引言

        發(fā)動(dòng)機(jī)燃燒室流動(dòng)是典型的多相流場(chǎng),在較大飛行過(guò)載條件下,凝相粒子存在更加復(fù)雜的運(yùn)動(dòng)過(guò)程、易造成粒子沉積和絕熱層燒蝕異常等現(xiàn)象。因此,研究和掌握燃燒室內(nèi)流場(chǎng)粒子流動(dòng)規(guī)律,對(duì)于帶有潛入式噴管的大型固體火箭發(fā)動(dòng)機(jī)設(shè)計(jì)具有重要意義。國(guó)內(nèi)外進(jìn)行了許多燃燒室凝相粒子運(yùn)動(dòng)規(guī)律的研究,王國(guó)輝[1]和李越森[2]等分別對(duì)過(guò)載狀態(tài)下固體火箭發(fā)動(dòng)機(jī)燃燒室內(nèi)二相流動(dòng)進(jìn)行了數(shù)值研究,獲得了過(guò)載對(duì)粒子運(yùn)動(dòng)的影響。魏超[3-4]探求固體火箭發(fā)動(dòng)機(jī)潛入噴管背壁區(qū)中粒子沉積的機(jī)理和預(yù)估方法,確定了燃燒室中粒子的直徑分布和兩種粒子捕獲判據(jù)。Andrew等[5]進(jìn)行了發(fā)動(dòng)機(jī)粒子沉積試驗(yàn)和數(shù)值計(jì)算,得到了粒子沉積分布規(guī)律。田維平等[6]總結(jié)分析了國(guó)內(nèi)外對(duì)飛行過(guò)載下固體火箭發(fā)動(dòng)機(jī)中出現(xiàn)的絕熱層燒蝕問(wèn)題的研究方法,詳細(xì)闡述了燃燒室粒子粒度參數(shù)確定方法、過(guò)載流場(chǎng)數(shù)值模擬方法及地面模擬過(guò)載試驗(yàn)方法等方面研究進(jìn)展。許團(tuán)委等[7]分析了戰(zhàn)術(shù)發(fā)動(dòng)機(jī)燃燒室內(nèi)凝相顆粒沖刷參數(shù)分布,討論了顆粒聚集狀態(tài)與絕熱層燒蝕之間的關(guān)系。陳福振等[8]針對(duì)φ315 mm發(fā)動(dòng)機(jī)工作狀態(tài)和結(jié)構(gòu)特點(diǎn),進(jìn)行了發(fā)動(dòng)機(jī)燃燒室內(nèi)過(guò)載條件下的氣粒兩相流動(dòng)數(shù)值模擬。夏勝勇[9]數(shù)值研究了Al2O3液滴尺寸、尺寸比及環(huán)境壓強(qiáng)對(duì)液滴碰撞結(jié)果的影響,結(jié)合理論分析,獲得了碰撞結(jié)果的影響規(guī)律。

        目前,對(duì)于實(shí)體發(fā)動(dòng)機(jī)粒子沉積規(guī)律以及定量計(jì)算,并沒有準(zhǔn)確且合理的計(jì)算方法,同時(shí)針對(duì)大型高能量發(fā)動(dòng)機(jī)的研究較少。

        本文針對(duì)小長(zhǎng)細(xì)比潛入式噴管高性能固體火箭發(fā)動(dòng)機(jī)燃燒室內(nèi)部?jī)上嗔鲌?chǎng)建立數(shù)值模型,基于10g軸向飛行過(guò)載,將數(shù)值計(jì)算結(jié)果與試驗(yàn)結(jié)果進(jìn)行了對(duì)比驗(yàn)證,并分析了不同軸向過(guò)載對(duì)固體火箭發(fā)動(dòng)機(jī)凝相粒子沉積的影響。

        1 數(shù)值模型

        基于標(biāo)準(zhǔn)k-ε模型建立數(shù)值模型。湍流動(dòng)能k和耗散率ε的輸運(yùn)方程分別為

        (1)

        (2)

        基于布拉修斯假設(shè),湍流動(dòng)能生成項(xiàng)Gk為

        (3)

        平均張力張量比率的系數(shù)S定義為

        (4)

        通過(guò)計(jì)算k和ε,湍流粘度μt為

        (5)

        在k方程中,耗散項(xiàng)是Yk=ρε。在ε方程中,生成項(xiàng)和耗散項(xiàng)為

        (6)

        本文模型常數(shù)選擇是Cμ=0.09,Clε=1.44,C2ε=1.92,σk=1.0,σε=1.3。

        2 數(shù)值計(jì)算

        2.1 計(jì)算區(qū)域和網(wǎng)格

        以某發(fā)動(dòng)機(jī)燃燒室為研究對(duì)象,劃分結(jié)構(gòu)化網(wǎng)格結(jié)構(gòu),如圖1所示。

        圖1中,網(wǎng)格數(shù)量為141×105個(gè),所有壁面處均進(jìn)行加密,喉部和后封頭沉積部位第一層網(wǎng)格高度為3 mm。劃分了4種不同數(shù)量的網(wǎng)格結(jié)構(gòu),壁面處網(wǎng)格基本一致,而在流動(dòng)核心區(qū)采用了不同的網(wǎng)格密度,網(wǎng)格數(shù)量分別為188×105、141×105、105×105、63×105。

        2.2 計(jì)算設(shè)置及方法

        把氣相作為連續(xù)相,凝相粒子作為離散相進(jìn)行處理,鑒于本文主要研究發(fā)動(dòng)機(jī)中氣相及凝相粒子在過(guò)載條件下的流動(dòng)特性,所以將復(fù)雜的兩相流動(dòng)模型進(jìn)行了適當(dāng)簡(jiǎn)化:不考慮氣相化學(xué)反應(yīng);不考慮凝相粒子的燃燒、蒸發(fā)、破碎過(guò)程。根據(jù)計(jì)算經(jīng)驗(yàn)和文獻(xiàn)驗(yàn)證[4-5],具體參數(shù)設(shè)置見表1。

        表1 計(jì)算參數(shù)設(shè)置

        入口采用壓力入口:8.2 MPa;出口為壓強(qiáng)出口:10 Pa。壁面為無(wú)滑移壁面,采用標(biāo)準(zhǔn)壁面函數(shù)。設(shè)置粒子的主要沉積位置為發(fā)動(dòng)機(jī)后封頭和噴管潛入段壁面。

        燃面和其他壁面設(shè)置為反射邊界條件,彈性恢復(fù)系數(shù)為

        (7)

        (8)

        其中,n、t分別為法向和切向速度分量,下標(biāo)1、2分別表示碰撞前后的量。本文恢復(fù)系數(shù)設(shè)置為0.9。

        沉積模型中顆粒的磨蝕和沉積速率分別定義為

        (9)

        (10)

        式中C(dp)為顆粒直徑的函數(shù);α為顆粒對(duì)壁面的沖擊角;f(a)為沖擊角的函數(shù);v為顆粒相對(duì)于壁面的速度;b(v)為此相對(duì)速度的函數(shù)。本文設(shè)置為C=1.8×10-9,f=1,b=0。

        采用標(biāo)準(zhǔn)k-ε模型,基于密度求解器進(jìn)行計(jì)算,二階迎風(fēng)格式進(jìn)行離散。首先進(jìn)行燃?xì)饬鲌?chǎng)穩(wěn)態(tài)計(jì)算,然后添加離散項(xiàng)進(jìn)行兩相流場(chǎng)計(jì)算,對(duì)后封頭等特殊位置進(jìn)行離散項(xiàng)流量監(jiān)測(cè),判斷計(jì)算收斂。

        本文研究流場(chǎng)中的離散相是為了計(jì)算粒子沉積情況。故在本文中,將隨流性較好的小粒徑煙塵粒子的質(zhì)量歸于連續(xù)相,不對(duì)其進(jìn)行計(jì)算分析,僅研究大粒徑顆粒在發(fā)動(dòng)機(jī)中的運(yùn)動(dòng)狀況。

        Rosin-Rammler分布假定在粒子直徑d與大于此直徑的粒子質(zhì)量分?jǐn)?shù)Yd之間存在指數(shù)關(guān)系:

        (11)

        式中的d為平均直徑(中位徑),n為分布指數(shù)。設(shè)大粒徑粒子占所有粒子質(zhì)量的20%。對(duì)于高性能推進(jìn)劑,目前沒有太多資料顯示粒徑分布。結(jié)合文獻(xiàn)[10]資料,并經(jīng)過(guò)多次計(jì)算修正,選取離散相粒子平均直徑Dm為45 μm,取分布指數(shù)為3.5,粒子直徑的分布曲線見圖2。計(jì)算中,將粒子直徑離散為6組。

        粒子密度對(duì)于計(jì)算粒子軌跡和粒子沉積非常關(guān)鍵,沉積物重量與粒子密度成比例[4]。粒子密度隨溫度變化很大,由以下經(jīng)驗(yàn)公式進(jìn)行計(jì)算:

        ρ=5 632-1.127T

        (12)

        式中T為燃燒室絕熱火焰溫度。

        推進(jìn)劑燃燒絕熱火焰溫度T取為3 600 K,由式(12)計(jì)算,ρ=1 575 kg/m3。

        軸向過(guò)載設(shè)置為1、3、5、8、10g,橫向過(guò)載均為0.5g進(jìn)行計(jì)算。

        3 計(jì)算算例及結(jié)果分析

        3.1 數(shù)值驗(yàn)證

        對(duì)不同數(shù)量的網(wǎng)格結(jié)構(gòu)進(jìn)行網(wǎng)格無(wú)關(guān)性驗(yàn)證,發(fā)現(xiàn)燃燒室中高溫高壓下,燃?xì)饬鲌?chǎng)對(duì)網(wǎng)格結(jié)構(gòu)并不十分敏感,而噴管流動(dòng)對(duì)網(wǎng)格分辨率要求更為嚴(yán)格。因此,在燃燒室壁面處進(jìn)行Y+驗(yàn)證,保證Y+值在30~200之間,并加密了噴管處網(wǎng)格。對(duì)于粒子運(yùn)動(dòng)的捕捉,63×105網(wǎng)格數(shù)量太少,導(dǎo)致結(jié)果與其他3種存在明顯差異,其他3種網(wǎng)格基本一致。為了獲得較好的計(jì)算精度,綜合分析后,選擇141×105網(wǎng)格進(jìn)行后續(xù)計(jì)算。

        圖3顯示的是軸向過(guò)載和沉積質(zhì)量的關(guān)系曲線。后封頭沉積質(zhì)量流率和軸向過(guò)載加速度成正比,在8g之前,隨著軸向加速度的增大而逐漸線性增大,8g和10g增長(zhǎng)速度減緩。在軸向過(guò)載為10g時(shí),沉積質(zhì)量流率達(dá)到了0.03 kg/(m2·s)。由發(fā)動(dòng)機(jī)工作中不同時(shí)刻的過(guò)載,得到了粒子沉積質(zhì)量流率隨時(shí)間的變化規(guī)律。

        將不同時(shí)刻沉積質(zhì)量流率和沉積位置暴露面積(包括后封頭和噴管潛入段)相乘,得到不同時(shí)刻沉積質(zhì)量,進(jìn)而計(jì)算出發(fā)動(dòng)機(jī)工作過(guò)程中的總沉積質(zhì)量規(guī)律,如圖4所示。對(duì)數(shù)值結(jié)果積分后,得到總沉積質(zhì)量為3.3 kg。發(fā)動(dòng)機(jī)飛行試驗(yàn)后,將收集到的全部殘?jiān)M(jìn)行分離稱重,通過(guò)元素分析法,得到Al2O3殘?jiān)|(zhì)量為3.02 kg。對(duì)于沉積質(zhì)量的計(jì)算,數(shù)值結(jié)果與實(shí)測(cè)值吻合較好,數(shù)值模型得到初步驗(yàn)證。由于缺少不同位置不同時(shí)刻粒子沉積量的實(shí)測(cè)值,因此暫無(wú)法做全面的驗(yàn)證分析。

        3.2 不同軸向過(guò)載沉積規(guī)律

        由圖5中不同過(guò)載(軸向1、5、10g)后封頭沉積率分布云圖可知,軸向過(guò)載越大,后封頭沉積質(zhì)量流率越大,而潛入段粒子沉積量反而變小。1g軸向過(guò)載時(shí),主要的粒子沉積區(qū)域位于噴管潛入段,潛入段沉積率遠(yuǎn)大于后封頭沉積率。在潛入段和后封頭過(guò)渡區(qū),同樣沉積較多的粒子,后封頭的大部分區(qū)域沉積并不明顯,僅在靠近潛入段處沉積較嚴(yán)重。隨著軸向過(guò)載增大,后封頭區(qū)域均出現(xiàn)粒子沉積現(xiàn)象,分布范圍較大。此時(shí),噴管潛入段沉積率有所減小,但仍大于后封頭。粒子主要沉積在潛入段防熱環(huán)部位,此處的沉積率遠(yuǎn)高于其他位置。因此,軸向過(guò)載增大,加重了后封頭區(qū)域的粒子沉積。

        圖6為后封頭燃?xì)饧傲W舆\(yùn)動(dòng)軌跡示意圖。圖6中,流線為燃?xì)膺\(yùn)動(dòng)軌跡,燃?xì)膺\(yùn)動(dòng)至防熱環(huán)位置處,并需要在此處做較大角度的偏轉(zhuǎn)繞過(guò)噴管潛入段進(jìn)入喉部。如圖6中箭頭所示,對(duì)于流動(dòng)性較差的凝相顆粒而言,在慣性力的作用下,部分凝相粒子較易沉積在噴管潛入段位置。防熱環(huán)形成的凹槽結(jié)構(gòu),更增大了粒子沉積的可能性。因此,在噴管潛入段的防熱環(huán)位置,凝相粒子沉積較為嚴(yán)重。

        圖7為不同過(guò)載燃燒室渦量分布圖。

        在軸向過(guò)載較小時(shí),由于橫向過(guò)載的影響,在發(fā)動(dòng)機(jī)Y軸負(fù)向燃面渦量較大。分析認(rèn)為,是由于Y軸正向燃面的粒子和Y軸負(fù)向粒子在此處相互影響,對(duì)燃?xì)饩哂袛_動(dòng)作用。隨著軸向過(guò)載加大,粒子加速噴出,此位置處的擾動(dòng)效果減弱,橫向過(guò)載對(duì)凝相粒子影響逐漸減小。受到橫向過(guò)載影響,Y軸負(fù)向粒子渦量值大于對(duì)應(yīng)位置處Y軸正向渦量值。隨著軸向過(guò)載的增大,燃燒室流場(chǎng)主要受到軸向過(guò)載的影響,但仍存在橫向過(guò)載對(duì)粒子運(yùn)動(dòng)的影響。

        圖8為后封頭和潛入段在不同軸向過(guò)載時(shí)平均沉積質(zhì)量流率的變化情況。由圖8可知,潛入段沉積質(zhì)量流率大于后封頭,后封頭和潛入段平均值分別為0.012 kg/(m2·s)和0.085 kg/(m2·s)。隨著過(guò)載增大,后封頭沉積質(zhì)量流率逐漸增大,而潛入段逐漸減小,圖9為不同過(guò)載下沉積平均粒徑分布曲線,圖10為粒子在燃燒室和噴管中的平均駐留時(shí)間。

        由圖9可知,后封頭沉積粒子的粒徑大于潛入段粒徑。對(duì)于后封頭,隨著過(guò)載增大,粒子平均粒徑逐漸減小。分析認(rèn)為,主要是由于小過(guò)載下,大粒子隨流性差,較易沉積于后封頭,此時(shí)平均粒徑較大。隨著軸向過(guò)載加大,小粒徑顆粒沉積逐漸增多。因此,出現(xiàn)平均粒徑減小現(xiàn)象。不同過(guò)載情況下,沉積粒子的平均粒徑均大于輸入的平均粒徑45 μm,說(shuō)明粒徑越大,越容易沉積,軸向過(guò)載對(duì)隨流性較差的大粒徑顆粒影響更大。

        由圖10可見,軸向過(guò)載增大,加大了粒子運(yùn)動(dòng)速度,減小了粒子在發(fā)動(dòng)機(jī)內(nèi)部的駐留時(shí)間。由圖10可知,燃燒室內(nèi)部粒子最小的駐留時(shí)間為0.13 s,而鋁凝團(tuán)粒子燃燒完全所需要的時(shí)間約為0.013~0.017 s。因此,數(shù)值計(jì)算結(jié)果表明,10g過(guò)載下的粒子駐留時(shí)間能滿足推進(jìn)劑充分燃燒。

        4 結(jié)論

        (1)隨著軸向過(guò)載變大,后封頭沉積質(zhì)量流率增大,而潛入段粒子沉積量反而變小,后封頭的沉積率小于潛入段。

        (2)受噴管潛入結(jié)構(gòu)的影響,凝相粒子主要沉積在潛入段防熱環(huán)部位,此處的沉積率遠(yuǎn)高于其他位置。

        (3)不同過(guò)載情況下,沉積粒子的平均粒徑均大于輸入的平均粒徑45 μm,說(shuō)明粒徑越大,越易沉積,后封頭沉積粒子的粒徑大于潛入段粒徑。

        (4)軸向過(guò)載增大,加大了粒子運(yùn)動(dòng)速度,減小了粒子在發(fā)動(dòng)機(jī)內(nèi)部的駐留時(shí)間。最小的平均駐留時(shí)間為0.13 s,能夠保證凝相粒子充分燃燒。

        [1] 王國(guó)輝,何國(guó)強(qiáng),劉佩進(jìn),等.過(guò)載狀態(tài)下固體火箭發(fā)動(dòng)機(jī)燃燒室內(nèi)二相流動(dòng)數(shù)值模擬[J].固體火箭技術(shù),2001,24(2):8-12.

        [2] 李越森,葉定友.高過(guò)載下固體發(fā)動(dòng)機(jī)內(nèi)Al2O3粒子運(yùn)動(dòng)狀況的數(shù)值模擬[J].固體火箭技術(shù),2008,31(1):24-27.

        [3] 魏超,侯曉.潛入噴管背壁區(qū)粒子沉積的機(jī)理分析與數(shù)值模擬[J].航空動(dòng)力學(xué)報(bào),2006,21(6):1109-1115.

        [4] 魏超.固體火箭發(fā)動(dòng)機(jī)潛入與非潛入噴管流場(chǎng)與粒子沉積研究[D].西安:航天動(dòng)力技術(shù)研究院,2005.

        [5] Andrew C Cortopassi,Eric Boyer and Kenneth K Kuo.A subscale solid rocket motor for characterization of submerged nozzle erosion[R].AIAA 2009-5172.

        [6] 田維平,許團(tuán)委,王建儒.過(guò)載下燃燒室粒子特性與絕熱層燒蝕研究進(jìn)展[J].固體火箭技術(shù),2015,38(1):30-36.

        [7] 許團(tuán)委,田維平,王建儒.中小過(guò)載下戰(zhàn)術(shù)發(fā)動(dòng)機(jī)內(nèi)流場(chǎng)數(shù)值模擬[J].推進(jìn)技術(shù),2015,36(4):532-539.

        [8] 陳福振,強(qiáng)洪夫,高巍然,等.固體火箭發(fā)動(dòng)機(jī)內(nèi)氣粒兩相流動(dòng)的SPH-FVM耦合方法數(shù)值模擬[J].推進(jìn)技術(shù),2015,36(2):175-185.

        [9] 夏勝勇.三氧化二鋁液滴碰撞機(jī)理及模型研究[D].西安:西北工業(yè)大學(xué),2015.

        [10] 李強(qiáng),甘曉松,劉佩進(jìn),等.大型固體發(fā)動(dòng)機(jī)潛入式噴管背壁區(qū)域熔渣沉積數(shù)值模擬[J].固體火箭技術(shù),2010,33(2):148-151.

        (編輯:崔賢彬)

        Simulation research for accretion characteristic of condensed phase particle in SRM chamber with flight overload

        LIU Chang-meng1,YU Zhen-yong2,LI Kan1,ZHAO Jin-ping1

        (1.The 41st Institute of the Fourth Academy of CSAC, Xi'an 710025,China; 2.The Fourth Academy of CSAC,Xi'an 710025,China)

        The partially accretion of condensed particles in combustion chamber seriously affects the performance of the rocket motor under the overloaded action of flight.By determining the particle parameters and establishing the numerical method of two-phase flow field in combustion chamber,the particle movement and accretion law was obtained with the different axial overloads.The accuracy of the calculation method was verified by comparing with the experimental data.The numerical results show that with the increase of the axial overload,the accretion quality of the particle is gradually increased in the rear dome and the entrance region of submerged nozzle.The size of accretion particles is greater than the average value and the particle with larger particle size is much easier to be deposited.The residence time of the particles is reduced with the increases of the axial overload in the SRM.It is greater than 0.13 s for the average residence time of the condensed phase particle.

        solid rocket motor;flight overload;particle accretion

        2016-03-21;

        2016-04-11。

        劉長(zhǎng)猛(1986—),男,博士,研究方向?yàn)榘l(fā)動(dòng)機(jī)總體設(shè)計(jì)。E-mail:liuchangmeng2126@126.com

        V435

        A

        1006-2793(2017)03-0302-05

        10.7673/j.issn.1006-2793.2017.03.006

        猜你喜歡
        發(fā)動(dòng)機(jī)
        元征X-431實(shí)測(cè):奔馳發(fā)動(dòng)機(jī)編程
        2015款寶馬525Li行駛中發(fā)動(dòng)機(jī)熄火
        2012年奔馳S600發(fā)動(dòng)機(jī)故障燈偶爾點(diǎn)亮
        發(fā)動(dòng)機(jī)空中起動(dòng)包線擴(kuò)展試飛組織與實(shí)施
        RD-180超級(jí)火箭發(fā)動(dòng)機(jī)的興衰
        太空探索(2016年8期)2016-07-10 09:21:58
        奔馳E200車發(fā)動(dòng)機(jī)故障燈常亮
        奔馳E260冷車時(shí)發(fā)動(dòng)機(jī)抖動(dòng)
        新一代MTU2000發(fā)動(dòng)機(jī)系列
        2013年車用發(fā)動(dòng)機(jī)排放控制回顧(下)
        2013年車用發(fā)動(dòng)機(jī)排放控制回顧(上)
        青青手机在线视频观看| 天躁夜夜躁狼狠躁| 亚洲色图视频在线| 亚洲av熟女天堂系列| 久久中文字幕一区二区| 国产精品毛片一区二区三区| 国产主播一区二区三区在线观看 | 国产国语按摩对白av在线观看 | 国产成人啪精品午夜网站| 国产精品视频免费一区二区三区| 久久国产精品婷婷激情| 欧美寡妇xxxx黑人猛交| 一级二级中文字幕在线视频| 风韵丰满妇啪啪区老老熟女杏吧| 自拍偷拍 视频一区二区| 欧美人妻日韩精品| 性欧美大战久久久久久久久| 免费一级黄色大片久久久| 在线观看一区二区三区国产| 乱色欧美激惰| 一本大道东京热无码| 玖玖资源站无码专区| 久久国产精品一区av瑜伽| 久久国产亚洲高清观看| 九九在线视频| 蜜桃av一区二区三区久久| 国产精品亚洲专区无码不卡| 51久久国产露脸精品国产| 囯产精品无码一区二区三区AV| 91精品国产九色综合久久香蕉 | 久久亚洲色www成人欧美| 久久青草国产精品一区| 日本免费看一区二区三区| 国产不卡精品一区二区三区| 三年在线观看免费大全下载| 无码一区二区三区不卡AV| 中文字幕人乱码中文字幕乱码在线| 777精品出轨人妻国产| 免费的成年私人影院网站| 国产精品麻豆A啊在线观看 | 人人爽人人爱|