Computing a HAP basis

From Polymath Wiki
Revision as of 14:44, 14 January 2010 by Harrison (talk | contribs) (New page: This is some code I wrote to compute an alternative basis (the HAP basis, or basis of homogeneous arithmetic progressions) for a given sequence. It prints the list of integers n such that ...)
(diff) ← Older revision | Latest revision (diff) | Newer revision → (diff)
Jump to navigationJump to search
The printable version is no longer supported and may have rendering errors. Please update your browser bookmarks and please use the default browser print function instead.

This is some code I wrote to compute an alternative basis (the HAP basis, or basis of homogeneous arithmetic progressions) for a given sequence. It prints the list of integers n such that {n, 2n, 3n, ...} is a basis vector to stdout.

It hopefully isn't that hard to use, although it's not incredibly well-commented. Replace the 1124-sequence I gave in the declaration of "sequence" with your sequence, and change "#define SEQSIZE foo" to match the length of your sequence. Compile, and run. 

/* Computes a HAP basis for a given sequence
 * Author: Harrison Brown
 * Date: 1/14/09
 * I hereby release this work into the public domain */

// This is a hacky algorithm that could almost certainly be made faster, but right now the sequences are small enough that it doesn't really matter.

#include <stdio.h>
#include <math.h>
#define SEQSIZE 1125 // size of the input sequence plus 1

