engineer_eden.c

/*
 ============================================================================
 Name        : engineer_edem.c
 Author      : 
 Version     :
 Copyright   : Yuriy Nazarov
 Compile     : gcc -std=c11 -O3 -fopenmp enginer_edem.c
 ============================================================================
 */

#include <omp.h>
#include <stdio.h>
#include <stdlib.h>
#include <math.h>
#include <string.h>
#include <stdint.h>

int check(uint32_t n){
	uint32_t s = 0;
	uint32_t i;
	for(i = 1; i*i <= n; i++){
		if (i-1 > s) return 0;
		if (0 == n % i) s += i;
	}
	i--;
	if(i*i)i--;
	for(; i >= 1; i--){
		if (0 == n % i){
			uint32_t j = n / i;
			if (j-1 > s) return 0;
			s += j;
		}
	}
	if (n-1 > s) return 0;
	return 1;
}

int preCheck5(uint32_t n){
	uint32_t i = 1;
	uint32_t s0 = 0, s1 = 0, s2 = 0, s3 = 0, s4 = 0;
	for(; i*i <= n; i++){
		if (i-1 > s0 || i-1 > s1 || i-1 > s2 || i-1 > s3 || i-1 > s4) return 0;
		if (0 == (n+0) % i) s0 = (s0>0xFFffFFff-i) ? 0xFFffFFff : s0+i;
		if (0 == (n+4) % i) s1 = (s0>0xFFffFFff-i) ? 0xFFffFFff : s1+i;
		if (0 == (n+8) % i) s2 = (s0>0xFFffFFff-i) ? 0xFFffFFff : s2+i;
		if (0 == (n+12) % i) s3 = (s0>0xFFffFFff-i) ? 0xFFffFFff : s3+i;
		if (0 == (n+16) % i) s4 = (s0>0xFFffFFff-i) ? 0xFFffFFff : s4+i;
	}
	return 1;
}

void eratosthenesOddSingleBlock(char* isComposite, uint32_t from, uint32_t to) {
	for (int i = 3; i*i <= to; i+=2) {
		if (i >= 3*3 && i % 3 == 0) continue;
		if (i >= 5*5 && i % 5 == 0) continue;
		if (i >= 7*7 && i % 7 == 0) continue;
		if (i >= 11*11 && i % 11 == 0) continue;
		if (i >= 13*13 && i % 13 == 0) continue;
		if (i >= 17*17 && i % 17 == 0) continue;
		if (i >= 19*19 && i % 19 == 0) continue;
		if (i >= 23*23 && i % 23 == 0) continue;
		int minJ = ((from+i-1)/i)*i;
		if (minJ < i*i) minJ = i*i;
		if ((minJ & 1) == 0) minJ += i;
		for (int j = minJ; j <= to; j += 2*i) {
			int index = j - from;
			isComposite[index/2] = 1;
		}
	}
}

#define prime(n) !isComposite[(n)/2]
void checkQuads(char *isComposite, uint32_t from, uint32_t to){
	for (uint64_t m = 1; from+m < to-18; m+=2){
		if (prime(m)    && !prime(m+2)  && !prime(m+4) &&
			prime(m+6)  && !prime(m+8)  && !prime(m+10) &&
			prime(m+12) && !prime(m+14) && !prime(m+16) &&
			prime(m+18) ){
			if(preCheck5(from+m+1) && check(from+m+1) && check(from+m+5) && check(from+m+9) && check(from+m+13) && check(from+m+17))
				printf("%ld %ld %ld %ld\n", from+m, from+m+6, from+m+12, from+m+18);
		}
	}
}

#define N 1200
#define M 1000000
int main(void) {
	#pragma omp parallel for schedule(static, 100)
	for(int i = 0; i < N; i++){
		char isComposite[M/2+10];
		memset(isComposite, 0, M/2+10);
		eratosthenesOddSingleBlock(isComposite, M*(i+0), M*(i+1)-1+20);
		checkQuads(isComposite, M*(i+0), M*(i+1)-1+20);
	}
}