1
0
mirror of https://github.com/golang/go synced 2024-11-21 23:44:39 -07:00

math: regularize build

This will be nicer to the automatic tools.
It requires a few more assembly stubs
but fewer Go files.

There are a few instances where it looks like
there are new blobs of code, but they are just
being copied out of deleted files.

There is no new code here.

Suppose you have a portable implementation for Sin
and a 386-specific assembly one.  The old way to
do this was to write three files

sin_decl.go
   func Sin(x float64) float64  // declaration only
sin_386.s
   assembly implementation

sin_port.go
   func Sin(x float64) float64 { ... }  // pure-Go impl

and then link in either sin_decl.go+sin_386.s or
just sin_port.go.  The Makefile actually did the magic
of linking in only the _port.go files for those without
assembly and only the _decl.go files for those with
assembly, or at least some of that magic.

The biggest problem with this, beyond being hard
to explain to the build system, is that once you do
explain it to the build system, godoc knows which
of sin_port.go or sin_decl.go are involved on a given
architecture, and it (correctly) ignores the other.
That means you have to put identical doc comments
in both files.

The new approach, which is more like what we did
in the later packages math/big and sync/atomic,
is to have

sin.go
   func Sin(x float64) float64  // decl only
   func sin(x float64) float64 {...}  // pure-Go impl

sin_386.s
   // assembly for Sin (ignores sin)
sin_amd64.s
   // assembly for Sin: jmp sin
sin_arm.s
   // assembly for Sin: jmp sin

Once we abandon Makefiles we can put all the assembly
stubs in one source file, so the number of files will
actually go down.

Chris asked whether the branches cost anything.
Given that they are branching to pure-Go implementations
that are not typically known for their speed, the single
direct branch is not going to be noticeable.  That is,
it's on the slow path.

An alternative would have been to preserve the old
"only write assembly files when there's an implementation"
and still have just one copy of the declaration of Sin
(and thus one doc comment) by doing:

sin.go
   func Sin(x float64) float64 { return sin(x) }

sin_decl.go
   func sin(x float64) float64 // declaration only
sin_386.s
   // assembly for sin

sin_port.go
   func sin(x float64) float64 { portable code }

In this version everyone would link in sin.go and
then either sin_decl.go+sin_386.s or sin_port.go.

This has an extra function call on all paths, including
the "fast path" to get to assembly, and it triples the
number of Go files involved compared to what I did
in this CL.  On the other hand you don't have to
write assembly stubs.  After starting down this path
I decided that the assembly stubs were the easier
approach.

As for generating the assembly stubs on the fly, much
of the goal here is to eliminate magic from the build
process, so that zero-configuration tools like goinstall
or the new go tool can handle this package.

R=golang-dev, r, cw, iant, r
CC=golang-dev
https://golang.org/cl/5488057
This commit is contained in:
Russ Cox 2011-12-13 15:20:12 -05:00
parent 6f975fbb31
commit dd8dc6f059
84 changed files with 763 additions and 655 deletions

View File

@ -6,44 +6,31 @@ include ../../Make.inc
TARG=math
OFILES_arm=\
sqrt_arm.$O\
OFILES_amd64=\
abs_amd64.$O\
dim_amd64.$O\
exp_amd64.$O\
hypot_amd64.$O\
log_amd64.$O\
sqrt_amd64.$O\
OFILES_386=\
abs_386.$O\
asin_386.$O\
atan_386.$O\
atan2_386.$O\
exp_386.$O\
exp2_386.$O\
expm1_386.$O\
floor_386.$O\
frexp_386.$O\
hypot_386.$O\
ldexp_386.$O\
log_386.$O\
log10_386.$O\
log1p_386.$O\
mod_386.$O\
modf_386.$O\
remainder_386.$O\
sin_386.$O\
sincos_386.$O\
sqrt_386.$O\
tan_386.$O\
OFILES=\
$(OFILES_$(GOARCH))
abs_$(GOARCH).$O\
asin_$(GOARCH).$O\
atan_$(GOARCH).$O\
atan2_$(GOARCH).$O\
dim_$(GOARCH).$O\
exp_$(GOARCH).$O\
exp2_$(GOARCH).$O\
expm1_$(GOARCH).$O\
floor_$(GOARCH).$O\
frexp_$(GOARCH).$O\
hypot_$(GOARCH).$O\
ldexp_$(GOARCH).$O\
log_$(GOARCH).$O\
log10_$(GOARCH).$O\
log1p_$(GOARCH).$O\
mod_$(GOARCH).$O\
modf_$(GOARCH).$O\
remainder_$(GOARCH).$O\
sin_$(GOARCH).$O\
sincos_$(GOARCH).$O\
sqrt_$(GOARCH).$O\
tan_$(GOARCH).$O\
ALLGOFILES=\
GOFILES=\
abs.go\
acosh.go\
asin.go\
@ -58,14 +45,11 @@ ALLGOFILES=\
dim.go\
erf.go\
exp.go\
exp_port.go\
exp2.go\
expm1.go\
floor.go\
frexp.go\
gamma.go\
hypot.go\
hypot_port.go\
j0.go\
j1.go\
jn.go\
@ -86,16 +70,8 @@ ALLGOFILES=\
sincos.go\
sinh.go\
sqrt.go\
sqrt_port.go\
tan.go\
tanh.go\
unsafe.go\
NOGOFILES=\
$(subst _$(GOARCH).$O,.go,$(OFILES_$(GOARCH)))
GOFILES=\
$(filter-out $(NOGOFILES),$(ALLGOFILES))\
$(subst .go,_decl.go,$(NOGOFILES))\
include ../../Make.pkg

View File

