prog.pl 1.0 KB

1234567891011121314151617181920212223242526272829303132333435363738394041424344454647484950515253545556575859606162
  1. #!/usr/bin/perl
  2. # a(n) is the index of the smallest n-gonal pyramidal number divisible by exactly n n-gonal pyramidal numbers.
  3. # https://oeis.org/A358059
  4. # Known terms:
  5. # 6, 7, 20, 79, 90, 203, 972, 3135, 374, 283815, 31824, 2232, 10240, 144584, 70784
  6. use 5.020;
  7. use warnings;
  8. use ntheory qw(:all);
  9. use experimental qw(signatures);
  10. sub pyramidal ($k, $r) {
  11. divint(vecprod($k, ($k+1), ($r-2)*$k + (5-$r)), 6);
  12. }
  13. sub a($n) {
  14. my %table;
  15. for(my $k = 1; ; ++$k) {
  16. my $t = pyramidal($k, $n);
  17. undef $table{$t};
  18. my $count = 0;
  19. foreach my $d (divisors($t)) {
  20. if (exists $table{$d}) {
  21. ++$count;
  22. last if ($count > $n);
  23. }
  24. }
  25. if ($count == $n) {
  26. return $k;
  27. }
  28. }
  29. }
  30. foreach my $n (3..100) {
  31. say "a($n) = ", a($n);
  32. }
  33. __END__
  34. a(3) = 6
  35. a(4) = 7
  36. a(5) = 20
  37. a(6) = 79
  38. a(7) = 90
  39. a(8) = 203
  40. a(9) = 972
  41. a(10) = 3135
  42. a(11) = 374
  43. a(12) = 283815
  44. a(13) = 31824
  45. a(14) = 2232
  46. a(15) = 10240
  47. a(16) = 144584
  48. a(17) = 70784