Spherical-harmonics expansion method is proposed to solve the quantum time-evolution equations for density matrices numerically in momentum space. This method facilitates efficient real-time simulations of quantum electron dynamics in continuum states through extension of our previous formalism developed for discrete states. Numerical accuracy and efficiency are demonstrated through two examples: (i) multiphoton ionization of a one-electron atom in intense laser field, and (ii) real-time dynamics of plasma oscillations in electron liquids. In case (i), coupled dynamics of density matrices for bound and continuum states reveals an enhancement of multiphoton ionization over the Keldysh approximation through resonant intermediate states. (C) 2015 Wiley Periodicals, Inc.