Lagrange Points

Program 1/2: RMlpoints.z80

The Program

This is a heavily-modified form of the Orbits program. The idea first occurred to me when I was writing that program, although I didn't take it any further at the time.

There are now 2 idealized bodies exerting a gravitational force on our spacecraft - the earth and the moon (lowercase - this isn't supposed to be a realistic representation). The moon follows a circular, anticlockwise, orbit. Its starting point is always (1, 0). In the program, it has a mass 1% that of the earth's. The time interval is set to 0.001.

The idea is that the craft should be placed at a Lagrange point with the same velocity as that of the L point. Then we leave the simulation to see how stable the craft's position is over time, relative to the L point.

When calculating the acceleration experienced by the spacecraft, we now have to add in the contribution from the moon and to do that we first need to calculate the position of the moon in its orbit.

That's not all though. Moon and earth are both orbiting about their centre of mass, the so-called barycentre. This is now the centre of our coordinate system. The earth's centre of mass is on the opposite side of the barycentre to the moon. If the moon's CoM is at one unit from the barycentre, then the earth's is at 0.01 units in the other direction. It's a small effect, but significant. The earth's starting point now is therefore always (-0.01, 0) and we have to calculate the position of the earth at each instant as well, before determining the pull due to the earth on the craft.

The craft can be placed at one of 4 L points:

(L3 located on the opposite side of the earth from the moon has been left for now).

Thrusters are not provided!

3 frames of reference (Centres of View) are provided:

Placing an object at L4

The image above shows what happens when a craft is placed at the L4 point, 60 degrees ahead of the moon, with orbital velocity equal to that of the moon, but rotated by 60 degrees. After one orbit, it is back at its starting point (within 1 percentage point, which I put down to accumulated error in the simulation). Its orbit is the same as the moon's and is indeed stable, at least over the interval shown.

Note: this program runs very slowly indeed and I've cranked up the emulation speed in Fuse by a (maximum) factor of 50 to have it run in a reasonable time. This is cheating of course and I really should consider re-coding in z80 Assembly language.

Placing an object far from L4

However, to confirm stability, we really need to compare it with some other points on the moon's orbit.

First, let's look at what happens when the craft is placed just 30 degrees beyond the position of the moon, with orbital velocity equal to that of the moon, but rotated by 30 degrees. (I hacked the code to achieve this - there is no L9!). The first image below shows the start and end points over one orbit do not meet - the orbit at this point is indeed not stable.

The second image below shows what happens when the craft is now placed initially at 90 degrees beyond the position of the moon, with orbital velocity equal to that of the moon, but rotated by 90 degrees. As it is further away from the moon, we might expect it to be more stable than the previous example and this does seem to be the case, although note that this time, the craft ends up beyond the moon's orbit. However, even though further from the moon than at the L4 point, the effect of the moon's gravity on the craft can be seen here - stability at L4 does indeed seem to be confirmed!

I would expect to see similar behaviour at the L5 point.

Interestingly, this article says that L4 and L5 are stable if the major body's mass is greater than about 25 times that of the smaller body's. That means for the points to be stable in the simulation, the moon's mass must be below about 0.04 times that of the earth's. I chose 0.01 (the actual figure for the moon is about 0.0123) so these points should be stable and that is the case as has been seen (for L4). I've tried a mass of 0.05, and now the L4 point is indeed unstable. However, I also see evidence of instability at 0.02! Perhaps changing the time interval to 0.0001 would reduce or eliminate the discrepancy?

Placing an object close to L4

I've added some variables to allow for simulating the placing of an object near an L point. These are hard-coded in the program, which therefore needs to be edited in order to change them:

(I may add prompts for these values in future).

Below are some screenshots of behaviour near L4, with dri=2 and various values for dvi. It seems that stable paths about L4 are possible for positive values of dri, if the velocity of the object with respect to the barycentre is the same as, or somewhat less than that of the L point. Note that in the 4th image, dvi =-2, the magnification factor is 8 times, compared with just two times in the other images: the object stays very close indeed to the L4 point at this distance and velocity.

dri=2, dvi=4, umag=1

dri=2, dvi=2, umag=2 (same angular speed as L point)

dri=2, dvi=0, umag=2 (same velocity as L point)

dri=2, dvi=-2, umag=8

dri=2, dvi=-4, umag=2

dri=2, dvi=-6, umag=2

That 4th case - dvi =-2, 8 times magnification - is again shown below, but this time in the rotating frame of the L4 point, and over a larger period of time, equivalent to about 3 orbits by the L4 point about the barycentre. Note the two large lobes.

dr=2, dv=-2, umag=8, CoV=3

Placing an object at L1, L2

The paths of objects placed at the L1 and L2 points (in this case at a point on the earth-moon line 0.151 before and beyond the moon respectively) are not stable, as shown below (each in the rotating frame of its L point). However, for L2, the path is shown magnified 4 times and we can see that the object stays close to the L point for about half an orbit of the L point about the barycentre.

Program 2/2: RMlpointsz3d.z80

Three dimensions

a work in progress...

So far, I've limited the program to looking at motion in a plane containing spacecraft, earth and moon. Apparently, interesting things happen in 3 dimensions, with semi-stable halo orbits about L1 and L2 and this is an attempt to write code to investigate those.

Moon and earth are still in the (x, y) plane and I have the object start there too, but now with a velocity with a z component (vzi, hard-coded near the start of the program). The existing code was extended easily enough to compute acceleration, position and velocity in the z direction. In addition, there is now a choice of view planes - (x, y), (y, z), ( z,x) - specified with the pv variable, values 1, 2 or 3 respectively.

To test the program, shown below are three views of the path of an object placed at L4 with vzi=0.2. All show the barycentre as the centre of view. The first shows the path in the (x, y) plane, the second in the (y, z) plane, and the third in the (z,x) plane. 

To squeeze in the details about the view plane, I've changed the text slightly to show view details (V) of cov/pv.

Interestingly, looking at the first and second views over the time of about one lunar orbit, we can see the path of the object is no longer closed. The third view doesn't show this, presumably as a result of perspective. Note too that the loop in the second image is wider than that in the second.

I've enhanced the program further still to allow for the creation of a stereoscopic pair of views, to make it easier to visualize the path. This is achieved by having an image rotated by a small angle anticlockwise about the vertical axis of the view plane. For simplicity, I have used a Pythagorean triple, so that sine and cosine functions don't need to be called. The triple is 21, 220, 221, for a sine of 21/221 (ie approximately 0.1), and cosine of 220/221.

However, the program will need some minor editing to implement this, as the line that does this is commented out by default and there is instead another, active, line where the angle is effectively 0, the sine being 0, the cosine 1.

As a result of this rotation, distant points move to the left and near points to the right and so this image corresponds, I think, to a view for the left eye.

The pair of pictures below again shows the path of an object at L4, with vzi=0.2, but this time I've chosen the L point as the centre of view and (z, x) as the view plane. IIRC, magnification was 2. The object starts at the beginning of the straight line segment. The second image has been rotated slightly as described above, and hence is intended for the left eye. The pair can be viewed, if you can see them side by side on a larger screen, using the cross-eyed method, which is easier to achieve than it may sound! The 3D effect works fairly well, I think.

All this extra processing does mean the program will run even more slowly now! Other than that, all looks OK so far, so at some point I will try to generate halo orbits for L1 and L2.