Runaway electrons (REs), generated during plasma disruptions in tokamaks, pose significant challenges due to the risk of causing damage to the first wall of a device. Understanding the interaction between REs and magnetohydrodynamic (MHD) instabilities is crucial for predicting a safe operation of large future tokamak devices in which RE generation will be drastically enhanced due to the high plasma current. In this work, we introduce a hybrid fluid–kinetic model within the three-dimensional nonlinear MHD code JOREK (Hoelzl et al. 2021 Nucl. Fusion, vol. 61, 065001; 2024 Nucl. Fusion, vol. 64, 112016), treating REs kinetically using a relativistic guiding-centre approach, while describing the background plasma by ansatz-based reduced MHD equations. At first, comprehensive benchmark studies are conducted regarding the two-dimensional equilibrium force balance with
$J_{total}= J_{RE}$, and the linear stability of three-dimensional tearing modes (TMs), verifying the accuracy of the model against analytical predictions and other numerical methods, e.g. the full-orbit approach in JOREK and the fluid model in M3D-C1. These benchmark studies build a solid foundation for applying our model to more complex nonlinear scenarios. In this respect, we confirm that the nonlinear saturation of TMs is significantly influenced by the presence of REs. Previous analytical studies (Helander et al. 2007 Phys. Plasmas vol. 14, 122102) suggest that in the case of small
$\varDelta ^\prime$, the saturation width of the magnetic island driven by REs is roughly 1.5 times larger than in the otherwise identical Ohmic current scenario. Our simulations are quantitatively in line with this prediction. Moreover, REs alter the energy evolution within the magnetic reconnection process and decouple the bulk plasma and magnetic fields. In summary, RE-driven magnetic reconnection leads to larger magnetic islands but weaker reconnection flows.