e_acosh.c - vx32 - Local 9vx git repository for patches.
 (HTM) git clone git://r-36.net/vx32
 (DIR) Log
 (DIR) Files
 (DIR) Refs
       ---
       e_acosh.c (1672B)
       ---
            1 /* @(#)e_acosh.c 5.1 93/09/24 */
            2 /*
            3  * ====================================================
            4  * Copyright (C) 1993 by Sun Microsystems, Inc. All rights reserved.
            5  *
            6  * Developed at SunPro, a Sun Microsystems, Inc. business.
            7  * Permission to use, copy, modify, and distribute this
            8  * software is freely granted, provided that this notice
            9  * is preserved.
           10  * ====================================================
           11  */
           12 
           13 #ifndef lint
           14 static char rcsid[] = "$FreeBSD: src/lib/msun/src/e_acosh.c,v 1.7 2002/05/28 17:03:12 alfred Exp $";
           15 #endif
           16 
           17 /* __ieee754_acosh(x)
           18  * Method :
           19  *        Based on
           20  *                acosh(x) = log [ x + sqrt(x*x-1) ]
           21  *        we have
           22  *                acosh(x) := log(x)+ln2,        if x is large; else
           23  *                acosh(x) := log(2x-1/(sqrt(x*x-1)+x)) if x>2; else
           24  *                acosh(x) := log1p(t+sqrt(2.0*t+t*t)); where t=x-1.
           25  *
           26  * Special cases:
           27  *        acosh(x) is NaN with signal if x<1.
           28  *        acosh(NaN) is NaN without signal.
           29  */
           30 
           31 #include "math.h"
           32 #include "math_private.h"
           33 
           34 static const double
           35 one        = 1.0,
           36 ln2        = 6.93147180559945286227e-01;  /* 0x3FE62E42, 0xFEFA39EF */
           37 
           38 double
           39 __ieee754_acosh(double x)
           40 {
           41         double t;
           42         int32_t hx;
           43         u_int32_t lx;
           44         EXTRACT_WORDS(hx,lx,x);
           45         if(hx<0x3ff00000) {                /* x < 1 */
           46             return (x-x)/(x-x);
           47         } else if(hx >=0x41b00000) {        /* x > 2**28 */
           48             if(hx >=0x7ff00000) {        /* x is inf of NaN */
           49                 return x+x;
           50             } else
           51                 return __ieee754_log(x)+ln2;        /* acosh(huge)=log(2x) */
           52         } else if(((hx-0x3ff00000)|lx)==0) {
           53             return 0.0;                        /* acosh(1) = 0 */
           54         } else if (hx > 0x40000000) {        /* 2**28 > x > 2 */
           55             t=x*x;
           56             return __ieee754_log(2.0*x-one/(x+__ieee754_sqrt(t-one)));
           57         } else {                        /* 1<x<2 */
           58             t = x-one;
           59             return log1p(t+__ieee754_sqrt(2.0*t+t*t));
           60         }
           61 }