How do the Math.divide & Math.sqrt algorithms work?

[garamor]

Hello, I have been trying to understand the Math.divide algorithm that was suggested in the videos series in unit 6.3 for a while. I have made some progress and walked myself through some example inputs and outputs, however it still mostly just seems like magic. Can someone help me really understand how it works? Here is my code:

function int divide(int x, int y) {
    var int q;
    if (y > x) {
        return 0;
    }
    let q = Math.divide(x, y + y)
    if (x - (q + q * y) < y) {
        return q + q;
    }
    else {
        return q + q + 1;
    }
}


I am also confused by the sqrt algorithm. if anyone has the time to help me understand this too, it would be greatly appreciated!:

sqrt(x):
    y = 0
    for j = n/2-1 ... 0 do
        if (y+2^j)^2 <= x then y = y + 2
    return y

[dolomiti7]

First of all your code is not correct or at least misleading:

1. divide: brackets missing
 if (x - ( ( q + q  ) * y) < y) {
it may still work without the brackets in Jack since there is no operator precedence, and depending on the parser the addition might still be evaluated as intended before the multiplication. However, this is risky since other parsers might not evaluate it in the same way. And for sure it is not intuitive without the brackets since mathematically (and in probably all other languages) q + q * y is evaluated as q + (q * y). I guess it is also easier to understand with the brackets:
- recursively shift left by 1 bit (y + y) until it is too big
- Afterwards you shift left the current result q (q + q), test if the (shifted) result still fits into x and depending on that set the lowest bit to 0 or 1 (either q+q or q+q+1)

Please note that this algorithm is at risk of returning wrong results due to an overflow for high numbers. If x>y>=16384 the y+y in the recursive division call will overflow to negative numbers and return wrong results.


2. sqrt: missing power at the end:
 if (y+(2^j))^2 <= x then y = y + (2 ^ j)

The algorithm is basically a binary search, trying out bit by bit.

Since the square root only needs half of the bits of the operand, you can start at n/2-1  (*):
       log2(x)=log2(sqrt(x)^2)=2*log2(sqrt(x))
 => log2(sqrt(x))/log2(x)=1/2 => half number of bits needed
For the Hack platform, the highest possible square root is 181 (sqrt(32767)),  which requires 8 bits (from bit 7 to bit 0) as expected. So the loop needs to go from (16/2)-1=7 to 0.
(*) a minor inaccuracy in the book: n/2 should be rounded up. Technically we are handling only 15-bit values for positive numbers, n/2-1 with n=15 would then lead to a loop from 6 to 0 which is wrong. Even for 15 bits we need to go from bit 7 to 0. In this case the round up could be achieved with (n+1)/2 - 1, which leads to 7 for both n=15 and n=16.

You try each bit from high to low (the order is important) and test if the square of the result would still be below or equal the operand. If yes, the respective bit must be 1 and remains in the result (y=y+2^j), if not you continue with the next lower bit:
2^7=128 -> is 128*128 <= x ?  -> if yes, keep bit 7, e.g. result is >=128
2^6=64 -> is (result + 64)*(result + 64) <= x ? if yes, keep bit 6,...

[garamor]

This has all been very helpful, and I am sorry for not replying sooner. I have gone on and am now working on the String library. When testing it I ran into an overflow error with my division. Do you have any clue how I could solve this?

[garamor]

I figured out how solve the division overflow error in cases where x is greater than 16384. Here is my code:

    function int divide(int x, int y) {
        var int q;

        if (y = 0) {
            do Sys.error(3);
            return 0;
        }

        let q = Math._div(Math.abs(x), Math.abs(y));

        if ((x < 0) = (y < 0)) {return q;}
        else                   {return -q;}
    }
    function int _div(int x, int y) {
        var int q, sum;
        if (x < y) {return 0;}
        if (y < 16384) {let q = Math._div(x, y + y);}
        else           {let q = 0;}
        let sum = q + q;
        if (x - ((sum) * y) < y) {return sum;}
        else                     {return sum + 1;}
    }

[dolomiti7]

It appears to me as if the algorithm catches most of the potential overflows, but have you tested the correctness with all possible inputs and corner cases? What about dividing -32768 by y (*)? In this case your Math.abs will overflow since there is no positive representation for -32768 in a 16-bit Two's complement. X will still be negative inside _div and therefore return 0. Or dividing by -32768: in this case x will never be smaller than y, y will be <16384 and your call to Math._div(x, y + y) will lead to a division by zero error, since (-32768)+(-32768)=0 in 16-bit...

As elegant as the recursive division algorithm is, I am not a big fan of it in terms of practicability. It takes quite some efforts to make it bulletproof in the implementation.

(*) In Jack -32768 cannot be entered as a literal, but needs to be constructed in the source code
i.e. (-32767-1)