1/* Copyright (C) 2011 by Valentin Ochs
2 *
3 * Permission is hereby granted, free of charge, to any person obtaining a copy
4 * of this software and associated documentation files (the "Software"), to
5 * deal in the Software without restriction, including without limitation the
6 * rights to use, copy, modify, merge, publish, distribute, sublicense, and/or
7 * sell copies of the Software, and to permit persons to whom the Software is
8 * furnished to do so, subject to the following conditions:
9 *
10 * The above copyright notice and this permission notice shall be included in
11 * all copies or substantial portions of the Software.
12 *
13 * THE SOFTWARE IS PROVIDED "AS IS", WITHOUT WARRANTY OF ANY KIND, EXPRESS OR
14 * IMPLIED, INCLUDING BUT NOT LIMITED TO THE WARRANTIES OF MERCHANTABILITY,
15 * FITNESS FOR A PARTICULAR PURPOSE AND NONINFRINGEMENT. IN NO EVENT SHALL THE
16 * AUTHORS OR COPYRIGHT HOLDERS BE LIABLE FOR ANY CLAIM, DAMAGES OR OTHER
17 * LIABILITY, WHETHER IN AN ACTION OF CONTRACT, TORT OR OTHERWISE, ARISING
18 * FROM, OUT OF OR IN CONNECTION WITH THE SOFTWARE OR THE USE OR OTHER DEALINGS
19 * IN THE SOFTWARE.
20 */
21
22/* Minor changes by Rich Felker for integration in musl, 2011-04-27. */
23
24#include <stdint.h>
25#include <stdlib.h>
26#include <string.h>
27
28#include "atomic.h"
29#define ntz(x) a_ctz_l((x))
30
31typedef int (*cmpfun)(const void *, const void *);
32
33static inline int pntz(size_t p[2]) {
34	int r = ntz(p[0] - 1);
35	if(r != 0 || (r = 8*sizeof(size_t) + ntz(p[1])) != 8*sizeof(size_t)) {
36		return r;
37	}
38	return 0;
39}
40
41static void cycle(size_t width, unsigned char* ar[], int n)
42{
43	unsigned char tmp[256];
44	size_t l;
45	int i;
46
47	if(n < 2) {
48		return;
49	}
50
51	ar[n] = tmp;
52	while(width) {
53		l = sizeof(tmp) < width ? sizeof(tmp) : width;
54		memcpy(ar[n], ar[0], l);
55		for(i = 0; i < n; i++) {
56			memcpy(ar[i], ar[i + 1], l);
57			ar[i] += l;
58		}
59		width -= l;
60	}
61}
62
63/* shl() and shr() need n > 0 */
64static inline void shl(size_t p[2], int n)
65{
66	if(n >= 8 * sizeof(size_t)) {
67		n -= 8 * sizeof(size_t);
68		p[1] = p[0];
69		p[0] = 0;
70	}
71	p[1] <<= n;
72	p[1] |= p[0] >> (sizeof(size_t) * 8 - n);
73	p[0] <<= n;
74}
75
76static inline void shr(size_t p[2], int n)
77{
78	if(n >= 8 * sizeof(size_t)) {
79		n -= 8 * sizeof(size_t);
80		p[0] = p[1];
81		p[1] = 0;
82	}
83	p[0] >>= n;
84	p[0] |= p[1] << (sizeof(size_t) * 8 - n);
85	p[1] >>= n;
86}
87
88static void sift(unsigned char *head, size_t width, cmpfun cmp, int pshift, size_t lp[])
89{
90	unsigned char *rt, *lf;
91	unsigned char *ar[14 * sizeof(size_t) + 1];
92	int i = 1;
93
94	ar[0] = head;
95	while(pshift > 1) {
96		rt = head - width;
97		lf = head - width - lp[pshift - 2];
98
99		if((*cmp)(ar[0], lf) >= 0 && (*cmp)(ar[0], rt) >= 0) {
100			break;
101		}
102		if((*cmp)(lf, rt) >= 0) {
103			ar[i++] = lf;
104			head = lf;
105			pshift -= 1;
106		} else {
107			ar[i++] = rt;
108			head = rt;
109			pshift -= 2;
110		}
111	}
112	cycle(width, ar, i);
113}
114
115static void trinkle(unsigned char *head, size_t width, cmpfun cmp, size_t pp[2], int pshift, int trusty, size_t lp[])
116{
117	unsigned char *stepson,
118	              *rt, *lf;
119	size_t p[2];
120	unsigned char *ar[14 * sizeof(size_t) + 1];
121	int i = 1;
122	int trail;
123
124	p[0] = pp[0];
125	p[1] = pp[1];
126
127	ar[0] = head;
128	while(p[0] != 1 || p[1] != 0) {
129		stepson = head - lp[pshift];
130		if((*cmp)(stepson, ar[0]) <= 0) {
131			break;
132		}
133		if(!trusty && pshift > 1) {
134			rt = head - width;
135			lf = head - width - lp[pshift - 2];
136			if((*cmp)(rt, stepson) >= 0 || (*cmp)(lf, stepson) >= 0) {
137				break;
138			}
139		}
140
141		ar[i++] = stepson;
142		head = stepson;
143		trail = pntz(p);
144		shr(p, trail);
145		pshift += trail;
146		trusty = 0;
147	}
148	if(!trusty) {
149		cycle(width, ar, i);
150		sift(head, width, cmp, pshift, lp);
151	}
152}
153
154void qsort(void *base, size_t nel, size_t width, cmpfun cmp)
155{
156	size_t lp[12*sizeof(size_t)];
157	size_t i, size = width * nel;
158	unsigned char *head, *high;
159	size_t p[2] = {1, 0};
160	int pshift = 1;
161	int trail;
162
163	if (!size) return;
164
165	head = base;
166	high = head + size - width;
167
168	/* Precompute Leonardo numbers, scaled by element width */
169	for(lp[0]=lp[1]=width, i=2; (lp[i]=lp[i-2]+lp[i-1]+width) < size; i++);
170
171	while(head < high) {
172		if((p[0] & 3) == 3) {
173			sift(head, width, cmp, pshift, lp);
174			shr(p, 2);
175			pshift += 2;
176		} else {
177			if(lp[pshift - 1] >= high - head) {
178				trinkle(head, width, cmp, p, pshift, 0, lp);
179			} else {
180				sift(head, width, cmp, pshift, lp);
181			}
182
183			if(pshift == 1) {
184				shl(p, 1);
185				pshift = 0;
186			} else {
187				shl(p, pshift - 1);
188				pshift = 1;
189			}
190		}
191
192		p[0] |= 1;
193		head += width;
194	}
195
196	trinkle(head, width, cmp, p, pshift, 0, lp);
197
198	while(pshift != 1 || p[0] != 1 || p[1] != 0) {
199		if(pshift <= 1) {
200			trail = pntz(p);
201			shr(p, trail);
202			pshift += trail;
203		} else {
204			shl(p, 2);
205			pshift -= 2;
206			p[0] ^= 7;
207			shr(p, 1);
208			trinkle(head - lp[pshift] - width, width, cmp, p, pshift + 1, 1, lp);
209			shl(p, 1);
210			p[0] |= 1;
211			trinkle(head - width, width, cmp, p, pshift, 1, lp);
212		}
213		head -= width;
214	}
215}
216