Several binary black holes (BBHs) have been observed using gravitational wave detectors. For the formation mechanism of BBHs, two main mechanisms, isolated binary evolution and dynamical formation in dense star clusters, have been suggested. Future observations are expected to provide more information about BBH distributions, and it will help us to distinguish the two formation mechanisms. For the star cluster channel, globular clusters have mainly been investigated. However, recent simulations have suggested that BBH formation in open clusters is not negligible. We estimate a local merger rate density of BBHs originated from open clusters using the results of our $N$-body simulations of open clusters with four different metallicities. We find that the merger rate per cluster is the highest for our 0.1 solar metallicity model. Assuming a cosmic star formation history and a metallicity evolution with dispersion, we estimate the local merger rate density of BBHs originated from open clusters to be $\sim 70~{\rm yr}^{-1} {\rm Gpc}^{-3}$. This value is comparable to the merger rate density expected from the first and second observation runs of LIGO and Virgo. In addition, we find that BBH mergers obtained from our simulations can reproduce the distribution of primary mass and mass ratio of merging BBHs estimated from the LIGO and Virgo observations.

We have devised fitting formulae for evolution tracks of massive stars with $8 \lesssim M/M_\odot \lesssim 160$ under extreme metal poor (EMP) environments for $\log (Z/Z_\odot) = -2, -4, -5, -6$, and $-8$, where $M_\odot$ and $Z_\odot$ are the solar mass and metallicity, respectively. Our fitting formulae are based on reference stellar models which we have newly obtained by simulating the time evolutions of EMP stars. Our fitting formulae take into account stars ending with blue supergiant (BSG) stars, and stars skipping Hertzsprung gap (HG) phases and blue loops, which are characteristics of massive EMP stars. In our fitting formulae, stars may remain BSG stars when they finish their core Helium burning (CHeB) phases. Our fitting formulae are in good agreement with our stellar evolution models. We can use these fitting formulae on the SSE, BSE, NBODY4, and NBODY6 codes, which are widely used for population synthesis calculations and star cluster simulations. These fitting formulae should be useful to make theoretical templates of binary black holes formed under EMP environments.

Black hole-main sequence star (BH-MS) binaries are one of the targets of the future data releases of the astrometric satellite {\it Gaia}. They are supposed to be formed in two main sites: a galactic field and star clusters. However, previous work has never predicted the number of BH-MS binaries originating in the latter site. In this paper, we estimate the number of BH-MS binaries formed in open clusters and detectable with {\it Gaia} based on the results of {\it N}-body simulations. By considering interstellar extinction in the Milky Way (MW) and observational constraints, we predict $\sim 10$ BH-MS binaries are observable. We also find that chemical abundance patterns of companion MSs will help us to identify the origin of the binaries as star clusters. Such MSs are not polluted by outflows of the BH progenitors, such as stellar winds and supernova ejecta. Chemical anomalies might be a good test to confirm the origin of binaries with relatively less massive MSs ($\lesssim 5M_{\odot}$), orbital periods ($\sim 1.5\;$year) and higher eccentricities ($e \gtrsim 0.1$).

We estimate formation rates of LB-1-like systems through dynamical interactions in the framework of the theory of stellar evolution before the discovery of the LB-1 system. The LB-1 system contains $\sim 70M_\odot$ black hole (BH), so-called pair instability (PI)-gap BH, and B-type star with solar metallicity, and has nearly zero eccentricity. The most efficient formation mechanism is as follows. In an open cluster, a naked helium (He) star (with $\sim 20M_\odot$) collides with a heavy main-sequence (MS) star (with $\sim 50M_\odot$) which has a B-type companion. The collision results in a binary consisting of the collision product and B-type star with a high eccentricity. The binary can be circularized through the dynamical tide with radiative damping of the collision-product envelope. Finally, the collision product collapses to a PI-gap BH, avoiding pulsational pair instability and pair instability supernovae because its He core is as massive as the pre-colliding naked He star. We find that the number of LB-1-like systems in the Milky Way galaxy is $\sim 0.01 (\rho_{\rm oc} / 10^4 M_\odot \mbox{pc}^{-3})$, where $\rho_{\rm oc}$ is the initial mass densities of open clusters. If we take into account LB-1-like systems with O-type companion stars, the number increases to $\sim 0.03 (\rho_{\rm oc} / 10^4 M_\odot \mbox{pc}^{-3})$. This mechanism can form LB-1-like systems at least 10 times more efficiently than the other mechanisms: captures of B-type stars by PI-gap BHs, stellar collisions between other type stars, and stellar mergers in hierarchical triple systems. We conclude that no dynamical mechanism can explain the presence of the LB-1 system.

