1
$\begingroup$

I am trying to implement (in Python for now) low thrust orbit propagation for spacecraft using universal variables. For a given central body with the gravitational parameter $\mu$ and an orbit with the semi-major axis $a$ and the initial position $\vec{r}_0$ and velocity $\vec{v}_0$ at $t = t_0$ the position for a given time is given by: $$ \vec{r} = \vec{r}_0 \ f(s) + \vec{v}_0 \ g(s) $$ and the velocity by $$ \vec{v} = \vec{r}_0\ \dot{f}(s) + \vec{v}_0\ \dot{g}(s) $$

Where $$ f(s) = 1-\left(\frac{\mu}{|\vec{r}_0|}\right) s^2 c_2(\alpha s^2) $$ $$ g(s) = t-t_0-\mu s^3c_3(\alpha s^2) $$

With $\alpha=\frac{\mu}{a}$, $c_n$ the n-th Stumpff function and $s$ the solution of $$ t-t_0=|\vec{r}_0| s \ c_1(\alpha s^2) + |\vec{r}_0| s^2 \ c_2(\alpha s^2) + \mu s^3 \ c_3(\alpha s^3) $$

I am calculating $s$ by using Newton's Method and everything works fine for scenarios without thrust. The orbit is elliptical and closed, parabolic for higher initial velocities, if I choose a $\vec{r}_0$, $\vec{v}_0$ and $t_0$ and propagate $t$ forward.

For cases with thrust where I would change the velocity in every iteration

$$ \vec{v}_{k} = \vec{v}_{k-1} + \vec{a}_{k} \cdot \Delta t $$

it is necessary to update $\vec{r}_0$ and $\vec{v}_0$ every time, propagate $\Delta t$ forward and repeat instead of choosing initial conditions and just propagating $t$ forward. In preparation for that I chose not to add any velocity so both approaches should yield the same results.

Here is my Problem: Even after the first iteration the differences are large and the error increases fast

blue: fixed initial conditions, propagating $t$; red: propagating initial conditions, fixed $\Delta t$

blue: fixed initial conditions, propagating $t$

red: propagating initial conditions, fixed $\Delta t$

I thought I could solve this by using a more advanced integration algorithm like Runge-Kutta and tried to calculate every iteration in as few steps as possible but I was not able to transform the equations so that I could use Runge-Kutta (since it is not an ODE) and reducing the steps did not help at all.

Can someone help to fix this or give any hints to why this error is so large? Thank you in advance!

|cite|improve this question
$\endgroup$
  • 1
    $\begingroup$ What do you mean by "reducing the steps did not help"? Reducing the number of steps or reducing the duration of the steps? The first thing you should to is change $\Delta t$. Reduce it by a factor of ten. The program will run more slowly, but it the problem may simply be that your $\Delta t$ is too large. $\endgroup$ – garyp Aug 17 '16 at 18:34
  • $\begingroup$ Sorry, I should have clarified this. By "reducing the steps" I meant that I tried avoiding numerical errors by only calculating for example the $\sin(\sqrt{x})$ of the Stumpff functions once even though it is needed in different places. $\endgroup$ – Ganymed_ Aug 17 '16 at 18:39
  • $\begingroup$ Avoiding recalculation might save time, but it won't change the numeric behavior of the code. It won't even save time if the environment is doing sufficient de-dpulication (I don't know anything about the internals of cpython (or any other implementation) so I won't even guess) because the environment might have made that optimization without telling you. $\endgroup$ – dmckee♦ Aug 17 '16 at 18:59
  • $\begingroup$ BTW, You might get a better response on Computational Science. $\endgroup$ – dmckee♦ Aug 17 '16 at 19:00
  • 4
    $\begingroup$ I'm voting to close this question as off-topic because it is about numerical computation, not the physics involved in the problem. $\endgroup$ – sammy gerbil Aug 18 '16 at 1:31