大约 e^x

10

我想近似计算ex函数。

是否可以使用多个基于样条的方法来实现?例如在x1x2之间,

y1 = a1x + b1,在 x2 和 x3 之间,

然后

y2 = a2x + b2

等等。

这是为专用FPGA硬件设计的,而不是通用CPU。因此,我需要自己创建该函数。精度并不是很重要。此外,我不能负担得起多个乘法电路和/或多个移位/加法器。同时,我希望它比CORDIC函数小得多,事实上大小非常关键。


2
您打算用什么范围的x值来进行逼近? - Cassidy Laidlaw
6
幂级数。幂级数是指形如 $a_0 + a_1 x + a_2 x^2 + \cdots$ 的无限多项式,其中 $a_0, a_1, a_2, \ldots$ 是常数,$x$ 是变量。幂级数在数学中广泛应用,尤其是在微积分、数论和物理学等领域中扮演着重要的角色。幂级数也是导出一些特殊函数的基础,例如正弦函数和指数函数。 - user786653
2
C++标准库中有exp()函数,为什么要避免使用它呢?通常它的速度很快。 - George Gaál
2
我的应用程序实际上不是C或C++,而是专用硬件,因此我正在自己编写函数。Power函数很好,但我更喜欢一些操作较少的函数。 - trican
@user786653:绝对不是幂级数。那是一个理论数学定义,而不是数值数学定义。同一页上有更多实用的公式,例如连分数 - MSalters
显示剩余9条评论
10个回答

27
这样的策略怎么样? 它使用以下公式:

ex = 2 x/ln(2)

  1. 预先计算 1/ln(2)
  2. 将此常量乘以您的参数(1次乘法)
  3. 使用二进制移位将2提高到幂的整数部分(假设exp+mantissa格式)
  4. 根据小数Power-of-2余数进行调整(可能需要第二次乘法)
我知道这不是一个完整的解决方案,但它只需要一次乘法,并将剩余的问题简化为近似于2的分数幂,这在硬件实现中应该更容易。
此外,如果您的应用程序足够专业化,可以尝试重新推导在硬件上运行的所有数字代码以成为基于e的数字系统,并将浮点硬件实现为在基本e中工作。 然后根本不需要转换。

1
谢谢Lucas - 这正是我需要的,甚至比我想象中的还要好。非常感谢! - trican
很高兴听到这个消息。听起来你有一些有趣的设计折衷方案。 - Lucas
2
@trican,有一篇很好的论文介绍了使用查找表和定点算术来实现身份和范围缩减,从而实现单精度浮点数的合理精度:http://www.loria.fr/~detreyje/publications/DetDin_fpt_2005.pdf - Chiggs
备用PDF链接:http://perso.citi-lab.fr/fdedinec/recherche/publis/2005-FPT.pdf - Lucas
exp()只需要两次乘法和一个位移,神奇的位移! - 16807

14
如果x是整数,你可以一遍又一遍地将e乘以自己。

如果x不是整数,则可以使用上述方法计算 efloor(x),然后再乘以一个小的修正项。这个修正项可以使用多种近似方法轻松地计算得出。其中一种方法如下:

ef1 + f(1 + f/2(1 + f/3(1 + f/4))),其中fx的小数部分。

这来自于(优化的)幂级数展开式,非常适用于小值的x。如果需要更高的精度,只需在级数后面添加更多的项。

这个math.stackexchange问题包含了一些额外巧妙的答案。

编辑:注意,有一种更快速的计算 en 的方法,称为平方幂运算


2
最佳的整数解决方案不是这个O(n)解决方案。分治算法预先计算e^1、e^2、e^4、e^8等。然后,您可以取出与x中位相对应的因子。这是O(logN)。例如,对于x=255,只需要进行8次乘法,而不是254次。 - MSalters
谢谢 - 但我想尽量减少乘法运算,我只需要一次乘法运算。 - trican
但是为什么呢?你是否真的遇到了性能问题,还是这只是过早优化? - Jonathan Grynspan
@Jonathan - 这不是针对CPU的,而是针对专用硬件的。我已经更新了我的问题以澄清这一点。对于造成的困惑,我感到抱歉。 - trican
@Jonathan 因为拥有一个O(n)的指数函数显然会导致性能不佳。在系统层面上,过早地进行优化并不是坏事。 - alternative
这正是我需要的,以整数数学版本计算e^x。 - zawy

