1 /*- 2 * Copyright (c) 1990 The Regents of the University of California. 3 * All rights reserved. 4 * 5 * %sccs.include.redist.c% 6 */ 7 8 #if defined(LIBC_SCCS) && !defined(lint) 9 static char sccsid[] = "@(#)radixsort.c 5.7 (Berkeley) 02/23/91"; 10 #endif /* LIBC_SCCS and not lint */ 11 12 #include <sys/types.h> 13 #include <limits.h> 14 #include <stdlib.h> 15 #include <stddef.h> 16 #include <string.h> 17 18 /* 19 * __rspartition is the cutoff point for a further partitioning instead 20 * of a shellsort. If it changes check __rsshell_increments. Both of 21 * these are exported, as the best values are data dependent. 22 */ 23 #define NPARTITION 40 24 int __rspartition = NPARTITION; 25 int __rsshell_increments[] = { 4, 1, 0, 0, 0, 0, 0, 0 }; 26 27 /* 28 * Stackp points to context structures, where each structure schedules a 29 * partitioning. Radixsort exits when the stack is empty. 30 * 31 * If the buckets are placed on the stack randomly, the worst case is when 32 * all the buckets but one contain (npartitions + 1) elements and the bucket 33 * pushed on the stack last contains the rest of the elements. In this case, 34 * stack growth is bounded by: 35 * 36 * limit = (nelements / (npartitions + 1)) - 1; 37 * 38 * This is a very large number, 52,377,648 for the maximum 32-bit signed int. 39 * 40 * By forcing the largest bucket to be pushed on the stack first, the worst 41 * case is when all but two buckets each contain (npartitions + 1) elements, 42 * with the remaining elements split equally between the first and last 43 * buckets pushed on the stack. In this case, stack growth is bounded when: 44 * 45 * for (partition_cnt = 0; nelements > npartitions; ++partition_cnt) 46 * nelements = 47 * (nelements - (npartitions + 1) * (nbuckets - 2)) / 2; 48 * The bound is: 49 * 50 * limit = partition_cnt * (nbuckets - 1); 51 * 52 * This is a much smaller number, 4590 for the maximum 32-bit signed int. 53 */ 54 #define NBUCKETS (UCHAR_MAX + 1) 55 56 typedef struct _stack { 57 const u_char **bot; 58 int indx, nmemb; 59 } CONTEXT; 60 61 #define STACKPUSH { \ 62 stackp->bot = p; \ 63 stackp->nmemb = nmemb; \ 64 stackp->indx = indx; \ 65 ++stackp; \ 66 } 67 #define STACKPOP { \ 68 if (stackp == stack) \ 69 break; \ 70 --stackp; \ 71 bot = stackp->bot; \ 72 nmemb = stackp->nmemb; \ 73 indx = stackp->indx; \ 74 } 75 76 /* 77 * A variant of MSD radix sorting; see Knuth Vol. 3, page 177, and 5.2.5, 78 * Ex. 10 and 12. Also, "Three Partition Refinement Algorithms, Paige 79 * and Tarjan, SIAM J. Comput. Vol. 16, No. 6, December 1987. 80 * 81 * This uses a simple sort as soon as a bucket crosses a cutoff point, 82 * rather than sorting the entire list after partitioning is finished. 83 * This should be an advantage. 84 * 85 * This is pure MSD instead of LSD of some number of MSD, switching to 86 * the simple sort as soon as possible. Takes linear time relative to 87 * the number of bytes in the strings. 88 */ 89 int 90 #if __STDC__ 91 radixsort(const u_char **l1, int nmemb, const u_char *tab, u_char endbyte) 92 #else 93 radixsort(l1, nmemb, tab, endbyte) 94 const u_char **l1; 95 register int nmemb; 96 const u_char *tab; 97 u_char endbyte; 98 #endif 99 { 100 register int i, indx, t1, t2; 101 register const u_char **l2; 102 register const u_char **p; 103 register const u_char **bot; 104 register const u_char *tr; 105 CONTEXT *stack, *stackp; 106 int c[NBUCKETS + 1], max; 107 u_char ltab[NBUCKETS]; 108 static void shellsort(); 109 110 if (nmemb <= 1) 111 return(0); 112 113 /* 114 * T1 is the constant part of the equation, the number of elements 115 * represented on the stack between the top and bottom entries. 116 * It doesn't get rounded as the divide by 2 rounds down (correct 117 * for a value being subtracted). T2, the nelem value, has to be 118 * rounded up before each divide because we want an upper bound; 119 * this could overflow if nmemb is the maximum int. 120 */ 121 t1 = ((__rspartition + 1) * (NBUCKETS - 2)) >> 1; 122 for (i = 0, t2 = nmemb; t2 > __rspartition; i += NBUCKETS - 1) 123 t2 = ((t2 + 1) >> 1) - t1; 124 if (i) { 125 if (!(stack = stackp = (CONTEXT *)malloc(i * sizeof(CONTEXT)))) 126 return(-1); 127 } else 128 stack = stackp = NULL; 129 130 /* 131 * There are two arrays, one provided by the user (l1), and the 132 * temporary one (l2). The data is sorted to the temporary stack, 133 * and then copied back. The speedup of using index to determine 134 * which stack the data is on and simply swapping stacks back and 135 * forth, thus avoiding the copy every iteration, turns out to not 136 * be any faster than the current implementation. 137 */ 138 if (!(l2 = (const u_char **)malloc(sizeof(u_char *) * nmemb))) 139 return(-1); 140 141 /* 142 * Tr references a table of sort weights; multiple entries may 143 * map to the same weight; EOS char must have the lowest weight. 144 */ 145 if (tab) 146 tr = tab; 147 else { 148 for (t1 = 0, t2 = endbyte; t1 < t2; ++t1) 149 ltab[t1] = t1 + 1; 150 ltab[t2] = 0; 151 for (t1 = endbyte + 1; t1 < NBUCKETS; ++t1) 152 ltab[t1] = t1; 153 tr = ltab; 154 } 155 156 /* First sort is entire stack */ 157 bot = l1; 158 indx = 0; 159 160 for (;;) { 161 /* Clear bucket count array */ 162 bzero((char *)c, sizeof(c)); 163 164 /* 165 * Compute number of items that sort to the same bucket 166 * for this index. 167 */ 168 for (p = bot, i = nmemb; --i >= 0;) 169 ++c[tr[(*p++)[indx]]]; 170 171 /* 172 * Sum the number of characters into c, dividing the temp 173 * stack into the right number of buckets for this bucket, 174 * this index. C contains the cumulative total of keys 175 * before and included in this bucket, and will later be 176 * used as an index to the bucket. c[NBUCKETS] contains 177 * the total number of elements, for determining how many 178 * elements the last bucket contains. At the same time 179 * find the largest bucket so it gets pushed first. 180 */ 181 for (i = max = t1 = 0, t2 = __rspartition; i <= NBUCKETS; ++i) { 182 if (c[i] > t2) { 183 t2 = c[i]; 184 max = i; 185 } 186 t1 = c[i] += t1; 187 } 188 189 /* 190 * Partition the elements into buckets; c decrements through 191 * the bucket, and ends up pointing to the first element of 192 * the bucket. 193 */ 194 for (i = nmemb; --i >= 0;) { 195 --p; 196 l2[--c[tr[(*p)[indx]]]] = *p; 197 } 198 199 /* Copy the partitioned elements back to user stack */ 200 bcopy(l2, bot, nmemb * sizeof(u_char *)); 201 202 ++indx; 203 /* 204 * Sort buckets as necessary; don't sort c[0], it's the 205 * EOS character bucket, and nothing can follow EOS. 206 */ 207 for (i = max; i; --i) { 208 if ((nmemb = c[i + 1] - (t1 = c[i])) < 2) 209 continue; 210 p = bot + t1; 211 if (nmemb > __rspartition) 212 STACKPUSH 213 else 214 shellsort(p, indx, nmemb, tr); 215 } 216 for (i = max + 1; i < NBUCKETS; ++i) { 217 if ((nmemb = c[i + 1] - (t1 = c[i])) < 2) 218 continue; 219 p = bot + t1; 220 if (nmemb > __rspartition) 221 STACKPUSH 222 else 223 shellsort(p, indx, nmemb, tr); 224 } 225 /* Break out when stack is empty */ 226 STACKPOP 227 } 228 229 free((char *)l2); 230 free((char *)stack); 231 return(0); 232 } 233 234 /* 235 * Shellsort (diminishing increment sort) from Data Structures and 236 * Algorithms, Aho, Hopcraft and Ullman, 1983 Edition, page 290; 237 * see also Knuth Vol. 3, page 84. The increments are selected from 238 * formula (8), page 95. Roughly O(N^3/2). 239 */ 240 static void 241 shellsort(p, indx, nmemb, tr) 242 register u_char **p, *tr; 243 register int indx, nmemb; 244 { 245 register u_char ch, *s1, *s2; 246 register int incr, *incrp, t1, t2; 247 248 for (incrp = __rsshell_increments; incr = *incrp++;) 249 for (t1 = incr; t1 < nmemb; ++t1) 250 for (t2 = t1 - incr; t2 >= 0;) { 251 s1 = p[t2] + indx; 252 s2 = p[t2 + incr] + indx; 253 while ((ch = tr[*s1++]) == tr[*s2] && ch) 254 ++s2; 255 if (ch > tr[*s2]) { 256 s1 = p[t2]; 257 p[t2] = p[t2 + incr]; 258 p[t2 + incr] = s1; 259 t2 -= incr; 260 } else 261 break; 262 } 263 } 264