From b03b26656da3af88409ca91986a44ff192cad1d7 Mon Sep 17 00:00:00 2001 From: Marshall Lochbaum Date: Sun, 2 Mar 2025 14:16:07 -0500 Subject: [PATCH 01/17] Fix usz=32 build --- src/singeli/src/mask.singeli | 4 ++-- src/singeli/src/scan.singeli | 7 ++++--- 2 files changed, 6 insertions(+), 5 deletions(-) diff --git a/src/singeli/src/mask.singeli b/src/singeli/src/mask.singeli index 5ff8abd2..6fa5c1d9 100644 --- a/src/singeli/src/mask.singeli +++ b/src/singeli/src/mask.singeli @@ -96,7 +96,7 @@ local def ml_exec{i, iter, vars0, bulk, M} = { # i0 - initial batch index; not used as begin because it's in a different scale compared to end def for_masked{bulk, i0}{vars,begin==0,end,iter} = { - l:u64 = end + l:u64 = promote{u64, end} m:u64 = l / bulk @for (i from i0 to m) ml_exec{i, iter, vars, bulk, mask_none} @@ -128,7 +128,7 @@ def for_masked_pos{bulk}{vars,begin==0,end:L,iter} = { # end is scalar element count # index given is a tuple of batch indexes to process def for_mu{bulk, unr, fromunr}{vars,begin==0,end,iter} = { - l:u64 = end + l:u64 = promote{u64, end} m:u64 = l / bulk if (unr==1) { diff --git a/src/singeli/src/scan.singeli b/src/singeli/src/scan.singeli index af3eec65..ac283d4f 100644 --- a/src/singeli/src/scan.singeli +++ b/src/singeli/src/scan.singeli @@ -358,10 +358,11 @@ def loose_mask_gen{V=[vl]T, l} = { # Slow, for ≠` only } def has_vecshift = hasarch{'AVX2'} or hasarch{'AARCH64'} def loose_mask_gen{V=[vl](u64), l if has_vecshift} = { + l64 := promote{u64, l} q := -make{V, 64*iota{vl}} # distance to next row boundary - def q_mod{} = { q+= V**l & -(q>>63) } - def q_mod{if hasarch{'SSE4.1'}} = { q = blend_top{q,q+V**l, q} } - o:u64 = width{V}; while (o>l) { o-=l; q_mod{} } + def q_mod{} = { q+= V**l64 & -(q>>63) } + def q_mod{if hasarch{'SSE4.1'}} = { q = blend_top{q,q+V**l64, q} } + o:u64 = width{V}; while (o>l64) { o-=l64; q_mod{} } {} => { m:= V**1 << q; if (not hasarch{'AVX2'}) m&= q < V**64 q-= V**o; q_mod{} From 90d0f3015af96db3c1df7d7440521f006f597a5a Mon Sep 17 00:00:00 2001 From: Marshall Lochbaum Date: Fri, 28 Feb 2025 13:09:53 -0500 Subject: [PATCH 02/17] Use rank-1 scan special code whenever cell size is 1 --- src/builtins/scan.c | 32 +++++++++++++++++--------------- 1 file changed, 17 insertions(+), 15 deletions(-) diff --git a/src/builtins/scan.c b/src/builtins/scan.c index 128aa726..5e63e8c5 100644 --- a/src/builtins/scan.c +++ b/src/builtins/scan.c @@ -70,7 +70,7 @@ B mul_c2(B, B, B); B scan_ne(B x, u64 p, u64 ia) { // consumes x u64* xp = bitany_ptr(x); - u64* rp; B r=m_bitarrv(&rp,ia); + u64* rp; B r=m_bitarrc(&rp,x); #if SINGELI si_scan_ne(p, xp, rp, BIT_N(ia)); #if USE_VALGRIND @@ -97,14 +97,14 @@ B scan_eq(B x, u64 ia) { // consumes x static B scan_or(B x, u64 ia) { // consumes x u64* xp = bitany_ptr(x); - u64* rp; B r=m_bitarrv(&rp,ia); + u64* rp; B r=m_bitarrc(&rp,x); usz n=BIT_N(ia); u64 xi; usz i=0; while (i,MIN,or ,asc) } #define MM2_ICASE(T,N,C,I) \ case el_##T : { \ if (wv!=(T)wv) { if (wv C 0) { r=C2(shape,m_f64(ia),w); break; } else wv=I; } \ - T* xp=T##any_ptr(x); T* rp; r=m_##T##arrv(&rp, ia); MINMAX_SCAN(T,N,C,wv); \ + T* xp=T##any_ptr(x); T* rp; r=m_##T##arrc(&rp, x); MINMAX_SCAN(T,N,C,wv); \ break; } #define MINMAX2(NAME,C,INIT,BIT,BI,ORD) \ i32 wv=0; if (q_i32(w)) { wv=o2fG(w); } else { x=taga(cpyF64Arr(x)); xe=el_f64; } \ @@ -195,7 +195,7 @@ SHOULD_INLINE B scan2_max_num(B w, B x, u8 xe, usz ia) { MINMAX2(max,>,MIN,or ,0 static B scan_lt(B x, u64 p, usz ia) { u64* xp = bitany_ptr(x); - u64* rp; B r=m_bitarrv(&rp,ia); usz n=BIT_N(ia); + u64* rp; B r=m_bitarrc(&rp,x); usz n=BIT_N(ia); u64 m = 0x5555555555555555; for (usz i=0; if; Arr* r = TI(x,slice)(x, 0, csz); return C2(shape, s, taga(r)); } - if (!(xr==1 && xe<=el_f64)) goto base; + if (xe > el_f64) goto base; + if (xr!=1 && arr_csz(x)!=1) goto base; if (xe==el_bit) switch (rtid) { default: goto base; case n_add: return scan_add_bool(x, ia); // + @@ -278,7 +279,7 @@ B scan_c1(Md1D* d, B x) { B f = d->f; if (!elInt(xe)) goto base; f64 x0 = o2fG(IGetU(x,0)); if (!q_fbit(x0)) goto base; - u64* rp; B r = m_bitarrv(&rp, ia); + u64* rp; B r = m_bitarrc(&rp, x); bool c = x0; rp[0] = c; if (xe==el_i8 ) { i8* xp=i8any_ptr (x); for (usz i=1; if; u8 rtid = RTID(f); if (rtid==n_rtack) { dec(w); return x; } if (rtid==n_ltack) return C2(shape, C1(fne, x), w); - if (!(xr==1 && elNum(xe) && xe<=el_f64)) goto base; + if (!(elNum(xe) && xe<=el_f64)) goto base; + if (xr!=1 && arr_csz(x)!=1) goto base; if (!isF64(w)) goto base; if (rtid==n_floor) return scan2_min_num(w, x, xe, ia); // ⌊ @@ -350,7 +352,7 @@ B scan_c2(Md1D* d, B w, B x) { B f = d->f; if (xe==el_bit) return scan_ne(x, -(u64)(wBit? o2bG(w) : 1&~*bitany_ptr(x)), ia); if (!wBit || !elInt(xe)) goto base; bool c = o2bG(w); - u64* rp; B r = m_bitarrv(&rp, ia); + u64* rp; B r = m_bitarrc(&rp, x); if (xe==el_i8 ) { i8* xp=i8any_ptr (x); for (usz i=0; i Date: Fri, 28 Feb 2025 14:24:58 -0500 Subject: [PATCH 03/17] =?UTF-8?q?Basic=20strided=20=E2=88=A7=E2=88=A8?= =?UTF-8?q?=E2=89=A0=3D=20scan=20implementations?= MIME-Version: 1.0 Content-Type: text/plain; charset=UTF-8 Content-Transfer-Encoding: 8bit --- src/builtins/scan.c | 19 ++++++++++++++++++- src/singeli/src/scan.singeli | 33 +++++++++++++++++++++++++++++++++ 2 files changed, 51 insertions(+), 1 deletion(-) diff --git a/src/builtins/scan.c b/src/builtins/scan.c index 5e63e8c5..d0418dea 100644 --- a/src/builtins/scan.c +++ b/src/builtins/scan.c @@ -258,7 +258,24 @@ B scan_c1(Md1D* d, B x) { B f = d->f; return C2(shape, s, taga(r)); } if (xe > el_f64) goto base; - if (xr!=1 && arr_csz(x)!=1) goto base; + if (xr != 1) { usz csz = arr_csz(x); if (csz != 1) { + #if SINGELI + i8 t = -1; bool neg = 0; + if (xe==el_bit) switch (rtid) { + CASE_N_OR: t=0; break; + CASE_N_AND: t=1; break; + case n_eq: neg=1; case n_ne: t=2; break; + } + if (t != -1) { + if (neg) x = bit_negate(x); + u64* rp; B r=m_bitarrc(&rp,x); + si_scan_bool_stride[t](bitany_ptr(x), rp, ia, csz); + if (neg) r = bit_negate(r); + decG(x); return r; + } + #endif + goto base; + }} if (xe==el_bit) switch (rtid) { default: goto base; case n_add: return scan_add_bool(x, ia); // + diff --git a/src/singeli/src/scan.singeli b/src/singeli/src/scan.singeli index ac283d4f..a3460b96 100644 --- a/src/singeli/src/scan.singeli +++ b/src/singeli/src/scan.singeli @@ -620,3 +620,36 @@ export{'si_scan_rows_and', scan_rows_andor{0}} export{'si_scan_rows_or', scan_rows_andor{1}} export{'si_scan_rows_ne', scan_rows_neq} export{'si_scan_rows_ltack', scan_rows_left} + + +# Strided boolean scans +fn scan_stride_bool_assoc{op}(x:*u64, r:*u64, nl:usz, l:usz) : void = { + assert{l > 1} + def {flip,opf} = if (same{op, &}) tup{~,|} else tup{{x}=>x,op} + nw:= cdiv{nl, 64} + if (l <= 64) { + c:u64 = 0 # carry l bits, no matter the alignment + @for (r, x over nw) { + c = opf{flip{x}, c >> (64-l)} + s:= l; while (s < 64) { c = opf{c, c<>(64-q) | p< Date: Fri, 28 Feb 2025 14:32:43 -0500 Subject: [PATCH 04/17] =?UTF-8?q?Use=20CLMUL-based=20=E2=89=A0`=20for=20po?= =?UTF-8?q?wer-of-two=20stride=20up=20to=2064?= MIME-Version: 1.0 Content-Type: text/plain; charset=UTF-8 Content-Transfer-Encoding: 8bit --- src/singeli/src/scan.singeli | 16 ++++++++++------ 1 file changed, 10 insertions(+), 6 deletions(-) diff --git a/src/singeli/src/scan.singeli b/src/singeli/src/scan.singeli index a3460b96..30cb3fce 100644 --- a/src/singeli/src/scan.singeli +++ b/src/singeli/src/scan.singeli @@ -106,13 +106,13 @@ def vec_prefix_byshift{op, sh} = { def scan_word_ne = prefix_byshift{^, <<} def scan_words_ne = vec_prefix_byshift{^, <<} -fn scan_neq{}(c:u64, x:*u64, r:*u64, nw:u64) : void = { +fn scan_neq{}(c:u64, x:*u64, r:*u64, nw:usz) : void = { @for (x, r over nw) { r = c ^ scan_word_ne{x} c = -(r>>63) # repeat sign bit } } -fn scan_neq{if has_simd}(c0:u64, x:*u64, r:*u64, nw:u64) : void = { +fn scan_neq{if has_simd}(c0:u64, x:*u64, r:*u64, nw:usz) : void = { def vl = arch_defvw / 64 def V = [vl]u64 c := V**c0 @@ -123,7 +123,7 @@ fn scan_neq{if has_simd}(c0:u64, x:*u64, r:*u64, nw:u64) : void = { c = broadcast_last{p} } } -fn clmul_scan_ne_any{if hasarch{'PCLMUL'}}(x:*void, r:*void, init:u64, words:u64, mark:u64) : void = { +fn clmul_scan_ne_any{if hasarch{'PCLMUL'}}(x:*void, r:*void, init:u64, words:usz, mark:u64) : void = { def V = [2]u64 m := V**mark def xor64{a, i, carry} = { # carry is 64-bit broadcasted current total @@ -144,10 +144,10 @@ fn clmul_scan_ne_any{if hasarch{'PCLMUL'}}(x:*void, r:*void, init:u64, words:u64 store{*u64~~(rv+e), clmul{load{V, *u64~~(xv+e), 1}, m, 0} ^ c, 1} } } -fn scan_neq{if hasarch{'PCLMUL'}}(init:u64, x:*u64, r:*u64, nw:u64) : void = { +fn scan_neq{if hasarch{'PCLMUL'}}(init:u64, x:*u64, r:*u64, nw:usz) : void = { clmul_scan_ne_any{}(*void~~x, *void~~r, init, nw, -(u64~~1)) } -fn scan_neq{if hasarch{'AVX512BW', 'VPCLMULQDQ', 'GFNI'}}(init:u64, x:*u64, r:*u64, nw:u64) : void = { +fn scan_neq{if hasarch{'AVX512BW', 'VPCLMULQDQ', 'GFNI'}}(init:u64, x:*u64, r:*u64, nw:usz) : void = { def V = [8]u64 def sse{a} = make{[2]u64, a, 0} carry := sse{init} @@ -561,7 +561,7 @@ fn scan_rows_neq(x:*u64, r:*u64, nl:usz, l:usz) : void = { c:u64 = 0 # carry while (1) { i+= l; ii := iw; iw = cdiv{i, 64} - scan_neq{}(c, x+ii, r+ii, promote{u64,iw-ii}) + scan_neq{}(c, x+ii, r+ii, iw-ii) if (i == nl) return{} s:= load{r, iw-1} q := i%64 @@ -628,6 +628,10 @@ fn scan_stride_bool_assoc{op}(x:*u64, r:*u64, nl:usz, l:usz) : void = { def {flip,opf} = if (same{op, &}) tup{~,|} else tup{{x}=>x,op} nw:= cdiv{nl, 64} if (l <= 64) { + if (same{op, ^} and hasarch{'PCLMUL'} and (l & (l-1)) == 0) { + clmul_scan_ne_any{}(*void~~x, *void~~r, 0, nw, aligned_spaced_mask{l}) + return{} + } c:u64 = 0 # carry l bits, no matter the alignment @for (r, x over nw) { c = opf{flip{x}, c >> (64-l)} From 360e4e8320bb47524a06bf0880d86cec5949fe87 Mon Sep 17 00:00:00 2001 From: Marshall Lochbaum Date: Fri, 28 Feb 2025 14:56:49 -0500 Subject: [PATCH 05/17] Save a word instead of writing to avoid stall for scan with 64>(64-q) | p< Date: Sat, 1 Mar 2025 08:42:36 -0500 Subject: [PATCH 06/17] =?UTF-8?q?Strided=20=E2=8C=8A=E2=8C=88=20scans=20wi?= =?UTF-8?q?th=20shuffles?= MIME-Version: 1.0 Content-Type: text/plain; charset=UTF-8 Content-Transfer-Encoding: 8bit --- src/builtins/scan.c | 7 ++++ src/singeli/src/scan.singeli | 67 ++++++++++++++++++++++++++++++++++++ 2 files changed, 74 insertions(+) diff --git a/src/builtins/scan.c b/src/builtins/scan.c index d0418dea..5de4cf58 100644 --- a/src/builtins/scan.c +++ b/src/builtins/scan.c @@ -273,6 +273,13 @@ B scan_c1(Md1D* d, B x) { B f = d->f; if (neg) r = bit_negate(r); decG(x); return r; } + if (rtid==n_floor | rtid==n_ceil) { + // boolean was handled as CASE_N_AND + B r; void* rp = m_tyarrc(&r, elWidth(xe), x, el2t(xe)); + void* xp = tyany_ptr(x); + si_scan_stride_minmax[4*(rtid==n_ceil) + xe-el_i8](xp, rp, ia, csz); + decG(x); return r; + } #endif goto base; }} diff --git a/src/singeli/src/scan.singeli b/src/singeli/src/scan.singeli index 9420957b..1d98c5bf 100644 --- a/src/singeli/src/scan.singeli +++ b/src/singeli/src/scan.singeli @@ -98,6 +98,73 @@ def shift_first{c:V=[l]_, p:V} = { else blend_first{c, rotate_right{p}} } + +# Strided scans +fn scan_stride_assoc{op, T}(xv:*void, rv:*void, ia:usz, l:usz) : void = { + def id = match (op) { {(min)} => maxvalue; {(max)} => minvalue } + def f = width{T}/8; def vl = 16/f + x:= *T~~xv; r:= *T~~rv + def has_shuf = hasarch{'SSSE3'} or hasarch{'AARCH64'} + if (has_shuf and T<=i32 and l= I**fq; def bl = blend_hom{..., bv} + def V = [vl]T + c:= V**id{T} + o:= l - q + if (l == 2*vl) { o = vl; bv = ~bv } + if (o == vl) { + p:= load{*V~~x}; store{*V~~r, 0, p} + @for_masked{vl} (x in tup{V, x+o}, r in tup{V, r+o} over ia-o) { + p = rot{p} + r = op{bl{c, p}, x} + c = p; p = r + } + } else { + @for_masked{vl} (x in tup{V, x+o}, r in tup{V, r+o}, p in tup{V, r} over ia-o) { + q:= rot{p} + r = op{bl{c, q}, x} + c = q + } + } + } else { + @for (r, x, p in r-l over _ from l to ia) r = op{p, x} + } + } +} +export_tab{'si_scan_stride_minmax', + flat_table{scan_stride_assoc, tup{min,max}, tup{i8,i16,i32,f64}} +} + + # xor scan def vec_prefix_byshift{op, sh} = { def pre{v:V, k} = if (k < elwidth{V}) pre{op{v, sh{v,k}}, 2*k} else v From 180c79e751174ecc3ef46346947f71dbb32b710e Mon Sep 17 00:00:00 2001 From: Marshall Lochbaum Date: Sat, 1 Mar 2025 20:41:58 -0500 Subject: [PATCH 07/17] Strided +`, widening and re-trying on overflow --- src/builtins/scan.c | 17 +++++++++++++++++ src/singeli/src/scan.singeli | 33 ++++++++++++++++++++++++++++----- 2 files changed, 45 insertions(+), 5 deletions(-) diff --git a/src/builtins/scan.c b/src/builtins/scan.c index 5de4cf58..471bd98d 100644 --- a/src/builtins/scan.c +++ b/src/builtins/scan.c @@ -280,6 +280,23 @@ B scan_c1(Md1D* d, B x) { B f = d->f; si_scan_stride_minmax[4*(rtid==n_ceil) + xe-el_i8](xp, rp, ia, csz); decG(x); return r; } + if (rtid==n_add) { + if (xe==el_bit) { x = toI8Any(x); xe=el_i8; } + restart: + B r; void* rp = m_tyarrc(&r, elWidth(xe), x, el2t(xe)); + void* xp = tyany_ptr(x); + bool done = si_scan_stride_add[xe-el_i8](xp, rp, ia, csz); + if (!done) { + decG(r); + switch (++xe) { + case el_i16: x = toI16Any(x); break; + case el_i32: x = toI32Any(x); break; + case el_f64: x = toF64Any(x); break; + } + goto restart; + } + decG(x); return r; + } #endif goto base; }} diff --git a/src/singeli/src/scan.singeli b/src/singeli/src/scan.singeli index 1d98c5bf..10c93e35 100644 --- a/src/singeli/src/scan.singeli +++ b/src/singeli/src/scan.singeli @@ -100,8 +100,11 @@ def shift_first{c:V=[l]_, p:V} = { # Strided scans -fn scan_stride_assoc{op, T}(xv:*void, rv:*void, ia:usz, l:usz) : void = { - def id = match (op) { {(min)} => maxvalue; {(max)} => minvalue } +fn scan_stride_assoc{op, T, Ret, check_over}(xv:*void, rv:*void, ia:usz, l:usz) : Ret = { + def id = match (op) { + {(min)} => maxvalue; {(max)} => minvalue + {(+)} => ({_}=>0) + } def f = width{T}/8; def vl = 16/f x:= *T~~xv; r:= *T~~rv def has_shuf = hasarch{'SSSE3'} or hasarch{'AARCH64'} @@ -112,8 +115,9 @@ fn scan_stride_assoc{op, T}(xv:*void, rv:*void, ia:usz, l:usz) : void = { spr:= I**16 - j + iv inds:= @collect (k) { v:= iv - (j &~ I~~(iv{}} +def check_add_over{w:V=[_]E, x:V, r:V} = { + o:= (if (not hasarch{'X86_64'} or width{E}<=16) any_hom{subs{r,w} != x} + else any_top{(w^r) & (x^r)}) + if (o) return{0} +} +def check_add_over{x, r} = check_add_over{r-x, x, r} export_tab{'si_scan_stride_minmax', flat_table{scan_stride_assoc, tup{min,max}, tup{i8,i16,i32,f64}} } +export_tab{'si_scan_stride_add', tup{ + ...each{scan_stride_assoc{+, ., u1, check_add_over}, tup{i8,i16,i32}}, + scan_stride_assoc{+, f64, u1, {..._}=>{}} +}} # xor scan From 4072bde8066e985653886d788fc8c40f4403f00b Mon Sep 17 00:00:00 2001 From: Marshall Lochbaum Date: Mon, 3 Mar 2025 07:48:48 -0500 Subject: [PATCH 08/17] Strided scan with AVX2 [8]i32 permute if possible --- src/singeli/src/scan.singeli | 70 +++++++++++++++++++++++------------- 1 file changed, 45 insertions(+), 25 deletions(-) diff --git a/src/singeli/src/scan.singeli b/src/singeli/src/scan.singeli index 10c93e35..e0daf81d 100644 --- a/src/singeli/src/scan.singeli +++ b/src/singeli/src/scan.singeli @@ -101,49 +101,69 @@ def shift_first{c:V=[l]_, p:V} = { # Strided scans fn scan_stride_assoc{op, T, Ret, check_over}(xv:*void, rv:*void, ia:usz, l:usz) : Ret = { + def minvalue{(f64)} = -1/0; def maxvalue{(f64)} = 1/0 def id = match (op) { {(min)} => maxvalue; {(max)} => minvalue {(+)} => ({_}=>0) } - def f = width{T}/8; def vl = 16/f x:= *T~~xv; r:= *T~~rv + # Architecture determination + # Use largest vector width with a full-width shuffle def has_shuf = hasarch{'SSSE3'} or hasarch{'AARCH64'} - if (has_shuf and T<=i32 and l=i32) [8]i32 else [16]i8 + def [il]IE = I; def selI = shuf{IE, ...} + def wT = width{T} + def f = wT/width{IE} + def vl = width{I}/wT + def V = [vl]T + if (has_shuf and l < vl) { + # Small stride: power-of-two shifts def small{k} = { - def I = [16]i8 - iv:= iota{I}; j:= I**cast_i{i8,l*f} - spr:= I**16 - j + iv - inds:= @collect (k) { + iv:= iota{I}; j:= I**cast_i{IE,l*f} + spr:= I**il - j + iv + def inds = @collect (k) { v:= iv - (j &~ I~~(iv= js); {x} => selI{x, v} & m } } - def V = [vl]T c:= V**id{T} @for_masked{vl} (x in tup{V, x}, r in tup{V, r} over ia) { - def sc{v, i} = op{shuf{i8, v, i}, v} - r = c = op{shuf{i8, c, spr}, fold{sc, x, inds}} + xs:= fold{{v, i} => op{i{v}, v}, x, inds} + r = c = op{shuf{IE, c, spr}, xs} check_over{x, r} # For +, infers other argument as r-x } } - if (f==1 and l<4) small{3} else small{if (f<=2) 2 else 1} - } else if (has_simd and T==f64 and l==vl) { - def V = [vl]f64 - p:= load{*V~~x}; store{*V~~r, 0, p} - @for (r in *V~~r, x in *V~~x over _ from 1 to ia/vl) { - r = p = op{p, x} + if (not (same{op,+} and V==[4]f64)) { + def max_k = lb{vl/2} + if (max_k<3 or l<4) small{max_k} else small{max_k-1} + } else { # Non-associative! + c:= V**0 + if (l==2) { + @for_masked{vl} (x in tup{V, x}, r in tup{V, r} over ia) { + a:= c + shuf{x, 0,1,0,1} + c = a + shuf{x, 2,3,2,3} + r = blend{a, c, 0,0,1,1} + } + } else { + assert{l==3} + @for_masked{vl} (x in tup{V, x}, r in tup{V, r} over ia) { + a:= shuf{c, 1,1,2,3} + blend{x, V**0, 0,1,1,1} + r = c = x + shuf{a, 1,2,3,0} + } + } } } else { + # Large stride: single shift, with saved register or memory def op_chk{p, x} = { r:= op{p, x}; check_over{p, x, r}; r } @for (r, x over l) r = x - if (has_shuf and T<=i32 and l<256/f) { - def I = [16]i8 - q:= l%vl; fq:= cast_i{i8, f*q} - def rot = shuf{i8, ., (iota{I} - I**fq) & I**15} + if (has_shuf and l<256/(wT/8)) { + def [il]IE = I + q:= l%vl; fq:= cast_i{IE, q*f} + def rot = shuf{IE, ., (iota{I} - I**fq) & I**(il-1)} bv:= iota{I} >= I**fq; def bl = blend_hom{..., bv} - def V = [vl]T c:= V**id{T} o:= l - q if (l == 2*vl) { o = vl; bv = ~bv } @@ -162,7 +182,7 @@ fn scan_stride_assoc{op, T, Ret, check_over}(xv:*void, rv:*void, ia:usz, l:usz) } } } else if (same{op, +} and T<=i32 and has_simd) { - def vl = arch_defvw/width{T}; def V = [vl]T + def vl = arch_defvw/wT; def V = [vl]T @for_masked{vl} (x in tup{V, x+l}, r in tup{V, r+l}, p in tup{V, r} over ia-l) { r = op_chk{p, x} } From 46300bcb51e5f606451ae2ffa201cf3a90d5dcbf Mon Sep 17 00:00:00 2001 From: Marshall Lochbaum Date: Mon, 3 Mar 2025 08:15:47 -0500 Subject: [PATCH 09/17] Use scalar case for small-stride +` on =vl)) { def vl = arch_defvw/wT; def V = [vl]T @for_masked{vl} (x in tup{V, x+l}, r in tup{V, r+l}, p in tup{V, r} over ia-l) { r = op_chk{p, x} } } else { - @for (r, x, p in r-l over _ from l to ia) r = op{p, x} + @for (r, x, p in r-l over _ from l to ia) r = op_chk{p, x} } } 1 } def scan_stride_assoc{op, T} = scan_stride_assoc{op, T, void, {..._}=>{}} +def check_add_over{w:T, x:T, r:T} = { if ((w^r) & (x^r) < 0) return{0} } def check_add_over{w:V=[_]E, x:V, r:V} = { o:= (if (not hasarch{'X86_64'} or width{E}<=16) any_hom{subs{r,w} != x} else any_top{(w^r) & (x^r)}) From 010e6a10d4d2b2312105ea9916f86bb573277629 Mon Sep 17 00:00:00 2001 From: Marshall Lochbaum Date: Mon, 3 Mar 2025 19:56:53 -0500 Subject: [PATCH 10/17] Implementation comments for new high-rank scans --- src/builtins/scan.c | 21 +++++++++++++++------ 1 file changed, 15 insertions(+), 6 deletions(-) diff --git a/src/builtins/scan.c b/src/builtins/scan.c index 471bd98d..3f1b2079 100644 --- a/src/builtins/scan.c +++ b/src/builtins/scan.c @@ -1,22 +1,31 @@ // Scan (`) // Empty 𝕩, and length 1 if no 𝕨: return 𝕩 -// Generic operand: +// Generic argument: // Constant: copy // ⊢ identity, ⊣ reshape 𝕨 or first cell -// Boolean operand, rank 1: +// Boolean argument, stride 1: // + AVX2 expansion (SHOULD have better generic, add SSE, NEON) // ∨⌈ ∧×⌊ search+copy, then memset (COULD vectorize search) // ≠ SWAR/SIMD shifts, CLMUL, VPCLMUL (SHOULD add NEON polynomial mul) // < SWAR // =≤≥>- in terms of ≠<∨∧+ with adjustments -// Arithmetic operand, rank 1: +// Numeric argument, stride 1: // ⌈⌊ Scalar, SIMD in log(vector width) steps // Check in 6-vector blocks to quickly write result if constant // + Overflow-checked scalar or AVX2 // Ad-hoc boolean-valued handling for ≠∨ -// SHOULD extend rank 1 special cases to cell bound 1 -// Higher-rank arithmetic, non-tiny cells: apply operand cell-wise -// SHOULD have dedicated high-rank scan optimizations +// Higher-rank arithmetic: +// Boolean ≠∨∧ and synonyms: SWAR; ⌊⌈+: SIMD with shuffle/permute +// Stride Date: Tue, 4 Mar 2025 04:04:04 +0200 Subject: [PATCH 11/17] fix old gcc --- src/builtins/scan.c | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/src/builtins/scan.c b/src/builtins/scan.c index 3f1b2079..864b2b91 100644 --- a/src/builtins/scan.c +++ b/src/builtins/scan.c @@ -291,7 +291,7 @@ B scan_c1(Md1D* d, B x) { B f = d->f; } if (rtid==n_add) { if (xe==el_bit) { x = toI8Any(x); xe=el_i8; } - restart: + restart:; B r; void* rp = m_tyarrc(&r, elWidth(xe), x, el2t(xe)); void* xp = tyany_ptr(x); bool done = si_scan_stride_add[xe-el_i8](xp, rp, ia, csz); From a11ed3d79a53e62c47e00a167ec136bddcda6103 Mon Sep 17 00:00:00 2001 From: Marshall Lochbaum Date: Tue, 4 Mar 2025 07:38:28 -0500 Subject: [PATCH 12/17] Mask out values past the end for strided +` overflow checks --- src/singeli/src/scan.singeli | 30 +++++++++++++++--------------- 1 file changed, 15 insertions(+), 15 deletions(-) diff --git a/src/singeli/src/scan.singeli b/src/singeli/src/scan.singeli index 47783b0d..9ffe1e97 100644 --- a/src/singeli/src/scan.singeli +++ b/src/singeli/src/scan.singeli @@ -130,10 +130,10 @@ fn scan_stride_assoc{op, T, Ret, check_over}(xv:*void, rv:*void, ia:usz, l:usz) else { m:= V~~(iv >= js); {x} => selI{x, v} & m } } c:= V**id{T} - @for_masked{vl} (x in tup{V, x}, r in tup{V, r} over ia) { + @for_masked{vl} (x in tup{V, x}, r in tup{V, r}, M in 'm' over ia) { xs:= fold{{v, i} => op{i{v}, v}, x, inds} r = c = op{shuf{IE, c, spr}, xs} - check_over{x, r} # For +, infers other argument as r-x + check_over{M, x, r} # For +, infers other argument as r-x } } if (not (same{op,+} and V==[4]f64)) { @@ -157,7 +157,7 @@ fn scan_stride_assoc{op, T, Ret, check_over}(xv:*void, rv:*void, ia:usz, l:usz) } } else { # Large stride: single shift, with saved register or memory - def op_chk{p, x} = { r:= op{p, x}; check_over{p, x, r}; r } + def op_chk{M, p, x} = { r:= op{p, x}; check_over{M, p, x, r}; r } @for (r, x over l) r = x if (has_shuf and l<256/(wT/8)) { def [il]IE = I @@ -169,37 +169,37 @@ fn scan_stride_assoc{op, T, Ret, check_over}(xv:*void, rv:*void, ia:usz, l:usz) if (l == 2*vl) { o = vl; bv = ~bv } if (o == vl) { p:= load{*V~~x}; store{*V~~r, 0, p} - @for_masked{vl} (x in tup{V, x+o}, r in tup{V, r+o} over ia-o) { + @for_masked{vl} (x in tup{V, x+o}, r in tup{V, r+o}, M in 'm' over ia-o) { p = rot{p} - r = op_chk{bl{c, p}, x} + r = op_chk{M, bl{c, p}, x} c = p; p = r } } else { - @for_masked{vl} (x in tup{V, x+o}, r in tup{V, r+o}, p in tup{V, r} over ia-o) { + @for_masked{vl} (x in tup{V, x+o}, r in tup{V, r+o}, p in tup{V, r}, M in 'm' over ia-o) { q:= rot{p} - r = op_chk{bl{c, q}, x} + r = op_chk{M, bl{c, q}, x} c = q } } } else if (same{op, +} and T<=i32 and has_simd and (has_shuf or l>=vl)) { def vl = arch_defvw/wT; def V = [vl]T - @for_masked{vl} (x in tup{V, x+l}, r in tup{V, r+l}, p in tup{V, r} over ia-l) { - r = op_chk{p, x} + @for_masked{vl} (x in tup{V, x+l}, r in tup{V, r+l}, p in tup{V, r}, M in 'm' over ia-l) { + r = op_chk{M, p, x} } } else { - @for (r, x, p in r-l over _ from l to ia) r = op_chk{p, x} + @for (r, x, p in r-l over _ from l to ia) r = op_chk{0, p, x} } } 1 } def scan_stride_assoc{op, T} = scan_stride_assoc{op, T, void, {..._}=>{}} -def check_add_over{w:T, x:T, r:T} = { if ((w^r) & (x^r) < 0) return{0} } -def check_add_over{w:V=[_]E, x:V, r:V} = { - o:= (if (not hasarch{'X86_64'} or width{E}<=16) any_hom{subs{r,w} != x} - else any_top{(w^r) & (x^r)}) +def check_add_over{_, w:T, x:T, r:T} = { if ((w^r) & (x^r) < 0) return{0} } +def check_add_over{M, w:V=[_]E, x:V, r:V} = { + o:= (if (not hasarch{'X86_64'} or width{E}<=16) any_hom{M, subs{r,w} != x} + else any_top{M, (w^r) & (x^r)}) if (o) return{0} } -def check_add_over{x, r} = check_add_over{r-x, x, r} +def check_add_over{M, x, r} = check_add_over{M, r-x, x, r} export_tab{'si_scan_stride_minmax', flat_table{scan_stride_assoc, tup{min,max}, tup{i8,i16,i32,f64}} } From d9a80ee83bba627c32f6544bffb3333e7b80cfd6 Mon Sep 17 00:00:00 2001 From: Marshall Lochbaum Date: Tue, 4 Mar 2025 10:57:28 -0500 Subject: [PATCH 13/17] Test bound and length instead of rank in scan where possible --- src/builtins/scan.c | 18 +++++++++--------- 1 file changed, 9 insertions(+), 9 deletions(-) diff --git a/src/builtins/scan.c b/src/builtins/scan.c index 864b2b91..27ff6a89 100644 --- a/src/builtins/scan.c +++ b/src/builtins/scan.c @@ -243,10 +243,9 @@ static B scan_plus(f64 r0, B x, u8 xe, usz ia) { extern B scan_arith(B f, B w, B x, usz* xsh); // from cells.c B scan_c1(Md1D* d, B x) { B f = d->f; - if (isAtm(x) || RNK(x)==0) thrM("𝔽`𝕩: 𝕩 cannot have rank 0"); - ur xr = RNK(x); - usz ia = IA(x); - if (*SH(x)<=1 || ia==0) return x; + if (isAtm(x)) { unit: thrM("𝔽`𝕩: 𝕩 cannot have rank 0"); } + usz ia = IA(x); if (ia <= 1) { if (ia==1 && RNK(x)==0) goto unit; return x; } + usz n = *SH(x); if (n <= 1) return x; if (RARE(!isFun(f))) { if (isMd(f)) thrM("Calling a modifier"); B xf = getFillR(x); @@ -267,8 +266,9 @@ B scan_c1(Md1D* d, B x) { B f = d->f; return C2(shape, s, taga(r)); } if (xe > el_f64) goto base; - if (xr != 1) { usz csz = arr_csz(x); if (csz != 1) { + if (ia != n) { // csz != 1 #if SINGELI + usz csz = arr_csz(x); i8 t = -1; bool neg = 0; if (xe==el_bit) switch (rtid) { CASE_N_OR: t=0; break; @@ -308,7 +308,7 @@ B scan_c1(Md1D* d, B x) { B f = d->f; } #endif goto base; - }} + } if (xe==el_bit) switch (rtid) { default: goto base; case n_add: return scan_add_bool(x, ia); // + @@ -340,7 +340,7 @@ B scan_c1(Md1D* d, B x) { B f = d->f; if (rtid==n_or) { x=num_squeezeChk(x); xe=TI(x,elType); if (xe==el_bit) return scan_or(x, ia); } } base:; - if (xr>1 && ia >= 6 * (u64)*SH(x) && isPervasiveDy(f)) return scan_arith(f, m_f64(0), x, SH(x)); + if (ia!=n && ia >= 6 * (u64)n && isPervasiveDy(f)) return scan_arith(f, m_f64(0), x, SH(x)); SLOW2("𝕎` 𝕩", f, x); B xf = getFillR(x); @@ -348,7 +348,7 @@ B scan_c1(Md1D* d, B x) { B f = d->f; SGet(x) FC2 fc2 = c2fn(f); - if (xr==1) { + if (ia == n) { r.a[0] = Get(x,0); for (usz i=1; if; if (rtid==n_rtack) { dec(w); return x; } if (rtid==n_ltack) return C2(shape, C1(fne, x), w); if (!(elNum(xe) && xe<=el_f64)) goto base; - if (xr!=1 && arr_csz(x)!=1) goto base; + if (xr!=1 && *SH(x)!=ia) goto base; if (!isF64(w)) goto base; if (rtid==n_floor) return scan2_min_num(w, x, xe, ia); // ⌊ From 26cd40dfd3a3d5407a91771f20439163f2f16235 Mon Sep 17 00:00:00 2001 From: dzaima Date: Tue, 4 Mar 2025 20:11:04 +0200 Subject: [PATCH 14/17] tests for single-element-cell scans --- test/cases/prims.bqn | 3 +++ 1 file changed, 3 insertions(+) diff --git a/test/cases/prims.bqn b/test/cases/prims.bqn index 934dd6f1..758dc5db 100644 --- a/test/cases/prims.bqn +++ b/test/cases/prims.bqn @@ -498,6 +498,9 @@ a←↕2 ⋄ ! "e" ≡ (↕10){b←a‿a‿a‿a‿a‿a‿a‿a‿a‿a ⋄ !"𝕨𝔽`𝕩: Shape of 𝕨 must match the cell of 𝕩 (2‿2 ≡ ≢𝕨, 3‿2‿3 ≡ ≢𝕩)" % (2‿2⥊1)+`↕3‿2‿3 !"𝕨𝔽`𝕩: Shape of 𝕨 must match the cell of 𝕩 (2‿2 ≡ ≢𝕨, 3‿3‿2 ≡ ≢𝕩)" % (2‿2⥊1)+`↕3‿3‿2 !"𝕨𝔽`𝕩: Shape of 𝕨 must match the cell of 𝕩 (⟨⟩ ≡ ≢𝕨, 3‿3 ≡ ≢𝕩)" % 2+`↕3‿3 +%USE eqvar ⋄ %USE k ⋄ ⟨+,-,×,÷,⋆,√,⌊,⌈,¬,∧,∨,<,>,≠,=,≤,≥,⊣,⊢⟩ {f𝕊arr: ! ( F _k`_k ⥊arr) ≡ ⥊_k F` _eqvar arr}⌜ {𝕩∾≍˘¨𝕩} ≍˘¨ ⟨⋈1, 1‿0‿1⟩ +%USE eqvar ⋄ %USE k ⋄ ⟨+,-,×,÷,⋆,√,⌊,⌈,¬,∧,∨,<,>,≠,=,≤,≥,⊣,⊢⟩ {f𝕊arr: ! (1 F _k`_k ⥊arr) ≡ ⥊_k (1¨⊏𝕩)⊸(F`) _eqvar arr}⌜ {𝕩∾≍˘¨𝕩} ≍˘¨ ⟨⋈1, 1‿0‿1⟩ + # ´ !"𝔽´𝕩: 𝕩 must be a list (⟨⟩ ≡ ≢𝕩)" % +´0 From 79fc8d5116230a186ee9255c16127fe706af1766 Mon Sep 17 00:00:00 2001 From: dzaima Date: Tue, 4 Mar 2025 23:18:30 +0200 Subject: [PATCH 15/17] use direct copy on known-different-eltype arrays --- src/builtins/scan.c | 10 +++++----- 1 file changed, 5 insertions(+), 5 deletions(-) diff --git a/src/builtins/scan.c b/src/builtins/scan.c index 27ff6a89..cbc3f74d 100644 --- a/src/builtins/scan.c +++ b/src/builtins/scan.c @@ -290,17 +290,17 @@ B scan_c1(Md1D* d, B x) { B f = d->f; decG(x); return r; } if (rtid==n_add) { - if (xe==el_bit) { x = toI8Any(x); xe=el_i8; } + if (xe==el_bit) { x = taga(cpyI8Arr(x)); xe=el_i8; } restart:; B r; void* rp = m_tyarrc(&r, elWidth(xe), x, el2t(xe)); void* xp = tyany_ptr(x); bool done = si_scan_stride_add[xe-el_i8](xp, rp, ia, csz); if (!done) { decG(r); - switch (++xe) { - case el_i16: x = toI16Any(x); break; - case el_i32: x = toI32Any(x); break; - case el_f64: x = toF64Any(x); break; + switch (++xe) { default: UD; + case el_i16: x = taga(cpyI16Arr(x)); break; + case el_i32: x = taga(cpyI32Arr(x)); break; + case el_f64: x = taga(cpyF64Arr(x)); break; } goto restart; } From 267ddc49875c24bb95d473ff31eb8e62ea752533 Mon Sep 17 00:00:00 2001 From: dzaima Date: Wed, 5 Mar 2025 01:33:35 +0200 Subject: [PATCH 16/17] clarifying comments --- src/singeli/src/scan.singeli | 7 ++++--- 1 file changed, 4 insertions(+), 3 deletions(-) diff --git a/src/singeli/src/scan.singeli b/src/singeli/src/scan.singeli index 9ffe1e97..b105d6c1 100644 --- a/src/singeli/src/scan.singeli +++ b/src/singeli/src/scan.singeli @@ -137,8 +137,8 @@ fn scan_stride_assoc{op, T, Ret, check_over}(xv:*void, rv:*void, ia:usz, l:usz) } } if (not (same{op,+} and V==[4]f64)) { - def max_k = lb{vl/2} - if (max_k<3 or l<4) small{max_k} else small{max_k-1} + def max_k = lb{vl/2} # Divide by two from assuming l≥2 + if (max_k<3 or l<4) small{max_k} else small{max_k-1} # l=2 and l=3 are the only cases needing the full max_k iterations; max_k<3 limits specialization to where it's significant } else { # Non-associative! c:= V**0 if (l==2) { @@ -160,6 +160,7 @@ fn scan_stride_assoc{op, T, Ret, check_over}(xv:*void, rv:*void, ia:usz, l:usz) def op_chk{M, p, x} = { r:= op{p, x}; check_over{M, p, x, r}; r } @for (r, x over l) r = x if (has_shuf and l<256/(wT/8)) { + # Make sure to load the previous row data at the same alignment to not hit bad store-to-load forwarding def [il]IE = I q:= l%vl; fq:= cast_i{IE, q*f} def rot = shuf{IE, ., (iota{I} - I**fq) & I**(il-1)} @@ -736,7 +737,7 @@ export{'si_scan_rows_ltack', scan_rows_left} # Strided boolean scans fn scan_stride_bool_assoc{op}(x:*u64, r:*u64, nl:usz, l:usz) : void = { assert{l > 1} - def {flip,opf} = if (same{op, &}) tup{~,|} else tup{{x}=>x,op} + def {flip,opf} = if (same{op, &}) tup{~,|} else tup{{x}=>x,op} # such that identity of opf is 0 nw:= cdiv{nl, 64} if (l <= 64) { if (same{op, ^} and hasarch{'PCLMUL'} and (l & (l-1)) == 0) { From 6d7983fd4805d79ec47d08d641f599142d33cc9f Mon Sep 17 00:00:00 2001 From: dzaima Date: Thu, 6 Mar 2025 01:29:46 +0200 Subject: [PATCH 17/17] test for scan overflow checking in tail --- test/cases/fuzz/scan.bqn | 2 ++ 1 file changed, 2 insertions(+) create mode 100644 test/cases/fuzz/scan.bqn diff --git a/test/cases/fuzz/scan.bqn b/test/cases/fuzz/scan.bqn new file mode 100644 index 00000000..c901725a --- /dev/null +++ b/test/cases/fuzz/scan.bqn @@ -0,0 +1,2 @@ +# tests scan on a matrix checking for overflow in tail elements +8‿16‿32 {e𝕊n‿w: m←2⋆e-2 ⋄ (m⊸+⌾(⟨¯1,w|n⟩⊸⊑) ⌽‿w⥊(n+w+1)⥊m) ≡ +` •internal.Squeeze ↑‿w⥊ (w⥊m)∾(n⥊0)∾m}⌜ 1+↕200‿65 \ No newline at end of file