#include
#include
using namespace std;
typedef unsigned long long int ulong;
//======================================================================
// This program generates a b-file for the lucky numbers (A000959)
// Run "lucky " to print elements 1 through
// It uses a "virtual sieve" and requires O() space
// It runs as fast of faster than an explicit sieving algorithm
int main(int argc, char* argv[])
{
// Obtain number of elements from command line
unsigned elems = argc <= 1 ? 10000 : atoi(argv[1]);
// Create a vector for our seive
vector lucky;
lucky.resize(elems);
// Set and print the first two elements explicitly
// Indexing from 0 simplifies the computation
if (elems >= 1)
{
lucky[0] = 1;
cout << "1 1" << endl;
}
if (elems >= 2)
{
lucky[1] = 3;
cout << "2 3" << endl;
}
// g is the largest index with lucky[g] <= n+1
unsigned g = 0;
// Compute the nth lucky number for 2 <= n <= elems
for (unsigned n = 2; n < elems; n++)
{
// Update g to largest index with lucky[g] <= n+1
if (lucky[g+1] <= n+1) g++;
// Now we are going to trace the position k of the nth
// lucky number backwards through the sieving process.
// k is the nth lucky number, so it is at position n
// after all the sieves.
ulong k = n;
// If lucky[i] > n+1, the sieve on lucky[i] does not alter
// the position of the nth lucky number, that is, does not
// alter k. So we need to run backwards through the sieves
// for which lucky[i] <= n+1. The last such sieve is the
// sieve for lucky[g], by definition of g.
// So, we run backwards through the sieves for lucky[g]
// down to the sieve for lucky[1] = 3.
for (unsigned i = g; i >= 1; i--)
{
// Here k is the position of the nth lucky number
// after the sieve on lucky[i]. Adjust the position
// prior to the sieve on lucky[i].
k = k*lucky[i]/(lucky[i]-1);
}
// Here k is the position of the nth lucky number prior to
// sieve on 3, that is, after the sieve on 2. Adjust the
// position prior to the sieve on 2.
k = 2*k;
// Here k is the position of the nth lucky number prior to
// the sieve on 2, that is, within the natural numbers
// (1, 2, 3, ...) indexed from 0. So the nth lucky number is
lucky[n] = k+1;
// Adjust n for 1-indexing and print our new value
cout << n+1 << " " << lucky[n] << endl;
}
// And we are done
return 0;
}