我正在解决一个问题,它涉及到快速排序10个数字(int32)。我的应用程序需要以尽可能快的速度对10个数字进行百万次排序。我正在对一个包含数十亿个元素的数据集进行采样,每次我都需要从中挑选10个数字(简化)并对它们进行排序(并从排序的10个元素列表中得出结论)。

目前我正在使用插入排序,但我想我可以为我的10个数字的特定问题实现一个非常快速的自定义排序算法,这将击败插入排序。

我该如何处理这个问题呢?


当前回答

那么展开的、无分支的选择排序呢?

#include <iostream>
#include <algorithm>
#include <random>

//return the index of the minimum element in array a
int min(const int * const a) {
  int m = a[0];
  int indx = 0;
  #define TEST(i) (m > a[i]) && (m = a[i], indx = i ); 
  //see http://stackoverflow.com/a/7074042/2140449
  TEST(1);
  TEST(2);
  TEST(3);
  TEST(4);
  TEST(5);
  TEST(6);
  TEST(7);
  TEST(8);
  TEST(9);
  #undef TEST
  return indx;
}

void sort( int * const a ){
  int work[10];
  int indx;
  #define GET(i) indx = min(a); work[i] = a[indx]; a[indx] = 2147483647; 
  //get the minimum, copy it to work and set it at max_int in a
  GET(0);
  GET(1);
  GET(2);
  GET(3);
  GET(4);
  GET(5);
  GET(6);
  GET(7);
  GET(8);
  GET(9);
  #undef GET
  #define COPY(i) a[i] = work[i];
  //copy back to a
  COPY(0);
  COPY(1);
  COPY(2);
  COPY(3);
  COPY(4);
  COPY(5);
  COPY(6);
  COPY(7);
  COPY(8);
  COPY(9);
  #undef COPY
}

int main() {
  //generating and printing a random array
  int a[10] = { 1,2,3,4,5,6,7,8,9,10 };
  std::random_device rd;
  std::mt19937 g(rd());
  std::shuffle( a, a+10, g);
  for (int i = 0; i < 10; i++) {
    std::cout << a[i] << ' ';
  }
  std::cout << std::endl;

  //sorting and printing again
  sort(a);
  for (int i = 0; i < 10; i++) {
    std::cout << a[i] << ' ';
  } 

  return 0;
}

http://coliru.stacked-crooked.com/a/71e18bc4f7fa18c6

唯一相关的行是前两行#define。

它使用两个列表,并将第一个列表重新检查十次,这将是一个糟糕的选择排序,但它避免了分支和变长循环,这可能会补偿现代处理器和如此小的数据集。


基准

我对排序网络进行了基准测试,我的代码似乎变慢了。然而,我试图删除展开和复制。运行以下代码:

#include <iostream>
#include <algorithm>
#include <random>
#include <chrono>

int min(const int * const a, int i) {
  int m = a[i];
  int indx = i++;
  for ( ; i<10; i++) 
    //see http://stackoverflow.com/a/7074042/2140449
    (m > a[i]) && (m = a[i], indx = i ); 
  return indx;
}

void sort( int * const a ){
  for (int i = 0; i<9; i++)
    std::swap(a[i], a[min(a,i)]); //search only forward
}


