sasin.S (2348B)
1| 2| sasin.sa 3.3 12/19/90 3| 4| Description: The entry point sAsin computes the inverse sine of 5| an input argument; sAsind does the same except for denormalized 6| input. 7| 8| Input: Double-extended number X in location pointed to 9| by address register a0. 10| 11| Output: The value arcsin(X) returned in floating-point register Fp0. 12| 13| Accuracy and Monotonicity: The returned result is within 3 ulps in 14| 64 significant bit, i.e. within 0.5001 ulp to 53 bits if the 15| result is subsequently rounded to double precision. The 16| result is provably monotonic in double precision. 17| 18| Speed: The program sASIN takes approximately 310 cycles. 19| 20| Algorithm: 21| 22| ASIN 23| 1. If |X| >= 1, go to 3. 24| 25| 2. (|X| < 1) Calculate asin(X) by 26| z := sqrt( [1-X][1+X] ) 27| asin(X) = atan( x / z ). 28| Exit. 29| 30| 3. If |X| > 1, go to 5. 31| 32| 4. (|X| = 1) sgn := sign(X), return asin(X) := sgn * Pi/2. Exit. 33| 34| 5. (|X| > 1) Generate an invalid operation by 0 * infinity. 35| Exit. 36| 37 38| Copyright (C) Motorola, Inc. 1990 39| All Rights Reserved 40| 41| For details on the license for this file, please see the 42| file, README, in this same directory. 43 44|SASIN idnt 2,1 | Motorola 040 Floating Point Software Package 45 46 |section 8 47 48PIBY2: .long 0x3FFF0000,0xC90FDAA2,0x2168C235,0x00000000 49 50 |xref t_operr 51 |xref t_frcinx 52 |xref t_extdnrm 53 |xref satan 54 55 .global sasind 56sasind: 57|--ASIN(X) = X FOR DENORMALIZED X 58 59 bra t_extdnrm 60 61 .global sasin 62sasin: 63 fmovex (%a0),%fp0 | ...LOAD INPUT 64 65 movel (%a0),%d0 66 movew 4(%a0),%d0 67 andil #0x7FFFFFFF,%d0 68 cmpil #0x3FFF8000,%d0 69 bges asinbig 70 71|--THIS IS THE USUAL CASE, |X| < 1 72|--ASIN(X) = ATAN( X / SQRT( (1-X)(1+X) ) ) 73 74 fmoves #0x3F800000,%fp1 75 fsubx %fp0,%fp1 | ...1-X 76 fmovemx %fp2-%fp2,-(%a7) 77 fmoves #0x3F800000,%fp2 78 faddx %fp0,%fp2 | ...1+X 79 fmulx %fp2,%fp1 | ...(1+X)(1-X) 80 fmovemx (%a7)+,%fp2-%fp2 81 fsqrtx %fp1 | ...SQRT([1-X][1+X]) 82 fdivx %fp1,%fp0 | ...X/SQRT([1-X][1+X]) 83 fmovemx %fp0-%fp0,(%a0) 84 bsr satan 85 bra t_frcinx 86 87asinbig: 88 fabsx %fp0 | ...|X| 89 fcmps #0x3F800000,%fp0 90 fbgt t_operr |cause an operr exception 91 92|--|X| = 1, ASIN(X) = +- PI/2. 93 94 fmovex PIBY2,%fp0 95 movel (%a0),%d0 96 andil #0x80000000,%d0 | ...SIGN BIT OF X 97 oril #0x3F800000,%d0 | ...+-1 IN SGL FORMAT 98 movel %d0,-(%sp) | ...push SIGN(X) IN SGL-FMT 99 fmovel %d1,%FPCR 100 fmuls (%sp)+,%fp0 101 bra t_frcinx 102 103 |end