where the |Φ k and E k are, respectively, the eigenvectors and the eigenvalues of the Hamiltonian in the Krylov subspace. The first vector of this subspace is the wavepacket at time t and other vectors are added until convergence (the contribution of the last |Φ k is smaller than 10 -12 ). Furthermore, to avoid numerical instabilities, a vector orthogonalization is enforced via a Gram-Schmidt procedure at each iteration k.