Galilean satellites: Evolutionary paths in deep resonance

ICARUS 70, 3 3 4 - 3 4 7 (1987)

Galilean Satellites" Evolutionary Paths in Deep Resonance RICHARD

GREENBERG

Lunar and Planetary Laboratory, University of Arizona, Tucson, Arizona 85721 R e c e i v e d S e p t e m b e r 8, 1986; r e v i s e d J a n u a r y 12, 1987

The Laplace resonance among the inner three Galilean satellites (mean motions na 3n2 + 2n3 = 0) has stable configurations in "deep resonance," i.e., where mean motions taken by pairs are in ratios very close to 2 : 1. The present satellite configuration, with the resonance variable ¢~ --- ,~1 - 3A2 + 2,~3 stable at 180 °, is unstable near this exact commensurability. But there is a continuous path of stable conditions branching from = 180 ° to higher and lower values of ~ and toward very deep resonance, according to a theory extended to third order in orbital eccentricity. This path provides a track for tidal evolution of the system. Thus, scenarios involving evolution (probably episodic) from deep resonance are viable, and eliminate the requirement by the alternative equilibrium hypothesis for rapid tidal dissipation in Jupiter. Evolution out from deep resonance is consistent with the free eccentricity of Ganymede, the free libration of ~ , and observational constraints on lo's secular acceleration. Also, the relatively large forced eccentricities in deep resonance may have controlled geophysical processes in the satellites by much greater tidal heating and global stress than at present. © 1987AcademicPress, Inc.

pairs. The mean motions of Io(1) and Europa(2) obey the relation nl - 2n2 = v where v is a quantity much smaller than either nj or n:. Current values are /l I --203.5°/day, rt 2 = 101.4°/day and v -~ 0.74°/ day. If v were precisely zero, Io would make two orbits for every one of Europa's so that the longitude of conjunction of the two satellites would be fixed. With nl - 2n2 actually not quite equal to zero, the longitude of conjunction moves at the slow rate v. The pair Europa(2) and Ganymede(3) are also near a 2:1 commensurability, with n2 - 2n3 equal to the same small value v. Thus conjunction of E u r o p a and G a n y m e d e also migrates at the same slow rate. These kinematic relationships yield periodic geometrical configurations which enhance the satellite's mutual gravitational perturbations with two major effects. First, the commensurability is maintained because the angle ~b between the two conjunction longitudes is stabilized at a constant value 180°. This result ensures that nl - 2n2 and n2 - 2n3 remain equal to the same small

I. I N T R O D U C T I O N

The Galilean satellites are locked into the Laplace orbital resonance, which is responsible for forcing Io and E u r o p a ' s significant orbital eccentricities (Greenberg 1976). The eccentricities in turn result in tidal energy dissipation within the satellites that is responsible for the dramatic thermal activity on Io (Peale et al. 1979), including both volcanism observed by Voyager (Smith et al. 1979) and the infrared heat flux measured from both Earth and spacecraft (McEwen et al. 1985). This forced tidal heating has probably played a dominant role in governing the geological and geophysical evolution o f E u r o p a as well (Cassen et al. 1982). The cause-and-effect loop is completed as the tidal dissipation, modulated by the satellite's geophysical properties, controls the long-term evolution of the orbital resonance. The resonance results from the nearly 2 : 1 commensurability of orbital periods (or equivalently of mean motions) taken by 334 0019-1035/87 $3.00 Copyright © 1987 by Academic Press, Inc. All fights of reproduction in any form reserved.

RESONANCE OF GALILEAN SATELLITES v, or equivalently that the Laplace relation nl - 3n2 + 2n3 = 0 is maintained. Second, orbital eccentricities are enhanced, promoting tidal heating with its important geological and geophysical consequences. The forced eccentricities are roughly inversely proportional to v. Tidal energy dissipation also affects the long-term orbital evolution of the satellites. Any tidal heat dissipated within a satellite must come from the satellite's orbital energy. Long-term evolution may also be affected by tides raised on Jupiter by the satellites. This process may add energy to the satellites' orbits at the expense of Jupiter's rotation. Because tidal effects are strongly dependent on the distance between planet and satellite, tides raised on or by Io, rather than other satellites, dominate the evolution. Yoder (1979) pointed out that tides raised on Io tend to drive the system out of resonance, by which we mean that the mean motion ratios become farther from the exact 2:1 ratio, the value of v increases, and the forced eccentricities become smaller. This conclusion followed from consideration of equations for variation of orbital elements. But it also follows from energy and angular momentum considerations: It is the only way evolution can go, given interactions that maintain the Laplace relation, a net orbital energy loss, and conservation of total orbital angular momentum. Yoder showed that in contrast tides raised on Jupiter, which add both orbital energy and angular momentum to the satellite system, tend to drive the system into " d e e p e r " resonance, that is toward smaller values of v and greater forced eccentricities. Yoder proposed that the system began out of resonance (v -> 0), but evolved into resonance due to tides on Jupiter. As v decreased, Io's forced eccentricity increased, so that tides on Io became increasingly important. Eventually, the decrease in v was halted by the contrary tendency of the tides on Io. Thus Yoder interpreted the present state of the system and present value of v as

335

representing an equilibrium between the competing effects. Yoder's proposed scenario thus has two component hypotheses: (a) resonance capture and (b) present equilibrium. The equilibrium is stable because an increase (or decrease) in v will decrease (or increase) Io's eccentricity and decrease (or increase) tidal heating in Io. This will drive v back toward the equilibrium value. The difficulty with Yoder's elegant theory is that the present (presumed equilibrium) value of v implies a specific ratio between the energy dissipation rates in Jupiter and in Io, and the high measured and inferred heat loss from Io implies a correspondingly high, and difficult to explain, tidal dissipation in Jupiter (Greenberg 1982). The original estimate of tidal dissipation in Io by Peale e t al. (1979) was 2 × 1013 W. Assuming Yoder's equilibrium, this value corresponds to Qj ~ 105, where Qj is the dissipation parameter (inversely proportional to dissipation rate) for Jupiter. In contrast, estimates at that time of physically plausible dissipation processes in Jupiter generally yielded Qj values >107 (see discussion in Greenberg 1982). This problem with the equilibrium model is exacerbated by the most recent evaluation and extrapolation of infrared mapping of Io (McEwen et al. 1985), which gives a heat loss rate of 8 × 1013 W --- 25% and a required value for equilibrium of Qj < 3 × 104. The only physical model that gives such rapid dissipation involves condensation of helium clouds deep in Jupiter's interior (Stevenson 1983). It actually allows Qj to be as small as -102. However, that model is highly speculative because the result depends sensitively on such unknowns as the amount of helium at the appropriate depth, the density of droplet nucleation sites, and the size of the helium droplets. The phenomenon, if it occurs at all, might be only transitory, ending when the helium rains out of the clouds. Also, any long-term value of Qj < l0 5 violates a lower limit noted by

336

RICHARD GREENBERG

Goldreich and Soter (1966) because it implies that Io is moving outward from Jupiter so fast that it must have emerged from the planet <4.5 byr ago. Thus to save the equilibrium hypothesis it appears necessary to accept a speculative interior model for Jupiter and assume that that dissipative mechanism only began to operate relatively recently. Yoder (1979) had noted that his resonance capture and equilibrium scenario would also explain the finite free libration amplitude of the Laplace relation (4~ oscillates about 180°) if capture occurred recently (<5 x 107 years). But there may be other ways to explain that free libration, either by asteroidal impact into one of the satellites (Greenberg 1980, 1982) or by a Sun-satellite resonance (Yoder and Peale 1981). The above discussion suggests that the system may not be in equilibrium. Two alternative evolutionary scenarios have been proposed (Greenberg 1980, 1982), each of which yields disequilibrium of the system with dominant dissipation in Io presently driving the system outward from deeper resonance. In one version, Io's interior alternately melts and refreezes in a feedback cycle with orbital variation. When Io melts, the higher amplitude of its tidal distortion allows greater heat production (as in Peale et al. 1979), so u can increase. The corresponding decrease in Io's orbital eccentricity then reduces the heating rate, so Io freezes again and the process can reverse and then repeat itself. At present the system may be at a point in the cycle where dissipation in Io dominates (/, > 0), rather than in the equilibrium state (b = 0) proposed by Yoder. This episodic heating scenario was originally speculative (Greenberg 1980, 1982), particularly in its geophysical assumptions, but recent modeling (Ojakangas and Stevenson 1984, 1986) shows that it may be reasonable. My other scenario for disequilibrium was that tides on Io may simply have always dominated the evolution, so that the system has monotonically evolved outward from

deep resonance (b > 0). An implication of that hypothesis is that the system may have started on the other side of exact commensurability (v < 0), evolved through exact commensurability (v = 0), and is now on its way out of resonance (toward v ~> 0). Both the episodic-heating scenario and this scenario with monotonically increasing v avoid the difficulties of Yoder's (/, = 0) equilibrium model. Yoder and Peale (1981) demonstrated an apparent serious difficulty, however, with any scenario that involves deep resonance in the past. Although analysis to first order in orbital eccentricities indicated that the Laplace relation is stable for small values of v as well as its present value (Greenberg 1981), Yoder and Peale's third-order analysis showed that 4~ is not stable at 180° for v less than about ~ of its present value. For such small v, eccentricities are sufficiently large that the higher-order terms can destabilize the configuration. In fact, this instability can also be shown with a second-order analysis (see Section II below). The system is stable for v < -0.5°/day, but at 4~ = 0° rather than 180° (Sinclair 1975, Yoder and Peale 1981). Thus, any evolution from negative v through deep resonance would have to jump the uncertain region in which 4~ is stable at neither 0° nor 180° (see Fig. 1). Similarly, the episodic-heating model would be suspect to the extent that it involves dipping into this "unstable" region. In this paper, I investigate possible behavior in deep resonance, and demonstrate that, while 4~ is not stable at either 0° or 180°, it is stable at intermediate values that vary with v. There are stable evolutionary paths in deep resonance, even where 4~ = 0 or 180° is unstable. Thus it may have been possible for the system to have evolved through deep resonance after all. One result to be demonstrated is that the forced eccentricities do not follow the same u dependence as for the equilibrium conditions at 4~ = 180°. I had originally hoped that the eccentricities might be smaller so

RESONANCE OF GALILEAN SATELLITES UNSTABLE EQ

360* STABLEEQ

~80~ UNSTABLE EQ

~

I /] "~i

0 STABLEEQ v<<0

II. V A R I A T I O N E Q U A T I O N S TO S E C O N D O R D E R IN E C C E N T R I C I T I E S

STABLEEQ )~-PRESENT VALUE

UNSTABLE EO V:0

337

V>>0

Flo. 1. Schematic of the v, ~ plane. The resonant system is in equilibrium for ~b = 0 ° and 180°, but Yoder and Peale (198 l) s h o w e d that these positions are stable only for v ~> 0?IS/day at ~b = 180° and v < -0?7/day at ~b = 0° (dark lines). Very near v = 0 (shaded area), forced eccentricities are too large to permit the type of analysis used in this paper. In Sections III and IV, stable paths are s h o w n to exist in deep resonance, qualitatively as s h o w n by the dashed lines.

that the analysis based on an expansion in e could be carried down to and beyond v - 0. In fact, h o w e v e r , the forced e ' s at these new stable equilibria are somewhat greater than on the ~b = 180° equilibria. Thus, there remains a region of uncertain behavior in very deep resonance. Nevertheless, I speculate that these new stable paths may provide a continuous connection from the stability at ~ = 0 ° for v < 0 to 4~ = 180° for v>0. In the following sections, I develop a formulation o f the equations for variation of orbital elements. Then I show how numerical integration starting at the u n s t a b l e equilibrium at small v and ~b = 180° can lead to oscillation about 4> ~ 135°. Stability at 135° is t h e n demonstrated analytically to second order in e. N e x t I map the stable paths in deep resonance by a numerical search for apparently stable conditions. And, finally, I discuss some implications of evolution from deep resonance that may fit the observed properties of the system, including the substantial free eccentricity of Ganymede, the free libration of th, the present secular acceleration of orbits, and the geological and geophysical states o f the satellites.

