1* $NetBSD: ssinh.sa,v 1.3 1994/10/26 07:50:05 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* ssinh.sa 3.1 12/10/90 35* 36* The entry point sSinh computes the hyperbolic sine of 37* an input argument; sSinhd 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 sinh(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 sSINH takes approximately 280 cycles. 51* 52* Algorithm: 53* 54* SINH 55* 1. If |X| > 16380 log2, go to 3. 56* 57* 2. (|X| <= 16380 log2) Sinh(X) is obtained by the formulae 58* y = |X|, sgn = sign(X), and z = expm1(Y), 59* sinh(X) = sgn*(1/2)*( z + z/(1+z) ). 60* Exit. 61* 62* 3. If |X| > 16480 log2, go to 5. 63* 64* 4. (16380 log2 < |X| <= 16480 log2) 65* sinh(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* sgn := sign(X) 70* sgnFact := sgn * 2**(16380) 71* Y' := Y - 16381 log2 72* sinh(X) := sgnFact * exp(Y'). 73* Exit. 74* 75* 5. (|X| > 16480 log2) sinh(X) must overflow. Return 76* sign(X)*Huge*Huge to generate overflow and an infinity with 77* the appropriate sign. Huge is the largest finite number in 78* extended format. Exit. 79* 80 81SSINH IDNT 2,1 Motorola 040 Floating Point Software Package 82 83 section 8 84 85T1 DC.L $40C62D38,$D3D64634 ... 16381 LOG2 LEAD 86T2 DC.L $3D6F90AE,$B1E75CC7 ... 16381 LOG2 TRAIL 87 88 xref t_frcinx 89 xref t_ovfl 90 xref t_extdnrm 91 xref setox 92 xref setoxm1 93 94 xdef ssinhd 95ssinhd: 96*--SINH(X) = X FOR DENORMALIZED X 97 98 bra t_extdnrm 99 100 xdef ssinh 101ssinh: 102 FMOVE.x (a0),FP0 ...LOAD INPUT 103 104 move.l (a0),d0 105 move.w 4(a0),d0 106 move.l d0,a1 save a copy of original (compacted) operand 107 AND.L #$7FFFFFFF,D0 108 CMP.L #$400CB167,D0 109 BGT.B SINHBIG 110 111*--THIS IS THE USUAL CASE, |X| < 16380 LOG2 112*--Y = |X|, Z = EXPM1(Y), SINH(X) = SIGN(X)*(1/2)*( Z + Z/(1+Z) ) 113 114 FABS.X FP0 ...Y = |X| 115 116 movem.l a1/d1,-(sp) 117 fmovem.x fp0,(a0) 118 clr.l d1 119 bsr setoxm1 ...FP0 IS Z = EXPM1(Y) 120 fmove.l #0,fpcr 121 movem.l (sp)+,a1/d1 122 123 FMOVE.X FP0,FP1 124 FADD.S #:3F800000,FP1 ...1+Z 125 FMOVE.X FP0,-(sp) 126 FDIV.X FP1,FP0 ...Z/(1+Z) 127 MOVE.L a1,d0 128 AND.L #$80000000,D0 129 OR.L #$3F000000,D0 130 FADD.X (sp)+,FP0 131 MOVE.L D0,-(sp) 132 133 fmove.l d1,fpcr 134 fmul.s (sp)+,fp0 ;last fp inst - possible exceptions set 135 136 bra t_frcinx 137 138SINHBIG: 139 cmp.l #$400CB2B3,D0 140 bgt t_ovfl 141 FABS.X FP0 142 FSUB.D T1(pc),FP0 ...(|X|-16381LOG2_LEAD) 143 clr.l -(sp) 144 move.l #$80000000,-(sp) 145 move.l a1,d0 146 AND.L #$80000000,D0 147 OR.L #$7FFB0000,D0 148 MOVE.L D0,-(sp) ...EXTENDED FMT 149 FSUB.D T2(pc),FP0 ...|X| - 16381 LOG2, ACCURATE 150 151 move.l d1,-(sp) 152 clr.l d1 153 fmovem.x fp0,(a0) 154 bsr setox 155 fmove.l (sp)+,fpcr 156 157 fmul.x (sp)+,fp0 ;possible exception 158 bra t_frcinx 159 160 end 161