# # "Tax the rat farms." - Lord Vetinari # # The following hash values are used: # sign : +,-,NaN,+inf,-inf # _d : denominator # _n : numeraotr (value = _n/_d) # _a : accuracy # _p : precision # _f : flags, used by MBR to flag parts of a rational as untouchable # You should not look at the innards of a BigRat - use the methods for this. package Math::BigRat; require 5.005_03; use strict; require Exporter; use Math::BigFloat; use vars qw($VERSION @ISA $upgrade $downgrade $accuracy $precision $round_mode $div_scale $_trap_nan $_trap_inf); @ISA = qw(Exporter Math::BigFloat); $VERSION = '0.13'; use overload; # inherit overload from Math::BigFloat BEGIN { *objectify = \&Math::BigInt::objectify; # inherit this from BigInt *AUTOLOAD = \&Math::BigFloat::AUTOLOAD; # can't inherit AUTOLOAD # we inherit these from BigFloat because currently it is not possible # that MBF has a different $MBI variable than we, because MBF also uses # Math::BigInt::config->('lib'); (there is always only one library loaded) *_e_add = \&Math::BigFloat::_e_add; *_e_sub = \&Math::BigFloat::_e_sub; } ############################################################################## # Global constants and flags. Access these only via the accessor methods! $accuracy = $precision = undef; $round_mode = 'even'; $div_scale = 40; $upgrade = undef; $downgrade = undef; # These are internally, and not to be used from the outside at all! $_trap_nan = 0; # are NaNs ok? set w/ config() $_trap_inf = 0; # are infs ok? set w/ config() # the package we are using for our private parts, defaults to: # Math::BigInt->config()->{lib} my $MBI = 'Math::BigInt::Calc'; my $nan = 'NaN'; my $class = 'Math::BigRat'; my $IMPORT = 0; sub isa { return 0 if $_[1] =~ /^Math::Big(Int|Float)/; # we aren't UNIVERSAL::isa(@_); } ############################################################################## sub _new_from_float { # turn a single float input into a rational number (like '0.1') my ($self,$f) = @_; return $self->bnan() if $f->is_nan(); return $self->binf($f->{sign}) if $f->{sign} =~ /^[+-]inf$/; $self->{_n} = $MBI->_copy( $f->{_m} ); # mantissa $self->{_d} = $MBI->_one(); $self->{sign} = $f->{sign} || '+'; if ($f->{_es} eq '-') { # something like Math::BigRat->new('0.1'); # 1 / 1 => 1/10 $MBI->_lsft ( $self->{_d}, $f->{_e} ,10); } else { # something like Math::BigRat->new('10'); # 1 / 1 => 10/1 $MBI->_lsft ( $self->{_n}, $f->{_e} ,10) unless $MBI->_is_zero($f->{_e}); } $self; } sub new { # create a Math::BigRat my $class = shift; my ($n,$d) = shift; my $self = { }; bless $self,$class; # input like (BigInt,BigInt) or (BigFloat,BigFloat) not handled yet if ((!defined $d) && (ref $n) && (!$n->isa('Math::BigRat'))) { if ($n->isa('Math::BigFloat')) { $self->_new_from_float($n); } if ($n->isa('Math::BigInt')) { # TODO: trap NaN, inf $self->{_n} = $MBI->_copy($n->{value}); # "mantissa" = $n $self->{_d} = $MBI->_one(); # d => 1 $self->{sign} = $n->{sign}; } if ($n->isa('Math::BigInt::Lite')) { # TODO: trap NaN, inf $self->{sign} = '+'; $self->{sign} = '-' if $$n < 0; $self->{_n} = $MBI->_new(abs($$n)); # "mantissa" = $n $self->{_d} = $MBI->_one(); # d => 1 } return $self->bnorm(); # normalize (120/1 => 12/10) } return $n->copy() if ref $n; # already a BigRat if (!defined $n) { $self->{_n} = $MBI->_zero(); # undef => 0 $self->{_d} = $MBI->_one(); $self->{sign} = '+'; return $self; } # string input with / delimiter if ($n =~ /\s*\/\s*/) { return $class->bnan() if $n =~ /\/.*\//; # 1/2/3 isn't valid return $class->bnan() if $n =~ /\/\s*$/; # 1/ isn't valid ($n,$d) = split (/\//,$n); # try as BigFloats first if (($n =~ /[\.eE]/) || ($d =~ /[\.eE]/)) { local $Math::BigFloat::accuracy = undef; local $Math::BigFloat::precision = undef; # one of them looks like a float my $nf = Math::BigFloat->new($n,undef,undef); $self->{sign} = '+'; return $self->bnan() if $nf->is_nan(); $self->{_n} = $MBI->_copy( $nf->{_m} ); # get mantissa # now correct $self->{_n} due to $n my $f = Math::BigFloat->new($d,undef,undef); return $self->bnan() if $f->is_nan(); $self->{_d} = $MBI->_copy( $f->{_m} ); # calculate the difference between nE and dE # XXX TODO: check that exponent() makes a copy to avoid copy() my $diff_e = $nf->exponent()->copy()->bsub( $f->exponent); if ($diff_e->is_negative()) { # < 0: mul d with it $MBI->_lsft( $self->{_d}, $MBI->_new( $diff_e->babs()), 10); } elsif (!$diff_e->is_zero()) { # > 0: mul n with it $MBI->_lsft( $self->{_n}, $MBI->_new( $diff_e), 10); } } else { # both d and n look like (big)ints $self->{sign} = '+'; # no sign => '+' $self->{_n} = undef; $self->{_d} = undef; if ($n =~ /^([+-]?)0*(\d+)\z/) # first part ok? { $self->{sign} = $1 || '+'; # no sign => '+' $self->{_n} = $MBI->_new($2 || 0); } if ($d =~ /^([+-]?)0*(\d+)\z/) # second part ok? { $self->{sign} =~ tr/+-/-+/ if ($1 || '') eq '-'; # negate if second part neg. $self->{_d} = $MBI->_new($2 || 0); } if (!defined $self->{_n} || !defined $self->{_d}) { $d = Math::BigInt->new($d,undef,undef) unless ref $d; $n = Math::BigInt->new($n,undef,undef) unless ref $n; if ($n->{sign} =~ /^[+-]$/ && $d->{sign} =~ /^[+-]$/) { # both parts are ok as integers (wierd things like ' 1e0' $self->{_n} = $MBI->_copy($n->{value}); $self->{_d} = $MBI->_copy($d->{value}); $self->{sign} = $n->{sign}; $self->{sign} =~ tr/+-/-+/ if $d->{sign} eq '-'; # -1/-2 => 1/2 return $self->bnorm(); } $self->{sign} = '+'; # a default sign return $self->bnan() if $n->is_nan() || $d->is_nan(); # handle inf cases: if ($n->is_inf() || $d->is_inf()) { if ($n->is_inf()) { return $self->bnan() if $d->is_inf(); # both are inf => NaN my $s = '+'; # '+inf/+123' or '-inf/-123' $s = '-' if substr($n->{sign},0,1) ne $d->{sign}; # +-inf/123 => +-inf return $self->binf($s); } # 123/inf => 0 return $self->bzero(); } } } return $self->bnorm(); } # simple string input if (($n =~ /[\.eE]/)) { # looks like a float, quacks like a float, so probably is a float $self->{sign} = 'NaN'; local $Math::BigFloat::accuracy = undef; local $Math::BigFloat::precision = undef; $self->_new_from_float(Math::BigFloat->new($n,undef,undef)); } else { # for simple forms, use $MBI directly if ($n =~ /^([+-]?)0*(\d+)\z/) { $self->{sign} = $1 || '+'; $self->{_n} = $MBI->_new($2 || 0); $self->{_d} = $MBI->_one(); } else { my $n = Math::BigInt->new($n,undef,undef); $self->{_n} = $MBI->_copy($n->{value}); $self->{_d} = $MBI->_one(); $self->{sign} = $n->{sign}; return $self->bnan() if $self->{sign} eq 'NaN'; return $self->binf($self->{sign}) if $self->{sign} =~ /^[+-]inf$/; } } $self->bnorm(); } 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 {}, $c; $self->{sign} = $x->{sign}; $self->{_d} = $MBI->_copy($x->{_d}); $self->{_n} = $MBI->_copy($x->{_n}); $self->{_a} = $x->{_a} if defined $x->{_a}; $self->{_p} = $x->{_p} if defined $x->{_p}; $self; } ############################################################################## sub config { # return (later set?) configuration data as hash ref my $class = shift || 'Math::BigFloat'; my $cfg = $class->SUPER::config(@_); # now we need only to override the ones that are different from our parent $cfg->{class} = $class; $cfg->{with} = $MBI; $cfg; } ############################################################################## sub bstr { my ($self,$x) = ref($_[0]) ? (undef,$_[0]) : objectify(1,@_); if ($x->{sign} !~ /^[+-]$/) # inf, NaN etc { my $s = $x->{sign}; $s =~ s/^\+//; # +inf => inf return $s; } my $s = ''; $s = $x->{sign} if $x->{sign} ne '+'; # '+3/2' => '3/2' return $s . $MBI->_str($x->{_n}) if $MBI->_is_one($x->{_d}); $s . $MBI->_str($x->{_n}) . '/' . $MBI->_str($x->{_d}); } sub bsstr { my ($self,$x) = ref($_[0]) ? (ref($_[0]),$_[0]) : objectify(1,@_); if ($x->{sign} !~ /^[+-]$/) # inf, NaN etc { my $s = $x->{sign}; $s =~ s/^\+//; # +inf => inf return $s; } my $s = ''; $s = $x->{sign} if $x->{sign} ne '+'; # +3 vs 3 $s . $MBI->_str($x->{_n}) . '/' . $MBI->_str($x->{_d}); } sub bnorm { # reduce the number to the shortest form my ($self,$x) = ref($_[0]) ? (ref($_[0]),$_[0]) : objectify(1,@_); # Both parts must be objects of whatever we are using today. # Second check because Calc.pm has ARRAY res as unblessed objects. if (ref($x->{_n}) ne $MBI && ref($x->{_n}) ne 'ARRAY') { require Carp; Carp::croak ("n is not $MBI but (".ref($x->{_n}).') in bnorm()'); } if (ref($x->{_d}) ne $MBI && ref($x->{_d}) ne 'ARRAY') { require Carp; Carp::croak ("d is not $MBI but (".ref($x->{_d}).') in bnorm()'); } # no normalize for NaN, inf etc. return $x if $x->{sign} !~ /^[+-]$/; # normalize zeros to 0/1 if ($MBI->_is_zero($x->{_n})) { $x->{sign} = '+'; # never leave a -0 $x->{_d} = $MBI->_one() unless $MBI->_is_one($x->{_d}); return $x; } return $x if $MBI->_is_one($x->{_d}); # no need to reduce # reduce other numbers my $gcd = $MBI->_copy($x->{_n}); $gcd = $MBI->_gcd($gcd,$x->{_d}); if (!$MBI->_is_one($gcd)) { $x->{_n} = $MBI->_div($x->{_n},$gcd); $x->{_d} = $MBI->_div($x->{_d},$gcd); } $x; } ############################################################################## # special values sub _bnan { # used by parent class bnan() to initialize number to NaN my $self = shift; if ($_trap_nan) { require Carp; my $class = ref($self); Carp::croak ("Tried to set $self to NaN in $class\::_bnan()"); } $self->{_n} = $MBI->_zero(); $self->{_d} = $MBI->_zero(); } sub _binf { # used by parent class bone() to initialize number to +inf/-inf my $self = shift; if ($_trap_inf) { require Carp; my $class = ref($self); Carp::croak ("Tried to set $self to inf in $class\::_binf()"); } $self->{_n} = $MBI->_zero(); $self->{_d} = $MBI->_zero(); } sub _bone { # used by parent class bone() to initialize number to +1/-1 my $self = shift; $self->{_n} = $MBI->_one(); $self->{_d} = $MBI->_one(); } sub _bzero { # used by parent class bzero() to initialize number to 0 my $self = shift; $self->{_n} = $MBI->_zero(); $self->{_d} = $MBI->_one(); } ############################################################################## # mul/add/div etc sub badd { # add two rational numbers # set up parameters my ($self,$x,$y,@r) = (ref($_[0]),@_); # objectify is costly, so avoid it if ((!ref($_[0])) || (ref($_[0]) ne ref($_[1]))) { ($self,$x,$y,@r) = objectify(2,@_); } # +inf + +inf => +inf, -inf + -inf => -inf return $x->binf(substr($x->{sign},0,1)) if $x->{sign} eq $y->{sign} && $x->{sign} =~ /^[+-]inf$/; # +inf + -inf or -inf + +inf => NaN return $x->bnan() if ($x->{sign} !~ /^[+-]$/ || $y->{sign} !~ /^[+-]$/); # 1 1 gcd(3,4) = 1 1*3 + 1*4 7 # - + - = --------- = -- # 4 3 4*3 12 # we do not compute the gcd() here, but simple do: # 5 7 5*3 + 7*4 41 # - + - = --------- = -- # 4 3 4*3 12 # and bnorm() will then take care of the rest $x->{_n} = $MBI->_mul( $x->{_n}, $y->{_d}); my $m = $MBI->_mul( $MBI->_copy( $y->{_n} ), $x->{_d} ); ($x->{_n}, $x->{sign}) = _e_add( $x->{_n}, $m, $x->{sign}, $y->{sign}); $x->{_d} = $MBI->_mul( $x->{_d}, $y->{_d}); # normalize and round $x->bnorm()->round(@r); } sub bsub { # subtract two rational numbers # set up parameters my ($self,$x,$y,@r) = (ref($_[0]),@_); # objectify is costly, so avoid it if ((!ref($_[0])) || (ref($_[0]) ne ref($_[1]))) { ($self,$x,$y,@r) = objectify(2,@_); } # flip sign of $x, call badd(), then flip sign of result $x->{sign} =~ tr/+-/-+/ unless $x->{sign} eq '+' && $MBI->_is_zero($x->{_n}); # not -0 $x->badd($y,@r); # does norm and round $x->{sign} =~ tr/+-/-+/ unless $x->{sign} eq '+' && $MBI->_is_zero($x->{_n}); # not -0 $x; } sub bmul { # multiply two rational numbers # set up parameters my ($self,$x,$y,@r) = (ref($_[0]),@_); # objectify is costly, so avoid it if ((!ref($_[0])) || (ref($_[0]) ne ref($_[1]))) { ($self,$x,$y,@r) = objectify(2,@_); } return $x->bnan() if ($x->{sign} eq 'NaN' || $y->{sign} eq 'NaN'); # 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() if ($x->{sign} =~ /^-/ && $y->{sign} =~ /^-/); return $x->binf('-'); } # x== 0 # also: or y == 1 or y == -1 return wantarray ? ($x,$self->bzero()) : $x if $x->is_zero(); # XXX TODO: # According to Knuth, this can be optimized by doing gcd twice (for d and n) # and reducing in one step. This would save us the bnorm() at the end. # 1 2 1 * 2 2 1 # - * - = ----- = - = - # 4 3 4 * 3 12 6 $x->{_n} = $MBI->_mul( $x->{_n}, $y->{_n}); $x->{_d} = $MBI->_mul( $x->{_d}, $y->{_d}); # compute new sign $x->{sign} = $x->{sign} eq $y->{sign} ? '+' : '-'; $x->bnorm()->round(@r); } sub bdiv { # (dividend: BRAT or num_str, divisor: BRAT or num_str) return # (BRAT,BRAT) (quo,rem) or BRAT (only rem) # set up parameters my ($self,$x,$y,@r) = (ref($_[0]),@_); # objectify is costly, so avoid it if ((!ref($_[0])) || (ref($_[0]) ne ref($_[1]))) { ($self,$x,$y,@r) = objectify(2,@_); } return $self->_div_inf($x,$y) if (($x->{sign} !~ /^[+-]$/) || ($y->{sign} !~ /^[+-]$/) || $y->is_zero()); # x== 0 # also: or y == 1 or y == -1 return wantarray ? ($x,$self->bzero()) : $x if $x->is_zero(); # XXX TODO: list context, upgrade # According to Knuth, this can be optimized by doing gcd twice (for d and n) # and reducing in one step. This would save us the bnorm() at the end. # 1 1 1 3 # - / - == - * - # 4 3 4 1 $x->{_n} = $MBI->_mul( $x->{_n}, $y->{_d}); $x->{_d} = $MBI->_mul( $x->{_d}, $y->{_n}); # compute new sign $x->{sign} = $x->{sign} eq $y->{sign} ? '+' : '-'; $x->bnorm()->round(@r); $x; } sub bmod { # compute "remainder" (in Perl way) of $x / $y # set up parameters my ($self,$x,$y,@r) = (ref($_[0]),@_); # objectify is costly, so avoid it if ((!ref($_[0])) || (ref($_[0]) ne ref($_[1]))) { ($self,$x,$y,@r) = objectify(2,@_); } return $self->_div_inf($x,$y) if (($x->{sign} !~ /^[+-]$/) || ($y->{sign} !~ /^[+-]$/) || $y->is_zero()); return $x if $x->is_zero(); # 0 / 7 = 0, mod 0 # compute $x - $y * floor($x/$y), keeping the sign of $x # copy x to u, make it positive and then do a normal division ($u/$y) my $u = bless { sign => '+' }, $self; $u->{_n} = $MBI->_mul( $MBI->_copy($x->{_n}), $y->{_d} ); $u->{_d} = $MBI->_mul( $MBI->_copy($x->{_d}), $y->{_n} ); # compute floor(u) if (! $MBI->_is_one($u->{_d})) { $u->{_n} = $MBI->_div($u->{_n},$u->{_d}); # 22/7 => 3/1 w/ truncate # no need to set $u->{_d} to 1, since below we set it to $y->{_d} anyway } # now compute $y * $u $u->{_d} = $MBI->_copy($y->{_d}); # 1 * $y->{_d}, see floor above $u->{_n} = $MBI->_mul($u->{_n},$y->{_n}); my $xsign = $x->{sign}; $x->{sign} = '+'; # remember sign and make x positive # compute $x - $u $x->bsub($u); $x->{sign} = $xsign; # put sign back $x->bnorm()->round(@r); } ############################################################################## # bdec/binc sub bdec { # decrement value (subtract 1) my ($self,$x,@r) = ref($_[0]) ? (ref($_[0]),@_) : objectify(1,@_); return $x if $x->{sign} !~ /^[+-]$/; # NaN, inf, -inf if ($x->{sign} eq '-') { $x->{_n} = $MBI->_add( $x->{_n}, $x->{_d}); # -5/2 => -7/2 } else { if ($MBI->_acmp($x->{_n},$x->{_d}) < 0) # n < d? { # 1/3 -- => -2/3 $x->{_n} = $MBI->_sub( $MBI->_copy($x->{_d}), $x->{_n}); $x->{sign} = '-'; } else { $x->{_n} = $MBI->_sub($x->{_n}, $x->{_d}); # 5/2 => 3/2 } } $x->bnorm()->round(@r); } sub binc { # increment value (add 1) my ($self,$x,@r) = ref($_[0]) ? (ref($_[0]),@_) : objectify(1,@_); return $x if $x->{sign} !~ /^[+-]$/; # NaN, inf, -inf if ($x->{sign} eq '-') { if ($MBI->_acmp($x->{_n},$x->{_d}) < 0) { # -1/3 ++ => 2/3 (overflow at 0) $x->{_n} = $MBI->_sub( $MBI->_copy($x->{_d}), $x->{_n}); $x->{sign} = '+'; } else { $x->{_n} = $MBI->_sub($x->{_n}, $x->{_d}); # -5/2 => -3/2 } } else { $x->{_n} = $MBI->_add($x->{_n},$x->{_d}); # 5/2 => 7/2 } $x->bnorm()->round(@r); } ############################################################################## # is_foo methods (the rest is inherited) sub is_int { # return true if arg (BRAT or num_str) is an integer my ($self,$x) = ref($_[0]) ? (undef,$_[0]) : objectify(1,@_); return 1 if ($x->{sign} =~ /^[+-]$/) && # NaN and +-inf aren't $MBI->_is_one($x->{_d}); # x/y && y != 1 => no integer 0; } sub is_zero { # return true if arg (BRAT or num_str) is zero my ($self,$x) = ref($_[0]) ? (undef,$_[0]) : objectify(1,@_); return 1 if $x->{sign} eq '+' && $MBI->_is_zero($x->{_n}); 0; } sub is_one { # return true if arg (BRAT or num_str) is +1 or -1 if signis given my ($self,$x) = ref($_[0]) ? (undef,$_[0]) : objectify(1,@_); my $sign = $_[2] || ''; $sign = '+' if $sign ne '-'; return 1 if ($x->{sign} eq $sign && $MBI->_is_one($x->{_n}) && $MBI->_is_one($x->{_d})); 0; } sub is_odd { # return true if arg (BFLOAT or num_str) is odd or false if even my ($self,$x) = ref($_[0]) ? (undef,$_[0]) : objectify(1,@_); return 1 if ($x->{sign} =~ /^[+-]$/) && # NaN & +-inf aren't ($MBI->_is_one($x->{_d}) && $MBI->_is_odd($x->{_n})); # x/2 is not, but 3/1 0; } sub is_even { # return true if arg (BINT or num_str) is even or false if odd my ($self,$x) = ref($_[0]) ? (undef,$_[0]) : objectify(1,@_); return 0 if $x->{sign} !~ /^[+-]$/; # NaN & +-inf aren't return 1 if ($MBI->_is_one($x->{_d}) # x/3 is never && $MBI->_is_even($x->{_n})); # but 4/1 is 0; } ############################################################################## # parts() and friends sub numerator { my ($self,$x) = ref($_[0]) ? (ref($_[0]),$_[0]) : objectify(1,@_); # NaN, inf, -inf return Math::BigInt->new($x->{sign}) if ($x->{sign} !~ /^[+-]$/); my $n = Math::BigInt->new($MBI->_str($x->{_n})); $n->{sign} = $x->{sign}; $n; } sub denominator { my ($self,$x) = ref($_[0]) ? (ref($_[0]),$_[0]) : objectify(1,@_); # NaN return Math::BigInt->new($x->{sign}) if $x->{sign} eq 'NaN'; # inf, -inf return Math::BigInt->bone() if $x->{sign} !~ /^[+-]$/; Math::BigInt->new($MBI->_str($x->{_d})); } sub parts { my ($self,$x) = ref($_[0]) ? (ref($_[0]),$_[0]) : objectify(1,@_); my $c = 'Math::BigInt'; return ($c->bnan(),$c->bnan()) if $x->{sign} eq 'NaN'; return ($c->binf(),$c->binf()) if $x->{sign} eq '+inf'; return ($c->binf('-'),$c->binf()) if $x->{sign} eq '-inf'; my $n = $c->new( $MBI->_str($x->{_n})); $n->{sign} = $x->{sign}; my $d = $c->new( $MBI->_str($x->{_d})); ($n,$d); } sub length { my ($self,$x) = ref($_[0]) ? (undef,$_[0]) : objectify(1,@_); return $nan unless $x->is_int(); $MBI->_len($x->{_n}); # length(-123/1) => length(123) } sub digit { my ($self,$x,$n) = ref($_[0]) ? (undef,$_[0],$_[1]) : objectify(1,@_); return $nan unless $x->is_int(); $MBI->_digit($x->{_n},$n || 0); # digit(-123/1,2) => digit(123,2) } ############################################################################## # special calc routines sub bceil { my ($self,$x) = ref($_[0]) ? (ref($_[0]),$_[0]) : objectify(1,@_); return $x if $x->{sign} !~ /^[+-]$/ || # not for NaN, inf $MBI->_is_one($x->{_d}); # 22/1 => 22, 0/1 => 0 $x->{_n} = $MBI->_div($x->{_n},$x->{_d}); # 22/7 => 3/1 w/ truncate $x->{_d} = $MBI->_one(); # d => 1 $x->{_n} = $MBI->_inc($x->{_n}) if $x->{sign} eq '+'; # +22/7 => 4/1 $x->{sign} = '+' if $MBI->_is_zero($x->{_n}); # -0 => 0 $x; } sub bfloor { my ($self,$x) = ref($_[0]) ? (ref($_[0]),$_[0]) : objectify(1,@_); return $x if $x->{sign} !~ /^[+-]$/ || # not for NaN, inf $MBI->_is_one($x->{_d}); # 22/1 => 22, 0/1 => 0 $x->{_n} = $MBI->_div($x->{_n},$x->{_d}); # 22/7 => 3/1 w/ truncate $x->{_d} = $MBI->_one(); # d => 1 $x->{_n} = $MBI->_inc($x->{_n}) if $x->{sign} eq '-'; # -22/7 => -4/1 $x; } sub bfac { my ($self,$x,@r) = ref($_[0]) ? (ref($_[0]),@_) : objectify(1,@_); # if $x is not an integer if (($x->{sign} ne '+') || (!$MBI->_is_one($x->{_d}))) { return $x->bnan(); } $x->{_n} = $MBI->_fac($x->{_n}); # since _d is 1, we don't need to reduce/norm the result $x->round(@r); } sub bpow { # power ($x ** $y) # set up parameters my ($self,$x,$y,@r) = (ref($_[0]),@_); # objectify is costly, so avoid it if ((!ref($_[0])) || (ref($_[0]) ne ref($_[1]))) { ($self,$x,$y,@r) = objectify(2,@_); } return $x if $x->{sign} =~ /^[+-]inf$/; # -inf/+inf ** x return $x->bnan() if $x->{sign} eq $nan || $y->{sign} eq $nan; return $x->bone(@r) if $y->is_zero(); return $x->round(@r) if $x->is_one() || $y->is_one(); if ($x->{sign} eq '-' && $MBI->_is_one($x->{_n}) && $MBI->_is_one($x->{_d})) { # if $x == -1 and odd/even y => +1/-1 return $y->is_odd() ? $x->round(@r) : $x->babs()->round(@r); # my Casio FX-5500L has a bug here: -1 ** 2 is -1, but -1 * -1 is 1; } # 1 ** -y => 1 / (1 ** |y|) # so do test for negative $y after above's clause return $x->round(@r) if $x->is_zero(); # 0**y => 0 (if not y <= 0) # shortcut y/1 (and/or x/1) if ($MBI->_is_one($y->{_d})) { # shortcut for x/1 and y/1 if ($MBI->_is_one($x->{_d})) { $x->{_n} = $MBI->_pow($x->{_n},$y->{_n}); # x/1 ** y/1 => (x ** y)/1 if ($y->{sign} eq '-') { # 0.2 ** -3 => 1/(0.2 ** 3) ($x->{_n},$x->{_d}) = ($x->{_d},$x->{_n}); # swap } # correct sign; + ** + => + if ($x->{sign} eq '-') { # - * - => +, - * - * - => - $x->{sign} = '+' if $MBI->_is_even($y->{_n}); } return $x->round(@r); } # x/z ** y/1 $x->{_n} = $MBI->_pow($x->{_n},$y->{_n}); # 5/2 ** y/1 => 5 ** y / 2 ** y $x->{_d} = $MBI->_pow($x->{_d},$y->{_n}); if ($y->{sign} eq '-') { # 0.2 ** -3 => 1/(0.2 ** 3) ($x->{_n},$x->{_d}) = ($x->{_d},$x->{_n}); # swap } # correct sign; + ** + => + if ($x->{sign} eq '-') { # - * - => +, - * - * - => - $x->{sign} = '+' if $MBI->_is_even($y->{_n}); } return $x->round(@r); } # regular calculation (this is wrong for d/e ** f/g) my $pow2 = $self->bone(); my $y1 = $MBI->_div ( $MBI->_copy($y->{_n}), $y->{_d}); my $two = $MBI->_two(); while (!$MBI->_is_one($y1)) { $pow2->bmul($x) if $MBI->_is_odd($y1); $MBI->_div($y1, $two); $x->bmul($x); } $x->bmul($pow2) unless $pow2->is_one(); # n ** -x => 1/n ** x ($x->{_d},$x->{_n}) = ($x->{_n},$x->{_d}) if $y->{sign} eq '-'; $x->bnorm()->round(@r); } sub blog { # set up parameters my ($self,$x,$y,@r) = (ref($_[0]),@_); # objectify is costly, so avoid it if ((!ref($_[0])) || (ref($_[0]) ne ref($_[1]))) { ($self,$x,$y,@r) = objectify(2,$class,@_); } # blog(1,Y) => 0 return $x->bzero() if $x->is_one() && $y->{sign} eq '+'; # $x <= 0 => NaN return $x->bnan() if $x->is_zero() || $x->{sign} ne '+' || $y->{sign} ne '+'; if ($x->is_int() && $y->is_int()) { return $self->new($x->as_number()->blog($y->as_number(),@r)); } # do it with floats $x->_new_from_float( $x->_as_float()->blog(Math::BigFloat->new("$y"),@r) ); } sub _float_from_part { my $x = shift; my $f = Math::BigFloat->bzero(); $f->{_m} = $MBI->_copy($x); $f->{_e} = $MBI->_zero(); $f; } sub _as_float { my $x = shift; local $Math::BigFloat::upgrade = undef; local $Math::BigFloat::accuracy = undef; local $Math::BigFloat::precision = undef; # 22/7 => 3.142857143.. my $a = $x->accuracy() || 0; if ($a != 0 || !$MBI->_is_one($x->{_d})) { # n/d return Math::BigFloat->new($x->{sign} . $MBI->_str($x->{_n}))->bdiv( $MBI->_str($x->{_d}), $x->accuracy()); } # just n Math::BigFloat->new($x->{sign} . $MBI->_str($x->{_n})); } sub broot { # set up parameters my ($self,$x,$y,@r) = (ref($_[0]),@_); # objectify is costly, so avoid it if ((!ref($_[0])) || (ref($_[0]) ne ref($_[1]))) { ($self,$x,$y,@r) = objectify(2,@_); } if ($x->is_int() && $y->is_int()) { return $self->new($x->as_number()->broot($y->as_number(),@r)); } # do it with floats $x->_new_from_float( $x->_as_float()->broot($y,@r) ); } sub bmodpow { # set up parameters my ($self,$x,$y,$m,@r) = (ref($_[0]),@_); # objectify is costly, so avoid it if ((!ref($_[0])) || (ref($_[0]) ne ref($_[1]))) { ($self,$x,$y,$m,@r) = objectify(3,@_); } # $x or $y or $m are NaN or +-inf => NaN return $x->bnan() if $x->{sign} !~ /^[+-]$/ || $y->{sign} !~ /^[+-]$/ || $m->{sign} !~ /^[+-]$/; if ($x->is_int() && $y->is_int() && $m->is_int()) { return $self->new($x->as_number()->bmodpow($y->as_number(),$m,@r)); } warn ("bmodpow() not fully implemented"); $x->bnan(); } sub bmodinv { # set up parameters my ($self,$x,$y,@r) = (ref($_[0]),@_); # objectify is costly, so avoid it if ((!ref($_[0])) || (ref($_[0]) ne ref($_[1]))) { ($self,$x,$y,@r) = objectify(2,@_); } # $x or $y are NaN or +-inf => NaN return $x->bnan() if $x->{sign} !~ /^[+-]$/ || $y->{sign} !~ /^[+-]$/; if ($x->is_int() && $y->is_int()) { return $self->new($x->as_number()->bmodinv($y->as_number(),@r)); } warn ("bmodinv() not fully implemented"); $x->bnan(); } sub bsqrt { my ($self,$x,@r) = ref($_[0]) ? (ref($_[0]),@_) : objectify(1,@_); return $x->bnan() if $x->{sign} !~ /^[+]/; # NaN, -inf or < 0 return $x if $x->{sign} eq '+inf'; # sqrt(inf) == inf return $x->round(@r) if $x->is_zero() || $x->is_one(); local $Math::BigFloat::upgrade = undef; local $Math::BigFloat::downgrade = undef; local $Math::BigFloat::precision = undef; local $Math::BigFloat::accuracy = undef; local $Math::BigInt::upgrade = undef; local $Math::BigInt::precision = undef; local $Math::BigInt::accuracy = undef; $x->{_n} = _float_from_part( $x->{_n} )->bsqrt(); $x->{_d} = _float_from_part( $x->{_d} )->bsqrt(); # XXX TODO: we probably can optimze this: # if sqrt(D) was not integer if ($x->{_d}->{_es} ne '+') { $x->{_n}->blsft($x->{_d}->exponent()->babs(),10); # 7.1/4.51 => 7.1/45.1 $x->{_d} = $MBI->_copy( $x->{_d}->{_m} ); # 7.1/45.1 => 71/45.1 } # if sqrt(N) was not integer if ($x->{_n}->{_es} ne '+') { $x->{_d}->blsft($x->{_n}->exponent()->babs(),10); # 71/45.1 => 710/45.1 $x->{_n} = $MBI->_copy( $x->{_n}->{_m} ); # 710/45.1 => 710/451 } # convert parts to $MBI again $x->{_n} = $MBI->_lsft( $MBI->_copy( $x->{_n}->{_m} ), $x->{_n}->{_e}, 10) if ref($x->{_n}) ne $MBI && ref($x->{_n}) ne 'ARRAY'; $x->{_d} = $MBI->_lsft( $MBI->_copy( $x->{_d}->{_m} ), $x->{_d}->{_e}, 10) if ref($x->{_d}) ne $MBI && ref($x->{_d}) ne 'ARRAY'; $x->bnorm()->round(@r); } sub blsft { my ($self,$x,$y,$b,@r) = objectify(3,@_); $b = 2 unless defined $b; $b = $self->new($b) unless ref ($b); $x->bmul( $b->copy()->bpow($y), @r); $x; } sub brsft { my ($self,$x,$y,$b,@r) = objectify(3,@_); $b = 2 unless defined $b; $b = $self->new($b) unless ref ($b); $x->bdiv( $b->copy()->bpow($y), @r); $x; } ############################################################################## # round sub round { $_[0]; } sub bround { $_[0]; } sub bfround { $_[0]; } ############################################################################## # comparing sub bcmp { # compare two signed numbers # set up parameters my ($self,$x,$y) = (ref($_[0]),@_); # objectify is costly, so avoid it if ((!ref($_[0])) || (ref($_[0]) ne ref($_[1]))) { ($self,$x,$y) = objectify(2,@_); } if (($x->{sign} !~ /^[+-]$/) || ($y->{sign} !~ /^[+-]$/)) { # handle +-inf and NaN return undef if (($x->{sign} eq $nan) || ($y->{sign} eq $nan)); return 0 if $x->{sign} eq $y->{sign} && $x->{sign} =~ /^[+-]inf$/; return +1 if $x->{sign} eq '+inf'; return -1 if $x->{sign} eq '-inf'; return -1 if $y->{sign} eq '+inf'; return +1; } # check sign for speed first return 1 if $x->{sign} eq '+' && $y->{sign} eq '-'; # does also 0 <=> -y return -1 if $x->{sign} eq '-' && $y->{sign} eq '+'; # does also -x <=> 0 # shortcut my $xz = $MBI->_is_zero($x->{_n}); my $yz = $MBI->_is_zero($y->{_n}); return 0 if $xz && $yz; # 0 <=> 0 return -1 if $xz && $y->{sign} eq '+'; # 0 <=> +y return 1 if $yz && $x->{sign} eq '+'; # +x <=> 0 my $t = $MBI->_mul( $MBI->_copy($x->{_n}), $y->{_d}); my $u = $MBI->_mul( $MBI->_copy($y->{_n}), $x->{_d}); my $cmp = $MBI->_acmp($t,$u); # signs are equal $cmp = -$cmp if $x->{sign} eq '-'; # both are '-' => reverse $cmp; } sub bacmp { # compare two numbers (as unsigned) # set up parameters my ($self,$x,$y) = (ref($_[0]),@_); # objectify is costly, so avoid it if ((!ref($_[0])) || (ref($_[0]) ne ref($_[1]))) { ($self,$x,$y) = objectify(2,$class,@_); } if (($x->{sign} !~ /^[+-]$/) || ($y->{sign} !~ /^[+-]$/)) { # handle +-inf and NaN return undef if (($x->{sign} eq $nan) || ($y->{sign} eq $nan)); return 0 if $x->{sign} =~ /^[+-]inf$/ && $y->{sign} =~ /^[+-]inf$/; return 1 if $x->{sign} =~ /^[+-]inf$/ && $y->{sign} !~ /^[+-]inf$/; return -1; } my $t = $MBI->_mul( $MBI->_copy($x->{_n}), $y->{_d}); my $u = $MBI->_mul( $MBI->_copy($y->{_n}), $x->{_d}); $MBI->_acmp($t,$u); # ignore signs } ############################################################################## # output conversation sub numify { # convert 17/8 => float (aka 2.125) my ($self,$x) = ref($_[0]) ? (ref($_[0]),$_[0]) : objectify(1,@_); return $x->bstr() if $x->{sign} !~ /^[+-]$/; # inf, NaN, etc # N/1 => N return $MBI->_num($x->{_n}) if $MBI->_is_one($x->{_d}); # N/D my $neg = 1; $neg = -1 if $x->{sign} ne '+'; $neg * $MBI->_num($x->{_n}) / $MBI->_num($x->{_d}); # return sign * N/D } sub as_number { my ($self,$x) = ref($_[0]) ? (undef,$_[0]) : objectify(1,@_); return Math::BigInt->new($x) if $x->{sign} !~ /^[+-]$/; # NaN, inf etc my $u = Math::BigInt->bzero(); $u->{sign} = $x->{sign}; $u->{value} = $MBI->_div( $MBI->_copy($x->{_n}), $x->{_d}); # 22/7 => 3 $u; } sub as_bin { my ($self,$x) = ref($_[0]) ? (undef,$_[0]) : objectify(1,@_); return $x unless $x->is_int(); my $s = $x->{sign}; $s = '' if $s eq '+'; $s . $MBI->_as_bin($x->{_n}); } sub as_hex { my ($self,$x) = ref($_[0]) ? (undef,$_[0]) : objectify(1,@_); return $x unless $x->is_int(); my $s = $x->{sign}; $s = '' if $s eq '+'; $s . $MBI->_as_hex($x->{_n}); } sub import { my $self = shift; my $l = scalar @_; my $lib = ''; my @a; $IMPORT++; for ( my $i = 0; $i < $l ; $i++) { # print "at $_[$i] (",$_[$i+1]||'undef',")\n"; if ( $_[$i] eq ':constant' ) { # this rest causes overlord er load to step in # print "overload @_\n"; overload::constant float => sub { $self->new(shift); }; } # elsif ($_[$i] eq 'upgrade') # { # # this causes upgrading # $upgrade = $_[$i+1]; # or undef to disable # $i++; # } elsif ($_[$i] eq 'downgrade') { # this causes downgrading $downgrade = $_[$i+1]; # or undef to disable $i++; } elsif ($_[$i] eq 'lib') { $lib = $_[$i+1] || ''; # default Calc $i++; } elsif ($_[$i] eq 'with') { $MBI = $_[$i+1] || 'Math::BigInt'; # default Math::BigInt $i++; } else { push @a, $_[$i]; } } # let use Math::BigInt lib => 'GMP'; use Math::BigRat; still work my $mbilib = eval { Math::BigInt->config()->{lib} }; if ((defined $mbilib) && ($MBI eq 'Math::BigInt')) { # MBI already loaded $MBI->import('lib',"$lib,$mbilib", 'objectify'); } else { # MBI not loaded, or not with "Math::BigInt" $lib .= ",$mbilib" if defined $mbilib; if ($] < 5.006) { # Perl < 5.6.0 dies with "out of memory!" when eval() and ':constant' is # used in the same script, or eval inside import(). my @parts = split /::/, $MBI; # Math::BigInt => Math BigInt my $file = pop @parts; $file .= '.pm'; # BigInt => BigInt.pm $file = File::Spec->catfile (@parts, $file); eval { require $file; $MBI->import( lib => '$lib', 'objectify' ); } } else { my $rc = "use $MBI lib => '$lib', 'objectify';"; eval $rc; } } if ($@) { require Carp; Carp::croak ("Couldn't load $MBI: $! $@"); } $MBI = Math::BigFloat->config()->{lib}; # any non :constant stuff is handled by our parent, Exporter # even if @_ is empty, to give it a chance $self->SUPER::import(@a); # for subclasses $self->export_to_level(1,$self,@a); # need this, too } 1; __END__ =head1 NAME Math::BigRat - arbitrarily big rational numbers =head1 SYNOPSIS use Math::BigRat; my $x = Math::BigRat->new('3/7'); $x += '5/9'; print $x->bstr(),"\n"; print $x ** 2,"\n"; my $y = Math::BigRat->new('inf'); print "$y ", ($y->is_inf ? 'is' : 'is not') , " infinity\n"; my $z = Math::BigRat->new(144); $z->bsqrt(); =head1 DESCRIPTION Math::BigRat complements Math::BigInt and Math::BigFloat by providing support for arbitrarily big rational numbers. =head2 MATH LIBRARY Math with the numbers is done (by default) by a module called Math::BigInt::Calc. This is equivalent to saying: use Math::BigRat lib => 'Calc'; You can change this by using: use Math::BigRat lib => 'BitVect'; The following would first try to find Math::BigInt::Foo, then Math::BigInt::Bar, and when this also fails, revert to Math::BigInt::Calc: use Math::BigRat lib => 'Foo,Math::BigInt::Bar'; Calc.pm uses as internal format an array of elements of some decimal base (usually 1e7, but this might be different for some systems) with the least significant digit first, while BitVect.pm uses a bit vector of base 2, most significant bit first. Other modules might use even different means of representing the numbers. See the respective module documentation for further details. Currently the following replacement libraries exist, search for them at CPAN: Math::BigInt::BitVect Math::BigInt::GMP Math::BigInt::Pari Math::BigInt::FastCalc =head1 METHODS Any methods not listed here are dervied from Math::BigFloat (or Math::BigInt), so make sure you check these two modules for further information. =head2 new() $x = Math::BigRat->new('1/3'); Create a new Math::BigRat object. Input can come in various forms: $x = Math::BigRat->new(123); # scalars $x = Math::BigRat->new('inf'); # infinity $x = Math::BigRat->new('123.3'); # float $x = Math::BigRat->new('1/3'); # simple string $x = Math::BigRat->new('1 / 3'); # spaced $x = Math::BigRat->new('1 / 0.1'); # w/ floats $x = Math::BigRat->new(Math::BigInt->new(3)); # BigInt $x = Math::BigRat->new(Math::BigFloat->new('3.1')); # BigFloat $x = Math::BigRat->new(Math::BigInt::Lite->new('2')); # BigLite =head2 numerator() $n = $x->numerator(); Returns a copy of the numerator (the part above the line) as signed BigInt. =head2 denominator() $d = $x->denominator(); Returns a copy of the denominator (the part under the line) as positive BigInt. =head2 parts() ($n,$d) = $x->parts(); Return a list consisting of (signed) numerator and (unsigned) denominator as BigInts. =head2 as_number() $x = Math::BigRat->new('13/7'); print $x->as_number(),"\n"; # '1' Returns a copy of the object as BigInt trunced it to integer. =head2 bfac() $x->bfac(); Calculates the factorial of $x. For instance: print Math::BigRat->new('3/1')->bfac(),"\n"; # 1*2*3 print Math::BigRat->new('5/1')->bfac(),"\n"; # 1*2*3*4*5 Works currently only for integers. =head2 blog() Is not yet implemented. =head2 bround()/round()/bfround() Are not yet implemented. =head2 bmod() use Math::BigRat; my $x = Math::BigRat->new('7/4'); my $y = Math::BigRat->new('4/3'); print $x->bmod($y); Set $x to the remainder of the division of $x by $y. =head2 is_one() print "$x is 1\n" if $x->is_one(); Return true if $x is exactly one, otherwise false. =head2 is_zero() print "$x is 0\n" if $x->is_zero(); Return true if $x is exactly zero, otherwise false. =head2 is_positive() print "$x is >= 0\n" if $x->is_positive(); Return true if $x is positive (greater than or equal to zero), otherwise false. Please note that '+inf' is also positive, while 'NaN' and '-inf' aren't. =head2 is_negative() print "$x is < 0\n" if $x->is_negative(); Return true if $x is negative (smaller than zero), otherwise false. Please note that '-inf' is also negative, while 'NaN' and '+inf' aren't. =head2 is_int() print "$x is an integer\n" if $x->is_int(); Return true if $x has a denominator of 1 (e.g. no fraction parts), otherwise false. Please note that '-inf', 'inf' and 'NaN' aren't integer. =head2 is_odd() print "$x is odd\n" if $x->is_odd(); Return true if $x is odd, otherwise false. =head2 is_even() print "$x is even\n" if $x->is_even(); Return true if $x is even, otherwise false. =head2 bceil() $x->bceil(); Set $x to the next bigger integer value (e.g. truncate the number to integer and then increment it by one). =head2 bfloor() $x->bfloor(); Truncate $x to an integer value. =head2 bsqrt() $x->bsqrt(); Calculate the square root of $x. =head2 config use Data::Dumper; print Dumper ( Math::BigRat->config() ); print Math::BigRat->config()->{lib},"\n"; Returns a hash containing the configuration, e.g. the version number, lib loaded etc. The following hash keys are currently filled in with the appropriate information. key RO/RW Description Example ============================================================ lib RO Name of the Math library Math::BigInt::Calc lib_version RO Version of 'lib' 0.30 class RO The class of config you just called Math::BigRat version RO version number of the class you used 0.10 upgrade RW To which class numbers are upgraded undef downgrade RW To which class numbers are downgraded undef precision RW Global precision undef accuracy RW Global accuracy undef round_mode RW Global round mode even div_scale RW Fallback acccuracy for div 40 trap_nan RW Trap creation of NaN (undef = no) undef trap_inf RW Trap creation of +inf/-inf (undef = no) undef By passing a reference to a hash you may set the configuration values. This works only for values that a marked with a C above, anything else is read-only. =head1 BUGS Some things are not yet implemented, or only implemented half-way: =over 2 =item inf handling (partial) =item NaN handling (partial) =item rounding (not implemented except for bceil/bfloor) =item $x ** $y where $y is not an integer =item bmod(), blog(), bmodinv() and bmodpow() (partial) =back =head1 LICENSE This program is free software; you may redistribute it and/or modify it under the same terms as Perl itself. =head1 SEE ALSO L and L as well as L, L and L. See L for a way to use Math::BigRat. The package at L may contain more documentation and examples as well as testcases. =head1 AUTHORS (C) by Tels L 2001, 2002, 2003, 2004. =cut