An algorithm for time propagation of the time-dependent Kohn-Sham equations is presented. The algorithm is based on dividing the Hamiltonian into small time steps and assuming that it is constant over these steps. This allows for the time-propagating Kohn-Sham wave function to be expanded in the instantaneous eigenstates of the Hamiltonian. The method is particularly efficient for basis sets which allow for a full diagonalization of the Hamiltonian matrix. One such basis is the linearized augmented plane waves. In this case we find it is sufficient to perform the evolution as a second-variational step alone, so long as sufficient number of first variational states are used. The algorithm is tested not just for nonmagnetic but also for fully non-collinear magnetic systems. We show that even for delicate properties, like the magnetization density, fairly large time-step sizes can be used demonstrating the stability and efficiency of the algorithm.