R精度/输出中的牛顿法

Question

所以，我应该编写代码来执行牛顿法，以指定精度（公差）计算任意数的平方根。

这是我的代码：

MySqrt <- function(x, eps = 1e-6, itmax = 100, verbose = TRUE) {
  GUESS <- 11
  myvector <- integer(0)
  i <- 1
  if (x < 0) {
    stop("Square root of negative value")
  }
  else {
    myvector[i] <- GUESS

    while (i <= itmax) {
      GUESS <- (GUESS + (x/GUESS)) * 0.5
      myvector[i+1] <- GUESS

      if (abs(GUESS-myvector[i]) < eps) {
        break()
      } 

      if (verbose) {
        cat("Iteration: ", formatC(i, width = 1), formatC(GUESS, digits = 10, width = 12),     "\n")
      }

      i <- i + 1

    }
  }    
  myvector[i]
}

eps 是公差。当我使用该函数计算 21 的平方根时，我得到了以下输出：

> MySqrt(21, eps = 1e-1, verbose = TRUE)
Iteration:  1  6.454545455 
Iteration:  2  4.854033291 
Iteration:  3   4.59016621 

但是，我不确定该函数是否会停止执行迭代。有人可以验证我的代码是否正确吗？将不胜感激！

Answer 1

您的代码几乎是正确的。它正在迭代正确的次数。唯一的错误是您不会在 break 语句之后增加 i，因此您不会返回最新的近似值。相反，您将返回前一个。

为了验证它是否在正确的时间停止，您可以将跟踪线向上移动到中断上方。您还可以将GUESS-myvector[i] 添加到跟踪中，这样只要差异变得足够小，您就可以看到它停止。如果您这样做并运行该函数，那么它在正确的时间停止的事实以及它返回错误值的事实将是显而易见的：

> MySqrt(21,eps=1e-1)
Iteration:  1 6.454545 -4.545455 
Iteration:  2 4.854033 -1.600512 
Iteration:  3 4.590166 -0.2638671 
Iteration:  4 4.582582 -0.007584239 
[1] 4.590166

虽然您的代码（几乎）是正确的，但它并不是以非常好的 R 风格编写的。例如，除非您想返回整个估计向量，否则没有理由需要将它们全部保留。此外，与其使用 while 循环，不如使用 for 循环更有意义。这是您功能的一种可能的改进版本：

MySqrt <- function(x, eps = 1e-6, itmax = 100, verbose = TRUE) {
  GUESS <- 11
  if (x < 0) {
    stop("Square root of negative value")
  }
  for(i in 1:itmax){
      nextGUESS <- (GUESS + (x/GUESS)) * 0.5
      if (verbose)
        cat("Iteration: ", i, nextGUESS, nextGUESS-GUESS, "\n")

      if (abs(GUESS-nextGUESS) < eps) 
        break

      GUESS<- nextGUESS
    }
  nextGUESS
}