CE_ArraySortImpl.h

Go to the documentation of this file.

/*
 * PROPRIETARY INFORMATION.  This software is proprietary to
 * Side Effects Software Inc., and is not to be reproduced,
 * transmitted, or disclosed in any way without written permission.
 *
 * NAME: CE_Array.h ( CE Library, C++)
 *
 * COMMENTS: UT_Array style class on GPU.
 */

#ifndef __CE_ArraySortImpl__
#define __CE_ArraySortImpl__

#include "CE_Array.h"

#include <UT/UT_Tracing.h>
#include <UT/UT_WorkBuffer.h>
#include <UT/UT_Swap.h>

// Appends a type name we use as part of type defines in array.cl
template <typename T> template <typename V>
void
CE_Array<T>::appendElemType(UT_WorkBuffer &wb)
{
    if constexpr (std::is_same<V, uint8>::value)
        wb.append("UCHAR");
    if constexpr (std::is_same<V, int8>::value)
        wb.append("CHAR");
    if constexpr (std::is_same<V, uint16>::value)
        wb.append("USHORT");
    if constexpr (std::is_same<V, int16>::value)
        wb.append("SHORT");
    if constexpr (std::is_same<V, int32>::value)
        wb.append("INT");
    if constexpr (std::is_same<V, uint32>::value)
        wb.append("UINT");
    if constexpr (std::is_same<V, int64>::value)
        wb.append("EXINT");
    if constexpr (std::is_same<V, uint64>::value)
        wb.append("ULONG");
    if constexpr (std::is_same<V, fpreal16>::value)
        wb.append("HALF");
    if constexpr (std::is_same<V, fpreal32>::value)
        wb.append("FLOAT");
    if constexpr (std::is_same<V, fpreal64>::value)
        wb.append("DOUBLE");
}

template<typename T> template 
<typename V>
void
CE_Array<T>::sortInternal(CE_Array<V> &vals, bool is_descending, int maxbits)
{
    utZoneScoped;
    exint nelem = CE_BufferDevice<T>::size();
    utZoneValue(nelem);
    bool has_vals = !vals.isEmpty();
    CE_Context *context = CE_Context::getContext();

    // Number of compute units and local memory size.
    cl::Device dev = context->getDevice();
    size_t ncompute = dev.getInfo<CL_DEVICE_MAX_COMPUTE_UNITS>();
    size_t nlocal = dev.getInfo<CL_DEVICE_LOCAL_MEM_SIZE>();

    // Bits per radix.
    constexpr size_t radixbits = 8;
    // Radix value.
    constexpr size_t radix = 1 << radixbits;
    // How many total bits in the type.
    size_t totalbits = 8 * sizeof(T);
    // Limit the number of bits considered if provided.
    // We still have to do it in multiples of radixbits.
    if (maxbits > 0)
    {
        totalbits = SYSmin(SYSroundUpToMultipleOf(size_t(maxbits), radixbits),
                           totalbits);
    }
    size_t npasses =  totalbits / radixbits;
    // Max number of elements each workitem can process.
    constexpr size_t maxelemperworkitem = 256;
    // Workitems per group.
    size_t nworkitems = 16;
    // Usually 1K per workitem at radixbits = 8.
    size_t localsize = sizeof(uint32) * radix * nworkitems;
    // Give a little extra headroom in local memory.
    while (localsize > nlocal * 1.2)
    {
        nworkitems /= 2;
        localsize = sizeof(uint32) * radix * nworkitems;
    }
    // Elements processed per group.
    size_t maxelempergroup = maxelemperworkitem * nworkitems;
    // Just a heuristic, this is the maximum possible number of groups.
    // We seem to get decreasing performance above 128.
    size_t ngroups = SYSmin(size_t(128), ncompute * 4);
    // Clamp groups from 1 to the min number necessary to process all elements.
    size_t mingroups = SYSroundUpToMultipleOf(size_t(nelem), maxelempergroup) / maxelempergroup;
    // size_t SYSclamp is ambigious on mac so use int version.
    ngroups = SYSclamp((int)ngroups, 1, (int)mingroups);

    // TODO - write a single workgroup/single pass sort for when ngroups == 1.
    cl::NDRange globalRange {ngroups * nworkitems};
    cl::NDRange localRange {nworkitems};
    const char *opt = is_descending ? "-D SORT_DESCENDING" : nullptr;
    UT_WorkBuffer wb;
    if (has_vals)
    {        
        // Append a value flag and type.
        wb.append(opt);
        wb.append(" -D HAS_values -D CEARRAY_VALUE_");
        appendElemType<V>(wb);
        opt = wb.buffer();
    }
    cl::Kernel khist = loadKernel("radix_sort_histogram", opt);
    cl::KernelFunctor hist = khist.bind(context->getQueue(), globalRange, localRange);
    cl::Kernel kreorder = loadKernel("radix_sort_reorder", opt);
    cl::KernelFunctor reorder = kreorder.bind(context->getQueue(), globalRange, localRange);

    // Each group takes nworkitems Kb at radixbits = 8.
    CE_UInt32Array histogram(radix * nworkitems * ngroups);
    CE_UInt32Array histosums(histogram.size());
    auto localarg = cl::__local(localsize);

    // Set up secondary key and value (if needed) buffers.
    CE_Array<T> dst_array(nelem);
    CE_Array<V> dst_vals(has_vals ? nelem: 0);

    for(size_t pass=0; pass < npasses; pass++)
    {
        const cl::Buffer &src = CE_BufferDevice<T>::buffer();
        // Histogram.
        hist(src, nelem, static_cast<int>(pass), histogram.buffer(), localarg);
        // Histogram prefix sum.
        histogram.prefixSum(histosums);
        // Sort, possibly with values.
        if (has_vals)
        {
            reorder(src, nelem, static_cast<int>(pass), histosums.buffer(),
                    dst_array.buffer(), vals.buffer(), dst_vals.buffer(), localarg);
            UTswap(vals, dst_vals);
        }
        else
        {
            reorder(src, nelem, static_cast<int>(pass),
                    histosums.buffer(), dst_array.buffer(), localarg);
        }
        // Swap buffers.
        UTswap(*this, dst_array);
    }
}

#endif