4
$\begingroup$

EDIT October 15, 2017: the problem in this question is cast as an optimal control problem and reformulated here.

The spacecraft I am designing is flying in an elliptical orbit (8 x 3000 km) around the Moon. At periapsis, the non-impulsive bi-propellant retromotor mounted upon the spacecraft will cancel the orbital velocity while at the same time increase the altitude from 8 to 28 km (the spacecraft will then start its free-fall phase to the surface, but that's not part of this problem). The Moon is assumed to be a perfect sphere with $r_{Moon}=1737.100$ $km$ and non-rotating.

The state vector consists of 7 components (position and velocity vectors, plus the mass). The initial state vector is defined as:

$$ \textbf{y}_0 = \begin{bmatrix}x \\ y \\ z \\ \dot{x} \\ \dot{y} \\ \dot{z} \\ m \end{bmatrix}_{t=0} = \begin{bmatrix}1745100 \enspace m\\ 0 \\ 0 \\ 0 \\ 2026.813 \enspace m/s\\ 0 \\ 70 \enspace kg \end{bmatrix} $$

The dynamics are given by:

$$ \ddot{\textbf{r}}=-\mu\frac{\textbf{r}}{r^{3}}+\frac{\textbf{T}}{m} $$

So:

$$ \dot{\textbf{y}} = \begin{bmatrix}\dot{x} \\ \dot{y} \\ \dot{z} \\ \ddot{x} \\ \ddot{y} \\ \ddot{z} \\ \dot{m} \end{bmatrix} \quad \textbf{f}(t,\textbf{y}) = \begin{bmatrix}y_4\\ y_5 \\ y_6 \\ -\mu\frac{y_1}{r^{3}}+\frac{T}{y_7}\cdot (k_x ?) \\ -\mu\frac{y_2}{r^{3}}+\frac{T}{y_7}\cdot (k_y?)\\ -\mu\frac{y_3}{r^{3}}+\frac{T}{y_7}\cdot (k_z?) \\ -\frac{T}{I_{sp}g_{0}} \end{bmatrix} $$

With $r=\sqrt{y_1^{2}+y_2^{2}+y_3^{2}}$ and constants:

$\mu = 4.90487E(+12) \enspace m^3/s^2$; $T = 200 \enspace N$; $I_{sp} = 300 \enspace s$; $g_{0} = 9.81 \enspace m/s^2$

The $k$'s in the set of equations should determine the thrust magnitude (and sign) in that particular direction ($x$, $y$ or $z$) at each integration step: $\textbf{k}$ is the unit vector of the thrust vector. The solution I'm seeking should be in the form of a vector set containing thrust vectors at each integration step of the breaking/climbing phase.

The final state vector should be:

$$ \textbf{y}_n = \begin{bmatrix}x_n \\ y_n \\ z_n \\ 0 \\ 0 \\ 0 \\ m_n \end{bmatrix} $$

Where $r_n=\sqrt{x_n^{2}+y_n^{2}+z_n^{2}}=1737100+28000=1765100 \enspace m$ and $m_n$ must be maximum.

Any ideas on how to find an optimal solution (read: optimal thrust-vectoring strategy, final mass $m_n$ should be maximum) for this boundary value problem?

Thanks a lot!

Improve this question
$\endgroup$
7
  • 1
    $\begingroup$ Haven't done this in several decades, but there are some general approaches to solving BVPs. One might be to choose a functional form for the thrust vector. I assume you'll use constant magnitude (maximum) thrust, so all you have to do is invent a function for the direction of the thrust within the orbital plane as a function of time; $\theta (t)$, right? It could be say a 3rd or 5th degree polynomial. I'd recommend some kind of Monte Carlo to find regions in coefficient space that roughly work, then you could use steepest descent in each. $\endgroup$
    – uhoh
    Sep 21, 2017 at 14:15
  • $\begingroup$ If your orbital plane is in the $x-y$ plane, then why not drop $z$ completely? $\endgroup$
    – uhoh
    Sep 21, 2017 at 14:17
  • $\begingroup$ Thanks for your reply! For this particular case you could drop the $z$ indeed. But I would like to create an optimizer that could also deal with orbits other than those in the $x-y$ plane. Unfortunately I do not have a lot of experience with optimization techniques. So what you suggest is to assume that the thrust angle follows a 3rd or 5th degree polynomial with unknown coefficients? These coefficient can then be guessed through Monte Carlo and steepest descent? $\endgroup$
    – woeterb
    Sep 22, 2017 at 10:57
  • 1
    $\begingroup$ That's right. You could make $\theta$ relative to the $\hat{y}$ direction or relative to the local horizontal, I'm not sure which one is easier. I can't say that this is necessarily the easiest way to solve it, or else I'd post it as an answer. I'll try to do a few tests in a day or two. You might start by finding a fixed thrust angle that gets you your goal. It may not be optimal, but then you could use it as a starting point and then add one new term at a time to the polynomial, re-optimizing after each step. One parameter means you don't have to use Monte Carlo to hunt 5D parameter space. $\endgroup$
    – uhoh
    Sep 22, 2017 at 11:51
  • $\begingroup$ Thanks for your efforts! In this situation I guess it's best to make $\theta$ relative to the azimuthal component of the velocity vector (i.e. to the local horizontal), since the direction of the thrust vector will depend on the direction of the velocity vector. Unfortunately, there's is no fixed thrust angle that gets me my goal. I do however get the behaviour I'd like to see in the trajectory: the orbital velocity is cancelled while a distance from the surface is maintained (although the final altitude is too large). I'm still seeking a way to determine the optimal thrust angle strategy. $\endgroup$
    – woeterb
    Sep 25, 2017 at 13:37

0

Reset to default

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service and acknowledge that you have read and understand our privacy policy and code of conduct.

Browse other questions tagged or ask your own question.