www.digitalmars.com         C & C++   DMDScript  

digitalmars.D - The Quick Mom Algorithm

reply Andrei Alexandrescu <SeeWebsiteForEmail erdani.org> writes:
http://dpaste.dzfl.pl/05a82699acc8

So over the past few days I've been in the zone working on a smooth 
implementation of "Median of Medians" 
(https://en.wikipedia.org/wiki/Median_of_medians). Its performance is 
much better compared to the straightforward implementation. However, in 
practice it's very hard to get it to beat simple heuristics (such as 
median of five, random pivot etc). These heuristics run `n` comparisons 
for an input of length `n`, whereas the MoM needs over twice as many 
comparisons. Also it does a bunch of swapping around the data. Overall I 
got it within 2.5x of the heuristic-based topN for most data sizes up to 
tens of millions.

While thinking of MoM and the core reasons of its being slow (adds nice 
structure to its input and then "forgets" most of it when recursing), I 
stumbled upon a different algorithm. It's much simpler, also 
deterministic and faster than MoM for many (most?) inputs. But it's not 
guaranteed to be linear. After having pounded at this for many hours, it 
is clear that I am in need of some serious due destruction. I call it a 
"quick median of medians" or in short "quick mom".

Consider the algorithm defined as follows over a range `r` of length 
`n`. It returns an index to an element likely to be in the second 
tertile of r. (A tertile is a third of the range. So the expectation is 
that the returned index x is such that e<=r[x] for at least one third of 
the range, and r[x]<=e again for at least one third of the range.)

0. If n<=3, compute median by rote and return its index.

1. Divide r in three equal adjacent subranges r0 = r[0 .. $/3], r1 = 
r[$/3 .. $*2/3], r2 = r[$*2/3 .. $].

2. Recurse to get the median of these three, call them m0, m1, m2.

3. Return the median of m0, m1, m2.

Note that no data has been written yet; we only have an estimate of a 
pivot. On random inputs it's expected that the median thus gotten is 
greater than 1/6 + 1/6 = 1/3 elements, and less than the same fraction.

The algorithm completes in linear time.

However, the pivot obtained is just an approximation. In the worst case 
it's possible that e.g. all recursions return the leftmost of the 
allowed index, so the bound deteriorates by 1/3 for each recursion 
depth, or generally to (1/3)^^log3(n).

Not good! That's where the quick mom pays attention. After it partitions 
the data using the pivot obtained by the quick mom method, the 
partitioning stage checks whether the resulting pivot falls within 
bounds. If not, it runs a precise selection method (such as proper 
median of medians - "thorough mom") to bring the pivot where it needs to 
be. The data patterns that cause the quick mom to fail systematically 
are rather odd but worst case is what it is.

Overall the partitioning either succeeds with the quick mom pivot in 
linear time, or does one more linear pass if "unlucky". So overall 
partition is linear. (Micro-optimization: not all data needs to be 
repartitioned, only that outside the not-so-good pivot.) Unlike the 
quick mom, the partition guarantees a pivot in the mid tertile.

After partitioning the classic quickselect algorithm may be implemented.

There's one more _really_ juicy detail. In fact, the quick mom may 
finish without looking at all data. Look at the implementation - there 
are two overloads of quickMom. The second gets the bounds and works as 
follows: if you've already computed two elements of a median, the third 
may be chosen only if in between the two. Otherwise, you only care 
whether it's larger or smaller than the other two. This allows the 
algorithm to finish certain recursion branches early.


Destroy!

Andrei
Jan 29
parent reply Ivan Kazmenko <gassa mail.ru> writes:
On Friday, 29 January 2016 at 22:47:44 UTC, Andrei Alexandrescu 
wrote:
 http://dpaste.dzfl.pl/05a82699acc8

 While thinking of MoM and the core reasons of its being slow 
 (adds nice structure to its input and then "forgets" most of it 
 when recursing), I stumbled upon a different algorithm. It's 
 much simpler, also deterministic and faster than MoM for many 
 (most?) inputs. But it's not guaranteed to be linear. After 
 having pounded at this for many hours, it is clear that I am in 
 need of some serious due destruction. I call it a "quick median 
 of medians" or in short "quick mom".
The code does not compile when I try to use topN, as medianOfUpTo is not present. Please include it. Anyway, for now, I just wanted to quickly check whether my previous attack on topN would have any effect here. I've replaced the line if (r.length <= 5) return medianOfUpTo!(5, lp)(r); with a quick fix: if (r.length == 2) { if (!less (r.front, r.back)) swap (r.front, r.back); return 0; } if (r.length == 1) return 0; Here is the code (just glued together the two snippets): http://dpaste.dzfl.pl/f648a600a11d Turns out like maybe it does have its effect (size = 5 million): building Element array: 2442 msecs checking Element array: 2432 msecs checking int array: 388 msecs checking random int array: 140 msecs checking increasing int array: 41 msecs checking decreasing int array: 42 msecs The int array built by the attack makes the program run slower than a random array, and the ratio increases with size.
Feb 02
parent Andrei Alexandrescu <SeeWebsiteForEmail erdani.org> writes:
On 02/02/2016 05:14 PM, Ivan Kazmenko wrote:
 On Friday, 29 January 2016 at 22:47:44 UTC, Andrei Alexandrescu wrote:
 http://dpaste.dzfl.pl/05a82699acc8

 While thinking of MoM and the core reasons of its being slow (adds
 nice structure to its input and then "forgets" most of it when
 recursing), I stumbled upon a different algorithm. It's much simpler,
 also deterministic and faster than MoM for many (most?) inputs. But
 it's not guaranteed to be linear. After having pounded at this for
 many hours, it is clear that I am in need of some serious due
 destruction. I call it a "quick median of medians" or in short "quick
 mom".
The code does not compile when I try to use topN, as medianOfUpTo is not present. Please include it.
Thanks for looking into this. I've been working a lot more on this and have an algorithm in testing that would be very interesting to develop an attack input for. Please give me to the end of the week and I'll send it along. Thx! -- Andrei
Feb 03