We study the emission properties of thermonuclear explosions in a helium white dwarf (WD) tidal disruption event (TDE). We consider a TDE where a 0.2 $M_{\odot}$ helium WD is disrupted by a $10^{2.5}\,M_{\odot}$ intermediate-mass black hole (IMBH). The helium WD is not only tidally disrupted but is also detonated by the tidal compression and by succeeding shocks. We focus on the emission powered by radioactive nuclei in the unbound TDE ejecta. We perform hydrodynamic simulations coupled with nuclear reactions, post-process detailed nucleosynthesis calculations, and then radiative transfer simulations. We thus derive multi-band light curves and spectra. The helium WD TDE shows rapid ($\Delta t_{1\mathrm{mag}}\simeq5\text{--}10$ days) and relatively faint ($L_{\mathrm{peak}}\simeq10^{42}\,\mathrm{erg}\,\mathrm{s}^{-1}$) light curves, because the ejecta mass and $^{56}$Ni mass are low ($0.12\,M_{\odot}$ and $0.03\,M_{\odot}$, respectively). The spectra show strong calcium and Fe-peak features and very weak silicon features, reflecting the peculiar elemental abundance. The key feature is the Doppler shift of the spectral lines up to $\simeq\pm12,000\,\mathrm{km}\,\mathrm{s}^{-1}$, depending on the viewing angle, due to the bulk motion of the ejecta. Our model matches well with two rapid and faint transients reported in Pursiainen et al. (2018). The particular model presented here does not match with observed SNe Iax, calcium-rich transients, or .Ia explosion candidates, either in the spectra or light curves. However, we expect a large variety of the observational signatures once a wide range of the WD/black hole masses and orbital parameters are considered. This study contributes to the search for WD TDEs with current and upcoming surveys, and to the identification of IMBHs as disrupters in the TDEs.

Using a prototype of the Tomo-e Gozen wide-field CMOS mosaic camera, we acquire wide-field optical images at a cadence of 2 Hz and search them for transient sources of duration 1.5 to 11.5 seconds. Over the course of eight nights, our survey encompasses the equivalent of roughly two days on one square degree, to a fluence equivalent to a limiting magnitude about $V = 15.6$ in a 1-second exposure. After examining by eye the candidates identified by a software pipeline, we find no sources which meet all our criteria. We compute upper limits to the rate of optical transients consistent with our survey, and compare those to the rates expected and observed for representative sources of ephemeral optical light.

We have studied double-detonation explosions in double-degenerate (DD) systems with different companion white dwarfs (WD) for modeling type Ia supernovae (SNe Ia) by means of high-resolution smoothed particle hydrodynamics (SPH) simulations. We have found that only the primary WDs explode in some of the DD systems, while the explosions of the primary WDs induce the explosions of the companion WDs in the other DD systems. The former case is so-called Dynamically-Driven Double-Degenerate Double-Detonation (D$^6$) explosion, or helium-ignited violent merger explosion. The supernova ejecta of the primary WDs strip materials from the companion WDs, whose mass is $\sim 10^{-3}M_\odot$. The stripped materials contain carbon and oxygen when the companion WDs are carbon-oxygen (CO) WDs with He shells $\lesssim 0.04M_\odot$. Since they contribute to low-velocity ejecta components as observationally interfered for iPTF14atg, D$^6$ explosions can be counterparts of sub-luminous SNe Ia. The stripped materials may contribute to low-velocity C seen in several SNe Ia. In the latter case, the companion WDs explode through He detonation if they are He WDs, and through double-detonation mechanism if they are CO WDs with He shells. We name these explosions "triple" and "quadruple" detonation (TD/QD) explosions after the number of detonations. The QD explosion may be counterparts of luminous SNe Ia, such as SN 1991T and SN 1999aa, since they yield a large amount of $^{56}$Ni, and their He-detonation products contribute to the early emissions accompanying such luminous SNe Ia. On the other hand, the TD explosion may not yield a sufficient amount of $^{56}$Ni to explain luminous SNe Ia.

We study tidal disruption events of rotating stars by a supermassive black hole in a galactic nucleus by using a smoothed-particle hydrodynamics (SPH) code. We compare mass infall rates of tidal-disruption debris of a non-rotating and of a rotating star when they come close to the supermassive black hole. Remarkably the mass distribution of debris bound to the black hole as a function of specific energy shows clear difference between rotating and non-rotating stars, even if the stellar rotation is far from the break-up limit. The debris of a star whose initial spin is parallel to the orbital angular momentum has a mass distribution which extends to lower energy than that of non-rotating star. The debris of a star with anti-parallel spin has a larger energy compared with a non-rotating counterpart. As a result, debris from a star with anti-parallel spin is bound more loosely to the black hole and the mass-infall rate rises later in time, while that of a star with a parallel spin is tightly bound and falls back to the black hole earlier. The different rising timescales of mass-infall rate may affect the early phase of flares due to the tidal disruptions. In the Appendix we study the disruptions by using a uniform-density ellipsoid model which approximately takes into account the effect of strong gravity of the black hole. We find the mass infall rate reaches its maximum earlier for strong gravity cases because the debris is trapped in a deeper potential well of the black hole.

