Python - 实现数值方程求解器 (Newton-Raphson)

Question

我警告你，这可能会造成混淆，而且我编写的代码与其说是完成的代码，不如说是思维导图。

我正在尝试实现 Newton-Raphson 方法来求解方程。我想不通的是如何写这个

Python 中的方程，用于根据上一个近似值 (xn) 计算下一个近似值 (xn+1)。我必须使用循环来越来越接近真实答案，并且当近似值之间的变化小于变量 h 时循环应该终止。

1. 如何编写方程式的代码？

2. 当近似值不再变化时如何终止循环？

   计算方程 f 在 x 点的导数，精度为 h(用于 solve() 的方程)

   def derivative(f, x, h):
       deriv = (1.0/(2*h))*(f(x+h)-f(x-h))
       return deriv
   

   数值方程求解器

   应该循环直到近似值之间的差异小于h

   def solve(f, x0, h):
       xn = x0
       prev = 0
   
       while ( approx - prev > h):
            xn = xn - (f(xn))/derivative(f, xn, h)
   
       return xn
   

Answer 1

这是 N-R 求解器的实现，它扩展了您在上面编写的内容(完整，有效)。我添加了一些额外的行来显示正在发生的事情......

def derivative(f, x, h):
      return (f(x+h) - f(x-h)) / (2.0*h)  # might want to return a small non-zero if ==0

def quadratic(x):
    return 2*x*x-5*x+1     # just a function to show it works

def solve(f, x0, h):
    lastX = x0
    nextX = lastX + 10* h  # "different than lastX so loop starts OK
    while (abs(lastX - nextX) > h):  # this is how you terminate the loop - note use of abs()
        newY = f(nextX)                     # just for debug... see what happens
        print "f(", nextX, ") = ", newY     # print out progress... again just debug
        lastX = nextX
        nextX = lastX - newY / derivative(f, lastX, h)  # update estimate using N-R
    return nextX

xFound = solve(quadratic, 5, 0.01)    # call the solver
print "solution: x = ", xFound        # print the result

输出:

f( 5.1 ) =  27.52
f( 3.31298701299 ) =  6.38683083151
f( 2.53900845771 ) =  1.19808560807
f( 2.30664271935 ) =  0.107987672721
f( 2.28109300639 ) =  0.00130557566462
solution: x =  2.28077645501

编辑 - 您还可以检查 newY 的值并在它“足够接近零”时停止 - 但通常你会继续这样做直到 x 发生变化是<=h (您可以在数值方法中争论 = 符号的值 - 我自己更喜欢更有力的 <，认为再迭代一次不会有坏处。)。