{'fraction_non_alphanumeric': 0.03486692348100157, 'fraction_numerical': 0.019085974435300298, 'mean_word_length': 4.27221324717286, 'pattern_counts': {'":': 0, '<': 17, '<?xml version=': 0, '>': 7, 'https://': 0, 'lorem ipsum': 0, 'www.': 0, 'xml': 0}, 'pii_count': 0, 'substrings_counts': 0, 'word_list_counts': {'cursed_substrings.json': 3, 'profanity_word_list.json': 0, 'sexual_word_list.json': 0, 'zh_pornsignals.json': 0}}
Migration of a moonlet in a ring of solid particles : Theory and application to Saturn's propellers

Hundred meter sized objects have been identified by the Cassini spacecraft in Saturn's A ring through the so-called"propeller"features they create in the ring. These moonlets should migrate, due to their gravitational interaction with the ring ; in fact, some orbital variation have been detected. The standard theory of type I migration of planets in protoplanetary disks can't be applied to the ring system, as it is pressureless. Thus, we compute the differential torque felt by a moonlet embedded in a two-dimensional disk of solid particles, with flat surface density profile, both analytically and numerically. We find that the corresponding migration rate is too small to explain the observed variations of the propeller's orbit in Saturn's A-ring. However, local density fluctuations (due to gravity wakes in the marginally gravitationally stable A-ring) may exert a stochastic torque on a moonlet. Our simulations show that this torque can be large enough to account for the observations, depending on the parameters of the rings. We find that on time scales of several years the migration of propellers is likely to be dominated by stochastic effects (while the former, non-stochastic migration dominates after ~ 10^{4-5} years). In that case, the migration rates provided by observations so far suggests that the surface density of the A ring should be of the order of 700 kg/m^2. The age of the propellers shouldn't exceed 1 to 100 million years, depending on the dominant migration regime.


Introduction
The theory of disk-satellite interactions has for the most part been developed after Voyager's encounter with Saturn.
The satellites that orbit beyond the outer edge of the rings perturb the dynamics of the particles composing the rings. This leads to an exchange of angular momentum between the rings and the satellites, and to the formation of density waves in the rings. Lin and Papaloizou (1979) and Tremaine (1979, 1980) have calculated, using two different methods, the total torque exerted by a satellite on a disk interior (or exterior) to its orbit. This torque is called the one-sided Lindblad torque, because it can be computed as the sum of the torques exerted at Lindblad resonances with the secondary body. In the lowest order, local approximation, the inner and outer torques are equal and opposite : the torque exerted by a satellite on a disk located inside its orbit is negative, with the same absolute value as the positive torque exerted on a disk located outside the orbit.
Reciprocally, a disk exerts a torque on the secondary body. When the strictly local approximation is relaxed, the inner and outer torques are not exactly equal and opposite (Ward 1986). Their sum, called the differential Lindblad torque, is generally negative. As a consequence, the orbital angular momentum of a body embedded in a disk decreases, and so does its semi-major axis (on circular Keplerian orbits, the orbital angular momentum is proportional to the square root of the semi-major axis). This is planetary migration of type I (Ward 1997). So far, this phenomenon has been mainly studied in the frame of planets embedded in protoplanetary gaseous disks (see Papaloizou et al. 2007, for a review).
The Cassini spacecraft has been orbiting the Saturnian ring system since 2004, offering the possibility to observe the coupled evolution of the ring system and the satellites. Due to short orbital timescales (1 year is equivalent to about 700 orbits of the A ring) it may be possible to observe the exchange of angular momentum between the two systems. One of the most striking discoveries of the Cassini spacecraft is the observation of propeller shaped features in the A ring (located between 122 000 and 137 000 km from Saturn), with longi-tudinal extent about 3 km (Tiscareno et al. 2006(Tiscareno et al. , 2008Sremčević et al. 2007). They are most probably caused by the presence of moonlets about hundred meters in size, embedded in the ring, and scattering ring particles (Spahn and Sremčević 2000). As they are embedded in the ring, these small bodies should exchange angular momentum with the ring, and migrate ). This migration could be detected by Cassini observations through the cumulative lag, or advance with time t of the orbital longitude φ induced by a small variation of the semi-major axis and the angular velocity Ω (δφ = δΩ × δt), offering for the first time the possibility to confront directly the planetary migration theory with observations, and to give insights and constrains on the physical properties of the rings and of the moonlets.
In this paper, we address the question of the theoretical migration rate of these propellers, using both numerical and analytical approaches. The theory is then confronted to observations. In Sect. 2, we review the standard theory of type I migration ; the differences between migration in protoplanetary disks and in Saturn's rings are explained, showing the need for a new calculation of the migration rate of embedded moonlets. This rate is given in Sect. 3 for an homogeneous, axisymmetric disk with a flat surface density profile, as a result of numerical computation in Sect. 3.1, and analytical calculation in Sect. 3.2. In Sect. 4, we consider the effect of density fluctuations in the rings, in particular the role of short-lived gravitating clumps, also called gravity wakes, which are known to be numerous in the A ring (Colwell et al. 2006). We then conclude in Sect. 5 on what our model tells us on the properties of the rings, given the observed migration rates.

