home *** CD-ROM | disk | FTP | other *** search
/ Game Zone - 1,000+ Games / GAMEZONE.BIN / Programs / PALM / Oh-One / src / e_acosh.c < prev    next >
C/C++ Source or Header  |  1997-08-16  |  2KB  |  72 lines
  1. // 15 August 1997, Rick Huebner:  Small changes made to adapt for MathLib
  2.  
  3. /* @(#)e_acosh.c 5.1 93/09/24 */
  4. /*
  5.  * ====================================================
  6.  * Copyright (C) 1993 by Sun Microsystems, Inc. All rights reserved.
  7.  *
  8.  * Developed at SunPro, a Sun Microsystems, Inc. business.
  9.  * Permission to use, copy, modify, and distribute this
  10.  * software is freely granted, provided that this notice 
  11.  * is preserved.
  12.  * ====================================================
  13.  */
  14.  
  15. #if defined(LIBM_SCCS) && !defined(lint)
  16. static char rcsid[] = "$NetBSD: e_acosh.c,v 1.9 1995/05/12 04:57:18 jtc Exp $";
  17. #endif
  18.  
  19. /* __ieee754_acosh(x)
  20.  * Method :
  21.  *    Based on 
  22.  *        acosh(x) = log [ x + sqrt(x*x-1) ]
  23.  *    we have
  24.  *        acosh(x) := log(x)+ln2,    if x is large; else
  25.  *        acosh(x) := log(2x-1/(sqrt(x*x-1)+x)) if x>2; else
  26.  *        acosh(x) := log1p(t+sqrt(2.0*t+t*t)); where t=x-1.
  27.  *
  28.  * Special cases:
  29.  *    acosh(x) is NaN with signal if x<1.
  30.  *    acosh(NaN) is NaN without signal.
  31.  */
  32.  
  33. #include "math.h"
  34. #include "math_private.h"
  35.  
  36. #ifdef __STDC__
  37. static const double 
  38. #else
  39. static double 
  40. #endif
  41. one    = 1.0,
  42. ln2    = 6.93147180559945286227e-01;  /* 0x3FE62E42, 0xFEFA39EF */
  43.  
  44. #ifdef __STDC__
  45.     double __ieee754_acosh(double x)
  46. #else
  47.     double __ieee754_acosh(x)
  48.     double x;
  49. #endif
  50. {    
  51.     double t;
  52.     int32_t hx;
  53.     u_int32_t lx;
  54.     EXTRACT_WORDS(hx,lx,x);
  55.     if(hx<0x3ff00000) {        /* x < 1 */
  56.         return (x-x)/(x-x);
  57.     } else if(hx >=0x41b00000) {    /* x > 2**28 */
  58.         if(hx >=0x7ff00000) {    /* x is inf of NaN */
  59.             return x+x;
  60.         } else 
  61.         return __ieee754_log(x)+ln2;    /* acosh(huge)=log(2x) */
  62.     } else if(((hx-0x3ff00000)|lx)==0) {
  63.         return 0.0;            /* acosh(1) = 0 */
  64.     } else if (hx > 0x40000000) {    /* 2**28 > x > 2 */
  65.         t=x*x;
  66.         return __ieee754_log(2.0*x-one/(x+__ieee754_sqrt(t-one)));
  67.     } else {            /* 1<x<2 */
  68.         t = x-one;
  69.         return __log1p(t+__ieee754_sqrt(2.0*t+t*t));
  70.     }
  71. }
  72.