  The parameterization of small-scale turbulent fluctuations in convective
systems and in the presence of strong stratification is a key issue for many
applied problems in oceanography, atmospheric science and planetology. In the
presence of stratification, one needs to cope with bulk turbulent fluctuations
and with inversion regions, where temperature, density -or both- develop highly
non-linear mean profiles due to the interactions between the turbulent boundary
layer and the unmixed -stable- flow above/below it. We present a second order
closure able to cope simultaneously with both bulk and boundary layer regions,
and we test it against high-resolution state-of-the-art 2D numerical
simulations in a convective and stratified belt for values of the Rayleigh
number, up to Ra = 10^9. Data are taken from a Rayleigh-Taylor system confined
by the existence of an adiabatic gradient.