We investigate metal pollution onto the surface of low-mass population III stars (Pop. III survivors) via interstellar objects floating in the Galactic interstellar medium. Only recently, Tanikawa et al. analytically estimated how much metal should collide to an orbiting Pop. III survivor encouraged by the recent discovery of 'Oumuamua and suggested that ISOs are the most dominant contributor of metal enrichment of Pop. III survivors. When we consider a distribution of interstellar objects in the Galactic disc, Pop. III survivors' orbits are significant properties to estimate the accretion rate of them though Tanikawa et al. assumed one modelled orbit. To take more realistic orbits into calculating the accretion rate, we use a high-resolution cosmological $N$-body simulation that resolves dark matter minihaloes. Pop. III survivors located at solar neighbourhood have a number of chances of ISO($> 100$ m) collisions, typically $5\times10^6$ times in the last $5$ Gyr, which is one order of magnitude greater than estimated in the previous study. When we assume a power-law parameter $\alpha$ of the ISO cumulative number density with size greater than $D$ as $n \propto D^{-\alpha}$, $0.80 \, M_{\odot}$ stars should be typically polluted [Fe/H]$\sim -2$ for the case of $\alpha=2.0$. Even in the cases of $0.70$ and $0.75 \, M_{\odot}$ stars, the typical surface metallicity are around [Fe/H]$=-6 \sim -5$. From the presence of stars with their [Fe/H], we can constrain on the lower limit of the power $\alpha$, as $\alpha \gtrsim 2.0$, which is consistent with $\alpha$ of km-size asteroids and comets in the solar system. Furthermore, we provide six candidates as the ISO-polluted Pop. III stars in the case of $\alpha \sim 2.5$. Metal-poor stars so far discovered are possible to be metal-free Pop. III stars on birth.

In order to investigate the formation rate of binary black holes (BBHs) in stellar clusters with a mass comparable to open clusters, we performed a series of direct $N$-body simulations of open clusters with a mass of $2.5\times10^3$ (Model A) and $10^4 M_{\odot}$ (Model B). Since such low-mass clusters would have been more populous than globular clusters when they were born, low-mass clusters are also candidates as the origin of BBHs which are the source of the gravitational waves. In model A, most of BBHs merged within 10 Gyr formed via dynamically formed main-sequence binary stars and stable and unstable mass transfer between them since open clusters collapse within the main-sequence life-time of massive stars. These binaries, therefore, have little eccentricities. The fraction of such binaries among all merging BBHs increases as the cluster mass decrease due to the shorter relaxation time. In our simulations, $4.0\times10^{-5}$ and $1.7\times10^{-5}$ BBHs per solar mass merged within 10 Gyr for models A and B, respectively. These values correspond to $\sim$ 20-50% of the number of mergers per solar mass originated from globular clusters with a mass of $10^5$-$10^6M_{\odot}$. Thus, the contribution of BBHs originated from open clusters is not negligible. The estimated mergers rate density in the local universe is about 0.3 yr$^{-1}$ Gpc$^{-3}$ assuming a cluster mass function with a power of $-$2.

We investigate the properties of young massive clusters (YMCs) in a galaxy-galaxy merger remnant by analyzing the data obtained by a gas rich major merger simulations in Matsui et al. 2012. We found that the YMCs are distributed at a few $\rm kpc$ and at $\sim 10~{\rm kpc}$ from the galactic center, in other words, there are two components of their distribution. The former are formed in filamentary and turbulent gas generated at a few $\rm kpc$ from the center because of galaxy encounters, and the latter are formed in tidal tails which are far from the center. The YMCs are much less concentrated than galaxy stars. The mass function of the YMCs is $dN/dM \propto M^{-2}$. Most of YMCs are formed from the second encounter to the final coalescence phase of the galactic cores, and their formation rate is especially high at final coalescence phase. Most of them consists of single stellar population in age, but YMCs with multi stellar populations in age are also formed. The multiple populations are produced by the following process: a YMC captures dense gas, and another generation stars form within the cluster. There are several YMCs formed in an isolated disk before the encounter of galaxies. These candidates contain stars with various age by capturing dense gas and forming stars. YMCs in a merger remnant, have various orbits, but large fraction of candidates have circular orbits.

