We present an implementation of the hierarchical tree algorithm on the individual timestep algorithm (the Hermite scheme) for collisional $N$-body simulations, running on GRAPE-9 system, a special-purpose hardware accelerator for gravitational many-body simulations. Such combination of the tree algorithm and the individual timestep algorithm was not easy on the previous GRAPE system mainly because its memory addressing scheme was limited only to sequential access to a full set of particle data. The present GRAPE-9 system has an indirect memory addressing unit and a particle memory large enough to store all particles data and also tree nodes data. The indirect memory addressing unit stores interaction lists for the tree algorithm, which is constructed on host computer, and, according to the interaction lists, force pipelines calculate only the interactions necessary. In our implementation, the interaction calculations are significantly reduced compared to direct $N^2$ summation in the original Hermite scheme. For example, we can archive about a factor 30 of speedup (equivalent to about 17 teraflops) against the Hermite scheme for a simulation of $N=10^6$ system, using hardware of a peak speed of 0.6 teraflops for the Hermite scheme.

We investigate effects of hardness of primordial binaries on whole evolution of star clusters by means of N-body simulations. Using newly developed code, GORILLA, we simulated eleven N=16384 clusters with primordial binaries whose binding energies are equal in each cluster in range of 1-300kT_0, where 1.5kT_0 is average stellar kinetic energy at the initial time. We found that, in both soft (< 3kT_0) and hard (> 300kT_0) limits, clusters experience deep core collapse. In the intermediate hardness (10-100kT_0), the core collapses halt halfway due to an energy releases of the primordial binaries. The core radii at the halt can be explained by their energy

We investigate the long-term dynamical evolution of star clusters in a steady tidal field produced by its parent galaxy. In this paper, we focus on the influence of mass profile of the parent galaxy. The previous studies were done with the simplification where the parent galaxy was expressed by point mass. We express different mass profiles of the parent galaxy by the tidal fields in which the ratios of the epicyclic frequency to the angular velocity are different. We compare the mass-loss timescale of star clusters whose tidal radii are identical but in parent galaxies with different mass profile, by means of orbits calculations in fixed cluster potential and N-body simulations. In this situation, a cluster rotates around the parent galaxy more rapidly as the parent galaxy has shallower mass profile. We found that the mass-loss timescale increase 20% and 50% for the cases that the mass density profile of the parent galaxies are proportional to R^-2 and R^-1.5 where R is the distance from the galaxy center, compared to the point-mass case, in moderately strong tidal field. Counterintuitively, a cluster which rotates around the parent galaxy more rapidly has a longer lifetime. The increase of lifetime is due to the fact that the fraction occupied by regular-like orbit increases in shallower profile. Finally, we derive an evaluation formula for the mass-loss timescale of clusters. Our formula can explain a property of the population of the observed galactic globular clusters that their half-mass radii become smaller as their distances from the galactic center become smaller.

In this paper, we describe the implementation and performance of GreeM, a massively parallel TreePM code for large-scale cosmological N-body simulations. GreeM uses a recursive multi-section algorithm for domain decomposition. The size of the domains are adjusted so that the total calculation time of the force becomes the same for all processes. The loss of performance due to non-optimal load balancing is around 4%, even for more than 10^3 CPU cores. GreeM runs efficiently on PC clusters and massively-parallel computers such as a Cray XT4. The measured calculation speed on Cray XT4 is 5 \times 10^4 particles per second per CPU core, for the case of an opening angle of \theta=0.5, if the number of particles per CPU core is larger than 10^6.

We analyzed the statistics of subhalo abundance of galaxy-sized and giant-galaxy-sized halos formed in a high-resolution cosmological simulation of a 46.5Mpc cube with the uniform mass resolution of $10^6 M_{\odot}$. We analyzed all halos with mass more than $1.5 \times 10^{12}M_{\odot}$ formed in this simulation box. The total number of halos was 125. We found that the subhalo abundance, measured by the number of subhalos with maximum rotation velocity larger than 10% of that of the parent halo, shows large halo-to-halo variations. The results of recent ultra-high-resolution runs fall within the variation of our samples. We found that the concentration parameter and the radius at the moment of the maximum expansion shows fairly tight correlation with the subhalo abundance. This correlation suggests that the variation of the subhalo abundance is at least partly due to the difference in the formation history. Halos formed earlier have smaller number of subhalos at present.

