64 / 32 => 32 Unsigned Divide

[Phil Pilgrim (PhiPi)]

I needed a division routine that divides a 64-bit unsigned dividend by a 32-bit unsigned divisor to return a 32-bit quotient and 32-bit remainder, leaving the original carry bit unchanged. The only limitation would have to be that the high long of the dividend be less than the divisor (i.e. no overflow). I couldn't find a routine like this in the forum or in the object exchange (although I'm sure that deSilva, who has a memory like a steel trap, can point to one immediately). Here's the code I came up with:




'Divide ra:rb by rx, leaving quotient in rb, remainder in ra, and rx unchanged.
'Precondition: ra < rx.

divabx mov cctr,#32 'Initialize loop counter.
:loop rcl rb,#1 wc 'Rotate quotient bit in from carry,
rcl ra,#1 wc ' and rotate dividend out to carry
if_c sub ra,rx 'rx < carry:ra if carry set, so just subtract, leaving carry set.
if_nc cmpsub ra,rx wc 'Otherwise, use cmpsub to do the dirty work.
djnz cctr,#:loop 'Back for more.
rcl rb,#1 wc 'Rotate last quotient bit into place, restoring original carry.
divabx_ret ret




It's pretty straightforward, except for the inconvenient fact that cmpsub doesn't treat the carry flag as a most-significant bit when doing its comparison; so that condition has to be tested separately. (The consequence of this behavior for division routines that use cmpsub exclusively is that the high bit of the divisor must be zero.) The original carry bit simply propagates through rb and back out again, so is preserved.

Exclusive of the call and return, this routine takes 8.2µs to execute with an 80MHz clock.

-Phil

[deSilva]

The standard 32/16 routine is this:


' Divide x[ 31..0] by y[ 15..0] (y[16] must be 0)
' on exit, quotient is in x[ 15..0] and remainder is in x[ 31..16]
'
divide shl y, #15 'get divisor into y[ 30..15]
mov t, #16 'ready for 16 quotient bits
:loop
cmpsub x, y wc 'if y =< x then subtract it, set C
rcl x, #1 'rotate c into quotient, shift dividend
djnz t, #:loop 'loop until done
divide_ret
ret 'quotient in x[ 15..0], rem. in x[ 31..16]


Which is essentially the same.
Note that they have not explicitely stated as Phil did, that "....the high long of the dividend be less than the divisor."

Now if we imagine we had only 16 bits words at our disposal how would that routine had to be modified?

We see immediately that we will already run into trouble with the first instruction (SHL #15). So, we first reduce our requirements to dividing 63 bit numbers by 31 bit numbers only http://forums.parallax.com/images/smilies/smile.gif Tis is not so grave, as signed numbers are just of that kind...



divide2 shl y, #16 ' get divisor into y[ 31..16]
shl x, #1 ' to compensate for this
mov t, #16 ' ready for 16 quotient bits
:loop


Now directly adapting it to 63/31 bit division would yield:



divide3 nop ' leave divisor in y
rcl xl,#1 wc ' to compensate for one bit position
rcl xh,#1
mov t, #32 ' ready for 32 quotient bits
:loop
cmpsub xh, y wc 'if yh =< x then subtract it, set C
rcl xl, #1 'rotate c into quotient...
rcl xh, #1 '... shift dividend
djnz t, #:loop 'loop until done

shr xh, #1 ' but the remainder had been shifted once too often
divide3_ret


It is obvious now, that we can save an instruction by CAREFULLY exchanging CMPSUB and the RCLs



divide4
mov t, #32 ' ready for 32 quotient bits
' but mind the very first carry!!
:loop
rcl xl, #1 wc 'rotate c into quotient...
rcl xh, #1 '... shift dividend
cmpsub xh, y wc 'if yh =< x then subtract it, set C
djnz t, #:loop 'loop until done

' now the remnainder had been shifted correctly, but not the result!
rcl xl, #1
divide4_ret



This looks like a straightforward transformation... It comes out very nicely withthe preserved carry flag. I think the special consideration Phil took for the CMPSUB overflow has to do with his strive for 64/32 rather than 63/31 division. I THINK it works, but I am not absolutely sure http://forums.parallax.com/images/smilies/smile.gif

Post Edited (deSilva) : 1/27/2008 12:17:10 AM GMT