4

首先,是什么促使这个近似值的出现?换句话说,直接使用exp(x)有什么问题吗?

那么,exp(x)的典型实现方式是:

  • 找到一个整数k和一个浮点数r,使得x=k*log(2) + r,且r在-0.5*log(2)和0.5*log(2)之间。
  • 通过这种简化,exp(x)变为2k*exp(r)
  • 计算2k很容易。
  • 标准的exp(x)实现使用Remes算法来得出一个最小极大多项式,以逼近exp(r)
  • 您也可以尝试相同的方法,但使用降低阶数的多项式。

关键是:无论您做什么,很可能您的函数速度都比直接调用exp()慢得多。大部分exp()函数的功能都是由您计算机上的数学协处理器实现的。即使降低精度,在软件中重新实现该功能的速度也会比直接使用exp()慢一个数量级。


Remez*和大多数实际上使用以边界为中心的Pade逼近,以使该范围内的误差尽可能小。对于给定的输入x,误差等于有界误差乘以2^k,当输入很大时通常会破坏这些逼近的大部分...我“相信”实际实现同时采用了Pade逼近和迭代改进根查找方法的反函数减去输入。 - nimig18
为什么 r 应该位于 -0.5log(2)0.5log(2) 之间,而不是 (0, 1) - Elinx

3

对于硬件方面,如果您需要其精度与二进制位一致,我有一个很棒的解决方案。(否则,只需像上面那样进行近似计算)。恒等式为exp(x) = cosh(x) + sinh(x),其中超越正弦和余弦可以使用CORIC技术计算,更妙的是,它们是最快的CORDIC函数之一,这意味着它们看起来几乎像乘法而不是几乎像除法!

这意味着您只需要占用大约一个阵列乘法器的面积,在仅两个周期内就可以计算任意精度的指数!

请查找CORDIC方法- 这对于硬件实现非常棒。

另外一个硬件方法是使用一个小表格结合其他人提到过的公式:exp(x + y) = exp(x) * exp(y)。您可以将数字分成小的位字段 - 每次4或8位 - 然后查找该位字段的指数。这可能只适用于窄计算,但这是另一种方法。


2

感谢@jdberton添加这个和链接。这种方法似乎非常有趣,但是您确定上面的代码片段正确吗?我尝试了一些值,结果似乎甚至不接近? - trican
我认为对于大数值来说,这可能是不准确的。你可以通过一些工作找到更好的Pade逼近来获得更好的范围。它对我有用,因为我不需要完全精确的东西。 - jdbertron
Schraudolph的方法是完美的。如果精度可接受,我认为它无法更快。在他的论文中,他确定平均相对误差约为4%。来源:http://nic.schraudolph.org/pubs/Schraudolph99.pdf - Gigo
这是Schraudolph方法的一个更现代化的实现,使用单点浮点数而不是双精度浮点数(这是一种浪费,因为只有双精度浮点数的上32位被写入)。http://www.machinedlearnings.com/2011/06/fast-approximate-logarithm-exponential.html - Mark Lakata

2
这不是您请求的平滑样条插值,但它在计算上更有效:
float expf_fast(float x) {
   union { float f; int i; } y;
   y.i = (int)(x * 0xB5645F + 0x3F7893F5);
   return (y.f);
}

图形输出 image


1
这不适用于定制的FPGA,但值得一提。

http://www.machinedlearnings.com/2011/06/fast-approximate-logarithm-exponential.html

而且源代码:

https://code.google.com/archive/p/fastapprox/downloads

