OpenCL 嵌套循环产生意外结果

OpenCL nested loop produces unexpected results

我对 OpenCL/C99 有点陌生,但不明白为什么下面的两个内核会给出不同的结果。 X 已初始化为零,但在每个外循环步骤都需要 "re-zeroing",否则会获得不正确的结果(见图表)。请注意,我在这里没有调用任何并行性;据我所知,这完全是连续的:

kernels.cl的内容:

#pragma OPENCL EXTENSION cl_khr_fp64: enable
#define __CL_ENABLE_EXCEPTIONS

__kernel void nested_sum_succeeds(
    int L,
    __global __read_only double* X,
    __global double* Y
)
{
    for (int k=0; k<=L; k++) {
        Y[k] = 0;
        for (int j=0; j<=k; j++) {
            Y[k] += X[j];
        }
    }
}

__kernel void nested_sum_fails(
    int L,
    __global __read_only double* X,
    __global double* Y
)
{
    for (int k=0; k<=L; k++) {
        // Y[k] = 0;
        for (int j=0; j<=k; j++) {
            Y[k] += X[j];
        }
    }
}

script.py的内容:

import numpy as np
import pyopencl as cl
import pyopencl.array as cl_array
import matplotlib.pyplot as plt

with open("./kernels.cl") as fp:
    prog_str = fp.read()
ctx = cl.create_some_context()
queue = cl.CommandQueue(ctx)
prog = cl.Program(ctx, prog_str).build()

L = 1000
X = np.linspace(0, 10, L)
X_dev = cl_array.to_device(queue, X)
Y_succeeds_dev = cl_array.to_device(queue, np.zeros(shape=X.shape, dtype=np.float64))
Y_fails_dev = cl_array.to_device(queue, np.zeros(shape=X.shape, dtype=np.float64))


nested_sum_succeeds = prog.nested_sum_succeeds
nested_sum_succeeds.set_scalar_arg_dtypes([
    np.int64,
    None,
    None,
])

nested_sum_succeeds(
    queue,
    (len(X),),
    None,
    L,
    X_dev.data,
    Y_succeeds_dev.data,
)


nested_sum_fails = prog.nested_sum_fails
nested_sum_fails.set_scalar_arg_dtypes([
    np.int64,
    None,
    None,
])

nested_sum_fails(
    queue,
    (len(X),),
    None,
    L,
    X_dev.data,
    Y_fails_dev.data,
)

np.allclose(Y_succeeds_dev.get(), Y_fails_dev.get()) #False

plt.ion()
plt.plot(Y_succeeds_dev.get())
plt.plot(Y_fails_dev.get())

结果:

Notice that I have not invoked any parallelism here;

是的,内核并行执行 运行 - 它计划在第一个维度中 运行 len(x) 工作项。将其更改为使用 1 个工作项进行处理后,一切正常:

nested_sum_succeeds(
    queue,
    (1,),
    None,
    L,
    X_dev.data,
    Y_succeeds_dev.data,
)

nested_sum_fails(
    queue,
    (1,),
    None,
    L,
    X_dev.data,
    Y_fails_dev.data,
)

然后 np.allclose(Y_succeeds_dev.get(), Y_fails_dev.get()) returns 正确。 您还可以从 nested_sum_succeeds 内核中删除该归零 Y[k] = 0;,因为它不需要。

此外,如果您想 运行 在其他设备中使用此内核,则需要进行一些小的修复,因为并非所有编译器都会接受第一个内核参数的类型在内核中 int 并且在 python 中被指定为 np.int64,它必须匹配内核中的内容,所以:

nested_sum_succeeds.set_scalar_arg_dtypes([
    np.int32,
    None,
    None,
])

nested_sum_fails.set_scalar_arg_dtypes([
    np.int32,
    None,
    None,
])

还有一件事适用于在其他设备上使用,我会删除 __read_only 访问限定符,它也不会在所有设备上编译。