Removed extra free

[pspp-builds.git] / src / algorithm.c

diff --git a/src/algorithm.c b/src/algorithm.c
index 8a569749e6ddfd85fc89d4a7a24367bd02e67b37..cfb1ba9fbf533ef55baefbeb06c594eaf9220611 100644 (file)
--- a/src/algorithm.c
+++ b/src/algorithm.c
@@ -14,8 +14,8 @@
 
    You should have received a copy of the GNU General Public License
    along with this program; if not, write to the Free Software
 
    You should have received a copy of the GNU General Public License
    along with this program; if not, write to the Free Software

-   Foundation, Inc., 59 Temple Place - Suite 330, Boston, MA
-   02111-1307, USA. */
+   Foundation, Inc., 51 Franklin Street, Fifth Floor, Boston, MA
+   02110-1301, USA. */

 
 /* Copyright (C) 2001 Free Software Foundation, Inc.
   
 
 /* Copyright (C) 2001 Free Software Foundation, Inc.
   


@@ -32,7 +32,7 @@
 
    You should have received a copy of the GNU General Public License along
    with this library; see the file COPYING.  If not, write to the Free
 
    You should have received a copy of the GNU General Public License along
    with this library; see the file COPYING.  If not, write to the Free

-   Software Foundation, 59 Temple Place - Suite 330, Boston, MA 02111-1307,
+   Software Foundation, 51 Franklin Street, Fifth Floor, Boston, MA 02110-1301,

    USA.
 
    As a special exception, you may use this file as part of a free software
    USA.
 
    As a special exception, you may use this file as part of a free software


@@ -86,28 +86,37 @@
 
    You should have received a copy of the GNU Lesser General Public
    License along with the GNU C Library; if not, write to the Free
 
    You should have received a copy of the GNU Lesser General Public
    License along with the GNU C Library; if not, write to the Free

-   Software Foundation, Inc., 59 Temple Place, Suite 330, Boston, MA
-   02111-1307 USA.  */
+   Software Foundation, Inc., 51 Franklin Street, Fifth Floor, Boston, MA
+   02110-1301 USA.  */

 
 #include <config.h>
 #include "algorithm.h"
 
 #include <config.h>
 #include "algorithm.h"

-#include <assert.h>
+#include <gsl/gsl_rng.h>

 #include <limits.h>
 #include <stdlib.h>
 #include <string.h>
 #include "alloc.h"
 #include <limits.h>
 #include <stdlib.h>
 #include <string.h>
 #include "alloc.h"

-#include "random.h"
+
+/* Some of the assertions in this file are very expensive.  We
+   don't use them by default. */
+#ifdef EXTRA_CHECKS
+#define expensive_assert(X) assert(X)
+#else
+#define expensive_assert(X) ((void) 0)
+#endif
+#include "error.h"

 \f
 /* Finds an element in ARRAY, which contains COUNT elements of
    SIZE bytes each, using COMPARE for comparisons.  Returns the
    first element in ARRAY that matches TARGET, or a null pointer
    on failure.  AUX is passed to each comparison as auxiliary
    data. */
 \f
 /* Finds an element in ARRAY, which contains COUNT elements of
    SIZE bytes each, using COMPARE for comparisons.  Returns the
    first element in ARRAY that matches TARGET, or a null pointer
    on failure.  AUX is passed to each comparison as auxiliary
    data. */

