Minimizing the numerical viscosity in smoothed particle hydrodynamics simulations of discs

Simulations using the smoothed particle hydrodynamics (SPH) technique typically include numerical viscosity to model shocks and maintain particle order on the kernel scale. This numerical viscosity is composed of linear and quadratic terms, with coefficients αSPH and βSPH, respectively. Setting these coefficients too high results in excessive numerical dissipation, whereas setting them too low may lead to unwanted effects such as particle penetration, which also leads to excess dissipation. In this study, we simulate accretion discs using the SPH code PHANTOM to investigate the effective disc viscosity arising from numerical viscosity. We model steady-state coplanar and circular discs with different values of αSPH and βSPH, from which we determine the coefficients that lead to minimum levels of numerical viscosity by maximizing the steady-state disc surface density for the same mass input rate. We find that, for planar and circular discs, the default values of the numerical viscosity parameters in the PHANTOM code can be too high particularly for the quadratic term. As higher values of the coefficients are required to adequately capture strong shocks in the flow, we suggest that the coefficient of the quadratic term should be time-dependent in a similar manner to the presently used ‘switches’ on the linear term. This can be simply achieved by setting βSPH to be a constant multiple of αSPH with αSPH determined by an appropriate switch, as previously advocated in the literature.