# Flying in Formation: The Orbital Dynamics of LISA’s Three Spacecraft

### A.  Introduction

The era of gravitational wave (GW) observational astronomy began on September 14, 2015, when LIGO’s (Laser Interferometer Gravitational Wave Observatory) two giant Michelson-like interferometers detected the gravitational radiation emitted as two black holes, with total mass of $\boldsymbol{70\, M_{\textup{Sun}}}$, collided and coalesced more than 1 billion light years from Earth. [Ref. 1]  This landmark discovery was followed by six others by LIGO and its sister interferometer VIRGO, located near Pisa, Italy.  In each case, the Earth-bound detectors captured a short burst of gravitational radiation emitted in the final seconds before two astronomical bodies — either black holes or neutron stars — collided and merged.  In those fleeting moments, the GW power radiated by the two colliding bodies was comparable to the light emitted by the entire visible universe!

Although GW astronomy is still in its infancy, scientists are already planning the next generation of GW observatories.  The huge LIGO detectors, with interferometer arms stretching for 4 km, are exquisitely sensitive, yet they are incapable of detecting other important astronomical events such as the merger of massive black holes, those with mass  $\boldsymbol{> 10^{4}\, M_{\mathrm{Sun}}}$.  The highest frequency of a GW emitted during the merger of two bodies is inversely proportional to the bodies’ total mass, and LIGO’s low frequency sensitivity is limited by its size and residual coupling to the low-pitched rumble of Earth’s seismic activity.  This is where LISA comes in.  Anticipated to launch in 2034, LISA (Laser Interferometer Space Antenna) is a much larger space-based GW observatory employing three identical spacecraft equipped with lasers, mirrors, telescopes and detectors to form a giant interferometer.  The spacecraft will be deployed about 0.35 AU from Earth and will fly in an equilateral triangular formation of side length $\boldsymbol{2.5\times 10^{6}\: \textup{km}}$ (Figure 1. )  Due to its enormous size and immunity to Earth’s vibrations, LISA will be sensitive to GW radiation at much lower frequencies than LIGO or VIRGO, thus complementing the Earth-bound instruments. [Ref. 2]

LISA is expected to gather data for 2-10 years.  To function for this length of time, its three spacecraft must maintain their separation and configuration with minimal assistance from on-board ion thrusters.  But simply injecting them into orbit with the desired spacing and orientation does not guarantee that they will remain so.  After injection, each spacecraft will move independently in its own elliptical orbit about the Sun, and unless the three orbits are chosen carefully, the configuration of the spacecraft will change unacceptably with time.

As an extreme example, consider the fleet of four spacecraft shown in Figure 2.  The “horizontal” and “vertical” spacecraft pairs represent two orthogonal arms of a GW interferometer, which are initially of length $\boldsymbol{L}$.  The four bodies lie in a plane that is inclined by angle φ to the plane of the ecliptic, and the center of the fleet moves about the Sun in a circular orbit of radius $\boldsymbol{R=1\: \textup{AU}}$.  Let $\boldsymbol{\phi =90^{\circ}}$ and assume $\boldsymbol{L\ll R}$.  (For LISA, $\boldsymbol{L/R=0.017}$.)  The distance of each spacecraft to the Sun is $\boldsymbol{r=\sqrt{R^{2}+L^{2}/4}\simeq R(1+L^{2}/8R^{2})}$, so to first order in $\boldsymbol{L/R}$$\boldsymbol{r=R}$ for each spacecraft.  Naïvely, it would seem that if all four spacecraft were to travel in circular orbits with period $\boldsymbol{T=1\: \textup{yr}}$, their configuration would be stable.  Clearly, if spacecraft C and D travel within the ecliptic plane, they would remain separated by $\boldsymbol{L}$ as they circle the Sun.  But the orbits of A and B are tilted from the ecliptic by the small angle $\boldsymbol{\delta =\pm L/2R}$, and after time $\boldsymbol{T/4=3\: \textup{mo}}$, they intersect and — in the worst case — the two spacecraft would collide!  A more physics-informed strategy is needed for a successful space-based interferometer. [Ref. 3]

