123456789101112131415161718192021222324252627282930313233343536373839404142434445464748495051525354555657585960616263646566676869707172737475767778798081828384858687888990919293949596979899100101102103104105106107108109110111112113114115116117118119120121122123124125126127128129130131132 |
- #!/usr/bin/perl
- # Daniel "Trizen" Șuteu
- # Date: 22 May 2020
- # https://github.com/trizen
- # Twos are all you need
- # https://projecteuler.net/problem=708
- # Solution by counting the number of k-almost primes <= n.
- # Let:
- # pi_k(n) = the number of k-almost primes <= n.
- # Then:
- # S(n) = Sum_{k=1..n} 2^bigomega(k)
- # = Sum_{k=1..floor(log_2(n))} 2^k * pi_k(n)
- # Runtime: 2 min, 23 sec
- use 5.020;
- use integer;
- use ntheory qw(:all);
- use experimental qw(signatures);
- my @pi_table;
- my $table_len = 1e5;
- {
- my $count = 0;
- foreach my $k(0..$table_len) {
- if (is_prime($k)) {
- ++$count;
- }
- $pi_table[$k] = $count;
- }
- }
- sub prime_count($n) {
- ($n <= $table_len) ? $pi_table[$n] : ntheory::prime_count($n);
- }
- sub k_prime_count ($n, $k, $i = 14) {
- if ($k == 1) {
- return [0, 4, 25, 168, 1229, 9592, 78498, 664579, 5761455, 50847534, 455052511, 4118054813, 37607912018, 346065536839, 3204941750802, 29844570422669, 279238341033925, 2623557157654233, 24739954287740860, 234057667276344607, 2220819602560918840, 21127269486018731928, 201467286689315906290, 1925320391606803968923, 18435599767349200867866, 176846309399143769411680, 1699246750872437141327603, 16352460426841680446427399]->[$i];
- }
- if ($k == 2) {
- return [0, 4, 34, 299, 2625, 23378, 210035, 1904324, 17427258, 160788536, 1493776443, 13959990342, 131126017178, 1237088048653, 11715902308080, 111329817298881, 1061057292827269, 10139482913717352, 97123037685177087, 932300026230174178, 8966605849641219022, 86389956293761485464]->[$i];
- }
- if ($k == 3) {
- return [0, 1, 22, 247, 2569, 25556, 250853, 2444359, 23727305, 229924367, 2227121996, 21578747909, 209214982913, 2030133769624, 19717814526785, 191693417109381, 1865380637252270, 18168907486812690, 177123437184971927, 1728190923820610000]->[$i];
- }
- if ($k == 4) {
- return [ 0, 0, 12, 149, 1712, 18744, 198062, 2050696, 20959322, 212385942, 2139236881, 21454599814, 214499908019, 2139634739326, 21306682904040, 211905511283590, 2105504493045818, 20905484578206982]->[$i];
- }
- if ($k == 5) {
- return [0, 0, 4, 76, 963, 11185, 124465, 1349779, 14371023, 150982388, 1570678136, 16218372618, 166497674684, 1701439985694, 17323079621014]->[$i];
- }
- if ($k == 6) {
- return [0, 0, 2, 37, 485, 5933, 68963, 774078, 8493366, 91683887, 977694273, 10327249593, 108264085934, 1128049914377, 11694704489580]->[$i];
- }
- if ($k == 7) {
- return [ 0, 0, 0, 14, 231, 2973, 35585, 409849, 4600247, 50678212, 550454756, 5913771637, 62981797962, 665997804082, 7001087934965]->[$i];
- }
- if ($k == 8) {
- return [ 0, 0, 0, 7, 105, 1418, 17572, 207207, 2367507, 26483012, 291646797, 3173159326, 34192782745, 365561221293, 3882841742380]->[$i];
- }
- if ($k == 9) {
- return [ 0, 0, 0, 2, 47, 671, 8491, 101787, 1180751, 13377156, 148930536, 1636170477, 17787688377, 191742524399, 2052389350029]->[$i];
- }
- if ($k == 10) {
- return [0, 0, 0, 0, 22, 306, 4016, 49163, 578154, 6618221, 74342563, 823164388, 9011965866, 97765974368, 1052666075366]->[$i];
- }
- if ($k == 11) {
- return [0, 0, 0, 0, 7, 138, 1878, 23448, 279286, 3230577, 36585097, 407818620, 4490844534, 48972151631, 529781669333]->[$i];
- }
- if ($k == 12) {
- return [0, 0, 0, 0, 3, 63, 865, 11068, 133862, 1563465, 17836903, 200051717, 2214357712, 24255601105, 263439785143]->[$i];
- }
- my $count = 0;
- sub ($m, $p, $r) {
- my $s = rootint(divint($n, $m), $r);
- if ($r == 2) {
- my $j = prime_count($p) - 2;
- forprimes {
- $count += (prime_count(divint($n, $m * $_)) - ++$j);
- } $p, $s;
- return;
- }
- for (my $q = $p ; $q <= $s ; $q = next_prime($q)) {
- __SUB__->($m * $q, $q, $r - 1);
- }
- }->(1, 2, $k);
- return $count;
- }
- sub S($i) {
- my $t = 1;
- my $n = powint(10, $i);
- foreach my $k(1..logint($n, 2)) {
- say "Computing: $k";
- $t += powint(2, $k) * k_prime_count($n, $k, $i);
- }
- return $t;
- }
- say "S(10^8) = ", S(8);
- say "S(10^14) = ", S(14);
|