Review of Type I migration and the differential Lindblad torque
In protoplanetary gaseous disks, the perturbation caused by a terrestrial planet leads to the formation of a one armed spiral density wave, leading the planet in the inner disk, and trailing behind the planet in the outer disk. This wave is pressure supported and generally called the wake, but it has nothing to do with the gravity wakes mentioned above : the latter are local features, while the planet wake spirals through the whole disk. The planet wake carries angular momentum, so that the angular momentum given by the planet to the disk is not deposited locally (e.g. Crida et al. 2006, Appendix C). Therefore, the disk profile is hardly modified in this linear regime. Still, the negative torque exerted by the outer disk on the planet through the wake is larger in absolute value than the positive torque from the inner disk. Without going into the details (for which the reader is refereed to Ward 1997), the main reason the outer disk wins over the inner disk lies in pressure effects : from the dispersion equation of a pressure supported wave, one finds that the location r L,m where the wave with azimuthal mode number m, corresponding to the mth Lindblad resonance with the planet, is launched, is not exactly the location of the resonance given by Kepler's laws. The shift is not symmetrical with respect to the planet position for inner and outer resonances, but favors the outer ones. As a consequence, the planet feels a negative total torque, called the differential Lindblad toque, and given by Tanaka et al. (2002) : where the index p refers to the planet, r p being the radius of its orbit and Ω p its angular velocity, q is the planet to primary mass ratio, Σ is the surface density of the disk in the neighborhood of the planetary orbit. Finally, h is the aspect ratio of the disk, being the ratio between its scale height and the distance to the central body, being proportional to the square root of the gas pressure. In a typical protoplanetary disk, h ≈ 0.05. The numerical coefficient C is given by C = 2.340 − 0.099ξ, where ξ is the index of the power law of the density profile : Σ ∝ r −ξ . This result is robust. In particular, the value of the negative torque is almost independent on the slope of the density profile ξ. This is due to the socalled pressure buffer : the resonances are shifted when the density gradient varies (Ward 1997). If the disk were pressureless, then the expression of C would be completely different. Also, the aspect ratio h in Eq.
(1) appears because r L,m doesn't converge towards r p when m tends to infinity, but towards r p (1 ± 2h/3), due to pressure effects. To sum up, the gas pressure plays a fundamental role in type I migration.
It should be mentioned for completeness that, in addition to the differential Lindblad torque discussed above, the horseshoe drag -exerted on the planet by the gas on horseshoe orbits around the planetary orbit -plays a significant role in type I migration (see e.g. Ward 1991;Masset 2001;Baruteau and Masset 2008;Kley and Crida 2008;Paardekooper and Papaloizou 2009;).
In contrast to gaseous protoplanetary disks, pressure effects in Saturn's rings are not important. The aspect ratio h is of the order of 10 −7 . The spiral density waves that are observed in the A ring are gravity supported, not pressure supported. Thus, the standard theory of type I migration does not apply. In particular, the spiral planet wake doesn't appear. The interaction of the moonlet responsible for the propeller structure with the disk is observed to take place within a few hundred kilometers. Resonances with m 10 3 are located within this distance and should play a significant role. However, in the standard type I migration the important resonances have m ∼ 1/h ∼ 10 7 . Therefore, Eq. (1) can't be directly applied to a moonlet in Saturn's rings. A new approach is needed, adapted to the two main characteristics of the problem : the fact that the rings are made of solid particles, and the fact that the interaction is taking place very close to the moonlet.

The ring-moonlet interaction
In this section, we compute the interaction between a moonlet and a ring test particle. In this analysis, the gravity of the other ring particles is neglected. This leads to the torque exerted on the moonlet by an initially unperturbed, homogeneous ring. In subsection 3.1, the computation is performed numerically. In subsection 3.2, it is derived analytically. Both results are in agreement, and a corresponding migration rate for the moonlet is given and discussed in subsection 3.3.
From now on, m denotes the mass of the moonlet (and not anymore the order of a resonance). The moonlet, has a circular orbit of radius r m around the central planet of mass M . The gravitational potential due to the moonlet is Ψ. The radial and azimuthal components of the equation of motion of a ring particle in two dimensional cylindrical polar coordinates (r, φ) are and The angular velocity of the moonlet is ω = GM/r m 3 . Let r 0 be the radius of the initially circular orbit of a test particle and Ω = GM/r 0 3 its angular velocity. We note b = r 0 − r m is the impact parameter, andb = b/r H the normalized impact parameter, where r H = r m m 3M 1/3 is the Hill radius of the moonlet.

