pyopencl 乘法问题

Question

这似乎是一个奇怪但非常基本的问题。我尝试在 pyopencl 中做一个简单的操作。下面给出的是代码，如果我将我的位置与 exp(-f_sum/sigma2)/sigma2 相乘，我得到 0（即使我的两个位置和 sigma 都有非零值）但是当我添加值时，我得到了正确的结果。

kernelsource = """  __kernel void forceFinder(
const int N,
const int dim,
const float sigma,
const float resistant,
__global float* datacl,
__constant float* poscl,
__global float* res

)
{
int i = get_global_id(0);
float f_sum ;
int k;
float sigma2 = sigma * sigma;
float tempo;
if (i < N ) {

    f_sum = 0;

    for (k = 0; k < dim; k++)
    {

        f_sum += pown((poscl[k] - datacl[i * dim + k]), 2);

    }




     for (k = 0; k < dim; k++)
    { 

        res[i * dim + k] = (datacl[i * dim + k] - poscl[k]) * exp(-f_sum/sigma2)/sigma2;

    }
  }

}
"""     

如果我用“+”替换最后一个循环中的“*”而不是“*”，我会得到输出

Answer 1

我试图通过编写一个完整的示例来了解您的内核在做什么（见下文）。虽然我无法完全理解发生了什么，但如果我在我的机器上运行下面的代码，我会收到一个 10x10 矩阵，所有条目的值为 -0.0024。 如果您需要进一步的帮助，请提供完整的示例或更多信息。

另外：你可以通过使用 exp(x)*exp(x)=exp(x+x) 来摆脱你的第一个 for 循环吗？

import pyopencl as cl
import numpy as np

kernelsource = """  
__kernel void forceFinder(  const int dim,
                            const float sigma,
                            __global float* datacl,
                            __constant float* poscl,
                            __global float* res ){
    int i = get_global_id(0);
    float f_sum = 0;
    float sigma2 = sigma * sigma;
    for (int k = 0; k < dim; k++){
        f_sum += pown((poscl[k] - datacl[i * dim + k]), 2);
    }

    for (int k = 0; k < dim; k++){ 
        res[i * dim + k] = (datacl[i * dim + k] - poscl[k]) * exp(-f_sum/sigma2)/sigma2;
    }
}
"""

device = cl.get_platforms()[0].get_devices()[0]
context = cl.Context([device])
program = cl.Program(context, kernelsource).build()
queue = cl.CommandQueue(context)

sigma = 20
dim = 10
N = 5
poscl_local = np.ones(dim).astype(np.float32) * 2.
datacl_local = np.ones((N,dim)).astype(np.float32)
res_local = np.zeros(datacl_local.shape).astype(np.float32)

poscl_buf = cl.Buffer(context, cl.mem_flags.READ_ONLY | cl.mem_flags.COPY_HOST_PTR, hostbuf=poscl_local)
datacl_buf = cl.Buffer(context, cl.mem_flags.READ_ONLY | cl.mem_flags.COPY_HOST_PTR, hostbuf=datacl_local)
res_buf = cl.Buffer(context, cl.mem_flags.WRITE_ONLY, res_local.nbytes)
program.forceFinder(queue,(N,), None, np.int32(dim), np.float32(sigma),datacl_buf,poscl_buf,res_buf)
cl.enqueue_copy(queue, res_local, res_buf)

print("result: {}".format(res_local))