  A self-consistent model is presented for performing steady-state fully
kinetic Particle-in-Cell simulations of magnetised plasma plumes. An
energy-based electron reflection prevents the numerical pump instability
associated with a typical open-outflow boundary, and is shown to be
sufficiently general that both the plume kinetics and plasma potential
demonstrate domain independence (within 4%). This is upheld by non-stationary
Robin-type boundary conditions on the Poisson's equation, coupled to a
capacitive circuit that allows physical evolution of the downstream potential
drop in the transient. The method has been validated against experiments,
providing results that fall within the uncertainty of measurements. Simulations
are then carried out to study collisional xenon discharges into axisymmetric
diverging magnetic nozzles. Particular discussion is given to the
identification of a potential well arising from charge separation at the edge
of the plume, the role of ion-neutral charge exchange, and a three-region
piecewise polytropic cooling regime for electrons. The polytropic index is
shown to depend on the degree of magnetisation. Specifically, in the region
near the thruster outlet, the plume is weakly-magnetised due to the cross-field
diffusion of electron-heavy particle collisions. Downstream, a
strongly-magnetised region of near-isothermal expansion occurs. Finally, in the
detached region, the polytropic index tends to that of a more adiabatic
unmagnetised case. With an increasing magnetic nozzle field strength, an
inferior limit is found to the average polytropic index of
$\bar{\gamma}_e\sim1.16$.