-void *find (const void *array, size_t count, size_t size,
-            const void *target,
-            algo_compare_func *compare, void *aux) 
+void *
+find (const void *array, size_t count, size_t size,
+      const void *target,
+      algo_compare_func *compare, void *aux) 

 {
 {

-  const unsigned char *element = array;
+  const char *element = array;

 
   while (count-- > 0) 
     {
 
   while (count-- > 0) 
     {


@@ -119,6 +128,51 @@ void *find (const void *array, size_t count, size_t size,
 
   return NULL;
 }
 
   return NULL;
 }

+
+/* Counts and return the number of elements in ARRAY, which
+   contains COUNT elements of SIZE bytes each, which are equal to
+   ELEMENT as compared with COMPARE.  AUX is passed as auxiliary
+   data to COMPARE. */
+size_t
+count_equal (const void *array, size_t count, size_t size,
+             const void *element,
+             algo_compare_func *compare, void *aux)
+{
+  const char *first = array;
+  size_t equal_cnt = 0;
+
+  while (count-- > 0) 
+    {
+      if (compare (element, first, aux) == 0)
+        equal_cnt++;
+      
+      first += size;
+    }
+
+  return equal_cnt;
+}
+
+/* Counts and return the number of elements in ARRAY, which
+   contains COUNT elements of SIZE bytes each, for which
+   PREDICATE returns nonzero.  AUX is passed as auxiliary data to
+   PREDICATE. */
+size_t
+count_if (const void *array, size_t count, size_t size,
+          algo_predicate_func *predicate, void *aux) 
+{
+  const char *first = array;
+  size_t nonzero_cnt = 0;
+
+  while (count-- > 0) 
+    {
+      if (predicate (first, aux) != 0)
+        nonzero_cnt++;
+      
+      first += size;
+    }
+
+  return nonzero_cnt;
+}

 \f
 /* Byte-wise swap two items of size SIZE. */
 #define SWAP(a, b, size)                        \
 \f
 /* Byte-wise swap two items of size SIZE. */
 #define SWAP(a, b, size)                        \


@@ -148,8 +202,12 @@ unique (void *array, size_t count, size_t size,
   for (;;) 
     {
       first += size;
   for (;;) 
     {
       first += size;

-      if (first >= last)
-        return count;
+      if (first >= last) 
+        {
+          assert (adjacent_find_equal (array, count,
+                                       size, compare, aux) == NULL);
+          return count; 
+        }

 
       if (compare (result, first, aux)) 
         {
 
       if (compare (result, first, aux)) 
         {


@@ -181,8 +239,9 @@ size_t
 partition (void *array, size_t count, size_t size,
            algo_predicate_func *predicate, void *aux) 
 {
 partition (void *array, size_t count, size_t size,
            algo_predicate_func *predicate, void *aux) 
 {

+  size_t nonzero_cnt = count;

   char *first = array;
   char *first = array;

-  char *last = first + count * size;
+  char *last = first + nonzero_cnt * size;

 
   for (;;)
     {
 
   for (;;)
     {


@@ -191,13 +250,13 @@ partition (void *array, size_t count, size_t size,
       for (;;) 
         {
           if (first == last)
       for (;;) 
         {
           if (first == last)

-            return count;
+            goto done;

           else if (!predicate (first, aux)) 
             break;
 
           first += size; 
         }
           else if (!predicate (first, aux)) 
             break;
 
           first += size; 
         }

-      count--;
+      nonzero_cnt--;

 
       /* Move LAST backward to point to last element that passes
          PREDICATE. */
 
       /* Move LAST backward to point to last element that passes
          PREDICATE. */


@@ -206,11 +265,11 @@ partition (void *array, size_t count, size_t size,
           last -= size;
 
           if (first == last)
           last -= size;
 
           if (first == last)

-            return count;
+            goto done;

           else if (predicate (last, aux)) 
             break;
           else
           else if (predicate (last, aux)) 
             break;
           else

-            count--;
+            nonzero_cnt--;

         }
       
       /* By swapping FIRST and LAST we extend the starting and
         }
       
       /* By swapping FIRST and LAST we extend the starting and


@@ -219,31 +278,32 @@ partition (void *array, size_t count, size_t size,
       SWAP (first, last, size);
       first += size;
     }
       SWAP (first, last, size);
       first += size;
     }

-}
-\f
-/* A algo_random_func that uses random.h. */
-unsigned
-algo_default_random (unsigned max, void *aux unused) 
-{
-  return rng_get_unsigned (pspp_rng ()) % max;
+
+ done:
+  assert (is_partitioned (array, count, size, nonzero_cnt, predicate, aux));
+  return nonzero_cnt; 

 }
 
 }
 

-/* Randomly reorders ARRAY, which contains COUNT elements of SIZE
-   bytes each.  Uses RANDOM as a source of random data, passing
-   AUX as the auxiliary data.  RANDOM may be null to use a
-   default random source. */
-void
-random_shuffle (void *array_, size_t count, size_t size,
-                algo_random_func *random, void *aux)
+/* Checks whether ARRAY, which contains COUNT elements of SIZE
+   bytes each, is partitioned such that PREDICATE returns nonzero
+   for the first NONZERO_CNT elements and zero for the remaining
+   elements.  AUX is passed as auxiliary data to PREDICATE. */
+int
+is_partitioned (const void *array, size_t count, size_t size,
+                size_t nonzero_cnt,
+                algo_predicate_func *predicate, void *aux) 

 {
 {

-  unsigned char *array = array_;
-  int i;
-
-  if (random == NULL)
-    random = algo_default_random;
-
-  for (i = 1; i < count; i++)
-    SWAP (array + i * size, array + random (i + 1, aux) * size, size);
+  const char *first = array;
+  size_t idx;
+
+  assert (nonzero_cnt <= count);
+  for (idx = 0; idx < nonzero_cnt; idx++)
+    if (predicate (first + idx * size, aux) == 0)
+      return 0;
+  for (idx = nonzero_cnt; idx < count; idx++)
+    if (predicate (first + idx * size, aux) != 0)
+      return 0;
+  return 1;

 }
 \f
 /* Copies the COUNT elements of SIZE bytes each from ARRAY to
 }
 \f
 /* Copies the COUNT elements of SIZE bytes each from ARRAY to


@@ -255,9 +315,10 @@ copy_if (const void *array, size_t count, size_t size,
          void *result,
          algo_predicate_func *predicate, void *aux) 
 {
          void *result,
          algo_predicate_func *predicate, void *aux) 
 {

-  const unsigned char *input = array;
-  const unsigned char *last = input + size * count;
-  unsigned char *output = result;
+  const char *input = array;
+  const char *last = input + size * count;
+  char *output = result;
+  size_t nonzero_cnt = 0;

   
   while (input < last)
     {
   
   while (input < last)
     {


@@ -265,14 +326,74 @@ copy_if (const void *array, size_t count, size_t size,
         {
           memcpy (output, input, size);
           output += size;
         {
           memcpy (output, input, size);
           output += size;

+          nonzero_cnt++;

         }
         }

-      else
-        count--;

 
       input += size;
     }
 
 
       input += size;
     }
 

-  return count;
+  assert (nonzero_cnt == count_if (array, count, size, predicate, aux));
+  assert (nonzero_cnt == count_if (result, nonzero_cnt, size, predicate, aux));
+
+  return nonzero_cnt;
+}
+
+/* Removes N elements starting at IDX from ARRAY, which consists
+   of COUNT elements of SIZE bytes each, by shifting the elements
+   following them, if any, into its position. */
+void
+remove_range (void *array_, size_t count, size_t size,
+              size_t idx, size_t n) 
+{
+  char *array = array_;
+  
+  assert (array != NULL);
+  assert (idx <= count);
+  assert (idx + n <= count);
+
+  if (idx + n < count)
+    memmove (array + idx * size, array + (idx + n) * size,
+             size * (count - idx - n));
+}
+
+/* Removes element IDX from ARRAY, which consists of COUNT
+   elements of SIZE bytes each, by shifting the elements
+   following it, if any, into its position. */
+void
+remove_element (void *array, size_t count, size_t size,
+                size_t idx) 
+{
+  remove_range (array, count, size, idx, 1);
+}
+
+/* Moves an element in ARRAY, which consists of COUNT elements of
+   SIZE bytes each, from OLD_IDX to NEW_IDX, shifting around
+   other elements as needed.  Runs in O(abs(OLD_IDX - NEW_IDX))
+   time. */
+void
+move_element (void *array_, size_t count, size_t size,
+              size_t old_idx, size_t new_idx) 
+{
+  assert (array_ != NULL || count == 0);
+  assert (old_idx < count);
+  assert (new_idx < count);
+  
+  if (old_idx != new_idx) 
+    {
+      char *array = array_;
+      char *element = xmalloc (size);
+      char *new = array + new_idx * size;
+      char *old = array + old_idx * size;
+
+      memcpy (element, old, size);
+      if (new < old)
+        memmove (new + size, new, (old_idx - new_idx) * size);
+      else
+        memmove (old, old + size, (new_idx - old_idx) * size);
+      memcpy (new, element, size);
+
+      free (element);
+    }

 }
 
 /* A predicate and its auxiliary data. */
 }
 
 /* A predicate and its auxiliary data. */


@@ -299,14 +420,14 @@ remove_equal (void *array, size_t count, size_t size,
               void *element,
               algo_compare_func *compare, void *aux) 
 {
               void *element,
               algo_compare_func *compare, void *aux) 
 {

-  unsigned char *first = array;
-  unsigned char *last = first + count * size;
-  unsigned char *result;
+  char *first = array;
+  char *last = first + count * size;
+  char *result;

 
   for (;;)
     {
       if (first >= last)
 
   for (;;)
     {
       if (first >= last)

-        return count;
+        goto done;

       if (compare (first, element, aux) == 0)
         break;
 
       if (compare (first, element, aux) == 0)
         break;
 


@@ -319,7 +440,7 @@ remove_equal (void *array, size_t count, size_t size,
     {
       first += size;
       if (first >= last)
     {
       first += size;
       if (first >= last)

-        return count;
+        goto done;

 
       if (compare (first, element, aux) == 0) 
         {
 
       if (compare (first, element, aux) == 0) 
         {


@@ -331,6 +452,8 @@ remove_equal (void *array, size_t count, size_t size,
       result += size;
     }
 
       result += size;
     }
 

+ done:
+  assert (count_equal (array, count, size, element, compare, aux) == 0);

   return count;
 }
 
   return count;
 }
 


@@ -365,14 +488,14 @@ binary_search (const void *array, size_t count, size_t size,
 
   if (count != 0) 
     {
 
   if (count != 0) 
     {

-      const unsigned char *first = array;
+      const char *first = array;

       int low = 0;
       int high = count - 1;
 
       while (low <= high) 
         {
           int middle = (low + high) / 2;
       int low = 0;
       int high = count - 1;
 
       while (low <= high) 
         {
           int middle = (low + high) / 2;

-          const unsigned char *element = first + middle * size;
+          const char *element = first + middle * size;

           int cmp = compare (value, element, aux);
 
           if (cmp > 0) 
           int cmp = compare (value, element, aux);
 
           if (cmp > 0) 


@@ -383,6 +506,8 @@ binary_search (const void *array, size_t count, size_t size,
             return (void *) element;
         }
     }
             return (void *) element;
         }
     }

+
+  expensive_assert (find (array, count, size, value, compare, aux) == NULL);

   return NULL;
 }
 \f
   return NULL;
 }
 \f


@@ -392,13 +517,13 @@ binary_search (const void *array, size_t count, size_t size,
    strcmp()-type result.  AUX is passed to COMPARE as auxiliary
    data. */
 int
    strcmp()-type result.  AUX is passed to COMPARE as auxiliary
    data. */
 int

-lexicographical_compare (const void *array1, size_t count1,
-                         const void *array2, size_t count2,
-                         size_t size,
-                         algo_compare_func *compare, void *aux) 
+lexicographical_compare_3way (const void *array1, size_t count1,
+                              const void *array2, size_t count2,
+                              size_t size,
+                              algo_compare_func *compare, void *aux) 

 {
 {

-  const unsigned char *first1 = array1;
-  const unsigned char *first2 = array2;
+  const char *first1 = array1;
+  const char *first2 = array2;

   size_t min_count = count1 < count2 ? count1 : count2;
 
   while (min_count > 0)
   size_t min_count = count1 < count2 ? count1 : count2;
 
   while (min_count > 0)


@@ -419,7 +544,6 @@ lexicographical_compare (const void *array1, size_t count1,
    Engineering a sort function; Jon Bentley and M. Douglas McIlroy;
    Software - Practice and Experience; Vol. 23 (11), 1249-1265, 1993.  */
 
    Engineering a sort function; Jon Bentley and M. Douglas McIlroy;
    Software - Practice and Experience; Vol. 23 (11), 1249-1265, 1993.  */
 

-#include <alloca.h>

 #include <limits.h>
 #include <stdlib.h>
 #include <string.h>
 #include <limits.h>
 #include <stdlib.h>
 #include <string.h>


@@ -471,21 +595,20 @@ typedef struct
       stack size is needed (actually O(1) in this case)!  */
 
 void
       stack size is needed (actually O(1) in this case)!  */
 
 void

-sort (void *const pbase, size_t total_elems, size_t size,
-      algo_compare_func *cmp, void *aux)
+sort (void *array, size_t count, size_t size,
+      algo_compare_func *compare, void *aux)

 {
 {

-  register char *base_ptr = (char *) pbase;
-
+  char *const first = array;

   const size_t max_thresh = MAX_THRESH * size;
 
   const size_t max_thresh = MAX_THRESH * size;
 

-  if (total_elems == 0)
+  if (count == 0)

     /* Avoid lossage with unsigned arithmetic below.  */
     return;
 
     /* Avoid lossage with unsigned arithmetic below.  */
     return;
 

-  if (total_elems > MAX_THRESH)
+  if (count > MAX_THRESH)

     {
     {

-      char *lo = base_ptr;
-      char *hi = &lo[size * (total_elems - 1)];
+      char *lo = first;
+      char *hi = &lo[size * (count - 1)];

       stack_node stack[STACK_SIZE];
       stack_node *top = stack + 1;
 
       stack_node stack[STACK_SIZE];
       stack_node *top = stack + 1;
 


@@ -502,13 +625,13 @@ sort (void *const pbase, size_t total_elems, size_t size,
 
          char *mid = lo + size * ((hi - lo) / size >> 1);
 
 
          char *mid = lo + size * ((hi - lo) / size >> 1);
 

-         if ((*cmp) ((void *) mid, (void *) lo, aux) < 0)
+         if (compare (mid, lo, aux) < 0)

            SWAP (mid, lo, size);
            SWAP (mid, lo, size);

-         if ((*cmp) ((void *) hi, (void *) mid, aux) < 0)
+         if (compare (hi, mid, aux) < 0)

            SWAP (mid, hi, size);
          else
            goto jump_over;
            SWAP (mid, hi, size);
          else
            goto jump_over;

-         if ((*cmp) ((void *) mid, (void *) lo, aux) < 0)
+         if (compare (mid, lo, aux) < 0)

            SWAP (mid, lo, size);
        jump_over:;
 
            SWAP (mid, lo, size);
        jump_over:;
 


@@ -520,10 +643,10 @@ sort (void *const pbase, size_t total_elems, size_t size,
             that this algorithm runs much faster than others. */
          do
            {
             that this algorithm runs much faster than others. */
          do
            {

-             while ((*cmp) ((void *) left_ptr, (void *) mid, aux) < 0)
+             while (compare (left_ptr, mid, aux) < 0)

                left_ptr += size;
 
                left_ptr += size;
 

-             while ((*cmp) ((void *) mid, (void *) right_ptr, aux) < 0)
+             while (compare (mid, right_ptr, aux) < 0)

                right_ptr -= size;
 
              if (left_ptr < right_ptr)
                right_ptr -= size;
 
              if (left_ptr < right_ptr)


@@ -577,18 +700,18 @@ sort (void *const pbase, size_t total_elems, size_t size,
         }
     }
 
         }
     }
 

-  /* Once the BASE_PTR array is partially sorted by quicksort the rest
+  /* Once the FIRST array is partially sorted by quicksort the rest

      is completely sorted using insertion sort, since this is efficient
      is completely sorted using insertion sort, since this is efficient

-     for partitions below MAX_THRESH size. BASE_PTR points to the beginning
+     for partitions below MAX_THRESH size. FIRST points to the beginning

      of the array to sort, and END_PTR points at the very last element in
      the array (*not* one beyond it!). */
 
 #define min(x, y) ((x) < (y) ? (x) : (y))
 
   {
      of the array to sort, and END_PTR points at the very last element in
      the array (*not* one beyond it!). */
 
 #define min(x, y) ((x) < (y) ? (x) : (y))
 
   {

-    char *const end_ptr = &base_ptr[size * (total_elems - 1)];
-    char *tmp_ptr = base_ptr;
-    char *thresh = min(end_ptr, base_ptr + max_thresh);
+    char *const end_ptr = &first[size * (count - 1)];
+    char *tmp_ptr = first;
+    char *thresh = min(end_ptr, first + max_thresh);

     register char *run_ptr;
 
     /* Find smallest element in first threshold and place it at the
     register char *run_ptr;
 
     /* Find smallest element in first threshold and place it at the


@@ -596,19 +719,19 @@ sort (void *const pbase, size_t total_elems, size_t size,
        and the operation speeds up insertion sort's inner loop. */
 
     for (run_ptr = tmp_ptr + size; run_ptr <= thresh; run_ptr += size)
        and the operation speeds up insertion sort's inner loop. */
 
     for (run_ptr = tmp_ptr + size; run_ptr <= thresh; run_ptr += size)

-      if ((*cmp) ((void *) run_ptr, (void *) tmp_ptr, aux) < 0)
+      if (compare (run_ptr, tmp_ptr, aux) < 0)

         tmp_ptr = run_ptr;
 
         tmp_ptr = run_ptr;
 

-    if (tmp_ptr != base_ptr)
-      SWAP (tmp_ptr, base_ptr, size);
+    if (tmp_ptr != first)
+      SWAP (tmp_ptr, first, size);

 
     /* Insertion sort, running from left-hand-side up to right-hand-side.  */
 
 
     /* Insertion sort, running from left-hand-side up to right-hand-side.  */
 

-    run_ptr = base_ptr + size;
+    run_ptr = first + size;

     while ((run_ptr += size) <= end_ptr)
       {
        tmp_ptr = run_ptr - size;
     while ((run_ptr += size) <= end_ptr)
       {
        tmp_ptr = run_ptr - size;

-       while ((*cmp) ((void *) run_ptr, (void *) tmp_ptr, aux) < 0)
+       while (compare (run_ptr, tmp_ptr, aux) < 0)

          tmp_ptr -= size;
 
        tmp_ptr += size;
          tmp_ptr -= size;
 
        tmp_ptr += size;


@@ -629,6 +752,25 @@ sort (void *const pbase, size_t total_elems, size_t size,
           }
       }
   }
           }
       }
   }

+
+  assert (is_sorted (array, count, size, compare, aux));
+}
+
+/* Tests whether ARRAY, which contains COUNT elements of SIZE
+   bytes each, is sorted in order according to COMPARE.  AUX is
+   passed to COMPARE as auxiliary data. */
+int
+is_sorted (const void *array, size_t count, size_t size,
+           algo_compare_func *compare, void *aux) 
+{
+  const char *first = array;
+  size_t idx;
+      
+  for (idx = 0; idx + 1 < count; idx++)
+    if (compare (first + idx * size, first + (idx + 1) * size, aux) > 0)
+      return 0; 
+  
+  return 1;

 }
 \f
 /* Computes the generalized set difference, ARRAY1 minus ARRAY2,
 }
 \f
 /* Computes the generalized set difference, ARRAY1 minus ARRAY2,


@@ -645,11 +787,11 @@ size_t set_difference (const void *array1, size_t count1,
                        void *result_,
                        algo_compare_func *compare, void *aux) 
 {
                        void *result_,
                        algo_compare_func *compare, void *aux) 
 {

-  const unsigned char *first1 = array1;
-  const unsigned char *last1 = first1 + count1 * size;
-  const unsigned char *first2 = array2;
-  const unsigned char *last2 = first2 + count2 * size;
-  unsigned char *result = result_;
+  const char *first1 = array1;
+  const char *last1 = first1 + count1 * size;
+  const char *first2 = array2;
+  const char *last2 = first2 + count2 * size;
+  char *result = result_;

   size_t result_count = 0;
   
   while (first1 != last1 && first2 != last2) 
   size_t result_count = 0;
   
   while (first1 != last1 && first2 != last2) 


@@ -692,8 +834,8 @@ void *
 adjacent_find_equal (const void *array, size_t count, size_t size,
                      algo_compare_func *compare, void *aux) 
 {
 adjacent_find_equal (const void *array, size_t count, size_t size,
                      algo_compare_func *compare, void *aux) 
 {

-  const unsigned char *first = array;
-  const unsigned char *last = first + count * size;
+  const char *first = array;
+  const char *last = first + count * size;

 
   while (first < last && first + size < last) 
     {
 
   while (first < last && first + size < last) 
     {


@@ -704,4 +846,142 @@ adjacent_find_equal (const void *array, size_t count, size_t size,
 
   return NULL;
 }
 
   return NULL;
 }

+\f
+/* ARRAY contains COUNT elements of SIZE bytes each.  Initially
+   the first COUNT - 1 elements of these form a heap, followed by
+   a single element not part of the heap.  This function adds the
+   final element, forming a heap of COUNT elements in ARRAY.
+   Uses COMPARE to compare elements, passing AUX as auxiliary
+   data. */
+void
+push_heap (void *array, size_t count, size_t size,
+           algo_compare_func *compare, void *aux) 
+{
+  char *first = array;
+  size_t i;
+  
+  expensive_assert (count < 1 || is_heap (array, count - 1,
+                                          size, compare, aux));
+  for (i = count; i > 1; i /= 2) 
+    {
+      char *parent = first + (i / 2 - 1) * size;
+      char *element = first + (i - 1) * size;
+      if (compare (parent, element, aux) < 0)
+        SWAP (parent, element, size);
+      else
+        break; 
+    }
+  expensive_assert (is_heap (array, count, size, compare, aux));
+}
+
+/* ARRAY contains COUNT elements of SIZE bytes each.  Initially
+   the children of ARRAY[idx - 1] are heaps, but ARRAY[idx - 1]
+   may be smaller than its children.  This function fixes that,
+   so that ARRAY[idx - 1] itself is a heap.  Uses COMPARE to
+   compare elements, passing AUX as auxiliary data. */
+static void
+heapify (void *array, size_t count, size_t size,
+         size_t idx,
+         algo_compare_func *compare, void *aux) 
+{
+  char *first = array;
+  
+  for (;;) 
+    {
+      size_t left = 2 * idx;
+      size_t right = left + 1;
+      size_t largest = idx;
+
+      if (left <= count
+          && compare (first + size * (left - 1),
+                      first + size * (idx - 1), aux) > 0)
+        largest = left;
+
+      if (right <= count
+          && compare (first + size * (right - 1),
+                      first + size * (largest - 1), aux) > 0)
+        largest = right;
+
+      if (largest == idx)
+        break;
+
+      SWAP (first + size * (idx - 1), first + size * (largest - 1), size);
+      idx = largest;
+    }
+}
+
+/* ARRAY contains COUNT elements of SIZE bytes each.  Initially
+   all COUNT elements form a heap.  This function moves the
+   largest element in the heap to the final position in ARRAY and
+   reforms a heap of the remaining COUNT - 1 elements at the
+   beginning of ARRAY.  Uses COMPARE to compare elements, passing
+   AUX as auxiliary data. */
+void
+pop_heap (void *array, size_t count, size_t size,
+          algo_compare_func *compare, void *aux) 
+{
+  char *first = array;
+
+  expensive_assert (is_heap (array, count, size, compare, aux));
+  SWAP (first, first + (count - 1) * size, size);
+  heapify (first, count - 1, size, 1, compare, aux);
+  expensive_assert (count < 1 || is_heap (array, count - 1,
+                                          size, compare, aux));
+}
+
+/* Turns ARRAY, which contains COUNT elements of SIZE bytes, into
+   a heap.  Uses COMPARE to compare elements, passing AUX as
+   auxiliary data. */
+void
+make_heap (void *array, size_t count, size_t size,
+           algo_compare_func *compare, void *aux) 
+{
+  size_t idx;
+  
+  for (idx = count / 2; idx >= 1; idx--)
+    heapify (array, count, size, idx, compare, aux);
+  expensive_assert (count < 1 || is_heap (array, count, size, compare, aux));
+}
+
+/* ARRAY contains COUNT elements of SIZE bytes each.  Initially
+   all COUNT elements form a heap.  This function turns the heap
+   into a fully sorted array.  Uses COMPARE to compare elements,
+   passing AUX as auxiliary data. */
+void
+sort_heap (void *array, size_t count, size_t size,
+           algo_compare_func *compare, void *aux) 
+{
+  char *first = array;
+  size_t idx;
+
+  expensive_assert (is_heap (array, count, size, compare, aux));
+  for (idx = count; idx >= 2; idx--)
+    {
+      SWAP (first, first + (idx - 1) * size, size);
+      heapify (array, idx - 1, size, 1, compare, aux);
+    }
+  expensive_assert (is_sorted (array, count, size, compare, aux));
+}
+
+/* ARRAY contains COUNT elements of SIZE bytes each.  This
+   function tests whether ARRAY is a heap and returns 1 if so, 0
+   otherwise.  Uses COMPARE to compare elements, passing AUX as
+   auxiliary data. */
+int
+is_heap (const void *array, size_t count, size_t size,
+         algo_compare_func *compare, void *aux) 
+{
+  const char *first = array;
+  size_t child;
+  
+  for (child = 2; child <= count; child++)
+    {
+      size_t parent = child / 2;
+      if (compare (first + (parent - 1) * size,
+                   first + (child - 1) * size, aux) < 0)
+        return 0;
+    }
+
+  return 1;
+}