OpenCL/C pow(x,0.5)!= sqrt(x)

Dsc*_*oni 6 c opencl pyopencl

我想我在一些非常奇怪的边界情况下,可能有双精度问题,我想知道,发生了什么.

在OpenCL内核中我使用:

#pragma OPENCL EXTENSION cl_khr_fp64 : enable

__private int k = 2; // I need k to be an int, because I want to use as a counter
__private double s = 18;
__private double a = 1;

a = a/(double)k; // just to show, that I make in-place typecasting of k
a = k+1;
k = (int)a; //to show that I store k in a double buffer in an intermediate-step
if ((k-1)==2)
{
//    k = 3;
    s = pow(s/(double)(k-1),0.5);
}
Run Code Online (Sandbox Code Playgroud)

这导致我s = 2.999[...]6 然而,如果我取消注释该k=3行,我得到(在我看来)正确的结果s = 3.这是为什么?

作为附带信息:当我这样做时,不会发生同样的行为

s = sqrt(s/(double)(k-1))
Run Code Online (Sandbox Code Playgroud)

下面是pyopencl的完整,最小内核和主机代码

内核(Minima.cl):

#pragma OPENCL EXTENSION cl_khr_fp64 : enable

__kernel void init_z(__global double * buffer)
{
    __private int x = get_global_id(0);
    __private int y = get_global_id(1);
    //w,h
    __private int w_y = get_global_size(1);
    __private int address = x*w_y+y;
    //h,w
    __private double init = 3.0;
    buffer[address]=init;
}

__kernel void root(__global double * buffer)
{
    __private int x = get_global_id(0);
    __private int y = get_global_id(1);
    //w,h
    __private int w_y = get_global_size(1);
    __private int address = x*w_y+y;
    //h,w
    __private double value = 18;
    __private int k;
    __private double out;
    k = (int) buffer[address];
  //k = 3;  If this line is uncommented, the result will be exact.
    out = pow(value/(double)(k-1), 0.5);
    buffer[address] = out;
}
Run Code Online (Sandbox Code Playgroud)

主办:

import pyopencl as cl
import numpy as np


platform = cl.get_platforms()[0]
devs = platform.get_devices()
device1 = devs[1]
h_buffer = np.empty((10,10)).astype(np.float64)
mf = cl.mem_flags
ctx = cl.Context([device1])
Queue1 = cl.CommandQueue(ctx,properties=cl.command_queue_properties.PROFILING_ENABLE)
Queue2 = cl.CommandQueue(ctx,properties=cl.command_queue_properties.PROFILING_ENABLE)
mf = cl.mem_flags
m_dic = {0:mf.READ_ONLY,1:mf.WRITE_ONLY,2:mf.READ_WRITE}


fi = open('Minimal.cl', 'r')
fstr = "".join(fi.readlines())
prg = cl.Program(ctx, fstr).build()
knl = prg.init_z
knl.set_scalar_arg_dtypes([None,])
knl_root = prg.root
knl_root.set_scalar_arg_dtypes([None,])



def f():
    d_buffer =  cl.Buffer(ctx,m_dic[2], int(10 * 10  * 8))
    knl.set_args(d_buffer)
    knl_root.set_args(d_buffer)
    a = cl.enqueue_nd_range_kernel(Queue2,knl,(10,10),None)
    b = cl.enqueue_nd_range_kernel(Queue2,knl_root,(10,10),None, wait_for = [a,])
    cl.enqueue_copy(Queue1,h_buffer,d_buffer,wait_for=[b,])
    return h_buffer
a = f()
a[0,0] # Getting the result on the host.
Run Code Online (Sandbox Code Playgroud)

编辑:由于一些更不清楚我再次更新这个问题.我明白,相同输入的值powsqrt不必相同.我的问题是,为什么pow显示SAME输入的不同输出,这取决于我从哪里得到它.

二进制文件位于pastebin上: k_explicitk_read

printf("a%\n", out)导致0x1.8p+1k=3线路和0x1.7ffffffffffffp+1当它被注释掉时.

Dsc*_*oni 0

因此,经过大量讨论,我的代码似乎没有问题。该行为似乎依赖于硬件,到目前为止只能在 Nvidia 硬件上重现。我仍然对“为什么”和“如何”发生感兴趣,但在这种情况下,问题得到了解答。