Radiative MHD simulations of photon bubbles in radiation-supported magnetized atmospheres of neutron stars with isotropic Thomson scattering

A major uncertainty in the structure and dynamics of magnetized, radiation pressure dominated neutron star accretion columns in X-ray pulsars and pulsating ultraluminous X-ray sources is that they are thought to be subject to the photon bubble instability. We present the results of two dimensional radiation relativistic magnetohydrodynamic simulations of a non-accreting, static atmosphere to study the development of this instability assuming isotropic Thomson scattering in the slow diffusion regime that is relevant to neutron star accretion columns. Photon bubbles generally grow faster toward shorter wavelengths, until a maximum growth rate is achieved at the radiation viscosity length scale, which is generally quite small and requires high numerical resolution to simulate. We confirm the consistency between our simulation results and linear theory in detail, and show that the nonlinear evolution inevitably leads to collapse of the atmosphere with the higher resolution simulation collapsing faster due to the presence of shorter length scale nonlinear structures. At least in static atmospheres with horizontally periodic boundary conditions, this resolution dependence may make simulations of the nonlinear dynamics of photon bubble instability in neutron star accretion columns challenging. It remains to be seen whether these difficulties will persist upon inclusion of an accretion flow through the top and magnetically-confined horizontal boundaries through which photons can escape. Our results here provide a foundation for such future work.