Numerical transport models based on the advection-dispersion equation (ADE) are built on the assumption that sub-grid cell transport is Fickian such that dispersive spreading around the average velocity is symmetric and without significant tailing on the front edge of a solute plume. However, anomalous diffusion in the form of super-diffusion due to preferential pathways in an aquifer has been observed in field data, challenging the assumption of Fickian dispersion at the local scale. This study develops a fully Lagrangian method to simulate sub-grid super-diffusion in a multi-dimensional regional-scale transport. The underlying concept is based on previous observations that solutions to space-fractional ADEs, which can describe super-diffusive dispersion, can be obtained by transforming solutions of classical ADEs. The transformations are equivalent to randomizing particle travel time or relative velocity for each model time step. Here, the time randomizing procedure known as subordination is applied to flow field output from MODFLOW simulations. Numerical tests check the applicability of the novel method in mapping regional-scale super-diffusive transport conditioned on local properties of multi-dimensional heterogeneous media.