In this section, I introduce the method of analysis by considering perturbations only through second order in eccentricities. The general approach and notation are similar to the first-order analysis by Greenberg (1981). The objective is to define the notation, to outline the way the variation equations are derived, to show the breakdown of stability at 4) -- 180° for small v, and to indicate the presence of other stable values of 4~- In subsequent sections, I extend the analysis to third order to map these deep resonance stability paths more precisely. The following notation is used: n, a, e, &, h, and e are the standard Keplerian elements mean motion, semimajor axis, eccentricity, longitude of pericenter, mean longitude, and mean longitude at epoch, respectively. Subscripts 1, 2, and 3 refer to Io, Europa, and G a n y m e d e , O/12 ~ a J a 2 , and a23 -= a z / a 3 . A satellite's mass m as a fraction of the planet's mass M is t~. Due to the commensurability of the mean motions, the following critical arguments are slow varying 011------2h2 -- h i - - &l 022 ~ 2h3 -- h2 012 ~"

2h2 -- hl

¢b2

- - ¢~2

023 -- 2h3 - h2 - ¢b3. (1) In the stable Laplace configuration the upper two O's are fixed at 0 ° and the lower two at 180° . In reality, the small oscillations about stability cause 023 to circulate through 360 °, while the other three librate about their stable values. In a subsequent section, I argue that this oscillation supports the hypothesis that the system was previously in deeper resonance. Note that the angle 4~ is given by = 022 -- 01Z = h 1 --

3h2 + 2h3.

I also define the quantities Vll = - - ( 2 n z -- n i -- ~ols) 1122 =

-(2n3 - n2

-- t~2s)

(2)

338 /'12 =

RICHARD GREENBERG -(2n2 - nl

-

/'23 =

¢~2s)

J¢11 = -/z2ot12nl{-4.39h12 + 3.4hN} + vlnhll

-(2n3 -- n2 -- rO3s) (3)

where the subscript s refers to the secular variation o f pericenter due to Jupiter's oblateness, as well as secular c o m p o n e n t s of satellites' mutual disturbing potentials. According to C h a o (1976), cojs = 0?16/day, o~2, = 0?04/day, and ~3s = 0.°01/day; these values include a s e c o n d - o r d e r correction due to the r e s o n a n c e terms. The disturbing function representing the perturbing potential at Io due to E u r o p a is ( B r o u w e r and C l e m e n c e 1961, Y o d e r and Peale 1981) R = ( G m 2 / a 2 ) { - 1 . 1 9 e l cos 0in + 0.43e2 cos 012 - 0.58ele2 cos(0jl - 012) + 1.70el2 cos 2011 - 4.97ele2 cos(011 + 012) + 3.59e 2 cos 2012} + O(e3). (4) F o r the disturbance at E u r o p a due to Io, m2 is replaced by m l . F o r the E u r o p a - G a n y m e d e interaction subscripts 1 and 2 are replaced by 2 and 3, respectively. Direct interactions b e t w e e n Io and G a n y m e d e are of third order in e. N o t e that the numerical coefficients in (4) are actually functions of Ot12 and are here evaluated for a12 = 0.63. T h e y are not v e r y sensitive to the small variations in al2 that a c c o m p a n y the behavior I will be discussing. The variation in orbital elements is obtained f r o m L a g r a n g e ' s variation equations for the four elements n, e, &, and e for each of the three satellites. The analysis is perf o r m e d assuming coplanar, equatorial motion, so these 12 equations completely describe the two-dimensional motion. The analysis is greatly simplified by replacing e ' s and O's b y h ' s and k's as defined below: hij =-- ej sin

Oij

kij ~ ei

cos

Oij.

(7) ]'/23 =

-/x2n3{0.43 - 5.55k22 + 7.18k23} -- /'23k23

J¢23 = tz2n3{-4.39h22

(9) ill2 = -/zln2{0.43 - 5.55knl + 7.18k12} - /z3ot23n2{-1.19 cos ~b - 5.55k23 cos ~b + 4.39h23 sin ~b + 3.4k22 cos ~b - 3.4h22 sin ~h} - v n k l z

(10)

--/xlnz{--4.39hll + 7.18h12} - /z3ot23n2{-1.19 sin ~b - 4.39h23 cos ~b - 5.55k23 sin ~b + 3.4h22 cos qb

k12 =

+ 3.4k23 sin q)} + v12h12 (11) J~22 = -/zln2{0.43 cos ~b - 5.55ktl cos ~b - 4.39hll sin ~b + 7.18k12 cos ~b + 7.18hn2 sin ~b} - /z3a23n2{- 1.19 - 5.55k23 + 3.4k22} -

/-'22k22.

(12)

F r o m the definitions of h ' s and k's we see that k22 is a function of h12, k12, and h22 ; a variation equation for k22 would be redundant. Also note that, if the highest-order terms in e (or equivalently in h and k) are deleted from these variation equations, the equations reduce to those shown in Greenberg (1981). The right-hand sides of these seven variation equations contain the variables n and a but these quantities vary only slowly in the b e h a v i o r of interest and are here a s s u m e d to be constant. The right sides also depend on ~b and on the four v's. N o w vll and v12 differ by only a constant value. Similarly v23 is a simple function of v22. Therefore we only need to k n o w the b e h a v i o r of oh, v12, and v22 in order to solve the h and k variation equations ((6)-(12)). By definition (~ =

hln = -/x2a12nl{-1.19 - 5.55k12 + 3.4kll} vllkll

+ 7.18h23} + lJ23h23

(5)

The variation of these elements then becomes

--

(8)

(6)

/)12 -- /'22"

(13)

The variation of the v's is simply iq2 = hi /)22 =

2h2

(14)

