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

    基于深度學(xué)習(xí)的高超聲速飛行器再入預(yù)測校正容錯制導(dǎo)

    2020-05-19 12:39:10余躍王宏倫
    兵工學(xué)報 2020年4期
    關(guān)鍵詞:深度

    余躍, 王宏倫

    (1.北京航空航天大學(xué) 自動化科學(xué)與電氣工程學(xué)院, 北京 100191; 2.北京航天自動控制研究所, 北京 100854;3.北京航空航天大學(xué) 飛行器控制一體化技術(shù)重點實驗室, 北京 100191)

    0 引言

    高超聲速飛行器因其飛行速度快、機動性能好、生存能力強、突防能力高以及遠(yuǎn)程精確打擊等優(yōu)點,具有重要的戰(zhàn)略意義和軍事、民用價值[1]。再入制導(dǎo)技術(shù)通過設(shè)計制導(dǎo)律使得在滿足路徑約束的前提下平穩(wěn)到達(dá)末端能量管理段,作為高超聲速飛行器的核心和關(guān)鍵技術(shù)之一,因再入環(huán)境的復(fù)雜性和不確定性以及再入過程中受熱流密度、動壓、過載的嚴(yán)格約束,成為各國航空航天領(lǐng)域競相研究的重點和難點問題[2]。

    再入制導(dǎo)算法一般分為標(biāo)稱軌跡制導(dǎo)與預(yù)測校正制導(dǎo)兩大類。前者是跟蹤預(yù)先設(shè)計好的符合各類約束的標(biāo)稱軌跡完成制導(dǎo)任務(wù),具有計算量小、實時性高的優(yōu)點,但需預(yù)先存儲標(biāo)稱軌跡信息,且對環(huán)境變化敏感,難以保證較高的制導(dǎo)精度;后者是對落點航程進(jìn)行預(yù)測并實時校正控制量的制導(dǎo)算法,相比于標(biāo)稱軌跡制導(dǎo),預(yù)測校正制導(dǎo)具有更高的落點精度,且對各種再入初始誤差和飛行過程中的干擾偏差具有更強的魯棒性,因此更加受到研究學(xué)者的青睞。

    預(yù)測校正制導(dǎo)算法最早由Schultz等[3]提出,該算法利用預(yù)測方程改變升阻比和滾轉(zhuǎn)角同時控制橫程和航程,但是制導(dǎo)律缺少軌跡阻尼。Lu[4]針對低升阻比飛行器進(jìn)行預(yù)測校正制導(dǎo)研究,在傾側(cè)角幅值參數(shù)化這一核心要素中設(shè)計了較大的傾側(cè)角終端值來增加算法魯棒性,同時提出通過改變初始傾側(cè)角幅值來減輕過載的策略。Xue等[5]針對中高升阻比飛行器提出全約束預(yù)測校正制導(dǎo)算法,利用修正的準(zhǔn)平衡滑翔條件,將速度- 高度空間的路徑約束轉(zhuǎn)化為與速度有關(guān)的傾側(cè)角幅值上界,提高了算法處理約束的能力,且保持了算法的簡易性和魯棒性。文獻(xiàn)[6]提出一種基于自適應(yīng)神經(jīng)模糊系統(tǒng)(ANFIS)的再入預(yù)測校正制導(dǎo)算法,設(shè)計并訓(xùn)練ANFIS控制器,將終端落點偏差信息輸入控制器快速求取校正后的傾側(cè)角幅值,減少了制導(dǎo)指令的解算時間。文獻(xiàn)[7]提出一種基于傾側(cè)角參數(shù)化的離線彈道優(yōu)化與在線預(yù)測校正相結(jié)合的再入制導(dǎo)算法,提升了預(yù)測校正算法的實時性,但是該制導(dǎo)算法較為復(fù)雜。文獻(xiàn)[8]針對傳統(tǒng)預(yù)測校正制導(dǎo)算法中的終端航程與高度不匹配問題,提出一種基于傾側(cè)角剖面的、嚴(yán)格約束終端高度的預(yù)測校正制導(dǎo)算法。文獻(xiàn)[9]采用二次函數(shù)模型進(jìn)行傾側(cè)角參數(shù)化設(shè)計,并利用考慮哥氏加速度項改進(jìn)的擬平衡滑翔條件,將再入走廊轉(zhuǎn)化為傾側(cè)角幅值約束,提出一種改進(jìn)的預(yù)測校正制導(dǎo)算法。

    從已有研究看,預(yù)測校正制導(dǎo)算法呈現(xiàn)出良好的發(fā)展態(tài)勢,但仍存在算法實時性不高、魯棒性不強的問題。此外,高超聲速飛行器經(jīng)常處于低溫真空、高輻射的惡劣環(huán)境中,其執(zhí)行機構(gòu)很容易產(chǎn)生故障,使得飛行器性能惡化甚至系統(tǒng)失穩(wěn),往往會造成經(jīng)濟(jì)上的重大損失甚至災(zāi)難性后果。此時,需要進(jìn)行姿態(tài)容錯控制,以保證系統(tǒng)的穩(wěn)定性和安全性。在一些情況下,姿態(tài)環(huán)的容錯控制可以使飛行器從故障中恢復(fù)。但是,當(dāng)發(fā)生較為嚴(yán)重的故障如多個舵面卡死時,原先的攻角指令可能因為飛行器不能旋轉(zhuǎn)配平而無法有效跟蹤。針對這種情況,姿態(tài)環(huán)的容錯控制必須結(jié)合制導(dǎo)環(huán)的容錯制導(dǎo),即根據(jù)故障自適應(yīng)調(diào)整制導(dǎo)指令,以獲得更大限度的容錯能力。

    關(guān)于高超聲速飛行器容錯制導(dǎo)的研究目前還非常少。文獻(xiàn)[10]設(shè)計了一種基于待飛最優(yōu)路徑法的容錯制導(dǎo)算法,主要包含外環(huán)制導(dǎo)重構(gòu)和在線軌跡重塑兩部分,其中在線軌跡重塑的步驟包括離線軌跡庫生成、離線軌跡庫建模與編碼、在線軌跡重塑。文獻(xiàn)[11]根據(jù)飛行器當(dāng)前飛行狀態(tài)和終端約束條件,使用勒讓德偽譜法進(jìn)行在線軌道重構(gòu),生成滿足各種軌道約束的最優(yōu)返回軌跡,并實時反饋更新當(dāng)前軌道控制量迎角和傾斜角,達(dá)到實時最優(yōu)閉環(huán)制導(dǎo)的目的。文獻(xiàn)[12]采用相鄰可行軌跡存在定理設(shè)計了容錯制導(dǎo)律,以解決再入段執(zhí)行器發(fā)生故障的軌跡重構(gòu)問題。

    從現(xiàn)有文獻(xiàn)看,目前關(guān)于容錯制導(dǎo)的算法主要存在兩方面問題:1)編碼復(fù)雜,需離線生成和存儲大量軌跡,對機載計算機存儲空間要求較高,如基于待飛最優(yōu)路徑法的容錯制導(dǎo)算法;2)算法實時性差,如基于故障下氣動參數(shù)估計的偽譜法。分析可知,容錯制導(dǎo)的一個關(guān)鍵問題是確定故障后可配平的飛行包線區(qū)域和升力系數(shù)、阻力系數(shù)。一旦求解出這些數(shù)據(jù),就可以借助現(xiàn)有成熟的再入制導(dǎo)技術(shù)實現(xiàn)故障下的容錯制導(dǎo)。

    預(yù)測校正制導(dǎo)算法對初始誤差不敏感,抗擾性、魯棒性和自適應(yīng)能力都很強,很適合解決飛行器故障條件下的容錯制導(dǎo)問題。然而,高超聲速飛行器具有快時變和強不確定性特性,必須進(jìn)一步提升傳統(tǒng)預(yù)測校正制導(dǎo)算法的實時性和魯棒性,以適應(yīng)飛行器容錯制導(dǎo)的需求。此外,還必須考慮飛行器故障時的特點,有針對性地設(shè)計容錯制導(dǎo)方案。

    針對傳統(tǒng)預(yù)測校正算法的實時性問題,可以從縱向制導(dǎo)律中的預(yù)測落點過程入手。預(yù)測落點的實質(zhì)為根據(jù)當(dāng)前狀態(tài)量和控制量,通過積分預(yù)測終點的過程,其中包含大量積分運算,會耗費大量的計算時間和成本。本質(zhì)上,根據(jù)狀態(tài)量和控制量求取終點是一個非線性映射,可以借助神經(jīng)網(wǎng)絡(luò)強大的逼近能力來擬合。然而傳統(tǒng)的人工神經(jīng)網(wǎng)絡(luò)只是一個淺層結(jié)構(gòu),參數(shù)訓(xùn)練速度慢且容易出現(xiàn)過擬合[13]。近年來深度學(xué)習(xí)算法因其緩解了傳統(tǒng)訓(xùn)練算法的局部最小性而引起廣泛關(guān)注,使得設(shè)計和訓(xùn)練深層神經(jīng)網(wǎng)絡(luò)成為可能[14]。

    考慮到故障和參數(shù)攝動都會對氣動參數(shù)變化量產(chǎn)生影響,如果深度神經(jīng)網(wǎng)絡(luò)沒有考慮氣動參數(shù)變化量的作用,只用標(biāo)稱氣動參數(shù)來訓(xùn)練,則飛行器故障時的制導(dǎo)精度必然受到影響。為了提高算法的魯棒性和容錯性,可將升力系數(shù)、阻力系數(shù)變化量作為深度神經(jīng)網(wǎng)絡(luò)的輸入,根據(jù)實際氣動參數(shù)(故障和氣動參數(shù)攝動后的氣動參數(shù))預(yù)測落點,提高算法的制導(dǎo)精度。

    基于上述分析,本文針對故障條件下高超聲速飛行器的容錯制導(dǎo)問題,借助預(yù)測校正制導(dǎo)算法抗擾性強、精度高的優(yōu)勢,提出一種基于深度學(xué)習(xí)的預(yù)測校正容錯制導(dǎo)算法。該容錯制導(dǎo)算法的主要創(chuàng)新點為:1)給出了故障條件下的可行攻角剖面和升力系數(shù)、阻力系數(shù)的求解方法,為后續(xù)設(shè)計縱向制導(dǎo)律打下了基礎(chǔ);2)利用深度神經(jīng)網(wǎng)絡(luò)代替縱向制導(dǎo)律中的積分運算預(yù)測落點的過程,大大提高了算法的實時性,同時深度神經(jīng)網(wǎng)絡(luò)的輸入?yún)?shù)中包含升力系數(shù)、阻力系數(shù)變化量,使得預(yù)測落點時考慮了故障和參數(shù)攝動的影響,算法魯棒性和容錯性大大提高;3)構(gòu)建了擴張狀態(tài)觀測器(ESO)實時估計升力系數(shù)、阻力系數(shù)變化量,作為深度神經(jīng)網(wǎng)絡(luò)的輸入。

    1 再入制導(dǎo)問題

    1.1 3自由度運動學(xué)方程

    不考慮地球自轉(zhuǎn)的影響,建立高超聲速飛行器3自由度無量綱運動方程[2,5]為

    (1)

    (2)

    K=0.5R0S/m,S為飛行器參考面積,m為飛行器質(zhì)量,ρ為大氣密度,C′L和C′D分別為故障和參數(shù)攝動下的升力系數(shù)、阻力系數(shù),

    (3)

    CL和CD分別為標(biāo)稱升力系數(shù)和阻力系數(shù),ΔCL和ΔCD分別為升力系數(shù)和阻力系數(shù)的百分比變化量,主要表征由氣動參數(shù)攝動和執(zhí)行機構(gòu)故障帶來的氣動參數(shù)變化。

    需要指出的是,在容錯制導(dǎo)中,描述舵面故障的一貫做法是將其轉(zhuǎn)換為升力系數(shù)、阻力系數(shù)變化量。這是因為如果直接考慮不同的舵面故障,由于有無數(shù)個可能的位置和多個舵面組合,處理起來特別繁瑣,難以工程應(yīng)用。從另一個角度看,舵面卡死會改變升力、阻力特性,同時也需要剩余健康舵面合理配置以實現(xiàn)飛行器的旋轉(zhuǎn)配平,這將進(jìn)一步改變升力、阻力特性。因此,舵面故障的影響可以用升力系數(shù)、阻力系數(shù)變化量來表征[10]。

    1.2 再入過程約束

    (4)

    (1/r-v2)/r-Lcosσ=0.

    (5)

    根據(jù)實際情況給定一個較小的傾側(cè)角下邊界約束值σEQ,(5)式可以轉(zhuǎn)化為如(6)式的擬平衡滑翔約束:

    (1/r-v2)/r-LcosσEQ≤0.

    (6)

    1.3 再入終端約束

    高超聲速飛行器3自由度無量綱運動方程一般以時間為自變量,然而終端時間的不確定性會給彈道積分帶來不便。再入終端高度和速度是已知的,故一般引入類似能量的變量:

    (7)

    再入飛行終端約束主要包括高度約束、速度約束和經(jīng)緯度約束,用公式表示為

    (8)

    式中:ef表示終端能量;rf、vf、θf和φf分別為終端高度、終端速度、終端經(jīng)度和終端緯度。

    2 基于深度學(xué)習(xí)的預(yù)測校正容錯制導(dǎo)

    預(yù)測校正制導(dǎo)包含縱向制導(dǎo)和側(cè)向制導(dǎo),縱向制導(dǎo)確定傾側(cè)角的幅值,側(cè)向制導(dǎo)確定傾側(cè)角的符號。縱向制導(dǎo)主要包括預(yù)測環(huán)節(jié)和校正環(huán)節(jié)。利用當(dāng)前狀態(tài)量和控制量對飛行器運動方程積分預(yù)測落點,得到當(dāng)前位置到預(yù)測落點的待飛航程信息,此為預(yù)測環(huán)節(jié)。利用當(dāng)前位置到預(yù)測落點和終端落點的待飛航程偏差信息對傾側(cè)角幅值進(jìn)行調(diào)節(jié),直至落點偏差在預(yù)定范圍內(nèi),然后輸出校正后的傾側(cè)角幅值,此為校正環(huán)節(jié)。側(cè)向制導(dǎo)采用基于航向角誤差走廊的側(cè)向制導(dǎo)邏輯,得到傾側(cè)角的符號。求得傾側(cè)角的幅值和符號后,結(jié)合預(yù)先設(shè)置的攻角剖面以及當(dāng)前狀態(tài)量,對運動學(xué)方程進(jìn)行數(shù)值積分直至下一個制導(dǎo)周期。不斷重復(fù)上述預(yù)測校正過程,直至飛行器到達(dá)預(yù)定范圍。

    上述傳統(tǒng)預(yù)測校正制導(dǎo)算法在縱向制導(dǎo)中伴隨著大量對運動方程積分的過程,由此帶來算法實時性問題。同時,高超聲速飛行器的一個顯著特征是具有強不確定性,即氣動參數(shù)存在大范圍攝動,并且執(zhí)行機構(gòu)故障也會給氣動參數(shù)帶來明顯變化,此時若仍對標(biāo)稱狀態(tài)(氣動參數(shù)未攝動)下的運動方程積分求預(yù)測落點,則會導(dǎo)致預(yù)測落點精確度不高,并最終影響預(yù)測校正制導(dǎo)精度。顯然,傳統(tǒng)預(yù)測校正制導(dǎo)算法的魯棒性和容錯性有待進(jìn)一步提高。

    針對實時性問題,本文利用離線訓(xùn)練好的深度神經(jīng)網(wǎng)絡(luò)代替?zhèn)鹘y(tǒng)算法中積分運算求預(yù)測落點的過程,輸入當(dāng)前狀態(tài)量、控制量以及氣動參數(shù)變化量即可實時預(yù)測落點,從而有效縮短制導(dǎo)指令解算時間,提高算法的實時性。

    針對魯棒性和容錯性問題,本文充分考慮氣動參數(shù)攝動和執(zhí)行機構(gòu)故障帶來的影響,將升力系數(shù)和阻力系數(shù)百分比變化量作為深度神經(jīng)網(wǎng)絡(luò)的輸入,從而可以根據(jù)氣動參數(shù)變化后的運動方程模型預(yù)測落點,獲得更高的制導(dǎo)精度。為獲取當(dāng)前時刻的氣動參數(shù)值,基于運動方程模型構(gòu)造擴張狀態(tài)觀測器估計升力系數(shù)、阻力系數(shù)變化量。基于深度學(xué)習(xí)的預(yù)測校正容錯制導(dǎo)系統(tǒng)原理框圖如圖1所示。

    圖1 基于深度學(xué)習(xí)的預(yù)測校正容錯制導(dǎo)系統(tǒng)原理框圖Fig.1 Block diagram of deep learning-based predictor-corrector fault-tolerant guidance system

    2.1 縱向制導(dǎo)律設(shè)計

    2.1.1 攻角剖面

    再入走廊是指再入過程中滿足所有約束條件的飛行區(qū)域,這些約束主要指熱流密度、動壓、過載和擬平衡滑翔約束。在計算再入走廊之前,需要提前制定好飛行中的攻角方案,一般采用常用的標(biāo)準(zhǔn)攻角方案[17]:

    (9)

    式中:Ma為馬赫數(shù)。但是,當(dāng)執(zhí)行機構(gòu)出現(xiàn)故障后,并非在所有攻角情況下飛行器都可以配平。為使飛行器可控,首先需要求取可配平的攻角剖面。

    故障條件下,利用控制分配算法配平飛行器。旋轉(zhuǎn)配平飛行器的必要條件是:所有控制面偏角產(chǎn)生的力矩系數(shù)與翼身組合體產(chǎn)生的力矩系數(shù)大小相等、方向相反,用公式表示為

    (10)

    式中:Cl,δδ、Cm,δδ、Cn,δδ分別為控制面偏角產(chǎn)生的滾轉(zhuǎn)力矩系數(shù)、俯仰力矩系數(shù)和偏航力矩系數(shù),δ為控制面偏角向量;Cm,α為翼身組合體俯仰力矩系數(shù)。

    控制分配器的首要任務(wù)是確定最優(yōu)控制面偏角使得控制不足最小,即(1)式中等式兩端的差最小,用公式表示為

    (11)

    (12)

    式中:W為權(quán)重向量。

    利用求得的舵偏角δ*,在攻角網(wǎng)格上生成翼身組合體和舵偏角產(chǎn)生的升力系數(shù)、阻力系數(shù),這些系數(shù)合起來就是總的升力系數(shù)、阻力系數(shù),表示為

    (13)

    式中:CL,α和CD,α分別為翼身組合體產(chǎn)生的升力系數(shù)和阻力系數(shù);CL,δ*δ*和CD,δ*δ*分別為舵偏角產(chǎn)生的升力系數(shù)和阻力系數(shù)。

    2.1.2 傾側(cè)角約束

    利用(11)式求得的攻角可行區(qū)域設(shè)計攻角方案,并將攻角方案代入再入過程約束(4)式和(6)式中,可得到滿足過程約束的高度- 速度剖面,即再入走廊。然而,在再入飛行過程中,如果在每一點上計算高度和速度的關(guān)系來驗證是否處于再入走廊,會導(dǎo)致巨大的計算量,通常利用擬平衡滑翔條件(5)式將再入走廊的約束轉(zhuǎn)換成傾側(cè)角幅值約束,使得飛行器在傾側(cè)角幅值約束范圍內(nèi)飛行就可以滿足再入過程約束。傾側(cè)角幅值限制表示為

    (14)

    (15)

    結(jié)合給定的傾側(cè)角下邊界約束值,可得傾側(cè)角約束為

    σEQ≤|σ|≤|σ|max.

    (16)

    2.1.3 傾側(cè)角幅值求解

    根據(jù)飛行特點,再入過程分為初始下降段和擬平衡滑翔段。初始下降段具有較高高度,氣動力較小,故采用常值傾側(cè)角|σ0|飛行。擬平衡滑翔段采用基于深度學(xué)習(xí)的預(yù)測校正容錯制導(dǎo)。具體說來,在每個制導(dǎo)周期內(nèi),給定初始傾側(cè)角的幅值|σi|,結(jié)合預(yù)設(shè)攻角剖面和當(dāng)前狀態(tài),由基于深度學(xué)習(xí)的神經(jīng)網(wǎng)絡(luò)算出當(dāng)前位置到預(yù)測落點的待飛航程,并結(jié)合終端落點計算待飛航程偏差。為使偏差為0 rad,采用常用的割線法求出傾側(cè)角幅值,并結(jié)合側(cè)向制導(dǎo)求解的傾側(cè)角符號,輸出制導(dǎo)指令。與傳統(tǒng)算法中對運動方程(1)式進(jìn)行積分來預(yù)測落點相比,本文算法可以避免大量的積分運算過程,從而大大減少在線計算時間,有力地提升預(yù)測校正制導(dǎo)算法的實時性。以下詳細(xì)介紹利用深層神經(jīng)網(wǎng)絡(luò)計算當(dāng)前位置到預(yù)測落點待飛航程的過程。

    2.1.3.1 深度神經(jīng)網(wǎng)絡(luò)

    深度學(xué)習(xí)的概念最早由多倫多大學(xué)Hinton等[18]于2006年提出,是指基于樣本數(shù)據(jù)通過一定的訓(xùn)練方法得到包含多個層級的深度網(wǎng)絡(luò)結(jié)構(gòu)的機器學(xué)習(xí)過程。傳統(tǒng)神經(jīng)網(wǎng)絡(luò)隨機初始化網(wǎng)絡(luò)中的權(quán)值,導(dǎo)致網(wǎng)絡(luò)很容易收斂到局部最小值,為解決這一問題,Hinton等[18]提出先使用無監(jiān)督預(yù)訓(xùn)練方法優(yōu)化網(wǎng)絡(luò)權(quán)值的初值、再進(jìn)行權(quán)值微調(diào)的方法,拉開了深度學(xué)習(xí)的序幕。

    深度學(xué)習(xí)所得到的深度網(wǎng)絡(luò)結(jié)構(gòu)包含大量單一元素(神經(jīng)元),每個神經(jīng)元與大量其他神經(jīng)元相連接,神經(jīng)元間的連接強度(權(quán)值)在學(xué)習(xí)過程中修改,并決定網(wǎng)絡(luò)的功能。通過深度學(xué)習(xí)得到的深度網(wǎng)絡(luò)結(jié)構(gòu)符合神經(jīng)網(wǎng)絡(luò)的特征,因此深度網(wǎng)絡(luò)就是深層次的神經(jīng)網(wǎng)絡(luò),即深度神經(jīng)網(wǎng)絡(luò)(DNN)[19]。

    2.1.3.2 預(yù)測環(huán)節(jié)的輸入輸出參數(shù)

    傳統(tǒng)預(yù)測校正算法預(yù)測環(huán)節(jié)為:根據(jù)當(dāng)前狀態(tài)變量和控制量,對(1)式進(jìn)行積分,得到預(yù)測落點。由此可得深度神經(jīng)網(wǎng)絡(luò)預(yù)測環(huán)節(jié)的輸入?yún)?shù)為地心距r、經(jīng)度θ、緯度φ、速度v、航跡傾角γ、航向角ψ、攻角α、傾側(cè)角σ、升力系數(shù)變化量ΔCL和阻力系數(shù)變化量ΔCD;輸出參數(shù)為當(dāng)前位置到預(yù)測落點的待飛航程sp.

    2.1.3.3 層數(shù)和節(jié)點數(shù)

    利用全連接神經(jīng)網(wǎng)絡(luò)來近似預(yù)測環(huán)節(jié)。關(guān)于隱含層層數(shù)和節(jié)點數(shù)的選取目前尚沒有相關(guān)理論支撐,隨著隱含層層數(shù)和節(jié)點數(shù)的增加,近似精度會提高,但同時也會增加計算量,因此在選擇層數(shù)和節(jié)點數(shù)時需要綜合考慮。本文選取5個隱含層,每個隱含層節(jié)點數(shù)為20,網(wǎng)絡(luò)結(jié)構(gòu)如圖2所示。

    圖2 深度神經(jīng)網(wǎng)絡(luò)結(jié)構(gòu)圖Fig.2 Structure diagram of deep neural network

    定義深度神經(jīng)網(wǎng)絡(luò)每層的變量為fk,1≤k≤7,輸入層變量f1=[r,θ,φ,v,γ,ψ,α,σ,ΔCL,ΔCD]T,輸出層變量f7=sp,選擇指數(shù)線性單元(ELU)為激活函數(shù):

    (17)

    則前向傳播計算過程為

    fk=fe(fk-1Wk-1+Bk-1),2≤k≤7,

    (18)

    式中:Wk-1表示fk-1和fk之間的權(quán)重矩陣;Bk-1表示fk的偏置向量。

    2.1.3.4 深度神經(jīng)網(wǎng)絡(luò)的訓(xùn)練

    深度神經(jīng)網(wǎng)絡(luò)的訓(xùn)練步驟為:

    1)樣本獲取。以方程(1)式中的狀態(tài)變量和攻角、傾側(cè)角以及升力、阻力系數(shù)變化量為輸入量,對微分方程(1)式進(jìn)行積分,得到預(yù)測落點,并計算當(dāng)前位置到預(yù)測落點的待飛航程。當(dāng)狀態(tài)變量、攻角和傾側(cè)角以及升力、阻力系數(shù)變化量按照表1不斷變化時,可以得到大量的輸入- 輸出數(shù)據(jù)對,用于深度神經(jīng)網(wǎng)絡(luò)的訓(xùn)練。

    表1 深度神經(jīng)網(wǎng)絡(luò)輸入量變化情況Tab.1 Input variables of deep neural network

    2)網(wǎng)絡(luò)訓(xùn)練。將1 865 676組輸入- 輸出數(shù)據(jù)對用深度神經(jīng)網(wǎng)絡(luò)擬合,圖3為其中一段100個樣本點深度神經(jīng)網(wǎng)絡(luò)的擬合值和通過積分運算得到的真實值的對比曲線。圖4和圖5分別為深度神經(jīng)網(wǎng)絡(luò)擬合的絕對誤差和相對誤差曲線,可以看出擬合相對誤差在10%以內(nèi)。表2為深度神經(jīng)網(wǎng)絡(luò)擬合誤差統(tǒng)計結(jié)果。從表2可知,擬合誤差均值為-1.541 9×10-4rad,標(biāo)準(zhǔn)差為0.036 3 rad. 由此可以看出,所訓(xùn)練的深度神經(jīng)網(wǎng)絡(luò)對落點預(yù)測環(huán)節(jié)具有很高的近似精度,具有表征該非線性映射的能力。

    圖3 深度神經(jīng)網(wǎng)絡(luò)擬合曲線Fig.3 Fitting curve of deep neural network

    圖4 深度神經(jīng)網(wǎng)絡(luò)擬合絕對誤差曲線Fig.4 Absolute fitting error curve of deep neural network

    圖5 深度神經(jīng)網(wǎng)絡(luò)擬合相對誤差曲線Fig.5 Relative fitting error curve of deep neural network

    2.1.3.5 割線法求幅值

    利用訓(xùn)練的深度神經(jīng)網(wǎng)絡(luò)求出待飛航程sp后,計算待飛航程偏差:

    表2 深度神經(jīng)網(wǎng)絡(luò)擬合誤差統(tǒng)計Tab.2 Fitting error statistics of deep neural network

    fi(|σi|)=sp-sf,

    (19)

    式中:sf為當(dāng)前位置到終端落點的待飛航程,

    sf=arccos (sinφsinφf+cosφcosφfcos (θf-θ)).

    (20)

    為使待飛航程偏差為0 rad,需求得方程fi(|σi|)=0 rad的零解,利用割線法可計算傾側(cè)角幅值為

    (21)

    2.2 側(cè)向制導(dǎo)律設(shè)計

    設(shè)計側(cè)向制導(dǎo)律主要是為了求取傾側(cè)角的符號。在初始階段,傾側(cè)角選取與再入點航向誤差相反的符號。在隨后的再入飛行過程中,為了滿足飛行過程中各種約束條件和預(yù)期的落點精度要求,需要在特定的時候使傾側(cè)角的符號反向。然而,過多的反轉(zhuǎn)次數(shù)會破壞飛行的自適應(yīng)性和魯棒性[10]。因此,一個好的反轉(zhuǎn)邏輯是使得末端航向角誤差在設(shè)定范圍內(nèi),同時具有較少的反轉(zhuǎn)次數(shù)。

    定義Φ為當(dāng)前位置到目標(biāo)點的視線方位角,其計算式為

    (22)

    由此可得航向角誤差為Δψ=ψ-Φ. 為使落點達(dá)到期望的精度范圍,航向角誤差需滿足:

    -Δψmax≤Δψ≤Δψmax,

    (23)

    式中:Δψmax為航向角誤差上界。航向角誤差走廊如圖6所示。

    圖6 航向角誤差走廊Fig.6 Heading angle error corridor

    傾側(cè)角反轉(zhuǎn)邏輯為:當(dāng)航向角誤差超出誤差上界Δψmax時,傾側(cè)角符號為負(fù);當(dāng)航向角誤差滿足誤差范圍(18)式時,傾側(cè)角符號保持不變;當(dāng)航向角誤差超出誤差下界-Δψmax時,傾側(cè)角符號為正。上述反轉(zhuǎn)邏輯用公式表示為

    (24)

    在實際飛行中,本節(jié)提出的深度神經(jīng)網(wǎng)絡(luò)是有一定適應(yīng)邊界的,這主要是因為訓(xùn)練深度神經(jīng)網(wǎng)絡(luò)的數(shù)據(jù)集存在一定的邊界。雖然深度神經(jīng)網(wǎng)絡(luò)具有一定的泛化能力,但是為了確保飛行安全,飛行的狀態(tài)量、控制量以及氣動參數(shù)變化量需要在該邊界范圍內(nèi)。具體地,高度范圍為[20 km, 80 km],速度范圍為[1 800 m/s,7 000 m/s]s,經(jīng)度范圍為[10°, 90 °],緯度范圍為[-20°, 30°],航跡傾角范圍為[-5°, 5°],航跡偏角范圍為[45°, 70°]°,攻角范圍為[10 °, 45°],傾側(cè)角范圍為[-70°, 70°],升力系數(shù)變化量為[-50%, 50%],阻力系數(shù)變化量為[-50%, 50%],當(dāng)飛行器的變量滿足上述范圍時,可采用本節(jié)所提基于深度學(xué)習(xí)的預(yù)測校正容錯制導(dǎo)算法。在一些情況下,飛行器的變量在上述范圍之外,可補充更大范圍的數(shù)據(jù)集對深度神經(jīng)網(wǎng)絡(luò)進(jìn)一步訓(xùn)練,以擴大本節(jié)容錯制導(dǎo)算法的適用范圍。

    關(guān)于基于深度學(xué)習(xí)的預(yù)測校正容錯制導(dǎo)算法的收斂性,現(xiàn)從兩方面分析:一是傳統(tǒng)預(yù)測校正制導(dǎo)算法的收斂性,如果利用傳統(tǒng)的預(yù)測校正算法沒法求出滿足再入終端約束的傾側(cè)角和攻角指令,則說明終端落點是不可達(dá)的,算法無法收斂,此時采用基于深度學(xué)習(xí)的預(yù)測校正容錯制導(dǎo)算法也無法收斂;二是深度神經(jīng)網(wǎng)絡(luò)近似非線性映射關(guān)系的收斂性,理論上,對于任意一個3層或者3層以上神經(jīng)網(wǎng)絡(luò),只要隱含層神經(jīng)元數(shù)目足夠多,該網(wǎng)絡(luò)就能以任意精度逼近一個非線性映射關(guān)系[20]。本文采用深度神經(jīng)網(wǎng)絡(luò),包含5個隱含層,每個隱含層節(jié)點數(shù)為20,能以很高精度(相對誤差小于10%)來近似“積分運算求預(yù)測落點”這個非線性映射關(guān)系。因此,只要終端落點是可達(dá)的,本文所提基于深度學(xué)習(xí)的預(yù)測校正容錯制導(dǎo)算法極大概率就是收斂的。

    3 氣動參數(shù)變化量估計

    為了在存在氣動參數(shù)攝動條件下準(zhǔn)確地預(yù)測落點,2.1.3節(jié)已經(jīng)構(gòu)造了輸入包含氣動參數(shù)變化量的深度神經(jīng)網(wǎng)絡(luò)。因此,要想預(yù)測落點,首先要考慮氣動參數(shù)變化量的獲取問題。借鑒自抗擾控制中利用ESO估計干擾的思想,本文在制導(dǎo)環(huán)構(gòu)建ESO對氣動參數(shù)變化量進(jìn)行估計,并將估計值輸入深度神經(jīng)網(wǎng)絡(luò)預(yù)測落點,輸出當(dāng)前位置到預(yù)測落點的待飛航程。

    結(jié)合飛行器3自由度運動方程(1)式,參考文獻(xiàn)[21-22]中ESO的構(gòu)造方法,構(gòu)造ESO估計阻力相關(guān)項:

    (25)

    式中:ev為速度的估計誤差;z1v為速度v的估計;z2v為阻力相關(guān)項-D的估計;β1v和β2v為觀測器增益。結(jié)合(3)式可得,阻力系數(shù)百分比變化量估計為

    (26)

    類似地,構(gòu)造ESO估計升力相關(guān)項:

    (27)

    式中:eγ為航跡傾角的估計誤差;z1γ為航跡傾角γ的估計;z2γ為升力相關(guān)項Lcosσ/v的估計;β1γ和β2γ為觀測器增益。結(jié)合(2)式可得升力系數(shù)百分比變化量估計為

    (28)

    可以證明,當(dāng)ESO增益取合適值時,升力系數(shù)百分比變化量ΔCL和阻力系數(shù)百分比變化量ΔCD能被ESO快速精確估計。

    4 仿真分析

    考慮文獻(xiàn)[2]中的再入初始狀態(tài),高度H0=80 km,速度v0=7 100 m/s,經(jīng)度θ0=10°,緯度φ0=-20°,航跡傾角γ0=-1°,航向角ψ0=45°,其中下標(biāo)0表示再入初始值。執(zhí)行機構(gòu)故障類型為:右內(nèi)側(cè)升降舵δ1、右外側(cè)升降舵δ2、左內(nèi)側(cè)升降舵δ3和左外側(cè)升降舵δ4均卡死在-20°. 終端約束為:經(jīng)度θf=90°,緯度φf=30°,高度Hf=20 km,速度vf=1 800 m/s. 要求終端位置誤差小于10 km,速度誤差小于100 m/s,高度誤差小于2 km,熱流密度約束為1.5 MW/m2,動壓約束為200 kPa,過載約束為4.5g. 為驗證本文所提容錯制導(dǎo)算法的可行性和優(yōu)越性,將文獻(xiàn)[5]中傳統(tǒng)的預(yù)測校正制導(dǎo)算法和基于深度學(xué)習(xí)的預(yù)測校正容錯制導(dǎo)算法仿真結(jié)果進(jìn)行對比。

    4.1 可行攻角剖面分析

    當(dāng)X-33飛行器的4個升降舵都卡死在-20°時,可供選擇的攻角范圍受到限制,如果攻角剖面選擇不合理,則飛行器會處于失控狀態(tài)。根據(jù)(11)式中的優(yōu)化原則,得到仿真樣例故障的配平不足曲線,如圖7所示。由圖7可知:當(dāng)攻角α∈[-10°,20°]∪[38.8°,50°]時,配平不足值JD=0,表明飛行器三軸力矩處于平衡狀態(tài),飛行器可旋轉(zhuǎn)配平;當(dāng)α∈(20°,38.8°)時,配平不足值JD>0,表明飛行器三軸力矩處于非平衡狀態(tài),飛行器不可旋轉(zhuǎn)配平。考慮故障條件下飛行器的可控性和安全性需求,本文選取故障條件下的攻角剖面為α=12°.

    圖7 故障下的配平不足曲線Fig.7 Trim deficiency curve in case of fault

    4.2 氣動參數(shù)變化量估計分析

    由(13)式分別求出標(biāo)稱條件下的升力系數(shù)、阻力系數(shù)CL和CD以及故障和氣動參數(shù)攝動下的升力系數(shù)、阻力系數(shù)C′L和C′D,并由(3)式求得升力系數(shù)、阻力系數(shù)變化量ΔCL和ΔCD,其曲線如圖8和圖9所示。由圖8和圖9可知,升力系數(shù)、阻力系數(shù)變化量大體處于合理范圍,但是升力系數(shù)變化量在零攻角處存在尖峰,這是因為標(biāo)稱升力系數(shù)在零攻角處值太小導(dǎo)致。在4.1節(jié)中,已經(jīng)選取攻角剖面為α=12°,則由圖8和圖9可得由于故障導(dǎo)致的升力、阻力系數(shù)百分比變化量分別為1.82%和16.91%.

    圖8 升力系數(shù)變化量Fig.8 Variation of lift coefficient

    圖9 阻力系數(shù)變化量Fig.9 Variation of drag coefficient

    除了故障會導(dǎo)致氣動參數(shù)變化之外,復(fù)雜的大氣環(huán)境、傳感器誤差以及建模誤差也會導(dǎo)致氣動參數(shù)存在攝動。為了更好地展示ESO對氣動參數(shù)變化量的估計,考慮氣動參數(shù)攝動在[-30%,30%]之內(nèi)變化,并用正弦函數(shù)30sin(1.5t)表征。在故障和氣動參數(shù)攝動的綜合作用下,總的升力系數(shù)、阻力系數(shù)變化量估計曲線如圖10和圖11所示。由圖10和圖11可見,ESO可對氣動參數(shù)變化量實時精準(zhǔn)估計,為將氣動參數(shù)變化量輸入深度神經(jīng)網(wǎng)絡(luò)、實現(xiàn)故障和氣動參數(shù)攝動下的飛行器落點精準(zhǔn)預(yù)測提供有力保障。

    圖10 升力系數(shù)變化量估計Fig.10 Estimation of lift coefficient variation

    圖11 阻力系數(shù)變化量估計Fig.11 Estimation of drag coefficient variation

    4.3 標(biāo)準(zhǔn)條件下制導(dǎo)仿真分析

    在前述標(biāo)準(zhǔn)初始狀態(tài)且不考慮過程擾動的情況下,分別采用傳統(tǒng)預(yù)測校正制導(dǎo)算法和基于深度學(xué)習(xí)的預(yù)測校正容錯制導(dǎo)算法進(jìn)行制導(dǎo)仿真。終端狀態(tài)誤差和算法運行時間如表3所示。

    從表3中可以看出,兩種預(yù)測校正制導(dǎo)算法均能達(dá)到終端誤差要求且誤差數(shù)值較小,在終端速度約束方面,兩種算法效果相當(dāng),在終端位置和高度約束方面,深度學(xué)習(xí)容錯制導(dǎo)算法優(yōu)于傳統(tǒng)制導(dǎo)算法。值得注意的是,深度學(xué)習(xí)容錯制導(dǎo)算法采用深度神經(jīng)網(wǎng)絡(luò)提前建立了狀態(tài)量、控制量以及氣動參數(shù)變化量與預(yù)測落點的非線性映射關(guān)系,避免了傳統(tǒng)制導(dǎo)算法中大量的積分運算過程,從而使得制導(dǎo)指令解算時間大大減少。傳統(tǒng)制導(dǎo)算法的制導(dǎo)指令解算時間為195 649 ms,而深度學(xué)習(xí)容錯制導(dǎo)算法為9 587 ms;具體到每個制導(dǎo)周期,傳統(tǒng)制導(dǎo)算法的制導(dǎo)指令平均解算時間為103.46 ms,而深度學(xué)習(xí)容錯制導(dǎo)算法為5.02 ms. 因此深度學(xué)習(xí)容錯制導(dǎo)算法實時性遠(yuǎn)遠(yuǎn)高于傳統(tǒng)制導(dǎo)算法。

    表3 終端狀態(tài)誤差和算法運行時間Tab.3 Terminal state error and operation time

    圖12為兩種預(yù)測校正制導(dǎo)算法下的三維軌跡曲線圖,可見兩條軌跡基本重合,均能從再入初始點跳躍滑翔至預(yù)定終端落點。圖13為兩種算法的傾側(cè)角- 航程距離曲線圖,其中s表示航程距離。從圖13中可以看出,傾側(cè)角在約束范圍內(nèi)變化且反轉(zhuǎn)次數(shù)較少。圖14~圖16分別為熱流密度、動壓和過載歷程曲線,橫虛線表示約束邊界,從中可以看出,再入飛行時過程約束均得到滿足。

    圖12 兩種預(yù)測校正制導(dǎo)算法的三維軌跡曲線Fig.12 Three-dimensional trajectories of two predictor-corrector guidance algorithms

    圖13 傾側(cè)角- 航程距離曲線Fig.13 Curves of bank angle-downrange

    圖14 熱流密度- 航程距離曲線Fig.14 Curves of heat flux density-downrange

    圖15 動壓- 航程距離曲線Fig.15 Curves of dynamic pressure-downrange

    圖16 過載- 航程距離曲線Fig.16 Curves of overload-downrange

    4.4 故障和干擾條件下制導(dǎo)仿真分析

    為驗證基于深度學(xué)習(xí)的預(yù)測校正容錯制導(dǎo)算法的容錯性和魯棒性,對存在初始狀態(tài)偏差、執(zhí)行機構(gòu)卡死故障以及氣動參數(shù)攝動的情況進(jìn)行仿真分析。偏差和攝動如表4所示。

    表4 初始狀態(tài)偏差和參數(shù)攝動Tab.4 Initial state error and parameter perturbation

    需要說明的是,表4中的升力系數(shù)、阻力系數(shù)變化量ΔCL和ΔCD僅僅是指由氣動參數(shù)攝動引起的變化,變量在區(qū)間[-30%,30%]內(nèi)服從高斯分布。此外,還有由故障引起的升力系數(shù)、阻力系數(shù)變化,分別為1.82%和16.91%,需要在仿真時一并考慮。

    在擾動條件下采用基于深度學(xué)習(xí)的預(yù)測校正容錯制導(dǎo)算法進(jìn)行100組Monte Carlo仿真,圖17~圖23給出了制導(dǎo)仿真曲線。由圖17~圖23可以看出,傾側(cè)角在約束范圍內(nèi)經(jīng)過5次左右反轉(zhuǎn)后,飛行器能在終端約束范圍內(nèi)達(dá)到目標(biāo)點,且再入過程中熱流密度、動壓和過載都滿足約束條件。

    圖17 三維軌跡圖Fig.17 Three-dimensional trajectory curve

    圖18 二維航跡圖Fig.18 Two-dimensional trajectory curve

    圖19 速度剖面圖Fig.19 Velocity curve

    圖20 傾側(cè)角剖面圖Fig.20 Bank angle curve

    圖21 熱流密度剖面圖Fig.21 Heat flux density curve

    圖22 動壓剖面圖Fig.22 Dynamic pressure curve

    圖23 過載剖面圖Fig.23 Overload curve

    圖24給出了故障和干擾條件下兩種制導(dǎo)算法100組Monte Carlo仿真的落點經(jīng)緯度散布圖。由圖24可見:采用基于深度學(xué)習(xí)的預(yù)測校正容錯制導(dǎo)算法時,飛行器落點位置誤差均小于5 km;采用傳統(tǒng)預(yù)測校正制導(dǎo)算法時,有相當(dāng)一部分落點位于目標(biāo)落點10 km之外。由此可見,在故障和干擾條件下,基于深度學(xué)習(xí)的預(yù)測校正容錯制導(dǎo)算法具有更高的制導(dǎo)精度,也證明了該算法具有更好的容錯性和魯棒性。

    圖24 故障和干擾條件下落點經(jīng)度和緯度散布圖Fig.24 Longitude and latitude dispersion under the conditions of fault and disturbance

    圖25 制導(dǎo)指令解算時間Fig.25 Calculating time of guidance command

    圖25給出了兩種制導(dǎo)算法制導(dǎo)指令解算時間的100組仿真數(shù)據(jù)。由圖25可以看出,基于深度學(xué)習(xí)的預(yù)測校正容錯制導(dǎo)算法在制導(dǎo)指令解算上遠(yuǎn)遠(yuǎn)快于傳統(tǒng)預(yù)測校正制導(dǎo)算法。表5進(jìn)一步給出了制導(dǎo)指令解算時間的統(tǒng)計數(shù)據(jù)。由表5可見:傳統(tǒng)預(yù)測校正制導(dǎo)算法仿真時間均值為199.32 s,標(biāo)準(zhǔn)差為4.44 s;基于深度學(xué)習(xí)的預(yù)測校正容錯制導(dǎo)算法的平均仿真時間僅為9.31 s,標(biāo)準(zhǔn)差為0.42 s. 以上結(jié)果表明本文所提算法能使傳統(tǒng)算法的制導(dǎo)指令解算時間減少95%,大大提高了傳統(tǒng)算法的實時性。這是因為提前離線生成深度神經(jīng)網(wǎng)絡(luò)來預(yù)測落點,有效避免了傳統(tǒng)算法中通過大量積分運算預(yù)測落點的過程。

    表5 制導(dǎo)指令解算時間統(tǒng)計數(shù)據(jù)Tab.5 Calculating time of guidance command s

    5 結(jié)論

    本文針對故障條件下高超聲速飛行器的容錯制導(dǎo)問題,借助預(yù)測校正制導(dǎo)算法抗擾性強、精度高的優(yōu)勢,提出一種基于深度學(xué)習(xí)的預(yù)測校正容錯制導(dǎo)算法。仿真結(jié)果表明:

    1)故障條件下通過控制分配算法能有效求解可行攻角剖面和升力系數(shù)、阻力系數(shù)。

    2)ESO可對氣動參數(shù)變化量實時精準(zhǔn)估計,為深度神經(jīng)網(wǎng)絡(luò)精確預(yù)測落點提供了基礎(chǔ)。

    3)在標(biāo)準(zhǔn)條件下,傳統(tǒng)預(yù)測校正制導(dǎo)算法與基于深度學(xué)習(xí)的預(yù)測校正容錯制導(dǎo)算法的制導(dǎo)精度相當(dāng)。

    4)在故障和干擾條件下,基于深度學(xué)習(xí)的預(yù)測校正容錯制導(dǎo)算法的制導(dǎo)精度大大優(yōu)于傳統(tǒng)預(yù)測校正制導(dǎo)算法,前者飛行器落點位置誤差均小于5 km,后者飛行器有相當(dāng)一部分落點位于目標(biāo)落點10 km之外。

    5)基于深度學(xué)習(xí)的預(yù)測校正容錯制導(dǎo)算法的實時性遠(yuǎn)遠(yuǎn)優(yōu)于傳統(tǒng)預(yù)測校正制導(dǎo)算法,能使傳統(tǒng)算法的制導(dǎo)指令解算時間減少95%.

    猜你喜歡
    深度
    深度理解不等關(guān)系
    四增四減 深度推進(jìn)
    深度理解一元一次方程
    深度觀察
    深度觀察
    深度觀察
    深度觀察
    芻議深度報道的深度與“文”度
    新聞傳播(2016年10期)2016-09-26 12:14:59
    提升深度報道量與質(zhì)
    新聞傳播(2015年10期)2015-07-18 11:05:40
    微小提議 深度思考
    狂野欧美激情性bbbbbb| 久久久久国产精品人妻一区二区| av卡一久久| 国产成人91sexporn| 午夜日本视频在线| 伊人亚洲综合成人网| 伊人亚洲综合成人网| 亚洲欧美一区二区三区黑人 | 亚洲欧美日韩另类电影网站| 色哟哟·www| 桃花免费在线播放| 亚洲欧洲国产日韩| 中文天堂在线官网| 国产高清不卡午夜福利| 性高湖久久久久久久久免费观看| 在线观看一区二区三区激情| 精品亚洲成a人片在线观看| 女的被弄到高潮叫床怎么办| 免费看不卡的av| 秋霞伦理黄片| 男人添女人高潮全过程视频| 亚洲无线观看免费| videosex国产| 免费av不卡在线播放| 人妻制服诱惑在线中文字幕| 岛国毛片在线播放| 国产熟女欧美一区二区| 精品国产一区二区三区久久久樱花| 少妇被粗大的猛进出69影院 | 精品国产一区二区久久| 日韩av在线免费看完整版不卡| 久久亚洲国产成人精品v| 亚洲高清免费不卡视频| 久久99热这里只频精品6学生| 少妇精品久久久久久久| 三上悠亚av全集在线观看| 看十八女毛片水多多多| 国产成人免费无遮挡视频| 精品国产露脸久久av麻豆| 久热久热在线精品观看| 男女免费视频国产| 亚洲欧美一区二区三区国产| 国产精品一国产av| 久久久a久久爽久久v久久| 日韩不卡一区二区三区视频在线| 亚洲国产精品国产精品| 欧美性感艳星| 天天操日日干夜夜撸| 欧美97在线视频| 国内精品宾馆在线| 七月丁香在线播放| 亚洲高清免费不卡视频| 国产欧美亚洲国产| 成年人午夜在线观看视频| 色吧在线观看| 免费人成在线观看视频色| 精品一区二区三区视频在线| 久久久久久伊人网av| 成人免费观看视频高清| 精品99又大又爽又粗少妇毛片| 久久综合国产亚洲精品| 美女大奶头黄色视频| 国产色婷婷99| 女性被躁到高潮视频| 国产欧美另类精品又又久久亚洲欧美| av.在线天堂| 狠狠精品人妻久久久久久综合| 久久午夜福利片| 菩萨蛮人人尽说江南好唐韦庄| 精品少妇黑人巨大在线播放| 91精品国产九色| 国产成人av激情在线播放 | 制服诱惑二区| 国产精品人妻久久久久久| 久热久热在线精品观看| 亚洲四区av| 久久午夜综合久久蜜桃| 国产成人91sexporn| 亚洲av日韩在线播放| 在线精品无人区一区二区三| 中文字幕制服av| 欧美xxxx性猛交bbbb| 国产精品久久久久久av不卡| 一本色道久久久久久精品综合| freevideosex欧美| 亚洲欧美中文字幕日韩二区| 精品一区二区三卡| 一个人看视频在线观看www免费| 国产成人精品福利久久| 日韩中文字幕视频在线看片| 精品久久久久久久久av| freevideosex欧美| 国产成人精品在线电影| av在线老鸭窝| 国产日韩欧美视频二区| 亚洲精品av麻豆狂野| 日韩中文字幕视频在线看片| 久久99蜜桃精品久久| 国产欧美亚洲国产| 日韩电影二区| 日日啪夜夜爽| 人妻少妇偷人精品九色| 国产成人免费无遮挡视频| 高清av免费在线| 日本色播在线视频| 纵有疾风起免费观看全集完整版| 日日摸夜夜添夜夜添av毛片| 丁香六月天网| 特大巨黑吊av在线直播| 曰老女人黄片| 免费观看无遮挡的男女| 国产精品女同一区二区软件| 国产国语露脸激情在线看| 永久网站在线| 中国美白少妇内射xxxbb| 老熟女久久久| 妹子高潮喷水视频| 老司机影院成人| 亚洲精品av麻豆狂野| 丰满乱子伦码专区| 69精品国产乱码久久久| 午夜福利在线观看免费完整高清在| 日本欧美视频一区| 狂野欧美白嫩少妇大欣赏| 天美传媒精品一区二区| 日韩制服骚丝袜av| 久久这里有精品视频免费| 纵有疾风起免费观看全集完整版| 大话2 男鬼变身卡| 日韩人妻高清精品专区| 国产午夜精品久久久久久一区二区三区| 国产av精品麻豆| 18在线观看网站| 美女国产视频在线观看| 久久影院123| 高清欧美精品videossex| 飞空精品影院首页| 少妇的逼水好多| 99国产综合亚洲精品| 黄片无遮挡物在线观看| 亚洲av福利一区| 一二三四中文在线观看免费高清| 狂野欧美白嫩少妇大欣赏| 中文乱码字字幕精品一区二区三区| 免费播放大片免费观看视频在线观看| 51国产日韩欧美| 人人妻人人添人人爽欧美一区卜| 日韩制服骚丝袜av| 久久久午夜欧美精品| 乱码一卡2卡4卡精品| 亚洲人与动物交配视频| 菩萨蛮人人尽说江南好唐韦庄| 亚洲人成77777在线视频| 新久久久久国产一级毛片| 亚洲精品中文字幕在线视频| 精品一区二区三区视频在线| 啦啦啦中文免费视频观看日本| 免费久久久久久久精品成人欧美视频 | 美女国产高潮福利片在线看| 中文字幕最新亚洲高清| 一边亲一边摸免费视频| 免费观看a级毛片全部| 这个男人来自地球电影免费观看 | 这个男人来自地球电影免费观看 | 国产黄色免费在线视频| 亚洲精品,欧美精品| 亚洲精品乱久久久久久| 国国产精品蜜臀av免费| 美女福利国产在线| 性色avwww在线观看| 丝袜脚勾引网站| 性色av一级| 午夜精品国产一区二区电影| 少妇精品久久久久久久| 亚洲人成网站在线观看播放| 国产精品.久久久| 国产av码专区亚洲av| 大话2 男鬼变身卡| 亚洲欧洲日产国产| 中国国产av一级| 中文字幕av电影在线播放| 99久久综合免费| 熟女人妻精品中文字幕| 亚洲国产成人一精品久久久| 久久久国产欧美日韩av| 欧美老熟妇乱子伦牲交| 女性生殖器流出的白浆| 人妻夜夜爽99麻豆av| 国产精品 国内视频| 满18在线观看网站| 日韩av在线免费看完整版不卡| 日韩伦理黄色片| 少妇人妻久久综合中文| av天堂久久9| 亚洲成人手机| 一级黄片播放器| 十八禁高潮呻吟视频| 日韩精品免费视频一区二区三区 | 婷婷色av中文字幕| 久久久久久久大尺度免费视频| 成人二区视频| 在线亚洲精品国产二区图片欧美 | 欧美日韩精品成人综合77777| 免费不卡的大黄色大毛片视频在线观看| 国产精品欧美亚洲77777| 毛片一级片免费看久久久久| 丰满饥渴人妻一区二区三| 少妇精品久久久久久久| 色5月婷婷丁香| 久久女婷五月综合色啪小说| av播播在线观看一区| 国产男女内射视频| 日本爱情动作片www.在线观看| 在现免费观看毛片| 丝袜在线中文字幕| av在线老鸭窝| 久久鲁丝午夜福利片| 国产成人精品一,二区| 日韩视频在线欧美| 久久久国产精品麻豆| 18在线观看网站| 国产一区亚洲一区在线观看| 一区二区三区四区激情视频| a级毛片黄视频| 国产亚洲精品第一综合不卡 | 两个人免费观看高清视频| 国产免费现黄频在线看| 两个人的视频大全免费| 国产国语露脸激情在线看| 视频中文字幕在线观看| 久久久久精品性色| 乱人伦中国视频| 欧美人与善性xxx| 狠狠婷婷综合久久久久久88av| 日本av免费视频播放| 啦啦啦在线观看免费高清www| 久久99热这里只频精品6学生| 成年人午夜在线观看视频| 丝瓜视频免费看黄片| 国产伦理片在线播放av一区| 欧美xxⅹ黑人| 99热国产这里只有精品6| 在线观看国产h片| 久久国内精品自在自线图片| 伊人久久国产一区二区| 十八禁网站网址无遮挡| 成人黄色视频免费在线看| 免费久久久久久久精品成人欧美视频 | 天天影视国产精品| 亚洲性久久影院| 永久免费av网站大全| 精品国产一区二区久久| 综合色丁香网| 一级黄片播放器| 久久婷婷青草| 三级国产精品欧美在线观看| 九九久久精品国产亚洲av麻豆| 日韩成人av中文字幕在线观看| 最近的中文字幕免费完整| 亚洲第一av免费看| 老司机影院毛片| 极品人妻少妇av视频| 一级a做视频免费观看| 国产老妇伦熟女老妇高清| 欧美成人午夜免费资源| 欧美+日韩+精品| 亚洲av国产av综合av卡| 亚洲一级一片aⅴ在线观看| 国产精品久久久久久av不卡| 国产精品久久久久久久久免| 国产毛片在线视频| 亚洲精品自拍成人| 亚洲国产av影院在线观看| 黄片播放在线免费| 国产 精品1| 国产日韩欧美在线精品| 狠狠婷婷综合久久久久久88av| 欧美人与性动交α欧美精品济南到 | 久久影院123| 三级国产精品片| kizo精华| 这个男人来自地球电影免费观看 | 天美传媒精品一区二区| 高清在线视频一区二区三区| 日韩中字成人| 9色porny在线观看| 国产无遮挡羞羞视频在线观看| 久久久久久久久久成人| 国产黄色免费在线视频| 亚洲美女视频黄频| 卡戴珊不雅视频在线播放| 中文字幕免费在线视频6| 国产成人午夜福利电影在线观看| 国产精品久久久久久久电影| 极品人妻少妇av视频| 精品少妇内射三级| 18禁动态无遮挡网站| 日韩成人伦理影院| 久久99精品国语久久久| 亚洲精品国产色婷婷电影| 在线 av 中文字幕| 亚洲欧美色中文字幕在线| 肉色欧美久久久久久久蜜桃| 十分钟在线观看高清视频www| av国产久精品久网站免费入址| 欧美成人午夜免费资源| 大香蕉久久网| 日韩不卡一区二区三区视频在线| 卡戴珊不雅视频在线播放| 国产亚洲欧美精品永久| 高清欧美精品videossex| 丁香六月天网| 国产熟女欧美一区二区| 中国三级夫妇交换| 精品一区二区三区视频在线| 亚洲精品久久成人aⅴ小说 | 母亲3免费完整高清在线观看 | 麻豆乱淫一区二区| 亚洲精品一二三| 国产精品偷伦视频观看了| 精品久久久久久久久av| 亚洲精品第二区| 国产精品国产av在线观看| 一级,二级,三级黄色视频| av天堂久久9| 久久精品久久久久久久性| a级片在线免费高清观看视频| 亚洲一区二区三区欧美精品| 狂野欧美激情性xxxx在线观看| 97精品久久久久久久久久精品| 午夜视频国产福利| 成人毛片a级毛片在线播放| 亚洲国产欧美在线一区| av电影中文网址| 免费少妇av软件| a级毛片在线看网站| 国产日韩一区二区三区精品不卡 | 一区二区日韩欧美中文字幕 | 亚洲精华国产精华液的使用体验| 久久免费观看电影| 九色成人免费人妻av| 国产 精品1| 80岁老熟妇乱子伦牲交| 亚洲精品久久久久久婷婷小说| 国产片内射在线| 99热这里只有是精品在线观看| kizo精华| 国产一区二区在线观看日韩| 王馨瑶露胸无遮挡在线观看| 九九爱精品视频在线观看| 国产成人精品在线电影| 亚州av有码| 精品人妻熟女av久视频| 80岁老熟妇乱子伦牲交| 蜜臀久久99精品久久宅男| 一区二区三区四区激情视频| 狠狠精品人妻久久久久久综合| 人体艺术视频欧美日本| 久久久久视频综合| 亚洲四区av| 国产又色又爽无遮挡免| 欧美bdsm另类| videossex国产| 欧美精品高潮呻吟av久久| 大话2 男鬼变身卡| av一本久久久久| 超碰97精品在线观看| 少妇高潮的动态图| 成人国产av品久久久| 亚洲精品美女久久av网站| 看非洲黑人一级黄片| 久久国产亚洲av麻豆专区| 黑人巨大精品欧美一区二区蜜桃 | 亚洲精品美女久久av网站| 亚洲欧洲国产日韩| 婷婷色综合www| 精品一品国产午夜福利视频| 国产午夜精品久久久久久一区二区三区| 中国国产av一级| 中文字幕制服av| 美女脱内裤让男人舔精品视频| 一级片'在线观看视频| 少妇人妻精品综合一区二区| 国产探花极品一区二区| 搡女人真爽免费视频火全软件| 一本久久精品| 免费看不卡的av| 欧美日韩在线观看h| 久久国产亚洲av麻豆专区| 久久久久久久久久久丰满| 国产av码专区亚洲av| 国产片内射在线| 99久久人妻综合| 免费观看av网站的网址| 日韩av免费高清视频| 国产视频内射| 国产免费现黄频在线看| 久久久久久久久大av| 国内精品宾馆在线| 久久av网站| 最黄视频免费看| 国产精品嫩草影院av在线观看| 成人午夜精彩视频在线观看| 免费观看a级毛片全部| 一级毛片我不卡| 少妇被粗大的猛进出69影院 | 熟女av电影| 少妇被粗大猛烈的视频| av电影中文网址| 欧美激情 高清一区二区三区| 亚洲国产av新网站| 一级毛片我不卡| 青春草视频在线免费观看| 婷婷色综合大香蕉| a级毛片黄视频| 亚洲国产精品999| 国产精品秋霞免费鲁丝片| 国产精品熟女久久久久浪| 免费看不卡的av| 国产国语露脸激情在线看| 亚洲精品日本国产第一区| 中国美白少妇内射xxxbb| 草草在线视频免费看| 五月天丁香电影| 麻豆成人av视频| 美女国产视频在线观看| 午夜激情福利司机影院| 日韩熟女老妇一区二区性免费视频| 国国产精品蜜臀av免费| 国产无遮挡羞羞视频在线观看| 亚洲成人av在线免费| 狂野欧美激情性xxxx在线观看| 黑人巨大精品欧美一区二区蜜桃 | 制服人妻中文乱码| 欧美精品亚洲一区二区| 欧美+日韩+精品| 成人国产av品久久久| 欧美三级亚洲精品| 免费久久久久久久精品成人欧美视频 | 久久99一区二区三区| 综合色丁香网| 一个人免费看片子| 亚洲美女黄色视频免费看| 最后的刺客免费高清国语| 亚洲国产精品999| 色吧在线观看| 丝袜在线中文字幕| 永久网站在线| 中文字幕av电影在线播放| h视频一区二区三区| 搡女人真爽免费视频火全软件| 曰老女人黄片| 国产成人精品福利久久| 亚洲内射少妇av| 在线观看www视频免费| 老司机影院毛片| 黄片播放在线免费| 久久久久久久亚洲中文字幕| 大香蕉久久成人网| 免费日韩欧美在线观看| 大香蕉久久网| 国产亚洲最大av| 成年av动漫网址| 亚洲精品乱码久久久久久按摩| 久久韩国三级中文字幕| 国产在线一区二区三区精| 国产av精品麻豆| 国产亚洲最大av| 亚洲人与动物交配视频| 亚洲精品一二三| xxx大片免费视频| 免费看av在线观看网站| 交换朋友夫妻互换小说| 在线观看免费高清a一片| 午夜福利,免费看| 一本色道久久久久久精品综合| 赤兔流量卡办理| 狠狠婷婷综合久久久久久88av| 国产精品一区二区在线不卡| av在线app专区| 欧美人与性动交α欧美精品济南到 | a级毛片黄视频| 久久女婷五月综合色啪小说| 18禁在线播放成人免费| 99久久精品一区二区三区| 久久久久久久久大av| 另类精品久久| 成人毛片a级毛片在线播放| 丝袜美足系列| 亚洲国产精品999| 亚洲内射少妇av| 亚洲精品乱码久久久v下载方式| av在线观看视频网站免费| 久久久久精品性色| 色网站视频免费| 日日爽夜夜爽网站| 亚洲久久久国产精品| 一区在线观看完整版| 亚洲国产精品专区欧美| 水蜜桃什么品种好| 在线看a的网站| 青青草视频在线视频观看| 在线观看一区二区三区激情| 久久精品久久久久久噜噜老黄| 亚洲人与动物交配视频| 中文字幕人妻丝袜制服| 王馨瑶露胸无遮挡在线观看| 色94色欧美一区二区| 成人无遮挡网站| 黑丝袜美女国产一区| 激情五月婷婷亚洲| 国产色婷婷99| 日本av免费视频播放| 免费黄网站久久成人精品| 3wmmmm亚洲av在线观看| 日韩一区二区三区影片| 午夜日本视频在线| 欧美日韩av久久| 亚洲国产精品一区二区三区在线| 在现免费观看毛片| 美女大奶头黄色视频| 男人爽女人下面视频在线观看| 亚洲欧美精品自产自拍| 最近手机中文字幕大全| 最近中文字幕2019免费版| 18禁在线无遮挡免费观看视频| 成人亚洲欧美一区二区av| 国产老妇伦熟女老妇高清| 亚洲av欧美aⅴ国产| 女性被躁到高潮视频| 国产免费一级a男人的天堂| 亚洲欧美清纯卡通| 成人亚洲欧美一区二区av| 啦啦啦视频在线资源免费观看| 国产在线一区二区三区精| 国产成人免费观看mmmm| 免费黄色在线免费观看| 精品久久久精品久久久| 男女无遮挡免费网站观看| 99热网站在线观看| 国产欧美日韩一区二区三区在线 | 中文字幕精品免费在线观看视频 | 在线天堂最新版资源| 亚洲,一卡二卡三卡| 国产精品99久久久久久久久| 免费大片18禁| 中文字幕人妻熟人妻熟丝袜美| 精品久久久噜噜| 一区在线观看完整版| 最近的中文字幕免费完整| 日韩av不卡免费在线播放| a级毛片黄视频| 边亲边吃奶的免费视频| 老司机影院成人| 男女免费视频国产| 国产在视频线精品| 成人亚洲欧美一区二区av| 男男h啪啪无遮挡| 久久久久久久久久久免费av| 久久99热这里只频精品6学生| 国产在线视频一区二区| 观看美女的网站| 亚洲久久久国产精品| 久久久久久久精品精品| 尾随美女入室| 亚洲情色 制服丝袜| 亚洲av免费高清在线观看| 亚洲人成77777在线视频| 精品少妇内射三级| 国产男女超爽视频在线观看| 男男h啪啪无遮挡| 成人亚洲欧美一区二区av| 国产精品久久久久久av不卡| 色94色欧美一区二区| 能在线免费看毛片的网站| 精品国产露脸久久av麻豆| 精品少妇内射三级| 91精品三级在线观看| 一级a做视频免费观看| 自线自在国产av| 国产精品 国内视频| 久久国产亚洲av麻豆专区| 国产淫语在线视频| 日韩伦理黄色片| 国产乱来视频区| 午夜福利在线观看免费完整高清在| 久久久久久久精品精品| 午夜免费观看性视频| 99热国产这里只有精品6| 日韩一区二区三区影片| 久久99精品国语久久久| 99热6这里只有精品| 免费高清在线观看视频在线观看| 在线观看www视频免费| 哪个播放器可以免费观看大片| 纯流量卡能插随身wifi吗| 亚洲av成人精品一二三区| 亚洲欧美色中文字幕在线| 大片电影免费在线观看免费| 午夜免费男女啪啪视频观看| 人人澡人人妻人| 国产成人精品福利久久| 最后的刺客免费高清国语| av视频免费观看在线观看| 只有这里有精品99| 97精品久久久久久久久久精品| 日本午夜av视频| 少妇精品久久久久久久| 国产高清国产精品国产三级| 成年人免费黄色播放视频| 大香蕉97超碰在线| 91久久精品电影网| 国产无遮挡羞羞视频在线观看| 亚洲av男天堂| 人妻 亚洲 视频|