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

橫向過載對固體火箭發動機推進劑點火建壓過程的影響

2021-11-01 09:08:10官典李世鵬劉筑王寧飛
兵工學報 2021年9期
關鍵詞:效應發動機模型

官典, 李世鵬, 劉筑, 王寧飛

(1.北京理工大學 宇航學院, 北京 100081; 2.中國運載火箭技術研究院 空間物理重點實驗室, 北京 100076)

0 引言

固體火箭發動機是一種結構緊湊、性能優良的動力裝置,其工作過程耦合了燃燒、傳熱、輻射、火焰傳播、壓力建立和凝聚相侵蝕等一系列復雜物理和化學現象[1]。有研究表明:在過載情況下(如旋轉產生的徑向過載或橫向機動產生的橫向過載),發動機內的上述物理或化學過程均會發生較大改變,致使真實飛行中發動機性能與地面靜止實驗結果相差很大,甚至會因為地面設計的參考不當造成發動機殉爆事故[1-6]。對大部分飛行試驗結果研究發現,發動機從初始點火到穩定燃燒過程極易出現故障。由于實驗難度大和數值描述復雜等原因,目前學術界對過載下在推進劑點火建壓過程中侵蝕燃燒與過載效應耦合規律的研究依然匱乏。為了預示這一規律并保證發動機大過載下正常開啟和工作,對橫向過載條件下推進劑點火建壓期間燃燒性能變化規律進行研究是基礎的且勢在必行的。

針對無過載時固體火箭發動機點火演變規律,國內外學者通過實驗與數值模擬等方法已經做了大量研究。各種測試手段(如壓力傳感器及其熱電偶測試[7]、鐳射多普勒速度場測試[8]以及高速攝影[9])能夠準確地捕捉實驗中的關鍵數據,反映部分實驗現象。但由于固體火箭點火是強瞬態過程,工程上很難測量和記錄豐富的瞬時數據。另一方面,數值仿真在研究多維下點火瞬態過程中已經逐漸得到成熟運用,如基于純流體假設的數值模型[10-11]、基于兩相流假設的數值模型[12-13]和流體- 固體耦合數值模型[14-15],在數值結果與實驗結果一致時,數值方法在研究發動機瞬態過程中具有特定優勢。

截至目前,已經有很多學者開展了過載條件下推進劑燃燒性能的研究。Northam[16]發現作用在燃面上的加速度矢量角度和推進劑成分對燃速變化影響較大。Ye等[17]指出燃面附近雷諾數可用于描述微觀燃燒模型中旋轉和壓力對氣相火焰偏轉角的影響。Caveny等[18]、Willoughby等[19]和Crowe[20]對過載情況下的推進劑燃燒規律都提出了半經驗模型。Greatrix[21]、Greatrix等[22]、Greatrix[23]在前人基礎上,引入燃面能量薄層和加速度質量通量系數,提出描述適用面更廣的過載條件下推進劑燃燒特性的模型。關于推進劑的侵蝕燃燒,已經有很多學者建立了多種推進劑燃燒侵蝕模型。其中包括:Lenoir等[24]基于中心高溫氣流對流傳熱理論提出了L-R模型,該模型廣泛應用于固體火箭發動機的設計中;Greatrix等[25]提出了優化型對流傳熱反饋模型,該模型對流動校正因子的依賴性較小;King[26]提出的火焰彎曲模型,從微觀擴散火焰傾斜加大燃面熱對流的方式解釋了侵蝕效應。在此基礎上,Guan等[27]采用數值模擬方法對高速自旋條件下小型發動機點火機理進行研究,得到了旋轉過載下點火規律與無過載情況不一致的理論描述。

本文針對橫向過載下推進劑點火建壓特點,以靜態點火建壓模型為基礎,將用戶自定義函數(UDF)首次引入橫向過載場模型、過載下推進劑燃燒模型以及侵蝕燃燒模型,對不同橫向過載下發動機內推進劑點火建壓規律進行研究。

1 物理模型和計算方法

1.1 發動機結構

為了檢驗侵蝕、傳熱和加質等模型在推進劑初始建壓模擬中的正確性,參考Peretz靜態點火實驗[28],建立與實驗一致的物理模型,如圖1所示。

圖1 內嵌雙面推進劑的發動機結構示意圖Fig.1 Structure diagram of rocket motor with two slices of propellant

