Modeling the multiwavelength spectra and images of hot accretion flows around black holes requires accurate treatment of key radiative processes within a general relativistic framework. To meet this challenge, we present Odyssey and Iliad, two GPU-accelerated general relativistic radiative transfer (GRRT) codes that simulate photon propagation and interaction in the curved spacetime of rotating (Kerr) black holes. The two codes implement complementary ray-tracing schemes: Odyssey adopts an observer-to-source approach, ideal for interpreting mm/sub-mm VLBI observations of supermassive black holes such as Sgr A* and M87, while Iliad employs a source-to-observer Monte Carlo method that naturally incorporates scattering, enabling exploration of X-ray spectra and high-energy emission. In this talk, I will introduce the astrophysical foundations of GRRT and showcase illustrative examples of black hole imaging and spectral modeling. Together, these tools offer a high-fidelity, GPU-powered framework for testing theoretical models of black hole accretion and emission against cutting-edge observations.

We describe a post-processing radiative transport code for computing the spectra, the coreshift, and the surface-brightness distribution of special relativistic jets with arbitrary optical thickness. The jet consists of an electron-positron pair plasma and an electron-proton normal plasma. Electrons and positrons are relativistic and composed of thermal and nonthermal components, while protons are non-relativistic and non-radiating. The fraction of a pair plasma, as well as the fraction of a nonthermal component can be arbitrarily chosen. Only the synchrotron process is considered for emission and absorption when the radiative-transfer equation is integrated along our lines of sight. We describe a suite of test problems, and confirm the frequency dependence of the coreshift in the Konigl jet model, when the plasma is composed of nonthermal component alone. We illustrate the capabilities of the code with model calculations, demonstrating that the jet will exhibit a limb-brightened structure in general when it is energized by the rotational energy of the black hole. It is also demonstrated that such limb-brightened jets show a ring-like structure in the brightness map when we observe the jet launching region nearly face-on. Applying the code to the M87 jet, and comparing with multi-wavelength observations, we successfully constrain the composition of the jet, jet-launching altitude, collimation index of the magnetic field, and the magnetic field strength at the jet base.

Growth of massive black holes (BHs) in the galactic centers are regulated by the environment. Modern cosmological galaxy-formation simulations suggest that supernova (SN) feedback evacuates the gas in galactic center, suppressing the BH growth until the host galaxies have grown sufficiently to develop a deep gravitational potential, leading to under-massive growth track relative to the local relationship. However, this scenario does not explain the over-massive nature of BHs observed at high redshift through JWST. In this work, we perform a suite of 3D high-resolution hydrodynamical simulations that investigate the properties of turbulent, multi-phase gas driven by individual SN explosions, and the dynamics of accreting gas onto a BH through its gravitational influence radius. We explore a broad parameter space of the BH mass (~ 10^4 - 10^7 Msun), density of the surrounding gas (~ 1-10^5 cm^-3), and frequency of explosions (given by star-formation timescale, tau ~ 10-10^4 Myr). When the density in the nucleus is as high as > 10^3 cm^-3 (tau / 10^2 Myr)^(-2), where the volume filling factor of SN bubbles within the BH influence radius is less than 0.1, the BH is fed at a high rate comparable to the Bondi accretion rate by dense cold gas formed between SN bubbles. This result, unlike most large-scale galaxy simulations that hardly resolve the nucleus, suggests that SN feedback is inefficient to expel the gas and prevent the BH from growing. These high-resolution simulations enable us to provide a physically motivated subgrid feedback model, which can be applied to large-scale simulations.

Relativistic jets from active galactic nuclei are observed to be collimated on the parsec scale. When the pressure between the jet and the ambient medium is mismatched, recollimation shocks and rarefaction shocks are formed (e.g., the HST-1 in the M87 jet). In this work, we perform 3D Relativistic MagnetoHydrodynamic (RMHD) simulations for over-pressured magnetized jets, showing that the current-driven kink instability grows after passing the recollimation shocks. To explore the polarized and time-domain features, we develop the special relativistic radiative transfer code RaptorP, which self-consistently incorporates synchrotron emission, absorption, and Faraday rotation. Cooperating with constraints for non-thermal electrons inspired by recent Particle-In-Cell (PIC) simulations, we reveal that the jet image differs at different frequencies. The turbulence is observed at low frequencies, while the inner structure appears at high frequencies, including bright knots and twisted kinks. We further compute the multi-wavelength spectrum energy distribution (SED) and catch the Quasi-Periodic Oscillation (QPO) signal in light curves.

Multiple nearby Active Galactic Nuclei have been reported as sources of high-energy neutrinos detected by the IceCube observatory. The results strongly suggest efficient proton acceleration to (sub‑)PeV energies, likely within Black Hole (BH) Coronae, given the lack of γ-ray counterparts. The acceleration mechanisms remain unconfirmed due to limited understanding of coronal environments and challenges in modeling hot, relativistic plasmas. While diffusive shock acceleration (DSA) has been proposed, a self-consistent treatment incorporating detailed kinetic plasma effects has been absent. In this study, we present the particle-in-cell (PIC) method as a first-principles approach to investigate particle acceleration by collisionless shocks. Using large-scale 1D simulations, we surveyed possible shock parameters inferred from multi-wavelength observations of BH Coronae, focusing on previously underexplored effects, such as initial temperature ratios between ions and electrons and trans-relativistic shock velocities. We found that collisionless shocks can accelerate protons to the energies required for observed neutrino spectra under a wide range of plasma conditions. The simulations also provide insights into the energy partition between protons and electrons, offering crucial constraints on coronal plasma conditions verified by observational counterparts in radiation.

Gravitational dynamics describes a wide range of structure of the Universe. Because gravity is scale free, gravitational structures of interest can form at small scales (e.g., binary stars with separation 1e-9 light years) to large scales (e.g., galaxy clusters 10^8 light years). Gravitational structure simulations of the Universe must be able to resolve a wide range of scales while remaining computationally efficient, which poses a formidable problem. Any efficient simulation of structure formation in the Universe must be adaptive in time and/or length.

Moreover, differential equations of gravity are approximated to excellent degree by a system of first order differential equations possessing symmetries of the continuous, differential, and discrete kind. Gravitational simulations, in some cases having to integrate billions of orbits, are prone to catastrophic error unless their underlying numerical methods contain knowledge of these symmetries. Such numerical methods are known as geometric integrators.

I present new geometric integration methods for gravitational dynamics that exploit the Hamiltonian and time-reversible nature of gravity. These methods have been implemented in several open source packages, including the popular Rebound package. They have been used for many purposes including constraining the nature of planets in exoplanetary systems like TRAPPIST-1.

Self-gravity plays a crucial role in the formation of cosmic structures, and many hydrodynamic simulations have incorporated its effects. The multigrid method is commonly used to solve the Poisson equation for self-gravity, but its parallel efficiency declines on massively parallel computers (> ~10⁴ cores), even with well-optimized codes. To address this, we have developed a new self-gravity solver based on the telegraph equation with a damping coefficient, κ. Unlike the elliptic Poisson equation, the telegraph equation is hyperbolic, making it significantly easier to parallelize. Through convergence tests, we determined the optimal non-dimensional damping coefficient (κ~2.5), and weak-scaling tests confirmed that our method maintains high parallel efficiency even on large-scale systems. When the simulation timestep is constrained by heating, cooling, or chemical processes rather than by the CFL condition, our method can outperform conventional solvers such as FFT and multigrid methods. This is because the gravitational phase velocity allowed by the CFL condition in such regimes can be much larger than the fluid velocity plus sound speed. Finally, we will present examples of astrophysical applications using this solver, demonstrating its practical effectiveness.

A key goal of the Milky Way galaxy simulation is to achieve sufficient resolution down to the level of individual stars. However, the scaling of massively parallel high-performance computing fails due to some small-scale, short-timescale phenomena, such as supernova explosions. We have developed a novel integration scheme for N-body/hydrodynamics simulations working with machine learning. This approach bypasses the short timesteps caused by supernova explosions using a surrogate model, thereby improving scalability. With this method, we reached 300 billion particles (1 Msun per star and gas particle) using 148,900 nodes, equivalent to 7,147,200 CPU cores of Fugaku supercomputer. We also performed a simulation with 1,000 GPU nodes of Miyabi supercomputer. To evaluate our new scheme, we compared the star formation rate and mass loading factor of our new simulation with those of a simulation that integrates with variable timesteps for a 1/100 Milky Way disk model and confirmed that the global features of the disk are equivalent. Our new scheme was 4 times faster than the conventional simulation, and this difference will become larger for Milky Way size galaxy simulations.

We introduce a robust and accurate method for solving the Poisson equation in spherical polar coordinates using a logarithmic radial grid. Our algorithm follows a hierarchical divide-and-conquer strategy, dividing the 2D grid after applying FFT into multiple grid levels. Starting from the smallest grid level, we compute the zero-boundary potential and then progressively combine them to obtain the global solution. The boundary potential is calculated using discrete Green’s functions. All these computations can be efficiently expressed as matrix-vector operations using precomputed kernel matrices, and the computational complexity of these calculations is O(N^3log⁡N), where N is the number of cells per dimension. The method is implemented in the FARGO3D MHD code, and its computational performance and second-order accuracy are verified through a series of numerical tests.

The particle acceleration process is intrinsically multi-scale in nature, which is expected to continue from microscopic all the way into global scales. However, PIC simulations, which requires resolving the microscopic plasma scales and accommodating the speed of light, are severely limited in their ability to simulate the acceleration process towards global scale, especially in the non-relativistic regime. We aim to develop numerical methods that enable the study of non-thermal particle acceleration on MHD scales while capturing essential physics. Here, we present the formulation, algorithm, and numerical tests of the ​MHD-PIC method with particles treated under the guiding center approximation, which we term the ​MHD-gPIC method, implemented in the ​Athena++​​ MHD code. The ​MHD-gPIC​ model consists of a thermal (cold) fluid component and non-thermal particles, where the latter are evolved via guiding center dynamics, including drift motion, with carefully evaluated momentum and energy source terms for particle backreaction. By incorporating a subgrid pitch-angle diffusion model and a particle injection recipe, we further study particle acceleration in non-relativistic magnetic reconnection, demonstrating the method’s capability in capturing key acceleration mechanisms.

