# metode numerik
f1 <- function(x,y){y/(2*x+1)}

# metode analitik
f2 <- function(x){sqrt(2*x+1)}
x0 <- 0
y0 <- 1
x <- x0
y <- y0

for(i in 1:100){
  y0 <- f2(x0+0.05)
  x0 <- x0+0.05
  x <- c(x, x0)
  y <- c(y, y0)
}
true <- data.frame(x=x, y=y)

1 Metode Heun

heun <- function(f, x0, y0, h, n, iter=1){
  x <- x0
  y <- y0
  
  for(i in 1:n){
    ypred0 <- f(x0,y0)
    ypred1 <- y0 + h*ypred0
    ypred2 <- f(x0+h,ypred1)
    ykor <- y0 + h*(ypred0+ypred2)/2
    if(iter!=1){
      for(i in 1:iter){
        ykor <- y0 + h*(ypred0+f(x0+h,ykor))/2
      }
    }
    y0 <- ykor
    x0 <- x0 + h
    x <- c(x, x0)
    y <- c(y, y0)
  }
  
  return(data.frame(x=x,y=y))
}
num <- heun(f1, x0=0, y0=1, h=0.05, n=100)

2 Metode Runge-Kutta Orde 4

rk4 <- function(f, x0, y0, h, n){
  x <- x0
  y <- y0
  
  for(i in 1:n){
    k1 <- f(x0,y0)
    k2 <- f(x0+0.5*h,y0+0.5*k1*h)
    k3 <- f(x0+0.5*h,y0+0.5*k2*h)
    k4 <- f(x0+h,y0+k3*h)
    y0 <- y0 + (1/6)*(k1+2*k2+2*k3+k4)*h
    x0 <- x0 + h
    x <- c(x, x0)
    y <- c(y, y0)
  }
  
  return(data.frame(x=x,y=y))
}
num <- rk4(f1, x0=0, y0=1, h=0.05, n=100)