On the Thermal Stability of Radiation Dominated Accretion Disks

We study the long-term thermal stability of radiation dominated disks in which the vertical structure is determined self-consistently by the balance of heating due to dissipation of MHD turbulence driven by the magneto-rotational instability (MRI), and cooling due to radiation emitted at the photosphere. The calculations adopt the local shearing box approximation, and utilize the recently developed radiation transfer module in the Athena MHD code based on a variable Eddington tensor rather than an assumed local closure. After saturation of the MRI, in many cases the disk maintains a steady vertical structure for many thermal times. However, in every case in which the box size in the horizontal directions is at least one pressure scale height, fluctuations associated with MRI turbulence and dynamo action in the disk eventually trigger a thermal runaway which causes the disk to either expand or contract until the calculation must be terminated. During runaway, the dependence of the heating and cooling rates on total pressure satisfy the simplest criterion for classical thermal instability. We identify several physical reasons why the thermal runaway observed in our simulations differ from the standard alpha disk model, for example the advection of radiation contributes a non-negligible fraction to the vertical energy flux at the largest radiation pressure, most of the dissipation does not happen in the disk mid-plane, and the change of dissipation scale height with mid-plane pressure is slower than the change of density scale height. We discuss how and why our results differ from those published previously. Such thermal runaway behavior might have important implications for interpreting temporal variability in observed systems, but fully global simulations are required to study the saturated state before detailed predictions can be made.