The Particle Mesh (PM) method is widely used in cosmological N-body simulations, often combined with the tree method. It approximately solves the gravitational potential on a uniform grid using the fast Fourier transform (FFT). Increasing the FFT parallelism (number of parallel axes) is considered optimal for accelerating PM calculations. However, depending on the FFT library employed, specific parallel data structure such as slabs, pencils, or cubes are required. On the other hand, due to the nature of gravity, the best load balance of the tree method is achieved by a highly non-uniform domain decomposition, which is significantly different from what is required by FFT. Resolving this mismatch requires additional calculation costs for data transformation and the associated communication, which could be a bottleneck in extreme scale simulations. Therefore, it is not clear whether increasing the FFT parallelism improves overall performance, and the performance evaluation taking this cost into account is necessary. In this study, we implement the PM method using multiple FFT libraries with different FFT parallelism and compare the overall performance to identify the optimal choice, using a large number of nodes of the Japanese flagship supercomputer Fugaku and particles.

Machine learning techniques offer new approaches for astronomers to efficiently explore and extract novel astrophysical information from large datasets. In this talk, I will show two examples of utilizing supervised and unsupervised learning techniques to explore the massive spectroscopic dataset provided by the Dark Energy Spectroscopic Instrument (DESI) project. First, with a supervised learning technique, we identify hidden Lyman alpha emitters (LAEs) at redshift>2 from millions of DESI spectra and construct a catalog of >10,000 Luminous LAEs, enabling follow-up studies of the physical properties of LAEs. Second, with an unsupervised learning technique, we show that high-dimensional DESI galaxy spectra can be described by a few coefficients, linking to galaxy physical properties, such as stellar mass, star-formation rate, line-ratio etc. This allows us to explore galaxy population in low-dimensional space, providing a machine-learning perspective on galaxy classification.

While galaxy formation is fundamentally governed by dark matter halo properties, the environment surrounding galaxies plays an equally crucial role in shaping their physical properties. Understanding how these environmental effects may vary with redshift provides important insights into galaxy formation and evolution processes. We construct an interpretable neural network framework to characterize the surrounding environment of galaxies and investigate the extent to which their properties are affected by neighboring galaxies in IllustrisTNG300 cosmological magnetohydrodynamic simulation data. Our models predict galaxy properties (stellar mass and star formation rate) given dark matter subhalo properties of both the host subhalo and surrounding galaxies, serving as an explainable, flexible galaxy-halo connection model. We present detailed results at z=0 and explore the potential for extending this analysis to higher redshifts (z=1, 2, 3). At z=0, we find that environmental influence follows a clear hierarchical pattern, with the nearest neighbor providing the dominant contribution. Environmental dependence varies significantly with galaxy type and mass: low-mass galaxies show 35-50% environmental contribution compared to just 8-30% for massive centrals, while satellite galaxies experience consistently stronger environmental effects than centrals. We will also present preliminary investigations into how these environmental dependencies may change at higher redshifts, exploring the framework’s applicability to z=1, 2, 3 snapshots and discussing the potential for understanding redshift evolution of environmental effects. These quantitative results offer guidance for constructing more sophisticated empirical and semi-analytic models of galaxy formation that explicitly include environmental dependence as a function of galaxy type and mass, with future potential for incorporating redshift-dependent effects.

Stellar spectroscopic analysis, though a mature and industrialized field supported by stellar-atmosphere radiative transfer theory, now faces new challenges. The abruptly increasing volume of observational data requires more efficient analysis methods, but recent data-driven methods suffer from unrealistic, non-physical features. We present Kurucz-a1, a physics-informed neural network (PINN) that emulates 1D stellar atmosphere models under Local Thermodynamic Equilibrium (LTE), addressing a this critical bottleneck in stellar spectroscopy. By incorporating hydrostatic equilibrium as a physical constraint during training, Kurucz-a1 creates a differentiable atmospheric structure solver that maintains physical consistency while achieving computational efficiency. Kurucz-a1 can achieve superior hydrostatic equilibrium and more consistent with the solar observed spectra compared to ATLAS-12 itself, demonstrating the advantages of modern optimization techniques. Combined with modern differentiable radiative transfer codes, this approach serves as the first step on data-driven optimization of universal physical parameters across diverse stellar populations-a capability essential for next-generation stellar astrophysics, and presenting as an ideal example for machine learning applications in astrophysics.

Mounting evidence has shown that solar eruptions generate two waves, i.e., a fast-mode MHD wave (or shock wave) followed by a slower non-wave component, which is sometimes called solar EIT waves. However, EIT wave is missing in some eruption events, and even when an EIT wave is present, the intensity distribution often exhibits strong anisotropy. This study is intended to unveil the formation condition of the solar EIT waves. With data-inspired MHD simulation, we propose that backward-inclined magnetic field lines favor the appearance of EIT waves. The more the field lines are forward-inclined, the weaker such wavelike fronts become.

Solar flares are the largest energy releasing events in the solar system, where the open magnetic field lines have reconnected and form the closed flare loops. During this process, magnetic reconnection proceeds rapidly with various plasma dynamics, such as shock waves and chromospheric evaporation. Superflares have been also observed for many stars and astrophysical jets, which have common features to solar flares. However, the plasma dynamics in the flare loops remain unclear owing to the lack of a comprehensive three-dimensional (3D) modeling and the corresponding observational evidences. Here we present a 3D high-resolution magnetohydrodynamics simulation of solar flares. In this study, it has been found that slow-mode shocks are recurrently formed after the collision between the post-reconnection downflows and evaporation flows in the flare loops, which subsequently propagate toward the chromosphere at speeds comparable to the evaporation flows. By comparing modeling results with spectrum observations of an M-class solar flare, we identified these slow shocks in evaporation flows. The profiles of the Fe XXI line exhibit two well-separated components on the arm of flare loops, a highly blueshifted component representing the hot evaporation flows, another small blueshift component accounting for the disrupted flows across the slow shocks. The detection of such a two-component blueshift between the flare loop-top and footpoints by spectrum observations indicates the existence of propagating slow shocks. These shocks remarkably facilitate energy release in flares and affect energy transport, suggesting a significant advancement in the standard flare model framework.

Understanding the complex and diverse physical processes that drive galaxy evolution requires the analysis of high-dimensional, multiwavelength observational datasets. Contemporary galaxy surveys routinely collect photometric and spectroscopic data across numerous bands, encompassing information related to stellar mass, star formation rate, dust attenuation, and gas content. While rich in information, the high dimensionality of such data presents substantial challenges for physical interpretation, direct comparison with theoretical simulations, and integration into machine-learning frameworks.

To address these challenges, we apply dimensionality-reduction techniques to galaxy photometry measured in 15 bands, ranging from the far-ultraviolet to the far-infrared. We use Principal Component Analysis (PCA) for linear dimensionality reduction and Uniform Manifold Approximation and Projection (UMAP) to capture nonlinear structure in the data. PCA reveals that the first three components alone account for 93 % of the total variance, and reconstruction-error analysis from UMAP similarly indicates that a three-dimensional representation is sufficient to preserve the essential structure of the dataset. Each of these three principal components can be physically interpreted as tracing stellar mass, star formation activity, and dust content, respectively.

This three-dimensional embedding defines a physically meaningful, data-driven “galaxy manifold.” In this presentation, we will explore how individual galaxies move along this manifold, providing a framework for interpreting galaxy evolutionary trajectories. We will also demonstrate how this manifold can serve as a powerful diagnostic tool and quantitative benchmark for evaluating the performance of numerical simulations and semi-analytic models of galaxy formation and evolution.

The angular momentum of molecular cloud cores plays a key role in the star formation process. However, the evolution of the angular momentum of molecular cloud cores formed in magnetized molecular filaments is still unclear. Observations show that the strength of the magnetic field in the filaments are relatively strong and that some of the cores are marginally supercritical. These observations suggests that it is important to investigate the physical properties of cores formed in strongly magnetized molecular filaments. We perform 3D magnetohydrodynamics simulations to reveal the effect of the strong magnetic field on the evolution of the angular momentum of molecular cloud cores formed through filament fragmentation. Our results show that the initial angular momentum of the core with mass of one solar mass is transferred during their formation phase. If the initial magnetic field is weak, the rotation direction of the cores is perpendicular to the filament longitudinal axis due to the effect of the filament geometry. On the other hand, if the magnetic field is strong, the signature of perpendicular rotation becomes less clear due to the strong magnetic braking. We also investigate the rotation direction with respect to the local magnetic field direction. Our results show that most of the cores experience the alignment of the rotation direction with the magnetic field, but the distribution of the angle between the rotation direction and magnetic field follows random distribution. This is because the angular momentum vector continues to oscillate even after the alignment. In addition to these results, we will also explain the impact of the non-ideal MHD effect on the core formation in the filament.

The orbital decay of binary systems is a key process in the formation of close binaries. We perform high-resolution three-dimensional magnetohydrodynamic (MHD) simulations of a binary system accreting gas from the surrounding environment. Our simulations reveal the presence of outflows and jets launched from both the circumstellar disks (CSDs) and the circumbinary disk (CBD). Magnetorotational instability (MRI) is also excited within the CBD. These magnetic processes efficiently extract orbital angular momentum from the binary, leading to orbital decay, whereas a purely hydrodynamical model shows orbital expansion. The decay rate reaches up to ~1% per orbital period, depending on the accretion rate onto each binary component. In the later stages, orbital decay slows down because the accumulation of magnetic flux suppresses accretion onto the stars; since both mass and angular momentum accretion drive the orbital evolution, the decay ceases as accretion halts. We discuss the possibility that magnetic fields can induce orbital decay and contribute to the formation of close binary systems.

The total gas mass in the galactic disk is about 10^9 M_sun and the star formation rate is approximately 3 M_sun / yr. Thus, the disk gas would be depleted within 1 Gyr, and star formation cannot continue over 1 Gyr. However, observations show star formation has continued at a nearly constant rate for about 10 Gyr. The mechanism that sustains long-term star formation remains unclear. Recent multi-wavelength observations have revealed ionized, metal-enriched gas in the galactic halo, extending over 100 kpc from the disk and possessing an estimated mass exceeding 10^10 M_sun. These findings suggest that gas inflows from the halo to the disk prevent gas depletion and sustain star formation for cosmic time. However, the detailed processes governing gas supply remain poorly understood. We investigate “thermal instability” as a key mechanism driving gas supply from the halo. Specifically, we analyze the stability of a gravitationally stratified halo gas at ~0.1 keV in thermal equilibrium, considering hydrodynamics, the galactic gravitational field, radiative cooling, and cosmic rays. Our results indicate that the halo gas is thermally unstable and that cold HI clouds can be formed Furthermore, we discuss the observability of these cold clouds as Intermediate- / High-Velocity Clouds, which are well-known HI clouds in the galactic halo with debated origins. Our talk will focus on three main topics: (1) the stability of equilibrium halo gas, (2) the observability of HI clouds as Intermediate- / High-Velocity Clouds, and (3) the contribution of thermal instability to the overall gas supply.

