3131#include <string.h>
3232
3333#define MIN (x ,y ) (((x)<(y)) ? (x) : (y))
34+ #define MEDIAN3 (a ,b ,c ) (((a)<(b)) ? \
35+ ((b)<(c) ? (b) : ((a)<(c) ? (c) : (a))) : \
36+ ((b)>(c) ? (b) : ((a)>(c) ? (c) : (a))))
3437
3538static void split (int64_t * I ,int64_t * V ,int64_t start ,int64_t len ,int64_t h )
3639{
37- int64_t i ,j ,k ,x ,tmp ,jj ,kk ;
40+ int64_t i ,j ,k ,x ,y , z , tmp ,jj ,kk ;
3841
3942 if (len < 16 ) {
4043 for (k = start ;k < start + len ;k += j ) {
@@ -55,7 +58,20 @@ static void split(int64_t *I,int64_t *V,int64_t start,int64_t len,int64_t h)
5558 return ;
5659 };
5760
58- x = V [I [start + len /2 ]+ h ];
61+ /* Select pivot, algorithm by Bentley & McIlroy */
62+ j = start + len /2 ;
63+ k = start + len - 1 ;
64+ x = V [I [j ]+ h ];
65+ y = V [I [start ]+ h ];
66+ z = V [I [k ]+ h ];
67+ if (len > 40 ) { /* Big array: Pseudomedian of 9 */
68+ tmp = len /8 ;
69+ x = MEDIAN3 (x ,V [I [j - tmp ]+ h ],V [I [j + tmp ]+ h ]);
70+ y = MEDIAN3 (y ,V [I [start + tmp ]+ h ],V [I [start + tmp + tmp ]+ h ]);
71+ z = MEDIAN3 (z ,V [I [k - tmp ]+ h ],V [I [k - tmp - tmp ]+ h ]);
72+ }; /* Else medium array: Pseudomedian of 3 */
73+ x = MEDIAN3 (x ,y ,z );
74+
5975 jj = 0 ;kk = 0 ;
6076 for (i = start ;i < start + len ;i ++ ) {
6177 if (V [I [i ]+ h ]< x ) jj ++ ;
0 commit comments