[Paul Curtis] over at Segger has an interesting series of blog posts about calculating division. This used to be a hotter topic, but nowadays many computers or computer languages have support for multiplication and division built-in. But some processors lack the instructions and a library to do it might be less than ideal. Knowing how to roll your own might allow you to optimize for speed or space. The current installment covers using Newton’s algorithm to do division.
Steve Martin had a famous bit about how to be a millionaire and never pay taxes. He started out by saying, “First… get a million dollar. Then…” This method is a bit like that since you first have to know how to multiply before you can divide. The basic premise is twofold: Newton’s method let you refine an estimate of a reciprocal by successive multiplications and then multiplying a number a reciprocal is the same as dividing. In other words, if we need to divide 34 by 6, you could rewrite 34/6 to 34 * 1/6 and the answer is the same.
Newton’s approximation for reciprocals lets you take a guess at the answer and then refine it through a series of multiplications. Each multiplication creates better accuracy. You can use this to perform a classic speed/space trade-off. For example, let’s just assume we want to find the reciprocal of a byte (presumably a fixed point byte). A look-up table of 256 elements would provide perfect accuracy and would be very fast. No more math is necessary. But what about 32 bits? Now the table is just too big. But you could look up, say, the first 8 bits of the 32-bit number. Or more. Or less. Depends on what’s important to you.
So now you have a poor estimate of your reciprocal. Sir Issac can make it better. For some number a, You take your estimate (x) and multiply them together. Subtract that number from 2 and you have a factor to multiply your old estimate by to get a new estimate. Skipping ahead, it is clear if your estimate was right, the multiplication would give you 1 which would not change the old estimate at all. If the estimate is off, you’ll get a scaling factor.
As a formula it looks like this:
So if you decide the reciprocal of 22 might be .02, the first pass will give you:
0.02*(2-22*0.02) = .0312 .0312*(2-22*.0312) = .0410 .0410*(2-22*.0410) = 0.0450
The right answer is a repeating decimal 0.0454545 and if you keep going, you’ll get there.
Of course, then you have to multiply one more time to do the division.
We liked that the post has a fixed-point implementation and then examines the resulting assembly code for ARM, RISC-V, and dsPIC30. Well worth a read.