H II regions formed by massive stars significantly affect the structure of the circumstellar medium, potentially quenching further star formation. Their expansion is strongly influenced by the density structure of the ambient gas and the spectra of radiation from multiple stars, making it essential to account for the multidimensional distribution of gas and stellar sources. In this study, we perform three-dimensional radiation hydrodynamics simulations to investigate radiative feedback from massive star formation. We analyze the effect of clustering of massive stars on the efficiency of radiative feedback. We discuss the implications of these findings for molecular cloud destruction and resulting star formation efficiency in realistic star-forming environments.

In the early universe, including the era of the formation of the first stars, direct observations of star-forming regions are extremely challenging. Therefore, theoretical studies using numerical simulations have been the primary tool for investigating their formation processes. Historically, magnetic fields in the early universe were expected to be extremely weak, and most star formation simulations neglected magnetic effects. However, recent studies suggest that magnetic fields in the early universe could have been strong enough to significantly influence the star formation process. In the formation sites of second-generation stars, which are formed after the supernova explosions of the first stars, the magnetic field is expected to be amplified and become dynamically important. Given this context, our study focuses on star formation under strong magnetic fields and low-metallicity environments. The strength of the magnetic field is characterized by the mass-to-flux ratio, μ₀=3, normalized by the critical value (2πG)⁻¹. To reproduce star-forming environments at different cosmic epochs, we constructed six models with different metallicities of 0, 10⁻⁶, 10⁻⁵, 10⁻⁴, 10⁻³, and 10⁻² Z_sun. To accurately calculate magnetic field dissipation, we performed three-dimensional non-ideal magnetohydrodynamics simulations including both Ohmic dissipation and ambipolar diffusion. As a result, magnetically driven interchange instability was observed in all models except the one with 10⁻⁴ Z_sun. Approximately 100 years after protostar formation, we measured the ratio of magnetic field strength at 10 au to that at 1 au from the protostar (B₁₀au / B₁au). In the models where interchange instability occurred, (B₁₀au / B₁au) ≈ 1, indicating that magnetic fields were effectively transported outward. In contrast, the model without the instability showed (B₁₀au / B₁au) ≈ 0.01. These results demonstrate that in strong-field, low-metallicity environments, magnetic interchange instability plays a crucial role in transporting magnetic fields from the region close to the protostar to the circumstellar region during the early star formation process.

Relativistic jets of radio galaxies are widely regarded as potential sources of ultra-high-energy cosmic rays (UHECRs). In this study, we investigate the acceleration of UHECRs within the jets of Fanaroff–Riley (FR) type I and II galaxies through numerical simulations. To accurately capture the nonlinear dynamics of ultra-relativistic jets, we employ a newly developed high-order relativistic hydrodynamic (RHD) code. Subsequently, we follow the transport and acceleration of cosmic rays using a Monte Carlo approach. Our results reveal that the most prominent acceleration occurs around the interface between the jet spine and backflow (cocoon) via the so-called relativistic shear acceleration. Additionally, though less significant, contributions arise from shocks and turbulence within the backflow region. Based on these simulations, we present model spectra of UHECRs accelerated in FR I and II jets, which can be used to interpret UHECR observations at Earth, particularly those originating from nearby radio galaxies.

Depending on the astrophysical source and its environment, the accretion flows can exhibit a variety of behaviors and characteristics in accordance with the type of solutions. To understand how the flow changes with different angular momentum, by changing the initial condition of the accretion torus minimally, we can simulate a steady, low-angular-momentum accretion flow around a Kerr black hole. We find that an accretion flow with an intermediate range of angular momentum differs significantly from high- or very-low-angular-momentum flows. With the increase in flow angular momentum, it develops a nonaxisymmetric nature. In a certain parameter range, we also find the standing shock feature in accretion flows. We will discuss the possibility of explaining the complex observation features of the supermassive black hole Sgr A* at our galactic center.

Stars play a central role in cosmic evolution, and understanding their formation is a fundamental question in astrophysics. Radiation hydrodynamics (RHD) simulations are now essential for studying this process, and recent advances in computational resources enable simulations across vast dynamic ranges that resolve scales from star-forming clouds (~1000 au) down to stellar interiors (<0.1 au). These simulations elucidate realistic protostellar evolution under the influence of global accretion flows. To perform such simulations over long timescales, numerical schemes that solve radiative transfer at low computational cost with high parallel efficiency are crucial. The M1 closure method combined with the reduced speed of light approximation (RSLA) is one such efficient scheme. However, it breaks down in extremely optically thick regions, such as protostellar interiors.

To overcome this limitation, we extend the M1 closure method with the RSLA and develop a new RHD scheme. This scheme maintains high parallel efficiency through explicit radiation transport and is applicable to optically thick regions, including protostellar interiors. In this presentation, we will outline the proposed scheme and demonstrate its effectiveness in revealing the 3D evolution of protostars under global gas accretion. Furthermore, as an application, we investigate the evolution of primordial protostars forming in the early universe under rapid accretion. Our simulation results reveal a 3D picture of protostellar evolution that differs from the conventional 1D picture.

Observations of accreting black hole (BH) systems, such as microquasars and supermassive black holes, often reveal a precessing jet with changing directions, indicating a misaligned accretion flow relative to the BH spin. The precession is commonly attributed to the Lense-Thirring (LT) effect, which arises from the BH’s rotation twisting the surrounding spacetime and accretion flow. In the strongly magnetized regime, which is preferred accretion flow conditions for M~87∗ and likely other jet-producing systems, the large-scale magnetic field can significantly influence the flow dynamics. Here, we perform large-scale three-dimensional general relativistic magnetohydrodynamic simulations of tilted accretion onto a rotating BH, and find a never-seen-before new retrograde precession. This precession arises from a magnetic torque on the disk generated by the poloidal magnetic field aligned with the BH’s rotation, opposing the LT torque. This finding highlights the unique property of highly magnetized accretion flows around BHs and provides a new interpretation of jet precession observed in many systems.

The galactic centre supermassive black hole, Sagittarius A* (Sgr A*), is powered by a hot, magnetised, and radiatively inefficient accretion flow. Because the plasma is collisionless, the electrons are not necessarily in thermal equilibrium. Accurate modelling of electron thermodynamics is therefore crucial in interpreting horizon-scale resolution observations, particularly those from the Event Horizon Telescope (EHT) . Current EHT simulation libraries rely on simplified, parametric electron models, which may not fully capture the complexity of the physical system and could bias inferred black hole properties. To address this limitation, we implemented a new, physically motivated electron temperature model into the general relativistic radiative transfer (GRRT) code ipole, and conducted an extensive numerical study. Our results show that this improved electron model significantly alters the optical depth of the accretion flow. Specifically, it concentrates emission closer to the black hole and produces a more asymmetric, optically thinner accretion disk. This leads to enhanced flux of the photon ring, with measurable effects in both polarization and variability. I will present the numerical details of this work and discuss its implications for future EHT modeling of Sgr A*.

Galaxy formation simulations are powerful tools for studying the formation of galaxies. Current galaxy formation simulations successfully reproduce galaxies’ observed properties in many aspects. Some of the simulations have sufficiently high resolution to resolve pc-scale phenomena. However, it is still challenging to resolve individual stars in galaxy formation simulations. In this talk, I will review the recent progress of galaxy formation simulations. I will show some star-by-star galaxy formation simulation projects, including our Simulations Resolving Individual Stars (SIRIUS) project. With SIRIUS simulations, we show that ~5% of extremely metal-poor stars are polluted by a single supernova. I will also show how to accelerate the computation of galaxy formation using deep learning.

Dwarf galaxies are thought of as the building blocks of large galaxies such as our Milky Way. This talk presents new high-resolution hydrodynamical simulations of dwarf galaxies and their intergalactic medium with the GIZMO code. Our simulations consider the key physical processes of galaxy evolution, such as gas cooling, chemistry, and stellar and black hole feedback. Unlike the previous work, the initial conditions of our simulations taking the dwarf galaxies of 2 − 5 × 1010 M⊙ from the realistic cosmology simulations, IllustrisTNG. We further increase the original resolution of IllustrisTNG by a factor of ∼ 100 via a particle splitting scheme. Our results show that the evolution of complex multiphase CGM and its metal content is sensitive to the redshift of dwarf galaxies. The accretion of CGM into dwarf galaxies plays a key role in providing 20% − 50% of the star-forming gas and replenishing 40% − 70% of the total mass in the galactic disk. Furthermore, the accretion history of supermassive black holes in the centers of high-z dwarf galaxies shows episodic patterns with high-accreting states close to ∼ 10% of the Eddington mass accretion rate, implying the rapid growth of supermassive black holes in the early universe, which may be revealed by the coming observations from the James Webb Space Telescope (JWST).

In the ΛCDM universe, structure formation is generally not a self-similar process, while some self-similarity remains in certain statistics which can greatly simplify our description and understanding of the cosmic structures.In this work, we show that the merger tree of dark matter halos is approximately self-similar by investigating the universality of the subhalo peak mass function (PMF) describing the mass distribution of progenitor halos. Using a set of cosmological simulations and identifying subhalos of different merger levels with hbt+, we verify that the level-1 subhalo PMF is close to universal across halo mass, redshift, and cosmology. This approximate self-similarity allows us to analytically derive the subhalo PMF for subhalos accreted at any level (i.e., for sub-sub…halos) through self-convolutions of the level-1 PMF, and the resulting model shows good agreement with simulation measurements. We further derive a number of analytical properties on the hierarchical origin of subhalos. We show that higher-level subhalos dominate at progressively lower peak mass in the PMF and are more likely to originate from major mergers than lower-level ones. At a given merger mass ratio, the subhalo accretion rates at each level track the growth rate of the host halo. At a fixed final mass ratio, however, subhalos of higher-level, higher-mass-ratio and in more massive haloes tend to be accreted more recently. Matching subhalo peak mass to galaxy mass, these results have direct implications on the hierarchical origin of satellite galaxies.

Dark matter subhaloes are highly nonlinear structures resulting from mergers among halos in the hierarchical universe and serve as a key link between galaxies and the underlying density field. Their formation and evolution involves complex physical processes and simulation studies of them can be affected by severe numerical effects. Despite this, many universal yet simple laws have been discovered describing the distributions of subhalos, including their initial and final mass functions and spatial distributions. Validating the fidelity of these laws against numerical effects and extending them to different cosmologies require physical understanding of their origins. In this talk, I will introduce an extremely simple idea to unify the understanding of these universal laws, assuming that subhalo progenitors are accreted unbiasedly from the surrounding environment. I will also introduce the extension of subhalo distribution laws to different cosmologies, different dark matter species and different levels of subhalos, and respond to recent challenges on artificial disruptions of subhalos according to our physical understanding of subhalo evolution.

