999精品在线视频,手机成人午夜在线视频,久久不卡国产精品无码,中日无码在线观看,成人av手机在线观看,日韩精品亚洲一区中文字幕,亚洲av无码人妻,四虎国产在线观看 ?

耗氧型惰化反應(yīng)器起燃特性

2021-03-27 04:48:34馮詩(shī)愚任童謝輝輝彭孝天王晨臣王洋洋潘俊
航空學(xué)報(bào) 2021年3期
關(guān)鍵詞:模型

馮詩(shī)愚, 任童, 謝輝輝, 彭孝天,王晨臣,王洋洋,潘俊

1. 南京航空航天大學(xué) 航空學(xué)院 飛行器環(huán)境控制與生命保障工業(yè)和信息化部重點(diǎn)實(shí)驗(yàn)室, 南京 210016 2. 航空工業(yè)江西洪都航空工業(yè)集團(tuán)有限責(zé)任公司,南昌 330024 3. 中國(guó)航空工業(yè)集團(tuán)有限公司 南京機(jī)電液壓工程研究中心航空機(jī)電系統(tǒng)綜合航空科技重點(diǎn)實(shí)驗(yàn)室,南京 211106

飛機(jī)油箱上部空間充滿易燃易爆的燃油蒸汽,加之常將燃油作為熱沉,增加了油箱的危險(xiǎn)性[1]。惰化技術(shù)將油箱上部空間氧含量降低至安全值內(nèi),減少燃爆危害,已廣泛應(yīng)用于各類軍民飛機(jī)上[2-4]。目前利用機(jī)載設(shè)備制造惰性氣體或者直接消耗燃油箱內(nèi)氧氣來惰化油箱成為惰化技術(shù)發(fā)展的新趨勢(shì),由此產(chǎn)生了分子篩式[5]、滲透膜式[6]、耗氧式[7]3類主流的機(jī)載惰化技術(shù)。2004年,美國(guó)Phyre公司率先應(yīng)用此技術(shù)為空軍研制了一種耗氧型惰化系統(tǒng)[8-10],相較于其他2種方法,該系統(tǒng)無需從發(fā)動(dòng)機(jī)引氣,除低壓氣泵外,無任何運(yùn)動(dòng)部件,具有功耗低、結(jié)構(gòu)緊湊,質(zhì)量輕等優(yōu)點(diǎn),目前中國(guó)也對(duì)其展開了研究[11]。

所謂耗氧型惰化技術(shù)是將油箱上部的氣相空間內(nèi)氣體與燃油蒸汽混合物引出后,通入一個(gè)催化燃燒反應(yīng)器中,在催化劑的作用下燃油蒸汽與氧氣發(fā)生化學(xué)反應(yīng)分解成水蒸氣與二氧化碳,反應(yīng)器出口的低氧含量的惰性氣體回充入燃油箱中,以此達(dá)到惰化目的。

在惰化系統(tǒng)啟動(dòng)過程中,需要對(duì)抽出的氣體進(jìn)行加熱,然后利用加熱后的氣體對(duì)冷態(tài)的催化床加熱,達(dá)到化學(xué)反應(yīng)所需的溫度,從而進(jìn)入工作狀態(tài),這一過程稱為反應(yīng)器的起燃。顯然反應(yīng)器起燃過程的快慢決定了耗氧型惰化系統(tǒng)的啟動(dòng)速度。由于耗氧型惰化系統(tǒng)重量的限制,對(duì)預(yù)熱氣體的溫度與速度都有限制,需要得到各種工作狀態(tài)對(duì)反應(yīng)器起燃特性的影響。因此反應(yīng)器起燃特性的研究對(duì)于耗氧型惰化系統(tǒng)的設(shè)計(jì)具有重要意義。