We study the hydrodynamics and nucleosynthesis in the double-detonation model of Type Ia supernovae (SNe~Ia) and the interaction between the ejecta and a surviving white dwarf (WD) companion in the double-degenerate scenario. We set up a binary star system with $1.0M_\odot$ and $0.6M_\odot$ carbon-oxygen (CO) WDs, where the primary WD consists of a CO core and helium (He) shell with $0.95$ and $0.05M_\odot$, respectively. We follow the evolution of the binary star system from the initiation of a He detonation, ignition and propagation of a CO detonation, and the interaction of SN ejecta with the companion WD. The companion (or surviving) WD gets a flung-away velocity of $\sim 1700$~km~s$^{-1}$, and captures $^{56}$Ni of $\sim 0.03M_\odot$, and He of $3 \times 10^{-4}M_\odot$. Such He can be detected on the surface of surviving WDs. The SN ejecta contains a `companion-origin stream', and unburned materials stripped from the companion WD ($\sim 3 \cdot 10^{-3}M_\odot$), although the stream compositions would depend on the He shell mass of the companion WD. The ejecta has also a velocity shift of $\sim 1000$~km~s$^{-1}$ due to the binary motion of the exploding primary WD. These features would be prominent in nebular-phase spectra of oxygen emission lines from the unburned materials like SN~2010lp and iPTF14atg, and of blue- or red-shifted Fe-group emission lines from the velocity shift like a part of sub-luminous SNe~Ia. We expect SN~Ia counterparts to the D$^6$ model would leave these fingerprints for SN~Ia observations.

We calculate accretion mass of interstellar objects (ISOs) like `Oumuamua onto low-mass population III stars (Pop.~III survivors), and estimate surface pollution of Pop.~III survivors. An ISO number density estimated from the discovery of `Oumuamua is so high ($\sim 0.2$~au$^{-3}$) that Pop.~III survivors have chances at colliding with ISOs $\gtrsim 10^5$ times per $1$~Gyr. `Oumuamua itself would be sublimated near Pop.~III survivors, since it has small size, $\sim 100$~m. However, ISOs with size $\gtrsim 3$~km would reach the Pop.~III survivor surfaces. Supposing an ISO cumulative number density with size larger than $D$ is $n \propto D^{-\alpha}$, Pop.~III survivors can accrete ISO mass $\gtrsim 10^{-16}M_\odot$, or ISO iron mass $\gtrsim 10^{-17}M_\odot$, if $\alpha < 4$. This iron mass is larger than the accretion mass of interstellar medium (ISM) by several orders of magnitude. Taking into account material mixing in a convection zone of Pop.~III survivors, we obtain their surface pollution is typically [Fe/H] $\lesssim -8$ in most cases, however the surface pollution of Pop.~III survivors with $0.8M_\odot$ can be [Fe/H] $\gtrsim -6$ because of the very shallow convective layer. If we apply to Pop.III survivors located at the Galactocentric distance of 8 kpc, the dependence of the metal pollustion is as follows. If $\alpha > 4$, Pop.~III survivors have no chance at colliding with ISOs with $D \gtrsim 3$~km, and keep metal-free. If $3 < \alpha < 4$, Pop.~III survivors would be most polluted by ISOs up to [Fe/H] $\sim -7$. If $\alpha < 3$ up to $D \sim 10$~km, Pop.~III survivors could hide in metal-poor stars so far discovered. Pop.~III survivors would be more polluted with decreasing the Galactocentric distance. Although the metal pollution depends on $\alpha$ and the Galactocentric distance, we first show the importance of ISOs for the metal pollution of Pop.~III survivors.