在同一橫向過載下,對稱燃面因過載相對矢量角的不同而出現燃燒差異,為了達到控制單一變量的目的,在模型驗證完后,將Peretz實驗裝置中雙面推進劑裝藥形式改為單面推進劑裝藥形式,原有的另一面推進劑用等厚度鋼材替代,如圖2所示。下文所有物性和尺寸參數均與Peretz實驗[28]統一。

圖2 內嵌單面推進劑的發動機結構示意圖Fig.2 Structure diagram of rocket motor with single slice of propellant

1.2 數值計算方法

發動機內流場使用基于k-ε雙方程湍流模型的Navier-Stokes方程求解[29]。基于Fluent軟件,使用UDF實現點火過程模擬。將過載場效應、侵蝕效應和過載燃速效應嵌入模型源項,表征橫向過載下發動機內物理現象。發動機內流場區域控制方程包括如下質量、動量、能量守恒方程和狀態方程:

(1)

僅考慮燃氣與推進劑表面之間對流傳熱。推進劑表面溫度Ts作為點火判定依據,由(2)式和(3)式耦合計算獲得。

(2)

式中:等號左項表示推進劑表面導熱狀態,右項表示對流換熱量;λp為推進劑導熱系數;Tp為推進劑溫度分布函數;n為代表燃面法向。推進劑內二維瞬態傳熱方程如下:

(3)

式中:Cp為推進劑比熱。

在本文中推進劑燃速rg主要受過載效應與侵蝕效應影響[21],如(4)式所示:

rg=r0+Δre+Δrb,

(4)

式中:r0為推進劑基礎燃速,r0=apn,a為燃速系數,n為燃速壓力指數;Δre為因侵蝕效應帶來的燃速增量,遵循L-R侵蝕燃速公式,參數依據文獻[28]給出(見表1):

表1 點火過程相關物性參數

Δre=ξ1hcexp(-βr0ρp/vrρr),

(5)

(6)

ξ1和β為比例系數,hc為xr處推進劑表面對流換熱系數,xr為距離裝藥前端的距離,vr為xr處氣體流速,ρr為xr處氣體密度,ξ2為比例系數,Pr為普朗特數,W為燃氣摩爾質量,Tsg=(Ts+Tf)/2,dh為通道濕潤周直徑,Ap為xr處通道面積,At為喉部面積;Δrb為橫向過載效應產生的燃速變化,由Greatrix等[22]和Greatrix[23]推進劑過載燃燒模型計算獲得:

(7)

Ti為推進劑初溫,λ為氣相熱導率,Ga為加速度的質量通量,如(8)式所示:

(8)

Ga0表示最大加速度質量通量,?d為增強定向角(即側向/縱向位移角),K為定向角修正因子,?為加速度偏角,δ0為燃燒區厚度,如(9)式所示:

(9)

ΔHp為凈放熱量。

由于本文采用內嵌單面推進劑的發動機結構,控制了橫向過載方向與推進劑表面矢量關系。圖3列舉了3種過載下過載場效應與過載質量通量Ga關系。圖3(a)所示發動機不受加速度作用;圖3(b)為推進劑所受加速度ay指出燃面,流場慣性加速a′y指入燃面(加速度a′y為正),此時Ga=|Ga0|;圖3(c)為當推進劑所受加速度ay指入燃面,流場慣性加速a′y指出燃面(加速度a′y為負),此時Ga=0[30-31].后續提到的加速度均指慣性加速度a′y.

圖3 過載場效應以及過載條件下的質量通量GaFig.3 Accelerative field effect and mass flux Ga in acceleration condition

依據(7)式、(8)式和(9)式,圖4給出了推進劑過載燃速(rb=Δrb+r0)隨過載和壓力的非線性變化規律。由圖4可見:當相對加速度指向燃面時a′y為正,燃速隨過載和壓力的增加而非線性增大;在加速度背離向燃面時a′y為負,燃速隨過載幅值增加而略微減小,隨著壓力增加而非線性增加;相比負向加速度,正向加速度對燃速幅值變化的影響更大。

圖4 推進劑燃速、表面過載與表面壓力的非線性關系Fig.4 Non-linear relationship between propellant burning rate and surface acceleration and surface pressure

固體火箭發動機在點火時,腔體內氣體高速流動,呈高雷諾數特征,并伴隨有壓縮波甚至沖擊波。由于k-ε湍流模型適合完全湍流流動,采用k-ε對上述N-S模型進行封閉是適用的。湍流動能k及湍流能量耗散率ε分別由以下兩式求得:

