In this paper, a τ-leap algorithm for stochastic simulation of chemical reaction was adopted to 13C labeling kinetics simulation of enzymatic reaction. The 13C labeling kinetics simulation of enzymatic reaction at the single-cell level was completed, and the kinetics simulations of enzymatic reaction in the multiple single-cells were carried out with the parameters from Escherichia coli cell. The numbers of isotopic molecules and the proportions of the numbers of isotopic molecules were calculated within a small time interval. The intergroup variances for multiple cell groups were calculated and compared. The simulation results show that the numbers of isotopic molecules and the proportions of the numbers of isotopic molecules of the metabolites are heterogeneous in the stochastic domain of the Michaelis–Menten mechanism at the single-cell level. With the same number of iterations, this heterogeneity gradually decreases with the increase in the initial numbers of isotopic molecules of the metabolites at the single-cell level. In addition, the intergroup heterogeneity of the numbers of isotopic molecules and the proportions of the numbers of isotopic molecules shows a decreasing trend with the increase in the number of simulated cells within the group. These results demonstrate the necessity and possibility for a complete 13C flux analysis at single-cell in future.