I'm trying to solve the Kepler problem in Cartesian coordinates, that is, I want to show that the trajectory is an ellipse using Cartesian coordinates instead of using polar coordinates, as is usually done. For simplicity, I will solve the one-body problem in which one of the two masses, with mass $M$, is assumed to be static, while only the other mass, $m$ is moving, and $M\gg m$. Assuming gravity is the only force acting on the particle of mass $m$, the equation of motion is $$m\ddot{\mathbf{r}}=-\frac{GMm}{|\mathbf{r}|^3}\mathbf{r}=-\frac{GMm}{|\mathbf{r}|^2}\hat{\mathbf{r}},$$ which is equivalent to the following equations for the $x$ and $y$ coordinates $$m\ddot{x}=-\frac{GMmx}{(x^2+y^2)^{\frac{3}{2}}}, \qquad m\ddot{y}=-\frac{GMmy} {(x^2+y^2)^{\frac{3}{2}}},$$ since the motion is confined to a plane. Now that I have these equations I don't know how to proceed and show that the trajectory is an ellipse. How does one show from these equations that the trajectory is a conic section of eccentricity $\epsilon$? Also, how does one derive Kepler's third law from this?

  • 5
    $\begingroup$ Oh my, why do you want to do this? Do you have a source that guarantees you are going to be able to solve this analytically? It wouldn't surprise me if you could not directly solve the equations in Cartesian coordinates with standard methods. You might have to use numerical methods to make progress. But also... what insight are you hoping to gain by solving the problem in Cartesian coordinates? Knowing that might help someone suggest an alternative way to approach the problem to get at what it is you are really after. $\endgroup$
    – Andrew
    Mar 31, 2022 at 1:31
  • 3
    $\begingroup$ @DonAl An important point is that the physics does not depend on the choice of coordinates, but the difficulty of the calculation very much can depend on the choice of coordinates. It is definitely possible to take a solvable problem, make a poor choice of coordinates, and end up with a problem that is not solvable with known methods except by undoing the transformation. I am not 100% sure that the Kepler problem can't be analytically solved in Cartesian coordinates, but it wouldn't surprise me. $\endgroup$
    – Andrew
    Mar 31, 2022 at 3:06
  • 2
    $\begingroup$ However, there are ways you can still check that the physics is invariant under changes of coordinates. For example, you could take the solution in polar coordinates, transform it to Cartesian coordinates, and check explicitly that it solves Newton's laws in Cartesian coordinates. $\endgroup$
    – Andrew
    Mar 31, 2022 at 3:07
  • 2
    $\begingroup$ To show that the trajectory is an ellipse, use reduce two body problem $\endgroup$
    – Eli
    Mar 31, 2022 at 13:57
  • 3
    $\begingroup$ If you're after an example where you can solve a 2D problem in two different sets coordinates, consider the 2D harmonic oscillator instead. In Cartesian coordinates, it's easy to show that the answer is an ellipse centered at the origin. But the problem can be solved in polar coordinates with the same result (though it's a bit more involved.) $\endgroup$ Mar 31, 2022 at 15:15

3 Answers 3


Absorb the dimensional GM into the units, to remove excuses for not recognizing the plane-geometry structure. Note the rotational and translational invariance to be used in fixing your coordinate system below. All vectors are then 2-vectors on that plane, $$\ddot{ x }=-{ { x}\over r^3},\qquad \ddot{ y }=-{ { y}\over r^3} ~. $$ It is then self-evident that $$ L=x\dot y -y\dot x $$ is a constant of the motion (in the z-direction, the only non vanishing one, of course), $\dot L=0$.

Moreover, the mass/normalized LRL 2-vector on that plane, $$ \vec e= L \begin{pmatrix} \dot y\\-\dot x\end{pmatrix}-{1\over r} \begin{pmatrix} x\\ y\end{pmatrix} $$ is also conserved, $\dot{\vec e }=0$. Its magnitude will turn out to be the eccentricity.

Dotting by $\vec r$ you have $$ \vec r\cdot \vec e =L^2 -r, \leadsto \\ r= L^2-\vec r\cdot \vec e , \leadsto \\ x^2+y^2= (L^2-\vec r\cdot \vec e )^2. $$ You may use rotational invariance to take $\vec e$ along the x-axis, $\vec e=-\epsilon \hat x$ to prettify your ellipse orientation, and work out the r.h.s. square, as a quadratic polynomial in x, with the obvious constants. Elementary algebra leads you to the ellipse your teacher taught you in terms of the constants Ξ΅ and L. That is, shifting the origin of xs to $L^2 \epsilon/(1-\epsilon^2)$ and taking $\epsilon ^2= 1-b^2/a^2$, you get $$ \frac{x^2}{a^2}+\frac{y^2}{b^2} = L^4 a^2/b^4. $$

  • 1
    $\begingroup$ +1, but I must note in passing the LRL trick only works because it's a $1/r^2$ force, so the Binet equation is still the wisest strategy for more general potentials. $\endgroup$
    – J.G.
    Mar 31, 2022 at 16:04
  • $\begingroup$ Of course, but even in polar coordinates, which the OP wishes to bypass, the LRL trick requires an 1/rr force. $\endgroup$ Mar 31, 2022 at 16:35
  • $\begingroup$ Well, some other special cases of the Binet equation are easily solved. The most physically relevant is adding a small $1/r^3$ perturbation to the force, which causes a precessing elliptical orbit. $\endgroup$
    – J.G.
    Mar 31, 2022 at 16:36
  • 1
    $\begingroup$ @Eli You appreciate the OP has chosen the angular momentum in the z-direction, so, then, r strictly in the x-y plane, no? $\endgroup$ Apr 4, 2022 at 14:27
  • 1
    $\begingroup$ This was not the OP's question. WP has plenty of general and elegant derivations of Kepler orbits and the OP is indicating he's not after those. The OP is strictly looking at the trajectory plane and is striving to precisely avoid /bypass polar coordinates. This is his desideratum, no? You want to go easy, just do Hamilton's hodograph, as described in detail there. But this is not the question here. $\endgroup$ Apr 4, 2022 at 15:37

I can't comment directly on posts, so I'll post a note here (for NinjaDarth). The other comment is right: the solution is too over-wrought and (even as noted in the reply) it can be simplified. So, I'll actually do it here, following the suggestions made in the reply. The result is much simpler and more direct.

We'll keep the notation: the Cartesian coordinates are wrapped up in 3D-vector form. The Kepler problem is framed as: $$ 𝐫'(t) = 𝐯(t), \hspace 1em 𝐫(0) = 𝐫_0, \\ 𝐯'(t) = -\frac{μ𝐫(t)}{r(t)^3}, \hspace 1em 𝐯(0) = 𝐯_0, $$ where $r = |𝐫|$, with a similar convention to denote magnitudes of vectors by using the light-face version of the bold-face letter. We'll change the clock to $G$ by writing $$\frac{dt}{dG} = r.$$

When this is applied to $𝐫$ and $r$, the results - after using some vector algebra - are: $$ \frac{d𝐫}{dG} = r𝐯, \hspace 1em \frac{d(r𝐯)}{dG} = z𝐫 + 𝐚, \\ \frac{dr}{dG} = 𝐫·𝐯, \hspace 1em \frac{d(𝐫·𝐯)}{dG} = zr + ΞΌ, $$ where $$z = v^2 - \frac{2ΞΌ}{r}, \hspace 1em 𝐚 = \frac{μ𝐫}{r} - 𝐯×(𝐫×𝐯),$$ both of which prove to be constant, and to be connected to the better-known constants of the Kepler problem: $$𝐑 = 𝐫×𝐯, \hspace 1em 𝐞 = \frac{𝐯×𝐑}{ΞΌ} - \frac{𝐫}{r}, \hspace 1em H = \frac{v^2}{2} - \frac{ΞΌ}{r}, $$ by $𝐚 = -μ𝐞$ and $z = 2H$.

The Solution:
We will use the functions $(C,S,D,T)$ of $G$ and $z$ given by $$\frac{d}{dG}(C,S,D,T) = (zS,C,S,D), \hspace 1em (C,S,D,T) = (1,0,0,0) @ G = 0.$$ Noting that the general solution to the problem $$a'(G) = b, \hspace 1em b'(G) = z a(G) + c, \hspace 1em a(0) = a_0, \hspace 1em b(0) = b_0,$$ is $$a = a_0 C + b_0 S + c D, \hspace 1em b = b_0 C + (c + a_0 z) S,$$ then it follows that: $$ 𝐫 = 𝐫_0 C + r_0 𝐯_0 S + 𝐚 D, \\ r = r_0 C + 𝐫_0·𝐯_0 S + μ D, \\ 𝐯 = \frac{r_0 𝐯_0 C + (𝐚 + z 𝐫_0) S}{r}, $$ and, upon integration of $dt/dG = r$ with $t(0) = 0$: $$t = r_0 S + 𝐫_0·𝐯_0 D + μ T.$$

To recover $t_0$ as the sixth Kepler constant (along with the five independent parameters that come out of $(𝐞, 𝐑, H)$), we could reset $G$ so that $G = 0$ when $r$ is at its minimum, at which time, we would have $𝐫_0·𝐯_0 = 0$; but it's not really necessary, since $𝐫_0$ and $𝐯_0$ already provide us with six independent parameters, in place of the six Kepler parameters.

Conic Section Equation:
Actual equations can written for $𝐫$ in vector form, which will show that it is a conic or a line, by using the identities: $$C = 1 + z D, \hspace 1em S^2 = 2 D + z D^2.$$ First: for the following, define $$𝐛 = 𝐚 + z 𝐫_0, \hspace 1em 𝐜 = 𝐛 Γ— 𝐯_0 = 𝐚 Γ— 𝐯_0 + z 𝐑.$$ Then, second: we have $$𝐫 - 𝐫_0 = r_0 𝐯_0 S + 𝐛 D,$$ $$\left(𝐫 - 𝐫_0\right) Γ— 𝐯_0 = \left(r_0 𝐯_0 S + 𝐛 D\right) Γ— 𝐯_0 = 𝐛 Γ— 𝐯_0 D = 𝐜 D,$$ $$𝐛 Γ— \left(𝐫 - 𝐫_0\right) = 𝐛 Γ— \left(r_0 𝐯_0 S + 𝐛 D\right) = r_0 𝐛 Γ— 𝐯_0 S = r_0 𝐜 S.$$ Then, from $S^2 = 2D + z D^2$, upon cancellation of the factor of $c^2$, follows: $$\left|𝐛 Γ— \left(𝐫 - 𝐫_0\right)\right|^2 = {r_0}^2 \left(2𝐜·\left(𝐫 - 𝐫_0\right) Γ— 𝐯_0 + \left|\left(𝐫 - 𝐫_0\right) Γ— 𝐯_0\right|^2 \right).$$ Otherwise, we have: $$𝐜 = 𝟎.$$

When expressed in terms of the constants $𝐞$, $𝐑$ and $H$, noting that: $$ \left(𝐫 - 𝐫_0\right)×𝐯_0 = 𝐫×𝐯_0 - 𝐑, \\ 𝐛 = 𝐚 + z 𝐫_0 = -ΞΌ 𝐞 + 2 H 𝐫_0, \\ 𝐜 = 𝐛 Γ— 𝐯_0 = 𝐚 Γ— 𝐯_0 + z 𝐑 = -ΞΌ 𝐞 Γ— 𝐯_0 + 2 H 𝐑, $$ we can write the following: $$ \left|(-ΞΌ 𝐞 + 2H 𝐫_0) Γ— \left(𝐫 - 𝐫_0\right)\right|^2 = {r_0}^2 \left(2 \left(-ΞΌ 𝐞 Γ— 𝐯_0 + 2 H 𝐑\right)Β·\left(𝐫 Γ— 𝐯_0 - 𝐑\right) + \left|𝐫 Γ— 𝐯_0 - 𝐑\right|^2 \right). $$

The case $𝐜 = 𝟎$ has to be handled separately. The condition entails that $$\left(𝐫 - 𝐫_0\right) Γ— 𝐯_0 = 𝐜 D = 𝟎,$$ and may only occur if either $𝐯_0 = 𝟎$ - which gives rise to linear motion $𝐫 = 𝐫_0 C + 𝐚 D = 𝐫_0 + 𝐛 D$ - or $𝐯_0 β‰  0$ - which also also gives rise to linear motion, with $𝐫 = 𝐫_0 + (r_0 S + 𝐛·𝐯_0/{v_0}^2 D) 𝐯_0$.


In Cartesian coordinates? Ok. We'll do a full treatment, here, that covers all the different orbit types (elliptical, circular, hyperbolic, recti-linear) in a seamless unification. To actually be able to do such a unification requires a fudge-parameter for wiggle room (which I'll call $Ξ›$, for inverse speed) - along with a specification of different conventions for normalization to more familiar forms. It also requires combining the trigonometric and hyberbolic functions (along with their parabolic limits) seamlessly. This is done with the following set of functions: $$e_n(z,G) = \frac{G^n}{n!} + z \frac{G^{n+2}}{(n+2)!} + z^2 \frac{G^{n+4}}{(n+4)!} + z^3 \frac{G^{n+6}}{(n+6)!} + β‹―,\hspace 1em (n = 0, 1, 2, 3, β‹―)$$ with the cases of interest being $$C = e_0(z,G), \hspace 1em S = e_1(z,G), \hspace 1em D = e_2(z,G), \hspace 1em T = e_3(z,G),$$ which is all you'll need. The angle $G$ combines the angle $E$ used in general treatments for the elliptic case with the angle $H$ used for hyperbolic orbits (which, however, I will refer to, below as $F$ - as $H$ is going to be used for the total energy). Also: it will be done primarily with 3-dimensional vector algebra. That includes the convention of using boldface for a vector (e.g. $𝐫$) and the lighter-face for its magnitude (e.g. $r = |𝐫|$). There's no conflict with Newton's constant $G$, since we're instead using a coefficient $ΞΌ$ that may, more properly, be called "Kepler's constant" for the given gravitating body. If the body is the Earth, then in terms of its mass $M$ and Newton's coefficient $G$, $ΞΌ = GM ≃ 4Γ—10¹⁴\ \mbox{metersΒ³/secondΒ²}$, which is the "gravitational charge" of the Earth.

