Ambipolar plasma thrusters are an appealing technology due to multiple system-related advantages, including propellant flexibility and the absence of electrodes or neutralizer. Understanding the plasma generation and acceleration mechanisms is key to improving the performance and capabilities of these thrusters. However, the source and plume regions inside are often simulated separately, and no self-consistent strategy exists which can couple these different simulations together. This paper introduces the MUlti-regime Plasma Equilibrium Transport Solver (MUPETS), a self-consistent coupled model integrating a fluid solver for the plasma dynamics in the source, which are collision-driven, with a kinetic Particle-In-Cell (PIC) code for the plasma dynamics in the magnetic nozzle, which involve expansion across a diverging magnetic field. The methodology begins by solving the plasma source with the classical Bohm condition at the thruster’s throat. The resulting plasma profiles (density, temperature, speed) are input into the PIC code for the magnetic nozzle. The PIC code calculates the plasma plume expansion and determines the electric field at the thruster’s throat. This electric field is then used as a boundary condition in the fluid code, where it replaces the Bohm assumption, and the fluid simulation is repeated. This iterative process continues until convergence. In comparing the MUPETS results with those for an experimental thruster, the plasma densities at the thruster’s throat differed by less than 2–5% between the fluid and PIC regions. The thrust predictions agreed with the experimental trend, and were kept well within the measurement’s uncertainty band. These results validate the effectiveness of the coupling strategy for enhancing plasma thruster simulation accuracy.