(10)

(11)

式中:ui為時均系數;μm為層流黏性系數;μt為湍流黏性系數,

(12)

表2 湍流模型中的常值參數[32]

本文點火模型的簡化借鑒文獻[28],并加入過載相關假設如下:

1)基于發動機結構,將物理模型簡化為二維平面模型。

2)使用固相點火準則作為點火判定條件,即壁面溫度Ts大于點火溫度視為燃面點燃;

3)不考慮化學反應,推進劑燃燒產物與點火器燃氣具有相同的Cf、W和λ,均視為可壓縮理想氣體。

4)考慮燃氣沖刷帶來的侵蝕效應。

5)考慮過載情況下,流場受慣性力作用。

6)過載燃速效應基于Greatrix模型。

1.3 計算模型和邊界條件

根據模型結構,本文建立Oxy平面二維數值模型(見圖5)。圖5中,沿著推進劑表面距離推進劑頭部0.05L、0.50L和0.95L位置設置A、B和C點,L為推進劑長度,Coupled壁面固相網格第1層高度1×10-6m,Coupled壁面氣相第1層網格高度5×10-5m. 由于氣流流速變化較大,燃燒室兩側加密。點火入口采用質量流量入口,質量通量為93.18 kg/(m2·s);出口采用壓力出口;流體- 固體耦合交界面采用Coupled壁面,其他壁面采用絕熱壁面,由于y+始終處于20~60之間,壁面函數采用標準壁面函數。初始溫度298 K;初始壓強101 325 Pa. 采用PISO算法對壓力和速度進行解耦。密度、動量和能量等方程采用2階迎風格式進行離散,瞬時項離散格式采用2階隱格式。為了滿足Courant-Friedrich-Levy(CFL)穩定條件,時間步長選取1×10-6s.

為了方便燃速研究,做如下定義:

Ab(i)=Δrb(i)/r0(i),

(13)

式中:Ab(i)為推進劑表面A點在i時刻由過載效應帶來的相對燃速增量,用以表示過載效應對推進劑燃燒的影響程度(Bb和Cb以此類推);Ae(i)=Δre(i)/r0(i),

(14)

式中:Ae(i)為推進劑表面A點在i時刻由侵蝕效應帶來的相對燃速增量,表示侵蝕效應對推進劑的影響程度(Be和Ce以此類推)。

βj(i)=Δrj(i)/(Δre(i)+Δrb(i)),

(15)

式中:βj(i)為燃速增量占比,表示不同i時刻過載或侵蝕效應在燃燒變化過程所占比重;j=e或b,代表侵蝕作用或過載作用。

1.4 數值驗證

數值仿真結果與Peretz的無過載固體火箭點火實驗數據對比如圖6所示。由圖6(a)可見,數值計算得到的瞬變壓力在點火滯后期和火焰傳播期與靜態數據十分貼近,峰值壓力誤差小于5%,表明該模型對瞬變壓力的預測具有較高精度。圖6(b)為不同時刻沿軸向壓力分布,可以發現該模型對發動軸線壓力分布預測準確性均較高。上述結果驗證了本文所建傳熱、侵蝕和加質模型的合理性。

圖6 Peretz等實驗結果[28]和本文模型計算結果對比Fig.6 Comparison of Peretz’s measured results[28] and numerical results calculated by the proposed model

為了排除網格對過載計算結果的影響,針對文中極端過載100g工況,對不同橫/縱節點分布(40×480,60×730和80×950)的3套結構化網格進行網格無關性分析。其中,粗網格(40×480)的Coupled壁面對應的氣相第1層網格高度為1×10-4m,中等數量網格(60×730)和精細網格(80×950)的Coupled壁面對應的氣相第1層網格高度為5×10-5m. 燃燒室內瞬變升壓結果如圖7所示。由圖7可以看出:3套網格對點火過程的預示差異不大,相較于精細網格的預測結果,粗網格預測壓力值在整個建壓過程均較小;中等密度網格與精細網格相比,趨勢十分接近,峰值壓力誤差小于2%. 因此在確保計算精度的前提條件下,本文擬在下文中所有算例中采用中等密度結構網格。從無過載點條件點火數值模型結果與實驗結果的比較,到添加過載效應模型后點火數值模型的網格無關性驗證,目的是保證本文模型預示結果的可信度。

