lib/Math/BigFloat.pm

   1 #!/usr/bin/perl -w
   2
   3 # The following hash values are internally used:
   4 #   _e: exponent (BigInt)
   5 #   _m: mantissa (absolute BigInt)
   6 # sign: +,-,"NaN" if not a number
   7 #   _a: accuracy
   8 #   _p: precision
   9 #   _f: flags, used to signal MBI not to touch our private parts
  10 # _cow: Copy-On-Write (NRY)
  11
  12 package Math::BigFloat;
  13
  14 $VERSION = '1.21';
  15 require 5.005;
  16 use Exporter;
  17 use Math::BigInt qw/objectify/;
  18 @ISA =       qw( Exporter Math::BigInt);
  19 # can not export bneg/babs since the are only in MBI
  20 @EXPORT_OK = qw(
  21                 bcmp
  22                 badd bmul bdiv bmod bnorm bsub
  23                 bgcd blcm bround bfround
  24                 bpow bnan bzero bfloor bceil
  25                 bacmp bstr binc bdec binf
  26                 is_odd is_even is_nan is_inf is_positive is_negative
  27                 is_zero is_one sign
  28                );
  29
  30 #@EXPORT = qw( );
  31 use strict;
  32 use vars qw/$AUTOLOAD $accuracy $precision $div_scale $rnd_mode/;
  33 my $class = "Math::BigFloat";
  34
  35 use overload
  36 '<=>'   =>      sub { $_[2] ?
  37                       ref($_[0])->bcmp($_[1],$_[0]) :
  38                       ref($_[0])->bcmp($_[0],$_[1])},
  39 'int'   =>      sub { $_[0]->as_number() },             # 'trunc' to bigint
  40 ;
  41
  42 ##############################################################################
  43 # global constants, flags and accessory
  44
  45 use constant MB_NEVER_ROUND => 0x0001;
  46
  47 # are NaNs ok?
  48 my $NaNOK=1;
  49 # constant for easier life
  50 my $nan = 'NaN';
  51
  52 # Rounding modes one of 'even', 'odd', '+inf', '-inf', 'zero' or 'trunc'
  53 $rnd_mode = 'even';
  54 $accuracy = undef;
  55 $precision = undef;
  56 $div_scale = 40;
  57
  58 # in case we call SUPER::->foo() and this wants to call modify()
  59 # sub modify () { 0; }
  60
  61 {
  62   # checks for AUTOLOAD
  63   my %methods = map { $_ => 1 }
  64    qw / fadd fsub fmul fdiv fround ffround fsqrt fmod fstr fsstr fpow fnorm
  65         fabs fneg fint fcmp fzero fnan finc fdec
  66       /;
  67
  68   sub method_valid { return exists $methods{$_[0]||''}; }
  69 }
  70
  71 ##############################################################################
  72 # constructors
  73
  74 sub new
  75   {
  76   # create a new BigFloat object from a string or another bigfloat object.
  77   # _e: exponent
  78   # _m: mantissa
  79   # sign  => sign (+/-), or "NaN"
  80
  81   my $class = shift;
  82
  83   my $wanted = shift; # avoid numify call by not using || here
  84   return $class->bzero() if !defined $wanted;      # default to 0
  85   return $wanted->copy() if ref($wanted) eq $class;
  86
  87   my $round = shift; $round = 0 if !defined $round; # no rounding as default
  88   my $self = {}; bless $self, $class;
  89   # shortcut for bigints and its subclasses
  90   if ((ref($wanted)) && (ref($wanted) ne $class))
  91     {
  92     $self->{_m} = $wanted->as_number();         # get us a bigint copy
  93     $self->{_e} = Math::BigInt->new(0);
  94     $self->{_m}->babs();
  95     $self->{sign} = $wanted->sign();
  96     return $self->bnorm();
  97     }
  98   # got string
  99   # handle '+inf', '-inf' first
 100   if ($wanted =~ /^[+-]inf$/)
 101     {
 102     $self->{_e} = Math::BigInt->new(0);
 103     $self->{_m} = Math::BigInt->new(0);
 104     $self->{sign} = $wanted;
 105     return $self->bnorm();
 106     }
 107   #print "new string '$wanted'\n";
 108   my ($mis,$miv,$mfv,$es,$ev) = Math::BigInt::_split(\$wanted);
 109   if (!ref $mis)
 110     {
 111     die "$wanted is not a number initialized to $class" if !$NaNOK;
 112     $self->{_e} = Math::BigInt->new(0);
 113     $self->{_m} = Math::BigInt->new(0);
 114     $self->{sign} = $nan;
 115     }
 116   else
 117     {
 118     # make integer from mantissa by adjusting exp, then convert to bigint
 119     $self->{_e} = Math::BigInt->new("$$es$$ev");        # exponent
 120     $self->{_m} = Math::BigInt->new("$$mis$$miv$$mfv"); # create mantissa
 121     # 3.123E0 = 3123E-3, and 3.123E-2 => 3123E-5
 122     $self->{_e} -= CORE::length($$mfv);
 123     $self->{sign} = $self->{_m}->sign(); $self->{_m}->babs();
 124     }
 125   #print "$wanted => $self->{sign} $self->{value}\n";
 126   $self->bnorm();       # first normalize
 127   # if any of the globals is set, round to them and thus store them insid $self
 128   $self->round($accuracy,$precision,$rnd_mode)
 129    if defined $accuracy || defined $precision;
 130   return $self;
 131   }
 132
 133 sub bnan
 134   {
 135   # create a bigfloat 'NaN', if given a BigFloat, set it to 'NaN'
 136   my $self = shift;
 137   $self = $class if !defined $self;
 138   if (!ref($self))
 139     {
 140     my $c = $self; $self = {}; bless $self, $c;
 141     }
 142   $self->{_m} = Math::BigInt->bzero();
 143   $self->{_e} = Math::BigInt->bzero();
 144   $self->{sign} = $nan;
 145   return $self;
 146   }
 147
 148 sub binf
 149   {
 150   # create a bigfloat '+-inf', if given a BigFloat, set it to '+-inf'
 151   my $self = shift;
 152   my $sign = shift; $sign = '+' if !defined $sign || $sign ne '-';
 153
 154   $self = $class if !defined $self;
 155   if (!ref($self))
 156     {
 157     my $c = $self; $self = {}; bless $self, $c;
 158     }
 159   $self->{_m} = Math::BigInt->bzero();
 160   $self->{_e} = Math::BigInt->bzero();
 161   $self->{sign} = $sign.'inf';
 162   return $self;
 163   }
 164
 165 sub bone
 166   {
 167   # create a bigfloat '+-1', if given a BigFloat, set it to '+-1'
 168   my $self = shift;
 169   my $sign = shift; $sign = '+' if !defined $sign || $sign ne '-';
 170
 171   $self = $class if !defined $self;
 172   if (!ref($self))
 173     {
 174     my $c = $self; $self = {}; bless $self, $c;
 175     }
 176   $self->{_m} = Math::BigInt->bone();
 177   $self->{_e} = Math::BigInt->bzero();
 178   $self->{sign} = $sign;
 179   return $self;
 180   }
 181
 182 sub bzero
 183   {
 184   # create a bigfloat '+0', if given a BigFloat, set it to 0
 185   my $self = shift;
 186   $self = $class if !defined $self;
 187   if (!ref($self))
 188     {
 189     my $c = $self; $self = {}; bless $self, $c;
 190     }
 191   $self->{_m} = Math::BigInt->bzero();
 192   $self->{_e} = Math::BigInt->bone();
 193   $self->{sign} = '+';
 194   return $self;
 195   }
 196
 197 ##############################################################################
 198 # string conversation
 199
 200 sub bstr
 201   {
 202   # (ref to BFLOAT or num_str ) return num_str
 203   # Convert number from internal format to (non-scientific) string format.
 204   # internal format is always normalized (no leading zeros, "-0" => "+0")
 205   my ($self,$x) = objectify(1,@_);
 206
 207   #die "Oups! e was $nan" if $x->{_e}->{sign} eq $nan;
 208   #die "Oups! m was $nan" if $x->{_m}->{sign} eq $nan;
 209   if ($x->{sign} !~ /^[+-]$/)
 210     {
 211     return $x->{sign} unless $x->{sign} eq '+inf';      # -inf, NaN
 212     return 'inf';                                       # +inf
 213     }
 214
 215   my $es = '0'; my $len = 1; my $cad = 0; my $dot = '.';
 216
 217   my $not_zero = !$x->is_zero();
 218   if ($not_zero)
 219     {
 220     $es = $x->{_m}->bstr();
 221     $len = CORE::length($es);
 222     if (!$x->{_e}->is_zero())
 223 #      {
 224 #      $es = $x->{sign}.$es if $x->{sign} eq '-';
 225 #      }
 226 #    else
 227       {
 228       if ($x->{_e}->sign() eq '-')
 229         {
 230         $dot = '';
 231         if ($x->{_e} <= -$len)
 232           {
 233           # print "style: 0.xxxx\n";
 234           my $r = $x->{_e}->copy(); $r->babs()->bsub( CORE::length($es) );
 235           $es = '0.'. ('0' x $r) . $es; $cad = -($len+$r);
 236           }
 237         else
 238           {
 239           # print "insert '.' at $x->{_e} in '$es'\n";
 240           substr($es,$x->{_e},0) = '.'; $cad = $x->{_e};
 241           }
 242         }
 243       else
 244         {
 245         # expand with zeros
 246         $es .= '0' x $x->{_e}; $len += $x->{_e}; $cad = 0;
 247         }
 248       }
 249     } # if not zero
 250   $es = $x->{sign}.$es if $x->{sign} eq '-';
 251   # if set accuracy or precision, pad with zeros
 252   if ((defined $x->{_a}) && ($not_zero))
 253     {
 254     # 123400 => 6, 0.1234 => 4, 0.001234 => 4
 255     my $zeros = $x->{_a} - $cad;                # cad == 0 => 12340
 256     $zeros = $x->{_a} - $len if $cad != $len;
 257     #print "acc padd $x->{_a} $zeros (len $len cad $cad)\n";
 258     $es .= $dot.'0' x $zeros if $zeros > 0;
 259     }
 260   elsif ($x->{_p} || 0 < 0)
 261     {
 262     # 123400 => 6, 0.1234 => 4, 0.001234 => 6
 263     my $zeros = -$x->{_p} + $cad;
 264     #print "pre padd $x->{_p} $zeros (len $len cad $cad)\n";
 265     $es .= $dot.'0' x $zeros if $zeros > 0;
 266     }
 267   return $es;
 268   }
 269
 270 sub bsstr
 271   {
 272   # (ref to BFLOAT or num_str ) return num_str
 273   # Convert number from internal format to scientific string format.
 274   # internal format is always normalized (no leading zeros, "-0E0" => "+0E0")
 275   my ($self,$x) = objectify(1,@_);
 276
 277   #die "Oups! e was $nan" if $x->{_e}->{sign} eq $nan;
 278   #die "Oups! m was $nan" if $x->{_m}->{sign} eq $nan;
 279   if ($x->{sign} !~ /^[+-]$/)
 280     {
 281     return $x->{sign} unless $x->{sign} eq '+inf';      # -inf, NaN
 282     return 'inf';                                       # +inf
 283     }
 284   my $sign = $x->{_e}->{sign}; $sign = '' if $sign eq '-';
 285   my $sep = 'e'.$sign;
 286   return $x->{_m}->bstr().$sep.$x->{_e}->bstr();
 287   }
 288
 289 sub numify
 290   {
 291   # Make a number from a BigFloat object
 292   # simple return string and let Perl's atoi()/atof() handle the rest
 293   my ($self,$x) = objectify(1,@_);
 294   return $x->bsstr();
 295   }
 296
 297 ##############################################################################
 298 # public stuff (usually prefixed with "b")
 299
 300 # really? Just for exporting them is not what I had in mind
 301 #sub babs
 302 #  {
 303 #  $class->SUPER::babs($class,@_);
 304 #  }
 305 #sub bneg
 306 #  {
 307 #  $class->SUPER::bneg($class,@_);
 308 #  }
 309
 310 # tels 2001-08-04
 311 # todo: this must be overwritten and return NaN for non-integer values
 312 # band(), bior(), bxor(), too
 313 #sub bnot
 314 #  {
 315 #  $class->SUPER::bnot($class,@_);
 316 #  }
 317
 318 sub bcmp
 319   {
 320   # Compares 2 values.  Returns one of undef, <0, =0, >0. (suitable for sort)
 321   # (BFLOAT or num_str, BFLOAT or num_str) return cond_code
 322   my ($self,$x,$y) = objectify(2,@_);
 323
 324   if (($x->{sign} !~ /^[+-]$/) || ($y->{sign} !~ /^[+-]$/))
 325     {
 326     # handle +-inf and NaN
 327     return undef if (($x->{sign} eq $nan) || ($y->{sign} eq $nan));
 328     return 0 if ($x->{sign} eq $y->{sign}) && ($x->{sign} =~ /^[+-]inf$/);
 329     return +1 if $x->{sign} eq '+inf';
 330     return -1 if $x->{sign} eq '-inf';
 331     return -1 if $y->{sign} eq '+inf';
 332     return +1 if $y->{sign} eq '-inf';
 333     }
 334
 335   # check sign for speed first
 336   return 1 if $x->{sign} eq '+' && $y->{sign} eq '-';   # does also 0 <=> -y
 337   return -1 if $x->{sign} eq '-' && $y->{sign} eq '+';  # does also -x <=> 0
 338
 339   # shortcut
 340   my $xz = $x->is_zero();
 341   my $yz = $y->is_zero();
 342   return 0 if $xz && $yz;                               # 0 <=> 0
 343   return -1 if $xz && $y->{sign} eq '+';                # 0 <=> +y
 344   return 1 if $yz && $x->{sign} eq '+';                 # +x <=> 0
 345
 346   # adjust so that exponents are equal
 347   my $lxm = $x->{_m}->length();
 348   my $lym = $y->{_m}->length();
 349   my $lx = $lxm + $x->{_e};
 350   my $ly = $lym + $y->{_e};
 351   # print "x $x y $y lx $lx ly $ly\n";
 352   my $l = $lx - $ly; $l = -$l if $x->{sign} eq '-';
 353   # print "$l $x->{sign}\n";
 354   return $l <=> 0 if $l != 0;
 355
 356   # lengths (corrected by exponent) are equal
 357   # so make mantissa euqal length by padding with zero (shift left)
 358   my $diff = $lxm - $lym;
 359   my $xm = $x->{_m};            # not yet copy it
 360   my $ym = $y->{_m};
 361   if ($diff > 0)
 362     {
 363     $ym = $y->{_m}->copy()->blsft($diff,10);
 364     }
 365   elsif ($diff < 0)
 366     {
 367     $xm = $x->{_m}->copy()->blsft(-$diff,10);
 368     }
 369   my $rc = $xm->bcmp($ym);
 370   $rc = -$rc if $x->{sign} eq '-';              # -124 < -123
 371   return $rc <=> 0;
 372   }
 373
 374 sub bacmp
 375   {
 376   # Compares 2 values, ignoring their signs.
 377   # Returns one of undef, <0, =0, >0. (suitable for sort)
 378   # (BFLOAT or num_str, BFLOAT or num_str) return cond_code
 379   my ($self,$x,$y) = objectify(2,@_);
 380   return undef if (($x->{sign} eq $nan) || ($y->{sign} eq $nan));
 381
 382   # signs are ignored, so check length
 383   # length(x) is length(m)+e aka length of non-fraction part
 384   # the longer one is bigger
 385   my $l = $x->length() - $y->length();
 386   #print "$l\n";
 387   return $l if $l != 0;
 388   #print "equal lengths\n";
 389
 390   # if both are equal long, make full compare
 391   # first compare only the mantissa
 392   # if mantissa are equal, compare fractions
 393
 394   return $x->{_m} <=> $y->{_m} || $x->{_e} <=> $y->{_e};
 395   }
 396
 397 sub badd
 398   {
 399   # add second arg (BFLOAT or string) to first (BFLOAT) (modifies first)
 400   # return result as BFLOAT
 401   my ($self,$x,$y,$a,$p,$r) = objectify(2,@_);
 402
 403   # inf and NaN handling
 404   if (($x->{sign} !~ /^[+-]$/) || ($y->{sign} !~ /^[+-]$/))
 405     {
 406     # NaN first
 407     return $x->bnan() if (($x->{sign} eq $nan) || ($y->{sign} eq $nan));
 408     # inf handline
 409     if (($x->{sign} =~ /^[+-]inf$/) && ($y->{sign} =~ /^[+-]inf$/))
 410       {
 411       # + and + => +, - and - => -, + and - => 0, - and + => 0
 412       return $x->bzero() if $x->{sign} ne $y->{sign};
 413       return $x;
 414       }
 415     # +-inf + something => +inf
 416     # something +-inf => +-inf
 417     $x->{sign} = $y->{sign}, return $x if $y->{sign} =~ /^[+-]inf$/;
 418     return $x;
 419     }
 420
 421   # speed: no add for 0+y or x+0
 422   return $x if $y->is_zero();                           # x+0
 423   if ($x->is_zero())                                    # 0+y
 424     {
 425     # make copy, clobbering up x (modify in place!)
 426     $x->{_e} = $y->{_e}->copy();
 427     $x->{_m} = $y->{_m}->copy();
 428     $x->{sign} = $y->{sign} || $nan;
 429     return $x->round($a,$p,$r,$y);
 430     }
 431
 432   # take lower of the two e's and adapt m1 to it to match m2
 433   my $e = $y->{_e}; $e = Math::BigInt::bzero() if !defined $e;  # if no BFLOAT
 434   $e = $e - $x->{_e};
 435   my $add = $y->{_m}->copy();
 436   if ($e < 0)
 437     {
 438     # print "e < 0\n";
 439     #print "\$x->{_m}: $x->{_m} ";
 440     #print "\$x->{_e}: $x->{_e}\n";
 441     my $e1 = $e->copy()->babs();
 442     $x->{_m} *= (10 ** $e1);
 443     $x->{_e} += $e;                     # need the sign of e
 444     #$x->{_m} += $y->{_m};
 445     #print "\$x->{_m}: $x->{_m} ";
 446     #print "\$x->{_e}: $x->{_e}\n";
 447     }
 448   elsif ($e > 0)
 449     {
 450     # print "e > 0\n";
 451     #print "\$x->{_m}: $x->{_m} \$y->{_m}: $y->{_m} \$e: $e ",ref($e),"\n";
 452     $add *= (10 ** $e);
 453     #$x->{_m} += $y->{_m} * (10 ** $e);
 454     #print "\$x->{_m}: $x->{_m}\n";
 455     }
 456   # else: both e are same, so leave them
 457   #print "badd $x->{sign}$x->{_m} +  $y->{sign}$add\n";
 458   # fiddle with signs
 459   $x->{_m}->{sign} = $x->{sign};
 460   $add->{sign} = $y->{sign};
 461   # finally do add/sub
 462   $x->{_m} += $add;
 463   # re-adjust signs
 464   $x->{sign} = $x->{_m}->{sign};
 465   $x->{_m}->{sign} = '+';
 466   #$x->bnorm();                         # delete trailing zeros
 467   return $x->round($a,$p,$r,$y);
 468   }
 469
 470 sub bsub
 471   {
 472   # (BigFloat or num_str, BigFloat or num_str) return BigFloat
 473   # subtract second arg from first, modify first
 474   my ($self,$x,$y) = objectify(2,@_);
 475
 476   $x->badd($y->bneg()); # badd does not leave internal zeros
 477   $y->bneg();           # refix y, assumes no one reads $y in between
 478   return $x;            # badd() already normalized and rounded
 479   }
 480
 481 sub binc
 482   {
 483   # increment arg by one
 484   my ($self,$x,$a,$p,$r) = objectify(1,@_);
 485   $x->badd($self->_one())->round($a,$p,$r);
 486   }
 487
 488 sub bdec
 489   {
 490   # decrement arg by one
 491   my ($self,$x,$a,$p,$r) = objectify(1,@_);
 492   $x->badd($self->_one('-'))->round($a,$p,$r);
 493   }
 494
 495 sub blcm
 496   {
 497   # (BINT or num_str, BINT or num_str) return BINT
 498   # does not modify arguments, but returns new object
 499   # Lowest Common Multiplicator
 500
 501   my ($self,@arg) = objectify(0,@_);
 502   my $x = $self->new(shift @arg);
 503   while (@arg) { $x = _lcm($x,shift @arg); }
 504   $x;
 505   }
 506
 507 sub bgcd
 508   {
 509   # (BINT or num_str, BINT or num_str) return BINT
 510   # does not modify arguments, but returns new object
 511   # GCD -- Euclids algorithm Knuth Vol 2 pg 296
 512
 513   my ($self,@arg) = objectify(0,@_);
 514   my $x = $self->new(shift @arg);
 515   while (@arg) { $x = _gcd($x,shift @arg); }
 516   $x;
 517   }
 518
 519 sub is_zero
 520   {
 521   # return true if arg (BINT or num_str) is zero (array '+', '0')
 522   my $x = shift; $x = $class->new($x) unless ref $x;
 523
 524   return 1 if $x->{sign} eq '+' && $x->{_m}->is_zero();
 525   return 0;
 526   }
 527
 528 sub is_one
 529   {
 530   # return true if arg (BINT or num_str) is +1 (array '+', '1')
 531   # or -1 if signis given
 532   my $x = shift; $x = $class->new($x) unless ref $x;
 533   #my ($self,$x) = objectify(1,@_);
 534   my $sign = $_[2] || '+';
 535   return ($x->{sign} eq $sign && $x->{_e}->is_zero() && $x->{_m}->is_one());
 536   }
 537
 538 sub is_odd
 539   {
 540   # return true if arg (BINT or num_str) is odd or false if even
 541   my $x = shift; $x = $class->new($x) unless ref $x;
 542   #my ($self,$x) = objectify(1,@_);
 543
 544   return 0 if $x->{sign} !~ /^[+-]$/;                   # NaN & +-inf aren't
 545   return ($x->{_e}->is_zero() && $x->{_m}->is_odd());
 546   }
 547
 548 sub is_even
 549   {
 550   # return true if arg (BINT or num_str) is even or false if odd
 551   my $x = shift; $x = $class->new($x) unless ref $x;
 552   #my ($self,$x) = objectify(1,@_);
 553
 554   return 0 if $x->{sign} !~ /^[+-]$/;                   # NaN & +-inf aren't
 555   return 1 if $x->{_m}->is_zero();                      # 0e1 is even
 556   return ($x->{_e}->is_zero() && $x->{_m}->is_even());  # 123.45 is never
 557   }
 558
 559 sub bmul
 560   {
 561   # multiply two numbers -- stolen from Knuth Vol 2 pg 233
 562   # (BINT or num_str, BINT or num_str) return BINT
 563   my ($self,$x,$y,$a,$p,$r) = objectify(2,@_);
 564
 565   # print "mbf bmul $x->{_m}e$x->{_e} $y->{_m}e$y->{_e}\n";
 566   return $x->bnan() if (($x->{sign} eq $nan) || ($y->{sign} eq $nan));
 567
 568   # handle result = 0
 569   return $x->bzero() if $x->is_zero() || $y->is_zero();
 570   # inf handling
 571   if (($x->{sign} =~ /^[+-]inf$/) || ($y->{sign} =~ /^[+-]inf$/))
 572     {
 573     # result will always be +-inf:
 574     # +inf * +/+inf => +inf, -inf * -/-inf => +inf
 575     # +inf * -/-inf => -inf, -inf * +/+inf => -inf
 576     return $x->binf() if ($x->{sign} =~ /^\+/ && $y->{sign} =~ /^\+/);
 577     return $x->binf() if ($x->{sign} =~ /^-/ && $y->{sign} =~ /^-/);
 578     return $x->binf('-');
 579     }
 580
 581   # aEb * cEd = (a*c)E(b+d)
 582   $x->{_m} = $x->{_m} * $y->{_m};
 583   #print "m: $x->{_m}\n";
 584   $x->{_e} = $x->{_e} + $y->{_e};
 585   #print "e: $x->{_m}\n";
 586   # adjust sign:
 587   $x->{sign} = $x->{sign} ne $y->{sign} ? '-' : '+';
 588   #print "s: $x->{sign}\n";
 589   $x->bnorm();
 590   return $x->round($a,$p,$r,$y);
 591   }
 592
 593 sub bdiv
 594   {
 595   # (dividend: BFLOAT or num_str, divisor: BFLOAT or num_str) return
 596   # (BFLOAT,BFLOAT) (quo,rem) or BINT (only rem)
 597   my ($self,$x,$y,$a,$p,$r) = objectify(2,@_);
 598
 599   # x / +-inf => 0, reminder x
 600   return wantarray ? ($x->bzero(),$x->copy()) : $x->bzero()
 601    if $y->{sign} =~ /^[+-]inf$/;
 602
 603   # NaN if x == NaN or y == NaN or x==y==0
 604   return wantarray ? ($x->bnan(),bnan()) : $x->bnan()
 605    if (($x->is_nan() || $y->is_nan()) ||
 606       ($x->is_zero() && $y->is_zero()));
 607
 608   # 5 / 0 => +inf, -6 / 0 => -inf
 609   return wantarray
 610    ? ($x->binf($x->{sign}),$self->bnan()) : $x->binf($x->{sign})
 611    if ($x->{sign} =~ /^[+-]$/ && $y->is_zero());
 612
 613   $y = $class->new($y) if ref($y) ne $class;            # promote bigints
 614
 615   # print "mbf bdiv $x ",ref($x)," ",$y," ",ref($y),"\n";
 616   # we need to limit the accuracy to protect against overflow
 617   my ($scale) = $x->_scale_a($accuracy,$rnd_mode,$a,$r);        # ignore $p
 618   my $fallback = 0;
 619   if (!defined $scale)
 620     {
 621     # simulate old behaviour
 622     $scale = $div_scale+1;      # one more for proper riund
 623     $a = $div_scale;            # and round to it
 624     $fallback = 1;              # to clear a/p afterwards
 625     }
 626   my $lx = $x->{_m}->length(); my $ly = $y->{_m}->length();
 627   $scale = $lx if $lx > $scale;
 628   $scale = $ly if $ly > $scale;
 629   #print "scale $scale $lx $ly\n";
 630   my $diff = $ly - $lx;
 631   $scale += $diff if $diff > 0;         # if lx << ly, but not if ly << lx!
 632
 633   return wantarray ? ($x,$self->bzero()) : $x if $x->is_zero();
 634
 635   $x->{sign} = $x->{sign} ne $y->sign() ? '-' : '+';
 636
 637   # check for / +-1 ( +/- 1E0)
 638   if ($y->is_one())
 639     {
 640     return wantarray ? ($x,$self->bzero()) : $x;
 641     }
 642
 643   # a * 10 ** b / c * 10 ** d => a/c * 10 ** (b-d)
 644   #print "self: $self x: $x ref(x) ", ref($x)," m: $x->{_m}\n";
 645   # my $scale_10 = 10 ** $scale; $x->{_m}->bmul($scale_10);
 646   $x->{_m}->blsft($scale,10);
 647   #print "m: $x->{_m} $y->{_m}\n";
 648   $x->{_m}->bdiv( $y->{_m} );   # a/c
 649   #print "m: $x->{_m}\n";
 650   #print "e: $x->{_e} $y->{_e}",$scale,"\n";
 651   $x->{_e}->bsub($y->{_e});     # b-d
 652   #print "e: $x->{_e}\n";
 653   $x->{_e}->bsub($scale);       # correct for 10**scale
 654   #print "after div: m: $x->{_m} e: $x->{_e}\n";
 655   $x->bnorm();                  # remove trailing 0's
 656   #print "after div: m: $x->{_m} e: $x->{_e}\n";
 657   $x->round($a,$p,$r);          # then round accordingly
 658   if ($fallback)
 659     {
 660     # clear a/p after round, since user did not request it
 661     $x->{_a} = undef;
 662     $x->{_p} = undef;
 663     }
 664
 665   if (wantarray)
 666     {
 667     my $rem = $x->copy();
 668     $rem->bmod($y,$a,$p,$r);
 669     if ($fallback)
 670       {
 671       # clear a/p after round, since user did not request it
 672       $x->{_a} = undef;
 673       $x->{_p} = undef;
 674       }
 675     return ($x,$rem);
 676     }
 677   return $x;
 678   }
 679
 680 sub bmod
 681   {
 682   # (dividend: BFLOAT or num_str, divisor: BFLOAT or num_str) return reminder
 683   my ($self,$x,$y,$a,$p,$r) = objectify(2,@_);
 684
 685   return $x->bnan() if ($x->{sign} eq $nan || $y->is_nan() || $y->is_zero());
 686   return $x->bzero() if $y->is_one();
 687
 688   # XXX tels: not done yet
 689   return $x->round($a,$p,$r,$y);
 690   }
 691
 692 sub bsqrt
 693   {
 694   # calculate square root; this should probably
 695   # use a different test to see whether the accuracy we want is...
 696   my ($self,$x,$a,$p,$r) = objectify(1,@_);
 697
 698   return $x->bnan() if $x->{sign} eq 'NaN' || $x->{sign} =~ /^-/; # <0, NaN
 699   return $x if $x->{sign} eq '+inf';                              # +inf
 700   return $x if $x->is_zero() || $x == 1;
 701
 702   # we need to limit the accuracy to protect against overflow
 703   my ($scale) = $x->_scale_a($accuracy,$rnd_mode,$a,$r);        # ignore $p
 704   my $fallback = 0;
 705   if (!defined $scale)
 706     {
 707     # simulate old behaviour
 708     $scale = $div_scale+1;      # one more for proper riund
 709     $a = $div_scale;            # and round to it
 710     $fallback = 1;              # to clear a/p afterwards
 711     }
 712   my $lx = $x->{_m}->length();
 713   $scale = $lx if $scale < $lx;
 714   my $e = Math::BigFloat->new("1E-$scale");     # make test variable
 715   return $x->bnan() if $e->sign() eq 'NaN';
 716
 717   # start with some reasonable guess
 718   #$x *= 10 ** ($len - $org->{_e}); $x /= 2;    # !?!?
 719   $lx = $lx+$x->{_e};
 720   $lx = 1 if $lx < 1;
 721   my $gs = Math::BigFloat->new('1'. ('0' x $lx));
 722
 723   # print "first guess: $gs (x $x) scale $scale\n";
 724
 725   my $diff = $e;
 726   my $y = $x->copy();
 727   my $two = Math::BigFloat->new(2);
 728   $x = Math::BigFloat->new($x) if ref($x) ne $class;    # promote BigInts
 729   # $scale = 2;
 730   while ($diff >= $e)
 731     {
 732     return $x->bnan() if $gs->is_zero();
 733     $r = $y->copy(); $r->bdiv($gs,$scale);
 734     $x = ($r + $gs);
 735     $x->bdiv($two,$scale);
 736     $diff = $x->copy()->bsub($gs)->babs();
 737     $gs = $x->copy();
 738     }
 739   $x->round($a,$p,$r);
 740   if ($fallback)
 741     {
 742     # clear a/p after round, since user did not request it
 743     $x->{_a} = undef;
 744     $x->{_p} = undef;
 745     }
 746   $x;
 747   }
 748
 749 sub bpow
 750   {
 751   # (BFLOAT or num_str, BFLOAT or num_str) return BFLOAT
 752   # compute power of two numbers, second arg is used as integer
 753   # modifies first argument
 754
 755   my ($self,$x,$y,$a,$p,$r) = objectify(2,@_);
 756
 757   return $x if $x->{sign} =~ /^[+-]inf$/;
 758   return $x->bnan() if $x->{sign} eq $nan || $y->{sign} eq $nan;
 759   return $x->bone() if $y->is_zero();
 760   return $x         if $x->is_one() || $y->is_one();
 761   my $y1 = $y->as_number();             # make bigint
 762   if ($x == -1)
 763     {
 764     # if $x == -1 and odd/even y => +1/-1  because +-1 ^ (+-1) => +-1
 765     return $y1->is_odd() ? $x : $x->babs(1);
 766     }
 767   return $x if $x->is_zero() && $y->{sign} eq '+'; # 0**y => 0 (if not y <= 0)
 768   # 0 ** -y => 1 / (0 ** y) => / 0! (1 / 0 => +inf)
 769   return $x->binf() if $x->is_zero() && $y->{sign} eq '-';
 770
 771   # calculate $x->{_m} ** $y and $x->{_e} * $y separately (faster)
 772   $y1->babs();
 773   $x->{_m}->bpow($y1);
 774   $x->{_e}->bmul($y1);
 775   $x->{sign} = $nan if $x->{_m}->{sign} eq $nan || $x->{_e}->{sign} eq $nan;
 776   $x->bnorm();
 777   if ($y->{sign} eq '-')
 778     {
 779     # modify $x in place!
 780     my $z = $x->copy(); $x->bzero()->binc();
 781     return $x->bdiv($z,$a,$p,$r);       # round in one go (might ignore y's A!)
 782     }
 783   return $x->round($a,$p,$r,$y);
 784   }
 785
 786 ###############################################################################
 787 # rounding functions
 788
 789 sub bfround
 790   {
 791   # precision: round to the $Nth digit left (+$n) or right (-$n) from the '.'
 792   # $n == 0 means round to integer
 793   # expects and returns normalized numbers!
 794   my $x = shift; $x = $class->new($x) unless ref $x;
 795
 796   return $x if $x->modify('bfround');
 797
 798   my ($scale,$mode) = $x->_scale_p($precision,$rnd_mode,@_);
 799   return $x if !defined $scale;                 # no-op
 800
 801   # never round a 0, +-inf, NaN
 802   return $x if $x->{sign} !~ /^[+-]$/ || $x->is_zero();
 803   # print "MBF bfround $x to scale $scale mode $mode\n";
 804
 805   if ($scale < 0)
 806     {
 807     # print "bfround scale $scale e $x->{_e}\n";
 808     # round right from the '.'
 809     return $x if $x->{_e} >= 0;                 # nothing to round
 810     $scale = -$scale;                           # positive for simplicity
 811     my $len = $x->{_m}->length();               # length of mantissa
 812     my $dad = -$x->{_e};                        # digits after dot
 813     my $zad = 0;                                # zeros after dot
 814     $zad = -$len-$x->{_e} if ($x->{_e} < -$len);# for 0.00..00xxx style
 815     # print "scale $scale dad $dad zad $zad len $len\n";
 816
 817     # number  bsstr   len zad dad
 818     # 0.123   123e-3    3   0 3
 819     # 0.0123  123e-4    3   1 4
 820     # 0.001   1e-3      1   2 3
 821     # 1.23    123e-2    3   0 2
 822     # 1.2345  12345e-4  5   0 4
 823
 824     # do not round after/right of the $dad
 825     return $x if $scale > $dad;                 # 0.123, scale >= 3 => exit
 826
 827      # round to zero if rounding inside the $zad, but not for last zero like:
 828      # 0.0065, scale -2, round last '0' with following '65' (scale == zad case)
 829      if ($scale < $zad)
 830       {
 831       return $x->bzero();
 832       }
 833     if ($scale == $zad)                 # for 0.006, scale -2 and trunc
 834       {
 835       $scale = -$len;
 836       }
 837     else
 838       {
 839       # adjust round-point to be inside mantissa
 840       if ($zad != 0)
 841         {
 842         $scale = $scale-$zad;
 843         }
 844       else
 845         {
 846         my $dbd = $len - $dad; $dbd = 0 if $dbd < 0;    # digits before dot
 847         $scale = $dbd+$scale;
 848         }
 849       }
 850     # print "round to $x->{_m} to $scale\n";
 851     }
 852   else
 853     {
 854     # 123 => 100 means length(123) = 3 - $scale (2) => 1
 855
 856     # calculate digits before dot
 857     my $dbt = $x->{_m}->length(); $dbt += $x->{_e} if $x->{_e}->sign() eq '-';
 858     if (($scale > $dbt) && ($dbt < 0))
 859       {
 860       # if not enough digits before dot, round to zero
 861       return $x->bzero();
 862       }
 863     if (($scale >= 0) && ($dbt == 0))
 864       {
 865       # 0.49->bfround(1): scale == 1, dbt == 0: => 0.0
 866       # 0.51->bfround(0): scale == 0, dbt == 0: => 1.0
 867       # 0.5->bfround(0):  scale == 0, dbt == 0: => 0
 868       # 0.05->bfround(0): scale == 0, dbt == 0: => 0
 869       # print "$scale $dbt $x->{_m}\n";
 870       $scale = -$x->{_m}->length();
 871       }
 872     elsif ($dbt > 0)
 873       {
 874       # correct by subtracting scale
 875       $scale = $dbt - $scale;
 876       }
 877     else
 878       {
 879       $scale = $x->{_m}->length() - $scale;
 880       }
 881     }
 882   # print "using $scale for $x->{_m} with '$mode'\n";
 883   # pass sign to bround for rounding modes '+inf' and '-inf'
 884   $x->{_m}->{sign} = $x->{sign};
 885   $x->{_m}->bround($scale,$mode);
 886   $x->{_m}->{sign} = '+';               # fix sign back
 887   $x->bnorm();
 888   }
 889
 890 sub bround
 891   {
 892   # accuracy: preserve $N digits, and overwrite the rest with 0's
 893   my $x = shift; $x = $class->new($x) unless ref $x;
 894   my ($scale,$mode) = $x->_scale_a($accuracy,$rnd_mode,@_);
 895   return $x if !defined $scale;                 # no-op
 896
 897   return $x if $x->modify('bround');
 898
 899   # print "bround $scale $mode\n";
 900   # 0 => return all digits, scale < 0 makes no sense
 901   return $x if ($scale <= 0);
 902   # never round a 0, +-inf, NaN
 903   return $x if $x->{sign} !~ /^[+-]$/ || $x->is_zero();
 904
 905   # if $e longer than $m, we have 0.0000xxxyyy style number, and must
 906   # subtract the delta from scale, to simulate keeping the zeros
 907   # -5 +5 => 1; -10 +5 => -4
 908   my $delta = $x->{_e} + $x->{_m}->length() + 1;
 909   # removed by tlr, since causes problems with fraction tests:
 910   # $scale += $delta if $delta < 0;
 911
 912   # if we should keep more digits than the mantissa has, do nothing
 913   return $x if $x->{_m}->length() <= $scale;
 914
 915   # pass sign to bround for '+inf' and '-inf' rounding modes
 916   $x->{_m}->{sign} = $x->{sign};
 917   $x->{_m}->bround($scale,$mode);       # round mantissa
 918   $x->{_m}->{sign} = '+';               # fix sign back
 919   $x->bnorm();                          # del trailing zeros gen. by bround()
 920   }
 921
 922 sub bfloor
 923   {
 924   # return integer less or equal then $x
 925   my ($self,$x,$a,$p,$r) = objectify(1,@_);
 926
 927   return $x if $x->modify('bfloor');
 928
 929   return $x if $x->{sign} !~ /^[+-]$/;  # nan, +inf, -inf
 930
 931   # if $x has digits after dot
 932   if ($x->{_e}->{sign} eq '-')
 933     {
 934     $x->{_m}->brsft(-$x->{_e},10);
 935     $x->{_e}->bzero();
 936     $x-- if $x->{sign} eq '-';
 937     }
 938   return $x->round($a,$p,$r);
 939   }
 940
 941 sub bceil
 942   {
 943   # return integer greater or equal then $x
 944   my ($self,$x,$a,$p,$r) = objectify(1,@_);
 945
 946   return $x if $x->modify('bceil');
 947   return $x if $x->{sign} !~ /^[+-]$/;  # nan, +inf, -inf
 948
 949   # if $x has digits after dot
 950   if ($x->{_e}->{sign} eq '-')
 951     {
 952     $x->{_m}->brsft(-$x->{_e},10);
 953     $x->{_e}->bzero();
 954     $x++ if $x->{sign} eq '+';
 955     }
 956   return $x->round($a,$p,$r);
 957   }
 958
 959 ###############################################################################
 960
 961 sub DESTROY
 962   {
 963   # going trough AUTOLOAD for every DESTROY is costly, so avoid it by empty sub
 964   }
 965
 966 sub AUTOLOAD
 967   {
 968   # make fxxx and bxxx work
 969   # my $self = $_[0];
 970   my $name = $AUTOLOAD;
 971
 972   $name =~ s/.*:://;    # split package
 973   #print "$name\n";
 974   if (!method_valid($name))
 975     {
 976     #no strict 'refs';
 977     ## try one level up
 978     #&{$class."::SUPER->$name"}(@_);
 979     # delayed load of Carp and avoid recursion
 980     require Carp;
 981     Carp::croak ("Can't call $class\-\>$name, not a valid method");
 982     }
 983   no strict 'refs';
 984   my $bname = $name; $bname =~ s/^f/b/;
 985   *{$class."\:\:$name"} = \&$bname;
 986   &$bname;      # uses @_
 987   }
 988
 989 sub exponent
 990   {
 991   # return a copy of the exponent
 992   my $self = shift;
 993   $self = $class->new($self) unless ref $self;
 994
 995   return bnan() if $self->is_nan();
 996   return $self->{_e}->copy();
 997   }
 998
 999 sub mantissa
