24 bit multiplication

[ACagliano]

Ok. I am particularly interested now in 2-byte multiplication and 4-byte square rooting. How would they be done?

[jacobly]

Code: Select All | Copy To Clipboard

1 2 3 4 5 6 7 8 9

// Multiply a times b temp = 0 repeat for each bit in a temp <<= 1 if (high bit of a set) temp += b a <<= 1 return temp

if a and b are 2 bytes, temp is 4 bytes, and you loop 16 times.

Spoiler for for code:

stolen from Axe
p_MulFull:
   ; Input in hl, result in cahl
   ld   c,h
   ld   a,l
   ld   hl,0   ;11
   ld   b,16   ;7
__MulFullNext:
   add   hl,hl   ;11
   rla      ;4
   rl   c   ;8
   jr   nc,__MulFullSkip   ;12/7
   add   hl,de   ;11
   adc   a,0   ;7
   jr   nc,__MulFullSkip
   inc   c
__MulFullSkip:
   djnz   __MulFullNext
   ret
__MulFullEnd:

Code: Select All | Copy To Clipboard

1 2 3 4 5 6 7 8 9 10 11 12 13 14

// Sqrt a temp = high byte of a a <<= 8 b = 0 repeat for every 2 bits in a test = b << 8 + 0x40 b <<= 1 if (temp >= test) temp -= test set low bit of b temp += high 2 bits of a a <<= 2 return b

If a is 4 bytes, then b and temp are 2 bytes, and you loop 16 times.

Spoiler for code:

stole my own routine from axe (and modified it)
p_Sqrt88:
   ; input in hlde, result in de
   ld   b,16
   ld   a,h
   ld   c,l
   push   de ; ld ixh,d
   pop   ix ; ld ixl,e
   ld   de,0
   ld   h,d
   ld   l,e
__Sqrt88Loop:
   sub   $40
   sbc   hl,de
   jr   nc,__Sqrt88Skip
   add   a,$40
   adc   hl,de
__Sqrt88Skip:
   ccf
   rl   e
   rl   d
   add   ix,ix
   rl   c
   rla
   adc   hl,hl
   add   ix,ix
   rl   c
   rla
   adc    hl,hl
   djnz   __Sqrt88Loop
   ret
__Sqrt88End:

[ACagliano]

Code: Select All | Copy To Clipboard

1 2 3 4 5 6 7 8 9

// Multiply a times b temp = 0 repeat for each bit in a temp <<= 1 if (high bit of a set) temp += b a <<= 1 return temp

if a and b are 2 bytes, temp is 4 bytes, and you loop 16 times.

Spoiler for for code:

stolen from Axe
p_MulFull:
   ; Input in hl, result in cahl
   ld   c,h
   ld   a,l
   ld   hl,0   ;11
   ld   b,16   ;7
__MulFullNext:
   add   hl,hl   ;11
   rla      ;4
   rl   c   ;8
   jr   nc,__MulFullSkip   ;12/7
   add   hl,de   ;11
   adc   a,0   ;7
   jr   nc,__MulFullSkip
   inc   c
__MulFullSkip:
   djnz   __MulFullNext
   ret
__MulFullEnd:

Code: Select All | Copy To Clipboard

1 2 3 4 5 6 7 8 9 10 11 12 13 14

// Sqrt a temp = high byte of a a <<= 8 b = 0 repeat for every 2 bits in a test = b << 8 + 0x40 b <<= 1 if (temp >= test) temp -= test set low bit of b temp += high 2 bits of a a <<= 2 return b

If a is 4 bytes, then b and temp are 2 bytes, and you loop 16 times.

Spoiler for code:

stole my own routine from axe (and modified it)
p_Sqrt88:
   ; input in hlde, result in de
   ld   b,16
   ld   a,h
   ld   c,l
   push   de ; ld ixh,d
   pop   ix ; ld ixl,e
   ld   de,0
   ld   h,d
   ld   l,e
__Sqrt88Loop:
   sub   $40
   sbc   hl,de
   jr   nc,__Sqrt88Skip
   add   a,$40
   adc   hl,de
__Sqrt88Skip:
   ccf
   rl   e
   rl   d
   add   ix,ix
   rl   c
   rla
   adc   hl,hl
   add   ix,ix
   rl   c
   rla
   adc    hl,hl
   djnz   __Sqrt88Loop
   ret
__Sqrt88End:

Kool. Thanks. But, shouldn't the first one have two inputs?

[Xeda112358]

So with two-byte multiplication, you can take advantage of the fact that add hl,hl is the same as shifting hl left. It even gives you the carry! So in this case:

Code: Select All | Copy To Clipboard

