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

基于滑移網格的傾轉旋翼機全機干擾流場研究

2021-04-08 15:16:41林沐陽招啟軍趙國慶
航空科學技術 2021年1期

林沐陽 招啟軍 趙國慶

摘要:為了對傾轉旋翼機全機干擾規律及機理進行研究,對某傾轉旋翼機的孤立機身流場特性、懸停及巡航狀態下全機干擾流場特性進行模擬分析。首先,分別生成圍繞機身和旋翼槳葉的貼體網格,然后通過網格搭接得到全機干擾流場計算網格;基于滑移網格方法,建立一套適用于傾轉旋翼機全機氣動干擾特性分析的數值方法,并以有試驗結果對比的ROBIN機身進行了算例驗證;采用建立的數值模擬方法對傾轉旋翼機的孤立機身流場、懸停及巡航狀態下全機干擾流場進行細致的模擬,重點分析了旋翼/機身間的氣動干擾特性。研究發現,懸停時,傾轉旋翼機會出現特有的“噴泉效應”現象,并導致機翼出現周期性下洗載荷,下洗載荷平均值約為旋翼總拉力的12.19%。巡航時,旋翼槳尖渦會向后脫落與機翼相互作用,導致傾轉旋翼機機身阻力出現周期性變化,同時機身阻力平均值也有所增加。

關鍵詞:傾轉旋翼機;全機干擾;噴泉效應;滑移網格;CFD方法

中圖分類號:V11.52文獻標識碼:ADOI:10.19452/j.issn1007-5453.2021.01.017

基金項目:國家自然科學基金(11872211);國家重點實驗室基金(61422200101);江蘇高校優勢學科建設工程資助項目

垂直起降飛行器因為有著對起飛、著陸場地要求較低的優點,在軍用和民用方面均有很大的應用空間[1-2]。傾轉旋翼機是一種兼具垂直起降和高速巡航優勢的新概念可變飛行模態飛行器:在垂直起降狀態,旋翼軸處于豎直位置并由旋翼提供升力;在高速巡航狀態,旋翼軸處于水平位置并由旋翼(螺旋槳)提供推進力,同時由機翼提供升力以克服機身重力[3]。

由于其獨特的設計,傾轉旋翼機的巡航速度遠大于常規直升機,具備滿足多種飛行任務需求的能力[4-5]。但也因此存在復雜的氣動干擾現象[6-8],各個部件間(如旋翼/機翼、旋翼/機身、旋翼/尾翼)會相互干擾,對其氣動和飛行性能產生巨大影響。因此,針對傾轉旋翼機典型飛行狀態流場開展數值模擬及全機氣動干擾分析研究,十分有利于對傾轉旋翼機氣動干擾機理規律的理解和掌握,從而有效減少該干擾現象。

目前針對傾轉旋翼機氣動干擾的數值模擬分析,主要是采用計算流體力學(CFD)方法求解黏性Navier-Stokes(NS)方程,對傾轉旋翼機在不同狀態的干擾流場進行分析。CFD方法具有精度高、能夠有效模擬流動細節等優勢,是目前對傾轉旋翼機全機干擾氣動特性分析的一種重要研究方法。2002年,Romander[9]采用CFD方法計算了不同巡航條件下的V-22傾轉旋翼的氣動特性,發現在某些狀態旋翼槳尖會出現激波和失速現象。2005年,Potsdam等[10]通過CFD方法對孤立旋翼及傾轉旋翼機在懸停狀態下的全機流場進行模擬分析,計算結果為改善傾轉旋翼的氣動性能提供了大量的流場細節。2009年,Lee-Rausch等[11]通過求解基于非結構運動嵌套網格的雷諾平均N-S方程對孤立傾轉旋翼在懸停狀態下的流場進行了細致模擬,并將計算的旋翼氣動性能和槳葉負載變化趨勢與試驗值進行了對比,驗證了方法的有效性。國內,劉全等[12]基于動量源方法和歐拉控制方程建立了一套用于傾轉旋翼機干擾流場計算的數值模擬方法,并對傾轉旋翼機旋翼/機身的組合流場在懸停和巡航狀態下的氣動特性進行計算分析。李鵬等[13]采用新型運動嵌套網格方法建立了一套用于求解傾轉旋翼機非定常流場的數值模擬方法,并基于所建立的方法對旋翼/機身氣動干擾現象、渦流動機理及相互作用機理開展了分析研究。但是上述研究大多以孤立旋翼或旋翼/機翼流場特性為研究對象,對于傾轉旋翼機的全機干擾流場特性分析研究較少。

