  High levels of deuterium fractionation of $\rm N_2H^+$ (i.e., $\rm
D_{frac}^{N_2H^+} \gtrsim 0.1$) are often observed in pre-stellar cores (PSCs)
and detection of $\rm N_2D^+$ is a promising method to identify elusive massive
PSCs. However, the physical and chemical conditions required to reach such high
levels of deuteration are still uncertain, as is the diagnostic utility of $\rm
N_2H^+$ and $\rm N_2D^+$ observations of PSCs. We perform 3D
magnetohydrodynamics simulations of a massive, turbulent, magnetised PSC,
coupled with a sophisticated deuteration astrochemical network. Although the
core has some magnetic/turbulent support, it collapses under gravity in about
one freefall time, which marks the end of the simulations. Our fiducial model
achieves relatively low $\rm D_{frac}^{N_2H^+} \sim 0.002$ during this time. We
then investigate effects of initial ortho-para ratio of $\rm H_2$ ($\rm
OPR^{H_2}$), temperature, cosmic ray (CR) ionization rate, CO and N-species
depletion factors and prior PSC chemical evolution. We find that high CR
ionization rates and high depletion factors allow the simulated $\rm
D_{frac}^{N_2H^+}$ and absolute abundances to match observational values within
one freefall time. For $\rm OPR^{H_2}$, while a lower initial value helps the
growth of $\rm D_{frac}^{N_2H^+}$, the spatial structure of deuteration is too
widespread compared to observed systems. For an example model with elevated CR
ionization rates and significant heavy element depletion, we then study the
kinematic and dynamic properties of the core as traced by its $\rm N_2D^+$
emission. The core, undergoing quite rapid collapse, exhibits disturbed
kinematics in its average velocity map. Still, because of magnetic support, the
core often appears kinematically sub-virial based on its $\rm N_2D^+$ velocity
dispersion.