圖7 橫向過載條件下網格無關性驗證Fig.7 Grid independence in accelerative field effect

2 結果分析

2.1 橫向過載大小對于點火升壓過程影響

圖8所示為不同橫向過載下點火壓力曲線圖。圖8中,Δp為壓力增量,p0g表示0g時的壓力。以0g過載壓力峰值pmax對各工況下壓力進行無量綱處理,從圖8中可見:橫向過載的出現對點火建壓及穩定燃燒整個過程均產生了明顯影響,過載導致的壓差隨時間增長而增大;正向過載作用下的內彈道壓力明顯高于無過載情況,在100g下,壓力增量最高達到117%;負向過載時的燃燒室壓力在點火過程中均小于無過載時,-100g時,壓力減小量在燃燒穩定時候達到最大,為無過載時的-10%.可見100g對壓力值的改變量幾乎是-100g時的12倍。由此可以推斷:在同一方向過載下,具有對稱燃面的固體發動機從點火初始就出現較大的燃燒差異,并隨著時間變化表現為累積效果,而這些變化在無橫向過載下是不明顯的。

圖8 不同橫向過載下點火壓力曲線圖Fig.8 Impact of different radial accelerations on ignition process

與50g相比,100g對燃燒室建壓過程影響相對較大。0g到50g的壓力增長幅度Δp/p0g為74%,大于從50g到100g的壓力增長幅度43%,可以看出:壓力增長與加速度現正相關,但壓力增速隨加速度的增加呈減緩趨勢。

圖9所示為不同橫向過載下點火過程中時刻點分布圖。對圖9中點火過程3個階段(點火滯后期、火焰傳播期和火焰填充階段)進行比較:-100~100g5種工況下點火滯后期時間分別為3.51 ms、3.24 ms、3.2 ms、3.05 ms和2.82 ms,可以發現正向過載下點火滯后期縮短。這主要是因為負向過載使得點火燃氣遠離推進劑表面,削弱對流換熱作用,導致推進劑表面在未點燃前熱通量下降,到達點火溫度時間延長。由于采用純氣體假設且推進劑首次點燃位置靠近點火入口,慣性作用對點火燃氣的偏向影響在推進劑頭部位置表現較弱,所以點火滯后變化量較小,隨著發動機點火燃氣顆粒、尺寸和通道寬度的增加,上述變化將會更加明顯;-100g~100g5種工況的火焰傳播時間分別為45.7 ms、45.5 ms、44.5 ms、35.2 ms和32.6 ms,表明正向過載加速了推進劑燃面火焰傳播速度,負向過載降低了推進劑燃面火焰傳播速度。從影響程度上看,正向過載對火焰傳播時間上的影響大于負向過載;-100~100g5種工況火焰填充時間分別為12.3 ms、11.1 ms、12.0 ms、14.6 ms和16.4 ms,可以發現正向過載引起更高的壓力峰值,延長加壓填充時間。由于正向過載對點火初期中火焰傳播時間的削減大于峰值壓力增加引起的填充時間的增加,整個點火過程表現為點火延遲縮短。同一個過載下,對稱面燃面在點火初始可能出現燃燒面積不對稱的現象,進而存在偏離設計彈道的可能。

圖9 不同橫向過載下點火過程中時刻點分布圖Fig.9 3 stages of inigition process under different radial accelerations

2.2 橫向過載對于推進劑燃燒過程影響

圖10所示為-100g、0g、50g和100g4種過載下A點、B點和C點由過載/侵蝕兩種效應分別引發的相對燃速增量,計算方式依據(13)式和(14)式所示。圖10可見:在推進劑前/中部,0g和-100g的Ae和Be相差較小,而推進劑后部Ce在-100g過載下整體略小于0g過載,這是因為負向過載致使燃速下降(Ab、Bb和Cb出現負值,其中,A點處出現最大減少量-13.5%),削弱上游燃氣加注強度,進而減弱了后端的侵蝕過程,這解釋了負向過載削弱侵蝕效應的現象;對比100g和0g的Ae,二者差異較小,表明正向過載的出現對推進劑頭部侵蝕燃燒的影響較小。對比100g和0g的Ce,二者差異較大,表明正向過載的出現對推進劑后段侵蝕效應的影響較為明顯;對比100g和50g,侵蝕燃速Be和Ce在推進劑中后段的影響明顯高于過載效應Bb和Cb,但在較高過載時(100g),推進劑前段過載效應Ab強于侵蝕效應Bb并起主導作用。橫向過載條件下,推進劑燃燒將發生變化,被改變程度由燃面所處位置決定,表現為:推進劑前段主要由于過載效應影響,后段主要由侵蝕效應影響;橫向過載效應對侵蝕效應的作用是通過影響推進劑頭部燃燒而加劇/削弱下游流速,進而加劇/削弱推進劑侵蝕效應實現的。

