  We investigate the thermalization of high-energy particles injected from the
perturbative decay of inflaton during the pre-thermal phase of reheating in
detail. In general, thermalization takes a relatively long time in a
low-temperature plasma; therefore, the instantaneous thermalization
approximation is not justified, even for the reheating of the Standard Model
(SM) sector. We consider a pure Yang-Mills (YM) theory as an approximation of
the SM sector or a possible dark sector, considering the
Landau-Pomeranchuk-Migdal effect, a quantum interference effect in a finite
temperature plasma. We perform the first numerical calculation to solve the
time evolution of the system, including the redshift due to the expansion of
the Universe, and show the details of the temperature evolution near the
maximum and the behavior of the quasi-attractors at later times. The maximal
temperature $T_\text{max}$ and time scale $t_\text{max}$ are determined
quantitatively, such as $T_\text{max} \simeq 0.05 \times (\Gamma_I
M_\text{Pl}^2/m_I^3)^{2/5} m_I$ and $t_\text{max} \simeq 2 \times 10^3 \times
(\Gamma_I M_\text{Pl}^2/m_I^3)^{-3/5} m_I^{-1}$ in the SM-like system, where
$m_I$ and $\Gamma_I$ are the mass and decay rate of inflaton. We also provide a
similar formula for pure $\operatorname*{SU}(N)$ and $\operatorname*{SO}(N)$ YM
theories for general values of $N$ and coupling constant $\alpha$, including
$T_\text{max} \propto \alpha^{4/5}$ and $t_\text{max} \propto N^{-2}
\alpha^{-16/5}$ behaviors and their numerical coefficients. The thermalization
occurs in a finite time scale, resulting in a lower maximal temperature of the
Universe after inflation than that under the instantaneous thermalization
approximation.
