This work concerns with the numerical approximation for the stochastic Lotka–Volterra model originally studied by Mao et al. (Stoch Process Appl 97(1):95–110, 2002). The natures of the model including multi-dimension, super-linearity of both the drift and diffusion coefficients and the positivity of the solution make most of the existing numerical methods fail. In particular, the super-linearity of the diffusion coefficient results in the explosion of the 1st moment of the analytical solution at a finite time. This becomes one of our main technical challenges. As a result, the convergence framework is to be set up under the θ th moment with 0 < θ < 1 . The idea developed in this paper will not only be able to cope with the stochastic Lotka–Volterra model but also work for a large class of multi-dimensional super-linear SDE models. [ABSTRACT FROM AUTHOR]