反應(yīng)器起燃過程是一個(gè)耦合化學(xué)反應(yīng)與傳質(zhì)傳熱的過程,利用CFD可對(duì)其過程求解,對(duì)反應(yīng)器內(nèi)復(fù)雜的傳質(zhì)傳熱過程進(jìn)行預(yù)測(cè)。早期的仿真大部分是一維模型,雖然對(duì)反應(yīng)器出口產(chǎn)物濃度有較好的預(yù)測(cè),但是無法得到反應(yīng)器軸向與徑向的參數(shù)分布。近年來,對(duì)反應(yīng)仿真主要有2種方法:一種是對(duì)催化床管道內(nèi)的控制方程直接求解,可得到反應(yīng)器內(nèi)微觀的傳熱傳質(zhì)情況,但該方法需要大量網(wǎng)格及復(fù)雜的邊界條件。另一種是將催化床簡(jiǎn)化為多孔介質(zhì),減少計(jì)算量,例如Wang等[12]將多孔介質(zhì)兩溫度模型應(yīng)用于甲烷與水蒸氣重整泡沫陶瓷型反應(yīng)器模型,得到了氣相與固相內(nèi)溫度云圖分布。Gao等[13]建立了一個(gè)耦合CO-DEO(Carbon monoxide-Diethyl Oxalate)反應(yīng)的雙溫度多孔介質(zhì)固定床反應(yīng)器模型,得到了軸向與徑向的組分濃度、氣固相溫度分布,并對(duì)模型進(jìn)行了驗(yàn)證,將仿真結(jié)果與徐艷等[14]早先進(jìn)行的實(shí)驗(yàn)數(shù)據(jù)對(duì)比,結(jié)果顯示多孔介質(zhì)兩溫度模型所測(cè)得的仿真結(jié)果與實(shí)驗(yàn)數(shù)據(jù)較為吻合。龍回龍等[15]采用多孔介質(zhì)兩溫度模型建立了耦合水氣變換反應(yīng)的反應(yīng)器模型,將仿真數(shù)據(jù)與雙溫度模型進(jìn)行對(duì)比,發(fā)現(xiàn)雙溫度模型更貼合實(shí)際情況。綜上所述,多孔介質(zhì)雙溫度模型在反應(yīng)器仿真中廣泛應(yīng)用且模型發(fā)展較為成熟,與實(shí)驗(yàn)結(jié)果貼合度高。

與化工反應(yīng)器不同的是,耗氧型惰化技術(shù)的目的是降低燃油箱的燃爆可能,反應(yīng)器設(shè)計(jì)中需考慮出口溫度等因素,追求轉(zhuǎn)化率則被放在次要位置上。由于保密原因,該技術(shù)的公開資料非常少,且國(guó)產(chǎn)燃油RP-3的碳?xì)湮餄舛扰c國(guó)外燃油也有較大差異。因此本文將采用多孔介質(zhì)兩溫度模型,選取RP-3燃油,對(duì)反應(yīng)器的起燃過程進(jìn)行仿真研究。

1 仿真模型及參數(shù)選取

1.1 雙溫度模型

反應(yīng)器結(jié)構(gòu)如圖1所示,x為軸向,r為徑向,中間為由150 mL催化劑顆粒填充的催化床,燃油蒸汽混合物通過加熱器后進(jìn)入此區(qū)域發(fā)生化學(xué)反應(yīng),并釋放熱量。催化床由多孔介質(zhì)模型描述,并認(rèn)為催化床是均勻的。

圖1 反應(yīng)器結(jié)構(gòu)圖Fig.1 Reactor structure diagram

選取燃油為國(guó)產(chǎn)大慶RP-3,其氫質(zhì)量分?jǐn)?shù)為14.46%,分子式可簡(jiǎn)化為C7.12H14.46,摩爾質(zhì)量為141 g/moL, 燃燒熱值為43.53 MJ/kg,反應(yīng)方程為

C7.12H14.46+10.74O2=7.12CO2+7.23H2O

(1)

化學(xué)反應(yīng)動(dòng)力學(xué)方程由實(shí)驗(yàn)測(cè)得[16]:

(2)

式中:RRP-3為RP-3的化學(xué)反應(yīng)速率;cRP-3為RP-3濃度;cO2為氧氣的濃度;φO2為氧氣摩爾分?jǐn)?shù);R為氣體摩爾常數(shù);T為溫度。

