# How should I set the turbulent prandtl number in a turbulent compressible flow simulation?

I'm working on simulating compressible fluid flow over a flat plate with turbulence. The energy equation for compressible flow with turbulence requires modeling a turbulent heat flux of the form

$$q_T=\frac{\mu_T c_p}{Pr_T}\frac{\partial T}{\partial x_j}$$

where $\mu_T, c_p,Pr_T$ are the turbulent viscosity, specific heat (constant pressure) and turbulent Prandtl number. Clearly, all three of these paramters vary throughout the domain and require a model to determine their values as a function of the flow variables (density $\rho$, velocity $u$, pressure $p$, and temperature $T$. The turbulent viscosity is estimated using the variables of the specifically chosen turbulence model. The specific heat $c_p$ is a function of temperature, and can be estimated by interpolating over the most recent temperature values.

However, I have no idea how to estimate the turbulent Prandtl number as a function of the temperature, velocity, density, or pressure. The only relevant discussion I've found on this issue is in Wilcox (2006)'s Turbulence Modeling for CFD (chapter 5, section 5.4.2, pg 250), where it states

A constant value for $Pr_T$ is often used and this is usually satisfactory for shock-free flows up to low supersonic speeds, provided the heat transfer rate is not too high. The most common values assumed for $Pr_T$ are 0.89 or 0.90, in the case of a boundary layer. Heat transfer predictions can usually be improved somewhat by letting $Pr_T$ vary through the boundary layer. Near the edge of a boundary layer and throughout a free shear layer, a value of the order of 0.5 is more appropriate for $Pr_T$.

While modeling $Pr_T$ as a constant is OK for some cases, I want to know if there is a standard model to prescribe its value as a function of the flow variables $T,u,\rho,$ and $p$. Preferably, this model should ascribe values near 0.9 within the boundary layer and 0.5 throughout a free shear layer, as described by Wilcox. How should I prescribe the turbulent Prandtl number throughout the domain as a function of the flow variables?