This is a live mirror of the Perl 5 development currently hosted at https://github.com/perl/perl5
Make also the bigintpm.t (like op/sprintf.t) be less demanding
[perl5.git] / lib / Math / BigInt.pm
CommitLineData
58cde26e
JH
1#!/usr/bin/perl -w
2
3# mark.biggar@TrustedSysLabs.com
4# eay@mincom.com is dead (math::BigInteger)
5# see: http://www.cypherspace.org/~adam/rsa/pureperl.html (contacted c. adam
6# on 2000/11/13 - but email is dead
7
8# todo:
9# - fully remove funky $# stuff (maybe)
10# - use integer; vs 1e7 as base
11# - speed issues (XS? Bit::Vector?)
12# - split out actual math code to Math::BigNumber
13
14# Qs: what exactly happens on numify of HUGE numbers? overflow?
15# $a = -$a is much slower (making copy of $a) than $a->bneg(), hm!?
16# (copy_on_write will help there, but that is not yet implemented)
17
18# The following hash values are used:
19# value: the internal array, base 100000
20# sign : +,-,NaN,+inf,-inf
21# _a : accuracy
22# _p : precision
23# _cow : copy on write: number of objects that share the data (NRY)
24# Internally the numbers are stored in an array with at least 1 element, no
25# leading zero parts (except the first) and in base 100000
26
27# USE_MUL: due to problems on certain os (os390, posix-bc) "* 1e-5" is used
28# instead of "/ 1e5" at some places, (marked with USE_MUL). But instead of
29# using the reverse only on problematic machines, I used it everytime to avoid
30# the costly comparisations. This _should_ work everywhere. Thanx Peter Prymmer
b4f14daa 31
58cde26e
JH
32package Math::BigInt;
33my $class = "Math::BigInt";
34
35$VERSION = 1.35;
36use Exporter;
37@ISA = qw( Exporter );
38@EXPORT_OK = qw( bneg babs bcmp badd bmul bdiv bmod bnorm bsub
39 bgcd blcm
40 bround
41 blsft brsft band bior bxor bnot bpow bnan bzero
42 bacmp bstr bsstr binc bdec bint binf bfloor bceil
43 is_odd is_even is_zero is_one is_nan is_inf sign
44 length as_number
45 trace objectify _swap
46 );
47
48#@EXPORT = qw( );
49use vars qw/$rnd_mode $accuracy $precision $div_scale/;
50use strict;
51
52# Inside overload, the first arg is always an object. If the original code had
53# it reversed (like $x = 2 * $y), then the third paramater indicates this
54# swapping. To make it work, we use a helper routine which not only reswaps the
55# params, but also makes a new object in this case. See _swap() for details,
56# especially the cases of operators with different classes.
57
58# For overloaded ops with only one argument we simple use $_[0]->copy() to
59# preserve the argument.
60
61# Thus inheritance of overload operators becomes possible and transparent for
62# our subclasses without the need to repeat the entire overload section there.
a0d0e21e 63
a5f75d66 64use overload
58cde26e
JH
65'=' => sub { $_[0]->copy(); },
66
67# '+' and '-' do not use _swap, since it is a triffle slower. If you want to
68# override _swap (if ever), then override overload of '+' and '-', too!
69# for sub it is a bit tricky to keep b: b-a => -a+b
70'-' => sub { my $c = $_[0]->copy; $_[2] ?
71 $c->bneg()->badd($_[1]) :
72 $c->bsub( $_[1]) },
73'+' => sub { $_[0]->copy()->badd($_[1]); },
74
75# some shortcuts for speed (assumes that reversed order of arguments is routed
76# to normal '+' and we thus can always modify first arg. If this is changed,
77# this breaks and must be adjusted.)
78'+=' => sub { $_[0]->badd($_[1]); },
79'-=' => sub { $_[0]->bsub($_[1]); },
80'*=' => sub { $_[0]->bmul($_[1]); },
81'/=' => sub { scalar $_[0]->bdiv($_[1]); },
82'**=' => sub { $_[0]->bpow($_[1]); },
83
84'<=>' => sub { $_[2] ?
85 $class->bcmp($_[1],$_[0]) :
86 $class->bcmp($_[0],$_[1])},
87'cmp' => sub {
88 $_[2] ?
89 $_[1] cmp $_[0]->bstr() :
90 $_[0]->bstr() cmp $_[1] },
91
92'int' => sub { $_[0]->copy(); },
93'neg' => sub { $_[0]->copy()->bneg(); },
94'abs' => sub { $_[0]->copy()->babs(); },
95'~' => sub { $_[0]->copy()->bnot(); },
96
97'*' => sub { my @a = ref($_[0])->_swap(@_); $a[0]->bmul($a[1]); },
98'/' => sub { my @a = ref($_[0])->_swap(@_);scalar $a[0]->bdiv($a[1]);},
99'%' => sub { my @a = ref($_[0])->_swap(@_); $a[0]->bmod($a[1]); },
100'**' => sub { my @a = ref($_[0])->_swap(@_); $a[0]->bpow($a[1]); },
101'<<' => sub { my @a = ref($_[0])->_swap(@_); $a[0]->blsft($a[1]); },
102'>>' => sub { my @a = ref($_[0])->_swap(@_); $a[0]->brsft($a[1]); },
103
104'&' => sub { my @a = ref($_[0])->_swap(@_); $a[0]->band($a[1]); },
105'|' => sub { my @a = ref($_[0])->_swap(@_); $a[0]->bior($a[1]); },
106'^' => sub { my @a = ref($_[0])->_swap(@_); $a[0]->bxor($a[1]); },
107
108# can modify arg of ++ and --, so avoid a new-copy for speed, but don't
109# use $_[0]->_one(), it modifies $_[0] to be 1!
110'++' => sub { $_[0]->binc() },
111'--' => sub { $_[0]->bdec() },
112
113# if overloaded, O(1) instead of O(N) and twice as fast for small numbers
114'bool' => sub {
115 # this kludge is needed for perl prior 5.6.0 since returning 0 here fails :-/
116 # v5.6.1 dumps on that: return !$_[0]->is_zero() || undef; :-(
117 my $t = !$_[0]->is_zero();
118 undef $t if $t == 0;
119 return $t;
120 },
a0d0e21e
LW
121
122qw(
58cde26e
JH
123"" bstr
1240+ numify), # Order of arguments unsignificant
a5f75d66 125;
a0d0e21e 126
58cde26e
JH
127##############################################################################
128# global constants, flags and accessory
129
130# are NaNs ok?
131my $NaNOK=1;
132# set to 1 for tracing
133my $trace = 0;
134# constants for easier life
135my $nan = 'NaN';
136my $BASE_LEN = 5;
137my $BASE = int("1e".$BASE_LEN); # var for trying to change it to 1e7
138my $RBASE = 1e-5; # see USE_MUL
139
140# Rounding modes one of 'even', 'odd', '+inf', '-inf', 'zero' or 'trunc'
141$rnd_mode = 'even';
142$accuracy = undef;
143$precision = undef;
144$div_scale = 40;
145
146sub round_mode
147 {
148 # make Class->round_mode() work
149 my $self = shift || $class;
150 # shift @_ if defined $_[0] && $_[0] eq $class;
151 if (defined $_[0])
152 {
153 my $m = shift;
154 die "Unknown round mode $m"
155 if $m !~ /^(even|odd|\+inf|\-inf|zero|trunc)$/;
156 $rnd_mode = $m; return;
157 }
158 return $rnd_mode;
159 }
160
161sub accuracy
162 {
163 # $x->accuracy($a); ref($x) a
164 # $x->accuracy(); ref($x);
165 # Class::accuracy(); # not supported
166 #print "MBI @_ ($class)\n";
167 my $x = shift;
168
169 die ("accuracy() needs reference to object as first parameter.")
170 if !ref $x;
171
172 if (@_ > 0)
173 {
174 $x->{_a} = shift;
175 $x->round() if defined $x->{_a};
176 }
177 return $x->{_a};
178 }
179
180sub precision
181 {
182 my $x = shift;
183
184 die ("precision() needs reference to object as first parameter.")
185 unless ref $x;
186
187 if (@_ > 0)
188 {
189 $x->{_p} = shift;
190 $x->round() if defined $x->{_p};
191 }
192 return $x->{_p};
193 }
194
195sub _scale_a
196 {
197 # select accuracy parameter based on precedence,
198 # used by bround() and bfround(), may return undef for scale (means no op)
199 my ($x,$s,$m,$scale,$mode) = @_;
200 $scale = $x->{_a} if !defined $scale;
201 $scale = $s if (!defined $scale);
202 $mode = $m if !defined $mode;
203 return ($scale,$mode);
204 }
205
206sub _scale_p
207 {
208 # select precision parameter based on precedence,
209 # used by bround() and bfround(), may return undef for scale (means no op)
210 my ($x,$s,$m,$scale,$mode) = @_;
211 $scale = $x->{_p} if !defined $scale;
212 $scale = $s if (!defined $scale);
213 $mode = $m if !defined $mode;
214 return ($scale,$mode);
215 }
216
217##############################################################################
218# constructors
219
220sub copy
221 {
222 my ($c,$x);
223 if (@_ > 1)
224 {
225 # if two arguments, the first one is the class to "swallow" subclasses
226 ($c,$x) = @_;
227 }
228 else
229 {
230 $x = shift;
231 $c = ref($x);
232 }
233 return unless ref($x); # only for objects
234
235 my $self = {}; bless $self,$c;
236 foreach my $k (keys %$x)
237 {
238 if (ref($x->{$k}) eq 'ARRAY')
239 {
240 $self->{$k} = [ @{$x->{$k}} ];
241 }
242 elsif (ref($x->{$k}) eq 'HASH')
243 {
244 # only one level deep!
245 foreach my $h (keys %{$x->{$k}})
246 {
247 $self->{$k}->{$h} = $x->{$k}->{$h};
248 }
249 }
250 elsif (ref($x->{$k}))
251 {
252 my $c = ref($x->{$k});
253 $self->{$k} = $c->new($x->{$k}); # no copy() due to deep rec
254 }
255 else
256 {
257 $self->{$k} = $x->{$k};
258 }
259 }
260 $self;
261 }
262
263sub new
264 {
265 # create a new BigInts object from a string or another bigint object.
266 # value => internal array representation
267 # sign => sign (+/-), or "NaN"
268
269 # the argument could be an object, so avoid ||, && etc on it, this would
270 # cause costly overloaded code to be called. The only allowed op are ref()
271 # and definend.
272
273 trace (@_);
274 my $class = shift;
275
276 my $wanted = shift; # avoid numify call by not using || here
277 return $class->bzero() if !defined $wanted; # default to 0
278 return $class->copy($wanted) if ref($wanted);
279
280 my $self = {}; bless $self, $class;
281 # handle '+inf', '-inf' first
282 if ($wanted =~ /^[+-]inf$/)
283 {
284 $self->{value} = [ 0 ];
285 $self->{sign} = $wanted;
286 return $self;
287 }
288 # split str in m mantissa, e exponent, i integer, f fraction, v value, s sign
289 my ($mis,$miv,$mfv,$es,$ev) = _split(\$wanted);
290 if (ref $mis && !ref $miv)
291 {
292 # _from_hex
293 $self->{value} = $mis->{value};
294 $self->{sign} = $mis->{sign};
295 return $self;
296 }
297 if (!ref $mis)
298 {
299 die "$wanted is not a number initialized to $class" if !$NaNOK;
300 #print "NaN 1\n";
301 $self->{value} = [ 0 ];
302 $self->{sign} = $nan;
303 return $self;
304 }
305 # make integer from mantissa by adjusting exp, then convert to bigint
306 $self->{sign} = $$mis; # store sign
307 $self->{value} = [ 0 ]; # for all the NaN cases
308 my $e = int("$$es$$ev"); # exponent (avoid recursion)
309 if ($e > 0)
310 {
311 my $diff = $e - CORE::length($$mfv);
312 if ($diff < 0) # Not integer
313 {
314 #print "NOI 1\n";
315 $self->{sign} = $nan;
316 }
317 else # diff >= 0
318 {
319 # adjust fraction and add it to value
320 # print "diff > 0 $$miv\n";
321 $$miv = $$miv . ($$mfv . '0' x $diff);
322 }
323 }
324 else
325 {
326 if ($$mfv ne '') # e <= 0
327 {
328 # fraction and negative/zero E => NOI
329 #print "NOI 2 \$\$mfv '$$mfv'\n";
330 $self->{sign} = $nan;
331 }
332 elsif ($e < 0)
333 {
334 # xE-y, and empty mfv
335 #print "xE-y\n";
336 $e = abs($e);
337 if ($$miv !~ s/0{$e}$//) # can strip so many zero's?
338 {
339 #print "NOI 3\n";
340 $self->{sign} = $nan;
341 }
342 }
343 }
344 $self->{sign} = '+' if $$miv eq '0'; # normalize -0 => +0
345 $self->_internal($miv) if $self->{sign} ne $nan; # as internal array
346 #print "$wanted => $self->{sign} $self->{value}->[0]\n";
347 # if any of the globals is set, round to them and thus store them insid $self
348 $self->round($accuracy,$precision,$rnd_mode)
349 if defined $accuracy || defined $precision;
350 return $self;
351 }
352
353# some shortcuts for easier life
354sub bint
355 {
356 # exportable version of new
357 trace(@_);
358 return $class->new(@_);
359 }
360
361sub bnan
362 {
363 # create a bigint 'NaN', if given a BigInt, set it to 'NaN'
b4f14daa 364 my $self = shift;
58cde26e
JH
365 $self = $class if !defined $self;
366 if (!ref($self))
367 {
368 my $c = $self; $self = {}; bless $self, $c;
369 }
370 return if $self->modify('bnan');
371 $self->{value} = [ 0 ];
372 $self->{sign} = $nan;
373 trace('NaN');
374 return $self;
b4f14daa 375 }
58cde26e
JH
376
377sub binf
378 {
379 # create a bigint '+-inf', if given a BigInt, set it to '+-inf'
380 # the sign is either '+', or if given, used from there
381 my $self = shift;
382 my $sign = shift; $sign = '+' if !defined $sign || $sign ne '-';
383 $self = $class if !defined $self;
384 if (!ref($self))
385 {
386 my $c = $self; $self = {}; bless $self, $c;
387 }
388 return if $self->modify('binf');
389 $self->{value} = [ 0 ];
390 $self->{sign} = $sign.'inf';
391 trace('inf');
392 return $self;
393 }
394
395sub bzero
396 {
397 # create a bigint '+0', if given a BigInt, set it to 0
398 my $self = shift;
399 $self = $class if !defined $self;
400 if (!ref($self))
401 {
402 my $c = $self; $self = {}; bless $self, $c;
403 }
404 return if $self->modify('bzero');
405 $self->{value} = [ 0 ];
406 $self->{sign} = '+';
407 trace('0');
408 return $self;
409 }
410
411##############################################################################
412# string conversation
413
414sub bsstr
415 {
416 # (ref to BFLOAT or num_str ) return num_str
417 # Convert number from internal format to scientific string format.
418 # internal format is always normalized (no leading zeros, "-0E0" => "+0E0")
419 trace(@_);
420 my ($self,$x) = objectify(1,@_);
421
422 return $x->{sign} if $x->{sign} !~ /^[+-]$/;
423 my ($m,$e) = $x->parts();
424 # can be only '+', so
425 my $sign = 'e+';
426 # MBF: my $s = $e->{sign}; $s = '' if $s eq '-'; my $sep = 'e'.$s;
427 return $m->bstr().$sign.$e->bstr();
428 }
429
430sub bstr
431 {
432 # (ref to BINT or num_str ) return num_str
433 # Convert number from internal base 100000 format to string format.
434 # internal format is always normalized (no leading zeros, "-0" => "+0")
435 trace(@_);
436 my $x = shift; $x = $class->new($x) unless ref $x;
437 # my ($self,$x) = objectify(1,@_);
438
439 return $x->{sign} if $x->{sign} !~ /^[+-]$/;
440 my $ar = $x->{value} || return $nan; # should not happen
441 my $es = "";
442 $es = $x->{sign} if $x->{sign} eq '-'; # get sign, but not '+'
443 my $l = scalar @$ar; # number of parts
444 return $nan if $l < 1; # should not happen
445 # handle first one different to strip leading zeros from it (there are no
446 # leading zero parts in internal representation)
447 $l --; $es .= $ar->[$l]; $l--;
448 # Interestingly, the pre-padd method uses more time
449 # the old grep variant takes longer (14 to 10 sec)
450 while ($l >= 0)
451 {
452 $es .= substr('0000'.$ar->[$l],-5); # fastest way I could think of
453 $l--;
a0d0e21e 454 }
58cde26e
JH
455 return $es;
456 }
457
458sub numify
459 {
460 # Make a number from a BigInt object
461 # old: simple return string and let Perl's atoi() handle the rest
462 # new: calc because it is faster than bstr()+atoi()
463 #trace (@_);
464 #my ($self,$x) = objectify(1,@_);
465 #return $x->bstr(); # ref($x);
466 my $x = shift; $x = $class->new($x) unless ref $x;
467
468 return $nan if $x->{sign} eq $nan;
469 my $fac = 1; $fac = -1 if $x->{sign} eq '-';
470 return $fac*$x->{value}->[0] if @{$x->{value}} == 1; # below $BASE
471 my $num = 0;
472 foreach (@{$x->{value}})
473 {
474 $num += $fac*$_; $fac *= $BASE;
475 }
476 return $num;
477 }
478
479##############################################################################
480# public stuff (usually prefixed with "b")
481
482sub sign
483 {
484 # return the sign of the number: +/-/NaN
485 my ($self,$x) = objectify(1,@_);
486 return $x->{sign};
487 }
488
489sub round
490 {
491 # After any operation or when calling round(), the result is rounded by
492 # regarding the A & P from arguments, local parameters, or globals.
493 # The result's A or P are set by the rounding, but not inspected beforehand
494 # (aka only the arguments enter into it). This works because the given
495 # 'first' argument is both the result and true first argument with unchanged
496 # A and P settings.
497 # This does not yet handle $x with A, and $y with P (which should be an
498 # error).
499 my $self = shift;
500 my $a = shift; # accuracy, if given by caller
501 my $p = shift; # precision, if given by caller
502 my $r = shift; # round_mode, if given by caller
503 my @args = @_; # all 'other' arguments (0 for unary, 1 for binary ops)
504
505 unshift @args,$self; # add 'first' argument
506
507 $self = new($self) unless ref($self); # if not object, make one
508
509 # find out class of argument to round
510 my $c = ref($args[0]);
511
512 # now pick $a or $p, but only if we have got "arguments"
513 if ((!defined $a) && (!defined $p) && (@args > 0))
514 {
515 foreach (@args)
516 {
517 # take the defined one, or if both defined, the one that is smaller
518 $a = $_->{_a} if (defined $_->{_a}) && (!defined $a || $_->{_a} < $a);
519 }
520 if (!defined $a) # if it still is not defined, take p
521 {
522 foreach (@args)
523 {
524 # take the defined one, or if both defined, the one that is smaller
525 $p = $_->{_p} if (defined $_->{_p}) && (!defined $p || $_->{_p} < $p);
1f45ae4a 526 }
58cde26e
JH
527 # if none defined, use globals (#2)
528 if (!defined $p)
529 {
530 no strict 'refs';
531 my $z = "$c\::accuracy"; $a = $$z;
532 if (!defined $a)
533 {
534 $z = "$c\::precision"; $p = $$z;
535 }
1f45ae4a 536 }
58cde26e
JH
537 } # endif !$a
538 } # endif !$a || !$P && args > 0
539 # for clearity, this is not merged at place (#2)
540 # now round, by calling fround or ffround:
541 if (defined $a)
542 {
543 $self->{_a} = $a; $self->bround($a,$r);
544 }
545 elsif (defined $p)
546 {
547 $self->{_p} = $p; $self->bfround($p,$r);
548 }
549 return $self->bnorm();
550 }
551
552sub bnorm
553 {
554 # (num_str or BINT) return BINT
555 # Normalize number -- no-op here
556 my $self = shift;
557
558 return $self;
559 }
560
561sub babs
562 {
563 # (BINT or num_str) return BINT
564 # make number absolute, or return absolute BINT from string
565 #my ($self,$x) = objectify(1,@_);
566 my $x = shift; $x = $class->new($x) unless ref $x;
567 return $x if $x->modify('babs');
568 # post-normalized abs for internal use (does nothing for NaN)
569 $x->{sign} =~ s/^-/+/;
570 $x;
571 }
572
573sub bneg
574 {
575 # (BINT or num_str) return BINT
576 # negate number or make a negated number from string
577 my ($self,$x,$a,$p,$r) = objectify(1,@_);
578 return $x if $x->modify('bneg');
579 # for +0 dont negate (to have always normalized)
580 return $x if $x->is_zero();
581 $x->{sign} =~ tr/+\-/-+/; # does nothing for NaN
582 # $x->round($a,$p,$r); # changing this makes $x - $y modify $y!!
583 $x;
584 }
585
586sub bcmp
587 {
588 # Compares 2 values. Returns one of undef, <0, =0, >0. (suitable for sort)
589 # (BINT or num_str, BINT or num_str) return cond_code
590 my ($self,$x,$y) = objectify(2,@_);
591 return undef if (($x->{sign} eq $nan) || ($y->{sign} eq $nan));
592 &cmp($x->{value},$y->{value},$x->{sign},$y->{sign}) <=> 0;
593 }
594
595sub bacmp
596 {
597 # Compares 2 values, ignoring their signs.
598 # Returns one of undef, <0, =0, >0. (suitable for sort)
599 # (BINT, BINT) return cond_code
600 my ($self,$x,$y) = objectify(2,@_);
601 return undef if (($x->{sign} eq $nan) || ($y->{sign} eq $nan));
602 acmp($x->{value},$y->{value}) <=> 0;
603 }
604
605sub badd
606 {
607 # add second arg (BINT or string) to first (BINT) (modifies first)
608 # return result as BINT
609 trace(@_);
610 my ($self,$x,$y,$a,$p,$r) = objectify(2,@_);
611
612 return $x if $x->modify('badd');
613 return $x->bnan() if (($x->{sign} eq $nan) || ($y->{sign} eq $nan));
614
615 # for round calls, make array
616 my @bn = ($a,$p,$r,$y);
617 # speed: no add for 0+y or x+0
618 return $x->bnorm(@bn) if $y->is_zero(); # x+0
619 if ($x->is_zero()) # 0+y
620 {
621 # make copy, clobbering up x
622 $x->{value} = [ @{$y->{value}} ];
623 $x->{sign} = $y->{sign} || $nan;
624 return $x->round(@bn);
625 }
626
627 # shortcuts
628 my $xv = $x->{value};
629 my $yv = $y->{value};
630 my ($sx, $sy) = ( $x->{sign}, $y->{sign} ); # get signs
631
632 if ($sx eq $sy)
633 {
634 add($xv,$yv); # if same sign, absolute add
635 $x->{sign} = $sx;
636 }
637 else
638 {
639 my $a = acmp ($yv,$xv); # absolute compare
640 if ($a > 0)
641 {
642 #print "swapped sub (a=$a)\n";
643 &sub($yv,$xv,1); # absolute sub w/ swapped params
644 $x->{sign} = $sy;
645 }
646 elsif ($a == 0)
647 {
648 # speedup, if equal, set result to 0
649 $x->{value} = [ 0 ];
650 $x->{sign} = '+';
651 }
652 else # a < 0
653 {
654 #print "unswapped sub (a=$a)\n";
655 &sub($xv, $yv); # absolute sub
656 $x->{sign} = $sx;
a0d0e21e 657 }
a0d0e21e 658 }
58cde26e
JH
659 return $x->round(@bn);
660 }
661
662sub bsub
663 {
664 # (BINT or num_str, BINT or num_str) return num_str
665 # subtract second arg from first, modify first
666 my ($self,$x,$y,$a,$p,$r) = objectify(2,@_);
667
668 trace(@_);
669 return $x if $x->modify('bsub');
670 $x->badd($y->bneg()); # badd does not leave internal zeros
671 $y->bneg(); # refix y, assumes no one reads $y in between
672 return $x->round($a,$p,$r,$y);
673 }
674
675sub binc
676 {
677 # increment arg by one
678 my ($self,$x,$a,$p,$r) = objectify(1,@_);
679 # my $x = shift; $x = $class->new($x) unless ref $x; my $self = ref($x);
680 trace(@_);
681 return $x if $x->modify('binc');
682 $x->badd($self->_one())->round($a,$p,$r);
683 }
684
685sub bdec
686 {
687 # decrement arg by one
688 my ($self,$x,$a,$p,$r) = objectify(1,@_);
689 trace(@_);
690 return $x if $x->modify('bdec');
691 $x->badd($self->_one('-'))->round($a,$p,$r);
692 }
693
694sub blcm
695 {
696 # (BINT or num_str, BINT or num_str) return BINT
697 # does not modify arguments, but returns new object
698 # Lowest Common Multiplicator
699 trace(@_);
700
701 my ($self,@arg) = objectify(0,@_);
702 my $x = $self->new(shift @arg);
703 while (@arg) { $x = _lcm($x,shift @arg); }
704 $x;
705 }
706
707sub bgcd
708 {
709 # (BINT or num_str, BINT or num_str) return BINT
710 # does not modify arguments, but returns new object
711 # GCD -- Euclids algorithm, variant C (Knuth Vol 3, pg 341 ff)
712 trace(@_);
713
714 my ($self,@arg) = objectify(0,@_);
715 my $x = $self->new(shift @arg);
716 while (@arg)
717 {
718 #$x = _gcd($x,shift @arg); last if $x->is_one(); # new fast, but is slower
719 $x = _gcd_old($x,shift @arg); last if $x->is_one(); # old, slow, but faster
720 }
721 $x;
722 }
723
724sub bmod
725 {
726 # modulus
727 # (BINT or num_str, BINT or num_str) return BINT
728 my ($self,$x,$y) = objectify(2,@_);
729
730 return $x if $x->modify('bmod');
731 (&bdiv($self,$x,$y))[1];
732 }
733
734sub bnot
735 {
736 # (num_str or BINT) return BINT
737 # represent ~x as twos-complement number
738 my ($self,$x) = objectify(1,@_);
739 return $x if $x->modify('bnot');
740 $x->bneg(); $x->bdec(); # was: bsub(-1,$x);, time it someday
741 $x;
742 }
743
744sub is_zero
745 {
746 # return true if arg (BINT or num_str) is zero (array '+', '0')
747 #my ($self,$x) = objectify(1,@_);
748 #trace(@_);
749 my $x = shift; $x = $class->new($x) unless ref $x;
750 return (@{$x->{value}} == 1) && ($x->{sign} eq '+')
751 && ($x->{value}->[0] == 0);
752 }
753
754sub is_nan
755 {
756 # return true if arg (BINT or num_str) is NaN
757 #my ($self,$x) = objectify(1,@_);
758 #trace(@_);
759 my $x = shift; $x = $class->new($x) unless ref $x;
760 return ($x->{sign} eq $nan);
761 }
762
763sub is_inf
764 {
765 # return true if arg (BINT or num_str) is +-inf
766 #my ($self,$x) = objectify(1,@_);
767 #trace(@_);
768 my $x = shift; $x = $class->new($x) unless ref $x;
769 my $sign = shift || '';
770
771 return $x->{sign} =~ /^[+-]inf/ if $sign eq '';
772 return $x->{sign} =~ /^[$sign]inf/;
773 }
774
775sub is_one
776 {
777 # return true if arg (BINT or num_str) is +1 (array '+', '1')
778 # or -1 if signis given
779 #my ($self,$x) = objectify(1,@_);
780 my $x = shift; $x = $class->new($x) unless ref $x;
781 my $sign = shift || '+'; #$_[2] || '+';
782 return (@{$x->{value}} == 1) && ($x->{sign} eq $sign)
783 && ($x->{value}->[0] == 1);
784 }
785
786sub is_odd
787 {
788 # return true when arg (BINT or num_str) is odd, false for even
789 my $x = shift; $x = $class->new($x) unless ref $x;
790 #my ($self,$x) = objectify(1,@_);
791 return (($x->{sign} ne $nan) && ($x->{value}->[0] & 1));
792 }
793
794sub is_even
795 {
796 # return true when arg (BINT or num_str) is even, false for odd
797 my $x = shift; $x = $class->new($x) unless ref $x;
798 #my ($self,$x) = objectify(1,@_);
799 return (($x->{sign} ne $nan) && (!($x->{value}->[0] & 1)));
800 }
801
802sub bmul
803 {
804 # multiply two numbers -- stolen from Knuth Vol 2 pg 233
805 # (BINT or num_str, BINT or num_str) return BINT
806 my ($self,$x,$y,$a,$p,$r) = objectify(2,@_);
807 #print "$self bmul $x ",ref($x)," $y ",ref($y),"\n";
808 trace(@_);
809 return $x if $x->modify('bmul');
810 return $x->bnan() if (($x->{sign} eq $nan) || ($y->{sign} eq $nan));
811
812 mul($x,$y); # do actual math
813 return $x->round($a,$p,$r,$y);
814 }
815
816sub bdiv
817 {
818 # (dividend: BINT or num_str, divisor: BINT or num_str) return
819 # (BINT,BINT) (quo,rem) or BINT (only rem)
820 trace(@_);
821 my ($self,$x,$y,$a,$p,$r) = objectify(2,@_);
822
823 return $x if $x->modify('bdiv');
824
825 # NaN?
826 return wantarray ? ($x->bnan(),bnan()) : $x->bnan()
827 if ($x->{sign} eq $nan || $y->{sign} eq $nan || $y->is_zero());
828
829 # 0 / something
830 return wantarray ? ($x,$self->bzero()) : $x if $x->is_zero();
831
832 # Is $x in the interval [0, $y) ?
833 my $cmp = acmp($x->{value},$y->{value});
834 if (($cmp < 0) and ($x->{sign} eq $y->{sign}))
835 {
836 return $x->bzero() unless wantarray;
837 my $t = $x->copy(); # make copy first, because $x->bzero() clobbers $x
838 return ($x->bzero(),$t);
839 }
840 elsif ($cmp == 0)
841 {
842 # shortcut, both are the same, so set to +/- 1
843 $x->_one( ($x->{sign} ne $y->{sign} ? '-' : '+') );
844 return $x unless wantarray;
845 return ($x,$self->bzero());
846 }
847
848 # calc new sign and in case $y == +/- 1, return $x
849 $x->{sign} = ($x->{sign} ne $y->{sign} ? '-' : '+');
850 # check for / +-1 (cant use $y->is_one due to '-'
851 if ((@{$y->{value}} == 1) && ($y->{value}->[0] == 1))
852 {
853 return wantarray ? ($x,$self->bzero()) : $x;
854 }
855
856 # call div here
857 my $rem = $self->bzero();
858 $rem->{sign} = $y->{sign};
859 ($x->{value},$rem->{value}) = div($x->{value},$y->{value});
860 # do not leave rest "-0";
861 $rem->{sign} = '+' if (@{$rem->{value}} == 1) && ($rem->{value}->[0] == 0);
862 if (($x->{sign} eq '-') and (!$rem->is_zero()))
863 {
864 $x->bdec();
865 }
866 $x->round($a,$p,$r,$y);
867 if (wantarray)
868 {
869 $rem->round($a,$p,$r,$x,$y);
870 return ($x,$y-$rem) if $x->{sign} eq '-'; # was $x,$rem
871 return ($x,$rem);
872 }
873 return $x;
874 }
875
876sub bpow
877 {
878 # (BINT or num_str, BINT or num_str) return BINT
879 # compute power of two numbers -- stolen from Knuth Vol 2 pg 233
880 # modifies first argument
881 #trace(@_);
882 my ($self,$x,$y,$a,$p,$r) = objectify(2,@_);
883
884 return $x if $x->modify('bpow');
885
886 return $x->bnan() if $x->{sign} eq $nan || $y->{sign} eq $nan;
887 return $x->_one() if $y->is_zero();
888 return $x if $x->is_one() || $y->is_one();
889 if ($x->{sign} eq '-' && @{$x->{value}} == 1 && $x->{value}->[0] == 1)
890 {
891 # if $x == -1 and odd/even y => +1/-1
892 return $y->is_odd() ? $x : $x->_set(1); # $x->babs() would work to
893 # my Casio FX-5500L has here a bug, -1 ** 2 is -1, but -1 * -1 is 1 LOL
894 }
895 # shortcut for $x ** 2
896 if ($y->{sign} eq '+' && @{$y->{value}} == 1 && $y->{value}->[0] == 2)
897 {
898 return $x->bmul($x)->bround($a,$p,$r);
899 }
900 # 1 ** -y => 1 / (1**y), so do test for negative $y after above's clause
901 return $x->bnan() if $y->{sign} eq '-';
902 return $x if $x->is_zero(); # 0**y => 0 (if not y <= 0)
903
904 # tels: 10**x is special (actually 100**x etc is special, too) but not here
905 #if ((@{$x->{value}} == 1) && ($x->{value}->[0] == 10))
906 # {
907 # # 10**2
908 # my $yi = int($y); my $yi5 = int($yi/5);
909 # $x->{value} = [];
910 # my $v = $x->{value};
911 # if ($yi5 > 0)
912 # {
913 # # $x->{value}->[$yi5-1] = 0; # pre-padd array (no use)
914 # for (my $i = 0; $i < $yi5; $i++)
915 # {
916 # $v->[$i] = 0;
917 # }
918 # }
919 # push @{$v}, int( '1'.'0' x ($yi % 5));
920 # if ($x->{sign} eq '-')
921 # {
922 # $x->{sign} = $y->is_odd() ? '-' : '+'; # -10**2 = 100, -10**3 = -1000
923 # }
924 # return $x;
925 # }
926
927 # based on the assumption that shifting in base 10 is fast, and that bpow()
928 # works faster if numbers are small: we count trailing zeros (this step is
929 # O(1)..O(N), but in case of O(N) we save much more time), stripping them
930 # out of the multiplication, and add $count * $y zeros afterwards:
931 # 300 ** 3 == 300*300*300 == 3*3*3 . '0' x 2 * 3 == 27 . '0' x 6
932 my $zeros = $x->_trailing_zeros();
933 if ($zeros > 0)
934 {
935 $x->brsft($zeros,10); # remove zeros
936 $x->bpow($y); # recursion (will not branch into here again)
937 $zeros = $y * $zeros; # real number of zeros to add
938 $x->blsft($zeros,10);
939 return $x;
940 }
941
942 my $pow2 = $self->_one();
943 my $y1 = $class->new($y);
944 my ($res);
945 while (!$y1->is_one())
946 {
947 #print "bpow: p2: $pow2 x: $x y: $y1 r: $res\n";
948 #print "len ",$x->length(),"\n";
949 ($y1,$res)=&bdiv($y1,2);
950 if (!$res->is_zero()) { &bmul($pow2,$x); }
951 if (!$y1->is_zero()) { &bmul($x,$x); }
952 }
953 #print "bpow: e p2: $pow2 x: $x y: $y1 r: $res\n";
954 &bmul($x,$pow2) if (!$pow2->is_one());
955 #print "bpow: e p2: $pow2 x: $x y: $y1 r: $res\n";
956 return $x->round($a,$p,$r);
957 }
958
959sub blsft
960 {
961 # (BINT or num_str, BINT or num_str) return BINT
962 # compute x << y, base n, y >= 0
963 my ($self,$x,$y,$n) = objectify(2,@_);
964
965 return $x if $x->modify('blsft');
966 return $x->bnan() if ($x->{sign} !~ /^[+-]$/ || $y->{sign} !~ /^[+-]$/);
967
968 $n = 2 if !defined $n; return $x if $n == 0;
969 return $x->bnan() if $n < 0 || $y->{sign} eq '-';
970 if ($n != 10)
971 {
972 $x->bmul( $self->bpow($n, $y) );
973 }
974 else
975 {
976 # shortcut (faster) for shifting by 10) since we are in base 10eX
977 # multiples of 5:
978 my $src = scalar @{$x->{value}}; # source
979 my $len = $y->numify(); # shift-len as normal int
980 my $rem = $len % 5; # reminder to shift
981 my $dst = $src + int($len/5); # destination
982
983 my $v = $x->{value}; # speed-up
984 my $vd; # further speedup
985 #print "src $src:",$v->[$src]||0," dst $dst:",$v->[$dst]||0," rem $rem\n";
986 $v->[$src] = 0; # avoid first ||0 for speed
987 while ($src >= 0)
988 {
989 $vd = $v->[$src]; $vd = '00000'.$vd;
990 #print "s $src d $dst '$vd' ";
991 $vd = substr($vd,-5+$rem,5-$rem);
992 #print "'$vd' ";
993 $vd .= $src > 0 ? substr('00000'.$v->[$src-1],-5,$rem) : '0' x $rem;
994 #print "'$vd' ";
995 $vd = substr($vd,-5,5) if length($vd) > 5;
996 #print "'$vd'\n";
997 $v->[$dst] = int($vd);
998 $dst--; $src--;
999 }
1000 # set lowest parts to 0
1001 while ($dst >= 0) { $v->[$dst--] = 0; }
1002 # fix spurios last zero element
1003 splice @$v,-1 if $v->[-1] == 0;
1004 #print "elems: "; my $i = 0;
1005 #foreach (reverse @$v) { print "$i $_ "; $i++; } print "\n";
1006 # old way: $x->bmul( $self->bpow($n, $y) );
1007 }
1008 return $x;
1009 }
1010
1011sub brsft
1012 {
1013 # (BINT or num_str, BINT or num_str) return BINT
1014 # compute x >> y, base n, y >= 0
1015 my ($self,$x,$y,$n) = objectify(2,@_);
1016
1017 return $x if $x->modify('brsft');
1018 return $x->bnan() if ($x->{sign} !~ /^[+-]$/ || $y->{sign} !~ /^[+-]$/);
1019
1020 $n = 2 if !defined $n; return $x->bnan() if $n <= 0 || $y->{sign} eq '-';
1021 if ($n != 10)
1022 {
1023 scalar bdiv($x, $self->bpow($n, $y));
1024 }
1025 else
1026 {
1027 # shortcut (faster) for shifting by 10)
1028 # multiples of 5:
1029 my $dst = 0; # destination
1030 my $src = $y->numify(); # as normal int
1031 my $rem = $src % 5; # reminder to shift
1032 $src = int($src / 5); # source
1033 my $len = scalar @{$x->{value}} - $src; # elems to go
1034 my $v = $x->{value}; # speed-up
1035 if ($rem == 0)
1036 {
1037 splice (@$v,0,$src); # even faster, 38.4 => 39.3
1038 }
1039 else
1040 {
1041 my $vd;
1042 $v->[scalar @$v] = 0; # avoid || 0 test inside loop
1043 while ($dst < $len)
1044 {
1045 $vd = '00000'.$v->[$src];
1046 #print "$dst $src '$vd' ";
1047 $vd = substr($vd,-5,5-$rem);
1048 #print "'$vd' ";
1049 $src++;
1050 $vd = substr('00000'.$v->[$src],-$rem,$rem) . $vd;
1051 #print "'$vd1' ";
1052 #print "'$vd'\n";
1053 $vd = substr($vd,-5,5) if length($vd) > 5;
1054 $v->[$dst] = int($vd);
1055 $dst++;
1056 }
1057 splice (@$v,$dst) if $dst > 0; # kill left-over array elems
1058 pop @$v if $v->[-1] == 0; # kill last element
1059 } # else rem == 0
1060 # old way: scalar bdiv($x, $self->bpow($n, $y));
1061 }
1062 return $x;
1063 }
1064
1065sub band
1066 {
1067 #(BINT or num_str, BINT or num_str) return BINT
1068 # compute x & y
1069 trace(@_);
1070 my ($self,$x,$y) = objectify(2,@_);
1071
1072 return $x if $x->modify('band');
1073
1074 return $x->bnan() if ($x->{sign} !~ /^[+-]$/ || $y->{sign} !~ /^[+-]$/);
1075 return $x->bzero() if $y->is_zero();
1076 my $r = $self->bzero(); my $m = new Math::BigInt 1; my ($xr,$yr);
1077 my $x10000 = new Math::BigInt (0x10000);
1078 my $y1 = copy(ref($x),$y); # make copy
1079 while (!$x->is_zero() && !$y1->is_zero())
1080 {
1081 ($x, $xr) = bdiv($x, $x10000);
1082 ($y1, $yr) = bdiv($y1, $x10000);
1083 $r->badd( bmul( new Math::BigInt ( $xr->numify() & $yr->numify()), $m ));
1084 $m->bmul($x10000);
1085 }
1086 $x = $r;
1087 }
1088
1089sub bior
1090 {
1091 #(BINT or num_str, BINT or num_str) return BINT
1092 # compute x | y
1093 trace(@_);
1094 my ($self,$x,$y) = objectify(2,@_);
1095
1096 return $x if $x->modify('bior');
1097
1098 return $x->bnan() if ($x->{sign} !~ /^[+-]$/ || $y->{sign} !~ /^[+-]$/);
1099 return $x if $y->is_zero();
1100 my $r = $self->bzero(); my $m = new Math::BigInt 1; my ($xr,$yr);
1101 my $x10000 = new Math::BigInt (0x10000);
1102 my $y1 = copy(ref($x),$y); # make copy
1103 while (!$x->is_zero() || !$y1->is_zero())
1104 {
1105 ($x, $xr) = bdiv($x,$x10000);
1106 ($y1, $yr) = bdiv($y1,$x10000);
1107 $r->badd( bmul( new Math::BigInt ( $xr->numify() | $yr->numify()), $m ));
1108 $m->bmul($x10000);
1109 }
1110 $x = $r;
1111 }
1112
1113sub bxor
1114 {
1115 #(BINT or num_str, BINT or num_str) return BINT
1116 # compute x ^ y
1117 my ($self,$x,$y) = objectify(2,@_);
1118
1119 return $x if $x->modify('bxor');
1120
1121 return $x->bnan() if ($x->{sign} eq $nan || $y->{sign} eq $nan);
1122 return $x if $y->is_zero();
1123 return $x->bzero() if $x == $y; # shortcut
1124 my $r = $self->bzero(); my $m = new Math::BigInt 1; my ($xr,$yr);
1125 my $x10000 = new Math::BigInt (0x10000);
1126 my $y1 = copy(ref($x),$y); # make copy
1127 while (!$x->is_zero() || !$y1->is_zero())
1128 {
1129 ($x, $xr) = bdiv($x, $x10000);
1130 ($y1, $yr) = bdiv($y1, $x10000);
1131 $r->badd( bmul( new Math::BigInt ( $xr->numify() ^ $yr->numify()), $m ));
1132 $m->bmul($x10000);
1133 }
1134 $x = $r;
1135 }
1136
1137sub length
1138 {
1139 my ($self,$x) = objectify(1,@_);
1140
1141 return (_digits($x->{value}), 0) if wantarray;
1142 _digits($x->{value});
1143 }
1144
1145sub digit
1146 {
1147 # return the nth digit, negative values count backward
1148 my $x = shift;
1149 my $n = shift || 0;
1150
1151 my $len = $x->length();
1152
1153 $n = $len+$n if $n < 0; # -1 last, -2 second-to-last
1154 $n = abs($n); # if negatives are to big
1155 $len--; $n = $len if $n > $len; # n to big?
1156
1157 my $elem = int($n / 5); # which array element
1158 my $digit = $n % 5; # which digit in this element
1159 $elem = '0000'.$x->{value}->[$elem]; # get element padded with 0's
1160 return substr($elem,-$digit-1,1);
1161 }
1162
1163sub _trailing_zeros
1164 {
1165 # return the amount of trailing zeros in $x
1166 my $x = shift;
1167 $x = $class->new($x) unless ref $x;
1168
1169 return 0 if $x->is_zero() || $x->is_nan();
1170 # check each array elem in _m for having 0 at end as long as elem == 0
1171 # Upon finding a elem != 0, stop
1172 my $zeros = 0; my $elem;
1173 foreach my $e (@{$x->{value}})
1174 {
1175 if ($e != 0)
1176 {
1177 $elem = "$e"; # preserve x
1178 $elem =~ s/.*?(0*$)/$1/; # strip anything not zero
1179 $zeros *= 5; # elems * 5
1180 $zeros += CORE::length($elem); # count trailing zeros
1181 last; # early out
1182 }
1183 $zeros ++; # real else branch: 50% slower!
1184 }
1185 return $zeros;
1186 }
1187
1188sub bsqrt
1189 {
1190 my ($self,$x) = objectify(1,@_);
1191
1192 return $x->bnan() if $x->{sign} =~ /\-|$nan/; # -x or NaN => NaN
1193 return $x->bzero() if $x->is_zero(); # 0 => 0
1194 return $x if $x == 1; # 1 => 1
1195
1196 my $y = $x->copy(); # give us one more digit accur.
1197 my $l = int($x->length()/2);
1198
1199 $x->bzero();
1200 $x->binc(); # keep ref($x), but modify it
1201 $x *= 10 ** $l;
1202
1203 # print "x: $y guess $x\n";
1204
1205 my $last = $self->bzero();
1206 while ($last != $x)
1207 {
1208 $last = $x;
1209 $x += $y / $x;
1210 $x /= 2;
1211 }
1212 return $x;
1213 }
1214
1215sub exponent
1216 {
1217 # return a copy of the exponent (here always 0, NaN or 1 for $m == 0)
1218 my ($self,$x) = objectify(1,@_);
1219
1220 return bnan() if $x->is_nan();
1221 my $e = $class->bzero();
1222 return $e->binc() if $x->is_zero();
1223 $e += $x->_trailing_zeros();
1224 return $e;
1225 }
1226
1227sub mantissa
1228 {
1229 # return a copy of the mantissa (here always $self)
1230 my ($self,$x) = objectify(1,@_);
1231
1232 return bnan() if $x->is_nan();
1233 my $m = $x->copy();
1234 # that's inefficient
1235 my $zeros = $m->_trailing_zeros();
1236 $m /= 10 ** $zeros if $zeros != 0;
1237 return $m;
1238 }
1239
1240sub parts
1241 {
1242 # return a copy of both the exponent and the mantissa (here 0 and self)
1243 my $self = shift;
1244 $self = $class->new($self) unless ref $self;
1245
1246 return ($self->mantissa(),$self->exponent());
1247 }
1248
1249##############################################################################
1250# rounding functions
1251
1252sub bfround
1253 {
1254 # precision: round to the $Nth digit left (+$n) or right (-$n) from the '.'
1255 # $n == 0 => round to integer
1256 my $x = shift; $x = $class->new($x) unless ref $x;
1257 my ($scale,$mode) = $x->_scale_p($precision,$rnd_mode,@_);
1258 return $x if !defined $scale; # no-op
1259
1260 # no-op for BigInts if $n <= 0
1261 return $x if $scale <= 0;
1262
1263 $x->bround( $x->length()-$scale, $mode);
1264 }
1265
1266sub _scan_for_nonzero
1267 {
1268 my $x = shift;
1269 my $pad = shift;
1270
1271 my $len = $x->length();
1272 return 0 if $len == 1; # '5' is trailed by invisible zeros
1273 my $follow = $pad - 1;
1274 return 0 if $follow > $len || $follow < 1;
1275 #print "checking $x $r\n";
1276 # old, slow way checking string for non-zero characters
1277 my $r = substr ("$x",-$follow);
1278 return 1 if $r =~ /[^0]/; return 0;
1279
1280 # faster way checking array contents; it is actually not faster (even in a
1281 # rounding-only-shoutout, so I leave the simpler code in)
1282 #my $rem = $follow % 5; my $div = $follow / 5; my $v = $x->{value};
1283 # pad with zeros and extract
1284 #print "last part : ",'00000'.$v->[$div]," $rem = '";
1285 #print substr('00000'.$v->[$div],-$rem,5),"'\n";
1286 #my $r1 = substr ('00000'.$v->[$div],-$rem,5);
1287 #print "$r1\n";
1288 #return 1 if $r1 =~ /[^0]/;
1289 #
1290 #for (my $j = $div-1; $j >= 0; $j --)
1291 # {
1292 # #print "part $v->[$j]\n";
1293 # return 1 if $v->[$j] != 0;
1294 # }
1295 #return 0;
1296 }
1297
1298sub fround
1299 {
1300 # to make life easier for switch between MBF and MBI (autoload fxxx()
1301 # like MBF does for bxxx()?)
1302 my $x = shift;
1303 return $x->bround(@_);
1304 }
1305
1306sub bround
1307 {
1308 # accuracy: +$n preserve $n digits from left,
1309 # -$n preserve $n digits from right (f.i. for 0.1234 style in MBF)
1310 # no-op for $n == 0
1311 # and overwrite the rest with 0's, return normalized number
1312 # do not return $x->bnorm(), but $x
1313 my $x = shift; $x = $class->new($x) unless ref $x;
1314 my ($scale,$mode) = $x->_scale_a($accuracy,$rnd_mode,@_);
1315 return $x if !defined $scale; # no-op
1316
1317 # print "MBI round: $x to $scale $mode\n";
1318 # -scale means what? tom? hullo? -$scale needed by MBF round, but what for?
1319 return $x if $x->is_nan() || $x->is_zero() || $scale == 0;
1320
1321 # we have fewer digits than we want to scale to
1322 my $len = $x->length();
1323 # print "$len $scale\n";
1324 return $x if $len < abs($scale);
1325
1326 # count of 0's to pad, from left (+) or right (-): 9 - +6 => 3, or |-6| => 6
1327 my ($pad,$digit_round,$digit_after);
1328 $pad = $len - $scale;
1329 $pad = abs($scale)+1 if $scale < 0;
1330 $digit_round = '0'; $digit_round = $x->digit($pad) if $pad < $len;
1331 $digit_after = '0'; $digit_after = $x->digit($pad-1) if $pad > 0;
1332 # print "r $x: pos:$pad l:$len s:$scale r:$digit_round a:$digit_after m: $mode\n";
1333
1334 # in case of 01234 we round down, for 6789 up, and only in case 5 we look
1335 # closer at the remaining digits of the original $x, remember decision
1336 my $round_up = 1; # default round up
1337 $round_up -- if
1338 ($mode eq 'trunc') || # trunc by round down
1339 ($digit_after =~ /[01234]/) || # round down anyway,
1340 # 6789 => round up
1341 ($digit_after eq '5') && # not 5000...0000
1342 ($x->_scan_for_nonzero($pad) == 0) &&
1343 (
1344 ($mode eq 'even') && ($digit_round =~ /[24680]/) ||
1345 ($mode eq 'odd') && ($digit_round =~ /[13579]/) ||
1346 ($mode eq '+inf') && ($x->{sign} eq '-') ||
1347 ($mode eq '-inf') && ($x->{sign} eq '+') ||
1348 ($mode eq 'zero') # round down if zero, sign adjusted below
1349 );
1350 # allow rounding one place left of mantissa
1351 #print "$pad $len $scale\n";
1352 # this is triggering warnings, and buggy for $scale < 0
1353 #if (-$scale != $len)
1354 {
1355 # split mantissa at $scale and then pad with zeros
1356 my $s5 = int($pad / 5);
1357 my $i = 0;
1358 while ($i < $s5)
1359 {
1360 $x->{value}->[$i++] = 0; # replace with 5 x 0
1361 }
1362 $x->{value}->[$s5] = '00000'.$x->{value}->[$s5]; # pad with 0
1363 my $rem = $pad % 5; # so much left over
1364 if ($rem > 0)
1365 {
1366 #print "remainder $rem\n";
1367 #print "elem $x->{value}->[$s5]\n";
1368 substr($x->{value}->[$s5],-$rem,$rem) = '0' x $rem; # stamp w/ '0'
1369 }
1370 $x->{value}->[$s5] = int ($x->{value}->[$s5]); # str '05' => int '5'
1371 }
1372 if ($round_up) # what gave test above?
1373 {
1374 $pad = $len if $scale < 0; # tlr: whack 0.51=>1.0
1375 # modify $x in place, undef, undef to avoid rounding
1376 $x->badd( Math::BigInt->new($x->{sign}.'1'.'0'x$pad),
1377 undef,undef );
1378 # str creation much faster than 10 ** something
1379 }
1380 $x;
1381 }
1382
1383sub bfloor
1384 {
1385 # return integer less or equal then number, since it is already integer,
1386 # always returns $self
1387 my ($self,$x,$a,$p,$r) = objectify(1,@_);
1388
1389 # not needed: return $x if $x->modify('bfloor');
1390
1391 return $x->round($a,$p,$r);
1392 }
1393
1394sub bceil
1395 {
1396 # return integer greater or equal then number, since it is already integer,
1397 # always returns $self
1398 my ($self,$x,$a,$p,$r) = objectify(1,@_);
1399
1400 # not needed: return $x if $x->modify('bceil');
1401
1402 return $x->round($a,$p,$r);
1403 }
1404
1405##############################################################################
1406# private stuff (internal use only)
1407
1408sub trace
1409 {
1410 # print out a number without using bstr (avoid deep recurse) for trace/debug
1411 return unless $trace;
1412
1413 my ($package,$file,$line,$sub) = caller(1);
1414 print "'$sub' called from '$package' line $line:\n ";
1415
1416 foreach my $x (@_)
1417 {
1418 if (!defined $x)
1419 {
1420 print "undef, "; next;
1421 }
1422 if (!ref($x))
1423 {
1424 print "'$x' "; next;
1425 }
1426 next if (ref($x) ne "HASH");
1427 print "$x->{sign} ";
1428 foreach (@{$x->{value}})
1429 {
1430 print "$_ ";
1431 }
1432 print ", ";
1433 }
1434 print "\n";
1435 }
1436
1437sub _set
1438 {
1439 # internal set routine to set X fast to an integer value < [+-]100000
1440 my $self = shift;
1441 my $wanted = shift || 0;
1442
1443 $self->{sign} = $nan, return if $wanted !~ /^[+-]?[0-9]+$/;
1444 $self->{sign} = '-'; $self->{sign} = '+' if $wanted >= 0;
1445 $self->{value} = [ abs($wanted) ];
1446 return $self;
1447 }
1448
1449sub _one
1450 {
1451 # internal speedup, set argument to 1, or create a +/- 1
1452 my $self = shift;
1453 my $x = $self->bzero(); $x->{value} = [ 1 ]; $x->{sign} = shift || '+'; $x;
1454 }
1455
1456sub _swap
1457 {
1458 # Overload will swap params if first one is no object ref so that the first
1459 # one is always an object ref. In this case, third param is true.
1460 # This routine is to overcome the effect of scalar,$object creating an object
1461 # of the class of this package, instead of the second param $object. This
1462 # happens inside overload, when the overload section of this package is
1463 # inherited by sub classes.
1464 # For overload cases (and this is used only there), we need to preserve the
1465 # args, hence the copy().
1466 # You can override this method in a subclass, the overload section will call
1467 # $object->_swap() to make sure it arrives at the proper subclass, with some
1468 # exceptions like '+' and '-'.
1469
1470 # object, (object|scalar) => preserve first and make copy
1471 # scalar, object => swapped, re-swap and create new from first
1472 # (using class of second object, not $class!!)
1473 my $self = shift; # for override in subclass
1474 #print "swap $self 0:$_[0] 1:$_[1] 2:$_[2]\n";
1475 if ($_[2])
1476 {
1477 my $c = ref ($_[0]) || $class; # fallback $class should not happen
1478 return ( $c->new($_[1]), $_[0] );
1479 }
1480 else
1481 {
1482 return ( $_[0]->copy(), $_[1] );
1483 }
1484 }
1485
1486sub objectify
1487 {
1488 # check for strings, if yes, return objects instead
1489
1490 # the first argument is number of args objectify() should look at it will
1491 # return $count+1 elements, the first will be a classname. This is because
1492 # overloaded '""' calls bstr($object,undef,undef) and this would result in
1493 # useless objects beeing created and thrown away. So we cannot simple loop
1494 # over @_. If the given count is 0, all arguments will be used.
1495
1496 # If the second arg is a ref, use it as class.
1497 # If not, try to use it as classname, unless undef, then use $class
1498 # (aka Math::BigInt). The latter shouldn't happen,though.
1499
1500 # caller: gives us:
1501 # $x->badd(1); => ref x, scalar y
1502 # Class->badd(1,2); => classname x (scalar), scalar x, scalar y
1503 # Class->badd( Class->(1),2); => classname x (scalar), ref x, scalar y
1504 # Math::BigInt::badd(1,2); => scalar x, scalar y
1505 # In the last case we check number of arguments to turn it silently into
1506 # $class,1,2. (We can not take '1' as class ;o)
1507 # badd($class,1) is not supported (it should, eventually, try to add undef)
1508 # currently it tries 'Math::BigInt' + 1, which will not work.
1509
1510 trace(@_);
1511 my $count = abs(shift || 0);
1512
1513 #print caller(),"\n";
1514
1515 my @a; # resulting array
1516 if (ref $_[0])
1517 {
1518 # okay, got object as first
1519 $a[0] = ref $_[0];
1520 }
1521 else
1522 {
1523 # nope, got 1,2 (Class->xxx(1) => Class,1 and not supported)
1524 $a[0] = $class;
1525 #print "@_\n"; sleep(1);
1526 $a[0] = shift if $_[0] =~ /^[A-Z].*::/; # classname as first?
1527 }
1528 #print caller(),"\n";
1529 # print "Now in objectify, my class is today $a[0]\n";
1530 my $k;
1531 if ($count == 0)
1532 {
1533 while (@_)
1534 {
1535 $k = shift;
1536 if (!ref($k))
1537 {
1538 $k = $a[0]->new($k);
1539 }
1540 elsif (ref($k) ne $a[0])
1541 {
1542 # foreign object, try to convert to integer
1543 $k->can('as_number') ? $k = $k->as_number() : $k = $a[0]->new($k);
e16b8f49 1544 }
58cde26e
JH
1545 push @a,$k;
1546 }
1547 }
1548 else
1549 {
1550 while ($count > 0)
1551 {
1552 #print "$count\n";
1553 $count--;
1554 $k = shift;
1555 if (!ref($k))
1556 {
1557 $k = $a[0]->new($k);
1558 }
1559 elsif (ref($k) ne $a[0])
1560 {
1561 # foreign object, try to convert to integer
1562 $k->can('as_number') ? $k = $k->as_number() : $k = $a[0]->new($k);
e16b8f49 1563 }
58cde26e
JH
1564 push @a,$k;
1565 }
1566 push @a,@_; # return other params, too
1567 }
1568 #my $i = 0;
1569 #foreach (@a)
1570 # {
1571 # print "o $i $a[0]\n" if $i == 0;
1572 # print "o $i ",ref($_),"\n" if $i != 0; $i++;
1573 # }
1574 #print "objectify done: would return ",scalar @a," values\n";
1575 #print caller(1),"\n" unless wantarray;
1576 die "$class objectify needs list context" unless wantarray;
1577 @a;
1578 }
1579
1580sub import
1581 {
1582 my $self = shift;
1583 #print "import $self @_\n";
1584 for ( my $i = 0; $i < @_ ; $i++ )
1585 {
1586 if ( $_[$i] eq ':constant' )
1587 {
1588 # this rest causes overlord er load to step in
1589 overload::constant integer => sub { $self->new(shift) };
1590 splice @_, $i, 1; last;
1591 }
1592 }
1593 # any non :constant stuff is handled by our parent, Exporter
1594 # even if @_ is empty, to give it a chance
1595 #$self->SUPER::import(@_); # does not work
1596 $self->export_to_level(1,$self,@_); # need this instead
1597 }
1598
1599sub _internal
1600 {
1601 # (ref to self, ref to string) return ref to num_array
1602 # Convert a number from string format to internal base 100000 format.
1603 # Assumes normalized value as input.
1604 my ($s,$d) = @_;
1605 my $il = CORE::length($$d)-1;
1606 # these leaves '00000' instead of int 0 and will be corrected after any op
1607 $s->{value} = [ reverse(unpack("a" . ($il%5+1) . ("a5" x ($il/5)), $$d)) ];
1608 $s;
1609 }
1610
1611sub _strip_zeros
1612 {
1613 # internal normalization function that strips leading zeros from the array
1614 # args: ref to array
1615 #trace(@_);
1616 my $s = shift;
1617
1618 my $cnt = scalar @$s; # get count of parts
1619 my $i = $cnt-1;
1620 #print "strip: cnt $cnt i $i\n";
1621 # '0', '3', '4', '0', '0',
1622 # 0 1 2 3 4
1623 # cnt = 5, i = 4
1624 # i = 4
1625 # i = 3
1626 # => fcnt = cnt - i (5-2 => 3, cnt => 5-1 = 4, throw away from 4th pos)
1627 # >= 1: skip first part (this can be zero)
1628 while ($i > 0) { last if $s->[$i] != 0; $i--; }
1629 $i++; splice @$s,$i if ($i < $cnt); # $i cant be 0
1630 return $s;
1631 }
1632
1633sub _from_hex
1634 {
1635 # convert a (ref to) big hex string to BigInt, return undef for error
1636 my $hs = shift;
1637
1638 my $x = Math::BigInt->bzero();
1639 return $x->bnan() if $$hs !~ /^[\-\+]?0x[0-9A-Fa-f]+$/;
1640
1641 my $mul = Math::BigInt->bzero(); $mul++;
1642 my $x65536 = Math::BigInt->new(65536);
1643
1644 my $len = CORE::length($$hs)-2; my $sign = '+';
1645 if ($$hs =~ /^\-/)
1646 {
1647 $sign = '-'; $len--;
1648 }
1649 $len = int($len/4); # 4-digit parts, w/o '0x'
1650 my $val; my $i = -4;
1651 while ($len >= 0)
1652 {
1653 $val = substr($$hs,$i,4);
1654 $val =~ s/^[\-\+]?0x// if $len == 0; # for last part only because
1655 $val = hex($val); # hex does not like wrong chars
1656 # print "$val ",substr($$hs,$i,4),"\n";
1657 $i -= 4; $len --;
1658 $x += $mul * $val if $val != 0;
1659 $mul *= $x65536 if $len >= 0; # skip last mul
1660 }
1661 $x->{sign} = $sign if !$x->is_zero();
1662 return $x;
1663 }
1664
1665sub _from_bin
1666 {
1667 # convert a (ref to) big binary string to BigInt, return undef for error
1668 my $bs = shift;
1669
1670 my $x = Math::BigInt->bzero();
1671 return $x->bnan() if $$bs !~ /^[\-\+]?0b[01]+$/;
1672
1673 my $mul = Math::BigInt->bzero(); $mul++;
1674 my $x256 = Math::BigInt->new(256);
1675
1676 my $len = CORE::length($$bs)-2; my $sign = '+';
1677 if ($$bs =~ /^\-/)
1678 {
1679 $sign = '-'; $len--;
1680 }
1681 $len = int($len/8); # 8-digit parts, w/o '0b'
1682 my $val; my $i = -8;
1683 while ($len >= 0)
1684 {
1685 $val = substr($$bs,$i,8);
1686 $val =~ s/^[\-\+]?0b// if $len == 0; # for last part only
1687 #$val = oct('0b'.$val); # does not work on Perl prior 5.6.0
1688 $val = ('0' x (8-CORE::length($val))).$val if CORE::length($val) < 8;
1689 $val = ord(pack('B8',$val));
1690 # print "$val ",substr($$bs,$i,16),"\n";
1691 $i -= 8; $len --;
1692 $x += $mul * $val if $val != 0;
1693 $mul *= $x256 if $len >= 0; # skip last mul
1694 }
1695 $x->{sign} = $sign if !$x->is_zero();
1696 return $x;
1697 }
1698
1699sub _split
1700 {
1701 # (ref to num_str) return num_str
1702 # internal, take apart a string and return the pieces
1703 my $x = shift;
1704
1705 # pre-parse input
1706 $$x =~ s/^\s+//g; # strip white space at front
1707 $$x =~ s/\s+$//g; # strip white space at end
1708 #$$x =~ s/\s+//g; # strip white space (no longer)
1709 return if $$x eq "";
1710
1711 return _from_hex($x) if $$x =~ /^[\-\+]?0x/; # hex string
1712 return _from_bin($x) if $$x =~ /^[\-\+]?0b/; # binary string
1713
1714 return if $$x !~ /^[\-\+]?\.?[0-9]/;
1715
1716 $$x =~ s/(\d)_(\d)/$1$2/g; # strip underscores between digits
1717 $$x =~ s/(\d)_(\d)/$1$2/g; # do twice for 1_2_3
1718
1719 # some possible inputs:
1720 # 2.1234 # 0.12 # 1 # 1E1 # 2.134E1 # 434E-10 # 1.02009E-2
1721 # .2 # 1_2_3.4_5_6 # 1.4E1_2_3 # 1e3 # +.2
1722
1723 #print "input: '$$x' ";
1724 my ($m,$e) = split /[Ee]/,$$x;
1725 $e = '0' if !defined $e || $e eq "";
1726 # print "m '$m' e '$e'\n";
1727 # sign,value for exponent,mantint,mantfrac
1728 my ($es,$ev,$mis,$miv,$mfv);
1729 # valid exponent?
1730 if ($e =~ /^([+-]?)0*(\d+)$/) # strip leading zeros
1731 {
1732 $es = $1; $ev = $2;
1733 #print "'$m' '$e' e: $es $ev ";
1734 # valid mantissa?
1735 return if $m eq '.' || $m eq '';
1736 my ($mi,$mf) = split /\./,$m;
1737 $mi = '0' if !defined $mi;
1738 $mi .= '0' if $mi =~ /^[\-\+]?$/;
1739 $mf = '0' if !defined $mf || $mf eq '';
1740 if ($mi =~ /^([+-]?)0*(\d+)$/) # strip leading zeros
1741 {
1742 $mis = $1||'+'; $miv = $2;
1743 #print "$mis $miv";
1744 # valid, existing fraction part of mantissa?
1745 return unless ($mf =~ /^(\d*?)0*$/); # strip trailing zeros
1746 $mfv = $1;
1747 #print " split: $mis $miv . $mfv E $es $ev\n";
1748 return (\$mis,\$miv,\$mfv,\$es,\$ev);
1749 }
1750 }
1751 return; # NaN, not a number
1752 }
1753
1754sub _digits
1755 {
1756 # computer number of digits in bigint, minus the sign
1757 # int() because add/sub leaves sometimes strings (like '00005') instead of
1758 # int ('5') in this place, causing length to fail
1759 my $cx = shift;
1760
1761 #print "len: ",(@$cx-1)*5+CORE::length(int($cx->[-1])),"\n";
1762 return (@$cx-1)*5+CORE::length(int($cx->[-1]));
1763 }
1764
1765sub as_number
1766 {
1767 # an object might be asked to return itself as bigint on certain overloaded
1768 # operations, this does exactly this, so that sub classes can simple inherit
1769 # it or override with their own integer conversion routine
1770 my $self = shift;
1771
1772 return $self->copy();
1773 }
1774
1775##############################################################################
1776# internal calculation routines
1777
1778sub acmp
1779 {
1780 # internal absolute post-normalized compare (ignore signs)
1781 # ref to array, ref to array, return <0, 0, >0
1782 # arrays must have at least on entry, this is not checked for
1783
1784 my ($cx, $cy) = @_;
1785
1786 #print "$cx $cy\n";
1787 my ($i,$a,$x,$y,$k);
1788 # calculate length based on digits, not parts
1789 $x = _digits($cx); $y = _digits($cy);
1790 # print "length: ",($x-$y),"\n";
1791 return $x-$y if ($x - $y); # if different in length
1792 #print "full compare\n";
1793 $i = 0; $a = 0;
1794 # first way takes 5.49 sec instead of 4.87, but has the early out advantage
1795 # so grep is slightly faster, but more unflexible. hm. $_ instead if $k
1796 # yields 5.6 instead of 5.5 sec huh?
1797 # manual way (abort if unequal, good for early ne)
1798 my $j = scalar @$cx - 1;
1799 while ($j >= 0)
1800 {
1801 # print "$cx->[$j] $cy->[$j] $a",$cx->[$j]-$cy->[$j],"\n";
1802 last if ($a = $cx->[$j] - $cy->[$j]); $j--;
1803 }
1804 return $a;
1805 # while it early aborts, it is even slower than the manual variant
1806 #grep { return $a if ($a = $_ - $cy->[$i++]); } @$cx;
1807 # grep way, go trough all (bad for early ne)
1808 #grep { $a = $_ - $cy->[$i++]; } @$cx;
1809 #return $a;
1810 }
1811
1812sub cmp
1813 {
1814 # post-normalized compare for internal use (honors signs)
1815 # ref to array, ref to array, return < 0, 0, >0
1816 my ($cx,$cy,$sx,$sy) = @_;
1817
1818 #return 0 if (is0($cx,$sx) && is0($cy,$sy));
1819
1820 if ($sx eq '+')
1821 {
1822 return 1 if $sy eq '-'; # 0 check handled above
1823 return acmp($cx,$cy);
1824 }
1825 else
1826 {
1827 # $sx eq '-'
1828 return -1 if ($sy eq '+');
1829 return acmp($cy,$cx);
1830 }
1831 return 0; # equal
1832 }
1833
1834sub add
1835 {
1836 # (ref to int_num_array, ref to int_num_array)
1837 # routine to add two base 1e5 numbers
1838 # stolen from Knuth Vol 2 Algorithm A pg 231
1839 # there are separate routines to add and sub as per Kunth pg 233
1840 # This routine clobbers up array x, but not y.
1841
1842 my ($x,$y) = @_;
1843
1844 # for each in Y, add Y to X and carry. If after that, something is left in
1845 # X, foreach in X add carry to X and then return X, carry
1846 # Trades one "$j++" for having to shift arrays, $j could be made integer
1847 # but this would impose a limit to number-length to 2**32.
1848 my $i; my $car = 0; my $j = 0;
1849 for $i (@$y)
1850 {
1851 $x->[$j] -= $BASE
1852 if $car = (($x->[$j] += $i + $car) >= $BASE) ? 1 : 0;
1853 $j++;
1854 }
1855 while ($car != 0)
1856 {
1857 $x->[$j] -= $BASE if $car = (($x->[$j] += $car) >= $BASE) ? 1 : 0; $j++;
1858 }
1859 }
1860
1861sub sub
1862 {
1863 # (ref to int_num_array, ref to int_num_array)
1864 # subtract base 1e5 numbers -- stolen from Knuth Vol 2 pg 232, $x > $y
1865 # subtract Y from X (X is always greater/equal!) by modifiyng x in place
1866 my ($sx,$sy,$s) = @_;
1867
1868 my $car = 0; my $i; my $j = 0;
1869 if (!$s)
1870 {
1871 #print "case 2\n";
1872 for $i (@$sx)
1873 {
1874 last unless defined $sy->[$j] || $car;
1875 #print "x: $i y: $sy->[$j] c: $car\n";
1876 $i += $BASE if $car = (($i -= ($sy->[$j] || 0) + $car) < 0); $j++;
1877 #print "x: $i y: $sy->[$j-1] c: $car\n";
1878 }
1879 # might leave leading zeros, so fix that
1880 _strip_zeros($sx);
1881 return $sx;
1882 }
1883 else
1884 {
1885 #print "case 1 (swap)\n";
1886 for $i (@$sx)
1887 {
1888 last unless defined $sy->[$j] || $car;
1889 #print "$sy->[$j] $i $car => $sx->[$j]\n";
1890 $sy->[$j] += $BASE
1891 if $car = (($sy->[$j] = $i-($sy->[$j]||0) - $car) < 0);
1892 #print "$sy->[$j] $i $car => $sy->[$j]\n";
1893 $j++;
1894 }
1895 # might leave leading zeros, so fix that
1896 _strip_zeros($sy);
1897 return $sy;
1898 }
1899 }
1900
1901sub mul
1902 {
1903 # (BINT, BINT) return nothing
1904 # multiply two numbers in internal representation
1905 # modifies first arg, second needs not be different from first
1906 my ($x,$y) = @_;
1907
1908 $x->{sign} = $x->{sign} ne $y->{sign} ? '-' : '+';
1909 my @prod = (); my ($prod,$car,$cty,$xi,$yi);
1910 my $xv = $x->{value};
1911 my $yv = $y->{value};
1912 # since multiplying $x with $x fails, make copy in this case
1913 $yv = [@$xv] if "$xv" eq "$yv";
1914 for $xi (@$xv)
1915 {
1916 $car = 0; $cty = 0;
1917 for $yi (@$yv)
1918 {
1919 $prod = $xi * $yi + ($prod[$cty] || 0) + $car;
1920 $prod[$cty++] =
1921 $prod - ($car = int($prod * 1e-5)) * $BASE; # see USE_MUL
1922 }
1923 $prod[$cty] += $car if $car; # need really to check for 0?
1924 $xi = shift @prod;
1925 }
1926 push @$xv, @prod;
1927 _strip_zeros($x->{value});
1928 # normalize (handled last to save check for $y->is_zero()
1929 $x->{sign} = '+' if @$xv == 1 && $xv->[0] == 0; # not is_zero due to '-'
1930 }
1931
1932sub div
1933 {
1934 # ref to array, ref to array, modify first array and return reminder if
1935 # in list context
1936 # does no longer handle sign
1937 my ($x,$yorg) = @_;
1938 my ($car,$bar,$prd,$dd,$xi,$yi,@q,$v2,$v1);
1939
1940 my (@d,$tmp,$q,$u2,$u1,$u0);
1941
1942 $car = $bar = $prd = 0;
1943
1944 my $y = [ @$yorg ];
1945 if (($dd = int($BASE/($y->[-1]+1))) != 1)
1946 {
1947 for $xi (@$x)
1948 {
1949 $xi = $xi * $dd + $car;
1950 $xi -= ($car = int($xi * $RBASE)) * $BASE; # see USE_MUL
1951 }
1952 push(@$x, $car); $car = 0;
1953 for $yi (@$y)
1954 {
1955 $yi = $yi * $dd + $car;
1956 $yi -= ($car = int($yi * $RBASE)) * $BASE; # see USE_MUL
1957 }
1958 }
1959 else
1960 {
1961 push(@$x, 0);
1962 }
1963 @q = (); ($v2,$v1) = @$y[-2,-1];
1964 $v2 = 0 unless $v2;
1965 while ($#$x > $#$y)
1966 {
1967 ($u2,$u1,$u0) = @$x[-3..-1];
1968 $u2 = 0 unless $u2;
1969 print "oups v1 is 0, u0: $u0 $y->[-2] $y->[-1] l ",scalar @$y,"\n"
1970 if $v1 == 0;
1971 $q = (($u0 == $v1) ? 99999 : int(($u0*$BASE+$u1)/$v1));
1972 --$q while ($v2*$q > ($u0*1e5+$u1-$q*$v1)*$BASE+$u2);
1973 if ($q)
1974 {
1975 ($car, $bar) = (0,0);
1976 for ($yi = 0, $xi = $#$x-$#$y-1; $yi <= $#$y; ++$yi,++$xi)
1977 {
1978 $prd = $q * $y->[$yi] + $car;
1979 $prd -= ($car = int($prd * $RBASE)) * $BASE; # see USE_MUL
1980 $x->[$xi] += 1e5 if ($bar = (($x->[$xi] -= $prd + $bar) < 0));
e16b8f49 1981 }
58cde26e
JH
1982 if ($x->[-1] < $car + $bar)
1983 {
1984 $car = 0; --$q;
1985 for ($yi = 0, $xi = $#$x-$#$y-1; $yi <= $#$y; ++$yi,++$xi)
1986 {
1987 $x->[$xi] -= 1e5
1988 if ($car = (($x->[$xi] += $y->[$yi] + $car) > $BASE));
1989 }
1990 }
1991 }
1992 pop(@$x); unshift(@q, $q);
e16b8f49 1993 }
58cde26e
JH
1994 if (wantarray)
1995 {
1996 @d = ();
1997 if ($dd != 1)
1998 {
1999 $car = 0;
2000 for $xi (reverse @$x)
2001 {
2002 $prd = $car * $BASE + $xi;
2003 $car = $prd - ($tmp = int($prd / $dd)) * $dd; # see USE_MUL
2004 unshift(@d, $tmp);
2005 }
2006 }
2007 else
2008 {
2009 @d = @$x;
2010 }
2011 @$x = @q;
2012 _strip_zeros($x);
2013 _strip_zeros(\@d);
2014 return ($x,\@d);
2015 }
2016 @$x = @q;
2017 _strip_zeros($x);
2018 return $x;
2019 }
2020
2021sub _lcm
2022 {
2023 # (BINT or num_str, BINT or num_str) return BINT
2024 # does modify first argument
2025 # LCM
2026
2027 my $x = shift; my $ty = shift;
2028 return $x->bnan() if ($x->{sign} eq $nan) || ($ty->{sign} eq $nan);
2029 return $x * $ty / bgcd($x,$ty);
2030 }
2031
2032sub _gcd_old
2033 {
2034 # (BINT or num_str, BINT or num_str) return BINT
2035 # does modify first arg
2036 # GCD -- Euclids algorithm E, Knuth Vol 2 pg 296
2037 trace(@_);
2038
2039 my $x = shift; my $ty = $class->new(shift); # preserve y
2040 return $x->bnan() if ($x->{sign} eq $nan) || ($ty->{sign} eq $nan);
2041
2042 while (!$ty->is_zero())
2043 {
2044 ($x, $ty) = ($ty,bmod($x,$ty));
2045 }
2046 $x;
2047 }
2048
2049sub _gcd
2050 {
2051 # (BINT or num_str, BINT or num_str) return BINT
2052 # does not modify arguments
2053 # GCD -- Euclids algorithm, variant L (Lehmer), Knuth Vol 3 pg 347 ff
2054 # unfortunately, it is slower and also seems buggy (the A=0, B=1, C=1, D=0
2055 # case..)
2056 trace(@_);
2057
2058 my $u = $class->new(shift); my $v = $class->new(shift); # preserve u,v
2059 return $u->bnan() if ($u->{sign} eq $nan) || ($v->{sign} eq $nan);
2060
2061 $u->babs(); $v->babs(); # Euclid is valid for |u| and |v|
2062
2063 my ($U,$V,$A,$B,$C,$D,$T,$Q); # single precision variables
2064 my ($t); # multiprecision variables
2065
2066 while ($v > $BASE)
2067 {
2068 #sleep 1;
2069 ($u,$v) = ($v,$u) if ($u < $v); # make sure that u >= v
2070 #print "gcd: $u $v\n";
2071 # step L1, initialize
2072 $A = 1; $B = 0; $C = 0; $D = 1;
2073 $U = $u->{value}->[-1]; # leading digits of u
2074 $V = $v->{value}->[-1]; # leading digits of v
2075
2076 # step L2, test quotient
2077 if (($V + $C != 0) && ($V + $D != 0)) # div by zero => go to L4
2078 {
2079 $Q = int(($U + $A)/($V + $C)); # quotient
2080 #print "L1 A=$A B=$B C=$C D=$D U=$U V=$V Q=$Q\n";
2081 # do L3? (div by zero => go to L4)
2082 while ($Q == int(($U + $B)/($V + $D)))
2083 {
2084 # step L3, emulate Euclid
2085 #print "L3a A=$A B=$B C=$C D=$D U=$U V=$V Q=$Q\n";
2086 $T = $A - $Q*$C; $A = $C; $C = $T;
2087 $T = $B - $Q*$D; $B = $D; $D = $T;
2088 $T = $U - $Q*$V; $U = $V; $V = $T;
2089 last if ($V + $D == 0) || ($V + $C == 0); # div by zero => L4
2090 $Q = int(($U + $A)/($V + $C)); # quotient for next test
2091 #print "L3b A=$A B=$B C=$C D=$D U=$U V=$V Q=$Q\n";
2092 }
2093 }
2094 # step L4, multiprecision step
2095 # was if ($B == 0)
2096 # in case A = 0, B = 1, C = 0 and D = 1, this case would simple swap u & v
2097 # and loop endless. Not sure why this happens, Knuth does not make a
2098 # remark about this special case. bug?
2099 if (($B == 0) || (($A == 0) && ($C == 1) && ($D == 0)))
2100 {
2101 #print "L4b1: u=$u v=$v\n";
2102 ($u,$v) = ($v,bmod($u,$v));
2103 #$t = $u % $v; $u = $v->copy(); $v = $t;
2104 #print "L4b12: u=$u v=$v\n";
2105 }
2106 else
2107 {
2108 #print "L4b: $u $v $A $B $C $D\n";
2109 $t = $A*$u + $B*$v; $v *= $D; $v += $C*$u; $u = $t;
2110 #print "L4b2: $u $v\n";
2111 }
2112 } # back to L1
e16b8f49 2113
58cde26e
JH
2114 return _gcd_old($u,$v) if $v != 1; # v too small
2115 return $v; # 1
2116 }
2117
2118###############################################################################
2119# this method return 0 if the object can be modified, or 1 for not
2120# We use a fast use constant statement here, to avoid costly calls. Subclasses
2121# may override it with special code (f.i. Math::BigInt::Constant does so)
2122
2123use constant modify => 0;
2124
2125#sub modify
2126# {
2127# my $self = shift;
2128# my $method = shift;
2129# print "original $self modify by $method\n";
2130# return 0; # $self;
2131# }
e16b8f49 2132
a0d0e21e 21331;
a5f75d66
AD
2134__END__
2135
2136=head1 NAME
2137
2138Math::BigInt - Arbitrary size integer math package
2139
2140=head1 SYNOPSIS
2141
2142 use Math::BigInt;
58cde26e
JH
2143
2144 # Number creation
2145 $x = Math::BigInt->new($str); # defaults to 0
2146 $nan = Math::BigInt->bnan(); # create a NotANumber
2147 $zero = Math::BigInt->bzero();# create a "+0"
2148
2149 # Testing
2150 $x->is_zero(); # return whether arg is zero or not
2151 $x->is_nan(); # return whether arg is NaN or not
2152 $x->is_one(); # return true if arg is +1
2153 $x->is_one('-'); # return true if arg is -1
2154 $x->is_odd(); # return true if odd, false for even
2155 $x->is_even(); # return true if even, false for odd
2156 $x->bcmp($y); # compare numbers (undef,<0,=0,>0)
2157 $x->bacmp($y); # compare absolutely (undef,<0,=0,>0)
2158 $x->sign(); # return the sign, either +,- or NaN
2159 $x->digit($n); # return the nth digit, counting from right
2160 $x->digit(-$n); # return the nth digit, counting from left
2161
2162 # The following all modify their first argument:
2163
2164 # set
2165 $x->bzero(); # set $x to 0
2166 $x->bnan(); # set $x to NaN
2167
2168 $x->bneg(); # negation
2169 $x->babs(); # absolute value
2170 $x->bnorm(); # normalize (no-op)
2171 $x->bnot(); # two's complement (bit wise not)
2172 $x->binc(); # increment x by 1
2173 $x->bdec(); # decrement x by 1
2174
2175 $x->badd($y); # addition (add $y to $x)
2176 $x->bsub($y); # subtraction (subtract $y from $x)
2177 $x->bmul($y); # multiplication (multiply $x by $y)
2178 $x->bdiv($y); # divide, set $x to quotient
2179 # return (quo,rem) or quo if scalar
2180
2181 $x->bmod($y); # modulus (x % y)
2182 $x->bpow($y); # power of arguments (x ** y)
2183 $x->blsft($y); # left shift
2184 $x->brsft($y); # right shift
2185 $x->blsft($y,$n); # left shift, by base $n (like 10)
2186 $x->brsft($y,$n); # right shift, by base $n (like 10)
2187
2188 $x->band($y); # bitwise and
2189 $x->bior($y); # bitwise inclusive or
2190 $x->bxor($y); # bitwise exclusive or
2191 $x->bnot(); # bitwise not (two's complement)
2192
2193 $x->bsqrt(); # calculate square-root
2194
2195 $x->round($A,$P,$round_mode); # round to accuracy or precision using mode $r
2196 $x->bround($N); # accuracy: preserve $N digits
2197 $x->bfround($N); # round to $Nth digit, no-op for BigInts
2198
2199 # The following do not modify their arguments in BigInt, but do in BigFloat:
2200 $x->bfloor(); # return integer less or equal than $x
2201 $x->bceil(); # return integer greater or equal than $x
2202
2203 # The following do not modify their arguments:
2204
2205 bgcd(@values); # greatest common divisor
2206 blcm(@values); # lowest common multiplicator
2207
2208 $x->bstr(); # normalized string
2209 $x->bsstr(); # normalized string in scientific notation
2210 $x->length(); # return number of digits in number
2211 ($x,$f) = $x->length(); # length of number and length of fraction part
2212
2213 $x->exponent(); # return exponent as BigInt
2214 $x->mantissa(); # return mantissa as BigInt
2215 $x->parts(); # return (mantissa,exponent) as BigInt
a5f75d66
AD
2216
2217=head1 DESCRIPTION
2218
58cde26e
JH
2219All operators (inlcuding basic math operations) are overloaded if you
2220declare your big integers as
a5f75d66 2221
58cde26e 2222 $i = new Math::BigInt '123_456_789_123_456_789';
a5f75d66 2223
58cde26e
JH
2224Operations with overloaded operators preserve the arguments which is
2225exactly what you expect.
a5f75d66
AD
2226
2227=over 2
2228
2229=item Canonical notation
2230
58cde26e 2231Big integer values are strings of the form C</^[+-]\d+$/> with leading
a5f75d66
AD
2232zeros suppressed.
2233
58cde26e
JH
2234 '-0' canonical value '-0', normalized '0'
2235 ' -123_123_123' canonical value '-123123123'
2236 '1_23_456_7890' canonical value '1234567890'
2237
a5f75d66
AD
2238=item Input
2239
58cde26e
JH
2240Input values to these routines may be either Math::BigInt objects or
2241strings of the form C</^\s*[+-]?[\d]+\.?[\d]*E?[+-]?[\d]*$/>.
2242
2243You can include one underscore between any two digits.
2244
2245This means integer values like 1.01E2 or even 1000E-2 are also accepted.
2246Non integer values result in NaN.
2247
2248Math::BigInt::new() defaults to 0, while Math::BigInt::new('') results
2249in 'NaN'.
2250
2251bnorm() on a BigInt object is now effectively a no-op, since the numbers
2252are always stored in normalized form. On a string, it creates a BigInt
2253object.
a5f75d66
AD
2254
2255=item Output
2256
58cde26e
JH
2257Output values are BigInt objects (normalized), except for bstr(), which
2258returns a string in normalized form.
2259Some routines (C<is_odd()>, C<is_even()>, C<is_zero()>, C<is_one()>,
2260C<is_nan()>) return true or false, while others (C<bcmp()>, C<bacmp()>)
2261return either undef, <0, 0 or >0 and are suited for sort.
a5f75d66
AD
2262
2263=back
2264
58cde26e 2265=head2 Rounding
a5f75d66 2266
58cde26e
JH
2267Only C<bround()> is defined for BigInts, for further details of rounding see
2268L<Math::BigFloat>.
2269
2270=over 2
a5f75d66 2271
5dc6f178 2272=item bfround ( +$scale )
a5f75d66 2273
5dc6f178 2274rounds to the $scale'th place left from the '.'
58cde26e 2275
5dc6f178
JH
2276=item bround ( +$scale )
2277
2278preserves accuracy to $scale significant digits counted from the left
2279and pads the number with zeros
2280
2281=item bround ( -$scale )
2282
2283preserves accuracy to $scale significant digits counted from the right
2284and pads the number with zeros.
58cde26e
JH
2285
2286=back
2287
2288C<bfround()> does nothing in case of negative C<$scale>. Both C<bround()> and
2289C<bfround()> are a no-ops for a scale of 0.
2290
2291All rounding functions take as a second parameter a rounding mode from one of
2292the following: 'even', 'odd', '+inf', '-inf', 'zero' or 'trunc'.
2293
2294The default is 'even'. By using C<< Math::BigInt->round_mode($rnd_mode); >>
2295you can get and set the default round mode for subsequent rounding.
2296
2297The second parameter to the round functions than overrides the default
2298temporarily.
2299
2300=head2 Internals
2301
2302Actual math is done in an internal format consisting of an array of
2303elements of base 100000 digits with the least significant digit first.
2304The sign C</^[+-]$/> is stored separately. The string 'NaN' is used to
2305represent the result when input arguments are not numbers, as well as
2306the result of dividing by zero.
2307
2308You sould neither care nor depend on the internal represantation, it might
2309change without notice. Use only method calls like C<< $x->sign(); >> instead
2310relying on the internal hash keys like in C<< $x->{sign}; >>.
2311
2312=head2 mantissa(), exponent() and parts()
2313
2314C<mantissa()> and C<exponent()> return the said parts of the BigInt such
2315that:
2316
2317 $m = $x->mantissa();
2318 $e = $x->exponent();
2319 $y = $m * ( 10 ** $e );
2320 print "ok\n" if $x == $y;
2321
2322C<($m,$e) = $x->parts()> is just a shortcut that gives you both of them in one
2323go. Both the returned mantissa and exponent do have a sign.
2324
2325Currently, for BigInts C<$e> will be always 0, except for NaN where it will be
2326NaN and for $x == 0, then it will be 1 (to be compatible with Math::BigFlaot's
2327internal representation of a zero as C<0E1>).
2328
2329C<$m> will always be a copy of the original number. The relation between $e
2330and $m might change in the future, but will be always equivalent in a
2331numerical sense, e.g. $m might get minized.
2332
2333=head1 EXAMPLES
2334
2335 use Math::BigInt qw(bstr bint);
2336 $x = bstr("1234") # string "1234"
2337 $x = "$x"; # same as bstr()
2338 $x = bneg("1234") # Bigint "-1234"
2339 $x = Math::BigInt->bneg("1234"); # Bigint "-1234"
2340 $x = Math::BigInt->babs("-12345"); # Bigint "12345"
2341 $x = Math::BigInt->bnorm("-0 00"); # BigInt "0"
2342 $x = bint(1) + bint(2); # BigInt "3"
2343 $x = bint(1) + "2"; # ditto (auto-BigIntify of "2")
2344 $x = bint(1); # BigInt "1"
2345 $x = $x + 5 / 2; # BigInt "3"
2346 $x = $x ** 3; # BigInt "27"
2347 $x *= 2; # BigInt "54"
2348 $x = new Math::BigInt; # BigInt "0"
2349 $x--; # BigInt "-1"
2350 $x = Math::BigInt->badd(4,5) # BigInt "9"
2351 $x = Math::BigInt::badd(4,5) # BigInt "9"
2352 print $x->bsstr(); # 9e+0
a5f75d66 2353
b3ac6de7
IZ
2354=head1 Autocreating constants
2355
58cde26e
JH
2356After C<use Math::BigInt ':constant'> all the B<integer> decimal constants
2357in the given scope are converted to C<Math::BigInt>. This conversion
b3ac6de7
IZ
2358happens at compile time.
2359
2360In particular
2361
58cde26e
JH
2362 perl -MMath::BigInt=:constant -e 'print 2**100,"\n"'
2363
2364prints the integer value of C<2**100>. Note that without conversion of
2365constants the expression 2**100 will be calculated as floating point
2366number.
2367
2368Please note that strings and floating point constants are not affected,
2369so that
2370
2371 use Math::BigInt qw/:constant/;
2372
2373 $x = 1234567890123456789012345678901234567890
2374 + 123456789123456789;
2375 $x = '1234567890123456789012345678901234567890'
2376 + '123456789123456789';
b3ac6de7 2377
58cde26e
JH
2378do both not work. You need a explicit Math::BigInt->new() around one of them.
2379
2380=head1 PERFORMANCE
2381
2382Using the form $x += $y; etc over $x = $x + $y is faster, since a copy of $x
2383must be made in the second case. For long numbers, the copy can eat up to 20%
2384of the work (in case of addition/subtraction, less for
2385multiplication/division). If $y is very small compared to $x, the form
2386$x += $y is MUCH faster than $x = $x + $y since making the copy of $x takes
2387more time then the actual addition.
2388
2389With a technic called copy-on-write the cost of copying with overload could
2390be minimized or even completely avoided. This is currently not implemented.
2391
2392The new version of this module is slower on new(), bstr() and numify(). Some
2393operations may be slower for small numbers, but are significantly faster for
2394big numbers. Other operations are now constant (O(1), like bneg(), babs()
2395etc), instead of O(N) and thus nearly always take much less time.
2396
2397For more benchmark results see http://bloodgate.com/perl/benchmarks.html
b3ac6de7 2398
a5f75d66
AD
2399=head1 BUGS
2400
58cde26e
JH
2401=over 2
2402
2403=item :constant and eval()
2404
2405Under Perl prior to 5.6.0 having an C<use Math::BigInt ':constant';> and
2406C<eval()> in your code will crash with "Out of memory". This is probably an
2407overload/exporter bug. You can workaround by not having C<eval()>
2408and ':constant' at the same time or upgrade your Perl.
2409
2410=back
2411
2412=head1 CAVEATS
2413
2414Some things might not work as you expect them. Below is documented what is
2415known to be troublesome:
2416
2417=over 1
2418
2419=item stringify, bstr(), bsstr() and 'cmp'
2420
2421Both stringify and bstr() now drop the leading '+'. The old code would return
2422'+3', the new returns '3'. This is to be consistent with Perl and to make
2423cmp (especially with overloading) to work as you expect. It also solves
2424problems with Test.pm, it's ok() uses 'eq' internally.
2425
2426Mark said, when asked about to drop the '+' altogether, or make only cmp work:
2427
2428 I agree (with the first alternative), don't add the '+' on positive
2429 numbers. It's not as important anymore with the new internal
2430 form for numbers. It made doing things like abs and neg easier,
2431 but those have to be done differently now anyway.
2432
2433So, the following examples will now work all as expected:
2434
2435 use Test;
2436 BEGIN { plan tests => 1 }
2437 use Math::BigInt;
2438
2439 my $x = new Math::BigInt 3*3;
2440 my $y = new Math::BigInt 3*3;
2441
2442 ok ($x,3*3);
2443 print "$x eq 9" if $x eq $y;
2444 print "$x eq 9" if $x eq '9';
2445 print "$x eq 9" if $x eq 3*3;
2446
2447Additionally, the following still works:
2448
2449 print "$x == 9" if $x == $y;
2450 print "$x == 9" if $x == 9;
2451 print "$x == 9" if $x == 3*3;
2452
2453There is now a C<bsstr()> method to get the string in scientific notation aka
2454C<1e+2> instead of C<100>. Be advised that overloaded 'eq' always uses bstr()
2455for comparisation, but Perl will represent some numbers as 100 and others
2456as 1e+308. If in doubt, convert both arguments to Math::BigInt before doing eq:
2457
2458 use Test;
2459 BEGIN { plan tests => 3 }
2460 use Math::BigInt;
2461
2462 $x = Math::BigInt->new('1e56'); $y = 1e56;
2463 ok ($x,$y); # will fail
2464 ok ($x->bsstr(),$y); # okay
2465 $y = Math::BigInt->new($y);
2466 ok ($x,$y); # okay
2467
2468=item int()
2469
2470C<int()> will return (at least for Perl v5.7.1 and up) another BigInt, not a
2471Perl scalar:
2472
2473 $x = Math::BigInt->new(123);
2474 $y = int($x); # BigInt 123
2475 $x = Math::BigFloat->new(123.45);
2476 $y = int($x); # BigInt 123
2477
2478In all Perl versions you can use C<as_number()> for the same effect:
2479
2480 $x = Math::BigFloat->new(123.45);
2481 $y = $x->as_number(); # BigInt 123
2482
2483This also works for other subclasses, like Math::String.
2484
2485=item bdiv
2486
2487The following will probably not do what you expect:
2488
2489 print $c->bdiv(10000),"\n";
2490
2491It prints both quotient and reminder since print calls C<bdiv()> in list
2492context. Also, C<bdiv()> will modify $c, so be carefull. You probably want
2493to use
2494
2495 print $c / 10000,"\n";
2496 print scalar $c->bdiv(10000),"\n"; # or if you want to modify $c
2497
2498instead.
2499
2500The quotient is always the greatest integer less than or equal to the
2501real-valued quotient of the two operands, and the remainder (when it is
2502nonzero) always has the same sign as the second operand; so, for
2503example,
2504
2505 1 / 4 => ( 0, 1)
2506 1 / -4 => (-1,-3)
2507 -3 / 4 => (-1, 1)
2508 -3 / -4 => ( 0,-3)
2509
2510As a consequence, the behavior of the operator % agrees with the
2511behavior of Perl's built-in % operator (as documented in the perlop
2512manpage), and the equation
2513
2514 $x == ($x / $y) * $y + ($x % $y)
2515
2516holds true for any $x and $y, which justifies calling the two return
2517values of bdiv() the quotient and remainder.
2518
2519Perl's 'use integer;' changes the behaviour of % and / for scalars, but will
2520not change BigInt's way to do things. This is because under 'use integer' Perl
2521will do what the underlying C thinks is right and this is different for each
2522system. If you need BigInt's behaving exactly like Perl's 'use integer', bug
2523the author to implement it ;)
2524
2525=item Modifying and =
2526
2527Beware of:
2528
2529 $x = Math::BigFloat->new(5);
2530 $y = $x;
2531
2532It will not do what you think, e.g. making a copy of $x. Instead it just makes
2533a second reference to the B<same> object and stores it in $y. Thus anything
2534that modifies $x will modify $y, and vice versa.
2535
2536 $x->bmul(2);
2537 print "$x, $y\n"; # prints '10, 10'
2538
2539If you want a true copy of $x, use:
2540
2541 $y = $x->copy();
2542
2543See also the documentation in L<overload> regarding C<=>.
2544
2545=item bpow
2546
2547C<bpow()> (and the rounding functions) now modifies the first argument and
2548return it, unlike the old code which left it alone and only returned the
2549result. This is to be consistent with C<badd()> etc. The first three will
2550modify $x, the last one won't:
2551
2552 print bpow($x,$i),"\n"; # modify $x
2553 print $x->bpow($i),"\n"; # ditto
2554 print $x **= $i,"\n"; # the same
2555 print $x ** $i,"\n"; # leave $x alone
2556
2557The form C<$x **= $y> is faster than C<$x = $x ** $y;>, though.
2558
2559=item Overloading -$x
2560
2561The following:
2562
2563 $x = -$x;
2564
2565is slower than
2566
2567 $x->bneg();
2568
2569since overload calls C<sub($x,0,1);> instead of C<neg($x)>. The first variant
2570needs to preserve $x since it does not know that it later will get overwritten.
2571This makes a copy of $x and takes O(N). But $x->bneg() is O(1).
2572
2573With Copy-On-Write, this issue will be gone. Stay tuned...
2574
2575=item Mixing different object types
2576
2577In Perl you will get a floating point value if you do one of the following:
2578
2579 $float = 5.0 + 2;
2580 $float = 2 + 5.0;
2581 $float = 5 / 2;
2582
2583With overloaded math, only the first two variants will result in a BigFloat:
2584
2585 use Math::BigInt;
2586 use Math::BigFloat;
2587
2588 $mbf = Math::BigFloat->new(5);
2589 $mbi2 = Math::BigInteger->new(5);
2590 $mbi = Math::BigInteger->new(2);
2591
2592 # what actually gets called:
2593 $float = $mbf + $mbi; # $mbf->badd()
2594 $float = $mbf / $mbi; # $mbf->bdiv()
2595 $integer = $mbi + $mbf; # $mbi->badd()
2596 $integer = $mbi2 / $mbi; # $mbi2->bdiv()
2597 $integer = $mbi2 / $mbf; # $mbi2->bdiv()
2598
2599This is because math with overloaded operators follows the first (dominating)
2600operand, this one's operation is called and returns thus the result. So,
2601Math::BigInt::bdiv() will always return a Math::BigInt, regardless whether
2602the result should be a Math::BigFloat or the second operant is one.
2603
2604To get a Math::BigFloat you either need to call the operation manually,
2605make sure the operands are already of the proper type or casted to that type
2606via Math::BigFloat->new():
2607
2608 $float = Math::BigFloat->new($mbi2) / $mbi; # = 2.5
2609
2610Beware of simple "casting" the entire expression, this would only convert
2611the already computed result:
2612
2613 $float = Math::BigFloat->new($mbi2 / $mbi); # = 2.0 thus wrong!
2614
2615Beware of the order of more complicated expressions like:
2616
2617 $integer = ($mbi2 + $mbi) / $mbf; # int / float => int
2618 $integer = $mbi2 / Math::BigFloat->new($mbi); # ditto
2619
2620If in doubt, break the expression into simpler terms, or cast all operands
2621to the desired resulting type.
2622
2623Scalar values are a bit different, since:
2624
2625 $float = 2 + $mbf;
2626 $float = $mbf + 2;
2627
2628will both result in the proper type due to the way the overloaded math works.
2629
2630This section also applies to other overloaded math packages, like Math::String.
2631
2632=item bsqrt()
2633
2634C<bsqrt()> works only good if the result is an big integer, e.g. the square
2635root of 144 is 12, but from 12 the square root is 3, regardless of rounding
2636mode.
2637
2638If you want a better approximation of the square root, then use:
2639
2640 $x = Math::BigFloat->new(12);
2641 $Math::BigFloat::precision = 0;
2642 Math::BigFloat->round_mode('even');
2643 print $x->copy->bsqrt(),"\n"; # 4
2644
2645 $Math::BigFloat::precision = 2;
2646 print $x->bsqrt(),"\n"; # 3.46
2647 print $x->bsqrt(3),"\n"; # 3.464
2648
2649=back
2650
2651=head1 LICENSE
2652
2653This program is free software; you may redistribute it and/or modify it under
2654the same terms as Perl itself.
a5f75d66 2655
58cde26e 2656=head1 AUTHORS
a5f75d66 2657
58cde26e
JH
2658Original code by Mark Biggar, overloaded interface by Ilya Zakharevich.
2659Completely rewritten by Tels http://bloodgate.com in late 2000, 2001.
a5f75d66
AD
2660
2661=cut