The three-dimensional flow transition is examined in the wake of a heated square cylinder subjected to horizontal cross-flow perpendicular to gravity utilizing a direct numerical simulation approach. The surface of the square cylinder is heated uniformly to an elevated temperature Tw, and the amount of excess temperature is represented as the over-heat ratio ε=(Tw−T∞)/T∞, where T∞ represents the surrounding temperature. The effects of large-scale heating on the transport properties and thermal straining of the fluid particles are captured using an in-house non-Oberbeck–Boussinesq compressible model. The compressible flow governing equations (in a body-fitted coordinate system) are solved using a variant of flux-based particle velocity upwind-modified+ (PVU-M+) technique [Ahmad et al., “On the formation and sustenance of the compressible vortex rings in starting axisymmetric jets: A phenomenological approach,” Phys. Fluids 32, 126114 (2020)]. In this investigation, all computations are conducted at a low Mach number (Ma = 0.1) and air (Prandtl number, Pr = 0.71) is used as the working fluid. As the heating level rises, the shape and wavelength of the vortical structure undergo significant alterations. At Re = 250, the mode-B transition with a shorter spanwise wavelength and the mode-D transition with a longer wavelength are observed, respectively, for heating levels ε=0.0−0.2 and ε=0.8−1.0. Furthermore, for heating levels in the range 0.4≤ε≤0.6, an intermediate wavelength of the mode-E transition is detected. The temporal variation of fluid properties such as the force coefficient (CL, CD) and the Nusselt number (Nu) are shown at various heating levels. In addition, surface vorticity is examined in order to comprehend the flow dynamics near the surface of a heated square cylinder.