在回答另一个Stack Overflow问题时,我偶然发现了一个有趣的子问题。对6个整数的数组进行排序的最快方法是什么?

因为问题层次很低:

我们不能假设库是可用的(而且调用本身也有开销),只有纯C 为了避免清空指令管道(这有非常高的成本),我们可能应该最小化分支、跳转和其他类型的控制流中断(比如隐藏在&&或||序列点后面的那些)。 空间是有限的,最小化寄存器和内存的使用是一个问题,理想情况下,就地排序可能是最好的。

实际上,这个问题是一种Golf,其目标不是最小化源长度,而是最小化执行时间。我称之为“Zening”代码,就像Michael Abrash在《Zen of code optimization》一书及其续集中所使用的那样。

至于为什么它有趣,有几个层面:

示例简单,易于理解和测量,不需要太多的C技能 它显示了对问题选择好的算法的影响,也显示了编译器和底层硬件的影响。

下面是我的参考(简单的,不是优化的)实现和测试集。

#include <stdio.h>

static __inline__ int sort6(int * d){

    char j, i, imin;
    int tmp;
    for (j = 0 ; j < 5 ; j++){
        imin = j;
        for (i = j + 1; i < 6 ; i++){
            if (d[i] < d[imin]){
                imin = i;
            }
        }
        tmp = d[j];
        d[j] = d[imin];
        d[imin] = tmp;
    }
}

static __inline__ unsigned long long rdtsc(void)
{
  unsigned long long int x;
     __asm__ volatile (".byte 0x0f, 0x31" : "=A" (x));
     return x;
}

int main(int argc, char ** argv){
    int i;
    int d[6][5] = {
        {1, 2, 3, 4, 5, 6},
        {6, 5, 4, 3, 2, 1},
        {100, 2, 300, 4, 500, 6},
        {100, 2, 3, 4, 500, 6},
        {1, 200, 3, 4, 5, 600},
        {1, 1, 2, 1, 2, 1}
    };

    unsigned long long cycles = rdtsc();
    for (i = 0; i < 6 ; i++){
        sort6(d[i]);
        /*
         * printf("d%d : %d %d %d %d %d %d\n", i,
         *  d[i][0], d[i][6], d[i][7],
         *  d[i][8], d[i][9], d[i][10]);
        */
    }
    cycles = rdtsc() - cycles;
    printf("Time is %d\n", (unsigned)cycles);
}

生的结果

随着变体的数量越来越多,我将它们都收集到一个测试套件中,可以在这里找到。在Kevin Stock的帮助下,实际使用的测试没有上面展示的那么简单。您可以在自己的环境中编译和执行它。我对不同目标架构/编译器上的行为很感兴趣。(好了,伙计们,把它放在答案里,我将+1一个新结果集的每个贡献者)。

一年前,我把答案给了Daniel Stutzbach(高尔夫),因为他是当时最快的解决方案(排序网络)的来源。

Linux 64位,gcc 4.6.1 64位,Intel Core 2 Duo E8400, -O2

Direct call to qsort library function : 689.38 Naive implementation (insertion sort) : 285.70 Insertion Sort (Daniel Stutzbach) : 142.12 Insertion Sort Unrolled : 125.47 Rank Order : 102.26 Rank Order with registers : 58.03 Sorting Networks (Daniel Stutzbach) : 111.68 Sorting Networks (Paul R) : 66.36 Sorting Networks 12 with Fast Swap : 58.86 Sorting Networks 12 reordered Swap : 53.74 Sorting Networks 12 reordered Simple Swap : 31.54 Reordered Sorting Network w/ fast swap : 31.54 Reordered Sorting Network w/ fast swap V2 : 33.63 Inlined Bubble Sort (Paolo Bonzini) : 48.85 Unrolled Insertion Sort (Paolo Bonzini) : 75.30

Linux 64位,gcc 4.6.1 64位,Intel Core 2 Duo E8400, -O1

