by

Powers of two: the solution

(This post solves the puzzle given in the previous post.)

Doing an empirical analysis

As we are interested in powers of two starting with a given set of digits, without caring about their overall magnitude, let’s plot a histogram of the values of

\displaystyle \frac{2^n}{10^m},

with n going from 1 to 100000 and m chosen to get values between 0.1 and 1.

By using a simple python script and matplotlib, we get the following distribution:

Distribution of the first digits (considered as a fraction between 0.10 and 0.99) of the first 100000 powers of two.

The distribution follows Benford’s law quite precisely

The previous plot, adding the prediction given by Benford's law.

making us expect a uniform distribution in log space

Distribution of the first digits (considered as a fraction between 0.10 and 0.99) of the first 100000 powers of two but grouped in logarithmic bins.

As the leading digits expressed as fractions cover the interval between 0.1 and 1 quite densely, we can expect that some power of two will start with the digits 31415926535897932384626433832795028841971, but this plot is quite far from being a proof.

Translating the problem

We are looking for a power of two 2^n starting with a given sequence of m digits, k. By expressing the digits sequence as a fraction between 0.1 and 1, k/10^m, we can represent the desired condition by the following inequality:

\displaystyle \frac{k}{10^m} \cdot 10^d \le 2^n < \frac{k+1}{10^m} \cdot 10^d.

As logarithms are monotonically increasing, we can apply them on both sides of the inequality:

\displaystyle \log_{10}\frac{k}{10^m} + d \le n \log_{10} 2 < \log_{10}\frac{k+1}{10^m} + d.

If we assume k+1 is not a power of ten, we can translate the inequality to the fractional parts of each member:

\displaystyle \left\{\log_{10}\frac{k}{10^m} + d\right\} \le \left\{n \log_{10} 2\right\} < \left\{\log_{10}\frac{k+1}{10^m} + d\right\}

\displaystyle \log_{10}\frac{k}{10^m} \le \left\{n \log_{10} 2\right\} < \log_{10}\frac{k+1}{10^m}.

Now the problem is reduced to see if we can find a multiple of \log_{10} 2 such that its fractional part falls between two given real numbers. In the next section we will see how this problem can be solved.

Solving the problem

By the Archimedean property of the real numbers, we know that there must be a natural number q such that

\displaystyle \frac{1}{q} < \log_{10}\frac{k+1}{10^m} - \log_{10}\frac{k}{10^m}.

Now, let’s analyze the distribution of \{n\log_{10}2\} for n in the range 1, …, q+1. As we have q+1 numbers inside the unit interval, we must have two numbers with separation less than or equal to 1/q (it’s easy to prove by contradiction). Then there are two exponents, r and s, such that

\displaystyle 0 \le \left\{r\log_{10}2\right\} - \left\{s\log_{10}2\right\} \le \frac{1}{q},

\displaystyle 0 < \left\{(r-s)\log_{10}2\right\} \le \frac{1}{q},

where we know that the inequality with zero is strict because \log_{10}2 is irrational. Then we can use the Archimedean property again to extend the inequality:

\displaystyle \frac{1}{t} < \left\{(r-s)\log_{10}2\right\} \le \frac{1}{q}.

It's intuitively easy to see (and not very difficult to prove) that we must have a number u such that

\displaystyle u\left\{(r-s)\log_{10}2\right\} < \log_{10}\frac{k}{10^m} < (u+1)\left\{(r-s)\log_{10}2\right\} < \log_{10}\frac{k+1}{10^m}

As all the involved numbers are in the range [0.1, 1), we can put the multipliers inside the fractional parts and we can also discard the first inequality because it's not important for our proof:

\displaystyle \log_{10}\frac{k}{10^m} < \left\{(u+1)(r-s)\log_{10}2\right\} < \log_{10}\frac{k+1}{10^m}.

Removing the fractional parts:

\displaystyle \log_{10}\frac{k}{10^m} + \lfloor(u+1)(r-s)\log_{10}2\rfloor < (u+1)(r-s)\log_{10}2

\displaystyle (u+1)(r-s)\log_{10}2 < \log_{10}\frac{k+1}{10^m} + \lfloor(u+1)(r-s)\log_{10}2\rfloor

Exponentiating:

\displaystyle \frac{k}{10^m} \cdot 10^{\lfloor(u+1)(r-s)\log_{10}2\rfloor} < 2^{(u+1)(r-s)}

