странное поведение функции powf ()

Question

Неожиданным образом powf производит странный вывод для нечетных базовых чисел, когда его тип равен int. Например, powf(-4,2) возвращает 16, но powf(-5,2) возвращает 24 !!!

После отслеживания корня неправильного вывода в длинных вычислениях я выяснил, что функция powf показывает странное поведение для нечетных чисел, когда тип вывода равен integer.

__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

EDIT:

После исследования в Matlab по предложению @Robert Crovella я обнаружил, что Command Window в Matlab показывает out4=25.000, в отличие от окна Variables, в котором отображается содержание out4 = 24.9999981.

Все должны быть очень осторожны, поскольку есть небольшая ошибка, связанная с выводом функции powf (24.9999981 вместо 25), которая может распространиться и стать проблемой при больших вычислениях

Answer 1

В качестве дополнения к ответу Роберта Кровеллы : CUDA является подмножеством C ++ и поэтому предоставляет перегруженные математические функции. В частности, он предлагает следующие четыре варианта pow():

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

Если вы проверите машинный код, сгенерированный для этих вариантов, с помощью cuobjdump --dump-sass, вы обнаружите, что используются четыре разные реализации. Как указал Роберт Кровелла, для конкретного случая возведения в квадрат лучше всего использовать умножение, но вы, безусловно, можете использовать 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

Answer 2

Я считаю, что это связано с неразумным использованием типов данных с feval.

Мне кажется, что feval преобразует тип возвращаемого значения в тот же тип, что и тип параметра. Это имеет смысл, поскольку возвращаемый тип извлекается из указателя на переданный аргумент для этого параметра.

обратите внимание, что powf принимает float параметров и возвращает float, а pow принимает double параметров и возвращает double. int количества не имеют отдельной функции (прототипа) в математическом API CUDA , поэтому, если вы используете их, они будут приводиться к типам с плавающей запятой и обратно.

Вот что я вижу в чистой 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 возвращает 24,999998 для (-5,2)
2. если мы конвертируем это в int, оно усекается до 24
3. если мы преобразуем это значение в double, а затем округлим до 3 десятичных разрядов, правильно округленный результат будет равен 25.000, как показано в выходных данных matlab

Предложения:

1. не делай этого
2. не использовать целочисленные типы с функциями с плавающей точкой (особенно приведение результата)
3. если вы хотите что-то возвести в квадрат, просто умножьте это на себя. Это определенно будет быстрее, чем использование powf(x, 2) и, возможно, будет более точным.

Если вы хотите знать, «почему CUDA powf(-5, 2) возвращает 24.999998?», Пожалуйста, задайте это в отдельном вопросе. Точность определена в руководстве по программированию , и я уверен, что это находится в пределах опубликованных границ ошибки.