We present a grid-based non-parametric approach to obtain a triaxial three-dimensional luminosity density from its surface brightness distribution. Triaxial deprojection is highly degenerate and our approach illustrates the resulting difficulties. Fortunately, for massive elliptical galaxies, many deprojections for a particular line of sight can be discarded, because their projection along other lines of sight does not resemble elliptical galaxies. The near-elliptical isophotes of these objects imply near ellipsoidal intrinsic shapes. In fact, deprojection is unique for densities distributed on ellipsoidal shells. The constrained non-parametric deprojection method we present here relaxes this constraint and assumes that the contours of the luminosity density are boxy/discy ellipsoids with radially varying axis ratios. With this approach we are able to reconstruct the intrinsic triaxial densities of our test models, including one drawn from an $N$-body simulation. The method also allows to compare the relative likelihood of deprojections at different viewing angles. We show that the viewing orientations of individual galaxies with nearly ellipsoidal isophotal shapes can be constrained from photometric data alone.

The coalescence of massive black hole binaries (BHBs) in galactic mergers is the primary source of gravitational waves (GWs) at low frequencies. Current estimates of GW detection rates for the Laser Interferometer Space Antenna and the Pulsar Timing Array vary by three orders of magnitude. To understand this variation, we simulate the merger of equal-mass, eccentric, galaxy pairs with central massive black holes and shallow inner density cusps. We model the formation and hardening of a central BHB using the Fast Multiple Method as a force solver, which features a $O(N)$ scaling with the number $N$ of particles and obtains results equivalent to direct-summation simulations. At $N \sim 5\times 10^5$, typical for contemporary studies, the eccentricity of the BHBs can vary significantly for different random realisations of the same initial condition, resulting in a substantial variation of the merger timescale. This scatter owes to the stochasticity of stellar encounters with the BHB and decreases with increasing $N$. We estimate that $N \sim 10^7$ within the stellar half-light radius suffices to reduce the scatter in the merger timescale to $\sim 10$\%. Our results suggest that at least some of the uncertainty in low-frequency GW rates owes to insufficient numerical resolution.

We simulate the passage through the Sun-Jupiter system of interstellar objects (ISOs) similar to 1I/`Oumuamua or 2I/Borisov. Capture of such objects is rare and overwhelmingly from low incoming speeds onto orbits akin to those of known long-period comets. This suggests that some of these comets could be of extra-solar origin, in particular inactive ones. Assuming ISOs follow the local stellar velocity distribution, we infer a volume capture rate of $0.051\,\mathrm{au}^3 \mathrm{yr}^{-1}$. Current estimates for orbital lifetimes and space densities then imply steady-state captured populations of $\sim10^2$ comets and $\sim10^5$ `Oumuamua-like rocks, of which 0.033% are within 6au at any time.

We report the discovery of a Milky Way satellite in the constellation of Antlia. The Antlia 2 dwarf galaxy is located behind the Galactic disc at a latitude of $b\sim 11^{\circ}$ and spans $1.26$ degrees, which corresponds to $\sim2.9$ kpc at its distance of 130 kpc. While similar in spatial extent to the Large Magellanic Cloud, Antlia 2 is orders of magnitude fainter at $M_V=-9$ mag, making it by far the lowest surface brightness system known (at $\sim31.9$ mag/arcsec$^2$), $\sim100$ times more diffuse than the so-called ultra diffuse galaxies. The satellite was identified using a combination of astrometry, photometry and variability data from \textit{Gaia} Data Release 2, and its nature confirmed with deep archival DECam imaging, which revealed a conspicuous BHB signal. We have also obtained follow-up spectroscopy using AAOmega on the AAT, identifying 159 member stars, and we used them to measure the dwarf's systemic velocity, $290.9\pm0.5$km/s, its velocity dispersion, $5.7\pm1.1$ km/s, and mean metallicity, [Fe/H]$=-1.4$. From these properties we conclude that Antlia 2 inhabits one of the least dense Dark Matter (DM) halos probed to date. Dynamical modelling and tidal-disruption simulations suggest that a combination of a cored DM profile and strong tidal stripping may explain the observed properties of this satellite. The origin of this core may be consistent with aggressive feedback, or may even require alternatives to cold dark matter (such as ultra-light bosons).

