-
Notifications
You must be signed in to change notification settings - Fork 34
/
modular_fibonacci_cassini.pl
executable file
·54 lines (36 loc) · 1.26 KB
/
modular_fibonacci_cassini.pl
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
#!/usr/bin/perl
# An efficient algorithm for computing the nth-Fibonacci number (mod m).
# Algorithm from:
# https://metacpan.org/source/KRYDE/Math-NumSeq-72/lib/Math/NumSeq/Fibonacci.pm
# See also:
# https://en.wikipedia.org/wiki/Fibonacci_number
use 5.020;
use warnings;
use experimental qw(signatures);
use Math::GMPz;
use Math::Prime::Util::GMP qw(consecutive_integer_lcm gcd);
sub fibmod ($n, $m) {
$n = Math::GMPz->new("$n");
$m = Math::GMPz->new("$m");
my ($f, $g, $a) = (0, 1, -2);
foreach my $bit (split(//, substr(Math::GMPz::Rmpz_get_str($n, 2), 1))) {
($g *= $g) %= $m;
($f *= $f) %= $m;
my $t = ($g << 2) - $f + $a;
$f += $g;
if ($bit) {
($f, $g, $a) = ($t - $f, $t, -2);
}
else {
($g, $a) = ($t - $f, 2);
}
}
return ($g % $m);
}
sub fibonacci_factorization ($n, $B = 10000) {
my $k = consecutive_integer_lcm($B); # lcm(1..B)
my $F = fibmod($k, $n); # Fibonacci(k) (mod n)
return gcd($F, $n);
}
say fibonacci_factorization("121095274043", 700); #=> 470783 (p+1 is 700-smooth)
say fibonacci_factorization("544812320889004864776853", 3000); #=> 333732865481 (p-1 is 3000-smooth)