A proper, working, stable optimisation for sort {$b cmp $a}
[p5sagit/p5-mst-13.2.git] / pp_sort.c
1 /*    pp_sort.c
2  *
3  *    Copyright (C) 1991, 1992, 1993, 1994, 1995, 1996, 1997, 1998, 1999,
4  *    2000, 2001, 2002, 2003, 2004, by Larry Wall and others
5  *
6  *    You may distribute under the terms of either the GNU General Public
7  *    License or the Artistic License, as specified in the README file.
8  *
9  */
10
11 /*
12  *   ...they shuffled back towards the rear of the line. 'No, not at the
13  *   rear!'  the slave-driver shouted. 'Three files up. And stay there...
14  */
15
16 #include "EXTERN.h"
17 #define PERL_IN_PP_SORT_C
18 #include "perl.h"
19
20 #if defined(UNDER_CE)
21 /* looks like 'small' is reserved word for WINCE (or somesuch)*/
22 #define small xsmall
23 #endif
24
25 static I32 sortcv(pTHX_ SV *a, SV *b);
26 static I32 sortcv_stacked(pTHX_ SV *a, SV *b);
27 static I32 sortcv_xsub(pTHX_ SV *a, SV *b);
28 static I32 sv_ncmp(pTHX_ SV *a, SV *b);
29 static I32 sv_i_ncmp(pTHX_ SV *a, SV *b);
30 static I32 amagic_ncmp(pTHX_ SV *a, SV *b);
31 static I32 amagic_i_ncmp(pTHX_ SV *a, SV *b);
32 static I32 amagic_cmp(pTHX_ SV *a, SV *b);
33 static I32 amagic_cmp_locale(pTHX_ SV *a, SV *b);
34
35 #define sv_cmp_static Perl_sv_cmp
36 #define sv_cmp_locale_static Perl_sv_cmp_locale
37
38 #define SORTHINTS(hintsv) \
39     (((hintsv) = GvSV(gv_fetchpv("sort::hints", GV_ADDMULTI, SVt_IV))), \
40     (SvIOK(hintsv) ? ((I32)SvIV(hintsv)) : 0))
41
42 #ifndef SMALLSORT
43 #define SMALLSORT (200)
44 #endif
45
46 /*
47  * The mergesort implementation is by Peter M. Mcilroy <pmcilroy@lucent.com>.
48  *
49  * The original code was written in conjunction with BSD Computer Software
50  * Research Group at University of California, Berkeley.
51  *
52  * See also: "Optimistic Merge Sort" (SODA '92)
53  *
54  * The integration to Perl is by John P. Linderman <jpl@research.att.com>.
55  *
56  * The code can be distributed under the same terms as Perl itself.
57  *
58  */
59
60
61 typedef char * aptr;            /* pointer for arithmetic on sizes */
62 typedef SV * gptr;              /* pointers in our lists */
63
64 /* Binary merge internal sort, with a few special mods
65 ** for the special perl environment it now finds itself in.
66 **
67 ** Things that were once options have been hotwired
68 ** to values suitable for this use.  In particular, we'll always
69 ** initialize looking for natural runs, we'll always produce stable
70 ** output, and we'll always do Peter McIlroy's binary merge.
71 */
72
73 /* Pointer types for arithmetic and storage and convenience casts */
74
75 #define APTR(P) ((aptr)(P))
76 #define GPTP(P) ((gptr *)(P))
77 #define GPPP(P) ((gptr **)(P))
78
79
80 /* byte offset from pointer P to (larger) pointer Q */
81 #define BYTEOFF(P, Q) (APTR(Q) - APTR(P))
82
83 #define PSIZE sizeof(gptr)
84
85 /* If PSIZE is power of 2, make PSHIFT that power, if that helps */
86
87 #ifdef  PSHIFT
88 #define PNELEM(P, Q)    (BYTEOFF(P,Q) >> (PSHIFT))
89 #define PNBYTE(N)       ((N) << (PSHIFT))
90 #define PINDEX(P, N)    (GPTP(APTR(P) + PNBYTE(N)))
91 #else
92 /* Leave optimization to compiler */
93 #define PNELEM(P, Q)    (GPTP(Q) - GPTP(P))
94 #define PNBYTE(N)       ((N) * (PSIZE))
95 #define PINDEX(P, N)    (GPTP(P) + (N))
96 #endif
97
98 /* Pointer into other corresponding to pointer into this */
99 #define POTHER(P, THIS, OTHER) GPTP(APTR(OTHER) + BYTEOFF(THIS,P))
100
101 #define FROMTOUPTO(src, dst, lim) do *dst++ = *src++; while(src<lim)
102
103
104 /* Runs are identified by a pointer in the auxilliary list.
105 ** The pointer is at the start of the list,
106 ** and it points to the start of the next list.
107 ** NEXT is used as an lvalue, too.
108 */
109
110 #define NEXT(P)         (*GPPP(P))
111
112
113 /* PTHRESH is the minimum number of pairs with the same sense to justify
114 ** checking for a run and extending it.  Note that PTHRESH counts PAIRS,
115 ** not just elements, so PTHRESH == 8 means a run of 16.
116 */
117
118 #define PTHRESH (8)
119
120 /* RTHRESH is the number of elements in a run that must compare low
121 ** to the low element from the opposing run before we justify
122 ** doing a binary rampup instead of single stepping.
123 ** In random input, N in a row low should only happen with
124 ** probability 2^(1-N), so we can risk that we are dealing
125 ** with orderly input without paying much when we aren't.
126 */
127
128 #define RTHRESH (6)
129
130
131 /*
132 ** Overview of algorithm and variables.
133 ** The array of elements at list1 will be organized into runs of length 2,
134 ** or runs of length >= 2 * PTHRESH.  We only try to form long runs when
135 ** PTHRESH adjacent pairs compare in the same way, suggesting overall order.
136 **
137 ** Unless otherwise specified, pair pointers address the first of two elements.
138 **
139 ** b and b+1 are a pair that compare with sense ``sense''.
140 ** b is the ``bottom'' of adjacent pairs that might form a longer run.
141 **
142 ** p2 parallels b in the list2 array, where runs are defined by
143 ** a pointer chain.
144 **
145 ** t represents the ``top'' of the adjacent pairs that might extend
146 ** the run beginning at b.  Usually, t addresses a pair
147 ** that compares with opposite sense from (b,b+1).
148 ** However, it may also address a singleton element at the end of list1,
149 ** or it may be equal to ``last'', the first element beyond list1.
150 **
151 ** r addresses the Nth pair following b.  If this would be beyond t,
152 ** we back it off to t.  Only when r is less than t do we consider the
153 ** run long enough to consider checking.
154 **
155 ** q addresses a pair such that the pairs at b through q already form a run.
156 ** Often, q will equal b, indicating we only are sure of the pair itself.
157 ** However, a search on the previous cycle may have revealed a longer run,
158 ** so q may be greater than b.
159 **
160 ** p is used to work back from a candidate r, trying to reach q,
161 ** which would mean b through r would be a run.  If we discover such a run,
162 ** we start q at r and try to push it further towards t.
163 ** If b through r is NOT a run, we detect the wrong order at (p-1,p).
164 ** In any event, after the check (if any), we have two main cases.
165 **
166 ** 1) Short run.  b <= q < p <= r <= t.
167 **      b through q is a run (perhaps trivial)
168 **      q through p are uninteresting pairs
169 **      p through r is a run
170 **
171 ** 2) Long run.  b < r <= q < t.
172 **      b through q is a run (of length >= 2 * PTHRESH)
173 **
174 ** Note that degenerate cases are not only possible, but likely.
175 ** For example, if the pair following b compares with opposite sense,
176 ** then b == q < p == r == t.
177 */
178
179
180 static IV
181 dynprep(pTHX_ gptr *list1, gptr *list2, size_t nmemb, SVCOMPARE_t cmp)
182 {
183     I32 sense;
184     register gptr *b, *p, *q, *t, *p2;
185     register gptr c, *last, *r;
186     gptr *savep;
187     IV runs = 0;
188
189     b = list1;
190     last = PINDEX(b, nmemb);
191     sense = (cmp(aTHX_ *b, *(b+1)) > 0);
192     for (p2 = list2; b < last; ) {
193         /* We just started, or just reversed sense.
194         ** Set t at end of pairs with the prevailing sense.
195         */
196         for (p = b+2, t = p; ++p < last; t = ++p) {
197             if ((cmp(aTHX_ *t, *p) > 0) != sense) break;
198         }
199         q = b;
200         /* Having laid out the playing field, look for long runs */
201         do {
202             p = r = b + (2 * PTHRESH);
203             if (r >= t) p = r = t;      /* too short to care about */
204             else {
205                 while (((cmp(aTHX_ *(p-1), *p) > 0) == sense) &&
206                        ((p -= 2) > q));
207                 if (p <= q) {
208                     /* b through r is a (long) run.
209                     ** Extend it as far as possible.
210                     */
211                     p = q = r;
212                     while (((p += 2) < t) &&
213                            ((cmp(aTHX_ *(p-1), *p) > 0) == sense)) q = p;
214                     r = p = q + 2;      /* no simple pairs, no after-run */
215                 }
216             }
217             if (q > b) {                /* run of greater than 2 at b */
218                 savep = p;
219                 p = q += 2;
220                 /* pick up singleton, if possible */
221                 if ((p == t) &&
222                     ((t + 1) == last) &&
223                     ((cmp(aTHX_ *(p-1), *p) > 0) == sense))
224                     savep = r = p = q = last;
225                 p2 = NEXT(p2) = p2 + (p - b); ++runs;
226                 if (sense) while (b < --p) {
227                     c = *b;
228                     *b++ = *p;
229                     *p = c;
230                 }
231                 p = savep;
232             }
233             while (q < p) {             /* simple pairs */
234                 p2 = NEXT(p2) = p2 + 2; ++runs;
235                 if (sense) {
236                     c = *q++;
237                     *(q-1) = *q;
238                     *q++ = c;
239                 } else q += 2;
240             }
241             if (((b = p) == t) && ((t+1) == last)) {
242                 NEXT(p2) = p2 + 1; ++runs;
243                 b++;
244             }
245             q = r;
246         } while (b < t);
247         sense = !sense;
248     }
249     return runs;
250 }
251
252
253 /* The original merge sort, in use since 5.7, was as fast as, or faster than,
254  * qsort on many platforms, but slower than qsort, conspicuously so,
255  * on others.  The most likely explanation was platform-specific
256  * differences in cache sizes and relative speeds.
257  *
258  * The quicksort divide-and-conquer algorithm guarantees that, as the
259  * problem is subdivided into smaller and smaller parts, the parts
260  * fit into smaller (and faster) caches.  So it doesn't matter how
261  * many levels of cache exist, quicksort will "find" them, and,
262  * as long as smaller is faster, take advanatge of them.
263  *
264  * By contrast, consider how the original mergesort algorithm worked.
265  * Suppose we have five runs (each typically of length 2 after dynprep).
266  * 
267  * pass               base                        aux
268  *  0              1 2 3 4 5
269  *  1                                           12 34 5
270  *  2                1234 5
271  *  3                                            12345
272  *  4                 12345
273  *
274  * Adjacent pairs are merged in "grand sweeps" through the input.
275  * This means, on pass 1, the records in runs 1 and 2 aren't revisited until
276  * runs 3 and 4 are merged and the runs from run 5 have been copied.
277  * The only cache that matters is one large enough to hold *all* the input.
278  * On some platforms, this may be many times slower than smaller caches.
279  *
280  * The following pseudo-code uses the same basic merge algorithm,
281  * but in a divide-and-conquer way.
282  *
283  * # merge $runs runs at offset $offset of list $list1 into $list2.
284  * # all unmerged runs ($runs == 1) originate in list $base.
285  * sub mgsort2 {
286  *     my ($offset, $runs, $base, $list1, $list2) = @_;
287  *
288  *     if ($runs == 1) {
289  *         if ($list1 is $base) copy run to $list2
290  *         return offset of end of list (or copy)
291  *     } else {
292  *         $off2 = mgsort2($offset, $runs-($runs/2), $base, $list2, $list1)
293  *         mgsort2($off2, $runs/2, $base, $list2, $list1)
294  *         merge the adjacent runs at $offset of $list1 into $list2
295  *         return the offset of the end of the merged runs
296  *     }
297  * }
298  * mgsort2(0, $runs, $base, $aux, $base);
299  *
300  * For our 5 runs, the tree of calls looks like 
301  *
302  *           5
303  *      3        2
304  *   2     1   1   1
305  * 1   1
306  *
307  * 1   2   3   4   5
308  *
309  * and the corresponding activity looks like
310  *
311  * copy runs 1 and 2 from base to aux
312  * merge runs 1 and 2 from aux to base
313  * (run 3 is where it belongs, no copy needed)
314  * merge runs 12 and 3 from base to aux
315  * (runs 4 and 5 are where they belong, no copy needed)
316  * merge runs 4 and 5 from base to aux
317  * merge runs 123 and 45 from aux to base
318  *
319  * Note that we merge runs 1 and 2 immediately after copying them,
320  * while they are still likely to be in fast cache.  Similarly,
321  * run 3 is merged with run 12 while it still may be lingering in cache.
322  * This implementation should therefore enjoy much of the cache-friendly
323  * behavior that quicksort does.  In addition, it does less copying
324  * than the original mergesort implementation (only runs 1 and 2 are copied)
325  * and the "balancing" of merges is better (merged runs comprise more nearly
326  * equal numbers of original runs).
327  *
328  * The actual cache-friendly implementation will use a pseudo-stack
329  * to avoid recursion, and will unroll processing of runs of length 2,
330  * but it is otherwise similar to the recursive implementation.
331  */
332
333 typedef struct {
334     IV  offset;         /* offset of 1st of 2 runs at this level */
335     IV  runs;           /* how many runs must be combined into 1 */
336 } off_runs;             /* pseudo-stack element */
337
338
339 static I32
340 cmp_desc(pTHX_ gptr a, gptr b)
341 {
342     return -PL_sort_RealCmp(aTHX_ a, b);
343 }
344
345 STATIC void
346 S_mergesortsv(pTHX_ gptr *base, size_t nmemb, SVCOMPARE_t cmp, U32 flags)
347 {
348     IV i, run, runs, offset;
349     I32 sense, level;
350     int iwhich;
351     register gptr *f1, *f2, *t, *b, *p, *tp2, *l1, *l2, *q;
352     gptr *aux, *list1, *list2;
353     gptr *p1;
354     gptr small[SMALLSORT];
355     gptr *which[3];
356     off_runs stack[60], *stackp;
357     SVCOMPARE_t savecmp;
358
359     if (nmemb <= 1) return;                     /* sorted trivially */
360
361     if (flags) {
362         savecmp = PL_sort_RealCmp;      /* Save current comparison routine, if any */
363         PL_sort_RealCmp = cmp;  /* Put comparison routine where cmp_desc can find it */
364         cmp = cmp_desc;
365     }
366
367     if (nmemb <= SMALLSORT) aux = small;        /* use stack for aux array */
368     else { New(799,aux,nmemb,gptr); }           /* allocate auxilliary array */
369     level = 0;
370     stackp = stack;
371     stackp->runs = dynprep(aTHX_ base, aux, nmemb, cmp);
372     stackp->offset = offset = 0;
373     which[0] = which[2] = base;
374     which[1] = aux;
375     for (;;) {
376         /* On levels where both runs have be constructed (stackp->runs == 0),
377          * merge them, and note the offset of their end, in case the offset
378          * is needed at the next level up.  Hop up a level, and,
379          * as long as stackp->runs is 0, keep merging.
380          */
381         if ((runs = stackp->runs) == 0) {
382             iwhich = level & 1;
383             list1 = which[iwhich];              /* area where runs are now */
384             list2 = which[++iwhich];            /* area for merged runs */
385             do {
386                 offset = stackp->offset;
387                 f1 = p1 = list1 + offset;               /* start of first run */
388                 p = tp2 = list2 + offset;       /* where merged run will go */
389                 t = NEXT(p);                    /* where first run ends */
390                 f2 = l1 = POTHER(t, list2, list1); /* ... on the other side */
391                 t = NEXT(t);                    /* where second runs ends */
392                 l2 = POTHER(t, list2, list1);   /* ... on the other side */
393                 offset = PNELEM(list2, t);
394                 while (f1 < l1 && f2 < l2) {
395                     /* If head 1 is larger than head 2, find ALL the elements
396                     ** in list 2 strictly less than head1, write them all,
397                     ** then head 1.  Then compare the new heads, and repeat,
398                     ** until one or both lists are exhausted.
399                     **
400                     ** In all comparisons (after establishing
401                     ** which head to merge) the item to merge
402                     ** (at pointer q) is the first operand of
403                     ** the comparison.  When we want to know
404                     ** if ``q is strictly less than the other'',
405                     ** we can't just do
406                     **    cmp(q, other) < 0
407                     ** because stability demands that we treat equality
408                     ** as high when q comes from l2, and as low when
409                     ** q was from l1.  So we ask the question by doing
410                     **    cmp(q, other) <= sense
411                     ** and make sense == 0 when equality should look low,
412                     ** and -1 when equality should look high.
413                     */
414
415
416                     if (cmp(aTHX_ *f1, *f2) <= 0) {
417                         q = f2; b = f1; t = l1;
418                         sense = -1;
419                     } else {
420                         q = f1; b = f2; t = l2;
421                         sense = 0;
422                     }
423
424
425                     /* ramp up
426                     **
427                     ** Leave t at something strictly
428                     ** greater than q (or at the end of the list),
429                     ** and b at something strictly less than q.
430                     */
431                     for (i = 1, run = 0 ;;) {
432                         if ((p = PINDEX(b, i)) >= t) {
433                             /* off the end */
434                             if (((p = PINDEX(t, -1)) > b) &&
435                                 (cmp(aTHX_ *q, *p) <= sense))
436                                  t = p;
437                             else b = p;
438                             break;
439                         } else if (cmp(aTHX_ *q, *p) <= sense) {
440                             t = p;
441                             break;
442                         } else b = p;
443                         if (++run >= RTHRESH) i += i;
444                     }
445
446
447                     /* q is known to follow b and must be inserted before t.
448                     ** Increment b, so the range of possibilities is [b,t).
449                     ** Round binary split down, to favor early appearance.
450                     ** Adjust b and t until q belongs just before t.
451                     */
452
453                     b++;
454                     while (b < t) {
455                         p = PINDEX(b, (PNELEM(b, t) - 1) / 2);
456                         if (cmp(aTHX_ *q, *p) <= sense) {
457                             t = p;
458                         } else b = p + 1;
459                     }
460
461
462                     /* Copy all the strictly low elements */
463
464                     if (q == f1) {
465                         FROMTOUPTO(f2, tp2, t);
466                         *tp2++ = *f1++;
467                     } else {
468                         FROMTOUPTO(f1, tp2, t);
469                         *tp2++ = *f2++;
470                     }
471                 }
472
473
474                 /* Run out remaining list */
475                 if (f1 == l1) {
476                        if (f2 < l2) FROMTOUPTO(f2, tp2, l2);
477                 } else              FROMTOUPTO(f1, tp2, l1);
478                 p1 = NEXT(p1) = POTHER(tp2, list2, list1);
479
480                 if (--level == 0) goto done;
481                 --stackp;
482                 t = list1; list1 = list2; list2 = t;    /* swap lists */
483             } while ((runs = stackp->runs) == 0);
484         }
485
486
487         stackp->runs = 0;               /* current run will finish level */
488         /* While there are more than 2 runs remaining,
489          * turn them into exactly 2 runs (at the "other" level),
490          * each made up of approximately half the runs.
491          * Stack the second half for later processing,
492          * and set about producing the first half now.
493          */
494         while (runs > 2) {
495             ++level;
496             ++stackp;
497             stackp->offset = offset;
498             runs -= stackp->runs = runs / 2;
499         }
500         /* We must construct a single run from 1 or 2 runs.
501          * All the original runs are in which[0] == base.
502          * The run we construct must end up in which[level&1].
503          */
504         iwhich = level & 1;
505         if (runs == 1) {
506             /* Constructing a single run from a single run.
507              * If it's where it belongs already, there's nothing to do.
508              * Otherwise, copy it to where it belongs.
509              * A run of 1 is either a singleton at level 0,
510              * or the second half of a split 3.  In neither event
511              * is it necessary to set offset.  It will be set by the merge
512              * that immediately follows.
513              */
514             if (iwhich) {       /* Belongs in aux, currently in base */
515                 f1 = b = PINDEX(base, offset);  /* where list starts */
516                 f2 = PINDEX(aux, offset);       /* where list goes */
517                 t = NEXT(f2);                   /* where list will end */
518                 offset = PNELEM(aux, t);        /* offset thereof */
519                 t = PINDEX(base, offset);       /* where it currently ends */
520                 FROMTOUPTO(f1, f2, t);          /* copy */
521                 NEXT(b) = t;                    /* set up parallel pointer */
522             } else if (level == 0) goto done;   /* single run at level 0 */
523         } else {
524             /* Constructing a single run from two runs.
525              * The merge code at the top will do that.
526              * We need only make sure the two runs are in the "other" array,
527              * so they'll end up in the correct array after the merge.
528              */
529             ++level;
530             ++stackp;
531             stackp->offset = offset;
532             stackp->runs = 0;   /* take care of both runs, trigger merge */
533             if (!iwhich) {      /* Merged runs belong in aux, copy 1st */
534                 f1 = b = PINDEX(base, offset);  /* where first run starts */
535                 f2 = PINDEX(aux, offset);       /* where it will be copied */
536                 t = NEXT(f2);                   /* where first run will end */
537                 offset = PNELEM(aux, t);        /* offset thereof */
538                 p = PINDEX(base, offset);       /* end of first run */
539                 t = NEXT(t);                    /* where second run will end */
540                 t = PINDEX(base, PNELEM(aux, t)); /* where it now ends */
541                 FROMTOUPTO(f1, f2, t);          /* copy both runs */
542                 NEXT(b) = p;                    /* paralled pointer for 1st */
543                 NEXT(p) = t;                    /* ... and for second */
544             }
545         }
546     }
547 done:
548     if (aux != small) Safefree(aux);    /* free iff allocated */
549     if (flags) {
550          PL_sort_RealCmp = savecmp;     /* Restore current comparison routine, if any */
551     }
552     return;
553 }
554
555 /*
556  * The quicksort implementation was derived from source code contributed
557  * by Tom Horsley.
558  *
559  * NOTE: this code was derived from Tom Horsley's qsort replacement
560  * and should not be confused with the original code.
561  */
562
563 /* Copyright (C) Tom Horsley, 1997. All rights reserved.
564
565    Permission granted to distribute under the same terms as perl which are
566    (briefly):
567
568     This program is free software; you can redistribute it and/or modify
569     it under the terms of either:
570
571         a) the GNU General Public License as published by the Free
572         Software Foundation; either version 1, or (at your option) any
573         later version, or
574
575         b) the "Artistic License" which comes with this Kit.
576
577    Details on the perl license can be found in the perl source code which
578    may be located via the www.perl.com web page.
579
580    This is the most wonderfulest possible qsort I can come up with (and
581    still be mostly portable) My (limited) tests indicate it consistently
582    does about 20% fewer calls to compare than does the qsort in the Visual
583    C++ library, other vendors may vary.
584
585    Some of the ideas in here can be found in "Algorithms" by Sedgewick,
586    others I invented myself (or more likely re-invented since they seemed
587    pretty obvious once I watched the algorithm operate for a while).
588
589    Most of this code was written while watching the Marlins sweep the Giants
590    in the 1997 National League Playoffs - no Braves fans allowed to use this
591    code (just kidding :-).
592
593    I realize that if I wanted to be true to the perl tradition, the only
594    comment in this file would be something like:
595
596    ...they shuffled back towards the rear of the line. 'No, not at the
597    rear!'  the slave-driver shouted. 'Three files up. And stay there...
598
599    However, I really needed to violate that tradition just so I could keep
600    track of what happens myself, not to mention some poor fool trying to
601    understand this years from now :-).
602 */
603
604 /* ********************************************************** Configuration */
605
606 #ifndef QSORT_ORDER_GUESS
607 #define QSORT_ORDER_GUESS 2     /* Select doubling version of the netBSD trick */
608 #endif
609
610 /* QSORT_MAX_STACK is the largest number of partitions that can be stacked up for
611    future processing - a good max upper bound is log base 2 of memory size
612    (32 on 32 bit machines, 64 on 64 bit machines, etc). In reality can
613    safely be smaller than that since the program is taking up some space and
614    most operating systems only let you grab some subset of contiguous
615    memory (not to mention that you are normally sorting data larger than
616    1 byte element size :-).
617 */
618 #ifndef QSORT_MAX_STACK
619 #define QSORT_MAX_STACK 32
620 #endif
621
622 /* QSORT_BREAK_EVEN is the size of the largest partition we should insertion sort.
623    Anything bigger and we use qsort. If you make this too small, the qsort
624    will probably break (or become less efficient), because it doesn't expect
625    the middle element of a partition to be the same as the right or left -
626    you have been warned).
627 */
628 #ifndef QSORT_BREAK_EVEN
629 #define QSORT_BREAK_EVEN 6
630 #endif
631
632 /* QSORT_PLAY_SAFE is the size of the largest partition we're willing
633    to go quadratic on.  We innoculate larger partitions against
634    quadratic behavior by shuffling them before sorting.  This is not
635    an absolute guarantee of non-quadratic behavior, but it would take
636    staggeringly bad luck to pick extreme elements as the pivot
637    from randomized data.
638 */
639 #ifndef QSORT_PLAY_SAFE
640 #define QSORT_PLAY_SAFE 255
641 #endif
642
643 /* ************************************************************* Data Types */
644
645 /* hold left and right index values of a partition waiting to be sorted (the
646    partition includes both left and right - right is NOT one past the end or
647    anything like that).
648 */
649 struct partition_stack_entry {
650    int left;
651    int right;
652 #ifdef QSORT_ORDER_GUESS
653    int qsort_break_even;
654 #endif
655 };
656
657 /* ******************************************************* Shorthand Macros */
658
659 /* Note that these macros will be used from inside the qsort function where
660    we happen to know that the variable 'elt_size' contains the size of an
661    array element and the variable 'temp' points to enough space to hold a
662    temp element and the variable 'array' points to the array being sorted
663    and 'compare' is the pointer to the compare routine.
664
665    Also note that there are very many highly architecture specific ways
666    these might be sped up, but this is simply the most generally portable
667    code I could think of.
668 */
669
670 /* Return < 0 == 0 or > 0 as the value of elt1 is < elt2, == elt2, > elt2
671 */
672 #define qsort_cmp(elt1, elt2) \
673    ((*compare)(aTHX_ array[elt1], array[elt2]))
674
675 #ifdef QSORT_ORDER_GUESS
676 #define QSORT_NOTICE_SWAP swapped++;
677 #else
678 #define QSORT_NOTICE_SWAP
679 #endif
680
681 /* swaps contents of array elements elt1, elt2.
682 */
683 #define qsort_swap(elt1, elt2) \
684    STMT_START { \
685       QSORT_NOTICE_SWAP \
686       temp = array[elt1]; \
687       array[elt1] = array[elt2]; \
688       array[elt2] = temp; \
689    } STMT_END
690
691 /* rotate contents of elt1, elt2, elt3 such that elt1 gets elt2, elt2 gets
692    elt3 and elt3 gets elt1.
693 */
694 #define qsort_rotate(elt1, elt2, elt3) \
695    STMT_START { \
696       QSORT_NOTICE_SWAP \
697       temp = array[elt1]; \
698       array[elt1] = array[elt2]; \
699       array[elt2] = array[elt3]; \
700       array[elt3] = temp; \
701    } STMT_END
702
703 /* ************************************************************ Debug stuff */
704
705 #ifdef QSORT_DEBUG
706
707 static void
708 break_here()
709 {
710    return; /* good place to set a breakpoint */
711 }
712
713 #define qsort_assert(t) (void)( (t) || (break_here(), 0) )
714
715 static void
716 doqsort_all_asserts(
717    void * array,
718    size_t num_elts,
719    size_t elt_size,
720    int (*compare)(const void * elt1, const void * elt2),
721    int pc_left, int pc_right, int u_left, int u_right)
722 {
723    int i;
724
725    qsort_assert(pc_left <= pc_right);
726    qsort_assert(u_right < pc_left);
727    qsort_assert(pc_right < u_left);
728    for (i = u_right + 1; i < pc_left; ++i) {
729       qsort_assert(qsort_cmp(i, pc_left) < 0);
730    }
731    for (i = pc_left; i < pc_right; ++i) {
732       qsort_assert(qsort_cmp(i, pc_right) == 0);
733    }
734    for (i = pc_right + 1; i < u_left; ++i) {
735       qsort_assert(qsort_cmp(pc_right, i) < 0);
736    }
737 }
738
739 #define qsort_all_asserts(PC_LEFT, PC_RIGHT, U_LEFT, U_RIGHT) \
740    doqsort_all_asserts(array, num_elts, elt_size, compare, \
741                  PC_LEFT, PC_RIGHT, U_LEFT, U_RIGHT)
742
743 #else
744
745 #define qsort_assert(t) ((void)0)
746
747 #define qsort_all_asserts(PC_LEFT, PC_RIGHT, U_LEFT, U_RIGHT) ((void)0)
748
749 #endif
750
751 /* ****************************************************************** qsort */
752
753 STATIC void /* the standard unstable (u) quicksort (qsort) */
754 S_qsortsvu(pTHX_ SV ** array, size_t num_elts, SVCOMPARE_t compare)
755 {
756    register SV * temp;
757
758    struct partition_stack_entry partition_stack[QSORT_MAX_STACK];
759    int next_stack_entry = 0;
760
761    int part_left;
762    int part_right;
763 #ifdef QSORT_ORDER_GUESS
764    int qsort_break_even;
765    int swapped;
766 #endif
767
768    /* Make sure we actually have work to do.
769    */
770    if (num_elts <= 1) {
771       return;
772    }
773
774    /* Innoculate large partitions against quadratic behavior */
775    if (num_elts > QSORT_PLAY_SAFE) {
776       register size_t n, j;
777       register SV **q;
778       for (n = num_elts, q = array; n > 1; ) {
779          j = (size_t)(n-- * Drand01());
780          temp = q[j];
781          q[j] = q[n];
782          q[n] = temp;
783       }
784    }
785
786    /* Setup the initial partition definition and fall into the sorting loop
787    */
788    part_left = 0;
789    part_right = (int)(num_elts - 1);
790 #ifdef QSORT_ORDER_GUESS
791    qsort_break_even = QSORT_BREAK_EVEN;
792 #else
793 #define qsort_break_even QSORT_BREAK_EVEN
794 #endif
795    for ( ; ; ) {
796       if ((part_right - part_left) >= qsort_break_even) {
797          /* OK, this is gonna get hairy, so lets try to document all the
798             concepts and abbreviations and variables and what they keep
799             track of:
800
801             pc: pivot chunk - the set of array elements we accumulate in the
802                 middle of the partition, all equal in value to the original
803                 pivot element selected. The pc is defined by:
804
805                 pc_left - the leftmost array index of the pc
806                 pc_right - the rightmost array index of the pc
807
808                 we start with pc_left == pc_right and only one element
809                 in the pivot chunk (but it can grow during the scan).
810
811             u:  uncompared elements - the set of elements in the partition
812                 we have not yet compared to the pivot value. There are two
813                 uncompared sets during the scan - one to the left of the pc
814                 and one to the right.
815
816                 u_right - the rightmost index of the left side's uncompared set
817                 u_left - the leftmost index of the right side's uncompared set
818
819                 The leftmost index of the left sides's uncompared set
820                 doesn't need its own variable because it is always defined
821                 by the leftmost edge of the whole partition (part_left). The
822                 same goes for the rightmost edge of the right partition
823                 (part_right).
824
825                 We know there are no uncompared elements on the left once we
826                 get u_right < part_left and no uncompared elements on the
827                 right once u_left > part_right. When both these conditions
828                 are met, we have completed the scan of the partition.
829
830                 Any elements which are between the pivot chunk and the
831                 uncompared elements should be less than the pivot value on
832                 the left side and greater than the pivot value on the right
833                 side (in fact, the goal of the whole algorithm is to arrange
834                 for that to be true and make the groups of less-than and
835                 greater-then elements into new partitions to sort again).
836
837             As you marvel at the complexity of the code and wonder why it
838             has to be so confusing. Consider some of the things this level
839             of confusion brings:
840
841             Once I do a compare, I squeeze every ounce of juice out of it. I
842             never do compare calls I don't have to do, and I certainly never
843             do redundant calls.
844
845             I also never swap any elements unless I can prove there is a
846             good reason. Many sort algorithms will swap a known value with
847             an uncompared value just to get things in the right place (or
848             avoid complexity :-), but that uncompared value, once it gets
849             compared, may then have to be swapped again. A lot of the
850             complexity of this code is due to the fact that it never swaps
851             anything except compared values, and it only swaps them when the
852             compare shows they are out of position.
853          */
854          int pc_left, pc_right;
855          int u_right, u_left;
856
857          int s;
858
859          pc_left = ((part_left + part_right) / 2);
860          pc_right = pc_left;
861          u_right = pc_left - 1;
862          u_left = pc_right + 1;
863
864          /* Qsort works best when the pivot value is also the median value
865             in the partition (unfortunately you can't find the median value
866             without first sorting :-), so to give the algorithm a helping
867             hand, we pick 3 elements and sort them and use the median value
868             of that tiny set as the pivot value.
869
870             Some versions of qsort like to use the left middle and right as
871             the 3 elements to sort so they can insure the ends of the
872             partition will contain values which will stop the scan in the
873             compare loop, but when you have to call an arbitrarily complex
874             routine to do a compare, its really better to just keep track of
875             array index values to know when you hit the edge of the
876             partition and avoid the extra compare. An even better reason to
877             avoid using a compare call is the fact that you can drop off the
878             edge of the array if someone foolishly provides you with an
879             unstable compare function that doesn't always provide consistent
880             results.
881
882             So, since it is simpler for us to compare the three adjacent
883             elements in the middle of the partition, those are the ones we
884             pick here (conveniently pointed at by u_right, pc_left, and
885             u_left). The values of the left, center, and right elements
886             are refered to as l c and r in the following comments.
887          */
888
889 #ifdef QSORT_ORDER_GUESS
890          swapped = 0;
891 #endif
892          s = qsort_cmp(u_right, pc_left);
893          if (s < 0) {
894             /* l < c */
895             s = qsort_cmp(pc_left, u_left);
896             /* if l < c, c < r - already in order - nothing to do */
897             if (s == 0) {
898                /* l < c, c == r - already in order, pc grows */
899                ++pc_right;
900                qsort_all_asserts(pc_left, pc_right, u_left + 1, u_right - 1);
901             } else if (s > 0) {
902                /* l < c, c > r - need to know more */
903                s = qsort_cmp(u_right, u_left);
904                if (s < 0) {
905                   /* l < c, c > r, l < r - swap c & r to get ordered */
906                   qsort_swap(pc_left, u_left);
907                   qsort_all_asserts(pc_left, pc_right, u_left + 1, u_right - 1);
908                } else if (s == 0) {
909                   /* l < c, c > r, l == r - swap c&r, grow pc */
910                   qsort_swap(pc_left, u_left);
911                   --pc_left;
912                   qsort_all_asserts(pc_left, pc_right, u_left + 1, u_right - 1);
913                } else {
914                   /* l < c, c > r, l > r - make lcr into rlc to get ordered */
915                   qsort_rotate(pc_left, u_right, u_left);
916                   qsort_all_asserts(pc_left, pc_right, u_left + 1, u_right - 1);
917                }
918             }
919          } else if (s == 0) {
920             /* l == c */
921             s = qsort_cmp(pc_left, u_left);
922             if (s < 0) {
923                /* l == c, c < r - already in order, grow pc */
924                --pc_left;
925                qsort_all_asserts(pc_left, pc_right, u_left + 1, u_right - 1);
926             } else if (s == 0) {
927                /* l == c, c == r - already in order, grow pc both ways */
928                --pc_left;
929                ++pc_right;
930                qsort_all_asserts(pc_left, pc_right, u_left + 1, u_right - 1);
931             } else {
932                /* l == c, c > r - swap l & r, grow pc */
933                qsort_swap(u_right, u_left);
934                ++pc_right;
935                qsort_all_asserts(pc_left, pc_right, u_left + 1, u_right - 1);
936             }
937          } else {
938             /* l > c */
939             s = qsort_cmp(pc_left, u_left);
940             if (s < 0) {
941                /* l > c, c < r - need to know more */
942                s = qsort_cmp(u_right, u_left);
943                if (s < 0) {
944                   /* l > c, c < r, l < r - swap l & c to get ordered */
945                   qsort_swap(u_right, pc_left);
946                   qsort_all_asserts(pc_left, pc_right, u_left + 1, u_right - 1);
947                } else if (s == 0) {
948                   /* l > c, c < r, l == r - swap l & c, grow pc */
949                   qsort_swap(u_right, pc_left);
950                   ++pc_right;
951                   qsort_all_asserts(pc_left, pc_right, u_left + 1, u_right - 1);
952                } else {
953                   /* l > c, c < r, l > r - rotate lcr into crl to order */
954                   qsort_rotate(u_right, pc_left, u_left);
955                   qsort_all_asserts(pc_left, pc_right, u_left + 1, u_right - 1);
956                }
957             } else if (s == 0) {
958                /* l > c, c == r - swap ends, grow pc */
959                qsort_swap(u_right, u_left);
960                --pc_left;
961                qsort_all_asserts(pc_left, pc_right, u_left + 1, u_right - 1);
962             } else {
963                /* l > c, c > r - swap ends to get in order */
964                qsort_swap(u_right, u_left);
965                qsort_all_asserts(pc_left, pc_right, u_left + 1, u_right - 1);
966             }
967          }
968          /* We now know the 3 middle elements have been compared and
969             arranged in the desired order, so we can shrink the uncompared
970             sets on both sides
971          */
972          --u_right;
973          ++u_left;
974          qsort_all_asserts(pc_left, pc_right, u_left, u_right);
975
976          /* The above massive nested if was the simple part :-). We now have
977             the middle 3 elements ordered and we need to scan through the
978             uncompared sets on either side, swapping elements that are on
979             the wrong side or simply shuffling equal elements around to get
980             all equal elements into the pivot chunk.
981          */
982
983          for ( ; ; ) {
984             int still_work_on_left;
985             int still_work_on_right;
986
987             /* Scan the uncompared values on the left. If I find a value
988                equal to the pivot value, move it over so it is adjacent to
989                the pivot chunk and expand the pivot chunk. If I find a value
990                less than the pivot value, then just leave it - its already
991                on the correct side of the partition. If I find a greater
992                value, then stop the scan.
993             */
994             while ((still_work_on_left = (u_right >= part_left))) {
995                s = qsort_cmp(u_right, pc_left);
996                if (s < 0) {
997                   --u_right;
998                } else if (s == 0) {
999                   --pc_left;
1000                   if (pc_left != u_right) {
1001                      qsort_swap(u_right, pc_left);
1002                   }
1003                   --u_right;
1004                } else {
1005                   break;
1006                }
1007                qsort_assert(u_right < pc_left);
1008                qsort_assert(pc_left <= pc_right);
1009                qsort_assert(qsort_cmp(u_right + 1, pc_left) <= 0);
1010                qsort_assert(qsort_cmp(pc_left, pc_right) == 0);
1011             }
1012
1013             /* Do a mirror image scan of uncompared values on the right
1014             */
1015             while ((still_work_on_right = (u_left <= part_right))) {
1016                s = qsort_cmp(pc_right, u_left);
1017                if (s < 0) {
1018                   ++u_left;
1019                } else if (s == 0) {
1020                   ++pc_right;
1021                   if (pc_right != u_left) {
1022                      qsort_swap(pc_right, u_left);
1023                   }
1024                   ++u_left;
1025                } else {
1026                   break;
1027                }
1028                qsort_assert(u_left > pc_right);
1029                qsort_assert(pc_left <= pc_right);
1030                qsort_assert(qsort_cmp(pc_right, u_left - 1) <= 0);
1031                qsort_assert(qsort_cmp(pc_left, pc_right) == 0);
1032             }
1033
1034             if (still_work_on_left) {
1035                /* I know I have a value on the left side which needs to be
1036                   on the right side, but I need to know more to decide
1037                   exactly the best thing to do with it.
1038                */
1039                if (still_work_on_right) {
1040                   /* I know I have values on both side which are out of
1041                      position. This is a big win because I kill two birds
1042                      with one swap (so to speak). I can advance the
1043                      uncompared pointers on both sides after swapping both
1044                      of them into the right place.
1045                   */
1046                   qsort_swap(u_right, u_left);
1047                   --u_right;
1048                   ++u_left;
1049                   qsort_all_asserts(pc_left, pc_right, u_left, u_right);
1050                } else {
1051                   /* I have an out of position value on the left, but the
1052                      right is fully scanned, so I "slide" the pivot chunk
1053                      and any less-than values left one to make room for the
1054                      greater value over on the right. If the out of position
1055                      value is immediately adjacent to the pivot chunk (there
1056                      are no less-than values), I can do that with a swap,
1057                      otherwise, I have to rotate one of the less than values
1058                      into the former position of the out of position value
1059                      and the right end of the pivot chunk into the left end
1060                      (got all that?).
1061                   */
1062                   --pc_left;
1063                   if (pc_left == u_right) {
1064                      qsort_swap(u_right, pc_right);
1065                      qsort_all_asserts(pc_left, pc_right-1, u_left, u_right-1);
1066                   } else {
1067                      qsort_rotate(u_right, pc_left, pc_right);
1068                      qsort_all_asserts(pc_left, pc_right-1, u_left, u_right-1);
1069                   }
1070                   --pc_right;
1071                   --u_right;
1072                }
1073             } else if (still_work_on_right) {
1074                /* Mirror image of complex case above: I have an out of
1075                   position value on the right, but the left is fully
1076                   scanned, so I need to shuffle things around to make room
1077                   for the right value on the left.
1078                */
1079                ++pc_right;
1080                if (pc_right == u_left) {
1081                   qsort_swap(u_left, pc_left);
1082                   qsort_all_asserts(pc_left+1, pc_right, u_left+1, u_right);
1083                } else {
1084                   qsort_rotate(pc_right, pc_left, u_left);
1085                   qsort_all_asserts(pc_left+1, pc_right, u_left+1, u_right);
1086                }
1087                ++pc_left;
1088                ++u_left;
1089             } else {
1090                /* No more scanning required on either side of partition,
1091                   break out of loop and figure out next set of partitions
1092                */
1093                break;
1094             }
1095          }
1096
1097          /* The elements in the pivot chunk are now in the right place. They
1098             will never move or be compared again. All I have to do is decide
1099             what to do with the stuff to the left and right of the pivot
1100             chunk.
1101
1102             Notes on the QSORT_ORDER_GUESS ifdef code:
1103
1104             1. If I just built these partitions without swapping any (or
1105                very many) elements, there is a chance that the elements are
1106                already ordered properly (being properly ordered will
1107                certainly result in no swapping, but the converse can't be
1108                proved :-).
1109
1110             2. A (properly written) insertion sort will run faster on
1111                already ordered data than qsort will.
1112
1113             3. Perhaps there is some way to make a good guess about
1114                switching to an insertion sort earlier than partition size 6
1115                (for instance - we could save the partition size on the stack
1116                and increase the size each time we find we didn't swap, thus
1117                switching to insertion sort earlier for partitions with a
1118                history of not swapping).
1119
1120             4. Naturally, if I just switch right away, it will make
1121                artificial benchmarks with pure ascending (or descending)
1122                data look really good, but is that a good reason in general?
1123                Hard to say...
1124          */
1125
1126 #ifdef QSORT_ORDER_GUESS
1127          if (swapped < 3) {
1128 #if QSORT_ORDER_GUESS == 1
1129             qsort_break_even = (part_right - part_left) + 1;
1130 #endif
1131 #if QSORT_ORDER_GUESS == 2
1132             qsort_break_even *= 2;
1133 #endif
1134 #if QSORT_ORDER_GUESS == 3
1135             int prev_break = qsort_break_even;
1136             qsort_break_even *= qsort_break_even;
1137             if (qsort_break_even < prev_break) {
1138                qsort_break_even = (part_right - part_left) + 1;
1139             }
1140 #endif
1141          } else {
1142             qsort_break_even = QSORT_BREAK_EVEN;
1143          }
1144 #endif
1145
1146          if (part_left < pc_left) {
1147             /* There are elements on the left which need more processing.
1148                Check the right as well before deciding what to do.
1149             */
1150             if (pc_right < part_right) {
1151                /* We have two partitions to be sorted. Stack the biggest one
1152                   and process the smallest one on the next iteration. This
1153                   minimizes the stack height by insuring that any additional
1154                   stack entries must come from the smallest partition which
1155                   (because it is smallest) will have the fewest
1156                   opportunities to generate additional stack entries.
1157                */
1158                if ((part_right - pc_right) > (pc_left - part_left)) {
1159                   /* stack the right partition, process the left */
1160                   partition_stack[next_stack_entry].left = pc_right + 1;
1161                   partition_stack[next_stack_entry].right = part_right;
1162 #ifdef QSORT_ORDER_GUESS
1163                   partition_stack[next_stack_entry].qsort_break_even = qsort_break_even;
1164 #endif
1165                   part_right = pc_left - 1;
1166                } else {
1167                   /* stack the left partition, process the right */
1168                   partition_stack[next_stack_entry].left = part_left;
1169                   partition_stack[next_stack_entry].right = pc_left - 1;
1170 #ifdef QSORT_ORDER_GUESS
1171                   partition_stack[next_stack_entry].qsort_break_even = qsort_break_even;
1172 #endif
1173                   part_left = pc_right + 1;
1174                }
1175                qsort_assert(next_stack_entry < QSORT_MAX_STACK);
1176                ++next_stack_entry;
1177             } else {
1178                /* The elements on the left are the only remaining elements
1179                   that need sorting, arrange for them to be processed as the
1180                   next partition.
1181                */
1182                part_right = pc_left - 1;
1183             }
1184          } else if (pc_right < part_right) {
1185             /* There is only one chunk on the right to be sorted, make it
1186                the new partition and loop back around.
1187             */
1188             part_left = pc_right + 1;
1189          } else {
1190             /* This whole partition wound up in the pivot chunk, so
1191                we need to get a new partition off the stack.
1192             */
1193             if (next_stack_entry == 0) {
1194                /* the stack is empty - we are done */
1195                break;
1196             }
1197             --next_stack_entry;
1198             part_left = partition_stack[next_stack_entry].left;
1199             part_right = partition_stack[next_stack_entry].right;
1200 #ifdef QSORT_ORDER_GUESS
1201             qsort_break_even = partition_stack[next_stack_entry].qsort_break_even;
1202 #endif
1203          }
1204       } else {
1205          /* This partition is too small to fool with qsort complexity, just
1206             do an ordinary insertion sort to minimize overhead.
1207          */
1208          int i;
1209          /* Assume 1st element is in right place already, and start checking
1210             at 2nd element to see where it should be inserted.
1211          */
1212          for (i = part_left + 1; i <= part_right; ++i) {
1213             int j;
1214             /* Scan (backwards - just in case 'i' is already in right place)
1215                through the elements already sorted to see if the ith element
1216                belongs ahead of one of them.
1217             */
1218             for (j = i - 1; j >= part_left; --j) {
1219                if (qsort_cmp(i, j) >= 0) {
1220                   /* i belongs right after j
1221                   */
1222                   break;
1223                }
1224             }
1225             ++j;
1226             if (j != i) {
1227                /* Looks like we really need to move some things
1228                */
1229                int k;
1230                temp = array[i];
1231                for (k = i - 1; k >= j; --k)
1232                   array[k + 1] = array[k];
1233                array[j] = temp;
1234             }
1235          }
1236
1237          /* That partition is now sorted, grab the next one, or get out
1238             of the loop if there aren't any more.
1239          */
1240
1241          if (next_stack_entry == 0) {
1242             /* the stack is empty - we are done */
1243             break;
1244          }
1245          --next_stack_entry;
1246          part_left = partition_stack[next_stack_entry].left;
1247          part_right = partition_stack[next_stack_entry].right;
1248 #ifdef QSORT_ORDER_GUESS
1249          qsort_break_even = partition_stack[next_stack_entry].qsort_break_even;
1250 #endif
1251       }
1252    }
1253
1254    /* Believe it or not, the array is sorted at this point! */
1255 }
1256
1257 /* Stabilize what is, presumably, an otherwise unstable sort method.
1258  * We do that by allocating (or having on hand) an array of pointers
1259  * that is the same size as the original array of elements to be sorted.
1260  * We initialize this parallel array with the addresses of the original
1261  * array elements.  This indirection can make you crazy.
1262  * Some pictures can help.  After initializing, we have
1263  *
1264  *  indir                  list1
1265  * +----+                 +----+
1266  * |    | --------------> |    | ------> first element to be sorted
1267  * +----+                 +----+
1268  * |    | --------------> |    | ------> second element to be sorted
1269  * +----+                 +----+
1270  * |    | --------------> |    | ------> third element to be sorted
1271  * +----+                 +----+
1272  *  ...
1273  * +----+                 +----+
1274  * |    | --------------> |    | ------> n-1st element to be sorted
1275  * +----+                 +----+
1276  * |    | --------------> |    | ------> n-th element to be sorted
1277  * +----+                 +----+
1278  *
1279  * During the sort phase, we leave the elements of list1 where they are,
1280  * and sort the pointers in the indirect array in the same order determined
1281  * by the original comparison routine on the elements pointed to.
1282  * Because we don't move the elements of list1 around through
1283  * this phase, we can break ties on elements that compare equal
1284  * using their address in the list1 array, ensuring stabilty.
1285  * This leaves us with something looking like
1286  *
1287  *  indir                  list1
1288  * +----+                 +----+
1289  * |    | --+       +---> |    | ------> first element to be sorted
1290  * +----+   |       |     +----+
1291  * |    | --|-------|---> |    | ------> second element to be sorted
1292  * +----+   |       |     +----+
1293  * |    | --|-------+ +-> |    | ------> third element to be sorted
1294  * +----+   |         |   +----+
1295  *  ...
1296  * +----+    | |   | |    +----+
1297  * |    | ---|-+   | +--> |    | ------> n-1st element to be sorted
1298  * +----+    |     |      +----+
1299  * |    | ---+     +----> |    | ------> n-th element to be sorted
1300  * +----+                 +----+
1301  *
1302  * where the i-th element of the indirect array points to the element
1303  * that should be i-th in the sorted array.  After the sort phase,
1304  * we have to put the elements of list1 into the places
1305  * dictated by the indirect array.
1306  */
1307
1308
1309 static I32
1310 cmpindir(pTHX_ gptr a, gptr b)
1311 {
1312     I32 sense;
1313     gptr *ap = (gptr *)a;
1314     gptr *bp = (gptr *)b;
1315
1316     if ((sense = PL_sort_RealCmp(aTHX_ *ap, *bp)) == 0)
1317          sense = (ap > bp) ? 1 : ((ap < bp) ? -1 : 0);
1318     return sense;
1319 }
1320
1321 static I32
1322 cmpindir_desc(pTHX_ gptr a, gptr b)
1323 {
1324     I32 sense;
1325     gptr *ap = (gptr *)a;
1326     gptr *bp = (gptr *)b;
1327
1328     /* Reverse the default */
1329     if ((sense = PL_sort_RealCmp(aTHX_ *ap, *bp)))
1330         return -sense;
1331     /* But don't reverse the stability test.  */
1332     return (ap > bp) ? 1 : ((ap < bp) ? -1 : 0);
1333
1334 }
1335
1336 STATIC void
1337 S_qsortsv(pTHX_ gptr *list1, size_t nmemb, SVCOMPARE_t cmp, U32 flags)
1338 {
1339     SV *hintsv;
1340
1341     if (SORTHINTS(hintsv) & HINT_SORT_STABLE) {
1342          register gptr **pp, *q;
1343          register size_t n, j, i;
1344          gptr *small[SMALLSORT], **indir, tmp;
1345          SVCOMPARE_t savecmp;
1346          if (nmemb <= 1) return;     /* sorted trivially */
1347
1348          /* Small arrays can use the stack, big ones must be allocated */
1349          if (nmemb <= SMALLSORT) indir = small;
1350          else { New(1799, indir, nmemb, gptr *); }
1351
1352          /* Copy pointers to original array elements into indirect array */
1353          for (n = nmemb, pp = indir, q = list1; n--; ) *pp++ = q++;
1354
1355          savecmp = PL_sort_RealCmp;     /* Save current comparison routine, if any */
1356          PL_sort_RealCmp = cmp; /* Put comparison routine where cmpindir can find it */
1357
1358          /* sort, with indirection */
1359          S_qsortsvu(aTHX_ (gptr *)indir, nmemb,
1360                     flags ? cmpindir_desc : cmpindir);
1361
1362          pp = indir;
1363          q = list1;
1364          for (n = nmemb; n--; ) {
1365               /* Assert A: all elements of q with index > n are already
1366                * in place.  This is vacuosly true at the start, and we
1367                * put element n where it belongs below (if it wasn't
1368                * already where it belonged). Assert B: we only move
1369                * elements that aren't where they belong,
1370                * so, by A, we never tamper with elements above n.
1371                */
1372               j = pp[n] - q;            /* This sets j so that q[j] is
1373                                          * at pp[n].  *pp[j] belongs in
1374                                          * q[j], by construction.
1375                                          */
1376               if (n != j) {             /* all's well if n == j */
1377                    tmp = q[j];          /* save what's in q[j] */
1378                    do {
1379                         q[j] = *pp[j];  /* put *pp[j] where it belongs */
1380                         i = pp[j] - q;  /* the index in q of the element
1381                                          * just moved */
1382                         pp[j] = q + j;  /* this is ok now */
1383                    } while ((j = i) != n);
1384                    /* There are only finitely many (nmemb) addresses
1385                     * in the pp array.
1386                     * So we must eventually revisit an index we saw before.
1387                     * Suppose the first revisited index is k != n.
1388                     * An index is visited because something else belongs there.
1389                     * If we visit k twice, then two different elements must
1390                     * belong in the same place, which cannot be.
1391                     * So j must get back to n, the loop terminates,
1392                     * and we put the saved element where it belongs.
1393                     */
1394                    q[n] = tmp;          /* put what belongs into
1395                                          * the n-th element */
1396               }
1397          }
1398
1399         /* free iff allocated */
1400          if (indir != small) { Safefree(indir); }
1401          /* restore prevailing comparison routine */
1402          PL_sort_RealCmp = savecmp;
1403     } else if (flags) {
1404          SVCOMPARE_t savecmp = PL_sort_RealCmp; /* Save current comparison routine, if any */
1405          PL_sort_RealCmp = cmp; /* Put comparison routine where cmp_desc can find it */
1406          cmp = cmp_desc;
1407          S_qsortsvu(aTHX_ list1, nmemb, cmp);
1408          /* restore prevailing comparison routine */
1409          PL_sort_RealCmp = savecmp;
1410     } else {
1411          S_qsortsvu(aTHX_ list1, nmemb, cmp);
1412     }
1413 }
1414
1415 /*
1416 =head1 Array Manipulation Functions
1417
1418 =for apidoc sortsv
1419
1420 Sort an array. Here is an example:
1421
1422     sortsv(AvARRAY(av), av_len(av)+1, Perl_sv_cmp_locale);
1423
1424 See lib/sort.pm for details about controlling the sorting algorithm.
1425
1426 =cut
1427 */
1428
1429 void
1430 Perl_sortsv(pTHX_ SV **array, size_t nmemb, SVCOMPARE_t cmp)
1431 {
1432     void (*sortsvp)(pTHX_ SV **array, size_t nmemb, SVCOMPARE_t cmp, U32 flags)
1433       = S_mergesortsv;
1434     SV *hintsv;
1435     I32 hints;
1436
1437     /*  Sun's Compiler (cc: WorkShop Compilers 4.2 30 Oct 1996 C 4.2) used 
1438         to miscompile this function under optimization -O.  If you get test 
1439         errors related to picking the correct sort() function, try recompiling 
1440         this file without optimiziation.  -- A.D.  4/2002.
1441     */
1442     hints = SORTHINTS(hintsv);
1443     if (hints & HINT_SORT_QUICKSORT) {
1444         sortsvp = S_qsortsv;
1445     }
1446     else {
1447         /* The default as of 5.8.0 is mergesort */
1448         sortsvp = S_mergesortsv;
1449     }
1450
1451     sortsvp(aTHX_ array, nmemb, cmp, 0);
1452 }
1453
1454
1455 void
1456 S_sortsv_desc(pTHX_ SV **array, size_t nmemb, SVCOMPARE_t cmp)
1457 {
1458     void (*sortsvp)(pTHX_ SV **array, size_t nmemb, SVCOMPARE_t cmp, U32 flags)
1459       = S_mergesortsv;
1460     SV *hintsv;
1461     I32 hints;
1462
1463     /*  Sun's Compiler (cc: WorkShop Compilers 4.2 30 Oct 1996 C 4.2) used 
1464         to miscompile this function under optimization -O.  If you get test 
1465         errors related to picking the correct sort() function, try recompiling 
1466         this file without optimiziation.  -- A.D.  4/2002.
1467     */
1468     hints = SORTHINTS(hintsv);
1469     if (hints & HINT_SORT_QUICKSORT) {
1470         sortsvp = S_qsortsv;
1471     }
1472     else {
1473         /* The default as of 5.8.0 is mergesort */
1474         sortsvp = S_mergesortsv;
1475     }
1476
1477     sortsvp(aTHX_ array, nmemb, cmp, 1);
1478 }
1479
1480 PP(pp_sort)
1481 {
1482     dSP; dMARK; dORIGMARK;
1483     register SV **p1 = ORIGMARK+1, **p2;
1484     register I32 max, i;
1485     AV* av = Nullav;
1486     HV *stash;
1487     GV *gv;
1488     CV *cv = 0;
1489     I32 gimme = GIMME;
1490     OP* nextop = PL_op->op_next;
1491     I32 overloading = 0;
1492     bool hasargs = FALSE;
1493     I32 is_xsub = 0;
1494     I32 sorting_av = 0;
1495     void (*sortsvp)(pTHX_ SV **array, size_t nmemb, SVCOMPARE_t cmp)
1496       = Perl_sortsv;
1497
1498     if (gimme != G_ARRAY) {
1499         SP = MARK;
1500         RETPUSHUNDEF;
1501     }
1502
1503     ENTER;
1504     SAVEVPTR(PL_sortcop);
1505     if (PL_op->op_flags & OPf_STACKED) {
1506         if (PL_op->op_flags & OPf_SPECIAL) {
1507             OP *kid = cLISTOP->op_first->op_sibling;    /* pass pushmark */
1508             kid = kUNOP->op_first;                      /* pass rv2gv */
1509             kid = kUNOP->op_first;                      /* pass leave */
1510             PL_sortcop = kid->op_next;
1511             stash = CopSTASH(PL_curcop);
1512         }
1513         else {
1514             cv = sv_2cv(*++MARK, &stash, &gv, 0);
1515             if (cv && SvPOK(cv)) {
1516                 STRLEN n_a;
1517                 char *proto = SvPV((SV*)cv, n_a);
1518                 if (proto && strEQ(proto, "$$")) {
1519                     hasargs = TRUE;
1520                 }
1521             }
1522             if (!(cv && CvROOT(cv))) {
1523                 if (cv && CvXSUB(cv)) {
1524                     is_xsub = 1;
1525                 }
1526                 else if (gv) {
1527                     SV *tmpstr = sv_newmortal();
1528                     gv_efullname3(tmpstr, gv, Nullch);
1529                     DIE(aTHX_ "Undefined sort subroutine \"%"SVf"\" called",
1530                         tmpstr);
1531                 }
1532                 else {
1533                     DIE(aTHX_ "Undefined subroutine in sort");
1534                 }
1535             }
1536
1537             if (is_xsub)
1538                 PL_sortcop = (OP*)cv;
1539             else {
1540                 PL_sortcop = CvSTART(cv);
1541                 SAVEVPTR(CvROOT(cv)->op_ppaddr);
1542                 CvROOT(cv)->op_ppaddr = PL_ppaddr[OP_NULL];
1543
1544                 PAD_SET_CUR(CvPADLIST(cv), 1);
1545             }
1546         }
1547     }
1548     else {
1549         PL_sortcop = Nullop;
1550         stash = CopSTASH(PL_curcop);
1551     }
1552
1553     /* optimiser converts "@a = sort @a" to "sort \@a";
1554      * in case of tied @a, pessimise: push (@a) onto stack, then assign
1555      * result back to @a at the end of this function */
1556     if (PL_op->op_private & OPpSORT_INPLACE) {
1557         assert( MARK+1 == SP && *SP && SvTYPE(*SP) == SVt_PVAV);
1558         (void)POPMARK; /* remove mark associated with ex-OP_AASSIGN */
1559         av = (AV*)(*SP);
1560         max = AvFILL(av) + 1;
1561         if (SvMAGICAL(av)) {
1562             MEXTEND(SP, max);
1563             p2 = SP;
1564             for (i=0; i < (U32)max; i++) {
1565                 SV **svp = av_fetch(av, i, FALSE);
1566                 *SP++ = (svp) ? *svp : Nullsv;
1567             }
1568         }
1569         else {
1570             p1 = p2 = AvARRAY(av);
1571             sorting_av = 1;
1572         }
1573     }
1574     else {
1575         p2 = MARK+1;
1576         max = SP - MARK;
1577    }
1578
1579     if (PL_op->op_private & OPpSORT_DESCEND) {
1580         sortsvp = S_sortsv_desc;
1581     }
1582
1583     /* shuffle stack down, removing optional initial cv (p1!=p2), plus any
1584      * nulls; also stringify any args */
1585     for (i=max; i > 0 ; i--) {
1586         if ((*p1 = *p2++)) {                    /* Weed out nulls. */
1587             SvTEMP_off(*p1);
1588             if (!PL_sortcop && !SvPOK(*p1)) {
1589                 STRLEN n_a;
1590                 if (SvAMAGIC(*p1))
1591                     overloading = 1;
1592                 else
1593                     (void)sv_2pv(*p1, &n_a);
1594             }
1595             p1++;
1596         }
1597         else
1598             max--;
1599     }
1600     if (sorting_av)
1601         AvFILLp(av) = max-1;
1602
1603     if (max > 1) {
1604         if (PL_sortcop) {
1605             PERL_CONTEXT *cx;
1606             SV** newsp;
1607             bool oldcatch = CATCH_GET;
1608
1609             SAVETMPS;
1610             SAVEOP();
1611
1612             CATCH_SET(TRUE);
1613             PUSHSTACKi(PERLSI_SORT);
1614             if (!hasargs && !is_xsub) {
1615                 if (PL_sortstash != stash || !PL_firstgv || !PL_secondgv) {
1616                     SAVESPTR(PL_firstgv);
1617                     SAVESPTR(PL_secondgv);
1618                     PL_firstgv = gv_fetchpv("a", TRUE, SVt_PV);
1619                     PL_secondgv = gv_fetchpv("b", TRUE, SVt_PV);
1620                     PL_sortstash = stash;
1621                 }
1622                 SAVESPTR(GvSV(PL_firstgv));
1623                 SAVESPTR(GvSV(PL_secondgv));
1624             }
1625
1626             PUSHBLOCK(cx, CXt_NULL, PL_stack_base);
1627             if (!(PL_op->op_flags & OPf_SPECIAL)) {
1628                 cx->cx_type = CXt_SUB;
1629                 cx->blk_gimme = G_SCALAR;
1630                 PUSHSUB(cx);
1631             }
1632             PL_sortcxix = cxstack_ix;
1633
1634             if (hasargs && !is_xsub) {
1635                 /* This is mostly copied from pp_entersub */
1636                 AV *av = (AV*)PAD_SVl(0);
1637
1638                 cx->blk_sub.savearray = GvAV(PL_defgv);
1639                 GvAV(PL_defgv) = (AV*)SvREFCNT_inc(av);
1640                 CX_CURPAD_SAVE(cx->blk_sub);
1641                 cx->blk_sub.argarray = av;
1642             }
1643            sortsvp(aTHX_ p1-max, max,
1644                    is_xsub ? sortcv_xsub : hasargs ? sortcv_stacked : sortcv);
1645
1646             POPBLOCK(cx,PL_curpm);
1647             PL_stack_sp = newsp;
1648             POPSTACK;
1649             CATCH_SET(oldcatch);
1650         }
1651         else {
1652             MEXTEND(SP, 20);    /* Can't afford stack realloc on signal. */
1653             sortsvp(aTHX_ sorting_av ? AvARRAY(av) : ORIGMARK+1, max,
1654                     (PL_op->op_private & OPpSORT_NUMERIC)
1655                         ? ( (PL_op->op_private & OPpSORT_INTEGER)
1656                             ? ( overloading ? amagic_i_ncmp : sv_i_ncmp)
1657                             : ( overloading ? amagic_ncmp : sv_ncmp))
1658                         : ( IN_LOCALE_RUNTIME
1659                             ? ( overloading
1660                                 ? amagic_cmp_locale
1661                                 : sv_cmp_locale_static)
1662                             : ( overloading ? amagic_cmp : sv_cmp_static)));
1663             if (PL_op->op_private & OPpSORT_REVERSE) {
1664                 SV **p = sorting_av ? AvARRAY(av) : ORIGMARK+1;
1665                 SV **q = p+max-1;
1666                 while (p < q) {
1667                     SV *tmp = *p;
1668                     *p++ = *q;
1669                     *q-- = tmp;
1670                 }
1671             }
1672         }
1673     }
1674     if (av && !sorting_av) {
1675         /* simulate pp_aassign of tied AV */
1676         SV *sv;
1677         SV** base, **didstore;
1678         for (base = ORIGMARK+1, i=0; i < max; i++) {
1679             sv = NEWSV(28,0);
1680             sv_setsv(sv, base[i]);
1681             base[i] = sv;
1682         }
1683         av_clear(av);
1684         av_extend(av, max);
1685         for (i=0; i < max; i++) {
1686             sv = base[i];
1687             didstore = av_store(av, i, sv);
1688             if (SvSMAGICAL(sv))
1689                 mg_set(sv);
1690             if (!didstore)
1691                 sv_2mortal(sv);
1692         }
1693     }
1694     LEAVE;
1695     PL_stack_sp = ORIGMARK + (sorting_av ? 0 : max);
1696     return nextop;
1697 }
1698
1699 static I32
1700 sortcv(pTHX_ SV *a, SV *b)
1701 {
1702     I32 oldsaveix = PL_savestack_ix;
1703     I32 oldscopeix = PL_scopestack_ix;
1704     I32 result;
1705     GvSV(PL_firstgv) = a;
1706     GvSV(PL_secondgv) = b;
1707     PL_stack_sp = PL_stack_base;
1708     PL_op = PL_sortcop;
1709     CALLRUNOPS(aTHX);
1710     if (PL_stack_sp != PL_stack_base + 1)
1711         Perl_croak(aTHX_ "Sort subroutine didn't return single value");
1712     if (!SvNIOKp(*PL_stack_sp))
1713         Perl_croak(aTHX_ "Sort subroutine didn't return a numeric value");
1714     result = SvIV(*PL_stack_sp);
1715     while (PL_scopestack_ix > oldscopeix) {
1716         LEAVE;
1717     }
1718     leave_scope(oldsaveix);
1719     return result;
1720 }
1721
1722 static I32
1723 sortcv_stacked(pTHX_ SV *a, SV *b)
1724 {
1725     I32 oldsaveix = PL_savestack_ix;
1726     I32 oldscopeix = PL_scopestack_ix;
1727     I32 result;
1728     AV *av;
1729
1730     av = GvAV(PL_defgv);
1731
1732     if (AvMAX(av) < 1) {
1733         SV** ary = AvALLOC(av);
1734         if (AvARRAY(av) != ary) {
1735             AvMAX(av) += AvARRAY(av) - AvALLOC(av);
1736             SvPVX(av) = (char*)ary;
1737         }
1738         if (AvMAX(av) < 1) {
1739             AvMAX(av) = 1;
1740             Renew(ary,2,SV*);
1741             SvPVX(av) = (char*)ary;
1742         }
1743     }
1744     AvFILLp(av) = 1;
1745
1746     AvARRAY(av)[0] = a;
1747     AvARRAY(av)[1] = b;
1748     PL_stack_sp = PL_stack_base;
1749     PL_op = PL_sortcop;
1750     CALLRUNOPS(aTHX);
1751     if (PL_stack_sp != PL_stack_base + 1)
1752         Perl_croak(aTHX_ "Sort subroutine didn't return single value");
1753     if (!SvNIOKp(*PL_stack_sp))
1754         Perl_croak(aTHX_ "Sort subroutine didn't return a numeric value");
1755     result = SvIV(*PL_stack_sp);
1756     while (PL_scopestack_ix > oldscopeix) {
1757         LEAVE;
1758     }
1759     leave_scope(oldsaveix);
1760     return result;
1761 }
1762
1763 static I32
1764 sortcv_xsub(pTHX_ SV *a, SV *b)
1765 {
1766     dSP;
1767     I32 oldsaveix = PL_savestack_ix;
1768     I32 oldscopeix = PL_scopestack_ix;
1769     I32 result;
1770     CV *cv=(CV*)PL_sortcop;
1771
1772     SP = PL_stack_base;
1773     PUSHMARK(SP);
1774     EXTEND(SP, 2);
1775     *++SP = a;
1776     *++SP = b;
1777     PUTBACK;
1778     (void)(*CvXSUB(cv))(aTHX_ cv);
1779     if (PL_stack_sp != PL_stack_base + 1)
1780         Perl_croak(aTHX_ "Sort subroutine didn't return single value");
1781     if (!SvNIOKp(*PL_stack_sp))
1782         Perl_croak(aTHX_ "Sort subroutine didn't return a numeric value");
1783     result = SvIV(*PL_stack_sp);
1784     while (PL_scopestack_ix > oldscopeix) {
1785         LEAVE;
1786     }
1787     leave_scope(oldsaveix);
1788     return result;
1789 }
1790
1791
1792 static I32
1793 sv_ncmp(pTHX_ SV *a, SV *b)
1794 {
1795     NV nv1 = SvNV(a);
1796     NV nv2 = SvNV(b);
1797     return nv1 < nv2 ? -1 : nv1 > nv2 ? 1 : 0;
1798 }
1799
1800 static I32
1801 sv_i_ncmp(pTHX_ SV *a, SV *b)
1802 {
1803     IV iv1 = SvIV(a);
1804     IV iv2 = SvIV(b);
1805     return iv1 < iv2 ? -1 : iv1 > iv2 ? 1 : 0;
1806 }
1807 #define tryCALL_AMAGICbin(left,right,meth,svp) STMT_START { \
1808           *svp = Nullsv;                                \
1809           if (PL_amagic_generation) { \
1810             if (SvAMAGIC(left)||SvAMAGIC(right))\
1811                 *svp = amagic_call(left, \
1812                                    right, \
1813                                    CAT2(meth,_amg), \
1814                                    0); \
1815           } \
1816         } STMT_END
1817
1818 static I32
1819 amagic_ncmp(pTHX_ register SV *a, register SV *b)
1820 {
1821     SV *tmpsv;
1822     tryCALL_AMAGICbin(a,b,ncmp,&tmpsv);
1823     if (tmpsv) {
1824         NV d;
1825  
1826         if (SvIOK(tmpsv)) {
1827             I32 i = SvIVX(tmpsv);
1828             if (i > 0)
1829                return 1;
1830             return i? -1 : 0;
1831         }
1832         d = SvNV(tmpsv);
1833         if (d > 0)
1834            return 1;
1835         return d? -1 : 0;
1836      }
1837      return sv_ncmp(aTHX_ a, b);
1838 }
1839
1840 static I32
1841 amagic_i_ncmp(pTHX_ register SV *a, register SV *b)
1842 {
1843     SV *tmpsv;
1844     tryCALL_AMAGICbin(a,b,ncmp,&tmpsv);
1845     if (tmpsv) {
1846         NV d;
1847
1848         if (SvIOK(tmpsv)) {
1849             I32 i = SvIVX(tmpsv);
1850             if (i > 0)
1851                return 1;
1852             return i? -1 : 0;
1853         }
1854         d = SvNV(tmpsv);
1855         if (d > 0)
1856            return 1;
1857         return d? -1 : 0;
1858     }
1859     return sv_i_ncmp(aTHX_ a, b);
1860 }
1861
1862 static I32
1863 amagic_cmp(pTHX_ register SV *str1, register SV *str2)
1864 {
1865     SV *tmpsv;
1866     tryCALL_AMAGICbin(str1,str2,scmp,&tmpsv);
1867     if (tmpsv) {
1868         NV d;
1869  
1870         if (SvIOK(tmpsv)) {
1871             I32 i = SvIVX(tmpsv);
1872             if (i > 0)
1873                return 1;
1874             return i? -1 : 0;
1875         }
1876         d = SvNV(tmpsv);
1877         if (d > 0)
1878            return 1;
1879         return d? -1 : 0;
1880     }
1881     return sv_cmp(str1, str2);
1882 }
1883
1884 static I32
1885 amagic_cmp_locale(pTHX_ register SV *str1, register SV *str2)
1886 {
1887     SV *tmpsv;
1888     tryCALL_AMAGICbin(str1,str2,scmp,&tmpsv);
1889     if (tmpsv) {
1890         NV d;
1891  
1892         if (SvIOK(tmpsv)) {
1893             I32 i = SvIVX(tmpsv);
1894             if (i > 0)
1895                return 1;
1896             return i? -1 : 0;
1897         }
1898         d = SvNV(tmpsv);
1899         if (d > 0)
1900            return 1;
1901         return d? -1 : 0;
1902     }
1903     return sv_cmp_locale(str1, str2);
1904 }