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

        ?

        改性C/C-ZrC-SiC復(fù)合材料熱化學燒蝕數(shù)值分析 ①

        2021-07-15 01:49:20梅宗書范學領(lǐng)石成英晉小超
        固體火箭技術(shù) 2021年3期
        關(guān)鍵詞:損失率改性復(fù)合材料

        梅宗書,范學領(lǐng),石成英,晉小超

        (1.火箭軍工程大學 核工程學院,西安 710025;2.西安交通大學 航天航空學院 機械結(jié)構(gòu)強度與振動國家重點實驗室,西安 710049)

        0 引言

        碳/碳(C/C)復(fù)合材料因其突出的耐燒蝕、高比強度、高比模量、低熱膨脹系數(shù)、耐熱沖擊等一系列優(yōu)異性能,廣泛應(yīng)用于航天、航空領(lǐng)域,如再入飛行器鼻錐、固體火箭噴管和航天飛機機翼前緣[1]。但C/C材料在高溫下易氧化,從而降低了其抗燒蝕性能。基體改性通過在碳基質(zhì)中加入改性抑制劑或密封劑封閉活性點來降低氧化速率從而提高抗燒蝕性能。改性材料根據(jù)引入相的不同,可分為SiC陶瓷、UHTCs、難熔金屬等,其中UHTCs顆粒摻雜物主要是難熔金屬Zr、Hf和Ta的碳化物及硼化物,研究較多的有ZrC、ZrB2、ZrB2-SiC、HfB2-SiC、ZrB2-ZrC-SiC等,具有較強的應(yīng)用前景,可在2200 ℃以上高溫對燒蝕結(jié)構(gòu)進行有效防護[2]。

        在數(shù)值計算方面,碳基復(fù)合材料主要從表面熱化學燒蝕和內(nèi)部材料熱響應(yīng)兩個方面,對其燒蝕表面氣體產(chǎn)物濃度、質(zhì)量燒蝕率、線燒蝕率等性能參數(shù)進行計算求解。YIN Tiantian等[3]根據(jù)C/C熱化學燒蝕在低、高溫下分別由反應(yīng)動力學速率和擴散速率控制的特點,研究了C/C復(fù)合材料前緣結(jié)構(gòu)在不均勻壓力、溫度分布下的表面退移率,同時探討了基體改性對材料燒蝕性能的影響。楊德軍等[4]、WANG Chen[5]根據(jù)表面能量平衡、質(zhì)量守恒及化學平衡方程構(gòu)建表面燒蝕模型,并通過虛擬失效、重新構(gòu)建網(wǎng)格的方法,實現(xiàn)了熱化學燒蝕、燒蝕表面退縮及溫度場耦合作用下C/C復(fù)合材料燒蝕性能變化規(guī)律有限元分析。有限速率表面燒蝕模型[6-8]是C/C、碳/酚醛復(fù)合材料表面熱化學燒蝕計算的又一重要而有效的途徑,建立了一種基于表面元素質(zhì)量守恒的有限速率燒蝕表面邊界條件公式群,可預(yù)測再入空間飛行器的熱罩表面的氣動加熱、化學物質(zhì)濃度和碳質(zhì)材料燒蝕速率。目前,改性C/C復(fù)合材料的研究主要是通過試驗的方法[9-11],數(shù)值計算大多針對未摻雜改性碳/碳復(fù)合材料而言[12-14],改性C/C復(fù)合材料的數(shù)值分析則較少[11,15],本文通過熱化學平衡常數(shù)和表面元素質(zhì)量守恒得到ZrC-SiC改性C/C復(fù)合材料燒蝕表面產(chǎn)物濃度,并構(gòu)建CCZS表面氧化模型。通過燒蝕試驗結(jié)果數(shù)據(jù)和熱-力-化多場耦合數(shù)值分析計算結(jié)果比較,對構(gòu)建的改性C/C復(fù)合材料氧化模型進行驗證。

        1 改性復(fù)合材料燒蝕行為數(shù)值分析

        利用流體動力學和結(jié)構(gòu)有限元數(shù)值分析方法,結(jié)合流-固耦合多場耦合方案,對CCZS復(fù)合材料的燒蝕行為和質(zhì)量燒蝕率進行數(shù)值分析,將多場耦合數(shù)值分析方法的范圍從C/C復(fù)合材料擴展到基體改性C/C復(fù)合材料。CCZS復(fù)合材料的燒蝕是一個多組分、多相態(tài)(固態(tài)、液態(tài)、氣態(tài))的復(fù)雜熱化學反應(yīng)過程,其最終產(chǎn)物的確定目前仍是基于化學平衡的研究方法。本章將基于熱化學平衡原則和最小自由能原理,對燒蝕體系中的各組分分壓和濃度進行理論推導(dǎo)和數(shù)值分析。

        1.1 熱化學表面氧化機制

        CCZS復(fù)合材料在高溫氣流作用下表現(xiàn)出復(fù)雜的熱化學燒蝕行為,涉及多種化學反應(yīng)和氣體擴散。根據(jù)最小自由能判斷,CCZS復(fù)合材料在氧化燒蝕過程中可能發(fā)生的化學反應(yīng)如式(1)~式(10)所示,熱化學燒蝕包括碳纖維、碳基體、碳化鋯和碳化硅的氧化反應(yīng)、碳氮反應(yīng)及碳的升華等。

        3C(s)=C3(g)

        (1)

        SiO2(s)=SiO2(1)

        (2)

        C(s)+O2(g)=CO2(g)

        (3)

        2C(s)+O2(g)=2CO(g)

        (4)

        2C(s)+N2(g)=2CN(g)

        (5)

        4C(s)+N2(g)=2C2N(g)

        (6)

        ZrC(s)+2O2(g)=ZrO2(s)+CO2(g)

        (7)

        SiC(s)+2O2(g)=SiO2(s)+CO2(g)

        (8)

        本文試驗設(shè)定了三種燒蝕溫度,分別為1400、1600、1750 ℃,在這三種不同燒蝕溫度下,溫度越高,發(fā)生的熱化學反應(yīng)越多,本節(jié)以燒蝕溫度為1750 ℃的工況進行組分濃度表達式的推導(dǎo),在此溫度下,式(8)中產(chǎn)生的固態(tài)SiO2會發(fā)生相變,呈液態(tài),會對碳材料的氧化反應(yīng)起到阻隔氧化性氣體的作用。但因其含量較少,在進行氧化模型構(gòu)建時,忽略SiO2熔化對氧化反應(yīng)的影響,即SiO2熔化后對氧化性氣體的阻隔作用忽略不計。CCZS復(fù)合材料在燒蝕過程中的表面熱化學反應(yīng),除了式(1)~式(8)中所列出的材料固體相與氣體的反應(yīng),還應(yīng)包括氣體自身的離解。

        O2(g)=2O(g)

        (9)

        N2(g)=2N(g)

        (10)

        1.2 熱化學反應(yīng)模型

        當反應(yīng)達到平衡時,根據(jù)化學平衡計算原理,邊界反應(yīng)層中產(chǎn)物分壓和反應(yīng)物分壓之比為一個以溫度為變量的平衡常數(shù),其表達式為

        (11)

        其中,pp,i、pr,i分別為式i(i=1~10)中產(chǎn)物分壓和反應(yīng)物分壓;α,β分別為式i中產(chǎn)物和反應(yīng)物的系數(shù)。式(1)~式(10)所對應(yīng)的平衡常數(shù)系數(shù)a、b經(jīng)過查表[16]后擬合可得表1。

        表1 平衡常數(shù)系數(shù)

        組元濃度Ci和組元分壓pi之間的關(guān)系式為

        (12)

        (13)

        根據(jù)元素質(zhì)量守恒和相容性條件,可得元素濃度和組元濃度之間的關(guān)系為

        (15)

        (16)

        (17)

        (18)

        式中f1為復(fù)合材料中ZrC質(zhì)量分數(shù);f2為復(fù)合材料中SiC質(zhì)量分數(shù);CO,e、CN,e分別為邊界層外緣氧元素和氮元素質(zhì)量分數(shù);B為無因次質(zhì)量損失率。

        聯(lián)立式(14)~式(16),可得

        (19)

        從而得到無因次質(zhì)量損失率為

        (20)

        1.3 產(chǎn)物濃度計算

        根據(jù)式(1)~式(10)得到的反應(yīng)平衡常數(shù)關(guān)系式和式(17)、式(18)、式(20)可得以下各組元濃度表達式:

        (21)

        (22)

        (23)

        (24)

        (25)

        (26)

        (27)

        (28)

        (29)

        (30)

        上列方程組需迭代求解,先給定B和T,由式(21)~式(30)可求出全部壁面組分質(zhì)量分數(shù),將算出的濃度代入式(20),可求出一個新的B值。若給定的B與求出的B之差滿足精度要求,則給定的B為所求的解;若不滿足,需要重新給出B,反復(fù)迭代,表面溫度T由表面能量平衡方程確定。經(jīng)過數(shù)值計算,得到CCZS復(fù)合材料表面燒蝕層的各組元濃度。

        2 CCZS試樣燒蝕試驗

        2.1 試樣及燒蝕平臺

        采用聚合物浸漬裂解(Precursor Infiltration and Pyrolysis,PIP)結(jié)合化學氣相沉積(Chemical Vapor Deposition,CVD)制備CCZS材料。為研究材料在三種試樣表面溫度下的抗氧化燒蝕性能,制備了3顆形狀及材料相同的錐體燒蝕試樣。編號1~3的試樣分別對應(yīng)表面溫度為1400、1600、1750 ℃的三種工況。由于制備工藝的限制,三顆燒蝕試樣從尺寸上有微小的差別,但這種外形尺寸上的微小差別對地面燒蝕流場和材料燒蝕行為的影響程度是有限的。因此,在做燒蝕試驗及試驗后的分析時,可認為三顆試樣是完全相同的。

        燒蝕體幾何尺寸和基本物理屬性見表2,試樣的密度約為(2.795±0.014) g/cm3,摻雜濃度(15±0.29)%,高度(11.2±0.04) mm,底部直徑(10±0.03) mm,錐角為42°±0.5°。由于三顆試樣的物理屬性差別是很小的,在后面的分析中,不考慮三顆試樣本身的差異對燒蝕性能的影響,而主要討論不同工況下試樣的燒蝕行為及性能的差異。在下文的數(shù)值模擬建模中,為了簡化建模過程,燒蝕體的物理參數(shù)設(shè)為密度2.795 g/cm3,高度11.2 mm,底部直徑10 mm,錐角42°。

        采用氧-丙烷熱沖擊試驗平臺進行地面燒蝕試驗,通過調(diào)節(jié)O2和C3H8的體積流量進行燒蝕溫度調(diào)節(jié)。圖1為燒蝕試驗設(shè)備及不同工況下的試樣表面溫度和氣體流量隨時間的變化曲線。圖1中,溫度在設(shè)定的目標溫度上下有微小的波動,溫度的穩(wěn)定區(qū)域可滿足燒蝕時長的需要。因此,工況是滿足試驗要求的。

        (a)Oxygen-propane thermal shock test platform (b)1400 ℃

        采用氧-丙烷火焰對1~3號試樣分別在1400、1600、1750 ℃的高溫下燒蝕90 s,火焰噴口離試樣駐點200 mm,火焰方向與試樣軸線平行。試樣燒蝕過程中,可觀察到在1750 ℃時試樣表面有亮斑,推斷為SiO2,推斷基于以下兩點:一是排除為ZrO2,因其熔點為2677 ℃,且具有低的蒸汽壓,SiO2的熔點為1670 ℃,且具有較高的蒸汽壓,在1750 ℃的氣流下,只可生成熔融的SiO2;二是亮斑的位置大多在錐體中間和尾端部位,應(yīng)是由于氣流的沖刷作用導(dǎo)致的,符合液態(tài)SiO2流動性好的特點。

        2.2 燒蝕形貌分析

        三種工況下,三顆試樣燒蝕后的宏觀形貌如圖2所示。根據(jù)圖2(a)、(b)可知,在較低溫度(≤1600 ℃)下,試樣的宏觀形貌變化并不明顯,但試樣表面形成蓬松易脫落的物質(zhì)。

        (a)1400 ℃ (b)1600 ℃ (c)1750 ℃

        通過掃描電子顯微鏡(SEM)和能譜儀(EDS),對試樣表面的微觀形貌進行觀察。研究發(fā)現(xiàn),表面蓬松物質(zhì)主要為ZrO2和SiO2呈顆粒形式的堆積物。根據(jù)圖2(c)可知,在高溫下試樣表面形成熔融SiO2,并在局部位置聚集成微小顆粒附著在試樣表面,燒蝕結(jié)束冷卻后在試樣表面呈現(xiàn)出白色顆粒。觀察圖2(c)試樣的微觀形貌,可觀察到明顯的白色線狀物質(zhì),即試樣表面形成的熔融的SiO2,由于SiC含量較少,所以沒有形成成片的SiO2保護膜,但仍對保護內(nèi)部碳質(zhì)材料和氧化形成的ZrO2起到了重要的作用。

        在1400、1600、1750 ℃燒蝕90 s工況下,CCZS試樣質(zhì)量損失速率較低,分別為0.838、0.875、0.881 mg/s,表明CCZS復(fù)合材料具有良好的抗燒蝕性能。根據(jù)化學反應(yīng)速率可知,熱化學燒蝕一般隨溫度的增加呈指數(shù)型增加,CCZS試樣的質(zhì)量損失速率雖然隨著溫度的增加提高,但相比于由1400 ℃到1600 ℃,由1600 ℃到1750 ℃時燒蝕速率增加并不明顯。根據(jù)上述微觀形貌的分析結(jié)果可推測,試樣表面的熔融SiO2能夠起到隔絕氧氣的作用,對內(nèi)部材料具有較好的保護作用,并使疏松的ZrO2不易被氣流沖走,且隨著溫度的升高,生成的SiO2量增多,從而降低了隨著溫度升高而引起的質(zhì)量損失。

        3 多物理場耦合分析

        3.1 燒蝕流場數(shù)值分析

        利用Fluent中的非預(yù)混燃燒模型表征氧氣與丙烷的燃燒反應(yīng),其中丙烷作為燃料流,氧氣作為氧化劑流。非預(yù)混燃燒計算中不使用有限速率化學反應(yīng)模型,而是用統(tǒng)一的混合物濃度作為未知變量進行求解。因此,無需計算代表組元生成或消失的源項。

        對燃氣發(fā)生器噴嘴和氧-丙烷焰燃燒狀態(tài)進行仿真,建立二維氧-丙烷焰熱流場物理模型,根據(jù)物理模型對稱,建模時只考慮其1/2模型,如圖3所示。模型尺寸為50 mm×520 mm,針對發(fā)生器的3層小孔,構(gòu)建3個1 mm混合氣體出口,每個出個間隔2 mm。沿著發(fā)生器噴口方向提供速度為1 mm/s的流動空氣進行替代。構(gòu)建燒蝕試樣壁面,在壁面、燃氣出口等區(qū)域進行網(wǎng)格加密處理。

        圖3 氧-丙烷熱流場計算域網(wǎng)格劃分

        對模型邊界條件進行設(shè)置,三個燃氣入口均采用質(zhì)量入口(Mass flux inlet),入口參數(shù)需要給定氣體在出口處的質(zhì)量流率和初始靜壓。對于燃燒過程在開放環(huán)境下的自由射流場,出口邊界選擇壓力出口邊界條件(Pressure outlet)。對于噴嘴部分,除了出口通氣孔外,其他部位均為壁面邊界條件。上邊界同樣選為出口壓力邊界,下邊界選為軸對稱(Axisymmetric)。

        對三種工況下的熱流場進行穩(wěn)態(tài)數(shù)值計算,流場溫度分布云圖如圖4(a)所示。由圖4(a)可知,隨著燃氣流量的增加,流場最高溫度并沒有明顯的升高,推測熱流場溫度應(yīng)與燃氣中氧氣和丙烷的混合比例有關(guān)。三種工況1400、1600、1750 ℃對應(yīng)的試樣駐點處的溫度分別為1471.78、1677.57、1931.13 ℃,流場對稱軸上的溫度變化曲線如圖4(b)所示。燒蝕試樣駐點處的溫度基本上與測量得到的主控溫度一致,說明流場模擬的準確性,可作為多場耦合數(shù)值分析的流場部分使用。

        (a)Cloud map of flow field temperature distribution at T=1750 ℃ (b)Variation curves of temperature

        3.2 燒蝕試樣建模與設(shè)置

        利用Abaqus對燒蝕試樣進行有限元建模分析,結(jié)構(gòu)的尺寸和坐標與Fluent中的一致。試樣各項尺寸參數(shù)根據(jù)2.1節(jié)中的設(shè)置,分別為密度2.795 g/cm3,高度11.2 mm,底部直徑10 mm,錐角42°,結(jié)構(gòu)建模如圖5所示。網(wǎng)格采用結(jié)構(gòu)化的四邊形網(wǎng)格和軸對稱傳熱-變形耦合單元CAX4T,單元數(shù)為4875,節(jié)點數(shù)為5061,在燒蝕最為劇烈的端頭部分,進行了網(wǎng)格加密處理。

        圖5 結(jié)構(gòu)有限元模型網(wǎng)格劃分

        (31)

        式中 下標w為固體表面;下標e為來流邊界。

        表3 材料熱力學性能

        (32)

        式中fr為機械剝蝕因子,一般取1或1.5,本文取為1.5。

        總的材料質(zhì)量損失率和線燒蝕率為

        (33)

        (34)

        將燒蝕速率嵌入到表面燒蝕退移子程序Umeshmotion中,對CCZS復(fù)合材料燒蝕試樣熱-力-化多場耦合模型進行計算。

        3.3 數(shù)值計算結(jié)果分析

        圖6為三種燒蝕溫度下燒蝕試樣內(nèi)部溫度分布,根據(jù)表面溫度可知,模型可較準確地模擬燒蝕試驗的溫度分布。其中,最大溫度處駐點的溫度分別為1384、1566、1734 ℃。得到的表面溫度均比試驗測得的主控溫度1400,1600,1750 ℃要低,推理可能是如下兩種因素導(dǎo)致:

        (1)主控溫度的測量采用的是紅外測溫儀,火焰可能會對紅外測溫儀造成影響,從而產(chǎn)生誤差;

        (2)流場能量項沒有考慮表面反應(yīng)熱,從而對流場熱量計算時偏低。

        (a)T=1400 ℃ (b)T=1600 ℃

        由于熱膨脹的作用,且熱膨脹導(dǎo)致的計算節(jié)點正位移要比氧化燒蝕導(dǎo)致的負位移要大,即膨脹速率大于燒蝕退移速率,計算得到的試樣尺寸隨著燒蝕時間逐漸變大,尤其是軸向的長度變化較大。隨著燒蝕溫度的升高,熱膨脹占主導(dǎo)作用的節(jié)點正位移逐漸提高,但單位溫度增長下提高的溫度幅度明顯降低。原因是隨著溫度的升高,熱膨脹率有減低趨勢;另外,溫度越高,氧化導(dǎo)致的負位移越高,抵消了部分熱膨脹位移。

        圖7為三種工況下,試樣駐點溫度隨時間的變化,線型數(shù)據(jù)表示通過數(shù)值計算得到的數(shù)據(jù),符號型數(shù)據(jù)代表試驗測得的試驗數(shù)據(jù)。

        圖7 駐點溫度變化曲線

        由圖7可知,在燒蝕的前10 s內(nèi),駐點溫度上升最快,之后隨著燒蝕時間增加,溫度上升速率逐漸降低,在燒蝕約40 s后,溫度已經(jīng)基本達到了主控溫度,之后略微有所上升。試驗數(shù)據(jù)和計算數(shù)據(jù)在燒蝕時長為10 s后吻合度較高,10 s之前由于溫度處于急速上升階段,存在一定的誤差。整體而言,試驗數(shù)據(jù)與計算數(shù)據(jù)變化趨勢一致,表明所構(gòu)建模型的有效性。

        圖8為燒蝕90 s時長中,試樣的軸向表面位置出的質(zhì)量損失率。由圖8可知,在前30 s內(nèi),質(zhì)量損失率增長較快,之后增長趨于平緩,燒蝕60 s后,質(zhì)量損失率逐漸趨于穩(wěn)定。穩(wěn)定后的三種工況下,表面溫度由低到高,試樣駐點的質(zhì)量損失率分別為2.41×10-3、2.83×10-3、3.23×10-3kg/(m2·s)。比較駐點處和尾部的質(zhì)量損失率可知,兩者之間的質(zhì)量損失率大小相差較小,氧化燒蝕在試樣空間位置的差異不是很明顯,原因應(yīng)是試樣尺寸較小,流場參數(shù)在試樣表面的幾乎相同。

        (a)T=1400 ℃ (b)T=1600 ℃

        將得到的單位面積的質(zhì)量損失率通過乘以試樣受氧化面積,得到單位時間試樣整體的質(zhì)量損失率,其單位變?yōu)閙g/s。圖9曲線數(shù)據(jù)為三種工況下,試樣整體質(zhì)量燒蝕率隨時間的變化。

        圖9中的矩形數(shù)據(jù)為燒蝕試驗得到的不同工況下燒蝕時長90 s下的平均質(zhì)量損失率??芍?,隨著溫度的升高,燒蝕質(zhì)量損失率增加趨勢變緩,這與上文中通過分析摻雜超高溫陶瓷顆粒燒蝕機理時得到的結(jié)論是一致的。同時,與試驗測量得到的平均燒蝕率相比,也說明計算數(shù)據(jù)與試驗數(shù)據(jù)得到的質(zhì)量損失速率較為一致,進一步證明了模型的有效性和可靠性。

        圖9 燒蝕質(zhì)量損失率

        4 結(jié)論

        利用熱化學平衡原則進行改性復(fù)合材料氧化燒蝕行為的表征,推導(dǎo)并驗證了基體改性C/C復(fù)合材料的氧化模型。仿真與試驗驗證了模型的有效性,為發(fā)展改性材料的氧化模型提供了理論參考。文中的理論方法是一種簡化近似計算方法,沒有對液相產(chǎn)物在阻隔氧氣擴散和降低覆蓋點位處材料剝蝕率等方面的影響,因此接下來需要對液相產(chǎn)物的物理尺寸及附著在材料表面的分布位置和稀疏度等分布屬性進行研究,從而獲得液相產(chǎn)物更為精確的防熱機制和表征方法。

        猜你喜歡
        損失率改性復(fù)合材料
        農(nóng)業(yè)農(nóng)村部印發(fā)《意見》提出到2025年農(nóng)產(chǎn)品加工環(huán)節(jié)損失率降到5%以下
        帶有治療函數(shù)及免疫損失率的SIRS流行病模型的動力學分析
        民機復(fù)合材料的適航鑒定
        復(fù)合材料無損檢測探討
        電子測試(2017年11期)2017-12-15 08:57:13
        12部使用一年后最廉價轉(zhuǎn)售車
        海外星云(2016年19期)2016-10-24 11:53:42
        2014~2015年冬季美國蜂群損失調(diào)查
        P(3,4HB)/PHBV共混改性及微生物降解研究
        中國塑料(2016年12期)2016-06-15 20:30:07
        ABS/改性高嶺土復(fù)合材料的制備與表征
        中國塑料(2015年11期)2015-10-14 01:14:14
        聚甲醛增強改性研究進展
        中國塑料(2015年9期)2015-10-14 01:12:17
        聚乳酸擴鏈改性及其擠出發(fā)泡的研究
        中國塑料(2015年4期)2015-10-14 01:09:19
        亚洲AV无码一区二区一二区色戒| 99精品久久99久久久久| 国产v片在线播放免费无码| 久久精品国产亚洲av大全| 日韩av无码成人无码免费| 亚洲视频天堂| 免费黄色福利| 国产aⅴ丝袜旗袍无码麻豆| 国产精品成人有码在线观看| 日韩在线精品免费观看| 日产一区日产2区日产| 蜜桃视频免费进入观看| 夹得好湿真拔不出来了动态图 | 国产黄色片在线观看| 国产一区二区三区最新视频| 亚洲一区二区三区视频免费看 | 蜜臀aⅴ国产精品久久久国产老师| 国产一区亚洲欧美成人| 性感人妻一区二区三区| 91久久大香伊蕉在人线国产| 一边摸一边抽搐一进一出口述| 未满十八勿入av网免费| 99精品热这里只有精品| 久草视频国产| 99re国产电影精品| 二区三区视频在线观看| 在线亚洲日本一区二区| 婷婷色国产精品视频二区| 国产成人精品999视频| a级毛片免费观看网站| 亚洲国产另类久久久精品小说 | 久久精品蜜桃亚洲av高清| 精品亚洲成a人在线观看| 天躁夜夜躁狼狠躁| 久久99久久99精品免视看国产成人| 亚洲国内精品一区二区在线| 亚洲一区二区三区激情在线观看 | 亚洲一区二区三区av资源 | 美女一区二区三区在线视频| 开心五月婷婷激情综合网| 人妻中文字幕无码系列|