Why is there a loop in this division as multiplication code?

Question

I got the js code below from an archive of hackers delight (view the source)

The code takes in a value (such as 7) and spits out a magic number to multiply with. Then you bitshift to get the results. I don't remember assembly or any math so I'm sure I'm wrong but I can't find the reason why I'm wrong

From my understanding you could get a magic number by writing ceil(1/divide * 1<<32) (or <<64 for 64bit values, but you'd need bigger ints). If you multiple an integer with imul you'd get the result in one register and the remainder in another. The result register is magically the correct result of a division with this magic number from my formula

I wrote some C++ code to show what I mean. However I only tested with the values below. It seems correct. The JS code has a loop and more and I was wondering, why? Am I missing something? What values can I use to get an incorrect result that the JS code would get correctly? I'm not very good at math so I didn't understand any of the comments

#include <cstdio>
#include <cassert>
int main(int argc, char *argv[])
{
    auto test_divisor = 7;
    auto test_value = 43;
    auto a = test_value*test_divisor;
    auto b = a-1; //One less test

    auto magic = (1ULL<<32)/test_divisor;
    if (((1ULL<<32)%test_divisor) != 0) {
        magic++; //Round up
    }
    auto answer1 = (a*magic) >> 32;
    auto answer2 = (b*magic) >> 32;
    assert(answer1 == test_value);
    assert(answer2 == test_value-1);
    printf("%lld %lld\n", answer1, answer2);
}

JS code from hackers delight

var two31 = 0x80000000
var two32 = 0x100000000
function magic_signed(d) { with(Math) {
    if (d >= two31) d = d - two32// Treat large positive as short for negative.
    var ad = abs(d)
    var t = two31 + (d >>> 31)
    var anc = t - 1 - t%ad       // Absolute value of nc.
    var p = 31                   // Init p.
    var q1 = floor(two31/anc)    // Init q1 = 2**p/|nc|.
    var r1 = two31 - q1*anc      // Init r1 = rem(2**p, |nc|).
    var q2 = floor(two31/ad)     // Init q2 = 2**p/|d|.
    var r2 = two31 - q2*ad       // Init r2 = rem(2**p, |d|).
    do {
        p = p + 1;
        q1 = 2*q1;                // Update q1 = 2**p/|nc|.
        r1 = 2*r1;                // Update r1 = rem(2**p, |nc|.
        if (r1 >= anc) {          // (Must be an unsigned
            q1 = q1 + 1;           // comparison here).
            r1 = r1 - anc;}
        q2 = 2*q2;                // Update q2 = 2**p/|d|.
        r2 = 2*r2;                // Update r2 = rem(2**p, |d|.
        if (r2 >= ad) {           // (Must be an unsigned
            q2 = q2 + 1;           // comparison here).
            r2 = r2 - ad;}
        var delta = ad - r2;
    } while (q1 < delta || (q1 == delta && r1 == 0))

    var mag = q2 + 1
    if (d < 0) mag = two32 - mag // Magic number and
    shift = p - 32               // shift amount to return.
    return mag
}}

Have you considered the impact of the code from HD using only 32-bit computation, while the replacement uses 64-bit computation? — njuffa, Jun 29 '22 at 04:09
is this some kind of reduction like Montgomery? why use floating operations for integer math? — Spektre, Jun 29 '22 at 06:53
Hacker's Delight does do some odd things. It appears to be looping until the quotient is below delta and it is a form of iteration. — QuentinUK, Jul 31 '22 at 11:18

Prem · Accepted Answer · 2022-07-31T11:01:55.583

In the C CODE:

auto magic = (1ULL<<32)/test_divisor;

We get Integer Value in magic because both (1ULL<<32) & test_divisor are Integers.
The Algorithms requires incrementing magic on certain conditions, which is the next conditional statement.

Now, multiplication also gives Integers:
auto answer1 = (a*magic) >> 32;
auto answer2 = (b*magic) >> 32;

C CODE is DONE !

In the JS CODE:

All Variables are var ; no Data types !
No Integer Division ; No Integer Multiplication !
Bitwise Operations are not easy and not suitable to use in this Algorithm. Numeric Data is via Number & BigInt which are not like "C Int" or "C Unsigned Long Long".

Hence the Algorithm is using loops to Iteratively add and compare whether "Division & Multiplication" has occurred to within the nearest Integer.

Both versions try to Implement the same Algorithm ; Both "should" give same answer, but JS Version is "buggy" & non-standard.
While there are many Issues with the JS version, I will highlight only 3:

(1) In the loop, while trying to get the best Power of 2, we have these two statements :

        p = p + 1;
        q1 = 2*q1;                // Update q1 = 2**p/|nc|.

It is basically incrementing a counter & multiplying a number by 2, which is a left shift in C++.
The C++ version will not require this rigmarole.

(2) The while Condition has 2 Equality comparisons on RHS of || :

while (q1 < delta || (q1 == delta && r1 == 0))

But both these will be false in floating Point Calculations [[ eg check "Math.sqrt(2)*Math.sqrt(0.5) == 1" : even though this must be true, it will almost always be false ]] hence the while Condition is basically the LHS of || , because RHS will always be false.

(3) The JS version returns only one variable mag but user is supposed to get (& use) even variable shift which is given by global variable access. Inconsistent & BAD !

Comparing , we see that the C Version is more Standard, but Point is to not use auto but use int64_t with known number of bits.

Hmm. I still don't understand the JS code but it sounds like what they are doing is the same thing? Would both code always produce the same answer? If not which of the two is correct? I rather not get incorrect results if I could help it — Cal, Jul 01 '22 at 03:37
Both try to Implement same Algorithm ; Both "should" give same answer, but JS Version is "buggy" , Eg the `while` has 2 Equality comparisons on RHS of || , but both will be false in floating Point Calculations [[ eg check "Math.sqrt(2)*Math.sqrt(0.5) == 1" : this will be false ]] hence the `while` Condition is basically LHS of ||, because RHS will always be false. The C Version is more Standard, but Point is to not use `auto` but use `int64_t` with known number of bits. Example of the JS Version loop :: AxB can be Calculated by M=A+A+A... where B is the loop Counter. C Version is better !! — Prem, Jul 01 '22 at 04:47
I'll give other people a few days to try to answer to but I'm likely to accept this — Cal, Jul 01 '22 at 05:04
I have updated my Answer with more Details !! Gentle reminder @Cal — Prem, Jul 31 '22 at 11:08
I completely forgot about this. I'm not on the site much. Accepted. — Cal, Aug 03 '22 at 23:56

score 0 · Answer 2 · answered Jun 29 '22 at 09:24

First I think ceil(1/divide * 1<<32) can, depending on the divide, have cases where the result is off by one. So you don't need a loop but sometimes you need a corrective factor.

Secondly the JS code seems to allow for other shifts than 32: shift = p - 32 // shift amount to return. But then it never returns that. So not sure what is going on there.

Why not implement the JS code in C++ as well and then run a loop over all int32_t and see if they give the same result? That shouldn't take too long.

And when you find a d where they differ you can then test a / d for all int32_t a using both magic numbers and compare a / d, a * m_ceil and a * m_js.

Why is there a loop in this division as multiplication code?

2 Answers2

In the C CODE:

In the JS CODE: