6th and 8th Order Hermite Integrator for N-body Simulations

We present sixth- and eighth-order Hermite integrators for astrophysical $N$-body simulations, which use the derivatives of accelerations up to second order ({\it snap}) and third order ({\it crackle}). These schemes do not require previous values for the corrector, and require only one previous value to construct the predictor. Thus, they are fairly easy to implemente. The additional cost of the calculation of the higher order derivatives is not very high. Even for the eighth-order scheme, the number of floating-point operations for force calculation is only about two times larger than that for traditional fourth-order Hermite scheme. The sixth order scheme is better than the traditional fourth order scheme for most cases. When the required accuracy is very high, the eighth-order one is the best. These high-order schemes have several practical advantages. For example, they allow a larger number of particles to be integrated in parallel than the fourth-order scheme does, resulting in higher execution efficiency in both general-purpose parallel computers and GRAPE systems.Comment: 21 pages, 6 figures, New Astronomy accepte

Majority-Carrier Mobilities in Undoped and \textit{n}-type Doped ZnO Epitaxial Layers

Transparent and conductive ZnO:Ga thin films are prepared by laser molecular-beam epitaxy. Their electron properties were investigated by the temperature-dependent Hall-effect technique. The 300-K carrier concentration and mobility were about $n_s \sim 10^{16}$ cm$^{-3}$ and 440 cm$^{2}$/Vs, respectively. In the experimental `mobility vs concentration' curve, unusual phenomenon was observed, i.e., mobilities at $n_s \sim 5\times$ 10$^{18}$ cm$^{-3}$ are significantly smaller than those at higher densities above $\sim 10^{20}$ cm$^{-3}$. Several types of scattering centers including ionized donors and oxygen traps are considered to account for the observed dependence of the Hall mobility on carrier concentration. The scattering mechanism is explained in terms of inter-grain potential barriers and charged impurities. A comparison between theoretical results and experimental data is made.Comment: 5 pages, 1 figure, conference on II-VI compounds, RevTe

N-body simulation for self-gravitating collisional systems with a new SIMD instruction set extension to the x86 architecture, Advanced Vector eXtensions

We present a high-performance N-body code for self-gravitating collisional systems accelerated with the aid of a new SIMD instruction set extension of the x86 architecture: Advanced Vector eXtensions (AVX), an enhanced version of the Streaming SIMD Extensions (SSE). With one processor core of Intel Core i7-2600 processor (8 MB cache and 3.40 GHz) based on Sandy Bridge micro-architecture, we implemented a fourth-order Hermite scheme with individual timestep scheme (Makino and Aarseth, 1992), and achieved the performance of 20 giga floating point number operations per second (GFLOPS) for double-precision accuracy, which is two times and five times higher than that of the previously developed code implemented with the SSE instructions (Nitadori et al., 2006b), and that of a code implemented without any explicit use of SIMD instructions with the same processor core, respectively. We have parallelized the code by using so-called NINJA scheme (Nitadori et al., 2006a), and achieved 90 GFLOPS for a system containing more than N = 8192 particles with 8 MPI processes on four cores. We expect to achieve about 10 tera FLOPS (TFLOPS) for a self-gravitating collisional system with N 105 on massively parallel systems with at most 800 cores with Sandy Bridge micro-architecture. This performance will be comparable to that of Graphic Processing Unit (GPU) cluster systems, such as the one with about 200 Tesla C1070 GPUs (Spurzem et al., 2010). This paper offers an alternative to collisional N-body simulations with GRAPEs and GPUs.Comment: 14 pages, 9 figures, 3 tables, accepted for publication in New Astronomy. The code is publicly available at http://code.google.com/p/phantom-grape

Performance analysis of direct N-body algorithms for astrophysical simulations on distributed systems

We discuss the performance of direct summation codes used in the simulation of astrophysical stellar systems on highly distributed architectures. These codes compute the gravitational interaction among stars in an exact way and have an O(N^2) scaling with the number of particles. They can be applied to a variety of astrophysical problems, like the evolution of star clusters, the dynamics of black holes, the formation of planetary systems, and cosmological simulations. The simulation of realistic star clusters with sufficiently high accuracy cannot be performed on a single workstation but may be possible on parallel computers or grids. We have implemented two parallel schemes for a direct N-body code and we study their performance on general purpose parallel computers and large computational grids. We present the results of timing analyzes conducted on the different architectures and compare them with the predictions from theoretical models. We conclude that the simulation of star clusters with up to a million particles will be possible on large distributed computers in the next decade. Simulating entire galaxies however will in addition require new hybrid methods to speedup the calculation.Comment: 22 pages, 8 figures, accepted for publication in Parallel Computin

On the Origin of Density Cusps in Elliptical Galaxies

We investigated the dynamical reaction of the central region of galaxies to a falling massive black hole by N-body simulations. As the initial galaxy model, we used an isothermal King model and placed a massive black hole at around the half-mass radius of the galaxy. We found that the central core of the galaxy is destroyed by the heating due to the black hole and that a very weak density cusp ($\rho \propto r^{-\alpha}$, with $\alpha \sim 0.5$) is formed around the black hole. This result is consistent with recent observations of large elliptical galaxies with Hubble Space Telescope. The velocity of the stars becomes tangentially anisotropic in the inner region, while in the outer region the stars have radially anisotropic velocity dispersion. The radius of the weak cusp region is larger for larger black hole mass. Our result naturally explains the formation of the weak cusp found in the previous simulations of galaxy merging, and implies that the weak cusp observed in large elliptical galaxies may be formed by the heating process by sinking black holes during merging events.Comment: 14 pages with 29 EPS figures; LaTeX; new results added; accepted for publication in Ap

Evolution of Star Clusters near the Galactic Center: Fully Self-consistent N-body Simulations

We have performed fully self-consistent $N$-body simulations of star clusters near the Galactic center (GC). Such simulations have not been performed because it is difficult to perform fast and accurate simulations of such systems using conventional methods. We used the Bridge code, which integrates the parent galaxy using the tree algorithm and the star cluster using the fourth-order Hermite scheme with individual timestep. The interaction between the parent galaxy and the star cluster is calculate with the tree algorithm. Therefore, the Bridge code can handle both the orbital and internal evolutions of star clusters correctly at the same time. We investigated the evolution of star clusters using the Bridge code and compared the results with previous studies. We found that 1) the inspiral timescale of the star clusters is shorter than that obtained with "traditional" simulations, in which the orbital evolution of star clusters is calculated analytically using the dynamical friction formula and 2) the core collapse of the star cluster increases the core density and help the cluster survive. The initial conditions of star clusters is not so severe as previously suggested.Comment: 19 pages, 19 figures, accepted for publication in Ap

Formation of Protoplanets from Massive Planetesimals in Binary Systems

More than half of stars reside in binary or multiple star systems and many planets have been found in binary systems. From theoretical point of view, however, whether or not the planetary formation proceeds in a binary system is a very complex problem, because secular perturbation from the companion star can easily stir up the eccentricity of the planetesimals and cause high-velocity, destructive collisions between planetesimals. Early stage of planetary formation process in binary systems has been studied by restricted three-body approach with gas drag and it is commonly accepted that accretion of planetesimals can proceed due to orbital phasing by gas drag. However, the gas drag becomes less effective as the planetesimals become massive. Therefore it is still uncertain whether the collision velocity remains small and planetary accretion can proceed, once the planetesimals become massive. We performed {\it N}-body simulations of planetary formation in binary systems starting from massive planetesimals whose size is about 100-500 km. We found that the eccentricity vectors of planetesimals quickly converge to the forced eccentricity due to the coupling of the perturbation of the companion and the mutual interaction of planetesimals if the initial disk model is sufficiently wide in radial distribution. This convergence decreases the collision velocity and as a result accretion can proceed much in the same way as in isolated systems. The basic processes of the planetary formation, such as runaway growth and oligarchic growth and final configuration of the protoplanets are essentially the same in binary systems and single star systems, at least in the late stage where the effect of gas drag is small.Comment: 26pages, 11 figures. ApJ accepte

Evolution of Massive Black Hole Binaries

We present the result of large-scale N-body simulations of the stellar-dynamical evolution of a massive black-hole binary at the center of a spherical galaxy. We focus on the dependence of the hardening rate on the relaxation timescale of the parent galaxy. A simple theoretical argument predicts that a binary black hole creates the ``loss cone'' around it. Once the loss cone is formed, the hardening rate is determined by the rate at which field stars diffuse into the loss cone. Therefore the hardening timescale becomes proportional to the relaxation timescale. Recent N-body simulations, however, have failed to confirm this theory and various explanations have been proposed. By performing simulations with sufficiently large N (up to $10^6$) for sufficiently long time, we found that the hardening rate does depend on N. Our result is consistent with the simple theoretical prediction that the hardening timescale is proportional to the relaxation timescale. This dependence implies that most massive black hole binaries are unlikely to merge within the Hubble time through interaction with field stars and gravitational wave radiation alone.Comment: Reviced version accepted for publication in ApJ. Scheduled to appear in the February 10, 2004 issu

Post-Collapse evolution of globular clusters

A number of globular clusters appear to have undergone core collapse, in the sense that their predicted collapse time is much shorter than their current age. Simulations using gas models and Fokker-Planck approximation have shown that the central density of a globular cluster after the collapse undergoes nonlinear oscillation with large amplitude (gravothermal oscillation). However, whether such an oscillation actually takes place in a real $N$-body system has remained unsolved, because an $N$-body simulation with a sufficiently high resolution would have required the computing resource of the order of several Gflops$\cdot$years. In the present paper, we report the result of such a simulation, performed on a dedicated special-purpose computer GRAPE-4. We simulated the evolution of isolated point-mass systems with up to 32,768 particles. The largest number of particles reported previously is 10,000. We confirmed that gravothermal oscillation takes place in an $N$-body system. The expansion phase shows all signatures that are considered as the evidences of the gravothermal nature of the oscillation. At the maximum expansion, the core radius is $\sim 1$\% of the half-mass radius for the run with 32,768 particles. The maximum core size $r_c$ depends on $N$, as $\propto N^{-1/3}$.Comment: To appear in Apj, Vol 470, Oct 20, 11 pages, 7 figures, Postscript (or GIF) version available at http://grape.c.u-tokyo.ac.jp/pub/people/makino/papers/osclong.ps or http://grape.c.u-tokyo.ac.jp/pub/people/makino/papers/osclong.htm