1package Math::BigInt::Lib;
2
3use 5.006001;
4use strict;
5use warnings;
6
7our $VERSION = '1.999837';
8$VERSION =~ tr/_//d;
9
10use Carp;
11
12use overload
13
14  # overload key: with_assign
15
16  '+'    => sub {
17                my $class = ref $_[0];
18                my $x = $class -> _copy($_[0]);
19                my $y = ref($_[1]) ? $_[1] : $class -> _new($_[1]);
20                return $class -> _add($x, $y);
21            },
22
23  '-'    => sub {
24                my $class = ref $_[0];
25                my ($x, $y);
26                if ($_[2]) {            # if swapped
27                    $y = $_[0];
28                    $x = ref($_[1]) ? $_[1] : $class -> _new($_[1]);
29                } else {
30                    $x = $class -> _copy($_[0]);
31                    $y = ref($_[1]) ? $_[1] : $class -> _new($_[1]);
32                }
33                return $class -> _sub($x, $y);
34            },
35
36  '*'    => sub {
37                my $class = ref $_[0];
38                my $x = $class -> _copy($_[0]);
39                my $y = ref($_[1]) ? $_[1] : $class -> _new($_[1]);
40                return $class -> _mul($x, $y);
41            },
42
43  '/'    => sub {
44                my $class = ref $_[0];
45                my ($x, $y);
46                if ($_[2]) {            # if swapped
47                    $y = $_[0];
48                    $x = ref($_[1]) ? $_[1] : $class -> _new($_[1]);
49                } else {
50                    $x = $class -> _copy($_[0]);
51                    $y = ref($_[1]) ? $_[1] : $class -> _new($_[1]);
52                }
53                return $class -> _div($x, $y);
54            },
55
56  '%'    => sub {
57                my $class = ref $_[0];
58                my ($x, $y);
59                if ($_[2]) {            # if swapped
60                    $y = $_[0];
61                    $x = ref($_[1]) ? $_[1] : $class -> _new($_[1]);
62                } else {
63                    $x = $class -> _copy($_[0]);
64                    $y = ref($_[1]) ? $_[1] : $class -> _new($_[1]);
65                }
66                return $class -> _mod($x, $y);
67            },
68
69  '**'   => sub {
70                my $class = ref $_[0];
71                my ($x, $y);
72                if ($_[2]) {            # if swapped
73                    $y = $_[0];
74                    $x = ref($_[1]) ? $_[1] : $class -> _new($_[1]);
75                } else {
76                    $x = $class -> _copy($_[0]);
77                    $y = ref($_[1]) ? $_[1] : $class -> _new($_[1]);
78                }
79                return $class -> _pow($x, $y);
80            },
81
82  '<<'   => sub {
83                my $class = ref $_[0];
84                my ($x, $y);
85                if ($_[2]) {            # if swapped
86                    $y = $class -> _num($_[0]);
87                    $x = ref($_[1]) ? $_[1] : $class -> _new($_[1]);
88                } else {
89                    $x = $_[0];
90                    $y = ref($_[1]) ? $class -> _num($_[1]) : $_[1];
91                }
92                return $class -> _lsft($x, $y);
93            },
94
95  '>>'   => sub {
96                my $class = ref $_[0];
97                my ($x, $y);
98                if ($_[2]) {            # if swapped
99                    $y = $_[0];
100                    $x = ref($_[1]) ? $_[1] : $class -> _new($_[1]);
101                } else {
102                    $x = $class -> _copy($_[0]);
103                    $y = ref($_[1]) ? $_[1] : $class -> _new($_[1]);
104                }
105                return $class -> _rsft($x, $y);
106            },
107
108  # overload key: num_comparison
109
110  '<'    => sub {
111                my $class = ref $_[0];
112                my ($x, $y);
113                if ($_[2]) {            # if swapped
114                    $y = $_[0];
115                    $x = ref($_[1]) ? $_[1] : $class -> _new($_[1]);
116                } else {
117                    $x = $class -> _copy($_[0]);
118                    $y = ref($_[1]) ? $_[1] : $class -> _new($_[1]);
119                }
120                return $class -> _acmp($x, $y) < 0;
121            },
122
123  '<='   => sub {
124                my $class = ref $_[0];
125                my ($x, $y);
126                if ($_[2]) {            # if swapped
127                    $y = $_[0];
128                    $x = ref($_[1]) ? $_[1] : $class -> _new($_[1]);
129                } else {
130                    $x = $class -> _copy($_[0]);
131                    $y = ref($_[1]) ? $_[1] : $class -> _new($_[1]);
132                }
133                return $class -> _acmp($x, $y) <= 0;
134            },
135
136  '>'    => sub {
137                my $class = ref $_[0];
138                my ($x, $y);
139                if ($_[2]) {            # if swapped
140                    $y = $_[0];
141                    $x = ref($_[1]) ? $_[1] : $class -> _new($_[1]);
142                } else {
143                    $x = $class -> _copy($_[0]);
144                    $y = ref($_[1]) ? $_[1] : $class -> _new($_[1]);
145                }
146                return $class -> _acmp($x, $y) > 0;
147            },
148
149  '>='   => sub {
150                my $class = ref $_[0];
151                my ($x, $y);
152                if ($_[2]) {            # if swapped
153                    $y = $_[0];
154                    $x = ref($_[1]) ? $_[1] : $class -> _new($_[1]);
155                } else {
156                    $x = $class -> _copy($_[0]);
157                    $y = ref($_[1]) ? $_[1] : $class -> _new($_[1]);
158                }
159                return $class -> _acmp($x, $y) >= 0;
160          },
161
162  '=='   => sub {
163                my $class = ref $_[0];
164                my $x = $class -> _copy($_[0]);
165                my $y = ref($_[1]) ? $_[1] : $class -> _new($_[1]);
166                return $class -> _acmp($x, $y) == 0;
167            },
168
169  '!='   => sub {
170                my $class = ref $_[0];
171                my $x = $class -> _copy($_[0]);
172                my $y = ref($_[1]) ? $_[1] : $class -> _new($_[1]);
173                return $class -> _acmp($x, $y) != 0;
174            },
175
176  # overload key: 3way_comparison
177
178  '<=>'  => sub {
179                my $class = ref $_[0];
180                my ($x, $y);
181                if ($_[2]) {            # if swapped
182                    $y = $_[0];
183                    $x = ref($_[1]) ? $_[1] : $class -> _new($_[1]);
184                } else {
185                    $x = $class -> _copy($_[0]);
186                    $y = ref($_[1]) ? $_[1] : $class -> _new($_[1]);
187                }
188                return $class -> _acmp($x, $y);
189            },
190
191  # overload key: binary
192
193  '&'    => sub {
194                my $class = ref $_[0];
195                my ($x, $y);
196                if ($_[2]) {            # if swapped
197                    $y = $_[0];
198                    $x = ref($_[1]) ? $_[1] : $class -> _new($_[1]);
199                } else {
200                    $x = $class -> _copy($_[0]);
201                    $y = ref($_[1]) ? $_[1] : $class -> _new($_[1]);
202                }
203                return $class -> _and($x, $y);
204            },
205
206  '|'    => sub {
207                my $class = ref $_[0];
208                my ($x, $y);
209                if ($_[2]) {            # if swapped
210                    $y = $_[0];
211                    $x = ref($_[1]) ? $_[1] : $class -> _new($_[1]);
212                } else {
213                    $x = $class -> _copy($_[0]);
214                    $y = ref($_[1]) ? $_[1] : $class -> _new($_[1]);
215                }
216                return $class -> _or($x, $y);
217            },
218
219  '^'    => sub {
220                my $class = ref $_[0];
221                my ($x, $y);
222                if ($_[2]) {            # if swapped
223                    $y = $_[0];
224                    $x = ref($_[1]) ? $_[1] : $class -> _new($_[1]);
225                } else {
226                    $x = $class -> _copy($_[0]);
227                    $y = ref($_[1]) ? $_[1] : $class -> _new($_[1]);
228                }
229                return $class -> _xor($x, $y);
230            },
231
232  # overload key: func
233
234  'abs'  => sub { $_[0] },
235
236  'sqrt' => sub {
237                my $class = ref $_[0];
238                return $class -> _sqrt($class -> _copy($_[0]));
239            },
240
241  'int'  => sub { $_[0] },
242
243  # overload key: conversion
244
245  'bool' => sub { ref($_[0]) -> _is_zero($_[0]) ? '' : 1; },
246
247  '""'   => sub { ref($_[0]) -> _str($_[0]); },
248
249  '0+'   => sub { ref($_[0]) -> _num($_[0]); },
250
251  '='    => sub { ref($_[0]) -> _copy($_[0]); },
252
253  ;
254
255sub _new {
256    croak "@{[(caller 0)[3]]} method not implemented";
257}
258
259sub _zero {
260    my $class = shift;
261    return $class -> _new("0");
262}
263
264sub _one {
265    my $class = shift;
266    return $class -> _new("1");
267}
268
269sub _two {
270    my $class = shift;
271    return $class -> _new("2");
272
273}
274sub _ten {
275    my $class = shift;
276    return $class -> _new("10");
277}
278
279sub _1ex {
280    my ($class, $exp) = @_;
281    $exp = $class -> _num($exp) if ref($exp);
282    return $class -> _new("1" . ("0" x $exp));
283}
284
285sub _copy {
286    my ($class, $x) = @_;
287    return $class -> _new($class -> _str($x));
288}
289
290# catch and throw away
291sub import { }
292
293##############################################################################
294# convert back to string and number
295
296sub _str {
297    # Convert number from internal base 1eN format to string format. Internal
298    # format is always normalized, i.e., no leading zeros.
299    croak "@{[(caller 0)[3]]} method not implemented";
300}
301
302sub _num {
303    my ($class, $x) = @_;
304    0 + $class -> _str($x);
305}
306
307##############################################################################
308# actual math code
309
310sub _add {
311    croak "@{[(caller 0)[3]]} method not implemented";
312}
313
314sub _sub {
315    croak "@{[(caller 0)[3]]} method not implemented";
316}
317
318sub _mul {
319    my ($class, $x, $y) = @_;
320    my $sum = $class -> _zero();
321    my $i   = $class -> _zero();
322    while ($class -> _acmp($i, $y) < 0) {
323        $sum = $class -> _add($sum, $x);
324        $i   = $class -> _inc($i);
325    }
326    return $sum;
327}
328
329sub _div {
330    my ($class, $x, $y) = @_;
331
332    croak "@{[(caller 0)[3]]} requires non-zero divisor"
333      if $class -> _is_zero($y);
334
335    my $r = $class -> _copy($x);
336    my $q = $class -> _zero();
337    while ($class -> _acmp($r, $y) >= 0) {
338        $q = $class -> _inc($q);
339        $r = $class -> _sub($r, $y);
340    }
341
342    return $q, $r if wantarray;
343    return $q;
344}
345
346sub _inc {
347    my ($class, $x) = @_;
348    $class -> _add($x, $class -> _one());
349}
350
351sub _dec {
352    my ($class, $x) = @_;
353    $class -> _sub($x, $class -> _one());
354}
355
356# Signed addition. If the flag is false, $xa might be modified, but not $ya. If
357# the false is true, $ya might be modified, but not $xa.
358
359sub _sadd {
360    my $class = shift;
361    my ($xa, $xs, $ya, $ys, $flag) = @_;
362    my ($za, $zs);
363
364    # If the signs are equal we can add them (-5 + -3 => -(5 + 3) => -8)
365
366    if ($xs eq $ys) {
367        if ($flag) {
368            $za = $class -> _add($ya, $xa);
369        } else {
370            $za = $class -> _add($xa, $ya);
371        }
372        $zs = $class -> _is_zero($za) ? '+' : $xs;
373        return $za, $zs;
374    }
375
376    my $acmp = $class -> _acmp($xa, $ya);       # abs(x) = abs(y)
377
378    if ($acmp == 0) {                           # x = -y or -x = y
379        $za = $class -> _zero();
380        $zs = '+';
381        return $za, $zs;
382    }
383
384    if ($acmp > 0) {                            # abs(x) > abs(y)
385        $za = $class -> _sub($xa, $ya, $flag);
386        $zs = $xs;
387    } else {                                    # abs(x) < abs(y)
388        $za = $class -> _sub($ya, $xa, !$flag);
389        $zs = $ys;
390    }
391    return $za, $zs;
392}
393
394# Signed subtraction. If the flag is false, $xa might be modified, but not $ya.
395# If the false is true, $ya might be modified, but not $xa.
396
397sub _ssub {
398    my $class = shift;
399    my ($xa, $xs, $ya, $ys, $flag) = @_;
400
401    # Swap sign of second operand and let _sadd() do the job.
402    $ys = $ys eq '+' ? '-' : '+';
403    $class -> _sadd($xa, $xs, $ya, $ys, $flag);
404}
405
406##############################################################################
407# testing
408
409sub _acmp {
410    # Compare two (absolute) values. Return -1, 0, or 1.
411    my ($class, $x, $y) = @_;
412    my $xstr = $class -> _str($x);
413    my $ystr = $class -> _str($y);
414
415    length($xstr) <=> length($ystr) || $xstr cmp $ystr;
416}
417
418sub _len {
419    my ($class, $x) = @_;
420    CORE::length($class -> _str($x));
421}
422
423sub _alen {
424    my ($class, $x) = @_;
425    $class -> _len($x);
426}
427
428sub _digit {
429    my ($class, $x, $n) = @_;
430    substr($class ->_str($x), -($n+1), 1);
431}
432
433sub _digitsum {
434    my ($class, $x) = @_;
435
436    my $len = $class -> _len($x);
437    my $sum = $class -> _zero();
438    for (my $i = 0 ; $i < $len ; ++$i) {
439        my $digit = $class -> _digit($x, $i);
440        $digit = $class -> _new($digit);
441        $sum = $class -> _add($sum, $digit);
442    }
443
444    return $sum;
445}
446
447sub _zeros {
448    my ($class, $x) = @_;
449    my $str = $class -> _str($x);
450    $str =~ /[^0](0*)\z/ ? CORE::length($1) : 0;
451}
452
453##############################################################################
454# _is_* routines
455
456sub _is_zero {
457    # return true if arg is zero
458    my ($class, $x) = @_;
459    $class -> _str($x) == 0;
460}
461
462sub _is_even {
463    # return true if arg is even
464    my ($class, $x) = @_;
465    substr($class -> _str($x), -1, 1) % 2 == 0;
466}
467
468sub _is_odd {
469    # return true if arg is odd
470    my ($class, $x) = @_;
471    substr($class -> _str($x), -1, 1) % 2 != 0;
472}
473
474sub _is_one {
475    # return true if arg is one
476    my ($class, $x) = @_;
477    $class -> _str($x) == 1;
478}
479
480sub _is_two {
481    # return true if arg is two
482    my ($class, $x) = @_;
483    $class -> _str($x) == 2;
484}
485
486sub _is_ten {
487    # return true if arg is ten
488    my ($class, $x) = @_;
489    $class -> _str($x) == 10;
490}
491
492###############################################################################
493# check routine to test internal state for corruptions
494
495sub _check {
496    # used by the test suite
497    my ($class, $x) = @_;
498    return "Input is undefined" unless defined $x;
499    return "$x is not a reference" unless ref($x);
500    return 0;
501}
502
503###############################################################################
504
505sub _mod {
506    # modulus
507    my ($class, $x, $y) = @_;
508
509    croak "@{[(caller 0)[3]]} requires non-zero second operand"
510      if $class -> _is_zero($y);
511
512    if ($class -> can('_div')) {
513        $x = $class -> _copy($x);
514        my ($q, $r) = $class -> _div($x, $y);
515        return $r;
516    } else {
517        my $r = $class -> _copy($x);
518        while ($class -> _acmp($r, $y) >= 0) {
519            $r = $class -> _sub($r, $y);
520        }
521        return $r;
522    }
523}
524
525##############################################################################
526# shifts
527
528sub _rsft {
529    my ($class, $x, $n, $b) = @_;
530    $b = $class -> _new($b) unless ref $b;
531    return scalar $class -> _div($x, $class -> _pow($class -> _copy($b), $n));
532}
533
534sub _lsft {
535    my ($class, $x, $n, $b) = @_;
536    $b = $class -> _new($b) unless ref $b;
537    return $class -> _mul($x, $class -> _pow($class -> _copy($b), $n));
538}
539
540sub _pow {
541    # power of $x to $y
542    my ($class, $x, $y) = @_;
543
544    if ($class -> _is_zero($y)) {
545        return $class -> _one();        # y == 0 => x => 1
546    }
547
548    if (($class -> _is_one($x)) ||      #    x == 1
549        ($class -> _is_one($y)))        # or y == 1
550    {
551        return $x;
552    }
553
554    if ($class -> _is_zero($x)) {
555        return $class -> _zero();       # 0 ** y => 0 (if not y <= 0)
556    }
557
558    my $pow2 = $class -> _one();
559
560    my $y_bin = $class -> _as_bin($y);
561    $y_bin =~ s/^0b//;
562    my $len = length($y_bin);
563
564    while (--$len > 0) {
565        $pow2 = $class -> _mul($pow2, $x) if substr($y_bin, $len, 1) eq '1';
566        $x = $class -> _mul($x, $x);
567    }
568
569    $x = $class -> _mul($x, $pow2);
570    return $x;
571}
572
573sub _nok {
574    # Return binomial coefficient (n over k).
575    my ($class, $n, $k) = @_;
576
577    # If k > n/2, or, equivalently, 2*k > n, compute nok(n, k) as
578    # nok(n, n-k), to minimize the number if iterations in the loop.
579
580    {
581        my $twok = $class -> _mul($class -> _two(), $class -> _copy($k));
582        if ($class -> _acmp($twok, $n) > 0) {
583            $k = $class -> _sub($class -> _copy($n), $k);
584        }
585    }
586
587    # Example:
588    #
589    # / 7 \       7!       1*2*3*4 * 5*6*7   5 * 6 * 7
590    # |   | = --------- =  --------------- = --------- = ((5 * 6) / 2 * 7) / 3
591    # \ 3 /   (7-3)! 3!    1*2*3*4 * 1*2*3   1 * 2 * 3
592    #
593    # Equivalently, _nok(11, 5) is computed as
594    #
595    # (((((((7 * 8) / 2) * 9) / 3) * 10) / 4) * 11) / 5
596
597    if ($class -> _is_zero($k)) {
598        return $class -> _one();
599    }
600
601    # Make a copy of the original n, in case the subclass modifies n in-place.
602
603    my $n_orig = $class -> _copy($n);
604
605    # n = 5, f = 6, d = 2 (cf. example above)
606
607    $n = $class -> _sub($n, $k);
608    $n = $class -> _inc($n);
609
610    my $f = $class -> _copy($n);
611    $f = $class -> _inc($f);
612
613    my $d = $class -> _two();
614
615    # while f <= n (the original n, that is) ...
616
617    while ($class -> _acmp($f, $n_orig) <= 0) {
618        $n = $class -> _mul($n, $f);
619        $n = $class -> _div($n, $d);
620        $f = $class -> _inc($f);
621        $d = $class -> _inc($d);
622    }
623
624    return $n;
625}
626
627#sub _fac {
628#    # factorial
629#    my ($class, $x) = @_;
630#
631#    my $two = $class -> _two();
632#
633#    if ($class -> _acmp($x, $two) < 0) {
634#        return $class -> _one();
635#    }
636#
637#    my $i = $class -> _copy($x);
638#    while ($class -> _acmp($i, $two) > 0) {
639#        $i = $class -> _dec($i);
640#        $x = $class -> _mul($x, $i);
641#    }
642#
643#    return $x;
644#}
645
646sub _fac {
647    # factorial
648    my ($class, $x) = @_;
649
650    # This is an implementation of the split recursive algorithm. See
651    # http://www.luschny.de/math/factorial/csharp/FactorialSplit.cs.html
652
653    my $p   = $class -> _one();
654    my $r   = $class -> _one();
655    my $two = $class -> _two();
656
657    my ($log2n) = $class -> _log_int($class -> _copy($x), $two);
658    my $h     = $class -> _zero();
659    my $shift = $class -> _zero();
660    my $k     = $class -> _one();
661
662    while ($class -> _acmp($h, $x)) {
663        $shift = $class -> _add($shift, $h);
664        $h = $class -> _rsft($class -> _copy($x), $log2n, $two);
665        $log2n = $class -> _dec($log2n) if !$class -> _is_zero($log2n);
666        my $high = $class -> _copy($h);
667        $high = $class -> _dec($high) if $class -> _is_even($h);
668        while ($class -> _acmp($k, $high)) {
669            $k = $class -> _add($k, $two);
670            $p = $class -> _mul($p, $k);
671        }
672        $r = $class -> _mul($r, $p);
673    }
674    return $class -> _lsft($r, $shift, $two);
675}
676
677sub _dfac {
678    # double factorial
679    my ($class, $x) = @_;
680
681    my $two = $class -> _two();
682
683    if ($class -> _acmp($x, $two) < 0) {
684        return $class -> _one();
685    }
686
687    my $i = $class -> _copy($x);
688    while ($class -> _acmp($i, $two) > 0) {
689        $i = $class -> _sub($i, $two);
690        $x = $class -> _mul($x, $i);
691    }
692
693    return $x;
694}
695
696sub _log_int {
697    # calculate integer log of $x to base $base
698    # calculate integer log of $x to base $base
699    # ref to array, ref to array - return ref to array
700    my ($class, $x, $base) = @_;
701
702    # X == 0 => NaN
703    return if $class -> _is_zero($x);
704
705    $base = $class -> _new(2)     unless defined($base);
706    $base = $class -> _new($base) unless ref($base);
707
708    # BASE 0 or 1 => NaN
709    return if $class -> _is_zero($base) || $class -> _is_one($base);
710
711    # X == 1 => 0 (is exact)
712    if ($class -> _is_one($x)) {
713        return $class -> _zero(), 1;
714    }
715
716    my $cmp = $class -> _acmp($x, $base);
717
718    # X == BASE => 1 (is exact)
719    if ($cmp == 0) {
720        return $class -> _one(), 1;
721    }
722
723    # 1 < X < BASE => 0 (is truncated)
724    if ($cmp < 0) {
725        return $class -> _zero(), 0;
726    }
727
728    my $y;
729
730    # log(x) / log(b) = log(xm * 10^xe) / log(bm * 10^be)
731    #                 = (log(xm) + xe*(log(10))) / (log(bm) + be*log(10))
732
733    {
734        my $x_str = $class -> _str($x);
735        my $b_str = $class -> _str($base);
736        my $xm    = "." . $x_str;
737        my $bm    = "." . $b_str;
738        my $xe    = length($x_str);
739        my $be    = length($b_str);
740        my $log10 = log(10);
741        my $guess = int((log($xm) + $xe * $log10) / (log($bm) + $be * $log10));
742        $y = $class -> _new($guess);
743    }
744
745    my $trial = $class -> _pow($class -> _copy($base), $y);
746    my $acmp  = $class -> _acmp($trial, $x);
747
748    # Did we get the exact result?
749
750    return $y, 1 if $acmp == 0;
751
752    # Too small?
753
754    while ($acmp < 0) {
755        $trial = $class -> _mul($trial, $base);
756        $y     = $class -> _inc($y);
757        $acmp  = $class -> _acmp($trial, $x);
758    }
759
760    # Too big?
761
762    while ($acmp > 0) {
763        $trial = $class -> _div($trial, $base);
764        $y     = $class -> _dec($y);
765        $acmp  = $class -> _acmp($trial, $x);
766    }
767
768    return $y, 1 if $acmp == 0;         # result is exact
769    return $y, 0;                       # result is too small
770}
771
772sub _sqrt {
773    # square-root of $y in place
774    my ($class, $y) = @_;
775
776    return $y if $class -> _is_zero($y);
777
778    my $y_str = $class -> _str($y);
779    my $y_len = length($y_str);
780
781    # Compute the guess $x.
782
783    my $xm;
784    my $xe;
785    if ($y_len % 2 == 0) {
786        $xm = sqrt("." . $y_str);
787        $xe = $y_len / 2;
788        $xm = sprintf "%.0f", int($xm * 1e15);
789        $xe -= 15;
790    } else {
791        $xm = sqrt(".0" . $y_str);
792        $xe = ($y_len + 1) / 2;
793        $xm = sprintf "%.0f", int($xm * 1e16);
794        $xe -= 16;
795    }
796
797    my $x;
798    if ($xe < 0) {
799        $x = substr $xm, 0, length($xm) + $xe;
800    } else {
801        $x = $xm . ("0" x $xe);
802    }
803
804    $x = $class -> _new($x);
805
806    # Newton's method for computing square root of y
807    #
808    # x(i+1) = x(i) - f(x(i)) / f'(x(i))
809    #        = x(i) - (x(i)^2 - y) / (2 * x(i))     # use if x(i)^2 > y
810    #        = x(i) + (y - x(i)^2) / (2 * x(i))     # use if x(i)^2 < y
811
812    # Determine if x, our guess, is too small, correct, or too large.
813
814    my $xsq = $class -> _mul($class -> _copy($x), $x);          # x(i)^2
815    my $acmp = $class -> _acmp($xsq, $y);                       # x(i)^2 <=> y
816
817    # Only assign a value to this variable if we will be using it.
818
819    my $two;
820    $two = $class -> _two() if $acmp != 0;
821
822    # If x is too small, do one iteration of Newton's method. Since the
823    # function f(x) = x^2 - y is concave and monotonically increasing, the next
824    # guess for x will either be correct or too large.
825
826    if ($acmp < 0) {
827
828        # x(i+1) = x(i) + (y - x(i)^2) / (2 * x(i))
829
830        my $numer = $class -> _sub($class -> _copy($y), $xsq);  # y - x(i)^2
831        my $denom = $class -> _mul($class -> _copy($two), $x);  # 2 * x(i)
832        my $delta = $class -> _div($numer, $denom);
833
834        unless ($class -> _is_zero($delta)) {
835            $x    = $class -> _add($x, $delta);
836            $xsq  = $class -> _mul($class -> _copy($x), $x);    # x(i)^2
837            $acmp = $class -> _acmp($xsq, $y);                  # x(i)^2 <=> y
838        }
839    }
840
841    # If our guess for x is too large, apply Newton's method repeatedly until
842    # we either have got the correct value, or the delta is zero.
843
844    while ($acmp > 0) {
845
846        # x(i+1) = x(i) - (x(i)^2 - y) / (2 * x(i))
847
848        my $numer = $class -> _sub($xsq, $y);                   # x(i)^2 - y
849        my $denom = $class -> _mul($class -> _copy($two), $x);  # 2 * x(i)
850        my $delta = $class -> _div($numer, $denom);
851        last if $class -> _is_zero($delta);
852
853        $x    = $class -> _sub($x, $delta);
854        $xsq  = $class -> _mul($class -> _copy($x), $x);        # x(i)^2
855        $acmp = $class -> _acmp($xsq, $y);                      # x(i)^2 <=> y
856    }
857
858    # When the delta is zero, our value for x might still be too large. We
859    # require that the outout is either exact or too small (i.e., rounded down
860    # to the nearest integer), so do a final check.
861
862    while ($acmp > 0) {
863        $x    = $class -> _dec($x);
864        $xsq  = $class -> _mul($class -> _copy($x), $x);        # x(i)^2
865        $acmp = $class -> _acmp($xsq, $y);                      # x(i)^2 <=> y
866    }
867
868    return $x;
869}
870
871sub _root {
872    my ($class, $y, $n) = @_;
873
874    return $y if $class -> _is_zero($y) || $class -> _is_one($y) ||
875                 $class -> _is_one($n);
876
877    # If y <= n, the result is always (truncated to) 1.
878
879    return $class -> _one() if $class -> _acmp($y, $n) <= 0;
880
881    # Compute the initial guess x of y^(1/n). When n is large, Newton's method
882    # converges slowly if the "guess" (initial value) is poor, so we need a
883    # good guess. It the guess is too small, the next guess will be too large,
884    # and from then on all guesses are too large.
885
886    my $DEBUG = 0;
887
888    # Split y into mantissa and exponent in base 10, so that
889    #
890    #   y = xm * 10^xe, where 0 < xm < 1 and xe is an integer
891
892    my $y_str  = $class -> _str($y);
893    my $ym = "." . $y_str;
894    my $ye = length($y_str);
895
896    # From this compute the approximate base 10 logarithm of y
897    #
898    #   log_10(y) = log_10(ym) + log_10(ye^10)
899    #             = log(ym)/log(10) + ye
900
901    my $log10y = log($ym) / log(10) + $ye;
902
903    # And from this compute the approximate base 10 logarithm of x, where
904    # x = y^(1/n)
905    #
906    #   log_10(x) = log_10(y)/n
907
908    my $log10x = $log10y / $class -> _num($n);
909
910    # From this compute xm and xe, the mantissa and exponent (in base 10) of x,
911    # where 1 < xm <= 10 and xe is an integer.
912
913    my $xe = int $log10x;
914    my $xm = 10 ** ($log10x - $xe);
915
916    # Scale the mantissa and exponent to increase the integer part of ym, which
917    # gives us better accuracy.
918
919    if ($DEBUG) {
920        print "\n";
921        print "y_str  = $y_str\n";
922        print "ym     = $ym\n";
923        print "ye     = $ye\n";
924        print "log10y = $log10y\n";
925        print "log10x = $log10x\n";
926        print "xm     = $xm\n";
927        print "xe     = $xe\n";
928    }
929
930    my $d = $xe < 15 ? $xe : 15;
931    $xm *= 10 ** $d;
932    $xe -= $d;
933
934    if ($DEBUG) {
935        print "\n";
936        print "xm     = $xm\n";
937        print "xe     = $xe\n";
938    }
939
940    # If the mantissa is not an integer, round up to nearest integer, and then
941    # convert the number to a string. It is important to always round up due to
942    # how Newton's method behaves in this case. If the initial guess is too
943    # small, the next guess will be too large, after which every succeeding
944    # guess converges the correct value from above. Now, if the initial guess
945    # is too small and n is large, the next guess will be much too large and
946    # require a large number of iterations to get close to the solution.
947    # Because of this, we are likely to find the solution faster if we make
948    # sure the initial guess is not too small.
949
950    my $xm_int = int($xm);
951    my $x_str = sprintf '%.0f', $xm > $xm_int ? $xm_int + 1 : $xm_int;
952    $x_str .= "0" x $xe;
953
954    my $x = $class -> _new($x_str);
955
956    if ($DEBUG) {
957        print "xm     = $xm\n";
958        print "xe     = $xe\n";
959        print "\n";
960        print "x_str  = $x_str (initial guess)\n";
961        print "\n";
962    }
963
964    # Use Newton's method for computing n'th root of y.
965    #
966    # x(i+1) = x(i) - f(x(i)) / f'(x(i))
967    #        = x(i) - (x(i)^n - y) / (n * x(i)^(n-1))   # use if x(i)^n > y
968    #        = x(i) + (y - x(i)^n) / (n * x(i)^(n-1))   # use if x(i)^n < y
969
970    # Determine if x, our guess, is too small, correct, or too large. Rather
971    # than computing x(i)^n and x(i)^(n-1) directly, compute x(i)^(n-1) and
972    # then the same value multiplied by x.
973
974    my $nm1     = $class -> _dec($class -> _copy($n));           # n-1
975    my $xpownm1 = $class -> _pow($class -> _copy($x), $nm1);     # x(i)^(n-1)
976    my $xpown   = $class -> _mul($class -> _copy($xpownm1), $x); # x(i)^n
977    my $acmp    = $class -> _acmp($xpown, $y);                   # x(i)^n <=> y
978
979    if ($DEBUG) {
980        print "\n";
981        print "x      = ", $class -> _str($x), "\n";
982        print "x^n    = ", $class -> _str($xpown), "\n";
983        print "y      = ", $class -> _str($y), "\n";
984        print "acmp   = $acmp\n";
985    }
986
987    # If x is too small, do one iteration of Newton's method. Since the
988    # function f(x) = x^n - y is concave and monotonically increasing, the next
989    # guess for x will either be correct or too large.
990
991    if ($acmp < 0) {
992
993        # x(i+1) = x(i) + (y - x(i)^n) / (n * x(i)^(n-1))
994
995        my $numer = $class -> _sub($class -> _copy($y), $xpown);    # y - x(i)^n
996        my $denom = $class -> _mul($class -> _copy($n), $xpownm1);  # n * x(i)^(n-1)
997        my $delta = $class -> _div($numer, $denom);
998
999        if ($DEBUG) {
1000            print "\n";
1001            print "numer  = ", $class -> _str($numer), "\n";
1002            print "denom  = ", $class -> _str($denom), "\n";
1003            print "delta  = ", $class -> _str($delta), "\n";
1004        }
1005
1006        unless ($class -> _is_zero($delta)) {
1007            $x       = $class -> _add($x, $delta);
1008            $xpownm1 = $class -> _pow($class -> _copy($x), $nm1);     # x(i)^(n-1)
1009            $xpown   = $class -> _mul($class -> _copy($xpownm1), $x); # x(i)^n
1010            $acmp    = $class -> _acmp($xpown, $y);                   # x(i)^n <=> y
1011
1012            if ($DEBUG) {
1013                print "\n";
1014                print "x      = ", $class -> _str($x), "\n";
1015                print "x^n    = ", $class -> _str($xpown), "\n";
1016                print "y      = ", $class -> _str($y), "\n";
1017                print "acmp   = $acmp\n";
1018            }
1019        }
1020    }
1021
1022    # If our guess for x is too large, apply Newton's method repeatedly until
1023    # we either have got the correct value, or the delta is zero.
1024
1025    while ($acmp > 0) {
1026
1027        # x(i+1) = x(i) - (x(i)^n - y) / (n * x(i)^(n-1))
1028
1029        my $numer = $class -> _sub($class -> _copy($xpown), $y);    # x(i)^n - y
1030        my $denom = $class -> _mul($class -> _copy($n), $xpownm1);  # n * x(i)^(n-1)
1031
1032        if ($DEBUG) {
1033            print "numer  = ", $class -> _str($numer), "\n";
1034            print "denom  = ", $class -> _str($denom), "\n";
1035        }
1036
1037        my $delta = $class -> _div($numer, $denom);
1038
1039        if ($DEBUG) {
1040            print "delta  = ", $class -> _str($delta), "\n";
1041        }
1042
1043        last if $class -> _is_zero($delta);
1044
1045        $x       = $class -> _sub($x, $delta);
1046        $xpownm1 = $class -> _pow($class -> _copy($x), $nm1);     # x(i)^(n-1)
1047        $xpown   = $class -> _mul($class -> _copy($xpownm1), $x); # x(i)^n
1048        $acmp    = $class -> _acmp($xpown, $y);                   # x(i)^n <=> y
1049
1050        if ($DEBUG) {
1051            print "\n";
1052            print "x      = ", $class -> _str($x), "\n";
1053            print "x^n    = ", $class -> _str($xpown), "\n";
1054            print "y      = ", $class -> _str($y), "\n";
1055            print "acmp   = $acmp\n";
1056        }
1057    }
1058
1059    # When the delta is zero, our value for x might still be too large. We
1060    # require that the outout is either exact or too small (i.e., rounded down
1061    # to the nearest integer), so do a final check.
1062
1063    while ($acmp > 0) {
1064        $x     = $class -> _dec($x);
1065        $xpown = $class -> _pow($class -> _copy($x), $n);     # x(i)^n
1066        $acmp  = $class -> _acmp($xpown, $y);                 # x(i)^n <=> y
1067    }
1068
1069    return $x;
1070}
1071
1072##############################################################################
1073# binary stuff
1074
1075sub _and {
1076    my ($class, $x, $y) = @_;
1077
1078    return $x if $class -> _acmp($x, $y) == 0;
1079
1080    my $m    = $class -> _one();
1081    my $mask = $class -> _new("32768");
1082
1083    my ($xr, $yr);                # remainders after division
1084
1085    my $xc = $class -> _copy($x);
1086    my $yc = $class -> _copy($y);
1087    my $z  = $class -> _zero();
1088
1089    until ($class -> _is_zero($xc) || $class -> _is_zero($yc)) {
1090        ($xc, $xr) = $class -> _div($xc, $mask);
1091        ($yc, $yr) = $class -> _div($yc, $mask);
1092        my $bits = $class -> _new($class -> _num($xr) & $class -> _num($yr));
1093        $z = $class -> _add($z, $class -> _mul($bits, $m));
1094        $m = $class -> _mul($m, $mask);
1095    }
1096
1097    return $z;
1098}
1099
1100sub _xor {
1101    my ($class, $x, $y) = @_;
1102
1103    return $class -> _zero() if $class -> _acmp($x, $y) == 0;
1104
1105    my $m    = $class -> _one();
1106    my $mask = $class -> _new("32768");
1107
1108    my ($xr, $yr);                # remainders after division
1109
1110    my $xc = $class -> _copy($x);
1111    my $yc = $class -> _copy($y);
1112    my $z  = $class -> _zero();
1113
1114    until ($class -> _is_zero($xc) || $class -> _is_zero($yc)) {
1115        ($xc, $xr) = $class -> _div($xc, $mask);
1116        ($yc, $yr) = $class -> _div($yc, $mask);
1117        my $bits = $class -> _new($class -> _num($xr) ^ $class -> _num($yr));
1118        $z = $class -> _add($z, $class -> _mul($bits, $m));
1119        $m = $class -> _mul($m, $mask);
1120    }
1121
1122    # The loop above stops when the smallest of the two numbers is exhausted.
1123    # The remainder of the longer one will survive bit-by-bit, so we simple
1124    # multiply-add it in.
1125
1126    $z = $class -> _add($z, $class -> _mul($xc, $m))
1127      unless $class -> _is_zero($xc);
1128    $z = $class -> _add($z, $class -> _mul($yc, $m))
1129      unless $class -> _is_zero($yc);
1130
1131    return $z;
1132}
1133
1134sub _or {
1135    my ($class, $x, $y) = @_;
1136
1137    return $x if $class -> _acmp($x, $y) == 0; # shortcut (see _and)
1138
1139    my $m    = $class -> _one();
1140    my $mask = $class -> _new("32768");
1141
1142    my ($xr, $yr);                # remainders after division
1143
1144    my $xc = $class -> _copy($x);
1145    my $yc = $class -> _copy($y);
1146    my $z  = $class -> _zero();
1147
1148    until ($class -> _is_zero($xc) || $class -> _is_zero($yc)) {
1149        ($xc, $xr) = $class -> _div($xc, $mask);
1150        ($yc, $yr) = $class -> _div($yc, $mask);
1151        my $bits = $class -> _new($class -> _num($xr) | $class -> _num($yr));
1152        $z = $class -> _add($z, $class -> _mul($bits, $m));
1153        $m = $class -> _mul($m, $mask);
1154    }
1155
1156    # The loop above stops when the smallest of the two numbers is exhausted.
1157    # The remainder of the longer one will survive bit-by-bit, so we simple
1158    # multiply-add it in.
1159
1160    $z = $class -> _add($z, $class -> _mul($xc, $m))
1161      unless $class -> _is_zero($xc);
1162    $z = $class -> _add($z, $class -> _mul($yc, $m))
1163      unless $class -> _is_zero($yc);
1164
1165    return $z;
1166}
1167
1168sub _sand {
1169    my ($class, $x, $sx, $y, $sy) = @_;
1170
1171    return ($class -> _zero(), '+')
1172      if $class -> _is_zero($x) || $class -> _is_zero($y);
1173
1174    my $sign = $sx eq '-' && $sy eq '-' ? '-' : '+';
1175
1176    my ($bx, $by);
1177
1178    if ($sx eq '-') {                   # if x is negative
1179        # two's complement: inc (dec unsigned value) and flip all "bits" in $bx
1180        $bx = $class -> _copy($x);
1181        $bx = $class -> _dec($bx);
1182        $bx = $class -> _as_hex($bx);
1183        $bx =~ s/^-?0x//;
1184        $bx =~ tr<0123456789abcdef>
1185                <\x0f\x0e\x0d\x0c\x0b\x0a\x09\x08\x07\x06\x05\x04\x03\x02\x01\x00>;
1186    } else {                            # if x is positive
1187        $bx = $class -> _as_hex($x);    # get binary representation
1188        $bx =~ s/^-?0x//;
1189        $bx =~ tr<fedcba9876543210>
1190                 <\x0f\x0e\x0d\x0c\x0b\x0a\x09\x08\x07\x06\x05\x04\x03\x02\x01\x00>;
1191    }
1192
1193    if ($sy eq '-') {                   # if y is negative
1194        # two's complement: inc (dec unsigned value) and flip all "bits" in $by
1195        $by = $class -> _copy($y);
1196        $by = $class -> _dec($by);
1197        $by = $class -> _as_hex($by);
1198        $by =~ s/^-?0x//;
1199        $by =~ tr<0123456789abcdef>
1200                <\x0f\x0e\x0d\x0c\x0b\x0a\x09\x08\x07\x06\x05\x04\x03\x02\x01\x00>;
1201    } else {
1202        $by = $class -> _as_hex($y);    # get binary representation
1203        $by =~ s/^-?0x//;
1204        $by =~ tr<fedcba9876543210>
1205                <\x0f\x0e\x0d\x0c\x0b\x0a\x09\x08\x07\x06\x05\x04\x03\x02\x01\x00>;
1206    }
1207
1208    # now we have bit-strings from X and Y, reverse them for padding
1209    $bx = reverse $bx;
1210    $by = reverse $by;
1211
1212    # padd the shorter string
1213    my $xx = "\x00"; $xx = "\x0f" if $sx eq '-';
1214    my $yy = "\x00"; $yy = "\x0f" if $sy eq '-';
1215    my $diff = CORE::length($bx) - CORE::length($by);
1216    if ($diff > 0) {
1217        # if $yy eq "\x00", we can cut $bx, otherwise we need to padd $by
1218        $by .= $yy x $diff;
1219    } elsif ($diff < 0) {
1220        # if $xx eq "\x00", we can cut $by, otherwise we need to padd $bx
1221        $bx .= $xx x abs($diff);
1222    }
1223
1224    # and the strings together
1225    my $r = $bx & $by;
1226
1227    # and reverse the result again
1228    $bx = reverse $r;
1229
1230    # One of $bx or $by was negative, so need to flip bits in the result. In both
1231    # cases (one or two of them negative, or both positive) we need to get the
1232    # characters back.
1233    if ($sign eq '-') {
1234        $bx =~ tr<\x0f\x0e\x0d\x0c\x0b\x0a\x09\x08\x07\x06\x05\x04\x03\x02\x01\x00>
1235                 <0123456789abcdef>;
1236    } else {
1237        $bx =~ tr<\x0f\x0e\x0d\x0c\x0b\x0a\x09\x08\x07\x06\x05\x04\x03\x02\x01\x00>
1238                 <fedcba9876543210>;
1239    }
1240
1241    # leading zeros will be stripped by _from_hex()
1242    $bx = '0x' . $bx;
1243    $bx = $class -> _from_hex($bx);
1244
1245    $bx = $class -> _inc($bx) if $sign eq '-';
1246
1247    # avoid negative zero
1248    $sign = '+' if $class -> _is_zero($bx);
1249
1250    return $bx, $sign;
1251}
1252
1253sub _sxor {
1254    my ($class, $x, $sx, $y, $sy) = @_;
1255
1256    return ($class -> _zero(), '+')
1257      if $class -> _is_zero($x) && $class -> _is_zero($y);
1258
1259    my $sign = $sx ne $sy ? '-' : '+';
1260
1261    my ($bx, $by);
1262
1263    if ($sx eq '-') {                   # if x is negative
1264        # two's complement: inc (dec unsigned value) and flip all "bits" in $bx
1265        $bx = $class -> _copy($x);
1266        $bx = $class -> _dec($bx);
1267        $bx = $class -> _as_hex($bx);
1268        $bx =~ s/^-?0x//;
1269        $bx =~ tr<0123456789abcdef>
1270                <\x0f\x0e\x0d\x0c\x0b\x0a\x09\x08\x07\x06\x05\x04\x03\x02\x01\x00>;
1271    } else {                            # if x is positive
1272        $bx = $class -> _as_hex($x);    # get binary representation
1273        $bx =~ s/^-?0x//;
1274        $bx =~ tr<fedcba9876543210>
1275                 <\x0f\x0e\x0d\x0c\x0b\x0a\x09\x08\x07\x06\x05\x04\x03\x02\x01\x00>;
1276    }
1277
1278    if ($sy eq '-') {                   # if y is negative
1279        # two's complement: inc (dec unsigned value) and flip all "bits" in $by
1280        $by = $class -> _copy($y);
1281        $by = $class -> _dec($by);
1282        $by = $class -> _as_hex($by);
1283        $by =~ s/^-?0x//;
1284        $by =~ tr<0123456789abcdef>
1285                <\x0f\x0e\x0d\x0c\x0b\x0a\x09\x08\x07\x06\x05\x04\x03\x02\x01\x00>;
1286    } else {
1287        $by = $class -> _as_hex($y);    # get binary representation
1288        $by =~ s/^-?0x//;
1289        $by =~ tr<fedcba9876543210>
1290                <\x0f\x0e\x0d\x0c\x0b\x0a\x09\x08\x07\x06\x05\x04\x03\x02\x01\x00>;
1291    }
1292
1293    # now we have bit-strings from X and Y, reverse them for padding
1294    $bx = reverse $bx;
1295    $by = reverse $by;
1296
1297    # padd the shorter string
1298    my $xx = "\x00"; $xx = "\x0f" if $sx eq '-';
1299    my $yy = "\x00"; $yy = "\x0f" if $sy eq '-';
1300    my $diff = CORE::length($bx) - CORE::length($by);
1301    if ($diff > 0) {
1302        # if $yy eq "\x00", we can cut $bx, otherwise we need to padd $by
1303        $by .= $yy x $diff;
1304    } elsif ($diff < 0) {
1305        # if $xx eq "\x00", we can cut $by, otherwise we need to padd $bx
1306        $bx .= $xx x abs($diff);
1307    }
1308
1309    # xor the strings together
1310    my $r = $bx ^ $by;
1311
1312    # and reverse the result again
1313    $bx = reverse $r;
1314
1315    # One of $bx or $by was negative, so need to flip bits in the result. In both
1316    # cases (one or two of them negative, or both positive) we need to get the
1317    # characters back.
1318    if ($sign eq '-') {
1319        $bx =~ tr<\x0f\x0e\x0d\x0c\x0b\x0a\x09\x08\x07\x06\x05\x04\x03\x02\x01\x00>
1320                 <0123456789abcdef>;
1321    } else {
1322        $bx =~ tr<\x0f\x0e\x0d\x0c\x0b\x0a\x09\x08\x07\x06\x05\x04\x03\x02\x01\x00>
1323                 <fedcba9876543210>;
1324    }
1325
1326    # leading zeros will be stripped by _from_hex()
1327    $bx = '0x' . $bx;
1328    $bx = $class -> _from_hex($bx);
1329
1330    $bx = $class -> _inc($bx) if $sign eq '-';
1331
1332    # avoid negative zero
1333    $sign = '+' if $class -> _is_zero($bx);
1334
1335    return $bx, $sign;
1336}
1337
1338sub _sor {
1339    my ($class, $x, $sx, $y, $sy) = @_;
1340
1341    return ($class -> _zero(), '+')
1342      if $class -> _is_zero($x) && $class -> _is_zero($y);
1343
1344    my $sign = $sx eq '-' || $sy eq '-' ? '-' : '+';
1345
1346    my ($bx, $by);
1347
1348    if ($sx eq '-') {                   # if x is negative
1349        # two's complement: inc (dec unsigned value) and flip all "bits" in $bx
1350        $bx = $class -> _copy($x);
1351        $bx = $class -> _dec($bx);
1352        $bx = $class -> _as_hex($bx);
1353        $bx =~ s/^-?0x//;
1354        $bx =~ tr<0123456789abcdef>
1355                <\x0f\x0e\x0d\x0c\x0b\x0a\x09\x08\x07\x06\x05\x04\x03\x02\x01\x00>;
1356    } else {                            # if x is positive
1357        $bx = $class -> _as_hex($x);     # get binary representation
1358        $bx =~ s/^-?0x//;
1359        $bx =~ tr<fedcba9876543210>
1360                 <\x0f\x0e\x0d\x0c\x0b\x0a\x09\x08\x07\x06\x05\x04\x03\x02\x01\x00>;
1361    }
1362
1363    if ($sy eq '-') {                   # if y is negative
1364        # two's complement: inc (dec unsigned value) and flip all "bits" in $by
1365        $by = $class -> _copy($y);
1366        $by = $class -> _dec($by);
1367        $by = $class -> _as_hex($by);
1368        $by =~ s/^-?0x//;
1369        $by =~ tr<0123456789abcdef>
1370                <\x0f\x0e\x0d\x0c\x0b\x0a\x09\x08\x07\x06\x05\x04\x03\x02\x01\x00>;
1371    } else {
1372        $by = $class -> _as_hex($y);     # get binary representation
1373        $by =~ s/^-?0x//;
1374        $by =~ tr<fedcba9876543210>
1375                <\x0f\x0e\x0d\x0c\x0b\x0a\x09\x08\x07\x06\x05\x04\x03\x02\x01\x00>;
1376    }
1377
1378    # now we have bit-strings from X and Y, reverse them for padding
1379    $bx = reverse $bx;
1380    $by = reverse $by;
1381
1382    # padd the shorter string
1383    my $xx = "\x00"; $xx = "\x0f" if $sx eq '-';
1384    my $yy = "\x00"; $yy = "\x0f" if $sy eq '-';
1385    my $diff = CORE::length($bx) - CORE::length($by);
1386    if ($diff > 0) {
1387        # if $yy eq "\x00", we can cut $bx, otherwise we need to padd $by
1388        $by .= $yy x $diff;
1389    } elsif ($diff < 0) {
1390        # if $xx eq "\x00", we can cut $by, otherwise we need to padd $bx
1391        $bx .= $xx x abs($diff);
1392    }
1393
1394    # or the strings together
1395    my $r = $bx | $by;
1396
1397    # and reverse the result again
1398    $bx = reverse $r;
1399
1400    # One of $bx or $by was negative, so need to flip bits in the result. In both
1401    # cases (one or two of them negative, or both positive) we need to get the
1402    # characters back.
1403    if ($sign eq '-') {
1404        $bx =~ tr<\x0f\x0e\x0d\x0c\x0b\x0a\x09\x08\x07\x06\x05\x04\x03\x02\x01\x00>
1405                 <0123456789abcdef>;
1406    } else {
1407        $bx =~ tr<\x0f\x0e\x0d\x0c\x0b\x0a\x09\x08\x07\x06\x05\x04\x03\x02\x01\x00>
1408                 <fedcba9876543210>;
1409    }
1410
1411    # leading zeros will be stripped by _from_hex()
1412    $bx = '0x' . $bx;
1413    $bx = $class -> _from_hex($bx);
1414
1415    $bx = $class -> _inc($bx) if $sign eq '-';
1416
1417    # avoid negative zero
1418    $sign = '+' if $class -> _is_zero($bx);
1419
1420    return $bx, $sign;
1421}
1422
1423sub _to_bin {
1424    # convert the number to a string of binary digits without prefix
1425    my ($class, $x) = @_;
1426    my $str    = '';
1427    my $tmp    = $class -> _copy($x);
1428    my $chunk = $class -> _new("16777216");     # 2^24 = 24 binary digits
1429    my $rem;
1430    until ($class -> _acmp($tmp, $chunk) < 0) {
1431        ($tmp, $rem) = $class -> _div($tmp, $chunk);
1432        $str = sprintf("%024b", $class -> _num($rem)) . $str;
1433    }
1434    unless ($class -> _is_zero($tmp)) {
1435        $str = sprintf("%b", $class -> _num($tmp)) . $str;
1436    }
1437    return length($str) ? $str : '0';
1438}
1439
1440sub _to_oct {
1441    # convert the number to a string of octal digits without prefix
1442    my ($class, $x) = @_;
1443    my $str    = '';
1444    my $tmp    = $class -> _copy($x);
1445    my $chunk = $class -> _new("16777216");     # 2^24 = 8 octal digits
1446    my $rem;
1447    until ($class -> _acmp($tmp, $chunk) < 0) {
1448        ($tmp, $rem) = $class -> _div($tmp, $chunk);
1449        $str = sprintf("%08o", $class -> _num($rem)) . $str;
1450    }
1451    unless ($class -> _is_zero($tmp)) {
1452        $str = sprintf("%o", $class -> _num($tmp)) . $str;
1453    }
1454    return length($str) ? $str : '0';
1455}
1456
1457sub _to_hex {
1458    # convert the number to a string of hexadecimal digits without prefix
1459    my ($class, $x) = @_;
1460    my $str    = '';
1461    my $tmp    = $class -> _copy($x);
1462    my $chunk = $class -> _new("16777216");     # 2^24 = 6 hexadecimal digits
1463    my $rem;
1464    until ($class -> _acmp($tmp, $chunk) < 0) {
1465        ($tmp, $rem) = $class -> _div($tmp, $chunk);
1466        $str = sprintf("%06x", $class -> _num($rem)) . $str;
1467    }
1468    unless ($class -> _is_zero($tmp)) {
1469        $str = sprintf("%x", $class -> _num($tmp)) . $str;
1470    }
1471    return length($str) ? $str : '0';
1472}
1473
1474sub _as_bin {
1475    # convert the number to a string of binary digits with prefix
1476    my ($class, $x) = @_;
1477    return '0b' . $class -> _to_bin($x);
1478}
1479
1480sub _as_oct {
1481    # convert the number to a string of octal digits with prefix
1482    my ($class, $x) = @_;
1483    return '0' . $class -> _to_oct($x);         # yes, 0 becomes "00"
1484}
1485
1486sub _as_hex {
1487    # convert the number to a string of hexadecimal digits with prefix
1488    my ($class, $x) = @_;
1489    return '0x' . $class -> _to_hex($x);
1490}
1491
1492sub _to_bytes {
1493    # convert the number to a string of bytes
1494    my ($class, $x) = @_;
1495    my $str    = '';
1496    my $tmp    = $class -> _copy($x);
1497    my $chunk = $class -> _new("65536");
1498    my $rem;
1499    until ($class -> _is_zero($tmp)) {
1500        ($tmp, $rem) = $class -> _div($tmp, $chunk);
1501        $str = pack('n', $class -> _num($rem)) . $str;
1502    }
1503    $str =~ s/^\0+//;
1504    return length($str) ? $str : "\x00";
1505}
1506
1507*_as_bytes = \&_to_bytes;
1508
1509sub _to_base {
1510    # convert the number to a string of digits in various bases
1511    my $class = shift;
1512    my $x     = shift;
1513    my $base  = shift;
1514    $base = $class -> _new($base) unless ref($base);
1515
1516    my $collseq;
1517    if (@_) {
1518        $collseq = shift;
1519        croak "The collation sequence must be a non-empty string"
1520          unless defined($collseq) && length($collseq);
1521    } else {
1522        if ($class -> _acmp($base, $class -> _new("94")) <= 0) {
1523            $collseq = '0123456789'                     #  48 ..  57
1524                     . 'ABCDEFGHIJKLMNOPQRSTUVWXYZ'     #  65 ..  90
1525                     . 'abcdefghijklmnopqrstuvwxyz'     #  97 .. 122
1526                     . '!"#$%&\'()*+,-./'               #  33 ..  47
1527                     . ':;<=>?@'                        #  58 ..  64
1528                     . '[\\]^_`'                        #  91 ..  96
1529                     . '{|}~';                          # 123 .. 126
1530        } else {
1531            croak "When base > 94, a collation sequence must be given";
1532        }
1533    }
1534
1535    my @collseq = split '', $collseq;
1536
1537    my $str   = '';
1538    my $tmp   = $class -> _copy($x);
1539    my $rem;
1540    until ($class -> _is_zero($tmp)) {
1541        ($tmp, $rem) = $class -> _div($tmp, $base);
1542        my $num = $class -> _num($rem);
1543        croak "no character to represent '$num' in collation sequence",
1544          " (collation sequence is too short)" if $num > $#collseq;
1545        my $chr = $collseq[$num];
1546        $str = $chr . $str;
1547    }
1548    return $collseq[0] unless length $str;
1549    return $str;
1550}
1551
1552sub _to_base_num {
1553    # Convert the number to an array of integers in any base.
1554    my ($class, $x, $base) = @_;
1555
1556    # Make sure the base is an object and >= 2.
1557    $base = $class -> _new($base) unless ref($base);
1558    my $two = $class -> _two();
1559    croak "base must be >= 2" unless $class -> _acmp($base, $two) >= 0;
1560
1561    my $out   = [];
1562    my $xcopy = $class -> _copy($x);
1563    my $rem;
1564
1565    # Do all except the last (most significant) element.
1566    until ($class -> _acmp($xcopy, $base) < 0) {
1567        ($xcopy, $rem) = $class -> _div($xcopy, $base);
1568        unshift @$out, $rem;
1569    }
1570
1571    # Do the last (most significant element).
1572    unless ($class -> _is_zero($xcopy)) {
1573        unshift @$out, $xcopy;
1574    }
1575
1576    # $out is empty if $x is zero.
1577    unshift @$out, $class -> _zero() unless @$out;
1578
1579    return $out;
1580}
1581
1582sub _from_hex {
1583    # Convert a string of hexadecimal digits to a number.
1584
1585    my ($class, $hex) = @_;
1586    $hex =~ s/^0[xX]//;
1587
1588    # Find the largest number of hexadecimal digits that we can safely use with
1589    # 32 bit integers. There are 4 bits pr hexadecimal digit, and we use only
1590    # 31 bits to play safe. This gives us int(31 / 4) = 7.
1591
1592    my $len = length $hex;
1593    my $rem = 1 + ($len - 1) % 7;
1594
1595    # Do the first chunk.
1596
1597    my $ret = $class -> _new(int hex substr $hex, 0, $rem);
1598    return $ret if $rem == $len;
1599
1600    # Do the remaining chunks, if any.
1601
1602    my $shift = $class -> _new(1 << (4 * 7));
1603    for (my $offset = $rem ; $offset < $len ; $offset += 7) {
1604        my $part = int hex substr $hex, $offset, 7;
1605        $ret = $class -> _mul($ret, $shift);
1606        $ret = $class -> _add($ret, $class -> _new($part));
1607    }
1608
1609    return $ret;
1610}
1611
1612sub _from_oct {
1613    # Convert a string of octal digits to a number.
1614
1615    my ($class, $oct) = @_;
1616
1617    # Find the largest number of octal digits that we can safely use with 32
1618    # bit integers. There are 3 bits pr octal digit, and we use only 31 bits to
1619    # play safe. This gives us int(31 / 3) = 10.
1620
1621    my $len = length $oct;
1622    my $rem = 1 + ($len - 1) % 10;
1623
1624    # Do the first chunk.
1625
1626    my $ret = $class -> _new(int oct substr $oct, 0, $rem);
1627    return $ret if $rem == $len;
1628
1629    # Do the remaining chunks, if any.
1630
1631    my $shift = $class -> _new(1 << (3 * 10));
1632    for (my $offset = $rem ; $offset < $len ; $offset += 10) {
1633        my $part = int oct substr $oct, $offset, 10;
1634        $ret = $class -> _mul($ret, $shift);
1635        $ret = $class -> _add($ret, $class -> _new($part));
1636    }
1637
1638    return $ret;
1639}
1640
1641sub _from_bin {
1642    # Convert a string of binary digits to a number.
1643
1644    my ($class, $bin) = @_;
1645    $bin =~ s/^0[bB]//;
1646
1647    # The largest number of binary digits that we can safely use with 32 bit
1648    # integers is 31. We use only 31 bits to play safe.
1649
1650    my $len = length $bin;
1651    my $rem = 1 + ($len - 1) % 31;
1652
1653    # Do the first chunk.
1654
1655    my $ret = $class -> _new(int oct '0b' . substr $bin, 0, $rem);
1656    return $ret if $rem == $len;
1657
1658    # Do the remaining chunks, if any.
1659
1660    my $shift = $class -> _new(1 << 31);
1661    for (my $offset = $rem ; $offset < $len ; $offset += 31) {
1662        my $part = int oct '0b' . substr $bin, $offset, 31;
1663        $ret = $class -> _mul($ret, $shift);
1664        $ret = $class -> _add($ret, $class -> _new($part));
1665    }
1666
1667    return $ret;
1668}
1669
1670sub _from_bytes {
1671    # convert string of bytes to a number
1672    my ($class, $str) = @_;
1673    my $x    = $class -> _zero();
1674    my $base = $class -> _new("256");
1675    my $n    = length($str);
1676    for (my $i = 0 ; $i < $n ; ++$i) {
1677        $x = $class -> _mul($x, $base);
1678        my $byteval = $class -> _new(unpack 'C', substr($str, $i, 1));
1679        $x = $class -> _add($x, $byteval);
1680    }
1681    return $x;
1682}
1683
1684sub _from_base {
1685    # convert a string to a decimal number
1686    my $class = shift;
1687    my $str   = shift;
1688    my $base  = shift;
1689    $base = $class -> _new($base) unless ref($base);
1690
1691    my $n = length($str);
1692    my $x = $class -> _zero();
1693
1694    my $collseq;
1695    if (@_) {
1696        $collseq = shift();
1697    } else {
1698        if ($class -> _acmp($base, $class -> _new("36")) <= 0) {
1699            $str = uc $str;
1700            $collseq = '0123456789' . 'ABCDEFGHIJKLMNOPQRSTUVWXYZ';
1701        } elsif ($class -> _acmp($base, $class -> _new("94")) <= 0) {
1702            $collseq = '0123456789'                     #  48 ..  57
1703                     . 'ABCDEFGHIJKLMNOPQRSTUVWXYZ'     #  65 ..  90
1704                     . 'abcdefghijklmnopqrstuvwxyz'     #  97 .. 122
1705                     . '!"#$%&\'()*+,-./'               #  33 ..  47
1706                     . ':;<=>?@'                        #  58 ..  64
1707                     . '[\\]^_`'                        #  91 ..  96
1708                     . '{|}~';                          # 123 .. 126
1709        } else {
1710            croak "When base > 94, a collation sequence must be given";
1711        }
1712        $collseq = substr $collseq, 0, $class -> _num($base);
1713    }
1714
1715    # Create a mapping from each character in the collation sequence to the
1716    # corresponding integer. Check for duplicates in the collation sequence.
1717
1718    my @collseq = split '', $collseq;
1719    my %collseq;
1720    for my $num (0 .. $#collseq) {
1721        my $chr = $collseq[$num];
1722        die "duplicate character '$chr' in collation sequence"
1723          if exists $collseq{$chr};
1724        $collseq{$chr} = $num;
1725    }
1726
1727    for (my $i = 0 ; $i < $n ; ++$i) {
1728        my $chr = substr($str, $i, 1);
1729        die "input character '$chr' does not exist in collation sequence"
1730          unless exists $collseq{$chr};
1731        $x = $class -> _mul($x, $base);
1732        my $num = $class -> _new($collseq{$chr});
1733        $x = $class -> _add($x, $num);
1734    }
1735
1736    return $x;
1737}
1738
1739sub _from_base_num {
1740    # Convert an array in the given base to a number.
1741    my ($class, $in, $base) = @_;
1742
1743    # Make sure the base is an object and >= 2.
1744    $base = $class -> _new($base) unless ref($base);
1745    my $two = $class -> _two();
1746    croak "base must be >= 2" unless $class -> _acmp($base, $two) >= 0;
1747
1748    # @$in = map { ref($_) ? $_ : $class -> _new($_) } @$in;
1749
1750    my $ele = $in -> [0];
1751
1752    $ele = $class -> _new($ele) unless ref($ele);
1753    my $x = $class -> _copy($ele);
1754
1755    for my $i (1 .. $#$in) {
1756        $x = $class -> _mul($x, $base);
1757        $ele = $in -> [$i];
1758        $ele = $class -> _new($ele) unless ref($ele);
1759        $x = $class -> _add($x, $ele);
1760    }
1761
1762    return $x;
1763}
1764
1765##############################################################################
1766# special modulus functions
1767
1768sub _modinv {
1769    # modular multiplicative inverse
1770    my ($class, $x, $y) = @_;
1771
1772    # modulo zero
1773    if ($class -> _is_zero($y)) {
1774        return;
1775    }
1776
1777    # modulo one
1778    if ($class -> _is_one($y)) {
1779        return ($class -> _zero(), '+');
1780    }
1781
1782    my $u = $class -> _zero();
1783    my $v = $class -> _one();
1784    my $a = $class -> _copy($y);
1785    my $b = $class -> _copy($x);
1786
1787    # Euclid's Algorithm for bgcd().
1788
1789    my $q;
1790    my $sign = 1;
1791    {
1792        ($a, $q, $b) = ($b, $class -> _div($a, $b));
1793        last if $class -> _is_zero($b);
1794
1795        my $vq = $class -> _mul($class -> _copy($v), $q);
1796        my $t = $class -> _add($vq, $u);
1797        $u = $v;
1798        $v = $t;
1799        $sign = -$sign;
1800        redo;
1801    }
1802
1803    # if the gcd is not 1, there exists no modular multiplicative inverse
1804    return unless $class -> _is_one($a);
1805
1806    ($v, $sign == 1 ? '+' : '-');
1807}
1808
1809sub _modpow {
1810    # modulus of power ($x ** $y) % $z
1811    my ($class, $num, $exp, $mod) = @_;
1812
1813    # a^b (mod 1) = 0 for all a and b
1814    if ($class -> _is_one($mod)) {
1815        return $class -> _zero();
1816    }
1817
1818    # 0^a (mod m) = 0 if m != 0, a != 0
1819    # 0^0 (mod m) = 1 if m != 0
1820    if ($class -> _is_zero($num)) {
1821        return $class -> _is_zero($exp) ? $class -> _one()
1822                                        : $class -> _zero();
1823    }
1824
1825    #  $num = $class -> _mod($num, $mod);   # this does not make it faster
1826
1827    my $acc = $class -> _copy($num);
1828    my $t   = $class -> _one();
1829
1830    my $expbin = $class -> _as_bin($exp);
1831    $expbin =~ s/^0b//;
1832    my $len = length($expbin);
1833
1834    while (--$len >= 0) {
1835        if (substr($expbin, $len, 1) eq '1') {
1836            $t = $class -> _mul($t, $acc);
1837            $t = $class -> _mod($t, $mod);
1838        }
1839        $acc = $class -> _mul($acc, $acc);
1840        $acc = $class -> _mod($acc, $mod);
1841    }
1842    return $t;
1843}
1844
1845sub _gcd {
1846    # Greatest common divisor.
1847
1848    my ($class, $x, $y) = @_;
1849
1850    # gcd(0, 0) = 0
1851    # gcd(0, a) = a, if a != 0
1852
1853    if ($class -> _acmp($x, $y) == 0) {
1854        return $class -> _copy($x);
1855    }
1856
1857    if ($class -> _is_zero($x)) {
1858        if ($class -> _is_zero($y)) {
1859            return $class -> _zero();
1860        } else {
1861            return $class -> _copy($y);
1862        }
1863    } else {
1864        if ($class -> _is_zero($y)) {
1865            return $class -> _copy($x);
1866        } else {
1867
1868            # Until $y is zero ...
1869
1870            $x = $class -> _copy($x);
1871            until ($class -> _is_zero($y)) {
1872
1873                # Compute remainder.
1874
1875                $x = $class -> _mod($x, $y);
1876
1877                # Swap $x and $y.
1878
1879                my $tmp = $x;
1880                $x = $class -> _copy($y);
1881                $y = $tmp;
1882            }
1883
1884            return $x;
1885        }
1886    }
1887}
1888
1889sub _lcm {
1890    # Least common multiple.
1891
1892    my ($class, $x, $y) = @_;
1893
1894    # lcm(0, x) = 0 for all x
1895
1896    return $class -> _zero()
1897      if ($class -> _is_zero($x) ||
1898          $class -> _is_zero($y));
1899
1900    my $gcd = $class -> _gcd($class -> _copy($x), $y);
1901    $x = $class -> _div($x, $gcd);
1902    $x = $class -> _mul($x, $y);
1903    return $x;
1904}
1905
1906sub _lucas {
1907    my ($class, $n) = @_;
1908
1909    $n = $class -> _num($n) if ref $n;
1910
1911    # In list context, use lucas(n) = lucas(n-1) + lucas(n-2)
1912
1913    if (wantarray) {
1914        my @y;
1915
1916        push @y, $class -> _two();
1917        return @y if $n == 0;
1918
1919        push @y, $class -> _one();
1920        return @y if $n == 1;
1921
1922        for (my $i = 2 ; $i <= $n ; ++ $i) {
1923            $y[$i] = $class -> _add($class -> _copy($y[$i - 1]), $y[$i - 2]);
1924        }
1925
1926        return @y;
1927    }
1928
1929    # In scalar context use that lucas(n) = fib(n-1) + fib(n+1).
1930    #
1931    # Remember that _fib() behaves differently in scalar context and list
1932    # context, so we must add scalar() to get the desired behaviour.
1933
1934    return $class -> _two() if $n == 0;
1935
1936    return $class -> _add(scalar($class -> _fib($n - 1)),
1937                          scalar($class -> _fib($n + 1)));
1938}
1939
1940sub _fib {
1941    my ($class, $n) = @_;
1942
1943    $n = $class -> _num($n) if ref $n;
1944
1945    # In list context, use fib(n) = fib(n-1) + fib(n-2)
1946
1947    if (wantarray) {
1948        my @y;
1949
1950        push @y, $class -> _zero();
1951        return @y if $n == 0;
1952
1953        push @y, $class -> _one();
1954        return @y if $n == 1;
1955
1956        for (my $i = 2 ; $i <= $n ; ++ $i) {
1957            $y[$i] = $class -> _add($class -> _copy($y[$i - 1]), $y[$i - 2]);
1958        }
1959
1960        return @y;
1961    }
1962
1963    # In scalar context use a fast algorithm that is much faster than the
1964    # recursive algorith used in list context.
1965
1966    my $cache = {};
1967    my $two = $class -> _two();
1968    my $fib;
1969
1970    $fib = sub {
1971        my $n = shift;
1972        return $class -> _zero() if $n <= 0;
1973        return $class -> _one()  if $n <= 2;
1974        return $cache -> {$n}    if exists $cache -> {$n};
1975
1976        my $k = int($n / 2);
1977        my $a = $fib -> ($k + 1);
1978        my $b = $fib -> ($k);
1979        my $y;
1980
1981        if ($n % 2 == 1) {
1982            # a*a + b*b
1983            $y = $class -> _add($class -> _mul($class -> _copy($a), $a),
1984                                $class -> _mul($class -> _copy($b), $b));
1985        } else {
1986            # (2*a - b)*b
1987            $y = $class -> _mul($class -> _sub($class -> _mul(
1988                   $class -> _copy($two), $a), $b), $b);
1989        }
1990
1991        $cache -> {$n} = $y;
1992        return $y;
1993    };
1994
1995    return $fib -> ($n);
1996}
1997
1998##############################################################################
1999##############################################################################
2000
20011;
2002
2003__END__
2004
2005=pod
2006
2007=head1 NAME
2008
2009Math::BigInt::Lib - virtual parent class for Math::BigInt libraries
2010
2011=head1 SYNOPSIS
2012
2013    # In the backend library for Math::BigInt et al.
2014
2015    package Math::BigInt::MyBackend;
2016
2017    use Math::BigInt::Lib;
2018    our @ISA = qw< Math::BigInt::Lib >;
2019
2020    sub _new { ... }
2021    sub _str { ... }
2022    sub _add { ... }
2023    str _sub { ... }
2024    ...
2025
2026    # In your main program.
2027
2028    use Math::BigInt lib => 'MyBackend';
2029
2030=head1 DESCRIPTION
2031
2032This module provides support for big integer calculations. It is not intended
2033to be used directly, but rather as a parent class for backend libraries used by
2034Math::BigInt, Math::BigFloat, Math::BigRat, and related modules.
2035
2036Other backend libraries include Math::BigInt::Calc, Math::BigInt::FastCalc,
2037Math::BigInt::GMP, and Math::BigInt::Pari.
2038
2039In order to allow for multiple big integer libraries, Math::BigInt was
2040rewritten to use a plug-in library for core math routines. Any module which
2041conforms to the API can be used by Math::BigInt by using this in your program:
2042
2043        use Math::BigInt lib => 'libname';
2044
2045'libname' is either the long name, like 'Math::BigInt::Pari', or only the short
2046version, like 'Pari'.
2047
2048=head2 General Notes
2049
2050A library only needs to deal with unsigned big integers. Testing of input
2051parameter validity is done by the caller, so there is no need to worry about
2052underflow (e.g., in C<_sub()> and C<_dec()>) or about division by zero (e.g.,
2053in C<_div()> and C<_mod()>)) or similar cases.
2054
2055Some libraries use methods that don't modify their argument, and some libraries
2056don't even use objects, but rather unblessed references. Because of this,
2057liberary methods are always called as class methods, not instance methods:
2058
2059    $x = Class -> method($x, $y);     # like this
2060    $x = $x -> method($y);            # not like this ...
2061    $x -> method($y);                 # ... or like this
2062
2063And with boolean methods
2064
2065    $bool = Class -> method($x, $y);  # like this
2066    $bool = $x -> method($y);         # not like this
2067
2068Return values are always objects, strings, Perl scalars, or true/false for
2069comparison routines.
2070
2071=head3 API version
2072
2073=over 4
2074
2075=item CLASS-E<gt>api_version()
2076
2077This method is no longer used and can be omitted. Methods that are not
2078implemented by a subclass will be inherited from this class.
2079
2080=back
2081
2082=head3 Constructors
2083
2084The following methods are mandatory: _new(), _str(), _add(), and _sub().
2085However, computations will be very slow without _mul() and _div().
2086
2087=over 4
2088
2089=item CLASS-E<gt>_new(STR)
2090
2091Convert a string representing an unsigned decimal number to an object
2092representing the same number. The input is normalized, i.e., it matches
2093C<^(0|[1-9]\d*)$>.
2094
2095=item CLASS-E<gt>_zero()
2096
2097Return an object representing the number zero.
2098
2099=item CLASS-E<gt>_one()
2100
2101Return an object representing the number one.
2102
2103=item CLASS-E<gt>_two()
2104
2105Return an object representing the number two.
2106
2107=item CLASS-E<gt>_ten()
2108
2109Return an object representing the number ten.
2110
2111=item CLASS-E<gt>_from_bin(STR)
2112
2113Return an object given a string representing a binary number. The input has a
2114'0b' prefix and matches the regular expression C<^0[bB](0|1[01]*)$>.
2115
2116=item CLASS-E<gt>_from_oct(STR)
2117
2118Return an object given a string representing an octal number. The input has a
2119'0' prefix and matches the regular expression C<^0[1-7]*$>.
2120
2121=item CLASS-E<gt>_from_hex(STR)
2122
2123Return an object given a string representing a hexadecimal number. The input
2124has a '0x' prefix and matches the regular expression
2125C<^0x(0|[1-9a-fA-F][\da-fA-F]*)$>.
2126
2127=item CLASS-E<gt>_from_bytes(STR)
2128
2129Returns an object given a byte string representing the number. The byte string
2130is in big endian byte order, so the two-byte input string "\x01\x00" should
2131give an output value representing the number 256.
2132
2133=item CLASS-E<gt>_from_base(STR, BASE, COLLSEQ)
2134
2135Returns an object given a string STR, a base BASE, and a collation sequence
2136COLLSEQ. Each character in STR represents a numerical value identical to the
2137character's position in COLLSEQ. All characters in STR must be present in
2138COLLSEQ.
2139
2140If BASE is less than or equal to 94, and a collation sequence is not specified,
2141the following default collation sequence is used. It contains of all the 94
2142printable ASCII characters except space/blank:
2143
2144    0123456789                  # ASCII  48 to  57
2145    ABCDEFGHIJKLMNOPQRSTUVWXYZ  # ASCII  65 to  90
2146    abcdefghijklmnopqrstuvwxyz  # ASCII  97 to 122
2147    !"#$%&'()*+,-./             # ASCII  33 to  47
2148    :;<=>?@                     # ASCII  58 to  64
2149    [\]^_`                      # ASCII  91 to  96
2150    {|}~                        # ASCII 123 to 126
2151
2152If the default collation sequence is used, and the BASE is less than or equal
2153to 36, the letter case in STR is ignored.
2154
2155For instance, with base 3 and collation sequence "-/|", the character "-"
2156represents 0, "/" represents 1, and "|" represents 2. So if STR is "/|-", the
2157output is 1 * 3**2 + 2 * 3**1 + 0 * 3**0 = 15.
2158
2159The following examples show standard binary, octal, decimal, and hexadecimal
2160conversion. All examples return 250.
2161
2162    $x = $class -> _from_base("11111010", 2)
2163    $x = $class -> _from_base("372", 8)
2164    $x = $class -> _from_base("250", 10)
2165    $x = $class -> _from_base("FA", 16)
2166
2167Some more examples, all returning 250:
2168
2169    $x = $class -> _from_base("100021", 3)
2170    $x = $class -> _from_base("3322", 4)
2171    $x = $class -> _from_base("2000", 5)
2172    $x = $class -> _from_base("caaa", 5, "abcde")
2173    $x = $class -> _from_base("42", 62)
2174    $x = $class -> _from_base("2!", 94)
2175
2176=item CLASS-E<gt>_from_base_num(ARRAY, BASE)
2177
2178Returns an object given an array of values and a base. This method is
2179equivalent to C<_from_base()>, but works on numbers in an array rather than
2180characters in a string. Unlike C<_from_base()>, all input values may be
2181arbitrarily large.
2182
2183    $x = $class -> _from_base_num([1, 1, 0, 1], 2)    # $x is 13
2184    $x = $class -> _from_base_num([3, 125, 39], 128)  # $x is 65191
2185
2186=back
2187
2188=head3 Mathematical functions
2189
2190=over 4
2191
2192=item CLASS-E<gt>_add(OBJ1, OBJ2)
2193
2194Addition. Returns the result of adding OBJ2 to OBJ1.
2195
2196=item CLASS-E<gt>_mul(OBJ1, OBJ2)
2197
2198Multiplication. Returns the result of multiplying OBJ2 and OBJ1.
2199
2200=item CLASS-E<gt>_div(OBJ1, OBJ2)
2201
2202Division. In scalar context, returns the quotient after dividing OBJ1 by OBJ2
2203and truncating the result to an integer. In list context, return the quotient
2204and the remainder.
2205
2206=item CLASS-E<gt>_sub(OBJ1, OBJ2, FLAG)
2207
2208=item CLASS-E<gt>_sub(OBJ1, OBJ2)
2209
2210Subtraction. Returns the result of subtracting OBJ2 by OBJ1. If C<flag> is false
2211or omitted, OBJ1 might be modified. If C<flag> is true, OBJ2 might be modified.
2212
2213=item CLASS-E<gt>_sadd(OBJ1, SIGN1, OBJ2, SIGN2)
2214
2215Signed addition. Returns the result of adding OBJ2 with sign SIGN2 to OBJ1 with
2216sign SIGN1.
2217
2218    ($obj3, $sign3) = $class -> _sadd($obj1, $sign1, $obj2, $sign2);
2219
2220=item CLASS-E<gt>_ssub(OBJ1, SIGN1, OBJ2, SIGN2)
2221
2222Signed subtraction. Returns the result of subtracting OBJ2 with sign SIGN2 to
2223OBJ1 with sign SIGN1.
2224
2225    ($obj3, $sign3) = $class -> _sadd($obj1, $sign1, $obj2, $sign2);
2226
2227=item CLASS-E<gt>_dec(OBJ)
2228
2229Returns the result after decrementing OBJ by one.
2230
2231=item CLASS-E<gt>_inc(OBJ)
2232
2233Returns the result after incrementing OBJ by one.
2234
2235=item CLASS-E<gt>_mod(OBJ1, OBJ2)
2236
2237Returns OBJ1 modulo OBJ2, i.e., the remainder after dividing OBJ1 by OBJ2.
2238
2239=item CLASS-E<gt>_sqrt(OBJ)
2240
2241Returns the square root of OBJ, truncated to an integer.
2242
2243=item CLASS-E<gt>_root(OBJ, N)
2244
2245Returns the Nth root of OBJ, truncated to an integer.
2246
2247=item CLASS-E<gt>_fac(OBJ)
2248
2249Returns the factorial of OBJ, i.e., the product of all positive integers up to
2250and including OBJ.
2251
2252=item CLASS-E<gt>_dfac(OBJ)
2253
2254Returns the double factorial of OBJ. If OBJ is an even integer, returns the
2255product of all positive, even integers up to and including OBJ, i.e.,
22562*4*6*...*OBJ. If OBJ is an odd integer, returns the product of all positive,
2257odd integers, i.e., 1*3*5*...*OBJ.
2258
2259=item CLASS-E<gt>_pow(OBJ1, OBJ2)
2260
2261Returns OBJ1 raised to the power of OBJ2. By convention, 0**0 = 1.
2262
2263=item CLASS-E<gt>_modinv(OBJ1, OBJ2)
2264
2265Returns the modular multiplicative inverse, i.e., return OBJ3 so that
2266
2267    (OBJ3 * OBJ1) % OBJ2 = 1 % OBJ2
2268
2269The result is returned as two arguments. If the modular multiplicative inverse
2270does not exist, both arguments are undefined. Otherwise, the arguments are a
2271number (object) and its sign ("+" or "-").
2272
2273The output value, with its sign, must either be a positive value in the range
22741,2,...,OBJ2-1 or the same value subtracted OBJ2. For instance, if the input
2275arguments are objects representing the numbers 7 and 5, the method must either
2276return an object representing the number 3 and a "+" sign, since (3*7) % 5 = 1
2277% 5, or an object representing the number 2 and a "-" sign, since (-2*7) % 5 = 1
2278% 5.
2279
2280=item CLASS-E<gt>_modpow(OBJ1, OBJ2, OBJ3)
2281
2282Returns the modular exponentiation, i.e., (OBJ1 ** OBJ2) % OBJ3.
2283
2284=item CLASS-E<gt>_rsft(OBJ, N, B)
2285
2286Returns the result after shifting OBJ N digits to thee right in base B. This is
2287equivalent to performing integer division by B**N and discarding the remainder,
2288except that it might be much faster.
2289
2290For instance, if the object $obj represents the hexadecimal number 0xabcde,
2291then C<_rsft($obj, 2, 16)> returns an object representing the number 0xabc. The
2292"remainer", 0xde, is discarded and not returned.
2293
2294=item CLASS-E<gt>_lsft(OBJ, N, B)
2295
2296Returns the result after shifting OBJ N digits to the left in base B. This is
2297equivalent to multiplying by B**N, except that it might be much faster.
2298
2299=item CLASS-E<gt>_log_int(OBJ, B)
2300
2301Returns the logarithm of OBJ to base BASE truncted to an integer. This method
2302has two output arguments, the OBJECT and a STATUS. The STATUS is Perl scalar;
2303it is 1 if OBJ is the exact result, 0 if the result was truncted to give OBJ,
2304and undef if it is unknown whether OBJ is the exact result.
2305
2306=item CLASS-E<gt>_gcd(OBJ1, OBJ2)
2307
2308Returns the greatest common divisor of OBJ1 and OBJ2.
2309
2310=item CLASS-E<gt>_lcm(OBJ1, OBJ2)
2311
2312Return the least common multiple of OBJ1 and OBJ2.
2313
2314=item CLASS-E<gt>_fib(OBJ)
2315
2316In scalar context, returns the nth Fibonacci number: _fib(0) returns 0, _fib(1)
2317returns 1, _fib(2) returns 1, _fib(3) returns 2 etc. In list context, returns
2318the Fibonacci numbers from F(0) to F(n): 0, 1, 1, 2, 3, 5, 8, 13, 21, 34, ...
2319
2320=item CLASS-E<gt>_lucas(OBJ)
2321
2322In scalar context, returns the nth Lucas number: _lucas(0) returns 2, _lucas(1)
2323returns 1, _lucas(2) returns 3, etc. In list context, returns the Lucas numbers
2324from L(0) to L(n): 2, 1, 3, 4, 7, 11, 18, 29,47, 76, ...
2325
2326=back
2327
2328=head3 Bitwise operators
2329
2330=over 4
2331
2332=item CLASS-E<gt>_and(OBJ1, OBJ2)
2333
2334Returns bitwise and.
2335
2336=item CLASS-E<gt>_or(OBJ1, OBJ2)
2337
2338Returns bitwise or.
2339
2340=item CLASS-E<gt>_xor(OBJ1, OBJ2)
2341
2342Returns bitwise exclusive or.
2343
2344=item CLASS-E<gt>_sand(OBJ1, OBJ2, SIGN1, SIGN2)
2345
2346Returns bitwise signed and.
2347
2348=item CLASS-E<gt>_sor(OBJ1, OBJ2, SIGN1, SIGN2)
2349
2350Returns bitwise signed or.
2351
2352=item CLASS-E<gt>_sxor(OBJ1, OBJ2, SIGN1, SIGN2)
2353
2354Returns bitwise signed exclusive or.
2355
2356=back
2357
2358=head3 Boolean operators
2359
2360=over 4
2361
2362=item CLASS-E<gt>_is_zero(OBJ)
2363
2364Returns a true value if OBJ is zero, and false value otherwise.
2365
2366=item CLASS-E<gt>_is_one(OBJ)
2367
2368Returns a true value if OBJ is one, and false value otherwise.
2369
2370=item CLASS-E<gt>_is_two(OBJ)
2371
2372Returns a true value if OBJ is two, and false value otherwise.
2373
2374=item CLASS-E<gt>_is_ten(OBJ)
2375
2376Returns a true value if OBJ is ten, and false value otherwise.
2377
2378=item CLASS-E<gt>_is_even(OBJ)
2379
2380Return a true value if OBJ is an even integer, and a false value otherwise.
2381
2382=item CLASS-E<gt>_is_odd(OBJ)
2383
2384Return a true value if OBJ is an even integer, and a false value otherwise.
2385
2386=item CLASS-E<gt>_acmp(OBJ1, OBJ2)
2387
2388Compare OBJ1 and OBJ2 and return -1, 0, or 1, if OBJ1 is numerically less than,
2389equal to, or larger than OBJ2, respectively.
2390
2391=back
2392
2393=head3 String conversion
2394
2395=over 4
2396
2397=item CLASS-E<gt>_str(OBJ)
2398
2399Returns a string representing OBJ in decimal notation. The returned string
2400should have no leading zeros, i.e., it should match C<^(0|[1-9]\d*)$>.
2401
2402=item CLASS-E<gt>_to_bin(OBJ)
2403
2404Returns the binary string representation of OBJ.
2405
2406=item CLASS-E<gt>_to_oct(OBJ)
2407
2408Returns the octal string representation of the number.
2409
2410=item CLASS-E<gt>_to_hex(OBJ)
2411
2412Returns the hexadecimal string representation of the number.
2413
2414=item CLASS-E<gt>_to_bytes(OBJ)
2415
2416Returns a byte string representation of OBJ. The byte string is in big endian
2417byte order, so if OBJ represents the number 256, the output should be the
2418two-byte string "\x01\x00".
2419
2420=item CLASS-E<gt>_to_base(OBJ, BASE, COLLSEQ)
2421
2422Returns a string representation of OBJ in base BASE with collation sequence
2423COLLSEQ.
2424
2425    $val = $class -> _new("210");
2426    $str = $class -> _to_base($val, 10, "xyz")  # $str is "zyx"
2427
2428    $val = $class -> _new("32");
2429    $str = $class -> _to_base($val, 2, "-|")  # $str is "|-----"
2430
2431See _from_base() for more information.
2432
2433=item CLASS-E<gt>_to_base_num(OBJ, BASE)
2434
2435Converts the given number to the given base. This method is equivalent to
2436C<_to_base()>, but returns numbers in an array rather than characters in a
2437string. In the output, the first element is the most significant. Unlike
2438C<_to_base()>, all input values may be arbitrarily large.
2439
2440    $x = $class -> _to_base_num(13, 2)        # $x is [1, 1, 0, 1]
2441    $x = $class -> _to_base_num(65191, 128)   # $x is [3, 125, 39]
2442
2443=item CLASS-E<gt>_as_bin(OBJ)
2444
2445Like C<_to_bin()> but with a '0b' prefix.
2446
2447=item CLASS-E<gt>_as_oct(OBJ)
2448
2449Like C<_to_oct()> but with a '0' prefix.
2450
2451=item CLASS-E<gt>_as_hex(OBJ)
2452
2453Like C<_to_hex()> but with a '0x' prefix.
2454
2455=item CLASS-E<gt>_as_bytes(OBJ)
2456
2457This is an alias to C<_to_bytes()>.
2458
2459=back
2460
2461=head3 Numeric conversion
2462
2463=over 4
2464
2465=item CLASS-E<gt>_num(OBJ)
2466
2467Returns a Perl scalar number representing the number OBJ as close as
2468possible. Since Perl scalars have limited precision, the returned value might
2469not be exactly the same as OBJ.
2470
2471=back
2472
2473=head3 Miscellaneous
2474
2475=over 4
2476
2477=item CLASS-E<gt>_copy(OBJ)
2478
2479Returns a true copy OBJ.
2480
2481=item CLASS-E<gt>_len(OBJ)
2482
2483Returns the number of the decimal digits in OBJ. The output is a Perl scalar.
2484
2485=item CLASS-E<gt>_zeros(OBJ)
2486
2487Returns the number of trailing decimal zeros. The output is a Perl scalar. The
2488number zero has no trailing decimal zeros.
2489
2490=item CLASS-E<gt>_digit(OBJ, N)
2491
2492Returns the Nth digit in OBJ as a Perl scalar. N is a Perl scalar, where zero
2493refers to the rightmost (least significant) digit, and negative values count
2494from the left (most significant digit). If $obj represents the number 123, then
2495
2496    CLASS->_digit($obj,  0)     # returns 3
2497    CLASS->_digit($obj,  1)     # returns 2
2498    CLASS->_digit($obj,  2)     # returns 1
2499    CLASS->_digit($obj, -1)     # returns 1
2500
2501=item CLASS-E<gt>_digitsum(OBJ)
2502
2503Returns the sum of the base 10 digits.
2504
2505=item CLASS-E<gt>_check(OBJ)
2506
2507Returns true if the object is invalid and false otherwise. Preferably, the true
2508value is a string describing the problem with the object. This is a check
2509routine to test the internal state of the object for corruption.
2510
2511=item CLASS-E<gt>_set(OBJ)
2512
2513xxx
2514
2515=back
2516
2517=head2 API version 2
2518
2519The following methods are required for an API version of 2 or greater.
2520
2521=head3 Constructors
2522
2523=over 4
2524
2525=item CLASS-E<gt>_1ex(N)
2526
2527Return an object representing the number 10**N where N E<gt>= 0 is a Perl
2528scalar.
2529
2530=back
2531
2532=head3 Mathematical functions
2533
2534=over 4
2535
2536=item CLASS-E<gt>_nok(OBJ1, OBJ2)
2537
2538Return the binomial coefficient OBJ1 over OBJ1.
2539
2540=back
2541
2542=head3 Miscellaneous
2543
2544=over 4
2545
2546=item CLASS-E<gt>_alen(OBJ)
2547
2548Return the approximate number of decimal digits of the object. The output is a
2549Perl scalar.
2550
2551=back
2552
2553=head1 WRAP YOUR OWN
2554
2555If you want to port your own favourite C library for big numbers to the
2556Math::BigInt interface, you can take any of the already existing modules as a
2557rough guideline. You should really wrap up the latest Math::BigInt and
2558Math::BigFloat testsuites with your module, and replace in them any of the
2559following:
2560
2561        use Math::BigInt;
2562
2563by this:
2564
2565        use Math::BigInt lib => 'yourlib';
2566
2567This way you ensure that your library really works 100% within Math::BigInt.
2568
2569=head1 BUGS
2570
2571Please report any bugs or feature requests to
2572C<bug-math-bigint at rt.cpan.org>, or through the web interface at
2573L<https://rt.cpan.org/Ticket/Create.html?Queue=Math-BigInt>
2574(requires login).
2575We will be notified, and then you'll automatically be notified of progress on
2576your bug as I make changes.
2577
2578=head1 SUPPORT
2579
2580You can find documentation for this module with the perldoc command.
2581
2582    perldoc Math::BigInt::Calc
2583
2584You can also look for information at:
2585
2586=over 4
2587
2588=item * RT: CPAN's request tracker
2589
2590L<https://rt.cpan.org/Public/Dist/Display.html?Name=Math-BigInt>
2591
2592=item * AnnoCPAN: Annotated CPAN documentation
2593
2594L<http://annocpan.org/dist/Math-BigInt>
2595
2596=item * CPAN Ratings
2597
2598L<https://cpanratings.perl.org/dist/Math-BigInt>
2599
2600=item * MetaCPAN
2601
2602L<https://metacpan.org/release/Math-BigInt>
2603
2604=item * CPAN Testers Matrix
2605
2606L<http://matrix.cpantesters.org/?dist=Math-BigInt>
2607
2608=item * The Bignum mailing list
2609
2610=over 4
2611
2612=item * Post to mailing list
2613
2614C<bignum at lists.scsys.co.uk>
2615
2616=item * View mailing list
2617
2618L<http://lists.scsys.co.uk/pipermail/bignum/>
2619
2620=item * Subscribe/Unsubscribe
2621
2622L<http://lists.scsys.co.uk/cgi-bin/mailman/listinfo/bignum>
2623
2624=back
2625
2626=back
2627
2628=head1 LICENSE
2629
2630This program is free software; you may redistribute it and/or modify it under
2631the same terms as Perl itself.
2632
2633=head1 AUTHOR
2634
2635Peter John Acklam, E<lt>pjacklam@gmail.comE<gt>
2636
2637Code and documentation based on the Math::BigInt::Calc module by Tels
2638E<lt>nospam-abuse@bloodgate.comE<gt>
2639
2640=head1 SEE ALSO
2641
2642L<Math::BigInt>, L<Math::BigInt::Calc>, L<Math::BigInt::GMP>,
2643L<Math::BigInt::FastCalc> and L<Math::BigInt::Pari>.
2644
2645=cut
2646