nmagic: a fast parallel implementation of a χ2-made-to-measure algorithm for modelling observational data

We describe a made-to-measure (M2M) algorithm for constructing N -particle models of stellar systems from observational data (χ²M2M), extending earlier ideas by Syer & Tremaine. The algorithm properly accounts for observational errors, is flexible, and can be applied to various systems and geometries. We implement this algorithm in a parallel code nmagic and carry out a sequence of tests to illustrate its power and performance. (i) We reconstruct an isotropic Hernquist model from density moments and projected kinematics and recover the correct differential energy distribution and intrinsic kinematics. (ii) We build a self-consistent oblate three-integral maximum rotator model and compare how the distribution function is recovered from integral field and slit kinematic data. (iii) We create a non-rotating and a figure rotating triaxial stellar particle model, reproduce the projected kinematics of the figure rotating system by a non-rotating system of the same intrinsic shape, and illustrate the signature of pattern rotation in this model. From these tests, we comment on the dependence of the results from χ²M2M on the initial model, the geometry, and the amount of available data.     

Orbits supporting bars within bars

High-resolution observations of the inner regions of barred disc galaxies have revealed many asymmetrical, small-scale central features, some of which are best described as secondary bars. Because orbital time-scales in the galaxy centre are short, secondary bars are likely to be dynamically decoupled from the main kiloparsec-scale bars. Here we show that regular orbits exist in such doubly barred potentials, and that they can support the bars in their motion. We find orbits in which particles remain on loops : closed curves which return to their original positions after two bars have come back to the same relative orientation. Stars trapped around stable loops could form the building blocks for a long-lived, doubly barred galaxy. Using the loop representation, we can find which orbits support the bars in their motion, and the constraints on the sizes and shapes of self-consistent double bars. In particular, it appears that a long-lived secondary bar may exist only when an inner Lindblad resonance is present in the primary bar, and that it would not extend beyond this resonance.     

Final stages of N-body star cluster encounters

We performed numerical simulations of star cluster encounters with Hernquist's treecode on a CRAY YMP-2E computer. We used different initial conditions (relative positions and velocities, cluster sizes, masses and concentration degrees) with the total number of particles per simulation ranging from 4608 to 20 480. Long-term interaction stages (up to 1 Gyr) when the pair coalesces into a single cluster are compared with isolated LMC clusters. Evidence is found that, when seen in a favourable plane, these resulting clusters show elliptical shapes as a result of the disruption of one of the companions. These elliptical shapes are essentially time-independent, but they do depend on the initial structural parameters of the pair components. We also analysed the fraction of stars that are ejected to the field by the interaction. We found that this fraction can be almost 50 per cent for the disrupted cluster. These simulations can represent a possible mechanism with which to explain the ellipticity observed in several star clusters in the Magellanic Clouds.     

Dynamical stability and evolution of the discs of Sc galaxies

We examine the local stability of galactic discs against axisymmetric density perturbations with special attention to the different dynamics of the stellar and gaseous components. In particular, the discs of the Milky Way and of NGC 6946 are studied. The Milky Way is shown to be stable, whereas the inner parts of NGC 6946, a typical Sc galaxy from the Kennicutt sample, are dynamically unstable. The ensuing dynamical evolution of the composite disc is studied by numerical simulations. The evolution is so fierce that the stellar disc heats up dynamically on a short time-scale to such a high degree, which seems to contradict the morphological appearance of the galaxy. The star formation rate required to cool the disc dynamically is estimated. Even if the star formation rate in NGC 6946 is at present high enough to meet this requirement, it is argued that the discs of Sc galaxies cannot sustain such a high star formation rate for extended periods.     

Self-consistent response of a galactic disc to vertical perturbations

We study the self-consistent, linear response of a galactic disc to vertical perturbations, as induced, say, by a tidal interaction. We calculate the self-gravitational potential corresponding to a non-axisymmetric, self-consistent density response of the disc using the Green's function method. The response potential is shown to oppose the perturbation potential because the self-gravity of the disc resists the imposed potential, and this resistance is stronger in the inner parts of a galactic disc. For the   m = 1  azimuthal wavenumber, the disc response opposes the imposed perturbation up to a radius that spans a range of 4–6 disc scalelengths, so that the disc shows a net warp only beyond this region. This physically explains the well known but so far unexplained observation that warps typically set in beyond this range of radii. We show that the inclusion of a dark matter halo in the calculation only marginally changes (by ∼10 per cent) the radius for the onset of warps. For perturbations with higher azimuthal wavenumbers, the net signature of the vertical perturbations can only be seen at larger radii – for example, beyond 7 exponential disc scalelengths for   m = 10  . Also, for the high- m cases, the magnitude of the negative disc response due to the disc self-gravity is much smaller. This is shown to result in corrugations of the mid-plane density, which explains the puzzling scalloping with   m = 10  detected in H  i in the outermost regions ∼30 kpc in the Galaxy.     

Stability of power-law discs — II. The global spiral modes

This paper reports on the in-plane normal modes in the self-consistent and the cut-out power-law discs. Although the cut-out discs are remarkably stable to bisymmetric perturbations, they are very susceptible to one-armed modes. For this harmonic, there is no inner Lindblad resonance, thus removing a powerful stabilizing influence. A physical mechanism for the generation of the one-armed instabilities is put forward. Incoming trailing waves are reflected as leading waves at the inner cut-out, thus completing the feedback for the swing-amplifier. Growing three-armed and four-armed modes occur only at very low temperatures. However, neutral m  = 3 and m  = 4 modes are possible at higher temperatures for some discs. The rotation curve index β has a marked effect on stability. For all azimuthal wavenumbers, any unstable modes persist to higher temperatures and grow more vigorously if the rotation curve is rising (β < 0) than if the rotation curve is falling (β > 0). If the central regions or outer parts of the disc are carved out more abruptly, any instabilities become more virulent. The self-consistent power-law discs possess a number of unusual stability properties. There is no natural time-scale in the self-consistent disc. If a mode is admitted at some pattern speed and growth rate, then it must be present at all pattern speeds and growth rates. Our analysis — although falling short of a complete proof — suggests that such a two-dimensional continuum of non-axisymmetric modes does not occur and that the self-consistent power-law discs admit no global non-axisymmetric normal modes whatsoever. Without reflecting boundaries or cut-outs, there is no resonant cavity and no possibility of unstable growing modes. The self-consistent power-law discs certainly admit equi-angular spirals as neutral modes, together with a one-dimensional continuum of growing axisymmetric modes.     

Triaxial haloes and particle dark matter detection

This paper presents the properties of a family of scale-free triaxial haloes. We adduce arguments to suggest that the velocity ellipsoids of such models are aligned in conical coordinates. We provide an algorithm to find the set of conically aligned velocity second moments that support a given density against the gravity field of the halo. The case of the logarithmic ellipsoidal model – the simplest triaxial generalization of the familiar isothermal sphere – is examined in detail. The velocity dispersions required to hold up the self-consistent model are analytic. The velocity distribution of the dark matter can be approximated as a triaxial Gaussian with semiaxes equal to the velocity dispersions.
There are roughly 20 experiments worldwide that are searching for evidence of scarce interactions between weakly interacting massive-particle dark matter (WIMP) and detector nuclei. The annual modulation signal, caused by the Earth's rotation around the Sun, is a crucial discriminant between WIMP events and the background. The greatest rate is in June, the least in December. We compute the differential detection rate for energy deposited by the rare WIMP–nucleus interactions in our logarithmic ellipsoidal halo models. Triaxiality and velocity anisotropy change the total rate by up to ∼40 per cent, and have a substantial effect on the amplitude of the annual modulation signal. The overall rate is greatest, but the amplitude of the modulation is weakest, in our radially anisotropic halo models. Even the sign of the signal can be changed. Restricting attention to low energy events, the models predict that the maximum rate occurs in December, and not in June.     

Vertical distribution of Galactic disc stars and gas constrained by a molecular cloud complex

We investigate the dynamical effects of a molecular cloud complex with a mass ∼ 10⁷ M_⊙ and a size ∼ a few 100 pc on the vertical distribution of stars and atomic hydrogen gas in a spiral galactic disc. Such massive complexes have now been observed in a number of spiral galaxies. The extended mass distribution in a complex, with an average mass density 6 times higher than the Oort limit, is shown to dominate the local gravitational field. This results in a significant redistribution or clustering of the surrounding disc components towards the mid-plane, with a resulting decrease in their vertical scaleheights.
The modified, self-consistent stellar density distribution is obtained by solving the combined Poisson equation and the force equation along the z -direction for an isothermal stellar disc on which the complex is imposed. The effect of the complex is strongest at its centre, where the stellar mid-plane density increases by a factor of 2.6 and the vertical scaleheight decreases by a factor of 3.4 compared with the undisturbed stellar disc. A surprising result is the large radial distance of ∼ 500 pc from the complex centre over which the complex influences the disc; this is due to the extended mass distribution in a complex. The complex has a comparable effect on the vertical distribution of the atomic hydrogen gas in the galactic disc. This 'pinching' or constraining effect should be detectable in the nearby spiral galaxies, as for example has been done for NGC 2403 by Sicking. Thus the gravitational field of a complex results in local corrugations of the stellar and H  i vertical scaleheights, and the galactic disc potential is highly non-uniform on scales of the intercomplex separation of ∼ 1 kpc.     

