Sunday, August 1, 2010

Celestial Mechanics (II): Energy Constants

We continue looking at celestial mechanics in this instalment, as we examine the role of energy constants in terms of a generic orbit defined by a limited set of orbital elements. (Diagram). First let's look at the problem left at the end of the previous blog ('Legendre Polynomials'). This had to do with the perturbation on Mercury's orbit by Pluto (generic approximation), then the value of 1/DELTA. We use the same functions given in the diagram accompanying the Legendre blog, except now the angle S = 60 degrees. Given the data, r3 = 40 AU for Pluto and r = 0.4 AU for Mercury, the generic term perturbation is easiest to obtain, viz: it will be of order (r/ r3) = 0.4/40 = 0.01. Meanwhile, computations using the same formula in the Legendre diagram yield: 1/DELTA = 0.025

It needs to be emphasized here that what we are working with are raw order of magnitude figures and data. To refine the quantities one must obviously use actual units (e.g. kg, m etc.) not simply ratios of units. But this is not difficult, and is left as a task for the enterprising and industrial reader- who can easily obtain standard masses and distances for the planets, say from an astronomy text, or ephemeris.

But in the case of the illustrated example in the blog text, for the perturbation of Earth by Jupiter, obviously then r = 1 AU = 1.49 x 10^11 m, while r3 =7.45 x 10^11 m and DELTA = 8.29 x 10^11 m. In planetary perturbation theory, the two quantities of interest are:

V = V(r, t) the 'third body' potential function, and

grad V (e.g. dV/dr, or more accurately the partial of V over r) which is the gradient of the potential function.

More exactly, one has:

V(r,t) = - Gm3[ 1/DELTA - r*r3/ (r3)^3]


grad V = Gm3 [(r - r3)/ DELTA^3 + r3/ (r3)^3]

where G = 6.7 x 10^-11 N-m^2/kg^2 is the Newtonian gravitational constant.

On working through with proper units, one finds:

V(r,t) = -1.93 x 10^5 J/kg and grad V = (9.64 x 10^-8 J/kg-m) for the perturbation effect of Jupiter on Earth. I'll leave it to industrious readers to work out the same for the perturbation of Pluto on Mercury.

Energy constants in celestial mechanics are very useful for quickly coming to terms with specific properties of an orbit such as shown in the accompanying sketch- designating a generic orbit in x-y-z space. In the diagram, w is the argument of the perihelion, U is the longitude of the ascending node, f is the true anomaly and i is the inclination of the orbit. The critical or key parameter here is h, the angular momentum vector for the orbiting system. It may be useful here to refer to the diagram (b) in the figure used for angular momentum vector (z) in an atomic system:

Note the vertical direction and compare it with the direction for h(z) in the accompanying diagram here.

Getting specific, assuming r and r' are r (radius vector) and d r/dt, respectively, the magnitude h, of the angular momentum vector is:

h = r x r’ =

(y z’ - z y’)

(z x’ - x z’) = (C1 C2 C3)

(x y’ - y z’)

so (r x r’) = (C1/ h, C2/ h, C3/h)

and inserting variables one finds:

C1/ h = sin U sin (i)

C2/ h = - cos U sin (i)

C3/h = cos(i)

Now since (i) is known (23.5 deg) and therefore cos(i) can be determined, then sin(i) can be as well.h can be determined, since: h = C3 / cos(i) = (GMm a (1 – e^2)^1/2 where all the constants are known (a = semi-major axis of orbit, e = eccentricity of orbit)

We also know: h = (C1^ 2 + C2^2 + C3^2)^1/2 = [C^2]^ 1/2

and we can take:

(C1/ h)/ (C2/ h) = sin U sin(i)/ [- cosU sin(i)]

or C1/ C2 = - tan U

We basically already know, from the above (and using some basic algebra), that:

k (const.) = (C1^2 + C2^2)^1/2

Also, U = W - w (difference between mean anomaly and argument of the perihelion) where W can be obtained from a table based on observations, and w can be obtained using a Fourier expansion of the mean anomaly M, e.g. w = M + (2e – e^3/ 4) sin M + 5 e^^2/4 sin 2M + ... etc.

Once U is known, C1 and C2 and C3 are known.

Let's work out the energy constant for planet Earth, given an eccentricity e = 0.016, and a = 1.0 AU = 1.49 x 10^11 m. Then:

(e^2 - 1) = 2a(e^2 - 1) C/u

and one can find to good approximation (since e ~ 0):

C = u/ 2a

where u = G(m1 + m2) (Using the same mass symbols as in the Legendre problem)

So: C = G(m1 + m2)/ 2a , and h = +/- {ua(1 - e^2)}^1/2
~ +/- [G(m1 + m2)a]^1/2
Problem: Estimate the energy constant, C, for Jupiter, and also its orbital angle of inclinaton, i, if a = 5.2 AU, and e = 0.048

No comments: