Reynolds-averaged turbulence models can result in large numerical error at flow discontinuities like shock waves. This is due to the nonconservative nature of the source terms in the governing equations. In this paper, the k-ϵ turbulence model is used to compute the canonical interaction of a normal shock with homogeneous isotropic turbulence. The characteristics of the nonconservative error is studied as a function of shock strength, grid resolution, and upstream turbulence level. The model equations are cast in an equivalent conservation form that gives physically consistent results at a shock wave. The predicted amplifications of turbulent kinetic energy and its dissipation rate match direct numerical simulation data and linear theory estimates over a range of Mach numbers.