In the following sections, we examine how the spacing and orientation of a fleet of spacecraft evolve with time, to understand the strategy underlying LISA’s design.  We begin by presenting an approximate treatment of motion in elliptical orbits, to derive expressions for the position of a spacecraft and the separation between spacecraft, as functions of time.  Using these expressions, we then carry out a straightforward numerical calculation of spacecraft position and separation vs. time, which will elucidate the orbital design adopted by LISA.  All of this can be accomplished using the elementary mathematics and physical principles found in the typical introductory calculus-level physics course.  Because of this and the current high interest in gravitational radiation, LISA offers instructors an ideal opportunity to illustrate basic orbital mechanics in an intriguing way.  The purpose of this post is to show one way to present this topic to undergraduates.

### B.  Position vs. Time in Elliptical Orbits

Each of LISA’s three spacecraft will revolve about the Sun in an elliptical orbit with period $\boldsymbol{T=1}\: \textup{yr}$, so each orbit must have a semi-major axis $\boldsymbol{a=1\: \textup{AU}}$.  If an orbit’s eccentricity $\boldsymbol{\epsilon \ll 1}$, as is the case for LISA, it can be approximated as a circle with the Sun offset from its center by a distance $\boldsymbol{\epsilon R}$, where $\boldsymbol{R}$ is the radius of the circle and also the semi-major axis of the ellipse ($\boldsymbol{R=a}$) (Figure 3).  [Ref. 4]

The position $\boldsymbol{\vec{r}(t)}$ of the orbiting body, relative to the Sun, is specified by its distance to the Sun $\boldsymbol{r}$ and the angle ν, called the true anomaly in astronomical literature.  Similarly, the position $\boldsymbol{\vec{R}(t)}$ of the body relative to the center of the circle is given by $\boldsymbol{R}$ and the angle θ, called the eccentric anomaly in the literature.  Let $\boldsymbol{t=0}$ and $\boldsymbol{\nu =\theta =0}$ when the body is at perihelion, its closest approach to the Sun.  The area swept out by $\boldsymbol{\vec{r}(t)}$ from perihelion to time $\boldsymbol{t}$ is equal to the circular sector swept out by $\boldsymbol{\vec{R}(t)}$ minus the shaded triangular area shown in Fig. 3.  Since the total area swept out by $\boldsymbol{\vec{r}(t)}$ in a full orbital period is $\boldsymbol{\pi R^{2}}$, Kepler’s 2nd (equal area) Law yields

$\LARGE\boldsymbol{\frac{t}{T}=\frac{\frac{1}{2}R^{2}\theta -\frac{1}{2}\epsilon R^{2}\mathrm{sin}\theta }{\pi R^{^{2}}}}$

or

$\LARGE \boldsymbol{\frac{2\pi t}{T}=\theta -\epsilon\, \mathrm{sin}\theta\equiv M}$,

where $\boldsymbol{M}$ is often called the mean anomaly.  (Note that if $\boldsymbol{\epsilon =0}$$\boldsymbol{M=\theta }$, so $\boldsymbol{M}$ is the angular position of a body moving in a circular orbit with period $\boldsymbol{T}$.)  This is a transcendential equation for $\boldsymbol{\theta (t)}$, and an approximate solution is easily found for $\boldsymbol{\epsilon \ll 1}$.  To first order in $\boldsymbol{\epsilon}$$\boldsymbol{\theta =M}$, so

$\large \boldsymbol{\theta(t) = M+\epsilon\, \mathrm{sin}M}$.

(1 )

Conversely, if the orbiting body is at aphelion at $\boldsymbol{t=0}$, then the Sun and the body are on opposite sides of the circle center at $\boldsymbol{t=0}$, and the shaded area in Fig. 3 must be added to the circular sector to find the area swept out by $\boldsymbol{\vec{r}(t)}$.  In this case,

$\large \boldsymbol{\theta (t)=M-\epsilon \, \mathrm{sin}M}$.

(2)

