# big integer division

links

* <https://skanthak.homepage.t-online.de/division.html>
* <https://news.ycombinator.com/item?id=26562819>
* <https://gmplib.org/~tege/division-paper.pdf>
* <https://github.com/Richard-Mace/huge-integer-class/blob/master/HugeInt.cpp> nice-looking well-commented c++ implementation
* <https://lists.libre-soc.org/pipermail/libre-soc-dev/2022-April/004739.html>

the most efficient division algorithm is probably Knuth's Algorithm D (with modifications from the exercises section of his book) which is O(n^2) and uses 2N-by-N-bit div/rem

an oversimplified version of the knuth algorithm d with 32-bit words is:
(TODO find original)

```
void div(uint32_t *n, uint32_t *d, uint32_t* q, int n_bytes, int d_bytes) {
    // assumes d[0] != 0, also n, d, and q have their most-significant-word in index 0
    int q_bytes = n_bytes - d_bytes;
    for(int i = 0; i < q_bytes / sizeof(n[0]); i++) {
        // calculate guess for quotient word
        q[i] = (((uint64_t)n[i] << 32) + n[i + 1]) / d[0];
        // n -= q[i] * d
        uint32_t carry = 0, carry2 = 0;
        for(int j = d_bytes / sizeof(d[0]) - 1; j >= 0; j--) {
            uint64_t v = (uint64_t)q[i] * d[j] + carry;
            carry = v >> 32;
            v = (uint32_t)v;
            v = n[i + j] - v + carry2;
            carry2 = v >> 32; // either ~0 or 0
            n[i + j] = v;
        }
        // fixup if carry2 != 0
    }
    // now remainder is in n and quotient is in q
}
```

```
On Sat, Apr 16, 2022, 22:06 Jacob Lifshay <programmerjake@gmail.com> wrote:
and a mrsubcarry (the one actually needed by bigint division):
# for big_c - big_a * word_b
result <- RC + ~(RA * RB) + CARRY # this expression is wrong, needs further thought
CARRY <- HIGH_HALF(result)
RT <- LOW_HALF(result)

turns out, after some checking with 4-bit words, afaict the correct algorithm for mrsubcarry is:
# for big_c - big_a * word_b
result <- RC + ~(RA * RB) + CARRY
result_high <- HIGH_HALF(result)
if CARRY <= 1 then # unsigned comparison
    result_high <- result_high + 1
end
CARRY <- result_high
RT <- LOW_HALF(result)

afaict, that'll make the following algorithm work:

so the inner loop in the bigint division algorithm would end up being (assuming n, d, and q all fit in registers):
li r3, 1 # carry in for subtraction
mtspr CARRY, r3 # init carry spr
setvl loop_count
sv.mrsubcarry rn.v, rd.v, rq.s, rn.v
```