powf() 函数的奇怪行为

strange behavior of powf() function

以一种意想不到的方式,powf 在其类型为 int 时为奇数基数生成了一个奇怪的输出。例如 powf(-4,2)returns 16powf(-5,2) returns24 !!!

在长时间的计算中追踪错误输出的根源后,我发现当输出类型为 integer.

时,powf 函数对奇数显示奇怪的行为
__global__ void intFoo( int* a) 
{
    *a = powf(*a, 2);
}
__global__ void doubleFoo( double* a) 
{
    *a = powf(*a, 2);
}

我可以在 Matlab 中调用这个内核(例如):

!nvcc -ptx test.cu 
k1 = parallel.gpu.CUDAKernel('test.ptx', 'test.cu', 'intFoo');
k2 = parallel.gpu.CUDAKernel('test.ptx', 'test.cu', 'doubleFoo');
out1 = feval(k1, -4)
out2 = feval(k1, -5)
out3 = feval(k2, -4)
out4 = feval(k2, -5)

结果:

out1 = 16
out2 = 24 //This hasn't to be 25 !!??
out3 = 16
out4 = 25.000

编辑:

根据@Robert Crovella 的建议在 Matlab 中进行调查后,我发现 Matlab 中的命令 Window 显示 out4=25.000 而变量 Window 显示 [=22 的内容=].

每个人都应该非常小心,因为有一个与 powf 函数(24.9999981 而不是 25)的输出相关的小错误可能会传播并成为大问题计算

我认为这是由于 feval.

对数据类型的不当使用造成的

在我看来,feval 将 return 类型转换为与参数类型相同的类型。这是有道理的,因为 return 类型是从指向该参数的传递参数的指针中提取的。

注意 powf 接受 float 参数和 return 一个 floatpow 接受 double 参数和 return 是 doubleint 数量在 the CUDA math API 中没有单独的函数(原型),因此如果您使用它们,它们将被转换为浮点类型。

这是我在纯 CUDA C++ 中看到的内容:

$ cat t32.cu
#include <math.h>
#include <stdio.h>

__global__ void Foo( int a, double b)
{
            float res = powf((float)a, 2);
            printf("powf_int: %d, %d, %f\n", a, (int)res, res);
            res = powf((float)b, 2);
            printf("powf_double: %f, %f, %f\n", b, (double)res, res);
            double dres = pow((double)a, 2);
            printf("pow_int: %d, %d, %f\n", a, (int)dres, dres);
            dres = pow((double)b, 2);
            printf("pow_double: %f, %f, %f\n", b, (double)dres, dres);
}

int main(){

        Foo<<<1,1>>>(-5, -5);
        cudaDeviceSynchronize();
}
$ nvcc -o t32 t32.cu
$ cuda-memcheck ./t32
========= CUDA-MEMCHECK
powf_int: -5, 24, 24.999998
powf_double: -5.000000, 24.999998, 24.999998
pow_int: -5, 25, 25.000000
pow_double: -5.000000, 25.000000, 25.000000
========= ERROR SUMMARY: 0 errors
$

注意:

  1. CUDA powf returns 24.999998 对于 (-5,2)
  2. 如果我们将其转换为 int,它会被截断为 24
  3. 如果我们将其转换为 double 然后四舍五入到小数点后 3 位,正确四舍五入的结果将是 25.000,正如您的 matlab 输出中显示的那样

建议:

  1. 不要这样做
  2. 不要将整数类型与浮点函数一起使用(尤其是转换结果)
  3. 如果你想对某物进行平方,只需将其与自身相乘即可。它肯定会比使用 powf(x, 2) 更快,也可能更准确。

如果您想知道 "why does CUDA powf(-5, 2) return 24.999998?",请在单独的问题中提问。准确度在 programming manual 中定义,我有理由相信这在发布的误差范围内。

作为 Robert Crovella 的附录:CUDA 是 C++ 的一个子集,因此提供了重载的数学函数。特别是它提供了 pow() 的以下四种变体:

float pow (float, int); 
double pow (double, int); 
float pow (float, float); 
double pow (double, double);

如果您使用 cuobjdump --dump-sass 检查为这些变体生成的机器代码,您会发现使用了四种不同的实现。正如 Robert Crovella 指出的那样,对于平方的特殊情况,最好只使用乘法,但如果你愿意,你当然可以使用 pow() ,如以下代码所示(为简洁起见省略了错误检查):

#include <cmath>
#include <cstdlib>
#include <cstdio>

__global__ void kernel (int ib, float fa, float fb, double da, double db)
{
    printf ("pow_float_int     = %15.8e\n", pow (fa, ib));
    printf ("pow_float_float   = %15.8e\n", pow (fa, fb));
    printf ("pow_double_int    = %23.16e\n", pow (da, ib));
    printf ("pow_double_double = %23.16e\n", pow (da, db));
}

int main (void)
{
    int ia = -5, ib = 2;
    float fa = ia, fb = ib;
    double da = ia, db = ib;

    kernel<<<1,1>>>(ib, fa, fb, da, db);
    cudaDeviceSynchronize();
    return EXIT_SUCCESS;
}

上述程序的输出应如下所示:

pow_float_int     =  2.50000000e+01
pow_float_float   =  2.49999981e+01
pow_double_int    =  2.5000000000000000e+01
pow_double_double =  2.5000000000000000e+01