Research

Numerical Astrophysics / Galaxies / Galaxy Clusters / Large-Scale Structure

Technical Research

Phantom-GRAPE

The gravitational interaction is a one of the most important physical processes in astrophysics, especially in the studies of the formation of stars, galaxies and the large-scale structure in the universe. Numerical studies of such self-gravitating objects are manifestly carried out with N-body simulations, in which mass distribution of dark matter and baryons are represented with a number of super-particles and compute their orbits by computing gravitational forces among them. However, such calculations of gravitational forces are quite demanding and time consuming. We develop a numerical library "Phantom-GRAPE" to compute the gravitational forces among many particles quite efficiently. In the "Phantom-GRAPE" library, we utilize the Single-Instruction-Multi-Data (SIMD) extensions of computer processors, which perform the same operations on multiple data simultaneously. Thus, the performance of the "Phantom-GRAPE" is much higer than that obatined with the conventional implementation. Following two papers are helpful to understand the detailed implementation of Phantom-GRAPE library.

The package of "Phantom-GRAPE" is available at Bitbucket repository .

As of December 2018, Phantom-GRAPE supports the AVX-512 instruction set, the lastest SIMD instruction set available on x86 microprocessors, which is briefly introduced in a paper in Research Note of the AAS.

High Order Advection Solver for Vlasov Simulation

We developed a set of new numerical schemes to solve one-dimensional advection equation

$\displaystyle \frac{\partial f}{\partial t}+c\frac{\partial f}{\partial x}=0$,
where $f=f(x,t)$ is a function of a spatial coordinate $x$ and time $t$, and $c$ is the advection speed with very high-order accuracy. Such numerical schemes are critical to solve the Vlasov equation
$\displaystyle \frac{\partial f}{\partial t} + \boldsymbol{v}\cdot\frac{\partial f}{\partial \boldsymbol{x}}+\boldsymbol{F}\cdot\frac{\partial f}{\partial \boldsymbol{p}} = 0$,
because we can numerically solve it by splitting it into six one-dimensional advection equations along each axis of the six dimensional phase space. Our new advection schemes have spatially fifth- and seventh-order accuracy and preserve the monotonicity and the positivity of numerical solutions. For demonstrations, we perform 2-dimensional advection for rigid rotation described by the equation
$\displaystyle \frac{\partial f}{\partial t} + (\boldsymbol{x}\times\boldsymbol{\omega})\cdot\frac{\partial f}{\partial \boldsymbol{x}}=0$,
where $\boldsymbol\omega$ is a constant rotation vector and parallel to the $z$-axis. We show numerical solutions obtained with our new schemes SL-MPP7 and SL-MPP5 schemes which are of spatially seventh- and fifth-order accuracy as well as the PFC schemes adopted in our previous study, where the initial distribution function $f(\boldsymbol{x},t=0)$ is taken from the logo of Joint Institute for Computational Fundamental Science (JICFuS). It can be clearly seen that the distribution function solved with the SL-MPP7 scheme preseves its initial shape almost unchanged even after eight rotations, while the numerical solution obtained with the PFC scheme is significantly smeared by strong nuerical diffusion. The difference between numerical solutions obtained with SL-MPP7 and SL-MPP5 schemes is also significant.

SL-MPP7

SL-MPP5

PFC

Numerical Radiation Transfer & Radiation Hydrodynamics

Radiation transfer is one of the most important processes in astrophysics since the radiation is the most efficient way to transfer energy in the universe. Hence, emission, absorption, and scattering of radiation in various sites of astrophysical objects such as stars, QSOs and galaxies is of great interest both theoretically and observationally. A numerical simulation of radiation transfer (RT) is a quite challenging issue in numerical astrophysics since it has to solve the distribution function of photons in six-dimensional phase space volume. We develop a simulation code to simulate the RT in the inter-galactic medium (IGM) based on the ray-tracing scheme, in which numbers of light rays are casted and the radiation transfer equation is solved along them. Since the radiation transfer equation along each light rays can be solved independently, it has a large degree of concurrency. We exploit the capability of processors with a highly parallel architechture such as graphics processing units (GPUs), which are able to perform numbers of threads in parallel, to accelerate RT simulations. Generally, ray-tracing calculations of RT are categorized into two types; one solves the transfer of radiation emitted by point-like radiation sources such as stars, QSOs, and the other handles the diffuse radiation emitted by spatially extended sources. In the former case, we have the light-rays originating from point-like sources and thus the computational costs is usually proportional to the number of point-like sources. The RT code we developed, "ARGOT" accelerates RT calculations with many radiating point-like sources by adopting the tree structure to handle radiating sources.