We perform $N$-body simulations of the early phases of open cluster evolution including a large population of planetesimals, initially arranged in Kuiper-belt like discs around each star. Using a new, 4th-order and time-reversible $N$-body code on Graphics Processing Units (GPUs), we evolve the whole system under the stellar gravity, i.e. treating planetesimals as test particles, and consider two types of initial cluster models, similar to IC348 and the Hyades, respectively. In both cases, planetesimals can be dynamically excited, transferred between stars or liberated to become free-floating (such as A/2017 U1 or 'Oumuamua) during the early cluster evolution. We find that planetesimals captured from another star are not necessarily dynamically distinct from those native to a star. After an encounter both native and captured planetesimals can exhibit aligned periastrons, qualitatively similar to that seen in the Solar system and commonly thought to be the signature of Planet 9. We discuss the implications of our results for both our Solar system and exoplanetary systems.

Gas clumps formed within massive gravitationally unstable circumstellar discs are potential seeds of gas giant planets, brown dwarfs and companion stars. Simulations show that competition between three processes -- migration, gas accretion and tidal disruption -- establishes what grows from a given seed. Here we investigate the robustness of numerical modelling of clump migration and accretion with the codes PHANTOM, GADGET, SPHINX, SEREN, GIZMO-MFM, SPHNG and FARGO. The test problem comprises a clump embedded in a massive disc at an initial separation of 120 AU. There is a general qualitative agreement between the codes, but the quantitative agreement in the planet migration rate ranges from $\sim 10$% to $\sim 50$%, depending on the numerical setup. We find that the artificial viscosity treatment and the sink particle prescription may account for much of the differences between the codes. In order to understand the wider implications of our work, we also attempt to reproduce the planet evolution tracks from our hydrodynamical simulations with prescriptions from three previous population synthesis studies. We find that the disagreement amongst the population synthesis models is far greater than that between our hydrodynamical simulations. The results of our code comparison project are therefore encouraging in that uncertainties in the given problem are probably dominated by the physics not yet included in the codes rather than by how hydrodynamics is modelled in them.

Linear perturbation is used to investigate the effect of gravitational softening on the retrieved two-armed spiral eigenmodes of razor-thin stellar discs. We explore four softening kernels with different degrees of gravity bias, and with/without compact support (compact in the sense that they yield exactly Newtonian forces outside the softening kernel). These kernels are applied to two disc galaxy models with well-known unsoftened unstable modes. We illustrate quantitatively the importance of a vanishing linear gravity bias to yield accurate frequency estimates of the unstable modes. As such, Plummer softening, while very popular amongst simulators, performs poorly in our tests. The best results, with excellent agreement between the softened and unsoftened mode properties, are obtained with softening kernels that have a reduced gravity bias, obtained by compensating for the sub-Newtonian forces at small interparticle distances with slightly super-Newtonian forces at radii near the softening length. We present examples of such kernels that, moreover, are analytically simple and computationally cheap. Finally, these results light the way to the construction of softening methods with even smaller gravity bias, although at the price of increasingly complex kernels.

Tidal debris from Galactic satellites generally forms one-dimensional elongated streams, since nearby Galactic orbits have almost identical frequency ratios. We show that the situation is different for orbits close to the Galactic disc, whose vertical frequency $\Omega_z$ is strongly amplitude dependent. As a consequence, stars stripped from a satellite obtain a range of values for $\Omega_z$ and hence of frequency ratios, and spread into two dimensions, forming a ribbon-like structure with vertical extent comparable to that of the progenitor orbit. In integrals-of-motion space, tidal ribbons are clumps, which offers the best chance of detection and allows the determination of the Galactic potential vertically across the disc.

We derive unbiased distance estimates for the Gaia-TGAS dataset by correcting for the bias due to the distance dependence of the selection function, which we measure directly from the data. From these distances and proper motions, we estimate the vertical and azimuthal velocities, $W$ and $V_\phi$, and angular momentum $L_z$ for stars in the Galactic centre and anti-centre directions. The resulting mean vertical motion $\overline{W}$ shows a linear increase with both $V_\phi$ and $L_z$ at $10 \sigma$ significance. Such a trend is expected from and consistent with the known Galactic warp. This signal extends to stars with guiding centre radii $R_g<R_0$, placing the onset of the warp at $R\lesssim7{\rm kpc}$. At equally high significance, we detect a previously unknown wave-like pattern of $\overline{W}$ over guiding centre $R_g$ with amplitude $\sim1{\rm kms}^{-1}$ and wavelength $\sim2.5{\rm kpc}$. This pattern is present in both the centre and anti-centre directions, consistent with a winding (corrugated) warp or bending wave, likely related to known features in the outer disc (TriAnd and Monoceros over-densities), and may be caused by the interaction with the Sgr dwarf galaxy $\sim1{\rm Gyr}$ ago. The only significant deviation from this simple fit is a stream-like feature near $R_g\sim9{\rm kpc}$ ($|L_z|\sim2150{\rm kpckms}^{-1}$).

The processes driving gas accretion on to supermassive black holes (SMBHs) are still poorly understood. Angular momentum conservation prevents gas within $\sim10\,$pc of the black hole from reaching radii $\sim10^{-3}\,$pc where viscous accretion becomes efficient. Here we present simulations of the collapse of a clumpy shell of swept-up isothermal gas, which is assumed to have formed as a result of feedback from a previous episode of AGN activity. The gas falls towards the SMBH forming clumps and streams, which intersect, collide, and often form a disc. These collisions promote partial cancellations of angular momenta, resulting in further infall and more collisions. This continued collisional cascade generates a tail of gas with sufficiently small angular momenta and provides a viable route for gas inflow to sub-parsec scales. The efficiency of this process hardly depends on details, such as gas temperature, initial virial ratio and power spectrum of the gas distribution, as long as it is not strongly rotating. Adding star formation to this picture might explain the near-simultaneous formation of the S-stars (from tidally disrupted binaries formed in plunging gas streams) and the sub-parsec young stellar disc around Sgr A$^{\!\star}$.

Dwarf spheroidal galaxies of the Local Group obey a relationship between the line-of-sight velocity dispersion and half-light radius, although there are a number of dwarfs that lie beneath this relation with suppressed velocity dispersion. The most discrepant of these (in the Milky Way) is the `feeble giant' Crater II. Using analytic arguments supported by controlled numerical simulations of tidally-stripped flattened two-component dwarf galaxies, we investigate interpretations of Crater II within standard galaxy formation theory. Heavy tidal disruption is necessary to explain the velocity-dispersion suppression which is plausible if the proper motion of Crater II is $(\mu_{\alpha*},\mu_\delta)=(-0.21\pm0.09,-0.24\pm0.09)\mathrm{mas\,yr}^{-1}$. Furthermore, we demonstrate that the velocity dispersion of tidally-disrupted systems is solely a function of the total mass loss even for weakly-embedded and flattened systems. The half-light radius evolution depends more sensitively on orbital phase and the properties of the dark matter profile. The half-light radius of weakly-embedded cusped systems rapidly decreases producing some tension with the Crater II observations. This tension is alleviated by cored dark matter profiles, in which the half-light radius can grow after tidal disruption. The evolution of flattened galaxies is characterised by two competing effects: tidal shocking makes the central regions rounder whilst tidal distortion produces a prolate tidally-locked outer envelope. After $\sim70$ per cent of the central mass is lost, tidal distortion becomes the dominant effect and the shape of the central regions of the galaxy tends to a universal prolate shape irrespective of the initial shape.

Computational efficiency demands discretised, hierarchically organised, and individually adaptive time-step sizes (known as the block-step scheme) for the time integration of N-body models. However, most existing N-body codes adapt individual step sizes in a way that violates time symmetry (and symplecticity), resulting in artificial secular dissipation (and often secular growth of energy errors). Using single-orbit integrations, I investigate various possibilities to reduce or eliminate irreversibility from the time stepping scheme. Significant improvements over the standard approach are possible at little extra effort. However, in order to reduce irreversible step-size changes to negligible amounts, such as suitable for long-term integrations of planetary systems, more computational effort is needed, while exact time reversibility appears elusive for discretised individual step sizes.

The symplectic Wisdom-Holman map revolutionized long-term integrations of planetary systems. There is freedom in such methods of how to split the Hamiltonian and which coordinate system to employ, and several options have been proposed in the literature. These choices lead to different integration errors, which we study analytically and numerically. The Wisdom-Holman method in Jacobi coordinates and the method of Hernandez, H16, compare favorably and avoid problems of some of the other maps, such as incorrect center-of-mass position or truncation errors even in the one-planet case. We use H16 to compute the evolution of Pluto's orbital elements over 500 million years in a new calculation.

We study analytically and experimentally certain symplectic and time-reversible N-body integrators which employ a Kepler solver for each pair-wise interaction, including the method of Hernandez & Bertschinger (2015). Owing to the Kepler solver, these methods treat close two-body interactions correctly, while close three-body encounters contribute to the truncation error at second order and above. The second-order errors can be corrected to obtain a fourth-order scheme with little computational overhead. We generalise this map to an integrator which employs a Kepler solver only for selected interactions and yet retains fourth-order accuracy without backward steps. In this case, however, two-body encounters not treated via a Kepler solver contribute to the truncation error.

Coalescing massive black hole binaries, formed during galaxy mergers, are expected to be a primary source of low frequency gravitational waves. Yet in isolated gas-free spherical stellar systems, the hardening of the binary stalls at parsec-scale separations owing to the inefficiency of relaxation-driven loss-cone refilling. Repopulation via collisionless orbit diffusion in triaxial systems is more efficient, but published simulation results are contradictory. While sustained hardening has been reported in simulations of galaxy mergers with $N \sim 10^6$ stars and in early simulations of rotating models, in isolated non-rotating triaxial models the hardening rate continues to fall with increasing N, a signature of spurious two-body relaxation. We present a novel approach for studying loss cone repopulation in galactic nuclei. Since loss cone repopulation in triaxial systems owes to orbit diffusion, it is a purely collisionless phenomenon and can be studied with an approximated force calculation technique, provided the force errors are well behaved and sufficiently small. We achieve this using an accurate fast multipole method and define a proxy for the hardening rate that depends only on stellar angular momenta. We find that the loss cone is efficiently replenished even in very mildly triaxial models (with axis ratios 1 : 0.9 : 0.8). Such triaxiality is unavoidable following galactic mergers and can drive binaries into the gravitational wave regime. We conclude that there is no final parsec problem.

Gamma-ray experiments seeking to detect evidence of dark matter annihilation in dwarf spheroidal galaxies require knowledge of the distribution of dark matter within these systems. We analyze the effects of flattening on the annihilation (J) and decay (D) factors of dwarf spheroidal galaxies with both analytic and numerical methods. Flattening has two consequences: first, there is a geometric effect as the squeezing (or stretching) of the dark matter distribution enhances (or diminishes) the J-factor; second, the line of sight velocity dispersion of stars must hold up the flattened baryonic component in the flattened dark matter halo. We provide analytic formulae and a simple numerical approach to estimate the correction to the J- and D-factors required over simple spherical modeling. The formulae are validated with a series of equilibrium models of flattened stellar distributions embedded in flattened dark-matter distributions. We compute corrections to the J- and D-factors for the Milky Way dwarf spheroidal galaxies under the assumption that they are all prolate or all oblate and find that the hierarchy of J-factors for the dwarf spheroidals is slightly altered (typical correction factors for an ellipticity of $0.4$ are $0.75$ for the oblate case and $1.6$ for the prolate case). We demonstrate that spherical estimates of the D-factors are very insensitive to the flattening and introduce uncertainties significantly less than the uncertainties in the D-factors from the other observables for all the dwarf spheroidals (for example, ${}^{+10\mathrm{\,percent}}_{-3\mathrm{\,percent}}$ for a typical ellipticity of $0.4$). We conclude by investigating the spread in correction factors produced by triaxial figures and provide uncertainties in the J-factors for the dwarf spheroidals using different physically-motivated assumptions for their intrinsic shape and axis alignments. (abridged)

We show that, provided the principal axes of the second velocity moment tensor of a stellar population are generally unequal and are oriented perpendicular to a set of orthogonal surfaces at each point, then those surfaces must be confocal quadric surfaces and the potential must be separable or Stackel. This is true under the mild assumption that the even part of the distribution function is invariant under time reversal $v_i \rightarrow -v_i$ of each velocity component. In particular, if the second velocity moment tensor is everywhere exactly aligned in spherical polar coordinates, then the potential must be of Stackel form (excepting degenerate cases where two or more of the semiaxes of ellipsoid are everywhere the same). The theorem also has consequences for alignment in cylindrical polar coordinates, which is used in the popular Jeans Anisotropic Models (JAM). We analyse data on the radial velocities and proper motions of a sample of $\sim 7400$ stars in the stellar halo of the Milky Way. We provide the distributions of the tilt angles or misalignments from the spherical polar coordinate systems. We show that in this sample the misalignment is always small (usually within $3^\circ$) for radii between 7 and 12 kpc. The velocity anisotropy is very radially biased ($\beta \approx 0.7$), and almost invariant across the volume in our study. Finally, we construct a triaxial stellar halo in a triaxial NFW dark matter halo using a made-to-measure method. Despite the triaxiality of the potential, the velocity ellipsoid of the stellar halo is nearly spherically aligned within $\sim6^\circ$ for large regions of space, particularly outside the scale radius of the stellar halo. We conclude that the second velocity moment ellipsoid can be close to spherically aligned for a much wider class of potentials than the strong constraints that arise from exact alignment might suggest.

We investigate the evolution of low mass (Md /Mb = 0.005) misaligned gaseous discs around eccentric supermassive black hole (SMBH) binaries. These are expected to form from randomly oriented accretion events onto a SMBH binary formed in a galaxy merger. When expanding the interaction terms between the binary and a circular ring to quadrupole order and averaging over the binary orbit, we expect four non-precessing disc orientations: aligned or counter-aligned with the binary, or polar orbits around the binary eccentricity vector with either sense of rotation. All other orientations precess around either of these, with the polar precession dominating for high eccentricity. These expectations are borne out by smoothed particle hydrodynamics simulations of initially misaligned viscous circumbinary discs, resulting in the formation of polar rings around highly eccentric binaries in contrast to the co-planar discs around circular binaries. Moreover, we observe disc tearing and violent interactions between differentially precessing rings in the disc significantly disrupting the disc structure and causing gas to fall onto the binary with little angular momentum. While accretion from a polar disc may not promote SMBH binary coalescence (solving the `final-parsec problem'), ejection of this infalling low-angular momentum material via gravitational slingshot is a possible mechanism to reduce the binary separation. Moreover, this process acts on dynamical rather than viscous time scales, and so is much faster.

The Kepler mission has recently discovered a number of exoplanetary systems, such as Kepler-11 and Kepler-32, in which ensembles of several planets are found in very closely packed orbits (often within a few percent of an AU of one another). These compact configurations present a challenge for traditional planet formation and migration scenarios. We present a dynamical study of the assembly of these systems, using an N-body method which incorporates a parametrized model of planet migration in a turbulent protoplanetary disc. We explore a wide parameter space, and find that under suitable conditions it is possible to form compact, close-packed planetary systems via traditional disc-driven migration. We find that simultaneous migration of multiple planets is a viable mechanism for the assembly of tightly-packed planetary systems, as long as the disc provides significant eccentricity damping and the level of turbulence in the disc is modest. We discuss the implications of our preferred parameters for the protoplanetary discs in which these systems formed, and comment on the occurrence and significance of mean-motion resonances in our simulations.

The approximate computation of all gravitational forces between $N$ interacting particles via the fast multipole method (FMM) can be made as accurate as direct summation, but requires less than $\mathcal{O}(N)$ operations. FMM groups particles into spatially bounded cells and uses cell-cell interactions to approximate the force at any position within the sink cell by a Taylor expansion obtained from the multipole expansion of the source cell. By employing a novel estimate for the errors incurred in this process, I minimise the computational effort required for a given accuracy and obtain a well-behaved distribution of force errors. For relative force errors of $\sim10^{-7}$, the computational costs exhibit an empirical scaling of $\propto N^{0.87}$. My implementation (running on a 16 core node) out-performs a GPU-based direct summation with comparable force errors for $N\gtrsim10^5$.

Non-axisymmetries in the Galactic potential (spiral arms and bar) induce kinematic groups such as the Hercules stream. Assuming that Hercules is caused by the effects of the Outer Lindblad Resonance of the Galactic bar, we model analytically its properties as a function of position in the Galaxy and its dependence on the bar's pattern speed and orientation. Using data from the RAVE survey we find that the azimuthal velocity of the Hercules structure decreases as a function of Galactocentric radius, in a manner consistent with our analytical model. This allows us to obtain new estimates of the parameters of the Milky Way's bar. The combined likelihood function of the bar's pattern speed and angle has its maximum for a pattern speed of $\Omega_b=(1.89 \pm 0.08) \times \Omega_0$ where $\Omega_0$ is the local circular frequency. Assuming a Solar radius of 8.05 kpc and a local circular velocity of 238 km/s, this corresponds to $\Omega_b=56\pm 2$ km/s/kpc. On the other hand, the bar's orientation $\phi_b$ cannot be constrained with the available data. In fact, the likelihood function shows that a tight correlation exists between the pattern speed and the orientation, implying that a better description of our best fit results is given by the linear relation $\Omega_b/\Omega_0=1.91 + 0.0044[\phi_b(\deg)-48]$, with standard deviation of 0.02. For example, for an angle of $\phi_b=30\deg$ the pattern speed is $54.0\pm0.5$ km/s/kpc. These results are not very sensitive to the other Galactic parameters such as the circular velocity curve or the peculiar motion of the Sun, and are robust to biases in distance.

We suggest a new picture of supermassive black hole (SMBH) growth in galaxy centers. Momentum-driven feedback from an accreting hole gives significant orbital energy but little angular momentum to the surrounding gas. Once central accretion drops, the feedback weakens and swept-up gas falls back towards the SMBH on near-parabolic orbits. These intersect near the black hole with partially opposed specific angular momenta, causing further infall and ultimately the formation of a small-scale accretion disk. The feeding rates into the disk typically exceed Eddington by factors of a few, growing the hole on the Salpeter timescale and stimulating further feedback. Natural consequences of this picture include (i) the formation and maintenance of a roughly toroidal distribution of obscuring matter near the hole; (ii) random orientations of successive accretion disk episodes; (iii) the possibility of rapid SMBH growth; (iv) tidal disruption of stars and close binaries formed from infalling gas, resulting in visible flares and ejection of hypervelocity stars; (v) super-solar abundances of the matter accreting on to the SMBH; and (vi) a lower central dark-matter density, and hence annihilation signal, than adiabatic SMBH growth implies. We also suggest a simple sub-grid recipe for implementing this process in numerical simulations.

Non-axisymmetric components, such as spirals and central bars, play a major role in shaping galactic discs. An important aspect of the disc secular evolution driven by these perturbers is the radial migration of stars. It has been suggested recently that migration can populate a thick-disc component from inner-disc stars with high vertical energies. Since this has never been demonstrated in simulations, we study in detail the effect of radial migration on the disc velocity dispersion and disc thickness, by separating simulated stars into migrators and non-migrators. We apply this method to three isolated barred Tree-SPH N-body galaxies with strong radial migration. Contrary to expectations, we find that as stellar samples migrate, on the average, their velocity dispersion change (by as much as 50%) in such a way as to approximately match the non-migrating population at the radius at which they arrive. We show that, in fact, migrators suppress heating in parts of the disc. To confirm the validity of our findings, we also apply our technique to three cosmological re-simulations, which use a completely different simulation scheme and, remarkably, find very similar results. We believe the inability of migration to thicken discs is a fundamental property of internal disc evolution, irrespective of the migration mechanism at work. We explain this with the approximate conservation of the (average) vertical and radial actions rather than the energy. This "action mixing" can be used to constrain the migration rate in the Milky Way: estimates of the average vertical action in observations for different populations of stars should reveal flattening with radius for older groups of stars.

We simulate the evolution of one-dimensional gravitating collisionless systems from non- equilibrium initial conditions, similar to the conditions that lead to the formation of dark- matter halos in three dimensions. As in the case of 3D halo formation we find that initially cold, nearly homogeneous particle distributions collapse to approach a final equilibrium state with a universal density profile. At small radii, this attractor exhibits a power-law behavior in density, {\rho}(x) \propto |x|^(-{\gamma}_crit), {\gamma}_crit \simeq 0.47, slightly but significantly shallower than the value {\gamma} = 1/2 suggested previously. This state develops from the initial conditions through a process of phase mixing and violent relaxation. This process preserves the energy ranks of particles. By warming the initial conditions, we illustrate a cross-over from this power-law final state to a final state containing a homogeneous core. We further show that inhomogeneous but cold power-law initial conditions, with initial exponent {\gamma}_i > {\gamma}_crit, do not evolve toward the attractor but reach a final state that retains their original power-law behavior in the interior of the profile, indicating a bifurcation in the final state as a function of the initial exponent. Our results rely on a high-fidelity event-driven simulation technique.

The numerical convergence of smoothed particle hydrodynamics (SPH) can be severely restricted by random force errors induced by particle disorder, especially in shear flows, which are ubiquitous in astrophysics. The increase in the number NH of neighbours when switching to more extended smoothing kernels at fixed resolution (using an appropriate definition for the SPH resolution scale) is insufficient to combat these errors. Consequently, trading resolution for better convergence is necessary, but for traditional smoothing kernels this option is limited by the pairing (or clumping) instability. Therefore, we investigate the suitability of the Wendland functions as smoothing kernels and compare them with the traditional B-splines. Linear stability analysis in three dimensions and test simulations demonstrate that the Wendland kernels avoid the pairing instability for all NH, despite having vanishing derivative at the origin (disproving traditional ideas about the origin of this instability; instead, we uncover a relation with the kernel Fourier transform and give an explanation in terms of the SPH density estimator). The Wendland kernels are computationally more convenient than the higher-order B-splines, allowing large NH and hence better numerical convergence (note that computational costs rise sub-linear with NH). Our analysis also shows that at low NH the quartic spline kernel with NH ~= 60 obtains much better convergence then the standard cubic spline.

Uniquely among the dwarf spheroidal (dSph) satellite galaxies of the Milky Way, Fornax hosts globular clusters. It remains a puzzle as to why dynamical friction has not yet dragged any of Fornax's five globular clusters to the centre, and also why there is no evidence that any similar star cluster has been in the past (for Fornax or any other dSph). We set up a suite of 2800 N-body simulations that sample the full range of globular-cluster orbits and mass models consistent with all existing observational constraints for Fornax. In agreement with previous work, we find that if Fornax has a large dark-matter core then its globular clusters remain close to their currently observed locations for long times. Furthermore, we find previously unreported behaviour for clusters that start inside the core region. These are pushed out of the core and gain orbital energy, a process we call 'dynamical buoyancy'. Thus a cored mass distribution in Fornax will naturally lead to a shell-like globular cluster distribution near the core radius, independent of the initial conditions. By contrast, CDM-type cusped mass distributions lead to the rapid infall of at least one cluster within \Delta t = 1-2Gyr, except when picking unlikely initial conditions for the cluster orbits (\sim 2% probability), and almost all clusters within \Delta t = 10Gyr. Alternatively, if Fornax has only a weakly cusped mass distribution, dynamical friction is much reduced. While over \Delta t = 10Gyr this still leads to the infall of 1-4 clusters from their present orbits, the infall of any cluster within \Delta t = 1-2Gyr is much less likely (with probability 0-70%, depending on \Delta t and the strength of the cusp). Such a solution to the timing problem requires that in the past the globular clusters were somewhat further from Fornax than today; they most likely did not form within Fornax, but were accreted.

Disc non-axisymmetrc components, such as spirals and central bars, are nowadays known to play an important role in shaping galactic discs. Here we use Tree-SPH N-body simulations to examine the effect of these perturbers on two aspects: the occurrence of multiple patterns in discs and the effects of radial migration on disc thickening. We find that, in addition to a central bar, multiple spiral patterns and lopsided modes develop in all models. Interaction among these asymmetric features results in a large scale stellar migration. However, we show that, despite the strong radial mixing, discs cannot be thickened sufficiently to match observed thick discs. We relate this to the adiabatic cooling as stars migrate radially outwards. We also find that the bulge contribution to a thick-disc component for an Sa-type galaxy at ~2.5 disc scale-lengths is less than 1% and zero in the case of a Milky Way-like, Sb-type. Our findings cast doubt on the plausibility of thick disc formation via stellar radial migration.

We consider the infall of a massive clump into a dark-matter halo as a simple and extreme model for the effect of baryonic physics (neglected in gravity-only simulations of large-scale structure formation) on the dark-matter. We find that such an infalling clump is extremely efficient in altering the structure of the halo and reducing its central density: a clump of 1% the mass of the halo can remove about twice its own mass from the inner halo and transform a cusp into a core or weaker cusp. If the clump is subsequently removed, mimicking a galactic wind, the central halo density is further reduced and the mass removed from the inner halo doubled. Lighter clumps are even more efficient: the ratio of removed mass to clump mass increases slightly towards smaller clump masses. This process is the more efficient the more radially anisotropic the initial dark-matter velocities. While such a clumpy infall may be somewhat unrealistic, it demonstrates that the baryons need to transfer only a small fraction of their initial energy to the dark matter via dynamical friction to explain the discrepancy between predicted dark-matter density profiles and those inferred from observations of dark-matter dominated galaxies.

We describe the astrophysical and numerical basis of N-body simulations, both of collisional stellar systems (dense star clusters and galactic centres) and collisionless stellar dynamics (galaxies and large-scale structure). We explain and discuss the state-of-the-art algorithms used for these quite different regimes, attempt to give a fair critique, and point out possible directions of future improvement and development. We briefly touch upon the history of N-body simulations and their most important results.

In smooth-particle hydrodynamics (SPH), artificial viscosity is necessary for the correct treatment of shocks, but often generates unwanted dissipation away from shocks. We present a novel method of controlling the amount of artificial viscosity, which uses the total time derivative of the velocity divergence as shock indicator and aims at completely eliminating viscosity away from shocks. We subject the new scheme to numerous tests and find that the method works at least as well as any previous technique in the strong-shock regime, but becomes virtually inviscid away from shocks, while still maintaining particle order. In particular sound waves or oscillations of gas spheres are hardly damped over many periods.

The chemodynamical evolution of spherical multi-component self-gravitating models for isolated dwarf galaxies is studied. We compare their evolution with and without feedback effects from star formation processes. We find that initially cuspy dark matter profiles flatten with time as a result of star formation, without any special tuning conditions. Thus the seemingly flattened profiles found in many dwarfs do not contradict the cuspy profiles predicted by cosmological models. We also calculate the chemical evolution of stars and gas, to permit comparisons with observational data.

We re-examine the stellar kinematics of the Solar neighbourhood in terms of the velocity of the Sun with respect to the local standard of rest. We show that the classical determination of its component V_sun in the direction of Galactic rotation via Stroemberg's relation is undermined by the metallicity gradient in the disc, which introduces a correlation between the colour of a group of stars and the radial gradients of its properties. Comparing the local stellar kinematics to a chemodynamical model which accounts for these effects, we obtain (U,V,W)_sun = (11.1 +/- 0.74, 12.24 +/- 0.47, 7.25 +/-0.37) km/s, with additional systematic uncertainties of ~ (1,2,0.5) km/s. In particular, V_sun is 7 km/s larger than previously estimated. The new values of solar motion are extremely insensitive to the metallicity gradient within the disc.

The made-to-measure N-body method (Syer & Tremaine 1996) slowly adapts the particle weights of an N-body model, whilst integrating the trajectories in an assumed static potential, until some constraints are satisfied, such as optimal fits to observational data. I propose a novel technique for this adaption procedure, which overcomes several limitations and shortcomings of the original method. The capability of the new technique is demonstrated by generating realistic N-body equilibrium models for dark-matter haloes with prescribed density profile, triaxial shape, and slowly outwardly growing radial velocity anisotropy

We present the second data release of the Radial Velocity Experiment (RAVE), an ambitious spectroscopic survey to measure radial velocities (RVs) and stellar atmosphere parameters of up to one million stars using the 6dF multi-object spectrograph on the 1.2-m UK Schmidt Telescope of the Anglo-Australian Observatory (AAO). It is obtaining medium resolution spectra (median R=7,500) in the Ca-triplet region (8,410--8,795 \AA) for southern hemisphere stars in the magnitude range 9<I<12. Following the first data release (Steinmetz et al. 2006) the current release doubles the sample of published RVs, now containing 51,829 RVs for 49,327 individual stars observed on 141 nights between April 11 2003 and March 31 2005. Comparison with external data sets shows that the new data collected since April 3 2004 show a standard deviation of 1.3 km/s, about twice better than for the first data release. For the first time this data release contains values of stellar parameters from 22,407 spectra of 21,121 individual stars. They were derived by a penalized \chi^2 method using an extensive grid of synthetic spectra calculated from the latest version of Kurucz models. From comparison with external data sets, our conservative estimates of errors of the stellar parameters (for a spectrum with S/N=40) are 400 K in temperature, 0.5 dex in gravity, and 0.2 dex in metallicity. We note however that the internal errors estimated from repeat RAVE observations of 822 stars are at least a factor 2 smaller. We demonstrate that the results show no systematic offsets if compared to values derived from photometry or complementary spectroscopic analyses. The data release includes proper motion and photometric measurements. It can be accessed via the RAVE webpage: http://www.rave-survey.org and through CDS.

We report new constraints on the local escape speed of our Galaxy. Our analysis is based on a sample of high velocity stars from the RAVE survey and two previously published datasets. We use cosmological simulations of disk galaxy formation to motivate our assumptions on the shape of the velocity distribution, allowing for a significantly more precise measurement of the escape velocity compared to previous studies. We find that the escape velocity lies within the range $498\kms < \ve < 608 \kms$ (90 per cent confidence), with a median likelihood of $544\kms$. The fact that $\ve^2$ is significantly greater than $2\vc^2$ (where $\vc=220\kms$ is the local circular velocity) implies that there must be a significant amount of mass exterior to the Solar circle, i.e. this convincingly demonstrates the presence of a dark halo in the Galaxy. For a simple isothermal halo, one can calculate that the minimum radial extent is $\sim58$ kpc. We use our constraints on $\ve$ to determine the mass of the Milky Way halo for three halo profiles. For example, an adiabatically contracted NFW halo model results in a virial mass of $1.42^{+1.14}_{-0.54}\times10^{12}M_\odot$ and virial radius of $305^{+66}_{-45}$ kpc (90 per cent confidence). For this model the circular velocity at the virial radius is $142^{+31}_{-21}\kms$. Although our halo masses are model dependent, we find that they are in good agreement with each other.

We present a general recipe for constructing N-body realizations of galaxies comprised of near-spherical and disc components. First, an exact spherical distribution function for the spheroids (halo & bulge) is determined, such that it is in equilibrium with the gravitational monopole of the disc components. Second, an N-body realisation of this model is adapted to the full disc potential by growing the latter adiabatically from its monopole. Finally, the disc is sampled with particles drawn from an appropriate distribution function, avoiding local-Maxwellian approximations. We performed test simulations and find that the halo and bulge radial density profile very closely match their target model, while they become slightly oblate due to the added disc gravity. Our findings suggest that vertical thickening of the initially thin disc is caused predominantly by spiral and bar instabilities, which also result in a radial re-distribution of matter, rather than scattering off interloping massive halo particles.

We perform collisionless N-body simulations of 1:1 galaxy mergers, using models which include a galaxy halo, disc and bulge, focusing on the behaviour of the halo component. The galaxy models are constructed without recourse to a Maxwellian approximation. We investigate the effect of varying the galaxies' orientation, their mutual orbit, and the initial velocity anisotropy or cusp strength of the haloes upon the remnant halo density profiles and shape, as well as on the kinematics. We observe that the halo density profile (determined as a spherical average, an approximation we find appropriate) is exceptionally robust in mergers, and that the velocity anisotropy of our remnant haloes is nearly independent of the orbits or initial anisotropy of the haloes. The remnants follow the halo anisotropy - local density slope (\beta-\gamma) relation suggested by Hansen & Moore (2006} in the inner parts of the halo, but \beta is systematically lower than this relation predicts in the outer parts. Remnant halo axis ratios are strongly dependent on the initial parameters of the haloes and on their orbits. We also find that the remnant haloes are significantly less spherical than those described in studies of simulations which include gas cooling.

We present the first data release of the Radial Velocity Experiment (RAVE), an ambitious spectroscopic survey to measure radial velocities and stellar atmosphere parameters (temperature, metallicity, surface gravity) of up to one million stars using the 6dF multi-object spectrograph on the 1.2-m UK Schmidt Telescope of the Anglo-Australian Observatory (AAO). The RAVE program started in 2003, obtaining medium resolution spectra (median R=7,500) in the Ca-triplet region ($\lambda\lambda$ 8,410--8,795 \AA) for southern hemisphere stars drawn from the Tycho-2 and SuperCOSMOS catalogs, in the magnitude range 9<I<12. The first data release is described in this paper and contains radial velocities for 24,748 individual stars (25,274 measurements when including re-observations). Those data were obtained on 67 nights between 11 April 2003 to 03 April 2004. The total sky coverage within this data release is $\sim$4,760 square degrees. The average signal to noise ratio of the observed spectra is 29.5, and 80% of the radial velocities have uncertainties better than 3.4 km/s. Combining internal errors and zero-point errors, the mode is found to be 2 km/s. Repeat observations are used to assess the stability of our radial velocity solution, resulting in a variance of 2.8 km/s. We demonstrate that the radial velocities derived for the first data set do not show any systematic trend with color or signal to noise. The RAVE radial velocities are complemented in the data release with proper motions from Starnet 2.0, Tycho-2 and SuperCOSMOS, in addition to photometric data from the major optical and infrared catalogs (Tycho-2, USNO-B, DENIS and 2MASS). The data release can be accessed via the RAVE webpage: http://www.rave-survey.org.

We re-analyse the velocity-dispersion profile sigma(r) at radii r>10 kpc in the Galactic stellar halo, recently derived by Battaglia et al. (2005), who concluded that, for constant velocity anisotropy of the tracers, these data rule out a flat circular-speed curve for the Milky Way. However, we demonstrate that if one makes the natural assumption that the tracer density is truncated at r_t >~ 160 kpc and falls off significantly more steeply than r^-3.5 at r >~ 80 kpc, then these data are consistent with a flat circular-speed curve and a constant velocity anisotropy comparable to that observed for halo stars in the Solar neighbourhood. We also consider a more detailed mass model with an exponential stellar disc and an extended non-isothermal dark-matter halo. In this two-component model, the Milky Way's virial radius and mass are r_vir = 200 kpc and M_vir = 1.5 x 10^12 M_sun. Still assuming the tracers' velocity anisotropy to be constant (at beta=0.5) we again find good agreement with the observed sigma(r), so long as the tracer density is truncated near the virial radius. These data by themselves do not allow to differentiate between different dark-halo or total-mass models for the Milky Way, nor between different velocity-anisotropy profiles for the tracers.

Low-mass X-ray binaries (LMXBs), which occur in old stellar populations, have velocities exceeding those of their parent distribution by at least 20 km/s. This makes them ideal probes for dark matter, in particular in dwarf spheroidals (dSph), where the LMXBs should penetrate well outside the visible galaxy. We argue that the most likely explanation of the observation of LMXBs in the Sculptor dSph by Maccarone et al (2005) is the presence of a dark matter halo of >~ 10^9 Msun, corresponding to a total-mass to light ratio of >~ 600 (M/L_V)_sun. In this case there should be an extended halo of LMXBs which may be observable.

It is believed that young massive stars orbiting Sgr A* in two stellar discs on scales of 0.1-0.2 parsecs were formed either farther out in the Galaxy and then quickly migrated inward, or in situ in a massive self-gravitating disc. Comparing N-body evolution of stellar orbits with observational constraints, we set upper limits on the masses of the two stellar systems. These masses turn out to be few times lower than the expected total stellar mass estimated from the observed young high-mass stellar population and the standard galactic IMF. If these stars were formed in situ, in a massive self-gravitating disc, our results suggest that the formation of low-mass stars was suppressed by a factor of at least a few, requiring a top-heavy initial mass function (IMF) for stars formed near sgr A*.

Angular momentum transfer from a rotating stellar bar has been proposed by Weinberg & Katz (2002) as a mechanism to destroy dark-matter cusps in a few rotation periods. The N-body simulations performed by these authors in support of their claim employed spherical harmonics for the force computation and were, as shown by Sellwood (2003), very sensitive to inclusion of asymmetric terms (odd l,m). In order to disentangle possible numerical artifacts due to the usage of spherical harmonics from genuine stellar dynamical effects, we performed similar experiments using a tree code and find that significant cusp destruction requires substantially more angular momentum than is realistically available. However, we find that the simplified model (a N-body halo torqued by a rotating bar pinned to the origin) undergoes an instability in which the cusp moves away from the origin. In presence of this off-centring, spherical density profiles centred on the origin display an apparent cusp-removal. We strongly suspect that it is this effect which Weinberg & Katz observed. When suppressing the artificial instability, cusp removal is very slow and requires much more angular momentum to be transferred to the halo than a realistic stellar bar possibly possesses.

We study the effect of a central mass concentration (CMC) on the secular evolution of a barred disc galaxy. Unlike previous studies, we use fully self-consistent 3D $N$-body simulations with live haloes, which are known to be important for bar evolution. The CMC is introduced gradually, to avoid transients. In all cases where the mass of the CMC is of the order of, or more than, a few per cent of the mass of the disc, the strength of the bar decreases noticeably. The amount of this decrease depends strongly on the bar type. For the same CMC, bars with exponential surface-density profile, which formed in a disk-dominated galaxy (MD-type bars), can be totally destroyed, while strong bars with a flat surface-density profile, whose evolution is largely due to the halo (MH-type bars), witness only a decrease of their strength. This decrease occurs simultaneously from both the innermost and outermost parts of the bar. The CMC has a stronger effect on the Fourier components of higher azimuthal wave number $m$, leading to fatter and/or less rectangular bars. Furthermore, the CMC changes the side-on outline from peanut-shaped to boxy or, for massive CMCs, to elliptical. Similarly, side-on initially boxy outlines can be destroyed. The CMC also influences the velocity dispersion profiles. Most of the decrease of the bar strength occurs while the mass of the CMC increases and it is accompanied by an increase of the pattern speed. In all our simulations, the mass of the CMC necessary in order to destroy the bar is at least several per cent of the mass of the disc. This argues that observed super-massive black holes are not likely to destroy pre-existing bars.

We investigate, using the spherical Jeans equation, self-gravitating dynamical equilibria satisfying a relation rho/sigma_r^3 propto r^-alpha, which holds for simulated dark-matter haloes over their whole resolved radial range. Considering first the case of velocity isotropy, we find that this problem has only one solution with realistic density profile, which occurs only for a critical value of alpha_crit = 35/18 ~= 1.94, which is consistent with the empirical value of 1.9+/-0.05. We extend our analysis in two ways: first we introduce a parameter epsilon to allow for a more general relation rho/\sigma_r^epsilon propto r^-alpha; second we consider velocity anisotropy. If we assume beta(r) := 1- sigma_theta^2 / sigma_r^2 to be linearly related to the logarithmic density slope gamma(r) := -dln(rho)/dln(r), which is in agreement with simulations, the problem remains analytically tractable and is equivalent to the simpler isotropic case: there exists only one physical solution, which occurs at a critical alpha value. Remarkably, this value of alpha and the density and velocity-dispersion profiles depend only on epsilon and the value beta_0 := beta(r=0), but not on the slope of the linear beta-gamma relation. For epsilon=3, alpha_crit = 35/18 - 2beta_0/9 and the resulting density profile is fully analytic (as are the velocity dispersion and circular speed) with an inner cusp rho propto r^{-(7+10beta_0)/9} and a very smooth transition to a steeper outer power-law asymptote. These models are in excellent agreement with the density, velocity-dispersion and anisotropy profiles of simulated dark-matter haloes over their full resolved radial range. If epsilon=3 is a universal constant, some scatter in beta_0 ~= 0 may account for some diversity in the density profiles. (ABRIDGED)

Collisionless stellar systems are driven towards equilibrium by mixing of phase-space elements. I show that the excess-mass function D(f)= int_{F(x,v)>f} (F(x,v)-f) d^3x d^3v (with F(x,v) the coarse-grained distribution function) always decreases on mixing. D(f) gives the excess mass from values of F(x,v))>f. This novel form of the mixing theorem extends the maximum phase-space density argument to all values of f. The excess-mass function can be computed from N-body simulations and is additive: the excess mass of a combination of non-overlapping systems is the sum of their individual D(f). I propose a novel interpretation for the coarse-grained distribution function, which avoids conceptual problems with the mixing theorem. As an example application, I show that for self-gravitating cusps (rho propto r^{-gamma} as r->0) the excess mass D propto f^{-2(3-gamma)/(6-gamma)} as f->oo, i.e. steeper cusps are less mixed than shallower ones, independent of the shape of surfaces of constant density or details of the distribution function (e.g. anisotropy). This property, together with the additivity of D(f) and the mixing theorem, implies that a merger remnant cannot have a cusp steeper than the steepest of its progenitors. Furthermore, I argue that the remnant's cusp should not be shallower either, implying that the steepest cusp always survives.

