Annotation of sys/arch/m68k/fpsp/stwotox.sa, Revision 1.1.1.1
1.1 nbrk 1: * $OpenBSD: stwotox.sa,v 1.2 1996/05/29 21:05:44 niklas Exp $
2: * $NetBSD: stwotox.sa,v 1.3 1994/10/26 07:50:15 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: * stwotox.sa 3.1 12/10/90
36: *
37: * stwotox --- 2**X
38: * stwotoxd --- 2**X for denormalized X
39: * stentox --- 10**X
40: * stentoxd --- 10**X for denormalized X
41: *
42: * Input: Double-extended number X in location pointed to
43: * by address register a0.
44: *
45: * Output: The function values are returned in Fp0.
46: *
47: * Accuracy and Monotonicity: The returned result is within 2 ulps in
48: * 64 significant bit, i.e. within 0.5001 ulp to 53 bits if the
49: * result is subsequently rounded to double precision. The
50: * result is provably monotonic in double precision.
51: *
52: * Speed: The program stwotox takes approximately 190 cycles and the
53: * program stentox takes approximately 200 cycles.
54: *
55: * Algorithm:
56: *
57: * twotox
58: * 1. If |X| > 16480, go to ExpBig.
59: *
60: * 2. If |X| < 2**(-70), go to ExpSm.
61: *
62: * 3. Decompose X as X = N/64 + r where |r| <= 1/128. Furthermore
63: * decompose N as
64: * N = 64(M + M') + j, j = 0,1,2,...,63.
65: *
66: * 4. Overwrite r := r * log2. Then
67: * 2**X = 2**(M') * 2**(M) * 2**(j/64) * exp(r).
68: * Go to expr to compute that expression.
69: *
70: * tentox
71: * 1. If |X| > 16480*log_10(2) (base 10 log of 2), go to ExpBig.
72: *
73: * 2. If |X| < 2**(-70), go to ExpSm.
74: *
75: * 3. Set y := X*log_2(10)*64 (base 2 log of 10). Set
76: * N := round-to-int(y). Decompose N as
77: * N = 64(M + M') + j, j = 0,1,2,...,63.
78: *
79: * 4. Define r as
80: * r := ((X - N*L1)-N*L2) * L10
81: * where L1, L2 are the leading and trailing parts of log_10(2)/64
82: * and L10 is the natural log of 10. Then
83: * 10**X = 2**(M') * 2**(M) * 2**(j/64) * exp(r).
84: * Go to expr to compute that expression.
85: *
86: * expr
87: * 1. Fetch 2**(j/64) from table as Fact1 and Fact2.
88: *
89: * 2. Overwrite Fact1 and Fact2 by
90: * Fact1 := 2**(M) * Fact1
91: * Fact2 := 2**(M) * Fact2
92: * Thus Fact1 + Fact2 = 2**(M) * 2**(j/64).
93: *
94: * 3. Calculate P where 1 + P approximates exp(r):
95: * P = r + r*r*(A1+r*(A2+...+r*A5)).
96: *
97: * 4. Let AdjFact := 2**(M'). Return
98: * AdjFact * ( Fact1 + ((Fact1*P) + Fact2) ).
99: * Exit.
100: *
101: * ExpBig
102: * 1. Generate overflow by Huge * Huge if X > 0; otherwise, generate
103: * underflow by Tiny * Tiny.
104: *
105: * ExpSm
106: * 1. Return 1 + X.
107: *
108:
109: STWOTOX IDNT 2,1 Motorola 040 Floating Point Software Package
110:
111: section 8
112:
113: include fpsp.h
114:
115: BOUNDS1 DC.L $3FB98000,$400D80C0 ... 2^(-70),16480
116: BOUNDS2 DC.L $3FB98000,$400B9B07 ... 2^(-70),16480 LOG2/LOG10
117:
118: L2TEN64 DC.L $406A934F,$0979A371 ... 64LOG10/LOG2
119: L10TWO1 DC.L $3F734413,$509F8000 ... LOG2/64LOG10
120:
121: L10TWO2 DC.L $BFCD0000,$C0219DC1,$DA994FD2,$00000000
122:
123: LOG10 DC.L $40000000,$935D8DDD,$AAA8AC17,$00000000
124:
125: LOG2 DC.L $3FFE0000,$B17217F7,$D1CF79AC,$00000000
126:
127: EXPA5 DC.L $3F56C16D,$6F7BD0B2
128: EXPA4 DC.L $3F811112,$302C712C
129: EXPA3 DC.L $3FA55555,$55554CC1
130: EXPA2 DC.L $3FC55555,$55554A54
131: EXPA1 DC.L $3FE00000,$00000000,$00000000,$00000000
132:
133: HUGE DC.L $7FFE0000,$FFFFFFFF,$FFFFFFFF,$00000000
134: TINY DC.L $00010000,$FFFFFFFF,$FFFFFFFF,$00000000
135:
136: EXPTBL
137: DC.L $3FFF0000,$80000000,$00000000,$3F738000
138: DC.L $3FFF0000,$8164D1F3,$BC030773,$3FBEF7CA
139: DC.L $3FFF0000,$82CD8698,$AC2BA1D7,$3FBDF8A9
140: DC.L $3FFF0000,$843A28C3,$ACDE4046,$3FBCD7C9
141: DC.L $3FFF0000,$85AAC367,$CC487B15,$BFBDE8DA
142: DC.L $3FFF0000,$871F6196,$9E8D1010,$3FBDE85C
143: DC.L $3FFF0000,$88980E80,$92DA8527,$3FBEBBF1
144: DC.L $3FFF0000,$8A14D575,$496EFD9A,$3FBB80CA
145: DC.L $3FFF0000,$8B95C1E3,$EA8BD6E7,$BFBA8373
146: DC.L $3FFF0000,$8D1ADF5B,$7E5BA9E6,$BFBE9670
147: DC.L $3FFF0000,$8EA4398B,$45CD53C0,$3FBDB700
148: DC.L $3FFF0000,$9031DC43,$1466B1DC,$3FBEEEB0
149: DC.L $3FFF0000,$91C3D373,$AB11C336,$3FBBFD6D
150: DC.L $3FFF0000,$935A2B2F,$13E6E92C,$BFBDB319
151: DC.L $3FFF0000,$94F4EFA8,$FEF70961,$3FBDBA2B
152: DC.L $3FFF0000,$96942D37,$20185A00,$3FBE91D5
153: DC.L $3FFF0000,$9837F051,$8DB8A96F,$3FBE8D5A
154: DC.L $3FFF0000,$99E04593,$20B7FA65,$BFBCDE7B
155: DC.L $3FFF0000,$9B8D39B9,$D54E5539,$BFBEBAAF
156: DC.L $3FFF0000,$9D3ED9A7,$2CFFB751,$BFBD86DA
157: DC.L $3FFF0000,$9EF53260,$91A111AE,$BFBEBEDD
158: DC.L $3FFF0000,$A0B0510F,$B9714FC2,$3FBCC96E
159: DC.L $3FFF0000,$A2704303,$0C496819,$BFBEC90B
160: DC.L $3FFF0000,$A43515AE,$09E6809E,$3FBBD1DB
161: DC.L $3FFF0000,$A5FED6A9,$B15138EA,$3FBCE5EB
162: DC.L $3FFF0000,$A7CD93B4,$E965356A,$BFBEC274
163: DC.L $3FFF0000,$A9A15AB4,$EA7C0EF8,$3FBEA83C
164: DC.L $3FFF0000,$AB7A39B5,$A93ED337,$3FBECB00
165: DC.L $3FFF0000,$AD583EEA,$42A14AC6,$3FBE9301
166: DC.L $3FFF0000,$AF3B78AD,$690A4375,$BFBD8367
167: DC.L $3FFF0000,$B123F581,$D2AC2590,$BFBEF05F
168: DC.L $3FFF0000,$B311C412,$A9112489,$3FBDFB3C
169: DC.L $3FFF0000,$B504F333,$F9DE6484,$3FBEB2FB
170: DC.L $3FFF0000,$B6FD91E3,$28D17791,$3FBAE2CB
171: DC.L $3FFF0000,$B8FBAF47,$62FB9EE9,$3FBCDC3C
172: DC.L $3FFF0000,$BAFF5AB2,$133E45FB,$3FBEE9AA
173: DC.L $3FFF0000,$BD08A39F,$580C36BF,$BFBEAEFD
174: DC.L $3FFF0000,$BF1799B6,$7A731083,$BFBCBF51
175: DC.L $3FFF0000,$C12C4CCA,$66709456,$3FBEF88A
176: DC.L $3FFF0000,$C346CCDA,$24976407,$3FBD83B2
177: DC.L $3FFF0000,$C5672A11,$5506DADD,$3FBDF8AB
178: DC.L $3FFF0000,$C78D74C8,$ABB9B15D,$BFBDFB17
179: DC.L $3FFF0000,$C9B9BD86,$6E2F27A3,$BFBEFE3C
180: DC.L $3FFF0000,$CBEC14FE,$F2727C5D,$BFBBB6F8
181: DC.L $3FFF0000,$CE248C15,$1F8480E4,$BFBCEE53
182: DC.L $3FFF0000,$D06333DA,$EF2B2595,$BFBDA4AE
183: DC.L $3FFF0000,$D2A81D91,$F12AE45A,$3FBC9124
184: DC.L $3FFF0000,$D4F35AAB,$CFEDFA1F,$3FBEB243
185: DC.L $3FFF0000,$D744FCCA,$D69D6AF4,$3FBDE69A
186: DC.L $3FFF0000,$D99D15C2,$78AFD7B6,$BFB8BC61
187: DC.L $3FFF0000,$DBFBB797,$DAF23755,$3FBDF610
188: DC.L $3FFF0000,$DE60F482,$5E0E9124,$BFBD8BE1
189: DC.L $3FFF0000,$E0CCDEEC,$2A94E111,$3FBACB12
190: DC.L $3FFF0000,$E33F8972,$BE8A5A51,$3FBB9BFE
191: DC.L $3FFF0000,$E5B906E7,$7C8348A8,$3FBCF2F4
192: DC.L $3FFF0000,$E8396A50,$3C4BDC68,$3FBEF22F
193: DC.L $3FFF0000,$EAC0C6E7,$DD24392F,$BFBDBF4A
194: DC.L $3FFF0000,$ED4F301E,$D9942B84,$3FBEC01A
195: DC.L $3FFF0000,$EFE4B99B,$DCDAF5CB,$3FBE8CAC
196: DC.L $3FFF0000,$F281773C,$59FFB13A,$BFBCBB3F
197: DC.L $3FFF0000,$F5257D15,$2486CC2C,$3FBEF73A
198: DC.L $3FFF0000,$F7D0DF73,$0AD13BB9,$BFB8B795
199: DC.L $3FFF0000,$FA83B2DB,$722A033A,$3FBEF84B
200: DC.L $3FFF0000,$FD3E0C0C,$F486C175,$BFBEF581
201:
202: N equ L_SCR1
203:
204: X equ FP_SCR1
205: XDCARE equ X+2
206: XFRAC equ X+4
207:
208: ADJFACT equ FP_SCR2
209:
210: FACT1 equ FP_SCR3
211: FACT1HI equ FACT1+4
212: FACT1LOW equ FACT1+8
213:
214: FACT2 equ FP_SCR4
215: FACT2HI equ FACT2+4
216: FACT2LOW equ FACT2+8
217:
218: xref t_unfl
219: xref t_ovfl
220: xref t_frcinx
221:
222: xdef stwotoxd
223: stwotoxd:
224: *--ENTRY POINT FOR 2**(X) FOR DENORMALIZED ARGUMENT
225:
226: fmove.l d1,fpcr ...set user's rounding mode/precision
227: Fmove.S #:3F800000,FP0 ...RETURN 1 + X
228: move.l (a0),d0
229: or.l #$00800001,d0
230: fadd.s d0,fp0
231: bra t_frcinx
232:
233: xdef stwotox
234: stwotox:
235: *--ENTRY POINT FOR 2**(X), HERE X IS FINITE, NON-ZERO, AND NOT NAN'S
236: FMOVEM.X (a0),FP0 ...LOAD INPUT, do not set cc's
237:
238: MOVE.L (A0),D0
239: MOVE.W 4(A0),D0
240: FMOVE.X FP0,X(a6)
241: ANDI.L #$7FFFFFFF,D0
242:
243: CMPI.L #$3FB98000,D0 ...|X| >= 2**(-70)?
244: BGE.B TWOOK1
245: BRA.W EXPBORS
246:
247: TWOOK1:
248: CMPI.L #$400D80C0,D0 ...|X| > 16480?
249: BLE.B TWOMAIN
250: BRA.W EXPBORS
251:
252:
253: TWOMAIN:
254: *--USUAL CASE, 2^(-70) <= |X| <= 16480
255:
256: FMOVE.X FP0,FP1
257: FMUL.S #:42800000,FP1 ...64 * X
258:
259: FMOVE.L FP1,N(a6) ...N = ROUND-TO-INT(64 X)
260: MOVE.L d2,-(sp)
261: LEA EXPTBL,a1 ...LOAD ADDRESS OF TABLE OF 2^(J/64)
262: FMOVE.L N(a6),FP1 ...N --> FLOATING FMT
263: MOVE.L N(a6),D0
264: MOVE.L D0,d2
265: ANDI.L #$3F,D0 ...D0 IS J
266: ASL.L #4,D0 ...DISPLACEMENT FOR 2^(J/64)
267: ADDA.L D0,a1 ...ADDRESS FOR 2^(J/64)
268: ASR.L #6,d2 ...d2 IS L, N = 64L + J
269: MOVE.L d2,D0
270: ASR.L #1,D0 ...D0 IS M
271: SUB.L D0,d2 ...d2 IS M', N = 64(M+M') + J
272: ADDI.L #$3FFF,d2
273: MOVE.W d2,ADJFACT(a6) ...ADJFACT IS 2^(M')
274: MOVE.L (sp)+,d2
275: *--SUMMARY: a1 IS ADDRESS FOR THE LEADING PORTION OF 2^(J/64),
276: *--D0 IS M WHERE N = 64(M+M') + J. NOTE THAT |M| <= 16140 BY DESIGN.
277: *--ADJFACT = 2^(M').
278: *--REGISTERS SAVED SO FAR ARE (IN ORDER) FPCR, D0, FP1, a1, AND FP2.
279:
280: FMUL.S #:3C800000,FP1 ...(1/64)*N
281: MOVE.L (a1)+,FACT1(a6)
282: MOVE.L (a1)+,FACT1HI(a6)
283: MOVE.L (a1)+,FACT1LOW(a6)
284: MOVE.W (a1)+,FACT2(a6)
285: clr.w FACT2+2(a6)
286:
287: FSUB.X FP1,FP0 ...X - (1/64)*INT(64 X)
288:
289: MOVE.W (a1)+,FACT2HI(a6)
290: clr.w FACT2HI+2(a6)
291: clr.l FACT2LOW(a6)
292: ADD.W D0,FACT1(a6)
293:
294: FMUL.X LOG2,FP0 ...FP0 IS R
295: ADD.W D0,FACT2(a6)
296:
297: BRA.W expr
298:
299: EXPBORS:
300: *--FPCR, D0 SAVED
301: CMPI.L #$3FFF8000,D0
302: BGT.B EXPBIG
303:
304: EXPSM:
305: *--|X| IS SMALL, RETURN 1 + X
306:
307: FMOVE.L d1,FPCR ;restore users exceptions
308: FADD.S #:3F800000,FP0 ...RETURN 1 + X
309:
310: bra t_frcinx
311:
312: EXPBIG:
313: *--|X| IS LARGE, GENERATE OVERFLOW IF X > 0; ELSE GENERATE UNDERFLOW
314: *--REGISTERS SAVE SO FAR ARE FPCR AND D0
315: MOVE.L X(a6),D0
316: TST.L D0
317: BLT.B EXPNEG
318:
319: bclr.b #7,(a0) ;t_ovfl expects positive value
320: bra t_ovfl
321:
322: EXPNEG:
323: bclr.b #7,(a0) ;t_unfl expects positive value
324: bra t_unfl
325:
326: xdef stentoxd
327: stentoxd:
328: *--ENTRY POINT FOR 10**(X) FOR DENORMALIZED ARGUMENT
329:
330: fmove.l d1,fpcr ...set user's rounding mode/precision
331: Fmove.S #:3F800000,FP0 ...RETURN 1 + X
332: move.l (a0),d0
333: or.l #$00800001,d0
334: fadd.s d0,fp0
335: bra t_frcinx
336:
337: xdef stentox
338: stentox:
339: *--ENTRY POINT FOR 10**(X), HERE X IS FINITE, NON-ZERO, AND NOT NAN'S
340: FMOVEM.X (a0),FP0 ...LOAD INPUT, do not set cc's
341:
342: MOVE.L (A0),D0
343: MOVE.W 4(A0),D0
344: FMOVE.X FP0,X(a6)
345: ANDI.L #$7FFFFFFF,D0
346:
347: CMPI.L #$3FB98000,D0 ...|X| >= 2**(-70)?
348: BGE.B TENOK1
349: BRA.W EXPBORS
350:
351: TENOK1:
352: CMPI.L #$400B9B07,D0 ...|X| <= 16480*log2/log10 ?
353: BLE.B TENMAIN
354: BRA.W EXPBORS
355:
356: TENMAIN:
357: *--USUAL CASE, 2^(-70) <= |X| <= 16480 LOG 2 / LOG 10
358:
359: FMOVE.X FP0,FP1
360: FMUL.D L2TEN64,FP1 ...X*64*LOG10/LOG2
361:
362: FMOVE.L FP1,N(a6) ...N=INT(X*64*LOG10/LOG2)
363: MOVE.L d2,-(sp)
364: LEA EXPTBL,a1 ...LOAD ADDRESS OF TABLE OF 2^(J/64)
365: FMOVE.L N(a6),FP1 ...N --> FLOATING FMT
366: MOVE.L N(a6),D0
367: MOVE.L D0,d2
368: ANDI.L #$3F,D0 ...D0 IS J
369: ASL.L #4,D0 ...DISPLACEMENT FOR 2^(J/64)
370: ADDA.L D0,a1 ...ADDRESS FOR 2^(J/64)
371: ASR.L #6,d2 ...d2 IS L, N = 64L + J
372: MOVE.L d2,D0
373: ASR.L #1,D0 ...D0 IS M
374: SUB.L D0,d2 ...d2 IS M', N = 64(M+M') + J
375: ADDI.L #$3FFF,d2
376: MOVE.W d2,ADJFACT(a6) ...ADJFACT IS 2^(M')
377: MOVE.L (sp)+,d2
378:
379: *--SUMMARY: a1 IS ADDRESS FOR THE LEADING PORTION OF 2^(J/64),
380: *--D0 IS M WHERE N = 64(M+M') + J. NOTE THAT |M| <= 16140 BY DESIGN.
381: *--ADJFACT = 2^(M').
382: *--REGISTERS SAVED SO FAR ARE (IN ORDER) FPCR, D0, FP1, a1, AND FP2.
383:
384: FMOVE.X FP1,FP2
385:
386: FMUL.D L10TWO1,FP1 ...N*(LOG2/64LOG10)_LEAD
387: MOVE.L (a1)+,FACT1(a6)
388:
389: FMUL.X L10TWO2,FP2 ...N*(LOG2/64LOG10)_TRAIL
390:
391: MOVE.L (a1)+,FACT1HI(a6)
392: MOVE.L (a1)+,FACT1LOW(a6)
393: FSUB.X FP1,FP0 ...X - N L_LEAD
394: MOVE.W (a1)+,FACT2(a6)
395:
396: FSUB.X FP2,FP0 ...X - N L_TRAIL
397:
398: clr.w FACT2+2(a6)
399: MOVE.W (a1)+,FACT2HI(a6)
400: clr.w FACT2HI+2(a6)
401: clr.l FACT2LOW(a6)
402:
403: FMUL.X LOG10,FP0 ...FP0 IS R
404:
405: ADD.W D0,FACT1(a6)
406: ADD.W D0,FACT2(a6)
407:
408: expr:
409: *--FPCR, FP2, FP3 ARE SAVED IN ORDER AS SHOWN.
410: *--ADJFACT CONTAINS 2**(M'), FACT1 + FACT2 = 2**(M) * 2**(J/64).
411: *--FP0 IS R. THE FOLLOWING CODE COMPUTES
412: *-- 2**(M'+M) * 2**(J/64) * EXP(R)
413:
414: FMOVE.X FP0,FP1
415: FMUL.X FP1,FP1 ...FP1 IS S = R*R
416:
417: FMOVE.D EXPA5,FP2 ...FP2 IS A5
418: FMOVE.D EXPA4,FP3 ...FP3 IS A4
419:
420: FMUL.X FP1,FP2 ...FP2 IS S*A5
421: FMUL.X FP1,FP3 ...FP3 IS S*A4
422:
423: FADD.D EXPA3,FP2 ...FP2 IS A3+S*A5
424: FADD.D EXPA2,FP3 ...FP3 IS A2+S*A4
425:
426: FMUL.X FP1,FP2 ...FP2 IS S*(A3+S*A5)
427: FMUL.X FP1,FP3 ...FP3 IS S*(A2+S*A4)
428:
429: FADD.D EXPA1,FP2 ...FP2 IS A1+S*(A3+S*A5)
430: FMUL.X FP0,FP3 ...FP3 IS R*S*(A2+S*A4)
431:
432: FMUL.X FP1,FP2 ...FP2 IS S*(A1+S*(A3+S*A5))
433: FADD.X FP3,FP0 ...FP0 IS R+R*S*(A2+S*A4)
434:
435: FADD.X FP2,FP0 ...FP0 IS EXP(R) - 1
436:
437:
438: *--FINAL RECONSTRUCTION PROCESS
439: *--EXP(X) = 2^M*2^(J/64) + 2^M*2^(J/64)*(EXP(R)-1) - (1 OR 0)
440:
441: FMUL.X FACT1(a6),FP0
442: FADD.X FACT2(a6),FP0
443: FADD.X FACT1(a6),FP0
444:
445: FMOVE.L d1,FPCR ;restore users exceptions
446: clr.w ADJFACT+2(a6)
447: move.l #$80000000,ADJFACT+4(a6)
448: clr.l ADJFACT+8(a6)
449: FMUL.X ADJFACT(a6),FP0 ...FINAL ADJUSTMENT
450:
451: bra t_frcinx
452:
453: end
CVSweb
![[BACK]](/icons/back.gif){kind=icon}
Return to stwotox.sa CVS log ![[TXT]](/icons/text.gif){kind=icon}
![[DIR]](/icons/dir.gif){kind=icon}
Up to [local] / sys / arch / m68k / fpsp