diff options
| author | Ian Lance Taylor <iant@golang.org> | 2018-09-24 21:46:21 +0000 |
|---|---|---|
| committer | Ian Lance Taylor <ian@gcc.gnu.org> | 2018-09-24 21:46:21 +0000 |
| commit | dd931d9b48647e898dc80927c532ae93cc09e192 (patch) | |
| tree | 71be2295cd79b8a182f6130611658db8628772d5 /libgo/go/math | |
| parent | 779d8a5ad09b01428726ea5a0e6c87bd9ac3c0e4 (diff) | |
| download | gcc-dd931d9b48647e898dc80927c532ae93cc09e192.zip gcc-dd931d9b48647e898dc80927c532ae93cc09e192.tar.gz gcc-dd931d9b48647e898dc80927c532ae93cc09e192.tar.bz2 | |
libgo: update to Go 1.11
Reviewed-on: https://go-review.googlesource.com/136435
gotools/:
* Makefile.am (mostlyclean-local): Run chmod on check-go-dir to
make sure it is writable.
(check-go-tools): Likewise.
(check-vet): Copy internal/objabi to check-vet-dir.
* Makefile.in: Rebuild.
From-SVN: r264546
Diffstat (limited to 'libgo/go/math')
35 files changed, 1039 insertions, 587 deletions
diff --git a/libgo/go/math/all_test.go b/libgo/go/math/all_test.go index 98437b0..00f2058 100644 --- a/libgo/go/math/all_test.go +++ b/libgo/go/math/all_test.go @@ -25,7 +25,7 @@ var vf = []float64{ } // The expected results below were computed by the high precision calculators -// at http://keisan.casio.com/. More exact input values (array vf[], above) +// at https://keisan.casio.com/. More exact input values (array vf[], above) // were obtained by printing them with "%.26f". The answers were calculated // to 26 digits (by using the "Digit number" drop-down control of each // calculator). @@ -946,6 +946,8 @@ var vferfSC = []float64{ 0, Inf(1), NaN(), + -1000, + 1000, } var erfSC = []float64{ -1, @@ -953,17 +955,23 @@ var erfSC = []float64{ 0, 1, NaN(), + -1, + 1, } var vferfcSC = []float64{ Inf(-1), Inf(1), NaN(), + -1000, + 1000, } var erfcSC = []float64{ 2, 0, NaN(), + 2, + 0, } var vferfinvSC = []float64{ @@ -1012,6 +1020,10 @@ var vfexpSC = []float64{ 1.48852223e+09, 1.4885222e+09, 1, + // near zero + 3.725290298461915e-09, + // denormal + -740, } var expSC = []float64{ 0, @@ -1023,6 +1035,8 @@ var expSC = []float64{ Inf(1), Inf(1), 2.718281828459045, + 1.0000000037252903, + 4.2e-322, } var vfexp2SC = []float64{ @@ -1033,6 +1047,10 @@ var vfexp2SC = []float64{ NaN(), // smallest float64 that overflows Exp2(x) 1024, + // near underflow + -1.07399999999999e+03, + // near zero + 3.725290298461915e-09, } var exp2SC = []float64{ 0, @@ -1041,6 +1059,8 @@ var exp2SC = []float64{ Inf(1), NaN(), Inf(1), + 5e-324, + 1.0000000025821745, } var vfexpm1SC = []float64{ @@ -1955,6 +1975,8 @@ var vfldexpBC = []fi{ {-1, -1075}, {1, 1024}, {-1, 1024}, + {1.0000000000000002, -1075}, + {1, -1075}, } var ldexpBC = []float64{ SmallestNonzeroFloat64, @@ -1965,6 +1987,8 @@ var ldexpBC = []float64{ Copysign(0, -1), Inf(1), Inf(-1), + SmallestNonzeroFloat64, + 0, } var logbBC = []float64{ @@ -2316,7 +2340,7 @@ func testExp2(t *testing.T, Exp2 func(float64) float64, name string) { } for i := 0; i < len(vfexp2SC); i++ { if f := Exp2(vfexp2SC[i]); !alike(exp2SC[i], f) { - t.Errorf("%s(%g) = %g, want %g", name, vfexpSC[i], f, expSC[i]) + t.Errorf("%s(%g) = %g, want %g", name, vfexp2SC[i], f, exp2SC[i]) } } for n := -1074; n < 1024; n++ { @@ -2419,6 +2443,10 @@ func TestMod(t *testing.T) { t.Errorf("Mod(%g, %g) = %g, want %g", vffmodSC[i][0], vffmodSC[i][1], f, fmodSC[i]) } } + // verify precision of result for extreme inputs + if f := Mod(5.9790119248836734e+200, 1.1258465975523544); 0.6447968302508578 != f { + t.Errorf("Remainder(5.9790119248836734e+200, 1.1258465975523544) = %g, want 0.6447968302508578", f) + } } func TestFrexp(t *testing.T) { @@ -2760,6 +2788,10 @@ func TestRemainder(t *testing.T) { t.Errorf("Remainder(%g, %g) = %g, want %g", vffmodSC[i][0], vffmodSC[i][1], f, fmodSC[i]) } } + // verify precision of result for extreme inputs + if f := Remainder(5.9790119248836734e+200, 1.1258465975523544); -0.4810497673014966 != f { + t.Errorf("Remainder(5.9790119248836734e+200, 1.1258465975523544) = %g, want -0.4810497673014966", f) + } } func TestRound(t *testing.T) { diff --git a/libgo/go/math/big/accuracy_string.go b/libgo/go/math/big/accuracy_string.go index 24ef7f1..1501ace 100644 --- a/libgo/go/math/big/accuracy_string.go +++ b/libgo/go/math/big/accuracy_string.go @@ -1,8 +1,8 @@ -// generated by stringer -type=Accuracy; DO NOT EDIT +// Code generated by "stringer -type=Accuracy"; DO NOT EDIT. package big -import "fmt" +import "strconv" const _Accuracy_name = "BelowExactAbove" @@ -10,8 +10,8 @@ var _Accuracy_index = [...]uint8{0, 5, 10, 15} func (i Accuracy) String() string { i -= -1 - if i < 0 || i+1 >= Accuracy(len(_Accuracy_index)) { - return fmt.Sprintf("Accuracy(%d)", i+-1) + if i < 0 || i >= Accuracy(len(_Accuracy_index)-1) { + return "Accuracy(" + strconv.FormatInt(int64(i+-1), 10) + ")" } return _Accuracy_name[_Accuracy_index[i]:_Accuracy_index[i+1]] } diff --git a/libgo/go/math/big/arith_amd64.go b/libgo/go/math/big/arith_amd64.go new file mode 100644 index 0000000..35164f3 --- /dev/null +++ b/libgo/go/math/big/arith_amd64.go @@ -0,0 +1,12 @@ +// Copyright 2017 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. + +// +build ignore_for_gccgo +// +build !math_big_pure_go + +package big + +import "internal/cpu" + +var support_adx = cpu.X86.HasADX && cpu.X86.HasBMI2 diff --git a/libgo/go/math/big/arith_test.go b/libgo/go/math/big/arith_test.go index 13b0436..cf386b3 100644 --- a/libgo/go/math/big/arith_test.go +++ b/libgo/go/math/big/arith_test.go @@ -142,6 +142,23 @@ func BenchmarkAddVV(b *testing.B) { } } +func BenchmarkSubVV(b *testing.B) { + for _, n := range benchSizes { + if isRaceBuilder && n > 1e3 { + continue + } + x := rndV(n) + y := rndV(n) + z := make([]Word, n) + b.Run(fmt.Sprint(n), func(b *testing.B) { + b.SetBytes(int64(n * _W)) + for i := 0; i < b.N; i++ { + subVV(z, x, y) + } + }) + } +} + type funVW func(z, x []Word, y Word) (c Word) type argVW struct { z, x nat @@ -255,6 +272,23 @@ func BenchmarkAddVW(b *testing.B) { } } +func BenchmarkSubVW(b *testing.B) { + for _, n := range benchSizes { + if isRaceBuilder && n > 1e3 { + continue + } + x := rndV(n) + y := rndW() + z := make([]Word, n) + b.Run(fmt.Sprint(n), func(b *testing.B) { + b.SetBytes(int64(n * _S)) + for i := 0; i < b.N; i++ { + subVW(z, x, y) + } + }) + } +} + type funVWW func(z, x []Word, y, r Word) (c Word) type argVWW struct { z, x nat diff --git a/libgo/go/math/big/calibrate_test.go b/libgo/go/math/big/calibrate_test.go index 2b96e74..4fa663f 100644 --- a/libgo/go/math/big/calibrate_test.go +++ b/libgo/go/math/big/calibrate_test.go @@ -28,24 +28,32 @@ import ( var calibrate = flag.Bool("calibrate", false, "run calibration test") +const ( + sqrModeMul = "mul(x, x)" + sqrModeBasic = "basicSqr(x)" + sqrModeKaratsuba = "karatsubaSqr(x)" +) + func TestCalibrate(t *testing.T) { - if *calibrate { - computeKaratsubaThresholds() - - // compute basicSqrThreshold where overhead becomes negligible - minSqr := computeSqrThreshold(10, 30, 1, 3) - // compute karatsubaSqrThreshold where karatsuba is faster - maxSqr := computeSqrThreshold(300, 500, 10, 3) - if minSqr != 0 { - fmt.Printf("found basicSqrThreshold = %d\n", minSqr) - } else { - fmt.Println("no basicSqrThreshold found") - } - if maxSqr != 0 { - fmt.Printf("found karatsubaSqrThreshold = %d\n", maxSqr) - } else { - fmt.Println("no karatsubaSqrThreshold found") - } + if !*calibrate { + return + } + + computeKaratsubaThresholds() + + // compute basicSqrThreshold where overhead becomes negligible + minSqr := computeSqrThreshold(10, 30, 1, 3, sqrModeMul, sqrModeBasic) + // compute karatsubaSqrThreshold where karatsuba is faster + maxSqr := computeSqrThreshold(200, 500, 10, 3, sqrModeBasic, sqrModeKaratsuba) + if minSqr != 0 { + fmt.Printf("found basicSqrThreshold = %d\n", minSqr) + } else { + fmt.Println("no basicSqrThreshold found") + } + if maxSqr != 0 { + fmt.Printf("found karatsubaSqrThreshold = %d\n", maxSqr) + } else { + fmt.Println("no karatsubaSqrThreshold found") } } @@ -109,16 +117,17 @@ func computeKaratsubaThresholds() { } } -func measureBasicSqr(words, nruns int, enable bool) time.Duration { +func measureSqr(words, nruns int, mode string) time.Duration { // more runs for better statistics initBasicSqr, initKaratsubaSqr := basicSqrThreshold, karatsubaSqrThreshold - if enable { - // set thresholds to use basicSqr at this number of words + switch mode { + case sqrModeMul: + basicSqrThreshold = words + 1 + case sqrModeBasic: basicSqrThreshold, karatsubaSqrThreshold = words-1, words+1 - } else { - // set thresholds to disable basicSqr for any number of words - basicSqrThreshold, karatsubaSqrThreshold = -1, -1 + case sqrModeKaratsuba: + karatsubaSqrThreshold = words - 1 } var testval int64 @@ -133,18 +142,18 @@ func measureBasicSqr(words, nruns int, enable bool) time.Duration { return time.Duration(testval) } -func computeSqrThreshold(from, to, step, nruns int) int { - fmt.Println("Calibrating thresholds for basicSqr via benchmarks of z.mul(x,x)") +func computeSqrThreshold(from, to, step, nruns int, lower, upper string) int { + fmt.Printf("Calibrating threshold between %s and %s\n", lower, upper) fmt.Printf("Looking for a timing difference for x between %d - %d words by %d step\n", from, to, step) var initPos bool var threshold int for i := from; i <= to; i += step { - baseline := measureBasicSqr(i, nruns, false) - testval := measureBasicSqr(i, nruns, true) + baseline := measureSqr(i, nruns, lower) + testval := measureSqr(i, nruns, upper) pos := baseline > testval delta := baseline - testval percent := delta * 100 / baseline - fmt.Printf("words = %3d deltaT = %10s (%4d%%) is basicSqr better: %v", i, delta, percent, pos) + fmt.Printf("words = %3d deltaT = %10s (%4d%%) is %s better: %v", i, delta, percent, upper, pos) if i == from { initPos = pos } diff --git a/libgo/go/math/big/float.go b/libgo/go/math/big/float.go index c042854..55b93c8 100644 --- a/libgo/go/math/big/float.go +++ b/libgo/go/math/big/float.go @@ -3,7 +3,7 @@ // license that can be found in the LICENSE file. // This file implements multi-precision floating-point numbers. -// Like in the GNU MPFR library (http://www.mpfr.org/), operands +// Like in the GNU MPFR library (https://www.mpfr.org/), operands // can be of mixed precision. Unlike MPFR, the rounding mode is // not specified with each operation, but with each operand. The // rounding mode of the result operand determines the rounding @@ -1429,8 +1429,6 @@ func (x *Float) ucmp(y *Float) int { // z's accuracy reports the result error relative to the exact (not rounded) // result. Add panics with ErrNaN if x and y are infinities with opposite // signs. The value of z is undefined in that case. -// -// BUG(gri) When rounding ToNegativeInf, the sign of Float values rounded to 0 is incorrect. func (z *Float) Add(x, y *Float) *Float { if debugFloat { x.validate() @@ -1466,6 +1464,9 @@ func (z *Float) Add(x, y *Float) *Float { z.usub(y, x) } } + if z.form == zero && z.mode == ToNegativeInf && z.acc == Exact { + z.neg = true + } return z } @@ -1530,6 +1531,9 @@ func (z *Float) Sub(x, y *Float) *Float { z.usub(y, x) } } + if z.form == zero && z.mode == ToNegativeInf && z.acc == Exact { + z.neg = true + } return z } diff --git a/libgo/go/math/big/float_test.go b/libgo/go/math/big/float_test.go index 5fd49bb..7d6bf03 100644 --- a/libgo/go/math/big/float_test.go +++ b/libgo/go/math/big/float_test.go @@ -1007,9 +1007,9 @@ func TestFloatFloat64(t *testing.T) { {"0x.fffffffffffffp-1022", smallestNormalFloat64 - math.SmallestNonzeroFloat64, Exact}, {"4503599627370495p-1074", smallestNormalFloat64 - math.SmallestNonzeroFloat64, Exact}, - // http://www.exploringbinary.com/php-hangs-on-numeric-value-2-2250738585072011e-308/ + // https://www.exploringbinary.com/php-hangs-on-numeric-value-2-2250738585072011e-308/ {"2.2250738585072011e-308", 2.225073858507201e-308, Below}, - // http://www.exploringbinary.com/java-hangs-when-converting-2-2250738585072012e-308/ + // https://www.exploringbinary.com/java-hangs-when-converting-2-2250738585072012e-308/ {"2.2250738585072012e-308", 2.2250738585072014e-308, Above}, } { for i := 0; i < 2; i++ { @@ -1258,6 +1258,31 @@ func TestFloatAdd(t *testing.T) { } } +// TestFloatAddRoundZero tests Float.Add/Sub rounding when the result is exactly zero. +// x + (-x) or x - x for non-zero x should be +0 in all cases except when +// the rounding mode is ToNegativeInf in which case it should be -0. +func TestFloatAddRoundZero(t *testing.T) { + for _, mode := range [...]RoundingMode{ToNearestEven, ToNearestAway, ToZero, AwayFromZero, ToPositiveInf, ToNegativeInf} { + x := NewFloat(5.0) + y := new(Float).Neg(x) + want := NewFloat(0.0) + if mode == ToNegativeInf { + want.Neg(want) + } + got := new(Float).SetMode(mode) + got.Add(x, y) + if got.Cmp(want) != 0 || got.neg != (mode == ToNegativeInf) { + t.Errorf("%s:\n\t %v\n\t+ %v\n\t= %v\n\twant %v", + mode, x, y, got, want) + } + got.Sub(x, x) + if got.Cmp(want) != 0 || got.neg != (mode == ToNegativeInf) { + t.Errorf("%v:\n\t %v\n\t- %v\n\t= %v\n\twant %v", + mode, x, x, got, want) + } + } +} + // TestFloatAdd32 tests that Float.Add/Sub of numbers with // 24bit mantissa behaves like float32 addition/subtraction // (excluding denormal numbers). @@ -1372,7 +1397,7 @@ func TestFloatMul(t *testing.T) { got.Mul(x, y) want := zbits.round(prec, mode) if got.Cmp(want) != 0 { - t.Errorf("i = %d, prec = %d, %s:\n\t %s %v\n\t* %s %v\n\t= %s\n\twant %s", + t.Errorf("i = %d, prec = %d, %s:\n\t %v %v\n\t* %v %v\n\t= %v\n\twant %v", i, prec, mode, x, xbits, y, ybits, got, want) } @@ -1382,7 +1407,7 @@ func TestFloatMul(t *testing.T) { got.Quo(z, x) want = ybits.round(prec, mode) if got.Cmp(want) != 0 { - t.Errorf("i = %d, prec = %d, %s:\n\t %s %v\n\t/ %s %v\n\t= %s\n\twant %s", + t.Errorf("i = %d, prec = %d, %s:\n\t %v %v\n\t/ %v %v\n\t= %v\n\twant %v", i, prec, mode, z, zbits, x, xbits, got, want) } } @@ -1465,13 +1490,13 @@ func TestIssue6866(t *testing.T) { z2.Sub(two, p) if z1.Cmp(z2) != 0 { - t.Fatalf("prec %d: got z1 = %s != z2 = %s; want z1 == z2\n", prec, z1, z2) + t.Fatalf("prec %d: got z1 = %v != z2 = %v; want z1 == z2\n", prec, z1, z2) } if z1.Sign() != 0 { - t.Errorf("prec %d: got z1 = %s; want 0", prec, z1) + t.Errorf("prec %d: got z1 = %v; want 0", prec, z1) } if z2.Sign() != 0 { - t.Errorf("prec %d: got z2 = %s; want 0", prec, z2) + t.Errorf("prec %d: got z2 = %v; want 0", prec, z2) } } } diff --git a/libgo/go/math/big/floatconv_test.go b/libgo/go/math/big/floatconv_test.go index 6d0f17d..269e265 100644 --- a/libgo/go/math/big/floatconv_test.go +++ b/libgo/go/math/big/floatconv_test.go @@ -286,9 +286,9 @@ func TestFloat64Text(t *testing.T) { {0.5, 'f', 0, "0"}, {1.5, 'f', 0, "2"}, - // http://www.exploringbinary.com/java-hangs-when-converting-2-2250738585072012e-308/ + // https://www.exploringbinary.com/java-hangs-when-converting-2-2250738585072012e-308/ {2.2250738585072012e-308, 'g', -1, "2.2250738585072014e-308"}, - // http://www.exploringbinary.com/php-hangs-on-numeric-value-2-2250738585072011e-308/ + // https://www.exploringbinary.com/php-hangs-on-numeric-value-2-2250738585072011e-308/ {2.2250738585072011e-308, 'g', -1, "2.225073858507201e-308"}, // Issue 2625. diff --git a/libgo/go/math/big/int.go b/libgo/go/math/big/int.go index 0eda9cd..47a288a 100644 --- a/libgo/go/math/big/int.go +++ b/libgo/go/math/big/int.go @@ -442,24 +442,28 @@ func (x *Int) BitLen() int { } // Exp sets z = x**y mod |m| (i.e. the sign of m is ignored), and returns z. -// If y <= 0, the result is 1 mod |m|; if m == nil or m == 0, z = x**y. +// If m == nil or m == 0, z = x**y unless y <= 0 then z = 1. // // Modular exponentation of inputs of a particular size is not a // cryptographically constant-time operation. func (z *Int) Exp(x, y, m *Int) *Int { // See Knuth, volume 2, section 4.6.3. - var yWords nat - if !y.neg { - yWords = y.abs + xWords := x.abs + if y.neg { + if m == nil || len(m.abs) == 0 { + return z.SetInt64(1) + } + // for y < 0: x**y mod m == (x**(-1))**|y| mod m + xWords = new(Int).ModInverse(x, m).abs } - // y >= 0 + yWords := y.abs var mWords nat if m != nil { mWords = m.abs // m.abs may be nil for m == 0 } - z.abs = z.abs.expNN(x.abs, yWords, mWords) + z.abs = z.abs.expNN(xWords, yWords, mWords) z.neg = len(z.abs) > 0 && x.neg && len(yWords) > 0 && yWords[0]&1 == 1 // 0 has no sign if z.neg && len(mWords) > 0 { // make modulus result positive @@ -485,162 +489,223 @@ func (z *Int) GCD(x, y, a, b *Int) *Int { } return z } - if x == nil && y == nil { - return z.lehmerGCD(a, b) - } - A := new(Int).Set(a) - B := new(Int).Set(b) - - X := new(Int) - lastX := new(Int).SetInt64(1) + return z.lehmerGCD(x, y, a, b) +} - q := new(Int) - temp := new(Int) +// lehmerSimulate attempts to simulate several Euclidean update steps +// using the leading digits of A and B. It returns u0, u1, v0, v1 +// such that A and B can be updated as: +// A = u0*A + v0*B +// B = u1*A + v1*B +// Requirements: A >= B and len(B.abs) >= 2 +// Since we are calculating with full words to avoid overflow, +// we use 'even' to track the sign of the cosequences. +// For even iterations: u0, v1 >= 0 && u1, v0 <= 0 +// For odd iterations: u0, v1 <= 0 && u1, v0 >= 0 +func lehmerSimulate(A, B *Int) (u0, u1, v0, v1 Word, even bool) { + // initialize the digits + var a1, a2, u2, v2 Word + + m := len(B.abs) // m >= 2 + n := len(A.abs) // n >= m >= 2 + + // extract the top Word of bits from A and B + h := nlz(A.abs[n-1]) + a1 = A.abs[n-1]<<h | A.abs[n-2]>>(_W-h) + // B may have implicit zero words in the high bits if the lengths differ + switch { + case n == m: + a2 = B.abs[n-1]<<h | B.abs[n-2]>>(_W-h) + case n == m+1: + a2 = B.abs[n-2] >> (_W - h) + default: + a2 = 0 + } - r := new(Int) - for len(B.abs) > 0 { - q, r = q.QuoRem(A, B, r) + // Since we are calculating with full words to avoid overflow, + // we use 'even' to track the sign of the cosequences. + // For even iterations: u0, v1 >= 0 && u1, v0 <= 0 + // For odd iterations: u0, v1 <= 0 && u1, v0 >= 0 + // The first iteration starts with k=1 (odd). + even = false + // variables to track the cosequences + u0, u1, u2 = 0, 1, 0 + v0, v1, v2 = 0, 0, 1 + + // Calculate the quotient and cosequences using Collins' stopping condition. + // Note that overflow of a Word is not possible when computing the remainder + // sequence and cosequences since the cosequence size is bounded by the input size. + // See section 4.2 of Jebelean for details. + for a2 >= v2 && a1-a2 >= v1+v2 { + q, r := a1/a2, a1%a2 + a1, a2 = a2, r + u0, u1, u2 = u1, u2, u1+q*u2 + v0, v1, v2 = v1, v2, v1+q*v2 + even = !even + } + return +} - A, B, r = B, r, A +// lehmerUpdate updates the inputs A and B such that: +// A = u0*A + v0*B +// B = u1*A + v1*B +// where the signs of u0, u1, v0, v1 are given by even +// For even == true: u0, v1 >= 0 && u1, v0 <= 0 +// For even == false: u0, v1 <= 0 && u1, v0 >= 0 +// q, r, s, t are temporary variables to avoid allocations in the multiplication +func lehmerUpdate(A, B, q, r, s, t *Int, u0, u1, v0, v1 Word, even bool) { + + t.abs = t.abs.setWord(u0) + s.abs = s.abs.setWord(v0) + t.neg = !even + s.neg = even + + t.Mul(A, t) + s.Mul(B, s) + + r.abs = r.abs.setWord(u1) + q.abs = q.abs.setWord(v1) + r.neg = even + q.neg = !even + + r.Mul(A, r) + q.Mul(B, q) + + A.Add(t, s) + B.Add(r, q) +} - temp.Set(X) - X.Mul(X, q) - X.Sub(lastX, X) - lastX.Set(temp) - } +// euclidUpdate performs a single step of the Euclidean GCD algorithm +// if extended is true, it also updates the cosequence Ua, Ub +func euclidUpdate(A, B, Ua, Ub, q, r, s, t *Int, extended bool) { + q, r = q.QuoRem(A, B, r) - if x != nil { - *x = *lastX - } + *A, *B, *r = *B, *r, *A - if y != nil { - // y = (z - a*x)/b - y.Mul(a, lastX) - y.Sub(A, y) - y.Div(y, b) + if extended { + // Ua, Ub = Ub, Ua - q*Ub + t.Set(Ub) + s.Mul(Ub, q) + Ub.Sub(Ua, s) + Ua.Set(t) } - - *z = *A - return z } // lehmerGCD sets z to the greatest common divisor of a and b, // which both must be > 0, and returns z. +// If x or y are not nil, their values are set such that z = a*x + b*y. // See Knuth, The Art of Computer Programming, Vol. 2, Section 4.5.2, Algorithm L. // This implementation uses the improved condition by Collins requiring only one // quotient and avoiding the possibility of single Word overflow. // See Jebelean, "Improving the multiprecision Euclidean algorithm", // Design and Implementation of Symbolic Computation Systems, pp 45-58. -func (z *Int) lehmerGCD(a, b *Int) *Int { - // ensure a >= b - if a.abs.cmp(b.abs) < 0 { - a, b = b, a - } +// The cosequences are updated according to Algorithm 10.45 from +// Cohen et al. "Handbook of Elliptic and Hyperelliptic Curve Cryptography" pp 192. +func (z *Int) lehmerGCD(x, y, a, b *Int) *Int { + var A, B, Ua, Ub *Int - // don't destroy incoming values of a and b - B := new(Int).Set(b) // must be set first in case b is an alias of z - A := z.Set(a) + A = new(Int).Set(a) + B = new(Int).Set(b) + + extended := x != nil || y != nil + + if extended { + // Ua (Ub) tracks how many times input a has been accumulated into A (B). + Ua = new(Int).SetInt64(1) + Ub = new(Int) + } // temp variables for multiprecision update - t := new(Int) + q := new(Int) r := new(Int) s := new(Int) - w := new(Int) + t := new(Int) + + // ensure A >= B + if A.abs.cmp(B.abs) < 0 { + A, B = B, A + Ub, Ua = Ua, Ub + } // loop invariant A >= B for len(B.abs) > 1 { - // initialize the digits - var a1, a2, u0, u1, u2, v0, v1, v2 Word - - m := len(B.abs) // m >= 2 - n := len(A.abs) // n >= m >= 2 - - // extract the top Word of bits from A and B - h := nlz(A.abs[n-1]) - a1 = (A.abs[n-1] << h) | (A.abs[n-2] >> (_W - h)) - // B may have implicit zero words in the high bits if the lengths differ - switch { - case n == m: - a2 = (B.abs[n-1] << h) | (B.abs[n-2] >> (_W - h)) - case n == m+1: - a2 = (B.abs[n-2] >> (_W - h)) - default: - a2 = 0 - } - - // Since we are calculating with full words to avoid overflow, - // we use 'even' to track the sign of the cosequences. - // For even iterations: u0, v1 >= 0 && u1, v0 <= 0 - // For odd iterations: u0, v1 <= 0 && u1, v0 >= 0 - // The first iteration starts with k=1 (odd). - even := false - // variables to track the cosequences - u0, u1, u2 = 0, 1, 0 - v0, v1, v2 = 0, 0, 1 - - // Calculate the quotient and cosequences using Collins' stopping condition. - // Note that overflow of a Word is not possible when computing the remainder - // sequence and cosequences since the cosequence size is bounded by the input size. - // See section 4.2 of Jebelean for details. - for a2 >= v2 && a1-a2 >= v1+v2 { - q := a1 / a2 - a1, a2 = a2, a1-q*a2 - u0, u1, u2 = u1, u2, u1+q*u2 - v0, v1, v2 = v1, v2, v1+q*v2 - even = !even - } + // Attempt to calculate in single-precision using leading words of A and B. + u0, u1, v0, v1, even := lehmerSimulate(A, B) - // multiprecision step + // multiprecision Step if v0 != 0 { - // simulate the effect of the single precision steps using the cosequences + // Simulate the effect of the single-precision steps using the cosequences. // A = u0*A + v0*B // B = u1*A + v1*B + lehmerUpdate(A, B, q, r, s, t, u0, u1, v0, v1, even) - t.abs = t.abs.setWord(u0) - s.abs = s.abs.setWord(v0) - t.neg = !even - s.neg = even - - t.Mul(A, t) - s.Mul(B, s) - - r.abs = r.abs.setWord(u1) - w.abs = w.abs.setWord(v1) - r.neg = even - w.neg = !even - - r.Mul(A, r) - w.Mul(B, w) - - A.Add(t, s) - B.Add(r, w) + if extended { + // Ua = u0*Ua + v0*Ub + // Ub = u1*Ua + v1*Ub + lehmerUpdate(Ua, Ub, q, r, s, t, u0, u1, v0, v1, even) + } } else { - // single-digit calculations failed to simluate any quotients - // do a standard Euclidean step - t.Rem(A, B) - A, B, t = B, t, A + // Single-digit calculations failed to simulate any quotients. + // Do a standard Euclidean step. + euclidUpdate(A, B, Ua, Ub, q, r, s, t, extended) } } if len(B.abs) > 0 { - // standard Euclidean algorithm base case for B a single Word + // extended Euclidean algorithm base case if B is a single Word if len(A.abs) > 1 { - // A is longer than a single Word - t.Rem(A, B) - A, B, t = B, t, A + // A is longer than a single Word, so one update is needed. + euclidUpdate(A, B, Ua, Ub, q, r, s, t, extended) } if len(B.abs) > 0 { - // A and B are both a single Word - a1, a2 := A.abs[0], B.abs[0] - for a2 != 0 { - a1, a2 = a2, a1%a2 + // A and B are both a single Word. + aWord, bWord := A.abs[0], B.abs[0] + if extended { + var ua, ub, va, vb Word + ua, ub = 1, 0 + va, vb = 0, 1 + even := true + for bWord != 0 { + q, r := aWord/bWord, aWord%bWord + aWord, bWord = bWord, r + ua, ub = ub, ua+q*ub + va, vb = vb, va+q*vb + even = !even + } + + t.abs = t.abs.setWord(ua) + s.abs = s.abs.setWord(va) + t.neg = !even + s.neg = even + + t.Mul(Ua, t) + s.Mul(Ub, s) + + Ua.Add(t, s) + } else { + for bWord != 0 { + aWord, bWord = bWord, aWord%bWord + } } - A.abs[0] = a1 + A.abs[0] = aWord } } + + if x != nil { + *x = *Ua + } + + if y != nil { + // y = (z - a*x)/b + y.Mul(a, Ua) + y.Sub(A, y) + y.Div(y, b) + } + *z = *A + return z } @@ -659,20 +724,33 @@ func (z *Int) Rand(rnd *rand.Rand, n *Int) *Int { } // ModInverse sets z to the multiplicative inverse of g in the ring ℤ/nℤ -// and returns z. If g and n are not relatively prime, the result is undefined. +// and returns z. If g and n are not relatively prime, g has no multiplicative +// inverse in the ring ℤ/nℤ. In this case, z is unchanged and the return value +// is nil. func (z *Int) ModInverse(g, n *Int) *Int { + // GCD expects parameters a and b to be > 0. + if n.neg { + var n2 Int + n = n2.Neg(n) + } if g.neg { - // GCD expects parameters a and b to be > 0. var g2 Int g = g2.Mod(g, n) } - var d Int - d.GCD(z, nil, g, n) - // x and y are such that g*x + n*y = d. Since g and n are - // relatively prime, d = 1. Taking that modulo n results in - // g*x = 1, therefore x is the inverse element. - if z.neg { - z.Add(z, n) + var d, x Int + d.GCD(&x, nil, g, n) + + // if and only if d==1, g and n are relatively prime + if d.Cmp(intOne) != 0 { + return nil + } + + // x and y are such that g*x + n*y = 1, therefore x is the inverse element, + // but it may be negative, so convert to the range 0 <= z < |n| + if x.neg { + z.Add(&x, n) + } else { + z.Set(&x) } return z } @@ -745,6 +823,30 @@ func (z *Int) modSqrt3Mod4Prime(x, p *Int) *Int { return z } +// modSqrt5Mod8 uses Atkin's observation that 2 is not a square mod p +// alpha == (2*a)^((p-5)/8) mod p +// beta == 2*a*alpha^2 mod p is a square root of -1 +// b == a*alpha*(beta-1) mod p is a square root of a +// to calculate the square root of any quadratic residue mod p quickly for 5 +// mod 8 primes. +func (z *Int) modSqrt5Mod8Prime(x, p *Int) *Int { + // p == 5 mod 8 implies p = e*8 + 5 + // e is the quotient and 5 the remainder on division by 8 + e := new(Int).Rsh(p, 3) // e = (p - 5) / 8 + tx := new(Int).Lsh(x, 1) // tx = 2*x + alpha := new(Int).Exp(tx, e, p) + beta := new(Int).Mul(alpha, alpha) + beta.Mod(beta, p) + beta.Mul(beta, tx) + beta.Mod(beta, p) + beta.Sub(beta, intOne) + beta.Mul(beta, x) + beta.Mod(beta, p) + beta.Mul(beta, alpha) + z.Mod(beta, p) + return z +} + // modSqrtTonelliShanks uses the Tonelli-Shanks algorithm to find the square // root of a quadratic residue modulo any prime. func (z *Int) modSqrtTonelliShanks(x, p *Int) *Int { @@ -811,12 +913,17 @@ func (z *Int) ModSqrt(x, p *Int) *Int { x = new(Int).Mod(x, p) } - // Check whether p is 3 mod 4, and if so, use the faster algorithm. - if len(p.abs) > 0 && p.abs[0]%4 == 3 { + switch { + case p.abs[0]%4 == 3: + // Check whether p is 3 mod 4, and if so, use the faster algorithm. return z.modSqrt3Mod4Prime(x, p) + case p.abs[0]%8 == 5: + // Check whether p is 5 mod 8, use Atkin's algorithm. + return z.modSqrt5Mod8Prime(x, p) + default: + // Otherwise, use Tonelli-Shanks. + return z.modSqrtTonelliShanks(x, p) } - // Otherwise, use Tonelli-Shanks. - return z.modSqrtTonelliShanks(x, p) } // Lsh sets z = x << n and returns z. diff --git a/libgo/go/math/big/int_test.go b/libgo/go/math/big/int_test.go index 270fec6..9930ed0 100644 --- a/libgo/go/math/big/int_test.go +++ b/libgo/go/math/big/int_test.go @@ -557,7 +557,7 @@ var expTests = []struct { {"0x8000000000000000", "3", "6719", "5447"}, {"0x8000000000000000", "1000", "6719", "1603"}, {"0x8000000000000000", "1000000", "6719", "3199"}, - {"0x8000000000000000", "-1000000", "6719", "1"}, + {"0x8000000000000000", "-1000000", "6719", "3663"}, // 3663 = ModInverse(3199, 6719) Issue #25865 {"0xffffffffffffffffffffffffffffffff", "0x12345678123456781234567812345678123456789", "0x01112222333344445555666677778889", "0x36168FA1DB3AAE6C8CE647E137F97A"}, @@ -675,21 +675,43 @@ func checkGcd(aBytes, bBytes []byte) bool { return x.Cmp(d) == 0 } -// euclidGCD is a reference implementation of Euclid's GCD -// algorithm for testing against optimized algorithms. +// euclidExtGCD is a reference implementation of Euclid's +// extended GCD algorithm for testing against optimized algorithms. // Requirements: a, b > 0 -func euclidGCD(a, b *Int) *Int { - +func euclidExtGCD(a, b *Int) (g, x, y *Int) { A := new(Int).Set(a) B := new(Int).Set(b) - t := new(Int) + // A = Ua*a + Va*b + // B = Ub*a + Vb*b + Ua := new(Int).SetInt64(1) + Va := new(Int) + + Ub := new(Int) + Vb := new(Int).SetInt64(1) + + q := new(Int) + temp := new(Int) + + r := new(Int) for len(B.abs) > 0 { - // A, B = B, A % B - t.Rem(A, B) - A, B, t = B, t, A + q, r = q.QuoRem(A, B, r) + + A, B, r = B, r, A + + // Ua, Ub = Ub, Ua-q*Ub + temp.Set(Ub) + Ub.Mul(Ub, q) + Ub.Sub(Ua, Ub) + Ua.Set(temp) + + // Va, Vb = Vb, Va-q*Vb + temp.Set(Vb) + Vb.Mul(Vb, q) + Vb.Sub(Va, Vb) + Va.Set(temp) } - return A + return A, Ua, Va } func checkLehmerGcd(aBytes, bBytes []byte) bool { @@ -700,12 +722,28 @@ func checkLehmerGcd(aBytes, bBytes []byte) bool { return true // can only test positive arguments } - d := new(Int).lehmerGCD(a, b) - d0 := euclidGCD(a, b) + d := new(Int).lehmerGCD(nil, nil, a, b) + d0, _, _ := euclidExtGCD(a, b) return d.Cmp(d0) == 0 } +func checkLehmerExtGcd(aBytes, bBytes []byte) bool { + a := new(Int).SetBytes(aBytes) + b := new(Int).SetBytes(bBytes) + x := new(Int) + y := new(Int) + + if a.Sign() <= 0 || b.Sign() <= 0 { + return true // can only test positive arguments + } + + d := new(Int).lehmerGCD(x, y, a, b) + d0, x0, y0 := euclidExtGCD(a, b) + + return d.Cmp(d0) == 0 && x.Cmp(x0) == 0 && y.Cmp(y0) == 0 +} + var gcdTests = []struct { d, x, y, a, b string }{ @@ -785,6 +823,10 @@ func TestGcd(t *testing.T) { if err := quick.Check(checkLehmerGcd, nil); err != nil { t.Error(err) } + + if err := quick.Check(checkLehmerExtGcd, nil); err != nil { + t.Error(err) + } } type intShiftTest struct { @@ -1318,7 +1360,7 @@ func BenchmarkModSqrt225_Tonelli(b *testing.B) { } } -func BenchmarkModSqrt224_3Mod4(b *testing.B) { +func BenchmarkModSqrt225_3Mod4(b *testing.B) { p := tri(225) x := new(Int).SetUint64(2) for i := 0; i < b.N; i++ { @@ -1327,27 +1369,25 @@ func BenchmarkModSqrt224_3Mod4(b *testing.B) { } } -func BenchmarkModSqrt5430_Tonelli(b *testing.B) { - if isRaceBuilder { - b.Skip("skipping on race builder") - } - p := tri(5430) - x := new(Int).SetUint64(2) +func BenchmarkModSqrt231_Tonelli(b *testing.B) { + p := tri(231) + p.Sub(p, intOne) + p.Sub(p, intOne) // tri(231) - 2 is a prime == 5 mod 8 + x := new(Int).SetUint64(7) for i := 0; i < b.N; i++ { - x.SetUint64(2) + x.SetUint64(7) x.modSqrtTonelliShanks(x, p) } } -func BenchmarkModSqrt5430_3Mod4(b *testing.B) { - if isRaceBuilder { - b.Skip("skipping on race builder") - } - p := tri(5430) - x := new(Int).SetUint64(2) +func BenchmarkModSqrt231_5Mod8(b *testing.B) { + p := tri(231) + p.Sub(p, intOne) + p.Sub(p, intOne) // tri(231) - 2 is a prime == 5 mod 8 + x := new(Int).SetUint64(7) for i := 0; i < b.N; i++ { - x.SetUint64(2) - x.modSqrt3Mod4Prime(x, p) + x.SetUint64(7) + x.modSqrt5Mod8Prime(x, p) } } @@ -1409,19 +1449,21 @@ var modInverseTests = []struct { {"1234567", "458948883992"}, {"239487239847", "2410312426921032588552076022197566074856950548502459942654116941958108831682612228890093858261341614673227141477904012196503648957050582631942730706805009223062734745341073406696246014589361659774041027169249453200378729434170325843778659198143763193776859869524088940195577346119843545301547043747207749969763750084308926339295559968882457872412993810129130294592999947926365264059284647209730384947211681434464714438488520940127459844288859336526896320919633919"}, {"-10", "13"}, // issue #16984 + {"10", "-13"}, + {"-17", "-13"}, } func TestModInverse(t *testing.T) { var element, modulus, gcd, inverse Int one := NewInt(1) - for i, test := range modInverseTests { + for _, test := range modInverseTests { (&element).SetString(test.element, 10) (&modulus).SetString(test.modulus, 10) (&inverse).ModInverse(&element, &modulus) (&inverse).Mul(&inverse, &element) (&inverse).Mod(&inverse, &modulus) if (&inverse).Cmp(one) != 0 { - t.Errorf("#%d: failed (e·e^(-1)=%s)", i, &inverse) + t.Errorf("ModInverse(%d,%d)*%d%%%d=%d, not 1", &element, &modulus, &element, &modulus, &inverse) } } // exhaustive test for small values @@ -1443,6 +1485,17 @@ func TestModInverse(t *testing.T) { } } +func BenchmarkModInverse(b *testing.B) { + p := new(Int).SetInt64(1) // Mersenne prime 2**1279 -1 + p.abs = p.abs.shl(p.abs, 1279) + p.Sub(p, intOne) + x := new(Int).Sub(p, intOne) + z := new(Int) + for i := 0; i < b.N; i++ { + z.ModInverse(x, p) + } +} + // testModSqrt is a helper for TestModSqrt, // which checks that ModSqrt can compute a square-root of elt^2. func testModSqrt(t *testing.T, elt, mod, sq, sqrt *Int) bool { diff --git a/libgo/go/math/big/nat.go b/libgo/go/math/big/nat.go index 3bb818f..a6f79ed 100644 --- a/libgo/go/math/big/nat.go +++ b/libgo/go/math/big/nat.go @@ -5,10 +5,16 @@ // This file implements unsigned multi-precision integers (natural // numbers). They are the building blocks for the implementation // of signed integers, rationals, and floating-point numbers. +// +// Caution: This implementation relies on the function "alias" +// which assumes that (nat) slice capacities are never +// changed (no 3-operand slice expressions). If that +// changes, alias needs to be updated for correctness. package big import ( + "encoding/binary" "math/bits" "math/rand" "sync" @@ -207,18 +213,17 @@ func (z nat) montgomery(x, y, m nat, k Word, n int) nat { if len(x) != n || len(y) != n || len(m) != n { panic("math/big: mismatched montgomery number lengths") } - z = z.make(n) + z = z.make(n * 2) z.clear() var c Word for i := 0; i < n; i++ { d := y[i] - c2 := addMulVVW(z, x, d) - t := z[0] * k - c3 := addMulVVW(z, m, t) - copy(z, z[1:]) + c2 := addMulVVW(z[i:n+i], x, d) + t := z[i] * k + c3 := addMulVVW(z[i:n+i], m, t) cx := c + c2 cy := cx + c3 - z[n-1] = cy + z[n+i] = cy if cx < c2 || cy < c3 { c = 1 } else { @@ -226,9 +231,11 @@ func (z nat) montgomery(x, y, m nat, k Word, n int) nat { } } if c != 0 { - subVV(z, z, m) + subVV(z[:n], z[n:], m) + } else { + copy(z[:n], z[n:]) } - return z + return z[:n] } // Fast version of z[0:n+n>>1].add(z[0:n+n>>1], x[0:n]) w/o bounds checks. @@ -351,6 +358,10 @@ func karatsuba(z, x, y nat) { } // alias reports whether x and y share the same base array. +// Note: alias assumes that the capacity of underlying arrays +// is never changed for nat values; i.e. that there are +// no 3-operand slice expressions in this code (or worse, +// reflect-based operations to the same effect). func alias(x, y nat) bool { return cap(x) > 0 && cap(y) > 0 && &x[0:cap(x)][cap(x)-1] == &y[0:cap(y)][cap(y)-1] } @@ -377,12 +388,12 @@ func max(x, y int) int { } // karatsubaLen computes an approximation to the maximum k <= n such that -// k = p<<i for a number p <= karatsubaThreshold and an i >= 0. Thus, the +// k = p<<i for a number p <= threshold and an i >= 0. Thus, the // result is the largest number that can be divided repeatedly by 2 before -// becoming about the value of karatsubaThreshold. -func karatsubaLen(n int) int { +// becoming about the value of threshold. +func karatsubaLen(n, threshold int) int { i := uint(0) - for n > karatsubaThreshold { + for n > threshold { n >>= 1 i++ } @@ -422,7 +433,7 @@ func (z nat) mul(x, y nat) nat { // y = yh*b + y0 (0 <= y0 < b) // b = 1<<(_W*k) ("base" of digits xi, yi) // - k := karatsubaLen(n) + k := karatsubaLen(n, karatsubaThreshold) // k <= n // multiply x0 and y0 via Karatsuba @@ -475,8 +486,8 @@ func (z nat) mul(x, y nat) nat { // basicSqr sets z = x*x and is asymptotically faster than basicMul // by about a factor of 2, but slower for small arguments due to overhead. -// Requirements: len(x) > 0, len(z) >= 2*len(x) -// The (non-normalized) result is placed in z[0 : 2 * len(x)]. +// Requirements: len(x) > 0, len(z) == 2*len(x) +// The (non-normalized) result is placed in z. func basicSqr(z, x nat) { n := len(x) t := make(nat, 2*n) // temporary variable to hold the products @@ -492,11 +503,47 @@ func basicSqr(z, x nat) { addVV(z, z, t) // combine the result } +// karatsubaSqr squares x and leaves the result in z. +// len(x) must be a power of 2 and len(z) >= 6*len(x). +// The (non-normalized) result is placed in z[0 : 2*len(x)]. +// +// The algorithm and the layout of z are the same as for karatsuba. +func karatsubaSqr(z, x nat) { + n := len(x) + + if n&1 != 0 || n < karatsubaSqrThreshold || n < 2 { + basicSqr(z[:2*n], x) + return + } + + n2 := n >> 1 + x1, x0 := x[n2:], x[0:n2] + + karatsubaSqr(z, x0) + karatsubaSqr(z[n:], x1) + + // s = sign(xd*yd) == -1 for xd != 0; s == 1 for xd == 0 + xd := z[2*n : 2*n+n2] + if subVV(xd, x1, x0) != 0 { + subVV(xd, x0, x1) + } + + p := z[n*3:] + karatsubaSqr(p, xd) + + r := z[n*4:] + copy(r, z[:n*2]) + + karatsubaAdd(z[n2:], r, n) + karatsubaAdd(z[n2:], r[n:], n) + karatsubaSub(z[n2:], p, n) // s == -1 for p != 0; s == 1 for p == 0 +} + // Operands that are shorter than basicSqrThreshold are squared using // "grade school" multiplication; for operands longer than karatsubaSqrThreshold -// the Karatsuba algorithm is used. +// we use the Karatsuba algorithm optimized for x == y. var basicSqrThreshold = 20 // computed by calibrate_test.go -var karatsubaSqrThreshold = 400 // computed by calibrate_test.go +var karatsubaSqrThreshold = 260 // computed by calibrate_test.go // z = x*x func (z nat) sqr(x nat) nat { @@ -514,18 +561,43 @@ func (z nat) sqr(x nat) nat { if alias(z, x) { z = nil // z is an alias for x - cannot reuse } - z = z.make(2 * n) if n < basicSqrThreshold { + z = z.make(2 * n) basicMul(z, x, x) return z.norm() } if n < karatsubaSqrThreshold { + z = z.make(2 * n) basicSqr(z, x) return z.norm() } - return z.mul(x, x) + // Use Karatsuba multiplication optimized for x == y. + // The algorithm and layout of z are the same as for mul. + + // z = (x1*b + x0)^2 = x1^2*b^2 + 2*x1*x0*b + x0^2 + + k := karatsubaLen(n, karatsubaSqrThreshold) + + x0 := x[0:k] + z = z.make(max(6*k, 2*n)) + karatsubaSqr(z, x0) // z = x0^2 + z = z[0 : 2*n] + z[2*k:].clear() + + if k < n { + var t nat + x0 := x0.norm() + x1 := x[k:] + t = t.mul(x0, x1) + addAt(z, t, k) + addAt(z, t, k) // z = 2*x1*x0*b + x0^2 + t = t.sqr(x1) + addAt(z, t, 2*k) // z = x1^2*b^2 + 2*x1*x0*b + x0^2 + } + + return z.norm() } // mulRange computes the product of all the unsigned integers in the @@ -718,8 +790,21 @@ func (x nat) trailingZeroBits() uint { return i*_W + uint(bits.TrailingZeros(uint(x[i]))) } +func same(x, y nat) bool { + return len(x) == len(y) && len(x) > 0 && &x[0] == &y[0] +} + // z = x << s func (z nat) shl(x nat, s uint) nat { + if s == 0 { + if same(z, x) { + return z + } + if !alias(z, x) { + return z.set(x) + } + } + m := len(x) if m == 0 { return z[:0] @@ -736,6 +821,15 @@ func (z nat) shl(x nat, s uint) nat { // z = x >> s func (z nat) shr(x nat, s uint) nat { + if s == 0 { + if same(z, x) { + return z + } + if !alias(z, x) { + return z.set(x) + } + } + m := len(x) n := m - int(s/_W) if n <= 0 { @@ -952,7 +1046,7 @@ func (z nat) expNN(x, y, m nat) nat { // x**1 mod m == x mod m if len(y) == 1 && y[0] == 1 && len(m) != 0 { - _, z = z.div(z, x, m) + _, z = nat(nil).div(z, x, m) return z } // y > 1 @@ -1125,7 +1219,7 @@ func (z nat) expNNMontgomery(x, y, m nat) nat { // RR = 2**(2*_W*len(m)) mod m RR := nat(nil).setWord(1) zz := nat(nil).shl(RR, uint(2*numWords*_W)) - _, RR = RR.div(RR, zz, m) + _, RR = nat(nil).div(RR, zz, m) if len(RR) < numWords { zz = zz.make(numWords) copy(zz, RR) @@ -1208,25 +1302,31 @@ func (z nat) bytes(buf []byte) (i int) { return } +// bigEndianWord returns the contents of buf interpreted as a big-endian encoded Word value. +func bigEndianWord(buf []byte) Word { + if _W == 64 { + return Word(binary.BigEndian.Uint64(buf)) + } + return Word(binary.BigEndian.Uint32(buf)) +} + // setBytes interprets buf as the bytes of a big-endian unsigned // integer, sets z to that value, and returns z. func (z nat) setBytes(buf []byte) nat { z = z.make((len(buf) + _S - 1) / _S) - k := 0 - s := uint(0) - var d Word - for i := len(buf); i > 0; i-- { - d |= Word(buf[i-1]) << s - if s += 8; s == _S*8 { - z[k] = d - k++ - s = 0 - d = 0 - } + i := len(buf) + for k := 0; i >= _S; k++ { + z[k] = bigEndianWord(buf[i-_S : i]) + i -= _S } - if k < len(z) { - z[k] = d + if i > 0 { + var d Word + for s := uint(0); i > 0; s += 8 { + d |= Word(buf[i-1]) << s + i-- + } + z[len(z)-1] = d } return z.norm() diff --git a/libgo/go/math/big/nat_test.go b/libgo/go/math/big/nat_test.go index c25cdf0..3c79495 100644 --- a/libgo/go/math/big/nat_test.go +++ b/libgo/go/math/big/nat_test.go @@ -267,6 +267,34 @@ func TestShiftRight(t *testing.T) { } } +func BenchmarkZeroShifts(b *testing.B) { + x := rndNat(800) + + b.Run("Shl", func(b *testing.B) { + for i := 0; i < b.N; i++ { + var z nat + z.shl(x, 0) + } + }) + b.Run("ShlSame", func(b *testing.B) { + for i := 0; i < b.N; i++ { + x.shl(x, 0) + } + }) + + b.Run("Shr", func(b *testing.B) { + for i := 0; i < b.N; i++ { + var z nat + z.shr(x, 0) + } + }) + b.Run("ShrSame", func(b *testing.B) { + for i := 0; i < b.N; i++ { + x.shr(x, 0) + } + }) +} + type modWTest struct { in string dividend string @@ -620,26 +648,26 @@ func TestSticky(t *testing.T) { } } -func testBasicSqr(t *testing.T, x nat) { +func testSqr(t *testing.T, x nat) { got := make(nat, 2*len(x)) want := make(nat, 2*len(x)) - basicSqr(got, x) - basicMul(want, x, x) + got = got.sqr(x) + want = want.mul(x, x) if got.cmp(want) != 0 { t.Errorf("basicSqr(%v), got %v, want %v", x, got, want) } } -func TestBasicSqr(t *testing.T) { +func TestSqr(t *testing.T) { for _, a := range prodNN { if a.x != nil { - testBasicSqr(t, a.x) + testSqr(t, a.x) } if a.y != nil { - testBasicSqr(t, a.y) + testSqr(t, a.y) } if a.z != nil { - testBasicSqr(t, a.z) + testSqr(t, a.z) } } } @@ -665,3 +693,22 @@ func BenchmarkNatSqr(b *testing.B) { }) } } + +func BenchmarkNatSetBytes(b *testing.B) { + const maxLength = 128 + lengths := []int{ + // No remainder: + 8, 24, maxLength, + // With remainder: + 7, 23, maxLength - 1, + } + n := make(nat, maxLength/_W) // ensure n doesn't need to grow during the test + buf := make([]byte, maxLength) + for _, l := range lengths { + b.Run(fmt.Sprint(l), func(b *testing.B) { + for i := 0; i < b.N; i++ { + n.setBytes(buf[:l]) + } + }) + } +} diff --git a/libgo/go/math/big/prime.go b/libgo/go/math/big/prime.go index 848affb..4c2c152 100644 --- a/libgo/go/math/big/prime.go +++ b/libgo/go/math/big/prime.go @@ -131,11 +131,11 @@ NextRandom: // // Baillie and Wagstaff, "Lucas Pseudoprimes", Mathematics of Computation 35(152), // October 1980, pp. 1391-1417, especially page 1401. -// http://www.ams.org/journals/mcom/1980-35-152/S0025-5718-1980-0583518-6/S0025-5718-1980-0583518-6.pdf +// https://www.ams.org/journals/mcom/1980-35-152/S0025-5718-1980-0583518-6/S0025-5718-1980-0583518-6.pdf // // Grantham, "Frobenius Pseudoprimes", Mathematics of Computation 70(234), // March 2000, pp. 873-891. -// http://www.ams.org/journals/mcom/2001-70-234/S0025-5718-00-01197-2/S0025-5718-00-01197-2.pdf +// https://www.ams.org/journals/mcom/2001-70-234/S0025-5718-00-01197-2/S0025-5718-00-01197-2.pdf // // Baillie, "Extra strong Lucas pseudoprimes", OEIS A217719, https://oeis.org/A217719. // diff --git a/libgo/go/math/big/prime_test.go b/libgo/go/math/big/prime_test.go index 7760519..bf50f34 100644 --- a/libgo/go/math/big/prime_test.go +++ b/libgo/go/math/big/prime_test.go @@ -29,13 +29,13 @@ var primes = []string{ "98920366548084643601728869055592650835572950932266967461790948584315647051443", "94560208308847015747498523884063394671606671904944666360068158221458669711639", - // http://primes.utm.edu/lists/small/small3.html + // https://primes.utm.edu/lists/small/small3.html "449417999055441493994709297093108513015373787049558499205492347871729927573118262811508386655998299074566974373711472560655026288668094291699357843464363003144674940345912431129144354948751003607115263071543163", "230975859993204150666423538988557839555560243929065415434980904258310530753006723857139742334640122533598517597674807096648905501653461687601339782814316124971547968912893214002992086353183070342498989426570593", "5521712099665906221540423207019333379125265462121169655563495403888449493493629943498064604536961775110765377745550377067893607246020694972959780839151452457728855382113555867743022746090187341871655890805971735385789993", "203956878356401977405765866929034577280193993314348263094772646453283062722701277632936616063144088173312372882677123879538709400158306567338328279154499698366071906766440037074217117805690872792848149112022286332144876183376326512083574821647933992961249917319836219304274280243803104015000563790123", - // ECC primes: http://tools.ietf.org/html/draft-ladd-safecurves-02 + // ECC primes: https://tools.ietf.org/html/draft-ladd-safecurves-02 "3618502788666131106986593281521497120414687020801267626233049500247285301239", // Curve1174: 2^251-9 "57896044618658097711785492504343953926634992332820282019728792003956564819949", // Curve25519: 2^255-19 "9850501549098619803069760025035903451269934817616361666987073351061430442874302652853566563721228910201656997576599", // E-382: 2^382-105 diff --git a/libgo/go/math/big/rat.go b/libgo/go/math/big/rat.go index b33fc69..46d58fc 100644 --- a/libgo/go/math/big/rat.go +++ b/libgo/go/math/big/rat.go @@ -422,7 +422,7 @@ func (z *Rat) norm() *Rat { neg := z.a.neg z.a.neg = false z.b.neg = false - if f := NewInt(0).lehmerGCD(&z.a, &z.b); f.Cmp(intOne) != 0 { + if f := NewInt(0).lehmerGCD(nil, nil, &z.a, &z.b); f.Cmp(intOne) != 0 { z.a.abs, _ = z.a.abs.div(nil, z.a.abs, f.abs) z.b.abs, _ = z.b.abs.div(nil, z.b.abs, f.abs) if z.b.abs.cmp(natOne) == 0 { diff --git a/libgo/go/math/big/ratconv_test.go b/libgo/go/math/big/ratconv_test.go index 56ac8d7..fe8b8b6 100644 --- a/libgo/go/math/big/ratconv_test.go +++ b/libgo/go/math/big/ratconv_test.go @@ -300,9 +300,9 @@ var float64inputs = []string{ // "1e-18446744073709551616", // "1e+18446744073709551616", - // http://www.exploringbinary.com/java-hangs-when-converting-2-2250738585072012e-308/ + // https://www.exploringbinary.com/java-hangs-when-converting-2-2250738585072012e-308/ "2.2250738585072012e-308", - // http://www.exploringbinary.com/php-hangs-on-numeric-value-2-2250738585072011e-308/ + // https://www.exploringbinary.com/php-hangs-on-numeric-value-2-2250738585072011e-308/ "2.2250738585072011e-308", // A very large number (initially wrongly parsed by the fast algorithm). diff --git a/libgo/go/math/big/roundingmode_string.go b/libgo/go/math/big/roundingmode_string.go index 05024b8..c7629eb 100644 --- a/libgo/go/math/big/roundingmode_string.go +++ b/libgo/go/math/big/roundingmode_string.go @@ -1,16 +1,16 @@ -// generated by stringer -type=RoundingMode; DO NOT EDIT +// Code generated by "stringer -type=RoundingMode"; DO NOT EDIT. package big -import "fmt" +import "strconv" const _RoundingMode_name = "ToNearestEvenToNearestAwayToZeroAwayFromZeroToNegativeInfToPositiveInf" var _RoundingMode_index = [...]uint8{0, 13, 26, 32, 44, 57, 70} func (i RoundingMode) String() string { - if i+1 >= RoundingMode(len(_RoundingMode_index)) { - return fmt.Sprintf("RoundingMode(%d)", i) + if i >= RoundingMode(len(_RoundingMode_index)-1) { + return "RoundingMode(" + strconv.FormatInt(int64(i), 10) + ")" } return _RoundingMode_name[_RoundingMode_index[i]:_RoundingMode_index[i+1]] } diff --git a/libgo/go/math/big/sqrt.go b/libgo/go/math/big/sqrt.go index 00433cf..b989649 100644 --- a/libgo/go/math/big/sqrt.go +++ b/libgo/go/math/big/sqrt.go @@ -128,18 +128,21 @@ func (z *Float) sqrtInverse(x *Float) { // g(t) = f(t)/f'(t) = -½t(1 - xt²) // and the next guess is given by // t2 = t - g(t) = ½t(3 - xt²) - u := new(Float) + u := newFloat(z.prec) + v := newFloat(z.prec) ng := func(t *Float) *Float { u.prec = t.prec + v.prec = t.prec u.Mul(t, t) // u = t² u.Mul(x, u) // = xt² - u.Sub(three, u) // = 3 - xt² - u.Mul(t, u) // = t(3 - xt²) + v.Sub(three, u) // v = 3 - xt² + u.Mul(t, v) // u = t(3 - xt²) return t.Mul(half, u) // = ½t(3 - xt²) } xf, _ := x.Float64() - sqi := NewFloat(1 / math.Sqrt(xf)) + sqi := newFloat(z.prec) + sqi.SetFloat64(1 / math.Sqrt(xf)) for prec := z.prec + 32; sqi.prec < prec; { sqi.prec *= 2 sqi = ng(sqi) @@ -149,3 +152,12 @@ func (z *Float) sqrtInverse(x *Float) { // x/√x = √x z.Mul(x, sqi) } + +// newFloat returns a new *Float with space for twice the given +// precision. +func newFloat(prec2 uint32) *Float { + z := new(Float) + // nat.make ensures the slice length is > 0 + z.mant = z.mant.make(int(prec2/_W) * 2) + return z +} diff --git a/libgo/go/math/bits/bits_test.go b/libgo/go/math/bits/bits_test.go index ba05210..5c34f6d 100644 --- a/libgo/go/math/bits/bits_test.go +++ b/libgo/go/math/bits/bits_test.go @@ -2,9 +2,10 @@ // Use of this source code is governed by a BSD-style // license that can be found in the LICENSE file. -package bits +package bits_test import ( + . "math/bits" "testing" "unsafe" ) @@ -83,7 +84,7 @@ func TestLeadingZeros(t *testing.T) { // Exported (global) variable serving as input for some // of the benchmarks to ensure side-effect free calls // are not optimized away. -var Input uint64 = deBruijn64 +var Input uint64 = DeBruijn64 // Exported (global) variable to store function results // during benchmarking to ensure side-effect free calls @@ -333,7 +334,7 @@ func BenchmarkOnesCount64(b *testing.B) { } func TestRotateLeft(t *testing.T) { - var m uint64 = deBruijn64 + var m uint64 = DeBruijn64 for k := uint(0); k < 128; k++ { x8 := uint8(m) diff --git a/libgo/go/math/bits/export_test.go b/libgo/go/math/bits/export_test.go new file mode 100644 index 0000000..8c6f933 --- /dev/null +++ b/libgo/go/math/bits/export_test.go @@ -0,0 +1,7 @@ +// Copyright 2018 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 bits + +const DeBruijn64 = deBruijn64 diff --git a/libgo/go/math/cmplx/cmath_test.go b/libgo/go/math/cmplx/cmath_test.go index 8d70562..80c3b33 100644 --- a/libgo/go/math/cmplx/cmath_test.go +++ b/libgo/go/math/cmplx/cmath_test.go @@ -39,7 +39,7 @@ var vc = []complex128{ } // The expected results below were computed by the high precision calculators -// at http://keisan.casio.com/. More exact input values (array vc[], above) +// at https://keisan.casio.com/. More exact input values (array vc[], above) // were obtained by printing them with "%.26f". The answers were calculated // to 26 digits (by using the "Digit number" drop-down control of each // calculator). diff --git a/libgo/go/math/erfinv.go b/libgo/go/math/erfinv.go index 21b5578..ee423d3 100644 --- a/libgo/go/math/erfinv.go +++ b/libgo/go/math/erfinv.go @@ -10,7 +10,7 @@ package math // This implementation is based on the rational approximation // of percentage points of normal distribution available from -// http://www.jstor.org/stable/2347330. +// https://www.jstor.org/stable/2347330. const ( // Coefficients for approximation to erf in |x| <= 0.85 diff --git a/libgo/go/math/example_test.go b/libgo/go/math/example_test.go index feaf9d8..a1f764b 100644 --- a/libgo/go/math/example_test.go +++ b/libgo/go/math/example_test.go @@ -89,3 +89,27 @@ func ExampleSqrt() { fmt.Printf("%.1f", c) // Output: 5.0 } + +func ExampleCeil() { + c := math.Ceil(1.49) + fmt.Printf("%.1f", c) + // Output: 2.0 +} + +func ExampleFloor() { + c := math.Floor(1.51) + fmt.Printf("%.1f", c) + // Output: 1.0 +} + +func ExamplePow() { + c := math.Pow(2, 3) + fmt.Printf("%.1f", c) + // Output: 8.0 +} + +func ExamplePow10() { + c := math.Pow10(2) + fmt.Printf("%.1f", c) + // Output: 100.0 +} diff --git a/libgo/go/math/floor_asm.go b/libgo/go/math/floor_asm.go deleted file mode 100644 index 761349a..0000000 --- a/libgo/go/math/floor_asm.go +++ /dev/null @@ -1,12 +0,0 @@ -// Copyright 2015 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. - -// +build ignore -// -build amd64 amd64p32 - -package math - -import "internal/cpu" - -var useSSE41 = cpu.X86.HasSSE41 diff --git a/libgo/go/math/hypot.go b/libgo/go/math/hypot.go index 59b9c74..844159a 100644 --- a/libgo/go/math/hypot.go +++ b/libgo/go/math/hypot.go @@ -28,12 +28,7 @@ func hypot(p, q float64) float64 { case IsNaN(p) || IsNaN(q): return NaN() } - if p < 0 { - p = -p - } - if q < 0 { - q = -q - } + p, q = Abs(p), Abs(q) if p < q { p, q = q, p } diff --git a/libgo/go/math/j0.go b/libgo/go/math/j0.go index fe26791..5523fc3 100644 --- a/libgo/go/math/j0.go +++ b/libgo/go/math/j0.go @@ -99,9 +99,7 @@ func J0(x float64) float64 { return 1 } - if x < 0 { - x = -x - } + x = Abs(x) if x >= 2 { s, c := Sincos(x) ss := s - c diff --git a/libgo/go/math/ldexp.go b/libgo/go/math/ldexp.go index e91a090c..8775cdb 100644 --- a/libgo/go/math/ldexp.go +++ b/libgo/go/math/ldexp.go @@ -37,7 +37,7 @@ func ldexp(frac float64, exp int) float64 { exp += e x := Float64bits(frac) exp += int(x>>shift)&mask - bias - if exp < -1074 { + if exp < -1075 { return Copysign(0, frac) // underflow } if exp > 1023 { // overflow @@ -48,8 +48,8 @@ func ldexp(frac float64, exp int) float64 { } var m float64 = 1 if exp < -1022 { // denormal - exp += 52 - m = 1.0 / (1 << 52) // 2**-52 + exp += 53 + m = 1.0 / (1 << 53) // 2**-53 } x &^= mask << shift x |= uint64(exp+bias) << shift diff --git a/libgo/go/math/lgamma.go b/libgo/go/math/lgamma.go index 19ac3ffa..7af5871 100644 --- a/libgo/go/math/lgamma.go +++ b/libgo/go/math/lgamma.go @@ -103,7 +103,7 @@ var _lgamA = [...]float64{ 4.48640949618915160150e-05, // 0x3F07858E90A45837 } var _lgamR = [...]float64{ - 1.0, // placeholder + 1.0, // placeholder 1.39200533467621045958e+00, // 0x3FF645A762C4AB74 7.21935547567138069525e-01, // 0x3FE71A1893D3DCDC 1.71933865632803078993e-01, // 0x3FC601EDCCFBDF27 diff --git a/libgo/go/math/rand/exp.go b/libgo/go/math/rand/exp.go index 4bc110f..5a8d946 100644 --- a/libgo/go/math/rand/exp.go +++ b/libgo/go/math/rand/exp.go @@ -13,7 +13,7 @@ import ( * * See "The Ziggurat Method for Generating Random Variables" * (Marsaglia & Tsang, 2000) - * http://www.jstatsoft.org/v05/i08/paper [pdf] + * https://www.jstatsoft.org/v05/i08/paper [pdf] */ const ( diff --git a/libgo/go/math/rand/normal.go b/libgo/go/math/rand/normal.go index ba4ea54..2c5a7aa 100644 --- a/libgo/go/math/rand/normal.go +++ b/libgo/go/math/rand/normal.go @@ -27,9 +27,9 @@ func absInt32(i int32) uint32 { return uint32(i) } -// NormFloat64 returns a normally distributed float64 in the range -// [-math.MaxFloat64, +math.MaxFloat64] with -// standard normal distribution (mean = 0, stddev = 1). +// NormFloat64 returns a normally distributed float64 in +// the range -math.MaxFloat64 through +math.MaxFloat64 inclusive, +// with standard normal distribution (mean = 0, stddev = 1). // To produce a different normal distribution, callers can // adjust the output using: // diff --git a/libgo/go/math/rand/rand.go b/libgo/go/math/rand/rand.go index 147c92f..04382e6 100644 --- a/libgo/go/math/rand/rand.go +++ b/libgo/go/math/rand/rand.go @@ -11,6 +11,9 @@ // The default Source is safe for concurrent use by multiple goroutines, but // Sources created by NewSource are not. // +// Mathematical interval notation such as [0, n) is used throughout the +// documentation for this package. +// // For random numbers suitable for security-sensitive work, see the crypto/rand // package. package rand diff --git a/libgo/go/math/rand/regress_test.go b/libgo/go/math/rand/regress_test.go index e31e6c5..1f30be8 100644 --- a/libgo/go/math/rand/regress_test.go +++ b/libgo/go/math/rand/regress_test.go @@ -121,219 +121,219 @@ func TestRegress(t *testing.T) { } var regressGolden = []interface{}{ - float64(4.668112973579268), // ExpFloat64() - float64(0.1601593871172866), // ExpFloat64() - float64(3.0465834105636), // ExpFloat64() - float64(0.06385839451671879), // ExpFloat64() - float64(1.8578917487258961), // ExpFloat64() - float64(0.784676123472182), // ExpFloat64() - float64(0.11225477361256932), // ExpFloat64() - float64(0.20173283329802255), // ExpFloat64() - float64(0.3468619496201105), // ExpFloat64() - float64(0.35601103454384536), // ExpFloat64() - float64(0.888376329507869), // ExpFloat64() - float64(1.4081362450365698), // ExpFloat64() - float64(1.0077753823151994), // ExpFloat64() - float64(0.23594100766227588), // ExpFloat64() - float64(2.777245612300007), // ExpFloat64() - float64(0.5202997830662377), // ExpFloat64() - float64(1.2842705247770294), // ExpFloat64() - float64(0.030307408362776206), // ExpFloat64() - float64(2.204156824853721), // ExpFloat64() - float64(2.09891923895058), // ExpFloat64() - float32(0.94519615), // Float32() - float32(0.24496509), // Float32() - float32(0.65595627), // Float32() - float32(0.05434384), // Float32() - float32(0.3675872), // Float32() - float32(0.28948045), // Float32() - float32(0.1924386), // Float32() - float32(0.65533215), // Float32() - float32(0.8971697), // Float32() - float32(0.16735445), // Float32() - float32(0.28858566), // Float32() - float32(0.9026048), // Float32() - float32(0.84978026), // Float32() - float32(0.2730468), // Float32() - float32(0.6090802), // Float32() - float32(0.253656), // Float32() - float32(0.7746542), // Float32() - float32(0.017480763), // Float32() - float32(0.78707397), // Float32() - float32(0.7993937), // Float32() - float64(0.9451961492941164), // Float64() - float64(0.24496508529377975), // Float64() - float64(0.6559562651954052), // Float64() - float64(0.05434383959970039), // Float64() - float64(0.36758720663245853), // Float64() - float64(0.2894804331565928), // Float64() - float64(0.19243860967493215), // Float64() - float64(0.6553321508148324), // Float64() - float64(0.897169713149801), // Float64() - float64(0.16735444255905835), // Float64() - float64(0.2885856518054551), // Float64() - float64(0.9026048462705047), // Float64() - float64(0.8497802817628735), // Float64() - float64(0.2730468047134829), // Float64() - float64(0.6090801919903561), // Float64() - float64(0.25365600644283687), // Float64() - float64(0.7746542391859803), // Float64() - float64(0.017480762156647272), // Float64() - float64(0.7870739563039942), // Float64() - float64(0.7993936979594545), // Float64() - int64(8717895732742165505), // Int() - int64(2259404117704393152), // Int() - int64(6050128673802995827), // Int() - int64(501233450539197794), // Int() - int64(3390393562759376202), // Int() - int64(2669985732393126063), // Int() - int64(1774932891286980153), // Int() - int64(6044372234677422456), // Int() - int64(8274930044578894929), // Int() - int64(1543572285742637646), // Int() - int64(2661732831099943416), // Int() - int64(8325060299420976708), // Int() - int64(7837839688282259259), // Int() - int64(2518412263346885298), // Int() - int64(5617773211005988520), // Int() - int64(2339563716805116249), // Int() - int64(7144924247938981575), // Int() - int64(161231572858529631), // Int() - int64(7259475919510918339), // Int() - int64(7373105480197164748), // Int() - int32(2029793274), // Int31() - int32(526058514), // Int31() - int32(1408655353), // Int31() - int32(116702506), // Int31() - int32(789387515), // Int31() - int32(621654496), // Int31() - int32(413258767), // Int31() - int32(1407315077), // Int31() - int32(1926657288), // Int31() - int32(359390928), // Int31() - int32(619732968), // Int31() - int32(1938329147), // Int31() - int32(1824889259), // Int31() - int32(586363548), // Int31() - int32(1307989752), // Int31() - int32(544722126), // Int31() - int32(1663557311), // Int31() - int32(37539650), // Int31() - int32(1690228450), // Int31() - int32(1716684894), // Int31() - int32(0), // Int31n(1) - int32(4), // Int31n(10) - int32(25), // Int31n(32) - int32(310570), // Int31n(1048576) - int32(857611), // Int31n(1048577) - int32(621654496), // Int31n(1000000000) - int32(413258767), // Int31n(1073741824) - int32(1407315077), // Int31n(2147483646) - int32(1926657288), // Int31n(2147483647) - int32(0), // Int31n(1) - int32(8), // Int31n(10) - int32(27), // Int31n(32) - int32(367019), // Int31n(1048576) - int32(209005), // Int31n(1048577) - int32(307989752), // Int31n(1000000000) - int32(544722126), // Int31n(1073741824) - int32(1663557311), // Int31n(2147483646) - int32(37539650), // Int31n(2147483647) - int32(0), // Int31n(1) - int32(4), // Int31n(10) - int64(8717895732742165505), // Int63() - int64(2259404117704393152), // Int63() - int64(6050128673802995827), // Int63() - int64(501233450539197794), // Int63() - int64(3390393562759376202), // Int63() - int64(2669985732393126063), // Int63() - int64(1774932891286980153), // Int63() - int64(6044372234677422456), // Int63() - int64(8274930044578894929), // Int63() - int64(1543572285742637646), // Int63() - int64(2661732831099943416), // Int63() - int64(8325060299420976708), // Int63() - int64(7837839688282259259), // Int63() - int64(2518412263346885298), // Int63() - int64(5617773211005988520), // Int63() - int64(2339563716805116249), // Int63() - int64(7144924247938981575), // Int63() - int64(161231572858529631), // Int63() - int64(7259475919510918339), // Int63() - int64(7373105480197164748), // Int63() - int64(0), // Int63n(1) - int64(2), // Int63n(10) - int64(19), // Int63n(32) - int64(959842), // Int63n(1048576) - int64(688912), // Int63n(1048577) - int64(393126063), // Int63n(1000000000) - int64(89212473), // Int63n(1073741824) - int64(834026388), // Int63n(2147483646) - int64(1577188963), // Int63n(2147483647) - int64(543572285742637646), // Int63n(1000000000000000000) - int64(355889821886249464), // Int63n(1152921504606846976) - int64(8325060299420976708), // Int63n(9223372036854775806) - int64(7837839688282259259), // Int63n(9223372036854775807) - int64(0), // Int63n(1) - int64(0), // Int63n(10) - int64(25), // Int63n(32) - int64(679623), // Int63n(1048576) - int64(882178), // Int63n(1048577) - int64(510918339), // Int63n(1000000000) - int64(782454476), // Int63n(1073741824) - int64(0), // Intn(1) - int64(4), // Intn(10) - int64(25), // Intn(32) - int64(310570), // Intn(1048576) - int64(857611), // Intn(1048577) - int64(621654496), // Intn(1000000000) - int64(413258767), // Intn(1073741824) - int64(1407315077), // Intn(2147483646) - int64(1926657288), // Intn(2147483647) - int64(543572285742637646), // Intn(1000000000000000000) - int64(355889821886249464), // Intn(1152921504606846976) - int64(8325060299420976708), // Intn(9223372036854775806) - int64(7837839688282259259), // Intn(9223372036854775807) - int64(0), // Intn(1) - int64(2), // Intn(10) - int64(14), // Intn(32) - int64(515775), // Intn(1048576) - int64(839455), // Intn(1048577) - int64(690228450), // Intn(1000000000) - int64(642943070), // Intn(1073741824) - float64(-0.28158587086436215), // NormFloat64() - float64(0.570933095808067), // NormFloat64() - float64(-1.6920196326157044), // NormFloat64() - float64(0.1996229111693099), // NormFloat64() - float64(1.9195199291234621), // NormFloat64() - float64(0.8954838794918353), // NormFloat64() - float64(0.41457072128813166), // NormFloat64() - float64(-0.48700161491544713), // NormFloat64() - float64(-0.1684059662402393), // NormFloat64() - float64(0.37056410998929545), // NormFloat64() - float64(1.0156889027029008), // NormFloat64() - float64(-0.5174422210625114), // NormFloat64() - float64(-0.5565834214413804), // NormFloat64() - float64(0.778320596648391), // NormFloat64() - float64(-1.8970718197702225), // NormFloat64() - float64(0.5229525761688676), // NormFloat64() - float64(-1.5515595563231523), // NormFloat64() - float64(0.0182029289376123), // NormFloat64() - float64(-0.6820951356608795), // NormFloat64() - float64(-0.5987943422687668), // NormFloat64() - []int{}, // Perm(0) - []int{0}, // Perm(1) - []int{0, 4, 1, 3, 2}, // Perm(5) - []int{3, 1, 0, 4, 7, 5, 2, 6}, // Perm(8) - []int{5, 0, 3, 6, 7, 4, 2, 1, 8}, // Perm(9) - []int{4, 5, 0, 2, 6, 9, 3, 1, 8, 7}, // Perm(10) + float64(4.668112973579268), // ExpFloat64() + float64(0.1601593871172866), // ExpFloat64() + float64(3.0465834105636), // ExpFloat64() + float64(0.06385839451671879), // ExpFloat64() + float64(1.8578917487258961), // ExpFloat64() + float64(0.784676123472182), // ExpFloat64() + float64(0.11225477361256932), // ExpFloat64() + float64(0.20173283329802255), // ExpFloat64() + float64(0.3468619496201105), // ExpFloat64() + float64(0.35601103454384536), // ExpFloat64() + float64(0.888376329507869), // ExpFloat64() + float64(1.4081362450365698), // ExpFloat64() + float64(1.0077753823151994), // ExpFloat64() + float64(0.23594100766227588), // ExpFloat64() + float64(2.777245612300007), // ExpFloat64() + float64(0.5202997830662377), // ExpFloat64() + float64(1.2842705247770294), // ExpFloat64() + float64(0.030307408362776206), // ExpFloat64() + float64(2.204156824853721), // ExpFloat64() + float64(2.09891923895058), // ExpFloat64() + float32(0.94519615), // Float32() + float32(0.24496509), // Float32() + float32(0.65595627), // Float32() + float32(0.05434384), // Float32() + float32(0.3675872), // Float32() + float32(0.28948045), // Float32() + float32(0.1924386), // Float32() + float32(0.65533215), // Float32() + float32(0.8971697), // Float32() + float32(0.16735445), // Float32() + float32(0.28858566), // Float32() + float32(0.9026048), // Float32() + float32(0.84978026), // Float32() + float32(0.2730468), // Float32() + float32(0.6090802), // Float32() + float32(0.253656), // Float32() + float32(0.7746542), // Float32() + float32(0.017480763), // Float32() + float32(0.78707397), // Float32() + float32(0.7993937), // Float32() + float64(0.9451961492941164), // Float64() + float64(0.24496508529377975), // Float64() + float64(0.6559562651954052), // Float64() + float64(0.05434383959970039), // Float64() + float64(0.36758720663245853), // Float64() + float64(0.2894804331565928), // Float64() + float64(0.19243860967493215), // Float64() + float64(0.6553321508148324), // Float64() + float64(0.897169713149801), // Float64() + float64(0.16735444255905835), // Float64() + float64(0.2885856518054551), // Float64() + float64(0.9026048462705047), // Float64() + float64(0.8497802817628735), // Float64() + float64(0.2730468047134829), // Float64() + float64(0.6090801919903561), // Float64() + float64(0.25365600644283687), // Float64() + float64(0.7746542391859803), // Float64() + float64(0.017480762156647272), // Float64() + float64(0.7870739563039942), // Float64() + float64(0.7993936979594545), // Float64() + int64(8717895732742165505), // Int() + int64(2259404117704393152), // Int() + int64(6050128673802995827), // Int() + int64(501233450539197794), // Int() + int64(3390393562759376202), // Int() + int64(2669985732393126063), // Int() + int64(1774932891286980153), // Int() + int64(6044372234677422456), // Int() + int64(8274930044578894929), // Int() + int64(1543572285742637646), // Int() + int64(2661732831099943416), // Int() + int64(8325060299420976708), // Int() + int64(7837839688282259259), // Int() + int64(2518412263346885298), // Int() + int64(5617773211005988520), // Int() + int64(2339563716805116249), // Int() + int64(7144924247938981575), // Int() + int64(161231572858529631), // Int() + int64(7259475919510918339), // Int() + int64(7373105480197164748), // Int() + int32(2029793274), // Int31() + int32(526058514), // Int31() + int32(1408655353), // Int31() + int32(116702506), // Int31() + int32(789387515), // Int31() + int32(621654496), // Int31() + int32(413258767), // Int31() + int32(1407315077), // Int31() + int32(1926657288), // Int31() + int32(359390928), // Int31() + int32(619732968), // Int31() + int32(1938329147), // Int31() + int32(1824889259), // Int31() + int32(586363548), // Int31() + int32(1307989752), // Int31() + int32(544722126), // Int31() + int32(1663557311), // Int31() + int32(37539650), // Int31() + int32(1690228450), // Int31() + int32(1716684894), // Int31() + int32(0), // Int31n(1) + int32(4), // Int31n(10) + int32(25), // Int31n(32) + int32(310570), // Int31n(1048576) + int32(857611), // Int31n(1048577) + int32(621654496), // Int31n(1000000000) + int32(413258767), // Int31n(1073741824) + int32(1407315077), // Int31n(2147483646) + int32(1926657288), // Int31n(2147483647) + int32(0), // Int31n(1) + int32(8), // Int31n(10) + int32(27), // Int31n(32) + int32(367019), // Int31n(1048576) + int32(209005), // Int31n(1048577) + int32(307989752), // Int31n(1000000000) + int32(544722126), // Int31n(1073741824) + int32(1663557311), // Int31n(2147483646) + int32(37539650), // Int31n(2147483647) + int32(0), // Int31n(1) + int32(4), // Int31n(10) + int64(8717895732742165505), // Int63() + int64(2259404117704393152), // Int63() + int64(6050128673802995827), // Int63() + int64(501233450539197794), // Int63() + int64(3390393562759376202), // Int63() + int64(2669985732393126063), // Int63() + int64(1774932891286980153), // Int63() + int64(6044372234677422456), // Int63() + int64(8274930044578894929), // Int63() + int64(1543572285742637646), // Int63() + int64(2661732831099943416), // Int63() + int64(8325060299420976708), // Int63() + int64(7837839688282259259), // Int63() + int64(2518412263346885298), // Int63() + int64(5617773211005988520), // Int63() + int64(2339563716805116249), // Int63() + int64(7144924247938981575), // Int63() + int64(161231572858529631), // Int63() + int64(7259475919510918339), // Int63() + int64(7373105480197164748), // Int63() + int64(0), // Int63n(1) + int64(2), // Int63n(10) + int64(19), // Int63n(32) + int64(959842), // Int63n(1048576) + int64(688912), // Int63n(1048577) + int64(393126063), // Int63n(1000000000) + int64(89212473), // Int63n(1073741824) + int64(834026388), // Int63n(2147483646) + int64(1577188963), // Int63n(2147483647) + int64(543572285742637646), // Int63n(1000000000000000000) + int64(355889821886249464), // Int63n(1152921504606846976) + int64(8325060299420976708), // Int63n(9223372036854775806) + int64(7837839688282259259), // Int63n(9223372036854775807) + int64(0), // Int63n(1) + int64(0), // Int63n(10) + int64(25), // Int63n(32) + int64(679623), // Int63n(1048576) + int64(882178), // Int63n(1048577) + int64(510918339), // Int63n(1000000000) + int64(782454476), // Int63n(1073741824) + int64(0), // Intn(1) + int64(4), // Intn(10) + int64(25), // Intn(32) + int64(310570), // Intn(1048576) + int64(857611), // Intn(1048577) + int64(621654496), // Intn(1000000000) + int64(413258767), // Intn(1073741824) + int64(1407315077), // Intn(2147483646) + int64(1926657288), // Intn(2147483647) + int64(543572285742637646), // Intn(1000000000000000000) + int64(355889821886249464), // Intn(1152921504606846976) + int64(8325060299420976708), // Intn(9223372036854775806) + int64(7837839688282259259), // Intn(9223372036854775807) + int64(0), // Intn(1) + int64(2), // Intn(10) + int64(14), // Intn(32) + int64(515775), // Intn(1048576) + int64(839455), // Intn(1048577) + int64(690228450), // Intn(1000000000) + int64(642943070), // Intn(1073741824) + float64(-0.28158587086436215), // NormFloat64() + float64(0.570933095808067), // NormFloat64() + float64(-1.6920196326157044), // NormFloat64() + float64(0.1996229111693099), // NormFloat64() + float64(1.9195199291234621), // NormFloat64() + float64(0.8954838794918353), // NormFloat64() + float64(0.41457072128813166), // NormFloat64() + float64(-0.48700161491544713), // NormFloat64() + float64(-0.1684059662402393), // NormFloat64() + float64(0.37056410998929545), // NormFloat64() + float64(1.0156889027029008), // NormFloat64() + float64(-0.5174422210625114), // NormFloat64() + float64(-0.5565834214413804), // NormFloat64() + float64(0.778320596648391), // NormFloat64() + float64(-1.8970718197702225), // NormFloat64() + float64(0.5229525761688676), // NormFloat64() + float64(-1.5515595563231523), // NormFloat64() + float64(0.0182029289376123), // NormFloat64() + float64(-0.6820951356608795), // NormFloat64() + float64(-0.5987943422687668), // NormFloat64() + []int{}, // Perm(0) + []int{0}, // Perm(1) + []int{0, 4, 1, 3, 2}, // Perm(5) + []int{3, 1, 0, 4, 7, 5, 2, 6}, // Perm(8) + []int{5, 0, 3, 6, 7, 4, 2, 1, 8}, // Perm(9) + []int{4, 5, 0, 2, 6, 9, 3, 1, 8, 7}, // Perm(10) []int{14, 2, 0, 8, 3, 5, 13, 12, 1, 4, 6, 7, 11, 9, 15, 10}, // Perm(16) - []int{}, // Perm(0) - []int{0}, // Perm(1) - []int{3, 0, 1, 2, 4}, // Perm(5) - []int{5, 1, 2, 0, 4, 7, 3, 6}, // Perm(8) - []int{4, 0, 6, 8, 1, 5, 2, 7, 3}, // Perm(9) - []int{8, 6, 1, 7, 5, 4, 3, 2, 9, 0}, // Perm(10) + []int{}, // Perm(0) + []int{0}, // Perm(1) + []int{3, 0, 1, 2, 4}, // Perm(5) + []int{5, 1, 2, 0, 4, 7, 3, 6}, // Perm(8) + []int{4, 0, 6, 8, 1, 5, 2, 7, 3}, // Perm(9) + []int{8, 6, 1, 7, 5, 4, 3, 2, 9, 0}, // Perm(10) []int{0, 3, 13, 2, 15, 4, 10, 1, 8, 14, 7, 6, 12, 9, 5, 11}, // Perm(16) []int{}, // Perm(0) []int{0}, // Perm(1) @@ -346,7 +346,7 @@ var regressGolden = []interface{}{ []byte{0x41, 0xd3, 0xff, 0x12, 0x4, 0x5b, 0x73, 0xc8}, // Read([0 0 0 0 0 0 0 0]) []byte{0x6e, 0x4f, 0xf9, 0x5f, 0xf6, 0x62, 0xa5, 0xee, 0xe8}, // Read([0 0 0 0 0 0 0 0 0]) []byte{0x2a, 0xbd, 0xf4, 0x4a, 0x2d, 0xb, 0x75, 0xfb, 0x18, 0xd}, // Read([0 0 0 0 0 0 0 0 0 0]) - []byte{0xaf}, // Read([0]) + []byte{0xaf}, // Read([0]) []byte{0x48, 0xa7, 0x9e, 0xe0, 0xb1, 0xd, 0x39}, // Read([0 0 0 0 0 0 0]) []byte{0x46, 0x51, 0x85, 0xf, 0xd4, 0xa1, 0x78, 0x89}, // Read([0 0 0 0 0 0 0 0]) []byte{0x2e, 0xe2, 0x85, 0xec, 0xe1, 0x51, 0x14, 0x55, 0x78}, // Read([0 0 0 0 0 0 0 0 0]) @@ -356,49 +356,49 @@ var regressGolden = []interface{}{ []byte{0xc6, 0xb1, 0xf8, 0x3b, 0x8e, 0x88, 0x3b, 0xbf}, // Read([0 0 0 0 0 0 0 0]) []byte{0x85, 0x7a, 0xab, 0x99, 0xc5, 0xb2, 0x52, 0xc7, 0x42}, // Read([0 0 0 0 0 0 0 0 0]) []byte{0x9c, 0x32, 0xf3, 0xa8, 0xae, 0xb7, 0x9e, 0xf8, 0x56, 0xf6}, // Read([0 0 0 0 0 0 0 0 0 0]) - []byte{0x59}, // Read([0]) + []byte{0x59}, // Read([0]) []byte{0xc1, 0x8f, 0xd, 0xce, 0xcc, 0x77, 0xc7}, // Read([0 0 0 0 0 0 0]) []byte{0x5e, 0x7a, 0x81, 0xbf, 0xde, 0x27, 0x5f, 0x67}, // Read([0 0 0 0 0 0 0 0]) []byte{0xcf, 0xe2, 0x42, 0xcf, 0x3c, 0xc3, 0x54, 0xf3, 0xed}, // Read([0 0 0 0 0 0 0 0 0]) []byte{0xe2, 0xd6, 0xbe, 0xcc, 0x4e, 0xa3, 0xae, 0x5e, 0x88, 0x52}, // Read([0 0 0 0 0 0 0 0 0 0]) - uint32(4059586549), // Uint32() - uint32(1052117029), // Uint32() - uint32(2817310706), // Uint32() - uint32(233405013), // Uint32() - uint32(1578775030), // Uint32() - uint32(1243308993), // Uint32() - uint32(826517535), // Uint32() - uint32(2814630155), // Uint32() - uint32(3853314576), // Uint32() - uint32(718781857), // Uint32() - uint32(1239465936), // Uint32() - uint32(3876658295), // Uint32() - uint32(3649778518), // Uint32() - uint32(1172727096), // Uint32() - uint32(2615979505), // Uint32() - uint32(1089444252), // Uint32() - uint32(3327114623), // Uint32() - uint32(75079301), // Uint32() - uint32(3380456901), // Uint32() - uint32(3433369789), // Uint32() - uint64(8717895732742165505), // Uint64() - uint64(2259404117704393152), // Uint64() - uint64(6050128673802995827), // Uint64() - uint64(9724605487393973602), // Uint64() - uint64(12613765599614152010), // Uint64() - uint64(11893357769247901871), // Uint64() - uint64(1774932891286980153), // Uint64() - uint64(15267744271532198264), // Uint64() - uint64(17498302081433670737), // Uint64() - uint64(1543572285742637646), // Uint64() - uint64(11885104867954719224), // Uint64() - uint64(17548432336275752516), // Uint64() - uint64(7837839688282259259), // Uint64() - uint64(2518412263346885298), // Uint64() - uint64(5617773211005988520), // Uint64() - uint64(11562935753659892057), // Uint64() - uint64(16368296284793757383), // Uint64() - uint64(161231572858529631), // Uint64() - uint64(16482847956365694147), // Uint64() - uint64(16596477517051940556), // Uint64() + uint32(4059586549), // Uint32() + uint32(1052117029), // Uint32() + uint32(2817310706), // Uint32() + uint32(233405013), // Uint32() + uint32(1578775030), // Uint32() + uint32(1243308993), // Uint32() + uint32(826517535), // Uint32() + uint32(2814630155), // Uint32() + uint32(3853314576), // Uint32() + uint32(718781857), // Uint32() + uint32(1239465936), // Uint32() + uint32(3876658295), // Uint32() + uint32(3649778518), // Uint32() + uint32(1172727096), // Uint32() + uint32(2615979505), // Uint32() + uint32(1089444252), // Uint32() + uint32(3327114623), // Uint32() + uint32(75079301), // Uint32() + uint32(3380456901), // Uint32() + uint32(3433369789), // Uint32() + uint64(8717895732742165505), // Uint64() + uint64(2259404117704393152), // Uint64() + uint64(6050128673802995827), // Uint64() + uint64(9724605487393973602), // Uint64() + uint64(12613765599614152010), // Uint64() + uint64(11893357769247901871), // Uint64() + uint64(1774932891286980153), // Uint64() + uint64(15267744271532198264), // Uint64() + uint64(17498302081433670737), // Uint64() + uint64(1543572285742637646), // Uint64() + uint64(11885104867954719224), // Uint64() + uint64(17548432336275752516), // Uint64() + uint64(7837839688282259259), // Uint64() + uint64(2518412263346885298), // Uint64() + uint64(5617773211005988520), // Uint64() + uint64(11562935753659892057), // Uint64() + uint64(16368296284793757383), // Uint64() + uint64(161231572858529631), // Uint64() + uint64(16482847956365694147), // Uint64() + uint64(16596477517051940556), // Uint64() } diff --git a/libgo/go/math/rand/rng.go b/libgo/go/math/rand/rng.go index f922417..f305df1 100644 --- a/libgo/go/math/rand/rng.go +++ b/libgo/go/math/rand/rng.go @@ -12,19 +12,16 @@ package rand */ const ( - _LEN = 607 - _TAP = 273 - _MAX = 1 << 63 - _MASK = _MAX - 1 - _A = 48271 - _M = (1 << 31) - 1 - _Q = 44488 - _R = 3399 + rngLen = 607 + rngTap = 273 + rngMax = 1 << 63 + rngMask = rngMax - 1 + int32max = (1 << 31) - 1 ) var ( - // Used for seeding. See gen_cooked.go for details. - rng_cooked [_LEN]int64 = [...]int64{ + // rngCooked used for seeding. See gen_cooked.go for details. + rngCooked [rngLen]int64 = [...]int64{ -4181792142133755926, -4576982950128230565, 1395769623340756751, 5333664234075297259, -6347679516498800754, 9033628115061424579, 7143218595135194537, 4812947590706362721, 7937252194349799378, 5307299880338848416, 8209348851763925077, -7107630437535961764, @@ -181,18 +178,24 @@ var ( ) type rngSource struct { - tap int // index into vec - feed int // index into vec - vec [_LEN]int64 // current feedback register + tap int // index into vec + feed int // index into vec + vec [rngLen]int64 // current feedback register } // seed rng x[n+1] = 48271 * x[n] mod (2**31 - 1) func seedrand(x int32) int32 { - hi := x / _Q - lo := x % _Q - x = _A*lo - _R*hi + const ( + A = 48271 + Q = 44488 + R = 3399 + ) + + hi := x / Q + lo := x % Q + x = A*lo - R*hi if x < 0 { - x += _M + x += int32max } return x } @@ -200,18 +203,18 @@ func seedrand(x int32) int32 { // Seed uses the provided seed value to initialize the generator to a deterministic state. func (rng *rngSource) Seed(seed int64) { rng.tap = 0 - rng.feed = _LEN - _TAP + rng.feed = rngLen - rngTap - seed = seed % _M + seed = seed % int32max if seed < 0 { - seed += _M + seed += int32max } if seed == 0 { seed = 89482311 } x := int32(seed) - for i := -20; i < _LEN; i++ { + for i := -20; i < rngLen; i++ { x = seedrand(x) if i >= 0 { var u int64 @@ -220,7 +223,7 @@ func (rng *rngSource) Seed(seed int64) { u ^= int64(x) << 20 x = seedrand(x) u ^= int64(x) - u ^= rng_cooked[i] + u ^= rngCooked[i] rng.vec[i] = u } } @@ -228,19 +231,19 @@ func (rng *rngSource) Seed(seed int64) { // Int63 returns a non-negative pseudo-random 63-bit integer as an int64. func (rng *rngSource) Int63() int64 { - return int64(rng.Uint64() & _MASK) + return int64(rng.Uint64() & rngMask) } // Uint64 returns a non-negative pseudo-random 64-bit integer as an uint64. func (rng *rngSource) Uint64() uint64 { rng.tap-- if rng.tap < 0 { - rng.tap += _LEN + rng.tap += rngLen } rng.feed-- if rng.feed < 0 { - rng.feed += _LEN + rng.feed += rngLen } x := rng.vec[rng.feed] + rng.vec[rng.tap] diff --git a/libgo/go/math/sin.go b/libgo/go/math/sin.go index 75579ab..e3166e2 100644 --- a/libgo/go/math/sin.go +++ b/libgo/go/math/sin.go @@ -137,9 +137,7 @@ func cos(x float64) float64 { // make argument positive sign := false - if x < 0 { - x = -x - } + x = Abs(x) j := int64(x * M4PI) // integer part of x/(Pi/4), as integer for tests on the phase angle y := float64(j) // integer part of x/(Pi/4), as float diff --git a/libgo/go/math/sinh.go b/libgo/go/math/sinh.go index 139b911..9dff71e 100644 --- a/libgo/go/math/sinh.go +++ b/libgo/go/math/sinh.go @@ -43,10 +43,11 @@ func Sinh(x float64) float64 { var temp float64 switch true { case x > 21: - temp = Exp(x) / 2 + temp = Exp(x) * 0.5 case x > 0.5: - temp = (Exp(x) - Exp(-x)) / 2 + ex := Exp(x) + temp = (ex - 1/ex) * 0.5 default: sq := x * x @@ -67,11 +68,10 @@ func Sinh(x float64) float64 { // Cosh(±Inf) = +Inf // Cosh(NaN) = NaN func Cosh(x float64) float64 { - if x < 0 { - x = -x - } + x = Abs(x) if x > 21 { - return Exp(x) / 2 + return Exp(x) * 0.5 } - return (Exp(x) + Exp(-x)) / 2 + ex := Exp(x) + return (ex + 1/ex) * 0.5 } |