Direct call to qsort library function : 705.93 Naive implementation (insertion sort) : 135.60 Insertion Sort (Daniel Stutzbach) : 142.11 Insertion Sort Unrolled : 126.75 Rank Order : 46.42 Rank Order with registers : 43.58 Sorting Networks (Daniel Stutzbach) : 115.57 Sorting Networks (Paul R) : 64.44 Sorting Networks 12 with Fast Swap : 61.98 Sorting Networks 12 reordered Swap : 54.67 Sorting Networks 12 reordered Simple Swap : 31.54 Reordered Sorting Network w/ fast swap : 31.24 Reordered Sorting Network w/ fast swap V2 : 33.07 Inlined Bubble Sort (Paolo Bonzini) : 45.79 Unrolled Insertion Sort (Paolo Bonzini) : 80.15

我包括了-O1和-O2的结果,因为令人惊讶的是,在一些程序中,O2的效率低于O1。我想知道什么具体的优化有这种效果?

对建议解决方案的评论

插入排序(丹尼尔·斯图茨巴赫)

正如预期的那样,最小化分支确实是一个好主意。

排序网络(丹尼尔·斯图茨巴赫)

比插入排序好。我想知道主要的效果是不是避免外部循环。我试着通过展开插入排序来检查,确实我们得到了大致相同的数字(代码在这里)。

排序网络(保罗R)

迄今为止最好的。我用来测试的实际代码在这里。目前还不知道为什么它的速度几乎是其他排序网络实现的两倍。参数传递?快速max ?

排序网络12 SWAP与快速交换

根据Daniel Stutzbach的建议,我将他的12交换排序网络与无分支快速交换相结合(代码在这里)。它确实更快,到目前为止最好的,只有很小的利润率(大约5%),因为可以使用更少的交换。

同样有趣的是,无分支交换似乎比在PPC架构上使用if的简单交换效率低得多(4倍)。

调用库qsort

To give another reference point I also tried as suggested to just call library qsort (code is here). As expected it is much slower : 10 to 30 times slower... as it became obvious with the new test suite, the main problem seems to be the initial load of the library after the first call, and it compares not so poorly with other version. It is just between 3 and 20 times slower on my Linux. On some architecture used for tests by others it seems even to be faster (I'm really surprised by that one, as library qsort use a more complex API).

等级次序

Rex Kerr proposed another completely different method : for each item of the array compute directly its final position. This is efficient because computing rank order do not need branch. The drawback of this method is that it takes three times the amount of memory of the array (one copy of array and variables to store rank orders). The performance results are very surprising (and interesting). On my reference architecture with 32 bits OS and Intel Core2 Quad E8300, cycle count was slightly below 1000 (like sorting networks with branching swap). But when compiled and executed on my 64 bits box (Intel Core2 Duo) it performed much better : it became the fastest so far. I finally found out the true reason. My 32bits box use gcc 4.4.1 and my 64bits box gcc 4.4.3 and the last one seems much better at optimizing this particular code (there was very little difference for other proposals).

更新:

正如上面公布的数字所示,这种效果在gcc的后续版本中仍然得到了增强,Rank Order的速度始终是其他任何替代版本的两倍。

用重新排序的交换对网络进行排序

The amazing efficiency of the Rex Kerr proposal with gcc 4.4.3 made me wonder : how could a program with 3 times as much memory usage be faster than branchless sorting networks? My hypothesis was that it had less dependencies of the kind read after write, allowing for better use of the superscalar instruction scheduler of the x86. That gave me an idea: reorder swaps to minimize read after write dependencies. More simply put: when you do SWAP(1, 2); SWAP(0, 2); you have to wait for the first swap to be finished before performing the second one because both access to a common memory cell. When you do SWAP(1, 2); SWAP(4, 5);the processor can execute both in parallel. I tried it and it works as expected, the sorting networks is running about 10% faster.

使用简单交换对网络进行排序

One year after the original post Steinar H. Gunderson suggested, that we should not try to outsmart the compiler and keep the swap code simple. It's indeed a good idea as the resulting code is about 40% faster! He also proposed a swap optimized by hand using x86 inline assembly code that can still spare some more cycles. The most surprising (it says volumes on programmer's psychology) is that one year ago none of used tried that version of swap. Code I used to test is here. Others suggested other ways to write a C fast swap, but it yields the same performances as the simple one with a decent compiler.

“最佳”代码如下:

static inline void sort6_sorting_network_simple_swap(int * d){
#define min(x, y) (x<y?x:y)
#define max(x, y) (x<y?y:x) 
#define SWAP(x,y) { const int a = min(d[x], d[y]); \
                    const int b = max(d[x], d[y]); \
                    d[x] = a; d[y] = b; }
    SWAP(1, 2);
    SWAP(4, 5);
    SWAP(0, 2);
    SWAP(3, 5);
    SWAP(0, 1);
    SWAP(3, 4);
    SWAP(1, 4);
    SWAP(0, 3);
    SWAP(2, 5);
    SWAP(1, 3);
    SWAP(2, 4);
    SWAP(2, 3);
#undef SWAP
#undef min
#undef max
}

如果我们相信我们的测试集(是的,它很差,它的唯一好处是简短,简单,易于理解我们所测量的内容),那么一个排序的结果代码的平均循环次数低于40个循环(执行6个测试)。这使得每次交换平均为4个周期。我称之为惊人的快。还有其他可能的改进吗?


当前回答

也许我来晚了,但至少我的贡献是一种新的方法。

The code really should be inlined even if inlined, there are too many branches the analysing part is basically O(N(N-1)) which seems OK for N=6 the code could be more effective if the cost of swap would be higher (irt the cost of compare) I trust on static functions being inlined. The method is related to rank-sort instead of ranks, the relative ranks (offsets) are used. the sum of the ranks is zero for every cycle in any permutation group. instead of SWAP()ing two elements, the cycles are chased, needing only one temp, and one (register->register) swap (new <- old).


更新:代码改动了一点,有些人用c++编译器编译C代码…

#include <stdio.h>

#if WANT_CHAR
typedef signed char Dif;
#else
typedef signed int Dif;
#endif

static int walksort (int *arr, int cnt);
static void countdifs (int *arr, Dif *dif, int cnt);
static void calcranks(int *arr, Dif *dif);

int wsort6(int *arr);

void do_print_a(char *msg, int *arr, unsigned cnt)
{
fprintf(stderr,"%s:", msg);
for (; cnt--; arr++) {
        fprintf(stderr, " %3d", *arr);
        }
fprintf(stderr,"\n");
}

void do_print_d(char *msg, Dif *arr, unsigned cnt)
{
fprintf(stderr,"%s:", msg);
for (; cnt--; arr++) {
        fprintf(stderr, " %3d", (int) *arr);
        }
fprintf(stderr,"\n");
}

static void inline countdifs (int *arr, Dif *dif, int cnt)
{
int top, bot;

for (top = 0; top < cnt; top++ ) {
        for (bot = 0; bot < top; bot++ ) {
                if (arr[top] < arr[bot]) { dif[top]--; dif[bot]++; }
                }
        }
return ;
}
        /* Copied from RexKerr ... */
static void inline calcranks(int *arr, Dif *dif){

dif[0] =     (arr[0]>arr[1])+(arr[0]>arr[2])+(arr[0]>arr[3])+(arr[0]>arr[4])+(arr[0]>arr[5]);
dif[1] = -1+ (arr[1]>=arr[0])+(arr[1]>arr[2])+(arr[1]>arr[3])+(arr[1]>arr[4])+(arr[1]>arr[5]);
dif[2] = -2+ (arr[2]>=arr[0])+(arr[2]>=arr[1])+(arr[2]>arr[3])+(arr[2]>arr[4])+(arr[2]>arr[5]);
dif[3] = -3+ (arr[3]>=arr[0])+(arr[3]>=arr[1])+(arr[3]>=arr[2])+(arr[3]>arr[4])+(arr[3]>arr[5]);
dif[4] = -4+ (arr[4]>=arr[0])+(arr[4]>=arr[1])+(arr[4]>=arr[2])+(arr[4]>=arr[3])+(arr[4]>arr[5]);
dif[5] = -(dif[0]+dif[1]+dif[2]+dif[3]+dif[4]);
}

