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

鴨式旋翼/機(jī)翼飛行器轉(zhuǎn)換末段氣動(dòng)特性

2011-03-15 12:37:22李毅波馬東立
關(guān)鍵詞:平尾

李毅波 馬東立

(北京航空航天大學(xué) 航空科學(xué)與工程學(xué)院,北京 100191)

牛凌宇

(北京系統(tǒng)工程研究所,北京 100101)

鴨式旋翼/機(jī)翼飛行器轉(zhuǎn)換末段氣動(dòng)特性

李毅波 馬東立

(北京航空航天大學(xué) 航空科學(xué)與工程學(xué)院,北京 100191)

牛凌宇

(北京系統(tǒng)工程研究所,北京 100101)

采用數(shù)值模擬方法研究鴨式旋翼/機(jī)翼(CRW,Canard Rotor/Wing)飛行器在轉(zhuǎn)換過(guò)程末段,旋翼轉(zhuǎn)速極低時(shí)全機(jī)氣動(dòng)特性變化規(guī)律及其產(chǎn)生原因.給出了旋翼旋轉(zhuǎn)一周時(shí),全機(jī)氣動(dòng)力、氣動(dòng)力矩、焦點(diǎn)位置變化規(guī)律,對(duì)此布局形式,轉(zhuǎn)換過(guò)程末段全機(jī)升力、阻力變化幅度可達(dá)10.7%,3.7%,焦點(diǎn)可移動(dòng)0.6m.研究顯示:旋翼處于前后不對(duì)稱流場(chǎng)及旋翼處于不同方位角時(shí)對(duì)機(jī)體的不對(duì)稱干擾是氣動(dòng)力與氣動(dòng)力矩變化原因,旋翼與平尾升力線斜率變化、旋翼自身焦點(diǎn)位置變化導(dǎo)致了全機(jī)焦點(diǎn)移動(dòng).

鴨式旋翼/機(jī)翼飛行器;氣動(dòng)特性;流場(chǎng)干擾;飛行模式轉(zhuǎn)換

鴨式旋翼/機(jī)翼(CRW,Canard Rotor/Wing)飛行器是一種新型垂直起降、高速巡航飛行器,其主要特征是具有一副可停轉(zhuǎn)旋翼,旋翼可高速旋轉(zhuǎn)為垂直起降提供拉力,又可鎖定成固定翼為高速巡航提供升力.CRW飛行器與高速直升機(jī)、傾轉(zhuǎn)旋翼飛行器相比,具有更寬的飛行包線,在軍民用領(lǐng)域的應(yīng)用前景十分廣泛[1],但由于CRW飛行器存在空氣動(dòng)力學(xué)/飛行動(dòng)力學(xué)/飛行控制系統(tǒng)的復(fù)雜耦合,其進(jìn)入實(shí)用化還需開(kāi)展大量研究.

美國(guó)從20世紀(jì)90年代初開(kāi)始對(duì)CRW飛行器展開(kāi)系統(tǒng)地理論與實(shí)驗(yàn)研究,1998—2005年進(jìn)行了X-50A“蜻蜓”項(xiàng)目,演示驗(yàn)證CRW飛行器關(guān)鍵技術(shù).韓國(guó)航空研究院對(duì)CRW飛行器也開(kāi)展了研究[2-3].近幾年來(lái),國(guó)內(nèi)也對(duì) CRW 飛行器相關(guān)技術(shù)進(jìn)行了研究,主要集中在噴氣旋翼的研究上[4-5],對(duì)轉(zhuǎn)換過(guò)程氣動(dòng)特性的研究尚處于空白.