Recent high-resolution simulations of the formation of dark-matter halos have shown that the distribution of subhalos is scale-free, in the sense that if scaled by the velocity dispersion of the parent halo, the velocity distribution function of galaxy-sized and cluster-sized halos are identical. For cluster-sized halos, simulation results agreed well with observations. Simulations, however, predicted far too many subhalos for galaxy-sized halos. Our galaxy has several tens of known dwarf galaxies. On the other hands, simulated dark-matter halos contain thousands of subhalos. We have performed simulation of a single large volume and measured the abundance of subhalos in all massive halos. We found that the variation of the subhalo abundance is very large, and those with largest number of subhalos correspond to simulated halos in previous studies. The subhalo abundance depends strongly on the local density of the background. Halos in high-density regions contain large number of subhalos. Our galaxy is in the low-density region. For our simulated halos in low-density regions, the number of subhalos is within a factor of three to that of our galaxy. We argue that the ``missing dwarf problem'' is not a real problem but caused by the biased selection of the initial conditions in previous studies, which were not appropriate for field galaxies.

We have developed PGPG (Pipeline Generator for Programmable GRAPE), a software which generates the low-level design of the pipeline processor and communication software for FPGA-based computing engines (FBCEs). An FBCE typically consists of one or multiple FPGA (Field-Programmable Gate Array) chips and local memory. Here, the term "Field-Programmable" means that one can rewrite the logic implemented to the chip after the hardware is completed, and therefore a single FBCE can be used for calculation of various functions, for example pipeline processors for gravity, SPH interaction, or image processing. The main problem with FBCEs is that the user need to develop the detailed hardware design for the processor to be implemented to FPGA chips. In addition, she or he has to write the control logic for the processor, communication and data conversion library on the host processor, and application program which uses the developed processor. These require detailed knowledge of hardware design, a hardware description language such as VHDL, the operating system and the application, and amount of human work is huge. A relatively simple design would require 1 person-year or more. The PGPG software generates all necessary design descriptions, except for the application software itself, from a high-level design description of the pipeline processor in the PGPG language. The PGPG language is a simple language, specialized to the description of pipeline processors. Thus, the design of pipeline processor in PGPG language is much easier than the traditional design. For real applications such as the pipeline for gravitational interaction, the pipeline processor generated by PGPG achieved the performance similar to that of hand-written code. In this paper we present a detailed description of PGPG version 1.0.

We present a novel approach to accelerate astrophysical hydrodynamical simulations. In astrophysical many-body simulations, GRAPE (GRAvity piPE) system has been widely used by many researchers. However, in the GRAPE systems, its function is completely fixed because specially developed LSI is used as a computing engine. Instead of using such LSI, we are developing a special purpose computing system using Field Programmable Gate Array (FPGA) chips as the computing engine. Together with our developed programming system, we have implemented computing pipelines for the Smoothed Particle Hydrodynamics (SPH) method on our PROGRAPE-3 system. The SPH pipelines running on PROGRAPE-3 system have the peak speed of 85 GFLOPS and in a realistic setup, the SPH calculation using one PROGRAPE-3 board is 5-10 times faster than the calculation on the host computer. Our results clearly shows for the first time that we can accelerate the speed of the SPH simulations of a simple astrophysical phenomena using considerable computing power offered by the hardware.

We present Particle-Particle-Particle-Mesh (PPPM) and Tree Particle-Mesh (TreePM) implementations on GRAPE-5 and GRAPE-6A systems, special-purpose hardware accelerators for gravitational many-body simulations. In our PPPM and TreePM implementations on GRAPE, the computational time is significantly reduced compared with the conventional implementations without GRAPE, especially under the strong particle clustering, and almost constant irrespective of the degree of particle clustering. We carry out the survey of two simulation parameters, the PM grid spacing and the opening parameter for the most optimal combination of force accuracy and computational speed. We also describe the parallelization of these implementations on a PC-GRAPE cluster, in which each node has one GRAPE board, and present the optimal configuration of simulation parameters for good parallel scalability.

In this paper, we describe the design and performance of GRAPE-6A, a special-purpose computer for gravitational many-body simulations. It was designed to be used with a PC cluster, in which each node has one GRAPE-6A. Such configuration is particularly effective in running parallel tree algorithm. Though the use of parallel tree algorithm was possible with the original GRAPE-6 hardware, it was not very cost-effective since a single GRAPE-6 board was still too fast and too expensive. Therefore, we designed GRAPE-6A as a single PCI card to minimize the reproduction cost and optimize the computing speed. The peak performance is 130 Gflops for one GRAPE-6A board and 3.1 Tflops for our 24 node cluster. We describe the implementation of the tree, TreePM and individual timestep algorithms on both a single GRAPE-6A system and GRAPE-6A cluster. Using the tree algorithm on our 16-node GRAPE-6A system, we can complete a collisionless simulation with 100 million particles (8000 steps) within 10 days.