化學(xué)反應(yīng)以UDF(User Defined Function)源項(xiàng)形式加入,以表示化學(xué)反應(yīng)中各組分的生成與消耗,并在多孔介質(zhì)區(qū)域通過UDS(User Defined Scalar)加入固相能量方程。由于催化反應(yīng)發(fā)生在催化劑顆粒表面,反應(yīng)熱量大部分由固相吸收,因此將反應(yīng)熱以能量源項(xiàng)形式加入的固相能量方程中,認(rèn)為化學(xué)反應(yīng)熱量全部由固相吸收,通過傳熱傳遞給氣相。催化劑顆粒密度ρs、催化劑比熱容cs、導(dǎo)熱系數(shù)кs、直徑dp、催化床孔隙率γ參數(shù)見表1。

表1 催化床中主要參數(shù)Table 1 Main parameters of reactor catalytic bed

1.2 控制方程與邊界條件

在多介質(zhì)區(qū)域中,存在氣固兩相,其內(nèi)部的控制方程如下所示。

氣相質(zhì)量守恒:

(3)

式中:v為氣體速度;ρg為氣相密度。

氣相動(dòng)量方程:

(4)

式中:τ為剪切應(yīng)力;p為氣相壓力;S為動(dòng)量源項(xiàng)。

(5)

式中:μg為氣體動(dòng)力黏度。

氣相狀態(tài)方程:

(6)

式中:Yi為氣相第i組分質(zhì)量分?jǐn)?shù);Mi為氣相第i組分摩爾質(zhì)量;Tg為氣體溫度;N為組分總數(shù)。

組分運(yùn)輸方程:

(7)

式中:Ri為第i組分的質(zhì)量源項(xiàng),即化學(xué)反應(yīng)中各組分的生成消耗量;Ji為氣相第i組分的質(zhì)量擴(kuò)散。

氣相能量方程

(8)

式中:hv氣固換熱系數(shù);Ts為固體溫度;Eg為氣體能量;κg為氣相導(dǎo)熱系數(shù);cp,g為氣相比定壓熱容;hi為氣相第i組分焓變。

由UDS添加的固相能量方程為

hv(Ts-Tg)+Φr

(9)

式中:Φr為反應(yīng)釋放的熱量,通過UDF添加在固相能量方程中。

氣固相傳熱系數(shù)[17]為

(10)

(11)

Pr=cp,gμg/κg

(12)

Re=ρg|v|Dbμg

(13)

式中:Db為催化劑管徑;Re為雷諾數(shù);Nu為努塞爾數(shù);Pr為普朗特?cái)?shù)。

Guardo等[18]使用Spalart-Allmaras、k-ω、標(biāo)準(zhǔn)k-ε、Realizablek-ε和RNGk-ε等5種模型來研究湍流模型對(duì)固定床壁面?zhèn)鳠岬挠绊?表明使用Spalart-Allmaras(S-A)模型可以更好地反映氣體在復(fù)雜床層間的流動(dòng)情況,即

(14)

在入口處邊界條件為

u=uin,v=0,T=Tin

(15)

式中:u和v分別為軸向和徑向速度;uin為反應(yīng)器入口氣體定常初始速度;Tin為氣體被加熱器預(yù)熱后的溫度。

在出口處邊界為

p=0

(16)

反應(yīng)器壁面熱量損失為

q=hc(to-tw)

(17)

式中:hc為表面?zhèn)鳠嵯禂?shù);to為反應(yīng)器外界氣體溫度,設(shè)置為30 ℃;tw為反應(yīng)器壁面溫度。

1.3 網(wǎng)格無關(guān)性驗(yàn)證

反應(yīng)器采用軸對(duì)稱模型,將三維模型簡(jiǎn)化為二維。為了證明網(wǎng)格的數(shù)量滿足仿真要求,分別采用4種不同數(shù)量網(wǎng)格來計(jì)算氣體通過反應(yīng)器產(chǎn)生的壓降。結(jié)果如圖2所示,75×150網(wǎng)格數(shù)量已滿足要求,因此采用此套數(shù)量網(wǎng)格進(jìn)行仿真。