CRW飛行器轉(zhuǎn)換過(guò)程氣動(dòng)特性具有強(qiáng)烈的非定常與非線性特點(diǎn),是CRW飛行器研究的重點(diǎn).轉(zhuǎn)換過(guò)程可分為前、中、末3個(gè)階段,轉(zhuǎn)換前段,鴨翼與平尾處于豎起狀態(tài),旋翼承擔(dān)絕大部分升力,全機(jī)氣動(dòng)特性主要受旋翼流場(chǎng)影響.轉(zhuǎn)換中段,鴨翼與平尾承擔(dān)的升力逐漸增大,旋翼升力逐漸減小,全機(jī)氣動(dòng)特性受旋翼流場(chǎng)與來(lái)流流場(chǎng)綜合影響,非定常氣動(dòng)力表現(xiàn)為與旋翼旋轉(zhuǎn)相應(yīng)的周期性高頻振蕩.轉(zhuǎn)換末段,鴨翼與平尾承擔(dān)全部升力,旋翼總距為0,轉(zhuǎn)速極低,此時(shí)全機(jī)氣動(dòng)特性主要受來(lái)流流場(chǎng)影響,與轉(zhuǎn)換前段與中段相比,非定常氣動(dòng)特性及其機(jī)理存在較大差別[6].研究轉(zhuǎn)換末段氣動(dòng)特性對(duì)CRW飛行器從直升機(jī)模式向固定翼模式的平穩(wěn)安全轉(zhuǎn)換具有重要意義.

本文用數(shù)值方法模擬CRW飛行器處于轉(zhuǎn)換過(guò)程末段,旋翼轉(zhuǎn)速極低,鴨翼與平尾產(chǎn)生全部升力情況下,全機(jī)氣動(dòng)力、氣動(dòng)力矩及焦點(diǎn)位置隨旋翼旋轉(zhuǎn)變化規(guī)律,并分析其變化原因.

1 分析模型與方法

1.1 幾何模型

CRW飛行器幾何外形與坐標(biāo)選取如圖1所示.采用三翼面布局形式,從前至后分別為鴨翼、旋翼/機(jī)翼、平尾.鴨翼與平尾在起降及轉(zhuǎn)換前段處于豎起狀態(tài),以減小旋翼拉力損失,轉(zhuǎn)換末段兩者承擔(dān)全部升力.旋翼為右旋旋翼,定義旋翼處于y軸正方向所在半平面時(shí)為右旋翼,與之相對(duì)的為左旋翼.為了防止旋翼停轉(zhuǎn)后須將一側(cè)機(jī)翼翻轉(zhuǎn)的情況,保證轉(zhuǎn)換過(guò)程平穩(wěn),旋翼采用無(wú)彎度橢圓翼型.

圖1 CRW飛行器

轉(zhuǎn)換過(guò)程末段,鴨翼襟翼偏轉(zhuǎn)角20°,全動(dòng)鴨翼偏轉(zhuǎn)角14°,平尾偏轉(zhuǎn)角10°,保證力及力矩平衡,見(jiàn)表1.來(lái)流迎角0°,來(lái)流馬赫數(shù) Ma=0.197,雷諾數(shù)Re=2.5×106,參考面積為旋翼面積.

表1 轉(zhuǎn)換末段翼面偏轉(zhuǎn)角度 (°)

1.2 數(shù)值計(jì)算方法

轉(zhuǎn)換過(guò)程末段,旋翼轉(zhuǎn)速極低,槳尖速度相對(duì)于來(lái)流速度很小,分析時(shí)可將旋翼固定于一系列方位角位置,研究每個(gè)方位角下全機(jī)氣動(dòng)特性.

采用有限體積法求解三維可壓雷諾平均N-S(Navier-Stokes)方程,三維任意坐標(biāo)系下的N-S方程[7]為

其中,Q為守恒矢量;E,F(xiàn),G為3個(gè)坐標(biāo)方向的對(duì)流通量;Ev,F(xiàn)v,Gv為3個(gè)坐標(biāo)方向的粘性通量.對(duì)流項(xiàng)采用一階迎風(fēng)格式,擴(kuò)散項(xiàng)采用中心差分格式,壓力和速度耦合采用SIMPLE算法.

