在 R 中使用 Halley 方法逼近根答案

【问题标题】：Approximation of a root using Halley's method in R在 R 中使用 Halley 方法逼近根
【发布时间】：2014-08-07 21:29:34
【问题描述】：

好的，所以我正在开发一个小型 R 程序，以便使用 Halley 方法进行近似。基本上我需要能够使用哈雷的方法将 59^(1/7) 的值近似到小数点后 9 位。

我对哈雷方法的一阶递归关系是这样的：

Xn+1 = Xn - [ f(Xn) / ( f'(Xn) - (f(Xn)f''(Xn)/2f'(Xn)) ) ]

到目前为止，这是我拥有的代码。

halleysMethodApprox = function(ftnH, x0, m0, k0, tol = 1e-8, max.iter=2) {
    x <- x0
    m <- m0
    k <- k0
    fx <- ftnH(x, m, k)
    iter <- 0
    b <- fx[1]/(fx[2] - (fx[1]*fx[3])/(2*fx[2]) )
    while( (abs(fx[1] - x) > tol) && (iter < max.iter) ) {
            # calculate X(n+1)
            b <- ( fx[1]/(fx[2] - ( (fx[1]*fx[3])/(2*fx[2]) ) ))
            x <- x - b
            fx <- ftnH(x, m-1, 0)
            iter <- iter + 1
            cat("At iteration", iter, "value of x is: ", x, "\n")
    }
    if( abs(x) > tol ) {
            cat("Algorithm failed to converge\n")
            return(NULL)
    } else {
            cat("Algorithm converged\n")
            return(x)
    }
}

这个函数用于生成一个包含 x 的函数及其导数的向量。

ftnH = function (x, m, k) {
        fx <- x^m - k
        dfx <- (m*x^(m-1))
        ddfx <- ((m-1)*m)*x^(m-2)
        return (c(fx, dfx, ddfx))
}


print(halleysMethodApprox(ftnH, 59, (1/7), 0))

我不太确定我应该如何使用上述哈雷方法的定义在数值上近似 59^(1/7)。

【问题讨论】：

我很不清楚你的问题是什么。您面临什么编程挑战？
我面临的编程挑战是重写一阶递归关系，使我能够近似函数 f(x) x^m - k 的值，其中 x = 59 和 m = 1/7

标签： r numerical-methods calculus

【解决方案1】：

如果你想计算x=59^(1/7)，你使用哈雷方法的函数是f(x) = x^7-59 = 0。您的代码中有一些拼写错误和小错误。这是一个更正的版本：

halleysMethodApprox = function(ftnH, x0, m0, k0, tol = 1e-8, max.iter=100) {
  # x0: starting estimate, eqn: x = k0^(1/m0) 
  x <- x0
  m <- m0
  k <- k0
  fx <- ftnH(x0, m0, k0)
  iter <- 0
  b <- fx[1]/(fx[2] - (fx[1]*fx[3])/(2*fx[2]) )
  while( (abs(fx[1]) > tol) && (iter < max.iter) ) {
    # calculate X(n+1)
    b <- ( fx[1]/(fx[2] - ( (fx[1]*fx[3])/(2*fx[2]) ) ))
    x <- x - b
    fx <- ftnH(x, m0, k0)
    iter <- iter + 1
    cat("At iteration", iter, "value of x is: ", x, "\n")
  }
  if(abs(fx[1]) > tol) {
    cat("Algorithm failed to converge\n")
    return(NULL)
  } else {
    cat("Algorithm converged\n")
    return(x)
  }
}

ftnH = function (x, m, k) {
  fx <- x^m - k
  dfx <- (m*x^(m-1))
  ddfx <- ((m-1)*m)*x^(m-2)
  return (c(fx, dfx, ddfx))
}

halleysMethodApprox(ftnH, 1, 7, 59) 
 # At iteration 1 value of x is:  1.320442 
 # At iteration 2 value of x is:  1.655396 
 # At iteration 3 value of x is:  1.78716 
 # At iteration 4 value of x is:  1.790519 
 # At iteration 5 value of x is:  1.790519 
 # Algorithm converged
 # [1] 1.790519

【讨论】：