Numerical simulations based on particle methods have been widely used in various fields including astrophysics. To date, simulation softwares have been developed by individual researchers or research groups in each field, with a huge amount of time and effort, even though numerical algorithms used are very similar. To improve the situation, we have developed a framework, called FDPS, which enables researchers to easily develop massively parallel particle simulation codes for arbitrary particle methods. Until version 3.0, FDPS have provided API only for C++ programing language. This limitation comes from the fact that FDPS is developed using the template feature in C++, which is essential to support arbitrary data types of particle. However, there are many researchers who use Fortran to develop their codes. Thus, the previous versions of FDPS require such people to invest much time to learn C++. This is inefficient. To cope with this problem, we newly developed a Fortran interface layer in FDPS, which provides API for Fortran. In order to support arbitrary data types of particle in Fortran, we design the Fortran interface layer as follows. Based on a given derived data type in Fortran representing particle, a Python script provided by us automatically generates a library that manipulates the C++ core part of FDPS. This library is seen as a Fortran module providing API of FDPS from the Fortran side and uses C programs internally to interoperate Fortran with C++. In this way, we have overcome several technical issues when emulating `template' in Fortran. By using the Fortran interface, users can develop all parts of their codes in Fortran. We show that the overhead of the Fortran interface part is sufficiently small and a code written in Fortran shows a performance practically identical to the one written in C++.

We run a suite of hydrodynamics simulations of tidal disruption events (TDEs) of a white dwarf (WD) by a black hole (BH) with a wide range of WD/BH masses and orbital parameters. We implement nuclear reactions to study nucleosynthesis and its dynamical effect through release of nuclear energy. The released nuclear energy effectively increases the fraction of unbound ejecta. This effect is weaker for a heavy WD with 1.2 $\mathrm{M}_{\odot}$, because the specific orbital energy distribution of the debris is predominantly determined by the tidal force, rather than by the explosive reactions. The elemental yield of a TDE depends critically on the initial composition of a WD, while the BH mass and the orbital parameters also affect the total amount of synthesized elements. Tanikawa et al. (2017) find that simulations of WD-BH TDEs with low resolution suffer from spurious heating and inaccurate nuclear reaction results. In order to examine the validity of our calculations, we compare the amounts of the synthesized elements with the upper limits of them derived in a way where we can avoid uncertainties due to low resolution. The results are largely consistent, and thus support our findings. We find particular TDEs where early self-intersection of a WD occurs during the first pericentre passage, promoting formation of an accretion disc. We expect that relativistic jets and/or winds would form in these cases because accretion rates would be super-Eddington. The WD-BH TDEs result in a variety of events depending on the WD/BH mass and pericentre radius of the orbit.

We demonstrate tidal detonation during a tidal disruption event (TDE) of a helium (He) white dwarf (WD) with $0.45M_\odot$ by an intermediate mass black hole (IMBH) by extremely high-resolution simulations. Tanikawa et al. (2017) have showed tidal detonation in previous studies results from unphysical heating due to low-resolution simulations, and such unphysical heating occurs in 3-dimensional (3D) smoothed particle hydrodynamics (SPH) simulations even with $10$ million SPH particles. In order to avoid such unphysical heating, we perform 3D SPH simulations up to $300$ million SPH particles, and 1D mesh simulations using flow structure in the 3D SPH simulations for 1D initial conditions. The 1D mesh simulations have higher resolution than the 3D SPH simulations. We show tidal detonation occurs, and confirm this result is perfectly converged with different space resolution in both 3D SPH and 1D mesh simulations. We find detonation waves independently arises in leading parts of the WD, and yield large amounts of $^{56}$Ni. Although detonation waves are not generated in trailing parts of the WD, the trailing parts receive detonation waves generated in the leading parts, and leave large amounts of Si group elements. Eventually, this He WD TDE would synthesize $^{56}$Ni of $0.30M_\odot$ and Si group elements of $0.08M_\odot$, and could be observed as a luminous thermonuclear transient comparable to type Ia supernovae.

We suggest "tidal double detonation": a new mechanism for a thermonuclear explosion of a white dwarf (WD) induced by a tidal disruption event (TDE). Tidal detonation is also a WD explosion induced by a TDE. In this case, helium (He) and carbon-oxygen (CO) detonation waves incinerate He~WD and CO~WD, respectively. On the other hand, for tidal double detonation, He detonation is first excited in the He shell of a CO~WD, and drives CO detonation in the CO core. We name this mechanism after the double detonation scenario in the context of type Ia supernovae. In this paper, we show tidal double detonation occurs in shallower encounter of a CO~WD with an intermediate mass black hole (IMBH) than simple tidal detonation, performing numerical simulations for CO~WDs with $0.60M_\odot$ with and without a He shell. We expect tidal double detonation spreads opportunity to WD~TDEs illuminating IMBHs.

Advanced LIGO achieved the first detection of the gravitational wave, which was from a merging binary black hole (BBH). In the near future, more merger events will be observed, and the mass distribution of them will become available. The mass distribution of merger events reflects the evolutionary path of BBHs: dynamical formation in dense star clusters or common envelope evolution from primordial binaries. In this paper, we estimate the detection rate of merging BBHs which dynamically formed in dense star clusters by combining the results of $N$-body simulations, modeling of globular clusters, and cosmic star-cluster formation history. We estimate that the merger rate density in the local universe within the redshift of 0.1 is 13--57 Gpc$^{-3}$ yr$^{-1}$. We find that the detection rate is 0.23--4.6 per year for the current sensitivity limit and that it would increase to 5.1--99 per year for the designed sensitivity which will be achieved in 2019. The distribution of merger rate density in the local universe as a function of redshifted chirp mass has a peak close to the low-mass end.The chirp mass function of the detected mergers, on the other hand, has a peak at the high-mass end, but is almost flat. This difference is simply because the detection range is larger for more massive BBHs.

We investigate nucleosynthesis in tidal disruption events (TDEs) of white dwarfs (WDs) by intermediate mass black holes (IMBHs). We consider various types of WDs with different masses and compositions by means of 3 dimensional (3D) smoothed particle hydrodynamics (SPH) simulations. We model these WDs with different numbers of SPH particles, $N$, from a few $10^4$ to a few $10^7$, in order to check mass resolution convergence, where SPH simulations with $N>10^7$ (or a space resolution of several $10^6$ cm) have unprecedentedly high resolution in this kind of simulations. We find that nuclear reactions become less active with increasing $N$, and that these nuclear reactions are excited by spurious heating due to low resolution. Moreover, we find no shock wave generation. In order to investigate the reason for the absence of a shock wave, we additionally perform 1 dimensional (1D) SPH and mesh-based simulations with a space resolution ranging from $10^4$ to $10^7$ cm, using characteristic flow structure extracted from the 3D SPH simulations. We find shock waves in these 1D high-resolution simulations. One of these shock waves triggers a detonation wave. However, we have to be careful of the fact that, if the shock wave emerged at a bit outer region, it could not trigger the detonation wave due to low density. Note that the 1D initial conditions lack accuracy to precisely determine where a shock wave emerges. We need to perform 3D simulations with $\lesssim 10^6$ cm space resolution in order to conclude that WD TDEs become optical transients powered by radioactive nuclei.

The giant impact (GI) is one of the most important hypotheses both in planetary science and geoscience, since it is related to the origin of the Moon and also the initial condition of the Earth. A number of numerical simulations have been done using the smoothed particle hydrodynamics (SPH) method. However, GI hypothesis is currently in a crisis. The "canonical" GI scenario failed to explain the identical isotope ratio between the Earth and the Moon. On the other hand, little has been known about the reliability of the result of GI simulations. In this paper, we discuss the effect of the resolution on the results of the GI simulations by varying the number of particles from $3 \times10^3$ to $10^8$. We found that the results does not converge, but shows oscillatory behaviour. We discuss the origin of this oscillatory behaviour.

We present the basic idea, implementation, measured performance and performance model of FDPS (Framework for developing particle simulators). FDPS is an application-development framework which helps the researchers to develop particle-based simulation programs for large-scale distributed-memory parallel supercomputers. A particle-based simulation program for distributed-memory parallel computers needs to perform domain decomposition, redistribution of particles, and gathering of particle information for interaction calculation. Also, even if distributed-memory parallel computers are not used, in order to reduce the amount of computation, algorithms such as Barnes-Hut tree method should be used for long-range interactions. For short-range interactions, some methods to limit the calculation to neighbor particles are necessary. FDPS provides all of these necessary functions for efficient parallel execution of particle-based simulations as "templates", which are independent of the actual data structure of particles and the functional form of the interaction. By using FDPS, researchers can write their programs with the amount of work necessary to write a simple, sequential and unoptimized program of O(N^2) calculation cost, and yet the program, once compiled with FDPS, will run efficiently on large-scale parallel supercomputers. A simple gravitational N-body program can be written in around 120 lines. We report the actual performance of these programs and the performance model. The weak scaling performance is very good, and almost linear speedup was obtained for up to the full system of K computer. The minimum calculation time per timestep is in the range of 30 ms (N=10^7) to 300 ms (N=10^9). These are currently limited by the time for the calculation of the domain decomposition and communication necessary for the interaction calculation. We discuss how we can overcome these bottlenecks.

Mergers of carbon-oxygen (CO) white dwarfs (WDs) are considered as one of the potential progenitors of type Ia supernovae (SNe Ia). Recent hydrodynamical simulations showed that the less massive (secondary) WD violently accretes onto the more massive (primary) one, carbon detonation occurs, the detonation wave propagates through the primary, and the primary finally explodes as a sub-Chandrasekhar mass SN Ia. Such an explosion mechanism is called the violent merger scenario. Based on the smoothed particle hydrodynamics (SPH) simulations of merging CO WDs, we derived more stringent critical mass ratio (qcr) leading to the violent merger scenario than the previous results. We conclude that this difference mainly comes from the differences in the initial condition, synchronously spinning of WDs or not. Using our new results, we estimated the brightness distribution of SNe Ia in the violent merger scenario and compared it with previous studies. We found that our new qcr does not significantly affect the brightness distribution. We present the direct outcome immediately following CO WD mergers for various primary masses and mass ratios. We also discussed the final fate of the central system of the bipolar planetary nebula Henize 2-428, which was recently suggested to be a double CO WD system whose total mass exceeds the Chandrasekhar-limiting mass, merging within the Hubble time. Even considering the uncertainties in the proposed binary parameters, we concluded that the final fate of this system is almost certainly a sub-Chandrasekhar mass SN Ia in the violent merger scenario.

Mergers of two carbon-oxygen (CO) white dwarfs (WDs) have been considered as progenitors of Type Ia supernovae (SNe Ia). Based on smoothed particle hydrodynamics (SPH) simulations, previous studies claimed that mergers of CO WDs lead to an SN Ia explosion either in the dynamical merger phase or stationary rotating merger remnant phase. However, the mass range of CO WDs that lead to an SN Ia has not been clearly identified yet. In the present work, we perform systematic SPH merger simulations for the WD masses ranging from $0.5~M_{\odot}$ to $1.1~M_{\odot}$ with higher resolutions than the previous systematic surveys and examine whether or not carbon burning occurs dynamically or quiescently in each phase. We further study the possibility of SN Ia explosion and estimate the mass range of CO WDs that lead to an SN Ia. We found that when the both WDs are massive, i.e., in the mass range of $0.9~M_{\odot} {\le} M_{1,2} {\le} 1.1~M_{\odot}$, they can explode as an SN Ia in the merger phase. On the other hand, when the more massive WD is in the range of $0.7~M_{\odot} {\le} M_{1} {\le} 0.9~M_{\odot}$ and the total mass exceeds $1.38~M_{\odot}$, they can finally explode in the stationary rotating merger remnant phase. We estimate the contribution of CO WD mergers to the entire SN Ia rate in our galaxy to be of ${\lt} 9\%$. So, it might be difficult to explain all galactic SNe Ia by CO WD mergers.

We perform smoothed particle hydrodynamics (SPH) simulations for merging binary carbon-oxygen (CO) white dwarfs (WDs) with masses of $1.1$ and $1.0$ $M_\odot$, until the merger remnant reaches a dynamically steady state. Using these results, we assess whether the binary could induce a thermonuclear explosion, and whether the explosion could be observed as a type Ia supernova (SN Ia). We investigate three explosion mechanisms: a helium-ignition following the dynamical merger (`helium-ignited violent merger model'), a carbon-ignition (`carbon-ignited violent merger model'), and an explosion following the formation of the Chandrasekhar mass WD (`Chandrasekhar mass model'). An explosion of the helium-ignited violent merger model is possible, while we predict that the resulting SN ejecta are highly asymmetric since its companion star is fully intact at the time of the explosion. The carbon-ignited violent merger model can also lead to an explosion. However, the envelope of the exploding WD spreads out to $\sim 0.1R_\odot$; it is much larger than that inferred for SN 2011fe ($< 0.1R_\odot $) while much smaller than that for SN 2014J ($\sim 1R_\odot$). For the particular combination of the WD masses studied in this work, the Chandrasekhar mass model is not successful to lead to an SN Ia explosion. Besides these assessments, we investigate the evolution of unbound materials ejected through the merging process (`merger ejecta'), assuming a case where the SN Ia explosion is not triggered by the helium- or carbon-ignition during the merger. The merger ejecta interact with the surrounding interstellar medium, and form a shell. The shell has a bolometric luminosity of more than $2 \times 10^{35}$ ergs$^{-1}$ lasting for $\sim 2 \times 10^4$ yr. If this is the case, Milky Way should harbor about $10$ such shells at any given time.

We perform N-body simulations on a multiple massive black hole (MBH) system in a host galaxy to derive the criteria for successive MBH merger. The calculations incorporate the dynamical friction by stars and general relativistic effects as pericentre shift and gravitational wave recoil. The orbits of MBHs are pursed down to ten Schwarzschild radii (~ 1 AU). As a result, it is shown that about a half of MBHs merge during 1 Gyr in a galaxy with mass $10^{11}M_{\odot}$ and stellar velocity dispersion 240 km/s, even if the recoil velocity is two times as high as the stellar velocity dispersion. The dynamical friction allows a binary MBH to interact frequently with other MBHs, and then the decay of the binary orbits leads to the merger through gravitational wave radiation, as shown by Tanikawa & Umemura (2011). We derive the MBH merger criteria for the masses, sizes, and luminosities of host galaxies. It is found that the successive MBH mergers are expected in bright galaxies, depending on redshifts. Furthermore, we find that the central stellar density is reduced by the sling-shot mechanism and that high-velocity stars with ~ 1000 km/s are generated intermittently in extremely radial orbits.

At the moment of deepest core collapse, a star cluster core contains less than ten stars. This small number makes the traditional treatment of hard binary formation, assuming a homogeneous background density, suspect. In a previous paper, we have found that indeed the conventional wisdom of binary formation, based on three-body encounters, is incorrect. Here we refine that insight, by further dissecting the subsequent steps leading to hard binary formation. For this purpose, we add some analysis tools in order to make the study less subjective. We find that the conventional treatment does remain valid for direct three-body scattering, but fails for resonant three-body scattering. Especially democratic resonance scattering, which forms an important part of the analytical theory of three-body binary formation, takes too much space and time to be approximated as being isolated, in the context of a cluster core around core collapse. We conclude that, while three-body encounters can be analytically approximated as isolated, subsequent strong perturbations typically occur whenever those encounters give rise to democratic resonances. We present analytical estimates postdicting our numerical results. If we only had been a bit more clever, we could have predicted this qualitative behaviour.

We have performed N-body simulations of globular clusters (GCs) in order to estimate a detection rate of mergers of Binary stellar-mass Black Holes (BBHs) by means of gravitational wave (GW) observatories. For our estimate, we have only considered mergers of BBHs which escape from GCs (BBH escapers). BBH escapers merge more quickly than BBHs inside GCs because of their small semi-major axes. N-body simulation can not deal with a GC with the number of stars N ~ 10^6 due to its high calculation cost. We have simulated dynamical evolution of small-N clusters (10^4 <~ N <~ 10^5), and have extrapolated our simulation results to large-N clusters. From our simulation results, we have found the following dependence of BBH properties on N. BBHs escape from a cluster at each two-body relaxation time at a rate proportional to N. Semi-major axes of BBH escapers are inversely proportional to N, if initial mass densities of clusters are fixed. Eccentricities, primary masses, and mass ratios of BBH escapers are independent of N. Using this dependence of BBH properties, we have artificially generated a population of BBH escapers from a GC with N ~ 10^6, and have estimated a detection rate of mergers of BBH escapers by next-generation GW observatories. We have assumed that all the GCs are formed 10 or 12Gyrs ago with their initial numbers of stars N_i=5 x 10^5 -- 2 x 10^6 and their initial stellar mass densities inside their half-mass radii \rho_h,i=6 x 10^3 -- 10^6M_sun pc^-3. Then, the detection rate of BBH escapers is 0.5 -- 20 yr^-1 for a BH retention fraction R_BH=0.5. A few BBH escapers are components of hierarchical triple systems, although we do not consider secular perturbation on such BBH escapers for our estimate. Our simulations have shown that BHs are still inside some of GCs at the present day. These BHs may marginally contribute to BBH detection.

(Abridged) We have developed a numerical software library for collisionless N-body simulations named "Phantom-GRAPE" which highly accelerates force calculations among particles by use of a new SIMD instruction set extension to the x86 architecture, AVX, an enhanced version of SSE. In our library, not only the Newton's forces, but also central forces with an arbitrary shape f(r), which has a finite cutoff radius r_cut (i.e. f(r)=0 at r>r_cut), can be quickly computed. Using an Intel Core i7--2600 processor, we measure the performance of our library for both the forces. In the case of Newton's forces, we achieve 2 x 10^9 interactions per second with 1 processor core, which is 20 times higher than the performance of an implementation without any explicit use of SIMD instructions, and 2 times than that with the SSE instructions. With 4 processor cores, we obtain the performance of 8 x 10^9 interactions per second. In the case of the arbitrarily shaped forces, we can calculate 1 x 10^9 and 4 x 10^9 interactions per second with 1 and 4 processor cores, respectively. The performance with 1 processor core is 6 times and 2 times higher than those of the implementations without any use of SIMD instructions and with the SSE instructions. These performances depend weakly on the number of particles. It is good contrast with the fact that the performance of force calculations accelerated by GPUs depends strongly on the number of particles. Substantially weak dependence of the performance on the number of particles is suitable to collisionless N-body simulations, since these simulations are usually performed with sophisticated N-body solvers such as Tree- and TreePM-methods combined with an individual timestep scheme. Collisionless N-body simulations accelerated with our library have significant advantage over those accelerated by GPUs, especially on massively parallel environments.

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.

The conventional wisdom for the formation of the first hard binary in core collapse is that three-body interactions of single stars form many soft binaries, most of which are quickly destroyed, but eventually one of them survives. We report on direct N-body simulations to test these ideas, for the first time. We find that both assumptions are often incorrect: 1) quite a few three-body interactions produce a hard binary from scratch; 2) and in many cases there are more than three bodies directly and simultaneously involved in the production of the first binary. The main reason for the discrepancies is that the core of a star cluster, at the first deep collapse, contains typically only five or so stars. Therefore, the homogeneous background assumption, which still would be reasonable for, say, 25 stars, utterly breaks down. There have been some speculations in this direction, but we demonstrate this result here explicitly, for the first time.