湍流模型采用Spalart-Allmaras模型,通過(guò)對(duì)自由剪切流、高雷諾數(shù)時(shí)的近壁區(qū)流動(dòng)、有限雷諾數(shù)時(shí)的近壁區(qū)流動(dòng)、包含層流區(qū)與轉(zhuǎn)捩的流動(dòng)等4種狀態(tài)的分析與建模,得到最終的方程[8]:

其中,封閉系數(shù) Cb1=0.135 5,σ =2/3;Cb2=0.622;K=0.41;Cw1=Cb1/K2+(1+Cb2)σ.

流動(dòng)區(qū)域?yàn)?0×10×8倍旋翼直徑,采用四面體網(wǎng)格對(duì)流場(chǎng)區(qū)域進(jìn)行劃分,用到的邊界類型包括壓力遠(yuǎn)場(chǎng)邊界以及無(wú)滑移絕熱固壁邊界.

將本文計(jì)算方法所得結(jié)果與文獻(xiàn)[6]中類似布局飛行器氣動(dòng)數(shù)據(jù)對(duì)比,見(jiàn)圖2,可見(jiàn)旋翼處于不同方位角時(shí),本文結(jié)果與文獻(xiàn)數(shù)據(jù)較為吻合.

圖2 全機(jī)升力系數(shù)C L隨方位角ψ的變化

2 轉(zhuǎn)換末段氣動(dòng)特性

2.1 氣動(dòng)力變化

旋翼轉(zhuǎn)至不同方位角時(shí),全機(jī)升力系數(shù)CL、阻力系數(shù)CD、側(cè)力系數(shù)CC的變化規(guī)律如圖3所示,氣動(dòng)力呈周期性變化.在1/2周期內(nèi)(ψ=0°~90°),全機(jī) CL單調(diào)下降,0°時(shí) CL最大,90°時(shí)最小;CD總體呈上升趨勢(shì),方位角增大,CD增大;CC在45°時(shí)最大,0°與90°時(shí)為0.

圖3 全機(jī)氣動(dòng)力隨方位角的變化

定義一個(gè)周期內(nèi)氣動(dòng)力系數(shù)平均值:

變化幅度:

相對(duì)變化幅度:

其中,Ct為 t時(shí)刻氣動(dòng)力系數(shù)瞬態(tài)值;Cn,max,Cn,min分別為一個(gè)周期內(nèi)氣動(dòng)力系數(shù)的最大值與最小值;T為周期.

根據(jù)上述定義,可知在一個(gè)周期內(nèi),全機(jī)升力系數(shù)變化幅度為0.12,相對(duì)變化10.7%,阻力系數(shù)變化幅度0.0087,相對(duì)變化3.7%,側(cè)力系數(shù)變化幅度0.026,如表2所示.

表2 氣動(dòng)力均值與幅值

全機(jī)升力系數(shù)變化主要由旋翼升力系數(shù)變化引起,如圖4所示.在方位角0°~90°半周期內(nèi),旋翼升力系數(shù)變化具有3個(gè)特征:①左右旋翼均產(chǎn)生負(fù)升力;②左右旋翼升力系數(shù)均隨方位角增大而減小;③相同時(shí)刻,右旋翼負(fù)升力絕對(duì)值大于左旋翼負(fù)升力絕對(duì)值.

圖4 部件升力系數(shù)隨方位角的變化

旋翼采用對(duì)稱橢圓翼型,總距為0時(shí)旋翼升力應(yīng)為0,但在轉(zhuǎn)換末段,旋翼實(shí)際產(chǎn)生周期性變化負(fù)升力,這一現(xiàn)象與原設(shè)想的旋翼完全卸載存在差別,與鴨翼、機(jī)身在旋翼旋轉(zhuǎn)平面誘導(dǎo)的不均勻流場(chǎng)有關(guān).

在旋翼旋轉(zhuǎn)平面上,速度分布左右對(duì)稱,前后不對(duì)稱,見(jiàn)圖5.旋轉(zhuǎn)平面大部分處于下洗流場(chǎng)中,只在機(jī)身前部狹窄區(qū)域內(nèi)存在上洗流場(chǎng).

