mag*_*852 3 python math calculus
from math import sqrt
a=1e-8
b=10
c=1e-8
x1 = ((-b)-sqrt((b**2)-(4*a*c)))/(2*a)
x2 = ((-b)+sqrt((b**2)-(4*a*c)))/(2*a)
print 'x1 = {}'.format(x1)
print 'x2 = {}'.format(x2)
print (4*a*c)
print (sqrt(b**2-4*a*c))
print b**2
print 2*a
Run Code Online (Sandbox Code Playgroud)
当我运行该程序时,它返回:
x1 = -1e+09
x2 = 0.0
4e-16
10.0
100.0
2e-08
Run Code Online (Sandbox Code Playgroud)
我需要的是x2等于-1e-9.
问题似乎与
sqrt((b**2)-(4*a*c))
Run Code Online (Sandbox Code Playgroud)
因为它给出10作为结果,显然因为4*(10 ^ -8)*(10 ^ -8)几乎等于0,并且被python认为是0.
这导致:
sqrt((b**2)-(4*a*c)) = sqrt(b**2) = sqrt(10**2) = 10
Run Code Online (Sandbox Code Playgroud)
任何帮助将不胜感激
使用十进制模块:
from decimal import Decimal
a = Decimal('1E-8')
b = 10
c = Decimal('1E-8')
x1 = ((-b)-((b**2)-(4*a*c)).sqrt())/(2*a)
x2 = ((-b)+((b**2)-(4*a*c)).sqrt())/(2*a)
print 'x1 = {}'.format(x1)
print 'x2 = {}'.format(x2)
Run Code Online (Sandbox Code Playgroud)
结果是
x1 = -999999999.999999999000000000
x2 = -1.0000000000E-9
Run Code Online (Sandbox Code Playgroud)
您不需要额外的精度来解决这个问题:Python float
已经具备足够的精度来完成这项工作.你只需要一个(稍微)聪明的算法.
你的问题源于减去两个几乎相等的计算值:对于b
正数和大数(与a
和相比c
),当你这样做时-b + sqrt(b*b-4*a*c)
,你最终会得到一个具有较大相对误差的结果.但请注意,这个问题只适用于两个根中的一个:在-b - sqrt(b*b-4*a*c)
,没有这样的问题.类似地,对于b
大和负,第一根是好的,但第二根可能会失去准确性.
解决方案是使用您现有的公式计算任何根没有取消问题,然后为另一个根使用不同的公式(实质上,使用您知道两个根的产品的事实c / a
) .那个公式是2c / (-b +/- sqrt(b*b-4*a*c))
.
这是一些示例代码.它用于math.copysign
选择不会导致取消错误的标志:
>>> from math import sqrt, copysign
>>> def quadratic_roots(a, b, c):
... discriminant = b*b - 4*a*c
... q = -b - copysign(sqrt(discriminant), b)
... root1 = q / (2*a)
... root2 = (2*c) / q
... return root1, root2
...
>>> quadratic_roots(a=1e-8, b=10, c=1e-8)
>>> (-1000000000.0, -1e-09)
Run Code Online (Sandbox Code Playgroud)
这涉及数值不稳定的最严重可能原因.在计算判别式时,还有第二个可能的原因,如果b*b
恰好非常接近4*a*c
.在这种情况下,可能会丢失一半正确的有效数字(因此每个根目录只能获得7-8个准确数字).获得全精度结果在这种情况下就需要使用计算扩展精度的判别.
关于失去重要性的维基百科文章包含对这个问题的有用讨论.
归档时间: |
|
查看次数: |
571 次 |
最近记录: |