Various arguments strongly suggest that the population of ultraluminous X-ray sources (ULXs: apparent X-ray luminosity > Eddington limit for 10 Msun ~10^39 erg/s) in nearby galaxies are mostly stellar-mass X-ray binaries in unusual evolutionary stages. However there are indications that the very brightest systems may be difficult to explain this way. Accordingly we consider the class of hyperluminous X-ray sources (HLXs) (i.e. those with apparent bolometric luminosities >~10^41 erg/s). Because this class is small (currently only the M82 object is a secure member) we do not need to invoke a new formation mechanism for its black holes. We explore instead the idea that HLXs may be the nuclei of satellite galaxies captured during hierarchical merging. The observed correlation between AGN and tidal interactions implies that HLX activity would switch on during passage through the host galaxy, close to pericentre. This suggests that HLXs should appear near the host galaxy, be associated with star formation, and thus possibly with ULXs.

In order to understand the extended massive tidal tails of the globular cluster Pal5, its very low mass and velocity dispersion, and its size, which is much larger than the theoretical tidal radius, we performed more than 1000 N-body simulations. Tidal shocks at disk crossings near perigalacticon dominate the evolution of extended low-concentration clusters, resulting in massive tidal tails and often in a quick destruction of the cluster. The overlarge size of Pal5 can be explained as the result of an expansion following the heating induced by the last strong disk shock ~150Myr ago. The models can reproduce the low observed velocity dispersion and the relative fractions of stars in the tails and between the inner and outer parts of the tails. The tidal tails of Pal5 show substantial structure not seen in our simulations. We argue that this structure is probably caused by Galactic substructure, such as giant molecular clouds, spiral arms, and dark-matter clumps. Clusters initially larger than their theoretical tidal limit remain so, because, after being shocked, they settle into a new equilibrium near apogalacticon, where they are unaffected by the perigalactic tidal field. This implies that, contrary to previous wisdom, globular clusters on eccentric orbits may well remain super-tidally limited and hence vulnerable to strong disk shocks, which dominate their evolution until destruction. Our simulations unambiguously predict the destruction of Pal5 at its next disk crossing in \~110Myr. This corresponds to only 1% of the cluster lifetime, suggesting that many more similar systems could once have populated the inner parts of the Milky Way, but have been transformed into debris streams by the Galactic tidal field. (Abridged)

