The advancement of nonadiabatic dynamics simulations relies on the parallel development of their two essential components: the way the nuclei are propagated and the electronic structure method used. These developments must also face the well-known and unavoidable accuracy vs. cost trade-off that leads to choosing one method over another. Regarding the nuclear propagation among different potential energy surfaces, the space of methods is spanned by trajectory-based mixed quantum-classical methods, such as trajectory surface hopping, on the one side3, and fully quantum wavefunction propagation, such as multiconfigurational time-dependent Hartree, on the other extreme of the spectrum. A family of methods that can be placed in the middle are Gaussian-based approaches. They rely on the expansion of the nuclear wavefunction in a basis of Gaussian functions leading to a trajectory-like approach while retaining a nuclear wavefunction in the formalism. The propagation of the trajectory basis functions (TBFs), whether they follow quantum, mixed quantum-classical, or classical trajectories, differentiates the three most popular Gaussian-based methods, variational multiconfigurational Gaussian, multiconfigurational Ehrenfest and full and Ab initio multiple spawning, 9,10 respectively. Despite their different formalisms, in the limit of a sufficiently large Gaussian basis and an exact evaluation of the matrix elements, all these methods will converge to an exact solution. Ab initio multiple spawning (AIMS) has risen to a popular choice of nonadiabatic dynamics methods for the description of the photochemistry of medium-sized molecules. It has successfully been used to model the excited state dynamics of important biomolecules, and to directly simulate timeresolved spectra. In particular, thanks to its extensions to the description of intersystem crossing, external fields, tunnelling and complex environments, AIMS has become a method of wide applicability that offers a good balance between accuracy and cost. The original most widely used implementation of AIMS, FMS90, has been integrated within the Molpro program package as well as TeraChem. Additionally, interfaces have been created for Columbus, GAMESS, and Mopac. There has been much recent interest in the development of broadly available software packages for nonadiabatic molecules dynamics simulations, following the example of the electronic structure community. In this spirit, the Levine group has developed PySpawn, an open-source software package for AIMS, which is designed to be highly extensible for integration with various electronic structure codes. Developed in Python, this program offers an intuitive interface for both software developers and end-users. It includes an interactive analysis module designed to facilitate the broader adoption of AIMS. Moreover, PySpawn employs a novel task-based implementation of the AIMS algorithm and utilizes modern structured output formats, enabling efficient deployment on shared high-performance computing resources.