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

    定后掠角密切錐乘波體的生成和設(shè)計(jì)方法

    2016-11-20 01:51:06段焰輝范召林吳文華
    航空學(xué)報(bào) 2016年10期
    關(guān)鍵詞:后掠角激波前緣

    段焰輝, 范召林, 吳文華

    中國(guó)空氣動(dòng)力研究與發(fā)展中心 計(jì)算空氣動(dòng)力研究所, 綿陽(yáng) 621000

    定后掠角密切錐乘波體的生成和設(shè)計(jì)方法

    段焰輝, 范召林, 吳文華*

    中國(guó)空氣動(dòng)力研究與發(fā)展中心 計(jì)算空氣動(dòng)力研究所, 綿陽(yáng) 621000

    對(duì)定后掠角密切錐乘波體(OCWRCAS)的生成方法和考慮黏性的設(shè)計(jì)方法進(jìn)行了研究。定后掠角乘波體的前緣具有特定的后掠角,能夠在上表面產(chǎn)生穩(wěn)定分離渦從而改善乘波體的氣動(dòng)性能。本文首先在傳統(tǒng)密切錐乘波體生成方法的基礎(chǔ)上給出了定后掠角密切錐乘波體的生成方法;從前緣后掠的幾何特征中提取了后掠角、激波角和前緣曲線程度等設(shè)計(jì)變量,并研究了設(shè)計(jì)變量的取值范圍;以遍歷設(shè)計(jì)空間的思路對(duì)兩類定后掠角密切錐乘波體進(jìn)行了設(shè)計(jì)分析,研究了升阻比、體積效率隨設(shè)計(jì)變量的變化規(guī)律,然后在設(shè)計(jì)空間內(nèi)進(jìn)行了多目標(biāo)尋優(yōu);最后使用計(jì)算流體力學(xué)方法對(duì)定后掠角乘波體的乘波特性和渦升力特性進(jìn)行了驗(yàn)證。結(jié)果表明,由本文生成方法得到的定后掠角密切錐乘波體具有明顯的乘波特性并且能夠在較高的升阻比時(shí)保證一定的體積效率;定后掠角前緣能夠在一定的迎角下在上表面產(chǎn)生穩(wěn)定的分離渦,產(chǎn)生渦升力。

    乘波體; 密切錐; 黏性; 后掠角; 渦升力

    乘波體經(jīng)過幾十年的發(fā)展,從早期的單一構(gòu)型逐漸發(fā)展為具有更多特點(diǎn)的復(fù)雜構(gòu)型,處理手段也愈加靈活。Nonweiler[1]在1959年提出了由二維流場(chǎng)構(gòu)造高超聲速飛行器的理論,并據(jù)此生成了“Λ”型乘波體,隨后Jones等[2]提出了基于軸對(duì)稱流場(chǎng)的乘波體生成理論,國(guó)內(nèi)也對(duì)錐導(dǎo)乘波體進(jìn)行了詳細(xì)研究[3-4],將軸對(duì)稱流場(chǎng)擴(kuò)展至更一般的三維流場(chǎng),就出現(xiàn)了源于帶迎角錐和橢圓錐體流場(chǎng)等非軸對(duì)稱流場(chǎng)的乘波體[5-6]。發(fā)動(dòng)機(jī)一體化設(shè)計(jì)需要下表面有較為均勻的流場(chǎng),源于楔型-錐型混合流場(chǎng)的乘波體[7]以及由反設(shè)計(jì)思路得到的密切錐或密切流場(chǎng)乘波體[8-10]具有“平坦”的下表面,符合上述設(shè)計(jì)條件。在諸多生成方法中,密切錐(或密切流場(chǎng))乘波體生成方法靈活,可以生成具有更多特性的乘波體,國(guó)外已將這種方法發(fā)展為成熟軟件,能夠生產(chǎn)具有復(fù)雜平面形狀的乘波體[11],國(guó)內(nèi)由賀旭照等[12-13]提出來的密切曲面錐乘波體也是基于這種反設(shè)計(jì)思路得到的。

    根據(jù)乘波體的原理,其氣動(dòng)性能主要來源于下表面,如果能有效利用乘波體的上表面將有可能改善乘波體的整體性能。傳統(tǒng)做法是將上表面整個(gè)或部分區(qū)域設(shè)計(jì)為膨脹面[14],若只有膨脹而無壓縮則必然降低乘波體的體積,這對(duì)體積效率本就很低的乘波體是非常不利的(若采用先壓縮后膨脹表面,雖然有可能提高體積,但壓縮區(qū)域必然增大乘波體阻力)。定后掠角乘波體(Waverider with Constant Angle of Sweepback, WRCAS),能夠通過具有固定后掠角的前緣在上表面產(chǎn)生穩(wěn)定分離渦從而提供更高的升力,而無需將上表面設(shè)計(jì)為犧牲體積的膨脹面;更重要的是這種乘波體的前緣后掠角能夠作為設(shè)計(jì)參數(shù)在設(shè)計(jì)階段予以控制。

    定后掠角乘波體由來已久,但在早期只是其他型乘波體的“附屬品”。由Nonweiler提出來的“Λ”型乘波體就是一種定后掠角乘波體,但是由于其升阻比不高,體積利用率低等問題,難以實(shí)用。Jones在對(duì)錐導(dǎo)乘波體進(jìn)行介紹時(shí),提到了一種流動(dòng)捕獲管(Flow Capture Tube, FCT)過圓錐頂點(diǎn)的乘波體,這種乘波體也屬于定后掠角乘波體,也因?yàn)閷?shí)用價(jià)值不高未受重視。真正以定后掠角前緣為研究重點(diǎn)的是洛克希德·馬丁公司的Rodi[15-17],其在文獻(xiàn)中對(duì)定后掠角密切錐乘波體(Osculating Cone Waverider with Constant Angle of Sweepback,OCWRCAS)和定后掠角密切流場(chǎng)乘波體(Osculating Flow Waverider with Constant Angle of Sweepback, OFWRCAS)進(jìn)行了研究,其研究重點(diǎn)在于此類乘波體的幾何特征和渦升力特性,卻并未說明怎樣生成定后掠角前緣對(duì)應(yīng)的乘波體下表面,而這也正是這種乘波體的重點(diǎn)之一,因?yàn)橹挥泻侠淼南卤砻娌拍艽_保乘波體的乘波特性。

    本文首先根據(jù)OCWRCAS的幾何特征及密切錐乘波體的生成方法,給出了定后掠角前緣對(duì)應(yīng)的乘波體下表面生成方法;之后提取了兩種典型定后掠角密切錐乘波體的設(shè)計(jì)變量,以遍歷設(shè)計(jì)空間的思路進(jìn)行了多目標(biāo)優(yōu)化設(shè)計(jì);最后使用CFD方法對(duì)此類乘波體的乘波特性和渦升力特性進(jìn)行了驗(yàn)證。

    1 OCWRCAS的生成及氣動(dòng)性能計(jì)算方法

    1.1 OCWRCAS的生成方法

    OCWRCAS是密切錐乘波體的一種,生成方法與其類似,因此本文從密切錐乘波體的生成方法入手,引出OCWRCAS的生成方法。密切錐乘波體的生成方法是一種反設(shè)計(jì)方法,這種方法由給定的激波形狀來確定乘波外形。所以在設(shè)計(jì)時(shí)除了要給定流動(dòng)捕獲管FCT,還需要給出激波形狀,即進(jìn)氣捕獲曲線(Inlet Capture Curve, ICC),該曲線與進(jìn)氣道進(jìn)口截面形狀有關(guān)。這兩條曲線必須一階導(dǎo)數(shù)連續(xù),ICC還要保證二階導(dǎo)數(shù)非負(fù)。對(duì)于任意的ICC,其對(duì)應(yīng)的流場(chǎng)一般都是三維的,直接計(jì)算這種流場(chǎng)計(jì)算量較大不利于快速生成乘波體。Sobieczky等[8]提出的相切錐理論大大簡(jiǎn)化了由激波形狀反向計(jì)算流場(chǎng)的過程,首先定義過ICC上任意一點(diǎn)并且垂直于這個(gè)點(diǎn)切線方向的平面為密切平面,然后使用一系列密切平面內(nèi)的錐型流近似這種三維流動(dòng),這樣就大大提高了計(jì)算效率。密切錐的頂點(diǎn)可由密切平面內(nèi)對(duì)應(yīng)點(diǎn)的曲率半徑和激波角確定,當(dāng)曲率半徑無窮大時(shí),當(dāng)做二維流場(chǎng)處理。密切錐確定之后,可根據(jù)錐型流理論求得解析解,再由FCT和激波形狀確定該密切面內(nèi)的前緣點(diǎn),并以該點(diǎn)作為起始點(diǎn)在密切平面內(nèi)做流線追蹤,得到一條流線。所有密切平面上的流線組成的流面即為密切錐乘波體的下表面。生成密切錐乘波體的兩個(gè)充分條件是:① 每個(gè)密切面內(nèi)流場(chǎng)的激波必須與激波形狀(ICC)一致;② 相鄰密切面內(nèi)的橫向流動(dòng)必須足夠小。要滿足這兩個(gè)條件,一般都限定每個(gè)密切面內(nèi)激波角相等。

    在滿足上述條件的前提下,要得到具有一定后掠角的直線前緣,必須保證前緣對(duì)應(yīng)的激波面和流動(dòng)捕獲管都是平面。鑒于此,定后掠角密切錐乘波體與一般密切錐乘波的不同之處在于:① ICC上使用一條直線段來生成帶有特定后掠角的前緣;② 對(duì)于ICC上曲率半徑無限大的情況要進(jìn)行判斷,如果是對(duì)應(yīng)發(fā)動(dòng)機(jī)進(jìn)口的直線段,使用二維流場(chǎng)計(jì)算,如果對(duì)應(yīng)后掠前緣則使用錐型流計(jì)算;③ 直線ICC的每個(gè)密切面上的密切錐頂點(diǎn)需要滿足一定的分布規(guī)律以使激波面為平面,本文取與直線激波平行的直線;④ FCT上與直線ICC對(duì)應(yīng)的部分應(yīng)為直線,為了使直線前緣處于水平面內(nèi),取水平直線段。

    圖1 兩種乘波體的密切錐示意圖Fig.1 Schematic diagram of osculating cone of two types of waverider

    圖1給出了密切錐乘波體(OCWR)和OCWRCAS的密切錐示意圖,實(shí)線分別表示FCT和ICC,虛線表示密切平面,“△”表示ICC上密切平面內(nèi)的點(diǎn),“□”表示對(duì)應(yīng)密切錐的頂點(diǎn)。OCWR的ICC比較任意,密切錐頂點(diǎn)由曲線當(dāng)?shù)氐那拾霃酱_定,由圖1(a)可以看出,其密切錐頂點(diǎn)分布與FCT無關(guān),靠近對(duì)稱面附近的水平直線曲率半徑為無限大,計(jì)算時(shí)使用二維方法。OCWRCAS的FCT和ICC如圖1(b)所示,OCWRCAS的ICC由圓弧(紅色部分)和直線段(藍(lán)色部分)組成,圓弧對(duì)應(yīng)曲線型頭部,直線對(duì)應(yīng)可控后掠前緣,其對(duì)應(yīng)密切錐頂點(diǎn)沿與其平行的直線分布。

    1.2 圓錐流場(chǎng)及黏性力計(jì)算方法

    由定后掠角密切錐乘波體的生成方法可知,生成該型乘波體需要錐型流場(chǎng)。本文使用Taylor-Maccoll方程[18]計(jì)算錐型流場(chǎng),該方程的無量綱形式可以表示為

    (1)

    (2)

    (3)

    定義無量綱速度為

    (4)

    式中:Vmax為氣體膨脹至絕對(duì)零度時(shí)的理想最大速度。該方程為常微分方程,可采用四階Runge-Kutta求解。

    乘波體的黏性阻力與波阻處于同一量級(jí),因此升阻比分析必須考慮黏性阻力。本文采用文獻(xiàn)[19]中的參考溫度方法,該方法比邊界層積分方法的精度稍低但計(jì)算速度卻要快得多,因此具有很大的實(shí)用價(jià)值。本文假設(shè)乘波體表面全部為湍流,并對(duì)不可壓平板的黏性力計(jì)算方法進(jìn)行了可壓縮修正,提高了計(jì)算精度。首先由式(5)計(jì)算參考溫度:

    (5)

    式中:T*為參考溫度;Tδ為邊界層外緣溫度;Maδ為邊界層外緣馬赫數(shù);TW為壁面溫度,可在計(jì)算時(shí)指定或者由氣動(dòng)熱計(jì)算方法計(jì)算。參考雷諾數(shù)為

    (6)

    式中:Vδ為邊界層外緣速度;μ*=fT*為參考黏性系數(shù),是參考溫度的函數(shù),可由撒特蘭公式求得;ρ*=gp*,T*為參考密度,可根據(jù)氣體狀態(tài)方程求得。附面層內(nèi)壓力梯度較小,所以取p*=pδ,pδ為邊界層外緣壓力。不可壓平板湍流黏性力系數(shù)為

    (7)

    可壓縮黏性力系數(shù)的修正因子為

    Fc=ρδ/ρ*

    (8)

    式中:ρδ為邊界層外緣密度。

    最終求得當(dāng)?shù)仞ば粤ο禂?shù)為

    (9)

    2 OCWRCAS的設(shè)計(jì)方法

    OCWRCAS的生成方法屬于密切錐乘波體的一種,為便于研究,本文將FCT設(shè)定為水平直線段,其外形就主要由ICC控制,因此生成此類乘波體的關(guān)鍵就是ICC曲線。一般密切錐乘波體的ICC只要滿足一階導(dǎo)數(shù)連續(xù),二階導(dǎo)數(shù)非負(fù)的條件即可,而OCWRCAS則必須有一條位于遠(yuǎn)離對(duì)稱面一側(cè)的直線段來生成具有一定后掠角的前緣。

    本文對(duì)其中兩種典型外形進(jìn)行研究,分別用OCWRCAS I和OCWRCAS II表示。OCWRCAS I頭部為尖頂點(diǎn),即整個(gè)前緣為后掠角可控的直線段前緣;OCWRCAS II頭部為曲線,除頭部外的部分前緣為后掠角可控的直線段前緣。而且ICC中的曲線都選擇圓弧,這樣減少了設(shè)計(jì)變量的個(gè)數(shù),便于將研究重點(diǎn)放在后掠角部分。通過簡(jiǎn)化,這兩種乘波體的設(shè)計(jì)變量都很少,因此文中采用遍歷的思路進(jìn)行優(yōu)化設(shè)計(jì)。

    乘波體的設(shè)計(jì)目標(biāo)一般都為升阻比最大,同時(shí)要求具有一定的體積效率。體積效率的定義方法有多種[20],本文使用

    τ=V2/3/Sp

    (10)

    式中:Sp為乘波體的平面面積。本文乘波體的設(shè)計(jì)狀態(tài)均為:30 km高度,來流馬赫數(shù)Ma∞=6,乘波體長(zhǎng)度取20 m。

    2.1 OCWRCAS I的設(shè)計(jì)

    OCWRCAS I的ICC由兩部分組成:圓弧段和直線段,并且圓弧的圓心O都位于FCT上。這類乘波體的俯視圖和后視圖(側(cè)視圖基本相同,此處不予畫出)如圖2所示,整個(gè)前緣為一條直線段VE,具有相同的后掠角,頭部為一個(gè)尖點(diǎn)。本文從其幾何特征入手,分析影響乘波體性能的設(shè)計(jì)變量,并提出了具有針對(duì)性的設(shè)計(jì)方法。

    圖2 OCWRCAS I的幾何特征Fig.2 Geometry characteristic of OCWRCAS I

    2.1.1 幾何特征分析

    如圖2所示,OCWRCAS I的ICC由兩部分組成:圓心在O點(diǎn)的圓弧RA和直線段AE,且兩者相切于A點(diǎn)。圓弧部分對(duì)應(yīng)頭部的頂點(diǎn),直線部分對(duì)應(yīng)具有某一后掠角的前緣。這樣的曲線由兩個(gè)變量控制,一個(gè)是圓弧的半徑r,一個(gè)是直線與FCT的夾角γ。定義機(jī)身長(zhǎng)度為l,機(jī)身寬度為s,激波角為β,后掠角為λ。

    r可由機(jī)身長(zhǎng)度和激波角確定,即

    r=OA=OR=l·tanβ

    γ可根據(jù)相切關(guān)系求得,即

    sinγ=OA/OE

    式中:OE為機(jī)身寬度,可由機(jī)身長(zhǎng)度和后掠角確定,即

    s=OE=l·tanλ

    最終可得[16]

    sinγ=tanλ·tanβ

    (11)

    與一般的密切錐乘波體相同,當(dāng)給定來流馬赫數(shù)和設(shè)計(jì)激波角后,乘波體下表面只和ICC有關(guān)。根據(jù)上述分析,對(duì)于OCWRCAS I,控制ICC的變量為l、β和λ。所以影響OCWRCAS I性能的變量有4個(gè):Ma∞、l、β和λ。Ma∞和l在設(shè)計(jì)前給定,所以O(shè)CWRCAS I最終的設(shè)計(jì)變量為激波角和后掠角,即β和λ。接下來分析β和λ的變化范圍。

    當(dāng)給定后掠角后,由式(10)可確定激波角上限為

    tanλ·tanβ<1

    β<90°-λ

    對(duì)于來流馬赫數(shù)Ma∞,馬赫波可以看做是該馬赫數(shù)下最弱的激波,因此可認(rèn)為馬赫角是最小的激波角,所以激波角的范圍為

    arcsin1/Ma∞<β<90°-λ

    (12)

    同時(shí)可得后掠角的上限為

    λ<90°-arcsin(1/Ma∞)

    (13)

    當(dāng)Ma∞=6時(shí),arcsin1/Ma∞=9.59°,此時(shí)λ<80.41°,取λ=70° 的乘波體為例,激波角的變化范圍為:9.59°<β<20°;若后掠角λ=65°,則激波角的變化范圍為:9.59°<β<35°。高超聲速飛行器的后掠角不能太小,因此本文取后掠角的下限為50°。

    通過上述分析,本文分別在后掠角為75°、70°、65°、60°、55°、50° 的情況下,對(duì)不同激波角下的乘波體進(jìn)行了分析,具體分配如表1所示,每個(gè)后掠角下激波角的變化間隔都為0.5°,共對(duì)210個(gè)乘波體進(jìn)行了分析。

    表1 OCWRCAS I設(shè)計(jì)變量分布Table 1 Distribution of design variables of OCWRCAS I (°)

    2.1.2 設(shè)計(jì)結(jié)果分析

    圖3 不同后掠角下升阻比和體積效率隨激波角變化曲線Fig.3 Lift to drag ratio and volumetric efficiency vs shock angle at different angles of sweepback

    不同后掠角下氣動(dòng)性能和體積效率隨激波角變化如圖3所示。圖3(a)給出了不同后掠角時(shí)升阻比隨激波角的變化曲線,可以看出,不同后掠角時(shí)升阻比CL/CD隨激波角的變化趨勢(shì)類似,都呈現(xiàn)先增大后減小的趨勢(shì)。升阻比的最大值都出現(xiàn)在11° 激波角時(shí),而且最大值隨著后掠角的增大而減小,但總體變化幅度不大。圖3(b)給出了不同后掠角時(shí)體積效率τ隨激波角的變化曲線,不同后掠角時(shí)體積效率都隨激波角的增大而增大,若給定激波角,體積效率隨后掠角的增大而增大。

    圖4給出了OCWRCAS I隨后掠角和激波角的變化規(guī)律,但并不利于從中選擇“最優(yōu)結(jié)果”。本文對(duì)圖3中的結(jié)果按升阻比隨體積效率的變化輸出,如圖4所示,這樣就產(chǎn)生了類似Pareto前沿線的多目標(biāo)優(yōu)化結(jié)果,該前沿線上的結(jié)果包含了升阻比最大、體積效率最大以及升阻比和體積效率適度的所有乘波體。由圖可知,后掠角為75°,激波角大于10.5°,所有乘波體都處于前沿線上,并且屬于升阻比較高部分。根據(jù)文獻(xiàn)[20],一般高超聲速飛行器的體積效率不應(yīng)低于0.2,本文以此為體積效率的下限,選擇升阻比最大的乘波體,圖中箭頭標(biāo)出了“最優(yōu)”的乘波體外形(Selected)。表2給出了升阻比最大、體積效率最大和“最優(yōu)”乘波體的性能參數(shù),此3種外形的三視圖如圖5所示?!白顑?yōu)”乘波體升阻比為6.58,體積效率為0.222 1,后掠角為75°,激波角為 12.0°;升阻比最大和體積效率最大的乘波體都為50° 后掠角乘波體,最大升阻比為7.53,最大體積效率為1.14。

    圖4 OCWRCAS I升阻比隨體積效率的變化規(guī)律 Fig.4 Lift to drag ratio vs volumetric efficiency for OCWRCAS I

    表2 OCWRCAS I優(yōu)化設(shè)計(jì)結(jié)果Table 2 Optimal design results of OCWRCAS I

    TypeCL/CDV2/3/Spλ/(°)β/(°)MaxCL/CD7.530.11355011.0MaxV2/3/Sp1.140.58615039.0Selected6.580.22217512.0

    圖5 OCWRCAS I 3種外形的三視圖Fig.5 Three views of three types of OCWRCAS I

    2.2 OCWRCAS II的設(shè)計(jì)

    OCWRCAS II是指頭部為光滑曲線的定后掠角密切錐乘波體。將圖2中圓心O的位置上移,其他方法不變,便可得到生成OCWRCAS II的ICC,如圖6所示。

    圖6 OCWRCAS II的幾何特征Fig.6 Geometry characteristic of OCWRCAS II

    2.2.1 幾何特征分析

    如圖6所示,ICC的直線段AE對(duì)應(yīng)乘波體的后掠前緣ME,這部分與OCWRCAS I的前緣相同,所以式(10)同樣適用,不同之處在于OCWRCAS I的O點(diǎn)與O1是重合的,而OCWRCAS II的O點(diǎn)與O1點(diǎn)存在一定距離,定義為Δr,有限的Δr確保了乘波體靠近頭部的前緣為曲線形狀,該值越大前緣曲線部分所占比例就越大。由Δr與其他已知變量可得圓弧半徑r,聯(lián)立:

    最終得

    (14)

    保持機(jī)身長(zhǎng)度和機(jī)身厚度不變,Δr越大,前緣曲線部分所占的寬度就越大,機(jī)身寬度也越大,這在一定程度上增大了乘波體的體積。但是Δr不能無限制地增大,圖7給出了Δr的變化范圍示意圖,保持激波角和后掠角不變,只增大Δr。圖7中Δr2對(duì)應(yīng)的ICC為RA2E2曲線,相比較于Δr1對(duì)應(yīng)的RA1E1曲線機(jī)身寬度增加,可以預(yù)見體積也有明顯增大。當(dāng)達(dá)到Δrmax時(shí),ICC為RE3整個(gè)一條圓弧,相應(yīng)乘波體為錐導(dǎo)乘波體的一部分,但是其前緣離圓錐軸線非常遠(yuǎn),體積效率很低。Δr作為控制變量,便于說明問題,但是不宜作為設(shè)計(jì)變量,因?yàn)椴煌げń堑摩變化范圍不同,難以進(jìn)行統(tǒng)一量化。本文使用前緣曲線程度變量Δw/s來替換Δr,其中Δw為圓弧對(duì)應(yīng)FCT上的寬度,很明顯該變量表征了前緣曲線部分所占比例的大小,Δw/s與Δr和其他已知量之間的關(guān)系為

    (15)

    通過幾何特征分析,OCWRCAS II相對(duì)于OCWRCAS I多了一個(gè)設(shè)計(jì)變量。對(duì)其設(shè)計(jì)時(shí),后掠角和激波角仍然沿用OCWRCAS I的設(shè)定,對(duì)于Δw/s,為0時(shí)生成的乘波體為OCWRCAS I,為1時(shí)生成的乘波體為錐導(dǎo)乘波體,因此本文選定其變化范圍為0.04~0.96,變化間隔為0.04,如表3所示,共對(duì)5 040個(gè)乘波體進(jìn)行了分析。

    圖7 Δr變化范圍示意圖Fig.7 Schematic diagram of variation scope of Δr

    表3 OCWRCAS II設(shè)計(jì)變量分布Table 3 Distribution of design variables of OCWRCAS II

    λ/(°)βmin/(°)βmax/(°)Δβ/(°)(Δw/s)min(Δw/s)maxΔ(Δw/s)5010.039.50.50.040.960.045510.034.50.50.040.960.046010.029.50.50.040.960.046510.024.50.50.040.960.047010.019.50.50.040.960.047510.014.50.50.040.960.04

    2.2.2 設(shè)計(jì)結(jié)果分析

    圖8 不同激波角下升阻比和體積效率隨Δw/s的變化曲線Fig.8 Lift to drag ratio and volumetric efficiency vs Δw/s at different shock angles

    在每個(gè)特定Δw/s下,OCWRCAS II隨激波角的變化趨勢(shì)與圖3(a)一致,并且都是在11° 激波角時(shí)出現(xiàn)升阻比最大的情況,本節(jié)不再給出OCWRCAS II隨激波角變化曲線,重點(diǎn)分析Δw/s大小對(duì)乘波體性能的影響。以75° 后掠角的設(shè)計(jì)結(jié)果為代表進(jìn)行研究,圖8給出了該型乘波體的升阻比和體積效率隨Δw/s的變化曲線,圖中曲線根據(jù)激波角的大小分為3類:小于11°,升阻比隨激波角的增大而增大,圖中用虛線畫出;大于11°,升阻比隨激波角的增大而減小,圖中用點(diǎn)線畫出;等于11°,圖中用實(shí)線表示。圖8(a)也給出了升阻比隨Δw/s的變化曲線,小于11° 的曲線升阻比隨著Δw/s的增大先增后減,減小幅度隨激波角的增大逐漸降低,11° 時(shí)的升阻比也為先增大后減小,但減小趨勢(shì)非常低,到大于11° 時(shí)就變?yōu)殡S著Δw/s的增大而增大。圖8(b)給出了體積效率隨Δw/s的變化曲線,其變化規(guī)律較升阻比相比要簡(jiǎn)單,都是隨著Δw/s的增大而減小。圖8中的結(jié)果說明Δw/s的增大對(duì)升阻比而言有有利情況也有不利情況,對(duì)體積效率則完全是不利的,因此要根據(jù)上述結(jié)論在設(shè)計(jì)空間中選擇“最優(yōu)”的結(jié)果很困難,還需要給出升阻比隨體積效率的變化趨勢(shì)。

    同樣給出OCWRCAS II的升阻比隨體積效率的變化規(guī)律,如圖9所示,得到類似Pareto最優(yōu)解的前沿線,占據(jù)前沿線最大比例的為后掠角為75° 的乘波體,這個(gè)結(jié)論與OCWRCAS I一致,說明對(duì)于定后掠角密切錐乘波體,后掠角最大的乘波體是升阻比和體積效率都較高的一類乘波體,具有最高的實(shí)用價(jià)值。

    如圖10所示,本節(jié)使用與OCWRCAS I相同的選擇標(biāo)準(zhǔn),對(duì)OCWRCAS II進(jìn)行分析。表4給出了升阻比最大、體積效率最大和本文選擇結(jié)果的性能參數(shù),升阻比最大的乘波體生成激波角為12°,后掠角為50°,Δw/s為最大值0.96,升阻比可達(dá)9.05,圖10(a)給出了升阻比最大的OCWRCAS II三視圖,從三視圖可以看出該型乘波體非常接近錐導(dǎo)乘波體,已經(jīng)沒有明顯的后掠前緣,而且厚度非常小,所以有很大的升阻比,但體積效率僅有0.083 9;體積效率最大的乘波體后掠角同樣為50°,生成激波角為最大激波角 39.5°,Δw/s為0.04,圖10(b)給出了體積效率最大的OCWRCAS II三視圖,該外形曲線型頭部非常小,具有與圓錐接近的外形,所以具有很大的體積效率,但升阻比只有1.15;本節(jié)選擇的“最優(yōu)”結(jié)果的生成激波角為12.5°,后掠角為75°,Δw/s為0.36,該型乘波體的三視圖如圖10(c)所示,其升阻比為7.15,優(yōu)于OCWRCAS I選擇的外形,但體積效率較OCWRCAS I要低。

    圖9 OCWRCAS II升阻比隨體積效率的變化規(guī)律Fig.9 Lift to drag ratio vs volumetric efficiency for OCWRCAS II

    圖10 OCWRCAS II 3種外形的三視圖Fig.10 Three views of three types of OCWRCAS II

    表4 OCWRCAS II優(yōu)化設(shè)計(jì)結(jié)果Table 4 Optimal design results of OCWRCAS II

    3 CFD驗(yàn)證

    3.1 乘波特性

    以O(shè)CWRCAS II中的“最優(yōu)”乘波體為例,使用CFD方法進(jìn)行黏性流場(chǎng)分析,驗(yàn)證乘波體的乘波氣動(dòng)特性。氣動(dòng)力結(jié)果比較如表5所示,可見CFD計(jì)算結(jié)果與理論結(jié)果是相近的。

    圖11給出乘波體下表面流場(chǎng)中的壓力等值線,說明該型乘波體激波整體附著在前緣,使流場(chǎng)在前緣處泄露非常少,具有很好的乘波特性。圖12 給出了乘波體底部處的壓力分布等值線以及設(shè)計(jì)該型乘波體時(shí)的ICC曲線,ICC即底部的理論激波形狀,由等值線反映出的激波形狀與ICC符合得很好,說明本文給出的設(shè)計(jì)方法是一種有效的密切錐乘波體設(shè)計(jì)方法。文獻(xiàn)[21]通過計(jì)算指出,雖然下表面的理論結(jié)果存在一定誤差,但由于數(shù)值差別較小而且面積不大,最終的升力系數(shù)和阻力系數(shù)仍然能與CFD計(jì)算結(jié)果保持一致,所以使用理論結(jié)果進(jìn)行設(shè)計(jì)是有效的。

    表5“最優(yōu)”O(jiān)CWRCASII的氣動(dòng)力結(jié)果比較

    Table5Comparisonofaerodynamicforceresultsof“optimal”O(jiān)CWRCASII

    MethodCL(CD)inviscous(CD)viscousCL/CDTheory0.02030.00160.00127.15CFD0.02090.00170.00146.74

    圖11 設(shè)計(jì)狀態(tài)時(shí)下表面流場(chǎng)的壓力等值線 Fig.11 Pressure contour at lower surface with design status

    圖12 理論激波形狀(”)與CFD計(jì)算結(jié)果(壓力等值線)比較Fig.12 Comparison of shock wave between theoretical shape (”) and CFD (pressure contour) computational results

    3.2 渦升力

    本節(jié)使用CFD方法分別對(duì)OCWRCAS II的“最優(yōu)”外形進(jìn)行了研究。計(jì)算狀態(tài)為:30 km,馬赫數(shù)6,全湍流。計(jì)算迎角α為0°、2°、4°、6°、8°、10°和12°。OCWRCAS II上表面的壓力分布系數(shù)Cp如圖13所示,可知該型乘波體隨著迎角增大,上表面壓力漸趨減小,由渦誘導(dǎo)的低壓區(qū)域逐漸增大。圖14給出了升力系數(shù)隨迎角變化曲線,在6° 迎角之后,升力曲線呈現(xiàn)明顯非線性增加趨勢(shì)。本節(jié)計(jì)算說明,定后掠角密切錐乘波體在大迎角時(shí)具有高升力特性。但是高超聲速時(shí),飛行器上表面的分離情況是非常復(fù)雜的,其生成、發(fā)展機(jī)理以及影響因素有待進(jìn)一步研究。

    圖13 OCWRCAS II “最優(yōu)”外形不同迎角時(shí)上表面壓力云圖Fig.13 Pressure contour at upper surface of “Optimal” OCWRCAS II shape with different angles of attack

    圖14 OCWRCAS II “最優(yōu)”外形升力系數(shù)隨迎角變化曲線 Fig.14 Lift coefficient of “Optimal” OCWRCAS II vs angle of attack

    4 結(jié) 論

    1) 本文給出的OCWRCAS定后掠角前緣對(duì)應(yīng)下表面的生成方法切實(shí)有效,得到的乘波體乘波特性明顯。

    2) 后掠角最大的定后掠角密切錐乘波體,絕大部分都位于“Pareto前沿線”上,即這類乘波體處于升阻比和體積效率均“較優(yōu)”的位置。

    3) 對(duì)于某一后掠角的定后掠角密切錐乘波體,升阻比隨激波角的增大呈現(xiàn)先增大后減小的趨勢(shì),體積效率隨激波角的增大逐漸增大。

    4) 對(duì)于OCWRCAS II,給定激波角和后掠角時(shí),升阻比隨前緣曲線程度的變化趨勢(shì)與激波角的大小有關(guān),體積效率則是隨前緣曲線程度的增大而增大。

    5) 通過對(duì)多個(gè)迎角的CFD計(jì)算,說明定后掠角密切錐乘波體在一定迎角時(shí),上表面會(huì)產(chǎn)生穩(wěn)定的分離,從而改善升力特性。

    6) 本文在對(duì)定后掠角密切錐乘波體進(jìn)行設(shè)計(jì)時(shí),組成ICC的曲線部分使用了圓弧,限制了設(shè)計(jì)空間,若采用一般曲線進(jìn)行優(yōu)化設(shè)計(jì)有望得到性能更優(yōu)的乘波體。

    [1] NONWEILER T R F. Aerodynamic problems of manned space vehicles[J]. Journal of the Royal Aeronautical Society, 1959, 63(585): 521-528.

    [2] JONES J G, MOORE K C, PIKE J, et al. A method for designing lifting configurations for high supersonic speeds using axisymmetric flow field[J]. Archive of Applied Mechanics, 1968, 37(1): 56-72.

    [3] 耿永兵, 劉宏, 姚文秀, 等. 錐形流乘波體優(yōu)化設(shè)計(jì)研究[J]. 航空學(xué)報(bào), 2006, 27(1): 23-28.

    GENG Y B, LIU H, YAO W X, et al. Viscous optimized design of waverider derived from cone flow [J]. Acta Aeronautica et Astronautica Sinica, 2006, 27(1): 23-28 (in Chinese).

    [4] 耿永兵, 劉宏, 雷麥芳, 等. 高升阻比乘波構(gòu)型優(yōu)化設(shè)計(jì)[J]. 力學(xué)學(xué)報(bào), 2006, 38(4): 540-546.

    GENG Y B, LIU H, LEI M F, et al. Optimized design of waverider with high lift over drag ratio[J]. Chinese Journal of Theoretical and Applied Mechanics, 2006, 38(4): 540-546 (in Chinese).

    [5] RASMUSSEN M L. Waverider configurations derived from inclined circular and elliptic cones[J]. Journal of Spacecraft and Rockets, 1980, 17(6): 537-545.

    [6] 樂貴高, 馬大為, 李自勇. 橢圓錐乘波體高超聲速流場(chǎng)數(shù)值計(jì)算[J]. 南京理工大學(xué)學(xué)報(bào)(自然科學(xué)版), 2006, 30(3): 257-260.

    LE G H, MA D W, LI Z Y. Computation of hypersonic flowfields for elliptic-cone-derived waverider[J]. Journal of Nanjing University of Science and Technology (Natural Science), 2006, 30(3): 257-260 (in Chinese).

    [7] TAKASHIMA N, LEWIS M J. Waverider configurations based on non-axisymmetric flow fields for engine-airframe integration: AIAA-1994-0380[R]. Reston: AIAA, 1994.

    [8] SOBIECZKY H, DOUGHERTY F C, JONES K. Hypersonic waverider design from given shock wave[C]//First International Waverider Symposium. Maryland: University of Maryland, 1990.

    [9] CENTER K, SOBIECZKY H, DOUGHERTY F C. Interactive design and analysis of hypersonic waverider geometries: AIAA-1991-1697[R]. Reston: AIAA, 1991.

    [10] KONTOGIANNIS K, SOBESTER A, TAYLOR N J. On the conceptual design of waverider forebody feometries: AIAA-2015-1009[R]. Reston: AIAA, 2015.

    [11] SZEMA K, LIU Z, MUNIPALLI R. An efficient GUI design tool for high-speed airbreathing propulsion integration: AIAA-2010-4362[R]. Reston: AIAA, 2010.

    [12] 賀旭照, 倪鴻禮. 密切曲面錐乘波體-設(shè)計(jì)方法與性能分析[J]. 力學(xué)學(xué)報(bào), 2011, 43(6): 1077-1082.

    HE X Z, NI H L. Osculating curved cone(OCC) waverider: Design methods and performance analysis[J]. Chinese Journal of Theoretical and Applied Mechanics, 2011, 43(6): 1077-1082 (in Chinese).

    [13] 賀旭照, 周正, 毛鵬飛, 等. 密切曲面內(nèi)錐乘波前體進(jìn)氣道設(shè)計(jì)和試驗(yàn)研究[J]. 實(shí)驗(yàn)流體力學(xué), 2014, 28(3): 39-44.

    HE X Z, ZHOU Z, MAO P F, et al. Design and experimental study of osculating inward turning cone waverider/inlet(OICWI)[J]. Journal of Experiments in Fluid Mechanics, 2014, 28(3): 39-44 (in Chinese).

    [14] BOWCUTT K G. Optimization of hypersonic waveriders derived from cone flows-including viscous effects[D]. Maryland: University of Maryland, 1986.

    [15] RODI P E. The osculating flowfield method of waverider geometry generation: AIAA-2005-0511[R]. Reston: AIAA, 2005.

    [16] RODI P E. Geometrical relationships for osculating cones and osculating flowfield waverider: AIAA-2011-1188[R]. Reston: AIAA, 2011.

    [17] RODI P E. Vortex lift waverider configurations: AIAA-2012-1238[R]. Reston: AIAA, 2012.

    [18] ANDERSON J D. Modern compressible flow[M]. 2nd ed. New York: McGraw-Hill Publishing Company, 1999: 294-307.

    [19] CORDA S, ANDERSON J. Viscous optimized hypersonic waveriders designed from axisymmetric flow fields: AIAA-1988-0369[R]. Reston: AIAA, 1988.

    [20] DIETER J, GOTTFRIED S, SIEGFRIED W. Basic research and technologies for two-stage-to-orbit vehicles[M]. Weinheim: WILEY-VCH Verlag Gmbh & Co. KGaA, 2005: 215-260.

    [21] CENTER K, SOBIECZKY H, DOUGHERTY F. Interactive design of hypersonic waverider geometries: AIAA-1991-1697[R].Reston: AIAA, 1991.

    段焰輝男, 博士。主要研究方向: 計(jì)算流體力學(xué), 氣動(dòng)外形優(yōu)化設(shè)計(jì)。

    E-mail: duanyanhui@foxmail.com

    吳文華男, 研究員。主要研究方向: 空氣動(dòng)力學(xué), 氣動(dòng)外形優(yōu)化設(shè)計(jì)。

    Tel.: 0816-2463132

    E-mail: 619677947@qq.com

    URL:www.cnki.net/kcms/detail/11.1929.V.20160226.1344.002.html

    Generationanddesignmethodsofosculatingconewaveriderwithconstantangleofsweepback

    DUANYanhui,FANZhaolin,WUWenhua*

    ComputationalAerodynamicsResearchInstitute,ChinaAerodynamicsResearchandDevelopmentCenter,Mianyang621000,China

    Inthispaper,thegenerationanddesignmethodsofosculatingconewaveriderwithconstantangleofsweepback(OCWRCAS)arestudied,andviscousisconsideredduringdesign.Stablevortexcanbegeneratedbytheleadingedgewithconstantangleofsweepback,whichwillimprovetheaerodynamicabilityofthewaverider.Firstly,thegenerationmethodofOCWRCASispresentedbasedonthegenerationmethodoftraditionalosculatingconewaverider.Thedesignvariablesofsweepbackangle,shockangleandthecurveshapeofheadareextractedbyanalyzingthegeometrycharacteroftheOCWRCAS,andthevariationtrendoflifttodragratioandvolumetricefficiencywiththesevariablesisalsostudied.Themulti-objectiveoptimalsolutionsarefoundedfromtwoclassicaltypesofOCWRCASbysearchingthetotaldesignspace.Finally,themethodofcomputationalfluiddynamicsisusedtoverifythecharacterofwaveridingandvortexlift.TheresultsshowthatOCWRCASwithgoodabilityofwaveridingandhighlifttodragratiokeepsrelationalvolumetricefficiency;vortexliftcanbegeneratedbytheleadingedgewithconstantangleofsweepbackatcertainangleofattack.

    waverider;osculatingcone;viscous;angleofsweepback;vortexlift

    2015-12-14;Revised2015-12-30;Accepted2016-01-18;Publishedonline2016-02-261344

    .Tel.:0816-2463132E-mail619677947@qq.com

    2015-12-14;退修日期2015-12-30;錄用日期2016-01-18; < class="emphasis_bold">網(wǎng)絡(luò)出版時(shí)間

    時(shí)間:2016-02-261344

    www.cnki.net/kcms/detail/11.1929.V.20160226.1344.002.html

    .Tel.:0816-2463132E-mail619677947@qq.com

    段焰輝, 范召林, 吳文華. 定后掠角密切錐乘波體的生成和設(shè)計(jì)方法J. 航空學(xué)報(bào),2016,37(10):3023-3034.DUANYH,FANZL,WUWH.GenerationanddesignmethodsofosculatingconewaveriderwithconstantangleofsweepbackJ.ActaAeronauticaetAstronauticaSinica,2016,37(10):3023-3034.

    http://hkxb.buaa.edu.cnhkxb@buaa.edu.cn

    10.7527/S1000-6893.2016.0024

    V211.5

    A

    1000-6893(2016)10-3023-12

    猜你喜歡
    后掠角激波前緣
    不同后掠角大展弦比復(fù)合材料機(jī)翼氣動(dòng)特性
    一種基于聚類分析的二維激波模式識(shí)別算法
    基于HIFiRE-2超燃發(fā)動(dòng)機(jī)內(nèi)流道的激波邊界層干擾分析
    一種飛機(jī)尾翼前緣除冰套安裝方式
    斜激波入射V形鈍前緣溢流口激波干擾研究
    適于可壓縮多尺度流動(dòng)的緊致型激波捕捉格式
    基于CFD的最優(yōu)變后掠規(guī)律研究
    不同后掠角三角翼的靜態(tài)地面效應(yīng)數(shù)值模擬
    深水沉積研究進(jìn)展及前緣問題
    亚洲成人手机| 国产精品免费大片| 日日爽夜夜爽网站| 丝袜在线中文字幕| 国产成人freesex在线| 交换朋友夫妻互换小说| 一级毛片 在线播放| 国产精品一区二区在线不卡| 中文乱码字字幕精品一区二区三区| 三上悠亚av全集在线观看 | 自拍偷自拍亚洲精品老妇| 精品久久久久久久久av| 中文乱码字字幕精品一区二区三区| 伊人久久国产一区二区| 国产精品不卡视频一区二区| 婷婷色麻豆天堂久久| 亚州av有码| 国模一区二区三区四区视频| 亚洲av免费高清在线观看| 高清欧美精品videossex| 国产精品99久久99久久久不卡 | 欧美最新免费一区二区三区| 视频区图区小说| 视频中文字幕在线观看| 欧美xxxx性猛交bbbb| 午夜激情福利司机影院| 黄色视频在线播放观看不卡| 麻豆精品久久久久久蜜桃| 欧美成人精品欧美一级黄| 国产亚洲一区二区精品| 亚洲精品国产av蜜桃| 国产在线一区二区三区精| 乱人伦中国视频| 国产免费又黄又爽又色| 内射极品少妇av片p| 成人亚洲欧美一区二区av| 精品国产一区二区三区久久久樱花| 久久久亚洲精品成人影院| 国产真实伦视频高清在线观看| av在线观看视频网站免费| 日本wwww免费看| 97在线视频观看| 亚洲色图综合在线观看| 国内精品宾馆在线| 免费看av在线观看网站| 免费少妇av软件| 国产黄片美女视频| 日本vs欧美在线观看视频 | 精品国产露脸久久av麻豆| 日本wwww免费看| 亚洲成人手机| 狠狠精品人妻久久久久久综合| 成人午夜精彩视频在线观看| freevideosex欧美| av国产精品久久久久影院| 三级国产精品片| kizo精华| 国产极品天堂在线| 嫩草影院新地址| 国产精品一二三区在线看| 少妇的逼水好多| 日韩av在线免费看完整版不卡| 自拍欧美九色日韩亚洲蝌蚪91 | 一级,二级,三级黄色视频| 亚洲欧洲日产国产| 国产精品久久久久久精品古装| 熟女av电影| 天天躁夜夜躁狠狠久久av| 欧美精品国产亚洲| 国产精品久久久久久av不卡| 中文字幕人妻熟人妻熟丝袜美| 精品国产一区二区三区久久久樱花| 我的女老师完整版在线观看| 简卡轻食公司| 精品视频人人做人人爽| 国产一区二区三区av在线| 9色porny在线观看| 纯流量卡能插随身wifi吗| 人妻制服诱惑在线中文字幕| 中文天堂在线官网| 色婷婷av一区二区三区视频| 久久人人爽人人爽人人片va| 精品久久国产蜜桃| 久久 成人 亚洲| 亚洲欧美日韩另类电影网站| 乱人伦中国视频| 尾随美女入室| 男女免费视频国产| 精华霜和精华液先用哪个| 有码 亚洲区| av卡一久久| 国产中年淑女户外野战色| 精品亚洲乱码少妇综合久久| 日韩欧美一区视频在线观看 | a级毛片免费高清观看在线播放| 在线观看免费高清a一片| 欧美日韩av久久| 如何舔出高潮| 久久狼人影院| .国产精品久久| 亚洲真实伦在线观看| 精品久久久久久电影网| 色哟哟·www| 一二三四中文在线观看免费高清| 国产精品伦人一区二区| 男人添女人高潮全过程视频| 国产毛片在线视频| 亚洲欧美一区二区三区黑人 | 国产免费又黄又爽又色| 亚洲欧美成人精品一区二区| 黄色一级大片看看| 自拍偷自拍亚洲精品老妇| 一级毛片 在线播放| 亚洲av电影在线观看一区二区三区| 国产免费视频播放在线视频| 下体分泌物呈黄色| 人妻少妇偷人精品九色| 在线亚洲精品国产二区图片欧美 | 简卡轻食公司| 亚洲精品成人av观看孕妇| 日韩视频在线欧美| 亚洲欧美日韩东京热| 亚洲经典国产精华液单| 大陆偷拍与自拍| 国产一区二区三区av在线| 多毛熟女@视频| 高清欧美精品videossex| 国产爽快片一区二区三区| 另类亚洲欧美激情| 天堂中文最新版在线下载| 国产一区二区三区av在线| 97在线视频观看| 男人和女人高潮做爰伦理| 久久久久久人妻| 在线播放无遮挡| 亚洲av欧美aⅴ国产| 亚洲精品色激情综合| 美女福利国产在线| 丝袜在线中文字幕| 麻豆成人av视频| 亚洲欧美精品专区久久| 99九九在线精品视频 | 一级毛片电影观看| 国产精品久久久久久久久免| 免费观看av网站的网址| 色视频在线一区二区三区| 亚洲美女搞黄在线观看| 久久免费观看电影| 男男h啪啪无遮挡| tube8黄色片| 寂寞人妻少妇视频99o| 国产一区二区在线观看av| 9色porny在线观看| 欧美日韩国产mv在线观看视频| 99久久综合免费| 最近2019中文字幕mv第一页| 一个人免费看片子| 日韩在线高清观看一区二区三区| 不卡视频在线观看欧美| 人妻夜夜爽99麻豆av| 成人黄色视频免费在线看| 欧美少妇被猛烈插入视频| 亚洲色图综合在线观看| 乱系列少妇在线播放| 国产精品久久久久久久久免| 观看免费一级毛片| 国产日韩欧美亚洲二区| 免费人成在线观看视频色| 黄色配什么色好看| 22中文网久久字幕| 中国三级夫妇交换| 中文字幕av电影在线播放| 国产69精品久久久久777片| a级一级毛片免费在线观看| 亚洲第一区二区三区不卡| 亚洲av成人精品一二三区| av国产精品久久久久影院| 啦啦啦啦在线视频资源| 97超碰精品成人国产| 久久久亚洲精品成人影院| 有码 亚洲区| a 毛片基地| 国产成人精品一,二区| 欧美三级亚洲精品| 爱豆传媒免费全集在线观看| 2021少妇久久久久久久久久久| 特大巨黑吊av在线直播| 纵有疾风起免费观看全集完整版| 黄色毛片三级朝国网站 | 一级毛片 在线播放| 在线免费观看不下载黄p国产| 亚洲中文av在线| 自拍欧美九色日韩亚洲蝌蚪91 | 丰满人妻一区二区三区视频av| 水蜜桃什么品种好| 国产白丝娇喘喷水9色精品| 国产国拍精品亚洲av在线观看| 啦啦啦啦在线视频资源| 国产熟女欧美一区二区| 色视频在线一区二区三区| 丝袜喷水一区| 日本与韩国留学比较| 欧美日韩国产mv在线观看视频| 国产精品国产av在线观看| 99久久精品热视频| 免费黄色在线免费观看| 人体艺术视频欧美日本| 最近中文字幕高清免费大全6| 在线观看一区二区三区激情| 亚洲av成人精品一二三区| 国产美女午夜福利| 精品国产露脸久久av麻豆| 久久影院123| 又爽又黄a免费视频| 亚洲精品自拍成人| 在线 av 中文字幕| 久久女婷五月综合色啪小说| 黄色怎么调成土黄色| 观看免费一级毛片| 少妇人妻 视频| 国产美女午夜福利| 午夜福利在线观看免费完整高清在| 一级二级三级毛片免费看| 丁香六月天网| 如日韩欧美国产精品一区二区三区 | 免费人妻精品一区二区三区视频| 日本欧美国产在线视频| 国产精品一区二区在线观看99| 伦理电影免费视频| 亚洲三级黄色毛片| 久久久久精品久久久久真实原创| 亚洲精品第二区| 黑丝袜美女国产一区| 偷拍熟女少妇极品色| 久久久国产欧美日韩av| 久久久久精品久久久久真实原创| 妹子高潮喷水视频| 久久久午夜欧美精品| 性色avwww在线观看| 亚洲国产精品999| 免费黄色在线免费观看| 男女免费视频国产| 亚洲va在线va天堂va国产| 伦精品一区二区三区| 午夜激情久久久久久久| 99久久人妻综合| 亚洲av国产av综合av卡| 99久国产av精品国产电影| 日本黄色日本黄色录像| 日韩av免费高清视频| 亚洲图色成人| 91aial.com中文字幕在线观看| 亚洲av免费高清在线观看| www.av在线官网国产| 男女国产视频网站| 精品一品国产午夜福利视频| 亚洲国产av新网站| 亚洲精品国产av蜜桃| 亚洲av国产av综合av卡| 国产精品一区www在线观看| 久久国产精品大桥未久av | 精品国产一区二区三区久久久樱花| 久久综合国产亚洲精品| 成人午夜精彩视频在线观看| 99久久精品热视频| 亚洲精品国产av蜜桃| 亚洲av国产av综合av卡| 黄色视频在线播放观看不卡| 你懂的网址亚洲精品在线观看| 午夜免费鲁丝| 久热这里只有精品99| 日本vs欧美在线观看视频 | 成人国产av品久久久| 久久精品熟女亚洲av麻豆精品| 亚洲av成人精品一二三区| 久久99热6这里只有精品| 两个人的视频大全免费| 国产免费一区二区三区四区乱码| 女性生殖器流出的白浆| 国精品久久久久久国模美| 91精品一卡2卡3卡4卡| 国产精品人妻久久久影院| 国产亚洲一区二区精品| 免费高清在线观看视频在线观看| 亚洲欧美一区二区三区国产| 国产成人精品久久久久久| 毛片一级片免费看久久久久| 欧美日本中文国产一区发布| 女人精品久久久久毛片| 国产淫片久久久久久久久| h视频一区二区三区| 亚洲第一区二区三区不卡| 亚洲怡红院男人天堂| av国产久精品久网站免费入址| 乱人伦中国视频| 国产一区二区三区综合在线观看 | 午夜老司机福利剧场| 亚洲图色成人| 欧美老熟妇乱子伦牲交| 男人舔奶头视频| 免费久久久久久久精品成人欧美视频 | 免费观看性生交大片5| 国产精品久久久久久久久免| 国产一区二区三区综合在线观看 | 精品人妻熟女毛片av久久网站| 日韩伦理黄色片| 大香蕉97超碰在线| 91精品一卡2卡3卡4卡| 亚洲欧洲国产日韩| 国产精品嫩草影院av在线观看| 久久毛片免费看一区二区三区| 日本午夜av视频| 在现免费观看毛片| 欧美老熟妇乱子伦牲交| 亚洲美女视频黄频| 精品人妻熟女av久视频| 成年人免费黄色播放视频 | 青春草视频在线免费观看| 美女国产视频在线观看| 久久精品国产鲁丝片午夜精品| 日本午夜av视频| 国产69精品久久久久777片| 亚洲精品自拍成人| 少妇熟女欧美另类| 91午夜精品亚洲一区二区三区| 我要看黄色一级片免费的| 久久毛片免费看一区二区三区| 国产黄频视频在线观看| 免费观看性生交大片5| 国产精品久久久久久久久免| 岛国毛片在线播放| 色吧在线观看| 2021少妇久久久久久久久久久| 久久这里有精品视频免费| 国产永久视频网站| 国产高清国产精品国产三级| 中文字幕亚洲精品专区| 少妇高潮的动态图| 3wmmmm亚洲av在线观看| 精品少妇内射三级| 乱码一卡2卡4卡精品| 日韩一本色道免费dvd| 亚洲欧美日韩卡通动漫| 亚洲国产成人一精品久久久| 波野结衣二区三区在线| 十分钟在线观看高清视频www | 如何舔出高潮| 色94色欧美一区二区| 国产亚洲91精品色在线| 超碰97精品在线观看| 夜夜看夜夜爽夜夜摸| 欧美性感艳星| av在线老鸭窝| 欧美日韩视频精品一区| 久久国内精品自在自线图片| 高清视频免费观看一区二区| 99热网站在线观看| 国产精品一区www在线观看| 成人免费观看视频高清| 少妇的逼好多水| 亚洲人成网站在线播| 亚洲av欧美aⅴ国产| 91精品国产国语对白视频| 国产午夜精品一二区理论片| www.色视频.com| 黑人巨大精品欧美一区二区蜜桃 | 久久99蜜桃精品久久| 建设人人有责人人尽责人人享有的| 日日啪夜夜撸| 国产亚洲一区二区精品| 精华霜和精华液先用哪个| 丝袜在线中文字幕| 日韩不卡一区二区三区视频在线| 十八禁高潮呻吟视频 | 精品视频人人做人人爽| 激情五月婷婷亚洲| 日日摸夜夜添夜夜添av毛片| 综合色丁香网| 人妻制服诱惑在线中文字幕| 少妇猛男粗大的猛烈进出视频| 日韩熟女老妇一区二区性免费视频| 午夜福利视频精品| 99精国产麻豆久久婷婷| 色婷婷久久久亚洲欧美| 亚洲一区二区三区欧美精品| 妹子高潮喷水视频| 欧美日韩视频高清一区二区三区二| 国产精品国产av在线观看| 噜噜噜噜噜久久久久久91| 在线播放无遮挡| 久久人人爽av亚洲精品天堂| 伦理电影大哥的女人| 91久久精品国产一区二区成人| 国产黄片美女视频| 久久国内精品自在自线图片| 国内少妇人妻偷人精品xxx网站| 国产成人aa在线观看| 亚洲精品日韩在线中文字幕| 亚洲av综合色区一区| 久久国产乱子免费精品| 老熟女久久久| 97超视频在线观看视频| 亚洲丝袜综合中文字幕| 伦理电影免费视频| 十八禁网站网址无遮挡 | 日本黄大片高清| 在线 av 中文字幕| 日日啪夜夜撸| 91精品国产九色| 免费播放大片免费观看视频在线观看| 91成人精品电影| 亚州av有码| 青春草视频在线免费观看| 亚洲欧洲日产国产| 麻豆成人av视频| 亚洲色图综合在线观看| 精品国产露脸久久av麻豆| 99久国产av精品国产电影| 高清午夜精品一区二区三区| 啦啦啦中文免费视频观看日本| 日韩欧美一区视频在线观看 | 热re99久久精品国产66热6| 精品酒店卫生间| √禁漫天堂资源中文www| 亚洲图色成人| 亚洲欧美精品自产自拍| 偷拍熟女少妇极品色| 伊人久久国产一区二区| 亚洲精品视频女| 啦啦啦中文免费视频观看日本| 久久狼人影院| 国产亚洲5aaaaa淫片| 91精品伊人久久大香线蕉| 欧美高清成人免费视频www| 菩萨蛮人人尽说江南好唐韦庄| 在线免费观看不下载黄p国产| 国产精品一区二区三区四区免费观看| 一本—道久久a久久精品蜜桃钙片| 水蜜桃什么品种好| 欧美区成人在线视频| 欧美性感艳星| 国产高清三级在线| 国产极品粉嫩免费观看在线 | 丰满乱子伦码专区| 国产高清三级在线| 精品一区二区免费观看| 欧美日韩亚洲高清精品| 这个男人来自地球电影免费观看 | 九九爱精品视频在线观看| 国产黄频视频在线观看| 亚洲不卡免费看| 国产精品秋霞免费鲁丝片| 久久久久久久国产电影| 亚州av有码| 91成人精品电影| 亚洲国产av新网站| 香蕉精品网在线| 国产黄色免费在线视频| 老司机影院成人| 女性生殖器流出的白浆| 一级毛片电影观看| 欧美精品高潮呻吟av久久| 草草在线视频免费看| 国产精品久久久久久精品古装| 两个人的视频大全免费| 啦啦啦中文免费视频观看日本| 精品一品国产午夜福利视频| 日韩伦理黄色片| 国产69精品久久久久777片| 一级毛片黄色毛片免费观看视频| 亚洲国产精品一区三区| 亚洲三级黄色毛片| 美女国产视频在线观看| av女优亚洲男人天堂| 看十八女毛片水多多多| 久久影院123| 久久狼人影院| 国产精品国产三级国产av玫瑰| 国产精品蜜桃在线观看| 久久热精品热| 欧美 亚洲 国产 日韩一| 少妇高潮的动态图| 日韩成人av中文字幕在线观看| 另类亚洲欧美激情| 国国产精品蜜臀av免费| 免费少妇av软件| 免费不卡的大黄色大毛片视频在线观看| 欧美成人午夜免费资源| 久久精品国产自在天天线| 欧美 日韩 精品 国产| 日本av免费视频播放| 在线观看www视频免费| 中文天堂在线官网| 女性被躁到高潮视频| 国产精品一区二区性色av| 人妻人人澡人人爽人人| 大陆偷拍与自拍| 午夜91福利影院| 人体艺术视频欧美日本| 能在线免费看毛片的网站| 三级国产精品欧美在线观看| 国产视频内射| 国精品久久久久久国模美| 成年美女黄网站色视频大全免费 | 91成人精品电影| 日韩强制内射视频| 夫妻性生交免费视频一级片| 少妇人妻精品综合一区二区| 国产精品秋霞免费鲁丝片| h视频一区二区三区| 久久精品久久精品一区二区三区| 国产成人精品久久久久久| 热re99久久精品国产66热6| 日日撸夜夜添| 中文资源天堂在线| 最近手机中文字幕大全| av在线观看视频网站免费| 日日啪夜夜爽| 在线观看国产h片| 亚洲成人av在线免费| 亚洲欧美日韩另类电影网站| 少妇人妻 视频| 99久久精品一区二区三区| 久久 成人 亚洲| 欧美精品一区二区大全| 18禁动态无遮挡网站| 大片电影免费在线观看免费| 黄色配什么色好看| av女优亚洲男人天堂| 五月伊人婷婷丁香| 日韩精品免费视频一区二区三区 | 日韩不卡一区二区三区视频在线| 日韩欧美一区视频在线观看 | 国产免费一级a男人的天堂| 免费看日本二区| 国产男人的电影天堂91| 99国产精品免费福利视频| 午夜老司机福利剧场| 亚洲真实伦在线观看| 在线观看人妻少妇| 综合色丁香网| 国产精品一区二区性色av| 嘟嘟电影网在线观看| 欧美精品国产亚洲| 精品酒店卫生间| 日韩 亚洲 欧美在线| 亚洲av成人精品一二三区| av视频免费观看在线观看| 亚洲成人手机| 日产精品乱码卡一卡2卡三| av在线app专区| 我要看日韩黄色一级片| 欧美区成人在线视频| 日本黄色片子视频| 777米奇影视久久| 欧美激情极品国产一区二区三区 | 男人舔奶头视频| 国产视频内射| 国产精品福利在线免费观看| 麻豆乱淫一区二区| 免费观看av网站的网址| 99久久精品国产国产毛片| 在线看a的网站| 成年人午夜在线观看视频| 中文天堂在线官网| 一本一本综合久久| 制服丝袜香蕉在线| 人妻人人澡人人爽人人| 色网站视频免费| 能在线免费看毛片的网站| 永久网站在线| av不卡在线播放| 熟女电影av网| 成人无遮挡网站| h日本视频在线播放| 久久韩国三级中文字幕| 性高湖久久久久久久久免费观看| 少妇猛男粗大的猛烈进出视频| 久久久久久久久久久免费av| 久久午夜福利片| 人体艺术视频欧美日本| 99久国产av精品国产电影| 交换朋友夫妻互换小说| 国产中年淑女户外野战色| 只有这里有精品99| 亚洲精品国产av蜜桃| 丰满人妻一区二区三区视频av| av专区在线播放| 免费大片18禁| av有码第一页| 日韩免费高清中文字幕av| 国产永久视频网站| 男人狂女人下面高潮的视频| 菩萨蛮人人尽说江南好唐韦庄| 亚洲三级黄色毛片| 久久影院123| av有码第一页| 日韩亚洲欧美综合| 最新的欧美精品一区二区| 中文字幕制服av| 一级毛片 在线播放| 日韩一区二区视频免费看| 在线观看美女被高潮喷水网站| 国产男女超爽视频在线观看| 国产精品久久久久久精品电影小说| 2022亚洲国产成人精品| 中文字幕免费在线视频6| 久久久久视频综合| 97超碰精品成人国产| 黄色毛片三级朝国网站 | 久久韩国三级中文字幕| 肉色欧美久久久久久久蜜桃| 国产美女午夜福利| 国产午夜精品久久久久久一区二区三区| 中文字幕制服av|