newton_optimize <- function(f, x0, tol = 1e-6, max_iter = 100) {
  iter <- 0
  while (iter < max_iter) {
    df <- grad(f, x0)  # Gradien
    ddf <- hessian(f, x0)  # Hessian
    dx <- solve(ddf, -df)  # Perbarui langkah
    x1 <- x0 + dx
    if (abs(f(x1) - f(x0)) < tol) return(list(x = x1, iter = iter))
    iter <- iter + 1
    x0 <- x1
  }
  return(NULL)
}
# Contoh penggunaan
x0 <- c(0, 0)
o1 <- newton_optimize(f, x0)

marquardt_optimize <- function(f, x0, tol = 1e-6, max_iter = 100, lambda = 1) {
  iter <- 0
  while (iter < max_iter) {
    df <- grad(f, x0)
    ddf <- hessian(f, x0)
    mu <- lambda * max(diag(ddf))  # Damping factor
    dx <- solve(ddf + mu * diag(length(x0)), -df)
    x1 <- x0 + dx
    if (abs(f(x1) - f(x0)) < tol) return(list(x = x1, iter = iter))
    iter <- iter + 1
    lambda <- ifelse(f(x1) < f(x0), lambda / 10, lambda * 10)
    x0 <- x1
  }
  return(NULL)
}
# Contoh penggunaan
o2 <- marquardt_optimize(f, x0)

dfp_optimize <- function(f, x0, B0, tol = 1e-6, max_iter = 100) {
  # Initialize variables
  iter <- 0
  n <- length(x0)
  x <- x0
  B <- B0
  
  # Loop until convergence or maximum number of iterations is reached
  while (iter < max_iter) {
    # Calculate the gradient using numerical differentiation
    g <- grad(f, x)
    
    # Calculate the search direction
    d <- -B %*% g
    
    # Perform line search to find optimal step size
    alpha <- optimize(function(alpha) f(x + alpha*d), interval = c(0, 1))$minimum
    x_new <- x + alpha * d
    
    # Calculate the new gradient using numerical differentiation
    g_new <- grad(f, x_new)
    
    # Calculate the update to the Hessian matrix
    s <- alpha * d
    y <- g_new - g
    B <- B + ((s %*% t(s)) / as.numeric(t(s) %*% y)) - ((B %*% y %*% t(y) %*% B) / as.numeric(t(y) %*% B %*% y))
    
    # Check for convergence
    if (abs(f(x_new) - f(x)) < tol) {
      return(list(x=x_new, iter=iter))
    }
    
    # Update the iteration counter and the guess
    iter <- iter + 1
    x <- x_new
  }
  
  # If the algorithm does not converge within the maximum number of iterations, return NULL
  return(NULL)
}

# DFP
B0<-diag(2)
o3 <- dfp_optimize(f, x0, B0)

# BFGS dengan optim
o4 <- optim(par = x0, fn = f, method = "BFGS", control = list(trace = TRUE))

gradient_descent <- function(obj_func, init_x, learning_rate = 0.01, maxIter = 1e5, tol = 1e-6) {
  x <- init_x
  gradx <- grad(obj_func, x)
  iter <- 0
  while (iter < maxIter && norm(as.matrix(gradx)) > tol) {
    gradx <- grad(obj_func, x)
    x <- x - learning_rate * gradx
    iter <- iter + 1
  }
  return(list(x = x, iter = iter))
}
# Contoh penggunaan
o5 <- gradient_descent(f, x0)