Equations 1 and 2 specify the position of the orbiting body as a function of time, relative to the center of the offset circle.  Both equations are needed below, as we will study the motion of two bodies that are either at perihelion or aphelion at $\boldsymbol{t=0}$.  Knowing $\boldsymbol{\theta (t)}$, it is straightforward to calculate the body’s position relative to the Sun.  Keep in mind that these results, and all of the calculations to follow, are accurate only to first order in $\boldsymbol{\epsilon }$.

### C.  Maintaining a Constant Separation between Spacecraft

Consider just the two “vertical” spacecraft A and B shown in Fig. 2.  Let the line $\boldsymbol{\overline{AB}}$ joining them be tilted by angle φ from the plane of the ecliptic, and suppose that the midpoint of $\boldsymbol{\overline{AB}}$ travels in a circular path about the Sun.  If $\boldsymbol{\phi < 90^{\circ}}$, then A and B are no longer equidistant from the Sun.  At the instant shown in Fig. 2, $\boldsymbol{r_{\mathrm{A,B}}=R\pm L\mathrm{cos}\phi /2}$.  To have a period of 1 yr, both spacecraft must move in elliptical orbits with semi-major axis $\boldsymbol{a=R=1\: \textrm{AU}}$.  If their orbital velocities are perpendicular to their position vectors at $\boldsymbol{t=0}$, A and B must be at aphelion and perihelion, respectively.  The aphelion and perihelion distances for an elliptical orbit are $\boldsymbol{r_{a,p}=a(1\pm \epsilon )}$, so $\boldsymbol{\epsilon =L\mathrm{cos}\phi /2R}$ for either orbit.  In addition, both orbits are inclined from the ecliptic plane by $\boldsymbol{\delta =\pm L\mathrm{sin}\phi /2R}$.  ($\boldsymbol{\left | \delta \right |< 0.5^{\circ}}$ for LISA).

Figure 4 shows the elliptical orbits of A and B as circles offset from the Sun by $\boldsymbol{\pm \epsilon R}$.  The two orbits still intersect, but now B passes through the intersection point ahead of A.  From Eqns. 1 and 2, when $\boldsymbol{M=\pi /2}$$\boldsymbol{\theta _{\mathrm{B}}-\pi /2=\pi /2-\theta _{\mathrm{A}}=\epsilon }$, so the separation s between the two spacecraft is

$\boldsymbol{s=2\epsilon R+R\, \mathrm{sin}(\theta _{\mathrm{B}}-\pi /2)+R\, \mathrm{sin}(\pi /2-\theta _{\mathrm{A}})\simeq 4\epsilon R=2L\mathrm{cos}\phi }$,

where we have used the small angle approximation $\boldsymbol{\mathrm{sin}\epsilon \simeq \epsilon }$ to simplify the expression.  Since the difference $\boldsymbol{\theta _{\mathrm{B}}-\theta _{\mathrm{A}}=2\epsilon\, \mathrm{sin}M}$, we expect $\boldsymbol{\left | s \right |}$ to oscillate between minima at $\boldsymbol{M=0,\pi }$ and maxima at $\boldsymbol{M=\pi /2,3\pi /2 }$ as the two spacecraft orbit the Sun.  This will be borne out in the following section.  But note that for $\boldsymbol{\phi =60^{\circ}}$, $\boldsymbol{s=L}$ at $\boldsymbol{M=0}$ and $\boldsymbol{M=\pi /2}$, suggesting that this may be the optimal angle of inclination for the plane of the space-based interferometer.

### D.  Numerical Calculation of Spacecraft Separation

Eqns. 1 and 2 allow us to calculate numerically the coordinates of spacecraft A and B, and find their separation, throughout their circumsolar orbits.  The calculations can be carried out using a common spreadsheet, such as Excel, which should make it easy for students to interpret and duplicate them.  First consider spacecraft A.  It is at aphelion at $\boldsymbol{t=0}$, and the plane of its orbit is tilted from the ecliptic plane by the small angle $\boldsymbol{\delta=L\, \mathrm{sin}\phi /2R}$.  Begin the spreadsheet by generating a column of values for the mean anomaly $\boldsymbol{0\leq M\leq 2\pi }$.  Next, use Eqn. 2 to calculate a corresponding column of values $\boldsymbol{\theta _{\mathrm{A}}=M-\epsilon\, \mathrm{sin}M}$, where $\boldsymbol{\epsilon =L\, \mathrm{cos}\phi /2R}$.  The coordinates of A, relative to the Sun, are then found using

