技術(shù)特征：

1.一種小天體附著探測下降軌跡優(yōu)化方法，其特征在于：采用內(nèi)球諧引力場模型估計目標(biāo)小天體附近的引力加速度，采用凸優(yōu)化算法解軌跡優(yōu)化方程，最終實現(xiàn)下降軌跡優(yōu)化。

2.如權(quán)利要求1所述的一種小天體附著探測下降軌跡優(yōu)化方法，其特征在于：所述凸優(yōu)化算法包括松弛、線性化、離散化和內(nèi)點法。

3.如權(quán)利要求1所述的一種小天體附著探測下降軌跡優(yōu)化方法，其特征在于：具體詳細步驟如下：

步驟一、由最小二乘法估計內(nèi)球諧引力場模型的球諧系數(shù)：

在小天體外部構(gòu)造內(nèi)布里淵球，所構(gòu)造內(nèi)布里淵球應(yīng)該與小天體表面的目標(biāo)著陸點相切，球心選取應(yīng)保證著陸器下降軌跡包含在內(nèi)布里淵球內(nèi)部；在內(nèi)布里淵球內(nèi)部任取N_data個點，由多面體模型計算N_data個點的引力加速度，并由所述引力加速度構(gòu)造3N_data×1維矩陣由最小二乘法求取內(nèi)球諧引力模型球諧系數(shù)；

$\[{\tilde{C}}_{n,m}^i\]=\left({\[Q_{n,m}^i\]}^{T}W{\[Q_{n,m}^i\]}^{-1}\right)\left\{{\[Q_{n,m}^i\]}^{T}W\left(\frac{\∂V}{\∂r}\right)\right\}---\left(1\right)$

其中，n和m分別表示勒讓德多項式的階次和冪次，為一個(n²+2n)×1維的矩陣，包含了所要求取的各階內(nèi)球諧系數(shù)

$\[{\tilde{C}}_{n,m}^i\]=\[C_{1,0}^i,C_{1,1}^i,S_{1,1}^i,C_{2,0}^i,C_{2,1}^i,S_{2,1}^i,C_{2,2}^i,S_{2,2}^i,C_{3,0}^i,...\]---\left(2\right)$

為引力加速度關(guān)于內(nèi)球諧系數(shù)的導(dǎo)數(shù)構(gòu)成的矩陣，維度為(n²+2n)×3N_data；W為一個3N_data維單位陣；

步驟二、構(gòu)造小天體附著探測下降軌跡優(yōu)化方程：

在小天體固連坐標(biāo)系下，著陸器滿足如下的動力學(xué)方程

