如何以有效的方式找到数字的立方根？我认为可以使用Newton-Raphson方法,但我不知道如何以编程方式猜测初始解决方案以最小化迭代次数.
鉴于超过"接受的答案"的"链接腐烂",我将给出一个更加独立的答案,重点关注快速获得适合超线性迭代的初始猜测的主题.
同性恋者(Wayback链接)的"调查" 为各种起始值/迭代组合提供了一些时序比较(包括Newton和Halley方法).它的参考文献是W. Kahan的作品,"计算一个真正的立方根 " 和K. Turkowski的作品,"计算立方根 ".
metamarist用这个片段更新了W. Kahan的DEC-VAX时代比特摆弄技术,该片段"假设32位整数"并依赖于IEEE 754格式的双精度"以5位精度生成初始估计":
inline double cbrt_5d(double d)
{
const unsigned int B1 = 715094163;
double t = 0.0;
unsigned int* pt = (unsigned int*) &t;
unsigned int* px = (unsigned int*) &d;
pt[1]=px[1]/3+B1;
return t;
}
K. Turkowski的代码通过传统的2次幂缩放提供稍微更高的精度("大约6位")float fr,然后在区间[0.125,1.0]上对其立方根进行二次近似:
/* Compute seed with a quadratic qpproximation */
fr = (-0.46946116F * fr + 1.072302F) * fr + 0.3812513F;/* 0.5<=fr<1 */
并随后恢复两个指数(调整为三分之一).该指数/尾数提取和恢复使用的数学库调用来frexp和ldexp.
与其他立方根"种子"近似值的比较
为了理解这些立方根近似,我们需要将它们与其他可能的形式进行比较.首先是判断的标准:我们考虑区间[1/8,1]的近似值,并且我们使用最佳(最小化最大值)相对误差.
也就是说,如果f(x)是拟议的近似值x^{1/3},我们会发现它的相对误差:
error_rel = max | f(x)/x^(1/3) - 1 | on [1/8,1]
最简单的近似当然是在间隔上使用单个常数,并且在这种情况下的最佳相对误差是通过拾取f_0(x) = sqrt(2)/2端点处的值的几何平均值来实现的.这给出了1.27位的相对精度,这是牛顿迭代的一个快速但又脏的起点.
更好的近似将是最好的一次多项式:
f_1(x) = 0.6042181313*x + 0.4531635984
这给出了4.12位的相对精度,这是一个很大的改进,但是缺少了Kahan和Turkowski各自方法所承诺的5-6位相对精度.但它只是在球场,并且只使用一次乘法(和一次加法).
最后,如果我们允许自己划分而不是乘法怎么办？事实证明,通过一个除法和两个"加法",我们可以得到最好的线性分数函数:
f_M(x) = 1.4774329094 - 0.8414323527/(x+0.7387320679)
它给出了7.265位的相对精度.
乍一看,这似乎是一种有吸引力的方法,但旧的经验法则是将FP分割的成本视为三次FP乘法(并且主要忽略加法和减法).然而,对于当前的FPU设计,这是不现实的.虽然增加/减少乘法的相对成本已降低,但在大多数情况下降低到两倍甚至相等,但除法成本并没有下降,但往往高达乘法成本的7-10倍.因此,我们必须吝啬我们的分工.