We present 2D magnetohydrodynamics numerical simulations of tearing-unstable current sheets coupled to a population of non-thermal test particles, in order to address the problem of numerical convergence with respect to grid resolution, numerical method, and physical resistivity. Numerical simulations are performed with the PLUTO code for astrophysical fluid dynamics through different combinations of Riemann solvers, reconstruction methods, and grid resolutions at various Lundquist numbers. The constrained transport method is employed to control the divergence-free condition of magnetic field. Our results indicate that the reconnection rate of the background tearing-unstable plasma converges only for finite values of the Lundquist number and for sufficiently large grid resolutions. In general, it is found that (for a second-order scheme) the minimum threshold for numerical convergence during the linear phases requires the number of computational zones covering the initial current sheet width to scale roughly as the square of the Lundquist number defined on the current sheet width. On the other hand, the process of particle acceleration is found to be nearly independent of the underlying numerical details inasmuch as the system becomes tearing-unstable and enters in its non-linear stages. In the limit of large Lundquist number, the ensuing power-law index quickly converge to p ≈ 1.7, consistently with the fast reconnection regime.