Heterogeneous computing is revolutionizing the landscape of astrophysical simulations, necessitating the development of well-optimized systems to fully harness the potential of modern computing capabilities. In this presentation, we will delve into KRATOS, a state-of-the-art high-performance heterogeneous simulation system developed by the speaker. KRATOS is distinguished by its meticulously crafted one-sided coupling mechanism that integrates seamlessly with geometrical and physical modules, culminating in a sophisticated infrastructure that offers unparalleled versatility and flexibility, which would be especially helpful in dealing with multi-scale astrophysical simulations. We will explore the applications of KRATOS in the realm of astrophysical simulations and touch upon its potential for other research domains, showcasing its adaptability and the breadth of its utility.

Hydrodynamic instabilities likely operate in protoplanetary disks. One candidate, convective overstability (COS) or radial convection, can develop in particular disk regions. The ensuing turbulence and flow structures are expected to affect planet formation. We present high-resolution simulations of the COS using local Boussinesq shearing box and global compressible cylindrical models with and without embedded dust grains. Under axisymmetry, we show that the COS generates zonal flows that trap dust; however, unsteady flows and dust feedback effects can limit dust concentrations. Global 3D simulations reveal that COS induces outward mass transport and large-scale vortex formation that significantly modifies the disk structure. Furthermore, these COS-produced vortices efficiently concentrate dust, particularly for larger grains at solar and supersolar metallicities. These findings highlight the COS’s role in shaping disk structures and dust distributions, offering insights into observed dust rings and asymmetries in protoplanetary disks, with implications for planet formation processes.

Turbulence in protoplanetary disks affects the evolution of dust and the formation of planetesimals. One possible source of turbulence is the vertical shear instability (VSI), which can develop when gas cools quickly. Interestingly, dust particles not only respond to VSI-driven turbulence but also influence the cooling and potentially regulate turbulence itself. As a result of this interaction, it remains unclear how turbulent structures and dust distributions form and persist. To investigate how the VSI sustains turbulence while balancing between dust settling and diffusion, we perform global two-dimensional two-fluid hydrodynamical simulations of an axisymmetric protoplanetary disk, including dust as a pressureless fluid. Our simulations capture the dynamic interplay between dust distribution, cooling rates, and gas turbulence. We find that VSI mixing, dust settling, and local dust cooling reach a quasi-equilibrium, forming a thick dust layer with a dimensionless vertical diffusion coefficient of approximately 10^{-3}. This balance also depends on the dust size and dust-to-gas mass ratio. Larger grains or lower mass ratios weaken turbulence, leading to dust settling. The condition for equilibrium is consistent with the prediction of the semi-analytic model presented by Fukuhara & Okuzumi (2024). Our results suggest that the VSI activity can influence dust growth efficiency in protoplanetary disks.

In planet formation, multiple characteristics, or components, of dust aggregates other than mass are highlighted. However, numerical simulations that calculate the coagulation between multi-component dust aggregates are 1. computationally high-cost if calculated directly, and 2. contradictory with observational results for porosity if calculated approximately (Okuzumi et al. 2012; Kataoka et al. 2013; Zhang et al. 2023). Thus, we need a numerical multi-component coagulation method with fewer approximations to compose a finer planet formation model.

We propose a novel multi-component dust coagulation algorithm using a tree algorithm from N-body gravity simulations. The distribution of dust aggregate mass and volume is represented by N^d bins, where d is the number of components and N is the number of bins per component. In the direct method, all interactions between these bins are calculated, ending up in the computational complexity of O(N^(2d)). On the other hand, in our tree method, we analogically assume mass ratio and volume ratio between coagulating aggregates as “distance”, and by averaging the mass and volume for the bins for smaller aggregates, we can reduce the number of interactions and the computational complexity to O(d N^d log N).

To evaluate the performance of our novel method, we compared the simulation results of the direct method and the tree method, and analytical solutions in terms of computational error and time needed to compute. When the number of components is d=1, the computational time was almost equivalent for the direct and the tree method. However, when the number of components is d=2, the tree method outperformed the direct method in computational time by orders of magnitude. The numerical error was worse for the tree method, but by adjusting the parameters, the error could be improved for the tree method.

With our novel tree-based coagulation algorithm, we are now able to numerically calculate the coagulations between multi-component aggregates with fewer approximations. We plan to implement the physical quantities and processes in protoplanetary disks and calculate the dust growth with porosity, revising the planet formation models.

Gravitational instability (GI) remains a promising yet debated pathway for giant planet formation in protoplanetary disks (PPDs), especially at wide orbital separations or around low-mass stars where core accretion faces significant challenges. However, two key uncertainties remain: under what conditions GI leads to disk fragmentation, and whether the resulting fragments can form gas giants without becoming overly massive. To address these questions, we conduct a suite of global three-dimensional radiation hydrodynamics (RHD) simulations of self-gravitating PPDs using the meshless finite-mass (MFM) method. Our simulations systematically vary disk mass and opacity, solving the radiation transport equations via the M1 closure scheme. We show that both increasing disk mass and lowering opacity promote fragmentation by enhancing radiative cooling. Non-fragmenting disks settle into a gravito-turbulent state with low-order spiral structures and effective angular momentum transport characterized by $\alpha \sim \beta_\mathrm{cool}^{-1}$. In fragmented runs, we identify and track gravitationally bound clumps and find that a subset survives as long-lived fragments. The initial masses of these fragments range from $\sim0.3–10,M_\mathrm{J}$, consistent with the formation of gas giants. Normalizing by an analytic mass estimate based on local disk properties yields a universal, dimensionless distribution across simulations, suggesting that fragment masses are governed by local gas conditions. Our results demonstrate that under favorable cooling conditions, GI can produce planet-mass fragments in realistic RHD simulations. This reinforces GI as a viable formation channel for gas giants and lays the groundwork for future studies on fragment survival, migration, and observational signatures.

Planetary systems often host debris disks. Recently direct imaging observation by JWST has detected a Saturn-mass giant planet and a debris disk, but in the future, it is expected that Earth-like planets and debris disks will be detected using infrared imaging. If the dust around a planet is dense and uniformly distributed, the thermal emission from the dust obscures the radiation from the planet. On the other hand, if a characteristic dust distribution formed by planetary gravity is clearly observed, it could serve as a clue suggesting the presence of a planet. Therefore, when detecting planets through direct imaging, it is necessary to accurately consider the dust distribution around the planet, i.e., the structure and radiation of the debris disk. We consider planetary systems with debris disks similar to the solar system, where the zodiacal dust distributes around the Earth. We investigated the dust distribution formed around an Earth-mass planet with circular orbit at 1 au. Dust particles drift to the planets from a planetesimal belt at 2-3 au by the Poynting-Robertson effect (P-R effect). Based on the orbital evolution of the planet and dust particles, we obtain the steady-state dust distribution around the planet via the following way. We divide the disk into radial and azimuthal meshes and calculate the time each dust particle spends in a mesh. By integrating the time, we obtained the number density distribution of dust within the disk. Next, we evaluated the infrared radiation flux ratio of the debris disk relative to the stellar radiation. In the planetesimal disk, destructive collisions happen for large dust particles prior to P-R drift, while small dust particles can drift into the region around the planet. Therefore, we estimated the dust size distribution based on the collision time in the planetesimal belt and the drift time due to the P-R effect. According to the dust size distribution, in bright debris disks, small dust particles are dominant, which produce structurally smooth disks. It is difficult to detect a planet embedded in such a smooth, bright disk. On the other hand, in dark debris disks, larger particles are dominant, which forms characteristic structures around a planet. Such disks may assist the detection of exoplanets.

Protoplanets can interact with their natal disks and generate gas and dust substructures such as gaps and rings. However, how these planet-induced substructures affect the disk temperature, and how that in turn influences the substructures, remains unclear. We aim to study disk substructures and the thermal structure self-consistently and explore their impact on volatile distribution. To this end, we perform iterative multi-fluid hydrodynamical and radiative transfer simulations of planet-disk interactions. We find that the temperature in a structured disk deviates significantly from that of a smooth disk due to giant planet formation. In particular, midplane temperatures in gaps can increase by tens of Kelvin, leading to volatile sublimation as well as radial shifts and multiplication of icelines. Comparing our multi-dust models with previous gas-only models, we find that the former produces slightly shallower gaps and temperatures about 10 K (~25%) higher. Furthermore, the temperature at dust rings formed by pressure bumps can drop by several Kelvin, creating volatile freeze-out regions. Nevertheless, the overall midplane ice distribution is not strongly sensitive to whether dust is included. We also investigate the effect of varying disk viscosity. Increasing $\alpha$ viscosity from $10^{-4}$ to $10^{-2}$ leads to a roughly 10 K (~25%) warmer midplane due to enhanced vertical dust mixing. However, higher viscosity suppresses gap opening and reduces the temperature enhancement within gaps. As a result, iceline locations do not follow a simple trend with viscosity. Finally, we propose an observational strategy using ALMA to test our predicted temperature changes within disk gaps.

Elastic dynamics plays a fundamental role in understanding the mechanical behavior of solid bodies, with applications ranging from planetary science to engineering development. Numerical methods, especially particle-based approaches, have become indispensable tools for simulating such phenomena. In this study, we develop a novel Lagrangian method that incorporates two kinds of forces: one is similar to pressure in Smoothed Particle Hydrodynamics (SPH), and the other describing attractive force. Previous extensions of SPH to elasticity solver typically rely on the time evolution of stress tensors, which not only increases computational cost but also leads to violations of angular momentum conservation due to non-central inter-particle forces. Our new approach resolves these issues by expressing the deviatoric stress using central forces, thereby eliminating the need for stress evolution equations and ensuring strict conservation of total angular momentum. This simplification dramatically reduces computational overhead while preserving physical accuracy. The proposed method is particularly effective in simulating dynamic elastic events such as impacts and deformations. Through numerical experiments, we demonstrate the robustness and efficiency of our approach.

Understanding planetesimal collisions is important for elucidating the fundamental physical mechanisms involved in planetary formation. Prior work has typically modeled planetesimals as monolithic rocky bodies with elastic material properties (e.g., Sugiura et al., 2018), an approach that has successfully reproduced the morphology of many observed asteroids. However, growing planetesimals are unlikely to have experienced thermal metamorphism. Therefore, modeling them as monolithic rocks may fail to capture their actual internal structure and mechanical behavior. Instead, it is necessary to assume a structure more similar to loosely bound dust aggregates, which better reflects their primordial, unprocessed nature. In this study, we implement the equation of state for dust aggregates proposed by Tatsuuma et al. (2019) into a Smoothed Particle Hydrodynamics (SPH) code to simulate planetesimal collisions. Each SPH particle is treated as a volume element of a dust aggregate. This approach allows us to reproduce the bulk stress response of dust aggregates without resolving individual dust particle interactions. Using this method, we carry out a systematic study of planetesimal collisions through numerical simulations. We employ the Friend-of-Friend (FoF) algorithm to identify collision fragments and quantitatively evaluate their masses. By varying parameters such as impact velocity, impact angle, and mass ratio, we investigate the dependence of collision outcomes on these parameters and the associated impact energy. These simulations provide new insights into how the internal structure of primordial planetesimals influences collisional outcomes, offering a more realistic framework for understanding the early stages of planetary accretion.