We investigate evolution of star clusters in external tidal field by means of $N$-body simulations. We followed seven sets of cluster models whose central concentration and strength of the tidal field are different. We found that the mass loss timescale due to escape of stars, $t_{mloss}$, and its dependence on the two-body relaxation timescale, $t_{rh,i}$, are determined by the strength of the tidal field. The logarithmic slope [$= dln(t_{mloss})/dln(t_{rh,i})$] approaches to near unity for the cluster models in weaker tidal field. The timescale and the dependence are almost independent of the central concentration for clusters in the tidal field of the same strength. In our results, the scaling found by Baumgardt (2001) can be seen only in the cluster models with moderately strong tidal field.

In this paper, we describe the architecture and performance of the GRAPE-6 system, a massively-parallel special-purpose computer for astrophysical $N$-body simulations. GRAPE-6 is the successor of GRAPE-4, which was completed in 1995 and achieved the theoretical peak speed of 1.08 Tflops. As was the case with GRAPE-4, the primary application of GRAPE-6 is simulation of collisional systems, though it can be used for collisionless systems. The main differences between GRAPE-4 and GRAPE-6 are (a) The processor chip of GRAPE-6 integrates 6 force-calculation pipelines, compared to one pipeline of GRAPE-4 (which needed 3 clock cycles to calculate one interaction), (b) the clock speed is increased from 32 to 90 MHz, and (c) the total number of processor chips is increased from 1728 to 2048. These improvements resulted in the peak speed of 64 Tflops. We also discuss the design of the successor of GRAPE-6.

We investigate the structure of the dark matter halo formed in the cold dark matter scenarios by N-body simulations with parallel treecode on GRAPE cluster systems. We simulated 8 halos with the mass of $4.4\times 10^{14}M_{\odot}$ to $1.6\times 10^{15}M_{\odot}$ in the SCDM and LCDM model using up to 30 million particles. With the resolution of our simulations, the density profile is reliable down to 0.2 percent of the virial radius. Our results show that the slope of inner cusp within 1 percent virial radius is shallower than -1.5, and the radius where the shallowing starts exhibits run-to-run variation, which means the innermost profile is not universal.

Challenging stellar dynamical problems, such as the study of gravothermal oscillations in star clusters, have in the past initiated the very successful building of GRAPE special purpose computers. It is discussed, that present day tasks such as the formation and evolution of galactic nuclei with one or more massive black holes and the coupled stellar and gas dynamical processes in the formation of nuclei and star clusters, demand a new kind of hybrid architecture, using both GRAPE and a reconfigurable logics board called RACE. For such a system we have developed first implementations and floating point performance studies in the case of the SPH algorithm (smoothed particle hydrodynamics), which will be of great advantage for SPH modelling and also for direct $N$-body simulations using the more efficient Ahmad-Cohen neighbour schemes.

We investigate the structure of the dark matter halo formed in three different cold dark matter scenarios. We performed N-body simulations of formation of 13 cluster-sized halos. In all runs, density cusps proportional to $r^{-1.5}$ developed at the center. This result was independent of the cosmological models we simulated. We could not reproduce the cusp shallower than $r^{-1.5}$, which was obtained in some of previous studies. We also found that in all runs the density structure evolves in a self-similar way, even in $\Omega\ne 1$ universes. These results show that the formation of structural form is a process decoupled from a background cosmology.

We report on the detailed analysis of the mean-pairwise peculiar velocity profile in high-resolution cosmological N-body simulations ($N=8.8\times 10^6$ particles in a sphere of $50 \sim 200$Mpc radius). In particular we examine the validity and limitations of the stable condition which states that the mean {\it physical} separation of particle pairs is constant on small scales. We find a significant time-variation (irregular oscillatory behavior) of the mean-pairwise peculiar velocity in nonlinear regimes. We argue that this behavior is not due to any numerical artifact, but a natural consequence of the continuous merging processes in the hierarchical clustering universe. While such a time-variation is significant in a relatively local patch of the universe, the global average over a huge spatial volume $\simgt$(200Mpc)$^3$ does not reveal any {\it systematic} departure from the stable condition. Thus we conclude that the mean-pairwise peculiar velocity is rather unstable statistics but still satisfies the stable condition when averaged over the cosmological volume.

