/usr/share/go-1.6/src/math/pow.go is in golang-1.6-src 1.6.1-0ubuntu1.
This file is owned by root:root, with mode 0o644.
The actual contents of the file can be viewed below.
1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60 61 62 63 64 65 66 67 68 69 70 71 72 73 74 75 76 77 78 79 80 81 82 83 84 85 86 87 88 89 90 91 92 93 94 95 96 97 98 99 100 101 102 103 104 105 106 107 108 109 110 111 112 113 114 115 116 117 118 119 120 121 122 123 124 125 126 127 128 129 130 131 132 133 134 135 136 137 | // Copyright 2009 The Go Authors. All rights reserved.
// Use of this source code is governed by a BSD-style
// license that can be found in the LICENSE file.
package math
func isOddInt(x float64) bool {
xi, xf := Modf(x)
return xf == 0 && int64(xi)&1 == 1
}
// Special cases taken from FreeBSD's /usr/src/lib/msun/src/e_pow.c
// updated by IEEE Std. 754-2008 "Section 9.2.1 Special values".
// Pow returns x**y, the base-x exponential of y.
//
// Special cases are (in order):
// Pow(x, ±0) = 1 for any x
// Pow(1, y) = 1 for any y
// Pow(x, 1) = x for any x
// Pow(NaN, y) = NaN
// Pow(x, NaN) = NaN
// Pow(±0, y) = ±Inf for y an odd integer < 0
// Pow(±0, -Inf) = +Inf
// Pow(±0, +Inf) = +0
// Pow(±0, y) = +Inf for finite y < 0 and not an odd integer
// Pow(±0, y) = ±0 for y an odd integer > 0
// Pow(±0, y) = +0 for finite y > 0 and not an odd integer
// Pow(-1, ±Inf) = 1
// Pow(x, +Inf) = +Inf for |x| > 1
// Pow(x, -Inf) = +0 for |x| > 1
// Pow(x, +Inf) = +0 for |x| < 1
// Pow(x, -Inf) = +Inf for |x| < 1
// Pow(+Inf, y) = +Inf for y > 0
// Pow(+Inf, y) = +0 for y < 0
// Pow(-Inf, y) = Pow(-0, -y)
// Pow(x, y) = NaN for finite x < 0 and finite non-integer y
func Pow(x, y float64) float64 {
switch {
case y == 0 || x == 1:
return 1
case y == 1:
return x
case y == 0.5:
return Sqrt(x)
case y == -0.5:
return 1 / Sqrt(x)
case IsNaN(x) || IsNaN(y):
return NaN()
case x == 0:
switch {
case y < 0:
if isOddInt(y) {
return Copysign(Inf(1), x)
}
return Inf(1)
case y > 0:
if isOddInt(y) {
return x
}
return 0
}
case IsInf(y, 0):
switch {
case x == -1:
return 1
case (Abs(x) < 1) == IsInf(y, 1):
return 0
default:
return Inf(1)
}
case IsInf(x, 0):
if IsInf(x, -1) {
return Pow(1/x, -y) // Pow(-0, -y)
}
switch {
case y < 0:
return 0
case y > 0:
return Inf(1)
}
}
absy := y
flip := false
if absy < 0 {
absy = -absy
flip = true
}
yi, yf := Modf(absy)
if yf != 0 && x < 0 {
return NaN()
}
if yi >= 1<<63 {
return Exp(y * Log(x))
}
// ans = a1 * 2**ae (= 1 for now).
a1 := 1.0
ae := 0
// ans *= x**yf
if yf != 0 {
if yf > 0.5 {
yf--
yi++
}
a1 = Exp(yf * Log(x))
}
// ans *= x**yi
// by multiplying in successive squarings
// of x according to bits of yi.
// accumulate powers of two into exp.
x1, xe := Frexp(x)
for i := int64(yi); i != 0; i >>= 1 {
if i&1 == 1 {
a1 *= x1
ae += xe
}
x1 *= x1
xe <<= 1
if x1 < .5 {
x1 += x1
xe--
}
}
// ans = a1*2**ae
// if flip { ans = 1 / ans }
// but in the opposite order
if flip {
a1 = 1 / a1
ae = -ae
}
return Ldexp(a1, ae)
}
|