The Last Days of Cassini: the Grand Finale

cassini_april27release-b
Figure 1:  Artist’s conception of Spacecraft Cassini threading the gap between Saturn’s cloud tops and its innermost ring.  The gap is only 3000 km wide.  (Image courtesy NASA/JPL)

A.  Introduction

A central premise of Physics from Planet Earth is that students will be more motivated to study physics if their coursework addresses important and newsworthy scientific events.  A good example is the Cassini-Huygens mission to Saturn, which ended dramatically on September 15, 2017 when the spacecraft dipped into the Saturnian atmosphere and self-destructed — its grand finale.  Cassini‘s fiery demise was carefully planned by NASA seven years ago (Ref. 1), when the spacecraft’s fuel supply began to run low, limiting opportunities for further exploration of Saturn and its surroundings.  Flight engineers chose to end the mission by crashing into Saturn to eliminate any possibility of the spacecraft contaminating Enceladus or Titan, two moons of Saturn that potentially harbor life.

Cassini was launched from Cape Canaveral 20 years ago.  After 4 gravity assists (2 from Venus, 1 from Earth, and 1 from Jupiter), it entered orbit about Saturn in July 2004, the first spacecraft ever to do so.  The ESA’s (European Space Agency) Huygens spaceprobe was dispatched from Cassini on Christmas day, 2004, and soft-landed on Titan three weeks later, the most distant landing ever achieved by a spacecraft from Earth.  Over the next 13 years, Cassini discovered 7 new moons of Saturn (bringing the total to 62) and conducted close flybys of 13 moons.  To visit them, it employed gravity assists, or “targeted flybys,” to alter its flight trajectory with minimal expenditure of fuel.  Most of these flybys — 127 in all — involved Titan, Saturn’s largest moon and the only one with sufficient mass to alter Cassini‘s trajectory significantly.  By taking advantage of gravity assists, NASA flight engineers were able to extend the mission by 9 years.  Their efforts were amply rewarded.  Together, Cassini and Huygens revolutionized our understanding of the Saturnian system.  In addition to providing stunning close-up images of the planet, its rings and moons, the two spacecraft discovered subsurface water oceans on Titan, Enceladus and Dione, surface oceans of methane (CH4) on Titan, and gaseous plumes  containing H2O, CO2, CH4 and H2 emanating from Enceladus, making the latter the most promising candidate in our solar system for hosting extraterrestrial life.  

Prior to its “grand finale” orbits, Cassini used a targeted flyby of Titan (designated T125) to enter a “ring-grazing” orbit, with periapse just outside Saturn’s F ring.  After completing 20 ring-grazing orbits, the spacecraft executed its final targeted flyby (T126) of Titan on April 22, 2017, and was inserted into a daring trajectory that — near its periapse — passed between the cloud tops of Saturn and its innermost ring, a gap only 3000 km wide.  Over the next 5 months, Cassini executed 22 grand finale orbits, until a high altitude flyby of Titan on September 12 caused the spacecraft to plunge into Saturn’s atmosphere three days later.

In this post, we focus on T126 to learn how Cassini was transferred from its last ring-grazing orbit to the first of its grand finale orbits.  This is done via a set of homework exercises and/or exam problems suitable for first year university physics students.  We gratefully acknowledge the generous help of Dr. Duane Roth, Cassini Navigation Team Chief, who kindly forwarded Reference 1 and provided detailed instructions for using the online software package WebGeoCalc (Ref. 2) to calculate the trajectory of Cassini relative to Saturn and Titan, both before and after T126.  We drew heavily WebGeoCalc’s calculations to write this post.  Thank you, Dr. Roth!

Cassini‘s grand finale was a headline-grabbing event, and perfectly timed for inclusion in the fall semester syllabus of an introductory mechanics course.  We hope your students will enjoy learning the physics behind the headlines.  Farewell, Cassini!  Right up to your final moments, you performed brilliantly!

 

B.  The T126 flyby

In the exercises below, the T126 flyby should be treated as an elastic “collision” between Cassini and Titan.  The following information will be useful for solving the problems.

Titan is in a near-circular orbit (ε = 0.0288) about Saturn, well-aligned to the planet’s equatorial plane, with period \boldsymbol{P_{T}=1.379\times 10^{6}} s (16.0 d).  Using Kepler’s 3rd law, its semi-major axis is found to be \boldsymbol{a_{T}=1.222\times 10^{6}} km.  Prior to the T126 flyby, Cassini was in a highly eccentric orbit (ε = 0.795) about Saturn with inclination about 60° relative to the equatorial plane, period \boldsymbol{P_{C}=6.205\times 10^{5}} s (7.2 d) and semi-major axis \boldsymbol{a_{C}=7.178\times 10^{5}} km.  The flyby occurred at a distance \boldsymbol{r=1.206\times 10^{6}} km from Saturn, not far from the spacecraft’s apoapse.  (See Figure 2.)  At the encounter point, the spacecraft crossed the plane of Titan’s orbit from below to above.  As a result of the flyby, Cassini entered an even more eccentric orbit (ε = 0.906) which, at periapsis, passed between the planet and its innermost ring.

Figure 2:  Intersecting orbits of Titan and Cassini.  The plane of the figure is the equatorial plane of Saturn.   Cassini’s orbit is tilted by about 60 degrees.  At the encounter, the spacecraft’s velocity vector points out of the page.  (Image credit:  NASA/JPL)

 

Exercises for students:

1. Titan is in a near-circular orbit about Saturn with period \boldsymbol{P_{T}=1.379\times 10^{6}} s.  The T126 flyby occurred when Titan was \boldsymbol{1.206\times 10^{6}} km from Saturn.  Find its speed at this time.  Use \boldsymbol{GM_{S}=3.793\times 10^{7}\:\mathrm{km^{3}/s^{2}}}.  (Ans: \boldsymbol{v_{T}=5.645\: \mathrm{km/s}})

2.  Just before the flyby, Cassini was \boldsymbol{1.19\times 10^{6}\: \mathrm{km}} from Saturn and moving with speed \boldsymbol{v_{C}=3.291\: \mathrm{km/s}}.  The eccentricity of its orbit was 0.795.  Show that, at periapsis, the spacecraft was just outside Saturn’s F ring, which has a radius of about \boldsymbol{1.40\times 10^{5}\: \mathrm{km}}.

3.  Titan’s orbit is within 0.3° of the plane of Saturn’s rings, whereas Cassini‘s orbit is tilted by about 60°.  As shown in this video,  Entering Final Orbits: Cassini’s Grand Finale (Ref 3), Cassini passed in front of Titan during T126.  (The encounter occurs 20 s from the beginning of the video.)  Let \boldsymbol{\vec{u}} and \boldsymbol{\vec{{u}'}} be the velocities of the spacecraft relative to Titan before and after the flyby.  Make a rough sketch of the velocity vectors \boldsymbol{\vec{v}_{T}}\boldsymbol{\vec{u}} and \boldsymbol{\vec{{u}'}}, and explain why the spacecraft’s speed decreased due to the flyby.  No calculations are necessary.

4.  In a reference frame whose x– axis is parallel to \boldsymbol{\vec{v}_{T}} , and whose z-axis is perpendicular to the moon’s orbit, the velocity of Titan just before the flyby may be expressed as \boldsymbol{\vec{v}_{T}=5.645\, \hat{i}\: \mathrm{km/s}}, and the velocity of Cassini as \boldsymbol{\vec{v}_{C}=1.161\, \hat{i}-0.701\, \hat{j}+2.998\, \hat{k}\: \mathrm{km/s}}.   Calculate the angle between the two vectors and also find the relative velocity vector \boldsymbol{\vec{u}=\vec{v}_{C}-\vec{v}_{T}}.  (Ans: 66°)

5.  As Cassini approached Titan, the angle between the relative velocity vector \boldsymbol{\vec{u}} and \boldsymbol{\vec{v}_{T}} was 145.5°.  The flyby deflected the direction of the relative velocity by 9.1°.  Show that this reduced the spacecraft’s speed relative to Saturn by about 850 m/s.

