I have implemented AREPO-RT, a novel radiation hydrodynamic (RHD) solver for the unstructured moving-mesh code AREPO . It solves the moment-based radiative transfer equations using the M1 closure relation. I achieve second-order accuracy by using a slope limited linear spatial extrapolation and a first order time prediction step to obtain the values of the primitive variables on both sides of the cell interface. A Harten-Lax-Van Leer flux function, suitably modified for moving meshes, is then used to solve the Riemann problem at the interface. The implementation is fully conservative and compatible with the individual timestepping scheme of AREPO. It incorporates atomic Hydrogen (H) and Helium (He) thermochemistry, which is used to couple the ultra-violet (UV) radiation field to the gas. Additionally, the infrared radia- tion is coupled to the gas under the assumption of local thermodynamic equilibrium between the gas and the dust. We successfully apply our code to a large number of test problems, in- cluding applications such as the expansion of H II regions, radiation pressure driven outflows and the levitation of optically thick layers of gas by trapped IR radiation. The new implemen- tation is suitable for studying various important astrophysical phenomena, such as the effect of radiative feedback in driving galactic scale outflows, radiation driven dusty winds in high redshift quasars, or simulating the reionisation history of the Universe in a self consistent manner.
Interstellar dust is an important component of the galactic ecosystem, playing a key role in multiple galaxy formation processes. We have been developing a novel numerical framework for the dynamics and size evolution of dust grains implemented in the moving-mesh hydrodynamics code AREPO suited for cosmological galaxy formation simulations. We employ a particle-based method for dust subject to dynamical forces including drag and gravity. The drag force is implemented using a second-order semi-implicit integrator and validated using several dust-hydrodynamical test problems. Each dust particle has a grain size distribution, describing the local abundance of grains of different sizes. The grain size distribution is discretised with a second-order piecewise linear method and evolves in time according to various dust physical processes, including accretion, sputtering, shattering, and coagulation.
The numerical modeling of the anisotropic diffusion equation is problematic and non-trivial. Widely used discretization approaches violate the second law of thermodynamics i.e., heat can flow from lower to higher temperatures. This accentuates temperature extrema causing numerical instabilities which can trigger unphysical temperature oscillations.
I implemented an extremum preserving anisotropic diffusion solver for the unstructured meshes of the moving mesh code AREPO (Kannan et al. 2016b). It relies on splitting the one sided facet fluxes into normal and oblique components. This is achieved by decomposing the gradient of temperature in the co-ordinate system defined by the cell center and its appropriate neighbors. The neighbors (2 in 2D and 3 in 3D) that form the new co-ordinate system are chosen such that the components of the vector in this system are all positive. The flux along the face normal will always be along the temperature gradient but the same cannot be said about the other oblique components. Therefore, the oblique flux is then non-linearly limited in such a way that the total flux is both locally conservative and also extremum preserving. The extremum preserving property of the scheme ensures that the second law of thermodynamics is not violated.
The values of the relevant variables at the cell faces are extrapolated from the mesh generating points using a a very simple yet robust interpolation scheme that works well for strong heterogeneous and highly anisotropic problems. The required discretization stencil is essentially small, consisting of just the point and its Delaunay connections. The numerical diffusivity was shown to be as low as ~1%, even at really low resolution and decreases at almost second order with increase in resolution. This ensures that our scheme has negligible numerical diffusivity in all practical applications.
We have also implemented a semi-implicit algorithm where, the non-linear flux terms that depend on the internal energy are integrated explicitly, while the other terms are integrated implicitly. This linear system is solved using HYPRE, which is a library for solving large, sparse linear systems of equations on massively parallel computers. This method is extremely fast as it requires only the solution of a linear system per timestep and is faster than the simple explicit scheme by a factor of 20.
As an example we show the diffusion on a hot patch of gas surrounded by a cold medium with circular magnetic fields
Another example if the Sedov-Taylor blast wave problem. This problem is a good test to determine the accuracy of coupling between two distinct physical processes: hydrodynamics and diffusion. We simulate three different configurations of the blast wave problem, the classical adiabatic blast wave test (left panel), the blast wave test with isotropic conduction (middle panel), and the blast wave test with anisotropic conduction (right panel) in 3D. The magnetic field in the anisotropic blast wave test points in the x-direction. The size of the box is (100 pc)^3.
The dynamics of a rapidly conducting stratified plasma differs from that of an adiabatic fluid. The temperature gradient and the local orientation of the magnetic field determine the plasma’s convective stability. When the temperature increases with height, the convective instability is known as the heat-flux-driven buoyancy instability (HBI). This instability saturates by reorienting the magnetic fields perpendicular to gravity. In the following visualization gravity points downwards. The size of the box is (0.1 cm)^2.
When the temperature decreases with height, the convective instability is known as the magneto thermal instability (MTI). This instability does not saturate and can drive sustained convection. In the following visualization gravity points downwards. The size of the box is (0.1 cm)^2. The HBI and MTI make the ICM unstable even if dS/dr > 0, which has important implications for the coupling between the injected feedback energy from the central blackhole and the ICM.