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

 幸运之星07812 发布于 2023-02-09 14:23

我警告你,这可能令人困惑,我写的代码更多的是思维导图而不是完成代码.

我正在尝试使用Newton-Raphson方法来求解方程.我无法弄清楚的是如何写这个

在此输入图像描述

Python中的等式,用于计算最后一次近似(xn)的下一个近似值(xn + 1).我必须使用循环,越来越接近真实的答案,并且当近似值之间的变化小于变量h时,循环应该终止.

    如何编写等式的代码?

    当近似值不再变化时,如何终止循环?

    计算方程f的导数,在点x,精度为h(这用于求解方程式())

    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
    

Floris.. 8

以下是NR解算器的实现,扩展了您上面写的内容(完整,有效).我添加了一些额外的行来显示正在发生的事情......

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(你可以争论的价值=符号的数值方法-我喜欢更加强调<自己,认为再多一次迭代不会受到伤害.)

1 个回答
  • 以下是NR解算器的实现,扩展了您上面写的内容(完整,有效).我添加了一些额外的行来显示正在发生的事情......

    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(你可以争论的价值=符号的数值方法-我喜欢更加强调<自己,认为再多一次迭代不会受到伤害.)

    2023-02-09 14:26 回答
撰写答案
今天,你开发时遇到什么问题呢?
立即提问
热门标签
PHP1.CN | 中国最专业的PHP中文社区 | PNG素材下载 | DevBox开发工具箱 | json解析格式化 |PHP资讯 | PHP教程 | 数据库技术 | 服务器技术 | 前端开发技术 | PHP框架 | 开发工具 | 在线工具
Copyright © 1998 - 2020 PHP1.CN. All Rights Reserved 京公网安备 11010802041100号 | 京ICP备19059560号-4 | PHP1.CN 第一PHP社区 版权所有