include "ConsoleWindow"

def tab 9

local fn dydx( x as double, y as double ) as double
end fn = x * sqr(y)

local fn exactY( x as long ) as double
end fn = ( x ^2 + 4 ) ^2 / 16

dim as long i
dim as double h, k1, k2, k3, k4, x, y, result

h = 0.1
y = 1
for i = 0 to 100
x = i * h
if x == int(x)
result = fn exactY( x )
print "y("; mid$( str$(x), 2, len(str$(x) )); ") = "; y, "Error = "; result - y
end if

k1 = h * fn dydx( x, y )
k2 = h * fn dydx( x + h / 2, y + k1 / 2 )
k3 = h * fn dydx( x + h / 2, y + k2 / 2 )
k4 = h * fn dydx( x + h,     y + k3     )

y = y + 1 / 6 * ( k1 + 2 * k2 + 2 * k3 + k4 )
next