home *** CD-ROM | disk | FTP | other *** search
/ GEMini Atari / GEMini_Atari_CD-ROM_Walnut_Creek_December_1993.iso / zip / gnu / mntlib16.lzh / MNTLIB16 / _MULSF3.S < prev    next >
Text File  |  1993-07-29  |  3KB  |  115 lines
  1. | single floating point multiplication routine
  2. |
  3. | written by Kai-Uwe Bloem (I5110401@dbstu1.bitnet).
  4. | Based on a 80x86 floating point packet from comp.os.minix, written by P.Housel
  5. |
  6. |
  7. | Revision 1.2, kub 01-90 :
  8. | added support for denormalized numbers
  9. |
  10. | Revision 1.1, kub 12-89 :
  11. | Created single float version for 68000. Code could be speed up by having
  12. | the accumulator in the 68000 register set ...
  13. |
  14. | Revision 1.0:
  15. | original 8088 code from P.S.Housel for double floats
  16.  
  17. BIAS4    =    0x7F-1
  18.  
  19.     .text
  20.     .even
  21.     .globl    __mulsf3, ___mulsf3
  22.  
  23. __mulsf3:
  24. ___mulsf3:
  25.     lea    sp@(4),a0
  26.     moveml    d2-d5,sp@-
  27.     moveml    a0@,d4/d5    | d4 = v, d5 = u
  28.     subw    #8,sp        | multiplication accumulator
  29.  
  30.     movel    d5,d0        | d0 = u.exp
  31.     swap    d0
  32.     movew    d0,d2        | d2 = u.sign
  33.     lsrw    #7,d0
  34.     andw    #0xff,d0    | kill sign bit
  35.  
  36.     movel    d4,d1        | d1 = v.exp
  37.     swap    d1
  38.     eorw    d1,d2        | d2 = u.sign ^ v.sign (in bit 31)
  39.     lsrw    #7,d1
  40.     andw    #0xff,d1    | kill sign bit
  41.  
  42.     andl    #0x7fffff,d5    | remove exponent from u.mantissa
  43.     tstw    d0        | check for zero exponent - no leading "1"
  44.     beq    0f
  45.     orl    #0x800000,d5    | restore implied leading "1"
  46.     bra    1f
  47. 0:    addw    #1,d0        | "normalize" exponent
  48. 1:    tstl    d5
  49.     beq    retz        | multiplying zero
  50.  
  51.     andl    #0x7fffff,d4    | remove exponent from v.mantissa
  52.     tstw    d1        | check for zero exponent - no leading "1"
  53.     beq    0f
  54.     orl    #0x800000,d4    | restore implied leading "1"
  55.     bra    1f
  56. 0:    addw    #1,d1        | "normalize" exponent
  57. 1:    tstl    d4
  58.     beq    retz        | multiply by zero
  59.  
  60.     addw    d1,d0        | add exponents,
  61.     subw    #BIAS4+16-8,d0    | remove excess bias, acnt for repositioning
  62.  
  63.     clrl    sp@        | initialize 64-bit product to zero
  64.     clrl    sp@(4)
  65.  
  66. | see Knuth, Seminumerical Algorithms, section 4.3. algorithm M
  67.  
  68.     movew    d4,d3
  69.     mulu    d5,d3        | mulitply with bigit from multiplier
  70.     movel    d3,sp@(4)    | store into result
  71.  
  72.     movel    d4,d3
  73.     swap    d3
  74.     mulu    d5,d3
  75.     addl    d3,sp@(2)    | add to result
  76.  
  77.     swap    d5        | [TOP 8 BITS SHOULD BE ZERO !]
  78.  
  79.     movew    d4,d3
  80.     mulu    d5,d3        | mulitply with bigit from multiplier
  81.     addl    d3,sp@(2)    | store into result (no carry can occur here)
  82.  
  83.     movel    d4,d3
  84.     swap    d3
  85.     mulu    d5,d3
  86.     addl    d3,sp@        | add to result
  87.                 | [TOP 16 BITS SHOULD BE ZERO !]
  88.     moveml    sp@(2),d4-d5    | get the 48 valid mantissa bits
  89.     clrw    d5        | (pad to 64)
  90. 2:
  91.     cmpl    #0x0000ffff,d4    | multiply (shift) until
  92.     bhi    3f        |  1 in upper 16 result bits
  93.     cmpw    #9,d0        | give up for denormalized numbers
  94.     ble    3f
  95.     swap    d4        | (we're getting here only when multiplying
  96.     swap    d5        |  with a denormalized number; there's an
  97.     movew    d5,d4        |  eventual loss of 4 bits in the rounding
  98.     clrw    d5        |  byte -- what a pity 8-)
  99.     subw    #16,d0        | decrement exponent
  100.     bra    2b
  101. 3:
  102.     movel    d5,d1        | get rounding bits
  103.     roll    #8,d1
  104.     movel    d1,d3        | see if sticky bit should be set
  105.     andl    #0xffffff00,d3
  106.     beq    4f
  107.     orb    #1,d1        | set "sticky bit" if any low-order set
  108. 4:    addw    #8,sp        | remove accumulator from stack
  109.     jmp    norm_sf        | (result in d4)
  110.  
  111. retz:    clrl    d0        | save zero as result
  112.     addw    #8,sp
  113.     moveml    sp@+,d2-d5
  114.     rts            | no normalizing neccessary
  115.