Thermally driven flows can be approximated as constant property flows only when temperature differences are relatively small. In that case, the Oberbeck-Boussinesq approximation holds and the set of governing equations is simplified. For larger temperature differences, the variation of the fluid properties with temperature cannot be ignored and in the case of gasses the governing equations take the low Mach form. This significantly complicates the numerical solution of variable property against constant property flows because of the emergence of a variable coefficient Poisson equation for the pressure. In the present study, a numerical methodology for the direct numerical simulation of thermally driven low-Mach flows is presented. In this framework, a pressure-splitting scheme is utilised to transform the variable coefficient Poisson equation for the pressure into a constant coefficient Poisson equation, improving the efficiency of the numerical solution. The consistent boundary conditions for the pressure are derived and presented. Furthermore, the proposed methodology is validated for the natural convection of air inside a differentially heated cavity, for a wide range of temperature differences, both within and outside the limits of applicability of the Oberbeck-Boussinesq approximation. Finally, to demonstrate the potential of this methodology, the three-dimensional natural convection of air inside a differentially heated cavity is simulated for a Rayleigh number Ra = 2.0 x 10(9) and for temperature differences Delta T = 50 K, 100 K, and 200 K. It is found that, with increasing temperature difference, the symmetry around the centre of the cavity that characterises the Oberbeck-Boussinesq solution is lost. In addition, the laminar-turbulent transition point on the heated and cooled walls changes position, moving to a more upstream position on the heated wall and a more downstream position on the cooled wall. (C) 2018 Elsevier Ltd. All rights reserved.