The present study focuses on the numerical analysis of turbulent natural convection within a square enclosure with differentially heated vertical walls, and this is achieved by employing a Large Eddy Simulation (LES) approach coupled with the Lattice Boltzmann Method (LBM). The lattice model used for this study is based on the D2Q9 lattice for both hydrodynamic and thermal fields, and buoyancy effects are modeled using the Boussinesq approximation. The results for this study are presented using contours of stream function, isotherms, and mean Nusselt number. At moderate values of Rayleigh number, i.e., for Ra ≤ 108, the thermal field within the core of the square cavity is found to be almost uniform, and this shows that heat transfer is due to conduction. However, for Ra = 109, significant distortions and oscillations within the isothermal contours show that turbulent convection is established within the square cavity.