  High Mach number shocks are ubiquitous in interstellar turbulence. The Pencil
Code is particularly well suited to the study of magnetohydrodynamics in weakly
compressible turbulence and the numerical investigation of dynamos because of
its high-order advection and time evolution algorithms. However, the high-order
algorithms and lack of Riemann solver to follow shocks make it less well suited
to handling high Mach number shocks, such as those produced by supernovae
(SNe). Here, we outline methods required to enable the code to efficiently and
accurately model SNe, using parameters that allow stable simulation of
SN-driven turbulence, in order to construct a physically realistic galactic
dynamo model. These include the resolution of shocks with artificial viscosity,
thermal conductivity, and mass diffusion; the correction of the mass diffusion
terms; and a novel generalization of the Courant condition to include all
source terms in the momentum and energy equations. We test our methods with the
numerical solution of the one-dimensional (1D) Riemann shock tube (Sod, J.
Comput. Phys. 1978, 27), also extended to a 1D adiabatic shock with parameters
and Mach number relevant to SN shock evolution, including shocks with radiative
losses. We extend our test with the three-dimensional (3D) numerical simulation
of individual SN remnant evolution for a range of ambient gas densities typical
of the interstellar medium and compare these to the analytical solutions of
Sedov-Taylor (adiabatic) and the snowplough and Cioffi, McKee and Bertschinger
(Astrophys. J. 1988, 334) results incorporating cooling and heating processes.
We show that our new timestep algorithm leads to linear rather than quadratic
resolution dependence as the strength of the artificial viscosity varies,
because of the corresponding change in the strength of interzone gradients.
