From: Jarkko Hietaniemi Date: Sat, 16 Jun 2001 14:59:38 +0000 (+0000) Subject: Math::BigInt 1.35 from Tels. X-Git-Url: http://git.shadowcat.co.uk/gitweb/gitweb.cgi?a=commitdiff_plain;h=58cde26e9566cfa700a6d18945c4428f22c6be85;p=p5sagit%2Fp5-mst-13.2.git Math::BigInt 1.35 from Tels. p4raw-id: //depot/perl@10628 --- diff --git a/MANIFEST b/MANIFEST index e5f39df..4cfeff8 100644 --- a/MANIFEST +++ b/MANIFEST @@ -1597,6 +1597,7 @@ t/lib/lc-currency.t See if Locale::Codes work t/lib/lc-language.t See if Locale::Codes work t/lib/lc-maketext.t See if Locale::Maketext works t/lib/lc-uk.t See if Locale::Codes work +t/lib/mbimbf.t BigInt/BigFloat accuracy, precicion and fallback, round_mode t/lib/md5-aaa.t See if Digest::MD5 extension works t/lib/md5-align.t See if Digest::MD5 extension works t/lib/md5-badf.t See if Digest::MD5 extension works diff --git a/lib/Math/BigFloat.pm b/lib/Math/BigFloat.pm index edde97a..04622ee 100644 --- a/lib/Math/BigFloat.pm +++ b/lib/Math/BigFloat.pm @@ -1,434 +1,1387 @@ +#!/usr/bin/perl -w + +# mark.biggar@TrustedSysLabs.com + +# The following hash values are internally used: +# _e: exponent (BigInt) +# _m: mantissa (absolute BigInt) +# sign: +,-,"NaN" if not a number +# _a: accuracy +# _p: precision +# _cow: Copy-On-Write (NRY) + package Math::BigFloat; -use Math::BigInt; +$VERSION = 1.15; +require 5.005; +use Exporter; +use Math::BigInt qw/trace objectify/; +@ISA = qw( Exporter Math::BigInt); +# can not export bneg/babs since the are only in MBI +@EXPORT_OK = qw( + bcmp + badd bmul bdiv bmod bnorm bsub + bgcd blcm bround bfround + bpow bnan bzero bfloor bceil + bacmp bstr binc bdec bint binf + is_odd is_even is_nan is_inf + is_zero is_one sign + ); -use Exporter; # just for use to be happy -@ISA = (Exporter); -$VERSION = '0.03'; +#@EXPORT = qw( ); +use strict; +use vars qw/$AUTOLOAD $accuracy $precision $div_scale $rnd_mode/; +my $class = "Math::BigFloat"; use overload -'+' => sub {new Math::BigFloat &fadd}, -'-' => sub {new Math::BigFloat - $_[2]? fsub($_[1],${$_[0]}) : fsub(${$_[0]},$_[1])}, -'<=>' => sub {$_[2]? fcmp($_[1],${$_[0]}) : fcmp(${$_[0]},$_[1])}, -'cmp' => sub {$_[2]? ($_[1] cmp ${$_[0]}) : (${$_[0]} cmp $_[1])}, -'*' => sub {new Math::BigFloat &fmul}, -'/' => sub {new Math::BigFloat - $_[2]? scalar fdiv($_[1],${$_[0]}) : - scalar fdiv(${$_[0]},$_[1])}, -'%' => sub {new Math::BigFloat - $_[2]? scalar fmod($_[1],${$_[0]}) : - scalar fmod(${$_[0]},$_[1])}, -'neg' => sub {new Math::BigFloat &fneg}, -'abs' => sub {new Math::BigFloat &fabs}, -'int' => sub {new Math::BigInt &f2int}, - -qw( -"" stringify -0+ numify) # Order of arguments unsignificant +'<=>' => sub { + $_[2] ? + $class->bcmp($_[1],$_[0]) : + $class->bcmp($_[0],$_[1])}, +'int' => sub { $_[0]->copy()->bround(0,'trunc'); }, ; -sub new { - my ($class) = shift; - my ($foo) = fnorm(shift); - bless \$foo, $class; +# are NaNs ok? +my $NaNOK=1; +# set to 1 for tracing +my $trace = 0; +# constant for easier life +my $nan = 'NaN'; +my $ten = Math::BigInt->new(10); # shortcut for speed + +# Rounding modes one of 'even', 'odd', '+inf', '-inf', 'zero' or 'trunc' +$rnd_mode = 'even'; +$accuracy = undef; +$precision = undef; +$div_scale = 40; + +{ + # checks for AUTOLOAD + my %methods = map { $_ => 1 } + qw / fadd fsub fmul fdiv fround ffround fsqrt fmod fstr fsstr fpow fnorm + fabs fneg fint fcmp fzero fnan finc fdec + /; + + sub method_valid { return exists $methods{$_[0]||''}; } } -sub numify { 0 + "${$_[0]}" } # Not needed, additional overhead - # comparing to direct compilation based on - # stringify -sub stringify { - my $n = ${$_[0]}; +############################################################################## +# constructors - my $minus = ($n =~ s/^([+-])// && $1 eq '-'); - $n =~ s/E//; +sub new + { + # create a new BigFloat object from a string or another bigfloat object. + # _e: exponent + # _m: mantissa + # sign => sign (+/-), or "NaN" - $n =~ s/([-+]\d+)$//; + trace (@_); + my $class = shift; + + my $wanted = shift; # avoid numify call by not using || here + return $class->bzero() if !defined $wanted; # default to 0 + return $wanted->copy() if ref($wanted) eq $class; - my $e = $1; - my $ln = length($n); + my $round = shift; $round = 0 if !defined $round; # no rounding as default + my $self = {}; bless $self, $class; + #shortcut for bigints + if (ref($wanted) eq 'Math::BigInt') + { + $self->{_m} = $wanted; + $self->{_e} = Math::BigInt->new(0); + $self->{_m}->babs(); + $self->{sign} = $wanted->sign(); + return $self; + } + # got string + # handle '+inf', '-inf' first + if ($wanted =~ /^[+-]inf$/) + { + $self->{_e} = Math::BigInt->new(0); + $self->{_m} = Math::BigInt->new(0); + $self->{sign} = $wanted; + return $self; + } + #print "new string '$wanted'\n"; + my ($mis,$miv,$mfv,$es,$ev) = Math::BigInt::_split(\$wanted); + if (!ref $mis) + { + die "$wanted is not a number initialized to $class" if !$NaNOK; + $self->{_e} = Math::BigInt->new(0); + $self->{_m} = Math::BigInt->new(0); + $self->{sign} = $nan; + } + else + { + # make integer from mantissa by adjusting exp, then convert to bigint + $self->{_e} = Math::BigInt->new("$$es$$ev"); # exponent + $self->{_m} = Math::BigInt->new("$$mis$$miv$$mfv"); # create mantissa + # 3.123E0 = 3123E-3, and 3.123E-2 => 3123E-5 + $self->{_e} -= CORE::length($$mfv); + $self->{sign} = $self->{_m}->sign(); $self->{_m}->babs(); + } + #print "$wanted => $self->{sign} $self->{value}->[0]\n"; + $self->bnorm(); # first normalize + # if any of the globals is set, round to them and thus store them insid $self + $self->round($accuracy,$precision,$rnd_mode) + if defined $accuracy || defined $precision; + return $self; + } - if ( defined $e ) +# some shortcuts for easier life +sub bfloat + { + # exportable version of new + trace(@_); + return $class->new(@_); + } + +sub bint + { + # exportable version of new + trace(@_); + return $class->new(@_,0)->bround(0,'trunc'); + } + +sub bnan + { + # create a bigfloat 'NaN', if given a BigFloat, set it to 'NaN' + my $self = shift; + $self = $class if !defined $self; + if (!ref($self)) { - if ($e > 0) { - $n .= "0" x $e . '.'; - } elsif (abs($e) < $ln) { - substr($n, $ln + $e, 0) = '.'; - } else { - $n = '.' . ("0" x (abs($e) - $ln)) . $n; - } + my $c = $self; $self = {}; bless $self, $c; } - $n = "-$n" if $minus; + $self->{_e} = new Math::BigInt 0; + $self->{_m} = new Math::BigInt 0; + $self->{sign} = $nan; + trace('NaN'); + return $self; + } - # 1 while $n =~ s/(.*\d)(\d\d\d)/$1,$2/; +sub binf + { + # create a bigfloat '+-inf', if given a BigFloat, set it to '+-inf' + my $self = shift; + my $sign = shift; $sign = '+' if !defined $sign || $sign ne '-'; - return $n; -} + $self = $class if !defined $self; + if (!ref($self)) + { + my $c = $self; $self = {}; bless $self, $c; + } + $self->{_e} = new Math::BigInt 0; + $self->{_m} = new Math::BigInt 0; + $self->{sign} = $sign.'inf'; + trace('inf'); + return $self; + } -sub import { - shift; - return unless @_; - if (@_ == 1 && $_[0] ne ':constant') { - if ($_[0] > 0) { - if ($VERSION < $_[0]) { - die __PACKAGE__.": $_[0] required--this is only version $VERSION"; +sub bzero + { + # create a bigfloat '+0', if given a BigFloat, set it to 0 + my $self = shift; + $self = $class if !defined $self; + if (!ref($self)) + { + my $c = $self; $self = {}; bless $self, $c; + } + $self->{_m} = new Math::BigInt 0; + $self->{_e} = new Math::BigInt 1; + $self->{sign} = '+'; + trace('0'); + return $self; + } + +############################################################################## +# string conversation + +sub bstr + { + # (ref to BFLOAT or num_str ) return num_str + # Convert number from internal format to (non-scientific) string format. + # internal format is always normalized (no leading zeros, "-0" => "+0") + trace(@_); + my ($self,$x) = objectify(1,@_); + + #return "Oups! e was $nan" if $x->{_e}->{sign} eq $nan; + #return "Oups! m was $nan" if $x->{_m}->{sign} eq $nan; + return $x->{sign} if $x->{sign} !~ /^[+-]$/; + return '0' if $x->is_zero(); + + my $es = $x->{_m}->bstr(); + if ($x->{_e}->is_zero()) + { + $es = $x->{sign}.$es if $x->{sign} eq '-'; + return $es; + } + + if ($x->{_e}->sign() eq '-') + { + if ($x->{_e} <= -CORE::length($es)) + { + # print "style: 0.xxxx\n"; + my $r = $x->{_e}->copy(); $r->babs()->bsub( CORE::length($es) ); + $es = '0.'. ('0' x $r) . $es; + } + else + { + # print "insert '.' at $x->{_e} in '$es'\n"; + substr($es,$x->{_e},0) = '.'; } - } else { - die __PACKAGE__.": unknown import: $_[0]"; } + else + { + # expand with zeros + $es .= '0' x $x->{_e}; + } + $es = $x->{sign}.$es if $x->{sign} eq '-'; + return $es; } - overload::constant float => sub {Math::BigFloat->new(shift)}; -} -$div_scale = 40; +sub bsstr + { + # (ref to BFLOAT or num_str ) return num_str + # Convert number from internal format to scientific string format. + # internal format is always normalized (no leading zeros, "-0E0" => "+0E0") + trace(@_); + my ($self,$x) = objectify(1,@_); -# Rounding modes one of 'even', 'odd', '+inf', '-inf', 'zero' or 'trunc'. + return "Oups! e was $nan" if $x->{_e}->{sign} eq $nan; + return "Oups! m was $nan" if $x->{_m}->{sign} eq $nan; + return $x->{sign} if $x->{sign} !~ /^[+-]$/; + my $sign = $x->{_e}->{sign}; $sign = '' if $sign eq '-'; + my $sep = 'e'.$sign; + return $x->{_m}->bstr().$sep.$x->{_e}->bstr(); + } + +sub numify + { + # Make a number from a BigFloat object + # simple return string and let Perl's atoi() handle the rest + trace (@_); + my ($self,$x) = objectify(1,@_); + return $x->bsstr(); + } -$rnd_mode = 'even'; +############################################################################## +# public stuff (usually prefixed with "b") + +# really? Just for exporting them is not what I had in mind +#sub babs +# { +# $class->SUPER::babs($class,@_); +# } +#sub bneg +# { +# $class->SUPER::bneg($class,@_); +# } +#sub bnot +# { +# $class->SUPER::bnot($class,@_); +# } + +sub bcmp + { + # Compares 2 values. Returns one of undef, <0, =0, >0. (suitable for sort) + # (BFLOAT or num_str, BFLOAT or num_str) return cond_code + my ($self,$x,$y) = objectify(2,@_); + return undef if (($x->{sign} eq $nan) || ($y->{sign} eq $nan)); + + # check sign + return 1 if $x->{sign} eq '+' && $y->{sign} eq '-'; + return -1 if $x->{sign} eq '-' && $y->{sign} eq '+'; # does also -x <=> 0 + + return 0 if $x->is_zero() && $y->is_zero(); # 0 <=> 0 + return -1 if $x->is_zero() && $y->{sign} eq '+'; # 0 <=> +y + return 1 if $y->is_zero() && $x->{sign} eq '+'; # +x <=> 0 + + # adjust so that exponents are equal + my $lx = $x->{_m}->length() + $x->{_e}; + my $ly = $y->{_m}->length() + $y->{_e}; + # print "x $x y $y lx $lx ly $ly\n"; + my $l = $lx - $ly; $l = -$l if $x->{sign} eq '-'; + # print "$l $x->{sign}\n"; + return $l if $l != 0; + + # lens are equal, so compare mantissa, if equal, compare exponents + # this assumes normaized numbers (no trailing zeros etc) + my $rc = $x->{_m} <=> $y->{_m} || $x->{_e} <=> $y->{_e}; + $rc = -$rc if $x->{sign} eq '-'; # -124 < -123 + return $rc; + } + +sub bacmp + { + # Compares 2 values, ignoring their signs. + # Returns one of undef, <0, =0, >0. (suitable for sort) + # (BFLOAT or num_str, BFLOAT or num_str) return cond_code + my ($self,$x,$y) = objectify(2,@_); + return undef if (($x->{sign} eq $nan) || ($y->{sign} eq $nan)); + + # signs are ignored, so check length + # length(x) is length(m)+e aka length of non-fraction part + # the longer one is bigger + my $l = $x->length() - $y->length(); + #print "$l\n"; + return $l if $l != 0; + #print "equal lengths\n"; + + # if both are equal long, make full compare + # first compare only the mantissa + # if mantissa are equal, compare fractions + + return $x->{_m} <=> $y->{_m} || $x->{_e} <=> $y->{_e}; + } -sub fadd; sub fsub; sub fmul; sub fdiv; -sub fneg; sub fabs; sub fcmp; -sub fround; sub ffround; -sub fnorm; sub fsqrt; - -# Convert a number to canonical string form. -# Takes something that looks like a number and converts it to -# the form /^[+-]\d+E[+-]\d+$/. -sub fnorm { #(string) return fnum_str - local($_) = @_; - s/\s+//g; # strip white space - no warnings; # $4 and $5 below might legitimately be undefined - if (/^([+-]?)(\d*)(\.(\d*))?([Ee]([+-]?\d+))?$/ && "$2$4" ne '') { - &norm(($1 ? "$1$2$4" : "+$2$4"),(($4 ne '') ? $6-length($4) : $6)); - } else { - 'NaN'; +sub badd + { + # add second arg (BFLOAT or string) to first (BFLOAT) (modifies first) + # return result as BFLOAT + trace(@_); + my ($self,$x,$y,$a,$p,$r) = objectify(2,@_); + + #print "add $x ",ref($x)," $y ",ref($y),"\n"; + return $x->bnan() if (($x->{sign} eq $nan) || ($y->{sign} eq $nan)); + + # speed: no add for 0+y or x+0 + return $x if $y->is_zero(); # x+0 + if ($x->is_zero()) # 0+y + { + # make copy, clobbering up x (modify in place!) + $x->{_e} = $y->{_e}->copy(); + $x->{_m} = $y->{_m}->copy(); + $x->{sign} = $y->{sign} || $nan; + return $x->round($a,$p,$r,$y); } -} + + # take lower of the two e's and adapt m1 to it to match m2 + my $e = $y->{_e}; $e = Math::BigInt::bzero() if !defined $e; # if no BFLOAT + $e = $e - $x->{_e}; + my $add = $y->{_m}->copy(); + if ($e < 0) + { + #print "e < 0\n"; + #print "\$x->{_m}: $x->{_m} "; + #print "\$x->{_e}: $x->{_e}\n"; + my $e1 = $e->copy()->babs(); + $x->{_m} *= (10 ** $e1); + $x->{_e} += $e; # need the sign of e + #$x->{_m} += $y->{_m}; + #print "\$x->{_m}: $x->{_m} "; + #print "\$x->{_e}: $x->{_e}\n"; + } + elsif ($e > 0) + { + #print "e > 0\n"; + #print "\$x->{_m}: $x->{_m} \$y->{_m}: $y->{_m} \$e: $e ",ref($e),"\n"; + $add *= (10 ** $e); + #$x->{_m} += $y->{_m} * (10 ** $e); + #print "\$x->{_m}: $x->{_m}\n"; + } + # else: both e are same, so leave them + #print "badd $x->{sign}$x->{_m} + $y->{sign}$add\n"; + # fiddle with signs + $x->{_m}->{sign} = $x->{sign}; + $add->{sign} = $y->{sign}; + # finally do add/sub + $x->{_m} += $add; + # re-adjust signs + $x->{sign} = $x->{_m}->{sign}; + $x->{_m}->{sign} = '+'; + return $x->round($a,$p,$r,$y); + } + +sub bsub + { + # (BINT or num_str, BINT or num_str) return num_str + # subtract second arg from first, modify first + my ($self,$x,$y) = objectify(2,@_); -# normalize number -- for internal use -sub norm { #(mantissa, exponent) return fnum_str - local($_, $exp) = @_; - $exp = 0 unless defined $exp; - if ($_ eq 'NaN') { - 'NaN'; - } else { - s/^([+-])0+/$1/; # strip leading zeros - if (length($_) == 1) { - '+0E+0'; - } else { - $exp += length($1) if (s/(0+)$//); # strip trailing zeros - sprintf("%sE%+ld", $_, $exp); - } + trace(@_); + $x->badd($y->bneg()); # badd does not leave internal zeros + $y->bneg(); # refix y, assumes no one reads $y in between + return $x; + } + +sub binc + { + # increment arg by one + my ($self,$x,$a,$p,$r) = objectify(1,@_); + trace(@_); + $x->badd($self->_one())->round($a,$p,$r); + } + +sub bdec + { + # decrement arg by one + my ($self,$x,$a,$p,$r) = objectify(1,@_); + trace(@_); + $x->badd($self->_one('-'))->round($a,$p,$r); + } + +sub blcm + { + # (BINT or num_str, BINT or num_str) return BINT + # does not modify arguments, but returns new object + # Lowest Common Multiplicator + trace(@_); + + my ($self,@arg) = objectify(0,@_); + my $x = $self->new(shift @arg); + while (@arg) { $x = _lcm($x,shift @arg); } + $x; + } + +sub bgcd + { + # (BINT or num_str, BINT or num_str) return BINT + # does not modify arguments, but returns new object + # GCD -- Euclids algorithm Knuth Vol 2 pg 296 + trace(@_); + + my ($self,@arg) = objectify(0,@_); + my $x = $self->new(shift @arg); + while (@arg) { $x = _gcd($x,shift @arg); } + $x; + } + +sub is_zero + { + # return true if arg (BINT or num_str) is zero (array '+', '0') + my $x = shift; $x = $class->new($x) unless ref $x; + #my ($self,$x) = objectify(1,@_); + trace(@_); + return ($x->{sign} ne $nan && $x->{_m}->is_zero()); + } + +sub is_one + { + # return true if arg (BINT or num_str) is +1 (array '+', '1') + # or -1 if signis given + my $x = shift; $x = $class->new($x) unless ref $x; + #my ($self,$x) = objectify(1,@_); + my $sign = $_[2] || '+'; + return ($x->{sign} eq $sign && $x->{_e}->is_zero() && $x->{_m}->is_one()); + } + +sub is_odd + { + # return true if arg (BINT or num_str) is odd or -1 if even + my $x = shift; $x = $class->new($x) unless ref $x; + #my ($self,$x) = objectify(1,@_); + return ($x->{sign} ne $nan && $x->{_e}->is_zero() && $x->{_m}->is_odd()); + } + +sub is_even + { + # return true if arg (BINT or num_str) is even or -1 if odd + my $x = shift; $x = $class->new($x) unless ref $x; + #my ($self,$x) = objectify(1,@_); + return 0 if $x->{sign} eq $nan; # NaN isn't + return 1 if $x->{_m}->is_zero(); # 0 is + return ($x->{_e}->is_zero() && $x->{_m}->is_even()); + } + +sub bmul + { + # multiply two numbers -- stolen from Knuth Vol 2 pg 233 + # (BINT or num_str, BINT or num_str) return BINT + my ($self,$x,$y,$a,$p,$r) = objectify(2,@_); + # trace(@_); + + #print "mul $x->{_m}e$x->{_e} $y->{_m}e$y->{_e}\n"; + return $x->bnan() if (($x->{sign} eq $nan) || ($y->{sign} eq $nan)); + + # print "$x $y\n"; + # aEb * cEd = (a*c)E(b+d) + $x->{_m} = $x->{_m} * $y->{_m}; + #print "m: $x->{_m}\n"; + $x->{_e} = $x->{_e} + $y->{_e}; + #print "e: $x->{_m}\n"; + # adjust sign: + $x->{sign} = $x->{sign} ne $y->{sign} ? '-' : '+'; + #print "s: $x->{sign}\n"; + return $x->round($a,$p,$r,$y); + } + +sub bdiv + { + # (dividend: BFLOAT or num_str, divisor: BFLOAT or num_str) return + # (BFLOAT,BFLOAT) (quo,rem) or BINT (only rem) + my ($self,$x,$y,$a,$p,$r) = objectify(2,@_); + + return wantarray ? ($x->bnan(),bnan()) : $x->bnan() + if ($x->{sign} eq $nan || $y->is_nan() || $y->is_zero()); + + # we need to limit the accuracy to protect against overflow + my ($scale,$mode) = $x->_scale_a($accuracy,$rnd_mode,$a,$r); # ignore $p + my $add = 1; # for proper rounding + my $fallback = 0; + if (!defined $scale) + { + $fallback = 1; $scale = $div_scale; # simulate old behaviour } -} + #print "div_scale $div_scale\n"; + my $lx = $x->{_m}->length(); + $scale = $lx if $lx > $scale; + my $ly = $y->{_m}->length(); + $scale = $ly if $ly > $scale; + #print "scale $scale $lx $ly\n"; + #$scale = $scale - $lx + $ly; + #print "scale $scale\n"; + $scale += $add; # calculate some more digits for proper rounding -# negation -sub fneg { #(fnum_str) return fnum_str - local($_) = fnorm($_[$[]); - vec($_,0,8) ^= ord('+') ^ ord('-') unless $_ eq '+0E+0'; # flip sign - s/^H/N/; - $_; -} + # print "bdiv $x $y scale $scale xl $lx yl $ly\n"; -# absolute value -sub fabs { #(fnum_str) return fnum_str - local($_) = fnorm($_[$[]); - s/^-/+/; # mash sign - $_; -} + return wantarray ? ($x,$self->bzero()) : $x if $x->is_zero(); + + $x->{sign} = $x->{sign} ne $y->sign() ? '-' : '+'; -# multiplication -sub fmul { #(fnum_str, fnum_str) return fnum_str - local($x,$y) = (fnorm($_[$[]),fnorm($_[$[+1])); - if ($x eq 'NaN' || $y eq 'NaN') { - 'NaN'; - } else { - local($xm,$xe) = split('E',$x); - local($ym,$ye) = split('E',$y); - &norm(Math::BigInt::bmul($xm,$ym),$xe+$ye); + # check for / +-1 ( +/- 1E0) + if ($y->is_one()) + { + return wantarray ? ($x,$self->bzero()) : $x; } -} -# addition -sub fadd { #(fnum_str, fnum_str) return fnum_str - local($x,$y) = (fnorm($_[$[]),fnorm($_[$[+1])); - if ($x eq 'NaN' || $y eq 'NaN') { - 'NaN'; - } else { - local($xm,$xe) = split('E',$x); - local($ym,$ye) = split('E',$y); - ($xm,$xe,$ym,$ye) = ($ym,$ye,$xm,$xe) if ($xe < $ye); - &norm(Math::BigInt::badd($ym,$xm.('0' x ($xe-$ye))),$ye); + # a * 10 ** b / c * 10 ** d => a/c * 10 ** (b-d) + #print "self: $self x: $x ref(x) ", ref($x)," m: $x->{_m}\n"; + # my $scale_10 = 10 ** $scale; $x->{_m}->bmul($scale_10); + $x->{_m}->blsft($scale,10); + #print "m: $x->{_m}\n"; + $x->{_m}->bdiv( $y->{_m} ); # a/c + #print "m: $x->{_m}\n"; + #print "e: $x->{_e} $y->{_e}",$scale,"\n"; + $x->{_e}->bsub($y->{_e}); # b-d + #print "e: $x->{_e}\n"; + $x->{_e}->bsub($scale); # correct for 10**scale + #print "e: $x->{_e}\n"; + $x->bnorm(); # remove trailing zeros + + # print "round $x to -$scale (-$add) mode $mode\n"; + #print "$x ",scalar ref($x), "=> $t",scalar ref($t),"\n"; + if ($fallback) + { + $scale -= $add; $x->round($scale,undef,$r); # round to less } -} + else + { + return $x->round($a,$p,$r,$y); + } + if (wantarray) + { + my $rem = $x->copy(); + $rem->bmod($y,$a,$p,$r); + return ($x,$rem->round($scale,undef,$r)) if $fallback; + return ($x,$rem->round($a,$p,$r,$y)); + } + return $x; + } -# subtraction -sub fsub { #(fnum_str, fnum_str) return fnum_str - fadd($_[$[],fneg($_[$[+1])); -} +sub bmod + { + # (dividend: BFLOAT or num_str, divisor: BFLOAT or num_str) return reminder + my ($self,$x,$y,$a,$p,$r) = objectify(2,@_); -# division -# args are dividend, divisor, scale (optional) -# result has at most max(scale, length(dividend), length(divisor)) digits -sub fdiv #(fnum_str, fnum_str[,scale]) return fnum_str -{ - local($x,$y,$scale) = (fnorm($_[$[]),fnorm($_[$[+1]),$_[$[+2]); - if ($x eq 'NaN' || $y eq 'NaN' || $y eq '+0E+0') { - 'NaN'; - } else { - local($xm,$xe) = split('E',$x); - local($ym,$ye) = split('E',$y); - $scale = $div_scale if (!$scale); - $scale = length($xm)-1 if (length($xm)-1 > $scale); - $scale = length($ym)-1 if (length($ym)-1 > $scale); - $scale = $scale + length($ym) - length($xm); - &norm(&round(Math::BigInt::bdiv($xm.('0' x $scale),$ym), - Math::BigInt::babs($ym)), - $xe-$ye-$scale); + return $x->bnan() if ($x->{sign} eq $nan || $y->is_nan() || $y->is_zero()); + return $x->bzero() if $y->is_one(); + + # XXX tels: not done yet + return $x->round($a,$p,$r,$y); + } + +sub bsqrt + { + # calculate square root + # this should use a different test to see wether the accuracy we want is... + my ($self,$x,$a,$p,$r) = objectify(1,@_); + + # we need to limit the accuracy to protect against overflow + my ($scale,$mode) = $x->_scale_a($accuracy,$rnd_mode,$a,$r); # ignore $p + $scale = $div_scale if (!defined $scale); # simulate old behaviour + # print "scale $scale\n"; + + return $x->bnan() if ($x->sign() eq '-') || ($x->sign() eq $nan); + return $x if $x->is_zero() || $x == 1; + + my $len = $x->{_m}->length(); + $scale = $len if $scale < $len; + print "scale $scale\n"; + $scale += 1; # because we need more than $scale to later round + my $e = Math::BigFloat->new("1E-$scale"); # make test variable + return $x->bnan() if $e->sign() eq 'NaN'; + + # print "$scale $e\n"; + + my $gs = Math::BigFloat->new(100); # first guess + my $org = $x->copy(); + + # start with some reasonable guess + #$x *= 10 ** ($len - $org->{_e}); + #$x /= 2; + #my $gs = Math::BigFloat->new(1); + # print "first guess: $gs (x $x)\n"; + + my $diff = $e; + my $y = $x->copy(); + my $two = Math::BigFloat->new(2); + $x = Math::BigFloat->new($x) if ref($x) ne $class; # promote BigInts + # $scale = 2; + while ($diff >= $e) + { + #sleep(1); + return $x->bnan() if $gs->is_zero(); + #my $r = $y / $gs; + #print "$y / $gs = ",$r," ref(\$r) ",ref($r),"\n"; + my $r = $y->copy(); $r->bdiv($gs,$scale); # $scale); + $x = ($r + $gs); + $x->bdiv($two,$scale); # $scale *= 2; + $diff = $x->copy()->bsub($gs)->babs(); + #print "gs: $gs x: $x \n"; + $gs = $x->copy(); + # print "$x $org $scale $gs\n"; + #$gs *= 2; + #$y = $org->copy(); + #$x += $y->bdiv($x, $scale); # need only $gs scale + # $y = $org->copy(); + #$x /= 2; + print "x $x diff $diff $e\n"; } -} + $x->bnorm($scale-1,undef,$mode); + } -# modular division -# args are dividend, divisor -sub fmod #(fnum_str, fnum_str) return fnum_str -{ - local($x,$y) = (fnorm($_[$[]),fnorm($_[$[+1])); - if ($x eq 'NaN' || $y eq 'NaN' || $y eq '+0E+0') { - 'NaN'; - } else { - local($xm,$xe) = split('E',$x); - local($ym,$ye) = split('E',$y); - if ( $xe < $ye ) - { - $ym .= ('0' x ($ye-$xe)); - } - else - { - $xm .= ('0' x ($xe-$ye)); - } - &norm(Math::BigInt::bmod($xm,$ym)); +sub _set + { + # set to a specific 'small' value, internal usage + my $x = shift; + my $v = shift||0; + + $x->{sign} = $nan, return if $v !~ /^[-+]?[0-9]+$/; + $x->{_m}->{value} = [abs($v)]; + $x->{_e}->{value} = [0]; + $x->{sign} = '+'; $x->{sign} = '-' if $v < 0; + return $x; + } + +sub bpow + { + # (BFLOAT or num_str, BFLOAT or num_str) return BFLOAT + # compute power of two numbers, second arg is used as integer + # modifies first argument + + my ($self,$x,$y,$a,$p,$r) = objectify(2,@_); + + return $x->bnan() if $x->{sign} eq $nan || $y->{sign} eq $nan; + return $x->_one() if $y->is_zero(); + return $x if $x->is_one() || $y->is_one(); + my $y1 = $y->as_number(); # make bigint + if ($x == -1) + { + # if $x == -1 and odd/even y => +1/-1 because +-1 ^ (+-1) => +-1 + return $y1->is_odd() ? $x : $x->_set(1); # $x->babs() would work to } -} -# round int $q based on fraction $r/$base using $rnd_mode -sub round { #(int_str, int_str, int_str) return int_str - local($q,$r,$base) = @_; - if ($q eq 'NaN' || $r eq 'NaN') { - 'NaN'; - } elsif ($rnd_mode eq 'trunc') { - $q; # just truncate - } else { - local($cmp) = Math::BigInt::bcmp(Math::BigInt::bmul($r,'+2'),$base); - if ( $cmp < 0 || - ($cmp == 0 && ( - ($rnd_mode eq 'zero' ) || - ($rnd_mode eq '-inf' && (substr($q,$[,1) eq '+')) || - ($rnd_mode eq '+inf' && (substr($q,$[,1) eq '-')) || - ($rnd_mode eq 'even' && $q =~ /[24680]$/ ) || - ($rnd_mode eq 'odd' && $q =~ /[13579]$/ ) ) - ) - ) { - $q; # round down - } else { - Math::BigInt::badd($q, ((substr($q,$[,1) eq '-') ? '-1' : '+1')); - # round up - } + return $x if $x->is_zero() && $y->{sign} eq '+'; # 0**y => 0 (if not y <= 0) + # 0 ** -y => 1 / (0 ** y) => / 0! + return $x->bnan() if $x->is_zero() && $y->{sign} eq '-'; + + # calculate $x->{_m} ** $y and $x->{_e} * $y separately (faster) + $y1->babs(); + $x->{_m}->bpow($y1); + $x->{_e}->bmul($y1); + $x->{sign} = $nan if $x->{_m}->{sign} eq $nan || $x->{_e}->{sign} eq $nan; + $x->bnorm(); + if ($y->{sign} eq '-') + { + # modify $x in place! + my $z = $x->copy(); $x->_set(1); + return $x->bdiv($z,$a,$p,$r); # round in one go (might ignore y's A!) } -} + return $x->round($a,$p,$r,$y); + } + +############################################################################### +# rounding functions + +sub bfround + { + # precision: round to the $Nth digit left (+$n) or right (-$n) from the '.' + # $n == 0 means round to integer + # expects and returns normalized numbers! + my $x = shift; $x = $class->new($x) unless ref $x; -# round the mantissa of $x to $scale digits -sub fround { #(fnum_str, scale) return fnum_str - local($x,$scale) = (fnorm($_[$[]),$_[$[+1]); - if ($x eq 'NaN' || $scale <= 0) { - $x; - } else { - local($xm,$xe) = split('E',$x); - if (length($xm)-1 <= $scale) { - $x; - } else { - &norm(&round(substr($xm,$[,$scale+1), - "+0".substr($xm,$[+$scale+1),"+1"."0" x length(substr($xm,$[+$scale+1))), - $xe+length($xm)-$scale-1); - } + return $x if $x->modify('bfround'); + + my ($scale,$mode) = $x->_scale_p($precision,$rnd_mode,@_); + return $x if !defined $scale; # no-op + + # print "MBF bfround $x to scale $scale mode $mode\n"; + return $x if $x->is_nan() or $x->is_zero(); + + if ($scale < 0) + { + # print "bfround scale $scale e $x->{_e}\n"; + # round right from the '.' + return $x if $x->{_e} >= 0; # nothing to round + $scale = -$scale; # positive for simplicity + my $len = $x->{_m}->length(); # length of mantissa + my $dad = -$x->{_e}; # digits after dot + my $zad = 0; # zeros after dot + $zad = -$len-$x->{_e} if ($x->{_e} < -$len);# for 0.00..00xxx style + # print "scale $scale dad $dad zad $zad len $len\n"; + + # number bsstr len zad dad + # 0.123 123e-3 3 0 3 + # 0.0123 123e-4 3 1 4 + # 0.001 1e-3 1 2 3 + # 1.23 123e-2 3 0 2 + # 1.2345 12345e-4 5 0 4 + + # do not round after/right of the $dad + return $x if $scale > $dad; # 0.123, scale >= 3 => exit + + # round to zero if rounding inside the $zad, but not for last zero like: + # 0.0065, scale -2, round last '0' with following '65' (scale == zad case) + if ($scale < $zad) + { + $x->{_m} = Math::BigInt->new(0); + $x->{_e} = Math::BigInt->new(1); + $x->{sign} = '+'; + return $x; + } + if ($scale == $zad) # for 0.006, scale -2 and trunc + { + $scale = -$len; + } + else + { + # adjust round-point to be inside mantissa + if ($zad != 0) + { + $scale = $scale-$zad; + } + else + { + my $dbd = $len - $dad; $dbd = 0 if $dbd < 0; # digits before dot + $scale = $dbd+$scale; + } + } + # print "round to $x->{_m} to $scale\n"; } -} + else + { + # 123 => 100 means length(123) = 3 - $scale (2) => 1 -# round $x at the 10 to the $scale digit place -sub ffround { #(fnum_str, scale) return fnum_str - local($x,$scale) = (fnorm($_[$[]),$_[$[+1]); - if ($x eq 'NaN') { - 'NaN'; - } else { - local($xm,$xe) = split('E',$x); - if ($xe >= $scale) { - $x; - } else { - $xe = length($xm)+$xe-$scale; - if ($xe < 1) { - '+0E+0'; - } elsif ($xe == 1) { - # The first substr preserves the sign, passing a non- - # normalized "-0" to &round when rounding -0.006 (for - # example), purely so &round won't lose the sign. - &norm(&round(substr($xm,$[,1).'0', - "+0".substr($xm,$[+1), - "+1"."0" x length(substr($xm,$[+1))), $scale); - } else { - &norm(&round(substr($xm,$[,$xe), - "+0".substr($xm,$[+$xe), - "+1"."0" x length(substr($xm,$[+$xe))), $scale); - } - } + # calculate digits before dot + my $dbt = $x->{_m}->length(); $dbt += $x->{_e} if $x->{_e}->sign() eq '-'; + if (($scale > $dbt) && ($dbt < 0)) + { + # if not enough digits before dot, round to zero + $x->{_m} = Math::BigInt->new(0); + $x->{_e} = Math::BigInt->new(1); + $x->{sign} = '+'; + return $x; + } + if (($scale >= 0) && ($dbt == 0)) + { + # 0.49->bfround(1): scale == 1, dbt == 0: => 0.0 + # 0.51->bfround(0): scale == 0, dbt == 0: => 1.0 + # 0.5->bfround(0): scale == 0, dbt == 0: => 0 + # 0.05->bfround(0): scale == 0, dbt == 0: => 0 + # print "$scale $dbt $x->{_m}\n"; + $scale = -$x->{_m}->length(); + } + elsif ($dbt > 0) + { + # correct by subtracting scale + $scale = $dbt - $scale; + } + else + { + $scale = $x->{_m}->length() - $scale; + } } -} + #print "using $scale for $x->{_m} with '$mode'\n"; + # pass sign to bround for '+inf' and '-inf' rounding modes + $x->{_m}->{sign} = $x->{sign}; + $x->{_m}->bround($scale,$mode); + $x->{_m}->{sign} = '+'; # fix sign back + $x->bnorm(); + } + +sub bround + { + # accuracy: preserve $N digits, and overwrite the rest with 0's + my $x = shift; $x = $class->new($x) unless ref $x; + my ($scale,$mode) = $x->_scale_a($accuracy,$rnd_mode,@_); + return $x if !defined $scale; # no-op + + return $x if $x->modify('bround'); + + # print "bround $scale $mode\n"; + # 0 => return all digits, scale < 0 makes no sense + return $x if ($scale <= 0); + return $x if $x->is_nan() or $x->is_zero(); # never round a 0 + + # if $e longer than $m, we have 0.0000xxxyyy style number, and must + # subtract the delta from scale, to simulate keeping the zeros + # -5 +5 => 1; -10 +5 => -4 + my $delta = $x->{_e} + $x->{_m}->length() + 1; + # removed by tlr, since causes problems with fraction tests: + # $scale += $delta if $delta < 0; + + # if we should keep more digits than the mantissa has, do nothing + return $x if $x->{_m}->length() <= $scale; -# Calculate the integer part of $x -sub f2int { #(fnum_str) return inum_str - local($x) = ${$_[$[]}; - if ($x eq 'NaN') { - die "Attempt to take int(NaN)"; - } else { - local($xm,$xe) = split('E',$x); - if ($xe >= 0) { - $xm . '0' x $xe; - } else { - $xe = length($xm)+$xe; - if ($xe <= 1) { - '+0'; - } else { - substr($xm,$[,$xe); - } - } + # pass sign to bround for '+inf' and '-inf' rounding modes + $x->{_m}->{sign} = $x->{sign}; + $x->{_m}->bround($scale,$mode); # round mantissa + $x->{_m}->{sign} = '+'; # fix sign back + return $x->bnorm(); # del trailing zeros gen. by bround() + } + +sub bfloor + { + # return integer less or equal then $x + my ($self,$x,$a,$p,$r) = objectify(1,@_); + + return $x if $x->modify('bfloor'); + + return $x if $x->{sign} !~ /^[+-]$/; # nan, +inf, -inf + + # if $x has digits after dot + if ($x->{_e}->{sign} eq '-') + { + $x->{_m}->brsft(-$x->{_e},10); + $x->{_e}->bzero(); + $x-- if $x->{sign} eq '-'; } -} + return $x->round($a,$p,$r); + } -# compare 2 values returns one of undef, <0, =0, >0 -# returns undef if either or both input value are not numbers -sub fcmp #(fnum_str, fnum_str) return cond_code -{ - local($x, $y) = (fnorm($_[$[]),fnorm($_[$[+1])); - if ($x eq "NaN" || $y eq "NaN") { - undef; - } else { - local($xm,$xe,$ym,$ye) = split('E', $x."E$y"); - if ($xm eq '+0' || $ym eq '+0') { - return $xm <=> $ym; - } - if ( $xe < $ye ) # adjust the exponents to be equal - { - $ym .= '0' x ($ye - $xe); - $ye = $xe; - } - elsif ( $ye < $xe ) # same here - { - $xm .= '0' x ($xe - $ye); - $xe = $ye; - } - return Math::BigInt::cmp($xm,$ym); +sub bceil + { + # return integer greater or equal then $x + my ($self,$x,$a,$p,$r) = objectify(1,@_); + + return $x if $x->modify('bceil'); + return $x if $x->{sign} !~ /^[+-]$/; # nan, +inf, -inf + + # if $x has digits after dot + if ($x->{_e}->{sign} eq '-') + { + $x->{_m}->brsft(-$x->{_e},10); + $x->{_e}->bzero(); + $x++ if $x->{sign} eq '+'; } -} + return $x->round($a,$p,$r); + } + +############################################################################### -# square root by Newtons method. -sub fsqrt { #(fnum_str[, scale]) return fnum_str - local($x, $scale) = (fnorm($_[$[]), $_[$[+1]); - if ($x eq 'NaN' || $x =~ /^-/) { - 'NaN'; - } elsif ($x eq '+0E+0') { - '+0E+0'; - } else { - local($xm, $xe) = split('E',$x); - $scale = $div_scale if (!$scale); - $scale = length($xm)-1 if ($scale < length($xm)-1); - local($gs, $guess) = (1, sprintf("1E%+d", (length($xm)+$xe-1)/2)); - while ($gs < 2*$scale) { - $guess = fmul(fadd($guess,fdiv($x,$guess,$gs*2)),".5"); - $gs *= 2; - } - new Math::BigFloat &fround($guess, $scale); +sub DESTROY + { + # going trough AUTOLOAD for every DESTROY is costly, so avoid it by empty sub + } + +sub AUTOLOAD + { + # make fxxx and bxxx work + # my $self = $_[0]; + my $name = $AUTOLOAD; + + $name =~ s/.*:://; # split package + #print "$name\n"; + if (!method_valid($name)) + { + #no strict 'refs'; + ## try one level up + #&{$class."::SUPER->$name"}(@_); + # delayed load of Carp and avoid recursion + require Carp; + Carp::croak ("Can't call $class\-\>$name, not a valid method"); } -} + no strict 'refs'; + my $bname = $name; $bname =~ s/^f/b/; + *{$class."\:\:$name"} = \&$bname; + &$bname; # uses @_ + } + +sub exponent + { + # return a copy of the exponent + my $self = shift; + $self = $class->new($self) unless ref $self; + + return bnan() if $self->is_nan(); + return $self->{_e}->copy(); + } + +sub mantissa + { + # return a copy of the mantissa + my $self = shift; + $self = $class->new($self) unless ref $self; + + return bnan() if $self->is_nan(); + my $m = $self->{_m}->copy(); # faster than going via bstr() + $m->bneg() if $self->{sign} eq '-'; + + return $m; + } + +sub parts + { + # return a copy of both the exponent and the mantissa + my $self = shift; + $self = $class->new($self) unless ref $self; + + return (bnan(),bnan()) if $self->is_nan(); + my $m = $self->{_m}->copy(); # faster than going via bstr() + $m->bneg() if $self->{sign} eq '-'; + return ($m,$self->{_e}->copy()); + } + +############################################################################## +# private stuff (internal use only) + +sub _one + { + # internal speedup, set argument to 1, or create a +/- 1 + # uses internal knowledge about MBI, thus (bad) + my $self = shift; + my $x = $self->bzero(); + $x->{_m}->{value} = [ 1 ]; $x->{_m}->{sign} = '+'; + $x->{_e}->{value} = [ 0 ]; $x->{_e}->{sign} = '+'; + $x->{sign} = shift || '+'; + return $x; + } + +sub import + { + my $self = shift; + #print "import $self\n"; + for ( my $i = 0; $i < @_ ; $i++ ) + { + if ( $_[$i] eq ':constant' ) + { + # this rest causes overlord er load to step in + # print "overload @_\n"; + overload::constant float => sub { $self->new(shift); }; + splice @_, $i, 1; last; + } + } + # any non :constant stuff is handled by our parent, Exporter + # even if @_ is empty, to give it a chance + #$self->SUPER::import(@_); # does not work (would call MBI) + $self->export_to_level(1,$self,@_); # need this instead + } + +sub bnorm + { + # adjust m and e so that m is smallest possible + # round number according to accuracy and precision settings + my $x = shift; + + return $x if $x->is_nan(); + + my $zeros = $x->{_m}->_trailing_zeros(); # correct for trailing zeros + if ($zeros != 0) + { + $x->{_m}->brsft($zeros,10); $x->{_e} += $zeros; + } + # for something like 0Ey, set y to 1 + $x->{_e}->bzero()->binc() if $x->{_m}->is_zero(); + return $x->SUPER::bnorm(@_); # call MBI bnorm for round + } + +############################################################################## +# internal calculation routines + +sub as_number + { + # return a bigint representation of this BigFloat number + my ($self,$x) = objectify(1,@_); + + my $z; + if ($x->{_e}->is_zero()) + { + $z = $x->{_m}->copy(); + $z->{sign} = $x->{sign}; + return $z; + } + if ($x->{_e} < 0) + { + $x->{_e}->babs(); + my $y = $x->{_m} / ($ten ** $x->{_e}); + $x->{_e}->bneg(); + $y->{sign} = $x->{sign}; + return $y; + } + $z = $x->{_m} * ($ten ** $x->{_e}); + $z->{sign} = $x->{sign}; + return $z; + } + +sub length + { + my $x = shift; $x = $class->new($x) unless ref $x; + + my $len = $x->{_m}->length(); + $len += $x->{_e} if $x->{_e}->sign() eq '+'; + if (wantarray()) + { + my $t = Math::BigInt::bzero(); + $t = $x->{_e}->copy()->babs() if $x->{_e}->sign() eq '-'; + return ($len,$t); + } + return $len; + } 1; __END__ =head1 NAME -Math::BigFloat - Arbitrary length float math package +Math::BigFloat - Arbitrary size floating point math package =head1 SYNOPSIS use Math::BigFloat; - $f = Math::BigFloat->new($string); - - $f->fadd(NSTR) return NSTR addition - $f->fsub(NSTR) return NSTR subtraction - $f->fmul(NSTR) return NSTR multiplication - $f->fdiv(NSTR[,SCALE]) returns NSTR division to SCALE places - $f->fmod(NSTR) returns NSTR modular remainder - $f->fneg() return NSTR negation - $f->fabs() return NSTR absolute value - $f->fcmp(NSTR) return CODE compare undef,<0,=0,>0 - $f->fround(SCALE) return NSTR round to SCALE digits - $f->ffround(SCALE) return NSTR round at SCALEth place - $f->fnorm() return (NSTR) normalize - $f->fsqrt([SCALE]) return NSTR sqrt to SCALE places + + # Number creation + $x = Math::BigInt->new($str); # defaults to 0 + $nan = Math::BigInt->bnan(); # create a NotANumber + $zero = Math::BigInt->bzero();# create a "+0" + + # Testing + $x->is_zero(); # return whether arg is zero or not + $x->is_one(); # return true if arg is +1 + $x->is_one('-'); # return true if arg is -1 + $x->is_odd(); # return true if odd, false for even + $x->is_even(); # return true if even, false for odd + $x->bcmp($y); # compare numbers (undef,<0,=0,>0) + $x->bacmp($y); # compare absolutely (undef,<0,=0,>0) + $x->sign(); # return the sign, either +,- or NaN + + # The following all modify their first argument: + + # set + $x->bzero(); # set $i to 0 + $x->bnan(); # set $i to NaN + + $x->bneg(); # negation + $x->babs(); # absolute value + $x->bnorm(); # normalize (no-op) + $x->bnot(); # two's complement (bit wise not) + $x->binc(); # increment x by 1 + $x->bdec(); # decrement x by 1 + + $x->badd($y); # addition (add $y to $x) + $x->bsub($y); # subtraction (subtract $y from $x) + $x->bmul($y); # multiplication (multiply $x by $y) + $x->bdiv($y); # divide, set $i to quotient + # return (quo,rem) or quo if scalar + + $x->bmod($y); # modulus + $x->bpow($y); # power of arguments (a**b) + $x->blsft($y); # left shift + $x->brsft($y); # right shift + # return (quo,rem) or quo if scalar + + $x->band($y); # bit-wise and + $x->bior($y); # bit-wise inclusive or + $x->bxor($y); # bit-wise exclusive or + $x->bnot(); # bit-wise not (two's complement) + + $x->bround($N); # accuracy: preserver $N digits + $x->bfround($N); # precision: round to the $Nth digit + + # The following do not modify their arguments: + + bgcd(@values); # greatest common divisor + blcm(@values); # lowest common multiplicator + + $x->bstr(); # return string + $x->bsstr(); # return string in scientific notation + + $x->exponent(); # return exponent as BigInt + $x->mantissa(); # return mantissa as BigInt + $x->parts(); # return (mantissa,exponent) as BigInt + + $x->length(); # number of digits (w/o sign and '.') + ($l,$f) = $x->length(); # number of digits, and length of fraction =head1 DESCRIPTION -All basic math operations are overloaded if you declare your big -floats as +All operators (inlcuding basic math operations) are overloaded if you +declare your big floating point numbers as - $float = new Math::BigFloat "2.123123123123123123123123123123123"; + $i = new Math::BigFloat '12_3.456_789_123_456_789E-2'; + +Operations with overloaded operators preserve the arguments, which is +exactly what you expect. + +=head2 Canonical notation + +Input to these routines are either BigFloat objects, or strings of the +following four forms: =over 2 -=item number format +=item * + +C -canonical strings have the form /[+-]\d+E[+-]\d+/ . Input values can -have embedded whitespace. +=item * -=item Error returns 'NaN' +C -An input parameter was "Not a Number" or divide by zero or sqrt of -negative number. +=item * -=item Division is computed to +C -C -digits by default. -Also used for default sqrt scale. +=item * -=item Rounding is performed +C -according to the value of -C<$Math::BigFloat::rnd_mode>: +=back + +all with optional leading and trailing zeros and/or spaces. Additonally, +numbers are allowed to have an underscore between any two digits. + +Empty strings as well as other illegal numbers results in 'NaN'. + +bnorm() on a BigFloat object is now effectively a no-op, since the numbers +are always stored in normalized form. On a string, it creates a BigFloat +object. + +=head2 Output + +Output values are BigFloat objects (normalized), except for bstr() and bsstr(). + +The string output will always have leading and trailing zeros stripped and drop +a plus sign. C will give you always the form with a decimal point, +while C (for scientific) gives you the scientific notation. + + Input bstr() bsstr() + '-0' '0' '0E1' + ' -123 123 123' '-123123123' '-123123123E0' + '00.0123' '0.0123' '123E-4' + '123.45E-2' '1.2345' '12345E-4' + '10E+3' '10000' '1E4' + +Some routines (C, C, C, C, +C) return true or false, while others (C, C) +return either undef, <0, 0 or >0 and are suited for sort. + +Actual math is done by using BigInts to represent the mantissa and exponent. +The sign C is stored separately. The string 'NaN' is used to +represent the result when input arguments are not numbers, as well as +the result of dividing by zero. + +=head2 C, C and C + +C and C return the said parts of the BigFloat +as BigInts such that: + + $m = $x->mantissa(); + $e = $x->exponent(); + $y = $m * ( 10 ** $e ); + print "ok\n" if $x == $y; + +C<< ($m,$e) = $x->parts(); >> is just a shortcut giving you both of them. + +A zero is represented and returned as C<0E1>, B C<0E0> (after Knuth). + +Currently the mantissa is reduced as much as possible, favouring higher +exponents over lower ones (e.g. returning 1e7 instead of 10e6 or 10000000e0). +This might change in the future, so do not depend on it. + +=head2 Accuracy vs. Precision + +See also: L. + +Math::BigFloat supports both precision and accuracy. (here should follow +a short description of both). - trunc truncate the value - zero round towards 0 - +inf round towards +infinity (round up) - -inf round towards -infinity (round down) - even round to the nearest, .5 to the even digit - odd round to the nearest, .5 to the odd digit +Precision: digits after the '.', laber, schwad +Accuracy: Significant digits blah blah -The default is C rounding. +Since things like sqrt(2) or 1/3 must presented with a limited precision lest +a operation consumes all resources, each operation produces no more than +C digits. + +In case the result of one operation has more precision than specified, +it is rounded. The rounding mode taken is either the default mode, or the one +supplied to the operation after the I: + + $x = Math::BigFloat->new(2); + Math::BigFloat::precision(5); # 5 digits max + $y = $x->copy()->bdiv(3); # will give 0.66666 + $y = $x->copy()->bdiv(3,6); # will give 0.666666 + $y = $x->copy()->bdiv(3,6,'odd'); # will give 0.666667 + Math::BigFloat::round_mode('zero'); + $y = $x->copy()->bdiv(3,6); # will give 0.666666 + +=head2 Rounding + +=over 2 + +=item ffround ( +$scale ) rounds to the $scale'th place left from the '.', counting from the dot. The first digit is numbered 1. + +=item ffround ( -$scale ) rounds to the $scale'th place right from the '.', counting from the dot + +=item ffround ( 0 ) rounds to an integer + +=item fround ( +$scale ) preserves accuracy to $scale digits from the left (aka significant digits) and paddes the rest with zeros. If the number is between 1 and -1, the significant digits count from the first non-zero after the '.' + +=item fround ( -$scale ) and fround ( 0 ) are a no-ops =back +All rounding functions take as a second parameter a rounding mode from one of +the following: 'even', 'odd', '+inf', '-inf', 'zero' or 'trunc'. + +The default rounding mode is 'even'. By using +C<< Math::BigFloat::round_mode($rnd_mode); >> you can get and set the default +mode for subsequent rounding. The usage of C<$Math::BigFloat::$rnd_mode> is +no longer supported. + The second parameter to the round functions then overrides the default +temporarily. + +The C<< as_number() >> function returns a BigInt from a Math::BigFloat. It uses +'trunc' as rounding mode to make it equivalent to: + + $x = 2.5; + $y = int($x) + 2; + +You can override this by passing the desired rounding mode as parameter to +C: + + $x = Math::BigFloat->new(2.5); + $y = $x->as_number('odd'); # $y = 3 + +=head1 EXAMPLES + + use Math::BigFloat qw(bstr bint); + # not ready yet + $x = bstr("1234") # string "1234" + $x = "$x"; # same as bstr() + $x = bneg("1234") # BigFloat "-1234" + $x = Math::BigFloat->bneg("1234"); # BigFloat "1234" + $x = Math::BigFloat->babs("-12345"); # BigFloat "12345" + $x = Math::BigFloat->bnorm("-0 00"); # BigFloat "0" + $x = bint(1) + bint(2); # BigFloat "3" + $x = bint(1) + "2"; # ditto (auto-BigFloatify of "2") + $x = bint(1); # BigFloat "1" + $x = $x + 5 / 2; # BigFloat "3" + $x = $x ** 3; # BigFloat "27" + $x *= 2; # BigFloat "54" + $x = new Math::BigFloat; # BigFloat "0" + $x--; # BigFloat "-1" + +=head1 Autocreating constants + +After C all the floating point constants +in the given scope are converted to C. This conversion +happens at compile time. + +In particular + + perl -MMath::BigFloat=:constant -e 'print 2E-100,"\n"' + +prints the value of C<2E-100>. Note that without conversion of +constants the expression 2E-100 will be calculated as normal floating point +number. + +=head1 PERFORMANCE + +Greatly enhanced ;o) +SectionNotReadyYet. + =head1 BUGS -The current version of this module is a preliminary version of the -real thing that is currently (as of perl5.002) under development. +=over 2 + +=item * + +The following does not work yet: + + $m = $x->mantissa(); + $e = $x->exponent(); + $y = $m * ( 10 ** $e ); + print "ok\n" if $x == $y; + +=item * + +There is no fmod() function yet. + +=back + +=head1 CAVEAT + +=over 1 + +=item stringify, bstr() + +Both stringify and bstr() now drop the leading '+'. The old code would return +'+1.23', the new returns '1.23'. See the documentation in L for +reasoning and details. + +=item bdiv + +The following will probably not do what you expect: + + print $c->bdiv(123.456),"\n"; + +It prints both quotient and reminder since print works in list context. Also, +bdiv() will modify $c, so be carefull. You probably want to use + + print $c / 123.456,"\n"; + print scalar $c->bdiv(123.456),"\n"; # or if you want to modify $c + +instead. + +=item Modifying and = + +Beware of: + + $x = Math::BigFloat->new(5); + $y = $x; + +It will not do what you think, e.g. making a copy of $x. Instead it just makes +a second reference to the B object and stores it in $y. Thus anything +that modifies $x will modify $y, and vice versa. + + $x->bmul(2); + print "$x, $y\n"; # prints '10, 10' + +If you want a true copy of $x, use: + + $y = $x->copy(); + +See also the documentation in L regarding C<=>. + +=item bpow + +C now modifies the first argument, unlike the old code which left +it alone and only returned the result. This is to be consistent with +C etc. The first will modify $x, the second one won't: + + print bpow($x,$i),"\n"; # modify $x + print $x->bpow($i),"\n"; # ditto + print $x ** $i,"\n"; # leave $x alone + +=back + +=head1 LICENSE -The printf subroutine does not use the value of -C<$Math::BigFloat::rnd_mode> when rounding values for printing. -Consequently, the way to print rounded values is -to specify the number of digits both as an -argument to C and in the C<%f> printf string, -as follows: +This program is free software; you may redistribute it and/or modify it under +the same terms as Perl itself. - printf "%.3f\n", $bigfloat->ffround(-3); +=head1 AUTHORS -=head1 AUTHOR +Mark Biggar, overloaded interface by Ilya Zakharevich. +Completely rewritten by Tels http://bloodgate.com in 2001. -Mark Biggar -Patches by John Peacock Apr 2001 =cut diff --git a/lib/Math/BigInt.pm b/lib/Math/BigInt.pm index 745e1c6..53d5b11 100644 --- a/lib/Math/BigInt.pm +++ b/lib/Math/BigInt.pm @@ -1,440 +1,2134 @@ -package Math::BigInt; -require Exporter; -@ISA = qw(Exporter); +#!/usr/bin/perl -w + +# mark.biggar@TrustedSysLabs.com +# eay@mincom.com is dead (math::BigInteger) +# see: http://www.cypherspace.org/~adam/rsa/pureperl.html (contacted c. adam +# on 2000/11/13 - but email is dead + +# todo: +# - fully remove funky $# stuff (maybe) +# - use integer; vs 1e7 as base +# - speed issues (XS? Bit::Vector?) +# - split out actual math code to Math::BigNumber + +# Qs: what exactly happens on numify of HUGE numbers? overflow? +# $a = -$a is much slower (making copy of $a) than $a->bneg(), hm!? +# (copy_on_write will help there, but that is not yet implemented) + +# The following hash values are used: +# value: the internal array, base 100000 +# sign : +,-,NaN,+inf,-inf +# _a : accuracy +# _p : precision +# _cow : copy on write: number of objects that share the data (NRY) +# Internally the numbers are stored in an array with at least 1 element, no +# leading zero parts (except the first) and in base 100000 + +# USE_MUL: due to problems on certain os (os390, posix-bc) "* 1e-5" is used +# instead of "/ 1e5" at some places, (marked with USE_MUL). But instead of +# using the reverse only on problematic machines, I used it everytime to avoid +# the costly comparisations. This _should_ work everywhere. Thanx Peter Prymmer -$VERSION='0.02'; +package Math::BigInt; +my $class = "Math::BigInt"; + +$VERSION = 1.35; +use Exporter; +@ISA = qw( Exporter ); +@EXPORT_OK = qw( bneg babs bcmp badd bmul bdiv bmod bnorm bsub + bgcd blcm + bround + blsft brsft band bior bxor bnot bpow bnan bzero + bacmp bstr bsstr binc bdec bint binf bfloor bceil + is_odd is_even is_zero is_one is_nan is_inf sign + length as_number + trace objectify _swap + ); + +#@EXPORT = qw( ); +use vars qw/$rnd_mode $accuracy $precision $div_scale/; +use strict; + +# Inside overload, the first arg is always an object. If the original code had +# it reversed (like $x = 2 * $y), then the third paramater indicates this +# swapping. To make it work, we use a helper routine which not only reswaps the +# params, but also makes a new object in this case. See _swap() for details, +# especially the cases of operators with different classes. + +# For overloaded ops with only one argument we simple use $_[0]->copy() to +# preserve the argument. + +# Thus inheritance of overload operators becomes possible and transparent for +# our subclasses without the need to repeat the entire overload section there. use overload -'+' => sub {new Math::BigInt &badd}, -'-' => sub {new Math::BigInt - $_[2]? bsub($_[1],${$_[0]}) : bsub(${$_[0]},$_[1])}, -'<=>' => sub {$_[2]? bcmp($_[1],${$_[0]}) : bcmp(${$_[0]},$_[1])}, -'cmp' => sub {$_[2]? ($_[1] cmp ${$_[0]}) : (${$_[0]} cmp $_[1])}, -'*' => sub {new Math::BigInt &bmul}, -'/' => sub {new Math::BigInt - $_[2]? scalar bdiv($_[1],${$_[0]}) : - scalar bdiv(${$_[0]},$_[1])}, -'%' => sub {new Math::BigInt - $_[2]? bmod($_[1],${$_[0]}) : bmod(${$_[0]},$_[1])}, -'**' => sub {new Math::BigInt - $_[2]? bpow($_[1],${$_[0]}) : bpow(${$_[0]},$_[1])}, -'neg' => sub {new Math::BigInt &bneg}, -'abs' => sub {new Math::BigInt &babs}, -'<<' => sub {new Math::BigInt - $_[2]? blsft($_[1],${$_[0]}) : blsft(${$_[0]},$_[1])}, -'>>' => sub {new Math::BigInt - $_[2]? brsft($_[1],${$_[0]}) : brsft(${$_[0]},$_[1])}, -'&' => sub {new Math::BigInt &band}, -'|' => sub {new Math::BigInt &bior}, -'^' => sub {new Math::BigInt &bxor}, -'~' => sub {new Math::BigInt &bnot}, -'int' => sub { shift }, +'=' => sub { $_[0]->copy(); }, + +# '+' and '-' do not use _swap, since it is a triffle slower. If you want to +# override _swap (if ever), then override overload of '+' and '-', too! +# for sub it is a bit tricky to keep b: b-a => -a+b +'-' => sub { my $c = $_[0]->copy; $_[2] ? + $c->bneg()->badd($_[1]) : + $c->bsub( $_[1]) }, +'+' => sub { $_[0]->copy()->badd($_[1]); }, + +# some shortcuts for speed (assumes that reversed order of arguments is routed +# to normal '+' and we thus can always modify first arg. If this is changed, +# this breaks and must be adjusted.) +'+=' => sub { $_[0]->badd($_[1]); }, +'-=' => sub { $_[0]->bsub($_[1]); }, +'*=' => sub { $_[0]->bmul($_[1]); }, +'/=' => sub { scalar $_[0]->bdiv($_[1]); }, +'**=' => sub { $_[0]->bpow($_[1]); }, + +'<=>' => sub { $_[2] ? + $class->bcmp($_[1],$_[0]) : + $class->bcmp($_[0],$_[1])}, +'cmp' => sub { + $_[2] ? + $_[1] cmp $_[0]->bstr() : + $_[0]->bstr() cmp $_[1] }, + +'int' => sub { $_[0]->copy(); }, +'neg' => sub { $_[0]->copy()->bneg(); }, +'abs' => sub { $_[0]->copy()->babs(); }, +'~' => sub { $_[0]->copy()->bnot(); }, + +'*' => sub { my @a = ref($_[0])->_swap(@_); $a[0]->bmul($a[1]); }, +'/' => sub { my @a = ref($_[0])->_swap(@_);scalar $a[0]->bdiv($a[1]);}, +'%' => sub { my @a = ref($_[0])->_swap(@_); $a[0]->bmod($a[1]); }, +'**' => sub { my @a = ref($_[0])->_swap(@_); $a[0]->bpow($a[1]); }, +'<<' => sub { my @a = ref($_[0])->_swap(@_); $a[0]->blsft($a[1]); }, +'>>' => sub { my @a = ref($_[0])->_swap(@_); $a[0]->brsft($a[1]); }, + +'&' => sub { my @a = ref($_[0])->_swap(@_); $a[0]->band($a[1]); }, +'|' => sub { my @a = ref($_[0])->_swap(@_); $a[0]->bior($a[1]); }, +'^' => sub { my @a = ref($_[0])->_swap(@_); $a[0]->bxor($a[1]); }, + +# can modify arg of ++ and --, so avoid a new-copy for speed, but don't +# use $_[0]->_one(), it modifies $_[0] to be 1! +'++' => sub { $_[0]->binc() }, +'--' => sub { $_[0]->bdec() }, + +# if overloaded, O(1) instead of O(N) and twice as fast for small numbers +'bool' => sub { + # this kludge is needed for perl prior 5.6.0 since returning 0 here fails :-/ + # v5.6.1 dumps on that: return !$_[0]->is_zero() || undef; :-( + my $t = !$_[0]->is_zero(); + undef $t if $t == 0; + return $t; + }, qw( -"" stringify -0+ numify) # Order of arguments unsignificant +"" bstr +0+ numify), # Order of arguments unsignificant ; -$NaNOK=1; - -sub new { - my($class) = shift; - my($foo) = bnorm(shift); - die "Not a number initialized to Math::BigInt" if !$NaNOK && $foo eq "NaN"; - bless \$foo, $class; -} -sub stringify { "${$_[0]}" } -sub numify { 0 + "${$_[0]}" } # Not needed, additional overhead - # comparing to direct compilation based on - # stringify -sub import { +############################################################################## +# global constants, flags and accessory + +# are NaNs ok? +my $NaNOK=1; +# set to 1 for tracing +my $trace = 0; +# constants for easier life +my $nan = 'NaN'; +my $BASE_LEN = 5; +my $BASE = int("1e".$BASE_LEN); # var for trying to change it to 1e7 +my $RBASE = 1e-5; # see USE_MUL + +# Rounding modes one of 'even', 'odd', '+inf', '-inf', 'zero' or 'trunc' +$rnd_mode = 'even'; +$accuracy = undef; +$precision = undef; +$div_scale = 40; + +sub round_mode + { + # make Class->round_mode() work + my $self = shift || $class; + # shift @_ if defined $_[0] && $_[0] eq $class; + if (defined $_[0]) + { + my $m = shift; + die "Unknown round mode $m" + if $m !~ /^(even|odd|\+inf|\-inf|zero|trunc)$/; + $rnd_mode = $m; return; + } + return $rnd_mode; + } + +sub accuracy + { + # $x->accuracy($a); ref($x) a + # $x->accuracy(); ref($x); + # Class::accuracy(); # not supported + #print "MBI @_ ($class)\n"; + my $x = shift; + + die ("accuracy() needs reference to object as first parameter.") + if !ref $x; + + if (@_ > 0) + { + $x->{_a} = shift; + $x->round() if defined $x->{_a}; + } + return $x->{_a}; + } + +sub precision + { + my $x = shift; + + die ("precision() needs reference to object as first parameter.") + unless ref $x; + + if (@_ > 0) + { + $x->{_p} = shift; + $x->round() if defined $x->{_p}; + } + return $x->{_p}; + } + +sub _scale_a + { + # select accuracy parameter based on precedence, + # used by bround() and bfround(), may return undef for scale (means no op) + my ($x,$s,$m,$scale,$mode) = @_; + $scale = $x->{_a} if !defined $scale; + $scale = $s if (!defined $scale); + $mode = $m if !defined $mode; + return ($scale,$mode); + } + +sub _scale_p + { + # select precision parameter based on precedence, + # used by bround() and bfround(), may return undef for scale (means no op) + my ($x,$s,$m,$scale,$mode) = @_; + $scale = $x->{_p} if !defined $scale; + $scale = $s if (!defined $scale); + $mode = $m if !defined $mode; + return ($scale,$mode); + } + +############################################################################## +# constructors + +sub copy + { + my ($c,$x); + if (@_ > 1) + { + # if two arguments, the first one is the class to "swallow" subclasses + ($c,$x) = @_; + } + else + { + $x = shift; + $c = ref($x); + } + return unless ref($x); # only for objects + + my $self = {}; bless $self,$c; + foreach my $k (keys %$x) + { + if (ref($x->{$k}) eq 'ARRAY') + { + $self->{$k} = [ @{$x->{$k}} ]; + } + elsif (ref($x->{$k}) eq 'HASH') + { + # only one level deep! + foreach my $h (keys %{$x->{$k}}) + { + $self->{$k}->{$h} = $x->{$k}->{$h}; + } + } + elsif (ref($x->{$k})) + { + my $c = ref($x->{$k}); + $self->{$k} = $c->new($x->{$k}); # no copy() due to deep rec + } + else + { + $self->{$k} = $x->{$k}; + } + } + $self; + } + +sub new + { + # create a new BigInts object from a string or another bigint object. + # value => internal array representation + # sign => sign (+/-), or "NaN" + + # the argument could be an object, so avoid ||, && etc on it, this would + # cause costly overloaded code to be called. The only allowed op are ref() + # and definend. + + trace (@_); + my $class = shift; + + my $wanted = shift; # avoid numify call by not using || here + return $class->bzero() if !defined $wanted; # default to 0 + return $class->copy($wanted) if ref($wanted); + + my $self = {}; bless $self, $class; + # handle '+inf', '-inf' first + if ($wanted =~ /^[+-]inf$/) + { + $self->{value} = [ 0 ]; + $self->{sign} = $wanted; + return $self; + } + # split str in m mantissa, e exponent, i integer, f fraction, v value, s sign + my ($mis,$miv,$mfv,$es,$ev) = _split(\$wanted); + if (ref $mis && !ref $miv) + { + # _from_hex + $self->{value} = $mis->{value}; + $self->{sign} = $mis->{sign}; + return $self; + } + if (!ref $mis) + { + die "$wanted is not a number initialized to $class" if !$NaNOK; + #print "NaN 1\n"; + $self->{value} = [ 0 ]; + $self->{sign} = $nan; + return $self; + } + # make integer from mantissa by adjusting exp, then convert to bigint + $self->{sign} = $$mis; # store sign + $self->{value} = [ 0 ]; # for all the NaN cases + my $e = int("$$es$$ev"); # exponent (avoid recursion) + if ($e > 0) + { + my $diff = $e - CORE::length($$mfv); + if ($diff < 0) # Not integer + { + #print "NOI 1\n"; + $self->{sign} = $nan; + } + else # diff >= 0 + { + # adjust fraction and add it to value + # print "diff > 0 $$miv\n"; + $$miv = $$miv . ($$mfv . '0' x $diff); + } + } + else + { + if ($$mfv ne '') # e <= 0 + { + # fraction and negative/zero E => NOI + #print "NOI 2 \$\$mfv '$$mfv'\n"; + $self->{sign} = $nan; + } + elsif ($e < 0) + { + # xE-y, and empty mfv + #print "xE-y\n"; + $e = abs($e); + if ($$miv !~ s/0{$e}$//) # can strip so many zero's? + { + #print "NOI 3\n"; + $self->{sign} = $nan; + } + } + } + $self->{sign} = '+' if $$miv eq '0'; # normalize -0 => +0 + $self->_internal($miv) if $self->{sign} ne $nan; # as internal array + #print "$wanted => $self->{sign} $self->{value}->[0]\n"; + # if any of the globals is set, round to them and thus store them insid $self + $self->round($accuracy,$precision,$rnd_mode) + if defined $accuracy || defined $precision; + return $self; + } + +# some shortcuts for easier life +sub bint + { + # exportable version of new + trace(@_); + return $class->new(@_); + } + +sub bnan + { + # create a bigint 'NaN', if given a BigInt, set it to 'NaN' my $self = shift; - return unless @_; - for ( my $i; $i < @_ ; $i++ ) { - if ( $_[$i] eq ':constant' ) { - overload::constant integer => sub {Math::BigInt->new(shift)}; - splice @_, $i, 1; - last; - } + $self = $class if !defined $self; + if (!ref($self)) + { + my $c = $self; $self = {}; bless $self, $c; + } + return if $self->modify('bnan'); + $self->{value} = [ 0 ]; + $self->{sign} = $nan; + trace('NaN'); + return $self; } - $self->SUPER::import(@_); -} - -$zero = 0; - -# overcome a floating point problem on certain osnames (posix-bc, os390) -BEGIN { - my $x = 100000.0; - my $use_mult = int($x*1e-5)*1e5 == $x ? 1 : 0; -} - -# normalize string form of number. Strip leading zeros. Strip any -# white space and add a sign, if missing. -# Strings that are not numbers result the value 'NaN'. - -sub bnorm { #(num_str) return num_str - local($_) = @_; - s/\s+//g; # strip white space - if (s/^([+-]?)0*(\d+)$/$1$2/) { # test if number - substr($_,$[,0) = '+' unless $1; # Add missing sign - s/^-0/+0/; - $_; - } else { - 'NaN'; - } -} - -# Convert a number from string format to internal base 100000 format. -# Assumes normalized value as input. -sub internal { #(num_str) return int_num_array - local($d) = @_; - ($is,$il) = (substr($d,$[,1),length($d)-2); - substr($d,$[,1) = ''; - ($is, reverse(unpack("a" . ($il%5+1) . ("a5" x ($il/5)), $d))); -} - -# Convert a number from internal base 100000 format to string format. -# This routine scribbles all over input array. -sub external { #(int_num_array) return num_str - $es = shift; - grep($_ > 9999 || ($_ = substr('0000'.$_,-5)), @_); # zero pad - &bnorm(join('', $es, reverse(@_))); # reverse concat and normalize -} - -# Negate input value. -sub bneg { #(num_str) return num_str - local($_) = &bnorm(@_); - return $_ if $_ eq '+0' or $_ eq 'NaN'; - vec($_,0,8) ^= ord('+') ^ ord('-'); - $_; -} - -# Returns the absolute value of the input. -sub babs { #(num_str) return num_str - &abs(&bnorm(@_)); -} - -sub abs { # post-normalized abs for internal use - local($_) = @_; - s/^-/+/; - $_; -} - -# Compares 2 values. Returns one of undef, <0, =0, >0. (suitable for sort) -sub bcmp { #(num_str, num_str) return cond_code - local($x,$y) = (&bnorm($_[$[]),&bnorm($_[$[+1])); - if ($x eq 'NaN') { - undef; - } elsif ($y eq 'NaN') { - undef; - } else { - &cmp($x,$y) <=> 0; - } -} - -sub cmp { # post-normalized compare for internal use - local($cx, $cy) = @_; - - return 0 if ($cx eq $cy); - - local($sx, $sy) = (substr($cx, 0, 1), substr($cy, 0, 1)); - local($ld); - - if ($sx eq '+') { - return 1 if ($sy eq '-' || $cy eq '+0'); - $ld = length($cx) - length($cy); - return $ld if ($ld); - return $cx cmp $cy; - } else { # $sx eq '-' - return -1 if ($sy eq '+'); - $ld = length($cy) - length($cx); - return $ld if ($ld); - return $cy cmp $cx; - } -} - -sub badd { #(num_str, num_str) return num_str - local(*x, *y); ($x, $y) = (&bnorm($_[$[]),&bnorm($_[$[+1])); - if ($x eq 'NaN') { - 'NaN'; - } elsif ($y eq 'NaN') { - 'NaN'; - } else { - @x = &internal($x); # convert to internal form - @y = &internal($y); - local($sx, $sy) = (shift @x, shift @y); # get signs - if ($sx eq $sy) { - &external($sx, &add(*x, *y)); # if same sign add - } else { - ($x, $y) = (&abs($x),&abs($y)); # make abs - if (&cmp($y,$x) > 0) { - &external($sy, &sub(*y, *x)); - } else { - &external($sx, &sub(*x, *y)); - } - } + +sub binf + { + # create a bigint '+-inf', if given a BigInt, set it to '+-inf' + # the sign is either '+', or if given, used from there + my $self = shift; + my $sign = shift; $sign = '+' if !defined $sign || $sign ne '-'; + $self = $class if !defined $self; + if (!ref($self)) + { + my $c = $self; $self = {}; bless $self, $c; + } + return if $self->modify('binf'); + $self->{value} = [ 0 ]; + $self->{sign} = $sign.'inf'; + trace('inf'); + return $self; + } + +sub bzero + { + # create a bigint '+0', if given a BigInt, set it to 0 + my $self = shift; + $self = $class if !defined $self; + if (!ref($self)) + { + my $c = $self; $self = {}; bless $self, $c; + } + return if $self->modify('bzero'); + $self->{value} = [ 0 ]; + $self->{sign} = '+'; + trace('0'); + return $self; + } + +############################################################################## +# string conversation + +sub bsstr + { + # (ref to BFLOAT or num_str ) return num_str + # Convert number from internal format to scientific string format. + # internal format is always normalized (no leading zeros, "-0E0" => "+0E0") + trace(@_); + my ($self,$x) = objectify(1,@_); + + return $x->{sign} if $x->{sign} !~ /^[+-]$/; + my ($m,$e) = $x->parts(); + # can be only '+', so + my $sign = 'e+'; + # MBF: my $s = $e->{sign}; $s = '' if $s eq '-'; my $sep = 'e'.$s; + return $m->bstr().$sign.$e->bstr(); + } + +sub bstr + { + # (ref to BINT or num_str ) return num_str + # Convert number from internal base 100000 format to string format. + # internal format is always normalized (no leading zeros, "-0" => "+0") + trace(@_); + my $x = shift; $x = $class->new($x) unless ref $x; + # my ($self,$x) = objectify(1,@_); + + return $x->{sign} if $x->{sign} !~ /^[+-]$/; + my $ar = $x->{value} || return $nan; # should not happen + my $es = ""; + $es = $x->{sign} if $x->{sign} eq '-'; # get sign, but not '+' + my $l = scalar @$ar; # number of parts + return $nan if $l < 1; # should not happen + # handle first one different to strip leading zeros from it (there are no + # leading zero parts in internal representation) + $l --; $es .= $ar->[$l]; $l--; + # Interestingly, the pre-padd method uses more time + # the old grep variant takes longer (14 to 10 sec) + while ($l >= 0) + { + $es .= substr('0000'.$ar->[$l],-5); # fastest way I could think of + $l--; } -} - -sub bsub { #(num_str, num_str) return num_str - &badd($_[$[],&bneg($_[$[+1])); -} - -# GCD -- Euclids algorithm Knuth Vol 2 pg 296 -sub bgcd { #(num_str, num_str) return num_str - local($x,$y) = (&bnorm($_[$[]),&bnorm($_[$[+1])); - if ($x eq 'NaN' || $y eq 'NaN') { - 'NaN'; - } else { - ($x, $y) = ($y,&bmod($x,$y)) while $y ne '+0'; - $x; - } -} - -# routine to add two base 1e5 numbers -# stolen from Knuth Vol 2 Algorithm A pg 231 -# there are separate routines to add and sub as per Kunth pg 233 -sub add { #(int_num_array, int_num_array) return int_num_array - local(*x, *y) = @_; - $car = 0; - for $x (@x) { - last unless @y || $car; - $x -= 1e5 if $car = (($x += (@y ? shift(@y) : 0) + $car) >= 1e5) ? 1 : 0; - } - for $y (@y) { - last unless $car; - $y -= 1e5 if $car = (($y += $car) >= 1e5) ? 1 : 0; - } - (@x, @y, $car); -} - -# subtract base 1e5 numbers -- stolen from Knuth Vol 2 pg 232, $x > $y -sub sub { #(int_num_array, int_num_array) return int_num_array - local(*sx, *sy) = @_; - $bar = 0; - for $sx (@sx) { - last unless @sy || $bar; - $sx += 1e5 if $bar = (($sx -= (@sy ? shift(@sy) : 0) + $bar) < 0); - } - @sx; -} - -# multiply two numbers -- stolen from Knuth Vol 2 pg 233 -sub bmul { #(num_str, num_str) return num_str - local(*x, *y); ($x, $y) = (&bnorm($_[$[]), &bnorm($_[$[+1])); - if ($x eq 'NaN') { - 'NaN'; - } elsif ($y eq 'NaN') { - 'NaN'; - } else { - @x = &internal($x); - @y = &internal($y); - &external(&mul(*x,*y)); - } -} - -# multiply two numbers in internal representation -# destroys the arguments, supposes that two arguments are different -sub mul { #(*int_num_array, *int_num_array) return int_num_array - local(*x, *y) = (shift, shift); - local($signr) = (shift @x ne shift @y) ? '-' : '+'; - @prod = (); - for $x (@x) { - ($car, $cty) = (0, $[); - for $y (@y) { - $prod = $x * $y + ($prod[$cty] || 0) + $car; - if ($use_mult) { - $prod[$cty++] = - $prod - ($car = int($prod * 1e-5)) * 1e5; + return $es; + } + +sub numify + { + # Make a number from a BigInt object + # old: simple return string and let Perl's atoi() handle the rest + # new: calc because it is faster than bstr()+atoi() + #trace (@_); + #my ($self,$x) = objectify(1,@_); + #return $x->bstr(); # ref($x); + my $x = shift; $x = $class->new($x) unless ref $x; + + return $nan if $x->{sign} eq $nan; + my $fac = 1; $fac = -1 if $x->{sign} eq '-'; + return $fac*$x->{value}->[0] if @{$x->{value}} == 1; # below $BASE + my $num = 0; + foreach (@{$x->{value}}) + { + $num += $fac*$_; $fac *= $BASE; + } + return $num; + } + +############################################################################## +# public stuff (usually prefixed with "b") + +sub sign + { + # return the sign of the number: +/-/NaN + my ($self,$x) = objectify(1,@_); + return $x->{sign}; + } + +sub round + { + # After any operation or when calling round(), the result is rounded by + # regarding the A & P from arguments, local parameters, or globals. + # The result's A or P are set by the rounding, but not inspected beforehand + # (aka only the arguments enter into it). This works because the given + # 'first' argument is both the result and true first argument with unchanged + # A and P settings. + # This does not yet handle $x with A, and $y with P (which should be an + # error). + my $self = shift; + my $a = shift; # accuracy, if given by caller + my $p = shift; # precision, if given by caller + my $r = shift; # round_mode, if given by caller + my @args = @_; # all 'other' arguments (0 for unary, 1 for binary ops) + + unshift @args,$self; # add 'first' argument + + $self = new($self) unless ref($self); # if not object, make one + + # find out class of argument to round + my $c = ref($args[0]); + + # now pick $a or $p, but only if we have got "arguments" + if ((!defined $a) && (!defined $p) && (@args > 0)) + { + foreach (@args) + { + # take the defined one, or if both defined, the one that is smaller + $a = $_->{_a} if (defined $_->{_a}) && (!defined $a || $_->{_a} < $a); + } + if (!defined $a) # if it still is not defined, take p + { + foreach (@args) + { + # take the defined one, or if both defined, the one that is smaller + $p = $_->{_p} if (defined $_->{_p}) && (!defined $p || $_->{_p} < $p); } - else { - $prod[$cty++] = - $prod - ($car = int($prod / 1e5)) * 1e5; + # if none defined, use globals (#2) + if (!defined $p) + { + no strict 'refs'; + my $z = "$c\::accuracy"; $a = $$z; + if (!defined $a) + { + $z = "$c\::precision"; $p = $$z; + } } + } # endif !$a + } # endif !$a || !$P && args > 0 + # for clearity, this is not merged at place (#2) + # now round, by calling fround or ffround: + if (defined $a) + { + $self->{_a} = $a; $self->bround($a,$r); + } + elsif (defined $p) + { + $self->{_p} = $p; $self->bfround($p,$r); + } + return $self->bnorm(); + } + +sub bnorm + { + # (num_str or BINT) return BINT + # Normalize number -- no-op here + my $self = shift; + + return $self; + } + +sub babs + { + # (BINT or num_str) return BINT + # make number absolute, or return absolute BINT from string + #my ($self,$x) = objectify(1,@_); + my $x = shift; $x = $class->new($x) unless ref $x; + return $x if $x->modify('babs'); + # post-normalized abs for internal use (does nothing for NaN) + $x->{sign} =~ s/^-/+/; + $x; + } + +sub bneg + { + # (BINT or num_str) return BINT + # negate number or make a negated number from string + my ($self,$x,$a,$p,$r) = objectify(1,@_); + return $x if $x->modify('bneg'); + # for +0 dont negate (to have always normalized) + return $x if $x->is_zero(); + $x->{sign} =~ tr/+\-/-+/; # does nothing for NaN + # $x->round($a,$p,$r); # changing this makes $x - $y modify $y!! + $x; + } + +sub bcmp + { + # Compares 2 values. Returns one of undef, <0, =0, >0. (suitable for sort) + # (BINT or num_str, BINT or num_str) return cond_code + my ($self,$x,$y) = objectify(2,@_); + return undef if (($x->{sign} eq $nan) || ($y->{sign} eq $nan)); + &cmp($x->{value},$y->{value},$x->{sign},$y->{sign}) <=> 0; + } + +sub bacmp + { + # Compares 2 values, ignoring their signs. + # Returns one of undef, <0, =0, >0. (suitable for sort) + # (BINT, BINT) return cond_code + my ($self,$x,$y) = objectify(2,@_); + return undef if (($x->{sign} eq $nan) || ($y->{sign} eq $nan)); + acmp($x->{value},$y->{value}) <=> 0; + } + +sub badd + { + # add second arg (BINT or string) to first (BINT) (modifies first) + # return result as BINT + trace(@_); + my ($self,$x,$y,$a,$p,$r) = objectify(2,@_); + + return $x if $x->modify('badd'); + return $x->bnan() if (($x->{sign} eq $nan) || ($y->{sign} eq $nan)); + + # for round calls, make array + my @bn = ($a,$p,$r,$y); + # speed: no add for 0+y or x+0 + return $x->bnorm(@bn) if $y->is_zero(); # x+0 + if ($x->is_zero()) # 0+y + { + # make copy, clobbering up x + $x->{value} = [ @{$y->{value}} ]; + $x->{sign} = $y->{sign} || $nan; + return $x->round(@bn); + } + + # shortcuts + my $xv = $x->{value}; + my $yv = $y->{value}; + my ($sx, $sy) = ( $x->{sign}, $y->{sign} ); # get signs + + if ($sx eq $sy) + { + add($xv,$yv); # if same sign, absolute add + $x->{sign} = $sx; + } + else + { + my $a = acmp ($yv,$xv); # absolute compare + if ($a > 0) + { + #print "swapped sub (a=$a)\n"; + &sub($yv,$xv,1); # absolute sub w/ swapped params + $x->{sign} = $sy; + } + elsif ($a == 0) + { + # speedup, if equal, set result to 0 + $x->{value} = [ 0 ]; + $x->{sign} = '+'; + } + else # a < 0 + { + #print "unswapped sub (a=$a)\n"; + &sub($xv, $yv); # absolute sub + $x->{sign} = $sx; } - $prod[$cty] += $car if $car; - $x = shift @prod; - } - ($signr, @x, @prod); -} - -# modulus -sub bmod { #(num_str, num_str) return num_str - (&bdiv(@_))[$[+1]; -} - -sub bdiv { #(dividend: num_str, divisor: num_str) return num_str - local (*x, *y); ($x, $y) = (&bnorm($_[$[]), &bnorm($_[$[+1])); - return wantarray ? ('NaN','NaN') : 'NaN' - if ($x eq 'NaN' || $y eq 'NaN' || $y eq '+0'); - return wantarray ? ('+0',$x) : '+0' if (&cmp(&abs($x),&abs($y)) < 0); - @x = &internal($x); @y = &internal($y); - $srem = $y[$[]; - $sr = (shift @x ne shift @y) ? '-' : '+'; - $car = $bar = $prd = 0; - if (($dd = int(1e5/($y[$#y]+1))) != 1) { - for $x (@x) { - $x = $x * $dd + $car; - if ($use_mult) { - $x -= ($car = int($x * 1e-5)) * 1e5; - } - else { - $x -= ($car = int($x / 1e5)) * 1e5; - } - } - push(@x, $car); $car = 0; - for $y (@y) { - $y = $y * $dd + $car; - if ($use_mult) { - $y -= ($car = int($y * 1e-5)) * 1e5; - } - else { - $y -= ($car = int($y / 1e5)) * 1e5; - } - } } - else { - push(@x, 0); - } - @q = (); ($v2,$v1) = @y[-2,-1]; - $v2 = 0 unless $v2; - while ($#x > $#y) { - ($u2,$u1,$u0) = @x[-3..-1]; - $u2 = 0 unless $u2; - $q = (($u0 == $v1) ? 99999 : int(($u0*1e5+$u1)/$v1)); - --$q while ($v2*$q > ($u0*1e5+$u1-$q*$v1)*1e5+$u2); - if ($q) { - ($car, $bar) = (0,0); - for ($y = $[, $x = $#x-$#y+$[-1; $y <= $#y; ++$y,++$x) { - $prd = $q * $y[$y] + $car; - if ($use_mult) { - $prd -= ($car = int($prd * 1e-5)) * 1e5; - } - else { - $prd -= ($car = int($prd / 1e5)) * 1e5; - } - $x[$x] += 1e5 if ($bar = (($x[$x] -= $prd + $bar) < 0)); - } - if ($x[$#x] < $car + $bar) { - $car = 0; --$q; - for ($y = $[, $x = $#x-$#y+$[-1; $y <= $#y; ++$y,++$x) { - $x[$x] -= 1e5 - if ($car = (($x[$x] += $y[$y] + $car) > 1e5)); - } - } - } - pop(@x); unshift(@q, $q); - } - if (wantarray) { - @d = (); - if ($dd != 1) { - $car = 0; - for $x (reverse @x) { - $prd = $car * 1e5 + $x; - $car = $prd - ($tmp = int($prd / $dd)) * $dd; - unshift(@d, $tmp); - } - } - else { - @d = @x; - } - (&external($sr, @q), &external($srem, @d, $zero)); - } else { - &external($sr, @q); - } -} - -# compute power of two numbers -- stolen from Knuth Vol 2 pg 233 -sub bpow { #(num_str, num_str) return num_str - local(*x, *y); ($x, $y) = (&bnorm($_[$[]), &bnorm($_[$[+1])); - if ($x eq 'NaN') { - 'NaN'; - } elsif ($y eq 'NaN') { - 'NaN'; - } elsif ($x eq '+1') { - '+1'; - } elsif ($x eq '-1') { - &bmod($x,2) ? '-1': '+1'; - } elsif ($y =~ /^-/) { - 'NaN'; - } elsif ($x eq '+0' && $y eq '+0') { - 'NaN'; - } else { - @x = &internal($x); - local(@pow2)=@x; - local(@pow)=&internal("+1"); - local($y1,$res,@tmp1,@tmp2)=(1); # need tmp to send to mul - while ($y ne '+0') { - ($y,$res)=&bdiv($y,2); - if ($res ne '+0') {@tmp=@pow2; @pow=&mul(*pow,*tmp);} - if ($y ne '+0') {@tmp=@pow2;@pow2=&mul(*pow2,*tmp);} - } - &external(@pow); - } -} - -# compute x << y, y >= 0 -sub blsft { #(num_str, num_str) return num_str - &bmul($_[$[], &bpow(2, $_[$[+1])); -} - -# compute x >> y, y >= 0 -sub brsft { #(num_str, num_str) return num_str - &bdiv($_[$[], &bpow(2, $_[$[+1])); -} - -# compute x & y -sub band { #(num_str, num_str) return num_str - local($x,$y,$r,$m,$xr,$yr) = (&bnorm($_[$[]),&bnorm($_[$[+1]),0,1); - if ($x eq 'NaN' || $y eq 'NaN') { - 'NaN'; - } else { - while ($x ne '+0' && $y ne '+0') { - ($x, $xr) = &bdiv($x, 0x10000); - ($y, $yr) = &bdiv($y, 0x10000); - $r = &badd(&bmul(int $xr & $yr, $m), $r); - $m = &bmul($m, 0x10000); + return $x->round(@bn); + } + +sub bsub + { + # (BINT or num_str, BINT or num_str) return num_str + # subtract second arg from first, modify first + my ($self,$x,$y,$a,$p,$r) = objectify(2,@_); + + trace(@_); + return $x if $x->modify('bsub'); + $x->badd($y->bneg()); # badd does not leave internal zeros + $y->bneg(); # refix y, assumes no one reads $y in between + return $x->round($a,$p,$r,$y); + } + +sub binc + { + # increment arg by one + my ($self,$x,$a,$p,$r) = objectify(1,@_); + # my $x = shift; $x = $class->new($x) unless ref $x; my $self = ref($x); + trace(@_); + return $x if $x->modify('binc'); + $x->badd($self->_one())->round($a,$p,$r); + } + +sub bdec + { + # decrement arg by one + my ($self,$x,$a,$p,$r) = objectify(1,@_); + trace(@_); + return $x if $x->modify('bdec'); + $x->badd($self->_one('-'))->round($a,$p,$r); + } + +sub blcm + { + # (BINT or num_str, BINT or num_str) return BINT + # does not modify arguments, but returns new object + # Lowest Common Multiplicator + trace(@_); + + my ($self,@arg) = objectify(0,@_); + my $x = $self->new(shift @arg); + while (@arg) { $x = _lcm($x,shift @arg); } + $x; + } + +sub bgcd + { + # (BINT or num_str, BINT or num_str) return BINT + # does not modify arguments, but returns new object + # GCD -- Euclids algorithm, variant C (Knuth Vol 3, pg 341 ff) + trace(@_); + + my ($self,@arg) = objectify(0,@_); + my $x = $self->new(shift @arg); + while (@arg) + { + #$x = _gcd($x,shift @arg); last if $x->is_one(); # new fast, but is slower + $x = _gcd_old($x,shift @arg); last if $x->is_one(); # old, slow, but faster + } + $x; + } + +sub bmod + { + # modulus + # (BINT or num_str, BINT or num_str) return BINT + my ($self,$x,$y) = objectify(2,@_); + + return $x if $x->modify('bmod'); + (&bdiv($self,$x,$y))[1]; + } + +sub bnot + { + # (num_str or BINT) return BINT + # represent ~x as twos-complement number + my ($self,$x) = objectify(1,@_); + return $x if $x->modify('bnot'); + $x->bneg(); $x->bdec(); # was: bsub(-1,$x);, time it someday + $x; + } + +sub is_zero + { + # return true if arg (BINT or num_str) is zero (array '+', '0') + #my ($self,$x) = objectify(1,@_); + #trace(@_); + my $x = shift; $x = $class->new($x) unless ref $x; + return (@{$x->{value}} == 1) && ($x->{sign} eq '+') + && ($x->{value}->[0] == 0); + } + +sub is_nan + { + # return true if arg (BINT or num_str) is NaN + #my ($self,$x) = objectify(1,@_); + #trace(@_); + my $x = shift; $x = $class->new($x) unless ref $x; + return ($x->{sign} eq $nan); + } + +sub is_inf + { + # return true if arg (BINT or num_str) is +-inf + #my ($self,$x) = objectify(1,@_); + #trace(@_); + my $x = shift; $x = $class->new($x) unless ref $x; + my $sign = shift || ''; + + return $x->{sign} =~ /^[+-]inf/ if $sign eq ''; + return $x->{sign} =~ /^[$sign]inf/; + } + +sub is_one + { + # return true if arg (BINT or num_str) is +1 (array '+', '1') + # or -1 if signis given + #my ($self,$x) = objectify(1,@_); + my $x = shift; $x = $class->new($x) unless ref $x; + my $sign = shift || '+'; #$_[2] || '+'; + return (@{$x->{value}} == 1) && ($x->{sign} eq $sign) + && ($x->{value}->[0] == 1); + } + +sub is_odd + { + # return true when arg (BINT or num_str) is odd, false for even + my $x = shift; $x = $class->new($x) unless ref $x; + #my ($self,$x) = objectify(1,@_); + return (($x->{sign} ne $nan) && ($x->{value}->[0] & 1)); + } + +sub is_even + { + # return true when arg (BINT or num_str) is even, false for odd + my $x = shift; $x = $class->new($x) unless ref $x; + #my ($self,$x) = objectify(1,@_); + return (($x->{sign} ne $nan) && (!($x->{value}->[0] & 1))); + } + +sub bmul + { + # multiply two numbers -- stolen from Knuth Vol 2 pg 233 + # (BINT or num_str, BINT or num_str) return BINT + my ($self,$x,$y,$a,$p,$r) = objectify(2,@_); + #print "$self bmul $x ",ref($x)," $y ",ref($y),"\n"; + trace(@_); + return $x if $x->modify('bmul'); + return $x->bnan() if (($x->{sign} eq $nan) || ($y->{sign} eq $nan)); + + mul($x,$y); # do actual math + return $x->round($a,$p,$r,$y); + } + +sub bdiv + { + # (dividend: BINT or num_str, divisor: BINT or num_str) return + # (BINT,BINT) (quo,rem) or BINT (only rem) + trace(@_); + my ($self,$x,$y,$a,$p,$r) = objectify(2,@_); + + return $x if $x->modify('bdiv'); + + # NaN? + return wantarray ? ($x->bnan(),bnan()) : $x->bnan() + if ($x->{sign} eq $nan || $y->{sign} eq $nan || $y->is_zero()); + + # 0 / something + return wantarray ? ($x,$self->bzero()) : $x if $x->is_zero(); + + # Is $x in the interval [0, $y) ? + my $cmp = acmp($x->{value},$y->{value}); + if (($cmp < 0) and ($x->{sign} eq $y->{sign})) + { + return $x->bzero() unless wantarray; + my $t = $x->copy(); # make copy first, because $x->bzero() clobbers $x + return ($x->bzero(),$t); + } + elsif ($cmp == 0) + { + # shortcut, both are the same, so set to +/- 1 + $x->_one( ($x->{sign} ne $y->{sign} ? '-' : '+') ); + return $x unless wantarray; + return ($x,$self->bzero()); + } + + # calc new sign and in case $y == +/- 1, return $x + $x->{sign} = ($x->{sign} ne $y->{sign} ? '-' : '+'); + # check for / +-1 (cant use $y->is_one due to '-' + if ((@{$y->{value}} == 1) && ($y->{value}->[0] == 1)) + { + return wantarray ? ($x,$self->bzero()) : $x; + } + + # call div here + my $rem = $self->bzero(); + $rem->{sign} = $y->{sign}; + ($x->{value},$rem->{value}) = div($x->{value},$y->{value}); + # do not leave rest "-0"; + $rem->{sign} = '+' if (@{$rem->{value}} == 1) && ($rem->{value}->[0] == 0); + if (($x->{sign} eq '-') and (!$rem->is_zero())) + { + $x->bdec(); + } + $x->round($a,$p,$r,$y); + if (wantarray) + { + $rem->round($a,$p,$r,$x,$y); + return ($x,$y-$rem) if $x->{sign} eq '-'; # was $x,$rem + return ($x,$rem); + } + return $x; + } + +sub bpow + { + # (BINT or num_str, BINT or num_str) return BINT + # compute power of two numbers -- stolen from Knuth Vol 2 pg 233 + # modifies first argument + #trace(@_); + my ($self,$x,$y,$a,$p,$r) = objectify(2,@_); + + return $x if $x->modify('bpow'); + + return $x->bnan() if $x->{sign} eq $nan || $y->{sign} eq $nan; + return $x->_one() if $y->is_zero(); + return $x if $x->is_one() || $y->is_one(); + if ($x->{sign} eq '-' && @{$x->{value}} == 1 && $x->{value}->[0] == 1) + { + # if $x == -1 and odd/even y => +1/-1 + return $y->is_odd() ? $x : $x->_set(1); # $x->babs() would work to + # my Casio FX-5500L has here a bug, -1 ** 2 is -1, but -1 * -1 is 1 LOL + } + # shortcut for $x ** 2 + if ($y->{sign} eq '+' && @{$y->{value}} == 1 && $y->{value}->[0] == 2) + { + return $x->bmul($x)->bround($a,$p,$r); + } + # 1 ** -y => 1 / (1**y), so do test for negative $y after above's clause + return $x->bnan() if $y->{sign} eq '-'; + return $x if $x->is_zero(); # 0**y => 0 (if not y <= 0) + + # tels: 10**x is special (actually 100**x etc is special, too) but not here + #if ((@{$x->{value}} == 1) && ($x->{value}->[0] == 10)) + # { + # # 10**2 + # my $yi = int($y); my $yi5 = int($yi/5); + # $x->{value} = []; + # my $v = $x->{value}; + # if ($yi5 > 0) + # { + # # $x->{value}->[$yi5-1] = 0; # pre-padd array (no use) + # for (my $i = 0; $i < $yi5; $i++) + # { + # $v->[$i] = 0; + # } + # } + # push @{$v}, int( '1'.'0' x ($yi % 5)); + # if ($x->{sign} eq '-') + # { + # $x->{sign} = $y->is_odd() ? '-' : '+'; # -10**2 = 100, -10**3 = -1000 + # } + # return $x; + # } + + # based on the assumption that shifting in base 10 is fast, and that bpow() + # works faster if numbers are small: we count trailing zeros (this step is + # O(1)..O(N), but in case of O(N) we save much more time), stripping them + # out of the multiplication, and add $count * $y zeros afterwards: + # 300 ** 3 == 300*300*300 == 3*3*3 . '0' x 2 * 3 == 27 . '0' x 6 + my $zeros = $x->_trailing_zeros(); + if ($zeros > 0) + { + $x->brsft($zeros,10); # remove zeros + $x->bpow($y); # recursion (will not branch into here again) + $zeros = $y * $zeros; # real number of zeros to add + $x->blsft($zeros,10); + return $x; + } + + my $pow2 = $self->_one(); + my $y1 = $class->new($y); + my ($res); + while (!$y1->is_one()) + { + #print "bpow: p2: $pow2 x: $x y: $y1 r: $res\n"; + #print "len ",$x->length(),"\n"; + ($y1,$res)=&bdiv($y1,2); + if (!$res->is_zero()) { &bmul($pow2,$x); } + if (!$y1->is_zero()) { &bmul($x,$x); } + } + #print "bpow: e p2: $pow2 x: $x y: $y1 r: $res\n"; + &bmul($x,$pow2) if (!$pow2->is_one()); + #print "bpow: e p2: $pow2 x: $x y: $y1 r: $res\n"; + return $x->round($a,$p,$r); + } + +sub blsft + { + # (BINT or num_str, BINT or num_str) return BINT + # compute x << y, base n, y >= 0 + my ($self,$x,$y,$n) = objectify(2,@_); + + return $x if $x->modify('blsft'); + return $x->bnan() if ($x->{sign} !~ /^[+-]$/ || $y->{sign} !~ /^[+-]$/); + + $n = 2 if !defined $n; return $x if $n == 0; + return $x->bnan() if $n < 0 || $y->{sign} eq '-'; + if ($n != 10) + { + $x->bmul( $self->bpow($n, $y) ); + } + else + { + # shortcut (faster) for shifting by 10) since we are in base 10eX + # multiples of 5: + my $src = scalar @{$x->{value}}; # source + my $len = $y->numify(); # shift-len as normal int + my $rem = $len % 5; # reminder to shift + my $dst = $src + int($len/5); # destination + + my $v = $x->{value}; # speed-up + my $vd; # further speedup + #print "src $src:",$v->[$src]||0," dst $dst:",$v->[$dst]||0," rem $rem\n"; + $v->[$src] = 0; # avoid first ||0 for speed + while ($src >= 0) + { + $vd = $v->[$src]; $vd = '00000'.$vd; + #print "s $src d $dst '$vd' "; + $vd = substr($vd,-5+$rem,5-$rem); + #print "'$vd' "; + $vd .= $src > 0 ? substr('00000'.$v->[$src-1],-5,$rem) : '0' x $rem; + #print "'$vd' "; + $vd = substr($vd,-5,5) if length($vd) > 5; + #print "'$vd'\n"; + $v->[$dst] = int($vd); + $dst--; $src--; + } + # set lowest parts to 0 + while ($dst >= 0) { $v->[$dst--] = 0; } + # fix spurios last zero element + splice @$v,-1 if $v->[-1] == 0; + #print "elems: "; my $i = 0; + #foreach (reverse @$v) { print "$i $_ "; $i++; } print "\n"; + # old way: $x->bmul( $self->bpow($n, $y) ); + } + return $x; + } + +sub brsft + { + # (BINT or num_str, BINT or num_str) return BINT + # compute x >> y, base n, y >= 0 + my ($self,$x,$y,$n) = objectify(2,@_); + + return $x if $x->modify('brsft'); + return $x->bnan() if ($x->{sign} !~ /^[+-]$/ || $y->{sign} !~ /^[+-]$/); + + $n = 2 if !defined $n; return $x->bnan() if $n <= 0 || $y->{sign} eq '-'; + if ($n != 10) + { + scalar bdiv($x, $self->bpow($n, $y)); + } + else + { + # shortcut (faster) for shifting by 10) + # multiples of 5: + my $dst = 0; # destination + my $src = $y->numify(); # as normal int + my $rem = $src % 5; # reminder to shift + $src = int($src / 5); # source + my $len = scalar @{$x->{value}} - $src; # elems to go + my $v = $x->{value}; # speed-up + if ($rem == 0) + { + splice (@$v,0,$src); # even faster, 38.4 => 39.3 + } + else + { + my $vd; + $v->[scalar @$v] = 0; # avoid || 0 test inside loop + while ($dst < $len) + { + $vd = '00000'.$v->[$src]; + #print "$dst $src '$vd' "; + $vd = substr($vd,-5,5-$rem); + #print "'$vd' "; + $src++; + $vd = substr('00000'.$v->[$src],-$rem,$rem) . $vd; + #print "'$vd1' "; + #print "'$vd'\n"; + $vd = substr($vd,-5,5) if length($vd) > 5; + $v->[$dst] = int($vd); + $dst++; + } + splice (@$v,$dst) if $dst > 0; # kill left-over array elems + pop @$v if $v->[-1] == 0; # kill last element + } # else rem == 0 + # old way: scalar bdiv($x, $self->bpow($n, $y)); + } + return $x; + } + +sub band + { + #(BINT or num_str, BINT or num_str) return BINT + # compute x & y + trace(@_); + my ($self,$x,$y) = objectify(2,@_); + + return $x if $x->modify('band'); + + return $x->bnan() if ($x->{sign} !~ /^[+-]$/ || $y->{sign} !~ /^[+-]$/); + return $x->bzero() if $y->is_zero(); + my $r = $self->bzero(); my $m = new Math::BigInt 1; my ($xr,$yr); + my $x10000 = new Math::BigInt (0x10000); + my $y1 = copy(ref($x),$y); # make copy + while (!$x->is_zero() && !$y1->is_zero()) + { + ($x, $xr) = bdiv($x, $x10000); + ($y1, $yr) = bdiv($y1, $x10000); + $r->badd( bmul( new Math::BigInt ( $xr->numify() & $yr->numify()), $m )); + $m->bmul($x10000); + } + $x = $r; + } + +sub bior + { + #(BINT or num_str, BINT or num_str) return BINT + # compute x | y + trace(@_); + my ($self,$x,$y) = objectify(2,@_); + + return $x if $x->modify('bior'); + + return $x->bnan() if ($x->{sign} !~ /^[+-]$/ || $y->{sign} !~ /^[+-]$/); + return $x if $y->is_zero(); + my $r = $self->bzero(); my $m = new Math::BigInt 1; my ($xr,$yr); + my $x10000 = new Math::BigInt (0x10000); + my $y1 = copy(ref($x),$y); # make copy + while (!$x->is_zero() || !$y1->is_zero()) + { + ($x, $xr) = bdiv($x,$x10000); + ($y1, $yr) = bdiv($y1,$x10000); + $r->badd( bmul( new Math::BigInt ( $xr->numify() | $yr->numify()), $m )); + $m->bmul($x10000); + } + $x = $r; + } + +sub bxor + { + #(BINT or num_str, BINT or num_str) return BINT + # compute x ^ y + my ($self,$x,$y) = objectify(2,@_); + + return $x if $x->modify('bxor'); + + return $x->bnan() if ($x->{sign} eq $nan || $y->{sign} eq $nan); + return $x if $y->is_zero(); + return $x->bzero() if $x == $y; # shortcut + my $r = $self->bzero(); my $m = new Math::BigInt 1; my ($xr,$yr); + my $x10000 = new Math::BigInt (0x10000); + my $y1 = copy(ref($x),$y); # make copy + while (!$x->is_zero() || !$y1->is_zero()) + { + ($x, $xr) = bdiv($x, $x10000); + ($y1, $yr) = bdiv($y1, $x10000); + $r->badd( bmul( new Math::BigInt ( $xr->numify() ^ $yr->numify()), $m )); + $m->bmul($x10000); + } + $x = $r; + } + +sub length + { + my ($self,$x) = objectify(1,@_); + + return (_digits($x->{value}), 0) if wantarray; + _digits($x->{value}); + } + +sub digit + { + # return the nth digit, negative values count backward + my $x = shift; + my $n = shift || 0; + + my $len = $x->length(); + + $n = $len+$n if $n < 0; # -1 last, -2 second-to-last + $n = abs($n); # if negatives are to big + $len--; $n = $len if $n > $len; # n to big? + + my $elem = int($n / 5); # which array element + my $digit = $n % 5; # which digit in this element + $elem = '0000'.$x->{value}->[$elem]; # get element padded with 0's + return substr($elem,-$digit-1,1); + } + +sub _trailing_zeros + { + # return the amount of trailing zeros in $x + my $x = shift; + $x = $class->new($x) unless ref $x; + + return 0 if $x->is_zero() || $x->is_nan(); + # check each array elem in _m for having 0 at end as long as elem == 0 + # Upon finding a elem != 0, stop + my $zeros = 0; my $elem; + foreach my $e (@{$x->{value}}) + { + if ($e != 0) + { + $elem = "$e"; # preserve x + $elem =~ s/.*?(0*$)/$1/; # strip anything not zero + $zeros *= 5; # elems * 5 + $zeros += CORE::length($elem); # count trailing zeros + last; # early out + } + $zeros ++; # real else branch: 50% slower! + } + return $zeros; + } + +sub bsqrt + { + my ($self,$x) = objectify(1,@_); + + return $x->bnan() if $x->{sign} =~ /\-|$nan/; # -x or NaN => NaN + return $x->bzero() if $x->is_zero(); # 0 => 0 + return $x if $x == 1; # 1 => 1 + + my $y = $x->copy(); # give us one more digit accur. + my $l = int($x->length()/2); + + $x->bzero(); + $x->binc(); # keep ref($x), but modify it + $x *= 10 ** $l; + + # print "x: $y guess $x\n"; + + my $last = $self->bzero(); + while ($last != $x) + { + $last = $x; + $x += $y / $x; + $x /= 2; + } + return $x; + } + +sub exponent + { + # return a copy of the exponent (here always 0, NaN or 1 for $m == 0) + my ($self,$x) = objectify(1,@_); + + return bnan() if $x->is_nan(); + my $e = $class->bzero(); + return $e->binc() if $x->is_zero(); + $e += $x->_trailing_zeros(); + return $e; + } + +sub mantissa + { + # return a copy of the mantissa (here always $self) + my ($self,$x) = objectify(1,@_); + + return bnan() if $x->is_nan(); + my $m = $x->copy(); + # that's inefficient + my $zeros = $m->_trailing_zeros(); + $m /= 10 ** $zeros if $zeros != 0; + return $m; + } + +sub parts + { + # return a copy of both the exponent and the mantissa (here 0 and self) + my $self = shift; + $self = $class->new($self) unless ref $self; + + return ($self->mantissa(),$self->exponent()); + } + +############################################################################## +# rounding functions + +sub bfround + { + # precision: round to the $Nth digit left (+$n) or right (-$n) from the '.' + # $n == 0 => round to integer + my $x = shift; $x = $class->new($x) unless ref $x; + my ($scale,$mode) = $x->_scale_p($precision,$rnd_mode,@_); + return $x if !defined $scale; # no-op + + # no-op for BigInts if $n <= 0 + return $x if $scale <= 0; + + $x->bround( $x->length()-$scale, $mode); + } + +sub _scan_for_nonzero + { + my $x = shift; + my $pad = shift; + + my $len = $x->length(); + return 0 if $len == 1; # '5' is trailed by invisible zeros + my $follow = $pad - 1; + return 0 if $follow > $len || $follow < 1; + #print "checking $x $r\n"; + # old, slow way checking string for non-zero characters + my $r = substr ("$x",-$follow); + return 1 if $r =~ /[^0]/; return 0; + + # faster way checking array contents; it is actually not faster (even in a + # rounding-only-shoutout, so I leave the simpler code in) + #my $rem = $follow % 5; my $div = $follow / 5; my $v = $x->{value}; + # pad with zeros and extract + #print "last part : ",'00000'.$v->[$div]," $rem = '"; + #print substr('00000'.$v->[$div],-$rem,5),"'\n"; + #my $r1 = substr ('00000'.$v->[$div],-$rem,5); + #print "$r1\n"; + #return 1 if $r1 =~ /[^0]/; + # + #for (my $j = $div-1; $j >= 0; $j --) + # { + # #print "part $v->[$j]\n"; + # return 1 if $v->[$j] != 0; + # } + #return 0; + } + +sub fround + { + # to make life easier for switch between MBF and MBI (autoload fxxx() + # like MBF does for bxxx()?) + my $x = shift; + return $x->bround(@_); + } + +sub bround + { + # accuracy: +$n preserve $n digits from left, + # -$n preserve $n digits from right (f.i. for 0.1234 style in MBF) + # no-op for $n == 0 + # and overwrite the rest with 0's, return normalized number + # do not return $x->bnorm(), but $x + my $x = shift; $x = $class->new($x) unless ref $x; + my ($scale,$mode) = $x->_scale_a($accuracy,$rnd_mode,@_); + return $x if !defined $scale; # no-op + + # print "MBI round: $x to $scale $mode\n"; + # -scale means what? tom? hullo? -$scale needed by MBF round, but what for? + return $x if $x->is_nan() || $x->is_zero() || $scale == 0; + + # we have fewer digits than we want to scale to + my $len = $x->length(); + # print "$len $scale\n"; + return $x if $len < abs($scale); + + # count of 0's to pad, from left (+) or right (-): 9 - +6 => 3, or |-6| => 6 + my ($pad,$digit_round,$digit_after); + $pad = $len - $scale; + $pad = abs($scale)+1 if $scale < 0; + $digit_round = '0'; $digit_round = $x->digit($pad) if $pad < $len; + $digit_after = '0'; $digit_after = $x->digit($pad-1) if $pad > 0; + # print "r $x: pos:$pad l:$len s:$scale r:$digit_round a:$digit_after m: $mode\n"; + + # in case of 01234 we round down, for 6789 up, and only in case 5 we look + # closer at the remaining digits of the original $x, remember decision + my $round_up = 1; # default round up + $round_up -- if + ($mode eq 'trunc') || # trunc by round down + ($digit_after =~ /[01234]/) || # round down anyway, + # 6789 => round up + ($digit_after eq '5') && # not 5000...0000 + ($x->_scan_for_nonzero($pad) == 0) && + ( + ($mode eq 'even') && ($digit_round =~ /[24680]/) || + ($mode eq 'odd') && ($digit_round =~ /[13579]/) || + ($mode eq '+inf') && ($x->{sign} eq '-') || + ($mode eq '-inf') && ($x->{sign} eq '+') || + ($mode eq 'zero') # round down if zero, sign adjusted below + ); + # allow rounding one place left of mantissa + #print "$pad $len $scale\n"; + # this is triggering warnings, and buggy for $scale < 0 + #if (-$scale != $len) + { + # split mantissa at $scale and then pad with zeros + my $s5 = int($pad / 5); + my $i = 0; + while ($i < $s5) + { + $x->{value}->[$i++] = 0; # replace with 5 x 0 + } + $x->{value}->[$s5] = '00000'.$x->{value}->[$s5]; # pad with 0 + my $rem = $pad % 5; # so much left over + if ($rem > 0) + { + #print "remainder $rem\n"; + #print "elem $x->{value}->[$s5]\n"; + substr($x->{value}->[$s5],-$rem,$rem) = '0' x $rem; # stamp w/ '0' + } + $x->{value}->[$s5] = int ($x->{value}->[$s5]); # str '05' => int '5' + } + if ($round_up) # what gave test above? + { + $pad = $len if $scale < 0; # tlr: whack 0.51=>1.0 + # modify $x in place, undef, undef to avoid rounding + $x->badd( Math::BigInt->new($x->{sign}.'1'.'0'x$pad), + undef,undef ); + # str creation much faster than 10 ** something + } + $x; + } + +sub bfloor + { + # return integer less or equal then number, since it is already integer, + # always returns $self + my ($self,$x,$a,$p,$r) = objectify(1,@_); + + # not needed: return $x if $x->modify('bfloor'); + + return $x->round($a,$p,$r); + } + +sub bceil + { + # return integer greater or equal then number, since it is already integer, + # always returns $self + my ($self,$x,$a,$p,$r) = objectify(1,@_); + + # not needed: return $x if $x->modify('bceil'); + + return $x->round($a,$p,$r); + } + +############################################################################## +# private stuff (internal use only) + +sub trace + { + # print out a number without using bstr (avoid deep recurse) for trace/debug + return unless $trace; + + my ($package,$file,$line,$sub) = caller(1); + print "'$sub' called from '$package' line $line:\n "; + + foreach my $x (@_) + { + if (!defined $x) + { + print "undef, "; next; + } + if (!ref($x)) + { + print "'$x' "; next; + } + next if (ref($x) ne "HASH"); + print "$x->{sign} "; + foreach (@{$x->{value}}) + { + print "$_ "; + } + print ", "; + } + print "\n"; + } + +sub _set + { + # internal set routine to set X fast to an integer value < [+-]100000 + my $self = shift; + my $wanted = shift || 0; + + $self->{sign} = $nan, return if $wanted !~ /^[+-]?[0-9]+$/; + $self->{sign} = '-'; $self->{sign} = '+' if $wanted >= 0; + $self->{value} = [ abs($wanted) ]; + return $self; + } + +sub _one + { + # internal speedup, set argument to 1, or create a +/- 1 + my $self = shift; + my $x = $self->bzero(); $x->{value} = [ 1 ]; $x->{sign} = shift || '+'; $x; + } + +sub _swap + { + # Overload will swap params if first one is no object ref so that the first + # one is always an object ref. In this case, third param is true. + # This routine is to overcome the effect of scalar,$object creating an object + # of the class of this package, instead of the second param $object. This + # happens inside overload, when the overload section of this package is + # inherited by sub classes. + # For overload cases (and this is used only there), we need to preserve the + # args, hence the copy(). + # You can override this method in a subclass, the overload section will call + # $object->_swap() to make sure it arrives at the proper subclass, with some + # exceptions like '+' and '-'. + + # object, (object|scalar) => preserve first and make copy + # scalar, object => swapped, re-swap and create new from first + # (using class of second object, not $class!!) + my $self = shift; # for override in subclass + #print "swap $self 0:$_[0] 1:$_[1] 2:$_[2]\n"; + if ($_[2]) + { + my $c = ref ($_[0]) || $class; # fallback $class should not happen + return ( $c->new($_[1]), $_[0] ); + } + else + { + return ( $_[0]->copy(), $_[1] ); + } + } + +sub objectify + { + # check for strings, if yes, return objects instead + + # the first argument is number of args objectify() should look at it will + # return $count+1 elements, the first will be a classname. This is because + # overloaded '""' calls bstr($object,undef,undef) and this would result in + # useless objects beeing created and thrown away. So we cannot simple loop + # over @_. If the given count is 0, all arguments will be used. + + # If the second arg is a ref, use it as class. + # If not, try to use it as classname, unless undef, then use $class + # (aka Math::BigInt). The latter shouldn't happen,though. + + # caller: gives us: + # $x->badd(1); => ref x, scalar y + # Class->badd(1,2); => classname x (scalar), scalar x, scalar y + # Class->badd( Class->(1),2); => classname x (scalar), ref x, scalar y + # Math::BigInt::badd(1,2); => scalar x, scalar y + # In the last case we check number of arguments to turn it silently into + # $class,1,2. (We can not take '1' as class ;o) + # badd($class,1) is not supported (it should, eventually, try to add undef) + # currently it tries 'Math::BigInt' + 1, which will not work. + + trace(@_); + my $count = abs(shift || 0); + + #print caller(),"\n"; + + my @a; # resulting array + if (ref $_[0]) + { + # okay, got object as first + $a[0] = ref $_[0]; + } + else + { + # nope, got 1,2 (Class->xxx(1) => Class,1 and not supported) + $a[0] = $class; + #print "@_\n"; sleep(1); + $a[0] = shift if $_[0] =~ /^[A-Z].*::/; # classname as first? + } + #print caller(),"\n"; + # print "Now in objectify, my class is today $a[0]\n"; + my $k; + if ($count == 0) + { + while (@_) + { + $k = shift; + if (!ref($k)) + { + $k = $a[0]->new($k); + } + elsif (ref($k) ne $a[0]) + { + # foreign object, try to convert to integer + $k->can('as_number') ? $k = $k->as_number() : $k = $a[0]->new($k); } - $r; - } -} - -# compute x | y -sub bior { #(num_str, num_str) return num_str - local($x,$y,$r,$m,$xr,$yr) = (&bnorm($_[$[]),&bnorm($_[$[+1]),0,1); - if ($x eq 'NaN' || $y eq 'NaN') { - 'NaN'; - } else { - while ($x ne '+0' || $y ne '+0') { - ($x, $xr) = &bdiv($x, 0x10000); - ($y, $yr) = &bdiv($y, 0x10000); - $r = &badd(&bmul(int $xr | $yr, $m), $r); - $m = &bmul($m, 0x10000); + push @a,$k; + } + } + else + { + while ($count > 0) + { + #print "$count\n"; + $count--; + $k = shift; + if (!ref($k)) + { + $k = $a[0]->new($k); + } + elsif (ref($k) ne $a[0]) + { + # foreign object, try to convert to integer + $k->can('as_number') ? $k = $k->as_number() : $k = $a[0]->new($k); } - $r; - } -} - -# compute x ^ y -sub bxor { #(num_str, num_str) return num_str - local($x,$y,$r,$m,$xr,$yr) = (&bnorm($_[$[]),&bnorm($_[$[+1]),0,1); - if ($x eq 'NaN' || $y eq 'NaN') { - 'NaN'; - } else { - while ($x ne '+0' || $y ne '+0') { - ($x, $xr) = &bdiv($x, 0x10000); - ($y, $yr) = &bdiv($y, 0x10000); - $r = &badd(&bmul(int $xr ^ $yr, $m), $r); - $m = &bmul($m, 0x10000); + push @a,$k; + } + push @a,@_; # return other params, too + } + #my $i = 0; + #foreach (@a) + # { + # print "o $i $a[0]\n" if $i == 0; + # print "o $i ",ref($_),"\n" if $i != 0; $i++; + # } + #print "objectify done: would return ",scalar @a," values\n"; + #print caller(1),"\n" unless wantarray; + die "$class objectify needs list context" unless wantarray; + @a; + } + +sub import + { + my $self = shift; + #print "import $self @_\n"; + for ( my $i = 0; $i < @_ ; $i++ ) + { + if ( $_[$i] eq ':constant' ) + { + # this rest causes overlord er load to step in + overload::constant integer => sub { $self->new(shift) }; + splice @_, $i, 1; last; + } + } + # any non :constant stuff is handled by our parent, Exporter + # even if @_ is empty, to give it a chance + #$self->SUPER::import(@_); # does not work + $self->export_to_level(1,$self,@_); # need this instead + } + +sub _internal + { + # (ref to self, ref to string) return ref to num_array + # Convert a number from string format to internal base 100000 format. + # Assumes normalized value as input. + my ($s,$d) = @_; + my $il = CORE::length($$d)-1; + # these leaves '00000' instead of int 0 and will be corrected after any op + $s->{value} = [ reverse(unpack("a" . ($il%5+1) . ("a5" x ($il/5)), $$d)) ]; + $s; + } + +sub _strip_zeros + { + # internal normalization function that strips leading zeros from the array + # args: ref to array + #trace(@_); + my $s = shift; + + my $cnt = scalar @$s; # get count of parts + my $i = $cnt-1; + #print "strip: cnt $cnt i $i\n"; + # '0', '3', '4', '0', '0', + # 0 1 2 3 4 + # cnt = 5, i = 4 + # i = 4 + # i = 3 + # => fcnt = cnt - i (5-2 => 3, cnt => 5-1 = 4, throw away from 4th pos) + # >= 1: skip first part (this can be zero) + while ($i > 0) { last if $s->[$i] != 0; $i--; } + $i++; splice @$s,$i if ($i < $cnt); # $i cant be 0 + return $s; + } + +sub _from_hex + { + # convert a (ref to) big hex string to BigInt, return undef for error + my $hs = shift; + + my $x = Math::BigInt->bzero(); + return $x->bnan() if $$hs !~ /^[\-\+]?0x[0-9A-Fa-f]+$/; + + my $mul = Math::BigInt->bzero(); $mul++; + my $x65536 = Math::BigInt->new(65536); + + my $len = CORE::length($$hs)-2; my $sign = '+'; + if ($$hs =~ /^\-/) + { + $sign = '-'; $len--; + } + $len = int($len/4); # 4-digit parts, w/o '0x' + my $val; my $i = -4; + while ($len >= 0) + { + $val = substr($$hs,$i,4); + $val =~ s/^[\-\+]?0x// if $len == 0; # for last part only because + $val = hex($val); # hex does not like wrong chars + # print "$val ",substr($$hs,$i,4),"\n"; + $i -= 4; $len --; + $x += $mul * $val if $val != 0; + $mul *= $x65536 if $len >= 0; # skip last mul + } + $x->{sign} = $sign if !$x->is_zero(); + return $x; + } + +sub _from_bin + { + # convert a (ref to) big binary string to BigInt, return undef for error + my $bs = shift; + + my $x = Math::BigInt->bzero(); + return $x->bnan() if $$bs !~ /^[\-\+]?0b[01]+$/; + + my $mul = Math::BigInt->bzero(); $mul++; + my $x256 = Math::BigInt->new(256); + + my $len = CORE::length($$bs)-2; my $sign = '+'; + if ($$bs =~ /^\-/) + { + $sign = '-'; $len--; + } + $len = int($len/8); # 8-digit parts, w/o '0b' + my $val; my $i = -8; + while ($len >= 0) + { + $val = substr($$bs,$i,8); + $val =~ s/^[\-\+]?0b// if $len == 0; # for last part only + #$val = oct('0b'.$val); # does not work on Perl prior 5.6.0 + $val = ('0' x (8-CORE::length($val))).$val if CORE::length($val) < 8; + $val = ord(pack('B8',$val)); + # print "$val ",substr($$bs,$i,16),"\n"; + $i -= 8; $len --; + $x += $mul * $val if $val != 0; + $mul *= $x256 if $len >= 0; # skip last mul + } + $x->{sign} = $sign if !$x->is_zero(); + return $x; + } + +sub _split + { + # (ref to num_str) return num_str + # internal, take apart a string and return the pieces + my $x = shift; + + # pre-parse input + $$x =~ s/^\s+//g; # strip white space at front + $$x =~ s/\s+$//g; # strip white space at end + #$$x =~ s/\s+//g; # strip white space (no longer) + return if $$x eq ""; + + return _from_hex($x) if $$x =~ /^[\-\+]?0x/; # hex string + return _from_bin($x) if $$x =~ /^[\-\+]?0b/; # binary string + + return if $$x !~ /^[\-\+]?\.?[0-9]/; + + $$x =~ s/(\d)_(\d)/$1$2/g; # strip underscores between digits + $$x =~ s/(\d)_(\d)/$1$2/g; # do twice for 1_2_3 + + # some possible inputs: + # 2.1234 # 0.12 # 1 # 1E1 # 2.134E1 # 434E-10 # 1.02009E-2 + # .2 # 1_2_3.4_5_6 # 1.4E1_2_3 # 1e3 # +.2 + + #print "input: '$$x' "; + my ($m,$e) = split /[Ee]/,$$x; + $e = '0' if !defined $e || $e eq ""; + # print "m '$m' e '$e'\n"; + # sign,value for exponent,mantint,mantfrac + my ($es,$ev,$mis,$miv,$mfv); + # valid exponent? + if ($e =~ /^([+-]?)0*(\d+)$/) # strip leading zeros + { + $es = $1; $ev = $2; + #print "'$m' '$e' e: $es $ev "; + # valid mantissa? + return if $m eq '.' || $m eq ''; + my ($mi,$mf) = split /\./,$m; + $mi = '0' if !defined $mi; + $mi .= '0' if $mi =~ /^[\-\+]?$/; + $mf = '0' if !defined $mf || $mf eq ''; + if ($mi =~ /^([+-]?)0*(\d+)$/) # strip leading zeros + { + $mis = $1||'+'; $miv = $2; + #print "$mis $miv"; + # valid, existing fraction part of mantissa? + return unless ($mf =~ /^(\d*?)0*$/); # strip trailing zeros + $mfv = $1; + #print " split: $mis $miv . $mfv E $es $ev\n"; + return (\$mis,\$miv,\$mfv,\$es,\$ev); + } + } + return; # NaN, not a number + } + +sub _digits + { + # computer number of digits in bigint, minus the sign + # int() because add/sub leaves sometimes strings (like '00005') instead of + # int ('5') in this place, causing length to fail + my $cx = shift; + + #print "len: ",(@$cx-1)*5+CORE::length(int($cx->[-1])),"\n"; + return (@$cx-1)*5+CORE::length(int($cx->[-1])); + } + +sub as_number + { + # an object might be asked to return itself as bigint on certain overloaded + # operations, this does exactly this, so that sub classes can simple inherit + # it or override with their own integer conversion routine + my $self = shift; + + return $self->copy(); + } + +############################################################################## +# internal calculation routines + +sub acmp + { + # internal absolute post-normalized compare (ignore signs) + # ref to array, ref to array, return <0, 0, >0 + # arrays must have at least on entry, this is not checked for + + my ($cx, $cy) = @_; + + #print "$cx $cy\n"; + my ($i,$a,$x,$y,$k); + # calculate length based on digits, not parts + $x = _digits($cx); $y = _digits($cy); + # print "length: ",($x-$y),"\n"; + return $x-$y if ($x - $y); # if different in length + #print "full compare\n"; + $i = 0; $a = 0; + # first way takes 5.49 sec instead of 4.87, but has the early out advantage + # so grep is slightly faster, but more unflexible. hm. $_ instead if $k + # yields 5.6 instead of 5.5 sec huh? + # manual way (abort if unequal, good for early ne) + my $j = scalar @$cx - 1; + while ($j >= 0) + { + # print "$cx->[$j] $cy->[$j] $a",$cx->[$j]-$cy->[$j],"\n"; + last if ($a = $cx->[$j] - $cy->[$j]); $j--; + } + return $a; + # while it early aborts, it is even slower than the manual variant + #grep { return $a if ($a = $_ - $cy->[$i++]); } @$cx; + # grep way, go trough all (bad for early ne) + #grep { $a = $_ - $cy->[$i++]; } @$cx; + #return $a; + } + +sub cmp + { + # post-normalized compare for internal use (honors signs) + # ref to array, ref to array, return < 0, 0, >0 + my ($cx,$cy,$sx,$sy) = @_; + + #return 0 if (is0($cx,$sx) && is0($cy,$sy)); + + if ($sx eq '+') + { + return 1 if $sy eq '-'; # 0 check handled above + return acmp($cx,$cy); + } + else + { + # $sx eq '-' + return -1 if ($sy eq '+'); + return acmp($cy,$cx); + } + return 0; # equal + } + +sub add + { + # (ref to int_num_array, ref to int_num_array) + # routine to add two base 1e5 numbers + # stolen from Knuth Vol 2 Algorithm A pg 231 + # there are separate routines to add and sub as per Kunth pg 233 + # This routine clobbers up array x, but not y. + + my ($x,$y) = @_; + + # for each in Y, add Y to X and carry. If after that, something is left in + # X, foreach in X add carry to X and then return X, carry + # Trades one "$j++" for having to shift arrays, $j could be made integer + # but this would impose a limit to number-length to 2**32. + my $i; my $car = 0; my $j = 0; + for $i (@$y) + { + $x->[$j] -= $BASE + if $car = (($x->[$j] += $i + $car) >= $BASE) ? 1 : 0; + $j++; + } + while ($car != 0) + { + $x->[$j] -= $BASE if $car = (($x->[$j] += $car) >= $BASE) ? 1 : 0; $j++; + } + } + +sub sub + { + # (ref to int_num_array, ref to int_num_array) + # subtract base 1e5 numbers -- stolen from Knuth Vol 2 pg 232, $x > $y + # subtract Y from X (X is always greater/equal!) by modifiyng x in place + my ($sx,$sy,$s) = @_; + + my $car = 0; my $i; my $j = 0; + if (!$s) + { + #print "case 2\n"; + for $i (@$sx) + { + last unless defined $sy->[$j] || $car; + #print "x: $i y: $sy->[$j] c: $car\n"; + $i += $BASE if $car = (($i -= ($sy->[$j] || 0) + $car) < 0); $j++; + #print "x: $i y: $sy->[$j-1] c: $car\n"; + } + # might leave leading zeros, so fix that + _strip_zeros($sx); + return $sx; + } + else + { + #print "case 1 (swap)\n"; + for $i (@$sx) + { + last unless defined $sy->[$j] || $car; + #print "$sy->[$j] $i $car => $sx->[$j]\n"; + $sy->[$j] += $BASE + if $car = (($sy->[$j] = $i-($sy->[$j]||0) - $car) < 0); + #print "$sy->[$j] $i $car => $sy->[$j]\n"; + $j++; + } + # might leave leading zeros, so fix that + _strip_zeros($sy); + return $sy; + } + } + +sub mul + { + # (BINT, BINT) return nothing + # multiply two numbers in internal representation + # modifies first arg, second needs not be different from first + my ($x,$y) = @_; + + $x->{sign} = $x->{sign} ne $y->{sign} ? '-' : '+'; + my @prod = (); my ($prod,$car,$cty,$xi,$yi); + my $xv = $x->{value}; + my $yv = $y->{value}; + # since multiplying $x with $x fails, make copy in this case + $yv = [@$xv] if "$xv" eq "$yv"; + for $xi (@$xv) + { + $car = 0; $cty = 0; + for $yi (@$yv) + { + $prod = $xi * $yi + ($prod[$cty] || 0) + $car; + $prod[$cty++] = + $prod - ($car = int($prod * 1e-5)) * $BASE; # see USE_MUL + } + $prod[$cty] += $car if $car; # need really to check for 0? + $xi = shift @prod; + } + push @$xv, @prod; + _strip_zeros($x->{value}); + # normalize (handled last to save check for $y->is_zero() + $x->{sign} = '+' if @$xv == 1 && $xv->[0] == 0; # not is_zero due to '-' + } + +sub div + { + # ref to array, ref to array, modify first array and return reminder if + # in list context + # does no longer handle sign + my ($x,$yorg) = @_; + my ($car,$bar,$prd,$dd,$xi,$yi,@q,$v2,$v1); + + my (@d,$tmp,$q,$u2,$u1,$u0); + + $car = $bar = $prd = 0; + + my $y = [ @$yorg ]; + if (($dd = int($BASE/($y->[-1]+1))) != 1) + { + for $xi (@$x) + { + $xi = $xi * $dd + $car; + $xi -= ($car = int($xi * $RBASE)) * $BASE; # see USE_MUL + } + push(@$x, $car); $car = 0; + for $yi (@$y) + { + $yi = $yi * $dd + $car; + $yi -= ($car = int($yi * $RBASE)) * $BASE; # see USE_MUL + } + } + else + { + push(@$x, 0); + } + @q = (); ($v2,$v1) = @$y[-2,-1]; + $v2 = 0 unless $v2; + while ($#$x > $#$y) + { + ($u2,$u1,$u0) = @$x[-3..-1]; + $u2 = 0 unless $u2; + print "oups v1 is 0, u0: $u0 $y->[-2] $y->[-1] l ",scalar @$y,"\n" + if $v1 == 0; + $q = (($u0 == $v1) ? 99999 : int(($u0*$BASE+$u1)/$v1)); + --$q while ($v2*$q > ($u0*1e5+$u1-$q*$v1)*$BASE+$u2); + if ($q) + { + ($car, $bar) = (0,0); + for ($yi = 0, $xi = $#$x-$#$y-1; $yi <= $#$y; ++$yi,++$xi) + { + $prd = $q * $y->[$yi] + $car; + $prd -= ($car = int($prd * $RBASE)) * $BASE; # see USE_MUL + $x->[$xi] += 1e5 if ($bar = (($x->[$xi] -= $prd + $bar) < 0)); } - $r; + if ($x->[-1] < $car + $bar) + { + $car = 0; --$q; + for ($yi = 0, $xi = $#$x-$#$y-1; $yi <= $#$y; ++$yi,++$xi) + { + $x->[$xi] -= 1e5 + if ($car = (($x->[$xi] += $y->[$yi] + $car) > $BASE)); + } + } + } + pop(@$x); unshift(@q, $q); } -} + if (wantarray) + { + @d = (); + if ($dd != 1) + { + $car = 0; + for $xi (reverse @$x) + { + $prd = $car * $BASE + $xi; + $car = $prd - ($tmp = int($prd / $dd)) * $dd; # see USE_MUL + unshift(@d, $tmp); + } + } + else + { + @d = @$x; + } + @$x = @q; + _strip_zeros($x); + _strip_zeros(\@d); + return ($x,\@d); + } + @$x = @q; + _strip_zeros($x); + return $x; + } + +sub _lcm + { + # (BINT or num_str, BINT or num_str) return BINT + # does modify first argument + # LCM + + my $x = shift; my $ty = shift; + return $x->bnan() if ($x->{sign} eq $nan) || ($ty->{sign} eq $nan); + return $x * $ty / bgcd($x,$ty); + } + +sub _gcd_old + { + # (BINT or num_str, BINT or num_str) return BINT + # does modify first arg + # GCD -- Euclids algorithm E, Knuth Vol 2 pg 296 + trace(@_); + + my $x = shift; my $ty = $class->new(shift); # preserve y + return $x->bnan() if ($x->{sign} eq $nan) || ($ty->{sign} eq $nan); + + while (!$ty->is_zero()) + { + ($x, $ty) = ($ty,bmod($x,$ty)); + } + $x; + } + +sub _gcd + { + # (BINT or num_str, BINT or num_str) return BINT + # does not modify arguments + # GCD -- Euclids algorithm, variant L (Lehmer), Knuth Vol 3 pg 347 ff + # unfortunately, it is slower and also seems buggy (the A=0, B=1, C=1, D=0 + # case..) + trace(@_); + + my $u = $class->new(shift); my $v = $class->new(shift); # preserve u,v + return $u->bnan() if ($u->{sign} eq $nan) || ($v->{sign} eq $nan); + + $u->babs(); $v->babs(); # Euclid is valid for |u| and |v| + + my ($U,$V,$A,$B,$C,$D,$T,$Q); # single precision variables + my ($t); # multiprecision variables + + while ($v > $BASE) + { + #sleep 1; + ($u,$v) = ($v,$u) if ($u < $v); # make sure that u >= v + #print "gcd: $u $v\n"; + # step L1, initialize + $A = 1; $B = 0; $C = 0; $D = 1; + $U = $u->{value}->[-1]; # leading digits of u + $V = $v->{value}->[-1]; # leading digits of v + + # step L2, test quotient + if (($V + $C != 0) && ($V + $D != 0)) # div by zero => go to L4 + { + $Q = int(($U + $A)/($V + $C)); # quotient + #print "L1 A=$A B=$B C=$C D=$D U=$U V=$V Q=$Q\n"; + # do L3? (div by zero => go to L4) + while ($Q == int(($U + $B)/($V + $D))) + { + # step L3, emulate Euclid + #print "L3a A=$A B=$B C=$C D=$D U=$U V=$V Q=$Q\n"; + $T = $A - $Q*$C; $A = $C; $C = $T; + $T = $B - $Q*$D; $B = $D; $D = $T; + $T = $U - $Q*$V; $U = $V; $V = $T; + last if ($V + $D == 0) || ($V + $C == 0); # div by zero => L4 + $Q = int(($U + $A)/($V + $C)); # quotient for next test + #print "L3b A=$A B=$B C=$C D=$D U=$U V=$V Q=$Q\n"; + } + } + # step L4, multiprecision step + # was if ($B == 0) + # in case A = 0, B = 1, C = 0 and D = 1, this case would simple swap u & v + # and loop endless. Not sure why this happens, Knuth does not make a + # remark about this special case. bug? + if (($B == 0) || (($A == 0) && ($C == 1) && ($D == 0))) + { + #print "L4b1: u=$u v=$v\n"; + ($u,$v) = ($v,bmod($u,$v)); + #$t = $u % $v; $u = $v->copy(); $v = $t; + #print "L4b12: u=$u v=$v\n"; + } + else + { + #print "L4b: $u $v $A $B $C $D\n"; + $t = $A*$u + $B*$v; $v *= $D; $v += $C*$u; $u = $t; + #print "L4b2: $u $v\n"; + } + } # back to L1 -# represent ~x as twos-complement number -sub bnot { #(num_str) return num_str - &bsub(-1,$_[$[]); -} + return _gcd_old($u,$v) if $v != 1; # v too small + return $v; # 1 + } + +############################################################################### +# this method return 0 if the object can be modified, or 1 for not +# We use a fast use constant statement here, to avoid costly calls. Subclasses +# may override it with special code (f.i. Math::BigInt::Constant does so) + +use constant modify => 0; + +#sub modify +# { +# my $self = shift; +# my $method = shift; +# print "original $self modify by $method\n"; +# return 0; # $self; +# } 1; __END__ @@ -446,84 +2140,514 @@ Math::BigInt - Arbitrary size integer math package =head1 SYNOPSIS use Math::BigInt; - $i = Math::BigInt->new($string); - - $i->bneg return BINT negation - $i->babs return BINT absolute value - $i->bcmp(BINT) return CODE compare numbers (undef,<0,=0,>0) - $i->badd(BINT) return BINT addition - $i->bsub(BINT) return BINT subtraction - $i->bmul(BINT) return BINT multiplication - $i->bdiv(BINT) return (BINT,BINT) division (quo,rem) just quo if scalar - $i->bmod(BINT) return BINT modulus - $i->bgcd(BINT) return BINT greatest common divisor - $i->bnorm return BINT normalization - $i->blsft(BINT) return BINT left shift - $i->brsft(BINT) return (BINT,BINT) right shift (quo,rem) just quo if scalar - $i->band(BINT) return BINT bit-wise and - $i->bior(BINT) return BINT bit-wise inclusive or - $i->bxor(BINT) return BINT bit-wise exclusive or - $i->bnot return BINT bit-wise not + + # Number creation + $x = Math::BigInt->new($str); # defaults to 0 + $nan = Math::BigInt->bnan(); # create a NotANumber + $zero = Math::BigInt->bzero();# create a "+0" + + # Testing + $x->is_zero(); # return whether arg is zero or not + $x->is_nan(); # return whether arg is NaN or not + $x->is_one(); # return true if arg is +1 + $x->is_one('-'); # return true if arg is -1 + $x->is_odd(); # return true if odd, false for even + $x->is_even(); # return true if even, false for odd + $x->bcmp($y); # compare numbers (undef,<0,=0,>0) + $x->bacmp($y); # compare absolutely (undef,<0,=0,>0) + $x->sign(); # return the sign, either +,- or NaN + $x->digit($n); # return the nth digit, counting from right + $x->digit(-$n); # return the nth digit, counting from left + + # The following all modify their first argument: + + # set + $x->bzero(); # set $x to 0 + $x->bnan(); # set $x to NaN + + $x->bneg(); # negation + $x->babs(); # absolute value + $x->bnorm(); # normalize (no-op) + $x->bnot(); # two's complement (bit wise not) + $x->binc(); # increment x by 1 + $x->bdec(); # decrement x by 1 + + $x->badd($y); # addition (add $y to $x) + $x->bsub($y); # subtraction (subtract $y from $x) + $x->bmul($y); # multiplication (multiply $x by $y) + $x->bdiv($y); # divide, set $x to quotient + # return (quo,rem) or quo if scalar + + $x->bmod($y); # modulus (x % y) + $x->bpow($y); # power of arguments (x ** y) + $x->blsft($y); # left shift + $x->brsft($y); # right shift + $x->blsft($y,$n); # left shift, by base $n (like 10) + $x->brsft($y,$n); # right shift, by base $n (like 10) + + $x->band($y); # bitwise and + $x->bior($y); # bitwise inclusive or + $x->bxor($y); # bitwise exclusive or + $x->bnot(); # bitwise not (two's complement) + + $x->bsqrt(); # calculate square-root + + $x->round($A,$P,$round_mode); # round to accuracy or precision using mode $r + $x->bround($N); # accuracy: preserve $N digits + $x->bfround($N); # round to $Nth digit, no-op for BigInts + + # The following do not modify their arguments in BigInt, but do in BigFloat: + $x->bfloor(); # return integer less or equal than $x + $x->bceil(); # return integer greater or equal than $x + + # The following do not modify their arguments: + + bgcd(@values); # greatest common divisor + blcm(@values); # lowest common multiplicator + + $x->bstr(); # normalized string + $x->bsstr(); # normalized string in scientific notation + $x->length(); # return number of digits in number + ($x,$f) = $x->length(); # length of number and length of fraction part + + $x->exponent(); # return exponent as BigInt + $x->mantissa(); # return mantissa as BigInt + $x->parts(); # return (mantissa,exponent) as BigInt =head1 DESCRIPTION -All basic math operations are overloaded if you declare your big -integers as +All operators (inlcuding basic math operations) are overloaded if you +declare your big integers as - $i = new Math::BigInt '123 456 789 123 456 789'; + $i = new Math::BigInt '123_456_789_123_456_789'; +Operations with overloaded operators preserve the arguments which is +exactly what you expect. =over 2 =item Canonical notation -Big integer value are strings of the form C with leading +Big integer values are strings of the form C with leading zeros suppressed. + '-0' canonical value '-0', normalized '0' + ' -123_123_123' canonical value '-123123123' + '1_23_456_7890' canonical value '1234567890' + =item Input -Input values to these routines may be strings of the form -C. +Input values to these routines may be either Math::BigInt objects or +strings of the form C. + +You can include one underscore between any two digits. + +This means integer values like 1.01E2 or even 1000E-2 are also accepted. +Non integer values result in NaN. + +Math::BigInt::new() defaults to 0, while Math::BigInt::new('') results +in 'NaN'. + +bnorm() on a BigInt object is now effectively a no-op, since the numbers +are always stored in normalized form. On a string, it creates a BigInt +object. =item Output -Output values always always in canonical form +Output values are BigInt objects (normalized), except for bstr(), which +returns a string in normalized form. +Some routines (C, C, C, C, +C) return true or false, while others (C, C) +return either undef, <0, 0 or >0 and are suited for sort. =back -Actual math is done in an internal format consisting of an array -whose first element is the sign (/^[+-]$/) and whose remaining -elements are base 100000 digits with the least significant digit first. -The string 'NaN' is used to represent the result when input arguments -are not numbers, as well as the result of dividing by zero. +=head2 Rounding -=head1 EXAMPLES +Only C is defined for BigInts, for further details of rounding see +L. + +=over 2 - '+0' canonical zero value - ' -123 123 123' canonical value '-123123123' - '1 23 456 7890' canonical value '+1234567890' +=item bfround ( +$scale ) rounds to the $scale'th place left from the '.' +=item bround ( +$scale ) preserves accuracy to $scale sighnificant digits counted from the left and paddes the number with zeros + +=item bround ( -$scale ) preserves accuracy to $scale significant digits counted from the right and paddes the number with zeros. + +=back + +C does nothing in case of negative C<$scale>. Both C and +C are a no-ops for a scale of 0. + +All rounding functions take as a second parameter a rounding mode from one of +the following: 'even', 'odd', '+inf', '-inf', 'zero' or 'trunc'. + +The default is 'even'. By using C<< Math::BigInt->round_mode($rnd_mode); >> +you can get and set the default round mode for subsequent rounding. + +The second parameter to the round functions than overrides the default +temporarily. + +=head2 Internals + +Actual math is done in an internal format consisting of an array of +elements of base 100000 digits with the least significant digit first. +The sign C is stored separately. The string 'NaN' is used to +represent the result when input arguments are not numbers, as well as +the result of dividing by zero. + +You sould neither care nor depend on the internal represantation, it might +change without notice. Use only method calls like C<< $x->sign(); >> instead +relying on the internal hash keys like in C<< $x->{sign}; >>. + +=head2 mantissa(), exponent() and parts() + +C and C return the said parts of the BigInt such +that: + + $m = $x->mantissa(); + $e = $x->exponent(); + $y = $m * ( 10 ** $e ); + print "ok\n" if $x == $y; + +C<($m,$e) = $x->parts()> is just a shortcut that gives you both of them in one +go. Both the returned mantissa and exponent do have a sign. + +Currently, for BigInts C<$e> will be always 0, except for NaN where it will be +NaN and for $x == 0, then it will be 1 (to be compatible with Math::BigFlaot's +internal representation of a zero as C<0E1>). + +C<$m> will always be a copy of the original number. The relation between $e +and $m might change in the future, but will be always equivalent in a +numerical sense, e.g. $m might get minized. + +=head1 EXAMPLES + + use Math::BigInt qw(bstr bint); + $x = bstr("1234") # string "1234" + $x = "$x"; # same as bstr() + $x = bneg("1234") # Bigint "-1234" + $x = Math::BigInt->bneg("1234"); # Bigint "-1234" + $x = Math::BigInt->babs("-12345"); # Bigint "12345" + $x = Math::BigInt->bnorm("-0 00"); # BigInt "0" + $x = bint(1) + bint(2); # BigInt "3" + $x = bint(1) + "2"; # ditto (auto-BigIntify of "2") + $x = bint(1); # BigInt "1" + $x = $x + 5 / 2; # BigInt "3" + $x = $x ** 3; # BigInt "27" + $x *= 2; # BigInt "54" + $x = new Math::BigInt; # BigInt "0" + $x--; # BigInt "-1" + $x = Math::BigInt->badd(4,5) # BigInt "9" + $x = Math::BigInt::badd(4,5) # BigInt "9" + print $x->bsstr(); # 9e+0 =head1 Autocreating constants -After C all the integer decimal constants -in the given scope are converted to C. This conversion +After C all the B decimal constants +in the given scope are converted to C. This conversion happens at compile time. In particular - perl -MMath::BigInt=:constant -e 'print 2**100' + perl -MMath::BigInt=:constant -e 'print 2**100,"\n"' + +prints the integer value of C<2**100>. Note that without conversion of +constants the expression 2**100 will be calculated as floating point +number. + +Please note that strings and floating point constants are not affected, +so that + + use Math::BigInt qw/:constant/; + + $x = 1234567890123456789012345678901234567890 + + 123456789123456789; + $x = '1234567890123456789012345678901234567890' + + '123456789123456789'; -print the integer value of C<2**100>. Note that without conversion of -constants the expression 2**100 will be calculated as floating point number. +do both not work. You need a explicit Math::BigInt->new() around one of them. + +=head1 PERFORMANCE + +Using the form $x += $y; etc over $x = $x + $y is faster, since a copy of $x +must be made in the second case. For long numbers, the copy can eat up to 20% +of the work (in case of addition/subtraction, less for +multiplication/division). If $y is very small compared to $x, the form +$x += $y is MUCH faster than $x = $x + $y since making the copy of $x takes +more time then the actual addition. + +With a technic called copy-on-write the cost of copying with overload could +be minimized or even completely avoided. This is currently not implemented. + +The new version of this module is slower on new(), bstr() and numify(). Some +operations may be slower for small numbers, but are significantly faster for +big numbers. Other operations are now constant (O(1), like bneg(), babs() +etc), instead of O(N) and thus nearly always take much less time. + +For more benchmark results see http://bloodgate.com/perl/benchmarks.html =head1 BUGS -The current version of this module is a preliminary version of the -real thing that is currently (as of perl5.002) under development. +=over 2 + +=item :constant and eval() + +Under Perl prior to 5.6.0 having an C and +C in your code will crash with "Out of memory". This is probably an +overload/exporter bug. You can workaround by not having C +and ':constant' at the same time or upgrade your Perl. + +=back + +=head1 CAVEATS + +Some things might not work as you expect them. Below is documented what is +known to be troublesome: + +=over 1 + +=item stringify, bstr(), bsstr() and 'cmp' + +Both stringify and bstr() now drop the leading '+'. The old code would return +'+3', the new returns '3'. This is to be consistent with Perl and to make +cmp (especially with overloading) to work as you expect. It also solves +problems with Test.pm, it's ok() uses 'eq' internally. + +Mark said, when asked about to drop the '+' altogether, or make only cmp work: + + I agree (with the first alternative), don't add the '+' on positive + numbers. It's not as important anymore with the new internal + form for numbers. It made doing things like abs and neg easier, + but those have to be done differently now anyway. + +So, the following examples will now work all as expected: + + use Test; + BEGIN { plan tests => 1 } + use Math::BigInt; + + my $x = new Math::BigInt 3*3; + my $y = new Math::BigInt 3*3; + + ok ($x,3*3); + print "$x eq 9" if $x eq $y; + print "$x eq 9" if $x eq '9'; + print "$x eq 9" if $x eq 3*3; + +Additionally, the following still works: + + print "$x == 9" if $x == $y; + print "$x == 9" if $x == 9; + print "$x == 9" if $x == 3*3; + +There is now a C method to get the string in scientific notation aka +C<1e+2> instead of C<100>. Be advised that overloaded 'eq' always uses bstr() +for comparisation, but Perl will represent some numbers as 100 and others +as 1e+308. If in doubt, convert both arguments to Math::BigInt before doing eq: + + use Test; + BEGIN { plan tests => 3 } + use Math::BigInt; + + $x = Math::BigInt->new('1e56'); $y = 1e56; + ok ($x,$y); # will fail + ok ($x->bsstr(),$y); # okay + $y = Math::BigInt->new($y); + ok ($x,$y); # okay + +=item int() + +C will return (at least for Perl v5.7.1 and up) another BigInt, not a +Perl scalar: + + $x = Math::BigInt->new(123); + $y = int($x); # BigInt 123 + $x = Math::BigFloat->new(123.45); + $y = int($x); # BigInt 123 + +In all Perl versions you can use C for the same effect: + + $x = Math::BigFloat->new(123.45); + $y = $x->as_number(); # BigInt 123 + +This also works for other subclasses, like Math::String. + +=item bdiv + +The following will probably not do what you expect: + + print $c->bdiv(10000),"\n"; + +It prints both quotient and reminder since print calls C in list +context. Also, C will modify $c, so be carefull. You probably want +to use + + print $c / 10000,"\n"; + print scalar $c->bdiv(10000),"\n"; # or if you want to modify $c + +instead. + +The quotient is always the greatest integer less than or equal to the +real-valued quotient of the two operands, and the remainder (when it is +nonzero) always has the same sign as the second operand; so, for +example, + + 1 / 4 => ( 0, 1) + 1 / -4 => (-1,-3) + -3 / 4 => (-1, 1) + -3 / -4 => ( 0,-3) + +As a consequence, the behavior of the operator % agrees with the +behavior of Perl's built-in % operator (as documented in the perlop +manpage), and the equation + + $x == ($x / $y) * $y + ($x % $y) + +holds true for any $x and $y, which justifies calling the two return +values of bdiv() the quotient and remainder. + +Perl's 'use integer;' changes the behaviour of % and / for scalars, but will +not change BigInt's way to do things. This is because under 'use integer' Perl +will do what the underlying C thinks is right and this is different for each +system. If you need BigInt's behaving exactly like Perl's 'use integer', bug +the author to implement it ;) + +=item Modifying and = + +Beware of: + + $x = Math::BigFloat->new(5); + $y = $x; + +It will not do what you think, e.g. making a copy of $x. Instead it just makes +a second reference to the B object and stores it in $y. Thus anything +that modifies $x will modify $y, and vice versa. + + $x->bmul(2); + print "$x, $y\n"; # prints '10, 10' + +If you want a true copy of $x, use: + + $y = $x->copy(); + +See also the documentation in L regarding C<=>. + +=item bpow + +C (and the rounding functions) now modifies the first argument and +return it, unlike the old code which left it alone and only returned the +result. This is to be consistent with C etc. The first three will +modify $x, the last one won't: + + print bpow($x,$i),"\n"; # modify $x + print $x->bpow($i),"\n"; # ditto + print $x **= $i,"\n"; # the same + print $x ** $i,"\n"; # leave $x alone + +The form C<$x **= $y> is faster than C<$x = $x ** $y;>, though. + +=item Overloading -$x + +The following: + + $x = -$x; + +is slower than + + $x->bneg(); + +since overload calls C instead of C. The first variant +needs to preserve $x since it does not know that it later will get overwritten. +This makes a copy of $x and takes O(N). But $x->bneg() is O(1). + +With Copy-On-Write, this issue will be gone. Stay tuned... + +=item Mixing different object types + +In Perl you will get a floating point value if you do one of the following: + + $float = 5.0 + 2; + $float = 2 + 5.0; + $float = 5 / 2; + +With overloaded math, only the first two variants will result in a BigFloat: + + use Math::BigInt; + use Math::BigFloat; + + $mbf = Math::BigFloat->new(5); + $mbi2 = Math::BigInteger->new(5); + $mbi = Math::BigInteger->new(2); + + # what actually gets called: + $float = $mbf + $mbi; # $mbf->badd() + $float = $mbf / $mbi; # $mbf->bdiv() + $integer = $mbi + $mbf; # $mbi->badd() + $integer = $mbi2 / $mbi; # $mbi2->bdiv() + $integer = $mbi2 / $mbf; # $mbi2->bdiv() + +This is because math with overloaded operators follows the first (dominating) +operand, this one's operation is called and returns thus the result. So, +Math::BigInt::bdiv() will always return a Math::BigInt, regardless whether +the result should be a Math::BigFloat or the second operant is one. + +To get a Math::BigFloat you either need to call the operation manually, +make sure the operands are already of the proper type or casted to that type +via Math::BigFloat->new(): + + $float = Math::BigFloat->new($mbi2) / $mbi; # = 2.5 + +Beware of simple "casting" the entire expression, this would only convert +the already computed result: + + $float = Math::BigFloat->new($mbi2 / $mbi); # = 2.0 thus wrong! + +Beware of the order of more complicated expressions like: + + $integer = ($mbi2 + $mbi) / $mbf; # int / float => int + $integer = $mbi2 / Math::BigFloat->new($mbi); # ditto + +If in doubt, break the expression into simpler terms, or cast all operands +to the desired resulting type. + +Scalar values are a bit different, since: + + $float = 2 + $mbf; + $float = $mbf + 2; + +will both result in the proper type due to the way the overloaded math works. + +This section also applies to other overloaded math packages, like Math::String. + +=item bsqrt() + +C works only good if the result is an big integer, e.g. the square +root of 144 is 12, but from 12 the square root is 3, regardless of rounding +mode. + +If you want a better approximation of the square root, then use: + + $x = Math::BigFloat->new(12); + $Math::BigFloat::precision = 0; + Math::BigFloat->round_mode('even'); + print $x->copy->bsqrt(),"\n"; # 4 + + $Math::BigFloat::precision = 2; + print $x->bsqrt(),"\n"; # 3.46 + print $x->bsqrt(3),"\n"; # 3.464 + +=back + +=head1 LICENSE + +This program is free software; you may redistribute it and/or modify it under +the same terms as Perl itself. -=head1 AUTHOR +=head1 AUTHORS -Mark Biggar, overloaded interface by Ilya Zakharevich. +Original code by Mark Biggar, overloaded interface by Ilya Zakharevich. +Completely rewritten by Tels http://bloodgate.com in late 2000, 2001. =cut diff --git a/t/lib/bigfltpm.t b/t/lib/bigfltpm.t index 8247e42..e8de58d 100755 --- a/t/lib/bigfltpm.t +++ b/t/lib/bigfltpm.t @@ -1,208 +1,294 @@ -#!./perl +#!/usr/bin/perl -w -BEGIN { - chdir 't' if -d 't'; - @INC = '../lib'; -} +use Test; +use strict; + +BEGIN + { + $| = 1; + unshift @INC, '../lib'; # for running manually + # chdir 't' if -d 't'; + plan tests => 514; + } use Math::BigFloat; +use Math::BigInt; -$test = 0; -$| = 1; -print "1..414\n"; -while () { - chomp; - if (s/^&//) { - $f = $_; - } elsif (/^\$.*/) { - eval "$_;"; - } else { - ++$test; - if (m|^(.*?):(/.+)$|) { - $ans = $2; - @args = split(/:/,$1,99); - } - else { - @args = split(/:/,$_,99); - $ans = pop(@args); - } - $try = "\$x = new Math::BigFloat \"$args[0]\";"; - if ($f eq "fnorm"){ - $try .= "\$x+0;"; - } elsif ($f eq "fneg") { - $try .= "-\$x;"; - } elsif ($f eq "fabs") { - $try .= "abs \$x;"; - } elsif ($f eq "fint") { - $try .= "int \$x;"; - } elsif ($f eq "fround") { - $try .= "0+\$x->fround($args[1]);"; - } elsif ($f eq "ffround") { - $try .= "0+\$x->ffround($args[1]);"; - } elsif ($f eq "fsqrt") { - $try .= "0+\$x->fsqrt;"; - } else { - $try .= "\$y = new Math::BigFloat \"$args[1]\";"; - if ($f eq "fcmp") { - $try .= "\$x <=> \$y;"; - } elsif ($f eq "fadd") { - $try .= "\$x + \$y;"; - } elsif ($f eq "fsub") { - $try .= "\$x - \$y;"; - } elsif ($f eq "fmul") { - $try .= "\$x * \$y;"; - } elsif ($f eq "fdiv") { - $try .= "\$x / \$y;"; - } elsif ($f eq "fmod") { - $try .= "\$x % \$y;"; - } else { warn "Unknown op"; } - } - #print ">>>",$try,"<<<\n"; - $ans1 = eval $try; - if ($ans =~ m|^/(.*)$|) { - my $pat = $1; - if ($ans1 =~ /$pat/) { - print "ok $test\n"; - } - else { - print "not ok $test\n"; - print "# '$try' expected: /$pat/ got: '$ans1'\n"; - } - } - else { +my ($x,$y,$f,@args,$ans,$try,$ans1,$ans1_str,$setup); +while () + { + chop; + $_ =~ s/#.*$//; # remove comments + $_ =~ s/\s+$//; # trailing spaces + next if /^$/; # skip empty lines & comments + if (s/^&//) + { + $f = $_; + } + elsif (/^\$/) + { + $setup = $_; $setup =~ s/^\$/\$Math::BigFloat::/; # rnd_mode, div_scale + # print "$setup\n"; + } + else + { + if (m|^(.*?):(/.+)$|) + { + $ans = $2; + @args = split(/:/,$1,99); + } + else + { + @args = split(/:/,$_,99); $ans = pop(@args); + } + $try = "\$x = new Math::BigFloat \"$args[0]\";"; + if ($f eq "fnorm") + { + $try .= "\$x;"; + } elsif ($f eq "binf") { + $try .= "\$x->binf('$args[1]');"; + } elsif ($f eq "bsstr") { + $try .= "\$x->bsstr();"; + } elsif ($f eq "_set") { + $try .= "\$x->_set('$args[1]'); \$x;"; + } elsif ($f eq "fneg") { + $try .= "-\$x;"; + } elsif ($f eq "bfloor") { + $try .= "\$x->bfloor();"; + } elsif ($f eq "bceil") { + $try .= "\$x->bceil();"; + } elsif ($f eq "is_zero") { + $try .= "\$x->is_zero()+0;"; + } elsif ($f eq "is_one") { + $try .= "\$x->is_one()+0;"; + } elsif ($f eq "is_odd") { + $try .= "\$x->is_odd()+0;"; + } elsif ($f eq "is_even") { + $try .= "\$x->is_even()+0;"; + } elsif ($f eq "as_number") { + $try .= "\$x->as_number();"; + } elsif ($f eq "fpow") { + $try .= "\$x ** $args[1];"; + } elsif ($f eq "fabs") { + $try .= "abs \$x;"; + }elsif ($f eq "fround") { + $try .= "$setup; \$x->fround($args[1]);"; + } elsif ($f eq "ffround") { + $try .= "$setup; \$x->ffround($args[1]);"; + } elsif ($f eq "fsqrt") { + $try .= "$setup; \$x->fsqrt();"; + } + else + { + $try .= "\$y = new Math::BigFloat \"$args[1]\";"; + if ($f eq "fcmp") { + $try .= "\$x <=> \$y;"; + } elsif ($f eq "fadd") { + $try .= "\$x + \$y;"; + } elsif ($f eq "fsub") { + $try .= "\$x - \$y;"; + } elsif ($f eq "fmul") { + $try .= "\$x * \$y;"; + } elsif ($f eq "fdiv") { + $try .= "$setup; \$x / \$y;"; + } elsif ($f eq "fmod") { + $try .= "\$x % \$y;"; + } else { warn "Unknown op '$f'"; } + } + $ans1 = eval $try; + if ($ans =~ m|^/(.*)$|) + { + my $pat = $1; + if ($ans1 =~ /$pat/) + { + ok (1,1); + } + else + { + print "# '$try' expected: /$pat/ got: '$ans1'\n" if !ok(1,0); + } + } + else + { + if ($ans eq "") + { + ok_undef ($ans1); + } + else + { + print "# Tried: '$try'\n" if !ok ($ans1, $ans); + } + } # end pattern or string + } + } # end while - $ans1_str = defined $ans1? "$ans1" : ""; - if ($ans1_str eq $ans) { #bug! - print "ok $test\n"; - } else { - print "not ok $test\n"; - print "# '$try' expected: '$ans' got: '$ans1'\n"; - } - } - } -} +# all done -{ - use Math::BigFloat ':constant'; +############################################################################### +# Perl 5.005 does not like ok ($x,undef) - $test++; - # print "# " . 2. * '1427247692705959881058285969449495136382746624' . "\n"; - print "not " - unless 2. * '1427247692705959881058285969449495136382746624' - == "2854495385411919762116571938898990272765493248."; - print "ok $test\n"; - $test++; - @a = (); - for ($i = 1.; $i < 10; $i++) { - push @a, $i; - } - print "not " unless "@a" eq "1. 2. 3. 4. 5. 6. 7. 8. 9."; - print "ok $test\n"; -} +sub ok_undef + { + my $x = shift; + ok (1,1) and return if !defined $x; + ok ($x,'undef'); + } + __END__ +&as_number +0:0 +1:1 +1.2:1 +2.345:2 +-2:-2 +-123.456:-123 +-200:-200 +&binf +1:+:+inf +2:-:-inf +3:abc:+inf +&bsstr ++inf:+inf +-inf:-inf +abc:NaN &fnorm ++inf:+inf +-inf:-inf ++infinity:NaN ++-inf:NaN abc:NaN 1 a:NaN 1bcd2:NaN 11111b:NaN +1z:NaN -1z:NaN -0:0. -+0:0. -+00:0. -+0 0 0:0. -000000 0000000 00000:0. --0:0. --0000:0. -+1:1. -+01:1. -+001:1. -+00000100000:100000. -123456789:123456789. --1:-1. --01:-1. --001:-1. --123456789:-123456789. --00000100000:-100000. +0:0 ++0:0 ++00:0 ++0_0_0:0 +000000_0000000_00000:0 +-0:0 +-0000:0 ++1:1 ++01:1 ++001:1 ++00000100000:100000 +123456789:123456789 +-1:-1 +-01:-1 +-001:-1 +-123456789:-123456789 +-00000100000:-100000 123.456a:NaN 123.456:123.456 -0.01:.01 -.002:.002 --0.0003:-.0003 --.0000000004:-.0000000004 -123456E2:12345600. +0.01:0.01 +.002:0.002 ++.2:0.2 +-0.0003:-0.0003 +-.0000000004:-0.0000000004 +123456E2:12345600 123456E-2:1234.56 --123456E2:-12345600. +-123456E2:-12345600 -123456E-2:-1234.56 -1e1:10. -2e-11:.00000000002 --3e111:-3000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000. --4e-1111:-.0000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000004 +1e1:10 +2e-11:0.00000000002 +-3e111:-3000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000 +-4e-1111:-0.0000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000004 +&fpow +2:2:4 +1:2:1 +1:3:1 +-1:2:1 +-1:3:-1 +123.456:2:15241.383936 +2:-2:0.25 +2:-3:0.125 +128:-2:0.00006103515625 &fneg abc:NaN -+0:0. -+1:-1. --1:1. -+123456789:-123456789. --123456789:123456789. ++0:0 ++1:-1 +-1:1 ++123456789:-123456789 +-123456789:123456789 +123.456789:-123.456789 -123456.789:123456.789 &fabs abc:NaN -+0:0. -+1:1. --1:1. -+123456789:123456789. --123456789:123456789. ++0:0 ++1:1 +-1:1 ++123456789:123456789 +-123456789:123456789 +123.456789:123.456789 -123456.789:123456.789 &fround -$Math::BigFloat::rnd_mode = 'trunc' +$rnd_mode = "trunc" +10123456789:5:10123000000 -10123456789:5:-10123000000 ++10123456789.123:5:10123000000 +-10123456789.123:5:-10123000000 +10123456789:9:10123456700 -10123456789:9:-10123456700 +101234500:6:101234000 -101234500:6:-101234000 -$Math::BigFloat::rnd_mode = 'zero' +$rnd_mode = "zero" +20123456789:5:20123000000 -20123456789:5:-20123000000 ++20123456789.123:5:20123000000 +-20123456789.123:5:-20123000000 +20123456789:9:20123456800 -20123456789:9:-20123456800 +201234500:6:201234000 -201234500:6:-201234000 -$Math::BigFloat::rnd_mode = '+inf' +$rnd_mode = "+inf" +30123456789:5:30123000000 -30123456789:5:-30123000000 ++30123456789.123:5:30123000000 +-30123456789.123:5:-30123000000 +30123456789:9:30123456800 -30123456789:9:-30123456800 +301234500:6:301235000 -301234500:6:-301234000 -$Math::BigFloat::rnd_mode = '-inf' +$rnd_mode = "-inf" +40123456789:5:40123000000 -40123456789:5:-40123000000 ++40123456789.123:5:40123000000 +-40123456789.123:5:-40123000000 +40123456789:9:40123456800 -40123456789:9:-40123456800 +401234500:6:401234000 -401234500:6:-401235000 -$Math::BigFloat::rnd_mode = 'odd' +$rnd_mode = "odd" +50123456789:5:50123000000 -50123456789:5:-50123000000 ++50123456789.123:5:50123000000 +-50123456789.123:5:-50123000000 +50123456789:9:50123456800 -50123456789:9:-50123456800 +501234500:6:501235000 -501234500:6:-501235000 -$Math::BigFloat::rnd_mode = 'even' +$rnd_mode = "even" +60123456789:5:60123000000 -60123456789:5:-60123000000 +60123456789:9:60123456800 -60123456789:9:-60123456800 +601234500:6:601234000 -601234500:6:-601234000 ++60123456789.0123:5:60123000000 +-60123456789.0123:5:-60123000000 &ffround -$Math::BigFloat::rnd_mode = 'trunc' +$rnd_mode = "trunc" +1.23:-1:1.2 ++1.234:-1:1.2 ++1.2345:-1:1.2 ++1.23:-2:1.23 ++1.234:-2:1.23 ++1.2345:-2:1.23 ++1.23:-3:1.23 ++1.234:-3:1.234 ++1.2345:-3:1.234 -1.23:-1:-1.2 +1.27:-1:1.2 -1.27:-1:-1.2 @@ -210,12 +296,22 @@ $Math::BigFloat::rnd_mode = 'trunc' -1.25:-1:-1.2 +1.35:-1:1.3 -1.35:-1:-1.3 +-0.0061234567890:-1:0 +-0.0061:-1:0 +-0.00612:-1:0 +-0.00612:-2:0 -0.006:-1:0 -0.006:-2:0 +-0.0006:-2:0 +-0.0006:-3:0 -0.0065:-3:/-0\.006|-6e-03 -0.0065:-4:/-0\.006(?:5|49{5}\d+)|-6\.5e-03 -0.0065:-5:/-0\.006(?:5|49{5}\d+)|-6\.5e-03 -$Math::BigFloat::rnd_mode = 'zero' +0.05:0:0 +0.5:0:0 +0.51:0:0 +0.41:0:0 +$rnd_mode = "zero" +2.23:-1:/2.2(?:0{5}\d+)? -2.23:-1:/-2.2(?:0{5}\d+)? +2.27:-1:/2.(?:3|29{5}\d+) @@ -229,7 +325,11 @@ $Math::BigFloat::rnd_mode = 'zero' -0.0065:-3:/-0\.006|-6e-03 -0.0065:-4:/-0\.006(?:5|49{5}\d+)|-6\.5e-03 -0.0065:-5:/-0\.006(?:5|49{5}\d+)|-6\.5e-03 -$Math::BigFloat::rnd_mode = '+inf' +0.05:0:0 +0.5:0:0 +0.51:0:1 +0.41:0:0 +$rnd_mode = "+inf" +3.23:-1:/3.2(?:0{5}\d+)? -3.23:-1:/-3.2(?:0{5}\d+)? +3.27:-1:/3.(?:3|29{5}\d+) @@ -243,7 +343,11 @@ $Math::BigFloat::rnd_mode = '+inf' -0.0065:-3:/-0\.006|-6e-03 -0.0065:-4:/-0\.006(?:5|49{5}\d+)|-6\.5e-03 -0.0065:-5:/-0\.006(?:5|49{5}\d+)|-6\.5e-03 -$Math::BigFloat::rnd_mode = '-inf' +0.05:0:0 +0.5:0:1 +0.51:0:1 +0.41:0:0 +$rnd_mode = "-inf" +4.23:-1:/4.2(?:0{5}\d+)? -4.23:-1:/-4.2(?:0{5}\d+)? +4.27:-1:/4.(?:3|29{5}\d+) @@ -257,7 +361,11 @@ $Math::BigFloat::rnd_mode = '-inf' -0.0065:-3:/-0\.007|-7e-03 -0.0065:-4:/-0\.006(?:5|49{5}\d+)|-6\.5e-03 -0.0065:-5:/-0\.006(?:5|49{5}\d+)|-6\.5e-03 -$Math::BigFloat::rnd_mode = 'odd' +0.05:0:0 +0.5:0:0 +0.51:0:1 +0.41:0:0 +$rnd_mode = "odd" +5.23:-1:/5.2(?:0{5}\d+)? -5.23:-1:/-5.2(?:0{5}\d+)? +5.27:-1:/5.(?:3|29{5}\d+) @@ -271,7 +379,11 @@ $Math::BigFloat::rnd_mode = 'odd' -0.0065:-3:/-0\.007|-7e-03 -0.0065:-4:/-0\.006(?:5|49{5}\d+)|-6\.5e-03 -0.0065:-5:/-0\.006(?:5|49{5}\d+)|-6\.5e-03 -$Math::BigFloat::rnd_mode = 'even' +0.05:0:0 +0.5:0:1 +0.51:0:1 +0.41:0:0 +$rnd_mode = "even" +6.23:-1:/6.2(?:0{5}\d+)? -6.23:-1:/-6.2(?:0{5}\d+)? +6.27:-1:/6.(?:3|29{5}\d+) @@ -282,9 +394,21 @@ $Math::BigFloat::rnd_mode = 'even' -6.35:-1:/-6.(?:4|39{5}\d+|29{8}\d+) -0.0065:-1:0 -0.0065:-2:/-0\.01|-1e-02 --0.0065:-3:/-0\.006|-7e-03|-6e-03 +-0.0065:-3:/-0\.006|-7e-03 -0.0065:-4:/-0\.006(?:5|49{5}\d+)|-6\.5e-03 -0.0065:-5:/-0\.006(?:5|49{5}\d+)|-6\.5e-03 +0.05:0:0 +0.5:0:0 +0.51:0:1 +0.41:0:0 +0.01234567:-3:0.012 +0.01234567:-4:0.0123 +0.01234567:-5:0.01235 +0.01234567:-6:0.012346 +0.01234567:-7:0.0123457 +0.01234567:-8:0.01234567 +0.01234567:-9:0.01234567 +0.01234567:-12:0.01234567 &fcmp abc:abc: abc:+0: @@ -312,140 +436,162 @@ abc:+0: +124:+123:1 -123:-124:1 -124:-123:-1 +0:0.01:-1 +0:0.0001:-1 +0:-0.0001:1 +0:-0.1:1 +0.1:0:1 +0.00001:0:1 +-0.0001:0:-1 +-0.1:0:-1 +0:0.0001234:-1 +0:-0.0001234:1 +0.0001234:0:1 +-0.0001234:0:-1 +0.0001:0.0005:-1 +0.0005:0.0001:1 +0.005:0.0001:1 +0.001:0.0005:1 +0.000001:0.0005:-2 # <0, but can't test this +0.00000123:0.0005:-2 # <0, but can't test this +0.00512:0.0001:1 +0.005:0.000112:1 +0.00123:0.0005:1 &fadd abc:abc:NaN abc:+0:NaN +0:abc:NaN -+0:+0:0. -+1:+0:1. -+0:+1:1. -+1:+1:2. --1:+0:-1. -+0:-1:-1. --1:-1:-2. --1:+1:0. -+1:-1:0. -+9:+1:10. -+99:+1:100. -+999:+1:1000. -+9999:+1:10000. -+99999:+1:100000. -+999999:+1:1000000. -+9999999:+1:10000000. -+99999999:+1:100000000. -+999999999:+1:1000000000. -+9999999999:+1:10000000000. -+99999999999:+1:100000000000. -+10:-1:9. -+100:-1:99. -+1000:-1:999. -+10000:-1:9999. -+100000:-1:99999. -+1000000:-1:999999. -+10000000:-1:9999999. -+100000000:-1:99999999. -+1000000000:-1:999999999. -+10000000000:-1:9999999999. -+123456789:+987654321:1111111110. --123456789:+987654321:864197532. --123456789:-987654321:-1111111110. -+123456789:-987654321:-864197532. ++0:+0:0 ++1:+0:1 ++0:+1:1 ++1:+1:2 +-1:+0:-1 ++0:-1:-1 +-1:-1:-2 +-1:+1:0 ++1:-1:0 ++9:+1:10 ++99:+1:100 ++999:+1:1000 ++9999:+1:10000 ++99999:+1:100000 ++999999:+1:1000000 ++9999999:+1:10000000 ++99999999:+1:100000000 ++999999999:+1:1000000000 ++9999999999:+1:10000000000 ++99999999999:+1:100000000000 ++10:-1:9 ++100:-1:99 ++1000:-1:999 ++10000:-1:9999 ++100000:-1:99999 ++1000000:-1:999999 ++10000000:-1:9999999 ++100000000:-1:99999999 ++1000000000:-1:999999999 ++10000000000:-1:9999999999 ++123456789:+987654321:1111111110 +-123456789:+987654321:864197532 +-123456789:-987654321:-1111111110 ++123456789:-987654321:-864197532 &fsub abc:abc:NaN abc:+0:NaN +0:abc:NaN -+0:+0:0. -+1:+0:1. -+0:+1:-1. -+1:+1:0. --1:+0:-1. -+0:-1:1. --1:-1:0. --1:+1:-2. -+1:-1:2. -+9:+1:8. -+99:+1:98. -+999:+1:998. -+9999:+1:9998. -+99999:+1:99998. -+999999:+1:999998. -+9999999:+1:9999998. -+99999999:+1:99999998. -+999999999:+1:999999998. -+9999999999:+1:9999999998. -+99999999999:+1:99999999998. -+10:-1:11. -+100:-1:101. -+1000:-1:1001. -+10000:-1:10001. -+100000:-1:100001. -+1000000:-1:1000001. -+10000000:-1:10000001. -+100000000:-1:100000001. -+1000000000:-1:1000000001. -+10000000000:-1:10000000001. -+123456789:+987654321:-864197532. --123456789:+987654321:-1111111110. --123456789:-987654321:864197532. -+123456789:-987654321:1111111110. ++0:+0:0 ++1:+0:1 ++0:+1:-1 ++1:+1:0 +-1:+0:-1 ++0:-1:1 +-1:-1:0 +-1:+1:-2 ++1:-1:2 ++9:+1:8 ++99:+1:98 ++999:+1:998 ++9999:+1:9998 ++99999:+1:99998 ++999999:+1:999998 ++9999999:+1:9999998 ++99999999:+1:99999998 ++999999999:+1:999999998 ++9999999999:+1:9999999998 ++99999999999:+1:99999999998 ++10:-1:11 ++100:-1:101 ++1000:-1:1001 ++10000:-1:10001 ++100000:-1:100001 ++1000000:-1:1000001 ++10000000:-1:10000001 ++100000000:-1:100000001 ++1000000000:-1:1000000001 ++10000000000:-1:10000000001 ++123456789:+987654321:-864197532 +-123456789:+987654321:-1111111110 +-123456789:-987654321:864197532 ++123456789:-987654321:1111111110 &fmul abc:abc:NaN abc:+0:NaN +0:abc:NaN -+0:+0:0. -+0:+1:0. -+1:+0:0. -+0:-1:0. --1:+0:0. -+123456789123456789:+0:0. -+0:+123456789123456789:0. --1:-1:1. --1:+1:-1. -+1:-1:-1. -+1:+1:1. -+2:+3:6. --2:+3:-6. -+2:-3:-6. --2:-3:6. -+111:+111:12321. -+10101:+10101:102030201. -+1001001:+1001001:1002003002001. -+100010001:+100010001:10002000300020001. -+10000100001:+10000100001:100002000030000200001. -+11111111111:+9:99999999999. -+22222222222:+9:199999999998. -+33333333333:+9:299999999997. -+44444444444:+9:399999999996. -+55555555555:+9:499999999995. -+66666666666:+9:599999999994. -+77777777777:+9:699999999993. -+88888888888:+9:799999999992. -+99999999999:+9:899999999991. ++0:+0:0 ++0:+1:0 ++1:+0:0 ++0:-1:0 +-1:+0:0 ++123456789123456789:+0:0 ++0:+123456789123456789:0 +-1:-1:1 +-1:+1:-1 ++1:-1:-1 ++1:+1:1 ++2:+3:6 +-2:+3:-6 ++2:-3:-6 +-2:-3:6 ++111:+111:12321 ++10101:+10101:102030201 ++1001001:+1001001:1002003002001 ++100010001:+100010001:10002000300020001 ++10000100001:+10000100001:100002000030000200001 ++11111111111:+9:99999999999 ++22222222222:+9:199999999998 ++33333333333:+9:299999999997 ++44444444444:+9:399999999996 ++55555555555:+9:499999999995 ++66666666666:+9:599999999994 ++77777777777:+9:699999999993 ++88888888888:+9:799999999992 ++99999999999:+9:899999999991 &fdiv +$div_scale = 40; $Math::BigFloat::rnd_mode = 'even' abc:abc:NaN abc:+1:abc:NaN +1:abc:NaN +0:+0:NaN -+0:+1:0. ++0:+1:0 +1:+0:NaN -+0:-1:0. ++0:-1:0 -1:+0:NaN -+1:+1:1. --1:-1:1. -+1:-1:-1. --1:+1:-1. -+1:+2:.5 -+2:+1:2. -+10:+5:2. -+100:+4:25. -+1000:+8:125. -+10000:+16:625. -+10000:-16:-625. -+999999999999:+9:111111111111. -+999999999999:+99:10101010101. -+999999999999:+999:1001001001. -+999999999999:+9999:100010001. -+999999999999999:+99999:10000100001. ++1:+1:1 +-1:-1:1 ++1:-1:-1 +-1:+1:-1 ++1:+2:0.5 ++2:+1:2 ++10:+5:2 ++100:+4:25 ++1000:+8:125 ++10000:+16:625 ++10000:-16:-625 ++999999999999:+9:111111111111 ++999999999999:+99:10101010101 ++999999999999:+999:1001001001 ++999999999999:+9999:100010001 ++999999999999999:+99999:10000100001 +1000000000:+9:111111111.1111111111111111111111111111111 +2000000000:+9:222222222.2222222222222222222222222222222 +3000000000:+9:333333333.3333333333333333333333333333333 @@ -454,12 +600,12 @@ abc:+1:abc:NaN +6000000000:+9:666666666.6666666666666666666666666666667 +7000000000:+9:777777777.7777777777777777777777777777778 +8000000000:+9:888888888.8888888888888888888888888888889 -+9000000000:+9:1000000000. ++9000000000:+9:1000000000 +35500000:+113:314159.2920353982300884955752212389380531 +71000000:+226:314159.2920353982300884955752212389380531 +106500000:+339:314159.2920353982300884955752212389380531 +1000000000:+3:333333333.3333333333333333333333333333333 -$Math::BigFloat::div_scale = 20 +$div_scale = 20 +1000000000:+9:111111111.11111111111 +2000000000:+9:222222222.22222222222 +3000000000:+9:333333333.33333333333 @@ -468,75 +614,95 @@ $Math::BigFloat::div_scale = 20 +6000000000:+9:666666666.66666666667 +7000000000:+9:777777777.77777777778 +8000000000:+9:888888888.88888888889 -+9000000000:+9:1000000000. -+35500000:+113:314159.292035398230088 -+71000000:+226:314159.292035398230088 ++9000000000:+9:1000000000 +# following two cases are the "old" behaviour, but are now (>v0.01) different +#+35500000:+113:314159.292035398230088 +#+71000000:+226:314159.292035398230088 ++35500000:+113:314159.29203539823009 ++71000000:+226:314159.29203539823009 +106500000:+339:314159.29203539823009 +1000000000:+3:333333333.33333333333 -$Math::BigFloat::div_scale = 40 -&fsqrt -+0:0 --1:/^(?i:0|\?|NaNQ?)$ --2:/^(?i:0|\?|NaNQ?)$ --16:/^(?i:0|\?|NaNQ?)$ --123.456:/^(?i:0|\?|NaNQ?)$ -+1:1. -+1.44:1.2 -+2:1.41421356237309504880168872420969807857 -+4:2. -+16:4. -+100:10. -+123.456:11.11107555549866648462149404118219234119 -+15241.383936:123.456 -&fint -+0:+0 -+1:+1 -+11111111111111111234:+11111111111111111234 --1:-1 --11111111111111111234:-11111111111111111234 -+0.3:+0 -+1.3:+1 -+23.3:+23 -+12345678901234567890:+12345678901234567890 -+12345678901234567.890:+12345678901234567 -+12345678901234567890E13:+123456789012345678900000000000000 -+12345678901234567.890E13:+123456789012345678900000000000 -+12345678901234567890E-3:+12345678901234567 -+12345678901234567.890E-3:+12345678901234 -+12345678901234567890E-13:+1234567 -+12345678901234567.890E-13:+1234 -+12345678901234567890E-17:+123 -+12345678901234567.890E-16:+1 -+12345678901234567.890E-17:+0 -+12345678901234567890E-19:+1 -+12345678901234567890E-20:+0 -+12345678901234567890E-21:+0 -+12345678901234567890E-225:+0 --0:+0 --0.3:+0 --1.3:-1 --23.3:-23 --12345678901234567890:-12345678901234567890 --12345678901234567.890:-12345678901234567 --12345678901234567890E13:-123456789012345678900000000000000 --12345678901234567.890E13:-123456789012345678900000000000 --12345678901234567890E-3:-12345678901234567 --12345678901234567.890E-3:-12345678901234 --12345678901234567890E-13:-1234567 --12345678901234567.890E-13:-1234 --12345678901234567890E-17:-123 --12345678901234567.890E-16:-1 --12345678901234567.890E-17:+0 --12345678901234567890E-19:-1 --12345678901234567890E-20:+0 --12345678901234567890E-21:+0 --12345678901234567890E-225:+0 +$div_scale = 1 +# div_scale will be 3 since $x has 3 digits ++124:+3:41.3 +# reset scale for further tests +$div_scale = 40 &fmod +0:0:NaN -+0:1:0. -+3:1:0. -+5:2:1. -+9:4:1. -+9:5:4. -+9000:56:40. -+56:9000:56. ++0:1:0 ++3:1:0 +#+5:2:1 +#+9:4:1 +#+9:5:4 +#+9000:56:40 +#+56:9000:56 +&fsqrt ++0:0 +-1:NaN +-2:NaN +-16:NaN +-123.45:NaN ++1:1 +#+1.44:1.2 +#+2:1.41421356237309504880168872420969807857 +#+4:2 +#+16:4 +#+100:10 +#+123.456:11.11107555549866648462149404118219234119 +#+15241.38393:123.456 +&is_odd +abc:0 +0:0 +-1:1 +-3:1 +1:1 +3:1 +1000001:1 +1000002:0 +2:0 +&is_even +abc:0 +0:1 +-1:0 +-3:0 +1:0 +3:0 +1000001:0 +1000002:1 +2:1 +&is_zero +NaNzero:0 +0:1 +-1:0 +1:0 +&is_one +0:0 +2:0 +1:1 +-1:0 +-2:0 +&_set +NaN:2:2 +2:abc:NaN +1:-1:-1 +2:1:1 +-2:0:0 +128:-2:-2 +&bfloor +0:0 +abc:NaN ++inf:+inf +-inf:-inf +1:1 +-51:-51 +-51.2:-52 +12.2:12 +&bceil +0:0 +abc:NaN ++inf:+inf +-inf:-inf +1:1 +-51:-51 +-51.2:-51 +12.2:13 diff --git a/t/lib/bigintpm.t b/t/lib/bigintpm.t index 6904c2d..f819104 100755 --- a/t/lib/bigintpm.t +++ b/t/lib/bigintpm.t @@ -1,114 +1,681 @@ -#!./perl +#!/usr/bin/perl -w -BEGIN { - chdir 't' if -d 't'; - @INC = '../lib'; -} +use strict; +use Test; + +BEGIN + { + $| = 1; + # chdir 't' if -d 't'; + unshift @INC, '../lib'; # for running manually + plan tests => 1190; + } + +############################################################################## +# for testing inheritance of _swap + +package Math::Foo; + +use Math::BigInt; +use vars qw/@ISA/; +@ISA = (qw/Math::BigInt/); + +use overload +# customized overload for sub, since original does not use swap there +'-' => sub { my @a = ref($_[0])->_swap(@_); + $a[0]->bsub($a[1])}; + +sub _swap + { + # a fake _swap, which reverses the params + my $self = shift; # for override in subclass + if ($_[2]) + { + my $c = ref ($_[0] ) || 'Math::Foo'; + return ( $_[0]->copy(), $_[1] ); + } + else + { + return ( Math::Foo->new($_[1]), $_[0] ); + } + } + +############################################################################## +package main; use Math::BigInt; -$test = 0; -$| = 1; -print "1..283\n"; -while () { - chop; - if (s/^&//) { - $f = $_; - } else { - ++$test; - @args = split(/:/,$_,99); - $ans = pop(@args); - $try = "\$x = new Math::BigInt \"$args[0]\";"; - if ($f eq "bnorm"){ - $try .= "\$x+0;"; - } elsif ($f eq "bneg") { - $try .= "-\$x;"; - } elsif ($f eq "babs") { - $try .= "abs \$x;"; - } elsif ($f eq "bint") { - $try .= "int \$x;"; - } else { - $try .= "\$y = new Math::BigInt \"$args[1]\";"; - if ($f eq "bcmp"){ - $try .= "\$x <=> \$y;"; - }elsif ($f eq "badd"){ - $try .= "\$x + \$y;"; - }elsif ($f eq "bsub"){ - $try .= "\$x - \$y;"; - }elsif ($f eq "bmul"){ - $try .= "\$x * \$y;"; - }elsif ($f eq "bdiv"){ - $try .= "\$x / \$y;"; - }elsif ($f eq "bmod"){ - $try .= "\$x % \$y;"; - }elsif ($f eq "bgcd"){ - $try .= "Math::BigInt::bgcd(\$x, \$y);"; - }elsif ($f eq "blsft"){ - $try .= "\$x << \$y;"; - }elsif ($f eq "brsft"){ - $try .= "\$x >> \$y;"; - }elsif ($f eq "band"){ - $try .= "\$x & \$y;"; - }elsif ($f eq "bior"){ - $try .= "\$x | \$y;"; - }elsif ($f eq "bxor"){ - $try .= "\$x ^ \$y;"; - }elsif ($f eq "bnot"){ - $try .= "~\$x;"; - } else { warn "Unknown op"; } - } - #print ">>>",$try,"<<<\n"; - $ans1 = eval $try; - if ("$ans1" eq $ans) { #bug! - print "ok $test\n"; - } else { - print "not ok $test\n"; - print "# '$try' expected: '$ans' got: '$ans1'\n"; - } - } -} - -{ - use Math::BigInt(0.02,':constant'); - - $test++; - print "not " - unless 2**150 eq "+1427247692705959881058285969449495136382746624"; - print "ok $test\n"; - $test++; - @a = (); - for ($i = 1; $i < 10; $i++) { - push @a, $i; +my (@args,$f,$try,$x,$y,$z,$a,$exp,$ans,$ans1,@a,$m,$e,$round_mode); + +while () + { + chop; + next if /^#/; # skip comments + if (s/^&//) + { + $f = $_; + } + elsif (/^\$/) + { + $round_mode = $_; + $round_mode =~ s/^\$/Math::BigInt->/; + # print "$round_mode\n"; + } + else + { + @args = split(/:/,$_,99); + $ans = pop(@args); + $try = "\$x = Math::BigInt->new(\"$args[0]\");"; + if ($f eq "bnorm"){ + # $try .= '$x+0;'; + } elsif ($f eq "_set") { + $try .= '$x->_set($args[1]); "$x";'; + } elsif ($f eq "is_zero") { + $try .= '$x->is_zero()+0;'; + } elsif ($f eq "is_one") { + $try .= '$x->is_one()+0;'; + } elsif ($f eq "is_odd") { + $try .= '$x->is_odd()+0;'; + } elsif ($f eq "is_even") { + $try .= '$x->is_even()+0;'; + } elsif ($f eq "binf") { + $try .= "\$x->binf('$args[1]');"; + } elsif ($f eq "bfloor") { + $try .= '$x->bfloor();'; + } elsif ($f eq "bceil") { + $try .= '$x->bceil();'; + } elsif ($f eq "is_inf") { + $try .= "\$x->is_inf('$args[1]')+0;"; + } elsif ($f eq "bsstr") { + $try .= '$x->bsstr();'; + } elsif ($f eq "bneg") { + $try .= '-$x;'; + } elsif ($f eq "babs") { + $try .= 'abs $x;'; + } elsif ($f eq "binc") { + $try .= '++$x;'; + } elsif ($f eq "bdec") { + $try .= '--$x;'; + }elsif ($f eq "bnot") { + $try .= '~$x;'; + }elsif ($f eq "bsqrt") { + $try .= '$x->bsqrt();'; + }elsif ($f eq "length") { + $try .= "\$x->length();"; + }elsif ($f eq "bround") { + $try .= "$round_mode; \$x->bround($args[1]);"; + }elsif ($f eq "exponent"){ + $try .= '$x = $x->exponent()->bstr();'; + }elsif ($f eq "mantissa"){ + $try .= '$x = $x->mantissa()->bstr();'; + }elsif ($f eq "parts"){ + $try .= "(\$m,\$e) = \$x->parts();"; + $try .= '$m = $m->bstr(); $m = "NaN" if !defined $m;'; + $try .= '$e = $e->bstr(); $e = "NaN" if !defined $e;'; + $try .= '"$m,$e";'; + } else { + $try .= "\$y = new Math::BigInt \"$args[1]\";"; + if ($f eq "bcmp"){ + $try .= '$x <=> $y;'; + }elsif ($f eq "bacmp"){ + $try .= '$x->bacmp($y);'; + }elsif ($f eq "badd"){ + $try .= "\$x + \$y;"; + }elsif ($f eq "bsub"){ + $try .= "\$x - \$y;"; + }elsif ($f eq "bmul"){ + $try .= "\$x * \$y;"; + }elsif ($f eq "bdiv"){ + $try .= "\$x / \$y;"; + }elsif ($f eq "bmod"){ + $try .= "\$x % \$y;"; + }elsif ($f eq "bgcd") + { + if (defined $args[2]) + { + $try .= " \$z = new Math::BigInt \"$args[2]\"; "; + } + $try .= "Math::BigInt::bgcd(\$x, \$y"; + $try .= ", \$z" if (defined $args[2]); + $try .= " );"; + } + elsif ($f eq "blcm") + { + if (defined $args[2]) + { + $try .= " \$z = new Math::BigInt \"$args[2]\"; "; + } + $try .= "Math::BigInt::blcm(\$x, \$y"; + $try .= ", \$z" if (defined $args[2]); + $try .= " );"; + }elsif ($f eq "blsft"){ + if (defined $args[2]) + { + $try .= "\$x->blsft(\$y,$args[2]);"; + } + else + { + $try .= "\$x << \$y;"; + } + }elsif ($f eq "brsft"){ + if (defined $args[2]) + { + $try .= "\$x->brsft(\$y,$args[2]);"; + } + else + { + $try .= "\$x >> \$y;"; + } + }elsif ($f eq "band"){ + $try .= "\$x & \$y;"; + }elsif ($f eq "bior"){ + $try .= "\$x | \$y;"; + }elsif ($f eq "bxor"){ + $try .= "\$x ^ \$y;"; + }elsif ($f eq "bpow"){ + $try .= "\$x ** \$y;"; + }elsif ($f eq "digit"){ + $try = "\$x = Math::BigInt->new(\"$args[0]\"); \$x->digit($args[1]);"; + } else { warn "Unknown op '$f'"; } + } + # print "trying $try\n"; + $ans1 = eval $try; + $ans =~ s/^[+]([0-9])/$1/; # remove leading '+' + if ($ans eq "") + { + ok_undef ($ans1); + } + else + { + #print "try: $try ans: $ans1 $ans\n"; + print "# Tried: '$try'\n" if !ok ($ans1, $ans); + } + # check internal state of number objects + is_valid($ans1) if ref $ans1; + } + } # endwhile data tests +close DATA; + +# test whether constant works or not +$try = "use Math::BigInt (1.31,'babs',':constant');"; +$try .= ' $x = 2**150; babs($x); $x = "$x";'; +$ans1 = eval $try; + +ok ( $ans1, "1427247692705959881058285969449495136382746624"); + +# test some more +@a = (); +for (my $i = 1; $i < 10; $i++) + { + push @a, $i; + } +ok "@a", "1 2 3 4 5 6 7 8 9"; + +# test whether selfmultiplication works correctly (result is 2**64) +$try = '$x = new Math::BigInt "+4294967296";'; +$try .= '$a = $x->bmul($x);'; +$ans1 = eval $try; +print "# Tried: '$try'\n" if !ok ($ans1, Math::BigInt->new(2) ** 64); + +# test whether op detroys args or not (should better not) + +$x = new Math::BigInt (3); +$y = new Math::BigInt (4); +$z = $x & $y; +ok ($x,3); +ok ($y,4); +ok ($z,0); +$z = $x | $y; +ok ($x,3); +ok ($y,4); +ok ($z,7); +$x = new Math::BigInt (1); +$y = new Math::BigInt (2); +$z = $x | $y; +ok ($x,1); +ok ($y,2); +ok ($z,3); + +$x = new Math::BigInt (5); +$y = new Math::BigInt (4); +$z = $x ^ $y; +ok ($x,5); +ok ($y,4); +ok ($z,1); + +$x = new Math::BigInt (-5); $y = -$x; +ok ($x, -5); + +$x = new Math::BigInt (-5); $y = abs($x); +ok ($x, -5); + +# check whether overloading cmp works +$try = "\$x = Math::BigInt->new(0);"; +$try .= "\$y = 10;"; +$try .= "'false' if \$x ne \$y;"; +$ans = eval $try; +print "# For '$try'\n" if (!ok "$ans" , "false" ); + +# we cant test for working cmpt with other objects here, we would need a dummy +# object with stringify overload for this. see Math::String tests + +############################################################################### +# check shortcuts +$try = "\$x = Math::BigInt->new(1); \$x += 9;"; +$try .= "'ok' if \$x == 10;"; +$ans = eval $try; +print "# For '$try'\n" if (!ok "$ans" , "ok" ); + +$try = "\$x = Math::BigInt->new(1); \$x -= 9;"; +$try .= "'ok' if \$x == -8;"; +$ans = eval $try; +print "# For '$try'\n" if (!ok "$ans" , "ok" ); + +$try = "\$x = Math::BigInt->new(1); \$x *= 9;"; +$try .= "'ok' if \$x == 9;"; +$ans = eval $try; +print "# For '$try'\n" if (!ok "$ans" , "ok" ); + +$try = "\$x = Math::BigInt->new(10); \$x /= 2;"; +$try .= "'ok' if \$x == 5;"; +$ans = eval $try; +print "# For '$try'\n" if (!ok "$ans" , "ok" ); + +############################################################################### +# check reversed order of arguments +$try = "\$x = Math::BigInt->new(10); \$x = 2 ** \$x;"; +$try .= "'ok' if \$x == 1024;"; $ans = eval $try; +print "# For '$try'\n" if (!ok "$ans" , "ok" ); + +$try = "\$x = Math::BigInt->new(10); \$x = 2 * \$x;"; +$try .= "'ok' if \$x == 20;"; $ans = eval $try; +print "# For '$try'\n" if (!ok "$ans" , "ok" ); + +$try = "\$x = Math::BigInt->new(10); \$x = 2 + \$x;"; +$try .= "'ok' if \$x == 12;"; $ans = eval $try; +print "# For '$try'\n" if (!ok "$ans" , "ok" ); + +$try = "\$x = Math::BigInt->new(10); \$x = 2 - \$x;"; +$try .= "'ok' if \$x == -8;"; $ans = eval $try; +print "# For '$try'\n" if (!ok "$ans" , "ok" ); + +$try = "\$x = Math::BigInt->new(10); \$x = 20 / \$x;"; +$try .= "'ok' if \$x == 2;"; $ans = eval $try; +print "# For '$try'\n" if (!ok "$ans" , "ok" ); + +############################################################################### +# check badd(4,5) form + +$try = "\$x = Math::BigInt::badd(4,5);"; +$try .= "'ok' if \$x == 9;"; +$ans = eval $try; +print "# For '$try'\n" if (!ok "$ans" , "ok" ); + +$try = "\$x = Math::BigInt->badd(4,5);"; +$try .= "'ok' if \$x == 9;"; +$ans = eval $try; +print "# For '$try'\n" if (!ok "$ans" , "ok" ); + +############################################################################### +# check proper length of internal arrays + +$x = Math::BigInt->new(99999); +ok ($x,99999); +ok (scalar @{$x->{value}}, 1); +$x += 1; +ok ($x,100000); +ok (scalar @{$x->{value}}, 2); +$x -= 1; +ok ($x,99999); +ok (scalar @{$x->{value}}, 1); + +############################################################################### +# check numify + +my $BASE = int(1e5); +$x = Math::BigInt->new($BASE-1); ok ($x->numify(),$BASE-1); +$x = Math::BigInt->new(-($BASE-1)); ok ($x->numify(),-($BASE-1)); +$x = Math::BigInt->new($BASE); ok ($x->numify(),$BASE); +$x = Math::BigInt->new(-$BASE); ok ($x->numify(),-$BASE); +$x = Math::BigInt->new( -($BASE*$BASE*1+$BASE*1+1) ); +ok($x->numify(),-($BASE*$BASE*1+$BASE*1+1)); + +############################################################################### +# test bug in _digits with length($c[-1]) where $c[-1] was "00001" instead of 1 + +$x = Math::BigInt->new(99998); $x++; $x++; $x++; $x++; +if ($x > 100000) { ok (1,1) } else { ok ("$x < 100000","$x > 100000"); } + +$x = Math::BigInt->new(100003); $x++; +$y = Math::BigInt->new(1000000); +if ($x < 1000000) { ok (1,1) } else { ok ("$x > 1000000","$x < 1000000"); } + +############################################################################### +# bug in sub where number with at least 6 trailing zeros after any op failed + +$x = Math::BigInt->new(123456); $z = Math::BigInt->new(10000); $z *= 10; +$x -= $z; +ok ($z, 100000); +ok ($x, 23456); + +############################################################################### +# bug with rest "-0" in div, causing further div()s to fail + +$x = Math::BigInt->new(-322056000); ($x,$y) = $x->bdiv('-12882240'); + +ok ($y,'0'); # not '-0' +is_valid($y); + +############################################################################### +# check undefs: NOT DONE YET + +############################################################################### +# bool + +$x = Math::BigInt->new(1); if ($x) { ok (1,1); } else { ok($x,'to be true') } +$x = Math::BigInt->new(0); if (!$x) { ok (1,1); } else { ok($x,'to be false') } + +############################################################################### +# objectify() + +@args = Math::BigInt::objectify(2,4,5); +ok (scalar @args,3); # 'Math::BigInt', 4, 5 +ok ($args[0],'Math::BigInt'); +ok ($args[1],4); +ok ($args[2],5); + +@args = Math::BigInt::objectify(0,4,5); +ok (scalar @args,3); # 'Math::BigInt', 4, 5 +ok ($args[0],'Math::BigInt'); +ok ($args[1],4); +ok ($args[2],5); + +@args = Math::BigInt::objectify(2,4,5); +ok (scalar @args,3); # 'Math::BigInt', 4, 5 +ok ($args[0],'Math::BigInt'); +ok ($args[1],4); +ok ($args[2],5); + +@args = Math::BigInt::objectify(2,4,5,6,7); +ok (scalar @args,5); # 'Math::BigInt', 4, 5, 6, 7 +ok ($args[0],'Math::BigInt'); +ok ($args[1],4); ok (ref($args[1]),$args[0]); +ok ($args[2],5); ok (ref($args[2]),$args[0]); +ok ($args[3],6); ok (ref($args[3]),''); +ok ($args[4],7); ok (ref($args[4]),''); + +@args = Math::BigInt::objectify(2,'Math::BigInt',4,5,6,7); +ok (scalar @args,5); # 'Math::BigInt', 4, 5, 6, 7 +ok ($args[0],'Math::BigInt'); +ok ($args[1],4); ok (ref($args[1]),$args[0]); +ok ($args[2],5); ok (ref($args[2]),$args[0]); +ok ($args[3],6); ok (ref($args[3]),''); +ok ($args[4],7); ok (ref($args[4]),''); + +############################################################################### +# test for flaoting-point input (other tests in bnorm() below) + +$z = 1050000000000000; # may be int on systems with 64bit? +$x = Math::BigInt->new($z); ok ($x->bsstr(),'105e+13'); # not 1.03e+15? +$z = 1e+129; # definitely a float +$x = Math::BigInt->new($z); ok ($x->bsstr(),$z); + +############################################################################### +# prime number tests, also test for **= and length() +# found on: http://www.utm.edu/research/primes/notes/by_year.html + +# ((2^148)-1)/17 +$x = Math::BigInt->new(2); $x **= 148; $x++; $x = $x / 17; +ok ($x,"20988936657440586486151264256610222593863921"); +ok ($x->length(),length "20988936657440586486151264256610222593863921"); + +# MM7 = 2^127-1 +$x = Math::BigInt->new(2); $x **= 127; $x--; +ok ($x,"170141183460469231731687303715884105727"); + +# I am afraid the following is not yet possible due to slowness +# Also, testing for 2 meg output is a bit hard ;) +#$x = new Math::BigInt(2); $x **= 6972593; $x--; + +# 593573509*2^332162+1 has exactly 100.000 digits +# takes over 16 mins and still not complete, so can not be done yet ;) +#$x = Math::BigInt->new(2); $x **= 332162; $x *= "593573509"; $x++; +#ok ($x->digits(),100000); + +############################################################################### +# inheritance and overriding of _swap + +$x = Math::Foo->new(5); +$x = $x - 8; # 8 - 5 instead of 5-8 +ok ($x,3); +ok (ref($x),'Math::Foo'); + +$x = Math::Foo->new(5); +$x = 8 - $x; # 5 - 8 instead of 8 - 5 +ok ($x,-3); +ok (ref($x),'Math::Foo'); + +############################################################################### +# all tests done + +# devel test, see whether valid catches errors +#$x = Math::BigInt->new(0); +#$x->{sign} = '-'; +#is_valid($x); # nok +# +#$x->{sign} = 'e'; +#is_valid($x); # nok +# +#$x->{value}->[0] = undef; +#is_valid($x); # nok +# +#$x->{value}->[0] = 1e6; +#is_valid($x); # nok +# +#$x->{value}->[0] = -2; +#is_valid($x); # nok +# +#$x->{sign} = '+'; +#is_valid($x); # ok + +############################################################################### +# Perl 5.005 does not like ok ($x,undef) + +sub ok_undef + { + my $x = shift; + + ok (1,1) and return if !defined $x; + ok ($x,'undef'); + } + +############################################################################### +# sub to check validity of a BigInt internally, to ensure that no op leaves a +# number object in an invalid state (f.i. "-0") + +sub is_valid + { + my $x = shift; + + my $error = ["",]; + + # ok as reference? + is_okay('ref($x)','Math::BigInt',ref($x),$error); + + # has ok sign? + is_okay('$x->{sign}',"'+', '-', '-inf', '+inf' or 'NaN'",$x->{sign},$error) + if $x->{sign} !~ /^(\+|-|\+inf|-inf|NaN)$/; + + # is not -0? + if (($x->{sign} eq '-') && (@{$x->{value}} == 1) && ($x->{value}->[0] == 0)) + { + is_okay("\$x ne '-0'","0",$x,$error); + } + # all parts are valid? + my $i = 0; my $j = scalar @{$x->{value}}; my $e; my $try; + while ($i < $j) + { + $e = $x->{value}->[$i]; $e = 'undef' unless defined $e; + $try = '=~ /^[\+]?[0-9]+\$/; '."($f, $x, $e)"; + last if $e !~ /^[+]?[0-9]+$/; + $try = ' < 0 || >= 1e5; '."($f, $x, $e)"; + last if $e <0 || $e >= 1e5; + # this test is disabled, since new/bnorm and certain ops (like early out + # in add/sub) are allowed/expected to leave '00000' in some elements + #$try = '=~ /^00+/; '."($f, $x, $e)"; + #last if $e =~ /^00+/; + $i++; + } + is_okay("\$x->{value}->[$i] $try","not $e",$e,$error) + if $i < $j; # trough all? + + # see whether errors crop up + $error->[1] = 'undef' unless defined $error->[1]; + if ($error->[0] ne "") + { + ok ($error->[1],$error->[2]); + print "# Tried: $error->[0]\n"; + } + else + { + ok (1,1); + } + } + +sub is_okay + { + my ($tried,$expected,$try,$error) = @_; + + return if $error->[0] ne ""; # error, no further testing + + @$error = ( $tried, $try, $expected ) if $try ne $expected; } - print "not " unless "@a" eq "+1 +2 +3 +4 +5 +6 +7 +8 +9"; - print "ok $test\n"; -} - + __END__ &bnorm +# binary input +0babc:NaN +0b123:NaN +0b0:0 +-0b0:0 +-0b1:-1 +0b0001:1 +0b001:1 +0b011:3 +0b101:5 +0b1000000000000000000000000000000:1073741824 +# hex input +-0x0:0 +0xabcdefgh:NaN +0x1234:4660 +0xabcdef:11259375 +-0xABCDEF:-11259375 +-0x1234:-4660 +0x12345678:305419896 +# inf input ++inf:+inf +-inf:-inf +0inf:NaN +# normal input +:NaN abc:NaN 1 a:NaN 1bcd2:NaN 11111b:NaN +1z:NaN -1z:NaN -0:+0 -+0:+0 -+00:+0 -+0 0 0:+0 -000000 0000000 00000:+0 --0:+0 --0000:+0 -+1:+1 -+01:+1 -+001:+1 -+00000100000:+100000 -123456789:+123456789 +0:0 ++0:0 ++00:0 ++000:0 +000000000000000000:0 +-0:0 +-0000:0 ++1:1 ++01:1 ++001:1 ++00000100000:100000 +123456789:123456789 -1:-1 -01:-1 -001:-1 -123456789:-123456789 -00000100000:-100000 +1_2_3:123 +_123:NaN +_123_:NaN +_123_:NaN +1__23:NaN +10000000000E-1_0:1 +1E2:100 +1E1:10 +1E0:1 +E1:NaN +E23:NaN +1.23E2:123 +1.23E1:NaN +1.23E-1:NaN +100E-1:10 +# floating point input +1.01E2:101 +1010E-1:101 +-1010E0:-1010 +-1010E1:-10100 +-1010E-2:NaN +-1.01E+1:NaN +-1.01E-1:NaN +&binf +1:+:+inf +2:-:-inf +3:abc:+inf +&is_inf ++inf::1 +-inf::1 +abc::0 +1::0 +NaN::0 +-1::0 ++inf:-:0 ++inf:+:1 +-inf:-:1 +-inf:+:0 +&blsft +abc:abc:NaN ++2:+2:+8 ++1:+32:+4294967296 ++1:+48:+281474976710656 ++8:-2:NaN +# excercise base 10 ++12345:4:10:123450000 +-1234:0:10:-1234 ++1234:0:10:+1234 ++2:2:10:200 ++12:2:10:1200 ++1234:-3:10:NaN +1234567890123:12:10:1234567890123000000000000 +&brsft +abc:abc:NaN ++8:+2:+2 ++4294967296:+32:+1 ++281474976710656:+48:+1 ++2:-2:NaN +# excercise base 10 +-1234:0:10:-1234 ++1234:0:10:+1234 ++200:2:10:2 ++1234:3:10:1 ++1234:2:10:12 ++1234:-3:10:NaN +310000:4:10:31 +12300000:5:10:123 +1230000000000:10:10:123 +09876123456789067890:12:10:9876123 +1234561234567890123:13:10:123456 +&bsstr +1e+34:1e+34 +123.456E3:123456e+0 +100:1e+2 +abc:NaN &bneg abd:NaN +0:+0 @@ -147,6 +714,29 @@ abc:+0: -123:-124:1 -124:-123:-1 +100:+5:1 +-123456789:+987654321:-1 ++123456789:-987654321:1 +-987654321:+123456789:-1 +&bacmp ++0:-0:0 ++0:+1:-1 +-1:+1:0 ++1:-1:0 +-1:+2:-1 ++2:-1:1 +-123456789:+987654321:-1 ++123456789:-987654321:-1 +-987654321:+123456789:1 +&binc +abc:NaN ++0:+1 ++1:+2 +-1:+0 +&bdec +abc:NaN ++0:-1 ++1:+0 +-1:-2 &badd abc:abc:NaN abc:+0:NaN @@ -256,6 +846,9 @@ abc:+0:NaN +77777777777:+9:+699999999993 +88888888888:+9:+799999999992 +99999999999:+9:+899999999991 ++25:+25:+625 ++12345:+12345:+152399025 ++99999:+11111:+1111088889 &bdiv abc:abc:NaN abc:+1:abc:NaN @@ -271,6 +864,7 @@ abc:+1:abc:NaN -1:+1:-1 +1:+2:+0 +2:+1:+2 ++1:+26:+0 +1000000000:+9:+111111111 +2000000000:+9:+222222222 +3000000000:+9:+333333333 @@ -293,6 +887,15 @@ abc:+1:abc:NaN +999999999999:+999:+1001001001 +999999999999:+9999:+100010001 +999999999999999:+99999:+10000100001 ++1111088889:+99999:+11111 +-5:-3:1 +4:3:1 +1:3:0 +-2:-3:0 +-2:3:-1 +1:-3:-1 +-5:3:-2 +4:-3:-2 &bmod abc:abc:NaN abc:+1:abc:NaN @@ -330,6 +933,17 @@ abc:+1:abc:NaN +999999999999:+999:+0 +999999999999:+9999:+0 +999999999999999:+99999:+0 +-9:+5:+1 ++9:-5:-1 +-9:-5:-4 +-5:3:1 +-2:3:1 +4:3:1 +1:3:1 +-5:-3:-2 +-2:-3:-2 +4:-3:-2 +1:-3:-2 &bgcd abc:abc:NaN abc:+0:NaN @@ -340,34 +954,41 @@ abc:+0:NaN +1:+1:+1 +2:+3:+1 +3:+2:+1 +-3:+2:+1 +100:+625:+25 +4096:+81:+1 -&blsft ++1034:+804:+2 ++27:+90:+56:+1 ++27:+90:+54:+9 +&blcm abc:abc:NaN -+2:+2:+8 -+1:+32:+4294967296 -+1:+48:+281474976710656 -+8:-2:NaN -&brsft -abc:abc:NaN -+8:+2:+2 -+4294967296:+32:+1 -+281474976710656:+48:+1 -+2:-2:NaN +abc:+0:NaN ++0:abc:NaN ++0:+0:NaN ++1:+0:+0 ++0:+1:+0 ++27:+90:+270 ++1034:+804:+415668 &band abc:abc:NaN +abc:0:NaN +0:abc:NaN +8:+2:+0 +281474976710656:+0:+0 +281474976710656:+1:+0 +281474976710656:+281474976710656:+281474976710656 &bior abc:abc:NaN +abc:0:NaN +0:abc:NaN +8:+2:+10 +281474976710656:+0:+281474976710656 +281474976710656:+1:+281474976710657 +281474976710656:+281474976710656:+281474976710656 &bxor abc:abc:NaN +abc:0:NaN +0:abc:NaN +8:+2:+10 +281474976710656:+0:+281474976710656 +281474976710656:+1:+281474976710657 @@ -377,9 +998,241 @@ abc:NaN +0:-1 +8:-9 +281474976710656:-281474976710657 -&bint -+0:+0 -+1:+1 -+11111111111111111234:+11111111111111111234 +&digit +0:0:0 +12:0:2 +12:1:1 +123:0:3 +123:1:2 +123:2:1 +123:-1:1 +123:-2:2 +123:-3:3 +123456:0:6 +123456:1:5 +123456:2:4 +123456:3:3 +123456:4:2 +123456:5:1 +123456:-1:1 +123456:-2:2 +123456:-3:3 +100000:-3:0 +100000:0:0 +100000:1:0 +&mantissa +abc:NaN +1e4:1 +2e0:2 +123:123 +-1:-1 +-2:-2 +&exponent +abc:NaN +1e4:4 +2e0:0 +123:0 +-1:0 +-2:0 +0:1 +&parts +abc:NaN,NaN +1e4:1,4 +2e0:2,0 +123:123,0 +-1:-1,0 +-2:-2,0 +0:0,1 +&bpow +0:0:1 +0:1:0 +0:2:0 +0:-1:NaN +0:-2:NaN +1:0:1 +1:1:1 +1:2:1 +1:3:1 +1:-1:1 +1:-2:1 +1:-3:1 +2:0:1 +2:1:2 +2:2:4 +2:3:8 +3:3:27 +2:-1:NaN +-2:-1:NaN +2:-2:NaN +-2:-2:NaN +# 1 ** -x => 1 / (1 ** x) +-1:0:1 +-2:0:1 +-1:1:-1 +-1:2:1 +-1:3:-1 +-1:4:1 +-1:5:-1 +-1:-1:-1 +-1:-2:1 +-1:-3:-1 +-1:-4:1 +10:2:100 +10:3:1000 +10:4:10000 +10:5:100000 +10:6:1000000 +10:7:10000000 +10:8:100000000 +10:9:1000000000 +10:20:100000000000000000000 +123456:2:15241383936 +&length +100:3 +10:2 +1:1 +0:1 +12345:5 +10000000000000000:17 +-123:3 +&bsqrt +144:12 +16:4 +4:2 +2:1 +12:3 +256:16 +100000000:10000 +4000000000000:2000000 +1:1 +0:0 +-2:NaN +Nan:NaN +&bround +$round_mode('trunc') +1234:0:1234 +1234:2:1200 +123456:4:123400 +123456:5:123450 +123456:6:123456 ++10123456789:5:+10123000000 +-10123456789:5:-10123000000 ++10123456789:9:+10123456700 +-10123456789:9:-10123456700 ++101234500:6:+101234000 +-101234500:6:-101234000 +#+101234500:-4:+101234000 +#-101234500:-4:-101234000 +$round_mode('zero') ++20123456789:5:+20123000000 +-20123456789:5:-20123000000 ++20123456789:9:+20123456800 +-20123456789:9:-20123456800 ++201234500:6:+201234000 +-201234500:6:-201234000 +#+201234500:-4:+201234000 +#-201234500:-4:-201234000 ++12345000:4:12340000 +-12345000:4:-12340000 +$round_mode('+inf') ++30123456789:5:+30123000000 +-30123456789:5:-30123000000 ++30123456789:9:+30123456800 +-30123456789:9:-30123456800 ++301234500:6:+301235000 +-301234500:6:-301234000 +#+301234500:-4:+301235000 +#-301234500:-4:-301234000 ++12345000:4:12350000 +-12345000:4:-12340000 +$round_mode('-inf') ++40123456789:5:+40123000000 +-40123456789:5:-40123000000 ++40123456789:9:+40123456800 +-40123456789:9:-40123456800 ++401234500:6:+401234000 ++401234500:6:+401234000 +#-401234500:-4:-401235000 +#-401234500:-4:-401235000 ++12345000:4:12340000 +-12345000:4:-12350000 +$round_mode('odd') ++50123456789:5:+50123000000 +-50123456789:5:-50123000000 ++50123456789:9:+50123456800 +-50123456789:9:-50123456800 ++501234500:6:+501235000 +-501234500:6:-501235000 +#+501234500:-4:+501235000 +#-501234500:-4:-501235000 ++12345000:4:12350000 +-12345000:4:-12350000 +$round_mode('even') ++60123456789:5:+60123000000 +-60123456789:5:-60123000000 ++60123456789:9:+60123456800 +-60123456789:9:-60123456800 ++601234500:6:+601234000 +-601234500:6:-601234000 +#+601234500:-4:+601234000 +#-601234500:-4:-601234000 +#-601234500:-9:0 +#-501234500:-9:0 +#-601234500:-8:0 +#-501234500:-8:0 ++1234567:7:1234567 ++1234567:6:1234570 ++12345000:4:12340000 +-12345000:4:-12340000 +&is_odd +abc:0 +0:0 +1:1 +3:1 +-1:1 +-3:1 +10000001:1 +10000002:0 +2:0 +&is_even +abc:0 +0:1 +1:0 +3:0 +-1:0 +-3:0 +10000001:0 +10000002:1 +2:1 +&is_zero +0:1 +NaNzero:0 +123:0 +-1:0 +1:0 +&_set +2:-1:-1 +-2:1:1 +NaN:2:2 +2:abc:NaN +&is_one +0:0 +1:1 +2:0 +-1:0 +-2:0 +# floor and ceil tests are pretty pointless in integer space...but play safe +&bfloor +0:0 -1:-1 --11111111111111111234:-11111111111111111234 +-2:-2 +2:2 +3:3 +abc:NaN +&bceil +0:0 +-1:-1 +-2:-2 +2:2 +3:3 +abc:NaN diff --git a/t/lib/mbimbf.t b/t/lib/mbimbf.t new file mode 100644 index 0000000..3948102 --- /dev/null +++ b/t/lib/mbimbf.t @@ -0,0 +1,214 @@ +#!/usr/bin/perl -w + +# test accuracy, precicion and fallback, round_mode + +use strict; +use Test; + +BEGIN + { + $| = 1; + # chdir 't' if -d 't'; + unshift @INC, '../lib'; # for running manually + plan tests => 103; + } + +use Math::BigInt; +use Math::BigFloat; + +my ($x,$y,$z,$u); + +############################################################################### +# test defaults and set/get + +ok_undef ($Math::BigInt::accuracy); +ok_undef ($Math::BigInt::precision); +ok ($Math::BigInt::div_scale,40); +ok (Math::BigInt::round_mode(),'even'); +ok ($Math::BigInt::rnd_mode,'even'); + +ok_undef ($Math::BigFloat::accuracy); +ok_undef ($Math::BigFloat::precision); +ok ($Math::BigFloat::div_scale,40); +ok ($Math::BigFloat::rnd_mode,'even'); + +# accuracy +foreach (qw/5 42 -1 0/) + { + ok ($Math::BigFloat::accuracy = $_,$_); + ok ($Math::BigInt::accuracy = $_,$_); + } +ok_undef ($Math::BigFloat::accuracy = undef); +ok_undef ($Math::BigInt::accuracy = undef); + +# precision +foreach (qw/5 42 -1 0/) + { + ok ($Math::BigFloat::precision = $_,$_); + ok ($Math::BigInt::precision = $_,$_); + } +ok_undef ($Math::BigFloat::precision = undef); +ok_undef ($Math::BigInt::precision = undef); + +# fallback +foreach (qw/5 42 1/) + { + ok ($Math::BigFloat::div_scale = $_,$_); + ok ($Math::BigInt::div_scale = $_,$_); + } +# illegal values are possible for fallback due to no accessor + +# round_mode +foreach (qw/odd even zero trunc +inf -inf/) + { + ok ($Math::BigFloat::rnd_mode = $_,$_); + ok ($Math::BigInt::rnd_mode = $_,$_); + } +$Math::BigFloat::rnd_mode = 4; +ok ($Math::BigFloat::rnd_mode,4); +ok ($Math::BigInt::rnd_mode,'-inf'); # from above + +$Math::BigInt::accuracy = undef; +$Math::BigInt::precision = undef; +# local copies +$x = Math::BigFloat->new(123.456); +ok_undef ($x->accuracy()); +ok ($x->accuracy(5),5); +ok_undef ($x->accuracy(undef),undef); +ok_undef ($x->precision()); +ok ($x->precision(5),5); +ok_undef ($x->precision(undef),undef); + +# see if MBF changes MBIs values +ok ($Math::BigInt::accuracy = 42,42); +ok ($Math::BigFloat::accuracy = 64,64); +ok ($Math::BigInt::accuracy,42); # should be still 42 +ok ($Math::BigFloat::accuracy,64); # should be still 64 + +############################################################################### +# see if creating a number under set A or P will round it + +$Math::BigInt::accuracy = 4; +$Math::BigInt::precision = 3; + +ok (Math::BigInt->new(123456),123500); # with A +$Math::BigInt::accuracy = undef; +ok (Math::BigInt->new(123456),123000); # with P + +$Math::BigFloat::accuracy = 4; +$Math::BigFloat::precision = -1; +$Math::BigInt::precision = undef; + +ok (Math::BigFloat->new(123.456),123.5); # with A +$Math::BigFloat::accuracy = undef; +ok (Math::BigFloat->new(123.456),123.5); # with P from MBF, not MBI! + +$Math::BigFloat::precision = undef; + +############################################################################### +# see if setting accuracy/precision actually rounds the number + +$x = Math::BigFloat->new(123.456); $x->accuracy(4); ok ($x,123.5); +$x = Math::BigFloat->new(123.456); $x->precision(-2); ok ($x,123.46); + +$x = Math::BigInt->new(123456); $x->accuracy(4); ok ($x,123500); +$x = Math::BigInt->new(123456); $x->precision(2); ok ($x,123500); + +############################################################################### +# test actual rounding via round() + +$x = Math::BigFloat->new(123.456); +ok ($x->copy()->round(5,2),123.46); +ok ($x->copy()->round(4,2),123.5); +ok ($x->copy()->round(undef,-2),123.46); +ok ($x->copy()->round(undef,2),100); + +$x = Math::BigFloat->new(123.45000); +ok ($x->copy()->round(undef,-1,'odd'),123.5); + +# see if rounding is 'sticky' +$x = Math::BigFloat->new(123.4567); +$y = $x->copy()->bround(); # no-op since nowhere A or P defined + +ok ($y,123.4567); +$y = $x->copy()->round(5,2); +ok ($y->accuracy(),5); +ok_undef ($y->precision()); # A has precedence, so P still unset +$y = $x->copy()->round(undef,2); +ok ($y->precision(),2); +ok_undef ($y->accuracy()); # P has precedence, so A still unset + +# does copy work? +$x = Math::BigFloat->new(123.456); $x->accuracy(4); $x->precision(2); +$z = $x->copy(); ok ($z->accuracy(),4); ok ($z->precision(),2); + +############################################################################### +# test wether operations round properly afterwards +# These tests are not complete, since they do not excercise every "return" +# statement in the op's. But heh, it's better than nothing... + +$x = Math::BigFloat->new(123.456); +$y = Math::BigFloat->new(654.321); +$x->{_a} = 5; # $x->accuracy(5) would round $x straightaway +$y->{_a} = 4; # $y->accuracy(4) would round $x straightaway + +$z = $x + $y; ok ($z,777.8); +$z = $y - $x; ok ($z,530.9); +$z = $y * $x; ok ($z,80780); +$z = $x ** 2; ok ($z,15241); +$z = $x * $x; ok ($z,15241); +# not yet: $z = -$x; ok ($z,-123.46); ok ($x,123.456); +$z = $x->copy(); $z->{_a} = 2; $z = $z / 2; ok ($z,62); +$x = Math::BigFloat->new(123456); $x->{_a} = 4; +$z = $x->copy; $z++; ok ($z,123500); + +$x = Math::BigInt->new(123456); +$y = Math::BigInt->new(654321); +$x->{_a} = 5; # $x->accuracy(5) would round $x straightaway +$y->{_a} = 4; # $y->accuracy(4) would round $x straightaway + +$z = $x + $y; ok ($z,777800); +$z = $y - $x; ok ($z,530900); +$z = $y * $x; ok ($z,80780000000); +$z = $x ** 2; ok ($z,15241000000); +# not yet: $z = -$x; ok ($z,-123460); ok ($x,123456); +$z = $x->copy; $z++; ok ($z,123460); +$z = $x->copy(); $z->{_a} = 2; $z = $z / 2; ok ($z,62000); + +############################################################################### +# test mixed arguments + +$x = Math::BigFloat->new(10); +$u = Math::BigFloat->new(2.5); +$y = Math::BigInt->new(2); + +$z = $x + $y; ok ($z,12); ok (ref($z),'Math::BigFloat'); +$z = $x / $y; ok ($z,5); ok (ref($z),'Math::BigFloat'); +$z = $u * $y; ok ($z,5); ok (ref($z),'Math::BigFloat'); + +$y = Math::BigInt->new(12345); +$z = $u->copy()->bmul($y,2,0,'odd'); ok ($z,31000); +$z = $u->copy()->bmul($y,3,0,'odd'); ok ($z,30900); +$z = $u->copy()->bmul($y,undef,0,'odd'); ok ($z,30863); +$z = $u->copy()->bmul($y,undef,1,'odd'); ok ($z,30860); +$z = $u->copy()->bmul($y,undef,-1,'odd'); ok ($z,30862.5); + +# breakage: +# $z = $y->copy()->bmul($u,2,0,'odd'); ok ($z,31000); +# $z = $y * $u; ok ($z,5); ok (ref($z),'Math::BigInt'); +# $z = $y + $x; ok ($z,12); ok (ref($z),'Math::BigInt'); +# $z = $y / $x; ok ($z,0); ok (ref($z),'Math::BigInt'); + +# all done + +############################################################################### +# Perl 5.005 does not like ok ($x,undef) + +sub ok_undef + { + my $x = shift; + + ok (1,1) and return if !defined $x; + ok ($x,'undef'); + } +