IBM Ponder This, March 2024
Solution to the IBM Ponder This problem for March 2024. Here is the problem statement, and here is the official solution.
X_1000 115192665 X_2024 117778830159
Approach
A prime sieve is initially populated by following Sieve of Eratosthenes. As an optimization, a bitset is used (instead of individual bytes) for 8-fold memory savings.
The whole sieve still took around 30 GB of RAM, which just about my system memory. An additional optimization could be populating the sieve contiguous segments at a time (instead of walking the whole memory space for each prime) to leverage RAM NUMA caching. (This is called the “Segmented Sieve” approach).
After that, a brute force approach is used, testing every integer starting from 1. The numbers seem to be all safely within the 64-bit range, so they can be efficiently computed by a modern 64-bit processor (without slow big integer libraries).
The code below uses a single core and takes around one hour of wall clock time on my machine to both compute the sieve and find the solutions for X_i, with i from 1 to 2024. The solution from the previous i-1 is the starting point for the brute force search for X_i to save on some recomputation.
Source code
#include <stdio.h>
#include <stdlib.h>
#include <string.h>
#include <stdint.h>
typedef int64_t i64;
typedef uint64_t u64;
// #define END 1000
// #define MAX_N 120000000LL
// #define SQRT_MAX_N 11000LL
#define END 2024
#define MAX_N 250000000000LL
#define SQRT_MAX_N 500000LL
#define BASE 64
u64 bitset[(MAX_N/BASE)+1];
#define bitset_test(X) (bitset[((u64)(X))/BASE] & ((u64)1 << (((u64)(X)) % BASE)))
#define bitset_set(X) do{ bitset[((u64)(X))/BASE] |= ((u64)1 << (((u64)(X)) % BASE)); }while(0)
#undef assert
#define assert(X) if(!(X)) { \
fprintf(stderr, "%s:%d: assert failed '%s'\n", __FILE__, __LINE__, #X); exit(1); }
static int is_prime(i64 x) {
assert(x < MAX_N);
return !bitset_test(x);
}
static int is_valid_fast(i64 n, i64 i) {
if (is_prime(i)) { return 0; }
i64 x = i;
for (i64 k = 1; k < n; k++) {
x = x + k;
if (is_prime(x)) { return 0; }
}
return 1;
}
static i64 solve_fast(i64 n, i64 hint) {
for (i64 i = hint; i < MAX_N; i++) {
if (is_valid_fast(n, i)) {
return i;
}
}
return -1;
}
static void populate_sieve() {
memset(bitset, 0, sizeof(bitset));
bitset_set(1);
for (i64 i = 2; i <= SQRT_MAX_N; i++) {
if (bitset_test(i)) { continue; }
for (i64 k = i*i; k <= MAX_N; k += i) {
bitset_set(k);
}
}
}
int main() {
fprintf(stderr, "Computing sieve up to %ld\n", MAX_N);
populate_sieve();
fprintf(stderr, "done\n");
i64 hint = 1;
for (int i = 1; i <= END; i++) {
i64 sol = solve_fast(i, hint);
if (hint != sol || i == END) {
printf("%d\t%ld\n", i, sol);
fflush(stdout);
}
hint = sol;
}
}