Math-Prime-Util
Math-Prime-Util copied to clipboard
[feature request] Image of sigma() function
Hi. I would like to suggest a new function, called inverse_sigma(n,k=1)
, which returns a list with all the solutions x
that satisfy sigma_k(x) = n
.
Similar in functionality to inverse_totient(n)
.
Algorithm in Perl (based on invphi.gp ver. 2.1 by Max Alekseyev):
use 5.020;
use strict;
use warnings;
use ntheory qw(:all);
use experimental qw(signatures);
sub dynamicPreimage ($N, $L) {
my %r = (1 => [1]);
foreach my $l (@$L) {
my %t;
foreach my $pair (@$l) {
my ($x, $y) = @$pair;
foreach my $d (divisors(divint($N, $x))) {
if (exists $r{$d}) {
push @{$t{mulint($x, $d)}}, map { mulint($_, $y) } @{$r{$d}};
}
}
}
while (my ($k, $v) = each %t) {
push @{$r{$k}}, @$v;
}
}
return if not exists $r{$N};
sort { $a <=> $b } @{$r{$N}};
}
sub cook_sigma ($N, $k) {
my %L;
foreach my $d (divisors($N)) {
next if ($d == 1);
foreach my $p (map { $_->[0] } factor_exp(subint($d, 1))) {
my $q = addint(mulint($d, subint(powint($p, $k), 1)), 1);
my $t = valuation($q, $p);
next if ($t <= $k or ($t % $k) or $q != powint($p, $t));
push @{$L{$p}}, [$d, powint($p, subint(divint($t, $k), 1))];
}
}
[values %L];
}
sub inverse_sigma ($N, $k = 1) {
($N == 1) ? (1) : dynamicPreimage($N, cook_sigma($N, $k));
}
say "solutions for sigma_1(x) = 5040: ", "[", join(", ", inverse_sigma(5040)), "]";
say "solutions for sigma_2(x) = 22100: ", "[", join(", ", inverse_sigma(22100, 2)), "]";
Extra: algorithm for computing the image of the usigma(n)
function (OEIS: A034448):
use 5.020;
use strict;
use warnings;
use ntheory qw(:all);
use experimental qw(signatures);
sub inverse_usigma ($n) {
my %r = (1 => [1]);
foreach my $d (divisors($n)) {
my $D = subint($d, 1);
is_prime_power($D) || next;
my %temp;
foreach my $f (divisors(divint($n, $d))) {
if (exists $r{$f}) {
push @{$temp{mulint($f, $d)}}, map { mulint($D, $_) }
grep { gcd($D, $_) == 1 } @{$r{$f}};
}
}
while (my ($key, $value) = each(%temp)) {
push @{$r{$key}}, @$value;
}
}
return if not exists $r{$n};
return sort { $a <=> $b } @{$r{$n}};
}
say "Solutions for usigma(x) = 5040: ", "[", join(', ', inverse_usigma(5040)), "]";