Dynamical evolution of galactic disks driven by interaction with a satellite

Dynamical evolution of galactic disks driven by interaction with satellite galaxies, particularly the problem of the disk warping and thickening is studied numerically. One of the main purpose of the study is to resolve the long standing problem of the origin of the disk warping. A possible cause of the warp is interaction with a satellite galaxy. In the case of the Milky Way, the LMC has been considered as the candidate. Some linear analysis have already given a positive result, but one had to wait for a fully self-consistent simulation as a proof. I have accomplished the numerical simulations with a million particles, by introducing a hybrid algorithm, SCF-TREE. Those simulations give us quantitative estimates for the Milky Way system. We have found an example in which large warp amplitudes are developed. We also found that the warp amplitudes depend on the halo distribution. Among our three models, the most massive and spherical halo is preferable for the observable warp excitation. This revised version was published online in August 2006 with corrections to the Cover Date.     

Invariant distributions and collisionless equilibria

This paper discusses the possibility of constructing time-independent solutions to the collisionless Boltzmann equation which depend on quantities other than global isolating integrals such as energy and angular momentum. The key point is that, at least in principle, a self-consistent equilibrium can be constructed from any set of time-independent phase-space building blocks which, when combined, generate the mass distribution associated with an assumed time-independent potential. This approach provides a way to justify Schwarzschild's method for the numerical construction of self-consistent equilibria with arbitrary time-independent potentials, generalizing thereby an approach developed by Vandervoort for integrable potentials. As a simple illustration, Schwarzschild's method is reformulated to allow for a straightforward computation of equilibria which depend only on one or two global integrals and no other quantities, as is reasonable, for example, for modelling axisymmetric configurations characterized by a non-integrable potential.     

The onset of warps in Spitzer observations of edge-on spiral galaxies

We analyse warps in the nearby edge-on spiral galaxies observed in the Spitzer /Infrared Array Camera (IRAC) 4.5-μm band. In our sample of 24 galaxies, we find evidence of warp in 14 galaxies. We estimate the observed onset radii for the warps in a subsample of 10 galaxies. The dark matter distribution in each of these galaxies are calculated using the mass distribution derived from the observed light distribution and the observed rotation curves. The theoretical predictions of the onset radii for the warps are then derived by applying a self-consistent linear response theory to the obtained mass models for six galaxies with rotation curves in the literature. By comparing the observed onset radii to the theoretical ones, we find that discs with constant thickness can not explain the observations; moderately flaring discs are needed. The required flaring is consistent with the observations. Our analysis shows that the onset of warp is not symmetric in our sample of galaxies. We define a new quantity called the onset-asymmetry index and study its dependence on galaxy properties. The onset asymmetries in warps tend to be larger in galaxies with smaller disc scalelengths. We also define and quantify the global asymmetry in the stellar light distribution, that we call the edge-on asymmetry in edge-on galaxies. It is shown that in most cases the onset asymmetry in warp is actually anticorrelated with the measured edge-on asymmetry in our sample of edge-on galaxies and this could plausibly indicate that the surrounding dark matter distribution is asymmetric.     

Bar-induced evolution of dark matter cusps

The evolution of a stellar bar transforms not only the galactic disc, but also the host dark matter halo. We present high-resolution, fully self-consistent N -body simulations that clearly demonstrate that dark matter halo central density cusps flatten as the bar torques the halo. This effect is independent of the bar formation mode and occurs even for rather short bars. The halo and bar evolution is mediated by resonant interactions between orbits in the halo and the bar pattern speed, as predicted by linear Hamiltonian perturbation theory. The bar lengthens and slows as it loses angular momentum, a process that occurs even in rather warm discs. We demonstrate that the bar and halo response can be critically underestimated for experiments that are unable to resolve the relevant resonant dynamics; this occurs when the phase space in the resonant region is undersampled or plagued by noise.     

Global m= 1 instabilities and lopsidedness in disc galaxies