We investigate the structure of the dark matter halo formed in the cold dark matter scenario using $N$-body simulations. We simulated 12 halos with the mass of $6.6\times 10^{11}M_{\odot}$ to $8.0\times 10^{14}M_{\odot}$. In almost all runs, the halos have density cusps proportional to $r^{-1.5}$ developed at the center, which is consistent with the results of recent high-resolution calculations. The density structure evolves in a self-similar way, and is universal in the sense that it is independent of the halo mass and initial random realization of density fluctuation. The density profile is in good agreement with the profile proposed by Moore et al. (1999), which has central slope proportional to $r^{-1.5}$ and outer slope proportional to $r^{-3}$. The halo grows through repeated accretion of diffuse smaller halos. We argue that the cusp is understood as a convergence slope for the accretion of tidally disrupted matter.

As an entry for the 1999 Gordon Bell price/performance prize, we report an astrophysical N-body simulation performed with a treecode on GRAPE-5 (Gravity Pipe 5) system, a special-purpose computer for astrophysical N-body simulations. The GRAPE-5 system has 32 pipeline processors specialized for the gravitational force calculation. Other operations, such as tree construction, tree traverse and time integration, are performed on a general purpose workstation. The total cost for the GRAPE-5 system is 40,900 dollars. We performed a cosmological N-body simulation with 2.1 million particles, which sustained a performance of 5.92 Gflops averaged over 8.37 hours. The price per performance obtained is 7.0 dollars per Mflops.

In this paper a cluster is modelled as a smooth potential (due to the cluster stars) plus the steady tidal field of the Galaxy. In this model there is a minimum energy below which stars cannot escape. Above this energy, however, the time scale on which stars escape varies with the orbital parameters of the star (mainly its energy) in a way which we attempt to quantify, with both theoretical arguments and computer simulations. Within the limitations of the model we show that the time scale is long enough to complicate the interpretation of full N-body simulations of clusters, and that stars above the escape energy may remain bound to the cluster for about a Hubble time.

We have developed a special-purpose computer for gravitational many-body simulations, GRAPE-5. GRAPE-5 is the successor of GRAPE-3. Both consist of eight custom pipeline chips (G5 chip and GRAPE chip). The difference between GRAPE-5 and GRAPE-3 are: (1) The G5 chip contains two pipelines operating at 80 MHz, while the GRAPE chip had one at 20 MHz. Thus, the calculation speed of the G5 chip and that of GRAPE-5 board are 8 times faster than that of GRAPE chip and GRAPE-3 board. (2) The GRAPE-5 board adopted PCI bus as the interface to the host computer instead of VME of GRAPE-3, resulting in the communication speed one order of magnitude faster. (3) In addition to the pure 1/r potential, the G5 chip can calculate forces with arbitrary cutoff functions, so that it can be applied to Ewald or P^3M methods. (4) The pairwise force calculated on GRAPE-5 is about 10 times more accurate than that on GRAPE-3. On one GRAPE-5 board, one timestep of 128k-body simulation with direct summation algorithm takes 14 seconds. With Barnes-Hut tree algorithm (theta = 0.75), one timestep of 10^6-body simulation can be done in 16 seconds.

We have developed PROGRAPE-1 (PROgrammable GRAPE-1), a programmable multi-purpose computer for many-body simulations. The main difference between PROGRAPE-1 and "traditional" GRAPE systems is that the former uses FPGA (Field Programmable Gate Array) chips as the processing elements, while the latter rely on the hardwired pipeline processor specialized to gravitational interactions. Since the logic implemented in FPGA chips can be reconfigured, we can use PROGRAPE-1 to calculate not only gravitational interactions but also other forms of interactions such as van der Waals force, hydrodynamical interactions in SPH calculation and so on. PROGRAPE-1 comprises two Altera EPF10K100 FPGA chips, each of which contains nominally 100,000 gates. To evaluate the programmability and performance of PROGRAPE-1, we implemented a pipeline for gravitational interaction similar to that of GRAPE-3. One pipeline fitted into a single FPGA chip, which operated at 16 MHz clock. Thus, for gravitational interaction, PROGRAPE-1 provided the speed of 0.96 Gflops-equivalent. PROGRAPE will prove to be useful for wide-range of particle-based simulations in which the calculation cost of interactions other than gravity is high, such as the evaluation of SPH interactions.

The next decade will be an exciting time for computational physicists. After 50 years of being forced to use standardized commercial equipment, it will finally become relatively straightforward to adapt one's computing tools to one's own needs. The breakthrough that opens this new era is the now wide-spread availability of programmable chips that allow virtually every computational scientist to design his or her own special-purpose computer.