We present an analysis of the asymmetries in the population of Galactic M-giant stars present in the 2MASS All Sky catalogue. Several large-scale asymmetries are detected, the most significant of which is a strong elliptical-shaped stellar over-density, close to the Galactic plane at (l=240, b=-8), in the constellation of Canis Major. A small grouping of globular clusters (NGC 1851, NGC 1904, NGC 2298, and NGC 2808), coincident in position and radial velocity, surround this structure, as do a number of open clusters. The population of M-giant stars in this over-density is similar in number to that in the core of the Sagittarius dwarf galaxy. We argue that this object is the likely dwarf galaxy progenitor of the ring-like structure that has recently been found at the edge of the Galactic disk. A numerical study of the tidal disruption of an accreted dwarf galaxy is presented. The simulated debris fits well the extant position, distance and velocity information on the ``Galactic Ring'', as well as that of the M-giant over-densities, suggesting that all these structures are the consequence of a single accretion event. The disrupted dwarf galaxy stream orbits close to the Galactic Plane, with a pericentre at approximately the Solar circle, an orbital eccentricity similar to that of stars in the Galactic thick disk, as well as a vertical scale height similar to that of the thick disk. This finding strongly suggests that the Canis Major dwarf galaxy is a building block of the Galactic thick disk, that the thick disk is continually growing, even up to the present time, and that thick disk globular clusters were accreted onto the Milky Way from dwarf galaxies in co-planar orbits.