ha - 2n3

(15)

-

RESONANCE OF GALILEAN SATELLITES where

339

p

i

0.2

0.4

i

hi = -3/z2cq2n~{- 1.19hll + 0.43h12 + 6.8kllh11 - 9.94klih12 - 9.94kl2hu + 14.36k12h12} O.E

h2 = -3/xln22{2.38hu - 0.86h12 - 13.6kllhll + 19.88kuh12 + 19.88k12hls - 28.7k12h12} - 3/z3ot23n2{- 1.19h22 + 0.43h23 + 6.8k22h22 - 9.94k22h23 - 9.94k23h22 + 14.36k23h23}

(Vdoy) 0.4

/

0.~

h3 = -3/z2nE{2.38h22 - 0.86h23 - 13.6k22h22 + 19.88k22h23 + 19.88k23h22 -- 28.7k23h23}. (16) Although the dynamical system is twelve-dimensional, with three bodies in planar motion, the nine equations ((6)-(11), (13)-(15)) form a closed set. Equation (12) is eliminated because h22 is a function of h12, k12, and ~b (see definitions and Fig. 2). The system is in equilibrium when the right sides o f all nine equations equal zero. One way to achieve that state is to have all four h's equal to 0. It follows by definition that 4) must be either 0 or 180°. Inspection shows that all k's are then 0. ~b is zero if Vl2 = v22 =- u. T h e n the three h equations ((6), (8), (10)) can be solved for the three values, k11, k12, and k23, that make the h's equal to zero. Thus for any value of u there is equilibrium at both ~ = 0 and 180°, and the

k22

e/h22'

hi2 ,k12

FIG. 2. By definition the vectors h12, k12, and h22 are of equal length e2 and are separated by the angle ~b -= 022 - 012.

01.6 °/day ls

FIG. 3. Eigenfrequencies for the system, evaluated to second order in e. Compare with the first-order analysis (Fig. 2 of Greenberg 1981) which does not show the instability for small v, and with the thirdorder analysis which more precisely defines the critical value of v (Yoder and Peale 1981). The instability was also identified by Wiesel (1981).

choice of u or ~b specifies the set of k values for equilibrium. The magnitudes of the k11, k~2, and k23 values are the forced eccentricities of the three satellites. For u > 0 and Ob = 180°, we find kll > 0, k12 < 0, and k23 < 0. In other words, 011 and 022 are in equilibrium at 0 °, and 012 and 023 are in equilibrium at 180° , in accord with the observed state of the system described in Section I. F o r Iv[ ~< 0?01/day, e2 = Ik121 ~> 0.1 and there is concern about validity of any analysis based on expansion in powers of eccentricity. Stability o f these equilibria is investigated by expanding about the equilibrium state to linearize the equations as was done in the first-order analysis (Greenberg 1982). Eigenfrequencies from that procedure are shown in Fig. 3. N o t e that in contrast to the first-order analysis (cf. Fig. 2 of Greenberg 1981), one f r e q u e n c y becomes imaginary for small v; i.e., the system goes unstable as found by Yoder and Peale (1981). The critical value of v found here is 0?26/day, larger than Y o d e r and Peale's value 0.°17/day, but I show in a subsequent section that their value is confirmed by extension to third or-

340

RICHARD GREENBERG

der in e. Similarly, one can confirm the stability at ~b = 0 ° for u < 0. III. S T A B I L I T Y

AT OTHER

VALUES

O F ~b

In order to understand the nature of the instability, I numerically integrated the nine equations of motion starting at an unstable equilibrium condition with ~b = 180° and v small. As expected, the system diverged exponentially from the initial condition. Then, however, it began to oscillate about a state with ~b ~ 135 °. Similar behavior had occurred in a numerical integration by a completely independent method by Wiesel (1981). One example of this type of behavior is shown in Fig. 4. In this case, the system was started at the unstable equilibrium at v = 0.15°/day and 4) = 180°. Although only the behavior of ~b is shown, r e m e m b e r that a full nine-dimensional system was integrated. Also, the oscillations about 135 ° continued for - 1 0 0 years although less time is shown in the figure. This remarkable numerical result led me to investigate analytically the possibility that there is an equilibrium state at th = 135° , which would represent an asymmetrical configuration. Careful evaluation of the mean values o f all nine variables in several test integrations showed that the center of oscillation appeared to have ~b = 135°, h22 = h23 = 0, 1'12 = /'22, and h12 = kl2. When these mean conditions were taken as initial conditions, numerical integration showed only very small oscillations about these values even for integration over several hundred

UNSTABLE

EQUILIBRIUMAT 180"

135°

50100 Tilde (DAYS)

10,000

FIG. 4. T h e i n s t a b i l i t y at ~b = 180 ° leads to o s c i l l a t i o n about ~ = 135° (cf. Wiesel 1981).

years, suggesting a stable equilibrium with ~b at 135 °. The equilibrium at ~b = 135° can be demonstrated analytically. First, I set v12 = u22 v, h23 = 0, and kjz = h~z. It follows that h22 = 0 and k22 = X/2kl2. With these choices, ~b and k23 are zero, according to (9) and (13). That leaves seven more variation equations to be satisfied with the appropriate choice of the remaining five unspecified variables, which turns out to be possible in this case because two equations are redundant. First, with h23 and h22 equal to zero, we see that/~J2 is just a multiple of ~'22, SO 1)22 = 0 is automatically satisfied if i'12 = 0. Also, algebraic manipulation shows that kl2 is a linear combination of the other equations. Thus only five independent variation equations need be satisfied by the five variables. The numerical solution is k~z = - 2 . 5 4 × 10 -2, kll --- +1.08 × 10 -2, k23 = - I . 8 4 × 10 -3, hli = 1.23 × 10 -3, and u = 0?18/day. Note that unlike the equilibrium at 4~ = 180 °, which exists o v e r a continuous range of v values, the equilibrium at ~b = 135° exists only at a single value of v. But numerical integration indicates that this new equilibrium is stable. Having found this first example of apparently nonsymmetrical (~b ~ 0 ° or 180°) stability, I numerically explored the neighborhood of that equilibrium condition by changing initial conditions slightly. The result was oscillations about a neighboring point in the nine-dimensional space. It too seemed to be a stable equilibrium point. Further numerical exploration in the same and opposite direction in nine-dimensional space revealed a continuous string of stable points. In the v, ~b plane (cf. Fig. 1) these points seemed to form loci (shown schematically by the dashed lines in Fig. 1) branching from the critical point on ~b = 180° where the system goes unstable, along two paths (symmetrical about ~b = 180 °) with ~b 180°, toward smaller v. Section IV describes results of a more precise computation of these paths.

