home *** CD-ROM | disk | FTP | other *** search
/ Atari FTP / ATARI_FTP_0693.zip / ATARI_FTP_0693 / Mint / mntlib32.zoo / _divsf3.cpp < prev    next >
C/C++ Source or Header  |  1993-06-17  |  6KB  |  272 lines

  1. |
  2. | single floating point divide routine
  3. |
  4. #ifndef __M68881__
  5.     .text
  6.     .even
  7.     .globl    __divsf3, ___divsf3
  8. #ifdef    ERROR_CHECK
  9. #include "errbase.h"
  10.     .globl    __infinitysf
  11.  
  12. LC0:
  13.     .ascii "floating point division by 0\12\15\0"
  14.     .even
  15. #endif    ERROR_CHECK
  16.  
  17. __divsf3:
  18. ___divsf3:
  19.  
  20. #ifdef    ERROR_CHECK
  21.     tstl    sp@(8)            | check if divisor is 0
  22.     bne    no_exception
  23.  
  24.     pea    pc@(LC0)
  25.     pea    Stderr
  26.     jbsr    _fprintf
  27.     addql    #8,sp
  28.  
  29.     moveq    #Erange,d0        | set _errno to ERANGE
  30.     Emove    d0,Errno
  31.     movel    __infinitysf,d0        | return signed infinity
  32.     btst    #31,sp@(4)        | transfer sign of dividend
  33.     beq    clear            | (mjr++)
  34.     bset    #31,d0            |
  35. clear:                    |
  36.     rts
  37.  
  38. no_exception:
  39. #endif    ERROR_CHECK
  40.  
  41. #ifndef sfp004
  42. |
  43. | written by Kai-Uwe Bloem (I5110401@dbstu1.bitnet).
  44. | Based on a 80x86 floating point packet from comp.os.minix, written by P.Housel
  45. | Revision 1.2.4 michal 05-93 (ntomczak@vm.ucs.ualberta.ca)
  46. |   + resynchro with errno codes
  47. |   + code smoothing
  48. |   + removed extra code 'retz:' and 'divz:'; these cases are
  49. |     going back through 'retzero:' and 'retinf:'
  50. |
  51. | patched by Olaf Flebbe (flebbe@tat.physik.uni-tuebingen.de)
  52. |
  53. | Revision 1.2.3 olaf 5-93
  54. |   + correct sign bug for retinf and retzero
  55. |
  56. | Revision 1.2.2 olaf 12-92
  57. |   + added support for NaN and Infinites
  58. |   + added support for -0
  59. |
  60. | Revision 1.2.1 olaf 11-92
  61. |   + prevent the tie rounding case if dividing is not exact.
  62. |      > paranoia now says: "Division appears to round correctly"
  63. |      ** requires _normsf Version 1.4.2 or later
  64. |
  65. | Revision 1.2, kub 01-90 :
  66. | added support for denormalized numbers
  67. |
  68. | Revision 1.1, kub 12-89 :
  69. | Created single float version for 68000
  70. |
  71. | Revision 1.0:
  72. | original 8088 code from P.S.Housel for double floats
  73.  
  74. BIAS4    =    0x7F-1
  75.  
  76.     lea    sp@(4),a0    | pointer to parameters u and v
  77.     moveml    d2-d5,sp@-    | save registers
  78.     moveml    a0@,d4/d5    | d4 = u, d5 = v
  79.  
  80.     movel    #0x7fffff,d3
  81.     movel    d4,d0        | d0 = u.exp
  82.     andl    d3,d4        | remove exponent from u.mantissa
  83.     swap    d0
  84.     movew    d0,d2        | d2 = u.sign
  85.  
  86.     movel    d5,d1        | d1 = v.exp
  87.     andl    d3,d5        | remove exponent from v.mantissa
  88.     swap    d1
  89.     eorw    d1,d2        | d2 = u.sign ^ v.sign (in bit 15)
  90.  
  91.     moveq    #15,d3
  92.     bclr    d3,d0        | kill sign bit
  93.     bclr    d3,d1        | kill sign bit
  94.     lsrw    #7,d0
  95.     lsrw    #7,d1
  96. |
  97. |
  98. |
  99.     moveq    #-1,d3
  100.     cmpb    d3,d0        | comparison with #0xff
  101.     beq    0f        | u == NaN || u== Inf
  102.     cmpb    d3,d1
  103.     beq    1f        | v == NaN || v == Inf
  104.     tstb    d0
  105.     bne    3f        | u not zero nor denorm
  106.     tstl    d4
  107.     beq    2f        | 0/ ?
  108.  
  109. 3:    tstw    d1
  110.     bne    nospec
  111.  
  112.     tstl    d5
  113.     bne    nospec
  114.     bra    retinf        | x/0 -> +/- Inf
  115.  
  116. 0:    tstl    d4        | u == NaN ?
  117.     bne    retnan        | NaN/ x
  118.     cmpb    d3,d1
  119.     beq    retnan        | Inf/Inf or Inf/NaN
  120. |    bra    retinf        | Inf/x | x != Inf && x != NaN
  121. |
  122. |    Return Infinity with correct sign
  123. |
  124. retinf:    movel    #0xff000000,d0
  125.     lslw    #1,d2
  126.     roxrl   #1,d0        | shift in high bit as given by d2
  127. return:    moveml    sp@+,d2-d5
  128.     rts
  129.  
  130. 1:    tstl    d5
  131.     bne    retnan        | x/NaN
  132. |    bra    retzero        | x/Inf -> +/- 0
  133. |
  134. |    Return correct signed zero
  135. |
  136. retzero:moveq    #0,d0        | zero destination
  137.     lslw    #1,d2        | set X bit accordingly
  138.     roxrl    #1,d0
  139.     bra    return
  140.  
  141. 2:    tstw    d1
  142.     bne    retzero        | 0/x ->+/- 0
  143.     tstl    d4
  144.     bne    retzero        | 0/x
  145. |    bra    retnan        | 0/0
  146. |
  147. |    Return NaN
  148. |
  149. retnan: movel    d3,d0        | d3 contains 0xffffffff
  150.     lsrl    #1,d0
  151.     bra    return
  152. |
  153. |    End of special handling
  154. |
  155. nospec:    moveq    #23,d3
  156.     bset    d3,d4        | restore implied leading "1"
  157.     tstw    d0        | check for zero exponent - no leading "1"
  158.     bne    1f
  159.     bclr    d3,d4        | remove it
  160.     addw    #1,d0        | "normalize" exponent
  161. 1:
  162.     tstl    d4
  163.     beq    retzero        | dividing zero
  164.  
  165.     bset    d3,d5        | restore implied leading "1"
  166.     tstw    d1        | check for zero exponent - no leading "1"
  167.     bne    1f
  168.     bclr    d3,d5        | remove it
  169.     addw    #1,d1        | "normalize" exponent
  170. 1:
  171. # ifndef ERROR_CHECK
  172.     tstl    d5
  173.     beq    retinf        | divide by zero
  174. # endif    ERROR_CHECK
  175.  
  176.     subw    d1,d0        | subtract exponents,
  177.     addw    #BIAS4-8+1,d0    | add bias back in, account for shift
  178.     addw    #34,d0        | add loop offset, +2 for extra rounding bits
  179.                 | for denormalized numbers (2 implied by dbra)
  180.     movew    #27,d1        | bit number for "implied" pos (+4 for rounding)
  181.     moveql    #-1,d3        | zero quotient (for speed a one''s complement)
  182.     subl    d5,d4        | initial subtraction, u = u - v
  183. 2:
  184.     btst    d1,d3        | divide until 1 in implied position
  185.     beq    5f
  186.  
  187.     addl    d4,d4
  188.     bcs    4f        | if carry is set, add, else subtract
  189.  
  190.     addxl    d3,d3        | shift quotient and set bit zero
  191.     subl    d5,d4        | subtract, u = u - v
  192.     dbra    d0,2b        | give up if result is denormalized
  193.     bra    5f
  194. 4:
  195.     addxl    d3,d3        | shift quotient and clear bit zero
  196.     addl    d5,d4        | add (restore), u = u + v
  197.     dbra    d0,2b        | give up if result is denormalized
  198. 5:    subqw    #2,d0        | remove rounding offset for denormalized nums
  199.     notl    d3        | invert quotient to get it right
  200.  
  201.     clrl    d1        | zero rounding bits
  202.     tstl     d4        | check for exact result
  203.     beq    1f
  204.     moveql    #-1,d1        | prevent tie case
  205. 1:    movel    d3,d4        | save quotient mantissa
  206.     jmp    norm_sf        | (registers on stack removed by norm_sf)
  207.  
  208. #if 0
  209. | this  is dead code right now - retzero and retinf used (mj)
  210. # ifndef ERROR_CHECK
  211. retz:    moveq    #0,d0        | zero destination
  212.     moveml    sp@+,d2-d5
  213.     rts            | no normalization needed
  214.  
  215. divz:    movel    __infinitysf,d0    | return infinity value
  216.     moveml    sp@+,d2-d5    | should really cause trap ?!?
  217.     btst    #31,sp@(4)    | transfer sign of dividend
  218.     beq    clear        | (mjr++)
  219.     bset    #31,d0        |
  220.     rts            |
  221. clear:                |
  222.     bclr    #31,d0        |
  223.     rts
  224. # endif    ERROR_CHECK
  225. #endif
  226.  
  227. #else
  228.  
  229. | single precision floating point stuff for Atari-gcc using the SFP004
  230. | or compatible boards with a memory mapped 68881
  231. | developed with gas
  232. |
  233. |  single floating point divide routine
  234. |
  235. | M. Ritzert (mjr at dmzrzu71)
  236. |            (ritzert@dfg.dbp.de)
  237. | 4.10.1990
  238. |
  239. | +_infinitysf returned instead of a NAN
  240. | the DOMAIN exception is not supported yet. In case of an exception
  241. | _errno is always set to ERANGE
  242.  
  243. | addresses of the 68881 data port. This choice is fastest when much data is
  244. | transferred between the two processors.
  245.  
  246. comm =     -6
  247. resp =    -16
  248. zahl =      0
  249.  
  250. | waiting loop ...
  251. |
  252. | wait:
  253. | ww:    cmpiw    #0x8900,a0@(resp)
  254. |     beq    ww
  255. | is coded directly by
  256. |    .long    0x0c688900, 0xfff067f8
  257.  
  258.     lea    0xfffffa50:w,a0
  259.     movew    #0x4400,a0@(comm)    | load first argument to fp0
  260.     cmpiw    #0x8900,a0@(resp)    | check
  261.     movel    sp@(4),a0@
  262.     movew    #0x4424,a0@(comm)
  263.     .long    0x0c688900, 0xfff067f8
  264.     movel    sp@(8),a0@
  265.     movew    #0x6400,a0@(comm)    | result to d0
  266.     .long    0x0c688900, 0xfff067f8
  267.     movel    a0@,d0            | REMARK: 0/0 returns a NAN
  268.     rts                | if ERROR_CHECK is disabled
  269.  
  270. #endif    sfp004
  271. #endif /* !__M68881__ */
  272.