N -body Models of Aggregation and Disruption

1. N-body Models of Aggregation and Disruption Derek C. Richardson University of Maryland

2. Overview • Introduction/the N-body problem. • Numerical method (pkdgrav). • Application: binary asteroids. • Non-idealized & strength models. • First results: “YORP” spinup of rubble piles & spin limits with strength.

3. Introduction • Many dynamical processes in the solar system can be modeled by gravity and collisions alone. E.g., • Reaccumulation after catastrophic disruption (collisional or rotational). • Planetary ring dynamics. • Planet formation. • Problems well suited to N-body code.

4. The N-body problem The orbit of any one planet depends on the combined motion of all the planets, not to mention the actions of all these on each other. To consider simultaneously all these causes of motion and to define these motions by exact laws allowing of convenient calculation exceeds, unless I am mistaken, the forces of the entire human intellect. — Isaac Newton, 1687.

5. The N-body problem Cost = N (N – 1) / 2 = O(N2)

6. Tree codes • Reduce computational cost by treating particles in groups.

7. Tree codes Replace many summations with single multipole expansion around center of mass.

8. Tree codes • Reduce computational cost by treating particles in groups. • Error controlled by opening angle criterion and order of expansion.

9. Tree codes Use multipole expansion if opening angle  < crit.  crit

10. Tree codes • Reduce computational cost by treating particles in groups. • Error controlled by opening angle criterion and order of expansion. • Particles organized into systematic hierarchical structure. • Ideally suited for recursive algorithms.

11. Tree codes E.g. Barnes & Hut (1986) two-dimensional tree. Cost = O(N log N)

12. Reducing cost further • Parallel methods: • Distribute work among Np processors. • N-body problem difficult—exploit tree. • Adaptive/hierarchical timestepping: • Focus work on most active particles. • Good object-oriented code structure. • Hard-core optimizations.

13. Integrating the equations of motion • Many techniques for solving coupled linear ordinary differential equations. • Most popular: • Runge-Kutta (explicit forward). • Bulirsch-Stoer (complex/expensive). • Leapfrog/symplectic methods. • Preserve phase space volume. • Timestep adaptability issues.

14. R1 R2 Collision detection • Particles collide when separation distance equals sum of radii.

15. Collision detection • Particles collide when separation distance equals sum of radii. • Two approaches: • Predict collisions before they occur. • Need neighbour-finding algorithm (tree!). • Detect collisions after they occur. • Detected by mutual overlap. • Adaptive timestepping essential.

16. Numerical method • Our group uses pkdgrav: • Parallel k-D tree code. • k-D: split along longest dimension. • Expand to hexadecapole order. • Second-order leapfrog integrator. • Hierarchical timestepping. • Collisions predicted before they occur. • Includes bouncing and sliding friction.

17. Parallelism in pkdgrav • master • controls overall flow • “mdl” • interface between pkdgrav and parallel primitives (e.g. mpi) • “pst” • loops over processors • “pkd” • loops over particles on one processor

18. Application: binary asteroids • Use N-body code to simulate: • Capture of collisional ejecta in Main Belt. • Michel et al., Durda et al.: collisions that make families also make satellites.

19. Application: binary asteroids Michel et al. 2001

20. Application: binary asteroids • Use N-body code to simulate: • Capture of collisional ejecta in Main Belt. • Michel et al., Durda et al.: collisions that make families also make satellites. • Rotational disruption of gravitational aggregates in near-Earth population. • Tidal disruption. • “YORP” thermal spin-up.

21. Application: binary asteroids

22. Tidal disruption vs. YORP • Tidal disruption makes binaries, but also destroys them quickly. • Binary NEA mean lifetime only ~ 1 Myr. • YORP thermal effect may form binaries through rotational disruption. • But, some internal strength/cohesion may be necessary to prevent material from just “dribbling” away (but that may be OK too!).

23. Forming binaries with YORP • Preliminary investigation: • Slowly spin up various rubble piles. • Find particles leak away from equator (no fission).

24. Forming binaries with YORP • Preliminary investigation: • Slowly spin up various rubble piles. • Find particles leak away from equator (no fission). Recoil: new mobility mechanism?

25. Forming binaries with YORP

26. Forming binaries with YORP • May need strength and/or irregular body shape to form binaries. • E.g., contact binary can separate.

27. Non-idealized models • Treating particles as idealized, rigid, independent spheres is convenient. • Components with different shapes may provide more realism. E.g., • Ellisoidal particles (Roig et al.) • Polyhedral (Korycansky & Asphaug). • We combine best of both worlds: allow spheres to “fuse” together…

