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

    自適應(yīng)設(shè)計空間擴展的高效代理模型氣動優(yōu)化設(shè)計方法

    2018-07-31 06:23:32王超高正紅張偉夏露黃江濤
    航空學(xué)報 2018年7期
    關(guān)鍵詞:加點代理樣本

    王超,高正紅,*,張偉,夏露,黃江濤

    1. 西北工業(yè)大學(xué) 航空學(xué)院,西安 710072 2. 中國空氣動力研究與發(fā)展中心 計算空氣動力研究所,綿陽 621000

    基于代理模型的優(yōu)化(Surrogate-Based Optimization,SBO)方法可以提高優(yōu)化設(shè)計效率,便于實現(xiàn)全局優(yōu)化搜索,因此在工程設(shè)計領(lǐng)域得到廣泛應(yīng)用[1-4]。對于確定的優(yōu)化問題,代理模型利用設(shè)計空間中有限的樣本,通過一定的模型假設(shè)構(gòu)建設(shè)計變量與目標(biāo)函數(shù)的近似關(guān)系,并對設(shè)計空間中的未知點進行預(yù)測。代理模型的預(yù)測精度是影響基于代理模型優(yōu)化效果的關(guān)鍵,工程設(shè)計中常采用具有描述強非線性問題能力的代理模型,如徑向基函數(shù)(RBF)[5]、Kriging[6]、支持向量回歸(SVR)[7]等模型,這類代理模型的預(yù)測精度與設(shè)計空間內(nèi)樣本的密度密切相關(guān)[8]。對于工程優(yōu)化設(shè)計問題,特別是在初步設(shè)計階段,通常需要選擇盡可能大的設(shè)計空間,以獲得具有全局最優(yōu)的設(shè)計結(jié)果。然而,設(shè)計空間的擴大將直接導(dǎo)致樣本在空間中的分布密度急劇下降,嚴(yán)重影響代理模型的精度。增加樣本數(shù)目可以改善代理模型的精度,但是由此會帶來巨大的計算代價,從而大幅度降低代理模型的優(yōu)化效率。因此,如何構(gòu)建盡可能包含全局最優(yōu)設(shè)計點的有效設(shè)計空間,以及如何選擇用于重構(gòu)代理模型的有效樣本是基于代理模型優(yōu)化方法需要解決的關(guān)鍵問題。

    Jones等[9]于1998年提出了基于Kriging模型的高效全局優(yōu)化(Efficient Global Optimization,EGO)方法。EGO方法將優(yōu)化搜索與代理模型加點重構(gòu)相結(jié)合,根據(jù)Kriging模型對目標(biāo)函數(shù)預(yù)測的均值和方差,將目標(biāo)函數(shù)改善的期望(Expected Improvement,EI)作為Kriging模型優(yōu)化的適應(yīng)值函數(shù)。由于EGO方法同時考慮了代理模型的預(yù)測值以及預(yù)測的不確定性,因此具有全局優(yōu)化能力[10]。Jones[11]對Kriging模型不同的加點函數(shù)進行對比研究時指出,EI加點準(zhǔn)則存在收斂緩慢,優(yōu)化效率低的問題。為了進一步提高EGO方法的全局性和高效性,Schonlau和Welch[12]將EI函數(shù)進行了推廣,引入了一個調(diào)節(jié)算法全局和局部搜索能力的參數(shù);Sasena[13]、韓忠華[14-15]等討論了不同的加點函數(shù)和約束處理方法;Soberster[16]、Ginsbourger[17-18]、Liu[19]等發(fā)展了多種并行加點策略。

    盡管各類改進的EGO方法提高了樣本的利用率,但是對于大設(shè)計空間內(nèi)的優(yōu)化問題,EGO方法仍然難以實現(xiàn)全局性與高效性的統(tǒng)一。例如,對于氣動外形優(yōu)化設(shè)計問題,尤其是機翼/翼型等自由曲線曲面類的外形設(shè)計,如果增大設(shè)計空間,優(yōu)化過程中會生成許多流場分離甚至無法計算收斂的“非翼型”外形,給樣本帶來許多“噪聲”;另一方面,大的設(shè)計空間會導(dǎo)致代理模型的精度急劇下降[20],代理優(yōu)化需要大量的加點和更新,從而導(dǎo)致無法在可接受的計算代價內(nèi)獲得理想的設(shè)計結(jié)果。因此,應(yīng)用EGO方法進行氣動優(yōu)化設(shè)計時,為了獲得全局最優(yōu)解,一般采用多輪優(yōu)化的策略,即以上一輪優(yōu)化的結(jié)果作為新的初始外形,重新劃定設(shè)計空間進行優(yōu)化。例如,Zhang等[21]對NACA0012翼型進行三輪優(yōu)化,達(dá)到了全局優(yōu)化的目的。但是,多輪優(yōu)化的方法會導(dǎo)致人工成本增加,優(yōu)化效率降低。

    針對EGO方法在氣動優(yōu)化設(shè)計中面臨的主要問題,本文對基于Kriging代理模型優(yōu)化的加點方法和設(shè)計空間的構(gòu)建問題進行了研究。首先,針對EGO方法本身的收斂性問題,本文采用先全局再局部的優(yōu)化思想,提出了一種混合加點方法,該方法通過引入EI閾值控制不同的加點準(zhǔn)則,以提高EGO方法在確定設(shè)計空間內(nèi)的收斂性。針對設(shè)計空間的構(gòu)建問題,本文對比了擴大設(shè)計變量范圍和多輪優(yōu)化兩種不同的設(shè)計空間構(gòu)建方法,分析了設(shè)計變量范圍對設(shè)計空間大小和樣本密度的影響,進而提出了自適應(yīng)設(shè)計空間擴展的高效代理模型優(yōu)化方法。自適應(yīng)設(shè)計空間擴展的方法在動態(tài)的設(shè)計空間中進行優(yōu)化搜索,只在有潛力的維度擴展設(shè)計空間邊界,通過構(gòu)建自適應(yīng)的設(shè)計空間高效地分配樣本。最后,通過AIAA氣動優(yōu)化設(shè)計討論組(Aerodynamic Design Optimization Discussion Group,ADDG)[22]發(fā)布的標(biāo)準(zhǔn)翼型氣動優(yōu)化設(shè)計算例對本文提出的自適應(yīng)設(shè)計空間擴展的優(yōu)化方法進行了驗證。

    1 基于Kriging模型的混合加點方法

    1.1 Kriging模型

    Kriging模型是一種基于統(tǒng)計理論的插值技術(shù),對于一定數(shù)目的樣本集X=[x1x2…xn]T及其響應(yīng)值y=[y1y2…yn]T,Kriging模型假設(shè)目標(biāo)函數(shù)值與設(shè)計變量之間的真實關(guān)系可以表示為

    y(x)=g(x)+z(x)

    (1)

    g(x)是回歸部分,是用來描述全局趨勢的函數(shù),一般由多項式線性疊加表示,即

    (2)

    式中:fi(x)為低階多項式基函數(shù),可為0階、一階或二階多項式;p為多項式基函數(shù)個數(shù);βi為多項式回歸系數(shù),一般通過廣義最小二乘法求得。

    z(x)是一個靜態(tài)的隨機過程,用來描述真實樣本與趨勢函數(shù)之間的偏差,滿足:

    (3)

    (4)

    (5)

    Kriging模型在未知點xnew的預(yù)測值表示為已知訓(xùn)練樣本的線性組合,即

    (6)

    式中:c為線性加權(quán)系數(shù)。通過預(yù)測方差最小的無偏估計[23]方法,可以求得Kriging模型的預(yù)測值為

    R-1(y-f(xnew)Tβ)

    (7)

    同時得到Kriging模型的預(yù)測方差為

    σ2(xnew)=σz2[1+ζT(FTR-1F)-1ζ-rTR-1r]

    (8)

    式中:ζ=FTR-1r-f,F=[f(x1)f(x2)…f(xn)]T,R為樣本點之間的相關(guān)函數(shù)R(xi,xj)組成的相關(guān)矩陣,r為未知點xnew與樣本點之間的相關(guān)函數(shù)R(xnew,xi)組成的相關(guān)矩陣。

    基于隨機過程的Kriging模型與其他代理模型最大的區(qū)別是它的輸出是一個隨機變量,可以同時給出預(yù)測的均值和方差。根據(jù)Kriging模型預(yù)測的均值和方差,通過一定的加點準(zhǔn)則[11]可以實現(xiàn)高效的優(yōu)化設(shè)計。

    1.2 EGO方法

    基于Kriging模型優(yōu)化最常用的是EGO方法[24-27],又稱EI加點準(zhǔn)則??紤]以下一般性優(yōu)化問題:

    (9)

    (10)

    假設(shè)當(dāng)前樣本中目標(biāo)函數(shù)的最小值為ymin,則Kriging模型在未知點處的預(yù)測值y相對于ymin的改善值為

    I(x)=max(ymin-y,0)

    (11)

    根據(jù)Kriging模型對未知點預(yù)測的概率分布,可以計算任意未知點相對于當(dāng)前最優(yōu)值改善的期望。圖1給出了目標(biāo)函數(shù)改善的幾何意義,圖中紅色部分為目標(biāo)函數(shù)改善的概率,它反映了目標(biāo)函數(shù)改善可能性大小。對目標(biāo)改善值I(x)

    求一階矩,可以得到目標(biāo)改善期望的函數(shù)為

    (12)

    式中:Φ為標(biāo)準(zhǔn)正態(tài)分布函數(shù);φ為標(biāo)準(zhǔn)正態(tài)分布的密度函數(shù)。

    綜上,對于式(9)中的原始優(yōu)化問題,EGO方法的優(yōu)化模型為

    (13)

    EGO方法將EI函數(shù)作為Kriging模型優(yōu)化的適應(yīng)值函數(shù),每次把搜索到EI函數(shù)的最大值點加入到樣本中,更新代理模型,直到滿足終止條件。由于EGO方法可以同時考慮Kriging模型的預(yù)測值及其不確定性,因此具有高效的全局優(yōu)化能力[9]。

    用二維Rastrigin函數(shù)對EGO方法進行測試, 其函數(shù)表達(dá)式為

    f(x,y)=x2+y2-10cos(2πx)-

    10cos(2πy)+20

    (14)

    圖2給出了二維Rastrigin函數(shù)圖,該函數(shù)在x,y∈[-1,1]的設(shè)計空間內(nèi),存在多個局部最優(yōu),其中全局最優(yōu)值為f(0,0)=0。

    采用拉丁超立方取樣(Latin Hypercubic Sampling, LHS)方法[28]生成20個樣本,構(gòu)建初始Kriging模型,然后采用EGO方法對Kriging模型進行優(yōu)化,優(yōu)化算法采用粒子群算法[29],加點50次。圖3給出了目標(biāo)函數(shù)fmin的收斂歷程,圖4給出了每次搜索到EI函數(shù)最大值EImax的收斂歷程。可見,在優(yōu)化后期,EI函數(shù)的最大值變得非常小,其量級已經(jīng)遠(yuǎn)小于目標(biāo)函數(shù)的最小值(如圖5所示),此時EI函數(shù)已經(jīng)無法對代理優(yōu)化進行有效地指引。圖6給出了初始樣本和加點樣本的位置,雖然EGO方法搜索到了最優(yōu)區(qū)域,但是由于EI函數(shù)逐漸失效的問題,導(dǎo)致EGO方法在優(yōu)化的后期基本變?yōu)殡S機搜索,優(yōu)化收斂十分緩慢。

    1.3 混合加點方法

    針對EGO方法中EI函數(shù)在加點后期失效的問題,本文提出混合(Hybrid)加點方法:當(dāng)EI函數(shù)的最大值小于設(shè)定的閾值CEI時,轉(zhuǎn)換為最小預(yù)測值(Minimum Prediction, MP)加點準(zhǔn)則。MP加點準(zhǔn)則直接將代理模型對目標(biāo)函數(shù)的預(yù)測值作為適應(yīng)值函數(shù),具有很強的局部挖掘能力,收斂速度快[11],因此可以彌補EI函數(shù)失效的問題。本文構(gòu)建的Hybrid加點方法的適應(yīng)值函數(shù)表示為

    Fitness=

    (15)

    式中:CEI=0.01ymin,即當(dāng)前樣本中目標(biāo)函數(shù)最小值改善的期望小于1%時,由EI加點準(zhǔn)則轉(zhuǎn)換為MP加點準(zhǔn)則。

    使用Hybrid加點方法對二維Rastrigin函數(shù)重新進行優(yōu)化,圖7給出了加點過程中最大EI值與最小目標(biāo)值的比值,可見當(dāng)加到第34個點時,EImax/ymin<0.01,此時EI準(zhǔn)則轉(zhuǎn)換為MP準(zhǔn)則。

    圖8給出了Hybrid加點過程中不同適應(yīng)值函數(shù)搜索到的最優(yōu)值,其中左軸代表EI函數(shù)最大值,右軸代表預(yù)測目標(biāo)函數(shù)的最小值。

    圖9給出了Hybrid加點方法與EGO方法中目標(biāo)函數(shù)的收斂歷程對比,可見Hybrid方法明顯優(yōu)于EGO方法,EI準(zhǔn)則向MP準(zhǔn)則的轉(zhuǎn)換大幅度加速了收斂進程。圖10給出了Hybrid方法的加點位置,可見Hybrid方法首先采用EI準(zhǔn)則,對設(shè)計空間的最優(yōu)區(qū)域進行全局探索,當(dāng)EI函數(shù)最大值小于CEI時,說明Kriging模型在設(shè)計空間最優(yōu)區(qū)域內(nèi)已經(jīng)具有較高預(yù)測精度,此時轉(zhuǎn)換為MP準(zhǔn)則。由圖10可見,Hybrid方法中的MP樣本點均在全局最優(yōu)值附近,可以對代理模型進行局部精細(xì)化挖掘,提高優(yōu)化的收斂速度。

    2 基于混合加點方法的翼型優(yōu)化設(shè)計

    分別采用EGO方法和Hybrid加點方法對RAE2822翼型進行氣動優(yōu)化設(shè)計。該算例為 ADODG[22]發(fā)布的第2個標(biāo)準(zhǔn)算例:跨聲速黏性條件下RAE2822翼型減阻優(yōu)化設(shè)計。翼型的設(shè)計狀態(tài)為:馬赫數(shù)Ma=0.734, 升力系數(shù)CL=0.824,雷諾數(shù)Re=6.5×105,優(yōu)化目標(biāo)為固定升力系數(shù)下阻力系數(shù)CD最小,約束條件為低頭力矩系數(shù)Cm≤0.92,同時翼型的面積A不減小。RAE2822翼型優(yōu)化的數(shù)學(xué)模型為

    (16)

    采用C型多塊結(jié)構(gòu)網(wǎng)格進行CFD計算,RAE2822翼型網(wǎng)格的拓?fù)浣Y(jié)構(gòu)如圖11所示,c為翼型弦長,網(wǎng)格的具體參數(shù)如表1所示。

    采用CST方法[30]對翼型進行參數(shù)化,翼型上下表面各12個設(shè)計變量,設(shè)計變量數(shù)目為24,給定每個設(shè)計變量的擾動范圍,翼型的設(shè)計空間如圖12所示。

    使用LHS方法在設(shè)計空間生成100個初始樣本,分別采用EGO和Hybrid方法添加200個樣本點,在Hybrid方法中,設(shè)定CEI=0.01CD_min,CD_min為樣本中的最小阻力系數(shù)。采用粒子群算法進行代理模型優(yōu)化搜索,圖13給出了不同優(yōu)化

    Table1ParametersofcomputationalgridofRAE2822airfoil

    參數(shù)數(shù)值遠(yuǎn)場大小50網(wǎng)格規(guī)模601×213物面第一層網(wǎng)格距離5×10-6物面法向網(wǎng)格增長率1.13翼型前緣網(wǎng)格距離0.001翼型后緣網(wǎng)格距離0.001

    方法中翼型阻力系數(shù)的收斂歷程,表2給出了不同方法優(yōu)化的翼型與初始翼型的計算結(jié)果對比,可見Hybrid加點方法的設(shè)計結(jié)果優(yōu)于EGO方法,在滿足約束的條件下,Hybrid方法優(yōu)化翼型的阻力系數(shù)比RAE2822翼型降低了91.8 counts。

    圖14給出了使用Hybrid方法優(yōu)化翼型與RAE2822翼型的外形對比,優(yōu)化翼型頭部變尖,最大厚度位置后移,后緣彎度增大。圖15給出了優(yōu)化翼型與RAE2822翼型的壓力系數(shù)Cp分布對比,圖16和圖17分別給出了優(yōu)化前后翼型流場的壓力云圖。可見,優(yōu)化翼型的前緣吸力峰增大,后加載增加,壓力恢復(fù)平緩,激波得到消除。

    Table2OptimizationresultsofEGOmethodandHybridinfillmethod

    方法CLCDCmARAE28220.8240.020 31-0.092 70.077 9EGO0.8240.011 25-0.091 90.077 9混合加點方法0.8240.011 13-0.091 90.077 9

    3 設(shè)計空間的擴展

    3.1 代理模型加點分析

    為了進一步探討設(shè)計空間對優(yōu)化結(jié)果的影響,首先對RAE2822翼型優(yōu)化算例中的設(shè)計變量進行歸一化處理。圖18通過平行坐標(biāo)系顯示了Hybrid方法中設(shè)計空間各維度的邊界以及在優(yōu)化過程中所有的加點樣本,圖19給出了所有加點樣本對應(yīng)的翼型外形。統(tǒng)計發(fā)現(xiàn),加點樣本中有98.5%的樣本的某一維(或多維)落在設(shè)計空間的邊界,其中最優(yōu)樣本也位于設(shè)計空間的邊界處。由此可見,第2節(jié)RAE2822翼型優(yōu)化算例中的設(shè)計空間不足以尋找全局最優(yōu),因此需要對翼型的設(shè)計空間進行擴展。

    3.2 擴展設(shè)計空間

    設(shè)計空間的構(gòu)建和樣本的選取是基于代理模型優(yōu)化的關(guān)鍵,為了分析不同的設(shè)計空間擴展方法和樣本數(shù)目對設(shè)計結(jié)果的影響,本文對比以下兩種設(shè)計空間擴展方法:

    方法1直接擴展空間。將原始設(shè)計空間每一維的變化范圍擴大一倍,擴展的翼型設(shè)計空間如圖20所示。該方法中,初始樣本數(shù)目為100,加點500次,樣本總數(shù)為原始空間優(yōu)化的2倍。

    方法2兩輪優(yōu)化。以原始設(shè)計空間優(yōu)化的翼型為初始,進行第2輪優(yōu)化,每一維設(shè)計變量的變化范圍與第1輪優(yōu)化保持一致,第2輪優(yōu)化翼型的設(shè)計空間如圖21所示。該方法中,初始樣本數(shù)目為100,加點200次,第2輪優(yōu)化使用的樣本總數(shù)與原始空間優(yōu)化相同。

    使用Opt1代表方法1的優(yōu)化結(jié)果,Opt2代表方法2的優(yōu)化結(jié)果,兩種方法優(yōu)化阻力系數(shù)的收斂對比如圖22所示,優(yōu)化的翼型外形以及壓力分布對比如圖23和圖24所示。可見,直接擴展空間優(yōu)化的翼型與原始空間優(yōu)化的翼型差異很大,翼型上表面仍然存在明顯的激波。表3給出了兩種設(shè)計空間擴展方法優(yōu)化翼型的計算結(jié)果對比,可以發(fā)現(xiàn),直接擴展空間的方法不如原始空間的設(shè)計結(jié)果,而第2輪優(yōu)化的結(jié)果優(yōu)于原始空間的設(shè)計結(jié)果。圖25給出了第2輪優(yōu)化中的加點樣本與設(shè)計空間邊界,可見最優(yōu)樣本落在設(shè)計空間內(nèi)部, 即第2輪優(yōu)化的設(shè)計空間包含了最優(yōu)點。

    該算例表明,設(shè)計空間的構(gòu)建方法對代理模型優(yōu)化設(shè)計結(jié)果有重要的影響,在相同的計算成本下,直接擴大設(shè)計空間的方法不如多輪優(yōu)化的方法。

    方法CFD調(diào)用次數(shù)CDCmA初始翼型0.020 31-0.092 70.077 9原始空間優(yōu)化3000.011 13-0.091 90.077 9第2輪優(yōu)化3000.011 06-0.091 80.077 9直接擴展空間6000.011 54-0.091 10.077 9

    為了進一步分析設(shè)計變量范圍對設(shè)計空間的影響,假設(shè)每一維設(shè)計變量的變化范圍均為R,維數(shù)為D,則設(shè)計空間的大小可以用空間的超體積表示為

    V=RD

    (17)

    設(shè)計空間超體積隨著變量范圍的變化率為

    (18)

    可見,設(shè)計空間的大小隨變量范圍呈D-1次多項式增長。對于RAE2822翼型優(yōu)化算例,如果每一維設(shè)計變量的范圍增加一倍,則設(shè)計空間擴大為原來的224倍。

    Kriging模型的精度與設(shè)計空間樣本的密度密切相關(guān)[20],直接擴大設(shè)計變量范圍會導(dǎo)致設(shè)計空間的暴增,從而帶來樣本嚴(yán)重的稀疏性,代理模型的精度大幅下降,如圖26所示。因此,如果直接擴展設(shè)計空間,代理優(yōu)化難以收斂,無法在可接受的計算代價內(nèi)獲得理想的設(shè)計結(jié)果。

    在兩輪優(yōu)化的方法中,以第1輪優(yōu)化的結(jié)果為初始外形,劃定與第1輪優(yōu)化同樣的變量范圍,則搜索空間只增加了一倍,所以第2輪優(yōu)化采用同樣的樣本得到了理想的設(shè)計結(jié)果。因此,在實際的翼型設(shè)計中,切實可行的方案是進行多輪優(yōu)化而不是一次性給定足夠大的變量范圍。

    4 自適應(yīng)設(shè)計空間擴展方法

    在基于代理模型的氣動優(yōu)化設(shè)計方法中,設(shè)計者一般難以直接構(gòu)建合理的設(shè)計空間,如果設(shè)計空間太大,將導(dǎo)致有限的樣本難以滿足代理模型的精度要求,從而無法獲得理想的設(shè)計結(jié)果;如果設(shè)計空間太小,則需要多輪優(yōu)化,但是多輪優(yōu)化不僅增加了人工成本,而且降低了設(shè)計效率。針對代理模型氣動優(yōu)化中設(shè)計空間的構(gòu)建問題,本文提出自適應(yīng)設(shè)計空間擴展的優(yōu)化方法,該方法不再固定設(shè)計空間的大小,設(shè)計空間在優(yōu)化過程中進行動態(tài)變化,只在有潛力的維度自適應(yīng)擴展設(shè)計變量范圍,通過構(gòu)建自適應(yīng)的設(shè)計空間,提高樣本的利用效率。

    在RAE2822翼型優(yōu)化算例中,加點樣本大部分都落到設(shè)計空間的邊界上,主要有兩方面的原因:一方面是因為設(shè)計空間太小,最優(yōu)解在設(shè)計空間之外,導(dǎo)致優(yōu)化搜索到了邊界;另一方面是因為代理模型在邊界處的精度太差,導(dǎo)致優(yōu)化搜索到了錯誤的最優(yōu)?;谏鲜鰞牲c原因,在進行設(shè)計空間擴展時,需要同時滿足邊界要求和精度要求。本文構(gòu)建的自適應(yīng)設(shè)計空間擴展的代理模型氣動優(yōu)化框架如圖27所示,具體步驟為

    Step1設(shè)定初始的設(shè)計空間,即設(shè)計變量的個數(shù)以及每個設(shè)計變量的變化范圍。

    Step2利用LHS方法在設(shè)計空間選取初始樣本,將初始樣本添加到樣本庫。

    Step3根據(jù)樣本庫確定設(shè)計空間的邊界Xmin和Xmax,構(gòu)建當(dāng)前的設(shè)計空間。其中Xmin為樣本庫中每一維的最小值組成的向量,Xmax為樣本庫中每一維的最大值組成的向量。

    Step4利用樣本庫中的樣本構(gòu)建代理模型。為了提高優(yōu)化效率,每添加10個樣本,重新訓(xùn)練Kriging模型的超參數(shù)。

    Step5根據(jù)代理模型的加點準(zhǔn)則,優(yōu)化搜索得到Xbest。

    Step6CFD評估Xbest,并將Xbest加入到樣本庫。

    Step7判斷Xbest是否滿足設(shè)計空間擴展條件。如果代理模型在Xbest處的預(yù)測值滿足精度要求,并且Xbest的某些維到達(dá)當(dāng)前設(shè)計空間的邊界,則生成自適應(yīng)樣本點Xadaptive,CFD評估Xadaptive,并將Xadaptive添加到樣本庫。

    Xadaptive由Xbest通過在某些維度擴展生成,其作用是擴大設(shè)計空間的邊界。Xadaptive的計算方式為

    (19)

    式中:coefficient為設(shè)計空間擴展系數(shù),本文取0.1。

    Step8判斷樣本總數(shù)N是否達(dá)到設(shè)定的最大樣本數(shù)目Nmax,如果N>Nmax,則終止程序,否則轉(zhuǎn)到Step 3。

    采用本文構(gòu)建的自適應(yīng)設(shè)計空間擴展的優(yōu)化方法對RAE2822翼型重新進行優(yōu)化設(shè)計,初始樣本數(shù)為100,加點200次。以初始設(shè)計空間為參考對樣本每一維設(shè)計變量進行歸一化處理,圖28給出了初始設(shè)計空間與最終設(shè)計空間中設(shè)計變量的邊界對比,對應(yīng)翼型外形邊界如圖29所示。由圖可見,初始設(shè)計空間在某些維度進行了自適應(yīng)地擴展。表4給出了不同方法的優(yōu)化結(jié)果,自適應(yīng)設(shè)計空間的優(yōu)化結(jié)果與第2輪優(yōu)化結(jié)果基本一致,該方法通過一次優(yōu)化便達(dá)到了多輪優(yōu)化的效果,因此大幅提高了氣動優(yōu)化的效率。

    方法CFD調(diào)用次數(shù)CDCmA初始翼型0.020 31-0.092 70.077 9原始空間優(yōu)化3000.011 13-0.091 90.077 9第2輪優(yōu)化3000.011 06-0.091 80.077 9直接擴展空間6000.011 54-0.091 10.077 9自適應(yīng)擴展空間3000.011 05-0.091 90.077 9

    圖30和圖31分別給出了使用不同方法優(yōu)化的翼型外形以及壓力分布對比,可見優(yōu)化的翼型外形以及壓力分布十分相似,說明該算例的“最優(yōu)”翼型集中在初始翼型附近的一個區(qū)域。

    為了對RAE2822翼型優(yōu)化結(jié)果的可靠性進行驗證,將本文的優(yōu)化結(jié)果與文獻(xiàn)[21]中的結(jié)果進行對比,如表5所示??梢园l(fā)現(xiàn),由于網(wǎng)格和CFD求解器的不同,不同方法對初始翼型和優(yōu)化翼型計算的氣動力系數(shù)存在一定的差異,因此難以直接衡量優(yōu)化結(jié)果的優(yōu)劣。從翼型的外形與壓力分布來看,文獻(xiàn)[21]中優(yōu)化翼型的外形以及壓力分布的形態(tài)與本文的結(jié)果十分相似,二者均得到一個典型的跨聲速無激波壓力分布,如圖32和圖33所示。

    與文獻(xiàn)[21]中基于代理模型優(yōu)化方法的計算代價進行對比,本文的CFD計算次數(shù)為300,文獻(xiàn)[21]采用多輪優(yōu)化的方式,CFD計算次數(shù)為592。由此可以看出,本文提出的自適應(yīng)設(shè)計空間擴展方法的優(yōu)化效率得到大幅提升。

    表5 與文獻(xiàn)[21]中RAE2822翼型的優(yōu)化結(jié)果對比Table 5 Comparison with optimization results of RAE2822 airfoil in Ref. [21]

    5 NACA0012翼型優(yōu)化設(shè)計

    為了進一步驗證自適應(yīng)設(shè)計空間擴展優(yōu)化方法的適用性,本節(jié)對NACA0012翼型進行優(yōu)化設(shè)計。該算例為ADODG[22]發(fā)布的第1個標(biāo)準(zhǔn)算例:跨聲速無黏條件下NACA0012翼型的減阻設(shè)計。從文獻(xiàn)[31-36]的結(jié)果可以看出,該算例的優(yōu)化結(jié)果與初始翼型相比,變形尺度非常大,可以更好地考驗自適應(yīng)設(shè)計空間擴展方法的能力。

    該算例的設(shè)計狀態(tài)為:Ma=0.85,迎角α=0°,優(yōu)化目標(biāo)為0°迎角下翼型阻力系數(shù)最小,約束條件為翼型每個位置處的厚度t不減小。NACA0012翼型優(yōu)化的數(shù)學(xué)模型為

    (20)

    由于該算例為0°迎角下對稱翼型的優(yōu)化,本文采用半模外形進行計算和優(yōu)化,全模翼型的阻力為半模外形的2倍。采用O型多塊結(jié)構(gòu)網(wǎng)格進行CFD計算,初始翼型網(wǎng)格的拓?fù)浣Y(jié)構(gòu)如圖34所示,網(wǎng)格參數(shù)的設(shè)置如表6所示。

    采用CST方法對翼型進行參數(shù)化,設(shè)計變量數(shù)目為24,給定每個設(shè)計變量的擾動范圍,翼型的初始設(shè)計空間如圖35所示。

    采用構(gòu)建的自適應(yīng)設(shè)計空間擴展的優(yōu)化方法對NACA0012翼型進行優(yōu)化,并與固定設(shè)計空間的優(yōu)化方法進行對比。初始樣本數(shù)目為100,加點200次。圖36給出了歸一化后初始設(shè)計空間與最終設(shè)計空間的邊界對比,對應(yīng)的翼型外形邊界如圖37所示,可見自適應(yīng)設(shè)計空間的方法在翼型的后緣大幅度地擴展了設(shè)計變量的范圍。

    Table6ParametersofNACA0012airfoilcomputationalgrid

    參數(shù)數(shù)值遠(yuǎn)場大小50網(wǎng)格規(guī)模201×201物面第1層網(wǎng)格距離5×10-4物面法向網(wǎng)格增長率1.05翼型前緣網(wǎng)格距離0.001翼型后緣網(wǎng)格距離0.001

    圖38給出了自適應(yīng)設(shè)計空間與固定設(shè)計空間方法優(yōu)化阻力系數(shù)的收斂歷程,在相同的計算花費下,自適應(yīng)設(shè)計空間的方法獲得了更優(yōu)的結(jié)果。表7給出了不同設(shè)計空間下翼型優(yōu)化結(jié)果對比,自適應(yīng)設(shè)計空間優(yōu)化的翼型阻力系數(shù)為55 counts,明顯小于固定設(shè)計空間的優(yōu)化結(jié)果。

    圖39和圖40分別給出了固定設(shè)計空間與自適應(yīng)設(shè)計空間優(yōu)化的翼型外形以及壓力分布對比??梢?,自適應(yīng)設(shè)計空間優(yōu)化的翼型前緣和后緣更加飽滿, 翼型前緣半徑增加,后緣變形尺度增大;從壓力分布來看, 自適應(yīng)設(shè)計空間優(yōu)化的翼型前緣的吸力峰更高,翼型后緣通過連續(xù)地膨脹和壓縮使壓力得到緩和恢復(fù),從而減弱了激波強度,降低了激波阻力。

    方法CLCDNACA001200.047 1固定設(shè)計空間00.017 5自適應(yīng)設(shè)計空間00.005 5

    從總壓損失的角度進一步分析優(yōu)化翼型阻力系數(shù)減小的原因。在可壓縮條件下,總壓與靜壓的關(guān)系以及總壓系數(shù)表示為

    (21)

    (22)

    式中:p0為總壓;p為靜壓;Ma為當(dāng)?shù)伛R赫數(shù);γ為氣體比熱比;p∞為遠(yuǎn)場靜壓;ρ∞為遠(yuǎn)場密度;v∞為自由來流速度。

    NACA0012翼型以及兩種不同方法優(yōu)化翼型的流場總壓系數(shù)云圖如圖41~圖43所示,截取翼型表面的總壓系數(shù)如圖44所示。由圖可見,NACA0012翼型激波最強,在激波位置處總壓下降最明顯。圖45給出了翼型后緣一倍弦長處(x/c=2)截面總壓系數(shù)的對比,可以發(fā)現(xiàn),該截面處自適應(yīng)設(shè)計空間優(yōu)化翼型的總壓相對來流總壓損失最小,即能量損失最小,因此該翼型的阻力系數(shù)最小。

    將本文的優(yōu)化結(jié)果與文獻(xiàn)[21]中結(jié)果進行對比,圖46和圖47分別給出了文獻(xiàn)[21]通過多輪優(yōu)化方法得到的翼型外形和壓力分布,該結(jié)果與本文的優(yōu)化結(jié)果具有很強的一致性。表8給出了兩種不同方法優(yōu)化的結(jié)果對比,可以發(fā)現(xiàn),本文的方法通過更少的CFD計算代價獲得了阻力更小的翼型。

    該算例進一步證實,自適應(yīng)設(shè)計空間擴展優(yōu)化方法對于大尺度變形的氣動優(yōu)化問題具有很好的適應(yīng)性。

    表8 與文獻(xiàn)[21]中NACA0012翼型的優(yōu)化結(jié)果對比

    6 結(jié) 論

    本文針對基于Kriging模型氣動優(yōu)化中的加點方法以及設(shè)計空間的構(gòu)建問題進行了研究,提出了混合加點方法和自適應(yīng)設(shè)計空間擴展的優(yōu)化方法,通過ADODG[22]標(biāo)準(zhǔn)翼型優(yōu)化算例對文中提出的方法進行了驗證。

    1) 文中構(gòu)建的混合加點方法通過引入EI閾值控制EI和MP加點準(zhǔn)則,利用先全局再局部的優(yōu)化思想,可以提高EGO算法在確定設(shè)計空間內(nèi)的收斂性。

    2) 設(shè)計空間的構(gòu)建方法對基于代理模型氣動優(yōu)化的結(jié)果有重要的影響,如果設(shè)計空間的尺度偏于保守,則搜索空間受到限制,從而失去全局優(yōu)化能力;如果設(shè)計空間的尺度偏于激進,將導(dǎo)致有限的樣本難以滿足代理模型的精度要求,從而無法獲得理想的設(shè)計結(jié)果。

    3) 本文構(gòu)建的自適應(yīng)設(shè)計空間擴展的方法打破了傳統(tǒng)代理模型優(yōu)化依賴人工經(jīng)驗給定設(shè)計空間的局限,通過構(gòu)建自適應(yīng)的設(shè)計空間,實現(xiàn)了樣本的高效配置,很好地兼顧了氣動優(yōu)化的高效性和全局性。通過ADODG標(biāo)準(zhǔn)翼型氣動優(yōu)化設(shè)計算例證實,自適應(yīng)設(shè)計空間擴展方法可以大幅提高氣動優(yōu)化效率。

    猜你喜歡
    加點代理樣本
    給地球加點綠
    用樣本估計總體復(fù)習(xí)點撥
    給電影加點特效
    代理圣誕老人
    代理手金寶 生意特別好
    推動醫(yī)改的“直銷樣本”
    泡腳可以加點藥
    隨機微分方程的樣本Lyapunov二次型估計
    復(fù)仇代理烏龜君
    為靦腆膽怯加點“料”,秀出你的不同凡響
    欧美少妇被猛烈插入视频| 免费高清在线观看日韩| 国产精品国产三级国产专区5o| 日日爽夜夜爽网站| 国产精品亚洲av一区麻豆| 麻豆国产av国片精品| 日本精品一区二区三区蜜桃| 国产精品国产三级国产专区5o| 777米奇影视久久| 中亚洲国语对白在线视频| 国产91精品成人一区二区三区 | 黑人欧美特级aaaaaa片| 午夜精品久久久久久毛片777| 99精国产麻豆久久婷婷| 国产主播在线观看一区二区| 美女国产高潮福利片在线看| 午夜精品久久久久久毛片777| 亚洲成av片中文字幕在线观看| 自线自在国产av| 岛国在线观看网站| 丰满人妻熟妇乱又伦精品不卡| 亚洲,欧美精品.| www.熟女人妻精品国产| 国产国语露脸激情在线看| av有码第一页| 美女午夜性视频免费| 91精品伊人久久大香线蕉| 视频区图区小说| 少妇裸体淫交视频免费看高清 | 欧美老熟妇乱子伦牲交| 亚洲少妇的诱惑av| 午夜91福利影院| 国产极品粉嫩免费观看在线| 亚洲视频免费观看视频| 亚洲精品一卡2卡三卡4卡5卡 | 法律面前人人平等表现在哪些方面 | 精品少妇黑人巨大在线播放| 免费在线观看影片大全网站| 两性午夜刺激爽爽歪歪视频在线观看 | 久久久精品94久久精品| 女警被强在线播放| 日韩熟女老妇一区二区性免费视频| 亚洲精品国产av蜜桃| 国产精品免费视频内射| 女性被躁到高潮视频| 久久久久久久精品精品| 伊人亚洲综合成人网| 欧美变态另类bdsm刘玥| 国产在线一区二区三区精| 我要看黄色一级片免费的| 91大片在线观看| 国产极品粉嫩免费观看在线| 国产精品久久久久久人妻精品电影 | 精品第一国产精品| 亚洲一卡2卡3卡4卡5卡精品中文| 日韩,欧美,国产一区二区三区| 18禁观看日本| 日韩 亚洲 欧美在线| 亚洲va日本ⅴa欧美va伊人久久 | 久久久精品国产亚洲av高清涩受| 久久精品国产亚洲av高清一级| 91精品三级在线观看| 男人添女人高潮全过程视频| 亚洲国产毛片av蜜桃av| 亚洲精品一卡2卡三卡4卡5卡 | 又黄又粗又硬又大视频| 亚洲国产欧美在线一区| 久久性视频一级片| 久久精品人人爽人人爽视色| 他把我摸到了高潮在线观看 | 亚洲欧美清纯卡通| 亚洲精品av麻豆狂野| 国产成人av教育| 90打野战视频偷拍视频| 精品国产一区二区三区四区第35| 天堂俺去俺来也www色官网| 一本大道久久a久久精品| 精品一区二区三区av网在线观看 | 在线天堂中文资源库| 免费少妇av软件| www.自偷自拍.com| 日韩一区二区三区影片| 大片免费播放器 马上看| 精品久久久久久电影网| 中亚洲国语对白在线视频| 成人亚洲精品一区在线观看| 久久精品国产亚洲av香蕉五月 | 男女下面插进去视频免费观看| 在线精品无人区一区二区三| 男女无遮挡免费网站观看| 黄网站色视频无遮挡免费观看| 免费在线观看日本一区| 久久精品国产亚洲av高清一级| 高清视频免费观看一区二区| 午夜影院在线不卡| 国产一级毛片在线| 秋霞在线观看毛片| 久久久久网色| 咕卡用的链子| 热re99久久国产66热| 欧美一级毛片孕妇| 亚洲美女黄色视频免费看| 肉色欧美久久久久久久蜜桃| 成人亚洲精品一区在线观看| 日本精品一区二区三区蜜桃| tube8黄色片| 免费在线观看影片大全网站| 久久久精品区二区三区| 高清欧美精品videossex| 中国国产av一级| 狂野欧美激情性xxxx| 久久精品亚洲av国产电影网| 狠狠狠狠99中文字幕| 午夜福利在线免费观看网站| 国产精品九九99| 一本一本久久a久久精品综合妖精| tube8黄色片| 国产亚洲av片在线观看秒播厂| 欧美精品高潮呻吟av久久| 下体分泌物呈黄色| 黄色a级毛片大全视频| 国产成人免费观看mmmm| 欧美日韩一级在线毛片| 一级片免费观看大全| 在线看a的网站| 免费高清在线观看日韩| 国产精品二区激情视频| 精品一品国产午夜福利视频| 欧美精品高潮呻吟av久久| 欧美日韩黄片免| 一级毛片精品| 丰满迷人的少妇在线观看| 欧美少妇被猛烈插入视频| 首页视频小说图片口味搜索| 精品久久久久久电影网| 99热国产这里只有精品6| 青草久久国产| av超薄肉色丝袜交足视频| 18在线观看网站| 满18在线观看网站| 国产精品九九99| 青春草视频在线免费观看| 国产在线观看jvid| 国产精品久久久av美女十八| 精品视频人人做人人爽| 国产又爽黄色视频| 欧美+亚洲+日韩+国产| 亚洲专区中文字幕在线| 狠狠精品人妻久久久久久综合| 动漫黄色视频在线观看| 婷婷成人精品国产| 精品少妇一区二区三区视频日本电影| 欧美激情久久久久久爽电影 | 丝瓜视频免费看黄片| 伊人久久大香线蕉亚洲五| 悠悠久久av| 悠悠久久av| 自线自在国产av| 久久精品国产a三级三级三级| 日韩制服丝袜自拍偷拍| 亚洲国产成人一精品久久久| 国产精品99久久99久久久不卡| 成年人午夜在线观看视频| 亚洲午夜精品一区,二区,三区| 老司机亚洲免费影院| 国产精品久久久久成人av| 在线十欧美十亚洲十日本专区| 欧美乱码精品一区二区三区| 亚洲欧美一区二区三区久久| 亚洲成人手机| 精品高清国产在线一区| 精品国内亚洲2022精品成人 | 黄片播放在线免费| 国产97色在线日韩免费| 自拍欧美九色日韩亚洲蝌蚪91| 国产精品亚洲av一区麻豆| 婷婷丁香在线五月| 亚洲精品国产一区二区精华液| 法律面前人人平等表现在哪些方面 | 国产成人av教育| 亚洲激情五月婷婷啪啪| 这个男人来自地球电影免费观看| 欧美少妇被猛烈插入视频| av天堂久久9| 精品少妇久久久久久888优播| av在线播放精品| 亚洲av成人一区二区三| 咕卡用的链子| 一边摸一边做爽爽视频免费| 97人妻天天添夜夜摸| 亚洲精品一卡2卡三卡4卡5卡 | 丝袜在线中文字幕| 每晚都被弄得嗷嗷叫到高潮| 国产欧美日韩一区二区三区在线| 黑人巨大精品欧美一区二区蜜桃| 精品亚洲乱码少妇综合久久| 免费黄频网站在线观看国产| 日本五十路高清| 一个人免费看片子| 国产xxxxx性猛交| 国产免费av片在线观看野外av| 精品一品国产午夜福利视频| 国产av一区二区精品久久| 在线永久观看黄色视频| 美女高潮喷水抽搐中文字幕| 国产三级黄色录像| 亚洲九九香蕉| 亚洲国产精品一区二区三区在线| 日韩人妻精品一区2区三区| 18在线观看网站| e午夜精品久久久久久久| 夜夜夜夜夜久久久久| 大片电影免费在线观看免费| 精品一区在线观看国产| 成人18禁高潮啪啪吃奶动态图| 窝窝影院91人妻| 精品久久久精品久久久| 黄色视频在线播放观看不卡| 国产亚洲一区二区精品| 久久久久久久大尺度免费视频| 午夜福利视频在线观看免费| 80岁老熟妇乱子伦牲交| 老汉色av国产亚洲站长工具| 国产精品.久久久| 成在线人永久免费视频| 国产伦人伦偷精品视频| 亚洲精品自拍成人| 亚洲人成电影免费在线| 十分钟在线观看高清视频www| 一级片'在线观看视频| 久久人人爽av亚洲精品天堂| 午夜91福利影院| 一边摸一边抽搐一进一出视频| 久久久精品区二区三区| 日韩制服丝袜自拍偷拍| 久久亚洲精品不卡| 1024香蕉在线观看| 国产一区二区三区综合在线观看| 国产免费现黄频在线看| 永久免费av网站大全| 在线观看免费午夜福利视频| 婷婷丁香在线五月| 亚洲欧美清纯卡通| 亚洲 国产 在线| 亚洲精品久久成人aⅴ小说| 亚洲精品美女久久av网站| 肉色欧美久久久久久久蜜桃| 久久性视频一级片| 99热网站在线观看| 欧美久久黑人一区二区| 黑人操中国人逼视频| 国产精品一区二区精品视频观看| 久久久精品区二区三区| 久久人人97超碰香蕉20202| 精品国产一区二区三区四区第35| 中文字幕人妻丝袜一区二区| 成在线人永久免费视频| 欧美亚洲 丝袜 人妻 在线| 久久精品久久久久久噜噜老黄| √禁漫天堂资源中文www| 91精品伊人久久大香线蕉| 亚洲午夜精品一区,二区,三区| 天堂俺去俺来也www色官网| www.自偷自拍.com| 国产淫语在线视频| 一本一本久久a久久精品综合妖精| 90打野战视频偷拍视频| 男女边摸边吃奶| 国产成人av教育| 久久久久精品国产欧美久久久 | a级毛片黄视频| 欧美亚洲日本最大视频资源| 搡老岳熟女国产| 久久久国产精品麻豆| 69av精品久久久久久 | 俄罗斯特黄特色一大片| 一级,二级,三级黄色视频| 免费观看a级毛片全部| 午夜久久久在线观看| 亚洲av成人一区二区三| 黑人巨大精品欧美一区二区mp4| 999久久久国产精品视频| 精品一品国产午夜福利视频| 国产免费现黄频在线看| 亚洲美女黄色视频免费看| 最新在线观看一区二区三区| 淫妇啪啪啪对白视频 | 国产无遮挡羞羞视频在线观看| 久久久欧美国产精品| 99国产综合亚洲精品| 制服诱惑二区| 午夜成年电影在线免费观看| 99九九在线精品视频| 久久久久久免费高清国产稀缺| 欧美精品高潮呻吟av久久| 五月开心婷婷网| 成年人黄色毛片网站| 日韩中文字幕视频在线看片| 国产日韩一区二区三区精品不卡| av电影中文网址| 亚洲欧美激情在线| 亚洲国产毛片av蜜桃av| 黄色视频不卡| 亚洲av国产av综合av卡| 欧美黄色淫秽网站| 久久久久久久精品精品| 日韩 亚洲 欧美在线| 新久久久久国产一级毛片| 成年人免费黄色播放视频| 亚洲美女黄色视频免费看| 欧美一级毛片孕妇| 国产1区2区3区精品| 婷婷成人精品国产| 精品一品国产午夜福利视频| 国产av精品麻豆| 亚洲人成电影免费在线| 国产精品一区二区在线不卡| 久久精品亚洲熟妇少妇任你| 国产在线一区二区三区精| 99精品欧美一区二区三区四区| 欧美精品一区二区大全| 国产一区二区在线观看av| 亚洲精品久久午夜乱码| 精品福利观看| 女人高潮潮喷娇喘18禁视频| 亚洲第一青青草原| 精品免费久久久久久久清纯 | 欧美精品亚洲一区二区| 在线观看人妻少妇| 日韩,欧美,国产一区二区三区| 超色免费av| 免费久久久久久久精品成人欧美视频| 丰满少妇做爰视频| 男人舔女人的私密视频| 亚洲精品美女久久久久99蜜臀| a级毛片在线看网站| 少妇的丰满在线观看| 脱女人内裤的视频| 国产亚洲一区二区精品| av视频免费观看在线观看| 51午夜福利影视在线观看| 99热网站在线观看| 午夜日韩欧美国产| 亚洲精品中文字幕在线视频| 久久亚洲国产成人精品v| 视频在线观看一区二区三区| 午夜激情久久久久久久| 免费在线观看视频国产中文字幕亚洲 | 国产欧美日韩一区二区三区在线| 国产精品国产av在线观看| 国产成人精品无人区| 亚洲成国产人片在线观看| 亚洲欧美一区二区三区久久| 欧美av亚洲av综合av国产av| 男女免费视频国产| 可以免费在线观看a视频的电影网站| 夫妻午夜视频| 美女国产高潮福利片在线看| 美女午夜性视频免费| 最近最新中文字幕大全免费视频| 亚洲第一欧美日韩一区二区三区 | 男女边摸边吃奶| 国产成人啪精品午夜网站| 免费在线观看视频国产中文字幕亚洲 | av欧美777| 大型av网站在线播放| av又黄又爽大尺度在线免费看| 王馨瑶露胸无遮挡在线观看| 午夜精品国产一区二区电影| 啦啦啦 在线观看视频| 成人黄色视频免费在线看| 久久综合国产亚洲精品| 国产无遮挡羞羞视频在线观看| 99精品欧美一区二区三区四区| 国产一区二区三区av在线| 免费av中文字幕在线| 色老头精品视频在线观看| 男女无遮挡免费网站观看| 久久人妻福利社区极品人妻图片| 久久九九热精品免费| 黑人欧美特级aaaaaa片| 婷婷色av中文字幕| 成年人免费黄色播放视频| 在线观看免费午夜福利视频| 一区二区三区四区激情视频| 久久久久久人人人人人| 女性生殖器流出的白浆| 国产成人av教育| 无限看片的www在线观看| 精品亚洲成国产av| 黄频高清免费视频| av又黄又爽大尺度在线免费看| 夜夜骑夜夜射夜夜干| 一个人免费看片子| 亚洲av日韩在线播放| 777米奇影视久久| 五月天丁香电影| 国产av一区二区精品久久| 色精品久久人妻99蜜桃| 国产高清videossex| 日韩 亚洲 欧美在线| 男人添女人高潮全过程视频| 最近最新免费中文字幕在线| 国产成人啪精品午夜网站| 97人妻天天添夜夜摸| 精品亚洲乱码少妇综合久久| 视频区欧美日本亚洲| 99精品久久久久人妻精品| 99国产精品99久久久久| 90打野战视频偷拍视频| 在线av久久热| 色婷婷久久久亚洲欧美| avwww免费| 久久精品亚洲av国产电影网| 日日夜夜操网爽| 乱人伦中国视频| 国产一区二区三区在线臀色熟女 | 午夜福利在线观看吧| 亚洲综合色网址| 一本大道久久a久久精品| 日韩视频在线欧美| 纵有疾风起免费观看全集完整版| 成年女人毛片免费观看观看9 | 少妇粗大呻吟视频| 久久久国产精品麻豆| 亚洲激情五月婷婷啪啪| 2018国产大陆天天弄谢| 97精品久久久久久久久久精品| 国产一卡二卡三卡精品| 欧美人与性动交α欧美精品济南到| 看免费av毛片| 午夜影院在线不卡| 天天操日日干夜夜撸| 亚洲国产精品一区三区| 999久久久精品免费观看国产| 国产又爽黄色视频| 永久免费av网站大全| 精品一区在线观看国产| 免费在线观看影片大全网站| 国产精品一区二区在线观看99| 亚洲成国产人片在线观看| 国产精品1区2区在线观看. | 人成视频在线观看免费观看| 久久人人爽av亚洲精品天堂| 99久久人妻综合| 男女之事视频高清在线观看| 午夜久久久在线观看| 曰老女人黄片| 精品乱码久久久久久99久播| 亚洲天堂av无毛| 亚洲精品国产av成人精品| 成年动漫av网址| 母亲3免费完整高清在线观看| 亚洲欧美日韩高清在线视频 | 老汉色∧v一级毛片| 国产有黄有色有爽视频| 夫妻午夜视频| 悠悠久久av| 狠狠狠狠99中文字幕| 91大片在线观看| 男人舔女人的私密视频| 欧美精品一区二区大全| 欧美大码av| 国产深夜福利视频在线观看| 狠狠狠狠99中文字幕| 中国国产av一级| 精品人妻1区二区| av网站免费在线观看视频| 欧美乱码精品一区二区三区| 欧美日本中文国产一区发布| 午夜激情久久久久久久| 精品人妻1区二区| 18在线观看网站| 男男h啪啪无遮挡| 亚洲国产看品久久| 最黄视频免费看| 老司机亚洲免费影院| 一本一本久久a久久精品综合妖精| 亚洲av男天堂| 男人爽女人下面视频在线观看| 十八禁网站免费在线| 女人精品久久久久毛片| 麻豆乱淫一区二区| 女人久久www免费人成看片| 国产精品欧美亚洲77777| 久久久久久久国产电影| 久久久久久久久久久久大奶| 丁香六月欧美| 男女无遮挡免费网站观看| 国产极品粉嫩免费观看在线| 久久久久网色| 亚洲精品美女久久av网站| 大码成人一级视频| 久久久久久人人人人人| 高清黄色对白视频在线免费看| 亚洲 欧美一区二区三区| 美女午夜性视频免费| 视频区图区小说| 久久久久国产精品人妻一区二区| 亚洲av日韩精品久久久久久密| 男女国产视频网站| 日韩欧美国产一区二区入口| 男女下面插进去视频免费观看| 一区二区三区四区激情视频| 在线观看免费视频网站a站| 99精品久久久久人妻精品| 精品福利观看| 国产男女超爽视频在线观看| 久久综合国产亚洲精品| 精品一区二区三卡| 国产成人一区二区三区免费视频网站| 精品少妇一区二区三区视频日本电影| 黄色视频,在线免费观看| 午夜免费鲁丝| 两个人免费观看高清视频| 亚洲av男天堂| 女人被躁到高潮嗷嗷叫费观| av超薄肉色丝袜交足视频| 欧美日韩国产mv在线观看视频| 久久天堂一区二区三区四区| 亚洲精品乱久久久久久| 国产成人欧美| 久久免费观看电影| 一本大道久久a久久精品| www.999成人在线观看| 色综合欧美亚洲国产小说| 成年av动漫网址| 欧美精品一区二区大全| 国产免费视频播放在线视频| 热99久久久久精品小说推荐| 狂野欧美激情性xxxx| 国产深夜福利视频在线观看| 最新的欧美精品一区二区| 老司机深夜福利视频在线观看 | 日本猛色少妇xxxxx猛交久久| 18禁观看日本| 亚洲人成电影观看| 99国产精品一区二区三区| 成人免费观看视频高清| 久久国产精品人妻蜜桃| 免费观看a级毛片全部| 97人妻天天添夜夜摸| 99精国产麻豆久久婷婷| 看免费av毛片| 考比视频在线观看| 精品欧美一区二区三区在线| 人人妻,人人澡人人爽秒播| 午夜精品久久久久久毛片777| 国产日韩欧美视频二区| 久久久久久久久久久久大奶| 两性午夜刺激爽爽歪歪视频在线观看 | 国产精品一区二区在线观看99| 我的亚洲天堂| 国产91精品成人一区二区三区 | 国产精品亚洲av一区麻豆| 好男人电影高清在线观看| 亚洲人成电影观看| 亚洲欧美精品综合一区二区三区| 午夜福利在线免费观看网站| 中文字幕av电影在线播放| 精品少妇黑人巨大在线播放| 欧美精品av麻豆av| 精品人妻在线不人妻| 亚洲成国产人片在线观看| 亚洲精品国产av蜜桃| 成年人黄色毛片网站| 精品少妇一区二区三区视频日本电影| 久久人妻福利社区极品人妻图片| 日日摸夜夜添夜夜添小说| 中文字幕人妻丝袜一区二区| 国产日韩一区二区三区精品不卡| 日韩欧美一区二区三区在线观看 | 啦啦啦视频在线资源免费观看| 人妻一区二区av| 青春草亚洲视频在线观看| 久久毛片免费看一区二区三区| 欧美日韩亚洲国产一区二区在线观看 | 婷婷色av中文字幕| 国产成人影院久久av| 老司机在亚洲福利影院| 亚洲欧美一区二区三区黑人| netflix在线观看网站| 啦啦啦 在线观看视频| 久久青草综合色| www.999成人在线观看| 人人澡人人妻人| 人人妻人人爽人人添夜夜欢视频| 中国国产av一级| 国产91精品成人一区二区三区 | 色综合欧美亚洲国产小说| 亚洲精品一区蜜桃| 日韩制服骚丝袜av| 国产精品影院久久| 极品人妻少妇av视频| 亚洲少妇的诱惑av| 女人精品久久久久毛片| 国产日韩欧美亚洲二区| 99久久精品国产亚洲精品| 精品人妻1区二区| 99久久国产精品久久久| 男人舔女人的私密视频| 黄色片一级片一级黄色片| 成年av动漫网址| 人人妻人人澡人人看| 999久久久国产精品视频| 人人妻人人澡人人看| 午夜视频精品福利| 中文字幕人妻熟女乱码| 好男人电影高清在线观看| 精品久久久久久电影网| 国产野战对白在线观看| 各种免费的搞黄视频| 国产成人系列免费观看| 午夜激情久久久久久久| 两个人免费观看高清视频|