Task05 Александр Исаенков

src/kernels/cl/radix_sort_01_local_counting.cl

@@ -5,14 +5,55 @@
 #include "helpers/rassert.cl"
 #include "../defines.h"
 
-__attribute__((reqd_work_group_size(1, 1, 1)))
+__attribute__((reqd_work_group_size(GROUP_SIZE, 1, 1)))
 __kernel void radix_sort_01_local_counting(
-    // это лишь шаблон! смело меняйте аргументы и используемые буфера! можете сделать даже больше кернелов, если это вызовет затруднения - смело спрашивайте в чате
-    // НЕ ПОДСТРАИВАЙТЕСЬ ПОД СИСТЕМУ! СВЕРНИТЕ С РЕЛЬС!! БУНТ!!! АНТИХАЙП!11!!1
     __global const uint* buffer1,
     __global       uint* buffer2,
     unsigned int a1,
     unsigned int a2)
 {
-    // TODO
+    const size_t elems_per_block = 2u * GROUP_SIZE;
+    const size_t group_id = get_group_id(0);
+    const size_t local_id = get_local_id(0);
+    const size_t base = group_id * elems_per_block;
+
+    uint private_hist[4];
+    private_hist[0] = 0u; 
+    private_hist[1] = 0u; 
+    private_hist[2] = 0u; 
+    private_hist[3] = 0u;
+
+    const size_t idx0 = base + local_id;
+    const size_t idx1 = base + GROUP_SIZE + local_id;
+    if (idx0 < a1) {
+        uint digit = (buffer1[idx0] >> a2) & 3u;
+        ++private_hist[digit];
+    }
+    if (idx1 < a1) {
+        uint digit = (buffer1[idx1] >> a2) & 3u;
+        ++private_hist[digit];
+    }
+
+    const size_t local_size = get_local_size(0);
+    __local uint reduction_buffer[GROUP_SIZE];
+    barrier(CLK_LOCAL_MEM_FENCE);
+    for (uint bits = 0; bits < 4; ++bits)
+    {
+        reduction_buffer[local_id] = private_hist[bits];
+        barrier(CLK_LOCAL_MEM_FENCE);
+
+        for (uint stride = local_size / 2; stride > 0; stride /= 2)
+        {
+            if (local_id < stride)
+            {
+                reduction_buffer[local_id] += reduction_buffer[local_id + stride];
+            }
+            barrier(CLK_LOCAL_MEM_FENCE);
+        }
+
+        if (local_id == 0u)
+        {
+            buffer2[4 * group_id + bits] = reduction_buffer[0];
+        }
+    }
 }

src/kernels/cl/radix_sort_02_global_prefixes_scan_sum_reduction.cl

@@ -5,13 +5,60 @@
 #include "helpers/rassert.cl"
 #include "../defines.h"
 
-__attribute__((reqd_work_group_size(1, 1, 1)))
+__attribute__((reqd_work_group_size(GROUP_SIZE, 1, 1)))
 __kernel void radix_sort_02_global_prefixes_scan_sum_reduction(
-    // это лишь шаблон! смело меняйте аргументы и используемые буфера! можете сделать даже больше кернелов, если это вызовет затруднения - смело спрашивайте в чате
-    // НЕ ПОДСТРАИВАЙТЕСЬ ПОД СИСТЕМУ! СВЕРНИТЕ С РЕЛЬС!! БУНТ!!! АНТИХАЙП!11!!1
     __global const uint* buffer1,
     __global       uint* buffer2,
-    unsigned int a1)
-{
-    // TODO
-}
+    unsigned int a1              
+){
+    const size_t n = a1;
+
+    const size_t local_id = get_local_id(0);
+    const size_t group_id = get_group_id(0);
+
+    const size_t elems_per_block = 2 * GROUP_SIZE;
+    const size_t base = group_id * elems_per_block;
+
+    const size_t idx0 = base + local_id;
+    const size_t idx1 = base + GROUP_SIZE + local_id;
+
+    uint v0 = 0;
+    uint v1 = 0;
+    if (idx0 < n) {
+        v0 = buffer1[idx0];
+    }
+    if (idx1 < n) {
+        v1 = buffer1[idx1];
+    }
+
+    __local uint s[2 * GROUP_SIZE];
+    s[local_id] = v0;
+    s[GROUP_SIZE + local_id] = v1;
+    barrier(CLK_LOCAL_MEM_FENCE);
+
+    // upsweep
+    for (size_t stride = 1; stride < elems_per_block; stride <<= 1) {
+        size_t index = ((local_id + 1) * stride * 2) - 1;
+        if (index < elems_per_block) {
+            s[index] += s[index - stride];
+        }
+        barrier(CLK_LOCAL_MEM_FENCE);
+    }
+
+    if (local_id == 0) {
+        buffer2[group_id] = s[elems_per_block - 1];
+        s[elems_per_block - 1] = 0;
+    }
+    barrier(CLK_LOCAL_MEM_FENCE);
+
+    // downsweep
+    for (size_t stride = elems_per_block >> 1; stride >= 1; stride >>= 1) {
+        size_t index = ((local_id + 1) * stride * 2) - 1;
+        if (index < elems_per_block) {
+            uint t = s[index - stride];
+            s[index - stride] = s[index];
+            s[index] += t;
+        }
+        barrier(CLK_LOCAL_MEM_FENCE);
+    }
+}

