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| 52| For details on the license for this file, please see the 53| file, README, in this same directory. 54 55|SSINH idnt 2,1 | Motorola 040 Floating Point Software Package 56 57 |section 8 58 59T1: .long 0x40C62D38,0xD3D64634 | ... 16381 LOG2 LEAD 60T2: .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 69ssinhd: 70|--SINH(X) = X FOR DENORMALIZED X 71 72 bra t_extdnrm 73 74 .global ssinh 75ssinh: 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 112SINHBIG: 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 135