# I use the function knn.reg to evaluate the numerical solution at a fixed grid of points
#
if( !require('FNN') ){
    install.packages('FNN', dependencies=TRUE, repos='http://cran.rstudio.com/')
}
library(FNN)

source('runge_kutta.R')

t0 = 0.0
y0 = 1
t1 = 2.0

# Where do we want to evaluate the solution:
#
x_grid = c( seq( t0, 0.5, by=0.1 ), c( 1.0, 1.5, 2.0 ) )

dy.dt = function(t, y){ 1 - t + 4*y }

h=0.2
res = runge_kutta(dy.dt, h=h, start=t0, y0=y0, end=t1)
rk_yhat_h1 = knn.reg( as.matrix( res$xs ), y=res$ys, test=as.matrix( x_grid ), k=1, algorithm='brute' )

h=0.1
res = runge_kutta(dy.dt, h=h, start=t0, y0=y0, end=t1)
rk_yhat_h2 = knn.reg( as.matrix( res$xs ), y=res$ys, test=as.matrix( x_grid ), k=1, algorithm='brute' )

h=0.05
res = runge_kutta(dy.dt, h=h, start=t0, y0=y0, end=t1)
rk_yhat_h3 = knn.reg( as.matrix( res$xs ), y=res$ys, test=as.matrix( x_grid ), k=1, algorithm='brute' )

A = rbind( rk_yhat_h1$pred, rk_yhat_h2$pred, rk_yhat_h3$pred )
rownames(A) = NULL
R = cbind( c( 0.2, 0.1, 0.05 ), A )
colnames(R) = c( 'h', sprintf( 'y(%.2f)', x_grid ) )
print(R)