The Problem:
First, define the equations of motion for the Kepler Problem in vector form as: $$\frac{d𝐫}{dt} = 𝐯, \hspace 1em \frac{d𝐯}{dt} = -\frac{μ𝐫}{r^3},$$ where $μ > 0$.

Second, note (as per many treatments of the problem) that the following are constants of motion: $$𝐑 = 𝐫×𝐯, \hspace 1em 𝐞 = \frac{𝐯×𝐑}{ΞΌ} - \frac{𝐫}{r}, \hspace 1em H = \frac{v^2}{2} - \frac{ΞΌ}{r}.$$ You can verify this, after the fact, by applying $d/dt$ to each of them, using the equations of motion. They satisfy the following identities: $$𝐑·𝐞 = 0, \hspace 1em 2H h^2 = ΞΌ^2 \left(e^2 - 1\right).$$ So, out of 6 constants of motion, 5 independent ones are included in this set. The 6th one is the time $t_0$ of nearest approach, which is not unique if the orbit is periodic.

Of these constants, $𝐑$ is generic to all central force problems, and a version of $H$ will appear in the cases of forces derivable from potentials. The potential, here, is $-μ/r$. The existence and form of the constant $𝐞$ is specific to the Kepler problem. I don't know of anything analogous to this for the generic case of the central force problem.

