上一篇blog翻譯了 Lars Blackmore(Lars Blackmore is principal rocket landing engineer at SpaceX)的文章,SpaceX 使用 CVXGEN 生成定制飛行代碼,實現超高速機載凸優化。利用地形相對導航實現了數十米量級的導航精度,著陸器在著陸過程中成像行星表面并將特征與機載地圖匹配來確定自己的位置。
今天來講解下Lars Blackmore 在MIT讀博時候的成果:Non-Convex Soft Landing Optimal Control
本文禁止轉載,如果對技術解讀有不同見解,可以郵箱反饋:fanzexuan135@163.com
非凸軟著陸最優控制問題的控制邊界和指向約束的無損凸化
1. 介紹
行星軟著陸是最優控制理論的基準問題之一,由于近期對探索太陽系行星(如火星)的關注增加而重新引起興趣。軟著陸問題在考慮所有相關約束的情況下,可以表述為一個有狀態和控制約束的有限時域最優控制問題。實時生成到行星表面規定位置的燃料最優路徑是一個具有挑戰性的問題,因為存在燃料、控制輸入和狀態的約束。求解這個約束問題的主要困難在于控制輸入上存在非凸約束,這是由于控制輸入幅值的非零下界以及其方向上的非凸約束。本文引入了對控制約束的一種凸化方法,并證明它是無損的;也就是說,軟著陸問題的最優解可以通過求解提出的問題的凸松弛來獲得。無損凸化使得能夠使用凸優化的內點法來獲得原始非凸最優控制問題的最優解。
1.1 符號說明
- R \mathbb{R} R: 實數集
- a.e. [a,b]: 條件幾乎處處成立的區間
- R n \mathbb{R}^n Rn: n維實向量空間
- ∥ v ∥ \|v\| ∥v∥: 向量 v v v的2范數
- 0 0 0: 零矩陣
- I I I: 單位矩陣
- e i e_i ei?: 第 i i i個元素為1,其余為0的列向量
- ( v 1 , … , v m ) (v_1,\dots,v_m) (v1?,…,vm?): 由向量 v 1 , … , v m v_1,\dots,v_m v1?,…,vm?組成的向量,即 ( v 1 , … , v m ) : = [ v 1 T v 2 T … v m T ] T (v_1,\dots,v_m):=\begin{bmatrix}v_1^T & v_2^T & \dots & v_m^T\end{bmatrix}^T (v1?,…,vm?):=[v1T??v2T??…?vmT??]T
- ? S \partial S ?S: 集合 S S S的邊界點集
- i n t S int\,S intS: 集合 S S S的內部
2. 推力指向約束下的行星著陸
行星軟著陸問題尋找推力(控制)剖面 T c T_c Tc?和伴隨的平移狀態軌跡 ( r , r ˙ ) (r,\dot{r}) (r,r˙),使著陸器從初始位置 r 0 r_0 r0?和速度 r ˙ 0 \dot{r}_0 r˙0?導航到行星表面規定的目標位置并以零速度相對于表面,即"軟著陸",同時最小化燃料消耗。該問題考慮具有恒定自轉速率(角速度)、均勻重力場和忽略空氣動力的行星。當從給定初始狀態無法到達目標點時,考慮精確著陸問題(或最小著陸誤差問題),目標是首先找到離目標最近的可達表面位置,然后獲得到該最近點的最小燃料狀態軌跡。我們提出了一種優先級優化方法,在統一的框架下處理這兩個問題,因此稱為行星軟著陸問題。
在這個問題中,有幾個狀態和控制約束。主要的狀態約束是位置向量上的下滑坡道約束和速度向量幅值的上界約束。下滑坡道約束如圖1所示,用于確保著陸器在到達目標前與地面保持安全距離。對于有大氣層的行星,需要對速度設置上限以避免超音速,因為在超音速下控制推力器可能不可靠。這兩個約束都是凸的,適合本文考慮的凸優化框架。控制約束卻具有挑戰性,因為它們定義了一個非凸的可行控制集合。我們有三個控制約束(見圖2):給定任意機動時間(飛行時間) t f t_f tf?,對所有 t ∈ [ 0 , t f ] t\in[0,t_f] t∈[0,tf?]:
- 推力的凸上界: ∥ T c ( t ) ∥ ≤ ρ 2 \|T_c(t)\|\le\rho_2 ∥Tc?(t)∥≤ρ2?
- 推力的非凸下界: ∥ T c ( t ) ∥ ≥ ρ 1 > 0 \|T_c(t)\|\ge\rho_1>0 ∥Tc?(t)∥≥ρ1?>0
- 推力指向約束: n ^ T T c ( t ) / ∥ T c ( t ) ∥ ≥ cos ? θ \hat{n}^TT_c(t)/\|T_c(t)\|\ge\cos\theta n^TTc?(t)/∥Tc?(t)∥≥cosθ,其中 ∥ n ^ ∥ = 1 \|\hat{n}\|=1 ∥n^∥=1是方向向量, 0 ≤ θ ≤ π 0\le\theta\le\pi 0≤θ≤π是相對于 n ^ \hat{n} n^給定方向的最大允許偏差角,當 θ ≤ π / 2 \theta\le\pi/2 θ≤π/2時它是凸的,當 θ > π / 2 \theta>\pi/2 θ>π/2時它是非凸的。
機載地形相對導航傳感器通常需要特定的觀測方向,這對飛行器姿態施加了約束。由于我們將飛行器建模為具有推力向量的點質量,所需的控制力通過將推力向量指向所需力的方向來施加。在這個框架下,我們可以通過簡單地限制推力向量可以指向的方向來對飛行器姿態施加約束。這也避免了將飛行器的姿態動力學合并到問題公式中,否則會顯著增加問題的復雜性。明確考慮姿態動力學并直接對姿態而不是推力方向施加指向約束可以成為未來研究的一部分,這可以受益于最近關于約束姿態控制的凸化結果[23]。
著陸器被建模為具有作為控制的推力矢量的集中質量,其動力學由下式描述:
x ˙ ( t ) = A ( ω ) x ( t ) + B ( g + T c ( t ) m ( t ) ) m ˙ ( t ) = ? α ∥ T c ( t ) ∥ \begin{aligned} \dot{x}(t) &= A(\omega)x(t)+B\left(g+\frac{T_c(t)}{m(t)}\right)\\ \dot{m}(t) &= -\alpha\|T_c(t)\| \end{aligned} x˙(t)m˙(t)?=A(ω)x(t)+B(g+m(t)Tc?(t)?)=?α∥Tc?(t)∥?
其中 x ( t ) = ( r ( t ) , r ˙ ( t ) ) : R + → R 6 x(t)=(r(t),\dot{r}(t)):R_+\to\mathbb{R}^6 x(t)=(r(t),r˙(t)):R+?→R6, m ( t ) : R + → R + m(t):R_+\to R_+ m(t):R+?→R+?是著陸器的質量,
A ( ω ) = [ 0 I ? S ( ω ) 2 ? 2 S ( ω ) ] , B = [ 0 I ] , S ( ω ) = [ 0 ? ω 3 ω 2 ω 3 0 ? ω 1 ? ω 2 ω 1 0 ] A(\omega)=\begin{bmatrix}0 & I\\ -S(\omega)^2 & -2S(\omega)\end{bmatrix},\quad B=\begin{bmatrix}0\\I\end{bmatrix},\quad S(\omega)=\begin{bmatrix} 0 & -\omega_3 & \omega_2\\ \omega_3 & 0 & -\omega_1\\ -\omega_2 & \omega_1 & 0 \end{bmatrix} A(ω)=[0?S(ω)2?I?2S(ω)?],B=[0I?],S(ω)= ?0ω3??ω2???ω3?0ω1??ω2??ω1?0? ?
ω = ( ω 1 , ω 2 , ω 3 ) ∈ R 3 \omega=(\omega_1,\omega_2,\omega_3)\in\mathbb{R}^3 ω=(ω1?,ω2?,ω3?)∈R3是行星的恒定角速度向量, g ∈ R 3 g\in\mathbb{R}^3 g∈R3是恒定重力向量, α > 0 \alpha>0 α>0是描述燃料消耗(質量損耗)率的常數。這里,矢量量的時間導數在固連于行星表面的參考系中表示,該參考系具有行星的角速率,我們還使用了火箭方程,它將燃料質量消耗率與施加的推力矢量聯系起來[24]。
我們使用著陸飛行器的集中質量剛體模型,其中平移動力學與旋轉(姿態)動力學解耦。這是一個在實踐中常用的假設,主要是因為姿態控制權限通常遠高于平移控制權限的帶寬。具體來說,為了使推力器指向正確的平移控制方向而需要的任何姿態機動都可以非常快速地完成,使得姿態和平移控制系統之間的相互作用最小化。因此,這是一個合理的假設,可以大大降低問題的復雜性。
給定約束、動力學方程和表面上的目標位置 ( 0 , q ) (0,q) (0,q),其中 q ∈ R 2 q\in\mathbb{R}^2 q∈R2表示零高度處目標的坐標,行星軟著陸問題可以表述為如下的優先級優化問題。
問題1(非凸最小著陸誤差問題):
min ? t f , T c ∥ E r ( t f ) ? q ∥ s.t. x ˙ ( t ) = A ( ω ) x ( t ) + B ( g + T c ( t ) m ) , m ˙ ( t ) = ? α ∥ T c ( t ) ∥ , ? t ∈ [ 0 , t f ] x ( t ) ∈ X , ? t ∈ [ 0 , t f ] 0 < ρ 1 ≤ ∥ T c ( t ) ∥ ≤ ρ 2 , n ^ T T c ( t ) ≥ ∥ T c ( t ) ∥ cos ? θ m ( 0 ) = m 0 , m ( t f ) ≥ m 0 ? m f > 0 r ( 0 ) = r 0 , r ˙ ( 0 ) = r ˙ 0 e 1 T r ( t f ) = 0 , r ˙ ( t f ) = 0 \begin{aligned} \min_{t_f,T_c}\quad & \|Er(t_f)-q\|\\ \text{s.t.}\quad & \dot{x}(t)=A(\omega)x(t)+B\left(g+\frac{T_c(t)}{m}\right),\quad \dot{m}(t)=-\alpha\|T_c(t)\|,\quad \forall t\in[0,t_f]\\ & x(t)\in X,\quad \forall t\in[0,t_f]\\ & 0<\rho_1\le\|T_c(t)\|\le\rho_2,\quad \hat{n}^TT_c(t)\ge\|T_c(t)\|\cos\theta\\ & m(0)=m_0,\quad m(t_f)\ge m_0-m_f > 0\\ & r(0)=r_0,\quad \dot{r}(0)=\dot{r}_0\\ & e_1^Tr(t_f)=0,\quad \dot{r}(t_f)=0 \end{aligned} tf?,Tc?min?s.t.?∥Er(tf?)?q∥x˙(t)=A(ω)x(t)+B(g+mTc?(t)?),m˙(t)=?α∥Tc?(t)∥,?t∈[0,tf?]x(t)∈X,?t∈[0,tf?]0<ρ1?≤∥Tc?(t)∥≤ρ2?,n^TTc?(t)≥∥Tc?(t)∥cosθm(0)=m0?,m(tf?)≥m0??mf?>0r(0)=r0?,r˙(0)=r˙0?e1T?r(tf?)=0,r˙(tf?)=0?
問題2(非凸最小燃料問題):
max ? t f , T c m ( t f ) ? m ( 0 ) = min ? t f , T c ∫ 0 t f α ∥ T c ( t ) ∥ d t s.t. 動力學方程和約束(4)-(9) ∥ E r ( t f ) ? q ∥ ≤ ∥ d P 1 ? ? q ∥ \begin{aligned} \max_{t_f,T_c}\quad & m(t_f)-m(0)=\min_{t_f,T_c}\int_0^{t_f}\alpha\|T_c(t)\|\,dt\\ \text{s.t.}\quad & \text{動力學方程和約束(4)-(9)}\\ & \|Er(t_f)-q\|\le\|d_{P1}^*-q\| \end{aligned} tf?,Tc?max?s.t.?m(tf?)?m(0)=tf?,Tc?min?∫0tf??α∥Tc?(t)∥dt動力學方程和約束(4)-(9)∥Er(tf?)?q∥≤∥dP1???q∥?
其中, d P 1 ? = E r P 1 ? ( t f ) ∈ R 2 d_{P1}^*=Er_{P1}^*(t_f)\in\mathbb{R}^2 dP1??=ErP1??(tf?)∈R2是通過求解問題1得到的最優成本對應的最終位置,即 d P 1 ? d_{P1}^* dP1??是離目標位置 q q q最近的可達表面點, m f m_f mf?是可用燃料, m 0 m_0 m0?是著陸器的初始質量,
E = [ e 2 T e 3 T ] = [ 0 1 0 0 0 1 ] . E=\begin{bmatrix}e_2^T \\ e_3^T\end{bmatrix}=\begin{bmatrix} 0 & 1 & 0\\ 0 & 0 & 1 \end{bmatrix}. E=[e2T?e3T??]=[00?10?01?].
我們用 X X X來定義飛行器可行位置和速度的集合,即下滑坡道約束和最大允許速度 V m a x V_{max} Vmax?的約束
X = { ( r , r ˙ ) ∈ R 6 : ∥ r ˙ ∥ ≤ V m a x , ∥ E ( r ? r ( t f ) ) ∥ ? c T ( r ? r ( t f ) ) ≤ 0 } X=\left\{(r,\dot{r})\in\mathbb{R}^6:\|\dot{r}\|\le V_{max},\quad \|E(r-r(t_f))\|-c^T(r-r(t_f))\le 0\right\} X={(r,r˙)∈R6:∥r˙∥≤Vmax?,∥E(r?r(tf?))∥?cT(r?r(tf?))≤0}
其中 c c c指定了一個頂點在 r ( t f ) r(t_f) r(tf?)處的可行錐
c : = e 1 tan ? γ g s , γ g s ∈ ( 0 , π / 2 ) . c:=\frac{e_1}{\tan\gamma_{gs}},\quad\gamma_{gs}\in(0,\pi/2). c:=tanγgs?e1??,γgs?∈(0,π/2).
這里, γ g s \gamma_{gs} γgs?是最小下滑坡道角,如圖1所示。下滑坡道約束(12)確保到目標的軌跡不能太淺或進入地下。 X X X是一個凸集,為了完整起見,我們給出 X X X內部的標準定義
i n t X : = { x ∈ X : ? ε > 0 使得? y ∈ X 如果? ∥ x ? y ∥ < ε } . int\,X:=\{x\in X:\exists\varepsilon>0\ \text{使得}\ y\in X\ \text{如果}\ \|x-y\|<\varepsilon\}. intX:={x∈X:?ε>0?使得?y∈X?如果?∥x?y∥<ε}.
X X X的邊界由 ? X : = { x ∈ X : x ? i n t X } \partial X:=\{x\in X:x\notin int\,X\} ?X:={x∈X:x∈/intX}給出。等式(7)定義了著陸器的初始質量,并確保使用的燃料不超過可用燃料。等式(8)定義了著陸器的初始位置和速度,而(9)約束了最終高度和速度。飛行時間 t f t_f tf?是一個優化變量,不是先驗固定的。
軟著陸問題的解通過求解問題1和問題2獲得。這種兩步優先級方法的動機是非常直觀的。本文考慮的行星著陸問題的主要目標是使飛行器盡可能接近給定目標著陸,即像問題1中那樣最小化著陸誤差。這個問題可能有多個最優解,因此我們有第二步,即像問題2中那樣找到消耗最少燃料的最小誤差解。
備注1: 在問題2中,最終位置的約束由不等式(11)給出。這個約束也可以是以下選擇之一
E r ( t f ) = d P 1 ? Er(t_f)=d_{P1}^* Er(tf?)=dP1??
∥ E r ( t f ) ? q ∥ = ∥ d P 1 ? ? q ∥ . \|Er(t_f)-q\|=\|d_{P1}^*-q\|. ∥Er(tf?)?q∥=∥dP1???q∥.
第三個選擇是由(11)給出的約束,它是第二個選擇的凸松弛,包括嚴格在這個圓內的終點位置(見圖3)。 顯然,通過這種進一步的松弛,我們永遠無法計算出最終位置嚴格在圓內的解,因為這個圓的半徑是可達到目標的最小距離。我們使用第三種選擇有兩個動機:(1)獲得具有盡可能接近期望目標的最終位置的最小燃料解;(2)使問題2的所有約束都是凸的。為了了解我們如何實現這兩個目標,令 F 1 F_1 F1?表示問題2的所有可行解,其中(11)被(15)替換。 F 2 F_2 F2?是相應的集合,其中(11)被(16)替換,而 F e F_e Fe?是問題2的集合。那么很容易注意到 F 1 ? F 2 ? F e F_1\subseteq F_2\subseteq F_e F1??F2??Fe?。此外,我們有 F 2 = F e F_2=F_e F2?=Fe?。這是因為不等式(11)不排除終點位置在圓上的解,而且我們不能進入圓內。現在,由于 F 1 ? F 2 F_1\subseteq F_2 F1??F2?,上述第二種選擇將產生使用最小燃料軌道并盡可能接近目標的解。所以最好的選擇是使用第二種選擇。但 ∥ E r ( t f ) ? q ∥ = ∥ d P 1 ? ? q ∥ \|Er(t_f)-q\|=\|d_{P1}^*-q\| ∥Er(tf?)?q∥=∥dP1???q∥是終點位置上的非凸約束。由于 F 2 = F e F_2=F_e F2?=Fe?,我們可以簡單地用凸不等式(11)替換這個非凸約束。因此,我們通過系統地優先考慮兩個成本,即可達距離和燃料,來實現兩個主要目標:(1)明確優先考慮兩個成本;(2)使兩個問題都凸化,以便可以通過內點法算法在多項式時間內求解到全局最優。注意,我們通過系統地優先考慮成本來實現這兩個目標。這種優化問題中的優先成本也稱為字典式目標規劃[25]。
定理1的關鍵挑戰是(6)中對推力幅值的下界 ρ 1 > 0 \rho_1>0 ρ1?>0,這意味著允許推力值的集合是非凸的(見圖2中的第一個插圖)。此外,當 ρ 1 = 0 \rho_1 = 0 ρ1?=0時,推力邊界約束是凸的;然而,由于(6)中的推力指向約束,控制約束仍然可能是非凸的,當 θ > π / 2 \theta>\pi/2 θ>π/2時它是非凸的(見圖2中的第二個和第三個插圖)。這些非凸控制約束阻止了直接使用凸優化技術來求解這個問題。此外,(4)中質量消耗 m ˙ ( t ) \dot{m}(t) m˙(t)的動力學定義了一個非線性微分方程,離散化后會導致非線性等式約束,這也是非凸的。
[1]中的關鍵結果包括對非凸推力邊界約束的松弛以及處理質量消耗動力學的方法,提供了問題2的松弛版本。這個松弛問題的最優解被證明也是問題2的最優解。然而,[1]中的凸化結果在存在任何推力指向約束時不成立,包括 θ ∈ [ 0 , π / 2 ] \theta\in[0,\pi/2] θ∈[0,π/2]的情況。本文將控制約束的凸化擴展到也包含推力指向約束的情況。在凸化具有非凸推力指向約束的問題時,我們發展了問題的幾何洞察,建立了與"正則系統"的聯系[18]。正則線性系統是在最優控制理論的背景下定義的,如果系統在可行控制集合的唯一點處最大化Hamilton函數,則稱該系統相對于一組可行控制是正則的。在可行控制集合是凸的情況下,系統正則意味著Hamilton函數在集合的極點處最大化[18]。我們的凸化結果具有類似的幾何解釋,因為它通過確保Hamilton函數在松弛可行控制集的投影的極點處最大化來建立無損凸化。然后證明這個集合包含在原始非凸可行控制集中,從而確立了我們可以通過求解其凸松弛來獲得原始非凸問題的最優解。
軟著陸問題的無損凸化的理論發展允許應用凸優化的內點法[10]-[12],[19],這些方法可以可靠地求解這些問題,并具有多項式時間收斂保證。此外,對快速實時凸優化領域的興趣激增[20]-[22]表明,內點法的計算速度提高了幾個數量級。這一領域的進步將顯著提高內點法的實時計算效率,從而使其能夠用于行星軟著陸。
1.2 符號列表
- R \mathbb{R} R是實數集。
- 如果條件在 [ a , b ] [a,b] [a,b]中不成立的點集是零測度集,則稱條件幾乎處處成立,記為a.e. [ a , b ] [a,b] [a,b]。
- R n \mathbb{R}^n Rn是n維實向量空間。
- ∥ v ∥ \|v\| ∥v∥是向量 v v v的2范數。
- 0 0 0是零矩陣。
- I I I是單位矩陣。
- e i e_i ei?是第 i i i個元素為1,其余為0的列向量。
- ( v 1 , v 2 , … , v m ) (v_1,v_2,\dots,v_m) (v1?,v2?,…,vm?)表示通過增廣向量 v 1 , … , v m v_1,\dots,v_m v1?,…,vm?得到的向量,使得 ( v 1 , v 2 , … , v m ) : = [ v 1 T v 2 T … s v m T ] T (v_1,v_2,\dots,v_m):=\begin{bmatrix}v_1^T & v_2^T & \dots & sv_m^T\end{bmatrix}^T (v1?,v2?,…,vm?):=[v1T??v2T??…?svmT??]T。
- ? S \partial S ?S表示集合 S S S的邊界點集, i n t S int\,S intS表示 S S S的內部。
3. 無損凸化
我們提出以下問題3和4作為問題1和2的凸松弛。
問題3(凸松弛最小著陸誤差問題):
min ? t f , T c , σ ∥ E r ( t f ) ? q ∥ s.t. (5),?(7),?(8),?(9) , x ˙ ( t ) = A ( ω ) x ( t ) + B ( g + T c ( t ) m ) , m ˙ ( t ) = ? α σ ( t ) , ? t ∈ [ 0 , t f ] ∥ T c ( t ) ∥ ≤ σ ( t ) , 0 < ρ 1 ≤ σ ( t ) ≤ ρ 2 n ^ T T c ( t ) ≥ cos ? θ σ ( t ) \begin{aligned} \min_{t_f,T_c,\sigma}\quad & \|Er(t_f)-q\|\\ \text{s.t.}\quad & \text{(5), (7), (8), (9)},\\ & \begin{aligned} \dot{x}(t) &= A(\omega)x(t)+B\left(g+\frac{T_c(t)}{m}\right),\\ \dot{m}(t) &= -\alpha\sigma(t),\quad\forall t\in[0,t_f] \end{aligned}\\ & \|T_c(t)\|\le \sigma(t),\quad 0<\rho_1\le\sigma(t)\le\rho_2\\ & \hat{n}^TT_c(t)\ge\cos\theta\,\sigma(t) \end{aligned} tf?,Tc?,σmin?s.t.?∥Er(tf?)?q∥(5),?(7),?(8),?(9),x˙(t)m˙(t)?=A(ω)x(t)+B(g+mTc?(t)?),=?ασ(t),?t∈[0,tf?]?∥Tc?(t)∥≤σ(t),0<ρ1?≤σ(t)≤ρ2?n^TTc?(t)≥cosθσ(t)?
問題4(凸松弛最小燃料問題):
min ? t f , T c , σ ∫ 0 t f σ ( t ) d t s.t. (5),?(7),?(8),?(9),?(17),?(18),?(19), ∥ E r ( t f ) ? q ∥ ≤ ∥ d P 3 ? ? q ∥ \begin{aligned} \min_{t_f,T_c,\sigma}\quad & \int_0^{t_f}\sigma(t)\,dt\\ \text{s.t.}\quad & \text{(5), (7), (8), (9), (17), (18), (19),}\\ & \|Er(t_f)-q\|\le\|d_{P3}^*-q\| \end{aligned} tf?,Tc?,σmin?s.t.?∫0tf??σ(t)dt(5),?(7),?(8),?(9),?(17),?(18),?(19),∥Er(tf?)?q∥≤∥dP3???q∥?
其中 d P 3 ? d_{P3}^* dP3??是問題3的最終最優位置。注意,問題1和2中非凸推力約束(6)已經被問題3和4中的凸約束(18)和(19)所取代。在[1]中,我們表明了對推力邊界的約束松弛(18)允許將問題4的離散時間形式表述為凸優化問題;加上凸推力指向約束(19),同樣的結論也成立。因此,我們不會在本文中討論這個問題的離散化,讀者可參考[1]。
定義1: F e F_e Fe?和 F f F_f Ff?分別表示問題1和2的可行解集,即如果 { t f , T c , x , m } \{t_f,T_c,x,m\} {tf?,Tc?,x,m}滿足問題1的所有狀態(5)、控制(6)、燃料(7)約束、動力學(4)和邊界條件(8)和(9),則 { t f , T c , x , m } ∈ F e \{t_f,T_c,x,m\}\in F_e {tf?,Tc?,x,m}∈Fe?,對 F f F_f Ff?也類似。 F e ? F_e^* Fe??和 F f ? F_f^* Ff??是相應的最優解集。 F r e F_{re} Fre?和 F r f F_{rf} Frf?是問題3和4的可行解 { t f , T c , σ , x , m } \{t_f,T_c,\sigma,x,m\} {tf?,Tc?,σ,x,m}的集合, F r e ? F_{re}^* Fre??和 F r f ? F_{rf}^* Frf??是最優解集。
引理1: 以下結論成立:
- F f ? F e F_f\subseteq F_e Ff??Fe?, F r f ? F r e F_{rf}\subseteq F_{re} Frf??Fre?, F f ? ? F e ? F_f^*\subseteq F_e^* Ff???Fe??, 且 F r f ? ? F r e ? F_{rf}^*\subseteq F_{re}^* Frf???Fre??;
- 如果 { t f , T c , σ , x , m } ∈ F r f \{t_f,T_c,\sigma,x,m\}\in F_{rf} {tf?,Tc?,σ,x,m}∈Frf?使得 ∥ T c ( t ) ∥ = σ ( t ) , ? [ 0 , t f ] \|T_c(t)\|=\sigma(t), \forall[0,t_f] ∥Tc?(t)∥=σ(t),?[0,tf?],則 { t f , T c , x , m } ∈ F f \{t_f,T_c,x,m\}\in F_f {tf?,Tc?,x,m}∈Ff?;
- 如果 { t f ? , T c ? , σ ? , x ? , m ? } ∈ F r f ? \{t_f^*,T_c^*,\sigma^*,x^*,m^*\}\in F_{rf}^* {tf??,Tc??,σ?,x?,m?}∈Frf??使得 ∥ T c ? ( t ) ∥ = σ ? ( t ) \|T_c^*(t)\|=\sigma^*(t) ∥Tc??(t)∥=σ?(t), a.e. [ 0 , t f ? ] [0,t_f^*] [0,tf??],則 { t f ? , T c ? , x ? , m ? } ∈ F f ? \{t_f^*,T_c^*,x^*,m^*\}\in F_f^* {tf??,Tc??,x?,m?}∈Ff??。
證明:
- 中的前兩個關系很容易證明。接下來的兩個關系來自于前兩個。考慮 F r f F_{rf} Frf?的一個子集 F  ̄ r f \overline{F}_{rf} Frf?,它由所有滿足 σ ( t ) = ∥ T c ( t ) ∥ , ? t ∈ [ 0 , t f ] \sigma(t)=\|T_c(t)\|,\forall t\in[0,t_f] σ(t)=∥Tc?(t)∥,?t∈[0,tf?]的 F r f F_{rf} Frf?的解組成。由于 { t f , T c , σ , x , m } ∈ F  ̄ r f \{t_f,T_c,\sigma,x,m\}\in\overline{F}_{rf} {tf?,Tc?,σ,x,m}∈Frf?,其中 ∥ T c ( t ) ∥ = σ ( t ) \|T_c(t)\|=\sigma(t) ∥Tc?(t)∥=σ(t),滿足問題2的所有動力學和約束條件, { t f , T c , x , m } ∈ F f \{t_f,T_c,x,m\}\in F_f {tf?,Tc?,x,m}∈Ff?。這證明了2)。因此,由于目標函數對于 ∥ T c ∥ = σ \|T_c\|=\sigma ∥Tc?∥=σ的問題4和問題2是相同的,定理的最后一個陳述成立。
為了利用引理1,我們需要計算問題4的最優解并檢查 σ ( t ) = ∥ T c ( t ) ∥ \sigma(t)=\|T_c(t)\| σ(t)=∥Tc?(t)∥是否 ? t ∈ [ 0 , t f ? ] \forall t\in[0,t_f^*] ?t∈[0,tf??]。但是在數值計算解之前,我們不知道這個條件是否會滿足。下面的定理建立了這個條件通常會滿足,但問題可能需要稍加修改,引入 ω ^ \hat{\omega} ω^和 ε \varepsilon ε。它還提供了對[1]和[2]中早期結果的推廣,以處理推力指向約束以及推力幅值的非零下界。因此,它建立了最小燃料軟著陸問題(問題2,從而問題1)中控制約束的無損凸化。
定理1: 考慮在問題4中用 ω ^ \hat{\omega} ω^替換 ω \omega ω,其中 ω ^ \hat{\omega} ω^定義如下:
ω ^ : = { ω 如果? S ( ω ) n ^ = 0 , N T S ( ω ) 2 n ^ = 0 ω + ε n ^ ⊥ 如果? S ( ω ) n ^ = 0 ω + ε n ^ 如果? S ( ω ) n ^ = 0 , N T S ( ω ) 2 n ^ ≠ 0 \hat{\omega}:=\begin{cases} \omega & \text{如果}\ S(\omega)\hat{n}=0,\, N^TS(\omega)^2\hat{n}=0\\ \omega+\varepsilon\hat{n}^\perp & \text{如果}\ S(\omega)\hat{n}=0\\ \omega+\varepsilon\hat{n} & \text{如果}\ S(\omega)\hat{n}=0,\, N^TS(\omega)^2\hat{n}\ne 0 \end{cases} ω^:=? ? ??ωω+εn^⊥ω+εn^?如果?S(ω)n^=0,NTS(ω)2n^=0如果?S(ω)n^=0如果?S(ω)n^=0,NTS(ω)2n^=0?
其中 n ^ ⊥ \hat{n}^\perp n^⊥是滿足 n ^ T n ^ ⊥ = 0 \hat{n}^T\hat{n}^\perp=0 n^Tn^⊥=0的單位向量, N ∈ R 3 × 2 N\in\mathbb{R}^{3\times 2} N∈R3×2的列張成 n ^ T \hat{n}^T n^T的零空間, ε > 0 \varepsilon>0 ε>0是一個(任意小的)實數。設 { t f ? , T c ? , σ ? , x ? , m ? } ∈ F r f ? \{t_f^*,T_c^*,\sigma^*,x^*,m^*\}\in F_{rf}^* {tf??,Tc??,σ?,x?,m?}∈Frf??使得相應的狀態軌跡 x ? ( t ) ∈ i n t X , ? t ∈ [ 0 , t f ? ] x^*(t)\in int\,X,\forall t\in[0,t_f^*] x?(t)∈intX,?t∈[0,tf??]。則有 { t f ? , T c ? , x ? , m ? } ∈ F f ? \{t_f^*,T_c^*,x^*,m^*\}\in F_f^* {tf??,Tc??,x?,m?}∈Ff??,其中的 ω \omega ω用 ω ^ \hat{\omega} ω^替換。
證明:
本證明使用了附錄中的引理2和3。
令 q ~ : = E r ? ( t f ? ) \tilde{q}:=Er^*(t_f^*) q~?:=Er?(tf??)。那么我們也可以考慮 { t f ? , T c ? , σ ? , x ? , m ? } \{t_f^*,T_c^*,\sigma^*,x^*,m^*\} {tf??,Tc??,σ?,x?,m?}作為問題4的最優燃料解,其中約束 ∥ E r ( t f ) ? q ∥ ≤ ∥ d P 1 ? ? q ∥ \|Er(t_f)-q\|\le\|d_{P1}^*-q\| ∥Er(tf?)?q∥≤∥dP1???q∥被 E r ( t f ) = q ~ Er(t_f)=\tilde{q} Er(tf?)=q~?代替。因此,在本證明中將使用這個版本的問題4,而不失一般性。
由于 x ? ( t ) ∈ i n t X x^*(t)\in int\,X x?(t)∈intX且對所有 t ∈ [ 0 , t f ] t\in[0,t_f] t∈[0,tf?]有 m ( t ) > m 0 ? m f m(t)>m_0-m_f m(t)>m0??mf?,根據最優控制的極大值原理[見[18,第V.3節]或[26,第1章]],存在常數 β ≤ 0 \beta\le 0 β≤0和絕對連續函數 λ : R + → R 6 \lambda:R_+\to\mathbb{R}^6 λ:R+?→R6和 η : R + → R \eta:R_+\to\mathbb{R} η:R+?→R,即余態向量,使得以下條件成立。
-
余態條件: ? t ∈ [ 0 , t f ? ] \forall t\in[0,t_f^*] ?t∈[0,tf??]
( β , λ ( t ) , η ( t ) ) ≠ 0 (\beta,\lambda(t),\eta(t))\ne 0 (β,λ(t),η(t))=0
λ ˙ ( t ) = ? A ( ω ^ ) T λ ( t ) \dot{\lambda}(t)=-A(\hat{\omega})^T\lambda(t) λ˙(t)=?A(ω^)Tλ(t)
η ˙ ( t ) = λ ( t ) T B T c ( t ) m ( t ) 2 . \dot{\eta}(t)=\frac{\lambda(t)^TBT_c(t)}{m(t)^2}. η˙?(t)=m(t)2λ(t)TBTc?(t)?. -
逐點極大值原理:
H ( φ ( t ) ) = M ( x ? ( t ) , m ? ( t ) , λ ( t ) , η ( t ) ) a.e.? t ∈ [ 0 , t f ? ] H(\varphi(t))=M(x^*(t),m^*(t),\lambda(t),\eta(t))\quad \text{a.e.}\ t\in[0,t_f^*] H(φ(t))=M(x?(t),m?(t),λ(t),η(t))a.e.?t∈[0,tf??]
其中 φ ( t ) = ( t , x ? ( t ) , m ? ( t ) , T c ? ( t ) , σ ? ( t ) , λ ( t ) , η ( t ) ) \varphi(t)=(t,x^*(t),m^*(t),T_c^*(t),\sigma^*(t),\lambda(t),\eta(t)) φ(t)=(t,x?(t),m?(t),Tc??(t),σ?(t),λ(t),η(t)), H H H是如下定義的Hamilton函數
H ( φ ) : = β + λ T ( A ( ω ^ ) x + B ( g + T c m ) ) ? α η H(\varphi):=\beta+\lambda^T\left(A(\hat{\omega})x+B\left(g+\frac{T_c}{m}\right)\right)-\alpha\eta H(φ):=β+λT(A(ω^)x+B(g+mTc??))?αη
并且,通過定義$V:={(T_c,\sigma)\in\mathbb{R}4:|T_c|\le\sigma,\rho_1\le\sigma\le\rho_2,\hat{n}TT_c\證明的后續部分如下:
≥ cos ? θ σ } \ge\cos\theta\,\sigma\} ≥cosθσ},有
M ( x ? , m ? , λ , η ) = max ? ( T c , σ ) ∈ V H ( φ ) . M(x^*,m^*,\lambda,\eta)=\max_{(T_c,\sigma)\in V}H(\varphi). M(x?,m?,λ,η)=(Tc?,σ)∈Vmax?H(φ).
- 橫截條件:
η ( t f ? ) = 0 和 H ( φ ( t f ? ) ) = 0. \eta(t_f^*)=0\quad\text{和}\quad H(\varphi(t_f^*))=0. η(tf??)=0和H(φ(tf??))=0.
以上最優性的必要條件1)和2)直接來自于極大值原理的陳述。但橫截條件需要進一步解釋。橫截條件意味著(見[18,第190頁,第V.3節]),對于松弛問題的最優解,由下式定義的向量 ψ \psi ψ
ψ : = ( H ( φ ( 0 ) ) , H ( φ ( t f ? ) ) , ? λ ( 0 ) , ? η ( 0 ) , λ ( t f ? ) , η ( t f ? ) ) \psi:=\left(H(\varphi(0)),H(\varphi(t_f^*)),-\lambda(0),-\eta(0),\lambda(t_f^*),\eta(t_f^*)\right) ψ:=(H(φ(0)),H(φ(tf??)),?λ(0),?η(0),λ(tf??),η(tf??))
必須正交于由初始和最終狀態的可行集合描述的流形,該流形由 ( 0 , t f , x 0 , m 0 , ( 0 , q ~ , 0 ) , m ( t f ? ) ) (0,t_f,x_0,m_0,(0,\tilde{q},0),m(t_f^*)) (0,tf?,x0?,m0?,(0,q~?,0),m(tf??))給出,為 s p a n { e 2 , e 14 } span\{e_2,e_{14}\} span{e2?,e14?}。上述結論來自這樣一個事實,即 t f t_f tf?和 m ( t f ) m(t_f) m(tf?)是邊界條件流形中唯一的自由變量。這就意味著 e 2 T ψ = 0 e_2^T\psi=0 e2T?ψ=0和 e 14 T ψ = 0 e_{14}^T\psi=0 e14T?ψ=0,即 H ( φ ( t f ? ) ) = 0 H(\varphi(t_f^*))=0 H(φ(tf??))=0和 η ( t f ? ) = 0 \eta(t_f^*)=0 η(tf??)=0。接下來我們證明
y ( t ) : = B T λ ( t ) = 0 a.e.? [ 0 , t f ? ] . y(t):=B^T\lambda(t)=0\quad\text{a.e.}\ [0,t_f^*]. y(t):=BTλ(t)=0a.e.?[0,tf??].
這將通過反證法來完成。假設條件(29)不成立。由于 y y y是系統(23)的輸出, y ( t ) = 0 , ? t ∈ [ 0 , t f ? ] y(t)=0,\forall t\in[0,t_f^*] y(t)=0,?t∈[0,tf??],或 y ( t ) = 0 y(t)=0 y(t)=0在 [ 0 , t f ? ] [0,t_f^*] [0,tf??]中只在可數個點處出現,這來自于引理2的第一個結論。假設 y ( t ) = 0 , ? t ∈ [ 0 , t f ? ] y(t)=0,\forall t\in[0,t_f^*] y(t)=0,?t∈[0,tf??]。注意到對 [ A ( ω ) , B ] [A(\omega),B] [A(ω),B]是可控的,這來自于 [ B A ( ω ^ ) B ] [B\ A(\hat{\omega})B] [B?A(ω^)B]是可逆矩陣這一事實。因此對 ( B T , ? A ( ω ^ ) T ) (B^T,-A(\hat{\omega})^T) (BT,?A(ω^)T)是可觀測的。因此, y ( t ) = 0 , ? t ∈ [ 0 , t f ? ] y(t)=0,\forall t\in[0,t_f^*] y(t)=0,?t∈[0,tf??]意味著 λ ( t ) = 0 , ? t ∈ [ 0 , t f ? ] \lambda(t)=0,\forall t\in[0,t_f^*] λ(t)=0,?t∈[0,tf??]。因此 η ˙ ( t ) = 0 , ? t ∈ [ 0 , t f ? ] \dot{\eta}(t)=0,\forall t\in[0,t_f^*] η˙?(t)=0,?t∈[0,tf??]。由于 η ( t f ? ) = 0 \eta(t_f^*)=0 η(tf??)=0,這就意味著 η ( t ) = 0 , ? t ∈ [ 0 , t f ? ] \eta(t)=0,\forall t\in[0,t_f^*] η(t)=0,?t∈[0,tf??]。這些都意味著 H ( φ ( t ) ) = β σ ( t ) H(\varphi(t))=\beta\sigma(t) H(φ(t))=βσ(t)。由于 H ( φ ( t f ? ) ) = 0 H(\varphi(t_f^*))=0 H(φ(tf??))=0且 σ ( t ) ≥ ρ 1 > 0 \sigma(t)\ge\rho_1>0 σ(t)≥ρ1?>0,這表明 β = 0 \beta=0 β=0。因此, ( β , λ ( t ) , η ( t ) ) = 0 , ? t ∈ [ 0 , t f ? ] (\beta,\lambda(t),\eta(t))=0,\forall t\in[0,t_f^*] (β,λ(t),η(t))=0,?t∈[0,tf??],這與上述必要條件1)矛盾。因此,在 [ 0 , t f ? ] [0,t_f^*] [0,tf??]中有可數個點處 y ( t ) = 0 y(t)=0 y(t)=0。由于可數集的測度為零,引理2的第二個結論意味著
y ( t ) = ? α ( t ) n ^ a.e.? [ 0 , t f ? ] , 其中 α ( t ) > 0. y(t)=-\alpha(t)\hat{n}\quad\text{a.e.}\ [0,t_f^*], \text{其中}\alpha(t)>0. y(t)=?α(t)n^a.e.?[0,tf??],其中α(t)>0.
由于條件(29)成立,a.e. [ 0 , t f ? ] [0,t_f^*] [0,tf??]使得 y ( t ) ≠ 0 y(t)\ne 0 y(t)=0,且對于給定的 σ ? ( t ) \sigma^*(t) σ?(t),最優控制推力 T c ? ( t ) T_c^*(t) Tc??(t)必須滿足
T c ? ( t ) = arg ? max ? ( T c , σ ? ( t ) ) ∈ V y ( t ) T T c = arg ? max ? T c ∈ U ( σ ? ) y ( t ) T T c T_c^*(t)=\arg\max_{(T_c,\sigma^*(t))\in V}y(t)^TT_c=\arg\max_{T_c\in U(\sigma^*)}y(t)^TT_c Tc??(t)=arg(Tc?,σ?(t))∈Vmax?y(t)TTc?=argTc?∈U(σ?)max?y(t)TTc?
其中 U ( σ ) : = { T c ∈ R 3 : ∥ T c ∥ ≤ σ , n ^ T T c ≥ cos ? θ σ } U(\sigma):=\{T_c\in\mathbb{R}^3:\|T_c\|\le\sigma,\hat{n}^TT_c\ge\cos\theta\,\sigma\} U(σ):={Tc?∈R3:∥Tc?∥≤σ,n^TTc?≥cosθσ}。此外,由于條件(30)成立,a.e. [ 0 , t f ? ] [0,t_f^*] [0,tf??]使得 y ( t ) ≠ 0 y(t)\ne 0 y(t)=0且 y ( t ) = ? α ( t ) n ^ y(t)=-\alpha(t)\hat{n} y(t)=?α(t)n^對某些 α ( t ) > 0 \alpha(t)>0 α(t)>0成立,(31)的最大化解必須在 U ( σ ? ) U(\sigma^*) U(σ?)的邊界點處,而且是 U ( σ ? ) U(\sigma^*) U(σ?)的極點,這來自于引理3。該引理還意味著 U ( σ ) U(\sigma) U(σ)的所有極點滿足 ∥ T c ∥ = σ \|T_c\|=\sigma ∥Tc?∥=σ,因此 ∥ T c ? ( t ) ∥ = σ ? ( t ) \|T_c^*(t)\|=\sigma^*(t) ∥Tc??(t)∥=σ?(t)
∥ T c ? ( t ) ∥ = σ ? ( t ) a.e.? [ 0 , t f ? ] \|T_c^*(t)\|=\sigma^*(t)\quad\text{a.e.}\ [0,t_f^*] ∥Tc??(t)∥=σ?(t)a.e.?[0,tf??]
這意味著松弛問題(4)的最優解滿足
0 < ρ 1 ≤ ∥ T c ? ( t ) ∥ ≤ ρ 2 , n ^ T T c ? ( t ) ≥ ∥ T c ? ( t ) ∥ cos ? θ a.e.? [ 0 , t f ? ] . 0<\rho_1\le\|T_c^*(t)\|\le\rho_2,\quad \hat{n}^TT_c^*(t)\ge\|T_c^*(t)\|\cos\theta\quad\text{a.e.}\ [0,t_f^*]. 0<ρ1?≤∥Tc??(t)∥≤ρ2?,n^TTc??(t)≥∥Tc??(t)∥cosθa.e.?[0,tf??].
因此, ( t f ? , T c ? , x ? , m ? ) ∈ F f (t_f^*,T_c^*,x^*,m^*)\in F_f (tf??,Tc??,x?,m?)∈Ff?。由于對任意 ( t f , T c , x , m ) ∈ F f (t_f,T_c,x,m)\in F_f (tf?,Tc?,x,m)∈Ff?,有 ( t f , T c , ∥ T c ∥ , x , m ) ∈ F r f (t_f,T_c,\|T_c\|,x,m)\in F_{rf} (tf?,Tc?,∥Tc?∥,x,m)∈Frf?,問題4的最優解具有不大于問題2最優成本的最優成本。這意味著根據引理1, ( t f ? , T c ? , x ? , m ? ) ∈ F f ? (t_f^*,T_c^*,x^*,m^*)\in F_f^* (tf??,Tc??,x?,m?)∈Ff??。
上述定理指出,我們可以通過求解其在問題4中的松弛來找到問題2的最優解,其中 ω \omega ω用 ω ^ \hat{\omega} ω^代替。顯然,對于 ω = ω ^ \omega=\hat{\omega} ω=ω^,我們可以通過求解其松弛來找到感興趣的原始問題的精確最優解。當 ω ^ ≠ ω \hat{\omega}\ne\omega ω^=ω時,我們可以找到一個問題的最優解,該問題可以通過簡單地選擇 ε > 0 \varepsilon>0 ε>0接近于零而變得任意接近于感興趣的問題。此外,當 θ = π \theta=\pi θ=π,即當沒有指向約束時,我們可以使用任何單位向量作為 n ^ \hat{n} n^,使得 ω = ω ^ \omega=\hat{\omega} ω=ω^。因此,我們可以找到感興趣的原始問題的精確最優解。
值得注意的是,我們的凸化結果依賴于對 [ A ( ω ^ ) , B ] [A(\hat{\omega}),B] [A(ω^),B]是可控的。因此,只要保持這個可控性,它仍然適用于系統參數變化的情況。這顯然不是一個限制性條件,因為可以合理地期望行星著陸飛行器的設計是可控的。
另一個重要的問題是對動態模型中未知參數變化或干擾力的不確定性的魯棒性。這里我們提供了一種計算要遵循的最優狀態軌跡和相應最優控制的方法。雖然本文沒有涉及,但必須考慮反饋控制作用來處理不確定性和干擾。我們可以提出幾種方法來實現這一點:(1)有一個跟蹤反饋控制器來跟蹤最優狀態軌跡;(2)隨著狀態知識(估計)在機動過程中更新而重新計算最優軌跡和控制;(3)跟蹤反饋和最優控制重計算的組合。所有行動都將利用當前最佳的狀態估計,并提供反饋作用以最小化不確定性和干擾的影響。魯棒反饋作用的設計將是我們當前論文的另一個未來擴展。
一個自然的問題是,這里的凸化方法是否可以擴展到其他控制約束。一些針對線性系統的新結果在[27]中提供,其中這種凸化過程應用于大類控制約束。然而,那里的凸化利用了最優控制位于松弛控制集的邊界點而不是極點這一事實。根據本文中獲得的關于松弛集極點的見解,可以進一步推廣這種凸化方法,這可以成為未來研究的主題。
3.1 變量替換
我們對推力向量和質量使用以下變量替換,以消除由 T c / m T_c/m Tc?/m引起的動力學中的非線性: σ : = m f m , u : = T c m , z : = ln ? m \sigma:=\frac{m_f}{m},u:=\frac{T_c}{m},z:=\ln m σ:=mmf??,u:=mTc??,z:=lnm。
質量損耗動力學然后可以重寫為
z ˙ = m ˙ ( t ) m ( t ) = ? α σ ( t ) . \dot{z}=\frac{\dot{m}(t)}{m(t)}=-\alpha\sigma(t). z˙=m(t)m˙(t)?=?ασ(t).
因此,變量替換產生了一組狀態動力學的線性方程。然而,控制約束不再是凸的。現在它們由下式給出
∥ u ( t ) ∥ ≤ σ ( t ) , n ^ T u ( t ) ≥ cos ? θ σ ( t ) , ? t ∈ [ 0 , t f ] \|u(t)\|\le\sigma(t),\quad \hat{n}^Tu(t)\ge\cos\theta\,\sigma(t),\quad \forall t\in[0,t_f] ∥u(t)∥≤σ(t),n^Tu(t)≥cosθσ(t),?t∈[0,tf?]
ρ 1 e ? z ( t ) ≤ σ ( t ) ≤ ρ 2 e ? z ( t ) , ? t ∈ [ 0 , t f ] . \rho_1e^{-z(t)}\le\sigma(t)\le\rho_2e^{-z(t)},\quad\forall t\in[0,t_f]. ρ1?e?z(t)≤σ(t)≤ρ2?e?z(t),?t∈[0,tf?].
如同在[1]中,我們對不等式(35)使用二階錐近似,它可以很容易地合并到SOCP解框架中,由下式給出
ρ 1 e ? z 0 [ 1 ? ( z ( t ) ? z 0 ( t ) ) + ( z ( t ) ? z 0 ( t ) ) 2 2 ] ≤ σ ( t ) ≤ ρ 2 e ? z 0 [ 1 ? ( z ( t ) ? z 0 ( t ) ) ] , ? t ∈ [ 0 , t f ] \rho_1e^{-z_0}\left[1-(z(t)-z_0(t))+\frac{(z(t)-z_0(t))^2}{2}\right]\le\sigma(t)\le\rho_2e^{-z_0}[1-(z(t)-z_0(t))],\quad\forall t\in[0,t_f] ρ1?e?z0?[1?(z(t)?z0?(t))+2(z(t)?z0?(t))2?]≤σ(t)≤ρ2?e?z0?[1?(z(t)?z0?(t))],?t∈[0,tf?]
其中 z 0 ( t ) = ln ? ( m 0 ? α ρ 2 t ) z_0(t)=\ln(m_0-\alpha\rho_2t) z0?(t)=ln(m0??αρ2?t), m 0 m_0 m0?是飛行器的初始質量。
在[1, 引理2]中,證明了由(36)給出的對(35)中不等式的近似總是產生松弛問題的可行解,并且通常對不等式的兩部分都非常準確,并推導出近似誤差的解析上界。由于本文的重點是控制約束的凸化,我們不會進一步詳細討論這個近似,讀者可參考[1]了解更多細節。
4. 數值例子
本節給出一個數值例子,以演示上一節提出的算法求解行星軟著陸問題。
行星軟著陸中的指向約束確保平移機動不需要飛行器的姿態超出期望的指向錐,這通常會導致性能下降。例如,隨著指向約束的收緊,由于推力指向能力的限制,所需的燃料和飛行時間通常會增加。這個結果將在下面的比較仿真中看到,這些仿真使用了具有以下特性的示例: m 0 = 2000 m_0=2000 m0?=2000 kg, m f = 300 m_f=300 mf?=300 kg, ρ 1 = 0.2 T m a x \rho_1=0.2T_{max} ρ1?=0.2Tmax?, ρ 2 = 0.8 T m a x \rho_2=0.8T_{max} ρ2?=0.8Tmax?, T m a x = 24000 T_{max}=24000 Tmax?=24000 N, α = 5 × 1 0 ? 4 \alpha=5\times 10^{-4} α=5×10?4 s/m,其中 T m a x T_{max} Tmax?是最大推力幅值。推力限制對應于20%和80%的最小和最大節流水平。使用如圖1所示的參考系來表示火星重力和旋轉向量,以及飛行器在動力下降點火時的初始狀態。仿真中飛行器的初始狀態由下式給出
r 0 = ( 2400 , 450 , ? 330 ) T m , r ˙ 0 = ( ? 10 , ? 40 , 10 ) T m/s r_0=(2400,450,-330)^T\text{ m},\quad \dot{r}_0=(-10,-40,10)^T\text{ m/s} r0?=(2400,450,?330)T?m,r˙0?=(?10,?40,10)T?m/s
目標著陸點在 q = ( 0 , 0 , 0 ) T q=(0,0,0)^T q=(0,0,0)T m處,即指導框架的原點。火星引力參數也在同一坐標系中表示如下: g = ( ? 3.71 , 0 , 0 ) T g=(-3.71,0,0)^T g=(?3.71,0,0)T m/s 2 ^2 2, ω = ( 2.53 × 1 0 ? 5 , 0 , 6.62 × 1 0 ? 5 ) T \omega=(2.53\times 10^{-5},0,6.62\times 10^{-5})^T ω=(2.53×10?5,0,6.62×10?5)T rad/s。由于 S ( ω ) n ^ = 0 S(\omega)\hat{n}=0 S(ω)n^=0且 N T S ( ω ) 2 n ^ = 0 N^TS(\omega)^2\hat{n}=0 NTS(ω)2n^=0,我們有 ω ^ = ω \hat{\omega}=\omega ω^=ω。進行了三次仿真,分別針對不同的指向約束限制:(1)無約束;(2) 90°約束;(3) 45°約束。這些仿真的結果如圖6所示。如圖所示,指向角相對于局部垂直方向,即將指向錐 n ^ \hat{n} n^向量與坐標框架的x軸對齊。姿態指向圖表明松弛問題的解確保滿足原問題規定的指向約束。節流圖顯示滿足推力邊界,這表明問題4中對推力幅值和指向的約束松弛對原問題[問題2]仍然有效。該圖進一步提供了一些關于約束收緊時性能權衡的見解。隨著指向限制的收緊,所需的飛行時間和燃料增加,如表I所總結的,這在位置軌跡中也可以看到。45°情況下沿y軸超調目標以滿足指向約束。有趣的是,90°約束路徑采取了更直接的路線。
下一個仿真(見圖7)給出了一個與 X X X邊界有有限次接觸時凸化仍然成立的情況。在這個例子中,我們有以下參數:
r 0 = ( 2400 , 3400 , 0 ) T m , r ˙ 0 = ( ? 40 , 45 , 0 ) T m/s r_0=(2400,3400,0)^T\text{ m},\dot{r}_0=(-40,45,0)^T\text{ m/s} r0?=(2400,3400,0)T?m,r˙0?=(?40,45,0)T?m/s
下滑坡道角度 γ g s = 3 0 ° \gamma_{gs}=30^\circ γgs?=30°, θ = 12 0 ° \theta=120^\circ θ=120°,最大速度為90 m/s。
5. 結論
本文提出了一種針對最優控制理論中的基準問題即軟著陸問題的推力指向約束的無損凸化方法。這擴展和統一了我們之前在該領域的工作,使得算法可以處理推力上下界、推力指向約束、位置和速度約束以及行星自轉。這種凸化使得行星軟著陸問題可以通過凸優化來最優地求解,因此適合機載實現,具有先驗的計算復雜度界和實時執行時間。
附錄
引理2: 考慮以下線性時不變系統:
λ ˙ ( t ) = ? A ( ω ^ ) T λ ( t ) \dot{\lambda}(t)=-A(\hat{\omega})^T\lambda(t) λ˙(t)=?A(ω^)Tλ(t)
y ( t ) = B T λ ( t ) y(t)=B^T\lambda(t) y(t)=BTλ(t)
其中 λ ( t ) ∈ R 6 \lambda(t)\in\mathbb{R}^6 λ(t)∈R6, y ( t ) ∈ R 3 y(t)\in\mathbb{R}^3 y(t)∈R3,且 A ( ω ^ ) A(\hat{\omega}) A(ω^)和 B B B由(2)給出。則以下條件成立;對任意有限區間 [ 0 , t f ] [0,t_f] [0,tf?]:
(i) y y y是解析函數, y ( t ) = 0 y(t)=0 y(t)=0要么在 [ 0 , t f ] [0,t_f] [0,tf?]上成立,要么只在可數個時刻成立;
(ii) 在 [ 0 , t f ] [0,t_f] [0,tf?]中有可數個時刻使得 y ( t ) = α ( t ) n ^ y(t)=\alpha(t)\hat{n} y(t)=α(t)n^對某些 α ( t ) > 0 \alpha(t)>0 α(t)>0成立。
證明:
λ \lambda λ,因此 y y y,是 t t t的解析函數,來自于[28,定理3,第213頁]。因此 y ( t ) = 0 , ? t ∈ [ 0 , t f ? ] y(t)=0,\forall t\in[0,t_f^*] y(t)=0,?t∈[0,tf??],或 y ( t ) = 0 y(t)=0 y(t)=0在 [ 0 , t f ? ] [0,t_f^*] [0,tf??]中只在可數個點處出現([28,命題4.1,第41頁])。這證明了(i)。
為證明條件(ii),假設存在一個區間 [ t 1 , t 2 ] ? [ 0 , t f ? ] [t_1,t_2]\subseteq[0,t_f^*] [t1?,t2?]?[0,tf??],使得 y ( t ) = ? α ( t ) n ^ , ? t ∈ [ t 1 , t 2 ] y(t)=-\alpha(t)\hat{n},\forall t\in[t_1,t_2] y(t)=?α(t)n^,?t∈[t1?,t2?],其中 α ( t ) > 0 \alpha(t)>0 α(t)>0。令 λ = ( λ 1 , λ 2 ) \lambda=(\lambda_1,\lambda_2) λ=(λ1?,λ2?),其中 λ 1 , 2 ∈ R 3 \lambda_{1,2}\in\mathbb{R}^3 λ1,2?∈R3, v 1 : = S ( ω ^ ) T n ^ v_1:=S(\hat{\omega})^T\hat{n} v1?:=S(ω^)Tn^和 v 2 : = S ( ω ^ ) T v 1 v_2:=S(\hat{\omega})^Tv_1 v2?:=S(ω^)Tv1?,該假設意味著以下動力學對 t ∈ [ t 1 , t 2 ] t\in[t_1,t_2] t∈[t1?,t2?]成立: λ ˙ 1 ( t ) = ? α ( t ) v 2 \dot{\lambda}_1(t)=-\alpha(t)v_2 λ˙1?(t)=?α(t)v2?和 λ ˙ 2 ( t ) = ? λ 1 ( t ) ? 2 α ( t ) v 1 \dot{\lambda}_2(t)=-\lambda_1(t)-2\alpha(t)v_1 λ˙2?(t)=?λ1?(t)?2α(t)v1?。
由于在問題4中對 ω ^ \hat{\omega} ω^的構造,我們有 v 1 ≠ 0 v_1\ne 0 v1?=0和 N T v 2 = 0 N^Tv_2=0 NTv2?=0。第一個不等式很容易證明。為了看到第二個,我們需要考慮(21)中的三種情況。首先考慮 S ( ω ) n ^ = 0 , N T S ( ω ) 2 n ^ = 0 S(\omega)\hat{n}=0,N^TS(\omega)^2\hat{n}=0 S(ω)n^=0,NTS(ω)2n^=0的情況。那么 ω = ω ^ \omega=\hat{\omega} ω=ω^且 N T v 2 = N T S ( ω ) 2 n ^ = 0 N^Tv_2=N^TS(\omega)^2\hat{n}=0 NTv2?=NTS(ω)2n^=0。第二,考慮 S ( ω ) n ^ = 0 S(\omega)\hat{n}=0 S(ω)n^=0的情況,即 ω = a n ^ \omega=a\hat{n} ω=an^對某些 a ∈ R a\in\mathbb{R} a∈R。那么 ω ^ = a n ^ + ε n ^ ⊥ \hat{\omega}=a\hat{n}+\varepsilon\hat{n}^\perp ω^=an^+εn^⊥。然后 v 2 = S ( ω ^ ) T n ^ = S ( ω ^ ) 2 n ^ = S ( a n ^ + ε n ^ ⊥ ) S ( ε n ^ ⊥ ) n ^ v_2=S(\hat{\omega})^T\hat{n}=S(\hat{\omega})^2\hat{n}=S(a\hat{n}+\varepsilon\hat{n}^\perp)S(\varepsilon\hat{n}^\perp)\hat{n} v2?=S(ω^)Tn^=S(ω^)2n^=S(an^+εn^⊥)S(εn^⊥)n^。由于 S ( ε n ^ ⊥ ) n ^ S(\varepsilon\hat{n}^\perp)\hat{n} S(εn^⊥)n^與 a n ^ a\hat{n} an^和 ε n ^ ⊥ \varepsilon\hat{n}^\perp εn^⊥都正交,它與它們的和正交,該和非零,因此 v 2 = 0 v_2=0 v2?=0且它在 n ^ T \hat{n}^T n^T的零空間中。因此 N T S ( ω ^ ) T n ^ = 0 N^TS(\hat{\omega})^T\hat{n}=0 NTS(ω^)Tn^=0。最后,考慮 S ( ω ) n ^ = 0 S(\omega)\hat{n}=0 S(ω)n^=0但 N T S ( ω ) T n ^ ≠ 0 N^TS(\omega)^T\hat{n}\ne 0 NTS(ω)Tn^=0的情況。那么我們有 S ( ω ^ ) 2 n ^ = S ( ω ) 2 n ^ + ε 2 S ( n ^ ) S ( ω ) n ^ S(\hat{\omega})^2\hat{n}=S(\omega)^2\hat{n}+\varepsilon^2 S(\hat{n})S(\omega)\hat{n} S(ω^)2n^=S(ω)2n^+ε2S(n^)S(ω)n^。這在 n ^ T \hat{n}^T n^T的零空間中有一個非零分量,因為 ε 2 S ( n ^ ) S ( ω ) n ^ \varepsilon^2 S(\hat{n})S(\omega)\hat{n} ε2S(n^)S(ω)n^非零且在 n ^ T \hat{n}^T n^T的零空間中,而 S ( ω ^ ) 2 n ^ S(\hat{\omega})^2\hat{n} S(ω^)2n^與 n ^ T \hat{n}^T n^T的零空間正交。因此,我們有 N T v 2 = 0 N^Tv_2=0 NTv2?=0。現在注意到 v 1 T n ^ = 0 v_1^T\hat{n}=0 v1T?n^=0。令 v 3 v_3 v3?為非零向量使得 v 3 T n ^ = v 3 T v 1 = 0 v_3^T\hat{n}=v_3^Tv_1=0 v3T?n^=v3T?v1?=0。因此 v 1 v_1 v1?和 v 3 v_3 v3?張成 n ^ \hat{n} n^的零空間。這意味著我們可以寫 v 2 = c 1 n ^ + c 2 v 3 + c 3 v 1 v_2=c_1\hat{n}+c_2v_3+c_3v_1 v2?=c1?n^+c2?v3?+c3?v1?對某些標量 c 1 c_1 c1?到 c 3 c_3 c3?。由于 v 2 T v 1 = 0 v_2^Tv_1=0 v2T?v1?=0和 N T v 2 = 0 N^Tv_2=0 NTv2?=0我們知道 c 3 = 0 c_3=0 c3?=0和 c 2 = 0 c_2=0 c2?=0。
觀察到
λ ˙ 2 ( t ) = λ 1 ( t 1 ) ? ∫ t 1 t α ( s ) d s v 2 ? 2 α ( t ) v 1 \dot{\lambda}_2(t)=\lambda_1(t_1)-\int_{t_1}^t\alpha(s)ds\,v_2-2\alpha(t)v_1 λ˙2?(t)=λ1?(t1?)?∫t1?t?α(s)dsv2??2α(t)v1?
這意味著
λ ˙ 2 ( t ) = λ 1 ( t 1 ) + g 1 ( t ) n ^ + g 2 ( t ) v 3 ? 2 α ( t ) v 1 \dot{\lambda}_2(t)=\lambda_1(t_1)+g_1(t)\hat{n}+g_2(t)v_3-2\alpha(t)v_1 λ˙2?(t)=λ1?(t1?)+g1?(t)n^+g2?(t)v3??2α(t)v1?
其中 g 1 ( t ) = c 1 ∫ t 1 t α ( s ) d s g_1(t)=c_1\int_{t_1}^t\alpha(s)ds g1?(t)=c1?∫t1?t?α(s)ds且 g 2 ( t ) = c 2 ∫ t 1 t α ( s ) d s = 0 g_2(t)=c_2\int_{t_1}^t\alpha(s)ds=0 g2?(t)=c2?∫t1?t?α(s)ds=0對 t ∈ ( t 1 , t 2 ) t\in(t_1,t_2) t∈(t1?,t2?)。由于 v 1 , v 3 v_1,v_3 v1?,v3?和 n ^ \hat{n} n^形成一組非零向量的正交集,如果 α \alpha α是關于時間的非常數函數,那么 N T λ ˙ 2 ( t ) = 0 N^T\dot{\lambda}_2(t)=0 NTλ˙2?(t)=0幾乎處處在 [ t 1 , t 2 ] [t_1,t_2] [t1?,t2?]上。如果 α \alpha α是常數,那么 g 2 g_2 g2?是關于時間的非常數線性函數(因為 c 2 ≠ 0 c_2\ne 0 c2?=0),因此 N T λ ˙ 2 ( t ) = 0 N^T\dot{\lambda}_2(t)=0 NTλ˙2?(t)=0幾乎處處在 [ t 1 , t 2 ] [t_1,t_2] [t1?,t2?]上。因此, N T λ ˙ 2 ( t ) = 0 N^T\dot{\lambda}_2(t)=0 NTλ˙2?(t)=0幾乎處處在 [ t 1 , t 2 ] [t_1,t_2] [t1?,t2?]上,這意味著存在 [ t 1 , t 2 ] [t_1,t_2] [t1?,t2?]中的一個開的有限區間使得 λ 2 \lambda_2 λ2?將離開由 n ^ \hat{n} n^定義的子空間,這與假設 λ 2 ( t ) = ? α ( t ) n ^ \lambda_2(t)=-\alpha(t)\hat{n} λ2?(t)=?α(t)n^在 [ t 1 , t 2 ] [t_1,t_2] [t1?,t2?]上成立相矛盾。因此不存在有限區間使得 y ( t ) = λ 2 ( t ) = ? α ( t ) n ^ y(t)=\lambda_2(t)=-\alpha(t)\hat{n} y(t)=λ2?(t)=?α(t)n^,即 κ ( t ) : = N T λ ˙ 2 ( t ) ≠ 0 \kappa(t):=N^T\dot{\lambda}_2(t)\ne 0 κ(t):=NTλ˙2?(t)=0。由于 κ \kappa κ是時間的解析函數(因為 λ \lambda λ是時間的解析函數),要么它有可數個零點,要么存在它為零的時間區間[28]。由于我們證明了后者是不可能的,那么它必須有可數個零點。這就完成了證明。
下面的引理在使用Pontryagin極大值原理表明最優控制發生在凸可行控制集的極點處起著重要作用。
引理3:
以下優化問題的最優解也是可行解集 U ( σ ) U(\sigma) U(σ)的極點: max ? T c y T T c \max_{T_c}y^TT_c maxTc??yTTc? 滿足 T c ∈ U ( σ ) T_c\in U(\sigma) Tc?∈U(σ),其中 U ( σ ) : = { T c : ∥ T c ∥ ≤ σ , n ^ T T c ≥ cos ? θ σ } U(\sigma):=\{T_c:\|T_c\|\le\sigma,\hat{n}^TT_c\ge\cos\theta\,\sigma\} U(σ):={Tc?:∥Tc?∥≤σ,n^TTc?≥cosθσ},且 y ≠ 0 y\ne 0 y=0和 y = ? α n ^ y=-\alpha\hat{n} y=?αn^對任意 α > 0 \alpha>0 α>0。因此最優解 T c ? T_c^* Tc??滿足 ∥ T c ? ∥ = σ \|T_c^*\|=\sigma ∥Tc??∥=σ。
證明:
由于優化問題的目標函數是線性的,因此是凸的,最大化問題導致 U U U的邊界 ? U \partial U ?U上的最優解[29]。 ? U \partial U ?U有一部分是半徑為 σ \sigma σ的球面,一部分是單位法向量為 n ^ \hat{n} n^的超平面。邊界的所有極點都在球面上:如果邊界上的任何點不在球面上,那么它在超平面上,在邊界上另外兩點之間的線段上。這意味著 U U U的極點滿足 ∥ T c ∥ = σ \|T_c\|=\sigma ∥Tc?∥=σ。
注意到 y = c 1 n ^ + c 2 n ^ ⊥ y=c_1\hat{n}+c_2\hat{n}^\perp y=c1?n^+c2?n^⊥對某些 c 2 ≠ 0 c_2\ne 0 c2?=0和與 n ^ \hat{n} n^正交的單位向量 n ^ ⊥ \hat{n}^\perp n^⊥。這意味著 T c ? = a 1 n ^ + a 2 n ^ ⊥ T_c^*=a_1\hat{n}+a_2\hat{n}^\perp Tc??=a1?n^+a2?n^⊥,其中 a 2 ≠ 0 a_2\ne 0 a2?=0且 a 1 2 + a 2 2 ≤ σ 2 a_1^2+a_2^2\le\sigma^2 a12?+a22?≤σ2。假設 n ^ T T c ? = cos ? θ σ \hat{n}^TT_c^*=\cos\theta\,\sigma n^TTc??=cosθσ,這意味著 a 1 = cos ? θ σ a_1=\cos\theta\,\sigma a1?=cosθσ。由于 y T T c ? = c 1 a 1 + c 2 a 2 = c 1 cos ? θ σ + c 2 a 2 y^TT_c^*=c_1a_1+c_2a_2=c_1\cos\theta\,\sigma+c_2a_2 yTTc??=c1?a1?+c2?a2?=c1?cosθσ+c2?a2?,如果 c 2 c_2 c2?和 a 2 a_2 a2?同號且 a 1 2 + a 2 2 = σ 2 a_1^2+a_2^2=\sigma^2 a12?+a22?=σ2,即 ∥ T c ? ∥ = σ \|T_c^*\|=\sigma ∥Tc??∥=σ,那么代價最大。
因此最優解必須是 U U U的極點。
總結
這篇文章研究了行星軟著陸問題,該問題可以表述為一個有狀態和控制約束的有限時域最優控制問題。求解這個問題的主要困難在于控制輸入上存在非凸約束。本文的主要貢獻是引入了對這些控制約束的凸化方法,并證明它是無損的,即可以通過求解凸化松弛問題得到原問題的最優解。
論文首先給出了行星軟著陸問題的數學描述,包括動力學方程、狀態約束(如下滑道約束)和控制約束(如推力幅值和指向約束)。然后提出了兩個優化問題的形式:最小著陸誤差問題和最小燃料問題。接下來引入了這兩個問題的凸松弛形式,其中非凸控制約束被凸約束代替。
論文的核心結果是定理1,它說明了在一定條件下,可以通過求解凸松弛問題得到原問題的最優解。證明過程用到了Pontryagin極大值原理,以及附錄中給出的兩個引理。引理2給出了余態方程的一些性質,引理3說明了最優控制發生在可行控制集的極點處。
為了方便求解,論文還引入了變量替換以消除動力學方程中的非線性,并使用了二階錐規劃對約束進行近似。最后,通過數值算例驗證了所提出方法的有效性。
我給出一個簡單的代碼示例,演示如何使用cvxpy庫在Python中求解凸化后的軟著陸問題。
import cvxpy as cp
import numpy as np# 問題參數
n = 100 # 時間步數
dt = 0.1 # 時間步長
m0 = 1000 # 初始質量
mf = 300 # 最終質量
g = 9.8 # 重力加速度
rho_1 = 1000 # 最小推力
rho_2 = 5000 # 最大推力
theta = np.pi/3 # 推力指向約束的半角
r0 = np.array([0, 0, 1000]) # 初始位置
v0 = np.array([0, 0, -50]) # 初始速度
rf = np.array([100, 100, 0]) # 目標位置# 決策變量
r = cp.Variable((n+1, 3)) # 位置
v = cp.Variable((n+1, 3)) # 速度
u = cp.Variable((n, 3)) # 推力加速度
sigma = cp.Variable(n) # 推力大小# 初始條件
constr = [r[0] == r0, v[0] == v0]# 動力學方程
for t in range(n):dr = v[t] * dtdv = (u[t] - g*np.array([0,0,1])) * dtconstr.append(r[t+1] == r[t] + dr)constr.append(v[t+1] == v[t] + dv) # 終端狀態約束
constr.append(r[-1] == rf)
constr.append(v[-1] == 0)# 控制約束
for t in range(n):constr.append(cp.norm(u[t]) <= sigma[t])constr.append(u[t] @ np.array([0,0,1]) >= cp.norm(u[t])*np.cos(theta))constr.append(rho_1/m0 <= sigma[t])constr.append(sigma[t] <= rho_2/m0)# 目標函數
objective = cp.Minimize(cp.sum(sigma)*dt)# 求解優化問題
prob = cp.Problem(objective, constr)
result = prob.solve()# 輸出結果
print(f"最優燃料消耗: {result*m0}")
-
首先定義了問題的參數,如時間步數、初始狀態、目標位置等。
-
然后定義了決策變量,包括各時刻的位置r、速度v、推力加速度u和推力大小sigma。
-
接下來是約束條件,包括初始條件、動力學方程約束、終端狀態約束和控制約束。其中控制約束對應于論文中的凸化松弛形式。
-
目標函數是最小化總推力的積分,對應最小燃料消耗。
-
最后,我們構建優化問題,調用求解器,輸出最優目標值。
這個示例模型與論文中的略有不同,例如使用了推力加速度而不是直接的推力作為控制變量,對質量的變化也進行了簡化。但基本的思路是一致的,即將原問題轉化為凸優化問題求解。
這篇論文中的參考文獻列表:
[1] B. A??kme?se and S. R. Ploen, “Convex programming approach to powered descent guidance for mars landing,” Journal of Guidance, Control, and Dynamics, vol. 30, no. 5, pp. 1353–1366, 2007.
[2] L. Blackmore, B. A??kme?se, and D. P. Scharf, “Minimum-landing-error powered-descent guidance for mars landing using convex optimization,” Journal of Guidance, Control, and Dynamics, vol. 33, no. 4, pp. 1161–1171, 2010.
[3] J. S. Meditch, “On the problem of optimal thrust programming for a lunar soft landing,” IEEE Transactions on Automatic Control, vol. 9, no. 4, pp. 477–484, 1964.
[4] A. Miele, “The calculus of variations in applied aerodynamics and flight mechanics,” in Optimization Techniques, G. Leitmann, Ed. Elsevier, 1962, pp. 99–170.
[5] A. R. Klumpp, “Apollo lunar descent guidance,” Automatica, vol. 10, no. 2, pp. 133–146, 1974.
[6] U. Topcu, J. Casoliva, and K. D. Mease, “Minimum-fuel powered descent for mars pinpoint landing,” Journal of Spacecraft and Rockets, vol. 44, no. 2, pp. 324–331, 2007.
[7] R. Sostaric and J. Rea, “Powered descent guidance methods for the moon and mars,” in AIAA Guidance, Navigation, and Control Conference and Exhibit, 2005, p. 6287.
[8] F. Najson and K. D. Mease, “A computationally non-expensive guidance algorithm for fuel efficient soft landing,” in AIAA Guidance, Navigation, and Control Conference and Exhibit, 2006, p. 6438.
[9] B. A. Steinfeldt, M. J. Grant, D. A. Matz, R. D. Braun, and G. H. Barton, “Guidance, navigation, and control system performance trades for mars pinpoint landing,” Journal of Spacecraft and Rockets, vol. 47, no. 1, pp. 188–198, 2010.
[10] S. Boyd and L. Vandenberghe, Convex Optimization. Cambridge University Press, 2004.
[11] Y. Nesterov and A. Nemirovskii, Interior-Point Polynomial Algorithms in Convex Programming. SIAM, 1994.
[12] J. Peng, C. Roos, and T. Terlaky, Self-Regularity: A New Paradigm for Primal-Dual Interior-Point Algorithms. Princeton University Press, 2002.
[13] J. F. Sturm, “Using SeDuMi 1.02, a MATLAB toolbox for optimization over symmetric cones,” Optimization methods and software, vol. 11, no. 1-4, pp. 625–653, 1999.
[14] C. D’Souza, “An optimal guidance law for planetary landing,” in AIAA Guidance, Navigation, and Control Conference, 1997, pp. 1376–1381.
[15] S. R. Ploen, B. Acikmese, and A. Wolf, “A comparison of powered descent guidance laws for mars pinpoint landing,” in AIAA/AAS Astrodynamics Specialist Conference and Exhibit, 2006, p. 6676.
[16] B. A. Steinfeldt, M. J. Grant, D. A. Matz, R. D. Braun, and G. H. Barton, “Guidance, navigation, and control technology system trades for mars pinpoint landing,” in AIAA Guidance, Navigation and Control Conference and Exhibit, 2008, p. 6216.
[17] J. M. Carson III, B. A??kme?se, and L. Blackmore, “Lossless convexification of powered-descent guidance with non-convex thrust bound and pointing constraints,” in 2011 American Control Conference, IEEE, 2011, pp. 2651–2656.
[18] L. D. Berkovitz, Optimal Control Theory, ser. Applied Mathematical Sciences. Springer New York, 1974, vol. 12.
[19] K.-C. Toh, M. J. Todd, and R. H. Tütüncü, “SDPT3—a matlab software package for semidefinite programming, version 1.3,” Optimization methods and software, vol. 11, no. 1-4, pp. 545–581, 1999.
[20] J. Mattingley and S. Boyd, “Real-time convex optimization in signal processing,” IEEE Signal Processing Magazine, vol. 27, no. 3, pp. 50–61, 2010.
[21] ——, “Automatic code generation for real-time convex optimization,” in Convex Optimization in Signal Processing and Communications, D. P. Palomar and Y. C. Eldar, Eds. Cambridge University Press, 2010.
[22] Y. Wang and S. Boyd, “Fast model predictive control using online optimization,” IEEE Transactions on Control Systems Technology, vol. 18, no. 2, pp. 267–278, 2009.
[23] Y. Kim and M. Mesbahi, “Quadratically constrained attitude control via semidefinite programming,” IEEE Transactions on Automatic Control, vol. 49, no. 5, pp. 731–735, 2004.
[24] G. P. Sutton and O. Biblarz, Rocket Propulsion Elements. John Wiley & Sons, 2010.
[25] M. Tamiz, D. Jones, and C. Romero, “Goal programming for decision making: An overview of the current state-of-the-art,” European journal of operational research, vol. 111, no. 3, pp. 569–581, 1998.
[26] L. S. Pontryagin, Mathematical Theory of Optimal Processes. Routledge, 2018.
[27] B. A??kme?se and L. Blackmore, “Lossless convexification of a class of optimal control problems with non-convex control constraints,” Automatica, vol. 47, no. 2, pp. 341–347, 2011.
[28] H. Cartan, Elementary Theory of Analytic Functions of One or Several Complex Variables. Courier Corporation, 1995.
[29] L. D. Berkovitz, Convexity and Optimization in R n R^n Rn. John Wiley & Sons, 2003.