home *** CD-ROM | disk | FTP | other *** search
/ The CDPD Public Domain Collection for CDTV 3 / CDPDIII.bin / pd / programming / gnuc / stdlib / rcs / radixsort.c,v < prev    next >
Encoding:
Text File  |  1992-07-04  |  9.1 KB  |  315 lines
  1. head    1.1;
  2. access;
  3. symbols
  4.     version39-41:1.1;
  5. locks;
  6. comment    @ * @;
  7.  
  8.  
  9. 1.1
  10. date    92.06.08.17.01.06;    author mwild;    state Exp;
  11. branches;
  12. next    ;
  13.  
  14.  
  15. desc
  16. @initial checkin
  17. @
  18.  
  19.  
  20. 1.1
  21. log
  22. @Initial revision
  23. @
  24. text
  25. @/*-
  26.  * Copyright (c) 1990 The Regents of the University of California.
  27.  * All rights reserved.
  28.  *
  29.  * Redistribution and use in source and binary forms, with or without
  30.  * modification, are permitted provided that the following conditions
  31.  * are met:
  32.  * 1. Redistributions of source code must retain the above copyright
  33.  *    notice, this list of conditions and the following disclaimer.
  34.  * 2. Redistributions in binary form must reproduce the above copyright
  35.  *    notice, this list of conditions and the following disclaimer in the
  36.  *    documentation and/or other materials provided with the distribution.
  37.  * 3. All advertising materials mentioning features or use of this software
  38.  *    must display the following acknowledgement:
  39.  *    This product includes software developed by the University of
  40.  *    California, Berkeley and its contributors.
  41.  * 4. Neither the name of the University nor the names of its contributors
  42.  *    may be used to endorse or promote products derived from this software
  43.  *    without specific prior written permission.
  44.  *
  45.  * THIS SOFTWARE IS PROVIDED BY THE REGENTS AND CONTRIBUTORS ``AS IS'' AND
  46.  * ANY EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT LIMITED TO, THE
  47.  * IMPLIED WARRANTIES OF MERCHANTABILITY AND FITNESS FOR A PARTICULAR PURPOSE
  48.  * ARE DISCLAIMED.  IN NO EVENT SHALL THE REGENTS OR CONTRIBUTORS BE LIABLE
  49.  * FOR ANY DIRECT, INDIRECT, INCIDENTAL, SPECIAL, EXEMPLARY, OR CONSEQUENTIAL
  50.  * DAMAGES (INCLUDING, BUT NOT LIMITED TO, PROCUREMENT OF SUBSTITUTE GOODS
  51.  * OR SERVICES; LOSS OF USE, DATA, OR PROFITS; OR BUSINESS INTERRUPTION)
  52.  * HOWEVER CAUSED AND ON ANY THEORY OF LIABILITY, WHETHER IN CONTRACT, STRICT
  53.  * LIABILITY, OR TORT (INCLUDING NEGLIGENCE OR OTHERWISE) ARISING IN ANY WAY
  54.  * OUT OF THE USE OF THIS SOFTWARE, EVEN IF ADVISED OF THE POSSIBILITY OF
  55.  * SUCH DAMAGE.
  56.  */
  57.  
  58. #if defined(LIBC_SCCS) && !defined(lint)
  59. static char sccsid[] = "@@(#)radixsort.c    5.7 (Berkeley) 2/23/91";
  60. #endif /* LIBC_SCCS and not lint */
  61.  
  62. #include <sys/types.h>
  63. #include <limits.h>
  64. #include <stdlib.h>
  65. #include <stddef.h>
  66. #include <string.h>
  67.  
  68. /*
  69.  * __rspartition is the cutoff point for a further partitioning instead
  70.  * of a shellsort.  If it changes check __rsshell_increments.  Both of
  71.  * these are exported, as the best values are data dependent.
  72.  */
  73. #define    NPARTITION    40
  74. int __rspartition = NPARTITION;
  75. int __rsshell_increments[] = { 4, 1, 0, 0, 0, 0, 0, 0 };
  76.  
  77. /*
  78.  * Stackp points to context structures, where each structure schedules a
  79.  * partitioning.  Radixsort exits when the stack is empty.
  80.  *
  81.  * If the buckets are placed on the stack randomly, the worst case is when
  82.  * all the buckets but one contain (npartitions + 1) elements and the bucket
  83.  * pushed on the stack last contains the rest of the elements.  In this case,
  84.  * stack growth is bounded by:
  85.  *
  86.  *    limit = (nelements / (npartitions + 1)) - 1;
  87.  *
  88.  * This is a very large number, 52,377,648 for the maximum 32-bit signed int.
  89.  *
  90.  * By forcing the largest bucket to be pushed on the stack first, the worst
  91.  * case is when all but two buckets each contain (npartitions + 1) elements,
  92.  * with the remaining elements split equally between the first and last
  93.  * buckets pushed on the stack.  In this case, stack growth is bounded when:
  94.  *
  95.  *    for (partition_cnt = 0; nelements > npartitions; ++partition_cnt)
  96.  *        nelements =
  97.  *            (nelements - (npartitions + 1) * (nbuckets - 2)) / 2;
  98.  * The bound is:
  99.  *
  100.  *    limit = partition_cnt * (nbuckets - 1);
  101.  *
  102.  * This is a much smaller number, 4590 for the maximum 32-bit signed int.
  103.  */
  104. #define    NBUCKETS    (UCHAR_MAX + 1)
  105.  
  106. typedef struct _stack {
  107.     const u_char **bot;
  108.     int indx, nmemb;
  109. } CONTEXT;
  110.  
  111. #define    STACKPUSH { \
  112.     stackp->bot = p; \
  113.     stackp->nmemb = nmemb; \
  114.     stackp->indx = indx; \
  115.     ++stackp; \
  116. }
  117. #define    STACKPOP { \
  118.     if (stackp == stack) \
  119.         break; \
  120.     --stackp; \
  121.     bot = stackp->bot; \
  122.     nmemb = stackp->nmemb; \
  123.     indx = stackp->indx; \
  124. }
  125.  
  126. /*
  127.  * A variant of MSD radix sorting; see Knuth Vol. 3, page 177, and 5.2.5,
  128.  * Ex. 10 and 12.  Also, "Three Partition Refinement Algorithms, Paige
  129.  * and Tarjan, SIAM J. Comput. Vol. 16, No. 6, December 1987.
  130.  *
  131.  * This uses a simple sort as soon as a bucket crosses a cutoff point,
  132.  * rather than sorting the entire list after partitioning is finished.
  133.  * This should be an advantage.
  134.  *
  135.  * This is pure MSD instead of LSD of some number of MSD, switching to
  136.  * the simple sort as soon as possible.  Takes linear time relative to
  137.  * the number of bytes in the strings.
  138.  */
  139. int
  140. #if __STDC__
  141. radixsort(const u_char **l1, int nmemb, const u_char *tab, u_char endbyte)
  142. #else
  143. radixsort(l1, nmemb, tab, endbyte)
  144.     const u_char **l1;
  145.     register int nmemb;
  146.     const u_char *tab;
  147.     u_char endbyte;
  148. #endif
  149. {
  150.     register int i, indx, t1, t2;
  151.     register const u_char **l2;
  152.     register const u_char **p;
  153.     register const u_char **bot;
  154.     register const u_char *tr;
  155.     CONTEXT *stack, *stackp;
  156.     int c[NBUCKETS + 1], max;
  157.     u_char ltab[NBUCKETS];
  158.     static void shellsort();
  159.  
  160.     if (nmemb <= 1)
  161.         return(0);
  162.  
  163.     /*
  164.      * T1 is the constant part of the equation, the number of elements
  165.      * represented on the stack between the top and bottom entries.
  166.      * It doesn't get rounded as the divide by 2 rounds down (correct
  167.      * for a value being subtracted).  T2, the nelem value, has to be
  168.      * rounded up before each divide because we want an upper bound;
  169.      * this could overflow if nmemb is the maximum int.
  170.      */
  171.     t1 = ((__rspartition + 1) * (NBUCKETS - 2)) >> 1;
  172.     for (i = 0, t2 = nmemb; t2 > __rspartition; i += NBUCKETS - 1)
  173.         t2 = ((t2 + 1) >> 1) - t1;
  174.     if (i) {
  175.         if (!(stack = stackp = (CONTEXT *)malloc(i * sizeof(CONTEXT))))
  176.             return(-1);
  177.     } else
  178.         stack = stackp = NULL;
  179.  
  180.     /*
  181.      * There are two arrays, one provided by the user (l1), and the
  182.      * temporary one (l2).  The data is sorted to the temporary stack,
  183.      * and then copied back.  The speedup of using index to determine
  184.      * which stack the data is on and simply swapping stacks back and
  185.      * forth, thus avoiding the copy every iteration, turns out to not
  186.      * be any faster than the current implementation.
  187.      */
  188.     if (!(l2 = (const u_char **)malloc(sizeof(u_char *) * nmemb)))
  189.         return(-1);
  190.  
  191.     /*
  192.      * Tr references a table of sort weights; multiple entries may
  193.      * map to the same weight; EOS char must have the lowest weight.
  194.      */
  195.     if (tab)
  196.         tr = tab;
  197.     else {
  198.         for (t1 = 0, t2 = endbyte; t1 < t2; ++t1)
  199.             ltab[t1] = t1 + 1;
  200.         ltab[t2] = 0;
  201.         for (t1 = endbyte + 1; t1 < NBUCKETS; ++t1)
  202.             ltab[t1] = t1;
  203.         tr = ltab;
  204.     }
  205.  
  206.     /* First sort is entire stack */
  207.     bot = l1;
  208.     indx = 0;
  209.  
  210.     for (;;) {
  211.         /* Clear bucket count array */
  212.         bzero((char *)c, sizeof(c));
  213.  
  214.         /*
  215.          * Compute number of items that sort to the same bucket
  216.          * for this index.
  217.          */
  218.         for (p = bot, i = nmemb; --i >= 0;)
  219.             ++c[tr[(*p++)[indx]]];
  220.  
  221.         /*
  222.          * Sum the number of characters into c, dividing the temp
  223.          * stack into the right number of buckets for this bucket,
  224.          * this index.  C contains the cumulative total of keys
  225.          * before and included in this bucket, and will later be
  226.          * used as an index to the bucket.  c[NBUCKETS] contains
  227.          * the total number of elements, for determining how many
  228.          * elements the last bucket contains.  At the same time
  229.          * find the largest bucket so it gets pushed first.
  230.          */
  231.         for (i = max = t1 = 0, t2 = __rspartition; i <= NBUCKETS; ++i) {
  232.             if (c[i] > t2) {
  233.                 t2 = c[i];
  234.                 max = i;
  235.             }
  236.             t1 = c[i] += t1;
  237.         }
  238.  
  239.         /*
  240.          * Partition the elements into buckets; c decrements through
  241.          * the bucket, and ends up pointing to the first element of
  242.          * the bucket.
  243.          */
  244.         for (i = nmemb; --i >= 0;) {
  245.             --p;
  246.             l2[--c[tr[(*p)[indx]]]] = *p;
  247.         }
  248.  
  249.         /* Copy the partitioned elements back to user stack */
  250.         bcopy(l2, bot, nmemb * sizeof(u_char *));
  251.  
  252.         ++indx;
  253.         /*
  254.          * Sort buckets as necessary; don't sort c[0], it's the
  255.          * EOS character bucket, and nothing can follow EOS.
  256.          */
  257.         for (i = max; i; --i) {
  258.             if ((nmemb = c[i + 1] - (t1 = c[i])) < 2)
  259.                 continue;
  260.             p = bot + t1;
  261.             if (nmemb > __rspartition)
  262.                 STACKPUSH
  263.             else
  264.                 shellsort(p, indx, nmemb, tr);
  265.         }
  266.         for (i = max + 1; i < NBUCKETS; ++i) {
  267.             if ((nmemb = c[i + 1] - (t1 = c[i])) < 2)
  268.                 continue;
  269.             p = bot + t1;
  270.             if (nmemb > __rspartition)
  271.                 STACKPUSH
  272.             else
  273.                 shellsort(p, indx, nmemb, tr);
  274.         }
  275.         /* Break out when stack is empty */
  276.         STACKPOP
  277.     }
  278.  
  279.     free((char *)l2);
  280.     free((char *)stack);
  281.     return(0);
  282. }
  283.  
  284. /*
  285.  * Shellsort (diminishing increment sort) from Data Structures and
  286.  * Algorithms, Aho, Hopcraft and Ullman, 1983 Edition, page 290;
  287.  * see also Knuth Vol. 3, page 84.  The increments are selected from
  288.  * formula (8), page 95.  Roughly O(N^3/2).
  289.  */
  290. static void
  291. shellsort(p, indx, nmemb, tr)
  292.     register u_char **p, *tr;
  293.     register int indx, nmemb;
  294. {
  295.     register u_char ch, *s1, *s2;
  296.     register int incr, *incrp, t1, t2;
  297.  
  298.     for (incrp = __rsshell_increments; incr = *incrp++;)
  299.         for (t1 = incr; t1 < nmemb; ++t1)
  300.             for (t2 = t1 - incr; t2 >= 0;) {
  301.                 s1 = p[t2] + indx;
  302.                 s2 = p[t2 + incr] + indx;
  303.                 while ((ch = tr[*s1++]) == tr[*s2] && ch)
  304.                     ++s2;
  305.                 if (ch > tr[*s2]) {
  306.                     s1 = p[t2];
  307.                     p[t2] = p[t2 + incr];
  308.                     p[t2 + incr] = s1;
  309.                     t2 -= incr;
  310.                 } else
  311.                     break;
  312.             }
  313. }
  314. @
  315.