treewide: Replace <name>_cnt by n_<name>s and <name>_cap by allocated_<name>.
[pspp] / src / libpspp / array.c
1 /* PSPP - a program for statistical analysis.
2    Copyright (C) 1997-9, 2000, 2010, 2011 Free Software Foundation, Inc.
3
4    This program is free software: you can redistribute it and/or modify
5    it under the terms of the GNU General Public License as published by
6    the Free Software Foundation, either version 3 of the License, or
7    (at your option) any later version.
8
9    This program is distributed in the hope that it will be useful,
10    but WITHOUT ANY WARRANTY; without even the implied warranty of
11    MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the
12    GNU General Public License for more details.
13
14    You should have received a copy of the GNU General Public License
15    along with this program.  If not, see <http://www.gnu.org/licenses/>. */
16
17 /* Copyright (C) 2001 Free Software Foundation, Inc.
18
19    This file is part of the GNU ISO C++ Library.  This library is free
20    software; you can redistribute it and/or modify it under the
21    terms of the GNU General Public License as published by the
22    Free Software Foundation; either version 2, or (at your option)
23    any later version.
24
25    This library is distributed in the hope that it will be useful,
26    but WITHOUT ANY WARRANTY; without even the implied warranty of
27    MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the
28    GNU General Public License for more details.
29
30    You should have received a copy of the GNU General Public License
31    along with this program.  If not, see <http://www.gnu.org/licenses/>.
32
33    As a special exception, you may use this file as part of a free software
34    library without restriction.  Specifically, if other files instantiate
35    templates or use macros or inline functions from this file, or you compile
36    this file and link it with other files to produce an executable, this
37    file does not by itself cause the resulting executable to be covered by
38    the GNU General Public License.  This exception does not however
39    invalidate any other reasons why the executable file might be covered by
40    the GNU General Public License. */
41
42 /*
43  *
44  * Copyright (c) 1994
45  * Hewlett-Packard Company
46  *
47  * Permission to use, copy, modify, distribute and sell this software
48  * and its documentation for any purpose is hereby granted without fee,
49  * provided that the above copyright notice appear in all copies and
50  * that both that copyright notice and this permission notice appear
51  * in supporting documentation.  Hewlett-Packard Company makes no
52  * representations about the suitability of this software for any
53  * purpose.  It is provided "as is" without express or implied warranty.
54  *
55  *
56  * Copyright (c) 1996
57  * Silicon Graphics Computer Systems, Inc.
58  *
59  * Permission to use, copy, modify, distribute and sell this software
60  * and its documentation for any purpose is hereby granted without fee,
61  * provided that the above copyright notice appear in all copies and
62  * that both that copyright notice and this permission notice appear
63  * in supporting documentation.  Silicon Graphics makes no
64  * representations about the suitability of this software for any
65  * purpose.  It is provided "as is" without express or implied warranty.
66  */
67
68 /* Copyright (C) 1991, 1992, 1996, 1997, 1999 Free Software Foundation, Inc.
69    This file is part of the GNU C Library.
70    Written by Douglas C. Schmidt (schmidt@ics.uci.edu).
71
72    The GNU C Library is free software; you can redistribute it and/or
73    modify it under the terms of the GNU Lesser General Public
74    License as published by the Free Software Foundation; either
75    version 2.1 of the License, or (at your option) any later version.
76
77    The GNU C Library is distributed in the hope that it will be useful,
78    but WITHOUT ANY WARRANTY; without even the implied warranty of
79    MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the GNU
80    Lesser General Public License for more details.
81
82    You should have received a copy of the GNU Lesser General Public
83    License along with the GNU C Library.  If not, see
84    <http://www.gnu.org/licenses/>. */
85
86 #include <config.h>
87
88 #include "array.h"
89
90 #include <limits.h>
91 #include <stdlib.h>
92 #include <stdint.h>
93 #include <string.h>
94 #include "libpspp/assertion.h"
95
96 #include "gl/xalloc.h"
97 #include "gl/minmax.h"
98 \f
99 /* Finds an element in ARRAY, which contains COUNT elements of
100    SIZE bytes each, using COMPARE for comparisons.  Returns the
101    first element in ARRAY that matches TARGET, or a null pointer
102    on failure.  AUX is passed to each comparison as auxiliary
103    data. */
104 void *
105 find (const void *array, size_t count, size_t size,
106       const void *target,
107       algo_compare_func *compare, const void *aux)
108 {
109   const char *element = array;
110
111   while (count-- > 0)
112     {
113       if (compare (target, element, aux) == 0)
114         return (void *) element;
115
116       element += size;
117     }
118
119   return NULL;
120 }
121
122 /* Counts and return the number of elements in ARRAY, which
123    contains COUNT elements of SIZE bytes each, which are equal to
124    ELEMENT as compared with COMPARE.  AUX is passed as auxiliary
125    data to COMPARE. */
126 size_t
127 count_equal (const void *array, size_t count, size_t size,
128              const void *element,
129              algo_compare_func *compare, const void *aux)
130 {
131   const char *first = array;
132   size_t n_equals = 0;
133
134   while (count-- > 0)
135     {
136       if (compare (element, first, aux) == 0)
137         n_equals++;
138
139       first += size;
140     }
141
142   return n_equals;
143 }
144
145 /* Counts and return the number of elements in ARRAY, which
146    contains COUNT elements of SIZE bytes each, for which
147    PREDICATE returns true.  AUX is passed as auxiliary data to
148    PREDICATE. */
149 size_t
150 count_if (const void *array, size_t count, size_t size,
151           algo_predicate_func *predicate, const void *aux)
152 {
153   const char *first = array;
154   size_t n_trues = 0;
155
156   while (count-- > 0)
157     {
158       if (predicate (first, aux) != 0)
159         n_trues++;
160
161       first += size;
162     }
163
164   return n_trues;
165 }
166 \f
167 /* Byte-wise swap objects A and B, each SIZE bytes. */
168 static void
169 swap (void *a_, void *b_, size_t size)
170 {
171   uint8_t *a = a_;
172   uint8_t *b = b_;
173   while (size-- > 0)
174     {
175       uint8_t tmp = *a;
176       *a++ = *b;
177       *b++ = tmp;
178     }
179 }
180
181 /* Makes the elements in ARRAY unique, by moving up duplicates,
182    and returns the new number of elements in the array.  Sorted
183    arrays only.  Arguments same as for sort() above. */
184 size_t
185 unique (void *array, size_t count, size_t size,
186         algo_compare_func *compare, const void *aux)
187 {
188   char *first = array;
189   char *last = first + size * count;
190   char *result = array;
191
192   for (;;)
193     {
194       first += size;
195       if (first >= last)
196         {
197           assert (adjacent_find_equal (array, count,
198                                        size, compare, aux) == NULL);
199           return count;
200         }
201
202       if (compare (result, first, aux))
203         {
204           result += size;
205           if (result != first)
206             memcpy (result, first, size);
207         }
208       else
209         count--;
210     }
211 }
212
213 /* Helper function that calls sort(), then unique(). */
214 size_t
215 sort_unique (void *array, size_t count, size_t size,
216              algo_compare_func *compare, const void *aux)
217 {
218   sort (array, count, size, compare, aux);
219   return unique (array, count, size, compare, aux);
220 }
221 \f
222 /* Reorders ARRAY, which contains COUNT elements of SIZE bytes
223    each, so that the elements for which PREDICATE returns true
224    precede those for which PREDICATE returns zero.  AUX is
225    passed to each predicate as auxiliary data.  Returns the
226    number of elements for which PREDICATE returns true.  Not
227    stable. */
228 size_t
229 partition (void *array, size_t count, size_t size,
230            algo_predicate_func *predicate, const void *aux)
231 {
232   size_t n_trues = count;
233   char *first = array;
234   char *last = first + n_trues * size;
235
236   for (;;)
237     {
238       /* Move FIRST forward to point to first element that fails
239          PREDICATE. */
240       for (;;)
241         {
242           if (first == last)
243             goto done;
244           else if (!predicate (first, aux))
245             break;
246
247           first += size;
248         }
249       n_trues--;
250
251       /* Move LAST backward to point to last element that passes
252          PREDICATE. */
253       for (;;)
254         {
255           last -= size;
256
257           if (first == last)
258             goto done;
259           else if (predicate (last, aux))
260             break;
261           else
262             n_trues--;
263         }
264
265       /* By swapping FIRST and LAST we extend the starting and
266          ending sequences that pass and fail, respectively,
267          PREDICATE. */
268       swap (first, last, size);
269       first += size;
270     }
271
272  done:
273   assert (is_partitioned (array, count, size, n_trues, predicate, aux));
274   return n_trues;
275 }
276
277 /* Checks whether ARRAY, which contains COUNT elements of SIZE
278    bytes each, is partitioned such that PREDICATE returns true
279    for the first TRUE_CNT elements and zero for the remaining
280    elements.  AUX is passed as auxiliary data to PREDICATE. */
281 bool
282 is_partitioned (const void *array, size_t count, size_t size,
283                 size_t n_trues,
284                 algo_predicate_func *predicate, const void *aux)
285 {
286   const char *first = array;
287   size_t idx;
288
289   assert (n_trues <= count);
290   for (idx = 0; idx < n_trues; idx++)
291     if (predicate (first + idx * size, aux) == 0)
292       return false;
293   for (idx = n_trues; idx < count; idx++)
294     if (predicate (first + idx * size, aux) != 0)
295       return false;
296   return true;
297 }
298 \f
299 /* Copies the COUNT elements of SIZE bytes each from ARRAY to
300    RESULT, except that elements for which PREDICATE is false are
301    not copied.  Returns the number of elements copied.  AUX is
302    passed to PREDICATE as auxiliary data.  */
303 size_t
304 copy_if (const void *array, size_t count, size_t size,
305          void *result,
306          algo_predicate_func *predicate, const void *aux)
307 {
308   const char *input = array;
309   const char *last = input + size * count;
310   char *output = result;
311   size_t n_nonzeros = 0;
312
313   while (input < last)
314     {
315       if (predicate (input, aux))
316         {
317           memcpy (output, input, size);
318           output += size;
319           n_nonzeros++;
320         }
321
322       input += size;
323     }
324
325   assert (n_nonzeros == count_if (array, count, size, predicate, aux));
326   assert (n_nonzeros == count_if (result, n_nonzeros, size, predicate, aux));
327
328   return n_nonzeros;
329 }
330
331 /* Removes N elements starting at IDX from ARRAY, which consists
332    of COUNT elements of SIZE bytes each, by shifting the elements
333    following them, if any, into its position. */
334 void
335 remove_range (void *array_, size_t count, size_t size,
336               size_t idx, size_t n)
337 {
338   char *array = array_;
339
340   assert (array != NULL);
341   assert (idx <= count);
342   assert (idx + n <= count);
343
344   if (idx + n < count)
345     memmove (array + idx * size, array + (idx + n) * size,
346              size * (count - idx - n));
347 }
348
349 /* Removes element IDX from ARRAY, which consists of COUNT
350    elements of SIZE bytes each, by shifting the elements
351    following it, if any, into its position. */
352 void
353 remove_element (void *array, size_t count, size_t size,
354                 size_t idx)
355 {
356   remove_range (array, count, size, idx, 1);
357 }
358
359 /* Makes room for N elements starting at IDX in ARRAY, which
360    initially consists of COUNT elements of SIZE bytes each, by
361    shifting elements IDX...COUNT (exclusive) to the right by N
362    positions. */
363 void
364 insert_range (void *array_, size_t count, size_t size,
365               size_t idx, size_t n)
366 {
367   char *array = array_;
368
369   assert (idx <= count);
370   memmove (array + (idx + n) * size, array + idx * size, (count - idx) * size);
371 }
372
373 /* Makes room for a new element at IDX in ARRAY, which initially
374    consists of COUNT elements of SIZE bytes each, by shifting
375    elements IDX...COUNT (exclusive) to the right by one
376    position. */
377 void
378 insert_element (void *array, size_t count, size_t size,
379                 size_t idx)
380 {
381   insert_range (array, count, size, idx, 1);
382 }
383
384 /* Moves an element in ARRAY, which consists of COUNT elements of
385    SIZE bytes each, from OLD_IDX to NEW_IDX, shifting around
386    other elements as needed.  Runs in O(abs(OLD_IDX - NEW_IDX))
387    time. */
388 void
389 move_element (void *array_, size_t count, size_t size,
390               size_t old_idx, size_t new_idx)
391 {
392   assert (array_ != NULL || count == 0);
393   assert (old_idx < count);
394   assert (new_idx < count);
395
396   if (old_idx != new_idx)
397     {
398       char *array = array_;
399       char *element = xmalloc (size);
400       char *new = array + new_idx * size;
401       char *old = array + old_idx * size;
402
403       memcpy (element, old, size);
404       if (new < old)
405         memmove (new + size, new, (old_idx - new_idx) * size);
406       else
407         memmove (old, old + size, (new_idx - old_idx) * size);
408       memcpy (new, element, size);
409
410       free (element);
411     }
412 }
413
414 /* Moves N elements in ARRAY starting at OLD_IDX, which consists
415    of COUNT elements of SIZE bytes each, so that they now start
416    at NEW_IDX, shifting around other elements as needed. */
417 void
418 move_range (void *array_, size_t count, size_t size,
419             size_t old_idx, size_t new_idx, size_t n)
420 {
421   assert (array_ != NULL || count == 0);
422   assert (n <= count);
423   assert (old_idx + n <= count);
424   assert (new_idx + n <= count);
425
426   if (old_idx != new_idx && n > 0)
427     {
428       char *array = array_;
429       char *range = xmalloc (size * n);
430       char *new = array + new_idx * size;
431       char *old = array + old_idx * size;
432
433       memcpy (range, old, size * n);
434       if (new < old)
435         memmove (new + size * n, new, (old_idx - new_idx) * size);
436       else
437         memmove (old, old + size * n, (new_idx - old_idx) * size);
438       memcpy (new, range, size * n);
439
440       free (range);
441     }
442 }
443
444 /* A predicate and its auxiliary data. */
445 struct pred_aux
446   {
447     algo_predicate_func *predicate;
448     const void *aux;
449   };
450
451 static bool
452 not (const void *data, const void *pred_aux_)
453 {
454   const struct pred_aux *pred_aux = pred_aux_;
455
456   return !pred_aux->predicate (data, pred_aux->aux);
457 }
458
459 /* Removes elements equal to ELEMENT from ARRAY, which consists
460    of COUNT elements of SIZE bytes each.  Returns the number of
461    remaining elements.  AUX is passed to COMPARE as auxiliary
462    data. */
463 size_t
464 remove_equal (void *array, size_t count, size_t size,
465               void *element,
466               algo_compare_func *compare, const void *aux)
467 {
468   char *first = array;
469   char *last = first + count * size;
470   char *result;
471
472   for (;;)
473     {
474       if (first >= last)
475         goto done;
476       if (compare (first, element, aux) == 0)
477         break;
478
479       first += size;
480     }
481
482   result = first;
483   count--;
484   for (;;)
485     {
486       first += size;
487       if (first >= last)
488         goto done;
489
490       if (compare (first, element, aux) == 0)
491         {
492           count--;
493           continue;
494         }
495
496       memcpy (result, first, size);
497       result += size;
498     }
499
500  done:
501   assert (count_equal (array, count, size, element, compare, aux) == 0);
502   return count;
503 }
504
505 /* Copies the COUNT elements of SIZE bytes each from ARRAY to
506    RESULT, except that elements for which PREDICATE is true are
507    not copied.  Returns the number of elements copied.  AUX is
508    passed to PREDICATE as auxiliary data.  */
509 size_t
510 remove_copy_if (const void *array, size_t count, size_t size,
511                 void *result,
512                 algo_predicate_func *predicate, const void *aux)
513 {
514   struct pred_aux pred_aux;
515   pred_aux.predicate = predicate;
516   pred_aux.aux = aux;
517   return copy_if (array, count, size, result, not, &pred_aux);
518 }
519 \f
520 /* Searches ARRAY, which contains COUNT of SIZE bytes each, using
521    a binary search.  Returns any element that equals VALUE, if
522    one exists, or a null pointer otherwise.  ARRAY must ordered
523    according to COMPARE.  AUX is passed to COMPARE as auxiliary
524    data. */
525 void *
526 binary_search (const void *array, size_t count, size_t size,
527                void *value,
528                algo_compare_func *compare, const void *aux)
529 {
530   assert (array != NULL || count == 0);
531   assert (count <= INT_MAX);
532   assert (compare != NULL);
533
534   if (count != 0)
535     {
536       const char *first = array;
537       int low = 0;
538       int high = count - 1;
539
540       while (low <= high)
541         {
542           int middle = (low + high) / 2;
543           const char *element = first + middle * size;
544           int cmp = compare (value, element, aux);
545
546           if (cmp > 0)
547             low = middle + 1;
548           else if (cmp < 0)
549             high = middle - 1;
550           else
551             return (void *) element;
552         }
553     }
554
555   expensive_assert (find (array, count, size, value, compare, aux) == NULL);
556   return NULL;
557 }
558 \f
559 /* Lexicographically compares ARRAY1, which contains COUNT1
560    elements of SIZE bytes each, to ARRAY2, which contains COUNT2
561    elements of SIZE bytes, according to COMPARE.  Returns a
562    strcmp()-type result.  AUX is passed to COMPARE as auxiliary
563    data. */
564 int
565 lexicographical_compare_3way (const void *array1, size_t count1,
566                               const void *array2, size_t count2,
567                               size_t size,
568                               algo_compare_func *compare, const void *aux)
569 {
570   const char *first1 = array1;
571   const char *first2 = array2;
572   size_t min_count = count1 < count2 ? count1 : count2;
573
574   while (min_count > 0)
575     {
576       int cmp = compare (first1, first2, aux);
577       if (cmp != 0)
578         return cmp;
579
580       first1 += size;
581       first2 += size;
582       min_count--;
583     }
584
585   return count1 < count2 ? -1 : count1 > count2;
586 }
587 \f
588 /* If you consider tuning this algorithm, you should consult first:
589    Engineering a sort function; Jon Bentley and M. Douglas McIlroy;
590    Software - Practice and Experience; Vol. 23 (11), 1249-1265, 1993.  */
591
592 #include <limits.h>
593 #include <stdlib.h>
594 #include <string.h>
595
596 /* Discontinue quicksort algorithm when partition gets below this size.
597    This particular magic number was chosen to work best on a Sun 4/260. */
598 #define MAX_THRESH 4
599
600 /* Stack node declarations used to store unfulfilled partition obligations. */
601 typedef struct
602   {
603     char *lo;
604     char *hi;
605   } stack_node;
606
607 /* The next 4 #defines implement a very fast in-line stack abstraction. */
608 /* The stack needs log (total_elements) entries (we could even subtract
609    log(MAX_THRESH)).  Since total_elements has type size_t, we get as
610    upper bound for log (total_elements):
611    bits per byte (CHAR_BIT) * sizeof(size_t).  */
612 #define STACK_SIZE      (CHAR_BIT * sizeof(size_t))
613 #define PUSH(low, high) ((void) ((top->lo = (low)), (top->hi = (high)), ++top))
614 #define POP(low, high)  ((void) (--top, (low = top->lo), (high = top->hi)))
615 #define STACK_NOT_EMPTY (stack < top)
616
617
618 /* Order size using quicksort.  This implementation incorporates
619    four optimizations discussed in Sedgewick:
620
621    1. Non-recursive, using an explicit stack of pointer that store the
622       next array partition to sort.  To save time, this maximum amount
623       of space required to store an array of SIZE_MAX is allocated on the
624       stack.  Assuming a 32-bit (64 bit) integer for size_t, this needs
625       only 32 * sizeof(stack_node) == 256 bytes (for 64 bit: 1024 bytes).
626       Pretty cheap, actually.
627
628    2. Chose the pivot element using a median-of-three decision tree.
629       This reduces the probability of selecting a bad pivot value and
630       eliminates certain extraneous comparisons.
631
632    3. Only quicksorts TOTAL_ELEMS / MAX_THRESH partitions, leaving
633       insertion sort to order the MAX_THRESH items within each partition.
634       This is a big win, since insertion sort is faster for small, mostly
635       sorted array segments.
636
637    4. The larger of the two sub-partitions is always pushed onto the
638       stack first, with the algorithm then concentrating on the
639       smaller partition.  This *guarantees* no more than log (total_elems)
640       stack size is needed (actually O(1) in this case)!  */
641
642 void
643 sort (void *array, size_t count, size_t size,
644       algo_compare_func *compare, const void *aux)
645 {
646   char *const first = array;
647   const size_t max_thresh = MAX_THRESH * size;
648
649   if (count == 0)
650     /* Avoid lossage with unsigned arithmetic below.  */
651     return;
652
653   if (count > MAX_THRESH)
654     {
655       char *lo = first;
656       char *hi = &lo[size * (count - 1)];
657       stack_node stack[STACK_SIZE];
658       stack_node *top = stack + 1;
659
660       while (STACK_NOT_EMPTY)
661         {
662           char *left_ptr;
663           char *right_ptr;
664
665           /* Select median value from among LO, MID, and HI. Rearrange
666              LO and HI so the three values are sorted. This lowers the
667              probability of picking a pathological pivot value and
668              skips a comparison for both the LEFT_PTR and RIGHT_PTR in
669              the while loops. */
670
671           char *mid = lo + size * ((hi - lo) / size >> 1);
672
673           if (compare (mid, lo, aux) < 0)
674             swap (mid, lo, size);
675           if (compare (hi, mid, aux) < 0)
676             swap (mid, hi, size);
677           else
678             goto jump_over;
679           if (compare (mid, lo, aux) < 0)
680             swap (mid, lo, size);
681         jump_over:;
682
683           left_ptr  = lo + size;
684           right_ptr = hi - size;
685
686           /* Here's the famous ``collapse the walls'' section of quicksort.
687              Gotta like those tight inner loops!  They are the main reason
688              that this algorithm runs much faster than others. */
689           do
690             {
691               while (compare (left_ptr, mid, aux) < 0)
692                 left_ptr += size;
693
694               while (compare (mid, right_ptr, aux) < 0)
695                 right_ptr -= size;
696
697               if (left_ptr < right_ptr)
698                 {
699                   swap (left_ptr, right_ptr, size);
700                   if (mid == left_ptr)
701                     mid = right_ptr;
702                   else if (mid == right_ptr)
703                     mid = left_ptr;
704                   left_ptr += size;
705                   right_ptr -= size;
706                 }
707               else if (left_ptr == right_ptr)
708                 {
709                   left_ptr += size;
710                   right_ptr -= size;
711                   break;
712                 }
713             }
714           while (left_ptr <= right_ptr);
715
716           /* Set up pointers for next iteration.  First determine whether
717              left and right partitions are below the threshold size.  If so,
718              ignore one or both.  Otherwise, push the larger partition's
719              bounds on the stack and continue sorting the smaller one. */
720
721           if ((size_t) (right_ptr - lo) <= max_thresh)
722             {
723               if ((size_t) (hi - left_ptr) <= max_thresh)
724                 /* Ignore both small partitions. */
725                 POP (lo, hi);
726               else
727                 /* Ignore small left partition. */
728                 lo = left_ptr;
729             }
730           else if ((size_t) (hi - left_ptr) <= max_thresh)
731             /* Ignore small right partition. */
732             hi = right_ptr;
733           else if ((right_ptr - lo) > (hi - left_ptr))
734             {
735               /* Push larger left partition indices. */
736               PUSH (lo, right_ptr);
737               lo = left_ptr;
738             }
739           else
740             {
741               /* Push larger right partition indices. */
742               PUSH (left_ptr, hi);
743               hi = right_ptr;
744             }
745         }
746     }
747
748   /* Once the FIRST array is partially sorted by quicksort the rest
749      is completely sorted using insertion sort, since this is efficient
750      for partitions below MAX_THRESH size. FIRST points to the beginning
751      of the array to sort, and END_PTR points at the very last element in
752      the array (*not* one beyond it!). */
753
754   {
755     char *const end_ptr = &first[size * (count - 1)];
756     char *tmp_ptr = first;
757     char *thresh = MIN (end_ptr, first + max_thresh);
758     register char *run_ptr;
759
760     /* Find smallest element in first threshold and place it at the
761        array's beginning.  This is the smallest array element,
762        and the operation speeds up insertion sort's inner loop. */
763
764     for (run_ptr = tmp_ptr + size; run_ptr <= thresh; run_ptr += size)
765       if (compare (run_ptr, tmp_ptr, aux) < 0)
766         tmp_ptr = run_ptr;
767
768     if (tmp_ptr != first)
769       swap (tmp_ptr, first, size);
770
771     /* Insertion sort, running from left-hand-side up to right-hand-side.  */
772
773     run_ptr = first + size;
774     while ((run_ptr += size) <= end_ptr)
775       {
776         tmp_ptr = run_ptr - size;
777         while (compare (run_ptr, tmp_ptr, aux) < 0)
778           tmp_ptr -= size;
779
780         tmp_ptr += size;
781         if (tmp_ptr != run_ptr)
782           {
783             char *trav;
784
785             trav = run_ptr + size;
786             while (--trav >= run_ptr)
787               {
788                 char c = *trav;
789                 char *hi, *lo;
790
791                 for (hi = lo = trav; (lo -= size) >= tmp_ptr; hi = lo)
792                   *hi = *lo;
793                 *hi = c;
794               }
795           }
796       }
797   }
798
799   assert (is_sorted (array, count, size, compare, aux));
800 }
801
802 /* Tests whether ARRAY, which contains COUNT elements of SIZE
803    bytes each, is sorted in order according to COMPARE.  AUX is
804    passed to COMPARE as auxiliary data. */
805 bool
806 is_sorted (const void *array, size_t count, size_t size,
807            algo_compare_func *compare, const void *aux)
808 {
809   const char *first = array;
810   size_t idx;
811
812   for (idx = 0; idx + 1 < count; idx++)
813     if (compare (first + idx * size, first + (idx + 1) * size, aux) > 0)
814       return false;
815
816   return true;
817 }
818 \f
819 /* Computes the generalized set difference, ARRAY1 minus ARRAY2,
820    into RESULT, and returns the number of elements written to
821    RESULT.  If a value appears M times in ARRAY1 and N times in
822    ARRAY2, then it will appear max(M - N, 0) in RESULT.  ARRAY1
823    and ARRAY2 must be sorted, and RESULT is sorted and stable.
824    ARRAY1 consists of COUNT1 elements, ARRAY2 of COUNT2 elements,
825    each SIZE bytes.  AUX is passed to COMPARE as auxiliary
826    data. */
827 size_t set_difference (const void *array1, size_t count1,
828                        const void *array2, size_t count2,
829                        size_t size,
830                        void *result_,
831                        algo_compare_func *compare, const void *aux)
832 {
833   const char *first1 = array1;
834   const char *last1 = first1 + count1 * size;
835   const char *first2 = array2;
836   const char *last2 = first2 + count2 * size;
837   char *result = result_;
838   size_t result_count = 0;
839
840   while (first1 != last1 && first2 != last2)
841     {
842       int cmp = compare (first1, first2, aux);
843       if (cmp < 0)
844         {
845           memcpy (result, first1, size);
846           first1 += size;
847           result += size;
848           result_count++;
849         }
850       else if (cmp > 0)
851         first2 += size;
852       else
853         {
854           first1 += size;
855           first2 += size;
856         }
857     }
858
859   while (first1 != last1)
860     {
861       memcpy (result, first1, size);
862       first1 += size;
863       result += size;
864       result_count++;
865     }
866
867   return result_count;
868 }
869 \f
870 /* Finds the first pair of adjacent equal elements in ARRAY,
871    which has COUNT elements of SIZE bytes.  Returns the first
872    element in ARRAY such that COMPARE returns zero when it and
873    its successor element are compared, or a null pointer if no
874    such element exists.  AUX is passed to COMPARE as auxiliary
875    data. */
876 void *
877 adjacent_find_equal (const void *array, size_t count, size_t size,
878                      algo_compare_func *compare, const void *aux)
879 {
880   const char *first = array;
881   const char *last = first + count * size;
882
883   while (first < last && first + size < last)
884     {
885       if (compare (first, first + size, aux) == 0)
886         return (void *) first;
887       first += size;
888     }
889
890   return NULL;
891 }
892 \f
893 /* ARRAY contains COUNT elements of SIZE bytes each.  Initially
894    the first COUNT - 1 elements of these form a heap, followed by
895    a single element not part of the heap.  This function adds the
896    final element, forming a heap of COUNT elements in ARRAY.
897    Uses COMPARE to compare elements, passing AUX as auxiliary
898    data. */
899 void
900 push_heap (void *array, size_t count, size_t size,
901            algo_compare_func *compare, const void *aux)
902 {
903   char *first = array;
904   size_t i;
905
906   expensive_assert (count < 1 || is_heap (array, count - 1,
907                                           size, compare, aux));
908   for (i = count; i > 1; i /= 2)
909     {
910       char *parent = first + (i / 2 - 1) * size;
911       char *element = first + (i - 1) * size;
912       if (compare (parent, element, aux) < 0)
913         swap (parent, element, size);
914       else
915         break;
916     }
917   expensive_assert (is_heap (array, count, size, compare, aux));
918 }
919
920 /* ARRAY contains COUNT elements of SIZE bytes each.  Initially
921    the children of ARRAY[idx - 1] are heaps, but ARRAY[idx - 1]
922    may be smaller than its children.  This function fixes that,
923    so that ARRAY[idx - 1] itself is a heap.  Uses COMPARE to
924    compare elements, passing AUX as auxiliary data. */
925 static void
926 heapify (void *array, size_t count, size_t size,
927          size_t idx,
928          algo_compare_func *compare, const void *aux)
929 {
930   char *first = array;
931
932   for (;;)
933     {
934       size_t left = 2 * idx;
935       size_t right = left + 1;
936       size_t largest = idx;
937
938       if (left <= count
939           && compare (first + size * (left - 1),
940                       first + size * (idx - 1), aux) > 0)
941         largest = left;
942
943       if (right <= count
944           && compare (first + size * (right - 1),
945                       first + size * (largest - 1), aux) > 0)
946         largest = right;
947
948       if (largest == idx)
949         break;
950
951       swap (first + size * (idx - 1), first + size * (largest - 1), size);
952       idx = largest;
953     }
954 }
955
956 /* ARRAY contains COUNT elements of SIZE bytes each.  Initially
957    all COUNT elements form a heap.  This function moves the
958    largest element in the heap to the final position in ARRAY and
959    reforms a heap of the remaining COUNT - 1 elements at the
960    beginning of ARRAY.  Uses COMPARE to compare elements, passing
961    AUX as auxiliary data. */
962 void
963 pop_heap (void *array, size_t count, size_t size,
964           algo_compare_func *compare, const void *aux)
965 {
966   char *first = array;
967
968   expensive_assert (is_heap (array, count, size, compare, aux));
969   swap (first, first + (count - 1) * size, size);
970   heapify (first, count - 1, size, 1, compare, aux);
971   expensive_assert (count < 1 || is_heap (array, count - 1,
972                                           size, compare, aux));
973 }
974
975 /* Turns ARRAY, which contains COUNT elements of SIZE bytes, into
976    a heap.  Uses COMPARE to compare elements, passing AUX as
977    auxiliary data. */
978 void
979 make_heap (void *array, size_t count, size_t size,
980            algo_compare_func *compare, const void *aux)
981 {
982   size_t idx;
983
984   for (idx = count / 2; idx >= 1; idx--)
985     heapify (array, count, size, idx, compare, aux);
986   expensive_assert (count < 1 || is_heap (array, count, size, compare, aux));
987 }
988
989 /* ARRAY contains COUNT elements of SIZE bytes each.  Initially
990    all COUNT elements form a heap.  This function turns the heap
991    into a fully sorted array.  Uses COMPARE to compare elements,
992    passing AUX as auxiliary data. */
993 void
994 sort_heap (void *array, size_t count, size_t size,
995            algo_compare_func *compare, const void *aux)
996 {
997   char *first = array;
998   size_t idx;
999
1000   expensive_assert (is_heap (array, count, size, compare, aux));
1001   for (idx = count; idx >= 2; idx--)
1002     {
1003       swap (first, first + (idx - 1) * size, size);
1004       heapify (array, idx - 1, size, 1, compare, aux);
1005     }
1006   expensive_assert (is_sorted (array, count, size, compare, aux));
1007 }
1008
1009 /* ARRAY contains COUNT elements of SIZE bytes each.  This
1010    function tests whether ARRAY is a heap and returns true if so,
1011    false otherwise.  Uses COMPARE to compare elements, passing
1012    AUX as auxiliary data. */
1013 bool
1014 is_heap (const void *array, size_t count, size_t size,
1015          algo_compare_func *compare, const void *aux)
1016 {
1017   const char *first = array;
1018   size_t child;
1019
1020   for (child = 2; child <= count; child++)
1021     {
1022       size_t parent = child / 2;
1023       if (compare (first + (parent - 1) * size,
1024                    first + (child - 1) * size, aux) < 0)
1025         return false;
1026     }
1027
1028   return true;
1029 }
1030
1031 /* Reverses the order of ARRAY, which contains COUNT elements of SIZE bytes
1032    each. */
1033 void
1034 reverse_array (void *array_, size_t count, size_t size)
1035 {
1036   uint8_t *array = array_;
1037   uint8_t *first = array;
1038   uint8_t *last = array + (count - 1) * size;
1039   for (size_t i = 0; i < count / 2; i++)
1040     {
1041       swap (first, last, size);
1042       first += size;
1043       last -= size;
1044     }
1045 }