The giant impact hypothesis is a widely supported scenario for the origin of the Moon(Hartmann & Davis 1975; Cameron & Ward 1976). According to this hypothesis, the Moon formed from the accretion of a circumterrestrial disk generated by the collision of a Mars-sized celestial body with the proto-Earth. The circumterrestrial disk, which forms within the Roche limit and cools and solidifies from hot gas, may consist of very small rocky particles. The circumterrestrial disk has been treated as a particle system and investigated using N-body simulations. Sasaki & Hosono (2018) modeled the circumterrestrial disk with particles ranging from 10,000 to 10,000,000 and investigated the effect of the numerical resolution on the luner accretion process. They concluded that the numerical resolution affects both the growth rate and the evolution of the disk structure. Therefore, to understand the numerial resolution on the luner accretion process is crucial.

In this presentation, we report the results of an investigation of the evolution of spiral structures in the Moon-forming disk and a stability analysis. We performed N-body simulations of disk with the number of particles ranging from 100,000 to 10,000,000. It was confirmed that leading spirals appear in the inner region of the disk only in high-resolution simulations, which is consistent with Sasaki & Hosono (2018). In comparison with the results of Salo et al. (2001), these leading spirals are suggested to form through viscous overstability. The leading spirals that form are non-steady structures, and radial structures develop within them. Furthermore, through the Toomre stability criterion analysis, it was found that gravitational instability intensifies along the leading spirals.

These leading spirals are generated by viscous overstability and, under the influence of gravitational instability, develop internal radial structures. Due to such complex structures in the moon-forming disk, the growth rate in the 10^7 simulation is slower compared to the 10^5 simulation.

We simulate asteroid collisions using the Material Point Method (MPM) to study the effects of impact velocity and angle on the shapes of the largest remnants. The MPM is a Lagrangian method that uses material points and background grids to handle collisions and deformation. Compared to the Smoothed Particle Hydrodynamics (SPH), it does not require neighbor searches and can easily handle boundaries and external forces. Based on previous SPH simulations, we simulate collisions between equal-mass asteroids with 50 km radii, varying the impact velocity from 50 to 400 m/s and the impact angle from 5 to 45 degrees. We use a GPU-accelerated MPM framework that incorporates self-gravity using a Barnes-Hut tree. In addition to the standard MPM, we use the Tillotson equation of state, initiate cracks based on a the Weibull distribution, and model the evolution of damage following Benz & Asphaug (1995). Damaged material is treated as Drucker–Prager granular media. These methods are widely used in modeling asteroid collisions. In our simulations, various shapes of remnants such as bilobed, spherical, flat, elongated, and hemispherical are generated depending on the collision conditions. These results help explain the irregular shapes observed in the main belt. We also compare the results of our MPM with those of previous studies using the SPH.

Meteorite impacts play a crucial role in the atmospheric evolution of young terrestrial planets. When a meteoroid enters the atmosphere, it generates high-temperature, high-pressure regions that disturb atmospheric flow. To study these complex dynamics, we performed high-resolution three-dimensional fluid simulations using Athena++. We focused on a meteoroid 10–100 meters in diameter, which typically breaks up mid-atmosphere, releasing energy and forming a chimney-like rarefied region along its trajectory (Shuvalov & Trubetskaya, 2007). These two regions, energy release and rarefied zone, significantly affect atmospheric structure. We modeled both regions within a background atmosphere based on Earth’s temperature profile. As a result, we found that expansion waves propagate more rapidly through the rarefied region, transporting lower atmospheric layers upward. To validate our model, we have simulated the 2013 Chelyabinsk meteor event using observational data (Popova et al., 2013). The simulated peak pressure distribution was consistent with the observed damage patterns, confirming the accuracy of our model. We have also investigated how energy partitioning (thermal vs. kinetic), energy release geometry, and atmospheric temperature structure influence the outcomes. This analysis allowed us to constrain key parameters governing the atmospheric response to meteoroid disruption.

The CompAS program at ASIAA is an incubator for next-generation codes that meet the demanding physics and numerical requirements of modern astrophysical computations. The developmental efforts underpin the dynamic simulations of physical processes encompassing star-forming processes and beyond, which impact a wide range of scientific areas. We build on the mission and vision established in CompAS, augmented by the science-driven physics development in CHARMS, utilizing the existing scientific and technical strengths. These flagship developments range from basic Newtonian physics surrounding newborn stars to general relativity surrounding black holes. Each in-house development addresses challenges that occur in real astrophysical systems. Future developments may enable the incorporation of the rapidly evolving trends in machine learning and artificial intelligence. We showcase some of the representative science and technical areas.

Recent observations discovered primordial Galactic stars that have very low metallicities, but their formation has not been understood clearly yet. One of formation scenarios suggests that these stars formed when the Milky Way was very young and have been intact since their formation. We elaborate this scenario further by using equilibrium fractions of primordial chemical network. In particular, we estimate the fraction of molecular hydrogen which would eventually evolve to form primordial stars. We also make use of cosmological simulations, e.g., Illustris TNG50, to figure out the conditions for the star forming regions within the Milky Way such as gas density and temperature. Radiation coming from the neighboring galaxy can also affect the formation of molecular hydrogen via photo-ionization and information on this background radiation source is also extracted from the cosmological simulations. By combining the cosmological simulations and the equilibrium fractions of molecular hydrogen, we propose a formation scenario for primordial Galactic stars by investigating the optimal conditions for star forming clouds in the young Milky Way.

Angular momentum transport by magnetic fields is important for formation and evolution of protoplanetary disks. The effects of magnetic fields are suppressed due to non-ideal magnetohydrodynamic (MHD) effects such as ambipolar diffusion and Ohmic dissipation, which depend on the degree of ionization. Cosmic rays (CRs) are the primary source of ionization in star-forming clouds, and their distribution is nonuniform as it is affected by gas density and magnetic fields. Therefore, CRs, magnetic fields, and gas interact with each other. In this work, we develop a new fully implicit cosmic ray transport module in Athena++ and perform three-dimensional simulations of disk formation from collapse of molecular cloud cores. Since CRs are strongly attenuated in the dense gas at the disk scale, distribution of magnetic fields is considerably altered compared to conventional models assuming a uniform ionization rate. While the total magnetic fluxes accreted onto the disks remain similar as the gas outside the disks remain sufficiently ionized and well coupled, the magnetic fields in the disks are less twisted due to the stronger non-ideal MHD effects. As a consequence, magnetic angular momentum transport is strongly suppressed at the disk scale, resulting in more gravitationally unstable disks with more prominent spiral arms. Our simulations demonstrate influence of non-uniform ionization resulting from CR transport and attenuation on the disk formation and evolution.

The collapse of singular magnetized toroids is a natural representation of an early phase in star formation, bridging the prestellar and protostellar phases of the collapse of molecular cloud cores. We revisit the collapse study of A. Allen et al., now with explicit nonideal MHD (ohmic diffusivity η) and higher resolution using a code able to cover a broader range of the magnetization parameter H_0. Galli–Shu equatorial pseudodisks form for all values of H_0 and η, and the asymptotic central mass growth rate is \dot{M}_{*}~(a^3/G)(1+H_0), where a is the isothermal sound speed, consistent with previous results and predictions. The explicit ohmic diffusivity makes the field line structure less radial than in previous work, connecting the pseudodisk more effectively to its surroundings. Matter can fall efficiently onto the pseudodisk surfaces, forming oblique shocks, where shock heating and large density gradients raise the possibility of rich astrochemistry. Pseudodisk size and structure are influenced by magnetic diffusivity. Force and velocity ratios were computed to explore the magnetic support within the pseudodisk and its induced reduction in infall velocity. Magnetic diffusivity was measured to control the strength of these effects and their location within the pseudodisk. The dependence of the field line configurations, pseudodisk structure, and velocity ratios on magnetic diffusivity has observable consequences for collapsing envelopes.

Observations of nearby star-forming regions have revealed that stars form in molecular cloud cores along the dense filaments via gravitational collapse. Some recent observations further suggest that gas inflow from the surrounding environment of prestellar cores can supply additional mass to the cores. These findings indicate that the external environment outside prestellar cores plays a crucial role in the star formation process, including protostellar mass growth. In our previous numerical study focusing on star formation inside molecular cloud cores, we performed sub-parsec scale three-dimensional numerical simulations to follow the evolution of cloud cores embedded in envelopes with different density profiles. This study showed that the mass accretion rate onto the protostar is temporarily enhanced in models with higher ambient density. Additionally, our results suggest that dense surroundings can significantly influence the evolution of circumstellar disks and the protostellar outflows. These observational and numerical findings suggest that, to understand the mechanism of stellar mass determination in nearby star-forming regions, it is necessary to trace the star formation process back to the formation process of molecular cloud cores. In particular, it is important to theoretically track the time evolution of molecular cloud cores as mass reservoirs for protostars under varying environmental conditions. In this study, we conducted parsec-scale numerical simulations using passive tracer particles to investigate the physical properties of star-forming cores as mass reservoirs for protostars. The simulations were performed under initial turbulence strengths characterized by root-mean-square Mach numbers of 2 and 10. Unlike the core identification algorithms used in observational studies, we identified and analyzed the star-forming cores by tracing the trajectories of tracer particles falling onto protostars. Our results show that the identified star-forming cores do not always coincide with the high-density regions identified as cores in observations. In clustered star-forming environments, the cores often fragment into clumpy structures due to selective accretion onto protostars. To quantify the spatial distribution of mass reservoirs, we defined a filling factor using the convex hulls of star-forming cores. We then investigated how the initial turbulence strength affects their physical properties, including the number of protostars, mass, size, and gravitational stability. These findings suggest that considering the broader core environment may be important for improving theoretical models of star formation.

