memory-management - 写入全局或本地内存会使内核执行时间增加10000％

Question

我有以下OpenCL内核:

kernel void ndft(
    global float *re, global float *im, int num_values,
    global float *spectrum_re, global float *spectrum_im,
    global float *spectrum_abs,
    global float *sin_array, global float *cos_array,
    float sqrt_num_values_reciprocal)
{
    // MATH MAGIC - DISREGARD FROM HERE -----------

    float x;
    float y;
    float sum_re = 0;
    float sum_im = 0;

    size_t thread_id = get_global_id(0);
    //size_t local_id = get_local_id(0);

    // num_values = 24 (live environment), 48 (test)
    for (int i = 0; i < num_values; i++)
    {
        x = cos_array[thread_id * num_values + i] * sqrt_num_values_reciprocal;
        y = sin_array[thread_id * num_values + i] * sqrt_num_values_reciprocal;
        sum_re = sum_re + re[i] * x + im[i] * y;
        sum_im = sum_im - re[i] * y + x * im[i];
    }

    // MATH MAGIC DONE ----------------------------

    //spectrum_re[thread_id] = sum_re;
    //spectrum_im[thread_id] = sum_im;
    //spectrum_abs[thread_id] = hypot(sum_re, sum_im);
    float asdf = hypot(sum_re, sum_im); // this is just a dummy calculation
}

这样，执行时间约为15 us(工作组大小= 567，14个工作组，总共7938个线程)。

但是，当然，我需要某种方式来检索操作的结果，这是最后几行的意思(注释掉了)。一旦执行了这些内存操作中的一个操作(并且 spectrum_X是 global，如示例中的 local，就没有关系)，内核的执行时间将增加到〜1.4到1.5 ms。

我以为执行时间的增加是某种固定的开销，所以我只会累积更多的数据，以使由于这种影响而造成的相对时间损失最小化。但是，当我将线程数量加倍时(即数据量的两倍)，执行时间也会加倍(至2.8〜3.0 ms)。

我发现，即使我只取消注释了其中一行，也有相同的执行时间，就像我取消注释了所有三行一样。即使添加 if (thread_id == 0)并运行它，我的执行时间也相同。但是，这样太慢了(我的应用程序的上限大约是30 us)。当我在CPU上以普通C代码运行它时，它的性能甚至提高了约5倍。

现在，我显然做错了事，但是我不确定从哪里开始寻找解决方案。

当我评论talonmies的答案时，我还进行了以下操作:

从上面的代码，我使最后四行看起来像

//spectrum_re[thread_id] = sum_re;
//spectrum_im[thread_id] = sum_im;
spectrum_abs[thread_id] = hypot(sum_re, sum_im);
//float asdf = hypot(sum_re, sum_im);

如预期的那样，执行时间约为1.8毫秒。
我的系统生成的汇编代码是:

//
// Generated by NVIDIA NVVM Compiler
// Compiler built on Tue Apr 03 12:42:39 2012 (1333449759)
// Driver 
//

.version 3.0
.target sm_21, texmode_independent
.address_size 32


