123456789101112131415161718192021222324252627282930313233343536373839404142434445464748495051525354555657585960616263646566676869707172 |
- #!/usr/bin/perl
- # Daniel "Trizen" Șuteu
- # Date: 18 July 2016
- # Edit: 23 October 2017
- # https://github.com/trizen
- # Compute the inverse of n-factorial.
- # The function is defined only for factorial numbers.
- # It may return non-sense for non-factorials.
- # See also:
- # https://oeis.org/A090368
- use 5.010;
- use strict;
- use warnings;
- use experimental qw(signatures);
- use ntheory qw(valuation factor factorial);
- sub factorial_prime_pow ($n, $p) {
- my $count = 0;
- my $ppow = $p;
- while ($ppow <= $n) {
- $count += int($n / $ppow);
- $ppow *= $p;
- }
- return $count;
- }
- sub p_adic_inverse ($p, $k) {
- my $n = $k * ($p - 1);
- while (factorial_prime_pow($n, $p) < $k) {
- $n -= $n % $p;
- $n += $p;
- }
- return $n;
- }
- sub inverse_of_factorial ($f) {
- return 1 if $f == 1;
- my $t = valuation($f, 2); # largest power of 2 in f
- my $z = p_adic_inverse(2, $t); # smallest number z such that 2^t divides z!
- my $d = (factor($z + 1))[-1]; # largest factor of z+1
- if (valuation($f, $d) != factorial_prime_pow($z + 1, $d)) {
- return $z;
- }
- return $z + 1;
- }
- foreach my $n (1 .. 30) {
- my $f = factorial($n);
- my $i = inverse_of_factorial($f);
- say "$i! = $f";
- if ($i != $n) {
- die "error: $i != $n";
- }
- }
|