因此,本文將結合滑移網格技術,建立適用于傾轉旋翼機全機干擾流場模擬及氣動特性分析的高精度數值方法,并以某傾轉旋翼機為分析對象,對其孤立機身流場和典型飛行狀態下全機干擾流場進行模擬,對比分析有無旋翼干擾情況下的流場特性和機身氣動載荷差異,對傾轉旋翼機在典型飛行狀態下的全機干擾機理進行初步探究。

1計算模型及方法

1.1幾何模型及網格劃分

基于CATIA軟件生成傾轉旋翼機幾何模型,表1給出了本文計算對象傾轉旋翼機的基本參數。

本文選取半徑為15L、高為20L的圓柱體作為流場計算域,其中L是機身長度。計算域分為遠場-機身區域和旋翼區域,包括機身網格、遠場網格、對稱面網格和交界面網格,整個計算域網格示意圖如圖1所示。由于傾轉旋翼全機干擾流場具有對稱性,為了減少網格量,僅對半機身/旋翼流場進行數值模擬。

為了提高對傾轉旋翼機流場數值模擬的精度,對主要氣動部件(旋翼、機翼、平尾、垂尾)及曲率較大的部件(短艙)網格進行加密處理。圖2為傾轉旋翼機機身表面網格,面網格總數為24萬,其中旋翼表面網格占比最大,約占面網格總數的70%。

為了有效模擬機身邊界層流動,在機身表面生成邊界層網格,第一層網格的高度為0.005mm,網格增長率為1.5,共生成15層邊界層網格。最終得到適用于全機干擾流場模擬的計算網格。

圖3給出了懸停狀態下的全機干擾網格結構圖。由于旋翼需要進行旋轉運動,因此結合滑移網格技術,網格分為運動網格和靜止網格,運動網格對應前文所述的旋翼區域,在數值模擬過程中,運動網格做旋轉運動,并通過滑移交界面進行運動網格與靜止網格的流場信息傳遞,需要保證旋翼區域和遠場-機身區域的交界面網格密度一致,防止出現交界面上的某個網格對應另一個交界面多個網格的情況。

1.2計算方法

本文通過求解雷諾平均的N-S方程進行全機干擾流場的計算,湍流模型選取k -ω湍流模型。在笛卡兒坐標下的N-S方程的形式如式(1)所示:

本文基于Fluent軟件,選用SIMPLE壓力速度耦合迭代方法,壓力離散采用二階精度。對于壁面邊界層,選擇自由剪切流,適用于存在逆壓力梯度時的流動和分離,邊界條件采用壓力遠場邊界條件。

1.3算例驗證

對本文所建立的CFD方法進行算例驗證計算,以驗證方法的有效性。對流場計算方法進行驗證,驗證算例對象選擇ROBIN機身,采用所建立的CFD方法對該機身在巡航狀態下的流場進行數值模擬,并將計算結果與試驗數據[14]進行對比。表2是驗證算例的基本參數。

圖4給出了對應ROBIN機身上的凸臺起始端和凸臺中段位置的機身橫截面表面壓力(壓強)系數分布,并與試驗值進行了對比。可以看出在各個截面位置的機身表面壓力系數分布的計算值與試驗值吻合較好,表明本文所建立的CFD方法能有效地對機身表面壓力及機身附近流場進行數值模擬。

2計算結果分析

選取傾轉旋翼機計算模型參數見表1,各個計算狀態中機身俯仰角、側滑角和偏航角均為0°,旋翼總距θ0= 10o,為了簡化計算,忽略旋翼槳轂,取旋翼根切為25%R。

2.1孤立機身流場氣動特性分析

開展傾轉旋翼機孤立機身氣動特性計算,獲得機身阻力以及機身流場特性,為進一步了解全機氣動干擾機理提供參考基礎。以巡航速度為自變量,對孤立機身在不同巡航速度下的氣動特性進行模擬,表3為計算結果。