RESONANCE OF GALILEAN SATELLITES None of these equilibria are really stable, strictly speaking. In the case of stability at ~b = 180° (e.g., the present condition of the system), the nine-dimensional system is really only in neutral equilibrium: An infinitesimal perturbation in the right direction (a small change in v's and k's) would put the system in a new equilibrium, rather than yield oscillation around the original equilibrium. All the equilibria on this branch are at ~b = 180°, so that value of ~b is stable, even though the system itself is strictly only in neutral equilibrium. It is this neutral equilibrium that allows the system to evolve along the ~b = 180° equilibrium track under the very small effects of tidal torques as discussed by Yoder (1979) and Greenberg (1980, 1982). Similarly, the equilibria along the ~b :# 180° track are also neutral, with small perturbations allowing the system to move to adjacent equilibria. I have speculated (Greenberg 1984) that this newly discovered set of equilibria may represent an evolutionary track that the system might have taken from the known stability at 6 = 0° (with v < 0) to the stability at ~b = 180° (with v > 0). The trend of the new stable branch seems likely to connect the 0 ° and 180° stability, However, before carrying this discussion further, it is necessary to extend the analysis to one order higher in e to confirm these general results. IV. EXTENSION

TO THIRD ORDER IN e

The analysis is extended to third order by using the full expansion of the disturbing function to that order. The algebra becomes extremely lengthy and complicated, and is not displayed here because in principle it is similar to the analysis to second order. The form of the disturbing function R used is identical to the long-period terms displayed by Yoder and Peale (1981, their Eq. (12)). Yoder and Peale neglected the perturbations at Io directly due to Ganymede and vice versa, although such perturbations do contain long-period terms of third order in e. They considered perturbations only of

341

adjacent satellites. I have included the direct Io-Ganymede terms in the work presented here. These terms do not qualitatively affect Yoder and Peale's results. Another third-order effect explicitly neglected by Yoder and Peale is the variation of e's, the mean longitudes at epoch. These contribute terms in the equations for h's and k's of the same order as other terms due to third-order parts of R. Again, although Yoder and Peale were formally inconsistent to neglect such terms, the effect is not large. I also include the terms kl - 3k2 + 2e3 in the equation for ~b (cf. Eq. (13)), to the same order as the b's were included in the h and k equations. Numerical integration of the higher-order versions of the nine governing equations confirmed the oscillations about points lying on a line with ~b :/: 180° or 0° in ~b, v space. However, a slow systematic drift along this line occurred on a time scale much slower than the oscillations. As discussed in the next section, such a systematic drift is to be expected if there is energy loss in the system, but constant angular momentum. With constant energy, this drift in the numerical integration corresponds to a systematic increase in the angular momentum (v is increasing and e is decreasing). The reason for the drift was that, in fact, the equations of motion even to third order in e (or to any finite power of e) do not necessarily conserve energy. Remember the expression for angular momentum includes the factor (1 - e2) 1/2, which corresponds to an infinite power series in e. Thus, the lack of fourth-order terms in my equations of motion allowed the slow artificial drift along the evolutionary path. In order to prevent this artificial drift, I simply modified the expressions for dnJdt and dn3/dt so that the equations for conservation of total energy and total angular momentum are exactly satisfied, given the unchanged expressions for dn2/dt, and h's and k's. The new formulae for dnJdt and dn3/dt are the same as the old ones to third order in e, but now the physical constraints of

342

RICHARD GREENBERG

energy and angular m o m e n t u m conservation are satisfied. Numerical integration was then used to search for apparent stable equilibrium points. The method involved selection of a set o f values for the nine variables expected to be near equilibrium. Numerical integration o v e r a period o f about 100 years would give oscillation of the elements about some set of mean values. The mean values were then taken as new initial conditions for numerical integration, with the general result that the amplitude of oscillations is considerably reduced. This process was iterated until a set of values was found at which the integration gave no change to a few significant figures. This was taken to be a stable equilibrium condition. A set of such equilibria was generated in this way, which mapped out a path of stability (Fig. 5). Strictly speaking, I have not proven stability on this path. Possibly longer-term behavior might go unstable. But in all the integrations for all the points studied (cumulative time of - 10,000 years), there was n e v e r any indication of instability for this class o f equilibrium condition. In Fig. 5, I show the equilibrium points in 4o, v~2 space. The region for q5 < 200 ° only is shown. Inspection o f the nine equations of motion shows that there is a path for ~b > 180° symmetrical about ~b = 180°. The curved branches with ~b ~ 180° are apparently stable, as is the branch to the right of the bifurcation point (v 12 = 0?175/day) on th = 180°. To the left of the bifurcation point the equilibrium at ~b = 180° is unstable, in accord with Y o d e r and Peale (1981). In the upper portion o f Fig. 5, I show the forced eccentricities (i.e., the equilibrium values of e) for each of the three satellites as a function of v12. The solid lines are for the stable equilibria and the dashed line for e2 is for the unstable equilibria. Unfortunately, the high values of e2 for Vl2 < 0?/day preclude application of this computational method in very deep resonance. The values of the h's, k's, and v's on the stable equilibrium paths are shown in Fig. 6 as functions

I

I

I

0.08 0.06 e

-,.<.___

0.04 0.02

e~ e3

0

200 ° 180 ° 160 °

+

140° 120 ° I00 °

UNSTABLE ) ~

