I am trying to solve numerically a problem in orbital motion, using Runge-Kutta 4 method with adaptive step size. Because energy is the most obvious theoretically conservative number in the motion, I want to dymically control the step size of the calculaion such that $m=\frac{\delta E_T}{E_0}$ (the relative error of the energy at the end of the calculation and the starting, analytical energy, the desired accuracy) is equal or smaller than some number.
I tried to follow this logic: first, I calculated the new state of the particle with some $dt$ (starting or given in previous iteration), use it to calculate the current energy and $s=\frac{\delta E_t}{E_0}$ (relative error in current step) and get a new timestep: $dt_{new} = dt(\frac{m}{s})^{0.2}$ (scaling the timestep as we have seen in class), and use the new timestep to calculate the state of the particle, and advance to next iteration. However, calculation attemps so far have been unsuccessful. I noticed that no matter how small initial $dt$ I give, the relative error eventually grows too large and in turn $dt$ collapses to extremly small numbers, and the calculation effectively becomes an infinite loop.
What did I do wrong? Thanks in advance
One observation is that you compare the local error over $dt$ to the (desired) global error over time $T$. This makes not much sense. What you want is that the local error is proportional to the global error, roughly $s:m = dt:T$. So you need to change the formula to $$ dt_{\rm new}=dt·\left(\frac{m·dt}{s}\right)^{0.2} ~~\text{ or fully }~~ dt_{\rm new}=dt·\left(\frac{m·dt}{s·T}\right)^{0.2}. $$ The exponent does only give a guess of the immediate optimal step size if the error in the energy has global order 5. However, as a regulation process it will tend towards an optimal step size over a small number of steps even if order and exponent are not that closely related.