src/kernels/cl/radix_sort_03_global_prefixes_scan_accumulation.cl

@@ -5,14 +5,25 @@
 #include "helpers/rassert.cl"
 #include "../defines.h"
 
-__attribute__((reqd_work_group_size(1, 1, 1)))
+__attribute__((reqd_work_group_size(GROUP_SIZE, 1, 1)))
 __kernel void radix_sort_03_global_prefixes_scan_accumulation(
-    // это лишь шаблон! смело меняйте аргументы и используемые буфера! можете сделать даже больше кернелов, если это вызовет затруднения - смело спрашивайте в чате
-    // НЕ ПОДСТРАИВАЙТЕСЬ ПОД СИСТЕМУ! СВЕРНИТЕ С РЕЛЬС!! БУНТ!!! АНТИХАЙП!11!!1
     __global const uint* buffer1,
     __global       uint* buffer2,
-    unsigned int a1,
-    unsigned int a2)
+    unsigned int a1,             
+    unsigned int a2)             
 {
-    // TODO
-}
+    const uint n    = a1;
+
+    const size_t gid = get_global_id(0);
+    if (gid >= n) {
+        return;
+    }
+
+    const size_t elems_per_block = 2 * GROUP_SIZE;
+    const size_t block_id = gid / elems_per_block;
+
+    if (block_id == 0) {
+        return;
+    }
+    buffer2[gid] += buffer1[block_id - 1];
+}

src/kernels/cl/radix_sort_04_scatter.cl

@@ -5,15 +5,71 @@
 #include "helpers/rassert.cl"
 #include "../defines.h"
 
