projects / p5sagit/p5-mst-13.2.git / blame
? search:
summary | shortlog | log | commit | commitdiff | tree
blob | history | HEAD
perhaps pthread_attr_t really needs to be set for the old api
[p5sagit/p5-mst-13.2.git] / lib / Math / BigFloat.pm
CommitLineData
58cde26e 1#!/usr/bin/perl -w
2
58cde26e 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
0716bf9b 9# _f: flags, used to signal MBI not to touch our private parts
58cde26e 10# _cow: Copy-On-Write (NRY)
11
a0d0e21e 12package Math::BigFloat;
13
bd05a461 14$VERSION = '1.21';
58cde26e 15require 5.005;
16use Exporter;
0716bf9b 17use Math::BigInt qw/objectify/;
58cde26e 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
574bacfe 25 bacmp bstr binc bdec binf
0716bf9b 26 is_odd is_even is_nan is_inf is_positive is_negative
58cde26e 27 is_zero is_one sign
28 );
a0d0e21e 29
58cde26e 30#@EXPORT = qw( );
31use strict;
32use vars qw/$AUTOLOAD $accuracy $precision $div_scale $rnd_mode/;
33my $class = "Math::BigFloat";
a0d0e21e 34
a5f75d66 35use overload
bd05a461 36'<=>' => sub { $_[2] ?
37 ref($_[0])->bcmp($_[1],$_[0]) :
38 ref($_[0])->bcmp($_[0],$_[1])},
0716bf9b 39'int' => sub { $_[0]->as_number() }, # 'trunc' to bigint
a5f75d66 40;
a0d0e21e 41
0716bf9b 42##############################################################################
43# global constants, flags and accessory
44
45use constant MB_NEVER_ROUND => 0x0001;
46
58cde26e 47# are NaNs ok?
48my $NaNOK=1;
58cde26e 49# constant for easier life
50my $nan = 'NaN';
58cde26e 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
574bacfe 58# in case we call SUPER::->foo() and this wants to call modify()
59# sub modify () { 0; }
60
58cde26e 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]||''}; }
a0d0e21e 69}
0e8b9368 70
58cde26e 71##############################################################################
72# constructors
a0d0e21e 73
58cde26e 74sub 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"
a0d0e21e 80
58cde26e 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;
a0d0e21e 86
58cde26e 87 my $round = shift; $round = 0 if !defined $round; # no rounding as default
88 my $self = {}; bless $self, $class;
b22b3e31 89 # shortcut for bigints and its subclasses
0716bf9b 90 if ((ref($wanted)) && (ref($wanted) ne $class))
58cde26e 91 {
0716bf9b 92 $self->{_m} = $wanted->as_number(); # get us a bigint copy
58cde26e 93 $self->{_e} = Math::BigInt->new(0);
94 $self->{_m}->babs();
95 $self->{sign} = $wanted->sign();
0716bf9b 96 return $self->bnorm();
58cde26e 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;
0716bf9b 105 return $self->bnorm();
58cde26e 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 }
0716bf9b 125 #print "$wanted => $self->{sign} $self->{value}\n";
58cde26e 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 }
a0d0e21e 132
58cde26e 133sub 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))
288d023a 139 {
58cde26e 140 my $c = $self; $self = {}; bless $self, $c;
a0d0e21e 141 }
574bacfe 142 $self->{_m} = Math::BigInt->bzero();
143 $self->{_e} = Math::BigInt->bzero();
58cde26e 144 $self->{sign} = $nan;
58cde26e 145 return $self;
146 }
a0d0e21e 147
58cde26e 148sub 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 '-';
a0d0e21e 153
58cde26e 154 $self = $class if !defined $self;
155 if (!ref($self))
156 {
157 my $c = $self; $self = {}; bless $self, $c;
158 }
574bacfe 159 $self->{_m} = Math::BigInt->bzero();
160 $self->{_e} = Math::BigInt->bzero();
58cde26e 161 $self->{sign} = $sign.'inf';
58cde26e 162 return $self;
163 }
a0d0e21e 164
574bacfe 165sub 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
58cde26e 182sub 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 }
574bacfe 191 $self->{_m} = Math::BigInt->bzero();
192 $self->{_e} = Math::BigInt->bone();
58cde26e 193 $self->{sign} = '+';
58cde26e 194 return $self;
195 }
196
197##############################################################################
198# string conversation
199
200sub 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")
58cde26e 205 my ($self,$x) = objectify(1,@_);
206
574bacfe 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} !~ /^[+-]$/)
58cde26e 210 {
574bacfe 211 return $x->{sign} unless $x->{sign} eq '+inf'; # -inf, NaN
212 return 'inf'; # +inf
58cde26e 213 }
214
574bacfe 215 my $es = '0'; my $len = 1; my $cad = 0; my $dot = '.';
216
217 my $not_zero = !$x->is_zero();
218 if ($not_zero)
58cde26e 219 {
574bacfe 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
58cde26e 227 {
574bacfe 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 }
82cf049f 248 }
574bacfe 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;
82cf049f 259 }
574bacfe 260 elsif ($x->{_p} || 0 < 0)
58cde26e 261 {
574bacfe 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;
58cde26e 266 }
58cde26e 267 return $es;
82cf049f 268 }
f216259d 269
58cde26e 270sub 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")
58cde26e 275 my ($self,$x) = objectify(1,@_);
a0d0e21e 276
574bacfe 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 }
58cde26e 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
289sub numify
290 {
291 # Make a number from a BigFloat object
574bacfe 292 # simple return string and let Perl's atoi()/atof() handle the rest
58cde26e 293 my ($self,$x) = objectify(1,@_);
294 return $x->bsstr();
295 }
a0d0e21e 296
58cde26e 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# }
574bacfe 309
310# tels 2001-08-04
311# todo: this must be overwritten and return NaN for non-integer values
312# band(), bior(), bxor(), too
58cde26e 313#sub bnot
314# {
315# $class->SUPER::bnot($class,@_);
316# }
317
318sub 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,@_);
58cde26e 323
0716bf9b 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
574bacfe 336 return 1 if $x->{sign} eq '+' && $y->{sign} eq '-'; # does also 0 <=> -y
58cde26e 337 return -1 if $x->{sign} eq '-' && $y->{sign} eq '+'; # does also -x <=> 0
338
574bacfe 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
58cde26e 345
346 # adjust so that exponents are equal
bd05a461 347 my $lxm = $x->{_m}->length();
348 my $lym = $y->{_m}->length();
349 my $lx = $lxm + $x->{_e};
350 my $ly = $lym + $y->{_e};
58cde26e 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";
bd05a461 354 return $l <=> 0 if $l != 0;
58cde26e 355
bd05a461 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);
58cde26e 370 $rc = -$rc if $x->{sign} eq '-'; # -124 < -123
bd05a461 371 return $rc <=> 0;
58cde26e 372 }
373
374sub 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 }
a0d0e21e 396
58cde26e 397sub badd
398 {
399 # add second arg (BFLOAT or string) to first (BFLOAT) (modifies first)
400 # return result as BFLOAT
58cde26e 401 my ($self,$x,$y,$a,$p,$r) = objectify(2,@_);
402
574bacfe 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
58cde26e 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);
a0d0e21e 430 }
58cde26e 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 {
0716bf9b 438 # print "e < 0\n";
58cde26e 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 {
0716bf9b 450 # print "e > 0\n";
58cde26e 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} = '+';
0716bf9b 466 #$x->bnorm(); # delete trailing zeros
58cde26e 467 return $x->round($a,$p,$r,$y);
468 }
469
470sub bsub
471 {
0716bf9b 472 # (BigFloat or num_str, BigFloat or num_str) return BigFloat
58cde26e 473 # subtract second arg from first, modify first
474 my ($self,$x,$y) = objectify(2,@_);
a0d0e21e 475
58cde26e 476 $x->badd($y->bneg()); # badd does not leave internal zeros
477 $y->bneg(); # refix y, assumes no one reads $y in between
0716bf9b 478 return $x; # badd() already normalized and rounded
58cde26e 479 }
480
481sub binc
482 {
483 # increment arg by one
484 my ($self,$x,$a,$p,$r) = objectify(1,@_);
58cde26e 485 $x->badd($self->_one())->round($a,$p,$r);
486 }
487
488sub bdec
489 {
490 # decrement arg by one
491 my ($self,$x,$a,$p,$r) = objectify(1,@_);
58cde26e 492 $x->badd($self->_one('-'))->round($a,$p,$r);
493 }
494
495sub 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
58cde26e 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
507sub 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
58cde26e 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
519sub 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;
574bacfe 523
524 return 1 if $x->{sign} eq '+' && $x->{_m}->is_zero();
525 return 0;
58cde26e 526 }
527
528sub 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
538sub is_odd
539 {
b22b3e31 540 # return true if arg (BINT or num_str) is odd or false if even
58cde26e 541 my $x = shift; $x = $class->new($x) unless ref $x;
542 #my ($self,$x) = objectify(1,@_);
0716bf9b 543
544 return 0 if $x->{sign} !~ /^[+-]$/; # NaN & +-inf aren't
545 return ($x->{_e}->is_zero() && $x->{_m}->is_odd());
58cde26e 546 }
547
548sub is_even
549 {
b22b3e31 550 # return true if arg (BINT or num_str) is even or false if odd
58cde26e 551 my $x = shift; $x = $class->new($x) unless ref $x;
552 #my ($self,$x) = objectify(1,@_);
0716bf9b 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
58cde26e 557 }
558
559sub 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,@_);
58cde26e 564
0716bf9b 565 # print "mbf bmul $x->{_m}e$x->{_e} $y->{_m}e$y->{_e}\n";
58cde26e 566 return $x->bnan() if (($x->{sign} eq $nan) || ($y->{sign} eq $nan));
567
574bacfe 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
58cde26e 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";
0716bf9b 589 $x->bnorm();
58cde26e 590 return $x->round($a,$p,$r,$y);
591 }
592
593sub 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
574bacfe 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
58cde26e 604 return wantarray ? ($x->bnan(),bnan()) : $x->bnan()
574bacfe 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());
0716bf9b 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";
58cde26e 616 # we need to limit the accuracy to protect against overflow
0716bf9b 617 my ($scale) = $x->_scale_a($accuracy,$rnd_mode,$a,$r); # ignore $p
574bacfe 618 my $fallback = 0;
58cde26e 619 if (!defined $scale)
620 {
0716bf9b 621 # simulate old behaviour
622 $scale = $div_scale+1; # one more for proper riund
574bacfe 623 $a = $div_scale; # and round to it
624 $fallback = 1; # to clear a/p afterwards
a0d0e21e 625 }
0716bf9b 626 my $lx = $x->{_m}->length(); my $ly = $y->{_m}->length();
58cde26e 627 $scale = $lx if $lx > $scale;
58cde26e 628 $scale = $ly if $ly > $scale;
629 #print "scale $scale $lx $ly\n";
0716bf9b 630 my $diff = $ly - $lx;
631 $scale += $diff if $diff > 0; # if lx << ly, but not if ly << lx!
a0d0e21e 632
58cde26e 633 return wantarray ? ($x,$self->bzero()) : $x if $x->is_zero();
634
635 $x->{sign} = $x->{sign} ne $y->sign() ? '-' : '+';
a0d0e21e 636
58cde26e 637 # check for / +-1 ( +/- 1E0)
638 if ($y->is_one())
639 {
640 return wantarray ? ($x,$self->bzero()) : $x;
a0d0e21e 641 }
a5f75d66 642
58cde26e 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);
0716bf9b 647 #print "m: $x->{_m} $y->{_m}\n";
58cde26e 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
0716bf9b 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
574bacfe 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 }
0716bf9b 664
58cde26e 665 if (wantarray)
666 {
667 my $rem = $x->copy();
668 $rem->bmod($y,$a,$p,$r);
574bacfe 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 }
0716bf9b 675 return ($x,$rem);
58cde26e 676 }
677 return $x;
678 }
a0d0e21e 679
58cde26e 680sub 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,@_);
a0d0e21e 684
58cde26e 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
692sub bsqrt
693 {
0716bf9b 694 # calculate square root; this should probably
695 # use a different test to see whether the accuracy we want is...
58cde26e 696 my ($self,$x,$a,$p,$r) = objectify(1,@_);
697
0716bf9b 698 return $x->bnan() if $x->{sign} eq 'NaN' || $x->{sign} =~ /^-/; # <0, NaN
699 return $x if $x->{sign} eq '+inf'; # +inf
58cde26e 700 return $x if $x->is_zero() || $x == 1;
701
0716bf9b 702 # we need to limit the accuracy to protect against overflow
703 my ($scale) = $x->_scale_a($accuracy,$rnd_mode,$a,$r); # ignore $p
574bacfe 704 my $fallback = 0;
0716bf9b 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
574bacfe 710 $fallback = 1; # to clear a/p afterwards
0716bf9b 711 }
712 my $lx = $x->{_m}->length();
713 $scale = $lx if $scale < $lx;
58cde26e 714 my $e = Math::BigFloat->new("1E-$scale"); # make test variable
715 return $x->bnan() if $e->sign() eq 'NaN';
716
58cde26e 717 # start with some reasonable guess
0716bf9b 718 #$x *= 10 ** ($len - $org->{_e}); $x /= 2; # !?!?
574bacfe 719 $lx = $lx+$x->{_e};
0716bf9b 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";
58cde26e 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 {
58cde26e 732 return $x->bnan() if $gs->is_zero();
0716bf9b 733 $r = $y->copy(); $r->bdiv($gs,$scale);
58cde26e 734 $x = ($r + $gs);
0716bf9b 735 $x->bdiv($two,$scale);
58cde26e 736 $diff = $x->copy()->bsub($gs)->babs();
58cde26e 737 $gs = $x->copy();
a0d0e21e 738 }
0716bf9b 739 $x->round($a,$p,$r);
574bacfe 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;
58cde26e 747 }
748
749sub 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
0716bf9b 757 return $x if $x->{sign} =~ /^[+-]inf$/;
58cde26e 758 return $x->bnan() if $x->{sign} eq $nan || $y->{sign} eq $nan;
574bacfe 759 return $x->bone() if $y->is_zero();
58cde26e 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
0716bf9b 765 return $y1->is_odd() ? $x : $x->babs(1);
288d023a 766 }
58cde26e 767 return $x if $x->is_zero() && $y->{sign} eq '+'; # 0**y => 0 (if not y <= 0)
574bacfe 768 # 0 ** -y => 1 / (0 ** y) => / 0! (1 / 0 => +inf)
769 return $x->binf() if $x->is_zero() && $y->{sign} eq '-';
58cde26e 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!
0716bf9b 780 my $z = $x->copy(); $x->bzero()->binc();
58cde26e 781 return $x->bdiv($z,$a,$p,$r); # round in one go (might ignore y's A!)
a0d0e21e 782 }
58cde26e 783 return $x->round($a,$p,$r,$y);
784 }
785
786###############################################################################
787# rounding functions
788
789sub 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;
a0d0e21e 795
58cde26e 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
574bacfe 801 # never round a 0, +-inf, NaN
802 return $x if $x->{sign} !~ /^[+-]$/ || $x->is_zero();
58cde26e 803 # print "MBF bfround $x to scale $scale mode $mode\n";
58cde26e 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 {
574bacfe 831 return $x->bzero();
58cde26e 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";
a0d0e21e 851 }
58cde26e 852 else
853 {
854 # 123 => 100 means length(123) = 3 - $scale (2) => 1
a5f75d66 855
58cde26e 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
574bacfe 861 return $x->bzero();
58cde26e 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 }
a0d0e21e 881 }
574bacfe 882 # print "using $scale for $x->{_m} with '$mode'\n";
883 # pass sign to bround for rounding modes '+inf' and '-inf'
58cde26e 884 $x->{_m}->{sign} = $x->{sign};
885 $x->{_m}->bround($scale,$mode);
886 $x->{_m}->{sign} = '+'; # fix sign back
887 $x->bnorm();
888 }
889
890sub 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);
574bacfe 902 # never round a 0, +-inf, NaN
903 return $x if $x->{sign} !~ /^[+-]$/ || $x->is_zero();
58cde26e 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;
f216259d 914
58cde26e 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
574bacfe 919 $x->bnorm(); # del trailing zeros gen. by bround()
58cde26e 920 }
921
922sub 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 '-';
f216259d 937 }
58cde26e 938 return $x->round($a,$p,$r);
939 }
288d023a 940
58cde26e 941sub 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 '+';
a0d0e21e 955 }
58cde26e 956 return $x->round($a,$p,$r);
957 }
958
959###############################################################################
a5f75d66 960
58cde26e 961sub DESTROY
962 {
963 # going trough AUTOLOAD for every DESTROY is costly, so avoid it by empty sub
964 }
965
966sub 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");
a0d0e21e 982 }
58cde26e 983 no strict 'refs';
984 my $bname = $name; $bname =~ s/^f/b/;
985 *{$class."\:\:$name"} = \&$bname;
986 &$bname; # uses @_
987 }
988
989sub 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
999sub 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
1012sub 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
1027sub _one
1028 {
1029 # internal speedup, set argument to 1, or create a +/- 1
0716bf9b 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);
58cde26e 1034 $x->{sign} = shift || '+';
1035 return $x;
1036 }
1037
1038sub 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
1058sub 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
0716bf9b 1064 return $x if $x->{sign} !~ /^[+-]$/; # inf, nan etc
58cde26e 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
574bacfe 1072 $x->{sign} = '+', $x->{_e}->bzero()->binc() if $x->{_m}->is_zero();
0716bf9b 1073 $x->{_m}->{_f} = MB_NEVER_ROUND;
1074 $x->{_e}->{_f} = MB_NEVER_ROUND;
1075 return $x; # MBI bnorm is no-op
58cde26e 1076 }
1077
1078##############################################################################
1079# internal calculation routines
1080
1081sub 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 }
0716bf9b 1093 $z = $x->{_m}->copy();
58cde26e 1094 if ($x->{_e} < 0)
1095 {
0716bf9b 1096 $z->brsft(-$x->{_e},10);
1097 }
1098 else
1099 {
1100 $z->blsft($x->{_e},10);
58cde26e 1101 }
58cde26e 1102 $z->{sign} = $x->{sign};
1103 return $z;
1104 }
1105
1106sub 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 }
a0d0e21e 1120
11211;
a5f75d66 1122__END__
1123
1124=head1 NAME
1125
58cde26e 1126Math::BigFloat - Arbitrary size floating point math package
a5f75d66 1127
1128=head1 SYNOPSIS
1129
a2008d6d 1130 use Math::BigFloat;
58cde26e 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
0716bf9b 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 '+')
58cde26e 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
a5f75d66 1198
1199=head1 DESCRIPTION
1200
58cde26e 1201All operators (inlcuding basic math operations) are overloaded if you
1202declare your big floating point numbers as
a5f75d66 1203
58cde26e 1204 $i = new Math::BigFloat '12_3.456_789_123_456_789E-2';
1205
1206Operations with overloaded operators preserve the arguments, which is
1207exactly what you expect.
1208
1209=head2 Canonical notation
1210
1211Input to these routines are either BigFloat objects, or strings of the
1212following four forms:
a5f75d66 1213
1214=over 2
1215
58cde26e 1216=item *
1217
1218C</^[+-]\d+$/>
a5f75d66 1219
58cde26e 1220=item *
a5f75d66 1221
58cde26e 1222C</^[+-]\d+\.\d*$/>
a5f75d66 1223
58cde26e 1224=item *
a5f75d66 1225
58cde26e 1226C</^[+-]\d+E[+-]?\d+$/>
a5f75d66 1227
58cde26e 1228=item *
a5f75d66 1229
58cde26e 1230C</^[+-]\d*\.\d+E[+-]?\d+$/>
5d7098d5 1231
58cde26e 1232=back
1233
1234all with optional leading and trailing zeros and/or spaces. Additonally,
1235numbers are allowed to have an underscore between any two digits.
1236
1237Empty strings as well as other illegal numbers results in 'NaN'.
1238
1239bnorm() on a BigFloat object is now effectively a no-op, since the numbers
1240are always stored in normalized form. On a string, it creates a BigFloat
1241object.
1242
1243=head2 Output
1244
1245Output values are BigFloat objects (normalized), except for bstr() and bsstr().
1246
1247The string output will always have leading and trailing zeros stripped and drop
1248a plus sign. C<bstr()> will give you always the form with a decimal point,
1249while 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
1258Some routines (C<is_odd()>, C<is_even()>, C<is_zero()>, C<is_one()>,
1259C<is_nan()>) return true or false, while others (C<bcmp()>, C<bacmp()>)
1260return either undef, <0, 0 or >0 and are suited for sort.
1261
1262Actual math is done by using BigInts to represent the mantissa and exponent.
1263The sign C</^[+-]$/> is stored separately. The string 'NaN' is used to
1264represent the result when input arguments are not numbers, as well as
1265the result of dividing by zero.
1266
1267=head2 C<mantissa()>, C<exponent()> and C<parts()>
1268
1269C<mantissa()> and C<exponent()> return the said parts of the BigFloat
1270as 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
1277C<< ($m,$e) = $x->parts(); >> is just a shortcut giving you both of them.
1278
1279A zero is represented and returned as C<0E1>, B<not> C<0E0> (after Knuth).
1280
1281Currently the mantissa is reduced as much as possible, favouring higher
1282exponents over lower ones (e.g. returning 1e7 instead of 10e6 or 10000000e0).
1283This might change in the future, so do not depend on it.
1284
1285=head2 Accuracy vs. Precision
1286
1287See also: L<Rounding|Rounding>.
1288
1289Math::BigFloat supports both precision and accuracy. (here should follow
1290a short description of both).
5d7098d5 1291
58cde26e 1292Precision: digits after the '.', laber, schwad
1293Accuracy: Significant digits blah blah
5d7098d5 1294
58cde26e 1295Since things like sqrt(2) or 1/3 must presented with a limited precision lest
1296a operation consumes all resources, each operation produces no more than
1297C<Math::BigFloat::precision()> digits.
1298
1299In case the result of one operation has more precision than specified,
1300it is rounded. The rounding mode taken is either the default mode, or the one
1301supplied 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
5dc6f178 1315=item ffround ( +$scale )
58cde26e 1316
0716bf9b 1317Rounds to the $scale'th place left from the '.', counting from the dot.
1318The first digit is numbered 1.
58cde26e 1319
5dc6f178 1320=item ffround ( -$scale )
58cde26e 1321
0716bf9b 1322Rounds to the $scale'th place right from the '.', counting from the dot.
58cde26e 1323
5dc6f178 1324=item ffround ( 0 )
1325
0716bf9b 1326Rounds to an integer.
5dc6f178 1327
1328=item fround ( +$scale )
1329
0716bf9b 1330Preserves accuracy to $scale digits from the left (aka significant digits)
1331and pads the rest with zeros. If the number is between 1 and -1, the
1332significant digits count from the first non-zero after the '.'
5dc6f178 1333
1334=item fround ( -$scale ) and fround ( 0 )
1335
0716bf9b 1336These are effetively no-ops.
5d7098d5 1337
a5f75d66 1338=back
1339
0716bf9b 1340All rounding functions take as a second parameter a rounding mode from one of
1341the following: 'even', 'odd', '+inf', '-inf', 'zero' or 'trunc'.
58cde26e 1342
1343The default rounding mode is 'even'. By using
0716bf9b 1344C<< Math::BigFloat::round_mode($rnd_mode); >> you can get and set the default
1345mode for subsequent rounding. The usage of C<$Math::BigFloat::$rnd_mode> is
1346no longer supported.
b22b3e31 1347The second parameter to the round functions then overrides the default
0716bf9b 1348temporarily.
58cde26e 1349
1350The 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
1356You can override this by passing the desired rounding mode as parameter to
1357C<as_number()>:
1358
1359 $x = Math::BigFloat->new(2.5);
1360 $y = $x->as_number('odd'); # $y = 3
1361
1362=head1 EXAMPLES
1363
58cde26e 1364 # not ready yet
58cde26e 1365
1366=head1 Autocreating constants
1367
1368After C<use Math::BigFloat ':constant'> all the floating point constants
1369in the given scope are converted to C<Math::BigFloat>. This conversion
1370happens at compile time.
1371
1372In particular
1373
1374 perl -MMath::BigFloat=:constant -e 'print 2E-100,"\n"'
1375
1376prints the value of C<2E-100>. Note that without conversion of
1377constants the expression 2E-100 will be calculated as normal floating point
1378number.
1379
a5f75d66 1380=head1 BUGS
1381
58cde26e 1382=over 2
1383
1384=item *
1385
1386The 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
1395There is no fmod() function yet.
1396
1397=back
1398
1399=head1 CAVEAT
1400
1401=over 1
1402
1403=item stringify, bstr()
1404
1405Both 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
1407reasoning and details.
1408
1409=item bdiv
1410
1411The following will probably not do what you expect:
1412
1413 print $c->bdiv(123.456),"\n";
1414
1415It prints both quotient and reminder since print works in list context. Also,
1416bdiv() 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
1421instead.
1422
1423=item Modifying and =
1424
1425Beware of:
1426
1427 $x = Math::BigFloat->new(5);
1428 $y = $x;
1429
1430It will not do what you think, e.g. making a copy of $x. Instead it just makes
1431a second reference to the B<same> object and stores it in $y. Thus anything
1432that modifies $x will modify $y, and vice versa.
1433
1434 $x->bmul(2);
1435 print "$x, $y\n"; # prints '10, 10'
1436
1437If you want a true copy of $x, use:
1438
1439 $y = $x->copy();
1440
1441See also the documentation in L<overload> regarding C<=>.
1442
1443=item bpow
1444
1445C<bpow()> now modifies the first argument, unlike the old code which left
1446it alone and only returned the result. This is to be consistent with
1447C<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
a5f75d66 1456
58cde26e 1457This program is free software; you may redistribute it and/or modify it under
1458the same terms as Perl itself.
5d7098d5 1459
58cde26e 1460=head1 AUTHORS
5d7098d5 1461
58cde26e 1462Mark Biggar, overloaded interface by Ilya Zakharevich.
1463Completely rewritten by Tels http://bloodgate.com in 2001.
a5f75d66 1464
a5f75d66 1465=cut
Unnamed repository; edit this file to name it for gitweb.