I am currently trying to get the Runge Kutta 4 Integration to work, but it outputs the following:
Runge Kutta value: inf
While it is supposed to give a value in the range of:
Runge Kutta value: 8.476271005220534e+16
I used the following code, but can't seem to output the right approximation.
endtime = 5
h = 0.01
def f(x, z):
return x*x*z
t = 0
y = 1
while t < endtime:
k1 = f(t, y)
k2 = f(t+(h/2), y+(k1/2))
k3 = f(t+(h/2), y+(k2/2))
k4 = f(t+h, y+k3)
y = y + h*(k1 + 2*k2 + 2*k3 + k4)/6
t = t + h
print("Runge Kutta value: " + str(y))
Is there anyone who knows where I made the mistake
your formulas are wrong because most of the time you're forgetting to multiply by h
and coefficients too:
(source: Wikipedia)
so:
endtime = 5
h = 0.01
def f(x, z):
return x*x*z
y = 1
count = 0
t = 0
for count in range(int(endtime/h)):
t = h*count
k1 = f(t, y)
k2 = f(t+(h/2), y+h*(k1/2))
k3 = f(t+(h/2), y+h*(k2/2))
k4 = f(t+h, y+h*k3)
y += h*(k1 + 2*k2 + 2*k3 + k4)/6
print("Runge Kutta value: " + str(y))
Also avoid floating point accumulation by computing the value of t
each time, instead of adding the fixed step, and trade the while
loop for a for
loop.
not sure it's perfect but I get a finite result now :)
Runge Kutta value: 1.245858162131811e+18