圖10 不同橫向過載下推進劑不同位置的相對 燃速增量分布散點圖Fig.10 Scatter plot of the non-dimensional burning-rate augmentation at different positions of propellant under different lateral accelerations

表3所示為不同橫向過載下侵蝕效應/過載效應去量綱燃速增量時變圖。表3從時間尺度表明:無論正向過載還是負向過載,燃速變化均出現明顯瞬變區和穩態區。將燃速變化率第1次出現小于5%時刻作為瞬態- 穩態轉折點,可以發現:過載燃速變化的轉折點連線近乎豎直(見表3內各過載效應圖中豎直虛線),表明過載燃速增量在推進劑表面各點出現穩定狀態的時間幾乎一致;而侵蝕效應引起的轉折點根據燃面所處位置靠后而推遲出現。這是因為由侵蝕(5)式和(6)式可知,侵蝕效應由壓力、速度和溫差等因素決定,在燃燒建立過程中,推進劑通道內流動特性相比壓力變化存在滯后,使得侵蝕效應瞬態- 穩態轉折點不統一;由橫向過載效應(7)式、(8)式和(9)式可知,橫向過載效應不受到縱向內流場影響,因此瞬態- 穩態轉折點出現時刻差異并不明顯。

表3 不同橫向過載下侵蝕效應/過載效應的元量綱燃速增量時變圖

通過比較圖11中100g和-100g時燃速增量占比曲線(計算方式依據(15)式)可以發現:燃速效應所占比重在推進劑燃燒建立過程中是瞬變的;在100g過載下,橫向過載效應在推進劑頭部的比重逐步增加,當燃燒達到穩定時,由過載效應帶來的燃速增量將占比90%,因此正向過載下發動機穩定工作段可以忽略由過載帶來的侵蝕效應;-100g過載下,侵蝕效應在初始到穩定各個階段占絕對主導,推進劑頭部受到的過載效應的影響在穩定燃燒階段出現一定的增加,但占比較少。

圖11 不同橫向過載下燃速增量占比曲線Fig.11 Ratio of burning-rate augmentation under different lateral accelerations

2.3 橫向過載對推進劑火焰傳播過程影響

圖12所示為火焰傳播速度曲線。由圖12可以發現:無論是否處于過載條件,火焰傳播速度均具有明顯的瞬變特性。隨著時間推移,火焰傳播速度不斷增加,在達到峰值過后迅速下降為0 m/s. 100g時,火焰傳播速度峰值達到77 m/s,是0g時峰值的2倍;-100g時火焰傳播速度峰值低于0過載時的傳播速度峰值,二者相差9%. 從火焰傳播時間上看,正向過載的出現較大程度地提前了峰值出現的時間,這是因為正向過載加劇了燃燒室上游推進劑的燃燒,點火燃氣與推進劑燃氣加強下游推進劑的對流換熱,而負向過載則延緩了峰值時間的出現,但是延緩程度較小。

圖12 不同橫向過載下火焰傳播速度曲線Fig.12 Flame-spread speed under different lateral accelerations

為了研究普適的規律,在圖13中將升壓速率、火焰傳播速度和燃面面積變化以各自條件下最大值作為歸一標準進行無量綱處理,可以發現:火焰傳播速度峰值時刻與推進劑表面首次全部點燃時刻以及升壓速率峰值點的時刻幾乎一致,這意味著無堵蓋敞口發動機點火研究中,可以使用實驗中獲得的壓力曲線以及升壓速率去分析和判定推進劑表面燃燒狀況,降低實驗觀察難度。

圖13 不同橫向過載下復合內彈道曲線圖Fig.13 Composite interior ballistic curves under different lateral overloads

3 結論

本文構建了單側推進劑橫向過載下的點火過程物理模型,研究和解釋了高橫向過載下固體火箭發動機過載效應/侵蝕效應耦合規律。通過分析了不同橫向過載下火箭點火升壓變化規律、燃燒變化耦合機制和推進劑表面火焰傳播特性。得到主要結論如下:

1)正向過載下點火壓力峰值升高,負過載下點火壓力峰值降低。在同一個過載下,內孔型燃面在點火過程將出現燃燒面積不對稱的可能,致使實際彈道偏離設計彈道。

2)橫向過載對推進劑的燃燒存在影響,影響程度由所處位置決定,主要表現為:正向過載下,推進劑前段燃燒主要由于過載效應影響,后段主要由侵蝕效應影響;負向過載下,過載效應對推進劑燃燒的影響在變化程度和持續時間上均較小。

3)正向過載加劇推進劑前段燃燒,加劇下游侵蝕效應,縮短火焰傳播時間;負向過載有削弱侵蝕效應作用,但削弱程度較小。

4)內彈道升壓速率可以用來分析推進劑表面燃燒狀況,降低實驗觀測難度。

實際上,發動機尺寸和結構之間的差異將影響上述物理現象的演變程度。本文采用無量綱分析消除了不同條件的幅值差異,所反映規律更客觀。使用單側平面推進劑對橫向過載下點火過程進行研究能清晰地反映出侵蝕和過載效應在點火過程的耦合關系,使得該規律具有更強普適性。推進劑在過載下點火過程的研究結論也為工作段過載環境下發動機內彈道研究提供更接近真實的初始條件。

猜你喜歡
效應發動機模型
一半模型
鈾對大型溞的急性毒性效應
懶馬效應
今日農業(2020年19期)2020-12-14 14:16:52
重要模型『一線三等角』
重尾非線性自回歸模型自加權M-估計的漸近分布
發動機空中起動包線擴展試飛組織與實施
應變效應及其應用
3D打印中的模型分割與打包
新一代MTU2000發動機系列
新型1.5L-Eco-Boost發動機
主站蜘蛛池模板: 无码aⅴ精品一区二区三区| 日韩AV手机在线观看蜜芽| 久久黄色毛片| 久久综合五月| 一级成人a毛片免费播放| 国产真实自在自线免费精品| 成人午夜视频在线| 亚洲啪啪网| 国产在线第二页| 国产一级毛片yw| 婷婷综合在线观看丁香| 三上悠亚在线精品二区| 欧美三级自拍| 免费Aⅴ片在线观看蜜芽Tⅴ| 97国内精品久久久久不卡| www.狠狠| 国产成人无码AV在线播放动漫| 特级毛片8级毛片免费观看| 欧美精品v欧洲精品| 亚洲久悠悠色悠在线播放| 高清码无在线看| 亚洲天天更新| 亚洲成A人V欧美综合| 在线国产三级| 欧美精品成人| 亚洲欧美激情小说另类| 欧美精品成人| 精品欧美视频| 99re经典视频在线| 欧美a级在线| 99精品影院| 2021国产精品自产拍在线观看| 天天视频在线91频| 欧美激情福利| 日本一区二区三区精品国产| 99久久精品视香蕉蕉| 国产91九色在线播放| 成人蜜桃网| 久久一级电影| 国产成人av大片在线播放| 一区二区三区毛片无码| www.91中文字幕| 国产一级做美女做受视频| 夜夜操国产| 久久伊人色| 在线精品亚洲一区二区古装| 亚洲一级毛片在线观| 日韩无码一二三区| 欧美日韩亚洲综合在线观看| 91精品专区国产盗摄| 人人妻人人澡人人爽欧美一区| 久久网综合| 青草精品视频| 中文字幕乱妇无码AV在线| 欧美一区二区三区不卡免费| 无码一区二区波多野结衣播放搜索| 欧美国产另类| AV无码无在线观看免费| 欧美一级在线播放| a网站在线观看| 亚洲香蕉在线| 午夜少妇精品视频小电影| 婷婷午夜天| 在线视频亚洲欧美| 久久久久国色AV免费观看性色| 成人福利在线免费观看| 久久天天躁夜夜躁狠狠| 亚洲美女一区| 久久精品人人做人人爽97| 国产精品乱偷免费视频| 国产精品久久自在自线观看| 国产女主播一区| 欧美在线天堂| 麻豆精品在线视频| 女同国产精品一区二区| 亚洲日本韩在线观看| 2021精品国产自在现线看| 亚洲性影院| 国产精品无码AV中文| 亚洲VA中文字幕| 亚洲色图在线观看| 亚洲男人在线天堂|