Fastest way to compute order of magnitude in x86 assembly

Question

The task is simple: write assembly that computes the order of magnitude of an integer using as few clock cycles as possible.

• Order of magnitude is defined as log10, not log2.
• The range of valid input is 0 to 1012, inclusive. Behaviour for input outside that range is undefined.
• Values should be rounded down to the nearest integer, with the exception that given input 0 the output should be 0. (You can consider this to be the best representation of negative infinity that's possible in unsigned integers).
• Must be x86 assembly.
• The integer must be a runtime value, not a static/inline integer. So we don't know what it is at compile time.
• Assume that you have an integer loaded into a register already. (But include setting the value in the register in the answer for clarity).
• Can't call any external libraries or functions.
• Free to use any of the available instructions in the Intel docs.
• No C.
• Any of the ~7 Intel Core architectures are acceptable (listed on page 10). Ideally Nehalem (Intel Core i7).

The winning answer is one that uses the fewest clock cycles possible. That is, it can have the most operations per second. Approximate clock cycle summaries are here: http://www.agner.org/optimize/instruction_tables.pdf. Calculating clock cycles can happen after the answer is posted.

Answer 1

7 cycles, constant time

Here's a solution based on my answer to this SO Question. It uses BSR to count how many bits are needed to hold the number. It looks up how many decimal digits needed to represent the largest number that many bits can hold. Then it subtracts 1 if the actual number is less than the nearest power of 10 with that many digits.

    .intel_syntax noprefix
    .globl  main
main:
    mov rdi, 1000000000000              #;your value here
    bsr rax, rdi
    movzx   eax, BYTE PTR maxdigits[1+rax]
    cmp rdi, QWORD PTR powers[0+eax*8]
    sbb al, 0
    ret
maxdigits:
    .byte   0
    .byte   0
    .byte   0
    .byte   0
    .byte   1
    .byte   1
    .byte   1
    .byte   2
    .byte   2
    .byte   2
    .byte   3
    .byte   3
    .byte   3
    .byte   3
    .byte   4
    .byte   4
    .byte   4
    .byte   5
    .byte   5
    .byte   5
    .byte   6
    .byte   6
    .byte   6
    .byte   6
    .byte   7
    .byte   7
    .byte   7
    .byte   8
    .byte   8
    .byte   8
    .byte   9
    .byte   9
    .byte   9
    .byte   9
    .byte   10
    .byte   10
    .byte   10
    .byte   11
    .byte   11
    .byte   11
    .byte   12
powers:
    .quad   0
    .quad   10
    .quad   100
    .quad   1000
    .quad   10000
    .quad   100000
    .quad   1000000
    .quad   10000000
    .quad   100000000
    .quad   1000000000
    .quad   10000000000
    .quad   100000000000
    .quad   1000000000000

Compiles on GCC 4.6.3 for ubuntu and returns the value in the exit code.

I'm not confident interpreting that cycles table for any modern processor, but using @DigitalTrauma's method, on the Nehalim processor, I get 7?

ins        uOps Latency
---        -    - 
BSR r,r  : 1    3
MOVZX r,m: 1    -
CMP m,r/i: 1    1 
SBB r,r/i: 2    2

Answer 2

Best case 8 cycles, Worst case 12 cycles

Since it is not clear in the question, I am basing this off the Ivy Bridge latencies.

The approach here is to use the bsr (bit scan reverse) instruction as a poor-man's log2(). The result is used as an index into a jump table which contains entries for bits 0 to 42. I am assuming that given that operation on 64bit data is implicitly required, then use of the bsr instruction is OK.

In the best case inputs, the jumptable entry can map the bsr result directly to the magnitude. e.g. For inputs in the range 32-63, the bsr result will be 5, which is mapped to a magnitude of 1. In this case, the instruction path is:

Instruction    Latency

bsrq                 3
jmp                  2
movl                 1
jmp                  2

total                8

In the worst case inputs, the bsr result will map to two possible magnitudes, so the jumptable entry does one additional cmp to check if the input is > 10ⁿ. E.g. for inputs in the range 64-127, the bsr result will be 6. The corresponding jumptable entry then checks if the input > 100 and sets the output magnitude accordingly.

In addition for the worst case path, we have an additional mov instruction to load a 64bit immediate value for use in the cmp, so the worst case instruction path is:

Instruction    Latency

bsrq                 3
jmp                  2
movabsq              1
cmpq                 1
ja                   2
movl                 1
jmp                  2

total               12

Here is the code:

    /* Input is loaded in %rdi */
    bsrq    %rdi, %rax
    jmp *jumptable(,%rax,8)
.m0:
    movl    $0, %ecx
    jmp .end
.m0_1:
    cmpq    $9, %rdi
    ja  .m1
    movl    $0, %ecx
    jmp .end
.m1:
    movl    $1, %ecx
    jmp .end
.m1_2:
    cmpq    $99, %rdi
    ja  .m2
    movl    $1, %ecx
    jmp .end
.m2:
    movl    $2, %ecx
    jmp .end
.m2_3:
    cmpq    $999, %rdi
    ja  .m3
    movl    $2, %ecx
    jmp .end
.m3:
    movl    $3, %ecx
    jmp .end
.m3_4:
    cmpq    $9999, %rdi
    ja  .m4
    movl    $3, %ecx
    jmp .end
.m4:
    movl    $4, %ecx
    jmp .end
.m4_5:
    cmpq    $99999, %rdi
    ja  .m5
    movl    $4, %ecx
    jmp .end
.m5:
    movl    $5, %ecx
    jmp .end
.m5_6:
    cmpq    $999999, %rdi
    ja  .m6
    movl    $5, %ecx
    jmp .end
.m6:
    movl    $6, %ecx
    jmp .end
.m6_7:
    cmpq    $9999999, %rdi
    ja  .m7
    movl    $6, %ecx
    jmp .end
.m7:
    movl    $7, %ecx
    jmp .end
.m7_8:
    cmpq    $99999999, %rdi
    ja  .m8
    movl    $7, %ecx
    jmp .end
.m8:
    movl    $8, %ecx
    jmp .end
.m8_9:
    cmpq    $999999999, %rdi
    ja  .m9
    movl    $8, %ecx
    jmp .end
.m9:
    movl    $9, %ecx
    jmp .end
.m9_10:
    movabsq $9999999999, %rax
    cmpq    %rax, %rdi
    ja  .m10
    movl    $9, %ecx
    jmp .end
.m10:
    movl    $10, %ecx
    jmp .end
.m10_11:
    movabsq $99999999999, %rax
    cmpq    %rax, %rdi
    ja  .m11
    movl    $10, %ecx
    jmp .end
.m11:
    movl    $11, %ecx
    jmp .end
.m11_12:
    movabsq $999999999999, %rax
    cmpq    %rax, %rdi
    ja  .m12
    movl    $11, %ecx
    jmp .end
.m12:
    movl    $12, %ecx
    jmp .end

jumptable:
    .quad   .m0
    .quad   .m0
    .quad   .m0
    .quad   .m0_1
    .quad   .m1
    .quad   .m1
    .quad   .m1_2
    .quad   .m2
    .quad   .m2
    .quad   .m2_3
    .quad   .m3
    .quad   .m3
    .quad   .m3
    .quad   .m3_4
    .quad   .m4
    .quad   .m4
    .quad   .m4_5
    .quad   .m5
    .quad   .m5
    .quad   .m5_6
    .quad   .m6
    .quad   .m6
    .quad   .m6
    .quad   .m6_7
    .quad   .m7
    .quad   .m7
    .quad   .m7_8
    .quad   .m8
    .quad   .m8
    .quad   .m8_9
    .quad   .m9
    .quad   .m9
    .quad   .m9
    .quad   .m9_10
    .quad   .m10
    .quad   .m10
    .quad   .m10_11
    .quad   .m11
    .quad   .m11
    .quad   .m11_12
    .quad   .m12
    .quad   .m12
    .quad   .m12

.end:
/* output is given in %ecx */

This was mostly generated from gcc assembler output for proof-of-concept C code I wrote. Note the C code uses a computable goto to implement the jump table. It also uses the __builtin_clzll() gcc builtin, which compiles to the bsr instruction (plus an xor).

---

I considered several solutions before arriving at this one:

• FYL2X to calculate the natural log, then FMUL by the necessary constant. This would presumably win this if it was an [tag:instruction:golf] contest. But FYL2X has as latency of 90-106 for Ivy bridge.

• Hard-coded binary search. This may actually be competitive - I'll leave it to someone else to implement :).

• Complete lookup table of results. This I'm sure is theoretically fastest, but would require a 1TB lookup table - not practical yet - maybe in a few years if Moore's Law continues to hold.