math/big: avoid unnecessary Newton iteration in Float.Sqrt

An initial draft of the Newton code for Float.Sqrt was structured like
this:

  for condition
    // do Newton iteration..
    prec *= 2

since prec, at the end of the loop, was double the precision used in
the last Newton iteration, the termination condition was set to
2*limit. The code was later rewritten in the form

  for condition
    prec *= 2
    // do Newton iteration..

but condition was not updated, and it's still 2*limit, which is about
double what we actually need, and is triggering the execution of an
additional, and unnecessary, Newton iteration.

This change adjusts the Newton termination condition to the (correct)
value of z.prec, plus 32 guard bits as a safety margin.

name                 old time/op    new time/op    delta
FloatSqrt/64-4          798ns ± 3%     802ns ± 3%     ~     (p=0.458 n=8+8)
FloatSqrt/128-4        1.65µs ± 1%    1.65µs ± 1%     ~     (p=0.290 n=8+8)
FloatSqrt/256-4        3.10µs ± 1%    2.10µs ± 0%  -32.32%  (p=0.000 n=8+7)
FloatSqrt/1000-4       8.83µs ± 1%    4.91µs ± 2%  -44.39%  (p=0.000 n=8+8)
FloatSqrt/10000-4       107µs ± 1%      40µs ± 1%  -62.68%  (p=0.000 n=8+8)
FloatSqrt/100000-4     2.91ms ± 1%    0.96ms ± 1%  -67.13%  (p=0.000 n=8+8)
FloatSqrt/1000000-4     240ms ± 1%      80ms ± 1%  -66.66%  (p=0.000 n=8+8)

name                 old alloc/op   new alloc/op   delta
FloatSqrt/64-4           416B ± 0%      416B ± 0%     ~     (all equal)
FloatSqrt/128-4          720B ± 0%      720B ± 0%     ~     (all equal)
FloatSqrt/256-4        1.34kB ± 0%    0.82kB ± 0%  -39.29%  (p=0.000 n=8+8)
FloatSqrt/1000-4       5.09kB ± 0%    2.50kB ± 0%  -50.94%  (p=0.000 n=8+8)
FloatSqrt/10000-4      45.9kB ± 0%    23.5kB ± 0%  -48.81%  (p=0.000 n=8+8)
FloatSqrt/100000-4      533kB ± 0%     251kB ± 0%  -52.90%  (p=0.000 n=8+8)
FloatSqrt/1000000-4    9.21MB ± 0%    4.61MB ± 0%  -49.98%  (p=0.000 n=8+8)

name                 old allocs/op  new allocs/op  delta
FloatSqrt/64-4           9.00 ± 0%      9.00 ± 0%     ~     (all equal)
FloatSqrt/128-4          13.0 ± 0%      13.0 ± 0%     ~     (all equal)
FloatSqrt/256-4          15.0 ± 0%      12.0 ± 0%  -20.00%  (p=0.000 n=8+8)
FloatSqrt/1000-4         24.0 ± 0%      19.0 ± 0%  -20.83%  (p=0.000 n=8+8)
FloatSqrt/10000-4        40.0 ± 0%      35.0 ± 0%  -12.50%  (p=0.000 n=8+8)
FloatSqrt/100000-4       66.0 ± 0%      55.0 ± 0%  -16.67%  (p=0.000 n=8+8)
FloatSqrt/1000000-4       143 ± 0%       122 ± 0%  -14.69%  (p=0.000 n=8+8)

Change-Id: I4868adb7f8960f2ca20e7792734c2e6211669fc0
Reviewed-on: https://go-review.googlesource.com/75010
Reviewed-by: Robert Griesemer <gri@golang.org>

diff --git a/src/math/big/sqrt.go b/src/math/big/sqrt.go
index d1deb776524930c708bf0be957c2caec02a06a4c..00433cfe7a704e2c022d7459e80e80cbceb6ed3c 100644 (file)
--- a/src/math/big/sqrt.go
+++ b/src/math/big/sqrt.go
@@ -140,7 +140,7 @@ func (z *Float) sqrtInverse(x *Float) {
 
        xf, _ := x.Float64()
        sqi := NewFloat(1 / math.Sqrt(xf))
-       for prec := 2 * z.prec; sqi.prec < prec; {
+       for prec := z.prec + 32; sqi.prec < prec; {
                sqi.prec *= 2
                sqi = ng(sqi)
        }