Skip to content

Instantly share code, notes, and snippets.

@Hermann-SW
Created January 3, 2026 21:21
Show Gist options
  • Select an option

  • Save Hermann-SW/459c5a4c59d7a810a4d7bdea4078db1f to your computer and use it in GitHub Desktop.

Select an option

Save Hermann-SW/459c5a4c59d7a810a4d7bdea4078db1f to your computer and use it in GitHub Desktop.
Determine pi(10^n) = pi_{k}(10^n) (for all k<=10^n)
// NOLINT(legal/copyright)
// g++ pi.cc -lgmp -lgmpxx -O3 -Wall -pedantic -Wextra -o pi
// (cpplinted and cppchecked)
//
#include <time.h>
#include <math.h>
#include <gmpxx.h>
#include <assert.h>
#include <iostream>
#include <cstdint>
int main(int argc, char *argv[]) {
mpz_class mx, a = 2;
if (argc < 2) {
std::cerr << "Format: ./pi exp\n";
exit(1);
}
uint64_t c = 0, u = atoi(argv[1]);
mpz_ui_pow_ui(mx.get_mpz_t(), 10, u);
do {
mpz_nextprime(a.get_mpz_t(), a.get_mpz_t());
++c;
}
while (a <= mx); // NOLINT
std::cerr << "pi(" << mx << ")=" << c << "\n";
return 0;
}
@Hermann-SW
Copy link
Author

Hermann-SW commented Jan 3, 2026

Above gist code is outperformed by simple sieve of Eratosthenes.
But GMP code does not need big memory, while sieve code does.

hermann@7950x:~/RSA_numbers_factored/c++$ time ./pi_ 8
5761455

real	0m0.735s
user	0m0.699s
sys	0m0.036s
hermann@7950x:~/RSA_numbers_factored/c++$ time ./pi_ 9
50847534

real	0m10.186s
user	0m9.879s
sys	0m0.306s
hermann@7950x:~/RSA_numbers_factored/c++$ time ./pi_ 10
455052511

real	1m51.612s
user	1m48.876s
sys	0m2.728s
hermann@7950x:~/RSA_numbers_factored/c++$ 
$ cat pi_.cc 
// NOLINT(legal/copyright)
// g++ pi_.cc -O3 -Wall -pedantic -Wextra -o pi_
// (cpplinted and cppchecked)
//
#include <string.h>

#include <iostream>
#include <cassert>

char *V;
typedef __int64_t num;

int main(int argc, char *argv[]) {
    assert(argc == 2);
    assert(atoi(argv[1]) >= 2);
    num mx = 1;
    for (int e = 1; e <= atoi(argv[1]); ++e) {
      mx *= 10;
    }

    V = new char[1 + mx];
    assert(V);
    explicit_bzero(V, 1+mx);

    num c = 0;
    for (num i = 2; i <= mx; ++i) {
      if (V[i] != 0)  continue;

      for (num j = 2*i; j <= mx; j += i)  V[j] = 1;
      ++c;
    }

    std::cout << c << "\n";

    return 0;
}
$

Sign up for free to join this conversation on GitHub. Already have an account? Sign in to comment