\displaystyle 2^{(u+1)(r-s)} < \frac{k+1}{10^m}\cdot 10^{\lfloor(u+1)(r-s)\log_{10}2\rfloor}

Conclusion

We have found that we can get a power of two starting with the sequence of digits we want, though we haven’t found any bounds on the required exponents. In fact, the powers of two follow Benford’s law and we can use that to estimate the magnitude of the required exponent:

\displaystyle P({\rm power\ of\ two\ starts\ with\ required\ digits}) \approx \log_{10}\left(1 + \frac{1}{3.2 \cdot 10^{40}}\right)

\displaystyle P({\rm power\ of\ two\ starts\ with\ required\ digits}) \approx \frac{1 + \frac{1}{3.2 \cdot 10^{40}} - 1}{\ln 10}

\displaystyle P({\rm power\ of\ two\ starts\ with\ required\ digits}) \approx 10^{-41}

So we can be quite confident that an exponent smaller than 10^{43} should meet the requirements, though it’s obviously impossible to find it by brute force. But it should be feasible to compute values of u, r and s and getting the first 100 digits of 2^{(u+1)(r-s)}.

by

Powers of two: a puzzle

While the answer to the previous post remains in its editing phase, this is a quick puzzle that cannot be brute forced:

Is there any power of two whose decimal representation starts with 31415926535897932384626433832795028841971?

My answer will be included in the next “quick post”.

by

Do we need to live with the uncertainty?

Introduction

In the previous post we have seen that the halting problem cannot be solved. But maybe the problem is just in our programming languages, so let’s see if we can improve them to disallow non-halting programs.

Bounded C

A program can fail to halt by two reasons (we are ignoring I/O, limited memory sizes and all that earthly issues :D ): it can get trapped in an infinite loop or we can enter an infinite recursion. So let’s take a language like C and remove these “troublesome” features from it by

  • only allowing iteration via the for (b = 0; b < a; b++) construct and
  • removing the possibility of doing either direct or indirect recursion.

At first sight we haven’t lost a lot. For example we can easily raise numbers to a power

int my_pow(int base, int exp)
{
    int p = 1, i;
    for (i = 0; i < exp; i++)
        p *= base;
    return p;
}

or compute the elements of the Fibonacci series.

int fib(int n)
{
    int a = 0, b = 1, i, t;
    for (i = 0; i < n; i++)
    {
        t = a + b;
        a = b;
        b = t;
    }
    return a;
}

But maybe we have lost something subtler…

Challenge

Can we make a function in “normal C” that cannot be implemented in our “bounded C”? The answer will be given in the following post… or in the comments :D

by

A really impossible problem

The problem

Previously we have solved a puzzle that claimed to be impossible but was, in fact, just a bit hard :-D But in this post we will see a problem that is really impossible to solve, the halting problem:

Given a description of an arbitrary computer program, decide whether the program finishes running or continues to run forever.

For some programs it’s easy to see if they halt or not:

def program_1():
    return 2 + 2

def program_2():
    while True: pass

def program_3():
    from itertools import count
    for i in count():
        if str(2**i)[0] == '7':
            break

def program_4():
    from itertools import count
    for i in count():
        if str(2**i)[-1] == '7':
            break

For others it’s not so simple:

def program_5():
    def z(m, n):
        if m == 0:
            return n + 1
        elif n == 0:
            return z(m - 1, 1)
        else:
            return z(m - 1, z(m, n - 1))
    z(4, 4)

And there are simple ones where it’s an open problem to know if they halt or not:

def program_6():
    from itertools import count
    for i in count(1):
        p = 2 ** i
        r = int(''.join(reversed(str(p))))
        for j in count(1):
            q = 5 ** j
            if q > r:
                break
            elif q == r:
                return

The proof

Let’s assume we have a procedure, represented by the function does_halt(), to determine if any given program will halt:

def does_halt(program):
    #
    # magic implementation
    #

We don’t care about the details of the does_halt() function. The only requirements are:

  • Its implementation should be finite.
  • It must accept any Python function that doesn’t take arguments as a program.
  • It must give an answer in an easily bounded time. For example, we should be able to say “if the input program has 1000 bytes, the function should return a boolean value indicating if it halts in no more than NNN seconds”.

The clever point of the proof, though it’s not a new kind of cleverness, is to build a program that “halts only if it doesn’t”:

def paradox():
    if does_halt(paradox):
        while True:
            pass

As the paradox() program is finite, including this small amount of code plus the code of does_halt(), does_halt() should return in a finite amount of time. If does_halt() returns False, indicating that paradox() does not halt, the program exits, finishing its execution in a finite time. On the other hand, if does_halt() returns True, indicating that paradox() does halt, it enters an infinite loop.

This proves that some of these propositions must be true:

  • does_halt() must not be finitely implementable.
  • does_halt() must reject some valid programs.
  • does_halt() execution time is unbounded.
by

Quines [remix]

As the quines in the previous post were criticized as boring and ordinary :-P , I did a remix:

#include <stdio.h>
#define s char s[]
s="#if 0\nimport json;r=json.dumps\nprint'#include <stdio.h>\\n#define s char s[]\\ns=%s;\\n%s'%(r(s),s)\n\"\"\" \"\n#elif 1\n#undef s\nint main(void)\n{\n  char *t = s;\n  printf(\"#include <stdio.h>\\n#define s char s[]\\ns=\\\"\");\n  while (*t)\n  {\n    if (*t == '\\n')\n      printf(\"\\\\n\");\n    else if (*t == '\"')\n      printf(\"\\\\\\\"\");\n    else if (*t == '\\\\')\n      printf(\"\\\\\\\\\");\n    else\n      printf(\"%c\", *t);\n    t++;\n  }\n  printf(\"\\\";\\n%s\\n\", s);\n  return 0;\n}\n#elif 0\n\" \"\"\"\n#endif";
#if 0
import json;r=json.dumps
print'#include <stdio.h>\n#define s char s[]\ns=%s;\n%s'%(r(s),s)
""" "
#elif 1
#undef s
int main(void)
{
  char *t = s;
  printf("#include <stdio.h>\n#define s char s[]\ns=\"");
  while (*t)
  {
    if (*t == '\n')
      printf("\\n");
    else if (*t == '"')
      printf("\\\"");
    else if (*t == '\\')
      printf("\\\\");
    else
      printf("%c", *t);
    t++;
  }
  printf("\";\n%s\n", s);
  return 0;
}
#elif 0
" """
#endif

You can check that it works on the web or by downloading the file and testing it:

$ python polyquine.c | diff polyquine.c -
$ gcc -ansi -pedantic -Wall polyquine.c -o polyquine && ./polyquine | diff polyquine.c -

I’m aware that there are some impressive examples out there, but I haven’t analyzed them to avoid spoiling the fun.

What other language should I add? Reader contributions are welcome!

by

Quines

Introduction

One interesting programming exercise is to write a quine, a program that outputs its own source code (excluding empty programs!). If we naively try to write it just by using the print statement we will get into an infinite regression:

print 'print \' print \\\'print \\\\\\\'...

The basic problem is that the instructions required to write other instructions lead to an infinite recursion. So how can we avoid this problem?

The trick

One general way to do it is the one used by Von Neumann to make his abstract self-replicating machines and by nature itself. The trick is to handle the same instructions in two different ways: as instructions to be followed and as data to be copied.

Writing a C quine

As we will need to access the same data in two different ways, we must define a variable to avoid duplicating it:

char s[] = "<data will go here>";

It’s easy to write the code that prints everything before the data,

int main(void)
{
  printf("char s[] = \"");

but now we need to print the data as represented inside the string. Fortunately, of the characters we use, only the newline, the quote and the backlash require escape sequences:

  char *t = s;
  while (*t)
  {
    if (*t == '\n')
      printf("\\n");
    else if (*t == '"')
      printf("\\\"");
    else if (*t == '\\')
      printf("\\\\");
    else
      printf("%c", *t);
    t++;
  }

Finally, we print the data as a normal string and exit main():

  printf("%s", s);
  return 0;
}

This program is finished with exception of filling in the data in variable s. As this is quite tedious to do by hand, we will do it using a program:

char s[] = "\";\nint main(void)\n{\n  printf(\"char s[] = \\\"\");\n  char *t = s;\n  while (*t)\n  {\n    if (*t == '\\n')\n      printf(\"\\\\n\");\n    else if (*t == '\"')\n      printf(\"\\\\\\\"\");\n    else if (*t == '\\\\')\n      printf(\"\\\\\\\\\");\n    else\n      printf(\"%c\", *t);\n    t++;\n  }\n  printf(\"%s\", s);\n  return 0;\n}\n";

