OpenCL 多重指数减少

Question

我有两个数组 - 一个源数组

a

和一个目标数组

b

。 我想通过从

a

索引它们来设置

b

中的某些值。为此，我有另一个索引数组

i

，其大小与

b

相同，但不一定与

a

相同，并为

b

中的每个元素分配

a

中的位置。如果

i

将两个元素分配到同一位置，则必须将它们相加。即

b = [1,2,3,4,7,8]
i = [0,0,1,3,3,1]
a = [1,2,3,4,5,6,7,8]
operation(b, i, a) = [3,11,3,11,6,7,8]

我更喜欢在 GPU 上执行此操作，并且因为某些现有内核会在具有结果工作大小（即

a

的大小）的内核中执行此操作。如果值是浮点，是否可以轻松做到这一点（因为 OpenCL 中没有原子浮点运算）？

Answer 1

从“i”值的角度来看效率不高（正如您猜测的那样，只有原子才容易使用它，并且只有当您可以将浮点数缩放为整数然后重新缩放回来时），它被称为“分散”模式。

相反，使用“聚集”模式遍历“a”索引的目标位置，也可能自动从缓存获取帮助。像这样：

• 当前线程将计算 a[k]，其中 k=thread-id
• 检查整个数组并计算（当所有数组都是 N 大小时，复杂度为 O(N x N)）

以下内核在具有 100 万个元素的 RTX4070 上需要 3.8 秒：

 //size of a,b,i: 1024*1024
 kernel void run(global float * a, global float * b, global int * i)
 {
            int id = get_global_id(0);
            const int m = 1024*1024;
            for(int k=0;k<m;k++)
            {
                    a[id] += (i[k] == id) * b[k];
            }
 }

因为它还没有优化而且还很暴力。要优化，您可以尝试其他方法，例如平铺（使用本地内存）、分层碰撞检查加速结构。

您还可以将浮点数重新缩放为整数并使用原子来仅获得 15 毫秒（RTX4070），但不能保证具有稳定的运行时间。暴力版本虽然速度慢，但运行时间恒定。

仅本地内存平铺就应该提供至少 20-30 倍的加速，因为主内存带宽约为 500GB/s，而本地内存带宽约为数十 TB/s。

本地内存平铺

// 1024 * 1024 workitems, 1024 workitems per group
kernel void run(global float * a, global float * b, global int * i)
{
    int id = get_global_id(0);
    int localId = id % 1024;
    const int m = 1024*1024;
    local int cacheI[1024];
    local float cacheB[1024];
    float val = a[id];
    for(int k=0;k<m;k+=1024)
    {
            cacheI[localId] = i[k + localId];
            cacheB[localId] = b[k + localId];
            barrier(CLK_LOCAL_MEM_FENCE);
            for(int p=0;p<1024;p++)                            
                val += (cacheI[p] == id) * cacheB[p];
            barrier(CLK_LOCAL_MEM_FENCE);
    }
    a[id] = val;
}

1M 元素为 337 毫秒（在 0-1M-1 之间随机初始化）或来自缓存 I 和缓存 B 读数的 23 TB/s 带宽。

23 TB/s 在 RTX4070 上已经足够好了，因此寄存器平铺没有多大帮助。但未来的一些其他架构可能会受益于寄存器平铺和本地内存平铺的组合。

编译器可能会重复使用寄存器来重复读取相同的cacheI和cacheB值，因为23TB/s高于假设的每个SM单元300GB/s的理论带宽（其中46个仅约13TB/s）。

仅使用 256k 工作项进行注册平铺 + 本地内存平铺

  kernel void run(global float * a, global float * b, global int * i)
    {
        int id = get_global_id(0);
        int localId = id % 1024;
        const int m = 1024*1024;
        local int cacheI[1024];
        local float cacheB[1024];
        const int id1 = id*4;
        const int id2 = id1 + 1;
        const int id3 = id1 + 2;
        const int id4 = id1 + 3;
        float val = a[id1];
        float val2 = a[id2];
        float val3 = a[id3];
        float val4 = a[id4];
        for(int k=0;k<m;k+=1024)
        {
                cacheI[localId] = i[k + localId];
                cacheB[localId] = b[k + localId];
                barrier(CLK_LOCAL_MEM_FENCE);
                for(int p=0;p<1024;p++)                            
                {
                    const int registerI = cacheI[p];
                    const float registerB = cacheB[p];
                    val += (registerI == id1) * registerB;
                    val2 += (registerI == id2) * registerB;
                    val3 += (registerI == id3) * registerB;
                    val4 += (registerI == id4) * registerB;
                }
                barrier(CLK_LOCAL_MEM_FENCE);
        }
        a[id1] = val;
        a[id2] = val2;
        a[id3] = val3;
        a[id4] = val4;
    }

348 毫秒。只需使用局部平铺就足够了。