We developed a PCI interface for GRAPE systems. GRAPE(GRAvity piPE) is a special-purpose computer for gravitational N-body simulations. A GRAPE system consists of GRAPE processor boards and a host computer. GRAPE processors perform the calculation of gravitational forces between particles. The host computer performs the rest of calculations. The newest of GRAPE machines, the GRAPE-4, achieved the peak performance of 1.08 Tflops. The GRAPE-4 system uses TURBOChannel for the interface to the host, which limits the selection of the host computer. The TURBOChannel bus is not supported by any of recent workstations. We developed a new host interface board which adopts the PCI bus instead of the TURBOChannel. PCI is an I/O bus standard developed by Intel. It has fairly high peak transfer speed, and is available on wide range of computers, from PCs to supercomputers. Thus, the new interface allows us to connect GRAPE-4 to a wide variety of host computers. In test runs with a Barnes-Hut treecode, we found that the performance of new system with PCI interface is 40% better than that of the original system.

Observed cusps with density profiles $\rho\propto r^{-1}$ or shallower, in the central regions of galaxies, cannot be reproduced in the standard Cold Dark Matter (CDM) picture of hierarchical clustering. Previous claims to the contrary were based on simulations with relatively few particles, and substantial softening. We present simulations with particle numbers an order of magnitude higher, and essentially no softening, and show that typical central density profiles are clearly steeper than $\rho\propto r^{-1}$. The observed shallower profiles may have formed through the smoothing effect of the spiral-in of central black holes in previous merger phases. In addition, we confirm the presence of a temperature inversion in the inner 5 kpc of massive galactic halos, and illustrate its formation as a natural result of the merging of unequal progenitors.

We report on resent N-body simulations of galaxy formation performed on the GRAPE-4 (GRAvity PipE) system, a special-purpose computer for astrophysical N-body simulations. We review the astrophysical motivation, the algorithm, the actual performance, and the price per performance. The performance obtained is 332 Gflops averaged over 185 hours for a simulation of a galaxy formation with 786,400 particles. The price per performance obtained is 4,600 dollars per Gflops. The configuration used for the simulation consists of 1,269 pipeline processors and has a peak speed of 663 Gflops.

The homogeneity of the cosmic microwave background radiation (CBR) is one of the most severe constraint for theories of the structure formation in the universe. We investigated the effect of the gravitational scattering (lensing) of galaxies, clusters of galaxies, and superclusters on the anisotropy of the CBR by numerical simulations. Although this effect was thought to be unimportant, we found that the gravitational scatterings by superclusters can significantly reduce the anisotropy of the CBR. We took into account the exponential growth of the distance between two rays due to multiple scatterings. The bending angle of each ray grows through the random walk process. On the other hand, difference between two rays grows exponentially while it is small. This exponential growth is caused by coherent scatterings that two rays suffer, and was neglected in the previous studies. The gravitational scattering by superclusters reduces the observed temperature anisotropy of the CBR at present time approximately by 40--60 $\%$ from that at the recombination time for angular scale up to a few degrees if the supercluster were formed at $z=$2--4.

We give an estimate of the relative error in the angular measurement of observations for high redshift objects induced by gravitational scatterings (lensing). Gunn (1967) concluded that the gravitational scatterings by galaxies induce the relative error of a few percent in the observations for objects at $z=1$. This estimate has been considered as a fundamental limitation of accuracy of the angular measurements in the observational cosmology. In multiple graviational scatterings, bending angle of single ray grows through the random work process. Gunn (1967) assumed that the difference of nearby rays also grows through the random walk process. However, distance between nearby photons grows exponentially because the two rays suffer coherent scatterings by the same scattering object. This exponential growth continues as long as the scattering is coherent. In the case of scattering by individual galaxies, the exponential growth continues until the angular distance reaches an arcminute or so. The relative error of the angular measurements under an arcminute due to the exponential growth is $\sim 30\%$ at $z=1$ and exceeds $100\%$ at $z=3$, in the case that the density parameter of galaxies is 0.2. The effects of clusters of galaxies or superclusters are more difficult to estimate accurately, but might be significant. In the case of supercluster the angular measurements up to a few degrees could be affected.

This paper is concerned with collisionless aspects of the early evolution of model star clusters. The effects of mass loss through stellar evolution and of a steady tidal field are modelled using $N$-body simulations. Our results (which depend on the assumed initial structure and the mass spectrum) agree qualitatively with those of Chernoff \& Weinberg (1990), who used a Fokker-Planck model with a spherically symmetric tidal cutoff. For those systems which are disrupted, the lifetime to disruption generally exceeds that found by Chernoff \& Weinberg, sometimes by as much as an order of magnitude. Because we do not model collisional effects correctly we cannot establish the fate of the survivors. In terms of theoretical interpretation, we find that tidal disruption must be understood as a loss of {\sl equilibrium}, and not a loss of {\sl stability}, as is sometimes stated.