Annotation of sys/arch/m68k/fpsp/stanh.sa, Revision 1.1.1.1
1.1 nbrk 1: * $OpenBSD: stanh.sa,v 1.2 1996/05/29 21:05:43 niklas Exp $
2: * $NetBSD: stanh.sa,v 1.3 1994/10/26 07:50:12 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: * stanh.sa 3.1 12/10/90
36: *
37: * The entry point sTanh computes the hyperbolic tangent of
38: * an input argument; sTanhd does the same except for denormalized
39: * input.
40: *
41: * Input: Double-extended number X in location pointed to
42: * by address register a0.
43: *
44: * Output: The value tanh(X) returned in floating-point register Fp0.
45: *
46: * Accuracy and Monotonicity: The returned result is within 3 ulps in
47: * 64 significant bit, i.e. within 0.5001 ulp to 53 bits if the
48: * result is subsequently rounded to double precision. The
49: * result is provably monotonic in double precision.
50: *
51: * Speed: The program stanh takes approximately 270 cycles.
52: *
53: * Algorithm:
54: *
55: * TANH
56: * 1. If |X| >= (5/2) log2 or |X| <= 2**(-40), go to 3.
57: *
58: * 2. (2**(-40) < |X| < (5/2) log2) Calculate tanh(X) by
59: * sgn := sign(X), y := 2|X|, z := expm1(Y), and
60: * tanh(X) = sgn*( z/(2+z) ).
61: * Exit.
62: *
63: * 3. (|X| <= 2**(-40) or |X| >= (5/2) log2). If |X| < 1,
64: * go to 7.
65: *
66: * 4. (|X| >= (5/2) log2) If |X| >= 50 log2, go to 6.
67: *
68: * 5. ((5/2) log2 <= |X| < 50 log2) Calculate tanh(X) by
69: * sgn := sign(X), y := 2|X|, z := exp(Y),
70: * tanh(X) = sgn - [ sgn*2/(1+z) ].
71: * Exit.
72: *
73: * 6. (|X| >= 50 log2) Tanh(X) = +-1 (round to nearest). Thus, we
74: * calculate Tanh(X) by
75: * sgn := sign(X), Tiny := 2**(-126),
76: * tanh(X) := sgn - sgn*Tiny.
77: * Exit.
78: *
79: * 7. (|X| < 2**(-40)). Tanh(X) = X. Exit.
80: *
81:
82: STANH IDNT 2,1 Motorola 040 Floating Point Software Package
83:
84: section 8
85:
86: include fpsp.h
87:
88: X equ FP_SCR5
89: XDCARE equ X+2
90: XFRAC equ X+4
91:
92: SGN equ L_SCR3
93:
94: V equ FP_SCR6
95:
96: BOUNDS1 DC.L $3FD78000,$3FFFDDCE ... 2^(-40), (5/2)LOG2
97:
98: xref t_frcinx
99: xref t_extdnrm
100: xref setox
101: xref setoxm1
102:
103: xdef stanhd
104: stanhd:
105: *--TANH(X) = X FOR DENORMALIZED X
106:
107: bra t_extdnrm
108:
109: xdef stanh
110: stanh:
111: FMOVE.X (a0),FP0 ...LOAD INPUT
112:
113: FMOVE.X FP0,X(a6)
114: move.l (a0),d0
115: move.w 4(a0),d0
116: MOVE.L D0,X(a6)
117: AND.L #$7FFFFFFF,D0
118: CMP2.L BOUNDS1(pc),D0 ...2**(-40) < |X| < (5/2)LOG2 ?
119: BCS.B TANHBORS
120:
121: *--THIS IS THE USUAL CASE
122: *--Y = 2|X|, Z = EXPM1(Y), TANH(X) = SIGN(X) * Z / (Z+2).
123:
124: MOVE.L X(a6),D0
125: MOVE.L D0,SGN(a6)
126: AND.L #$7FFF0000,D0
127: ADD.L #$00010000,D0 ...EXPONENT OF 2|X|
128: MOVE.L D0,X(a6)
129: AND.L #$80000000,SGN(a6)
130: FMOVE.X X(a6),FP0 ...FP0 IS Y = 2|X|
131:
132: move.l d1,-(a7)
133: clr.l d1
134: fmovem.x fp0,(a0)
135: bsr setoxm1 ...FP0 IS Z = EXPM1(Y)
136: move.l (a7)+,d1
137:
138: FMOVE.X FP0,FP1
139: FADD.S #:40000000,FP1 ...Z+2
140: MOVE.L SGN(a6),D0
141: FMOVE.X FP1,V(a6)
142: EOR.L D0,V(a6)
143:
144: FMOVE.L d1,FPCR ;restore users exceptions
145: FDIV.X V(a6),FP0
146: bra t_frcinx
147:
148: TANHBORS:
149: CMP.L #$3FFF8000,D0
150: BLT.W TANHSM
151:
152: CMP.L #$40048AA1,D0
153: BGT.W TANHHUGE
154:
155: *-- (5/2) LOG2 < |X| < 50 LOG2,
156: *--TANH(X) = 1 - (2/[EXP(2X)+1]). LET Y = 2|X|, SGN = SIGN(X),
157: *--TANH(X) = SGN - SGN*2/[EXP(Y)+1].
158:
159: MOVE.L X(a6),D0
160: MOVE.L D0,SGN(a6)
161: AND.L #$7FFF0000,D0
162: ADD.L #$00010000,D0 ...EXPO OF 2|X|
163: MOVE.L D0,X(a6) ...Y = 2|X|
164: AND.L #$80000000,SGN(a6)
165: MOVE.L SGN(a6),D0
166: FMOVE.X X(a6),FP0 ...Y = 2|X|
167:
168: move.l d1,-(a7)
169: clr.l d1
170: fmovem.x fp0,(a0)
171: bsr setox ...FP0 IS EXP(Y)
172: move.l (a7)+,d1
173: move.l SGN(a6),d0
174: FADD.S #:3F800000,FP0 ...EXP(Y)+1
175:
176: EOR.L #$C0000000,D0 ...-SIGN(X)*2
177: FMOVE.S d0,FP1 ...-SIGN(X)*2 IN SGL FMT
178: FDIV.X FP0,FP1 ...-SIGN(X)2 / [EXP(Y)+1 ]
179:
180: MOVE.L SGN(a6),D0
181: OR.L #$3F800000,D0 ...SGN
182: FMOVE.S d0,FP0 ...SGN IN SGL FMT
183:
184: FMOVE.L d1,FPCR ;restore users exceptions
185: FADD.X fp1,FP0
186:
187: bra t_frcinx
188:
189: TANHSM:
190: CLR.W XDCARE(a6)
191:
192: FMOVE.L d1,FPCR ;restore users exceptions
193: FMOVE.X X(a6),FP0 ;last inst - possible exception set
194:
195: bra t_frcinx
196:
197: TANHHUGE:
198: *---RETURN SGN(X) - SGN(X)EPS
199: MOVE.L X(a6),D0
200: AND.L #$80000000,D0
201: OR.L #$3F800000,D0
202: FMOVE.S d0,FP0
203: AND.L #$80000000,D0
204: EOR.L #$80800000,D0 ...-SIGN(X)*EPS
205:
206: FMOVE.L d1,FPCR ;restore users exceptions
207: FADD.S d0,FP0
208:
209: bra t_frcinx
210:
211: end
CVSweb