# _p: precision
# _f: flags, used to signal MBI not to touch our private parts
-$VERSION = '1.39';
+$VERSION = '1.40';
require 5.005;
use Exporter;
@ISA = qw(Exporter Math::BigInt);
# Convert number from internal format to (non-scientific) string format.
# internal format is always normalized (no leading zeros, "-0" => "+0")
my ($self,$x) = ref($_[0]) ? (ref($_[0]),$_[0]) : objectify(1,@_);
- #my $x = shift; my $class = ref($x) || $x;
- #$x = $class->new(shift) unless ref($x);
if ($x->{sign} !~ /^[+-]$/)
{
return $x->{sign} unless $x->{sign} eq '+inf'; # -inf, NaN
return 'inf'; # +inf
}
-
+
my $es = '0'; my $len = 1; my $cad = 0; my $dot = '.';
- my $not_zero = ! $x->is_zero();
+ # $x is zero?
+ my $not_zero = !($x->{sign} eq '+' && $x->{_m}->is_zero());
if ($not_zero)
{
$es = $x->{_m}->bstr();
$len = CORE::length($es);
- if (!$x->{_e}->is_zero())
+ my $e = $x->{_e}->numify();
+ if ($e < 0)
{
- if ($x->{_e}->sign() eq '-')
+ $dot = '';
+ # if _e is bigger than a scalar, the following will blow your memory
+ if ($e <= -$len)
{
- $dot = '';
- if ($x->{_e} <= -$len)
- {
- #print "style: 0.xxxx\n";
- my $r = $x->{_e}->copy(); $r->babs()->bsub( CORE::length($es) );
- $es = '0.'. ('0' x $r) . $es; $cad = -($len+$r);
- }
- else
- {
- #print "insert '.' at $x->{_e} in '$es'\n";
- substr($es,$x->{_e},0) = '.'; $cad = $x->{_e};
- }
+ #print "style: 0.xxxx\n";
+ my $r = abs($e) - $len;
+ $es = '0.'. ('0' x $r) . $es; $cad = -($len+$r);
}
else
{
- # expand with zeros
- $es .= '0' x $x->{_e}; $len += $x->{_e}; $cad = 0;
+ #print "insert '.' at $e in '$es'\n";
+ substr($es,$e,0) = '.'; $cad = $x->{_e};
}
}
+ elsif ($e > 0)
+ {
+ # expand with zeros
+ $es .= '0' x $e; $len += $e; $cad = 0;
+ }
} # if not zero
- $es = $x->{sign}.$es if $x->{sign} eq '-';
- # if set accuracy or precision, pad with zeros
+ $es = '-'.$es if $x->{sign} eq '-';
+ # if set accuracy or precision, pad with zeros on the right side
if ((defined $x->{_a}) && ($not_zero))
{
# 123400 => 6, 0.1234 => 4, 0.001234 => 4
$zeros = $x->{_a} - $len if $cad != $len;
$es .= $dot.'0' x $zeros if $zeros > 0;
}
- elsif ($x->{_p} || 0 < 0)
+ elsif ((($x->{_p} || 0) < 0))
{
# 123400 => 6, 0.1234 => 4, 0.001234 => 6
my $zeros = -$x->{_p} + $cad;
sub broot
{
# calculate $y'th root of $x
-
+
# set up parameters
my ($self,$x,$y,@r) = (ref($_[0]),@_);
+
+ $y = $self->new(2) unless defined $y;
+
# objectify is costly, so avoid it
- if ((!ref($_[0])) || (ref($_[0]) ne ref($_[1])))
+ if ((!ref($x)) || (ref($x) ne ref($y)))
{
- ($self,$x,$y,@r) = objectify(2,@_);
+ ($self,$x,$y,@r) = $self->objectify(2,@_);
}
return $x if $x->modify('broot');
return $x->round(@r)
if $x->is_zero() || $x->is_one() || $x->is_inf() || $y->is_one();
- return $upgrade->broot($x,@r) if defined $upgrade;
+ return $upgrade->new($x)->broot($upgrade->new($y),@r) if defined $upgrade;
if ($CALC->can('_root'))
{
# since we take at least a cubic root, and only 8 ** 1/3 >= 2 (==2):
return $x->bone('+',@r) if $x < 8; # $x=2..7 => 1
- my $org = $x->copy();
- my $l = int($x->length()/$y->numify());
-
- $x->bone(); # keep ref($x), but modify it
- $x->blsft($l,10) if $l != 0; # first guess: 1.('0' x (l/$y))
+ my $num = $x->numify();
- my $last = $self->bzero();
- my $lastlast = $self->bzero();
- #my $lastlast = $x+$y;
- my $divider = $self->new(2);
- my $up = $y-1;
- #print "start $org divider $divider up $up\n";
- while ($last != $x && $lastlast != $x)
+ if ($num <= 1000000)
{
- #print "at $x ($last $lastlast)\n";
- $lastlast = $last; $last = $x->copy();
- #print "at $x ($last ",($org / ($x ** $up)),"\n";
- $x->badd($org / ($x ** 2));
- $x->bdiv($divider);
+ $x = $self->new( int($num ** (1 / $y->numify()) ));
+ return $x->round(@r);
+ }
+
+ # if $n is a power of two, we can repeatedly take sqrt($X) and find the
+ # proper result, because sqrt(sqrt($x)) == root($x,4)
+ # See Calc.pm for more details
+ my $b = $y->as_bin();
+ if ($b =~ /0b1(0+)/)
+ {
+ my $count = CORE::length($1); # 0b100 => len('00') => 2
+ my $cnt = $count; # counter for loop
+ my $shift = $self->new(6);
+ $x->blsft($shift); # add some zeros (even amount)
+ while ($cnt-- > 0)
+ {
+ # 'inflate' $X by adding more zeros
+ $x->blsft($shift);
+ # calculate sqrt($x), $x is now a bit too big, again. In the next
+ # round we make even bigger, again.
+ $x->bsqrt($x);
+ }
+ # $x is still to big, so truncate result
+ $x->brsft($shift);
}
- #print $x ** $y," org ",$org,"\n";
- # correct overshot
- while ($x ** $y < $org)
+ else
{
- #print "correcting $x to ";
- $x->binc();
- #print "$x ( $x ** $y == ",$x ** $y,")\n";
+ # Should compute a guess of the result (by rule of thumb), then improve it
+ # via Newton's method or something similiar.
+ # XXX TODO
+ warn ('broot() not fully implemented in BigInt.');
}
- $x->round(@r);
+ return $x->round(@r);
}
sub exponent
}
my $up = ${"$a[0]::upgrade"};
- #print "Now in objectify, my class is today $a[0]\n";
+ #print "Now in objectify, my class is today $a[0], count = $count\n";
if ($count == 0)
{
while (@_)
div_scale Fallback acccuracy for div
40
-The following values can be set by passing config a reference to a hash:
+The following values can be set by passing C<config()> a reference to a hash:
trap_inf trap_nan
upgrade downgrade precision accuracy round_mode div_scale
The return true when the argument satisfies the condition. C<NaN>, C<+inf>,
C<-inf> are not integers and are neither odd nor even.
+In BigInt, all numbers except C<NaN>, C<+inf> and C<-inf> are integers.
+
=head2 bcmp
$x->bcmp($y);
# determine how many digits fit into an integer and can be safely added
# together plus carry w/o causing an overflow
- # this below detects 15 on a 64 bit system, because after that it becomes
- # 1e16 and not 1000000 :/ I can make it detect 18, but then I get a lot of
- # test failures. Ugh! (Tomake detect 18: uncomment lines marked with *)
use integer;
- my $bi = 5; # approx. 16 bit
- $num = int('9' x $bi);
- # $num = 99999; # *
- # while ( ($num+$num+1) eq '1' . '9' x $bi) # *
- while ( int($num+$num+1) eq '1' . '9' x $bi)
- {
- $bi++; $num = int('9' x $bi);
- # $bi++; $num *= 10; $num += 9; # *
- }
- $bi--; # back off one step
+
+ ############################################################################
+ # the next block is no longer important
+
+ ## this below detects 15 on a 64 bit system, because after that it becomes
+ ## 1e16 and not 1000000 :/ I can make it detect 18, but then I get a lot of
+ ## test failures. Ugh! (Tomake detect 18: uncomment lines marked with *)
+
+ #my $bi = 5; # approx. 16 bit
+ #$num = int('9' x $bi);
+ ## $num = 99999; # *
+ ## while ( ($num+$num+1) eq '1' . '9' x $bi) # *
+ #while ( int($num+$num+1) eq '1' . '9' x $bi)
+ # {
+ # $bi++; $num = int('9' x $bi);
+ # # $bi++; $num *= 10; $num += 9; # *
+ # }
+ #$bi--; # back off one step
# by setting them equal, we ignore the findings and use the default
# one-size-fits-all approach from former versions
- $bi = $e; # XXX, this should work always
+ my $bi = $e; # XXX, this should work always
__PACKAGE__->_base_len($e,$bi); # set and store
# find out how many bits _and, _or and _xor can take (old default = 16)
# I don't think anybody has yet 128 bit scalars, so let's play safe.
local $^W = 0; # don't warn about 'nonportable number'
- $AND_BITS = 15; $XOR_BITS = 15; $OR_BITS = 15;
+ $AND_BITS = 15; $XOR_BITS = 15; $OR_BITS = 15;
# find max bits, we will not go higher than numberofbits that fit into $BASE
# to make _and etc simpler (and faster for smaller, slower for large numbers)
sub _root
{
# take n'th root of $x in place (n >= 3)
- # Compute a guess of the result (by rule of thumb), then improve it via
- # Newton's method.
my ($c,$x,$n) = @_;
if (scalar @$x == 1)
return $x;
}
- # XXX TODO
-
+ # X is more than one element
+ # if $n is a power of two, we can repeatedly take sqrt($X) and find the
+ # proper result, because sqrt(sqrt($x)) == root($x,4)
+ my $b = _as_bin($c,$n);
+ if ($$b =~ /0b1(0+)/)
+ {
+ my $count = CORE::length($1); # 0b100 => len('00') => 2
+ my $cnt = $count; # counter for loop
+ unshift (@$x, 0); # add one element, together with one
+ # more below in the loop this makes 2
+ while ($cnt-- > 0)
+ {
+ # 'inflate' $X by adding one element, basically computing
+ # $x * $BASE * $BASE. This gives us more $BASE_LEN digits for result
+ # since len(sqrt($X)) approx == len($x) / 2.
+ unshift (@$x, 0);
+ # calculate sqrt($x), $x is now one element to big, again. In the next
+ # round we make that two, again.
+ _sqrt($c,$x);
+ }
+ # $x is now one element to big, so truncate result by removing it
+ splice (@$x,0,1);
+ }
+ else
+ {
+ # Should compute a guess of the result (by rule of thumb), then improve it
+ # via Newton's method or something similiar.
+ # XXX TODO
+ warn ('_root() not fully implemented in Calc.');
+ }
$x;
}
}
print "# INC = @INC\n";
- plan tests => 2668;
+ plan tests => 2684;
}
use Math::BigInt lib => 'BareCalc';
8:3:2
-8:3:NaN
# fourths root
-#16:4:2
-#81:4:3
+16:4:2
+81:4:3
+# 2 ** 32
+18446744073709551616:4:65536
+18446744073709551616:8:256
+18446744073709551616:16:16
+18446744073709551616:32:4
+18446744073709551616:64:2
+18446744073709551616:128:1
&bsqrt
145:12
144:12
my $location = $0; $location =~ s/bigintpm.t//;
unshift @INC, $location; # to locate the testing files
chdir 't' if -d 't';
- plan tests => 2668;
+ plan tests => 2684;
}
use Math::BigInt;
# we create the numbers from "patterns", e.g. get a random number and a
# random count and string them together. This means things like
# "100000999999999999911122222222" are much more likely. If we just strung
- # together digits, we would end up with "1272398823211223" etc.
+ # together digits, we would end up with "1272398823211223" etc. It also means
+ # that we get more frequently equal numbers or other special cases.
while (length($As) < $la) { $As .= int(rand(100)) x int(rand(16)); }
while (length($Bs) < $lb) { $Bs .= int(rand(100)) x int(rand(16)); }
print "# seed $seed, ". join(' ',Math::BigInt::Calc->_base_len()),"\n".
"# tried $ADB * $B + $two*$AMB - $AMB\n"
unless ok ($ADB*$B+$two*$AMB-$AMB,$As);
- print "\$ADB * \$B / \$B = ", $ADB * $B / $B, " != $ADB (\$B=$B)\n"
+ print "# seed: $seed, \$ADB * \$B / \$B = ", $ADB * $B / $B, " != $ADB (\$B=$B)\n"
unless ok ($ADB*$B/$B,$ADB);
# swap 'em and try this, too
# $X = ($B/$A)*$A + $B % $A;
unless ok ($ADB*$A+$two*$AMB-$AMB,$Bs);
# print " +$two * $AMB = ",$ADB * $A + $two * $AMB,"\n";
# print " -$AMB = ",$ADB * $A + $two * $AMB - $AMB,"\n";
- print "\$ADB * \$A / \$A = ", $ADB * $A / $A, " != $ADB (\$A=$A)\n"
+ print "# seed $seed, \$ADB * \$A / \$A = ", $ADB * $A / $A, " != $ADB (\$A=$A)\n"
unless ok ($ADB*$A/$A,$ADB);
}
}
print "# INC = @INC\n";
- plan tests => 2668
+ plan tests => 2684
+ 5; # +5 own tests
}
-#include this file into another for subclass testing
+# include this file into another for subclass testing
# This file is nearly identical to bigintpm.t, except that certain results
# are _requird_ to be different due to "upgrading" or "promoting" to BigFloat.
# The reverse is not true, any unmarked results can be either BigInt or
# BigFloat, depending on how good the internal optimization is (e.g. it
-# is usually desirable to have 2 ** 2 return an BigInt, not an BigFloat).
+# is usually desirable to have 2 ** 2 return a BigInt, not a BigFloat).
-# Results that are required to be BigFloat are marked with an "^" at the end.
+# Results that are required to be BigFloat are marked with C<^> at the end.
# Please note that the testcount goes up by two for each extra result marked
# with ^, since then we test whether it has the proper class and that it left
$try .= '$x <=> $y;';
} elsif ($f eq "bround") {
$try .= "$round_mode; \$x->bround(\$y);";
+ } elsif ($f eq "broot") {
+ $try .= "\$x->broot(\$y);";
} elsif ($f eq "bacmp"){
$try .= '$x->bacmp($y);';
} elsif ($f eq "badd"){
10000000000000000:17
-123:3
215960156869840440586892398248:30
+# broot always upgrades
+&broot
+144:2:12^
+123:2:11.09053650640941716205160010260993291846^
+# bsqrt always upgrades
&bsqrt
145:12.04159457879229548012824103037860805243^
144:12^
}
print "# INC = @INC\n";
- plan tests => 2074
+ plan tests => 2082
+ 2; # our own tests
}
projects / p5sagit/p5-mst-13.2.git / commitdiff