圖2 網(wǎng)格數(shù)量與反應(yīng)器壓降關(guān)系Fig.2 Relationship between number of grids and reactor pressure drop

2 仿真結(jié)果與分析

2.1 操作條件對(duì)反應(yīng)的影響

首先需要了解反應(yīng)器進(jìn)口氣體溫度、速度、RP-3和O2摩爾分?jǐn)?shù)對(duì)反應(yīng)器轉(zhuǎn)化率的影響。圖3 給出了氧氣摩爾分?jǐn)?shù)為21%,氣體速度為1 m/s 時(shí),3種不同RP-3摩爾分?jǐn)?shù)下, RP-3轉(zhuǎn)化率隨反應(yīng)溫度的變化關(guān)系。從圖中可見轉(zhuǎn)化率隨著反應(yīng)溫度增加而增加,反應(yīng)溫度在460 K時(shí),轉(zhuǎn)化率增長(zhǎng)速率明顯增加,且不同RP-3摩爾分?jǐn)?shù)下的轉(zhuǎn)化率出現(xiàn)差異。因此,需要將催化床預(yù)熱到至少460 K左右,才能使反應(yīng)器有較高的效率。

固定反應(yīng)溫度為460 K,RP-3摩爾分?jǐn)?shù)為0.3%,圖4得到了不同氣體速度(vf)下RP-3轉(zhuǎn)化率隨O2摩爾分?jǐn)?shù)的變化關(guān)系。可以發(fā)現(xiàn)當(dāng)氧氣摩爾分?jǐn)?shù)低于12%時(shí),氧氣摩爾分?jǐn)?shù)增加促進(jìn)了化學(xué)反應(yīng),RP-3轉(zhuǎn)化率隨氧氣摩爾分?jǐn)?shù)增加而增加;而氧氣摩爾分?jǐn)?shù)高于12%后,氧氣相對(duì)RP-3過量,RP-3轉(zhuǎn)化率幾乎不變。鑒于耗氧型惰化系統(tǒng)啟動(dòng)時(shí),燃油箱內(nèi)氧氣摩爾分?jǐn)?shù)較高,因此可忽略氧氣摩爾濃度對(duì)反應(yīng)器起燃特性的影響。此外,當(dāng)氧氣過量時(shí),氣體速度為1.5 m/s條件下,RP-3轉(zhuǎn)化率大約為20%,氣體速度為1 m/s時(shí),轉(zhuǎn)化率大約為30%,氣體速度為0.5 m/s 時(shí),轉(zhuǎn)化率大約為50%,由于減小氣體流速后,造成反應(yīng)物在催化床內(nèi)停留時(shí)間較長(zhǎng)從而促進(jìn)了反應(yīng)進(jìn)行,因此RP-3轉(zhuǎn)化率受氣體速度影響較大,與氣體速度接近反比關(guān)系,隨著氣體速度減小而增加。

圖3 不同RP-3摩爾分?jǐn)?shù)下RP-3轉(zhuǎn)化率隨溫度的變化Fig.3 Variation of RP-3 conversion rate with reaction temperature under different RP-3 mole fractions

圖4 不同氣體速度下RP-3轉(zhuǎn)化率隨氧氣摩爾分?jǐn)?shù)的變化Fig.4 Variation of RP-3 conversion rate with oxygen mole fraction at different gas velocities

2.2 催化床內(nèi)的溫度變化

設(shè)置反應(yīng)器進(jìn)口氧氣摩爾分?jǐn)?shù)為21%,RP-3為0.3%,氣體溫度為408 K,催化床初始溫度為303 K,得到了在3個(gè)不同反應(yīng)時(shí)間下,反應(yīng)器軸向溫度與x=0.075 m處反應(yīng)器徑向溫度變化,如圖5、圖6所示。從圖5可見,在100 s時(shí)間內(nèi),受到氣體的加熱催化床前部分溫度最先上升;200 s 后,由于催化床前段溫度已接近氣體的溫度,溫度無法提升,而催化床后半段溫度快速上升;300 s后,催化床后半段由于化學(xué)反應(yīng)加劇溫度進(jìn)一步上升,最高達(dá)到550 K。