Numerical computation of the ring moonlet interaction
In this section, the numerical integration of the above equations of motion is performed, in order to find the trajectories of ring particles in the presence of a perturbing moonlet in the frame corotating with the moonlet. To measure the tiny asymmetry between the inner and the outer part of the ring, the full equations are integrated, without linearization or simplification. A Bulirsch-Stoer algorithm (Press et al. 1992) is used, and a Taylor expansion is performed in the code when necessary to subtract accurately large numbers, in order to achieve machine double precision (10 −16 ). We have checked that the Jacobi constant is conserved to this precision along the trajectories. Examples of obtained trajectories are given in Fig. 1.
It is well known within the framework of the restricted 3-body problem that if |b| is small enough, the test particle has a horseshoe shaped orbit, while if |b| is larger than ∼ 2.5, the test particle is circulating, and scattered into an eccentric orbit. This can be seen in Fig. 1. We perform many nu-  merical integrations with various b in the case of a moonlet of mass m = 3 × 10 −12 M , starting the particle at an azimuth |φ 0 | = 3000 r H /r m = 0.3 (where the moonlet is at φ = 0). This angle is large enough so that at this location the influence of the moonlet is negligible and the orbital parameters of the test particle are not disturbed, as will be checked later. We find that the horseshoe regime occurs forb < 1.8, and the scattered regime occurs forb > 2.5. For 1.774 <b < 2.503, however, the trajectory approaches the center of the moonlet to within a distance smaller than 0.95 r H . In that case, if one assumes the moonlet is a point mass, the test particle eventually leaves the Hill sphere, either on a horseshoe or a circulating trajectory, but the outcome changes several times with increasingb. In the case we are concerned about here, the moonlet most likely almost fills its Roche lobe, and therefore we stop the integration of the trajectory as soon as the distance between the test particle and the moonlet is less than 0.95 r H , assuming a collision. The specific orbital angular momentum J = r 2 (dφ/dt) of the test particles is computed along the trajectories. Angular momentum is exchanged during the close encounter with the moonlet. For b 2.503, the test particle is scattered onto an eccentric orbit of larger angular momentum than initially, which results in a gain in angular momentum. The variation of orbital angular momentum along the trajectory is shown in the bottom panel of Fig. 2 for the caseb = 3, where the top panel is the trajectory. The difference in angular momentum between the initial circular orbit at φ 0 = 0.3 sgn(b) and the end of the integration, when |φ| = 0.3 again, is noted ∆J. In the figure, only the interval −0.05 < φ < 0.05 is displayed, for convenience. Most of the exchange of angular momentum occurs when |φ| < 0.01. Figure 3 shows |∆J| (top thick curve) as a function ofb, in units of the specific angular momentum of the moonlet J m = r 2 m ω. For 0 <b 1.774, the horseshoe trajectory corresponds to a U-turn towards the central planet, and to a loss of angular momentum for the test particle. More precisely, as for circular orbits J ∝ r 1/2 , one expects for such a U-turn ∆J/J = 1 2 ∆r r = −b/r m ; this is indeed the case forb < 1.3. In the case where the test particle collides with the moonlet, we assume that it gives all its orbital angular momentum to the moon- . This also appears in Fig. 3. The opposite holds for b < 0.
Computing ∆J as a function of b to numerical precision enables us to also compute the difference between the inner and outer disk : . This quantity is small with respect to ∆J(b), but nonetheless well determined and converged in our simulations : is constant after the encounter to a precision better than 0.5% for all |φ| > 0.02. This validates our choice of φ 0 . In Fig. 3, the bottom thick dashed curve shows δJ in the same scale as |∆J|. We see that δJ > 0 for all b > 0 and that δJ ≪ ∆J, with for circulating trajectories, and for horseshoe orbits. In the following subsection, the empirically found Eq. (4) is derived analytically and justified.

Analytic model for the ring moonlet interaction
In this section, we consider only circulating trajectories. Developing to second order the exchange of angular momentum during an encounter with the moonlet ∆J, we can find the asymmetry δJ.

Solution for the perturbed moonlet orbit
Let us start again from Eq. (2) and (3). The ring particle is assumed to be on an unperturbed circular orbit of radius r 0 = r m + b. It orbits with angular velocity Ω = GM/r 0 3 such that φ = Ωt, where without loss of generality we have defined the origin of time t = 0 to be when the particle is at φ = 0. Under the perturbation induced by Ψ, the particle moves to r = r 0 + x, and φ = Ωt+y/r m , where x and y are assumed to be small. Linearizing Eq. (2) and (3) about the circular orbit state, we obtain equations for x and y in the form -Top panel : trajectory of a test particle with impact parameterb = 3 ; the motion of the particle is toward negative φ. Bottom panel : variation of the specific orbital angular momentum J of the same particle along its trajectory. Here the subscript 0 denotes evaluation on the unperturbed particle orbit.
We suppose that the perturbation is induced by a moonlet of mass m that is on a circular orbit of radius r m and has an angular velocity ω. Then its azimuthal coordinate φ m = ωt + φ m,0 , with φ m,0 being a constant. The perturbing potential Ψ = Gm becomes a function of time through substituting Using this in Eq. (6) and integrating with respect to time, we obtain which when combined with Eq. (5) gives an equation for x in the form