圖5 旋翼旋轉(zhuǎn)平面速度分布(無(wú)旋翼)

將旋翼升力系數(shù)表示為方位角為ψ的函數(shù):

其中,α為來(lái)流迎角,轉(zhuǎn)換末段α=0°;εψ為鴨翼、機(jī)身誘導(dǎo)的下洗角;α0為旋翼零升迎角,旋翼采用零彎度對(duì)稱橢圓翼型,α0=0°.故可得

可見(jiàn)旋翼升力系數(shù)在不同方位角下均為負(fù).

在1/2周期內(nèi),方位角增大,旋翼升力系數(shù)減小,原因?yàn)?①ψ=0°時(shí),旋翼部分處于機(jī)身上洗流中,ψ增大,上洗逐漸變?yōu)橄孪矗蚁孪唇侵饾u增大;②ψ由0°增至90°,相當(dāng)于旋翼由斜掠翼變?yōu)槠街币恚€斜率增大.由式(4)可得旋翼升力隨ψ增大而減小.

旋翼處于前后不對(duì)稱流場(chǎng)中導(dǎo)致左右旋翼升力系數(shù)不對(duì)稱.在ψ=0°~90°范圍內(nèi),左旋翼(后行)處于前部流場(chǎng),受機(jī)身上洗影響,局部迎角與局部速度較大,升力較大;右旋翼(前行)處于后部流場(chǎng),受機(jī)身對(duì)氣流的阻滯與鴨翼下洗影響,升力較小.

旋翼在方位角由0°轉(zhuǎn)到90°過(guò)程中,迎風(fēng)面積增大,同時(shí)旋翼與機(jī)身干擾增大,全機(jī)阻力系數(shù)呈上升趨勢(shì).

側(cè)力系數(shù)在ψ=0°~90°范圍內(nèi)先增大后減小,ψ=45°時(shí)最大.機(jī)身與垂尾是全機(jī)側(cè)力變化的主要來(lái)源,見(jiàn)圖6.機(jī)身側(cè)力的變化主要由旋翼處于不同方位角時(shí)對(duì)機(jī)身的干擾不同引起,在ψ =0°與90°時(shí),旋翼-機(jī)身組合體左右對(duì)稱,側(cè)力為0.ψ≠0°與90°時(shí),旋翼-機(jī)身組合體左右不對(duì)稱,機(jī)身產(chǎn)生不對(duì)稱側(cè)力.

圖6 部件側(cè)力系數(shù)隨方位角的變化

2.2 氣動(dòng)力矩變化

圖7表示全機(jī)力矩系數(shù)隨方位角變化,在1/2周期內(nèi),俯仰力矩系數(shù)單調(diào)減小,ψ=90°時(shí)為0,此時(shí)全機(jī)處于配平狀態(tài);滾轉(zhuǎn)力矩系數(shù)與偏航力矩系數(shù)絕對(duì)值先增大后減小,在0°與90°位置時(shí)為0,45°位置左右達(dá)到峰值.

圖7 全機(jī)力矩系數(shù)隨方位角的變化

俯仰力矩系數(shù)的變化主要由旋翼與平尾引起.圖8顯示了ψ=0°時(shí)旋翼對(duì)稱面壓力分布,表現(xiàn)為:較大正升力-負(fù)升力-較小正升力,總體效果抬頭力矩.ψ=90°時(shí),旋翼壓力中心對(duì)參考點(diǎn)的力臂長(zhǎng)度較短,對(duì)全機(jī)俯仰力矩的貢獻(xiàn)可忽略.因此旋翼從0°轉(zhuǎn)到90°,抬頭力矩逐漸減小直至0.

平尾俯仰力矩系數(shù)可表示為