1 2 3 4 5 6 7 8 9 10 11 12 13

ld hl,0      ld a,16 MultLoop:      add hl,hl      ;shifts hl left      rl e \ rl d    ;shifts de left and if hl overflowed, it overflows into de      jr nc,$+6      ;if the bit in DE is o, skip this chunk        add hl,bc    ;add bc to hl (think of this as the first number)        jr nc,$+3    ;overflow into de          inc de      dec a      jr nz,MultLoop      ret

That will multiply DE times BC and return the result in DEHL. I will see if I can port a square root routine for 32-bit...

EDIT: changed inc e to inc de

[jacobly]

Code: Select All | Copy To Clipboard

1 2 3 4 5 6 7 8 9

// Multiply a times b temp = 0 repeat for each bit in a temp <<= 1 if (high bit of a set) temp += b a <<= 1 return temp

if a and b are 2 bytes, temp is 4 bytes, and you loop 16 times.

Spoiler for for code:

stolen from Axe
p_MulFull:
   ; Input in hl and de, result in cahl
   ld   c,h
   ld   a,l
   ld   hl,0   ;11
   ld   b,16   ;7
__MulFullNext:
   add   hl,hl   ;11
   rla      ;4
   rl   c   ;8
   jr   nc,__MulFullSkip   ;12/7
   add   hl,de   ;11
   adc   a,0   ;7
   jr   nc,__MulFullSkip
   inc   c
__MulFullSkip:
   djnz   __MulFullNext
   ret
__MulFullEnd:

Code: Select All | Copy To Clipboard

1 2 3 4 5 6 7 8 9 10 11 12 13 14

// Sqrt a temp = high byte of a a <<= 8 b = 0 repeat for every 2 bits in a test = b << 8 + 0x40 b <<= 1 if (temp >= test) temp -= test set low bit of b temp += high 2 bits of a a <<= 2 return b

If a is 4 bytes, then b and temp are 2 bytes, and you loop 16 times.

Spoiler for code:

stole my own routine from axe (and modified it)
p_Sqrt88:
   ; input in hlde, result in de
   ld   b,16
   ld   a,h
   ld   c,l
   push   de ; ld ixh,d
   pop   ix ; ld ixl,e
   ld   de,0
   ld   h,d
   ld   l,e
__Sqrt88Loop:
   sub   $40
   sbc   hl,de
   jr   nc,__Sqrt88Skip
   add   a,$40
   adc   hl,de
__Sqrt88Skip:
   ccf
   rl   e
   rl   d
   add   ix,ix
   rl   c
   rla
   adc   hl,hl
   add   ix,ix
   rl   c
   rla
   adc    hl,hl
   djnz   __Sqrt88Loop
   ret
__Sqrt88End:

Kool. Thanks. But, shouldn't the first one have two inputs?

Of course, hl and de, isn't that what I said

[FloppusMaximus]

My first multiplication routine takes 2746 - 4570 cycles, the second takes 1680 - 2880 cycles.

Oh boy, optimization time

