/* References:
* saffah_cc_v41_agg1 (this account), https://duck.ac/submission/121607:
* copied verbatim as the base -- the previous board leader and this file's parent:
* the inherited KG/CDQ engine with the [SORT2] two-element leaf-bucket ordering;
* every earlier citation carried by that file is kept unchanged below.
* saffah_codex_6a_agg3, https://duck.ac/submission/121774:
* the three- and four-event compare-exchange networks used below for buckets that
* hold three or four events ([CE34]) are taken from that submission's leaf sort; it
* pairs them with a compacted list of multi-event buckets, whereas here they sit in
* the mask-driven visitor instead and no such list is built. Public sources
* displayed no separate software license notices.
* Idea [SORTB]: each leaf ends by ordering its own run with a counting sort whose key
* is E_KA(v) >> (LFSH + LFD); the pass that repaired the intra-bucket order used to
* walk all LB bucket slots directly, and at n = 1e6 LB = G_LIM = 3908 while a leaf
* holds at most KG_BASE = 3072 events, so over 85% of the slots are empty or hold a
* single event and the `if (hi - lo > 1)` guard sitting on them is a data-dependent
* branch that flips unpredictably. The per-bucket sizes cnt[b] - cnt[b-1] are now
* formed with one 256-bit subtract off a padded array (cnt[-1] == 0 makes the load
* legal), compared against 1 and folded into a 32-bit mask four vector pairs at a
* time, so the discovery loop takes no data-dependent branch and only the buckets
* that really hold two or more events are ever visited. A bucket with exactly two
* events is ordered by the [SORT2] branch-free compare-exchange, three and four by
* the [CE34] networks, and larger ones by insertion sort as before.
* [LDU2] the 32-byte loads of that pass are written as one inline-asm `vmovdqu`,
* because gcc 9.3 splits a 256-bit load whose alignment it cannot prove into
* vmovdqu xmm + vinserti128 (2 uops).
* Bit-identical by construction: the bucket assignment, the prefix sums, the scatter
* and every repair are the same operations on the same ranges, in the same bucket
* order, so the leaf's sorted run -- and therefore the whole merge tree -- is unchanged.
* Purpose: Official experiment removing the last data-dependent branch from the leaf
* bucket pass and making every small-bucket repair branch-free; no change to the
* engine, to any key, or to the merge order.
*/
/* References:
* saffah_codex_6a_agg3, https://duck.ac/submission/120748:
* copied verbatim as the base (the previous board leader: compile-time NEED
* specialization, sentinel merge, AVX2 cross re-key, prefix-byte leaf counters)
* and its inherited citations; the earlier citations in this file are kept.
* saffah_cc_v41_agg1, https://duck.ac/submission/121023 (this account):
* the immediately preceding best file = #120748 verbatim + [CLR0] + [PFCT].
* Public sources displayed no separate software license notices.
* Idea [SORT2]: the per-leaf bucket pass that orders a leaf's events for its parent's
* merge ends with an insertion sort whose inner `while (q >= lo && v < out[q])` is a
* data-dependent branch; with LEAFNB 4096 the buckets hold ~0.5 events on average, so
* the great majority of multi-event buckets hold EXACTLY two and the mispredict rate of
* that branch is ~50%. Special-case hi-lo == 2 with a straight two-way cmov swap
* (a < b ? ... : ...), which is order-identical and equally stable. Local phase
* instrumentation (rdtsc, per-depth, n=1e6 mode10) attributes 5.6% of the whole run to
* that insertion-sort phase; this change takes it 5.44e8 -> 4.77e8 ticks (-12.3%, two
* identical runs, both directions reproduced).
* Purpose: Official experiment removing branch mispredictions from the leaf ordering
* pass. Bit-identical to #121023/#120748: n=1e6 mode10 checksum 15581398097 / FNV
* 12542427962008566018 unchanged; brute force exact for k=6 n=1000/3000/20000/40000
* x mode 0/1/3 and k=4/5 n=5000/60000; judge gcc 9.3 -O2 -std=c++17 compiles clean and
* its asm keeps vinserti128=80 / vmovdqu-xmm=82 / vpermq=0 exactly as the base.
*/
/* References:
* saffah_codex_6a_agg3, https://duck.ac/submission/120748:
* copied verbatim as the base -- the current board leader: compile-time NEED
* specialization of solve_t/solve_b, sentinel merge loop, AVX2 vectorised cross
* re-key, prefix-byte leaf counters and its inherited citations.
* saffah_cc_v41_agg1, https://duck.ac/submission/120565 (and the chain it cites):
* inherited engine and the external-leader lineage carried by that submission.
* Public sources displayed no separate software license notices.
* Idea: [CLR0] ctr_solve's dirt-list clear only has to bring the 32-byte CT0 leaf
* block back to zero -- leaf_add is exactly undone and the counter must end empty --
* so two 16-byte zero stores replace load+subtract+store (written as two 16-byte
* stores because gcc 9.3, the judge's compiler, splits a 256-bit store intrinsic
* into vmovdqu xmm + vinserti128). [PFCT] drop the two prefetches aimed at CT1/CT2
* in the same forward loop: those arrays are 64 KB / 4 KB and stay L1-resident, only
* CT0 (1 MB) needs the lead.
* Purpose: Official experiment removing per-insert counter-clear work and pure-uop
* prefetches in the fused innermost-dimension counter. Verified bit-identical to
* #120748 (n=1e6 mode10 checksum 15581398097 / FNV 12542427962008566018, plus seeds
* 1/2/3 at n=300000) and brute-force exact for k=4/5/6, n=1000..40000, modes 0/1/3.
*/
/* References:
* saffah_codex_6a_agg3, https://duck.ac/submission/120649:
* retained the accepted sentinel/vector-rekey engine and inherited citations.
* saffah_cc_v41_agg1, https://duck.ac/submission/120534:
* adapted making the output-required property a compile-time recursion parameter.
* Public sources displayed no separate software license notices.
* Idea: Specialize recursive merge functions for whether their parent needs a
* sorted output run. This removes per-event runtime condition checks and lets
* the compiler optimize distinct write-producing and count-only bodies.
* Purpose: Official experiment reducing invariant branches in recursive merges.
*/
/* References:
* saffah_codex_6a_agg3, https://duck.ac/submission/120630:
* retained the vectorized re-key engine and inherited citations.
* saffah_cc_v41_agg1, https://duck.ac/submission/108101:
* adapted its sentinel-based single-boundary merge-loop technique.
* Public source displayed no separate software license notice.
* Idea: Reserve one extra left-run slot containing an infinite sentinel.
* The main merge can then run to the right endpoint with no left bound check,
* absorbing the former right tail. Left-tail processing and event order stay exact.
* Purpose: Official experiment removing one merge-loop boundary dependency.
*/
/* References:
* saffah_codex_6a_agg3, https://duck.ac/submission/120600:
* retained its accepted engine with all inherited citations.
* Public sources displayed no separate software license notices.
* Idea: Re-key four cross events at a time with AVX2 gathers and packed 64-bit
* masks/shifts, retaining the same ahead-of-use coordinate prefetch and scalar tail.
* Purpose: Official experiment reducing hot-loop instruction and memory work.
*/
/* References:
* saffah_cc_v41_agg1, https://duck.ac/submission/120565:
* retained the current external leader and its inherited citations.
* saffah_codex_6a_agg3, https://duck.ac/submission/120558:
* ported guarded byte-prefix counters and the co-located upper-level storage.
* Public source displayed no separate license notice.
* Idea: Use a scalar prefix-byte lookup instead of a horizontal leaf reduction.
* Keep the safe per-block overflow guard and one-base counter layout, preserving
* the leader's separately prefetched event re-keying pass.
* Purpose: Official experiment transferring improvements from related dimensions.
*/
// ===== REFERENCES =====
// [1] duck.ac 用户 saffah_codex_6a_agg3, 提交 #120433 <https://duck.ac/submission/120433>
// 用途：**逐字复制**该提交的全部实现 —— 本题现役最快的 KG / CDQ 引擎：事件序列分治
// （u64 内联域 [ID:20][KB:20][TY:1][KA:20]）、位置二分树递归到 KG_BASE=3072 的
// 4-查询组 AVX2 叶（MF_LEAFG=2 / MF_LEAFLP=128）、D=KT-3 处融合的四层前缀计数器
// [10ez6-L4]（增量掩码 SFX16/SFX32/SFX32W + 批量清零 ctr_reset_hi）。抄的是整个
// 文件；本发只在其上做下述三处改动（同一族：cross 重定键的搬运方式）。
// （该件自身头段声明其来源：CDQ 引擎承自 saffah_cc_v41_260924 #117669
// <https://duck.ac/submission/117669>；四层计数器/增量掩码/批量清零承自本账号
// saffah_cc_v41_agg1 #112070 <https://duck.ac/submission/112070> —— 一并转引。）
// [2] 本账号, 提交 #117444 <https://duck.ac/submission/117444>（`x12aa_sub_c1.cpp`）
// 用途：本发改动①②的出处 —— [C1]（cross 重定键外提为一趟独立扫描）是该件的唯一
// 变量；该件与 #117426（= #117412 的逐字复刻，同批锚）之差即 [C1] 的判题机定价：
// 3974.468325 → 3682.330464 ms = **−292.138 ms（−7.35%）**。
// [3] 本账号, 提交 #120513 <https://duck.ac/submission/120513>
// 用途：如实登记的**负结果**（本发的前提）：把 [C1] 单独叠到对手件 #120433 上
// = 3652.958984 ms，比对手原件的 3598.223054 **慢 +54.727930 ms**。⇒ 该负号的
// 机理分析见下（合并环 KG_PF 前瞻在 [C1] 下变成纯浪费），改动②③即由此而来。
// ======================
// ===== 思路 =====
// 承接 #120513 的负结果：对手件上 [C1] 之所以由正转负，是因为它的合并环里有一对
// `__builtin_prefetch(&XN[E_ID(retA[ia+64])])`（KG_PF=64）——在原设计里它覆盖的是
// **合并环自己**那次 4 MB 坐标行 gather；[C1] 把 gather 搬走之后，这两条前瞻仍然照发，
// 但被前瞻的行要等到**整个合并环跑完**、第二趟扫描时才被读，届时早已被挤出缓存
// ⇒ 前瞻变成纯浪费（发 uop + 拉一轮注定被淘汰的 DRAM 流量），把 [C1] 的收益吃掉还倒亏。
//
// 本发三处改动（自包含、彼此同族，均只动"重定键这一趟"的取数方式）：
// ① [C1] 合并环（solve_b 主循环）只写原始事件 `cross[cx] = v`；两个尾巴环同理
// （不再逐事件算重定键表达式、不再在环内发 4 MB 随机 gather）。
// ② [PF0] `KG_PF 64 → 0`：删掉合并环里那两条已成空转的前瞻（**保留它们必负**，见 (3)）。
// ③ [PFR64] 在第二趟（重定键扫描）里发前瞻：`if (j + 64 < cx)
// __builtin_prefetch(&XN[E_ID(cross[j+64])], 0, 1);` —— 这一趟是纯 gather 循环，
// 前瞻距离 64 恰好覆盖它自己的取数（对手在合并环用的也是 64，同量级已验证）。
// 逐位等价：① 施用的表达式与顺序与原来**逐字相同**（只是延后到合并完成之后）；
// ②③ 只是前瞻提示（不改语义）。⇒ 输出与对手件逐位相同（实测见下）。
//
// 闸门（本机 g++ 15.2 -O2 -std=c++17；驱动 = 题内长期沿用的 base.cpp 数据生成器）：
// * 与对手件 #120433 **逐位一致**：n=1e6 × seed=1/2/3 × mode10（判题机真实工作点尺度）
// checksum + FNV hash 三对全同（seed1 = 15581398097，与生产格在册值相同）；
// * 暴力对拍：k=4/5/6 × n=1000/3000 × mode=0..4 = 30/30 OK；k=6 × n=20000/40000 ×
// mode=0/1/3 = 6/6 OK（n=40000 已覆盖到 D=3 层的合并与融合计数器）；
// * 正控（#120513 那一轮做的，同一趟代码）：把第二趟的 KB 置 0 ⇒ 小格 17635 mismatches、
// 生产 checksum 15247643 ≠ 15581398097 ⇒ 闸门非空转、第二趟确实在生成码执行路径上；
// * 本地 A/B（同 harness、交错、taskset 钉核、min-of-5）：五个 arm 五轮**全轮同号**
// `对手件 6152.2 ↔ 本件(pf64) 5341.3 ms`（−13.2%）；前瞻距离 16/32/64/96/128 的
// 差别在本机噪声内（本机 L3 远大于判题机 6 MB ⇒ 幅度不可外推，只作符号判据）。
// ================
/* References:
* saffah_cc_v41_260924, https://duck.ac/submission/117669:
* retained its complete CDQ engine and inherited citations.
* saffah_cc_v41_agg1, https://duck.ac/submission/112070:
* ported the four-level cumulative counter, increment masks and bulk reset.
* No separate software license notices appeared on the public submissions.
* Idea: Transfer the accepted four-upper-level hierarchy into this dimension/size
* engine, adjusting counter prefetch indices and reset bounds consistently.
* Queries use one fewer upper-level load; inserts use one fewer update.
* Purpose: Official experiment of a shared counter improvement across related tasks.
*/
// ARMVAR KG_BASE=2048 MF_CTRPF=32
#define KG_BASE 3072
#define MF_CTRPF 32
#ifndef MF_CTRPF
#define MF_CTRPF 0
#endif
#ifndef MF_CTR5
#define MF_CTR5 1
#endif
#ifndef MF_LEAFLP
#define MF_LEAFLP 128
#endif
#ifndef MF_LEAFG
#define MF_LEAFG 2
#endif
#pragma GCC optimize("O3","unroll-loops")
#pragma GCC target("avx2,bmi,bmi2,popcnt,lzcnt")
// duck.ac 1010 (4D) / 1011 (5D) -- "count for every point the number of points that are
// smaller in *every* dimension" (n = 300000, coords in [0,n), strict comparisons).
//
// Algorithm: D&C over dimensions on an event sequence, keys carried inline.
// * 2n events: one query event + one insert event per point, keyed by dim 0
// (a query event is stored *before* an insert event of the same key, which encodes
// the strict-inequality rule).
// * solve(seq,m,d): seq is ordered by dim d and every event carries the coordinate of
// dim d+1 in bits 32..63 of its u64. Split by position, recurse both halves, then
// merge the two halves by the carried key (a plain unsigned `<` on the packed u64,
// so equal keys keep query-before-insert order); the same pass collects the cross
// events (left inserts + right queries), re-keys them to dim d+2 and recurses one
// dimension deeper. Ordering by dim d with that tie rule makes "insert before query"
// exactly "all dims so far strictly smaller", and equal coordinates are never counted.
// * d == K-1: every insert-before-query pair now qualifies, so one O(m) running counter
// finishes the innermost dimension (no log factor).
// * m <= 32: direct O(m^2) counting plus an insertion sort.
//
// Local timings (i5-6500 @ 3.2GHz, best of 3, n = 300000 uniform random):
// 4D 1.57 s 5D 3.85 s
// (Board #1: 508.596 ms / 745.358 ms -- this is still 3x / 5x behind.)
// k-dimensional dominance counting: out[i] = #{j : all dims of j < dims of i}
//
// Event-sequence CDQ with keys carried inline.
// event (u64) = (key << 32) | (ins << 31) | id ins = 1 for an "insert" event
// key = coordinate of dim (level+1): the field used by this level's merge.
// The packing makes a plain unsigned comparison put, for equal keys, the query
// event first -- exactly the strict-inequality rule (equal coords must not count).
//
// solve(seq, m, d): seq is ordered by dim d (tie rule above), elements carry the
// key for dim d+1; count pairs (insert before query) with dims d+1..K-1 strictly
// smaller; return `out` = the same events merged by the carried key (= dim d+1).
// Cross = A-inserts + B-queries re-keyed to dim d+2, recursed at level d+1.
// d == K-1 : every (insert before query) pair qualifies -> running counter, O(m).
#include <cstring>
#include <immintrin.h>
#pragma GCC push_options
#pragma GCC target("avx2,popcnt,bmi2")
#include <cstdint>
typedef unsigned u32;
typedef uint64_t u64;
typedef unsigned char u8;
typedef unsigned short u16;
#ifndef KG_COND
#define KG_COND 0
#endif
#ifndef KG_PF
#define KG_PF 0
#endif
#ifndef KG_BASE
#define KG_BASE 4096
#endif
#ifndef LEAFNB
#define LEAFNB 4096
#endif
#define KG_MAXK 6
#define KG_ARENA (160u << 20)
namespace kg {
// ------------ fused innermost-dimension counter (prefix count over [0,n)) -------------
// AVX2. Levels: u8 leaf /32 values, u16 per 32 and per 256, u32 per 2048/16384/131072.
// Every u16 block count must stay <= 32767 (madd uses signed lanes); the guard (max
// multiplicity of a last-dim coordinate <= 127) bounds 32*127=4064 and 256*127=32512.
#define CT_MAXV 1048580u
struct alignas(64) PrefixCounterStorage {
alignas(64) u8 CT0[CT_MAXV];
alignas(64) u16 CT1[CT_MAXV / 32 + 64];
alignas(64) u16 CT2[CT_MAXV / 512 + 64];
alignas(64) u32 CT3[CT_MAXV / 8192 + 64];
alignas(64) u32 CT4[CT_MAXV / 131072 + 16];
};
static PrefixCounterStorage prefixct;
#define CT0 (prefixct.CT0)
#define CT1 (prefixct.CT1)
#define CT2 (prefixct.CT2)
#define CT3 (prefixct.CT3)
#define CT4 (prefixct.CT4)
/* [10ez6-L4] four upper layers instead of five: bits 5-8 / 9-12 / 13-16 / 17-19 (leaf = bits 0-4).
Same suffix-ladder semantics, one fewer scalar load per query and one fewer RMW per insert. */
// [SFX] suffix-counter increment masks. Row r holds a 1 in every lane strictly
// above r, so an insert landing in block r adds one to each suffix position of
// r's 8-block group that must start counting it.
static const u16 INC16[16][16] __attribute__((aligned(64))) = { {0,1,1,1,1,1,1,1,1,1,1,1,1,1,1,1},{0,0,1,1,1,1,1,1,1,1,1,1,1,1,1,1},{0,0,0,1,1,1,1,1,1,1,1,1,1,1,1,1},{0,0,0,0,1,1,1,1,1,1,1,1,1,1,1,1},{0,0,0,0,0,1,1,1,1,1,1,1,1,1,1,1},{0,0,0,0,0,0,1,1,1,1,1,1,1,1,1,1},{0,0,0,0,0,0,0,1,1,1,1,1,1,1,1,1},{0,0,0,0,0,0,0,0,1,1,1,1,1,1,1,1},{0,0,0,0,0,0,0,0,0,1,1,1,1,1,1,1},{0,0,0,0,0,0,0,0,0,0,1,1,1,1,1,1},{0,0,0,0,0,0,0,0,0,0,0,1,1,1,1,1},{0,0,0,0,0,0,0,0,0,0,0,0,1,1,1,1},{0,0,0,0,0,0,0,0,0,0,0,0,0,1,1,1},{0,0,0,0,0,0,0,0,0,0,0,0,0,0,1,1},{0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,1},{0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0} };
static const u32 INC32W[16][16] __attribute__((aligned(64))) = { {0,1,1,1,1,1,1,1,1,1,1,1,1,1,1,1},{0,0,1,1,1,1,1,1,1,1,1,1,1,1,1,1},{0,0,0,1,1,1,1,1,1,1,1,1,1,1,1,1},{0,0,0,0,1,1,1,1,1,1,1,1,1,1,1,1},{0,0,0,0,0,1,1,1,1,1,1,1,1,1,1,1},{0,0,0,0,0,0,1,1,1,1,1,1,1,1,1,1},{0,0,0,0,0,0,0,1,1,1,1,1,1,1,1,1},{0,0,0,0,0,0,0,0,1,1,1,1,1,1,1,1},{0,0,0,0,0,0,0,0,0,1,1,1,1,1,1,1},{0,0,0,0,0,0,0,0,0,0,1,1,1,1,1,1},{0,0,0,0,0,0,0,0,0,0,0,1,1,1,1,1},{0,0,0,0,0,0,0,0,0,0,0,0,1,1,1,1},{0,0,0,0,0,0,0,0,0,0,0,0,0,1,1,1},{0,0,0,0,0,0,0,0,0,0,0,0,0,0,1,1},{0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,1},{0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0} };
static const u32 INC32[8][8] __attribute__((aligned(32))) = {
{0,1,1,1,1,1,1,1},{0,0,1,1,1,1,1,1},{0,0,0,1,1,1,1,1},{0,0,0,0,1,1,1,1},
{0,0,0,0,0,1,1,1},{0,0,0,0,0,0,1,1},{0,0,0,0,0,0,0,1},{0,0,0,0,0,0,0,0}};
#define SFX16(p, r) do { __m256i _a = _mm256_load_si256((const __m256i *)(p)); \
__m256i _b = _mm256_load_si256((const __m256i *)(((const u16 *)INC16) + 16 * (size_t)(r))); \
_mm256_store_si256((__m256i *)(p), _mm256_add_epi16(_a, _b)); } while (0)
#define SFX32(p, r) do { __m256i _a = _mm256_load_si256((const __m256i *)(p)); \
__m256i _b = _mm256_load_si256((const __m256i *)(INC32 + (r))); \
_mm256_store_si256((__m256i *)(p), _mm256_add_epi32(_a, _b)); } while (0)
#define SFX32W(p, r) do { __m256i _a = _mm256_load_si256((const __m256i *)(p)); \
__m256i _c = _mm256_load_si256((const __m256i *)(((const u32 *)INC32W) + 16 * (size_t)(r))); \
_mm256_store_si256((__m256i *)(p), _mm256_add_epi32(_a, _c)); \
__m256i _d = _mm256_load_si256((const __m256i *)(p) + 1); \
__m256i _e = _mm256_load_si256((const __m256i *)(((const u32 *)INC32W) + 16 * (size_t)(r)) + 1); \
_mm256_store_si256((__m256i *)(p) + 1, _mm256_add_epi32(_d, _e)); } while (0)
static int FUSE_OK;
static inline u32 su8_32(const u8 *p, int cnt) {
__m256i v = _mm256_loadu_si256((const __m256i *)p);
__m256i m = _mm256_cmpgt_epi8(_mm256_set1_epi8((char)cnt),
_mm256_setr_epi8(0, 1, 2, 3, 4, 5, 6, 7, 8, 9, 10, 11, 12, 13, 14, 15,
16, 17, 18, 19, 20, 21, 22, 23, 24, 25, 26, 27, 28, 29, 30, 31));
__m256i x = _mm256_and_si256(v, m);
__m256i s = _mm256_sad_epu8(x, _mm256_setzero_si256());
__m128i t = _mm_add_epi64(_mm256_castsi256_si128(s), _mm256_extracti128_si256(s, 1));
return (u32)((unsigned long long)_mm_cvtsi128_si64(t) + (unsigned long long)_mm_extract_epi64(t, 1));
}
static inline u32 su16_8(const u16 *p, int cnt) {
__m128i v = _mm_loadu_si128((const __m128i *)p);
__m128i x = _mm_and_si128(v, _mm_cmpgt_epi16(_mm_set1_epi16((short)cnt),
_mm_setr_epi16(0, 1, 2, 3, 4, 5, 6, 7)));
__m128i s = _mm_madd_epi16(x, _mm_set1_epi16(1));
s = _mm_add_epi32(s, _mm_shuffle_epi32(s, _MM_SHUFFLE(1, 0, 3, 2)));
s = _mm_add_epi32(s, _mm_shuffle_epi32(s, _MM_SHUFFLE(2, 3, 0, 1)));
return (u32)_mm_cvtsi128_si32(s);
}
static inline u32 su32_8(const u32 *p, int cnt) {
__m256i v = _mm256_loadu_si256((const __m256i *)p);
__m256i x = _mm256_and_si256(v, _mm256_cmpgt_epi32(_mm256_set1_epi32(cnt),
_mm256_setr_epi32(0, 1, 2, 3, 4, 5, 6, 7)));
__m128i t = _mm_add_epi32(_mm256_castsi256_si128(x), _mm256_extracti128_si256(x, 1));
t = _mm_add_epi32(t, _mm_shuffle_epi32(t, _MM_SHUFFLE(1, 0, 3, 2)));
t = _mm_add_epi32(t, _mm_shuffle_epi32(t, _MM_SHUFFLE(2, 3, 0, 1)));
return (u32)_mm_cvtsi128_si32(t);
}
static const u8 INC8[32][32] __attribute__((aligned(32))) = {
{0,1,1,1,1,1,1,1,1,1,1,1,1,1,1,1,1,1,1,1,1,1,1,1,1,1,1,1,1,1,1,1},
{0,0,1,1,1,1,1,1,1,1,1,1,1,1,1,1,1,1,1,1,1,1,1,1,1,1,1,1,1,1,1,1},
{0,0,0,1,1,1,1,1,1,1,1,1,1,1,1,1,1,1,1,1,1,1,1,1,1,1,1,1,1,1,1,1},
{0,0,0,0,1,1,1,1,1,1,1,1,1,1,1,1,1,1,1,1,1,1,1,1,1,1,1,1,1,1,1,1},
{0,0,0,0,0,1,1,1,1,1,1,1,1,1,1,1,1,1,1,1,1,1,1,1,1,1,1,1,1,1,1,1},
{0,0,0,0,0,0,1,1,1,1,1,1,1,1,1,1,1,1,1,1,1,1,1,1,1,1,1,1,1,1,1,1},
{0,0,0,0,0,0,0,1,1,1,1,1,1,1,1,1,1,1,1,1,1,1,1,1,1,1,1,1,1,1,1,1},
{0,0,0,0,0,0,0,0,1,1,1,1,1,1,1,1,1,1,1,1,1,1,1,1,1,1,1,1,1,1,1,1},
{0,0,0,0,0,0,0,0,0,1,1,1,1,1,1,1,1,1,1,1,1,1,1,1,1,1,1,1,1,1,1,1},
{0,0,0,0,0,0,0,0,0,0,1,1,1,1,1,1,1,1,1,1,1,1,1,1,1,1,1,1,1,1,1,1},
{0,0,0,0,0,0,0,0,0,0,0,1,1,1,1,1,1,1,1,1,1,1,1,1,1,1,1,1,1,1,1,1},
{0,0,0,0,0,0,0,0,0,0,0,0,1,1,1,1,1,1,1,1,1,1,1,1,1,1,1,1,1,1,1,1},
{0,0,0,0,0,0,0,0,0,0,0,0,0,1,1,1,1,1,1,1,1,1,1,1,1,1,1,1,1,1,1,1},
{0,0,0,0,0,0,0,0,0,0,0,0,0,0,1,1,1,1,1,1,1,1,1,1,1,1,1,1,1,1,1,1},
{0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,1,1,1,1,1,1,1,1,1,1,1,1,1,1,1,1,1},
{0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,1,1,1,1,1,1,1,1,1,1,1,1,1,1,1,1},
{0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,1,1,1,1,1,1,1,1,1,1,1,1,1,1,1},
{0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,1,1,1,1,1,1,1,1,1,1,1,1,1,1},
{0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,1,1,1,1,1,1,1,1,1,1,1,1,1},
{0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,1,1,1,1,1,1,1,1,1,1,1,1},
{0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,1,1,1,1,1,1,1,1,1,1,1},
{0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,1,1,1,1,1,1,1,1,1,1},
{0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,1,1,1,1,1,1,1,1,1},
{0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,1,1,1,1,1,1,1,1},
{0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,1,1,1,1,1,1,1},
{0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,1,1,1,1,1,1},
{0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,1,1,1,1,1},
{0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,1,1,1,1},
{0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,1,1,1},
{0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,1,1},
{0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,1},
{0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0}
};
static inline void leaf_add(u32 v) {
__m256i *p=(__m256i *)(CT0+(v&~31u));
_mm256_store_si256(p,_mm256_add_epi8(_mm256_load_si256(p),_mm256_load_si256((const __m256i *)INC8[v&31])));
}
static inline void leaf_sub(u32 v) {
__m256i *p=(__m256i *)(CT0+(v&~31u));
_mm256_store_si256(p,_mm256_sub_epi8(_mm256_load_si256(p),_mm256_load_si256((const __m256i *)INC8[v&31])));
}
static inline void ctr_ins(u32 v) {
leaf_add(v);
SFX16(CT1 + ((v >> 9) << 4), (v >> 5) & 15);
SFX16(CT2 + ((v >> 13) << 4), (v >> 9) & 15);
SFX32W(CT3 + ((v >> 17) << 4), (v >> 13) & 15);
SFX32(CT4, (v >> 17) & 7);
}
static inline void ctr_clr(u32 v) {
/* [CLR0] the block must end at zero anyway (leaf_add is exactly undone), so a plain
32-byte zero store is enough -- no load, no subtract. Written as two 16-byte stores
so gcc 9.3 (the judge's compiler) emits exactly two stores instead of a split 256-bit
store (vmovdqu xmm + vinserti128). */
u8 *p = CT0 + (v & ~31u);
_mm_storeu_si128((__m128i *)p, _mm_setzero_si128());
_mm_storeu_si128((__m128i *)(p + 16), _mm_setzero_si128());
(void)0;
}
// [CTRR] live bound of the coordinate values seen by ctr_* (set in run_t)
static u32 CTN;
static inline void ctr_reset_hi(void) {
memset(CT1, 0, sizeof(u16) * (((size_t)CTN >> 5) + 32));
memset(CT2, 0, sizeof(u16) * (((size_t)CTN >> 9) + 32));
memset(CT3, 0, sizeof(u32) * (((size_t)CTN >> 13) + 32));
memset(CT4, 0, sizeof(CT4));
}
static inline u32 ctr_qry(u32 v) {
u32 s = CT4[v >> 17];
s += CT3[v >> 13];
s += CT2[v >> 9];
s += CT1[v >> 5];
s += CT0[v];
return s;
}
static int K;
#ifndef LFAUTO
#define LFAUTO 0
#endif
static int G_LIM;
static int LFSH;
static int LFD;
static const u32 *XD[KG_MAXK];
static u32 *OUT;
static u64 *AB, *ABtop;
// ===== v2 PACKED EVENT =====================================================
// v = [KA:20][TY:1][KB:20][ID:20] (bits 43..62, 42, 22..41, 0..19)
// KA = the compare key for THIS level (dim D+1), KB = the key for the NEXT
// level (dim D+2), carried so the cross needs no re-key gather.
// A plain u64 compare gives (KA, TY, KB, ID) -- KA primary, TY breaks equal
// KA with query(TY=0) before insert(TY=1) = the strict-inequality rule.
// Requires id < 2^20 and every coordinate < 2^20 (true: n = 1e6 < 2^20).
#define E_ID(v) ((u32)((v) & 0xFFFFFu))
#define E_KB(v) ((u32)(((v) >> 22) & 0xFFFFFu))
#define E_TY(v) ((u32)(((v) >> 42) & 1u))
#define E_KA(v) ((u32)((v) >> 43))
static inline u64 mk2(u32 ka, u32 ty, u32 kb, u32 id) {
return ((u64)ka << 43) | ((u64)ty << 42) | ((u64)kb << 22) | (u64)id;
}
static inline int eid(u64 v) { return (int)(v & 0x7FFFFFFFu); }
static inline int lfbucket(u64 v) { return (int)(((u32)(v >> 32)) >> LFSH); }
static inline int eq_(u64 v) { return (int)((v >> 31) & 1); }
static inline int keep_test(int takeA, int role) { return 1 - (takeA ^ role); }
static inline u64 mk(u32 key, int q, int id) {
return ((u64)key << 32) | ((u64)(1 - q) << 31) | (u64)id;
}
// resolve a fusion-level cross: events carry the last dim (K-1) coordinate in bits 32..63
// and are ordered by dim K-2; a left-insert before a right-query with a strictly smaller
// last-dim value is a hit.
static void ctr_solve(u64 *cross, int cx, u32 *dirt) {
int nd = 0;
for (int i = 0; i < cx; i++) {
#if MF_CTRPF > 0
if (i + MF_CTRPF < cx) {
u32 vf = E_KA(cross[i + MF_CTRPF]);
__builtin_prefetch(&CT0[vf], 1, 0); /* [PFCT] CT1/CT2 live in L1; only CT0 (1 MB) needs the lead */
}
#endif
u64 v = cross[i];
u32 val = E_KA(v); // v2: at the fuse level the counter key is dim K-1, and a
// cross event has just had KB promoted into the KA slot
if (E_TY(v)) { ctr_ins(val); dirt[nd++] = val; }
else OUT[E_ID(v)] += ctr_qry(val);
}
for (int i = 0; i < nd; i++) {
#if MF_CTRPF > 0
if (i + MF_CTRPF < nd) {
u32 vf = dirt[i + MF_CTRPF];
__builtin_prefetch(&CT0[vf], 1, 1);
}
#endif
ctr_clr(dirt[i]);
}
ctr_reset_hi();
}
// small base: compact insert list + SIMD comparison, then insertion sort by the key
static const u32 VMASK[8][8] = {
{0,0,0,0,0,0,0,0},
{~0u,0,0,0,0,0,0,0},
{~0u,~0u,0,0,0,0,0,0},
{~0u,~0u,~0u,0,0,0,0,0},
{~0u,~0u,~0u,~0u,0,0,0,0},
{~0u,~0u,~0u,~0u,~0u,0,0,0},
{~0u,~0u,~0u,~0u,~0u,~0u,0,0},
{~0u,~0u,~0u,~0u,~0u,~0u,~0u,0}};
// [LDU2] one 32-byte unaligned load; gcc 9.3 (judge) otherwise splits a 256-bit
// load whose alignment it cannot prove into vmovdqu xmm + vinserti128 (2 uops).
static inline __m256i ldu256(const void *p) {
__m256i v;
__asm__("vmovdqu %1, %0" : "=x"(v) : "m"(*(const __m256i *)p));
return v;
}
static inline u32 hsum256(__m256i v) {
__m128i t = _mm_add_epi32(_mm256_castsi256_si128(v), _mm256_extracti128_si256(v, 1));
t = _mm_add_epi32(t, _mm_shuffle_epi32(t, _MM_SHUFFLE(1, 0, 3, 2)));
t = _mm_add_epi32(t, _mm_shuffle_epi32(t, _MM_SHUFFLE(2, 3, 0, 1)));
return (u32)_mm_cvtsi128_si32(t);
}
template<int ND>
static void base_case_t(u64 *seq, int m, int d, u64 *out, int needOut) {
const int nd = ND;
__attribute__((aligned(32)))
u32 IK[KG_BASE + 8], IC0[KG_BASE + 8], IC1[KG_BASE + 8], IC2[KG_BASE + 8], IC3[KG_BASE + 8];
u32 QK[KG_BASE + 8], QC0[KG_BASE + 8], QC1[KG_BASE + 8], QC2[KG_BASE + 8], QC3[KG_BASE + 8];
u32 QID[KG_BASE + 8], QN[KG_BASE + 8];
const u32 *X1 = nd > 1 ? XD[d + 3] : 0, *X2 = nd > 2 ? XD[d + 4] : 0,
*X3 = nd > 3 ? XD[d + 5] : 0;
int ni = 0, nq = 0;
for (int t = 0; t < m; t++) {
#if MF_LEAFLP > 0
if (t + MF_LEAFLP < m) {
int idf = (int)E_ID(seq[t + MF_LEAFLP]);
if (nd > 1) __builtin_prefetch(&X1[idf], 0, 1);
if (nd > 2) __builtin_prefetch(&X2[idf], 0, 1);
if (nd > 3) __builtin_prefetch(&X3[idf], 0, 1);
}
#endif
u64 v = seq[t];
int id = (int)E_ID(v);
u32 k = E_KA(v);
u32 c0 = 0, c1 = 0, c2 = 0, c3 = 0;
if (nd > 0) c0 = E_KB(v); // dim d+2 is IN the event -- no gather
if (nd > 1) c1 = X1[id];
if (nd > 2) c2 = X2[id];
if (nd > 3) c3 = X3[id];
int isins = (int)E_TY(v);
IK[ni] = k; IC0[ni] = c0; IC1[ni] = c1; IC2[ni] = c2, IC3[ni] = c3; ni += isins;
QK[nq] = k; QC0[nq] = c0; QC1[nq] = c1; QC2[nq] = c2, QC3[nq] = c3; QID[nq] = (u32)id; QN[nq] = (u32)ni; nq += 1 - isins;
}
#if MF_LEAFG >= 2
/* c18/4-query group: one insert-block load feeds 4 chains */
#define MIN4(a0,a1,a2,a3) ((((a0) < (a1) ? (a0) : (a1))) < (((a2) < (a3) ? (a2) : (a3))) ? (((a0) < (a1) ? (a0) : (a1))) : (((a2) < (a3) ? (a2) : (a3))))
int q = 0;
for (; q + 4 <= nq; q += 4) {
const int n = MIN4((int)QN[q], (int)QN[q + 1], (int)QN[q + 2], (int)QN[q + 3]);
const int nA = (int)QN[q], nB = (int)QN[q + 1];
__m256i kvA = _mm256_set1_epi32((int)QK[q]), kvB = _mm256_set1_epi32((int)QK[q + 1]);
__m256i v0A = _mm256_set1_epi32((int)QC0[q]), v0B = _mm256_set1_epi32((int)QC0[q + 1]);
__m256i v1A = _mm256_set1_epi32((int)QC1[q]), v1B = _mm256_set1_epi32((int)QC1[q + 1]);
__m256i v2A = _mm256_set1_epi32((int)QC2[q]), v2B = _mm256_set1_epi32((int)QC2[q + 1]);
__m256i v3A = _mm256_set1_epi32((int)QC3[q]), v3B = _mm256_set1_epi32((int)QC3[q + 1]);
__m256i accA = _mm256_setzero_si256(), accB = _mm256_setzero_si256();
const int nC = (int)QN[q + 2], nD = (int)QN[q + 3];
__m256i kvC = _mm256_set1_epi32((int)QK[q + 2]), kvD = _mm256_set1_epi32((int)QK[q + 3]);
__m256i v0C = _mm256_set1_epi32((int)QC0[q + 2]), v0D = _mm256_set1_epi32((int)QC0[q + 3]);
__m256i v1C = _mm256_set1_epi32((int)QC1[q + 2]), v1D = _mm256_set1_epi32((int)QC1[q + 3]);
__m256i v2C = _mm256_set1_epi32((int)QC2[q + 2]), v2D = _mm256_set1_epi32((int)QC2[q + 3]);
__m256i v3C = _mm256_set1_epi32((int)QC3[q + 2]), v3D = _mm256_set1_epi32((int)QC3[q + 3]);
__m256i accC = _mm256_setzero_si256(), accD = _mm256_setzero_si256();
int j = 0;
if (nd >= 3) {
for (; j + 8 <= n; j += 8) {
__m256i IKv = _mm256_load_si256((const __m256i *)(IK + j));
__m256i C0v = _mm256_load_si256((const __m256i *)(IC0 + j));
__m256i C1v = _mm256_load_si256((const __m256i *)(IC1 + j));
__m256i C2v = _mm256_load_si256((const __m256i *)(IC2 + j));
__m256i mA = _mm256_and_si256(_mm256_cmpgt_epi32(kvA, IKv), _mm256_cmpgt_epi32(v0A, C0v));
__m256i mC = _mm256_and_si256(_mm256_cmpgt_epi32(kvC, IKv), _mm256_cmpgt_epi32(v0C, C0v));
__m256i mB = _mm256_and_si256(_mm256_cmpgt_epi32(kvB, IKv), _mm256_cmpgt_epi32(v0B, C0v));
__m256i mD = _mm256_and_si256(_mm256_cmpgt_epi32(kvD, IKv), _mm256_cmpgt_epi32(v0D, C0v));
mA = _mm256_and_si256(mA, _mm256_cmpgt_epi32(v1A, C1v));
mC = _mm256_and_si256(mC, _mm256_cmpgt_epi32(v1C, C1v));
mB = _mm256_and_si256(mB, _mm256_cmpgt_epi32(v1B, C1v));
mD = _mm256_and_si256(mD, _mm256_cmpgt_epi32(v1D, C1v));
mA = _mm256_and_si256(mA, _mm256_cmpgt_epi32(v2A, C2v));
mC = _mm256_and_si256(mC, _mm256_cmpgt_epi32(v2C, C2v));
mB = _mm256_and_si256(mB, _mm256_cmpgt_epi32(v2B, C2v));
mD = _mm256_and_si256(mD, _mm256_cmpgt_epi32(v2D, C2v));
if (nd > 3) { /* m161i key 3 */
__m256i C3v = _mm256_load_si256((const __m256i *)(IC3 + j));
mA = _mm256_and_si256(mA, _mm256_cmpgt_epi32(v3A, C3v));
mC = _mm256_and_si256(mC, _mm256_cmpgt_epi32(v3C, C3v));
mB = _mm256_and_si256(mB, _mm256_cmpgt_epi32(v3B, C3v));
mD = _mm256_and_si256(mD, _mm256_cmpgt_epi32(v3D, C3v));
}
accA = _mm256_sub_epi32(accA, mA);
accC = _mm256_sub_epi32(accC, mC);
accB = _mm256_sub_epi32(accB, mB);
accD = _mm256_sub_epi32(accD, mD);
}
} else {
for (; j + 8 <= n; j += 8) {
__m256i IKv = _mm256_load_si256((const __m256i *)(IK + j));
__m256i mA = _mm256_cmpgt_epi32(kvA, IKv);
__m256i mC = _mm256_cmpgt_epi32(kvC, IKv);
__m256i mB = _mm256_cmpgt_epi32(kvB, IKv);
__m256i mD = _mm256_cmpgt_epi32(kvD, IKv);
if (nd > 0) {
__m256i C0v = _mm256_load_si256((const __m256i *)(IC0 + j));
mA = _mm256_and_si256(mA, _mm256_cmpgt_epi32(v0A, C0v));
mC = _mm256_and_si256(mC, _mm256_cmpgt_epi32(v0C, C0v));
mB = _mm256_and_si256(mB, _mm256_cmpgt_epi32(v0B, C0v));
mD = _mm256_and_si256(mD, _mm256_cmpgt_epi32(v0D, C0v));
}
if (nd > 1) {
__m256i C1v = _mm256_load_si256((const __m256i *)(IC1 + j));
mA = _mm256_and_si256(mA, _mm256_cmpgt_epi32(v1A, C1v));
mC = _mm256_and_si256(mC, _mm256_cmpgt_epi32(v1C, C1v));
mB = _mm256_and_si256(mB, _mm256_cmpgt_epi32(v1B, C1v));
mD = _mm256_and_si256(mD, _mm256_cmpgt_epi32(v1D, C1v));
}
accA = _mm256_sub_epi32(accA, mA);
accC = _mm256_sub_epi32(accC, mC);
accB = _mm256_sub_epi32(accB, mB);
accD = _mm256_sub_epi32(accD, mD);
}
}
{
int jj = j;
__m256i acc = accA;
for (; jj + 8 <= nA; jj += 8) {
__m256i msk = _mm256_cmpgt_epi32(kvA, _mm256_load_si256((const __m256i *)(IK + jj)));
if (nd > 0) msk = _mm256_and_si256(msk, _mm256_cmpgt_epi32(v0A, _mm256_load_si256((const __m256i *)(IC0 + jj))));
if (nd > 1) msk = _mm256_and_si256(msk, _mm256_cmpgt_epi32(v1A, _mm256_load_si256((const __m256i *)(IC1 + jj))));
if (nd > 2) msk = _mm256_and_si256(msk, _mm256_cmpgt_epi32(v2A, _mm256_load_si256((const __m256i *)(IC2 + jj))));
if (nd > 3) msk = _mm256_and_si256(msk, _mm256_cmpgt_epi32(v3A, _mm256_load_si256((const __m256i *)(IC3 + jj))));
acc = _mm256_sub_epi32(acc, msk);
}
if (jj < nA) {
int n = nA;
__m256i msk = _mm256_cmpgt_epi32(kvA, _mm256_load_si256((const __m256i *)(IK + jj)));
if (nd > 0) msk = _mm256_and_si256(msk, _mm256_cmpgt_epi32(v0A, _mm256_load_si256((const __m256i *)(IC0 + jj))));
if (nd > 1) msk = _mm256_and_si256(msk, _mm256_cmpgt_epi32(v1A, _mm256_load_si256((const __m256i *)(IC1 + jj))));
if (nd > 2) msk = _mm256_and_si256(msk, _mm256_cmpgt_epi32(v2A, _mm256_load_si256((const __m256i *)(IC2 + jj))));
if (nd > 3) msk = _mm256_and_si256(msk, _mm256_cmpgt_epi32(v3A, _mm256_load_si256((const __m256i *)(IC3 + jj))));
msk = _mm256_and_si256(msk, _mm256_loadu_si256((const __m256i *)VMASK[n - jj]));
acc = _mm256_sub_epi32(acc, msk);
}
OUT[QID[q]] += hsum256(acc);
}
{
int jj = j;
__m256i acc = accB;
for (; jj + 8 <= nB; jj += 8) {
__m256i msk = _mm256_cmpgt_epi32(kvB, _mm256_load_si256((const __m256i *)(IK + jj)));
if (nd > 0) msk = _mm256_and_si256(msk, _mm256_cmpgt_epi32(v0B, _mm256_load_si256((const __m256i *)(IC0 + jj))));
if (nd > 1) msk = _mm256_and_si256(msk, _mm256_cmpgt_epi32(v1B, _mm256_load_si256((const __m256i *)(IC1 + jj))));
if (nd > 2) msk = _mm256_and_si256(msk, _mm256_cmpgt_epi32(v2B, _mm256_load_si256((const __m256i *)(IC2 + jj))));
if (nd > 3) msk = _mm256_and_si256(msk, _mm256_cmpgt_epi32(v3B, _mm256_load_si256((const __m256i *)(IC3 + jj))));
acc = _mm256_sub_epi32(acc, msk);
}
if (jj < nB) {
int n = nB;
__m256i msk = _mm256_cmpgt_epi32(kvB, _mm256_load_si256((const __m256i *)(IK + jj)));
if (nd > 0) msk = _mm256_and_si256(msk, _mm256_cmpgt_epi32(v0B, _mm256_load_si256((const __m256i *)(IC0 + jj))));
if (nd > 1) msk = _mm256_and_si256(msk, _mm256_cmpgt_epi32(v1B, _mm256_load_si256((const __m256i *)(IC1 + jj))));
if (nd > 2) msk = _mm256_and_si256(msk, _mm256_cmpgt_epi32(v2B, _mm256_load_si256((const __m256i *)(IC2 + jj))));
if (nd > 3) msk = _mm256_and_si256(msk, _mm256_cmpgt_epi32(v3B, _mm256_load_si256((const __m256i *)(IC3 + jj))));
msk = _mm256_and_si256(msk, _mm256_loadu_si256((const __m256i *)VMASK[n - jj]));
acc = _mm256_sub_epi32(acc, msk);
}
OUT[QID[q + 1]] += hsum256(acc);
}
{
int jj = j;
__m256i acc = accC;
for (; jj + 8 <= nC; jj += 8) {
__m256i msk = _mm256_cmpgt_epi32(kvC, _mm256_load_si256((const __m256i *)(IK + jj)));
if (nd > 0) msk = _mm256_and_si256(msk, _mm256_cmpgt_epi32(v0C, _mm256_load_si256((const __m256i *)(IC0 + jj))));
if (nd > 1) msk = _mm256_and_si256(msk, _mm256_cmpgt_epi32(v1C, _mm256_load_si256((const __m256i *)(IC1 + jj))));
if (nd > 2) msk = _mm256_and_si256(msk, _mm256_cmpgt_epi32(v2C, _mm256_load_si256((const __m256i *)(IC2 + jj))));
if (nd > 3) msk = _mm256_and_si256(msk, _mm256_cmpgt_epi32(v3C, _mm256_load_si256((const __m256i *)(IC3 + jj))));
acc = _mm256_sub_epi32(acc, msk);
}
if (jj < nC) {
int n = nC;
__m256i msk = _mm256_cmpgt_epi32(kvC, _mm256_load_si256((const __m256i *)(IK + jj)));
if (nd > 0) msk = _mm256_and_si256(msk, _mm256_cmpgt_epi32(v0C, _mm256_load_si256((const __m256i *)(IC0 + jj))));
if (nd > 1) msk = _mm256_and_si256(msk, _mm256_cmpgt_epi32(v1C, _mm256_load_si256((const __m256i *)(IC1 + jj))));
if (nd > 2) msk = _mm256_and_si256(msk, _mm256_cmpgt_epi32(v2C, _mm256_load_si256((const __m256i *)(IC2 + jj))));
if (nd > 3) msk = _mm256_and_si256(msk, _mm256_cmpgt_epi32(v3C, _mm256_load_si256((const __m256i *)(IC3 + jj))));
msk = _mm256_and_si256(msk, _mm256_loadu_si256((const __m256i *)VMASK[n - jj]));
acc = _mm256_sub_epi32(acc, msk);
}
OUT[QID[q + 2]] += hsum256(acc);
}
{
int jj = j;
__m256i acc = accD;
for (; jj + 8 <= nD; jj += 8) {
__m256i msk = _mm256_cmpgt_epi32(kvD, _mm256_load_si256((const __m256i *)(IK + jj)));
if (nd > 0) msk = _mm256_and_si256(msk, _mm256_cmpgt_epi32(v0D, _mm256_load_si256((const __m256i *)(IC0 + jj))));
if (nd > 1) msk = _mm256_and_si256(msk, _mm256_cmpgt_epi32(v1D, _mm256_load_si256((const __m256i *)(IC1 + jj))));
if (nd > 2) msk = _mm256_and_si256(msk, _mm256_cmpgt_epi32(v2D, _mm256_load_si256((const __m256i *)(IC2 + jj))));
if (nd > 3) msk = _mm256_and_si256(msk, _mm256_cmpgt_epi32(v3D, _mm256_load_si256((const __m256i *)(IC3 + jj))));
acc = _mm256_sub_epi32(acc, msk);
}
if (jj < nD) {
int n = nD;
__m256i msk = _mm256_cmpgt_epi32(kvD, _mm256_load_si256((const __m256i *)(IK + jj)));
if (nd > 0) msk = _mm256_and_si256(msk, _mm256_cmpgt_epi32(v0D, _mm256_load_si256((const __m256i *)(IC0 + jj))));
if (nd > 1) msk = _mm256_and_si256(msk, _mm256_cmpgt_epi32(v1D, _mm256_load_si256((const __m256i *)(IC1 + jj))));
if (nd > 2) msk = _mm256_and_si256(msk, _mm256_cmpgt_epi32(v2D, _mm256_load_si256((const __m256i *)(IC2 + jj))));
if (nd > 3) msk = _mm256_and_si256(msk, _mm256_cmpgt_epi32(v3D, _mm256_load_si256((const __m256i *)(IC3 + jj))));
msk = _mm256_and_si256(msk, _mm256_loadu_si256((const __m256i *)VMASK[n - jj]));
acc = _mm256_sub_epi32(acc, msk);
}
OUT[QID[q + 3]] += hsum256(acc);
}
}
for (; q + 2 <= nq; q += 2) {
const int n = ((int)QN[q] < (int)QN[q + 1]) ? (int)QN[q] : (int)QN[q + 1];
const int nA = (int)QN[q], nB = (int)QN[q + 1];
__m256i kvA = _mm256_set1_epi32((int)QK[q]), kvB = _mm256_set1_epi32((int)QK[q + 1]);
__m256i v0A = _mm256_set1_epi32((int)QC0[q]), v0B = _mm256_set1_epi32((int)QC0[q + 1]);
__m256i v1A = _mm256_set1_epi32((int)QC1[q]), v1B = _mm256_set1_epi32((int)QC1[q + 1]);
__m256i v2A = _mm256_set1_epi32((int)QC2[q]), v2B = _mm256_set1_epi32((int)QC2[q + 1]);
__m256i v3A = _mm256_set1_epi32((int)QC3[q]), v3B = _mm256_set1_epi32((int)QC3[q + 1]);
__m256i accA = _mm256_setzero_si256(), accB = _mm256_setzero_si256();
int j = 0;
if (nd >= 3) {
for (; j + 8 <= n; j += 8) {
__m256i IKv = _mm256_load_si256((const __m256i *)(IK + j));
__m256i C0v = _mm256_load_si256((const __m256i *)(IC0 + j));
__m256i C1v = _mm256_load_si256((const __m256i *)(IC1 + j));
__m256i C2v = _mm256_load_si256((const __m256i *)(IC2 + j));
__m256i mA = _mm256_and_si256(_mm256_cmpgt_epi32(kvA, IKv), _mm256_cmpgt_epi32(v0A, C0v));
__m256i mB = _mm256_and_si256(_mm256_cmpgt_epi32(kvB, IKv), _mm256_cmpgt_epi32(v0B, C0v));
mA = _mm256_and_si256(mA, _mm256_cmpgt_epi32(v1A, C1v));
mB = _mm256_and_si256(mB, _mm256_cmpgt_epi32(v1B, C1v));
mA = _mm256_and_si256(mA, _mm256_cmpgt_epi32(v2A, C2v));
mB = _mm256_and_si256(mB, _mm256_cmpgt_epi32(v2B, C2v));
if (nd > 3) { /* m161i key 3 */
__m256i C3v = _mm256_load_si256((const __m256i *)(IC3 + j));
mA = _mm256_and_si256(mA, _mm256_cmpgt_epi32(v3A, C3v));
mB = _mm256_and_si256(mB, _mm256_cmpgt_epi32(v3B, C3v));
}
accA = _mm256_sub_epi32(accA, mA);
accB = _mm256_sub_epi32(accB, mB);
}
} else {
for (; j + 8 <= n; j += 8) {
__m256i IKv = _mm256_load_si256((const __m256i *)(IK + j));
__m256i mA = _mm256_cmpgt_epi32(kvA, IKv);
__m256i mB = _mm256_cmpgt_epi32(kvB, IKv);
if (nd > 0) {
__m256i C0v = _mm256_load_si256((const __m256i *)(IC0 + j));
mA = _mm256_and_si256(mA, _mm256_cmpgt_epi32(v0A, C0v));
mB = _mm256_and_si256(mB, _mm256_cmpgt_epi32(v0B, C0v));
}
if (nd > 1) {
__m256i C1v = _mm256_load_si256((const __m256i *)(IC1 + j));
mA = _mm256_and_si256(mA, _mm256_cmpgt_epi32(v1A, C1v));
mB = _mm256_and_si256(mB, _mm256_cmpgt_epi32(v1B, C1v));
}
accA = _mm256_sub_epi32(accA, mA);
accB = _mm256_sub_epi32(accB, mB);
}
}
// per-query tails: each member finishes [j, QN[.]) on its own
{
int jj = j;
__m256i acc = accA;
for (; jj + 8 <= nA; jj += 8) {
__m256i msk = _mm256_cmpgt_epi32(kvA, _mm256_load_si256((const __m256i *)(IK + jj)));
if (nd > 0) msk = _mm256_and_si256(msk, _mm256_cmpgt_epi32(v0A, _mm256_load_si256((const __m256i *)(IC0 + jj))));
if (nd > 1) msk = _mm256_and_si256(msk, _mm256_cmpgt_epi32(v1A, _mm256_load_si256((const __m256i *)(IC1 + jj))));
if (nd > 2) msk = _mm256_and_si256(msk, _mm256_cmpgt_epi32(v2A, _mm256_load_si256((const __m256i *)(IC2 + jj))));
if (nd > 3) msk = _mm256_and_si256(msk, _mm256_cmpgt_epi32(v3A, _mm256_load_si256((const __m256i *)(IC3 + jj))));
acc = _mm256_sub_epi32(acc, msk);
}
if (jj < nA) {
int n = nA;
__m256i msk = _mm256_cmpgt_epi32(kvA, _mm256_load_si256((const __m256i *)(IK + jj)));
if (nd > 0) msk = _mm256_and_si256(msk, _mm256_cmpgt_epi32(v0A, _mm256_load_si256((const __m256i *)(IC0 + jj))));
if (nd > 1) msk = _mm256_and_si256(msk, _mm256_cmpgt_epi32(v1A, _mm256_load_si256((const __m256i *)(IC1 + jj))));
if (nd > 2) msk = _mm256_and_si256(msk, _mm256_cmpgt_epi32(v2A, _mm256_load_si256((const __m256i *)(IC2 + jj))));
if (nd > 3) msk = _mm256_and_si256(msk, _mm256_cmpgt_epi32(v3A, _mm256_load_si256((const __m256i *)(IC3 + jj))));
msk = _mm256_and_si256(msk, _mm256_loadu_si256((const __m256i *)VMASK[n - jj]));
acc = _mm256_sub_epi32(acc, msk);
}
OUT[QID[q]] += hsum256(acc);
}
{
int jj = j;
__m256i acc = accB;
for (; jj + 8 <= nB; jj += 8) {
__m256i msk = _mm256_cmpgt_epi32(kvB, _mm256_load_si256((const __m256i *)(IK + jj)));
if (nd > 0) msk = _mm256_and_si256(msk, _mm256_cmpgt_epi32(v0B, _mm256_load_si256((const __m256i *)(IC0 + jj))));
if (nd > 1) msk = _mm256_and_si256(msk, _mm256_cmpgt_epi32(v1B, _mm256_load_si256((const __m256i *)(IC1 + jj))));
if (nd > 2) msk = _mm256_and_si256(msk, _mm256_cmpgt_epi32(v2B, _mm256_load_si256((const __m256i *)(IC2 + jj))));
if (nd > 3) msk = _mm256_and_si256(msk, _mm256_cmpgt_epi32(v3B, _mm256_load_si256((const __m256i *)(IC3 + jj))));
acc = _mm256_sub_epi32(acc, msk);
}
if (jj < nB) {
int n = nB;
__m256i msk = _mm256_cmpgt_epi32(kvB, _mm256_load_si256((const __m256i *)(IK + jj)));
if (nd > 0) msk = _mm256_and_si256(msk, _mm256_cmpgt_epi32(v0B, _mm256_load_si256((const __m256i *)(IC0 + jj))));
if (nd > 1) msk = _mm256_and_si256(msk, _mm256_cmpgt_epi32(v1B, _mm256_load_si256((const __m256i *)(IC1 + jj))));
if (nd > 2) msk = _mm256_and_si256(msk, _mm256_cmpgt_epi32(v2B, _mm256_load_si256((const __m256i *)(IC2 + jj))));
if (nd > 3) msk = _mm256_and_si256(msk, _mm256_cmpgt_epi32(v3B, _mm256_load_si256((const __m256i *)(IC3 + jj))));
msk = _mm256_and_si256(msk, _mm256_loadu_si256((const __m256i *)VMASK[n - jj]));
acc = _mm256_sub_epi32(acc, msk);
}
OUT[QID[q + 1]] += hsum256(acc);
}
}
#endif
for (; q < nq; q++) {
const int n = (int)QN[q];
__m256i kv = _mm256_set1_epi32((int)QK[q]);
__m256i v0 = _mm256_set1_epi32((int)QC0[q]);
__m256i v1 = _mm256_set1_epi32((int)QC1[q]);
__m256i v2 = _mm256_set1_epi32((int)QC2[q]);
__m256i v3 = _mm256_set1_epi32((int)QC3[q]);
__m256i acc = _mm256_setzero_si256();
int j = 0;
switch (nd) {
case 0:
for (; j + 8 <= n; j += 8)
acc = _mm256_sub_epi32(acc, _mm256_cmpgt_epi32(kv, _mm256_load_si256((const __m256i *)(IK + j))));
break;
case 1:
for (; j + 8 <= n; j += 8) {
__m256i msk = _mm256_cmpgt_epi32(kv, _mm256_load_si256((const __m256i *)(IK + j)));
msk = _mm256_and_si256(msk, _mm256_cmpgt_epi32(v0, _mm256_load_si256((const __m256i *)(IC0 + j))));
acc = _mm256_sub_epi32(acc, msk);
}
break;
case 2:
for (; j + 8 <= n; j += 8) {
__m256i msk = _mm256_cmpgt_epi32(kv, _mm256_load_si256((const __m256i *)(IK + j)));
msk = _mm256_and_si256(msk, _mm256_cmpgt_epi32(v0, _mm256_load_si256((const __m256i *)(IC0 + j))));
msk = _mm256_and_si256(msk, _mm256_cmpgt_epi32(v1, _mm256_load_si256((const __m256i *)(IC1 + j))));
acc = _mm256_sub_epi32(acc, msk);
}
break;
default:
for (; j + 8 <= n; j += 8) {
__m256i msk = _mm256_cmpgt_epi32(kv, _mm256_load_si256((const __m256i *)(IK + j)));
msk = _mm256_and_si256(msk, _mm256_cmpgt_epi32(v0, _mm256_load_si256((const __m256i *)(IC0 + j))));
msk = _mm256_and_si256(msk, _mm256_cmpgt_epi32(v1, _mm256_load_si256((const __m256i *)(IC1 + j))));
msk = _mm256_and_si256(msk, _mm256_cmpgt_epi32(v2, _mm256_load_si256((const __m256i *)(IC2 + j))));
if (nd > 3) msk = _mm256_and_si256(msk, _mm256_cmpgt_epi32(v3, _mm256_load_si256((const __m256i *)(IC3 + j))));
acc = _mm256_sub_epi32(acc, msk);
}
break;
}
if (j < n) { // partial block, masked by a table lookup
__m256i msk = _mm256_cmpgt_epi32(kv, _mm256_load_si256((const __m256i *)(IK + j)));
if (nd > 0) msk = _mm256_and_si256(msk, _mm256_cmpgt_epi32(v0, _mm256_load_si256((const __m256i *)(IC0 + j))));
if (nd > 1) msk = _mm256_and_si256(msk, _mm256_cmpgt_epi32(v1, _mm256_load_si256((const __m256i *)(IC1 + j))));
if (nd > 2) msk = _mm256_and_si256(msk, _mm256_cmpgt_epi32(v2, _mm256_load_si256((const __m256i *)(IC2 + j))));
if (nd > 3) msk = _mm256_and_si256(msk, _mm256_cmpgt_epi32(v3, _mm256_load_si256((const __m256i *)(IC3 + j))));
msk = _mm256_and_si256(msk, _mm256_loadu_si256((const __m256i *)VMASK[n - j]));
acc = _mm256_sub_epi32(acc, msk);
}
u32 c = hsum256(acc);
OUT[QID[q]] += c;
}
if (!needOut) return;
#if LEAFNB
{
const int NB = (LFAUTO && m <= LEAFNB) ? (m < 8 ? 8 : m) : LEAFNB;
const int LB = (G_LIM < NB) ? G_LIM : NB; // used-bucket bound
/* [SORTB] one extra leading slot so that cnt[-1] == 0 is a legal read: the
bucket-b discovery pass below can then vectorise over every bucket. */
u32 cntb[LEAFNB + 3];
u32 *cnt = cntb + 1;
cnt[-1] = 0;
for (int b = 0; b < LB; b++) cnt[b] = 0;
for (int t = 0; t < m; t++) cnt[(int)(E_KA(seq[t]) >> (LFSH + LFD)) + 1]++;
for (int b = 0; b + 1 < LB; b++) cnt[b + 1] += cnt[b];
for (int t = 0; t < m; t++) { u64 v = seq[t]; out[cnt[(int)(E_KA(v) >> (LFSH + LFD))]++] = v; }
/* [SORTB] buckets holding 0 or 1 events need no work, and with LEAFNB >= m the
great majority are of that kind; a 256-bit subtract/compare/movemask finds the
multi-event buckets eight at a time, so the scan loop takes no data-dependent
branch at all -- only the 1/8-odd that really hold two or more are revisited. */
{
const __m256i one1 = _mm256_set1_epi32(1);
int b = 0;
/* [SORTB3] thirty-two buckets per iteration: four vector pairs build one
32-bit mask, so the loop takes one empty-mask exit per 32 slots. */
for (; b + 32 <= LB - 1; b += 32) {
unsigned mask = 0;
for (int h = 0; h < 4; h++) {
__m256i zz = _mm256_sub_epi32(ldu256((const __m256i *)(cnt + b + 8 * h)),
ldu256((const __m256i *)(cnt + b + 8 * h - 1)));
mask |= ((unsigned)_mm256_movemask_ps(_mm256_castsi256_ps(_mm256_cmpgt_epi32(zz, one1)))) << (8 * h);
}
while (mask) {
int k = __builtin_ctz(mask);
mask &= mask - 1;
int bb = b + k;
int lo = (int)cnt[bb - 1];
int hi = (int)cnt[bb];
if (hi - lo == 2) {
u64 a = out[lo], c = out[lo + 1];
int lt = (a < c);
out[lo] = lt ? a : c;
out[lo + 1] = lt ? c : a;
} else if (hi - lo <= 4) {
/* [CE34] fixed compare-exchange networks for three and four events */
u64 a = out[lo], b = out[lo + 1], c = out[lo + 2];
if (hi - lo == 3) {
u64 x = a, y = b; int lt = x < y; a = lt ? x : y; b = lt ? y : x;
x = b; y = c; lt = x < y; b = lt ? x : y; c = lt ? y : x;
x = a; y = b; lt = x < y; a = lt ? x : y; b = lt ? y : x;
out[lo] = a; out[lo + 1] = b; out[lo + 2] = c;
} else {
u64 d = out[lo + 3], x, y; int lt;
x = a; y = b; lt = x < y; a = lt ? x : y; b = lt ? y : x;
x = c; y = d; lt = x < y; c = lt ? x : y; d = lt ? y : x;
x = a; y = c; lt = x < y; a = lt ? x : y; c = lt ? y : x;
x = b; y = d; lt = x < y; b = lt ? x : y; d = lt ? y : x;
x = b; y = c; lt = x < y; b = lt ? x : y; c = lt ? y : x;
out[lo] = a; out[lo + 1] = b; out[lo + 2] = c; out[lo + 3] = d;
}
} else
for (int t = lo + 1; t < hi; t++) {
u64 v = out[t];
int q = t - 1;
while (q >= lo && v < out[q]) { out[q + 1] = out[q]; q--; }
out[q + 1] = v;
}
}
}
/* [SORTB2] sixteen buckets per iteration: one mask (and one empty-mask exit)
per sixteen slots instead of per eight. */
for (; b + 16 <= LB - 1; b += 16) {
__m256i z0 = _mm256_sub_epi32(ldu256((const __m256i *)(cnt + b)),
ldu256((const __m256i *)(cnt + b - 1)));
__m256i z1 = _mm256_sub_epi32(ldu256((const __m256i *)(cnt + b + 8)),
ldu256((const __m256i *)(cnt + b + 7)));
unsigned mask = (unsigned)_mm256_movemask_ps(_mm256_castsi256_ps(_mm256_cmpgt_epi32(z0, one1)))
| ((unsigned)_mm256_movemask_ps(_mm256_castsi256_ps(_mm256_cmpgt_epi32(z1, one1))) << 8);
while (mask) {
int k = __builtin_ctz(mask);
mask &= mask - 1;
int bb = b + k;
int lo = (int)cnt[bb - 1];
int hi = (int)cnt[bb];
if (hi - lo == 2) {
u64 a = out[lo], c = out[lo + 1];
int lt = (a < c);
out[lo] = lt ? a : c;
out[lo + 1] = lt ? c : a;
} else
for (int t = lo + 1; t < hi; t++) {
u64 v = out[t];
int q = t - 1;
while (q >= lo && v < out[q]) { out[q + 1] = out[q]; q--; }
out[q + 1] = v;
}
}
}
for (; b + 8 <= LB - 1; b += 8) {
__m256i cur = ldu256((const __m256i *)(cnt + b));
__m256i pre = ldu256((const __m256i *)(cnt + b - 1));
__m256i sz = _mm256_sub_epi32(cur, pre);
unsigned mask = (unsigned)_mm256_movemask_ps(_mm256_castsi256_ps(_mm256_cmpgt_epi32(sz, one1)));
while (mask) {
int k = __builtin_ctz(mask);
mask &= mask - 1;
int bb = b + k;
int lo = (int)cnt[bb - 1];
int hi = (int)cnt[bb];
if (hi - lo == 2) {
/* [SORT2] two straight comparisons, no data-dependent branch */
u64 a = out[lo], c = out[lo + 1];
int lt = (a < c);
out[lo] = lt ? a : c;
out[lo + 1] = lt ? c : a;
} else
for (int t = lo + 1; t < hi; t++) {
u64 v = out[t];
int q = t - 1;
while (q >= lo && v < out[q]) { out[q + 1] = out[q]; q--; }
out[q + 1] = v;
}
}
}
for (; b + 1 < LB; b++) {
int lo = (int)cnt[b - 1];
int hi = (int)cnt[b];
if (hi - lo > 1) {
if (hi - lo == 2) {
u64 a = out[lo], c = out[lo + 1];
int lt = (a < c);
out[lo] = lt ? a : c;
out[lo + 1] = lt ? c : a;
} else
for (int t = lo + 1; t < hi; t++) {
u64 v = out[t];
int q = t - 1;
while (q >= lo && v < out[q]) { out[q + 1] = out[q]; q--; }
out[q + 1] = v;
}
}
}
}
}
#else
for (int t = 1; t < m; t++) {
u64 v = seq[t];
int q = t - 1;
while (q >= 0 && v < seq[q]) { seq[q + 1] = seq[q]; q--; }
seq[q + 1] = v;
}
memcpy(out, seq, sizeof(u64) * (size_t)m);
#endif
}
static void base_case(u64 *seq, int m, int d, u64 *out, int needOut) {
switch (K - (d + 2)) {
case 0: base_case_t<0>(seq, m, d, out, needOut); break;
case 1: base_case_t<1>(seq, m, d, out, needOut); break;
case 2: base_case_t<2>(seq, m, d, out, needOut); break;
case 3: base_case_t<3>(seq, m, d, out, needOut); break;
case 4: base_case_t<4>(seq, m, d, out, needOut); break;
default: base_case_t<4>(seq, m, d, out, needOut); break;
}
}
// ---- ND specialisation dispatcher (the fleet's `DN` recipe) ----
static void base_case(u64 *seq, int m, int d, u64 *out, int needOut);
// ================= depth-templated solve: D and KT are compile-time =================
template <bool B> struct TB {};
template <int V> struct NCL { static const int v = (V < 0) ? 0 : ((V > 4) ? 4 : V); };
template <int KT, int D, bool NEED> static inline void solve_t(u64 *seq, int m, u64 *out);
// ---- terminal level (D == KT-1): every insert-before-query pair qualifies ----
template <int KT, int D, bool NEED>
static void solve_b(u64 *seq, int m, u64 *out, TB<true>) {
int cnt = 0;
for (int t = 0; t < m; t++) {
u64 v = seq[t];
if (E_TY(v)) cnt++; // insert
else OUT[E_ID(v)] += (u32)cnt; // query
}
}
// ---- general level (D < KT-1) ----
template <int KT, int D, bool NEED>
static void solve_b(u64 *seq, int m, u64 *out, TB<false>) {
if (m <= KG_BASE) {
base_case_t<NCL< KT - (D + 2) >::v>(seq, m, D, out, NEED);
return;
}
int h = m >> 1, r = m - h;
u64 *save = ABtop;
u64 *retA = ABtop; ABtop += h+1;
u64 *retB = ABtop; ABtop += r;
solve_t<KT, D, true>(seq, h, retA);
solve_t<KT, D, true>(seq + h, r, retB);
u64 *cross = ABtop; ABtop += m + 1;
// v2: PROMOTE moves the carried KB up into the KA slot (a new level's key);
// GATHER fetches the next KB (dim D+3) -- the ONLY remaining per-level gather,
// and it exists only at D = 0 for KT = 4.
const int PROM = (D + 2 <= KT - 1);
const int GAT = (D + 3 <= KT - 1);
const u32 *XN = XD[(D + 3 <= KG_MAXK - 1) ? (D + 3) : 0];
retA[h]=~0ULL;
int ia = 0, ib = 0, cx = 0, t = 0;
while (ib < r) {
#if KG_PF
if (ia + KG_PF < h) __builtin_prefetch(&XN[E_ID(retA[ia + KG_PF])], 0, 1); /* c18: DENSE, distance 64 */
if (ib + KG_PF < r) __builtin_prefetch(&XN[E_ID(retB[ib + KG_PF])], 0, 1);
#endif
u64 ea = retA[ia], eb = retB[ib];
int takeA = (ea < eb);
ia += takeA; ib += 1 - takeA;
u64 v = takeA ? ea : eb;
if (NEED) out[t] = v;
t++;
int role = (int)E_TY(v);
#if KG_COND
if (keep_test(takeA, role)) { cross[cx] = v; cx++; }
#else
cross[cx] = v;
cx += keep_test(takeA, role);
#endif
}
while (ia < h) {
u64 v = retA[ia++];
if (NEED) out[t] = v;
t++;
if (E_TY(v)) cross[cx++] = v;
}
if (GAT) { // ----- C1: cross re-key, second pass -----
int j=0;
const __m256i idmask=_mm256_set1_epi64x(0xfffffULL);
const __m256i kmask=_mm256_set1_epi64x(0xfffffULL<<43);
const __m256i keepmask=_mm256_set1_epi64x((1ULL<<42)|0xfffffULL);
for(;j+4<=cx;j+=4){
if(j+67<cx)for(int z=0;z<4;z++)__builtin_prefetch(&XN[E_ID(cross[j+64+z])],0,1);
__m256i v=_mm256_loadu_si256((const __m256i*)(cross+j));
__m256i ix=_mm256_and_si256(v,idmask);
__m128i c=_mm256_i64gather_epi32((const int*)XN,ix,4);
__m256i kb=_mm256_slli_epi64(_mm256_cvtepu32_epi64(c),22);
__m256i ka=_mm256_and_si256(_mm256_slli_epi64(v,21),kmask);
__m256i res=_mm256_or_si256(_mm256_or_si256(ka,kb),_mm256_and_si256(v,keepmask));
_mm256_storeu_si256((__m256i*)(cross+j),res);
}
for (; j < cx; j++) {
if (j + 64 < cx) __builtin_prefetch(&XN[E_ID(cross[j + 64])], 0, 1);
u64 v = cross[j];
cross[j] = ((v << 21) & (0xFFFFFull << 43)) | (v & (1ull << 42)) | (v & 0xFFFFFull) |
((u64)XN[E_ID(v)] << 22);
}
} else if (PROM) {
for (int j = 0; j < cx; j++) {
u64 v = cross[j];
cross[j] = ((v << 21) & (0xFFFFFull << 43)) | (v & (1ull << 42)) | (v & 0xFFFFFull);
}
}
if (FUSE_OK && D == KT - 3) {
u32 *dirt = (u32 *)ABtop; ABtop += (cx + 3) / 2 + 1;
ctr_solve(cross, cx, dirt);
} else {
u64 *cres = ABtop; ABtop += cx;
solve_t<KT, D + 1, false>(cross, cx, cres);
}
ABtop = save;
}
template <int KT, int D, bool NEED>
static inline void solve_t(u64 *seq, int m, u64 *out) {
if (m <= 1) return;
solve_b<KT, D, NEED>(seq, m, out, TB<(D >= KT - 1)>());
}
// ====================================================================================
static u32 *g_ord, *g_x;
static void sort0(int n) {
const u32 *X0 = XD[0];
u32 *cnt = g_x;
memset(cnt, 0, sizeof(u32) * (size_t)(n + 1));
for (int i = 0; i < n; i++) cnt[X0[i] + 1]++;
for (int i = 0; i < n; i++) cnt[i + 1] += cnt[i];
for (int i = 0; i < n; i++) g_ord[cnt[X0[i]]++] = (u32)i;
}
template <int KT>
static void run_t(int n, const u32 **x, u32 *out) {
LFSH = 0; LFD = 0; CTN = (u32)n;
for (int d = 0; d < KT; d++) XD[d] = x[d];
OUT = out;
ABtop = AB;
if (n <= 1) return;
FUSE_OK = (KT >= 4 && n > 2 && (u32)n < CT_MAXV);
if(FUSE_OK){
const u32*XL=XD[KT-1];size_t count=((size_t)n+31)/32;
memset(CT1,0,(count+8)*sizeof(u16));
for(int i=0;i<n;i++)if(++CT1[XL[i]>>5]>255){FUSE_OK=0;break;}
memset(CT1,0,(count+8)*sizeof(u16));
}
{
int bit = 0;
while ((n >> bit) >= LEAFNB && bit < 30) bit++;
LFSH = bit; LFD = 0;
}
G_LIM = (int)(((u32)(n - 1)) >> (LFSH + LFD)) + 2;
sort0(n);
const u32 *X0 = XD[0], *X1 = XD[1];
u64 *ev = AB; // top-level event array: 2n entries
int m = 0, g = 0;
while (g < n) {
int g1 = g;
u32 v = X0[g_ord[g]];
while (g1 < n && X0[g_ord[g1]] == v) g1++;
const u32 *X2 = XD[2];
for (int t = g; t < g1; t++) { int q = (int)g_ord[t]; ev[m++] = mk2(X1[q], 0, X2[q], (u32)q); }
for (int t = g; t < g1; t++) { int q = (int)g_ord[t]; ev[m++] = mk2(X1[q], 1, X2[q], (u32)q); }
g = g1;
}
ABtop = AB + m;
u64 *res = ABtop; ABtop += m;
solve_t<KT, 0, false>(ev, m, res);
}
} // namespace kg
static u64 g_arena[KG_ARENA];
static u32 g_ord_st[1048576], g_x_st[1048576];
static int g_ready = 0;
static void setup() {
if (g_ready) return;
g_ready = 1;
kg::AB = g_arena;
kg::g_ord = g_ord_st;
kg::g_x = g_x_st;
}
void count_4d(int n, const unsigned *x[4], unsigned *out) {
setup();
memset(out, 0, sizeof(unsigned) * (size_t)n);
kg::run_t<4>(n, (const unsigned **)x, out);
}
void count_5d(int n, const unsigned *x[5], unsigned *out) {
setup();
memset(out, 0, sizeof(unsigned) * (size_t)n);
kg::run_t<5>(n, (const unsigned **)x, out);
}
void count_6d(int n, const unsigned *x[6], unsigned *out) {
setup();
memset(out, 0, sizeof(unsigned) * (size_t)n);
kg::run_t<6>(n, (const unsigned **)x, out);
}
#pragma GCC pop_options
| Compilation | N/A | N/A | Compile OK | Score: N/A | 显示更多 |
| Testcase #1 | 3.191 s | 98 MB + 608 KB | Accepted | Score: 100 | 显示更多 |