This question was also posted at the Space Exploration StackExchange, but I thought it might have a better chance of receiving an answer here.
For an in-plane non-impulsive orbital maneuver, I'd like to find the thrust-direction history $\beta(t)$ to minimize the time needed to transfer a spacecraft from a specified initial state to a specified terminal state, subject to the following state equations:
$$\dot{\textbf{x}}(\textbf{x},\beta,t)=\begin{bmatrix}\dot{r} & \dot{V}_r & \dot{V}_{\theta}\end{bmatrix}$$ $$\dot{r}=v_r$$ $$\dot{v}_r=\frac{v_{\theta}^2}{r}-\frac{\mu}{r^2}+\frac{T\sin{\beta}}{m_0-|\dot{m}|t}$$ $$\dot{v}_{\theta}=-\frac{v_{r}v_{\theta}}{r}+\frac{T\cos{\beta}}{m_0-|\dot{m}|t}$$
where $\mu$, $T$, $m_0$ and $\dot{m}$ are constants, with state variables:
$$\textbf{x}(t)=\begin{bmatrix}r & V_r & V_{\theta}\end{bmatrix}$$
and initial and terminal conditions:
$$\begin{align*} \textbf{x}(t_0) &=\begin{bmatrix}1745100 & 0 & 2027\end{bmatrix}\\ \textbf{x}(t_f) &=\begin{bmatrix}1765100 & 0 & 0\end{bmatrix}\\ \end{align*}$$
i.e. we deal with a two-point boundary value problem where a spacecraft makes a transfer from an initial position to a higher position, while decreasing its initial velocity to zero terminal velocity. The problem is cast as an optimal control problem for which I have formulated the following minimum-time cost functional:
$$J=\min_{\beta}\int_{t_0}^{t_f}L(\textbf{x},\beta,t)dt=\min_{\beta}\int_{t_0}^{t_f}(1)dt$$
where $L$ denotes the Lagrangian. The Hamiltonian is, therefore,
$$\begin{align*} H(\textbf{x},\beta,t) &= L(\textbf{x},\beta,t)+\lambda^T\dot{\textbf{x}}(\textbf{x},\beta,t) \\ &= 1+\lambda_{r}\dot{r}+\lambda_{v_r}\dot{v}_r+\lambda_{v_{\theta}}\dot{v_{\theta}} \\ &= 1+\lambda_{r}\cdot v_r+\lambda_{v_r}\cdot (\frac{v_{\theta}^2}{r}-\frac{\mu}{r^2}+\frac{T\sin{\beta}}{m_0-|\dot{m}|t})+\lambda_{v_{\theta}}\cdot (-\frac{v_{r}v_{\theta}}{r}+\frac{T\cos{\beta}}{m_0-|\dot{m}|t}) \end{align*}$$
where $\mathbf{\lambda}(t)$ is the Lagrange multiplier vector, whose elements are the costate variables. The $\lambda$-dynamics that capture the behavior of the Lagrange multipliers, the costate equations, are given by:
$$\begin{align*} \dot{\lambda}_{r} =-\frac{\delta H}{\delta r} &= -\lambda_{v_r} \cdot (-\frac{v_{\theta}^2}{r^2}+\frac{2\mu}{r^3}) - \lambda_{v_{\theta}} \cdot (\frac{v_{r}v_{\theta}}{r^2})\\ \dot{\lambda}_{v_r} =-\frac{\delta H}{\delta v_r} &= -\lambda_{r} + \lambda_{v_{\theta}} \cdot (\frac{v_{\theta}}{r})\\ \dot{\lambda}_{v_{\theta}} =-\frac{\delta H}{\delta v_{\theta}} &= -\lambda_{v_r} \cdot (\frac{2v_{\theta}}{r}) + \lambda_{v_{\theta}} \cdot (\frac{v_r}{r})\\ \end{align*}$$
The partial of $H$ with respect to $\beta$ must equal zero, so
$$\frac{\delta H}{\delta \beta}=\lambda_{v_r} \cdot (\frac{T\cos{\beta}}{m_0-|\dot{m}|t}) - \lambda_{v_{\theta}} \cdot (\frac{T\sin{\beta}}{m_0-|\dot{m}|t})=0$$
leading to the following control law:
$$\tan{\beta}=\frac{-\lambda_{v_r}}{-\lambda_{v_{\theta}}}$$
or,
$$\sin{\beta}=\frac{-\lambda_{v_r}}{\sqrt{\lambda_{v_r}^2+\lambda_{v_{\theta}}^2}} \qquad \cos{\beta}=\frac{-\lambda_{v_{\theta}}}{\sqrt{\lambda_{v_r}^2+\lambda_{v_{\theta}}^2}}$$
Note: the signs in the ratio are negative because we are looking for a minimum.
This is the first time I'm working with an optimal control problem, and I have some trouble understanding how to proceed from here and solve the problem with MATLAB. Is there perhaps someone out here who knows how to do this and could explain what the next steps would be?
Thanks a lot!