static int walksort (int *arr, int cnt)
{
int idx, src,dst, nswap;

Dif difs[cnt];

#if WANT_REXK
calcranks(arr, difs);
#else
for (idx=0; idx < cnt; idx++) difs[idx] =0;
countdifs(arr, difs, cnt);
#endif
calcranks(arr, difs);

#define DUMP_IT 0
#if DUMP_IT
do_print_d("ISteps ", difs, cnt);
#endif

nswap = 0;
for (idx=0; idx < cnt; idx++) {
        int newval;
        int step,cyc;
        if ( !difs[idx] ) continue;
        newval = arr[idx];
        cyc = 0;
        src = idx;
        do      {
                int oldval;
                step = difs[src];
                difs[src] =0;
                dst = src + step;
                cyc += step ;
                if(dst == idx+1)idx=dst;
                oldval = arr[dst];
#if (DUMP_IT&1)
                fprintf(stderr, "[Nswap=%d] Cyc=%d Step=%2d Idx=%d  Old=%2d New=%2d #### Src=%d Dst=%d[%2d]->%2d <-- %d\n##\n"
                        , nswap, cyc, step, idx, oldval, newval
                        , src, dst, difs[dst], arr[dst]
                        , newval  );
                do_print_a("Array ", arr, cnt);
                do_print_d("Steps ", difs, cnt);
#endif

                arr[dst] = newval;
                newval = oldval;
                nswap++;
                src = dst;
                } while( cyc);
        }

return nswap;
}
/*************/
int wsort6(int *arr)
{
return walksort(arr, 6);
}

其他回答

//Bruteforce compute unrolled count dumbsort(min to 0-index)
void bcudc_sort6(int* a)
{
    int t[6] = {0};
    int r1,r2;

    r1=0;
    r1 += (a[0] > a[1]);
    r1 += (a[0] > a[2]);
    r1 += (a[0] > a[3]);
    r1 += (a[0] > a[4]);
    r1 += (a[0] > a[5]);
    while(t[r1]){r1++;}
    t[r1] = a[0];

    r2=0;
    r2 += (a[1] > a[0]);
    r2 += (a[1] > a[2]);
    r2 += (a[1] > a[3]);
    r2 += (a[1] > a[4]);
    r2 += (a[1] > a[5]);
    while(t[r2]){r2++;} 
    t[r2] = a[1];

    r1=0;
    r1 += (a[2] > a[0]);
    r1 += (a[2] > a[1]);
    r1 += (a[2] > a[3]);
    r1 += (a[2] > a[4]);
    r1 += (a[2] > a[5]);
    while(t[r1]){r1++;}
    t[r1] = a[2];

    r2=0;
    r2 += (a[3] > a[0]);
    r2 += (a[3] > a[1]);
    r2 += (a[3] > a[2]);
    r2 += (a[3] > a[4]);
    r2 += (a[3] > a[5]);
    while(t[r2]){r2++;} 
    t[r2] = a[3];

    r1=0;
    r1 += (a[4] > a[0]);
    r1 += (a[4] > a[1]);
    r1 += (a[4] > a[2]);
    r1 += (a[4] > a[3]);
    r1 += (a[4] > a[5]);
    while(t[r1]){r1++;}
    t[r1] = a[4];

    r2=0;
    r2 += (a[5] > a[0]);
    r2 += (a[5] > a[1]);
    r2 += (a[5] > a[2]);
    r2 += (a[5] > a[3]);
    r2 += (a[5] > a[4]);
    while(t[r2]){r2++;} 
    t[r2] = a[5];

    a[0]=t[0];
    a[1]=t[1];
    a[2]=t[2];
    a[3]=t[3];
    a[4]=t[4];
    a[5]=t[5];
}

static __inline__ void sort6(int* a)
{
    #define wire(x,y); t = a[x] ^ a[y] ^ ( (a[x] ^ a[y]) & -(a[x] < a[y]) ); a[x] = a[x] ^ t; a[y] = a[y] ^ t;
    register int t;

    wire( 0, 1); wire( 2, 3); wire( 4, 5);
    wire( 3, 5); wire( 0, 2); wire( 1, 4);
    wire( 4, 5); wire( 2, 3); wire( 0, 1); 
    wire( 3, 4); wire( 1, 2); 
    wire( 2, 3);

    #undef wire
}

下面是一个使用排序网络的实现:

inline void Sort2(int *p0, int *p1)
{
    const int temp = min(*p0, *p1);
    *p1 = max(*p0, *p1);
    *p0 = temp;
}

inline void Sort3(int *p0, int *p1, int *p2)
{
    Sort2(p0, p1);
    Sort2(p1, p2);
    Sort2(p0, p1);
}

