xref: /netbsd/sys/arch/m68k/fpsp/scosh.sa (revision bf9ec67e)
1*	$NetBSD: scosh.sa,v 1.2 1994/10/26 07:49:39 cgd Exp $
2
3*	MOTOROLA MICROPROCESSOR & MEMORY TECHNOLOGY GROUP
4*	M68000 Hi-Performance Microprocessor Division
5*	M68040 Software Package
6*
7*	M68040 Software Package Copyright (c) 1993, 1994 Motorola Inc.
8*	All rights reserved.
9*
10*	THE SOFTWARE is provided on an "AS IS" basis and without warranty.
11*	To the maximum extent permitted by applicable law,
12*	MOTOROLA DISCLAIMS ALL WARRANTIES WHETHER EXPRESS OR IMPLIED,
13*	INCLUDING IMPLIED WARRANTIES OF MERCHANTABILITY OR FITNESS FOR A
14*	PARTICULAR PURPOSE and any warranty against infringement with
15*	regard to the SOFTWARE (INCLUDING ANY MODIFIED VERSIONS THEREOF)
16*	and any accompanying written materials.
17*
18*	To the maximum extent permitted by applicable law,
19*	IN NO EVENT SHALL MOTOROLA BE LIABLE FOR ANY DAMAGES WHATSOEVER
20*	(INCLUDING WITHOUT LIMITATION, DAMAGES FOR LOSS OF BUSINESS
21*	PROFITS, BUSINESS INTERRUPTION, LOSS OF BUSINESS INFORMATION, OR
22*	OTHER PECUNIARY LOSS) ARISING OF THE USE OR INABILITY TO USE THE
23*	SOFTWARE.  Motorola assumes no responsibility for the maintenance
24*	and support of the SOFTWARE.
25*
26*	You are hereby granted a copyright license to use, modify, and
27*	distribute the SOFTWARE so long as this entire notice is retained
28*	without alteration in any modified and/or redistributed versions,
29*	and that such modified versions are clearly identified as such.
30*	No licenses are granted by implication, estoppel or otherwise
31*	under any patents or trademarks of Motorola, Inc.
32
33*
34*	scosh.sa 3.1 12/10/90
35*
36*	The entry point sCosh computes the hyperbolic cosine of
37*	an input argument; sCoshd does the same except for denormalized
38*	input.
39*
40*	Input: Double-extended number X in location pointed to
41*		by address register a0.
42*
43*	Output: The value cosh(X) returned in floating-point register Fp0.
44*
45*	Accuracy and Monotonicity: The returned result is within 3 ulps in
46*		64 significant bit, i.e. within 0.5001 ulp to 53 bits if the
47*		result is subsequently rounded to double precision. The
48*		result is provably monotonic in double precision.
49*
50*	Speed: The program sCOSH takes approximately 250 cycles.
51*
52*	Algorithm:
53*
54*	COSH
55*	1. If |X| > 16380 log2, go to 3.
56*
57*	2. (|X| <= 16380 log2) Cosh(X) is obtained by the formulae
58*		y = |X|, z = exp(Y), and
59*		cosh(X) = (1/2)*( z + 1/z ).
60*		Exit.
61*
62*	3. (|X| > 16380 log2). If |X| > 16480 log2, go to 5.
63*
64*	4. (16380 log2 < |X| <= 16480 log2)
65*		cosh(X) = sign(X) * exp(|X|)/2.
66*		However, invoking exp(|X|) may cause premature overflow.
67*		Thus, we calculate sinh(X) as follows:
68*		Y	:= |X|
69*		Fact	:=	2**(16380)
70*		Y'	:= Y - 16381 log2
71*		cosh(X) := Fact * exp(Y').
72*		Exit.
73*
74*	5. (|X| > 16480 log2) sinh(X) must overflow. Return
75*		Huge*Huge to generate overflow and an infinity with
76*		the appropriate sign. Huge is the largest finite number in
77*		extended format. Exit.
78*
79
80SCOSH	IDNT	2,1 Motorola 040 Floating Point Software Package
81
82	section	8
83
84	xref	t_ovfl
85	xref	t_frcinx
86	xref	setox
87
88T1	DC.L $40C62D38,$D3D64634 ... 16381 LOG2 LEAD
89T2	DC.L $3D6F90AE,$B1E75CC7 ... 16381 LOG2 TRAIL
90
91TWO16380 DC.L $7FFB0000,$80000000,$00000000,$00000000
92
93	xdef	scoshd
94scoshd:
95*--COSH(X) = 1 FOR DENORMALIZED X
96
97	FMOVE.S		#:3F800000,FP0
98
99	FMOVE.L		d1,FPCR
100	FADD.S		#:00800000,FP0
101	bra		t_frcinx
102
103	xdef	scosh
104scosh:
105	FMOVE.X		(a0),FP0	...LOAD INPUT
106
107	move.l		(a0),d0
108	move.w		4(a0),d0
109	ANDI.L		#$7FFFFFFF,d0
110	CMPI.L		#$400CB167,d0
111	BGT.B		COSHBIG
112
113*--THIS IS THE USUAL CASE, |X| < 16380 LOG2
114*--COSH(X) = (1/2) * ( EXP(X) + 1/EXP(X) )
115
116	FABS.X		FP0		...|X|
117
118	move.l		d1,-(sp)
119	clr.l		d1
120	fmovem.x	fp0,(a0)	;pass parameter to setox
121	bsr		setox		...FP0 IS EXP(|X|)
122	FMUL.S		#:3F000000,FP0	...(1/2)EXP(|X|)
123	move.l		(sp)+,d1
124
125	FMOVE.S		#:3E800000,FP1	...(1/4)
126	FDIV.X		FP0,FP1	 	...1/(2 EXP(|X|))
127
128	FMOVE.L		d1,FPCR
129	FADD.X		fp1,FP0
130
131	bra		t_frcinx
132
133COSHBIG:
134	CMPI.L		#$400CB2B3,d0
135	BGT.B		COSHHUGE
136
137	FABS.X		FP0
138	FSUB.D		T1(pc),FP0		...(|X|-16381LOG2_LEAD)
139	FSUB.D		T2(pc),FP0		...|X| - 16381 LOG2, ACCURATE
140
141	move.l		d1,-(sp)
142	clr.l		d1
143	fmovem.x	fp0,(a0)
144	bsr		setox
145	fmove.l		(sp)+,fpcr
146
147	FMUL.X		TWO16380(pc),FP0
148	bra		t_frcinx
149
150COSHHUGE:
151	fmove.l		#0,fpsr		;clr N bit if set by source
152	bclr.b		#7,(a0)		;always return positive value
153	fmovem.x	(a0),fp0
154	bra		t_ovfl
155
156	end
157