libgo: update to go1.7rc3

    
    Reviewed-on: https://go-review.googlesource.com/25150

From-SVN: r238662

44 files changed, 759 insertions, 441 deletions

diff --git a/libgo/go/math/acosh.go b/libgo/go/math/acosh.go
index e394008..dce21b2 100644
--- a/libgo/go/math/acosh.go
+++ b/libgo/go/math/acosh.go
@@ -6,7 +6,7 @@ package math
 
 // The original C code, the long comment, and the constants
 // below are from FreeBSD's /usr/src/lib/msun/src/e_acosh.c
-// and came with this notice.  The go code is a simplified
+// and came with this notice. The go code is a simplified
 // version of the original C.
 //
 // ====================================================
diff --git a/libgo/go/math/all_test.go b/libgo/go/math/all_test.go
index 968a7b1..d9ea1fd 100644
--- a/libgo/go/math/all_test.go
+++ b/libgo/go/math/all_test.go
@@ -1225,12 +1225,6 @@ var hypotSC = []float64{
 	NaN(),
 }
 
-var vfilogbSC = []float64{
-	Inf(-1),
-	0,
-	Inf(1),
-	NaN(),
-}
 var ilogbSC = []int{
 	MaxInt32,
 	MinInt32,
@@ -1756,7 +1750,6 @@ func tolerance(a, b, e float64) bool {
 	}
 	return d < e
 }
-func kindaclose(a, b float64) bool { return tolerance(a, b, 1e-8) }
 func close(a, b float64) bool      { return tolerance(a, b, 1e-14) }
 func veryclose(a, b float64) bool  { return tolerance(a, b, 4e-16) }
 func soclose(a, b, e float64) bool { return tolerance(a, b, e) }
diff --git a/libgo/go/math/asinh.go b/libgo/go/math/asinh.go
index ff2de02..3b793b0 100644
--- a/libgo/go/math/asinh.go
+++ b/libgo/go/math/asinh.go
@@ -6,7 +6,7 @@ package math
 
 // The original C code, the long comment, and the constants
 // below are from FreeBSD's /usr/src/lib/msun/src/s_asinh.c
-// and came with this notice.  The go code is a simplified
+// and came with this notice. The go code is a simplified
 // version of the original C.
 //
 // ====================================================
diff --git a/libgo/go/math/atanh.go b/libgo/go/math/atanh.go
index 113d5c1..d59a847 100644
--- a/libgo/go/math/atanh.go
+++ b/libgo/go/math/atanh.go
@@ -6,7 +6,7 @@ package math
 
 // The original C code, the long comment, and the constants
 // below are from FreeBSD's /usr/src/lib/msun/src/e_atanh.c
-// and came with this notice.  The go code is a simplified
+// and came with this notice. The go code is a simplified
 // version of the original C.
 //
 // ====================================================
diff --git a/libgo/go/math/big/arith_decl.go b/libgo/go/math/big/arith_decl.go
index 1707aa4..5433b6d 100644
--- a/libgo/go/math/big/arith_decl.go
+++ b/libgo/go/math/big/arith_decl.go
@@ -1,4 +1,4 @@
-// Copyright 2010 The Go Authors.  All rights reserved.
+// 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.
 
diff --git a/libgo/go/math/big/arith_decl_pure.go b/libgo/go/math/big/arith_decl_pure.go
index e760a38..21775dd 100644
--- a/libgo/go/math/big/arith_decl_pure.go
+++ b/libgo/go/math/big/arith_decl_pure.go
@@ -1,4 +1,4 @@
-// Copyright 2015 The Go Authors.  All rights reserved.
+// 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.
 
diff --git a/libgo/go/math/big/arith_test.go b/libgo/go/math/big/arith_test.go
index f46a494..75862b4 100644
--- a/libgo/go/math/big/arith_test.go
+++ b/libgo/go/math/big/arith_test.go
@@ -5,6 +5,7 @@
 package big
 
 import (
+	"fmt"
 	"math/rand"
 	"testing"
 )
@@ -118,28 +119,22 @@ func rndV(n int) []Word {
 	return v
 }
 
-func benchmarkFunVV(b *testing.B, f funVV, n int) {
-	x := rndV(n)
-	y := rndV(n)
-	z := make([]Word, n)
-	b.SetBytes(int64(n * _W))
-	b.ResetTimer()
-	for i := 0; i < b.N; i++ {
-		f(z, x, y)
+var benchSizes = []int{1, 2, 3, 4, 5, 1e1, 1e2, 1e3, 1e4, 1e5}
+
+func BenchmarkAddVV(b *testing.B) {
+	for _, n := range benchSizes {
+		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++ {
+				addVV(z, x, y)
+			}
+		})
 	}
 }
 
