2

我知道这已经在 C/C++ 中解决了,但是我对这些语言不够满意,无法将其转换为 python。我正在尝试在 python中创建它。我能来的最接近的是:

#This is meant to work for functions of the form x^a + b = 0

def secant(base, exp=2, it=20):
    def f(x):
        return x**exp - base
    x1 = base / float(exp**2)
    xnm1 = x1 - 5
    xnm2 = x1 + 5
    xn = 0
    for n in range(it):
        q = (xnm1-xnm2)/float(f(xnm1)-f(xnm2))
        xn = xnm1 - (f(xnm1)*q)
        xnm1, xnm2 = xn, xnm1
    return xn

print secant(2, 2)

这将返回错误:

Traceback (most recent call last):
  File "/Users/Joe/Desktop/secant.py", line 16, in <module>
    print secant(2, 2)
  File "/Users/Joe/Desktop/secant.py", line 11, in secant
    q = (xnm1-xnm2)/float(f(xnm1)-f(xnm2))
ZeroDivisionError: float division by zero

但是,我能够编写基于此代码的牛顿方法。如果有帮助,这里是:

def newton(base, exp=2, it=20):
    def f(x):
        return x**exp - base
    def df(x):
        return exp*(x**(exp-1))
    x1 = base / float(exp**2)
    xnp = x1
    xn = 0
    for n in range(it):
        xn = xnp - ((f(xnp)/df(xnp)))
        xnp = xn
    return xn

以下方法在 20 次迭代后给出了 12 位精度的答案。任何帮助,将不胜感激。

4

2 回答 2

2

您得到除以零误差,因为该算法在 Python 浮点数的精度范围内收敛到答案。除了迭代最大次数(以避免无限循环)之外,您还应该检查最后两个猜测是否“足够接近”。

于 2012-10-11T02:59:20.430 回答
2

重要的是浮点数的轮次。在前一种情况下,经过一些迭代,和之间的差异f(xnm1)很小f(xnm2),以至于浮点数无法表示它,因此将其四舍五入为零,然后抛出错误。

在后一种中,我们直接使用导数计算该点的斜率,虽然它很小,但它是非零的(但在接下来的迭代中它几乎保持相同的值),所以它不会抛出错误。

这是一个例子,和IEEE 浮点数

于 2012-10-11T03:16:38.553 回答