Commit | Line | Data |
---|---|---|
34849bd9 KB |
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) | |
f0a345ab | 9 | static char sccsid[] = "@(#)radixsort.c 5.7 (Berkeley) %G%"; |
34849bd9 KB |
10 | #endif /* LIBC_SCCS and not lint */ |
11 | ||
12 | #include <sys/types.h> | |
34849bd9 KB |
13 | #include <limits.h> |
14 | #include <stdlib.h> | |
15 | #include <stddef.h> | |
f0a345ab | 16 | #include <string.h> |
34849bd9 | 17 | |
34849bd9 | 18 | /* |
34849bd9 KB |
19 | * __rspartition is the cutoff point for a further partitioning instead |
20 | * of a shellsort. If it changes check __rsshell_increments. Both of | |
33f8130f | 21 | * these are exported, as the best values are data dependent. |
34849bd9 KB |
22 | */ |
23 | #define NPARTITION 40 | |
24 | int __rspartition = NPARTITION; | |
25 | int __rsshell_increments[] = { 4, 1, 0, 0, 0, 0, 0, 0 }; | |
34849bd9 KB |
26 | |
27 | /* | |
dfad239c KB |
28 | * Stackp points to context structures, where each structure schedules a |
29 | * partitioning. Radixsort exits when the stack is empty. | |
34849bd9 | 30 | * |
a5baf581 | 31 | * If the buckets are placed on the stack randomly, the worst case is when |
e061e641 | 32 | * all the buckets but one contain (npartitions + 1) elements and the bucket |
a5baf581 KB |
33 | * pushed on the stack last contains the rest of the elements. In this case, |
34 | * stack growth is bounded by: | |
dfad239c | 35 | * |
e061e641 KB |
36 | * limit = (nelements / (npartitions + 1)) - 1; |
37 | * | |
38 | * This is a very large number, 52,377,648 for the maximum 32-bit signed int. | |
dfad239c | 39 | * |
e061e641 KB |
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: | |
dfad239c | 44 | * |
a5baf581 | 45 | * for (partition_cnt = 0; nelements > npartitions; ++partition_cnt) |
dfad239c KB |
46 | * nelements = |
47 | * (nelements - (npartitions + 1) * (nbuckets - 2)) / 2; | |
48 | * The bound is: | |
49 | * | |
50 | * limit = partition_cnt * (nbuckets - 1); | |
a5baf581 | 51 | * |
e061e641 | 52 | * This is a much smaller number, 4590 for the maximum 32-bit signed int. |
34849bd9 | 53 | */ |
e061e641 KB |
54 | #define NBUCKETS (UCHAR_MAX + 1) |
55 | ||
34849bd9 | 56 | typedef struct _stack { |
f0a345ab | 57 | const u_char **bot; |
34849bd9 KB |
58 | int indx, nmemb; |
59 | } CONTEXT; | |
60 | ||
34849bd9 | 61 | #define STACKPUSH { \ |
34849bd9 KB |
62 | stackp->bot = p; \ |
63 | stackp->nmemb = nmemb; \ | |
64 | stackp->indx = indx; \ | |
65 | ++stackp; \ | |
66 | } | |
34849bd9 KB |
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, | |
dfad239c KB |
78 | * Ex. 10 and 12. Also, "Three Partition Refinement Algorithms, Paige |
79 | * and Tarjan, SIAM J. Comput. Vol. 16, No. 6, December 1987. | |
34849bd9 | 80 | * |
dfad239c KB |
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. | |
34849bd9 | 84 | * |
dfad239c KB |
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. | |
34849bd9 | 88 | */ |
f0a345ab DS |
89 | int |
90 | #if __STDC__ | |
91 | radixsort(const u_char **l1, int nmemb, const u_char *tab, u_char endbyte) | |
92 | #else | |
34849bd9 | 93 | radixsort(l1, nmemb, tab, endbyte) |
f0a345ab | 94 | const u_char **l1; |
34849bd9 | 95 | register int nmemb; |
f0a345ab DS |
96 | const u_char *tab; |
97 | u_char endbyte; | |
98 | #endif | |
34849bd9 KB |
99 | { |
100 | register int i, indx, t1, t2; | |
f0a345ab DS |
101 | register const u_char **l2; |
102 | register const u_char **p; | |
103 | register const u_char **bot; | |
104 | register const u_char *tr; | |
dfad239c | 105 | CONTEXT *stack, *stackp; |
e061e641 KB |
106 | int c[NBUCKETS + 1], max; |
107 | u_char ltab[NBUCKETS]; | |
33f8130f | 108 | static void shellsort(); |
34849bd9 KB |
109 | |
110 | if (nmemb <= 1) | |
111 | return(0); | |
112 | ||
a5baf581 | 113 | /* |
dfad239c KB |
114 | * T1 is the constant part of the equation, the number of elements |
115 | * represented on the stack between the top and bottom entries. | |
a5baf581 KB |
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 | */ | |
e061e641 KB |
121 | t1 = ((__rspartition + 1) * (NBUCKETS - 2)) >> 1; |
122 | for (i = 0, t2 = nmemb; t2 > __rspartition; i += NBUCKETS - 1) | |
5fa271fa | 123 | t2 = ((t2 + 1) >> 1) - t1; |
dfad239c | 124 | if (i) { |
a5baf581 | 125 | if (!(stack = stackp = (CONTEXT *)malloc(i * sizeof(CONTEXT)))) |
dfad239c KB |
126 | return(-1); |
127 | } else | |
128 | stack = stackp = NULL; | |
129 | ||
34849bd9 | 130 | /* |
dfad239c | 131 | * There are two arrays, one provided by the user (l1), and the |
34849bd9 KB |
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 | */ | |
f0a345ab | 138 | if (!(l2 = (const u_char **)malloc(sizeof(u_char *) * nmemb))) |
34849bd9 KB |
139 | return(-1); |
140 | ||
34849bd9 | 141 | /* |
dfad239c | 142 | * Tr references a table of sort weights; multiple entries may |
34849bd9 KB |
143 | * map to the same weight; EOS char must have the lowest weight. |
144 | */ | |
145 | if (tab) | |
146 | tr = tab; | |
147 | else { | |
34849bd9 | 148 | for (t1 = 0, t2 = endbyte; t1 < t2; ++t1) |
f0a345ab DS |
149 | ltab[t1] = t1 + 1; |
150 | ltab[t2] = 0; | |
e061e641 | 151 | for (t1 = endbyte + 1; t1 < NBUCKETS; ++t1) |
f0a345ab DS |
152 | ltab[t1] = t1; |
153 | tr = ltab; | |
34849bd9 KB |
154 | } |
155 | ||
dfad239c | 156 | /* First sort is entire stack */ |
34849bd9 KB |
157 | bot = l1; |
158 | indx = 0; | |
159 | ||
160 | for (;;) { | |
dfad239c | 161 | /* Clear bucket count array */ |
34849bd9 KB |
162 | bzero((char *)c, sizeof(c)); |
163 | ||
164 | /* | |
dfad239c | 165 | * Compute number of items that sort to the same bucket |
34849bd9 KB |
166 | * for this index. |
167 | */ | |
33f8130f | 168 | for (p = bot, i = nmemb; --i >= 0;) |
34849bd9 KB |
169 | ++c[tr[(*p++)[indx]]]; |
170 | ||
171 | /* | |
dfad239c | 172 | * Sum the number of characters into c, dividing the temp |
34849bd9 KB |
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 | |
e061e641 | 176 | * used as an index to the bucket. c[NBUCKETS] contains |
34849bd9 | 177 | * the total number of elements, for determining how many |
dfad239c | 178 | * elements the last bucket contains. At the same time |
e061e641 | 179 | * find the largest bucket so it gets pushed first. |
34849bd9 | 180 | */ |
e061e641 KB |
181 | for (i = max = t1 = 0, t2 = __rspartition; i <= NBUCKETS; ++i) { |
182 | if (c[i] > t2) { | |
183 | t2 = c[i]; | |
dfad239c KB |
184 | max = i; |
185 | } | |
e061e641 | 186 | t1 = c[i] += t1; |
dfad239c | 187 | } |
34849bd9 KB |
188 | |
189 | /* | |
e061e641 KB |
190 | * Partition the elements into buckets; c decrements through |
191 | * the bucket, and ends up pointing to the first element of | |
192 | * the bucket. | |
34849bd9 | 193 | */ |
33f8130f | 194 | for (i = nmemb; --i >= 0;) { |
34849bd9 KB |
195 | --p; |
196 | l2[--c[tr[(*p)[indx]]]] = *p; | |
197 | } | |
198 | ||
dfad239c | 199 | /* Copy the partitioned elements back to user stack */ |
34849bd9 KB |
200 | bcopy(l2, bot, nmemb * sizeof(u_char *)); |
201 | ||
202 | ++indx; | |
203 | /* | |
dfad239c | 204 | * Sort buckets as necessary; don't sort c[0], it's the |
34849bd9 KB |
205 | * EOS character bucket, and nothing can follow EOS. |
206 | */ | |
dfad239c KB |
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 | |
33f8130f | 214 | shellsort(p, indx, nmemb, tr); |
dfad239c | 215 | } |
e061e641 | 216 | for (i = max + 1; i < NBUCKETS; ++i) { |
34849bd9 KB |
217 | if ((nmemb = c[i + 1] - (t1 = c[i])) < 2) |
218 | continue; | |
219 | p = bot + t1; | |
220 | if (nmemb > __rspartition) | |
221 | STACKPUSH | |
222 | else | |
33f8130f | 223 | shellsort(p, indx, nmemb, tr); |
34849bd9 | 224 | } |
dfad239c | 225 | /* Break out when stack is empty */ |
34849bd9 KB |
226 | STACKPOP |
227 | } | |
228 | ||
229 | free((char *)l2); | |
230 | free((char *)stack); | |
34849bd9 KB |
231 | return(0); |
232 | } | |
33f8130f KB |
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 | } |