Fast modulo 10 in c

Question

I am looking for a fast modulo 10 algorithm because I need to speed up my program which does many modulo operations in cycles.

I have checked out this page which compares some alternatives. As far as I understand it correctly, T3 was the fastest of all. My question is, how would x % y look like using T3 technique?

I copied T3 technique here for simplicity in case the link gets down.

for (int x = 0; x < max; x++)
{
        if (y > (threshold - 1))
        {
               y = 0; //reset
               total += x;
        }
        y += 1;
}

Regarding to comments, if this is not really faster then regular mod, I am looking for at least 2 times faster modulo than using %. I have seen many examples with use power of two, but since 10 is not, how can I get it to work?

Edit:

For my program, let's say I have 2 for cycles where n=1 000 000 and m=1000.

Looks like this:

for (i = 1; i <= n; i++) {
        D[(i%10)*m] = i;
        for (j = 1; j <= m; j++) {
           ...
        }
}

Answer 1

Here's the fastest modulo-10 function you can write:

unsigned mod10(unsigned x)
{
    return x % 10;
}

And here's what it looks like once compiled:

movsxd rax, edi
imul rcx, rax, 1717986919
mov rdx, rcx
shr rdx, 63
sar rcx, 34
add ecx, edx
add ecx, ecx
lea ecx, [rcx + 4*rcx]
sub eax, ecx
ret

Note the lack of division/modulus instructions, the mysterious constants, the use of an instruction which was originally intended for complex array indexing, etc. Needless to say, the compiler knows a lot of tricks to make your program as fast as possible. You'll rarely beat it on tasks like this.

Answer 2

You likely can't beat the compiler.

Debug build

//     int foo = x % 10;
010341C5  mov         eax,dword ptr [x]  
010341C8  cdq  
010341C9  mov         ecx,0Ah  
010341CE  idiv        eax,ecx  
010341D0  mov         dword ptr [foo],edx  

Retail build (doing some ninja math there...)

//    int foo = x % 10;
00BD100E  mov         eax,66666667h  
00BD1013  imul        esi  
00BD1015  sar         edx,2  
00BD1018  mov         ecx,edx  
00BD101A  shr         ecx,1Fh  
00BD101D  add         ecx,edx  
00BD101F  lea         eax,[ecx+ecx*4]  
00BD1022  add         eax,eax  
00BD1024  sub         esi,eax

Answer 3

The code isn’t a direct substitute for modulo, it substitutes modulo in that situation. You can write your own mod by analogy (for a, b > 0):

int mod(int a, int b) {
    while (a >= b) a -= b;
    return a;
}

… but whether that’s faster than % is highly questionable.

Answer 4

I come across this discussion and while for uint64_t the best way to perform a mod 10 operation is indeed through the usage of the compiler on my standard laptop. However for unt128_t on my recent ubuntu linux i get, for the routine:

for (int i = 0; i < 1000000000; i++)
{
  uint128_t x = n + i;
  s += x % 10;
}

The timing:

   Executed in   21,74 secs   fish           external 
   usr time      21,73 secs  420,00 micros   21,73 secs 
   sys time       0,00 secs  237,00 micros    0,00 secs 

This is very different from the result I got from using uint64_t instead. So one could expect to do something clever here (and I bet that in future versions of gcc they would implement some form of the following trick). We can take advantage of the rules,

(a+b) mod 10   = (a mod 10 + b mod 10) mod 10

And

(ab) mode 10 = ((a mode 10)*(b mod 10)) mod 10

To produce the code,

for (int i = 0; i < 1000000000; i++)
{
  uint128_t x = n + i;
  uint64_t  a = (uint64_t)(x >> 64);
  uint64_t  b = (uint64_t)(x & (~0UL));
  
  s += ((a%10)*2*((1UL<<63)%10) + (b%10))%10;
}

This benchmarks at speed,

Executed in    3,55 secs   fish           external 
usr time       3,55 secs  409,00 micros    3,55 secs 
sys time       0,00 secs  233,00 micros    0,00 secse here

A nice 5x speedup for the modulo 10 operation. Note that 10 is not magic here apart from the compiler possible being extra smart about 10 for 64 bit unsigned integers. Similar trick can be done for integer division by 10 where we note that we always can write a number x as x = a 10 + b, where a = x/10 and b = x%10, then again we can study x1*x2 and x1+x2 to deduce similar rules for the integer division of 128 bit integers taking advantage of the fast versions of the 64 bit ones. If one do the work we can produce the following code,

inline uint128_t div10q(uint128_t x)
{
  uint64_t   x1  = (uint64_t)(x >> 64);
  uint128_t  x2  = ((unt128_t)1) << 64;
  uint64_t   x3  = (uint64_t)(x & (~0UL));

  uint64_t   b1   = x1%10;
  uint128_t  y1  = x1/10;

  uint64_t   b2  = x2%10;
  uint128_t  y2  = x2/10;

  uint128_t yy1 = y1*y2*10+b1*y2+b2*y1 + (b1*b2)/10;
  uint64_t  bb1 = (b1*b2)%10;

  uint64_t  bb2 = x3 % 10;
  uint128_t yy2 = x3 / 10;

  return yy1+yy2+(bb1+bb2)/10;
}

