Commit | Line | Data |
---|---|---|
1da177e4 LT |
1 | | |
2 | | sacos.sa 3.3 12/19/90 | |
3 | | | |
4 | | Description: The entry point sAcos computes the inverse cosine of | |
5 | | an input argument; sAcosd 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 arccos(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 sCOS takes approximately 310 cycles. | |
19 | | | |
20 | | Algorithm: | |
21 | | | |
22 | | ACOS | |
23 | | 1. If |X| >= 1, go to 3. | |
24 | | | |
25 | | 2. (|X| < 1) Calculate acos(X) by | |
26 | | z := (1-X) / (1+X) | |
27 | | acos(X) = 2 * atan( sqrt(z) ). | |
28 | | Exit. | |
29 | | | |
30 | | 3. If |X| > 1, go to 5. | |
31 | | | |
32 | | 4. (|X| = 1) If X > 0, return 0. Otherwise, return Pi. 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 | | THIS IS UNPUBLISHED PROPRIETARY SOURCE CODE OF MOTOROLA | |
42 | | The copyright notice above does not evidence any | |
43 | | actual or intended publication of such source code. | |
44 | ||
45 | |SACOS idnt 2,1 | Motorola 040 Floating Point Software Package | |
46 | ||
47 | |section 8 | |
48 | ||
49 | PI: .long 0x40000000,0xC90FDAA2,0x2168C235,0x00000000 | |
50 | PIBY2: .long 0x3FFF0000,0xC90FDAA2,0x2168C235,0x00000000 | |
51 | ||
52 | |xref t_operr | |
53 | |xref t_frcinx | |
54 | |xref satan | |
55 | ||
56 | .global sacosd | |
57 | sacosd: | |
58 | |--ACOS(X) = PI/2 FOR DENORMALIZED X | |
59 | fmovel %d1,%fpcr | ...load user's rounding mode/precision | |
60 | fmovex PIBY2,%fp0 | |
61 | bra t_frcinx | |
62 | ||
63 | .global sacos | |
64 | sacos: | |
65 | fmovex (%a0),%fp0 | ...LOAD INPUT | |
66 | ||
67 | movel (%a0),%d0 | ...pack exponent with upper 16 fraction | |
68 | movew 4(%a0),%d0 | |
69 | andil #0x7FFFFFFF,%d0 | |
70 | cmpil #0x3FFF8000,%d0 | |
71 | bges ACOSBIG | |
72 | ||
73 | |--THIS IS THE USUAL CASE, |X| < 1 | |
74 | |--ACOS(X) = 2 * ATAN( SQRT( (1-X)/(1+X) ) ) | |
75 | ||
76 | fmoves #0x3F800000,%fp1 | |
77 | faddx %fp0,%fp1 | ...1+X | |
78 | fnegx %fp0 | ... -X | |
79 | fadds #0x3F800000,%fp0 | ...1-X | |
80 | fdivx %fp1,%fp0 | ...(1-X)/(1+X) | |
81 | fsqrtx %fp0 | ...SQRT((1-X)/(1+X)) | |
82 | fmovemx %fp0-%fp0,(%a0) | ...overwrite input | |
83 | movel %d1,-(%sp) |save original users fpcr | |
84 | clrl %d1 | |
85 | bsr satan | ...ATAN(SQRT([1-X]/[1+X])) | |
86 | fmovel (%sp)+,%fpcr |restore users exceptions | |
87 | faddx %fp0,%fp0 | ...2 * ATAN( STUFF ) | |
88 | bra t_frcinx | |
89 | ||
90 | ACOSBIG: | |
91 | fabsx %fp0 | |
92 | fcmps #0x3F800000,%fp0 | |
93 | fbgt t_operr |cause an operr exception | |
94 | ||
95 | |--|X| = 1, ACOS(X) = 0 OR PI | |
96 | movel (%a0),%d0 | ...pack exponent with upper 16 fraction | |
97 | movew 4(%a0),%d0 | |
98 | cmpl #0,%d0 |D0 has original exponent+fraction | |
99 | bgts ACOSP1 | |
100 | ||
101 | |--X = -1 | |
102 | |Returns PI and inexact exception | |
103 | fmovex PI,%fp0 | |
104 | fmovel %d1,%FPCR | |
105 | fadds #0x00800000,%fp0 |cause an inexact exception to be put | |
106 | | ;into the 040 - will not trap until next | |
107 | | ;fp inst. | |
108 | bra t_frcinx | |
109 | ||
110 | ACOSP1: | |
111 | fmovel %d1,%FPCR | |
112 | fmoves #0x00000000,%fp0 | |
113 | rts |Facos ; of +1 is exact | |
114 | ||
115 | |end |