Third, set up an orthonormal frame $(𝐒, 𝐣, 𝐀)$, such that $𝐞 = e 𝐒$, $𝐑 = h 𝐀$ and $𝐣 = 𝐀 Γ— 𝐒$. Then, the following are true:

  • the cases $𝐞 = 𝟎$ and $𝐑 = 𝟎$ cannot occur together;
  • if $𝐞 = 𝟎$, the orbit will turn out to be a circle, and the vectors $𝐒$ and $𝐣$ can be chosen as any two mutually orthogonal vectors perpendicular to $𝐑$; otherwise, the vector $𝐒$ points in the direction of the point of nearest approach;
  • if $𝐑 = 𝟎$, the orbit is in a line - a collision orbit - the point of "nearest approach" is the center of force, itself, which $𝐒$ points toward, and the vectors $𝐣$ and $𝐀$ can be chosen as any two mutually orthogonal vectors perpendicular to $𝐞$; otherwise, the orbit is in a plane perpendicular to $𝐑$.
Thus, without loss of generality, we may write: $$𝐫 = x 𝐒 + y 𝐣, \hspace 1em 𝐯 = xΜ‡ 𝐒 + yΜ‡ 𝐣, \hspace 1em r = \sqrt{x^2 + y^2} β‰₯ 0,$$ and upon substitution into the expressions for the constants of motion, we have $$ 𝐑 = h 𝐀 = (x yΜ‡ - y xΜ‡) 𝐀 β‡’ h = x yΜ‡ - y xΜ‡,\\ 𝐞 = e 𝐒 = \left(\frac{yΜ‡ h}{ΞΌ} - \frac{x}{r}\right) 𝐒 - \left(\frac{xΜ‡ h}{ΞΌ} + \frac{y}{r}\right) 𝐣 β‡’ xΜ‡ = -\frac{ΞΌ}{h} \frac{y}{r}, \hspace 1em yΜ‡ = \frac{ΞΌ}{h} \left(\frac{x}{r} + e\right). $$ Thus $$h = \frac{ΞΌ}{h} \left(x\left(\frac{x}{r} + e\right) + y \frac{y}{r}\right) = \frac{ΞΌ}{h} (r + xe) β‡’ r = \frac{h^2}{ΞΌ} - x e.$$

