src/math/acosh.go

   1 // Copyright 2010 The Go Authors. All rights reserved.
   2 // Use of this source code is governed by a BSD-style
   3 // license that can be found in the LICENSE file.
   4
   5 package math
   6
   7 // The original C code, the long comment, and the constants
   8 // below are from FreeBSD's /usr/src/lib/msun/src/e_acosh.c
   9 // and came with this notice. The go code is a simplified
  10 // version of the original C.
  11 //
  12 // ====================================================
  13 // Copyright (C) 1993 by Sun Microsystems, Inc. All rights reserved.
  14 //
  15 // Developed at SunPro, a Sun Microsystems, Inc. business.
  16 // Permission to use, copy, modify, and distribute this
  17 // software is freely granted, provided that this notice
  18 // is preserved.
  19 // ====================================================
  20 //
  21 //
  22 // __ieee754_acosh(x)
  23 // Method :
  24 //      Based on
  25 //              acosh(x) = log [ x + sqrt(x*x-1) ]
  26 //      we have
  27 //              acosh(x) := log(x)+ln2, if x is large; else
  28 //              acosh(x) := log(2x-1/(sqrt(x*x-1)+x)) if x>2; else
  29 //              acosh(x) := log1p(t+sqrt(2.0*t+t*t)); where t=x-1.
  30 //
  31 // Special cases:
  32 //      acosh(x) is NaN with signal if x<1.
  33 //      acosh(NaN) is NaN without signal.
  34 //
  35
  36 // Acosh returns the inverse hyperbolic cosine of x.
  37 //
  38 // Special cases are:
  39 //
  40 //      Acosh(+Inf) = +Inf
  41 //      Acosh(x) = NaN if x < 1
  42 //      Acosh(NaN) = NaN
  43 func Acosh(x float64) float64 {
  44         if haveArchAcosh {
  45                 return archAcosh(x)
  46         }
  47         return acosh(x)
  48 }
  49
  50 func acosh(x float64) float64 {
  51         const Large = 1 << 28 // 2**28
  52         // first case is special case
  53         switch {
  54         case x < 1 || IsNaN(x):
  55                 return NaN()
  56         case x == 1:
  57                 return 0
  58         case x >= Large:
  59                 return Log(x) + Ln2 // x > 2**28
  60         case x > 2:
  61                 return Log(2*x - 1/(x+Sqrt(x*x-1))) // 2**28 > x > 2
  62         }
  63         t := x - 1
  64         return Log1p(t + Sqrt(2*t+t*t)) // 2 >= x > 1
  65 }