xref: /openbsd-src/lib/libc/stdlib/qsort.c (revision 4520aa4f05703eeee22781f0c24f3cdcd1421431)
1 /*	$OpenBSD: qsort.c,v 1.16 2017/05/20 12:48:56 millert Exp $ */
2 /*-
3  * Copyright (c) 1992, 1993
4  *	The Regents of the University of California.  All rights reserved.
5  *
6  * Redistribution and use in source and binary forms, with or without
7  * modification, are permitted provided that the following conditions
8  * are met:
9  * 1. Redistributions of source code must retain the above copyright
10  *    notice, this list of conditions and the following disclaimer.
11  * 2. Redistributions in binary form must reproduce the above copyright
12  *    notice, this list of conditions and the following disclaimer in the
13  *    documentation and/or other materials provided with the distribution.
14  * 3. Neither the name of the University nor the names of its contributors
15  *    may be used to endorse or promote products derived from this software
16  *    without specific prior written permission.
17  *
18  * THIS SOFTWARE IS PROVIDED BY THE REGENTS AND CONTRIBUTORS ``AS IS'' AND
19  * ANY EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT LIMITED TO, THE
20  * IMPLIED WARRANTIES OF MERCHANTABILITY AND FITNESS FOR A PARTICULAR PURPOSE
21  * ARE DISCLAIMED.  IN NO EVENT SHALL THE REGENTS OR CONTRIBUTORS BE LIABLE
22  * FOR ANY DIRECT, INDIRECT, INCIDENTAL, SPECIAL, EXEMPLARY, OR CONSEQUENTIAL
23  * DAMAGES (INCLUDING, BUT NOT LIMITED TO, PROCUREMENT OF SUBSTITUTE GOODS
24  * OR SERVICES; LOSS OF USE, DATA, OR PROFITS; OR BUSINESS INTERRUPTION)
25  * HOWEVER CAUSED AND ON ANY THEORY OF LIABILITY, WHETHER IN CONTRACT, STRICT
26  * LIABILITY, OR TORT (INCLUDING NEGLIGENCE OR OTHERWISE) ARISING IN ANY WAY
27  * OUT OF THE USE OF THIS SOFTWARE, EVEN IF ADVISED OF THE POSSIBILITY OF
28  * SUCH DAMAGE.
29  */
30 
31 #include <sys/types.h>
32 #include <stdlib.h>
33 
34 static __inline char	*med3(char *, char *, char *, int (*)(const void *, const void *));
35 static __inline void	 swapfunc(char *, char *, size_t, int);
36 
37 #define min(a, b)	(a) < (b) ? a : b
38 
39 /*
40  * Qsort routine from Bentley & McIlroy's "Engineering a Sort Function".
41  *
42  * This version differs from Bentley & McIlroy in the following ways:
43  *   1. The partition value is swapped into a[0] instead of being
44  *	stored out of line.
45  *
46  *   2. It uses David Musser's introsort algorithm to fall back to
47  *	heapsort(3) when the recursion depth reaches 2*lg(n + 1).
48  *	This avoids quicksort's quadratic behavior for pathological
49  *	input without appreciably changing the average run time.
50  *
51  *   3. Tail recursion is eliminated when sorting the larger of two
52  *	subpartitions to save stack space.
53  */
54 #define swapcode(TYPE, parmi, parmj, n) { 		\
55 	size_t i = (n) / sizeof (TYPE); 		\
56 	TYPE *pi = (TYPE *) (parmi); 			\
57 	TYPE *pj = (TYPE *) (parmj); 			\
58 	do { 						\
59 		TYPE	t = *pi;			\
60 		*pi++ = *pj;				\
61 		*pj++ = t;				\
62         } while (--i > 0);				\
63 }
64 
65 #define SWAPINIT(a, es) swaptype = ((char *)a - (char *)0) % sizeof(long) || \
66 	es % sizeof(long) ? 2 : es == sizeof(long)? 0 : 1;
67 
68 static __inline void
69 swapfunc(char *a, char *b, size_t n, int swaptype)
70 {
71 	if (swaptype <= 1)
72 		swapcode(long, a, b, n)
73 	else
74 		swapcode(char, a, b, n)
75 }
76 
77 #define swap(a, b)					\
78 	if (swaptype == 0) {				\
79 		long t = *(long *)(a);			\
80 		*(long *)(a) = *(long *)(b);		\
81 		*(long *)(b) = t;			\
82 	} else						\
83 		swapfunc(a, b, es, swaptype)
84 
85 #define vecswap(a, b, n) 	if ((n) > 0) swapfunc(a, b, n, swaptype)
86 
87 static __inline char *
88 med3(char *a, char *b, char *c, int (*cmp)(const void *, const void *))
89 {
90 	return cmp(a, b) < 0 ?
91 	       (cmp(b, c) < 0 ? b : (cmp(a, c) < 0 ? c : a ))
92               :(cmp(b, c) > 0 ? b : (cmp(a, c) < 0 ? a : c ));
93 }
94 
95 static void
96 introsort(char *a, size_t n, size_t es, size_t maxdepth,
97     int (*cmp)(const void *, const void *))
98 {
99 	char *pa, *pb, *pc, *pd, *pl, *pm, *pn;
100 	int cmp_result, swaptype;
101 	size_t r, s;
102 
103 loop:	if (maxdepth == 0) {
104 		if (heapsort(a, n, es, cmp) == 0)
105 			return;
106 	}
107 	maxdepth--;
108 	SWAPINIT(a, es);
109 	if (n < 7) {
110 		for (pm = a + es; pm < a + n * es; pm += es)
111 			for (pl = pm; pl > a && cmp(pl - es, pl) > 0;
112 			     pl -= es)
113 				swap(pl, pl - es);
114 		return;
115 	}
116 	pm = a + (n / 2) * es;
117 	if (n > 7) {
118 		pl = a;
119 		pn = a + (n - 1) * es;
120 		if (n > 40) {
121 			s = (n / 8) * es;
122 			pl = med3(pl, pl + s, pl + 2 * s, cmp);
123 			pm = med3(pm - s, pm, pm + s, cmp);
124 			pn = med3(pn - 2 * s, pn - s, pn, cmp);
125 		}
126 		pm = med3(pl, pm, pn, cmp);
127 	}
128 	swap(a, pm);
129 	pa = pb = a + es;
130 	pc = pd = a + (n - 1) * es;
131 	for (;;) {
132 		while (pb <= pc && (cmp_result = cmp(pb, a)) <= 0) {
133 			if (cmp_result == 0) {
134 				swap(pa, pb);
135 				pa += es;
136 			}
137 			pb += es;
138 		}
139 		while (pb <= pc && (cmp_result = cmp(pc, a)) >= 0) {
140 			if (cmp_result == 0) {
141 				swap(pc, pd);
142 				pd -= es;
143 			}
144 			pc -= es;
145 		}
146 		if (pb > pc)
147 			break;
148 		swap(pb, pc);
149 		pb += es;
150 		pc -= es;
151 	}
152 
153 	pn = a + n * es;
154 	r = min(pa - a, pb - pa);
155 	vecswap(a, pb - r, r);
156 	r = min(pd - pc, pn - pd - es);
157 	vecswap(pb, pn - r, r);
158 	/*
159 	 * To save stack space we sort the smaller side of the partition first
160 	 * using recursion and eliminate tail recursion for the larger side.
161 	 */
162 	r = pb - pa;
163 	s = pd - pc;
164 	if (r < s) {
165 		/* Recurse for 1st side, iterate for 2nd side. */
166 		if (s > es) {
167 			if (r > es)
168 				introsort(a, r / es, es, maxdepth, cmp);
169 			a = pn - s;
170 			n = s / es;
171 			goto loop;
172 		}
173 	} else {
174 		/* Recurse for 2nd side, iterate for 1st side. */
175 		if (r > es) {
176 			if (s > es)
177 				introsort(pn - s, s / es, es, maxdepth, cmp);
178 			n = r / es;
179 			goto loop;
180 		}
181 	}
182 }
183 
184 void
185 qsort(void *a, size_t n, size_t es, int (*cmp)(const void *, const void *))
186 {
187 	size_t i, maxdepth = 0;
188 
189 	/* Approximate 2*ceil(lg(n + 1)) */
190 	for (i = n; i > 0; i >>= 1)
191 		maxdepth++;
192 	maxdepth *= 2;
193 
194 	introsort(a, n, es, maxdepth, cmp);
195 }
196 
197 DEF_STRONG(qsort);
198