其中,CMy0,h為平尾零升力矩系數(shù);CL,h(ψ)為平尾升力系數(shù),是旋翼方位角的函數(shù);ah為平尾力臂,由于力臂較長(zhǎng),故平尾升力系數(shù)產(chǎn)生微小變化可引起力矩系數(shù)較大變化.由圖4可知,ψ=0°時(shí)平尾升力較小,產(chǎn)生低頭力矩較小,ψ=90°時(shí)升力較大,產(chǎn)生低頭力矩較大.

滾轉(zhuǎn)力矩系數(shù)變化同樣與旋翼升力系數(shù)分布及力臂變化有關(guān).定義旋翼升力差為

旋翼升力差隨方位角的變化如圖9所示,差值在ψ=45°時(shí)達(dá)到最大.由于左右旋翼力臂相等,因此右旋翼產(chǎn)生的負(fù)滾轉(zhuǎn)力矩大于左旋翼產(chǎn)生的正滾轉(zhuǎn)力矩,總滾轉(zhuǎn)力矩為負(fù).

圖9 左右旋翼升力差隨方位角的變化

偏航力矩系數(shù)的變化與機(jī)身、垂尾側(cè)力系數(shù)的變化有關(guān).機(jī)身側(cè)力在ψ=45°達(dá)到最大,且前半機(jī)身產(chǎn)生負(fù)側(cè)力,后半機(jī)身產(chǎn)生正側(cè)力,兩者均產(chǎn)生正偏航力矩.垂尾側(cè)力對(duì)全機(jī)偏航力矩的影響與機(jī)身相反.

2.3 焦點(diǎn)位置變化

轉(zhuǎn)換末段旋翼處于不同方位角時(shí),旋翼與全機(jī)焦點(diǎn)位置發(fā)生較大變化,全機(jī)焦點(diǎn)在ψ=0°時(shí)處于后限位置,ψ=90°時(shí)處于前限位置,見(jiàn)圖10、圖11.

全機(jī)焦點(diǎn)位置變化主要受旋翼、平尾焦點(diǎn)位置及升力線斜率影響,可用XfCLα表征翼面對(duì)焦點(diǎn)位置的拉動(dòng)能力,XfCLα越大,對(duì)焦點(diǎn)后拉能力越強(qiáng).圖12表示ψ=0°~90°時(shí)旋翼、平尾對(duì)焦點(diǎn)的拉動(dòng)能力變化.平尾流場(chǎng)受鴨翼、旋翼、機(jī)身流場(chǎng)影響,旋翼處于不同方位角時(shí),對(duì)平尾的影響存在差別,從而使平尾升力線斜率產(chǎn)生微小變化.

圖10 全機(jī)焦點(diǎn)位置變化

圖11 ψ=0°與90°時(shí)旋翼與全機(jī)焦點(diǎn)位置

圖12 旋翼、平尾焦點(diǎn)位置與升力線斜率變化

ψ=0°時(shí),盡管旋翼焦點(diǎn)處于前限位置,但升力線斜率最低,故旋翼將全機(jī)焦點(diǎn)向前拉動(dòng)作用仍最小.此時(shí)平尾升力線斜率最大,將全機(jī)焦點(diǎn)拉至最后.ψ=90°時(shí),旋翼升力線斜率最大,將焦點(diǎn)向前拉動(dòng)作用最強(qiáng),平尾升力線斜率最小,將焦點(diǎn)向后拉動(dòng)作用最弱,此時(shí)焦點(diǎn)處于前限位置.可見(jiàn)旋翼、平尾升力線斜率變化是全機(jī)焦點(diǎn)位置變化的主要原因.

3 結(jié)束語(yǔ)

1)鴨式旋翼/機(jī)翼飛行器轉(zhuǎn)換過(guò)程末段,全機(jī)氣動(dòng)力及氣動(dòng)力矩隨旋翼方位角均發(fā)生周期性變化,ψ=90°時(shí)(旋翼垂直于機(jī)體軸線)升力最小、阻力最大、抬頭力矩最小,ψ=45°時(shí)側(cè)力、滾轉(zhuǎn)力矩、偏航力矩最大(小).其中升力系數(shù)變化幅度10.7%,阻力系數(shù)變化幅度3.7%.

