1 
/* * Mode: C; tabwidth: 8; indenttabsmode: nil; cbasicoffset: 4 * 
2 
* 
3 
* ***** BEGIN LICENSE BLOCK ***** 
4 
* Version: MPL 1.1/GPL 2.0/LGPL 2.1 
5 
* 
6 
* The contents of this file are subject to the Mozilla Public License Version 
7 
* 1.1 (the "License"); you may not use this file except in compliance with 
8 
* the License. You may obtain a copy of the License at 
9 
* http://www.mozilla.org/MPL/ 
10 
* 
11 
* Software distributed under the License is distributed on an "AS IS" basis, 
12 
* WITHOUT WARRANTY OF ANY KIND, either express or implied. See the License 
13 
* for the specific language governing rights and limitations under the 
14 
* License. 
15 
* 
16 
* The Original Code is Mozilla Communicator client code, released 
17 
* March 31, 1998. 
18 
* 
19 
* The Initial Developer of the Original Code is 
20 
* Sun Microsystems, Inc. 
21 
* Portions created by the Initial Developer are Copyright (C) 1998 
22 
* the Initial Developer. All Rights Reserved. 
23 
* 
24 
* Contributor(s): 
25 
* 
26 
* Alternatively, the contents of this file may be used under the terms of 
27 
* either of the GNU General Public License Version 2 or later (the "GPL"), 
28 
* or the GNU Lesser General Public License Version 2.1 or later (the "LGPL"), 
29 
* in which case the provisions of the GPL or the LGPL are applicable instead 
30 
* of those above. If you wish to allow use of your version of this file only 
31 
* under the terms of either the GPL or the LGPL, and not to allow others to 
32 
* use your version of this file under the terms of the MPL, indicate your 
33 
* decision by deleting the provisions above and replace them with the notice 
34 
* and other provisions required by the GPL or the LGPL. If you do not delete 
35 
* the provisions above, a recipient may use your version of this file under 
36 
* the terms of any one of the MPL, the GPL or the LGPL. 
37 
* 
38 
* ***** END LICENSE BLOCK ***** */ 
39 

40 
/* @(#)e_cosh.c 1.3 95/01/18 */ 
41 
/* 
42 
* ==================================================== 
43 
* Copyright (C) 1993 by Sun Microsystems, Inc. All rights reserved. 
44 
* 
45 
* Developed at SunSoft, a Sun Microsystems, Inc. business. 
46 
* Permission to use, copy, modify, and distribute this 
47 
* software is freely granted, provided that this notice 
48 
* is preserved. 
49 
* ==================================================== 
50 
*/ 
51 

52 
/* __ieee754_cosh(x) 
53 
* Method : 
54 
* mathematically cosh(x) if defined to be (exp(x)+exp(x))/2 
55 
* 1. Replace x by x (cosh(x) = cosh(x)). 
56 
* 2. 
57 
* [ exp(x)  1 ]^2 
58 
* 0 <= x <= ln2/2 : cosh(x) := 1 +  
59 
* 2*exp(x) 
60 
* 
61 
* exp(x) + 1/exp(x) 
62 
* ln2/2 <= x <= 22 : cosh(x) :=  
63 
* 2 
64 
* 22 <= x <= lnovft : cosh(x) := exp(x)/2 
65 
* lnovft <= x <= ln2ovft: cosh(x) := exp(x/2)/2 * exp(x/2) 
66 
* ln2ovft < x : cosh(x) := huge*huge (overflow) 
67 
* 
68 
* Special cases: 
69 
* cosh(x) is x if x is +INF, INF, or NaN. 
70 
* only cosh(0)=1 is exact for finite x. 
71 
*/ 
72 

73 
#include "fdlibm.h" 
74 

75 
#ifdef _WIN32 
76 
#define huge myhuge 
77 
#endif 
78 

79 
#ifdef __STDC__ 
80 
static const double one = 1.0, half=0.5, really_big = 1.0e300; 
81 
#else 
82 
static double one = 1.0, half=0.5, really_big = 1.0e300; 
83 
#endif 
84 

85 
#ifdef __STDC__ 
86 
double __ieee754_cosh(double x) 
87 
#else 
88 
double __ieee754_cosh(x) 
89 
double x; 
90 
#endif 
91 
{ 
92 
fd_twoints u; 
93 
double t,w; 
94 
int ix; 
95 
unsigned lx; 
96 

97 
/* High word of x. */ 
98 
u.d = x; 
99 
ix = __HI(u); 
100 
ix &= 0x7fffffff; 
101 

102 
/* x is INF or NaN */ 
103 
if(ix>=0x7ff00000) return x*x; 
104 

105 
/* x in [0,0.5*ln2], return 1+expm1(x)^2/(2*exp(x)) */ 
106 
if(ix<0x3fd62e43) { 
107 
t = fd_expm1(fd_fabs(x)); 
108 
w = one+t; 
109 
if (ix<0x3c800000) return w; /* cosh(tiny) = 1 */ 
110 
return one+(t*t)/(w+w); 
111 
} 
112 

113 
/* x in [0.5*ln2,22], return (exp(x)+1/exp(x)/2; */ 
114 
if (ix < 0x40360000) { 
115 
t = __ieee754_exp(fd_fabs(x)); 
116 
return half*t+half/t; 
117 
} 
118 

119 
/* x in [22, log(maxdouble)] return half*exp(x) */ 
120 
if (ix < 0x40862E42) return half*__ieee754_exp(fd_fabs(x)); 
121 

122 
/* x in [log(maxdouble), overflowthresold] */ 
123 
lx = *( (((*(unsigned*)&one)>>29)) + (unsigned*)&x); 
124 
if (ix<0x408633CE  
125 
(ix==0x408633ce)&&(lx<=(unsigned)0x8fb9f87d)) { 
126 
w = __ieee754_exp(half*fd_fabs(x)); 
127 
t = half*w; 
128 
return t*w; 
129 
} 
130 

131 
/* x > overflowthresold, cosh(x) overflow */ 
132 
return really_big*really_big; 
133 
} 