由圖6可見,由于催化床外壁有熱量損失,所以催化床溫度在靠近壁面處與軸線處相差較大。因此在起燃過程中床層內(nèi)部溫差較大。

圖5 催化床軸向溫度變化Fig.5 Axial temperature variation of catalytic bed

圖6 x=0.075 m處催化床徑向溫度變化Fig.6 Radial temperature variation of catalytic bed at x= 0.075 m

2.3 起燃條件

耗氧型惰化系統(tǒng)啟動(dòng)后需要對(duì)燃油箱抽出的氣體進(jìn)行加熱,加熱后的氣體通入反應(yīng)器的催化床中,與催化床進(jìn)行熱量交換,使催化床達(dá)到適宜的催化溫度。由于對(duì)系統(tǒng)中加熱器的重量有限制,無法提供足夠大功率的加熱器。在滿足系統(tǒng)使用要求的前提下,需要加熱后的氣體溫度盡可能低,且氣體速度盡量要小。

如前所述,催化床溫度應(yīng)在460 K以上才固定氣體速度為1 m/s,氧氣摩爾分?jǐn)?shù)為21%,RP-3摩爾分?jǐn)?shù)為0.3%,得到了不同預(yù)熱氣體溫度下催化床的平均溫度隨時(shí)間變化,如圖7所示。從圖中可見,當(dāng)氣體預(yù)熱溫度在398 K與403 K條件下,催化床的溫度從環(huán)境溫度(303 K)升高至預(yù)熱氣體的溫度后,床層溫度無法進(jìn)一步提升;而氣體預(yù)熱溫度在408 K與413 K條件下,催化床平均溫度隨時(shí)間上升后,超越了氣體預(yù)熱溫度達(dá)到了更高。當(dāng)氣體預(yù)熱溫度在398 K與403 K時(shí),氣體無法使催化床加熱至所需溫度。因此,在一定操作參數(shù)下,存在一個(gè)起燃所需最低溫度,氣體溫度必須大于此溫度時(shí),才能使反應(yīng)器成功起燃,下文將此溫度稱為起燃溫度。

如前所述,催化床溫度應(yīng)在460 K以上才有較高效率,所以本文將反應(yīng)器催化床平均溫度在5 min內(nèi)達(dá)到460 K以上作為成功起燃的標(biāo)準(zhǔn)。

圖7 氣體溫度對(duì)催化床平均溫度的影響Fig.7 Effect of gas temperature on average temperature of catalytic bed

由圖4可見,氣體速度減小促進(jìn)了化學(xué)反應(yīng)進(jìn)行。從而增加了反應(yīng)放熱,有助于縮短起燃時(shí)間。但氣體速度減少會(huì)減少催化床與氣體交換的熱量。因此需進(jìn)一步探討氣體速度大小對(duì)起燃過程的影響。固定氣體溫度為408 K,RP-3摩爾分?jǐn)?shù)為0.3%,氧氣摩爾分?jǐn)?shù)為21%。圖8給出了在不同氣體速度下催化床平均溫度隨時(shí)間的變化關(guān)系,圖9給出了不同氣體速度下反應(yīng)器出口RP-3摩爾分?jǐn)?shù)的變化。由圖8可知,在氣體速度為0.2 m/s時(shí),催化床平均溫度增加始終十分緩慢,在330 s時(shí)刻僅達(dá)到了350 K;而氣體速度在0.6 m/s與1.0 m/s條件下,催化床溫度隨時(shí)間增加呈指數(shù)型增長(zhǎng);氣體速度在1.4 m/s條件下,催化床溫度在前180 s增加速度均大于其他條件,但是隨后趨于不變。圖9可知,在氣體速度為0.2 m/s 與1.4 m/s條件下RP-3摩爾分?jǐn)?shù)在反應(yīng)器預(yù)熱過程中減小很少,而氣體速度為0.6 m/s與1.0 m/s時(shí),RP-3摩爾分?jǐn)?shù)隨時(shí)間急劇下降。