void sortNet10(int * const data) {  // ten-input sorting network by Waksman, 1969
    int swap;
    if (data[0] > data[5]) { swap = data[0]; data[0] = data[5]; data[5] = swap; }
    if (data[1] > data[6]) { swap = data[1]; data[1] = data[6]; data[6] = swap; }
    if (data[2] > data[7]) { swap = data[2]; data[2] = data[7]; data[7] = swap; }
    if (data[3] > data[8]) { swap = data[3]; data[3] = data[8]; data[8] = swap; }
    if (data[4] > data[9]) { swap = data[4]; data[4] = data[9]; data[9] = swap; }
    if (data[0] > data[3]) { swap = data[0]; data[0] = data[3]; data[3] = swap; }
    if (data[5] > data[8]) { swap = data[5]; data[5] = data[8]; data[8] = swap; }
    if (data[1] > data[4]) { swap = data[1]; data[1] = data[4]; data[4] = swap; }
    if (data[6] > data[9]) { swap = data[6]; data[6] = data[9]; data[9] = swap; }
    if (data[0] > data[2]) { swap = data[0]; data[0] = data[2]; data[2] = swap; }
    if (data[3] > data[6]) { swap = data[3]; data[3] = data[6]; data[6] = swap; }
    if (data[7] > data[9]) { swap = data[7]; data[7] = data[9]; data[9] = swap; }
    if (data[0] > data[1]) { swap = data[0]; data[0] = data[1]; data[1] = swap; }
    if (data[2] > data[4]) { swap = data[2]; data[2] = data[4]; data[4] = swap; }
    if (data[5] > data[7]) { swap = data[5]; data[5] = data[7]; data[7] = swap; }
    if (data[8] > data[9]) { swap = data[8]; data[8] = data[9]; data[9] = swap; }
    if (data[1] > data[2]) { swap = data[1]; data[1] = data[2]; data[2] = swap; }
    if (data[3] > data[5]) { swap = data[3]; data[3] = data[5]; data[5] = swap; }
    if (data[4] > data[6]) { swap = data[4]; data[4] = data[6]; data[6] = swap; }
    if (data[7] > data[8]) { swap = data[7]; data[7] = data[8]; data[8] = swap; }
    if (data[1] > data[3]) { swap = data[1]; data[1] = data[3]; data[3] = swap; }
    if (data[4] > data[7]) { swap = data[4]; data[4] = data[7]; data[7] = swap; }
    if (data[2] > data[5]) { swap = data[2]; data[2] = data[5]; data[5] = swap; }
    if (data[6] > data[8]) { swap = data[6]; data[6] = data[8]; data[8] = swap; }
    if (data[2] > data[3]) { swap = data[2]; data[2] = data[3]; data[3] = swap; }
    if (data[4] > data[5]) { swap = data[4]; data[4] = data[5]; data[5] = swap; }
    if (data[6] > data[7]) { swap = data[6]; data[6] = data[7]; data[7] = swap; }
    if (data[3] > data[4]) { swap = data[3]; data[3] = data[4]; data[4] = swap; }
    if (data[5] > data[6]) { swap = data[5]; data[5] = data[6]; data[6] = swap; }
}


std::chrono::duration<double> benchmark( void(*func)(int * const), const int seed ) {
  std::mt19937 g(seed);
  int a[10] = {10,11,12,13,14,15,16,17,18,19};
  std::chrono::high_resolution_clock::time_point t1, t2; 
  t1 = std::chrono::high_resolution_clock::now();
  for (long i = 0; i < 1e7; i++) {
    std::shuffle( a, a+10, g);
    func(a);
  }
  t2 = std::chrono::high_resolution_clock::now();
  return std::chrono::duration_cast<std::chrono::duration<double>>(t2 - t1);
}

int main() {
  std::random_device rd;
  for (int i = 0; i < 10; i++) {
    const int seed = rd();
    std::cout << "seed = " << seed << std::endl;
    std::cout << "sortNet10: " << benchmark(sortNet10, seed).count() << std::endl;
    std::cout << "sort:      " << benchmark(sort,      seed).count() << std::endl;
  }
  return 0;
}

与排序网络相比,我一直在为无分支选择排序获得更好的结果。

$ gcc -v
gcc version 5.2.0 (GCC) 
$ g++ -std=c++11 -Ofast sort.cpp && ./a.out
seed = -1727396418
sortNet10: 2.24137
sort:      2.21828
seed = 2003959850
sortNet10: 2.23914
sort:      2.21641
seed = 1994540383
sortNet10: 2.23782
sort:      2.21778
seed = 1258259982
sortNet10: 2.25199
sort:      2.21801
seed = 1821086932
sortNet10: 2.25535
sort:      2.2173
seed = 412262735
sortNet10: 2.24489
sort:      2.21776
seed = 1059795817
sortNet10: 2.29226
sort:      2.21777
seed = -188551272
sortNet10: 2.23803
sort:      2.22996
seed = 1043757247
sortNet10: 2.2503
sort:      2.23604
seed = -268332483
sortNet10: 2.24455
sort:      2.24304

其他回答

当您处理这个固定大小时,请查看排序网络。这些算法有固定的运行时间,并且独立于它们的输入。对于您的用例,您没有某些排序算法所具有的这种开销。

二进制排序就是这种网络的一种实现。这个方法在CPU上使用len(n) <= 32时效果最好。对于更大的输入,你可以考虑使用GPU。

顺便说一下,比较排序算法的一个好页面是这个(尽管它缺少二进制排序):

排序算法动画

(根据@HelloWorld的建议,研究排序网络。)

似乎29个比较/交换网络是进行10个输入排序的最快方法。在这个例子中,我使用了Waksman在1969年发现的JavaScript网络,它应该直接转换成C语言,因为它只是一个if语句、比较和交换的列表。

