RTAI API: base/math/e_cosh.c Source File

00001 /* @(#)e_cosh.c 5.1 93/09/24 */ 00002 /* 00003 * ==================================================== 00004 * Copyright (C) 1993 by Sun Microsystems, Inc. All rights reserved. 00005 * 00006 * Developed at SunPro, a Sun Microsystems, Inc. business. 00007 * Permission to use, copy, modify, and distribute this 00008 * software is freely granted, provided that this notice 00009 * is preserved. 00010 * ==================================================== 00011 */ 00012 00013 #if defined(LIBM_SCCS) && !defined(lint) 00014 static char rcsid[] = "$NetBSD: e_cosh.c,v 1.7 1995/05/10 20:44:58 jtc Exp $"; 00015 #endif 00016 00017 /* __ieee754_cosh(x) 00018 * Method : 00019 * mathematically cosh(x) if defined to be (exp(x)+exp(-x))/2 00020 * 1. Replace x by |x| (cosh(x) = cosh(-x)). 00021 * 2. 00022 * [ exp(x) - 1 ]^2 00023 * 0 <= x <= ln2/2 : cosh(x) := 1 + ------------------- 00024 * 2*exp(x) 00025 * 00026 * exp(x) + 1/exp(x) 00027 * ln2/2 <= x <= 22 : cosh(x) := ------------------- 00028 * 2 00029 * 22 <= x <= lnovft : cosh(x) := exp(x)/2 00030 * lnovft <= x <= ln2ovft: cosh(x) := exp(x/2)/2 * exp(x/2) 00031 * ln2ovft < x : cosh(x) := huge*huge (overflow) 00032 * 00033 * Special cases: 00034 * cosh(x) is |x| if x is +INF, -INF, or NaN. 00035 * only cosh(0)=1 is exact for finite x. 00036 */ 00037 00038 #include "math.h" 00039 #include "mathP.h" 00040 00041 #ifdef __STDC__ 00042 static const double one = 1.0, half=0.5, huge = 1.0e300; 00043 #else 00044 static double one = 1.0, half=0.5, huge = 1.0e300; 00045 #endif 00046 00047 #ifdef __STDC__ 00048 double __ieee754_cosh(double x) 00049 #else 00050 double __ieee754_cosh(x) 00051 double x; 00052 #endif 00053 { 00054 double t,w; 00055 int32_t ix; 00056 u_int32_t lx; 00057 00058 /* High word of |x|. */ 00059 GET_HIGH_WORD(ix,x); 00060 ix &= 0x7fffffff; 00061 00062 /* x is INF or NaN */ 00063 if(ix>=0x7ff00000) return x*x; 00064 00065 /* |x| in [0,0.5*ln2], return 1+expm1(|x|)^2/(2*exp(|x|)) */ 00066 if(ix<0x3fd62e43) { 00067 t = expm1(fabs(x)); 00068 w = one+t; 00069 if (ix<0x3c800000) return w; /* cosh(tiny) = 1 */ 00070 return one+(t*t)/(w+w); 00071 } 00072 00073 /* |x| in [0.5*ln2,22], return (exp(|x|)+1/exp(|x|)/2; */ 00074 if (ix < 0x40360000) { 00075 t = __ieee754_exp(fabs(x)); 00076 return half*t+half/t; 00077 } 00078 00079 /* |x| in [22, log(maxdouble)] return half*exp(|x|) */ 00080 if (ix < 0x40862E42) return half*__ieee754_exp(fabs(x)); 00081 00082 /* |x| in [log(maxdouble), overflowthresold] */ 00083 GET_LOW_WORD(lx,x); 00084 if (ix<0x408633CE || 00085 ((ix==0x408633ce)&&(lx<=(u_int32_t)0x8fb9f87d))) { 00086 w = __ieee754_exp(half*fabs(x)); 00087 t = half*w; 00088 return t*w; 00089 } 00090 00091 /* |x| > overflowthresold, cosh(x) overflow */ 00092 return huge*huge; 00093 }