projects / p5sagit/p5-mst-13.2.git / blame
? search:
summary | shortlog | log | commit | commitdiff | tree
blob | history | HEAD
Upgrade to Math::BigInt 1.51.
[p5sagit/p5-mst-13.2.git] / lib / Math / BigFloat.pm
CommitLineData
58cde26e 1# The following hash values are internally used:
2# _e: exponent (BigInt)
3# _m: mantissa (absolute BigInt)
4# sign: +,-,"NaN" if not a number
5# _a: accuracy
6# _p: precision
0716bf9b 7# _f: flags, used to signal MBI not to touch our private parts
58cde26e 8
a0d0e21e 9package Math::BigFloat;
10
b3abae2a 11$VERSION = '1.28';
58cde26e 12require 5.005;
13use Exporter;
0716bf9b 14use Math::BigInt qw/objectify/;
58cde26e 15@ISA = qw( Exporter Math::BigInt);
394e6ffb 16
58cde26e 17use strict;
027dc388 18use vars qw/$AUTOLOAD $accuracy $precision $div_scale $round_mode $rnd_mode/;
b3abae2a 19use vars qw/$upgrade $downgrade/;
58cde26e 20my $class = "Math::BigFloat";
a0d0e21e 21
a5f75d66 22use overload
bd05a461 23'<=>' => sub { $_[2] ?
24 ref($_[0])->bcmp($_[1],$_[0]) :
25 ref($_[0])->bcmp($_[0],$_[1])},
0716bf9b 26'int' => sub { $_[0]->as_number() }, # 'trunc' to bigint
b3abae2a 27'log' => sub { $_[0]->blog() },
a5f75d66 28;
a0d0e21e 29
0716bf9b 30##############################################################################
31# global constants, flags and accessory
32
33use constant MB_NEVER_ROUND => 0x0001;
34
58cde26e 35# are NaNs ok?
36my $NaNOK=1;
58cde26e 37# constant for easier life
38my $nan = 'NaN';
58cde26e 39
ee15d750 40# class constants, use Class->constant_name() to access
41$round_mode = 'even'; # one of 'even', 'odd', '+inf', '-inf', 'zero' or 'trunc'
42$accuracy = undef;
43$precision = undef;
44$div_scale = 40;
58cde26e 45
b3abae2a 46$upgrade = undef;
47$downgrade = undef;
48
027dc388 49##############################################################################
50# the old code had $rnd_mode, so we need to support it, too
51
52$rnd_mode = 'even';
53sub TIESCALAR { my ($class) = @_; bless \$round_mode, $class; }
54sub FETCH { return $round_mode; }
55sub STORE { $rnd_mode = $_[0]->round_mode($_[1]); }
56
57BEGIN { tie $rnd_mode, 'Math::BigFloat'; }
58
59##############################################################################
60
574bacfe 61# in case we call SUPER::->foo() and this wants to call modify()
62# sub modify () { 0; }
63
58cde26e 64{
ee15d750 65 # valid method aliases for AUTOLOAD
58cde26e 66 my %methods = map { $_ => 1 }
67 qw / fadd fsub fmul fdiv fround ffround fsqrt fmod fstr fsstr fpow fnorm
b3abae2a 68 fint facmp fcmp fzero fnan finf finc fdec flog ffac
61f5c3f5 69 fceil ffloor frsft flsft fone flog
ee15d750 70 /;
61f5c3f5 71 # valid method's that can be hand-ed up (for AUTOLOAD)
ee15d750 72 my %hand_ups = map { $_ => 1 }
73 qw / is_nan is_inf is_negative is_positive
394e6ffb 74 accuracy precision div_scale round_mode fneg fabs babs fnot
58cde26e 75 /;
76
ee15d750 77 sub method_alias { return exists $methods{$_[0]||''}; }
78 sub method_hand_up { return exists $hand_ups{$_[0]||''}; }
a0d0e21e 79}
0e8b9368 80
58cde26e 81##############################################################################
82# constructors
a0d0e21e 83
58cde26e 84sub new
85 {
86 # create a new BigFloat object from a string or another bigfloat object.
87 # _e: exponent
88 # _m: mantissa
89 # sign => sign (+/-), or "NaN"
a0d0e21e 90
61f5c3f5 91 my ($class,$wanted,@r) = @_;
b3abae2a 92
61f5c3f5 93 # avoid numify-calls by not using || on $wanted!
94 return $class->bzero() if !defined $wanted; # default to 0
95 return $wanted->copy() if UNIVERSAL::isa($wanted,'Math::BigFloat');
a0d0e21e 96
58cde26e 97 my $self = {}; bless $self, $class;
b22b3e31 98 # shortcut for bigints and its subclasses
0716bf9b 99 if ((ref($wanted)) && (ref($wanted) ne $class))
58cde26e 100 {
0716bf9b 101 $self->{_m} = $wanted->as_number(); # get us a bigint copy
027dc388 102 $self->{_e} = Math::BigInt->bzero();
58cde26e 103 $self->{_m}->babs();
104 $self->{sign} = $wanted->sign();
0716bf9b 105 return $self->bnorm();
58cde26e 106 }
107 # got string
108 # handle '+inf', '-inf' first
ee15d750 109 if ($wanted =~ /^[+-]?inf$/)
58cde26e 110 {
027dc388 111 $self->{_e} = Math::BigInt->bzero();
112 $self->{_m} = Math::BigInt->bzero();
58cde26e 113 $self->{sign} = $wanted;
ee15d750 114 $self->{sign} = '+inf' if $self->{sign} eq 'inf';
0716bf9b 115 return $self->bnorm();
58cde26e 116 }
117 #print "new string '$wanted'\n";
118 my ($mis,$miv,$mfv,$es,$ev) = Math::BigInt::_split(\$wanted);
119 if (!ref $mis)
120 {
121 die "$wanted is not a number initialized to $class" if !$NaNOK;
027dc388 122 $self->{_e} = Math::BigInt->bzero();
123 $self->{_m} = Math::BigInt->bzero();
58cde26e 124 $self->{sign} = $nan;
125 }
126 else
127 {
128 # make integer from mantissa by adjusting exp, then convert to bigint
61f5c3f5 129 # undef,undef to signal MBI that we don't need no bloody rounding
130 $self->{_e} = Math::BigInt->new("$$es$$ev",undef,undef); # exponent
131 $self->{_m} = Math::BigInt->new("$$miv$$mfv",undef,undef); # create mant.
58cde26e 132 # 3.123E0 = 3123E-3, and 3.123E-2 => 3123E-5
027dc388 133 $self->{_e} -= CORE::length($$mfv) if CORE::length($$mfv) != 0;
134 $self->{sign} = $$mis;
58cde26e 135 }
b3abae2a 136 # print "mbf new $self->{sign} $self->{_m} e $self->{_e}\n";
61f5c3f5 137 $self->bnorm()->round(@r); # first normalize, then round
58cde26e 138 }
a0d0e21e 139
58cde26e 140sub bnan
141 {
142 # create a bigfloat 'NaN', if given a BigFloat, set it to 'NaN'
143 my $self = shift;
144 $self = $class if !defined $self;
145 if (!ref($self))
288d023a 146 {
58cde26e 147 my $c = $self; $self = {}; bless $self, $c;
a0d0e21e 148 }
574bacfe 149 $self->{_m} = Math::BigInt->bzero();
150 $self->{_e} = Math::BigInt->bzero();
58cde26e 151 $self->{sign} = $nan;
61f5c3f5 152 $self->{_a} = undef; $self->{_p} = undef;
153 $self;
58cde26e 154 }
a0d0e21e 155
58cde26e 156sub binf
157 {
158 # create a bigfloat '+-inf', if given a BigFloat, set it to '+-inf'
159 my $self = shift;
160 my $sign = shift; $sign = '+' if !defined $sign || $sign ne '-';
a0d0e21e 161
58cde26e 162 $self = $class if !defined $self;
163 if (!ref($self))
164 {
165 my $c = $self; $self = {}; bless $self, $c;
166 }
574bacfe 167 $self->{_m} = Math::BigInt->bzero();
168 $self->{_e} = Math::BigInt->bzero();
58cde26e 169 $self->{sign} = $sign.'inf';
61f5c3f5 170 $self->{_a} = undef; $self->{_p} = undef;
171 $self;
58cde26e 172 }
a0d0e21e 173
574bacfe 174sub bone
175 {
176 # create a bigfloat '+-1', if given a BigFloat, set it to '+-1'
177 my $self = shift;
178 my $sign = shift; $sign = '+' if !defined $sign || $sign ne '-';
179
180 $self = $class if !defined $self;
181 if (!ref($self))
182 {
183 my $c = $self; $self = {}; bless $self, $c;
184 }
185 $self->{_m} = Math::BigInt->bone();
186 $self->{_e} = Math::BigInt->bzero();
187 $self->{sign} = $sign;
61f5c3f5 188 if (@_ > 0)
189 {
190 $self->{_a} = $_[0]
191 if (defined $self->{_a} && defined $_[0] && $_[0] > $self->{_a});
192 $self->{_p} = $_[1]
193 if (defined $self->{_p} && defined $_[1] && $_[1] < $self->{_p});
194 }
574bacfe 195 return $self;
196 }
197
58cde26e 198sub bzero
199 {
200 # create a bigfloat '+0', if given a BigFloat, set it to 0
201 my $self = shift;
202 $self = $class if !defined $self;
203 if (!ref($self))
204 {
205 my $c = $self; $self = {}; bless $self, $c;
206 }
574bacfe 207 $self->{_m} = Math::BigInt->bzero();
208 $self->{_e} = Math::BigInt->bone();
58cde26e 209 $self->{sign} = '+';
61f5c3f5 210 if (@_ > 0)
211 {
212 $self->{_a} = $_[0]
213 if (defined $self->{_a} && defined $_[0] && $_[0] > $self->{_a});
214 $self->{_p} = $_[1]
215 if (defined $self->{_p} && defined $_[1] && $_[1] < $self->{_p});
216 }
58cde26e 217 return $self;
218 }
219
220##############################################################################
221# string conversation
222
223sub bstr
224 {
225 # (ref to BFLOAT or num_str ) return num_str
226 # Convert number from internal format to (non-scientific) string format.
227 # internal format is always normalized (no leading zeros, "-0" => "+0")
ee15d750 228 my ($self,$x) = ref($_[0]) ? (ref($_[0]),$_[0]) : objectify(1,@_);
229 #my $x = shift; my $class = ref($x) || $x;
230 #$x = $class->new(shift) unless ref($x);
58cde26e 231
574bacfe 232 #die "Oups! e was $nan" if $x->{_e}->{sign} eq $nan;
233 #die "Oups! m was $nan" if $x->{_m}->{sign} eq $nan;
234 if ($x->{sign} !~ /^[+-]$/)
58cde26e 235 {
574bacfe 236 return $x->{sign} unless $x->{sign} eq '+inf'; # -inf, NaN
237 return 'inf'; # +inf
58cde26e 238 }
239
574bacfe 240 my $es = '0'; my $len = 1; my $cad = 0; my $dot = '.';
241
242 my $not_zero = !$x->is_zero();
243 if ($not_zero)
58cde26e 244 {
574bacfe 245 $es = $x->{_m}->bstr();
246 $len = CORE::length($es);
247 if (!$x->{_e}->is_zero())
58cde26e 248 {
574bacfe 249 if ($x->{_e}->sign() eq '-')
250 {
251 $dot = '';
252 if ($x->{_e} <= -$len)
253 {
254 # print "style: 0.xxxx\n";
255 my $r = $x->{_e}->copy(); $r->babs()->bsub( CORE::length($es) );
256 $es = '0.'. ('0' x $r) . $es; $cad = -($len+$r);
257 }
258 else
259 {
260 # print "insert '.' at $x->{_e} in '$es'\n";
261 substr($es,$x->{_e},0) = '.'; $cad = $x->{_e};
262 }
263 }
264 else
265 {
266 # expand with zeros
267 $es .= '0' x $x->{_e}; $len += $x->{_e}; $cad = 0;
268 }
82cf049f 269 }
574bacfe 270 } # if not zero
271 $es = $x->{sign}.$es if $x->{sign} eq '-';
272 # if set accuracy or precision, pad with zeros
273 if ((defined $x->{_a}) && ($not_zero))
274 {
275 # 123400 => 6, 0.1234 => 4, 0.001234 => 4
276 my $zeros = $x->{_a} - $cad; # cad == 0 => 12340
277 $zeros = $x->{_a} - $len if $cad != $len;
574bacfe 278 $es .= $dot.'0' x $zeros if $zeros > 0;
82cf049f 279 }
574bacfe 280 elsif ($x->{_p} || 0 < 0)
58cde26e 281 {
574bacfe 282 # 123400 => 6, 0.1234 => 4, 0.001234 => 6
283 my $zeros = -$x->{_p} + $cad;
574bacfe 284 $es .= $dot.'0' x $zeros if $zeros > 0;
58cde26e 285 }
58cde26e 286 return $es;
82cf049f 287 }
f216259d 288
58cde26e 289sub bsstr
290 {
291 # (ref to BFLOAT or num_str ) return num_str
292 # Convert number from internal format to scientific string format.
293 # internal format is always normalized (no leading zeros, "-0E0" => "+0E0")
ee15d750 294 my ($self,$x) = ref($_[0]) ? (ref($_[0]),$_[0]) : objectify(1,@_);
295 #my $x = shift; my $class = ref($x) || $x;
296 #$x = $class->new(shift) unless ref($x);
a0d0e21e 297
574bacfe 298 #die "Oups! e was $nan" if $x->{_e}->{sign} eq $nan;
299 #die "Oups! m was $nan" if $x->{_m}->{sign} eq $nan;
300 if ($x->{sign} !~ /^[+-]$/)
301 {
302 return $x->{sign} unless $x->{sign} eq '+inf'; # -inf, NaN
303 return 'inf'; # +inf
304 }
58cde26e 305 my $sign = $x->{_e}->{sign}; $sign = '' if $sign eq '-';
306 my $sep = 'e'.$sign;
307 return $x->{_m}->bstr().$sep.$x->{_e}->bstr();
308 }
309
310sub numify
311 {
312 # Make a number from a BigFloat object
574bacfe 313 # simple return string and let Perl's atoi()/atof() handle the rest
ee15d750 314 my ($self,$x) = ref($_[0]) ? (ref($_[0]),$_[0]) : objectify(1,@_);
58cde26e 315 return $x->bsstr();
316 }
a0d0e21e 317
58cde26e 318##############################################################################
319# public stuff (usually prefixed with "b")
320
574bacfe 321# tels 2001-08-04
322# todo: this must be overwritten and return NaN for non-integer values
323# band(), bior(), bxor(), too
58cde26e 324#sub bnot
325# {
326# $class->SUPER::bnot($class,@_);
327# }
328
329sub bcmp
330 {
331 # Compares 2 values. Returns one of undef, <0, =0, >0. (suitable for sort)
332 # (BFLOAT or num_str, BFLOAT or num_str) return cond_code
333 my ($self,$x,$y) = objectify(2,@_);
58cde26e 334
0716bf9b 335 if (($x->{sign} !~ /^[+-]$/) || ($y->{sign} !~ /^[+-]$/))
336 {
337 # handle +-inf and NaN
338 return undef if (($x->{sign} eq $nan) || ($y->{sign} eq $nan));
339 return 0 if ($x->{sign} eq $y->{sign}) && ($x->{sign} =~ /^[+-]inf$/);
340 return +1 if $x->{sign} eq '+inf';
341 return -1 if $x->{sign} eq '-inf';
342 return -1 if $y->{sign} eq '+inf';
b3abae2a 343 return +1;
0716bf9b 344 }
345
346 # check sign for speed first
574bacfe 347 return 1 if $x->{sign} eq '+' && $y->{sign} eq '-'; # does also 0 <=> -y
58cde26e 348 return -1 if $x->{sign} eq '-' && $y->{sign} eq '+'; # does also -x <=> 0
349
574bacfe 350 # shortcut
351 my $xz = $x->is_zero();
352 my $yz = $y->is_zero();
353 return 0 if $xz && $yz; # 0 <=> 0
354 return -1 if $xz && $y->{sign} eq '+'; # 0 <=> +y
355 return 1 if $yz && $x->{sign} eq '+'; # +x <=> 0
58cde26e 356
357 # adjust so that exponents are equal
bd05a461 358 my $lxm = $x->{_m}->length();
359 my $lym = $y->{_m}->length();
360 my $lx = $lxm + $x->{_e};
361 my $ly = $lym + $y->{_e};
b3abae2a 362 my $l = $lx - $ly; $l->bneg() if $x->{sign} eq '-';
bd05a461 363 return $l <=> 0 if $l != 0;
58cde26e 364
bd05a461 365 # lengths (corrected by exponent) are equal
366 # so make mantissa euqal length by padding with zero (shift left)
367 my $diff = $lxm - $lym;
368 my $xm = $x->{_m}; # not yet copy it
369 my $ym = $y->{_m};
370 if ($diff > 0)
371 {
372 $ym = $y->{_m}->copy()->blsft($diff,10);
373 }
374 elsif ($diff < 0)
375 {
376 $xm = $x->{_m}->copy()->blsft(-$diff,10);
377 }
378 my $rc = $xm->bcmp($ym);
58cde26e 379 $rc = -$rc if $x->{sign} eq '-'; # -124 < -123
b3abae2a 380 $rc <=> 0;
58cde26e 381 }
382
383sub bacmp
384 {
385 # Compares 2 values, ignoring their signs.
386 # Returns one of undef, <0, =0, >0. (suitable for sort)
387 # (BFLOAT or num_str, BFLOAT or num_str) return cond_code
388 my ($self,$x,$y) = objectify(2,@_);
ee15d750 389
390 # handle +-inf and NaN's
abcfbf51 391 if ($x->{sign} !~ /^[+-]$/ || $y->{sign} !~ /^[+-]$/)
ee15d750 392 {
393 return undef if (($x->{sign} eq $nan) || ($y->{sign} eq $nan));
394 return 0 if ($x->is_inf() && $y->is_inf());
395 return 1 if ($x->is_inf() && !$y->is_inf());
b3abae2a 396 return -1;
ee15d750 397 }
398
399 # shortcut
400 my $xz = $x->is_zero();
401 my $yz = $y->is_zero();
402 return 0 if $xz && $yz; # 0 <=> 0
403 return -1 if $xz && !$yz; # 0 <=> +y
404 return 1 if $yz && !$xz; # +x <=> 0
405
406 # adjust so that exponents are equal
407 my $lxm = $x->{_m}->length();
408 my $lym = $y->{_m}->length();
409 my $lx = $lxm + $x->{_e};
410 my $ly = $lym + $y->{_e};
394e6ffb 411 my $l = $lx - $ly;
ee15d750 412 return $l <=> 0 if $l != 0;
58cde26e 413
ee15d750 414 # lengths (corrected by exponent) are equal
394e6ffb 415 # so make mantissa equal-length by padding with zero (shift left)
ee15d750 416 my $diff = $lxm - $lym;
417 my $xm = $x->{_m}; # not yet copy it
418 my $ym = $y->{_m};
419 if ($diff > 0)
420 {
421 $ym = $y->{_m}->copy()->blsft($diff,10);
422 }
423 elsif ($diff < 0)
424 {
425 $xm = $x->{_m}->copy()->blsft(-$diff,10);
426 }
b3abae2a 427 $xm->bcmp($ym) <=> 0;
58cde26e 428 }
a0d0e21e 429
58cde26e 430sub badd
431 {
432 # add second arg (BFLOAT or string) to first (BFLOAT) (modifies first)
433 # return result as BFLOAT
58cde26e 434 my ($self,$x,$y,$a,$p,$r) = objectify(2,@_);
435
b3abae2a 436 #print "mbf badd $x $y\n";
574bacfe 437 # inf and NaN handling
438 if (($x->{sign} !~ /^[+-]$/) || ($y->{sign} !~ /^[+-]$/))
439 {
440 # NaN first
441 return $x->bnan() if (($x->{sign} eq $nan) || ($y->{sign} eq $nan));
442 # inf handline
443 if (($x->{sign} =~ /^[+-]inf$/) && ($y->{sign} =~ /^[+-]inf$/))
444 {
445 # + and + => +, - and - => -, + and - => 0, - and + => 0
446 return $x->bzero() if $x->{sign} ne $y->{sign};
447 return $x;
448 }
449 # +-inf + something => +inf
450 # something +-inf => +-inf
451 $x->{sign} = $y->{sign}, return $x if $y->{sign} =~ /^[+-]inf$/;
452 return $x;
453 }
454
58cde26e 455 # speed: no add for 0+y or x+0
456 return $x if $y->is_zero(); # x+0
457 if ($x->is_zero()) # 0+y
458 {
459 # make copy, clobbering up x (modify in place!)
460 $x->{_e} = $y->{_e}->copy();
461 $x->{_m} = $y->{_m}->copy();
462 $x->{sign} = $y->{sign} || $nan;
463 return $x->round($a,$p,$r,$y);
a0d0e21e 464 }
58cde26e 465
466 # take lower of the two e's and adapt m1 to it to match m2
467 my $e = $y->{_e}; $e = Math::BigInt::bzero() if !defined $e; # if no BFLOAT
468 $e = $e - $x->{_e};
469 my $add = $y->{_m}->copy();
470 if ($e < 0)
471 {
58cde26e 472 my $e1 = $e->copy()->babs();
473 $x->{_m} *= (10 ** $e1);
474 $x->{_e} += $e; # need the sign of e
58cde26e 475 }
476 elsif ($e > 0)
477 {
58cde26e 478 $add *= (10 ** $e);
58cde26e 479 }
61f5c3f5 480 # else: both e are the same, so just leave them
481 $x->{_m}->{sign} = $x->{sign}; # fiddle with signs
58cde26e 482 $add->{sign} = $y->{sign};
61f5c3f5 483 $x->{_m} += $add; # finally do add/sub
484 $x->{sign} = $x->{_m}->{sign}; # re-adjust signs
485 $x->{_m}->{sign} = '+'; # mantissa always positiv
486 # delete trailing zeros, then round
487 return $x->bnorm()->round($a,$p,$r,$y);
58cde26e 488 }
489
490sub bsub
491 {
0716bf9b 492 # (BigFloat or num_str, BigFloat or num_str) return BigFloat
58cde26e 493 # subtract second arg from first, modify first
e745a66c 494 my ($self,$x,$y,$a,$p,$r) = objectify(2,@_);
a0d0e21e 495
e745a66c 496 if (!$y->is_zero()) # don't need to do anything if $y is 0
497 {
498 $y->{sign} =~ tr/+\-/-+/; # does nothing for NaN
499 $x->badd($y,$a,$p,$r); # badd does not leave internal zeros
500 $y->{sign} =~ tr/+\-/-+/; # refix $y (does nothing for NaN)
501 }
502 $x; # already rounded by badd()
58cde26e 503 }
504
505sub binc
506 {
507 # increment arg by one
ee15d750 508 my ($self,$x,$a,$p,$r) = ref($_[0]) ? (ref($_[0]),@_) : objectify(1,@_);
e745a66c 509
510 if ($x->{_e}->sign() eq '-')
511 {
512 return $x->badd($self->bone(),$a,$p,$r); # digits after dot
513 }
514
515 if (!$x->{_e}->is_zero())
516 {
517 $x->{_m}->blsft($x->{_e},10); # 1e2 => 100
518 $x->{_e}->bzero();
519 }
520 # now $x->{_e} == 0
521 if ($x->{sign} eq '+')
522 {
523 $x->{_m}->binc();
524 return $x->bnorm()->bround($a,$p,$r);
525 }
526 elsif ($x->{sign} eq '-')
527 {
528 $x->{_m}->bdec();
529 $x->{sign} = '+' if $x->{_m}->is_zero(); # -1 +1 => -0 => +0
530 return $x->bnorm()->bround($a,$p,$r);
531 }
532 # inf, nan handling etc
533 $x->badd($self->__one(),$a,$p,$r); # does round
58cde26e 534 }
535
536sub bdec
537 {
538 # decrement arg by one
ee15d750 539 my ($self,$x,$a,$p,$r) = ref($_[0]) ? (ref($_[0]),@_) : objectify(1,@_);
e745a66c 540
541 if ($x->{_e}->sign() eq '-')
542 {
543 return $x->badd($self->bone('-'),$a,$p,$r); # digits after dot
544 }
545
546 if (!$x->{_e}->is_zero())
547 {
548 $x->{_m}->blsft($x->{_e},10); # 1e2 => 100
549 $x->{_e}->bzero();
550 }
551 # now $x->{_e} == 0
552 my $zero = $x->is_zero();
553 # <= 0
554 if (($x->{sign} eq '-') || $zero)
555 {
556 $x->{_m}->binc();
557 $x->{sign} = '-' if $zero; # 0 => 1 => -1
558 $x->{sign} = '+' if $x->{_m}->is_zero(); # -1 +1 => -0 => +0
559 return $x->bnorm()->round($a,$p,$r);
560 }
561 # > 0
562 elsif ($x->{sign} eq '+')
563 {
564 $x->{_m}->bdec();
565 return $x->bnorm()->round($a,$p,$r);
566 }
567 # inf, nan handling etc
568 $x->badd($self->bone('-'),$a,$p,$r); # does round
58cde26e 569 }
570
61f5c3f5 571sub blog
572 {
b3abae2a 573 my ($self,$x,$base,$a,$p,$r) = ref($_[0]) ? (ref($_[0]),@_) : objectify(2,@_);
61f5c3f5 574
575 # http://www.efunda.com/math/taylor_series/logarithmic.cfm?search_string=log
576
577 # u = x-1, v = x +1
578 # _ _
579 # taylor: | u 1 u^3 1 u^5 |
580 # ln (x) = 2 | --- + - * --- + - * --- + ... | x > 0
b3abae2a 581 # |_ v 3 v^3 5 v^5 _|
61f5c3f5 582
b3abae2a 583 # we need to limit the accuracy to protect against overflow
584 my $fallback = 0;
585 my $scale = 0;
586 my @params = $x->_find_round_parameters($a,$p,$r);
61f5c3f5 587
b3abae2a 588 # no rounding at all, so must use fallback
589 if (scalar @params == 1)
590 {
591 # simulate old behaviour
592 $params[1] = $self->div_scale(); # and round to it as accuracy
593 $scale = $params[1]+4; # at least four more for proper round
594 $params[3] = $r; # round mode by caller or undef
595 $fallback = 1; # to clear a/p afterwards
596 }
597 else
598 {
599 # the 4 below is empirical, and there might be cases where it is not
600 # enough...
601 $scale = abs($params[1] || $params[2]) + 4; # take whatever is defined
602 }
61f5c3f5 603
b3abae2a 604 return $x->bzero(@params) if $x->is_one();
605 return $x->bnan() if $x->{sign} ne '+' || $x->is_zero();
606 #return $x->bone('+',@params) if $x->bcmp($base) == 0;
61f5c3f5 607
b3abae2a 608 # when user set globals, they would interfere with our calculation, so
609 # disable then and later re-enable them
610 no strict 'refs';
611 my $abr = "$self\::accuracy"; my $ab = $$abr; $$abr = undef;
612 my $pbr = "$self\::precision"; my $pb = $$pbr; $$pbr = undef;
613 # we also need to disable any set A or P on $x (_find_round_parameters took
614 # them already into account), since these would interfere, too
615 delete $x->{_a}; delete $x->{_p};
616 # need to disable $upgrade in BigInt, to aoid deep recursion
617 local $Math::BigInt::upgrade = undef;
618
619 my $v = $x->copy(); $v->binc(); # v = x+1
620 $x->bdec(); my $u = $x->copy(); # u = x-1; x = x-1
621
622 $x->bdiv($v,$scale); # first term: u/v
623
624 my $below = $v->copy();
625 my $over = $u->copy();
626 $u *= $u; $v *= $v; # u^2, v^2
627 $below->bmul($v); # u^3, v^3
628 $over->bmul($u);
61f5c3f5 629 my $factor = $self->new(3); my $two = $self->new(2);
630
631 my $diff = $self->bone();
b3abae2a 632 my $limit = $self->new("1E-". ($scale-1)); my $last;
61f5c3f5 633 # print "diff $diff limit $limit\n";
b3abae2a 634 while ($diff->bcmp($limit) > 0)
61f5c3f5 635 {
b3abae2a 636 #print "$x $over $below $factor\n";
61f5c3f5 637 $diff = $x->copy()->bsub($last)->babs();
b3abae2a 638 #print "diff $diff $limit\n";
61f5c3f5 639 $last = $x->copy();
b3abae2a 640 $x += $over->copy()->bdiv($below->copy()->bmul($factor),$scale);
61f5c3f5 641 $over *= $u; $below *= $v; $factor->badd($two);
642 }
643 $x->bmul($two);
b3abae2a 644
645 # shortcut to not run trough _find_round_parameters again
646 if (defined $params[1])
647 {
648 $x->bround($params[1],$params[3]); # then round accordingly
649 }
650 else
651 {
652 $x->bfround($params[2],$params[3]); # then round accordingly
653 }
654 if ($fallback)
655 {
656 # clear a/p after round, since user did not request it
657 $x->{_a} = undef; $x->{_p} = undef;
658 }
659 # restore globals
660 $$abr = $ab; $$pbr = $pb;
661
662 $x;
61f5c3f5 663 }
664
58cde26e 665sub blcm
666 {
ee15d750 667 # (BFLOAT or num_str, BFLOAT or num_str) return BFLOAT
58cde26e 668 # does not modify arguments, but returns new object
669 # Lowest Common Multiplicator
58cde26e 670
671 my ($self,@arg) = objectify(0,@_);
672 my $x = $self->new(shift @arg);
673 while (@arg) { $x = _lcm($x,shift @arg); }
674 $x;
675 }
676
677sub bgcd
678 {
ee15d750 679 # (BFLOAT or num_str, BFLOAT or num_str) return BINT
58cde26e 680 # does not modify arguments, but returns new object
681 # GCD -- Euclids algorithm Knuth Vol 2 pg 296
58cde26e 682
683 my ($self,@arg) = objectify(0,@_);
684 my $x = $self->new(shift @arg);
685 while (@arg) { $x = _gcd($x,shift @arg); }
686 $x;
687 }
688
b3abae2a 689###############################################################################
690# is_foo methods (is_negative, is_positive are inherited from BigInt)
691
692sub is_int
693 {
694 # return true if arg (BFLOAT or num_str) is an integer
695 my ($self,$x) = ref($_[0]) ? (ref($_[0]),$_[0]) : objectify(1,@_);
696
697 return 1 if ($x->{sign} =~ /^[+-]$/) && # NaN and +-inf aren't
698 $x->{_e}->{sign} eq '+'; # 1e-1 => no integer
699 0;
700 }
701
58cde26e 702sub is_zero
703 {
b3abae2a 704 # return true if arg (BFLOAT or num_str) is zero
ee15d750 705 my ($self,$x) = ref($_[0]) ? (ref($_[0]),$_[0]) : objectify(1,@_);
574bacfe 706
707 return 1 if $x->{sign} eq '+' && $x->{_m}->is_zero();
b3abae2a 708 0;
58cde26e 709 }
710
711sub is_one
712 {
b3abae2a 713 # return true if arg (BFLOAT or num_str) is +1 or -1 if signis given
ee15d750 714 my ($self,$x) = ref($_[0]) ? (ref($_[0]),$_[0]) : objectify(1,@_);
715
716 my $sign = shift || ''; $sign = '+' if $sign ne '-';
717 return 1
718 if ($x->{sign} eq $sign && $x->{_e}->is_zero() && $x->{_m}->is_one());
b3abae2a 719 0;
58cde26e 720 }
721
722sub is_odd
723 {
ee15d750 724 # return true if arg (BFLOAT or num_str) is odd or false if even
725 my ($self,$x) = ref($_[0]) ? (ref($_[0]),$_[0]) : objectify(1,@_);
0716bf9b 726
b3abae2a 727 return 1 if ($x->{sign} =~ /^[+-]$/) && # NaN & +-inf aren't
728 ($x->{_e}->is_zero() && $x->{_m}->is_odd());
729 0;
58cde26e 730 }
731
732sub is_even
733 {
b22b3e31 734 # return true if arg (BINT or num_str) is even or false if odd
ee15d750 735 my ($self,$x) = ref($_[0]) ? (ref($_[0]),$_[0]) : objectify(1,@_);
0716bf9b 736
737 return 0 if $x->{sign} !~ /^[+-]$/; # NaN & +-inf aren't
b3abae2a 738# return 1 if $x->{_m}->is_zero(); # 0e1 is even
739 return 1 if ($x->{_e}->{sign} eq '+' # 123.45 is never
740 && $x->{_m}->is_even()); # but 1200 is
741 0;
58cde26e 742 }
743
744sub bmul
745 {
746 # multiply two numbers -- stolen from Knuth Vol 2 pg 233
747 # (BINT or num_str, BINT or num_str) return BINT
748 my ($self,$x,$y,$a,$p,$r) = objectify(2,@_);
58cde26e 749
0716bf9b 750 # print "mbf bmul $x->{_m}e$x->{_e} $y->{_m}e$y->{_e}\n";
58cde26e 751 return $x->bnan() if (($x->{sign} eq $nan) || ($y->{sign} eq $nan));
752
574bacfe 753 # handle result = 0
754 return $x->bzero() if $x->is_zero() || $y->is_zero();
755 # inf handling
756 if (($x->{sign} =~ /^[+-]inf$/) || ($y->{sign} =~ /^[+-]inf$/))
757 {
758 # result will always be +-inf:
759 # +inf * +/+inf => +inf, -inf * -/-inf => +inf
760 # +inf * -/-inf => -inf, -inf * +/+inf => -inf
761 return $x->binf() if ($x->{sign} =~ /^\+/ && $y->{sign} =~ /^\+/);
762 return $x->binf() if ($x->{sign} =~ /^-/ && $y->{sign} =~ /^-/);
763 return $x->binf('-');
764 }
765
58cde26e 766 # aEb * cEd = (a*c)E(b+d)
394e6ffb 767 $x->{_m}->bmul($y->{_m});
768 $x->{_e}->badd($y->{_e});
58cde26e 769 # adjust sign:
770 $x->{sign} = $x->{sign} ne $y->{sign} ? '-' : '+';
394e6ffb 771 return $x->bnorm()->round($a,$p,$r,$y);
58cde26e 772 }
773
774sub bdiv
775 {
776 # (dividend: BFLOAT or num_str, divisor: BFLOAT or num_str) return
777 # (BFLOAT,BFLOAT) (quo,rem) or BINT (only rem)
778 my ($self,$x,$y,$a,$p,$r) = objectify(2,@_);
779
574bacfe 780 # x / +-inf => 0, reminder x
781 return wantarray ? ($x->bzero(),$x->copy()) : $x->bzero()
782 if $y->{sign} =~ /^[+-]inf$/;
783
784 # NaN if x == NaN or y == NaN or x==y==0
58cde26e 785 return wantarray ? ($x->bnan(),bnan()) : $x->bnan()
574bacfe 786 if (($x->is_nan() || $y->is_nan()) ||
787 ($x->is_zero() && $y->is_zero()));
788
789 # 5 / 0 => +inf, -6 / 0 => -inf
790 return wantarray
791 ? ($x->binf($x->{sign}),$self->bnan()) : $x->binf($x->{sign})
792 if ($x->{sign} =~ /^[+-]$/ && $y->is_zero());
0716bf9b 793
394e6ffb 794 # x== 0 or y == 1 or y == -1
795 return wantarray ? ($x,$self->bzero()) : $x if $x->is_zero();
0716bf9b 796
58cde26e 797 # we need to limit the accuracy to protect against overflow
574bacfe 798 my $fallback = 0;
ee15d750 799 my $scale = 0;
ee15d750 800 my @params = $x->_find_round_parameters($a,$p,$r,$y);
801
802 # no rounding at all, so must use fallback
803 if (scalar @params == 1)
58cde26e 804 {
0716bf9b 805 # simulate old behaviour
ee15d750 806 $params[1] = $self->div_scale(); # and round to it as accuracy
61f5c3f5 807 $scale = $params[1]+4; # at least four more for proper round
ee15d750 808 $params[3] = $r; # round mode by caller or undef
809 $fallback = 1; # to clear a/p afterwards
810 }
811 else
812 {
813 # the 4 below is empirical, and there might be cases where it is not
814 # enough...
815 $scale = abs($params[1] || $params[2]) + 4; # take whatever is defined
a0d0e21e 816 }
0716bf9b 817 my $lx = $x->{_m}->length(); my $ly = $y->{_m}->length();
58cde26e 818 $scale = $lx if $lx > $scale;
58cde26e 819 $scale = $ly if $ly > $scale;
0716bf9b 820 my $diff = $ly - $lx;
821 $scale += $diff if $diff > 0; # if lx << ly, but not if ly << lx!
b3abae2a 822
823 # make copy of $x in case of list context for later reminder calculation
824 my $rem;
825 if (wantarray && !$y->is_one())
826 {
827 $rem = $x->copy();
828 }
a0d0e21e 829
58cde26e 830 $x->{sign} = $x->{sign} ne $y->sign() ? '-' : '+';
a0d0e21e 831
58cde26e 832 # check for / +-1 ( +/- 1E0)
394e6ffb 833 if (!$y->is_one())
58cde26e 834 {
394e6ffb 835 # promote BigInts and it's subclasses (except when already a BigFloat)
836 $y = $self->new($y) unless $y->isa('Math::BigFloat');
837
838 # calculate the result to $scale digits and then round it
839 # a * 10 ** b / c * 10 ** d => a/c * 10 ** (b-d)
840 $x->{_m}->blsft($scale,10);
841 $x->{_m}->bdiv( $y->{_m} ); # a/c
842 $x->{_e}->bsub( $y->{_e} ); # b-d
843 $x->{_e}->bsub($scale); # correct for 10**scale
844 $x->bnorm(); # remove trailing 0's
a0d0e21e 845 }
a5f75d66 846
ee15d750 847 # shortcut to not run trough _find_round_parameters again
848 if (defined $params[1])
849 {
61f5c3f5 850 $x->bround($params[1],$params[3]); # then round accordingly
ee15d750 851 }
852 else
853 {
854 $x->bfround($params[2],$params[3]); # then round accordingly
855 }
574bacfe 856 if ($fallback)
857 {
858 # clear a/p after round, since user did not request it
ee15d750 859 $x->{_a} = undef; $x->{_p} = undef;
574bacfe 860 }
0716bf9b 861
58cde26e 862 if (wantarray)
863 {
394e6ffb 864 if (!$y->is_one())
865 {
b3abae2a 866 $rem->bmod($y,$params[1],$params[2],$params[3]); # copy already done
394e6ffb 867 }
868 else
869 {
870 $rem = $self->bzero();
871 }
574bacfe 872 if ($fallback)
873 {
874 # clear a/p after round, since user did not request it
ee15d750 875 $rem->{_a} = undef; $rem->{_p} = undef;
574bacfe 876 }
0716bf9b 877 return ($x,$rem);
58cde26e 878 }
879 return $x;
880 }
a0d0e21e 881
58cde26e 882sub bmod
883 {
884 # (dividend: BFLOAT or num_str, divisor: BFLOAT or num_str) return reminder
885 my ($self,$x,$y,$a,$p,$r) = objectify(2,@_);
a0d0e21e 886
61f5c3f5 887 if (($x->{sign} !~ /^[+-]$/) || ($y->{sign} !~ /^[+-]$/))
888 {
889 my ($d,$re) = $self->SUPER::_div_inf($x,$y);
890 return $re->round($a,$p,$r,$y);
891 }
892 return $x->bnan() if $x->is_zero() && $y->is_zero();
893 return $x if $y->is_zero();
894 return $x->bnan() if $x->is_nan() || $y->is_nan();
895 return $x->bzero() if $y->is_one() || $x->is_zero();
58cde26e 896
61f5c3f5 897 # inf handling is missing here
898
899 my $cmp = $x->bacmp($y); # equal or $x < $y?
900 return $x->bzero($a,$p) if $cmp == 0; # $x == $y => result 0
901
902 # only $y of the operands negative?
903 my $neg = 0; $neg = 1 if $x->{sign} ne $y->{sign};
904
905 $x->{sign} = $y->{sign}; # calc sign first
906 return $x->round($a,$p,$r) if $cmp < 0 && $neg == 0; # $x < $y => result $x
907
908 my $ym = $y->{_m}->copy();
909
910 # 2e1 => 20
911 $ym->blsft($y->{_e},10) if $y->{_e}->{sign} eq '+' && !$y->{_e}->is_zero();
912
913 # if $y has digits after dot
914 my $shifty = 0; # correct _e of $x by this
915 if ($y->{_e}->{sign} eq '-') # has digits after dot
916 {
917 # 123 % 2.5 => 1230 % 25 => 5 => 0.5
918 $shifty = $y->{_e}->copy()->babs(); # no more digits after dot
919 $x->blsft($shifty,10); # 123 => 1230, $y->{_m} is already 25
920 }
921 # $ym is now mantissa of $y based on exponent 0
b3abae2a 922
61f5c3f5 923 my $shiftx = 0; # correct _e of $x by this
924 if ($x->{_e}->{sign} eq '-') # has digits after dot
925 {
926 # 123.4 % 20 => 1234 % 200
927 $shiftx = $x->{_e}->copy()->babs(); # no more digits after dot
928 $ym->blsft($shiftx,10);
929 }
930 # 123e1 % 20 => 1230 % 20
931 if ($x->{_e}->{sign} eq '+' && !$x->{_e}->is_zero())
932 {
933 $x->{_m}->blsft($x->{_e},10);
934 }
935 $x->{_e} = Math::BigInt->bzero() unless $x->{_e}->is_zero();
936
937 $x->{_e}->bsub($shiftx) if $shiftx != 0;
938 $x->{_e}->bsub($shifty) if $shifty != 0;
939
940 # now mantissas are equalized, exponent of $x is adjusted, so calc result
b3abae2a 941# $ym->{sign} = '-' if $neg; # bmod() will make the correction for us
942
61f5c3f5 943 $x->{_m}->bmod($ym);
944
945 $x->{sign} = '+' if $x->{_m}->is_zero(); # fix sign for -0
946 $x->bnorm();
947
948 if ($neg != 0) # one of them negative => correct in place
949 {
950 my $r = $y - $x;
951 $x->{_m} = $r->{_m};
952 $x->{_e} = $r->{_e};
953 $x->{sign} = '+' if $x->{_m}->is_zero(); # fix sign for -0
954 $x->bnorm();
955 }
956
957 $x->round($a,$p,$r,$y); # round and return
58cde26e 958 }
959
960sub bsqrt
961 {
0716bf9b 962 # calculate square root; this should probably
963 # use a different test to see whether the accuracy we want is...
ee15d750 964 my ($self,$x,$a,$p,$r) = ref($_[0]) ? (ref($_[0]),@_) : objectify(1,@_);
58cde26e 965
0716bf9b 966 return $x->bnan() if $x->{sign} eq 'NaN' || $x->{sign} =~ /^-/; # <0, NaN
967 return $x if $x->{sign} eq '+inf'; # +inf
394e6ffb 968 return $x if $x->is_zero() || $x->is_one();
58cde26e 969
61f5c3f5 970 # we need to limit the accuracy to protect against overflow
574bacfe 971 my $fallback = 0;
61f5c3f5 972 my $scale = 0;
973 my @params = $x->_find_round_parameters($a,$p,$r);
974
975 # no rounding at all, so must use fallback
976 if (scalar @params == 1)
0716bf9b 977 {
978 # simulate old behaviour
61f5c3f5 979 $params[1] = $self->div_scale(); # and round to it as accuracy
980 $scale = $params[1]+4; # at least four more for proper round
981 $params[3] = $r; # round mode by caller or undef
ee15d750 982 $fallback = 1; # to clear a/p afterwards
0716bf9b 983 }
61f5c3f5 984 else
985 {
986 # the 4 below is empirical, and there might be cases where it is not
987 # enough...
988 $scale = abs($params[1] || $params[2]) + 4; # take whatever is defined
989 }
990
991 # when user set globals, they would interfere with our calculation, so
992 # disable then and later re-enable them
993 no strict 'refs';
994 my $abr = "$self\::accuracy"; my $ab = $$abr; $$abr = undef;
b3abae2a 995 my $pbr = "$self\::precision"; my $pb = $$pbr; $$pbr = undef;
61f5c3f5 996 # we also need to disable any set A or P on $x (_find_round_parameters took
997 # them already into account), since these would interfere, too
998 delete $x->{_a}; delete $x->{_p};
b3abae2a 999 # need to disable $upgrade in BigInt, to aoid deep recursion
1000 local $Math::BigInt::upgrade = undef;
61f5c3f5 1001
394e6ffb 1002 my $xas = $x->as_number();
1003 my $gs = $xas->copy()->bsqrt(); # some guess
b3abae2a 1004
394e6ffb 1005 if (($x->{_e}->{sign} ne '-') # guess can't be accurate if there are
1006 # digits after the dot
b3abae2a 1007 && ($xas->bacmp($gs * $gs) == 0)) # guess hit the nail on the head?
394e6ffb 1008 {
1009 # exact result
61f5c3f5 1010 $x->{_m} = $gs; $x->{_e} = Math::BigInt->bzero(); $x->bnorm();
1011 # shortcut to not run trough _find_round_parameters again
1012 if (defined $params[1])
1013 {
1014 $x->bround($params[1],$params[3]); # then round accordingly
1015 }
1016 else
1017 {
1018 $x->bfround($params[2],$params[3]); # then round accordingly
1019 }
1020 if ($fallback)
1021 {
1022 # clear a/p after round, since user did not request it
1023 $x->{_a} = undef; $x->{_p} = undef;
1024 }
b3abae2a 1025 ${"$self\::accuracy"} = $ab; ${"$self\::precision"} = $pb;
61f5c3f5 1026 return $x;
394e6ffb 1027 }
61f5c3f5 1028 $gs = $self->new( $gs ); # BigInt to BigFloat
394e6ffb 1029
0716bf9b 1030 my $lx = $x->{_m}->length();
1031 $scale = $lx if $scale < $lx;
394e6ffb 1032 my $e = $self->new("1E-$scale"); # make test variable
b3abae2a 1033# return $x->bnan() if $e->sign() eq 'NaN';
58cde26e 1034
58cde26e 1035 my $y = $x->copy();
394e6ffb 1036 my $two = $self->new(2);
61f5c3f5 1037 my $diff = $e;
ee15d750 1038 # promote BigInts and it's subclasses (except when already a BigFloat)
1039 $y = $self->new($y) unless $y->isa('Math::BigFloat');
61f5c3f5 1040
ee15d750 1041 my $rem;
58cde26e 1042 while ($diff >= $e)
1043 {
61f5c3f5 1044 $rem = $y->copy()->bdiv($gs,$scale)->badd($gs)->bdiv($two,$scale);
1045 $diff = $rem->copy()->bsub($gs)->babs();
1046 $gs = $rem->copy();
a0d0e21e 1047 }
61f5c3f5 1048 # copy over to modify $x
1049 $x->{_m} = $rem->{_m}; $x->{_e} = $rem->{_e};
1050
1051 # shortcut to not run trough _find_round_parameters again
1052 if (defined $params[1])
1053 {
1054 $x->bround($params[1],$params[3]); # then round accordingly
1055 }
1056 else
1057 {
1058 $x->bfround($params[2],$params[3]); # then round accordingly
1059 }
574bacfe 1060 if ($fallback)
1061 {
1062 # clear a/p after round, since user did not request it
ee15d750 1063 $x->{_a} = undef; $x->{_p} = undef;
574bacfe 1064 }
61f5c3f5 1065 # restore globals
b3abae2a 1066 $$abr = $ab; $$pbr = $pb;
574bacfe 1067 $x;
58cde26e 1068 }
1069
b3abae2a 1070sub bfac
1071 {
1072 # (BINT or num_str, BINT or num_str) return BINT
1073 # compute factorial numbers
1074 # modifies first argument
1075 my ($self,$x,@r) = objectify(1,@_);
1076
1077 return $x->bnan() if $x->{sign} ne '+'; # inf, NnN, <0 etc => NaN
1078 return $x->bone(@r) if $x->is_zero() || $x->is_one(); # 0 or 1 => 1
1079
1080 return $x->bnan() if $x->{_e}->{sign} ne '+'; # digits after dot?
1081
1082 # use BigInt's bfac() for faster calc
1083 $x->{_m}->blsft($x->{_e},10); # un-norm m
1084 $x->{_e}->bzero(); # norm $x again
1085 $x->{_m}->bfac(); # factorial
1086 $x->bnorm();
1087
1088 #my $n = $x->copy();
1089 #$x->bone();
1090 #my $f = $self->new(2);
1091 #while ($f->bacmp($n) < 0)
1092 # {
1093 # $x->bmul($f); $f->binc();
1094 # }
1095 #$x->bmul($f); # last step
1096 $x->round(@r); # round
1097 }
1098
58cde26e 1099sub bpow
1100 {
1101 # (BFLOAT or num_str, BFLOAT or num_str) return BFLOAT
1102 # compute power of two numbers, second arg is used as integer
1103 # modifies first argument
1104
1105 my ($self,$x,$y,$a,$p,$r) = objectify(2,@_);
1106
0716bf9b 1107 return $x if $x->{sign} =~ /^[+-]inf$/;
58cde26e 1108 return $x->bnan() if $x->{sign} eq $nan || $y->{sign} eq $nan;
574bacfe 1109 return $x->bone() if $y->is_zero();
58cde26e 1110 return $x if $x->is_one() || $y->is_one();
ee15d750 1111 my $y1 = $y->as_number(); # make bigint (trunc)
394e6ffb 1112 # if ($x == -1)
1113 if ($x->{sign} eq '-' && $x->{_m}->is_one() && $x->{_e}->is_zero())
58cde26e 1114 {
1115 # if $x == -1 and odd/even y => +1/-1 because +-1 ^ (+-1) => +-1
0716bf9b 1116 return $y1->is_odd() ? $x : $x->babs(1);
288d023a 1117 }
58cde26e 1118 return $x if $x->is_zero() && $y->{sign} eq '+'; # 0**y => 0 (if not y <= 0)
574bacfe 1119 # 0 ** -y => 1 / (0 ** y) => / 0! (1 / 0 => +inf)
1120 return $x->binf() if $x->is_zero() && $y->{sign} eq '-';
58cde26e 1121
1122 # calculate $x->{_m} ** $y and $x->{_e} * $y separately (faster)
1123 $y1->babs();
1124 $x->{_m}->bpow($y1);
1125 $x->{_e}->bmul($y1);
1126 $x->{sign} = $nan if $x->{_m}->{sign} eq $nan || $x->{_e}->{sign} eq $nan;
1127 $x->bnorm();
1128 if ($y->{sign} eq '-')
1129 {
1130 # modify $x in place!
0716bf9b 1131 my $z = $x->copy(); $x->bzero()->binc();
58cde26e 1132 return $x->bdiv($z,$a,$p,$r); # round in one go (might ignore y's A!)
a0d0e21e 1133 }
58cde26e 1134 return $x->round($a,$p,$r,$y);
1135 }
1136
1137###############################################################################
1138# rounding functions
1139
1140sub bfround
1141 {
1142 # precision: round to the $Nth digit left (+$n) or right (-$n) from the '.'
1143 # $n == 0 means round to integer
1144 # expects and returns normalized numbers!
ee15d750 1145 my $x = shift; my $self = ref($x) || $x; $x = $self->new(shift) if !ref($x);
a0d0e21e 1146
58cde26e 1147 return $x if $x->modify('bfround');
1148
ee15d750 1149 my ($scale,$mode) = $x->_scale_p($self->precision(),$self->round_mode(),@_);
58cde26e 1150 return $x if !defined $scale; # no-op
1151
574bacfe 1152 # never round a 0, +-inf, NaN
61f5c3f5 1153 if ($x->is_zero())
1154 {
1155 $x->{_p} = $scale if !defined $x->{_p} || $x->{_p} < $scale; # -3 < -2
1156 return $x;
1157 }
1158 return $x if $x->{sign} !~ /^[+-]$/;
58cde26e 1159 # print "MBF bfround $x to scale $scale mode $mode\n";
58cde26e 1160
ee15d750 1161 # don't round if x already has lower precision
1162 return $x if (defined $x->{_p} && $x->{_p} < 0 && $scale < $x->{_p});
1163
1164 $x->{_p} = $scale; # remember round in any case
1165 $x->{_a} = undef; # and clear A
58cde26e 1166 if ($scale < 0)
1167 {
1168 # print "bfround scale $scale e $x->{_e}\n";
1169 # round right from the '.'
1170 return $x if $x->{_e} >= 0; # nothing to round
1171 $scale = -$scale; # positive for simplicity
1172 my $len = $x->{_m}->length(); # length of mantissa
1173 my $dad = -$x->{_e}; # digits after dot
1174 my $zad = 0; # zeros after dot
1175 $zad = -$len-$x->{_e} if ($x->{_e} < -$len);# for 0.00..00xxx style
ee15d750 1176 #print "scale $scale dad $dad zad $zad len $len\n";
58cde26e 1177
1178 # number bsstr len zad dad
1179 # 0.123 123e-3 3 0 3
1180 # 0.0123 123e-4 3 1 4
1181 # 0.001 1e-3 1 2 3
1182 # 1.23 123e-2 3 0 2
1183 # 1.2345 12345e-4 5 0 4
1184
1185 # do not round after/right of the $dad
1186 return $x if $scale > $dad; # 0.123, scale >= 3 => exit
1187
ee15d750 1188 # round to zero if rounding inside the $zad, but not for last zero like:
1189 # 0.0065, scale -2, round last '0' with following '65' (scale == zad case)
1190 return $x->bzero() if $scale < $zad;
1191 if ($scale == $zad) # for 0.006, scale -3 and trunc
58cde26e 1192 {
b3abae2a 1193 $scale = -$len;
58cde26e 1194 }
1195 else
1196 {
1197 # adjust round-point to be inside mantissa
1198 if ($zad != 0)
1199 {
1200 $scale = $scale-$zad;
1201 }
1202 else
1203 {
1204 my $dbd = $len - $dad; $dbd = 0 if $dbd < 0; # digits before dot
1205 $scale = $dbd+$scale;
1206 }
1207 }
1208 # print "round to $x->{_m} to $scale\n";
a0d0e21e 1209 }
58cde26e 1210 else
1211 {
1212 # 123 => 100 means length(123) = 3 - $scale (2) => 1
a5f75d66 1213
b3abae2a 1214 my $dbt = $x->{_m}->length();
1215 # digits before dot
1216 my $dbd = $dbt + $x->{_e};
1217 # should be the same, so treat it as this
1218 $scale = 1 if $scale == 0;
1219 # shortcut if already integer
1220 return $x if $scale == 1 && $dbt <= $dbd;
1221 # maximum digits before dot
1222 ++$dbd;
1223
1224 if ($scale > $dbd)
1225 {
1226 # not enough digits before dot, so round to zero
1227 return $x->bzero;
1228 }
1229 elsif ( $scale == $dbd )
1230 {
1231 # maximum
1232 $scale = -$dbt;
1233 }
58cde26e 1234 else
b3abae2a 1235 {
1236 $scale = $dbd - $scale;
1237 }
1238
a0d0e21e 1239 }
574bacfe 1240 # print "using $scale for $x->{_m} with '$mode'\n";
1241 # pass sign to bround for rounding modes '+inf' and '-inf'
58cde26e 1242 $x->{_m}->{sign} = $x->{sign};
1243 $x->{_m}->bround($scale,$mode);
1244 $x->{_m}->{sign} = '+'; # fix sign back
1245 $x->bnorm();
1246 }
1247
1248sub bround
1249 {
1250 # accuracy: preserve $N digits, and overwrite the rest with 0's
ee15d750 1251 my $x = shift; my $self = ref($x) || $x; $x = $self->new(shift) if !ref($x);
1252
1253 die ('bround() needs positive accuracy') if ($_[0] || 0) < 0;
58cde26e 1254
ee15d750 1255 my ($scale,$mode) = $x->_scale_a($self->accuracy(),$self->round_mode(),@_);
1256 return $x if !defined $scale; # no-op
61f5c3f5 1257
58cde26e 1258 return $x if $x->modify('bround');
61f5c3f5 1259
ee15d750 1260 # scale is now either $x->{_a}, $accuracy, or the user parameter
1261 # test whether $x already has lower accuracy, do nothing in this case
1262 # but do round if the accuracy is the same, since a math operation might
1263 # want to round a number with A=5 to 5 digits afterwards again
1264 return $x if defined $_[0] && defined $x->{_a} && $x->{_a} < $_[0];
58cde26e 1265
61f5c3f5 1266 # scale < 0 makes no sense
1267 # never round a +-inf, NaN
1268 return $x if ($scale < 0) || $x->{sign} !~ /^[+-]$/;
58cde26e 1269
61f5c3f5 1270 # 1: $scale == 0 => keep all digits
1271 # 2: never round a 0
1272 # 3: if we should keep more digits than the mantissa has, do nothing
1273 if ($scale == 0 || $x->is_zero() || $x->{_m}->length() <= $scale)
1274 {
1275 $x->{_a} = $scale if !defined $x->{_a} || $x->{_a} > $scale;
1276 return $x;
1277 }
f216259d 1278
58cde26e 1279 # pass sign to bround for '+inf' and '-inf' rounding modes
1280 $x->{_m}->{sign} = $x->{sign};
1281 $x->{_m}->bround($scale,$mode); # round mantissa
1282 $x->{_m}->{sign} = '+'; # fix sign back
61f5c3f5 1283 # $x->{_m}->{_a} = undef; $x->{_m}->{_p} = undef;
ee15d750 1284 $x->{_a} = $scale; # remember rounding
1285 $x->{_p} = undef; # and clear P
574bacfe 1286 $x->bnorm(); # del trailing zeros gen. by bround()
58cde26e 1287 }
1288
1289sub bfloor
1290 {
1291 # return integer less or equal then $x
ee15d750 1292 my ($self,$x,$a,$p,$r) = ref($_[0]) ? (ref($_[0]),@_) : objectify(1,@_);
58cde26e 1293
1294 return $x if $x->modify('bfloor');
1295
1296 return $x if $x->{sign} !~ /^[+-]$/; # nan, +inf, -inf
1297
1298 # if $x has digits after dot
1299 if ($x->{_e}->{sign} eq '-')
1300 {
1301 $x->{_m}->brsft(-$x->{_e},10);
1302 $x->{_e}->bzero();
1303 $x-- if $x->{sign} eq '-';
f216259d 1304 }
61f5c3f5 1305 $x->round($a,$p,$r);
58cde26e 1306 }
288d023a 1307
58cde26e 1308sub bceil
1309 {
1310 # return integer greater or equal then $x
ee15d750 1311 my ($self,$x,$a,$p,$r) = ref($_[0]) ? (ref($_[0]),@_) : objectify(1,@_);
58cde26e 1312
1313 return $x if $x->modify('bceil');
1314 return $x if $x->{sign} !~ /^[+-]$/; # nan, +inf, -inf
1315
1316 # if $x has digits after dot
1317 if ($x->{_e}->{sign} eq '-')
1318 {
1319 $x->{_m}->brsft(-$x->{_e},10);
1320 $x->{_e}->bzero();
1321 $x++ if $x->{sign} eq '+';
a0d0e21e 1322 }
61f5c3f5 1323 $x->round($a,$p,$r);
58cde26e 1324 }
1325
394e6ffb 1326sub brsft
1327 {
1328 # shift right by $y (divide by power of 2)
1329 my ($self,$x,$y,$n,$a,$p,$r) = objectify(2,@_);
1330
1331 return $x if $x->modify('brsft');
1332 return $x if $x->{sign} !~ /^[+-]$/; # nan, +inf, -inf
1333
1334 $n = 2 if !defined $n; $n = Math::BigFloat->new($n);
1335 $x->bdiv($n ** $y,$a,$p,$r,$y);
1336 }
1337
1338sub blsft
1339 {
1340 # shift right by $y (divide by power of 2)
1341 my ($self,$x,$y,$n,$a,$p,$r) = objectify(2,@_);
1342
1343 return $x if $x->modify('brsft');
1344 return $x if $x->{sign} !~ /^[+-]$/; # nan, +inf, -inf
1345
1346 $n = 2 if !defined $n; $n = Math::BigFloat->new($n);
1347 $x->bmul($n ** $y,$a,$p,$r,$y);
1348 }
1349
58cde26e 1350###############################################################################
a5f75d66 1351
58cde26e 1352sub DESTROY
1353 {
ee15d750 1354 # going through AUTOLOAD for every DESTROY is costly, so avoid it by empty sub
58cde26e 1355 }
1356
1357sub AUTOLOAD
1358 {
b3abae2a 1359 # make fxxx and bxxx both work by selectively mapping fxxx() to MBF::bxxx()
1360 # or falling back to MBI::bxxx()
58cde26e 1361 my $name = $AUTOLOAD;
1362
1363 $name =~ s/.*:://; # split package
ee15d750 1364 no strict 'refs';
1365 if (!method_alias($name))
58cde26e 1366 {
ee15d750 1367 if (!defined $name)
1368 {
1369 # delayed load of Carp and avoid recursion
1370 require Carp;
1371 Carp::croak ("Can't call a method without name");
1372 }
ee15d750 1373 if (!method_hand_up($name))
1374 {
1375 # delayed load of Carp and avoid recursion
1376 require Carp;
1377 Carp::croak ("Can't call $class\-\>$name, not a valid method");
1378 }
1379 # try one level up, but subst. bxxx() for fxxx() since MBI only got bxxx()
1380 $name =~ s/^f/b/;
1381 return &{'Math::BigInt'."::$name"}(@_);
a0d0e21e 1382 }
58cde26e 1383 my $bname = $name; $bname =~ s/^f/b/;
b3abae2a 1384 *{$class."::$name"} = \&$bname;
58cde26e 1385 &$bname; # uses @_
1386 }
1387
1388sub exponent
1389 {
1390 # return a copy of the exponent
ee15d750 1391 my ($self,$x) = ref($_[0]) ? (ref($_[0]),$_[0]) : objectify(1,@_);
58cde26e 1392
ee15d750 1393 if ($x->{sign} !~ /^[+-]$/)
1394 {
1395 my $s = $x->{sign}; $s =~ s/^[+-]//;
1396 return $self->new($s); # -inf, +inf => +inf
1397 }
1398 return $x->{_e}->copy();
58cde26e 1399 }
1400
1401sub mantissa
1402 {
1403 # return a copy of the mantissa
ee15d750 1404 my ($self,$x) = ref($_[0]) ? (ref($_[0]),$_[0]) : objectify(1,@_);
58cde26e 1405
ee15d750 1406 if ($x->{sign} !~ /^[+-]$/)
1407 {
1408 my $s = $x->{sign}; $s =~ s/^[+]//;
1409 return $self->new($s); # -inf, +inf => +inf
1410 }
1411 my $m = $x->{_m}->copy(); # faster than going via bstr()
1412 $m->bneg() if $x->{sign} eq '-';
58cde26e 1413
61f5c3f5 1414 $m;
58cde26e 1415 }
1416
1417sub parts
1418 {
1419 # return a copy of both the exponent and the mantissa
ee15d750 1420 my ($self,$x) = ref($_[0]) ? (ref($_[0]),$_[0]) : objectify(1,@_);
58cde26e 1421
ee15d750 1422 if ($x->{sign} !~ /^[+-]$/)
1423 {
1424 my $s = $x->{sign}; $s =~ s/^[+]//; my $se = $s; $se =~ s/^[-]//;
1425 return ($self->new($s),$self->new($se)); # +inf => inf and -inf,+inf => inf
1426 }
1427 my $m = $x->{_m}->copy(); # faster than going via bstr()
1428 $m->bneg() if $x->{sign} eq '-';
1429 return ($m,$x->{_e}->copy());
58cde26e 1430 }
1431
1432##############################################################################
1433# private stuff (internal use only)
1434
58cde26e 1435sub import
1436 {
1437 my $self = shift;
b3abae2a 1438 my $l = scalar @_; my $j = 0; my @a = @_;
1439 for ( my $i = 0; $i < $l ; $i++, $j++)
58cde26e 1440 {
1441 if ( $_[$i] eq ':constant' )
1442 {
1443 # this rest causes overlord er load to step in
1444 # print "overload @_\n";
1445 overload::constant float => sub { $self->new(shift); };
b3abae2a 1446 splice @a, $j, 1; $j--;
1447 }
1448 elsif ($_[$i] eq 'upgrade')
1449 {
1450 # this causes upgrading
1451 $upgrade = $_[$i+1]; # or undef to disable
1452 my $s = 2; $s = 1 if @a-$j < 2; # avoid "can not modify non-existant..."
1453 splice @a, $j, $s; $j -= $s;
58cde26e 1454 }
1455 }
1456 # any non :constant stuff is handled by our parent, Exporter
1457 # even if @_ is empty, to give it a chance
b3abae2a 1458 $self->SUPER::import(@a); # for subclasses
1459 $self->export_to_level(1,$self,@a); # need this, too
58cde26e 1460 }
1461
1462sub bnorm
1463 {
1464 # adjust m and e so that m is smallest possible
1465 # round number according to accuracy and precision settings
ee15d750 1466 my ($self,$x) = ref($_[0]) ? (ref($_[0]),$_[0]) : objectify(1,@_);
58cde26e 1467
0716bf9b 1468 return $x if $x->{sign} !~ /^[+-]$/; # inf, nan etc
58cde26e 1469
1470 my $zeros = $x->{_m}->_trailing_zeros(); # correct for trailing zeros
1471 if ($zeros != 0)
1472 {
1473 $x->{_m}->brsft($zeros,10); $x->{_e} += $zeros;
1474 }
ee15d750 1475 # for something like 0Ey, set y to 1, and -0 => +0
1476 $x->{sign} = '+', $x->{_e}->bone() if $x->{_m}->is_zero();
1477 # this is to prevent automatically rounding when MBI's globals are set
0716bf9b 1478 $x->{_m}->{_f} = MB_NEVER_ROUND;
1479 $x->{_e}->{_f} = MB_NEVER_ROUND;
ee15d750 1480 # 'forget' that mantissa was rounded via MBI::bround() in MBF's bfround()
1481 $x->{_m}->{_a} = undef; $x->{_e}->{_a} = undef;
1482 $x->{_m}->{_p} = undef; $x->{_e}->{_p} = undef;
61f5c3f5 1483 $x; # MBI bnorm is no-op, so dont call it
1484 }
58cde26e 1485
1486##############################################################################
1487# internal calculation routines
1488
1489sub as_number
1490 {
394e6ffb 1491 # return copy as a bigint representation of this BigFloat number
1492 my ($self,$x) = ref($_[0]) ? (ref($_[0]),$_[0]) : objectify(1,@_);
58cde26e 1493
1494 my $z;
1495 if ($x->{_e}->is_zero())
1496 {
1497 $z = $x->{_m}->copy();
1498 $z->{sign} = $x->{sign};
1499 return $z;
1500 }
0716bf9b 1501 $z = $x->{_m}->copy();
58cde26e 1502 if ($x->{_e} < 0)
1503 {
0716bf9b 1504 $z->brsft(-$x->{_e},10);
1505 }
1506 else
1507 {
1508 $z->blsft($x->{_e},10);
58cde26e 1509 }
58cde26e 1510 $z->{sign} = $x->{sign};
61f5c3f5 1511 $z;
58cde26e 1512 }
1513
1514sub length
1515 {
ee15d750 1516 my $x = shift;
1517 my $class = ref($x) || $x;
1518 $x = $class->new(shift) unless ref($x);
58cde26e 1519
ee15d750 1520 return 1 if $x->{_m}->is_zero();
58cde26e 1521 my $len = $x->{_m}->length();
1522 $len += $x->{_e} if $x->{_e}->sign() eq '+';
1523 if (wantarray())
1524 {
1525 my $t = Math::BigInt::bzero();
1526 $t = $x->{_e}->copy()->babs() if $x->{_e}->sign() eq '-';
1527 return ($len,$t);
1528 }
61f5c3f5 1529 $len;
58cde26e 1530 }
a0d0e21e 1531
15321;
a5f75d66 1533__END__
1534
1535=head1 NAME
1536
58cde26e 1537Math::BigFloat - Arbitrary size floating point math package
a5f75d66 1538
1539=head1 SYNOPSIS
1540
a2008d6d 1541 use Math::BigFloat;
58cde26e 1542
b3abae2a 1543 # Number creation
1544 $x = Math::BigFloat->new($str); # defaults to 0
1545 $nan = Math::BigFloat->bnan(); # create a NotANumber
1546 $zero = Math::BigFloat->bzero(); # create a +0
1547 $inf = Math::BigFloat->binf(); # create a +inf
1548 $inf = Math::BigFloat->binf('-'); # create a -inf
1549 $one = Math::BigFloat->bone(); # create a +1
1550 $one = Math::BigFloat->bone('-'); # create a -1
58cde26e 1551
1552 # Testing
b3abae2a 1553 $x->is_zero(); # true if arg is +0
1554 $x->is_nan(); # true if arg is NaN
0716bf9b 1555 $x->is_one(); # true if arg is +1
1556 $x->is_one('-'); # true if arg is -1
1557 $x->is_odd(); # true if odd, false for even
1558 $x->is_even(); # true if even, false for odd
1559 $x->is_positive(); # true if >= 0
1560 $x->is_negative(); # true if < 0
b3abae2a 1561 $x->is_inf(sign); # true if +inf, or -inf (default is '+')
1562
58cde26e 1563 $x->bcmp($y); # compare numbers (undef,<0,=0,>0)
1564 $x->bacmp($y); # compare absolutely (undef,<0,=0,>0)
1565 $x->sign(); # return the sign, either +,- or NaN
b3abae2a 1566 $x->digit($n); # return the nth digit, counting from right
1567 $x->digit(-$n); # return the nth digit, counting from left
58cde26e 1568
1569 # The following all modify their first argument:
b3abae2a 1570
58cde26e 1571 # set
1572 $x->bzero(); # set $i to 0
1573 $x->bnan(); # set $i to NaN
b3abae2a 1574 $x->bone(); # set $x to +1
1575 $x->bone('-'); # set $x to -1
1576 $x->binf(); # set $x to inf
1577 $x->binf('-'); # set $x to -inf
58cde26e 1578
1579 $x->bneg(); # negation
1580 $x->babs(); # absolute value
1581 $x->bnorm(); # normalize (no-op)
1582 $x->bnot(); # two's complement (bit wise not)
1583 $x->binc(); # increment x by 1
1584 $x->bdec(); # decrement x by 1
1585
1586 $x->badd($y); # addition (add $y to $x)
1587 $x->bsub($y); # subtraction (subtract $y from $x)
1588 $x->bmul($y); # multiplication (multiply $x by $y)
1589 $x->bdiv($y); # divide, set $i to quotient
1590 # return (quo,rem) or quo if scalar
1591
1592 $x->bmod($y); # modulus
1593 $x->bpow($y); # power of arguments (a**b)
1594 $x->blsft($y); # left shift
1595 $x->brsft($y); # right shift
1596 # return (quo,rem) or quo if scalar
1597
61f5c3f5 1598 $x->blog($base); # logarithm of $x, base defaults to e
1599 # (other bases than e not supported yet)
1600
58cde26e 1601 $x->band($y); # bit-wise and
1602 $x->bior($y); # bit-wise inclusive or
1603 $x->bxor($y); # bit-wise exclusive or
1604 $x->bnot(); # bit-wise not (two's complement)
b3abae2a 1605
1606 $x->bsqrt(); # calculate square-root
1607 $x->bfac(); # factorial of $x (1*2*3*4*..$x)
1608
58cde26e 1609 $x->bround($N); # accuracy: preserver $N digits
1610 $x->bfround($N); # precision: round to the $Nth digit
1611
1612 # The following do not modify their arguments:
58cde26e 1613 bgcd(@values); # greatest common divisor
1614 blcm(@values); # lowest common multiplicator
1615
1616 $x->bstr(); # return string
1617 $x->bsstr(); # return string in scientific notation
b3abae2a 1618
1619 $x->bfloor(); # return integer less or equal than $x
1620 $x->bceil(); # return integer greater or equal than $x
1621
58cde26e 1622 $x->exponent(); # return exponent as BigInt
1623 $x->mantissa(); # return mantissa as BigInt
1624 $x->parts(); # return (mantissa,exponent) as BigInt
1625
1626 $x->length(); # number of digits (w/o sign and '.')
1627 ($l,$f) = $x->length(); # number of digits, and length of fraction
a5f75d66 1628
1629=head1 DESCRIPTION
1630
58cde26e 1631All operators (inlcuding basic math operations) are overloaded if you
1632declare your big floating point numbers as
a5f75d66 1633
58cde26e 1634 $i = new Math::BigFloat '12_3.456_789_123_456_789E-2';
1635
1636Operations with overloaded operators preserve the arguments, which is
1637exactly what you expect.
1638
1639=head2 Canonical notation
1640
1641Input to these routines are either BigFloat objects, or strings of the
1642following four forms:
a5f75d66 1643
1644=over 2
1645
58cde26e 1646=item *
1647
1648C</^[+-]\d+$/>
a5f75d66 1649
58cde26e 1650=item *
a5f75d66 1651
58cde26e 1652C</^[+-]\d+\.\d*$/>
a5f75d66 1653
58cde26e 1654=item *
a5f75d66 1655
58cde26e 1656C</^[+-]\d+E[+-]?\d+$/>
a5f75d66 1657
58cde26e 1658=item *
a5f75d66 1659
58cde26e 1660C</^[+-]\d*\.\d+E[+-]?\d+$/>
5d7098d5 1661
58cde26e 1662=back
1663
1664all with optional leading and trailing zeros and/or spaces. Additonally,
1665numbers are allowed to have an underscore between any two digits.
1666
1667Empty strings as well as other illegal numbers results in 'NaN'.
1668
1669bnorm() on a BigFloat object is now effectively a no-op, since the numbers
1670are always stored in normalized form. On a string, it creates a BigFloat
1671object.
1672
1673=head2 Output
1674
1675Output values are BigFloat objects (normalized), except for bstr() and bsstr().
1676
1677The string output will always have leading and trailing zeros stripped and drop
1678a plus sign. C<bstr()> will give you always the form with a decimal point,
1679while C<bsstr()> (for scientific) gives you the scientific notation.
1680
1681 Input bstr() bsstr()
1682 '-0' '0' '0E1'
1683 ' -123 123 123' '-123123123' '-123123123E0'
1684 '00.0123' '0.0123' '123E-4'
1685 '123.45E-2' '1.2345' '12345E-4'
1686 '10E+3' '10000' '1E4'
1687
1688Some routines (C<is_odd()>, C<is_even()>, C<is_zero()>, C<is_one()>,
1689C<is_nan()>) return true or false, while others (C<bcmp()>, C<bacmp()>)
1690return either undef, <0, 0 or >0 and are suited for sort.
1691
1692Actual math is done by using BigInts to represent the mantissa and exponent.
1693The sign C</^[+-]$/> is stored separately. The string 'NaN' is used to
1694represent the result when input arguments are not numbers, as well as
1695the result of dividing by zero.
1696
1697=head2 C<mantissa()>, C<exponent()> and C<parts()>
1698
1699C<mantissa()> and C<exponent()> return the said parts of the BigFloat
1700as BigInts such that:
1701
1702 $m = $x->mantissa();
1703 $e = $x->exponent();
1704 $y = $m * ( 10 ** $e );
1705 print "ok\n" if $x == $y;
1706
1707C<< ($m,$e) = $x->parts(); >> is just a shortcut giving you both of them.
1708
1709A zero is represented and returned as C<0E1>, B<not> C<0E0> (after Knuth).
1710
1711Currently the mantissa is reduced as much as possible, favouring higher
1712exponents over lower ones (e.g. returning 1e7 instead of 10e6 or 10000000e0).
1713This might change in the future, so do not depend on it.
1714
1715=head2 Accuracy vs. Precision
1716
1717See also: L<Rounding|Rounding>.
1718
027dc388 1719Math::BigFloat supports both precision and accuracy. For a full documentation,
1720examples and tips on these topics please see the large section in
1721L<Math::BigInt>.
5d7098d5 1722
58cde26e 1723Since things like sqrt(2) or 1/3 must presented with a limited precision lest
1724a operation consumes all resources, each operation produces no more than
1725C<Math::BigFloat::precision()> digits.
1726
1727In case the result of one operation has more precision than specified,
1728it is rounded. The rounding mode taken is either the default mode, or the one
1729supplied to the operation after the I<scale>:
1730
1731 $x = Math::BigFloat->new(2);
1732 Math::BigFloat::precision(5); # 5 digits max
1733 $y = $x->copy()->bdiv(3); # will give 0.66666
1734 $y = $x->copy()->bdiv(3,6); # will give 0.666666
1735 $y = $x->copy()->bdiv(3,6,'odd'); # will give 0.666667
1736 Math::BigFloat::round_mode('zero');
1737 $y = $x->copy()->bdiv(3,6); # will give 0.666666
1738
1739=head2 Rounding
1740
1741=over 2
1742
5dc6f178 1743=item ffround ( +$scale )
58cde26e 1744
0716bf9b 1745Rounds to the $scale'th place left from the '.', counting from the dot.
1746The first digit is numbered 1.
58cde26e 1747
5dc6f178 1748=item ffround ( -$scale )
58cde26e 1749
0716bf9b 1750Rounds to the $scale'th place right from the '.', counting from the dot.
58cde26e 1751
5dc6f178 1752=item ffround ( 0 )
1753
0716bf9b 1754Rounds to an integer.
5dc6f178 1755
1756=item fround ( +$scale )
1757
0716bf9b 1758Preserves accuracy to $scale digits from the left (aka significant digits)
1759and pads the rest with zeros. If the number is between 1 and -1, the
1760significant digits count from the first non-zero after the '.'
5dc6f178 1761
1762=item fround ( -$scale ) and fround ( 0 )
1763
0716bf9b 1764These are effetively no-ops.
5d7098d5 1765
a5f75d66 1766=back
1767
0716bf9b 1768All rounding functions take as a second parameter a rounding mode from one of
1769the following: 'even', 'odd', '+inf', '-inf', 'zero' or 'trunc'.
58cde26e 1770
1771The default rounding mode is 'even'. By using
ee15d750 1772C<< Math::BigFloat::round_mode($round_mode); >> you can get and set the default
1773mode for subsequent rounding. The usage of C<$Math::BigFloat::$round_mode> is
0716bf9b 1774no longer supported.
b22b3e31 1775The second parameter to the round functions then overrides the default
0716bf9b 1776temporarily.
58cde26e 1777
1778The C<< as_number() >> function returns a BigInt from a Math::BigFloat. It uses
1779'trunc' as rounding mode to make it equivalent to:
1780
1781 $x = 2.5;
1782 $y = int($x) + 2;
1783
1784You can override this by passing the desired rounding mode as parameter to
1785C<as_number()>:
1786
1787 $x = Math::BigFloat->new(2.5);
1788 $y = $x->as_number('odd'); # $y = 3
1789
1790=head1 EXAMPLES
1791
58cde26e 1792 # not ready yet
58cde26e 1793
1794=head1 Autocreating constants
1795
1796After C<use Math::BigFloat ':constant'> all the floating point constants
1797in the given scope are converted to C<Math::BigFloat>. This conversion
1798happens at compile time.
1799
1800In particular
1801
1802 perl -MMath::BigFloat=:constant -e 'print 2E-100,"\n"'
1803
1804prints the value of C<2E-100>. Note that without conversion of
1805constants the expression 2E-100 will be calculated as normal floating point
1806number.
1807
a5f75d66 1808=head1 BUGS
1809
58cde26e 1810=over 2
1811
1812=item *
1813
1814The following does not work yet:
1815
1816 $m = $x->mantissa();
1817 $e = $x->exponent();
1818 $y = $m * ( 10 ** $e );
1819 print "ok\n" if $x == $y;
1820
1821=item *
1822
1823There is no fmod() function yet.
1824
1825=back
1826
1827=head1 CAVEAT
1828
1829=over 1
1830
1831=item stringify, bstr()
1832
1833Both stringify and bstr() now drop the leading '+'. The old code would return
1834'+1.23', the new returns '1.23'. See the documentation in L<Math::BigInt> for
1835reasoning and details.
1836
1837=item bdiv
1838
1839The following will probably not do what you expect:
1840
1841 print $c->bdiv(123.456),"\n";
1842
1843It prints both quotient and reminder since print works in list context. Also,
1844bdiv() will modify $c, so be carefull. You probably want to use
1845
1846 print $c / 123.456,"\n";
1847 print scalar $c->bdiv(123.456),"\n"; # or if you want to modify $c
1848
1849instead.
1850
1851=item Modifying and =
1852
1853Beware of:
1854
1855 $x = Math::BigFloat->new(5);
1856 $y = $x;
1857
1858It will not do what you think, e.g. making a copy of $x. Instead it just makes
1859a second reference to the B<same> object and stores it in $y. Thus anything
1860that modifies $x will modify $y, and vice versa.
1861
1862 $x->bmul(2);
1863 print "$x, $y\n"; # prints '10, 10'
1864
1865If you want a true copy of $x, use:
1866
1867 $y = $x->copy();
1868
1869See also the documentation in L<overload> regarding C<=>.
1870
1871=item bpow
1872
1873C<bpow()> now modifies the first argument, unlike the old code which left
1874it alone and only returned the result. This is to be consistent with
1875C<badd()> etc. The first will modify $x, the second one won't:
1876
1877 print bpow($x,$i),"\n"; # modify $x
1878 print $x->bpow($i),"\n"; # ditto
1879 print $x ** $i,"\n"; # leave $x alone
1880
1881=back
1882
1883=head1 LICENSE
a5f75d66 1884
58cde26e 1885This program is free software; you may redistribute it and/or modify it under
1886the same terms as Perl itself.
5d7098d5 1887
58cde26e 1888=head1 AUTHORS
5d7098d5 1889
58cde26e 1890Mark Biggar, overloaded interface by Ilya Zakharevich.
1891Completely rewritten by Tels http://bloodgate.com in 2001.
a5f75d66 1892
a5f75d66 1893=cut
Unnamed repository; edit this file to name it for gitweb.