-__attribute__((reqd_work_group_size(1, 1, 1)))
+__attribute__((reqd_work_group_size(GROUP_SIZE, 1, 1)))
 __kernel void radix_sort_04_scatter(
-    // это лишь шаблон! смело меняйте аргументы и используемые буфера! можете сделать даже больше кернелов, если это вызовет затруднения - смело спрашивайте в чате
-    // НЕ ПОДСТРАИВАЙТЕСЬ ПОД СИСТЕМУ! СВЕРНИТЕ С РЕЛЬС!! БУНТ!!! АНТИХАЙП!11!!1
     __global const uint* buffer1,
     __global const uint* buffer2,
-                   uint* buffer3,
+    __global       uint* buffer3,
     unsigned int a1,
     unsigned int a2)
 {
-    // TODO
+    const size_t elems_per_block = 2u * GROUP_SIZE;
+    const size_t group_id = get_group_id(0);
+    const size_t local_id = get_local_id(0);
+    const size_t base_idx = group_id * elems_per_block;
+
+    const size_t idx0 = base_idx + local_id;
+    const size_t idx1 = base_idx + GROUP_SIZE + local_id;
+    const uint v0_valid = (idx0 < a1);
+    const uint v1_valid = (idx1 < a1);
+    const uint v0 = v0_valid ? buffer1[idx0] : 0u;
+    const uint v1 = v1_valid ? buffer1[idx1] : 0u;
+    const uint d0 = (v0 >> a2) & 3u;
+    const uint d1 = (v1 >> a2) & 3u;
+
+    uint base[4];
+    base[0] = buffer2[4 * group_id + 0];
+    base[1] = buffer2[4 * group_id + 1];
+    base[2] = buffer2[4 * group_id + 2];
+    base[3] = buffer2[4 * group_id + 3];
+
+    __local uint s[2 * GROUP_SIZE];
+
+    for (uint b = 0u; b < 4u; ++b) {
+        s[local_id] = (v0_valid && d0 == b) ? 1u : 0u;
+        s[local_id + GROUP_SIZE] = (v1_valid && d1 == b) ? 1u : 0u;
+        barrier(CLK_LOCAL_MEM_FENCE);
+
+        // upsweep
+        for (size_t stride = 1; stride < elems_per_block; stride <<= 1) {
+            size_t index = ((local_id + 1) * stride * 2) - 1;
+            if (index < elems_per_block) {
+                s[index] += s[index - stride];
+            }
+            barrier(CLK_LOCAL_MEM_FENCE);
+        }
+        if (local_id == 0) {
+            s[elems_per_block - 1] = 0;
+        }
+        barrier(CLK_LOCAL_MEM_FENCE);
+
+        // downsweep
+        for (size_t stride = elems_per_block >> 1; stride >= 1; stride >>= 1) {
+            size_t index = ((local_id + 1) * stride * 2) - 1;
+            if (index < elems_per_block) {
+                uint t = s[index - stride];
+                s[index - stride] = s[index];
+                s[index] += t;
+            }
+            barrier(CLK_LOCAL_MEM_FENCE);
+        }
+
+        if (v0_valid && d0 == b) {
+            buffer3[base[b] + s[local_id]] = v0;
+        }
+        if (v1_valid && d1 == b) {
+            buffer3[base[b] + s[local_id + GROUP_SIZE]] = v1;
+        }
+        barrier(CLK_LOCAL_MEM_FENCE);
+    }
 }

src/main_radix_sort.cpp

@@ -50,6 +50,7 @@ void run(int argc, char** argv)
 
     FastRandom r;
 
+    //int n = 10000;
     int n = 100*1000*1000; // TODO при отладке используйте минимальное n (например n=5 или n=10) при котором воспроизводится бага
     int max_value = std::numeric_limits<int>::max(); // TODO при отладке используйте минимальное max_value (например max_value=8) при котором воспроизводится бага
     std::vector<unsigned int> as(n, 0);
@@ -85,53 +86,59 @@ void run(int argc, char** argv)
         std::cout << "CPU std::sort effective RAM bandwidth: " << memory_size_gb / t.elapsed() << " GB/s (" << n / 1000 / 1000 / t.elapsed() << " uint millions/s)" << std::endl;
     }
 
+    const size_t elems_per_block = 2u * GROUP_SIZE;
+    const size_t num_blocks = (n + elems_per_block - 1) / elems_per_block;
+
     // Аллоцируем буферы в VRAM
     gpu::gpu_mem_32u input_gpu(n);
-    gpu::gpu_mem_32u buffer1_gpu(n), buffer2_gpu(n), buffer3_gpu(n), buffer4_gpu(n); // TODO это просто шаблонка, можете переименовать эти буферы, сделать другого размера/типа, удалить часть, добавить новые
+    gpu::gpu_mem_32u sort_src_gpu(n), sort_dst_gpu(n);
+    gpu::gpu_mem_32u hist_gpu(num_blocks * 4u);
+    gpu::gpu_mem_32u prefix_gpu(num_blocks * 4u);
     gpu::gpu_mem_32u buffer_output_gpu(n);
 
     // Прогружаем входные данные по PCI-E шине: CPU RAM -> GPU VRAM
     input_gpu.writeN(as.data(), n);