6.  True or false:  the change in the speed of the spacecraft is equal to \boldsymbol{\left | {\vec{v}}'_{C}-\vec{v}_{C} \right |}, i.e., the absolute value of the change in the velocity of the spacecraft.

7.  At its closest approach to Titan, Cassini was only 980 km above the surface of the moon, or 3554 km from its center.  What was its speed relative to Titan at this time?  (Ans:  \boldsymbol{u_{max}= 5.880 \: \mathrm{km/s}})

8.  The unfueled mass of Cassini was 2125 kg, and the exhaust speed of its MMH/NTO fuel mixture is 2.8 km/s.  What is the minimum mass of fuel that would have been needed to change the spacecraft’s speed (relative to Saturn) by 850 m/s?  Assume that the maneuver totally depleted the fuel supply.  (Ans: 750 kg)

9.  In its new trajectory (after T126), Cassini made 22 full orbits with period \boldsymbol{5.565\times 10^{5}\: \mathrm{s}} before a long distance encounter with Titan caused it to plunge into the Saturnian atmosphere.  Why 22 orbits?  (Hint: compare the periods of Cassini and Titan.)

The following two exercises can be completed after studying Chapter 12 of PPE or similar material.

10.  Prove that in the reference frame of Titan, the spacecraft was deflected by 9.1° due to T126.  (This is the angle between \boldsymbol{{\vec{u}}'} and \boldsymbol{\vec{u}}.  

11.  On September 12, Cassini came within \boldsymbol{1.2\times 10^{5}\: \mathrm{km}} of Titan, and its relative speed at closest approach was 5.2 km/s.  Estimate the angle by which the spacecraft was deflected (in Titan’s reference frame), and the maximum change of speed relative to Saturn.  This change in speed was enough to ensure that the spacecraft entered Saturn’s atmosphere at periapse.

References

1.  Brent Buffington et al, IAC-10-C1.9.2, 61st International Astronomical Congress, Prague, CZ (2010)

2.  http://wgc.jpl.nasa.gov:8080/webgeocalc/#OrbitalElements

3.  Entering Final Orbits: Cassini‘s Grand Finale.  Video is available courtesy of NASA/JPL at https://saturn.jpl.nasa.gov/resources/7508/?category=videos .  

Advertisements

The Range of North Korean Ballistic Missiles

14dc-missile-1-master768
Figure 1:  According to North Korean sources, this image shows the launch of a Hwasong-14 missile, the same model as the one launched on July 28, 2017, which is believed to be capable of reaching the continental United States.

A. Introduction

North Korea’s missile launches of July 2017  demonstrate that it has developed ballistic missiles (Figure 1) with intercontinental range (ICBMs).  This, plus its parallel development of nuclear weapons, has added a grave new threat to world order.  This is a grim topic, but because it has such immediate importance — and encompasses fundamental physical principles — it deserves to be discussed in our introductory physics classes.

The missile launch of July 4 flew for 37 minutes (Ref. 1), and the one launched on July 28 was aloft for 47 minutes (Ref. 2).  Based on this information, analysts in the United States and Japan estimate the maximum range of the former to be 6,700 km, and the latter to be 10,400 km, placing Los Angeles, Denver, Chicago, New York, Boston and other North American population centers within targeting range of Pyongyang.  The purpose of this post is to show, using mathematics and concepts appropriate to the introductory calculus-based mechanics course, how to estimate missile range from the (observable) flight time.  Following a treatment of Kepler’s Laws, it should be possible to present this material in a single class period.

B.  Determining the Burnout Speed of a Ballistic Missile

To simplify the discussion, assume the missile altitude at burnout (when rocket fuel is exhausted) is negligible compared to the Earth’s radius, neglect atmospheric drag, and also neglect the Earth’s rotation.  North Korea uses highly “lofted” trajectories to test its ICBMs; that is, they are launched at much steeper angles (approaching 90º) than they would be to achieve maximum range.  In this section, we will assume that the missile is launched vertically, and calculate its flight time \boldsymbol{t_{f}} and maximum altitude \boldsymbol{h_{max}} as functions of its burnout speed \boldsymbol{v_{0}}.  For small \boldsymbol{v_{0}} and \boldsymbol{h_{max}}, this is an easy task, as the gravitational acceleration \boldsymbol{g} is roughly constant, and the familiar equations of kinematics are applicable: \boldsymbol{t_{f}=2v_{0}/g} and \boldsymbol{h_{max}=\frac{1}{8}gt_{f}^{2}}.  (Note: \boldsymbol{t_{f}} is twice the time from launch to maximum altitude.)  But a lofted ICBM will reach an altitude comparable to the Earth’s radius \boldsymbol{r_{E}}, so the gravitational acceleration varies significantly during the missile’s flight.  This complicates the analysis considerably.

Let \boldsymbol{G} be the gravitational constant, \boldsymbol{M} be the mass of the Earth, and \boldsymbol{r} be the instantaneous distance of the missile from the center of the planet.  Launched vertically, the speed of the missile at maximum altitude (\boldsymbol{r=r_{max}}) is zero, so by energy conservation,

\large \boldsymbol{\frac{1}{2}v_{0}^{2}-\frac{GM}{r_{E}}=0-\frac{GM}{r_{max}}}

or

\large \boldsymbol{r_{max}=\frac{r_{E}}{1-r_{E}v_{0}^{2}/2GM}=\frac{r_{E}}{1-u_{0}^{2}}}

(1)

where \boldsymbol{u_{0}=v_{0}/v_{esc}} and \boldsymbol{v_{esc}=\sqrt{2GM/r_{E}}} is the escape velocity from the planet’s surface.

Our first task is to find the flight time \boldsymbol{t_{f}} in terms of the burnout speed \boldsymbol{u_{0}}.  One way to proceed is to integrate the equation of motion numerically.  The most straightforward numerical technique is the second-order Taylor approximation (Ref. 3)in which the stepwise equations for position and speed are reminiscent of the kinematic equations for projectile motion:

\large \boldsymbol{r_{n+1}=r_{n}+v_{n}\Delta t+\frac{1}{2}a_{n}\Delta t^{2}}

and

\large \boldsymbol{v_{n+1}=v_{n}+\frac{1}{2}\left ( a_{n}+a_{n+1} \right )\Delta t}

where the acceleration \boldsymbol{a_{n}=-GM/r_{n}^{2}} and \boldsymbol{\Delta t} is the time step.  A simple Excel spreadsheet is adequate for the integration.  Figure 2 shows the computed altitude \boldsymbol{h(t)=r(t)-r_{E}} for a projectile that has reached its maximum altitude at \boldsymbol{t=0}, for the three values \boldsymbol{u_{0}} = 0.5, 0.6 and 0.7.  The intersection of each curve with the x-axis (\boldsymbol{h=0}) marks half the total flight time of the missile.

New Fig 2
Figure 2:  Altitude h(t) for a missile vs. its normalized burnout speed.  The missile is at its peak altitude at time t = 0.

Some instructors may prefer an analytic expression for \boldsymbol{h(t)}.  A power series solution for \boldsymbol{h(t)}, due to Foong (Ref. 4), is outlined in the Appendix and is displayed in Fig. 2 along with the numerical solution.  Figure 3 shows the total flight time \boldsymbol{t_{f}} as a function of \boldsymbol{u_{0}}, calculated numerically as well as by the power series approximation.  Knowing \boldsymbol{t_{f}}, one may find \boldsymbol{u_{0}} directly from Fig. 3.

New Fig 3
Figure 3:  Flight time of ICBM as a function of burnout speed, calculated    numerically (red) and by using the power series expansion (blue) described in the Appendix.  The black line is calculated for the simple case of constant acceleration, to illustrate that it converges with the real solution as u0 → 0.

C.  Finding the Maximum Range of the Missile

When the burnout speed is low (\boldsymbol{u_{0}\ll 1}), the maximum altitude and range of the missile are small, and its trajectory is adequately described as a parabola over a flat Earth.  At the higher speeds needed for intercontinental targeting, the parabolic approximation is no longer adequate, and the trajectory must be treated properly as a segment of an ellipse having one focal point at the center of the planet.  (This follows from Kepler’s Laws.) Changing the launch angle changes the eccentricity \boldsymbol{\varepsilon } of the ellipse  — which affects the range of the missile — but not the length of its semi-major axis \boldsymbol{a}.  The goal of this section is to find the optimum launch angle (the one that yields the maximum range), and to calculate that range as a function of the burnout speed \boldsymbol{u_{0}}.  This can be done easily by appealing to the reflection property of an ellipse (Ref. 5).

In Figure 4, line \boldsymbol{\overline{mn}} is tangent to the ellipse at point P.  The reflection, or focal property states that \boldsymbol{\overline{F_{1}P}} and \boldsymbol{\overline{F_{2}P}} make equal angles with the tangent line.  So at point P, angle φ = φ′ , and at point P′, angle γ = γ′.  If the ellipse were a reflective surface, any light ray emanating from one focus would, after reflection, pass through the other focus.  For example, a parabola is an ellipse with one focus at infinity.  As is well known, light rays originating at infinity (parallel rays) converge at the focal point of a parabolic mirror.

New Fig 4
Figure 4:  The reflection property of an ellipse.  Angles γ=γ’ and φ=φ’.  Any ray emanating from F1 will be reflected into F2, and vice versa.

Next we review the energetics of orbital motion governed by an inverse square central force.  If the eccentricity ε of the elliptical orbit is zero, the trajectory is circular with radius r.  The orbital speed is found by equating the centripetal acceleration to the gravitational force per unit mass, \boldsymbol{v^{2}/r = GM/r^{2}}, so

\large \boldsymbol{\frac{E}{m}=\frac{1}{2}v^{2}-\frac{GM}{r}=-\frac{GM}{2r}}

where is the total mechanical energy of the missile and m is its mass.  When \boldsymbol{\varepsilon \neq 0}, a similar relation holds (Ref. 6): the energy per unit mass is fixed — not by the orbital radius r but by the length of the semi-major axis a:

\large \boldsymbol{\frac{E}{m}=\frac{1}{2}v^{2}-\frac{GM}{r}=-\frac{GM}{2a}}.

An ICBM’s energy and semi-major axis are set by its burnout speed \boldsymbol{v_{0}}, which we have assumed is reached when the missile is close to the Earth’s surface (\boldsymbol{r\simeq r_{E}}):

\large \boldsymbol{\frac{E}{m}=\frac{1}{2}v_{0}^{2}-\frac{GM}{r_{E}}=-\frac{GM}{2a}}.

Using \boldsymbol{v_{esc}=\sqrt{2GM/r_{E}}} and \boldsymbol{u_{0}=v_{0}/v_{esc}}, we obtain

\large \boldsymbol{\frac{2a}{r_{E}}=\frac{1}{1-u_{0}^{2}}}.

(3)

Figure 5a shows the missile trajectory as an ellipse with focus \boldsymbol{F_{1}} at the Earth’s center (the center of force).  The launch point P and the target point Q both lie on the planet’s surface, and the actual missile trajectory is the segment of the ellipse lying outside of the planet (\boldsymbol{r\geq r_{E}}).  The range of the missile is \boldsymbol{R=2r_{E}\theta }, where θ is the angle between \boldsymbol{\overline{F_{1}P} } (or \boldsymbol{\overline{F_{1}Q} }) and the major axis of the ellipse.  The location of the second focus \boldsymbol{F_{2}} depends on \boldsymbol{\varepsilon }, which is governed by the missile launch angle.  For all points P (and, by symmetry,  Q), \boldsymbol{\overline{F_{1}P}+\overline{F_{2}P}=2a}, regardless of the location of \boldsymbol{F_{2}} along the major axis.  Since \boldsymbol{\overline{F_{1}P}=r_{E}}, then \boldsymbol{\overline{F_{2}P}=2a-r_{E}}.  By inspection, the angle θ is maximized when \boldsymbol{\overline{F_{2}P}} is perpendicular to the major axis.  Therefore,

\large \boldsymbol{sin\theta _{max}=\frac{2a-r_{E}}{r_{E}}=\frac{2a}{r_{E}}-1=\frac{u_{0}^{2}}{1-u{_{0}^{2}}}}

and

\large \boldsymbol{R_{max}=2r_{E}\theta _{max}} .

(4)

Besides calculating the maximum range, it is interesting to find the launch angle associated with it.  In Figure 5b, let \boldsymbol{\overline{mn}} once again be tangent to the ellipse at the launch point P.  Using the reflection property discussed above, angle β = β′, so \boldsymbol{\alpha +2\beta =\pi }.  (Note: β is the launch angle measured with respect to the local vertical direction at P.)   Since  \boldsymbol{F_{1}PF_{2}} is a right triangle, angle \boldsymbol{\alpha =\pi /2 -\theta _{max}}.  Therefore,

\large \boldsymbol{\beta =\frac{\pi }{2}-\left ( \frac{\pi }{4}-\frac{\theta _{max}}{2} \right )=\frac{\pi }{4}+\frac{\theta _{max}}{2}} .

In the literature, the launch angle is more commonly measured relative to the local horizontal, i.e., \boldsymbol{\varphi =\pi /2 - \beta }, or

\large \boldsymbol{\varphi =\frac{\pi }{4}-\frac{\theta _{max}}{2}} .

(5)

Note that as \boldsymbol{u_{0}\rightarrow 0}\boldsymbol{\theta _{max}\rightarrow 0} and \boldsymbol{\varphi \rightarrow \pi /4}, or 45°, just as we would expect from basic kinematics.  If \boldsymbol{u_{0}=1/\sqrt{2}} (the condition for low Earth orbit), \boldsymbol{\theta _{max}=\pi /2} and \boldsymbol{\varphi =0}, i.e., the missile should be launched parallel to the Earth’s surface.

New Fig 5
Figure 5: (a) The maximum range of the missile is achieved when P is  farthest from the major axis of the ellipse.  This occurs when F2P is  perpendicular to the major axis.  (b) From the reflection property, the launch angle β = β’, so α = π – 2β.  When β is chosen for maximum range, α = π/2 – θmax, allowing β to be calculated.

According to US Pacific Command, the Hwasong-14 missile launched by North Korea on July 4 was aloft for 37 minutes.  From Figure 2, this indicates that \boldsymbol{u_{0}=0.57}, or \boldsymbol{v_{0}=6.4} km/s.  Eqn. 1 then tells us that the missile reached a maximum altitude \boldsymbol{h_{max}\simeq 3000} km during its test flight.  From Eqn. 4, \boldsymbol{sin\theta _{max}=0.48}, or \boldsymbol{\theta _{max}=0.50} (radians), so the maximum range is \boldsymbol{R_{max}=2r_{E}\theta _{max}\simeq 6400} km.  The launch angle associated with this range is \boldsymbol{\varphi =\pi /4 - 0.25 = 0.54} rad, or 31°.  These results are in good agreement with those of Dr. David Wright, Senior Scientist at the Union of Concerned Scientists, who estimated (Ref. 7) the maximum altitude to be “more than 2800 km” and the maximum range to be “roughly 6700 km.”

D.  Questions for Students

1. The July 4 missile flew for 37 minutes and splashed down about 950 km from its launch site.  Show that the launch angle must have been close to 90°.  (Ans: \boldsymbol{\bar{v_{\theta }}=950/(37\times 60)=0.43} km/s, while \boldsymbol{v_{0}=6.4} km/s.  \boldsymbol{cos\varphi =0.43/6.4 = 0.067}, or \boldsymbol{\varphi \simeq 86}°)

2.   The July 28 missile flew for 47 minutes.  Find the maximum altitude during the test flight, the maximum range  \boldsymbol{R_{max}}, and the optimum launch angle \boldsymbol{\varphi }.  (Ans: \boldsymbol{h_{max}\simeq 4200} km, \boldsymbol{R_{max}\simeq 9200} km, and \boldsymbol{\varphi \simeq 24}°)

3.  Your calculation of \boldsymbol{R_{max}} above assumed that burnout occurred at zero altitude.  Realistically, burnout occurs at an altitude of several hundred kilometers.  The “official” estimate of the range of the July 28 missile was 10,400 km.  Using the launch angle calculated above, estimate the burnout altitude of the missile.  (Ans: about 250 km)

References

  1. Choe Sang-Hun, “U.S. Confirms North Korea Fires an Intercontinental Ballistic Missile,” New York Times, July 4, 2017
  2. Choe Sang-Hun and Eileen Sullivan, “North Korea Launches Ballistic Missile, the Pentagon Says,” New York Times, July 28, 2017
  3. Robert W. Stanley, “Numerical methods in mechanics,” Am. J. Phys. 52(6), 499-507 (1984)
  4. S. K. Foong, “From Moon-falls to motions under inverse square laws,” Eur. J. Phys. 29, 987-1003 (2008)
  5. For a simple proof of the reflection property, see Physics from Planet Earth, J. C. Amato and E. J. Galvez, Taylor & Francis/CRC Press (2015) p. 471
  6. ibid, section 8.2, “Kepler’s Laws”
  7.  David Wright, Union of Concerned Scientists blog, July 3, 2017:  http://allthingsnuclear.org/dwright/north-korea-appears-to-launch-missile-with-6700-km-range

Appendix

The series approximation shown in Figures 2 and 3 is derived by Foong in Ref. 4.  If the missile is at its peak altitude at time \boldsymbol{t=0}, then \boldsymbol{r(t)=r(-t)}; i.e., \boldsymbol{r(t)} is an even function of time, and may be expressed as a power series:

\large \boldsymbol{r(t)=r_{max}+c_{2}t^{2}+c_{4}t^{4}+\cdots } .

The coefficients \boldsymbol{c_{2n} } are chosen to satisfy the equation of motion

\large \boldsymbol{\frac{d^{2}r}{dt^{2}}=-\frac{GM}{r^{2}}} .

Following Foong, we introduce dimensionless variables \boldsymbol{x=r/r_{max}} and \boldsymbol{\tau =t/T}, where \boldsymbol{T=\sqrt{r_{max}^{3}/GM}}.  We then rewrite the power series as \boldsymbol{x=1+{c_{2}}'\tau ^{2}+{c_{4}}'\tau ^{4}+\cdots }, which must satisfy the equation of motion expressed in dimensionless units:

\large \boldsymbol{1+x^{2}\frac{d^{2}x}{d\tau ^{2}}=0}

or

\large \boldsymbol{1+(1+{c_{2}}'\tau ^{2}+{c_{4}}'\tau ^{4}+\cdots )^{2}(2{c_{2}}'+12{c_{4}}'\tau ^{2}+\cdots )=0} .

Expanding the parentheses and multiplying,

\large \boldsymbol{\left ( 1+2{c_{2}}' \right )+ \left ( 4{c_{2}}'^{2} + 12{c_{4}}' \right )\tau ^{2} +\cdots =0} .

For the equality to hold over the entire flight time of the missile, the coefficient of each power of τ must equal zero, so \boldsymbol{{c_{2}}'=-1/2} and \boldsymbol{{c_{4}}'=-1/12}.  See Ref. 4 for the values of higher order coefficients.  Returning to the physical variables r and t, our power series solution is

\large \boldsymbol{r=r_{max}-\frac{1}{2}\left ( \frac{GM}{r_{max}^{2}} \right )t^{2}-\frac{1}{12}\frac{1}{r_{max}^{3}}\left ( \frac{GM}{r_{max}^{2}} \right )^{2}t^{4}+\cdots }

where \boldsymbol{r_{max}} is given as a function of \boldsymbol{u_{0}} in Eqn. 3.  The truncated series containing just the terms shown above is plotted in Fig. 2 for \boldsymbol{u_{0}} = 0.5, 0.6 and 0.7.  The agreement with the numerical solution is quite good for \boldsymbol{u_{0}< 0.5}.

 

 

 

The Trappist-1 Solar System: Seven Earth-sized Planets Potentially Harboring Life

Trappist-1 system vs. Jupiter
Figure 1:  The red dwarf Trappist-1 and its seven planets, compared to the Sun and its solar system (above), and also to Jupiter and its moons (below).  Note that the horizontal dimension displays the period of each orbiting body.  Image credit: Nature

A. Introduction

The recent discovery (Ref. 1) of seven exoplanets orbiting the red dwarf star Trappist-1 is an exciting and potentially key development in the search for extraterrestrial life in the universe.  All of the seven planets are Earth-sized, and at least three lie in the temperate “Goldilocks” zone, neither too near nor too far from their host star for surface oceans to exist.  Given the mass and temperature of these planets, their atmospheres are likely to contain nitrogen, oxygen, water vapor, and carbon dioxide, just like on Earth.  As an added bonus, the small size of the host star will help astronomers determine the chemical composition of each planet’s atmosphere.

The goal of this post is to show students how scientists use elementary physics and trigonometry to derive many of the physical properties of this and other extrasolar systems.  We present straightforward explanations appropriate to an introductory physics course, and offer exercises suitable for a homework assignment or recitation class to help students appreciate the significance of the Trappist-1 discovery.  We hope these materials will be useful to you, and that your students will enjoy and benefit from them.

B.  Background

In February 2017, an international team of astronomers, led by Michaël Gillon of Belgium, announced the discovery (Ref. 1) of seven Earth-sized exoplanets orbiting the red dwarf Trappist-1, so named for the instrument in Chile (TRAPPIST: Transiting Planets and Planetesimals Small Telescope) first employed in their work.  Observations continued for one year using an array of ground-based telescopes located in Chile, Morocco, Hawaii, Spain, and South Africa, plus the orbiting Spitzer Space Telescope.  The exoplanet discoveries were made by the transit photometry technique:  when a planet eclipses its host star, i.e., as it passes between the star and Earth, it blocks a tiny fraction of starlight from reaching Earth.  The resulting light curve (flux vs. time, see Fig. 2) and its period P (time between successive transits), plus estimates of the star’s mass \boldsymbol{M_{\ast }} and radius \boldsymbol{R_{\ast }}, suffice to determine the planet’s radius \boldsymbol{R_{P}} and the radius of its orbit  \boldsymbol{a}.  This information can then be used to estimate the surface temperature \boldsymbol{T} and the potential for life on the planet.

light curves
Figure 2:  Light curves from Trappist-1 planets, taken using NASA’s Spitzer Space Telescope.  Each curve shows a dip in brightness as a planet eclipses the star, and each is labeled with that planet’s orbital period (days).  Image credit: M. Gillon et al/ESO

Astronomers determine stellar masses by studying binary stars.  For each binary system, they measure the distance from Earth by parallax, and measure the apparent brightness of each star to calculate its luminosity \boldsymbol{L_{\ast }}.  Mass is found by studying the orbital motion of the binary, as described in PPE, Section 8.11.  By compiling data on many binaries, an empirical mass-luminosity relation \boldsymbol{L_{\ast }=f(M_{\ast })} can be constructed.  A mass-radius relation can be constructed in a similar way, using data from eclipsing binary stars and long baseline interferometry.  For now, let’s accept the mass and radius of Trappist-1 given in Ref. 1: \boldsymbol{M_{\ast }/M_{\bigodot }=0.0802\pm .0073}\boldsymbol{R_{\ast }/R_{\bigodot }=0.117\pm 0.0036}, where the subscript \boldsymbol{\bigodot } indicates the Sun.  Table 1 is adopted from Ref. 1, and lists many of the properties of the Trappist-1 exosystem.

Data Table
Table 1:  Parameters of the Trappist-1 system, extracted from Table 1 of Ref. 1

1.   A commonly quoted mass-luminosity relation, suitable for some small-mass main sequence stars, is \boldsymbol{L_{\ast }/L_{\bigodot }=(M_{\ast }/M_{\bigodot })^{2.62}}, where the subscript \boldsymbol{\bigodot } indicates the Sun.  The measured luminosity of Trappist-1 is \boldsymbol{L_{\ast }/L_{\bigodot }=5.24\times 10^{-4}}.  Show that this leads to a value of \boldsymbol{M_{\ast }} that agrees with the value quoted in Ref. 1.

C.  Kepler’s Laws and Orbital Radii

Assume that the mass of each planet is much less than that of the host star, \boldsymbol{M_{p}\ll M_{\ast }}, and that each planet is in a near-circular orbit of radius \boldsymbol{a\gg R_{\ast }}, just like the planets of our own solar system.  Kepler’s 3rd law is

\LARGE \boldsymbol{\frac{P^{2}}{a^{3}}=\frac{4\pi ^{2}}{GM_{\ast }}}

(1)

where P is the orbital period.  After calculating \boldsymbol{M_{\ast }} from its luminosity, and measuring P , we can use Eqn. 1 to derive a.

2.  Astronomers have recently (February 2017) discovered seven exoplanets orbiting the star Trappist-1, an ultracool red dwarf about 12 pc from Earth.  This is an exciting discovery because at least three  planets are in the “temperate zone” of the host star, with surface temperatures favorable for life.  The table below lists properties of six bodies that orbit either Trappist-1 (\boldsymbol{M_{\ast }\simeq 0.08M_{\bigodot }}) or one other central body with mass \boldsymbol{M\simeq 0.001M_{\bigodot }}.  Fill in the missing numbers.  (Hint: Use ratios rather than converting to meters, seconds and kilograms.)  Can you identify M and the bodies in orbit about it?

Question 2

D.  Analyzing the Light Curve

Fig. 2 contains a great deal of information that can be used to further characterize the Trappist-1 system.  Even if \boldsymbol{M_{\ast }} and \boldsymbol{R_{\ast }} are not known a priori, the period P and the shape of the transiting planet’s light curve provide enough information to determine \boldsymbol{R_{p}}, the orbital radius a and inclination i, and the density of the host star \boldsymbol{\rho _{\ast }}.  Combining the star’s density with the appropriate mass-radius relation for red dwarfs yields \boldsymbol{M_{\ast }} and \boldsymbol{R_{\ast }}.  This strategy uses Kepler’s 3rd law plus simple kinematics and trigonometry, making it amenable to first-year physics students.  It is described in detail by Seager and Mallén-Ornelas in Ref. 2, and we will follow their treatment closely.

Figure 3
Figure 3:  (a) A sideways view of a planetary transit.  If the inclination i = 90 deg, the orbit is viewed edge-wise by obervers on Earth.  The impact parameter b < 1 for a transit to occur.  (b) A view of the transit from Earth.  The time between the first and last “contact” between the stellar and planetary discs is tt.  The planet is fully within the stellar disc for time tf.  Delta-F is the change in flux reaching Earth due to the eclipsing planet.  Adapted from Fig. 1 of Ref. 2

As usual, assume that the exoplanet’s mass is much less than the star’s, and that it is in a circular orbit of radius \boldsymbol{a\gg R_{\ast }}.  Fig. 3a is a “side view” of the orbit, perpendicular to the line of sight from Earth to the star.  The inclination angle is defined so that if \boldsymbol{i=90\: \textup{deg}}, observers on Earth view the orbit edge-on.  The impact parameter b is defined by \boldsymbol{bR_{\ast }=acosi}.  Clearly, \boldsymbol{b<1} for a transit to occur, and since \boldsymbol{a\gg R_{\ast }}, this implies that \boldsymbol{i\simeq 90\; \textup{deg}}.

From Fig. 3b, the amount of light blocked by the planet is proportional to its area, so the fractional change in flux reaching Earth during a transit is

\LARGE \boldsymbol{\frac{\Delta F}{F}=\frac{R_{p}^{2}}{R_{\ast }^{2}}}.

(2)

In practice, the straight line segments of the light curve shown in Fig. 3b are rounded by limb-darkening (as seen in Fig. 2): less light reaches us from areas near the edge of the stellar disc than from equal areas near the center of the disc.  For simplicity, we will ignore this effect.

3.  Fig. 2 displays the light curves of the seven planets orbiting Trappist-1.  Which planet has the largest radius?  (Ans: Trappist-1b)

Since \boldsymbol{a\gg R_{\ast }}, we may approximate the trajectory of the planet as a straight line as it moves across the face of the star.  For a circular orbit, the planet’s speed is \boldsymbol{v=2\pi a/P}.  By the Pythagorean theorem, the distance traveled by the planet in time \boldsymbol{t_{t}} is \boldsymbol{d_{t}=vt_{t}=2\sqrt{(R_{\ast }+R_{p})^{2}-b^{2}R_{\ast }^{2}}}, so

\LARGE \boldsymbol{t_{t}=\frac{P}{\pi a}\sqrt{(R_{\ast }+R_{p})^{2}-b^{2}R_{\ast }^{2}}}.

(3)

4.  Derive the corresponding equation for \large \boldsymbol{t_{f}}
(Ans: \large \boldsymbol{t_{t}=\frac{P}{\pi a}\sqrt{(R_{\ast }-R_{p})^{2}-b^{2}R_{\ast }^{2}}})

If \boldsymbol{R_{\ast }} and \boldsymbol{M_{\ast }} are not known a priori, then the only information available is the light curve and its period.  Nevertheless, it is still possible to proceed.  From Eqn. 3 and the corresponding expression for  \boldsymbol{t_{f}} (Question 4), it is easy to show that

\LARGE \boldsymbol{(t_{t}^{2}-t_{f}^{2})^{\frac{1}{2}}=\frac{2P}{\pi a}R_{\ast }\left ( \frac{\Delta F}{F} \right )^{\frac{1}{4}}},

where we have used Eqn. 2 to express \boldsymbol{R_{p}/R_{\ast }=\sqrt{\Delta F/F}}.  Rearranging the above expression,

\LARGE\boldsymbol{\frac{a}{R_{\ast }}=\frac{2P}{\pi }\frac{(\Delta F/F)^{\frac{1}{4}}}{(t_{t}^{2}-t_{f}^{2})^{\frac{1}{2}}}}

or

\LARGE \boldsymbol{\frac{a^{3}}{R_{\ast }^{3}}=\frac{8P^{3}}{\pi^{3} }\frac{(\Delta F/F)^{\frac{3}{4}}}{(t_{t}^{2}-t_{f}^{2})^{\frac{3}{2}}}}.

(4)

Combining Eqn. 4 with \boldsymbol{M_{\ast }=4\pi ^{2}a^{3}/GP^{2}} (Kepler’s 3rd law) yields an expression for the stellar “density”:

\LARGE \boldsymbol{\rho _{\ast }=\frac{M_{\ast }}{R_{\ast }^{3}}=\frac{32P}{\pi G}\frac{(\Delta F/F)^{\frac{3}{4}}}{(t_{t}^{2}-t_{f}^{2})^{\frac{3}{2}}}}.

It is convenient to re-express this in terms of the solar “density”

\LARGE \boldsymbol{\frac{\rho _{\ast }}{\rho _{\bigodot }}=\frac{32R_{\bigodot }^{3}}{\pi GM_{\bigodot }}\: \frac{(\Delta F/F)^{\frac{3}{4}}}{(t_{t}^{2}-t_{f}^{2})^{\frac{3}{2}}}\: P}

(5)

where the pre-factor \boldsymbol{32R_{\bigodot }^{3}/\pi GM_{\bigodot }=3.46\times 10^{-3}\: {\textrm{day}}^{2}}.  This is the key result of this section.

For classroom purposes, we can crudely express the mass-radius relation as \boldsymbol{R_{\ast }=M_{\ast }^{x}}. Choosing \boldsymbol{x=0.85}, we can now calculate the mass and radius of Trappist-1 in solar mass units:

\LARGE\boldsymbol{\frac{M_{\ast }}{M_{\bigodot }}=\frac{\rho _{\ast }}{\rho _{\bigodot }}\frac{R_{\ast }^{3}}{R_{\bigodot }^{3}}=\frac{\rho _{\ast }}{\rho _{\bigodot }}\left ( \frac{M_{\ast }}{M_{\bigodot }} \right )^{2.55}}

or

\LARGE \boldsymbol{\frac{M_{\ast }}{M_{\bigodot }}=\left ( \frac{\rho _{\ast }}{\rho _{\bigodot }} \right )^{-0.65}\; and\; \; \; \frac{R_{\ast }}{R_{\bigodot }}=\left ( \frac{\rho_{\ast } }{\rho_{\bigodot } } \right )^{-0.55}},

(6)

where \boldsymbol{\rho _{\bigodot }\equiv M_{\bigodot }/R_{\bigodot }^{3}=5.90\times 10^{3}\: \mathrm{kg/m}^{3}}.  The validity of these results depends sensitively on the mass-radius relation chosen.  Admittedly, the exponent \boldsymbol{x=0.85} was chosen to yield results consistent with those of Ref. 1.  In the literature, \boldsymbol{x=0.80} is more commonly quoted.

5.  According to Ref. 1, the transit duration of Trappist-1d is \boldsymbol{t_{t}=49.1\: \textrm{min}}, the period is \boldsymbol{P=4.05\: \textrm{days}}, and the transit depth is \boldsymbol{\Delta F/F=0.0037}.

a.  Use the values of \boldsymbol{M_{\ast }} and \boldsymbol{R_{\ast }} given in section B to compute \boldsymbol{\rho _{\ast }/\rho _{\bigodot }} and then find \boldsymbol{t_{f}}.  Check your answer against the light curve in Fig. 2.  (Ans: \boldsymbol{t_{f}\approx 0.88t_{t}=43.3\: \textrm{min}})

b.  Find the radius of the planet’s orbit.  (Ans: \boldsymbol{a=2.14\times 10^{-2}\: \textup{\textrm{AU}}})

c.  Find the radius of the planet.  Express \boldsymbol{R_{p}} as a fraction of \boldsymbol{R_{Earth}}.  (Ans: \boldsymbol{R_{p}=0.78\, R_{Earth}})

6.  Researchers at Trump University recently reported the discovery of an eighth planet orbiting Trappist-1.  In a series of tweets, team leaders provided a summary of their observations taken with the privately funded Andrew Jackson Telescope in Palm Beach, FL:  \boldsymbol{P=15.0\pm 0.2\: \textrm{min}}\boldsymbol{\Delta F/F=0.005\pm 0.0001}\boldsymbol{t_{t}=86.4\pm 0.2\: \textrm{min}}, and \boldsymbol{t_{f}=55.4\pm 0.2\: \textrm{min}}.  Their measurements have not yet been confirmed.  If you were a program director at the National Science Foundation, would you approve a grant to the Trump exoplanet program?

E.  Too Hot, Too Cold, or Just Right?

A blackbody is an object that absorbs all of the radiation incident on its surface.  (Since it reflects no light, it appears black.)  But a blackbody is not just an ideal absorber of radiant energy, it is an ideal emitter of energy as well.  The power S emitted by a blackbody per unit surface area is given by the Stefan-Boltzmann Law,

\LARGE \boldsymbol{S=\sigma _{B}T^{4}},

(7)

where T is the surface temperature in kelvin, and \boldsymbol{\sigma _{B}=5.67\times 10^{-8}\: \textrm{W/m}^{2}\textrm{-K}^{4}} is called the Stefan-Boltzmann constant.  The Sun and a cloudless Earth behave approximately as blackbodies.  The power emitted by the Sun, its luminosity, is \boldsymbol{L_{\bigodot }=4\pi \sigma _{B}R_{\bigodot }^{2}T_{\bigodot }^{4}=3.83\times 10^{26}\: \textrm{W}}.  At a distance r from the Sun, this power is spread evenly over a sphere of area \boldsymbol{4\pi r^{2}}.  A spherical planet of radius \boldsymbol{R_{P}} presents a cross-sectional area \boldsymbol{\pi R_{P}^{2}} to the Sun, so if it acts as a blackbody, it intercepts and absorbs energy at a rate \boldsymbol{W_{abs}=L_{\bigodot }R_{P}^{2}/4r^{2}}.  This energy warms the planet’s surface to a temperature , causing it to emit energy at a rate \boldsymbol{W_{emit}=4\pi R_{P}^{2}\cdot \sigma _{B}T^{4}}. When \boldsymbol{W_{emit}=W_{abs}}, the temperature reaches a steady-state value given by

\LARGE \boldsymbol{T=\left ( \frac{L_{\bigodot }}{16\pi \sigma _{B}r^{2}} \right )^{\frac{1}{4}}}.

(8)

7.  The Earth is 1 AU (\boldsymbol{1.5\times 10^{11}\: \textrm{m}}) from the Sun.  What is the steady-state temperature on Earth, according to Eqn. 8?  (Ans: 278 K.  Keep in mind, though, that this ignores the greenhouse effect, which raises the temperature significantly.)

To apply Eqn. 8 to an extrasolar system, simply replace \boldsymbol{L_{\bigodot }} by the luminosity of the host star \boldsymbol{L_{\ast }}.

8.  The period of Trappist-1d is 4.05 days.  Use the star’s mass \boldsymbol{M_{\ast }/M_{\bigodot }=0.08} and luminosity \boldsymbol{L_{\ast }/L_{\bigodot }=5.24\times 10^{-4}} to calculate the steady-state surface temperature on the planet.  Based only on your answer, what is the likelihood of finding liquid water on this planet?  (Ans: T ≈ 300 K.  However, if the planet has an Earth-like atmosphere, the greenhouse effect would raise its temperature too high to harbor liquid water on its surface.  This is mentioned briefly in Ref. 1.)

9.  What is the likelihood of finding intelligent life on Trappist-1i, the planet claimed by researchers at Trump University?

References

1.   Michaël Gillon et al, Nature 542, 456 (2017)

2.  S. Seager and G. Mallén-Ornelas, Ap. J. 585, 1038 (2003)

Work, Energy and the “Satellite Drag Paradox”

iss_nasa02-from-ss-endeavor
Figure 1: A view of the International Space Station taken from Space Shuttle Endeavor.  Every few months, the ISS’s orbit must be boosted to maintain an altitude of about 400 km.  Image credit: NASA

A.  Introduction

Here is a lovely application of energy conservation and the “work-energy theorem” that will surely intrigue your first-year physics students.  A satellite in low Earth orbit experiences a drag force \boldsymbol{\vec{F}_{d}} as it passes through the upper  reaches of the planet’s atmosphere.  (See PPE, Section 6.10, p. 188 ff.)  The drag force is antiparallel to the satellite’s velocity \boldsymbol{\vec{v}}, so it robs the spacecraft of energy and causes it to spiral inward toward the Earth’s surface.  However, assuming a gradual decay and a near-circular trajectory, it is easy to show (see below) that the speed of the spacecraft must increase as the orbital radius decays.  It seems unphysical that by applying a drag (friction) force to a body, its speed will increase!

This “paradox” has been around since the 1950’s, but to my knowledge it is not treated in introductory physics textbooks.  One of the clearest and most complete solutions was given by B. D. Mills in 1959, and we will use his approach in the following paragraphs.  As usual, the presentation will be simplified as much as possible to make it tractable for introductory students.

B.  The Satellite Drag Paradox

Let the satellite be moving on a low Earth near-circular trajectory with a slowly decreasing orbital radius \boldsymbol{r(t)} due to the drag force \boldsymbol{\vec{F}_{d}}.  Mills (Ref. 1) shows that for a weak drag force, the tangential speed of the spacecraft is at all times nearly equal to that of a body moving in a perfectly circular orbit of the same instantaneous radius.  (See also Ref. 5.)  Let \boldsymbol{M} and \boldsymbol{m} denote the mass of the Earth and satellite, respectively.  From Newton’s Second Law,  \boldsymbol{mv^{2}/r=GMm/r}, the tangential speed \boldsymbol{v_{c}} is

\large \boldsymbol{v_{c}=\sqrt{GM/r}}

(1)

and it increases as the orbital radius decreases.  The kinetic energy of a body \large \boldsymbol{m} in circular orbit is found using Eqn. 1:

\large\boldsymbol{K=\frac{1}{2}mv_{c}^{2}=\frac{1}{2}m\frac{GM}{r}=-\frac{1}{2}U},

(2)

where \boldsymbol{U=-GMm/r} is the potential energy, so the total energy of the orbiting body is

\large \boldsymbol{E=K+U=\frac{1}{2}U=-\frac{GMm}{2r}}.

(3)

These results should be familiar to students.

1.  Using Eqns. 2 and 3, verify that as \boldsymbol{r} decreases, \boldsymbol{K} increases even though  \boldsymbol{U} and \boldsymbol{E} decrease.

If the spacecraft’s orbit is initially elliptical (as opposed to circular), it will soon circularize due to the atmospheric drag force.  The drag force increases dramatically at lower altitudes, so it will be greatest when the spacecraft passes through perigee. Imagine that the drag force produces an instantaneous negative impulse at perigee, with change in velocity \boldsymbol{\Delta v} antiparallel to the spacecraft’s velocity.  This reduces the total energy of the spacecraft, shrinking its semi-major axis without changing its perigee position \boldsymbol{r_{p}}.  As a result, the apogee distance \boldsymbol{r_{a}} shrinks, and as \boldsymbol{r_{a}\rightarrow r_{p}}, the orbit becomes more circular.  So even if the initial orbit is significantly elliptical, atmospheric drag will soon circularize it, and the analysis presented below will be applicable.

For a spacecraft experiencing a drag force \boldsymbol{\vec{F}_{d}}, the rate of change of \boldsymbol{E} is equal to \boldsymbol{\vec{F}_{d}\cdot \vec{v}}, where \boldsymbol{\vec{v}} is the total velocity of the spacecraft, including the inspiral motion.  This follows because \boldsymbol{\vec{F}_{d}} is the only external force acting on the Earth-satellite system.  Taking the time derivative of Eqn. 3, and noting that \boldsymbol{\vec{F}_{d}} is antiparallel to the velocity \boldsymbol{\vec{v}}, we find

\large \boldsymbol{\frac{dE}{dt}=\frac{GMm}{2r^{2}}\frac{dr}{dt}=-F_{d}v}

or

\large \boldsymbol{\frac{dr}{dt}=-\frac{2r^{2}}{GMm}F_{d}v}.

Figure 2 illustrates the forces acting on the satellite, and the components of its velocity: \boldsymbol{\vec{v}=v_{c}\hat{i}+\frac{dr}{dt}\hat{j}}.  For satellites passing through the thin upper atmosphere of a planet, \boldsymbol{\left | \frac{dr}{dt} \right |\ll v_{c}}, so we can approximate \boldsymbol{\vec{F}_{d}} to be antiparallel to  \boldsymbol{\vec{v}_{c}}, as shown.  Thus, the above result may be written as

\large \boldsymbol{\frac{dr}{dt}=-\frac{2r^{2}}{GMm}F_{d}v_{c}}.

(4)

presentation3
Figure 2:  The orbiting spacecraft experiences a drag force Fd as well as a gravitational force Fg.  Its velocity has a small radial component dr/dt in addition to the tangential speed vc.  

Using the work-energy theorem \boldsymbol{\Delta K=\vec{F}\cdot \Delta \vec{r}}, where \boldsymbol{\vec{F}} is the total force acting on the body, we obtain

\large \boldsymbol{\frac{dK}{dt}=\vec{F}\cdot \vec{v}=F_{g}\frac{dr}{dt}-F_{d}v_{c}}

(5)

where \boldsymbol{F_{g}=-GMm/r^{2}} is the gravitational force exerted on the spacecraft by the Earth.  Eqn. 5 reveals the source of the apparent paradox: while the drag force decreases \boldsymbol{K}, the satellite gains speed as it falls toward Earth.  To evaluate \boldsymbol{dK/dt}, substitute Eqn. 4 into Eqn. 5:

\large \boldsymbol{\frac{dK}{dt}=-\frac{GMm}{r^{2}}\left ( -\frac{2r^{2}}{GMm}F_{d}v_{c} \right )-F_{d}v_{c}=F_{d}v_{c}}.

(6)

Surprisingly, the drag term appears in Eqn. 6 without a minus sign!  As described by Mills, it is “as if the air drag force, reversed, were pushing the satellite.”

2.   In your own words, resolve the paradox noted by Mills.

C.  Application: The International Space Station (ISS)

The International Space Station (Fig. 1) is maintained in a near-circular orbit of approximate altitude 400 km.  Figure 3 (Ref. 2) is a plot of the ISS’s mean altitude vs. time, and shows how the altitude is periodically boosted (the near-vertical line segments) and how it decays between boosts.  For an on-board demonstration of what ISS astronauts experience during a “reboost,” see Ref. 4.

orbitheightplot
Figure 3:  Mean altitude of the ISS vs. time.  The vertical lines indicate orbital boosts lasting only a few minutes.  Image credit:  Chris Peat, http://www.heavens-above.com

In November 2016, the ISS was boosted to a mean altitude of 406.5 km, and then fell steadily to 404.5 km over the subsequent three month interval (92 d).

3.  What is the orbital speed of the ISS at an altitude of about 400 km?  (Ans: 7.66 km/s)

4.  From November 2016 to February 2017, the altitude decreased by 2 km.  Verify that \boldsymbol{\left | dr/dt \right |\ll v_{c}}, as asserted above.

5.  Calculate the change in the orbital speed \boldsymbol{v_{c}} during that three month period.  (Hint:  \boldsymbol{\Delta v/v_{c}=(?)\Delta r/r}.  Ans: \boldsymbol{\Delta v_{c}=1.13\: \textup{m/s}})

Now let’s calculate how much rocket fuel must be consumed to boost the ISS back to its original altitude of 406.5 km.  To keep things simple, let’s assume that the boost is done in a single short burn, as shown in Figure 4.  A single burn will change a circular orbit to an elliptical one, with perigee \boldsymbol{r_{p}=R_{E}+404.5 \: \textup{km}} at the site of the burn.  If we wish the semi-major axis to be \boldsymbol{a=R_{E}+406.5\: \textup{km}}, then the apogee must be \boldsymbol{r_{a}=R_{E}+408.5\: \textup{km}}, located 180º from the site of the burn.  This results in an eccentricity \boldsymbol{\epsilon \approx 3\times 10^{-4}}, which we will ignore.  (To improve the circularity of the boosted orbit, a two burn strategy is sometimes used.  See Ref. 3.)

single-burn
Figure 4:  A single burn introduces a slight eccentricity, with the altitude at perigee slightly less that the altitude at apogee.

The change in the ISS’s orbital energy due to the burn is

\large \boldsymbol{E_{f}-E_{i}=\frac{1}{2}m(v_{i}+\Delta v)^{2}-\frac{1}{2}mv_{i}^{2}\simeq mv_{i}\Delta v},

where \boldsymbol{v_{i}} is the speed before the burn, and \boldsymbol{v_{i}+\Delta v} is the speed afterwards.  We can solve for \boldsymbol{\Delta v} using Eqn. 3, or, since the altitude change is so small, set \boldsymbol{\Delta U=mg\Delta a}, where \boldsymbol{g} is the gravitational acceleration at the altitude of the ISS, and \boldsymbol{\Delta a=406.5-404.5=2.0\: \textup{km}} is the change in the mean altitude.  From Eqn. 3, \boldsymbol{\Delta E=\frac{1}{2}\Delta U}.

6.  What is the change of velocity needed to boost the altitude of the ISS by 2.0 km?  (Ans: 1.13 km/s)

ISS orbit boosts are carried out using thrusters on the cargo spacecraft that supply the ISS, such as the Russian Progress modules.  (In the past, the US Space Shuttles were used.)  The Progress thrusters use a fuel (UDMH + NTO) with exhaust velocity 2.7 km/s.

7.  The mass of the ISS is \boldsymbol{4.2\times 10^{5}\: \textup{kg}}.  Calculate the fuel mass needed to boost its altitude by 2.0 km.  (Ans: 176 kg)

8.  Calculate the magnitude of the drag force F_{d} acting on the ISS during November 2016 to February 2017.  You might be surprised by your answer!  (Ans: .06 N)

9.  In ref. 4, astronaut Jeff Williams discusses what happens during a typical orbit boost.  For the maneuver shown in the video, \boldsymbol{\Delta v=2.7\: \textup{m/s}}, and the acceleration was \boldsymbol{1.85\: cm/s^{2}}.  What was the duration of the burn?  (Ans: 146 s)

10.  Does your answer to Question 7 include the fuel required to overcome the drag that lowered the orbit in the first place?  Why is this a small correction?  (Ans: the duration of the reboost is only a few minutes, whereas the orbital period is about 90 minutes.  It took many orbits for the drag force to lower the ISS’s altitude.)

11.  There are at least two minor misstatements in the video of Ref. 4.  Can you spot them?  (Ans: At 1:09, Williams states that because of atmospheric drag, “over a period of time we slow down and our altitude over the Earth decreases.  At 3:04 he states “We’re in weightlessness right now, and there’s no acceleration …”  Actually, he and the ISS are undergoing centripetal acceleration of about \boldsymbol{8.66\: m/s^{2}} towards the Earth.)

Acknowledgment

After writing a draft of this post, I asked Dr. Philip Blanco – a valued contributor to this blog – to review it for correctness and to offer suggestions for improvement.  He not only did that, but also supplied a numerical simulation of the trajectory of an ISS-like spacecraft immersed in a thin atmosphere.  His solution, shown below (Figure 5), shows that the spacecraft indeed follows a near-circular orbit right up to the point where a catastrophe occurs, whence it plunges toward Earth.  It would be interesting to calculate when that catastrophe occurs.

My thanks to Dr. Blanco for his assistance.  Of course, I bear the full responsibility for any remaining errors.

orbit-decay-simulation
Figure 5:  An ISS-like spacecraft starts out in a circular orbit of altitude 150 km.  Its trajectory remains near-circular as it loses altitude, until just before it plunges to the Earth.  Courtesy:  Dr. Philip Blanco

References

 1.  Blake D. Mills, Am. J. Phys. 27, 115 (1959)

2.   Chris Peat, http://www.heavens-above.com/IssHeight.aspx

3.  Robert Frost, https://www.quora.com/How-does-the-International-Space-Station-maintain-its-orbit-and-what-propellant-does-it-use

4.  https://www.youtube.com/watch?v=sI8ldDyr3G0, courtesy NASA Johnson Space Center

5.  Leon Blitzer, Am. J. Phys. 39, 882 (1971)

Rosetta’s Final Act: Free-Fall to the Surface of Comet 67P

rosetta-comet67p-768x432
Figure 1:  Artist’s conception of Rosetta approaching comet 67P in September 2016.  Credit:  ESA/Rosetta

The European Space Agency’s historic Rosetta mission is over.  On September 30, 2016, the spacecraft was allowed to fall ballistically (unpowered) to the surface of comet 67P/Churyumov-Gerasimenko, capturing numerous high-resolution photos of the comet’s surface during its 14 hour descent.  The final image, taken from an altitude of 20 meters, is shown in Figure 2.

Comets are presumed to be the primordial building blocks of the solar system, so by observing 67P at close range, Rosetta was investigating the origins of our solar system.

rosetta_s_last_image
Figure 2:  Rosetta‘s last image of the comet surface, taken from an altitude of 20 m, shows an area of about 1 square meter.  Credit:  ESA/Rosetta

During its 12 year lifespan, Rosetta became the first spacecraft to orbit a comet and observe the chemical processes taking place as it passed through perihelion in August 2015, when the rate of gas and dust evolution peaked.  In November 2014, Rosetta was the first spacecraft to deploy a lander, Philae, to the surface of a comet.  Although Philae functioned only briefly, the two allied space probes made numerous discoveries about the chemical makeup and geology of comets.

Comet 67P has a period of 6.45 years and an aphelion distance of 5.7 AU.  At that distance, Rosetta‘s onboard heating system would have been incapable of keeping the spacecraft warm enough for its instruments to hibernate safely.  Consequently, flight directors opted to end the mission by crash-landing the spacecraft onto the comet, gathering as much data as possible while its instruments were still functioning.  Today, 67P is 4.3 AU from the Sun.  As it heads toward aphelion in November, 2018, it carries the remains of Rosetta and Philae, the two pioneering spacecraft that doggedly chased it for over a decade.

In this blog post, we present a set of student-ready exercises based on Rosetta‘s final maneuver and descent to the surface of 67P.  These exercises were inspired by Dr. Philip Blanco’s recent letter (Ref. 1) to the editor of The Physics Teacher.  Both Philip and I are deeply indebted to Pablo Muñoz, flight dynamics engineer at the ESA’s European Space Operations Centre in Darmstadt, Germany.  Pablo supplied indispensible data on the position and velocity of the spacecraft before and after its final thruster burn, plus the mass of Rosetta and the amount of fuel consumed during that final burn.  Without his input, none of these exercises could have been written.

To make the exercises suitable for introductory students, I have modified the details slightly, as indicated below.   To fully appreciate the complexity of Rosetta‘s dance with its cometary partner, follow the links in Ref. 4 and 5 to the beautiful ESA videos of the mission.  Have fun!

A.  Problems

In these exercises, \boldsymbol{ GM=666.15\mathrm{ m^{3}/s^{2}}}, where
\boldsymbol{ M} is the mass of the comet (\boldsymbol{ M\simeq 1.0\times 10^{13}} \textup{ kg}), and model the comet as a homogeneous sphere of radius 2.00 km.

Prior to its plunge to the comet surface, Rosetta was placed in an unpowered elliptical orbit with apocenter (farthest distance to the comet’s CM) \boldsymbol{r_{a}=23.44}\textup{ km} and pericenter \boldsymbol{r_{p}=13.58}\textup{ km}.  This is shown nicely in Ref. 4.

1.  Calculate the specific energy \boldsymbol{ E/m} (where \boldsymbol{m} is Rosetta‘s mass) and the orbital period of the spacecraft.  (Ans: \boldsymbol{E/m=-1.800\times 10^{-2}}\textup{ J/kg}\boldsymbol{ T=6.13\times 10^{5}}\textup{ s}\simeq 1\textup{ week})

2.  What was the speed \boldsymbol{ v_{0}} of the spacecraft at apocenter? Compare this to the speed of a spacecraft in low Earth orbit.  Why the big difference?  (Ans: \boldsymbol{ v_{0}=0.1444}\textup{ m/s}\boldsymbol{M_{67P}/M_{Earth}\sim 10^{-12}})

The final maneuver that sent Rosetta diving to the comet surface occurred at a true anomaly of 208º, i.e., 28º beyond apocenter (0º corresponds to pericenter).  (See Ref. 4.)  To streamline the following questions, assume instead that the maneuver took place exactly at apocenter.  This alters the results modestly.

3.   Immediately after the burn, Rosetta‘s velocity was \boldsymbol{ \vec{v}{_{0}}'=-0.0211\hat{i}-0.3255\hat{j}}\textup{ m/s}, where \boldsymbol{ \hat{i}} and \boldsymbol{ \hat{j}} refer to the radial and azimuthal directions (Ref. 2) at apocenter.  See Figure 3.  Show that the burn inserted the spacecraft into an unbounded hyperbolic trajectory.  (Ans: \boldsymbol{ {E}'/m=2.478\times 10^{-2}}\textup{ J/kg}>0)

figure-3v2
Figure 3:  Spacecraft velocity immediately after its final burn at apocenter

4.  Note that \boldsymbol{ {\vec{v}}'_{0}} is not perfectly radial, i.e., directed toward the comet’s CM.  Convince someone that the spacecraft will still strike the comet.  This question and the solution below were suggested by Philip Blanco.          (Ans:  Compare the angle subtended by the comet at apocenter to the angle between the velocity vector \boldsymbol{{\vec{v}}'_{0}} and the radial direction.)

5.  With what speed \boldsymbol{ v} will the spacecraft strike the surface?  Compare this to the surface escape speed.  (Ans: \boldsymbol{ v=0.846}\textup{ m/s}\boldsymbol{ v_{esc}=0.816}\textup{ m/s})

6.  Find the spacecraft’s velocity \boldsymbol{ \vec{v}=v_{x}\hat{i}+\vec{v}_{y}\hat{j}} just before it strikes the surface.  Hint: first conserve angular momentum.  (Ans: \boldsymbol{ \vec{v}=-0.809\hat{i}-0.247\hat{j}}\textup{ m/s})

7.  Calculate the magnitude of the change in velocity caused by the thruster burn.  The mass of the spacecraft at this time was 1422.5 kg, and the burn consumed 0.2 kg of rocket fuel.  What was the effective exhaust velocity of the fuel (a mixture of monomethyl hydrazine (MMH) and nitrogen tetroxide (\large \mathrm{N_{2}O_{4}} or NOX))?  (Ans: \boldsymbol{ \Delta v_{0}=.365}\textup{ m/s}\boldsymbol{ v_{ex}=2.6}\textup{ km/s}.  This is substantially less than the “published” value of 3.35 km/s (Ref. 3) for MMH/NOX because a fuel leak forced flight engineers to operate the thrusters at lower than optimal pressure.)

8.  In September 2016, comet 67P has a rotation period of 12.055 hours.  This is less than its period in 2014 (12.4 hours), when Philae was released, due to outgassing from the comet as it passed through perihelion.   Assume that the comet’s rotation axis was perpendicular to the plane of the spacecraft’s motion (Figure 3).  Rosetta‘s velocity \boldsymbol{ \vec{v}_{0}{}'} was selected so that its impact velocity relative to the rotating surface was nearly vertical, i.e., perpendicular to the impact surface.  (Hence, an observer standing on the surface at the impact point would have seen the spacecraft descending vertically.)  Show that this is roughly consistent with your calculated value of \boldsymbol{v_{x}} in Question 6.  (Note: because the comet has such an irregular shape, the local vertical direction was significantly different from the radial direction.) (Ans:\boldsymbol{ v_{rot}=-\omega R=-29.0\textup{ cm/s} \approx v_{x}=-24.7\textup{ cm/s}})

The following question is suitable for users of Physics from Planet Earth who have studied Section 12.9: Application: the “Slingshot Effect” Revisited, or who have seen a comparable treatment elsewhere.

9.  Suppose Rosetta‘s final maneuver had been incorrectly executed, and the spacecraft just missed striking the comet, instead grazing its surface as it flew by.

a.  From Question 3, we know that Rosetta‘s speed after the burn was 0.3262 m/s at apocenter.  In the comet’s reference frame, what will be its speed after it has flown past the comet and is very far from it?  (Ans: \boldsymbol{ v_{\infty }=0.2226\textup{ m/s}})

b.  By what angle Θ will the spacecraft be deflected after its close fly-by of the comet?  (Ans:  \boldsymbol{ \tan (\Theta /2)=1.769}\boldsymbol{ \Theta =121^{\circ}}.  See Eqn. 12.12 of PPE.)

 

B.  Sphere of Influence of Comet 67P

Here is an additional “space-themed” topic that is quite suitable for discussion in an introductory mechanics course.  As far as I know, it is not discussed in any currently available textbook (including PPE).  I thank Philip Blanco for suggesting it.

So far, we have ignored the influence of the Sun on the motion of Rosetta about comet 67P.  We might justify this by invoking Einstein’s Principle of Equivalence (see Section 7.5 of PPE) which states that in a free-fall reference frame, the local effects of gravity disappear.  Comet 67P’s heliocentric orbit is free-fall motion about the Sun, so in the comet’s reference frame, we can ignore the Sun’s gravitational force on the comet (which is obvious).  However, there is a small but important complication when we consider the motion of the spacecraft in the same reference frame.  Since the Sun’s gravitational force varies inversely with the square of the distance, the spacecraft’s acceleration toward the Sun varies as the distance between it and the Sun changes along its orbit.  Let \boldsymbol{ R} be the distance from the Sun to the comet, and \boldsymbol{ r} be the distance from the comet’s CM to the spacecraft.  See Figure 4.  When the spacecraft is at position a(b) in the figure, its distance from the Sun is a maximum (minimum), and its acceleration due to the Sun is

\LARGE \boldsymbol{a=\frac{GM_{Sun}}{(R\pm r)^{^{2}}}\simeq \frac{GM_{Sun}}{R^{2}}\left ( 1\mp \frac{2r}{R} \right )}.

(1)

figure-4
Figure 4:  R is the distance from the Sun to the comet CM, and r is the distance from the comet CM to the spacecraft.  At a and b, the craft is farthest and nearest to the Sun.

The acceleration of the comet toward the Sun is \boldsymbol{ {a}'=GM_{Sun}/R^{2}}, so the maximum difference between the spacecraft’s acceleration and that of the comet is \boldsymbol{ \Delta a_{max}=2GM_{Sun}\, r/R^{3}}.  How close is the spacecraft’s motion to an unperturbed Keplerian orbit about the comet?  The answer depends on the ratio between \boldsymbol{ \Delta a_{max}} and the centripetal acceleration of the spacecraft toward the comet, \boldsymbol{ a_{c}=GM/r^{2}}, where \boldsymbol{ M} is the comet mass.  The smaller this ratio is, the smaller the influence of the Sun, and the better the motion can be described as an unperturbed Keplerian orbit about the comet.

Now let’s ask the opposite question.  Imagine the spacecraft is far enough from the comet that its motion is better described as a heliocentric orbit perturbed slightly by the comet.  How close is its motion to a pure heliocentric orbit?  The answer depends on the ratio between the acceleration \boldsymbol{ a_{c}} toward the comet and the acceleration toward the Sun \boldsymbol{ {a}'}:

\LARGE\boldsymbol{ \frac{a_{c}}{{a}'}=\frac{M}{M_{Sun}}\frac{R^{2}}{r^{2}}}.

(2)

The smaller this ratio is, the smaller the influence of the comet, and the better the motion can be described as a Keplerian orbit about the Sun.  The “sphere of influence” of the comet is the region within which the motion of the spacecraft is only weakly perturbed by the Sun.  By convention, its radius \boldsymbol{ r_{soi}} is found by equating the ratios in Eqns. 1 and 2, and solving for \boldsymbol{ r\textup{ }\: (=r_{soi})}:

\LARGE\boldsymbol{\frac{2M_{Sun}}{M}\frac{r^{3}}{R^{3}}=\frac{M}{M_{Sun}}\frac{R^{2}}{r^{2}}}

or

\LARGE\boldsymbol{ \frac{r^{5}}{R^{5}}=\frac{1}{2}\left ( \frac{M}{M_{Sun}} \right )^{2}}.

Ignoring the factor \boldsymbol{ (1/2)^{1/5}=0.87\approx 1}, we obtain the radius of the sphere of influence:

\LARGE\boldsymbol{ r_{soi}=R\left ( \frac{M}{M_{Sun}} \right )^{\frac{2}{5}}}.

(3)

10.  a.  Philae was released from Rosetta when comet 67P was about 3 AU from the Sun, from an initial distance 22.5 km from the comet’s CM.  Show that the spacecraft was well within the comet’s sphere of influence at that time.  (Ans:  \boldsymbol{ r_{soi}=54\textup{ km}})

      b.  At perihelion (August 13, 2015), 67P was 1.24 AU from the Sun.  To escape possible damage from the comet’s efflux, Rosetta was removed to a distance of 330 km.  Was it in an unpowered orbit about 67P at that time?  (Ans: No, \boldsymbol{ r_{soi}=22\textup{ km}}.  See Ref. 5, fast-forward to August 2015.)

C. References

1a.  “Modeling Rosetta’s final descent,” Philip Blanco, The Physics Teacher 54, 516 (2016)

2a.  If students are familiar with polar coordinates, substitute \boldsymbol{ \hat{r}} and \boldsymbol{ \hat{\theta}} for \boldsymbol{ \hat{i}} and \boldsymbol{ \hat{j}}.

3.   https://en.wikipedia.org/wiki/Liquid_rocket_propellant

4.  https://www.youtube.com/watch?v=tD4c8evE4YA

5.  https://www.youtube.com/watch?v=dc-ICdwX5I0

Gravitational Radiation 3: GW151226, Encore from LIGO

gw151226-image
Figure 1: A filtered, reconstructed image of the gravitational wave detected by LIGO on December 26, 2015 (UTC).  The signal was generated from the merger of two black holes, of approximately 14.2 and 7.7 solar masses, roughly 1.4 billion light years from  Earth.  Figure taken from Abbott et al, Ref. 3.

LIGO (Laser Interferometer Gravitational-Wave Observatory) consists of two gargantuan Michelson-like interferometers, with arms 4 km long, located 3000 km apart in Livingston LA and Hanford WA.  Its first observing session began on September 12, 2015.  Incredibly, just two days later, it captured the fleeting signal (designated GW150914) produced by the merger of two inspiraling black holes more than 1 billion light years from Earth.  (Ref. 1)  That groundbreaking discovery, announced on February 11, 2016, ended a 60 year race to directly detect gravitational waves, and came 100 years after the phenomenon was first predicted by Albert Einstein.

On Christmas Day (in the USA), LIGO scientists detected a second event, GW151226 (Ref. 3), proving that the earlier discovery was not a fluke, and that black hole mergers are frequent events.  (Ref. 4,5)  Since January, LIGO’s sensitivity has been improved, and it will begin a new observing run in Autumn, 2016.  This time, the two US detectors will be joined by VIRGO, a similar instrument with arms 3 km long, located near Pisa, Italy.  The addition of VIRGO will greatly assist in localizing – by triangulation – the sources of gravitational radiation.  In the near future, when LIGO and VIRGO reach their full design sensitivities, they may be capable of detecting one or more black hole collisions daily!  Indeed, the 21st Century seems destined to become the era of gravitational wave astronomy.

In this brief post, we will use the mathematical treatment outlined in Ref. 2 to analyze the published data given in Ref. 3 for GW151226.  Our educational goal is to add to the collection of homework exercises related to gravitational waves that are suitable for first year physics students.  In this way, we hope to inspire instructors to include this exciting topic in their mechanics syllabi.

 

Analyzing GW151226

In Ref. 2, the chirp mass of the binary black hole system was defined as \large \mathfrak{M}=\frac{(m_{1}m_{2})^{3/5}}{(m_{1}+m_{2})^{1/5}}, where \large m_{1} and \large m_{2} are the masses of the black holes.  This quantity is important because it can be found directly from the time dependence of the gravitational wave’s (GW) frequency (Eqn. 10 of Ref. 2):

\large \mathfrak{M}=\frac{c^{3}}{G}\left ( \frac{5}{96}\: \pi ^{-8/3}f^{-11/3}\: \frac{df}{dt} \right ) .

(1)

Following Ref. 2, define \large A=\frac{c^{5}}{G^{5/3}}\: \frac{5}{96}\: \pi ^{-8/3}=5.45\times 10^{56} (SI units), and integrate Eqn. 1 over the time interval \large \Delta t=t_{2}-t_{1} to obtain

\large \mathfrak{M}^{5/3}\Delta t=A\int_{t_{1}}^{t^{_{2}}}f^{-11/3}df=-\frac{3}{8}Af^{-8/3}\mid _{t_{1}}^{t_{2}} .

(2)

1.  Referring to Fig. 1 above, at 0.80 s before the two black holes coalesce, the frequency of the GW was 39.2 Hz.  Later, 0.40 s before the merger, the frequency was 50.3 Hz. (Ref. 6)  Calculate the chirp mass and express your answer in terms of \large M_{Sun}.  (Ans: \large 9.7\: M_{Sun}.)

2.  Show that the total mass \large M=m_{1}+m_{2} of the binary system before coalescence was at least \large 22\: M_{Sun}.  (Hint: Let \large m_{2}=\alpha m_{1}.  Express \large M\large \mathfrak{M}, and \large M/\mathfrak{M} in terms of \large m_{1} and \large \alpha.  Show that the ratio is a minimum when \large \alpha=1.)

3.  The GW frequency at the moment of coalescence of the two black holes was 420 Hz.  (See Fig. 1.)  Recall from Ref. 2 that the orbital frequency of the binary system is half the frequency of the GW.   Use the total mass given in Ref. 3, \large M=22\: M_{Sun}, to find the distance between the bodies just before they merged.  (Ans: 120 km)

4.  Assume that one black hole was twice as massive as the other, i.e., \large m_{2}=2m_{1}, which is about what the LIGO team concluded.  If the orbits were circular, what were their radii just before coalescence?  What were their speeds? Ignore relativity.   (Ans: 40 and 80 km; \large v_{1}=0.35\, c)

References:

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

2.  this blog:  “Gravitational Radiation 2: The Chirp Heard Round the World”

3.  B. P. Abbott et al, “Observation of Gravitational Waves from a 22-Solar-Mass Binary Black Hole Coalescence,”, Phys. Rev. Lett. 116, 241103 (2016)

4.  A. Cho, “LIGO Detects Another Black Hole Crash,” Science 352, 1374, 17 June 2016

5.  D. Castelvecchi, “LIGO Sees a Second Black Hole Crash,” Nature 534, 448, 23 June 2016

6.  We thank Jonah Kanner at the LIGO Open Science Center for providing these numbers, which were calculated using the tutorial on the Center’s webpage:  https://losc.ligo.org.  The numbers obtained from the tutorial differ slightly from those reported in Ref. 3 because the analysis is not identical to the one used for the paper.  (In Ref. 3, \large \mathfrak{M}=8.9\pm 0.3\: M_{Sun}.)  The LIGO Open Science Center is a service of LIGO Laboratory and the LIGO Scientific Collaboration.    LIGO is funded by the U.S. National Science Foundation.

Proxima Centauri b: An Earth-like Neighbor in our Galactic backyard

proxima-centauri-b (II)
Figure 1:  Artist’s conception of sunset on Proxima Centauri b.  The binary star Alpha Centauri is also visible.  Inset:  Proxima Centauri is the Sun’s nearest neighbor, 1.294 pc (4.224 ly) away.  Image credit:  Pale Red Dot/ESO

In August (2016), astronomers at the European Southern Observatory in Chile announced the discovery of an Earth-sized exoplanet orbiting the red dwarf Proxima Centauri, the nearest star to our Sun (Ref. 1, 2).  The planet, called Proxima Centauri b, was detected by the radial velocity method:  the orbiting planet causes the star to execute a much smaller orbit of its own (as required by momentum conservation), and light from the star is Doppler-shifted due to its radial motion relative to Earth.  The discovery of Proxima Centauri b is exciting because its orbit lies within the “Goldilocks,” or temperate, zone of the star, where the planet’s surface temperature would allow water to exist in liquid form – a likely prerequisite for life.  Red dwarfs are the most common stars in our galaxy, so if this particular red dwarf harbors a life-supporting planet, it is likely that the Milky Way, which contains over 100 billion stars, is brimming  with life.

In this post, we will use the data reported in Reference 1, plus the strategy presented in Section 8.12 of PPE, to calculate the mass and orbital radius of the exoplanet Proxima Centauri b.  Our goal is to help students share the excitement of this new discovery through their understanding of introductory mechanics.

 

A.  Calculating the exoplanet’s mass and orbital radius

In the following, we will assume \large m\ll M, where \large m and \large M are the masses of the planet and star, respectively.  We will also assume that their orbits are circular, which is consistent with observations, and begin by treating the simplest case where the orbits are viewed edge-on from Earth.  (Later, we’ll relax this restriction.)  Radial velocity measurements of Proxima Centauri collected over the past 16 years are compiled in Figure 2.  Careful analysis indicates that the period \large T of the star’s orbit is 11.186 ± .002 d, and its orbital speed \large v is 1.38 ± 0.02 m/s.

figure-2
Figure 2:  Radial velocity measurements of Proxima Centauri, compiled over 16 years of observations.  Note the sinusoidal shape of the best-fit curve.

Kepler’s 3rd law states that

\large \frac{T^{2}}{r^{3}}=\frac{4\pi ^{2}}{GM},

where \large r is the orbital radius of the planet.  Observations of the star’s luminosity and color indicate that its mass \large M=0.120\pm .015\, M_{Sun}.

 

 

 

1.  Find the planetary orbit radius \large r (in AU).  Hint: use  \large T=11.186/365.25=.03063\, \mathrm{yr} and \large r^{3}\propto MT^{2} to compare this orbit to Earth’s orbit about the Sun.  (Ans: \large 0.048 \: \mathrm{AU} = 7.24\times 10^{9}\: \mathrm{m}.

2.  The orbital radius of the star \large R can be found from its period and radial velocity \large v\large vT=2\pi R.  Using  \ 1\: \mathrm{d}=8.64\times 10^{4}\: \boldsymbol{\mathrm{s}}, find \large R.  (Ans: \large 2.12\times 10^{5}\: \mathrm{m})

3.  The planet and star co-orbit their center of mass.  Find the mass of the exoplanet.  (Ans:   \large m=MR/r=2.93\times 10^{-5}\, M=6.99\times 10^{24}\, \mathrm{kg}=1.18\, m_{Earth})

If the orbit is not viewed edge-on, but is inclined to the line of sight (LOS), the measured radial speed of the star is less than the actual orbital speed \large v used in the above analysis:  \large v_{meas}=vsin(i), where i, the angle of inclination, is unmeasurable.  See Figure 3.  In this case, \large 2\pi R=Tv=Tv_{meas}/sin(i)).

figure-3
Figure 3:  In the general case, the orbits of the star and planet are not viewed edge-on (i = 90 deg) from Earth.  The measured speed is less than the actual speed by a factor of sin(i).

4.  Prove that, in the general case, the radial velocity measurements are only sufficient to determine \large msin(i).  Is the exoplanet mass you found in Question 3 the maximum, or minimum, mass of the planet?

5.  Because red dwarfs are much cooler than the Sun, their temperate zones are much closer to them than the Earth is to the Sun.  In addition, red dwarfs are less massive than the Sun.  Explain why both of these factors favor the detection of small habitable exoplanets.  (See Ref. 2 for a brief discussion.)

6.  Proxima Centauri subtends an angle of 1.02 ± .08 mas (milliarcseconds) when viewed from Earth.  (See PPE, Section 1.7)  Find its radius \large R_{s}.  (Ans: \large 1.0\times 10^{8}\: \mathrm{m}=0.14\, R_{\mathrm{Sun}})

If Proxima centauri b transits, or passes in front of, its star, as viewed from Earth, it will intercept some of the starlight and allow its radius to be determined.  Additionally, if the planet has an atmosphere, it will absorb certain wavelengths of light preferentially, enabling astronomers to detect the presence of atmospheric gases such as water, oxygen, carbon dioxide, and methane.  Methane (CH4), in particular, is often associated with the presence of life.

The probability of transit is found from the geometric argument outlined in Figure 4.  In the figure, the horizontal dotted line is the LOS from Earth, and the vertical line represents the plane perpendicular to the LOS passing through the star’s center.  For an orbit to be transiting, its normal (shown as a red line for each of the two orbits depicted) must lie within \large \pm \theta _{max}=\pm R_{s}/r of the vertical plane, where \large R_{s} is the star’s radius.  But the normal to the orbit does not need to lie in the plane of the figure: any transiting orbit rotated by any angle (0 – 2π) about the LOS would remain a transiting orbit.  The total solid angle available for transiting orbits is therefore \large 2\pi \cdot 2\theta _{max}=4\pi R_{s}/r, whereas the total solid angle associated with all possible orientations is \large 4\pi.  Therefore, the probability of a transiting orbit is \large R_{s}/r.  (Check: if \large r=R_{s}, then the probability equals 1.)  For a fuller discussion, see  http://kepler.nasa.gov/Science/about/characteristicsOfTransits/ .

figure-4
Figure 4:  In a transiting orbit, the planet passes in front of the star, as seen from Earth.  For this to happen, the radius of the planet’s orbit r must be less than the stellar radius Rs.

7.  How probable is it that Proxima Centauri b transits its star?  (Ans: 1.5%)

B.  How Warm is Proxima Centauri b?

The following questions are only suitable for students who have studied blackbody radiation in a course of thermal physics.

8.  The luminosity of a star (its total radiated power) is proportional to its surface area and the fourth power of its surface temperature: \large L_{star}\propto R_{star}^{2}T_{star}^{4}.  Proxima Centauri’s surface temperature is about 3050 K, whereas the Sun’s is 5780 K.  Show that \large L_{Proxima}\simeq 0.0015\, L_{Sun}.  This indicates that the temperate zone lies much closer than 1 AU to the star.

9.  The fraction of the star’s radiation that is intercepted by a planet is equal to the solid angle subtended by the planet divided by 4π:  \large \pi R_{p}^{2}/4\pi r^{2}.  The absorbed radiation warms the planet, which reaches a steady state temperature \large T_{p} when the blackbody radiation emitted by the planet equals the radiation absorbed from the star:
\large 4\pi R_{p}^{2}\sigma _{B}T_{p}^{4}=L_{Proxima}\pi R_{p}^{2}/4\pi r^{2}, where  \large \sigma _{B} is called the Stefan-Boltzmann constant.  (Its value is not needed to answer this question.)  Use this information to compare the surface temperature of Proxima Centauri b to that of the Earth, assuming both planets are perfect black bodies (or have the same albedo).  (Ans:   \large T_{p}/T_{Earth}=0.90)

References

1.    G.  Anglada-Escudé et al, Nature 536, 437 (25 August, 2016)

2.   A. P. Hatzes, Nature 536, 408 (25 August, 2016)