.entry ndft(
    .param .u32 .ptr .global .align 4 ndft_param_0,
    .param .u32 .ptr .global .align 4 ndft_param_1,
    .param .u32 ndft_param_2,
    .param .u32 .ptr .global .align 4 ndft_param_3,
    .param .u32 .ptr .global .align 4 ndft_param_4,
    .param .u32 .ptr .global .align 4 ndft_param_5,
    .param .u32 .ptr .global .align 4 ndft_param_6,
    .param .u32 .ptr .global .align 4 ndft_param_7,
    .param .f32 ndft_param_8
)
{
    .reg .f32   %f;
    .reg .pred  %p;
    .reg .s32   %r;


    ld.param.u32    %r3, [ndft_param_2];
    // inline asm
    mov.u32     %r18, %envreg3;
    // inline asm
    // inline asm
    mov.u32     %r19, %ntid.x;
    // inline asm
    // inline asm
    mov.u32     %r20, %ctaid.x;
    // inline asm
    // inline asm
    mov.u32     %r21, %tid.x;
    // inline asm
    add.s32     %r22, %r21, %r18;
    mad.lo.s32  %r11, %r20, %r19, %r22;
    setp.gt.s32     %p1, %r3, 0;
    @%p1 bra    BB0_2;

    mov.f32     %f46, 0f00000000;
    mov.f32     %f45, %f46;
    bra.uni     BB0_4;

BB0_2:
    ld.param.u32    %r38, [ndft_param_2];
    mul.lo.s32  %r27, %r38, %r11;
    shl.b32     %r28, %r27, 2;
    ld.param.u32    %r40, [ndft_param_6];
    add.s32     %r12, %r40, %r28;
    ld.param.u32    %r41, [ndft_param_7];
    add.s32     %r13, %r41, %r28;
    mov.f32     %f46, 0f00000000;
    mov.f32     %f45, %f46;
    mov.u32     %r43, 0;
    mov.u32     %r42, %r43;

BB0_3:
    add.s32     %r29, %r13, %r42;
    ld.global.f32   %f18, [%r29];
    ld.param.f32    %f44, [ndft_param_8];
    mul.f32     %f19, %f18, %f44;
    add.s32     %r30, %r12, %r42;
    ld.global.f32   %f20, [%r30];
    mul.f32     %f21, %f20, %f44;
    ld.param.u32    %r35, [ndft_param_0];
    add.s32     %r31, %r35, %r42;
    ld.global.f32   %f22, [%r31];
    fma.rn.f32  %f23, %f22, %f19, %f46;
    ld.param.u32    %r36, [ndft_param_1];
    add.s32     %r32, %r36, %r42;
    ld.global.f32   %f24, [%r32];
    fma.rn.f32  %f46, %f24, %f21, %f23;
    neg.f32     %f25, %f22;
    fma.rn.f32  %f26, %f25, %f21, %f45;
    fma.rn.f32  %f45, %f24, %f19, %f26;
    add.s32     %r42, %r42, 4;
    add.s32     %r43, %r43, 1;
    ld.param.u32    %r37, [ndft_param_2];
    setp.lt.s32     %p2, %r43, %r37;
    @%p2 bra    BB0_3;

BB0_4:
    // inline asm
    abs.f32     %f27, %f46;
    // inline asm
    // inline asm
    abs.f32     %f29, %f45;
    // inline asm
    setp.gt.f32     %p3, %f27, %f29;
    selp.f32    %f8, %f29, %f27, %p3;
    selp.f32    %f32, %f27, %f29, %p3;
    // inline asm
    abs.f32     %f31, %f32;
    // inline asm
    setp.gt.f32     %p4, %f31, 0f7E800000;
    mov.f32     %f47, %f32;
    @%p4 bra    BB0_6;

    mov.f32     %f48, %f8;
    bra.uni     BB0_7;

BB0_6:
    mov.f32     %f33, 0f3E800000;
    mul.rn.f32  %f10, %f8, %f33;
    mul.rn.f32  %f47, %f32, %f33;
    mov.f32     %f48, %f10;

BB0_7:
    mov.f32     %f13, %f48;
    // inline asm
    div.approx.f32  %f34, %f13, %f47;
    // inline asm
    mul.rn.f32  %f39, %f34, %f34;
    add.f32     %f38, %f39, 0f3F800000;
    // inline asm
    sqrt.approx.f32     %f37, %f38;      // <-- this is part of hypot()
    // inline asm
    mul.rn.f32  %f40, %f32, %f37;
    add.f32     %f41, %f32, %f8;
    setp.eq.f32     %p5, %f32, 0f00000000;
    selp.f32    %f42, %f41, %f40, %p5;
    setp.eq.f32     %p6, %f32, 0f7F800000;
    setp.eq.f32     %p7, %f8, 0f7F800000;
    or.pred     %p8, %p6, %p7;
    selp.f32    %f43, 0f7F800000, %f42, %p8;
    shl.b32     %r33, %r11, 2;
    ld.param.u32    %r39, [ndft_param_5];
    add.s32     %r34, %r39, %r33;
    st.global.f32   [%r34], %f43;    // <-- stores the hypot's result in spectrum_abs
    ret;
}

确实，我所有的计算操作都在那里-大量的添加/混合以及 sqrt函数的 hypot。从上面的asm代码中，我删除了第二行:
st.global.f32 [%r34], %f43;
这是实际将数据存储在全局数组 spectrum_abs中的行。然后，我使用 clCreateProgramWithBinary并使用修改后的asm代码文件作为输入。执行时间减少到20 us。

Answer 1

我想您会看到编译器优化的效果。

NVIDIA编译器非常积极地消除“死代码”，该“死代码”不直接参与对全局内存的写入。因此，在您的内核中，如果您不编写sum_re或sum_im，则编译器将优化整个计算循环(可能还有其他所有内容)，并将您的内核留给一个空内核，其中只包含no-op。您所看到的15微秒执行时间主要只是内核启动开销，而没有太多其他开销。当取消注释全局存储器写操作时，编译器将所有计算代码保留在原位，您会看到代码的真正执行时间。

