After playing a game called "Kerbal Space Program" I got interested in orbital mechanics and started messing with simplified calculations to determine $\Delta v$ requirements. In which I compared two trajectories to get into orbit/land from orbit to see how much it would matter if you would burn vertical or horizontal (Hohmann) at the surface. Graphs of some results can be seen here.
However after this I also wanted to know how much time it would take to perform these maneuvers. I assumed that changes in velocity are instantaneous, so this would be equal to the time between the two changes in velocity. For the Hohmann-like-transfer it is simple, since you travel from periapsis tot apoapsis, which takes half of the orbital period (due to symmetry). But for the other transfer, which applies a vertical change in velocity at the surface, it is very difficult (if there is a sidereal rotational velocity). I tried solving it using Kepler's second law, given this formula:
$$ r=\frac{a(1-e^2)}{1+e\cos{\theta}} \tag{1} $$
Since the area underneath the path from $\theta_0$ to $\theta_1$ can be calculated like this:
$$ A_{\theta_0,\theta_1}=\int_{\theta_0}^{\theta_1}{\frac{1}{2}r^2d\theta} \tag{2} $$
Which can be used to determine the amount of time it would take, since the total area of the ellipse would take one orbital period $T=2\pi\sqrt{\frac{a^3}{\mu}}$:
$$ t_{\theta_0,\theta_1}=T\frac{\int_{\theta_0}^{\theta_1}{\frac{1}{2}r^2d\theta}}{\int_{0}^{2\pi}{\frac{1}{2}r^2d\theta}}=T\frac{\int_{\theta_0}^{\theta_1}{\frac{1}{2}r^2d\theta}}{\pi a^2\sqrt{1-e^2}} \tag{3} $$
I am not that good at solving integrals, so I both tried to solve it with MATLAB and WolframAlpha, which gave me this as the result:
$$ \int{\frac{1}{2}r^2d\theta}=a^2\sqrt{1-e^2}\left(\tan^{-1}\left(\frac{\sqrt{1-e^2}\tan{\frac{\theta}{2}}}{1+e}\right)-\frac{e\sqrt{1-e^2}\sin{\theta}}{2(1+e\cos{\theta})}\right) \tag{4} $$
So:
$$ t_{\theta_0,\theta_1}=\sqrt{\frac{a^3}{\mu}}\left[2\tan^{-1}\left(\frac{\sqrt{1-e^2}\tan{\frac{\theta}{2}}}{1+e}\right)-\frac{e\sqrt{1-e^2}\sin{\theta}}{1+e\cos{\theta}}\right]^{\theta_1}_{\theta_0} \tag{5} $$
However when I tried to determine the integral with symbolic parameters and bounderies, MATLAB returned "Warning: Explicit integral could not be found." probably because there is not a general solution for both closed and open orbits. But if I use this formula and filled in the boundaries myself, it seemed that I got wrong results. However when I let MATLAB calculate the integral for given values of $a$, $e$, $\theta_0$ and $\theta_1$ I do get results which seem correct.
So is this a good approach to find a time dependency of an (elliptical) orbit. If so is there an (continuous) equation derivable from this. And what might be a better approach?
EDIT 1: The unexpected results form this formula might be caused by the fact that $\tan^{-1}$ returns results from $-\frac{\pi}{2}$ to $\frac{\pi}{2}$. Since when I calculated the time from $\theta=0$ to $\theta=\pi$ I got results which agree with the conservation of energy and angular momentum. For this is used the following equations to determine the time derivatives:
$$ \omega=\frac{\delta\theta}{\delta t}=\frac{\delta \theta / \delta \theta}{\delta t / \delta \theta}=\sqrt{\frac{\mu}{a^3\left(1-e^2\right)^3}}\left(1+e\cos{\theta}\right)^2 \tag{6} $$
$$ \dot{r}=\frac{\delta r}{\delta t}=\frac{\delta r / \delta\theta}{\delta t / \delta\theta}=\omega\frac{\delta r}{\delta\theta}=\omega\frac{ae(1-e^2)\sin{\theta}}{\left(1+e\cos{\theta}\right)^2} \tag{7} $$
Which yields that the specific orbital energy $\epsilon$ and specific angular momentum $h$ are constant when substituting theses therms into there expressions:
$$ \epsilon=\frac{v^2}{2}-\frac{\mu}{r}=\frac{\dot{r}^2+\omega^2r^2}{2}-\frac{\mu}{r}=-\frac{\mu}{2a} \tag{8} $$
$$ h=\omega r^2=\sqrt{a\mu\left(1-e^2\right)} \tag{9} $$
Edit 2: To show that equation $(5)$ is equivalent to Kepler's equation without relying on Wolfram Alpha or Mathematica to do trigonometric simplifications, you can multiply both sides by $\sqrt{\frac{\mu}{a^3}}$, which turns the left hand side into the mean anomaly. The first term can easilly be shown to be equal to the eccentric anomaly, because the relation between them is defined as
$$ \tan \frac{\theta}{2} = \sqrt{\frac{1+e}{1-e}} \tan\frac{E}{2}, \tag{10a} $$
from which follows
$$ \tan\frac{E}{2} = \sqrt{\frac{1-e}{1+e}} \sqrt{\frac{1+e}{1+e}} \tan \frac{\theta}{2} = \frac{\sqrt{1-e^2}}{1+e} \tan \frac{\theta}{2}, \tag{10b} $$
$$ E = 2 \tan^{-1}\left(\frac{\sqrt{1-e^2}}{1+e} \tan \frac{\theta}{2}\right). \tag{10c} $$
For the second term it is a little more work to show that it is equal to minus $e\sin E$. The minus and $e$ can be factored out in both terms, such that we need to prove that
$$ \sin E = \frac{\sqrt{1-e^2}\sin{\theta}}{1+e\cos{\theta}}. \tag{11} $$
The following trigonometric identities will be used:
- $\tan x = \frac{\sin x}{\cos x}$
- $\sin 2x = \frac{2 \tan x}{1 + \tan^2 x}$
- $\sin^2 x + \cos^2 x = 1$
- $\sin 2x = 2 \sin x \cos x$
- $\cos 2x = \cos^2 x - \sin^2 x$
We have an expression for $\tan\frac{E}{2}$ as a function of $\theta$ from equation $\text{(10a)}$, this will initially be represented by $\alpha$
$$ \sin E = \frac{2 \tan \frac{E}{2}}{1 + \tan^2 \frac{E}{2}} = \frac{2 \alpha}{1 + \alpha^2}. \tag{12a} $$
Substituting equation $\text{(10b)}$ for $\alpha$ into the last equation yields
$$ \sin E = \frac{2 \frac{\sqrt{1-e^2}}{1+e} \tan \frac{\theta}{2}}{1 + \left(\frac{\sqrt{1-e^2}}{1+e} \tan \frac{\theta}{2}\right)^2} = \frac{2 \frac{\sqrt{1-e^2}}{1+e} \frac{\sin \frac{\theta}{2}}{\cos \frac{\theta}{2}}}{1 + \frac{1-e^2}{(1+e)^2} \frac{\sin^2\!\frac{\theta}{2}}{\cos^2\!\frac{\theta}{2}}} = \frac{2 \sqrt{1-e^2} \sin \frac{\theta}{2} \cos \frac{\theta}{2}}{(1+e) \cos^2\!\frac{\theta}{2} + (1-e) \sin^2\!\frac{\theta}{2}} = \frac{\sqrt{1-e^2} \sin \theta}{\cos^2\!\frac{\theta}{2} + \sin^2\!\frac{\theta}{2} + e \left(\cos^2\!\frac{\theta}{2} - \sin^2\!\frac{\theta}{2}\right)} = \frac{\sqrt{1-e^2} \sin \theta}{1 + e \cos\theta}. \tag{12b} $$
So this shows that equation $\text{(11)}$ is true and equation $\text{(5)}$ is equivalent to Kepler's equation.
Your equation is equivalent to Kepler's equation, just in terms of true anomaly $\theta$ instead of eccentric anomaly $E$. Read more about it here.
$\tan \frac{\theta}{2} = \sqrt{\frac{1 + e}{1 - e}} \tan \frac{E}{2} \\ E = 2 \tan^{-1} \left( \sqrt{\frac{1 - e}{1 + e}} \tan \frac{\theta}{2} \right) = 2 \tan^{-1} \left( \frac{\sqrt{1 - e^2}}{1 + e} \tan \frac{\theta}{2} \right) \\ n t = M = E - e \sin E \\ t = \frac{P}{2 \pi} ( E - e \sin E )$
Having Wolfram Alpha or Mathematica do the trig simplification, we get (moving $i$ inside the square root)
$\sin E = \frac{\sqrt{1-e^2} \sin \theta}{1 + e \cos \theta}$