進一步對機身受力進行無量綱處理,得到機身各個氣動力系數隨巡航速度變化的曲線圖,結果如圖5所示。在巡航速度為40~240km/h的范圍內,阻力系數隨速度變化較小,基本保持在0.0372,升力系數在速度較小時隨巡航速度的增大而減小,當巡航速度繼續增大,升力系數基本維持不變。

圖6給出了孤立機身在巡航(160km/h)狀態下的等渦量圖。可以看出氣流在機翼后端、尾翼后端及機翼與機身連接處后方的旋渦脫落較明顯。

2.2巡航狀態全機干擾流場氣動特性分析

圖7為在巡航速度為160km/h狀態下計入旋翼影響后的機身阻力隨旋翼旋轉方位角變化曲線圖,由于旋翼/機身全機干擾流場的非定常特性,機身阻力出現非定常波動并且不完全是周期波動,這表明巡航狀態下,旋翼/機翼流場相互干擾作用機理更為復雜。

對不同巡航速度下計入旋翼影響的全機干擾流場進行數值模擬,計算得到的機身氣動力(見表4),表中數值為旋翼旋轉一周時機身氣動力的平均值。

圖8給出了計入旋翼影響前后的機身阻力隨巡航速度變化曲線圖,可以看出在計入旋翼影響后機身阻力增加,這主要是由于旋翼對機身的干擾引起的。

圖9為考慮旋翼干擾前后的機身升阻比隨巡航速度變化曲線圖。在考慮旋翼干擾前,機身升阻比隨巡航速度增加而略微減小,但基本在2.6左右。考慮旋翼干擾后,機身升阻比隨巡航速度增加而增加。

進一步分析旋翼在傾轉旋翼機巡航狀態下旋翼下洗流對機翼的影響,圖10給出了巡航速度為160km/h時孤立機身流場和旋翼/機身干擾流場中機翼不同展向位置剖面的壓強系數分布對比。可以看出在計入旋翼影響后,位于旋翼內側后方的機翼剖面上表面的負峰值增大,機翼升力增大,這是由于旋翼對后方氣流起到加速的作用;機翼翼根處由于旋翼的影響,上表面峰值增大,升力增加。

圖11給出在巡航速度160km/h下,傾轉旋翼機全機干擾流場的渦量等值面圖,對比圖6發現,旋翼脫落的槳尖渦會隨著來流向后運動,最終與機翼發生相互干擾,在圖11的俯視圖中可以明顯看出旋翼/機身相互干擾作用的區域,這是與孤立機身流場的明顯區別,同時也是導致兩者機翼壓強分布存在差異的主要原因。

圖12給出了計入旋翼影響前后的機身表面對稱截面的壓強系數分布對比,可以看出計入旋翼影響前后的機身表面壓強系數差別較小,僅在機尾處存在略微差別,這主要是由于兩者垂尾處的面網格劃分存在一定的差異。

進一步分析全機干擾流場與孤立機身流場的尾翼氣動特性的差異。本文計算狀態中偏航角均為0°,垂尾的側向力和升力幾乎為0,因此僅對在計入旋翼影響前后的垂尾阻力進行對比分析,結果如圖13所示。可以看出,計入旋翼影響后,垂尾阻力變化不大,這主要是由于傾轉旋翼機采用的是T形尾,垂尾處于機身對稱截面上,受到旋翼下洗流影響較小,這與前文的分析一致。

圖14給出了計入旋翼影響前后的平尾升/阻力對比圖。可以看出,計入旋翼影響后平尾阻力變化不大,在巡航速度較大時略微減小;但是平尾升力出現明顯減小,這表明旋翼對平尾氣動特性產生干擾作用。

圖15為巡航速度160km/h下傾轉旋翼機全機干擾流場的流線圖,左側為空間三維流線,右側為垂直于X軸且經過點(1.3m,0m,0m)的截面上的面流線圖,由左側三維流線可以看出,當氣流經過機翼上表面時發生明顯的偏折,同時動壓增大靜壓減小,表明氣體加速運動;由右側可以看出旋翼后方氣流流動較為復雜。

2.3懸停狀態全機干擾流場氣動特性分析

