Reusable Launch Vehicle Ascent

The reusable launch vehicle ascent optimal control problem is stated as follows. The trajectory is divided into four phases: phase 1 runs from staging to main engine cutoff, phases 2 and 4 are powered orbital maneuvering phases, and phase 3 is a coast phase. Maximize the final mass, equivalently minimize \[J=-\frac{m^{(4)}(t_f^{(4)})}{10^5}.\] In each phase, the state is \[\mathbf{x}^{(i)}=(r,\theta,\phi,v,\gamma,\psi,m),\qquad i=1,2,3,4,\] where \(r\) is geocentric radius, \(\theta\) is longitude, \(\phi\) is geocentric latitude, \(v\) is Earth-relative speed, \(\gamma\) is flight-path angle, \(\psi\) is heading angle, and \(m\) is mass. In phases 1, 2, and 4, the control is \[\mathbf{u}^{(i)}=(\alpha,w_1,w_2).\] The differential equations in phases 1, 2, and 4 are \[\begin{array}{lcl}\dot r &=& v\sin\gamma, \\ \dot\theta &=& \displaystyle \frac{v\cos\gamma\cos\psi}{r\cos\phi}, \\ \dot\phi &=& \displaystyle \frac{v\cos\gamma\sin\psi}{r}, \\ \dot v &=& \displaystyle \frac{F_T}{m}-g\sin\gamma+\omega^2 r\cos\phi\left(\sin\gamma\cos\phi-\cos\gamma\sin\phi\sin\psi\right), \\ \dot\gamma &=& \displaystyle \frac{1}{v}\left[\frac{F_N w_1}{m}-\left(g-\frac{v^2}{r}\right)\cos\gamma+2\omega v\cos\phi\cos\psi+\omega^2 r\cos\phi\left(\cos\gamma\cos\phi+\sin\gamma\sin\phi\sin\psi\right)\right], \\ \dot\psi &=& \displaystyle \frac{1}{v}\left[\frac{F_N w_2}{m\cos\gamma}-\frac{v^2\cos\gamma\cos\psi\tan\phi}{r}+2\omega v\left(\tan\gamma\cos\phi\sin\psi-\sin\phi\right)-\frac{\omega^2 r\sin\phi\cos\phi\cos\psi}{\cos\gamma}\right], \\ \dot m &=& \displaystyle -\frac{T}{I_{sp}g_0},\end{array}\] where \[g=\frac{\mu}{r^2}.\] In the coast phase, phase 3, the differential equations are obtained by setting \[F_T=0,\qquad F_N=0,\qquad \dot m=0.\] The powered phases satisfy the path constraint \[w_1^2+w_2^2=1,\qquad i=1,2,4.\] The phase linkage event constraints are \[\begin{array}{lcl}\mathbf{x}^{(2)}(t_0^{(2)})-\mathbf{x}^{(1)}(t_f^{(1)}) &=& (0,0,0,0,0,0,-m_{\mathrm{ET,empty}}), \\ t_0^{(2)}-t_f^{(1)} &=& 0, \\ \mathbf{x}^{(3)}(t_0^{(3)})-\mathbf{x}^{(2)}(t_f^{(2)}) &=& \mathbf{0}, \\ t_0^{(3)}-t_f^{(2)} &=& 0, \\ \mathbf{x}^{(4)}(t_0^{(4)})-\mathbf{x}^{(3)}(t_f^{(3)}) &=& \mathbf{0}, \\ t_0^{(4)}-t_f^{(3)} &=& 0.\end{array}\] The initial boundary conditions at staging are \[\begin{array}{lcl} t_0^{(1)} &=& 126.1~\mathrm{s}, \\ h^{(1)}(t_0^{(1)}) &=& 46447.588~\mathrm{m}, \\ \theta^{(1)}(t_0^{(1)}) &=& -120.74~\mathrm{deg}, \\ \phi^{(1)}(t_0^{(1)}) &=& 34.1~\mathrm{deg}, \\ v^{(1)}(t_0^{(1)}) &=& 1384.7~\mathrm{m}\,\mathrm{s}^{-1}, \\ \gamma^{(1)}(t_0^{(1)}) &=& 26.409~\mathrm{deg}, \\ \psi^{(1)}(t_0^{(1)}) &=& 259~\mathrm{deg}, \\ m^{(1)}(t_0^{(1)}) &=& 675702.8~\mathrm{kg}.\end{array}\] The main-engine-cutoff event constraints are \[h^{(1)}(t_f^{(1)})=105564.1~\mathrm{m},\qquad v_I^{(1)}(t_f^{(1)})=7734.0~\mathrm{m}\,\mathrm{s}^{-1},\qquad \gamma_I^{(1)}(t_f^{(1)})=0.65~\mathrm{deg},\qquad i^{(1)}(t_f^{(1)})=98~\mathrm{deg}.\] The final-orbit event constraints are \[h^{(4)}(t_f^{(4)})=203720~\mathrm{m},\qquad v_I^{(4)}(t_f^{(4)})=\sqrt{\frac{\mu}{r^{(4)}(t_f^{(4)})}},\qquad \gamma_I^{(4)}(t_f^{(4)})=0,\qquad i^{(4)}(t_f^{(4)})=98~\mathrm{deg}.\] The inertial event quantities are computed from \[\mathbf{r}=(r\cos\phi\cos\theta,\;r\cos\phi\sin\theta,\;r\sin\phi),\qquad \theta=\theta_{\mathrm{lon}}+\omega(t-t_0^{(1)}),\] \[\mathbf{v}_I=\mathbf{v}_{\mathrm{rel}}+\boldsymbol{\omega}\times\mathbf{r},\qquad \gamma_I=\sin^{-1}\left(\frac{\mathbf{r}\cdot\mathbf{v}_I}{\|\mathbf{r}\|\,\|\mathbf{v}_I\|}\right),\qquad i=\cos^{-1}\left(\frac{h_3}{\|\mathbf{h}\|}\right),\qquad \mathbf{h}=\mathbf{r}\times\mathbf{v}_I.\] The solution to the reusable launch vehicle ascent problem using GPOPS-II is shown in the figures below.