因此，您可能应该问的真正问题是如何优化内核，以将其执行时间从当前所需的1.5毫秒减少到(非常雄心勃勃的)30微秒目标。

尽管对原始答案表示怀疑，但这是一个完整的repro案例，它支持断言这是与编译器相关的效果:

#include <iostream>
#include <OpenCL/opencl.h>

size_t source_size;
const char * source_str =
"kernel void ndft(                                                                     \n" \
"    global float *re, global float *im, int num_values,                               \n" \
"    global float *spectrum_re, global float *spectrum_im,                             \n" \
"    global float *spectrum_abs,                                                       \n" \
"    global float *sin_array, global float *cos_array,                                 \n" \
"    float sqrt_num_values_reciprocal)                                                 \n" \
"{                                                                                     \n" \
"    // MATH MAGIC - DISREGARD FROM HERE -----------                                   \n" \
"                                                                                      \n" \
"    float x;                                                                          \n" \
"    float y;                                                                          \n" \
"    float sum_re = 0;                                                                 \n" \
"    float sum_im = 0;                                                                 \n" \
"                                                                                      \n" \
"    size_t thread_id = get_global_id(0);                                              \n" \
"                                                                                      \n" \
"    for (int i = 0; i < num_values; i++)                                              \n" \
"    {                                                                                 \n" \
"        x = cos_array[thread_id * num_values + i] * sqrt_num_values_reciprocal;       \n" \
"        y = sin_array[thread_id * num_values + i] * sqrt_num_values_reciprocal;       \n" \
"        sum_re += re[i] * x + im[i] * y;                                              \n" \
"        sum_im -= re[i] * y + x * im[i];                                              \n" \
"    }                                                                                 \n" \
"                                                                                      \n" \
"    // MATH MAGIC DONE ----------------------------                                   \n" \
"                                                                                      \n" \
"    //spectrum_re[thread_id] = sum_re;                                                \n" \
"    //spectrum_im[thread_id] = sum_im;                                                \n" \
"    //spectrum_abs[thread_id] = hypot(sum_re, sum_im);                                \n" \
"}                                                                                     \n";

int main(void)
{
    int err;

    cl_device_id device_id;  
    clGetDeviceIDs(NULL, CL_DEVICE_TYPE_GPU, 1, &device_id, NULL);
    cl_context context = clCreateContext(0, 1, &device_id, NULL, NULL, &err);
    cl_program program = clCreateProgramWithSource(context, 1, (const char **)&source_str, (const size_t *)&source_size, &err);

    err = clBuildProgram(program, 1, &device_id, NULL, NULL, NULL);

    cl_uint program_num_devices;
    clGetProgramInfo(program, CL_PROGRAM_NUM_DEVICES, sizeof(cl_uint), &program_num_devices, NULL);

    size_t * binaries_sizes = new size_t[program_num_devices];
    clGetProgramInfo( program, CL_PROGRAM_BINARY_SIZES, program_num_devices*sizeof(size_t), binaries_sizes, NULL);

    char **binaries = new char*[program_num_devices];
    for (size_t i = 0; i < program_num_devices; i++)
        binaries[i] = new char[binaries_sizes[i]+1];

    clGetProgramInfo(program, CL_PROGRAM_BINARIES, program_num_devices*sizeof(size_t), binaries, NULL);
    for (size_t i = 0; i < program_num_devices; i++)
    {
        binaries[i][binaries_sizes[i]] = '\0';
        std::cout << "Program " << i << ":" << std::endl;
        std::cout << binaries[i];
    }
    return 0;
}

编译并运行后，它将从OpenCL运行时发出以下PTX代码:

Program 0:
bplist00?^clBinaryDriver\clBinaryData_clBinaryVersionWCLH 1.0O!.version 1.5
.target sm_12
.target texmode_independent

.reg .b32  r<126>; /* define r0..125 */ 
.reg .b64  x<126>; /* define r0..125 */ 
.reg .b32  f<128>; /* define f0..127 */ 
.reg .pred p<32>; /* define p0..31 */ 
.reg .u32  sp;

.reg .b8   wb0,wb1,wb2,wb3; /* 8-bit write buffer */ 
.reg .b16  ws0,ws1,ws2,ws3; /* 16-bit write buffer */ 
.reg .b32  tb0,tb1,tb2,tb3; /* read tex buffer */ 
.reg .b64  vl0,vl1; /* 64-bit vector buffer */ 
.reg .b16  cvt16_0,cvt16_1; /* tmps for conversions */ 