$\large \boldsymbol{x_{\mathrm{A}}(t)=R(\mathrm{cos}\theta _{\mathrm{A}}+\epsilon )\mathrm{cos}\delta }$

$\large \boldsymbol{y_{\mathrm{A}}(t)=R\, \mathrm{sin}\theta _{\mathrm{A}}}$

$\large \boldsymbol{z_{\mathrm{A}}(t)=R(\mathrm{cos}\theta _{\mathrm{A}}+\epsilon )\mathrm{sin}\delta }$.

(3)

Spacecraft B is at perihelion at $\boldsymbol{t=0}$, so $\boldsymbol{\theta _{\mathrm{B}}(t)=M+\epsilon \, \mathrm{sin}M}$.  Its orbit in inclined by $\boldsymbol{-\delta }$, so its coordinates are given by

$\large \boldsymbol{x_{\mathrm{B}}(t)=R(\mathrm{cos\mathrm{\theta _{\mathrm{B}}}}-\epsilon )\mathrm{cos}(-\delta )}$

$\large \boldsymbol{y_{\mathrm{B}}(t)=R\, \mathrm{sin}\theta _{\mathrm{B}}}$

$\large \boldsymbol{z_{\mathrm{B}}(t)=R(\mathrm{cos\mathrm{\theta _{\mathrm{B}}}}-\epsilon )\mathrm{sin}(-\delta )}$.

(4)

After compiling the positions of A and B vs. time, their separation can be computed: $\boldsymbol{s(t)=\sqrt{(x_{\mathrm{A}}-x_{\mathrm{B}})^{2}+(y_{\mathrm{A}}-y_{\mathrm{B}})^{2}+\cdots }}$, and is shown in Figure 5 for tilt angles $\boldsymbol{\phi =0^{\circ},15^{\circ},\cdots ,90^{\circ}}$.  As anticipated, $\boldsymbol{s(t)=L}$ independent of time for $\boldsymbol{\phi =60^{\circ}}$.

