數(shù)學(xué)建模獲獎(jiǎng)?wù)撐?嫦娥三號軟著陸軌道設(shè)計(jì)與控制策略

嫦娥三號軟著陸軌道設(shè)計(jì)與控制策略摘要嫦娥三號的核心任務(wù)是實(shí)施高可靠高安全的月面軟著陸,要求著陸器必需具備自主障礙識別與規(guī)避能力。本文將通過研究月球軟著陸主減速段、快速調(diào)整段、避障段和緩速下降段的飛行動力學(xué)模型，在動力學(xué)方程的基礎(chǔ)上對各階段制導(dǎo)律進(jìn)行優(yōu)化設(shè)計(jì)后，結(jié)合工程實(shí)際確立嫦娥三號軟著陸軌道的下降軌跡和控制策略。針對第一問，本文基于開普勒定律求解著陸準(zhǔn)備軌道的運(yùn)動方程，得到近月點(diǎn)與遠(yuǎn)月點(diǎn)的速度。由于著陸點(diǎn)的位置是由近月點(diǎn)的位置和軟著陸軌跡共同決定的。因此本文通過求解嫦娥三號在主減速段和粗避障段的經(jīng)過的橫向距離，利用著陸點(diǎn)的位置近似確定近月點(diǎn)的位置。實(shí)驗(yàn)結(jié)果表明，此近似求解的近月點(diǎn)位置誤差在10%以內(nèi)。同時(shí)，本文所建立的基于飛行器動力學(xué)模型的主減速段燃料最優(yōu)制導(dǎo)律，通過蟻群算法的參數(shù)優(yōu)化，也取得誤差5%以內(nèi)的良好結(jié)果。本文還將證明主減速段的燃料最優(yōu)制導(dǎo)律等效于時(shí)間最優(yōu)制導(dǎo)律。針對第二問，本文通過查閱相關(guān)資料，給出了軟著陸各個(gè)階段的動力學(xué)方程和制導(dǎo)律。主減速段是耗時(shí)最長，燃料消耗最多的階段，因此采用了燃料最優(yōu)制導(dǎo)律?？焖僬{(diào)整階段是過渡階段，擬采用推力線性變化制導(dǎo)律。粗避障段的主要目的是避開大的隕石坑，由于著陸器此時(shí)接近垂直下降，因而采用最優(yōu)開關(guān)制導(dǎo)律。為實(shí)現(xiàn)粗避障，本文采用對光學(xué)圖像進(jìn)行分格，并且判斷各格為障礙區(qū)（0）或非障礙區(qū)（1），最終通過在0-1圖上尋找具有安全半徑的連片非障礙區(qū)，找到合適的著陸區(qū)域，取得了良好的實(shí)驗(yàn)效果。精避障段則要構(gòu)建三維高程圖和平均坡面，計(jì)算平均坡度，估算出障礙高度，最終找出安全著陸區(qū)。緩速下降段主要考慮的是著陸的安全性，本文將利用變推力制導(dǎo)律和垂直軟著陸模型進(jìn)行分析。針對第三問，本文將在模型的分析與檢驗(yàn)中給出相應(yīng)的誤差分析，實(shí)驗(yàn)結(jié)果表明，本文所提出的模型與嫦娥三號實(shí)際結(jié)果相符合，具有實(shí)用性和準(zhǔn)確性。關(guān)鍵詞：制導(dǎo)律蟻群算法數(shù)字高程圖嫦娥三號軟著陸1.問題的重述1.1背景實(shí)施在月球表面的軟著陸是月球勘探計(jì)劃的重要一步。所謂月球軟著陸,是指著陸器在制動系統(tǒng)作用下以很小的速度準(zhǔn)確降落到月面指定區(qū)域,以保證試驗(yàn)設(shè)備和宇航員的安全。由于月球表面沒有大氣，因此整個(gè)軟著陸過程需要制導(dǎo)發(fā)動機(jī)的控制。正確地建立著陸器的飛行動力學(xué)模型，并在此基礎(chǔ)對著陸器的下降軌跡和制導(dǎo)律進(jìn)行優(yōu)化設(shè)計(jì)，對于著陸器飛行程序的設(shè)計(jì)和燃料資源消耗等方面具有重要的指導(dǎo)意義。1.2需要解決的問題根據(jù)題目信息和附件所示內(nèi)容，分析嫦娥三號軟著陸過程中的6個(gè)階段，研究以下的問題：（1）根據(jù)附件1與附件2研究嫦娥三號的著陸準(zhǔn)備軌道，確定該橢圓軌道上近月點(diǎn)與遠(yuǎn)月點(diǎn)的位置，以及嫦娥三號相應(yīng)的速度大小與方向。（2）建立合理的著陸器飛行動力學(xué)方程，根據(jù)軟著陸過程中6個(gè)階段不同的特點(diǎn)與要求，結(jié)合附件1與附件2的內(nèi)容，確立嫦娥三號的下降軌跡和各個(gè)階段的最優(yōu)控制策略。（3）根據(jù)問題一和問題二建立的數(shù)學(xué)模型，對模型進(jìn)行相應(yīng)的誤差分析和敏感性分析。2.模型假設(shè)（1）假設(shè)諸如月球引力非球項(xiàng)、日月引力攝動等影響因素均可忽略不計(jì)（2）月球表面附近沒有大氣，所以在飛行器的動力學(xué)模型中沒有大氣阻力項(xiàng)（3）嫦娥三號在軟著陸全過程不發(fā)生意外碰撞和震蕩3.符號說明符號符號說明嫦娥三號進(jìn)行霍曼轉(zhuǎn)移的初始時(shí)刻嫦娥三號到達(dá)近月點(diǎn)的時(shí)刻嫦娥三號在主減速段的終止時(shí)刻著陸器的月心距月球半徑遠(yuǎn)月點(diǎn)與月面的距離著陸器在近月點(diǎn)的月心距近月點(diǎn)與月面的月心距月球的引力常數(shù)在時(shí)刻的橫向速度在時(shí)刻的徑向速度在時(shí)刻的橫向速度在時(shí)刻的徑向速度橫向速度徑向速度著陸器的月心距著陸器的極角推力加速度發(fā)動機(jī)的推力大小著陸器的初始質(zhì)量推力方向角的系數(shù)，待優(yōu)化參量燃料消耗率燃料消耗質(zhì)量比沖終止時(shí)刻合速度城市和的距離t時(shí)刻螞蟻k由城市轉(zhuǎn)移到城市的概率信息量表示螞蟻k下一步允許走過的城市集合信息量對螞蟻的影響力大小由城市i轉(zhuǎn)移到城市j的期望程度終止時(shí)刻，待優(yōu)化的時(shí)間參量的作用信息增量接近段的水平距離主減速段的水平距離近月點(diǎn)與著陸點(diǎn)的距離粗避障段的初始時(shí)刻粗避障段的終止時(shí)刻月面高度下降速度與垂直方向的夾角，稱為下降角月球的重力加速度制動力的開關(guān)量等效推力加速度的斜率近月點(diǎn)與著陸點(diǎn)的水平距離相對誤差近月點(diǎn)速度相對誤差主減速仰幅角的相對誤差主減速段時(shí)間相對誤差單格內(nèi)障礙點(diǎn)數(shù)單格內(nèi)障礙點(diǎn)閾值過暗灰度臨界值過亮灰度臨界值注：未列出的符號及重復(fù)的符號以出現(xiàn)處為準(zhǔn)。4.問題一的建模與求解4.1問題分析根據(jù)題目信息，著陸準(zhǔn)備軌道為近月點(diǎn)15km，遠(yuǎn)月點(diǎn)100km的橢圓軌道，根據(jù)開普勒定律，可以建立橢圓軌道上的運(yùn)動方程，求解出近月點(diǎn)與遠(yuǎn)月點(diǎn)的速度的大小與方向。著陸點(diǎn)的位置是由近月點(diǎn)在月心坐標(biāo)系中的位置和軟著陸形態(tài)共同決定的。因此可以通過求解軟著陸的下降軌跡，結(jié)合著陸點(diǎn)的位置，確定近月點(diǎn)在月心坐標(biāo)系中的切向距離。此外，利用月球的半徑參數(shù)，可以求得近月點(diǎn)的徑向距離，進(jìn)而確定近月點(diǎn)的位置。相應(yīng)地，遠(yuǎn)月點(diǎn)的位置可以結(jié)合近月點(diǎn)的位置，利用橢圓的軌跡方程進(jìn)行確定。在本問題中，軟著陸的下降軌跡的求解是一個(gè)難點(diǎn)，考慮到著陸器在近月點(diǎn)處離月面距離較高，因此可以在不引起較大的誤差的情況下，建立簡單的模型進(jìn)行求解。根據(jù)相關(guān)參考資料[1]有：在接近段中，著陸器采用接近直線下降方式逐步接近著陸區(qū)，利用附件1與2提供的接近段的高度信息，可以求解出對應(yīng)的水平距離。而主減速段的軌跡可以建立飛行動力學(xué)模型進(jìn)行求解。考慮到主減速段是軟著陸過程用時(shí)最長、推進(jìn)劑消耗最多的任務(wù)段。該段的主要任務(wù)是消除較大的動力下降段初始水平速度(約1.7km/s),因此推進(jìn)劑消耗優(yōu)化是該段制導(dǎo)律的主要設(shè)計(jì)目標(biāo)。因此，求解基于動力學(xué)模型的燃料消耗最低制導(dǎo)律方程，可以確定主減速段的水平距離。綜合上述求解出的距離信息，就可以在一定的誤差允許范圍內(nèi)給出近月點(diǎn)的位置，進(jìn)而得到遠(yuǎn)月點(diǎn)的位置?？傮w而言，問題一的流程圖如下所示：圖4-1 問題一的流程圖4.2模型的建立4.2.1著陸準(zhǔn)備軌道的運(yùn)動學(xué)方程假設(shè)嫦娥三號在時(shí)刻開始進(jìn)行霍曼轉(zhuǎn)移，在時(shí)刻到達(dá)近月點(diǎn)，則在慣性坐標(biāo)系中，以月心為原點(diǎn)的極坐標(biāo)形式的著陸準(zhǔn)備軌道的運(yùn)動學(xué)方程[2]為：在時(shí)刻： （4-1）其中的為著陸器的月心距，為月球半徑，可取值為月球的平均半徑1737.013km，為遠(yuǎn)月點(diǎn)與月面的距離，取值為100km，為著陸器在近月點(diǎn)的月心距，為近月點(diǎn)與月面的距離，取值為15km，為月球的引力常數(shù)，其值為：，為初始橫向速度，即著陸器從遠(yuǎn)月點(diǎn)開始運(yùn)動的速度，為初始徑向速度，其值為0.在時(shí)刻： （4-2）其中的為初始橫向速度，即著陸器從近月點(diǎn)開始運(yùn)動的速度，為初始徑向速度，其值為主減速段的動力學(xué)方程在慣性坐標(biāo)系中，以月心為原點(diǎn)的極坐標(biāo)形式受控飛行器動力學(xué)方程為[2]： （4-3）其中為月球引力常數(shù)，為徑向速度，為橫向速度，為著陸器的月心距，為著陸器的極角。其中為推力加速度，其值可以由下式計(jì)算： （4-4）其中為發(fā)動機(jī)推力的大小，查閱相關(guān)資料[1]可得：主減速段的主要任務(wù)是消除較大的動力下降段初始水平速度，所以發(fā)動機(jī)一直工作在最大推力模式即7500N。為著陸器的初始質(zhì)量，其值為2400kg，為比沖（即單位質(zhì)量的推進(jìn)劑產(chǎn)生的推力），其值為2940m/s。為推力方向角，即推力方向與當(dāng)?shù)厮骄€的夾角，查閱相關(guān)資料[3]可得： （4-5）其中,,,為待優(yōu)化的參量。4.2.3基于蟻群算法[3]的主減速段燃料最優(yōu)制導(dǎo)律燃料最優(yōu)制導(dǎo)律設(shè)主減速段的初始時(shí)刻為，終止時(shí)刻為，對于推力幅值恒定的飛行器，燃料消耗最省的性能指標(biāo)為： （4-6）其中為燃料消耗率，為了使消耗的燃料質(zhì)量最小，由于初始時(shí)刻，發(fā)動機(jī)推力和比沖均為定值，因此只能使主減速段的終止時(shí)刻最小。所以，對于主減速段而言，燃料最優(yōu)制導(dǎo)律等同于時(shí)間最優(yōu)制導(dǎo)律。查閱相關(guān)資料[1]可得，主減速末期的俯仰角為，因此可列出如下等式： （4-7）其中為終止時(shí)刻的合速度，為57m/s.因此，只要確定了的函數(shù)表達(dá)式或函數(shù)圖像，就可以用解析法或作圖法求得。而由式（4-3）所確定，式（4—3）的唯一參量是，又由參量,,,所確定，因此優(yōu)化的問題就轉(zhuǎn)換為優(yōu)化參量,,,的問題。根據(jù)式(4-5)，參量,,,的優(yōu)化可以用多個(gè)離散值和值來多階擬合求得。其中： （4-8）用函數(shù)逼近法進(jìn)行參數(shù)化的相關(guān)參數(shù)設(shè)置為：將軌跡離散化化為9段，一共10個(gè)離散點(diǎn)，相應(yīng)時(shí)刻的值可以由來確定： （4-9）因此待優(yōu)化參數(shù)共11個(gè)，即10個(gè)推力方向角和1個(gè)終端時(shí)刻。對于終端時(shí)刻，根據(jù)齊奧爾科夫斯基公式和軟著陸初始條件，可由下式估計(jì)： （4-10）式中和分別表示著陸器的終端速度和初始速度，經(jīng)計(jì)算確定搜索范圍為(單位秒)： （4-11）我們將用蟻群算法來找到這11個(gè)參數(shù)。蟻群算法[3]設(shè)有n個(gè)城市,(i,j=1,2,…,n)表示城市i和j間的距離,表示在t時(shí)刻城市i和j之間的信息量。我們以此來模擬實(shí)際螞蟻的分泌物。設(shè)共有m只螞蟻,用表示在t時(shí)刻螞蟻k由城市i轉(zhuǎn)移到城市j的概率,則 (4-12)其中,表示螞蟻k下一步允許走過的城市的集合,表示路徑上的信息量對螞蟻選擇路徑所起的作用大小,為由城市i轉(zhuǎn)移到城市j的期望程度,例如,可以?。健１硎镜淖饔?當(dāng)時(shí),算法就是傳統(tǒng)的貪心算法;而當(dāng)時(shí),就成了純粹的正反饋的啟發(fā)式算法.經(jīng)過n個(gè)時(shí)刻,螞蟻可走完所有的城市,完成一次循環(huán).每只螞蟻所走過的路徑就是一個(gè)解.此時(shí),要根據(jù)下面公式對各路徑上的信息量作更新: (4-13)其中,表示信息量隨時(shí)間的推移而衰減的程度.信息增量可表示為 (4-14)表示螞蟻k在本次循環(huán)中在城市i和j之間留下的信息量。下面給出蟻群算法應(yīng)用過程的具體步驟。(1)設(shè)置初始參數(shù)，包括螞蟻數(shù)N，循環(huán)次數(shù)K，揮發(fā)系數(shù)Rho，調(diào)節(jié)系數(shù)Lambda、Q，所有路徑信息素量初值Tau。，螞蟻初始位置等。(2)根據(jù)式(4-8)計(jì)算每只螞蟻的轉(zhuǎn)移概率，然后依據(jù)賭輪原則為每只螞蟻選擇下一個(gè)路徑。重復(fù)上述操作直至所有螞蟻均完成一次循環(huán)。(3)將每只螞蟻的路徑解碼為優(yōu)化參數(shù)值和和1個(gè)終端時(shí)刻，用這10個(gè)數(shù)據(jù)擬合求出，從而得到參量的表達(dá)式，根據(jù)式(4-6)計(jì)算目標(biāo)函數(shù)值求得,對比當(dāng)前最優(yōu)值，判斷找出當(dāng)前最好的路徑。然后根據(jù)式(4-9)、(4-10)更新相關(guān)路徑上的信息素。(4)判斷是否滿足終止條件，不滿足則重復(fù)(1)至(3)步過程；滿足則結(jié)束計(jì)算輸出結(jié)果。本節(jié)模型綜述綜上所述，4.2.3節(jié)的模型可以總結(jié)為下式所述的規(guī)劃問題：目標(biāo)函數(shù)：約束條件： （4-15）4.2.4遠(yuǎn)月點(diǎn)的位置嫦娥三號的運(yùn)行軌道近似一個(gè)橢圓，月球則是該橢圓的焦點(diǎn)之一，若設(shè)橢圓半長軸為a，半短軸為b，焦點(diǎn)為c。則焦點(diǎn)（即月球圓心）距離橢圓軌道（即嫦娥三號的繞月軌道）長軸兩端的距離分別為：。也就是遠(yuǎn)地點(diǎn)和近地點(diǎn)的高度（未減去月球半徑）。根據(jù)橢圓方程： （4-16）可以求得近月點(diǎn)和遠(yuǎn)月點(diǎn)的位置關(guān)系。4.3模型的求解4.3.1近月點(diǎn)和遠(yuǎn)月點(diǎn)速度利用式（4-1）可求得遠(yuǎn)月點(diǎn)的速度為： 利用式（4-2）可求得近月點(diǎn)的速度為：4.3.2近月點(diǎn)的位置接近段的水平距離查閱相關(guān)資料[1]可得，處于接近段的著陸器采用接近直線下降方式逐步接近著陸區(qū)，因此接近段的水平距離可以近似等于接近段的垂直距離。主減速段的水平距離距離可以根據(jù)主減速段最終時(shí)刻的最優(yōu)參量，以及的表達(dá)式積分求得，如下所示： （4-17）根據(jù)蟻群算法的優(yōu)化的優(yōu)化結(jié)果，求得: 從而得到的表達(dá)式為： （4-18）運(yùn)行結(jié)果圖如下：注：時(shí)間大于交點(diǎn)處時(shí)的線段不在主減速階段內(nèi)，故無效圖4-2 的函數(shù)圖像圖4-3 的函數(shù)圖像由圖可得，近月點(diǎn)與著陸點(diǎn)的距離在一定的誤差允許范圍內(nèi)，可以將近月點(diǎn)與著陸點(diǎn)的距離做如下近似計(jì)算：近月點(diǎn)的速度方向速度方向?yàn)?，即主減速段的俯仰角為5.問題二的建模與求解5.1問題分析嫦娥三號軟著陸過程大致可以劃分為6個(gè)階段，著陸準(zhǔn)備階段、主減速段、快速調(diào)整段、粗避障段、精避障段、緩速下降段。問題一中已經(jīng)確立了著陸準(zhǔn)備階段和主減速段的運(yùn)動方程和動力學(xué)模型，主減速段的主要目的是通過推力發(fā)動機(jī)抵消近月點(diǎn)處較大的初始速度，所以燃料消耗最優(yōu)是這一階段的主要優(yōu)化目標(biāo)。問題一中已給出了相應(yīng)的數(shù)學(xué)模型。對于快速調(diào)整段，其介于主減速段和粗避障段之間，主要目的是緩沖作用，即使著陸器平緩地從主減速段過渡到接近段。查閱相關(guān)資料[1]，為了滿足主發(fā)動機(jī)推力和著陸器姿態(tài)勻速過渡的要求,我們利用推力大小和方向線性變化的制導(dǎo)律,其中制導(dǎo)參數(shù)利用實(shí)際的主減速段末端狀態(tài)和接近段初始狀態(tài)約束確定。對于粗避障段，主要任務(wù)是要求避開大的隕石坑，實(shí)現(xiàn)在設(shè)計(jì)著陸點(diǎn)上方100m處懸停，并初步確定落月地點(diǎn)。查閱相關(guān)資料[4],為了滿足制導(dǎo)目標(biāo)的位置、速度、姿態(tài)以及初始高度和速度等多項(xiàng)約束，我們構(gòu)建平面月球模型，并利用最優(yōu)開關(guān)制導(dǎo)律模擬。在粗避障過程中，利用月球巖石坑的圖像特征，采用基于光學(xué)圖像的成像，通過對圖像進(jìn)行分析和處理，判斷出過亮和過暗障礙，識別出當(dāng)前下降區(qū)域上的大障礙，從而確定安全區(qū)域。對于精避障段，主要目的是在粗避障選取的較安全區(qū)域內(nèi)進(jìn)行精確的障礙檢測。由于精避障時(shí)已極接近月球表面，務(wù)必識別并剔除危及安全的小尺度障礙，確保落點(diǎn)安全。為此需要進(jìn)行高分辨率的三維成像，查閱相關(guān)資料，運(yùn)用基于三維圖像的精避障識別和安全著陸區(qū)選取算法處理三維圖像數(shù)據(jù)。對于緩速下降段，主要考慮的是著陸的安全性，查閱相關(guān)資料[5],我們利用變推力制導(dǎo)律和垂直軟著陸模型進(jìn)行分析。問題二的流程圖如下所示： 圖5-1 問題二的流程圖5.2模型的建立5.2.1快速調(diào)整階段為了滿足主發(fā)動機(jī)推力和著陸器姿態(tài)勻速過渡的要求，建立快速調(diào)整階段的推力模型如下所示： （5-1）其中為推力的大小，為推力方向角，即推力方向與當(dāng)?shù)厮骄€的夾角。假設(shè)快速調(diào)整階段的初始時(shí)刻為，終止時(shí)刻為，本階段的動力學(xué)方程可以將式（5-1）、式（4-3）與式（4-4）聯(lián)立，結(jié)合下述初值條件： （5-2）來構(gòu)造。其中,,與之前的定義相同。5.2.2粗避障段最優(yōu)開關(guān)制導(dǎo)律假設(shè)粗避障段的初始時(shí)刻為，終止時(shí)刻為，基于平面月球模型[6]的運(yùn)動方程如下所示： （5-3） 其中為月面高度，為下降速度與垂直方向的夾角,稱為下降角，為月球的重力加速度，為制動力的開關(guān)量，其定義如下所示： （5-4）式（5-4）表明粗避障段的著陸器將先做自由落體運(yùn)動，然后依靠發(fā)動機(jī)制動下降到接近段終端。為了簡化過程，在一定的誤差允許范圍內(nèi)，引入以下兩個(gè)簡化條件：（1）簡化條件1：由于著陸器在本階段接近垂直下降，因此可以將視為小量，因而有： （5-5）（2）簡化條件2：根據(jù)相關(guān)資料[1]，本階段歷時(shí)時(shí)間較短，故引入質(zhì)量不變假設(shè)： （5-6）其中下標(biāo)表示粗避障段初始狀態(tài)，下標(biāo)表示粗避障段終止?fàn)顟B(tài)。最終，我們得到了簡化后的運(yùn)動方程如下所示： （5-7）基于光學(xué)圖像的粗障礙檢測算法通過觀察，月球表面的巖石坑圖像特征如下：1）巖坑表面石頭具有高亮度；2）巖坑存在明顯的陰影區(qū)3）亮區(qū)與陰影之間存在強(qiáng)對比度。這些亮暗區(qū)都是我們剔除的目標(biāo)。（1）光學(xué)圖像的分格與障礙區(qū)的識別將光學(xué)圖像分格為個(gè)等大方格，A題提供的TIF圖為分辨率，將其灰度值分格為，每格像素為。然后對每格的灰度值進(jìn)行統(tǒng)計(jì)。每格共有個(gè)灰度，而灰度過大則表示過亮，即巖坑表面石頭；過小表示巖坑，為暗區(qū)。當(dāng)一個(gè)格的灰度過大和過小的點(diǎn)過多，這個(gè)格將被認(rèn)定為障礙區(qū)，不能降落。因此，需要對灰度過大和過小的點(diǎn)進(jìn)行統(tǒng)計(jì)，統(tǒng)計(jì)方法如下：設(shè)障礙點(diǎn)數(shù)為初值為0，障礙點(diǎn)數(shù)閾值為，過暗灰度臨界值為,過亮灰度臨界值為，對于一格內(nèi)個(gè)灰度點(diǎn)，若： 則表示該灰度點(diǎn)是障礙點(diǎn)，統(tǒng)計(jì)值加一。當(dāng)一格統(tǒng)計(jì)完成時(shí)，若：即障礙點(diǎn)數(shù)超過閾值，此時(shí)判斷該格為障礙區(qū)。為方便起見，標(biāo)記障礙區(qū)所在格為0，非障礙區(qū)為1。（2）采用螺旋搜索算法尋找安全著陸區(qū)通過螺旋搜索算法（圖5-2），從上步得到的的表格中間開始搜索非障礙區(qū)。當(dāng)找到非障礙區(qū)后，為其設(shè)置一個(gè)安全半徑。繼續(xù)搜索，若在安全半徑內(nèi)無障礙區(qū)，則將該半徑內(nèi)區(qū)域設(shè)置為安全著陸區(qū)。圖5-2 螺旋搜索算法示意圖[7]5.2.3精避障段基于三維圖像的精避障識別和安全著陸區(qū)選取[7]設(shè)計(jì)的精障礙識別和安全著陸區(qū)選取算法為：1）數(shù)據(jù)預(yù)處理，構(gòu)建測量坐標(biāo)系下的地形三維高程圖；2）平均坡面構(gòu)建，采用最小二乘法擬合一定單元區(qū)域的平均坡面；3）平均坡度計(jì)算，根據(jù)平均坡面計(jì)算該區(qū)域的平均坡度；4）障礙高度計(jì)算，根據(jù)平均坡面估算該區(qū)域內(nèi)每個(gè)單元格的障礙高度；5）安全著陸區(qū)選取，采用從著陸區(qū)中心開始順時(shí)針陸炫前進(jìn)搜索的方法（圖），直至找到符合著陸要求的區(qū)域位置。圖5-3 中心螺旋式安全著陸區(qū)搜索方法[4]5.2.4緩速下降段垂直軟著陸模型[6]接近段結(jié)束后,著陸器距離月面很近,且下降速度與當(dāng)?shù)厮矫鎺缀醮怪薄２捎梅侄慰刂品绞?著陸器依次經(jīng)過懸停、勻加速、勻減速和關(guān)機(jī)降落幾個(gè)過程,各段加速度分別如下所示： （5-8）等效連續(xù)變推力制導(dǎo)律[6]假設(shè)等效推力的加速度是線性變化的，取向上為正方向，有： （5-9）其中,表示等效推力加速度的斜率,表示推力逐漸增大,表示推力逐漸減??；表示等效推力加速度相對于月心重力加速度；的初始偏移量,表示初始推力加速度大于,表示初始推力加速度小于。結(jié)合關(guān)機(jī)條件： （5-10）對式（5-9）等式兩邊連續(xù)積分兩次，得： （5-11）其中，下標(biāo)和分別表示著陸段初始和關(guān)機(jī)時(shí)刻狀態(tài),下標(biāo)表示著陸段,下標(biāo)表示關(guān)機(jī)時(shí)刻。5.3模型的求解5.3.1粗避障光學(xué)圖像的分格與障礙區(qū)的識別通過MATLAB的imread函數(shù)，讀取粗避障光學(xué)圖像的灰度圖?？紤]到該灰度圖有個(gè)像素，用MATLAB進(jìn)行障礙區(qū)的判斷會相當(dāng)慢，于是將該灰度圖錄入C++下進(jìn)行模型的求解。用C++進(jìn)行分格和識別障礙區(qū)時(shí)，經(jīng)過多次測試設(shè)定各值：閾值過暗灰度臨界值過亮灰度臨界值得到障礙區(qū)的0-1圖如下，通過與原TIF圖比較可驗(yàn)證其正確性:圖5-4 障礙區(qū)0-1圖與原tif圖的比較顯然，0區(qū)對應(yīng)原圖的巖石坑區(qū)，1區(qū)對應(yīng)原圖較平坦的區(qū)域，1區(qū)才是適合著陸的安全區(qū)。計(jì)算結(jié)果和實(shí)際TIF圖表示的地形非常吻合。圓圈所包括的區(qū)域，是最后得出粗避障過程后適合安全降落的區(qū)域。5.3.2精避障三維高程圖圖5-5 精避障三維高程圖6.模型的分析與檢驗(yàn)6.1 著陸準(zhǔn)備段與主減速段模型的分析與檢驗(yàn)根據(jù)相關(guān)資料[5],著陸器動力下降初始狀態(tài)為:經(jīng)度19.0464°W,緯度28.9989°N,速度1695.7m/s,俯仰姿態(tài)角85°，主減速階段經(jīng)歷時(shí)間為487s。根據(jù)經(jīng)緯度距離計(jì)算公式，計(jì)算出實(shí)際近月點(diǎn)與著陸點(diǎn)的水平距離為： 解得：所以計(jì)算結(jié)果誤差為： 其中，，，，依次為近月點(diǎn)與著陸點(diǎn)水平距離相對誤差、近月點(diǎn)速度相對誤差、主減速段仰幅角相對誤差、主減速段時(shí)間相對誤差。結(jié)論：由上述數(shù)據(jù)可知，與后三項(xiàng)的值相比，距離誤差較大。從近月點(diǎn)與著陸點(diǎn)距離計(jì)算過程中可以知道，為了簡化計(jì)算，我們只考慮了主減速段和粗避障段產(chǎn)生的水平距離，而忽略了其他階段產(chǎn)生的水平距離，因而誤差相對較大，但是仍然控制在10%以內(nèi)，說明在一定誤差允許范圍內(nèi)，此結(jié)果具有一定的參考性。而后三項(xiàng)的相對誤差比較小，說明他們都與嫦娥三號實(shí)際飛行得到的結(jié)果比較符合。有力地論證了著陸準(zhǔn)備階段和主減速段模型的準(zhǔn)確性與可靠性。6.2粗避障光學(xué)圖像分格與障礙區(qū)識別的分析與檢驗(yàn)粗避障的模型通過模型的求解檢驗(yàn)了其正確性，在著陸器飛行時(shí)，若著陸器能及時(shí)計(jì)算出安全區(qū)域，該模型還是具有一定可行性的。該粗避障模型計(jì)算量偏大，對著陸器計(jì)算能力無疑是一大考驗(yàn)。該粗避障模型計(jì)算結(jié)果與實(shí)際地形誤差較少，若要進(jìn)一步減少誤差，提高準(zhǔn)確性，則需提高分格的格數(shù)，同時(shí)還要再度優(yōu)化過暗灰度臨界值、過亮灰度臨界值和閾值。7.模型的創(chuàng)新點(diǎn)與不足之處7.1模型的創(chuàng)新點(diǎn)（1）將月球軟著陸的軌跡優(yōu)化問題轉(zhuǎn)化為參數(shù)優(yōu)化問題，將求燃料最優(yōu)問題轉(zhuǎn)化為求時(shí)間最優(yōu)問題，使問題的求解變得簡單。（2）使用蟻群算法求解復(fù)雜的角度變化參數(shù)優(yōu)化問題，具有很好的試驗(yàn)效果。（3）利用光學(xué)圖像的分格和設(shè)定閾值來尋找連片非障礙區(qū)作為安全著陸區(qū)域。7.2模型的不足之處蟻群算法比較適合處理離散的優(yōu)化問題，從而限制了其在連續(xù)優(yōu)化問題中的應(yīng)用，本模型中將角度變化離散化，存在較大的誤差。8.參考資料[1]張洪華,關(guān)軼峰,黃翔宇,等.嫦娥三號著陸器動力下降的制導(dǎo)導(dǎo)航與控制.中國科學(xué):技術(shù)科學(xué),2014,44:377–384[2]王劼，李俊峰，崔乃剛，劉暾，登月飛行器軟著陸發(fā)動機(jī)持續(xù)工作初始點(diǎn)選擇研究，工程力學(xué)，1000-4750，(2003)06-0145-04[3]段佳佳，徐世杰，朱建豐，基于蟻群算法的月球軟著陸軌跡優(yōu)化，宇航學(xué)報(bào)，第29卷第2期：1000-1328，(2008)02.0476.06[4]KlumppAR.Apollolunardescentguidance.Automatica,1974,10:133–146[5]張洪華,梁俊,黃翔宇,等.嫦娥三號自主避障軟著陸控制技術(shù).中國科學(xué):技術(shù)科學(xué),2014,44:559–568[6]王鵬基，張熇，曲廣吉，月球軟著陸飛行動力學(xué)和制導(dǎo)控制建模與仿真，中國科學(xué)E輯:技術(shù)科學(xué)2009年，第39卷第3期:521~527[7]張洪華,梁俊,黃翔宇,等.嫦娥三號自主避障軟著陸控制技術(shù).中國科學(xué):技術(shù)科學(xué),2014,44:559–5689.附錄MATLAB源程序function[BESTX,BESTY,ALLX,ALLY]=ACOUCP(K,N,Rho,Q,Lambda,LB,UB)%%ACOUCP.m%在命令行中輸入LB=[00000000000];UB=[11111111111];ACOUCP(50,30,0.95,1,0.5,LB,UB)%%輸入?yún)?shù)列表%K迭代次數(shù)%N蟻群規(guī)模%Rho信息素蒸發(fā)系數(shù)，取值0～1之間，推薦取值0.7～0.95%Q信息素增加強(qiáng)度，大于0，推薦取值1左右%Lambda螞蟻爬行速度，取值0～1之間，推薦取值0.1～0.5%LB決策變量的下界，M×1的向量%UB決策變量的上界，M×1的向量%%輸出參數(shù)列表%BESTXK×1細(xì)胞結(jié)構(gòu)，每一個(gè)元素是M×1向量，記錄每一代的最優(yōu)螞蟻%BESTYK×1矩陣，記錄每一代的最優(yōu)螞蟻的評價(jià)函數(shù)值%ALLXK×1細(xì)胞結(jié)構(gòu)，每一個(gè)元素是M×N矩陣，記錄每一代螞蟻的位置%ALLYK×N矩陣，記錄每一代螞蟻的評價(jià)函數(shù)值%%參考設(shè)置%[BESTX,BESTY,ALLX,ALLY]=ACOUCP(50,30,0.95,1,0.5,LB,UB)%%第一步：初始化M=length(LB);%決策變量的個(gè)數(shù)%蟻群位置初始化X=zeros(M,N);fori=1:Mx=unifrnd(LB(i),UB(i),1,N);%unidrnd()離散分布X(i,:)=x;end%輸出變量初始化ALLX=cell(K,1);%細(xì)胞結(jié)構(gòu)，每一個(gè)元素是M×N矩陣，記錄每一代的個(gè)體ALLY=zeros(K,N);%K×N矩陣，記錄每一代評價(jià)函數(shù)值BESTX=cell(K,1);%細(xì)胞結(jié)構(gòu)，每一個(gè)元素是M×1向量，記錄每一代的最優(yōu)個(gè)體BESTY=zeros(K,1);%K×1矩陣，記錄每一代的最優(yōu)個(gè)體的評價(jià)函數(shù)值k=1;%迭代計(jì)數(shù)器初始化Tau=ones(1,N);%信息素初始化Y=zeros(1,N);%適應(yīng)值初始化%%第二步：迭代過程whilek<=KYY=zeros(1,N);forn=1:Nx=X(:,n);YY(n)=FIT(x);endmaxYY=max(YY);temppos=find(YY==maxYY);POS=temppos(1);%螞蟻隨機(jī)探路forn=1:Nifn~=POSx=X(:,n);Fx=FIT(x);mx=GaussMutation(x);Fmx=FIT(mx);ifFmx<FxX(:,n)=mx;Y(n)=Fmx;elseifrand>1-(1/(sqrt(k)))X(:,n)=mx;Y(n)=Fmx;elseX(:,n)=x;Y(n)=Fx;endendendforn=1:Nifn~=POSx=X(:,n);Fx=FIT(x);mx=GaussMutation(x);Fmx=FIT(mx);ifFmx<FxY(n)=Fmx;elseifrand>1-(1/(sqrt(k)))X(:,n)=mx;Y(n)=Fmx;elseX(:,n)=x;Y(n)=Fx;endendend%朝信息素最大的地方移動forn=1:Nifn~=POSx=X(:,n);r=(K+k)/(K+K);p=randperm(N);t=ceil(r*N);pos=p(1:t);TempTau=Tau(pos);maxTempTau=max(TempTau);pos2=find(TempTau==maxTempTau);pos3=pos(pos2(1));x2=X(:,pos3(1));x3=(1-Lambda)*x+Lambda*x2;Fx=FIT(x);Fx3=FIT(mx);ifFx3<FxX(:,n)=x3;Y(n)=Fx3;elseifrand>1-(1/(sqrt(k)))X(:,n)=x3;Y(n)=Fx3;elseX(:,n)=x;Y(n)=Fx;endendend%更新信息素并記錄Tau=Tau*(1-Rho);maxY=max(Y);minY=min(Y);DeltaTau=(maxY-Y)/(maxY-minY);Tau=Tau+Q*DeltaTau;ALLX{k}=X;ALLY(k,:)=Y;minY=min(Y);pos4=find(Y==minY);BESTX{k}=X(:,pos4(1));BESTY(k)=minY;disp(k);k=k+1;end%從蟻群算法中得到：tf和十個(gè)角度值angle(單位：度)%計(jì)算各個(gè)時(shí)刻functiontf1=angle(M,tf)fori=1:10td(i)=0+(i-1)*(tf-0)/9;end%三次性項(xiàng)擬合globalA;A=polyfit(td,M,3)aL=1737.013*10^3;%月球半徑%求解飛行器動力學(xué)方程[T,Y]=ode45('dip',[0750],[148.165,aL+15*10^3,1693.53,0]);%A=[-6.1321*10^(-8)8.0647*10^(-5)0.01635]%figure%plot(T,Y(:,3),'b-',0:750,51.66,'r-')%從數(shù)值解上選取坐標(biāo)(tf,51.66)tf1=interp1(Y(:,3),T,51.66,'nearest','extrap');%[tf1,Vr]=ginput;functiondy=dip(t,y)u=0.4902802627*10^13;%月球引力常數(shù)globalA;%飛行器動力學(xué)方程%y1是徑向速度，y2是月心距，y3是切向速度,y4是橫向距離%[Vr,Vo,r,Q]=dsolve('DVr=-u/(r^2)+(Vo^2)/r+7500/(2400-7500*t/2940)*sin(a0+a1*t+a2*t^2+a3*t^3)','DVo=-Vr*Vo/r+7500/(2400-7500*t/2940)*cos(a0+a1*t+a2*t^2+a3*t^3)','Dr=Vr','DQ=Vo/r','(Vr(0))^2+(Vo(0))^2=2890000','r(0)=aL+15*10^3','t');dy=zeros(4,1);dy(1)=-u/(y(2))^2+(y(3))^2/y(2)-7500/(2400-7500*t/2940)*sin((A(4)+A(3)*t+A(2)*t^2+A(1)*t^3)*pi/180);dy(2)=y(1);dy(3)=-y(1)*y(3)/y(2)-7500/(2400-7500*t/2940)*cos((A(4)+A(3)*t+A(2)*t^2+A(1)*t^3)*pi/180);dy(4)=y(3);end%使優(yōu)化參數(shù)反歸一化，代入目標(biāo)檢測函數(shù)中functionk=FIT(x)M=[x(1),x(2),x(3),x(4),x(5),x(6),x(7),x(8),x(9),x(10)];M=90*M;M=sort(M);tf2=angle(M,300+x(11)*200);k=tf2;%高斯變異functionmx=GaussMutation(x)l=normrnd(0,1,11,1);k=1/max(l);mx=x.*(1+k*l);mx=mx-floor(mx);aL=1737.013*10^3;%月球半徑globalA;%A=[-6.1321*10^(-8)8.0647*10^(-5)0.01635];A=[-6.4716*10^(-7)5.7154*10^-40.00582.2066][T,Y]=ode45('dip',[0750],[148.165,aL+15*10^3,1693.53,0]);figureplot(T,Y(:,3),'b-',0:1:750,51.66,'-');holdonfigureplot(T,Y(:,4),'b-');C++源程序：//粗避障分格、設(shè)置障礙區(qū)（0）和非障礙區(qū)（1）//先用MATLAB讀取TIF圖的灰度矩陣保存至TXT文件，再用C++通過讀取TXT從而讀取灰度矩陣#include<iostream>#include<fstream>usingnamespacestd;intmain(){ifstreamfile1("d:\\b.txt",ios::in);short*x;x=newshort[2300*2300];charc;for(inti=0;i<2300*2300;i++){file1>>x[i]>>c;}file1.close();inty[23*23];intthreshold=700,sum;shortDark=84;shortLight=105;intt0,t1;for(intxvec=0;xvec<23;xvec++){for(intyvec=0;yvec<23;yvec++){sum=0;t0=2300*100*xvec; t1=100*yvec; for(inti=0;i<100;i++){for(intj=0;j<100;j++) {if((x[t0+2300*i+t1+j]<Dark)||(x[t0+2300*i+t1+j]>Light)) sum++;}}if(sum<threshold)y[xvec*23+yvec]=1; elsey[xvec*23+yvec]=0;}}ofstreamfile2("d:\\testc.txt",ios::out);for(intxvec=0;xvec<23;xvec++){for(intyvec=0;yvec<23;yvec++){file2<<y[xvec*23+yvec]<<'';}file2<<'\n';}file2.close();}============================================================================基于C8051F單片機(jī)直流電動機(jī)反饋控制系統(tǒng)的設(shè)計(jì)與研究基于單片機(jī)的嵌入式Web服務(wù)器的研究MOTOROLA單片機(jī)MC68HC（8）05PV8/A內(nèi)嵌EEPROM的工藝和制程方法及對良率的影響研究基于模糊控制的電阻釬焊單片機(jī)溫度控制系統(tǒng)的研制基于MCS-51系列單片機(jī)的通用控制模塊的研究基于單片機(jī)實(shí)現(xiàn)的供暖系統(tǒng)最佳啟停自校正（STR）調(diào)節(jié)器單片機(jī)控制的二級倒立擺系統(tǒng)的研究基于增強(qiáng)型51系列單片機(jī)的TCP/IP協(xié)議棧的實(shí)現(xiàn)基于單片機(jī)的蓄電池自動監(jiān)測系統(tǒng)基于32位嵌入式單片機(jī)系統(tǒng)的圖像采集與處理技術(shù)的研究基于單片機(jī)的作物營養(yǎng)診斷專家系統(tǒng)的研究基于單片機(jī)的交流伺服電機(jī)運(yùn)動控制系統(tǒng)研究與開發(fā)基于單片機(jī)的泵管內(nèi)壁硬度測試儀的研制基于單片機(jī)的自動找平控制系統(tǒng)研究基于C8051F040單片機(jī)的嵌入式系統(tǒng)開發(fā)基于單片機(jī)的液壓動力系統(tǒng)狀態(tài)監(jiān)測儀開發(fā)模糊Smith智能控制方法的研究及其單片機(jī)實(shí)現(xiàn)一種基于單片機(jī)的軸快流CO〈,2〉激光器的手持控制面板的研制基于雙單片機(jī)沖床數(shù)控系統(tǒng)的研究基于CYGNAL單片機(jī)的在線間歇式濁度儀的研制基于單片機(jī)的噴油泵試驗(yàn)臺控制器的研制基于單片機(jī)的軟起動器的研究和設(shè)計(jì)基于單片機(jī)控制的高速快走絲電火花線切割機(jī)床短循環(huán)走絲方式研究基于單片機(jī)的機(jī)電產(chǎn)品控制系統(tǒng)開發(fā)基于PIC單片機(jī)的智能手機(jī)充電器基于單片機(jī)的實(shí)時(shí)內(nèi)核設(shè)計(jì)及其應(yīng)用研究基于單片機(jī)的遠(yuǎn)程抄表系統(tǒng)的設(shè)計(jì)與研究基于單片機(jī)的煙氣二氧化硫濃度檢測儀的研制基于微型光譜儀的單片機(jī)系統(tǒng)單片機(jī)系統(tǒng)軟件構(gòu)件開發(fā)的技術(shù)研究基于單片機(jī)的液體點(diǎn)滴速度自動檢測儀的研制基于單片機(jī)系統(tǒng)的多功能溫度測量儀的研制基于PIC單片機(jī)的電能采集終端的設(shè)計(jì)和應(yīng)用基于單片機(jī)的光纖光柵解調(diào)儀的研制氣壓式線性摩擦焊機(jī)單片機(jī)控制系統(tǒng)的研制基于單片機(jī)的數(shù)字磁通門傳感器基于單片機(jī)的旋轉(zhuǎn)變壓器-數(shù)字轉(zhuǎn)換器的研究基于單片機(jī)的光纖Bragg光柵解調(diào)系統(tǒng)的研究單片機(jī)控制的便攜式多功能乳腺治療儀的研制基于C8051F020單片機(jī)的多生理信號檢測儀基于單片機(jī)的電機(jī)運(yùn)動控制系統(tǒng)設(shè)計(jì)Pico專用單片機(jī)核的可測性設(shè)計(jì)研究基于MCS-51單片機(jī)的熱量計(jì)基于雙單片機(jī)的智能遙測微型氣象站MCS-51單片機(jī)構(gòu)建機(jī)器人的實(shí)踐研究基于單片機(jī)的輪軌力檢測基于單片機(jī)的GPS定位儀的研究與實(shí)現(xiàn)基于單片機(jī)的電液伺服控制系統(tǒng)用于單片機(jī)系統(tǒng)的MMC卡文件系統(tǒng)研制基于單片機(jī)的時(shí)控和計(jì)數(shù)系統(tǒng)性能優(yōu)化的研究基于單片機(jī)和CPLD的粗光柵位移測量系統(tǒng)研究單片機(jī)控制的后備式方波UPS提升高職學(xué)生單片機(jī)應(yīng)用能力的探究基于單片機(jī)控制的自動低頻減載裝置研究基于單片機(jī)控制的水下焊接電源的研究基于單片機(jī)的多通道數(shù)據(jù)采集系統(tǒng)基于uPSD3234單片機(jī)的氚表面污染測量儀的研制基于單片機(jī)的紅外測油儀的研究96系列單片機(jī)仿真器研究與設(shè)計(jì)基于單片機(jī)的單晶金剛石刀具刃磨設(shè)備的數(shù)控改造基于單片機(jī)的溫度智能控制系統(tǒng)的設(shè)計(jì)與實(shí)現(xiàn)基于MSP430單片機(jī)的電梯門機(jī)控制器的研制基于單片機(jī)的氣體測漏儀的研究基于三菱M16C/6N系列單片機(jī)的CAN/USB協(xié)議轉(zhuǎn)換器基于單片機(jī)和DSP的變壓器油色譜在線監(jiān)測技術(shù)研究基于單片機(jī)的膛壁溫度報(bào)警系統(tǒng)設(shè)計(jì)基于AVR單片機(jī)的低壓無功補(bǔ)償控制器的設(shè)計(jì)基于單片機(jī)船舶電力推進(jìn)電機(jī)監(jiān)測系統(tǒng)基于單片機(jī)網(wǎng)絡(luò)的振動信號的采集系統(tǒng)基于單片機(jī)的大容量數(shù)據(jù)存儲技術(shù)的應(yīng)用研究基于單片機(jī)的疊圖機(jī)研究與教學(xué)方法實(shí)踐基于單片機(jī)嵌入式Web服務(wù)器技術(shù)的研究及實(shí)現(xiàn)基于AT89S52單片機(jī)的通用數(shù)據(jù)采集系統(tǒng)基于單片機(jī)的多道脈沖幅度分析儀研究機(jī)器人旋轉(zhuǎn)電弧傳感角焊縫跟蹤單片機(jī)控制系統(tǒng)基于單片機(jī)的控制系統(tǒng)在PLC虛擬教學(xué)實(shí)驗(yàn)中的應(yīng)用研究基于單片機(jī)系統(tǒng)的網(wǎng)絡(luò)通信研究與應(yīng)用基于PIC16F877單片機(jī)的莫爾斯碼自動譯碼系統(tǒng)設(shè)計(jì)與研究基于單片機(jī)的模糊控制器在工業(yè)電阻爐上的應(yīng)用研究基于雙單片機(jī)沖床數(shù)控系統(tǒng)的研究與開發(fā)基于Cygnal單片機(jī)的μC/OS-Ⅱ的研究基于單片機(jī)的一體化智能差示掃描量熱儀系統(tǒng)研究基于TCP/IP協(xié)議的單片機(jī)與Internet互聯(lián)的研究與實(shí)現(xiàn)變頻調(diào)速液壓電梯單片機(jī)控制器的研究基