1000   {
1001   # return a copy of the mantissa
1002   my $self = shift;
1003   $self = $class->new($self) unless ref $self;
1004
1005   return bnan() if $self->is_nan();
1006   my $m = $self->{_m}->copy();  # faster than going via bstr()
1007   $m->bneg() if $self->{sign} eq '-';
1008
1009   return $m;
1010   }
1011
1012 sub parts
1013   {
1014   # return a copy of both the exponent and the mantissa
1015   my $self = shift;
1016   $self = $class->new($self) unless ref $self;
1017
1018   return (bnan(),bnan()) if $self->is_nan();
1019   my $m = $self->{_m}->copy();  # faster than going via bstr()
1020   $m->bneg() if $self->{sign} eq '-';
1021   return ($m,$self->{_e}->copy());
1022   }
1023
1024 ##############################################################################
1025 # private stuff (internal use only)
1026
1027 sub _one
1028   {
1029   # internal speedup, set argument to 1, or create a +/- 1
1030   my $self = shift; $self = ref($self) if ref($self);
1031   my $x = {}; bless $x, $self;
1032   $x->{_m} = Math::BigInt->new(1);
1033   $x->{_e} = Math::BigInt->new(0);
1034   $x->{sign} = shift || '+';
1035   return $x;
1036   }
1037
1038 sub import
1039   {
1040   my $self = shift;
1041   #print "import $self\n";
1042   for ( my $i = 0; $i < @_ ; $i++ )
1043     {
1044     if ( $_[$i] eq ':constant' )
1045       {
1046       # this rest causes overlord er load to step in
1047       # print "overload @_\n";
1048       overload::constant float => sub { $self->new(shift); };
1049       splice @_, $i, 1; last;
1050       }
1051     }
1052   # any non :constant stuff is handled by our parent, Exporter
1053   # even if @_ is empty, to give it a chance
1054   #$self->SUPER::import(@_);            # does not work (would call MBI)
1055   $self->export_to_level(1,$self,@_);   # need this instead
1056   }
1057
1058 sub bnorm
1059   {
1060   # adjust m and e so that m is smallest possible
1061   # round number according to accuracy and precision settings
1062   my $x = shift;
1063
1064   return $x if $x->{sign} !~ /^[+-]$/;          # inf, nan etc
1065
1066   my $zeros = $x->{_m}->_trailing_zeros();      # correct for trailing zeros
1067   if ($zeros != 0)
1068     {
1069     $x->{_m}->brsft($zeros,10); $x->{_e} += $zeros;
1070     }
1071   # for something like 0Ey, set y to 1
1072   $x->{sign} = '+', $x->{_e}->bzero()->binc() if $x->{_m}->is_zero();
1073   $x->{_m}->{_f} = MB_NEVER_ROUND;
1074   $x->{_e}->{_f} = MB_NEVER_ROUND;
1075   return $x;                                    # MBI bnorm is no-op
1076   }
1077
1078 ##############################################################################
1079 # internal calculation routines
1080
1081 sub as_number
1082   {
1083   # return a bigint representation of this BigFloat number
1084   my ($self,$x) = objectify(1,@_);
1085
1086   my $z;
1087   if ($x->{_e}->is_zero())
1088     {
1089     $z = $x->{_m}->copy();
1090     $z->{sign} = $x->{sign};
1091     return $z;
1092     }
1093   $z = $x->{_m}->copy();
1094   if ($x->{_e} < 0)
1095     {
1096     $z->brsft(-$x->{_e},10);
1097     }
1098   else
1099     {
1100     $z->blsft($x->{_e},10);
1101     }
1102   $z->{sign} = $x->{sign};
1103   return $z;
1104   }
1105
1106 sub length
1107   {
1108   my $x = shift; $x = $class->new($x) unless ref $x;
1109
1110   my $len = $x->{_m}->length();
1111   $len += $x->{_e} if $x->{_e}->sign() eq '+';
1112   if (wantarray())
1113     {
1114     my $t = Math::BigInt::bzero();
1115     $t = $x->{_e}->copy()->babs() if $x->{_e}->sign() eq '-';
1116     return ($len,$t);
1117     }
1118   return $len;
1119   }
1120
1121 1;
1122 __END__
1123
1124 =head1 NAME
1125
1126 Math::BigFloat - Arbitrary size floating point math package
1127
1128 =head1 SYNOPSIS
1129
1130   use Math::BigFloat;
1131
1132   # Number creation
1133   $x = Math::BigInt->new($str); # defaults to 0
1134   $nan  = Math::BigInt->bnan(); # create a NotANumber
1135   $zero = Math::BigInt->bzero();# create a "+0"
1136
1137   # Testing
1138   $x->is_zero();                # return whether arg is zero or not
1139   $x->is_nan();                 # return whether arg is NaN or not
1140   $x->is_one();                 # true if arg is +1
1141   $x->is_one('-');              # true if arg is -1
1142   $x->is_odd();                 # true if odd, false for even
1143   $x->is_even();                # true if even, false for odd
1144   $x->is_positive();            # true if >= 0
1145   $x->is_negative();            # true if <  0
1146   $x->is_inf(sign)              # true if +inf or -inf (sign default '+')
1147   $x->bcmp($y);                 # compare numbers (undef,<0,=0,>0)
1148   $x->bacmp($y);                # compare absolutely (undef,<0,=0,>0)
1149   $x->sign();                   # return the sign, either +,- or NaN
1150
1151   # The following all modify their first argument:
1152
1153   # set
1154   $x->bzero();                  # set $i to 0
1155   $x->bnan();                   # set $i to NaN
1156
1157   $x->bneg();                   # negation
1158   $x->babs();                   # absolute value
1159   $x->bnorm();                  # normalize (no-op)
1160   $x->bnot();                   # two's complement (bit wise not)
1161   $x->binc();                   # increment x by 1
1162   $x->bdec();                   # decrement x by 1
1163
1164   $x->badd($y);                 # addition (add $y to $x)
1165   $x->bsub($y);                 # subtraction (subtract $y from $x)
1166   $x->bmul($y);                 # multiplication (multiply $x by $y)
1167   $x->bdiv($y);                 # divide, set $i to quotient
1168                                 # return (quo,rem) or quo if scalar
1169
1170   $x->bmod($y);                 # modulus
1171   $x->bpow($y);                 # power of arguments (a**b)
1172   $x->blsft($y);                # left shift
1173   $x->brsft($y);                # right shift
1174                                 # return (quo,rem) or quo if scalar
1175
1176   $x->band($y);                 # bit-wise and
1177   $x->bior($y);                 # bit-wise inclusive or
1178   $x->bxor($y);                 # bit-wise exclusive or
1179   $x->bnot();                   # bit-wise not (two's complement)
1180
1181   $x->bround($N);               # accuracy: preserver $N digits
1182   $x->bfround($N);              # precision: round to the $Nth digit
1183
1184   # The following do not modify their arguments:
1185
1186   bgcd(@values);                # greatest common divisor
1187   blcm(@values);                # lowest common multiplicator
1188
1189   $x->bstr();                   # return string
1190   $x->bsstr();                  # return string in scientific notation
1191
1192   $x->exponent();               # return exponent as BigInt
1193   $x->mantissa();               # return mantissa as BigInt
1194   $x->parts();                  # return (mantissa,exponent) as BigInt
1195
1196   $x->length();                 # number of digits (w/o sign and '.')
1197   ($l,$f) = $x->length();       # number of digits, and length of fraction
1198
1199 =head1 DESCRIPTION
1200
1201 All operators (inlcuding basic math operations) are overloaded if you
1202 declare your big floating point numbers as
1203
1204   $i = new Math::BigFloat '12_3.456_789_123_456_789E-2';
1205
1206 Operations with overloaded operators preserve the arguments, which is
1207 exactly what you expect.
1208
1209 =head2 Canonical notation
1210
1211 Input to these routines are either BigFloat objects, or strings of the
1212 following four forms:
1213
1214 =over 2
1215
1216 =item *
1217
1218 C</^[+-]\d+$/>
1219
1220 =item *
1221
1222 C</^[+-]\d+\.\d*$/>
1223
1224 =item *
1225
1226 C</^[+-]\d+E[+-]?\d+$/>
1227
1228 =item *
1229
1230 C</^[+-]\d*\.\d+E[+-]?\d+$/>
1231
1232 =back
1233
1234 all with optional leading and trailing zeros and/or spaces. Additonally,
1235 numbers are allowed to have an underscore between any two digits.
1236
1237 Empty strings as well as other illegal numbers results in 'NaN'.
1238
1239 bnorm() on a BigFloat object is now effectively a no-op, since the numbers
1240 are always stored in normalized form. On a string, it creates a BigFloat
1241 object.
1242
1243 =head2 Output
1244
1245 Output values are BigFloat objects (normalized), except for bstr() and bsstr().
1246
1247 The string output will always have leading and trailing zeros stripped and drop
1248 a plus sign. C<bstr()> will give you always the form with a decimal point,
1249 while C<bsstr()> (for scientific) gives you the scientific notation.
1250
1251         Input                   bstr()          bsstr()
1252         '-0'                    '0'             '0E1'
1253         '  -123 123 123'        '-123123123'    '-123123123E0'
1254         '00.0123'               '0.0123'        '123E-4'
1255         '123.45E-2'             '1.2345'        '12345E-4'
1256         '10E+3'                 '10000'         '1E4'
1257
1258 Some routines (C<is_odd()>, C<is_even()>, C<is_zero()>, C<is_one()>,
1259 C<is_nan()>) return true or false, while others (C<bcmp()>, C<bacmp()>)
1260 return either undef, <0, 0 or >0 and are suited for sort.
1261
1262 Actual math is done by using BigInts to represent the mantissa and exponent.
1263 The sign C</^[+-]$/> is stored separately. The string 'NaN' is used to
1264 represent the result when input arguments are not numbers, as well as
1265 the result of dividing by zero.
1266
1267 =head2 C<mantissa()>, C<exponent()> and C<parts()>
1268
1269 C<mantissa()> and C<exponent()> return the said parts of the BigFloat
1270 as BigInts such that:
1271
1272         $m = $x->mantissa();
1273         $e = $x->exponent();
1274         $y = $m * ( 10 ** $e );
1275         print "ok\n" if $x == $y;
1276
1277 C<< ($m,$e) = $x->parts(); >> is just a shortcut giving you both of them.
1278
1279 A zero is represented and returned as C<0E1>, B<not> C<0E0> (after Knuth).
1280
1281 Currently the mantissa is reduced as much as possible, favouring higher
1282 exponents over lower ones (e.g. returning 1e7 instead of 10e6 or 10000000e0).
1283 This might change in the future, so do not depend on it.
1284
1285 =head2 Accuracy vs. Precision
1286
1287 See also: L<Rounding|Rounding>.
1288
1289 Math::BigFloat supports both precision and accuracy. (here should follow
1290 a short description of both).
1291
1292 Precision: digits after the '.', laber, schwad
1293 Accuracy: Significant digits blah blah
1294
1295 Since things like sqrt(2) or 1/3 must presented with a limited precision lest
1296 a operation consumes all resources, each operation produces no more than
1297 C<Math::BigFloat::precision()> digits.
1298
1299 In case the result of one operation has more precision than specified,
1300 it is rounded. The rounding mode taken is either the default mode, or the one
1301 supplied to the operation after the I<scale>:
1302
1303         $x = Math::BigFloat->new(2);
1304         Math::BigFloat::precision(5);           # 5 digits max
1305         $y = $x->copy()->bdiv(3);               # will give 0.66666
1306         $y = $x->copy()->bdiv(3,6);             # will give 0.666666
1307         $y = $x->copy()->bdiv(3,6,'odd');       # will give 0.666667
1308         Math::BigFloat::round_mode('zero');
1309         $y = $x->copy()->bdiv(3,6);             # will give 0.666666
1310
1311 =head2 Rounding
1312
1313 =over 2
1314
1315 =item ffround ( +$scale )
1316
1317 Rounds to the $scale'th place left from the '.', counting from the dot.
1318 The first digit is numbered 1.
1319
1320 =item ffround ( -$scale )
1321
1322 Rounds to the $scale'th place right from the '.', counting from the dot.
1323
1324 =item ffround ( 0 )
1325
1326 Rounds to an integer.
1327
1328 =item fround  ( +$scale )
1329
1330 Preserves accuracy to $scale digits from the left (aka significant digits)
1331 and pads the rest with zeros. If the number is between 1 and -1, the
1332 significant digits count from the first non-zero after the '.'
1333
1334 =item fround  ( -$scale ) and fround ( 0 )
1335
1336 These are effetively no-ops.
1337
1338 =back
1339
1340 All rounding functions take as a second parameter a rounding mode from one of
1341 the following: 'even', 'odd', '+inf', '-inf', 'zero' or 'trunc'.
1342
1343 The default rounding mode is 'even'. By using
1344 C<< Math::BigFloat::round_mode($rnd_mode); >> you can get and set the default
1345 mode for subsequent rounding. The usage of C<$Math::BigFloat::$rnd_mode> is
1346 no longer supported.
1347 The second parameter to the round functions then overrides the default
1348 temporarily.
1349
1350 The C<< as_number() >> function returns a BigInt from a Math::BigFloat. It uses
1351 'trunc' as rounding mode to make it equivalent to:
1352
1353         $x = 2.5;
1354         $y = int($x) + 2;
1355
1356 You can override this by passing the desired rounding mode as parameter to
1357 C<as_number()>:
1358
1359         $x = Math::BigFloat->new(2.5);
1360         $y = $x->as_number('odd');      # $y = 3
1361
1362 =head1 EXAMPLES
1363
1364   # not ready yet
1365
1366 =head1 Autocreating constants
1367
1368 After C<use Math::BigFloat ':constant'> all the floating point constants
1369 in the given scope are converted to C<Math::BigFloat>. This conversion
1370 happens at compile time.
1371
1372 In particular
1373
1374   perl -MMath::BigFloat=:constant -e 'print 2E-100,"\n"'
1375
1376 prints the value of C<2E-100>.  Note that without conversion of
1377 constants the expression 2E-100 will be calculated as normal floating point
1378 number.
1379
1380 =head1 BUGS
1381
1382 =over 2
1383
1384 =item *
1385
1386 The following does not work yet:
1387
1388         $m = $x->mantissa();
1389         $e = $x->exponent();
1390         $y = $m * ( 10 ** $e );
1391         print "ok\n" if $x == $y;
1392
1393 =item *
1394
1395 There is no fmod() function yet.
1396
1397 =back
1398
1399 =head1 CAVEAT
1400
1401 =over 1
1402
1403 =item stringify, bstr()
1404
1405 Both stringify and bstr() now drop the leading '+'. The old code would return
1406 '+1.23', the new returns '1.23'. See the documentation in L<Math::BigInt> for
1407 reasoning and details.
1408
1409 =item bdiv
1410
1411 The following will probably not do what you expect:
1412
1413         print $c->bdiv(123.456),"\n";
1414
1415 It prints both quotient and reminder since print works in list context. Also,
1416 bdiv() will modify $c, so be carefull. You probably want to use
1417
1418         print $c / 123.456,"\n";
1419         print scalar $c->bdiv(123.456),"\n";  # or if you want to modify $c
1420
1421 instead.
1422
1423 =item Modifying and =
1424
1425 Beware of:
1426
1427         $x = Math::BigFloat->new(5);
1428         $y = $x;
1429
1430 It will not do what you think, e.g. making a copy of $x. Instead it just makes
1431 a second reference to the B<same> object and stores it in $y. Thus anything
1432 that modifies $x will modify $y, and vice versa.
1433
1434         $x->bmul(2);
1435         print "$x, $y\n";       # prints '10, 10'
1436
1437 If you want a true copy of $x, use:
1438
1439         $y = $x->copy();
1440
1441 See also the documentation in L<overload> regarding C<=>.
1442
1443 =item bpow
1444
1445 C<bpow()> now modifies the first argument, unlike the old code which left
1446 it alone and only returned the result. This is to be consistent with
1447 C<badd()> etc. The first will modify $x, the second one won't:
1448
1449         print bpow($x,$i),"\n";         # modify $x
1450         print $x->bpow($i),"\n";        # ditto
1451         print $x ** $i,"\n";            # leave $x alone
1452
1453 =back
1454
1455 =head1 LICENSE
1456
1457 This program is free software; you may redistribute it and/or modify it under
1458 the same terms as Perl itself.
1459
1460 =head1 AUTHORS
1461
1462 Mark Biggar, overloaded interface by Ilya Zakharevich.
1463 Completely rewritten by Tels http://bloodgate.com in 2001.
1464
1465 =cut