[v6,4/6] stdlib: Implement introsort for qsort (BZ 19305)
Checks
Context |
Check |
Description |
redhat-pt-bot/TryBot-apply_patch |
success
|
Patch applied to master at the time it was sent
|
linaro-tcwg-bot/tcwg_glibc_build--master-arm |
success
|
Testing passed
|
linaro-tcwg-bot/tcwg_glibc_build--master-aarch64 |
success
|
Testing passed
|
linaro-tcwg-bot/tcwg_glibc_check--master-arm |
success
|
Testing passed
|
linaro-tcwg-bot/tcwg_glibc_check--master-aarch64 |
success
|
Testing passed
|
Commit Message
This patch makes the quicksort implementation to acts as introsort, to
avoid worse-case performance (and thus making it O(nlog n)). It switch
to heapsort when the depth level reaches 2*log2(total elements). The
heapsort is a textbook implementation.
Checked on x86_64-linux-gnu.
---
stdlib/qsort.c | 85 ++++++++++++++++++++++++++++++++++++++++++++++----
1 file changed, 79 insertions(+), 6 deletions(-)
Comments
On Tue, Jul 18, 2023 at 9:19 AM Adhemerval Zanella via Libc-alpha
<libc-alpha@sourceware.org> wrote:
>
> This patch makes the quicksort implementation to acts as introsort, to
> avoid worse-case performance (and thus making it O(nlog n)). It switch
> to heapsort when the depth level reaches 2*log2(total elements). The
> heapsort is a textbook implementation.
>
> Checked on x86_64-linux-gnu.
> ---
> stdlib/qsort.c | 85 ++++++++++++++++++++++++++++++++++++++++++++++----
> 1 file changed, 79 insertions(+), 6 deletions(-)
>
> diff --git a/stdlib/qsort.c b/stdlib/qsort.c
> index 640896a598..054c900b02 100644
> --- a/stdlib/qsort.c
> +++ b/stdlib/qsort.c
> @@ -119,6 +119,7 @@ typedef struct
> {
> char *lo;
> char *hi;
> + size_t depth;
You don't actually need to track depth in the "stack"
when you pop hi/lo you can take log2(hi - lo) to compute
what depth counter "should" be i.e
size_t loop_cnt = WORD_SIZE - clz(hi - lo);
then if --loop_cnt == 0 -> introsort that chunk.
> } stack_node;
>
> /* The stack needs log (total_elements) entries (we could even subtract
> @@ -128,22 +129,83 @@ typedef struct
> enum { STACK_SIZE = CHAR_BIT * sizeof (size_t) };
>
> static inline stack_node *
> -push (stack_node *top, char *lo, char *hi)
> +push (stack_node *top, char *lo, char *hi, size_t depth)
> {
> top->lo = lo;
> top->hi = hi;
> + top->depth = depth;
> return ++top;
> }
>
> static inline stack_node *
> -pop (stack_node *top, char **lo, char **hi)
> +pop (stack_node *top, char **lo, char **hi, size_t *depth)
> {
> --top;
> *lo = top->lo;
> *hi = top->hi;
> + *depth = top->depth;
> return top;
> }
>
> +/* A fast, small, non-recursive O(nlog n) heapsort, adapted from Linux
> + lib/sort.c. Used on introsort implementation as a fallback routine with
> + worst-case performance of O(nlog n) and worst-case space complexity of
> + O(1). */
> +
> +static inline void
> +siftdown (void *base, size_t size, size_t k, size_t n,
> + enum swap_type_t swap_type, __compar_d_fn_t cmp, void *arg)
> +{
> + while (k <= n / 2)
> + {
> + size_t j = 2 * k;
> + if (j < n && cmp (base + (j * size), base + ((j + 1) * size), arg) < 0)
> + j++;
> +
> + if (cmp (base + (k * size), base + (j * size), arg) >= 0)
> + break;
> +
> + do_swap (base + (size * j), base + (k * size), size, swap_type);
> + k = j;
> + }
> +}
> +
> +static inline void
> +heapify (void *base, size_t size, size_t n, enum swap_type_t swap_type,
> + __compar_d_fn_t cmp, void *arg)
> +{
> + size_t k = n / 2;
> + while (1)
> + {
> + siftdown (base, size, k, n, swap_type, cmp, arg);
> + if (k-- == 0)
> + break;
> + }
> +}
> +
> +static void
> +heapsort_r (void *base, void *end, size_t size, enum swap_type_t swap_type,
> + __compar_d_fn_t cmp, void *arg)
> +{
> + const size_t count = ((uintptr_t) end - (uintptr_t) base) / size;
> +
> + if (count < 2)
> + return;
> +
> + size_t n = count - 1;
> +
> + /* Build the binary heap, largest value at the base[0]. */
> + heapify (base, size, n, swap_type, cmp, arg);
> +
> + /* On each iteration base[0:n] is the binary heap, while base[n:count]
> + is sorted. */
> + while (n > 0)
> + {
> + do_swap (base, base + (n * size), size, swap_type);
> + n--;
> + siftdown (base, size, 0, n, swap_type, cmp, arg);
> + }
> +}
>
> /* Order size using quicksort. This implementation incorporates
> four optimizations discussed in Sedgewick:
> @@ -229,7 +291,7 @@ _quicksort (void *const pbase, size_t total_elems, size_t size,
>
> const size_t max_thresh = MAX_THRESH * size;
>
> - if (total_elems == 0)
> + if (total_elems <= 1)
> /* Avoid lossage with unsigned arithmetic below. */
> return;
>
> @@ -241,6 +303,10 @@ _quicksort (void *const pbase, size_t total_elems, size_t size,
> else
> swap_type = SWAP_BYTES;
>
> + /* Maximum depth before quicksort switches to heapsort. */
> + size_t depth = 2 * (sizeof (size_t) * CHAR_BIT - 1
> + - __builtin_clzl (total_elems));
> +
> if (total_elems > MAX_THRESH)
> {
> char *lo = base_ptr;
> @@ -250,6 +316,13 @@ _quicksort (void *const pbase, size_t total_elems, size_t size,
>
> while (stack < top)
> {
> + if (depth == 0)
> + {
> + heapsort_r (lo, hi, size, swap_type, cmp, arg);
> + top = pop (top, &lo, &hi, &depth);
> + continue;
> + }
> +
> char *left_ptr;
> char *right_ptr;
>
> @@ -313,7 +386,7 @@ _quicksort (void *const pbase, size_t total_elems, size_t size,
> {
> if ((size_t) (hi - left_ptr) <= max_thresh)
> /* Ignore both small partitions. */
> - top = pop (top, &lo, &hi);
> + top = pop (top, &lo, &hi, &depth);
> else
> /* Ignore small left partition. */
> lo = left_ptr;
> @@ -324,13 +397,13 @@ _quicksort (void *const pbase, size_t total_elems, size_t size,
> else if ((right_ptr - lo) > (hi - left_ptr))
> {
> /* Push larger left partition indices. */
> - top = push (top, lo, right_ptr);
> + top = push (top, lo, right_ptr, depth - 1);
> lo = left_ptr;
> }
> else
> {
> /* Push larger right partition indices. */
> - top = push (top, left_ptr, hi);
> + top = push (top, left_ptr, hi, depth - 1);
> hi = right_ptr;
> }
> }
> --
> 2.34.1
>
On 18/07/23 14:37, Noah Goldstein wrote:
> On Tue, Jul 18, 2023 at 9:19 AM Adhemerval Zanella via Libc-alpha
> <libc-alpha@sourceware.org> wrote:
>>
>> This patch makes the quicksort implementation to acts as introsort, to
>> avoid worse-case performance (and thus making it O(nlog n)). It switch
>> to heapsort when the depth level reaches 2*log2(total elements). The
>> heapsort is a textbook implementation.
>>
>> Checked on x86_64-linux-gnu.
>> ---
>> stdlib/qsort.c | 85 ++++++++++++++++++++++++++++++++++++++++++++++----
>> 1 file changed, 79 insertions(+), 6 deletions(-)
>>
>> diff --git a/stdlib/qsort.c b/stdlib/qsort.c
>> index 640896a598..054c900b02 100644
>> --- a/stdlib/qsort.c
>> +++ b/stdlib/qsort.c
>> @@ -119,6 +119,7 @@ typedef struct
>> {
>> char *lo;
>> char *hi;
>> + size_t depth;
>
> You don't actually need to track depth in the "stack"
> when you pop hi/lo you can take log2(hi - lo) to compute
> what depth counter "should" be i.e
I don't think we can take the assumption quicksort will always iterate the
input array as a binary tree, so we can't really use log2(hi - lo) as the
depth.
However I found another issue, which a change to use 'stack_node *top = stack + 1'
removed the initial depth initialization.
>
> size_t loop_cnt = WORD_SIZE - clz(hi - lo);
>
> then if --loop_cnt == 0 -> introsort that chunk
On Tue, Jul 18, 2023 at 2:08 PM Adhemerval Zanella Netto
<adhemerval.zanella@linaro.org> wrote:
>
>
>
> On 18/07/23 14:37, Noah Goldstein wrote:
> > On Tue, Jul 18, 2023 at 9:19 AM Adhemerval Zanella via Libc-alpha
> > <libc-alpha@sourceware.org> wrote:
> >>
> >> This patch makes the quicksort implementation to acts as introsort, to
> >> avoid worse-case performance (and thus making it O(nlog n)). It switch
> >> to heapsort when the depth level reaches 2*log2(total elements). The
> >> heapsort is a textbook implementation.
> >>
> >> Checked on x86_64-linux-gnu.
> >> ---
> >> stdlib/qsort.c | 85 ++++++++++++++++++++++++++++++++++++++++++++++----
> >> 1 file changed, 79 insertions(+), 6 deletions(-)
> >>
> >> diff --git a/stdlib/qsort.c b/stdlib/qsort.c
> >> index 640896a598..054c900b02 100644
> >> --- a/stdlib/qsort.c
> >> +++ b/stdlib/qsort.c
> >> @@ -119,6 +119,7 @@ typedef struct
> >> {
> >> char *lo;
> >> char *hi;
> >> + size_t depth;
> >
> > You don't actually need to track depth in the "stack"
> > when you pop hi/lo you can take log2(hi - lo) to compute
> > what depth counter "should" be i.e
>
> I don't think we can take the assumption quicksort will always iterate the
> input array as a binary tree, so we can't really use log2(hi - lo) as the
> depth.
> However I found another issue, which a change to use 'stack_node *top = stack + 1'
> removed the initial depth initialization.
>
Bah, sorry missed your reply!
What i mean is not that you take `log2(hi - lo)` as the current depth,
but whenever
you pop a new set of bounds you set `max_depth = log2(hi - lo)`. Then
in the loop if
`--max_depth <= 0` we are in pessimistic state -> do introsort.
Think this works to avoid the O(N^2) state and saves a lot of state tracking.
> >
> > size_t loop_cnt = WORD_SIZE - clz(hi - lo);
> >
> > then if --loop_cnt == 0 -> introsort that chunk
On 29/08/23 04:45, Noah Goldstein wrote:
> On Tue, Jul 18, 2023 at 2:08 PM Adhemerval Zanella Netto
> <adhemerval.zanella@linaro.org> wrote:
>>
>>
>>
>> On 18/07/23 14:37, Noah Goldstein wrote:
>>> On Tue, Jul 18, 2023 at 9:19 AM Adhemerval Zanella via Libc-alpha
>>> <libc-alpha@sourceware.org> wrote:
>>>>
>>>> This patch makes the quicksort implementation to acts as introsort, to
>>>> avoid worse-case performance (and thus making it O(nlog n)). It switch
>>>> to heapsort when the depth level reaches 2*log2(total elements). The
>>>> heapsort is a textbook implementation.
>>>>
>>>> Checked on x86_64-linux-gnu.
>>>> ---
>>>> stdlib/qsort.c | 85 ++++++++++++++++++++++++++++++++++++++++++++++----
>>>> 1 file changed, 79 insertions(+), 6 deletions(-)
>>>>
>>>> diff --git a/stdlib/qsort.c b/stdlib/qsort.c
>>>> index 640896a598..054c900b02 100644
>>>> --- a/stdlib/qsort.c
>>>> +++ b/stdlib/qsort.c
>>>> @@ -119,6 +119,7 @@ typedef struct
>>>> {
>>>> char *lo;
>>>> char *hi;
>>>> + size_t depth;
>>>
>>> You don't actually need to track depth in the "stack"
>>> when you pop hi/lo you can take log2(hi - lo) to compute
>>> what depth counter "should" be i.e
>>
>> I don't think we can take the assumption quicksort will always iterate the
>> input array as a binary tree, so we can't really use log2(hi - lo) as the
>> depth.
>> However I found another issue, which a change to use 'stack_node *top = stack + 1'
>> removed the initial depth initialization.
>>
>
> Bah, sorry missed your reply!
>
> What i mean is not that you take `log2(hi - lo)` as the current depth,
> but whenever
> you pop a new set of bounds you set `max_depth = log2(hi - lo)`. Then
> in the loop if
> `--max_depth <= 0` we are in pessimistic state -> do introsort.
>
> Think this works to avoid the O(N^2) state and saves a lot of state tracking.
The problem with this approach is not all architectures has a clz instruction,
so log2 will be a libgcc function call. And due how partition are pushed on
the stack, it is guarantee that no more than log (total_elems) state track is
required (O(1)). To track the state it is an extra of STACK_SIZE size_t state,
so I am not sure if this will be a large improvement even archs with clz
instructions.
On Thu, Aug 31, 2023 at 7:13 AM Adhemerval Zanella Netto
<adhemerval.zanella@linaro.org> wrote:
>
>
>
> On 29/08/23 04:45, Noah Goldstein wrote:
> > On Tue, Jul 18, 2023 at 2:08 PM Adhemerval Zanella Netto
> > <adhemerval.zanella@linaro.org> wrote:
> >>
> >>
> >>
> >> On 18/07/23 14:37, Noah Goldstein wrote:
> >>> On Tue, Jul 18, 2023 at 9:19 AM Adhemerval Zanella via Libc-alpha
> >>> <libc-alpha@sourceware.org> wrote:
> >>>>
> >>>> This patch makes the quicksort implementation to acts as introsort, to
> >>>> avoid worse-case performance (and thus making it O(nlog n)). It switch
> >>>> to heapsort when the depth level reaches 2*log2(total elements). The
> >>>> heapsort is a textbook implementation.
> >>>>
> >>>> Checked on x86_64-linux-gnu.
> >>>> ---
> >>>> stdlib/qsort.c | 85 ++++++++++++++++++++++++++++++++++++++++++++++----
> >>>> 1 file changed, 79 insertions(+), 6 deletions(-)
> >>>>
> >>>> diff --git a/stdlib/qsort.c b/stdlib/qsort.c
> >>>> index 640896a598..054c900b02 100644
> >>>> --- a/stdlib/qsort.c
> >>>> +++ b/stdlib/qsort.c
> >>>> @@ -119,6 +119,7 @@ typedef struct
> >>>> {
> >>>> char *lo;
> >>>> char *hi;
> >>>> + size_t depth;
> >>>
> >>> You don't actually need to track depth in the "stack"
> >>> when you pop hi/lo you can take log2(hi - lo) to compute
> >>> what depth counter "should" be i.e
> >>
> >> I don't think we can take the assumption quicksort will always iterate the
> >> input array as a binary tree, so we can't really use log2(hi - lo) as the
> >> depth.
> >> However I found another issue, which a change to use 'stack_node *top = stack + 1'
> >> removed the initial depth initialization.
> >>
> >
> > Bah, sorry missed your reply!
> >
> > What i mean is not that you take `log2(hi - lo)` as the current depth,
> > but whenever
> > you pop a new set of bounds you set `max_depth = log2(hi - lo)`. Then
> > in the loop if
> > `--max_depth <= 0` we are in pessimistic state -> do introsort.
> >
> > Think this works to avoid the O(N^2) state and saves a lot of state tracking.
>
> The problem with this approach is not all architectures has a clz instruction,
> so log2 will be a libgcc function call. And due how partition are pushed on
> the stack, it is guarantee that no more than log (total_elems) state track is
> required (O(1)). To track the state it is an extra of STACK_SIZE size_t state,
> so I am not sure if this will be a large improvement even archs with clz
> instructions.
absolutely not necessary, just a thought.
@@ -119,6 +119,7 @@ typedef struct
{
char *lo;
char *hi;
+ size_t depth;
} stack_node;
/* The stack needs log (total_elements) entries (we could even subtract
@@ -128,22 +129,83 @@ typedef struct
enum { STACK_SIZE = CHAR_BIT * sizeof (size_t) };
static inline stack_node *
-push (stack_node *top, char *lo, char *hi)
+push (stack_node *top, char *lo, char *hi, size_t depth)
{
top->lo = lo;
top->hi = hi;
+ top->depth = depth;
return ++top;
}
static inline stack_node *
-pop (stack_node *top, char **lo, char **hi)
+pop (stack_node *top, char **lo, char **hi, size_t *depth)
{
--top;
*lo = top->lo;
*hi = top->hi;
+ *depth = top->depth;
return top;
}
+/* A fast, small, non-recursive O(nlog n) heapsort, adapted from Linux
+ lib/sort.c. Used on introsort implementation as a fallback routine with
+ worst-case performance of O(nlog n) and worst-case space complexity of
+ O(1). */
+
+static inline void
+siftdown (void *base, size_t size, size_t k, size_t n,
+ enum swap_type_t swap_type, __compar_d_fn_t cmp, void *arg)
+{
+ while (k <= n / 2)
+ {
+ size_t j = 2 * k;
+ if (j < n && cmp (base + (j * size), base + ((j + 1) * size), arg) < 0)
+ j++;
+
+ if (cmp (base + (k * size), base + (j * size), arg) >= 0)
+ break;
+
+ do_swap (base + (size * j), base + (k * size), size, swap_type);
+ k = j;
+ }
+}
+
+static inline void
+heapify (void *base, size_t size, size_t n, enum swap_type_t swap_type,
+ __compar_d_fn_t cmp, void *arg)
+{
+ size_t k = n / 2;
+ while (1)
+ {
+ siftdown (base, size, k, n, swap_type, cmp, arg);
+ if (k-- == 0)
+ break;
+ }
+}
+
+static void
+heapsort_r (void *base, void *end, size_t size, enum swap_type_t swap_type,
+ __compar_d_fn_t cmp, void *arg)
+{
+ const size_t count = ((uintptr_t) end - (uintptr_t) base) / size;
+
+ if (count < 2)
+ return;
+
+ size_t n = count - 1;
+
+ /* Build the binary heap, largest value at the base[0]. */
+ heapify (base, size, n, swap_type, cmp, arg);
+
+ /* On each iteration base[0:n] is the binary heap, while base[n:count]
+ is sorted. */
+ while (n > 0)
+ {
+ do_swap (base, base + (n * size), size, swap_type);
+ n--;
+ siftdown (base, size, 0, n, swap_type, cmp, arg);
+ }
+}
/* Order size using quicksort. This implementation incorporates
four optimizations discussed in Sedgewick:
@@ -229,7 +291,7 @@ _quicksort (void *const pbase, size_t total_elems, size_t size,
const size_t max_thresh = MAX_THRESH * size;
- if (total_elems == 0)
+ if (total_elems <= 1)
/* Avoid lossage with unsigned arithmetic below. */
return;
@@ -241,6 +303,10 @@ _quicksort (void *const pbase, size_t total_elems, size_t size,
else
swap_type = SWAP_BYTES;
+ /* Maximum depth before quicksort switches to heapsort. */
+ size_t depth = 2 * (sizeof (size_t) * CHAR_BIT - 1
+ - __builtin_clzl (total_elems));
+
if (total_elems > MAX_THRESH)
{
char *lo = base_ptr;
@@ -250,6 +316,13 @@ _quicksort (void *const pbase, size_t total_elems, size_t size,
while (stack < top)
{
+ if (depth == 0)
+ {
+ heapsort_r (lo, hi, size, swap_type, cmp, arg);
+ top = pop (top, &lo, &hi, &depth);
+ continue;
+ }
+
char *left_ptr;
char *right_ptr;
@@ -313,7 +386,7 @@ _quicksort (void *const pbase, size_t total_elems, size_t size,
{
if ((size_t) (hi - left_ptr) <= max_thresh)
/* Ignore both small partitions. */
- top = pop (top, &lo, &hi);
+ top = pop (top, &lo, &hi, &depth);
else
/* Ignore small left partition. */
lo = left_ptr;
@@ -324,13 +397,13 @@ _quicksort (void *const pbase, size_t total_elems, size_t size,
else if ((right_ptr - lo) > (hi - left_ptr))
{
/* Push larger left partition indices. */
- top = push (top, lo, right_ptr);
+ top = push (top, lo, right_ptr, depth - 1);
lo = left_ptr;
}
else
{
/* Push larger right partition indices. */
- top = push (top, left_ptr, hi);
+ top = push (top, left_ptr, hi, depth - 1);
hi = right_ptr;
}
}