For additional insight, it is useful to switch to a primed coordinate system rotating uniformly with angular speed $\boldsymbol{\Omega =2\pi /T}$ (so $\boldsymbol{\Omega t=M}$).  The $\boldsymbol{{x}'}$ axis points from the Sun through the midpoint of line $\boldsymbol{\overline{AB}}$, the $\boldsymbol{{y}'}$ axis points in the direction of increasing $\boldsymbol{\theta }$, and the $\boldsymbol{{z}'}$ axis is perpendicular to the plane of the ecliptic.  Transforming to the rotating coordinate system [Ref. 5],

$\large \boldsymbol{{x}'=x\, \mathrm{cos}M+y\, \mathrm{sin}M}$

$\large \boldsymbol{{y}'=-x\, \mathrm{sin}M+y\, \mathrm{cos}M}$

$\large \boldsymbol{{z}'=z}$.

(5)

Furthermore, by subtracting $\boldsymbol{R}$ from $\boldsymbol{{x}'}$, we can place the origin of the rotating system at the midpoint of line $\boldsymbol{\overline{AB}}$, which follows a circular path on the ecliptic.  Figure 6 shows the result for spacecraft A for $\boldsymbol{\phi =60^{\circ}}$.  Note that $\boldsymbol{{x}'-R}$,  $\boldsymbol{{y}'}$ and $\boldsymbol{{z}'}$ vary sinusoidally with period 1 yr, and that the value of $\boldsymbol{\sqrt{({x}'-R)^{2}+{y}'^{2}+{z}'^{2}}}$ is constant and equal to $\boldsymbol{L/2}$.  This indicates that as A orbits the Sun, it also executes a circular “epicycle” of radius $\boldsymbol{L/2}$ about the midpoint of $\boldsymbol{\overline{AB}}$, moving clockwise in the plane inclined by φ from the ecliptic.  Spacecraft B moves similarly, but occupies a position diametrically opposed to A on the same epicycle, so the separation between A and B is equal to $\boldsymbol{L}$ at all times.  If $\boldsymbol{\phi \neq 60^{\circ}}$, A and B still execute one epicycle per year, but it is elliptical rather than circular, with semi-axes $\boldsymbol{L/2}$ and $\boldsymbol{L\, \mathrm{cos}\phi }$.  In this case, the separation of the two spacecraft varies with time.

### E.  Maintaining the Shape of the Constellation

Imagine a fleet of n spacecraft arrayed about a central point that moves in a circular path about the Sun.  How should we choose their orbits so that their spacing and configuration are stable?  Generalizing from the preceding discussion, we want all of the spacecraft to move synchronously about a circular epicycle when viewed in the rotating coordinate system.  This is guaranteed if the fleet is embedded in a plane tilted from the ecliptic plane by $\boldsymbol{\phi =60^{\circ}}$.  If the spacecraft occupy the vertices of a regular polygon (equilateral triangle, square, pentagon, etc.), then each will pass through aphelion (or perihelion) at intervals $\boldsymbol{\Delta t=T/n}$ or $\boldsymbol{\Delta M=2\pi /n}$.  For example, in the preceding example of two spacecraft, B and A arrive at aphelion six months apart ($\boldsymbol{\Delta M=\pi }$).  Equally important, B’s orbit is the same as A’s, but rotated by 180° about the z-axis (perpendicular to the ecliptic plane).  So for n spacecraft, all orbits should have exactly the same shape and be inclined from the ecliptic by the same angle δ, but be rotated with respect to one another by $\boldsymbol{2\pi /n}$.

To illustrate, let us return to the four spacecraft shown in Fig. 2, and assume that $\boldsymbol{\phi =60^{\circ}}$.  Rather than following the ecliptic (as assumed in Section A), the orbits of C and D should have the same shape and tilt as A’s, but should be rotated about the z-axis by $\boldsymbol{\pm 90^{\circ}}$.  This is shown in Figure 7.  If A passes through aphelion at time $\boldsymbol{M=0}$, then C will pass through its aphelion at $\boldsymbol{M=\pi /2}$, while D has already done so at $\boldsymbol{M=-\pi /2}$.  At the instant shown in Fig. 2 ($\boldsymbol{M=0}$), C is moving upward through the ecliptic plane toward its aphelion, whereas D is moving downward.  Let’s check this by calculating the x and y positions of C and D at $\boldsymbol{M=0}$, which should be $\boldsymbol{x=R}$$\boldsymbol{y=\pm L/2}$.  From Eqn. 2, during the time interval $\boldsymbol{0< M< \pi /2}$, each spacecraft moves through an angle $\boldsymbol{\Delta \theta =\pi /2-\epsilon }$, measured from the center of its (offset) circular orbit.  Ignoring the small tilt of the orbits, the x and y coordinates (see Fig. 7) of C and D at $\boldsymbol{M=0}$ are, to first order in $\boldsymbol{\epsilon }$$\boldsymbol{x=R\, \mathrm{cos}\epsilon \simeq R }$ and $\boldsymbol{y=\pm (\epsilon R+R\,\mathrm{sin}\epsilon )\simeq \pm 2R\epsilon=\pm L/2 }$, as expected. [Ref. 7]

We conclude that to maintain its size and triangular configuration, LISA’s three spacecraft must lie in a plane tilted by 60° from the ecliptic plane, and must follow identical orbits rotated by 120°.  In fact, any number of spacecraft lying in this plane in any configuration will move synchronously if their orbits have the same semi-major axis a, eccentricity ε, and inclination δ, and in addition, are suitable rotated with respect to one another.  One should remember, though, that these results are valid only to first order in ε (or $\boldsymbol{L/R}$).  Exact solution of the equations of motion for LISA shows that the spacecraft separation will still vary by ± 1%. [Ref. 7]

### F.  Questions for Students

1. Draw a figure similar to Fig. 3, and derive Eqn. 2 for the position of a body that passes through aphelion at $\boldsymbol{t=0}$.

2.  A body is in an elliptical orbit about the Sun with period 1 yr and eccentricity  $\boldsymbol{\epsilon =0.10}$.  (a) What is the approximate eccentric anomaly θ two months after perihelion?  (b) Find the true anomaly ν and the distance r to the Sun at this time.  (Ans: (a) $\boldsymbol{M=60^{\circ}}$, so $\boldsymbol{\theta =65^{\circ}}$ (b) $\boldsymbol{\nu =70.4^{\circ}}$$\boldsymbol{r=0.96\: \textrm{AU}}$)

3. Use Fig. 1 to show that the eccentricity of LISA’s orbits is $\boldsymbol{\epsilon =L/(2R\sqrt{3})\simeq 0.005}$.

4. How can an L-shaped (as opposed to an equilateral triangle) fleet of three spacecraft be placed in an Earth-trailing orbit such that the spacing and configuration of the spacecraft remain stable?  (Ans: $\boldsymbol{a=1\: \textrm{AU}}$ for all three spacecraft,  which lie in a 60° plane; let $\boldsymbol{\epsilon =0}$ for the center craft, which follows the ecliptic around the Sun; $\boldsymbol{\epsilon =L/4R}$  for the other two, whose orbits are tilted by $\boldsymbol{\delta =\sqrt{3}L/4R}$ and are rotated by 90° with respect to each other.)

5. A fleet of five spacecraft is to be inserted in a circular orbit between Mars and Jupiter, 2.5 AU from the Sun.  The spacecraft are to remain in a regular pentagonal formation with minimal use of on-board rocket fuel.  (a) Explain how this can be done.  (b) One of the spacecraft passes through perihelion at $\boldsymbol{t=0}$.  How long before another passes through perihelion?  (Ans: $\boldsymbol{T=4\: \textrm{yr}}$, so $\boldsymbol{\Delta t=T/5=9.6\: \textrm{mo}}$.)

6. Suppose the plane of the fleet were inclined from the ecliptic by $\boldsymbol{\phi =120^{\circ}}$ (e.g., in Fig. 2, A is now closer to the Sun than B).  Would the configuration of the spacecraft be stable?

### Acknowledgments

We thank Jonathan Levine and Francesco Meschia for suggesting numerous ways to improve the clarity of the physics in this post.  The author is responsible for any remaining errors.

### References

1.  B. P. Abbott et al, “Observation of Gravity Waves from a Binary Black Hole Merger,”  Phys. Rev. Lett. 116, 061102 (2016)

2.  “LISA: Laser Interferometer Space Antenna: A cornerstone Mission for the Observation of Gravitational Waves,” ESA-SCI (2000) 11 July, 2000

3.  A nice treatment of a similar problem, by similar methods, is given by E. I. Butikov, “Relative Motion of Orbiting Bodies,” Am. J. Phys. 69, 63-67 (2001)

4.  This approximation was first used by ancient Greek astronomers to describe the motion of the Moon and Sun within the Ptolemaic (geocentric) system of the world.  See A. Berry, “A Short History of Astronomy,” (1898), available from Dover Publications.

5.  For an easy derivation of these transformation equations, see J. C. Amato and E. J. Galvez, Physics from Planet Earth, Chapter 2, problem 2.18 (Taylor & Francis/CRC Press, 2015).

6.  This topic is treated on a higher mathematical level by S. V. Dhurandhar et al, “Fundamentals of the LISA Stable Flight Formation,” Class. Quantum Grav. 22, 481-487 (2005).

7.  LISA Phase 0 CDF Study – Internal Final Presentation, ESTEC (European Space Research and Technology Centre) Conference, Noordwijk, Netherlands, March 8 – May 5, 2017,  slide 7, p. 110  http://sci.esa.int/lisa/59336-lisa-internal-phase-0-cdf-study-final-presentation/#