<p>As presence in cislunar space expands, new methods for characterizing the capabilities of low-thrust spacecraft in the regions beyond geosynchronous orbit are becoming increasingly important. Spacecraft reachability analysis, which refers to the determination of the set of states that a spacecraft can travel to or come from using its propulsion system over a given amount of time, is fundamental to understanding spacecraft trajectories and their neighboring phase-space regions in this new dynamical topography. Reachability analysis offers multifaceted utility to the space situational awareness community for tracking, detection, maneuver reconstruction, and evaluation of collision probabilities. However, the computation of reachable sets is a particularly formidable problem in astrodynamics and few methods exist to accurately and rapidly compute reachable sets for low-thrust spacecraft in strongly nonlinear regimes. Furthermore, computing reachable sets in the xGEO (beyond geosynchronous) regime presents unique challenges due to the highly sensitive and chaotic dynamics environment. In this paper, we present a unique technique for propagating the reachable sets of low-thrust spacecraft in the cislunar environment that leverages state-transition tensors and set-based computing techniques. In contrast to the prevalently used sample-based reachability techniques, which represent reachable sets as a collection of trajectories, this paper leverages a continuous set representation based on polynomial zonotopes, which offers additional benefits, including several set operations of interest and intersection checking.</p>