Solution of the linearized equations
To solve Eq. (10), we note that the perturbing potential Eq. (7) evaluated on the unperturbed orbits is a periodic function of time with period 2π/|ω−Ω| = 2π/β. Thus we should look for a periodic response. In order to do this we have to introduce a small frictional term into Eq. (10) to enable transients to decay and a net torque on the moonlet to be set up. When the frictional term is small it is expected that the resulting torque should not depend on it (e.g. Goldreich and Tremaine 1980). Hence we add a frictional term γ(dx/dt) to the left hand side of Eq. (10), where γ/Ω is a small constant parameter so that it now reads (11) As the potential is periodic in time we can adopt a Fourier series of the form where it is implicit that the real parts of such complex expressions is to be taken, and The periodic solution of Eq. (11) is now readily written down as We may write this in terms of a Green's function defined through Then the solution for x may be written Note that as the orbit of a ring particle relative to the planet is periodic, the solution given by Eq. (16) includes the effects of infinite numbers of repeating encounters. However, we wish to consider the case when dissipative effects, although weak, are strong enough to recircularize orbits between encounters in which case they will be independent of each other. This condition requires that γ/|ω−Ω| = γ/β ≫ 1. This is equivalent to requiring that the damping time scale be short compared to the relative orbital period between moonlet and ring particle. On account of the length scale of the encounters of interest being comparable to the Hill radius of the moonlet, this is much longer than the orbital period itself, so that we may adopt the ordering In order to make use of the above ordering we write down the form of the Green's function derived in the appendix (see Eq. (A7) ) valid for 0 < t < 2π/β.
where ω γ = Ω 2 − γ 2 /4. The function is defined elsewhere through its periodicity with period 2π/β. Making use of the inequality Eq. (17) we may replace the Green's function Eq. (18) by the simple expression Then the solution Eq. (16) gives

where
and (22) To make use of the above expressions, we consider the situation when the ring particle has a close encounter with the moonlet at time t = 0, thus we take φ m.0 = 0 (we note that a non zero φ m,0 can be dealt with by rotating the coordinate system and shifting the origin of time). The source term S is then expected to be highly peaked around t ′ = 0, and almost all of the contributions to the above integrals will occur for |t ′ | <∼ 2π/Ω. Furthermore, during this dynamical interaction, dissipation will be negligible. Thus, if we are interested in times after the main interaction, but before significant dissipation takes place, we may set γ = 0 and extend the limits of the integration to ±∞. However, in practice one may have to apply a cut off to the potential at large distances from the moonlet in order to do that (see below). But this should not matter if the important interaction occurs when the moonlet and ring particle are close.
Then we simply have and Thus A and B are constants representing epicyclic oscillation amplitudes induced after the close approach of the ring particle to the moonlet. We remark that the approximations made in obtaining Eq. (23) and Eq. (24) relate to how dissipation is treated. There has been no assumption that the particle trajectories are symmetric on opposite sides of the moonlet so that curvature effects remain fully incorporated during particle moonlet encounters. When dissipation is negligible during the encounter, then immediately afterward an epicyclic oscillation is established. The assumption that dissipation circularizes orbits between encounters implies that we should consider approaching ring particles to be on circular orbits. The above discussion indicates that errors associated with this assumption are exponentially small.

Angular momentum transfer
For the set up considered here, symmetry considerations imply that S(t) is an even function of time (see below), so that B = 0. The generation of the epicyclic oscillation is associated with an angular momentum transfer between the moonlet and particle. To find this we firstly note that where a f and e are the post encounter semi-major axis and eccentricity of the particle. We also note that the Jacobi constant implies that the change of the particle orbital energy and angular momentum are related by ∆E = GM [1/(2r 0 ) − 1/(2a f )] = ω ∆J. This can be used to eliminate a f in Eq. (25) after which ∆J may be found correct to second order in e ≡ A/r 0 with the result that ∆J = Ω 2 A 2 /[2(ω − Ω)]. This in turn may be simply determined after evaluating A. Note that ∆J < 0 for particles interior to the moonlet which have Ω > ω and conversely ∆J > 0 for particles orbiting exterior to the moonlet.

Development of the perturbing potential
We now consider We begin by recalling that In order to evaluate the Fourier transform as specified by Eq. (23), which was derived under the assumption that the interaction occurs only near closest approach, we must truncate the potential at large |t|. As the encounter takes place over a time ≪ 1/β, this can be achieved by replacing cos(βt) in Eq. (27) by 1 − β 2 t 2 /2. Note that a dimensionless estimate of the error involved is of order (β/ω) 2 ∼ (r H /r m ) 2 , where r H is the Hill radius of the moonlet. This is small enough that the leading order asymmetry in the angular momentum transferred to orbits with the same impact parameter on either side of the disk can be estimated.
As the first stage in evaluating the Fourier transform of S specified in Eq. (23) that gives the epicyclic amplitude, we evaluate dφ Ω .
(28) This can also be expressed as where ξ 0 = (Ω|r 0 − r m |)/(β √ r 0 r m ), and K j denotes the modified Bessel function of the second kind of order j.