That compiles to a similar speedup of 5X using the optimization -O3 in gcc.

Answer 5

This will work for (multiword) values larger than the machineword (but assuming a binary computer ...):

---

#include <stdio.h>

unsigned long mod10(unsigned long val)
{
unsigned res=0;

res =val &0xf;
while (res>=10) { res -= 10; }

for(val >>= 4; val; val >>= 4){
        res += 6 * (val&0xf);
        while (res >= 10) { res -= 10; }
        }

return res;
}

int main (int argc, char **argv)
{
unsigned long val;
unsigned res;

sscanf(argv[1], "%lu", &val);

res = mod10(val);
printf("%lu -->%u\n", val,res);

return 0;
}

---

UPDATE: With some extra effort, you could get the algoritm free of multiplications, and with the proper amount of optimisation we can even get the recursive call inlined:

---

static unsigned long mod10_1(unsigned long val)
{
unsigned char res=0; //just to show that we don't need a big accumulator

res =val &0xf; // res can never be > 15
if (res>=10) { res -= 10; }

for(val >>= 4; val; val >>= 4){
        res += (val&0xf)<<2 | (val&0xf) <<1;
        res= mod10_1(res); // the recursive call
        }

return res;
}

---

And the result for mod10_1 appears to be mul/div free and almost without branches:

---

mod10_1:
.LFB25:
    .cfi_startproc
    movl    %edi, %eax
    andl    $15, %eax
    leal    -10(%rax), %edx
    cmpb    $10, %al
    cmovnb  %edx, %eax
    movq    %rdi, %rdx
    shrq    $4, %rdx
    testq   %rdx, %rdx
    je      .L12
    pushq   %r12
    .cfi_def_cfa_offset 16
    .cfi_offset 12, -16
    pushq   %rbp
    .cfi_def_cfa_offset 24
    .cfi_offset 6, -24
    pushq   %rbx
    .cfi_def_cfa_offset 32
    .cfi_offset 3, -32
.L4:
    movl    %edx, %ecx
    andl    $15, %ecx
    leal    (%rcx,%rcx,2), %ecx
    leal    (%rax,%rcx,2), %eax
    movl    %eax, %ecx
    movzbl  %al, %esi
    andl    $15, %ecx
    leal    -10(%rcx), %r9d
    cmpb    $9, %cl
    cmovbe  %ecx, %r9d
    shrq    $4, %rsi
    leal    (%rsi,%rsi,2), %ecx
    leal    (%r9,%rcx,2), %ecx
    movl    %ecx, %edi
    movzbl  %cl, %ecx
    andl    $15, %edi
    testq   %rsi, %rsi
    setne   %r10b
    cmpb    $9, %dil
    leal    -10(%rdi), %eax
    seta    %sil
    testb   %r10b, %sil
    cmove   %edi, %eax
    shrq    $4, %rcx
    andl    $1, %r10d
    leal    (%rcx,%rcx,2), %r8d
    movl    %r10d, %r11d
    leal    (%rax,%r8,2), %r8d
    movl    %r8d, %edi
    andl    $15, %edi
    testq   %rcx, %rcx
    setne   %sil
    leal    -10(%rdi), %ecx
    andl    %esi, %r11d
    cmpb    $9, %dil
    seta    %bl
    testb   %r11b, %bl
    cmovne  %ecx, %edi
    andl    $1, %r11d
    andl    $240, %r8d
    leal    6(%rdi), %ebx
    setne   %cl
    movl    %r11d, %r8d
    andl    %ecx, %r8d
    leal    -4(%rdi), %ebp
    cmpb    $9, %bl
    seta    %r12b
    testb   %r8b, %r12b
    cmovne  %ebp, %ebx
    andl    $1, %r8d
    cmovne  %ebx, %edi
    xorl    $1, %ecx
    andl    %r11d, %ecx
    orb     %r8b, %cl
    cmovne  %edi, %eax
    xorl    $1, %esi
    andl    %r10d, %esi
    orb     %sil, %cl
    cmove   %r9d, %eax
    shrq    $4, %rdx
    testq   %rdx, %rdx
    jne     .L4
    popq    %rbx
    .cfi_restore 3
    .cfi_def_cfa_offset 24
    popq    %rbp
    .cfi_restore 6
    .cfi_def_cfa_offset 16
    movzbl  %al, %eax
    popq    %r12
    .cfi_restore 12
    .cfi_def_cfa_offset 8
    ret
.L12:
    movzbl  %al, %eax
    ret
    .cfi_endproc
.LFE25:
    .size   mod10_1, .-mod10_1
    .p2align 4,,15
    .globl  mod10
    .type   mod10, @function