A new algorithm for numerical simulation of Langevin equations
  Formulated is a new systematic method for obtaining higher order corrections in numerical simulation of stochastic differential equations (SDEs), i.e., Langevin equations. Random walk step algorithms within a given order of finite $\Delta t$, are obtained so as to reproduce within that order a corresponding transition density of the Fokker-Planck equations, in the weak Taylor approximation scheme. A great advantage of our method is its straightforwardness such that direct perturbative calculations produce the algorithm as an end result, so that the procedure is tractable by computer. Examples in general form for curved space cases as well as flat space cases are given in some order of approximations. Simulations are performed for specific examples of U(1) system and SU(2) systems, respectively. 
H. Nakajima and S. Furui