Turbulent heat fluxes are usually modeled by algebraic relations containing the turbulent Prandtl number. The underlying hypothesis of an analogy between turbulent momentum and heat transport is, however, not valid in complex flow situations. Transport equations are then needed to predict the heat fluxes. Two such equations are implemented in a code for statistically two-dimensional flows with a Reynolds stress model. The behaviour of this SOC model is tested in cases of compressible boundary layers along adiabatic walls at different Mach numbers and for hypersonic boundary layers along cooled walls. Finally, shock/turbulent boundary layer interactions are computed and discussed.