New appromixation for integer division by 63 or 127

[nikolaos....@gmail.com]

Dear all,

I believe I have invented a new approximation for quotient calculation that works for x/63 and x/127.

The formula is as follows and it is divisionless and multiplierless:

y = (((x>>n)+x+((1<<n)+1))>>n)-1;

Use n=6 for 63, and n=7 for 127. 1<<n is the strength-reduced form for calculating 2^n.

My formula does not require any double-word arithmetic (like with multiplying by multiplicative inverse and then performing adjustment steps). It does not seem to work for other dividers, but a systematic fix may be possible.

Are there any known references on this?

Best regards
Nikolaos Kavvadias
http://www.nkavvadias.com
http://github.com/nkkav/

[nikolaos....@gmail.com]

I have updated my claims for the algorithm. It works for all tested constant divisors (n=2 to n=16), and given the following assumption:

the range for x is [0,2^(2*n)-2], i.e. for n = 6, the range indeed is
0:4094.

And a test program:

[code]
#include <stdio.h>
#include <stdlib.h>
#include <math.h>

/* main:
*/
int main(void)
{
int qapprox, qexact;
int i, j, k;

for (i = 2; i < 8; i++) {
for (j = 1; j < (1<<i)*(1<<i)-1; j++) {
k = (1<<i) - 1;
qapprox = (((j>>i)+j+((1<<i)+1))>>i)-1;
qexact = j / k;
if (qapprox != qexact) {
fprintf(stderr, "qapprox = (%d/%d) = %d\tqexact = (%d/%d) = %d\n",
j, k, qapprox, j, k, qexact);
}
}
}
return 0;
}
[/code]

[Terje Mathisen]

nikolaos....@gmail.com wrote:
> Dear all,
>
> I believe I have invented a new approximation for quotient
> calculation that works for x/63 and x/127.
>
> The formula is as follows and it is divisionless and multiplierless:
>
> y = (((x>>n)+x+((1<<n)+1))>>n)-1;
>
> Use n=6 for 63, and n=7 for 127. 1<<n is the strength-reduced form
> for calculating 2^n.

Nothing new here.:-)

Instead of dividing by 63 you are multiplying by 65, and since
(x+1)(x-1) = x^2 - 1, the error is manageable.

Effectively you are doing reciprocal mul with a multiplier with very few
significant bits.

This also means that your formula only works for a limited range of
inputs (i.e. number of bits in x), and it fails badly when (for a 32-bit
number) x is so close to 2^32 that adding x/64 to it causes a carry.

BTW, I use a similar trick in my data calculations, when I need to
divide a small range of inputs by 100 in order to convert a year number
to century:

// y100 = y / 100; // Reciprocal mul
y100 = y * 41 >> 12; // Faster version with very short reciprocal

The multiplier has just three bits set (101001) so it can be calculated
with a pair of LEAs:

lea ebx, [eax+eax*4] ;; y*5
lea ebx, [eax+ebx*8] ;; y*41
shr ebx,12

The result is exact for all possible inputs in the 0..400 range which is
what I need at this point.

Terje

--
- <Terje.Mathisen at tmsw.no>
"almost all programming can be viewed as an exercise in caching"

[nikolaos....@gmail.com]

Hi Terje,

indeed, my formula:

- will actually work for any divisor of the form 2^n-1
- can be simplified to (x + (x>>n) + 1)>>n
- works for a limited range: x should be in the closed interval: [0, 2^(2n)-2]

I understand this is not new, but it has never been properly documented.
I have spent quite some time tracking a textbook (or paper) on this. Implicitly, yes, it's there, it can be inferred; but no direct reference to it.

BTW pkhuong contributed a nice number-theoretical proof for the formula regarding a somewhat larger range here:
http://www.reddit.com/r/programming/comments/2j18nt/multiplierless_hack_for_constant_division_by_2n1/

I also have a blog article on the topic:
- http://nkavvadias.com/blog/2014/10/11/multless-kdiv/

Further, it is useful and even more portable to have the formula working in any high-level language, I prefer it compared to x86 or ARM assembly trickery :). Your div-by-100 using leas for scaled multiplication is nice and you can do a similar thing with ARM's addressing modes.

I have written a tool for generating optimum code for ANSI C and "generic assembly". It is named kdiv and you can find it here: http://github.com/nkkav/kdiv
It is based on Hacker's Delight' ideas and code for calculating the magic number (multiplicative inverse)

There is also kmul, its sister (or brother?) project here: http://github.com/nkkav/kmul . This is based on Preston Briggs' code for Bernstein's algorithm on multiplication by integer constants.

Best regards
Nikolaos Kavvadias
http://www.nkavvadias.com