0

I was using the 4th order Runge Kutta method to solve differential equations of the duffing oscillator with NumPy arrays, but I had received an Error.

RuntimeWarning: overflow encountered in double_scalars

Does anyone know the possible sources of error that may cause the overflow error and ways in which I would be able to try to resolve it?

t = np.linspace(0,1,steps)
    # start at t=0 
    h = t[1] - t[0] 
    xn = np.zeros(steps)
    xn[0] = 1
    vn = np.zeros(steps)
    vn[0] = 1
    for n in range(1,steps):
        k1 = vn[n-1]
        l1 = g * np.cos(w * t[n-1]) - d * vn[n-1] - a * xn[n-1] - b * xn[n-1] **3
        k2 = vn[n-1] + h*k1/2
        l2 = g * np.cos(w * (t[n-1] + h/2)) - d * (vn[n-1] + h*k1/2) - a * (xn[n-1] + h*l1/2) - b * (xn[n-1] + h*l1/2)**3

        k3 = vn[n-1] + h*k2/2
        l3 = g * np.cos(w * (t[n-1] + h/2)) - d * (vn[n-1] + h*k2/2) - a * (xn[n-1] + h*l2/2) - b * (xn[n-1] + h*l2/2)**3

        k4 = vn[n-1]*k3
        l4 = g * np.cos(w * (t[n-1] + h)) - d * (vn[n-1] + h*k3) - a * (xn[n-1] + h*l3) - b * (xn[n-1] + h*l3)**3
        vn[n] = vn[n-1] + h/6 * (k1 + 2* k2 + 2 * k3+ k4)
        xn[n] = xn[n-1] + h/6 * (l1 + 2* l2 + 2* l3 + l4)

Here is my code for reference should anyone needs it. vn represents v_numerical A link to Wikipedia for reference of the formula: https://en.wikipedia.org/wiki/Runge%E2%80%93Kutta_methods I don't think it is very important to understand the formula, some possible solutions that I could try resolve the problem would be helpful.

James
  • 27
  • 5

1 Answers1

0

The k are the updates for x, the l the updates for v, as v=dx/dt. You are changing the ODE into a different system of two half-coupled first-order equations.

Also note the typo or deletion error in k4.

Lutz Lehmann
  • 25,219
  • 2
  • 22
  • 51
  • thank you for the help, I have made the changes to my code, however, the overflow error still persists, do you know how to resolve it? I was planning to do an Poincare map, but keep running into issues with overflow – James Apr 27 '21 at 14:38
  • After some debugging, I realised that some of my x_n terms reach 10^100, which was what that had caused the error – James Apr 27 '21 at 14:45