-#!/usr/bin/perl -w
+package Math::BigFloat;
+
+#
+# Mike grinned. 'Two down, infinity to go' - Mike Nostrus in Before and After
+#
# The following hash values are internally used:
# _e: exponent (BigInt)
# _p: precision
# _f: flags, used to signal MBI not to touch our private parts
-package Math::BigFloat;
-
-$VERSION = '1.27';
+$VERSION = '1.30';
require 5.005;
use Exporter;
use Math::BigInt qw/objectify/;
use strict;
use vars qw/$AUTOLOAD $accuracy $precision $div_scale $round_mode $rnd_mode/;
+use vars qw/$upgrade $downgrade/;
my $class = "Math::BigFloat";
use overload
$precision = undef;
$div_scale = 40;
+$upgrade = undef;
+$downgrade = undef;
+
##############################################################################
# the old code had $rnd_mode, so we need to support it, too
# valid method aliases for AUTOLOAD
my %methods = map { $_ => 1 }
qw / fadd fsub fmul fdiv fround ffround fsqrt fmod fstr fsstr fpow fnorm
- fint facmp fcmp fzero fnan finf finc fdec flog
+ fint facmp fcmp fzero fnan finf finc fdec flog ffac
fceil ffloor frsft flsft fone flog
/;
# valid method's that can be hand-ed up (for AUTOLOAD)
my %hand_ups = map { $_ => 1 }
qw / is_nan is_inf is_negative is_positive
accuracy precision div_scale round_mode fneg fabs babs fnot
+ objectify upgrade downgrade
+ bone binf bnan bzero
/;
sub method_alias { return exists $methods{$_[0]||''}; }
# sign => sign (+/-), or "NaN"
my ($class,$wanted,@r) = @_;
-
+
# avoid numify-calls by not using || on $wanted!
return $class->bzero() if !defined $wanted; # default to 0
return $wanted->copy() if UNIVERSAL::isa($wanted,'Math::BigFloat');
# handle '+inf', '-inf' first
if ($wanted =~ /^[+-]?inf$/)
{
+ return $downgrade->new($wanted) if $downgrade;
+
$self->{_e} = Math::BigInt->bzero();
$self->{_m} = Math::BigInt->bzero();
$self->{sign} = $wanted;
if (!ref $mis)
{
die "$wanted is not a number initialized to $class" if !$NaNOK;
+
+ return $downgrade->bnan() if $downgrade;
+
$self->{_e} = Math::BigInt->bzero();
$self->{_m} = Math::BigInt->bzero();
$self->{sign} = $nan;
$self->{_e} -= CORE::length($$mfv) if CORE::length($$mfv) != 0;
$self->{sign} = $$mis;
}
- # print "mbf new ",join(' ',@r),"\n";
+ # if downgrade, inf, NaN or integers go down
+
+ if ($downgrade && $self->{_e}->{sign} eq '+')
+ {
+# print "downgrading $$miv$$mfv"."E$$es$$ev";
+ if ($self->{_e}->is_zero())
+ {
+ $self->{_m}->{sign} = $$mis; # negative if wanted
+ return $downgrade->new($self->{_m});
+ }
+ return $downgrade->new("$$mis$$miv$$mfv"."E$$es$$ev");
+ }
+
+ # print "mbf new $self->{sign} $self->{_m} e $self->{_e}\n";
$self->bnorm()->round(@r); # first normalize, then round
}
-sub bnan
+sub _bnan
{
- # create a bigfloat 'NaN', if given a BigFloat, set it to 'NaN'
+ # used by parent class bone() to initialize number to 1
my $self = shift;
- $self = $class if !defined $self;
- if (!ref($self))
- {
- my $c = $self; $self = {}; bless $self, $c;
- }
$self->{_m} = Math::BigInt->bzero();
$self->{_e} = Math::BigInt->bzero();
- $self->{sign} = $nan;
- $self->{_a} = undef; $self->{_p} = undef;
- $self;
}
-sub binf
+sub _binf
{
- # create a bigfloat '+-inf', if given a BigFloat, set it to '+-inf'
+ # used by parent class bone() to initialize number to 1
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;
- }
$self->{_m} = Math::BigInt->bzero();
$self->{_e} = Math::BigInt->bzero();
- $self->{sign} = $sign.'inf';
- $self->{_a} = undef; $self->{_p} = undef;
- $self;
}
-sub bone
+sub _bone
{
- # create a bigfloat '+-1', if given a BigFloat, set it to '+-1'
+ # used by parent class bone() to initialize number to 1
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;
- }
$self->{_m} = Math::BigInt->bone();
$self->{_e} = Math::BigInt->bzero();
- $self->{sign} = $sign;
- if (@_ > 0)
- {
- $self->{_a} = $_[0]
- if (defined $self->{_a} && defined $_[0] && $_[0] > $self->{_a});
- $self->{_p} = $_[1]
- if (defined $self->{_p} && defined $_[1] && $_[1] < $self->{_p});
- }
- return $self;
}
-sub bzero
+sub _bzero
{
- # create a bigfloat '+0', if given a BigFloat, set it to 0
+ # used by parent class bone() to initialize number to 1
my $self = shift;
- $self = $class if !defined $self;
- if (!ref($self))
- {
- my $c = $self; $self = {}; bless $self, $c;
- }
$self->{_m} = Math::BigInt->bzero();
$self->{_e} = Math::BigInt->bone();
- $self->{sign} = '+';
- if (@_ > 0)
- {
- $self->{_a} = $_[0]
- if (defined $self->{_a} && defined $_[0] && $_[0] > $self->{_a});
- $self->{_p} = $_[1]
- if (defined $self->{_p} && defined $_[1] && $_[1] < $self->{_p});
- }
- return $self;
}
##############################################################################
my $es = '0'; my $len = 1; my $cad = 0; my $dot = '.';
- my $not_zero = !$x->is_zero();
+ my $not_zero = ! $x->is_zero();
if ($not_zero)
{
$es = $x->{_m}->bstr();
return +1 if $x->{sign} eq '+inf';
return -1 if $x->{sign} eq '-inf';
return -1 if $y->{sign} eq '+inf';
- return +1 if $y->{sign} eq '-inf';
+ return +1;
}
# check sign for speed first
# adjust so that exponents are equal
my $lxm = $x->{_m}->length();
my $lym = $y->{_m}->length();
- my $lx = $lxm + $x->{_e};
- my $ly = $lym + $y->{_e};
- # print "x $x y $y lx $lx ly $ly\n";
+ # the numify somewhat limits our length, but makes it much faster
+ my $lx = $lxm + $x->{_e}->numify();
+ my $ly = $lym + $y->{_e}->numify();
my $l = $lx - $ly; $l = -$l if $x->{sign} eq '-';
- # print "$l $x->{sign}\n";
return $l <=> 0 if $l != 0;
# lengths (corrected by exponent) are equal
- # so make mantissa euqal length by padding with zero (shift left)
+ # so make mantissa equal length by padding with zero (shift left)
my $diff = $lxm - $lym;
my $xm = $x->{_m}; # not yet copy it
my $ym = $y->{_m};
{
$xm = $x->{_m}->copy()->blsft(-$diff,10);
}
- my $rc = $xm->bcmp($ym);
+ my $rc = $xm->bacmp($ym);
$rc = -$rc if $x->{sign} eq '-'; # -124 < -123
- return $rc <=> 0;
+ $rc <=> 0;
}
sub bacmp
return undef if (($x->{sign} eq $nan) || ($y->{sign} eq $nan));
return 0 if ($x->is_inf() && $y->is_inf());
return 1 if ($x->is_inf() && !$y->is_inf());
- return -1 if (!$x->is_inf() && $y->is_inf());
+ return -1;
}
# shortcut
# adjust so that exponents are equal
my $lxm = $x->{_m}->length();
my $lym = $y->{_m}->length();
- my $lx = $lxm + $x->{_e};
- my $ly = $lym + $y->{_e};
- # print "x $x y $y lx $lx ly $ly\n";
+ # the numify somewhat limits our length, but makes it much faster
+ my $lx = $lxm + $x->{_e}->numify();
+ my $ly = $lym + $y->{_e}->numify();
my $l = $lx - $ly;
- # print "$l $x->{sign}\n";
return $l <=> 0 if $l != 0;
# lengths (corrected by exponent) are equal
{
$xm = $x->{_m}->copy()->blsft(-$diff,10);
}
- my $rc = $xm->bcmp($ym);
- return $rc <=> 0;
+ $xm->bacmp($ym) <=> 0;
}
sub badd
# return result as BFLOAT
my ($self,$x,$y,$a,$p,$r) = objectify(2,@_);
+ #print "mbf badd $x $y\n";
# inf and NaN handling
if (($x->{sign} !~ /^[+-]$/) || ($y->{sign} !~ /^[+-]$/))
{
# NaN first
return $x->bnan() if (($x->{sign} eq $nan) || ($y->{sign} eq $nan));
- # inf handline
+ # inf handling
if (($x->{sign} =~ /^[+-]inf$/) && ($y->{sign} =~ /^[+-]inf$/))
{
- # + and + => +, - and - => -, + and - => 0, - and + => 0
- return $x->bzero() if $x->{sign} ne $y->{sign};
- return $x;
+ # +inf++inf or -inf+-inf => same, rest is NaN
+ return $x if $x->{sign} eq $y->{sign};
+ return $x->bnan();
}
# +-inf + something => +inf
# something +-inf => +-inf
}
# speed: no add for 0+y or x+0
- return $x if $y->is_zero(); # x+0
+ return $x->bround($a,$p,$r) if $y->is_zero(); # x+0
if ($x->is_zero()) # 0+y
{
# make copy, clobbering up x (modify in place!)
}
# 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 $e = $y->{_e};
+ $e = Math::BigInt::bzero() if !defined $e; # if no BFLOAT ?
+ $e = $e->copy(); # make copy (didn't do it yet)
+ $e->bsub($x->{_e});
my $add = $y->{_m}->copy();
- if ($e < 0)
+# if ($e < 0) # < 0
+ if ($e->{sign} eq '-') # < 0
{
my $e1 = $e->copy()->babs();
- $x->{_m} *= (10 ** $e1);
+ #$x->{_m} *= (10 ** $e1);
+ $x->{_m}->blsft($e1,10);
$x->{_e} += $e; # need the sign of e
}
- elsif ($e > 0)
+# if ($e > 0) # > 0
+ elsif (!$e->is_zero()) # > 0
{
- $add *= (10 ** $e);
+ #$add *= (10 ** $e);
+ $add->blsft($e,10);
}
# else: both e are the same, so just leave them
$x->{_m}->{sign} = $x->{sign}; # fiddle with signs
# subtract second arg from first, modify first
my ($self,$x,$y,$a,$p,$r) = objectify(2,@_);
- if (!$y->is_zero()) # don't need to do anything if $y is 0
+ if ($y->is_zero()) # still round for not adding zero
{
- $y->{sign} =~ tr/+\-/-+/; # does nothing for NaN
- $x->badd($y,$a,$p,$r); # badd does not leave internal zeros
- $y->{sign} =~ tr/+\-/-+/; # refix $y (does nothing for NaN)
+ return $x->round($a,$p,$r);
}
+
+ $y->{sign} =~ tr/+\-/-+/; # does nothing for NaN
+ $x->badd($y,$a,$p,$r); # badd does not leave internal zeros
+ $y->{sign} =~ tr/+\-/-+/; # refix $y (does nothing for NaN)
$x; # already rounded by badd()
}
sub blog
{
- my ($self,$x,$base,@r) = ref($_[0]) ? (ref($_[0]),@_) : objectify(2,@_);
+ my ($self,$x,$base,$a,$p,$r) = ref($_[0]) ? (ref($_[0]),@_) : objectify(2,@_);
# http://www.efunda.com/math/taylor_series/logarithmic.cfm?search_string=log
# _ _
# taylor: | u 1 u^3 1 u^5 |
# ln (x) = 2 | --- + - * --- + - * --- + ... | x > 0
- # |_ v 3 v 5 v _|
+ # |_ v 3 v^3 5 v^5 _|
- return $x->bzero(@r) if $x->is_one();
- return $x->bone(@r) if $x->bcmp($base) == 0;
+ # we need to limit the accuracy to protect against overflow
+ my $fallback = 0;
+ my $scale = 0;
+ my @params = $x->_find_round_parameters($a,$p,$r);
- my $d = $r[0] || $self->accuracy() || 40;
- $d += 2; # 2 more for rounding
-
- my $u = $x->copy(); $u->bdec();
- my $v = $x->copy(); $v->binc();
+ # no rounding at all, so must use fallback
+ if (scalar @params == 1)
+ {
+ # simulate old behaviour
+ $params[1] = $self->div_scale(); # and round to it as accuracy
+ $scale = $params[1]+4; # at least four more for proper round
+ $params[3] = $r; # round mode by caller or undef
+ $fallback = 1; # to clear a/p afterwards
+ }
+ else
+ {
+ # the 4 below is empirical, and there might be cases where it is not
+ # enough...
+ $scale = abs($params[1] || $params[2]) + 4; # take whatever is defined
+ }
+
+ return $x->bzero(@params) if $x->is_one();
+ return $x->bnan() if $x->{sign} ne '+' || $x->is_zero();
+ #return $x->bone('+',@params) if $x->bcmp($base) == 0;
+
+ # when user set globals, they would interfere with our calculation, so
+ # disable then and later re-enable them
+ no strict 'refs';
+ my $abr = "$self\::accuracy"; my $ab = $$abr; $$abr = undef;
+ my $pbr = "$self\::precision"; my $pb = $$pbr; $$pbr = undef;
+ # we also need to disable any set A or P on $x (_find_round_parameters took
+ # them already into account), since these would interfere, too
+ delete $x->{_a}; delete $x->{_p};
+ # need to disable $upgrade in BigInt, to aoid deep recursion
+ local $Math::BigInt::upgrade = undef;
+
+ my $v = $x->copy(); $v->binc(); # v = x+1
+ $x->bdec(); my $u = $x->copy(); # u = x-1; x = x-1
- $x->bdec()->bdiv($v,$d); # first term: u/v
+ $x->bdiv($v,$scale); # first term: u/v
- $u *= $u; $v *= $v;
- my $below = $v->copy()->bmul($v);
- my $over = $u->copy()->bmul($u);
+ my $below = $v->copy();
+ my $over = $u->copy();
+ $u *= $u; $v *= $v; # u^2, v^2
+ $below->bmul($v); # u^3, v^3
+ $over->bmul($u);
my $factor = $self->new(3); my $two = $self->new(2);
my $diff = $self->bone();
- my $limit = $self->new("1E-". ($d-1)); my $last;
+ my $limit = $self->new("1E-". ($scale-1)); my $last;
# print "diff $diff limit $limit\n";
- while ($diff > $limit)
+ while ($diff->bcmp($limit) > 0)
{
- print "$x $over $below $factor\n";
+ #print "$x $over $below $factor\n";
$diff = $x->copy()->bsub($last)->babs();
- print "diff $diff $limit\n";
+ #print "diff $diff $limit\n";
$last = $x->copy();
- $x += $over->copy()->bdiv($below->copy()->bmul($factor),$d);
+ $x += $over->copy()->bdiv($below->copy()->bmul($factor),$scale);
$over *= $u; $below *= $v; $factor->badd($two);
}
$x->bmul($two);
- return $x->round(@r);
+
+ # shortcut to not run trough _find_round_parameters again
+ if (defined $params[1])
+ {
+ $x->bround($params[1],$params[3]); # then round accordingly
+ }
+ else
+ {
+ $x->bfround($params[2],$params[3]); # then round accordingly
+ }
+ if ($fallback)
+ {
+ # clear a/p after round, since user did not request it
+ $x->{_a} = undef; $x->{_p} = undef;
+ }
+ # restore globals
+ $$abr = $ab; $$pbr = $pb;
+
+ $x;
}
sub blcm
$x;
}
+###############################################################################
+# is_foo methods (is_negative, is_positive are inherited from BigInt)
+
+sub is_int
+ {
+ # return true if arg (BFLOAT or num_str) is an integer
+ my ($self,$x) = ref($_[0]) ? (ref($_[0]),$_[0]) : objectify(1,@_);
+
+ return 1 if ($x->{sign} =~ /^[+-]$/) && # NaN and +-inf aren't
+ $x->{_e}->{sign} eq '+'; # 1e-1 => no integer
+ 0;
+ }
+
sub is_zero
{
- # return true if arg (BFLOAT or num_str) is zero (array '+', '0')
+ # return true if arg (BFLOAT or num_str) is zero
my ($self,$x) = ref($_[0]) ? (ref($_[0]),$_[0]) : objectify(1,@_);
return 1 if $x->{sign} eq '+' && $x->{_m}->is_zero();
- return 0;
+ 0;
}
sub is_one
{
- # return true if arg (BFLOAT or num_str) is +1 (array '+', '1')
- # or -1 if signis given
+ # return true if arg (BFLOAT or num_str) is +1 or -1 if signis given
my ($self,$x) = ref($_[0]) ? (ref($_[0]),$_[0]) : objectify(1,@_);
my $sign = shift || ''; $sign = '+' if $sign ne '-';
return 1
if ($x->{sign} eq $sign && $x->{_e}->is_zero() && $x->{_m}->is_one());
- return 0;
+ 0;
}
sub is_odd
# return true if arg (BFLOAT or num_str) is odd or false if even
my ($self,$x) = ref($_[0]) ? (ref($_[0]),$_[0]) : objectify(1,@_);
- return 0 if $x->{sign} !~ /^[+-]$/; # NaN & +-inf aren't
- return 1 if ($x->{_e}->is_zero() && $x->{_m}->is_odd());
- return 0;
+ return 1 if ($x->{sign} =~ /^[+-]$/) && # NaN & +-inf aren't
+ ($x->{_e}->is_zero() && $x->{_m}->is_odd());
+ 0;
}
sub is_even
my ($self,$x) = ref($_[0]) ? (ref($_[0]),$_[0]) : objectify(1,@_);
return 0 if $x->{sign} !~ /^[+-]$/; # NaN & +-inf aren't
- return 1 if $x->{_m}->is_zero(); # 0e1 is even
- return 1 if ($x->{_e}->is_zero() && $x->{_m}->is_even()); # 123.45 is never
- return 0;
+ return 1 if ($x->{_e}->{sign} eq '+' # 123.45 is never
+ && $x->{_m}->is_even()); # but 1200 is
+ 0;
}
sub bmul
# (BINT or num_str, BINT or num_str) return BINT
my ($self,$x,$y,$a,$p,$r) = objectify(2,@_);
- # print "mbf bmul $x->{_m}e$x->{_e} $y->{_m}e$y->{_e}\n";
return $x->bnan() if (($x->{sign} eq $nan) || ($y->{sign} eq $nan));
- # handle result = 0
- return $x->bzero() if $x->is_zero() || $y->is_zero();
# inf handling
if (($x->{sign} =~ /^[+-]inf$/) || ($y->{sign} =~ /^[+-]inf$/))
{
+ return $x->bnan() if $x->is_zero() || $y->is_zero();
# result will always be +-inf:
# +inf * +/+inf => +inf, -inf * -/-inf => +inf
# +inf * -/-inf => -inf, -inf * +/+inf => -inf
return $x->binf() if ($x->{sign} =~ /^-/ && $y->{sign} =~ /^-/);
return $x->binf('-');
}
+ # handle result = 0
+ return $x->bzero() if $x->is_zero() || $y->is_zero();
# aEb * cEd = (a*c)E(b+d)
$x->{_m}->bmul($y->{_m});
# (BFLOAT,BFLOAT) (quo,rem) or BINT (only rem)
my ($self,$x,$y,$a,$p,$r) = objectify(2,@_);
- # x / +-inf => 0, reminder x
- return wantarray ? ($x->bzero(),$x->copy()) : $x->bzero()
- if $y->{sign} =~ /^[+-]inf$/;
-
- # NaN if x == NaN or y == NaN or x==y==0
- return wantarray ? ($x->bnan(),bnan()) : $x->bnan()
- if (($x->is_nan() || $y->is_nan()) ||
- ($x->is_zero() && $y->is_zero()));
+ return $self->_div_inf($x,$y)
+ if (($x->{sign} !~ /^[+-]$/) || ($y->{sign} !~ /^[+-]$/) || $y->is_zero());
- # 5 / 0 => +inf, -6 / 0 => -inf
- return wantarray
- ? ($x->binf($x->{sign}),$self->bnan()) : $x->binf($x->{sign})
- if ($x->{sign} =~ /^[+-]$/ && $y->is_zero());
-
- # x== 0 or y == 1 or y == -1
+ # x== 0 # also: or y == 1 or y == -1
return wantarray ? ($x,$self->bzero()) : $x if $x->is_zero();
+ # upgrade
+ return $upgrade->bdiv($x,$y,$a,$p,$r) if defined $upgrade;
+
# we need to limit the accuracy to protect against overflow
my $fallback = 0;
my $scale = 0;
$scale = $ly if $ly > $scale;
my $diff = $ly - $lx;
$scale += $diff if $diff > 0; # if lx << ly, but not if ly << lx!
+
+ # make copy of $x in case of list context for later reminder calculation
+ my $rem;
+ if (wantarray && !$y->is_one())
+ {
+ $rem = $x->copy();
+ }
$x->{sign} = $x->{sign} ne $y->sign() ? '-' : '+';
if (wantarray)
{
- my $rem;
if (!$y->is_one())
{
- $rem = $x->copy();
- $rem->bmod($y,$params[1],$params[2],$params[3]);
+ $rem->bmod($y,$params[1],$params[2],$params[3]); # copy already done
}
else
{
$x->blsft($shifty,10); # 123 => 1230, $y->{_m} is already 25
}
# $ym is now mantissa of $y based on exponent 0
-
+
my $shiftx = 0; # correct _e of $x by this
if ($x->{_e}->{sign} eq '-') # has digits after dot
{
$x->{_e}->bsub($shifty) if $shifty != 0;
# now mantissas are equalized, exponent of $x is adjusted, so calc result
+# $ym->{sign} = '-' if $neg; # bmod() will make the correction for us
+
$x->{_m}->bmod($ym);
$x->{sign} = '+' if $x->{_m}->is_zero(); # fix sign for -0
# disable then and later re-enable them
no strict 'refs';
my $abr = "$self\::accuracy"; my $ab = $$abr; $$abr = undef;
- $abr = "$self\::precision"; my $pb = $$abr; $$abr = undef;
+ my $pbr = "$self\::precision"; my $pb = $$pbr; $$pbr = undef;
# we also need to disable any set A or P on $x (_find_round_parameters took
# them already into account), since these would interfere, too
delete $x->{_a}; delete $x->{_p};
+ # need to disable $upgrade in BigInt, to aoid deep recursion
+ local $Math::BigInt::upgrade = undef;
my $xas = $x->as_number();
my $gs = $xas->copy()->bsqrt(); # some guess
+
if (($x->{_e}->{sign} ne '-') # guess can't be accurate if there are
# digits after the dot
- && ($xas->bcmp($gs * $gs) == 0)) # guess hit the nail on the head?
+ && ($xas->bacmp($gs * $gs) == 0)) # guess hit the nail on the head?
{
# exact result
$x->{_m} = $gs; $x->{_e} = Math::BigInt->bzero(); $x->bnorm();
# clear a/p after round, since user did not request it
$x->{_a} = undef; $x->{_p} = undef;
}
+ ${"$self\::accuracy"} = $ab; ${"$self\::precision"} = $pb;
return $x;
}
$gs = $self->new( $gs ); # BigInt to BigFloat
my $lx = $x->{_m}->length();
$scale = $lx if $scale < $lx;
my $e = $self->new("1E-$scale"); # make test variable
- return $x->bnan() if $e->sign() eq 'NaN';
+# return $x->bnan() if $e->sign() eq 'NaN';
my $y = $x->copy();
my $two = $self->new(2);
$y = $self->new($y) unless $y->isa('Math::BigFloat');
my $rem;
-# my $steps = 0;
while ($diff >= $e)
{
-# return $x->bnan() if $gs->is_zero();
-
$rem = $y->copy()->bdiv($gs,$scale)->badd($gs)->bdiv($two,$scale);
$diff = $rem->copy()->bsub($gs)->babs();
$gs = $rem->copy();
-# $steps++;
}
-# print "steps $steps\n";
# copy over to modify $x
$x->{_m} = $rem->{_m}; $x->{_e} = $rem->{_e};
$x->{_a} = undef; $x->{_p} = undef;
}
# restore globals
- ${"$self\::accuracy"} = $ab; ${"$self\::precision"} = $pb;
+ $$abr = $ab; $$pbr = $pb;
$x;
}
+sub bfac
+ {
+ # (BFLOAT or num_str, BFLOAT or num_str) return BFLOAT
+ # compute factorial numbers
+ # modifies first argument
+ my ($self,$x,@r) = objectify(1,@_);
+
+ return $x->bnan()
+ if (($x->{sign} ne '+') || # inf, NaN, <0 etc => NaN
+ ($x->{_e}->{sign} ne '+')); # digits after dot?
+
+ return $x->bone(@r) if $x->is_zero() || $x->is_one(); # 0 or 1 => 1
+
+ # use BigInt's bfac() for faster calc
+ $x->{_m}->blsft($x->{_e},10); # un-norm m
+ $x->{_e}->bzero(); # norm $x again
+ $x->{_m}->bfac(); # factorial
+ $x->bnorm()->round(@r);
+ }
+
sub bpow
{
# (BFLOAT or num_str, BFLOAT or num_str) return BFLOAT
# if $x == -1 and odd/even y => +1/-1 because +-1 ^ (+-1) => +-1
return $y1->is_odd() ? $x : $x->babs(1);
}
- return $x if $x->is_zero() && $y->{sign} eq '+'; # 0**y => 0 (if not y <= 0)
- # 0 ** -y => 1 / (0 ** y) => / 0! (1 / 0 => +inf)
- return $x->binf() if $x->is_zero() && $y->{sign} eq '-';
+ if ($x->is_zero())
+ {
+ return $x if $y->{sign} eq '+'; # 0**y => 0 (if not y <= 0)
+ # 0 ** -y => 1 / (0 ** y) => / 0! (1 / 0 => +inf)
+ $x->binf();
+ }
# calculate $x->{_m} ** $y and $x->{_e} * $y separately (faster)
$y1->babs();
my $z = $x->copy(); $x->bzero()->binc();
return $x->bdiv($z,$a,$p,$r); # round in one go (might ignore y's A!)
}
- return $x->round($a,$p,$r,$y);
+ $x->round($a,$p,$r,$y);
}
###############################################################################
return $x->bzero() if $scale < $zad;
if ($scale == $zad) # for 0.006, scale -3 and trunc
{
- $scale = -$len-1;
+ $scale = -$len;
}
else
{
{
# 123 => 100 means length(123) = 3 - $scale (2) => 1
- # calculate digits before dot
- my $dbt = $x->{_m}->length(); $dbt += $x->{_e} if $x->{_e}->sign() eq '-';
- # if not enough digits before dot, round to zero
- return $x->bzero() if ($scale > $dbt) && ($dbt < 0);
- # scale always >= 0 here
- if ($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;
- }
+ my $dbt = $x->{_m}->length();
+ # digits before dot
+ my $dbd = $dbt + $x->{_e};
+ # should be the same, so treat it as this
+ $scale = 1 if $scale == 0;
+ # shortcut if already integer
+ return $x if $scale == 1 && $dbt <= $dbd;
+ # maximum digits before dot
+ ++$dbd;
+
+ if ($scale > $dbd)
+ {
+ # not enough digits before dot, so round to zero
+ return $x->bzero;
+ }
+ elsif ( $scale == $dbd )
+ {
+ # maximum
+ $scale = -$dbt;
+ }
else
- {
- $scale = $x->{_m}->length() - $scale;
- }
+ {
+ $scale = $dbd - $scale;
+ }
+
}
# print "using $scale for $x->{_m} with '$mode'\n";
# pass sign to bround for rounding modes '+inf' and '-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 '-';
+ #$x->{_m}->brsft(-$x->{_e},10);
+ #$x->{_e}->bzero();
+ #$x-- if $x->{sign} eq '-';
+
+ $x->{_e}->{sign} = '+'; # negate e
+ $x->{_m}->brsft($x->{_e},10); # cut off digits after dot
+ $x->{_e}->bzero(); # trunc/norm
+ $x->{_m}->binc() if $x->{sign} eq '-'; # decrement if negative
}
$x->round($a,$p,$r);
}
# 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 '+';
+ #$x->{_m}->brsft(-$x->{_e},10);
+ #$x->{_e}->bzero();
+ #$x++ if $x->{sign} eq '+';
+
+ $x->{_e}->{sign} = '+'; # negate e
+ $x->{_m}->brsft($x->{_e},10); # cut off digits after dot
+ $x->{_e}->bzero(); # trunc/norm
+ $x->{_m}->binc() if $x->{sign} eq '+'; # decrement if negative
}
$x->round($a,$p,$r);
}
sub AUTOLOAD
{
- # make fxxx and bxxx work
- # my $self = $_[0];
+ # make fxxx and bxxx both work by selectively mapping fxxx() to MBF::bxxx()
+ # or falling back to MBI::bxxx()
my $name = $AUTOLOAD;
$name =~ s/.*:://; # split package
- #print "$name\n";
no strict 'refs';
if (!method_alias($name))
{
return &{'Math::BigInt'."::$name"}(@_);
}
my $bname = $name; $bname =~ s/^f/b/;
- *{$class."\:\:$name"} = \&$bname;
+ *{$class."::$name"} = \&$bname;
&$bname; # uses @_
}
sub import
{
my $self = shift;
- for ( my $i = 0; $i < @_ ; $i++ )
+ my $l = scalar @_; my $j = 0; my @a = @_;
+ for ( my $i = 0; $i < $l ; $i++, $j++)
{
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;
+ splice @a, $j, 1; $j--;
+ }
+ elsif ($_[$i] eq 'upgrade')
+ {
+ # this causes upgrading
+ $upgrade = $_[$i+1]; # or undef to disable
+ my $s = 2; $s = 1 if @a-$j < 2; # avoid "can not modify non-existant..."
+ splice @a, $j, $s; $j -= $s;
+ }
+ elsif ($_[$i] eq 'downgrade')
+ {
+ # this causes downgrading
+ $downgrade = $_[$i+1]; # or undef to disable
+ my $s = 2; $s = 1 if @a-$j < 2; # avoid "can not modify non-existant..."
+ splice @a, $j, $s; $j -= $s;
}
}
# any non :constant stuff is handled by our parent, Exporter
# even if @_ is empty, to give it a chance
- $self->SUPER::import(@_); # for subclasses
- $self->export_to_level(1,$self,@_); # need this, too
+ $self->SUPER::import(@a); # for subclasses
+ $self->export_to_level(1,$self,@a); # need this, too
}
sub bnorm
return $x if $x->{sign} !~ /^[+-]$/; # inf, nan etc
- 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, and -0 => +0
- $x->{sign} = '+', $x->{_e}->bone() if $x->{_m}->is_zero();
+# if (!$x->{_m}->is_odd())
+# {
+ my $zeros = $x->{_m}->_trailing_zeros(); # correct for trailing zeros
+ if ($zeros != 0)
+ {
+ $x->{_m}->brsft($zeros,10); $x->{_e}->badd($zeros);
+ }
+ # for something like 0Ey, set y to 1, and -0 => +0
+ $x->{sign} = '+', $x->{_e}->bone() if $x->{_m}->is_zero();
+# }
# this is to prevent automatically rounding when MBI's globals are set
$x->{_m}->{_f} = MB_NEVER_ROUND;
$x->{_e}->{_f} = MB_NEVER_ROUND;
# return copy as a bigint representation of this BigFloat number
my ($self,$x) = ref($_[0]) ? (ref($_[0]),$_[0]) : objectify(1,@_);
- my $z;
- if ($x->{_e}->is_zero())
+ my $z = $x->{_m}->copy();
+ if ($x->{_e}->{sign} eq '-') # < 0
{
- $z = $x->{_m}->copy();
- $z->{sign} = $x->{sign};
- return $z;
- }
- $z = $x->{_m}->copy();
- if ($x->{_e} < 0)
- {
- $z->brsft(-$x->{_e},10);
+ $x->{_e}->{sign} = '+'; # flip
+ $z->brsft($x->{_e},10);
+ $x->{_e}->{sign} = '-'; # flip back
}
- else
+ elsif (!$x->{_e}->is_zero()) # > 0
{
$z->blsft($x->{_e},10);
}
use Math::BigFloat;
- # Number creation
- $x = Math::BigInt->new($str); # defaults to 0
- $nan = Math::BigInt->bnan(); # create a NotANumber
- $zero = Math::BigInt->bzero();# create a "+0"
+ # Number creation
+ $x = Math::BigFloat->new($str); # defaults to 0
+ $nan = Math::BigFloat->bnan(); # create a NotANumber
+ $zero = Math::BigFloat->bzero(); # create a +0
+ $inf = Math::BigFloat->binf(); # create a +inf
+ $inf = Math::BigFloat->binf('-'); # create a -inf
+ $one = Math::BigFloat->bone(); # create a +1
+ $one = Math::BigFloat->bone('-'); # create a -1
# Testing
- $x->is_zero(); # return whether arg is zero or not
- $x->is_nan(); # return whether arg is NaN or not
+ $x->is_zero(); # true if arg is +0
+ $x->is_nan(); # true if arg is NaN
$x->is_one(); # true if arg is +1
$x->is_one('-'); # true if arg is -1
$x->is_odd(); # true if odd, false for even
$x->is_even(); # true if even, false for odd
$x->is_positive(); # true if >= 0
$x->is_negative(); # true if < 0
- $x->is_inf(sign) # true if +inf or -inf (sign default '+')
+ $x->is_inf(sign); # true if +inf, or -inf (default is '+')
+
$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 $i to 0
$x->bnan(); # set $i to NaN
+ $x->bone(); # set $x to +1
+ $x->bone('-'); # set $x to -1
+ $x->binf(); # set $x to inf
+ $x->binf('-'); # set $x to -inf
$x->bneg(); # negation
$x->babs(); # absolute value
$x->bior($y); # bit-wise inclusive or
$x->bxor($y); # bit-wise exclusive or
$x->bnot(); # bit-wise not (two's complement)
-
+
+ $x->bsqrt(); # calculate square-root
+ $x->bfac(); # factorial of $x (1*2*3*4*..$x)
+
$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->bfloor(); # return integer less or equal than $x
+ $x->bceil(); # return integer greater or equal than $x
+
$x->exponent(); # return exponent as BigInt
$x->mantissa(); # return mantissa as BigInt
$x->parts(); # return (mantissa,exponent) as BigInt
projects / p5sagit/p5-mst-13.2.git / blobdiff