Using wide-field photometric data from the Sloan Digital Sky Survey (SDSS) we recently showed that the Galactic globular cluster Palomar 5 is in the process of being tidally disrupted. Its tidal tails were initially detected in a 2.5 degree wide band along the celestial equator. A new analysis of SDSS data for a larger field now reveals that the tails of Pal 5 have a much larger spatial extent and can be traced over an arc of 10 deg across the sky, corresponding to a projected length of 4 kpc at the distance of the cluster. The number of former cluster stars found in the tails adds up to about 1.2 times the number of stars in the cluster. The radial profile of stellar surface density in the tails follows approximately a power law r^gamma with -1.5 < gamma < -1.2. The stream of debris from Pal 5 is significantly curved, which demonstrates its acceleration by the Galactic potential. The cluster is presently near the apocenter but has repeatedly undergone disk crossings in the inner part of the Galaxy leading to strong tidal shocks. Our results suggest that the observed debris originates mostly from mass loss within the last 2 Gyrs. The cluster is likely to be destroyed after the next disk crossing, which will happen in about 100 Myr. (abridged)

We study, using direct orbit integrations, the kinematic response of the outer stellar disk to the presence of a central bar, as in the Milky-Way. We find that the bar's outer Lindblad resonance (OLR) causes significant perturbations of the velocity moments. With increasing velocity dispersion, the radius of these perturbations is shifted outwards, beyond the nominal position of the OLR, but also the disk becomes less responsive. If we follow Dehnen (2000) in assuming that the OLR occurs just inside the Solar circle and that the Sun lags the bar major axis by ~20 degrees, we find (1) no significant radial motion of the local standard of rest (LSR), (2) a vertex deviation of \~10 degrees and (3) a lower ratio sigma_2/sigma_1 of the principal components of the velocity- dispersion tensor than for an unperturbed disk. All of these are actually consistent with the observations of the Solar-neighbourhood kinematics. Thus it seems that at least the lowest-order deviations of the local-disk kinematics from simple expectations based on axisymmetric equilibrium can be attributed entirely to the influence of the Galactic bar.