Using highly-accurate $N$-body simulations, we explore the evolution of multiple massive black holes (hereafter, MBHs) in a primordial galaxy that is composed of stars and MBHs. The evolution is pursed with a fourth-order Hermite scheme, where not only three-body interaction of MBHs but also dynamical friction by stars are incorporated. Initially, ten MBHs with equal mass of $10^7M_\odot$ are set in a host galaxy with $10^{11}M_\odot$. It is found that 4 - 6 MBHs merge successively within 1 Gyr, emitting gravitational wave radiation. The key process for the successive merger of MBHs is the dynamical friction by field stars, which enhances three-body interactions of MBHs when they enter the central regions of the galaxy. The heaviest MBH always composes a close binary at the galactic center, which shrinks owing to the angular momentum transfer by the third MBH and eventually merges. The angular momentum transfer by the third MBH is due to the sling-shot mechanism. We find that the secular Kozai mechanism does not work for a binary to merge if we include the relativistic pericenter shift. The simulations show that a multiple MBH system can produce a heavier MBH at the galactic center purely through $N$-body process. This merger path can be of great significance for the growth of MBHs in a primordial galaxy. The merger of multiple MBHs may be a potential source of gravitational waves for the Laser Interferometer Space Antenna (LISA) and pulsar timing.

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.

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.