Star formation is the driver of galaxy evolution. Understanding the properties of the cold, star-forming gas in the interstellar medium (ISM) is therefore of crucial importance. In this talk, I will present recent developments in high-resolution (sub-parsec) hydrodynamical simulations of the stellar feedback-regulated ISM and their predictions for chemical properties and line emissions, leveraging a novel hybrid method for ISM chemistry that captures the non-equilibrium effects of molecular hydrogen (H2). I will show that steady-state chemistry significantly over-predicts the abundances of H2 but not carbon monoxide (CO), leading to a reduced conversion factor (X_CO), especially at low metallicities where the H2 formation time becomes much longer than the dynamical time. On parsec scales, X_CO varies by orders of magnitude from place to place, primarily driven by the transition from atomic carbon to CO. Finally, I will present simulations combining ISM chemistry and dust evolution (sputtering and dust growth) and show how dust growth helps explain the observed CO luminosity in the low-metallicity WLM dwarf galaxy. Our results have important implications for galaxies observed in the early universe by JWST.

Coronal mass ejections (CMEs) are among the most spectacular eruptive phenomena in the solar system and are the primary driver of hazardous space weather. Current data-driven magnetohydrodynamic (MHD) simulations at the active-region scale have demonstrated preliminary success in reproducing the initiation of solar eruptions, i.e., the origin of CMEs. However, these simulations cannot track the evolution and propagation of CMEs as they form and enter the solar wind. On the other hand, global simulations spanning the corona to interplanetary space typically assume non-equilibrium magnetic flux ropes as the initial CME structures, neglecting the pre-eruptive magnetic field and triggering processes. To achieve a more self-consistent description of the complete process of CME origin and propagation, we developed a high-precision, global data-driven model covering the solar surface to 1.5 AU, capable of tracking the entire chain of CME processes—from their genesis in the source region, through near-Sun evolution, to interplanetary propagation and arrival at Earth’s orbit. The simulation provides critical insights into the 3D structure, evolution mechanisms, propagation processes, and key factors determining the geoeffectiveness of CMEs.

The dynamical processes driving common envelope evolution—including core inspiral, envelope ejection, and the emergence of luminous red novae—remain poorly understood despite their critical role in binary stellar evolution. Using 1D radiation-hydrodynamic simulations, we investigate how radiation forces and recombination contribute to envelope unbinding. We perform a comprehensive parameter space survey of a 1D radiation hydrodynamic ejecta model to study their impact on the fraction of the unbound material and light curves. The parameters include the radiation energy to gas energy ratio $\beta$, equation of state (EoS), and speed of the ejecta. By comparing simulation results across different models, we can conclude that radiation force is important to the unbinding of the envelope. Our results demonstrate that radiation forces dominate in a layer below the recombination front, where high opacity and luminosity maximize their influence. A realistic EoS enhances the pressure gradient, further promoting envelope ejection. Both $\beta$ and EoS significantly alter light curve shapes, providing observable diagnostics for these processes. These findings advance our understanding of the mechanisms governing common envelope dynamics and their observational signatures.

The acceleration process remains a significant challenge in solar physics. Solar energetic particle (SEP) events are broadly categorized into impulsive events—rich in electrons and heavy ions like helium, carbon, oxygen, and iron, with enhanced isotopic ratios such as ³He to ⁴He—and gradual events, dominated by protons and higher-energy particles. These classifications stem from their distinct temporal profiles and acceleration mechanisms. Observational data highlight key characteristics of SEP groups:i) Multi-component Composition: The solar atmosphere primarily comprises hydrogen (71.5%) and helium (27%), with heavier elements making up just a small fraction of the total mass. Each element exhibits a unique energy spectrum, which helps pinpoint the underlying acceleration mechanism. This study focuses on hydrogen and helium, the primary components modeled. ii) Variable Elemental Abundances: The solar atmosphere consists of plasma with varying chemical elements or ion concentrations across its layers (chromosphere, transition region, or corona). Abundances are typically expressed as ratios relative to hydrogen, the Sun’s most abundant element. For example, helium’s abundance is roughly 10% of hydrogen’s by number. These abundances differ across solar regions (e.g., active regions vs. quiet Sun) or features like coronal loops and prominences, often due to processes like fractionation. Measurements are typically reported as number density (particles per unit volume), energy flux (e.g., particles·cm⁻²·s⁻¹·sr⁻¹·MeV⁻¹), or flux ratios at specific energies.iii) Isotopic Diversity: Isotopes in the solar atmosphere are variants of elements with the same number of protons but different neutron counts, leading to distinct atomic masses (e.g., carbon-12 vs. carbon-13). Isotopic abundance, the relative proportion of an isotope compared to others of the same element, is critical for understanding the Sun’s composition, nuclear processes, and spectroscopic signatures. For instance, hydrogen in the solar atmosphere is predominantly ¹H (protium), with trace amounts of ²H (deuterium). .iv) Varying Ionization Levels: The solar atmosphere’s plasma displays significant spatial variation in ionization, ranging from partially to fully ionized states, primarily driven by local temperature conditions. v) Diverse Charge States: The ionization state of a particle, indicated by the number of electrons lost (e.g., O⁶⁺ or Fe¹⁰⁺), provides insights into SEP acceleration in the heliosphere. The Fe/O ratio, for example, reflects the relative abundance of iron to oxygen. Enrichment of heavy ions or ³He in impulsive events points to selective acceleration mechanisms, such as resonant wave-particle interactions in flare regions. Analyzing these species’ abundances, charge states, and energy distributions sheds light on the dynamics of SEP acceleration. The multi-component, abundance, and isotopic properties are critical for understanding the acceleration mechanisms of solar energetic particles (SEPs). These are typically characterized using elemental ratios (e.g., Fe/O for iron to oxygen, He/H for helium to hydrogen, C/O for carbon to oxygen), isotopic abundance ratios (e.g., ³He/⁴He), and element abundances (e.g., electrons, protons, helium, carbon, and other metals within the simulation domain). This study focuses on impulsive SEP events driven by fully coupled hydrodynamics-magnetohydrodynamics-kinetic continuous scale (FC-HDK-CS) 3D flare-loop considerable temporal-spatial turbulence magnetic reconnection (LTSTMR). These events, characterized by current sheet thicknesses relative to characteristic electron lengths (e.g., Larmor radius for low-β or electron inertial length for high-β, with ratios on the order of 10¹⁰–10¹¹, and evolution times to electron cyclotron time ratios on the same order), are rich in electrons and heavy ions and associated with type III radio bursts. We developed the Multi-Component Abundance-Isotope (M-CAI) physical and mathematical model, implemented a relativistic algorithm, and enhanced the Relativistic Hybrid Particle-in-Cell Lattice Boltzmann Method (RHPIC-LBM) code. These advancements provide novel tools for analyzing Langmuir turbulence acceleration (LTA) driven by nonlinear resonant wave-particle interactions.

Population III (or Pop III) stars are thought to end their lives as energetic pair-instability supernovae (PISNe) or hypernovae, which emit copious X-ray radiation. These X-rays ionize the intergalactic medium and the gas in minihalos, enhancing the formation of molecular hydrogen (H_2) via the H- channel. The increased H_2 abundance promotes further Pop III star formation, which in turn leads to more PISNe/hypernovae, reinforcing the X-ray feedback loop. However, previous simulations of Pop III star and galaxy formation have not adequately accounted for this X-ray feedback loop. In this talk, I propose an “on-the-fly” calculation of the X-ray background as a means to estimate the X-ray feedback loop, and I present its impact on the expected number of Pop III stars.

Although the influence of evaporation, stellar evolution, binary stars, and smooth external potential due to the galaxy on star cluster evolution has been intensively investigated, less attention has been paid on another possibly important mechanism: interaction with giant molecular clouds. In this work, we utilise the results of detailed SPH simulations of the interstellar medium in a Milky Way-like galaxy as a background potential for open cluster evolution. Star cluster dynamics is calculated by the code nbody6 including stellar evolution. We find that the life-time of open clusters is reduced by a factor of 1.5 to 2 when molecular clouds are included at the surface density of the Galaxy, and by a factor of 3 to 5 at the surface density 10x higher than that of the Galaxy. The presence of molecular clouds causes elevated virial ratio (T/|U|) of the clusters, which can be higher than 1 for a substantial fraction of the life-time of the clusters. The stars are weakly unbound, but it still takes them tens of Myr to escape from the clusters. After the first stronger encounter, the clusters get scattered vertically out of the plane, which extends their life-time by diminishing further encounters. Surprisingly, the life-time of more massive clusters is reduced more by the presence of molecular clouds than that of lower mass clusters.

We investigate the formation and quenching of massive galaxies using high-resolution cosmological zoom-in hydrodynamical simulations performed with a modified GADGET-3 code. Our simulations incorporate both supernova (SN) and active galactic nucleus (AGN) feedback, allowing us to probe their relative roles in regulating star formation across cosmic time. We find that mechanical AGN feedback—implemented as high-velocity winds—efficiently quenches in-situ star formation in massive galaxies, while SN feedback strongly influences the growth of their ex-situ stellar components by regulating the star formation in low-mass progenitors. Comparison with observations, including stellar metallicity distributions in galaxy outskirts, satellite galaxy abundances, and central stellar core sizes, provides stringent tests for our feedback models. These results underscore the importance of realistic feedback prescriptions in galaxy formation simulations and highlight how simulations can bridge the gap between theory and observation.

The center of our Milky Way galaxy hosts a series of energetic outbursts, including the well-known Fermi and eROSITA bubbles, Galactic center chimneys, the inner 15-pc Sgr A lobes. Are they long-lasting or fast evolving explosive events? What causes these structures? The Fermi and eROSITA bubbles may correspond to typical galactic feedback processes occurring in our own Galaxy in the near past. Galactic feedback is one central unsolved problem in contemporary astronomy, and the Fermi and eROSITA bubbles are also galactic-scale accelerators of cosmic rays, whose origin remains a century-long mystery. In this talk, I will talk about our double-episode jet model of the eROSITA and Fermi bubbles. I will also present our new TDE jet model for the origin of the Sgr A lobes.