Finally, adopt the trick - standard to the treatment of the Kepler Problem - of replacing the time $t$ with a clock, $G$, more suited to the problem: $$\frac{dt}{dG} = Ξ› r$$ where $Ξ› > 0$ is an arbitrary constant. Then, noting that $$r \frac{dr}{dt} = 𝐫·\frac{d𝐫}{dt} = 𝐫·𝐯 β‡’ \frac{dr}{dt} = \frac{𝐫·𝐯}{r}$$ and $$\frac{d}{dG} = Ξ› r \frac{d}{dt}$$ we have: $$\frac{d^2r}{dG^2} = Ξ›r \frac{d}{dt} \left(Ξ›r \frac{dr}{dt}\right) = Ξ›^2 r \frac{d}{dt} (𝐫·𝐯) = Ξ›^2 r \left(𝐯·𝐯 - \frac{𝐫·μ𝐫}{r^3}\right) = Ξ›^2 \left(r v^2 - ΞΌ\right) = Ξ›^2 (2Hr + ΞΌ) = zr + Ξ²$$ where we set: $$z = 2 Ξ›^2 H, \hspace 1em Ξ² = Ξ›^2 ΞΌ > 0.$$ In terms of Ξ², noting the assumption $Ξ› > 0$, we have $$Ξ› = \sqrt{\frac{Ξ²}{ΞΌ}}.$$ Since $r β‰₯ 0$, then the solution for $r$ must be bounded below, and we'll set $G = 0$ to be the time where $r$ is at the minimum, which we'll call $Ξ±$: $$G = 0 β‡’ r = Ξ±, \hspace 1em \frac{dr}{dG} = 0.$$ More precisely, we set it to the most recent time of a local minimum, or else to the current time if $r$ is locally constant.