@ -9,7 +9,9 @@ package math
// Special cases are:
// Abs(±Inf) = +Inf
// Abs(NaN) = NaN
func Abs(x float64) float64 {
func Abs(x float64) float64
func abs(x float64) float64 {
switch {
case x < 0:
return -x

View File

@ -1,7 +1,6 @@
// Copyright 2010 The Go Authors. All rights reserved.
// Copyright 2011 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 Log(x float64) float64
TEXT ·Abs(SB),7,$0
B ·abs(SB)

View File

@ -2062,6 +2062,34 @@ func TestHypot(t *testing.T) {
}
}
func TestHypotSqrtGo(t *testing.T) {
for i := 0; i < len(vf); i++ {
a := Abs(1e200 * tanh[i] * Sqrt(2))
if f := HypotSqrtGo(1e200*tanh[i], 1e200*tanh[i]); !veryclose(a, f) {
t.Errorf("HypotSqrtGo(%g, %g) = %g, want %g", 1e200*tanh[i], 1e200*tanh[i], f, a)
}
}
for i := 0; i < len(vfhypotSC); i++ {
if f := HypotSqrtGo(vfhypotSC[i][0], vfhypotSC[i][1]); !alike(hypotSC[i], f) {
t.Errorf("HypotSqrtGo(%g, %g) = %g, want %g", vfhypotSC[i][0], vfhypotSC[i][1], f, hypotSC[i])
}
}
}
func TestHypotNoSqrtGo(t *testing.T) {
for i := 0; i < len(vf); i++ {
a := Abs(1e200 * tanh[i] * Sqrt(2))
if f := HypotNoSqrtGo(1e200*tanh[i], 1e200*tanh[i]); !veryclose(a, f) {
t.Errorf("HypotNoSqrtGo(%g, %g) = %g, want %g", 1e200*tanh[i], 1e200*tanh[i], f, a)
}
}
for i := 0; i < len(vfhypotSC); i++ {
if f := HypotNoSqrtGo(vfhypotSC[i][0], vfhypotSC[i][1]); !alike(hypotSC[i], f) {
t.Errorf("HypotNoSqrtGo(%g, %g) = %g, want %g", vfhypotSC[i][0], vfhypotSC[i][1], f, hypotSC[i])
}
}
}
func TestIlogb(t *testing.T) {
for i := 0; i < len(vf); i++ {
a := frexp[i].i - 1 // adjust because fr in the interval [½, 1)
@ -2713,9 +2741,15 @@ func BenchmarkHypot(b *testing.B) {
}
}
func BenchmarkHypotGo(b *testing.B) {
func BenchmarkHypotNoSqrtGo(b *testing.B) {
for i := 0; i < b.N; i++ {
HypotGo(3, 4)
HypotNoSqrtGo(3, 4)
}
}
func BenchmarkHypotSqrtGo(b *testing.B) {
for i := 0; i < b.N; i++ {
HypotSqrtGo(3, 4)
}
}

View File

@ -16,7 +16,9 @@ package math
// Special cases are:
// Asin(±0) = ±0
// Asin(x) = NaN if x < -1 or x > 1
func Asin(x float64) float64 {
func Asin(x float64) float64
func asin(x float64) float64 {
if x == 0 {
return x // special case
}
@ -46,4 +48,8 @@ func Asin(x float64) float64 {
//
// Special case is:
// Acos(x) = NaN if x < -1 or x > 1
func Acos(x float64) float64 { return Pi/2 - Asin(x) }
func Acos(x float64) float64
func acos(x float64) float64 {
return Pi/2 - Asin(x)
}

View File

@ -0,0 +1,9 @@
// Copyright 2011 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.
TEXT ·Asin(SB),7,$0
JMP ·asin(SB)
TEXT ·Acos(SB),7,$0
JMP ·acos(SB)

9
src/pkg/math/asin_arm.s Normal file
View File

@ -0,0 +1,9 @@
// Copyright 2011 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.
TEXT ·Asin(SB),7,$0
B ·asin(SB)
TEXT ·Acos(SB),7,$0
B ·acos(SB)

View File

@ -1,8 +0,0 @@
// Copyright 2010 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 Acos(x float64) float64
func Asin(x float64) float64

View File

@ -51,7 +51,9 @@ func satan(arg float64) float64 {
// Special cases are:
// Atan(±0) = ±0
// Atan(±Inf) = ±Pi/2
func Atan(x float64) float64 {
func Atan(x float64) float64
func atan(x float64) float64 {
if x == 0 {
return x
}

View File

@ -26,7 +26,9 @@ package math
// Atan2(y<0, -Inf) = -Pi
// Atan2(+Inf, x) = +Pi/2
// Atan2(-Inf, x) = -Pi/2
func Atan2(y, x float64) float64 {
func Atan2(y, x float64) float64
func atan2(y, x float64) float64 {
// TODO(rsc): Remove manual inlining of IsNaN, IsInf
// when compiler does it for us
// special cases

View File

@ -1,7 +1,6 @@
// Copyright 2010 The Go Authors. All rights reserved.
// Copyright 2011 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 Exp(x float64) float64
TEXT ·Atan2(SB),7,$0
JMP ·atan2(SB)

View File

@ -1,7 +1,6 @@
// Copyright 2010 The Go Authors. All rights reserved.
// Copyright 2011 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 Tan(x float64) float64
TEXT ·Atan2(SB),7,$0
B ·atan2(SB)

View File

@ -0,0 +1,6 @@
// Copyright 2011 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.
TEXT ·Atan(SB),7,$0
JMP ·atan(SB)

6
src/pkg/math/atan_arm.s Normal file
View File

@ -0,0 +1,6 @@
// Copyright 2011 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.
TEXT ·Atan(SB),7,$0
B ·atan(SB)

View File

@ -10,8 +10,10 @@ package math
// Dim(+Inf, +Inf) = NaN
// Dim(-Inf, -Inf) = NaN
// Dim(x, NaN) = Dim(NaN, x) = NaN
func Dim(x, y float64) float64 {
return Max(x-y, 0)
func Dim(x, y float64) float64
func dim(x, y float64) float64 {
return max(x-y, 0)
}
// Max returns the larger of x or y.
@ -21,7 +23,9 @@ func Dim(x, y float64) float64 {
// Max(x, NaN) = Max(NaN, x) = NaN
// Max(+0, ±0) = Max(±0, +0) = +0
// Max(-0, -0) = -0
func Max(x, y float64) float64 {
func Max(x, y float64) float64
func max(x, y float64) float64 {
// TODO(rsc): Remove manual inlining of IsNaN, IsInf
// when compiler does it for us
// special cases
@ -48,7 +52,9 @@ func Max(x, y float64) float64 {
// Min(x, -Inf) = Min(-Inf, x) = -Inf
// Min(x, NaN) = Min(NaN, x) = NaN
// Min(-0, ±0) = Min(±0, -0) = -0
func Min(x, y float64) float64 {
func Min(x, y float64) float64
func min(x, y float64) float64 {
// TODO(rsc): Remove manual inlining of IsNaN, IsInf
// when compiler does it for us
// special cases

12
src/pkg/math/dim_386.s Normal file
View File

@ -0,0 +1,12 @@
// Copyright 2011 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.
TEXT ·Dim(SB),7,$0
JMP ·dim(SB)
TEXT ·Max(SB),7,$0
JMP ·max(SB)
TEXT ·Min(SB),7,$0
JMP ·min(SB)

12
src/pkg/math/dim_arm.s Normal file
View File

@ -0,0 +1,12 @@
// Copyright 2011 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.
TEXT ·Dim(SB),7,$0
B ·dim(SB)
TEXT ·Min(SB),7,$0
B ·min(SB)
TEXT ·Max(SB),7,$0
B ·max(SB)

View File

@ -1,9 +0,0 @@
// Copyright 2010 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 Dim(x, y float64) float64
func Max(x, y float64) float64
func Min(x, y float64) float64

View File

@ -11,4 +11,185 @@ package math
// Exp(NaN) = NaN
// Very large values overflow to 0 or +Inf.
// Very small values underflow to 1.
func Exp(x float64) float64 { return expGo(x) }
func Exp(x float64) float64
// The original C code, the long comment, and the constants
// below are from FreeBSD's /usr/src/lib/msun/src/e_exp.c
// and came with this notice. The go code is a simplified
// version of the original C.
//
// ====================================================
// Copyright (C) 2004 by Sun Microsystems, Inc. All rights reserved.
//
// Permission to use, copy, modify, and distribute this
// software is freely granted, provided that this notice
// is preserved.
// ====================================================
//
//
// exp(x)
// Returns the exponential of x.
//
// Method
// 1. Argument reduction:
// Reduce x to an r so that |r| <= 0.5*ln2 ~ 0.34658.
// Given x, find r and integer k such that
//
// x = k*ln2 + r, |r| <= 0.5*ln2.
//
// Here r will be represented as r = hi-lo for better
// accuracy.
//
// 2. Approximation of exp(r) by a special rational function on
// the interval [0,0.34658]:
// Write
// R(r**2) = r*(exp(r)+1)/(exp(r)-1) = 2 + r*r/6 - r**4/360 + ...
// We use a special Remes algorithm on [0,0.34658] to generate
// a polynomial of degree 5 to approximate R. The maximum error
// of this polynomial approximation is bounded by 2**-59. In
// other words,
// R(z) ~ 2.0 + P1*z + P2*z**2 + P3*z**3 + P4*z**4 + P5*z**5
// (where z=r*r, and the values of P1 to P5 are listed below)
// and
// | 5 | -59
// | 2.0+P1*z+...+P5*z - R(z) | <= 2
// | |
// The computation of exp(r) thus becomes
// 2*r
// exp(r) = 1 + -------
// R - r
// r*R1(r)
// = 1 + r + ----------- (for better accuracy)
// 2 - R1(r)
// where
// 2 4 10
// R1(r) = r - (P1*r + P2*r + ... + P5*r ).
//
// 3. Scale back to obtain exp(x):
// From step 1, we have
// exp(x) = 2**k * exp(r)
//
// Special cases:
// exp(INF) is INF, exp(NaN) is NaN;
// exp(-INF) is 0, and
// for finite argument, only exp(0)=1 is exact.
//
// Accuracy:
// according to an error analysis, the error is always less than
// 1 ulp (unit in the last place).
//
// Misc. info.
// For IEEE double
// if x > 7.09782712893383973096e+02 then exp(x) overflow
// if x < -7.45133219101941108420e+02 then exp(x) underflow
//
// Constants:
// The hexadecimal values are the intended ones for the following
// constants. The decimal values may be used, provided that the
// compiler will convert from decimal to binary accurately enough
// to produce the hexadecimal values shown.
func exp(x float64) float64 {
const (
Ln2Hi = 6.93147180369123816490e-01
Ln2Lo = 1.90821492927058770002e-10
Log2e = 1.44269504088896338700e+00
Overflow = 7.09782712893383973096e+02
Underflow = -7.45133219101941108420e+02
NearZero = 1.0 / (1 << 28) // 2**-28
)
// TODO(rsc): Remove manual inlining of IsNaN, IsInf
// when compiler does it for us
// special cases
switch {
case x != x || x > MaxFloat64: // IsNaN(x) || IsInf(x, 1):
return x
case x < -MaxFloat64: // IsInf(x, -1):
return 0
case x > Overflow:
return Inf(1)
case x < Underflow:
return 0
case -NearZero < x && x < NearZero:
return 1 + x
}
// reduce; computed as r = hi - lo for extra precision.
var k int
switch {
case x < 0:
k = int(Log2e*x - 0.5)
case x > 0:
k = int(Log2e*x + 0.5)
}
hi := x - float64(k)*Ln2Hi
lo := float64(k) * Ln2Lo
// compute
return expmulti(hi, lo, k)
}
// Exp2 returns 2**x, the base-2 exponential of x.
//
// Special cases are the same as Exp.
func Exp2(x float64) float64
func exp2(x float64) float64 {
const (
Ln2Hi = 6.93147180369123816490e-01
Ln2Lo = 1.90821492927058770002e-10
Overflow = 1.0239999999999999e+03
Underflow = -1.0740e+03
)
// TODO: remove manual inlining of IsNaN and IsInf
// when compiler does it for us
// special cases
switch {
case x != x || x > MaxFloat64: // IsNaN(x) || IsInf(x, 1):
return x
case x < -MaxFloat64: // IsInf(x, -1):
return 0
case x > Overflow:
return Inf(1)
case x < Underflow:
return 0
}
// argument reduction; x = r×lg(e) + k with |r| ≤ ln(2)/2.
// computed as r = hi - lo for extra precision.
var k int
switch {
case x > 0:
k = int(x + 0.5)
case x < 0:
k = int(x - 0.5)
}
t := x - float64(k)
hi := t * Ln2Hi
lo := -t * Ln2Lo
// compute
return expmulti(hi, lo, k)
}
// exp1 returns e**r × 2**k where r = hi - lo and |r| ≤ ln(2)/2.
func expmulti(hi, lo float64, k int) float64 {
const (
P1 = 1.66666666666666019037e-01 /* 0x3FC55555; 0x5555553E */
P2 = -2.77777777770155933842e-03 /* 0xBF66C16C; 0x16BEBD93 */
P3 = 6.61375632143793436117e-05 /* 0x3F11566A; 0xAF25DE2C */
P4 = -1.65339022054652515390e-06 /* 0xBEBBBD41; 0xC5D26BF1 */
P5 = 4.13813679705723846039e-08 /* 0x3E663769; 0x72BEA4D0 */
)
r := hi - lo
t := r * r
c := r - t*(P1+t*(P2+t*(P3+t*(P4+t*P5))))
y := 1 - ((lo - (r*c)/(2-c)) - hi)
// TODO(rsc): make sure Ldexp can handle boundary k
return Ldexp(y, k)
}

View File

@ -1,10 +0,0 @@
// Copyright 2010 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
// Exp2 returns 2**x, the base-2 exponential of x.
//
// Special cases are the same as Exp.
func Exp2(x float64) float64 { return exp2Go(x) }

View File

@ -0,0 +1,6 @@
// Copyright 2011 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.
TEXT ·Exp2(SB),7,$0
JMP ·exp2(SB)

6
src/pkg/math/exp2_arm.s Normal file
View File

@ -0,0 +1,6 @@
// Copyright 2011 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.
TEXT ·Exp2(SB),7,$0
B ·exp2(SB)

6
src/pkg/math/exp_arm.s Normal file
View File

@ -0,0 +1,6 @@
// Copyright 2011 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.
TEXT ·Exp(SB),7,$0
B ·exp(SB)

View File

@ -1,191 +0,0 @@
// 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
// The original C code, the long comment, and the constants
// below are from FreeBSD's /usr/src/lib/msun/src/e_exp.c
// and came with this notice. The go code is a simplified
// version of the original C.
//
// ====================================================
// Copyright (C) 2004 by Sun Microsystems, Inc. All rights reserved.
//
// Permission to use, copy, modify, and distribute this
// software is freely granted, provided that this notice
// is preserved.
// ====================================================
//
//
// exp(x)
// Returns the exponential of x.
//
// Method
// 1. Argument reduction:
// Reduce x to an r so that |r| <= 0.5*ln2 ~ 0.34658.
// Given x, find r and integer k such that
//
// x = k*ln2 + r, |r| <= 0.5*ln2.
//
// Here r will be represented as r = hi-lo for better
// accuracy.
//
// 2. Approximation of exp(r) by a special rational function on
// the interval [0,0.34658]:
// Write
// R(r**2) = r*(exp(r)+1)/(exp(r)-1) = 2 + r*r/6 - r**4/360 + ...
// We use a special Remes algorithm on [0,0.34658] to generate
// a polynomial of degree 5 to approximate R. The maximum error
// of this polynomial approximation is bounded by 2**-59. In
// other words,
// R(z) ~ 2.0 + P1*z + P2*z**2 + P3*z**3 + P4*z**4 + P5*z**5
// (where z=r*r, and the values of P1 to P5 are listed below)
// and
// | 5 | -59
// | 2.0+P1*z+...+P5*z - R(z) | <= 2
// | |
// The computation of exp(r) thus becomes
// 2*r
// exp(r) = 1 + -------
// R - r
// r*R1(r)
// = 1 + r + ----------- (for better accuracy)
// 2 - R1(r)
// where
// 2 4 10
// R1(r) = r - (P1*r + P2*r + ... + P5*r ).
//
// 3. Scale back to obtain exp(x):
// From step 1, we have
// exp(x) = 2**k * exp(r)
//
// Special cases:
// exp(INF) is INF, exp(NaN) is NaN;
// exp(-INF) is 0, and
// for finite argument, only exp(0)=1 is exact.
//
// Accuracy:
// according to an error analysis, the error is always less than
// 1 ulp (unit in the last place).
//
// Misc. info.
// For IEEE double
// if x > 7.09782712893383973096e+02 then exp(x) overflow
// if x < -7.45133219101941108420e+02 then exp(x) underflow
//
// Constants:
// The hexadecimal values are the intended ones for the following
// constants. The decimal values may be used, provided that the
// compiler will convert from decimal to binary accurately enough
// to produce the hexadecimal values shown.
// Exp returns e**x, the base-e exponential of x.
//
// Special cases are:
// Exp(+Inf) = +Inf
// Exp(NaN) = NaN
// Very large values overflow to 0 or +Inf.
// Very small values underflow to 1.
func expGo(x float64) float64 {
const (
Ln2Hi = 6.93147180369123816490e-01
Ln2Lo = 1.90821492927058770002e-10
Log2e = 1.44269504088896338700e+00
Overflow = 7.09782712893383973096e+02
Underflow = -7.45133219101941108420e+02
NearZero = 1.0 / (1 << 28) // 2**-28
)
// TODO(rsc): Remove manual inlining of IsNaN, IsInf
// when compiler does it for us
// special cases
switch {
case x != x || x > MaxFloat64: // IsNaN(x) || IsInf(x, 1):
return x
case x < -MaxFloat64: // IsInf(x, -1):
return 0
case x > Overflow:
return Inf(1)
case x < Underflow:
return 0
case -NearZero < x && x < NearZero:
return 1 + x
}
// reduce; computed as r = hi - lo for extra precision.
var k int
switch {
case x < 0:
k = int(Log2e*x - 0.5)
case x > 0:
k = int(Log2e*x + 0.5)
}
hi := x - float64(k)*Ln2Hi
lo := float64(k) * Ln2Lo
// compute
return exp(hi, lo, k)
}
// Exp2 returns 2**x, the base-2 exponential of x.
//
// Special cases are the same as Exp.
func exp2Go(x float64) float64 {
const (
Ln2Hi = 6.93147180369123816490e-01
Ln2Lo = 1.90821492927058770002e-10
Overflow = 1.0239999999999999e+03
Underflow = -1.0740e+03
)
// TODO: remove manual inlining of IsNaN and IsInf
// when compiler does it for us
// special cases
switch {
case x != x || x > MaxFloat64: // IsNaN(x) || IsInf(x, 1):
return x
case x < -MaxFloat64: // IsInf(x, -1):
return 0
case x > Overflow:
return Inf(1)
case x < Underflow:
return 0
}
// argument reduction; x = r×lg(e) + k with |r| ≤ ln(2)/2.
// computed as r = hi - lo for extra precision.
var k int
switch {
case x > 0:
k = int(x + 0.5)
case x < 0:
k = int(x - 0.5)
}
t := x - float64(k)
hi := t * Ln2Hi
lo := -t * Ln2Lo
// compute
return exp(hi, lo, k)
}
// exp returns e**r × 2**k where r = hi - lo and |r| ≤ ln(2)/2.
func exp(hi, lo float64, k int) float64 {
const (
P1 = 1.66666666666666019037e-01 /* 0x3FC55555; 0x5555553E */
P2 = -2.77777777770155933842e-03 /* 0xBF66C16C; 0x16BEBD93 */
P3 = 6.61375632143793436117e-05 /* 0x3F11566A; 0xAF25DE2C */
P4 = -1.65339022054652515390e-06 /* 0xBEBBBD41; 0xC5D26BF1 */
P5 = 4.13813679705723846039e-08 /* 0x3E663769; 0x72BEA4D0 */
)
r := hi - lo
t := r * r
c := r - t*(P1+t*(P2+t*(P3+t*(P4+t*P5))))
y := 1 - ((lo - (r*c)/(2-c)) - hi)
// TODO(rsc): make sure Ldexp can handle boundary k
return Ldexp(y, k)
}

View File

@ -1,10 +0,0 @@
// Copyright 2010 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
// Make expGo and exp2Go available for testing.
func ExpGo(x float64) float64 { return expGo(x) }
func Exp2Go(x float64) float64 { return exp2Go(x) }

View File

@ -121,7 +121,9 @@ package math
// Expm1(-Inf) = -1
// Expm1(NaN) = NaN
// Very large values overflow to -1 or +Inf.
func Expm1(x float64) float64 {
func Expm1(x float64) float64
func expm1(x float64) float64 {
const (
Othreshold = 7.09782712893383973096e+02 // 0x40862E42FEFA39EF
Ln2X56 = 3.88162421113569373274e+01 // 0x4043687a9f1af2b1

View File

@ -0,0 +1,6 @@
// Copyright 2011 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.
TEXT ·Expm1(SB),7,$0
JMP ·expm1(SB)

6
src/pkg/math/expm1_arm.s Normal file
View File

@ -0,0 +1,6 @@
// Copyright 2011 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.
TEXT ·Expm1(SB),7,$0
B ·expm1(SB)

View File

@ -1,7 +0,0 @@
// Copyright 2010 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 Expm1(x float64) float64

View File

@ -0,0 +1,12 @@
// Copyright 2011 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
// Export internal functions for testing.
var ExpGo = exp
var Exp2Go = exp2
var HypotSqrtGo = hypotSqrt
var HypotNoSqrtGo = hypotNoSqrt
var SqrtGo = sqrt

View File

@ -10,7 +10,9 @@ package math
// Floor(±0) = ±0
// Floor(±Inf) = ±Inf
// Floor(NaN) = NaN
func Floor(x float64) float64 {
func Floor(x float64) float64
func floor(x float64) float64 {
// TODO(rsc): Remove manual inlining of IsNaN, IsInf
// when compiler does it for us
if x == 0 || x != x || x > MaxFloat64 || x < -MaxFloat64 { // x == 0 || IsNaN(x) || IsInf(x, 0)
@ -33,7 +35,11 @@ func Floor(x float64) float64 {
// Ceil(±0) = ±0
// Ceil(±Inf) = ±Inf
// Ceil(NaN) = NaN
func Ceil(x float64) float64 { return -Floor(-x) }
func Ceil(x float64) float64
func ceil(x float64) float64 {
return -Floor(-x)
}
// Trunc returns the integer value of x.
//
@ -41,7 +47,9 @@ func Ceil(x float64) float64 { return -Floor(-x) }
// Trunc(±0) = ±0
// Trunc(±Inf) = ±Inf
// Trunc(NaN) = NaN
func Trunc(x float64) float64 {
func Trunc(x float64) float64
func trunc(x float64) float64 {
// TODO(rsc): Remove manual inlining of IsNaN, IsInf
// when compiler does it for us
if x == 0 || x != x || x > MaxFloat64 || x < -MaxFloat64 { // x == 0 || IsNaN(x) || IsInf(x, 0)

View File

@ -0,0 +1,12 @@
// Copyright 2011 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.
TEXT ·Floor(SB),7,$0
JMP ·floor(SB)
TEXT ·Ceil(SB),7,$0
JMP ·ceil(SB)
TEXT ·Trunc(SB),7,$0
JMP ·trunc(SB)

12
src/pkg/math/floor_arm.s Normal file
View File

@ -0,0 +1,12 @@
// Copyright 2011 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.
TEXT ·Floor(SB),7,$0
B ·floor(SB)
TEXT ·Ceil(SB),7,$0
B ·ceil(SB)
TEXT ·Trunc(SB),7,$0
B ·trunc(SB)

View File

@ -1,9 +0,0 @@
// Copyright 2010 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 Ceil(x float64) float64
func Floor(x float64) float64
func Trunc(x float64) float64

View File

@ -13,7 +13,9 @@ package math
// Frexp(±0) = ±0, 0
// Frexp(±Inf) = ±Inf, 0
// Frexp(NaN) = NaN, 0
func Frexp(f float64) (frac float64, exp int) {
func Frexp(f float64) (frac float64, exp int)
func frexp(f float64) (frac float64, exp int) {
// TODO(rsc): Remove manual inlining of IsNaN, IsInf
// when compiler does it for us
// special cases

View File

@ -0,0 +1,6 @@
// Copyright 2011 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.
TEXT ·Frexp(SB),7,$0
JMP ·frexp(SB)

6
src/pkg/math/frexp_arm.s Normal file
View File

@ -0,0 +1,6 @@
// Copyright 2011 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.
TEXT ·Frexp(SB),7,$0
B ·frexp(SB)

View File

@ -1,7 +0,0 @@
// Copyright 2010 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 Frexp(x float64) (f float64, e int)

View File

@ -14,7 +14,9 @@ package math
// Special cases are:
// Hypot(p, q) = +Inf if p or q is infinite
// Hypot(p, q) = NaN if p or q is NaN
func Hypot(p, q float64) float64 {
func Hypot(p, q float64) float64
func hypotSqrt(p, q float64) float64 {
// TODO(rsc): Remove manual inlining of IsNaN, IsInf
// when compiler does it for us
// special cases
@ -39,3 +41,46 @@ func Hypot(p, q float64) float64 {
q = q / p
return p * Sqrt(1+q*q)
}
func hypotNoSqrt(p, q float64) float64 {
// TODO(rsc): Remove manual inlining of IsNaN, IsInf
// when compiler does it for us
// special cases
switch {
case p < -MaxFloat64 || p > MaxFloat64 || q < -MaxFloat64 || q > MaxFloat64: // IsInf(p, 0) || IsInf(q, 0):
return Inf(1)
case p != p || q != q: // IsNaN(p) || IsNaN(q):
return NaN()
}
if p < 0 {
p = -p
}
if q < 0 {
q = -q
}
if p < q {
p, q = q, p
}
if p == 0 {
return 0
}
pfac := p
q = q / p
r := q
p = 1
for {
r = r * r
s := r + 4
if s == 4 {
return p * pfac
}
r = r / s
p = p + 2*r*p
q = q * r
r = q / p
}
panic("unreachable")
}

View File

@ -1,7 +1,6 @@
// Copyright 2010 The Go Authors. All rights reserved.
// Copyright 2011 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 Atan(x float64) float64
TEXT ·Hypot(SB),7,$0
B ·hypotNoSqrt(SB)

View File

@ -1,7 +0,0 @@
// Copyright 2010 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 Hypot(x, y float64) float64

View File

@ -1,63 +0,0 @@
// Copyright 2009-2010 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
/*
Hypot -- sqrt(p*p + q*q), but overflows only if the result does.
See:
Cleve Moler and Donald Morrison,
Replacing Square Roots by Pythagorean Sums
IBM Journal of Research and Development,
Vol. 27, Number 6, pp. 577-581, Nov. 1983
*/
// Hypot computes Sqrt(p*p + q*q), taking care to avoid
// unnecessary overflow and underflow.
//
// Special cases are:
// Hypot(p, q) = +Inf if p or q is infinite
// Hypot(p, q) = NaN if p or q is NaN
func hypotGo(p, q float64) float64 {
// TODO(rsc): Remove manual inlining of IsNaN, IsInf
// when compiler does it for us
// special cases
switch {
case p < -MaxFloat64 || p > MaxFloat64 || q < -MaxFloat64 || q > MaxFloat64: // IsInf(p, 0) || IsInf(q, 0):
return Inf(1)
case p != p || q != q: // IsNaN(p) || IsNaN(q):
return NaN()
}
if p < 0 {
p = -p
}
if q < 0 {
q = -q
}
if p < q {
p, q = q, p
}
if p == 0 {
return 0
}
pfac := p
q = q / p
r := q
p = 1
for {
r = r * r
s := r + 4
if s == 4 {
return p * pfac
}
r = r / s
p = p + 2*r*p
q = q * r
r = q / p
}
panic("unreachable")
}

View File

@ -1,9 +0,0 @@
// Copyright 2010 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
// Make hypotGo available for testing.
func HypotGo(x, y float64) float64 { return hypotGo(x, y) }

View File

@ -11,7 +11,9 @@ package math
// Ldexp(±0, exp) = ±0
// Ldexp(±Inf, exp) = ±Inf
// Ldexp(NaN, exp) = NaN
func Ldexp(frac float64, exp int) float64 {
func Ldexp(frac float64, exp int) float64
func ldexp(frac float64, exp int) float64 {
// TODO(rsc): Remove manual inlining of IsNaN, IsInf
// when compiler does it for us
// special cases

View File

@ -0,0 +1,6 @@
// Copyright 2011 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.
TEXT ·Ldexp(SB),7,$0
JMP ·ldexp(SB)

6
src/pkg/math/ldexp_arm.s Normal file
View File

@ -0,0 +1,6 @@
// Copyright 2011 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.
TEXT ·Ldexp(SB),7,$0
B ·ldexp(SB)

View File

@ -1,7 +0,0 @@
// Copyright 2010 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 Ldexp(f float64, e int) float64

View File

@ -77,7 +77,9 @@ package math
// Log(0) = -Inf
// Log(x < 0) = NaN
// Log(NaN) = NaN
func Log(x float64) float64 {
func Log(x float64) float64
func log(x float64) float64 {
const (
Ln2Hi = 6.93147180369123816490e-01 /* 3fe62e42 fee00000 */
Ln2Lo = 1.90821492927058770002e-10 /* 3dea39ef 35793c76 */

View File

@ -6,8 +6,16 @@ package math
// Log10 returns the decimal logarithm of x.
// The special cases are the same as for Log.
func Log10(x float64) float64 { return Log(x) * (1 / Ln10) }
func Log10(x float64) float64
func log10(x float64) float64 {
return Log(x) * (1 / Ln10)
}
// Log2 returns the binary logarithm of x.
// The special cases are the same as for Log.
func Log2(x float64) float64 { return Log(x) * (1 / Ln2) }
func Log2(x float64) float64
func log2(x float64) float64 {
return Log(x) * (1 / Ln2)
}

View File

@ -0,0 +1,9 @@
// Copyright 2011 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.
TEXT ·Log10(SB),7,$0
JMP ·log10(SB)
TEXT ·Log2(SB),7,$0
JMP ·log2(SB)

9
src/pkg/math/log10_arm.s Normal file
View File

@ -0,0 +1,9 @@
// Copyright 2011 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.
TEXT ·Log10(SB),7,$0
B ·log10(SB)
TEXT ·Log2(SB),7,$0
B ·log2(SB)

View File

@ -1,8 +0,0 @@
// Copyright 2010 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 Log10(x float64) float64
func Log2(x float64) float64

View File

@ -92,7 +92,9 @@ package math
// Log1p(-1) = -Inf
// Log1p(x < -1) = NaN
// Log1p(NaN) = NaN
func Log1p(x float64) float64 {
func Log1p(x float64) float64
func log1p(x float64) float64 {
const (
Sqrt2M1 = 4.142135623730950488017e-01 // Sqrt(2)-1 = 0x3fda827999fcef34
Sqrt2HalfM1 = -2.928932188134524755992e-01 // Sqrt(2)/2-1 = 0xbfd2bec333018866

View File

@ -0,0 +1,6 @@
// Copyright 2011 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.
TEXT ·Log1p(SB),7,$0
JMP ·log1p(SB)

6
src/pkg/math/log1p_arm.s Normal file
View File

@ -0,0 +1,6 @@
// Copyright 2011 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.
TEXT ·Log1p(SB),7,$0
B ·log1p(SB)

View File

@ -1,7 +0,0 @@
// Copyright 2010 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 Log1p(x float64) float64

6
src/pkg/math/log_arm.s Normal file
View File

@ -0,0 +1,6 @@
// Copyright 2011 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.
TEXT ·Log(SB),7,$0
B ·log(SB)

View File

@ -18,7 +18,9 @@ package math
// Mod(x, 0) = NaN
// Mod(x, ±Inf) = x
// Mod(x, NaN) = NaN
func Mod(x, y float64) float64 {
func Mod(x, y float64) float64
func mod(x, y float64) float64 {
// TODO(rsc): Remove manual inlining of IsNaN, IsInf
// when compiler does it for us.
if y == 0 || x > MaxFloat64 || x < -MaxFloat64 || x != x || y != y { // y == 0 || IsInf(x, 0) || IsNaN(x) || IsNan(y)

6
src/pkg/math/mod_amd64.s Normal file
View File

@ -0,0 +1,6 @@
// Copyright 2011 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.
TEXT ·Mod(SB),7,$0
JMP ·mod(SB)

6
src/pkg/math/mod_arm.s Normal file
View File

@ -0,0 +1,6 @@
// Copyright 2011 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.
TEXT ·Mod(SB),7,$0
B ·mod(SB)

View File

@ -1,7 +0,0 @@
// Copyright 2010 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 Mod(x, y float64) float64

View File

@ -10,7 +10,9 @@ package math
// Special cases are:
// Modf(±Inf) = ±Inf, NaN
// Modf(NaN) = NaN, NaN
func Modf(f float64) (int float64, frac float64) {
func Modf(f float64) (int float64, frac float64)
func modf(f float64) (int float64, frac float64) {
if f < 1 {
if f < 0 {
int, frac = Modf(-f)

View File

@ -0,0 +1,6 @@
// Copyright 2011 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.
TEXT ·Modf(SB),7,$0
JMP ·modf(SB)

6
src/pkg/math/modf_arm.s Normal file
View File

@ -0,0 +1,6 @@
// Copyright 2011 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.
TEXT ·Modf(SB),7,$0
B ·modf(SB)

View File

@ -1,7 +0,0 @@
// Copyright 2010 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 Modf(f float64) (int float64, frac float64)

View File

@ -34,7 +34,9 @@ package math
// Remainder(x, 0) = NaN
// Remainder(x, ±Inf) = x
// Remainder(x, NaN) = NaN
func Remainder(x, y float64) float64 {
func Remainder(x, y float64) float64
func remainder(x, y float64) float64 {
const (
Tiny = 4.45014771701440276618e-308 // 0x0020000000000000
HalfMax = MaxFloat64 / 2

View File

@ -1,7 +1,6 @@
// Copyright 2010 The Go Authors. All rights reserved.
// Copyright 2011 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 Atan2(y, x float64) float64
TEXT ·Remainder(SB),7,$0
JMP ·remainder(SB)

View File

@ -0,0 +1,6 @@
// Copyright 2011 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.
TEXT ·Remainder(SB),7,$0
B ·remainder(SB)

View File

@ -1,7 +0,0 @@
// Copyright 2010 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 Remainder(x, y float64) float64

View File

@ -113,7 +113,9 @@ var _cos = [...]float64{
// Special cases are:
// Cos(±Inf) = NaN
// Cos(NaN) = NaN
func Cos(x float64) float64 {
func Cos(x float64) float64
func cos(x float64) float64 {
const (
PI4A = 7.85398125648498535156E-1 // 0x3fe921fb40000000, Pi/4 split into three parts
PI4B = 3.77489470793079817668E-8 // 0x3e64442d00000000,
@ -170,7 +172,9 @@ func Cos(x float64) float64 {
// Sin(±0) = ±0
// Sin(±Inf) = NaN
// Sin(NaN) = NaN
func Sin(x float64) float64 {
func Sin(x float64) float64
func sin(x float64) float64 {
const (
PI4A = 7.85398125648498535156E-1 // 0x3fe921fb40000000, Pi/4 split into three parts
PI4B = 3.77489470793079817668E-8 // 0x3e64442d00000000,

View File

@ -2,6 +2,8 @@
// Use of this source code is governed by a BSD-style
// license that can be found in the LICENSE file.
package math
TEXT ·Sin(SB),7,$0
JMP ·sin(SB)
func Exp2(x float64) float64
TEXT ·Cos(SB),7,$0
JMP ·cos(SB)

View File

@ -2,6 +2,8 @@
// Use of this source code is governed by a BSD-style
// license that can be found in the LICENSE file.
package math
TEXT ·Sin(SB),7,$0
B ·sin(SB)
func Abs(x float64) float64
TEXT ·Cos(SB),7,$0
B ·cos(SB)

View File

@ -1,8 +0,0 @@
// 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 Cos(x float64) float64
func Sin(x float64) float64

View File

@ -12,7 +12,9 @@ package math
// Sincos(±0) = ±0, 1
// Sincos(±Inf) = NaN, NaN
// Sincos(NaN) = NaN, NaN
func Sincos(x float64) (sin, cos float64) {
func Sincos(x float64) (sin, cos float64)
func sincos(x float64) (sin, cos float64) {
const (
PI4A = 7.85398125648498535156E-1 // 0x3fe921fb40000000, Pi/4 split into three parts
PI4B = 3.77489470793079817668E-8 // 0x3e64442d00000000,

View File

@ -0,0 +1,6 @@
// Copyright 2011 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.
TEXT ·Sincos(SB),7,$0
B ·sincos(SB)

View File

@ -1,7 +0,0 @@
// Copyright 2010 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 Sincos(x float64) (sin, cos float64)

View File

@ -11,4 +11,142 @@ package math
// Sqrt(±0) = ±0
// Sqrt(x < 0) = NaN
// Sqrt(NaN) = NaN
func Sqrt(x float64) float64 { return sqrtGo(x) }
func Sqrt(x float64) float64
// The original C code and the long comment below are
// from FreeBSD's /usr/src/lib/msun/src/e_sqrt.c and
// came with this notice. The go code is a simplified
// version of the original C.
//
// ====================================================
// Copyright (C) 1993 by Sun Microsystems, Inc. All rights reserved.
//
// Developed at SunPro, a Sun Microsystems, Inc. business.
// Permission to use, copy, modify, and distribute this
// software is freely granted, provided that this notice
// is preserved.
// ====================================================
//
// __ieee754_sqrt(x)
// Return correctly rounded sqrt.
// -----------------------------------------
// | Use the hardware sqrt if you have one |
// -----------------------------------------
// Method:
// Bit by bit method using integer arithmetic. (Slow, but portable)
// 1. Normalization
// Scale x to y in [1,4) with even powers of 2:
// find an integer k such that 1 <= (y=x*2**(2k)) < 4, then
// sqrt(x) = 2**k * sqrt(y)
// 2. Bit by bit computation
// Let q = sqrt(y) truncated to i bit after binary point (q = 1),
// i 0
// i+1 2
// s = 2*q , and y = 2 * ( y - q ). (1)
// i i i i
//
// To compute q from q , one checks whether
// i+1 i
//
// -(i+1) 2
// (q + 2 ) <= y. (2)
// i
// -(i+1)
// If (2) is false, then q = q ; otherwise q = q + 2 .
// i+1 i i+1 i
//
// With some algebraic manipulation, it is not difficult to see
// that (2) is equivalent to
// -(i+1)
// s + 2 <= y (3)
// i i
//
// The advantage of (3) is that s and y can be computed by
// i i
// the following recurrence formula:
// if (3) is false
//
// s = s , y = y ; (4)
// i+1 i i+1 i
//
// otherwise,
// -i -(i+1)
// s = s + 2 , y = y - s - 2 (5)
// i+1 i i+1 i i
//
// One may easily use induction to prove (4) and (5).
// Note. Since the left hand side of (3) contain only i+2 bits,
// it does not necessary to do a full (53-bit) comparison
// in (3).
// 3. Final rounding
// After generating the 53 bits result, we compute one more bit.
// Together with the remainder, we can decide whether the
// result is exact, bigger than 1/2ulp, or less than 1/2ulp
// (it will never equal to 1/2ulp).
// The rounding mode can be detected by checking whether
// huge + tiny is equal to huge, and whether huge - tiny is
// equal to huge for some floating point number "huge" and "tiny".
//
//
// Notes: Rounding mode detection omitted. The constants "mask", "shift",
// and "bias" are found in src/pkg/math/bits.go
// Sqrt returns the square root of x.
//
// Special cases are:
// Sqrt(+Inf) = +Inf
// Sqrt(±0) = ±0
// Sqrt(x < 0) = NaN
// Sqrt(NaN) = NaN
func sqrt(x float64) float64 {
// special cases
// TODO(rsc): Remove manual inlining of IsNaN, IsInf
// when compiler does it for us
switch {
case x == 0 || x != x || x > MaxFloat64: // x == 0 || IsNaN(x) || IsInf(x, 1):
return x
case x < 0:
return NaN()
}
ix := Float64bits(x)
// normalize x
exp := int((ix >> shift) & mask)
if exp == 0 { // subnormal x
for ix&1<<shift == 0 {
ix <<= 1
exp--
}
exp++
}
exp -= bias // unbias exponent
ix &^= mask << shift
ix |= 1 << shift
if exp&1 == 1 { // odd exp, double x to make it even
ix <<= 1
}
exp >>= 1 // exp = exp/2, exponent of square root
// generate sqrt(x) bit by bit
ix <<= 1
var q, s uint64 // q = sqrt(x)
r := uint64(1 << (shift + 1)) // r = moving bit from MSB to LSB
for r != 0 {
t := s + r
if t <= ix {
s = t + r
ix -= t
q += r
}
ix <<= 1
r >>= 1
}
// final rounding
if ix != 0 { // remainder, result not exact
q += q & 1 // round according to extra bit
}
ix = q>>1 + uint64(exp-1+bias)<<shift // significand + biased exponent
return Float64frombits(ix)
}
func sqrtC(f float64, r *float64) {
*r = sqrt(f)
}

View File

@ -1,7 +0,0 @@
// 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 Sqrt(x float64) float64

View File

@ -1,147 +0,0 @@
// 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
/*
Floating-point square root.
*/
// The original C code and the long comment below are
// from FreeBSD's /usr/src/lib/msun/src/e_sqrt.c and
// came with this notice. The go code is a simplified
// version of the original C.
//
// ====================================================
// Copyright (C) 1993 by Sun Microsystems, Inc. All rights reserved.
//
// Developed at SunPro, a Sun Microsystems, Inc. business.
// Permission to use, copy, modify, and distribute this
// software is freely granted, provided that this notice
// is preserved.
// ====================================================
//
// __ieee754_sqrt(x)
// Return correctly rounded sqrt.
// -----------------------------------------
// | Use the hardware sqrt if you have one |
// -----------------------------------------
// Method:
// Bit by bit method using integer arithmetic. (Slow, but portable)
// 1. Normalization
// Scale x to y in [1,4) with even powers of 2:
// find an integer k such that 1 <= (y=x*2**(2k)) < 4, then
// sqrt(x) = 2**k * sqrt(y)
// 2. Bit by bit computation
// Let q = sqrt(y) truncated to i bit after binary point (q = 1),
// i 0
// i+1 2
// s = 2*q , and y = 2 * ( y - q ). (1)
// i i i i
//
// To compute q from q , one checks whether
// i+1 i
//
// -(i+1) 2
// (q + 2 ) <= y. (2)
// i
// -(i+1)
// If (2) is false, then q = q ; otherwise q = q + 2 .
// i+1 i i+1 i
//
// With some algebraic manipulation, it is not difficult to see
// that (2) is equivalent to
// -(i+1)
// s + 2 <= y (3)
// i i
//
// The advantage of (3) is that s and y can be computed by
// i i
// the following recurrence formula:
// if (3) is false
//
// s = s , y = y ; (4)
// i+1 i i+1 i
//
// otherwise,
// -i -(i+1)
// s = s + 2 , y = y - s - 2 (5)
// i+1 i i+1 i i
//
// One may easily use induction to prove (4) and (5).
// Note. Since the left hand side of (3) contain only i+2 bits,
// it does not necessary to do a full (53-bit) comparison
// in (3).
// 3. Final rounding
// After generating the 53 bits result, we compute one more bit.
// Together with the remainder, we can decide whether the
// result is exact, bigger than 1/2ulp, or less than 1/2ulp
// (it will never equal to 1/2ulp).
// The rounding mode can be detected by checking whether
// huge + tiny is equal to huge, and whether huge - tiny is
// equal to huge for some floating point number "huge" and "tiny".
//
//
// Notes: Rounding mode detection omitted. The constants "mask", "shift",
// and "bias" are found in src/pkg/math/bits.go
// Sqrt returns the square root of x.
//
// Special cases are:
// Sqrt(+Inf) = +Inf
// Sqrt(±0) = ±0
// Sqrt(x < 0) = NaN
// Sqrt(NaN) = NaN
func sqrtGo(x float64) float64 {
// special cases
// TODO(rsc): Remove manual inlining of IsNaN, IsInf
// when compiler does it for us
switch {
case x == 0 || x != x || x > MaxFloat64: // x == 0 || IsNaN(x) || IsInf(x, 1):
return x
case x < 0:
return NaN()
}
ix := Float64bits(x)
// normalize x
exp := int((ix >> shift) & mask)
if exp == 0 { // subnormal x
for ix&1<<shift == 0 {
ix <<= 1
exp--
}
exp++
}
exp -= bias // unbias exponent
ix &^= mask << shift
ix |= 1 << shift
if exp&1 == 1 { // odd exp, double x to make it even
ix <<= 1
}
exp >>= 1 // exp = exp/2, exponent of square root
// generate sqrt(x) bit by bit
ix <<= 1
var q, s uint64 // q = sqrt(x)
r := uint64(1 << (shift + 1)) // r = moving bit from MSB to LSB
for r != 0 {
t := s + r
if t <= ix {
s = t + r
ix -= t
q += r
}
ix <<= 1
r >>= 1
}
// final rounding
if ix != 0 { // remainder, result not exact
q += q & 1 // round according to extra bit
}
ix = q>>1 + uint64(exp-1+bias)<<shift // significand + biased exponent
return Float64frombits(ix)
}
func sqrtGoC(f float64, r *float64) {
*r = sqrtGo(f)
}

View File

@ -1,9 +0,0 @@
// Copyright 2010 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
// Make sqrtGo available for testing.
func SqrtGo(x float64) float64 { return sqrtGo(x) }

View File

@ -79,7 +79,9 @@ var _tanQ = [...]float64{
// Tan(±0) = ±0
// Tan(±Inf) = NaN
// Tan(NaN) = NaN
func Tan(x float64) float64 {
func Tan(x float64) float64
func tan(x float64) float64 {
const (
PI4A = 7.85398125648498535156E-1 // 0x3fe921fb40000000, Pi/4 split into three parts
PI4B = 3.77489470793079817668E-8 // 0x3e64442d00000000,

6
src/pkg/math/tan_amd64.s Normal file
View File

@ -0,0 +1,6 @@
// Copyright 2011 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.
TEXT ·Tan(SB),7,$0
JMP ·tan(SB)

6
src/pkg/math/tan_arm.s Normal file
View File

@ -0,0 +1,6 @@
// Copyright 2011 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.
TEXT ·Tan(SB),7,$0
B ·tan(SB)

View File

@ -15,7 +15,7 @@
#define FLAGS_V (1 << 28)
void runtime·abort(void);
void math·sqrtGoC(uint64, uint64*);
void math·sqrtC(uint64, uint64*);
static uint32 trace = 0;
@ -359,7 +359,7 @@ stage3: // regd, regm are 4bit variables
break;
case 0xeeb10bc0: // D[regd] = sqrt D[regm]
math·sqrtGoC(getd(regm), &uval);
math·sqrtC(getd(regm), &uval);
putd(regd, uval);
if(trace)