-func BenchmarkAddVV_1(b *testing.B)   { benchmarkFunVV(b, addVV, 1) }
-func BenchmarkAddVV_2(b *testing.B)   { benchmarkFunVV(b, addVV, 2) }
-func BenchmarkAddVV_3(b *testing.B)   { benchmarkFunVV(b, addVV, 3) }
-func BenchmarkAddVV_4(b *testing.B)   { benchmarkFunVV(b, addVV, 4) }
-func BenchmarkAddVV_5(b *testing.B)   { benchmarkFunVV(b, addVV, 5) }
-func BenchmarkAddVV_1e1(b *testing.B) { benchmarkFunVV(b, addVV, 1e1) }
-func BenchmarkAddVV_1e2(b *testing.B) { benchmarkFunVV(b, addVV, 1e2) }
-func BenchmarkAddVV_1e3(b *testing.B) { benchmarkFunVV(b, addVV, 1e3) }
-func BenchmarkAddVV_1e4(b *testing.B) { benchmarkFunVV(b, addVV, 1e4) }
-func BenchmarkAddVV_1e5(b *testing.B) { benchmarkFunVV(b, addVV, 1e5) }
-
 type funVW func(z, x []Word, y Word) (c Word)
 type argVW struct {
 	z, x nat
@@ -158,21 +153,6 @@ var sumVW = []argVW{
 	{nat{585}, nat{314}, 271, 0},
 }
 
-var prodVW = []argVW{
-	{},
-	{nat{0}, nat{0}, 0, 0},
-	{nat{0}, nat{_M}, 0, 0},
-	{nat{0}, nat{0}, _M, 0},
-	{nat{1}, nat{1}, 1, 0},
-	{nat{22793}, nat{991}, 23, 0},
-	{nat{0, 0, 0, 22793}, nat{0, 0, 0, 991}, 23, 0},
-	{nat{0, 0, 0, 0}, nat{7893475, 7395495, 798547395, 68943}, 0, 0},
-	{nat{0, 0, 0, 0}, nat{0, 0, 0, 0}, 894375984, 0},
-	{nat{_M << 1 & _M}, nat{_M}, 1 << 1, _M >> (_W - 1)},
-	{nat{_M << 7 & _M}, nat{_M}, 1 << 7, _M >> (_W - 7)},
-	{nat{_M << 7 & _M, _M, _M, _M}, nat{_M, _M, _M, _M}, 1 << 7, _M >> (_W - 7)},
-}
-
 var lshVW = []argVW{
 	{},
 	{nat{0}, nat{0}, 0, 0},
@@ -251,28 +231,20 @@ func TestFunVW(t *testing.T) {
 	}
 }
 
-func benchmarkFunVW(b *testing.B, f funVW, n int) {
-	x := rndV(n)
-	y := rndW()
-	z := make([]Word, n)
-	b.SetBytes(int64(n * _S))
-	b.ResetTimer()
-	for i := 0; i < b.N; i++ {
-		f(z, x, y)
+func BenchmarkAddVW(b *testing.B) {
+	for _, n := range benchSizes {
+		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++ {
+				addVW(z, x, y)
+			}
+		})
 	}
 }
 
-func BenchmarkAddVW_1(b *testing.B)   { benchmarkFunVW(b, addVW, 1) }
-func BenchmarkAddVW_2(b *testing.B)   { benchmarkFunVW(b, addVW, 2) }
-func BenchmarkAddVW_3(b *testing.B)   { benchmarkFunVW(b, addVW, 3) }
-func BenchmarkAddVW_4(b *testing.B)   { benchmarkFunVW(b, addVW, 4) }
-func BenchmarkAddVW_5(b *testing.B)   { benchmarkFunVW(b, addVW, 5) }
-func BenchmarkAddVW_1e1(b *testing.B) { benchmarkFunVW(b, addVW, 1e1) }
-func BenchmarkAddVW_1e2(b *testing.B) { benchmarkFunVW(b, addVW, 1e2) }
-func BenchmarkAddVW_1e3(b *testing.B) { benchmarkFunVW(b, addVW, 1e3) }
-func BenchmarkAddVW_1e4(b *testing.B) { benchmarkFunVW(b, addVW, 1e4) }
-func BenchmarkAddVW_1e5(b *testing.B) { benchmarkFunVW(b, addVW, 1e5) }
-
 type funVWW func(z, x []Word, y, r Word) (c Word)
 type argVWW struct {
 	z, x nat
@@ -397,28 +369,20 @@ func TestMulAddWWW(t *testing.T) {
 	}
 }
 
-func benchmarkAddMulVVW(b *testing.B, n int) {
-	x := rndV(n)
-	y := rndW()
-	z := make([]Word, n)
-	b.SetBytes(int64(n * _W))
-	b.ResetTimer()
-	for i := 0; i < b.N; i++ {
-		addMulVVW(z, x, y)
+func BenchmarkAddMulVVW(b *testing.B) {
+	for _, n := range benchSizes {
+		x := rndV(n)
+		y := rndW()
+		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++ {
+				addMulVVW(z, x, y)
+			}
+		})
 	}
 }
 
-func BenchmarkAddMulVVW_1(b *testing.B)   { benchmarkAddMulVVW(b, 1) }
-func BenchmarkAddMulVVW_2(b *testing.B)   { benchmarkAddMulVVW(b, 2) }
-func BenchmarkAddMulVVW_3(b *testing.B)   { benchmarkAddMulVVW(b, 3) }
-func BenchmarkAddMulVVW_4(b *testing.B)   { benchmarkAddMulVVW(b, 4) }
-func BenchmarkAddMulVVW_5(b *testing.B)   { benchmarkAddMulVVW(b, 5) }
-func BenchmarkAddMulVVW_1e1(b *testing.B) { benchmarkAddMulVVW(b, 1e1) }
-func BenchmarkAddMulVVW_1e2(b *testing.B) { benchmarkAddMulVVW(b, 1e2) }
-func BenchmarkAddMulVVW_1e3(b *testing.B) { benchmarkAddMulVVW(b, 1e3) }
-func BenchmarkAddMulVVW_1e4(b *testing.B) { benchmarkAddMulVVW(b, 1e4) }
-func BenchmarkAddMulVVW_1e5(b *testing.B) { benchmarkAddMulVVW(b, 1e5) }
-
 func testWordBitLen(t *testing.T, fname string, f func(Word) int) {
 	for i := 0; i <= _W; i++ {
 		x := Word(1) << uint(i-1) // i == 0 => x == 0
@@ -435,23 +399,15 @@ func TestWordBitLen(t *testing.T) {
 }
 
 // runs b.N iterations of bitLen called on a Word containing (1 << nbits)-1.
-func benchmarkBitLenN(b *testing.B, nbits uint) {
-	testword := Word((uint64(1) << nbits) - 1)
-	for i := 0; i < b.N; i++ {
-		bitLen(testword)
+func BenchmarkBitLen(b *testing.B) {
+	// Individual bitLen tests. Numbers chosen to examine both sides
+	// of powers-of-two boundaries.
+	for _, nbits := range []uint{0, 1, 2, 3, 4, 5, 8, 9, 16, 17, 31} {
+		testword := Word((uint64(1) << nbits) - 1)
+		b.Run(fmt.Sprint(nbits), func(b *testing.B) {
+			for i := 0; i < b.N; i++ {
+				bitLen(testword)
+			}
+		})
 	}
 }
-
-// Individual bitLen tests.  Numbers chosen to examine both sides
-// of powers-of-two boundaries.
-func BenchmarkBitLen0(b *testing.B)  { benchmarkBitLenN(b, 0) }
-func BenchmarkBitLen1(b *testing.B)  { benchmarkBitLenN(b, 1) }
-func BenchmarkBitLen2(b *testing.B)  { benchmarkBitLenN(b, 2) }
-func BenchmarkBitLen3(b *testing.B)  { benchmarkBitLenN(b, 3) }
-func BenchmarkBitLen4(b *testing.B)  { benchmarkBitLenN(b, 4) }
-func BenchmarkBitLen5(b *testing.B)  { benchmarkBitLenN(b, 5) }
-func BenchmarkBitLen8(b *testing.B)  { benchmarkBitLenN(b, 8) }
-func BenchmarkBitLen9(b *testing.B)  { benchmarkBitLenN(b, 9) }
-func BenchmarkBitLen16(b *testing.B) { benchmarkBitLenN(b, 16) }
-func BenchmarkBitLen17(b *testing.B) { benchmarkBitLenN(b, 17) }
-func BenchmarkBitLen31(b *testing.B) { benchmarkBitLenN(b, 31) }
diff --git a/libgo/go/math/big/float.go b/libgo/go/math/big/float.go
index b1c748c..7a9c2b3 100644
--- a/libgo/go/math/big/float.go
+++ b/libgo/go/math/big/float.go
@@ -392,15 +392,13 @@ func (z *Float) round(sbit uint) {
 	// m > 0 implies z.prec > 0 (checked by validate)
 
 	m := uint32(len(z.mant)) // present mantissa length in words
-	bits := m * _W           // present mantissa bits
+	bits := m * _W           // present mantissa bits; bits > 0
 	if bits <= z.prec {
 		// mantissa fits => nothing to do
 		return
 	}
 	// bits > z.prec
 
-	n := (z.prec + (_W - 1)) / _W // mantissa length in words for desired precision
-
 	// Rounding is based on two bits: the rounding bit (rbit) and the
 	// sticky bit (sbit). The rbit is the bit immediately before the
 	// z.prec leading mantissa bits (the "0.5"). The sbit is set if any
@@ -415,111 +413,77 @@ func (z *Float) round(sbit uint) {
 
 	// bits > z.prec: mantissa too large => round
 	r := uint(bits - z.prec - 1) // rounding bit position; r >= 0
-	rbit := z.mant.bit(r)        // rounding bit
+	rbit := z.mant.bit(r) & 1    // rounding bit; be safe and ensure it's a single bit
 	if sbit == 0 {
+		// TODO(gri) if rbit != 0 we don't need to compute sbit for some rounding modes (optimization)
 		sbit = z.mant.sticky(r)
 	}
-	if debugFloat && sbit&^1 != 0 {
-		panic(fmt.Sprintf("invalid sbit %#x", sbit))
-	}
-
-	// convert ToXInf rounding modes
-	mode := z.mode
-	switch mode {
-	case ToNegativeInf:
-		mode = ToZero
-		if z.neg {
-			mode = AwayFromZero
-		}
-	case ToPositiveInf:
-		mode = AwayFromZero
-		if z.neg {
-			mode = ToZero
-		}
-	}
+	sbit &= 1 // be safe and ensure it's a single bit
 
 	// cut off extra words
+	n := (z.prec + (_W - 1)) / _W // mantissa length in words for desired precision
 	if m > n {
 		copy(z.mant, z.mant[m-n:]) // move n last words to front
 		z.mant = z.mant[:n]
 	}
 
-	// determine number of trailing zero bits t
-	t := n*_W - z.prec // 0 <= t < _W
-	lsb := Word(1) << t
-
-	// make rounding decision
-	// TODO(gri) This can be simplified (see Bits.round in bits_test.go).
-	switch mode {
-	case ToZero:
-		// nothing to do
-	case ToNearestEven, ToNearestAway:
-		if rbit == 0 {
-			// rounding bits == 0b0x
-			mode = ToZero
-		} else if sbit == 1 {
-			// rounding bits == 0b11
-			mode = AwayFromZero
-		}
-	case AwayFromZero:
-		if rbit|sbit == 0 {
-			mode = ToZero
-		}
-	default:
-		// ToXInf modes have been converted to ToZero or AwayFromZero
-		panic("unreachable")
-	}
-
-	// round and determine accuracy
-	switch mode {
-	case ToZero:
-		if rbit|sbit != 0 {
-			z.acc = Below
+	// determine number of trailing zero bits (ntz) and compute lsb mask of mantissa's least-significant word
+	ntz := n*_W - z.prec // 0 <= ntz < _W
+	lsb := Word(1) << ntz
+
+	// round if result is inexact
+	if rbit|sbit != 0 {
+		// Make rounding decision: The result mantissa is truncated ("rounded down")
+		// by default. Decide if we need to increment, or "round up", the (unsigned)
+		// mantissa.
+		inc := false
+		switch z.mode {
+		case ToNegativeInf:
+			inc = z.neg
+		case ToZero:
+			// nothing to do
+		case ToNearestEven:
+			inc = rbit != 0 && (sbit != 0 || z.mant[0]&lsb != 0)
+		case ToNearestAway:
+			inc = rbit != 0
+		case AwayFromZero:
+			inc = true
+		case ToPositiveInf:
+			inc = !z.neg
+		default:
+			panic("unreachable")
 		}
 
-	case ToNearestEven, ToNearestAway:
-		if debugFloat && rbit != 1 {
-			panic("internal error in rounding")
-		}
-		if mode == ToNearestEven && sbit == 0 && z.mant[0]&lsb == 0 {
-			z.acc = Below
-			break
-		}
-		// mode == ToNearestAway || sbit == 1 || z.mant[0]&lsb != 0
-		fallthrough
-
-	case AwayFromZero:
-		// add 1 to mantissa
-		if addVW(z.mant, z.mant, lsb) != 0 {
-			// overflow => shift mantissa right by 1 and add msb
-			shrVU(z.mant, z.mant, 1)
-			z.mant[n-1] |= 1 << (_W - 1)
-			// adjust exponent
-			if z.exp < MaxExp {
+		// A positive result (!z.neg) is Above the exact result if we increment,
+		// and it's Below if we truncate (Exact results require no rounding).
+		// For a negative result (z.neg) it is exactly the opposite.
+		z.acc = makeAcc(inc != z.neg)
+
+		if inc {
+			// add 1 to mantissa
+			if addVW(z.mant, z.mant, lsb) != 0 {
+				// mantissa overflow => adjust exponent
+				if z.exp >= MaxExp {
+					// exponent overflow
+					z.form = inf
+					return
+				}
 				z.exp++
-			} else {
-				// exponent overflow
-				z.acc = makeAcc(!z.neg)
-				z.form = inf
-				return
+				// adjust mantissa: divide by 2 to compensate for exponent adjustment
+				shrVU(z.mant, z.mant, 1)
+				// set msb == carry == 1 from the mantissa overflow above
+				const msb = 1 << (_W - 1)
+				z.mant[n-1] |= msb
 			}
 		}
-		z.acc = Above
 	}
 
 	// zero out trailing bits in least-significant word
 	z.mant[0] &^= lsb - 1
 
-	// update accuracy
-	if z.acc != Exact && z.neg {
-		z.acc = -z.acc
-	}
-
 	if debugFloat {
 		z.validate()
 	}
-
-	return
 }
 
 func (z *Float) setBits64(neg bool, x uint64) *Float {
@@ -874,21 +838,43 @@ func (x *Float) Float32() (float32, Accuracy) {
 			emax  = bias              //   127  largest unbiased exponent (normal)
 		)
 
-		// Float mantissa m is 0.5 <= m < 1.0; compute exponent for floatxx mantissa.
-		e := x.exp - 1 // exponent for mantissa m with 1.0 <= m < 2.0
-		p := mbits + 1 // precision of normal float
+		// Float mantissa m is 0.5 <= m < 1.0; compute exponent e for float32 mantissa.
+		e := x.exp - 1 // exponent for normal mantissa m with 1.0 <= m < 2.0
 
-		// If the exponent is too small, we may have a denormal number
-		// in which case we have fewer mantissa bits available: reduce
-		// precision accordingly.
+		// Compute precision p for float32 mantissa.
+		// If the exponent is too small, we have a denormal number before
+		// rounding and fewer than p mantissa bits of precision available
+		// (the exponent remains fixed but the mantissa gets shifted right).
+		p := mbits + 1 // precision of normal float
 		if e < emin {
-			p -= emin - int(e)
-			// Make sure we have at least 1 bit so that we don't
-			// lose numbers rounded up to the smallest denormal.
-			if p < 1 {
-				p = 1
+			// recompute precision
+			p = mbits + 1 - emin + int(e)
+			// If p == 0, the mantissa of x is shifted so much to the right
+			// that its msb falls immediately to the right of the float32
+			// mantissa space. In other words, if the smallest denormal is
+			// considered "1.0", for p == 0, the mantissa value m is >= 0.5.
+			// If m > 0.5, it is rounded up to 1.0; i.e., the smallest denormal.
+			// If m == 0.5, it is rounded down to even, i.e., 0.0.
+			// If p < 0, the mantissa value m is <= "0.25" which is never rounded up.
+			if p < 0 /* m <= 0.25 */ || p == 0 && x.mant.sticky(uint(len(x.mant))*_W-1) == 0 /* m == 0.5 */ {
+				// underflow to ±0
+				if x.neg {
+					var z float32
+					return -z, Above
+				}
+				return 0.0, Below
+			}
+			// otherwise, round up
+			// We handle p == 0 explicitly because it's easy and because
+			// Float.round doesn't support rounding to 0 bits of precision.
+			if p == 0 {
+				if x.neg {
+					return -math.SmallestNonzeroFloat32, Below
+				}
+				return math.SmallestNonzeroFloat32, Above
 			}
 		}
+		// p > 0
 
 		// round
 		var r Float
@@ -898,12 +884,8 @@ func (x *Float) Float32() (float32, Accuracy) {
 
 		// Rounding may have caused r to overflow to ±Inf
 		// (rounding never causes underflows to 0).
-		if r.form == inf {
-			e = emax + 1 // cause overflow below
-		}
-
-		// If the exponent is too large, overflow to ±Inf.
-		if e > emax {
+		// If the exponent is too large, also overflow to ±Inf.
+		if r.form == inf || e > emax {
 			// overflow
 			if x.neg {
 				return float32(math.Inf(-1)), Below
@@ -921,17 +903,12 @@ func (x *Float) Float32() (float32, Accuracy) {
 		// Rounding may have caused a denormal number to
 		// become normal. Check again.
 		if e < emin {
-			// denormal number
-			if e < dmin {
-				// underflow to ±0
-				if x.neg {
-					var z float32
-					return -z, Above
-				}
-				return 0.0, Below
-			}
-			// bexp = 0
-			mant = msb32(r.mant) >> (fbits - r.prec)
+			// denormal number: recompute precision
+			// Since rounding may have at best increased precision
+			// and we have eliminated p <= 0 early, we know p > 0.
+			// bexp == 0 for denormals
+			p = mbits + 1 - emin + int(e)
+			mant = msb32(r.mant) >> uint(fbits-p)
 		} else {
 			// normal number: emin <= e <= emax
 			bexp = uint32(e+bias) << mbits
@@ -981,21 +958,43 @@ func (x *Float) Float64() (float64, Accuracy) {
 			emax  = bias              //  1023  largest unbiased exponent (normal)
 		)
 
-		// Float mantissa m is 0.5 <= m < 1.0; compute exponent for floatxx mantissa.
-		e := x.exp - 1 // exponent for mantissa m with 1.0 <= m < 2.0
-		p := mbits + 1 // precision of normal float
+		// Float mantissa m is 0.5 <= m < 1.0; compute exponent e for float64 mantissa.
+		e := x.exp - 1 // exponent for normal mantissa m with 1.0 <= m < 2.0
 
-		// If the exponent is too small, we may have a denormal number
-		// in which case we have fewer mantissa bits available: reduce
-		// precision accordingly.
+		// Compute precision p for float64 mantissa.
+		// If the exponent is too small, we have a denormal number before
+		// rounding and fewer than p mantissa bits of precision available
+		// (the exponent remains fixed but the mantissa gets shifted right).
+		p := mbits + 1 // precision of normal float
 		if e < emin {
-			p -= emin - int(e)
-			// Make sure we have at least 1 bit so that we don't
-			// lose numbers rounded up to the smallest denormal.
-			if p < 1 {
-				p = 1
+			// recompute precision
+			p = mbits + 1 - emin + int(e)
+			// If p == 0, the mantissa of x is shifted so much to the right
+			// that its msb falls immediately to the right of the float64
+			// mantissa space. In other words, if the smallest denormal is
+			// considered "1.0", for p == 0, the mantissa value m is >= 0.5.
+			// If m > 0.5, it is rounded up to 1.0; i.e., the smallest denormal.
+			// If m == 0.5, it is rounded down to even, i.e., 0.0.
+			// If p < 0, the mantissa value m is <= "0.25" which is never rounded up.
+			if p < 0 /* m <= 0.25 */ || p == 0 && x.mant.sticky(uint(len(x.mant))*_W-1) == 0 /* m == 0.5 */ {
+				// underflow to ±0
+				if x.neg {
+					var z float64
+					return -z, Above
+				}
+				return 0.0, Below
+			}
+			// otherwise, round up
+			// We handle p == 0 explicitly because it's easy and because
+			// Float.round doesn't support rounding to 0 bits of precision.
+			if p == 0 {
+				if x.neg {
+					return -math.SmallestNonzeroFloat64, Below
+				}
+				return math.SmallestNonzeroFloat64, Above
 			}
 		}
+		// p > 0
 
 		// round
 		var r Float
@@ -1005,12 +1004,8 @@ func (x *Float) Float64() (float64, Accuracy) {
 
 		// Rounding may have caused r to overflow to ±Inf
 		// (rounding never causes underflows to 0).
-		if r.form == inf {
-			e = emax + 1 // cause overflow below
-		}
-
-		// If the exponent is too large, overflow to ±Inf.
-		if e > emax {
+		// If the exponent is too large, also overflow to ±Inf.
+		if r.form == inf || e > emax {
 			// overflow
 			if x.neg {
 				return math.Inf(-1), Below
@@ -1028,17 +1023,12 @@ func (x *Float) Float64() (float64, Accuracy) {
 		// Rounding may have caused a denormal number to
 		// become normal. Check again.
 		if e < emin {
-			// denormal number
-			if e < dmin {
-				// underflow to ±0
-				if x.neg {
-					var z float64
-					return -z, Above
-				}
-				return 0.0, Below
-			}
-			// bexp = 0
-			mant = msb64(r.mant) >> (fbits - r.prec)
+			// denormal number: recompute precision
+			// Since rounding may have at best increased precision
+			// and we have eliminated p <= 0 early, we know p > 0.
+			// bexp == 0 for denormals
+			p = mbits + 1 - emin + int(e)
+			mant = msb64(r.mant) >> uint(fbits-p)
 		} else {
 			// normal number: emin <= e <= emax
 			bexp = uint64(e+bias) << mbits
@@ -1427,7 +1417,7 @@ func (z *Float) Add(x, y *Float) *Float {
 	}
 
 	if x.form == finite && y.form == finite {
-		// x + y (commom case)
+		// x + y (common case)
 		z.neg = x.neg
 		if x.neg == y.neg {
 			// x + y == x + y
diff --git a/libgo/go/math/big/float_test.go b/libgo/go/math/big/float_test.go
index d3b214b..464619b 100644
--- a/libgo/go/math/big/float_test.go
+++ b/libgo/go/math/big/float_test.go
@@ -829,7 +829,7 @@ func TestFloatFloat32(t *testing.T) {
 	}{
 		{"0", 0, Exact},
 
-		// underflow
+		// underflow to zero
 		{"1e-1000", 0, Below},
 		{"0x0.000002p-127", 0, Below},
 		{"0x.0000010p-126", 0, Below},
@@ -843,6 +843,46 @@ func TestFloatFloat32(t *testing.T) {
 		{"1p-149", math.SmallestNonzeroFloat32, Exact},
 		{"0x.fffffep-126", math.Float32frombits(0x7fffff), Exact}, // largest denormal
 
+		// special denormal cases (see issues 14553, 14651)
+		{"0x0.0000001p-126", math.Float32frombits(0x00000000), Below}, // underflow to zero
+		{"0x0.0000008p-126", math.Float32frombits(0x00000000), Below}, // underflow to zero
+		{"0x0.0000010p-126", math.Float32frombits(0x00000000), Below}, // rounded down to even
+		{"0x0.0000011p-126", math.Float32frombits(0x00000001), Above}, // rounded up to smallest denormal
+		{"0x0.0000018p-126", math.Float32frombits(0x00000001), Above}, // rounded up to smallest denormal
+
+		{"0x1.0000000p-149", math.Float32frombits(0x00000001), Exact}, // smallest denormal
+		{"0x0.0000020p-126", math.Float32frombits(0x00000001), Exact}, // smallest denormal
+		{"0x0.fffffe0p-126", math.Float32frombits(0x007fffff), Exact}, // largest denormal
+		{"0x1.0000000p-126", math.Float32frombits(0x00800000), Exact}, // smallest normal
+
+		{"0x0.8p-149", math.Float32frombits(0x000000000), Below}, // rounded down to even
+		{"0x0.9p-149", math.Float32frombits(0x000000001), Above}, // rounded up to smallest denormal
+		{"0x0.ap-149", math.Float32frombits(0x000000001), Above}, // rounded up to smallest denormal
+		{"0x0.bp-149", math.Float32frombits(0x000000001), Above}, // rounded up to smallest denormal
+		{"0x0.cp-149", math.Float32frombits(0x000000001), Above}, // rounded up to smallest denormal
+
+		{"0x1.0p-149", math.Float32frombits(0x000000001), Exact}, // smallest denormal
+		{"0x1.7p-149", math.Float32frombits(0x000000001), Below},
+		{"0x1.8p-149", math.Float32frombits(0x000000002), Above},
+		{"0x1.9p-149", math.Float32frombits(0x000000002), Above},
+
+		{"0x2.0p-149", math.Float32frombits(0x000000002), Exact},
+		{"0x2.8p-149", math.Float32frombits(0x000000002), Below}, // rounded down to even
+		{"0x2.9p-149", math.Float32frombits(0x000000003), Above},
+
+		{"0x3.0p-149", math.Float32frombits(0x000000003), Exact},
+		{"0x3.7p-149", math.Float32frombits(0x000000003), Below},
+		{"0x3.8p-149", math.Float32frombits(0x000000004), Above}, // rounded up to even
+
+		{"0x4.0p-149", math.Float32frombits(0x000000004), Exact},
+		{"0x4.8p-149", math.Float32frombits(0x000000004), Below}, // rounded down to even
+		{"0x4.9p-149", math.Float32frombits(0x000000005), Above},
+
+		// specific case from issue 14553
+		{"0x7.7p-149", math.Float32frombits(0x000000007), Below},
+		{"0x7.8p-149", math.Float32frombits(0x000000008), Above},
+		{"0x7.9p-149", math.Float32frombits(0x000000008), Above},
+
 		// normals
 		{"0x.ffffffp-126", math.Float32frombits(0x00800000), Above}, // rounded up to smallest normal
 		{"1p-126", math.Float32frombits(0x00800000), Exact},         // smallest normal
@@ -881,7 +921,7 @@ func TestFloatFloat32(t *testing.T) {
 			x := makeFloat(tx)
 			out, acc := x.Float32()
 			if !alike32(out, tout) || acc != tacc {
-				t.Errorf("%s: got %g (%#x, %s); want %g (%#x, %s)", tx, out, math.Float32bits(out), acc, test.out, math.Float32bits(test.out), tacc)
+				t.Errorf("%s: got %g (%#08x, %s); want %g (%#08x, %s)", tx, out, math.Float32bits(out), acc, test.out, math.Float32bits(test.out), tacc)
 			}
 
 			// test that x.SetFloat64(float64(f)).Float32() == f
@@ -903,18 +943,48 @@ func TestFloatFloat64(t *testing.T) {
 	}{
 		{"0", 0, Exact},
 
-		// underflow
+		// underflow to zero
 		{"1e-1000", 0, Below},
 		{"0x0.0000000000001p-1023", 0, Below},
 		{"0x0.00000000000008p-1022", 0, Below},
 
 		// denormals
 		{"0x0.0000000000000cp-1022", math.SmallestNonzeroFloat64, Above}, // rounded up to smallest denormal
-		{"0x0.0000000000001p-1022", math.SmallestNonzeroFloat64, Exact},  // smallest denormal
+		{"0x0.00000000000010p-1022", math.SmallestNonzeroFloat64, Exact}, // smallest denormal
 		{"0x.8p-1073", math.SmallestNonzeroFloat64, Exact},
 		{"1p-1074", math.SmallestNonzeroFloat64, Exact},
 		{"0x.fffffffffffffp-1022", math.Float64frombits(0x000fffffffffffff), Exact}, // largest denormal
 
+		// special denormal cases (see issues 14553, 14651)
+		{"0x0.00000000000001p-1022", math.Float64frombits(0x00000000000000000), Below}, // underflow to zero
+		{"0x0.00000000000004p-1022", math.Float64frombits(0x00000000000000000), Below}, // underflow to zero
+		{"0x0.00000000000008p-1022", math.Float64frombits(0x00000000000000000), Below}, // rounded down to even
+		{"0x0.00000000000009p-1022", math.Float64frombits(0x00000000000000001), Above}, // rounded up to smallest denormal
+		{"0x0.0000000000000ap-1022", math.Float64frombits(0x00000000000000001), Above}, // rounded up to smallest denormal
+
+		{"0x0.8p-1074", math.Float64frombits(0x00000000000000000), Below}, // rounded down to even
+		{"0x0.9p-1074", math.Float64frombits(0x00000000000000001), Above}, // rounded up to smallest denormal
+		{"0x0.ap-1074", math.Float64frombits(0x00000000000000001), Above}, // rounded up to smallest denormal
+		{"0x0.bp-1074", math.Float64frombits(0x00000000000000001), Above}, // rounded up to smallest denormal
+		{"0x0.cp-1074", math.Float64frombits(0x00000000000000001), Above}, // rounded up to smallest denormal
+
+		{"0x1.0p-1074", math.Float64frombits(0x00000000000000001), Exact},
+		{"0x1.7p-1074", math.Float64frombits(0x00000000000000001), Below},
+		{"0x1.8p-1074", math.Float64frombits(0x00000000000000002), Above},
+		{"0x1.9p-1074", math.Float64frombits(0x00000000000000002), Above},
+
+		{"0x2.0p-1074", math.Float64frombits(0x00000000000000002), Exact},
+		{"0x2.8p-1074", math.Float64frombits(0x00000000000000002), Below}, // rounded down to even
+		{"0x2.9p-1074", math.Float64frombits(0x00000000000000003), Above},
+
+		{"0x3.0p-1074", math.Float64frombits(0x00000000000000003), Exact},
+		{"0x3.7p-1074", math.Float64frombits(0x00000000000000003), Below},
+		{"0x3.8p-1074", math.Float64frombits(0x00000000000000004), Above}, // rounded up to even
+
+		{"0x4.0p-1074", math.Float64frombits(0x00000000000000004), Exact},
+		{"0x4.8p-1074", math.Float64frombits(0x00000000000000004), Below}, // rounded down to even
+		{"0x4.9p-1074", math.Float64frombits(0x00000000000000005), Above},
+
 		// normals
 		{"0x.fffffffffffff8p-1022", math.Float64frombits(0x0010000000000000), Above}, // rounded up to smallest normal
 		{"1p-1022", math.Float64frombits(0x0010000000000000), Exact},                 // smallest normal
@@ -958,7 +1028,7 @@ func TestFloatFloat64(t *testing.T) {
 			x := makeFloat(tx)
 			out, acc := x.Float64()
 			if !alike64(out, tout) || acc != tacc {
-				t.Errorf("%s: got %g (%#x, %s); want %g (%#x, %s)", tx, out, math.Float64bits(out), acc, test.out, math.Float64bits(test.out), tacc)
+				t.Errorf("%s: got %g (%#016x, %s); want %g (%#016x, %s)", tx, out, math.Float64bits(out), acc, test.out, math.Float64bits(test.out), tacc)
 			}
 
 			// test that x.SetFloat64(f).Float64() == f
diff --git a/libgo/go/math/big/floatconv.go b/libgo/go/math/big/floatconv.go
index 37d5c06..a884df6 100644
--- a/libgo/go/math/big/floatconv.go
+++ b/libgo/go/math/big/floatconv.go
@@ -85,7 +85,7 @@ func (z *Float) scan(r io.ByteScanner, base int) (f *Float, b int, err error) {
 	if fcount < 0 {
 		// The mantissa has a "decimal" point ddd.dddd; and
 		// -fcount is the number of digits to the right of '.'.
-		// Adjust relevant exponent accodingly.
+		// Adjust relevant exponent accordingly.
 		d := int64(fcount)
 		switch b {
 		case 10:
diff --git a/libgo/go/math/big/floatconv_test.go b/libgo/go/math/big/floatconv_test.go
index b6f9993..b2a1ab0 100644
--- a/libgo/go/math/big/floatconv_test.go
+++ b/libgo/go/math/big/floatconv_test.go
@@ -290,6 +290,11 @@ func TestFloat64Text(t *testing.T) {
 		// Issue 2625.
 		{383260575764816448, 'f', 0, "383260575764816448"},
 		{383260575764816448, 'g', -1, "3.8326057576481645e+17"},
+
+		// Issue 15918.
+		{1, 'f', -10, "1"},
+		{1, 'f', -11, "1"},
+		{1, 'f', -12, "1"},
 	} {
 		// The test cases are from the strconv package which tests float64 values.
 		// When formatting values with prec = -1 (shortest representation),
diff --git a/libgo/go/math/big/floatmarsh.go b/libgo/go/math/big/floatmarsh.go
index 44987ee..3725d4b 100644
--- a/libgo/go/math/big/floatmarsh.go
+++ b/libgo/go/math/big/floatmarsh.go
@@ -6,7 +6,94 @@
 
 package big
 
-import "fmt"
+import (
+	"encoding/binary"
+	"fmt"
+)
+
+// Gob codec version. Permits backward-compatible changes to the encoding.
+const floatGobVersion byte = 1
+
+// GobEncode implements the gob.GobEncoder interface.
+// The Float value and all its attributes (precision,
+// rounding mode, accuracy) are marshalled.
+func (x *Float) GobEncode() ([]byte, error) {
+	if x == nil {
+		return nil, nil
+	}
+
+	// determine max. space (bytes) required for encoding
+	sz := 1 + 1 + 4 // version + mode|acc|form|neg (3+2+2+1bit) + prec
+	n := 0          // number of mantissa words
+	if x.form == finite {
+		// add space for mantissa and exponent
+		n = int((x.prec + (_W - 1)) / _W) // required mantissa length in words for given precision
+		// actual mantissa slice could be shorter (trailing 0's) or longer (unused bits):
+		// - if shorter, only encode the words present
+		// - if longer, cut off unused words when encoding in bytes
+		//   (in practice, this should never happen since rounding
+		//   takes care of it, but be safe and do it always)
+		if len(x.mant) < n {
+			n = len(x.mant)
+		}
+		// len(x.mant) >= n
+		sz += 4 + n*_S // exp + mant
+	}
+	buf := make([]byte, sz)
+
+	buf[0] = floatGobVersion
+	b := byte(x.mode&7)<<5 | byte((x.acc+1)&3)<<3 | byte(x.form&3)<<1
+	if x.neg {
+		b |= 1
+	}
+	buf[1] = b
+	binary.BigEndian.PutUint32(buf[2:], x.prec)
+
+	if x.form == finite {
+		binary.BigEndian.PutUint32(buf[6:], uint32(x.exp))
+		x.mant[len(x.mant)-n:].bytes(buf[10:]) // cut off unused trailing words
+	}
+
+	return buf, nil
+}
+
+// GobDecode implements the gob.GobDecoder interface.
+// The result is rounded per the precision and rounding mode of
+// z unless z's precision is 0, in which case z is set exactly
+// to the decoded value.
+func (z *Float) GobDecode(buf []byte) error {
+	if len(buf) == 0 {
+		// Other side sent a nil or default value.
+		*z = Float{}
+		return nil
+	}
+
+	if buf[0] != floatGobVersion {
+		return fmt.Errorf("Float.GobDecode: encoding version %d not supported", buf[0])
+	}
+
+	oldPrec := z.prec
+	oldMode := z.mode
+
+	b := buf[1]
+	z.mode = RoundingMode((b >> 5) & 7)
+	z.acc = Accuracy((b>>3)&3) - 1
+	z.form = form((b >> 1) & 3)
+	z.neg = b&1 != 0
+	z.prec = binary.BigEndian.Uint32(buf[2:])
+
+	if z.form == finite {
+		z.exp = int32(binary.BigEndian.Uint32(buf[6:]))
+		z.mant = z.mant.setBytes(buf[10:])
+	}
+
+	if oldPrec != 0 {
+		z.mode = oldMode
+		z.SetPrec(uint(oldPrec))
+	}
+
+	return nil
+}
 
 // MarshalText implements the encoding.TextMarshaler interface.
 // Only the Float value is marshaled (in full precision), other
diff --git a/libgo/go/math/big/floatmarsh_test.go b/libgo/go/math/big/floatmarsh_test.go
index d7ef2fc..5bd906d 100644
--- a/libgo/go/math/big/floatmarsh_test.go
+++ b/libgo/go/math/big/floatmarsh_test.go
@@ -5,7 +5,10 @@
 package big
 
 import (
+	"bytes"
+	"encoding/gob"
 	"encoding/json"
+	"io"
 	"testing"
 )
 
@@ -23,6 +26,85 @@ var floatVals = []string{
 	"Inf",
 }
 
+func TestFloatGobEncoding(t *testing.T) {
+	var medium bytes.Buffer
+	enc := gob.NewEncoder(&medium)
+	dec := gob.NewDecoder(&medium)
+	for _, test := range floatVals {
+		for _, sign := range []string{"", "+", "-"} {
+			for _, prec := range []uint{0, 1, 2, 10, 53, 64, 100, 1000} {
+				for _, mode := range []RoundingMode{ToNearestEven, ToNearestAway, ToZero, AwayFromZero, ToNegativeInf, ToPositiveInf} {
+					medium.Reset() // empty buffer for each test case (in case of failures)
+					x := sign + test
+
+					var tx Float
+					_, _, err := tx.SetPrec(prec).SetMode(mode).Parse(x, 0)
+					if err != nil {
+						t.Errorf("parsing of %s (%dbits, %v) failed (invalid test case): %v", x, prec, mode, err)
+						continue
+					}
+
+					// If tx was set to prec == 0, tx.Parse(x, 0) assumes precision 64. Correct it.
+					if prec == 0 {
+						tx.SetPrec(0)
+					}
+
+					if err := enc.Encode(&tx); err != nil {
+						t.Errorf("encoding of %v (%dbits, %v) failed: %v", &tx, prec, mode, err)
+						continue
+					}
+
+					var rx Float
+					if err := dec.Decode(&rx); err != nil {
+						t.Errorf("decoding of %v (%dbits, %v) failed: %v", &tx, prec, mode, err)
+						continue
+					}
+
+					if rx.Cmp(&tx) != 0 {
+						t.Errorf("transmission of %s failed: got %s want %s", x, rx.String(), tx.String())
+						continue
+					}
+
+					if rx.Prec() != prec {
+						t.Errorf("transmission of %s's prec failed: got %d want %d", x, rx.Prec(), prec)
+					}
+
+					if rx.Mode() != mode {
+						t.Errorf("transmission of %s's mode failed: got %s want %s", x, rx.Mode(), mode)
+					}
+
+					if rx.Acc() != tx.Acc() {
+						t.Errorf("transmission of %s's accuracy failed: got %s want %s", x, rx.Acc(), tx.Acc())
+					}
+				}
+			}
+		}
+	}
+}
+
+func TestFloatCorruptGob(t *testing.T) {
+	var buf bytes.Buffer
+	tx := NewFloat(4 / 3).SetPrec(1000).SetMode(ToPositiveInf)
+	if err := gob.NewEncoder(&buf).Encode(tx); err != nil {
+		t.Fatal(err)
+	}
+	b := buf.Bytes()
+
+	var rx Float
+	if err := gob.NewDecoder(bytes.NewReader(b)).Decode(&rx); err != nil {
+		t.Fatal(err)
+	}
+
+	if err := gob.NewDecoder(bytes.NewReader(b[:10])).Decode(&rx); err != io.ErrUnexpectedEOF {
+		t.Errorf("got %v want EOF", err)
+	}
+
+	b[1] = 0
+	if err := gob.NewDecoder(bytes.NewReader(b)).Decode(&rx); err == nil {
+		t.Fatal("got nil want version error")
+	}
+}
+
 func TestFloatJSONEncoding(t *testing.T) {
 	for _, test := range floatVals {
 		for _, sign := range []string{"", "+", "-"} {
diff --git a/libgo/go/math/big/ftoa.go b/libgo/go/math/big/ftoa.go
index c5cdb5e..57b16e1 100644
--- a/libgo/go/math/big/ftoa.go
+++ b/libgo/go/math/big/ftoa.go
@@ -41,8 +41,11 @@ import (
 // x.Prec() mantissa bits.
 // The prec value is ignored for the 'b' or 'p' format.
 func (x *Float) Text(format byte, prec int) string {
-	const extra = 10 // TODO(gri) determine a good/better value here
-	return string(x.Append(make([]byte, 0, prec+extra), format, prec))
+	cap := 10 // TODO(gri) determine a good/better value here
+	if prec > 0 {
+		cap += prec
+	}
+	return string(x.Append(make([]byte, 0, cap), format, prec))
 }
 
 // String formats x like x.Text('g', 10).
@@ -333,9 +336,9 @@ func (x *Float) fmtB(buf []byte) []byte {
 	return strconv.AppendInt(buf, e, 10)
 }
 
-// fmtP appends the string of x in the format 0x." mantissa "p" exponent
-// with a hexadecimal mantissa and a binary exponent, or 0" if x is zero,
-// ad returns the extended buffer.
+// fmtP appends the string of x in the format "0x." mantissa "p" exponent
+// with a hexadecimal mantissa and a binary exponent, or "0" if x is zero,
+// and returns the extended buffer.
 // The mantissa is normalized such that 0.5 <= 0.mantissa < 1.0.
 // The sign of x is ignored, and x must not be an Inf.
 func (x *Float) fmtP(buf []byte) []byte {
@@ -374,12 +377,11 @@ func min(x, y int) int {
 }
 
 // Format implements fmt.Formatter. It accepts all the regular
-// formats for floating-point numbers ('e', 'E', 'f', 'F', 'g',
-// 'G') as well as 'b', 'p', and 'v'. See (*Float).Text for the
-// interpretation of 'b' and 'p'. The 'v' format is handled like
-// 'g'.
+// formats for floating-point numbers ('b', 'e', 'E', 'f', 'F',
+// 'g', 'G') as well as 'p' and 'v'. See (*Float).Text for the
+// interpretation of 'p'. The 'v' format is handled like 'g'.
 // Format also supports specification of the minimum precision
-// in digits, the output field width, as well as the format verbs
+// in digits, the output field width, as well as the format flags
 // '+' and ' ' for sign control, '0' for space or zero padding,
 // and '-' for left or right justification. See the fmt package
 // for details.
diff --git a/libgo/go/math/big/gcd_test.go b/libgo/go/math/big/gcd_test.go
index c0b9f58..a929bf5 100644
--- a/libgo/go/math/big/gcd_test.go
+++ b/libgo/go/math/big/gcd_test.go
@@ -20,13 +20,27 @@ func randInt(r *rand.Rand, size uint) *Int {
 }
 
 func runGCD(b *testing.B, aSize, bSize uint) {
+	b.Run("WithoutXY", func(b *testing.B) {
+		runGCDExt(b, aSize, bSize, false)
+	})
+	b.Run("WithXY", func(b *testing.B) {
+		runGCDExt(b, aSize, bSize, true)
+	})
+}
+
+func runGCDExt(b *testing.B, aSize, bSize uint, calcXY bool) {
 	b.StopTimer()
 	var r = rand.New(rand.NewSource(1234))
 	aa := randInt(r, aSize)
 	bb := randInt(r, bSize)
+	var x, y *Int
+	if calcXY {
+		x = new(Int)
+		y = new(Int)
+	}
 	b.StartTimer()
 	for i := 0; i < b.N; i++ {
-		new(Int).GCD(nil, nil, aa, bb)
+		new(Int).GCD(x, y, aa, bb)
 	}
 }
 
diff --git a/libgo/go/math/big/int.go b/libgo/go/math/big/int.go
index 67ab704..f2a75d1 100644
--- a/libgo/go/math/big/int.go
+++ b/libgo/go/math/big/int.go
@@ -459,11 +459,11 @@ func (z *Int) GCD(x, y, a, b *Int) *Int {
 	q := new(Int)
 	temp := new(Int)
 
+	r := new(Int)
 	for len(B.abs) > 0 {
-		r := new(Int)
 		q, r = q.QuoRem(A, B, r)
 
-		A, B = B, r
+		A, B, r = B, r, A
 
 		temp.Set(X)
 		X.Mul(X, q)
diff --git a/libgo/go/math/big/intconv.go b/libgo/go/math/big/intconv.go
index 56a75f8..daf674ae 100644
--- a/libgo/go/math/big/intconv.go
+++ b/libgo/go/math/big/intconv.go
@@ -52,16 +52,16 @@ func writeMultiple(s fmt.State, text string, count int) {
 	}
 }
 
-// Format is a support routine for fmt.Formatter. It accepts
-// the formats 'b' (binary), 'o' (octal), 'd' (decimal), 'x'
-// (lowercase hexadecimal), and 'X' (uppercase hexadecimal).
+// Format implements fmt.Formatter. It accepts the formats
+// 'b' (binary), 'o' (octal), 'd' (decimal), 'x' (lowercase
+// hexadecimal), and 'X' (uppercase hexadecimal).
 // Also supported are the full suite of package fmt's format
-// verbs for integral types, including '+', '-', and ' '
-// for sign control, '#' for leading zero in octal and for
-// hexadecimal, a leading "0x" or "0X" for "%#x" and "%#X"
-// respectively, specification of minimum digits precision,
-// output field width, space or zero padding, and left or
-// right justification.
+// flags for integral types, including '+' and ' ' for sign
+// control, '#' for leading zero in octal and for hexadecimal,
+// a leading "0x" or "0X" for "%#x" and "%#X" respectively,
+// specification of minimum digits precision, output field
+// width, space or zero padding, and '-' for left or right
+// justification.
 //
 func (x *Int) Format(s fmt.State, ch rune) {
 	// determine base
diff --git a/libgo/go/math/big/nat.go b/libgo/go/math/big/nat.go
index 79cf6e0..2e65d2a 100644
--- a/libgo/go/math/big/nat.go
+++ b/libgo/go/math/big/nat.go
@@ -8,7 +8,10 @@
 
 package big
 
-import "math/rand"
+import (
+	"math/rand"
+	"sync"
+)
 
 // An unsigned integer x of the form
 //
@@ -539,6 +542,21 @@ func (z nat) div(z2, u, v nat) (q, r nat) {
 	return
 }
 
+// getNat returns a nat of len n. The contents may not be zero.
+func getNat(n int) nat {
+	var z nat
+	if v := natPool.Get(); v != nil {
+		z = v.(nat)
+	}
+	return z.make(n)
+}
+
+func putNat(x nat) {
+	natPool.Put(x)
+}
+
+var natPool sync.Pool
+
 // q = (uIn-r)/v, with 0 <= r < y
 // Uses z as storage for q, and u as storage for r if possible.
 // See Knuth, Volume 2, section 4.3.1, Algorithm D.
@@ -557,7 +575,7 @@ func (z nat) divLarge(u, uIn, v nat) (q, r nat) {
 	}
 	q = z.make(m + 1)
 
-	qhatv := make(nat, n+1)
+	qhatv := getNat(n + 1)
 	if alias(u, uIn) || alias(u, v) {
 		u = nil // u is an alias for uIn or v - cannot reuse
 	}
@@ -565,10 +583,11 @@ func (z nat) divLarge(u, uIn, v nat) (q, r nat) {
 	u.clear() // TODO(gri) no need to clear if we allocated a new u
 
 	// D1.
+	var v1 nat
 	shift := nlz(v[n-1])
 	if shift > 0 {
 		// do not modify v, it may be used by another goroutine simultaneously
-		v1 := make(nat, n)
+		v1 = getNat(n)
 		shlVU(v1, v, shift)
 		v = v1
 	}
@@ -609,6 +628,10 @@ func (z nat) divLarge(u, uIn, v nat) (q, r nat) {
 
 		q[j] = qhat
 	}
+	if v1 != nil {
+		putNat(v1)
+	}
+	putNat(qhatv)
 
 	q = q.norm()
 	shrVU(u, u, shift)
@@ -647,7 +670,7 @@ func trailingZeroBits(x Word) uint {
 	// x & -x leaves only the right-most bit set in the word. Let k be the
 	// index of that bit. Since only a single bit is set, the value is two
 	// to the power of k. Multiplying by a power of two is equivalent to
-	// left shifting, in this case by k bits.  The de Bruijn constant is
+	// left shifting, in this case by k bits. The de Bruijn constant is
 	// such that all six bit, consecutive substrings are distinct.
 	// Therefore, if we have a left shifted version of this constant we can
 	// find by how many bits it was shifted by looking at which six bit
@@ -1018,7 +1041,7 @@ func (z nat) expNNWindowed(x, y, m nat) nat {
 		for j := 0; j < _W; j += n {
 			if i != len(y)-1 || j != 0 {
 				// Unrolled loop for significant performance
-				// gain.  Use go test -bench=".*" in crypto/rsa
+				// gain. Use go test -bench=".*" in crypto/rsa
 				// to check performance before making changes.
 				zz = zz.mul(z, z)
 				zz, z = z, zz
diff --git a/libgo/go/math/big/nat_test.go b/libgo/go/math/big/nat_test.go
index 563ccb3..ebb2985 100644
--- a/libgo/go/math/big/nat_test.go
+++ b/libgo/go/math/big/nat_test.go
@@ -5,6 +5,7 @@
 package big
 
 import (
+	"fmt"
 	"runtime"
 	"strings"
 	"testing"
@@ -509,24 +510,20 @@ func TestExpNN(t *testing.T) {
 	}
 }
 
-func ExpHelper(b *testing.B, x, y Word) {
-	var z nat
-	for i := 0; i < b.N; i++ {
-		z.expWW(x, y)
+func BenchmarkExp3Power(b *testing.B) {
+	const x = 3
+	for _, y := range []Word{
+		0x10, 0x40, 0x100, 0x400, 0x1000, 0x4000, 0x10000, 0x40000, 0x100000, 0x400000,
+	} {
+		b.Run(fmt.Sprintf("%#x", y), func(b *testing.B) {
+			var z nat
+			for i := 0; i < b.N; i++ {
+				z.expWW(x, y)
+			}
+		})
 	}
 }
 
-func BenchmarkExp3Power0x10(b *testing.B)     { ExpHelper(b, 3, 0x10) }
-func BenchmarkExp3Power0x40(b *testing.B)     { ExpHelper(b, 3, 0x40) }
-func BenchmarkExp3Power0x100(b *testing.B)    { ExpHelper(b, 3, 0x100) }
-func BenchmarkExp3Power0x400(b *testing.B)    { ExpHelper(b, 3, 0x400) }
-func BenchmarkExp3Power0x1000(b *testing.B)   { ExpHelper(b, 3, 0x1000) }
-func BenchmarkExp3Power0x4000(b *testing.B)   { ExpHelper(b, 3, 0x4000) }
-func BenchmarkExp3Power0x10000(b *testing.B)  { ExpHelper(b, 3, 0x10000) }
-func BenchmarkExp3Power0x40000(b *testing.B)  { ExpHelper(b, 3, 0x40000) }
-func BenchmarkExp3Power0x100000(b *testing.B) { ExpHelper(b, 3, 0x100000) }
-func BenchmarkExp3Power0x400000(b *testing.B) { ExpHelper(b, 3, 0x400000) }
-
 func fibo(n int) nat {
 	switch n {
 	case 0:
diff --git a/libgo/go/math/big/natconv.go b/libgo/go/math/big/natconv.go
index d2ce667..4454784 100644
--- a/libgo/go/math/big/natconv.go
+++ b/libgo/go/math/big/natconv.go
@@ -302,7 +302,7 @@ func (x nat) itoa(neg bool, base int) []byte {
 		}
 
 	} else {
-		bb, ndigits := maxPow(Word(b))
+		bb, ndigits := maxPow(b)
 
 		// construct table of successive squares of bb*leafSize to use in subdivisions
 		// result (table != nil) <=> (len(x) > leafSize > 0)
@@ -391,7 +391,7 @@ func (q nat) convertWords(s []byte, b Word, ndigits int, bb Word, table []diviso
 				// this appears to be faster for BenchmarkString10000Base10
 				// and smaller strings (but a bit slower for larger ones)
 				t := r / 10
-				s[i] = '0' + byte(r-t<<3-t-t) // TODO(gri) replace w/ t*10 once compiler produces better code
+				s[i] = '0' + byte(r-t*10)
 				r = t
 			}
 		}
diff --git a/libgo/go/math/big/natconv_test.go b/libgo/go/math/big/natconv_test.go
index 028e5a8..79901d1 100644
--- a/libgo/go/math/big/natconv_test.go
+++ b/libgo/go/math/big/natconv_test.go
@@ -6,6 +6,7 @@ package big
 
 import (
 	"bytes"
+	"fmt"
 	"io"
 	"strings"
 	"testing"
@@ -273,102 +274,58 @@ func BenchmarkStringPiParallel(b *testing.B) {
 	})
 }
 
-func BenchmarkScan10Base2(b *testing.B)     { ScanHelper(b, 2, 10, 10) }
-func BenchmarkScan100Base2(b *testing.B)    { ScanHelper(b, 2, 10, 100) }
-func BenchmarkScan1000Base2(b *testing.B)   { ScanHelper(b, 2, 10, 1000) }
-func BenchmarkScan10000Base2(b *testing.B)  { ScanHelper(b, 2, 10, 10000) }
-func BenchmarkScan100000Base2(b *testing.B) { ScanHelper(b, 2, 10, 100000) }
-
-func BenchmarkScan10Base8(b *testing.B)     { ScanHelper(b, 8, 10, 10) }
-func BenchmarkScan100Base8(b *testing.B)    { ScanHelper(b, 8, 10, 100) }
-func BenchmarkScan1000Base8(b *testing.B)   { ScanHelper(b, 8, 10, 1000) }
-func BenchmarkScan10000Base8(b *testing.B)  { ScanHelper(b, 8, 10, 10000) }
-func BenchmarkScan100000Base8(b *testing.B) { ScanHelper(b, 8, 10, 100000) }
-
-func BenchmarkScan10Base10(b *testing.B)     { ScanHelper(b, 10, 10, 10) }
-func BenchmarkScan100Base10(b *testing.B)    { ScanHelper(b, 10, 10, 100) }
-func BenchmarkScan1000Base10(b *testing.B)   { ScanHelper(b, 10, 10, 1000) }
-func BenchmarkScan10000Base10(b *testing.B)  { ScanHelper(b, 10, 10, 10000) }
-func BenchmarkScan100000Base10(b *testing.B) { ScanHelper(b, 10, 10, 100000) }
-
-func BenchmarkScan10Base16(b *testing.B)     { ScanHelper(b, 16, 10, 10) }
-func BenchmarkScan100Base16(b *testing.B)    { ScanHelper(b, 16, 10, 100) }
-func BenchmarkScan1000Base16(b *testing.B)   { ScanHelper(b, 16, 10, 1000) }
-func BenchmarkScan10000Base16(b *testing.B)  { ScanHelper(b, 16, 10, 10000) }
-func BenchmarkScan100000Base16(b *testing.B) { ScanHelper(b, 16, 10, 100000) }
-
-func ScanHelper(b *testing.B, base int, x, y Word) {
-	b.StopTimer()
-	var z nat
-	z = z.expWW(x, y)
-
-	s := z.utoa(base)
-	if t := itoa(z, base); !bytes.Equal(s, t) {
-		b.Fatalf("scanning: got %s; want %s", s, t)
+func BenchmarkScan(b *testing.B) {
+	const x = 10
+	for _, base := range []int{2, 8, 10, 16} {
+		for _, y := range []Word{10, 100, 1000, 10000, 100000} {
+			b.Run(fmt.Sprintf("%d/Base%d", y, base), func(b *testing.B) {
+				b.StopTimer()
+				var z nat
+				z = z.expWW(x, y)
+
+				s := z.utoa(base)
+				if t := itoa(z, base); !bytes.Equal(s, t) {
+					b.Fatalf("scanning: got %s; want %s", s, t)
+				}
+				b.StartTimer()
+
+				for i := 0; i < b.N; i++ {
+					z.scan(bytes.NewReader(s), base, false)
+				}
+			})
+		}
 	}
-	b.StartTimer()
+}
 
-	for i := 0; i < b.N; i++ {
-		z.scan(bytes.NewReader(s), base, false)
+func BenchmarkString(b *testing.B) {
+	const x = 10
+	for _, base := range []int{2, 8, 10, 16} {
+		for _, y := range []Word{10, 100, 1000, 10000, 100000} {
+			b.Run(fmt.Sprintf("%d/Base%d", y, base), func(b *testing.B) {
+				b.StopTimer()
+				var z nat
+				z = z.expWW(x, y)
+				z.utoa(base) // warm divisor cache
+				b.StartTimer()
+
+				for i := 0; i < b.N; i++ {
+					_ = z.utoa(base)
+				}
+			})
+		}
 	}
 }
 
-func BenchmarkString10Base2(b *testing.B)     { StringHelper(b, 2, 10, 10) }
-func BenchmarkString100Base2(b *testing.B)    { StringHelper(b, 2, 10, 100) }
-func BenchmarkString1000Base2(b *testing.B)   { StringHelper(b, 2, 10, 1000) }
-func BenchmarkString10000Base2(b *testing.B)  { StringHelper(b, 2, 10, 10000) }
-func BenchmarkString100000Base2(b *testing.B) { StringHelper(b, 2, 10, 100000) }
-
-func BenchmarkString10Base8(b *testing.B)     { StringHelper(b, 8, 10, 10) }
-func BenchmarkString100Base8(b *testing.B)    { StringHelper(b, 8, 10, 100) }
-func BenchmarkString1000Base8(b *testing.B)   { StringHelper(b, 8, 10, 1000) }
-func BenchmarkString10000Base8(b *testing.B)  { StringHelper(b, 8, 10, 10000) }
-func BenchmarkString100000Base8(b *testing.B) { StringHelper(b, 8, 10, 100000) }
-
-func BenchmarkString10Base10(b *testing.B)     { StringHelper(b, 10, 10, 10) }
-func BenchmarkString100Base10(b *testing.B)    { StringHelper(b, 10, 10, 100) }
-func BenchmarkString1000Base10(b *testing.B)   { StringHelper(b, 10, 10, 1000) }
-func BenchmarkString10000Base10(b *testing.B)  { StringHelper(b, 10, 10, 10000) }
-func BenchmarkString100000Base10(b *testing.B) { StringHelper(b, 10, 10, 100000) }
-
-func BenchmarkString10Base16(b *testing.B)     { StringHelper(b, 16, 10, 10) }
-func BenchmarkString100Base16(b *testing.B)    { StringHelper(b, 16, 10, 100) }
-func BenchmarkString1000Base16(b *testing.B)   { StringHelper(b, 16, 10, 1000) }
-func BenchmarkString10000Base16(b *testing.B)  { StringHelper(b, 16, 10, 10000) }
-func BenchmarkString100000Base16(b *testing.B) { StringHelper(b, 16, 10, 100000) }
-
-func StringHelper(b *testing.B, base int, x, y Word) {
-	b.StopTimer()
-	var z nat
-	z = z.expWW(x, y)
-	z.utoa(base) // warm divisor cache
-	b.StartTimer()
-
-	for i := 0; i < b.N; i++ {
-		_ = z.utoa(base)
+func BenchmarkLeafSize(b *testing.B) {
+	for n := 0; n <= 16; n++ {
+		b.Run(fmt.Sprint(n), func(b *testing.B) { LeafSizeHelper(b, 10, n) })
+	}
+	// Try some large lengths
+	for _, n := range []int{32, 64} {
+		b.Run(fmt.Sprint(n), func(b *testing.B) { LeafSizeHelper(b, 10, n) })
 	}
 }
 
-func BenchmarkLeafSize0(b *testing.B)  { LeafSizeHelper(b, 10, 0) } // test without splitting
-func BenchmarkLeafSize1(b *testing.B)  { LeafSizeHelper(b, 10, 1) }
-func BenchmarkLeafSize2(b *testing.B)  { LeafSizeHelper(b, 10, 2) }
-func BenchmarkLeafSize3(b *testing.B)  { LeafSizeHelper(b, 10, 3) }
-func BenchmarkLeafSize4(b *testing.B)  { LeafSizeHelper(b, 10, 4) }
-func BenchmarkLeafSize5(b *testing.B)  { LeafSizeHelper(b, 10, 5) }
-func BenchmarkLeafSize6(b *testing.B)  { LeafSizeHelper(b, 10, 6) }
-func BenchmarkLeafSize7(b *testing.B)  { LeafSizeHelper(b, 10, 7) }
-func BenchmarkLeafSize8(b *testing.B)  { LeafSizeHelper(b, 10, 8) }
-func BenchmarkLeafSize9(b *testing.B)  { LeafSizeHelper(b, 10, 9) }
-func BenchmarkLeafSize10(b *testing.B) { LeafSizeHelper(b, 10, 10) }
-func BenchmarkLeafSize11(b *testing.B) { LeafSizeHelper(b, 10, 11) }
-func BenchmarkLeafSize12(b *testing.B) { LeafSizeHelper(b, 10, 12) }
-func BenchmarkLeafSize13(b *testing.B) { LeafSizeHelper(b, 10, 13) }
-func BenchmarkLeafSize14(b *testing.B) { LeafSizeHelper(b, 10, 14) }
-func BenchmarkLeafSize15(b *testing.B) { LeafSizeHelper(b, 10, 15) }
-func BenchmarkLeafSize16(b *testing.B) { LeafSizeHelper(b, 10, 16) }
-func BenchmarkLeafSize32(b *testing.B) { LeafSizeHelper(b, 10, 32) } // try some large lengths
-func BenchmarkLeafSize64(b *testing.B) { LeafSizeHelper(b, 10, 64) }
-
 func LeafSizeHelper(b *testing.B, base, size int) {
 	b.StopTimer()
 	originalLeafSize := leafSize
diff --git a/libgo/go/math/big/rat.go b/libgo/go/math/big/rat.go
index 2cd9ed0..56ce33d 100644
--- a/libgo/go/math/big/rat.go
+++ b/libgo/go/math/big/rat.go
@@ -63,7 +63,7 @@ func (z *Rat) SetFloat64(f float64) *Rat {
 
 // quotToFloat32 returns the non-negative float32 value
 // nearest to the quotient a/b, using round-to-even in
-// halfway cases.  It does not mutate its arguments.
+// halfway cases. It does not mutate its arguments.
 // Preconditions: b is non-zero; a and b have no common factors.
 func quotToFloat32(a, b nat) (f float32, exact bool) {
 	const (
@@ -161,7 +161,7 @@ func quotToFloat32(a, b nat) (f float32, exact bool) {
 
 // quotToFloat64 returns the non-negative float64 value
 // nearest to the quotient a/b, using round-to-even in
-// halfway cases.  It does not mutate its arguments.
+// halfway cases. It does not mutate its arguments.
 // Preconditions: b is non-zero; a and b have no common factors.
 func quotToFloat64(a, b nat) (f float64, exact bool) {
 	const (
diff --git a/libgo/go/math/big/ratconv.go b/libgo/go/math/big/ratconv.go
index 4566ff4..ef2b675 100644
--- a/libgo/go/math/big/ratconv.go
+++ b/libgo/go/math/big/ratconv.go
@@ -15,7 +15,7 @@ import (
 )
 
 func ratTok(ch rune) bool {
-	return strings.IndexRune("+-/0123456789.eE", ch) >= 0
+	return strings.ContainsRune("+-/0123456789.eE", ch)
 }
 
 // Scan is a support routine for fmt.Scanner. It accepts the formats
@@ -25,7 +25,7 @@ func (z *Rat) Scan(s fmt.ScanState, ch rune) error {
 	if err != nil {
 		return err
 	}
-	if strings.IndexRune("efgEFGv", ch) < 0 {
+	if !strings.ContainsRune("efgEFGv", ch) {
 		return errors.New("Rat.Scan: invalid verb")
 	}
 	if _, ok := z.SetString(string(tok)); !ok {
@@ -88,6 +88,12 @@ func (z *Rat) SetString(s string) (*Rat, bool) {
 		return nil, false
 	}
 
+	// special-case 0 (see also issue #16176)
+	if len(z.a.abs) == 0 {
+		return z, true
+	}
+	// len(z.a.abs) > 0
+
 	// correct exponent
 	if ecorr < 0 {
 		exp += int64(ecorr)
@@ -178,7 +184,7 @@ func scanExponent(r io.ByteScanner, binExpOk bool) (exp int64, base int, err err
 			}
 			break // i > 0
 		}
-		digits = append(digits, byte(ch))
+		digits = append(digits, ch)
 	}
 	// i > 0 => we have at least one digit
 
diff --git a/libgo/go/math/big/ratconv_test.go b/libgo/go/math/big/ratconv_test.go
index da2fdab..35ad6cc 100644
--- a/libgo/go/math/big/ratconv_test.go
+++ b/libgo/go/math/big/ratconv_test.go
@@ -48,6 +48,7 @@ var setStringTests = []StringTest{
 	{"53/70893980658822810696", "53/70893980658822810696", true},
 	{"106/141787961317645621392", "53/70893980658822810696", true},
 	{"204211327800791583.81095", "4084226556015831676219/20000", true},
+	{"0e9999999999", "0", true}, // issue #16176
 	{in: "1/0"},
 }
 
@@ -137,7 +138,7 @@ func TestFloatString(t *testing.T) {
 	}
 }
 
-// Test inputs to Rat.SetString.  The prefix "long:" causes the test
+// Test inputs to Rat.SetString. The prefix "long:" causes the test
 // to be skipped in --test.short mode.  (The threshold is about 500us.)
 var float64inputs = []string{
 	// Constants plundered from strconv/testfp.txt.
diff --git a/libgo/go/math/cmplx/cmath_test.go b/libgo/go/math/cmplx/cmath_test.go
index 18d9be8..d904be8 100644
--- a/libgo/go/math/cmplx/cmath_test.go
+++ b/libgo/go/math/cmplx/cmath_test.go
@@ -9,6 +9,9 @@ import (
 	"testing"
 )
 
+// The higher-precision values in vc26 were used to derive the
+// input arguments vc (see also comment below). For reference
+// only (do not delete).
 var vc26 = []complex128{
 	(4.97901192488367350108546816 + 7.73887247457810456552351752i),
 	(7.73887247457810456552351752 - 0.27688005719200159404635997i),
@@ -21,6 +24,7 @@ var vc26 = []complex128{
 	(1.82530809168085506044576505 - 8.68592476857560136238589621i),
 	(-8.68592476857560136238589621 + 4.97901192488367350108546816i),
 }
+
 var vc = []complex128{
 	(4.9790119248836735e+00 + 7.7388724745781045e+00i),
 	(7.7388724745781045e+00 - 2.7688005719200159e-01i),
@@ -448,8 +452,7 @@ func tolerance(a, b, e float64) bool {
 	}
 	return d < e
 }
-func soclose(a, b, e float64) bool { return tolerance(a, b, e) }
-func veryclose(a, b float64) bool  { return tolerance(a, b, 4e-16) }
+func veryclose(a, b float64) bool { return tolerance(a, b, 4e-16) }
 func alike(a, b float64) bool {
 	switch {
 	case a != a && b != b: // math.IsNaN(a) && math.IsNaN(b):
diff --git a/libgo/go/math/cmplx/sqrt.go b/libgo/go/math/cmplx/sqrt.go
index 276be07..72f81e9 100644
--- a/libgo/go/math/cmplx/sqrt.go
+++ b/libgo/go/math/cmplx/sqrt.go
@@ -43,7 +43,7 @@ import "math"
 // Cancelation error in r-x or r+x is avoided by using the
 // identity  2 Re w Im w  =  y.
 //
-// Note that -w is also a square root of z.  The root chosen
+// Note that -w is also a square root of z. The root chosen
 // is always in the right half plane and Im w has the same sign as y.
 //
 // ACCURACY:
diff --git a/libgo/go/math/erf.go b/libgo/go/math/erf.go
index 4cd80f8..8ddd5f9 100644
--- a/libgo/go/math/erf.go
+++ b/libgo/go/math/erf.go
@@ -10,7 +10,7 @@ package math
 
 // The original C code and the long comment below are
 // from FreeBSD's /usr/src/lib/msun/src/s_erf.c and
-// came with this notice.  The go code is a simplified
+// came with this notice. The go code is a simplified
 // version of the original C.
 //
 // ====================================================
diff --git a/libgo/go/math/exp.go b/libgo/go/math/exp.go
index 51330c2..cbc955d 100644
--- a/libgo/go/math/exp.go
+++ b/libgo/go/math/exp.go
@@ -21,7 +21,7 @@ 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
+// and came with this notice. The go code is a simplified
 // version of the original C.
 //
 // ====================================================
diff --git a/libgo/go/math/expm1.go b/libgo/go/math/expm1.go
index b4dcdc5..5a96218 100644
--- a/libgo/go/math/expm1.go
+++ b/libgo/go/math/expm1.go
@@ -6,7 +6,7 @@ package math
 
 // The original C code, the long comment, and the constants
 // below are from FreeBSD's /usr/src/lib/msun/src/s_expm1.c
-// and came with this notice.  The go code is a simplified
+// and came with this notice. The go code is a simplified
 // version of the original C.
 //
 // ====================================================
diff --git a/libgo/go/math/export_test.go b/libgo/go/math/export_test.go
index 02992d7..368308e 100644
--- a/libgo/go/math/export_test.go
+++ b/libgo/go/math/export_test.go
@@ -1,4 +1,4 @@
-// Copyright 2011 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.
 
diff --git a/libgo/go/math/gamma.go b/libgo/go/math/gamma.go
index 164f54f..841ec11 100644
--- a/libgo/go/math/gamma.go
+++ b/libgo/go/math/gamma.go
@@ -21,7 +21,7 @@ package math
 //
 // DESCRIPTION:
 //
-// Returns gamma function of the argument.  The result is
+// Returns gamma function of the argument. The result is
 // correctly signed, and the sign (+1 or -1) is also
 // returned in a global (extern) variable named signgam.
 // This variable is also filled in by the logarithmic gamma
diff --git a/libgo/go/math/j0.go b/libgo/go/math/j0.go
index de77388..cbef7aa 100644
--- a/libgo/go/math/j0.go
+++ b/libgo/go/math/j0.go
@@ -10,7 +10,7 @@ package math
 
 // The original C code and the long comment below are
 // from FreeBSD's /usr/src/lib/msun/src/e_j0.c and
-// came with this notice.  The go code is a simplified
+// came with this notice. The go code is a simplified
 // version of the original C.
 //
 // ====================================================
diff --git a/libgo/go/math/j1.go b/libgo/go/math/j1.go
index c537a72..d359d90 100644
--- a/libgo/go/math/j1.go
+++ b/libgo/go/math/j1.go
@@ -10,7 +10,7 @@ package math
 
 // The original C code and the long comment below are
 // from FreeBSD's /usr/src/lib/msun/src/e_j1.c and
-// came with this notice.  The go code is a simplified
+// came with this notice. The go code is a simplified
 // version of the original C.
 //
 // ====================================================
diff --git a/libgo/go/math/jn.go b/libgo/go/math/jn.go
index ffb8a00..721112f 100644
--- a/libgo/go/math/jn.go
+++ b/libgo/go/math/jn.go
@@ -10,7 +10,7 @@ package math
 
 // The original C code and the long comment below are
 // from FreeBSD's /usr/src/lib/msun/src/e_jn.c and
-// came with this notice.  The go code is a simplified
+// came with this notice. The go code is a simplified
 // version of the original C.
 //
 // ====================================================
diff --git a/libgo/go/math/lgamma.go b/libgo/go/math/lgamma.go
index 6a02c41..19ac3ff 100644
--- a/libgo/go/math/lgamma.go
+++ b/libgo/go/math/lgamma.go
@@ -10,7 +10,7 @@ package math
 
 // The original C code and the long comment below are
 // from FreeBSD's /usr/src/lib/msun/src/e_lgamma_r.c and
-// came with this notice.  The go code is a simplified
+// came with this notice. The go code is a simplified
 // version of the original C.
 //
 // ====================================================
diff --git a/libgo/go/math/log.go b/libgo/go/math/log.go
index 60b5755..cf242e8 100644
--- a/libgo/go/math/log.go
+++ b/libgo/go/math/log.go
@@ -10,7 +10,7 @@ package math
 
 // The original C code, the long comment, and the constants
 // below are from FreeBSD's /usr/src/lib/msun/src/e_log.c
-// and came with this notice.  The go code is a simpler
+// and came with this notice. The go code is a simpler
 // version of the original C.
 //
 // ====================================================
diff --git a/libgo/go/math/log1p.go b/libgo/go/math/log1p.go
index c8daaaa..6834cae 100644
--- a/libgo/go/math/log1p.go
+++ b/libgo/go/math/log1p.go
@@ -6,7 +6,7 @@ package math
 
 // The original C code, the long comment, and the constants
 // below are from FreeBSD's /usr/src/lib/msun/src/s_log1p.c
-// and came with this notice.  The go code is a simplified
+// and came with this notice. The go code is a simplified
 // version of the original C.
 //
 // ====================================================
diff --git a/libgo/go/math/modf.go b/libgo/go/math/modf.go
index 5d2f489..b59d4b7 100644
--- a/libgo/go/math/modf.go
+++ b/libgo/go/math/modf.go
@@ -5,7 +5,7 @@
 package math
 
 // Modf returns integer and fractional floating-point numbers
-// that sum to f.  Both values have the same sign as f.
+// that sum to f. Both values have the same sign as f.
 //
 // Special cases are:
 //	Modf(±Inf) = ±Inf, NaN
diff --git a/libgo/go/math/rand/race_test.go b/libgo/go/math/rand/race_test.go
new file mode 100644
index 0000000..48f6c29
--- /dev/null
+++ b/libgo/go/math/rand/race_test.go
@@ -0,0 +1,47 @@
+// Copyright 2016 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 rand
+
+import (
+	"sync"
+	"testing"
+)
+
+// TestConcurrent exercises the rand API concurrently, triggering situations
+// where the race detector is likely to detect issues.
+func TestConcurrent(t *testing.T) {
+	const (
+		numRoutines = 10
+		numCycles   = 10
+	)
+	var wg sync.WaitGroup
+	defer wg.Wait()
+	wg.Add(numRoutines)
+	for i := 0; i < numRoutines; i++ {
+		go func(i int) {
+			defer wg.Done()
+			buf := make([]byte, 997)
+			for j := 0; j < numCycles; j++ {
+				var seed int64
+				seed += int64(ExpFloat64())
+				seed += int64(Float32())
+				seed += int64(Float64())
+				seed += int64(Intn(Int()))
+				seed += int64(Int31n(Int31()))
+				seed += int64(Int63n(Int63()))
+				seed += int64(NormFloat64())
+				seed += int64(Uint32())
+				for _, p := range Perm(10) {
+					seed += int64(p)
+				}
+				Read(buf)
+				for _, b := range buf {
+					seed += int64(b)
+				}
+				Seed(int64(i*j) * seed)
+			}
+		}(i)
+	}
+}
diff --git a/libgo/go/math/rand/rand.go b/libgo/go/math/rand/rand.go
index d693bfb..dd8d43c 100644
--- a/libgo/go/math/rand/rand.go
+++ b/libgo/go/math/rand/rand.go
@@ -33,14 +33,32 @@ func NewSource(seed int64) Source {
 // A Rand is a source of random numbers.
 type Rand struct {
 	src Source
+
+	// readVal contains remainder of 63-bit integer used for bytes
+	// generation during most recent Read call.
+	// It is saved so next Read call can start where the previous
+	// one finished.
+	readVal int64
+	// readPos indicates the number of low-order bytes of readVal
+	// that are still valid.
+	readPos int8
 }
 
 // New returns a new Rand that uses random values from src
 // to generate other random values.
-func New(src Source) *Rand { return &Rand{src} }
+func New(src Source) *Rand { return &Rand{src: src} }
 
 // Seed uses the provided seed value to initialize the generator to a deterministic state.
-func (r *Rand) Seed(seed int64) { r.src.Seed(seed) }
+// Seed should not be called concurrently with any other Rand method.
+func (r *Rand) Seed(seed int64) {
+	if lk, ok := r.src.(*lockedSource); ok {
+		lk.seedPos(seed, &r.readPos)
+		return
+	}
+
+	r.src.Seed(seed)
+	r.readPos = 0
+}
 
 // Int63 returns a non-negative pseudo-random 63-bit integer as an int64.
 func (r *Rand) Int63() int64 { return r.src.Int63() }
@@ -160,15 +178,29 @@ func (r *Rand) Perm(n int) []int {
 
 // Read generates len(p) random bytes and writes them into p. It
 // always returns len(p) and a nil error.
+// Read should not be called concurrently with any other Rand method.
 func (r *Rand) Read(p []byte) (n int, err error) {
-	for i := 0; i < len(p); i += 7 {
-		val := r.src.Int63()
-		for j := 0; i+j < len(p) && j < 7; j++ {
-			p[i+j] = byte(val)
-			val >>= 8
+	if lk, ok := r.src.(*lockedSource); ok {
+		return lk.read(p, &r.readVal, &r.readPos)
+	}
+	return read(p, r.Int63, &r.readVal, &r.readPos)
+}
+
+func read(p []byte, int63 func() int64, readVal *int64, readPos *int8) (n int, err error) {
+	pos := *readPos
+	val := *readVal
+	for n = 0; n < len(p); n++ {
+		if pos == 0 {
+			val = int63()
+			pos = 7
 		}
+		p[n] = byte(val)
+		val >>= 8
+		pos--
 	}
-	return len(p), nil
+	*readPos = pos
+	*readVal = val
+	return
 }
 
 /*
@@ -179,7 +211,9 @@ var globalRand = New(&lockedSource{src: NewSource(1)})
 
 // Seed uses the provided seed value to initialize the default Source to a
 // deterministic state. If Seed is not called, the generator behaves as
-// if seeded by Seed(1).
+// if seeded by Seed(1). Seed values that have the same remainder when
+// divided by 2^31-1 generate the same pseudo-random sequence.
+// Seed, unlike the Rand.Seed method, is safe for concurrent use.
 func Seed(seed int64) { globalRand.Seed(seed) }
 
 // Int63 returns a non-negative pseudo-random 63-bit integer as an int64
@@ -226,6 +260,7 @@ func Perm(n int) []int { return globalRand.Perm(n) }
 
 // Read generates len(p) random bytes from the default Source and
 // writes them into p. It always returns len(p) and a nil error.
+// Read, unlike the Rand.Read method, is safe for concurrent use.
 func Read(p []byte) (n int, err error) { return globalRand.Read(p) }
 
 // NormFloat64 returns a normally distributed float64 in the range
@@ -266,3 +301,19 @@ func (r *lockedSource) Seed(seed int64) {
 	r.src.Seed(seed)
 	r.lk.Unlock()
 }
+
+// seedPos implements Seed for a lockedSource without a race condiiton.
+func (r *lockedSource) seedPos(seed int64, readPos *int8) {
+	r.lk.Lock()
+	r.src.Seed(seed)
+	*readPos = 0
+	r.lk.Unlock()
+}
+
+// read implements Read for a lockedSource without a race condition.
+func (r *lockedSource) read(p []byte, readVal *int64, readPos *int8) (n int, err error) {
+	r.lk.Lock()
+	n, err = read(p, r.src.Int63, readVal, readPos)
+	r.lk.Unlock()
+	return
+}
diff --git a/libgo/go/math/rand/rand_test.go b/libgo/go/math/rand/rand_test.go
index 8d68335..6f31279 100644
--- a/libgo/go/math/rand/rand_test.go
+++ b/libgo/go/math/rand/rand_test.go
@@ -5,13 +5,16 @@
 package rand
 
 import (
+	"bytes"
 	"errors"
 	"fmt"
 	"internal/testenv"
+	"io"
 	"math"
 	"os"
 	"runtime"
 	"testing"
+	"testing/iotest"
 )
 
 const (
@@ -373,7 +376,7 @@ func testReadUniformity(t *testing.T, n int, seed int64) {
 	checkSampleDistribution(t, samples, expected)
 }
 
-func TestRead(t *testing.T) {
+func TestReadUniformity(t *testing.T) {
 	testBufferSizes := []int{
 		2, 4, 7, 64, 1024, 1 << 16, 1 << 20,
 	}
@@ -394,7 +397,42 @@ func TestReadEmpty(t *testing.T) {
 	if n != 0 {
 		t.Errorf("Read into empty buffer returned unexpected n of %d", n)
 	}
+}
+
+func TestReadByOneByte(t *testing.T) {
+	r := New(NewSource(1))
+	b1 := make([]byte, 100)
+	_, err := io.ReadFull(iotest.OneByteReader(r), b1)
+	if err != nil {
+		t.Errorf("read by one byte: %v", err)
+	}
+	r = New(NewSource(1))
+	b2 := make([]byte, 100)
+	_, err = r.Read(b2)
+	if err != nil {
+		t.Errorf("read: %v", err)
+	}
+	if !bytes.Equal(b1, b2) {
+		t.Errorf("read by one byte vs single read:\n%x\n%x", b1, b2)
+	}
+}
 
+func TestReadSeedReset(t *testing.T) {
+	r := New(NewSource(42))
+	b1 := make([]byte, 128)
+	_, err := r.Read(b1)
+	if err != nil {
+		t.Errorf("read: %v", err)
+	}
+	r.Seed(42)
+	b2 := make([]byte, 128)
+	_, err = r.Read(b2)
+	if err != nil {
+		t.Errorf("read: %v", err)
+	}
+	if !bytes.Equal(b1, b2) {
+		t.Errorf("mismatch after re-seed:\n%x\n%x", b1, b2)
+	}
 }
 
 // Benchmarks
diff --git a/libgo/go/math/rand/regress_test.go b/libgo/go/math/rand/regress_test.go
index 9ae5357..4dd965c 100644
--- a/libgo/go/math/rand/regress_test.go
+++ b/libgo/go/math/rand/regress_test.go
@@ -1,4 +1,4 @@
-// Copyright 2014 The Go Authors.  All rights reserved.
+// Copyright 2014 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.
 
@@ -342,25 +342,25 @@ var regressGolden = []interface{}{
 	[]int{8, 7, 5, 3, 4, 6, 0, 1, 2},    // Perm(9)
 	[]int{1, 0, 2, 5, 7, 6, 9, 8, 3, 4}, // Perm(10)
 	[]byte{0x1},                         // Read([0])
-	[]byte{0xc0, 0x41, 0xd3, 0xff, 0x12, 0x4, 0x5b},                   // Read([0 0 0 0 0 0 0])
-	[]byte{0x73, 0xc8, 0x6e, 0x4f, 0xf9, 0x5f, 0xf6, 0x62},            // Read([0 0 0 0 0 0 0 0])
-	[]byte{0x4a, 0x2d, 0xb, 0x75, 0xfb, 0x18, 0xd, 0xaf, 0x48},        // Read([0 0 0 0 0 0 0 0 0])
-	[]byte{0x39, 0x46, 0x51, 0x85, 0xf, 0xd4, 0xa1, 0x78, 0x89, 0x2e}, // Read([0 0 0 0 0 0 0 0 0 0])
-	[]byte{0x51}, // Read([0])
-	[]byte{0x4e, 0xe2, 0xd3, 0xd0, 0xd0, 0xde, 0x6b},                   // Read([0 0 0 0 0 0 0])
-	[]byte{0xf8, 0xf9, 0xb4, 0x4c, 0xe8, 0x5f, 0xf0, 0x44},             // Read([0 0 0 0 0 0 0 0])
-	[]byte{0x3b, 0xbf, 0x85, 0x7a, 0xab, 0x99, 0xc5, 0xb2, 0x52},       // Read([0 0 0 0 0 0 0 0 0])
-	[]byte{0xa8, 0xae, 0xb7, 0x9e, 0xf8, 0x56, 0xf6, 0x59, 0xc1, 0x8f}, // Read([0 0 0 0 0 0 0 0 0 0])
-	[]byte{0xc7}, // Read([0])
-	[]byte{0x5f, 0x67, 0xcf, 0xe2, 0x42, 0xcf, 0x3c},                   // Read([0 0 0 0 0 0 0])
-	[]byte{0xc3, 0x54, 0xf3, 0xed, 0xe2, 0xd6, 0xbe, 0xcc},             // Read([0 0 0 0 0 0 0 0])
-	[]byte{0x6a, 0x9f, 0x4a, 0x57, 0x8b, 0xcb, 0x9e, 0xf2, 0xd4},       // Read([0 0 0 0 0 0 0 0 0])
-	[]byte{0x6d, 0x29, 0x97, 0x61, 0xea, 0x9e, 0x4f, 0x5a, 0xa6, 0xae}, // Read([0 0 0 0 0 0 0 0 0 0])
-	[]byte{0xaa}, // Read([0])
-	[]byte{0x20, 0xef, 0xcd, 0x6c, 0xea, 0x84, 0xb6},                   // Read([0 0 0 0 0 0 0])
-	[]byte{0x92, 0x5e, 0x60, 0x7b, 0xe0, 0x63, 0x71, 0x6f},             // Read([0 0 0 0 0 0 0 0])
-	[]byte{0x4, 0x5c, 0x3f, 0x0, 0xf, 0x8a, 0x79, 0x6b, 0xce},          // Read([0 0 0 0 0 0 0 0 0])
-	[]byte{0xaa, 0xca, 0xee, 0xdf, 0xad, 0x5b, 0x50, 0x66, 0x64, 0xe8}, // Read([0 0 0 0 0 0 0 0 0 0])
+	[]byte{0x94, 0xfd, 0xc2, 0xfa, 0x2f, 0xfc, 0xc0},                 // Read([0 0 0 0 0 0 0])
+	[]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{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])
+	[]byte{0x8, 0x75, 0xd6, 0x4e, 0xe2, 0xd3, 0xd0, 0xd0, 0xde, 0x6b}, // Read([0 0 0 0 0 0 0 0 0 0])
+	[]byte{0xf8}, // Read([0])
+	[]byte{0xf9, 0xb4, 0x4c, 0xe8, 0x5f, 0xf0, 0x44},                   // Read([0 0 0 0 0 0 0])
+	[]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{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()
diff --git a/libgo/go/math/remainder.go b/libgo/go/math/remainder.go
index 98bb04d..c6a4c7d 100644
--- a/libgo/go/math/remainder.go
+++ b/libgo/go/math/remainder.go
@@ -6,7 +6,7 @@ package math
 
 // The original C code and the comment below are from
 // FreeBSD's /usr/src/lib/msun/src/e_remainder.c and came
-// with this notice.  The go code is a simplified version of
+// with this notice. The go code is a simplified version of
 // the original C.
 //
 // ====================================================
diff --git a/libgo/go/math/sqrt.go b/libgo/go/math/sqrt.go
index 86c0452..c65c32e 100644
--- a/libgo/go/math/sqrt.go
+++ b/libgo/go/math/sqrt.go
@@ -13,7 +13,7 @@ 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
+// came with this notice. The go code is a simplified
 // version of the original C.
 //
 // ====================================================
@@ -86,7 +86,7 @@ func Sqrt(x float64) float64 {
 //      equal to huge for some floating point number "huge" and "tiny".
 //
 //
-// Notes:  Rounding mode detection omitted.  The constants "mask", "shift",
+// Notes:  Rounding mode detection omitted. The constants "mask", "shift",
 // and "bias" are found in src/math/bits.go
 
 // Sqrt returns the square root of x.
@@ -148,7 +148,3 @@ func sqrt(x float64) float64 {
 	ix = q>>1 + uint64(exp-1+bias)<<shift // significand + biased exponent
 	return Float64frombits(ix)
 }
-
-func sqrtC(f float64, r *float64) {
-	*r = sqrt(f)
-}