圖8 氣體速度對(duì)催化床平均溫度的影響Fig.8 Effect of gas velocity on average temperature of catalytic bed

圖9 氣體速度對(duì)出口RP-3摩爾分?jǐn)?shù)的影響Fig.9 Effect of gas velocity on RP-3 mole fraction at reactor outlet

可見催化床的熱量升高主要是由反應(yīng)熱量驅(qū)動(dòng)的。當(dāng)氣體速度過低時(shí),催化床溫度難以上升從而驅(qū)動(dòng)化學(xué)反應(yīng);當(dāng)氣體速度過高又抑制了化學(xué)反應(yīng)的發(fā)生,造成催化床溫度無法進(jìn)一步上升。因此在一定條件下,有利于起燃的氣體速度具有一定范圍,不能過高也不能過低。

圖10得到了不同RP-3摩爾分?jǐn)?shù)下催化床平均溫度隨時(shí)間的變化。從圖10中可見,RP-3摩爾分?jǐn)?shù)越大,催化床溫度上升越迅速,因此增大RP-3摩爾分?jǐn)?shù)能大幅度縮短起燃時(shí)間。由圖3可知相同溫度下,不同RP-3摩爾分?jǐn)?shù)工況下的RP-3轉(zhuǎn)化率接近,但化學(xué)反應(yīng)所釋放的熱量與RP-3摩爾分?jǐn)?shù)成正比。

固定RP-3摩爾分?jǐn)?shù)為0.2%,圖11給出了反應(yīng)器起燃溫度隨氣體速度的變化關(guān)系。起燃溫度隨氣體速度先減小后增加,最低值大約在氣體速度為0.75 m/s時(shí)。結(jié)合圖8可知,在催化床溫度較低時(shí),催化床的熱量來源主要由熱氣體提供,而氣體速度較低時(shí)無法快速加熱催化床以促進(jìn)化學(xué)反應(yīng),因此需要較高的起燃溫度;而氣體速度過大又抑制了化學(xué)反應(yīng)的發(fā)生,所以對(duì)應(yīng)起燃溫度較高。

圖10 RP-3摩爾分?jǐn)?shù)對(duì)催化床平均溫度的影響Fig.10 Effect of RP-3 mole fraction on average temperature of catalytic bed

固定氣流速度為0.4 m/s,得到起燃所需最低氣溫隨RP-3摩爾分?jǐn)?shù)的變化,如圖12所示。發(fā)現(xiàn)起燃溫度隨著RP-3摩爾摩爾分?jǐn)?shù)增加而逐漸減小,最終趨于平緩。毫無疑問,增加RP-3濃度后,反應(yīng)所釋放熱量增加,導(dǎo)致起燃溫度變低。

圖11 起燃溫度隨氣體速度變化Fig.11 Variation of light-off temperature with gas velocity

圖12 起燃溫度隨RP-3摩爾分?jǐn)?shù)變化Fig.12 Variation of light-off temperature with RP-3 mole fraction

3 結(jié) 論

利用CFD方法建立了一個(gè)耦合化學(xué)反應(yīng)的多孔介質(zhì)二維擬均相反應(yīng)器模型,通過此模型給出了耗氧型反應(yīng)器起燃過程內(nèi)部催化床溫度變化過程并研究了不同參數(shù)對(duì)耗氧型反應(yīng)器起燃反應(yīng)器過程的影響,研究結(jié)果顯示:

1) 氧氣濃度高于12%以后,對(duì)RP-3的轉(zhuǎn)化率影響不大,在實(shí)際工況下,可忽略氧氣濃度的影響。

2) 氣體溫度越高,起燃時(shí)間越短。且存在一個(gè)起燃溫度,低于此溫度反應(yīng)器無法起燃。RP-3濃度的增加會(huì)降低起燃溫度。