2)旋翼處于鴨翼、機(jī)身誘導(dǎo)的前后不對(duì)稱流場(chǎng)與旋翼轉(zhuǎn)至不同方位角時(shí)對(duì)機(jī)身、尾翼的左右不對(duì)稱干擾導(dǎo)致了全機(jī)氣動(dòng)力與氣動(dòng)力矩周期性變化.

3)轉(zhuǎn)換過(guò)程末段,旋翼焦點(diǎn)位置變化、旋翼升力線斜率變化、平尾升力線斜率變化導(dǎo)致全機(jī)焦點(diǎn)位置變化幅度較大,移動(dòng)量為0.6m,全機(jī)縱向靜穩(wěn)定性變化明顯.

References)

[1] Mitchell CA,Vogel B J.The canard rotorwing(CRW)aircraft:a new way to fly[R].AIAA 2003-2517,2003

[2] Kim C,Lee J.Hub drag-reduction strategies of canard rotor/wing unmanned aerial vehicles[J].Journal of Aircraft,2007,44(4):1388-1391

[3] Lee JW,Jeon K S,Kim M,et al.Rotor design for the performance optimization of canard rotor/wing aircraft[J].Lecture Notes in Computer Science,2006,3984:932 -941

[4]曹旭,詹浩,鄧陽(yáng)平.尖部噴氣驅(qū)動(dòng)旋轉(zhuǎn)機(jī)翼地面試驗(yàn)裝置的研制[J].實(shí)驗(yàn)流體力學(xué),2007,21(2):89 -92 Cao Xu,Zhan Hao,Deng Yangping.Developmentof experimental equipment for reaction driving rotor wing[J].Journal of Experiments in Fluid Mechanics,2007,21(2):89 -92(in Chinese)

[5]詹浩,鄧陽(yáng)平,高正紅.橢圓翼型低速氣動(dòng)特性研究[J].航空計(jì)算技術(shù),2008,38(3):25 -27 Zhan Hao,Deng Yangping,Gao Zhenghong.Investigation on aerodynamics performance of elliptic airfoil at low speed[J].Aeronautical Computing Technique,2008,38(3):25 - 27(in Chinese)

[6] Pandya SA,AftosmisM J.Computation of external aerodynamics for a canard rotor/wing aircraft[R].AIAA 2001-0997,2001

[7]陸志良.空氣動(dòng)力學(xué)[M].北京:北京航空航天大學(xué)出版社,2009:256-260 Lu Zhiliang.Aerodynamics[M].Beijing:Beijing University of Aeronautics and Astronautics Press,2009:256 - 260(in Chinese)

[8]朱自強(qiáng),吳子牛,李津,等.應(yīng)用計(jì)算流體力學(xué)[M].北京:北京航空航天出版社,1998:136-137 Zhu Ziqiang,Wu Ziniu,Li Jin,et al.Application of computational fluid dynamics[M].Beijing:Beijing University of Aeronautics and Astronautics Press,1998:136 -137(in Chinese)

(編 輯:李 晶)

Aerodynamic characteristic of canard rotor/wing aircraft in conversion

Li Yibo Ma Dongli

(School of Aeronautic Science and Engineering,Beijing University of Aeronautics and Astronautics,Beijing 100191,China)

Niu Lingyu

(Beijing System Engineering Institute,Beijing 100101,China)

The aerodynamic characteristics and mechanism of canard rotor/wing(CRW)aircraft during conversion from rotary to fixed-wing flightwas numerically investigated.The variation of forces,moments and aerodynamic center with respect to rotor position are presented,the amplitude of lift,drag and aerodynamic center for this configuration in conversion can reach 10.7%,3.7%and 0.6m separately.The investigation shows that the cause of forces and moments variation is the asymmetry flow field in rotary plane and asymmetry interference between rotor and fuselage,the motion of aerodynamic center can be explained by the motion of rotor aerodynamic center and the variation of lift curve slope of rotor and horizontal tail.