Since $dr/dG = -e dx/dG$ and $e β‰₯ 0$, then $$G = 0 β‡’ \frac{dx}{dG} = 0 β‡’ 0 = Ξ› r xΜ‡ = -\frac{ΞΌΞ›}{h} y β‡’ y = 0.$$ (In the case $h = 0$, the orbit is in a line, aligned with $𝐒$, so $y = 0$ already also holds there at $G = 0$ and at all times, so there's no loss of generality here.)

Finally, noting that $$h^2 = |𝐫×𝐯|^2 = r^2 v^2 - (𝐫·𝐯)^2,$$ the condition that $r$ be minimized when $G = 0$ yields the following: $$G = 0 β‡’ 0 ≀ \frac{d^2r}{dt^2} = \frac{d}{dt} \left(\frac{𝐫·𝐯}{r}\right) = \frac{𝐯·𝐯}{r} + \frac{𝐫·-μ𝐫/r^3}{r} - \frac{𝐫·𝐯}{r^2} \frac{𝐫·𝐯}{r} = \frac{\left(r^2 v^2 - (𝐫·𝐯)^2\right)}{r^3} - \frac{ΞΌ}{r^2} = \frac{h^2}{r^3} - \frac{ΞΌ}{r^2}$$ or $$h^2 β‰₯ ΞΌ Ξ±.$$ From this, in turn, noting the identity $𝐫·𝐯×𝐑 = 𝐫×𝐯·𝐑 = h^2$, it follows that $$G = 0 β‡’ x e = 𝐫·𝐞 = 𝐫·\left(\frac{𝐯×𝐑}{ΞΌ} - \frac{𝐫}{r}\right) = \frac{h^2}{ΞΌ} - r β‰₯ 0,$$ so that $x > 0$ at $G = 0$. (In the case $e = 0$, as seen here, the orbit is a circle of radius $r = h^2/ΞΌ$, so one can set $G = 0$ at arbitrary times and orient $𝐒$ and $𝐣$ arbitrarily in the plane perpendicular to $𝐀$ to make $x > 0$ and $y = 0$ at $G = 0$ true by fiat, so again: no loss of generality here.)

The Solution:
Define the following functions of $G$: $$ C = 1 + z \frac{G^2}{2!} + z^2 \frac{G^4}{4!} + z^3 \frac{G^6}{6!} + β‹―, \\ S = G + z \frac{G^3}{3!} + z^2 \frac{G^5}{5!} + z^3 \frac{G^7}{7!} + β‹―, \\ D = \frac{G^2}{2!} + z \frac{G^4}{4!} + z^2 \frac{G^6}{6!} + β‹―, \\ T = \frac{G^3}{3!} + z \frac{G^5}{5!} + z^2 \frac{G^7}{7!} + β‹―. $$ They satisfy the following boundary value problem: $$ G = 0 β‡’ (C, S, D, T) = (1, 0, 0, 0), \\ \frac{βˆ‚}{βˆ‚G} (C, S, D, T) = (z S, C, S, D) $$ as well as the following identities: $$ C = 1 + z D, \hspace 1em S = G + z T,\\ C^2 - z S^2 = 1, \hspace 1em S^2 - z D^2 = 2D, \hspace 1em S^2 = (C + 1) D $$

Then, the following boundary problem is set up for $r$, in terms of $G$: $$ \frac{d^2r}{dG^2} = z r + Ξ²,\\ G = 0 β‡’ r = Ξ±, \hspace 1em \frac{dr}{dG} = 0 $$ which has the general solution $$r = r_0 + r_1 S + r_2 D.$$ Substituting into the differential equation and boundary conditions yields the following: $$r_0 = Ξ±, \hspace 1em r_1 = 0, \hspace 1em r_2 = z r_0 + Ξ²,$$ thus $$r = Ξ± (1 + z D) + Ξ² D = Ξ± C + Ξ² D.$$ From this, we also have, for the following boundary value problem: $$ \frac{dt}{dG} = Ξ› r = Ξ› (Ξ± C + Ξ² D),\\ G = 0 β‡’ t = t_0, $$ the following solution: $$t = t_0 + Ξ› (Ξ± S + Ξ² T) = t_0 + \sqrt{\frac{Ξ²}{ΞΌ}} (Ξ± S + Ξ² T).$$

For $𝐫$, we have the following differential equation: $$\frac{d^2𝐫}{dG^2} = Ξ› r \frac{d}{dt} \left(Ξ› r \frac{d𝐫}{dt}\right) = Ξ›^2 r \frac{d}{dt} (r 𝐯) = Λ² \left(r \frac{dr}{dt} 𝐯 + r^2 \frac{d𝐯}{dt}\right) = Ξ›^2 \left(𝐫·𝐯 𝐯 - \frac{ΞΌ 𝐫}{r}\right)$$ With the cross-product identity $$𝐯 Γ— 𝐑 = 𝐯 Γ— (𝐫 Γ— 𝐯) = v^2 𝐫 - 𝐫·𝐯 𝐯$$ and the definition of $𝐞$, we can write $$𝐫·𝐯 𝐯 = v^2 𝐫 - 𝐯 Γ— 𝐑 = v^2 𝐫 - ΞΌ \left(𝐞 + \frac{𝐫}{r}\right).$$ Thus $$\frac{d^2𝐫}{dG^2} = Ξ›^2 \left(v^2 𝐫 - ΞΌ \left(𝐞 + \frac{𝐫}{r}\right) - \frac{ΞΌ 𝐫}{r}\right) = Ξ›^2 \left(\left(v^2 - \frac{2ΞΌ}{r}\right) 𝐫 - ΞΌ 𝐞\right) = 2 Ξ›^2 H 𝐫 - Ξ›^2 ΞΌ 𝐞 = z 𝐫 - Ξ² 𝐞.$$ Since $r dr/dG = 𝐫·d𝐫/dG$, then we have, at $G = 0$, $𝐫·d𝐫/dG = 0$, which leads to the following boundary value problem: $$ \frac{d^2𝐫}{dG^2} = z 𝐫 - Ξ² 𝐞,\\ G = 0 β‡’ 𝐫·\frac{d𝐫}{dG} = 0. $$

This leads to the following boundary value problem for the coordinates, individually: $$ \frac{d^2x}{dG^2} = x z - Ξ² e, \hspace 1em \frac{d^2y}{dG^2} = y z,\\ G = 0 β‡’ x β‰₯ 0, \hspace 1em \frac{dx}{dG} = 0, \hspace 1em y = 0, \hspace 1em r = Ξ± β‡’ x = Ξ±. $$ Defining $$Ξ³ = \frac{dy}{dG} @G = 0,$$ and substituting into the general solution: $$x = x_0 + x_1 S + x_2 D, \hspace 1em y = y_0 C + y_1 S$$ yields $$x_2 = Ξ± z - Ξ² e, \hspace 1em y_0 = 0, \hspace 1em x_1 = 0, \hspace 1em x_0 = Ξ±, \hspace 1em y_1 = Ξ³$$ thus: $$x = Ξ± + (Ξ± z - Ξ² e) D = Ξ± (1 + z D) - Ξ² e D = Ξ± C - Ξ² e D, \hspace 1em y = Ξ³ S.$$ Matching to the condition that $x^2 + y^2 = r^2$ yields the following conditions on the parameters: $$(Ξ± C - Ξ² e D)^2 + (Ξ³ S)^2 = (Ξ± C + Ξ² D)^2$$ or, using the identities $C = 1 + z D$ and $S^2 = 2 D + z D^2$: $$Ξ³^2 = Ξ±Ξ²(e + 1), \hspace 1em Ξ±z = Ξ²(e - 1).$$ Using the second of these identities, we may rewrite $$x = Ξ± C - Ξ² e D = Ξ± + (Ξ± z - Ξ² e) D = Ξ± - Ξ² D.$$

So, the general solution is: $$ 𝐫 = (Ξ± - Ξ² D) 𝐒 + Ξ³ S 𝐣, \hspace 1em r = Ξ± C + Ξ² D, \\ t = t_0 + Ξ› (Ξ± S + Ξ² T), \hspace 1em 𝐯 = \frac{-Ξ² S 𝐒 + Ξ³ C 𝐣}{Ξ› r}, $$ where $$Ξ› = \sqrt{\frac{Ξ²}{ΞΌ}}; \hspace 1em Ξ±, Ξ³ β‰₯ 0, \hspace 1em Ξ² > 0, \hspace 1em Ξ³^2 = Ξ±Ξ²(e + 1), \hspace 1em Ξ±z = Ξ²(e - 1).$$ with the constants of motion reducing to: $$ 𝐑 = 𝐫×𝐯 = (Ξ± 𝐒)Γ—\frac{Ξ³ 𝐣}{Ξ› Ξ±} = \frac{Ξ³}{Ξ›} 𝐀 = h 𝐀 β‡’ h = \frac{Ξ³}{Ξ›}, \\ 𝐞 = \frac{𝐯×𝐑}{ΞΌ} - \frac{𝐫}{r} = \frac{Ξ³ 𝐣}{Ξ› Ξ±}Γ—\frac{Ξ³ 𝐀}{Ξ›} \frac{1}{ΞΌ} - 𝐒 = \left(\frac{Ξ³^2}{Ξ›^2 Ξ± ΞΌ} - 1\right) 𝐒 = \left(\frac{Ξ³^2}{Ξ± Ξ²} - 1\right) 𝐒 = e 𝐒 β‡’ e = \frac{Ξ³^2}{Ξ± Ξ²} - 1, \\ H = \frac{v^2}{2} - \frac{ΞΌ}{r} = \frac{1}{2} \frac{Ξ³^2}{(Ξ› Ξ±)^2} - \frac{ΞΌ}{Ξ±} = \frac{1}{2} \frac{\left(Ξ³^2 - 2Ξ±Ξ²\right)}{(Ξ› Ξ±)^2} = \frac{1}{2} \frac{Ξ±Ξ²(e - 1)}{Ξ›^2 Ξ±^2} = \frac{ΞΌ z}{2 Ξ²} $$

Normalizations and Reductions:
The different normalizations and reductions are described (up to a change in notation) in a response to the following query here: https://space.stackexchange.com/questions/20085/calculate-true-anomaly-at-future-point-in-time-with-hyperbolic-orbits.

The normalizations include the following cases:

  • $h > 0 β‡’ Ξ± = Ξ² (e + 1) = Ξ³, \hspace 1em z = \frac{e - 1}{e + 1}$,
  • $h > 0 β‡’ Ξ² = Ξ± (e + 1) = Ξ³, \hspace 1em z = e^2 - 1$,
  • $h > 0 β‡’ Ξ± = \frac{Ξ³}{\sqrt{e + 1}} = Ξ², \hspace 1em z = e - 1$,
  • $z > 0 β‡’ z = +1 β‡’ Ξ± = Ξ² (e - 1) = A, \hspace 1em Ξ³ = Ξ± \sqrt{e^2 - 1} = B$,
  • $z < 0 β‡’ z = -1 β‡’ Ξ± = Ξ² (1 - e) = a, \hspace 1em Ξ³ = Ξ± \sqrt{1 - e^2} = b$.
For recti-linear orbits, $Ξ± = 0 = Ξ³$, and $e = 1$, so that only $Ξ²$ and $z$ survive.

Irrespective of normalizations, more familiar forms associated with the elliptical/circular and hyperbolic cases can be recovered as follows: which also includes the recti-linear cases when $e = 1$, respectively with $b = 0$ or $B = 0$.

Elliptical Case ($z < 0$): major semi-axis $a$, minor semi-axis $b = a \sqrt{1 - e^2}$ $$ (G, C, S, D, T) = \left(\frac{E}{\sqrt{-z}}, \cos E, \frac{\sin E}{\sqrt{-z}}, \frac{\cos E - 1}{z}, \frac{\sin E - E}{z\sqrt{-z}}\right), \\ (Ξ±, Ξ², Ξ³) = (a(1 - e), -az, b \sqrt{-z}), \\ 𝐫 = a (\cos E - e) 𝐒 + b \sin E 𝐣, \hspace 1em r = a (1 - e \cos E), \\ t = t_0 + \sqrt{\frac{a^3}{ΞΌ}} (E - e \sin E), \hspace 1em 𝐯 = \frac{b \cos E 𝐣 - a \sin E 𝐒}{r} \sqrt{\frac{ΞΌ}{a}}, \\ h = b \sqrt{\frac{ΞΌ}{a}}, \hspace 1em H = -\frac{ΞΌ}{2a}, \\ b β‰  0 β‡’ \left(\frac{x}{a} + e\right)^2 + \left(\frac{y}{b}\right)^2 = 1. $$

Hyperbolic Case ($z > 0$): major semi-axis $A$, minor semi-axis $B = A \sqrt{e^2 - 1}$
$$ (G, C, S, D, T) = \left(\frac{F}{\sqrt{z}}, \cosh F, \frac{\sinh F}{\sqrt{z}}, \frac{\cosh F - 1}{z}, \frac{\sinh F - F}{z \sqrt{z}}\right), \\ (Ξ±, Ξ², Ξ³) = (A(e - 1), A z, B \sqrt{z}), \\ 𝐫 = A (e - \cosh F) 𝐒 + B \sinh F 𝐣, \hspace 1em r = A (e \cosh F - 1), \\ t = t_0 + \sqrt{\frac{A^3}{ΞΌ}} (e \sinh F - F), \hspace 1em 𝐯 = \frac{B \cosh F 𝐣 - A \sinh F 𝐒}{r} \sqrt{\frac{ΞΌ}{A}}, \\ h = B \sqrt{\frac{ΞΌ}{A}}, \hspace 1em H = \frac{ΞΌ}{2A}, \\ B β‰  0 β‡’ \left(\frac{x}{A} - e\right)^2 - \left(\frac{y}{B}\right)^2 = 1. $$

Finally, we have the

Parabolic Case ($z = 0$) (including the recti-linear case, when $Ξ± = 0 = Ξ³$):
$$ (G, C, S, D, T) = \left(G, 1, G, \frac{G^2}{2}, \frac{G^3}{6}\right), \\ Ξ³ = \sqrt{2Ξ±Ξ²}, \hspace 1em e = 1, \\ 𝐫 = \left(Ξ± - Ξ² \frac{G^2}{2}\right) 𝐒 + Ξ³ G 𝐣, \hspace 1em r = Ξ± + Ξ² \frac{G^2}{2}, \\ t = t_0 + \sqrt{\frac{Ξ²}{ΞΌ}} \left(Ξ± G + Ξ² \frac{G^3}{6}\right), \hspace 1em 𝐯 = \frac{-Ξ² G 𝐒 + Ξ³ 𝐣}{r} \sqrt{\frac{ΞΌ}{Ξ²}}, \\ h = \sqrt{2ΞΌΞ±}, \hspace 1em H = 0, \\ Ξ± β‰  0 β‡’ \frac{x}{Ξ±} + \left(\frac{y}{2Ξ±}\right)^2 = 1. $$

  • 2
    $\begingroup$ This answer seems to me to be pretty long, hard to follow and overcomplicating things, it is not as simple and straightforward as the other answer. $\endgroup$
    – Don Al
    Nov 26, 2022 at 1:22
  • 1
    $\begingroup$ The other answers are also wrong: in general, the orbit is not an ellipse. Extra conditions are required to get an ellipse (namely those spelled out in my reply), otherwise it can be a line, parabola or hyperbola. The solution posted is what results from directly solving the problem (which is what we're seeking) and provides the general solution - for all orbit types. It is as simple as you can make the general solution, other than (say) just setting Ξ› = 1 and using dt/dG = r or writing it directly in terms of initial values (𝐫₀, 𝐯₀) = (𝐫,𝐯) @ time t = 0. $\endgroup$
    – NinjaDarth
    Nov 28, 2022 at 3:16

Your Answer

By clicking β€œPost Your Answer”, you agree to our terms of service and acknowledge you have read our privacy policy.