3) 有助于起燃的氣體速度有一定范圍。氣體速度過小,反應(yīng)熱量低難以驅(qū)動(dòng)起燃;氣體速度過大,使化學(xué)反應(yīng)進(jìn)行不夠徹底同樣對(duì)起燃不利。

此外,本文研究?jī)H單獨(dú)研究了反應(yīng)器,來探討反應(yīng)器對(duì)惰化系統(tǒng)的影響,反應(yīng)器進(jìn)口參數(shù)皆考慮定常,后續(xù)研究應(yīng)當(dāng)考慮系統(tǒng)惰化過程中氧氣摩爾分?jǐn)?shù)、RP-3摩爾分?jǐn)?shù)的變化,為燃油箱惰化系統(tǒng)的設(shè)計(jì)提供豐富的理論基礎(chǔ)。

猜你喜歡
模型
一半模型
一種去中心化的域名服務(wù)本地化模型
適用于BDS-3 PPP的隨機(jī)模型
提煉模型 突破難點(diǎn)
函數(shù)模型及應(yīng)用
p150Glued在帕金森病模型中的表達(dá)及分布
函數(shù)模型及應(yīng)用
重要模型『一線三等角』
重尾非線性自回歸模型自加權(quán)M-估計(jì)的漸近分布
3D打印中的模型分割與打包
主站蜘蛛池模板: 亚洲另类色| 欧美日韩在线第一页| 国产精品美女免费视频大全| 免费人成网站在线观看欧美| 992tv国产人成在线观看| 综合久久久久久久综合网| 精品国产免费人成在线观看| 91毛片网| 久久久91人妻无码精品蜜桃HD| 国产91线观看| 一本久道热中字伊人| 在线欧美a| 国产永久免费视频m3u8| 国产成人调教在线视频| 日本黄色不卡视频| 国产主播一区二区三区| 国产情侣一区| 国产乱肥老妇精品视频| 国产99视频在线| 国产精品成人AⅤ在线一二三四| 青青操国产| 国产不卡在线看| 国产精品熟女亚洲AV麻豆| 国产一区二区丝袜高跟鞋| 无码一区二区三区视频在线播放| 国产精品19p| 欧美色99| 夜色爽爽影院18禁妓女影院| 一本一道波多野结衣一区二区| 国产欧美视频综合二区| 亚洲免费成人网| 亚洲性视频网站| 超级碰免费视频91| 国产在线观看91精品亚瑟| 国产成人综合久久精品尤物| 91蝌蚪视频在线观看| 992tv国产人成在线观看| 波多野结衣第一页| 她的性爱视频| 丁香亚洲综合五月天婷婷| 狠狠色成人综合首页| 国产激情国语对白普通话| 亚洲开心婷婷中文字幕| av在线手机播放| 国产欧美精品一区aⅴ影院| 国产亚洲视频播放9000| 日韩a在线观看免费观看| 亚洲无码高清免费视频亚洲| 亚洲天堂视频在线播放| 亚洲精品大秀视频| 尤物精品视频一区二区三区| 激情国产精品一区| 99热这里只有精品国产99| 人妻出轨无码中文一区二区| 日韩123欧美字幕| 手机精品福利在线观看| 欧美成人日韩| 乱色熟女综合一区二区| 国产成人精品男人的天堂| 久久精品女人天堂aaa| 成人a免费α片在线视频网站| 欧美第二区| 色婷婷亚洲十月十月色天| 国产精品私拍在线爆乳| 欧美三级自拍| 国产大全韩国亚洲一区二区三区| 熟妇无码人妻| 国产精品冒白浆免费视频| 精品自窥自偷在线看| 国产成人亚洲无码淙合青草| 国产凹凸一区在线观看视频| 国产在线欧美| 欧美亚洲欧美区| 色精品视频| 日本高清免费不卡视频| 欧美a在线| 久久精品这里只有精99品| 国模沟沟一区二区三区| 一级毛片免费播放视频| Jizz国产色系免费| 国产人人乐人人爱| 成人免费视频一区二区三区|