Julia 中的 exp 函数评估为 0

exp function in Julia evaluating to 0

我想在 Julia 中计算和绘制波函数的概率密度。我写了一小段 Julia 代码来评估以下函数:

P = A^2 e^{-(\sqrt(Cm)/\hbar)x^2}

Julia(不完整)代码是:

set_bigfloat_precision(100)
A = 10
C = 5
m = BigFloat(9.10938356e-31)
ℏ = BigFloat(1.054571800e-34)
t = exp(-(sqrt(C * m) / ℏ))

我计算 t 的最后一行给出 0.000000000000...。我也尝试设置 BigFloat 的精度。没有运气!我究竟做错了什么?帮助表示赞赏。

Chris Rackauckas 在评论中指出您输入的公式有误。我认为无论如何回答这个问题都足够有趣

让我们分解一下,看看我们在筹集什么:

A = 10
C = 5
m = BigFloat(9.10938356e-31)
h = BigFloat(1.054571800e-34)

z = -sqrt(C * m)/h
t = exp(z)

所以 z =-2.0237336022083455711032042949257e+19 所以非常粗略z=-2e19)

大致t=exp(-2e19)(即t=1/((e^(2*10^19))) 这是一个非常小的数字。

考虑一下 exp(big"-1e+10") = 9.278...e-4342944820exp(big"-1e+18") = 2.233...e-434294481903251828

是的,朱莉娅说: exp(big"-2e+19) = 0.0000

exp(big"-2e+19) 是一个很小的数字。 我希望这能让我们了解情况。数量很少。


所以 julia 依赖于 MPFR 的 BigFloats 你可以试试MPFR online。精度为 8192,exp(-2e10)=0 结果一样。

现在,我们关心的不是精度。 而是指数的范围。

MPFR 使用类似于 IEEE 风格的浮点数,其中精度是尾数的长度,然后你有一个指数。 2^exponent * mantissa

所以指数的范围是有限制的。

参见:MPFR docs:

Function: mpfr_exp_t mpfr_get_emin (void) Function: mpfr_exp_t mpfr_get_emax (void)

Return the (current) smallest and largest exponents allowed for a floating-point variable. The smallest positive value of a floating-point variable is one half times 2 raised to the smallest exponent and the largest value has the form (1 - epsilon) times 2 raised to the largest exponent, where epsilon depends on the precision of the considered variable.

现在 julia 确实将这些设置为相当默认的 MPFR 编译允许的最大范围。我一直在挖掘 MPFR 来源,试图找到它的设置位置,但找不到。我相信这与 Int64 可以容纳的最大故障有关。

Base.MPFR.get_emin() = -4611686018427387903 =typemin(Int64)>>1 + 1

您可以调整它,但只能向上调整。

所以无论如何

0.5*big"2.0"^(Base.MPFR.get_emin()) = 8.5096913117408361391297879096205e-1388255822130839284

但是 0.5*big"2.0"^(Base.MPFR.get_emin()-1) = 0.00000000000...


现在我们知道

exp(x) = 2^(log(2,e)*x)

所以我们可以exp(z) = 2^(log(2,e)*z)
log(2,e)*z = -29196304319863382016
Base.MPFR.get_emin() = -4611686018427387903

因此,由于指数(大约 -2.9e19)小于允许的最小指数(大约 -4.3e17)。 发生下溢。

因此你的答案是为什么你得到零。

可能(或不可能)用 Int128 指数重新编译 MPFR,但 julia 没有。

也许 julia 应该抛出一个下溢异常。 鼓励免费在 Julia Bug Tracker 上报告该问题。