圖16為懸停狀態下垂直于X軸且經過旋翼旋轉軸的截面上的面流線圖,對該截面中氣流的運動情況進行分析。可以看出,在傾轉旋翼機懸停時,旋翼流場呈現漏斗狀,旋翼上方的氣流受到吸引加速向下運動,位于旋翼外側的氣流由于沒有機翼的阻擋會形成脫落的槳尖渦,而位于內側的氣流在經過旋翼繼續向下運動時,由于機翼的阻擋不能無障礙地向下流動,導致機翼左右側均會出現氣流沿機翼展向流動,兩股氣流最終在機身對稱面相遇后被卷起向上運動,隨著氣流繼續向上運動受到旋翼的吸引作用,最終在機翼上方形成一個循環流動,形成“噴泉效應”現象。

圖17為相同截面的渦量圖,旋翼槳尖及槳根均存在渦量較大的區域,分別對應槳尖渦和槳根渦,后者是由于計算模型簡化,對旋翼進行根切產生的。可以看到內側旋翼脫落的槳尖渦向下運動,并由于“噴泉效應”的影響,在機翼上方出現渦量值較大的區域。此外機翼端部及機翼/機身連接處也由于下洗流的影響使得流動變得復雜,導致該區域渦量較大。

進一步給出該截面的壓力云圖,如圖18所示,可以看出左右旋翼內側與機翼之間均存在較大的高壓區域,對應誘導速度云圖,該區域下洗流減速,氣流動壓減小,靜壓增大,這也是導致機翼產生下洗載荷的主要原因,此外由于機翼的阻擋作用在機翼下方還存在較大的負壓區域,該位置氣流分離現象顯著。

圖19是旋翼旋轉一周時旋翼及全機升力的變化曲線圖,可以看出升力呈現周期性變化,在旋翼旋轉一圈的時間里升力變化曲線存在三個峰值,這是因為槳葉片數為三片。旋翼的平均升力為1173.95N,全機的平均升力為1032.03N,全機升力小于旋翼升力,這主要是由于旋翼下洗流打到機翼,產生額外的下洗載荷。

提取出機身的下洗載荷變化如圖20所示,可以看出下洗載荷也呈現周期性變化,下洗載荷平均值為143.10N,約為旋翼升力的12.19%。

圖21為懸停狀態下計入旋翼影響后機翼上下洗載荷沿翼展方向的分布圖。由圖可知曲線上存在兩個極大值,第一個極大值出現在機翼與機身的連接處(翼展位置0.3m處),位于該位置的機翼剖面下表面存在因氣流分離產生的負壓區,使得機翼載荷較大。第二個極大值出現在旋翼槳尖下方機翼(翼展位置0.8m處),這是由于旋翼槳尖附近的氣流向下流動速度大,并撞擊到機翼表面,導致位于該位置的機翼剖面上表面存在高壓區,機翼載荷較大。

圖22給出了傾轉旋翼機在懸停狀態下全機干擾流場渦量值為160的渦量等值面圖。由圖可以看出,旋翼槳尖渦撞擊在機翼上表面,發生相互干擾;在內側旋翼槳尖到機身對稱面之間的機翼上方可以看到氣流存在橫向流動現象,兩側的氣流最終在機身對稱面相遇,機翼上方在機身對稱面位置可以看到凸起的區域,這是機翼兩側氣流碰撞形成“噴泉效應”現象導致的,同時部分氣流沿機身上表面分別向機頭和機尾流動。

3結束語

通過研究,可以得出以下結論:

(1)本文結合滑移網格方法建立了一套能有效模擬傾轉旋翼全機干擾流場CFD方法,并運用所建立的CFD方法對某傾轉旋翼機的孤立機身流場、巡航及懸停的全機干擾流場進行模擬分析。

(2)懸停時,由于旋翼/機翼的相互干擾,傾轉旋翼機會出現特有的“噴泉效應”現象,并導致機身出現周期性下洗載荷,其平均值約為旋翼總拉力平均值的12.19%。

(3)巡航時,旋翼槳尖渦會向后脫落與機翼相互作用,導致傾轉旋翼機機身阻力出現周期性變化,同時平均值也有所增加。

本文下一步工作重點將考慮對傾轉旋翼機過渡狀態下的全機干擾流場進行計算分析,并基于本文初步總結的全機干擾規律,開展傾轉旋翼機氣動布局設計優化。

參考文獻

[1]李耐和.美軍下一次戰爭需要的10項關鍵技術[J].新時代國防, 2011(12): 1-8. Li Naihe. The US army ten key techniques needed for next warfare[J]. New Time Defence, 2011(12): 1-8. (in Chinese)

