How to solve the stiff equation in this Restricted Three Body Problem numerically?

I've come across a stiff equation in solving the Circular Restricted Three Body Problem. [An object is moving considering the effect of the gravitational forces caused by two gravitational sources fixed in a 2D Space.]

The equations are these:

$x''=-\frac{GM_1 (x-x_1)}{\sqrt{(x-x_1)^2+y^2}^3}--\frac{GM_2 (x-x_2)}{\sqrt{(x-x_2)^2+y^2}^3}$

$y''=-\frac{GM_1 y}{\sqrt{(x-x_1)^2+y^2}^3}--\frac{GM_2 y}{\sqrt{(x-x_2)^2+y^2}^3}$

Neither Euler Method or Runge Kutta will work as the property near $(x_1, 0)$ or $(x_1, 0)$ is not good. The derivatives change too fast. The simulation can't be solved out right. The object is too easy to hit on the gravitational source.

How can I fix this?

Thank you!

-
What do you mean by "Neither Euler Method or Runge Kutta will work"? Are these your own implementations of these methods? Are you using a fixed time step size? I believe that both of these methods should work just fine even close to the singularities if you choose an adaptive time step size. –  Wolfgang Bangerth Jan 16 at 14:58
Have you tried an implicit method, such as Backward Euler? –  Paul May 25 at 1:52

migrated from physics.stackexchange.comJan 16 at 11:23

This question came from our site for active researchers, academics and students of physics.