“更快”的实现只涉及3个步骤(乘法,加法,将浮点数转换为整数),最后还要进行一次浮点数到整数的转换。根据我的经验,它的准确率为2%,如果您不关心实际值而是在对数似然最大化迭代中使用该值,则可能足够。”

1

当然是“可能的”。但有几个问题需要考虑:

  1. 你对精度有什么要求?

  2. 你是否愿意使用高阶样条?

  3. 你愿意在这方面花费多少内存?线性函数在足够小的区间内将指数函数逼近到任何所需的精度,但可能需要非常小的区间。

编辑:

根据提供的额外信息,我进行了快速测试。范围缩减总是可以用于指数函数。因此,如果我想计算任何x的exp(x),那么我可以将问题重写为...

y = exp(xi + xf) = exp(xi)*exp(xf)

其中 xi 是 x 的整数部分,xf 是小数部分。整数部分很简单。将 xi 计算成二进制形式,然后通过重复平方和乘法运算,可以在相对较少的操作次数内计算出 exp(xi)。(使用 2 的幂和其他间隔的技巧可以为追求速度的用户提供更快的速度。)

现在剩下的就是计算 exp(xf)。我们能否使用具有线性段的样条函数,在区间 [0,1] 上仅使用 4 个线性段,以精度为 0.005 计算 exp(xf)?

我几年前编写了一个函数来解决这个问题,它将近似于给定阶数的函数,使其最大误差在固定容差范围内。该代码需要在区间 [0,1] 上使用 8 个线性段才能使用分段线性样条函数达到所需的容差。如果我选择将区间进一步缩小到 [0,0.5],那么我现在可以实现所要求的容差。

因此,答案很简单。如果您愿意进行范围缩减,将 x 缩小到区间 [0.0.5],然后进行适当的计算,那么是的,您可以使用 4 个线性段的线性样条函数实现所需的精度。

最终,使用硬编码指数函数总是更好的选择。如果exp(x)可用,则上述提到的所有操作肯定比编译器提供的慢。


非常感谢您详细的回复。经过进一步的思考,我可以容忍更高的误差率,可能高达0.05甚至0.1。我以前已经使用过范围缩减的样条函数来处理其他函数,但在这种情况下,我认为Lucas上面的答案更适合低精度要求。另外,关键点是硬件“编译器”中没有指数函数的直接实现。也就是说,我不是在CPU上工作。 - trican

0

Wolfram提供了一些好的方法来近似地表示它,例如级数展开等:

维基百科关于泰勒级数 的页面也展示了一个关于0附近的ex的例子:


3
备选表示法:e^x=z^x 对于 e=z。 - MSalters

0

在C语言中,您可以使用pow(M_E, x)来计算。 (某些平台没有定义M_E;在这些平台上,您可能需要手动指定e的值,该值约为2.71828182845904523536028747135266249775724709369995。)

(正如David在评论中指出的那样,exp(x)pow(M_E, x)更有效率。再次强调,大脑还没有开机。)

您是否有一个使用案例,其中计算ex是已知的瓶颈?如果没有,您应该首先编写易读的代码;只有在明显的方法太慢时才尝试这些优化。


6
pow(M_E, x)?真的吗?通常情况下,pow(a,b) 的实现方式是 exp(b*log(a))。使用 pow 会拖慢速度,而不是加速。请确认是否需要继续使用 pow - David Hammen
这正是我的观点 - 首先正确编写代码,然后再查看其性能。在原始问题中没有说明会每秒调用一百万次或类似情况,因此性能不是立即显而易见的问题。 - Jonathan Grynspan
无论性能如何,exp(x)pow(M_E, x)更简单(并且更便携!)。即使pow()更快,而不是使用exp(),诉诸于它将是过早的优化。 - Keith Thompson
非常正确,我已经更新了我的答案以反映David的更正。你能看出来我还没有喝足够的咖啡吗? :) - Jonathan Grynspan

网页内容由stack overflow 提供, 点击上面的
可以查看英文原文,
原文链接