-    // Советую занулить (или еще лучше - заполнить какой-то уникальной константой, например 255) все буферы
-    // В некоторых случаях это ускоряет отладку, но обратите внимание, что fill реализован через копию множества нулей по PCI-E, то есть он очень медленный
-    // Если вам нужно занулять буферы в процессе вычислений - используйте кернел который это сделает (см. кернел fill_buffer_with_zeros)
-    buffer1_gpu.fill(255);
-    buffer2_gpu.fill(255);
-    buffer3_gpu.fill(255);
-    buffer4_gpu.fill(255);
-    buffer_output_gpu.fill(255);
+    input_gpu.copyToN(sort_src_gpu, n);
+
+    std::vector<unsigned int> prefix_cpu(num_blocks * 4u);
 
     // Запускаем кернел (несколько раз и с замером времени выполнения)
     std::vector<double> times;
-    for (int iter = 0; iter < 10; ++iter) { // TODO при отладке запускайте одну итерацию
+    for (int iter = 0; iter < 10; ++iter) {
         timer t;
-
-        // Запускаем кернел, с указанием размера рабочего пространства и передачей всех аргументов
-        // Если хотите - можете удалить ветвление здесь и оставить только тот код который соответствует вашему выбору API
         if (context.type() == gpu::Context::TypeOpenCL) {
-            // TODO
-            throw std::runtime_error(CODE_IS_NOT_IMPLEMENTED);
-            // ocl_fillBufferWithZeros.exec();
-            // ocl_radixSort01LocalCounting.exec();
-            // ocl_radixSort02GlobalPrefixesScanSumReduction.exec();
-            // ocl_radixSort03GlobalPrefixesScanAccumulation.exec();
-            // ocl_radixSort04Scatter.exec();
-        } else if (context.type() == gpu::Context::TypeCUDA) {
-            // TODO
-            throw std::runtime_error(CODE_IS_NOT_IMPLEMENTED);
-            // cuda::fill_buffer_with_zeros();
-            // cuda::radix_sort_01_local_counting();
-            // cuda::radix_sort_02_global_prefixes_scan_sum_reduction();
-            // cuda::radix_sort_03_global_prefixes_scan_accumulation();
-            // cuda::radix_sort_04_scatter();
-        } else if (context.type() == gpu::Context::TypeVulkan) {
-            // TODO
-            throw std::runtime_error(CODE_IS_NOT_IMPLEMENTED);
-            // vk_fillBufferWithZeros.exec();
-            // vk_radixSort01LocalCounting.exec();
-            // vk_radixSort02GlobalPrefixesScanSumReduction.exec();
-            // vk_radixSort03GlobalPrefixesScanAccumulation.exec();
-            // vk_radixSort04Scatter.exec();
+            input_gpu.copyToN(sort_src_gpu, n);
+            gpu::gpu_mem_32u* values_src = &sort_src_gpu;
+            gpu::gpu_mem_32u* values_dst = &sort_dst_gpu;
+
+            for (unsigned int bit_offset = 0u; bit_offset < 32u; bit_offset += 2u) {
+                const size_t global_size = num_blocks * GROUP_SIZE;
+                gpu::WorkSize ws(static_cast<size_t>(GROUP_SIZE), global_size);
+                ocl_radixSort01LocalCounting.exec(ws, *values_src, hist_gpu, static_cast<unsigned int>(n), bit_offset);
+
+                std::vector<unsigned int> hist_cpu = hist_gpu.readVector();
+                unsigned int cumulative[4] = {0};
+                for (unsigned int bits = 0u; bits < 4u; ++bits) {
+                    unsigned int sum = 0;
+                    for (size_t i = 0; i < num_blocks; ++i) {
+                        prefix_cpu[i * 4u + bits] = sum;
+                        sum += hist_cpu[i * 4u + bits];
+                    }
+                    unsigned int total_b = sum;
+                    if (bits + 1u < 4u)
+                        cumulative[bits + 1u] = cumulative[bits] + total_b;
+                }
+                for (size_t i = 0; i < num_blocks; ++i)
+                    for (unsigned int bits = 0u; bits < 4u; ++bits)
+                        prefix_cpu[i * 4u + bits] += cumulative[bits];
+                prefix_gpu.writeN(prefix_cpu.data(), num_blocks * 4u);
+
+                ocl_radixSort04Scatter.exec(ws, *values_src, prefix_gpu, *values_dst, static_cast<unsigned int>(n), bit_offset);
+
+                std::swap(values_src, values_dst);
+            }
+
+            values_src->copyToN(buffer_output_gpu, n);
         } else {
             rassert(false, 4531412341, context.type());
         }