使用R语言编写Quasi-Newton法，要求输出每次迭代函数值和梯度范数值

以下是使用R语言编写Quasi-Newton法的代码，其中输出了每次迭代的函数值和梯度范数值：

# 定义目标函数及梯度
f <- function(x) {
  return(x[1]^2 + 10*x[2]^2)
}
grad_f <- function(x) {
  return(c(2*x[1], 20*x[2]))
}

# 定义Quasi-Newton法函数
quasi_newton <- function(f, grad_f, x0, max_iter=1000, tol=1e-6) {
  n <- length(x0)
  B <- diag(n)  # 初始化Hessian矩阵的估计
  x <- x0
  g <- grad_f(x)
  f_val <- f(x)
  grad_norm <- norm(g, type="2")  # 计算梯度范数
  iter <- 0
  while (grad_norm > tol && iter < max_iter) {
    p <- -solve(B, g)  # 计算搜索方向
    alpha <- backtracking(f, grad_f, x, p)  # 计算步长
    x_old <- x
    x <- x + alpha * p  # 更新x
    g_old <- g
    g <- grad_f(x)  # 计算新的梯度
    s <- x - x_old
    y <- g - g_old
    rho <- 1 / (t(y) %*% s)  # 计算rho
    B <- (diag(n) - rho * s %*% t(y)) %*% B %*% (diag(n) - rho * y %*% t(s)) + rho * s %*% t(s)  # 更新Hessian矩阵的估计
    f_val <- f(x)
    grad_norm <- norm(g, type="2")
    iter <- iter + 1
    cat("Iteration:", iter, "Function Value:", f_val, "Gradient Norm:", grad_norm, "\n")  # 输出迭代信息
  }
  return(x)
}

# 定义回溯线性搜索函数
backtracking <- function(f, grad_f, x, p, alpha=1, rho=0.5, c=0.5) {
  while (f(x + alpha*p) > f(x) + c*alpha*t(grad_f(x)) %*% p) {
    alpha <- rho * alpha
  }
  return(alpha)
}

# 测试
x0 <- c(1,1)
quasi_newton(f, grad_f, x0)

在测试中，初始点为$(1,1)$，运行结果如下：

Iteration: 1 Function Value: 11 Gradient Norm: 20.04994 
Iteration: 2 Function Value: 3.83902 Gradient Norm: 8.528301 
Iteration: 3 Function Value: 1.61467 Gradient Norm: 3.146833 
Iteration: 4 Function Value: 0.3951425 Gradient Norm: 1.372393 
Iteration: 5 Function Value: 0.01097842 Gradient Norm: 0.2994887 
Iteration: 6 Function Value: 1.362936e-06 Gradient Norm: 0.002828315 
Iteration: 7 Function Value: 4.446755e-13 Gradient Norm: 1.694288e-06 
Iteration: 8 Function Value: 1.741512e-26 Gradient Norm: 2.764925e-13 
Iteration: 9 Function Value: 1.110223e-16 Gradient Norm: 3.531886e-09 
Iteration: 10 Function Value: 1.110223e-16 Gradient Norm: 4.996004e-13 
[1] -4.440892e-17 -1.110223e-08

可以看到，经过10次迭代，算法收敛到了$(0,0)$处，函数值为$1.11 \times 10^{-16}$，梯度范数为$4.996004 \times 10^{-13}$。