Putting it all together, we get a program that prints its own source:

char s[] = "\";\nint main(void)\n{\n  printf(\"char s[] = \\\"\");\n  char *t = s;\n  while (*t)\n  {\n    if (*t == '\\n')\n      printf(\"\\\\n\");\n    else if (*t == '\"')\n      printf(\"\\\\\\\"\");\n    else if (*t == '\\\\')\n      printf(\"\\\\\\\\\");\n    else\n      printf(\"%c\", *t);\n    t++;\n  }\n  printf(\"%s\", s);\n  return 0;\n}\n";
int main(void)
{
  printf("char s[] = \"");
  char *t = s;
  while (*t)
  {
    if (*t == '\n')
      printf("\\n");
    else if (*t == '"')
      printf("\\\"");
    else if (*t == '\\')
      printf("\\\\");
    else
      printf("%c", *t);
    t++;
  }
  printf("%s", s);
  return 0;
}

But, though it works in GCC, it’s not a strict C99 file:

cc1: warnings being treated as errors
prog.c: In function ‘main’:
prog.c:4: error: implicit declaration of function ‘printf’
prog.c:4: error: incompatible implicit declaration of built-in function ‘printf’

This can be fixed by including the stdio.h header:

#include <stdio.h>
char s[] = "\";\nint main(void)\n{\n  printf(\"#include <stdio.h>\\nchar s[] = \\\"\");\n  char *t = s;\n  while (*t)\n  {\n    if (*t == '\\n')\n      printf(\"\\\\n\");\n    else if (*t == '\"')\n      printf(\"\\\\\\\"\");\n    else if (*t == '\\\\')\n      printf(\"\\\\\\\\\");\n    else\n      printf(\"%c\", *t);\n    t++;\n  }\n  printf(\"%s\", s);\n  return 0;\n}\n";
int main(void)
{
  printf("#include <stdio.h>\nchar s[] = \"");
  char *t = s;
  while (*t)
  {
    if (*t == '\n')
      printf("\\n");
    else if (*t == '"')
      printf("\\\"");
    else if (*t == '\\')
      printf("\\\\");
    else
      printf("%c", *t);
    t++;
  }
  printf("%s", s);
  return 0;
}

Doing a Python quine

Doing a Python quine is much easier, because the built-in function repr() gives a string representation of an object, allowing us to skip most of the code of the C quine. Using it we can get this natural three line quine:

s = "print 's = %s' % repr(s)\nprint s"
print 's = %s' % repr(s)
print s

If we want to do code golfing, we can:

  • Remove unnecesary spaces and line breaks.
  • Printing the representation of the data and the data using only one print statement.

This give us a much shorter quine:

s="print's=%s;%s'%(repr(s),s)";print's=%s;%s'%(repr(s),s)

though still longer than the record one.

by

The fall of the slinky I

The video

The video shows how a spring suspended from one of its ends reacts when its dropped. It can be observed that the lower end “doesn’t know what happened” until a wave propagates to it. In this post we will make a computer simulation of its behavior, to see if we can reproduce the phenomenon, and in the next we will apply a more analytic approach to the same problem.

Discretization

As we are studying the actions of gravity and inertia, we need to model the slinky as a massive spring. In (macroscopic) reality the mass and elasticity are continuously distributed throughout the slinky (macroscopically speaking) but, for the purpose of simulating its behavior with a computer model, we will represent the object as a series of masses connected with massless ideal springs:

Slinky discretization, showing how the discretized element properties relate to the original ones.

If we divide a slinky of total mass M in N small masses, each of them will have the value m = M/N. As the spring constants of series springs add as the resistances of parallel resistors, if we have a slinky with overall spring constant K divided into N-1 smaller springs, each of them will have a bigger spring constant, k = K(N-1) and (obviously) a smaller unstressed length, l_0 = L/(N-1).

Writing the simulation

Now it’s just a question of applying Newton’s laws of motion and the ideal spring equation to get a system of ordinary differential equations:

\displaystyle \frac{d^2x_i}{dt^2} = g + \frac{k}{m}\left[(x_{i+1}-x_i-l_0)H[i]-(x_i-x_{i-1}-l_0)H[N-1-i]\right]

where x_i are the coordinates of the masses (with i going from 0 to N-1), g is the acceleration of gravity, k is the spring constant of each massless spring, m is the value of each small mass, l_0 is the unstressed length of each massless spring and H[n] is the discrete Heaviside step function (used to avoid depending on undefined values).

