Abstract
In this article a method of numerical solution of the Schrödinger equation is proposed. The approach corrects the Ehrenfest approximation by using several trajectories/configurations with their amplitudes coupled within and across configurations, thus making the method formally exact. Accurate results are obtained for the spin-boson model with up to 2000 bath modes treated on fully quantum level without approximations.