integrate Math::BigInt-1.63

Subject: Re: [perl #16997] Math::BigFloat hang on bsqrt [ANNOUNCE v1.63
Message-Id: <200209082022.g88KMGY20194@crypt.org>

p4raw-id: //depot/perl@17880

7 files changed, 67 insertions, 41 deletions

diff --git a/lib/Math/BigFloat.pm b/lib/Math/BigFloat.pm
index f58aaa7fe5..4e93a2ffb7 100644
--- a/lib/Math/BigFloat.pm
+++ b/lib/Math/BigFloat.pm
@@ -12,12 +12,10 @@ package Math::BigFloat;
 #   _p: precision
 #   _f: flags, used to signal MBI not to touch our private parts
 
-$VERSION = '1.37';
+$VERSION = '1.38';
 require 5.005;
 use Exporter;
-use File::Spec;
-# use Math::BigInt;
-@ISA =       qw( Exporter Math::BigInt);
+@ISA =       qw(Exporter Math::BigInt);
 
 use strict;
 use vars qw/$AUTOLOAD $accuracy $precision $div_scale $round_mode $rnd_mode/;
@@ -1074,7 +1072,8 @@ sub bsqrt
   my @params = $x->_find_round_parameters($a,$p,$r);
 
   # no rounding at all, so must use fallback
-  if (scalar @params == 1)
+  if ((scalar @params == 1) ||
+      (!defined($params[1] || $params[2])))
     {
     # simulate old behaviour
     $params[1] = $self->div_scale();	# and round to it as accuracy
@@ -1086,7 +1085,7 @@ sub bsqrt
     {
     # the 4 below is empirical, and there might be cases where it is not
     # enough...
-    $scale = abs($params[1] || $params[2]) + 4;	# take whatever is defined
+    $scale = abs($params[1] || $params[2]) + 4; # take whatever is defined
     }
 
   # when user set globals, they would interfere with our calculation, so
@@ -1103,7 +1102,6 @@ sub bsqrt
   my $xas = $x->as_number();
   my $gs = $xas->copy()->bsqrt();	# some guess
 
-#  print "guess $gs\n";
   if (($x->{_e}->{sign} ne '-')		# guess can't be accurate if there are
 					# digits after the dot
    && ($xas->bacmp($gs * $gs) == 0))	# guess hit the nail on the head?
@@ -1128,29 +1126,37 @@ sub bsqrt
     ${"$self\::accuracy"} = $ab; ${"$self\::precision"} = $pb;
     return $x;
     }
-  $gs = $self->new( $gs );		# BigInt to BigFloat
-
-  my $lx = $x->{_m}->length();
-  $scale = $lx if $scale < $lx;
-  my $e = $self->new("1E-$scale");	# make test variable
-
-  my $y = $x->copy();
-  my $two = $self->new(2);
-  my $diff = $e;
-  # promote BigInts and it's subclasses (except when already a BigFloat)
-  $y = $self->new($y) unless $y->isa('Math::BigFloat'); 
+ 
+  # sqrt(2) = 1.4 because sqrt(2*100) = 1.4*10; so we can increase the accuracy
+  # of the result by multipyling the input by 100 and then divide the integer
+  # result of sqrt(input) by 10. Rounding afterwards returns the real result.
+  # this will transform 123.456 (in $x) into 123456 (in $y1)
+  my $y1 = $x->{_m}->copy();
+  # We now make sure that $y1 has the same odd or even number of digits than
+  # $x had. So when _e of $x is odd, we must shift $y1 by one digit left,
+  # because we always must multiply by steps of 100 (sqrt(100) is 10) and not
+  # steps of 10. The length of $x does not count, since an even or odd number
+  # of digits before the dot is not changed by adding an even number of digits
+  # after the dot (the result is still odd or even digits long).
+  $y1->bmul(10) if $x->{_e}->is_odd();
+  # now calculate how many digits the result of sqrt(y1) would have
+  my $digits = int($y1->length() / 2);
+  # but we need at least $scale digits, so calculate how many are missing
+  my $shift = $scale - $digits;
+  # that should never happen (we take care of integer guesses above)
+  # $shift = 0 if $shift < 0; 
+  # multiply in steps of 100, by shifting left two times the "missing" digits
+  $y1->blsft($shift*2,10);
+  # now take the square root and truncate to integer
+  $y1->bsqrt();
+  # By "shifting" $y1 right (by creating a negative _e) we calculate the final
+  # result, which is than later rounded to the desired scale.
+  $x->{_m} = $y1;
+  # gs->length() is the number of digits before the dot. Since gs is always
+  # truncated (9.99 => 9), it is always right (if gs was rounded, it would be
+  # '10' and thus gs->length() == 2, which would be wrong).
+  $x->{_e} = $MBI->new(- $y1->length() + $gs->length());
 
-  my $rem;
-  while ($diff->bacmp($e) >= 0)
-    {
-    $rem = $y->copy()->bdiv($gs,$scale);
-    $rem = $y->copy()->bdiv($gs,$scale)->badd($gs)->bdiv($two,$scale);
-    $diff = $rem->copy()->bsub($gs);
-    $gs = $rem->copy();
-    }
-  # copy over to modify $x
-  $x->{_m} = $rem->{_m}; $x->{_e} = $rem->{_e};
-  
   # shortcut to not run trough _find_round_parameters again
   if (defined $params[1])
     {
diff --git a/lib/Math/BigInt.pm b/lib/Math/BigInt.pm
index eef5b3cb67..478a1e7045 100644
--- a/lib/Math/BigInt.pm
+++ b/lib/Math/BigInt.pm
@@ -18,7 +18,7 @@ package Math::BigInt;
 my $class = "Math::BigInt";
 require 5.005;
 
-$VERSION = '1.62';
+$VERSION = '1.63';
 use Exporter;
 @ISA =       qw( Exporter );
 @EXPORT_OK = qw( objectify _swap bgcd blcm); 
@@ -2035,11 +2035,11 @@ sub bsqrt
   my $lastlast = $x+$two;
   while ($last != $x && $lastlast != $x)
     {
-    $lastlast = $last; $last = $x; 
-    $x += $y / $x; 
-    $x /= $two;
+    $lastlast = $last; $last = $x->copy(); 
+    $x->badd($y / $x); 
+    $x->bdiv($two);
     }
-  $x-- if $x * $x > $y;				# overshot?
+  $x->bdec() if $x * $x > $y;				# overshot?
   $x->round(@r);
   }
 
diff --git a/lib/Math/BigInt/t/bare_mif.t b/lib/Math/BigInt/t/bare_mif.t
index b38532a797..15a1448e0a 100644
--- a/lib/Math/BigInt/t/bare_mif.t
+++ b/lib/Math/BigInt/t/bare_mif.t
@@ -28,7 +28,7 @@ BEGIN
     }
   print "# INC = @INC\n";
 