This second order system can be reduced to system of first order ODEs in the usual way. Integrating it using RK4, we get the following Python code:

def get_xdot(x):
    sk = K * (NX - 1)
    sl = L / (NX - 1)
    mm = M / NX
    xdot = [x[NX + i] if i < NX else G
            for i in range(2 * NX)]
    for i in range(NX - 1):
        a = sk * (x[i + 1] - x[i] - sl) / mm
        xdot[NX + i] += a
        xdot[NX + i + 1] -= a
    return xdot

def rk4_step(x, dt):
    k1 = get_xdot(x)
    k2 = get_xdot([x[i] + dt * k1[i] / 2.0 for i in range(len(x))])
    k3 = get_xdot([x[i] + dt * k2[i] / 2.0 for i in range(len(x))])
    k4 = get_xdot([x[i] + dt * k3[i] for i in range(len(x))])
    return [x[i] + dt * (k1[i] + 2.0 * (k2[i] + k3[i]) + k4[i]) / 6.0
            for i in range(len(x))]

Now we need to define the initial conditions. If we just start with the masses separated by a distance l_0 and at rest, we won’t match the initial conditions of the slinky, because it was being stretched by the action of gravity. It’s not very difficult to compute initial conditions that leave the system at rest if it’s being held:

def initial_x():
    sl0 = L / (NX - 1)
    mm = M / NX
    sk = K * (NX - 1)
    w = M - mm
    x = [0.0 for i in range(2 * NX)]
    for i in range(1, NX):
        x[i] = x[i - 1] + sl0 + w / sk
        w -= mm
    return x

The remaining code is just plumbing and matplotlib presentation code. The whole program can be seen at the repository.

Running the simulation

If we run the simulation with the parameters

NT = 1000 # number of timesteps
NX = 40   # number of masses
T = 1.0   # simulation duration

L = 0.5   # slinky length
K = 1.0   # slinky overall spring constant
M = 1.0   # slinky mass
G = 1.0   # gravitational acceleration

we get the following results:

A simulation where the springs are too soft, giving some negative spring lengths (and consequent overlapping) near t = 1.

In this plot the gray lines represent the trajectory of the small masses and the black lines the trajectory of the slinky’s ends.

Clearly the springs are too soft and we are getting unphysical results, as we spring lengths go negative when t nears 1. To fix that, let’s run the simulations with a greater spring constant, K = 5:

A simulation where the springs have a physically reasonable constant, giving an intriguing behavior to its ends.

Now we get a more reasonable result, showing a phenomenon that is more similar to the one observed in the video: the bottom of the slinky remains in place while the top begins to fall. Now we can check if the slinky remains stationary when held by m_0:

def get_xdot(x):
    sk = K * (NX - 1)
    sl = L / (NX - 1)
    mm = M / NX
    xdot = [x[NX + i] if i < NX else G
            for i in range(2 * NX)]
    for i in range(NX - 1):
        a = sk * (x[i + 1] - x[i] - sl) / mm
        xdot[NX + i] += a
        xdot[NX + i + 1] -= a
    xdot[NX] = 0.0 # holding the slinky
    return xdot

Simulation to check if the slinky remains stationary when held from its upper end.

Conclusions

This simulation validates the main point of the original video: the lower end “doesn’t know” the upper end was released until a compression wave reaches it, at t \approx 0.45 in our simulation. But the detailed behavior differs, as the slinky only shows a compression wave once it reaches the nonlinear regime (when is no more space between the spires).

In the next post we will show an analysis of this nonlinear behavior and the analytical solution to the idealized slinky drop that was numerically simulated in this post.

by

Top k – Solution

[Answer to the problem in Top k.]

This was the problem:

Given a list of n numbers, get the top k ones in O(n) time. (Partial credit for o(n log n) solutions.)

n log k solution

The easiest way to solve this problem in o(n log n) time is to scan the n elements, keeping track of the greatest k elements seen up to the current element. Using a heap to store the top k elements, we get the following algorithm:

#define SWAP(type, a, b) {type swap_temp = a; a = b; b = swap_temp;}