/ I 0.1

/~-

STABLE

Bp/FUN.RCAT'O"

I 0.2

I 0.5 °/day

1/12 FIG. 5. Stable equilibrium p a t h s in tb, /)12 space. T h e s e are s y m m e t r i c a l a b o u t tb = 0 ° a n d 180 °. C o m p a r e with the s c h e m a t i c v e r s i o n in Fig. 1. Also s h o w n are the c o r r e s p o n d i n g forced eccentricity values. F o r e:, the d a s h e d line s h o w s the value on the unstable equilibrium with tb = 180 ° and small u.

of u~2 for the 4) = 180° case and as functions of 4~ for the + 4:180 ° branch. Note the continuity of values approaching the bifurcation point from both directions. Note also that on the ~b 4:180 ° track, I could not actually find any points closer to the bifurcation point than at ~b = 176 °, probably because the "potential well" so close to bifurcation must be very shallow; there is an unstable equilibrium (~b = 180 °) very near. Numerical integrations in this area tended to oscillate across ~b = 180°. My trial and error search method simply has small probability of catching a stable equilibrium in this region. The conclusion of this analysis is that

RESONANCE OF GALILEAN SATELLITES

343

// I

I

I

I

I

I

i

I

,/

ii

i

::::: :2

-3.2 -3.0

a.~~h, ~z

:28

JO

"2.6

k23 \

~

1

0/ 2.4

"2.4

i

",,e 2.6

i

i 2.8

~, \ ,

:;; J3 ,,,

~,

3.0

3.2

rr

4

vj 2 (163rodldoy)

FIG. 6. Values of h's, k's, and u's on the stable equilibrium paths, shown as functions of u12 for the ~b = 180° case, and as functions of~b for the 4~ ¢ 180° branch. Units are 10-2 for klt, - 10-2 for hi2 and k12, - 1 0 -3 for k23, and 10 -4 for h u . F o r example, k23 = -0.0023 on the left side of the figure.

there does seem to be a stable path in deep resonance (u~2 < 0?175/day). The stable path involves ~b :~ 180° and as such involves asymmetrical conjunction configurations, an unusual celestial mechanical result. The path seems to trend from the bifurcation point on 4) = 180° toward the stable configurations with 4) = 0 at /-:12 < 0 , but such a continuous connection is uncertain due to the b r e a k d o w n of the m e t h o d of analysis at u12 -~ 0?l/day. In the next section we discuss possible evolution along this path. V. E V O L U T I O N F R O M D E E P R E S O N A N C E

The presence of stable paths along which the system can evolve given small perturbations such as tidal effects opens the possibility that the system has evolved through deep resonance. Or in the case o f the episodic hypothesis, it is conceivable that as v periodically dips down to low values, the system passes through the bifurcation point and onto one o f the 4) ~ 180° tracks. In this section, I consider the rates o f evolution and the implications for the history of the system. Yoder (1979) showed how the rate of variation of u depends on tidal dissipation

