fcn yp(t,y) { t * y.sqrt() }
fcn exact(t){ u:=0.25*t*t + 1.0; u*u }

fcn rk4_step([(y,t)],h){
   k1:=h * yp(t,y);
   k2:=h * yp(t + 0.5*h, y + 0.5*k1);
   k3:=h * yp(t + 0.5*h, y + 0.5*k2);
   k4:=h * yp(t + h, y + k3);
   T(y + (k1+k4)/6.0 + (k2+k3)/3.0, t + h);
}

fcn loop(h,n,[(y,t)]){
   if(n % 10 == 1)
      print("t = %f,\ty = %f,\terr = %g\n".fmt(t,y,(y - exact(t)).abs()));
   if(n < 102) return(loop(h,(n+1),rk4_step(T(y,t),h))) //tail recursion
}