newton_method <- function(f, df, x_0,
                          stop_criteria = .Machine$double.eps,
                          max_iter = 10^3,
                          left_b = -Inf, right_b = Inf){

  # criteria for starting point x_0
  if (left_b > x_0 | right_b < x_0){
    warning("your starting point is out of bound,
            please try point in the bounds")
    return(NA)
  }

  # variable initialization
  iters <- 0
  x_1 <- 0
  distance <- Inf

  # do while until
  # 1) x does not move
  # 2) exceed the maximun iteration
  # 3) the function value become close to zero
  while ((distance > stop_criteria) &
         (iters < max_iter) &
         (abs(f(x_0)) > stop_criteria)){
    # newton method
    x_1 <- x_0 - ( f(x_0) / df(x_0) )

    # out of bound case
    if (x_1 <= left_b) {
      x_1 <- left_b
    }
    if (right_b <= x_1){
      x_1 <- right_b
    }
    distance <- abs(x_0 - x_1)
    iters <- iters + 1
    x_0 <- x_1
  }

  # return the final point
  if (x_0 == left_b | x_0 == right_b){
    warning("\n", "Final point is on the boundary", "\n")
  }
  if (iters > max_iter){
    warning("\n", "Cannot find the root during the iteration.", "\n")
  }
  cat("Final value of function:", f(x_0), "\n")
}

g <- function(x){ - sin(x)}
out<-sapply(c(1,4,8), newton_method, f = cos, df = g)
sol_plot<-data.frame(x_star=out,y_star=cos(out))


library(ggplot2)
x_grid<-seq(0,10,length.out = 200)
df_plot<-data.frame(x=x_grid,y=cos(x_grid),zero=0)

ggplot(data=df_plot)+geom_line(aes(x=x,y=y),colour="blue")+
  geom_line(aes(x=x,y=zero),colour="black")+ theme_bw()+
  geom_point(data=sol_plot,aes(x=x_star,y=y_star),shape=4,size=5)