PySCES is a highly parallelized PYTHON program to perform nonadiabatic molecular dynamics simulation for large–scale open quantum systems through the Linearized Semiclassical Initial Value Representation (LSC–IVR). Every dynamical variable in a simulation is propagated in phase space as prescribed by LSC-IVR. The electronic structure variables necessary to integrate the equations of motion are provided "on-the-fly" through the interface established between the code and external electronic structure softwares, specifically TeraChem and GAMESS. The current implementation is capable of computing an electronic population correlation function through three different population estimators. For more details, we refer readers to earlier studies: LSC–IVR• W. H. Miller, J. Phys. Chem. A 105, 2942– 2955 (2001)• Q. Shi and E. Geva, J. Chem. Phys. 118, 8173 (2003) On-the-fly implementation of LSC–IVR and different population estimators• K. Miyazaki and N. Ananth, J. Chem. Phys., 159, 124110 (2023)