We introduce a simple ansatz for the wavefunction of a many-body system based
on coupled forward and backward-propagating semiclassical trajectories. This
method is primarily aimed at, but not limited to, treating nonequilibrium
dynamics in electron-phonon systems. The time-evolution of the system is
obtained from the Euler-Lagrange variational principle, and we show that this
ansatz yields Ehrenfest mean field theory in the limit that the forward and
backward trajectories are orthogonal, and in the limit that they coalesce. We
investigate accuracy and performance of this method by simulating electronic
relaxation in the spin-boson model and the Holstein model. Although this method
involves only pairs of semiclassical trajectories, it shows a substantial
improvement over mean field theory, capturing quantum coherence of nuclear
dynamics as well as electron-nuclear correlations. This improvement is
particularly evident in nonadiabatic systems, where the accuracy of this
coupled trajectory method extends well beyond the perturbative electron-phonon
coupling regime. This approach thus provides an attractive route forward to the
ab-initio description of relaxation processes, such as thermalization, in
condensed phase systems