The best I have so far is somewhere around 1800 cycles average (I'm too lazy to work out the exact probabilities at the moment, and not counting memory delays) using a squaring table and undocumented IX instructions.  Input is BDE and CHL, output is BCDEAL.  This routine works by expanding the formula 2xy = x²+y²-|x-y|², summed over each of the 9 pairs of bytes in the input.

(I'm not saying this is practical - unless you really have thousands of 24-bit multiplications to perform, you don't need this kind of speed.  This is just for fun.)

Code: Select All | Copy To Clipboard

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60 61 62 63 64 65 66 67 68 69 70 71 72 73 74 75 76 77 78 79 80 81 82 83 84 85 86 87 88 89 90 91 92 93 94 95 96 97 98 99 100 101 102 103 104 105 106 107 108 109 110 111 112 113 114 115 116 117 118 119 120 121 122 123 124 125 126 127 128 129 130 131 132 133 134 135 136 137 138 139 140 141 142 143 144 145 146 147 148 149 150 151 152 153 154 155 156 157 158 159 160 161 162 163 164 165 166 167 168 169 170 171 172 173 174 175 176 177 178 179 180 181 182 183 184 185 186 187 188 189 190 191 192 193 194 195 196 197 198 199 200 201 202 203 204 205 206 207 208 209 210 211 212 213 214 215 216 217 218 219 220 221 222

SUBFIRST .macro src1, src2, hdest, ldest exx ld a, src1 sub src2 jr nc, $ + 4 neg exx ld l, a ld a, ldest sub (hl) ld ldest, a inc h ld a, hdest sbc a, (hl) ld hdest, a   .endm SUBNEXT .macro src1, src2, hdest, ldest dec h ex af, af' exx ld a, src1 sub src2 jr nc, $ + 4 neg exx ld l, a ex af, af' ld a, ldest sbc a, (hl) ld ldest, a inc h ld a, hdest sbc a, (hl) ld hdest, a   .endm BDE_times_CHL_sqrdiff_v3: ld a, d exx ld h, high(sqrtab) ld l, a ld e, (hl) inc h ld d, (hl) ; DE = d² exx ld a, b exx ld l, a ld b, (hl) dec h ld c, (hl) ; BC = b² exx ld a, e exx ld l, a ld a, (hl) inc h ld h, (hl) ld l, a ; HL = e² call BC_DE_HL_times_10101 push bc push hl   push de    exx    ld a, h    exx    ld h, high(sqrtab)    ld l, a    ld e, (hl)    inc h    ld d, (hl) ; DE = h²    exx    ld a, c    exx    ld l, a    ld b, (hl)    dec h    ld c, (hl) ; BC = c²    exx    ld a, l    exx    ld l, a    ld a, (hl)    inc h    ld h, (hl)    ld l, a ; HL = l²    call BC_DE_HL_times_10101    pop ix   add ix, de   pop de adc hl, de ex de, hl pop hl adc hl, bc ld b, h ld c, l ; BCDEIX = total push af ld h, high(sqrtab) SUBFIRST e, l, ixh, ixl SUBNEXT  d, h, d, e SUBNEXT  b, c, b, c jp nc, BDE_times_CHL_sqrdiff_v3_nc1 pop af ccf push af BDE_times_CHL_sqrdiff_v3_nc1: inc b dec h SUBFIRST e, h, e, ixh SUBNEXT  d, c, c, d jr nc, BDE_times_CHL_sqrdiff_v3_nc2 dec b jp nz, BDE_times_CHL_sqrdiff_v3_nc2 pop af ccf push af BDE_times_CHL_sqrdiff_v3_nc2: dec h SUBFIRST d, l, e, ixh SUBNEXT  b, h, c, d jr nc, BDE_times_CHL_sqrdiff_v3_nc3 dec b jp nz, BDE_times_CHL_sqrdiff_v3_nc3 pop af ccf push af BDE_times_CHL_sqrdiff_v3_nc3: inc c dec h SUBFIRST b, l, d, e jr nc, BDE_times_CHL_sqrdiff_v3_nc4 dec c jp nz, BDE_times_CHL_sqrdiff_v3_nc4 dec b jp nz, BDE_times_CHL_sqrdiff_v3_nc4 pop af ccf push af BDE_times_CHL_sqrdiff_v3_nc4: dec h SUBFIRST e, c, d, e pop hl jr nc, BDE_times_CHL_sqrdiff_v3_nc5 dec c jp nz, BDE_times_CHL_sqrdiff_v3_nc5 dec b jp nz, BDE_times_CHL_sqrdiff_v3_nc5 inc l BDE_times_CHL_sqrdiff_v3_nc5: dec b dec c rr l rr b rr c rr d rr e ld a, ixl ld l, a ld a, ixh rra rr l ret BC_DE_HL_times_10101: push bc ld a, h ex af, af' sub a ld c, a ld b, l add hl, bc adc a, a ld b, e add hl, bc adc a, c ; AHL = [ L+H+E L ] pop bc push hl push bc   ld c, a   ld b, 0   ex af, af'   ld h, a   add hl, bc ; no way this can carry (initial HL is a square)   ld c, a   ld b, e   sub a   add hl, bc   adc a, a ; AHL(SP+2) = [ H+E L+H L+H+E L ]   add hl, de   adc a, 0 ; AHL(SP+2) = [ H+E+D L+H+E L+H+E L ]   pop bc add hl, bc adc a, 0 ; AHL(SP) = [ H+E+D+B L+H+E+C L+H+E L ] ld e, d ld d, c add hl, de adc a, b jr nc, BC_DE_HL_times_10101_nc1 inc b ; BAHL(SP) = [ B B H+E+D+C+B L+H+E+D+C L+H+E L ] BC_DE_HL_times_10101_nc1: add a, e jr nc, BC_DE_HL_times_10101_nc2 inc b ; BAHL(SP) = [ B D+B H+E+D+C+B L+H+E+D+C L+H+E L ] BC_DE_HL_times_10101_nc2: pop de add a, c ld c, a ret nc inc b ; BCHLDE = [ B D+C+B H+E+D+C+B L+H+E+D+C L+H+E L ] ret

To get back to the topic somewhat, ACagliano, it sounds like you're more interested in squaring than in general multiplication.  Squaring can be considerably faster, especially if you use a lookup table (e.g., my best 16-bit squaring routine is around 170 cycles, versus around 800 for general multiplication.)

[cerzus69]

I do have a 24-bit floating-point multiplication routine

saved 2 bytes, 1149 cycles saved by using iy too (and in a way compatible with TIOS, imagine that)

Code: Select All | Copy To Clipboard

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23

; hldebc = hlc * bde ld (iy+asm_Flag1),b xor a ld ix,0 ld b,24 Loop: add ix,ix rla rl c adc hl,hl jr nc,Next add ix,de adc a,(iy+asm_Flag1) rl c jr nc,Next inc hl Next: djnz Loop ld e,a ld d,c push ix ; ld c,ixl pop bc ; ld b,ixh

Jacobly, are you sure this works because I've been going through the code and it seems to me like the second 'rl c' should instead be 'add carry flag to c'. 2 times 'rl c' per loop seems wrong to me. Could you explain please? Because I've tried it as well in wabbitemu, taking the different input in account, and it is still not doing the right thing.

[ACagliano]

Yeah, all I need is 16-bit subtraction (which 'sub' supports, I think), 16-bit squaring, 32-bit addition, then 32-bit square rooting (or will I need to go up to 40-bit?).

[Xeda112358]

16-bit subtraction

Code: Select All | Copy To Clipboard

1 2 3

or a     ;to make sure the c flag is reset. Not always necessary if you know the c flag will be reset sbc hl,bc  ;you can do sbc hl,de also.

32-bit addition (you mean two 32-bit inputs?)

Code: Select All | Copy To Clipboard

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18

;Inputs: ;     HLBC is one of the 32-bit inputs ;     DE points to the other 32-bit input in RAM ;Outputs: ;     HLBC is the 32-bit result ;     DE is incremented 3 times ;     A=H ;     c flag is set if there is an overflow      ld a,(de) \ inc de      add a,c \ ld c,a      ld a,(de) \ inc de      adc a,b \ ld b,a      ld a,(de) \ inc de      adc a,l \ ld l,a      ld a,(de)      adc a,h \ ld h,a      ret

Squaring and square rooting... I will think on it

Also, I am working on a mini math library that will include RAM based math (so all the values will be in RAM). It seems like a few of these commands will need to rely on some memory. If they do, I suggest using the OP registers (11 bytes of RAM each).

[jacobly]

I do have a 24-bit floating-point multiplication routine

saved 2 bytes, 1149 cycles saved by using iy too (and in a way compatible with TIOS, imagine that)

Code: Select All | Copy To Clipboard

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23

; hldebc = hlc * bde ld (iy+asm_Flag1),b xor a ld ix,0 ld b,24 Loop: add ix,ix rla rl c adc hl,hl jr nc,Next add ix,de adc a,(iy+asm_Flag1) rl c jr nc,Next inc hl Next: djnz Loop ld e,a ld d,c push ix ; ld c,ixl pop bc ; ld b,ixh

Jacobly, are you sure this works because I've been going through the code and it seems to me like the second 'rl c' should instead be 'add carry flag to c'. 2 times 'rl c' per loop seems wrong to me. Could you explain please? Because I've tried it as well in wabbitemu, taking the different input in account, and it is still not doing the right thing.

That's strange. My test program must not have been working right, because when I went back and changed it a bit, it suddenly started telling me that the second routine doesn't work.
Anyway, my new test program seems to agree with this change.

Code: Select All | Copy To Clipboard

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24

; hldebc = hlc * bde ld (iy+asm_Flag1),b xor a ld ix,0 ld b,24 Loop: add ix,ix rla rl c adc hl,hl jr nc,Next add ix,de adc a,(iy+asm_Flag1) jr nc,Next inc c jr nz,Next inc hl Next: djnz Loop ld e,a ld d,c push ix ; ld c,ixl pop bc ; ld b,ixh

[cerzus69]

That's strange. My test program must not have been working right, because when I went back and changed it a bit, it suddenly started telling me that the second routine doesn't work.
Anyway, my new test program seems to agree with this change.

Code: Select All | Copy To Clipboard

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24

; hldebc = hlc * bde ld (iy+asm_Flag1),b xor a ld ix,0 ld b,24 Loop: add ix,ix rla rl c adc hl,hl jr nc,Next add ix,de adc a,(iy+asm_Flag1) jr nc,Next inc c jr nz,Next inc hl Next: djnz Loop ld e,a ld d,c push ix ; ld c,ixl pop bc ; ld b,ixh

Cool, thanks a lot, indeed it works now!