.const .align 1 .b8 ndft_gid_base[52];
.local .align 16 .b8 ndft_stack[8];
.entry ndft(
    .param.b32 ndft_0  /* re */,
    .param.b32 ndft_1  /* im */,
    .param.b32 ndft_2  /* num_values */,
    .param.b32 ndft_3  /* spectrum_re */,
    .param.b32 ndft_4  /* spectrum_im */,
    .param.b32 ndft_5  /* spectrum_abs */,
    .param.b32 ndft_6  /* sin_array */,
    .param.b32 ndft_7  /* cos_array */,
    .param.f32 ndft_8  /* sqrt_num_values_reciprocal */
) {
    mov.u32 sp, ndft_stack;
    mov.u32 r0, 4294967295;
    ld.param.u32 r1, [ndft_2 + 0];
LBB1_1:
    add.u32 r0, r0, 1;
    setp.lt.s32 p0, r0, r1;
    @p0 bra LBB1_1;
LBB1_2:
    ret;
}

IE。一个不包含任何计算循环的内核 stub 。当内核的最后三行中的三个全局内存写入未注释时，它发出以下消息:

Program 0:
S.version 1.5inaryDriver\clBinaryData_clBinaryVersionWCLH 1.0O
.target sm_12
.target texmode_independent

.reg .b32  r<126>; /* define r0..125 */ 
.reg .b64  x<126>; /* define r0..125 */ 
.reg .b32  f<128>; /* define f0..127 */ 
.reg .pred p<32>; /* define p0..31 */ 
.reg .u32  sp;

.reg .b8   wb0,wb1,wb2,wb3; /* 8-bit write buffer */ 
.reg .b16  ws0,ws1,ws2,ws3; /* 16-bit write buffer */ 
.reg .b32  tb0,tb1,tb2,tb3; /* read tex buffer */ 
.reg .b64  vl0,vl1; /* 64-bit vector buffer */ 
.reg .b16  cvt16_0,cvt16_1; /* tmps for conversions */ 