in Io and Jupiter under present conditions (moderately large u). Yoder and Peale (1981) generalized this expression by using energy and conservation o f total angular m o m e n t u m L to include the significant variation of eccentricities in deep resonance. H o w e v e r , since they c o n s e r v e d angular m o m e n t u m , their formula cannot be used to include Jovian tides. Their formula can be further generalized by manipulation of the equation for total energy and L to the form d E = ~ E ~ ( - O . 6 6 d e 2 - 0.44de22 + 0.84du/n2 -

1.38dL/(

1.74de]

GV'-G--M~alm l ) )

(17)

where E is the total orbital energy of the system and Ei is the orbital energy of a single satellite: Ei =

x~

--~n~il~

~2~2 i t~ i .

(18)

(There is an interesting numerical curiosity that the mass and semimajor axis ratios are such that E1 = E2 + E3 and 2E~ ~ E3 .) Nevertheless, as long as tides raised on Io dominate the evolution, a reasonably accurate rate relation is obtained by setting d L / d t = 0. Given the now known relation between the e's and u we can relate the rate o f energy dissipation in Io to 1:,. At the present

344

RICHARD GREENBERG

rate (7 x 1013 W -¢- 25%, McEwen et al. 1985), the system would have evolved out from deep resonance in only --10 7 years (Yoder and Peale 1981). Thus the scenario of simple passage through and out of resonance suffers from two major weaknesses: (1) That such shortlived evolution took place recently is rather improbable. (As discussed in the Introduction, the equilibrium hypothesis may have similar difficulties.) (2) While this paper shows some possible phenomena in deep resonance, the dynamics of passage through exact resonance (v = 0) are not yet understood. Both of these difficulties can be avoided, while still retaining the disequilibrium hypothesis with the system presently evolving outward from resonance, if we adopt the episodic scenario described in the Introduction. Dissipation in Io alternately turns on and off, so that rapid recent evolution outward from deep resonance could be interpreted as only part of a periodic cycle. Recall that the advantage of this idea over the equilibrium hypothesis is that it does not require improbably small values for Q j. On the other hand, the episodic model cannot work for arbitrarily large values of Qj. If Q~ is much greater than 107, the phase of the cycle with Io turned off would last for the age of the solar system--hardly episodic. With Qj = 1 0 7 that phase would last --10 9 years; the probability that we would observe Io during its 107-year hot phase would be only 1%. Perhaps Qj ~ l06 provides an acceptable compromise, being not too far from most Jovian models, while giving a reasonable probability (10%) in the context of the episodic model for a presently active Io with evolution outward from resonance. Rapid evolution out of deep resonance provides a mechanism for stimulating the present free eccentricity of Ganymede. Suppose the evolution has progressed from deep resonance (v < 0.°175/day), through the bifurcation point to its present state u ~0.°74/day. Equation (17) above shows that continuous energy dissipation corresponds

to continuous variation of v. Now consider the th :# 180° path in Fig. 5. As it approaches the bifurcation point, it becomes perpendicular to the ~b = 180° line. Thus for continuous variation of v there must be a discontinuous jump in ~b. This jump must stimulate oscillations in the nine-dimensional system, including the h's and k's, as well as ~b. These oscillations are highly coupled deep in resonance. In other words, there would have been an induced free eccentricity superimposed on the forced eccentricity, and this would have happened ~ 1 0 7 years ago. Tides raised on a satellite generally tend to damp orbital eccentricity (Goldreich 1963). For the Galilean satellites, the orbital resonance keeps eccentricities up to the forced value, while tides damp the free (or oscillatory) component of the eccentricity. The rate of damping (Goldreich and Soter 1966) is b/e ~ - l O ( k / Q ) ( R / a ) 5 t x - J n ,

(19)

where k is the tidal Love number (inversely proportional to rigidity of the material) and R is the satellite's radius. For Io with its rapid dissipation, Q / k is effectively < 10, so the time scale r for damping e is < 104 years. Its free eccentricity must be damped shortly after stimulation. For Europa, assuming the reasonable value Q / k ~ 500 (Cassen et al. 1982), r ~ 107 years, most of the free e2 (whether primordial or induced at the bifurcation point) would be damped, consistent with the low present value. For Ganymede, with the rigidity of ice and Q ~ 100, free eccentricity would damp in ~"- 108 years. Even with rigidity an order of magnitude greater than ice (like rock), 7 would be less than the age of the solar system. Any primordial free eccentricity would probably have decayed away by now. But if the system has passed through the bifurcation point on its way out from deep resonance, the free eccentricity thus s t i m u l a t e d - 1 0 7 years ago would still remain. Thus, a substantial free eccentricity and corresponding circulation of 023 would

RESONANCE OF GALILEAN SATELLITES be expected at present. That state is exactly what is observed as described in the Introduction. For Ganymede, the free eccentricity is more than twice the present forced component; its value is very close to (and may be a fossil of) the forced value at the bifurcation point (Fig. 5). Similarly, the stimulation of oscillations at the bifurcation point could explain the present free libration of ~b as well. However, as discussed in the Introduction, there could be other ways to explain it. In any case, the present free oscillations in the system are consistent with evolution from deep resonance. The competing equilibrium and disequilibrium hypotheses make distinct predictions regarding the present rates of orbital acceleration. Greenberg et al. (1986) combined various assumptions about the ratio of dissipation in Io to that in Jupiter, with equations for rates of total orbital energy and angular momentum dissipation in the system, to obtain equations relating orbital acceleration to dissipation in Io. If Io were not involved in a resonance, and if dissipation in Jupiter were negligible, energy dissipation in the satellite at rate E would be related to orbital acceleration by

E/E1

= ~hl/nl.

(20)

Note that both/~ and El are negative. Because Io shares orbital energy and angular momentum with the other satellites via the Laplace resonance, this expression must be modified slightly to

i~/El = 0.56hJnl.

(21)

In that case, with negligible dissipation in Jupiter, the system would be evolving rapidly outward from deep resonance. Now if tidal dissipation in Jupiter is sufficiently great to give the equilibrium hypothesized by Yoder (1979), we have the very different relation

i~/El = -- 1.54hl/nl

(22)

345

(Greenberg et al. 1986) where E still represents dissipation in Io. For the observed value of/~ (McEwen et al., 1985), the assumption of little dissipation in Jupiter (disequilibrium with evolution outward from deePer resonance) gives hl/nl -- 3 X 10-l° year -I (from Eq. (21) above). This value is in good agreement with determinations of hi from the historical observational record by DeSitter (1928) and by Goldstein and Jacobs (1986). However, Lieske (1986) believes that there are important timing errors in those analyses. Lieske finds the value hi~n1 = ( - 7 . 4 - 8.7) × 10-12 year -1 consistent with the negative sign predicted by the equilibrium hypothesis (Eq. (22) above) but too slow by more than an order of magnitude. Lieske's value thus implies that dissipation in Jupiter is too slow to maintain Yoder's proposed equilibrium model. Thus, either determination of hi (Goldstein and Jacobs' or Lieske's) combined with the measured Io heat flux indicates disequilibrium, with evolution outward from deep resonance. The logical situation is summarized schematically in Fig. 7. The equilibrium hypothesis requires a specific ratio between the tidal dissipation rates (diagonal solid line). The rate of dissipation in Io is known, perhaps, from the heat flux. The intersection of these requirements predicts a dissipation rate in Jupiter. The latter rate is faster than produced by most interior models, with the exception of Stevenson's (1983) helium clouds. The value of hi found by Goldstein and Jacobs (1986) is consistent with very little dissipation in Jupiter, and hence evolution outward from resonance. The value of hi found by Lieske (1986) requires much more dissipation than given by most interior models, but distinctly less than is consistent with the equilibrium model. Neither determination of h~ would be consistent with equilibrium. It is evident from inspection of Fig. 7 that Lieske's determination of n~ could only be consistent with the equilibrium hypothesis if the true energy dissipation rate in Io time

346

RICHARD GREENBERG

ENERGY DISSIPATION IN ,JUPITER EXPRESSED

AS WQj

.

.

.

.

.

.

.

.

.

.

.

.

.

.

.

.

MOST JUPITER, MODELS

io"w

iol2w ,oI3w ~

ioISw

ENERGY DISSIPATION RATE OBSERVEDRATE IN I0

FIG. 7. Constraints on tidal dissipation in Io and in Jupiter, from Yoder's (1979) equilibrium model, from infrared observations (McEwen et al. 1985), and from measurements of Io's orbital acceleration (differing values from Goldstein and Jacobs 1986, and from Lieske 1986). If the system lies below the line for the equilibrium (/, = 0) model, it is evolving outward from resonance (v > 0).

averaged over the 300-year orbital observation interval were at least an order of magnitude less than the value inferred from the thermal flux. That would require great geophysical changes in Io in a remarkably short time. One argument in favor of the equilibrium model is the determination by Lieske (1986) of ]i,/nl] < 0.5 × 10-]l year -j. On the other hand, if the infrared-observed Io dissipation rate is relevant and Lieske's ht is correct, then total energy and angular momentum constraints require i,/nj ~ 4 × 10 -1] year -t , in which case the system is out of equilibrium evolving outward from deeper resonance. Finally, I note geological and geophysical properties of the satellites that seem consistent with substantially greater heating rates in the past than at present. This evidence also supports the disequilibrium hypothesis with evolution outward from deeper resonance, because in deeper resonance the forced eccentricities, and hence tidal heating rates, for all three satellites would have been enhanced. First, Ganymede with its grooved terrain has clearly undergone

