1

I am trying to simulate a rocket launch by solving the equations of motion for a rocket launched from earth. The system of equations is as follows: $$ \frac{d^{2}y}{dt^{2}} = \frac{T}{m(t)} - \frac{GM}{y^{2}} \quad , \quad \frac{dm}{dt} = \frac{-T}{I_{sp}g_0} \quad , \quad y(0)=y_0 \ , \ v(0)=0$$

Where $m(t)$ is the time-dependant mass of the rocket, $y_0$ is the radius of the earth, $g_{0}$ is the gravitational constant at $y_{0}$ and $T,I_{sp}$ are rocket parameters to be chosen freely in my model. I have implemented explicit Euler, velocity Verlet and RK4 (all correctly i believe) and they all seem to work fine, even for large time steps.

My issue is that i have been asked to motivate my choice of integrator and show some form of accuracy metric for this decision. In previous problems of the same kind, there was an exact solution to compare to or, in the case of velocity Verlet, you could use energy conservation in for example a harmonic pendulum to judge accuracy. In this case, I am the one to model the rocket freely so there is no exact solution and i assume the energy is not conserved because of the fuel burning so i can't use that either. How can i choose one of these and show its accuracy?

I would appreciate any suggestions.

Lere1224
  • 43
  • 1
    The easiest way is to compare to a second solution with the same method at half or double step. See also "Richardson extrapolation". See for an example of this strategy https://math.stackexchange.com/questions/3058387/empirical-error-proof-runge-kutta-algorithm-when-not-knowing-exact-solution/3058957#3058957 – Lutz Lehmann Jan 18 '23 at 14:14
  • @LutzLehmann Does halving the step size and showing that it converges count as proving the accuracy? When decreasing the step size, i of course see convergence towards a specific value but there is no way of knowing whether that convergence value is the right one. – Lere1224 Jan 18 '23 at 17:13
  • 1
    You can now compare the errors of the different methods. Or you can compute a sequence of errors for the same method but different step sizes and check if they behave according to the order of the method. The finiteness of floating-point numbers imposes a lower bound on the step sizes down to where the error follows the order. – Lutz Lehmann Jan 18 '23 at 19:05
  • @LutzLehmann Sorry but i don't quite understand which error? I don't know what to calculate the error to since there is no "real" solution. Is it the error to the convergence value? – Lere1224 Jan 18 '23 at 21:47
  • 1
    Step doubling gives you an error estimate which is quite good. It can happen that the error oscillates and the interval endpoint is close to a root of this error function. Thus take the maximum or $L^1$ norm or .. over the whole interval. – Lutz Lehmann Jan 19 '23 at 06:55
  • @LutzLehmann Thank you for your answers, i think i get what you mean now. Can i do this for the global error i.e the final altitude of the rocket or should i do it locally for the first step? Perhaps it doesn't matter? – Lere1224 Jan 19 '23 at 08:52
  • 1
    This is all for a fixed integration interval, if that is to the final altitude or some small initial segment. But never for the first step, as that depends on the step size, is not a fixed interval. With a small segment you get errors that are closer to the theoretical bounds. Over longer intervals the theoretical bounds that always estimate the worst case, can be overly pessimistic – Lutz Lehmann Jan 19 '23 at 09:02

0 Answers0