| Linus Torvalds | 1da177e | 2005-04-16 15:20:36 -0700 | [diff] [blame] | 1 | | | 
|  | 2 | |	ssinh.sa 3.1 12/10/90 | 
|  | 3 | | | 
|  | 4 | |       The entry point sSinh computes the hyperbolic sine of | 
|  | 5 | |       an input argument; sSinhd 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 sinh(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 sSINH takes approximately 280 cycles. | 
|  | 19 | | | 
|  | 20 | |       Algorithm: | 
|  | 21 | | | 
|  | 22 | |       SINH | 
|  | 23 | |       1. If |X| > 16380 log2, go to 3. | 
|  | 24 | | | 
|  | 25 | |       2. (|X| <= 16380 log2) Sinh(X) is obtained by the formulae | 
|  | 26 | |               y = |X|, sgn = sign(X), and z = expm1(Y), | 
|  | 27 | |               sinh(X) = sgn*(1/2)*( z + z/(1+z) ). | 
|  | 28 | |          Exit. | 
|  | 29 | | | 
|  | 30 | |       3. If |X| > 16480 log2, go to 5. | 
|  | 31 | | | 
|  | 32 | |       4. (16380 log2 < |X| <= 16480 log2) | 
|  | 33 | |               sinh(X) = sign(X) * exp(|X|)/2. | 
|  | 34 | |          However, invoking exp(|X|) may cause premature overflow. | 
|  | 35 | |          Thus, we calculate sinh(X) as follows: | 
|  | 36 | |             Y       := |X| | 
|  | 37 | |             sgn     := sign(X) | 
|  | 38 | |             sgnFact := sgn * 2**(16380) | 
|  | 39 | |             Y'      := Y - 16381 log2 | 
|  | 40 | |             sinh(X) := sgnFact * exp(Y'). | 
|  | 41 | |          Exit. | 
|  | 42 | | | 
|  | 43 | |       5. (|X| > 16480 log2) sinh(X) must overflow. Return | 
|  | 44 | |          sign(X)*Huge*Huge to generate overflow and an infinity with | 
|  | 45 | |          the appropriate sign. Huge is the largest finite number in | 
|  | 46 | |          extended format. Exit. | 
|  | 47 | | | 
|  | 48 |  | 
|  | 49 | |		Copyright (C) Motorola, Inc. 1990 | 
|  | 50 | |			All Rights Reserved | 
|  | 51 | | | 
| Matt Waddel | e00d82d | 2006-02-11 17:55:48 -0800 | [diff] [blame] | 52 | |       For details on the license for this file, please see the | 
|  | 53 | |       file, README, in this same directory. | 
| Linus Torvalds | 1da177e | 2005-04-16 15:20:36 -0700 | [diff] [blame] | 54 |  | 
|  | 55 | |SSINH	idnt	2,1 | Motorola 040 Floating Point Software Package | 
|  | 56 |  | 
|  | 57 | |section	8 | 
|  | 58 |  | 
|  | 59 | T1:	.long 0x40C62D38,0xD3D64634 | ... 16381 LOG2 LEAD | 
|  | 60 | T2:	.long 0x3D6F90AE,0xB1E75CC7 | ... 16381 LOG2 TRAIL | 
|  | 61 |  | 
|  | 62 | |xref	t_frcinx | 
|  | 63 | |xref	t_ovfl | 
|  | 64 | |xref	t_extdnrm | 
|  | 65 | |xref	setox | 
|  | 66 | |xref	setoxm1 | 
|  | 67 |  | 
|  | 68 | .global	ssinhd | 
|  | 69 | ssinhd: | 
|  | 70 | |--SINH(X) = X FOR DENORMALIZED X | 
|  | 71 |  | 
|  | 72 | bra	t_extdnrm | 
|  | 73 |  | 
|  | 74 | .global	ssinh | 
|  | 75 | ssinh: | 
|  | 76 | fmovex	(%a0),%fp0	| ...LOAD INPUT | 
|  | 77 |  | 
|  | 78 | movel	(%a0),%d0 | 
|  | 79 | movew	4(%a0),%d0 | 
|  | 80 | movel	%d0,%a1		| save a copy of original (compacted) operand | 
|  | 81 | andl	#0x7FFFFFFF,%d0 | 
|  | 82 | cmpl	#0x400CB167,%d0 | 
|  | 83 | bgts	SINHBIG | 
|  | 84 |  | 
|  | 85 | |--THIS IS THE USUAL CASE, |X| < 16380 LOG2 | 
|  | 86 | |--Y = |X|, Z = EXPM1(Y), SINH(X) = SIGN(X)*(1/2)*( Z + Z/(1+Z) ) | 
|  | 87 |  | 
|  | 88 | fabsx	%fp0		| ...Y = |X| | 
|  | 89 |  | 
|  | 90 | moveml	%a1/%d1,-(%sp) | 
|  | 91 | fmovemx %fp0-%fp0,(%a0) | 
|  | 92 | clrl	%d1 | 
|  | 93 | bsr	setoxm1		| ...FP0 IS Z = EXPM1(Y) | 
|  | 94 | fmovel	#0,%fpcr | 
|  | 95 | moveml	(%sp)+,%a1/%d1 | 
|  | 96 |  | 
|  | 97 | fmovex	%fp0,%fp1 | 
|  | 98 | fadds	#0x3F800000,%fp1	| ...1+Z | 
|  | 99 | fmovex	%fp0,-(%sp) | 
|  | 100 | fdivx	%fp1,%fp0		| ...Z/(1+Z) | 
|  | 101 | movel	%a1,%d0 | 
|  | 102 | andl	#0x80000000,%d0 | 
|  | 103 | orl	#0x3F000000,%d0 | 
|  | 104 | faddx	(%sp)+,%fp0 | 
|  | 105 | movel	%d0,-(%sp) | 
|  | 106 |  | 
|  | 107 | fmovel	%d1,%fpcr | 
|  | 108 | fmuls	(%sp)+,%fp0	|last fp inst - possible exceptions set | 
|  | 109 |  | 
|  | 110 | bra	t_frcinx | 
|  | 111 |  | 
|  | 112 | SINHBIG: | 
|  | 113 | cmpl	#0x400CB2B3,%d0 | 
|  | 114 | bgt	t_ovfl | 
|  | 115 | fabsx	%fp0 | 
|  | 116 | fsubd	T1(%pc),%fp0	| ...(|X|-16381LOG2_LEAD) | 
|  | 117 | movel	#0,-(%sp) | 
|  | 118 | movel	#0x80000000,-(%sp) | 
|  | 119 | movel	%a1,%d0 | 
|  | 120 | andl	#0x80000000,%d0 | 
|  | 121 | orl	#0x7FFB0000,%d0 | 
|  | 122 | movel	%d0,-(%sp)	| ...EXTENDED FMT | 
|  | 123 | fsubd	T2(%pc),%fp0	| ...|X| - 16381 LOG2, ACCURATE | 
|  | 124 |  | 
|  | 125 | movel	%d1,-(%sp) | 
|  | 126 | clrl	%d1 | 
|  | 127 | fmovemx %fp0-%fp0,(%a0) | 
|  | 128 | bsr	setox | 
|  | 129 | fmovel	(%sp)+,%fpcr | 
|  | 130 |  | 
|  | 131 | fmulx	(%sp)+,%fp0	|possible exception | 
|  | 132 | bra	t_frcinx | 
|  | 133 |  | 
|  | 134 | |end |