1 | 1 |
2 | ssin.sa 3.3 7/29/91 2 | ssin.sa 3.3 7/29/91
3 | 3 |
4 | The entry point sSIN computes the sine 4 | The entry point sSIN computes the sine of an input argument
5 | sCOS computes the cosine, and sSINCOS 5 | sCOS computes the cosine, and sSINCOS computes both. The
6 | corresponding entry points with a "d" 6 | corresponding entry points with a "d" computes the same
7 | corresponding function values for deno 7 | corresponding function values for denormalized inputs.
8 | 8 |
9 | Input: Double-extended number X in loc 9 | Input: Double-extended number X in location pointed to
10 | by address register a0. 10 | by address register a0.
11 | 11 |
12 | Output: The function value sin(X) or c 12 | Output: The function value sin(X) or cos(X) returned in Fp0 if SIN or
13 | COS is requested. Otherwise, f 13 | COS is requested. Otherwise, for SINCOS, sin(X) is returned
14 | in Fp0, and cos(X) is returned 14 | in Fp0, and cos(X) is returned in Fp1.
15 | 15 |
16 | Modifies: Fp0 for SIN or COS; both Fp0 16 | Modifies: Fp0 for SIN or COS; both Fp0 and Fp1 for SINCOS.
17 | 17 |
18 | Accuracy and Monotonicity: The returne 18 | Accuracy and Monotonicity: The returned result is within 1 ulp in
19 | 64 significant bit, i.e. withi 19 | 64 significant bit, i.e. within 0.5001 ulp to 53 bits if the
20 | result is subsequently rounded 20 | result is subsequently rounded to double precision. The
21 | result is provably monotonic i 21 | result is provably monotonic in double precision.
22 | 22 |
23 | Speed: The programs sSIN and sCOS take 23 | Speed: The programs sSIN and sCOS take approximately 150 cycles for
24 | input argument X such that |X| 24 | input argument X such that |X| < 15Pi, which is the usual
25 | situation. The speed for sSINC 25 | situation. The speed for sSINCOS is approximately 190 cycles.
26 | 26 |
27 | Algorithm: 27 | Algorithm:
28 | 28 |
29 | SIN and COS: 29 | SIN and COS:
30 | 1. If SIN is invoked, set AdjN := 0; o 30 | 1. If SIN is invoked, set AdjN := 0; otherwise, set AdjN := 1.
31 | 31 |
32 | 2. If |X| >= 15Pi or |X| < 2**(-40), g 32 | 2. If |X| >= 15Pi or |X| < 2**(-40), go to 7.
33 | 33 |
34 | 3. Decompose X as X = N(Pi/2) + r wher 34 | 3. Decompose X as X = N(Pi/2) + r where |r| <= Pi/4. Let
35 | k = N mod 4, so in particular, 35 | k = N mod 4, so in particular, k = 0,1,2,or 3. Overwrite
36 | k by k := k + AdjN. 36 | k by k := k + AdjN.
37 | 37 |
38 | 4. If k is even, go to 6. 38 | 4. If k is even, go to 6.
39 | 39 |
40 | 5. (k is odd) Set j := (k-1)/2, sgn := 40 | 5. (k is odd) Set j := (k-1)/2, sgn := (-1)**j. Return sgn*cos(r)
41 | where cos(r) is approximated b 41 | where cos(r) is approximated by an even polynomial in r,
42 | 1 + r*r*(B1+s*(B2+ ... + s*B8) 42 | 1 + r*r*(B1+s*(B2+ ... + s*B8)), s = r*r.
43 | Exit. 43 | Exit.
44 | 44 |
45 | 6. (k is even) Set j := k/2, sgn := (- 45 | 6. (k is even) Set j := k/2, sgn := (-1)**j. Return sgn*sin(r)
46 | where sin(r) is approximated b 46 | where sin(r) is approximated by an odd polynomial in r
47 | r + r*s*(A1+s*(A2+ ... + s*A7) 47 | r + r*s*(A1+s*(A2+ ... + s*A7)), s = r*r.
48 | Exit. 48 | Exit.
49 | 49 |
50 | 7. If |X| > 1, go to 9. 50 | 7. If |X| > 1, go to 9.
51 | 51 |
52 | 8. (|X|<2**(-40)) If SIN is invoked, r 52 | 8. (|X|<2**(-40)) If SIN is invoked, return X; otherwise return 1.
53 | 53 |
54 | 9. Overwrite X by X := X rem 2Pi. Now 54 | 9. Overwrite X by X := X rem 2Pi. Now that |X| <= Pi, go back to 3.
55 | 55 |
56 | SINCOS: 56 | SINCOS:
57 | 1. If |X| >= 15Pi or |X| < 2**(-40), g 57 | 1. If |X| >= 15Pi or |X| < 2**(-40), go to 6.
58 | 58 |
59 | 2. Decompose X as X = N(Pi/2) + r wher 59 | 2. Decompose X as X = N(Pi/2) + r where |r| <= Pi/4. Let
60 | k = N mod 4, so in particular, 60 | k = N mod 4, so in particular, k = 0,1,2,or 3.
61 | 61 |
62 | 3. If k is even, go to 5. 62 | 3. If k is even, go to 5.
63 | 63 |
64 | 4. (k is odd) Set j1 := (k-1)/2, j2 := 64 | 4. (k is odd) Set j1 := (k-1)/2, j2 := j1 (EOR) (k mod 2), i.e.
65 | j1 exclusive or with the l.s.b 65 | j1 exclusive or with the l.s.b. of k.
66 | sgn1 := (-1)**j1, sgn2 := (-1) 66 | sgn1 := (-1)**j1, sgn2 := (-1)**j2.
67 | SIN(X) = sgn1 * cos(r) and COS 67 | SIN(X) = sgn1 * cos(r) and COS(X) = sgn2*sin(r) where
68 | sin(r) and cos(r) are computed 68 | sin(r) and cos(r) are computed as odd and even polynomials
69 | in r, respectively. Exit 69 | in r, respectively. Exit
70 | 70 |
71 | 5. (k is even) Set j1 := k/2, sgn1 := 71 | 5. (k is even) Set j1 := k/2, sgn1 := (-1)**j1.
72 | SIN(X) = sgn1 * sin(r) and COS 72 | SIN(X) = sgn1 * sin(r) and COS(X) = sgn1*cos(r) where
73 | sin(r) and cos(r) are computed 73 | sin(r) and cos(r) are computed as odd and even polynomials
74 | in r, respectively. Exit 74 | in r, respectively. Exit
75 | 75 |
76 | 6. If |X| > 1, go to 8. 76 | 6. If |X| > 1, go to 8.
77 | 77 |
78 | 7. (|X|<2**(-40)) SIN(X) = X and COS(X 78 | 7. (|X|<2**(-40)) SIN(X) = X and COS(X) = 1. Exit.
79 | 79 |
80 | 8. Overwrite X by X := X rem 2Pi. Now 80 | 8. Overwrite X by X := X rem 2Pi. Now that |X| <= Pi, go back to 2.
81 | 81 |
82 82
83 | Copyright (C) Motorola, Inc. 1 83 | Copyright (C) Motorola, Inc. 1990
84 | All Rights Reserved 84 | All Rights Reserved
85 | 85 |
86 | For details on the license for this fi 86 | For details on the license for this file, please see the
87 | file, README, in this same directory. 87 | file, README, in this same directory.
88 88
89 |SSIN idnt 2,1 | Motorola 040 Floating Po 89 |SSIN idnt 2,1 | Motorola 040 Floating Point Software Package
90 90
91 |section 8 91 |section 8
92 92
93 #include "fpsp.h" 93 #include "fpsp.h"
94 94
95 BOUNDS1: .long 0x3FD78000,0x4004BC7E 95 BOUNDS1: .long 0x3FD78000,0x4004BC7E
96 TWOBYPI: .long 0x3FE45F30,0x6DC9C883 96 TWOBYPI: .long 0x3FE45F30,0x6DC9C883
97 97
98 SINA7: .long 0xBD6AAA77,0xCCC994F5 98 SINA7: .long 0xBD6AAA77,0xCCC994F5
99 SINA6: .long 0x3DE61209,0x7AAE8DA1 99 SINA6: .long 0x3DE61209,0x7AAE8DA1
100 100
101 SINA5: .long 0xBE5AE645,0x2A118AE4 101 SINA5: .long 0xBE5AE645,0x2A118AE4
102 SINA4: .long 0x3EC71DE3,0xA5341531 102 SINA4: .long 0x3EC71DE3,0xA5341531
103 103
104 SINA3: .long 0xBF2A01A0,0x1A018B59,0x00000000 104 SINA3: .long 0xBF2A01A0,0x1A018B59,0x00000000,0x00000000
105 105
106 SINA2: .long 0x3FF80000,0x88888888,0x888859AF 106 SINA2: .long 0x3FF80000,0x88888888,0x888859AF,0x00000000
107 107
108 SINA1: .long 0xBFFC0000,0xAAAAAAAA,0xAAAAAA99 108 SINA1: .long 0xBFFC0000,0xAAAAAAAA,0xAAAAAA99,0x00000000
109 109
110 COSB8: .long 0x3D2AC4D0,0xD6011EE3 110 COSB8: .long 0x3D2AC4D0,0xD6011EE3
111 COSB7: .long 0xBDA9396F,0x9F45AC19 111 COSB7: .long 0xBDA9396F,0x9F45AC19
112 112
113 COSB6: .long 0x3E21EED9,0x0612C972 113 COSB6: .long 0x3E21EED9,0x0612C972
114 COSB5: .long 0xBE927E4F,0xB79D9FCF 114 COSB5: .long 0xBE927E4F,0xB79D9FCF
115 115
116 COSB4: .long 0x3EFA01A0,0x1A01D423,0x00000000 116 COSB4: .long 0x3EFA01A0,0x1A01D423,0x00000000,0x00000000
117 117
118 COSB3: .long 0xBFF50000,0xB60B60B6,0x0B61D438 118 COSB3: .long 0xBFF50000,0xB60B60B6,0x0B61D438,0x00000000
119 119
120 COSB2: .long 0x3FFA0000,0xAAAAAAAA,0xAAAAAB5E 120 COSB2: .long 0x3FFA0000,0xAAAAAAAA,0xAAAAAB5E
121 COSB1: .long 0xBF000000 121 COSB1: .long 0xBF000000
122 122
123 INVTWOPI: .long 0x3FFC0000,0xA2F9836E,0x4E4415 123 INVTWOPI: .long 0x3FFC0000,0xA2F9836E,0x4E44152A
124 124
125 TWOPI1: .long 0x40010000,0xC90FDAA2,0x00000000 125 TWOPI1: .long 0x40010000,0xC90FDAA2,0x00000000,0x00000000
126 TWOPI2: .long 0x3FDF0000,0x85A308D4,0x00000000 126 TWOPI2: .long 0x3FDF0000,0x85A308D4,0x00000000,0x00000000
127 127
128 |xref PITBL 128 |xref PITBL
129 129
130 .set INARG,FP_SCR4 130 .set INARG,FP_SCR4
131 131
132 .set X,FP_SCR5 132 .set X,FP_SCR5
133 .set XDCARE,X+2 133 .set XDCARE,X+2
134 .set XFRAC,X+4 134 .set XFRAC,X+4
135 135
136 .set RPRIME,FP_SCR1 136 .set RPRIME,FP_SCR1
137 .set SPRIME,FP_SCR2 137 .set SPRIME,FP_SCR2
138 138
139 .set POSNEG1,L_SCR1 139 .set POSNEG1,L_SCR1
140 .set TWOTO63,L_SCR1 140 .set TWOTO63,L_SCR1
141 141
142 .set ENDFLAG,L_SCR2 142 .set ENDFLAG,L_SCR2
143 .set N,L_SCR2 143 .set N,L_SCR2
144 144
145 .set ADJN,L_SCR3 145 .set ADJN,L_SCR3
146 146
147 | xref t_frcinx 147 | xref t_frcinx
148 |xref t_extdnrm 148 |xref t_extdnrm
149 |xref sto_cos 149 |xref sto_cos
150 150
151 .global ssind 151 .global ssind
152 ssind: 152 ssind:
153 |--SIN(X) = X FOR DENORMALIZED X 153 |--SIN(X) = X FOR DENORMALIZED X
154 bra t_extdnrm 154 bra t_extdnrm
155 155
156 .global scosd 156 .global scosd
157 scosd: 157 scosd:
158 |--COS(X) = 1 FOR DENORMALIZED X 158 |--COS(X) = 1 FOR DENORMALIZED X
159 159
160 fmoves #0x3F800000,%fp0 160 fmoves #0x3F800000,%fp0
161 | 161 |
162 | 9D25B Fix: Sometimes the previous fmov 162 | 9D25B Fix: Sometimes the previous fmove.s sets fpsr bits
163 | 163 |
164 fmovel #0,%fpsr 164 fmovel #0,%fpsr
165 | 165 |
166 bra t_frcinx 166 bra t_frcinx
167 167
168 .global ssin 168 .global ssin
169 ssin: 169 ssin:
170 |--SET ADJN TO 0 170 |--SET ADJN TO 0
171 movel #0,ADJN(%a6) 171 movel #0,ADJN(%a6)
172 bras SINBGN 172 bras SINBGN
173 173
174 .global scos 174 .global scos
175 scos: 175 scos:
176 |--SET ADJN TO 1 176 |--SET ADJN TO 1
177 movel #1,ADJN(%a6) 177 movel #1,ADJN(%a6)
178 178
179 SINBGN: 179 SINBGN:
180 |--SAVE FPCR, FP1. CHECK IF |X| IS TOO SMALL O 180 |--SAVE FPCR, FP1. CHECK IF |X| IS TOO SMALL OR LARGE
181 181
182 fmovex (%a0),%fp0 | ...L 182 fmovex (%a0),%fp0 | ...LOAD INPUT
183 183
184 movel (%a0),%d0 184 movel (%a0),%d0
185 movew 4(%a0),%d0 185 movew 4(%a0),%d0
186 fmovex %fp0,X(%a6) 186 fmovex %fp0,X(%a6)
187 andil #0x7FFFFFFF,%d0 187 andil #0x7FFFFFFF,%d0 | ...COMPACTIFY X
188 188
189 cmpil #0x3FD78000,%d0 189 cmpil #0x3FD78000,%d0 | ...|X| >= 2**(-40)?
190 bges SOK1 190 bges SOK1
191 bra SINSM 191 bra SINSM
192 192
193 SOK1: 193 SOK1:
194 cmpil #0x4004BC7E,%d0 194 cmpil #0x4004BC7E,%d0 | ...|X| < 15 PI?
195 blts SINMAIN 195 blts SINMAIN
196 bra REDUCEX 196 bra REDUCEX
197 197
198 SINMAIN: 198 SINMAIN:
199 |--THIS IS THE USUAL CASE, |X| <= 15 PI. 199 |--THIS IS THE USUAL CASE, |X| <= 15 PI.
200 |--THE ARGUMENT REDUCTION IS DONE BY TABLE LOO 200 |--THE ARGUMENT REDUCTION IS DONE BY TABLE LOOK UP.
201 fmovex %fp0,%fp1 201 fmovex %fp0,%fp1
202 fmuld TWOBYPI,%fp1 | ...X 202 fmuld TWOBYPI,%fp1 | ...X*2/PI
203 203
204 |--HIDE THE NEXT THREE INSTRUCTIONS 204 |--HIDE THE NEXT THREE INSTRUCTIONS
205 lea PITBL+0x200,%a1 | ...T 205 lea PITBL+0x200,%a1 | ...TABLE OF N*PI/2, N = -32,...,32
206 206
207 207
208 |--FP1 IS NOW READY 208 |--FP1 IS NOW READY
209 fmovel %fp1,N(%a6) 209 fmovel %fp1,N(%a6) | ...CONVERT TO INTEGER
210 210
211 movel N(%a6),%d0 211 movel N(%a6),%d0
212 asll #4,%d0 212 asll #4,%d0
213 addal %d0,%a1 | ...A1 IS THE 213 addal %d0,%a1 | ...A1 IS THE ADDRESS OF N*PIBY2
214 | ...WHICH IS IN 214 | ...WHICH IS IN TWO PIECES Y1 & Y2
215 215
216 fsubx (%a1)+,%fp0 | ...X 216 fsubx (%a1)+,%fp0 | ...X-Y1
217 |--HIDE THE NEXT ONE 217 |--HIDE THE NEXT ONE
218 fsubs (%a1),%fp0 | ...F 218 fsubs (%a1),%fp0 | ...FP0 IS R = (X-Y1)-Y2
219 219
220 SINCONT: 220 SINCONT:
221 |--continuation from REDUCEX 221 |--continuation from REDUCEX
222 222
223 |--GET N+ADJN AND SEE IF SIN(R) OR COS(R) IS N 223 |--GET N+ADJN AND SEE IF SIN(R) OR COS(R) IS NEEDED
224 movel N(%a6),%d0 224 movel N(%a6),%d0
225 addl ADJN(%a6),%d0 | ...S 225 addl ADJN(%a6),%d0 | ...SEE IF D0 IS ODD OR EVEN
226 rorl #1,%d0 | ...D0 WAS OD 226 rorl #1,%d0 | ...D0 WAS ODD IFF D0 IS NEGATIVE
227 cmpil #0,%d0 227 cmpil #0,%d0
228 blt COSPOLY 228 blt COSPOLY
229 229
230 SINPOLY: 230 SINPOLY:
231 |--LET J BE THE LEAST SIG. BIT OF D0, LET SGN 231 |--LET J BE THE LEAST SIG. BIT OF D0, LET SGN := (-1)**J.
232 |--THEN WE RETURN SGN*SIN(R). SGN*SIN(R) 232 |--THEN WE RETURN SGN*SIN(R). SGN*SIN(R) IS COMPUTED BY
233 |--R' + R'*S*(A1 + S(A2 + S(A3 + S(A4 + ... + 233 |--R' + R'*S*(A1 + S(A2 + S(A3 + S(A4 + ... + SA7)))), WHERE
234 |--R' = SGN*R, S=R*R. THIS CAN BE REWRITTEN AS 234 |--R' = SGN*R, S=R*R. THIS CAN BE REWRITTEN AS
235 |--R' + R'*S*( [A1+T(A3+T(A5+TA7))] + [S(A2+T( 235 |--R' + R'*S*( [A1+T(A3+T(A5+TA7))] + [S(A2+T(A4+TA6))])
236 |--WHERE T=S*S. 236 |--WHERE T=S*S.
237 |--NOTE THAT A3 THROUGH A7 ARE STORED IN DOUBL 237 |--NOTE THAT A3 THROUGH A7 ARE STORED IN DOUBLE PRECISION
238 |--WHILE A1 AND A2 ARE IN DOUBLE-EXTENDED FORM 238 |--WHILE A1 AND A2 ARE IN DOUBLE-EXTENDED FORMAT.
239 fmovex %fp0,X(%a6) | ...X 239 fmovex %fp0,X(%a6) | ...X IS R
240 fmulx %fp0,%fp0 | ...F 240 fmulx %fp0,%fp0 | ...FP0 IS S
241 |---HIDE THE NEXT TWO WHILE WAITING FOR FP0 241 |---HIDE THE NEXT TWO WHILE WAITING FOR FP0
242 fmoved SINA7,%fp3 242 fmoved SINA7,%fp3
243 fmoved SINA6,%fp2 243 fmoved SINA6,%fp2
244 |--FP0 IS NOW READY 244 |--FP0 IS NOW READY
245 fmovex %fp0,%fp1 245 fmovex %fp0,%fp1
246 fmulx %fp1,%fp1 | ...F 246 fmulx %fp1,%fp1 | ...FP1 IS T
247 |--HIDE THE NEXT TWO WHILE WAITING FOR FP1 247 |--HIDE THE NEXT TWO WHILE WAITING FOR FP1
248 248
249 rorl #1,%d0 249 rorl #1,%d0
250 andil #0x80000000,%d0 250 andil #0x80000000,%d0
251 | ...LEAST SIG. 251 | ...LEAST SIG. BIT OF D0 IN SIGN POSITION
252 eorl %d0,X(%a6) | ...X 252 eorl %d0,X(%a6) | ...X IS NOW R'= SGN*R
253 253
254 fmulx %fp1,%fp3 | ...T 254 fmulx %fp1,%fp3 | ...TA7
255 fmulx %fp1,%fp2 | ...T 255 fmulx %fp1,%fp2 | ...TA6
256 256
257 faddd SINA5,%fp3 | ...A5+TA7 257 faddd SINA5,%fp3 | ...A5+TA7
258 faddd SINA4,%fp2 | ...A4+TA6 258 faddd SINA4,%fp2 | ...A4+TA6
259 259
260 fmulx %fp1,%fp3 | ...T 260 fmulx %fp1,%fp3 | ...T(A5+TA7)
261 fmulx %fp1,%fp2 | ...T 261 fmulx %fp1,%fp2 | ...T(A4+TA6)
262 262
263 faddd SINA3,%fp3 | ...A3+T(A 263 faddd SINA3,%fp3 | ...A3+T(A5+TA7)
264 faddx SINA2,%fp2 | ...A2+T(A 264 faddx SINA2,%fp2 | ...A2+T(A4+TA6)
265 265
266 fmulx %fp3,%fp1 | ...T 266 fmulx %fp3,%fp1 | ...T(A3+T(A5+TA7))
267 267
268 fmulx %fp0,%fp2 | ...S 268 fmulx %fp0,%fp2 | ...S(A2+T(A4+TA6))
269 faddx SINA1,%fp1 | ...A1+T(A 269 faddx SINA1,%fp1 | ...A1+T(A3+T(A5+TA7))
270 fmulx X(%a6),%fp0 | ...R 270 fmulx X(%a6),%fp0 | ...R'*S
271 271
272 faddx %fp2,%fp1 | ...[ 272 faddx %fp2,%fp1 | ...[A1+T(A3+T(A5+TA7))]+[S(A2+T(A4+TA6))]
273 |--FP3 RELEASED, RESTORE NOW AND TAKE SOME ADV 273 |--FP3 RELEASED, RESTORE NOW AND TAKE SOME ADVANTAGE OF HIDING
274 |--FP2 RELEASED, RESTORE NOW AND TAKE FULL ADV 274 |--FP2 RELEASED, RESTORE NOW AND TAKE FULL ADVANTAGE OF HIDING
275 275
276 276
277 fmulx %fp1,%fp0 277 fmulx %fp1,%fp0 | ...SIN(R')-R'
278 |--FP1 RELEASED. 278 |--FP1 RELEASED.
279 279
280 fmovel %d1,%FPCR 280 fmovel %d1,%FPCR |restore users exceptions
281 faddx X(%a6),%fp0 281 faddx X(%a6),%fp0 |last inst - possible exception set
282 bra t_frcinx 282 bra t_frcinx
283 283
284 284
285 COSPOLY: 285 COSPOLY:
286 |--LET J BE THE LEAST SIG. BIT OF D0, LET SGN 286 |--LET J BE THE LEAST SIG. BIT OF D0, LET SGN := (-1)**J.
287 |--THEN WE RETURN SGN*COS(R). SGN*COS(R) 287 |--THEN WE RETURN SGN*COS(R). SGN*COS(R) IS COMPUTED BY
288 |--SGN + S'*(B1 + S(B2 + S(B3 + S(B4 + ... + S 288 |--SGN + S'*(B1 + S(B2 + S(B3 + S(B4 + ... + SB8)))), WHERE
289 |--S=R*R AND S'=SGN*S. THIS CAN BE REWRITTEN A 289 |--S=R*R AND S'=SGN*S. THIS CAN BE REWRITTEN AS
290 |--SGN + S'*([B1+T(B3+T(B5+TB7))] + [S(B2+T(B4 290 |--SGN + S'*([B1+T(B3+T(B5+TB7))] + [S(B2+T(B4+T(B6+TB8)))])
291 |--WHERE T=S*S. 291 |--WHERE T=S*S.
292 |--NOTE THAT B4 THROUGH B8 ARE STORED IN DOUBL 292 |--NOTE THAT B4 THROUGH B8 ARE STORED IN DOUBLE PRECISION
293 |--WHILE B2 AND B3 ARE IN DOUBLE-EXTENDED FORM 293 |--WHILE B2 AND B3 ARE IN DOUBLE-EXTENDED FORMAT, B1 IS -1/2
294 |--AND IS THEREFORE STORED AS SINGLE PRECISION 294 |--AND IS THEREFORE STORED AS SINGLE PRECISION.
295 295
296 fmulx %fp0,%fp0 | ...F 296 fmulx %fp0,%fp0 | ...FP0 IS S
297 |---HIDE THE NEXT TWO WHILE WAITING FOR FP0 297 |---HIDE THE NEXT TWO WHILE WAITING FOR FP0
298 fmoved COSB8,%fp2 298 fmoved COSB8,%fp2
299 fmoved COSB7,%fp3 299 fmoved COSB7,%fp3
300 |--FP0 IS NOW READY 300 |--FP0 IS NOW READY
301 fmovex %fp0,%fp1 301 fmovex %fp0,%fp1
302 fmulx %fp1,%fp1 | ...F 302 fmulx %fp1,%fp1 | ...FP1 IS T
303 |--HIDE THE NEXT TWO WHILE WAITING FOR FP1 303 |--HIDE THE NEXT TWO WHILE WAITING FOR FP1
304 fmovex %fp0,X(%a6) | ...X 304 fmovex %fp0,X(%a6) | ...X IS S
305 rorl #1,%d0 305 rorl #1,%d0
306 andil #0x80000000,%d0 306 andil #0x80000000,%d0
307 | ...LEAST SIG. BIT OF D 307 | ...LEAST SIG. BIT OF D0 IN SIGN POSITION
308 308
309 fmulx %fp1,%fp2 | ...T 309 fmulx %fp1,%fp2 | ...TB8
310 |--HIDE THE NEXT TWO WHILE WAITING FOR THE XU 310 |--HIDE THE NEXT TWO WHILE WAITING FOR THE XU
311 eorl %d0,X(%a6) | ...X 311 eorl %d0,X(%a6) | ...X IS NOW S'= SGN*S
312 andil #0x80000000,%d0 312 andil #0x80000000,%d0
313 313
314 fmulx %fp1,%fp3 | ...T 314 fmulx %fp1,%fp3 | ...TB7
315 |--HIDE THE NEXT TWO WHILE WAITING FOR THE XU 315 |--HIDE THE NEXT TWO WHILE WAITING FOR THE XU
316 oril #0x3F800000,%d0 | ...D 316 oril #0x3F800000,%d0 | ...D0 IS SGN IN SINGLE
317 movel %d0,POSNEG1(%a6) 317 movel %d0,POSNEG1(%a6)
318 318
319 faddd COSB6,%fp2 | ...B6+TB8 319 faddd COSB6,%fp2 | ...B6+TB8
320 faddd COSB5,%fp3 | ...B5+TB7 320 faddd COSB5,%fp3 | ...B5+TB7
321 321
322 fmulx %fp1,%fp2 | ...T 322 fmulx %fp1,%fp2 | ...T(B6+TB8)
323 fmulx %fp1,%fp3 | ...T 323 fmulx %fp1,%fp3 | ...T(B5+TB7)
324 324
325 faddd COSB4,%fp2 | ...B4+T(B 325 faddd COSB4,%fp2 | ...B4+T(B6+TB8)
326 faddx COSB3,%fp3 | ...B3+T(B 326 faddx COSB3,%fp3 | ...B3+T(B5+TB7)
327 327
328 fmulx %fp1,%fp2 | ...T 328 fmulx %fp1,%fp2 | ...T(B4+T(B6+TB8))
329 fmulx %fp3,%fp1 | ...T 329 fmulx %fp3,%fp1 | ...T(B3+T(B5+TB7))
330 330
331 faddx COSB2,%fp2 | ...B2+T(B 331 faddx COSB2,%fp2 | ...B2+T(B4+T(B6+TB8))
332 fadds COSB1,%fp1 | ...B1+T(B 332 fadds COSB1,%fp1 | ...B1+T(B3+T(B5+TB7))
333 333
334 fmulx %fp2,%fp0 | ...S 334 fmulx %fp2,%fp0 | ...S(B2+T(B4+T(B6+TB8)))
335 |--FP3 RELEASED, RESTORE NOW AND TAKE SOME ADV 335 |--FP3 RELEASED, RESTORE NOW AND TAKE SOME ADVANTAGE OF HIDING
336 |--FP2 RELEASED. 336 |--FP2 RELEASED.
337 337
338 338
339 faddx %fp1,%fp0 339 faddx %fp1,%fp0
340 |--FP1 RELEASED 340 |--FP1 RELEASED
341 341
342 fmulx X(%a6),%fp0 342 fmulx X(%a6),%fp0
343 343
344 fmovel %d1,%FPCR 344 fmovel %d1,%FPCR |restore users exceptions
345 fadds POSNEG1(%a6),%fp0 345 fadds POSNEG1(%a6),%fp0 |last inst - possible exception set
346 bra t_frcinx 346 bra t_frcinx
347 347
348 348
349 SINBORS: 349 SINBORS:
350 |--IF |X| > 15PI, WE USE THE GENERAL ARGUMENT 350 |--IF |X| > 15PI, WE USE THE GENERAL ARGUMENT REDUCTION.
351 |--IF |X| < 2**(-40), RETURN X OR 1. 351 |--IF |X| < 2**(-40), RETURN X OR 1.
352 cmpil #0x3FFF8000,%d0 352 cmpil #0x3FFF8000,%d0
353 bgts REDUCEX 353 bgts REDUCEX
354 354
355 355
356 SINSM: 356 SINSM:
357 movel ADJN(%a6),%d0 357 movel ADJN(%a6),%d0
358 cmpil #0,%d0 358 cmpil #0,%d0
359 bgts COSTINY 359 bgts COSTINY
360 360
361 SINTINY: 361 SINTINY:
362 movew #0x0000,XDCARE(%a6) 362 movew #0x0000,XDCARE(%a6) | ...JUST IN CASE
363 fmovel %d1,%FPCR 363 fmovel %d1,%FPCR |restore users exceptions
364 fmovex X(%a6),%fp0 364 fmovex X(%a6),%fp0 |last inst - possible exception set
365 bra t_frcinx 365 bra t_frcinx
366 366
367 367
368 COSTINY: 368 COSTINY:
369 fmoves #0x3F800000,%fp0 369 fmoves #0x3F800000,%fp0
370 370
371 fmovel %d1,%FPCR 371 fmovel %d1,%FPCR |restore users exceptions
372 fsubs #0x00800000,%fp0 372 fsubs #0x00800000,%fp0 |last inst - possible exception set
373 bra t_frcinx 373 bra t_frcinx
374 374
375 375
376 REDUCEX: 376 REDUCEX:
377 |--WHEN REDUCEX IS USED, THE CODE WILL INEVITA 377 |--WHEN REDUCEX IS USED, THE CODE WILL INEVITABLY BE SLOW.
378 |--THIS REDUCTION METHOD, HOWEVER, IS MUCH FAS 378 |--THIS REDUCTION METHOD, HOWEVER, IS MUCH FASTER THAN USING
379 |--THE REMAINDER INSTRUCTION WHICH IS NOW IN S 379 |--THE REMAINDER INSTRUCTION WHICH IS NOW IN SOFTWARE.
380 380
381 fmovemx %fp2-%fp5,-(%a7) | ...s 381 fmovemx %fp2-%fp5,-(%a7) | ...save FP2 through FP5
382 movel %d2,-(%a7) 382 movel %d2,-(%a7)
383 fmoves #0x00000000,%fp1 383 fmoves #0x00000000,%fp1
384 |--If compact form of abs(arg) in d0=$7ffeffff 384 |--If compact form of abs(arg) in d0=$7ffeffff, argument is so large that
385 |--there is a danger of unwanted overflow in f 385 |--there is a danger of unwanted overflow in first LOOP iteration. In this
386 |--case, reduce argument by one remainder step 386 |--case, reduce argument by one remainder step to make subsequent reduction
387 |--safe. 387 |--safe.
388 cmpil #0x7ffeffff,%d0 |is ar 388 cmpil #0x7ffeffff,%d0 |is argument dangerously large?
389 bnes LOOP 389 bnes LOOP
390 movel #0x7ffe0000,FP_SCR2(%a6) 390 movel #0x7ffe0000,FP_SCR2(%a6) |yes
391 | ;creat 391 | ;create 2**16383*PI/2
392 movel #0xc90fdaa2,FP_SCR2+4(%a6) 392 movel #0xc90fdaa2,FP_SCR2+4(%a6)
393 clrl FP_SCR2+8(%a6) 393 clrl FP_SCR2+8(%a6)
394 ftstx %fp0 |test 394 ftstx %fp0 |test sign of argument
395 movel #0x7fdc0000,FP_SCR3(%a6) 395 movel #0x7fdc0000,FP_SCR3(%a6) |create low half of 2**16383*
396 | ;PI/2 396 | ;PI/2 at FP_SCR3
397 movel #0x85a308d3,FP_SCR3+4(%a6) 397 movel #0x85a308d3,FP_SCR3+4(%a6)
398 clrl FP_SCR3+8(%a6) 398 clrl FP_SCR3+8(%a6)
399 fblt red_neg 399 fblt red_neg
400 orw #0x8000,FP_SCR2(%a6) |posit 400 orw #0x8000,FP_SCR2(%a6) |positive arg
401 orw #0x8000,FP_SCR3(%a6) 401 orw #0x8000,FP_SCR3(%a6)
402 red_neg: 402 red_neg:
403 faddx FP_SCR2(%a6),%fp0 403 faddx FP_SCR2(%a6),%fp0 |high part of reduction is exact
404 fmovex %fp0,%fp1 |save 404 fmovex %fp0,%fp1 |save high result in fp1
405 faddx FP_SCR3(%a6),%fp0 405 faddx FP_SCR3(%a6),%fp0 |low part of reduction
406 fsubx %fp0,%fp1 406 fsubx %fp0,%fp1 |determine low component of result
407 faddx FP_SCR3(%a6),%fp1 407 faddx FP_SCR3(%a6),%fp1 |fp0/fp1 are reduced argument.
408 408
409 |--ON ENTRY, FP0 IS X, ON RETURN, FP0 IS X REM 409 |--ON ENTRY, FP0 IS X, ON RETURN, FP0 IS X REM PI/2, |X| <= PI/4.
410 |--integer quotient will be stored in N 410 |--integer quotient will be stored in N
411 |--Intermediate remainder is 66-bit long; (R,r 411 |--Intermediate remainder is 66-bit long; (R,r) in (FP0,FP1)
412 412
413 LOOP: 413 LOOP:
414 fmovex %fp0,INARG(%a6) | ...+ 414 fmovex %fp0,INARG(%a6) | ...+-2**K * F, 1 <= F < 2
415 movew INARG(%a6),%d0 415 movew INARG(%a6),%d0
416 movel %d0,%a1 | ...s 416 movel %d0,%a1 | ...save a copy of D0
417 andil #0x00007FFF,%d0 417 andil #0x00007FFF,%d0
418 subil #0x00003FFF,%d0 | ...D 418 subil #0x00003FFF,%d0 | ...D0 IS K
419 cmpil #28,%d0 419 cmpil #28,%d0
420 bles LASTLOOP 420 bles LASTLOOP
421 CONTLOOP: 421 CONTLOOP:
422 subil #27,%d0 | ...D0 IS L 422 subil #27,%d0 | ...D0 IS L := K-27
423 movel #0,ENDFLAG(%a6) 423 movel #0,ENDFLAG(%a6)
424 bras WORK 424 bras WORK
425 LASTLOOP: 425 LASTLOOP:
426 clrl %d0 | ...D 426 clrl %d0 | ...D0 IS L := 0
427 movel #1,ENDFLAG(%a6) 427 movel #1,ENDFLAG(%a6)
428 428
429 WORK: 429 WORK:
430 |--FIND THE REMAINDER OF (R,r) W.R.T. 2**L * 430 |--FIND THE REMAINDER OF (R,r) W.R.T. 2**L * (PI/2). L IS SO CHOSEN
431 |--THAT INT( X * (2/PI) / 2**(L) ) < 2**29. 431 |--THAT INT( X * (2/PI) / 2**(L) ) < 2**29.
432 432
433 |--CREATE 2**(-L) * (2/PI), SIGN(INARG)*2**(63 433 |--CREATE 2**(-L) * (2/PI), SIGN(INARG)*2**(63),
434 |--2**L * (PIby2_1), 2**L * (PIby2_2) 434 |--2**L * (PIby2_1), 2**L * (PIby2_2)
435 435
436 movel #0x00003FFE,%d2 | ...B 436 movel #0x00003FFE,%d2 | ...BIASED EXPO OF 2/PI
437 subl %d0,%d2 | ...B 437 subl %d0,%d2 | ...BIASED EXPO OF 2**(-L)*(2/PI)
438 438
439 movel #0xA2F9836E,FP_SCR1+4( 439 movel #0xA2F9836E,FP_SCR1+4(%a6)
440 movel #0x4E44152A,FP_SCR1+8( 440 movel #0x4E44152A,FP_SCR1+8(%a6)
441 movew %d2,FP_SCR1(%a6) 441 movew %d2,FP_SCR1(%a6) | ...FP_SCR1 is 2**(-L)*(2/PI)
442 442
443 fmovex %fp0,%fp2 443 fmovex %fp0,%fp2
444 fmulx FP_SCR1(%a6),%fp2 444 fmulx FP_SCR1(%a6),%fp2
445 |--WE MUST NOW FIND INT(FP2). SINCE WE NEED TH 445 |--WE MUST NOW FIND INT(FP2). SINCE WE NEED THIS VALUE IN
446 |--FLOATING POINT FORMAT, THE TWO FMOVE'S 446 |--FLOATING POINT FORMAT, THE TWO FMOVE'S FMOVE.L FP <--> N
447 |--WILL BE TOO INEFFICIENT. THE WAY AROUND IT 447 |--WILL BE TOO INEFFICIENT. THE WAY AROUND IT IS THAT
448 |--(SIGN(INARG)*2**63 + FP2) - SIGN(IN 448 |--(SIGN(INARG)*2**63 + FP2) - SIGN(INARG)*2**63 WILL GIVE
449 |--US THE DESIRED VALUE IN FLOATING POINT. 449 |--US THE DESIRED VALUE IN FLOATING POINT.
450 450
451 |--HIDE SIX CYCLES OF INSTRUCTION 451 |--HIDE SIX CYCLES OF INSTRUCTION
452 movel %a1,%d2 452 movel %a1,%d2
453 swap %d2 453 swap %d2
454 andil #0x80000000,%d2 454 andil #0x80000000,%d2
455 oril #0x5F000000,%d2 | ...D 455 oril #0x5F000000,%d2 | ...D2 IS SIGN(INARG)*2**63 IN SGL
456 movel %d2,TWOTO63(%a6) 456 movel %d2,TWOTO63(%a6)
457 457
458 movel %d0,%d2 458 movel %d0,%d2
459 addil #0x00003FFF,%d2 | ...B 459 addil #0x00003FFF,%d2 | ...BIASED EXPO OF 2**L * (PI/2)
460 460
461 |--FP2 IS READY 461 |--FP2 IS READY
462 fadds TWOTO63(%a6),%fp2 462 fadds TWOTO63(%a6),%fp2 | ...THE FRACTIONAL PART OF FP1 IS ROUNDED
463 463
464 |--HIDE 4 CYCLES OF INSTRUCTION; creating 2**( 464 |--HIDE 4 CYCLES OF INSTRUCTION; creating 2**(L)*Piby2_1 and 2**(L)*Piby2_2
465 movew %d2,FP_SCR2(%a6) 465 movew %d2,FP_SCR2(%a6)
466 clrw FP_SCR2+2(%a6) 466 clrw FP_SCR2+2(%a6)
467 movel #0xC90FDAA2,FP_SCR2+4( 467 movel #0xC90FDAA2,FP_SCR2+4(%a6)
468 clrl FP_SCR2+8(%a6) 468 clrl FP_SCR2+8(%a6) | ...FP_SCR2 is 2**(L) * Piby2_1
469 469
470 |--FP2 IS READY 470 |--FP2 IS READY
471 fsubs TWOTO63(%a6),%fp2 471 fsubs TWOTO63(%a6),%fp2 | ...FP2 is N
472 472
473 addil #0x00003FDD,%d0 473 addil #0x00003FDD,%d0
474 movew %d0,FP_SCR3(%a6) 474 movew %d0,FP_SCR3(%a6)
475 clrw FP_SCR3+2(%a6) 475 clrw FP_SCR3+2(%a6)
476 movel #0x85A308D3,FP_SCR3+4( 476 movel #0x85A308D3,FP_SCR3+4(%a6)
477 clrl FP_SCR3+8(%a6) 477 clrl FP_SCR3+8(%a6) | ...FP_SCR3 is 2**(L) * Piby2_2
478 478
479 movel ENDFLAG(%a6),%d0 479 movel ENDFLAG(%a6),%d0
480 480
481 |--We are now ready to perform (R+r) - N*P1 - 481 |--We are now ready to perform (R+r) - N*P1 - N*P2, P1 = 2**(L) * Piby2_1 and
482 |--P2 = 2**(L) * Piby2_2 482 |--P2 = 2**(L) * Piby2_2
483 fmovex %fp2,%fp4 483 fmovex %fp2,%fp4
484 fmulx FP_SCR2(%a6),%fp4 484 fmulx FP_SCR2(%a6),%fp4 | ...W = N*P1
485 fmovex %fp2,%fp5 485 fmovex %fp2,%fp5
486 fmulx FP_SCR3(%a6),%fp5 486 fmulx FP_SCR3(%a6),%fp5 | ...w = N*P2
487 fmovex %fp4,%fp3 487 fmovex %fp4,%fp3
488 |--we want P+p = W+w but |p| <= half ulp of 488 |--we want P+p = W+w but |p| <= half ulp of P
489 |--Then, we need to compute A := R-P and a 489 |--Then, we need to compute A := R-P and a := r-p
490 faddx %fp5,%fp3 490 faddx %fp5,%fp3 | ...FP3 is P
491 fsubx %fp3,%fp4 491 fsubx %fp3,%fp4 | ...W-P
492 492
493 fsubx %fp3,%fp0 493 fsubx %fp3,%fp0 | ...FP0 is A := R - P
494 faddx %fp5,%fp4 494 faddx %fp5,%fp4 | ...FP4 is p = (W-P)+w
495 495
496 fmovex %fp0,%fp3 496 fmovex %fp0,%fp3 | ...FP3 A
497 fsubx %fp4,%fp1 497 fsubx %fp4,%fp1 | ...FP1 is a := r - p
498 498
499 |--Now we need to normalize (A,a) to "new (R, 499 |--Now we need to normalize (A,a) to "new (R,r)" where R+r = A+a but
500 |--|r| <= half ulp of R. 500 |--|r| <= half ulp of R.
501 faddx %fp1,%fp0 501 faddx %fp1,%fp0 | ...FP0 is R := A+a
502 |--No need to calculate r if this is the last 502 |--No need to calculate r if this is the last loop
503 cmpil #0,%d0 503 cmpil #0,%d0
504 bgt RESTORE 504 bgt RESTORE
505 505
506 |--Need to calculate r 506 |--Need to calculate r
507 fsubx %fp0,%fp3 507 fsubx %fp0,%fp3 | ...A-R
508 faddx %fp3,%fp1 508 faddx %fp3,%fp1 | ...FP1 is r := (A-R)+a
509 bra LOOP 509 bra LOOP
510 510
511 RESTORE: 511 RESTORE:
512 fmovel %fp2,N(%a6) 512 fmovel %fp2,N(%a6)
513 movel (%a7)+,%d2 513 movel (%a7)+,%d2
514 fmovemx (%a7)+,%fp2-%fp5 514 fmovemx (%a7)+,%fp2-%fp5
515 515
516 516
517 movel ADJN(%a6),%d0 517 movel ADJN(%a6),%d0
518 cmpil #4,%d0 518 cmpil #4,%d0
519 519
520 blt SINCONT 520 blt SINCONT
521 bras SCCONT 521 bras SCCONT
522 522
523 .global ssincosd 523 .global ssincosd
524 ssincosd: 524 ssincosd:
525 |--SIN AND COS OF X FOR DENORMALIZED X 525 |--SIN AND COS OF X FOR DENORMALIZED X
526 526
527 fmoves #0x3F800000,%fp1 527 fmoves #0x3F800000,%fp1
528 bsr sto_cos |store 528 bsr sto_cos |store cosine result
529 bra t_extdnrm 529 bra t_extdnrm
530 530
531 .global ssincos 531 .global ssincos
532 ssincos: 532 ssincos:
533 |--SET ADJN TO 4 533 |--SET ADJN TO 4
534 movel #4,ADJN(%a6) 534 movel #4,ADJN(%a6)
535 535
536 fmovex (%a0),%fp0 | ...L 536 fmovex (%a0),%fp0 | ...LOAD INPUT
537 537
538 movel (%a0),%d0 538 movel (%a0),%d0
539 movew 4(%a0),%d0 539 movew 4(%a0),%d0
540 fmovex %fp0,X(%a6) 540 fmovex %fp0,X(%a6)
541 andil #0x7FFFFFFF,%d0 541 andil #0x7FFFFFFF,%d0 | ...COMPACTIFY X
542 542
543 cmpil #0x3FD78000,%d0 543 cmpil #0x3FD78000,%d0 | ...|X| >= 2**(-40)?
544 bges SCOK1 544 bges SCOK1
545 bra SCSM 545 bra SCSM
546 546
547 SCOK1: 547 SCOK1:
548 cmpil #0x4004BC7E,%d0 548 cmpil #0x4004BC7E,%d0 | ...|X| < 15 PI?
549 blts SCMAIN 549 blts SCMAIN
550 bra REDUCEX 550 bra REDUCEX
551 551
552 552
553 SCMAIN: 553 SCMAIN:
554 |--THIS IS THE USUAL CASE, |X| <= 15 PI. 554 |--THIS IS THE USUAL CASE, |X| <= 15 PI.
555 |--THE ARGUMENT REDUCTION IS DONE BY TABLE LOO 555 |--THE ARGUMENT REDUCTION IS DONE BY TABLE LOOK UP.
556 fmovex %fp0,%fp1 556 fmovex %fp0,%fp1
557 fmuld TWOBYPI,%fp1 | ...X 557 fmuld TWOBYPI,%fp1 | ...X*2/PI
558 558
559 |--HIDE THE NEXT THREE INSTRUCTIONS 559 |--HIDE THE NEXT THREE INSTRUCTIONS
560 lea PITBL+0x200,%a1 | ...T 560 lea PITBL+0x200,%a1 | ...TABLE OF N*PI/2, N = -32,...,32
561 561
562 562
563 |--FP1 IS NOW READY 563 |--FP1 IS NOW READY
564 fmovel %fp1,N(%a6) 564 fmovel %fp1,N(%a6) | ...CONVERT TO INTEGER
565 565
566 movel N(%a6),%d0 566 movel N(%a6),%d0
567 asll #4,%d0 567 asll #4,%d0
568 addal %d0,%a1 | ...A 568 addal %d0,%a1 | ...ADDRESS OF N*PIBY2, IN Y1, Y2
569 569
570 fsubx (%a1)+,%fp0 | ...X 570 fsubx (%a1)+,%fp0 | ...X-Y1
571 fsubs (%a1),%fp0 | ...F 571 fsubs (%a1),%fp0 | ...FP0 IS R = (X-Y1)-Y2
572 572
573 SCCONT: 573 SCCONT:
574 |--continuation point from REDUCEX 574 |--continuation point from REDUCEX
575 575
576 |--HIDE THE NEXT TWO 576 |--HIDE THE NEXT TWO
577 movel N(%a6),%d0 577 movel N(%a6),%d0
578 rorl #1,%d0 578 rorl #1,%d0
579 579
580 cmpil #0,%d0 | ...D 580 cmpil #0,%d0 | ...D0 < 0 IFF N IS ODD
581 bge NEVEN 581 bge NEVEN
582 582
583 NODD: 583 NODD:
584 |--REGISTERS SAVED SO FAR: D0, A0, FP2. 584 |--REGISTERS SAVED SO FAR: D0, A0, FP2.
585 585
586 fmovex %fp0,RPRIME(%a6) 586 fmovex %fp0,RPRIME(%a6)
587 fmulx %fp0,%fp0 | ... 587 fmulx %fp0,%fp0 | ...FP0 IS S = R*R
588 fmoved SINA7,%fp1 | ...A 588 fmoved SINA7,%fp1 | ...A7
589 fmoved COSB8,%fp2 | ...B 589 fmoved COSB8,%fp2 | ...B8
590 fmulx %fp0,%fp1 | ... 590 fmulx %fp0,%fp1 | ...SA7
591 movel %d2,-(%a7) 591 movel %d2,-(%a7)
592 movel %d0,%d2 592 movel %d0,%d2
593 fmulx %fp0,%fp2 | ... 593 fmulx %fp0,%fp2 | ...SB8
594 rorl #1,%d2 594 rorl #1,%d2
595 andil #0x80000000,%d2 595 andil #0x80000000,%d2
596 596
597 faddd SINA6,%fp1 | ...A 597 faddd SINA6,%fp1 | ...A6+SA7
598 eorl %d0,%d2 598 eorl %d0,%d2
599 andil #0x80000000,%d2 599 andil #0x80000000,%d2
600 faddd COSB7,%fp2 | ...B 600 faddd COSB7,%fp2 | ...B7+SB8
601 601
602 fmulx %fp0,%fp1 | ... 602 fmulx %fp0,%fp1 | ...S(A6+SA7)
603 eorl %d2,RPRIME(%a6) 603 eorl %d2,RPRIME(%a6)
604 movel (%a7)+,%d2 604 movel (%a7)+,%d2
605 fmulx %fp0,%fp2 | ... 605 fmulx %fp0,%fp2 | ...S(B7+SB8)
606 rorl #1,%d0 606 rorl #1,%d0
607 andil #0x80000000,%d0 607 andil #0x80000000,%d0
608 608
609 faddd SINA5,%fp1 | ...A 609 faddd SINA5,%fp1 | ...A5+S(A6+SA7)
610 movel #0x3F800000,POSNEG1(%a 610 movel #0x3F800000,POSNEG1(%a6)
611 eorl %d0,POSNEG1(%a6) 611 eorl %d0,POSNEG1(%a6)
612 faddd COSB6,%fp2 | ...B 612 faddd COSB6,%fp2 | ...B6+S(B7+SB8)
613 613
614 fmulx %fp0,%fp1 | ... 614 fmulx %fp0,%fp1 | ...S(A5+S(A6+SA7))
615 fmulx %fp0,%fp2 | ... 615 fmulx %fp0,%fp2 | ...S(B6+S(B7+SB8))
616 fmovex %fp0,SPRIME(%a6) 616 fmovex %fp0,SPRIME(%a6)
617 617
618 faddd SINA4,%fp1 | ...A 618 faddd SINA4,%fp1 | ...A4+S(A5+S(A6+SA7))
619 eorl %d0,SPRIME(%a6) 619 eorl %d0,SPRIME(%a6)
620 faddd COSB5,%fp2 | ...B 620 faddd COSB5,%fp2 | ...B5+S(B6+S(B7+SB8))
621 621
622 fmulx %fp0,%fp1 | ... 622 fmulx %fp0,%fp1 | ...S(A4+...)
623 fmulx %fp0,%fp2 | ... 623 fmulx %fp0,%fp2 | ...S(B5+...)
624 624
625 faddd SINA3,%fp1 | ...A 625 faddd SINA3,%fp1 | ...A3+S(A4+...)
626 faddd COSB4,%fp2 | ...B 626 faddd COSB4,%fp2 | ...B4+S(B5+...)
627 627
628 fmulx %fp0,%fp1 | ... 628 fmulx %fp0,%fp1 | ...S(A3+...)
629 fmulx %fp0,%fp2 | ... 629 fmulx %fp0,%fp2 | ...S(B4+...)
630 630
631 faddx SINA2,%fp1 | ...A 631 faddx SINA2,%fp1 | ...A2+S(A3+...)
632 faddx COSB3,%fp2 | ...B 632 faddx COSB3,%fp2 | ...B3+S(B4+...)
633 633
634 fmulx %fp0,%fp1 | ... 634 fmulx %fp0,%fp1 | ...S(A2+...)
635 fmulx %fp0,%fp2 | ... 635 fmulx %fp0,%fp2 | ...S(B3+...)
636 636
637 faddx SINA1,%fp1 | ...A 637 faddx SINA1,%fp1 | ...A1+S(A2+...)
638 faddx COSB2,%fp2 | ...B 638 faddx COSB2,%fp2 | ...B2+S(B3+...)
639 639
640 fmulx %fp0,%fp1 | ... 640 fmulx %fp0,%fp1 | ...S(A1+...)
641 fmulx %fp2,%fp0 | ... 641 fmulx %fp2,%fp0 | ...S(B2+...)
642 642
643 643
644 644
645 fmulx RPRIME(%a6),%fp1 645 fmulx RPRIME(%a6),%fp1 | ...R'S(A1+...)
646 fadds COSB1,%fp0 | ...B 646 fadds COSB1,%fp0 | ...B1+S(B2...)
647 fmulx SPRIME(%a6),%fp0 647 fmulx SPRIME(%a6),%fp0 | ...S'(B1+S(B2+...))
648 648
649 movel %d1,-(%sp) |resto 649 movel %d1,-(%sp) |restore users mode & precision
650 andil #0xff,%d1 650 andil #0xff,%d1 |mask off all exceptions
651 fmovel %d1,%FPCR 651 fmovel %d1,%FPCR
652 faddx RPRIME(%a6),%fp1 652 faddx RPRIME(%a6),%fp1 | ...COS(X)
653 bsr sto_cos |store 653 bsr sto_cos |store cosine result
654 fmovel (%sp)+,%FPCR |resto 654 fmovel (%sp)+,%FPCR |restore users exceptions
655 fadds POSNEG1(%a6),%fp0 655 fadds POSNEG1(%a6),%fp0 | ...SIN(X)
656 656
657 bra t_frcinx 657 bra t_frcinx
658 658
659 659
660 NEVEN: 660 NEVEN:
661 |--REGISTERS SAVED SO FAR: FP2. 661 |--REGISTERS SAVED SO FAR: FP2.
662 662
663 fmovex %fp0,RPRIME(%a6) 663 fmovex %fp0,RPRIME(%a6)
664 fmulx %fp0,%fp0 | ... 664 fmulx %fp0,%fp0 | ...FP0 IS S = R*R
665 fmoved COSB8,%fp1 665 fmoved COSB8,%fp1 | ...B8
666 fmoved SINA7,%fp2 666 fmoved SINA7,%fp2 | ...A7
667 fmulx %fp0,%fp1 | ... 667 fmulx %fp0,%fp1 | ...SB8
668 fmovex %fp0,SPRIME(%a6) 668 fmovex %fp0,SPRIME(%a6)
669 fmulx %fp0,%fp2 | ... 669 fmulx %fp0,%fp2 | ...SA7
670 rorl #1,%d0 670 rorl #1,%d0
671 andil #0x80000000,%d0 671 andil #0x80000000,%d0
672 faddd COSB7,%fp1 | ...B 672 faddd COSB7,%fp1 | ...B7+SB8
673 faddd SINA6,%fp2 | ...A 673 faddd SINA6,%fp2 | ...A6+SA7
674 eorl %d0,RPRIME(%a6) 674 eorl %d0,RPRIME(%a6)
675 eorl %d0,SPRIME(%a6) 675 eorl %d0,SPRIME(%a6)
676 fmulx %fp0,%fp1 | ... 676 fmulx %fp0,%fp1 | ...S(B7+SB8)
677 oril #0x3F800000,%d0 677 oril #0x3F800000,%d0
678 movel %d0,POSNEG1(%a6) 678 movel %d0,POSNEG1(%a6)
679 fmulx %fp0,%fp2 | ... 679 fmulx %fp0,%fp2 | ...S(A6+SA7)
680 680
681 faddd COSB6,%fp1 | ...B 681 faddd COSB6,%fp1 | ...B6+S(B7+SB8)
682 faddd SINA5,%fp2 | ...A 682 faddd SINA5,%fp2 | ...A5+S(A6+SA7)
683 683
684 fmulx %fp0,%fp1 | ... 684 fmulx %fp0,%fp1 | ...S(B6+S(B7+SB8))
685 fmulx %fp0,%fp2 | ... 685 fmulx %fp0,%fp2 | ...S(A5+S(A6+SA7))
686 686
687 faddd COSB5,%fp1 | ...B 687 faddd COSB5,%fp1 | ...B5+S(B6+S(B7+SB8))
688 faddd SINA4,%fp2 | ...A 688 faddd SINA4,%fp2 | ...A4+S(A5+S(A6+SA7))
689 689
690 fmulx %fp0,%fp1 | ... 690 fmulx %fp0,%fp1 | ...S(B5+...)
691 fmulx %fp0,%fp2 | ... 691 fmulx %fp0,%fp2 | ...S(A4+...)
692 692
693 faddd COSB4,%fp1 | ...B 693 faddd COSB4,%fp1 | ...B4+S(B5+...)
694 faddd SINA3,%fp2 | ...A 694 faddd SINA3,%fp2 | ...A3+S(A4+...)
695 695
696 fmulx %fp0,%fp1 | ... 696 fmulx %fp0,%fp1 | ...S(B4+...)
697 fmulx %fp0,%fp2 | ... 697 fmulx %fp0,%fp2 | ...S(A3+...)
698 698
699 faddx COSB3,%fp1 | ...B 699 faddx COSB3,%fp1 | ...B3+S(B4+...)
700 faddx SINA2,%fp2 | ...A 700 faddx SINA2,%fp2 | ...A2+S(A3+...)
701 701
702 fmulx %fp0,%fp1 | ... 702 fmulx %fp0,%fp1 | ...S(B3+...)
703 fmulx %fp0,%fp2 | ... 703 fmulx %fp0,%fp2 | ...S(A2+...)
704 704
705 faddx COSB2,%fp1 | ...B 705 faddx COSB2,%fp1 | ...B2+S(B3+...)
706 faddx SINA1,%fp2 | ...A 706 faddx SINA1,%fp2 | ...A1+S(A2+...)
707 707
708 fmulx %fp0,%fp1 | ... 708 fmulx %fp0,%fp1 | ...S(B2+...)
709 fmulx %fp2,%fp0 | ... 709 fmulx %fp2,%fp0 | ...s(a1+...)
710 710
711 711
712 712
713 fadds COSB1,%fp1 | ...B 713 fadds COSB1,%fp1 | ...B1+S(B2...)
714 fmulx RPRIME(%a6),%fp0 714 fmulx RPRIME(%a6),%fp0 | ...R'S(A1+...)
715 fmulx SPRIME(%a6),%fp1 715 fmulx SPRIME(%a6),%fp1 | ...S'(B1+S(B2+...))
716 716
717 movel %d1,-(%sp) |save 717 movel %d1,-(%sp) |save users mode & precision
718 andil #0xff,%d1 718 andil #0xff,%d1 |mask off all exceptions
719 fmovel %d1,%FPCR 719 fmovel %d1,%FPCR
720 fadds POSNEG1(%a6),%fp1 720 fadds POSNEG1(%a6),%fp1 | ...COS(X)
721 bsr sto_cos |store 721 bsr sto_cos |store cosine result
722 fmovel (%sp)+,%FPCR |resto 722 fmovel (%sp)+,%FPCR |restore users exceptions
723 faddx RPRIME(%a6),%fp0 723 faddx RPRIME(%a6),%fp0 | ...SIN(X)
724 724
725 bra t_frcinx 725 bra t_frcinx
726 726
727 SCBORS: 727 SCBORS:
728 cmpil #0x3FFF8000,%d0 728 cmpil #0x3FFF8000,%d0
729 bgt REDUCEX 729 bgt REDUCEX
730 730
731 731
732 SCSM: 732 SCSM:
733 movew #0x0000,XDCARE(%a6) 733 movew #0x0000,XDCARE(%a6)
734 fmoves #0x3F800000,%fp1 734 fmoves #0x3F800000,%fp1
735 735
736 movel %d1,-(%sp) |save 736 movel %d1,-(%sp) |save users mode & precision
737 andil #0xff,%d1 737 andil #0xff,%d1 |mask off all exceptions
738 fmovel %d1,%FPCR 738 fmovel %d1,%FPCR
739 fsubs #0x00800000,%fp1 739 fsubs #0x00800000,%fp1
740 bsr sto_cos |store 740 bsr sto_cos |store cosine result
741 fmovel (%sp)+,%FPCR |resto 741 fmovel (%sp)+,%FPCR |restore users exceptions
742 fmovex X(%a6),%fp0 742 fmovex X(%a6),%fp0
743 bra t_frcinx 743 bra t_frcinx
744 744
745 |end 745 |end
Linux® is a registered trademark of Linus Torvalds in the United States and other countries.
TOMOYO® is a registered trademark of NTT DATA CORPORATION.