Lopsidedness is common in spiral galaxies. Often, there is no obvious external cause, such as an interaction with a nearby galaxy, for such features. Alternatively, the lopsidedness may have an internal cause, such as a dynamical instability. In order to explore this idea, we have developed a computer code that searches for self-consistent perturbations in razor-thin disc galaxies and performed a thorough mode-analysis of a suite of dynamical models for disc galaxies embedded in an inert dark matter halo with varying amounts of rotation and radial anisotropy.
Models with two equal-mass counter-rotating discs and fully rotating models both show growing lopsided modes. For the counter-rotating models, this is the well-known counter-rotating instability, becoming weaker as the net rotation increases. The m = 1 mode of the maximally rotating models, on the other hand, becomes stronger with increasing net rotation. This rotating m = 1 mode is reminiscent of the eccentricity instability in near-Keplerian discs.
To unravel the physical origin of these two different m = 1 instabilities, we studied the individual stellar orbits in the perturbed potential and found that the presence of the perturbation gives rise to a very rich orbital behaviour. In the linear regime, both instabilities are supported by aligned loop orbits. In the non-linear regime, other orbit families exist that can help support the modes. In terms of density waves, the counter-rotating m = 1 mode is due to a purely growing Jeans-type instability. The rotating m = 1 mode, on the other hand, grows as a result of the swing amplifier working inside the resonance cavity that extends from the disc centre out to the radius where non-rotating waves are stabilized by the model's outwardly rising Q profile.     

An MHD gadget for cosmological simulations

Various radio observations have shown that the hot atmospheres of galaxy clusters are magnetized. However, our understanding of the origin of these magnetic fields, their implications on structure formation and their interplay with the dynamics of the cluster atmosphere, especially in the centres of galaxy clusters, is still very limited. In preparation for the upcoming new generation of radio telescopes (like Expanded Very Large Array, Low Wavelength Array, Low Frequency Array and Square Kilometer Array), a huge effort is being made to learn more about cosmological magnetic fields from the observational perspective. Here we present the implementation of magnetohydrodynamics (MHD) in the cosmological smoothed particle hydrodynamics (SPH) code gadget . We discuss the details of the implementation and various schemes to suppress numerical instabilities as well as regularization schemes, in the context of cosmological simulations. The performance of the SPH–MHD code is demonstrated in various one- and two-dimensional test problems, which we performed with a fully, three-dimensional set-up to test the code under realistic circumstances. Comparing solutions obtained using athena , we find excellent agreement with our SPH–MHD implementation. Finally, we apply our SPH–MHD implementation to galaxy cluster formation within a large, cosmological box. Performing a resolution study we demonstrate the robustness of the predicted shape of the magnetic field profiles in galaxy clusters, which is in good agreement with previous studies.     

Gas-driven evolution of stellar orbits in barred galaxies

We carry out a detailed orbit analysis of gravitational potentials selected at different times from an evolving self-consistent model galaxy consisting of a two-component disc (stars+gas) and a live halo. The results are compared with a pure stellar model, subject to nearly identical initial conditions, which are chosen so as to make the models develop a large-scale stellar bar. The bars are also subject to hose-pipe (buckling) instability which modifies the vertical structure of the disc. The diverging morphological evolution of both models is explained in terms of gas radial inflow, the resulting change in the gravitational potential at smaller radii, and the subsequent modification of the main families of orbits, both in and out of the disc plane.   We find that dynamical instabilities become milder in the presence of the gas component, and that the stability of planar and 3D stellar orbits is strongly affected by the related changes in the potential — both are destabilized, with the gas accumulation at the centre. This is reflected in the overall lower amplitude of the bar mode and in the substantial weakening of the bar, which appears to be a gradual process. The vertical buckling of the bar is much less pronounced and the characteristic peanut shape of the galactic bulge almost disappears when there is a substantial gas inflow towards the centre. Milder instability results in a smaller bulge, the basic parameters of which are in agreement with observations. We also find that the overall evolution in the model with a gas component is accelerated because of the larger central mass concentration and the resulting decrease in the characteristic dynamical time.     

The impact of accretion disc winds on the X-ray spectrum of AGN – I. xscort