.const .align 1 .b8 ndft_gid_base[52];
.local .align 16 .b8 ndft_stack[8];
.entry ndft(
    .param.b32 ndft_0  /* re */,
    .param.b32 ndft_1  /* im */,
    .param.b32 ndft_2  /* num_values */,
    .param.b32 ndft_3  /* spectrum_re */,
    .param.b32 ndft_4  /* spectrum_im */,
    .param.b32 ndft_5  /* spectrum_abs */,
    .param.b32 ndft_6  /* sin_array */,
    .param.b32 ndft_7  /* cos_array */,
    .param.f32 ndft_8  /* sqrt_num_values_reciprocal */
) {
    mov.u32 sp, ndft_stack;
    cvt.u32.u16 r0, %tid.x;
    cvt.u32.u16 r1, %ntid.x;
    cvt.u32.u16 r2, %ctaid.x;
    mad24.lo.u32 r0, r2, r1, r0;
    mov.u32 r1, 0;
    shl.b32 r2, r1, 2;
    mov.u32 r3, ndft_gid_base;
    add.u32 r2, r2, r3;
    ld.const.u32 r2, [r2 + 40];
    add.u32 r0, r0, r2;
    ld.param.u32 r2, [ndft_2 + 0];
    mul.lo.u32 r3, r0, r2;
    shl.b32 r3, r3, 2;
    mov.f32 f0, 0f00000000 /* 0.000000e+00 */;
    ld.param.f32 f1, [ndft_8 + 0];
    ld.param.u32 r4, [ndft_7 + 0];
    ld.param.u32 r5, [ndft_6 + 0];
    ld.param.u32 r6, [ndft_5 + 0];
    ld.param.u32 r7, [ndft_4 + 0];
    ld.param.u32 r8, [ndft_3 + 0];
    ld.param.u32 r9, [ndft_1 + 0];
    ld.param.u32 r10, [ndft_0 + 0];
    mov.u32 r11, r1;
    mov.f32 f2, f0;
LBB1_1:
    setp.ge.s32 p0, r11, r2;
    @!p0 bra    LBB1_7;
LBB1_2:
    shl.b32 r1, r0, 2;
    add.u32 r2, r8, r1;
    st.global.f32 [r2+0], f0;
    add.u32 r1, r7, r1;
    st.global.f32 [r1+0], f2;
    abs.f32 f1, f2;
    abs.f32 f0, f0;
    setp.gt.f32 p0, f0, f1;
    selp.f32 f2, f0, f1, p0;
    abs.f32 f3, f2;
    mov.f32 f4, 0f7E800000 /* 8.507059e+37 */;
    setp.gt.f32 p1, f3, f4;
    selp.f32 f0, f1, f0, p0;
    shl.b32 r0, r0, 2;
    add.u32 r0, r6, r0;
    @!p1 bra    LBB1_8;
LBB1_3:
    mul.rn.f32 f3, f2, 0f3E800000 /* 2.500000e-01 */;
    mul.rn.f32 f1, f0, 0f3E800000 /* 2.500000e-01 */;
LBB1_4:
    mov.f32 f4, 0f00000000 /* 0.000000e+00 */;
    setp.eq.f32 p0, f2, f4;
    @!p0 bra    LBB1_9;
LBB1_5:
    add.f32 f1, f2, f0;
LBB1_6:
    mov.f32 f3, 0f7F800000 /* inf */;
    setp.eq.f32 p0, f0, f3;
    setp.eq.f32 p1, f2, f3;
    or.pred p0, p1, p0;
    selp.f32 f0, f3, f1, p0;
    st.global.f32 [r0+0], f0;
    ret;
LBB1_7:
    add.u32 r12, r3, r1;
    add.u32 r13, r4, r12;
    ld.global.f32 f3, [r13+0];
    mul.rn.f32 f3, f3, f1;
    add.u32 r13, r9, r1;
    ld.global.f32 f4, [r13+0];
    mul.rn.f32 f5, f3, f4;
    add.u32 r12, r5, r12;
    ld.global.f32 f6, [r12+0];
    mul.rn.f32 f6, f6, f1;
    add.u32 r12, r10, r1;
    ld.global.f32 f7, [r12+0];
    mul.rn.f32 f8, f7, f6;
    add.f32 f5, f8, f5;
    sub.f32 f2, f2, f5;
    mul.rn.f32 f4, f4, f6;
    mul.rn.f32 f3, f7, f3;
    add.f32 f3, f3, f4;
    add.f32 f0, f0, f3;
    add.u32 r11, r11, 1;
    add.u32 r1, r1, 4;
    bra LBB1_1;
LBB1_8:
    mov.f32 f1, f0;
    mov.f32 f3, f2;
    bra LBB1_4;
LBB1_9:
    div.approx.f32 f1, f1, f3;
    mul.rn.f32 f1, f1, f1;
    add.f32 f1, f1, 0f3F800000 /* 1.000000e+00 */;
    sqrt.approx.ftz.f32 f1, f1;
    mul.rn.f32 f1, f2, f1;
    bra LBB1_6;
}

我认为这是无可辩驳的证据，证明是编译器优化导致了运行时差异，并且仅取决于内核代码中是否包含内存写操作。

我猜最后一个问题就变成了为什么它这么慢(不管关于它是否是由编译器优化引起的争论)。您所看到的1.5毫秒运行时是对代码性能的真实反射(reflect)，真正的问题是原因。根据我对内核代码的阅读，答案似乎在于内存访问模式，这对于GPU来说是非常可怕的。在计算循环内，您有两个全局内存读取，步幅非常大，如下所示:

x = cos_array[thread_id * num_values + i] * sqrt_num_values_reciprocal;

根据您的代码中的注释 num_values是24或48。这意味着内存读取不可能合并，并且Fermi GPU上的L1缓存也无济于事。这将对内存带宽利用率产生巨大的负面影响，并使代码非常缓慢。如果您坚持输入数据的排序，那么一种更快的解决方案是使用warp来计算一个输出(因此将warp范围内的约简化为最终和)。这会将读取步幅从24或48减少到1，并合并从那两个大输入数组读取的全局内存读取。

在循环内部，还重复读取了 re和 im的24或48个元素的全局内存:

    sum_re += re[i] * x + im[i] * y;
    sum_im -= re[i] * y + x * im[i];

这是不必要的，并且浪费了大量的全局内存带宽或缓存效率(GPU的寄存器不足，无法让编译器将每个数组的全部内容保存在寄存器中)。最好让每个工作组一次将这两个数组读入 __local内存数组，并在计算循环内使用本地内存副本。如果让每个工作组进行多次计算而不是一次计算，则可以节省大量的全局内存带宽并分摊初始读取，直到几乎可用为止。