lib/Math/BigRat.pm

   1
   2 #
   3 # "Tax the rat farms." - Lord Vetinari
   4 #
   5
   6 # The following hash values are used:
   7 #   sign : +,-,NaN,+inf,-inf
   8 #   _d   : denominator
   9 #   _n   : numeraotr (value = _n/_d)
  10 #   _a   : accuracy
  11 #   _p   : precision
  12 # You should not look at the innards of a BigRat - use the methods for this.
  13
  14 package Math::BigRat;
  15
  16 # anythig older is untested, and unlikely to work
  17 use 5.006;
  18 use strict;
  19
  20 use Math::BigFloat;
  21 use vars qw($VERSION @ISA $upgrade $downgrade
  22             $accuracy $precision $round_mode $div_scale $_trap_nan $_trap_inf);
  23
  24 @ISA = qw(Math::BigFloat);
  25
  26 $VERSION = '0.21';
  27
  28 use overload;                   # inherit overload from Math::BigFloat
  29
  30 BEGIN
  31   {
  32   *objectify = \&Math::BigInt::objectify;       # inherit this from BigInt
  33   *AUTOLOAD = \&Math::BigFloat::AUTOLOAD;       # can't inherit AUTOLOAD
  34   # we inherit these from BigFloat because currently it is not possible
  35   # that MBF has a different $MBI variable than we, because MBF also uses
  36   # Math::BigInt::config->('lib'); (there is always only one library loaded)
  37   *_e_add = \&Math::BigFloat::_e_add;
  38   *_e_sub = \&Math::BigFloat::_e_sub;
  39   *as_int = \&as_number;
  40   *is_pos = \&is_positive;
  41   *is_neg = \&is_negative;
  42   }
  43
  44 ##############################################################################
  45 # Global constants and flags. Access these only via the accessor methods!
  46
  47 $accuracy = $precision = undef;
  48 $round_mode = 'even';
  49 $div_scale = 40;
  50 $upgrade = undef;
  51 $downgrade = undef;
  52
  53 # These are internally, and not to be used from the outside at all!
  54
  55 $_trap_nan = 0;                         # are NaNs ok? set w/ config()
  56 $_trap_inf = 0;                         # are infs ok? set w/ config()
  57
  58 # the package we are using for our private parts, defaults to:
  59 # Math::BigInt->config()->{lib}
  60 my $MBI = 'Math::BigInt::Calc';
  61
  62 my $nan = 'NaN';
  63 my $class = 'Math::BigRat';
  64
  65 sub isa
  66   {
  67   return 0 if $_[1] =~ /^Math::Big(Int|Float)/;         # we aren't
  68   UNIVERSAL::isa(@_);
  69   }
  70
  71 ##############################################################################
  72
  73 sub _new_from_float
  74   {
  75   # turn a single float input into a rational number (like '0.1')
  76   my ($self,$f) = @_;
  77
  78   return $self->bnan() if $f->is_nan();
  79   return $self->binf($f->{sign}) if $f->{sign} =~ /^[+-]inf$/;
  80
  81   $self->{_n} = $MBI->_copy( $f->{_m} );        # mantissa
  82   $self->{_d} = $MBI->_one();
  83   $self->{sign} = $f->{sign} || '+';
  84   if ($f->{_es} eq '-')
  85     {
  86     # something like Math::BigRat->new('0.1');
  87     # 1 / 1 => 1/10
  88     $MBI->_lsft ( $self->{_d}, $f->{_e} ,10);
  89     }
  90   else
  91     {
  92     # something like Math::BigRat->new('10');
  93     # 1 / 1 => 10/1
  94     $MBI->_lsft ( $self->{_n}, $f->{_e} ,10) unless
  95       $MBI->_is_zero($f->{_e});
  96     }
  97   $self;
  98   }
  99
 100 sub new
 101   {
 102   # create a Math::BigRat
 103   my $class = shift;
 104
 105   my ($n,$d) = @_;
 106
 107   my $self = { }; bless $self,$class;
 108
 109   # input like (BigInt) or (BigFloat):
 110   if ((!defined $d) && (ref $n) && (!$n->isa('Math::BigRat')))
 111     {
 112     if ($n->isa('Math::BigFloat'))
 113       {
 114       $self->_new_from_float($n);
 115       }
 116     if ($n->isa('Math::BigInt'))
 117       {
 118       # TODO: trap NaN, inf
 119       $self->{_n} = $MBI->_copy($n->{value});           # "mantissa" = N
 120       $self->{_d} = $MBI->_one();                       # d => 1
 121       $self->{sign} = $n->{sign};
 122       }
 123     if ($n->isa('Math::BigInt::Lite'))
 124       {
 125       # TODO: trap NaN, inf
 126       $self->{sign} = '+'; $self->{sign} = '-' if $$n < 0;
 127       $self->{_n} = $MBI->_new(abs($$n));               # "mantissa" = N
 128       $self->{_d} = $MBI->_one();                       # d => 1
 129       }
 130     return $self->bnorm();                              # normalize (120/1 => 12/10)
 131     }
 132
 133   # input like (BigInt,BigInt) or (BigLite,BigLite):
 134   if (ref($d) && ref($n))
 135     {
 136     # do N first (for $self->{sign}):
 137     if ($n->isa('Math::BigInt'))
 138       {
 139       # TODO: trap NaN, inf
 140       $self->{_n} = $MBI->_copy($n->{value});           # "mantissa" = N
 141       $self->{sign} = $n->{sign};
 142       }
 143     elsif ($n->isa('Math::BigInt::Lite'))
 144       {
 145       # TODO: trap NaN, inf
 146       $self->{sign} = '+'; $self->{sign} = '-' if $$n < 0;
 147       $self->{_n} = $MBI->_new(abs($$n));               # "mantissa" = $n
 148       }
 149     else
 150       {
 151       require Carp;
 152       Carp::croak(ref($n) . " is not a recognized object format for Math::BigRat->new");
 153       }
 154     # now D:
 155     if ($d->isa('Math::BigInt'))
 156       {
 157       # TODO: trap NaN, inf
 158       $self->{_d} = $MBI->_copy($d->{value});           # "mantissa" = D
 159       # +/+ or -/- => +, +/- or -/+ => -
 160       $self->{sign} = $d->{sign} ne $self->{sign} ? '-' : '+';
 161       }
 162     elsif ($d->isa('Math::BigInt::Lite'))
 163       {
 164       # TODO: trap NaN, inf
 165       $self->{_d} = $MBI->_new(abs($$d));               # "mantissa" = D
 166       my $ds = '+'; $ds = '-' if $$d < 0;
 167       # +/+ or -/- => +, +/- or -/+ => -
 168       $self->{sign} = $ds ne $self->{sign} ? '-' : '+';
 169       }
 170     else
 171       {
 172       require Carp;
 173       Carp::croak(ref($d) . " is not a recognized object format for Math::BigRat->new");
 174       }
 175     return $self->bnorm();                              # normalize (120/1 => 12/10)
 176     }
 177   return $n->copy() if ref $n;                          # already a BigRat
 178
 179   if (!defined $n)
 180     {
 181     $self->{_n} = $MBI->_zero();                        # undef => 0
 182     $self->{_d} = $MBI->_one();
 183     $self->{sign} = '+';
 184     return $self;
 185     }
 186
 187   # string input with / delimiter
 188   if ($n =~ /\s*\/\s*/)
 189     {
 190     return $class->bnan() if $n =~ /\/.*\//;    # 1/2/3 isn't valid
 191     return $class->bnan() if $n =~ /\/\s*$/;    # 1/ isn't valid
 192     ($n,$d) = split (/\//,$n);
 193     # try as BigFloats first
 194     if (($n =~ /[\.eE]/) || ($d =~ /[\.eE]/))
 195       {
 196       local $Math::BigFloat::accuracy = undef;
 197       local $Math::BigFloat::precision = undef;
 198
 199       # one of them looks like a float
 200       my $nf = Math::BigFloat->new($n,undef,undef);
 201       $self->{sign} = '+';
 202       return $self->bnan() if $nf->is_nan();
 203
 204       $self->{_n} = $MBI->_copy( $nf->{_m} );   # get mantissa
 205
 206       # now correct $self->{_n} due to $n
 207       my $f = Math::BigFloat->new($d,undef,undef);
 208       return $self->bnan() if $f->is_nan();
 209       $self->{_d} = $MBI->_copy( $f->{_m} );
 210
 211       # calculate the difference between nE and dE
 212       my $diff_e = $nf->exponent()->bsub( $f->exponent);
 213       if ($diff_e->is_negative())
 214         {
 215         # < 0: mul d with it
 216         $MBI->_lsft( $self->{_d}, $MBI->_new( $diff_e->babs()), 10);
 217         }
 218       elsif (!$diff_e->is_zero())
 219         {
 220         # > 0: mul n with it
 221         $MBI->_lsft( $self->{_n}, $MBI->_new( $diff_e), 10);
 222         }
 223       }
 224     else
 225       {
 226       # both d and n look like (big)ints
 227
 228       $self->{sign} = '+';                                      # no sign => '+'
 229       $self->{_n} = undef;
 230       $self->{_d} = undef;
 231       if ($n =~ /^([+-]?)0*([0-9]+)\z/)                         # first part ok?
 232         {
 233         $self->{sign} = $1 || '+';                              # no sign => '+'
 234         $self->{_n} = $MBI->_new($2 || 0);
 235         }
 236
 237       if ($d =~ /^([+-]?)0*([0-9]+)\z/)                         # second part ok?
 238         {
 239         $self->{sign} =~ tr/+-/-+/ if ($1 || '') eq '-';        # negate if second part neg.
 240         $self->{_d} = $MBI->_new($2 || 0);
 241         }
 242
 243       if (!defined $self->{_n} || !defined $self->{_d})
 244         {
 245         $d = Math::BigInt->new($d,undef,undef) unless ref $d;
 246         $n = Math::BigInt->new($n,undef,undef) unless ref $n;
 247
 248         if ($n->{sign} =~ /^[+-]$/ && $d->{sign} =~ /^[+-]$/)
 249           {
 250           # both parts are ok as integers (wierd things like ' 1e0'
 251           $self->{_n} = $MBI->_copy($n->{value});
 252           $self->{_d} = $MBI->_copy($d->{value});
 253           $self->{sign} = $n->{sign};
 254           $self->{sign} =~ tr/+-/-+/ if $d->{sign} eq '-';      # -1/-2 => 1/2
 255           return $self->bnorm();
 256           }
 257
 258         $self->{sign} = '+';                                    # a default sign
 259         return $self->bnan() if $n->is_nan() || $d->is_nan();
 260
 261         # handle inf cases:
 262         if ($n->is_inf() || $d->is_inf())
 263           {
 264           if ($n->is_inf())
 265             {
 266             return $self->bnan() if $d->is_inf();               # both are inf => NaN
 267             my $s = '+';                # '+inf/+123' or '-inf/-123'
 268             $s = '-' if substr($n->{sign},0,1) ne $d->{sign};
 269             # +-inf/123 => +-inf
 270             return $self->binf($s);
 271             }
 272           # 123/inf => 0
 273           return $self->bzero();
 274           }
 275         }
 276       }
 277
 278     return $self->bnorm();
 279     }
 280
 281   # simple string input
 282   if (($n =~ /[\.eE]/))
 283     {
 284     # looks like a float, quacks like a float, so probably is a float
 285     $self->{sign} = 'NaN';
 286     local $Math::BigFloat::accuracy = undef;
 287     local $Math::BigFloat::precision = undef;
 288     $self->_new_from_float(Math::BigFloat->new($n,undef,undef));
 289     }
 290   else
 291     {
 292     # for simple forms, use $MBI directly
 293     if ($n =~ /^([+-]?)0*([0-9]+)\z/)
 294       {
 295       $self->{sign} = $1 || '+';
 296       $self->{_n} = $MBI->_new($2 || 0);
 297       $self->{_d} = $MBI->_one();
 298       }
 299     else
 300       {
 301       my $n = Math::BigInt->new($n,undef,undef);
 302       $self->{_n} = $MBI->_copy($n->{value});
 303       $self->{_d} = $MBI->_one();
 304       $self->{sign} = $n->{sign};
 305       return $self->bnan() if $self->{sign} eq 'NaN';
 306       return $self->binf($self->{sign}) if $self->{sign} =~ /^[+-]inf$/;
 307       }
 308     }
 309   $self->bnorm();
 310   }
 311
 312 sub copy
 313   {
 314   # if two arguments, the first one is the class to "swallow" subclasses
 315   my ($c,$x) = @_;
 316
 317   if (scalar @_ == 1)
 318     {
 319     $x = $_[0];
 320     $c = ref($x);
 321     }
 322   return unless ref($x); # only for objects
 323
 324   my $self = bless {}, $c;
 325
 326   $self->{sign} = $x->{sign};
 327   $self->{_d} = $MBI->_copy($x->{_d});
 328   $self->{_n} = $MBI->_copy($x->{_n});
 329   $self->{_a} = $x->{_a} if defined $x->{_a};
 330   $self->{_p} = $x->{_p} if defined $x->{_p};
 331   $self;
 332   }
 333
 334 ##############################################################################
 335
 336 sub config
 337   {
 338   # return (later set?) configuration data as hash ref
 339   my $class = shift || 'Math::BigRat';
 340
 341   if (@_ == 1 && ref($_[0]) ne 'HASH')
 342     {
 343     my $cfg = $class->SUPER::config();
 344     return $cfg->{$_[0]};
 345     }
 346
 347   my $cfg = $class->SUPER::config(@_);
 348
 349   # now we need only to override the ones that are different from our parent
 350   $cfg->{class} = $class;
 351   $cfg->{with} = $MBI;
 352   $cfg;
 353   }
 354
 355 ##############################################################################
 356
 357 sub bstr
 358   {
 359   my ($self,$x) = ref($_[0]) ? (undef,$_[0]) : objectify(1,@_);
 360
 361   if ($x->{sign} !~ /^[+-]$/)           # inf, NaN etc
 362     {
 363     my $s = $x->{sign}; $s =~ s/^\+//;  # +inf => inf
 364     return $s;
 365     }
 366
 367   my $s = ''; $s = $x->{sign} if $x->{sign} ne '+';     # '+3/2' => '3/2'
 368
 369   return $s . $MBI->_str($x->{_n}) if $MBI->_is_one($x->{_d});
 370   $s . $MBI->_str($x->{_n}) . '/' . $MBI->_str($x->{_d});
 371   }
 372
 373 sub bsstr
 374   {
 375   my ($self,$x) = ref($_[0]) ? (undef,$_[0]) : objectify(1,@_);
 376
 377   if ($x->{sign} !~ /^[+-]$/)           # inf, NaN etc
 378     {
 379     my $s = $x->{sign}; $s =~ s/^\+//;  # +inf => inf
 380     return $s;
 381     }
 382
 383   my $s = ''; $s = $x->{sign} if $x->{sign} ne '+';     # +3 vs 3
 384   $s . $MBI->_str($x->{_n}) . '/' . $MBI->_str($x->{_d});
 385   }
 386
 387 sub bnorm
 388   {
 389   # reduce the number to the shortest form
 390   my ($self,$x) = ref($_[0]) ? (undef,$_[0]) : objectify(1,@_);
 391
 392   # Both parts must be objects of whatever we are using today.
 393   if ( my $c = $MBI->_check($x->{_n}) )
 394     {
 395     require Carp; Carp::croak ("n did not pass the self-check ($c) in bnorm()");
 396     }
 397   if ( my $c = $MBI->_check($x->{_d}) )
 398     {
 399     require Carp; Carp::croak ("d did not pass the self-check ($c) in bnorm()");
 400     }
 401
 402   # no normalize for NaN, inf etc.
 403   return $x if $x->{sign} !~ /^[+-]$/;
 404
 405   # normalize zeros to 0/1
 406   if ($MBI->_is_zero($x->{_n}))
 407     {
 408     $x->{sign} = '+';                                   # never leave a -0
 409     $x->{_d} = $MBI->_one() unless $MBI->_is_one($x->{_d});
 410     return $x;
 411     }
 412
 413   return $x if $MBI->_is_one($x->{_d});                 # no need to reduce
 414
 415   # reduce other numbers
 416   my $gcd = $MBI->_copy($x->{_n});
 417   $gcd = $MBI->_gcd($gcd,$x->{_d});
 418
 419   if (!$MBI->_is_one($gcd))
 420     {
 421     $x->{_n} = $MBI->_div($x->{_n},$gcd);
 422     $x->{_d} = $MBI->_div($x->{_d},$gcd);
 423     }
 424   $x;
 425   }
 426
 427 ##############################################################################
 428 # sign manipulation
 429
 430 sub bneg
 431   {
 432   # (BRAT or num_str) return BRAT
 433   # negate number or make a negated number from string
 434   my ($self,$x) = ref($_[0]) ? (undef,$_[0]) : objectify(1,@_);
 435
 436   return $x if $x->modify('bneg');
 437
 438   # for +0 dont negate (to have always normalized +0). Does nothing for 'NaN'
 439   $x->{sign} =~ tr/+-/-+/ unless ($x->{sign} eq '+' && $MBI->_is_zero($x->{_n}));
 440   $x;
 441   }
 442
 443 ##############################################################################
 444 # special values
 445
 446 sub _bnan
 447   {
 448   # used by parent class bnan() to initialize number to NaN
 449   my $self = shift;
 450
 451   if ($_trap_nan)
 452     {
 453     require Carp;
 454     my $class = ref($self);
 455     # "$self" below will stringify the object, this blows up if $self is a
 456     # partial object (happens under trap_nan), so fix it beforehand
 457     $self->{_d} = $MBI->_zero() unless defined $self->{_d};
 458     $self->{_n} = $MBI->_zero() unless defined $self->{_n};
 459     Carp::croak ("Tried to set $self to NaN in $class\::_bnan()");
 460     }
 461   $self->{_n} = $MBI->_zero();
 462   $self->{_d} = $MBI->_zero();
 463   }
 464
 465 sub _binf
 466   {
 467   # used by parent class bone() to initialize number to +inf/-inf
 468   my $self = shift;
 469
 470   if ($_trap_inf)
 471     {
 472     require Carp;
 473     my $class = ref($self);
 474     # "$self" below will stringify the object, this blows up if $self is a
 475     # partial object (happens under trap_nan), so fix it beforehand
 476     $self->{_d} = $MBI->_zero() unless defined $self->{_d};
 477     $self->{_n} = $MBI->_zero() unless defined $self->{_n};
 478     Carp::croak ("Tried to set $self to inf in $class\::_binf()");
 479     }
 480   $self->{_n} = $MBI->_zero();
 481   $self->{_d} = $MBI->_zero();
 482   }
 483
 484 sub _bone
 485   {
 486   # used by parent class bone() to initialize number to +1/-1
 487   my $self = shift;
 488   $self->{_n} = $MBI->_one();
 489   $self->{_d} = $MBI->_one();
 490   }
 491
 492 sub _bzero
 493   {
 494   # used by parent class bzero() to initialize number to 0
 495   my $self = shift;
 496   $self->{_n} = $MBI->_zero();
 497   $self->{_d} = $MBI->_one();
 498   }
 499
 500 ##############################################################################
 501 # mul/add/div etc
 502
 503 sub badd
 504   {
 505   # add two rational numbers
 506
 507   # set up parameters
 508   my ($self,$x,$y,@r) = (ref($_[0]),@_);
 509   # objectify is costly, so avoid it
 510   if ((!ref($_[0])) || (ref($_[0]) ne ref($_[1])))
 511     {
 512     ($self,$x,$y,@r) = objectify(2,@_);
 513     }
 514
 515   # +inf + +inf => +inf,  -inf + -inf => -inf
 516   return $x->binf(substr($x->{sign},0,1))
 517     if $x->{sign} eq $y->{sign} && $x->{sign} =~ /^[+-]inf$/;
 518
 519   # +inf + -inf or -inf + +inf => NaN
 520   return $x->bnan() if ($x->{sign} !~ /^[+-]$/ || $y->{sign} !~ /^[+-]$/);
 521
 522   #  1   1    gcd(3,4) = 1    1*3 + 1*4    7
 523   #  - + -                  = --------- = --
 524   #  4   3                      4*3       12
 525
 526   # we do not compute the gcd() here, but simple do:
 527   #  5   7    5*3 + 7*4   43
 528   #  - + -  = --------- = --
 529   #  4   3       4*3      12
 530
 531   # and bnorm() will then take care of the rest
 532
 533   # 5 * 3
 534   $x->{_n} = $MBI->_mul( $x->{_n}, $y->{_d});
 535
 536   # 7 * 4
 537   my $m = $MBI->_mul( $MBI->_copy( $y->{_n} ), $x->{_d} );
 538
 539   # 5 * 3 + 7 * 4
 540   ($x->{_n}, $x->{sign}) = _e_add( $x->{_n}, $m, $x->{sign}, $y->{sign});
 541
 542   # 4 * 3
 543   $x->{_d} = $MBI->_mul( $x->{_d}, $y->{_d});
 544
 545   # normalize result, and possible round
 546   $x->bnorm()->round(@r);
 547   }
 548
 549 sub bsub
 550   {
 551   # subtract two rational numbers
 552
 553   # set up parameters
 554   my ($self,$x,$y,@r) = (ref($_[0]),@_);
 555   # objectify is costly, so avoid it
 556   if ((!ref($_[0])) || (ref($_[0]) ne ref($_[1])))
 557     {
 558     ($self,$x,$y,@r) = objectify(2,@_);
 559     }
 560
 561   # flip sign of $x, call badd(), then flip sign of result
 562   $x->{sign} =~ tr/+-/-+/
 563     unless $x->{sign} eq '+' && $MBI->_is_zero($x->{_n});       # not -0
 564   $x->badd($y,@r);                              # does norm and round
 565   $x->{sign} =~ tr/+-/-+/
 566     unless $x->{sign} eq '+' && $MBI->_is_zero($x->{_n});       # not -0
 567   $x;
 568   }
 569
 570 sub bmul
 571   {
 572   # multiply two rational numbers
 573
 574   # set up parameters
 575   my ($self,$x,$y,@r) = (ref($_[0]),@_);
 576   # objectify is costly, so avoid it
 577   if ((!ref($_[0])) || (ref($_[0]) ne ref($_[1])))
 578     {
 579     ($self,$x,$y,@r) = objectify(2,@_);
 580     }
 581
 582   return $x->bnan() if ($x->{sign} eq 'NaN' || $y->{sign} eq 'NaN');
 583
 584   # inf handling
 585   if (($x->{sign} =~ /^[+-]inf$/) || ($y->{sign} =~ /^[+-]inf$/))
 586     {
 587     return $x->bnan() if $x->is_zero() || $y->is_zero();
 588     # result will always be +-inf:
 589     # +inf * +/+inf => +inf, -inf * -/-inf => +inf
 590     # +inf * -/-inf => -inf, -inf * +/+inf => -inf
 591     return $x->binf() if ($x->{sign} =~ /^\+/ && $y->{sign} =~ /^\+/);
 592     return $x->binf() if ($x->{sign} =~ /^-/ && $y->{sign} =~ /^-/);
 593     return $x->binf('-');
 594     }
 595
 596   # x== 0 # also: or y == 1 or y == -1
 597   return wantarray ? ($x,$self->bzero()) : $x if $x->is_zero();
 598
 599   # XXX TODO:
 600   # According to Knuth, this can be optimized by doing gcd twice (for d and n)
 601   # and reducing in one step. This would save us the bnorm() at the end.
 602
 603   #  1   2    1 * 2    2    1
 604   #  - * - =  -----  = -  = -
 605   #  4   3    4 * 3    12   6
 606
 607   $x->{_n} = $MBI->_mul( $x->{_n}, $y->{_n});
 608   $x->{_d} = $MBI->_mul( $x->{_d}, $y->{_d});
 609
 610   # compute new sign
 611   $x->{sign} = $x->{sign} eq $y->{sign} ? '+' : '-';
 612
 613   $x->bnorm()->round(@r);
 614   }
 615
 616 sub bdiv
 617   {
 618   # (dividend: BRAT or num_str, divisor: BRAT or num_str) return
 619   # (BRAT,BRAT) (quo,rem) or BRAT (only rem)
 620
 621   # set up parameters
 622   my ($self,$x,$y,@r) = (ref($_[0]),@_);
 623   # objectify is costly, so avoid it
 624   if ((!ref($_[0])) || (ref($_[0]) ne ref($_[1])))
 625     {
 626     ($self,$x,$y,@r) = objectify(2,@_);
 627     }
 628
 629   return $self->_div_inf($x,$y)
 630    if (($x->{sign} !~ /^[+-]$/) || ($y->{sign} !~ /^[+-]$/) || $y->is_zero());
 631
 632   # x== 0 # also: or y == 1 or y == -1
 633   return wantarray ? ($x,$self->bzero()) : $x if $x->is_zero();
 634
 635   # XXX TODO: list context, upgrade
 636   # According to Knuth, this can be optimized by doing gcd twice (for d and n)
 637   # and reducing in one step. This would save us the bnorm() at the end.
 638
 639   # 1     1    1   3
 640   # -  /  - == - * -
 641   # 4     3    4   1
 642
 643   $x->{_n} = $MBI->_mul( $x->{_n}, $y->{_d});
 644   $x->{_d} = $MBI->_mul( $x->{_d}, $y->{_n});
 645
 646   # compute new sign
 647   $x->{sign} = $x->{sign} eq $y->{sign} ? '+' : '-';
 648
 649   $x->bnorm()->round(@r);
 650   $x;
 651   }
 652
 653 sub bmod
 654   {
 655   # compute "remainder" (in Perl way) of $x / $y
 656
 657   # set up parameters
 658   my ($self,$x,$y,@r) = (ref($_[0]),@_);
 659   # objectify is costly, so avoid it
 660   if ((!ref($_[0])) || (ref($_[0]) ne ref($_[1])))
 661     {
 662     ($self,$x,$y,@r) = objectify(2,@_);
 663     }
 664
 665   return $self->_div_inf($x,$y)
 666    if (($x->{sign} !~ /^[+-]$/) || ($y->{sign} !~ /^[+-]$/) || $y->is_zero());
 667
 668   return $x if $x->is_zero();           # 0 / 7 = 0, mod 0
 669
 670   # compute $x - $y * floor($x/$y), keeping the sign of $x
 671
 672   # copy x to u, make it positive and then do a normal division ($u/$y)
 673   my $u = bless { sign => '+' }, $self;
 674   $u->{_n} = $MBI->_mul( $MBI->_copy($x->{_n}), $y->{_d} );
 675   $u->{_d} = $MBI->_mul( $MBI->_copy($x->{_d}), $y->{_n} );
 676
 677   # compute floor(u)
 678   if (! $MBI->_is_one($u->{_d}))
 679     {
 680     $u->{_n} = $MBI->_div($u->{_n},$u->{_d});   # 22/7 => 3/1 w/ truncate
 681     # no need to set $u->{_d} to 1, since below we set it to $y->{_d} anyway
 682     }
 683
 684   # now compute $y * $u
 685   $u->{_d} = $MBI->_copy($y->{_d});             # 1 * $y->{_d}, see floor above
 686   $u->{_n} = $MBI->_mul($u->{_n},$y->{_n});
 687
 688   my $xsign = $x->{sign}; $x->{sign} = '+';     # remember sign and make x positive
 689   # compute $x - $u
 690   $x->bsub($u);
 691   $x->{sign} = $xsign;                          # put sign back
 692
 693   $x->bnorm()->round(@r);
 694   }
 695
 696 ##############################################################################
 697 # bdec/binc
 698
 699 sub bdec
 700   {
 701   # decrement value (subtract 1)
 702   my ($self,$x,@r) = ref($_[0]) ? (ref($_[0]),@_) : objectify(1,@_);
 703
 704   return $x if $x->{sign} !~ /^[+-]$/;  # NaN, inf, -inf
 705
 706   if ($x->{sign} eq '-')
 707     {
 708     $x->{_n} = $MBI->_add( $x->{_n}, $x->{_d});         # -5/2 => -7/2
 709     }
 710   else
 711     {
 712     if ($MBI->_acmp($x->{_n},$x->{_d}) < 0)             # n < d?
 713       {
 714       # 1/3 -- => -2/3
 715       $x->{_n} = $MBI->_sub( $MBI->_copy($x->{_d}), $x->{_n});
 716       $x->{sign} = '-';
 717       }
 718     else
 719       {
 720       $x->{_n} = $MBI->_sub($x->{_n}, $x->{_d});        # 5/2 => 3/2
 721       }
 722     }
 723   $x->bnorm()->round(@r);
 724   }
 725
 726 sub binc
 727   {
 728   # increment value (add 1)
 729   my ($self,$x,@r) = ref($_[0]) ? (ref($_[0]),@_) : objectify(1,@_);
 730
 731   return $x if $x->{sign} !~ /^[+-]$/;  # NaN, inf, -inf
 732
 733   if ($x->{sign} eq '-')
 734     {
 735     if ($MBI->_acmp($x->{_n},$x->{_d}) < 0)
 736       {
 737       # -1/3 ++ => 2/3 (overflow at 0)
 738       $x->{_n} = $MBI->_sub( $MBI->_copy($x->{_d}), $x->{_n});
 739       $x->{sign} = '+';
 740       }
 741     else
 742       {
 743       $x->{_n} = $MBI->_sub($x->{_n}, $x->{_d});        # -5/2 => -3/2
 744       }
 745     }
 746   else
 747     {
 748     $x->{_n} = $MBI->_add($x->{_n},$x->{_d});           # 5/2 => 7/2
 749     }
 750   $x->bnorm()->round(@r);
 751   }
 752
 753 ##############################################################################
 754 # is_foo methods (the rest is inherited)
 755
 756 sub is_int
 757   {
 758   # return true if arg (BRAT or num_str) is an integer
 759   my ($self,$x) = ref($_[0]) ? (undef,$_[0]) : objectify(1,@_);
 760
 761   return 1 if ($x->{sign} =~ /^[+-]$/) &&       # NaN and +-inf aren't
 762     $MBI->_is_one($x->{_d});                    # x/y && y != 1 => no integer
 763   0;
 764   }
 765
 766 sub is_zero
 767   {
 768   # return true if arg (BRAT or num_str) is zero
 769   my ($self,$x) = ref($_[0]) ? (undef,$_[0]) : objectify(1,@_);
 770
 771   return 1 if $x->{sign} eq '+' && $MBI->_is_zero($x->{_n});
 772   0;
 773   }
 774
 775 sub is_one
 776   {
 777   # return true if arg (BRAT or num_str) is +1 or -1 if signis given
 778   my ($self,$x) = ref($_[0]) ? (undef,$_[0]) : objectify(1,@_);
 779
 780   my $sign = $_[2] || ''; $sign = '+' if $sign ne '-';
 781   return 1
 782    if ($x->{sign} eq $sign && $MBI->_is_one($x->{_n}) && $MBI->_is_one($x->{_d}));
 783   0;
 784   }
 785
 786 sub is_odd
 787   {
 788   # return true if arg (BFLOAT or num_str) is odd or false if even
 789   my ($self,$x) = ref($_[0]) ? (undef,$_[0]) : objectify(1,@_);
 790
 791   return 1 if ($x->{sign} =~ /^[+-]$/) &&               # NaN & +-inf aren't
 792     ($MBI->_is_one($x->{_d}) && $MBI->_is_odd($x->{_n})); # x/2 is not, but 3/1
 793   0;
 794   }
 795
 796 sub is_even
 797   {
 798   # return true if arg (BINT or num_str) is even or false if odd
 799   my ($self,$x) = ref($_[0]) ? (undef,$_[0]) : objectify(1,@_);
 800
 801   return 0 if $x->{sign} !~ /^[+-]$/;                   # NaN & +-inf aren't
 802   return 1 if ($MBI->_is_one($x->{_d})                  # x/3 is never
 803      && $MBI->_is_even($x->{_n}));                      # but 4/1 is
 804   0;
 805   }
 806
 807 ##############################################################################
 808 # parts() and friends
 809
 810 sub numerator
 811   {
 812   my ($self,$x) = ref($_[0]) ? (ref($_[0]),$_[0]) : objectify(1,@_);
 813
 814   # NaN, inf, -inf
 815   return Math::BigInt->new($x->{sign}) if ($x->{sign} !~ /^[+-]$/);
 816
 817   my $n = Math::BigInt->new($MBI->_str($x->{_n})); $n->{sign} = $x->{sign};
 818   $n;
 819   }
 820
 821 sub denominator
 822   {
 823   my ($self,$x) = ref($_[0]) ? (ref($_[0]),$_[0]) : objectify(1,@_);
 824
 825   # NaN
 826   return Math::BigInt->new($x->{sign}) if $x->{sign} eq 'NaN';
 827   # inf, -inf
 828   return Math::BigInt->bone() if $x->{sign} !~ /^[+-]$/;
 829
 830   Math::BigInt->new($MBI->_str($x->{_d}));
 831   }
 832
 833 sub parts
 834   {
 835   my ($self,$x) = ref($_[0]) ? (ref($_[0]),$_[0]) : objectify(1,@_);
 836
 837   my $c = 'Math::BigInt';
 838
 839   return ($c->bnan(),$c->bnan()) if $x->{sign} eq 'NaN';
 840   return ($c->binf(),$c->binf()) if $x->{sign} eq '+inf';
 841   return ($c->binf('-'),$c->binf()) if $x->{sign} eq '-inf';
 842
 843   my $n = $c->new( $MBI->_str($x->{_n}));
 844   $n->{sign} = $x->{sign};
 845   my $d = $c->new( $MBI->_str($x->{_d}));
 846   ($n,$d);
 847   }
 848
 849 sub length
 850   {
 851   my ($self,$x) = ref($_[0]) ? (undef,$_[0]) : objectify(1,@_);
 852
 853   return $nan unless $x->is_int();
 854   $MBI->_len($x->{_n});                         # length(-123/1) => length(123)
 855   }
 856
 857 sub digit
 858   {
 859   my ($self,$x,$n) = ref($_[0]) ? (undef,$_[0],$_[1]) : objectify(1,@_);
 860
 861   return $nan unless $x->is_int();
 862   $MBI->_digit($x->{_n},$n || 0);               # digit(-123/1,2) => digit(123,2)
 863   }
 864
 865 ##############################################################################
 866 # special calc routines
 867
 868 sub bceil
 869   {
 870   my ($self,$x) = ref($_[0]) ? (ref($_[0]),$_[0]) : objectify(1,@_);
 871
 872   return $x if $x->{sign} !~ /^[+-]$/ ||        # not for NaN, inf
 873             $MBI->_is_one($x->{_d});            # 22/1 => 22, 0/1 => 0
 874
 875   $x->{_n} = $MBI->_div($x->{_n},$x->{_d});     # 22/7 => 3/1 w/ truncate
 876   $x->{_d} = $MBI->_one();                      # d => 1
 877   $x->{_n} = $MBI->_inc($x->{_n})
 878     if $x->{sign} eq '+';                       # +22/7 => 4/1
 879   $x->{sign} = '+' if $MBI->_is_zero($x->{_n}); # -0 => 0
 880   $x;
 881   }
 882
 883 sub bfloor
 884   {
 885   my ($self,$x) = ref($_[0]) ? (ref($_[0]),$_[0]) : objectify(1,@_);
 886
 887   return $x if $x->{sign} !~ /^[+-]$/ ||        # not for NaN, inf
 888             $MBI->_is_one($x->{_d});            # 22/1 => 22, 0/1 => 0
 889
 890   $x->{_n} = $MBI->_div($x->{_n},$x->{_d});     # 22/7 => 3/1 w/ truncate
 891   $x->{_d} = $MBI->_one();                      # d => 1
 892   $x->{_n} = $MBI->_inc($x->{_n})
 893     if $x->{sign} eq '-';                       # -22/7 => -4/1
 894   $x;
 895   }
 896
 897 sub bfac
 898   {
 899   my ($self,$x,@r) = ref($_[0]) ? (ref($_[0]),@_) : objectify(1,@_);
 900
 901   # if $x is not an integer
 902   if (($x->{sign} ne '+') || (!$MBI->_is_one($x->{_d})))
 903     {
 904     return $x->bnan();
 905     }
 906
 907   $x->{_n} = $MBI->_fac($x->{_n});
 908   # since _d is 1, we don't need to reduce/norm the result
 909   $x->round(@r);
 910   }
 911
 912 sub bpow
 913   {
 914   # power ($x ** $y)
 915
 916   # set up parameters
 917   my ($self,$x,$y,@r) = (ref($_[0]),@_);
 918   # objectify is costly, so avoid it
 919   if ((!ref($_[0])) || (ref($_[0]) ne ref($_[1])))
 920     {
 921     ($self,$x,$y,@r) = objectify(2,@_);
 922     }
 923
 924   return $x if $x->{sign} =~ /^[+-]inf$/;       # -inf/+inf ** x
 925   return $x->bnan() if $x->{sign} eq $nan || $y->{sign} eq $nan;
 926   return $x->bone(@r) if $y->is_zero();
 927   return $x->round(@r) if $x->is_one() || $y->is_one();
 928
 929   if ($x->{sign} eq '-' && $MBI->_is_one($x->{_n}) && $MBI->_is_one($x->{_d}))
 930     {
 931     # if $x == -1 and odd/even y => +1/-1
 932     return $y->is_odd() ? $x->round(@r) : $x->babs()->round(@r);
 933     # my Casio FX-5500L has a bug here: -1 ** 2 is -1, but -1 * -1 is 1;
 934     }
 935   # 1 ** -y => 1 / (1 ** |y|)
 936   # so do test for negative $y after above's clause
 937
 938   return $x->round(@r) if $x->is_zero();  # 0**y => 0 (if not y <= 0)
 939
 940   # shortcut y/1 (and/or x/1)
 941   if ($MBI->_is_one($y->{_d}))
 942     {
 943     # shortcut for x/1 and y/1
 944     if ($MBI->_is_one($x->{_d}))
 945       {
 946       $x->{_n} = $MBI->_pow($x->{_n},$y->{_n});         # x/1 ** y/1 => (x ** y)/1
 947       if ($y->{sign} eq '-')
 948         {
 949         # 0.2 ** -3 => 1/(0.2 ** 3)
 950         ($x->{_n},$x->{_d}) = ($x->{_d},$x->{_n});      # swap
 951         }
 952       # correct sign; + ** + => +
 953       if ($x->{sign} eq '-')
 954         {
 955         # - * - => +, - * - * - => -
 956         $x->{sign} = '+' if $MBI->_is_even($y->{_n});
 957         }
 958       return $x->round(@r);
 959       }
 960     # x/z ** y/1
 961     $x->{_n} = $MBI->_pow($x->{_n},$y->{_n});           # 5/2 ** y/1 => 5 ** y / 2 ** y
 962     $x->{_d} = $MBI->_pow($x->{_d},$y->{_n});
 963     if ($y->{sign} eq '-')
 964       {
 965       # 0.2 ** -3 => 1/(0.2 ** 3)
 966       ($x->{_n},$x->{_d}) = ($x->{_d},$x->{_n});        # swap
 967       }
 968     # correct sign; + ** + => +
 969     if ($x->{sign} eq '-')
 970       {
 971       # - * - => +, - * - * - => -
 972       $x->{sign} = '+' if $MBI->_is_even($y->{_n});
 973       }
 974     return $x->round(@r);
 975     }
 976
 977   # regular calculation (this is wrong for d/e ** f/g)
 978   my $pow2 = $self->bone();
 979   my $y1 = $MBI->_div ( $MBI->_copy($y->{_n}), $y->{_d});
 980   my $two = $MBI->_two();
 981
 982   while (!$MBI->_is_one($y1))
 983     {
 984     $pow2->bmul($x) if $MBI->_is_odd($y1);
 985     $MBI->_div($y1, $two);
 986     $x->bmul($x);
 987     }
 988   $x->bmul($pow2) unless $pow2->is_one();
 989   # n ** -x => 1/n ** x
 990   ($x->{_d},$x->{_n}) = ($x->{_n},$x->{_d}) if $y->{sign} eq '-';
 991   $x->bnorm()->round(@r);
 992   }
 993
 994 sub blog
 995   {
 996   # set up parameters
 997   my ($self,$x,$y,@r) = (ref($_[0]),@_);
 998
 999   # objectify is costly, so avoid it
1000   if ((!ref($_[0])) || (ref($_[0]) ne ref($_[1])))
1001     {
1002     ($self,$x,$y,@r) = objectify(2,$class,@_);
1003     }
1004
1005   # blog(1,Y) => 0
1006   return $x->bzero() if $x->is_one() && $y->{sign} eq '+';
1007
1008   # $x <= 0 => NaN
1009   return $x->bnan() if $x->is_zero() || $x->{sign} ne '+' || $y->{sign} ne '+';
1010
1011   if ($x->is_int() && $y->is_int())
1012     {
1013     return $self->new($x->as_number()->blog($y->as_number(),@r));
1014     }
1015
1016   # do it with floats
1017   $x->_new_from_float( $x->_as_float()->blog(Math::BigFloat->new("$y"),@r) );
1018   }
1019
1020 sub bexp
1021   {
1022   # set up parameters
1023   my ($self,$x,$y,$a,$p,$r) = (ref($_[0]),@_);
1024
1025   # objectify is costly, so avoid it
1026   if ((!ref($_[0])) || (ref($_[0]) ne ref($_[1])))
1027     {
1028     ($self,$x,$y,$a,$p,$r) = objectify(2,$class,@_);
1029     }
1030
1031   return $x->binf() if $x->{sign} eq '+inf';
1032   return $x->bzero() if $x->{sign} eq '-inf';
1033
1034   # we need to limit the accuracy to protect against overflow
1035   my $fallback = 0;
1036   my ($scale,@params);
1037   ($x,@params) = $x->_find_round_parameters($a,$p,$r);
1038
1039   # also takes care of the "error in _find_round_parameters?" case
1040   return $x if $x->{sign} eq 'NaN';
1041
1042   # no rounding at all, so must use fallback
1043   if (scalar @params == 0)
1044     {
1045     # simulate old behaviour
1046     $params[0] = $self->div_scale();    # and round to it as accuracy
1047     $params[1] = undef;                 # P = undef
1048     $scale = $params[0]+4;              # at least four more for proper round
1049     $params[2] = $r;                    # round mode by caller or undef
1050     $fallback = 1;                      # to clear a/p afterwards
1051     }
1052   else
1053     {
1054     # the 4 below is empirical, and there might be cases where it's not enough...
1055     $scale = abs($params[0] || $params[1]) + 4; # take whatever is defined
1056     }
1057
1058   return $x->bone(@params) if $x->is_zero();
1059
1060   # See the comments in Math::BigFloat on how this algorithm works.
1061   # Basically we calculate A and B (where B is faculty(N)) so that A/B = e
1062
1063   my $x_org = $x->copy();
1064   if ($scale <= 75)
1065     {
1066     # set $x directly from a cached string form
1067     $x->{_n} = $MBI->_new("90933395208605785401971970164779391644753259799242");
1068     $x->{_d} = $MBI->_new("33452526613163807108170062053440751665152000000000");
1069     $x->{sign} = '+';
1070     }
1071   else
1072     {
1073     # compute A and B so that e = A / B.
1074
1075     # After some terms we end up with this, so we use it as a starting point:
1076     my $A = $MBI->_new("90933395208605785401971970164779391644753259799242");
1077     my $F = $MBI->_new(42); my $step = 42;
1078
1079     # Compute how many steps we need to take to get $A and $B sufficiently big
1080     my $steps = Math::BigFloat::_len_to_steps($scale - 4);
1081 #    print STDERR "# Doing $steps steps for ", $scale-4, " digits\n";
1082     while ($step++ <= $steps)
1083       {
1084       # calculate $a * $f + 1
1085       $A = $MBI->_mul($A, $F);
1086       $A = $MBI->_inc($A);
1087       # increment f
1088       $F = $MBI->_inc($F);
1089       }
1090     # compute $B as factorial of $steps (this is faster than doing it manually)
1091     my $B = $MBI->_fac($MBI->_new($steps));
1092
1093 #  print "A ", $MBI->_str($A), "\nB ", $MBI->_str($B), "\n";
1094
1095     $x->{_n} = $A;
1096     $x->{_d} = $B;
1097     $x->{sign} = '+';
1098     }
1099
1100   # $x contains now an estimate of e, with some surplus digits, so we can round
1101   if (!$x_org->is_one())
1102     {
1103     # raise $x to the wanted power and round it in one step:
1104     $x->bpow($x_org, @params);
1105     }
1106   else
1107     {
1108     # else just round the already computed result
1109     delete $x->{_a}; delete $x->{_p};
1110     # shortcut to not run through _find_round_parameters again
1111     if (defined $params[0])
1112       {
1113       $x->bround($params[0],$params[2]);                # then round accordingly
1114       }
1115     else
1116       {
1117       $x->bfround($params[1],$params[2]);               # then round accordingly
1118       }
1119     }
1120   if ($fallback)
1121     {
1122     # clear a/p after round, since user did not request it
1123     delete $x->{_a}; delete $x->{_p};
1124     }
1125
1126   $x;
1127   }
1128
1129 sub bnok
1130   {
1131   # set up parameters
1132   my ($self,$x,$y,@r) = (ref($_[0]),@_);
1133
1134   # objectify is costly, so avoid it
1135   if ((!ref($_[0])) || (ref($_[0]) ne ref($_[1])))
1136     {
1137     ($self,$x,$y,@r) = objectify(2,$class,@_);
1138     }
1139
1140   # do it with floats
1141   $x->_new_from_float( $x->_as_float()->bnok(Math::BigFloat->new("$y"),@r) );
1142   }
1143
1144 sub _float_from_part
1145   {
1146   my $x = shift;
1147
1148   my $f = Math::BigFloat->bzero();
1149   $f->{_m} = $MBI->_copy($x);
1150   $f->{_e} = $MBI->_zero();
1151
1152   $f;
1153   }
1154
1155 sub _as_float
1156   {
1157   my $x = shift;
1158
1159   local $Math::BigFloat::upgrade = undef;
1160   local $Math::BigFloat::accuracy = undef;
1161   local $Math::BigFloat::precision = undef;
1162   # 22/7 => 3.142857143..
1163
1164   my $a = $x->accuracy() || 0;
1165   if ($a != 0 || !$MBI->_is_one($x->{_d}))
1166     {
1167     # n/d
1168     return Math::BigFloat->new($x->{sign} . $MBI->_str($x->{_n}))->bdiv( $MBI->_str($x->{_d}), $x->accuracy());
1169     }
1170   # just n
1171   Math::BigFloat->new($x->{sign} . $MBI->_str($x->{_n}));
1172   }
1173
1174 sub broot
1175   {
1176   # set up parameters
1177   my ($self,$x,$y,@r) = (ref($_[0]),@_);
1178   # objectify is costly, so avoid it
1179   if ((!ref($_[0])) || (ref($_[0]) ne ref($_[1])))
1180     {
1181     ($self,$x,$y,@r) = objectify(2,@_);
1182     }
1183
1184   if ($x->is_int() && $y->is_int())
1185     {
1186     return $self->new($x->as_number()->broot($y->as_number(),@r));
1187     }
1188
1189   # do it with floats
1190   $x->_new_from_float( $x->_as_float()->broot($y,@r) );
1191   }
1192
1193 sub bmodpow
1194   {
1195   # set up parameters
1196   my ($self,$x,$y,$m,@r) = (ref($_[0]),@_);
1197   # objectify is costly, so avoid it
1198   if ((!ref($_[0])) || (ref($_[0]) ne ref($_[1])))
1199     {
1200     ($self,$x,$y,$m,@r) = objectify(3,@_);
1201     }
1202
1203   # $x or $y or $m are NaN or +-inf => NaN
1204   return $x->bnan()
1205    if $x->{sign} !~ /^[+-]$/ || $y->{sign} !~ /^[+-]$/ ||
1206    $m->{sign} !~ /^[+-]$/;
1207
1208   if ($x->is_int() && $y->is_int() && $m->is_int())
1209     {
1210     return $self->new($x->as_number()->bmodpow($y->as_number(),$m,@r));
1211     }
1212
1213   warn ("bmodpow() not fully implemented");
1214   $x->bnan();
1215   }
1216
1217 sub bmodinv
1218   {
1219   # set up parameters
1220   my ($self,$x,$y,@r) = (ref($_[0]),@_);
1221   # objectify is costly, so avoid it
1222   if ((!ref($_[0])) || (ref($_[0]) ne ref($_[1])))
1223     {
1224     ($self,$x,$y,@r) = objectify(2,@_);
1225     }
1226
1227   # $x or $y are NaN or +-inf => NaN
1228   return $x->bnan()
1229    if $x->{sign} !~ /^[+-]$/ || $y->{sign} !~ /^[+-]$/;
1230
1231   if ($x->is_int() && $y->is_int())
1232     {
1233     return $self->new($x->as_number()->bmodinv($y->as_number(),@r));
1234     }
1235
1236   warn ("bmodinv() not fully implemented");
1237   $x->bnan();
1238   }
1239
1240 sub bsqrt
1241   {
1242   my ($self,$x,@r) = ref($_[0]) ? (ref($_[0]),@_) : objectify(1,@_);
1243
1244   return $x->bnan() if $x->{sign} !~ /^[+]/;    # NaN, -inf or < 0
1245   return $x if $x->{sign} eq '+inf';            # sqrt(inf) == inf
1246   return $x->round(@r) if $x->is_zero() || $x->is_one();
1247
1248   local $Math::BigFloat::upgrade = undef;
1249   local $Math::BigFloat::downgrade = undef;
1250   local $Math::BigFloat::precision = undef;
1251   local $Math::BigFloat::accuracy = undef;
1252   local $Math::BigInt::upgrade = undef;
1253   local $Math::BigInt::precision = undef;
1254   local $Math::BigInt::accuracy = undef;
1255
1256   $x->{_n} = _float_from_part( $x->{_n} )->bsqrt();
1257   $x->{_d} = _float_from_part( $x->{_d} )->bsqrt();
1258
1259   # XXX TODO: we probably can optimze this:
1260
1261   # if sqrt(D) was not integer
1262   if ($x->{_d}->{_es} ne '+')
1263     {
1264     $x->{_n}->blsft($x->{_d}->exponent()->babs(),10);   # 7.1/4.51 => 7.1/45.1
1265     $x->{_d} = $MBI->_copy( $x->{_d}->{_m} );           # 7.1/45.1 => 71/45.1
1266     }
1267   # if sqrt(N) was not integer
1268   if ($x->{_n}->{_es} ne '+')
1269     {
1270     $x->{_d}->blsft($x->{_n}->exponent()->babs(),10);   # 71/45.1 => 710/45.1
1271     $x->{_n} = $MBI->_copy( $x->{_n}->{_m} );           # 710/45.1 => 710/451
1272     }
1273
1274   # convert parts to $MBI again
1275   $x->{_n} = $MBI->_lsft( $MBI->_copy( $x->{_n}->{_m} ), $x->{_n}->{_e}, 10)
1276     if ref($x->{_n}) ne $MBI && ref($x->{_n}) ne 'ARRAY';
1277   $x->{_d} = $MBI->_lsft( $MBI->_copy( $x->{_d}->{_m} ), $x->{_d}->{_e}, 10)
1278     if ref($x->{_d}) ne $MBI && ref($x->{_d}) ne 'ARRAY';
1279
1280   $x->bnorm()->round(@r);
1281   }
1282
1283 sub blsft
1284   {
1285   my ($self,$x,$y,$b,@r) = objectify(3,@_);
1286
1287   $b = 2 unless defined $b;
1288   $b = $self->new($b) unless ref ($b);
1289   $x->bmul( $b->copy()->bpow($y), @r);
1290   $x;
1291   }
1292
1293 sub brsft
1294   {
1295   my ($self,$x,$y,$b,@r) = objectify(3,@_);
1296
1297   $b = 2 unless defined $b;
1298   $b = $self->new($b) unless ref ($b);
1299   $x->bdiv( $b->copy()->bpow($y), @r);
1300   $x;
1301   }
1302
1303 ##############################################################################
1304 # round
1305
1306 sub round
1307   {
1308   $_[0];
1309   }
1310
1311 sub bround
1312   {
1313   $_[0];
1314   }
1315
1316 sub bfround
1317   {
1318   $_[0];
1319   }
1320
1321 ##############################################################################
1322 # comparing
1323
1324 sub bcmp
1325   {
1326   # compare two signed numbers
1327
1328   # set up parameters
1329   my ($self,$x,$y) = (ref($_[0]),@_);
1330   # objectify is costly, so avoid it
1331   if ((!ref($_[0])) || (ref($_[0]) ne ref($_[1])))
1332     {
1333     ($self,$x,$y) = objectify(2,@_);
1334     }
1335
1336   if (($x->{sign} !~ /^[+-]$/) || ($y->{sign} !~ /^[+-]$/))
1337     {
1338     # handle +-inf and NaN
1339     return undef if (($x->{sign} eq $nan) || ($y->{sign} eq $nan));
1340     return 0 if $x->{sign} eq $y->{sign} && $x->{sign} =~ /^[+-]inf$/;
1341     return +1 if $x->{sign} eq '+inf';
1342     return -1 if $x->{sign} eq '-inf';
1343     return -1 if $y->{sign} eq '+inf';
1344     return +1;
1345     }
1346   # check sign for speed first
1347   return 1 if $x->{sign} eq '+' && $y->{sign} eq '-';   # does also 0 <=> -y
1348   return -1 if $x->{sign} eq '-' && $y->{sign} eq '+';  # does also -x <=> 0
1349
1350   # shortcut
1351   my $xz = $MBI->_is_zero($x->{_n});
1352   my $yz = $MBI->_is_zero($y->{_n});
1353   return 0 if $xz && $yz;                               # 0 <=> 0
1354   return -1 if $xz && $y->{sign} eq '+';                # 0 <=> +y
1355   return 1 if $yz && $x->{sign} eq '+';                 # +x <=> 0
1356
1357   my $t = $MBI->_mul( $MBI->_copy($x->{_n}), $y->{_d});
1358   my $u = $MBI->_mul( $MBI->_copy($y->{_n}), $x->{_d});
1359
1360   my $cmp = $MBI->_acmp($t,$u);                         # signs are equal
1361   $cmp = -$cmp if $x->{sign} eq '-';                    # both are '-' => reverse
1362   $cmp;
1363   }
1364
1365 sub bacmp
1366   {
1367   # compare two numbers (as unsigned)
1368
1369   # set up parameters
1370   my ($self,$x,$y) = (ref($_[0]),@_);
1371   # objectify is costly, so avoid it
1372   if ((!ref($_[0])) || (ref($_[0]) ne ref($_[1])))
1373     {
1374     ($self,$x,$y) = objectify(2,$class,@_);
1375     }
1376
1377   if (($x->{sign} !~ /^[+-]$/) || ($y->{sign} !~ /^[+-]$/))
1378     {
1379     # handle +-inf and NaN
1380     return undef if (($x->{sign} eq $nan) || ($y->{sign} eq $nan));
1381     return 0 if $x->{sign} =~ /^[+-]inf$/ && $y->{sign} =~ /^[+-]inf$/;
1382     return 1 if $x->{sign} =~ /^[+-]inf$/ && $y->{sign} !~ /^[+-]inf$/;
1383     return -1;
1384     }
1385
1386   my $t = $MBI->_mul( $MBI->_copy($x->{_n}), $y->{_d});
1387   my $u = $MBI->_mul( $MBI->_copy($y->{_n}), $x->{_d});
1388   $MBI->_acmp($t,$u);                                   # ignore signs
1389   }
1390
1391 ##############################################################################
1392 # output conversation
1393
1394 sub numify
1395   {
1396   # convert 17/8 => float (aka 2.125)
1397   my ($self,$x) = ref($_[0]) ? (undef,$_[0]) : objectify(1,@_);
1398
1399   return $x->bstr() if $x->{sign} !~ /^[+-]$/;  # inf, NaN, etc
1400
1401   # N/1 => N
1402   my $neg = ''; $neg = '-' if $x->{sign} eq '-';
1403   return $neg . $MBI->_num($x->{_n}) if $MBI->_is_one($x->{_d});
1404
1405   $x->_as_float()->numify() + 0.0;
1406   }
1407
1408 sub as_number
1409   {
1410   my ($self,$x) = ref($_[0]) ? (undef,$_[0]) : objectify(1,@_);
1411
1412   # NaN, inf etc
1413   return Math::BigInt->new($x->{sign}) if $x->{sign} !~ /^[+-]$/;
1414
1415   my $u = Math::BigInt->bzero();
1416   $u->{sign} = $x->{sign};
1417   $u->{value} = $MBI->_div( $MBI->_copy($x->{_n}), $x->{_d});   # 22/7 => 3
1418   $u;
1419   }
1420
1421 sub as_bin
1422   {
1423   my ($self,$x) = ref($_[0]) ? (undef,$_[0]) : objectify(1,@_);
1424
1425   return $x unless $x->is_int();
1426
1427   my $s = $x->{sign}; $s = '' if $s eq '+';
1428   $s . $MBI->_as_bin($x->{_n});
1429   }
1430
1431 sub as_hex
1432   {
1433   my ($self,$x) = ref($_[0]) ? (undef,$_[0]) : objectify(1,@_);
1434
1435   return $x unless $x->is_int();
1436
1437   my $s = $x->{sign}; $s = '' if $s eq '+';
1438   $s . $MBI->_as_hex($x->{_n});
1439   }
1440
1441 sub as_oct
1442   {
1443   my ($self,$x) = ref($_[0]) ? (undef,$_[0]) : objectify(1,@_);
1444
1445   return $x unless $x->is_int();
1446
1447   my $s = $x->{sign}; $s = '' if $s eq '+';
1448   $s . $MBI->_as_oct($x->{_n});
1449   }
1450
1451 ##############################################################################
1452
1453 sub from_hex
1454   {
1455   my $class = shift;
1456
1457   $class->new(@_);
1458   }
1459
1460 sub from_bin
1461   {
1462   my $class = shift;
1463
1464   $class->new(@_);
1465   }
1466
1467 sub from_oct
1468   {
1469   my $class = shift;
1470
1471   my @parts;
1472   for my $c (@_)
1473     {
1474     push @parts, Math::BigInt->from_oct($c);
1475     }
1476   $class->new ( @parts );
1477   }
1478
1479 ##############################################################################
1480 # import
1481
1482 sub import
1483   {
1484   my $self = shift;
1485   my $l = scalar @_;
1486   my $lib = ''; my @a;
1487   my $try = 'try';
1488
1489   for ( my $i = 0; $i < $l ; $i++)
1490     {
1491     if ( $_[$i] eq ':constant' )
1492       {
1493       # this rest causes overlord er load to step in
1494       overload::constant float => sub { $self->new(shift); };
1495       }
1496 #    elsif ($_[$i] eq 'upgrade')
1497 #      {
1498 #     # this causes upgrading
1499 #      $upgrade = $_[$i+1];             # or undef to disable
1500 #      $i++;
1501 #      }
1502     elsif ($_[$i] eq 'downgrade')
1503       {
1504       # this causes downgrading
1505       $downgrade = $_[$i+1];            # or undef to disable
1506       $i++;
1507       }
1508     elsif ($_[$i] =~ /^(lib|try|only)\z/)
1509       {
1510       $lib = $_[$i+1] || '';            # default Calc
1511       $try = $1;                        # lib, try or only
1512       $i++;
1513       }
1514     elsif ($_[$i] eq 'with')
1515       {
1516       # this argument is no longer used
1517       #$MBI = $_[$i+1] || 'Math::BigInt::Calc'; # default Math::BigInt::Calc
1518       $i++;
1519       }
1520     else
1521       {
1522       push @a, $_[$i];
1523       }
1524     }
1525   require Math::BigInt;
1526
1527   # let use Math::BigInt lib => 'GMP'; use Math::BigRat; still have GMP
1528   if ($lib ne '')
1529     {
1530     my @c = split /\s*,\s*/, $lib;
1531     foreach (@c)
1532       {
1533       $_ =~ tr/a-zA-Z0-9://cd;                    # limit to sane characters
1534       }
1535     $lib = join(",", @c);
1536     }
1537   my @import = ('objectify');
1538   push @import, $try => $lib if $lib ne '';
1539
1540   # MBI already loaded, so feed it our lib arguments
1541   Math::BigInt->import( @import );
1542
1543   $MBI = Math::BigFloat->config()->{lib};
1544
1545   # register us with MBI to get notified of future lib changes
1546   Math::BigInt::_register_callback( $self, sub { $MBI = $_[0]; } );
1547
1548   # any non :constant stuff is handled by our parent, Exporter (loaded
1549   # by Math::BigFloat, even if @_ is empty, to give it a chance
1550   $self->SUPER::import(@a);             # for subclasses
1551   $self->export_to_level(1,$self,@a);   # need this, too
1552   }
1553
1554 1;
1555
1556 __END__
1557
1558 =head1 NAME
1559
1560 Math::BigRat - Arbitrary big rational numbers
1561
1562 =head1 SYNOPSIS
1563
1564         use Math::BigRat;
1565
1566         my $x = Math::BigRat->new('3/7'); $x += '5/9';
1567
1568         print $x->bstr(),"\n";
1569         print $x ** 2,"\n";
1570
1571         my $y = Math::BigRat->new('inf');
1572         print "$y ", ($y->is_inf ? 'is' : 'is not') , " infinity\n";
1573
1574         my $z = Math::BigRat->new(144); $z->bsqrt();
1575
1576 =head1 DESCRIPTION
1577
1578 Math::BigRat complements Math::BigInt and Math::BigFloat by providing support
1579 for arbitrary big rational numbers.
1580
1581 =head2 MATH LIBRARY
1582
1583 You can change the underlying module that does the low-level
1584 math operations by using:
1585
1586         use Math::BigRat try => 'GMP';
1587
1588 Note: This needs Math::BigInt::GMP installed.
1589
1590 The following would first try to find Math::BigInt::Foo, then
1591 Math::BigInt::Bar, and when this also fails, revert to Math::BigInt::Calc:
1592
1593         use Math::BigRat try => 'Foo,Math::BigInt::Bar';
1594
1595 If you want to get warned when the fallback occurs, replace "try" with
1596 "lib":
1597
1598         use Math::BigRat lib => 'Foo,Math::BigInt::Bar';
1599
1600 If you want the code to die instead, replace "try" with
1601 "only":
1602
1603         use Math::BigRat only => 'Foo,Math::BigInt::Bar';
1604
1605 =head1 METHODS
1606
1607 Any methods not listed here are derived from Math::BigFloat (or
1608 Math::BigInt), so make sure you check these two modules for further
1609 information.
1610
1611 =head2 new()
1612
1613         $x = Math::BigRat->new('1/3');
1614
1615 Create a new Math::BigRat object. Input can come in various forms:
1616
1617         $x = Math::BigRat->new(123);                            # scalars
1618         $x = Math::BigRat->new('inf');                          # infinity
1619         $x = Math::BigRat->new('123.3');                        # float
1620         $x = Math::BigRat->new('1/3');                          # simple string
1621         $x = Math::BigRat->new('1 / 3');                        # spaced
1622         $x = Math::BigRat->new('1 / 0.1');                      # w/ floats
1623         $x = Math::BigRat->new(Math::BigInt->new(3));           # BigInt
1624         $x = Math::BigRat->new(Math::BigFloat->new('3.1'));     # BigFloat
1625         $x = Math::BigRat->new(Math::BigInt::Lite->new('2'));   # BigLite
1626
1627         # You can also give D and N as different objects:
1628         $x = Math::BigRat->new(
1629                 Math::BigInt->new(-123),
1630                 Math::BigInt->new(7),
1631                 );                      # => -123/7
1632
1633 =head2 numerator()
1634
1635         $n = $x->numerator();
1636
1637 Returns a copy of the numerator (the part above the line) as signed BigInt.
1638
1639 =head2 denominator()
1640
1641         $d = $x->denominator();
1642
1643 Returns a copy of the denominator (the part under the line) as positive BigInt.
1644
1645 =head2 parts()
1646
1647         ($n,$d) = $x->parts();
1648
1649 Return a list consisting of (signed) numerator and (unsigned) denominator as
1650 BigInts.
1651
1652 =head2 numify()
1653
1654         my $y = $x->numify();
1655
1656 Returns the object as a scalar. This will lose some data if the object
1657 cannot be represented by a normal Perl scalar (integer or float), so
1658 use as_int() instead.
1659
1660 This routine is automatically used whenever a scalar is required:
1661
1662         my $x = Math::BigRat->new('3/1');
1663         @array = (1,2,3);
1664         $y = $array[$x];                # set $y to 3
1665
1666 =head2 as_int()/as_number()
1667
1668         $x = Math::BigRat->new('13/7');
1669         print $x->as_int(),"\n";                # '1'
1670
1671 Returns a copy of the object as BigInt, truncated to an integer.
1672
1673 C<as_number()> is an alias for C<as_int()>.
1674
1675 =head2 as_hex()
1676
1677         $x = Math::BigRat->new('13');
1678         print $x->as_hex(),"\n";                # '0xd'
1679
1680 Returns the BigRat as hexadecimal string. Works only for integers.
1681
1682 =head2 as_bin()
1683
1684         $x = Math::BigRat->new('13');
1685         print $x->as_bin(),"\n";                # '0x1101'
1686
1687 Returns the BigRat as binary string. Works only for integers.
1688
1689 =head2 as_oct()
1690
1691         $x = Math::BigRat->new('13');
1692         print $x->as_oct(),"\n";                # '015'
1693
1694 Returns the BigRat as octal string. Works only for integers.
1695
1696 =head2 from_hex()/from_bin()/from_oct()
1697
1698         my $h = Math::BigRat->from_hex('0x10');
1699         my $b = Math::BigRat->from_bin('0b10000000');
1700         my $o = Math::BigRat->from_oct('020');
1701
1702 Create a BigRat from an hexadecimal, binary or octal number
1703 in string form.
1704
1705 =head2 length()
1706
1707         $len = $x->length();
1708
1709 Return the length of $x in digitis for integer values.
1710
1711 =head2 digit()
1712
1713         print Math::BigRat->new('123/1')->digit(1);     # 1
1714         print Math::BigRat->new('123/1')->digit(-1);    # 3
1715
1716 Return the N'ths digit from X when X is an integer value.
1717
1718 =head2 bnorm()
1719
1720         $x->bnorm();
1721
1722 Reduce the number to the shortest form. This routine is called
1723 automatically whenever it is needed.
1724
1725 =head2 bfac()
1726
1727         $x->bfac();
1728
1729 Calculates the factorial of $x. For instance:
1730
1731         print Math::BigRat->new('3/1')->bfac(),"\n";    # 1*2*3
1732         print Math::BigRat->new('5/1')->bfac(),"\n";    # 1*2*3*4*5
1733
1734 Works currently only for integers.
1735
1736 =head2 bround()/round()/bfround()
1737
1738 Are not yet implemented.
1739
1740 =head2 bmod()
1741
1742         use Math::BigRat;
1743         my $x = Math::BigRat->new('7/4');
1744         my $y = Math::BigRat->new('4/3');
1745         print $x->bmod($y);
1746
1747 Set $x to the remainder of the division of $x by $y.
1748
1749 =head2 bneg()
1750
1751         $x->bneg();
1752
1753 Used to negate the object in-place.
1754
1755 =head2 is_one()
1756
1757         print "$x is 1\n" if $x->is_one();
1758
1759 Return true if $x is exactly one, otherwise false.
1760
1761 =head2 is_zero()
1762
1763         print "$x is 0\n" if $x->is_zero();
1764
1765 Return true if $x is exactly zero, otherwise false.
1766
1767 =head2 is_pos()/is_positive()
1768
1769         print "$x is >= 0\n" if $x->is_positive();
1770
1771 Return true if $x is positive (greater than or equal to zero), otherwise
1772 false. Please note that '+inf' is also positive, while 'NaN' and '-inf' aren't.
1773
1774 C<is_positive()> is an alias for C<is_pos()>.
1775
1776 =head2 is_neg()/is_negative()
1777
1778         print "$x is < 0\n" if $x->is_negative();
1779
1780 Return true if $x is negative (smaller than zero), otherwise false. Please
1781 note that '-inf' is also negative, while 'NaN' and '+inf' aren't.
1782
1783 C<is_negative()> is an alias for C<is_neg()>.
1784
1785 =head2 is_int()
1786
1787         print "$x is an integer\n" if $x->is_int();
1788
1789 Return true if $x has a denominator of 1 (e.g. no fraction parts), otherwise
1790 false. Please note that '-inf', 'inf' and 'NaN' aren't integer.
1791
1792 =head2 is_odd()
1793
1794         print "$x is odd\n" if $x->is_odd();
1795
1796 Return true if $x is odd, otherwise false.
1797
1798 =head2 is_even()
1799
1800         print "$x is even\n" if $x->is_even();
1801
1802 Return true if $x is even, otherwise false.
1803
1804 =head2 bceil()
1805
1806         $x->bceil();
1807
1808 Set $x to the next bigger integer value (e.g. truncate the number to integer
1809 and then increment it by one).
1810
1811 =head2 bfloor()
1812
1813         $x->bfloor();
1814
1815 Truncate $x to an integer value.
1816
1817 =head2 bsqrt()
1818
1819         $x->bsqrt();
1820
1821 Calculate the square root of $x.
1822
1823 =head2 broot()
1824
1825         $x->broot($n);
1826
1827 Calculate the N'th root of $x.
1828
1829 =head2 badd()/bmul()/bsub()/bdiv()/bdec()/binc()
1830
1831 Please see the documentation in L<Math::BigInt>.
1832
1833 =head2 copy()
1834
1835         my $z = $x->copy();
1836
1837 Makes a deep copy of the object.
1838
1839 Please see the documentation in L<Math::BigInt> for further details.
1840
1841 =head2 bstr()/bsstr()
1842
1843         my $x = Math::BigInt->new('8/4');
1844         print $x->bstr(),"\n";                  # prints 1/2
1845         print $x->bsstr(),"\n";                 # prints 1/2
1846
1847 Return a string representating this object.
1848
1849 =head2 bacmp()/bcmp()
1850
1851 Used to compare numbers.
1852
1853 Please see the documentation in L<Math::BigInt> for further details.
1854
1855 =head2 blsft()/brsft()
1856
1857 Used to shift numbers left/right.
1858
1859 Please see the documentation in L<Math::BigInt> for further details.
1860
1861 =head2 bpow()
1862
1863         $x->bpow($y);
1864
1865 Compute $x ** $y.
1866
1867 Please see the documentation in L<Math::BigInt> for further details.
1868
1869 =head2 bexp()
1870
1871         $x->bexp($accuracy);            # calculate e ** X
1872
1873 Calculates two integers A and B so that A/B is equal to C<e ** $x>, where C<e> is
1874 Euler's number.
1875
1876 This method was added in v0.20 of Math::BigRat (May 2007).
1877
1878 See also L<blog()>.
1879
1880 =head2 bnok()
1881
1882         $x->bnok($y);              # x over y (binomial coefficient n over k)
1883
1884 Calculates the binomial coefficient n over k, also called the "choose"
1885 function. The result is equivalent to:
1886
1887         ( n )      n!
1888         | - |  = -------
1889         ( k )    k!(n-k)!
1890
1891 This method was added in v0.20 of Math::BigRat (May 2007).
1892
1893 =head2 config()
1894
1895         use Data::Dumper;
1896
1897         print Dumper ( Math::BigRat->config() );
1898         print Math::BigRat->config()->{lib},"\n";
1899
1900 Returns a hash containing the configuration, e.g. the version number, lib
1901 loaded etc. The following hash keys are currently filled in with the
1902 appropriate information.
1903
1904         key             RO/RW   Description
1905                                 Example
1906         ============================================================
1907         lib             RO      Name of the Math library
1908                                 Math::BigInt::Calc
1909         lib_version     RO      Version of 'lib'
1910                                 0.30
1911         class           RO      The class of config you just called
1912                                 Math::BigRat
1913         version         RO      version number of the class you used
1914                                 0.10
1915         upgrade         RW      To which class numbers are upgraded
1916                                 undef
1917         downgrade       RW      To which class numbers are downgraded
1918                                 undef
1919         precision       RW      Global precision
1920                                 undef
1921         accuracy        RW      Global accuracy
1922                                 undef
1923         round_mode      RW      Global round mode
1924                                 even
1925         div_scale       RW      Fallback accuracy for div
1926                                 40
1927         trap_nan        RW      Trap creation of NaN (undef = no)
1928                                 undef
1929         trap_inf        RW      Trap creation of +inf/-inf (undef = no)
1930                                 undef
1931
1932 By passing a reference to a hash you may set the configuration values. This
1933 works only for values that a marked with a C<RW> above, anything else is
1934 read-only.
1935
1936 =head1 BUGS
1937
1938 Some things are not yet implemented, or only implemented half-way:
1939
1940 =over 2
1941
1942 =item inf handling (partial)
1943
1944 =item NaN handling (partial)
1945
1946 =item rounding (not implemented except for bceil/bfloor)
1947
1948 =item $x ** $y where $y is not an integer
1949
1950 =item bmod(), blog(), bmodinv() and bmodpow() (partial)
1951
1952 =back
1953
1954 =head1 LICENSE
1955
1956 This program is free software; you may redistribute it and/or modify it under
1957 the same terms as Perl itself.
1958
1959 =head1 SEE ALSO
1960
1961 L<Math::BigFloat> and L<Math::Big> as well as L<Math::BigInt::BitVect>,
1962 L<Math::BigInt::Pari> and  L<Math::BigInt::GMP>.
1963
1964 See L<http://search.cpan.org/search?dist=bignum> for a way to use
1965 Math::BigRat.
1966
1967 The package at L<http://search.cpan.org/search?dist=Math%3A%3ABigRat>
1968 may contain more documentation and examples as well as testcases.
1969
1970 =head1 AUTHORS
1971
1972 (C) by Tels L<http://bloodgate.com/> 2001 - 2007.
1973
1974 =cut