function sortNet10(data) { // ten-input sorting network by Waksman, 1969 var swap; if (data[0] > data[5]) { swap = data[0]; data[0] = data[5]; data[5] = swap; } if (data[1] > data[6]) { swap = data[1]; data[1] = data[6]; data[6] = swap; } if (data[2] > data[7]) { swap = data[2]; data[2] = data[7]; data[7] = swap; } if (data[3] > data[8]) { swap = data[3]; data[3] = data[8]; data[8] = swap; } if (data[4] > data[9]) { swap = data[4]; data[4] = data[9]; data[9] = swap; } if (data[0] > data[3]) { swap = data[0]; data[0] = data[3]; data[3] = swap; } if (data[5] > data[8]) { swap = data[5]; data[5] = data[8]; data[8] = swap; } if (data[1] > data[4]) { swap = data[1]; data[1] = data[4]; data[4] = swap; } if (data[6] > data[9]) { swap = data[6]; data[6] = data[9]; data[9] = swap; } if (data[0] > data[2]) { swap = data[0]; data[0] = data[2]; data[2] = swap; } if (data[3] > data[6]) { swap = data[3]; data[3] = data[6]; data[6] = swap; } if (data[7] > data[9]) { swap = data[7]; data[7] = data[9]; data[9] = swap; } if (data[0] > data[1]) { swap = data[0]; data[0] = data[1]; data[1] = swap; } if (data[2] > data[4]) { swap = data[2]; data[2] = data[4]; data[4] = swap; } if (data[5] > data[7]) { swap = data[5]; data[5] = data[7]; data[7] = swap; } if (data[8] > data[9]) { swap = data[8]; data[8] = data[9]; data[9] = swap; } if (data[1] > data[2]) { swap = data[1]; data[1] = data[2]; data[2] = swap; } if (data[3] > data[5]) { swap = data[3]; data[3] = data[5]; data[5] = swap; } if (data[4] > data[6]) { swap = data[4]; data[4] = data[6]; data[6] = swap; } if (data[7] > data[8]) { swap = data[7]; data[7] = data[8]; data[8] = swap; } if (data[1] > data[3]) { swap = data[1]; data[1] = data[3]; data[3] = swap; } if (data[4] > data[7]) { swap = data[4]; data[4] = data[7]; data[7] = swap; } if (data[2] > data[5]) { swap = data[2]; data[2] = data[5]; data[5] = swap; } if (data[6] > data[8]) { swap = data[6]; data[6] = data[8]; data[8] = swap; } if (data[2] > data[3]) { swap = data[2]; data[2] = data[3]; data[3] = swap; } if (data[4] > data[5]) { swap = data[4]; data[4] = data[5]; data[5] = swap; } if (data[6] > data[7]) { swap = data[6]; data[6] = data[7]; data[7] = swap; } if (data[3] > data[4]) { swap = data[3]; data[3] = data[4]; data[4] = swap; } if (data[5] > data[6]) { swap = data[5]; data[5] = data[6]; data[6] = swap; } return(data); } document.write(sortNet10([5,7,1,8,4,3,6,9,2,0]));

这里是网络的图形表示,分为独立的阶段。

为了利用并行处理的优势,可以将5-4-3-3 - 4-4-2 -3-2分组改为4-4-4-2 -4-4-3-2分组。

以下是运行在10个CUDA线程上的CUDA内核(秩排序算法),在42毫秒内对1000个数组排序1000次,每次排序42纳秒或每次排序~70个周期(1.7 GHz):

inline
__device__ int findOrder(const int mask, const int data0)
{
    const int order1 = data0>__shfl_sync(mask,data0,0);
    const int order2 = data0>__shfl_sync(mask,data0,1);
    const int order3 = data0>__shfl_sync(mask,data0,2);
    const int order4 = data0>__shfl_sync(mask,data0,3);
    const int order5 = data0>__shfl_sync(mask,data0,4);
    const int order6 = data0>__shfl_sync(mask,data0,5);
    const int order7 = data0>__shfl_sync(mask,data0,6);
    const int order8 = data0>__shfl_sync(mask,data0,7);
    const int order9 = data0>__shfl_sync(mask,data0,8);
    const int order10 = data0>__shfl_sync(mask,data0,9);
    return order1 + order2 + order3 + order4 + order5 + order6 + order7 + order8 + order9 + order10;
}

// launch this with 10 CUDA threads in 1 block in 1 grid
// sorts 10 elements using only SIMT registers
__global__ void rankSort(int * __restrict__ buffer)
{    
    const int id  = threadIdx.x;

    // enable first 10 lanes of warp for shuffling 
    const int mask = __activemask();

    __shared__ int data[10000];

    for(int i=0;i<1000;i++)
    {
        data[id+i*10] = buffer[id+i*10];
    }
    __syncwarp();
    // benchmark 1000 times
    for(int k=0;k<1000;k++)
    {

        // sort 1000 arrays
        for(int j=0;j<1000;j+=5)
        {
            // sort 5 arrays at once to hide latency
            const int data1 = data[id+j*10];
            const int data2 = data[id+(j+1)*10];
            const int data3 = data[id+(j+2)*10];
            const int data4 = data[id+(j+3)*10];
            const int data5 = data[id+(j+4)*10];

            const int order1 = findOrder(mask,data1);        
            const int order2 = findOrder(mask,data2);
            const int order3 = findOrder(mask,data3);
            const int order4 = findOrder(mask,data4);
            const int order5 = findOrder(mask,data5);

            data[order1+j*10]=data1;         
            data[order2+(j+1)*10]=data2;           
            data[order3+(j+2)*10]=data3;
            data[order4+(j+3)*10]=data4;
            data[order5+(j+4)*10]=data5;
        }

    }
    __syncwarp();
    for(int i=0;i<1000;i++)
    {
        buffer[id+i*10] = data[id+i*10];
    }
}   

由于所有10个线程都在同一个SIMT单元上处理,它类似于运行在矢量寄存器上的AVX512优化版本,但除了更多的寄存器空间以隐藏更多的延迟之外。此外,SIMT单元是32宽的,因此它可以运行线性时间复杂度直到32个元素。

该算法假设元素是唯一的。对于重复的数据,需要一个额外的缩减步骤来将重复的顺序值解压缩为10个元素。首先,它计算每个元素的秩,然后直接将它们复制到数组中作为索引的秩。顺序值需要蛮力(O(N x N))次比较,为了减少延迟,使用warp-shuffles从(向量寄存器的)不同的warp-lanes收集数据。

插入排序平均需要29,6次与排序10个输入的比较,最佳情况为9,最差情况为45(给定输入的顺序是相反的)。

{9,6,1} shellsort平均需要25.5次比较来排序10个输入。最好的情况是14次比较,最坏的情况是34次,对反向输入排序需要22次。

因此,使用shellsort而不是插入排序可以减少14%的平均情况。尽管最佳情况增加了56%,但最坏情况减少了24%,这对于控制最坏情况性能非常重要的应用程序来说意义重大。反之则减少51%。

因为你似乎对插入排序很熟悉,你可以将算法实现为{9,6}的排序网络,然后在那之后附加插入排序({1}):

i[0] with i[9]    // {9}

i[0] with i[6]    // {6}
i[1] with i[7]    // {6}
i[2] with i[8]    // {6}
i[3] with i[9]    // {6}

i[0 ... 9]        // insertion sort

这个问题并没有说这是某种基于web的应用程序。有一件事引起了我的注意:

我正在对一个包含数十亿个元素的数据集进行采样,每次我都需要从中挑选10个数字(简化)并对它们进行排序(并从排序的10个元素列表中得出结论)。

As a software and hardware engineer this absolutely screams FPGA to me. I don't know what kind of conclusions you need to draw from the sorted set of numbers or where the data comes from, but I know it would be almost trivial to process somewhere between one hundred million and a billion of these "sort-and-analyze" operations per second. I've done FPGA-assisted DNA sequencing work in the past. It is nearly impossible to beat the massive processing power of FPGAs when the problem is well suited for that type of a solution.

在某种程度上,唯一的限制因素是将数据铲入FPGA的速度有多快,以及取出数据的速度有多快。

As a point of reference, I designed a high performance real-time image processor that received 32 bit RGB image data at a rate of about 300 million pixels per second. The data streamed through FIR filters, matrix multipliers, lookup tables, spatial edge detection blocks and a number of other operations before coming out the other end. All of this on a relatively small Xilinx Virtex2 FPGA with internal clocking spanning from about 33 MHz to, if I remember correctly, 400 MHz. Oh, yes, it also had a DDR2 controller implementation and ran two banks of DDR2 memory.

当工作在数百MHz时,FPGA可以在每次时钟转换中输出10个32位数字。当数据填满处理管道时,操作开始时会有短暂的延迟。在此之后,您应该能够在每个时钟获得一个结果。如果可以通过复制排序和分析管道使处理并行化,则会更多。原则上,解决方案几乎是微不足道的。

关键在于:如果应用程序不是pc绑定的,并且数据流和处理与FPGA解决方案“兼容”(无论是独立的还是作为机器中的协处理器卡),那么无论使用哪种算法,用任何语言编写的软件都无法击败可达到的性能水平。

我只是快速搜索了一下,找到了一篇可能对你有用的论文。看起来要追溯到2012年。在今天(甚至在过去),您可以在性能上做得更好。下面就是:

fpga上的排序网络