[BACK]Return to satan.sa CVS log [TXT][DIR] Up to [local] / sys / arch / m68k / fpsp

Annotation of sys/arch/m68k/fpsp/satan.sa, Revision 1.1.1.1

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

CVSweb