canard rotor/wing(CRW)aircraft;aerodynamic characteristic;flow interference;flight mode conversion

V 275+.3

A

1001-5965(2011)03-0311-05

2010-07-14

武器裝備預(yù)研重點(diǎn)基金資助項(xiàng)目(9140A25010208HK0101)

李毅波(1983-),男,湖南婁底人,博士生,liyibo1983@tom.com.

猜你喜歡
平尾
基于流場(chǎng)顯示技術(shù)的冰污染平尾失速試飛技術(shù)
直升機(jī)平尾電磁散射特性研究
雙層平尾對(duì)旋翼/平尾干擾的抑制機(jī)理研究
民用飛機(jī)平尾載荷的不確定性及全局靈敏度分析
直升機(jī)平尾載荷標(biāo)定工裝夾具設(shè)計(jì)
某型飛機(jī)全動(dòng)平尾安裝結(jié)構(gòu)優(yōu)化設(shè)計(jì)
教練機(jī)(2018年3期)2018-11-29 06:52:24
平尾損傷計(jì)算中特征載荷的算法研究
日本27歲囚犯越獄只因“感到無(wú)聊”
看世界(2018年10期)2018-05-25 04:53:40
全動(dòng)式水平尾翼
大飛機(jī)(2018年1期)2018-05-14 15:59:08
某型飛機(jī)平尾活動(dòng)間隙的定量檢查與控制
主站蜘蛛池模板: 日韩精品免费一线在线观看| 亚洲成人精品久久| 操国产美女| 亚洲一区二区约美女探花| 深爱婷婷激情网| 人妻无码AⅤ中文字| 黄色网址手机国内免费在线观看| 亚洲精品色AV无码看| 91无码国产视频| 91 九色视频丝袜| 久久 午夜福利 张柏芝| 久久性妇女精品免费| 亚洲最大综合网| 激情乱人伦| 最新亚洲人成无码网站欣赏网 | m男亚洲一区中文字幕| 国产视频欧美| 天天躁日日躁狠狠躁中文字幕| 欧美日韩成人在线观看| 国产精品免费露脸视频| 日韩成人午夜| 国产jizzjizz视频| 国产在线视频导航| 精品国产网| 欧美午夜小视频| 国产91在线|日本| 精品国产欧美精品v| 五月天香蕉视频国产亚| 激情五月婷婷综合网| 亚洲综合色婷婷| 美女免费黄网站| 国产精品无码AV片在线观看播放| 亚洲无码视频图片| 91无码国产视频| 在线欧美a| 欧美日韩另类国产| 毛片在线播放网址| 草草影院国产第一页| 国内精品手机在线观看视频| 成人国产一区二区三区| 青青国产成人免费精品视频| 少妇精品久久久一区二区三区| 免费国产好深啊好涨好硬视频| 日韩欧美视频第一区在线观看| yjizz视频最新网站在线| 真人免费一级毛片一区二区| 色一情一乱一伦一区二区三区小说| 九色视频线上播放| 日本一区二区三区精品国产| 在线国产毛片| 国产福利拍拍拍| 99国产精品一区二区| 亚洲黄色高清| 日本AⅤ精品一区二区三区日| 这里只有精品在线播放| 91精品小视频| 国产丰满大乳无码免费播放| 福利在线一区| 麻豆精品国产自产在线| 无码一区18禁| 婷婷色丁香综合激情| 在线观看亚洲成人| 米奇精品一区二区三区| 91网在线| 色老头综合网| 亚洲二区视频| 亚洲第一成人在线| 永久免费精品视频| 综合色区亚洲熟妇在线| 波多野结衣一二三| 国产情精品嫩草影院88av| 日韩高清一区 | 国产肉感大码AV无码| 国产欧美日本在线观看| 日韩久久精品无码aV| 精品亚洲国产成人AV| 五月天福利视频| 亚洲精品无码高潮喷水A| h网址在线观看| 国产原创第一页在线观看| 欧美啪啪精品| 日韩毛片在线视频|