123456789101112131415161718192021222324252627282930313233343536373839404142434445464748495051525354555657585960616263646566676869707172737475767778798081828384858687888990 |
- #!/usr/bin/perl
- # Daniel "Trizen" Șuteu
- # Date: 06 March 2019
- # https://github.com/trizen
- # Generalized algorithm for generating numbers that are smooth over a set A of primes, bellow a given limit.
- use 5.020;
- use warnings;
- use experimental qw(signatures);
- use Math::GMPz;
- use ntheory qw(:all);
- sub check_valuation ($n, $p) {
- ($n*$p)%2 == 0 or return 0;
- #~ if ($p > 13) {
- #~ return ( ($n % $p) != 0);
- #~ }
- 1;
- }
- sub smooth_numbers ($limit, $primes) {
- my @h = (1);
- foreach my $p (@$primes) {
- say "Prime: $p";
- foreach my $n (@h) {
- if ($n * $p <= $limit and check_valuation($n, $p)) {
- push @h, $n * $p;
- }
- }
- }
- return \@h;
- }
- #
- # Example for finding numbers `m` such that:
- # sigma(m) * phi(m) = n^k
- # for some `n` and `k`, with `n > 1` and `k > 1`.
- #
- # See also: https://oeis.org/A306724
- #
- my $t = 282669887501;
- my $base = 5;
- sub isok ($n) {
- #is_power(Math::GMPz->new(divisor_sum($n)) * euler_phi($n));
- my $k = $n+1;
- if ($k < $t) {
- return 0;
- }
- powmod($base, $k - 1, Math::GMPz->new($t) * $k) == 1;
- }
- my $h = smooth_numbers(2402694043751, primes(63));
- say "\nFound: ", scalar(@$h), " terms";
- my @list;
- foreach my $n (@$h) {
- my $p = isok($n);
- #if ($p >= 8) {
- if ($p) {
- say "Found: ", $n+1, " -> ", join(' * ', map { "$_->[0]^$_->[1]" } factor_exp($n)), ' -> ', is_prime($n+1) ? 'prime' : 'NOT PRIME';
- #push @{$table{$p}}, $n;
- #say "Found: ", $n+1;
- push @list, $n+1;
- }
- }
- say vecmin(@list);
- __END__
- a(29) = 282669887501
|