#ifndef THIER_SORT2_H #define THIER_SORT2_H #include #include "qsort/schwab_sort.h" #include "threepass.h" /* A non-implace tricky float-hackz based bucket sort variant. Uses schwabsort or threepass! */ #ifdef _MSC_VER #define KM_PREFETCH(x) #define LIKELY(x) #define UNLIKELY(x) #define KM_NOINLINE __declspec(noinline) #define KM_ALWAYS_INLINE __forceinline #else #define KM_PREFETCH(x) __builtin_prefetch(x) #define LIKELY(x) __builtin_expect((x),1) /* alternative: #define LIKELY(x) __builtin_expect(((x) != 0),1) */ #define UNLIKELY(x) __builtin_expect((x),0) #define KM_NOINLINE __attribute__ ((noinline)) #define KM_ALWAYS_INLINE __attribute__ ((always_inline)) #endif /* Float and unsigned32 reinterpreter */ union th2_fu { float f; uint32_t u; }; typedef union th2_fu th2_fu; /** Tells from the key which bucket it is in. */ static inline uint32_t witch_bucket(uint32_t key) { /* This is hackz to misuse int->float converter HEAVILY and IEE for bucketing */ /* https://en.wikipedia.org/wiki/Single-precision_floating-point_format */ /* Old Hungarian ASM trick I know from Tomcat/Abaddon mailing list and prog.hu */ /* "A következő nagyon fontos gondolat, hogy két lebegőpontos számot, */ /* amennyiben pozitív, "simán" is összehasonlíthatunk" */ /* See: https://prog.hu/cikkek/100239/fpu-gems */ /* This approach uses 12 bits from a 32 bit float to map onto a byte bucket index */ th2_fu as; as.f = (float) key; uint32_t witch_base = (key <= 2) ? 0 : (as.u >> 23) - 128; // 0, [127..159] -> [0..31] // return witch_base * 8 + ((as.u >> (23 - 3)) & 7); // 0..255 return witch_base * 128 + ((as.u >> (23 - 7)) & 127); // 0..2047 /* Alternative (but I measure it being worse): return (as.u >> 23); */ } /** * Sort the array using the temporary array of the same size with fast bucket sort thiersort. * * @param arr The array to sort, will contain result afterwards * @param temparr The temporary array with same size * @param n Number of elements in arr and temparr * @param rstate Create with sch_rand_state rstate = schwab_rand_state(junk_uint32_t); */ static inline void thiersort2(uint32_t *arr, uint32_t *temparr, int n, sch_rand_state *rstate) { int bucket[4096]; /* Inclusive */ int bucket_end[4096]; /* Not inclusive */ /* Check if need to sort at all - needed for invariants later */ if(n < 2) { return; } /* Count */ #pragma GCC unroll 64 for(int i = 0; i < 4096; ++i) { bucket[i] = 0; } #pragma GCC unroll 64 for(int i = 0; i < n; ++i) { ++bucket[witch_bucket(arr[i])]; } /* Prefix sum (like in Magyarsort) */ uint32_t prev = 0; #pragma GCC unroll 4 for (int i = 0; i < 4096; i++) { bucket[i] += prev; prev = bucket[i]; } /* Save end-offsets */ #pragma GCC unroll 64 for(int i = 0; i < 4096; ++i) { bucket_end[i] = bucket[i]; } /* arr -> temparr */ /* Move to the buckets */ /* Rem.: This also changes bucket[i] so they will point to bucket beginnings */ #pragma GCC unroll 64 for(int i = 0; i < n; ++i) { uint32_t num = arr[i]; uint32_t witch = witch_bucket(num); int offset = (--bucket[witch]); temparr[offset] = num; } /* temparr -> arr each bucket and sort them in-place */ #pragma GCC unroll 64 for(int b = 0; b < 4096; ++b) { int begin = bucket[b]; int end = bucket_end[b]; /* Ensure exists */ if(begin >= end) { continue; } /* We make a three-way FAST quicksort partitioning with first elem pivot: */ /* Basically a Lomuto-like unidirectional partitioning for pivot and two-way for small/big */ int smalli = begin; int biggi = end - 1; /* always exists */ int i = begin; uint32_t pivot = temparr[i]; #pragma GCC unroll 4 for(int j = begin + 1; j < end; ++j) { if(LIKELY(temparr[j] != pivot)) { /* Branchless partitioning */ /* copy to left */ /* Rem.: Because we overwrite by copy we can simplify this line */ /* arr[smalli] = (temparr[j] < pivot) ? temparr[j] : arr[smalli]; */ arr[smalli] = temparr[j]; smalli += (temparr[j] < pivot); /* copy to right */ /* Rem.: Because we overwrite by copy we can simplify this line */ /* arr[biggi] = (temparr[j] > pivot) ? temparr[j] : arr[biggi]; */ arr[biggi] = temparr[j]; biggi -= (temparr[j] > pivot); } else { /* swap to front partition */ ++i; uint32_t tmp = temparr[i]; temparr[i] = temparr[j]; temparr[j] = tmp; } } /* Copy the mid elements back */ int target = smalli; #pragma GCC unroll 64 for(int j = begin; j < i + 1; ++j) { arr[target++] = temparr[j]; } /* Call schwabsort - only to [begin..smalli) and (biggie..end) */ #ifdef USE_SCHWAB schwab_sort(arr, begin, smalli - 1, rstate); schwab_sort(arr, biggi + 1, end - 1, rstate); #else threepass(&arr[begin], smalli - begin); threepass(&arr[biggi + 1], end - (biggi + 1)); #endif /* USE_SCHWAB */ } } #endif /* THIER_SORT2_H */