The Oort constants describe the local variations of the stellar streaming field. Classically, they are determined from stellar proper motions. We discuss problems arising in this procedure. A large, hitherto overlooked, source of systematic error arises from longitudinal variations of the mean stellar parallax. Together with the solar reflex motion, these variations create contributions to the longitudinal proper motions mu_l*(l) that are indistinguishable from the Oort Constants at the 20% level. Fortunately, we can correct for this "mode mixing" using the latitudinal proper motions mu_b(l). We use ~10^6 stars from the ACT/Tycho-2 catalogs with proper motion errors of ~ 3 mas/yr. We find significant deviations from expectations based on a smooth axisymmetric equilibrium disk, in particular non-zero C for old red giant stars. We also find variations of the Oort constants with the asymmetric drift of the sub-sample. These correlations are different in nature than those expected for an axisymmetric Galaxy. The most reliable tracers for the ``true'' Oort constants are red giants, which are old enough to be in equilibrium and distant enough to be unaffected by possible local anomalies. For these stars we find, A ~ 16, B ~- 17, A - B ~ 33, and C ~- 10 km/s/kpc with internal errors of about 1-2 and external error of perhaps the same order. These values are consistent with our knowledge of the Milky Way (flat rotation curve and Omega=A - B ~ 28 +/- 2). Based on observations made with the ESA Hipparcos astrometry satellite. (Abridged)

We detected extended, curved stellar tidal tails emanating from the sparse, disrupting halo globular cluster Pal 5, which cover 10 degrees on the sky. These streams allow us to infer the orbit of Pal 5 and to ultimately constrain the Galactic potential at its location.

Sloan Digital Sky Survey data for the field of the globular cluster Pal 5 reveal the existence of a long massive stream of tidal debris spanning an arc of 10 degrees on the sky. Pal 5 thus provides an outstanding example for tidal disruption of globular clusters in the Milky Way. Radial velocities from VLT spectra show that Pal 5 has an extremely low velocity dispersion, in accordance with the very low mass derived from its total luminosity.

Wide-field photometric data from the Sloan Digital Sky Survey have recently revealed that the Galactic globular cluster Palomar 5 is in the process of being tidally disrupted (Odenkirchen et al. 2001). Here we investigate the kinematics of this sparse remote star cluster using high resolution spectra from the Very Large Telescope (VLT). Twenty candidate cluster giants located within 6 arcmin of the cluster center have been observed with the UV-Visual Echelle Spectrograph (UVES) on VLT-UT2. The spectra provide radial velocities with a typical accuracy of 0.15 km/s. We find that the sample contains 17 certain cluster members with very coherent kinematics, two unrelated field dwarfs, and one giant with a deviant velocity, which is most likely a cluster binary showing fast orbital motion. From the confirmed members we determine the heliocentric velocity of the cluster as -58.7 +- 0.2 km/s. The total line-of-sight velocity dispersion of the cluster stars is 1.1 +- 0.2 km/s (all members) or 0.9 +- 0.2 km/s (stars on the red giant branch only). This is the lowest velocity dispersion that has so far been measured for a stellar system classified as a globular cluster. The shape of the velocity distribution suggests that there is a significant contribution from orbital motions of binaries and that the dynamical part of the velocity dispersion is therefore still substantially smaller than the total dispersion. ... (abridged)

A novel code for the approximate computation of long-range forces between N mutually interacting bodies is presented. The code is based on a hierarchical tree of cubic cells and features mutual cell-cell interactions which are calculated via a Cartesian Taylor expansion in a symmetric way, such that total momentum is conserved. The code benefits from an improved and simple multipole acceptance criterion that reduces the force error and the computational effort. For N>10^4, the computational costs are found empirically to rise sub-linear with N. For applications in stellar dynamics, this is the first competitive code with complexity O(N), it is faster than the standard tree code by a factor ten or more.