28. Non-idealized models

29. Strength model • Colliding particles/aggregates can: • Stick on contact; • Bounce; • Liberate particle(s) from aggregate(s). • Outcome currently parameterized by impact speed.

30. Strength model • In addition, bonded aggregates can have a size-dependent bulk tensile and/or shear strength. • Particles experiencing stress (relative to center of mass) in excess of strength are liberated. • Global model (no fractures/cracks).

31. Strength model For a demo of the new strength model in action, see Patrick’s presentation!

32. Testing strength: spin limits • One way to test the strength model is to compare with analytical predictions of global failure (e.g. Holsapple). • Found good match for cohesionless models (Richardson et al. 2005). • Science motivation: spin-up past critical limit could make binaries (e.g. YORP).

33. Spin limits: preliminary results Work in progress!

38. Summary • N-body methods allow modeling of complex phenomena involving gravity & collisions. • Examples include post-disruption gravitational reaccumulation to form binaries & families. • Binaries: more work needed to assess YORP (including survivability against BYORP!). • New pkdgrav strength model provides added realism/complexity, but needs fracture model.

39. Extra Slides…

40. What is YORP? • Yarkovsky-O'Keefe-Radzievskii-Paddack effect. • Irregular bodies reflect/re-radiate solar photons in different directions: net torque  spin-up/down.

41. Results: Many Binaries Close approach distance q Encounter speed v∞ • High rates of production for: • Low q. • Low v∞. • Rapid spin. • Large elongation. Spin period P Elongation ε

42. Retrograde Orbital Properties Semimajor axis a (50% > 10 Rp) Eccentricity e (97% > 0.1) • High eccentricity. • Range of semi-major axis. • Binary orbit aligned more with approach orbit than progenitor spin. • Retrograde orbits possible. Spin-orbit angle Inclination I

43. Physical Properties Size ratio • Size ratio peaks at 0.1–0.2 (10–5:1). • Obliquities: • Primary spin aligned with binary orbit. • Wide range of secondary spin axes. • Spin Periods: • Primary has narrow range (3.5  6.0 h). • Secondary has wide range (4.0  20+ h). Obliquities Spin periods

44. Evolutionary Effects • Mutual tides damp eccentricity in ~ 1–10 My. • Repeated encounters may strip binary. • NEA population refreshed by MBAs (some of which may be binary). • Thermal effects (YORP) important?

45. Steady-state (Monte Carlo) Model • We know… • Binary production efficiency from tidal disruption (Walsh & Richardson 2006); • Planetary encounter circumstances (Bottke et al. 1994); • Distribution of NEA lifetimes (Gladman et al. 2000); • Shape and spin of source bodies (Harris et al. 2005); • Tidal evolution effects (Weidenschilling et al. 1989); • Effects of binary encounters with planets (Bottke & Melosh 1996; this work); • Small binary MBAs formed in collisional simulations (Durda et al. 2004).

46. Steady-state (Monte Carlo) Model • In one timestep… • All asteroids in the simulation are tested for: • End of lifetime (median ~ 10 Myr); • Close planetary encounter < 3REarth (one every ~3 Myr). • All binaries are tested for: • End of lifetime; • Close planetary encounter < 24REarth: explicit 3-body integration. • If neither happen, the binary is tidally evolved. • Removed NEAs/binaries are immediately replaced. • “Fresh” asteroids take spin/shape characteristics of MBAs, with a variable percentage being binaries. • MBA binaries have characteristics determined from the Durda et al. 2004 simulations.

47. Steady-state Results For 2000 asteroids: • Find ~2% binary fraction. • Binary NEA mean lifetime ~ 1 Myr. • 93% of removed binaries destroyed by planetary encounters. • MBA initial binary percentage has little effect (mean lifetime ~0.32 Myr).

48. Steady-state Results • The resultant steady-state binaries… • Have slightly larger semi-major axes than observed; Observed Steady-state

49. Steady-state Results • The resultant steady-state binaries… • Have slightly larger semi-major axes than observed; • Mostly have low eccentricities (< 0.2), consistent with observations. Eccentricity

50. A Word About Rubble Piles • Rubble piles are low-tensile-strength, medium-porosity gravitational aggregates. • In simulations, rubble piles consist of perfectly smooth spheres; some dissipation. • Used in a variety of contexts: planetesimal collisions, tidal disruption, spin-up. • How do they differ from perfect fluids?