greater endogenic evolution than the otherwise similar Callisto. Resonant tidal heating would have contributed to this difference, although the greater radioactive heating of Ganymede also may have been important (Cassen e t al. 1982). Second, the very smooth surface of Europa suggests significant tidal heating. Past residence in deep resonance would have allowed melting of a thick water layer, which would have tended to level the thin ice layer above it. Once such melting occurred, it would tend to promote increased dissipation due to the higher amplitude of tidal distortion, and would thus sustain a liquid water layer even as forced e2 decreased to its present value. If the system evolved from outside resonance (smaller e2), as in Yoder's hypothesis, melting would probably have been unable to get underway (Cassen e t al. 1982). Third, the removal of sulfur, which is so much in evidence near the surface of Io, from FeS in the interior likely occurred early in the geophysical history, before core formation which would have carried sulfur into the depths. Such removal would have required an early heat source. The extreme tidal heating in deep resonance could have served that purpose. The disequilibrium hypothesis seems to fit the relevant observable characteristics of the Galilean satellite system. It avoids some of the problems with the equilibrium model, especially the need for very rapid tidal energy dissipation in Jupiter. It does require recent passage out of deep resonance. The stable evolutionary paths in deep resonance described in this paper eliminate a major objection to such models. There remains uncertainty about behavior in very deep resonance. These difficulties are avoided in the episodic model, with present evolution outward from deep resonance representing one phase of periodic evolution. Although the episodic model contradicts the equilibrium hypothesis, it does not necessarily preclude earlier capture into resonance.

RESONANCE

OF GALILEAN

ACKNOWLEDGMENTS I thank S. J. Peale and J. Henrard for valuable reviews of this manuscript. This research is supported by Grant NAGW-944 from the NASA Planetary Geology and Geophysics Program. REFERENCES BROUWER, D., AND CLEMENCE, G. M. 1961. Methods o f Celestial Mechanics. Academic Press, New York. CASSEN, P. M., S. J. PEALE, AND R. T. REYNOLDS 1982. Structure and thermal evolution of the Galilean satellites. In Satellites o f Jupiter (D. Morrison, Ed.), pp. 93-128. Univ. of Arizona Press, Tucson. CHAO, C.-C. 1976. A General Perturbation Method and Its Application to the Motion o f the Four Massive Satellites o f Jupiter. Doctoral dissertation, UCLA. DESITTER, W. 1928. Leiden Annals, Vol. 16, Part 2. GOLDREICH, P. 1963. On the eccentricity of satellite orbits in the solar system. Mon. Not. R. Astron. Soc. 126, 257-268. GOLDREICH, P., AND S. SOTER 1966. Q in the solar system. Icarus 5, 375-389. GOLDSTEIN, S. J., AND K. C. JACOBS 1986. The contraction of Io's orbit. Astron. J. 92, 199-202. GREENaERG, R. 1976. The motions of satellites and asteroids. In Jupiter (T. Gehreis, Ed.), pp. 122-132. Univ. of Arizona Press, Tucson. GREENaERG, R. 1980. Orbital evolution of the Galilean satellites. Presentation at the Satellites o f Jupiter Conference, Kona, Hawaii. GREENBERG, R. 1981. Tidal evolution of the Galilean satellites: A linearized theory. Icarus 46, 415-423. GREENBERG, R. 1982. Orbital evolution of the Galilean satellites. In Satellites o f Jupiter (D. Morrison, Ed.), pp. 65-92. Univ. of Arizona Press, Tucson.

SATELLITES

347

GREENBERG, R. 1984. Stability of the Laplace relation in deep resonance. Bull. Amer. Astron. Soc. 16, 687. GREENaEgG, R., S. J. GOLDSTEIN, AND K. C. JACOBS 1986. Orbital acceleration and the energy budget in the Galilean satellite system. Nature (London) 323, 789-791. LIESKE, J. H. 1986. Galilean Satellite Evolution: Observational Evidence for Secular Changes in Mean Motions. Preprint. McEWEN, A. S., D. L. MATSON, T. V. JOHNSON, AND L. A. SODERBLOM 1985. Volcanic hot spots on Iv: Correlation with low-albedo calderas. J. Geophys. Res. 90, 12,345-12,379. OJAKANGAS, G. w., AND D. J. STEVENSON 1984. Episodic heat flow in tidally heated satellites. Bull. Amer. Astron. Soc. 16, 661. OJAKANGAS, G. W., AND D. J. STEVENSON 1986. Episodic volcanism of tidally heated satellites with application to Io. Icarus 66, 341-358. PEALE, S. J., P. CASSEN, AND R. T. REYNOLDS 1979. Melting of Io by tidal dissipation. Science 203, 892894, SINCLAIR, A. T. 1975. The orbital resonance amongst the Galilean satellites of Jupiter. Mon. Not. R. Astron. Soc. 171, 59-72. SMITH, B. A., AND THE VOYAGER IMAGING TEAM 1979. The Jupiter system through the eyes of Voyager 1. Science 204, 951-972. STEVENSON, D. J. 1983. The origin and evolution of the great resonance in the Jovian satellite system. J. Geophys. Res. 88, 2445-2455. WmSEL, W. E. 1981. The origin and evolution of the great resonance in the Jovian satellite system. Astron. J. 86, 611-618. YODER, C. F. 1979. How tidal heating in Io drives the Galilean resonance locks. Nature 279, 767. YODER. C. F., AND S. J. PEALE 1981. The tides of Io. Icarus 47, 1-35.