We investigate numerically a recent BGK-type model for a multi-component mixture of monatomic gases, undergoing a reversible bimolecular chemical reaction. The model replaces each collisional term of the Boltzmann equation with a relaxation term, thereby describing separately the effects of the mechanical processes and the chemical reaction. Additionally, the model exhibits consistency properties. The correct entropy production is ensured when auxiliary temperatures in the chemical contributions share a common value. We assume isotropic distributions and perform numerical simulations for the macroscopic fields to appraise how the dynamics push the mixture toward thermalization and chemical equilibrium. We show that the hypothesis on the equalization of fictitious species temperatures is justifiable to ensure the monotonicity of the classical H-Boltzmann functional. Simulations show that, when initial temperatures are far from equilibrium, the relaxation toward equilibrium occurs at a later stage and the classical H-Boltzmann functional is not monotone during the initial transient.