Skip to content

Commit

Permalink
- Added the Number inverse_psi_len(n) method.
Browse files Browse the repository at this point in the history 
Returns the number of solutions to Dedekind's psi function: psi(x) = n.
Equivalent to `n.inverse_psi.len`, but much faster.

- Added the Number `inverse_sigma_len(n,k=1)` method.

Returns the number of solutions to the sigma sum of divisors function: sigma_k(x) = n.
Equivalent to `n.inverse_sigma(k).len`, but much faster.

- Optimized `inverse_phi_len(n)` to use the new Number "_dynamic_length" private function.
  • Loading branch information
@trizen
trizen committed Jan 15, 2021
1 parent 9d0d882 commit 91b851f
Show file tree
Hide file tree
Showing 3 changed files with 116 additions and 14 deletions.
103 changes: 89 additions & 14 deletions lib/Sidef/Types/Number/Number.pm
View file
Original file line number Diff line number Diff line change
Expand Up @@ -13544,6 +13544,8 @@ package Sidef::Types::Number::Number {
sub _dynamic_preimage {
my ($N, $L, %opt) = @_;

# Based on invphi.gp ver. 2.1 by Max Alekseyev.

my %R = (1 => [$ONE]);
my $u = Math::GMPz::Rmpz_init();

Expand Down Expand Up @@ -13592,6 +13594,45 @@ package Sidef::Types::Number::Number {
$R{$N} // [];
}

sub _dynamic_length {
my ($N, $L) = @_;

# Based on invphi.gp ver. 2.1 by Max Alekseyev.

my %R = (1 => 1);
my $u = Math::GMPz::Rmpz_init();

foreach my $l (@$L) {
my %t;

foreach my $pair (@$l) {

my $x = $pair->[0];
Math::GMPz::Rmpz_divexact($u, $N, $x);

foreach my $d (_divisors($u)) {
if (exists $R{$d}) {

($d < ULONG_MAX)
? Math::GMPz::Rmpz_mul_ui($u, $x, $d)
: do {
Math::GMPz::Rmpz_set_str($u, $d, 10);
Math::GMPz::Rmpz_mul($u, $u, $x);
};

$t{Math::GMPz::Rmpz_get_str($u, 10)} += $R{$d};
}
}
}

while (my ($key, $value) = each %t) {
$R{$key} += $value;
}
}

$R{$N} // 0;
}

sub _cook_euler_phi {
my ($N) = @_;

Expand Down Expand Up @@ -13662,6 +13703,30 @@ package Sidef::Types::Number::Number {

*inverse_phi = \&inverse_totient;

sub inverse_totient_len {
my ($n) = @_;

$n = _any2mpz($$n) // return ZERO;

if (Math::GMPz::Rmpz_sgn($n) <= 0) {
return ONE if !Math::GMPz::Rmpz_sgn($n);
return ZERO;
}

my $r;

if (HAS_PRIME_UTIL and Math::GMPz::Rmpz_fits_ulong_p($n)) {
$r = scalar Math::Prime::Util::inverse_totient(Math::GMPz::Rmpz_get_ui($n));
}
else {
$r = _dynamic_length($n, _cook_euler_phi($n));
}

($r < ULONG_MAX) ? __PACKAGE__->_set_uint($r) : __PACKAGE__->_set_str('int', $r);
}

*inverse_phi_len = \&inverse_totient_len;

sub _cook_dedekind_psi {
my ($N, $k) = @_;

Expand Down Expand Up @@ -13718,24 +13783,21 @@ package Sidef::Types::Number::Number {

*inverse_psi = \&inverse_dedekind_psi;

sub inverse_totient_len {
sub inverse_dedekind_psi_len {
my ($n) = @_;

my $z = _any2mpz($$n) // return ZERO;
$n = _any2mpz($$n) // return ZERO;

if (Math::GMPz::Rmpz_sgn($z) <= 0) {
return ONE if !Math::GMPz::Rmpz_sgn($z);
if (Math::GMPz::Rmpz_sgn($n) <= 0) {
return ONE if !Math::GMPz::Rmpz_sgn($n);
return ZERO;
}

if (HAS_PRIME_UTIL) {
return Sidef::Types::Number::Number->_set_uint(scalar Math::Prime::Util::inverse_totient($z));
}

$n->inverse_totient->len;
my $r = _dynamic_length($n, _cook_dedekind_psi($n));
($r < ULONG_MAX) ? __PACKAGE__->_set_uint($r) : __PACKAGE__->_set_str('int', $r);
}

*inverse_phi_len = \&inverse_totient_len;
*inverse_psi_len = \&inverse_dedekind_psi_len;

sub _cook_usigma {
my ($N) = @_;
Expand Down Expand Up @@ -13816,6 +13878,8 @@ package Sidef::Types::Number::Number {
sub _cook_sigma {
my ($N, $k) = @_;

# Based on invphi.gp ver. 2.1 by Max Alekseyev.

my $q = Math::GMPz::Rmpz_init();
my $s = Math::GMPz::Rmpz_init();
my $u = Math::GMPz::Rmpz_init();
Expand Down Expand Up @@ -13879,14 +13943,25 @@ package Sidef::Types::Number::Number {
return Sidef::Types::Array::Array->new;
}

if (Math::GMPz::Rmpz_cmp_ui($n, 1) == 0) {
return Sidef::Types::Array::Array->new(ONE);
}

my $result = _dynamic_preimage($n, _cook_sigma($n, $k));
Sidef::Types::Array::Array->new([map { bless \$_ } sort { Math::GMPz::Rmpz_cmp($a, $b) } @$result]);
}

sub inverse_sigma_len {
my ($n, $k) = @_;

$n = _any2mpz($$n) // return ZERO;
$k = defined($k) ? do { _valid(\$k); _any2ui($$k) // return ZERO } : 1;

if (Math::GMPz::Rmpz_sgn($n) <= 0) {
return ONE if !Math::GMPz::Rmpz_sgn($n);
return ZERO;
}

my $r = _dynamic_length($n, _cook_sigma($n, $k));
($r < ULONG_MAX) ? __PACKAGE__->_set_uint($r) : __PACKAGE__->_set_str('int', $r);
}

sub jordan_totient {
my ($n, $k) = @_;
$k //= ONE;
Expand Down
18 changes: 18 additions & 0 deletions lib/Sidef/Types/Number/Number.pod
View file
Original file line number Diff line number Diff line change
Expand Up @@ -2162,6 +2162,16 @@ Aliases: I<inverse_dedekind_psi>

=cut

=head2 inverse_psi_len

Number.inverse_psi_len() -> I<Obj>

Return the

Aliases: I<inverse_dedekind_psi_len>

=cut

=head2 inverse_sigma

Number.inverse_sigma() -> I<Obj>
Expand All @@ -2170,6 +2180,14 @@ Return the

=cut

=head2 inverse_sigma_len

Number.inverse_sigma_len() -> I<Obj>

Return the

=cut

=head2 inverse_uphi

Number.inverse_uphi() -> I<Obj>
Expand Down
9 changes: 9 additions & 0 deletions scripts/Tests/number_methods.sf
View file
Original file line number Diff line number Diff line change
Expand Up @@ -39,6 +39,15 @@ assert_eq(20.of { .uphi }, 20.of { _ ? .factor_map { |p,e| p**e - 1 }.prod : 0
assert_eq(1..100 -> grep { .is_smooth(5) }, 1..100 -> grep { .gpf <= 5 })
assert_eq(2..100 -> grep { .is_rough(5) }, 2..100 -> grep { .lpf >= 5 })

assert_eq(50.of { .inverse_phi.len }, 50.of { .inverse_phi_len })
assert_eq(50.of { .inverse_psi.len }, 50.of { .inverse_psi_len })
assert_eq(50.of { .inverse_sigma.len }, 50.of { .inverse_sigma_len })
assert_eq(50.of { .inverse_sigma(2).len }, 50.of { .inverse_sigma_len(2) })

assert_eq(inverse_phi(2**64).len, inverse_phi_len(2**64))
assert_eq(inverse_psi(2**64).len, inverse_psi_len(2**64))
assert_eq(inverse_sigma(2**64).len, inverse_sigma_len(2**64))

do {
var a = Math.smooth_numbers(2,3,5,7) # 7-smooth numbers
var b = Math.smooth_numbers(2,5,7) # 7-smooth numbers not divisible by 3
Expand Down

0 comments on commit 91b851f

Please sign in to comment.