To constrain cosmological models beyond ACDM, the development of the Euclid analysis pipeline requires simulations that capture the non-linear phenomenology of such models. We present an overview of numerical methods and N-body simulation codes developed to study the non-linear regime of structure formation in alternative dark energy and modified gravity theories. We review a variety of numerical techniques and approximations employed in cosmological N-body simulations to model the complex phenomenology of scenarios beyond ACDM. This includes discussions on solving non-linear field equations, accounting for fifth forces, and implementing screening mechanisms. Furthermore, we conduct a code comparison exercise to assess the reliability and convergence of different simulation codes across a range of models. Our analysis demonstrates a high degree of agreement among the outputs of different simulation codes, typically within 2% for the predicted modification of the matter power spectrum and within 4% for the predicted modification of the halo mass function, although some approximations degrade accuracy a bit further. This provides confidence in current numerical methods of modelling cosmic structure formation beyond ACDM. We highlight recent advances made in simulating the non-linear scales of structure formation, which are essential for leveraging the full scientific potential of the forthcoming observational data from the Euclid mission.