We investigate the formation of multiple spiral structures, which are typically transient and recurrent, in disk–halo systems using a GPU $N$-body code without provoking perturbations. We explore how adjusting numerical resolution, varying local stability, and adopting a live halo influence the initial development and subsequent growth of spiral arms. To study the emerging spiral modes, we compute Fourier components from $m=6$ to $m=1$ and measure their RMS sum values. In a marginally unstable disk of our low-resolution models (with $N_{\star}=5\times10^6$ and $N_{\rm{DM}}=1.1\times10^7$), faint spiral arms appear earlier within the first 0.5 Gyr due to numerical noise compared to high-resolution models with ten times more particles. Increasing the resolution delays the formation of non-axisymmetric features, though the time delay is significantly shorter in more unstable disks. We observe an `inverse cascading’ in the Fourier modes: initially, higher-$m$ spiral modes form and decay repeatedly, then progressively lower-$m$ modes emerge, featuring a delayed onset with a radially inward-drifting epicenter. This sequence reflects the transient nature of swing amplification, where short-wavelength modes grow first, and subsequently larger-scale features dominate as the disk evolves. Our results converge on this behavior for both low- and high-resolution runs when the star-to-dark-matter mass ratio is below 1:10. In contrast, it is absent in models with a fixed potential halo or a halo containing only $1.14\times10^6$ particles, regardless of the softening length, and those models do not exhibit a net-like pattern in the redistribution of stellar angular momentum around each epicenter. Our findings indicate that dynamical friction or impedance from the halo to the spiral arms plays a key role in forming, growing, and repeatedly regenerating multi-armed spiral structures in harmonic order, even though the total angular momentum exchange is only fractional.

Faint dwarf galaxies in the vicinity of the Milky Way have extremely low metallicities and are considered ideal systems for studying the nature of the first stars. In this study, we combine the semi-analytic model A-SLOTH (Ancient Stars and Local Observables by Tracing Halos) with high-resolution cosmological simulations to follow the formation and chemical evolution of Milky Way dwarf galaxies from the early universe to the present. We focus on the relationship between their star formation histories and metallicity distributions. In this talk, we present our results and discuss how the first stars contribute to the chemical evolution of dwarf galaxies. In the near future, the Galactic archaeology program using the Prime Focus Spectrograph on the Subaru Telescope is expected to significantly improve measurements of metallicities and radial velocities for a number of stars within dwarf galaxies. This study provides a theoretical basis for interpreting such upcoming data.

Recent observations of nearby spiral galaxies by JWST revealed that spiral arms are full of bubbles. These bubbles are supposed to be originated from supernova remnants. Recent theories suggest that they contribute to the formation of molecular clouds in the interstellar medium. However, the connection between the observed bubbles and supernova remnants has not yet been quantitatively discussed from a theoretical perspective. If this connection is confirmed, it would significantly clarify our understanding of interstellar medium and the paradigm of star formation. In this work, we investigate whether supernova remnants can quantitatively account for the observed bubbles. First, we simulate the evolution of supernova remnants and estimate the size of bubbles created by the primary shock front. During these simulations, we encountered numerical challenges due to the artificial broadening of shock fronts and to resolve the cooling length. To overcome this problem, we developed a technique called the Retarded Cooling Method, which delays the cooling process at shock surfaces where radiative cooling should not occur physically. Next, we statistically construct a histogram of bubble sizes in the galactic disk. The histogram is then compared with observational data to examine the physical conditions of spiral galaxies.

The merger of two neutron stars can form a system composed of a central object (either a neutron star or black hole) and a centrifugally supported disk. Inside the disk, magnetorotational instability generates a turbulent state, which then induces an effective viscosity. The viscous angular momentum transport and heating can evolve the system and trigger mass ejection from the disk on a timescale of seconds. The post-merger mass ejection contributes to the total ejecta in addition to the violent merger phase and to shaping the abundance pattern of heavy nuclei produced via the r-process. In this talk, I will present the results of our numerical simulations of such systems and their implications.

We performed a second-long numerical relativity radiation MHD simulation for a binary neutron star merger in which a black hole surrounded by an accretion torus is promptly formed after the merger. We find that a magneto-sphere dominated by an aligned global magnetic field penetrating the black hole developes thanks to the alpha-Omega dynamo activated in the accretion torus. As a result, a collimated Poynting-flux dominated jet is launched which may drive a low luminosity short gamma-ray bursts.

We report the first numerical demonstration of ringdown gravitational waves from non-merging equal mass black hole encounters. Through numerical relativity simulations of close fly-by events, we identified ringdown signals produced by dynamical tidal deformations during hyperbolic scattering. The extracted ringdown frequencies agree remarkably well with quasi-normal modes of individual black holes, revealing that strong gravitational interactions can generate detectable ringdown radiation without merger. This finding opens new possibilities for gravitational wave astronomy in dense stellar environments where such encounters may occur.

Neutron stars exhibit sudden changes of its rotational velocity, known as “pulsar glitches”. It has been believed that glitches are mainly caused by superfluid neutron vortices in the inner crust of neutron stars. However, importance of contributions of the outer core has been recently discussed, and further microscopic investigations of quantum vortices and fluxtubes in the outer core of neutron stars are highly desired. In this study, we investigate the interaction between quantum vortices of $^3P_2$ superfluid neutrons and fluxtubes of $^1S_0$ superconducting protons in the outer core of neutron stars, based on a successful bosonic theory of superfluid, the Gross-Pitaevskii equation (GPE). In this talk, we will discuss how the interaction of the $^3P_2$ superfluid vortices and the proton fluxtubes under a magnetic field in the outer core of neutron star affect the structure of the vortices in the neutron star.

We numerically calculate the axisymmetric, steady-state, force-free magnetic field configuration around a rotating black hole. Following many previous studies, we solve the general relativistic Grad–Shafranov equation using the relaxation method. As boundary and initial conditions, we adopt a jet-like configuration based on a combined magnetic field, which is constructed from a superposition of a uniform vertical field (the Wald solution) and a Blandford–Znajek split-monopole field. The resulting force-free solution resembles a parabolic configuration near the black hole and approaches a uniform vertical field at large distances. Due to these characteristics, the solution possesses only a single singularity (the inner light surface). This allows us to determine the relation among the phi component of the vector potential, the angular velocity of magnetic field lines, and the toroidal magnetic field component from the regularity condition at infinity. We perform numerical calculations under various conditions.

To shed light on the nature of dark matter and dark energy, the next-generation cosmological surveys, i.e., Stage-IV surveys, will cover the widest areas with deep imaging. The high-fidelity observational data will enable the measurement of cosmological statistics at unprecedented precision. In the statistical analysis, simulations are critical to estimate possible systematic uncertainties and cosmic variances. To this end, the simulations are required to cover the large volume of Stage-IV surveys and a large number of independent realizations is also important to robustly estimate the covariance matrix of cosmological statistics. Furthermore, the cross-correlations of different observables, e.g., galaxy-galaxy lensing, are employed to extract more information from the Stage-IV survey data. For the cross-correlations, multiple observables must be modeled in simulations, in a correlated manner. In this talk, I will introduce the recent campaign to construct a suite of multi-wavelength cosmological simulations conducted with the supercomputer Fugaku, Japanese flagship supercomputer. We first produce the large suite of dark-matter only simulations and post-process them to generate various mock observations: weak lensing, cosmic microwave background, and galaxy distributions. This simulation suite will be ideal for measuring covariance matrices of statistics such as auto- and cross-correlations.

The upcoming generation of large-scale surveys—such as CSST, DESI, EUCLID, and SKA—will map billions of galaxies across cosmic time, offering unprecedented opportunities to probe dark matter, dark energy, and galaxy formation. To fully exploit these surveys, we need realistic cosmological mock catalogs that bridge the gap between theory and observation. In this talk, I will present our efforts to construct such mocks, combining large-volume N-body simulations with semi-analytical galaxy formation models to capture both the large-scale structure and small-scale physics. I will particularly discuss the development of emission-line galaxy models and physically motivated AGN SEDs, which are crucial for connecting simulations with spectroscopic and photometric observations. These mock datasets provide a unified framework for interpreting galaxy populations from high redshift to the present day, testing cosmological models, and preparing for the challenges of next-generation surveys.

Observations favor cosmological models with a time-varying dark energy component. But how does dynamical dark energy (DDE) influence the growth of structure in an expanding Universe? We investigate this question using high-resolution N-body simulations based on a DDE cosmology constrained by first-year DESI data (DESIY1+DDE), characterized by a 4% lower Hubble constant (H0) and 10% higher matter density (Omega0) than the Planck-2018 LCDM model. We examine the impact on the matter power spectrum, halo abundances, clustering, and Baryonic Acoustic Oscillations (BAO). We find that DESIY1+DDE exhibits a 10% excess in power at small scales and a 15% suppression at large scales, driven primarily by its higher Omega0. This trend is reflected in the halo mass function: DESIY1+DDE predicts up to 70% more massive halos at z=2 and a 40% excess at z=0.3. Clustering analysis reveals a 3.71% shift of the BAO peak towards smaller scales in DESIY1+DDE, consistent with its reduced sound horizon compared to Planck18. Measurements of the BAO dilation parameter, using halo samples with DESI-like tracer number densities across 0<z<1.5, agree with the expected DESIY1+DDE-to-Planck18 sound horizon ratio. After accounting for cosmology-dependent distances, the simulation-based observational dilation parameter closely matches DESI Y1 data. We find that the impact of DDE is severely limited by current observational constraints, which strongly favor cosmological models – whether including DDE or not – with a tightly constrained parameter Omega0h^2 ~ 0.143, within 1-2% uncertainty. Indeed, our results demonstrate that variations in cosmological parameters, particularly Omega0, have a greater influence on structure formation than the DDE component alone.

Modern cosmological simulations demand unprecedented computational resources to accurately model complex astrophysical phenomena across multiple scales. In this talk, I introduce RAMSES-DARWIN, a state-of-the-art numerical framework specifically designed for upcoming heterogeneous computing infrastructures in Korea, solving the computational bottlenecks encountered in large-scale cosmological hydrodynamic simulations. Our approach combines innovative algorithmic strategies with cutting-edge parallel computing architectures to address the substantial computational demands of physical modules involved in star formation and its feedback processes, galaxy and cluster formation, and evolution of large-scale structures. RAMSES-DARWIN emphasizes hybrid computing paradigms that synergistically exploit both CPU and accelerator (e.g., GPU) technologies. By implementing adaptive load balancing and memory-efficient data structures, the framework delivers significant performance improvements while preserving numerical accuracy. Additionally, it features dynamic resource allocation algorithms that optimize workload distribution based on hardware capabilities, alongside specialized numerical solvers tailored for heterogeneous computing environments. The presentation will demonstrate performance benchmarks across various hardware systems, showing how these developments enable new research possibilities in computational astrophysics. I will also discuss implementation strategies for next-generation exascale computing systems and explore their implications for advancing our understanding of cosmic evolution through high-fidelity hydrodynamic simulations.

Standard cosmology suggests that the first galaxies formed between redshifts z ≈ 10–20. Recent advances with the James Webb Space Telescope (JWST) have enabled the observation of galaxies in the early universe at z > 10. Interestingly, these galaxies often exhibit very compact and clumpy internal structures, highlighting the need for revised or extended models of early galaxy formation. In this study, we perform high-resolution cosmological simulations targeting a first galaxy, resolving gas temperatures down to T ≈ 2000 K and gas number densities up to n_H ≈ 10^5 cm^-3, to investigate gas dynamics and star formation during early galaxy assembly. We find that a compact gas clump (Σ_gas ~ 10^3 Msun/pc^2) forms at the halo center through Lyman-Werner radiation and atomic cooling, followed by a starburst triggered by metal-line cooling. The resulting star cluster has a stellar mass of ~10^6 Msun, a high stellar surface density (Σ_* ≈ 10^3–10^4 Msun/pc^2), and an average metallicity of -2.5 < log (Z/Zsun) < -2—similar to that of observed globular clusters. Moreover, we find that more massive first galaxies form through the hierarchical mergers of such compact star clusters.

Tidal disruption events (TDEs), in which a star is disrupted by a supermassive black hole (SMBH), have recently been observed to produce radio emission. Most radio-emitting TDEs lack on-axis relativistic jets (non-jetted TDEs), which provide a new observational window into gas dynamics and mass outflows near SMBHs. This talk overviews the observational properties of non-jetted TDEs and introduces a one-dimensional thin-shell model based on Hayasaki & Yamazaki (2023, 2025). The model incorporates time-dependent mass injection and SMBH gravity, showing transitions of the radio-emitting shell among gravity-, energy-, and momentum-dominated phases. In the gravity-dominated phase, a universal relation connecting the shell radius and velocity enables SMBH mass estimation from radio observations. This talk also discusses how this theoretical framework can be tested using numerical simulations and observational data, aiming to bridge theory, simulation, and observation.

To estimate the black hole spin parameter and the inclination angle from images, we develop a neural network based on E(2)-Equivariant Steerable Convolutional Neural Networks (CNN) architecture, which is a kind of CNN. We generate approximately $10^4$ images using a semi-analytic model that depends only on the spin parameter and the inclination angle, and train the network. As a result, for the validation data, we successfully developed a neural network capable of recovering the spin parameter and the inclination angle with higher accuracy than in previous works. Specifically, the standard deviations of the maximum absolute error between the true and predicted values of the spin parameter and the inclination angle are approximately $0.005$ and $0.62^\circ$. Additionally, by comparing the quantitative performance of our network with the conventional CNN used in previous works, we find that our network achieves superior performance. We further assume the Black Hole Explorer observation and investigate the performance for the blurred images. As a result, it is found that the spin parameter and inclination angle of rapidly spinning black holes (dimensionless spin parameter is approximately greater than $0.5$) can be recovered with especially higher accuracy. In this case, the standard deviations of the maximum absolute error between the true and predicted values of the spin parameter and the inclination angle are approximately $0.04$ and $3.0^\circ$. Even for the case of slowly rotating black holes, we confirm that the prediction accuracy exceeds that of previous studies. Finally, we also investigate the performance for the images of accretion disks with the finite scale height. It is found that the absolute error between the true and estimated values is approximately smaller than $0.1$ for the spin parameter, and $10^\circ$ for the inclination angle. In addition, when these images are blurred with the BHEX resolution, the performance of our network remains high. Our result demonstrates the feasibility of independently and accurately estimating the spin parameter and inclination angle from images obtained by high-resolution observations, such as those expected from Black Hole Explorer.

Understanding gas outflows in active galactic nuclei (AGNs) is key to revealing how supermassive black holes interact with their surroundings. We investigate the dynamical structures of both dusty and dust-free gas on sub-parsec scales around a black hole using the two-dimensional radiation hydrodynamics simulations. Our simulations show that radiation pressure drives time-dependent, multi-shell outflows. These outflows exhibit non-spherical structures shaped by shocks, inflowing gas, and variations in the dust sublimation radius, which changes with both angle and time. We find that an Eddington ratio above 10^{-3} significantly increases the outflow velocity and expands the sublimation region, while the overall temperature and density distributions remain largely unaffected. By combining numerical and analytical approaches, we provide a clearer picture of how radiation-driven processes determine the terminal velocity in AGNs.

In this presentation, we analyze the acceleration of a stationary high-temperature gas into relativistic velocity using numerical simulations. The Special Relativistic Godunov Smoothed Particle Hydrodynamics (SRGSPH) method [Kitajima et al. 2025, revised] is employed to model the fluid as a collection of discrete particles (SPH particles), allowing the fluid motion to be tracked through the interactions between particles and their environment.

By leveraging the advantages of the SPH method—such as its ability to handle vacuum boundaries and large deformations—we investigate the dynamics of fluids expanding into vacuum, a challenge for conventional grid-based methods. This study demonstrates the utility of SRGSPH in simulating jet formation from a high-temperature source and elucidates its acceleration mechanism into a vacuum.

Our findings offer valuable insights into the driving processes of relativistic jets, with implications for phenomena such as those observed in active galactic nuclei and gamma-ray bursts.

We investigate the time evolution of sub-Keplerian transonic accretion flow onto a non-rotating black hole using axisymmetric viscous hydrodynamics simulations. We simulate the shocked accretion flow by using boundary values from semi-analytical analysis. First, we stabilise the accretion disk in the presence of cooling. Cooling drives the shock inward by reducing the temperature gradient term and forms an almost stable shocked disk. Then, by choosing suitable viscosity parameters, we explored the dynamics of the disk configuration. We find that turbulence can develop in the post-shock region, with the intensity of turbulence increasing as the viscosity parameter increases. Viscosity tends to push the shock surface outward by redistributing the angular momentum, ultimately influencing the stability and dynamics of the shock. We compute the disk luminosity for various configurations. The oscillations of the post-shock disk lead to quasi-periodic oscillations in the synthetic light curve. These results are similar to some complex observational features of black hole accretion physics.

We investigate the formation and morphology of galaxies and their morphological evolution using the Horizon Run 5 cosmological simulation. We measure asymmetry and morphology of stellar mass component of the galaxies with stellar mass M> 2×109 M⊙ to classify them into disk, spheroid, and irregular types. We find that the initial morphology of the galaxies in the cosmic morning is dominantly disk type with the Sersic index < 1.5. The fraction of disk-type galaxies is ~2/3 and that of irregular or spheroid type is ~1/6. Irregular or spheroidal morphology is incidental and transient. The fractions are roughly independent of redshift and also of stellar mass up to 1010M⊙. The fraction of spheridal and irregular types increases for more massive galaxies, but at a fixed mass the fraction of disks monotonically increases . We use the same morphological parameters to classify the high-redshift galaxies in the JWST observations and confirm the dominance of disk galaxies at redshifts > 1 and upto 6 with a very good quantitative agreement in morphology fractions and their redshift evolution. Some other results from the HR5 simulation are reviewed too.

In the Λ-CDM model, gas and dark matter (DM) mix together within DM haloes. Gas cools down and condenses into the centre of haloes, forming new galaxies. It is natural to assume that gas and DM share identical specific AM (sAM, e.g. Mo et al., 1998). However, modern cosmological simulations have challenged this assumption. Even in the non-radiative hydrodynamical simulations, removing the effects from radiative cooling and feedback, gas on average has 30%-40% more sAM than DM (Chen et al., 2003; Sharma & Steinmetz, 2005; Zjupa & Springel, 2017). This work reproduces this result by analysing ∼50,000 well-resolved DM haloes in a non-radiative simulation. This can be pinned down to the excess AM of gas in the inner halo, which hints that DM may transfer AM to gas. We uncover the leading driver for this AM difference through a series of control simulations of a collapsing ellipsoidal top-hat. These runs reveal that the pressurized inner gas shells collapse more slowly, causing the DM ellipsoid to spin ahead of the gas ellipsoid. The rising torque generally transfers AM from the DM to the gas. The amount of AM transferred via this mode depends on the initial spin, the initial axes ratios, and the collapse factor. These quantities can be combined in a single dimensionless parameter, which robustly predicts the AM transfer of the ellipsoidal collapse. This model can robustly explain the average AM excess in controlled and cosmological simulations. Additionally, we find this AM transfer mode happens in major mergers of DM haloes, which leads to the gas-to-DM spin ratio change considerably during major mergers. The novelty of this study is that it identifies a mechanism through which AM is exchanged between gas and DM during halo formation and major mergers, and further research is required to extend our understanding to cosmological hydrodynamical simulations with full galaxy formation physics.

We investigate the role of intermediate-mass black holes (IMBHs) in shaping the morphology of dwarf galaxies in the early universe. Using a suite of high-resolution zoom-in cosmological simulations at redshift ~2, we analyze how different black hole seed masses, formation times, and feedback models influence the evolution of central dwarf galaxies residing in halos of similar mass. Our findings demonstrate that AGN feedback, particularly the strength of black hole-driven winds, plays a critical role in regulating gas content, star formation, and galactic structure. Galaxies experiencing stronger feedback tend to develop flatter morphologies, lower stellar masses, and distinct central concentrations, often characterized by low Sersic indices and intermediate rotational support. We also discuss how commonly used low-redshift morphological diagnostics may misrepresent high-redshift systems, especially when viewed through synthetic JWST imaging, where resolution effects can inflate perceived galaxy sizes. These results offer new insights into the connection between IMBH growth and the emergence of galaxy structure in the early universe.

In cosmological simulations with the AMR code RAMSES, we observe unphysical “jumps” in global quantities like the cosmic star formation rate density (cSFRD). We found these jumps are caused by the “holdback” refinement method, which is used to maintain a constant physical resolution. This method can trigger massive, simultaneous refinement events, which artificially boosts local gas densities and triggers spurious star formation.

To fix this problem, we are testing a solution with three main parts. First, we modified the refinement strategy itself. Our new hybrid method uses a standard quasi-Lagrangian refinement for large cells, but applies a stricter, sub-Lagrangian rule once a cell’s size (Δx) drops below our target resolution (Δr). This makes it harder for already small cells to refine further and prevents the large, simultaneous refinement events.

Second, to be consistent with our new refinement scheme, we apply a pressure floor based on the Truelove criterion. This prevents artificial fragmentation by making sure the Jeans length is always properly resolved. Finally, we are working on smoothing the density source for the Poisson solver. This smoothing will also be tied to our target resolution to keep the gravity calculation stable.

Our initial tests with the first two modifications already show that the artificial jumps in the cSFRD are significantly reduced. At the meeting, we plan to share our latest results from this approach. We hope to receive technical feedback, especially on our ongoing work with the gravity solver.

We investigate the depletion radius of dark matter haloes by analyzing stacked mass flow rate (MFR) profiles derived from a large N-body cosmological simulation, focusing on its dependence on various halo properties. We find that the MFR profiles are nearly self-similar across a wide range of halo masses. The depletion radius is primarily determined by the conventional mass accretion rate, defined as the logarithmic growth rate of the virial mass with respect to the scale factor, consistent with similar previous findings. By analyzing the phase-space structure of haloes, we find that this dependence originates from the underestimated halo boundary definition and the crossing of non-smooth accretion material through this boundary. In contrast, the depletion radius associated with smooth accretion component shows little dependence on the conventional accretion rate.