The paper presents a novel pressure-corrected formulation of the immersed boundary method (IBM) for the simulation of fully compressible non-Boussinesq natural convection flows. The formulation incorporated into the pressure-based fractional step approach facilitates simulation of the flows in the presence of an immersed body characterized by a complex geometry. Here, we first present extensive grid independence and verification studies addressing incompressible pressure-driven flow in an extended channel and non-Boussinesq natural convection flow in a differentially heated cavity. Next, the steady-state non-Boussinesq natural convection flow developing in the presence of hot cylinders of various diameters placed within a cold square cavity is thoroughly investigated. The obtained results are presented and analyzed in terms of the spatial distribution of path lines and temperature fields and of heat flux values typical of the hot cylinder and the cold cavity surfaces. Flow characteristics of multiple steady-state solutions discovered for several configurations are presented and discussed in detail.