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:


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:



$$\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!


There are many ways to deal with this problem. You may use a multiple shooting method or a sort of orthogonal collocation method (see here for a nice overview).

In any case, you would need a relatively good initial guess at the adjoint variables, which is not easy at all. I suggest that you first try to solve the problem as a direct optimization problem by minimizing $J$ while satisfying the end-point constraints. Then you can try to recover the adjoint variables from the suboptimal control $\bar\beta(t)$ and use them as the initial guess.

In any case, this might require quite a lot of programming so, maybe you can think of using some existing software (see, e.g., here or here). Some of the packages are free of charge when used for academic purposes.

  • 1
    $\begingroup$ Thanks for your suggestions Dmitry. I'm starting to understand that although solving the problem may seem straightforward, it isn't at all. Yesterday I found a trajectory optimization library for MATLAB, OptimTraj. I solved the problem above with direct collocation, using both the trapezoid method and Hermite-Simpson method. $\endgroup$ – woeterb Oct 20 '17 at 14:39

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service, privacy policy and cookie policy

Not the answer you're looking for? Browse other questions tagged or ask your own question.