TY - JOUR
AU1 - Potter, William, J.
AU2 - Balbus, Steven, A.
AB - Abstract Understanding what determines the strength of MHD turbulence in accretion discs is a question of fundamental theoretical and observational importance. In this work, we investigate whether the dependence of the turbulent accretion disc stress (α) on the magnetic Prandtl number (Pm) is sufficiently sensitive to induce thermal–viscous instability using 3D MHD simulations. We first investigate whether the α–Pm dependence, found by many previous authors, has a physical or numerical origin by conducting a suite of local shearing-box simulations. We find that a definite α–Pm dependence persists when simultaneously increasing numerical resolution and decreasing the absolute values of both the (microscopic) viscous and resistive dissipation coefficients. This points to a physical origin of the α–Pm dependence. Using a further set of simulations which include realistic turbulent heating and radiative cooling, by giving Pm a realistic physical dependence on the plasma temperature and density, we demonstrate that the α–Pm dependence is sufficiently strong to lead to a local instability. We show that the instability manifests itself as an unstable limit cycle by mapping the local thermal-equilibrium curve of the disc. This is the first self-consistent MHD simulation demonstrating the Pm instability from first principles. This result is important because a physical Pm instability could lead to the global propagation of heating and cooling fronts and a transition between disc states on time-scales compatible with the observed hard/soft state transitions in black hole binaries. accretion, accretion discs, black hole physics, instabilities, magnetic reconnection, MHD, turbulence 1 INTRODUCTION Accretion discs surrounding compact objects exhibit complex and dramatic cyclic changes in X-ray luminosity on a variety of time-scales (e.g. Fender, Belloni & Gallo 2004; Remillard & McClintock 2006; Done, Gierliński & Kubota 2007). Black hole binaries (hereafter ‘BHBs’) are observed to undergo flaring and quiescent cycles on month-long time-scales in which the X-ray luminosity of the disc changes by orders of magnitude. This long time-scale cyclic flaring is thought to be caused by a thermal–viscous instability due to the sharp increase in disc opacity when hydrogen becomes ionized. This is known as the disc instability mechanism (DIM) (Faulkner, Lin & Papaloizou 1983; Hameury et al. 1998; Lasota 2001). In the last decade, cyclic changes in the disc luminosity and X-ray spectrum have been observed in which the disc cycles between a radiatively efficient high soft state, resembling a classic composite blackbody and a radiatively inefficient low hard state, with a power law non-thermal hard X-ray spectrum. These changes in disc state are of particular interest because they are intimately linked to the production of a relativistic jet in the system; in the hard state a radio jet is observed, whilst in the soft state no jet is observed. However, the physical mechanism responsible for these state changes is not yet understood. In this paper we present simulations of a promising new thermal–viscous instability which may be involved with the hard/soft changes in disc state. In a previous paper (Potter & Balbus 2014) (hereafter PB), we put forth arguments for the existence of a new type of thermal–viscous disc instability triggered when the disc stress (usually parametrized by α, Shakura & Sunyaev 1973) depends sensitively on the physical properties of the disc plasma, especially temperature (see also the instability calculations in Takahashi & Masada 2011). Simulations both of accretion disc turbulence (Fromang et al. 2007; Lesur & Longaretti 2007; Simon & Hawley 2009; Meheut et al. 2015), and of driven turbulent dynamos (Schekochihin et al. 2004) had shown that when it is near unity, the magnetic Prandtl number (the ratio of the microscopic viscosity to resistivity, Pm = ν/η), can have a significant effect on the strength and the maintenance of turbulent fluctuations. A possible explanation for the dependence of the disc α parameter on Pm is that the rate of small scale magnetic reconnection in turbulent flow is sensitive to the strength of the viscous stress in a reconnecting layer when Pm is near unity (Balbus & Henri 2008). Magnetic and kinetic energy are extracted from the differential rotation of the disc by the magnetorotational instability (MRI), predominantly on spatial scales set by the disc scaleheight Balbus & Hawley (1998). Once the MRI is fully non-linear, the gas hosts a turbulent cascade, transferring energy downwards to ever smaller length scales, until finally it is dissipated as heat. In the Pm > 1 regime, the viscosity is larger than the resistivity on a given scale, and so the viscous dissipation is correspondingly larger than resistive dissipation; the opposite is true when Pm < 1. Both the growth rate of the MRI and the magnitude of the turbulent disc stress (most of which is magnetic) are related to the RMS magnetic field strength of the disc. The saturated magnetic field strength in turn depends upon the rate of magnetic reconnection in the plasma. When Pm > 1, the viscous dissipation scale exceeds the resistive scale and so the viscosity damps out velocity fluctuations on the resistive scale. This leads to a lower rate of reconnection via dissipation, a build-up of the magnetic energy, and an increased value of α. This is because the field dissipation process will generally require large velocity gradients over the resistive lengthscale in order to bring misaligned magnetic field lines together. If these large gradients produce correspondingly large viscous stresses, the dissipation will be inhibited. If on the other hand Pm < 1, the resistive length scale is larger than the viscous length scale, and velocity gradients on this resistive scale will not produce important dynamical stress. Magnetic dissipation unencumbered by viscosity will proceed apace, the cascade decreasing the saturated magnetic field strength, and with it the α stress. The disc radius at which the transition between Pm > 1 and Pm < 1 occurs is r ∼ 100rs in a typical BHB accretion disc (PB). The dependence of α on Pm is typically modelled as a power law, α ∝ Pmn. PB showed that within the classic α formalism, when 6/13 < n > 10/3, the disc is unstable. Many (though not all) MHD simulations of the dependence of α on Pm had found n ∼ 0.5 − 1, suggesting that the instability might plausibly be present in astrophysical discs (Fromang et al. 2007; Lesur & Longaretti 2007; Simon & Hawley 2009). Using an idealized 1D dynamic thin disc approximation, PB demonstrated that the instability does indeed lead to the formation of a local unstable limit cycle. Dynamical heating and cooling fronts move throughout the disc, as in classical dwarf nova modelling (e.g. Hameury et al. 1998). These preliminary results suggest two further avenues of exploration, which we follow here. The first is to investigate with some care the dependence of the disc stress upon the magnetic Prandtl number in a variety of initial conditions using 3D MHD isothermal shearing-box simulations. It is obviously critical to establish that the Pm dependence is physical, and not a numerical artefact. The second challenge is to conduct simulations of the disc thermal equilibrium curve using more realistic cooling, including a temperature and density dependent Pm. Does the suspected instability actually occur in a first-principle MHD simulation? The paper is organized as follows. In Section 2, the numerical setup of the MHD simulations is explained. This includes details of the cooling function and the temperature and density dependence of Pm. Section 3 contains the results of our isothermal local shearing-box simulations. These simulations study the effect on α, when the resistivity, viscosity, numerical resolution and initial net magnetic field configuration are changed. The purpose of these simulations is to test whether the α–Pm dependence is a numerical or physical effect and if the α–Pm dependence is sufficiently strong to trigger the Pm instability. In Section 4, the local thermal-equilibrium curve of the disc is calculated using a set of MHD simulations which include realistic turbulent heating and radiative cooling, and in which Pm has a realistic physical temperature and density dependence. This allows us to determine whether the Pm instability occurs in a self-consistent MHD simulation. In Section 5, the main results and conclusions of the paper are summarized. 2 NUMERICAL SETUP We carry out this investigation with the PLUTO MHD code (Mignone et al. 2007), which is widely-used and publicly available. The dissipative MHD equations are solved in the local shearing-box approximation using a Godunov-type finite volume scheme with explicit dissipation terms (Bodo et al. 2008; Hawley, Gammie & Balbus 1995; Mignone et al. 2012; Bodo et al. 2014). In conservation form, the equations are: \begin{equation} {\frac{\mathrm{\partial} \rho }{\mathrm{\partial} t}+\nabla \cdot (\rho {\boldsymbol v})=0,} \end{equation} (1) \begin{eqnarray} &&\!\!\!\!&&{\frac{\mathrm{\partial} (\rho {\boldsymbol v})}{\mathrm{\partial} t}+\nabla \cdot (\rho {\boldsymbol {vv}}-{\boldsymbol {BB}})+\nabla P_{\mathrm{tot}}} \nonumber\\ &&\quad=\rho {\boldsymbol g}_{s} - 2\Omega _{0}{\boldsymbol {\hat{z}}}\times \rho {\boldsymbol v}+\nabla \cdot \Pi , \end{eqnarray} (2) \begin{eqnarray} \frac{\mathrm{\partial} E}{\mathrm{\partial} t}&+&\nabla \cdot [(E+P_{\mathrm{tot}}){\boldsymbol v}-({\boldsymbol {v\cdot B}}){\boldsymbol B}]\nonumber \\ \quad &&=\, \rho {\boldsymbol v}\cdot {\boldsymbol g}_{s}-\nabla \cdot [(\eta \cdot {\boldsymbol J})\times {\boldsymbol B}]+\nabla \cdot ({\boldsymbol v}\cdot \Pi )-\Lambda , \end{eqnarray} (3) \begin{equation} \frac{\mathrm{\partial} {\boldsymbol B}}{\mathrm{\partial} t}-\nabla \times ({\boldsymbol v}\times {\boldsymbol B})=-\nabla \times (\eta \cdot {\boldsymbol J}), \end{equation} (4) \begin{equation} {\boldsymbol J}=\nabla \times {\boldsymbol B}, \qquad P_{\mathrm{tot}}=P+\frac{B^{2}}{2}, \end{equation} (5) \begin{equation} E=\rho e+\frac{\rho v^{2}}{2}+\frac{B^{2}}{2}, \qquad {\boldsymbol g}_{s}=\Omega _{0}^{2}(2qx{\boldsymbol {\hat{x}}}-z{\boldsymbol {\hat{z}}}), \end{equation} (6) \begin{equation} \Pi _{ij}=\nu \left[\frac{\mathrm{\partial} v_{i}}{\mathrm{\partial} x_{j}}+\frac{\mathrm{\partial} v_{j}}{\mathrm{\partial} x_{i}}-\frac{2}{3}\nabla \cdot {\boldsymbol v}\delta _{ij}\right], \end{equation} (7) where ρ is the mass density, |${\boldsymbol v}$| the fluid velocity, |${\boldsymbol B}$| the magnetic field, Ptot the total pressure (thermal plus magnetic), P the thermal pressure, |${\boldsymbol g}_{\text{s}}$| the effective gravity in the shearing-box approximation, Ω0 is the Keplerian angular velocity at the centre of the shearing box, Π the viscous stress tensor, E the energy density, e the internal energy per unit mass, η the resistivity, ν the shear viscosity (we assume η and ν are diagonal and isotropic tensors, i.e. scalars), |${\boldsymbol J}$| the electric current density, q is a local measure of the differential rotation, d lnΩ/d ln  R(q = 3/2 for a Keplerian disc) and Λ the rate of energy loss per unit volume due to radiative cooling. 2.1 Cooling We implement an effective cooling function in the simulation assuming that the opacity is dominated by electron scattering. This is likely to be the case in the inner regions of the disc, where we expect the Pm = 1 transition to occur (PB). Following Faulkner et al. 1983; Latter & Papaloizou 2012, we use a bulk cooling function to stand in for the diffusive cooling of the form \begin{equation} \Lambda =\frac{\sigma T_{\text{e}}^{4}}{H}, \end{equation} (8) where H is the disc scaleheight H = cs/Ω, cs is the thermal sound speed i.e. |$P=\rho c_{\text{s}}^{2}$|⁠, Ω is the Keplerian angular velocity and Te is the effective surface temperature given by \begin{equation} T_{\text{e}}^{4}=\frac{4T^{4}}{3\Sigma \kappa }. \end{equation} (9) Here T is the temperature of the simulated central disc plasma, κ = 0.4 cm2g−1 is the electron scattering opacity and Σ = ρH is the disc surface density. The cooling is then calculated using the temperature and density in each simulation cell. This differs slightly from the implementation in Latter & Papaloizou (2012), which used spatially averaged quantities and bulk cooling. We have tested both cooling methods to confirm that this choice does not affect our results. 2.2 Instability criterion and Pm In PB, we derived a thermal–viscous instability criterion for a standard, radiatively efficient thin disc in which the α parameter is a variable depending on ρ and Tc. This is \begin{equation} \frac{\mathrm{\partial} \ln \alpha }{\mathrm{\partial} \ln \Sigma }+\frac{1}{4}\frac{\mathrm{\partial} \ln \tau }{\mathrm{\partial} \ln \Sigma } + 1<0, \end{equation} (10) where τ = Σκ is the optical depth of the disc. Assuming a power-law dependence for α upon the magnetic Prandtl number, α ∝ Pmn, a parametrization suggested by simulations (Lesur & Longaretti 2007; Simon & Hawley 2009), the disc is unstable if 6/13 < n > 10/3 (PB, equations 46 and 47). Here, α is determined by the usual formula \begin{equation} \alpha P_{\mathrm{tot}}=\langle \rho (\delta v_{r} \delta v_{\phi }- v_{Ar}v_{A\phi })\rangle , \end{equation} (11) where δvi is the residual velocity after subtraction of the local Keplerian circular velocity (i.e. δvϕ = vϕ − ΩR, δvr = vr), vAi is the Alfvén velocity in the i-direction, |$v^{2}_{Ai}=B_{i}^{2}$| and the angle brackets denote a spatial average. For a typical BHB accretion disc environment, the radiative viscosity dominates the Coulomb viscosity, and one finds: \begin{eqnarray} \eta &=&\frac{5.6\times 10^{11}\ln \Lambda _{eH}}{T^{3/2}} \mathrm{cm}^{2}\mathrm{s}^{-1}, \nonumber\\ \nu _{\mathrm{Rad}}&=&\frac{6.7\times 10^{-26}T^{4}}{\kappa \rho ^{2}} \mathrm{cm}^{2}\mathrm{s}^{-1}, \end{eqnarray} (12) \begin{equation} {\rm Pm}\simeq \frac{\nu _{\mathrm{Rad}}}{\eta }=1.9\times 10^{-38}\frac{T^{11/2}}{\kappa \rho ^{2}}, \end{equation} (13) where κ is the opacity of the plasma and ln ΛeH is the electron–proton Coulomb collision factor, which we take to have a value of |$\sqrt{40}\simeq 6.3$| (see PB section 2.1.3 for more detail). 3 HOW DOES α DEPEND ON DISSIPATION AND RESOLUTION? In this section, we wish to address two questions important to the simulations: (i) to what extent do the viscous and resistive dissipation coefficients determine the strength of the saturated disc turbulence quantified by α and (ii) to what extent are these effects physical or numerical in origin? Recent results from several groups show that α depends on a whole array of physical and numerical parameters, such as the height of the simulating box, stratification, the presence of net magnetic fields, numerical resolution and convection, etc. (e.g. Hawley et al. 1995; Fromang et al. 2007; Simon, Beckwith & Armitage 2012; Hirose et al. 2014; Ryan et al. 2017). Here, we are interested in isolating the effect of dissipation coefficients. To avoid conflating this with other complications, we initially choose the simplest isothermal, unstratified local shearing box. We wish to study the effect of dissipation coefficients, numerical resolution and the initial magnetic field configuration. In these isothermal simulations, we maintain fixed dissipation coefficients throughout an individual run. (We shall later allow the viscosity and resistivity to become time-dependent functions of temperature and density.) In Figs 1–3 and Table 1, we summarize the results of the isothermal simulations for a variety of values of resistivity, viscosity, net magnetic field and resolution. The resistivity and viscosity can be expressed in terms of the dimensionless Reynolds number, Re, and magnetic Reynolds number Rm, given by \begin{equation} \mathrm{Rm}=\frac{c_{\text{s}}H}{\eta }, \qquad \mathrm{Re}=\frac{c_{\text{s}}H}{\nu }. \end{equation} (14) The magnetic Prandtl number is then \begin{equation} {\rm Pm}=\frac{\mathrm{Rm}}{\mathrm{Re}}. \end{equation} (15) Figure 1. Open in new tabDownload slide Average turbulent disc stress α for different values of Re and Rm with a net B-field in the azimuthal direction βϕ = 100. Large values of the resistivity have a more pronounced effect than either the viscosity or Pm (⁠|$\rm {Rm}/\rm {Re}$|⁠), which suggests that a resistivity Rm < 3200 affects the macroscale linear growth rate of the MRI. Values of Rm in excess of 12 800 appear to exhibit convergence in these low resolution runs (IsoBphi1, Table 1), essentially becoming insignificant compared to the numerical grid resistivity. Figure 1. Open in new tabDownload slide Average turbulent disc stress α for different values of Re and Rm with a net B-field in the azimuthal direction βϕ = 100. Large values of the resistivity have a more pronounced effect than either the viscosity or Pm (⁠|$\rm {Rm}/\rm {Re}$|⁠), which suggests that a resistivity Rm < 3200 affects the macroscale linear growth rate of the MRI. Values of Rm in excess of 12 800 appear to exhibit convergence in these low resolution runs (IsoBphi1, Table 1), essentially becoming insignificant compared to the numerical grid resistivity. Figure 2. Open in new tabDownload slide Average turbulent disc stress α for different values of the magnetic Prandtl number and net magnetic fields at fixed Rm = 12 800. Larger net magnetic field strengths lead to larger values of α in the low Pm regime. The dashed line shows the instability threshold (10), α–Pm dependencies with larger gradients than the dashed line are unstable. We expect the instability to occur for low net Bϕ fields with βϕ < 1000. The simulations of Pm < 1 show convergence to a constant value as the viscosity becomes dominated by the numerical grid viscosity. Simulation parameters are given in Table 1. Figure 2. Open in new tabDownload slide Average turbulent disc stress α for different values of the magnetic Prandtl number and net magnetic fields at fixed Rm = 12 800. Larger net magnetic field strengths lead to larger values of α in the low Pm regime. The dashed line shows the instability threshold (10), α–Pm dependencies with larger gradients than the dashed line are unstable. We expect the instability to occur for low net Bϕ fields with βϕ < 1000. The simulations of Pm < 1 show convergence to a constant value as the viscosity becomes dominated by the numerical grid viscosity. Simulation parameters are given in Table 1. Figure 3. Open in new tabDownload slide The effect of changing numerical resolution and the initial fixed value of Rm on the average turbulent disc stress α for a fixed range of values of the magnetic Prandtl number. All these simulations have a net Bϕ field with βϕ = 104 and Nz is the number of grid zones per disc scaleheight, H, in the z-direction. The α–Pm dependence is strong enough in all of these simulations at Pm > 1 to lead to a thermal instability (10), as indicated by possessing a gradient greater than the dotted lines. The results strongly support the hypothesis that there is a real physical α–Pm dependence since this dependence is maintained when simultaneously increasing the simulation resolution and decreasing the resistivity and viscosity (see the text for discussion). Figure 3. Open in new tabDownload slide The effect of changing numerical resolution and the initial fixed value of Rm on the average turbulent disc stress α for a fixed range of values of the magnetic Prandtl number. All these simulations have a net Bϕ field with βϕ = 104 and Nz is the number of grid zones per disc scaleheight, H, in the z-direction. The α–Pm dependence is strong enough in all of these simulations at Pm > 1 to lead to a thermal instability (10), as indicated by possessing a gradient greater than the dotted lines. The results strongly support the hypothesis that there is a real physical α–Pm dependence since this dependence is maintained when simultaneously increasing the simulation resolution and decreasing the resistivity and viscosity (see the text for discussion). Table 1. Table showing the simulation parameters used in this work, where βi = (cs/vAi)2 is the plasma beta of the initial net magnetic field in the i-direction and κ is the electron scattering opacity used in the cooling function in cgs units. Simulation set . βϕ . βz . Pm . Rm . Resolution (x,y,z) . E.O.S . ρ/ρ0 . κ . Pminitial . IsoBphi1 100 n.a. 1/2, 1, 2, 4 Variable 64 × 100 × 64 Isothermal n.a. n.a. n.a. IsoBphi2 100 n.a. 1/16, 1/4, 1, 4, 16 12 800 64 × 100 × 64 Isothermal n.a. n.a. n.a. IsoBphi3 1000 n.a. 1/16, 1/4, 1, 4, 16 12 800 64 × 100 × 64 Isothermal n.a. n.a. n.a. IsoBphi4 10 000 n.a. 1/16, 1/4, 1, 4, 16 12 800 64 × 100 × 64 Isothermal n.a. n.a. n.a. IsoBphi5 10 000 n.a. 1/16, 1/4, 1, 4, 16 12 800 128 × 200 × 128 Isothermal n.a. n.a. n.a. IsoBphi6 10 000 n.a. 1/16, 1/4, 1, 4, 16 25 600 128 × 200 × 128 Isothermal n.a. n.a. n.a. IsoBphi7 10 000 n.a. 1/16, 1/4, 1, 4, 16 51 200 256 × 400 × 256 Isothermal n.a. n.a. n.a. IsoBz1 n.a. 100 1/16, 1/4, 1, 4, 16 12 800 64 × 100 × 64 Isothermal n.a. n.a. n.a. IsoBz2 n.a. 1000 1/16, 1/4, 1, 4, 16 12 800 64 × 100 × 64 Isothermal n.a. n.a. n.a. IsoBz3 n.a. 10 000 1/16, 1/4, 1, 4, 16 12 800 64 × 100 × 64 Isothermal n.a. n.a. n.a. IsononetB n.a. n.a. 1/16, 1/4, 1, 4, 16 12 800 64 × 100 × 64 Isothermal n.a. n.a. n.a. IdealBphi1 10 000 n.a. variable (1/16-16) 25 600 128 × 200 × 128 Ideal 0.6|$\dot{7}$|⁠, 0.7|$\dot{3}$|⁠, 0.8, 0.9|$\dot{3}$|⁠, 1.0, 1.0|$\dot{6}$| 0.4 16 IdealBphi2 10 000 n.a. variable (1/16-16) 25 600 128 × 200 × 128 Ideal 0.6|$\dot{7}$|⁠, 0.7|$\dot{3}$|⁠, 0.8, 0.9|$\dot{3}$|⁠, 1.0, 1.0|$\dot{6}$| 0.4 1/16 Simulation set . βϕ . βz . Pm . Rm . Resolution (x,y,z) . E.O.S . ρ/ρ0 . κ . Pminitial . IsoBphi1 100 n.a. 1/2, 1, 2, 4 Variable 64 × 100 × 64 Isothermal n.a. n.a. n.a. IsoBphi2 100 n.a. 1/16, 1/4, 1, 4, 16 12 800 64 × 100 × 64 Isothermal n.a. n.a. n.a. IsoBphi3 1000 n.a. 1/16, 1/4, 1, 4, 16 12 800 64 × 100 × 64 Isothermal n.a. n.a. n.a. IsoBphi4 10 000 n.a. 1/16, 1/4, 1, 4, 16 12 800 64 × 100 × 64 Isothermal n.a. n.a. n.a. IsoBphi5 10 000 n.a. 1/16, 1/4, 1, 4, 16 12 800 128 × 200 × 128 Isothermal n.a. n.a. n.a. IsoBphi6 10 000 n.a. 1/16, 1/4, 1, 4, 16 25 600 128 × 200 × 128 Isothermal n.a. n.a. n.a. IsoBphi7 10 000 n.a. 1/16, 1/4, 1, 4, 16 51 200 256 × 400 × 256 Isothermal n.a. n.a. n.a. IsoBz1 n.a. 100 1/16, 1/4, 1, 4, 16 12 800 64 × 100 × 64 Isothermal n.a. n.a. n.a. IsoBz2 n.a. 1000 1/16, 1/4, 1, 4, 16 12 800 64 × 100 × 64 Isothermal n.a. n.a. n.a. IsoBz3 n.a. 10 000 1/16, 1/4, 1, 4, 16 12 800 64 × 100 × 64 Isothermal n.a. n.a. n.a. IsononetB n.a. n.a. 1/16, 1/4, 1, 4, 16 12 800 64 × 100 × 64 Isothermal n.a. n.a. n.a. IdealBphi1 10 000 n.a. variable (1/16-16) 25 600 128 × 200 × 128 Ideal 0.6|$\dot{7}$|⁠, 0.7|$\dot{3}$|⁠, 0.8, 0.9|$\dot{3}$|⁠, 1.0, 1.0|$\dot{6}$| 0.4 16 IdealBphi2 10 000 n.a. variable (1/16-16) 25 600 128 × 200 × 128 Ideal 0.6|$\dot{7}$|⁠, 0.7|$\dot{3}$|⁠, 0.8, 0.9|$\dot{3}$|⁠, 1.0, 1.0|$\dot{6}$| 0.4 1/16 Open in new tab Table 1. Table showing the simulation parameters used in this work, where βi = (cs/vAi)2 is the plasma beta of the initial net magnetic field in the i-direction and κ is the electron scattering opacity used in the cooling function in cgs units. Simulation set . βϕ . βz . Pm . Rm . Resolution (x,y,z) . E.O.S . ρ/ρ0 . κ . Pminitial . IsoBphi1 100 n.a. 1/2, 1, 2, 4 Variable 64 × 100 × 64 Isothermal n.a. n.a. n.a. IsoBphi2 100 n.a. 1/16, 1/4, 1, 4, 16 12 800 64 × 100 × 64 Isothermal n.a. n.a. n.a. IsoBphi3 1000 n.a. 1/16, 1/4, 1, 4, 16 12 800 64 × 100 × 64 Isothermal n.a. n.a. n.a. IsoBphi4 10 000 n.a. 1/16, 1/4, 1, 4, 16 12 800 64 × 100 × 64 Isothermal n.a. n.a. n.a. IsoBphi5 10 000 n.a. 1/16, 1/4, 1, 4, 16 12 800 128 × 200 × 128 Isothermal n.a. n.a. n.a. IsoBphi6 10 000 n.a. 1/16, 1/4, 1, 4, 16 25 600 128 × 200 × 128 Isothermal n.a. n.a. n.a. IsoBphi7 10 000 n.a. 1/16, 1/4, 1, 4, 16 51 200 256 × 400 × 256 Isothermal n.a. n.a. n.a. IsoBz1 n.a. 100 1/16, 1/4, 1, 4, 16 12 800 64 × 100 × 64 Isothermal n.a. n.a. n.a. IsoBz2 n.a. 1000 1/16, 1/4, 1, 4, 16 12 800 64 × 100 × 64 Isothermal n.a. n.a. n.a. IsoBz3 n.a. 10 000 1/16, 1/4, 1, 4, 16 12 800 64 × 100 × 64 Isothermal n.a. n.a. n.a. IsononetB n.a. n.a. 1/16, 1/4, 1, 4, 16 12 800 64 × 100 × 64 Isothermal n.a. n.a. n.a. IdealBphi1 10 000 n.a. variable (1/16-16) 25 600 128 × 200 × 128 Ideal 0.6|$\dot{7}$|⁠, 0.7|$\dot{3}$|⁠, 0.8, 0.9|$\dot{3}$|⁠, 1.0, 1.0|$\dot{6}$| 0.4 16 IdealBphi2 10 000 n.a. variable (1/16-16) 25 600 128 × 200 × 128 Ideal 0.6|$\dot{7}$|⁠, 0.7|$\dot{3}$|⁠, 0.8, 0.9|$\dot{3}$|⁠, 1.0, 1.0|$\dot{6}$| 0.4 1/16 Simulation set . βϕ . βz . Pm . Rm . Resolution (x,y,z) . E.O.S . ρ/ρ0 . κ . Pminitial . IsoBphi1 100 n.a. 1/2, 1, 2, 4 Variable 64 × 100 × 64 Isothermal n.a. n.a. n.a. IsoBphi2 100 n.a. 1/16, 1/4, 1, 4, 16 12 800 64 × 100 × 64 Isothermal n.a. n.a. n.a. IsoBphi3 1000 n.a. 1/16, 1/4, 1, 4, 16 12 800 64 × 100 × 64 Isothermal n.a. n.a. n.a. IsoBphi4 10 000 n.a. 1/16, 1/4, 1, 4, 16 12 800 64 × 100 × 64 Isothermal n.a. n.a. n.a. IsoBphi5 10 000 n.a. 1/16, 1/4, 1, 4, 16 12 800 128 × 200 × 128 Isothermal n.a. n.a. n.a. IsoBphi6 10 000 n.a. 1/16, 1/4, 1, 4, 16 25 600 128 × 200 × 128 Isothermal n.a. n.a. n.a. IsoBphi7 10 000 n.a. 1/16, 1/4, 1, 4, 16 51 200 256 × 400 × 256 Isothermal n.a. n.a. n.a. IsoBz1 n.a. 100 1/16, 1/4, 1, 4, 16 12 800 64 × 100 × 64 Isothermal n.a. n.a. n.a. IsoBz2 n.a. 1000 1/16, 1/4, 1, 4, 16 12 800 64 × 100 × 64 Isothermal n.a. n.a. n.a. IsoBz3 n.a. 10 000 1/16, 1/4, 1, 4, 16 12 800 64 × 100 × 64 Isothermal n.a. n.a. n.a. IsononetB n.a. n.a. 1/16, 1/4, 1, 4, 16 12 800 64 × 100 × 64 Isothermal n.a. n.a. n.a. IdealBphi1 10 000 n.a. variable (1/16-16) 25 600 128 × 200 × 128 Ideal 0.6|$\dot{7}$|⁠, 0.7|$\dot{3}$|⁠, 0.8, 0.9|$\dot{3}$|⁠, 1.0, 1.0|$\dot{6}$| 0.4 16 IdealBphi2 10 000 n.a. variable (1/16-16) 25 600 128 × 200 × 128 Ideal 0.6|$\dot{7}$|⁠, 0.7|$\dot{3}$|⁠, 0.8, 0.9|$\dot{3}$|⁠, 1.0, 1.0|$\dot{6}$| 0.4 1/16 Open in new tab 3.1 Initial conditions Explicit parameters of our numerical simulations are provided in Table 1. All simulations start with a zero net field in the z-direction, Bz = B0sin (2πx/H) and plasma beta, β0 = 100, where |$\beta _{0}=2 P/B_{0}^{2}$| is the ratio of thermal to magnetic pressure. We shall subsequently quantify the strength of any net magnetic field as |$\beta _{i}=2 P/B^{2}_{i}$|⁠, the ratio of thermal pressure to the magnetic pressure of the net field in the i-direction. The dimensionless orbital frequency and thermal sound speed are chosen to be Ω = cs = 0.001, so H = 1, and box lengths are Lx = H, Ly = 4H and Lz = H. Simulations were initialized with random pressure perturbations of maximum amplitude 1 per cent, as in Simon & Hawley (2009). 3.2 Isothermal simulations of the α–Pm dependence We first investigate disc turbulence for a variety of values of Rm and Re. Our findings are in basic agreement with previous investigations, which found that larger values of viscosity increase the value of α, whilst large resistivities decrease α (Fig. 1). There is clearly a complex interaction between the viscosity, resistivity, net B-field and resolution. Large values of the viscosity and resistivity, Re < 800, Rm < 800, result in sufficiently high dissipation on large lengthscales to prevent sustainable turbulence. These values are not representative of what would be expected in realistic astrophysical disc environments, and we avoid this regime. On the other hand, excessively small values of the dissipation coefficients (i.e. large Reynolds numbers) will produce no measurable effect in the simulations, since the relevant lengthscales on which the viscosity and resistivity become important will be below the grid resolution (and unavoidable numerical dissipation) of the simulations. Moreover, whilst an excessively large viscosity reduces the linear growth rate of the MRI, a large resistivity eliminates the linear MRI. We focus therefore on the intermediate asymptotic regime in which the resistivity is held constant at the lowest marginally resolvable value, and allow only the viscosity to vary in order to study the effect of altering Pm. 3.2.1 Numerical approach When the viscous and resistive lengthscales are well below the grid scale in a simulation, changing these dissipation coefficients will clearly have no effect on MRI turbulence. In this regime, there can be no α–Pm effect. If it is found that changing Pm does affect α, then we must be in a regime in which the dissipation lengthscales are partially or fully resolved. Our simulations formally span a large range, 1/16 ≥ Pm ≥ 16. Simulations for which Pm ⪅ 1 do not fully resolve the viscous lengthscale. This can be seen in Figs 2 and 3, where the α–Pm dependence flattens at low Pm. This flattening occurs because the viscous dissipation falls below the numerical grid dissipation in the Pm ≪ 1 regime, and thus has a diminishing effect on the turbulence. In general, it is expected that the α–Pm effect should be most pronounced close to Pm ∼ 1, since this is the transition between dominant viscous dissipation and dominant resistive dissipation. Simulations by Meheut et al. (2015) investigated the small Pm regime and found the α–Pm dependence continues to decrease down to Pm ≈ 1/4, before levelling off to a constant α for Pm ⪅ 0.1. These simulations used βz = 103 and Rm =400 (a relatively large resistivity) and required very powerful computational resources in order to fully resolve the viscous lengthscale even at moderately small Pm. In this paper, we have chosen to work in the opposite regime and probe the α–Pm dependence in the regime Pm ≳ 1. The main reason for this choice is so that we are able to use lower values of the resistivity, ensuring that the resistivity is not so large that it significantly affects the macroscopic linear growth rate of the MRI. The intent of our simulations is twofold: first to establish whether the α–Pm effect is physical or numerical in origin and secondly to test whether the α–Pm dependence is sufficiently sensitive to induce the Pm disc instability. The values of viscosity and resistivity are necessarily unrealistically large in all current simulations due to limited resolution and computational resources. To establish that the α–Pm effect is physical and not an artefact of these unrealistically large dissipation coefficients, we use the lowest possible viscosity and resistivity that are still adequately resolved in simulations with a broad range of numerical resolutions. The critical point is that if an α–Pm dependence remains unaffected even when the simulation resolution is increased and the absolute values of the viscosity and resistivity are decreased, then the α–Pm dependence is not a numerical artefact. This approach has the advantage that the dissipation coefficients are small to start and decrease further at higher numerical resolution. This increases the scale separation between the large-scale fastest growing MRI modes and the small-scale dissipation lengthscales at higher resolution, allowing us to test whether the α–Pm dependence is physical or numerical. 3.2.2 Results The results of varying Pm (with fixed Rm) are shown in Figs 2 and 3 and Table 1. The results, as expected, show an increase in the strength of the saturated MRI turbulence for values of Pm > 1. At small values of Pm, α tends to a constant value determined by the initial net B-field in the simulation, Rm, and the resolution of the simulation. The presence of an imposed initial net magnetic field is expected to increase the strength of saturated disc turbulence because it provides a sort of ‘backbone’ magnetic field which cannot be dissipated and is thus a permanent source of magnetic field. The initial net magnetic field is crucial in determining the minimum value of α at small Pm. Net Bϕ fields have a strong impact on the α value at Pm ≲ 1, with larger Bϕ fields increasing the value of α. Net Bz fields increase α at all values of Pm (the strongly enhancing effect of a net Bz field has been known for some time; Hawley et al. 1995). To establish the dependence of α on Pm as a physical, rather than a numerical, effect, we need to show that variations in α exist which depend only upon the ratio of dissipation coefficients (i.e. Pm) and not on their absolute values. This is done by comparing the results from simulations covering the same range of Pm with the same initial conditions, but as the numerical resolution of the simulation is increased, the absolute values of viscosity and resistivity are simultaneously decreased. Results are shown in Fig. 3, where it can be seen that using the same initial setup, the simulations at higher numerical resolution (which at a given value of Pm have substantially decreased values of both viscous and resistive dissipation coefficients) obtain remarkably similar results to the lower resolution simulations. If the Pm-α effect were due to artificially large dissipation lengthscales and insufficient scale separation, we would expect the effect to decrease with increased resolution and decreased values of the dissipation coefficients. Fig. 3 shows this is clearly not the case. In fact, the Pm-α effect is slightly more pronounced for Pm < 1 at higher resolution, even when the absolute dissipation coefficients are smaller. This is a clear indication that whatever residual numerical effects may be present, there appears to be a distinct and genuine physical Pm effect at work (Fig. 4). Figure 4. Open in new tabDownload slide Examples of the magnetic field structure in the case of (a): a low Pm equilibrium solution, (b): a high Pm equilibrium solution. The radial magnetic field component is shown in an x–z slice at the centre of the box for simulations (IsoBphi7) with Pm = 1/4 and Pm = 16, respectively. The magnetic field is dominated by larger scale structure in (b) with Pm = 16, than in (a) with Pm = 1/4. This is due to the larger viscosity at larger Pm, which damps small-scale turbulent velocities. This affects the structure of the magnetic field, even though the resistivity is the same in the two simulations. This is evidence in favour of the α–Pm dependence being caused, in part, by a lower reconnection rate at high Pm. The dimensionless magnetic field strength on the colour bar can be compared to the dimensionless isothermal sound speed cs = 0.001. Figure 4. Open in new tabDownload slide Examples of the magnetic field structure in the case of (a): a low Pm equilibrium solution, (b): a high Pm equilibrium solution. The radial magnetic field component is shown in an x–z slice at the centre of the box for simulations (IsoBphi7) with Pm = 1/4 and Pm = 16, respectively. The magnetic field is dominated by larger scale structure in (b) with Pm = 16, than in (a) with Pm = 1/4. This is due to the larger viscosity at larger Pm, which damps small-scale turbulent velocities. This affects the structure of the magnetic field, even though the resistivity is the same in the two simulations. This is evidence in favour of the α–Pm dependence being caused, in part, by a lower reconnection rate at high Pm. The dimensionless magnetic field strength on the colour bar can be compared to the dimensionless isothermal sound speed cs = 0.001. The formal instability criterion is ∂ln α/∂ln Pm > 6/13 (PB). From Fig. 3 this is satisfied in the Pm > 1 regime for the high resolution simulations with small mean magnetic field strengths. 4 THERMAL EQUILIBRIUM CURVE Let us now investigate whether the dependence of α on Pm leads to the thermal–viscous instability outlined in PB. To test this, the local thermal-equilibrium curve of the disc is simulated. It is necessary to include both turbulent heating and radiative cooling (the gas is no longer isothermal), as well as a variable magnetic Prandtl number which has the correct dependence on temperature and density calculated in (13). To minimize the direct effect of a large resistivity on the linear MRI growth rate, we adopt a constant low value for the resistivity (Rm ≃ 25600) and include the temperature and density dependence of Pm in the viscosity alone. The explicit resistivity and viscosity are given by \begin{equation} \eta =\frac{c_{\text{s}0}H_{0}}{25\,600}, \qquad \nu =\frac{c_{\text{s}0}H_{0}}{25\,600}\left(\frac{T}{T_{0}}\right)^{11/2} \left(\frac{\rho }{\rho _{0}}\right)^{-2}. \end{equation} (16) Instability, if present, will manifest itself as two sets of overlapping quasi-stable solutions (stable on the thermal time-scale but evolving on the mass accretion time-scale, hence the term thermal–viscous instability.) To find such solutions, we carry out two simulations with nearly identical initial conditions, the sole difference being that one simulation starts with a fixed high Pm = 16 and one with a fixed low Pm = 1/16. The MRI is allowed to grow, become fully turbulent and attain a dynamical-thermal equilibrium in the two cases. From Fig. 2, the high Pm simulation is expected to have stronger turbulent heating due to the larger α value and so to reach a higher equilibrium temperature than the low Pm simulation. The physical size of both low and high Pm simulations is set to be equal to the equilibrium disc scaleheight in the hotter, high Pm simulation i.e. Lz = H (Pm = 16). This ensures that the scalelength of the maximally growing MRI modes are captured in both sets of simulations; this is important for correctly determining α as the pressure changes (e.g. Sano et al. 2004; Ross, Latter & Guilet 2016). The analytic Pm instability criterion is found to be rather insensitive to the assumed linear scaling of the disc stress with pressure. The close agreement between simulations and analytic calculation evinced in Fig. 5 is further confirmation of this insensitivity. Choosing both simulations to have the same physical size has the advantage that we do not introduce an intrinsic bias towards finding two sets of solutions from two different initial conditions. This means that the simulations will not tend to find an instability if none actually exists, and in the regime in which only one stable solution exists, both simulations will converge to the same final equilibrium. The disadvantage is that the simulation size (Lx, Ly, Lz) of the cooler low Pm solution will be slightly larger than the disc scaleheight. This is because the disc scaleheight of the hotter, high Pm solution is slightly larger than that of the cooler, low Pm solution, i.e. H (Pm = 16) ≈ 1.7H (Pm = 1/16) from Fig. 5. Figure 5. Open in new tabDownload slide Local thermal equilibrium curve showing the two overlapping sets of stable solutions for low and high Pm. The simulations show the characteristic unstable S-curve associated with a thermal–viscous instability. The black lines show the expected thermal equilibrium curve obtained by solving (17), with solid indicating stable and dashed indicating unstable solutions. This demonstrates the existence of the Pm instability in a self-consistent set of MHD simulations. Figure 5. Open in new tabDownload slide Local thermal equilibrium curve showing the two overlapping sets of stable solutions for low and high Pm. The simulations show the characteristic unstable S-curve associated with a thermal–viscous instability. The black lines show the expected thermal equilibrium curve obtained by solving (17), with solid indicating stable and dashed indicating unstable solutions. This demonstrates the existence of the Pm instability in a self-consistent set of MHD simulations. After equilibrium has been achieved, Pm is freed and its value is calculated via equation (16) for each grid cell. Thus, Pm becomes a function of the local fluid temperature and density. The two simulations are then evolved until they relax to a new quasi-stable thermal equilibrium, if this in fact differs from the initial equilibrium. If the Pm instability is present in the simulation, this should manifest itself as two stable thermal equlibria at the same surface density Σ, producing the characteristic unstable S-curve, as in the case of the hydrogen ionization instability in dwarf novae (e.g. Lasota 2001; Frank, King & Raine 2002 and PB). Fig. 5 shows the local thermal equilibrium curve of the instability calculated using a set of simulations with a weak net magnetic field in the azimuthal direction (βϕ = 104). Two stable thermal equilibria exist at the same density in the unstable range of Pm (Pm ≈ 1 − 16). This is the first self-consistent MHD simulation demonstrating the Pm instability from first principles. It shows that the dependence of α on Pm in a full MHD simulation is sufficiently sensitive to induce a thermal–viscous instability in the disc plasma. The black curve in the figure is the expected theoretical thermal equilibrium curve calculated by balancing turbulent heating and radiative cooling. This is in close agreement with the results of the simulation. The black curve is calculated from the thermal equilibrium equation \begin{equation} \frac{\sigma T_{\text{e}}^{4}}{H}=\frac{3}{2}\Omega \alpha\, ({\rm Pm}) P_{\mathrm{tot}}, \end{equation} (17) where the LHS corresponds to radiative cooling and the RHS to turbulent heating (see e.g. Balbus & Hawley 1998). The value of alpha is determined as a function of Pm for the black curve by using the α–Pm dependence found in simulation IsoBphi6, which is the equivalent isothermal simulation. Due to the finite resolution of the simulation, turbulent fluctuations in the spatially averaged values of α and T are much larger than would be expected in an actual disc (in part, because of Poisson noise from the finite number of grid points in the simulation) and so this precludes accurately simulating the solutions close to the edges or corners of the solid S-curve. Close to the edges of the upper and lower branches of solutions small temporal fluctuations in temperature are sufficient to move the average temperature of the plasma over the dashed line of solutions which are unstable to perturbations. The dashed line represents the watershed between states which will experience runaway heating (above the dashed line) or cooling (below the dashed line) until they reach an equilibrium solution on the solid curve. Examples of the evolution of thermally stable and unstable simulations are shown in Fig. 6. The precise properties of the thermal-equilibrium curve are clearly dependent on the initial net magnetic field, which is likely to vary between different astrophysical systems. Figure 6. Open in new tabDownload slide The spatially averaged temperature as a function of time for two sets of simulations shown in Fig. 5. The red and green curves show the two stable thermal equilibria which exist at high and low Pm at a density of |$\rho =0.9\dot{3}\rho _{0}$|⁠. At a density of |$0.\dot{6}\rho _{0}$| only the low Pm thermal equilibrium solution exists and so both the blue and pink curve converge on the low Pm solution (after the initial high fixed Pm value is freed at ≈80Ω−1 for the simulation shown in pink). Figure 6. Open in new tabDownload slide The spatially averaged temperature as a function of time for two sets of simulations shown in Fig. 5. The red and green curves show the two stable thermal equilibria which exist at high and low Pm at a density of |$\rho =0.9\dot{3}\rho _{0}$|⁠. At a density of |$0.\dot{6}\rho _{0}$| only the low Pm thermal equilibrium solution exists and so both the blue and pink curve converge on the low Pm solution (after the initial high fixed Pm value is freed at ≈80Ω−1 for the simulation shown in pink). The limited resolution of our simulations prevents us from capturing the entire expected range of variation of α with Pm. The α–Pm dependence continues to change steeply down to Pm ≈ 1/4 in the high-resolution simulations by Meheut et al. (2015). This means that higher resolution simulations able to resolve a larger range of the α–Pm variation would almost certainly enhance the instability, strengthening the results of this paper. It is therefore important to conduct further simulations of the Pm disc instability at higher resolution in order to precisely determine the local thermal-equilibrium curve and physical effects induced by this instability. 5 CONCLUSION In this paper, we address two important questions: (i) to what extent does the turbulent disc stress depend on the magnetic Prandtl number, and is this effect numerical or physical? and (ii) is the α–Pm dependence sufficiently sensitive to trigger a thermal–viscous instability in the disc? In the first section of the paper, we address (i) by conducting a suite of isothermal 3D MHD local shearing-box simulations to investigate the dependence of the turbulent disc stress α on the viscous and resistive dissipation coefficients. In agreement with previous studies we find that α depends on the magnetic Prandtl number (Pm), the ratio of viscosity to resistivity. In the regime 1 < Pm < 16, α increases with Pm and the dependence is sufficiently sensitive to induce an important new thermal–viscous instability, the magnetic Prandtl number disc instability (see PB). We investigate whether the α–Pm dependence is physical or numerical in origin by conducting a suite of simulations covering the same range in Pm, but at different numerical resolutions and with different absolute values of viscosity and resistivity. If the α–Pm dependence were numerical and caused by artificially large dissipation coefficients or a lack of scale separation in the turbulent cascade, the dependence should decrease as the simulation resolution is increased and the dissipation coefficients are decreased. In fact, the α–Pm dependence persists as the simulation resolution is increased and the absolute values of the dissipation coefficients are decreased. This shows that the Pm-α effect is not a numerical artefact caused by artificially large dissipation coefficients or a lack of scale separation. It is firm evidence that the α–Pm dependence is physical in origin. To investigate whether the α–Pm dependence is sufficiently sensitive to induce the magnetic Prandtl number disc instability we conducted a further set of MHD simulations which include realistic turbulent heating and radiative cooling to map out the local thermal-equilibrium curve of the disc. In these simulations, Pm is no longer held fixed; instead it is given a realistic physical dependence on the plasma density and temperature. Significantly, we find that the α–Pm dependence is sufficiently strong to trigger the Pm instability using physical conditions appropriate for an accretion disc. The thermal-equilibrium curve is found to have a characteristic S-shape forming an unstable limit cycle. In the unstable region corresponding to the regime 1 < Pm < 16, the instability manifests itself as two sets of overlapping stable solutions with different equilibrium temperatures but the same surface density. The two sets of stable solutions correspond to a hotter, high α branch of solutions with large Pm and a cooler, low α branch of solutions corresponding to Pm ≲ 1 as predicted in PB. This is the first self-consistent MHD simulation demonstrating the existence of the Pm disc instability from first principles. It was shown in PB that the local Pm instability leads to the production of global heating and cooling fronts propagating through the disc, resulting in changes to the disc state on time-scales substantially shorter than those of the DIM. This makes the magnetic Prandtl number disc instability a good potential candidate to explain the hard/soft changes in disc state observed on week long time-scales. A detailed comparison of the expected time-scale and spectral properties of the global instability will be the subject of future work. Acknowledgements WJP is supported by a Junior Research Fellowship from University College, University of Oxford. SAB acknowledges support from the Royal Society in the form of a Wolfson Research Merit Award. We acknowledge support from STFC in the form of a Consolidated Grant (ST/N000919/1) to Oxford Astrophysics. The simulations presented here were run on the DiRAC Complexity Cluster at Leicester and on Berg, the DiRAC facility jointly funded by STFC, the Large Facilities Capital Fund of BIS and the University of Oxford. WJP is particularly grateful to Julien Faure for important numerical advice and to Sebastien Fromang for hosting a visit to the CEA Saclay. We would like to thank Julien Faure, Sebastien Fromang, Charles Gammie, Rob Fender, John Miller, Jean-Pierre Lasota and Alexander Schekochihin for helpful conversations and suggestions. REFERENCES Balbus S. A. , Hawley J. F., 1998 , Rev. Mod. Phys. , 70 , 1 Crossref Search ADS Balbus S. A. , Henri P., 2008 , ApJ , 674 , 408 Crossref Search ADS Bodo G. , Mignone A., Cattaneo F., Rossi P., Ferrari A., 2008 , A&A , 487 , 1 Crossref Search ADS Bodo G. , Cattaneo F., Mignone A., Rossi P., 2014 , ApJL , 787 , L13 Crossref Search ADS Done C. , Gierliński M., Kubota A., 2007 , A&ARv , 15 , 1 Crossref Search ADS Faulkner J. , Lin D. N. C., Papaloizou J., 1983 , MNRAS , 205 , 359 Crossref Search ADS Fender R. P. , Belloni T. M., Gallo E., 2004 , MNRAS , 355 , 1105 Crossref Search ADS Frank J. King A. Raine D. J., 2002 , Accretion Power in Astrophysics : Third Edition Google Scholar Crossref Search ADS Google Scholar Google Preview WorldCat COPAC Fromang S. , Papaloizou J., Lesur G., Heinemann T., 2007 , A&A , 476 , 1123 Crossref Search ADS Hameury J.-M. , Menou K., Dubus G., Lasota J.-P., Hure J.-M., 1998 , MNRAS , 298 , 1048 Crossref Search ADS Hawley J. F. , Gammie C. F., Balbus S. A., 1995 , ApJ , 440 , 742 Crossref Search ADS Hirose S. , Blaes O., Krolik J. H., Coleman M. S. B., Sano T., 2014 , ApJ , 787 , 1 Crossref Search ADS Lasota J.-P. , 2001 , New Astron. Rev. , 45 , 449 Crossref Search ADS Latter H. N. , Papaloizou J. C. B., 2012 , MNRAS , 426 , 1107 Crossref Search ADS Lesur G. , Longaretti P.-Y., 2007 , MNRAS , 378 , 1471 Crossref Search ADS Meheut H. , Fromang S., Lesur G., Joos M., Longaretti P.-Y., 2015 , A&A , 579 , A117 Crossref Search ADS Mignone A. , Bodo G., Massaglia S., Matsakos T., Tesileanu O., Zanni C., Ferrari A., 2007 , ApJS , 170 , 228 Crossref Search ADS Mignone A. , Flock M., Stute M., Kolb S. M., Muscianisi G., 2012 , A&A , 545 , A152 Crossref Search ADS Potter W. J. , Balbus S. A., 2014 , MNRAS , 441 , 681 Crossref Search ADS Remillard R. A. , McClintock J. E., 2006 , ARA&A , 44 , 49 Crossref Search ADS Ross J. , Latter H. N., Guilet J., 2016 , MNRAS , 455 , 526 Crossref Search ADS Ryan B. R. Gammie C. F. Fromang S. Kestener P., 2017 , preprints (arXive-prints) Sano T. , Inutsuka S.-i., Turner N. J., Stone J. M., 2004 , ApJ , 605 , 321 Crossref Search ADS Schekochihin A. A. , Cowley S. C., Taylor S. F., Maron J. L., McWilliams J. C., 2004 , ApJ , 612 , 276 Crossref Search ADS Shakura N. I. , Sunyaev R. A., 1973 , A&A , 24 , 337 Simon J. B. , Hawley J. F., 2009 , ApJ , 707 , 833 Crossref Search ADS Simon J. B. , Beckwith K., Armitage P. J., 2012 , MNRAS , 422 , 2685 Crossref Search ADS Takahashi H. R. , Masada Y., 2011 , ApJ , 727 , 106 Crossref Search ADS © 2017 The Authors Published by Oxford University Press on behalf of the Royal Astronomical Society
TI - Demonstration of a magnetic Prandtl number disc instability from first principles
JF - Monthly Notices of the Royal Astronomical Society
DO - 10.1093/mnras/stx2055
DA - 2017-12-11
UR - https://www.deepdyve.com/lp/oxford-university-press/demonstration-of-a-magnetic-prandtl-number-disc-instability-from-first-xI0cK4Oe15
SP - 3021
VL - 472
IS - 3
DP - DeepDyve
ER -