This paper proposes a fractional step method for the calculation of compressible Navier–Stokes equations. The purpose of this study is to develop a robust and efficient numerical method for the simulation of low Mach number flows in which the poorly distributed eigenvalues usually result in the numerical difficulties. The method takes advantage of the pressure-based and the density-based methods. It thus accelerates the numerical convergence by adjusting the eigenvalues of the Jacobian matrix. In order to control the numerical instability and the spurious wave reflections at the numerical boundaries, a characteristic boundary condition is formulated. This numerical method as well as the boundary condition treatment extend the traditional fractional step method to the simulation of compressible Navier–Stokes equations. The performance and the accuracy of the method to calculate the flows with large amplitude of acoustic waves and strong heat transfer have been demonstrated by the direct numerical simulation of three distinct cases: a one-dimensional Euler equation with large amplitude of acoustic waves, an adiabatic turbulent boundary layer, and a turbulent boundary layer with heated wall.