The accretion disc in active galactic nucleus (AGN) is expected to produce strong outflows, in particular an ultraviolet (UV)-line-driven wind. Several observed spectral features, including the soft X-ray excess, have been associated with the accretion disc wind. However, current spectral models of the X-ray spectrum of AGN observed through an accretion disc wind, known to provide a good fit to the observed X-ray data, are ad hoc in their treatment of the outflow velocity and density of the wind material. In order to address these limitations we adopt a numerical computational method that links a series of radiative transfer calculations, incorporating the effect of a global velocity field in a self-consistent manner { xstar Simulation Chain for Outflows with Radiative Transfer ( xscort )}. We present a series of example spectra from the xscort code that allow us to examine the shape of AGN X-ray spectra seen through a smooth wind with terminal velocity of 0.3 c , as appropriate for a UV-line-driven wind. We calculate spectra for a range of different acceleration laws, density distributions, total column densities and ionization parameters, but all these have sharp features that contrast strongly with both the previous 'smeared absorption' models, and with the observed smoothness of the soft X-ray excess. This rules out absorption in a radiatively driven accretion disc wind as the origin of the soft X-ray excess, though a larger terminal velocity, possibly associated with material in a magnetically driven outflow/jet, may allow outflow models to recover a smooth excess.     

Modelling gaseous and stellar kinematics in the disc galaxies NGC 772, 3898 and 7782

We present V -band surface photometry and major-axis kinematics of stars and ionized gas of three early-type spiral galaxies, namely NGC 772, 3898 and 7782. For each galaxy we present a self-consistent Jeans model for the stellar kinematics, adopting the light distribution of bulge and disc derived by means of a two-dimensional parametric photometric decomposition. This allows us to investigate the presence of non-circular gas motions, and derive the mass distribution of luminous and dark matter in these objects.
NGC 772 and 7782 have apparently normal kinematics with the ionized gas tracing the gravitational equilibrium circular speed. This is not true in the innermost region (| r |≲8 arcsec) of NGC 3898, where the ionized gas is rotating more slowly than the circular velocity predicted by dynamical modelling. This phenomenon is common in the bulge-dominated galaxies for which dynamical modelling enables us to make the direct comparison between the gas velocity and the circular speed, and it poses questions about the reliability of galaxy mass distributions derived by the direct decomposition of the observed ionized-gas rotation curve into the contributions of luminous and dark matter.     

Multimass schemes for collisionless N-body simulations

We present a general scheme for constructing Monte Carlo realizations of equilibrium, collisionless galaxy models with known distribution function (DF) f ₀. Our method uses importance sampling to find the sampling DF f _s that minimizes the mean-square formal errors in a given set of projections of the DF f ₀. The result is a multimass N -body realization of the galaxy model in which 'interesting' regions of phase space are densely populated by lots of low-mass particles, increasing the effective N there, and less interesting regions by fewer, higher mass particles.
As a simple application, we consider the case of minimizing the shot noise in estimates of the acceleration field for an N -body model of a spherical Hernquist model. Models constructed using our scheme easily yield a factor of ∼100 reduction in the variance at the central acceleration field when compared to a traditional equal-mass model with the same number of particles. When evolving both models with a real N -body code, the diffusion coefficients in our model are reduced by a similar factor. Therefore, for certain types of problems, our scheme is a practical method for reducing the two-body relaxation effects, thereby bringing the N -body simulations closer to the collisionless ideal.     

Density–potential pairs for spherical stellar systems with Sérsic light profiles and (optional) power-law cores

Popular models for describing the luminosity-density profiles of dynamically hot stellar systems (e.g. Jaffe, Hernquist, Dehnen) were constructed with the desire to match the deprojected form of an   R ^1/4  light profile. Real galaxies, however, are now known to have a range of different light-profile shapes that scale with mass. Consequently, although highly useful, the above models have implicit limitations, and this is illustrated here through their application to a number of real galaxy density profiles. On the other hand, the analytical density profile given by Prugniel & Simien closely matches the deprojected form of Sérsic   R ^{1/ n}   light profiles – including deprojected exponential light profiles. It is thus applicable for describing bulges in spiral galaxies, dwarf elliptical galaxies, and both ordinary and giant elliptical galaxies. Moreover, the observed Sérsic quantities define the parameters of the density model. Here we provide simple equations, in terms of elementary and special functions, for the gravitational potential and force associated with this density profile. Furthermore, to match galaxies with partially depleted cores, and better explore the supermassive black hole/galaxy connection, we have added a power-law core to this density profile and derived similar expressions for the potential and force of this hybrid profile. Expressions for the mass and velocity dispersion, assuming isotropy, are also given. These spherical models may also prove appropriate for describing the dark matter distribution in haloes formed from ΛCDM cosmological simulations.