[2]劉凱,葉賦晨.垂直起降飛行器的發展動態和趨勢分析[J].航空工程進展, 2015, 6(2): 127-138. Liu Kai, Ye Fuchen. Review and analysis of recent development for VTOL vehicles[J]. Advances in Aeronautical Science and Engineering, 2015, 6(2): 127-138. (in Chinese)

[3]王偉,段卓毅,周林.傾轉旋翼機設計特點及難點淺析[J].航空科學技術, 2015, 26(3): 1-4. Wang Wei, Duan Zhuoyi, Zhou Lin. Brief analysis on the design features and difficulties of tiltrotor[J]. Aeronautical Science & Technology, 2015, 26(3): 1-4. (in Chinese)

[4]陳恒,左曉陽,張玉琢.傾轉旋翼飛機技術發展研究[J].飛行力學, 2007(3): 5-8. Chen Heng, Zuo Xiaoyang, Zhang Yuzhuo. Tiltrotor aircraft key technology developing research[J]. Flight Dynamics, 2007(3): 5-8. (in Chinese)

[5]尹昱康,王紅州,蔡恒欲,等.傾轉旋翼機的發展[J].兵器裝備工程學報, 2018, 39(2): 111-113. Yin Yukang, Wang Hongzhou, Cai Hengyu, et al. Development and technical characteristics of tiltrotor aircraft[J]. Journal of Ordnance Equipment Engineering, 2018, 39(2): 111-113. (in Chinese)

[6]吳希明,仲唯貴,陳平劍.傾轉旋翼機氣動設計技術[J].航空科學技術, 2012,21(4): 17-24. Wu Ximing, Zhong Weigui, Chen Pingjian. Aerodynamic design technology of tilt-rotor aircraft[J]. Aeronautical Science&Technology, 2012,21(4): 17-24. (in Chinese)

[7]薛蒙,孫強.傾轉旋翼機軍事需求與關鍵技術分析[J].直升機技術, 2020(1): 47-49. Xue Meng, Sun Qiang. Tiltrotor military requirement and critical technology analysis[J]. Helicopter Technique, 2020(1): 47-49. (in Chinese)

[8]徐敏.傾轉旋翼機的發展與關鍵技術綜述[J].直升機技術, 2003(2): 40-44. Xu Min. Review of the development and key technologies of tiltrotor aircraft[J]. Helicopter Technique, 2003(2): 40-44. (in Chinese)

[9]Romander E A. Computational simulation of propellers in cruise[R]. ICAS Paper,2002.

[10]Potsdam M A,Strawn R C. CFD simulations of tiltrotor configurations in hover[J]. Journal of the American Helicopter Society,2005,50(1):82-94.

[11]Lee-Rausch E M,Biedron R T. Simulation of an isolated tiltrotor in hover with an unstructured overset-grid RANS solver[C]// The 65th Annual Forum of the American Helicopter Society. TX,2009:27-29.

[12]劉全.懸停和前飛狀態傾轉旋翼機流場的數值分析[D].南京:南京航空航天大學, 2009. Liu Quan. Numerical analysis of flow field of tilt-rotor aircraft in hover and cruise[D]. Nanjing: Nanjing University of Aeronautics andAstronautics, 2009. (in Chinese)

[13]李鵬.傾轉旋翼機非定常氣動特性分析及氣動設計研究[D].南京:南京航空航天大學, 2016. Li Peng. Analysis of unsteady aerodynamic characteristics and aerodynamic design of tilt-rotor aircraft[D]. Nanjing: Nanjing University ofAeronautics andAstronautics, 2016. (in Chinese)

[14]McKee J W,Naeseth R L. Experimental investigation of the drag of flat plates and cylinders in the slipstream of a hovering rotor[R]. NACATN 4239,1958.(責任編輯余培紅)

作者簡介

林沐陽(1996-)男,碩士。主要研究方向:直升機空氣動力學、全機氣動干擾。

招啟軍(1977-)男,博士,教授。主要研究方向:直升機計算流體力學、直升機空氣動力學、氣動聲學等。

Tel:025-84893753

E-mail:zhaoqijun@nuaa.edu.cn

趙國慶(1984-)男,博士,講師。主要研究方向:直升機空氣動力學、計算流體力學。

