숫자를 추가 할 때 파이썬/부동 소수점 단위의 부동 소수점 정밀도가

Question

numpy와 함께 사용되는 숫자가 너무 적어서 몇 가지 문제가 있습니다. float64 정밀도가 손실되면 함수에 수레를 더할 때 사실과 수치 적 통합을 통해 내 문제를 추적하는 데 몇 주가 걸렸습니다. 합계가 아닌 제품으로 수학적으로 동일한 계산을 수행하면 값이 올바르게됩니다. 당신이 볼 수 있듯이

from matplotlib.pyplot import * 
from numpy import vectorize, arange 
import math 

def func_product(x): 
    return math.exp(-x)/(1+math.exp(x)) 

def func_sum(x): 
    return math.exp(-x)-1/(1+math.exp(x)) 

#mathematically, both functions are the same 

vecfunc_sum = vectorize(func_sum) 
vecfunc_product = vectorize(func_product) 

x = arange(0.,300.,1.) 
y_sum = vecfunc_sum(x) 
y_product = vecfunc_product(x) 

plot(x,y_sum, 'k.-', label='sum') 
plot(x,y_product,'r--',label='product') 

yscale('symlog', linthreshy=1E-256) 
legend(loc='lower right') 
show() 

가 0에 흩어져있다 매우 낮은 또는 합계 값이 동안 정확히 제로입니다 : 여기

은 코드 샘플과 결과의 플롯이다 곱해진 값은 괜찮습니다 ...

누군가가 도와 줄 수 있습니까? 고마워요!

Answer 1

부동 소수점 정밀도는 반올림 오류로 인해 더하기/빼기에 매우 민감합니다. 결국 1+exp(x)이 커져서 exp (x)에 1을 더하면 exp (x)와 같은 결과가 나온다. exp(x) == 1e16 주위 어딘가에 배정 밀도 :

>>> (1e16 + 1) == (1e16) 
True 
>>> (1e15 + 1) == (1e15) 
False 

주 math.log(1e16) 것을는 약 37입니다 - 상황이 플롯에 미쳐 곳 정도된다.

당신은 같은 문제가있을 수 있지만 다른 비늘 : 당신의 정권에서 점의 대다수

>>> (1e-16 + 1.) == (1.) 
True 
>>> (1e-15 + 1.) == (1.) 
False 

, 당신의 func_product 실제로 계산된다

exp(-x)/exp(x) == exp(-2*x) 

가 어떤 이유 그래프의 기울기가 -2가됩니다.

다른 극단에 가지고 가서, 다른 버전 (적어도 약) 계산하는 것 :

exp(-x) - 1./exp(x) 

약 문제는 당신 func_sum이 numerically unstable 때문이라는 것이다

exp(-x) - exp(-x) 

Answer 2

인 그것은 매우 가까운 두 값 사이의 뺄셈을 포함합니다.

math.exp(200).hex() 
0x1.73f60ea79f5b9p+288 

(math.exp(200) + 1).hex() 
0x1.73f60ea79f5b9p+288 

(1/(math.exp(200) + 1)).hex() 
0x1.6061812054cfap-289 

math.exp(-200).hex() 
0x1.6061812054cfap-289 

: math.exp(200)-1을 추가하는 효과가 없기 때문에 그것은 64 비트 부동 소수점 정밀도 외부에 있으므로 func_sum(200) 계산에

예를 들어, math.exp(-200) 및 1/(1+math.exp(200))은 동일한 값을 가질 이것은 func_sum(200)이 0을주는 이유를 설명하지만, x 축에서 벗어나는 점은 무엇입니까? 이는 또한 부동 소수점 부정확도 때문에 발생합니다. math.exp(-x)이 1/math.exp(x)과 같지 않은 경우가 종종 발생합니다. 이상적으로 math.exp(x)은 e^x에 가장 가까운 부동 소수점 값이고 1/math.exp(x)은 math.exp(x)에 의해 계산 된 부동 소수점 수의 역수에 가장 가까운 부동 소수점 값이며 반드시 e^-x이 아니어야합니다.실제로, math.exp(-100) 및 1/(1+math.exp(100))는 마지막 단위에서 차이가 매우 가까운 사실입니다 (있는 경우)

math.exp(-100).hex() 
0x1.a8c1f14e2af5dp-145 

(1/math.exp(100)).hex() 
0x1.a8c1f14e2af5cp-145 

(1/(1+math.exp(100))).hex() 
0x1.a8c1f14e2af5cp-145 

func_sum(100).hex() 
0x1.0000000000000p-197 

은 그래서 당신이 실제로 계산 한 것은 math.exp(-x) 및 1/math.exp(x) 사이의 차이입니다. math.pow(2, -52) * math.exp(-x) 함수 행을 추적하여 양수 값 func_sum (52는 64 비트 부동 소수점의 유효 숫자의 크기 임)을 통과하는지 확인할 수 있습니다.

Answer 3

이것은 catastrophic cancellation의 예입니다. 계산이 틀려서 첫 번째 점에서

살펴 보자, 때 x = 36.0

In [42]: np.exp(-x) 
Out[42]: 2.3195228302435691e-16 

In [43]: - 1/(1+np.exp(x)) 
Out[43]: -2.3195228302435691e-16 

In [44]: np.exp(-x) - 1/(1+np.exp(x)) 
Out[44]: 0.0 

는 치명적인 취소를 피할 수 있도록, 거의 같은 숫자를 뺄하지 않습니다 func_product를 사용하여 계산. 당신이 math.expnp.exp로 변경하는 경우 그런데

---

는, 당신은 (느린입니다) np.vectorize 제거 할 수 있습니다`np.float64`가 충분하지 않습니다

def func_product(x): 
    return np.exp(-x)/(1+np.exp(x)) 

def func_sum(x): 
    return np.exp(-x)-1/(1+np.exp(x)) 

y_sum = func_sum_sum(x) 
y_product = func_product_product(x)