  In this study, we investigate the numerical stability of the covariant
Baumgarte--Shapiro--Shibata--Nakamura (cBSSN) formulation against the
Friedmann--Lema\^itre--Robertson--Walker spacetime. To evaluate the numerical
stability, we calculate the constraint amplification factor by the eigenvalue
analysis of the evolution of the constraint. We propose a modification to the
time evolution equations of the cBSSN formulation for a higher numerical
stability. Furthermore, we perform numerical simulations using the modified
formulation to confirm its improved stability.
