// ===== 引用区 =====
// [1] 原作者账号: saffah_cc_v41_agg1
// 原提交 URL: https://duck.ac/submission/130654
// 所用内容: 直接完整复制其计算引擎正文(为满足代码长度限制删除散文注释;代码逐字节等价,
// 以同一输入在评测机路径驱动下校验输出 checksum 与未删注释版相同),跨题用于本题。
// 许可证遵守: duck.ac 公开提交均未附独立许可证声明;按项目 INSTRUCTIONS 要求在同一评测站内
// 逐项注明原作者账号与原提交 URL 后使用;本发不新增第三方代码。
// [2] 上述提交正文内继承的全链引用(原注释所列提交号,逐一列明):
// #103361, #103367, #106719, #106747, #106762, #108184, #108329, #108359, #108399, #108401, #108502, #109438, #111338, #121252, #122788, #122903, #123047, #123142, #126143, #126941, #127524, #128321, #128473, #129085, #129575, #130301, #130561, #130562
// ===== 思路区 =====
// 本次提交的目的(校准/基线): 逐字采用本题当前 T 件引擎正文(去散文注释,代码逐字节等价;
// 本地同输入 checksum 与未删注释版相符),以本账号提交复现读数,建立本账号在本题的官方基线,
// 供后续单变量实验对比。算法、数据、常数均未改,本次不声称算法新思路。
// ============================================================
#line 397
#define MF_NOSDMP 1
#pragma GCC optimize("O3,unroll-loops,no-strict-aliasing,schedule-insns,web")
#pragma GCC optimize("ira-loop-pressure")
#pragma GCC target("avx2,bmi,bmi2,popcnt,lzcnt,tune=skylake")
#undef __code_model_small__
#undef __MMX_WITH_SSE__
#define __tune_skylake__ 1
#define __code_model_small__ 1
#define __SSE3__ 1
#define __SSSE3__ 1
#define __SSE4_1__ 1
#define __SSE4_2__ 1
#define __AVX__ 1
#define __AVX2__ 1
#define __BMI__ 1
#define __BMI2__ 1
#define __LZCNT__ 1
#define __CRC32__ 1
#define __POPCNT__ 1
#define __XSAVE__ 1
#define __MMX_WITH_SSE__ 1
#define MF_QCOEF 140
#define MF_NOFILTB 1
#ifndef MF_SPL_ALIGN
#define MF_SPL_ALIGN 1011712
#endif
#define MF_SPL ((size_t)MF_SPL_ALIGN)
#include <sys/mman.h>
#ifndef MADV_HUGEPAGE
#define MADV_HUGEPAGE 14
#endif
#include <immintrin.h>
#include <type_traits>
#pragma GCC optimize("tracer")
#include <type_traits>
#ifdef LOCAL_TEST
#include <cstdio>
#include <cstdlib>
#include <ctime>
#else
#ifdef __cplusplus
extern "C" void *malloc(unsigned long);
#else
void *malloc(unsigned long);
#endif
#endif
typedef unsigned u32;
typedef unsigned char u8;
typedef unsigned short u16;
typedef unsigned long long u64;
#define LDU256(p) ({ __m256i _v256; __asm__("vmovdqu %1, %0" : "=x"(_v256) : "m"(*(const __m256i *)(p))); _v256; })
#ifndef MAXN
#define MAXN 1000005
#endif
#define MAXK 8
#define STRIDE 8
#define PSTRIDE 16u
#ifndef MF_APF
#define MF_APF 0
#endif
#ifndef MF_QCOEF
#define MF_QCOEF 140
#endif
#ifndef MF_KSPAD
#define MF_KSPAD 1
#endif
#ifndef MF_FILTAUTO
#define MF_FILTAUTO 1
#endif
#ifndef MF_NT
#define MF_NT 1
#endif
#ifndef MF_FR_NOSORT
#define MF_FR_NOSORT 0
#endif
#ifndef MF_FR_NOCAND
#define MF_FR_NOCAND 0
#endif
#ifndef MF_BUILD_NOPOOL
#define MF_BUILD_NOPOOL 0
#endif
#ifndef MF_FR_NOREAD
#define MF_FR_NOREAD 0
#endif
#ifndef MF_NOSDMP
#define MF_NOSDMP 0
#endif
#ifndef MF_FB2
#define MF_FB2 0
#endif
#ifndef MF_FB3
#define MF_FB3 1
#endif
#ifndef MF_PLEN
#define MF_PLEN 1
#endif
#ifndef MF_BUDGETMB
#define MF_BUDGETMB 1800
#endif
#ifndef MF_NOFILTB
#define MF_NOFILTB 0
#endif
#ifndef MF_PFB
#define MF_PFB 6
#define MF_PFSC 32
#define MF_PFSC2 32
#endif
#ifndef MF_ASCR
#define MF_ASCR 1
#endif
#ifndef MF_WTRIM
#define MF_WTRIM 0
#endif
#define AS_WORDS 20480
#ifndef MF_PFA
#define MF_PFA 2
#endif
#ifndef MF_FR_NOMETA
#define MF_FR_NOMETA 0
#endif
#ifndef MAXQ
#define MAXQ 4096
#endif
extern "C" int madvise(void *addr, size_t len, int advice) noexcept;
static u32 s_ord[MAXK][MAXN];
static u32 s_frord[MAXK][MAXN];
static u32 s_gt[MAXK][MAXN];
static u32 s_pm[(size_t)MAXN * PSTRIDE] __attribute__((aligned(64)));
#define E6Z_BDR (((size_t)(MAXN) + 31u) & ~(size_t)31u)
static u16 s_bdT[(size_t)MAXK * E6Z_BDR];
static u16 s_bd[(size_t)MAXN * MAXK];
static u8 s_g[MAXN] __attribute__((aligned(64)));
static u32 s_qord[MAXN], s_tmp[MAXN], s_cntv[MAXN];
static u16 s_bkt[MAXN];
static u16 s_bks[MAXN];
static u32 s_nid[MAXN];
static u32 s_pos[MAXK][MAXN];
static u32 s_posp[(size_t)MAXK * MAXN];
#ifndef MF_INCA
#define MF_INCA 1
#endif
#ifndef MF_FRPF
#define MF_FRPF 32
#endif
#ifndef MF_INCA_PFD
#define MF_INCA_PFD 32
#endif
static int g_inca = 0;
static u32 *s_ho = 0;
static size_t s_ho_cap = 0;
static u64 s_facc[(size_t)(MAXK - 1) * (((MAXN) + 63) / 64 + 16)];
static u32 s_bnde[(size_t)MAXK * (MAXQ + 2)];
static u64 *s_fbm = 0;
static u64 s_fbuf[8192] __attribute__((aligned(64)));
static __m256i khi8v, khi8hv;
static u32 g_Qc = 32;
static u32 g_fbm_W = 0;
static int g_use_filt = 0;
static int g_lexsort = 0;
static int g_nofilt = 0;
static u32 s_fpos[(MAXK + 1) * MAXQ];
static u32 s_bval[(MAXK + 1) * MAXQ];
static u64 *s_bs = 0;
static u64 s_acc[MAXN / 64 + 256] __attribute__((aligned(64)));
static u32 *s_sdmp[MAXK];
static u64 s_ascr[AS_WORDS] __attribute__((aligned(64)));
static u32 s_maxr[(size_t)MAXK * MAXK * MAXQ];
static u8 *s_sig = 0;
static u8 b14af_csig[(size_t)MAXN*8] __attribute__((aligned(64)));
static int g_sig_on = 0;
static int g_sig_merged = 0;
static u32 g_sigsh = 0;
static u32 s_pofs[(size_t)MAXK * MAXQ];
static u32 s_plen[(size_t)MAXK * MAXQ];
static int g_plen_on = 0;
static u32 g_W, g_Q, g_B;
static u32 solve_generation;
static u32 g_KS;
static int g_nosort = 0;
static int g_nosdmp = 0;
static int g_skip_and = 0, g_skip_fringe = 0;
#define FPOS(d,b) s_fpos[(size_t)(d) * MAXQ + (b)]
#define BVAL(d,b) s_bval[(size_t)(d) * MAXQ + (b)]
static void *pool_alloc(size_t bytes) {
void *raw = malloc(bytes + (2u << 20));
if (!raw) return 0;
void *p = (void *)(((size_t)raw + (2u << 20) - 1) & ~(size_t)((2u << 20) - 1));
{ void *m = p; size_t mb = bytes & ~(size_t)((2u << 20) - 1);
if (mb) madvise(m, mb, MADV_HUGEPAGE); }
return p;
}
__attribute__((target("avx2,popcnt")))
static inline u32 hsum256(__m256i v) {
__m128i lo = _mm256_castsi256_si128(v);
__m128i hi = _mm256_extracti128_si256(v, 1);
__m128i s = _mm_add_epi64(lo, hi);
return (u32)((u64)_mm_cvtsi128_si64(s) + (u64)_mm_extract_epi64(s, 1));
}
template<int KK>
__attribute__((target("avx2,popcnt")))
static u32 and_popcount_range(const u64 *const *bs, u32 Wr, u32 rem, const u64 *const *bsNext) {
const __m256i lm = _mm256_set1_epi8(0x0f);
const __m256i lk = _mm256_setr_epi8(0,1,1,2,1,2,2,3,1,2,2,3,2,3,3,4,
0,1,1,2,1,2,2,3,1,2,2,3,2,3,3,4);
const __m256i z = _mm256_setzero_si256();
__m256i tot = z, acc0 = z, acc1 = z;
u32 w = 0, cnt = 0, total = 0;
u32 wlim = (Wr > 96u) ? (Wr - 96u) : 0u;
for (; w + 8 <= wlim; w += 8) {
if (KK > 1) _mm_prefetch((const char *)(bs[1] + w + 96), _MM_HINT_T0);
if (KK > 2) _mm_prefetch((const char *)(bs[2] + w + 96), _MM_HINT_T0);
if (KK > 3) _mm_prefetch((const char *)(bs[3] + w + 96), _MM_HINT_T0);
if (KK > 4) _mm_prefetch((const char *)(bs[4] + w + 96), _MM_HINT_T0);
if (KK > 5) _mm_prefetch((const char *)(bs[5] + w + 96), _MM_HINT_T0);
if (KK > 6) _mm_prefetch((const char *)(bs[6] + w + 96), _MM_HINT_T0);
if (KK > 7) _mm_prefetch((const char *)(bs[7] + w + 96), _MM_HINT_T0);
__m256i a = LDU256((bs[0] + w));
for (int d = 1; d < KK; d++) a = _mm256_and_si256(a, LDU256((bs[d] + w)));
__m256i lo0 = _mm256_and_si256(a, lm);
__m256i hi0 = _mm256_and_si256(_mm256_srli_epi16(a, 4), lm);
acc0 = _mm256_add_epi8(acc0, _mm256_add_epi8(_mm256_shuffle_epi8(lk, lo0), _mm256_shuffle_epi8(lk, hi0)));
__m256i b = LDU256((bs[0] + w + 4));
for (int d = 1; d < KK; d++) b = _mm256_and_si256(b, LDU256((bs[d] + w + 4)));
__m256i lo1 = _mm256_and_si256(b, lm);
__m256i hi1 = _mm256_and_si256(_mm256_srli_epi16(b, 4), lm);
acc1 = _mm256_add_epi8(acc1, _mm256_add_epi8(_mm256_shuffle_epi8(lk, lo1), _mm256_shuffle_epi8(lk, hi1)));
if (++cnt == 31) {
tot = _mm256_add_epi64(tot, _mm256_sad_epu8(acc0, z));
tot = _mm256_add_epi64(tot, _mm256_sad_epu8(acc1, z));
acc0 = acc1 = z; cnt = 0;
}
}
if (bsNext) {
u32 off = 0;
for (; w + 8 <= Wr; w += 8, off += 8) {
if (KK > 1) _mm_prefetch((const char *)(bsNext[1] + off), _MM_HINT_T0);
if (KK > 2) _mm_prefetch((const char *)(bsNext[2] + off), _MM_HINT_T0);
if (KK > 3) _mm_prefetch((const char *)(bsNext[3] + off), _MM_HINT_T0);
if (KK > 4) _mm_prefetch((const char *)(bsNext[4] + off), _MM_HINT_T0);
if (KK > 5) _mm_prefetch((const char *)(bsNext[5] + off), _MM_HINT_T0);
if (KK > 6) _mm_prefetch((const char *)(bsNext[6] + off), _MM_HINT_T0);
if (KK > 7) _mm_prefetch((const char *)(bsNext[7] + off), _MM_HINT_T0);
__m256i a = LDU256((bs[0] + w));
for (int d = 1; d < KK; d++) a = _mm256_and_si256(a, LDU256((bs[d] + w)));
__m256i lo0 = _mm256_and_si256(a, lm);
__m256i hi0 = _mm256_and_si256(_mm256_srli_epi16(a, 4), lm);
acc0 = _mm256_add_epi8(acc0, _mm256_add_epi8(_mm256_shuffle_epi8(lk, lo0), _mm256_shuffle_epi8(lk, hi0)));
__m256i b = LDU256((bs[0] + w + 4));
for (int d = 1; d < KK; d++) b = _mm256_and_si256(b, LDU256((bs[d] + w + 4)));
__m256i lo1 = _mm256_and_si256(b, lm);
__m256i hi1 = _mm256_and_si256(_mm256_srli_epi16(b, 4), lm);
acc1 = _mm256_add_epi8(acc1, _mm256_add_epi8(_mm256_shuffle_epi8(lk, lo1), _mm256_shuffle_epi8(lk, hi1)));
if (++cnt == 31) {
tot = _mm256_add_epi64(tot, _mm256_sad_epu8(acc0, z));
tot = _mm256_add_epi64(tot, _mm256_sad_epu8(acc1, z));
acc0 = acc1 = z; cnt = 0;
}
}
}
for (; w + 8 <= Wr; w += 8) {
__m256i a = LDU256((bs[0] + w));
for (int d = 1; d < KK; d++) a = _mm256_and_si256(a, LDU256((bs[d] + w)));
__m256i lo0 = _mm256_and_si256(a, lm);
__m256i hi0 = _mm256_and_si256(_mm256_srli_epi16(a, 4), lm);
acc0 = _mm256_add_epi8(acc0, _mm256_add_epi8(_mm256_shuffle_epi8(lk, lo0), _mm256_shuffle_epi8(lk, hi0)));
__m256i b = LDU256((bs[0] + w + 4));
for (int d = 1; d < KK; d++) b = _mm256_and_si256(b, LDU256((bs[d] + w + 4)));
__m256i lo1 = _mm256_and_si256(b, lm);
__m256i hi1 = _mm256_and_si256(_mm256_srli_epi16(b, 4), lm);
acc1 = _mm256_add_epi8(acc1, _mm256_add_epi8(_mm256_shuffle_epi8(lk, lo1), _mm256_shuffle_epi8(lk, hi1)));
if (++cnt == 31) {
tot = _mm256_add_epi64(tot, _mm256_sad_epu8(acc0, z));
tot = _mm256_add_epi64(tot, _mm256_sad_epu8(acc1, z));
acc0 = acc1 = z; cnt = 0;
}
}
for (; w + 4 <= Wr; w += 4) {
__m256i a = LDU256((bs[0] + w));
for (int d = 1; d < KK; d++) a = _mm256_and_si256(a, LDU256((bs[d] + w)));
__m256i lo2 = _mm256_and_si256(a, lm);
__m256i hi2 = _mm256_and_si256(_mm256_srli_epi16(a, 4), lm);
acc0 = _mm256_add_epi8(acc0, _mm256_add_epi8(_mm256_shuffle_epi8(lk, lo2), _mm256_shuffle_epi8(lk, hi2)));
}
tot = _mm256_add_epi64(tot, _mm256_sad_epu8(acc0, z));
tot = _mm256_add_epi64(tot, _mm256_sad_epu8(acc1, z));
total = hsum256(tot);
for (; w < Wr; w++) {
u64 v = bs[0][w];
for (int d = 1; d < KK; d++) v &= bs[d][w];
total += (u32)__builtin_popcountll(v);
}
if (rem) {
u64 v = bs[0][Wr];
for (int d = 1; d < KK; d++) v &= bs[d][Wr];
total += (u32)__builtin_popcountll(v & ((1ull << rem) - 1));
}
return total;
}
static const u32 g_r2j_kmaskT[32]={0u,1u,3u,7u,15u,31u,63u,127u,255u,511u,1023u,2047u,
4095u,8191u,16383u,32767u,65535u,131071u,262143u,524287u,1048575u,2097151u,
4194303u,8388607u,16777215u,33554431u,67108863u,134217727u,268435455u,
536870911u,1073741823u,2147483647u};
#define TIE_CAP 524288u
#define TIE_EW 12
alignas(64) static u32 tie_data[(size_t)TIE_CAP*TIE_EW],tie_points[TIE_CAP];
static u32 tie_count;
static void flush_ties(u32 *out){
u32 n=tie_count;tie_count=0;
const __m256i all8=_mm256_set1_epi32(-1);
for(u32 b=0;b<n;b+=256){
u32 end=b+256<n?b+256:n;
for(u32 j=b;j<end;j++){
const u32 *e=tie_data+(size_t)j*TIE_EW;
u32 p=s_ord[e[2]][e[1]];tie_points[j]=p;
_mm_prefetch((const char*)(s_pm+(size_t)p*PSTRIDE),_MM_HINT_T0);
}
for(u32 j=b;j<end;j++){
const u32 *e=tie_data+(size_t)j*TIE_EW;
__m256i v=LDU256(s_pm+(size_t)tie_points[j]*PSTRIDE);
__m256i q=_mm256_loadu_si256((const __m256i*)(e+3));
out[e[0]]+=(u32)_mm256_testc_si256(_mm256_cmpgt_epi32(q,v),all8);
}
}
}
template<int KK>
__attribute__((target("avx2,popcnt")))
static void query_group(u32 N, u32 *out, u32 ds, u32 q0, u32 q1, u32 Q, u32 W) {
const int k = KK;
const int K1 = KK - 1;
const u32 ASTRL = (KK <= 9) ? 8u : 16u;
#if MF_ASCR
const u32 d1 = (ds == 0) ? 1u : 0u;
const u32 e1 = (d1 < ds) ? d1 : d1 - 1;
{
const u32 *const fo1 = s_frord[d1];
for (u32 qi = q0; qi < q1; qi++) s_qord[qi] = fo1[qi] & 0xFFFFFu;
}
#else
const u32 d1 = 0, e1 = 0;
#endif
if (g_lexsort) {
u32 *cntv_ = s_cntv;
for (int dd = k - 1; dd >= 0; dd--) {
if (dd == ds) continue;
for (u32 v = 0; v <= Q; v++) cntv_[v] = 0;
for (u32 qi = q0; qi < q1; qi++) cntv_[s_bdT[(size_t)dd * E6Z_BDR + s_qord[qi]] + 1]++;
{ u32 ac = q0; for (u32 v = 0; v <= Q; v++) { u32 c = cntv_[v]; cntv_[v] = ac; ac += c; } }
for (u32 qi = q0; qi < q1; qi++) { u32 id = s_qord[qi]; s_tmp[cntv_[s_bdT[(size_t)dd * E6Z_BDR + id] + 1]++] = id; }
for (u32 qi = q0; qi < q1; qi++) s_qord[qi] = s_tmp[qi];
}
}
#ifdef PF_DIST
for (u32 qi = q0; qi < q0 + PF_DIST && qi + PF_DIST < q1; qi++) {
u32 j = s_qord[qi + PF_DIST];
const u16 *bdn = s_bd + (size_t)j * MAXK;
for (int e = 0; e < K1; e++) {
int d = (e < ds) ? e : e + 1;
const char *pp = (const char *)(s_bs + ((size_t)e * (Q + 1) + bdn[d]) * W);
_mm_prefetch(pp, _MM_HINT_T0);
_mm_prefetch(pp + 64, _MM_HINT_T0);
_mm_prefetch(pp + 128, _MM_HINT_T0);
}
}
#endif
#if MF_ASCR
{
#if MF_INCA
const u32 *a13_ho = 0;
u32 a13_app = 0;
if (g_inca) {
a13_ho = s_ho;
for (u32 w = 0; w < AS_WORDS; w++) s_ascr[w] = 0;
}
#endif
u32 qi = q0;
while (qi < q1) {
const u32 *const fo1b = s_frord[d1];
u32 b1 = fo1b[qi] >> 20;
u32 qe = qi + 1;
while (qe < q1 && (fo1b[qe] >> 20) == b1) qe++;
const u64 *bp1 = g_plen_on ? (s_bs + s_pofs[(size_t)e1 * MAXQ + b1])
: (s_bs + ((size_t)e1 * (Q + 1) + b1) * W);
#if MF_INCA
if (a13_ho) {
u32 pto = (b1 < Q) ? FPOS(d1, b1) : (u32)N;
const u32 *hp = a13_ho;
#if MF_INCA_PFD > 0
while (a13_app < pto) {
if (a13_app + MF_INCA_PFD < pto)
_mm_prefetch((const char *)(s_ascr + (hp[a13_app + MF_INCA_PFD] >> 6)), _MM_HINT_T0);
u32 q = hp[a13_app]; s_ascr[q >> 6] |= 1ull << (q & 63); a13_app++; }
#else
while (a13_app < pto) { u32 q = hp[a13_app]; s_ascr[q >> 6] |= 1ull << (q & 63); a13_app++; }
#endif
bp1 = s_ascr;
} else
#endif
{
u32 mx = 1;
for (u32 t = qi; t < qe; t++) { u32 rr = s_gt[ds][s_qord[t]]; if (rr > mx) mx = rr; }
u32 nwd = (mx >> 6) + 1;
if (nwd <= AS_WORDS) {
const u64 *src = bp1;
for (u32 w = 0; w < nwd; w++) s_ascr[w] = src[w];
bp1 = s_ascr;
}
}
const u64 *bp2v[2][MAXK];
u32 Rcur = s_gt[ds][s_qord[qi]];
{
const u16 *bd0 = s_bd + (size_t)s_qord[qi] * MAXK;
for (int e = 0; e < K1; e++) {
if (e == (int)e1) bp2v[qi & 1u][e] = bp1;
else { int d = (e < ds) ? e : e + 1;
bp2v[qi & 1u][e] = g_plen_on ? (s_bs + s_pofs[(size_t)e * MAXQ + bd0[d]])
: (s_bs + ((size_t)e * (Q + 1) + bd0[d]) * W); }
}
}
for (u32 qq = qi; qq < qe; qq++) {
u32 i = s_qord[qq];
#if MF_PFA > 0
{
u32 qj = qq + MF_PFA;
if (qj < qe) {
u32 j = s_qord[qj];
_mm_prefetch((const char *)(s_bd + (size_t)j * MAXK), _MM_HINT_T0);
_mm_prefetch((const char *)(s_gt[ds] + j), _MM_HINT_T0);
__asm__ __volatile__("prefetchw %0" :: "m"(*(const char *)(out + j)));
}
}
#endif
if (qq + 1 < qe) {
const u16 *bdn = s_bd + (size_t)s_qord[qq + 1] * MAXK;
const u64 **bn = bp2v[(qq + 1) & 1u];
for (int e = 0; e < K1; e++) {
if (e == (int)e1) bn[e] = bp1;
else { int d = (e < ds) ? e : e + 1;
bn[e] = g_plen_on ? (s_bs + s_pofs[(size_t)e * MAXQ + bdn[d]])
: (s_bs + ((size_t)e * (Q + 1) + bdn[d]) * W); }
}
}
const u32 R = Rcur;
if (qq + 1 < qe) Rcur = s_gt[ds][s_qord[qq + 1]];
u32 cnt = 0;
if (R && !g_skip_and) {
cnt = and_popcount_range<KK - 1>(bp2v[qq & 1u], R >> 6, R & 63u,
(qq + 1 < qe) ? bp2v[(qq + 1) & 1u] : 0);
}
out[i] = cnt;
}
qi = qe;
}
}
#else
for (u32 qi = q0; qi < q1; qi++) {
u32 i = s_qord[qi];
const u16 *bd = s_bd + (size_t)i * MAXK;
#if MF_PFA > 0
{
u32 qj = qi + MF_PFA;
if (qj < q1) {
u32 j = s_qord[qj];
_mm_prefetch((const char *)(s_bd + (size_t)j * MAXK), _MM_HINT_T0);
_mm_prefetch((const char *)(s_gt[ds] + j), _MM_HINT_T0);
_mm_prefetch((const char *)(out + j), _MM_HINT_T0);
}
}
#endif
#ifdef PF_DIST
if (qi + PF_DIST < q1) {
u32 j = s_qord[qi + PF_DIST];
const u16 *bdn = s_bd + (size_t)j * MAXK;
for (int e = 0; e < K1; e++) {
int d = (e < ds) ? e : e + 1;
const char *pp = (const char *)(s_bs + ((size_t)e * (Q + 1) + bdn[d]) * W);
_mm_prefetch(pp, _MM_HINT_T0);
_mm_prefetch(pp + 64, _MM_HINT_T0);
_mm_prefetch(pp + 128, _MM_HINT_T0);
}
}
#endif
const u32 R = s_gt[ds][i];
u32 cnt = 0;
if (R && !g_skip_and) {
const u64 *bp[MAXK];
for (int e = 0; e < K1; e++) {
int d = (e < ds) ? e : e + 1;
bp[e] = g_plen_on ? (s_bs + s_pofs[(size_t)e * MAXQ + bd[d]])
: (s_bs + ((size_t)e * (Q + 1) + bd[d]) * W);
#ifdef DBG_PLEN
if (g_plen_on) {
u32 need = (R >> 6) + ((R & 63) ? 1 : 0);
if (need > s_plen[(size_t)e * MAXQ + bd[d]])
printf("DBG ds=%d e=%d d=%d b=%u R=%u need=%u len=%u\n", ds, e, d, bd[d], R, need, s_plen[(size_t)e * MAXQ + bd[d]]);
}
#endif
}
cnt = and_popcount_range<KK - 1>(bp, R >> 6, R & 63u, 0);
}
out[i] = cnt;
}
#endif
if (g_skip_fringe) return;
static u32 *avpool = 0; static size_t avcap = 0;
{
size_t needav = (size_t)N * STRIDE * 4 + 256;
static u32 avN = 0;
if (needav > avcap) { avpool = (u32 *)pool_alloc(needav); avcap = avpool ? needav : 0; avN = 0; }
if ((int)PSTRIDE == STRIDE && avpool && avN != solve_generation) {
for (u32 ii = 0; ii < N; ii++) {
const u16 *bdi = s_bd + (size_t)ii * MAXK;
u32 *av = avpool + (size_t)ii * ASTRL;
for (int e = 0; e < k && e < (int)ASTRL; e++) av[e] = BVAL(e, bdi[e]);
for (int e = k; e < (int)ASTRL; e++) av[e] = 0x7FFFFFFFu;
}
avN = solve_generation;
}
}
auto fringe_dim = [&](auto dim) {
constexpr int d=decltype(dim)::value;
if(d==ds) return;
const u32 *fo = s_frord[d];
u32 exmask = 0;
for (int e = 0; e < d; e++) if (e != ds) exmask |= 1u << e;
u32 hiA[8], hiAh[8];
for (int e = 0; e < 8; e++) hiA[e] = ((exmask >> e) & 1u) ? 0u : 0x7FFFFFFFu;
for (int e = 8; e < 16; e++) hiAh[e - 8] = ((exmask >> e) & 1u) ? 0u : 0x7FFFFFFFu;
const __m256i notexlo = LDU256(hiA);
const __m256i notexloh = LDU256(hiAh);
const __m256i all8 = _mm256_set1_epi32(-1);
const u8 *const sbase = s_sig + (size_t)d * (k-1) * MF_SPL;
const u32 *const odd = s_ord[d];
for (u32 qi = q0; qi < q1; qi++) {
u32 i = fo[qi] & 0xFFFFFu;
#if MF_PFB > 0
{
u32 qj = qi + MF_PFB;
if (qj < q1) {
u32 j = fo[qj] & 0xFFFFFu;
_mm_prefetch((const char *)(s_pm + (size_t)j * PSTRIDE), _MM_HINT_T0);
_mm_prefetch((const char *)(s_gt[d] + j), _MM_HINT_T0);
__asm__ __volatile__("prefetchw %0" :: "m"(*(const char *)(out + j)));
}
}
#endif
const u32 lo = FPOS(d, fo[qi] >> 20);
const u32 hi = s_gt[d][i];
if (lo >= hi) continue;
const u32 *rec = s_pm + (size_t)i * PSTRIDE;
const u32 *Avp = avpool + (size_t)i * ASTRL;
__m256i qr = _mm256_or_si256(_mm256_load_si256((const __m256i *)rec), khi8v);
__m256i Avm = _mm256_or_si256((PSTRIDE > STRIDE) ? LDU256(rec + 8) : LDU256(Avp), notexlo);
__m256i qT = _mm256_min_epi32(qr, Avm);
__m256i qTh;
if (k > 9) {
__m256i qrh = _mm256_or_si256(LDU256((rec + 8)), khi8hv);
__m256i Avmh = _mm256_or_si256(LDU256((Avp + 8)), notexloh);
qTh = _mm256_min_epi32(qrh, Avmh);
} else if (k > 8) {
qTh = _mm256_or_si256(LDU256((rec + 8)), khi8hv);
} else {
qTh = all8;
}
u32 add = 0;
const u32 *c0 = s_sdmp[d] ? (s_sdmp[d] + (size_t)lo * g_KS) : 0;
if (g_use_filt && KK == 9) {
const u32 w0 = lo >> 6;
const u32 nw = ((hi - 1) >> 6) - w0 + 1;
const u64 wcs = (u64)N / g_Qc + 1;
const u32 FW = g_fbm_W;
const u64 *rows[8];
int nr = 0;
for (int e = 0; e < 9; e++) {
if (e == d) continue;
u32 cell = (u32)((u64)rec[e] / wcs) + 1;
if (cell > g_Qc) cell = g_Qc;
rows[nr++] = s_fbm + ((size_t)(d * 9 + e) * (g_Qc + 1) + cell) * FW + w0;
}
const u64 *r0 = rows[0];
const u64 *r1 = rows[1];
const u64 *r2 = rows[2];
const u64 *r3 = rows[3];
const u64 *r4 = rows[4];
const u64 *r5 = rows[5];
const u64 *r6 = rows[6];
const u64 *r7 = rows[7];
const u32 remh = hi & 63u;
{
u64 vv = r0[0] & r1[0] & r2[0] & r3[0] & r4[0] & r5[0] & r6[0] & r7[0];
vv &= ~0ull << (lo & 63);
if (nw == 1u && remh) vv &= (1ull << remh) - 1;
u32 word = w0;
#if MF_FR_NOCAND
vv = 0;
#endif
while (vv) {
int t = __builtin_ctzll(vv);
vv &= vv - 1;
u32 j = (word << 6) + (u32)t;
const u32 *cp = c0 ? (c0 + (size_t)(j - lo) * g_KS) : (s_pm + (size_t)s_ord[d][j] * PSTRIDE);
__m256i pv = LDU256(cp);
__m256i pvh = LDU256(cp + 8);
u32 ok = (u32)_mm256_testc_si256(_mm256_cmpgt_epi32(qT, pv), all8);
ok &= (u32)_mm256_testc_si256(_mm256_cmpgt_epi32(qTh, pvh), all8);
add += ok;
}
}
for (u32 w = 1; w + 1 < nw; w++) {
u64 vv = r0[w] & r1[w] & r2[w] & r3[w] & r4[w] & r5[w] & r6[w] & r7[w];
u32 word = w0 + w;
#if MF_FR_NOCAND
vv = 0;
#endif
while (vv) {
int t = __builtin_ctzll(vv);
vv &= vv - 1;
u32 j = (word << 6) + (u32)t;
const u32 *cp = c0 ? (c0 + (size_t)(j - lo) * g_KS) : (s_pm + (size_t)s_ord[d][j] * PSTRIDE);
__m256i pv = LDU256(cp);
__m256i pvh = LDU256(cp + 8);
u32 ok = (u32)_mm256_testc_si256(_mm256_cmpgt_epi32(qT, pv), all8);
ok &= (u32)_mm256_testc_si256(_mm256_cmpgt_epi32(qTh, pvh), all8);
add += ok;
}
}
if (nw > 1u) {
u32 w = nw - 1u;
u64 vv = r0[w] & r1[w] & r2[w] & r3[w] & r4[w] & r5[w] & r6[w] & r7[w];
if (remh) vv &= (1ull << remh) - 1;
u32 word = w0 + w;
#if MF_FR_NOCAND
vv = 0;
#endif
while (vv) {
int t = __builtin_ctzll(vv);
vv &= vv - 1;
u32 j = (word << 6) + (u32)t;
const u32 *cp = c0 ? (c0 + (size_t)(j - lo) * g_KS) : (s_pm + (size_t)s_ord[d][j] * PSTRIDE);
__m256i pv = LDU256(cp);
__m256i pvh = LDU256(cp + 8);
u32 ok = (u32)_mm256_testc_si256(_mm256_cmpgt_epi32(qT, pv), all8);
ok &= (u32)_mm256_testc_si256(_mm256_cmpgt_epi32(qTh, pvh), all8);
add += ok;
}
}
out[i] += add;
continue;
}
if (g_use_filt) {
const u32 w0 = lo >> 6;
const u32 nw = ((hi - 1) >> 6) - w0 + 1;
u64 *fb = s_fbuf;
const u64 wcs = (u64)N / g_Qc + 1;
const u32 FW = g_fbm_W;
int first = 1;
#if MF_FR_NOREAD
first = 0;
for (u32 w = 0; w < nw; w++) fb[w] = 0;
#endif
for (int e = 0; e < k; e++) {
if (e == d) continue;
u32 c = (u32)((u64)rec[e] / wcs) + 1;
if (c > g_Qc) c = g_Qc;
const u64 *row = s_fbm + ((size_t)(d * k + e) * (g_Qc + 1) + c) * FW + w0;
if (first) { for (u32 w = 0; w < nw; w++) fb[w] = row[w]; first = 0; }
else { for (u32 w = 0; w < nw; w++) fb[w] &= row[w]; }
}
fb[0] &= ~0ull << (lo & 63);
u32 r2 = hi & 63u;
if (r2) fb[nw - 1] &= (1ull << r2) - 1;
for (u32 w = 0; w < nw; w++) {
u64 vv = fb[w];
#if MF_FR_NOCAND
vv = 0;
#endif
while (vv) {
int t = __builtin_ctzll(vv);
vv &= vv - 1;
u32 j = ((w0 + w) << 6) + (u32)t;
const u32 *cp = c0 ? (c0 + (size_t)(j - lo) * g_KS) : (s_pm + (size_t)s_ord[d][j] * PSTRIDE);
__m256i pv = LDU256(cp);
u32 ok = (u32)_mm256_testc_si256(_mm256_cmpgt_epi32(qT, pv), all8);
if (k > 8) {
__m256i pvh = LDU256((cp + 8));
ok &= (u32)_mm256_testc_si256(_mm256_cmpgt_epi32(qTh, pvh), all8);
}
add += ok;
}
}
out[i] += add;
continue;
}
if (g_sig_on && (hi - lo) >= 2) {
__m256i qsb[16];
{
const __m256i zero8 = _mm256_setzero_si256();
const __m256i shv = _mm256_set1_epi32((int)g_sigsh);
const __m256i bias8 = _mm256_set1_epi8((char)0x80);
const __m256i sv_lo = _mm256_srlv_epi32(qT, shv);
const __m128i u16=_mm_packus_epi32(_mm256_castsi256_si128(sv_lo),_mm256_extracti128_si256(sv_lo,1));
const __m128i u8s=_mm_packus_epi16(u16,u16);
const __m256i packed=_mm256_xor_si256(_mm256_broadcastsi128_si256(u8s),bias8);
for (int e = 0; e < k; e++) {
if(e==(int)d)continue;
qsb[e-(e>(int)d)]=_mm256_shuffle_epi8(packed,_mm256_set1_epi8((char)e));
}
}
u32 aj = lo, aacc = 0, aexact = 0;
for (; aj + 32 <= hi; aj += 32) {
__m256i ac=_mm256_subs_epi8(LDU256(sbase+aj),qsb[0]);
for(int e=1;e<k-1;e++){
__m256i cv=LDU256(sbase+(size_t)e*MF_SPL+aj);
ac=_mm256_max_epi8(ac,_mm256_subs_epi8(cv,qsb[e]));
}
u32 decm = (u32)_mm256_movemask_epi8(ac);
u32 amb = (u32)_mm256_movemask_epi8(_mm256_cmpeq_epi8(ac, _mm256_setzero_si256()));
aacc += (u32)__builtin_popcount(decm);
while (amb) {
int t = __builtin_ctz(amb); amb &= amb - 1;
if(KK==8 && !c0 && tie_count<TIE_CAP){
u32 *e=tie_data+(size_t)tie_count++*TIE_EW;
e[0]=i;e[1]=aj+(u32)t;e[2]=d;
_mm256_storeu_si256((__m256i*)(e+3),qT);
continue;
}
const u32 *cp = c0 ? (c0 + (size_t)(aj + (u32)t - lo) * g_KS) : (s_pm + (size_t)odd[aj + (u32)t] * PSTRIDE);
__m256i pv = LDU256(cp);
u32 okl = (u32)_mm256_testc_si256(_mm256_cmpgt_epi32(qT, pv), all8);
if (k > 8) {
__m256i pvh = LDU256((cp + 8));
okl &= (u32)_mm256_testc_si256(_mm256_cmpgt_epi32(qTh, pvh), all8);
}
aacc += okl;
}
}
{
u32 r = hi - aj;
if (r) {
u32 base = aj, lk = 0u;
u32 kmask = g_r2j_kmaskT[r];
__m256i ac=_mm256_subs_epi8(LDU256(sbase+base),qsb[0]);
for(int e=1;e<k-1;e++){
__m256i cv=LDU256(sbase+(size_t)e*MF_SPL+base);
ac=_mm256_max_epi8(ac,_mm256_subs_epi8(cv,qsb[e]));
}
u32 decm = ((u32)_mm256_movemask_epi8(ac)) & kmask;
u32 amb = ((u32)_mm256_movemask_epi8(_mm256_cmpeq_epi8(ac, _mm256_setzero_si256()))) & kmask;
aacc += (u32)__builtin_popcount(decm);
while (amb) {
int t = __builtin_ctz(amb); amb &= amb - 1;
if(KK==8 && !c0 && tie_count<TIE_CAP){
u32 *e=tie_data+(size_t)tie_count++*TIE_EW;
e[0]=i;e[1]=base+(u32)t;e[2]=d;
_mm256_storeu_si256((__m256i*)(e+3),qT);
continue;
}
const u32 *cp = c0 ? (c0 + (size_t)(base + (u32)t - lo) * g_KS) : (s_pm + (size_t)odd[base + (u32)t] * PSTRIDE);
__m256i pv = LDU256(cp);
u32 okl = (u32)_mm256_testc_si256(_mm256_cmpgt_epi32(qT, pv), all8);
if (k > 8) {
__m256i pvh = LDU256((cp + 8));
okl &= (u32)_mm256_testc_si256(_mm256_cmpgt_epi32(qTh, pvh), all8);
}
aacc += okl;
}
}
}
add = aacc;
} else if (c0) {
if (k <= 8) {
const u32 *c = c0;
u32 a0 = 0, a1 = 0, a2 = 0, a3 = 0;
u32 j = lo;
for (; j + 4 <= hi; j += 4) {
__m256i p0 = LDU256((c));
__m256i p1 = LDU256((c + g_KS));
__m256i p2 = LDU256((c + 2 * g_KS));
__m256i p3 = LDU256((c + 3 * g_KS));
a0 += (u32)_mm256_testc_si256(_mm256_cmpgt_epi32(qT, p0), all8);
a1 += (u32)_mm256_testc_si256(_mm256_cmpgt_epi32(qT, p1), all8);
a2 += (u32)_mm256_testc_si256(_mm256_cmpgt_epi32(qT, p2), all8);
a3 += (u32)_mm256_testc_si256(_mm256_cmpgt_epi32(qT, p3), all8);
c += 4 * g_KS;
}
for (; j < hi; j++, c += g_KS) {
__m256i pv = LDU256(c);
a0 += (u32)_mm256_testc_si256(_mm256_cmpgt_epi32(qT, pv), all8);
}
add = a0 + a1 + a2 + a3;
} else {
const u32 *c = c0;
u32 a0 = 0, a1 = 0, a2 = 0, a3 = 0;
u32 j = lo;
for (; j + 4 <= hi; j += 4) {
__m256i p0 = LDU256((c));
__m256i p1 = LDU256((c + g_KS));
__m256i p2 = LDU256((c + 2 * g_KS));
__m256i p3 = LDU256((c + 3 * g_KS));
__m256i h0 = LDU256((c + 8));
__m256i h1 = LDU256((c + g_KS + 8));
__m256i h2 = LDU256((c + 2 * g_KS + 8));
__m256i h3 = LDU256((c + 3 * g_KS + 8));
a0 += (u32)_mm256_testc_si256(_mm256_cmpgt_epi32(qT, p0), all8) & (u32)_mm256_testc_si256(_mm256_cmpgt_epi32(qTh, h0), all8);
a1 += (u32)_mm256_testc_si256(_mm256_cmpgt_epi32(qT, p1), all8) & (u32)_mm256_testc_si256(_mm256_cmpgt_epi32(qTh, h1), all8);
a2 += (u32)_mm256_testc_si256(_mm256_cmpgt_epi32(qT, p2), all8) & (u32)_mm256_testc_si256(_mm256_cmpgt_epi32(qTh, h2), all8);
a3 += (u32)_mm256_testc_si256(_mm256_cmpgt_epi32(qT, p3), all8) & (u32)_mm256_testc_si256(_mm256_cmpgt_epi32(qTh, h3), all8);
c += 4 * g_KS;
}
for (; j < hi; j++, c += g_KS) {
__m256i pv = LDU256(c);
__m256i pvh = LDU256((c + 8));
a0 += (u32)_mm256_testc_si256(_mm256_cmpgt_epi32(qT, pv), all8) &
(u32)_mm256_testc_si256(_mm256_cmpgt_epi32(qTh, pvh), all8);
}
add = a0 + a1 + a2 + a3;
}
} else {
for (u32 j = lo; j < hi; j++) {
const u32 *c = s_pm + (size_t)odd[j] * PSTRIDE;
__m256i pv = _mm256_load_si256((const __m256i *)c);
u32 ok = (u32)_mm256_testc_si256(_mm256_cmpgt_epi32(qT, pv), all8);
if (k > 8) {
__m256i pvh = _mm256_load_si256((const __m256i *)(c + 8));
ok &= (u32)_mm256_testc_si256(_mm256_cmpgt_epi32(qTh, pvh), all8);
}
add += ok;
}
}
out[i] += add;
}
};
if constexpr(KK>0) fringe_dim(std::integral_constant<int,0>{});
if constexpr(KK>1) fringe_dim(std::integral_constant<int,1>{});
if constexpr(KK>2) fringe_dim(std::integral_constant<int,2>{});
if constexpr(KK>3) fringe_dim(std::integral_constant<int,3>{});
if constexpr(KK>4) fringe_dim(std::integral_constant<int,4>{});
if constexpr(KK>5) fringe_dim(std::integral_constant<int,5>{});
if constexpr(KK>6) fringe_dim(std::integral_constant<int,6>{});
if constexpr(KK>7) fringe_dim(std::integral_constant<int,7>{});
if constexpr(KK>8) fringe_dim(std::integral_constant<int,8>{});
if constexpr(KK>9) fringe_dim(std::integral_constant<int,9>{});
}
__attribute__((target("avx2")))
__attribute__((aligned(1024))) static void solve_mf(u32 N, const unsigned **x, u32 *out, int k) {
++solve_generation;
const u32 W = ((N + 63) / 64 + 7u) & ~7u;
g_W = W;
#ifdef MF_SKIP_AND
g_skip_and = 1;
#endif
#ifdef MF_FILT
g_use_filt = 1;
#endif
#ifdef MF_NOFILT
g_use_filt = 0;
#endif
#ifdef MF_Q
g_Q = MF_Q; g_B = (u32)((N + (MF_Q) - 1) / (MF_Q));
#endif
#ifdef MF_FILT
g_use_filt = 1;
#endif
#ifdef MF_NOFILT
g_use_filt = 0;
#endif
#ifdef MF_SKIP_FRINGE
g_skip_fringe = 1;
#endif
if (!g_Q) {
u64 r = (u64)N, s2 = 0;
while (s2 * s2 < r) s2++;
u64 overhead = (u64)(5 * k + 21) * N * 4;
#if !MF_NOFILTB
overhead += (u64)k * k * 4 * N;
#endif
overhead += (u64)k * k * ((u64)N + 64);
overhead += (u64)N * ((k <= 8) ? 8u : (u32)k) * 4 * k;
u64 budget = ((u64)MF_BUDGETMB << 20);
if (budget > overhead + ((u64)64 << 20)) budget -= overhead; else budget = (u64)64 << 20;
#if MF_PLEN
static const u32 kfill[7] = {72, 67, 62, 58, 54, 50, 46};
u64 perfill = ((u64)(k - 1) * W * 8) * kfill[k - 4] / 100;
#else
u64 perfill = (u64)(k - 1) * W * 8;
#endif
u64 Qmem = budget / (perfill ? perfill : 1);
if (Qmem < 8) Qmem = 8;
static const u32 kscale[7] = {82, 91, 100, 108, 115, 122, 129};
u32 q = (u32)((s2 * MF_QCOEF + 99) / 100);
q = (u32)(((u64)q * kscale[k - 4] + 99) / 100);
if (k == 9 && N >= 900000) q = 1400;
if ((u64)q > Qmem) q = (u32)Qmem;
if (q > N) q = N;
if (q > (u32)(MAXQ - 1)) q = (u32)(MAXQ - 1);
g_Q = q;
g_B = (N + q - 1) / q;
}
const u32 Q = g_Q, B = g_B;
const u64 bmag = ((1ull << 40) + B - 1) / B;
#if MF_PLEN
g_plen_on = 1;
#endif
#if MF_WTRIM
u32 Wpool = W;
#else
const u32 Wpool = W;
#endif
#if MF_FILTAUTO
if (!g_use_filt && B >= 400) g_use_filt = 1;
#endif
if (g_use_filt && (u64)k * k * (g_Qc + 1) * W * 8 > ((u64)700 << 20)) g_use_filt = 0;
#if !MF_FILTAUTO
g_use_filt = 0;
#endif
{
int lo8[8], hi8[8];
for (int e = 0; e < 8; e++) lo8[e] = (e < k) ? 0 : 0x7FFFFFFF;
for (int e = 0; e < 8; e++) hi8[e] = (8 + e < k) ? 0 : 0x7FFFFFFF;
khi8v = LDU256(lo8);
khi8hv = LDU256(hi8);
}
const int K1 = k - 1;
for (int d = 0; d < k; d++) {
const unsigned *xd = x[d];
u32 *od = s_ord[d];
static u8 hist8[MAXN+1] __attribute__((aligned(64)));
for(u32 i=0;i<=N;i++)hist8[i]=0;
for(u32 i=0;i<N;i++){
if(i+64<N)_mm_prefetch((const char*)(xd+i+64),_MM_HINT_NTA);
u32 v=xd[i];hist8[v<N?v:0]++;
}
u32 sum=0;
u32 hv=0;
for(;hv+8<=N+1;hv+=8){
__m256i c=_mm256_cvtepu8_epi32(_mm_loadl_epi64((const __m128i*)(hist8+hv)));
__m256i p=_mm256_add_epi32(c,_mm256_slli_si256(c,4));
p=_mm256_add_epi32(p,_mm256_slli_si256(p,8));
__m128i lo=_mm256_castsi256_si128(p),hi=_mm256_extracti128_si256(p,1);
hi=_mm_add_epi32(hi,_mm_shuffle_epi32(lo,255));
__m256i inc=_mm256_inserti128_si256(_mm256_castsi128_si256(lo),hi,1);
__m256i ex=_mm256_add_epi32(_mm256_sub_epi32(inc,c),_mm256_set1_epi32(sum));
_mm256_storeu_si256((__m256i*)(s_tmp+hv),ex);
sum+=(u32)_mm_extract_epi32(hi,3);
}
for(;hv<=N;hv++){u32 c=hist8[hv];s_tmp[hv]=sum;sum+=c;}
bool compact=(sum==N);
if(!compact){
for(u32 i=0;i<=N;i++)s_cntv[i]=0;
for(u32 i=0;i<N;i++){u32 v=xd[i];s_cntv[v<N?v:0]++;}
sum=0;for(u32 i=0;i<=N;i++){u32 c=s_cntv[i];s_cntv[i]=sum;s_tmp[i]=sum;sum+=c;}
}
u32 *gt=s_gt[d],*pp=s_pos[d];
auto fill=[&](auto tag){
for(u32 i=0;i<N;i++){
if(i+MF_PFSC2<N){u32 vp=xd[i+MF_PFSC2];if(vp>=N)vp=0;_mm_prefetch((const char*)(s_tmp+vp),_MM_HINT_T0);}
if(i+64<N)_mm_prefetch((const char*)(xd+i+64),_MM_HINT_NTA);
u32 v=xd[i];if(v>=N)v=0;
u32 st,q,num;
if constexpr(decltype(tag)::value){u32 w=s_tmp[v];st=w&0xfffffu;q=st+(w>>20);s_tmp[v]=w+(1u<<20);num=(s_tmp[v+1]&0xfffffu)-1u;}
else{st=s_cntv[v];q=s_tmp[v]++;num=s_cntv[v+1]-1u;}
_mm_stream_si32((int*)(gt+i),(int)st);_mm_stream_si32((int*)(pp+i),(int)q);
#if MF_FB2
s_posp[(size_t)i*MAXK+d]=q;
#endif
if(num>N)num=N;
u32 bb=(u32)(((u64)num*bmag)>>40);
s_bdT[(size_t)d*E6Z_BDR+i]=bb>Q-1?Q-1:bb;
}
};
if(compact)fill(std::true_type{});else fill(std::false_type{});
for(u32 i=0;i<N;i++){
#if MF_PFSC2 > 0
if(i+64<N)_mm_prefetch((const char*)(pp+i+64),_MM_HINT_NTA);
if(i+MF_PFSC2<N)_mm_prefetch((const char*)(od+pp[i+MF_PFSC2]),_MM_HINT_T0);
#endif
od[pp[i]]=i;
}
BVAL(d, 0) = 0; FPOS(d, 0) = 0;
for (u32 b = 1; b < Q; b++) {
u32 pos = b * B < N ? b * B : N - 1;
u32 v = xd[od[pos]];
BVAL(d, b) = v;
u32 gs = pos;
while (gs > 0 && xd[od[gs - 1]] == v) gs--;
FPOS(d, b) = gs;
}
BVAL(d, Q) = 0xFFFFFFFFu; FPOS(d, Q) = 0;
}
_mm_sfence();
for (u32 i0 = 0; i0 < N; i0 += 64) {
u32 i1 = i0 + 64; if (i1 > N) i1 = N;
for (u32 i = i0; i < i1; i++) {
u16 *dst = s_bd + (size_t)i * MAXK;
for (int dd = 0; dd < k; dd++) dst[dd] = s_bdT[(size_t)dd * E6Z_BDR + i];
}
}
{ u32 nn_ = N - 1; u32 sh_ = 0; while (nn_ > 255u) { nn_ >>= 1; sh_++; } g_sigsh = sh_; }
const int b14af_fill_sig = (k == STRIDE) && ((size_t)k * k * MF_SPL <= ((size_t)220 << 20)) && ((size_t)N + 32 <= MF_SPL);
if (k == STRIDE) {
for (u32 i = 0; i < N; i++) {
u32 *p = s_pm + (size_t)i * PSTRIDE;
__m256i pv = _mm256_setr_epi32((int)x[0][i],(int)x[1][i],(int)x[2][i],(int)x[3][i],
(int)x[4][i],(int)x[5][i],(int)x[6][i],(int)x[7][i]);
_mm256_stream_si256((__m256i *)p, pv);
if ((int)PSTRIDE > (int)STRIDE) {
const u16 *bdi = s_bd + (size_t)i * MAXK;
const __m256i avv = _mm256_setr_epi32(
(int)BVAL(0,bdi[0]),(int)BVAL(1,bdi[1]),(int)BVAL(2,bdi[2]),(int)BVAL(3,bdi[3]),
(int)BVAL(4,bdi[4]),(int)BVAL(5,bdi[5]),(int)BVAL(6,bdi[6]),(int)BVAL(7,bdi[7]));
_mm256_stream_si256((__m256i *)(p + 8), avv);
}
if (b14af_fill_sig) { u8 *cs_ = b14af_csig + (size_t)i * 8;
for (int e = 0; e < 8; e++) cs_[e] = (u8)((x[e][i] >> g_sigsh) ^ 0x80u); }
}
_mm_sfence();
} else {
for (u32 i = 0; i < N; i++) {
u32 *p = s_pm + (size_t)i * PSTRIDE;
for (int d = 0; d < k; d++) p[d] = x[d][i];
for (int d = k; d < STRIDE; d++) p[d] = 0;
}
}
{
u32 nn = N - 1; u32 sh = 0;
while (nn > 255u) { nn >>= 1; sh++; }
g_sigsh = sh;
size_t need = (size_t)k * (k-1) * MF_SPL;
static u8 *sigpool = 0; static size_t sigcap = 0;
if (need > ((size_t)220 << 20) || (size_t)N + 32 > MF_SPL) { sigpool = 0; sigcap = 0; }
else if (need > sigcap) { sigpool = (u8 *)pool_alloc(need); sigcap = sigpool ? need : 0; }
if (sigpool) {
s_sig = sigpool;
g_sig_on = 1;
g_sig_merged = 0;
}
}
#if MF_NOFILTB
g_use_filt = 0;
#endif
#if MF_FB3
if (g_use_filt) {
u32 Qc = g_Qc;
u64 wc = (u64)N / Qc + 1;
for (int e = 0; e < k; e++) {
const unsigned *xe = x[e];
const u32 *pe = s_ord[e];
u32 ptr = 0;
u32 *bn = s_bnde + (size_t)e * (MAXQ + 2);
for (u32 c = 0; c <= Qc; c++) {
u64 v = (u64)c * wc;
if (v > 0xFFFFFFFFull) v = 0xFFFFFFFFull;
while (v && ptr < N && (u64)xe[pe[ptr]] < v) ptr++;
bn[c] = ptr;
}
}
}
#endif
if (g_use_filt) {
u32 Qc = g_Qc;
u64 wc = (u64)N / Qc + 1;
size_t need2 = (size_t)k * k * (Qc + 1) * W * sizeof(u64);
static u64 *fpool = 0; static size_t fcap = 0;
if (need2 > fcap) { fpool = (u64 *)pool_alloc(need2 + 4096); fcap = need2; }
if (fpool) {
s_fbm = fpool; g_fbm_W = W;
#if MF_FB2
for (int d = 0; d < k; d++)
for (int e = 0; e < k; e++) {
if (e == d) continue;
u64 *b0 = s_fbm + (size_t)(d * k + e) * (Qc + 1) * W;
for (u32 w = 0; w < W; w++) b0[w] = 0;
}
for (int e = 0; e < k; e++) {
const unsigned *xe = x[e];
const u32 *pe = s_ord[e];
int dmap[MAXK], km1 = 0;
for (int d = 0; d < k; d++) if (d != e) dmap[km1++] = d;
for (int i = 0; i < km1; i++) { u64 *a = s_facc + (size_t)i * W; for (u32 w = 0; w < W; w++) a[w] = 0; }
u32 ptr = 0;
for (u32 c = 1; c <= Qc; c++) {
u64 v = (u64)c * wc;
if (v > 0xFFFFFFFFull) v = 0xFFFFFFFFull;
while (ptr < N && (u64)xe[pe[ptr]] < v) {
const u32 *pv = s_posp + (size_t)pe[ptr] * MAXK;
for (int i = 0; i < km1; i++) {
u32 q = pv[dmap[i]];
s_facc[(size_t)i * W + (q >> 6)] |= 1ull << (q & 63);
}
ptr++;
}
for (int i = 0; i < km1; i++) {
int d = dmap[i];
u64 *dst = s_fbm + ((size_t)(d * k + e) * (Qc + 1) + c) * W;
const u64 *src = s_facc + (size_t)i * W;
#if MF_NT
{ u32 w = 0;
for (; w + 4 <= W; w += 4) _mm256_stream_si256((__m256i *)(dst + w), _mm256_load_si256((const __m256i *)(src + w)));
for (; w < W; w++) dst[w] = src[w]; }
#else
for (u32 w = 0; w < W; w++) dst[w] = src[w];
#endif
}
}
}
#else
for (int d = 0; d < k; d++) {
for (int e = 0; e < k; e++) {
if (e == d) continue;
const unsigned *xe = x[e];
const u32 *pe = s_ord[e];
const u32 *pd = s_pos[d];
u64 *base = s_fbm + (size_t)(d * k + e) * (Qc + 1) * W;
for (u32 w = 0; w < W; w++) base[w] = 0;
u64 *cur = s_acc;
for (u32 w = 0; w < W; w++) cur[w] = 0;
u32 ptr = 0;
for (u32 c = 1; c <= Qc; c++) {
#if MF_FB3
u32 pend = s_bnde[(size_t)e * (MAXQ + 2) + c];
#else
u64 v = (u64)c * wc;
if (v > 0xFFFFFFFFull) v = 0xFFFFFFFFull;
u32 pend = ptr;
{ const unsigned *xe_ = xe; const u32 *pe_ = pe;
while (pend < N && (u64)xe_[pe_[pend]] < v) pend++; }
#endif
while (ptr < pend) { u32 q = pd[pe[ptr]]; cur[q >> 6] |= 1ull << (q & 63); ptr++; }
u64 *dst = base + (size_t)c * W;
#if MF_NT
{ u32 w = 0;
for (; w + 4 <= W; w += 4) _mm256_stream_si256((__m256i *)(dst + w), _mm256_load_si256((const __m256i *)(cur + w)));
for (; w < W; w++) dst[w] = cur[w]; }
#else
for (u32 w = 0; w < W; w++) dst[w] = cur[w];
#endif
}
}
}
#endif
} else { g_use_filt = 0; }
_mm_sfence();
}
{
size_t mrn = (size_t)MAXK * MAXK * MAXQ;
for (size_t z = 0; z < mrn; z++) s_maxr[z] = 0;
}
{
static u32 y2_dsb[1024], y2_bst[1024];
for (u32 i0 = 0; i0 < N; i0 += 1024) {
u32 i1 = i0 + 1024; if (i1 > N) i1 = N;
for (u32 i = i0; i < i1; i++) {
u32 best = 0xFFFFFFFFu; int ds = 0;
for (int d = 0; d < k; d++) {
u32 g = s_gt[d][i];
if (g < best) { best = g; ds = d; }
}
y2_dsb[i - i0] = (u32)ds;
y2_bst[i - i0] = best;
s_g[i] = (u8)ds;
}
if (g_plen_on) {
for (u32 i = i0; i < i1; i++) {
u32 ds = y2_dsb[i - i0];
u32 bb = y2_bst[i - i0];
const u16 *bdv = s_bd + (size_t)i * MAXK;
u32 *mr = s_maxr + ((size_t)ds * MAXK) * MAXQ;
for (int d = 0; d < k; d++) {
if (d == ds) continue;
u32 b = bdv[d];
{ u32 *mp = &mr[(size_t)d * MAXQ + b]; u32 v = *mp; *mp = (bb > v) ? bb : v; }
}
}
}
}
}
{
for (u32 v = 0; v <= (u32)k; v++) s_cntv[v] = 0;
for (u32 i = 0; i < N; i++) s_cntv[s_g[i] + 1]++;
{ u32 ac = 0; for (u32 v = 0; v <= (u32)k; v++) { u32 c = s_cntv[v]; s_cntv[v] = ac; ac += c; } }
for (u32 i = 0; i < N; i++) { s_tmp[s_cntv[s_g[i] + 1]++] = i; }
for (u32 i = 0; i < N; i++) s_qord[i] = s_tmp[i];
}
u32 gstart[MAXK + 2];
{ u32 c = 0; for (int d = 0; d <= k; d++) { gstart[d] = c; while (c < N && (int)s_g[s_qord[c]] == d) c++; } gstart[k + 1] = c; }
for (int d = 0; d < k; d++) {
u32 pos[MAXK];
for (int ds = 0; ds < k; ds++) pos[ds] = gstart[ds];
u32 bcap_ = (N - 1) / B; if (bcap_ > Q - 1) bcap_ = Q - 1;
u32 bb_ = 0;
u32 rnext = (bcap_ >= 1u) ? FPOS(d, 1u) : 0xFFFFFFFFu;
for (u32 r = 0; r < N; r++) {
if (r >= rnext) { do { bb_++; rnext = (bb_ < bcap_) ? FPOS(d, bb_ + 1u) : 0xFFFFFFFFu; } while (r >= rnext); }
#if MF_FRPF > 0
if (r + MF_FRPF < N) _mm_prefetch((const char *)(s_g + s_ord[d][r + MF_FRPF]), _MM_HINT_T0);
#endif
u32 id = s_ord[d][r];
u32 dst = pos[s_g[id]]++;
s_frord[d][dst] = id | (bb_ << 20);
}
}
#if MF_NOSDMP
g_nosdmp = 1;
#endif
if (!g_nosdmp) {
u32 KS = k;
#if MF_KSPAD
if (k <= 8) KS = 8;
#endif
size_t per = (size_t)N * KS * 4;
size_t lim = (size_t)150 << 20;
if (k >= 9) lim = (size_t)470 << 20;
if (per * k <= lim) {
static u32 *pool = 0; static size_t cap = 0;
if (per * k > cap) { pool = (u32 *)pool_alloc(per * k + 256); cap = per * k; }
if (pool) {
g_KS = KS;
for (int d = 0; d < k; d++) {
s_sdmp[d] = pool + (size_t)d * N * g_KS;
const u32 *od = s_ord[d];
u32 *dst = s_sdmp[d];
const u32 KSx = g_KS;
for (u32 j = 0; j < N; j++) {
#if MF_PFB > 0
if (j + MF_PFB < N) _mm_prefetch((const char *)(s_pm + (size_t)od[j + MF_PFB] * PSTRIDE), _MM_HINT_T0);
#endif
const u32 *row = s_pm + (size_t)od[j] * PSTRIDE;
u32 *dp = dst + (size_t)j * KSx;
for (u32 e = 0; e < KSx; e++) dp[e] = (e < (u32)k) ? row[e] : 0u;
if (s_sig) {
u8 *sb = s_sig + (size_t)d * (k-1) * MF_SPL;
for(int e=0;e<k;e++)if(e!=d){u32 z_=row[e]>>g_sigsh;sb[(size_t)(e-(e>d))*MF_SPL+j]=(u8)(z_^0x80u);}
}
}
if (s_sig) g_sig_merged = 1;
}
}
}
}
if (s_sig && !g_sig_merged && k==8) {
const u8 *const cached_signature = b14af_csig;
for(int d=0;d<k;d++) {
const u32 *od=s_ord[d];
u8 *base=s_sig+(size_t)d*(k-1)*MF_SPL;
u32 j=0;
auto emit=[&](auto tag){constexpr int OMIT=decltype(tag)::value;
alignas(64) u8 pbuf[8][64]; u32 jc = 0, jg_ = 0;
for(;j+8<=N;j+=8) {
if (jc == 0) jg_ = j;
if(j+128<N)for(int e=0;e<8;e++)_mm_prefetch((const char *)(cached_signature+(size_t)od[j+128+e]*8),_MM_HINT_T0);
__m128i r0=_mm_loadl_epi64((const __m128i *)(cached_signature+(size_t)od[j+0]*8));
__m128i r1=_mm_loadl_epi64((const __m128i *)(cached_signature+(size_t)od[j+1]*8));
__m128i r2=_mm_loadl_epi64((const __m128i *)(cached_signature+(size_t)od[j+2]*8));
__m128i r3=_mm_loadl_epi64((const __m128i *)(cached_signature+(size_t)od[j+3]*8));
__m128i r4=_mm_loadl_epi64((const __m128i *)(cached_signature+(size_t)od[j+4]*8));
__m128i r5=_mm_loadl_epi64((const __m128i *)(cached_signature+(size_t)od[j+5]*8));
__m128i r6=_mm_loadl_epi64((const __m128i *)(cached_signature+(size_t)od[j+6]*8));
__m128i r7=_mm_loadl_epi64((const __m128i *)(cached_signature+(size_t)od[j+7]*8));
__m128i t0=_mm_unpacklo_epi8(r0,r1);
__m128i t1=_mm_unpacklo_epi8(r2,r3);
__m128i t2=_mm_unpacklo_epi8(r4,r5);
__m128i t3=_mm_unpacklo_epi8(r6,r7);
__m128i u0=_mm_unpacklo_epi16(t0,t1);
__m128i u1=_mm_unpackhi_epi16(t0,t1);
__m128i u2=_mm_unpacklo_epi16(t2,t3);
__m128i u3=_mm_unpackhi_epi16(t2,t3);
__m128i v0=_mm_unpacklo_epi32(u0,u2);
__m128i v1=_mm_unpackhi_epi32(u0,u2);
__m128i v2=_mm_unpacklo_epi32(u1,u3);
__m128i v3=_mm_unpackhi_epi32(u1,u3);
__m128i h0=_mm_unpackhi_epi8(r0,r1);
__m128i h1=_mm_unpackhi_epi8(r2,r3);
__m128i h2=_mm_unpackhi_epi8(r4,r5);
__m128i h3=_mm_unpackhi_epi8(r6,r7);
__m128i hu0=_mm_unpacklo_epi16(h0,h1),hu1=_mm_unpacklo_epi16(h2,h3);
__m128i v4=_mm_unpacklo_epi32(hu0,hu1);
if constexpr(OMIT!=0)_mm_storel_epi64((__m128i*)(pbuf[(0-(0>OMIT))]+jc*8),v0);
if constexpr(OMIT!=1)_mm_storel_epi64((__m128i*)(pbuf[(1-(1>OMIT))]+jc*8),_mm_srli_si128(v0,8));
if constexpr(OMIT!=2)_mm_storel_epi64((__m128i*)(pbuf[(2-(2>OMIT))]+jc*8),v1);
if constexpr(OMIT!=3)_mm_storel_epi64((__m128i*)(pbuf[(3-(3>OMIT))]+jc*8),_mm_srli_si128(v1,8));
if constexpr(OMIT!=4)_mm_storel_epi64((__m128i*)(pbuf[(4-(4>OMIT))]+jc*8),v2);
if constexpr(OMIT!=5)_mm_storel_epi64((__m128i*)(pbuf[(5-(5>OMIT))]+jc*8),_mm_srli_si128(v2,8));
if constexpr(OMIT!=6)_mm_storel_epi64((__m128i*)(pbuf[(6-(6>OMIT))]+jc*8),v3);
if constexpr(OMIT!=7)_mm_storel_epi64((__m128i*)(pbuf[(7-(7>OMIT))]+jc*8),_mm_srli_si128(v3,8));
if (++jc == 8) {
jc = 0;
#define R2J_FLUSH_PLANE(E,BUFIDX) if constexpr(OMIT!=E) { __m256i _a=LDU256(pbuf[(BUFIDX)]), _b=LDU256(pbuf[(BUFIDX)]+32); \
_mm256_stream_si256((__m256i*)(base+(size_t)(E-(E>OMIT))*MF_SPL+jg_),_a); \
_mm256_stream_si256((__m256i*)(base+(size_t)(E-(E>OMIT))*MF_SPL+jg_+32),_b); }
R2J_FLUSH_PLANE(0,0-(0>OMIT)) R2J_FLUSH_PLANE(1,1-(1>OMIT)) R2J_FLUSH_PLANE(2,2-(2>OMIT)) R2J_FLUSH_PLANE(3,3-(3>OMIT))
R2J_FLUSH_PLANE(4,4-(4>OMIT)) R2J_FLUSH_PLANE(5,5-(5>OMIT)) R2J_FLUSH_PLANE(6,6-(6>OMIT)) R2J_FLUSH_PLANE(7,7-(7>OMIT))
#undef R2J_FLUSH_PLANE
}
}
for (u32 t2 = 0; t2 < jc; t2++) {
#define R2J_RES_PLANE(E,BUFIDX) if constexpr(OMIT!=E) _mm_storel_epi64((__m128i*)(base+(size_t)(E-(E>OMIT))*MF_SPL+(j-8*(jc-t2))+0),_mm_loadl_epi64((const __m128i*)(pbuf[(BUFIDX)]+t2*8)));
R2J_RES_PLANE(0,0-(0>OMIT)) R2J_RES_PLANE(1,1-(1>OMIT)) R2J_RES_PLANE(2,2-(2>OMIT)) R2J_RES_PLANE(3,3-(3>OMIT))
R2J_RES_PLANE(4,4-(4>OMIT)) R2J_RES_PLANE(5,5-(5>OMIT)) R2J_RES_PLANE(6,6-(6>OMIT)) R2J_RES_PLANE(7,7-(7>OMIT))
#undef R2J_RES_PLANE
}
};
switch(d){case 0:emit(std::integral_constant<int,0>{});break;case 1:emit(std::integral_constant<int,1>{});break;case 2:emit(std::integral_constant<int,2>{});break;case 3:emit(std::integral_constant<int,3>{});break;case 4:emit(std::integral_constant<int,4>{});break;case 5:emit(std::integral_constant<int,5>{});break;case 6:emit(std::integral_constant<int,6>{});break;case 7:emit(std::integral_constant<int,7>{});break;}
for(;j<N;j++) {
const u8 *row=cached_signature+(size_t)od[j]*8;
for(int e=0;e<k;e++)if(e!=d)base[(size_t)(e-(e>d))*MF_SPL+j]=row[e];
}
}
_mm_sfence();
g_sig_merged=1;
}
if (s_sig && !g_sig_merged) {
for (int d = 0; d < k; d++) {
const u32 *od = s_ord[d];
u8 *base = s_sig + (size_t)d * (k-1) * MF_SPL;
for (u32 j = 0; j < N; j++) {
#if MF_PFB > 0
if (j + MF_PFB < N) _mm_prefetch((const char *)(s_pm + (size_t)od[j + MF_PFB] * PSTRIDE), _MM_HINT_T0);
#endif
const u32 *row = s_pm + (size_t)od[j] * PSTRIDE;
for(int e=0;e<k;e++)if(e!=d){u32 z_=row[e]>>g_sigsh;base[(size_t)(e-(e>d))*MF_SPL+j]=(u8)(z_^0x80u);}
}
}
}
u64 poolwords = (u64)K1 * (Q + 1) * W;
#if MF_PLEN
if (g_plen_on) {
u64 best = 0;
for (int dsx = 0; dsx < k; dsx++) {
u64 tot = 0;
for (int e = 0; e < K1; e++) {
int d = (e < dsx) ? e : e + 1;
const u32 *mr = s_maxr + ((size_t)dsx * MAXK + d) * MAXQ;
for (u32 b = 0; b <= Q; b++) {
u32 RR = mr[b];
u32 len = ((RR >> 6) + 8u) & ~7u;
if (len > W) len = W;
tot += len;
}
}
if (tot > best) best = tot;
}
if (best && best < poolwords) poolwords = best;
}
#endif
size_t need = (size_t)poolwords * sizeof(u64);
static u64 *pool = 0; static size_t cap2 = 0;
if (need > cap2) { pool = (u64 *)pool_alloc(need + 4096); cap2 = need; }
s_bs = pool;
#if MF_INCA
{
size_t needho = (size_t)N * sizeof(u32) + 4096;
if (needho > s_ho_cap) { s_ho = (u32 *)pool_alloc(needho); s_ho_cap = s_ho ? needho : 0; }
g_inca = (s_ho && (u64)N <= (u64)AS_WORDS * 64ull) ? 1 : 0;
}
#endif
for (int ds = 0; ds < k; ds++) {
u32 q0 = gstart[ds], q1 = gstart[ds + 1];
if (q0 >= q1) continue;
#if MF_WTRIM
{
u32 mx = 1;
for (u32 qi = q0; qi < q1; qi++) { u32 rr = s_gt[ds][s_qord[qi]]; if (rr > mx) mx = rr; }
u32 wp = ((mx >> 6) + 8u) & ~7u;
Wpool = (wp < W) ? wp : W;
}
#endif
#if MF_PLEN
if (g_plen_on) {
u64 tot = 0;
for (int e = 0; e < K1; e++) {
int d = (e < ds) ? e : e + 1;
const u32 *mr = s_maxr + ((size_t)ds * MAXK + d) * MAXQ;
u32 *po = s_pofs + (size_t)e * MAXQ;
u32 *pl = s_plen + (size_t)e * MAXQ;
for (u32 b = 0; b <= Q; b++) {
u32 RR = mr[b];
u32 len = ((RR >> 6) + 8u) & ~7u;
if (len > Wpool) len = Wpool;
po[b] = (u32)tot; pl[b] = len; tot += len;
}
}
if (tot > (u64)K1 * (Q + 1) * W) { g_plen_on = 0; }
}
#endif
const u32 *const spp = s_pos[ds];
for (int e = 0; e < K1; e++) {
int d = (e < ds) ? e : e + 1;
const unsigned *xd = x[d];
const u32 *od = s_ord[d];
#if MF_INCA
if (g_inca && d == ((ds == 0) ? 1 : 0)) {
u32 *dstho = s_ho;
const u32 *hp = spp;
for (u32 j = 0; j < N; j++) {
if (j + 64 < N) _mm_prefetch((const char *)(od + j + 64), _MM_HINT_NTA);
_mm_stream_si32((int *)&dstho[j], (int)hp[od[j]]); }
_mm_sfence();
continue;
}
#endif
u64 *base = g_plen_on ? s_bs : (s_bs + (size_t)e * (Q + 1) * Wpool);
u64 *cur = s_acc;
for (u32 w = 0; w < Wpool; w++) cur[w] = 0;
if (g_plen_on) {
u64 *d0 = s_bs + s_pofs[(size_t)e * MAXQ];
u32 l0 = s_plen[(size_t)e * MAXQ];
for (u32 w = 0; w < l0; w++) d0[w] = 0;
}
u32 ptr = 0;
#if MF_BUILD_NOPOOL
while (ptr < N) ptr++;
#else
for (u32 b = 1; b <= Q; b++) {
u32 pend = (b < Q) ? FPOS(d, b) : (u32)N;
while (ptr < pend) {
u32 nx = ptr + 32u; if (nx >= N) nx = ptr;
_mm_prefetch((const char *)&spp[od[nx]], _MM_HINT_T0);
_mm_prefetch((const char *)(od + ptr + 64), _MM_HINT_NTA);
u32 q = spp[od[ptr]]; cur[q >> 6] |= 1ull << (q & 63); ptr++; }
u64 *dst = g_plen_on ? (base + s_pofs[(size_t)e * MAXQ + b])
: (base + (size_t)b * Wpool);
u32 wlen = g_plen_on ? s_plen[(size_t)e * MAXQ + b] : Wpool;
#if MF_NT
{
u32 w = 0;
for (; w + 4 <= wlen; w += 4)
_mm256_stream_si256((__m256i *)(dst + w), _mm256_load_si256((const __m256i *)(cur + w)));
for (; w < wlen; w++) dst[w] = cur[w];
}
#else
for (u32 w = 0; w < wlen; w++) dst[w] = cur[w];
#endif
}
_mm_sfence();
#endif
}
tie_count=0;
query_group<8>(N, out, ds, q0, q1, Q, Wpool);
flush_ties(out);
}
}
static void hud_pages(void *p, size_t n) {
if (!p || !n) return;
size_t a = (size_t)p & ~(size_t)4095;
madvise((void *)a, (size_t)(((size_t)p + n + 4095) & ~(size_t)4095) - a, 14 );
}
static void hud_all(void) {
char *lo = (char *)&s_ord, *hi = (char *)&s_ord;
char *c[17] = { (char *)&s_ord, (char *)&s_frord, (char *)&s_gt, (char *)&s_pm, (char *)&s_bdT,
(char *)&s_bd, (char *)&s_g, (char *)&s_qord, (char *)&s_tmp, (char *)&s_cntv,
(char *)&s_bkt, (char *)&s_bks, (char *)&s_nid, (char *)&s_pos, (char *)&s_posp,
(char *)&s_facc, (char *)&s_bnde };
for (int i = 0; i < 17; i++) { if (c[i] < lo) lo = c[i]; if (c[i] > hi) hi = c[i]; }
hi += 8u << 20;
madvise((void *)((size_t)lo & ~(size_t)4095), (size_t)((size_t)hi - ((size_t)lo & ~(size_t)4095)), 14);
hud_pages(s_fbm, (size_t)g_fbm_W * 8 * 1024);
hud_pages(s_ho, s_ho_cap * 4);
}
void count_8d(int n, const unsigned *x[8], unsigned *out) { hud_all(); solve_mf((u32)n, (const unsigned **)x, out, 8); }
#ifdef LOCAL_TEST
static void brute_nd(int n, const u32 **x, u32 *out, int k) {
for (int i = 0; i < n; i++) {
u32 c = 0;
for (int j = 0; j < n; j++) {
int ok = 1;
for (int d = 0; d < k; d++) if (!(x[d][j] < x[d][i])) { ok = 0; break; }
c += ok;
}
out[i] = c;
}
}
int main(int argc, char **argv) {
int k = (argc > 1) ? atoi(argv[1]) : 6;
int n = (argc > 2) ? atoi(argv[2]) : 200;
int seed = (argc > 3) ? atoi(argv[3]) : 1;
int mode = (argc > 4) ? atoi(argv[4]) : 0;
int qq = (argc > 5 && atoi(argv[5]) > 0) ? atoi(argv[5]) : 0;
g_Q = qq ? (u32)qq : (u32)(n / 250 + 1);
g_B = (u32)((n + g_Q - 1) / g_Q);
if (argc > 6) g_nosort = atoi(argv[6]);
if (argc > 7) g_nosdmp = atoi(argv[7]);
if (argc > 8) g_skip_and = atoi(argv[8]);
if (argc > 9) g_skip_fringe = atoi(argv[9]);
if (argc > 10) g_use_filt = atoi(argv[10]);
if (argc > 11) g_lexsort = atoi(argv[11]);
const int dm = (mode >= 10) ? mode - 10 : mode;
srand(seed);
static u32 *xs[10];
u32 *out = (u32 *)malloc(sizeof(u32) * n);
u32 *ref = (u32 *)malloc(sizeof(u32) * n);
for (int d = 0; d < k; d++) {
xs[d] = (u32 *)malloc(sizeof(u32) * n);
for (int i = 0; i < n; i++) {
if (dm == 0) xs[d][i] = (u32)(rand() % n);
else if (dm == 1) xs[d][i] = (u32)(rand() % 3);
else if (dm == 2) xs[d][i] = 0;
else if (dm == 3) xs[d][i] = (u32)(n - 1 - (rand() % n));
else if (dm == 4) xs[d][i] = (u32)(i);
else xs[d][i] = (u32)(rand() % n);
}
}
if (mode >= 10) {
struct timespec t0, t1;
clock_gettime(CLOCK_MONOTONIC, &t0);
switch (k) {
case 4: count_4d(n, (const unsigned **)xs, out); break;
case 5: count_5d(n, (const unsigned **)xs, out); break;
case 6: count_6d(n, (const unsigned **)xs, out); break;
case 7: count_7d(n, (const unsigned **)xs, out); break;
case 8: count_8d(n, (const unsigned **)xs, out); break;
case 9: count_9d(n, (const unsigned **)xs, out); break;
default: count_10d(n, (const unsigned **)xs, out); break;
}
clock_gettime(CLOCK_MONOTONIC, &t1);
double ms = (t1.tv_sec - t0.tv_sec) * 1e3 + (t1.tv_nsec - t0.tv_nsec) / 1e6;
u64 s = 0; for (int i = 0; i < n; i++) s += out[i];
printf("k=%d n=%d Q=%u TIME %.1f ms checksum %llu\n", k, n, g_Q, ms, s);
return 0;
}
switch (k) {
case 4: count_4d(n, (const unsigned **)xs, out); break;
case 5: count_5d(n, (const unsigned **)xs, out); break;
case 6: count_6d(n, (const unsigned **)xs, out); break;
case 7: count_7d(n, (const unsigned **)xs, out); break;
case 8: count_8d(n, (const unsigned **)xs, out); break;
case 9: count_9d(n, (const unsigned **)xs, out); break;
default: count_10d(n, (const unsigned **)xs, out); break;
}
brute_nd(n, (const u32 **)xs, ref, k);
int bad = 0;
for (int i = 0; i < n; i++) if (out[i] != ref[i]) { if (bad < 5) printf("MISMATCH i=%d got=%u ref=%u\n", i, out[i], ref[i]); bad++; }
printf("k=%d n=%d mode=%d Q=%u %s (%d mismatches)\n", k, n, mode, g_Q, bad ? "FAIL" : "OK", bad);
return bad != 0;
}
#endif
| Compilation | N/A | N/A | Compile OK | Score: N/A | 显示更多 |
| Testcase #1 | 4.64 s | 720 MB + 612 KB | Accepted | Score: 100 | 显示更多 |