We have investigated the spatial extent and structure of the Draco dwarf spheroidal galaxy using deep wide-field multicolor CCD photometry from the Sloan Digital Sky Survey (SDSS). Our study covers an area of 27 square degrees around the center of the Draco dwarf. We show that the spatial distribution of Draco's red giants, red horizontal branch stars and subgiants down to i=21.7 mag does not provide evidence for the existence of tidally induced tails or a halo of unbound stars. The radial profile can be fit by King models as well as by a generalized exponential. The core radius and the limiting (or tidal) radius along the major axis are 7.7' and 40.1', respectively, making Draco 40% larger than previously measured. Down to our magnitude limit tidal effects can only exist at a level of 10^-3 of the central surface density of Draco or below. The regular structure of Draco found from the new data argues against it being a portion of an unbound tidal stream and lends support to the assumption of dynamical equilibrium. We estimate Draco's total mass to be 2.2 - 3.5 times 10^7 solar masses. We obtain an overall mass-to-light ratio of 146+-42 or 92+-28 solar masses depending on the details of the mass and luminosity estimates. In summary, our results strengthen the case for a strongly dark matter dominated, bound stellar system. (Abstract strongly abridged).

We report the discovery of two well-defined tidal tails emerging from the sparse remote globular cluster Palomar 5. These tails stretch out symmetrically to both sides of the cluster in the direction of constant Galactic latitude and subtend an angle of 2.6 degrees on the sky. The tails have been detected in commissioning data of the Sloan Digital Sky Survey (SDSS), providing deep five-color photometry in a 2.5 degrees wide band along the equator. The stars in the tails make up a substantial part (~1/3) of the current total population of cluster stars in the magnitude interval 19.5 < i* < 22.0. This reveals that the cluster is subject to heavy mass loss. The orientation of the tails provides an important key for the determination of the cluster's Galactic orbit.

In N-body simulations of collisionless stellar systems, the forces are softened to reduce the shot noise. Softening modifies gravity at r=|x-y| smaller than softening length epsilon and the softened forces are increasingly biased for ever larger epsilon. There is, thus, some optimum between reducing the fluctuations and introducing a bias. Here, analytical relations are derived for the amplitudes of the bias and the fluctuations in the limit of small epsilon and large N. It is shown that the fluctuations of the force are generated locally, in contrast to the variations of the potential, which originate from noise in the whole system. Based on the asymptotic relations and using numerical experiments, I study the dependence of the resulting force error on N, epsilon, and on the functional form by which Newtonian gravity is replaced. The Plummer softening, where each body is replaced by a Plummer sphere of scale radius epsilon, yields significantly larger force errors than do methods in which the bodies are replaced by density kernels of finite extent. I also give special kernels, which reduce the errors even further. These kernels largely compensate the errors made with too small inter-particle forces at r<epsilon by exceeding Newtonian forces at r epsilon. Additionally, the possibilities of locally adapting epsilon and of using unequal weights for the bodies are investigated. These various techniques allow, without increasing N, to reduce the rms force error by a factor 2 compared to Plummer softening with constant epsilon. The results of this study are directly relevant to N-body simulations using direct summation techniques or the tree code. (abridged)

High-quality Hipparcos data for a complete sample of nearly 12000 main-sequence and subgiant stars, together with Padua isochrones, are used to constrain the star-formation history of the solar neigbourhood and the processes that stochastically accelerate disk stars. The velocity dispersion of a coeval group of stars is found to increase with time from ~8 kms at birth as t^{0.33}. In the fits, the slope of the IMF near 1 Msun proves to be degenerate with the rate at which the star-formation rate declines. If the slope of the IMF is to lie near Salpeter's value, -2.35, the star-formation rate has to be very nearly constant. The age of the solar neighbourhood is found to be 11.2+/-0.75 Gyr with remarkably little sensitivity to variations in the assumed metallicity distribution of old disk stars. This age is only a Gyr younger than the age of the oldest globular clusters when the same isochrones and distance scale are employed. It is compatible with current indications of the redshift of luminous galaxy formation only if there is a large cosmological constant. A younger age is formally excluded because it provides a poor fit to the number density of red stars. Since this density is subject to a significantly uncertain selection function, ages as low as 9 Gyr are plausible even though they lie outside our formal error bars.

The tree code for the approximate evaluation of gravitational forces is extended and substantially accelerated by including mutual cell-cell interactions. These are computed by a Taylor series in Cartesian coordinates and in a completely symmetric fashion, such that Newton's third law is satisfied by construction and hence momentum exactly conserved. The computational effort is further reduced by exploiting the mutual symmetry of the interactions. For typical astrophysical problems with N=10^5 and at the same level of accuracy, the new code is about four times faster than the tree code. For large N, the computational costs are found to scale almost linearly with N, which can also be supported by a theoretical argument, and the advantage over the tree code increases with ever larger N.

Hydro-dynamical modeling of the inner Galaxy suggest that the radius of the outer Lindblad resonance (OLR) of the Galactic bar lies in the vicinity of the Sun. How does this resonance affect the distribution function in the outer parts of a barred disk, and can we identify any effect of the resonance in the velocity distribution f(v) actually observed in the solar neighborhood? To answer these questions, detailed simulations of f(v) in the outer parts of an exponential stellar disks with nearly flat rotation curves and a rotating central bar have been performed. For a model resembling the old stellar disk, the OLR causes a distinct feature in f(v) over a significant fraction of the outer disk. For positions <2kpc outside the OLR radius and at bar angles of \~10-70 degrees, f(v) inhibits a bi-modality between the low-velocity stars moving like the local standard of rest (LSR) and a secondary mode of stars predominantly moving outward and rotating more slowly than the LSR. Such a bi-modality is indeed present in f(v) inferred from the Hipparcos data for late-type stars in the solar neighborhood. If one interpretes this observed bi-modality as induced by the OLR -- and there are hardly any viable alternatives -- then one is forced to deduce that the OLR radius is slightly smaller than Ro. Moreover, by a quantitative comparison of the observed with the simulated distributions one finds that the pattern speed of the bar is 1.85+/-0.15 times the local circular frequency, where the error is dominated by the uncertainty in bar angle and local circular speed. Also other, less prominent but still significant, features in the observed f(v) resemble properties of the simulated velocity distributions, in particular a ripple caused by orbits trapped in the outer 1:1 resonance.

Most late-type stars in the solar neighborhood have velocities similar to the local standard of rest (LSR), but there is a clearly separated secondary component corresponding to a slower rotation and a mean outward motion. Detailed simulations of the response of a stellar disk to a central bar show that such a bi-modality is expected from outer-Lindblad resonant scattering. When constraining the run of the rotation curve by the proper motion of Sgr A* and the terminal gas velocities, the value observed for the rotation velocity separating the two components results in a value of (53+/-3)km/s/kpc for the pattern speed of the bar, only weakly dependent on the precise values for Ro and bar angle phi.

Already slightly eccentric orbits, such as those occupied by many old stars in the Galactic disk, are not well approximated by Lindblad's epicycle theory. Here, alternative approximations for flat orbits in axisymmetric stellar systems are derived and compared to results from numeric integrations. All of these approximations are more accurate than Lindblad's classical theory. I also present approximate, but canonical, maps from ordinary phase-space coordinates to a set of action-angle variables. Unfortunately, the most accurate orbit approximation leads to non-analytical R(t). However, from this approximation simple and yet very accurate estimates can be derived for the peri- and apo-centers, frequencies, and actions integrals of galactic orbits, even for high eccentricities. Moreover, further approximating this approximation allows for an analytical R(t) and still an accurate approximation to galactic orbits, even with high eccentricities.