int main(void){
	const char sequence[SEQSIZE] = {0, +1, -1, +1, +1, -1, -1, -1, -1, +1, +1, -1, +1, +1, +1, -1, -1, +1, -1, +1, +1, -1, -1, -1, +1, -1, +1, +1, -1, +1, -1,
-1, +1, +1, -1, +1, -1, -1, +1, -1, -1, +1, +1, -1, +1, +1, -1, -1, -1, +1, +1, +1, -1, -1, +1, -1, +1, -1, -1, +1, +1,
+1, +1, -1, -1, +1, -1, -1, +1, +1, -1, -1, +1, +1, +1, -1, -1, +1, +1, -1, +1, -1, -1, +1, -1, -1, +1, +1, -1, +1, -1,
-1, +1, -1, +1, +1, +1, -1, -1, +1, -1, +1, -1, -1, +1, +1, -1, -1, -1, +1, +1, -1, +1, +1, +1, -1, -1, -1, +1, -1, -1,
+1, -1, +1, +1, +1, +1, -1, -1, -1, -1, +1, +1, -1, +1, +1, +1, -1, -1, +1, -1, -1, +1, -1, -1, +1, -1, +1, -1, +1, +1,
-1, +1, +1, +1, -1, -1, -1, -1, +1, +1, +1, +1, -1, -1, -1, +1, +1, -1, +1, -1, -1, +1, -1, +1, +1, -1, -1, +1, -1, +1,
+1, -1, +1, -1, +1, -1, +1, +1, -1, -1, -1, +1, +1, -1, +1, +1, -1, -1, +1, -1, -1, +1, -1, -1, +1, +1, +1, +1, -1, +1,
-1, -1, +1, -1, -1, +1, +1, -1, -1, +1, -1, +1, +1, -1, -1, -1, +1, -1, +1, +1, -1, +1, +1, +1, -1, -1, +1, +1, -1, -1,
+1, -1, -1, +1, -1, +1, +1, -1, -1, +1, -1, +1, +1, -1, +1, +1, -1, -1, +1, -1, -1, +1, -1, +1, +1, -1, -1, +1, +1, -1,
+1, -1, +1, -1, -1, +1, -1, +1, +1, +1, -1, -1, +1, -1, +1, +1, -1, -1, -1, -1, +1, +1, -1, -1, +1, -1, +1, +1, -1, +1,
+1, -1, -1, +1, -1, +1, +1, -1, -1, +1, +1, -1, +1, -1, -1, +1, -1, +1, -1, -1, +1, +1, -1, -1, +1, -1, +1, +1, +1, -1,
+1, -1, -1, +1, -1, +1, -1, -1, +1, +1, +1, +1, -1, -1, -1, +1, +1, -1, +1, -1, -1, +1, -1, +1, +1, -1, -1, +1, -1, +1,
+1, -1, +1, +1, -1, -1, +1, -1, -1, +1, -1, +1, +1, -1, -1, +1, +1, -1, -1, -1, +1, +1, -1, -1, +1, -1, +1, +1, -1, +1,
+1, -1, -1, +1, -1, -1, +1, +1, +1, +1, -1, -1, -1, -1, +1, +1, -1, +1, +1, -1, +1, +1, -1, -1, +1, -1, -1, +1, +1, -1,
-1, -1, +1, +1, -1, +1, +1, -1, -1, +1, -1, +1, +1, -1, +1, -1, -1, -1, +1, -1, +1, +1, -1, +1, +1, -1, -1, +1, -1, -1,
+1, -1, +1, +1, -1, -1, +1, +1, -1, +1, +1, +1, -1, -1, -1, +1, -1, +1, +1, -1, +1, +1, -1, -1, +1, -1, -1, -1, +1, +1,
+1, -1, -1, +1, -1, +1, -1, -1, +1, +1, +1, -1, +1, -1, +1, +1, -1, +1, -1, -1, -1, +1, +1, -1, +1, -1, +1, +1, -1, -1,
+1, -1, -1, +1, -1, +1, +1, -1, -1, +1, -1, +1, -1, +1, +1, -1, +1, -1, +1, -1, -1, +1, -1, +1, +1, -1, -1, +1, -1, -1,
+1, -1, +1, +1, +1, -1, -1, -1, +1, +1, -1, +1, +1, +1, -1, +1, -1, -1, +1, -1, +1, +1, -1, -1, +1, -1, +1, +1, -1, +1,
-1, -1, -1, +1, -1, +1, +1, -1, +1, +1, -1, -1, +1, -1, -1, -1, +1, +1, +1, -1, -1, +1, -1, +1, +1, -1, -1, +1, +1, -1,
+1, -1, +1, +1, -1, +1, -1, -1, -1, +1, -1, -1, +1, -1, +1, +1, -1, -1, +1, -1, +1, +1, -1, +1, +1, -1, -1, +1, -1, +1,
+1, -1, +1, +1, -1, -1, +1, -1, -1, +1, -1, +1, +1, -1, -1, +1, -1, -1, +1, -1, +1, +1, -1, +1, -1, -1, +1, +1, -1, +1,
+1, -1, -1, +1, -1, -1, +1, -1, +1, +1, -1, -1, +1, -1, +1, +1, -1, +1, +1, -1, -1, +1, -1, +1, +1, -1, -1, +1, +1, -1,
+1, -1, -1, +1, -1, +1, +1, -1, -1, +1, -1, -1, +1, -1, +1, +1, -1, -1, +1, -1, +1, +1, -1, +1, +1, -1, +1, +1, -1, -1,
+1, -1, +1, +1, -1, -1, +1, -1, -1, +1, -1, +1, +1, -1, -1, +1, -1, +1, -1, -1, +1, +1, +1, -1, +1, -1, -1, +1, -1, +1,
+1, -1, -1, +1, -1, +1, -1, -1, +1, +1, +1, -1, -1, -1, +1, +1, -1, +1, +1, -1, -1, +1, +1, -1, +1, -1, +1, +1, -1, -1,
+1, -1, -1, +1, -1, -1, +1, -1, +1, +1, -1, +1, -1, -1, +1, +1, +1, -1, +1, -1, -1, +1, -1, +1, +1, -1, -1, +1, -1, +1,
-1, -1, +1, +1, -1, -1, +1, -1, +1, +1, -1, +1, +1, -1, -1, +1, +1, -1, +1, -1, -1, +1, -1, -1, +1, -1, +1, +1, -1, +1,
+1, -1, -1, +1, -1, +1, +1, -1, +1, +1, -1, -1, +1, -1, -1, +1, +1, +1, -1, -1, -1, +1, -1, -1, +1, -1, +1, +1, -1, -1,
+1, +1, +1, +1, -1, +1, -1, -1, +1, +1, -1, -1, +1, -1, +1, +1, -1, -1, +1, -1, -1, +1, -1, +1, +1, -1, -1, +1, -1, +1,
+1, -1, +1, +1, -1, -1, +1, -1, -1, +1, -1, +1, +1, -1, -1, -1, +1, +1, -1, -1, +1, +1, -1, -1, +1, -1, +1, +1, -1, +1,
+1, -1, -1, +1, -1, -1, +1, -1, +1, +1, +1, -1, +1, -1, -1, +1, -1, +1, +1, -1, -1, +1, -1, +1, +1, -1, +1, -1, +1, -1,
-1, +1, -1, +1, -1, +1, +1, -1, -1, +1, -1, +1, -1, -1, +1, +1, +1, -1, +1, -1, -1, -1, +1, +1, +1, -1, -1, +1, -1, -1,
+1, -1, +1, +1, -1, -1, +1, +1, +1, +1, -1, +1, -1, -1, -1, +1, -1, +1, +1, -1, +1, +1, -1, -1, +1, -1, -1, +1, +1, +1,
-1, -1, -1, +1, -1, -1, +1, -1, +1, +1, +1, -1, +1, -1, +1, +1, -1, +1, -1, -1, -1, +1, -1, -1, +1, +1, +1, +1, -1, -1,
+1, -1, +1, +1, -1, -1, +1, -1, +1, +1, -1, +1, -1, -1, +1, -1, +1, +1, -1, -1, -1, +1, -1, -1, +1, +1, -1, +1, -1, +1,
+1, -1, +1, -1, -1, -1, +1, -1, +1, +1, +1, +1, -1, +1, -1, -1, -1, -1, +1, -1, +1, +1, -1, +1, +1, -1, -1, +1, +1, -1,
+1, +1, -1, +1, -1, -1, +1, -1, +1, +1, -1, +1, -1, -1};
	char basis[SEQSIZE];
	int i, j;
	int fb;		// "factor bound"
	int resetFB;	// will act as a boolean
	char basisVector, initialize;
	basis[0] = 0;
	initialize = sequence[1];
	basis[1] = initialize;
	fb = 1;		// fb will be the floor of sqrt(n) throughout
	resetFB = 0;
	for(i = 2; i < SEQSIZE; i++){
		basisVector = initialize * sequence[i];
		if(resetFB){
			fb++;
			resetFB = 0;
		}		
		if((fb+1) * (fb+1) == i){
			basisVector *= basis[fb+1];	//we don't want to double count square roots, so handle them separately
			resetFB = 1;
		}
		for(j = 2; j <= fb; j++){
			if(i % j == 0) basisVector *= (basis[j] * basis[i/j]);
		}
		basis[i] = basisVector;
		if(basisVector == -1){
			printf("%d ", i);
		}
	}
	printf("\n");
	return 0;
}
Retrieved from "https://michaelnielsen.org/polymath/index.php?title=Computing_a_HAP_basis&oldid=2725"