  In this paper, we are interested in constructing a scheme solving
compressible Navier--Stokes equations, with desired properties including high
order spatial accuracy, conservation, and positivity-preserving of density and
internal energy under a standard hyperbolic type CFL constraint on the time
step size, e.g., $\Delta t=\mathcal O(\Delta x)$. Strang splitting is used to
approximate convection and diffusion operators separately. For the convection
part, i.e., the compressible Euler equation, the high order accurate
postivity-preserving Runge--Kutta discontinuous Galerkin method can be used.
For the diffusion part, the equation of internal energy instead of the total
energy is considered, and a first order semi-implicit time discretization is
used for the ease of achieving positivity. A suitable interior penalty
discontinuous Galerkin method for the stress tensor can ensure the conservation
of momentum and total energy for any high order polynomial basis. In
particular, positivity can be proven with $\Delta t=\mathcal{O}(\Delta x)$ if
the Laplacian operator of internal energy is approximated by the $\mathbb{Q}^k$
spectral element method with $k=1,2,3$. So the full scheme with $\mathbb{Q}^k$
($k=1,2,3$) basis is conservative and positivity-preserving with $\Delta
t=\mathcal{O}(\Delta x)$, which is robust for demanding problems such as
solutions with low density and low pressure induced by high-speed shock
diffraction. Even though the full scheme is only first order accurate in time,
numerical tests indicate that higher order polynomial basis produces much
better numerical solutions, e.g., better resolution for capturing the roll-ups
during shock reflection.