inline void Sort4(int *p0, int *p1, int *p2, int *p3)
{
    Sort2(p0, p1);
    Sort2(p2, p3);
    Sort2(p0, p2);  
    Sort2(p1, p3);  
    Sort2(p1, p2);  
}

inline void Sort6(int *p0, int *p1, int *p2, int *p3, int *p4, int *p5)
{
    Sort3(p0, p1, p2);
    Sort3(p3, p4, p5);
    Sort2(p0, p3);  
    Sort2(p2, p5);  
    Sort4(p1, p2, p3, p4);  
}

你真的需要非常高效的无分支的min和max实现,因为这实际上就是这段代码归结为的-一个min和max操作序列(总共13个)。我把这个留给读者做练习。

注意,这种实现很容易实现矢量化(例如SIMD -大多数SIMD isa都有矢量最小/最大指令),也很容易实现GPU(例如CUDA -无分支,不存在扭曲发散等问题)。

参见:对非常小的列表进行排序的快速算法实现

如果插入排序在这里是合理的竞争,我建议尝试shell排序。我担心6个元素可能太少了,不足以跻身最佳之列,但它可能值得一试。

示例代码,未测试,未调试等。您希望调优inc = 4和inc -= 3序列以找到最优序列(例如,尝试inc = 2, inc -= 1)。

static __inline__ int sort6(int * d) {
    char j, i;
    int tmp;
    for (inc = 4; inc > 0; inc -= 3) {
        for (i = inc; i < 5; i++) {
            tmp = a[i];
            j = i;
            while (j >= inc && a[j - inc] > tmp) {
                a[j] = a[j - inc];
                j -= inc;
            }
            a[j] = tmp;
        }
    }
}

我不认为这个会赢,但是如果有人发了一个关于排序10个元素的问题,谁知道呢……

根据维基百科,这甚至可以与排序网络相结合: 普拉特,V(1979)。贝壳排序和排序网络(计算机科学杰出论文)。花环。ISBN 0-824-04406-1

虽然我真的很喜欢交换宏提供:

#define min(x, y) (y ^ ((x ^ y) & -(x < y)))
#define max(x, y) (x ^ ((x ^ y) & -(x < y)))
#define SWAP(x,y) { int tmp = min(d[x], d[y]); d[y] = max(d[x], d[y]); d[x] = tmp; }

我看到了一个改进(一个好的编译器可能会做到):

#define SWAP(x,y) { int tmp = ((x ^ y) & -(y < x)); y ^= tmp; x ^= tmp; }

我们注意到min和max是如何工作的,并显式地提取公共子表达式。这完全消除了min和max宏。

几天前,我无意中从谷歌中发现了这个问题,因为我还需要快速排序一个由6个整数组成的固定长度数组。然而,在我的情况下,我的整数只有8位(而不是32位),我没有严格的要求只使用c。我想我无论如何都会分享我的发现,以防他们可能对某人有帮助……

我在程序集中实现了一个网络排序的变体,它使用SSE尽可能地向量化比较和交换操作。需要六次“传递”才能对数组进行完全排序。我使用了一种新颖的机制,直接将PCMPGTB(向量化比较)的结果转换为PSHUFB(向量化交换)的洗牌参数,只使用PADDB(向量化添加),在某些情况下还使用PAND(位与)指令。

这种方法也有产生真正无分支函数的副作用。没有任何跳跃指令。

这个实现似乎比目前在问题(“用简单交换排序网络12”)中被标记为最快选项的实现快38%左右。在测试期间,我修改了该实现以使用char数组元素,以使比较公平。

我应该指出,这种方法可以应用于任何大小不超过16个元素的数组。我希望在更大的数组中,相对于替代方案的速度优势会越来越大。

代码是用MASM编写的,适用于带有SSSE3的x86_64处理器。该函数使用“new”Windows x64调用约定。在这儿……

PUBLIC simd_sort_6

.DATA

ALIGN 16

