A question was asked over at math.SE (here) about whether or not there are infinitely many superpalindromes. I'm going to rephrase the definition to a more suitable one for coding purposes:
Definition: A superpalindrome is a product of primes p(1) * p(2) * ... * p(k), where p(i)<=p(i+1) for 1<=i<=k-1, such that p(1) * p(2) * ... * p(r) is a palindrome for all 1 <= r <= k.
A natural question to ask, is whether or not there's an infinite number of superpalindromes with all of its prime factors <=N. Thus, we can use a depth-first search.
Here is my implementation in C using GMP, but I'm wondering if there is some ways to give non-trivial run-time improvements.
#include <stdio.h>
#include <gmp.h>
// search for superpalindromes containing prime factors in {2,3,...,q}
// where q is the MAX_NR_PRIMES-th prime
#define MAX_NR_PRIMES 20000
#define MAX_SEARCH_DEPTH 100
#define MAX_NR_DIGITS 10000
// start backtracking at the (STARTING_PRIME_ID+1)-th prime
#define STARTING_PRIME_ID 0
mpz_t list_of_small_primes[MAX_NR_PRIMES];
mpz_t n_new[MAX_SEARCH_DEPTH];
int max_depth_reached=0;
mpz_t nr_superpalindromes_found[MAX_SEARCH_DEPTH];
// checks if n is a palindrome
// idea "borrowed" from http://gmplib.org/list-archives/gmp-discuss/2012-February/004876.html
int is_palindrome(mpz_t n) {
char m[MAX_NR_DIGITS];
int len=gmp_sprintf(m,"%Zd",n);
for(int i=0;i<len;i++) {
if(m[i]!=m[len-i-1]) return 0;
}
return 1;
}
// depth-first search for superpalindromes
// searches for a prime p:=list_of_small_primes[i_new], with i_new>=i
// such that n*p is a palindrome; continues searching if one is found
int extend_superpalindrome_backtracking_algorithm(mpz_t n,int i,int depth) {
for(int i_new=i;i_new<MAX_NR_PRIMES;i_new++) {
// n_new[depth]:=n*p
mpz_mul(n_new[depth],n,list_of_small_primes[i_new]);
if(is_palindrome(n_new[depth])) {
// increment count of number of superpalindromes with depth+1 prime factors
mpz_add_ui(nr_superpalindromes_found[depth],nr_superpalindromes_found[depth],1);
// print out the first superpalindrome found with >depth+1 prime factors
if(depth>max_depth_reached) {
max_depth_reached=depth;
gmp_printf("superpalindrome found: %Zd at depth %d: ",n_new[depth],depth+1);
mpz_t primes[MAX_SEARCH_DEPTH];
for(int i=0;i<MAX_SEARCH_DEPTH;i++) mpz_init(primes[i]);
mpz_set(primes[0],n_new[0]);
for(int i=1;i<=depth;i++) mpz_divexact(primes[i],n_new[i],n_new[i-1]);
for(int i=0;i<=depth;i++) gmp_printf("%Zd ",primes[i]);
printf("\n");
for(int i=0;i<MAX_SEARCH_DEPTH;i++) mpz_clear(primes[i]);
}
// continue the depth-first search if superpalindrome found
extend_superpalindrome_backtracking_algorithm(n_new[depth],i_new,depth+1);
}
}
}
int main() {
mpz_t n,p;
mpz_init_set_ui(n,1);
mpz_init_set(p,list_of_small_primes[STARTING_PRIME_ID]);
for(int i=0;i<MAX_SEARCH_DEPTH;i++) mpz_init(n_new[i]);
for(int i=0;i<MAX_SEARCH_DEPTH;i++) mpz_init(nr_superpalindromes_found[i]);
for(int i=0;i<MAX_NR_PRIMES;i++) mpz_init(list_of_small_primes[i]);
// pre-compute small primes
mpz_set_ui(list_of_small_primes[0],2);
for(int i=1;i<MAX_NR_PRIMES;i++) mpz_nextprime(list_of_small_primes[i],list_of_small_primes[i-1]);
extend_superpalindrome_backtracking_algorithm(n,STARTING_PRIME_ID,0);
// output results
gmp_printf("\nfound all superpalindromes with prime factors in {2,3,...,%Zd}\n",list_of_small_primes[MAX_NR_PRIMES-1]);
mpz_t total_nr_superpalindromes;
mpz_init_set_ui(total_nr_superpalindromes,0);
for(int i=0;i<=max_depth_reached;i++) {
mpz_add(total_nr_superpalindromes,total_nr_superpalindromes,nr_superpalindromes_found[i]);
gmp_printf("nr superpalindromes with %d prime factors: %Zd\n",i+1,nr_superpalindromes_found[i]);
}
gmp_printf("total nr superpalindromes: %Zd\n",total_nr_superpalindromes);
mpz_clear(total_nr_superpalindromes);
for(int i=0;i<MAX_SEARCH_DEPTH;i++) mpz_clear(n_new[i]);
for(int i=0;i<MAX_NR_PRIMES;i++) mpz_clear(list_of_small_primes[i]);
mpz_clear(n);
mpz_clear(p);
return 0;
}
Here's the output:
superpalindrome found: 4 at depth 2: 2 2
superpalindrome found: 8 at depth 3: 2 2 2
superpalindrome found: 88 at depth 4: 2 2 2 11
superpalindrome found: 2552 at depth 5: 2 2 2 11 29
superpalindrome found: 257752 at depth 6: 2 2 2 11 29 101
superpalindrome found: 67788776 at depth 7: 2 2 2 11 29 101 263
superpalindrome found: 616267762616 at depth 8: 2 2 2 11 29 101 263 9091
superpalindrome found: 6101667117661016 at depth 9: 2 2 2 11 29 101 263 9091 9901
superpalindrome found: 20302629368699686392620302 at depth 10: 2 11 11 11 101 9091 9091 9091 9901 10151
superpalindrome found: 1001004004006006004004001001 at depth 11: 7 11 13 101 101 101 101 9901 9901 9901 9901
found all superpalindromes with prime factors in {2,3,...,224737}
nr superpalindromes with 1 prime factors: 113
nr superpalindromes with 2 prime factors: 428
nr superpalindromes with 3 prime factors: 1022
nr superpalindromes with 4 prime factors: 1539
nr superpalindromes with 5 prime factors: 1603
nr superpalindromes with 6 prime factors: 1137
nr superpalindromes with 7 prime factors: 565
nr superpalindromes with 8 prime factors: 217
nr superpalindromes with 9 prime factors: 50
nr superpalindromes with 10 prime factors: 13
nr superpalindromes with 11 prime factors: 1
total nr superpalindromes: 6688
Currently, I'm mostly concerned about the efficiency of is_palindrome(mpz_t n) since it feels like a "two-pass" method: (a) copy the number to a string, (b) check if the string is a palindrome. But please highlight any other areas I might not be concerned about but should be.
Most of the numbers encountered by is_palindrome(mpz_t n) will not be palindromes, but my attempt at checking only the first and last digits was thwarted by mpz_sizeinbase not giving the exact number of digits in base 10 (and it added too much overhead to add a "check if the number of digits is correct" function).