Distribution functions (DFs) for dynamically warm thin stellar disks residing in arbitrary axisymmetric potentials are presented which approximately reproduce pre-described surface-density and velocity-dispersion profiles. The functional form of the DFs is obtained by `warming-up' a model made entirely of circular orbits. This can be done in various ways giving different functional forms for the DF. In the best case, the DF reproduces the pre-described profiles to within a few per cent for a typical case reminiscent to the old stellar disk in the Milky Way. This match may be improved to about one per cent or better by a simple iterative method. An algorithm is given to draw phase-space points randomly from the DFs for the purpose of, e.g., N-body simulations. All the relevant computer programs are available from the author.

An extension of Schwarzschild's galaxy-building technique is presented that, for the first time, enables one to build Schwarzschild models with known distribution functions (DFs). The new extension makes it possible to combine a DF that depends only on classical integrals with orbits that respect non-classical integrals. With such a combination, Schwarzschild's orbits are used only to represent the difference between the true galaxy DF and an approximating classical DF. The new method is used to construct a dynamical model of the inner Galaxy. The model is based on an orbit library that contains 22168 regular orbits. The model aims to reproduce the three-dimensional mass density of Binney, Gerhard & Spergel (1997), which was obtained through deprojection of the COBE surface photometry, and to reproduce the observed kinematics in three windows - namely Baade's Window and two off-axis fields. The model fits essentially all the available data within the innermost 3 kpc. The axis ratio and the morphology of the projected density contours of the COBE bar are recovered to good accuracy within corotation. The kinematic quantities - the line-of-sight streaming velocity and velocity dispersion, as well as the proper motions when available - are recovered, not merely for the fitted fields, but also for three new fields. The dynamical model deviates most from the input density close to the Galactic plane just outside corotation, where the deprojection of the surface photometry is suspect. The dynamical model does not reproduce the kinematics at the most distant window, where disk contamination may be severe.

From the tangential velocities of stars all over the sky, one can, in a statistical way, infer their 3D velocity distribution. An application to Hipparcos data reveals rich structure in the planar stellar motions: there are several moving groups, increasing in number with stellar type but decreasing in importance. A distinct group of outward moving stars with low rotation velocities might be associated with the Galactic bar: growing a bar in a smooth model for the Galactic disk results in such a group, provided the Sun is outside the outer Lindblad resonance (OLR) and the orientation angle with the bar is in the first quadrant. The vertical motions show less structure, but the mean <v_z> increases for large rotational velocities, which can be nicely explained by a warp of the outer stellar disk.

(abridged) The velocity distribution f(v) of nearby stars is estimated, via a maximum- likelihood algorithm, from the positions and tangential velocities of a kinematically unbiased sample of 14369 stars observed by the HIPPARCOS satellite. f(v) shows rich structure in the radial and azimuthal motions, v_R and v_phi, but not in the vertical velocity, v_z: there are four prominent and many smaller maxima, many of which correspond to well known moving groups. While samples of early-type stars are dominated by these maxima, also up to 25% of red main-sequence stars are associated with them. These moving groups are responsible for the vertex deviation measured even for samples of late-type stars; they appear more frequently for ever redder samples; and as a whole they follow an asymmetric-drift relation, in the sense that those only present in red samples predominantly have large |v_R| and lag in v_phi w.r.t. the local standard of rest (LSR). The question arise, how these old moving groups got on their eccentric orbits. A plausible mechanism, known from solar system dynamics, which is able to manage a shift in orbit space involves locking into an orbital resonance. Apart from these moving groups, there is a smooth background distribution, akin to Schwarzschild's ellipsoidal model, with axis ratio of about 1:0.6:0.35 in v_R, v_phi, and v_z. The contours are aligned with the $v_r$ direction, but not w.r.t. the v_phi and v_z axes: the mean v_z increases for stars rotating faster than the LSR. This effect can be explained by the stellar warp of the Galactic disk. If this explanation is correct, the warp's inner edge must not be within the solar circle, while its pattern rotates with frequency of about 13 km/s/kpc or more retrograde w.r.t. the stellar orbits.

I review recent studies about the gravitational potential and stellar dynamics of M87 in particular, and the dynamics of the stars in the presence of a super-massive central black hole, in general. At large radii, investigations of both the X-ray emitting gas and the velocity distribution of globular clusters indicate the presence of large amounts of non-luminous matter, possibly belonging to the inner parts of the Virgo cluster. At small radii, there is no evidence from the stellar kinematics, at most a hint, for the existence of a central point mass, whereas the gas dynamics reveal the presence of a highly concentrated mass in the centre of M87, possibly a super-massive black hole (BH). Given the existence of such a central mass, the stellar kinematics indicate a strong tangential anisotropy of the stellar motion inside a few arcseconds. The implications of this result for the evolution and formation history of M87 and its central BH are discussed. I also discuss in more general terms the structural changes that a highly concentrated central mass can induce in its parent galaxy.

(shortened) From a kinematically unbiased subsample of the Hipparcos catalogue we have redetermined as a function of colour the kinematics of main-sequence stars. The stars' mean heliocentric velocity nicely follows the asymmetric drift relation, except for stars blueward of B-V=0.1. Extrapolating to zero dispersion yields for the velocity of the Sun w.r.t. the LSR in km/s: U_0=10.00+/-0.36 (radially inwards), V_0=5.23+/-0.62 (in direction of galactic rotation), and W_0=7.17+/-0.38 (vertically upwards). A plot of velocity dispersion vs. colour beautifully shows Parenago's discontinuity: the dispersion is constant for B-V>0.62 and decreases towards bluer colour. We determine the velocity-dispersion tensor sigma^2_ij as function of B-V. The mixed moments involving vertical motion are zero within the errors, while sigma^2_xy is non-zero at about (10km/s)^2 independent of colour. The resulting vertex deviations are about 20 deg for early-type stars and 10+/-4 deg for old-disc stars. The persistence of the vertex deviation to late-type stars implies that the Galactic potential is significantly non-axisymmetric at the solar radius. If spiral arms are responsible for this, they cannot be tightly wound. Except for stars bluer than B-V=0.1 the ratios of the principal velocity dispersions are 2.2 : 1.4 :1, while the absolute values increase with colour from sigma_1=20km/s at B-V=0.2 to sigma_1=38km/s at Parenago's discontinuity and beyond. These ratios imply significant heating of the disc by spiral structure and that R_0/R_d=3 to 3.5, where R_d is the scale length of the disc.

A parameterized model of the mass distribution within the Milky Way is fitted to the available observational constraints. The most important single parameter is the ratio of the scale length R_d* of the stellar disk to R0. The disk and bulge dominate v_c(R) at R<R0 only for R_d*/R0< 0.3. Since the only knowledge we have of the halo derives from studies like the present one, we allow it to contribute to the density at all radii. When allowed this freedom, however, the halo causes changes in assumptions relating to R << R0 to affect profoundly the structure of the best-fitting model at R >> R0. For example, changing the disk slightly from an exponential surface-density profile significantly changes the form of v_c(R) at R >> R0, where the disk makes a negligible contribution to v_c. Moreover, minor changes in the constraints can cause the halo to develop a deep hole at its centre that is not physically plausible. These problems call into question the proposition that flat rotation curves arise because galaxies have physically distinct halos rather than outwards-increasing mass-to-light ratios. The mass distribution of the Galaxy and the relative importance of its various components will remain very uncertain until more observational data can be used to constrain mass models. Data that constrain the Galactic force field at z > R and at R > R0 are especially important.

We present the kinematics of the ionized gas over the inner 140 arcsec (10 kpc) of NGC 1068 from observations with the HIFI Fabry-Perot interferometer. There is clear evidence for density wave streaming and bar-driven streaming across the field, with bi-symmetric arms that penetrate to within 200 pc of the nucleus. CO maps show linear structures along (although slightly offset from) the bar consistent with a strong shock. Along the spiral arms which encircle the bar, the HII regions lie downstream of the CO gas in the rest frame of the bar, as do the dust lanes, only if the gas outruns the stellar bar. As a first step towards understanding the details of the gas kinematics, and attempting to determine the mass inflow rate towards the nucleus, we build a mass model for the central disk constrained by near-infrared images. We plan to use this model as gravitational background potential for hydro-dynamical simulations of the gas response to the bar. Comparing these with the data presented should enable us to constrain various quantities such as pattern speed, stellar mass-to-light ratio, central mass concentration, and gas fueling rate.

A straightforward determination of the circular-speed curve vc(R) of the Milky Way suggests that near the Sun, vc starts to rise approximately linearly with R. If this result were correct, the Galactic mass density would have to be independent of radius at R ~> R0. We show that the apparent linear rise in v_c arises naturally if the true circular-speed curve is about constant or gently falling at R0 < R ~< 2 R0, but most tracers that appear to be at R ~> 1.25 R0 are actually concentrated into a ring of radius ~1.6 R0.

The only way to map the Galaxy's gravitational potential $\Phi({\bf x})$ and the distribution of matter that produces it is by modelling the dynamics of stars and gas. Observations of the kinematics of gas provide key information about gradients of $\Phi$ within the plane, but little information about the structure of $\Phi$ out of the plane. Traditional Galaxy models {\em assume}, for each of the Galaxy's components, arbitrary flattenings, which together with the components' relative masses yield the model's equipotentials. However, the Galaxy's isopotential surfaces should be {\em determined\/} directly from the motions of stars that move far from the plane. Moreover, from the kinematics of samples of such stars that have well defined selection criteria, one should be able not only to map $\Phi$ at all positions, but to determine the distribution function $f_i({\bf x},{\bf v})$ of each stellar population $i$ studied. These distribution functions will contain a wealth of information relevant to the formation and evolution of the Galaxy. An approach to fitting a wide class of dynamical models to the very heterogeneous body of available data is described and illustrated.

A technique for the construction of axisymmetric distribution functions for individual galaxies is presented. It starts from the observed surface bright- ness distribution, which is deprojected to gain the axisymmetric luminosity density, from which follows the stars' gravitational potential. After adding dark mass components, such as a central black hole, the two-integral distribu- tion function (2I-DF) f(E,Lz), which depends only on the classical integrals of motion in an axisymmetric potential, is constructed using the Richardson- Lucy algorithm. This algorithm proved to be very efficient in finding f(E,Lz) provided the integral equation to be solved has been properly modified. Once the 2I-\df\ is constructed, its kinematics can be computed and compared with those observed. Many discrepancies may be remedied by altering the assumed inclination angle, mass-to-light ratio, dark components, and odd part of the 2I-DF. Remaining discrepancies may indicate, that the distribution function depends on the non-classical third integral, or is non-axisymmetric. The method has been applied to the nearby elliptical galaxy M32. A 2I-DF with ~55 degrees inclination and a central black hole (or other compact dark mass inside ~1pc) of 1.6-2*10^6 Msun fits the high-spatial-resolution kinema- tic data of van der Marel et al. remarkably well. 2I-DFs with a significantly less or more massive central dark mass or with edge-on inclination can be ruled out for M32. Predictions are made for HST-observations: spectroscopy using its smallest square aperture of 0.09"*0.09" should yield a non-gaussian central velocity profile with broad wings, true and gaussian-fit velocity dispersion of 150-170km/s and 120-130km/s, respectively.