pass1_shuffle   OWORD   0F0E0D0C0B0A09080706040503010200h
pass1_add       OWORD   0F0E0D0C0B0A09080706050503020200h
pass2_shuffle   OWORD   0F0E0D0C0B0A09080706030405000102h
pass2_and       OWORD   00000000000000000000FE00FEFE00FEh
pass2_add       OWORD   0F0E0D0C0B0A09080706050405020102h
pass3_shuffle   OWORD   0F0E0D0C0B0A09080706020304050001h
pass3_and       OWORD   00000000000000000000FDFFFFFDFFFFh
pass3_add       OWORD   0F0E0D0C0B0A09080706050404050101h
pass4_shuffle   OWORD   0F0E0D0C0B0A09080706050100020403h
pass4_and       OWORD   0000000000000000000000FDFD00FDFDh
pass4_add       OWORD   0F0E0D0C0B0A09080706050403020403h
pass5_shuffle   OWORD   0F0E0D0C0B0A09080706050201040300h
pass5_and       OWORD 0000000000000000000000FEFEFEFE00h
pass5_add       OWORD   0F0E0D0C0B0A09080706050403040300h
pass6_shuffle   OWORD   0F0E0D0C0B0A09080706050402030100h
pass6_add       OWORD   0F0E0D0C0B0A09080706050403030100h

.CODE

simd_sort_6 PROC FRAME

    .endprolog

    ; pxor xmm4, xmm4
    ; pinsrd xmm4, dword ptr [rcx], 0
    ; pinsrb xmm4, byte ptr [rcx + 4], 4
    ; pinsrb xmm4, byte ptr [rcx + 5], 5
    ; The benchmarked 38% faster mentioned in the text was with the above slower sequence that tied up the shuffle port longer.  Same on extract
    ; avoiding pins/extrb also means we don't need SSE 4.1, but SSSE3 CPUs without SSE4.1 (e.g. Conroe/Merom) have slow pshufb.
    movd    xmm4, dword ptr [rcx]
    pinsrw  xmm4,  word ptr [rcx + 4], 2  ; word 2 = bytes 4 and 5


    movdqa xmm5, xmm4
    pshufb xmm5, oword ptr [pass1_shuffle]
    pcmpgtb xmm5, xmm4
    paddb xmm5, oword ptr [pass1_add]
    pshufb xmm4, xmm5

    movdqa xmm5, xmm4
    pshufb xmm5, oword ptr [pass2_shuffle]
    pcmpgtb xmm5, xmm4
    pand xmm5, oword ptr [pass2_and]
    paddb xmm5, oword ptr [pass2_add]
    pshufb xmm4, xmm5

    movdqa xmm5, xmm4
    pshufb xmm5, oword ptr [pass3_shuffle]
    pcmpgtb xmm5, xmm4
    pand xmm5, oword ptr [pass3_and]
    paddb xmm5, oword ptr [pass3_add]
    pshufb xmm4, xmm5

    movdqa xmm5, xmm4
    pshufb xmm5, oword ptr [pass4_shuffle]
    pcmpgtb xmm5, xmm4
    pand xmm5, oword ptr [pass4_and]
    paddb xmm5, oword ptr [pass4_add]
    pshufb xmm4, xmm5

    movdqa xmm5, xmm4
    pshufb xmm5, oword ptr [pass5_shuffle]
    pcmpgtb xmm5, xmm4
    pand xmm5, oword ptr [pass5_and]
    paddb xmm5, oword ptr [pass5_add]
    pshufb xmm4, xmm5

    movdqa xmm5, xmm4
    pshufb xmm5, oword ptr [pass6_shuffle]
    pcmpgtb xmm5, xmm4
    paddb xmm5, oword ptr [pass6_add]
    pshufb xmm4, xmm5

    ;pextrd dword ptr [rcx], xmm4, 0    ; benchmarked with this
    ;pextrb byte ptr [rcx + 4], xmm4, 4 ; slower version
    ;pextrb byte ptr [rcx + 5], xmm4, 5
    movd   dword ptr [rcx], xmm4
    pextrw  word ptr [rcx + 4], xmm4, 2  ; x86 is little-endian, so this is the right order

    ret

simd_sort_6 ENDP

END

您可以将其编译为可执行对象,并将其链接到您的C项目中。有关如何在Visual Studio中执行此操作的说明,您可以阅读这篇文章。你可以使用下面的C原型从你的C代码中调用这个函数:

void simd_sort_6(char *values);