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