File indexing completed on 2025-05-11 08:23:48
0001 #include "fpsp-namespace.h"
0002 //
0003 //
0004 // sasin.sa 3.3 12/19/90
0005 //
0006 // Description: The entry point sAsin computes the inverse sine of
0007 // an input argument; sAsind does the same except for denormalized
0008 // input.
0009 //
0010 // Input: Double-extended number X in location pointed to
0011 // by address register a0.
0012 //
0013 // Output: The value arcsin(X) returned in floating-point register Fp0.
0014 //
0015 // Accuracy and Monotonicity: The returned result is within 3 ulps in
0016 // 64 significant bit, i.e. within 0.5001 ulp to 53 bits if the
0017 // result is subsequently rounded to double precision. The
0018 // result is provably monotonic in double precision.
0019 //
0020 // Speed: The program sASIN takes approximately 310 cycles.
0021 //
0022 // Algorithm:
0023 //
0024 // ASIN
0025 // 1. If |X| >= 1, go to 3.
0026 //
0027 // 2. (|X| < 1) Calculate asin(X) by
0028 // z := sqrt( [1-X][1+X] )
0029 // asin(X) = atan( x / z ).
0030 // Exit.
0031 //
0032 // 3. If |X| > 1, go to 5.
0033 //
0034 // 4. (|X| = 1) sgn := sign(X), return asin(X) := sgn * Pi/2. Exit.
0035 //
0036 // 5. (|X| > 1) Generate an invalid operation by 0 * infinity.
0037 // Exit.
0038 //
0039
0040 // Copyright (C) Motorola, Inc. 1990
0041 // All Rights Reserved
0042 //
0043 // THIS IS UNPUBLISHED PROPRIETARY SOURCE CODE OF MOTOROLA
0044 // The copyright notice above does not evidence any
0045 // actual or intended publication of such source code.
0046
0047 //SASIN idnt 2,1 | Motorola 040 Floating Point Software Package
0048
0049 |section 8
0050
0051 PIBY2: .long 0x3FFF0000,0xC90FDAA2,0x2168C235,0x00000000
0052
0053 |xref t_operr
0054 |xref t_frcinx
0055 |xref t_extdnrm
0056 |xref satan
0057
0058 .global sasind
0059 sasind:
0060 //--ASIN(X) = X FOR DENORMALIZED X
0061
0062 bra t_extdnrm
0063
0064 .global sasin
0065 sasin:
0066 fmovex (%a0),%fp0 // ...LOAD INPUT
0067
0068 movel (%a0),%d0
0069 movew 4(%a0),%d0
0070 andil #0x7FFFFFFF,%d0
0071 cmpil #0x3FFF8000,%d0
0072 bges asinbig
0073
0074 //--THIS IS THE USUAL CASE, |X| < 1
0075 //--ASIN(X) = ATAN( X / SQRT( (1-X)(1+X) ) )
0076
0077 fmoves #0x3F800000,%fp1
0078 fsubx %fp0,%fp1 // ...1-X
0079 fmovemx %fp2-%fp2,-(%a7)
0080 fmoves #0x3F800000,%fp2
0081 faddx %fp0,%fp2 // ...1+X
0082 fmulx %fp2,%fp1 // ...(1+X)(1-X)
0083 fmovemx (%a7)+,%fp2-%fp2
0084 fsqrtx %fp1 // ...SQRT([1-X][1+X])
0085 fdivx %fp1,%fp0 // ...X/SQRT([1-X][1+X])
0086 fmovemx %fp0-%fp0,(%a0)
0087 bsr satan
0088 bra t_frcinx
0089
0090 asinbig:
0091 fabsx %fp0 // ...|X|
0092 fcmps #0x3F800000,%fp0
0093 fbgt t_operr //cause an operr exception
0094
0095 //--|X| = 1, ASIN(X) = +- PI/2.
0096
0097 fmovex PIBY2,%fp0
0098 movel (%a0),%d0
0099 andil #0x80000000,%d0 // ...SIGN BIT OF X
0100 oril #0x3F800000,%d0 // ...+-1 IN SGL FORMAT
0101 movel %d0,-(%sp) // ...push SIGN(X) IN SGL-FMT
0102 fmovel %d1,%FPCR
0103 fmuls (%sp)+,%fp0
0104 bra t_frcinx
0105
0106 |end
