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

    一種低速風(fēng)洞虛擬飛行試驗(yàn)裝置的建模與仿真

    2017-11-01 06:02:52郭林亮祝明紅鐘誠(chéng)文
    關(guān)鍵詞:試驗(yàn)裝置迎角風(fēng)洞

    郭林亮, 祝明紅, 傅 澔, 孔 鵬, 鐘誠(chéng)文

    (1. 西北工業(yè)大學(xué) 航空學(xué)院, 陜西 西安 710072; 2. 中國(guó)空氣動(dòng)力研究與發(fā)展中心 低速空氣動(dòng)力研究所, 四川 綿陽(yáng) 621000)

    一種低速風(fēng)洞虛擬飛行試驗(yàn)裝置的建模與仿真

    郭林亮1,*, 祝明紅2, 傅 澔2, 孔 鵬2, 鐘誠(chéng)文1

    (1. 西北工業(yè)大學(xué) 航空學(xué)院, 陜西 西安 710072; 2. 中國(guó)空氣動(dòng)力研究與發(fā)展中心 低速空氣動(dòng)力研究所, 四川 綿陽(yáng) 621000)

    為實(shí)現(xiàn)大角度范圍、多自由度的機(jī)動(dòng)動(dòng)作模擬,研發(fā)一種低速風(fēng)洞三自由度動(dòng)態(tài)試驗(yàn)支撐機(jī)構(gòu),可模擬繞速度矢滾轉(zhuǎn)機(jī)動(dòng)動(dòng)作以及失速偏離、尾旋等危險(xiǎn)飛行狀態(tài)。該機(jī)構(gòu)通過(guò)兩自由度轉(zhuǎn)臺(tái)和旋轉(zhuǎn)曲桿的組合運(yùn)動(dòng)模擬飛機(jī)模型的三軸姿態(tài)變化?;诙囿w動(dòng)力學(xué)理論,采用拉格朗日乘子法推導(dǎo)出該機(jī)構(gòu)曲桿-飛機(jī)模型的動(dòng)力學(xué)數(shù)學(xué)模型;模型中考慮了機(jī)構(gòu)與試驗(yàn)?zāi)P偷募s束關(guān)系、機(jī)構(gòu)摩擦力矩的影響。仿真結(jié)果表明:采用該支撐機(jī)構(gòu),飛機(jī)模型可在水平風(fēng)洞中實(shí)現(xiàn)繞速度矢量滾轉(zhuǎn)等典型機(jī)動(dòng)動(dòng)作;曲桿和試驗(yàn)?zāi)P偷臐L轉(zhuǎn)運(yùn)動(dòng)基本同步;曲桿主要影響速率響應(yīng)的動(dòng)態(tài)過(guò)程,摩擦力矩對(duì)速率的動(dòng)態(tài)過(guò)程和穩(wěn)態(tài)值有一定影響。以上數(shù)學(xué)建模和仿真驗(yàn)證可為風(fēng)洞試驗(yàn)提供理論依據(jù)。

    風(fēng)洞虛擬飛行;動(dòng)力學(xué)相似;三自由度動(dòng)態(tài)試驗(yàn);多體動(dòng)力學(xué);飛行仿真

    0 引 言

    目前飛行力學(xué)分析和飛行控制系統(tǒng)設(shè)計(jì)時(shí)主要采用線性疊加方法建立氣動(dòng)力模型,構(gòu)建氣動(dòng)力模型的數(shù)據(jù)主要通過(guò)靜態(tài)測(cè)力試驗(yàn)和小振幅動(dòng)導(dǎo)數(shù)試驗(yàn)獲得。在線性迎角范圍內(nèi),動(dòng)導(dǎo)數(shù)與迎角、振幅及頻率等關(guān)系不大;而在大迎角區(qū)域,迎角、振幅及頻率對(duì)動(dòng)導(dǎo)數(shù)的非線性影響明顯加大[1]?,F(xiàn)有設(shè)計(jì)流程沒(méi)有考慮這些非線性氣動(dòng)力,因此需要通過(guò)原型機(jī)飛行試驗(yàn)對(duì)氣動(dòng)力和飛行控制系統(tǒng)的性能進(jìn)行驗(yàn)證。目前獲取動(dòng)態(tài)氣動(dòng)力的試驗(yàn)裝置存在一定的局限性:動(dòng)導(dǎo)數(shù)裝置僅能在單個(gè)自由度上強(qiáng)迫振動(dòng),其幅值小、頻率調(diào)節(jié)范圍有限;旋轉(zhuǎn)天平裝置僅能實(shí)現(xiàn)繞速度矢作錐形運(yùn)動(dòng)。而研究大迎角飛行出現(xiàn)的氣動(dòng)力問(wèn)題,如機(jī)翼?yè)u滾、橫向偏離、上仰等,就需要把試驗(yàn)裝置的某些自由度放開(kāi)來(lái)復(fù)現(xiàn)某些復(fù)雜飛行現(xiàn)象[2]。比如自由搖滾、自由偏航等試驗(yàn)裝置較好地解決了搖晃、偏離等非線性飛行現(xiàn)象的模擬和研究。大迎角氣動(dòng)問(wèn)題具有非線性、縱橫向耦合強(qiáng)等特點(diǎn),因此為了盡量逼真地模擬飛行器飛行機(jī)動(dòng)動(dòng)作,新的多自由度試驗(yàn)裝置的研究一直在進(jìn)行之中。近年來(lái)發(fā)展的一些兩自由度裝置模擬了兩自由度耦合振蕩運(yùn)動(dòng)、旋轉(zhuǎn)流場(chǎng)下的振蕩運(yùn)動(dòng)狀態(tài)[3],多自由度機(jī)動(dòng)模擬裝置在風(fēng)洞中實(shí)現(xiàn)了對(duì)大飛行包線機(jī)動(dòng)動(dòng)作的模擬和多自由度耦合運(yùn)動(dòng)的研究[4-5]。

    另外,近年來(lái)發(fā)展的風(fēng)洞虛擬飛行試驗(yàn)技術(shù)[6-8]不僅要求多自由度耦合、模型自由轉(zhuǎn)動(dòng)、流場(chǎng)環(huán)境逼真;還要求在模型上集成飛行控制設(shè)備,使模型運(yùn)動(dòng)具備實(shí)時(shí)可控的能力[9-10]。虛擬飛行試驗(yàn)技術(shù)中模型支撐機(jī)構(gòu)設(shè)計(jì)是其中的關(guān)鍵技術(shù)之一。支撐機(jī)構(gòu)要求氣動(dòng)干擾小、摩擦影響小、角度范圍大、模擬更逼真。美國(guó)Magill等采用6根張線和球鉸的組合方式支撐BOA導(dǎo)彈進(jìn)行試驗(yàn),結(jié)果與飛行試驗(yàn)結(jié)果吻合較好[11-12]。英國(guó)Lowenberg等發(fā)展了單自由度、兩自由度、三自由度和五自由度等機(jī)動(dòng)試驗(yàn)裝置[13-15],開(kāi)展了支桿-飛機(jī)多體動(dòng)力學(xué)建模、非線性氣動(dòng)力、參數(shù)辨識(shí)等方面的研究,獲得了試驗(yàn)裝置的動(dòng)力學(xué)模型描述,以及Hawk模型和M2370模型三軸的靜、動(dòng)導(dǎo)數(shù)數(shù)據(jù)[16-17]。俄羅斯Sohi采用腹撐三自由度機(jī)構(gòu)在水平風(fēng)洞研究穩(wěn)定尾旋性能,與立式風(fēng)洞試驗(yàn)結(jié)果較為吻合[18]。TsAGI采用背撐三自由度機(jī)構(gòu)開(kāi)展了機(jī)翼?yè)u滾問(wèn)題的研究,通過(guò)H∞控制技術(shù)抑制了試驗(yàn)?zāi)P偷臋C(jī)翼?yè)u滾現(xiàn)象[19-20]。中國(guó)空氣動(dòng)力研究與發(fā)展中心在2.4 m跨聲速風(fēng)洞采用縱向單吊臂的支撐方式開(kāi)展了某導(dǎo)彈虛擬飛行試驗(yàn),驗(yàn)證了俯仰/滾轉(zhuǎn)耦合運(yùn)動(dòng)的解耦控制方法[21]。

    綜合分析以上動(dòng)態(tài)試驗(yàn)的支撐方式可以看出,采用直支桿腹撐方式的三自由度裝置迎角模擬范圍、滾轉(zhuǎn)運(yùn)動(dòng)范圍有限,難以滿足大迎角氣動(dòng)/運(yùn)動(dòng)耦合特性的研究。本文提出了一種三自由度裝置,能夠?qū)崿F(xiàn)大角度范圍的機(jī)動(dòng)動(dòng)作模擬。針對(duì)此支撐裝置,采用拉格朗日乘子法導(dǎo)出該虛擬飛行試驗(yàn)裝置的動(dòng)力學(xué)數(shù)學(xué)模型,并進(jìn)行典型激勵(lì)下縱向、橫航向動(dòng)力學(xué)響應(yīng)以及大迎角混合操縱響應(yīng)的非線性仿真與分析。

    1 試驗(yàn)裝置

    1.1裝置介紹

    圖1所示的三自由度虛擬飛行試驗(yàn)裝置通過(guò)置于模型內(nèi)部的萬(wàn)向鉸實(shí)現(xiàn)飛機(jī)模型的俯仰和偏航運(yùn)動(dòng),通過(guò)與旋轉(zhuǎn)曲桿一起轉(zhuǎn)動(dòng)實(shí)現(xiàn)飛機(jī)模型的滾轉(zhuǎn)運(yùn)動(dòng)。旋轉(zhuǎn)曲桿經(jīng)過(guò)三次預(yù)彎,一是為了擴(kuò)展模型的運(yùn)動(dòng)范圍,二是為了減小曲桿慣量對(duì)系統(tǒng)的影響。旋轉(zhuǎn)曲桿與模型腹部連接的部分進(jìn)行了第三次預(yù)彎,其偏斜角即該部分曲桿軸線(圖1中虛線所示的偏航軸)與水平面的夾角為45°,這樣可以減小大迎角狀態(tài)下模型的腹部開(kāi)孔。旋轉(zhuǎn)曲桿上方設(shè)計(jì)了具有流線外形的配重,其高度可調(diào)節(jié),可使曲桿重心在曲桿的旋轉(zhuǎn)軸上以保持動(dòng)平衡。試驗(yàn)?zāi)P偷娜S姿態(tài)及角速率可以通過(guò)機(jī)載的航姿參考系統(tǒng)進(jìn)行測(cè)量,迎角、側(cè)滑角可利用風(fēng)洞來(lái)流方向不變的特點(diǎn)由三軸姿態(tài)推算導(dǎo)出。模型的升降舵、副翼、方向舵等舵面與舵機(jī)連接采用平行四邊形傳動(dòng)方式;平行于舵面轉(zhuǎn)軸在機(jī)身或機(jī)翼內(nèi)部布置舵機(jī)轉(zhuǎn)軸,舵面和舵機(jī)轉(zhuǎn)軸上設(shè)計(jì)搖臂,通過(guò)推桿連接搖臂推動(dòng)舵面轉(zhuǎn)動(dòng)。該裝置能實(shí)現(xiàn)接近90°迎角的極限飛行動(dòng)作模擬,可開(kāi)展失速偏離、尾旋初始階段和穩(wěn)定尾旋的研究。該裝置的主要的幾何尺寸見(jiàn)圖2 。

    1.2技術(shù)分析

    飛機(jī)的三自由度角運(yùn)動(dòng)是快變化過(guò)程,相對(duì)慢變化的線運(yùn)動(dòng),其與飛行安全的關(guān)系更加緊密。現(xiàn)目前立式風(fēng)洞尾旋試驗(yàn)、水平風(fēng)洞自由飛試驗(yàn)中模型的線位移十分有限,其實(shí)質(zhì)也主要為三自由度角運(yùn)動(dòng)。因此,三自由度動(dòng)態(tài)試驗(yàn)裝置可以捕獲飛機(jī)三軸轉(zhuǎn)動(dòng)運(yùn)動(dòng)的主要特征,以達(dá)到開(kāi)展氣動(dòng)/運(yùn)動(dòng)耦合、多軸運(yùn)動(dòng)耦合研究的目的。

    采用常規(guī)直支桿腹撐或背撐方式實(shí)現(xiàn)大角度范圍的運(yùn)動(dòng)模擬,需要在模型表面開(kāi)設(shè)較大的孔,這樣將嚴(yán)重影響模型的氣動(dòng)特性。而采取翼型主支桿和預(yù)彎曲桿的方式,既大幅減小了對(duì)后方流場(chǎng)的干擾,又可以很好控制對(duì)模型外形的破壞。當(dāng)然,由于支撐的存在會(huì)對(duì)模型的氣動(dòng)力產(chǎn)生一定的干擾,比如支撐會(huì)在模型區(qū)誘導(dǎo)一定的上洗氣流,從而對(duì)模型Cm0的產(chǎn)生一個(gè)平移量等。支架干擾的貢獻(xiàn)量可以通過(guò)CFD計(jì)算模擬或風(fēng)洞支架干擾試驗(yàn)等手段獲取,并在氣動(dòng)力建模中予以考慮。

    繞速度矢橫滾是飛機(jī)最基本的機(jī)動(dòng)動(dòng)作之一。常規(guī)腹撐或背撐方式無(wú)法實(shí)現(xiàn)360°連續(xù)滾轉(zhuǎn),但利用風(fēng)洞的氣流方向保持不變的特點(diǎn),通過(guò)一個(gè)小巧的曲桿與模型同步滾轉(zhuǎn),就可以巧妙地解決繞速度矢旋轉(zhuǎn)的問(wèn)題。另外,置于模型內(nèi)部的萬(wàn)向鉸可以實(shí)現(xiàn)模型迎角、側(cè)滑角的連續(xù)變化。這樣,在可控的風(fēng)洞流場(chǎng)環(huán)境下就實(shí)現(xiàn)了物理意義明晰的三自由度運(yùn)動(dòng)。

    2 數(shù)學(xué)模型

    2.1前提和假設(shè)

    風(fēng)洞試驗(yàn)時(shí),模型質(zhì)心位于風(fēng)洞試驗(yàn)段中心保持不變,但模型的姿態(tài)可繞三軸轉(zhuǎn)動(dòng),因此,在動(dòng)力學(xué)建模時(shí)可作如下假設(shè):

    (1) 假設(shè)飛機(jī)模型及試驗(yàn)裝置均為剛體,不存在變形情況;

    (2) 假定曲桿所受氣動(dòng)力為零,曲桿的旋轉(zhuǎn)軸平行于風(fēng)洞來(lái)流方向;

    (3) 模型的轉(zhuǎn)動(dòng)中心為其重心,且模型和曲桿的重心在曲桿的旋轉(zhuǎn)軸上;

    (4) 試驗(yàn)?zāi)P偷膽T量和旋轉(zhuǎn)速率較小,可忽略陀螺力矩的影響。

    2.2軸系定義

    在系統(tǒng)動(dòng)力學(xué)分析時(shí)需涉及三種坐標(biāo)系(如圖3所示),分別為:

    地軸系oxgygzg,簡(jiǎn)稱Sg,其ozg沿鉛垂方向向下,oxg在水平面內(nèi),并與風(fēng)洞軸線平行,與來(lái)流方向相反,oyg與平面oxgzg垂直,指向右。

    模型體軸系oxayaza,簡(jiǎn)稱Sa,原點(diǎn)位于模型重心,oxa軸在飛行器對(duì)稱平面內(nèi),平行于機(jī)身軸線或機(jī)翼的平均氣動(dòng)弦線,指向前;oza位于對(duì)稱面內(nèi),垂直于oxa,指向下;oya垂直于對(duì)稱平面,指向右。

    曲桿體軸系oxryrzr,簡(jiǎn)稱Sr,原點(diǎn)位于曲桿重心,軸系定義與Sa類似。

    用φr、θr、ψr分別表示曲桿滾轉(zhuǎn)角、俯仰角、偏航角,導(dǎo)出Sr與Sg的轉(zhuǎn)換矩陣Tgr、Trg。

    2.3基于絕對(duì)坐標(biāo)方法的動(dòng)力學(xué)方程

    2.3.1 絕對(duì)坐標(biāo)

    利用絕對(duì)坐標(biāo)方法進(jìn)行動(dòng)力學(xué)分析,首先假設(shè)系統(tǒng)內(nèi)所有剛體為無(wú)約束的自由剛體,以各剛體的質(zhì)心笛卡爾坐標(biāo)和繞質(zhì)心轉(zhuǎn)動(dòng)的角度坐標(biāo)或歐拉參數(shù)作為系統(tǒng)的絕對(duì)坐標(biāo),對(duì)各剛體建立無(wú)約束的動(dòng)力學(xué)方程,再利用拉格朗日乘子法將其與系統(tǒng)鉸約束方程聯(lián)立,構(gòu)成完整的動(dòng)力學(xué)方程。

    系統(tǒng)的絕對(duì)坐標(biāo)為:

    剛體基于自身體軸的角速度向量為:

    ωi=(piqiri)T,i=r,a

    2.3.2 動(dòng)力學(xué)方程的一般形式

    確定坐標(biāo)后,根據(jù)拉格朗日乘子法,可以給出受約束的動(dòng)力學(xué)方程一般形式如下[22]:

    其中,A是與無(wú)約束動(dòng)力方程有關(guān)的系數(shù)矩陣,B是該方程組右側(cè)向量,Φq為約束方程Φ(q)的Jacobian矩陣,λ為拉格朗日乘子,ζ為加速度約束方程組右側(cè)向量。下面將給出該動(dòng)力學(xué)方程的推導(dǎo)過(guò)程。

    2.3.3 無(wú)約束的動(dòng)力學(xué)方程

    一般情況下,剛體Ri繞質(zhì)心轉(zhuǎn)動(dòng)的運(yùn)動(dòng)方程為如下形式:

    式中ωabi為剛體Ri的絕對(duì)角速度,其可以表示為:

    式中ω0為多體系統(tǒng)零剛體的角速度,系統(tǒng)選定支撐座為零剛體,因此ω0=0,式(5)可簡(jiǎn)化為:

    將上式代入式(4)后得到:

    此處引入反對(duì)稱陣的概念:設(shè)有向量a=(a1a2a3)T,則a的反對(duì)稱陣為:

    則式(7)可以整理為:

    2.3.4 約束條件

    因此,本系統(tǒng)的幾何約束方程組Φ(q)為:

    上述約束條件為不顯含時(shí)間的定常約束,上式對(duì)t求導(dǎo)即可得到速度約束方程:

    式中Φq為約束方程Φ(q)對(duì)坐標(biāo)q的雅可比矩陣。

    速度約束方程再次對(duì)時(shí)間求導(dǎo),得到加速度約束方程:

    2.3.5 拉格朗日乘子法

    引入3個(gè)拉格朗日乘子λ=(λ1λ2λ3)T,將上述約束方程與相同標(biāo)號(hào)的拉格朗日乘子相乘,加入動(dòng)力學(xué)方程中,得到下式:

    與加速度約束方程聯(lián)立,得到受約束的動(dòng)力學(xué)方程:

    2.4力和力矩

    2.4.1 氣動(dòng)力模型

    本文基于常規(guī)風(fēng)洞試驗(yàn)得到的基本氣動(dòng)力、舵效及動(dòng)導(dǎo)數(shù)等數(shù)據(jù)建立了如下的氣動(dòng)力模型,各系數(shù)均在體軸系下描述,氣動(dòng)力模型中支架干擾的貢獻(xiàn)量可通過(guò)后續(xù)補(bǔ)充風(fēng)洞試驗(yàn)獲得。該飛機(jī)模型舵面全零狀態(tài)下的基本氣動(dòng)力矩見(jiàn)圖4所示。

    (18)

    式中V為來(lái)流速度,δe、δr、δf、δa分別表示模型升降舵、方向舵、襟翼和副翼角度,各氣動(dòng)系數(shù)按照基本量疊加增量的方式計(jì)算。以俯仰力矩系數(shù)Cm為例,按照下式進(jìn)行計(jì)算:

    Cm=Cmbasic(α,δf)+ΔCmβ(α,β,δf)+ΔCmsup(α,β)+

    2.4.2 摩擦力矩計(jì)算

    式中,σ0為剛性系數(shù),N·m/rad;σ1為阻尼系數(shù),σ2為粘性系數(shù),N·m·s/rad;Z為剛毛平均形變,Ts為靜摩擦力矩,Tc為庫(kù)侖摩擦力矩,N·m;參數(shù)αf的引入使該模型符合Stribeck效應(yīng),單位為s/m;參數(shù)νd用以減弱高轉(zhuǎn)速下剛毛偏轉(zhuǎn)的影響,并且保證了摩擦力矩模型的鈍性,單位為rad/s。

    當(dāng)系統(tǒng)處于穩(wěn)態(tài)時(shí),dZj/dt=0,此時(shí)摩擦力矩模型可以簡(jiǎn)化為式(23)。圖5給出了該裝置俯仰運(yùn)動(dòng)方向摩擦特性的地面測(cè)試結(jié)果。

    (23)

    3 仿真與結(jié)果分析

    3.1仿真程序

    仿真程序在Matlab環(huán)境下編寫,采用Nelder-Mead算法在給定初始狀態(tài)下配平飛機(jī);導(dǎo)出約束方程的Jacobian矩陣的解析形式,與無(wú)約束的動(dòng)力學(xué)方程聯(lián)立構(gòu)成動(dòng)力學(xué)模型,然后采用無(wú)違約算法進(jìn)行求解[25-27];采用多維線性插值方法獲取氣動(dòng)力,根據(jù)LuGre模型解算摩擦力矩。圖6給出了計(jì)算程序的相關(guān)模塊和主要流程。

    為了分析曲桿、摩擦力矩等因素對(duì)飛機(jī)模型響應(yīng)特性的影響規(guī)律,選取三種工況進(jìn)行對(duì)比研究,即Ⅰ.僅試驗(yàn)?zāi)P偷娜杂啥冗\(yùn)動(dòng)(無(wú)摩擦), Ⅱ.試驗(yàn)?zāi)P秃颓鷹U的三自由度運(yùn)動(dòng)(無(wú)摩擦), Ⅲ.試驗(yàn)?zāi)P秃颓鷹U的三自由度運(yùn)動(dòng)(有摩擦)。通過(guò)Ⅰ和Ⅱ的對(duì)比,可以獲得曲桿對(duì)系統(tǒng)的影響規(guī)律;通過(guò)Ⅱ和Ⅲ的對(duì)比,可以獲得摩擦對(duì)系統(tǒng)的影響規(guī)律。

    3.2算法驗(yàn)證

    Adams是美國(guó)MSC公司開(kāi)發(fā)的可用于開(kāi)展多體動(dòng)力學(xué)仿真的成熟軟件。為驗(yàn)證前述動(dòng)力學(xué)模型與算法的正確性,利用Adams-Matlab聯(lián)合仿真功能進(jìn)行確認(rèn)。采用較為常用的斜拉桿操縱進(jìn)行對(duì)比(升降舵-10°,副翼-3°),兩種方法的對(duì)比驗(yàn)證結(jié)果見(jiàn)圖7,圖中Rig為本文算法的結(jié)果。不同算法下,縱向響應(yīng)的一致性很好;橫航向響應(yīng)基本吻合,曲線形態(tài)相似,只是其動(dòng)態(tài)振蕩過(guò)程略有差異??偟膩?lái)說(shuō),不同算法的差異可以接受,達(dá)到了算法驗(yàn)證的目的。后文在此基礎(chǔ)上利用本文的方法對(duì)試驗(yàn)裝置進(jìn)行詳細(xì)分析。

    3.3典型操縱

    為了便于分析,針對(duì)縱、橫、航向分別采用典型的階躍信號(hào)激勵(lì)飛機(jī)的動(dòng)態(tài)響應(yīng);縱向、橫向及航向的操縱信號(hào)分別為1°升降舵階躍、1°副翼階躍和1°方向舵階躍信號(hào)。同時(shí)給出一組大迎角混合操縱的響應(yīng)結(jié)果,各舵面輸入信號(hào)為:在1 s時(shí)輸入-30°升降舵階躍,4 s時(shí)輸入6°方向舵階躍,12 s時(shí)輸入10°副翼階躍。飛機(jī)模型的機(jī)翼后掠角約30°,其主要參數(shù)見(jiàn)表1。仿真初始狀態(tài)的主要參數(shù)見(jiàn)表2。圖8~圖11給出了相應(yīng)的動(dòng)力學(xué)響應(yīng)仿真結(jié)果。圖中標(biāo)注3-DOF為基于飛機(jī)模型三自由度轉(zhuǎn)動(dòng)方程無(wú)桿狀態(tài)的仿真結(jié)果,Rig為飛機(jī)模型在虛擬飛行試驗(yàn)裝置上的仿真結(jié)果,Rig+Friction為在該裝置基礎(chǔ)上增加摩擦力影響的仿真結(jié)果。

    表1 模型主要參數(shù)Table 1 Major parameters of the test model

    表2 初始狀態(tài)主要參數(shù)Table 2 Main initial state parameters

    3.4縱向階躍操縱響應(yīng)分析

    升降舵操縱的響應(yīng)結(jié)果如圖8所示,有桿、無(wú)桿仿真結(jié)果吻合較好。從理論分析,純縱向運(yùn)動(dòng)中曲桿不會(huì)對(duì)俯仰運(yùn)動(dòng)的響應(yīng)產(chǎn)生影響,因此模型與機(jī)構(gòu)均無(wú)滾轉(zhuǎn)、偏航等橫航向運(yùn)動(dòng)。

    但摩擦力矩對(duì)俯仰角、俯仰速率均有一定影響。對(duì)俯仰速率的影響主要體現(xiàn)在阻尼作用和死區(qū)現(xiàn)象。阻尼作用體現(xiàn)在俯仰速率峰值變小,在第一、二個(gè)峰值分別比無(wú)摩擦仿真結(jié)果減少0.11°/s和0.25°/s;死區(qū)現(xiàn)象表現(xiàn)為俯仰速率曲線在第二次過(guò)零時(shí)(約1.8 s)不再振蕩衰減,而是自此時(shí)以后一直為零。與之相應(yīng)的,俯仰角響應(yīng)峰值略有減小,并在1.8 s后不再變化,提前0.7s進(jìn)入穩(wěn)定狀態(tài),其穩(wěn)態(tài)值與無(wú)摩擦狀態(tài)相比約有0.02°的微小差量。

    3.5橫向階躍操縱響應(yīng)分析

    副翼操縱的響應(yīng)結(jié)果如圖9所示,曲桿使模型滾轉(zhuǎn)速率、偏航速率響應(yīng)變慢,速率響應(yīng)振蕩稍平緩,穩(wěn)態(tài)值差異較??;曲桿帶來(lái)模型滾轉(zhuǎn)角響應(yīng)差量約0.8°(仿真結(jié)束時(shí))。模型側(cè)滑角、偏航角響應(yīng)差異不大。曲桿滾轉(zhuǎn)角與模型滾轉(zhuǎn)角的響應(yīng)基本同步一致,最大差量約0.3°。

    摩擦力矩對(duì)速率響應(yīng)的動(dòng)態(tài)過(guò)程和穩(wěn)態(tài)值均有一定影響。模型滾轉(zhuǎn)速率峰值較無(wú)摩擦結(jié)果小0.4°/s,穩(wěn)態(tài)值存在0.4°/s的差量;偏航速率峰值較無(wú)摩擦結(jié)果相差約0.13°/s,穩(wěn)態(tài)值差量為0.04°/s。由于模型滾轉(zhuǎn)角由速率積分決定,有無(wú)摩擦兩種狀態(tài)滾轉(zhuǎn)角響應(yīng)出現(xiàn)明顯的斜率差異,滾轉(zhuǎn)角最大差量達(dá)3.4°;模型偏航角最大相差約0.18°,側(cè)滑角穩(wěn)態(tài)值相差約0.05°。

    3.6航向階躍操縱響應(yīng)分析

    方向舵操縱的響應(yīng)結(jié)果如圖10所示,三種狀態(tài)下總體響應(yīng)趨勢(shì)一致,但也存在一定差異。曲桿使偏航速率振蕩峰值先減小后增加,峰值時(shí)間略有推遲,穩(wěn)態(tài)值變化不大;與之對(duì)應(yīng),曲桿使側(cè)滑角、偏航角響應(yīng)的振蕩幅值也略增大。此種差異主要來(lái)自于萬(wàn)向鉸偏航旋轉(zhuǎn)軸的偏斜,曲桿慣量的影響次之。

    曲桿使模型滾轉(zhuǎn)速率峰值減小,峰值時(shí)間稍推遲,穩(wěn)態(tài)值基本不變;但操縱初始階段響應(yīng)曲線形態(tài)差異明顯。在操縱后0.2 s內(nèi)無(wú)桿狀態(tài)下模型先正滾轉(zhuǎn)后減速變?yōu)樨?fù)滾轉(zhuǎn),而有桿狀態(tài)下模型直接表現(xiàn)為負(fù)滾轉(zhuǎn)。模型滾轉(zhuǎn)角響應(yīng)的差異與速率差異吻合,二者在仿真結(jié)束時(shí)存在0.26°的差量。另一方面模型滾轉(zhuǎn)角與曲桿滾轉(zhuǎn)角在操縱初始時(shí)運(yùn)動(dòng)趨勢(shì)相反,曲桿滾轉(zhuǎn)角正滾轉(zhuǎn)最大達(dá)0.8°。結(jié)合加速度曲線分析發(fā)現(xiàn),曲桿的存在使?jié)L轉(zhuǎn)運(yùn)動(dòng)對(duì)應(yīng)慣量增大,導(dǎo)致滾轉(zhuǎn)模態(tài)時(shí)間常數(shù)變大,使得滾轉(zhuǎn)角速率響應(yīng)變慢,加上萬(wàn)向鉸偏航旋轉(zhuǎn)軸偏斜的影響,導(dǎo)致了滾轉(zhuǎn)速率在1~5 s的較大差異;而萬(wàn)向鉸偏航旋轉(zhuǎn)軸的偏斜使偏航運(yùn)動(dòng)分解為繞該軸的旋轉(zhuǎn)和繞速度矢的旋轉(zhuǎn),后者造成了曲桿與模型在初始響應(yīng)時(shí)滾轉(zhuǎn)方向不一致。

    摩擦力矩的影響主要表現(xiàn)為阻尼作用,模型的偏航速率和滾轉(zhuǎn)速率收斂更快,其中偏航速率在3.6~3.8 s時(shí)呈現(xiàn)了死區(qū)現(xiàn)象,偏航速率穩(wěn)態(tài)值差量約0.04°/s,使偏航角斜率略有差異,在仿真結(jié)束10 s時(shí)差量約0.3°,側(cè)滑角差量約0.05°;滾轉(zhuǎn)速率穩(wěn)態(tài)值與無(wú)摩擦結(jié)果存在約0.4°/s的差量,使?jié)L轉(zhuǎn)角響應(yīng)的斜率呈現(xiàn)明顯差別,仿真結(jié)束時(shí)滾轉(zhuǎn)角相差約3.2°。

    3.7大迎角混合操縱響應(yīng)分析

    如圖11所示,輸入升降舵信號(hào)后,模型快速抬頭,迎角在1.3 s時(shí)達(dá)到峰值,其后迅速收斂穩(wěn)定在28.6°;輸入方向舵信號(hào)后,模型4.6 s時(shí)達(dá)到側(cè)滑角峰值后略有收斂但振蕩劇烈,均值在8°左右,同時(shí)產(chǎn)生約-49°/s的滾轉(zhuǎn)速率和-31°/s的偏航速率,模型開(kāi)始繞速度矢連續(xù)滾轉(zhuǎn);輸入副翼信號(hào)后,滾轉(zhuǎn)速率和偏航速率繼續(xù)負(fù)增長(zhǎng),即模型繞速度矢加速滾轉(zhuǎn)。

    與前面小迎角狀態(tài)相比,模型的響應(yīng)特性有所不同。主要體現(xiàn)在:一是縱向和橫航向之間的耦合更為嚴(yán)重。操縱方向舵后,模型迎角、俯仰速率出現(xiàn)小幅振蕩,且振蕩形態(tài)與側(cè)滑角類似;這主要是由于模型隨側(cè)滑變大在俯仰方向會(huì)產(chǎn)生一定的上仰力矩(見(jiàn)圖4、圖11)。二是典型模態(tài)的響應(yīng)特性和操縱性發(fā)生變化??v向運(yùn)動(dòng)方面,小迎角下單位舵偏產(chǎn)生的迎角穩(wěn)態(tài)響應(yīng)為1.1°,而大迎角狀態(tài)下為0.73°,與升降舵效率隨迎角增大而降低的規(guī)律吻合。航向運(yùn)動(dòng)方面,小迎角下單位舵偏的側(cè)滑穩(wěn)態(tài)響應(yīng)為1.4°,大迎角狀態(tài)下為1.2°,量值相當(dāng);但大迎角狀態(tài)下模型荷蘭滾模態(tài)響應(yīng)特性惡化、振蕩加劇,這是由于隨迎角增大航向穩(wěn)定性降低、動(dòng)態(tài)阻尼下降所致。橫向運(yùn)動(dòng)方面,小迎角下單位舵偏滾轉(zhuǎn)速率的穩(wěn)態(tài)響應(yīng)為-8.7°/s,而大迎角狀態(tài)下為-5.4°/s,與副翼舵效、滾轉(zhuǎn)阻尼隨迎角和側(cè)滑角增加而顯著降低的規(guī)律相符。

    大迎角狀態(tài)下,機(jī)構(gòu)摩擦對(duì)響應(yīng)結(jié)果的穩(wěn)態(tài)值有微小影響,其影響規(guī)律與小迎角狀態(tài)響應(yīng)結(jié)果一致。從圖中還發(fā)現(xiàn),大迎角狀態(tài)下曲桿對(duì)模型航向操縱響應(yīng)的影響明顯增大,其差異仍來(lái)自于萬(wàn)向鉸偏航旋轉(zhuǎn)軸的偏斜以及曲桿慣量,但與小迎角狀態(tài)不同的是,此時(shí)曲桿慣量為主要影響因素,偏航旋轉(zhuǎn)軸的偏斜次之。這是因?yàn)椋浩叫D(zhuǎn)軸的偏斜主要影響偏航運(yùn)動(dòng),曲桿慣量主要影響滾轉(zhuǎn)運(yùn)動(dòng),并且隨著模型迎角增大,曲桿慣量對(duì)于滾轉(zhuǎn)運(yùn)動(dòng)的影響加大。而在航向操縱響應(yīng)中偏航運(yùn)動(dòng)和滾轉(zhuǎn)運(yùn)動(dòng)又是相互耦合的,最終就導(dǎo)致了側(cè)滑響應(yīng)的明顯差異。目前設(shè)計(jì)的曲桿慣量約占飛機(jī)滾轉(zhuǎn)慣量的20%,小迎角狀態(tài)下差異不大;但在大迎角狀態(tài),側(cè)滑響應(yīng)峰值相差約30%。因此曲桿慣量最好在20%以內(nèi)甚至更小,以保證試驗(yàn)?zāi)M的合理性和相似性。

    4 結(jié) 論

    本文基于絕對(duì)坐標(biāo)方法建立了一種虛擬飛行試驗(yàn)裝置的動(dòng)力學(xué)模型;基于該模型所獲得的典型操縱響應(yīng)結(jié)果趨勢(shì)合理,量值正確,表明該動(dòng)力學(xué)模型可從理論上有效指導(dǎo)該類試驗(yàn)裝置的設(shè)計(jì)與實(shí)現(xiàn)。通過(guò)仿真分析得到:

    1) 曲桿對(duì)試驗(yàn)?zāi)P蜐L轉(zhuǎn)角速率響應(yīng)變慢;對(duì)偏航運(yùn)動(dòng)的振蕩過(guò)程有一定的增強(qiáng)效果,該影響隨模型迎角增大而變大;但曲桿對(duì)穩(wěn)態(tài)響應(yīng)影響不大。因此建議旋轉(zhuǎn)曲桿應(yīng)盡量采用輕質(zhì)材料制作,其慣量最好控制在飛機(jī)滾轉(zhuǎn)慣量的20%以內(nèi),以降低對(duì)模型操縱響應(yīng)的影響。

    2) 摩擦力矩對(duì)短周期、荷蘭滾等模態(tài)的動(dòng)態(tài)過(guò)程和各特征量的穩(wěn)態(tài)值有一定影響,在低旋轉(zhuǎn)速度下可能產(chǎn)生死區(qū)現(xiàn)象。因此該裝置應(yīng)采用摩擦貢獻(xiàn)較小的精密軸承,以盡量減小摩擦力矩對(duì)模型操縱響應(yīng)的影響。

    本文提出的裝置尚未開(kāi)展風(fēng)洞試驗(yàn);后續(xù)將對(duì)現(xiàn)有試驗(yàn)裝置和模型進(jìn)行適應(yīng)性改造,以開(kāi)展大迎角風(fēng)洞虛擬飛行試驗(yàn)進(jìn)一步驗(yàn)證數(shù)學(xué)模型和算法。

    [1]Vicroy D D, Loeser T D, Schutte A. SACCON forced oscillation tests at DNW-NWB and NASA Langley 14x22-foot tunnel, AIAA-2010-4394[R]. Illinois: AIAA, 2010.

    [2]Liu W, Yang X L, Zhang H X, et al. A review on investigations of wing rock problems under high angles of attack[J]. Advances in Mechanics, 2008, 38(2): 214-228. (in Chinese)

    劉偉, 楊小亮, 張涵信, 等. 大攻角運(yùn)動(dòng)時(shí)的機(jī)翼?yè)u滾問(wèn)題研究綜述[J]. 力學(xué)進(jìn)展, 2008, 38(2): 214-228.

    [3]Bu C, Du X Q, Huang L J, et al. Investigation of unsteady aerodynamic characteristics for the large amplitude rolling under rotary flow field[J]. Journal of Experiments in Fluid Mechanics, 2008, 22(1): 46-54. (in Chinese)

    卜忱, 杜希奇, 黃麗婧, 等. 旋轉(zhuǎn)流場(chǎng)下飛機(jī)大幅滾轉(zhuǎn)振蕩時(shí)的動(dòng)態(tài)橫向氣動(dòng)特性實(shí)驗(yàn)研究[J]. 實(shí)驗(yàn)流體力學(xué), 2008, 22(1): 46-54.

    [4]Huebner A R, Bergmann A, Loeser T. Experimental and numerical investigations of unsteady force and pressure distributions of moving transport aircraft configurations, AIAA-2009-0091[R]. Florida: AIAA, 2009.

    [5]Xie Z J, Sun X Y, Sun H S, et al. Mechanism design and dynamics analysis of high speed parallel robot for dynamic test in low speed wind tunnel[J]. Acta Aeronautica et Astronautica Sinica, 2013, 34(3): 487-494. (in Chinese)

    謝志江, 孫小勇, 孫海生, 等. 低速風(fēng)洞動(dòng)態(tài)試驗(yàn)的高速并聯(lián)機(jī)構(gòu)設(shè)計(jì)及動(dòng)力學(xué)分析[J]. 航空學(xué)報(bào), 2013, 47(5): 487-494.

    [6]Huang M, Wang Z W. A review of wind tunnel based virtual flight testing techniques for evaluation of flight control systems[J]. International Journal of Aerospace Engineer, 2015, 2015(1): 1-22.

    [7]Sidoryuk M E. Robust control design to suppress wing rock motion of a wind-tunnel aircraft model in 3DOF gimbals. On possibility of critical flight regime study in wind tunnels using three-degree-of-freedom gimbals[J]. TsAGI Science Journal, 2014, 45(8): 977-992.

    [8]Stenfelt G, Ringertz U. Yaw control of a tailless aircraft configuration[J]. Journal of Aircraft, 2010, 47(5): 1807-1810.

    [9]Strub G, Theodoulisy S, Gassman V, et al. Pitch axis control for a guided projectile in a wind tunnel-based hardware-in-the-loop setup, AIAA-2015-0153[R]. Florida, AIAA, 2015.

    [10]Ignatyev D I, Zaripov K G. Wind tunnel tests for validation of control algorithms at high angles of attack using autonomous aircraft model mounted in 3DOF gimbals, AIAA-2016-3106[R]. Washington: AIAA, 2016.

    [11]Lawrence F C, Mills B H. Status update of the AEDC virtual flight testing development program, AIAA-2002-0168[R]. Reston: AIAA, 2002.

    [12]Magill J C, Cataldi P, Morency J R, et al. Demonstration of a wire suspension for wind-tunnel virtual flight testing[J]. Journal of Spacecraft and Rockets, 2009, 46(3): 624-633.

    [13]Gatto A. Application of a pendulum support test rig for aircraft stability derivative estimation[J]. Journal of Aircraft, 2006, 46(3): 927-934.

    [14]Gatto A, Lowenberg M H. Evaluation of a three-degree-of-freedom test rig for stability derivative estimation[J]. Journal of Aircraft, 2006, 43(6): 1747-1762.

    [15]Pattinson J, Lowenberg M H, Goman M G. A multi-degree-of-freedom rig for the wind tunnel determination of dynamic data, AIAA-2009-5727[R]. Chicago: AIAA, 2009.

    [16]Pattinson J, Lowenberg M H, Goman M G. Multi-degree-of-freedom wind-tunnel maneuver rig for dynamic simulation and aerodynamic model identification[J]. Journal of Aircraft, 2013, 50(2): 551-566.

    [17]Araujo-Estrada S A, Lowenberg M H, Neild S, et al. Evaluation of aircraft model upset behavior using wind tunnel maneuver rig, AIAA-2015-0750[R]. Florida: AIAA, 2015.

    [18]Sohi N P. Modeling of spin modes of supersonic aircraft in horizontal wind tunnel[R]. Yokohama, Japan, ICAS, 2004.

    [19]Grishin I, Khrabrov A, Kolinko A, et al. Wind tunnel investigation of critical flight regimes using dynamically scaled actively controlled model in 3 DOF gimbal[R]. St. Petersburg, Russia, ICAS, 2014.

    [20]Khrabrov A N, Sidoryuk M E, Kolesnikov E N, et al. On possibility of critical flight regime study in wind tunnels using three-degree-of-freedom gimbals[J]. TsAGI Science Journal, 2014, 45(8): 825-39.

    [21]Zhao Z L, Wu J Q, Li H, et al. Investigation of virtual flight testing technique based on 2.4 m transonic wind tunnel[J]. Acta Aeronautica et Astronautica Sinica, 2016, 37(2): 504-512. (in Chinese)

    趙忠良, 吳軍強(qiáng), 李浩, 等. 2.4 m跨聲速風(fēng)洞虛擬飛行試驗(yàn)技術(shù)初步研究[J]. 航空學(xué)報(bào), 2016, 37(2): 504-512.

    [22]Liu Y Z, Pan Z K, Ge X S. Dynamics of multibody systems[M]. Beijing: High Education Press, 2014: 173-215. (in Chinese)

    劉延柱, 潘振寬, 戈新生. 多體系統(tǒng)動(dòng)力學(xué)[M]. 北京: 高等教育出版社, 2014: 173-215.

    [23]Wit C C, Olsson H, Astrom K J, et al. A new model for control systems with friction[J]. IEEE Transactions on Automatic Control, 1995, 40(3): 419-425.

    [24]Kermani M, Patel R, Moallem M. Friction identification in robotic manipulators: case studies[C]//Proceedings of the 2005 IEEE Conference on Control Applications. Toronto: IEEE, 2005: 1170-1175.

    [25]Qi Z H, Xu Y S, Fang H Q. Study on redundant constraints in multibody systems[J]. Chinese Journal of Theoretical and Applied Mechanics, 2011, 43(2): 390-399. (in Chinese)

    齊朝暉, 許永生, 方慧青. 多體系統(tǒng)中的冗余約束[J]. 力學(xué)學(xué)報(bào), 2011, 43(2): 390-399.

    [26]Ma X T, Zhai Y B, Luo S Q. Numerical method of multibody dynamics based onθ1method[J]. Chinese Journal of Theoretical and Applied Mechanics, 2011, 43(5): 931-938. (in Chinese)

    馬秀騰, 翟彥博, 羅書強(qiáng). 基于θ1方法的多體動(dòng)力學(xué)數(shù)值算法研究[J]. 力學(xué)學(xué)報(bào), 2011, 43(5): 931-938.

    [27]Ma X T, Zhai Y B, Luo S Q. New generalized-method for over determined motion equations in multibody dynamics[J]. Acta Aeronautica et Astronautica Sinica, 2012, 44(5): 948-952. (in Chinese)

    馬秀騰, 翟彥博, 羅書強(qiáng). 多體動(dòng)力學(xué)超定運(yùn)動(dòng)方程廣義-求解新算法[J]. 力學(xué)學(xué)報(bào), 2012, 44(5): 948-952.

    Modelingandsimulationforalowspeedwindtunnelvirtualflighttestrig

    GUO Linliang1,*, ZHU Minghong2, FU Hao2, KONG Peng2, ZHONG Chengwen1

    (1.SchoolofAeronautics,NorthwesternPolytechnicalUniversity,Xi’an710072,China; 2.LowSpeedAerodynamicsInstituteofChinaAerodynamicsResearchandDevelopmentCenter,Mianyang621000,China)

    In order to realize the maneuver simulation of large angle and multi-degree-of-freedom, a three degree-of-freedom dynamic test rig, which enables simulation of the velocity vector roll and uncontrolled flight condition such as stall/departure/spin, is developed to solve the aerodynamic/motion coupling issues. The rotational motion of the three axis of the aircraft model is simulated by the combined motion of the 2-DOF gimbal. Based on multi-body dynamics theory, the dynamical equation of the curved rod and aircraft model with terms accounting for friction in the gimbals is obtained with Lagrange multiplier method. The simulation results show that, the rig can provide the ability to realize the velocity vector roll, and the rotation of the curved rod is synchronous with the rolling motion of the aircraft model. The curved rod has an effect on the dynamic response of rotational rates, while the frictional torque makes a difference on both the dynamic response and steady state value of rotational motions. The mathematical model and simulation validation provides a guide for wind tunnel test.

    wind tunnel virtual flight test; dynamically similarity; 3-DOF dynamical test; multi-body dynamics; flight simulation

    V212.1

    A

    10.7638/kqdlxxb-2017.0164

    0258-1825(2017)05-0708-10

    2017-08-11;

    2017-09-28

    國(guó)家重點(diǎn)基礎(chǔ)研究發(fā)展計(jì)劃資助(2015CB755800)

    郭林亮*(1982-),男,湖北人,博士研究生,研究方向:實(shí)驗(yàn)空氣動(dòng)力學(xué),飛行動(dòng)力學(xué)與控制. E-mail: guolinliangliang@163.com

    郭林亮, 祝明紅, 傅澔, 等. 一種低速風(fēng)洞虛擬飛行試驗(yàn)裝置的建模與仿真[J]. 空氣動(dòng)力學(xué)學(xué)報(bào), 2017, 35(5): 708-717, 726.

    10.7638/kqdlxxb-2017.0164 GUO L L, ZHU M H, FU H. Modeling and simulation for a low speed wind tunnel virtual flight test rig[J]. Acta Aerodynamica Sinica, 2017, 35(5): 708-717, 726.

    猜你喜歡
    試驗(yàn)裝置迎角風(fēng)洞
    連續(xù)變迎角試驗(yàn)數(shù)據(jù)自適應(yīng)分段擬合濾波方法
    斑頭雁進(jìn)風(fēng)洞
    秋千動(dòng)載性能試驗(yàn)裝置的研制
    黃風(fēng)洞貂鼠精
    自行車車閘的試驗(yàn)裝置的概述
    基于NI cRIO平臺(tái)的脈沖燃燒風(fēng)洞控制系統(tǒng)設(shè)計(jì)
    自行車前叉組件的疲勞試驗(yàn)裝置的專利分布
    2MV陡前沿沖擊試驗(yàn)裝置同步技術(shù)研究
    失速保護(hù)系統(tǒng)迎角零向跳變研究
    科技傳播(2014年4期)2014-12-02 01:59:42
    飛行器風(fēng)洞模型的快速制造技術(shù)
    久久精品国产a三级三级三级| 三级经典国产精品| 欧美老熟妇乱子伦牲交| 亚洲中文av在线| 午夜福利在线观看免费完整高清在| 亚洲欧美成人精品一区二区| 国产精品一区二区在线观看99| 在线天堂最新版资源| 亚洲精品亚洲一区二区| 黄片无遮挡物在线观看| 少妇丰满av| 五月玫瑰六月丁香| videossex国产| 久久99热6这里只有精品| 高清毛片免费看| 丝袜喷水一区| 亚洲国产日韩一区二区| 天天躁夜夜躁狠狠久久av| 日韩不卡一区二区三区视频在线| 欧美另类一区| 日韩免费高清中文字幕av| 亚洲图色成人| 观看免费一级毛片| 免费少妇av软件| 80岁老熟妇乱子伦牲交| 午夜免费男女啪啪视频观看| 免费人成在线观看视频色| 国产精品伦人一区二区| 在线 av 中文字幕| 人体艺术视频欧美日本| 欧美成人午夜免费资源| 秋霞在线观看毛片| 免费观看无遮挡的男女| 观看美女的网站| 国产免费又黄又爽又色| 黑人巨大精品欧美一区二区蜜桃 | 99久久精品一区二区三区| 美女中出高潮动态图| 日韩欧美一区视频在线观看 | 中国国产av一级| 好男人视频免费观看在线| 国产精品国产三级国产专区5o| 新久久久久国产一级毛片| 免费人成在线观看视频色| 成人亚洲精品一区在线观看| 久久99蜜桃精品久久| 亚洲国产欧美在线一区| 国产精品久久久久久av不卡| 香蕉精品网在线| av专区在线播放| 99视频精品全部免费 在线| 少妇人妻一区二区三区视频| 亚洲av.av天堂| 久久影院123| 在现免费观看毛片| 久久精品国产a三级三级三级| 欧美人与善性xxx| 天天躁夜夜躁狠狠久久av| 亚洲av男天堂| 午夜91福利影院| 久久精品国产亚洲网站| 精品久久久噜噜| 日韩熟女老妇一区二区性免费视频| 亚洲av国产av综合av卡| xxx大片免费视频| 久久6这里有精品| 51国产日韩欧美| 丁香六月天网| 一级毛片黄色毛片免费观看视频| 国产成人91sexporn| 中文字幕av电影在线播放| 国产免费视频播放在线视频| 大香蕉97超碰在线| 国产淫片久久久久久久久| 青春草国产在线视频| 国产欧美日韩一区二区三区在线 | 综合色丁香网| 一级a做视频免费观看| 2018国产大陆天天弄谢| 日韩大片免费观看网站| 日日撸夜夜添| 国产成人精品婷婷| 亚州av有码| 国产成人91sexporn| 精品一区二区免费观看| 99久久中文字幕三级久久日本| 国产精品秋霞免费鲁丝片| 美女主播在线视频| 3wmmmm亚洲av在线观看| 黑人高潮一二区| 亚洲欧洲精品一区二区精品久久久 | 国产一区有黄有色的免费视频| 五月开心婷婷网| 国产成人精品久久久久久| 人人妻人人澡人人看| 97在线视频观看| 国产免费一区二区三区四区乱码| 午夜91福利影院| 校园人妻丝袜中文字幕| 精品久久国产蜜桃| 婷婷色麻豆天堂久久| 久久热精品热| freevideosex欧美| 久久99精品国语久久久| 2022亚洲国产成人精品| 亚洲第一区二区三区不卡| 国内少妇人妻偷人精品xxx网站| 久久免费观看电影| 国产极品天堂在线| 在线观看免费日韩欧美大片 | 尾随美女入室| 国产精品一区二区在线观看99| 欧美区成人在线视频| 亚洲真实伦在线观看| 国产一区二区在线观看av| 中文字幕av电影在线播放| 丝袜在线中文字幕| 女人久久www免费人成看片| 久久久午夜欧美精品| 国产欧美另类精品又又久久亚洲欧美| 激情五月婷婷亚洲| 日韩,欧美,国产一区二区三区| 国产色爽女视频免费观看| 日韩电影二区| 亚洲av综合色区一区| 国产日韩欧美视频二区| 久久精品国产a三级三级三级| 精品一区二区三卡| 亚洲国产日韩一区二区| 少妇高潮的动态图| 老司机影院毛片| 亚洲av成人精品一区久久| 99国产精品免费福利视频| h日本视频在线播放| 夫妻性生交免费视频一级片| 欧美日韩综合久久久久久| av福利片在线观看| 天天操日日干夜夜撸| 最新的欧美精品一区二区| 欧美日韩精品成人综合77777| 桃花免费在线播放| 成人无遮挡网站| 老熟女久久久| 天堂中文最新版在线下载| 精品久久久久久电影网| 国产精品熟女久久久久浪| 国产日韩欧美亚洲二区| 性高湖久久久久久久久免费观看| 美女脱内裤让男人舔精品视频| 精品国产一区二区三区久久久樱花| 亚洲经典国产精华液单| 大陆偷拍与自拍| 亚洲欧美精品自产自拍| 亚洲欧美一区二区三区黑人 | 下体分泌物呈黄色| 国产精品99久久99久久久不卡 | 又爽又黄a免费视频| 国产亚洲av片在线观看秒播厂| 一级毛片aaaaaa免费看小| 欧美97在线视频| 国产免费又黄又爽又色| 国产又色又爽无遮挡免| 国产欧美亚洲国产| 国产精品久久久久久久电影| 久久毛片免费看一区二区三区| 我的女老师完整版在线观看| 免费观看无遮挡的男女| 久久综合国产亚洲精品| 丰满乱子伦码专区| 久久久久人妻精品一区果冻| 我要看日韩黄色一级片| 18禁在线无遮挡免费观看视频| 我要看黄色一级片免费的| 久久久久精品性色| 亚洲精品第二区| 狠狠精品人妻久久久久久综合| 中文字幕亚洲精品专区| 丝瓜视频免费看黄片| 老女人水多毛片| av在线老鸭窝| 国产一区有黄有色的免费视频| 久久精品熟女亚洲av麻豆精品| 日韩精品免费视频一区二区三区 | 熟女电影av网| 女人久久www免费人成看片| 日韩一本色道免费dvd| av免费在线看不卡| 六月丁香七月| 久久精品国产亚洲av天美| 久久久久国产网址| 免费播放大片免费观看视频在线观看| 久久 成人 亚洲| 三级国产精品片| av视频免费观看在线观看| 久久99一区二区三区| 婷婷色综合大香蕉| 国产成人a∨麻豆精品| 美女视频免费永久观看网站| 黄色配什么色好看| 国产一区二区在线观看av| 亚洲,一卡二卡三卡| 天堂俺去俺来也www色官网| 少妇被粗大的猛进出69影院 | 日产精品乱码卡一卡2卡三| 亚洲国产av新网站| 成人无遮挡网站| 亚洲图色成人| 精品人妻熟女毛片av久久网站| 免费人妻精品一区二区三区视频| 国产男人的电影天堂91| 亚洲国产色片| 国产午夜精品一二区理论片| 肉色欧美久久久久久久蜜桃| 天天操日日干夜夜撸| 成年人午夜在线观看视频| 日本猛色少妇xxxxx猛交久久| 国产高清国产精品国产三级| 日日摸夜夜添夜夜爱| 午夜激情福利司机影院| 亚洲av成人精品一区久久| 国产在线男女| 蜜桃久久精品国产亚洲av| 99久久人妻综合| 十八禁高潮呻吟视频 | 亚洲精品亚洲一区二区| 边亲边吃奶的免费视频| 国产 一区精品| 一本色道久久久久久精品综合| 日韩三级伦理在线观看| 美女中出高潮动态图| 欧美成人精品欧美一级黄| 亚洲av在线观看美女高潮| 99久久精品一区二区三区| 人妻 亚洲 视频| 赤兔流量卡办理| 边亲边吃奶的免费视频| 久久人人爽人人爽人人片va| 99久久精品国产国产毛片| 免费少妇av软件| 伦精品一区二区三区| 国产真实伦视频高清在线观看| 看免费成人av毛片| 色网站视频免费| 男人舔奶头视频| 日本vs欧美在线观看视频 | 色吧在线观看| 观看美女的网站| 久久精品夜色国产| 男女无遮挡免费网站观看| 97在线人人人人妻| 九九爱精品视频在线观看| 欧美变态另类bdsm刘玥| 国产免费一级a男人的天堂| 亚洲国产av新网站| av在线播放精品| 最近2019中文字幕mv第一页| 乱人伦中国视频| 国产毛片在线视频| 天天躁夜夜躁狠狠久久av| 啦啦啦在线观看免费高清www| 热re99久久国产66热| 在线观看www视频免费| 成人午夜精彩视频在线观看| 免费观看在线日韩| 99久久精品热视频| 久久99热这里只频精品6学生| 亚洲av免费高清在线观看| 人妻少妇偷人精品九色| 午夜精品国产一区二区电影| 另类亚洲欧美激情| 色94色欧美一区二区| 亚洲激情五月婷婷啪啪| 一级毛片 在线播放| 天堂8中文在线网| 卡戴珊不雅视频在线播放| 中文乱码字字幕精品一区二区三区| 又爽又黄a免费视频| 国产精品99久久久久久久久| 91午夜精品亚洲一区二区三区| 亚洲一级一片aⅴ在线观看| 日韩强制内射视频| 99久久精品国产国产毛片| 国产精品久久久久久精品古装| 男女无遮挡免费网站观看| 欧美日韩一区二区视频在线观看视频在线| 国产女主播在线喷水免费视频网站| 中文字幕精品免费在线观看视频 | 能在线免费看毛片的网站| 麻豆乱淫一区二区| 国产高清三级在线| 在线播放无遮挡| 在线天堂最新版资源| 婷婷色综合大香蕉| 丰满迷人的少妇在线观看| 一本一本综合久久| 午夜精品国产一区二区电影| 亚洲丝袜综合中文字幕| 日本-黄色视频高清免费观看| 免费观看a级毛片全部| 伦理电影免费视频| 丰满少妇做爰视频| 国产精品久久久久久精品古装| 欧美高清成人免费视频www| 九草在线视频观看| 亚洲一区二区三区欧美精品| 免费播放大片免费观看视频在线观看| 久久影院123| 在线免费观看不下载黄p国产| 丝袜喷水一区| 国产男女超爽视频在线观看| 黄片无遮挡物在线观看| 日日爽夜夜爽网站| av在线app专区| 久久亚洲国产成人精品v| 欧美亚洲 丝袜 人妻 在线| 久久ye,这里只有精品| 国模一区二区三区四区视频| 男人添女人高潮全过程视频| 街头女战士在线观看网站| 日日爽夜夜爽网站| av在线app专区| 亚洲av欧美aⅴ国产| 国产精品女同一区二区软件| 久久久久久人妻| 亚洲电影在线观看av| 精品久久国产蜜桃| 亚洲欧洲日产国产| 久久97久久精品| 久久久a久久爽久久v久久| 亚洲精品自拍成人| av国产精品久久久久影院| 亚洲av在线观看美女高潮| 日韩欧美 国产精品| 国产黄片视频在线免费观看| 狠狠精品人妻久久久久久综合| 高清黄色对白视频在线免费看 | 99久久人妻综合| 日韩人妻高清精品专区| 国精品久久久久久国模美| 国产色爽女视频免费观看| 日韩欧美 国产精品| 日日爽夜夜爽网站| 在线观看一区二区三区激情| 国产精品国产三级专区第一集| 精品卡一卡二卡四卡免费| 日韩一区二区视频免费看| 亚洲熟女精品中文字幕| 99久久精品国产国产毛片| 国产又色又爽无遮挡免| 国产极品粉嫩免费观看在线 | 欧美日韩精品成人综合77777| 丰满少妇做爰视频| 日韩电影二区| 免费观看a级毛片全部| 国产亚洲精品久久久com| 熟女av电影| 又大又黄又爽视频免费| 丝袜脚勾引网站| 黄色配什么色好看| 日本猛色少妇xxxxx猛交久久| 亚洲人成网站在线播| 欧美精品人与动牲交sv欧美| 久久精品久久精品一区二区三区| 新久久久久国产一级毛片| 91午夜精品亚洲一区二区三区| 人人妻人人添人人爽欧美一区卜| 十八禁高潮呻吟视频 | 久久久久人妻精品一区果冻| 亚洲人与动物交配视频| 国产探花极品一区二区| 黄片无遮挡物在线观看| 欧美日韩国产mv在线观看视频| av视频免费观看在线观看| 又大又黄又爽视频免费| 亚洲精品日本国产第一区| 亚洲在久久综合| 久久久久久伊人网av| 99热这里只有是精品在线观看| 国产在线视频一区二区| 我的女老师完整版在线观看| 少妇人妻精品综合一区二区| 亚洲av二区三区四区| 亚洲av福利一区| 香蕉精品网在线| 在线观看av片永久免费下载| 视频中文字幕在线观看| 亚洲精品国产av蜜桃| 99视频精品全部免费 在线| a级毛片免费高清观看在线播放| 成人二区视频| 亚洲成人av在线免费| 男女国产视频网站| 国产精品久久久久久久久免| 欧美日韩一区二区视频在线观看视频在线| 美女xxoo啪啪120秒动态图| 久久精品国产亚洲av天美| 国产精品成人在线| 国内精品宾馆在线| 亚洲性久久影院| 久久人人爽人人爽人人片va| 纯流量卡能插随身wifi吗| 国产女主播在线喷水免费视频网站| 综合色丁香网| 97在线人人人人妻| 只有这里有精品99| 国产在视频线精品| 成人综合一区亚洲| 国产熟女午夜一区二区三区 | 成年人午夜在线观看视频| www.av在线官网国产| 在线观看三级黄色| 啦啦啦啦在线视频资源| 国产免费福利视频在线观看| 日本av手机在线免费观看| 精品久久久噜噜| 免费播放大片免费观看视频在线观看| 最近2019中文字幕mv第一页| 看免费成人av毛片| av一本久久久久| 亚洲精品一区蜜桃| 婷婷色麻豆天堂久久| 一区二区三区乱码不卡18| 秋霞伦理黄片| 亚洲av.av天堂| 最近的中文字幕免费完整| 老司机影院毛片| 一级二级三级毛片免费看| 熟妇人妻不卡中文字幕| 久久热精品热| 精品少妇内射三级| 99久久精品一区二区三区| 美女国产视频在线观看| 国产 精品1| 女性生殖器流出的白浆| 国产精品一区二区在线不卡| 国产一区二区三区av在线| 婷婷色av中文字幕| 三级经典国产精品| 天天操日日干夜夜撸| 精华霜和精华液先用哪个| 男女国产视频网站| 建设人人有责人人尽责人人享有的| 亚洲精品456在线播放app| 午夜免费男女啪啪视频观看| 亚洲精品视频女| 亚洲情色 制服丝袜| 99热6这里只有精品| 国产成人精品福利久久| 大片电影免费在线观看免费| 亚洲精品中文字幕在线视频 | 成人亚洲欧美一区二区av| 日韩大片免费观看网站| 最近手机中文字幕大全| 午夜影院在线不卡| 高清黄色对白视频在线免费看 | 中文字幕制服av| 日韩,欧美,国产一区二区三区| 免费播放大片免费观看视频在线观看| 一二三四中文在线观看免费高清| 性色av一级| 国产欧美亚洲国产| 日韩不卡一区二区三区视频在线| 欧美+日韩+精品| 人人妻人人澡人人爽人人夜夜| 亚洲av男天堂| 国产成人aa在线观看| 午夜福利,免费看| 十分钟在线观看高清视频www | 午夜精品国产一区二区电影| 亚洲国产精品一区二区三区在线| 噜噜噜噜噜久久久久久91| 精品亚洲成a人片在线观看| 亚洲av日韩在线播放| 国产成人免费观看mmmm| 欧美日本中文国产一区发布| 亚洲欧美精品专区久久| 性色avwww在线观看| 国产伦理片在线播放av一区| 黑人高潮一二区| 欧美 日韩 精品 国产| 国内少妇人妻偷人精品xxx网站| 亚洲高清免费不卡视频| 久久精品熟女亚洲av麻豆精品| 欧美日韩精品成人综合77777| 国产免费一区二区三区四区乱码| 尾随美女入室| 又粗又硬又长又爽又黄的视频| 亚洲精品国产av蜜桃| 日韩欧美 国产精品| 啦啦啦中文免费视频观看日本| 亚洲图色成人| 欧美精品人与动牲交sv欧美| 嫩草影院入口| 十八禁网站网址无遮挡 | a 毛片基地| 夜夜看夜夜爽夜夜摸| 青春草亚洲视频在线观看| 久热久热在线精品观看| 夫妻性生交免费视频一级片| 精品久久国产蜜桃| 美女视频免费永久观看网站| 国产熟女欧美一区二区| 亚洲色图综合在线观看| 在线观看免费视频网站a站| 最新中文字幕久久久久| 免费看av在线观看网站| 交换朋友夫妻互换小说| 精品国产国语对白av| 日韩电影二区| 婷婷色综合大香蕉| 在线 av 中文字幕| 在线精品无人区一区二区三| 亚州av有码| 特大巨黑吊av在线直播| 美女大奶头黄色视频| 国产成人a∨麻豆精品| 久久久精品免费免费高清| 99久久精品热视频| 五月天丁香电影| 国产探花极品一区二区| 国产国拍精品亚洲av在线观看| 如何舔出高潮| 两个人的视频大全免费| av卡一久久| 内射极品少妇av片p| 亚洲精品乱久久久久久| 国产精品久久久久久精品古装| 国产成人精品无人区| av国产精品久久久久影院| 久久婷婷青草| 搡女人真爽免费视频火全软件| 狂野欧美激情性bbbbbb| 大陆偷拍与自拍| 国产综合精华液| 精品久久久噜噜| 国产 精品1| 免费少妇av软件| 久久99一区二区三区| 中国国产av一级| 国产色婷婷99| 少妇被粗大猛烈的视频| 熟女人妻精品中文字幕| 看免费成人av毛片| a级毛片免费高清观看在线播放| 如何舔出高潮| 男女边摸边吃奶| 高清欧美精品videossex| 女性生殖器流出的白浆| av一本久久久久| 国产免费又黄又爽又色| 国产伦精品一区二区三区视频9| 特大巨黑吊av在线直播| 亚洲精品日韩av片在线观看| 国产精品一区www在线观看| 黄色一级大片看看| 99九九线精品视频在线观看视频| 国产亚洲av片在线观看秒播厂| 婷婷色综合大香蕉| 一级av片app| 日产精品乱码卡一卡2卡三| 简卡轻食公司| 亚洲欧美一区二区三区黑人 | 国产男人的电影天堂91| 亚洲av在线观看美女高潮| 欧美日韩精品成人综合77777| 国产视频内射| freevideosex欧美| 人妻少妇偷人精品九色| av网站免费在线观看视频| 一本久久精品| 插逼视频在线观看| 国产黄片视频在线免费观看| 国产色爽女视频免费观看| 精品国产乱码久久久久久小说| 99久国产av精品国产电影| h日本视频在线播放| 天美传媒精品一区二区| 少妇被粗大的猛进出69影院 | 国产片特级美女逼逼视频| 午夜福利在线观看免费完整高清在| 观看美女的网站| 男女无遮挡免费网站观看| 七月丁香在线播放| 久久久久久久久久久丰满| 欧美bdsm另类| 啦啦啦啦在线视频资源| 亚洲激情五月婷婷啪啪| 春色校园在线视频观看| 精品少妇黑人巨大在线播放| 青春草亚洲视频在线观看| 亚洲久久久国产精品| 观看美女的网站| 日本黄色日本黄色录像| 国产黄片视频在线免费观看| 人人妻人人爽人人添夜夜欢视频 | 精品国产国语对白av| 久久久国产精品麻豆| 中文字幕av电影在线播放| 久久女婷五月综合色啪小说| 精品国产国语对白av| 免费观看av网站的网址| 亚洲va在线va天堂va国产| 国产男女超爽视频在线观看| 精品少妇久久久久久888优播| 久久久久久久久久成人| 午夜免费鲁丝| 亚洲欧美日韩东京热| 热re99久久国产66热| 久久久久国产网址| 久久韩国三级中文字幕| 久久久久久久久久人人人人人人| 亚洲精品国产av成人精品| 99精国产麻豆久久婷婷| 毛片一级片免费看久久久久| 人人妻人人澡人人爽人人夜夜|