In relativistic jets formed by active galactic nuclei charged particles seem to be accelerated to high energies far away from the central black hole. Diffusive shock acceleration, also known as second order Fermi acceleration, is a possible mechanism explaining the observed particle spectra. In the test-particle regime the transport equations describing the plasma are Fokker-Planck type equations which makes it possible to calculate numerical solutions using stochastic differential equations. Advantages and disadvantages of this Monte-Carlo method are discussed and it is applied to the case of diffusive shock acceleration. The accuracy of the results depends on the choice of the numerical integrator used and a number of schemes are tested and compared. Basic integrators like the CauchyEuler scheme fail to predict the acceleration accurately in a scenario with steep shock gradients. It is shown that a semi-implicit second-order scheme produces more precise output in this case. Future applications in combination with magnetohydrodynamic simulations are discussed.