Research on the Whole Aircraft Interference Flow Field of Tilt-rotor Aircraft Based on Sliding Grid

Lin Muyang,Zhao Qijun*,Zhao Guoqing

National Key Laboratory of Rotorcraft Aeromechanics,Nanjing University of Aeronautics and Astronautics,Nanjing 210016,China

Abstract: In order to study the law and mechanism of the whole aircraft interference of the tilt-rotor aircraft, the flow field characteristics of the isolated fuselage and the whole aircraft interference in hovering and cruising of a certain type of tilt-rotor UAV were simulated and analyzed. First, the grids around the fuselage and rotor blades were generated respectively, and then the calculation grid of the whole aircraft interference flow field was obtained through grid overlap. Then, based on the sliding grid method, a set of numerical methods suitable for the analysis of the aerodynamic interference characteristics of the tilt-rotor aircraft was established, and the method was verified by the example of ROBIN fuselage with experimental results. Then, the numerical simulation method was used to simulate the flow field of the isolated fuselage and the whole aircraft interference in hovering and cruising of the tilt-rotor aircraft. The aerodynamic interference characteristics between the rotor and the fuselage were mainly analyzed. The analysis shows that: When hovering, there is a unique "fountain effect" phenomenon in tilt-rotor aircraft flow field, and the phenomenon cause periodic downwash loads on the wings. The average downwash load is approximately 12.19% of the total rotor pull. When cruising, the tip vortex of the rotor blade will fall off backwards and interact with the wing, causing periodic changes in the fuselage drag and an increase in the average value of the fuselage drag.

Key Words: tilt-rotor aircraft; whole aircraft interference; fountain effect; sliding grid; CFD method

主站蜘蛛池模板: 日本爱爱精品一区二区| 精品久久久久久中文字幕女 | 国产精品夜夜嗨视频免费视频| 欧美人人干| 欧美一级高清视频在线播放| 国产成人综合久久精品尤物| 青草国产在线视频| 国产青青草视频| 欧美午夜在线观看| 久久国产高潮流白浆免费观看| 91欧美亚洲国产五月天| 久久九九热视频| 国产乱子伦视频在线播放| 国产成人无码综合亚洲日韩不卡| 综1合AV在线播放| 热热久久狠狠偷偷色男同| 97一区二区在线播放| 亚洲男人的天堂在线观看| 亚洲国产日韩视频观看| 亚洲一级毛片在线观播放| 国产99免费视频| 成人在线观看一区| 国产18页| 综合久久五月天| 在线亚洲精品自拍| AV在线天堂进入| 熟女成人国产精品视频| 91极品美女高潮叫床在线观看| 丁香五月婷婷激情基地| 国产福利一区在线| 国产自在线播放| 成人年鲁鲁在线观看视频| 就去吻亚洲精品国产欧美| 国产高清无码第一十页在线观看| 成人国产精品网站在线看| 亚洲一区网站| 欧美激情综合一区二区| 色久综合在线| 亚洲一级毛片在线观| 亚洲精品日产精品乱码不卡| 日韩精品高清自在线| 成人永久免费A∨一级在线播放| 国产成人综合久久精品下载| 日韩色图区| 日韩a级毛片| 伊人中文网| 成年人视频一区二区| 国产全黄a一级毛片| 欧美日在线观看| 美女高潮全身流白浆福利区| av一区二区无码在线| 国产精品久线在线观看| 日本91视频| 日韩精品成人网页视频在线| 99中文字幕亚洲一区二区| 国产浮力第一页永久地址| 在线看片国产| 草逼视频国产| 久久国产免费观看| 久久久噜噜噜久久中文字幕色伊伊| 久久毛片基地| 欧美国产日本高清不卡| 精品小视频在线观看| 五月天久久综合| 国产呦视频免费视频在线观看| 色成人亚洲| 91麻豆国产视频| 91免费片| 国产免费久久精品99re丫丫一| 色135综合网| 欧美成人午夜视频免看| 亚洲男人天堂2018| 免费一级全黄少妇性色生活片| 九九视频免费看| 亚洲第一天堂无码专区| 国产情侣一区二区三区| 91色爱欧美精品www| 国产区福利小视频在线观看尤物| 尤物精品视频一区二区三区| 午夜日本永久乱码免费播放片| 中文无码伦av中文字幕| 国产极品美女在线观看|