// References:
// - saffah_codex_6s_agg2, https://duck.ac/submission/105858: copied our accepted sorter and preserved its inherited citations below. Its public source has no separate license notice.
// Approach: On the initial histogram pass, one eight-byte load and BMI2 extraction yield the first five characters of four adjacent suffixes. Increment their same two-letter radix buckets from four shifted packed keys, while leaving the scatter and all later stages unchanged.
// Purpose: Measure whether shared character extraction reduces the official histogram time enough to improve the judged result.
// References:
// [1] saffah_codex_6s_agg2, https://duck.ac/submission/105598 . Directly reused our accepted suffix sorter, four-key MSD scatter, and locality-3 output-position prefetch.
// [2] saffah_cc_v41_agg1, https://duck.ac/submission/105404 . Its earlier sorter and ancestor references are inherited through [1] and retained below.
// The public sources contain no explicit license notice. Exact authors, addresses, and reused code are credited here and in the inherited comments before code.
// Approach:
// Prefetch the future radix output position 104 rather than 112 input positions ahead. This alters the overlap between the bucket-count dependency and the write-allocate request without altering any record or output byte.
// Purpose:
// Measure the nearby prefetch distance on official 1006e6 judging after the 112-position variant remained above the current strict target.
// References:
// [1] saffah_cc_v41_agg1, duck.ac #105404, https://duck.ac/submission/105404.
// Copied its accepted suffix sorter and inherited its citations. No
// separate license notice was attached to its public source.
// [2] saffah_codex_6s_agg2, duck.ac #105445, https://duck.ac/submission/105445.
// Directly reused our rolling-key, split-scatter and LCP-guard base.
// Approach:
// Prefetch the future radix output position 112 input positions ahead with
// locality hint 3; this changes the write-allocate timing and cache level.
// Purpose:
// Measure whether this distinct memory-scheduling choice meets the judge target.
#pragma GCC optimize("O3","unroll-loops","rename-registers")
#pragma GCC target("arch=skylake,avx2,bmi,bmi2,popcnt,lzcnt")
#define JUDGE 1
// 1006e6 v3 (lane e_variant_lane): SA + LCP for a uniformly-random lowercase string.
// stage 1 build u64 records [key30(c0..c5, 5 bits each) | idx20] (BMI2 pext, no gathers)
// stage 2 3 LSD radix passes, 10 bits each, 1024 buckets (L1 counters)
// stage 3 one walk: emit SA digits, derive LCP from packed-key xor (clz), refine equal-key runs
// stage 4 emit the LCP line (byte array -> digits)
// Fallback to SA-IS when a 6-char-key run exceeds RUNMAX (non-random data).
#include <cstdint>
#include <cstring>
#include <immintrin.h>
#include <immintrin.h>
#define MAXN 1000005
#define PADB 64
static unsigned char str[MAXN + PADB];
__attribute__((aligned(64))) static uint64_t RA[MAXN], RB[MAXN];
static uint32_t CNT[1024];
/* second radix pass's histogram, tallied during the FIRST pass's scatter (copied from
the rival's #102685). A histogram over a multiset does not care which permutation
of it you scan, so counting over `cur[]` is exactly counting over `nx[]`. */
static uint32_t CNT2[1024];
static uint32_t START[1025];
#define SCRSZ 65536
static uint64_t SCR[SCRSZ];
static unsigned char LCPB[MAXN + 8];
static uint32_t D4[10000];
static unsigned char L10[64];
/* LANE l1006e6b, copied from #100017: `(clz(x)-2)/5` is a multiply-high chain ON the
critical path of the four LCPB byte stores in the walk; clz of a non-zero 30-bit key
is 2..31, so the quotient is a 32-entry table (and clz(0) = 32 maps to 6, which is what
the old `(__builtin_clz(0)-2)/5` produced too). */
static const unsigned char D5B[64] = {0,0,0,0,0,0,0,1,1,1,1,1,2,2,2,2,2,3,3,3,3,3,4,4,4,4,4,5,5,5,5,5,6,6,6,6,6,6,6,6,6,6,6,6,6,6,6,6,6,6,6,6,6,6,6,6,6,6,6,6,6,6,6,6};
static unsigned char CLA[64];
static uint64_t PAIR[256];
static const unsigned PW[12] = {0, 1, 10, 100, 1000, 10000, 100000, 1000000,
10000000, 100000000, 1000000000u, 0};
/* AVX2 exclusive prefix sum over the 1024 radix counters. Directly copied from the
rival's #99905, whose comment reads: the scalar form was a serial `s += c` chain of
1024 dependent adds plus 1024 single-word stores, i.e. ~1024 cycles per call, and it
is called twice per non-empty bucket (~676 of them at n=1e6). `acc - original`
recovers the exclusive prefix, so the only shuffles left are the two vpermd broadcasts
of lane 3 and lane 7. */
static inline void prefix_excl864(uint32_t *C) {
__m256i carry = _mm256_setzero_si256();
const __m256i msk = _mm256_setr_epi32(0, 0, 0, 0, -1, -1, -1, -1);
const __m256i id3 = _mm256_set1_epi32(3), id7 = _mm256_set1_epi32(7);
for (int d = 0; d < 864; d += 8) {
__m256i o = _mm256_loadu_si256((const __m256i *)(C + d));
__m256i v = o;
v = _mm256_add_epi32(v, _mm256_slli_si256(v, 4));
v = _mm256_add_epi32(v, _mm256_slli_si256(v, 8));
v = _mm256_add_epi32(v, _mm256_and_si256(_mm256_permutevar8x32_epi32(v, id3), msk));
__m256i acc = _mm256_add_epi32(v, carry);
_mm256_storeu_si256((__m256i *)(C + d), _mm256_sub_epi32(acc, o));
carry = _mm256_permutevar8x32_epi32(acc, id7);
}
}
static inline uint64_t pext64(uint64_t v, uint64_t m) {
uint64_t r; __asm__("pext %2,%1,%0" : "=r"(r) : "r"(v), "r"(m)); return r;
}
#define M5X8 0x1F1F1F1F1F1F1F1FULL
#define REC_KEY(r) ((uint32_t)((r) >> 20))
#define REC_IDX(r) ((uint32_t)((r) & 0xFFFFFu))
static inline uint32_t key_at(int i) {
uint64_t v; memcpy(&v, str + i, 8);
return (uint32_t)(pext64(__builtin_bswap64(v), M5X8) >> 10);
}
/* [x6 cbd] the MSD bucket of suffix i is exactly bits 20..29 of key_at(i), i.e.
c1 | c0<<5 with cj = str[i+j] & 31 -- obtainable from a single two-byte load with no
bswap and no pext. Used by the COUNT pass only; the scatter still needs the full
30-bit key for the record it writes. */
static inline uint32_t bk_dig(int i) {
uint16_t v; memcpy(&v, str + i, 2);
return ((uint32_t)(v >> 8) & 31u) | (((uint32_t)v & 31u) << 5);
}
static inline int len10c(uint32_t v) { return CLA[31 - __builtin_clz(v | 1)]; }
static inline int len10(uint32_t v) {
int b = 31 - __builtin_clz(v | 1);
int l = L10[b];
l -= (v < PW[l]);
return l;
}
// pack v (<= 9999999) as digits followed by a space in an 8-byte word (left-aligned, low byte first)
/* LANE pragma_lane: TWO-LEVEL CHUNK WRITER. The old form built the number from a 4+4 digit split
and then LENGTH-ARITHMETIC'd it into place: 2 magic divisions, 6 compares, a cmov, and TWO variable
shifts (`>>= (8-l)*8` and `' ' << (8*l)`). This form looks the digits up in two 8 KB L1-resident
tables of PRE-SHIFTED, pre-spaced 8-byte words (`T4a[A]` = A's digits + space + its length in byte
7; `T3b[B]` = B as 3 zero-padded digits + space) and places the low part with ONE shift by
8*len(A). Byte-identical: the same digits in the same bytes, and the bytes ABOVE the space are
always overwritten by the next store (the last number of the line still goes through
`putu_exact`). MEASURED in a faithful mikro: -19.3 % instructions at unchanged IPC (3.29 -> 3.24),
i.e. -18.1 % cycles, on a loop that IS throughput-bound. */
static uint64_t T3x[1024], T4a[1024], T4d[1024]; static uint32_t T3b[1024];
static unsigned char LENA[1024];
static void build_t4(void) {
for (int v = 0; v < 1024; v++) {
char t[8] = {0,0,0,0,0,0,0,0}; int l;
if (v < 10) { t[0] = (char)('0' + v); t[1] = ' '; l = 2; }
else if (v < 100) { t[0] = (char)('0' + v / 10); t[1] = (char)('0' + v % 10); t[2] = ' '; l = 3; }
else { t[0] = (char)('0' + v / 100); t[1] = (char)('0' + (v / 10) % 10); t[2] = (char)('0' + v % 10); t[3] = ' '; l = 4; }
uint64_t u = 0; memcpy(&u, t, 8);
T3x[v] = u | ((uint64_t)(l - 1) << 56); /* for the A == 0 case */
char z[8] = {(char)('0' + (v / 100) % 10), (char)('0' + (v / 10) % 10), (char)('0' + v % 10), ' ', 0,0,0,0};
uint64_t uz = 0; memcpy(&uz, z, 8); T3b[v] = (uint32_t)uz; /* 3 zero-padded digits + space */
char w[8] = {0,0,0,0,0,0,0,0}; int la;
if (v < 10) { w[0] = (char)('0' + v); w[1] = ' '; la = 1; }
else if (v < 100) { w[0] = (char)('0' + v / 10); w[1] = (char)('0' + v % 10); w[2] = ' '; la = 2; }
else if (v < 1000) { w[0] = (char)('0' + v / 100); w[1] = (char)('0' + (v / 10) % 10); w[2] = (char)('0' + v % 10); w[3] = ' '; la = 3; }
else { w[0] = '1'; w[1] = '0'; w[2] = '0'; w[3] = '0'; w[4] = ' '; la = 4; }
uint64_t uw = 0; memcpy(&uw, w, 8); T4a[v] = uw | ((uint64_t)la << 56); LENA[v] = (unsigned char)la; T4d[v] = uw;
}
}
static inline uint64_t pack_num(uint32_t v, unsigned &adv) {
/* [k6b magic] `v / 1000` written out as its own unsigned-division magic
(274877907 / 2^38 == 1/1000, exact for every uint32_t). With the source
written as `v / 1000`, gcc 9.3 emits the quotient TWICE (`mov rcx,r13;
shr rcx,32; shr rcx,6` next to `shr r13,38`) because it needs it as both a
32-bit and a 64-bit value; writing the magic by hand gives it one value
number and one shift. Judge-visible effect: 27 fewer instructions in the
4-wide walk loop per 4 elements (148 -> 121 at gcc 9.3 -O2). */
uint64_t q = (uint64_t)v * 274877907ULL;
uint64_t A = q >> 38;
uint32_t B = (uint32_t)(v - (uint32_t)A * 1000);
uint64_t ta = T4d[A], tb = (uint64_t)T3b[B];
if (!A) { uint64_t tn = T3x[B]; adv = (unsigned)((tn >> 56) + 1); return tn; }
unsigned la = LENA[A];
uint64_t u = ta | (tb << (la << 3));
adv = la + 4;
return u;
}
/* [e6 c5] `pack_num` + store in one call, so the packed word has to live for exactly
one table lookup instead of four. The 4-wide caller used to keep u0..u3 alive across
all four lookups and gcc spilled every one of them (`movq %rcx,128(%rsp)` ...
`movl 136(%rsp),%ebx` ... `movq %rax,(%o)`), i.e. three memory ops per element where
one is needed. Byte-identical output: the same word, stored at the same `o`, and `o`
advances by the same total, so the overlapping-store trick is unchanged. */
static inline char *pack_emit(char *o, uint32_t v) {
unsigned adv;
uint64_t u = pack_num(v, adv);
memcpy(o, &u, 8);
return o + adv;
}
// write v (<= 9999999) followed by a space; 8-byte store, advances by len+1
// byte-exact writer (no overrun) for the final element of each line
static inline char *putu_exact(char *p, unsigned v) {
char tmp[8];
int i = 8;
while (v >= 100) { unsigned q = v / 100, r = v - q * 100; tmp[--i] = '0' + r % 10; tmp[--i] = '0' + r / 10; v = q; }
if (v >= 10) { tmp[--i] = (char)('0' + v % 10); tmp[--i] = (char)('0' + v / 10); }
else tmp[--i] = (char)('0' + v);
memcpy(p, tmp + i, 8 - i);
return p + (8 - i);
}
static inline void build_D4_if_needed() {
static int ready = 0;
if (ready) return;
for (int i = 0; i < 10000; i++) {
char t[4];
t[0] = (char)('0' + i / 1000); t[1] = (char)('0' + (i / 100) % 10);
t[2] = (char)('0' + (i / 10) % 10); t[3] = (char)('0' + i % 10);
memcpy(&D4[i], t, 4);
}
ready = 1;
}
static inline char *wr(char *p, uint32_t v) {
build_D4_if_needed();
uint32_t hi = v / 10000, lo = v % 10000;
uint64_t u = (uint64_t)D4[hi] | ((uint64_t)D4[lo] << 32);
int l = len10(v);
u >>= (8 - l) * 8;
u |= (uint64_t)' ' << (8 * l);
memcpy(p, &u, 8);
return p + l + 1;
}
// small-value writer for the LCP line (v < 1000 typical)
static uint64_t T3[1024];
static void build_T3() {
for (int v = 0; v < 1024; v++) {
char t[8]; int l;
if (v < 10) { t[0] = (char)('0' + v); t[1] = ' '; l = 2; }
else if (v < 100) { t[0] = (char)('0' + v / 10); t[1] = (char)('0' + v % 10); t[2] = ' '; l = 3; }
else { t[0] = (char)('0' + v / 100); t[1] = (char)('0' + (v / 10) % 10); t[2] = (char)('0' + v % 10); t[3] = ' '; l = 4; }
uint64_t u = 0; memcpy(&u, t, 8);
T3[v] = u; // left-aligned: low byte = first char
}
}
static inline char *wrs(char *p, uint32_t v) {
if (v >= 1024) return wr(p, v); // rare: long LCP (non-random data)
uint64_t u = T3[v];
memcpy(p, &u, 8);
return p + ((v < 10) ? 2 : (v < 100) ? 3 : 4);
}
// ---------------- SA-IS fallback (unchanged behaviour) ----------------
#define TPOOL 24
static int g_bkt[MAXN + 2];
static int g_cntc[MAXN + 2];
static unsigned char t_pool[TPOOL][(MAXN + 7) / 8];
static inline int tget(const unsigned char *t, int i) { return (t[i >> 3] >> (i & 7)) & 1; }
static inline void tset(unsigned char *t, int i, int b) {
if (b) t[i >> 3] |= (unsigned char)(1 << (i & 7));
else t[i >> 3] &= (unsigned char)~(1 << (i & 7));
}
static inline int isLMS(const unsigned char *t, int i) { return i > 0 && tget(t, i) && !tget(t, i - 1); }
template <typename T> static void getCounts(const T *s, int *C, int n, int K) {
for (int i = 0; i <= K; i++) C[i] = 0;
for (int i = 0; i < n; i++) C[(int)s[i]]++;
}
static void getBuckets(const int *C, int *B, int K, bool end) {
int sum = 0;
if (end) { for (int i = 0; i <= K; i++) { sum += C[i]; B[i] = sum; } }
else { for (int i = 0; i <= K; i++) { B[i] = sum; sum += C[i]; } }
}
template <typename T> static void induceSAl(const unsigned char *t, int *SA, const T *s, int *bkt, int n, int K) {
getBuckets(g_cntc, bkt, K, false);
for (int i = 0; i < n; i++) { int j = SA[i] - 1; if (j >= 0 && !tget(t, j)) SA[bkt[(int)s[j]]++] = j; }
}
template <typename T> static void induceSAs(const unsigned char *t, int *SA, const T *s, int *bkt, int n, int K) {
getBuckets(g_cntc, bkt, K, true);
for (int i = n - 1; i >= 0; i--) { int j = SA[i] - 1; if (j >= 0 && tget(t, j)) SA[--bkt[(int)s[j]]] = j; }
}
template <typename T> static void sais(const T *s, int *SA, int n, int K, int depth) {
unsigned char *t = t_pool[depth];
int *bkt = g_bkt;
int i, j;
tset(t, n - 2, 0); tset(t, n - 1, 1);
for (i = n - 3; i >= 0; i--) tset(t, i, (s[i] < s[i + 1] || (s[i] == s[i + 1] && tget(t, i + 1))));
getCounts(s, g_cntc, n, K);
getBuckets(g_cntc, bkt, K, true);
for (i = 0; i < n; i++) SA[i] = -1;
for (i = 1; i < n; i++) if (isLMS(t, i)) SA[--bkt[(int)s[i]]] = i;
induceSAl(t, SA, s, bkt, n, K);
induceSAs(t, SA, s, bkt, n, K);
int n1 = 0;
for (i = 0; i < n; i++) if (isLMS(t, SA[i])) SA[n1++] = SA[i];
for (i = n1; i < n; i++) SA[i] = -1;
int name = 0, prev = -1;
for (i = 0; i < n1; i++) {
int pos = SA[i];
bool diff = false;
int dmax = n - (pos > prev ? pos : prev);
for (int d = 0; d < dmax; d++)
if (prev == -1 || s[pos + d] != s[prev + d] || tget(t, pos + d) != tget(t, prev + d)) { diff = true; break; }
else if (d > 0 && (isLMS(t, pos + d) || isLMS(t, prev + d))) break;
if (diff) { name++; prev = pos; }
pos = (pos % 2 == 0) ? pos / 2 : (pos - 1) / 2;
SA[n1 + pos] = name - 1;
}
for (i = n - 1, j = n - 1; i >= n1; i--) if (SA[i] >= 0) SA[j--] = SA[i];
int *SA1 = SA, *s1 = SA + n - n1;
if (name < n1) sais<int>(s1, SA1, n1, name - 1, depth + 1);
else for (i = 0; i < n1; i++) SA1[s1[i]] = i;
getCounts(s, g_cntc, n, K);
getBuckets(g_cntc, bkt, K, true);
for (i = 1, j = 0; i < n; i++) if (isLMS(t, i)) s1[j++] = i;
for (i = 0; i < n1; i++) SA1[i] = s1[SA1[i]];
for (i = n1; i < n; i++) SA[i] = -1;
for (i = n1 - 1; i >= 0; i--) { j = SA[i]; SA[i] = -1; SA[--bkt[(int)s[j]]] = j; }
induceSAl(t, SA, s, bkt, n, K);
induceSAs(t, SA, s, bkt, n, K);
}
static int g_sa[MAXN];
// ---------------- helpers: full suffix compare / lcp ----------------
static inline bool suf_less(uint32_t i, uint32_t j) {
const unsigned char *p = str + i, *q = str + j;
for (;;) {
uint64_t a; memcpy(&a, p, 8); a = __builtin_bswap64(a);
uint64_t b; memcpy(&b, q, 8); b = __builtin_bswap64(b);
if (a != b) return a < b;
p += 8; q += 8;
}
}
static inline int lcp_len(uint32_t i, uint32_t j) {
const unsigned char *p = str + i, *q = str + j;
int h = 0;
for (;;) {
uint64_t a; memcpy(&a, p, 8); a = __builtin_bswap64(a);
uint64_t b; memcpy(&b, q, 8); b = __builtin_bswap64(b);
if (a != b) return h + ((int)(__builtin_clzll(a ^ b) >> 3));
h += 8; p += 8; q += 8;
}
}
#define RUNMAX 48
static int g_fallback = 0;
static char *g_sa_out_start;
static long g_osal;
#ifdef PHASE
static unsigned long long P[16];
#define MARK(k) P[k] = __builtin_ia32_rdtsc()
#else
#define MARK(k)
#endif
static int solve_string(const unsigned char *q, long sn, char **opo) {
int n = 0;
while (n < sn && (q[n] == ' ' || q[n] == '\n' || q[n] == '\r' || q[n] == '\t')) { n++; }
{ /* [x6 parse] was: byte-at-a-time store with four comparisons per byte. Locate the
token end with 32-byte SIMD compares, then copy the whole token in one memcpy.
Same token: the loop stops at the first byte in {0, ' ', '\n', '\r', '\t'}. */
long st = n, p = n;
int done = 0;
const __m256i sp8 = _mm256_set1_epi8(' '), nl8 = _mm256_set1_epi8('\n');
const __m256i cr8 = _mm256_set1_epi8('\r'), tb8 = _mm256_set1_epi8('\t');
const __m256i zz8 = _mm256_setzero_si256();
while (p + 32 <= sn) {
__m256i v; memcpy(&v, q + p, 32);
__m256i hit = _mm256_or_si256(
_mm256_or_si256(_mm256_cmpeq_epi8(v, sp8), _mm256_cmpeq_epi8(v, nl8)),
_mm256_or_si256(_mm256_cmpeq_epi8(v, cr8),
_mm256_or_si256(_mm256_cmpeq_epi8(v, tb8), _mm256_cmpeq_epi8(v, zz8))));
unsigned msk = (unsigned)_mm256_movemask_epi8(hit);
if (msk) { p += (long)__builtin_ctz(msk); done = 1; break; }
p += 32;
}
if (!done)
while (p < sn && q[p] != 0 && q[p] != ' ' && q[p] != '\n' && q[p] != '\r' && q[p] != '\t') p++;
n = (int)(p - st);
if (n > 0) memcpy(str, q + st, (size_t)n);
}
if (n <= 0) return 0;
memset(str + n, 0, PADB);
// L10[b] = decimal digits of the largest value with bit-length b+1 (= digits(2^(b+1)-1))
for (int b = 0; b < 32; b++) {
unsigned long long hi = (b == 31) ? 0xFFFFFFFFull : ((1ull << (b + 1)) - 1);
int l = 1; unsigned long long pp = 1;
while (pp * 10 <= hi) { pp *= 10; l++; }
L10[b] = (unsigned char)l;
}
for (int d = 0; d < 256; d++) {
uint64_t u = (d < 10) ? ((uint64_t)('0' + d) | ((uint64_t)' ' << 8)) : 0;
PAIR[d] = u;
}
build_T3();
build_t4();
MARK(0);
// ---- build records + MSD count (fused), then MSD scatter by chars c0,c1 ----
memset(CNT, 0, sizeof(CNT));
int ci = 0;
for (; ci + 4 <= n; ci += 4) {
uint64_t v; memcpy(&v, str + ci, 8);
uint64_t p = pext64(__builtin_bswap64(v), M5X8);
CNT[(p >> 30) & 1023]++;
CNT[(p >> 25) & 1023]++;
CNT[(p >> 20) & 1023]++;
CNT[(p >> 15) & 1023]++;
}
for (; ci < n; ci++) CNT[bk_dig(ci)]++;
{
uint32_t sum = 0;
for (int d = 0; d < 1024; d++) { uint32_t c = CNT[d]; START[d] = sum; CNT[d] = sum; sum += c; }
START[1024] = sum;
}
MARK(1);
uint32_t k0 = key_at(0), k1 = key_at(1), k2 = key_at(2), k3 = key_at(3);
int i = 0;
for (; i + 108 <= n; i += 4) {
const uint32_t future_bucket = bk_dig(i + 104);
__builtin_prefetch(RB + CNT[future_bucket], 1, 3);
// Four independent rolling keys retain the exact c0..c5 key for suffixes i..i+3.
RB[CNT[k0 >> 20]++] = ((uint64_t)k0 << 20) | (uint64_t)(uint32_t)(i + 1);
RB[CNT[k1 >> 20]++] = ((uint64_t)k1 << 20) | (uint64_t)(uint32_t)(i + 2);
RB[CNT[k2 >> 20]++] = ((uint64_t)k2 << 20) | (uint64_t)(uint32_t)(i + 3);
RB[CNT[k3 >> 20]++] = ((uint64_t)k3 << 20) | (uint64_t)(uint32_t)(i + 4);
uint64_t v; memcpy(&v, str + i + 6, 8);
uint64_t p = pext64(__builtin_bswap64(v), M5X8);
k0 = ((k0 << 20) & 0x3FFFFFFFu) | (uint32_t)((p >> 20) & 0xFFFFFu);
k1 = ((k1 << 20) & 0x3FFFFFFFu) | (uint32_t)((p >> 15) & 0xFFFFFu);
k2 = ((k2 << 20) & 0x3FFFFFFFu) | (uint32_t)((p >> 10) & 0xFFFFFu);
k3 = ((k3 << 20) & 0x3FFFFFFFu) | (uint32_t)((p >> 5) & 0xFFFFFu);
}
for (; i + 4 <= n; i += 4) {
// Four independent rolling keys retain the exact c0..c5 key for suffixes i..i+3.
RB[CNT[k0 >> 20]++] = ((uint64_t)k0 << 20) | (uint64_t)(uint32_t)(i + 1);
RB[CNT[k1 >> 20]++] = ((uint64_t)k1 << 20) | (uint64_t)(uint32_t)(i + 2);
RB[CNT[k2 >> 20]++] = ((uint64_t)k2 << 20) | (uint64_t)(uint32_t)(i + 3);
RB[CNT[k3 >> 20]++] = ((uint64_t)k3 << 20) | (uint64_t)(uint32_t)(i + 4);
uint64_t v; memcpy(&v, str + i + 6, 8);
uint64_t p = pext64(__builtin_bswap64(v), M5X8);
k0 = ((k0 << 20) & 0x3FFFFFFFu) | (uint32_t)((p >> 20) & 0xFFFFFu);
k1 = ((k1 << 20) & 0x3FFFFFFFu) | (uint32_t)((p >> 15) & 0xFFFFFu);
k2 = ((k2 << 20) & 0x3FFFFFFFu) | (uint32_t)((p >> 10) & 0xFFFFFu);
k3 = ((k3 << 20) & 0x3FFFFFFFu) | (uint32_t)((p >> 5) & 0xFFFFFu);
}
for (; i < n; i++) {
uint32_t k = key_at(i);
RB[CNT[k >> 20]++] = ((uint64_t)k << 20) | (uint64_t)(uint32_t)(i + 1);
}
MARK(2);
{ // oversized bucket => non-random data: bail out to SA-IS
int bad = 0;
for (int d = 0; d < 1024; d++) if (START[d + 1] - START[d] > SCRSZ) { bad = 1; break; }
if (bad) {
g_fallback = 1;
sais<unsigned char>(str, g_sa, n + 1, 256, 0);
for (int x = 0; x < n; x++) g_sa[x] = g_sa[x + 1];
char *o = *opo;
for (int x = 0; x < n - 1; x++) o = wr(o, (uint32_t)g_sa[x] + 1);
o = putu_exact(o, (uint32_t)g_sa[n - 1] + 1);
*o++ = '\n';
if (n > 1) {
for (int x = 1; x < n - 1; x++) o = wrs(o, (uint32_t)lcp_len((uint32_t)g_sa[x - 1], (uint32_t)g_sa[x]));
o = putu_exact(o, (uint32_t)lcp_len((uint32_t)g_sa[n - 2], (uint32_t)g_sa[n - 1]));
}
*o++ = '\n';
*opo = o;
return 1;
}
}
MARK(3);
// ---- per-bucket LSD by (c4,c5) then (c2,c3): fully cache resident ----
for (int bk = 0; bk < 864; bk++) {
uint32_t lo = START[bk], hi = START[bk + 1];
uint32_t m = hi - lo;
if (m < 2) continue;
uint64_t *cur = RB + lo, *nx = SCR;
{
int sh = 20;
memset(CNT, 0, sizeof(uint32_t) * 864);
for (uint32_t i = 0; i < m; i++) CNT[(uint32_t)(cur[i] >> sh) & 1023]++;
prefix_excl864(CNT);
memset(CNT2, 0, sizeof(uint32_t) * 864);
for (uint32_t i = 0; i < m; i++) {
uint64_t record = cur[i];
nx[CNT[(uint32_t)(record >> sh) & 1023]++] = record;
CNT2[(uint32_t)(record >> 30) & 1023]++;
}
}
{
prefix_excl864(CNT2);
for (uint32_t i = 0; i < m; i++) {
uint64_t record = nx[i];
cur[CNT2[(uint32_t)(record >> 30) & 1023]++] = record;
}
}
}
MARK(4);
// ---- single fused walk: run refinement + lcp + emit SA line + LCP line ----
char *o = *opo;
{
const char *ostart = o;
/* [k6b nox] `kprev` starts at ~0u (no 30-bit key can equal it), so the two
`if (x)` guards on the fast path are dead: at x == 0 the carry test
`k0 == kprev` is false and `LCPB[0]` gets `D5B[clz(~0u ^ k0)] == D5B[0]
== 0`, which is exactly the 0 the line above writes. LCPB[0] is never
read (the LCP line is emitted from index 1). */
uint32_t kprev = 0xFFFFFFFFu;
int x = 0;
LCPB[0] = 0;
/* LANE pragma_lane: THE FAST PATH USED TO EXIT ON THE FIRST REPEATED 6-GRAM AND NEVER RESUME.
`if (k1==k0 || ...) break;` fires as soon as two ADJACENT sorted suffixes share a 30-bit
key, i.e. the first repeated 6-gram -- and a random 1e6 string contains ~1600 of them
(n^2/2/26^6), the first at x=732. So the 4-at-a-time path processed 732 of 1e6 elements
and the per-element run loop did the other 999 268. Restructure: the fast loop is now the
OUTER loop and the run loop handles EXACTLY ONE equal-key run before control returns to it.
Byte-exact: same key, same `kprev` carry, same `g_fallback` bail (the `break` inside the
plain block exits the outer `while (x < n)`), same RUNMAX semantics. */
while (x < n) {
while (x + 5 <= n) { // sequential guards reduce live ranges in GCC 9.3
uint32_t k0 = (uint32_t)(RB[x] >> 20);
if (k0 == kprev) break; // run straddling the previous batch boundary
uint64_t r0 = RB[x], r1 = RB[x + 1], r2 = RB[x + 2], r3 = RB[x + 3];
uint32_t k1 = (uint32_t)(r1 >> 20);
if (k1 == k0) break; // run inside or just past the batch
uint32_t k2 = (uint32_t)(r2 >> 20);
if (k2 == k1) break;
uint32_t k3 = (uint32_t)(r3 >> 20);
if (k3 == k2) break;
if ((uint32_t)(RB[x + 4] >> 20) == k3) break;
LCPB[x] = (unsigned char)(D5B[(unsigned)__builtin_clz(kprev ^ k0)]);
LCPB[x + 1] = (unsigned char)(D5B[__builtin_clz(k0 ^ k1)]);
LCPB[x + 2] = (unsigned char)(D5B[__builtin_clz(k1 ^ k2)]);
LCPB[x + 3] = (unsigned char)(D5B[__builtin_clz(k2 ^ k3)]);
kprev = k3;
o = pack_emit(o, (uint32_t)(r0 & 0xFFFFFu));
o = pack_emit(o, (uint32_t)(r1 & 0xFFFFFu));
o = pack_emit(o, (uint32_t)(r2 & 0xFFFFFu));
o = pack_emit(o, (uint32_t)(r3 & 0xFFFFFu));
x += 4;
}
if (g_fallback) break;
{ // ONE equal-key run, then back to the fast path
uint32_t k = (uint32_t)(RB[x] >> 20);
int j = x + 1;
while (j < n && (uint32_t)(RB[j] >> 20) == k) j++;
if (j - x > RUNMAX) { g_fallback = 1; o = (char *)ostart; break; }
if (j - x > 1) {
for (int y = x + 1; y < j; y++) {
uint64_t r = RB[y];
int z = y - 1;
while (z >= x && suf_less((uint32_t)(r & 0xFFFFFu) - 1, (uint32_t)(RB[z] & 0xFFFFFu) - 1)) { RB[z + 1] = RB[z]; z--; }
RB[z + 1] = r;
}
}
for (int y = x; y < j; y++) {
uint32_t id = (uint32_t)(RB[y] & 0xFFFFFu);
uint32_t ky = (uint32_t)(RB[y] >> 20);
if (y != 0) LCPB[y] = (ky != kprev) ? (unsigned char)((__builtin_clz(kprev ^ ky) - 2) / 5)
: (unsigned char)lcp_len((uint32_t)(RB[y - 1] & 0xFFFFFu) - 1, id - 1);
kprev = ky;
if (y == n - 1) o = putu_exact(o, id);
else {
unsigned width;
uint64_t packed = pack_num(id, width);
memcpy(o, &packed, 8); o += width;
}
}
x = j;
}
}
if (!g_fallback) {
MARK(5);
*o++ = '\n';
if (n > 1) {
int y = 1;
/* [x6 lcp] eight LCP bytes -> sixteen output bytes per iteration. PAIR[d] for
d < 10 is exactly the u16 ('0'+d) | (' '<<8), i.e. byte-identical to
`_mm_unpacklo_epi8(_mm_add_epi8(d, '0'), ' ')`, and PAIR[d] == 0 for d >= 10,
which is what makes the scalar loop bail out. Any byte >= 9 (unsigned test,
so >= 128 included) therefore falls back to the original scalar loop. */
const __m128i nine = _mm_set1_epi8(9);
const __m256i digit_space = _mm256_set1_epi16(0x2030);
while (y + 16 <= n - 1) {
__m128i v = _mm_loadu_si128((const __m128i *)(LCPB + y));
if (_mm_movemask_epi8(_mm_cmpeq_epi8(_mm_min_epu8(v, nine), nine))) break;
__m256i words = _mm256_cvtepu8_epi16(v);
_mm256_storeu_si256((__m256i *)o, _mm256_add_epi16(words, digit_space));
o += 32; y += 16;
}
for (; y < n - 1; y++) o = wrs(o, LCPB[y]);
o = putu_exact(o, LCPB[n - 1]);
}
*o++ = '\n';
}
}
if (g_fallback) {
sais<unsigned char>(str, g_sa, n + 1, 256, 0);
for (int x = 0; x < n; x++) g_sa[x] = g_sa[x + 1];
char *o = *opo;
for (int x = 0; x < n - 1; x++) o = wr(o, (uint32_t)g_sa[x] + 1);
o = putu_exact(o, (uint32_t)g_sa[n - 1] + 1);
*o++ = '\n';
if (n > 1) {
for (int x = 1; x < n - 1; x++) o = wrs(o, (uint32_t)lcp_len((uint32_t)g_sa[x - 1], (uint32_t)g_sa[x]));
o = putu_exact(o, (uint32_t)lcp_len((uint32_t)g_sa[n - 2], (uint32_t)g_sa[n - 1]));
}
*o++ = '\n';
*opo = o;
return 1;
}
*opo = o;
return 1;
}
#ifdef BENCH
#include <cstdio>
#include <cstdlib>
#define OUTMAX (12u << 20)
static char OUTB[OUTMAX];
int main(int argc, char **argv) {
int n = (argc > 1) ? atoi(argv[1]) : 1000000;
unsigned long long seed = 123456789ULL;
for (int i = 0; i < n; i++) {
seed ^= seed << 13; seed ^= seed >> 7; seed ^= seed << 17;
str[i] = (unsigned char)('a' + (seed % 26));
}
str[n] = 0;
char *o = OUTB;
unsigned long long t0 = __builtin_ia32_rdtsc();
solve_string(str, n, &o);
unsigned long long t1 = __builtin_ia32_rdtsc();
unsigned long long oh = 1469598103934665603ULL;
for (long i = 0; i < (long)(o - OUTB); i++) oh = (oh ^ (unsigned char)OUTB[i]) * 1099511628211ULL;
// same checksum recipe as problems/1006e6/sol_key.cpp (0-based sa_, lcp_ 1..n-1)
unsigned long long sum = 0;
if (!g_fallback) {
for (int i = 0; i < n; i++) sum = sum * 131 + (unsigned)REC_IDX(RB[i]);
for (int i = 1; i < n; i++) sum = sum * 131 + (unsigned)LCPB[i];
} else {
for (int i = 0; i < n; i++) sum = sum * 131 + (unsigned)g_sa[i];
for (int i = 1; i < n; i++) sum = sum * 131 + (unsigned)sum; // fallback: not comparable
}
printf("n=%d total=%llu fb=%d olen=%ld ohash=%llu chk=%llu\n", n, t1 - t0, g_fallback, (long)(o - OUTB), oh, sum);
#ifdef PHASE
{ unsigned long long prev = t0; const char *nm[6]={"parse","build","radix","scanruns","emit+LCPB","lcpline"};
for (int k = 0; k < 6; k++) { printf("%s=%llu ", nm[k], (P[k]-prev)/1000); prev = P[k]; }
printf("perbucket=%llu walk=%llu afterwalk=%llu\n",(P[3]-P[2])/1000,(P[4]-P[3])/1000,(t1-P[5])/1000); }
#endif
if (argc > 2 && argv[2][0] == 'd') { fwrite(OUTB, 1, (size_t)(o - OUTB), stdout); }
return 0;
}
#elif defined(JUDGE)
struct DI {
unsigned long abi;
const char *s; unsigned long sn;
char *o; unsigned long ol; unsigned long os;
char *e; unsigned long el; unsigned long es;
const char *IB; unsigned long IBl;
char *OB; unsigned long OBl;
unsigned long tsc;
} __attribute__((packed));
extern "C" void __libc_start_main(void *m, int argc, char **argv) {
(void)m;
unsigned long *p = (unsigned long *)(argv + argc + 1);
while (*p) p++;
p++;
DI *d = 0;
for (int i = 0; i < 32 && p[0]; i++, p += 2)
if (p[0] == 0x6b637564UL) { d = (DI *)p[1]; break; }
if (d) {
char *o = d->o;
if (solve_string((const unsigned char *)d->s, (long)d->sn, &o)) d->os = (unsigned long)(o - d->o);
}
__asm__ volatile("syscall" ::"a"(60), "D"(0) : "rcx", "r11", "memory");
for (;;);
}
int main() { return 0; }
#endif