NormaJ: An efficient numerical code for geodynamic modeling of tectonic processes

Geodynamic numerical modelling is a common approach for investigating the mechanical evolution of large-scale tectonic features such as accretionary wedges, foreland fold-and-thrust belts, and continental rifting. Several numerical techniques have been applied in the field of geodynamics, such as the finite difference method, the finite volume method, the finite element method, and the spectral method. The long-term mechanical behaviours of large-scale tectonic processes are usually represented by the implementation of visco-elasto-plastic constitutive models to mimic the rheological mechanisms of geomaterials. However, such often nonlinear rock behaviours result in computational challenges in geodynamic numerical modeling, especially when considering the geological time scale. Here, we aim at developing and presenting an efficient numerical code (NormaJ) that acknowledges visco-elasto-plastic rheology to study large-scale tectonic processes. We apply the finite difference method with a fully staggered Eulerian grid strategy to solve nonlinear partial differential equations, integrating the marker-in-cell technique to track deformation. The code was programmed in Julia, which provides a better computational performance and data organization than Matlab and Python, while it is more user-friendly than C and Fortran. Thus, Julia is particularly suitable for large-scale geodynamics modeling. We have benchmarked our new code NormaJ with the existing fully-vectorized Matlab code Norma (https://github.com/Norma-VEP) using a baseline example, demonstrating an increase in time efficiency by a factor of about two. Such improvement indicates a good prospect for future applications and development. NormaJ will be openly available, providing the geodynamics community an easy-to-run and accessible tool to solve tectonic evolution problems.