static void minheap_pushdown(int *h, size_t hs, size_t i)
{
    size_t j = 0;
    if (2 * i + 2 < hs)
        j = (h[2 * i + 1] < h[2 * i + 2]) ? 2 * i + 1 : 2 * i + 2;
    else if (2 * i + 1 < hs)
        j = 2 * i + 1;
    if (j != 0 && h[j] < h[i])
    {
        SWAP(int, h[i], h[j]);
        minheap_pushdown(h, hs, j);
    }
}

static void minheap_raise(int *h, size_t i)
{
    if (i == 0)
        return;
    if (h[i] < h[(i - 1) / 2])
    {
        SWAP(int, h[i], h[(i - 1) / 2]);
        minheap_raise(h, (i - 1) / 2);
    }
}

void top_k_nlogk_cp(int *top_k, size_t k, const int *l, size_t n)
{
    size_t i;

    for (i = 0; i < k; i++)
    {
        top_k[i] = l[i];
        minheap_raise(top_k, i);
    }

    for (i = k; i < n; i++)
    {
        if (l[i] > top_k[0])
        {
            top_k[0] = l[i];
            minheap_pushdown(top_k, k, 0);
        }
    }
}

As minheap_pushdown() and minheap_raise() have O(log k) cost, the total cost of this solution will be O(n log k). This is better than O(n log n) if k grows slower than any power of n (if k is constant or O(log n), for example).

The previous solution is enough to get partial credit, but we can do better :-D

Linear time solution

Quicksort is based on a linear time partition operation. This operation starts with an array and an element of this array, called pivot, and ends with a partially sorted array, having the pivot in his sorted position, smaller elements before it and bigger ones after it. Then doing a single partition operation would suffice… if we could guarantee that the chosen pivot is the (n-k+1)-th element of the sorted array :-D

We cannot do that, but we can do something that is asymptotically equivalent: doing a quicksort-style recursion, though only in the partition that contains the desired pivot. The result is this (converting the tail recursion to iteration):

int *top_k_n_ip(int *l, size_t n, size_t k)
{
    size_t lo, hi, pos, i, j;
    int pivot;
   
    lo = 0;
    hi = n;
    pos = n - k;

    while (hi - lo > 1)
    {
        i = lo + 1;
        j = hi - 1;
        pivot = l[lo];

        while (1)
        {
            while (i < hi && l[i] <= pivot)
                i++;
            while (j > lo && l[j] > pivot)
                j--;
            if (i > j)
                break;
            SWAP(int, l[i], l[j]);
        }

        SWAP(int, l[lo], l[j]);

        if (j < pos)
            lo = j + 1;
        else if (j > pos)
            hi = j;
        else
            break;
    }

    return l + n - k;
}

Now we need to check that this gives us a linear-time algorithm. The execution time of the body of the while (1) loop is bounded from above by C times hi – lo, where C is a constant. If we assume that the chosen pivot divides the l [ lo .. hi ] array section in two roughly equal parts (this assumption can be made rigorous by choosing the pivot carefully), we get

T(n) \leq A + (B + C \cdot n) + (B + C \cdot n/2) + \hdots + (B + C \cdot 1)

where A, B and C are constants.

Summing all the terms containing B and taking out all the C factors we get:

T(n) \leq A + \lceil\log_2 n\rceil B + C (1 + \hdots + n/2 + n)

T(n) \leq A + \lceil\log_2 n\rceil B + 3 C n

T(n) = O(n)

As we know that any correct algorithm must employ \Omega(n) time (it needs to inspect all elements), we also know that this algorithm has \Theta(n) time complexity. Now let’s check that against reality :-D

Benchmarks

We will benchmark the following algorithm implementations:

Execution times for all the algorithms when k is fixed at 30.

As k has a fixed value in the previous plot, the O(n log k) algorithm is effectively an O(n) algorithm. In fact, it performs better than the truly linear algorithms, probably due to cache effects.

Another thing that can be seen in these plots is the difficulty of seeing the difference between n and n log n in a log-log plot, but that is something to be expected.

Execution times for all the algorithms when k is 30% of n.

Here the qualitative behavior of the O(n log k) solution is the expected one. The C implementation of the O(n) algorithm looks surprisingly good, with a performance competitive to the C++ implementation included as part of the standard library.

by

Top k

Given a list of n numbers, get the top k ones in O(n) time.

This problem is too easy to solve in O(n log n) time

def top_k_n_log_n(l, k):
    sl = sorted(l)
    return sl[-k:]

but o(n log n) solutions can be presented for partial credit. :-P