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 Erastosthenes. 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;
	}
}