arch/m68k-all/m680x0/fpsp/satan.sa

   1 *       $NetBSD: satan.sa,v 1.3 1994/10/26 07:49:31 cgd Exp $
   2
   3 *       MOTOROLA MICROPROCESSOR & MEMORY TECHNOLOGY GROUP
   4 *       M68000 Hi-Performance Microprocessor Division
   5 *       M68040 Software Package
   6 *
   7 *       M68040 Software Package Copyright (c) 1993, 1994 Motorola Inc.
   8 *       All rights reserved.
   9 *
  10 *       THE SOFTWARE is provided on an "AS IS" basis and without warranty.
  11 *       To the maximum extent permitted by applicable law,
  12 *       MOTOROLA DISCLAIMS ALL WARRANTIES WHETHER EXPRESS OR IMPLIED,
  13 *       INCLUDING IMPLIED WARRANTIES OF MERCHANTABILITY OR FITNESS FOR A
  14 *       PARTICULAR PURPOSE and any warranty against infringement with
  15 *       regard to the SOFTWARE (INCLUDING ANY MODIFIED VERSIONS THEREOF)
  16 *       and any accompanying written materials.
  17 *
  18 *       To the maximum extent permitted by applicable law,
  19 *       IN NO EVENT SHALL MOTOROLA BE LIABLE FOR ANY DAMAGES WHATSOEVER
  20 *       (INCLUDING WITHOUT LIMITATION, DAMAGES FOR LOSS OF BUSINESS
  21 *       PROFITS, BUSINESS INTERRUPTION, LOSS OF BUSINESS INFORMATION, OR
  22 *       OTHER PECUNIARY LOSS) ARISING OF THE USE OR INABILITY TO USE THE
  23 *       SOFTWARE.  Motorola assumes no responsibility for the maintenance
  24 *       and support of the SOFTWARE.
  25 *
  26 *       You are hereby granted a copyright license to use, modify, and
  27 *       distribute the SOFTWARE so long as this entire notice is retained
  28 *       without alteration in any modified and/or redistributed versions,
  29 *       and that such modified versions are clearly identified as such.
  30 *       No licenses are granted by implication, estoppel or otherwise
  31 *       under any patents or trademarks of Motorola, Inc.
  32
  33 *
  34 *       satan.sa 3.3 12/19/90
  35 *
  36 *       The entry point satan computes the arctagent of an
  37 *       input value. satand does the same except the input value is a
  38 *       denormalized number.
  39 *
  40 *       Input: Double-extended value in memory location pointed to by address
  41 *               register a0.
  42 *
  43 *       Output: Arctan(X) returned in floating-point register Fp0.
  44 *
  45 *       Accuracy and Monotonicity: The returned result is within 2 ulps in
  46 *               64 significant bit, i.e. within 0.5001 ulp to 53 bits if the
  47 *               result is subsequently rounded to double precision. The
  48 *               result is provably monotonic in double precision.
  49 *
  50 *       Speed: The program satan takes approximately 160 cycles for input
  51 *               argument X such that 1/16 < |X| < 16. For the other arguments,
  52 *               the program will run no worse than 10% slower.
  53 *
  54 *       Algorithm:
  55 *       Step 1. If |X| >= 16 or |X| < 1/16, go to Step 5.
  56 *
  57 *       Step 2. Let X = sgn * 2**k * 1.xxxxxxxx...x. Note that k = -4, -3,..., or 3.
  58 *               Define F = sgn * 2**k * 1.xxxx1, i.e. the first 5 significant bits
  59 *               of X with a bit-1 attached at the 6-th bit position. Define u
  60 *               to be u = (X-F) / (1 + X*F).
  61 *
  62 *       Step 3. Approximate arctan(u) by a polynomial poly.
  63 *
  64 *       Step 4. Return arctan(F) + poly, arctan(F) is fetched from a table of values
  65 *               calculated beforehand. Exit.
  66 *
  67 *       Step 5. If |X| >= 16, go to Step 7.
  68 *
  69 *       Step 6. Approximate arctan(X) by an odd polynomial in X. Exit.
  70 *
  71 *       Step 7. Define X' = -1/X. Approximate arctan(X') by an odd polynomial in X'.
  72 *               Arctan(X) = sign(X)*Pi/2 + arctan(X'). Exit.
  73 *
  74
  75 satan   IDNT    2,1 Motorola 040 Floating Point Software Package
  76
  77         section 8
  78
  79         include fpsp.h
  80
  81 BOUNDS1 DC.L $3FFB8000,$4002FFFF
  82
  83 ONE     DC.L $3F800000
  84
  85         DC.L $00000000
  86
  87 ATANA3  DC.L $BFF6687E,$314987D8
  88 ATANA2  DC.L $4002AC69,$34A26DB3
  89
  90 ATANA1  DC.L $BFC2476F,$4E1DA28E
  91 ATANB6  DC.L $3FB34444,$7F876989
  92
  93 ATANB5  DC.L $BFB744EE,$7FAF45DB
  94 ATANB4  DC.L $3FBC71C6,$46940220
  95
  96 ATANB3  DC.L $BFC24924,$921872F9
  97 ATANB2  DC.L $3FC99999,$99998FA9
  98
  99 ATANB1  DC.L $BFD55555,$55555555
 100 ATANC5  DC.L $BFB70BF3,$98539E6A
 101
 102 ATANC4  DC.L $3FBC7187,$962D1D7D
 103 ATANC3  DC.L $BFC24924,$827107B8
 104
 105 ATANC2  DC.L $3FC99999,$9996263E
 106 ATANC1  DC.L $BFD55555,$55555536
 107
 108 PPIBY2  DC.L $3FFF0000,$C90FDAA2,$2168C235,$00000000
 109 NPIBY2  DC.L $BFFF0000,$C90FDAA2,$2168C235,$00000000
 110 PTINY   DC.L $00010000,$80000000,$00000000,$00000000
 111 NTINY   DC.L $80010000,$80000000,$00000000,$00000000
 112
 113 ATANTBL:
 114         DC.L    $3FFB0000,$83D152C5,$060B7A51,$00000000
 115         DC.L    $3FFB0000,$8BC85445,$65498B8B,$00000000
 116         DC.L    $3FFB0000,$93BE4060,$17626B0D,$00000000
 117         DC.L    $3FFB0000,$9BB3078D,$35AEC202,$00000000
 118         DC.L    $3FFB0000,$A3A69A52,$5DDCE7DE,$00000000
 119         DC.L    $3FFB0000,$AB98E943,$62765619,$00000000
 120         DC.L    $3FFB0000,$B389E502,$F9C59862,$00000000
 121         DC.L    $3FFB0000,$BB797E43,$6B09E6FB,$00000000
 122         DC.L    $3FFB0000,$C367A5C7,$39E5F446,$00000000
 123         DC.L    $3FFB0000,$CB544C61,$CFF7D5C6,$00000000
 124         DC.L    $3FFB0000,$D33F62F8,$2488533E,$00000000
 125         DC.L    $3FFB0000,$DB28DA81,$62404C77,$00000000
 126         DC.L    $3FFB0000,$E310A407,$8AD34F18,$00000000
 127         DC.L    $3FFB0000,$EAF6B0A8,$188EE1EB,$00000000
 128         DC.L    $3FFB0000,$F2DAF194,$9DBE79D5,$00000000
 129         DC.L    $3FFB0000,$FABD5813,$61D47E3E,$00000000
 130         DC.L    $3FFC0000,$8346AC21,$0959ECC4,$00000000
 131         DC.L    $3FFC0000,$8B232A08,$304282D8,$00000000
 132         DC.L    $3FFC0000,$92FB70B8,$D29AE2F9,$00000000
 133         DC.L    $3FFC0000,$9ACF476F,$5CCD1CB4,$00000000
 134         DC.L    $3FFC0000,$A29E7630,$4954F23F,$00000000
 135         DC.L    $3FFC0000,$AA68C5D0,$8AB85230,$00000000
 136         DC.L    $3FFC0000,$B22DFFFD,$9D539F83,$00000000
 137         DC.L    $3FFC0000,$B9EDEF45,$3E900EA5,$00000000
 138         DC.L    $3FFC0000,$C1A85F1C,$C75E3EA5,$00000000
 139         DC.L    $3FFC0000,$C95D1BE8,$28138DE6,$00000000
 140         DC.L    $3FFC0000,$D10BF300,$840D2DE4,$00000000
 141         DC.L    $3FFC0000,$D8B4B2BA,$6BC05E7A,$00000000
 142         DC.L    $3FFC0000,$E0572A6B,$B42335F6,$00000000
 143         DC.L    $3FFC0000,$E7F32A70,$EA9CAA8F,$00000000
 144         DC.L    $3FFC0000,$EF888432,$64ECEFAA,$00000000
 145         DC.L    $3FFC0000,$F7170A28,$ECC06666,$00000000
 146         DC.L    $3FFD0000,$812FD288,$332DAD32,$00000000
 147         DC.L    $3FFD0000,$88A8D1B1,$218E4D64,$00000000
 148         DC.L    $3FFD0000,$9012AB3F,$23E4AEE8,$00000000
 149         DC.L    $3FFD0000,$976CC3D4,$11E7F1B9,$00000000
 150         DC.L    $3FFD0000,$9EB68949,$3889A227,$00000000
 151         DC.L    $3FFD0000,$A5EF72C3,$4487361B,$00000000
 152         DC.L    $3FFD0000,$AD1700BA,$F07A7227,$00000000
 153         DC.L    $3FFD0000,$B42CBCFA,$FD37EFB7,$00000000
 154         DC.L    $3FFD0000,$BB303A94,$0BA80F89,$00000000
 155         DC.L    $3FFD0000,$C22115C6,$FCAEBBAF,$00000000
 156         DC.L    $3FFD0000,$C8FEF3E6,$86331221,$00000000
 157         DC.L    $3FFD0000,$CFC98330,$B4000C70,$00000000
 158         DC.L    $3FFD0000,$D6807AA1,$102C5BF9,$00000000
 159         DC.L    $3FFD0000,$DD2399BC,$31252AA3,$00000000
 160         DC.L    $3FFD0000,$E3B2A855,$6B8FC517,$00000000
 161         DC.L    $3FFD0000,$EA2D764F,$64315989,$00000000
 162         DC.L    $3FFD0000,$F3BF5BF8,$BAD1A21D,$00000000
 163         DC.L    $3FFE0000,$801CE39E,$0D205C9A,$00000000
 164         DC.L    $3FFE0000,$8630A2DA,$DA1ED066,$00000000
 165         DC.L    $3FFE0000,$8C1AD445,$F3E09B8C,$00000000
 166         DC.L    $3FFE0000,$91DB8F16,$64F350E2,$00000000
 167         DC.L    $3FFE0000,$97731420,$365E538C,$00000000
 168         DC.L    $3FFE0000,$9CE1C8E6,$A0B8CDBA,$00000000
 169         DC.L    $3FFE0000,$A22832DB,$CADAAE09,$00000000
 170         DC.L    $3FFE0000,$A746F2DD,$B7602294,$00000000
 171         DC.L    $3FFE0000,$AC3EC0FB,$997DD6A2,$00000000
 172         DC.L    $3FFE0000,$B110688A,$EBDC6F6A,$00000000
 173         DC.L    $3FFE0000,$B5BCC490,$59ECC4B0,$00000000
 174         DC.L    $3FFE0000,$BA44BC7D,$D470782F,$00000000
 175         DC.L    $3FFE0000,$BEA94144,$FD049AAC,$00000000
 176         DC.L    $3FFE0000,$C2EB4ABB,$661628B6,$00000000
 177         DC.L    $3FFE0000,$C70BD54C,$E602EE14,$00000000
 178         DC.L    $3FFE0000,$CD000549,$ADEC7159,$00000000
 179         DC.L    $3FFE0000,$D48457D2,$D8EA4EA3,$00000000
 180         DC.L    $3FFE0000,$DB948DA7,$12DECE3B,$00000000
 181         DC.L    $3FFE0000,$E23855F9,$69E8096A,$00000000
 182         DC.L    $3FFE0000,$E8771129,$C4353259,$00000000
 183         DC.L    $3FFE0000,$EE57C16E,$0D379C0D,$00000000
 184         DC.L    $3FFE0000,$F3E10211,$A87C3779,$00000000
 185         DC.L    $3FFE0000,$F919039D,$758B8D41,$00000000
 186         DC.L    $3FFE0000,$FE058B8F,$64935FB3,$00000000
 187         DC.L    $3FFF0000,$8155FB49,$7B685D04,$00000000
 188         DC.L    $3FFF0000,$83889E35,$49D108E1,$00000000
 189         DC.L    $3FFF0000,$859CFA76,$511D724B,$00000000
 190         DC.L    $3FFF0000,$87952ECF,$FF8131E7,$00000000
 191         DC.L    $3FFF0000,$89732FD1,$9557641B,$00000000
 192         DC.L    $3FFF0000,$8B38CAD1,$01932A35,$00000000
 193         DC.L    $3FFF0000,$8CE7A8D8,$301EE6B5,$00000000
 194         DC.L    $3FFF0000,$8F46A39E,$2EAE5281,$00000000
 195         DC.L    $3FFF0000,$922DA7D7,$91888487,$00000000
 196         DC.L    $3FFF0000,$94D19FCB,$DEDF5241,$00000000
 197         DC.L    $3FFF0000,$973AB944,$19D2A08B,$00000000
 198         DC.L    $3FFF0000,$996FF00E,$08E10B96,$00000000
 199         DC.L    $3FFF0000,$9B773F95,$12321DA7,$00000000
 200         DC.L    $3FFF0000,$9D55CC32,$0F935624,$00000000
 201         DC.L    $3FFF0000,$9F100575,$006CC571,$00000000
 202         DC.L    $3FFF0000,$A0A9C290,$D97CC06C,$00000000
 203         DC.L    $3FFF0000,$A22659EB,$EBC0630A,$00000000
 204         DC.L    $3FFF0000,$A388B4AF,$F6EF0EC9,$00000000
 205         DC.L    $3FFF0000,$A4D35F10,$61D292C4,$00000000
 206         DC.L    $3FFF0000,$A60895DC,$FBE3187E,$00000000
 207         DC.L    $3FFF0000,$A72A51DC,$7367BEAC,$00000000
 208         DC.L    $3FFF0000,$A83A5153,$0956168F,$00000000
 209         DC.L    $3FFF0000,$A93A2007,$7539546E,$00000000
 210         DC.L    $3FFF0000,$AA9E7245,$023B2605,$00000000
 211         DC.L    $3FFF0000,$AC4C84BA,$6FE4D58F,$00000000
 212         DC.L    $3FFF0000,$ADCE4A4A,$606B9712,$00000000
 213         DC.L    $3FFF0000,$AF2A2DCD,$8D263C9C,$00000000
 214         DC.L    $3FFF0000,$B0656F81,$F22265C7,$00000000
 215         DC.L    $3FFF0000,$B1846515,$0F71496A,$00000000
 216         DC.L    $3FFF0000,$B28AAA15,$6F9ADA35,$00000000
 217         DC.L    $3FFF0000,$B37B44FF,$3766B895,$00000000
 218         DC.L    $3FFF0000,$B458C3DC,$E9630433,$00000000
 219         DC.L    $3FFF0000,$B525529D,$562246BD,$00000000
 220         DC.L    $3FFF0000,$B5E2CCA9,$5F9D88CC,$00000000
 221         DC.L    $3FFF0000,$B692CADA,$7ACA1ADA,$00000000
 222         DC.L    $3FFF0000,$B736AEA7,$A6925838,$00000000
 223         DC.L    $3FFF0000,$B7CFAB28,$7E9F7B36,$00000000
 224         DC.L    $3FFF0000,$B85ECC66,$CB219835,$00000000
 225         DC.L    $3FFF0000,$B8E4FD5A,$20A593DA,$00000000
 226         DC.L    $3FFF0000,$B99F41F6,$4AFF9BB5,$00000000
 227         DC.L    $3FFF0000,$BA7F1E17,$842BBE7B,$00000000
 228         DC.L    $3FFF0000,$BB471285,$7637E17D,$00000000
 229         DC.L    $3FFF0000,$BBFABE8A,$4788DF6F,$00000000
 230         DC.L    $3FFF0000,$BC9D0FAD,$2B689D79,$00000000
 231         DC.L    $3FFF0000,$BD306A39,$471ECD86,$00000000
 232         DC.L    $3FFF0000,$BDB6C731,$856AF18A,$00000000
 233         DC.L    $3FFF0000,$BE31CAC5,$02E80D70,$00000000
 234         DC.L    $3FFF0000,$BEA2D55C,$E33194E2,$00000000
 235         DC.L    $3FFF0000,$BF0B10B7,$C03128F0,$00000000
 236         DC.L    $3FFF0000,$BF6B7A18,$DACB778D,$00000000
 237         DC.L    $3FFF0000,$BFC4EA46,$63FA18F6,$00000000
 238         DC.L    $3FFF0000,$C0181BDE,$8B89A454,$00000000
 239         DC.L    $3FFF0000,$C065B066,$CFBF6439,$00000000
 240         DC.L    $3FFF0000,$C0AE345F,$56340AE6,$00000000
 241         DC.L    $3FFF0000,$C0F22291,$9CB9E6A7,$00000000
 242
 243 X       equ     FP_SCR1
 244 XDCARE  equ     X+2
 245 XFRAC   equ     X+4
 246 XFRACLO equ     X+8
 247
 248 ATANF   equ     FP_SCR2
 249 ATANFHI equ     ATANF+4
 250 ATANFLO equ     ATANF+8
 251
 252
 253         xref    t_frcinx
 254         xref    t_extdnrm
 255
 256         xdef    satand
 257 satand:
 258 *--ENTRY POINT FOR ATAN(X) FOR DENORMALIZED ARGUMENT
 259
 260         bra             t_extdnrm
 261
 262         xdef    satan
 263 satan:
 264 *--ENTRY POINT FOR ATAN(X), HERE X IS FINITE, NON-ZERO, AND NOT NAN'S
 265
 266         FMOVE.X         (A0),FP0        ...LOAD INPUT
 267
 268         MOVE.L          (A0),D0
 269         MOVE.W          4(A0),D0
 270         FMOVE.X         FP0,X(a6)
 271         ANDI.L          #$7FFFFFFF,D0
 272
 273         CMPI.L          #$3FFB8000,D0           ...|X| >= 1/16?
 274         BGE.B           ATANOK1
 275         BRA.W           ATANSM
 276
 277 ATANOK1:
 278         CMPI.L          #$4002FFFF,D0           ...|X| < 16 ?
 279         BLE.B           ATANMAIN
 280         BRA.W           ATANBIG
 281
 282
 283 *--THE MOST LIKELY CASE, |X| IN [1/16, 16). WE USE TABLE TECHNIQUE
 284 *--THE IDEA IS ATAN(X) = ATAN(F) + ATAN( [X-F] / [1+XF] ).
 285 *--SO IF F IS CHOSEN TO BE CLOSE TO X AND ATAN(F) IS STORED IN
 286 *--A TABLE, ALL WE NEED IS TO APPROXIMATE ATAN(U) WHERE
 287 *--U = (X-F)/(1+XF) IS SMALL (REMEMBER F IS CLOSE TO X). IT IS
 288 *--TRUE THAT A DIVIDE IS NOW NEEDED, BUT THE APPROXIMATION FOR
 289 *--ATAN(U) IS A VERY SHORT POLYNOMIAL AND THE INDEXING TO
 290 *--FETCH F AND SAVING OF REGISTERS CAN BE ALL HIDED UNDER THE
 291 *--DIVIDE. IN THE END THIS METHOD IS MUCH FASTER THAN A TRADITIONAL
 292 *--ONE. NOTE ALSO THAT THE TRADITIONAL SCHEME THAT APPROXIMATE
 293 *--ATAN(X) DIRECTLY WILL NEED TO USE A RATIONAL APPROXIMATION
 294 *--(DIVISION NEEDED) ANYWAY BECAUSE A POLYNOMIAL APPROXIMATION
 295 *--WILL INVOLVE A VERY LONG POLYNOMIAL.
 296
 297 *--NOW WE SEE X AS +-2^K * 1.BBBBBBB....B <- 1. + 63 BITS
 298 *--WE CHOSE F TO BE +-2^K * 1.BBBB1
 299 *--THAT IS IT MATCHES THE EXPONENT AND FIRST 5 BITS OF X, THE
 300 *--SIXTH BITS IS SET TO BE 1. SINCE K = -4, -3, ..., 3, THERE
 301 *--ARE ONLY 8 TIMES 16 = 2^7 = 128 |F|'S. SINCE ATAN(-|F|) IS
 302 *-- -ATAN(|F|), WE NEED TO STORE ONLY ATAN(|F|).
 303
 304 ATANMAIN:
 305
 306         CLR.W           XDCARE(a6)              ...CLEAN UP X JUST IN CASE
 307         ANDI.L          #$F8000000,XFRAC(a6)    ...FIRST 5 BITS
 308         ORI.L           #$04000000,XFRAC(a6)    ...SET 6-TH BIT TO 1
 309         CLR.L           XFRACLO(a6)             ...LOCATION OF X IS NOW F
 310
 311         FMOVE.X         FP0,FP1                 ...FP1 IS X
 312         FMUL.X          X(a6),FP1               ...FP1 IS X*F, NOTE THAT X*F > 0
 313         FSUB.X          X(a6),FP0               ...FP0 IS X-F
 314         FADD.S          #:3F800000,FP1          ...FP1 IS 1 + X*F
 315         FDIV.X          FP1,FP0                 ...FP0 IS U = (X-F)/(1+X*F)
 316
 317 *--WHILE THE DIVISION IS TAKING ITS TIME, WE FETCH ATAN(|F|)
 318 *--CREATE ATAN(F) AND STORE IT IN ATANF, AND
 319 *--SAVE REGISTERS FP2.
 320
 321         MOVE.L          d2,-(a7)        ...SAVE d2 TEMPORARILY
 322         MOVE.L          d0,d2           ...THE EXPO AND 16 BITS OF X
 323         ANDI.L          #$00007800,d0   ...4 VARYING BITS OF F'S FRACTION
 324         ANDI.L          #$7FFF0000,d2   ...EXPONENT OF F
 325         SUBI.L          #$3FFB0000,d2   ...K+4
 326         ASR.L           #1,d2
 327         ADD.L           d2,d0           ...THE 7 BITS IDENTIFYING F
 328         ASR.L           #7,d0           ...INDEX INTO TBL OF ATAN(|F|)
 329         LEA             ATANTBL,a1
 330         ADDA.L          d0,a1           ...ADDRESS OF ATAN(|F|)
 331         MOVE.L          (a1)+,ATANF(a6)
 332         MOVE.L          (a1)+,ATANFHI(a6)
 333         MOVE.L          (a1)+,ATANFLO(a6)       ...ATANF IS NOW ATAN(|F|)
 334         MOVE.L          X(a6),d0                ...LOAD SIGN AND EXPO. AGAIN
 335         ANDI.L          #$80000000,d0   ...SIGN(F)
 336         OR.L            d0,ATANF(a6)    ...ATANF IS NOW SIGN(F)*ATAN(|F|)
 337         MOVE.L          (a7)+,d2        ...RESTORE d2
 338
 339 *--THAT'S ALL I HAVE TO DO FOR NOW,
 340 *--BUT ALAS, THE DIVIDE IS STILL CRANKING!
 341
 342 *--U IN FP0, WE ARE NOW READY TO COMPUTE ATAN(U) AS
 343 *--U + A1*U*V*(A2 + V*(A3 + V)), V = U*U
 344 *--THE POLYNOMIAL MAY LOOK STRANGE, BUT IS NEVERTHELESS CORRECT.
 345 *--THE NATURAL FORM IS U + U*V*(A1 + V*(A2 + V*A3))
 346 *--WHAT WE HAVE HERE IS MERELY  A1 = A3, A2 = A1/A3, A3 = A2/A3.
 347 *--THE REASON FOR THIS REARRANGEMENT IS TO MAKE THE INDEPENDENT
 348 *--PARTS A1*U*V AND (A2 + ... STUFF) MORE LOAD-BALANCED
 349
 350
 351         FMOVE.X         FP0,FP1
 352         FMUL.X          FP1,FP1
 353         FMOVE.D         ATANA3,FP2
 354         FADD.X          FP1,FP2         ...A3+V
 355         FMUL.X          FP1,FP2         ...V*(A3+V)
 356         FMUL.X          FP0,FP1         ...U*V
 357         FADD.D          ATANA2,FP2      ...A2+V*(A3+V)
 358         FMUL.D          ATANA1,FP1      ...A1*U*V
 359         FMUL.X          FP2,FP1         ...A1*U*V*(A2+V*(A3+V))
 360
 361         FADD.X          FP1,FP0         ...ATAN(U), FP1 RELEASED
 362         FMOVE.L         d1,FPCR         ;restore users exceptions
 363         FADD.X          ATANF(a6),FP0   ...ATAN(X)
 364         bra             t_frcinx
 365
 366 ATANBORS:
 367 *--|X| IS IN d0 IN COMPACT FORM. FP1, d0 SAVED.
 368 *--FP0 IS X AND |X| <= 1/16 OR |X| >= 16.
 369         CMPI.L          #$3FFF8000,d0
 370         BGT.W           ATANBIG ...I.E. |X| >= 16
 371
 372 ATANSM:
 373 *--|X| <= 1/16
 374 *--IF |X| < 2^(-40), RETURN X AS ANSWER. OTHERWISE, APPROXIMATE
 375 *--ATAN(X) BY X + X*Y*(B1+Y*(B2+Y*(B3+Y*(B4+Y*(B5+Y*B6)))))
 376 *--WHICH IS X + X*Y*( [B1+Z*(B3+Z*B5)] + [Y*(B2+Z*(B4+Z*B6)] )
 377 *--WHERE Y = X*X, AND Z = Y*Y.
 378
 379         CMPI.L          #$3FD78000,d0
 380         BLT.W           ATANTINY
 381 *--COMPUTE POLYNOMIAL
 382         FMUL.X          FP0,FP0 ...FP0 IS Y = X*X
 383
 384
 385         CLR.W           XDCARE(a6)
 386
 387         FMOVE.X         FP0,FP1
 388         FMUL.X          FP1,FP1         ...FP1 IS Z = Y*Y
 389
 390         FMOVE.D         ATANB6,FP2
 391         FMOVE.D         ATANB5,FP3
 392
 393         FMUL.X          FP1,FP2         ...Z*B6
 394         FMUL.X          FP1,FP3         ...Z*B5
 395
 396         FADD.D          ATANB4,FP2      ...B4+Z*B6
 397         FADD.D          ATANB3,FP3      ...B3+Z*B5
 398
 399         FMUL.X          FP1,FP2         ...Z*(B4+Z*B6)
 400         FMUL.X          FP3,FP1         ...Z*(B3+Z*B5)
 401
 402         FADD.D          ATANB2,FP2      ...B2+Z*(B4+Z*B6)
 403         FADD.D          ATANB1,FP1      ...B1+Z*(B3+Z*B5)
 404
 405         FMUL.X          FP0,FP2         ...Y*(B2+Z*(B4+Z*B6))
 406         FMUL.X          X(a6),FP0               ...X*Y
 407
 408         FADD.X          FP2,FP1         ...[B1+Z*(B3+Z*B5)]+[Y*(B2+Z*(B4+Z*B6))]
 409
 410
 411         FMUL.X          FP1,FP0 ...X*Y*([B1+Z*(B3+Z*B5)]+[Y*(B2+Z*(B4+Z*B6))])
 412
 413         FMOVE.L         d1,FPCR         ;restore users exceptions
 414         FADD.X          X(a6),FP0
 415
 416         bra             t_frcinx
 417
 418 ATANTINY:
 419 *--|X| < 2^(-40), ATAN(X) = X
 420         CLR.W           XDCARE(a6)
 421
 422         FMOVE.L         d1,FPCR         ;restore users exceptions
 423         FMOVE.X         X(a6),FP0       ;last inst - possible exception set
 424
 425         bra             t_frcinx
 426
 427 ATANBIG:
 428 *--IF |X| > 2^(100), RETURN     SIGN(X)*(PI/2 - TINY). OTHERWISE,
 429 *--RETURN SIGN(X)*PI/2 + ATAN(-1/X).
 430         CMPI.L          #$40638000,d0
 431         BGT.W           ATANHUGE
 432
 433 *--APPROXIMATE ATAN(-1/X) BY
 434 *--X'+X'*Y*(C1+Y*(C2+Y*(C3+Y*(C4+Y*C5)))), X' = -1/X, Y = X'*X'
 435 *--THIS CAN BE RE-WRITTEN AS
 436 *--X'+X'*Y*( [C1+Z*(C3+Z*C5)] + [Y*(C2+Z*C4)] ), Z = Y*Y.
 437
 438         FMOVE.S         #:BF800000,FP1  ...LOAD -1
 439         FDIV.X          FP0,FP1         ...FP1 IS -1/X
 440
 441
 442 *--DIVIDE IS STILL CRANKING
 443
 444         FMOVE.X         FP1,FP0         ...FP0 IS X'
 445         FMUL.X          FP0,FP0         ...FP0 IS Y = X'*X'
 446         FMOVE.X         FP1,X(a6)               ...X IS REALLY X'
 447
 448         FMOVE.X         FP0,FP1
 449         FMUL.X          FP1,FP1         ...FP1 IS Z = Y*Y
 450
 451         FMOVE.D         ATANC5,FP3
 452         FMOVE.D         ATANC4,FP2
 453
 454         FMUL.X          FP1,FP3         ...Z*C5
 455         FMUL.X          FP1,FP2         ...Z*B4
 456
 457         FADD.D          ATANC3,FP3      ...C3+Z*C5
 458         FADD.D          ATANC2,FP2      ...C2+Z*C4
 459
 460         FMUL.X          FP3,FP1         ...Z*(C3+Z*C5), FP3 RELEASED
 461         FMUL.X          FP0,FP2         ...Y*(C2+Z*C4)
 462
 463         FADD.D          ATANC1,FP1      ...C1+Z*(C3+Z*C5)
 464         FMUL.X          X(a6),FP0               ...X'*Y
 465
 466         FADD.X          FP2,FP1         ...[Y*(C2+Z*C4)]+[C1+Z*(C3+Z*C5)]
 467
 468
 469         FMUL.X          FP1,FP0         ...X'*Y*([B1+Z*(B3+Z*B5)]
 470 *                                       ...     +[Y*(B2+Z*(B4+Z*B6))])
 471         FADD.X          X(a6),FP0
 472
 473         FMOVE.L         d1,FPCR         ;restore users exceptions
 474
 475         btst.b          #7,(a0)
 476         beq.b           pos_big
 477
 478 neg_big:
 479         FADD.X          NPIBY2,FP0
 480         bra             t_frcinx
 481
 482 pos_big:
 483         FADD.X          PPIBY2,FP0
 484         bra             t_frcinx
 485
 486 ATANHUGE:
 487 *--RETURN SIGN(X)*(PIBY2 - TINY) = SIGN(X)*PIBY2 - SIGN(X)*TINY
 488         btst.b          #7,(a0)
 489         beq.b           pos_huge
 490
 491 neg_huge:
 492         FMOVE.X         NPIBY2,fp0
 493         fmove.l         d1,fpcr
 494         fsub.x          NTINY,fp0
 495         bra             t_frcinx
 496
 497 pos_huge:
 498         FMOVE.X         PPIBY2,fp0
 499         fmove.l         d1,fpcr
 500         fsub.x          PTINY,fp0
 501         bra             t_frcinx
 502
 503         end