This study presents the development of a TKE-l parameterization of the diffusion coefficients for the representation of turbulent diffusion in neutral and stable conditions in large-scale atmospheric models. The parameterization has been carefully designed to be completely tunable in the sense that all adjustable parameters have been clearly identified and their number minimized as much as possible to help the calibration and to thoroughly assess the parametric sensitivity. We choose a mixing length formulation that depends on both static stability and wind shear to cover the different regimes of stable boundary layers. We follow a heuristic approach for expressing the stability functions and turbulent Prandlt number in order to guarantee the versatility of the scheme and its applicability for planetary atmospheres composed of an ideal and perfect gas such as that of Earth and Mars. Particular attention has also been paid to the numerical stability at typical time steps used in General Circulation Models. Test, parametric sensitivity assessment and preliminary tuning are performed on single-column idealized simulations of the weakly stable boundary layer. The robustness and versatility of the scheme are also assessed through its implementation in the LMDZ General Circulation Model and the Mars Planetary Climate Model and by running simulations of the Antarctic and Martian nocturnal boundary layers.