  The global coupling of few-level quantum systems ("spins") to a discrete set
of bosonic modes is a key ingredient for many applications in quantum science,
including large-scale entanglement generation, quantum simulation of the
dynamics of long-range interacting spin models, and hybrid platforms for force
and spin sensing. We present a general numerical framework for treating the
out-of-equilibrium dynamics of such models based on matrix product states. Our
approach applies for generic spin-boson systems: it treats any spatial and
operator dependence of the two-body spin-boson coupling and places no
restrictions on relative energy scales. We show that the full counting
statistics of collective spin measurements and infidelity of quantum simulation
due to spin-boson entanglement, both of which are difficult to obtain by other
techniques, are readily calculable in our approach. We benchmark our method
using a recently developed exact solution for a particular spin-boson coupling
relevant to trapped ion quantum simulators. Finally, we show how decoherence
can be incorporated within our framework using the method of quantum
trajectories, and study the dynamics of an open-system spin-boson model with
spatially non-uniform spin-boson coupling relevant for trapped atomic ion
crystals in the presence of molecular ion impurities.
