Python牛顿·拉夫森;计算精度
Python newton raphson; precision in calculations
我在python中写了这个函数:
def f2(x):
return (5.0*x + log1p(x) - 10000.0)
def dfdx2(x):
return (5.0-(1.0/x))
def newtonRaphson2(f, dfdx, x, tol):
x0 = x
for i in range(1, 2000):
if f(x) == 0.0:
return x
if dfdx(x) == 0.0:
print(dfdx(x))
break
x = x - (f(x) / dfdx(x))
#print(x)
Er = abs(x0-x)/abs(x0)
if Er <= tol:
return x
print(Er)
x0 = x
return x
然后我这样执行:
task2 = newtonRaphson2(f2, dfdx2, 1, 0.000001);
print(task2)
对于输出检查 Er 在它之前输出的最终精度为 4.245665128208564e-05 returns x.
X 在 1998.479871524306 处返回,这是一个很好的估计,但我最好至少将其降低到 1.0e-06。将 tol 变量更改为 1.0e-08 似乎没有任何作用。
我猜也许将每个变量都变成双精度是个更好的主意,但我仍然不知道为什么我的代码停在它该做的地方。我对 python 也不是那么稳定,这就是我问的原因。我已经编写了其中一个有效的方法,但它用于一个简单得多的方程式。
正如@alex-dubrovsky提到的,意味着收敛的计算应该用条件循环来实现,即:
while True:
if f(x) == 0.0:
return x
if dfdx(x) == 0.0:
print(dfdx(x))
break
x = x - (f(x) / dfdx(x))
#print(x)
Er = abs(x0-x)/abs(x0)
if Er <= tol:
return x
print(Er)
x0 = x
使用这种方法,您总是面临无限循环的风险,但这或多或少没问题,正如算法所暗示的那样"running until converged"
你的代码工作正常,一旦你正确缩进并添加 from math import log1p
。只需将 print(Er)
行放在计算后的右边即可查看其最终值。 Er
达到 ~10^-9
。这对我有用:
from math import log1p
def f2(x):
return (5.0*x + log1p(x) - 10000.0)
def dfdx2(x):
return (5.0-(1.0/x))
def newtonRaphson2(f, dfdx, x, tol):
x0 = x
for i in range(1, 2000):
if f(x) == 0.0:
return x
if dfdx(x) == 0.0:
print(dfdx(x))
break
x = x - (f(x) / dfdx(x))
#print(x)
Er = abs(x0-x)/abs(x0)
print('Er = {}'.format(Er))
if Er <= tol:
return x
x0 = x
return x
x = newtonRaphson2(f2, dfdx2, 1, 0.0000001)
print 'X = {}'.format(x)
输出为:
Er = 2498.5767132
Er = 0.200506616666
Er = 4.24566512821e-05
Er = 8.49642413214e-09
X = 1998.47987152
考虑在此处使用 while
。 Newton-Raphson 算法通常收敛得非常快,因此在大多数情况下 运行 不需要多次迭代。
这给出了相同的结果:
from math import log1p
def f2(x):
return (5.0*x + log1p(x) - 10000.0)
def dfdx2(x):
return (5.0-(1.0/x))
def newtonRaphson2(f, dfdx, x, tol):
x0 = x
Er = 1
while Er >= tol:
if f(x) == 0.0:
return x
if dfdx(x) == 0.0:
print(dfdx(x))
break
x = x - (f(x) / dfdx(x))
#print(x)
Er = abs(x0-x)/abs(x0)
print('Er = {}'.format(Er))
x0 = x
return x
x = newtonRaphson2(f2, dfdx2, 1, 0.0000001)
print 'X = {}'.format(x)
我在python中写了这个函数:
def f2(x):
return (5.0*x + log1p(x) - 10000.0)
def dfdx2(x):
return (5.0-(1.0/x))
def newtonRaphson2(f, dfdx, x, tol):
x0 = x
for i in range(1, 2000):
if f(x) == 0.0:
return x
if dfdx(x) == 0.0:
print(dfdx(x))
break
x = x - (f(x) / dfdx(x))
#print(x)
Er = abs(x0-x)/abs(x0)
if Er <= tol:
return x
print(Er)
x0 = x
return x
然后我这样执行:
task2 = newtonRaphson2(f2, dfdx2, 1, 0.000001);
print(task2)
对于输出检查 Er 在它之前输出的最终精度为 4.245665128208564e-05 returns x.
X 在 1998.479871524306 处返回,这是一个很好的估计,但我最好至少将其降低到 1.0e-06。将 tol 变量更改为 1.0e-08 似乎没有任何作用。
我猜也许将每个变量都变成双精度是个更好的主意,但我仍然不知道为什么我的代码停在它该做的地方。我对 python 也不是那么稳定,这就是我问的原因。我已经编写了其中一个有效的方法,但它用于一个简单得多的方程式。
正如@alex-dubrovsky提到的,意味着收敛的计算应该用条件循环来实现,即:
while True:
if f(x) == 0.0:
return x
if dfdx(x) == 0.0:
print(dfdx(x))
break
x = x - (f(x) / dfdx(x))
#print(x)
Er = abs(x0-x)/abs(x0)
if Er <= tol:
return x
print(Er)
x0 = x
使用这种方法,您总是面临无限循环的风险,但这或多或少没问题,正如算法所暗示的那样"running until converged"
你的代码工作正常,一旦你正确缩进并添加 from math import log1p
。只需将 print(Er)
行放在计算后的右边即可查看其最终值。 Er
达到 ~10^-9
。这对我有用:
from math import log1p
def f2(x):
return (5.0*x + log1p(x) - 10000.0)
def dfdx2(x):
return (5.0-(1.0/x))
def newtonRaphson2(f, dfdx, x, tol):
x0 = x
for i in range(1, 2000):
if f(x) == 0.0:
return x
if dfdx(x) == 0.0:
print(dfdx(x))
break
x = x - (f(x) / dfdx(x))
#print(x)
Er = abs(x0-x)/abs(x0)
print('Er = {}'.format(Er))
if Er <= tol:
return x
x0 = x
return x
x = newtonRaphson2(f2, dfdx2, 1, 0.0000001)
print 'X = {}'.format(x)
输出为:
Er = 2498.5767132
Er = 0.200506616666
Er = 4.24566512821e-05
Er = 8.49642413214e-09
X = 1998.47987152
考虑在此处使用 while
。 Newton-Raphson 算法通常收敛得非常快,因此在大多数情况下 运行 不需要多次迭代。
这给出了相同的结果:
from math import log1p
def f2(x):
return (5.0*x + log1p(x) - 10000.0)
def dfdx2(x):
return (5.0-(1.0/x))
def newtonRaphson2(f, dfdx, x, tol):
x0 = x
Er = 1
while Er >= tol:
if f(x) == 0.0:
return x
if dfdx(x) == 0.0:
print(dfdx(x))
break
x = x - (f(x) / dfdx(x))
#print(x)
Er = abs(x0-x)/abs(x0)
print('Er = {}'.format(Er))
x0 = x
return x
x = newtonRaphson2(f2, dfdx2, 1, 0.0000001)
print 'X = {}'.format(x)