Metallic nanoparticles support surface plasmon resonances after illumination and such collective excitations (involving many high-lying electronic states) can lead to many interesting dynamical phenomena, such as plasmon-mediated reduction and dissociation and plasmon-induced charge transfer. Near a metal surface, the interaction of a molecular system with a continium of states in a metal can also lead to chemicurrents, unusual vibrational relaxation, and other profoundly non-adiabatic effects. In order to model such essential (but complicated, many-body) phenomena, accurate and efficient tools for simulating coupled electron-nuclei dynamics are sorely needed, especially methods that can accommodate a large number of electronic excited states and many nuclear degrees of freedom. Simulating dozens or hundreds of electronic states with ab initio potentials remains a key challenge for physical chemists. Now, if we seek a molecular description of dynamics near a metal surface, practical considerations dictate that we must employ a semiclassical approach, one whereby we compute the electronic wavefunction using high-level quantum theory but the nuclear wavepacket is modeled by an ensemble of classical trajectories. Over the past two decades, beyond Ehrenfest dynamics, a few semiclassical methods have been proposed to treat dynamics near a metal surface based on Tully's fewest-switches surface hopping (FSSH) algorithm, including independent electron surface hopping (IESH), surface hopping with electronic friction (FSSH-ER), and a broadened classical master equation (BCME). Each of these methods requires propagating trajectories on one surface, with hops between potential energy surfaces to account for electronic relaxation. And in principle, two of these methods (IESH and FSSH-ER) can be combined with ab initio electronic structure calculations so as to simulate an arbitrary number of electronic states in the presence of many nuclear degrees of freedom.