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

    跨聲速運(yùn)輸機(jī)機(jī)翼氣動(dòng)/結(jié)構(gòu)優(yōu)化平臺(tái)AeroStruct的發(fā)展及應(yīng)用

    2022-05-29 08:01:27張科施凌圣博韓忠華
    航空科學(xué)技術(shù) 2022年4期
    關(guān)鍵詞:氣動(dòng)結(jié)構(gòu)設(shè)計(jì)

    張科施 凌圣博 韓忠華

    摘要:跨聲速運(yùn)輸類(lèi)飛機(jī)由于機(jī)翼展弦比大,使得氣動(dòng)加載下機(jī)翼變形量大,氣動(dòng)/結(jié)構(gòu)耦合效應(yīng)明顯,因而開(kāi)展耦合優(yōu)化設(shè)計(jì)對(duì)提升飛機(jī)綜合性能具有重要意義。本文研發(fā)了一種對(duì)該類(lèi)機(jī)翼進(jìn)行變可信度氣動(dòng)/結(jié)構(gòu)耦合優(yōu)化設(shè)計(jì)的自主軟件平臺(tái)(AeroStruct),氣動(dòng)數(shù)值模擬可選用N-S/Euler/全速勢(shì)方程,結(jié)構(gòu)數(shù)值模擬集成了結(jié)構(gòu)有限元分析平臺(tái),氣動(dòng)/結(jié)構(gòu)耦合分析采用徑向基函數(shù)插值方法進(jìn)行數(shù)據(jù)傳遞,在自主開(kāi)發(fā)的SurroOpt代理優(yōu)化軟件基礎(chǔ)上集成分析模塊,實(shí)現(xiàn)了考慮靜氣動(dòng)彈性效應(yīng)的機(jī)翼氣動(dòng)/結(jié)構(gòu)耦合優(yōu)化設(shè)計(jì)。對(duì)uCRM-9標(biāo)模機(jī)翼進(jìn)行了包含76個(gè)設(shè)計(jì)變量的減重優(yōu)化設(shè)計(jì),最優(yōu)方案比基準(zhǔn)機(jī)翼結(jié)構(gòu)質(zhì)量減小約13%,優(yōu)化效果顯著。

    關(guān)鍵詞:機(jī)翼設(shè)計(jì);多學(xué)科優(yōu)化設(shè)計(jì);氣動(dòng)/結(jié)構(gòu)耦合數(shù)值模擬;氣動(dòng)設(shè)計(jì);結(jié)構(gòu)設(shè)計(jì)

    中圖分類(lèi)號(hào):V221文獻(xiàn)標(biāo)識(shí)碼:ADOI:10.19452/j.issn1007-5453.2022.04.007

    基金項(xiàng)目:航空科學(xué)基金(2019ZA053004)

    先進(jìn)機(jī)翼氣動(dòng)/結(jié)構(gòu)耦合優(yōu)化設(shè)計(jì)技術(shù),采用高可信度數(shù)值模擬并充分考慮氣動(dòng)與結(jié)構(gòu)學(xué)科間的強(qiáng)耦合關(guān)系,是未來(lái)航空飛行器提升綜合性能的一項(xiàng)核心關(guān)鍵技術(shù)。尤其對(duì)具有大變形特征的跨聲速運(yùn)輸類(lèi)飛機(jī)機(jī)翼設(shè)計(jì)來(lái)說(shuō),氣動(dòng)/結(jié)構(gòu)耦合多學(xué)科優(yōu)化設(shè)計(jì)將成為未來(lái)必然選擇。

    機(jī)翼氣動(dòng)/結(jié)構(gòu)耦合優(yōu)化的研究始于20世紀(jì)70年代[1]。21世紀(jì)以來(lái),隨著多學(xué)科優(yōu)化設(shè)計(jì)(MDO)技術(shù)和高性能計(jì)算機(jī)并行計(jì)算能力的不斷進(jìn)步,在機(jī)翼氣動(dòng)/結(jié)構(gòu)耦合優(yōu)化研究領(lǐng)域,高可信度耦合數(shù)值模擬、大規(guī)模設(shè)計(jì)變量?jī)?yōu)化、高效全局優(yōu)化算法等方面的研究不斷取得突破。

    首先,高可信度氣動(dòng)/結(jié)構(gòu)耦合數(shù)值模擬研究不斷深入。高可信度耦合數(shù)值模擬是關(guān)系機(jī)翼氣動(dòng)/結(jié)構(gòu)耦合優(yōu)化設(shè)計(jì)技術(shù)是否具有工程實(shí)用性的一項(xiàng)核心關(guān)鍵技術(shù)。國(guó)內(nèi)外研究機(jī)構(gòu)在高可信度的“Navier-Stokes方程+有限元分析”數(shù)值模擬[2-6]以及高/中可信度的非對(duì)稱(chēng)數(shù)值模擬[7-10]方面取得重要進(jìn)展,并在優(yōu)化中考慮學(xué)科間相互影響——靜氣動(dòng)彈性[2-10],進(jìn)而考慮顫振[11]。其次,高效Adjoint敏度分析方法的發(fā)展,使得“維數(shù)災(zāi)難”問(wèn)題得到緩解。機(jī)翼氣動(dòng)/結(jié)構(gòu)優(yōu)化涉及上千個(gè)機(jī)翼外形和結(jié)構(gòu)尺寸變量。由于數(shù)值模擬非常耗時(shí),使得“維數(shù)災(zāi)難”一直以來(lái)都是困擾機(jī)翼氣動(dòng)/結(jié)構(gòu)耦合優(yōu)化的一個(gè)關(guān)鍵問(wèn)題。氣動(dòng)/結(jié)構(gòu)耦合Adjoint敏度分析方法[12]的提出大大減少了梯度求解的時(shí)間,為解決機(jī)翼氣動(dòng)/結(jié)構(gòu)優(yōu)化的“維數(shù)災(zāi)難”問(wèn)題提供了一條有效的途徑,使采用梯度優(yōu)化算法的機(jī)翼設(shè)計(jì)的變量維數(shù)達(dá)到上千量級(jí)[13-15]。但是,梯度優(yōu)化屬于局部?jī)?yōu)化,容易陷入局部最優(yōu),優(yōu)化結(jié)果強(qiáng)烈依賴(lài)于初始方案。此外,高效全局優(yōu)化算法得到高度重視和快速發(fā)展。隨著人們對(duì)全局優(yōu)化、離散/連續(xù)變量混合優(yōu)化(如復(fù)合材料機(jī)翼結(jié)構(gòu)設(shè)計(jì))需求的不斷增加[16],啟發(fā)式算法開(kāi)始用于機(jī)翼氣動(dòng)/結(jié)構(gòu)耦合優(yōu)化[17]。然而,該方法計(jì)算量巨大,因此,既能實(shí)現(xiàn)全局優(yōu)化又大大降低計(jì)算量的基于代理模型的優(yōu)化算法(簡(jiǎn)稱(chēng)代理優(yōu)化)[18-20]應(yīng)運(yùn)而生,并在機(jī)翼氣動(dòng)/結(jié)構(gòu)耦合優(yōu)化中應(yīng)用日益廣泛[8-9,21-24],已能實(shí)現(xiàn)上百個(gè)變量的機(jī)翼優(yōu)化設(shè)計(jì)。

    雖然機(jī)翼氣動(dòng)/結(jié)構(gòu)耦合優(yōu)化技術(shù)的研究已取得巨大突破,但要全面走向工程應(yīng)用,促進(jìn)我國(guó)飛機(jī)設(shè)計(jì)技術(shù)的進(jìn)步,僅僅依靠通用的優(yōu)化軟件/程序是不夠的,還需要發(fā)展專(zhuān)門(mén)的機(jī)翼優(yōu)化設(shè)計(jì)工具。由此帶來(lái)的機(jī)翼優(yōu)化參數(shù)化、數(shù)據(jù)傳遞等應(yīng)用技術(shù)的突破,還將大大降低其他布局機(jī)翼氣動(dòng)/結(jié)構(gòu)優(yōu)化技術(shù)的研究成本。本文基于變可信度氣動(dòng)數(shù)值模擬和結(jié)構(gòu)有限元數(shù)值模擬方法,考慮靜氣動(dòng)彈性效應(yīng),結(jié)合SurroOpt代理優(yōu)化軟件,旨在研發(fā)一個(gè)可實(shí)現(xiàn)跨聲速運(yùn)輸機(jī)機(jī)翼氣動(dòng)/結(jié)構(gòu)全局優(yōu)化設(shè)計(jì)的平臺(tái),命名為AeroStruct。

    1機(jī)翼氣動(dòng)/結(jié)構(gòu)數(shù)值模擬方法

    1.1氣動(dòng)數(shù)值模擬

    氣動(dòng)數(shù)值模擬可通過(guò)求解考慮邊界層黏性修正的全速勢(shì)方程、Euler方程或N-S方程實(shí)現(xiàn),分別依托于兩款自主發(fā)展的程序來(lái)完成,即WBCAero和PMNS3D。

    WBCAero是一款求解全速勢(shì)方程的程序,同時(shí)考慮了機(jī)身及機(jī)翼表面的黏性邊界層影響。它集成了自動(dòng)網(wǎng)格劃分功能,可以快速且較準(zhǔn)確地預(yù)測(cè)機(jī)翼的升阻力系數(shù)?;谟叙?無(wú)黏迭代的方法可以比較準(zhǔn)確地模擬接近分離的邊界層流動(dòng),從而有效預(yù)測(cè)機(jī)翼的升阻特性。實(shí)踐證明,對(duì)巡航狀態(tài)的運(yùn)輸機(jī)機(jī)翼氣動(dòng)分析,該方法的計(jì)算結(jié)果與NS方程計(jì)算結(jié)果比較接近,而計(jì)算效率遠(yuǎn)高于N-S方程,是一種適用于運(yùn)輸機(jī)機(jī)翼氣動(dòng)/結(jié)構(gòu)耦合優(yōu)化設(shè)計(jì)的快速方法。PMNS3D程序[25]可求解任意三維定常、非定常、可壓縮Euler方程及N-S方程。它包含豐富的湍流模型,如零方程的BL模型、一方程的SA模型及改進(jìn)型、兩方程的k -ω模型等。與全速勢(shì)方程相比,PMNS3D提供了更加精確的氣動(dòng)數(shù)值模擬方法,但由于網(wǎng)格量大且求解過(guò)程復(fù)雜,因此計(jì)算效率較低,適用于優(yōu)化過(guò)程中進(jìn)行少量高可信度氣動(dòng)/結(jié)構(gòu)耦合數(shù)值模擬。不同可信度的氣動(dòng)數(shù)值模擬程序可滿足快速設(shè)計(jì)、高可信度設(shè)計(jì)及變可信度設(shè)計(jì)的不同需求。

    1.2結(jié)構(gòu)數(shù)值模擬

    采用有限元方法,基于小變形假設(shè)對(duì)機(jī)翼進(jìn)行線性靜力分析,獲得機(jī)翼在氣動(dòng)載荷作用下的位移。為在工程設(shè)計(jì)中滿足與不同結(jié)構(gòu)求解模型搭建接口的需要,本文氣動(dòng)/結(jié)構(gòu)耦合數(shù)值模擬系統(tǒng)的結(jié)構(gòu)求解器可選用商業(yè)軟件Nastran或ANSYS。

    ANSYS求解器可通過(guò)其獨(dú)有的APDL語(yǔ)言方便、快捷地實(shí)現(xiàn)有限元分析的整個(gè)流程,包括幾何模型的建立、網(wǎng)格劃分、邊界條件的施加、求解及結(jié)果后處理等。封裝好的APDL命令流對(duì)各類(lèi)型的機(jī)翼均具有較好的適用性,能夠根據(jù)優(yōu)化設(shè)計(jì)需求參數(shù)化地改變機(jī)翼部件厚度、翼型、平面形狀、結(jié)構(gòu)布局等。本文在進(jìn)行優(yōu)化設(shè)計(jì)時(shí),多采用ANSYS作為結(jié)構(gòu)求解器。Nastran是一種在工程領(lǐng)域廣泛采用的功能強(qiáng)大的結(jié)構(gòu)分析軟件。為滿足不同對(duì)象的分析需求,本文也將其集成到機(jī)翼氣動(dòng)/結(jié)構(gòu)耦合優(yōu)化設(shè)計(jì)平臺(tái)中。對(duì)于需要精細(xì)化設(shè)計(jì)的對(duì)象而言,應(yīng)手動(dòng)建立精細(xì)的有限元模型,然后采用Nastran進(jìn)行結(jié)構(gòu)有限元分析。

    本文建立結(jié)構(gòu)有限元模型所用的幾何模型由機(jī)翼氣動(dòng)/結(jié)構(gòu)一體化參數(shù)化程序生成。幾何模型的形式分三種:關(guān)鍵點(diǎn)模型、線框模型和完整的面模型。針對(duì)不同類(lèi)型的結(jié)構(gòu)模型,本文發(fā)展的程序均提供了相應(yīng)接口,可以高效地實(shí)現(xiàn)機(jī)翼結(jié)構(gòu)靜力分析,獲得機(jī)翼在氣動(dòng)載荷作用下的節(jié)點(diǎn)位移和節(jié)點(diǎn)應(yīng)力等,使程序具有更廣泛的適用性。

    1.3氣動(dòng)/結(jié)構(gòu)數(shù)據(jù)傳遞

    氣動(dòng)/結(jié)構(gòu)耦合分析基于“弱耦合方法”進(jìn)行?!叭躐詈戏椒ā钡膬?yōu)點(diǎn)在于靜氣彈分析過(guò)程中可有效利用現(xiàn)有軟件/程序,進(jìn)行相對(duì)獨(dú)立的學(xué)科分析。

    由于氣動(dòng)、結(jié)構(gòu)學(xué)科所采用的分析方法和計(jì)算網(wǎng)格不同,因此需要通過(guò)某種數(shù)據(jù)近似/插值方式,進(jìn)行學(xué)科間的數(shù)據(jù)(氣動(dòng)分布力/結(jié)構(gòu)節(jié)點(diǎn)變形)轉(zhuǎn)換和傳遞。本文可通過(guò)兩種方式來(lái)實(shí)現(xiàn)氣動(dòng)與結(jié)構(gòu)學(xué)科間的數(shù)據(jù)傳遞:移置力矩的就近加載法(nnFM)、徑向基函數(shù)插值方法(RBFs)。

    移置力矩的就近加載法是基于力系數(shù)守恒原理,將計(jì)算流體力學(xué)(CFD)數(shù)值模擬所獲得的氣動(dòng)力傳遞至計(jì)算結(jié)構(gòu)力學(xué)(CSD)網(wǎng)格節(jié)點(diǎn)。由于通常CFD計(jì)算網(wǎng)格相較于結(jié)構(gòu)有限元分析網(wǎng)格更為密集,因此以結(jié)構(gòu)有限元網(wǎng)格節(jié)點(diǎn)為基礎(chǔ),在所有CFD網(wǎng)格點(diǎn)中搜索與之最近的氣動(dòng)網(wǎng)格節(jié)點(diǎn),將該節(jié)點(diǎn)三個(gè)方向的氣動(dòng)力移置到結(jié)構(gòu)有限元模型網(wǎng)格節(jié)點(diǎn),并施加相應(yīng)的移置力矩(見(jiàn)圖1),從而在力系數(shù)守恒的前提下實(shí)現(xiàn)氣動(dòng)力的傳遞。該方法簡(jiǎn)單高效,可快速獲得施加在各結(jié)構(gòu)節(jié)點(diǎn)的氣動(dòng)力信息,但存在一定程度的誤差。

    RBF方法[26]是一種插值類(lèi)方法,廣泛應(yīng)用于數(shù)據(jù)傳遞。該方法首先通過(guò)氣動(dòng)節(jié)點(diǎn)信息和結(jié)構(gòu)節(jié)點(diǎn)信息建立學(xué)科間的插值模型,再根據(jù)初始域內(nèi)需要傳遞的信息(如氣動(dòng)力、節(jié)點(diǎn)位移等)求得目標(biāo)域內(nèi)對(duì)應(yīng)的響應(yīng)值。該方法誤差小,可實(shí)現(xiàn)高精度數(shù)據(jù)傳遞。RBF插值方法基于能量守恒原理,認(rèn)為在氣動(dòng)與結(jié)構(gòu)的耦合界面上流體載荷(即氣動(dòng)力)與固體力在界面位移上所做的虛功相等,即

    1.4氣動(dòng)/結(jié)構(gòu)耦合數(shù)值模擬算例驗(yàn)證

    本節(jié)基于標(biāo)準(zhǔn)CRM和HIRENASD機(jī)翼靜氣彈研究模型,對(duì)所發(fā)展的氣動(dòng)/結(jié)構(gòu)耦合數(shù)值模擬程序進(jìn)行了驗(yàn)證。

    CRM模型[27]由阻力預(yù)測(cè)研討會(huì)(DPW)提出,其布局形式由美國(guó)國(guó)家航空航天局(NASA)的亞聲速固定翼空氣動(dòng)力技術(shù)研究室和DPW共同研發(fā),詳細(xì)的氣動(dòng)設(shè)計(jì)則由波音公司完成。它最初作為標(biāo)準(zhǔn)氣動(dòng)研究模型被提出,后逐漸被學(xué)者發(fā)展成為靜氣彈研究模型。此構(gòu)型包含機(jī)身、機(jī)翼、尾翼和發(fā)動(dòng)機(jī)短艙,本文所用的驗(yàn)證模型為其不含配平的風(fēng)洞試驗(yàn)?zāi)P停瑑H包含機(jī)身和機(jī)翼。CRM的風(fēng)洞試驗(yàn)結(jié)果來(lái)源于歐洲跨聲速風(fēng)洞(ETW),結(jié)構(gòu)有限元模型來(lái)源于DPW。CRM模型的CFD計(jì)算網(wǎng)格和結(jié)構(gòu)有限元網(wǎng)格如圖2、圖3所示。

    由于CRM風(fēng)洞試驗(yàn)?zāi)P蛣偠却?,靜氣彈變形量極小,通過(guò)本文的氣動(dòng)/結(jié)構(gòu)耦合數(shù)值模擬計(jì)算,在迭代三次后分析即達(dá)到收斂條件。靜氣彈收斂機(jī)翼壓力分布計(jì)算與風(fēng)洞試驗(yàn)對(duì)比如圖4所示。收斂后翼尖前緣的變形量為14.9009mm,與試驗(yàn)值(15.45545mm)之間的相對(duì)誤差為3.72159%。表明本文對(duì)CRM模型的靜氣彈分析具有較高的精度。

    HIRENASD模型[28-29]是由NASA的氣動(dòng)彈性預(yù)測(cè)研討會(huì)(Aeroelastic Prediction Workshop)提出的供學(xué)術(shù)界研究的標(biāo)準(zhǔn)靜氣彈預(yù)測(cè)模型。HIRENASD模型的風(fēng)洞試驗(yàn)由德國(guó)亞琛工業(yè)大學(xué)在歐洲跨聲速風(fēng)洞(ETW)完成,擁有詳細(xì)的試驗(yàn)數(shù)據(jù)。該模型包含機(jī)翼和機(jī)身兩部分,機(jī)翼為三段式無(wú)扭轉(zhuǎn)無(wú)上反機(jī)翼,機(jī)身為簡(jiǎn)化的修型鼓包,初始構(gòu)型為型架外形,HIRENASD模型的CFD計(jì)算網(wǎng)格和結(jié)構(gòu)有限元網(wǎng)格如圖5、圖6所示。

    通過(guò)本文的氣動(dòng)/結(jié)構(gòu)耦合數(shù)值模擬計(jì)算,經(jīng)過(guò)4次迭代后,靜氣彈分析達(dá)到收斂,獲得了變形后的機(jī)翼外形。靜氣彈收斂機(jī)翼壓力分布計(jì)算與風(fēng)洞試驗(yàn)對(duì)比如圖7所示。位于機(jī)翼展向0.95相對(duì)位置處的前緣變形為12.8mm,與試驗(yàn)值(12.5mm)之間的相對(duì)誤差僅為2.34%。表明本文對(duì)HIRENAD模型的靜氣彈分析具有較高的精度。

    基于以上兩種標(biāo)模機(jī)翼的靜氣彈數(shù)值模擬驗(yàn)證表明,本文所發(fā)展的機(jī)翼氣動(dòng)/結(jié)構(gòu)耦合數(shù)值模擬程序可行且精度較高。關(guān)于數(shù)值模擬驗(yàn)證的詳細(xì)論述見(jiàn)參考文獻(xiàn)[30]。

    2 AeroStruct:基于代理模型的跨聲速機(jī)翼氣動(dòng)/結(jié)構(gòu)耦合優(yōu)化平臺(tái)

    2.1 SurroOpt代理優(yōu)化工具箱

    SurroOpt代理優(yōu)化軟件是西北工業(yè)大學(xué)韓忠華教授[18, 31-32]團(tuán)隊(duì)開(kāi)發(fā)的基于代理模型的通用優(yōu)化工具,它可以求解任意單目標(biāo)、加權(quán)多目標(biāo)、Pareto多目標(biāo)的無(wú)約束、多約束優(yōu)化問(wèn)題。該工具箱包含豐富的試驗(yàn)設(shè)計(jì)方法,如拉丁超立方(LHS)、均勻設(shè)計(jì)(UD)、蒙特卡羅抽樣(MC);包含豐富的代理模型,如多項(xiàng)式二次響應(yīng)面(PRSM)、Kriging模型、梯度增強(qiáng)Kriging模型(GEK)、分層Kriging模型(HK)、徑向基函數(shù)(RBF)等;包含多種加點(diǎn)準(zhǔn)則,如最小化代理模型預(yù)測(cè)準(zhǔn)則(MSP)、改善期望準(zhǔn)則(EI)、改善概率準(zhǔn)則(PI)、均方差準(zhǔn)則(MSE)和置信下界準(zhǔn)則(LCB)等;包含多種子優(yōu)化方法,如Hooke-Jeeves模式搜索、擬牛頓梯度優(yōu)化、序列二次規(guī)劃(SQP)、遺傳算法等。SurroOpt經(jīng)大量數(shù)值算例驗(yàn)證,已成功應(yīng)用于多個(gè)實(shí)際工程優(yōu)化設(shè)計(jì)問(wèn)題中。

    2.2機(jī)翼模型參數(shù)化方法

    在機(jī)翼氣動(dòng)/結(jié)構(gòu)耦合優(yōu)化設(shè)計(jì)過(guò)程中,設(shè)計(jì)對(duì)象是機(jī)翼的氣動(dòng)外形、結(jié)構(gòu)布局及厚度參數(shù),因此需要將機(jī)翼的氣動(dòng)外形和結(jié)構(gòu)布局通過(guò)參數(shù)表示。本文采用了一體化的氣動(dòng)/結(jié)構(gòu)參數(shù)化,一次計(jì)算同時(shí)獲得參數(shù)化的機(jī)翼氣動(dòng)外形和結(jié)構(gòu)布局。AeroStruct平臺(tái)參數(shù)化流程如圖8所示:分為氣動(dòng)外形參數(shù)化和結(jié)構(gòu)參數(shù)化兩大模塊,基于同一輸入文件,各自提取所需信息進(jìn)行計(jì)算,生成相應(yīng)的外形數(shù)據(jù)文件和結(jié)構(gòu)線框模型。詳細(xì)過(guò)程如下。

    (1)氣動(dòng)外形參數(shù)化采用型函數(shù)/類(lèi)函數(shù)變換方法(CST)或PARSEC方法實(shí)現(xiàn)。在對(duì)三維機(jī)翼進(jìn)行參數(shù)化的過(guò)程中,首先通過(guò)二維參數(shù)化方法將機(jī)翼展向各剖面參數(shù)化,再根據(jù)機(jī)翼平面外形參數(shù)(展弦比、后掠角等)參數(shù)化整個(gè)機(jī)翼,最終獲得由數(shù)據(jù)集描述的機(jī)翼氣動(dòng)外形。

    (2)結(jié)構(gòu)布局參數(shù)化可通過(guò)ANSYS的APDL語(yǔ)言編程,讀入外形數(shù)據(jù)點(diǎn)及結(jié)構(gòu)尺寸信息直接建模。也可以基于CATIA的VBA二次開(kāi)發(fā)實(shí)現(xiàn)精細(xì)化建模:基于參數(shù)化機(jī)翼氣動(dòng)外形,根據(jù)梁、肋等結(jié)構(gòu)部件的位置、角度等信息生成結(jié)構(gòu)部件的外輪廓;為了保證數(shù)據(jù)傳輸?shù)木?,參?shù)化后獲得的是幾何文件iges格式的機(jī)翼結(jié)構(gòu)線框模型。

    2.3基于代理模型的跨聲速機(jī)翼氣動(dòng)/結(jié)構(gòu)優(yōu)化設(shè)計(jì)平臺(tái)——AeroStruct

    基于SurroOpt代理優(yōu)化軟件及所開(kāi)發(fā)的機(jī)翼氣動(dòng)/結(jié)構(gòu)耦合數(shù)值模擬系統(tǒng),搭建機(jī)翼氣動(dòng)/結(jié)構(gòu)耦合全局優(yōu)化設(shè)計(jì)平臺(tái)。如圖9所示,右側(cè)為機(jī)翼氣動(dòng)/結(jié)構(gòu)耦合數(shù)值模擬的詳細(xì)流程,這里不再贅述;左側(cè)為完整的代理優(yōu)化流程,分為以下5步:(1)建立優(yōu)化問(wèn)題的數(shù)學(xué)模型,確定優(yōu)化目標(biāo)和約束,確定設(shè)計(jì)變量及其上下界;(2)通過(guò)試驗(yàn)設(shè)計(jì)(DOE)方法選取適當(dāng)數(shù)量的初始樣本點(diǎn),通過(guò)氣動(dòng)/結(jié)構(gòu)耦合數(shù)值模擬程序計(jì)算初始樣本點(diǎn)對(duì)應(yīng)的響應(yīng)值,響應(yīng)值包括機(jī)翼升力、阻力、力矩及結(jié)構(gòu)質(zhì)量、單元應(yīng)力、最大位移等;(3)根據(jù)初始樣本點(diǎn)建立目標(biāo)函數(shù)和約束函數(shù)的初始代理模型;(4)在當(dāng)前代理模型上進(jìn)行子優(yōu)化,通過(guò)加點(diǎn)準(zhǔn)則指導(dǎo)加入新的樣本點(diǎn),并由數(shù)值模擬程序計(jì)算響應(yīng)值,加入到樣本點(diǎn)數(shù)據(jù)庫(kù)中,用于更新代理模型;(5)在未達(dá)到優(yōu)化收斂條件時(shí)重復(fù)以上過(guò)程,達(dá)到最大樣本點(diǎn)個(gè)數(shù)或代理模型預(yù)測(cè)誤差條件后終止迭代,優(yōu)化結(jié)束。

    3 uCRM-9標(biāo)模機(jī)翼氣動(dòng)/結(jié)構(gòu)耦合優(yōu)化設(shè)計(jì)

    3.1優(yōu)化問(wèn)題描述

    uCRM-9(undeflected Common Research Model-9)模型是Brooks等[33]在CRM機(jī)翼基礎(chǔ)上設(shè)計(jì)的未變形機(jī)翼模型,即機(jī)翼型架構(gòu)型。uCRM-9機(jī)翼扭轉(zhuǎn)和上反特征復(fù)雜,與波音777客機(jī)的機(jī)翼相當(dāng),適用于機(jī)翼氣動(dòng)/結(jié)構(gòu)耦合數(shù)值模擬研究及優(yōu)化設(shè)計(jì)研究。uCRM-9模型包含機(jī)身、機(jī)翼和平尾,其構(gòu)型如圖10所示。

    3.2優(yōu)化設(shè)計(jì)結(jié)果及分析

    優(yōu)化采用拉丁超立方采樣獲得設(shè)計(jì)空間內(nèi)的152個(gè)初始樣本點(diǎn),建立Kriging代理模型,利用EI+PI+LCB+MSP加點(diǎn)準(zhǔn)則在每輪優(yōu)化迭代中并行加入7個(gè)樣本點(diǎn)(其中EI兩個(gè)點(diǎn)、PI兩個(gè)點(diǎn)、LCB兩個(gè)點(diǎn)、MSP一個(gè)點(diǎn)),優(yōu)化迭代終止條件為優(yōu)化加點(diǎn)達(dá)到最大樣本點(diǎn)個(gè)數(shù)1000。

    uCRM-9機(jī)翼氣動(dòng)/結(jié)構(gòu)耦合優(yōu)化迭代收斂歷程如圖11所示。機(jī)翼結(jié)構(gòu)質(zhì)量在第一輪尋優(yōu)后迅速下降至最優(yōu)值附近,而后緩慢降低,直至達(dá)到最大迭代步數(shù),代理模型在最優(yōu)點(diǎn)處預(yù)測(cè)誤差為0.00001063%。初始機(jī)翼結(jié)構(gòu)質(zhì)量10944kg,優(yōu)化后為9501kg,減重約13%,優(yōu)化效果顯著。

    圖12為優(yōu)化后機(jī)翼上下翼面蒙皮厚度的展向分布及與優(yōu)化前厚度分布對(duì)比,其中縱坐標(biāo)為厚度,橫坐標(biāo)代表由翼根至翼尖劃分的38個(gè)設(shè)計(jì)區(qū)域。優(yōu)化設(shè)計(jì)前上下蒙皮厚度由翼根至翼尖均勻減小,優(yōu)化后翼根到翼尖上下蒙皮厚度均先緩慢增加再減小??梢园l(fā)現(xiàn):通過(guò)優(yōu)化設(shè)計(jì)減掉的機(jī)翼結(jié)構(gòu)質(zhì)量大部分由翼根位置貢獻(xiàn),因此翼根位置蒙皮厚度遠(yuǎn)低于初始設(shè)計(jì)厚度。這是因?yàn)椋弘m然靠近翼根處機(jī)翼所受彎矩最大,但由于翼根剖面翼型相對(duì)厚度大,翼盒的截面積大,因而應(yīng)力水平較低;而靠近內(nèi)外翼轉(zhuǎn)折(Kink)處的翼盒截面積相對(duì)翼根小,應(yīng)力水平較高,故蒙皮最厚。在大部分設(shè)計(jì)區(qū)域內(nèi),機(jī)翼下蒙皮厚度整體小于上蒙皮厚度,是由于上蒙皮主要承受壓應(yīng)力和切應(yīng)力,更容易失穩(wěn),符合真實(shí)機(jī)翼特征。

    圖13為uCRM-9機(jī)翼最優(yōu)方案的節(jié)點(diǎn)應(yīng)力云圖(左)及壓力云圖(右),圖14為機(jī)翼優(yōu)化前后節(jié)點(diǎn)應(yīng)力云圖,圖15為優(yōu)化前后y方向節(jié)點(diǎn)位移云圖。通過(guò)以上計(jì)算結(jié)果可以看出:對(duì)于基準(zhǔn)機(jī)翼,由于翼根到翼梢蒙皮厚度采用簡(jiǎn)單的線性遞減方案,應(yīng)力集中在機(jī)翼Kink偏外翼段的位置;且整體剛度過(guò)大,導(dǎo)致翼尖撓度過(guò)小。而對(duì)于優(yōu)化后的機(jī)翼,上下蒙皮厚度分布更加合理(由翼根到翼尖呈非線性分布),內(nèi)翼段蒙皮厚度較基準(zhǔn)方案減小,主承力區(qū)蒙皮厚度有所增加。因此優(yōu)化后機(jī)翼應(yīng)力分布更加均勻,最大應(yīng)力位于靠近機(jī)翼機(jī)身連接處后緣位置。由于真實(shí)機(jī)翼翼根位置通常布置加強(qiáng)肋,因此最大應(yīng)力位于此處是可以接受的。優(yōu)化后的機(jī)翼在滿足翼尖最大位移約束和升力系數(shù)約束的前提下,獲得了更合理的蒙皮厚度分布,減輕了結(jié)構(gòu)質(zhì)量。

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

    本文基于變可信度機(jī)翼氣動(dòng)/結(jié)構(gòu)耦合數(shù)值模擬程序,發(fā)展了基于代理模型的機(jī)翼氣動(dòng)/結(jié)構(gòu)高效全局耦合優(yōu)化設(shè)計(jì)平臺(tái)AeroStruct。對(duì)uCRM-9機(jī)翼包含76個(gè)設(shè)計(jì)變量的氣動(dòng)/結(jié)構(gòu)進(jìn)行優(yōu)化設(shè)計(jì),最優(yōu)方案相比于基準(zhǔn)機(jī)翼結(jié)構(gòu)質(zhì)量減小約13%,優(yōu)化后的機(jī)翼上下翼面蒙皮厚度分布合理,符合真實(shí)機(jī)翼特征,驗(yàn)證了所發(fā)展的優(yōu)化設(shè)計(jì)平臺(tái)的正確性和有效性。本文研發(fā)的機(jī)翼氣動(dòng)/結(jié)構(gòu)優(yōu)化設(shè)計(jì)系統(tǒng),為未來(lái)飛行器機(jī)翼優(yōu)化設(shè)計(jì)和新型機(jī)翼設(shè)計(jì)技術(shù)探索提供了重要手段。

    參考文獻(xiàn)

    [1]Haftka R T. Optimization of flexible wing structures subject to strength and induced drag constraints[J]. AIAA Journal,1977,14(8):1106-1977.

    [2]Martins J R R A,Alonso J J,Reuther J J. High-fidelity aerostructural design optimization of a supersonic business jet[J]. Journal ofAircraft,2004,41(3):523-530.

    [3]Kenway G K W,Martins J R R A. Multi-point high-fidelity aerostructural optimization of a transport aircraft configuration[J]. Journal ofAircraft,2014,51(1):144-160.

    [4]Ronzheimer A,Natterer F J,Brezillon J. Aircraft wing optimization using high fidelity closely coupled CFD and CSM methods[R].AIAA2010-9078,2010.

    [5]Abu-Zurayk M,Brezillon J. Development and evaluation of adjoint approach for aeroelastic wing optimization[J]. New Results in Numerical and Experimental Fluid Mechanics,2013,121:59-66.

    [6]Ghazlane I,Carrier G,Dumont A,et al. Aerostructural adjoint method for flexible wing optimization[R]. AIAA 2012-1924,2012.

    [7]Chittick I R,Martins J R R A. An asymmetric suboptimization approach to aerostructural optimization[J]. Optimization and Engineering,2009,10(1):133-152.

    [8]Zhang K S,Han Z H,Gao Z J,et al. Constraint aggregation for large number of constraints in wing surrogate-based optimization[J]. Structural and Multidisciplinary Optimization,2019,59(2):421-438.

    [9]左英桃,王曉鵬,陳云,等.一種高效的CFD/CSD耦合飛行器多學(xué)科優(yōu)化設(shè)計(jì)方法[J].航空動(dòng)力學(xué)報(bào), 2014, 29(12): 2898-2904. Zuo Yingtao, Wang Xiaopeng, Chen Yun, et al. An efficient method for multidisciplinary design optimization of aircraft based on CFD/CSD coupling[J]. Journal of Aerospace Power, 2014, 29(12): 2898-2904. (in Chinese)

    [10]胡婕,王如華,王穩(wěn)江,等.客機(jī)機(jī)翼氣動(dòng)/結(jié)構(gòu)多學(xué)科優(yōu)化方法[J].南京航空航天大學(xué)學(xué)報(bào), 2012, 44(4): 458-463. Hu Jie, Wang Ruhua, Wang Wenjiang, et al. Multidisciplinary optimizationoftransportwingaerodynamic/structural integrateddesign[J].JournalofNanjingUniversityof Aeronautics &Astronautics, 2012, 44(4): 458-463. (in Chinese)

    [11]Jonsson E,Kenway G K W,Martins J R RA,et al. Developmentofflutterconstraintsforhigh-fidelity aerostructural optimization[R]AIAA2017-4455,2017.

    [12]Martins J R R A,Alonso J J,Reuther J J. A coupled-adjoint sensitivity analysis method for high-fidelity aero-structural design[J]. Optimization and Engineering,2005,6(1):33-62.

    [13]Brooks T R,Kenway G K W,Martins J R R A. Undeflected common research model(uCRM):an aerostructural model for the study of high aspect ratio transport aircraft wings[R]. AIAA 2017-4456,2017.

    [14]黃江濤,劉剛,高正紅,等.飛行器多學(xué)科耦合伴隨體系的現(xiàn)狀與發(fā)展趨勢(shì)綜述[J].航空學(xué)報(bào), 2020, 41(5): 6-29. Huang Jiangtao, Liu Gang, Gao Zhenghong, et al. Current situation and development trend of multidisciplinary coupled adjoint system for aircraft[J]. Acta Aeronautica et Astronautica Sinica, 2020, 41(5): 6-29. (in Chinese)

    [15]雷銳午,白俊強(qiáng),許丹陽(yáng),等.基于耦合伴隨方法的串/并行氣動(dòng)結(jié)構(gòu)優(yōu)化設(shè)計(jì)對(duì)比[J].航空動(dòng)力學(xué)報(bào), 2019, 34(5): 1036-1049. Lei Ruiwu, Bai Junqiang, Xu Danyang, et al. Speciality assessmentofsequentialandconcurrentaerostructural optimization based on coupled adjoint technique[J]. Journal of Aerospace Power, 2019, 34(5): 1036-1049. (in Chinese)

    [16]趙占文,蘇雁飛,崔勇江.變剛度復(fù)合材料層壓板的軸壓穩(wěn)定性優(yōu)化研究[J].航空科學(xué)技術(shù), 2021, 32(8): 57-61. ZhaoZhanwen,SuYanfei,CuiYongjiang.Buckling optimization study on variable stiffness composite laminate undercompressionload[J].AeronauticalScience& Technology, 2021, 32(8): 57-61. (in Chinese)

    [17]宋倩,萬(wàn)志強(qiáng).飛翼式客機(jī)機(jī)翼氣動(dòng)/結(jié)構(gòu)綜合優(yōu)化方法研究[J].民用飛機(jī)設(shè)計(jì)與研究, 2018(4): 6-14. Song Qian, Wan Zhiqiang. Aerodynamic /structural integrated design for the wing of flying wing airliner[J]. Civil Aircraft Design and Research, 2018 (4): 6-14. (in Chinese)

    [18]韓忠華. Kriging模型及代理優(yōu)化算法研究新進(jìn)展[J].航空學(xué)報(bào), 2016, 37(11): 3197-3225. Han Zhonghua. Kriging surrogate model and its application to design optimization: a review of recent progress[J]. Acta Aeronautica et Astronautica Sinica, 2016, 37(11): 3197-3225.(in Chinese)

    [19]Yondo R,Andres E,Valero E. A review on design of experiments and surrogate models in aircraft real-time and many-query aerodynamic analyses[J]. Progress in Aerospace Sciences,2018,96:23-61.

    [20]張陽(yáng),韓忠華,柳斐,等.高超聲速飛行器寬速域翼型高效多目標(biāo)優(yōu)化設(shè)計(jì)方法研究[J].航空科學(xué)技術(shù), 2020, 31(11): 14-24. Zhang Yang, Han Zhonghua, Liu Fei, et al. Efficient multiobjective shape optimization method of hypersonic wide-machnumber-range airfoil[J]. Aeronautical Science & Technology, 2020, 31(11): 14-24.(in Chinese)

    [21]張科施,韓忠華,李為吉,等.一種考慮氣動(dòng)彈性的運(yùn)輸機(jī)機(jī)翼多學(xué)科優(yōu)化方法[J].空氣動(dòng)力學(xué)報(bào), 2008, 26(1): 1-7. Zhang Keshi, Han Zhonghua, Li Weiji, et al. A method of coupled aerodynamic/structural integration optimization for transport-wing design[J]. Acta Aerodynamica Sinica, 2008, 26(1): 1-7. (in Chinese)

    [22]Klimmek T,Schulze M,Abu-Zurayk M,et al. CPACS-MONA:an independent and in high fidelity based MDO tasks integrated process for the structural and aeroelastic design for aircraft configurations[C]//International Forum on Aeroelasticity and Structural Dynamics,IFASD 2019,2019.

    [23]Goertz S,Ilic C,Jepson J,et al. Multi-level MDO of a longrange transport aircraft using a distributed analysis framework[R].AIAA2017-4326,2017.

    [24]Wunderlich T,Daehne S,Heinrich L,et al. Multidisciplinary optimization of an NLF forward swept wing in combination with aeroelastic tailoring using CFRP[J]. CEAS Aeronautical Journal,2017,8:673-690.

    [25]Han Z H,He F,Song W P,et al. A Preconditioned multigrid methodforefficientsimulationofthree-dimensional compressible and incompressible flows[J]. Chinese Journal of Aeronautics,2007,20(4):289-296.

    [26]Armin B,Holger W. Multivariate interpolation for fluidstructure-interaction problems using radial basis functions[J]. Aerospace Science and Technology,2001,5(2):125-134.

    [27]Vassberg J,Dehaan M A,Rivers S M,et al. Development of a common research model for applied CFD validation studies[R].AIAA2008-6919,2008.

    [28]Reimer L,Braun C,Chen B H,et al. Computational aeroelastic design and analysis of the HIRENASD wind tunnel model and tests[C]// International Forum on Aeroelasticity and Structural Dynamics,IF-077,Stockholm,2007.

    [29]Wieseman C,Chwalowski P,Heeg J,et al. Structural dynamicsmodelingofHIRENASDinsupportofthe aeroelastic prediction workshop[C]// 54th AIAA/ASME/ASCE/ AHS /ASC Structures,Structure Dynamics,and Materials Conference,Boston,Massachusetts,2013.

    [30]Ling S B,Zhang K S,Han Z H,et al. Evaluation and validation of aeroelastic analysis program for the“AeroStruct”wing optimization platform[C]//32nd Congress of the International Council of the Aeronautical Sciences,Shanghai,2021.

    [31]Han Z H. SurroOpt:a generic surrogate -based optimization code for aerodynamic and multidisciplinary design[C]//30th Congress of the International Council of the Aeronautical Sciences,2016.

    [32]韓忠華,許晨舟,喬建領(lǐng),等.基于代理模型的高效全局氣動(dòng)優(yōu)化設(shè)計(jì)方法研究進(jìn)展[J].航空學(xué)報(bào), 2020, 41(5): 30-70. Han Zhonghua, Xu Chenzhou, Qiao Jianling, et al. Recent progress of efficient global aerodynamic shape optimization usingsurrogate-basedapproach[J].ActaAeronauticaet Astronautica Sinica, 2020, 41(5): 30-70. (in Chinese)

    [33]Brooks T R,Kenway G K W,Martins J R R A. Benchmark aerostructural models for the study of transonic aircraft wings[J].AIAAJournal,2018,56(7):2840-2855.

    [34]Kim Y,Lee D H,Kim Y,et al. Multidisciplinary design optimization of supersonic fighter wing using response surface methodology[R].AIAA2002-5408,2002.

    [35]Wakayama S,Kroo I. The challenge and promise of blendwing-body optimization[R].AIAA98-4736,1998.

    Development and Application of AeroStruct,an Aerodynamic/Structural Optimization Platform for Transonic Transport Aircraft Wings

    Zhang Keshi1,Ling Shengbo1,Han Zhonghua1,2

    1. Northwestern Polytechnical University,Xi’an 710072,China

    2. National Key laboratory of Science and Technology on Aerodynamic Design and Research,Xi’an 710072,China

    Abstract: As their aspect-ratios are high, the transonic aircraft wings exhibit geometrically nonlinear behavior so that it cannot be neglected. The advances in Multidisciplinary Design Optimization (MDO) have enabled optimizations of aircraft wings using high-fidelity simulations of coupled aerodynamic and structural behavior. This work aims to establish an aero/structural design optimization platform for the transonic wings, which is meaningful for promoting aircraft design technology. In the "AeroStruct" platform, the codes of Navier-Stokes (N-S), Euler, full-potential equations are available for the aerodynamic simulations, NASTRAN and ANSYS platforms are possible for the structural finite element analysis, and the Radius Basis Function (RBF) is realized for interpolation of fluid/structural interaction. The in-house code "SurroOpt", as a tool of surrogate-based optimizations, is used to perform wing aero/ structural global optimization. Our platform is preliminarily validated by an uCRM-9 wing optimization problem with 76 design variables, in which the wing structural weight is reduced by 13% with the limitation of some aero/structural constraints.

    Key Words:wing design;multidisciplinary design optimization;CFD/CSD coupling numerical simulation; aerodynamic design; structural design

    猜你喜歡
    氣動(dòng)結(jié)構(gòu)設(shè)計(jì)
    高層建筑連體結(jié)構(gòu)設(shè)計(jì)
    中寰氣動(dòng)執(zhí)行機(jī)構(gòu)
    中寰氣動(dòng)執(zhí)行機(jī)構(gòu)
    超限高層建筑結(jié)構(gòu)設(shè)計(jì)與優(yōu)化思考
    結(jié)構(gòu)設(shè)計(jì)優(yōu)化在房屋建筑結(jié)構(gòu)設(shè)計(jì)中的應(yīng)用
    基于NACA0030的波紋狀翼型氣動(dòng)特性探索
    基于反饋線性化的RLV氣動(dòng)控制一體化設(shè)計(jì)
    一種新型衛(wèi)生防疫車(chē)結(jié)構(gòu)設(shè)計(jì)
    BIM結(jié)構(gòu)設(shè)計(jì)應(yīng)用
    某排水虹吸井結(jié)構(gòu)設(shè)計(jì)
    亚洲久久久久久中文字幕| 特级一级黄色大片| 国产精品永久免费网站| 亚洲欧美日韩高清专用| 亚洲欧美清纯卡通| 丝袜美腿在线中文| 91在线精品国自产拍蜜月| 乱人视频在线观看| 国产精品三级大全| 国产av不卡久久| 日韩欧美三级三区| 中文精品一卡2卡3卡4更新| 国产一区二区在线av高清观看| 国产精品人妻久久久影院| 九色成人免费人妻av| 中文精品一卡2卡3卡4更新| av天堂在线播放| 天堂中文最新版在线下载 | 日韩高清综合在线| 六月丁香七月| 毛片一级片免费看久久久久| 色播亚洲综合网| 美女内射精品一级片tv| 精品久久久久久久久亚洲| 联通29元200g的流量卡| 啦啦啦韩国在线观看视频| 看片在线看免费视频| 欧美一级a爱片免费观看看| 国产成人a区在线观看| 乱码一卡2卡4卡精品| 国内久久婷婷六月综合欲色啪| 成人二区视频| 久久久色成人| 能在线免费观看的黄片| 亚洲成人久久爱视频| 亚洲欧美中文字幕日韩二区| 国产人妻一区二区三区在| 最近2019中文字幕mv第一页| 美女 人体艺术 gogo| 韩国av在线不卡| 干丝袜人妻中文字幕| 精品久久久久久成人av| 中出人妻视频一区二区| 欧美色视频一区免费| 午夜激情欧美在线| 欧美潮喷喷水| 欧美区成人在线视频| 欧美又色又爽又黄视频| a级毛片免费高清观看在线播放| 在现免费观看毛片| 久久精品国产鲁丝片午夜精品| 欧美日本视频| 日韩人妻高清精品专区| 亚洲欧洲日产国产| 一个人观看的视频www高清免费观看| 99热6这里只有精品| 美女被艹到高潮喷水动态| 欧美区成人在线视频| 五月伊人婷婷丁香| 一夜夜www| 婷婷色av中文字幕| 不卡视频在线观看欧美| 五月伊人婷婷丁香| 偷拍熟女少妇极品色| 精品一区二区三区视频在线| 午夜精品一区二区三区免费看| 亚洲av二区三区四区| 在线观看一区二区三区| 日本av手机在线免费观看| 久久久久性生活片| 尾随美女入室| av国产免费在线观看| 中文精品一卡2卡3卡4更新| 亚洲五月天丁香| 精品日产1卡2卡| 日韩成人伦理影院| 久久精品国产自在天天线| 欧美一区二区亚洲| 国产单亲对白刺激| 久久久久久久久久久免费av| 热99在线观看视频| 深爱激情五月婷婷| 麻豆精品久久久久久蜜桃| 国产成人午夜福利电影在线观看| 亚洲激情五月婷婷啪啪| 舔av片在线| 麻豆精品久久久久久蜜桃| 国产午夜精品一二区理论片| 神马国产精品三级电影在线观看| 日韩欧美三级三区| 99热6这里只有精品| 精品久久久久久久人妻蜜臀av| 日韩,欧美,国产一区二区三区 | 久久久a久久爽久久v久久| 成年版毛片免费区| 亚洲美女视频黄频| 99久国产av精品国产电影| 黄色欧美视频在线观看| 欧美一区二区精品小视频在线| 爱豆传媒免费全集在线观看| 一区二区三区四区激情视频 | 国产精品.久久久| 成人三级黄色视频| 美女xxoo啪啪120秒动态图| 51国产日韩欧美| 国产精品人妻久久久久久| 日韩一本色道免费dvd| 男女边吃奶边做爰视频| 老司机影院成人| 日韩视频在线欧美| 2022亚洲国产成人精品| 午夜亚洲福利在线播放| 久久精品国产鲁丝片午夜精品| 国产伦理片在线播放av一区 | 丰满人妻一区二区三区视频av| 亚洲欧美日韩卡通动漫| 又黄又爽又刺激的免费视频.| 成人亚洲欧美一区二区av| 精品一区二区免费观看| 一级毛片久久久久久久久女| 欧美3d第一页| 国产黄片美女视频| 哪里可以看免费的av片| 简卡轻食公司| 亚洲无线观看免费| 26uuu在线亚洲综合色| 色播亚洲综合网| 成年av动漫网址| 不卡视频在线观看欧美| 深夜精品福利| 欧美xxxx性猛交bbbb| 亚洲av成人精品一区久久| 国产一级毛片七仙女欲春2| 黑人高潮一二区| 男人的好看免费观看在线视频| 啦啦啦观看免费观看视频高清| 免费观看在线日韩| 久久精品综合一区二区三区| 桃色一区二区三区在线观看| 一个人看视频在线观看www免费| 日产精品乱码卡一卡2卡三| 亚洲无线在线观看| 亚洲不卡免费看| av天堂中文字幕网| 精品久久久久久久久av| 国产精品蜜桃在线观看 | 99久久中文字幕三级久久日本| 亚洲最大成人手机在线| av专区在线播放| 此物有八面人人有两片| av在线天堂中文字幕| 91精品国产九色| 婷婷亚洲欧美| 精品午夜福利在线看| 国产成人aa在线观看| 久久久久久久久中文| 日日啪夜夜撸| 99久久无色码亚洲精品果冻| 精品久久久久久成人av| 亚洲一区二区三区色噜噜| 长腿黑丝高跟| 美女高潮的动态| 亚洲,欧美,日韩| 97人妻精品一区二区三区麻豆| 嫩草影院入口| 美女被艹到高潮喷水动态| 亚洲精品自拍成人| 免费人成视频x8x8入口观看| 成人一区二区视频在线观看| 久久国产乱子免费精品| 国产精品久久久久久精品电影小说 | 亚洲欧美日韩东京热| 日韩人妻高清精品专区| 亚洲国产精品成人综合色| 亚洲欧美精品综合久久99| 中国美白少妇内射xxxbb| 村上凉子中文字幕在线| 六月丁香七月| 久久亚洲国产成人精品v| 国产一区二区三区av在线 | 国产亚洲精品久久久久久毛片| 国产私拍福利视频在线观看| 日本欧美国产在线视频| 欧美激情久久久久久爽电影| 国产高清视频在线观看网站| 色视频www国产| 亚洲色图av天堂| 久久久久久久久大av| 丝袜美腿在线中文| 国产精品蜜桃在线观看 | 国产精品久久久久久久电影| 爱豆传媒免费全集在线观看| 欧美丝袜亚洲另类| 免费无遮挡裸体视频| 丰满的人妻完整版| 三级毛片av免费| 伦理电影大哥的女人| 久久久久久久久大av| 在现免费观看毛片| 白带黄色成豆腐渣| 亚洲欧洲日产国产| 一级黄片播放器| 看十八女毛片水多多多| 国产av不卡久久| 青春草国产在线视频 | 久久久久久久久久黄片| 国产中年淑女户外野战色| 国内精品一区二区在线观看| 国产男人的电影天堂91| av女优亚洲男人天堂| 91狼人影院| 日韩精品有码人妻一区| 在现免费观看毛片| 波野结衣二区三区在线| 国产一区二区激情短视频| 国产精品久久久久久精品电影小说 | 日韩欧美在线乱码| 男人舔女人下体高潮全视频| 亚洲国产精品成人久久小说 | 国产精品三级大全| 成人二区视频| 精品一区二区三区人妻视频| 校园人妻丝袜中文字幕| 国产综合懂色| 自拍偷自拍亚洲精品老妇| 久久6这里有精品| 国产乱人偷精品视频| 亚洲精品亚洲一区二区| 97热精品久久久久久| 国内精品久久久久精免费| 男人舔奶头视频| 狠狠狠狠99中文字幕| 国产大屁股一区二区在线视频| 精品一区二区三区人妻视频| 欧美xxxx黑人xx丫x性爽| 99热精品在线国产| 久久精品国产自在天天线| 精品人妻熟女av久视频| 国语自产精品视频在线第100页| 精品国产三级普通话版| 国产真实乱freesex| 国产探花极品一区二区| 美女黄网站色视频| 午夜福利成人在线免费观看| 日本黄大片高清| 97在线视频观看| 精品日产1卡2卡| 亚洲一区二区三区色噜噜| 成人午夜高清在线视频| 成人欧美大片| 少妇被粗大猛烈的视频| 最近2019中文字幕mv第一页| 欧美日韩综合久久久久久| 精品久久久久久久久av| 欧美成人精品欧美一级黄| h日本视频在线播放| 久久99蜜桃精品久久| 两个人视频免费观看高清| 秋霞在线观看毛片| 久久久欧美国产精品| 日韩精品有码人妻一区| 韩国av在线不卡| 色综合色国产| 丰满乱子伦码专区| 国内久久婷婷六月综合欲色啪| 精品不卡国产一区二区三区| 欧美日韩乱码在线| 国产精品麻豆人妻色哟哟久久 | 日韩成人伦理影院| 丰满人妻一区二区三区视频av| 极品教师在线视频| 午夜福利高清视频| av卡一久久| 18+在线观看网站| 99久久精品热视频| 色视频www国产| 边亲边吃奶的免费视频| 性插视频无遮挡在线免费观看| 国产精品永久免费网站| 最好的美女福利视频网| 毛片女人毛片| 高清午夜精品一区二区三区 | 黄色配什么色好看| 噜噜噜噜噜久久久久久91| 中文字幕免费在线视频6| 精品人妻偷拍中文字幕| 最后的刺客免费高清国语| 亚洲国产日韩欧美精品在线观看| 亚洲av免费在线观看| 麻豆国产av国片精品| 欧美+日韩+精品| 亚洲国产欧美人成| 中文字幕人妻熟人妻熟丝袜美| 亚洲精品日韩av片在线观看| 麻豆成人午夜福利视频| 九九在线视频观看精品| 卡戴珊不雅视频在线播放| 国产伦精品一区二区三区四那| 久久精品国产自在天天线| 日本欧美国产在线视频| 三级经典国产精品| 超碰av人人做人人爽久久| 免费黄网站久久成人精品| 中文字幕av成人在线电影| 成人无遮挡网站| 亚洲av中文字字幕乱码综合| 亚州av有码| 中出人妻视频一区二区| 2022亚洲国产成人精品| 简卡轻食公司| 亚洲欧美精品自产自拍| 精品久久久久久久久av| 亚洲av不卡在线观看| 欧美+亚洲+日韩+国产| 18禁在线无遮挡免费观看视频| 久久精品国产亚洲av涩爱 | a级一级毛片免费在线观看| 国产高清激情床上av| 日韩亚洲欧美综合| 人妻夜夜爽99麻豆av| 日韩av在线大香蕉| 免费一级毛片在线播放高清视频| 国产伦精品一区二区三区视频9| 久久草成人影院| 三级国产精品欧美在线观看| 卡戴珊不雅视频在线播放| 久久精品国产亚洲av涩爱 | 特级一级黄色大片| 国产一区二区在线av高清观看| 国内久久婷婷六月综合欲色啪| www日本黄色视频网| 国产精品一区www在线观看| 亚洲av二区三区四区| 一区二区三区四区激情视频 | 亚洲欧美成人综合另类久久久 | 国产白丝娇喘喷水9色精品| 中文精品一卡2卡3卡4更新| 欧美高清性xxxxhd video| 波多野结衣高清无吗| 一本久久精品| 一级黄色大片毛片| 99久久无色码亚洲精品果冻| www.av在线官网国产| 欧美性猛交黑人性爽| 99久久人妻综合| 少妇被粗大猛烈的视频| 国产69精品久久久久777片| 午夜福利视频1000在线观看| 成年av动漫网址| 蜜臀久久99精品久久宅男| 亚洲最大成人手机在线| 成年版毛片免费区| 观看美女的网站| 日韩一区二区视频免费看| 亚洲av成人av| 级片在线观看| 神马国产精品三级电影在线观看| 少妇人妻精品综合一区二区 | 一边亲一边摸免费视频| 十八禁国产超污无遮挡网站| 久久亚洲精品不卡| 日本黄色片子视频| 亚洲精品久久久久久婷婷小说 | 一区二区三区免费毛片| 伊人久久精品亚洲午夜| 少妇的逼好多水| 国产综合懂色| 男女做爰动态图高潮gif福利片| 如何舔出高潮| 91午夜精品亚洲一区二区三区| 性插视频无遮挡在线免费观看| 日韩,欧美,国产一区二区三区 | a级毛片a级免费在线| 亚洲国产精品国产精品| 少妇被粗大猛烈的视频| 蜜桃亚洲精品一区二区三区| 22中文网久久字幕| 国产亚洲av嫩草精品影院| 波野结衣二区三区在线| 欧美在线一区亚洲| 国产精品人妻久久久久久| 91久久精品国产一区二区三区| 丝袜喷水一区| 九色成人免费人妻av| 国产精品麻豆人妻色哟哟久久 | 日本免费一区二区三区高清不卡| 在现免费观看毛片| 国产高清不卡午夜福利| 97热精品久久久久久| 色综合色国产| 最近中文字幕高清免费大全6| 欧美日本视频| 成人漫画全彩无遮挡| 久久久久免费精品人妻一区二区| 中国国产av一级| 国产伦精品一区二区三区四那| 精品久久久久久久久av| 欧美潮喷喷水| 哪里可以看免费的av片| 成人午夜精彩视频在线观看| 久久久精品94久久精品| 免费人成在线观看视频色| 高清日韩中文字幕在线| 高清在线视频一区二区三区 | 最好的美女福利视频网| 嘟嘟电影网在线观看| 国产极品天堂在线| 人妻少妇偷人精品九色| 日日啪夜夜撸| 国产精品一二三区在线看| 欧美日本亚洲视频在线播放| 少妇熟女欧美另类| 亚洲精品久久久久久婷婷小说 | 亚洲国产精品成人久久小说 | 日韩亚洲欧美综合| 少妇高潮的动态图| videossex国产| 波野结衣二区三区在线| 日韩av在线大香蕉| 成人亚洲欧美一区二区av| 日韩av在线大香蕉| 国产精品乱码一区二三区的特点| 男女视频在线观看网站免费| 老熟妇乱子伦视频在线观看| 成人无遮挡网站| 免费av不卡在线播放| 搡老妇女老女人老熟妇| 亚洲成人av在线免费| 毛片一级片免费看久久久久| 在线天堂最新版资源| 嫩草影院入口| 亚洲人成网站在线播放欧美日韩| 麻豆成人午夜福利视频| 在线观看av片永久免费下载| 精品无人区乱码1区二区| 亚洲精品日韩在线中文字幕 | 国产极品天堂在线| 99热6这里只有精品| 亚洲国产精品sss在线观看| 一级av片app| 神马国产精品三级电影在线观看| 又爽又黄无遮挡网站| 亚洲欧美精品专区久久| 亚洲国产精品久久男人天堂| 精品国产三级普通话版| 国产在视频线在精品| 久久精品影院6| 99热精品在线国产| 国产精品女同一区二区软件| 国内精品一区二区在线观看| 国产精品精品国产色婷婷| 不卡视频在线观看欧美| 丝袜美腿在线中文| 亚洲第一区二区三区不卡| 亚洲欧美日韩高清专用| 国产乱人视频| 日本黄色片子视频| 亚洲精品成人久久久久久| 最近最新中文字幕大全电影3| 少妇丰满av| 丝袜喷水一区| 激情 狠狠 欧美| 久99久视频精品免费| 免费看a级黄色片| 校园春色视频在线观看| 亚洲国产欧美在线一区| 亚洲成av人片在线播放无| 99国产精品一区二区蜜桃av| 亚洲欧美精品自产自拍| 99热这里只有是精品在线观看| 少妇丰满av| 99热这里只有是精品50| 狂野欧美激情性xxxx在线观看| 成人高潮视频无遮挡免费网站| 美女内射精品一级片tv| 悠悠久久av| 狂野欧美白嫩少妇大欣赏| av免费观看日本| 亚洲一区高清亚洲精品| 欧美日韩精品成人综合77777| 国产综合懂色| 九色成人免费人妻av| 国产精品99久久久久久久久| 亚洲七黄色美女视频| 成人特级黄色片久久久久久久| 日本欧美国产在线视频| 在线观看66精品国产| 日韩亚洲欧美综合| 22中文网久久字幕| 热99re8久久精品国产| 黄片无遮挡物在线观看| 韩国av在线不卡| 亚洲精品乱码久久久久久按摩| 菩萨蛮人人尽说江南好唐韦庄 | 波野结衣二区三区在线| 午夜亚洲福利在线播放| 网址你懂的国产日韩在线| 国产精品嫩草影院av在线观看| 成年版毛片免费区| 免费观看的影片在线观看| 日本五十路高清| 久久国内精品自在自线图片| av在线蜜桃| 国产成人aa在线观看| av.在线天堂| 国产成人aa在线观看| 美女脱内裤让男人舔精品视频 | 蜜桃久久精品国产亚洲av| 五月玫瑰六月丁香| 免费av观看视频| 亚洲av一区综合| 国产熟女欧美一区二区| 日韩亚洲欧美综合| 91av网一区二区| 99久久成人亚洲精品观看| 日本黄大片高清| 国产69精品久久久久777片| 26uuu在线亚洲综合色| 男人狂女人下面高潮的视频| 99久久精品热视频| 国内精品久久久久精免费| 久久久精品欧美日韩精品| 久久鲁丝午夜福利片| 变态另类丝袜制服| 亚洲性久久影院| 国产黄片视频在线免费观看| 特级一级黄色大片| 午夜精品在线福利| 国产精品久久久久久亚洲av鲁大| 99久国产av精品| www.色视频.com| 十八禁国产超污无遮挡网站| 色综合站精品国产| 久久久久免费精品人妻一区二区| 国产精品女同一区二区软件| 欧美最黄视频在线播放免费| 免费大片18禁| 99视频精品全部免费 在线| 国产亚洲精品久久久com| 久久久久久九九精品二区国产| 嘟嘟电影网在线观看| 久久久久久伊人网av| 久久亚洲国产成人精品v| 欧美在线一区亚洲| 亚洲最大成人中文| 毛片女人毛片| 精品不卡国产一区二区三区| 亚洲人与动物交配视频| 久久久午夜欧美精品| 欧美另类亚洲清纯唯美| 亚洲欧洲国产日韩| 啦啦啦韩国在线观看视频| 乱人视频在线观看| 精品熟女少妇av免费看| 三级毛片av免费| 少妇被粗大猛烈的视频| 国产成人a区在线观看| 日本爱情动作片www.在线观看| 欧美区成人在线视频| a级毛色黄片| 国内揄拍国产精品人妻在线| 男女那种视频在线观看| av在线播放精品| 精品不卡国产一区二区三区| 亚洲精品久久国产高清桃花| 91午夜精品亚洲一区二区三区| 99久久精品国产国产毛片| 国产精品福利在线免费观看| 99久久成人亚洲精品观看| 一区二区三区高清视频在线| 亚洲国产精品久久男人天堂| 国产成人福利小说| 亚洲最大成人中文| 国内精品宾馆在线| 男人舔女人下体高潮全视频| 青春草视频在线免费观看| 69av精品久久久久久| 全区人妻精品视频| 免费在线观看成人毛片| 国产精品av视频在线免费观看| 日本三级黄在线观看| 有码 亚洲区| 人妻制服诱惑在线中文字幕| 亚洲成人久久爱视频| 国产不卡一卡二| 赤兔流量卡办理| 国产精品久久久久久av不卡| 一区二区三区高清视频在线| 欧美极品一区二区三区四区| 高清日韩中文字幕在线| 国产亚洲精品av在线| 欧美成人一区二区免费高清观看| 亚洲真实伦在线观看| 别揉我奶头 嗯啊视频| 美女xxoo啪啪120秒动态图| 日本成人三级电影网站| 久久九九热精品免费| 久久99蜜桃精品久久| 国产精品久久久久久久久免| 少妇被粗大猛烈的视频| 欧美一级a爱片免费观看看| 亚州av有码| 日韩一区二区视频免费看| or卡值多少钱| 在线播放无遮挡| 给我免费播放毛片高清在线观看| 国产精品福利在线免费观看| 久久久久久久久久久丰满| 亚洲最大成人av| av.在线天堂| 国产精品美女特级片免费视频播放器| 亚洲精品成人久久久久久| 成人综合一区亚洲| 亚洲五月天丁香| 亚洲精品乱码久久久久久按摩|