Our quantum molecular dynamics simulations are performed by propagating vibronic wavepackets and classical-like trajectories. For the former, we will focus on the use of different ideas to define the time-independent grid, i.e., direct product versus Smolyak grids, in order to be able to include efficiently at a reasonable computational cost the secondary modes representing the environment. The calculations are performed with the ElVibRot code . For the latter, we will test and compare Tully surface hopping (TSH), an independent-trajectory method, to its coupled-trajectory formulation derived from the exact factorization of the electron-nuclear wavefunction , namely the coupled-trajectory Tully surface hopping (CT-TSH) algorithm. In particular, we will present an in-depth analysis of the computational cost of independent-trajectory and coupled-trajectory schemes as the number of degrees of freedom increases. The trajectory-based methods employed in this work are implemented in the G-CTMQC code . Both ElVibRot and G-CTMQC are interfaced with QuantumModelLib , a library of analytical model potentials which includes the model system-bath Hamiltonian used in this work.