• <tr id="yyy80"></tr>
  • <sup id="yyy80"></sup>
  • <tfoot id="yyy80"><noscript id="yyy80"></noscript></tfoot>
  • 99热精品在线国产_美女午夜性视频免费_国产精品国产高清国产av_av欧美777_自拍偷自拍亚洲精品老妇_亚洲熟女精品中文字幕_www日本黄色视频网_国产精品野战在线观看 ?

    超臨界自然層流機(jī)翼設(shè)計(jì)及基于TSP技術(shù)的邊界層轉(zhuǎn)捩風(fēng)洞試驗(yàn)

    2019-04-22 10:44:34張彥軍段卓毅雷武濤白俊強(qiáng)徐家寬
    航空學(xué)報(bào) 2019年4期
    關(guān)鍵詞:層流風(fēng)洞試驗(yàn)馬赫數(shù)

    張彥軍,段卓毅,雷武濤,白俊強(qiáng),徐家寬

    1. 航空工業(yè)第一飛機(jī)設(shè)計(jì)研究院,西安 710089 2. 西北工業(yè)大學(xué) 航空學(xué)院,西安 710072

    隨著環(huán)境保護(hù)形勢(shì)的日益嚴(yán)峻,國(guó)際航空運(yùn)輸協(xié)會(huì)提出了航空工業(yè)減少排放物和降低噪聲的新要求。而能夠?qū)崿F(xiàn)這一目標(biāo)的相關(guān)關(guān)鍵技術(shù)當(dāng)中,氣動(dòng)減阻技術(shù)成為空氣動(dòng)力學(xué)設(shè)計(jì)者重點(diǎn)研究的對(duì)象。隨著航空工業(yè)設(shè)計(jì)技術(shù)和制造工藝的突飛猛進(jìn),層流流動(dòng)設(shè)計(jì)逐步成為可能。對(duì)于民用客機(jī)而言,層流機(jī)翼設(shè)計(jì)技術(shù)(機(jī)翼、垂尾、平尾)和層流短艙設(shè)計(jì)技術(shù)可以降低摩擦阻力30%左右,即提高巡航效率15%左右。氣動(dòng)收益非常明顯,提升氣動(dòng)性能的同時(shí)減少燃油消耗和污染物排放以及降低噪聲。

    在20世紀(jì)80年代,美國(guó)和歐盟在自然層流機(jī)翼和混合層流機(jī)翼方面進(jìn)行了大量的理論研究、風(fēng)洞試驗(yàn)驗(yàn)證和飛行試驗(yàn)驗(yàn)證[1-3]。比較有代表性的是:Boeing公司在B757上的層流減阻研究,歐盟在SAAB 2000飛機(jī)的機(jī)翼上進(jìn)行的層流控制飛行試驗(yàn)驗(yàn)證和在A320垂尾上進(jìn)行的飛行試驗(yàn)驗(yàn)證。Falcon50飛機(jī)飛行試驗(yàn)是進(jìn)行混合層流設(shè)計(jì)項(xiàng)目研究的一部分,目的是在未來(lái)商用飛機(jī)的飛行馬赫數(shù)、雷諾數(shù)和后掠角范圍內(nèi)研究層流控制的可行性。在國(guó)外進(jìn)行的層流設(shè)計(jì)研究和案例當(dāng)中,德國(guó)宇航研究院(DLR)早在20世紀(jì)80年代就開(kāi)始進(jìn)行了前掠翼布局自然層流飛機(jī)(Forward Swept Wing-Natural Laminar Flow:FSW-NLF)的研究[4]。隨著現(xiàn)代工業(yè)的發(fā)展,滿足層流流動(dòng)對(duì)表面波紋度、光潔度等加工要求的機(jī)體可以實(shí)現(xiàn),自然層流設(shè)計(jì)也終于應(yīng)用到了工程上面。Honda Jet輕型公務(wù)機(jī)采用自然層流機(jī)身頭部和自然層流機(jī)翼[5-6],于2003年成功首飛并達(dá)到預(yù)期的設(shè)計(jì)目標(biāo)和要求。2006年,意大利Piaggio Aero Industries公司的研究人員與意大利宇航研究院合作提出了一種跨聲速自然層流超臨界機(jī)翼設(shè)計(jì)方法,并且基于此制作了相應(yīng)的風(fēng)洞模型,命名為UW-5006自然層流機(jī)翼,進(jìn)行了風(fēng)洞試驗(yàn)[7-8]。2010年,波音公司和美國(guó)國(guó)家航空航天局(NASA)進(jìn)行了高雷諾數(shù)后掠機(jī)翼自由轉(zhuǎn)捩風(fēng)洞試驗(yàn),分別在美國(guó)NTF(National Transonic Facility)和歐洲ETW(European Transonic Windtunnel)風(fēng)洞進(jìn)行溫敏涂層轉(zhuǎn)捩測(cè)量試驗(yàn),研究不同雷諾數(shù)下TS (Tollmien-Schlichting)波和橫流(Crossflow)波主導(dǎo)的轉(zhuǎn)捩,并研究了模型表面加工粗糙度對(duì)橫流駐波誘導(dǎo)轉(zhuǎn)捩的影響[9]。

    中國(guó)對(duì)超臨界自然層流設(shè)計(jì)的研究仍然處于比較初步的階段。西北工業(yè)大學(xué)的喬志德[10]研究了自然層流超臨界翼型的設(shè)計(jì)方法,解決了維持層流所需的有一定順壓梯度壓力分布形態(tài)和無(wú)激波超臨界翼型的屋頂狀壓力分布要求的矛盾,為自然層流機(jī)翼的設(shè)計(jì)奠定了基礎(chǔ)。北京航空航天大學(xué)的額日其太等[11]針對(duì)層流控制在飛機(jī)減阻、外表面紅外隱身方法以及抑制氣動(dòng)熱的生成等方面進(jìn)行了研究,分析和試驗(yàn)結(jié)果證明:前緣吸氣具有很好的層流控制效果。孫智偉[12]、黃江濤[13-14]等進(jìn)行了超臨界翼型和機(jī)翼的優(yōu)化設(shè)計(jì)研究,而針對(duì)超臨界自然層流翼型,西北工業(yè)大學(xué)的喬志德等[15-16]進(jìn)行了較為詳細(xì)的設(shè)計(jì)思想、設(shè)計(jì)方法和風(fēng)洞試驗(yàn)研究,清華大學(xué)的張宇飛等[17]針對(duì)超臨界自然層流翼型和機(jī)翼進(jìn)行了優(yōu)化設(shè)計(jì)策略的研究。西北工業(yè)大學(xué)的韓忠華等[18]采用代理模型,對(duì)自然層流機(jī)翼進(jìn)行優(yōu)化設(shè)計(jì)研究。

    在高精度邊界層轉(zhuǎn)捩預(yù)測(cè)方法方面,近些年來(lái)得到了長(zhǎng)足的發(fā)展。其中,Langtry等[19-22]提出的基于經(jīng)驗(yàn)關(guān)系式的輸運(yùn)方程轉(zhuǎn)捩模式在機(jī)械流動(dòng)和航空流動(dòng)中應(yīng)用廣泛。另一種基于穩(wěn)定性理論分析的半經(jīng)驗(yàn)轉(zhuǎn)捩預(yù)測(cè)方法是eN方法,該方法主要使用線性穩(wěn)定性理論,描述小擾動(dòng)行波——TS波的振幅沿邊界層流向的線性放大階段,并根據(jù)經(jīng)驗(yàn)選定判定轉(zhuǎn)捩發(fā)生的方法因子臨界N值,從而預(yù)測(cè)低湍流度下的各類(lèi)擾動(dòng)波主導(dǎo)的轉(zhuǎn)捩現(xiàn)象。eN方法最早在20世紀(jì)中期由Smith[23]和van Ingen[24]等發(fā)展而來(lái),隨后Gleyzes[25]和Drela[26]等進(jìn)一步提出了近似包絡(luò)方法。近似包絡(luò)方法通過(guò)采用線性穩(wěn)定性方法分析得到F-S(Falkner-Skan)速度型及其對(duì)應(yīng)的擾動(dòng)放大因子n與動(dòng)量厚度雷諾數(shù)的曲線,并將其用數(shù)學(xué)描述,得到對(duì)應(yīng)不同速度型的擾動(dòng)放大因子包絡(luò)線,將其作為轉(zhuǎn)捩判斷的數(shù)據(jù)庫(kù)。2013年,Coder和Maughmer[27-28]基于前述數(shù)據(jù)庫(kù)里的n因子與形狀因子和動(dòng)量損失厚度雷諾數(shù)的關(guān)系,構(gòu)造出了流向擾動(dòng)放大因子的輸運(yùn)方程,與Menterk-ωSST(Shear Stress Transport)湍流模式[29]耦合形成基于線性穩(wěn)定性理論的湍流轉(zhuǎn)捩模式。該方法所有變量均能夠當(dāng)?shù)厍蠼?,與現(xiàn)代CFD大規(guī)模并行求解兼容,且具有高精度的穩(wěn)定性分析基礎(chǔ)。徐家寬和白俊強(qiáng)[30]使用標(biāo)量輸運(yùn)方程的形式實(shí)現(xiàn)了包絡(luò)近似方法中放大因子的當(dāng)?shù)鼗蠼?,?shí)現(xiàn)了自然轉(zhuǎn)捩和分離泡轉(zhuǎn)捩的建模。

    在邊界層轉(zhuǎn)捩試驗(yàn)研究方面,中航工業(yè)氣動(dòng)院的尚金奎等[31]對(duì)溫度敏感材料涂層(Temperature Sensitive Paint, TSP)轉(zhuǎn)捩預(yù)測(cè)試驗(yàn)技術(shù)進(jìn)行了研究,采用TSP技術(shù)對(duì)某民機(jī)半模進(jìn)行試驗(yàn),預(yù)測(cè)轉(zhuǎn)捩位置,并通過(guò)與紅外試驗(yàn)技術(shù)預(yù)測(cè)結(jié)果進(jìn)行對(duì)比,驗(yàn)證了TSP試驗(yàn)方法的精度;北京大學(xué)的朱一丁等[32]采用瑞利散射流動(dòng)顯示、高頻動(dòng)態(tài)壓力傳感器以及粒子圖像測(cè)速等方法,在北京大學(xué)高超馬赫風(fēng)洞中開(kāi)展試驗(yàn),對(duì)高超聲速邊界層轉(zhuǎn)捩及湍流產(chǎn)生機(jī)理進(jìn)行了研究。

    當(dāng)前,國(guó)內(nèi)針對(duì)跨聲速超臨界自然層流機(jī)翼在高雷諾數(shù)下的邊界層轉(zhuǎn)捩試驗(yàn)研究非常罕見(jiàn),面對(duì)未來(lái)綠色高效飛行器的設(shè)計(jì)需求,這一方面的研究急需進(jìn)行和完善。本文對(duì)超臨界自然層流翼型和機(jī)翼進(jìn)行了設(shè)計(jì),并應(yīng)用高精度轉(zhuǎn)捩預(yù)測(cè)方法進(jìn)行氣動(dòng)特性評(píng)估,隨后加工制造了高質(zhì)量的風(fēng)洞試驗(yàn)?zāi)P筒⑦M(jìn)行了精細(xì)的高雷諾數(shù)邊界層轉(zhuǎn)捩風(fēng)洞試驗(yàn)驗(yàn)證,與高精度轉(zhuǎn)捩數(shù)值模擬結(jié)果進(jìn)行對(duì)比分析,得到超臨界自然層流機(jī)翼的邊界層轉(zhuǎn)捩特性,預(yù)期對(duì)該類(lèi)型機(jī)翼的研究和發(fā)展起到一定的推動(dòng)作用。

    1 超臨界自然層流機(jī)翼

    1.1 超臨界自然層流機(jī)翼的設(shè)計(jì)

    自然層流機(jī)翼需要在合適的氣動(dòng)布局和設(shè)計(jì)約束下才能發(fā)揮最佳的減阻效果,如適當(dāng)?shù)娘w行雷諾數(shù)、較小的機(jī)翼前緣后掠角、機(jī)翼上最好不要安裝發(fā)動(dòng)機(jī)等,因此研究背景飛機(jī)最終選定為尾吊布局噴氣式飛機(jī),設(shè)計(jì)目標(biāo)為支線客機(jī)和公務(wù)機(jī),其氣動(dòng)布局三面圖,如圖1所示,機(jī)翼平面形狀和展向參數(shù)分布如圖2所示。機(jī)翼的具體參數(shù)見(jiàn)表1。

    需要指出的是全機(jī)升力系數(shù)為0.4,考慮配平損失后對(duì)翼身組合體構(gòu)型設(shè)計(jì)升力系數(shù)定為0.42。 設(shè)計(jì)馬赫數(shù)為0.75;前緣后掠角為17.5°,屬于小后掠機(jī)翼范疇;設(shè)計(jì)飛行雷諾數(shù)為1.8×107左右。但考慮到風(fēng)洞試驗(yàn)技術(shù)、研究經(jīng)費(fèi)等原因,可在風(fēng)洞試驗(yàn)中驗(yàn)證的最大雷諾數(shù)在1×107左右,根據(jù)相關(guān)文獻(xiàn)中大量的穩(wěn)定性分析和飛行試驗(yàn)數(shù)據(jù)(如圖3所示)可知,在該后掠角和雷諾數(shù)組合狀態(tài)下,TS波失穩(wěn)主導(dǎo)轉(zhuǎn)捩,尚未出現(xiàn)橫流不穩(wěn)定性轉(zhuǎn)捩。且試驗(yàn)風(fēng)洞為低湍流度風(fēng)洞,暫不考慮橫流行波失穩(wěn)。針對(duì)橫流駐波主導(dǎo)的失穩(wěn),其擾動(dòng)源主要是壁面粗糙度。自然層流機(jī)翼風(fēng)洞試驗(yàn)?zāi)P蛻?yīng)比普通測(cè)力測(cè)壓的試驗(yàn)?zāi)P途哂懈叩墓鉂嵍?,普通試?yàn)?zāi)P蜋C(jī)翼表面粗糙度為Ra=0.8 μm,自然層流機(jī)翼風(fēng)洞試驗(yàn)?zāi)P蜋C(jī)翼表面粗糙度應(yīng)達(dá)到Ra=0.4 μm,具有比較高的橫流駐波失穩(wěn)臨界雷諾數(shù),不容易發(fā)生該類(lèi)型的轉(zhuǎn)捩。因此,在進(jìn)行機(jī)翼設(shè)計(jì)時(shí)暫不考慮橫流不穩(wěn)定性轉(zhuǎn)捩,從而采用先進(jìn)行基本翼型設(shè)計(jì),再進(jìn)行三維機(jī)翼設(shè)計(jì)的策略。

    圖1 背景飛機(jī)三面圖Fig.1 Plane three-view layout

    圖2 機(jī)翼平面形狀和參數(shù)分布Fig.2 Plane shape and parameters distribution of wing

    表1 機(jī)翼形狀具體參數(shù)Table 1 Detailed parameters of wing

    圖3 飛行試驗(yàn)和穩(wěn)定性分析結(jié)果總結(jié)而來(lái)的主導(dǎo) 失穩(wěn)類(lèi)型與前緣后掠角、雷諾數(shù)之間的關(guān)系Fig.3 Relations among dominated instability mode, leading edge swept angle and Reynolds number from the results of flight test and stability analysis

    機(jī)翼飛行雷諾數(shù)較高,飛行馬赫數(shù)較高,翼面上具有60%~70%弦長(zhǎng)的維持自然層流所需的有一定順壓梯度壓力分布形態(tài)是不現(xiàn)實(shí)的,因此在基本翼型設(shè)計(jì)時(shí),下翼面壓力分布順壓范圍定在50%左右,上翼面45%左右,以弱激波結(jié)束上翼面壓力分布順壓形態(tài),如圖4所示(圖中Cp為壓力系數(shù),c為參考弦長(zhǎng),x為弦向坐標(biāo));翼型具有適度的后加載,有利于減小低頭力矩,保證翼型后部的厚度。

    整個(gè)機(jī)翼由4個(gè)翼型控制剖面進(jìn)行三維構(gòu)造,翼根、拐折、70%展長(zhǎng)位置和翼梢4個(gè)設(shè)計(jì)翼型。因機(jī)翼當(dāng)?shù)乩字Z數(shù)內(nèi)翼大,外翼小,內(nèi)翼剖面最大厚度位置相對(duì)基本翼型前移,外翼后移,內(nèi)翼上翼面順壓梯度相對(duì)弦長(zhǎng)范圍減小,外翼增加。內(nèi)翼相對(duì)厚度大,外翼相對(duì)厚度小,外翼相對(duì)內(nèi)翼幾何負(fù)扭轉(zhuǎn),與一般機(jī)翼設(shè)計(jì)規(guī)律一致。最終設(shè)計(jì)所得4個(gè)剖面翼型如圖5所示,y為垂直于弦向的坐標(biāo),η為展向位置與展長(zhǎng)的比值。最終所得展向相對(duì)厚度分布及扭轉(zhuǎn)角分布如圖6所示,其中T為厚度。

    圖4 基本翼型設(shè)計(jì)狀態(tài)壓力系數(shù)分布Fig.4 Distributions of base airfoil pressure coefficients on design point

    圖5 機(jī)翼翼根、拐折、η=0.7處和翼梢的翼型Fig.5 Airfoils at root, kink, η=0.7 and tip of wing

    圖6 展向相對(duì)厚度分布和扭轉(zhuǎn)角分布Fig.6 Distribution of relative thickness and twist angle in spanwise direction

    1.2 邊界層轉(zhuǎn)捩預(yù)測(cè)方法

    包絡(luò)法中求解擾動(dòng)放大因子時(shí),n可以定義為

    (1)

    式中:s0和s分別表示沿流向積分的起始點(diǎn)和當(dāng)前位置;n的值取決于當(dāng)?shù)剡吔鐚有螤钜蜃右约皠?dòng)量損失厚度雷諾數(shù)Reθ,如果當(dāng)?shù)剡吔鐚有螤钜蜃右约皠?dòng)量厚度能夠合理地進(jìn)行當(dāng)?shù)鼗?,流?chǎng)中任意一點(diǎn)擾動(dòng)因子的當(dāng)?shù)卦鲩L(zhǎng)就可以求出。Coder和Maughmer[27]通過(guò)分析邊界層相似性解,構(gòu)建了合理的計(jì)算當(dāng)?shù)匦螤钜蜃拥墓?,使用輸運(yùn)方程對(duì)放大因子進(jìn)行求解:

    (2)

    (3)

    式中:σf=1.0;Pγ和Eγ分別為產(chǎn)生源項(xiàng)和破壞源項(xiàng);間歇因子γ與Menterk-ω剪切應(yīng)力輸運(yùn)(SST)湍流模式的耦合方式和各源項(xiàng)的詳細(xì)計(jì)算公式見(jiàn)文獻(xiàn)[27]。關(guān)于該方法的可靠性校核驗(yàn)證見(jiàn)文獻(xiàn)[27-28],本文不再贅述。

    1.3 氣動(dòng)特性評(píng)估

    CFD求解過(guò)程中,采用格心格式有限體積法求解可壓縮Navier-Stokes方程,無(wú)黏通量通過(guò)Roe的通量差分分裂FDS(Flux Difference Splitting)格式求解,黏性通量采用中心差分格式進(jìn)行離散,時(shí)間推進(jìn)采用近似因子分解(Approximate Factorization)方法。使用多重網(wǎng)格和網(wǎng)格序列技術(shù)加速求解的收斂。程序通過(guò)基于MPI(Message Passing Interface)的分布式并行策略提高計(jì)算速度。

    數(shù)值模擬采用的網(wǎng)格分布如圖7所示。氣動(dòng)特性分析分別使用自由轉(zhuǎn)捩和全湍流計(jì)算,在馬赫數(shù)為0.75、雷諾數(shù)為1×107的工況下,圖8給出了計(jì)算所得升阻力系數(shù)極曲線,CL為升力系數(shù),CD為阻力系數(shù)。自然層流設(shè)計(jì)帶來(lái)的減阻效果非常明顯,翼身組合體的阻力減小在30 counts(1 count=1.0×10-4)左右,如果加上層流短艙等的貢獻(xiàn),減阻量將更加可觀,由此可見(jiàn)自然層流設(shè)計(jì)的巨大潛力和可觀收益。

    圖7 計(jì)算網(wǎng)格分布Fig.7 Distribution of computational mesh

    圖9給出了自由轉(zhuǎn)捩和全湍流工況下的翼身組合體構(gòu)型的阻力發(fā)散曲線,分別評(píng)估了定升力系數(shù)0.38、0.42和0.46 3個(gè)升力狀態(tài)。無(wú)論全湍流還是自由轉(zhuǎn)捩工況,Ma=0.77與設(shè)計(jì)點(diǎn)Ma=0.75阻力系數(shù)變化不超過(guò)20 counts,滿足馬赫數(shù)增加0.02、阻力系數(shù)增加不超過(guò)20 counts的要求,因而阻力發(fā)散特性良好。

    圖8 自由轉(zhuǎn)捩和全湍流工況下翼身組合體的 升阻力系數(shù)極曲線Fig.8 Curves between lift coefficient and drag coefficient at the free transition and fully turbulent condition

    圖9 自由轉(zhuǎn)捩和全湍流工況下翼身組合體的 阻力發(fā)散曲線Fig.9 Curves of drag divergence between lift coefficient and drag coefficient at the free transition and fully turbulent condition

    下面將詳細(xì)評(píng)估該機(jī)翼在設(shè)計(jì)點(diǎn)附近的壓力分布和邊界層轉(zhuǎn)捩特性。

    首先,在設(shè)計(jì)點(diǎn)附近進(jìn)行了馬赫數(shù)擾動(dòng)變化的壓力分布分析,定升力系數(shù)為0.42,Ma=0.74, 0.75,0.76,計(jì)算所得不同展向位置的壓力分布對(duì)比如圖10所示。由圖可知,馬赫數(shù)的小幅度變化對(duì)機(jī)翼幾乎所有展向位置的壓力分布均有較為明顯的影響,尤其是上表面,其主要是由于激波位置的前后移動(dòng)所致,馬赫數(shù)越大,激波位置越靠后,下表面則變化幅度非常小。

    然后,圖11給出了Ma=0.75,CL=0.38,0.42 ,0.46時(shí),計(jì)算所得不同展向位置的壓力分布對(duì)比。由圖可知,不同升力系數(shù)直接影響的是飛行迎角,與馬赫數(shù)變化產(chǎn)生的效應(yīng)類(lèi)似,不同升力系數(shù)下的壓力分布差異主要集中在上表面激波位置附近,升力系數(shù)越大,對(duì)應(yīng)迎角越大,激波位置后移,但是后移程度不及馬赫數(shù)變化產(chǎn)生的偏移量。同樣地,下表面壓力分布受影響很微弱。

    最后,在設(shè)計(jì)點(diǎn)(Re=1×107,Ma=0.75,CL=0.42)狀態(tài),機(jī)翼上下表面摩擦力系數(shù)Cf分布如圖12所示。分析云圖可知,除了機(jī)翼和機(jī)身結(jié)合部的干擾所致轉(zhuǎn)捩,其余部分均與壓力分布形態(tài)對(duì)應(yīng)良好。從內(nèi)翼段到外翼段,上下表面的轉(zhuǎn)捩位置均出現(xiàn)在逆壓梯度出現(xiàn)的壓力恢復(fù)區(qū)域,尤其是上表面在激波出現(xiàn)的位置附近,該狀態(tài)下并無(wú)激波誘導(dǎo)附面層分離泡轉(zhuǎn)捩出現(xiàn)。這也與設(shè)計(jì)的目標(biāo)一致:即在順壓梯度保證TS波的抑制發(fā)展,逆壓梯度區(qū)TS波則會(huì)快速增長(zhǎng),誘發(fā)轉(zhuǎn)捩,由此獲得層流設(shè)計(jì)。

    經(jīng)過(guò)精細(xì)的氣動(dòng)設(shè)計(jì)和高精度CFD驗(yàn)證之后,將設(shè)計(jì)構(gòu)型加工成風(fēng)洞試驗(yàn)?zāi)P?,進(jìn)行風(fēng)洞試驗(yàn)驗(yàn)證。關(guān)于該機(jī)翼其他工況下的轉(zhuǎn)捩特性分析和驗(yàn)證將在后續(xù)章節(jié)與風(fēng)洞試驗(yàn)結(jié)果一起進(jìn)行。

    圖10 Ma=0.74,0.75,0.76時(shí)機(jī)翼不同展向位置的壓力系數(shù)分布對(duì)比Fig.10 Comparison of pressure coefficient distribution at different spanwise sections of wing at Ma =0.74, 0.75, 0.76

    圖11 CL=0.38, 0.42, 0.46時(shí)機(jī)翼不同展向位置的壓力系數(shù)分布對(duì)比Fig.11 Comparison of pressure coefficient distribution at different spanwise sections of wing at CL=0.38, 0.42, 0.46

    圖12 設(shè)計(jì)點(diǎn)機(jī)翼表面摩擦力系數(shù)云圖(Re=1×107,Ma=0.75,CL=0.42)Fig.12 Contour of skin friction coefficient on the wing surface at design point (Re=1×107, Ma=0.75, CL=0.42)

    2 風(fēng)洞試驗(yàn)設(shè)施和測(cè)量技術(shù)

    2.1 荷蘭HST風(fēng)洞和模型加工情況

    風(fēng)洞試驗(yàn)在荷蘭阿姆斯特丹的DNW(German-Dutch Wind tunnels)風(fēng)洞群內(nèi)的HST(High Speed wind Tunnel)跨聲速風(fēng)洞進(jìn)行[34-35]。需要指出的是,HST風(fēng)洞自1950年前后建成運(yùn)營(yíng)以來(lái),進(jìn)行了各種各樣的民機(jī)和軍機(jī)的風(fēng)洞試驗(yàn),取得了非常高的試驗(yàn)數(shù)據(jù)精度和風(fēng)洞品質(zhì)。該風(fēng)洞的宏觀構(gòu)造圖如圖13所示,它是一款可變密度的回流式風(fēng)洞,風(fēng)洞滯止壓強(qiáng)范圍為20~390 kPa,風(fēng)洞試驗(yàn)段尺寸1.8 m×2.0 m,馬赫數(shù)覆蓋范圍為0.1~1.3,雷諾數(shù)上限可達(dá)1×107,流場(chǎng)品質(zhì)較高,尤其適合高雷諾數(shù)跨聲速風(fēng)洞試驗(yàn)。

    風(fēng)洞試驗(yàn)?zāi)P褪紫劝凑毡?中的詳細(xì)參數(shù)對(duì)設(shè)計(jì)構(gòu)型進(jìn)行了1∶10.4的三維等比例縮放,精密加工之后,接著對(duì)設(shè)計(jì)構(gòu)型使用結(jié)構(gòu)有限元進(jìn)行了強(qiáng)度和剛度校核,計(jì)算出安全系數(shù)滿足HST風(fēng)洞要求。機(jī)翼外露部分表面粗糙度為Ra=0.4 μm, 機(jī)身外露部分表面粗糙度為Ra=0.8 μm,其余部分粗糙度為Ra=1.6 μm。因?yàn)槭前肽T囼?yàn),所以需在機(jī)身與風(fēng)洞洞壁之間添加附面層隔板,外置天平進(jìn)行測(cè)力,翼身組合體試驗(yàn)構(gòu)型在風(fēng)洞試驗(yàn)段的安裝情況如圖14所示。

    圖13 DNW-HST跨聲速風(fēng)洞Fig.13 DNW-HST transonic wind tunnel

    表2 風(fēng)洞試驗(yàn)翼身組合體具體參數(shù)

    圖14 風(fēng)洞試驗(yàn)段中翼身組合體構(gòu)型的 半模試驗(yàn)構(gòu)型Fig.14 Overview of wing-body combination half-model in wind tunnel test section

    2.2 TSP技術(shù)和測(cè)量細(xì)節(jié)

    TSP技術(shù)主要利用光學(xué)技術(shù)實(shí)現(xiàn)風(fēng)洞模型表面溫度的測(cè)量[36]。具體操作為:首先將溫敏材料均勻涂于機(jī)翼表面,確保氣動(dòng)外形不受影響。然后打開(kāi)風(fēng)洞制冷裝置,對(duì)風(fēng)洞的氣流和風(fēng)洞試驗(yàn)機(jī)翼進(jìn)行冷卻。準(zhǔn)備進(jìn)行轉(zhuǎn)捩測(cè)量時(shí),關(guān)閉制冷裝置,此時(shí)吹入風(fēng)洞的氣流溫度相比于物面溫度較高,在機(jī)翼表面將會(huì)進(jìn)行較為強(qiáng)烈的熱傳遞現(xiàn)象。湍流邊界層的熱傳導(dǎo)效率較高,而層流邊界層則較低,因此會(huì)在機(jī)翼表面出現(xiàn)明顯的溫度差,通過(guò)光學(xué)技術(shù)對(duì)物面溫度進(jìn)行拍照識(shí)別,從而判定層流-湍流區(qū)域。

    TSP涂層厚度為150~200 μm,3個(gè)標(biāo)志帶,展向25%、55%和85%展長(zhǎng)處,每10%弦長(zhǎng)一個(gè)標(biāo)志點(diǎn)。此外,轉(zhuǎn)捩帶位于距離前緣7%弦長(zhǎng)處,如圖15 所示。測(cè)壓孔徑在機(jī)翼表面為?=0.2 mm,保證測(cè)壓孔軸線與當(dāng)?shù)匦兔娣ň€方向一致,偏差小于3′。測(cè)壓孔周?chē)鷽](méi)有毛刺、雜質(zhì)、倒角和凹凸不平。測(cè)壓管選取外徑1.0 mm,內(nèi)徑0.7 mm的不銹鋼管,使用前按照要求會(huì)進(jìn)行氣密性檢查。測(cè)壓孔布置在展向35%和49%展長(zhǎng)處,如圖16所示。每個(gè)剖面上下表面各12個(gè)測(cè)壓孔,監(jiān)測(cè)10%~80%弦長(zhǎng)區(qū)間內(nèi)的離散壓力分布。

    圖15 流向等間距標(biāo)志孔和固定轉(zhuǎn)捩帶布置Fig.15 Uniformly spaced marked holes in streamwise direction and distribution of fixed transition tripping dots

    圖16 兩個(gè)測(cè)壓剖面的位置Fig.16 Overview of two cross-sections with pressure taps

    3 風(fēng)洞試驗(yàn)結(jié)果和數(shù)值模擬分析

    風(fēng)洞試驗(yàn)全部試驗(yàn)工況涵蓋了:馬赫數(shù)Ma=0.70,0.75,0.78,0.80,雷諾數(shù)Re=6×106,8×106, 9×106,10×106,迎角為1°和2°。但是隨著風(fēng)洞試驗(yàn)的進(jìn)行,有很多工況下的轉(zhuǎn)捩測(cè)量由于氣流污染物、物面污染物等因素的影響導(dǎo)致測(cè)量不是非常穩(wěn)定。表3給出了最終風(fēng)洞試驗(yàn)結(jié)果中流場(chǎng)品質(zhì)和轉(zhuǎn)捩測(cè)量效果均處于高水平高質(zhì)量的試驗(yàn)工況。下文將結(jié)合風(fēng)洞試驗(yàn)結(jié)果和數(shù)值模擬結(jié)果,詳細(xì)探討和研究這些變化的參數(shù)對(duì)跨聲速自然層流機(jī)翼在設(shè)計(jì)點(diǎn)附近的邊界層轉(zhuǎn)捩特性的影響。

    表3 風(fēng)洞試驗(yàn)高質(zhì)量測(cè)量工況

    3.1 同狀態(tài)不同車(chē)次試驗(yàn)結(jié)果

    首先對(duì)風(fēng)洞試驗(yàn)結(jié)果的重復(fù)性試驗(yàn)精度進(jìn)行了驗(yàn)證,圖17展示了Ma=0.75、Re=6×106、迎角α=2°工況下,不同車(chē)次的轉(zhuǎn)捩分布。由內(nèi)外翼段的轉(zhuǎn)捩線分布可知,雖然后續(xù)車(chē)次物面被些許污染物污染,但是整體轉(zhuǎn)捩位置變化非常小,證實(shí)了測(cè)量試驗(yàn)的高精度和合理性。

    圖17 不同車(chē)次的TSP測(cè)量所得機(jī)翼 表面層流-湍流區(qū)域分布Fig.17 Distribution of laminar-turbulent region on wing surface measured by TSP technique at different test numbers

    3.2 變馬赫數(shù)分析

    雷諾數(shù)Re=8×106、迎角α=1°,4個(gè)試驗(yàn)馬赫數(shù)Ma=0.70,0.75,0.78,0.80工況下的TSP技術(shù)拍攝層流-湍流分布如圖18所示。由試驗(yàn)結(jié)果可知,設(shè)計(jì)雷諾數(shù)和迎角1°固定,馬赫數(shù)從0.70 逐漸增大到0.80,機(jī)翼下表面的層流范圍幾乎不變,而上表面的邊界層轉(zhuǎn)捩特性則會(huì)經(jīng)歷一個(gè)比較復(fù)雜的變化過(guò)程。馬赫數(shù)從0.70增大到0.75時(shí),上表面外翼段層流區(qū)縮短,內(nèi)翼段略有增長(zhǎng)。馬赫數(shù)繼續(xù)增大到0.78和0.80,上翼面的層流區(qū)都急劇增加,且馬赫數(shù)0.78和 0.80 工況下的上表面層流區(qū)范圍差異很小。

    為了分析其轉(zhuǎn)捩特性變化的原因,本文對(duì)試驗(yàn)構(gòu)型進(jìn)行了高精度轉(zhuǎn)捩數(shù)值模擬,如圖19所示。截取不同馬赫數(shù)下機(jī)翼上表面和下表面的典型展向位置的壓力分布,對(duì)比分析可知在機(jī)翼下表面,隨著馬赫數(shù)逐漸增大,壓力分布形態(tài)變化很微弱,逆壓梯度起始點(diǎn)略有前移,因此轉(zhuǎn)捩位置略微前移,轉(zhuǎn)捩形式均為T(mén)S波急劇失穩(wěn)產(chǎn)生的自然轉(zhuǎn)捩。在機(jī)翼上表面,馬赫數(shù)為0.70時(shí),機(jī)翼上表面還未有明顯激波出現(xiàn),轉(zhuǎn)捩發(fā)生于弱逆壓梯度的發(fā)展過(guò)程中,TS波逐漸失穩(wěn)形成自然轉(zhuǎn)捩;馬赫數(shù)為0.75時(shí),機(jī)翼上表面形成較為明顯的激波,轉(zhuǎn)捩也發(fā)生在較強(qiáng)逆壓梯度的激波形成區(qū)域,轉(zhuǎn)捩形式依然為自然轉(zhuǎn)捩;馬赫數(shù)為0.78時(shí),轉(zhuǎn)捩模式受機(jī)翼的影響,在內(nèi)翼段預(yù)測(cè)所得轉(zhuǎn)捩位置略微靠前,在展向中部和外部區(qū)域轉(zhuǎn)捩位置預(yù)測(cè)均與試驗(yàn)數(shù)據(jù)吻合較好。該工況下明顯的特征就是在馬赫數(shù)為0.78時(shí),機(jī)翼上表面順壓梯度區(qū)非常長(zhǎng),可以達(dá)到75%左右,因此上表面層流區(qū)顯著增長(zhǎng)。由此帶來(lái)的缺點(diǎn)是過(guò)長(zhǎng)的較強(qiáng)順壓會(huì)形成很強(qiáng)的壓力恢復(fù)導(dǎo)致強(qiáng)激波誘導(dǎo)附面層分離,因此該工況下機(jī)翼上表面出現(xiàn)了激波誘導(dǎo)附面層分離引起的轉(zhuǎn)捩。

    圖18 TSP測(cè)量所得不同馬赫數(shù)下 機(jī)翼表面層流-湍流分布Fig.18 Distribution of laminar-turbulent region on wing surface measured by TSP technique at different Mach numbers

    圖19 轉(zhuǎn)捩計(jì)算所得不同馬赫數(shù)下機(jī)翼 表面層流-湍流分布Fig.19 Distribution of laminar-turbulent region on wing surface measured by transition calculations at different Mach numbers

    3.3 變雷諾數(shù)分析

    馬赫數(shù)Ma=0.75、飛行迎角α=1°、3個(gè)試驗(yàn)雷諾數(shù)Re=6×106,8×106,10×106工況下的TSP技術(shù)拍攝層流-湍流分布如圖20所示。由圖可知,隨著雷諾數(shù)的增加,機(jī)翼上下表面轉(zhuǎn)捩位置隨著雷諾數(shù)的增加略微前移,但是變化幅度非常小。圖21給出了6×106、8×106和10×1063個(gè)雷諾數(shù)工況下機(jī)翼上下表面的轉(zhuǎn)捩數(shù)值模擬結(jié)果。由壓力分布形態(tài)可知,隨著雷諾數(shù)的增加,激波位置受到細(xì)微的影響,但幅度很小,機(jī)翼上表面轉(zhuǎn)捩位置整體隨著雷諾數(shù)增加略微前移,而下表面轉(zhuǎn)捩位置則與逆壓梯度起始點(diǎn)一致,并未隨著雷諾數(shù)增加而明顯變化。

    圖20 TSP測(cè)量所得不同雷諾數(shù)下機(jī)翼 表面層流-湍流分布(α=1°)Fig.20 Distribution of laminar-turbulent region on wing surface measured by TSP technique at different Reynolds numbers (α=1°)

    馬赫數(shù)Ma=0.75,飛行迎角α=2°,3個(gè)試驗(yàn)雷諾數(shù)Re=6×106,8×106,9×106工況下的TSP技術(shù)拍攝層流-湍流分布如圖22所示。由圖可知,隨著雷諾數(shù)的變化,機(jī)翼上下表面轉(zhuǎn)捩位置依然變化非常微弱。

    圖23給出了6×106、8×106和9×1063個(gè)雷諾數(shù)工況下機(jī)翼上下表面的轉(zhuǎn)捩數(shù)值模擬結(jié)果。截取內(nèi)中外翼3個(gè)典型展向位置的壓力分布,可見(jiàn)在該跨聲速狀態(tài)下,機(jī)翼上下表面的壓力分布和激波形態(tài)隨著雷諾數(shù)的增加變化微小,下表面的轉(zhuǎn)捩預(yù)測(cè)位置與試驗(yàn)數(shù)據(jù)吻合良好,均發(fā)生在逆壓梯度起始位置;在上表面,轉(zhuǎn)捩模式預(yù)測(cè)的內(nèi)翼段轉(zhuǎn)捩位置會(huì)隨著雷諾數(shù)的增加有些許前移,計(jì)算所得TS波在接近但還未抵達(dá)激波位置時(shí)就發(fā)生了轉(zhuǎn)捩,但試驗(yàn)結(jié)果顯示在內(nèi)翼段的層流區(qū)并未受到雷諾數(shù)的影響。而在機(jī)翼的中部和外部,數(shù)值模擬所得轉(zhuǎn)捩位置均與激波起始位置和試驗(yàn)測(cè)量數(shù)據(jù)保持一致。但是在2°迎角工況下,機(jī)翼上表面中外翼段的強(qiáng)激波將會(huì)誘導(dǎo)附面層分離觸發(fā)轉(zhuǎn)捩。這種激波誘導(dǎo)附面層分離現(xiàn)象在Re=6×106時(shí)最為明顯,隨著雷諾數(shù)的增加,分離泡現(xiàn)象會(huì)逐漸減弱。

    圖21 轉(zhuǎn)捩計(jì)算所得不同雷諾數(shù)下機(jī)翼 表面層流-湍流分布(α=1°)Fig.21 Distribution of laminar-turbulent region on wing surface measured by transition calculations at different Reynolds numbers (α=1°)

    圖22 TSP測(cè)量所得不同雷諾數(shù)下機(jī)翼 表面層流-湍流分布(α=2°)Fig.22 Distribution of laminar-turbulent region on wing surface measured by TSP technique at different Reynolds numbers (α=2°)

    圖23 轉(zhuǎn)捩計(jì)算所得不同雷諾數(shù)下機(jī)翼 表面層流-湍流分布(α=2°)Fig.23 Distribution of laminar-turbulent region on wing surface measured by transition calculations at different Reynolds numbers (α=2°)

    在試驗(yàn)中展向35%和49%展長(zhǎng)位置處的測(cè)壓探針?biāo)鶞y(cè)得壓力分布數(shù)據(jù)與CFD數(shù)值模擬的結(jié)果對(duì)比如圖24所示。由圖可知,無(wú)論轉(zhuǎn)捩預(yù)測(cè)還是試驗(yàn)數(shù)據(jù)均顯示該工況下壓力分布對(duì)雷諾數(shù)變化并不敏感。順壓梯度的保持使得層流區(qū)域得以維系,隨著強(qiáng)激波的出現(xiàn),會(huì)產(chǎn)生激波誘導(dǎo)附面層分離現(xiàn)象,從而觸發(fā)轉(zhuǎn)捩。

    3.4 變迎角分析

    馬赫數(shù)Ma=0.75,飛行雷諾數(shù)Re=8×106,2個(gè)飛行迎角α=1°、α=2°工況下的TSP技術(shù)拍攝層流-湍流分布如圖25所示。由圖可知,迎角由1°增加到2°,機(jī)翼上表面層流區(qū)域有所增長(zhǎng),而下表面則變化非常小。與之對(duì)應(yīng)的轉(zhuǎn)捩數(shù)值模擬結(jié)果如圖26所示,由圖可知迎角變化最明顯的影響就是激波位置和強(qiáng)度。迎角由1°增加到2°時(shí),機(jī)翼下表面壓力分布形態(tài)非常接近,轉(zhuǎn)捩位置也變化很小;機(jī)翼上表面激波位置后移,層流區(qū)也隨之增長(zhǎng),但中外翼段轉(zhuǎn)捩形態(tài)由自然轉(zhuǎn)捩變?yōu)閺?qiáng)激波誘導(dǎo)附面層分離泡轉(zhuǎn)捩。

    圖24 轉(zhuǎn)捩計(jì)算所得不同雷諾數(shù)下機(jī)翼表面壓力系數(shù)分布與試驗(yàn)數(shù)據(jù)的對(duì)比(α=2°,Ma=0.75)Fig.24 Comparison of pressure coefficient distribution on wing surface between transition calculations and measured data at different Reynolds numbers (α=2°,Ma=0.75)

    本文重點(diǎn)關(guān)注超臨界層流機(jī)翼的邊界層轉(zhuǎn)捩特性,對(duì)于阻力測(cè)力結(jié)果則只做定性說(shuō)明。風(fēng)洞半模測(cè)力試驗(yàn)誤差比較大,因此阻力系數(shù)的絕對(duì)數(shù)值分析價(jià)值不高。但是其所反應(yīng)的趨勢(shì),與CFD數(shù)值模擬一致:以固定轉(zhuǎn)捩為例(轉(zhuǎn)捩帶位于距離前緣7%弦長(zhǎng)位置處),馬赫數(shù)Ma=0.75, 迎角α=1°,雷諾數(shù)從6×106增加到10×106,翼身組合體的阻力以近似線性的關(guān)系在減小。分析其原因在于隨著雷諾數(shù)的增加,翼身組合體的當(dāng)量厚度(實(shí)際厚度+附面層厚度)減小,壓差阻力減小,其中摩擦阻力的減小量非常小。

    圖25 TSP測(cè)量所得不同迎角下機(jī)翼 表面層流-湍流分布Fig.25 Distribution of laminar-turbulent region on wing surface measured by TSP technique at different angles of attack

    圖26 轉(zhuǎn)捩計(jì)算所得不同迎角下機(jī)翼 表面層流-湍流分布Fig.26 Distribution of laminar-turbulent region on wing surface measured by transition calculations at different angles of attack

    4 結(jié) 論

    通過(guò)超臨界自然層流機(jī)翼的氣動(dòng)設(shè)計(jì)和風(fēng)洞試驗(yàn),以及轉(zhuǎn)捩數(shù)值模擬情況進(jìn)行詳細(xì)的介紹和對(duì)比分析,本文對(duì)于超臨界自然層流機(jī)翼的邊界層特性和設(shè)計(jì)理念得出以下幾點(diǎn)結(jié)論:

    1) 對(duì)于超臨界自然層流機(jī)翼,在一定范圍內(nèi)(雷諾數(shù)<1×107,升力系數(shù)<0.5),雷諾數(shù)的增加會(huì)使轉(zhuǎn)捩位置略微前移,但雷諾數(shù)并不是主導(dǎo)自然轉(zhuǎn)捩的關(guān)鍵因素;其他因素不變,在一定范圍內(nèi)(雷諾數(shù)<1×107,升力系數(shù)<0.5),雷諾數(shù)增加,飛行器總阻力減小。

    2) 對(duì)于超臨界自然層流機(jī)翼,在1×107雷諾數(shù)量級(jí),馬赫數(shù)和飛行迎角是主導(dǎo)邊界層轉(zhuǎn)捩的主導(dǎo)因素,因?yàn)檫@2個(gè)因素將直接影響壓力分布形態(tài),能夠改變順壓梯度區(qū)和激波位置以及強(qiáng)度,從而決定轉(zhuǎn)捩位置和轉(zhuǎn)捩類(lèi)型。順壓梯度區(qū)過(guò)長(zhǎng),則層流區(qū)域增加,但會(huì)形成強(qiáng)激波誘導(dǎo)附面層分離泡轉(zhuǎn)捩;順壓梯度區(qū)過(guò)短,則層流區(qū)域縮短,自然轉(zhuǎn)捩會(huì)在弱逆壓梯度區(qū)域形成。

    3) 基于擾動(dòng)放大因子的轉(zhuǎn)捩模式對(duì)本文構(gòu)型的轉(zhuǎn)捩預(yù)測(cè)基本與試驗(yàn)數(shù)據(jù)吻合良好,為超臨界自然翼型和機(jī)翼的設(shè)計(jì)提供了可靠的計(jì)算分析工具。

    4) 對(duì)于20°以下的后掠角,1×107量級(jí)的超臨界自然層流機(jī)翼的設(shè)計(jì)規(guī)律就是對(duì)于順壓梯度的設(shè)計(jì)必須和激波位置匹配,過(guò)強(qiáng)的順壓梯度會(huì)導(dǎo)致強(qiáng)激波誘導(dǎo)附面層分離,過(guò)弱的順壓梯度則很難維持充足的層流區(qū)域。對(duì)更大后掠、更高雷諾數(shù)的機(jī)翼需要考慮橫流不穩(wěn)定性轉(zhuǎn)捩的影響,需要直接進(jìn)行三維穩(wěn)定性分析、三維氣動(dòng)設(shè)計(jì)以及流動(dòng)控制。

    5) 在流場(chǎng)品質(zhì)不佳的環(huán)境下,機(jī)翼表面非常容易被污染,層流很難維持;隨著時(shí)間的推進(jìn),層流區(qū)域?qū)?huì)被污染物強(qiáng)制轉(zhuǎn)捩成湍流。因此,層流設(shè)計(jì)的發(fā)展需要突破工業(yè)加工精度和如何保持機(jī)翼表面光潔等技術(shù)難題。

    6) 該試驗(yàn)?zāi)P涂梢宰鳛檫吔鐚愚D(zhuǎn)捩研究者的驗(yàn)證模型,其包含跨聲速可壓縮邊界層高雷諾數(shù)工況的TS波和激波誘導(dǎo)附面層轉(zhuǎn)捩現(xiàn)象。

    猜你喜歡
    層流風(fēng)洞試驗(yàn)馬赫數(shù)
    高馬赫數(shù)激波作用下單模界面的Richtmyer-Meshkov不穩(wěn)定性數(shù)值模擬
    爆炸與沖擊(2024年7期)2024-11-01 00:00:00
    層流輥道電機(jī)IP56防護(hù)等級(jí)結(jié)構(gòu)設(shè)計(jì)
    摻氫對(duì)二甲醚層流燃燒特性的影響
    一維非等熵可壓縮微極流體的低馬赫數(shù)極限
    載荷分布對(duì)可控?cái)U(kuò)散葉型性能的影響
    層流切應(yīng)力誘導(dǎo)microRNA-101下調(diào)EZH2抑制血管新生
    低風(fēng)壓架空導(dǎo)線的風(fēng)洞試驗(yàn)
    電線電纜(2017年5期)2017-10-18 00:52:03
    滾轉(zhuǎn)機(jī)動(dòng)載荷減緩風(fēng)洞試驗(yàn)
    遮擋條件下超高層建筑風(fēng)洞試驗(yàn)研究
    重慶建筑(2014年12期)2014-07-24 14:00:32
    X80鋼層流冷卻溫度場(chǎng)的有限元模擬
    丁香六月天网| 香蕉国产在线看| 9热在线视频观看99| 肉色欧美久久久久久久蜜桃| 国产免费视频播放在线视频| 国产又色又爽无遮挡免| 欧美中文综合在线视频| 男女无遮挡免费网站观看| 日日撸夜夜添| 涩涩av久久男人的天堂| 99热国产这里只有精品6| 亚洲熟女毛片儿| 日韩av在线免费看完整版不卡| 中文乱码字字幕精品一区二区三区| 日韩熟女老妇一区二区性免费视频| 精品一区二区三区av网在线观看 | 国产伦理片在线播放av一区| 男女边吃奶边做爰视频| 搡老岳熟女国产| 十分钟在线观看高清视频www| 18禁观看日本| 国产精品一国产av| 建设人人有责人人尽责人人享有的| 超色免费av| 日韩一区二区三区影片| 国产黄色免费在线视频| 九色亚洲精品在线播放| 老鸭窝网址在线观看| 美女国产高潮福利片在线看| 精品一区在线观看国产| 免费久久久久久久精品成人欧美视频| 免费看不卡的av| av天堂久久9| 国产成人午夜福利电影在线观看| 久久久精品94久久精品| 久久狼人影院| 中文精品一卡2卡3卡4更新| 国产成人系列免费观看| 美女扒开内裤让男人捅视频| 日韩精品免费视频一区二区三区| 精品一区在线观看国产| 久久久精品区二区三区| 国产在线免费精品| 啦啦啦视频在线资源免费观看| 久久久国产精品麻豆| 中文字幕精品免费在线观看视频| 亚洲精品aⅴ在线观看| 水蜜桃什么品种好| 亚洲熟女精品中文字幕| 十八禁高潮呻吟视频| 久久久久久久大尺度免费视频| 在线观看免费午夜福利视频| 亚洲av在线观看美女高潮| 狂野欧美激情性xxxx| 亚洲人成网站在线观看播放| 啦啦啦中文免费视频观看日本| 国产成人精品在线电影| 中文字幕精品免费在线观看视频| 亚洲成av片中文字幕在线观看| 69精品国产乱码久久久| 美女福利国产在线| 色网站视频免费| 亚洲精品美女久久av网站| 国产日韩欧美视频二区| 一级毛片电影观看| 国产激情久久老熟女| 国产精品秋霞免费鲁丝片| 久久久久国产精品人妻一区二区| 免费日韩欧美在线观看| 欧美精品亚洲一区二区| 韩国高清视频一区二区三区| 老司机靠b影院| 国产精品久久久久久精品古装| 国产精品嫩草影院av在线观看| 一级a爱视频在线免费观看| 亚洲三区欧美一区| 亚洲三区欧美一区| 国产一区二区 视频在线| 精品人妻一区二区三区麻豆| av国产精品久久久久影院| 国产精品欧美亚洲77777| 搡老岳熟女国产| 欧美日韩亚洲国产一区二区在线观看 | 最新在线观看一区二区三区 | 亚洲精品国产区一区二| 一本大道久久a久久精品| av在线播放精品| 一级毛片 在线播放| 成人亚洲精品一区在线观看| 国产亚洲精品第一综合不卡| 精品久久久精品久久久| 午夜av观看不卡| 狠狠精品人妻久久久久久综合| 国产欧美日韩综合在线一区二区| av线在线观看网站| e午夜精品久久久久久久| 一级爰片在线观看| 成人国产av品久久久| 日日摸夜夜添夜夜爱| 久久鲁丝午夜福利片| 卡戴珊不雅视频在线播放| 日韩熟女老妇一区二区性免费视频| 91精品伊人久久大香线蕉| 夫妻性生交免费视频一级片| 丰满少妇做爰视频| 亚洲七黄色美女视频| 激情五月婷婷亚洲| 嫩草影视91久久| 菩萨蛮人人尽说江南好唐韦庄| 中文字幕亚洲精品专区| 国产欧美日韩综合在线一区二区| 如何舔出高潮| 日韩一卡2卡3卡4卡2021年| 午夜激情av网站| 久久久久久久久久久久大奶| 水蜜桃什么品种好| 老司机影院毛片| 咕卡用的链子| 免费观看人在逋| 丝袜脚勾引网站| 亚洲一区二区三区欧美精品| 国产精品免费视频内射| 91aial.com中文字幕在线观看| 飞空精品影院首页| 国产精品国产三级专区第一集| 成年人午夜在线观看视频| 男女边摸边吃奶| 麻豆精品久久久久久蜜桃| 一区在线观看完整版| 久久人人97超碰香蕉20202| 18禁裸乳无遮挡动漫免费视频| 欧美变态另类bdsm刘玥| 欧美另类一区| av视频免费观看在线观看| 两个人看的免费小视频| 岛国毛片在线播放| 国产精品熟女久久久久浪| 婷婷色麻豆天堂久久| 蜜桃国产av成人99| 国产成人免费观看mmmm| 精品人妻在线不人妻| 黄色怎么调成土黄色| 精品免费久久久久久久清纯 | 国产精品 国内视频| 天天躁夜夜躁狠狠躁躁| 狂野欧美激情性xxxx| 宅男免费午夜| 午夜日韩欧美国产| 亚洲av中文av极速乱| 18禁动态无遮挡网站| 亚洲av男天堂| 少妇被粗大的猛进出69影院| 亚洲,一卡二卡三卡| 成人毛片60女人毛片免费| 曰老女人黄片| 99精品久久久久人妻精品| 高清在线视频一区二区三区| 中国三级夫妇交换| 1024香蕉在线观看| 天天躁日日躁夜夜躁夜夜| 国产黄频视频在线观看| 免费高清在线观看日韩| 老鸭窝网址在线观看| 在线天堂中文资源库| 男人添女人高潮全过程视频| 丝袜美足系列| 啦啦啦在线观看免费高清www| 国产亚洲av片在线观看秒播厂| 美女脱内裤让男人舔精品视频| 国产黄色视频一区二区在线观看| 哪个播放器可以免费观看大片| 成年动漫av网址| 卡戴珊不雅视频在线播放| 国产精品国产三级专区第一集| 精品亚洲乱码少妇综合久久| 国产高清国产精品国产三级| 日本av免费视频播放| 青草久久国产| 大陆偷拍与自拍| 久久亚洲国产成人精品v| 亚洲伊人色综图| 亚洲三区欧美一区| 国产老妇伦熟女老妇高清| 性少妇av在线| www.av在线官网国产| 国产有黄有色有爽视频| 亚洲欧美日韩另类电影网站| 日韩一卡2卡3卡4卡2021年| 人妻一区二区av| 777久久人妻少妇嫩草av网站| 最黄视频免费看| 亚洲av电影在线进入| 国产毛片在线视频| 美女视频免费永久观看网站| 国产精品人妻久久久影院| av在线app专区| 黑人欧美特级aaaaaa片| 国产成人精品久久久久久| 丝袜在线中文字幕| av天堂久久9| 欧美黄色片欧美黄色片| 国产av国产精品国产| 人体艺术视频欧美日本| 午夜91福利影院| 国产av码专区亚洲av| 纯流量卡能插随身wifi吗| 国产片内射在线| 精品国产乱码久久久久久男人| 热re99久久精品国产66热6| 亚洲国产欧美在线一区| 欧美激情高清一区二区三区 | 一区二区三区乱码不卡18| 久久精品国产亚洲av高清一级| 看免费成人av毛片| 伦理电影大哥的女人| 婷婷色av中文字幕| 街头女战士在线观看网站| 悠悠久久av| 亚洲国产av影院在线观看| 欧美精品高潮呻吟av久久| 美国免费a级毛片| 久久久久精品久久久久真实原创| 熟女av电影| 国产黄色视频一区二区在线观看| 2018国产大陆天天弄谢| 亚洲av电影在线进入| 国产成人欧美| 国产精品女同一区二区软件| 男的添女的下面高潮视频| 国产日韩欧美在线精品| 欧美国产精品va在线观看不卡| 亚洲伊人久久精品综合| 人人妻人人澡人人爽人人夜夜| 纯流量卡能插随身wifi吗| 黄片播放在线免费| 一级,二级,三级黄色视频| 亚洲成人免费av在线播放| av在线播放精品| 看非洲黑人一级黄片| 国产精品.久久久| 女人被躁到高潮嗷嗷叫费观| 两性夫妻黄色片| 精品国产一区二区三区四区第35| xxx大片免费视频| 精品免费久久久久久久清纯 | 天堂8中文在线网| 高清黄色对白视频在线免费看| 国产精品麻豆人妻色哟哟久久| 亚洲综合精品二区| 视频在线观看一区二区三区| 国产午夜精品一二区理论片| 国产精品一区二区精品视频观看| 国产激情久久老熟女| 最新的欧美精品一区二区| 欧美黑人精品巨大| 国产乱来视频区| 丰满迷人的少妇在线观看| 秋霞伦理黄片| 久久女婷五月综合色啪小说| 国产成人91sexporn| 成人手机av| 中文字幕另类日韩欧美亚洲嫩草| 亚洲欧美色中文字幕在线| 九色亚洲精品在线播放| 久久精品国产综合久久久| 韩国高清视频一区二区三区| 一边摸一边做爽爽视频免费| 午夜激情av网站| av女优亚洲男人天堂| 国产精品嫩草影院av在线观看| 人妻人人澡人人爽人人| 午夜激情久久久久久久| 亚洲人成电影观看| 性高湖久久久久久久久免费观看| 侵犯人妻中文字幕一二三四区| 免费不卡黄色视频| 国产精品久久久久久人妻精品电影 | 狠狠婷婷综合久久久久久88av| 免费观看人在逋| 操出白浆在线播放| 国产99久久九九免费精品| 免费黄色在线免费观看| 久久国产精品男人的天堂亚洲| 伦理电影大哥的女人| 欧美在线一区亚洲| 精品久久久久久电影网| av国产精品久久久久影院| 老司机深夜福利视频在线观看 | 亚洲国产欧美日韩在线播放| 久久人人97超碰香蕉20202| 热99久久久久精品小说推荐| 精品人妻熟女毛片av久久网站| 日韩一本色道免费dvd| 亚洲成国产人片在线观看| 男男h啪啪无遮挡| 青春草国产在线视频| 男女无遮挡免费网站观看| 香蕉丝袜av| 午夜福利免费观看在线| 亚洲精品日韩在线中文字幕| 制服诱惑二区| 观看美女的网站| 国产成人免费无遮挡视频| 少妇 在线观看| 欧美av亚洲av综合av国产av | 天堂8中文在线网| 国产毛片在线视频| 在线看a的网站| 精品国产一区二区久久| 青青草视频在线视频观看| 精品卡一卡二卡四卡免费| 99久国产av精品国产电影| 久久精品国产a三级三级三级| 一个人免费看片子| 在线观看一区二区三区激情| 午夜福利在线免费观看网站| 日韩一区二区视频免费看| 热99久久久久精品小说推荐| 欧美最新免费一区二区三区| 日本av免费视频播放| 中文字幕人妻熟女乱码| 精品一区在线观看国产| 不卡av一区二区三区| 97在线人人人人妻| av有码第一页| 免费女性裸体啪啪无遮挡网站| 又大又黄又爽视频免费| 久久97久久精品| 亚洲精品国产av蜜桃| 黄片小视频在线播放| 电影成人av| 成人18禁高潮啪啪吃奶动态图| 久热爱精品视频在线9| 精品国产一区二区三区四区第35| 性高湖久久久久久久久免费观看| av又黄又爽大尺度在线免费看| 精品久久久精品久久久| 少妇猛男粗大的猛烈进出视频| 9热在线视频观看99| 97在线人人人人妻| 黄色一级大片看看| 99久久精品国产亚洲精品| 久久精品久久久久久久性| 欧美xxⅹ黑人| 色播在线永久视频| 考比视频在线观看| 亚洲伊人久久精品综合| 丰满饥渴人妻一区二区三| 满18在线观看网站| 久久精品国产综合久久久| 国产精品国产av在线观看| 高清黄色对白视频在线免费看| 久久ye,这里只有精品| 欧美精品高潮呻吟av久久| 18禁裸乳无遮挡动漫免费视频| 老司机靠b影院| 亚洲欧美一区二区三区黑人| 一区二区三区四区激情视频| 国产免费又黄又爽又色| 日韩制服骚丝袜av| 侵犯人妻中文字幕一二三四区| 国产又色又爽无遮挡免| 日韩av在线免费看完整版不卡| 久久久久视频综合| 成人黄色视频免费在线看| svipshipincom国产片| 精品少妇一区二区三区视频日本电影 | 亚洲成人av在线免费| 精品一区二区三区四区五区乱码 | 日本色播在线视频| 黑人猛操日本美女一级片| 在现免费观看毛片| 少妇精品久久久久久久| 在线天堂最新版资源| 精品国产一区二区久久| e午夜精品久久久久久久| 好男人视频免费观看在线| 在线观看免费视频网站a站| 纵有疾风起免费观看全集完整版| www.熟女人妻精品国产| 中文字幕亚洲精品专区| 中文字幕另类日韩欧美亚洲嫩草| 日本猛色少妇xxxxx猛交久久| 欧美亚洲日本最大视频资源| 人体艺术视频欧美日本| 啦啦啦在线观看免费高清www| 最近2019中文字幕mv第一页| 国产一区有黄有色的免费视频| 久久久久久久久久久久大奶| 在线观看免费高清a一片| 18禁观看日本| 伦理电影大哥的女人| 成年女人毛片免费观看观看9 | 免费在线观看黄色视频的| 曰老女人黄片| 国产成人精品久久久久久| 亚洲精品国产色婷婷电影| 女人被躁到高潮嗷嗷叫费观| 亚洲国产欧美一区二区综合| 毛片一级片免费看久久久久| 欧美日本中文国产一区发布| 青青草视频在线视频观看| 女的被弄到高潮叫床怎么办| 中文字幕人妻熟女乱码| 婷婷成人精品国产| 亚洲精品,欧美精品| 热99国产精品久久久久久7| 天天添夜夜摸| 啦啦啦视频在线资源免费观看| 亚洲视频免费观看视频| svipshipincom国产片| 亚洲av综合色区一区| 美女主播在线视频| 久久精品久久久久久久性| a级毛片在线看网站| 亚洲精品国产av成人精品| 日韩熟女老妇一区二区性免费视频| 青青草视频在线视频观看| 国产精品久久久人人做人人爽| 精品午夜福利在线看| 日韩 亚洲 欧美在线| 欧美激情 高清一区二区三区| 中文精品一卡2卡3卡4更新| 97在线人人人人妻| 国产成人a∨麻豆精品| 午夜免费观看性视频| 亚洲成国产人片在线观看| 久久久久久人妻| 青草久久国产| 狂野欧美激情性bbbbbb| 精品国产露脸久久av麻豆| 中文字幕高清在线视频| 亚洲,欧美精品.| 久久久精品94久久精品| 久久狼人影院| av又黄又爽大尺度在线免费看| 亚洲图色成人| 黄片无遮挡物在线观看| 麻豆乱淫一区二区| 日韩精品免费视频一区二区三区| 成人黄色视频免费在线看| 亚洲三区欧美一区| 欧美人与善性xxx| √禁漫天堂资源中文www| 亚洲国产精品国产精品| 少妇人妻 视频| 久久久久久人人人人人| 大香蕉久久网| 国产精品人妻久久久影院| 亚洲第一av免费看| 各种免费的搞黄视频| 国产欧美日韩一区二区三区在线| 国产高清不卡午夜福利| 菩萨蛮人人尽说江南好唐韦庄| 欧美日韩亚洲高清精品| 中文欧美无线码| 国产老妇伦熟女老妇高清| 成年美女黄网站色视频大全免费| 色播在线永久视频| 大香蕉久久网| 啦啦啦在线免费观看视频4| 深夜精品福利| 亚洲情色 制服丝袜| 久热这里只有精品99| 又粗又硬又长又爽又黄的视频| 成年av动漫网址| 69精品国产乱码久久久| 中文字幕制服av| 亚洲精品aⅴ在线观看| 免费观看性生交大片5| 亚洲av成人不卡在线观看播放网 | av国产久精品久网站免费入址| 午夜精品国产一区二区电影| 岛国毛片在线播放| 七月丁香在线播放| 日韩一区二区三区影片| 亚洲视频免费观看视频| 国产一区二区三区av在线| 天天躁日日躁夜夜躁夜夜| 久久精品aⅴ一区二区三区四区| 国产精品无大码| 亚洲欧美一区二区三区久久| av在线app专区| 中文天堂在线官网| 亚洲情色 制服丝袜| 欧美日韩精品网址| 国产一卡二卡三卡精品 | 18在线观看网站| 一级毛片我不卡| 亚洲,欧美精品.| 日韩大片免费观看网站| 国产精品免费视频内射| 久久久久网色| 十八禁网站网址无遮挡| 亚洲精品国产一区二区精华液| 狠狠精品人妻久久久久久综合| 国产精品 国内视频| 日日摸夜夜添夜夜爱| 欧美97在线视频| 99香蕉大伊视频| 三上悠亚av全集在线观看| 日本一区二区免费在线视频| 99久久人妻综合| 国产成人91sexporn| 丁香六月天网| 久久这里只有精品19| 十八禁高潮呻吟视频| 中国三级夫妇交换| 少妇人妻精品综合一区二区| 九草在线视频观看| 超碰成人久久| 天天操日日干夜夜撸| 最近2019中文字幕mv第一页| 亚洲av欧美aⅴ国产| 蜜桃国产av成人99| 国产一区亚洲一区在线观看| 中文字幕av电影在线播放| 最近的中文字幕免费完整| 少妇被粗大猛烈的视频| 国产淫语在线视频| 黄色毛片三级朝国网站| 亚洲欧美精品自产自拍| 欧美亚洲日本最大视频资源| 亚洲国产欧美在线一区| 日韩,欧美,国产一区二区三区| 热99国产精品久久久久久7| 如日韩欧美国产精品一区二区三区| 午夜福利视频在线观看免费| 久久综合国产亚洲精品| 国产伦人伦偷精品视频| 国产一级毛片在线| 久久久国产一区二区| 在线天堂中文资源库| 国产老妇伦熟女老妇高清| 美女主播在线视频| 精品一品国产午夜福利视频| 亚洲,一卡二卡三卡| 交换朋友夫妻互换小说| 自线自在国产av| 欧美少妇被猛烈插入视频| 人人妻人人澡人人爽人人夜夜| 最近的中文字幕免费完整| 日韩制服骚丝袜av| 亚洲第一青青草原| 国产精品久久久久成人av| 欧美日韩成人在线一区二区| 无遮挡黄片免费观看| 欧美日韩视频高清一区二区三区二| 久久人妻熟女aⅴ| 欧美亚洲日本最大视频资源| av网站在线播放免费| 国产成人精品福利久久| a 毛片基地| 国产精品秋霞免费鲁丝片| 亚洲视频免费观看视频| 日韩免费高清中文字幕av| 91aial.com中文字幕在线观看| 免费少妇av软件| 国精品久久久久久国模美| 久久久久久久久久久久大奶| 亚洲av在线观看美女高潮| 制服人妻中文乱码| 国产熟女欧美一区二区| 欧美亚洲日本最大视频资源| 熟女少妇亚洲综合色aaa.| 亚洲av日韩精品久久久久久密 | 国产黄色视频一区二区在线观看| 日韩制服骚丝袜av| 人妻一区二区av| 老司机影院成人| 久久久久精品久久久久真实原创| 午夜免费男女啪啪视频观看| 久久精品国产综合久久久| 国产精品成人在线| 精品国产一区二区三区四区第35| 美女福利国产在线| 免费人妻精品一区二区三区视频| 国产日韩欧美在线精品| 国产免费视频播放在线视频| 这个男人来自地球电影免费观看 | 成人漫画全彩无遮挡| 你懂的网址亚洲精品在线观看| 999精品在线视频| 欧美中文综合在线视频| 色网站视频免费| 超碰97精品在线观看| 美女扒开内裤让男人捅视频| 国产男女超爽视频在线观看| 热re99久久精品国产66热6| 欧美日韩亚洲国产一区二区在线观看 | 国产精品国产av在线观看| 最近手机中文字幕大全| 中文字幕人妻熟女乱码| 久久久欧美国产精品| 欧美 亚洲 国产 日韩一| 精品国产乱码久久久久久小说| 日韩欧美一区视频在线观看| 亚洲婷婷狠狠爱综合网| 国产精品无大码| 午夜免费观看性视频| 国产在视频线精品| 少妇被粗大猛烈的视频| 亚洲精华国产精华液的使用体验| 在现免费观看毛片| 亚洲国产av影院在线观看| 男女床上黄色一级片免费看| 欧美日韩av久久| 女人精品久久久久毛片| 啦啦啦 在线观看视频| 国产亚洲av高清不卡| 777米奇影视久久| 欧美最新免费一区二区三区| 999精品在线视频| 激情五月婷婷亚洲| 亚洲欧美精品自产自拍|