-1
def second(y, z, y0, z0, t0, tf, dt):
N = (tf-t0)/dt
  y.append(y0)
  for i in range(int(N)):
    t1 = t0 + dt
    y1 = y0 + z0*dt
    zk1 = z(t0, y0, z0)
    zk2 = z(t0+1/2*dt, y0+1/2*dt, z0+1/2*zk1*dt)
    z1 = z0 + zk2*dt
    y.append(y1)
    t0 = t1
    y0 = y1
    z0 = z1
  return y

I am getting something like this: [1, 1, -8.0, -26.0, -269.0]. Can you check if I understand how to apply the second order runge-kutta method?

  • Could you disentangle the situation a little bit. You are using `z` for the velocity and acceleration function, it is not really clear what method you are trying to use. Is it just the explicit midpoint method with reduction of the linear steps? – Lutz Lehmann Apr 05 '21 at 09:49
  • @LutzLehmann it is the second-order differential equation, which is the y''+4y^2=0. So searching the internet I found out that to solve the 2nd order ODE we need to the y'=z and z'= -4y^2 – Ulanbek Kurmaniyazov Apr 05 '21 at 09:53
  • I’m voting to close this question because this question is more likely to get an answer on [Code Review](https://codereview.stackexchange.com/) – j3ff Apr 05 '21 at 13:47

1 Answers1

0

You want to implement the explicit midpoint method

k1 = f(t,u)
k2 = f(t+dt/2, u+k1*dt/2)
u_next = u+k2*dt

for a second order equation y''=a(t,y,y'). To that end you set z=y' and u=[y,z]. In components you get

k1y = z
k1z = a(t,y,z)
k2y = z + k1z*dt/2
k2z = a(t+dt/2, y+k1y*dt/2, z+k1z*dt/2)

y_next = y + k2y*dt
z_next = z+k2z*dt

Then to use the special structure of the first order system of a second order equation, you eliminate k1y,k2y which gives

k1z = a(t,y,z)
k2z = a(t+dt/2, y+z*dt/2, z+k1z*dt/2)
y_next = y + z*dt + k1z*dt^2/2
z_next = y + k2z*dt

Now, compared to your code, there are several differences.

Lutz Lehmann
  • 25,219
  • 2
  • 22
  • 51