This work presents a Streamline-Upwind Petrov–Galerkin (SUPG) framework with finite elements discretization for the prediction of non-ionized hypersonic flows in thermal non-equilibrium. The formulation is enhanced with a residual-based discontinuity-capturing (DC) operator. The numerical framework solves the set of Navier–Stokes equations for the compressible reacting flows with an additional equation for the vibrational–electronic energy conservation. The set of equations, which involves one continuity equation for each chemical species considered, is solved in the pressure-primitive variables. The numerical framework capabilities are assessed through several benchmark cases, including validation of the chemical model and its coupling with the two-temperature model in the context of critical flow features such as shock–shock and shock wave–boundary layer interactions. The accuracy of the results is assessed through the comparison with the numerical and experimental data available in the literature for all the test cases presented. The numerical results demonstrate the suitability of the formulation in predicting non-ionized reacting hypersonic flows in thermal non-equilibrium.