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

    求解中途飛越燃料最優(yōu)轉(zhuǎn)移軌道的同倫方法

    2017-05-18 07:39:18泮斌峰
    宇航學(xué)報 2017年4期
    關(guān)鍵詞:內(nèi)點性能指標(biāo)小行星

    潘 迅,泮斌峰,唐 碩

    (1. 西北工業(yè)大學(xué)航天學(xué)院,西安710072;2. 航天飛行動力學(xué)技術(shù)重點實驗室,西安 710072)

    求解中途飛越燃料最優(yōu)轉(zhuǎn)移軌道的同倫方法

    潘 迅1,2,泮斌峰1,2,唐 碩1

    (1. 西北工業(yè)大學(xué)航天學(xué)院,西安710072;2. 航天飛行動力學(xué)技術(shù)重點實驗室,西安 710072)

    提出一種新的同倫方法,用于求解深空探測中對其他天體進(jìn)行中途飛越的小推力燃料最優(yōu)轉(zhuǎn)移軌道,克服由于其存在內(nèi)點約束及不連續(xù)Bang-Bang控制所導(dǎo)致的數(shù)值優(yōu)化方法的求解困難。該同倫方法將同倫參數(shù)同時嵌入到性能指標(biāo)和內(nèi)點約束方程中,將容易求解的無內(nèi)點約束且控制量連續(xù)變化的最優(yōu)控制問題作為初始問題,求解一系列同倫參數(shù)遞增所對應(yīng)的同倫迭代子問題,直到得到原問題的解。該方法能夠有效地解決中途飛越所導(dǎo)致的優(yōu)化變量增加、求解難度增大等難題,能夠快速、穩(wěn)定地求解考慮中途飛越的小推力燃料最優(yōu)轉(zhuǎn)移軌道。最后,以地球到火星交會并中途飛越小行星和地球到木星交會并中途飛越火星兩個任務(wù)為例進(jìn)行數(shù)值仿真驗證該同倫方法在求解中途飛越的燃料最優(yōu)問題中的有效性和優(yōu)越性。

    中途飛越;最優(yōu)控制;燃料最優(yōu);小推力;同倫方法;深空探測

    0 引 言

    深空探測是21世紀(jì)航天領(lǐng)域的一大熱點。從經(jīng)濟(jì)性角度出發(fā),多任務(wù)、多目標(biāo)深空探測能夠在節(jié)省經(jīng)費(fèi)的同時獲得更多的科學(xué)回報。與單目標(biāo)探測的航天任務(wù)相比,多目標(biāo)探測任務(wù)對轉(zhuǎn)移軌道設(shè)計提出了新的要求。

    對天體的探測方式包括飛越、交會或撞擊目標(biāo)天體。對于多目標(biāo)探測任務(wù),對其中部分天體進(jìn)行飛越是較為合理的探測方式。文獻(xiàn)[1]對脈沖推進(jìn)的航天任務(wù),提出基于牛頓算法的多次中途飛越任務(wù)的轉(zhuǎn)移軌道優(yōu)化方法,利用該方法能使Galileo探測器完成多達(dá)11次中途飛越,從而對多個目標(biāo)天體進(jìn)行探測;文獻(xiàn)[2]分析Deep Space探測器對Hartley完成兩次飛越后,通過對剩余燃料進(jìn)行合理規(guī)劃,對小行星2002 GT進(jìn)行飛越的任務(wù)設(shè)計;文獻(xiàn)[3]研究航天器通過脈沖推進(jìn)對火星進(jìn)行飛越,然后跟隨火星運(yùn)行半圈,最后通過脈沖飛越火星返回地球。我國的嫦娥2號任務(wù)在完成對月探測后,也對小行星進(jìn)行了飛越探測。文獻(xiàn)[4]針對嫦娥二號從日地系L2點Lissajous軌道到Toutatis飛越的轉(zhuǎn)移軌道進(jìn)行設(shè)計研究;文獻(xiàn)[5]研究嫦娥二號在實現(xiàn)對4179小行星飛越之后,利用剩余燃料選取小行星1997XF11和2005VS作為飛越探測目標(biāo)的任務(wù)設(shè)計。

    隨著發(fā)動機(jī)技術(shù)的發(fā)展,小推力發(fā)動機(jī)具有比沖大、控制精度高等優(yōu)點,能極大地提高探測器的有效載荷的比重,因此利用小推力發(fā)動機(jī)進(jìn)行深空探測將是未來的發(fā)展趨勢。然而小推力發(fā)動機(jī)推力小,持續(xù)時間長,燃料最優(yōu)轉(zhuǎn)移軌道分為推進(jìn)段和滑行段,在發(fā)動機(jī)開關(guān)機(jī)點推力大小存在階躍變化,因此在求解過程中存在較大困難。文獻(xiàn)[6]研究了基于不變流形理論和Gauss偽譜法的日-地系Halo軌道到日-火系Halo軌道的小推力轉(zhuǎn)移;文獻(xiàn)[7]綜述了當(dāng)前軌道優(yōu)化設(shè)計領(lǐng)域的常用方法,針對連續(xù)小推力作用下的深空探測任務(wù)軌道優(yōu)化方法,對比和分析各方法的優(yōu)缺點;文獻(xiàn)[8]采用六次多項式形狀逼近策略,并結(jié)合改進(jìn)微分進(jìn)化算法進(jìn)行小推力行星借力轉(zhuǎn)移軌道的初始設(shè)計;文獻(xiàn)[9]結(jié)合我國即將開展的小行星探測任務(wù),給出多目標(biāo)、多任務(wù)小行星探測的一種推薦方案,分析得到6個可能的目標(biāo)小行星,并對可行方案進(jìn)行電推進(jìn)轉(zhuǎn)移軌道設(shè)計;文獻(xiàn)[10]推導(dǎo)了包含時變內(nèi)點約束的分塊偽譜法,并對地球出發(fā)-小行星飛越-地球返回的小推力轉(zhuǎn)移軌道進(jìn)行設(shè)計;文獻(xiàn)[11] 通過間接法對包含引力輔助的燃料最優(yōu)小推力轉(zhuǎn)移軌道進(jìn)行研究,利用粒子群算法選取初值,利用同倫算法和開關(guān)函數(shù)檢測法克服Bang-Bang控制的困難。

    對于考慮中途飛越的航天任務(wù)的轉(zhuǎn)移軌道,由于存在內(nèi)點約束,其在優(yōu)化過程中引入多個對應(yīng)于約束條件的拉格朗日乘子,增加了打靶所需求解的變量,從而導(dǎo)致求解難度增加。雖然文獻(xiàn)[8-11]能優(yōu)化得到滿足要求的轉(zhuǎn)移軌道,但在優(yōu)化過程中優(yōu)化變量較多,需要結(jié)合微分進(jìn)化算法,粒子群法等其他數(shù)值算法,計算過程較為復(fù)雜。本文提出一種新的同倫方法,先求解不包含內(nèi)點約束的轉(zhuǎn)移軌道,在該軌道上選取合適位置作為初始中途飛越點,從而在同倫初始時刻消除內(nèi)點約束所產(chǎn)生的影響;然后通過對飛越點位置進(jìn)行同倫,進(jìn)而得到對目標(biāo)天體進(jìn)行中途飛越的轉(zhuǎn)移軌道;同時對性能指標(biāo)進(jìn)行同倫,從而克服燃料最優(yōu)問題Bang-Bang控制所導(dǎo)致的求解困難。

    1 中途飛越的最優(yōu)化問題

    在進(jìn)行航天器軌道初步設(shè)計時,不考慮其他天體引力、太陽光壓以及其他攝動力,即假設(shè)航天器在中心引力場中運(yùn)動,僅受到中心天體的引力作用和發(fā)動機(jī)推力作用,則航天器的動力學(xué)方程為

    (1)

    為方便計算以及保證計算精度,通常需要對動力學(xué)方程中的變量進(jìn)行無量綱化。在日心黃道坐標(biāo)系中,長度量、時間量和質(zhì)量分別以AU(1AU=149597870700m)、a(1a=365.25d)、航天器初始質(zhì)量m0進(jìn)行無量綱化,太陽引力系數(shù)無量綱化后為39.476926AU3/a,同時需要對Tmax、Isp和g0進(jìn)行無量綱化[9]。

    對于中途飛越問題,除了在初始時刻和終端時刻存在約束外,在航天器進(jìn)行中途飛越時刻也存在約束,即內(nèi)點約束

    r(tm)-rO(tm)=0

    (2)

    式中:tm為中途飛越時刻,rO表示中途飛越目標(biāo)天體的位置矢量。

    將端點約束和內(nèi)點約束表示成p維的等式約束

    φ(x(t0),x(tf),x(tm),t0,tf,tm)=0

    (3)

    式中:x=[r,v,m]T表示狀態(tài)量,t0和tf分別表示初始時刻和結(jié)束時刻。

    對軌道設(shè)計問題,優(yōu)化的性能指標(biāo)通常分為三種:時間最優(yōu),能量最優(yōu)和燃料最優(yōu),都只與推力大小u有關(guān),與推力方向α無關(guān),因此,可將性能指標(biāo)表示為

    (4)

    根據(jù)龐德里亞金極小值原理,通過引入拉格朗日乘子構(gòu)造廣義指標(biāo)Jm和哈密爾頓函數(shù)H,可得

    (5)

    (6)

    對于控制無約束的情況,廣義指標(biāo)的一階變分為零是性能指標(biāo)取得極值的一階必要條件。通過對廣義指標(biāo)求一階變分,可得協(xié)態(tài)方程為

    (7)

    根據(jù)最優(yōu)控制理論,初始狀態(tài)和終端狀態(tài)的橫截條件為[12]

    (8)

    內(nèi)點約束時刻前后瞬時狀態(tài)約束對應(yīng)的橫截條件為

    (9)

    內(nèi)點約束時刻的靜態(tài)條件為

    (10)

    式中:vO為中途飛越目標(biāo)天體的速度矢量。

    根據(jù)最優(yōu)控制理論,可推導(dǎo)得到最優(yōu)推力方向為

    (11)

    在推導(dǎo)推力大小u的最優(yōu)控制律時,需要考慮不同性能指標(biāo)的影響。在此,分別對性能指標(biāo)為能量最優(yōu)和燃料最優(yōu)時的推力大小最優(yōu)控制律進(jìn)行推導(dǎo):

    1)當(dāng)性能指標(biāo)為能量最優(yōu)時,其最小化性能指標(biāo)可以表示成

    (12)

    將其代入到哈密爾頓函數(shù),并對u求偏導(dǎo),可得

    (13)

    為使哈密爾頓函數(shù)取得極小值,可得推力大小的最優(yōu)控制律為

    (14)

    式中:S為開關(guān)函數(shù),其表達(dá)式為

    (15)

    2)當(dāng)性能指標(biāo)為燃料最優(yōu)時,其最小化性能指標(biāo)可以表示成

    (16)

    同理可得其最優(yōu)控制律為

    (17)

    當(dāng)S=0時,最優(yōu)推力大小不能確定,此時存在奇異控制。由于開關(guān)函數(shù)S連續(xù)變化,不會在某一時間區(qū)間內(nèi)恒等于零,只會在有限孤立點處為零,因此可以不考慮S=0的情況。

    推導(dǎo)得到推力方向和推力大小的最優(yōu)控制律后,則將軌道優(yōu)化問題轉(zhuǎn)換為兩點邊值問題。對于無飛越的軌道轉(zhuǎn)移,其終端約束為航天器位置速度與目標(biāo)天體的位置速度一致,且不存在內(nèi)點約束,則該打靶問題為7維問題,打靶變量為初始時刻的協(xié)態(tài)變量[λr(t0),λv(t0),λm(t0)],打靶約束條件為

    (18)

    對于存在中途飛越的軌道轉(zhuǎn)移問題,在中途飛越時刻,航天器位置與飛越天體位置一致,即每個中途飛越對應(yīng)三個約束條件,同時由于中途飛越時刻tm未知,則該問題為十一維的打靶問題,打靶變量為[λr(t0),λv(t0),λm(t0),χ,tm],打靶方程為式(2),式(18)和式(10),共11個等式約束,式(9)用于協(xié)態(tài)變量λr在tm時刻的更新。通過對該打靶問題進(jìn)行求解,則可得到滿足條件的轉(zhuǎn)移軌道。

    2 同倫方法

    對于燃料最優(yōu)轉(zhuǎn)移軌道,從式(17)可知,其推力值u為Bang-Bang控制。由于其推力大小存在突變,利用間接法對其進(jìn)行求解過程中存在很大困難,因此通常先求解得到能量最優(yōu)轉(zhuǎn)移軌道,通過對性能指標(biāo)進(jìn)行能量最優(yōu)到燃料最優(yōu)的同倫迭代,最終得到燃料最優(yōu)轉(zhuǎn)移軌道[11]。但在考慮中途飛越問題時,由于內(nèi)點約束所帶來的在飛越時刻協(xié)態(tài)變量和哈密爾頓函數(shù)的不連續(xù)等問題,上述同倫方法則難以適用。本文針對這一問題,提出一種新的同倫方法進(jìn)行求解。

    2.1 中途飛越內(nèi)點約束方程的同倫

    考慮初始時刻和終端時刻的狀態(tài)量都給定的問題,端點約束對狀態(tài)量的偏導(dǎo)數(shù)為零,因此不需要對端點約束對應(yīng)的拉格朗日乘子進(jìn)行求解。中途飛越的轉(zhuǎn)移軌道優(yōu)化問題只增加了三維的內(nèi)點約束,同時增加了三維的拉格朗日乘子χ。由上文的推導(dǎo)過程中可知,在打靶過程中,與乘子χ有關(guān)的式并不作為打靶方程的約束,而是在tm時刻對協(xié)態(tài)變量λr進(jìn)行更新??紤]到若能在減小χ初值猜測難度的情況下,找到另外一條中間軌道作為初始的中途飛越軌道,則可對中途飛越軌道進(jìn)行平滑過渡,使航天器最終實現(xiàn)對目標(biāo)天體的軌道進(jìn)行飛越。

    將航天器需要進(jìn)行中途飛越的目標(biāo)天體軌道用經(jīng)典軌道六要素Ef=[af,ef,if,ωf,Ωf,Mf]來表示,選取初始飛越軌道為Ei=[ai,ei,ii,ωi,Ωi,Mi],軌道要素中的平近點角M都為初始時刻t0對應(yīng)的平近點角。具體的處理方法為:

    2.2 性能指標(biāo)的同倫

    從式中可以看出,當(dāng)性能指標(biāo)為能量最優(yōu)時,其最優(yōu)推力大小是連續(xù)變化的,此時較容易求解。因此,針對燃料最優(yōu)的Bang-Bang控制問題,Bertrand等在2002年提出了一種平滑技術(shù),即所謂的同倫方法,通過構(gòu)造新的性能指標(biāo),降低了求解難度[13]。

    引入同倫參數(shù)q2,構(gòu)造新的性能指標(biāo)為

    (19)

    當(dāng)q2:0→1過程中,性能指標(biāo)從能量最優(yōu)轉(zhuǎn)變?yōu)槿剂献顑?yōu)。根據(jù)最優(yōu)性一階必要條件,推力大小u的控制律為

    (20)

    式中:ε=1-q2,當(dāng)q2≠1時,u為連續(xù)控制力,且q2趨近于1時,即可得到燃料最優(yōu)的轉(zhuǎn)移軌道。

    在得到能量最優(yōu)問題的轉(zhuǎn)移軌道之后,同倫參數(shù)q2按照一定步長進(jìn)行迭代,并將當(dāng)前得到的解作為下一步迭代的初值進(jìn)行計算,并最終得到燃料最優(yōu)的轉(zhuǎn)移軌道。

    3 數(shù)值仿真

    本節(jié)通過兩個算例進(jìn)行數(shù)值仿真,驗證本文所提出的同倫方法的有效性。

    算例1.考慮航天器從地球出發(fā),與火星進(jìn)行交會,并對小行星2005UK5中途飛越的軌道轉(zhuǎn)移問題。航天器從地球出發(fā)時與地球具有相同的日心位置和速度,到達(dá)火星時與火星具有相同的日心位置和速度。小行星2005UK5的數(shù)據(jù)來源于小天體動力學(xué)網(wǎng)http://newton.dm.unipi.it/neodys/。算例相關(guān)參數(shù)如表1所示。

    表1 地球到火星交會并中途飛越小行星的算例的參數(shù)

    先通過打靶求解無中途飛越時的從地球到火星的能量最優(yōu)轉(zhuǎn)移軌道,根據(jù)航天器在轉(zhuǎn)移過程中與小行星的相對距離進(jìn)行初始飛越點的選取。選取轉(zhuǎn)移軌道上與小行星距離最近的時刻作為初始飛越時刻,轉(zhuǎn)移軌道上的該點作為初始飛越點,然后分兩步進(jìn)行同倫迭代求解:1.保持性能指標(biāo)為能量最優(yōu),對中途飛越點進(jìn)行同倫,即q2=0,q1:0→1;2.得到對小行星進(jìn)行飛越的能量最優(yōu)轉(zhuǎn)移軌道后,再對性能指標(biāo)進(jìn)行同倫,即q1=1,q2:0→1,從而得到航天器對小行星2005UK5進(jìn)行中途飛越的燃料最優(yōu)轉(zhuǎn)移軌道。

    最終得到的轉(zhuǎn)移軌道如圖1所示,飛越時刻為294.9028天,航天器最終剩余質(zhì)量為1649.41kg。圖1中的黑色虛線和紅色實線分別表示同倫初始時刻的轉(zhuǎn)移軌道和同倫結(jié)束時刻的燃料最優(yōu)轉(zhuǎn)移軌道,同倫結(jié)束時刻對小行星2005UK5實現(xiàn)了中途飛越,飛越位置已在圖中標(biāo)出。圖2為優(yōu)化變量在對中途飛越位置同倫過程中隨同倫參數(shù)q1的變化曲線,圖3為優(yōu)化變量在對性能指標(biāo)的同倫過程中隨同倫參數(shù)q2的變化曲線,從圖中可以看出,同倫參數(shù)連續(xù)緩慢變化,同倫過程順利進(jìn)行。圖4為航天器在轉(zhuǎn)移過程中開關(guān)函數(shù)和推力隨時間的變化曲線,對于q1=q2=0和q1=1,q2=0的兩條轉(zhuǎn)移軌道,其性能指標(biāo)同為能量最優(yōu),兩條軌道形狀較為接近,開關(guān)函數(shù)和推力大小變化曲線無明顯差異;對于q1=q2=1的燃料最優(yōu)轉(zhuǎn)移軌道,其推力大小變化為Bang-Bang控制,發(fā)動機(jī)共開機(jī)4次。

    描述量數(shù)值初始時刻(UTC)2021-11-160:0:0.0飛行時間ttof/d2201.0天初始位置/AU[0.5876420, 0.7954627, -3.845203×10-5]初始速度[-5.155764,3.707833, -3.191945×10-4]末端位置[-5.204974,1.495369,0.1102444]末端速度[-0.7936872,-2.523063,2.823220×10-2]發(fā)動機(jī)比沖Isp/s6000發(fā)動機(jī)最大推力Tmax/N2.26航天器初始質(zhì)量m0/kg20000.0火星軌道根數(shù)[a,e,i,Ω,ω,f]J2000日心黃道坐標(biāo)系中,參考?xì)v元[1.52363312AU,0.09327933,1.84785414°,49.48935357°,286.6709081°,328.887552°]2024-3-200:0:0.0(坐標(biāo)時)

    值得注意的是,在同倫過程中,可設(shè)置兩個同倫參數(shù)相等,即q1=q2,此時則可通過一次同倫過程得到所需的解。在此算例中,令q=q1=q2,同時對中途飛越位置和性能指標(biāo)進(jìn)行同倫,通過一次同倫得到中途飛越火星并與木星交會的燃料最優(yōu)轉(zhuǎn)移軌道。

    最終得到的航天器燃料最優(yōu)轉(zhuǎn)移軌道如圖5所示,航天器對火星進(jìn)行中途飛越的時刻為732.8713天,最終剩余質(zhì)量為14657.05kg,其中黑色虛線為同倫初始時刻的轉(zhuǎn)移軌道,該軌道為能量最優(yōu)轉(zhuǎn)移軌道,對其進(jìn)行求解時不考慮中途飛越內(nèi)點約束,然后在該軌道上進(jìn)行選取合適的初始中途飛越位置;紅色實線表示同倫結(jié)束時的轉(zhuǎn)移軌道,該軌道為對火星實現(xiàn)中途飛越的燃料最優(yōu)轉(zhuǎn)移軌道。圖6為航天器進(jìn)行中途飛越位置隨同倫參數(shù)的變化曲線,從圖中可以看出,當(dāng)同倫初始時刻時,航天器并未經(jīng)過火星,而是對虛擬的目標(biāo)進(jìn)行飛越;隨著同倫過程的進(jìn)行,最終實現(xiàn)對火星的飛越。圖7為優(yōu)化變量隨同倫參數(shù)的變化曲線,從圖中可以看出,各優(yōu)化變量變化平緩,不存在突變或跳躍的情況。圖8為當(dāng)同倫參數(shù)q分別取值為0,0.4,0.7,0.9和1.0時的開關(guān)函數(shù)和推力大小隨時間的變化曲線,從圖中可以看出,對于燃料最優(yōu)轉(zhuǎn)移軌道,其最優(yōu)推力為Bang-Bang控制,在轉(zhuǎn)移過程中發(fā)動機(jī)開機(jī)4次。

    4 結(jié) 論

    本文針對航天器需要對目標(biāo)天體進(jìn)行交會,并對其他天體進(jìn)行中途飛越的連續(xù)小推力燃料最優(yōu)轉(zhuǎn)移軌道進(jìn)行研究。本文建立了含中途飛越約束的動力學(xué)模型,推導(dǎo)得到了最優(yōu)性條件,提出的同倫方法包含對中途飛越軌道的同倫和對性能指標(biāo)的同倫。通過數(shù)值仿真表明,利用該同倫方法能有效地克服由中途飛越約束所引起的困難,并得到滿足任務(wù)需求,發(fā)動機(jī)多次開關(guān)機(jī)的燃料最優(yōu)轉(zhuǎn)移軌道。后續(xù)工作中,可對包含引力輔助等更復(fù)雜的內(nèi)點約束的小推力轉(zhuǎn)移軌道進(jìn)行研究。

    [1]D'AmarioL,ByrnesD,StanfordR.Anewmethodforoptimizingmultiple-flybytrajectories[J].JournalofGuidance,Control,andDynamics, 1981, 4(5): 591-596.

    [2]GrebowD,BhaskaranS,ChesleyS.Search&selectionforfutureflybytargetsfortheDI/EPOXIspacecraft[C].AIAA/AASAstrodynamicsSpecialistConference,Minnesota,USA.August13-16, 2012.

    [3]JesickM.Marsdouble-flybyfreereturns[J].JournalofSpacecraftandRockets, 2015, 52(5): 1348-1360.

    [4] 喬棟, 黃江川, 崔平遠(yuǎn), 等. 嫦娥二號衛(wèi)星飛越Toutatis小行星轉(zhuǎn)移軌道設(shè)計[J]. 中國科學(xué): 技術(shù)科學(xué), 2013, 43(5): 487-492. [QiaoDong,HuangJiang-chuan,CuiPing-yuan,etal.TrajectorydesignofCE-2flybyToutatisasteroidmission[J].ScientiaSinicaTechnologica, 2013, 43(5): 487-492.]

    [5] 劉磊, 吳偉仁, 唐歌實, 等. “嫦娥二號” 后續(xù)小行星飛越探測任務(wù)設(shè)計[J]. 國防科技大學(xué)學(xué)報, 2014, 36(2): 13-17. [LiuLei,WuWei-ren,TangGe-shi,etal.Designofanasteroidflying-bymissionforCHANG’E-2 [J].JournalofNationalUniversityofDefenseTechnology, 2014, 36(2): 13-17.]

    [6] 曹喜濱, 張相宇, 王峰. 采用Gauss偽譜法的小推力日-火Halo軌道轉(zhuǎn)移優(yōu)化設(shè)計[J]. 宇航學(xué)報, 2013, 34(8):1047-1054. [CaoXi-bin,ZhangXiang-yu,WangFeng.Optimizationoflow-thrusttransfertrajectoryfortheSun-Marshaloorbitbasedongausspseudospectralmethod[J].JournalofAstronautics, 2013, 34(8):1047-1054. ]

    [7] 李俊峰, 蔣方華. 連續(xù)小推力航天器的深空探測軌道優(yōu)化方法綜述[J].力學(xué)與實踐,2011,33(3):1-6. [LiJun-feng,JiangFang-hua.Surveyoflow-thrusttrajectoryoptimizationmethodsfordeepspaceexploration[J].MechanicsinEngineering, 2011, 33(3):1-6.]

    [8] 尚海濱, 崔平遠(yuǎn), 徐瑞,等. 結(jié)合行星借力飛行技術(shù)的小推力轉(zhuǎn)移軌道初始設(shè)計[J]. 宇航學(xué)報, 2011, 32(1):29-38. [ShangHai-bin,CuiPing-yuan,XuRui,etal.Preliminarydesignoflow-thrusttransfertrajectorywithplanetaryswing-bys[J].JournalofAstronautics, 2011, 32(1):29-38.]

    [9] 陳楊, 寶音賀西, 李俊峰. 我國小行星探測目標(biāo)分析與電推進(jìn)軌道設(shè)計[J]. 中國科學(xué) 物理學(xué) 力學(xué) 天文學(xué), 2011, 41(9): 1104-1111. [ChenYang,Bao-yinHe-xi,LiJun-feng.Targetanalysisandlow-thrusttrajectorydesignofChineseasteroidexplorationmission[J].ScienceChinaPhysics,Mechanics&Astronomy, 2011, 41(9): 1104-1111.]

    [10] 郭鐵丁. 深空探測小推力軌跡優(yōu)化的間接法與偽譜法研究[D].北京: 清華大學(xué), 2012. [GuoTie-ding.Studyofindirectandpseudospecrtralmethodsforlowthrusttrajectoryoptimizationindeepspaceexploration[D].Beijing:TsinghuaUniversity, 2012.]

    [11]JiangF,BaoyinH,LiJ.Practicaltechniquesforlow-thrusttrajectoryoptimizationwithhomotopicapproach[J].JournalofGuidance,Control,andDynamics, 2012, 35(1): 245-258.

    [12] 李俊峰, 寶音賀西, 蔣方華, 等. 深空探測動力學(xué)與控制[M]. 北京: 清華大學(xué)出版社, 2014.[LiJun-feng,Bao-yinHe-xi,JiangFang-hua,etal.Dynamicsandcontrolofinterplanetaryflight[M].Beijing:TsinghuaUniversityPress, 2014.]

    [13]BertrandR,EpenoyR.Newsmoothingtechniquesforsolvingbang-bangoptimalcontrolproblemsnumericalresultsandstatisticalinterpretation[J].OptimalControlApplicationsandMethods, 2002, 23(4): 171-197.

    通信地址:陜西省西安市碑林區(qū)友誼西路127號西北工業(yè)大學(xué)251信箱(710072)

    電話:15191433469

    E-mail:xpan2012@gmail.com

    (編輯:張宇平)

    Homotopy Method for Fuel-Optimal Trajectory Design in Flyby Mission

    PAN Xun1,2, PAN Bin-feng1,2, TANG Shuo1

    (1. School of Astronautics, Northwestern Polytechnical University, Xi’an 710072, China; 2. National Key Laboratory of Aerospace Flight Dynamics, Xi’an 710072, China)

    A new homotopy method is proposed to solve the fuel-optimal transfer trajectory with flyby of other celestial bodies in the deep space exploration, overcoming the difficulty due to the interior point constraints and discontinuous structure of the Bang-Bang control. The homotopic parameter is embedded in both the performance index and interior point constraint equations. By starting from a related and easier energy-optimal problem without inner constraints, and solving a series of iterative subproblems with homotopic parameter increased, the solution of the original problem is obtained when the parameter reaches one. With the approach, the difficulty of the optimization variables increased caused by the flyby inner constraints can be resolved, and the low-thrust fuel-optimal transfers can be completed effectively. At last, two examples about the Jupiter rendezvous with the Mars flyby and the Mars rendezvous with an asteroid flyby are given to substantiate the effectiveness and optimality of the homotopy method.

    Flyby; Optimal control; Fuel optimal; Low thrust; Homotopy method; Deep space exploration

    2016-07-20;

    2016-12-19

    V448.2

    A

    1000-1328(2017)04-0393-08

    10.3873/j.issn.1000-1328.2017.04.009

    潘 迅(1990-),男,博士生,主要研究方向為飛行動力學(xué)與控制,深空探測軌道優(yōu)化設(shè)計。

    猜你喜歡
    內(nèi)點性能指標(biāo)小行星
    NASA宣布成功撞擊小行星
    軍事文摘(2022年24期)2023-01-05 03:38:22
    我國發(fā)現(xiàn)2022年首顆近地小行星
    瀝青膠結(jié)料基本高溫性能指標(biāo)相關(guān)性研究
    石油瀝青(2021年1期)2021-04-13 01:31:08
    基于罰函數(shù)內(nèi)點法的泄露積分型回聲狀態(tài)網(wǎng)的參數(shù)優(yōu)化
    儲熱水箱分層性能指標(biāo)的研究進(jìn)展
    小行星:往左走
    太空探索(2016年1期)2016-07-12 09:55:54
    WebGIS關(guān)鍵性能指標(biāo)測試技術(shù)研究
    基于內(nèi)點方法的DSD算法與列生成算法
    一個新的求解半正定規(guī)劃問題的原始對偶內(nèi)點算法
    “隼鳥”2再探小行星
    太空探索(2014年11期)2014-07-12 15:17:00
    狠狠婷婷综合久久久久久88av| 国产免费视频播放在线视频| 国产亚洲av高清不卡| 黄色 视频免费看| h视频一区二区三区| 一级毛片黄色毛片免费观看视频| 精品久久久久久电影网| 天天躁日日躁夜夜躁夜夜| 国产一级毛片在线| 黄色视频不卡| 精品一区二区三区四区五区乱码 | 丝袜脚勾引网站| 欧美性长视频在线观看| 国产一级毛片在线| 久久久国产欧美日韩av| 在线精品无人区一区二区三| 大香蕉久久网| 久久精品成人免费网站| 两人在一起打扑克的视频| 黑人欧美特级aaaaaa片| 不卡av一区二区三区| 欧美 日韩 精品 国产| 国产欧美日韩精品亚洲av| 看免费成人av毛片| 大话2 男鬼变身卡| 欧美精品一区二区免费开放| 又粗又硬又长又爽又黄的视频| 天堂8中文在线网| 亚洲一区中文字幕在线| 黄色 视频免费看| 亚洲第一av免费看| 亚洲国产日韩一区二区| 国产成人精品久久久久久| 中文欧美无线码| 色婷婷av一区二区三区视频| 在线观看免费高清a一片| 每晚都被弄得嗷嗷叫到高潮| 久久久国产精品麻豆| 中文字幕人妻熟女乱码| 我的亚洲天堂| 91字幕亚洲| 19禁男女啪啪无遮挡网站| 观看av在线不卡| 黄色怎么调成土黄色| 性色av一级| 欧美激情 高清一区二区三区| 亚洲欧美一区二区三区黑人| 国产精品三级大全| 国产精品一二三区在线看| 18禁国产床啪视频网站| 免费在线观看黄色视频的| 看免费av毛片| 久久精品国产a三级三级三级| 国产亚洲欧美在线一区二区| 色婷婷av一区二区三区视频| 曰老女人黄片| 人妻 亚洲 视频| 亚洲黑人精品在线| 国产日韩一区二区三区精品不卡| 久久久久久久国产电影| 欧美日韩亚洲综合一区二区三区_| 亚洲伊人久久精品综合| 欧美精品啪啪一区二区三区 | h视频一区二区三区| 电影成人av| 一级毛片电影观看| 老汉色∧v一级毛片| 精品少妇黑人巨大在线播放| 别揉我奶头~嗯~啊~动态视频 | 中文字幕制服av| 岛国毛片在线播放| 操美女的视频在线观看| 国产片特级美女逼逼视频| 欧美97在线视频| 亚洲国产精品一区二区三区在线| 伦理电影免费视频| 老司机靠b影院| 国产在线视频一区二区| 少妇人妻久久综合中文| 一级毛片我不卡| 国产日韩欧美亚洲二区| 亚洲欧美精品自产自拍| 啦啦啦在线观看免费高清www| 在线 av 中文字幕| 我的亚洲天堂| 亚洲欧美激情在线| 黄片小视频在线播放| 欧美亚洲 丝袜 人妻 在线| 黑人欧美特级aaaaaa片| av线在线观看网站| 最近最新中文字幕大全免费视频 | 亚洲国产精品999| 一区二区三区精品91| 亚洲情色 制服丝袜| 免费看av在线观看网站| 久久人妻福利社区极品人妻图片 | 欧美激情 高清一区二区三区| 成人手机av| 久久久久久久精品精品| www.自偷自拍.com| 国产成人精品久久二区二区91| 国产精品 欧美亚洲| 国产视频一区二区在线看| 亚洲男人天堂网一区| 免费观看a级毛片全部| 啦啦啦视频在线资源免费观看| 91字幕亚洲| 国产熟女欧美一区二区| 丰满人妻熟妇乱又伦精品不卡| 亚洲精品成人av观看孕妇| 欧美人与善性xxx| 色精品久久人妻99蜜桃| 在线av久久热| 久久ye,这里只有精品| 午夜免费鲁丝| 久久影院123| 午夜日韩欧美国产| 一边摸一边抽搐一进一出视频| 飞空精品影院首页| 日韩制服丝袜自拍偷拍| av视频免费观看在线观看| 99国产精品免费福利视频| 欧美xxⅹ黑人| 91国产中文字幕| 成年女人毛片免费观看观看9 | 色视频在线一区二区三区| a 毛片基地| 人人妻,人人澡人人爽秒播 | 欧美日韩福利视频一区二区| 中文字幕高清在线视频| 久久狼人影院| 亚洲色图 男人天堂 中文字幕| 91精品国产国语对白视频| 精品亚洲乱码少妇综合久久| 99国产精品一区二区蜜桃av | 国产免费现黄频在线看| 午夜久久久在线观看| 亚洲精品第二区| 高清视频免费观看一区二区| 国产精品亚洲av一区麻豆| 亚洲九九香蕉| 狂野欧美激情性xxxx| 欧美日韩精品网址| 老汉色av国产亚洲站长工具| 大片电影免费在线观看免费| 国产精品久久久久成人av| 欧美成人午夜精品| 国产精品欧美亚洲77777| 国产成人a∨麻豆精品| 亚洲专区国产一区二区| 国产亚洲av片在线观看秒播厂| 久久鲁丝午夜福利片| 18在线观看网站| 又粗又硬又长又爽又黄的视频| 免费在线观看视频国产中文字幕亚洲 | 欧美黄色片欧美黄色片| 丝袜美腿诱惑在线| 女警被强在线播放| 国产高清videossex| 一二三四在线观看免费中文在| 国产成人欧美| 久久亚洲精品不卡| 亚洲七黄色美女视频| 亚洲国产av影院在线观看| 啦啦啦视频在线资源免费观看| 国产在线一区二区三区精| 国产熟女欧美一区二区| 亚洲精品国产一区二区精华液| 亚洲精品一卡2卡三卡4卡5卡 | 日韩 欧美 亚洲 中文字幕| 欧美激情 高清一区二区三区| 观看av在线不卡| 精品久久蜜臀av无| 国产片特级美女逼逼视频| 日韩视频在线欧美| 你懂的网址亚洲精品在线观看| 中文精品一卡2卡3卡4更新| 欧美 亚洲 国产 日韩一| www.自偷自拍.com| 亚洲人成77777在线视频| 操美女的视频在线观看| 亚洲国产欧美网| 无遮挡黄片免费观看| 亚洲欧美一区二区三区黑人| 国产极品粉嫩免费观看在线| 大片免费播放器 马上看| 亚洲精品国产av蜜桃| 国产精品一区二区在线观看99| 七月丁香在线播放| 一级毛片我不卡| av在线老鸭窝| 美女大奶头黄色视频| 国产精品一区二区免费欧美 | 悠悠久久av| 亚洲中文av在线| 十八禁人妻一区二区| 日本猛色少妇xxxxx猛交久久| 97精品久久久久久久久久精品| 亚洲av美国av| 国产国语露脸激情在线看| 国精品久久久久久国模美| 亚洲欧美激情在线| 国产精品麻豆人妻色哟哟久久| 欧美黄色片欧美黄色片| 天天躁夜夜躁狠狠久久av| 又粗又硬又长又爽又黄的视频| 99久久99久久久精品蜜桃| 精品国产超薄肉色丝袜足j| 一区二区日韩欧美中文字幕| 亚洲av日韩在线播放| 日韩一卡2卡3卡4卡2021年| 午夜日韩欧美国产| 国产欧美日韩一区二区三区在线| 在线观看免费视频网站a站| 国产免费福利视频在线观看| 99久久精品国产亚洲精品| 亚洲欧美精品综合一区二区三区| 99国产综合亚洲精品| 国产精品一区二区在线不卡| 巨乳人妻的诱惑在线观看| 制服诱惑二区| 国产精品久久久久成人av| 十八禁网站网址无遮挡| 国产人伦9x9x在线观看| 巨乳人妻的诱惑在线观看| 麻豆国产av国片精品| 欧美性长视频在线观看| 色网站视频免费| 欧美大码av| 午夜福利视频在线观看免费| 国产成人一区二区三区免费视频网站 | svipshipincom国产片| 夜夜骑夜夜射夜夜干| 精品视频人人做人人爽| 制服诱惑二区| 50天的宝宝边吃奶边哭怎么回事| 建设人人有责人人尽责人人享有的| 在线观看一区二区三区激情| 国产精品av久久久久免费| 色婷婷av一区二区三区视频| 三上悠亚av全集在线观看| 美女大奶头黄色视频| 久久精品国产综合久久久| kizo精华| 一级片'在线观看视频| 久久久久久久国产电影| 亚洲色图 男人天堂 中文字幕| 国产老妇伦熟女老妇高清| 亚洲,欧美,日韩| 91麻豆av在线| 极品少妇高潮喷水抽搐| 青草久久国产| 一二三四在线观看免费中文在| 精品久久久久久久毛片微露脸 | 9191精品国产免费久久| 成年人免费黄色播放视频| 欧美少妇被猛烈插入视频| 极品少妇高潮喷水抽搐| 一区二区av电影网| 亚洲国产欧美日韩在线播放| 十八禁人妻一区二区| 欧美性长视频在线观看| 99国产综合亚洲精品| 久久久久久久精品精品| 亚洲,欧美精品.| 久久久久国产一级毛片高清牌| 欧美成人精品欧美一级黄| 亚洲av欧美aⅴ国产| 丰满少妇做爰视频| 午夜福利在线免费观看网站| av国产精品久久久久影院| 另类亚洲欧美激情| 亚洲国产中文字幕在线视频| 国产野战对白在线观看| 菩萨蛮人人尽说江南好唐韦庄| 亚洲av综合色区一区| 视频区欧美日本亚洲| 精品一区二区三区四区五区乱码 | 成人午夜精彩视频在线观看| 亚洲国产av新网站| 亚洲专区中文字幕在线| 99香蕉大伊视频| 久久久久网色| 免费观看av网站的网址| 午夜福利影视在线免费观看| 色网站视频免费| 亚洲色图综合在线观看| 国产成人欧美| 亚洲精品日韩在线中文字幕| 无限看片的www在线观看| 99国产精品免费福利视频| 久久久精品94久久精品| 欧美黑人精品巨大| 首页视频小说图片口味搜索 | 欧美日韩亚洲国产一区二区在线观看 | av有码第一页| 啦啦啦中文免费视频观看日本| 久久毛片免费看一区二区三区| 国产精品人妻久久久影院| 18禁观看日本| 别揉我奶头~嗯~啊~动态视频 | 99re6热这里在线精品视频| 久久久精品94久久精品| 考比视频在线观看| 男女边吃奶边做爰视频| 国产xxxxx性猛交| 国产成人免费观看mmmm| 极品少妇高潮喷水抽搐| 亚洲av男天堂| 精品亚洲成a人片在线观看| 99香蕉大伊视频| 性色av一级| 亚洲一区二区三区欧美精品| 欧美日韩av久久| 免费观看人在逋| 午夜福利视频精品| 亚洲成色77777| 亚洲 欧美一区二区三区| 国产一卡二卡三卡精品| 制服人妻中文乱码| 国产午夜精品一二区理论片| 国产片内射在线| 美女中出高潮动态图| 女人精品久久久久毛片| 午夜福利视频在线观看免费| 各种免费的搞黄视频| av天堂在线播放| 在线观看免费午夜福利视频| 欧美日韩亚洲高清精品| 69精品国产乱码久久久| 欧美 日韩 精品 国产| 国产色视频综合| 国产精品国产三级国产专区5o| 电影成人av| 成年人午夜在线观看视频| 黄频高清免费视频| 日本欧美国产在线视频| 永久免费av网站大全| 欧美性长视频在线观看| 高清黄色对白视频在线免费看| 中文字幕精品免费在线观看视频| 国产男女超爽视频在线观看| 美国免费a级毛片| 午夜福利在线免费观看网站| 少妇猛男粗大的猛烈进出视频| www.999成人在线观看| 黑人欧美特级aaaaaa片| 老汉色av国产亚洲站长工具| 啦啦啦在线免费观看视频4| 一二三四社区在线视频社区8| 成人影院久久| 人人妻人人爽人人添夜夜欢视频| 亚洲欧美日韩高清在线视频 | 国产精品久久久av美女十八| 成人三级做爰电影| 亚洲国产看品久久| 黑人猛操日本美女一级片| 精品视频人人做人人爽| 男女高潮啪啪啪动态图| 免费看不卡的av| av一本久久久久| 另类亚洲欧美激情| 久久综合国产亚洲精品| 成年人免费黄色播放视频| 国产99久久九九免费精品| 午夜福利,免费看| 国产精品一区二区免费欧美 | 欧美日韩视频精品一区| 亚洲精品国产av蜜桃| 99热国产这里只有精品6| 女人被躁到高潮嗷嗷叫费观| 操出白浆在线播放| av天堂在线播放| 操美女的视频在线观看| 午夜影院在线不卡| 高清av免费在线| avwww免费| 男女边吃奶边做爰视频| 亚洲第一青青草原| 黄片播放在线免费| 19禁男女啪啪无遮挡网站| 高清不卡的av网站| 波野结衣二区三区在线| 我要看黄色一级片免费的| 欧美xxⅹ黑人| 亚洲欧美精品综合一区二区三区| 只有这里有精品99| 欧美日韩亚洲高清精品| 菩萨蛮人人尽说江南好唐韦庄| 亚洲少妇的诱惑av| 国产一区亚洲一区在线观看| 一区二区日韩欧美中文字幕| 国产一区二区 视频在线| 国产精品二区激情视频| 在线观看免费高清a一片| 国产成人精品在线电影| 亚洲自偷自拍图片 自拍| 性色av乱码一区二区三区2| 超碰成人久久| 久久天躁狠狠躁夜夜2o2o | 免费在线观看黄色视频的| 自线自在国产av| 日本wwww免费看| av天堂在线播放| 成人国语在线视频| 午夜免费观看性视频| 久久久久精品国产欧美久久久 | 国产成人精品在线电影| av不卡在线播放| 大话2 男鬼变身卡| 精品少妇一区二区三区视频日本电影| 亚洲五月婷婷丁香| 亚洲国产成人一精品久久久| 亚洲av男天堂| 日韩精品免费视频一区二区三区| 久久影院123| 一二三四社区在线视频社区8| 中文字幕av电影在线播放| 国产精品一区二区在线观看99| av电影中文网址| 国产福利在线免费观看视频| 日韩一区二区三区影片| 在线观看免费午夜福利视频| 日本91视频免费播放| 国产精品人妻久久久影院| 午夜福利乱码中文字幕| 久久久久久久久免费视频了| 色94色欧美一区二区| 女人爽到高潮嗷嗷叫在线视频| 国产亚洲精品久久久久5区| 手机成人av网站| 丝袜人妻中文字幕| 三上悠亚av全集在线观看| 99久久精品国产亚洲精品| 国产日韩欧美在线精品| 国产人伦9x9x在线观看| 九草在线视频观看| 夜夜骑夜夜射夜夜干| 99国产精品一区二区蜜桃av | 久久久精品国产亚洲av高清涩受| 一个人免费看片子| 黑人猛操日本美女一级片| 国产成人一区二区在线| 国产成人a∨麻豆精品| 中文字幕人妻丝袜制服| 777久久人妻少妇嫩草av网站| 一区二区三区精品91| 18禁观看日本| 婷婷丁香在线五月| 午夜免费男女啪啪视频观看| 丁香六月欧美| 国产成人精品久久久久久| 少妇 在线观看| 欧美精品一区二区免费开放| 午夜免费成人在线视频| 美女扒开内裤让男人捅视频| 人妻一区二区av| 啦啦啦视频在线资源免费观看| 国产有黄有色有爽视频| 久久久久国产精品人妻一区二区| 如日韩欧美国产精品一区二区三区| 国产精品国产av在线观看| 中国国产av一级| 大香蕉久久网| 高清欧美精品videossex| 色精品久久人妻99蜜桃| 久久久久久人人人人人| 久久久久久亚洲精品国产蜜桃av| 日韩欧美一区视频在线观看| 色视频在线一区二区三区| 国产男女超爽视频在线观看| 黄色视频在线播放观看不卡| 国产熟女午夜一区二区三区| 免费黄频网站在线观看国产| av一本久久久久| 19禁男女啪啪无遮挡网站| 亚洲欧美清纯卡通| 国产精品国产三级专区第一集| 久久久久久久国产电影| 又粗又硬又长又爽又黄的视频| 熟女av电影| 日日摸夜夜添夜夜爱| 少妇裸体淫交视频免费看高清 | 每晚都被弄得嗷嗷叫到高潮| 亚洲精品乱久久久久久| 亚洲第一青青草原| 亚洲欧美一区二区三区久久| 亚洲精品自拍成人| www.自偷自拍.com| 国产亚洲午夜精品一区二区久久| 亚洲精品美女久久av网站| 最近手机中文字幕大全| 十八禁人妻一区二区| 嫩草影视91久久| 国产男人的电影天堂91| 精品少妇黑人巨大在线播放| 高清欧美精品videossex| 18在线观看网站| 日韩制服丝袜自拍偷拍| 国产亚洲午夜精品一区二区久久| 大香蕉久久网| 亚洲欧美中文字幕日韩二区| 后天国语完整版免费观看| 人人妻,人人澡人人爽秒播 | 一个人免费看片子| 两个人看的免费小视频| 叶爱在线成人免费视频播放| 久久国产精品影院| bbb黄色大片| 美女中出高潮动态图| 校园人妻丝袜中文字幕| 亚洲欧美精品自产自拍| 性色av一级| 成人国产av品久久久| 欧美精品亚洲一区二区| 精品国产乱码久久久久久小说| 日本欧美视频一区| 人人澡人人妻人| 久久中文字幕一级| 国产亚洲一区二区精品| 欧美黄色淫秽网站| 欧美老熟妇乱子伦牲交| 亚洲欧美一区二区三区久久| 汤姆久久久久久久影院中文字幕| 十八禁高潮呻吟视频| 国产视频首页在线观看| 久久国产精品大桥未久av| 亚洲专区国产一区二区| 国产麻豆69| 国产午夜精品一二区理论片| 好男人视频免费观看在线| 母亲3免费完整高清在线观看| 悠悠久久av| 成人三级做爰电影| 国产欧美亚洲国产| 男人爽女人下面视频在线观看| 日本午夜av视频| 日韩 欧美 亚洲 中文字幕| xxxhd国产人妻xxx| 天天操日日干夜夜撸| 每晚都被弄得嗷嗷叫到高潮| 高清av免费在线| 欧美激情极品国产一区二区三区| a 毛片基地| 亚洲欧美成人综合另类久久久| 人妻 亚洲 视频| 大片免费播放器 马上看| 久久精品成人免费网站| av在线app专区| 男女国产视频网站| 国产精品.久久久| 日韩 亚洲 欧美在线| 亚洲av成人精品一二三区| 亚洲人成网站在线观看播放| 午夜日韩欧美国产| 国产精品国产三级专区第一集| 精品亚洲成a人片在线观看| 两个人看的免费小视频| 亚洲欧美精品自产自拍| 亚洲国产成人一精品久久久| 多毛熟女@视频| 免费高清在线观看日韩| 天堂俺去俺来也www色官网| 一本色道久久久久久精品综合| 在线精品无人区一区二区三| 成年人午夜在线观看视频| 99热国产这里只有精品6| 成人午夜精彩视频在线观看| av国产久精品久网站免费入址| 欧美日韩成人在线一区二区| 亚洲欧美色中文字幕在线| 一区二区三区四区激情视频| 一级片免费观看大全| 晚上一个人看的免费电影| 亚洲欧洲日产国产| 色网站视频免费| 另类精品久久| 日韩av免费高清视频| 90打野战视频偷拍视频| 亚洲中文字幕日韩| 精品欧美一区二区三区在线| 五月开心婷婷网| 秋霞在线观看毛片| 日韩中文字幕视频在线看片| 国产三级黄色录像| 99热全是精品| 大码成人一级视频| 日韩熟女老妇一区二区性免费视频| 亚洲欧美日韩另类电影网站| 亚洲精品乱久久久久久| 岛国毛片在线播放| 国产黄色视频一区二区在线观看| 精品久久久久久电影网| 一本一本久久a久久精品综合妖精| 欧美精品啪啪一区二区三区 | 美女扒开内裤让男人捅视频| 午夜福利在线免费观看网站| 亚洲精品在线美女| 国产激情久久老熟女| 国产av精品麻豆| 蜜桃国产av成人99| 国产精品三级大全| 久久国产精品大桥未久av| 婷婷成人精品国产| 多毛熟女@视频| 香蕉丝袜av| 国产欧美日韩一区二区三区在线| 亚洲精品国产色婷婷电影| 你懂的网址亚洲精品在线观看| 久久久久久久大尺度免费视频| 伦理电影免费视频| 男女免费视频国产| 老熟女久久久| 国产精品一区二区免费欧美 | 国产片内射在线|