Compute fast log base 2 ceiling in python

Question

for given x < 10^15, quickly and accurately determine the maximum integer p such that 2^p <= x

Here are some things I've tried:

First I tried this but it's not accurate for large numbers:

>>> from math import log
>>> x = 2**3
>>> x
8
>>> p = int(log(x, 2))
>>> 2**p == x
True
>>> x = 2**50
>>> p = int(log(x, 2))
>>> 2**p == x #not accurate for large numbers?
False

I could try something like:

p = 1
i = 1
while True:
    if i * 2 > n:
        break
    i *= 2
    p += 1
    not_p = n - p

Which would take up to 50 operations if p was 50

I could pre-compute all the powers of 2 up until 2^50, and use binary search to find p. This would take around log(50) operations but seems a bit excessive and ugly?

I found this thread for C based solutions: Compute fast log base 2 ceiling

However It seems a bit ugly and I wasn't exactly sure how to convert it to python.

You may want to read docs for [bitwise shift](http://docs.python.org/reference/expressions.html#shifting-operations) and [and/or/xor](http://docs.python.org/reference/expressions.html#binary-bitwise-operations). This can help you translate the C code to Python. — mtrw, Oct 28 '12 at 02:55
Perhaps you meant **floor** instead of ceiling? Because p is the max integer such that `2^p <= x` then `p == floor(log(x,2))`. — Bob Stein, Apr 01 '15 at 16:14
...In the [other question you linked to](http://stackoverflow.com/questions/3272424/compute-fast-log-base-2-ceiling) they really do mean ceiling, as in `q == ceil(log(x,2))` — Bob Stein, Apr 01 '15 at 17:24

score 36 · Accepted Answer · edited Apr 18 '23 at 01:57

In Python >= 2.7, you can use the .bit_length() method of integers:

def brute(x):
    # determine max p such that 2^p <= x
    p = 0
    while 2**p <= x:
        p += 1
    return p-1

def easy(x):
    return x.bit_length() - 1

which gives

>>> brute(0), brute(2**3-1), brute(2**3)
(-1, 2, 3)
>>> easy(0), easy(2**3-1), easy(2**3)
(-1, 2, 3)
>>> brute(2**50-1), brute(2**50), brute(2**50+1)
(49, 50, 50)
>>> easy(2**50-1), easy(2**50), easy(2**50+1)
(49, 50, 50)
>>> 
>>> all(brute(n) == easy(n) for n in range(10**6))
True
>>> nums = (max(2**x+d, 0) for x in range(200) for d in range(-50, 50))
>>> all(brute(n) == easy(n) for n in nums)
True

This computes floor(log2(n)), as evidenced by `easy(7) == 2` and `easy(8) == 3` — qwr, Jan 05 '23 at 14:04

Alex Mykyta · Answer 2 · 2022-02-20T17:52:16.257

Be careful! The accepted answer returns floor(log(n, 2)), NOT ceil(log(n, 2)) like the title of the question implies!

If you came here for a clog2 implementation, do this:

def clog2(x):
    """Ceiling of log2"""
    if x <= 0:
        raise ValueError("domain error")
    return (x-1).bit_length()

And for completeness:

def flog2(x):
    """Floor of log2"""
    if x <= 0:
        raise ValueError("domain error")
    return x.bit_length() - 1

score 5 · Answer 3 · edited May 23 '17 at 12:10

You specify in comments your x is an integer, but for anyone coming here where their x is already a float, then math.frexp() would be pretty fast at extracting log base 2:

log2_slow = int(floor(log(x, 2)))
log2_fast = frexp(x)[1]-1

The C function that frexp() calls just grabs and tweaks the exponent. Some more 'splainin:

The subscript[1] is because frexp() returns a tuple (significand, exponent).
The subtract-1 accounts for the significand being in the range [0.5,1.0). For example 2⁵⁰ is stored as 0.5x2⁵¹.
The floor() is because you specified 2^p <= x, so p == floor(log(x,2)).

(Derived from another answer.)

score 3 · Answer 4 · answered Oct 28 '12 at 02:35

You could try the log2 function from numpy, which appears to work for powers up to 2^62:

>>> 2**np.log2(2**50) == 2**50
True
>>> 2**np.log2(2**62) == 2**62
True

Above that (at least for me) it fails due to the limtiations of numpy's internal number types, but that will handle data in the range you say you're dealing with.

Russell Borogove · Answer 5 · 2012-10-28T02:43:29.003

2

Works for me, Python 2.6.5 (CPython) on OSX 10.7:

>>> x = 2**50
>>> x
1125899906842624L
>>> p = int(log(x,2))
>>> p
50
>>> 2**p == x
True

It continues to work at least for exponents up to 1e9, by which time it starts to take quite a while to do the math. What are you actually getting for x and p in your test? What version of Python, on what OS, are you running?

edited Oct 28 '12 at 02:43

answered Oct 28 '12 at 02:35

Russell Borogove

18,516
4
43
50

Here's an [example of the failure](http://repl.it/DGb) -- you have to click the "run session" button. – Brian Cain Oct 28 '12 at 02:39
1

I get 49.999999 instead of 50.0 for p. Maybe it's because my PC is 32 bit – Rusty Rob Oct 28 '12 at 02:40
I wouldn't expect 32 vs 64 bit to affect math.log. Strange. – Russell Borogove Oct 28 '12 at 02:42

Brian Cain · Answer 6 · 2012-10-28T02:47:54.443

With respect to "not accurate for large numbers" your challenge here is that the floating point representation is indeed not as precise as you need it to be (49.999999999993 != 50.0). A great reference is "What Every Computer Scientist Should Know About Floating-Point Arithmetic."

The good news is that the transformation of the C routine is very straightforward:

def getpos(value):
    if (value == 0):
        return -1
    pos = 0
    if (value & (value - 1)):
        pos = 1
    if (value & 0xFFFFFFFF00000000):
        pos += 32
        value = value >> 32
    if (value & 0x00000000FFFF0000):
        pos += 16
        value = value >> 16
    if (value & 0x000000000000FF00):
        pos += 8
        value = value >> 8
    if (value & 0x00000000000000F0):
        pos += 4
        value = value >> 4
    if (value & 0x000000000000000C):
        pos += 2
        value = value >> 2
    if (value & 0x0000000000000002):
        pos += 1
        value = value >> 1
    return pos

Another alternative is that you could round to the nearest integer, instead of truncating:

   log(x,2)
=> 49.999999999999993
   round(log(x,2),1)
=> 50.0

score 0 · Answer 7 · answered Oct 31 '16 at 23:42

I needed to calculate the upper bound power of two (to figure out how many bytes of entropy was needed to generate a random number in a given range using the modulus operator).

From a rough experiment I think the calculation below gives the minimum integer p such that val < 2^p

It's probably about as fast as you can get, and uses exclusively bitwise integer arithmetic.

def log2_approx(val):
    from math import floor
    val = floor(val)
    approx = 0
    while val != 0:
        val &= ~ (1<<approx)
        approx += 1
    return approx

Your slightly different value would be calculated for a given n by

log2_approx(n) - 1

...maybe. But in any case, the bitwise arithmetic could give you a clue how to do this fast.

Compute fast log base 2 ceiling in python

7 Answers7

Linked