This paper presents a Particle-In-Cell code designed for the simulation of water-fuelled Hall Effect Thrusters, including two different propellants: water vapour and oxygen (the latter being intended for water electrolysis propulsion with oxygen supplied through the anode and hydrogen through the cathode). The reactive model is structured in two stages, encompassing not only the initial reactions of water vapour and oxygen molecules (1st stage) but also the reactions of the diatomic and monoatomic products resulting from them (2nd stage). Specifically, the model accounts for 45 reactions in the case of water vapour and 18 reactions in the case of oxygen, including the most relevant excitation events for each species. The Particle-In-Cell code uses a combination of a 0-dimensional model with a 2-dimensional model. The 0-dimensional model provides initial neutral and electron densities, as well as the most significant reactions, to facilitate the convergence of the 2-dimensional model without the computational burden of starting a simulation from scratch. The 0-dimensional model reveals that the reactions considered within the 2nd stage are crucial for the plasma species composition of the discharge. The oxygen plasmas consist mainly of O+ and O2+ ions in a similar proportion, while double and negative ions do not play a significant role. Neutrals (O2 and O) also show similar distributions, depending on the thruster’s operating conditions. Water vapour plasmas are dominated by OH+ , H2O+ , H+ , and O+ ions, with other species such as H2+ , O++ , and negative ions being negligible. The neutral population is predominantly composed of monoatomic H particles. The O-dimensional model also demonstrates that all ionisation fractions of the species follow an exponentially increasing trend with the electron temperature. Finally, the 2-dimensional model provides additional insight into the plasma evolution, electron temperature, power losses coming from the reactive model and equilibrium points of the system.