-  plan tests => 661
+  plan tests => 669
     + 1;		# our onw tests
   }
 
diff --git a/lib/Math/BigInt/t/constant.t b/lib/Math/BigInt/t/constant.t
index 2e42a6155a..4e5a17eaa7 100644
--- a/lib/Math/BigInt/t/constant.t
+++ b/lib/Math/BigInt/t/constant.t
@@ -16,13 +16,12 @@ BEGIN
     }
   } 
 
-no warnings 'portable';	# this must wait until after the version check
 use Math::BigInt ':constant';
 
 ok (2 ** 255,'57896044618658097711785492504343953926634992332820282019728792003956564819968');
 
 {
-  local $^W = 0;	# protect against "non-portable" warnings
+  no warnings 'portable';	# protect against "non-portable" warnings
 # hexadecimal constants
 ok (0x123456789012345678901234567890,
     Math::BigInt->new('0x123456789012345678901234567890'));
diff --git a/lib/Math/BigInt/t/mbimbf.inc b/lib/Math/BigInt/t/mbimbf.inc
index 77f63803c1..a00183f0ec 100644
--- a/lib/Math/BigInt/t/mbimbf.inc
+++ b/lib/Math/BigInt/t/mbimbf.inc
@@ -91,8 +91,10 @@ foreach my $class ($mbi,$mbf)
   ok (${"$mbf\::round_mode"},'zero');
   ok (${"$mbi\::round_mode"},'-inf');	# from above
 
+  # reset for further tests
   ${"$mbi\::accuracy"} = undef;
   ${"$mbi\::precision"} = undef;
+  ${"$mbf\::div_scale"} = 40;
 }
 
 # local copies
@@ -319,6 +321,25 @@ $mbf->round_mode('even');
 $x = $mbf->new('740.7')->fdiv('6',4,undef,'zero'); ok ($x,'123.4');
 
 ###############################################################################
+# test that bop(0) does the same than bop(undef)
+
+$x = $mbf->new('1234567890');
+ok ($x->copy()->bsqrt(0),$x->copy()->bsqrt(undef));
+ok ($x->copy->bsqrt(0),'35136.41828644462161665823116758077037159');
+
+ok_undef ($x->{_a});
+
+# test that bsqrt() modifies $x and does not just return something else
+# (especially under BareCalc)
+$z = $x->bsqrt();
+ok ($z,$x); ok ($x,'35136.41828644462161665823116758077037159');
+
+$x = $mbf->new('1.234567890123456789');
+ok ($x->copy()->bpow('0.5',0),$x->copy()->bpow('0.5',undef));
+ok ($x->copy()->bpow('0.5',0),$x->copy()->bsqrt(undef));
+ok ($x->copy()->bpow('2',0),'1.524157875323883675019051998750190521');
+
+###############################################################################
 # test (also under Bare) that bfac() rounds at last step
 
 ok ($mbi->new(12)->bfac(),'479001600');
diff --git a/lib/Math/BigInt/t/mbimbf.t b/lib/Math/BigInt/t/mbimbf.t
index fcc9554b79..67645ef244 100644
--- a/lib/Math/BigInt/t/mbimbf.t
+++ b/lib/Math/BigInt/t/mbimbf.t
@@ -31,12 +31,12 @@ BEGIN
     }
   print "# INC = @INC\n";
 
-  plan tests => 661 
+  plan tests => 669 
     + 16;		# own tests
   }
 
-use Math::BigInt 1.62;
-use Math::BigFloat 1.37;
+use Math::BigInt 1.63;
+use Math::BigFloat 1.38;
 
 use vars qw/$mbi $mbf/;
 
diff --git a/lib/Math/BigInt/t/sub_mif.t b/lib/Math/BigInt/t/sub_mif.t
index 8e5656bbf0..365fe617d7 100644
--- a/lib/Math/BigInt/t/sub_mif.t
+++ b/lib/Math/BigInt/t/sub_mif.t
@@ -28,7 +28,7 @@ BEGIN
     }
   print "# INC = @INC\n";
 
-  plan tests => 661;
+  plan tests => 669;
   }
 
 use Math::BigInt::Subclass;