$\left\{\begin{array}{c}\stackrel{\·}{r}=v\\ \stackrel{\·}{v}=-2\mathrm{\ω}\×v-\mathrm{\ω}\×(\mathrm{\ω}\×r)+\▿V+\frac{T}{m}\\ \stackrel{\·}{m}=-\frac{\left|\right|T\left|\right|}{I_{sp}{g}_e}\end{array}\right.---\left(3\right)$

其中，r＝[x,y,z]^T和分別為固連坐標(biāo)系下著陸器的位置矢量和速度矢量；ω＝[0,0,ω]^T為小天體自旋角速度矢量；T＝[T_x,T_y,T_z]^T為著陸器推力矢量；m_e為著陸器的質(zhì)量；I_sp為發(fā)動機比沖；g_e＝9.807為地球引力加速度常數(shù)；為引力加速度矢量，由以下公式計算得到：

$\begin{array}{c}\frac{\∂V}{\∂x}=\frac{GM}{{\left(R\right)}^2}\underset{n=1}{\overset{\mathrm{\∞}}{\Σ}}\underset{m=0}{\overset{n}{\Σ}}\left[\begin{array}{c}{C}_{n,m}^i\\ S_{n,m}^i\end{array}\right]\·\[-\frac{1}{2}(1+{\mathrm{\δ}}_{0,m})b_{n-1,m+1}^i+\frac{1}{2}\frac{(n+m)!}{(n+m-2)!}{b}_{n-1,m-1}^i\]\\ \frac{\∂V}{\∂y}=\frac{GM}{{\left(R\right)}^2}\underset{n=1}{\overset{\mathrm{\∞}}{\Σ}}\underset{m=0}{\overset{n}{\Σ}}\left[\begin{array}{c}{S}_{n,m}^i\\ -C_{n,m}^i\end{array}\right]\·\[\frac{1}{2}(1+{\mathrm{\δ}}_{0,m})b_{n-1,m+1}^i+\frac{1}{2}\frac{(n+m)!}{(n+m-2)!}{b}_{n-1,m-1}^i\]\\ \frac{\∂V}{\∂z}=\frac{GM}{{\left(R\right)}^2}\underset{n=1}{\overset{\mathrm{\∞}}{\Σ}}\underset{m=0}{\overset{n}{\Σ}}\left[\begin{array}{c}{C}_{n,m}^i\\ S_{n,m}^i\end{array}\right]\·\[(n+m)b_{n-1,m}^i\]\end{array}---\left(16\right)$

其中，G為萬有引力常數(shù)，M為目標(biāo)小天體質(zhì)量，R為步驟一中內(nèi)布里淵球半徑，δ_0,m為Kronecker delta(克勞內(nèi)克)函數(shù)(當(dāng)m＝0時，δ_0,m＝1)；為步驟一所求的內(nèi)球諧系數(shù)；為步驟一中內(nèi)球諧引力場模型基函數(shù)，滿足如下遞推關(guān)系；

$\left\{\begin{array}{c}{b}_{0,0}^i=\left[\begin{array}{c}1\\ 0\end{array}\right]\\ b_{n,n}^i=(2n-1)\left[\begin{array}{cc}x/R& -y/R\\ y/R& x/R\end{array}\right]b_{n-1,n-1}^i\\ b_{n,m}^i=\frac{2n-1}{n-m}\frac{z}{R^i}{b}_{n-1,m}^i-\frac{n+m-1}{n-m}{\left(\frac{\sqrt{x^2+y^2+z^2}}{R}\right)}^2{b}_{n-2,m}^i\\ b_{n,n-1}^i=(2n-1)\frac{z}{R}{b}_{n-1,n-1}^i\end{array}\right.---\left(5\right)$

著陸器滿足如下的邊界條件

$\begin{array}{c}r\left(0\right)=r_0,v\left(t_0\right)=v_0,m_e\left(0\right)=m_0\\ r\left(t_f\right)=r_f,v\left(t_f\right)=v_f\end{array}---\left(6\right)$

其中，r₀，v₀和m₀分別為起始位置處著陸器的位置、速度和質(zhì)量；t_f為終端時間，r_f為目標(biāo)著陸點，v_f為終端速度約束，對于軟著陸問題，v_f＝0；

著陸器推力矢量滿足如下約束條件

0≤||T||≤T_max (7)

燃料最優(yōu)性能指標(biāo)表示如下

$J={\∫}_0^{t_f}\left|\right|T\left|\right|dt---\left(8\right)$

公式(3)、(6)、(7)、(8)構(gòu)成了下降軌跡優(yōu)化方程；

步驟三、向步驟二所得的下降軌跡優(yōu)化方程中引入松弛變量Γ，對所述下降軌跡優(yōu)化方程進行松弛：

引入松弛變量Γ替代下降軌跡優(yōu)化方程中的||T||，則松弛后的軌跡優(yōu)化方程為：

$\begin{array}{c}\begin{array}{cc}\min & {\∫}_0^{t_f}\mathrm{\Γ}dt\end{array}\\ \begin{array}{ccc}subjecto& \stackrel{\·}{m}=-\frac{\mathrm{\Γ}}{I_{sp}{g}_e}& \stackrel{\·}{r}=v\end{array}\\ \stackrel{\·}{v}=-2\mathrm{\ω}\×v-\mathrm{\ω}\×(\mathrm{\ω}\×r)+\▿V+\frac{T}{m}\\ \begin{array}{ccc}r\left(0\right)=r_0& v\left(0\right)=v_0& m\left(0\right)=m_0\end{array}\\ \begin{array}{cc}r\left(t_f\right)=r_f& v\left(t_f\right)=0\end{array}\\ \begin{array}{cc}0\≤\mathrm{\Γ}\≤T_{\max }& \sqrt{T_x^2+T_y^2+T_z^2}\≤\mathrm{\Γ}\end{array}\end{array}---\left(9\right)$

步驟四、對步驟三所得方程線性化處理：

定義新變量如下

$\begin{array}{ccc}\mathrm{\σ}\stackrel{\mathrm{\Δ}}{=}\frac{\mathrm{\Γ}}{m}& u\stackrel{\mathrm{\Δ}}{=}\frac{T}{m}& p\stackrel{\mathrm{\Δ}}{=}lnm\end{array}---\left(10\right)$

將以上變量代入步驟三的優(yōu)化方程中，得到線性化的優(yōu)化方程為

$\begin{array}{c}\begin{array}{cc}\min & J={\∫}_0^{t_f}\mathrm{\σ}dt\end{array}\\ \begin{array}{ccc}subjecto& \stackrel{\·}{p}=-\frac{\mathrm{\σ}}{I_{sp}{g}_e}& \stackrel{\·}{r}=v\end{array}\\ \stackrel{\·}{v}=-2\mathrm{\ω}\×v-\mathrm{\ω}\×(\mathrm{\ω}\×r)+\▿V+u\\ \begin{array}{ccc}r\left(0\right)=r_0& v\left(0\right)=v_0& p\left(0\right)=\ln m_0\end{array}\\ \begin{array}{cc}r\left(t_f\right)=r_f& v\left(t_f\right)=0\end{array}\\ \begin{array}{cc}\left|\right|u\left|\right|\≤\mathrm{\σ}& 0\≤\mathrm{\σ}\≤T_{\max }{e}^{-p_0}\[1-(p-p_0)\]\end{array}\end{array}---\left(11\right)$

其中，t表示時間變量；

步驟五、對步驟四所得的優(yōu)化方程進行離散化處理：

將時間區(qū)間[0,t_f]等分成N份，得到對步驟四的優(yōu)化方程進行離散化處理，經(jīng)過離散化以后，軌跡優(yōu)化方程轉(zhuǎn)換為參數(shù)優(yōu)化方程，參數(shù)優(yōu)化方程的表達式如下：

$\begin{array}{c}\min \underset{u_k,{\mathrm{\σ}}_k}{J}=\underset{k=0}{\overset{N}{\Σ}}{\mathrm{\σ}}_k\\ subjectto\left\{\begin{array}{c}{X}_0={\[r_0^T,v_0^T,p_0\]}^T\\ {\mathrm{MX}}_N={\[r_f^T,0_{1\×3}\]}^T\\ \begin{array}{cc}\left|\right|u_k\left|\right|\≤{\mathrm{\σ}}_k& k=0,\mathrm{...},N\end{array}\\ \begin{array}{cc}0\≤{\mathrm{\σ}}_k\≤T_{\max }{e}^{-p_0\left(t_k\right)}\[1-(p_k-p_0(t_k\left)\right)\]& k=0,\mathrm{...},N\end{array}\end{array}\right.\end{array}---\left(12\right)$

其中，M＝[I₆,0_6×1]，t_k表示第k個時間節(jié)點，u_k和σ_k分別表示t_k時刻u和σ的取值；

步驟六、采用內(nèi)點法對步驟五中的參數(shù)優(yōu)化方程進行迭代求解，具體求解過程如下：

1)令引力加速度為一常值▽V₀，采用內(nèi)點法求解步驟五中的參數(shù)優(yōu)化方程，得到一條軌跡；

2)將步驟1)所得到的軌跡作為參考軌跡，計算參考軌跡中每個節(jié)點(即k＝0,…,N)處的引力加速度，并帶入步驟五中的參數(shù)優(yōu)化方程中，采用內(nèi)點法進行求解，得到一條新的軌跡，將得到的新軌跡作為下一次迭代的參考軌跡；

3)當(dāng)?shù)玫降能壽E收斂，則得到最優(yōu)解，即得到探測最優(yōu)下降軌跡。

4.如權(quán)利要求3所述的一種小天體附著探測下降軌跡優(yōu)化方法，其特征在于：步驟五所述離散化處理方法采用顯式四階龍格庫塔積分公式法。