(指数规模)互补误差函数通常由erfcx ,在数学上定义为erfcx(x):= ex2 erfc(x)。 它经常发生在物理学和化学的扩散问题中。 虽然一些数学环境(如MATLAB和GNU Octave)提供了此功能,但它仅在提供erf()和erfc()的C标准数学库中不存在。
虽然直接基于数学定义可以实现自己的erfcx() ,但它只能在有限的输入域上工作,因为在正半平面erfc()下溢中等幅度的参数,而exp()溢出时,正如在这个问题中所指出的那样。
为了与C一起使用,可以调整一些erfcx()开源实现,例如Faadeeva包中的开源实现,正如对此问题的回应所指出的那样。 但是,这些实现通常不会为给定的浮点格式提供完整的准确性。 例如,使用232个测试向量进行的测试显示,由Faadeeva软件包提供的erfcx()的最大误差在正半平面为8.41 ulps,在负半平面为511.68 ulps。
精确实现的合理界限是4个ulps,对应于英特尔矢量数学库的LA配置文件中的数学函数的准确界限,我发现这对于非平凡的数学函数实现来说是合理的界限,准确性和良好的性能。
如何准确地实现erfcx()和相应的单精度版本erfcxf() ,同时只使用C标准数学库,并且不需要外部库? 我们可以假设C的float double类型映射为IEEE 754-2008 binary32和binary64浮点类型。 可以假设硬件支持融合乘加操作(FMA),因为目前所有主流处理器架构都支持这一功能。
迄今为止我发现的erfcx()实现的最佳方法是基于以下文章:
MM Shepherd和JG Laframboise,“(1 + 2 x)exp(x2)erfc x在0≤x
本文提出了巧妙的转换,将经过缩放的互补误差函数映射为可以直接进行多项式逼近的紧束缚函数。 为了表现,我已经尝试过变化的变化,但所有这些都对准确性有负面影响。 在变换(x - K)/(x + K)中常数K的选择与核心近似的精度之间存在着非明显的关系。 我凭经验确定了与文章不同的“最佳”值。
核心近似和中间结果的变换返回到erfcx结果会产生额外的舍入误差。 为了减轻他们对准确性的影响,我们需要采取补偿步骤,我在之前关于erfcf问题和回答中详细列出了这些步骤。 FMA的可用性极大地简化了这项任务。
生成的单精度代码如下所示:/*
* Based on: M. M. Shepherd and J. G. Laframboise, "Chebyshev Approximation of
* (1+2x)exp(x^2)erfc x in 0 <= x < INF", Mathematics of Computation, Vol. 36,
* No. 153, January 1981, pp. 249-253.
*
*/
float my_erfcxf (float x)
{
float a, d, e, m, p, q, r, s,