The simulation of nonlinear self-focusing phenomenon using ray-tracing method can macroscopically provide an intuitive picture of the propagation of light in a self-focusing material, without adopting paraxial approximation or self-similar hypothesis. In this paper, propagation of light is sampled by discrete slices along a certain direction. Thus nonlinear propagation is turned into the combination of optical modulation of the refractive index on separate slices and linear propagation between each two adjacent slices. On each slice, after calculating the flux, we use a novel algorithm to suppress the quantized errors. For the linear propagating process, Adams method is adopted to solve the ray equations, which solve the problem that the widely used Runge-Kutta method cannot be used in simulation of light in nonlinear materials. The simulation results reveal that there are several foci along the propagating axis and the location of the first focus becomes closer to the incident plane as the power of light goes up. Furthermore, because the program traces real rays, it is possible to reach the non-paraxial region and reveal the phenomenon of ring-structure flux distributions caused by self-focusing. This is significant for the safety of high-power laser systems. Some commercial optical design and simulation software are also based on ray-tracing methods. Thus the systems including both nonlinear and linear materials are possible to simulate, which can guide people to set up the corresponding experimental systems.