Pulsed plasma thrusters (PPTs) and other electric propulsion devices often operate in flow regimes that range from continuum to rarefied, posing significant challenges to conventional fluid or kinetic modeling approaches. In particular, accurately capturing magnetohydrodynamic (MHD) effects in plasma flows is crucial for improving thruster performance and design. The Discrete Unified Gas Kinetic Scheme (DUGKS), initially formulated for low-speed isothermal flows and later extended to compressible regimes , provides a promising framework for bridging these scales. By integrating MHD physics into DUGKS, we aim to develop a unified solver that can simulate plasma flows across the entire Knudsen number spectrum—an important capability for the electric propulsion community. Our work extends an existing OpenFOAM-based DUGKS implementation (dugksFoam) by adding the Lorentz force and magnetic induction equations. Specifically, we incorporate a force term into the Boltzmann–BGK equation whose second-order moment recovers the Maxwell stress tensor, thus allowing us to capture electromagnetic interactions in the flow. The induction equation for the magnetic field is then solved using DUGKS-computed macroscopic flow variables. This new solver, mhddugksFoam, includes models for PPT applications: an RLC circuit to represent the discharge process, a PTFE two-phase ablation model for propellant characterization, a Saha equation-based thermochemical plasma model, and simplified transport-coefficient formulations. These components collectively enable simulations of one-fluid, one-temperature plasma behavior under both continuum and rarefied conditions. Initial validation studies have been conducted on canonical MHD benchmarks—such as the Brio–Wu shock tube, MHD Rayleigh flow, and the Orszag–Tang vortex—demonstrating the solver’s accuracy in the continuum regime. Early tests on a simplified PPT configuration show the solver can capture key plasma discharge characteristics and plume expansion trends. We present preliminary PPT simulations, detailing performance metrics, ablation rates, and plume behavior across the Knudsen range. These findings are expected to provide new insights into PPT operation and pave the way for more robust, high-fidelity simulations in electric propulsion applications.