Total angular momentum exchange
We may now use the above expression together with Eq. (26) to evaluate the epicyclic amplitude Eq. (23) (noting that the radial derivative is with respect to r 0 with other quantities held fixed) so obtaining The associated angular momentum exchanged is then given by In a strictly local approximation under which the inner and outer sides are symmetric, the contributions from orbits equidistant from the moonlet would cancel, leaving the net result to be determined by the surface density profile. However, although we have assumed the interactions are local, we did not assume symmetry between the exterior and interior orbits. Accordingly we evaluate the difference in the magnitude of ∆J evaluated from orbits equidistant from the moonlet : r 0 = r m ±b. The leading order contribution to ∆J is symmetric in b. The lowest order contribution is antisymmetric and accordingly leads to cancellation between the two sides. We make use of the where α = 3 4 + (6K 1 (2/3) + 3K 0 (2/3)) (4K 0 (2/3) + 2K 1 (2/3) = 2.46 . (33) The first order term of Eq. (32) was already given by Goldreich and Tremaine (1980). It is plotted as a straight green dashed line in Fig. 3. Our expansion to second order enables us to go further, and to give the expression of the magnitude of the asymmetry between the two sides of the disk : It is such that for an orbit with a given impact parameter, the angular momentum exchanged in the outer disk is the larger. In the case studied numerically, we had r H = 10 −4 , so that |b|/r m = 10 −4b . Then, Eq. (34) remarkably agrees with the numerical fit Eq. (4). Fig. 3 displays 4.92 × 10 −4b ∆J.

The light blue dot-dashed curve in
In the context of the above, we note that approximations made in obtaining equation Eq. (31) such as effectively starting and truncating the interaction at some finite though large distance from the moonlet could conceivably lead to changes comparable to those given by Eq. (34). However, such changes are again approximately symmetric for trajectories on both sides of the moonlet and thus approximately cancel so we do not expect such effects to significantly alter Eq. (34).

Migration rate and discussion
If the surface density of ring particles is Σ, the total rate of angular momentum transferred to the moonlet is where the integral is taken over the disk. The particles exterior to the moonlet contribute negatively while those interior contribute positively. The cumulative torque exerted by the moonlet on the region of the ring located within a distance b to its orbit reads then : The normalized cumulative torque is plotted in Fig. 4. The proportionality to Σ is obvious ; that T c ∝ (m/M ) 4/3 is numerically verified for 3 × 10 −15 m/M 3 × 10 −9 , and has been already found analytically by Ward (1991) for the horseshoe drag in a similar context.
Most of the total torque comes from scattered, circulating particles, in particular the ones with smallest impact parameterb ≈ 2.5. This makes the total torque sensitive to the physical size of the moonlet (taken as 0.95 r H here), as some particles colliding with the moonlet could be circulating if it were smaller.
The role of the horseshoe drag appears to be non negligible, amounting to ∼ 4.1(Σ/M r −2 m )(m/M ) 4/3 M r 2 m ω 2 . The expression of Ward (1991) for the torque arising from material executing horseshoe turns, called the horseshoe drag, is for a Keplerian disk with flat density profile : where w is the half-width of the horseshoe region. In our case, w = 1.774 r H , which gives T HS = 2.6 (Σ/M r −2 m )(m/M ) 4/3 M r 2 m ω 2 . The agreement is good because Ward's analysis is based only on geometrical effects and angular momentum variation in a Keplerian disk, without any pressure effect. Therefore, it also applies in Saturn's ring. We remark that taking w = 2r H in Eq. (37) gives a perfect match with what we find numerically for the total horseshoe drag.
In conclusion, from Fig. 4, the total torque felt by a moonlet of mass m on a circular orbit of radius r m around a planet of mass M can be written as Note that to get the same dependency of the type I torque in the parameters of the system, one has to assume h ∝ r H /r m in Eq. (1) ; however, in a protoplanetary disk, h is fixed and independent of the mass of the secondary body, so that this proportionality would not be justified. The torque is related to the migration speed through T = 0.5 m r m Ω(dr m /dt). Hence we deduce that The migration rate is here proportional to the mass of the moonlet to the power 1/3, in contrast to standard type I migration where dr m /dt ∝ m. A numerical application to the case of an m = 10 −18 M Saturn = 5.68 × 10 8 kg moonlet in orbit in the A-ring of density Σ = 400 kg m −2 at r m = 130 000 km from Saturn gives dr m /dt = −0.23 m yr −1 . Increasing the mass by two orders of magnitude to correspond to a radius of ∼ 200 m speeds up the migration rate by a factor of only ∼ 4.5 to ∼ −1 m yr −1 .
After time t, a migrating propeller will be shifted longitudinally with respect to a corresponding non migrating one by a distance r m δφ = 3Ω |dr m /dt| t 2 /4. For the above parameters, this gives 713 [t/(1 year)] 2 m. A shift of this magnitude is potentially detectable on a timescale of a year to a few years (Porco et al. 2004 and note also Burns et al. 2009). Actually, migration of propellers has already been detected (Burns et al. 2009;Tiscareno et al. 2010). During one time period of nearly a year, a particular propeller has been seen moving outward at a rate of ∼ 110 m yr −1 ; and during a later similar time period, the same propeller has been seen moving inward at a rate of ∼ 40 m yr −1 (Tiscareno et al. 2010, and personal communication).
These observations are not compatible with the above theory. But we recall that the process of migration of a moonlet described above assumed a smooth particle disk with constant surface density. Here we note that there are features and mechanisms that might produce a significantly faster migration rate, possibly in both directions inward and outward, and non constant in time. One can first think of a radial density gradient : as there is no pressure buffer here, this would directly affect the balance between the torques from the inner and outer parts of the ring. This would also affect the torque from the horseshoe region, which could turn positive. However, if the migration is governed by the gradient of some quantity, it seems likely that the moonlet would have approached an extremum in that quantity, and thus should have attained a migration rate comparable to that estimated for a constant surface density.
Another possibility resulting in the moonlet migrating faster than what the previous calculation indicates, and possibly outward, is a runaway migration in a planetesimals disk (Ida et al. 2000;Levison et al. 2007, for a review), similar to the type III migration of planets in protoplanetary disks (Masset and Papaloizou 2003). In this regime, the migration of the moonlet in the disk leads to a positive feedback on its migration rate, because of the material of the inner (resp. outer) disk making horseshoe U-turns to the outer (resp. inner) disk. This speeds up the migration, possibly leading to a runaway. However, this leads inevitably to an asymmetry in the horseshoe region, while the propeller structures observed are rather symmetrical.
Finally, the A-ring of Saturn is not homogeneous. It is close to gravitational instability, which should lead to the formation of gravity wakes and density fluctuations. The effect of these density fluctuations on the moonlet is studied in next section.

The role of density fluctuations and resulting stochastic migration
The analytic calculations and numerical simulations in the previous chapters assume an inflow of particles on circular orbits only perturbed by the nearby moonlet. However, we know that Saturn's A ring is marginally gravitationally stable (Daisaka et al. 2001). The Toomre Q parameter (Toomre 1964), which is a measure of the importance of self-gravity, is expected to be of the order of 2 ∼ 7, indicating that the ring particles' mutual gravity is indeed a strong effect. It leads to the regular formation and dispersion of gravity wakes, which are local density enhancements elongated in parallel directions by the Keplerian shear. Those over-densities give rise to stochastic forces which act on the embedded moonlet.
A very similar effect is expected to occur in protoplanetary disks. These disks are thought to be turbulent due to the magneto-rotational instability (MRI, Balbus and Hawley 1991). The turbulent fluctuations create over-densities which interact gravitationally with embedded small mass planets. The stochastic forces make the planet undergo a random walk. An analytic model of this random walk has been derived by Rein and Papaloizou (2009). In the following, we apply this model to moonlets embedded in Saturn's rings. To do that, we need to get an estimate of the amplitude of the stochastic forces.

Numerical calculations
We perform three-dimensional simulations of ring particles, in a shearing box, similarly to Salo (1995). The simulations are done in a local cube of size H with shear periodic boundary conditions, and the origin of the box is fixed at a semi major axis of a = 130 000 km. A BH tree code (Barnes and Hut 1986) is used to calculate the self-gravity between ring particles and resolve inelastic collisions. Collisions between particles are resolved using the instantaneous collision model and a velocity dependent coefficient of restitution given by Bridges et al. (1984) : where v is the impact speed projected on the vector joining the centers of the two particles. The code is described in more detail in .
The size of ring particles is not well constrained. Therefore and to be able to scale to different locations in Saturn's rings, we perform multiple simulations. For a given simulation, all the particles are assumed to be spherical and have the same size (or radius), which varies from simulation to simulation from 0.52 to 13 meters. The simulation parameters are listed in Table 1. The nomenclature and physical parameters are, for easy comparison, the same as in Lewis and Stewart (2009), as our simulations are similar to theirs.
The moonlet is not taken into account in the simulations. We measure the specific gravitational forcef (or acceleration) felt by a passive testparticle sitting at the origin. We calculate the force in two different ways, in order to avoid the singularity at the origin and to account for the physical size of the moonlet. In the first case, we use a cut off at the moonlet's radius d and exclude all particles within that radius from the force calculation. In the second case we use a smoothed gravitational force per unit mass in the form wherer is the vector linking the origin to the particle and m part is the mass of the particle. The smoothing length d is set equal to the moonlet's size. In a self-consistent simulation, one should include the moonlet with it's real physical size. However, this goes beyond the scope of this paper and will be considered in future work (Rein and Papaloizou 2010). Our purpose here is to estimate the underlying stochastic fluctuations in the migration rate that occur independently of the moonlet. This procedure is reasonable as long as the moonlet is in a steady state, namely if it doesn't accumulate or lose a large amount of mass over one orbit. In all our simulations we assume a moonlet size of d = 200 m.  Table 1: Simulation parameters. The first column identifies the simulation, following the convention of Lewis and Stewart (2009). The second and third column give the size (or radius) of the particles and their density. The fourth and fifth column give the surface density and the size of the computational domain, respectively. The last column lists the number of particles.

Results
We measure the amplitude and the correlation time of the stochastic forces in all simulations. The results are listed in Table 2. We also plot the time evolution of the azimuthal (y) force component in Fig. 5. Whereas the correlation time in all simulations is almost the same, the amplitude of the stochastic fluctuations varies by almost a factor of 10 3 . The forces in the vertical direction are negligible and not presented here. The diffusion coefficient, being a measure of the strength of stochastic forces, is defined as D = 2 f 2 τ (see Rein and Papaloizou 2009), where f 2 1/2 and τ are the root mean square value and the approximate correlation time of the specific stochastic force in one direction.
The change in semi major axis a due to the effect of stochastic forces with diffusion coefficient D after time t is given by where ω is the mean motion of the moonlet (Rein and Papaloizou 2009). Note that in this regime, the acceleration of the moonlet doesn't depend on its mass (as can be seen in Eq. (41)). Therefore, the migration rate and the diffusion coefficient are independent of the mass of the moonlet ; this might be an observational indication for this migration regime.
Assuming an initial semi-major axis of a = 130 000 km and D ∼ 10 −17 m 2 s −3 as found in simulation S1, one can calculate the expected difference in semi major axis after one orbit due to stochastic forces which turns out to be ∆a = 0.01m. For simulation S3, assuming D ∼ 10 −15 m 2 s −3 , one finds ∆a = 0.11m. For the simulation L2, assuming D ∼ 10 −11 m 2 s −3 , one finds ∆a = 10.5 m. These translate to random walks with standard deviation given as a function of time by ∆a = 0.27 t/(1 year) m, ∆a = 2.7 t/(1 year) m, and ∆a = 270 t/(1 year) m respectively.

Discussion
From the above, it can be seen that increasing the surface density by factors 4 − 5 changes the migration rate by two orders of magnitude. Thus, the results show clearly that the surface density Σ is much more important than in the regular, type I like migration, presented in Sect. 3. This can be easily understood with a toy model. The critical unstable wavelength λ scales linearly with Σ (Toomre 1964). If we assume a fixed moonlet size, the ratio of moonlet size to λ therefore changes with Σ. In the limit where λ is much smaller than the moonlet radius, the stochastic forces are negligible as the density distribution is approximately homogeneous on the relevant scales. This is the case in simulation S1. In the other limit where λ is larger than the moonlet, the moonlet undergoes a random walk that is similar to that of individual ring particles, as seen in simulation L2.
The range in migration rates found in simulations shows that over time scales of several years, the migration of a moonlet of mass m ∼ 10 −16 M Saturn may be dominated by a random walk in some situations (eg. those in simulations S3 and L2, the latter carried out with particles of radius 13 m). However, in regions of the rings where the surface density is small (eg. simulation S1), the moonlet may be in a regular, nonstochastic migration regime. In that case, the model from Sect. 3 can be applied.
This dependence offers an exciting possibility to constrain the nature of the ring particles and the physical processes occurring in the rings by  Table 2: Simulation results. The first column gives the name of the simulation, as defined in Table 1. The second and third columns give the correlation time in the x (radial) and y (azimuthal) direction, respectively. The fourth and fifth columns list the diffusion coefficients. The sixth column is the Toomre Q parameter, as measured in the simulation.
measuring the migration of moonlets. But note that because regular migration gives a decrease in the semi-major axis that is linear in time, provided it continues to operate, it will always ultimately dominate the behavior for large time because the spreading of the semi-major axis associated with stochastic migration increases only as the square root of time.
At the present day, the number of observed migration rates is not sufficient to draw a statistically significant conclusion. However, the fact that the migration varies in rate and direction clearly favors the stochastic migration model presented in this section. Considering a migration rate of |∆a| = 100 m in t = 1 year in Eq. (42), one finds D = 1.4 × 10 −12 m 2 s −3 . Taking |∆a| = 40 m in t = 1 year gives D = 2.2 × 10 −13 m 2 s −3 . This is in the range obtained in the simulations, and tends to favor the case of simulation L2, and Σ ∼ 700 kg m −2 in the A ring.

Conclusion
In this paper we have calculated the differential torque exerted on a moonlet by the outer and the inner disk with a smooth, flat surface density profile. We performed both an accurate numerical integration and a second order analytical calculation. These approaches were found to be in excellent agreement where their domains of validity overlap. The migration rate found in this case is proportional to the mass of the moonlet to the power 1/3. It is about −1 m yr −1 for a 200 m radius moonlet in the A-ring. This is way too low to explain the observed migration of the propellers in Saturn's rings. Nonetheless, density fluctuations in the rings, due to their proximity to gravitational instability, can lead to stochastic torques on a moonlet, that may dominate on the timescale of the Cassini mission, to an extent that depends mainly on the surface density of the rings. These stochastic torques may account for the observations.
The possibility that the migration of propellers is induced by stochastic processes rather than by a regular type I like migration is therefore very exciting : this may help to infer the local surface density, and therefore the size of the ring particles. Indeed, both quantities are linked through the optical depth, which is observationally well constrained. Our estimate of Σ ∼ 700 kg m −2 is equivalent to ∼ 10 m size particles in the A ring for a single-sized population. This is in good agreement with available estimates from stellar occultations. For the A ring, Voyager occultations (Zebker et al. 1985) find that the radius r a of particles follows : 0.1 m (assumed) < r a < 11 m, with a power index ∼ −3. For the 28 Sgr occultation (French and Nicholson 2000), a range 1m < r a < 20m is found, with a power index between −2.7 and −3, and an effective size (the single average size accounting for the fluctuations in photon count) being about 7m (Cuzzi et al. 2009).
Fortunately the Cassini mission has been extended to 2017. In the meantime, numerous observations of the propellers will hopefully give a clear picture of their orbital evolution for a period of 10 years, representing about 10 000 orbits. This will allow us to test the hypothesis presented in this paper, and in particular to check whether the propellers really are in stochastic migration, whereas first results seem to favor the stochastic hypothesis.
Among the questions that still need to be addressed is why all the propellers seem to be gathered in the A ring, in places apparently devoid of density waves ? Indeed the propellers seem gathered in a couple of narrow radial bands of only about 1000 km width (Tiscareno et al. 2008). This is especially surprising since the A ring is densely populated by numerous density waves launched by the nearby small moons (Atlas, Prometheus, Pandora, Janus, Epimetheus). Is there a systematic mechanism that would eject the propellers away from density waves ? Or does the present location of propellers just reflect the initial location of the parent body, assuming that the propeller population comprises the fragments resulting from the destruction of an ancient moon orbiting within the rings ?
If the moonlets really undergo stochastic migration, then Eq. (42) may strongly constrain the age of the propellers, which can't be larger than the time needed to diffuse over ∆a > 1000 km. Unfortunately, as long as D is unknown Eq. (42) doesn't provide any useful information. However, considering that ∆a is proportional to the square root of the time, and assuming that a moonlet migrates about 100 m in 1 year, one finds that ∆a = 1000 km for t = 100 million years. Assuming |∆a| = 40 m in 1 year, we find that it requires ∼ 625 million years to diffuse over 1000 km. Note that for 200 m radius moonlets, the spreading due to stochastic migration equates to the contraction of the semi-major axes occurring as a result of smooth, regular migration after ∼ 10 4 yr, and then ∆a ≈ 10 km. Thus it would take about a million years to migrate through ∆a = 1000 km in this case, largely through the action of the nonstochastic, regular migration process, if that can be assumed to operate smoothly and simultaneously with the stochastic migration process. A 100 m radius moonlet would migrate through only 500 km in the same period, so that the smooth, regular migration process spreads a population of moonlets of various sizes over 1000 km in one to two million years. These could be the times since the catastrophic disruption of a small moon orbiting at 130 000 km from Saturn, that was broken into smaller moonlets by a meteoritic impact. On the other hand, an estimate of the lifetime of a Pan size moon (∼ 14 km in radius) against the today's cometary flux is provided by Dones et al. (2009) and gives a range between 100 Myr and 16 Gyr, depending on the size distribution of impactors. Therefore, the recent occurrence of such an event, about 4 to 4.5 billion years after solar system formation, is possible. Note also that an age of about 100 Myr is coherent with some estimates of Saturn's ring age despite of the lack of fully satisfactory explanation for their origin (see Charnoz et al. 2009a, for a review).
We see that the question of propeller's migration is inextricably linked to the issue of the origin of Saturn's moons embedded in the rings, which is still a mystery. Porco et al. (2007) and Charnoz et al. (2007) have jointly proposed that small moons embedded in the rings could be aggregates of material on an initial shard denser than ice. When destroyed by meteoritic bombardment, these could release dense chunks of material that could explain the origin of the propellers. However, the origin of Saturn's ring system is still a matter of debate (Harris 1984;Charnoz et al. 2009b). Knowledge of the age of the propellers could provide important constraints on the age of the main ring system and its embedded moons, as there are strong indications that these could have about the same age, provided these moonlets hide a dense shard (Charnoz et al. 2007;Porco et al. 2007). Understanding the migration rate of the propellers is therefore an important piece of this puzzle.
To perform the summation we use the general result that if for a general periodic function with period 2π/β, and b(n) being defined as an integrable function, we set then g(τ ) = 2π n=∞ n=−∞ F (τ + 2πn/β).