From 1e143644b6672c44fed30f821c2fe8902509795c Mon Sep 17 00:00:00 2001 From: Marshall Lochbaum Date: Thu, 8 Sep 2022 20:55:45 -0400 Subject: [PATCH 01/19] =?UTF-8?q?Filter=20out=20=C2=AF1s=20in=20Group=20wi?= =?UTF-8?q?th=20/=20if=20there=20are=20enough=20of=20them?= MIME-Version: 1.0 Content-Type: text/plain; charset=UTF-8 Content-Transfer-Encoding: 8bit --- src/builtins/sfns.c | 23 ++++++++++++++++------- 1 file changed, 16 insertions(+), 7 deletions(-) diff --git a/src/builtins/sfns.c b/src/builtins/sfns.c index 8e1dd6ab..d525ea0f 100644 --- a/src/builtins/sfns.c +++ b/src/builtins/sfns.c @@ -1292,6 +1292,7 @@ B shifta_c2(B t, B w, B x) { return qWithFill(mut_fcd(r, x), f); } +extern B ne_c2(B, B, B); extern B rt_group; B group_c2(B t, B w, B x) { if (isArr(w)&isArr(x) && RNK(w)==1 && RNK(x)==1 && depth(w)==1) { @@ -1303,19 +1304,27 @@ B group_c2(B t, B w, B x) { if (we!=el_i32) w = taga(cpyI32Arr(w)); i32* wp = i32any_ptr(w); i64 ria = wia==xia? 0 : wp[xia]; - if (ria<-1) thrM("โŠ”: ๐•จ can't contain elements less than ยฏ1"); + bool bad = ria < -1; usz neg=0; ria--; - for (usz i = 0; i < xia; i++) if (wp[i]>ria) ria = wp[i]; + for (usz i = 0; i < xia; i++) { + i32 n = wp[i]; + if (n>ria) ria = n; + bad |= n < -1; + neg += n == -1; + } + if (bad) thrM("โŠ”: ๐•จ can't contain elements less than ยฏ1"); if (ria > (i64)(USZ_MAX-1)) thrOOM(); + if (xia>32 && neg>xia/4+xia/8) { + if (wia>xia) w = take_c2(m_f64(0), m_f64(xia), w); + B m = ne_c2(m_f64(0), m_f64(-1), inc(w)); + w = slash_c2(m_f64(0), inc(m), w); if (TI(w,elType)!=el_i32) w = taga(cpyI32Arr(w)); wp = i32any_ptr(w); + x = slash_c2(m_f64(0), m, x); xia = IA(x); + } ria++; TALLOC(i32, lenO, ria+1); i32* len = lenO+1; TALLOC(i32, pos, ria); for (usz i = 0; i < ria; i++) len[i] = pos[i] = 0; - for (usz i = 0; i < xia; i++) { - i32 n = wp[i]; - if (n<-1) thrM("โŠ”: ๐•จ can't contain elements less than ยฏ1"); - len[n]++; // overallocation makes this safe after n<-1 check - } + for (usz i = 0; i < xia; i++) len[wp[i]]++; // overallocation makes this safe after n<-1 check Arr* r = arr_shVec(m_fillarrp(ria)); fillarr_setFill(r, m_f64(0)); B* rp = fillarr_ptr(r); From 569dcacac1aa942568befa40520d70f6bb8c884f Mon Sep 17 00:00:00 2001 From: Marshall Lochbaum Date: Thu, 8 Sep 2022 21:52:34 -0400 Subject: [PATCH 02/19] =?UTF-8?q?Do=20Group's=20range=20check=20and=20?= =?UTF-8?q?=C2=AF1=20counting=20on=201-=20or=202-byte=20ints?= MIME-Version: 1.0 Content-Type: text/plain; charset=UTF-8 Content-Transfer-Encoding: 8bit Doesn't vectorize, but at least it speeds up filtering later --- src/builtins/sfns.c | 36 ++++++++++++++++++++++-------------- 1 file changed, 22 insertions(+), 14 deletions(-) diff --git a/src/builtins/sfns.c b/src/builtins/sfns.c index d525ea0f..417eab4f 100644 --- a/src/builtins/sfns.c +++ b/src/builtins/sfns.c @@ -1301,26 +1301,34 @@ B group_c2(B t, B w, B x) { if (wia-xia > 1) thrF("โŠ”: โ‰ ๐•จ must be either โ‰ ๐•ฉ or one bigger (%sโ‰กโ‰ ๐•จ, %sโ‰กโ‰ ๐•ฉ)", wia, xia); u8 we = TI(w,elType); if (elInt(we)) { - if (we!=el_i32) w = taga(cpyI32Arr(w)); - i32* wp = i32any_ptr(w); - i64 ria = wia==xia? 0 : wp[xia]; - bool bad = ria < -1; usz neg=0; - ria--; - for (usz i = 0; i < xia; i++) { - i32 n = wp[i]; - if (n>ria) ria = n; - bad |= n < -1; - neg += n == -1; - } + if (we==el_bit) w = taga(cpyI8Arr(w)); + i64 ria = 0; + bool bad = false; + usz neg = 0; + void *wp0 = tyany_ptr(w); + #define CASE(T) case el_##T: { \ + T max = -1; \ + for (usz i = 0; i < xia; i++) { \ + T n = ((T*)wp0)[i]; \ + if (n>max) max = n; \ + bad |= n < -1; \ + neg += n == -1; \ + } \ + if (wia>xia) { ria=((T*)wp0)[xia]; bad|=ria<-1; } \ + i64 m=(i64)max+1; if (m>ria) ria=m; \ + break; } + switch (we) { default:UD; case el_bit: CASE(i8) CASE(i16) CASE(i32) } + #undef CASE if (bad) thrM("โŠ”: ๐•จ can't contain elements less than ยฏ1"); - if (ria > (i64)(USZ_MAX-1)) thrOOM(); + if (ria > (i64)(USZ_MAX)) thrOOM(); if (xia>32 && neg>xia/4+xia/8) { if (wia>xia) w = take_c2(m_f64(0), m_f64(xia), w); B m = ne_c2(m_f64(0), m_f64(-1), inc(w)); - w = slash_c2(m_f64(0), inc(m), w); if (TI(w,elType)!=el_i32) w = taga(cpyI32Arr(w)); wp = i32any_ptr(w); + w = slash_c2(m_f64(0), inc(m), w); x = slash_c2(m_f64(0), m, x); xia = IA(x); } - ria++; + if (TI(w,elType)!=el_i32) w = taga(cpyI32Arr(w)); + i32* wp = i32any_ptr(w); TALLOC(i32, lenO, ria+1); i32* len = lenO+1; TALLOC(i32, pos, ria); for (usz i = 0; i < ria; i++) len[i] = pos[i] = 0; From ef2dee9a4dcc45cd57153160a7811657dcf9012d Mon Sep 17 00:00:00 2001 From: Marshall Lochbaum Date: Fri, 9 Sep 2022 13:21:08 -0400 Subject: [PATCH 03/19] =?UTF-8?q?Special=20case=20data=20movement=20for=20?= =?UTF-8?q?=E2=8A=94=20on=20sorted=20=F0=9D=95=A8?= MIME-Version: 1.0 Content-Type: text/plain; charset=UTF-8 Content-Transfer-Encoding: 8bit --- src/builtins/sfns.c | 25 ++++++++++++++++++++++--- 1 file changed, 22 insertions(+), 3 deletions(-) diff --git a/src/builtins/sfns.c b/src/builtins/sfns.c index 417eab4f..43d946ea 100644 --- a/src/builtins/sfns.c +++ b/src/builtins/sfns.c @@ -1303,16 +1303,18 @@ B group_c2(B t, B w, B x) { if (elInt(we)) { if (we==el_bit) w = taga(cpyI8Arr(w)); i64 ria = 0; - bool bad = false; + bool bad = false, sort = true; usz neg = 0; void *wp0 = tyany_ptr(w); #define CASE(T) case el_##T: { \ - T max = -1; \ + T max = -1, prev = -1; \ for (usz i = 0; i < xia; i++) { \ T n = ((T*)wp0)[i]; \ if (n>max) max = n; \ bad |= n < -1; \ neg += n == -1; \ + sort &= prev <= n; \ + prev = n; \ } \ if (wia>xia) { ria=((T*)wp0)[xia]; bad|=ria<-1; } \ i64 m=(i64)max+1; if (m>ria) ria=m; \ @@ -1326,6 +1328,7 @@ B group_c2(B t, B w, B x) { B m = ne_c2(m_f64(0), m_f64(-1), inc(w)); w = slash_c2(m_f64(0), inc(m), w); x = slash_c2(m_f64(0), m, x); xia = IA(x); + neg = 0; } if (TI(w,elType)!=el_i32) w = taga(cpyI32Arr(w)); i32* wp = i32any_ptr(w); @@ -1347,9 +1350,25 @@ B group_c2(B t, B w, B x) { case el_i16: case el_c16: case el_i32: case el_c32: case el_f64: { u8 width = elWidth(xe); - for (usz i = 0; i < ria; i++) m_tyarrv(rp+i, width, len[i], el2t(xe)); void* xp = tyany_ptr(x); + B z = taga(rf); + u8 xt = el2t(xe); + if (sort) { + for (usz j=0, i=neg*width; j=0) ((u8* )tyarr_ptr(rp[n]))[pos[n]++] = ((u8* )xp)[i]; } break; case 2: for (usz i = 0; i < xia; i++) { i32 n = wp[i]; if (n>=0) ((u16*)tyarr_ptr(rp[n]))[pos[n]++] = ((u16*)xp)[i]; } break; From d8abc4ef7e0de7ab8edcf3e02410436c82e256a9 Mon Sep 17 00:00:00 2001 From: Marshall Lochbaum Date: Fri, 9 Sep 2022 18:11:12 -0400 Subject: [PATCH 04/19] =?UTF-8?q?Separate=20handling=20for=20=E2=8A=94=20w?= =?UTF-8?q?hen=20=F0=9D=95=A8=20changes=20value=20rarely?= MIME-Version: 1.0 Content-Type: text/plain; charset=UTF-8 Content-Transfer-Encoding: 8bit --- src/builtins/sfns.c | 78 +++++++++++++++++++++++++++++++++++++-------- 1 file changed, 64 insertions(+), 14 deletions(-) diff --git a/src/builtins/sfns.c b/src/builtins/sfns.c index 43d946ea..0a8109a6 100644 --- a/src/builtins/sfns.c +++ b/src/builtins/sfns.c @@ -1293,6 +1293,7 @@ B shifta_c2(B t, B w, B x) { } extern B ne_c2(B, B, B); +extern B select_c2(B, B, B); extern B rt_group; B group_c2(B t, B w, B x) { if (isArr(w)&isArr(x) && RNK(w)==1 && RNK(x)==1 && depth(w)==1) { @@ -1304,7 +1305,7 @@ B group_c2(B t, B w, B x) { if (we==el_bit) w = taga(cpyI8Arr(w)); i64 ria = 0; bool bad = false, sort = true; - usz neg = 0; + usz neg = 0, change = 0; void *wp0 = tyany_ptr(w); #define CASE(T) case el_##T: { \ T max = -1, prev = -1; \ @@ -1314,6 +1315,7 @@ B group_c2(B t, B w, B x) { bad |= n < -1; \ neg += n == -1; \ sort &= prev <= n; \ + change += prev != n; \ prev = n; \ } \ if (wia>xia) { ria=((T*)wp0)[xia]; bad|=ria<-1; } \ @@ -1323,6 +1325,65 @@ B group_c2(B t, B w, B x) { #undef CASE if (bad) thrM("โŠ”: ๐•จ can't contain elements less than ยฏ1"); if (ria > (i64)(USZ_MAX)) thrOOM(); + + Arr* r = arr_shVec(m_fillarrp(ria)); fillarr_setFill(r, m_f64(0)); + B* rp = fillarr_ptr(r); + for (usz i = 0; i < ria; i++) rp[i] = m_f64(0); // don't break if allocation errors + B xf = getFillQ(x); + + Arr* rf = arr_shVec(m_fillarrp(0)); fillarr_setFill(rf, m_f64(0)); + B z = taga(rf); + fillarr_setFill(r, z); + + TALLOC(i32, pos, 2*ria+1); i32* len = pos+ria+1; + // Both cases needed to make sure wia>0 for ip[wia-1] below + if (ria==0) goto intvec_ret; + if (neg==xia) { + for (usz i = 0; i < ria; i++) rp[i] = inc(z); + goto intvec_ret; + } + + u8 xe = TI(x,elType); + u8 width = elWidth(xe); + u64 xw; + if (xia>64 && (xw=(u64)xia*width)<=I32_MAX && changexia) w = C2(take, m_f64(xia), w); + B c = C2(ne, C2(drop, m_f64(-1), inc(w)), + C2(drop, m_f64( 1), inc(w))); + B ind = C1(slash, C2(join, m_f64(-1!=IGetU(w,0).f), c)); + w = C2(select, inc(ind), w); + #undef C1 + #undef C2 + if (TI(ind,elType)!=el_i32) ind = taga(cpyI32Arr(ind)); + if (TI(w ,elType)!=el_i32) w = taga(cpyI32Arr(w )); + wia = IA(ind); + + i32* ip = i32any_ptr(ind); + i32* wp = i32any_ptr(w); + usz i0 = ip[0]; + for (usz i=0; i32 && neg>xia/4+xia/8) { if (wia>xia) w = take_c2(m_f64(0), m_f64(xia), w); B m = ne_c2(m_f64(0), m_f64(-1), inc(w)); @@ -1332,26 +1393,14 @@ B group_c2(B t, B w, B x) { } if (TI(w,elType)!=el_i32) w = taga(cpyI32Arr(w)); i32* wp = i32any_ptr(w); - TALLOC(i32, lenO, ria+1); i32* len = lenO+1; - TALLOC(i32, pos, ria); for (usz i = 0; i < ria; i++) len[i] = pos[i] = 0; for (usz i = 0; i < xia; i++) len[wp[i]]++; // overallocation makes this safe after n<-1 check - Arr* r = arr_shVec(m_fillarrp(ria)); fillarr_setFill(r, m_f64(0)); - B* rp = fillarr_ptr(r); - for (usz i = 0; i < ria; i++) rp[i] = m_f64(0); // don't break if allocation errors - B xf = getFillQ(x); - - Arr* rf = arr_shVec(m_fillarrp(0)); fillarr_setFill(rf, m_f64(0)); - fillarr_setFill(r, taga(rf)); - u8 xe = TI(x,elType); switch (xe) { default: UD; case el_i8: case el_c8: case el_i16: case el_c16: case el_i32: case el_c32: case el_f64: { - u8 width = elWidth(xe); void* xp = tyany_ptr(x); - B z = taga(rf); u8 xt = el2t(xe); if (sort) { for (usz j=0, i=neg*width; j Date: Sat, 10 Sep 2022 08:57:41 -0400 Subject: [PATCH 05/19] Fast path for / with empty result, and faster non-Singeli Where --- src/builtins/sfns.c | 9 ++++++--- 1 file changed, 6 insertions(+), 3 deletions(-) diff --git a/src/builtins/sfns.c b/src/builtins/sfns.c index 0a8109a6..1079734e 100644 --- a/src/builtins/sfns.c +++ b/src/builtins/sfns.c @@ -429,6 +429,7 @@ B slash_c1(B t, B x) { if (RARE(isAtm(x)) || RARE(RNK(x)!=1)) thrF("/: Argument must have rank 1 (%H โ‰ก โ‰ข๐•ฉ)", x); u64 s = usum(x); if (s>=USZ_MAX) thrOOM(); + if (s==0) { decG(x); return emptyIVec(); } usz xia = IA(x); if (RARE(xia>=I32_MAX)) { usz xia = IA(x); @@ -471,9 +472,11 @@ B slash_c1(B t, B x) { if (xe==el_bit) { u64* xp = bitarr_ptr(x); while (xia>0 && !bitp_get(xp,xia-1)) xia--; - for (u64 i = 0; i < xia; i++) { - *rp = i; - rp+= bitp_get(xp, i); + u8* x8 = (u8*)xp; + u8 q=xia%8; if (q) x8[xia/8] &= (1<>=1; } } } else if (xe==el_i8) { i8* xp = i8any_ptr(x); From bb3ecdca6ab2dc765e9e2908ade9b0102ed2c369 Mon Sep 17 00:00:00 2001 From: Marshall Lochbaum Date: Sat, 10 Sep 2022 09:12:16 -0400 Subject: [PATCH 06/19] Sparse Where --- src/builtins/sfns.c | 31 +++++++++++++++++++++++-------- src/h.h | 1 + 2 files changed, 24 insertions(+), 8 deletions(-) diff --git a/src/builtins/sfns.c b/src/builtins/sfns.c index 1079734e..297e228b 100644 --- a/src/builtins/sfns.c +++ b/src/builtins/sfns.c @@ -444,9 +444,23 @@ B slash_c1(B t, B x) { } B r; u8 xe = TI(x,elType); - #if SINGELI && defined(__BMI2__) if (xe==el_bit) { u64* xp = bitarr_ptr(x); + // Sparse method with CTZ + #if SINGELI && defined(__BMI2__) + if (xia>128 && s < xia/8+(xia<=32768?xia/4:0)) { + #else + if (xia<=128 || s < xia/2+(xia<=32768?xia/4:0)) { + #endif + usz q=xia%64; if (q) xp[xia/64] &= ((u64)1<0 && !bitp_get(xp,xia-1)) xia--; u8* x8 = (u8*)xp; u8 q=xia%8; if (q) x8[xia/8] &= (1<>=1; } } - } else if (xe==el_i8) { + } + #endif + } else { + i32* rp; r = m_i32arrv(&rp, s); + if (xe==el_i8) { i8* xp = i8any_ptr(x); while (xia>0 && !xp[xia-1]) xia--; for (u64 i = 0; i < xia; i++) { diff --git a/src/h.h b/src/h.h index 4afb5a87..7fb6509e 100644 --- a/src/h.h +++ b/src/h.h @@ -130,6 +130,7 @@ typedef double f64; #define NORETURN __attribute__((noreturn)) #define AUTO __auto_type #define CLZ(X) __builtin_clzll(X) +#define CTZ(X) __builtin_ctzll(X) #define POPC(X) __builtin_popcountll(X) #define LIKELY(X) __builtin_expect(X,1) #define RARE(X) __builtin_expect(X,0) From 017fcaea2acb684f115db5996cd12a087d847452 Mon Sep 17 00:00:00 2001 From: Marshall Lochbaum Date: Sat, 10 Sep 2022 09:58:42 -0400 Subject: [PATCH 07/19] Move slash functions to a separate file --- makefile | 2 +- src/builtins/sfns.c | 506 +----------------------------------------- src/builtins/slash.c | 507 +++++++++++++++++++++++++++++++++++++++++++ src/opt/single.c | 1 + 4 files changed, 512 insertions(+), 504 deletions(-) create mode 100644 src/builtins/slash.c diff --git a/makefile b/makefile index 7d39a8e1..a70156b7 100644 --- a/makefile +++ b/makefile @@ -198,7 +198,7 @@ ${bd}/%.o: src/jit/%.c @echo $< | cut -c 5- @$(CC_INC) $@.d -o $@ -c $< -builtins: ${addprefix ${bd}/, arithm.o arithd.o cmp.o sfns.o squeeze.o select.o sort.o selfsearch.o md1.o md2.o fns.o sysfn.o internal.o inverse.o} +builtins: ${addprefix ${bd}/, arithm.o arithd.o cmp.o sfns.o squeeze.o select.o slash.o sort.o selfsearch.o md1.o md2.o fns.o sysfn.o internal.o inverse.o} ${bd}/%.o: src/builtins/%.c @echo $< | cut -c 5- @$(CC_INC) $@.d -o $@ -c $< diff --git a/src/builtins/sfns.c b/src/builtins/sfns.c index 297e228b..13dc8910 100644 --- a/src/builtins/sfns.c +++ b/src/builtins/sfns.c @@ -323,509 +323,6 @@ B pick_c2(B t, B w, B x) { return r; } -#ifdef __BMI2__ - #include - - #if USE_VALGRIND - #define DBG_VG_SLASH 0 - u64 loadMask(u64* p, u64 unk, u64 exp, u64 i, u64 pos) { - // #if DBG_VG_SLASH - // if (pos==0) printf("index %2ld, got %016lx\n", i, p[i]); - // #endif - if (pos==0) return ~(p[i]^exp); - u64 res = loadMask(p, unk, exp, i, pos<<1); - if (unk&pos) res&= loadMask(p, unk, exp, i|pos, pos<<1); - return res; - } - NOINLINE u64 vg_load64(u64* p, u64 i) { - u64 unk = ~vg_getDefined_u64(i); - u64 res = p[vg_withDefined_u64(i, ~0ULL)]; // result value will always be the proper indexing operation - - i32 undefCount = POPC(unk); - if (undefCount>0) { - if (undefCount>8) err("too many unknown bits in index of vg_load64"); - res = vg_withDefined_u64(res, loadMask(p, unk, res, i & ~unk, 1)); - } - #if DBG_VG_SLASH - vg_printDefined_u64("idx", i); - vg_printDefined_u64("res", res); - #endif - return res; - } - NOINLINE u64 vg_pext_u64(u64 src, u64 mask) { - u64 maskD = vg_getDefined_u64(mask); - u64 r = vg_undef_u64(0); - i32 ri = 0; - u64 undefMask = 0; - for (i32 i = 0; i < 64; i++) { - u64 c = 1ull<>i)&1) { - r|= (c&1) << i; - c>>= 1; - } - } - #if DBG_VG_SLASH - printf("pdep:\n"); - vg_printDefined_u64("src", src); - vg_printDefined_u64("msk", mask); - vg_printDefined_u64("res", r); - vg_printDefined_u64("exp", _pdep_u64(src, mask)); - #endif - return r; - } - NOINLINE u64 rand_popc64(u64 x) { - u64 def = vg_getDefined_u64(x); - if (def==~0ULL) return POPC(x); - i32 min = POPC(x & def); - i32 diff = POPC(~def); - i32 res = min + vgRand64Range(diff); - #if DBG_VG_SLASH - printf("popc:\n"); - vg_printDefined_u64("x", x); - printf("popc in %d-%d; res: %d\n", min, min+diff, res); - #endif - return res; - } - #define _pext_u32 vg_pext_u64 - #define _pext_u64 vg_pext_u64 - #define _pdep_u32 vg_pdep_u64 - #define _pdep_u64 vg_pdep_u64 - #else - #define vg_load64(p, i) p[i] - #define rand_popc64(X) POPC(X) - #endif - - void storeu_u64(u64* p, u64 v) { memcpy(p, &v, 8); } - u64 loadu_u64(u64* p) { u64 v; memcpy(&v, p, 8); return v; } - #if SINGELI - #pragma GCC diagnostic push - #pragma GCC diagnostic ignored "-Wunused-variable" - #include "../singeli/gen/slash.c" - #pragma GCC diagnostic pop - #endif -#endif - -extern B rt_slash; -B slash_c1(B t, B x) { - if (RARE(isAtm(x)) || RARE(RNK(x)!=1)) thrF("/: Argument must have rank 1 (%H โ‰ก โ‰ข๐•ฉ)", x); - u64 s = usum(x); - if (s>=USZ_MAX) thrOOM(); - if (s==0) { decG(x); return emptyIVec(); } - usz xia = IA(x); - if (RARE(xia>=I32_MAX)) { - usz xia = IA(x); - SGetU(x) - f64* rp; B r = m_f64arrv(&rp, s); usz ri = 0; - for (usz i = 0; i < xia; i++) { - usz c = o2s(GetU(x, i)); - for (usz j = 0; j < c; j++) rp[ri++] = i; - } - decG(x); - return r; - } - B r; - u8 xe = TI(x,elType); - if (xe==el_bit) { - u64* xp = bitarr_ptr(x); - // Sparse method with CTZ - #if SINGELI && defined(__BMI2__) - if (xia>128 && s < xia/8+(xia<=32768?xia/4:0)) { - #else - if (xia<=128 || s < xia/2+(xia<=32768?xia/4:0)) { - #endif - usz q=xia%64; if (q) xp[xia/64] &= ((u64)1<0 && !bitp_get(xp,xia-1)) xia--; - u8* x8 = (u8*)xp; - u8 q=xia%8; if (q) x8[xia/8] &= (1<>=1; } - } - } - #endif - } else { - i32* rp; r = m_i32arrv(&rp, s); - if (xe==el_i8) { - i8* xp = i8any_ptr(x); - while (xia>0 && !xp[xia-1]) xia--; - for (u64 i = 0; i < xia; i++) { - i32 c = xp[i]; - if (LIKELY(c==0 || c==1)) { - *rp = i; - rp+= c; - } else { - for (i32 j = 0; j < c; j++) *rp++ = i; - } - } - } else if (xe==el_i32) { - i32* xp = i32any_ptr(x); - while (xia>0 && !xp[xia-1]) xia--; - for (u64 i = 0; i < xia; i++) { - i32 c = xp[i]; - if (LIKELY(c==0 || c==1)) { - *rp = i; - rp+= c; - } else { - for (i32 j = 0; j < c; j++) *rp++ = i; - } - } - } else { - SLOW1("/๐•ฉ", x); - SGetU(x) - for (u64 i = 0; i < xia; i++) { - usz c = o2s(GetU(x, i)); - for (u64 j = 0; j < c; j++) *rp++ = i; - } - } - } - decG(x); - return r; -} - -B slash_c2(B t, B w, B x) { - if (isArr(x) && RNK(x)==1 && isArr(w) && RNK(w)==1 && depth(w)==1) { - usz wia = IA(w); - usz xia = IA(x); - if (RARE(wia!=xia)) { - if (wia==0) { decG(w); return x; } - thrF("/: Lengths of components of ๐•จ must match ๐•ฉ (%s โ‰  %s)", wia, xia); - } - B xf = getFillQ(x); - - usz ri = 0; - if (TI(w,elType)==el_bit) { - B r; - u64* wp = bitarr_ptr(w); - u8 xe = TI(x,elType); - #ifdef __BMI2__ - usz wsum = bit_sum(wp, wia); - switch(xe) { - case el_bit: { - u64* xp = bitarr_ptr(x); u64* rp; r = m_bitarrv(&rp,wsum+128); a(r)->ia = wsum; - u64 cw = 0; // current word - u64 ro = 0; // offset in word where next bit should be written; never 64 - for (usz i=0; i=64) { - *(rp++) = cw; - cw = ro? v>>(64-ro) : 0; - } - ro = ro2&63; - } - if (ro) *rp = cw; - goto bit_ret; - } - #if SINGELI - case el_i8: case el_c8: { i8* xp=tyany_ptr(x); i8* rp=m_tyarrvO(&r,1,wsum,el2t(xe), 8); bmipopc_2slash8 (wp, xp, rp, wia); goto bit_ret; } - case el_i16:case el_c16: { i16* xp=tyany_ptr(x); i16* rp=m_tyarrvO(&r,2,wsum,el2t(xe), 16); bmipopc_2slash16(wp, xp, rp, wia); goto bit_ret; } - case el_i32:case el_c32: { - i32* xp=tyany_ptr(x); i32* rp=m_tyarrv(&r,4,wsum,el2t(xe)); - usz b = 1<<7; - TALLOC(i8, buf, b); - i32* rq=rp; i32* end=xp+xia-b; - while (xp < end) { - bmipopc_1slash8(wp, buf, b); - usz bs = bit_sum(wp, b); - for (usz j=0; j0 && !bitp_get(wp,wia-1)) wia--; - #ifndef __BMI2__ - usz wsum = bit_sum(wp, wia); - #endif - if (wsum==0) { decG(w); decG(x); return q_N(xf)? emptyHVec() : isF64(xf)? emptyIVec() : isC32(xf)? emptyCVec() : m_emptyFVec(xf); } - switch(xe) { default: UD; - - #if !SINGELI - case el_i8: case el_c8: { i8* xp=tyany_ptr(x); i8* rp=m_tyarrv(&r,1,wsum,el2t(xe)); for (usz i=0; i0 && !wp[wia-1]) wia--; \ - i64 wsum = 0; \ - u32 or = 0; \ - for (usz i = 0; i < wia; i++) { \ - wsum+= wp[i]; \ - or|= (u32)wp[i]; \ - } \ - if (or>>SIGN) thrM("/: ๐•จ must consist of natural numbers"); \ - if (TI(x,elType)==el_bit) { \ - u64* xp = bitarr_ptr(x); u64 ri=0; \ - u64* rp; B r = m_bitarrv(&rp, wsum); \ - if (or<2) for (usz i = 0; i < wia; i++) { \ - bitp_set(rp, ri, bitp_get(xp,i)); \ - ri+= wp[i]; \ - } else for (usz i = 0; i < wia; i++) { \ - WT cw = wp[i]; bool cx = bitp_get(xp,i); \ - for (i64 j = 0; j < cw; j++) bitp_set(rp, ri++, cx); \ - } \ - decG(w); decG(x); return r; \ - } \ - CASE(WT,i8) CASE(WT,i16) CASE(WT,i32) CASE(WT,f64) \ - SLOW2("๐•จ/๐•ฉ", w, x); \ - M_HARR(r, wsum) SGetU(x) \ - for (usz i = 0; i < wia; i++) { \ - i32 cw = wp[i]; if (cw==0) continue; \ - B cx = incBy(GetU(x, i), cw); \ - for (i64 j = 0; j < cw; j++) HARR_ADDA(r, cx);\ - } \ - decG(w); decG(x); \ - return withFill(HARR_FV(r), xf); \ - } - if (TI(w,elType)==el_i8 ) TYPED(i8,7); - if (TI(w,elType)==el_i32) TYPED(i32,31); - #undef TYPED - #undef CASE - SLOW2("๐•จ/๐•ฉ", w, x); - u64 ria = usum(w); - if (ria>=USZ_MAX) thrOOM(); - M_HARR(r, ria) SGetU(w) SGetU(x) - for (usz i = 0; i < wia; i++) { - usz c = o2s(GetU(w, i)); - if (c) { - B cx = incBy(GetU(x, i), c); - for (usz j = 0; RARE(j < c); j++) HARR_ADDA(r, cx); - } - } - decG(w); decG(x); - return withFill(HARR_FV(r), xf); - } - if (isArr(x) && RNK(x)==1 && q_i32(w)) { - usz xia = IA(x); - i32 wv = o2i(w); - if (wv<=0) { - if (wv<0) thrM("/: ๐•จ cannot be negative"); - return taga(arr_shVec(TI(x,slice)(x, 0, 0))); - } - if (TI(x,elType)==el_i32) { - i32* xp = i32any_ptr(x); - i32* rp; B r = m_i32arrv(&rp, xia*wv); - for (usz i = 0; i < xia; i++) { - for (i64 j = 0; j < wv; j++) *rp++ = xp[i]; - } - decG(x); - return r; - } else { - SLOW2("๐•จ/๐•ฉ", w, x); - B xf = getFillQ(x); - HArr_p r = m_harrUv(xia*wv); - SGetU(x) - for (usz i = 0; i < xia; i++) { - B cx = incBy(GetU(x, i), wv); - for (i64 j = 0; j < wv; j++) *r.a++ = cx; - } - decG(x); - return withFill(r.b, xf); - } - } - return c2(rt_slash, w, x); -} - - -B slash_im(B t, B x) { - if (!isArr(x) || RNK(x)!=1) thrM("/โผ: Argument must be an array"); - u8 xe = TI(x,elType); - usz xia = IA(x); - if (xia==0) { decG(x); return emptyIVec(); } - switch(xe) { default: UD; - case el_bit: { - usz sum = bit_sum(bitarr_ptr(x), xia); - usz ria = 1 + (sum>0); - f64* rp; B r = m_f64arrv(&rp, ria); - rp[sum>0] = sum; rp[0] = xia - sum; - decG(x); return num_squeeze(r); - } -#define CASE_SMALL(N) \ - case el_i##N: { \ - i##N* xp = i##N##any_ptr(x); \ - usz m=1<xp[a-1]) a++; \ - max=xp[a-1]; \ - if (a==xia) { /* Sorted unique argument */ \ - usz ria = max + 1; \ - u64* rp; r = m_bitarrv(&rp, ria); \ - for (usz i=0; imax) max=c; } \ - if ((i##N)max<0) thrM("/โผ: Argument cannot contain negative numbers"); \ - usz ria = max+1; \ - i##N* rp; r = m_i##N##arrv(&rp, ria); for (usz i=0; im/2) thrM("/โผ: Argument cannot contain negative numbers"); \ - i32* rp; r = m_i32arrv(&rp, ria); for (usz i=0; imax?c:max; if (c<0) thrM("/โผ: Argument cannot contain negative numbers"); } - usz ria = max+1; - if (i==xia) { - u64* rp; r = m_bitarrv(&rp, ria); for (usz i=0; imax?c:max; if (c<0) thrM("/โผ: Argument cannot contain negative numbers"); } - usz ria = max+1; if (ria==0) thrOOM(); - if (i==xia) { - u64* rp; r = m_bitarrv(&rp, ria); for (usz i=0; ia; } - usz i,j; B r; i64 max=-1; - for (i = 0; i < xia; i++) { i64 c=o2i64(xp[i]); if (c<=max) break; max=c; } - for (j = i; j < xia; j++) { i64 c=o2i64(xp[j]); max=c>max?c:max; if (c<0) thrM("/โผ: Argument cannot contain negative numbers"); } - if (max > USZ_MAX-1) thrOOM(); - usz ria = max+1; - if (i==xia) { - u64* rp; r = m_bitarrv(&rp, ria); for (usz i=0; iuc1 = pick_uc1; c(BFn,bi_reverse)->uc1 = reverse_uc1; diff --git a/src/builtins/slash.c b/src/builtins/slash.c new file mode 100644 index 00000000..4518689a --- /dev/null +++ b/src/builtins/slash.c @@ -0,0 +1,507 @@ +#include "../core.h" +#include "../utils/mut.h" +#include "../utils/talloc.h" +#include "../builtins.h" + +#ifdef __BMI2__ + #include + + #if USE_VALGRIND + #define DBG_VG_SLASH 0 + u64 loadMask(u64* p, u64 unk, u64 exp, u64 i, u64 pos) { + // #if DBG_VG_SLASH + // if (pos==0) printf("index %2ld, got %016lx\n", i, p[i]); + // #endif + if (pos==0) return ~(p[i]^exp); + u64 res = loadMask(p, unk, exp, i, pos<<1); + if (unk&pos) res&= loadMask(p, unk, exp, i|pos, pos<<1); + return res; + } + NOINLINE u64 vg_load64(u64* p, u64 i) { + u64 unk = ~vg_getDefined_u64(i); + u64 res = p[vg_withDefined_u64(i, ~0ULL)]; // result value will always be the proper indexing operation + + i32 undefCount = POPC(unk); + if (undefCount>0) { + if (undefCount>8) err("too many unknown bits in index of vg_load64"); + res = vg_withDefined_u64(res, loadMask(p, unk, res, i & ~unk, 1)); + } + #if DBG_VG_SLASH + vg_printDefined_u64("idx", i); + vg_printDefined_u64("res", res); + #endif + return res; + } + NOINLINE u64 vg_pext_u64(u64 src, u64 mask) { + u64 maskD = vg_getDefined_u64(mask); + u64 r = vg_undef_u64(0); + i32 ri = 0; + u64 undefMask = 0; + for (i32 i = 0; i < 64; i++) { + u64 c = 1ull<>i)&1) { + r|= (c&1) << i; + c>>= 1; + } + } + #if DBG_VG_SLASH + printf("pdep:\n"); + vg_printDefined_u64("src", src); + vg_printDefined_u64("msk", mask); + vg_printDefined_u64("res", r); + vg_printDefined_u64("exp", _pdep_u64(src, mask)); + #endif + return r; + } + NOINLINE u64 rand_popc64(u64 x) { + u64 def = vg_getDefined_u64(x); + if (def==~0ULL) return POPC(x); + i32 min = POPC(x & def); + i32 diff = POPC(~def); + i32 res = min + vgRand64Range(diff); + #if DBG_VG_SLASH + printf("popc:\n"); + vg_printDefined_u64("x", x); + printf("popc in %d-%d; res: %d\n", min, min+diff, res); + #endif + return res; + } + #define _pext_u32 vg_pext_u64 + #define _pext_u64 vg_pext_u64 + #define _pdep_u32 vg_pdep_u64 + #define _pdep_u64 vg_pdep_u64 + #else + #define vg_load64(p, i) p[i] + #define rand_popc64(X) POPC(X) + #endif + + void storeu_u64(u64* p, u64 v) { memcpy(p, &v, 8); } + u64 loadu_u64(u64* p) { u64 v; memcpy(&v, p, 8); return v; } + #if SINGELI + #pragma GCC diagnostic push + #pragma GCC diagnostic ignored "-Wunused-variable" + #include "../singeli/gen/slash.c" + #pragma GCC diagnostic pop + #endif +#endif + +extern B rt_slash; +B slash_c1(B t, B x) { + if (RARE(isAtm(x)) || RARE(RNK(x)!=1)) thrF("/: Argument must have rank 1 (%H โ‰ก โ‰ข๐•ฉ)", x); + u64 s = usum(x); + if (s>=USZ_MAX) thrOOM(); + if (s==0) { decG(x); return emptyIVec(); } + usz xia = IA(x); + if (RARE(xia>=I32_MAX)) { + usz xia = IA(x); + SGetU(x) + f64* rp; B r = m_f64arrv(&rp, s); usz ri = 0; + for (usz i = 0; i < xia; i++) { + usz c = o2s(GetU(x, i)); + for (usz j = 0; j < c; j++) rp[ri++] = i; + } + decG(x); + return r; + } + B r; + u8 xe = TI(x,elType); + if (xe==el_bit) { + u64* xp = bitarr_ptr(x); + // Sparse method with CTZ + #if SINGELI && defined(__BMI2__) + if (xia>128 && s < xia/8+(xia<=32768?xia/4:0)) { + #else + if (xia<=128 || s < xia/2+(xia<=32768?xia/4:0)) { + #endif + usz q=xia%64; if (q) xp[xia/64] &= ((u64)1<0 && !bitp_get(xp,xia-1)) xia--; + u8* x8 = (u8*)xp; + u8 q=xia%8; if (q) x8[xia/8] &= (1<>=1; } + } + } + #endif + } else { + i32* rp; r = m_i32arrv(&rp, s); + if (xe==el_i8) { + i8* xp = i8any_ptr(x); + while (xia>0 && !xp[xia-1]) xia--; + for (u64 i = 0; i < xia; i++) { + i32 c = xp[i]; + if (LIKELY(c==0 || c==1)) { + *rp = i; + rp+= c; + } else { + for (i32 j = 0; j < c; j++) *rp++ = i; + } + } + } else if (xe==el_i32) { + i32* xp = i32any_ptr(x); + while (xia>0 && !xp[xia-1]) xia--; + for (u64 i = 0; i < xia; i++) { + i32 c = xp[i]; + if (LIKELY(c==0 || c==1)) { + *rp = i; + rp+= c; + } else { + for (i32 j = 0; j < c; j++) *rp++ = i; + } + } + } else { + SLOW1("/๐•ฉ", x); + SGetU(x) + for (u64 i = 0; i < xia; i++) { + usz c = o2s(GetU(x, i)); + for (u64 j = 0; j < c; j++) *rp++ = i; + } + } + } + decG(x); + return r; +} + +B slash_c2(B t, B w, B x) { + if (isArr(x) && RNK(x)==1 && isArr(w) && RNK(w)==1 && depth(w)==1) { + usz wia = IA(w); + usz xia = IA(x); + if (RARE(wia!=xia)) { + if (wia==0) { decG(w); return x; } + thrF("/: Lengths of components of ๐•จ must match ๐•ฉ (%s โ‰  %s)", wia, xia); + } + B xf = getFillQ(x); + + usz ri = 0; + if (TI(w,elType)==el_bit) { + B r; + u64* wp = bitarr_ptr(w); + u8 xe = TI(x,elType); + #ifdef __BMI2__ + usz wsum = bit_sum(wp, wia); + switch(xe) { + case el_bit: { + u64* xp = bitarr_ptr(x); u64* rp; r = m_bitarrv(&rp,wsum+128); a(r)->ia = wsum; + u64 cw = 0; // current word + u64 ro = 0; // offset in word where next bit should be written; never 64 + for (usz i=0; i=64) { + *(rp++) = cw; + cw = ro? v>>(64-ro) : 0; + } + ro = ro2&63; + } + if (ro) *rp = cw; + goto bit_ret; + } + #if SINGELI + case el_i8: case el_c8: { i8* xp=tyany_ptr(x); i8* rp=m_tyarrvO(&r,1,wsum,el2t(xe), 8); bmipopc_2slash8 (wp, xp, rp, wia); goto bit_ret; } + case el_i16:case el_c16: { i16* xp=tyany_ptr(x); i16* rp=m_tyarrvO(&r,2,wsum,el2t(xe), 16); bmipopc_2slash16(wp, xp, rp, wia); goto bit_ret; } + case el_i32:case el_c32: { + i32* xp=tyany_ptr(x); i32* rp=m_tyarrv(&r,4,wsum,el2t(xe)); + usz b = 1<<7; + TALLOC(i8, buf, b); + i32* rq=rp; i32* end=xp+xia-b; + while (xp < end) { + bmipopc_1slash8(wp, buf, b); + usz bs = bit_sum(wp, b); + for (usz j=0; j0 && !bitp_get(wp,wia-1)) wia--; + #ifndef __BMI2__ + usz wsum = bit_sum(wp, wia); + #endif + if (wsum==0) { decG(w); decG(x); return q_N(xf)? emptyHVec() : isF64(xf)? emptyIVec() : isC32(xf)? emptyCVec() : m_emptyFVec(xf); } + switch(xe) { default: UD; + + #if !SINGELI + case el_i8: case el_c8: { i8* xp=tyany_ptr(x); i8* rp=m_tyarrv(&r,1,wsum,el2t(xe)); for (usz i=0; i0 && !wp[wia-1]) wia--; \ + i64 wsum = 0; \ + u32 or = 0; \ + for (usz i = 0; i < wia; i++) { \ + wsum+= wp[i]; \ + or|= (u32)wp[i]; \ + } \ + if (or>>SIGN) thrM("/: ๐•จ must consist of natural numbers"); \ + if (TI(x,elType)==el_bit) { \ + u64* xp = bitarr_ptr(x); u64 ri=0; \ + u64* rp; B r = m_bitarrv(&rp, wsum); \ + if (or<2) for (usz i = 0; i < wia; i++) { \ + bitp_set(rp, ri, bitp_get(xp,i)); \ + ri+= wp[i]; \ + } else for (usz i = 0; i < wia; i++) { \ + WT cw = wp[i]; bool cx = bitp_get(xp,i); \ + for (i64 j = 0; j < cw; j++) bitp_set(rp, ri++, cx); \ + } \ + decG(w); decG(x); return r; \ + } \ + CASE(WT,i8) CASE(WT,i16) CASE(WT,i32) CASE(WT,f64) \ + SLOW2("๐•จ/๐•ฉ", w, x); \ + M_HARR(r, wsum) SGetU(x) \ + for (usz i = 0; i < wia; i++) { \ + i32 cw = wp[i]; if (cw==0) continue; \ + B cx = incBy(GetU(x, i), cw); \ + for (i64 j = 0; j < cw; j++) HARR_ADDA(r, cx);\ + } \ + decG(w); decG(x); \ + return withFill(HARR_FV(r), xf); \ + } + if (TI(w,elType)==el_i8 ) TYPED(i8,7); + if (TI(w,elType)==el_i32) TYPED(i32,31); + #undef TYPED + #undef CASE + SLOW2("๐•จ/๐•ฉ", w, x); + u64 ria = usum(w); + if (ria>=USZ_MAX) thrOOM(); + M_HARR(r, ria) SGetU(w) SGetU(x) + for (usz i = 0; i < wia; i++) { + usz c = o2s(GetU(w, i)); + if (c) { + B cx = incBy(GetU(x, i), c); + for (usz j = 0; RARE(j < c); j++) HARR_ADDA(r, cx); + } + } + decG(w); decG(x); + return withFill(HARR_FV(r), xf); + } + if (isArr(x) && RNK(x)==1 && q_i32(w)) { + usz xia = IA(x); + i32 wv = o2i(w); + if (wv<=0) { + if (wv<0) thrM("/: ๐•จ cannot be negative"); + return taga(arr_shVec(TI(x,slice)(x, 0, 0))); + } + if (TI(x,elType)==el_i32) { + i32* xp = i32any_ptr(x); + i32* rp; B r = m_i32arrv(&rp, xia*wv); + for (usz i = 0; i < xia; i++) { + for (i64 j = 0; j < wv; j++) *rp++ = xp[i]; + } + decG(x); + return r; + } else { + SLOW2("๐•จ/๐•ฉ", w, x); + B xf = getFillQ(x); + HArr_p r = m_harrUv(xia*wv); + SGetU(x) + for (usz i = 0; i < xia; i++) { + B cx = incBy(GetU(x, i), wv); + for (i64 j = 0; j < wv; j++) *r.a++ = cx; + } + decG(x); + return withFill(r.b, xf); + } + } + return c2(rt_slash, w, x); +} + + +B slash_im(B t, B x) { + if (!isArr(x) || RNK(x)!=1) thrM("/โผ: Argument must be an array"); + u8 xe = TI(x,elType); + usz xia = IA(x); + if (xia==0) { decG(x); return emptyIVec(); } + switch(xe) { default: UD; + case el_bit: { + usz sum = bit_sum(bitarr_ptr(x), xia); + usz ria = 1 + (sum>0); + f64* rp; B r = m_f64arrv(&rp, ria); + rp[sum>0] = sum; rp[0] = xia - sum; + decG(x); return num_squeeze(r); + } +#define CASE_SMALL(N) \ + case el_i##N: { \ + i##N* xp = i##N##any_ptr(x); \ + usz m=1<xp[a-1]) a++; \ + max=xp[a-1]; \ + if (a==xia) { /* Sorted unique argument */ \ + usz ria = max + 1; \ + u64* rp; r = m_bitarrv(&rp, ria); \ + for (usz i=0; imax) max=c; } \ + if ((i##N)max<0) thrM("/โผ: Argument cannot contain negative numbers"); \ + usz ria = max+1; \ + i##N* rp; r = m_i##N##arrv(&rp, ria); for (usz i=0; im/2) thrM("/โผ: Argument cannot contain negative numbers"); \ + i32* rp; r = m_i32arrv(&rp, ria); for (usz i=0; imax?c:max; if (c<0) thrM("/โผ: Argument cannot contain negative numbers"); } + usz ria = max+1; + if (i==xia) { + u64* rp; r = m_bitarrv(&rp, ria); for (usz i=0; imax?c:max; if (c<0) thrM("/โผ: Argument cannot contain negative numbers"); } + usz ria = max+1; if (ria==0) thrOOM(); + if (i==xia) { + u64* rp; r = m_bitarrv(&rp, ria); for (usz i=0; ia; } + usz i,j; B r; i64 max=-1; + for (i = 0; i < xia; i++) { i64 c=o2i64(xp[i]); if (c<=max) break; max=c; } + for (j = i; j < xia; j++) { i64 c=o2i64(xp[j]); max=c>max?c:max; if (c<0) thrM("/โผ: Argument cannot contain negative numbers"); } + if (max > USZ_MAX-1) thrOOM(); + usz ria = max+1; + if (i==xia) { + u64* rp; r = m_bitarrv(&rp, ria); for (usz i=0; i Date: Sat, 10 Sep 2022 13:13:41 -0400 Subject: [PATCH 08/19] Branchless i32 sparse Where --- src/builtins/slash.c | 43 +++++++++++++++++++++++++++++++++---------- 1 file changed, 33 insertions(+), 10 deletions(-) diff --git a/src/builtins/slash.c b/src/builtins/slash.c index 4518689a..4d9baa95 100644 --- a/src/builtins/slash.c +++ b/src/builtins/slash.c @@ -104,6 +104,11 @@ #endif #endif +#define WHERE_SPARSE(X,R,S) do { \ + for (usz i=0, j=0; j128 && s < xia/8+(xia<=32768?xia/4:0)) { + if (xia>128 && xia<=32768 && s < xia/16) { #else if (xia<=128 || s < xia/2+(xia<=32768?xia/4:0)) { #endif usz q=xia%64; if (q) xp[xia/64] &= ((u64)1<= b/8+b/16) { + bmipopc_1slash16(xp, buf, b); + for (usz j=0; j= b/256) { + for (usz j=0; j>=24; + p=(u32)u&(top-1); rq[j]+=(3*top)|p; j+=POPC(p); u>>=24; + p=(u32)u ; rq[j]+=(3*top)|p; j+=POPC(p); + } + u64 t=((u64)i<<21)-2*top; + for (usz j=0; j>24) + CTZ((u32)t); + t &= t-1; + } + } else { + for (usz ii=i/64, j=0; j Date: Sat, 10 Sep 2022 22:12:29 -0400 Subject: [PATCH 09/19] Use blocked method for i32 Where with or without Singeli --- src/builtins/slash.c | 67 +++++++++++++++++++++++++++----------------- 1 file changed, 42 insertions(+), 25 deletions(-) diff --git a/src/builtins/slash.c b/src/builtins/slash.c index 4d9baa95..a92f7878 100644 --- a/src/builtins/slash.c +++ b/src/builtins/slash.c @@ -131,34 +131,40 @@ B slash_c1(B t, B x) { u8 xe = TI(x,elType); if (xe==el_bit) { u64* xp = bitarr_ptr(x); - // Sparse method with CTZ - #if SINGELI && defined(__BMI2__) - if (xia>128 && xia<=32768 && s < xia/16) { - #else - if (xia<=128 || s < xia/2+(xia<=32768?xia/4:0)) { - #endif + if (xia > 32768) { usz q=xia%64; if (q) xp[xia/64] &= ((u64)1<xia-i) { + b = xia-i; + bs = s-(rq-rp); + } else { + bs = bit_sum(xp,b); + } + #if SINGELI && defined(__BMI2__) if (bs >= b/8+b/16) { bmipopc_1slash16(xp, buf, b); for (usz j=0; j= b/256) { + #else + if (bs >= b/2) { + for (usz ii=0; ii<(b+7)/8; ii++) { + u8 v = ((u8*)xp)[ii]; + i32* rs=rq; + for (usz k=0; k<8; k++) { *rs=i+8*ii+k; rs+=v&1; v>>=1; } + } + #endif + } else if (bs >= b/256) { // Branchless sparse for (usz j=0; j>=24; p=(u32)u&(top-1); rq[j]+=(3*top)|p; j+=POPC(p); u>>=24; @@ -170,7 +176,7 @@ B slash_c1(B t, B x) { rq[j] = 8*(t>>24) + CTZ((u32)t); t &= t-1; } - } else { + } else { // Branched very sparse for (usz ii=i/64, j=0; j128 && s < xia/16) { + #else + if (xia<=128 || s < xia/2) { + #endif + usz q=xia%64; if (q) xp[xia/64] &= ((u64)1<0 && !bitp_get(xp,xia-1)) xia--; + i32* rp = m_tyarrvO(&r, 4, s, t_i32arr, 4); u8* x8 = (u8*)xp; u8 q=xia%8; if (q) x8[xia/8] &= (1< Date: Sun, 11 Sep 2022 08:31:20 -0400 Subject: [PATCH 10/19] Separate where into its own function and split different types completely --- src/builtins/slash.c | 171 +++++++++++++++++++++++-------------------- 1 file changed, 90 insertions(+), 81 deletions(-) diff --git a/src/builtins/slash.c b/src/builtins/slash.c index a92f7878..6bb45fe4 100644 --- a/src/builtins/slash.c +++ b/src/builtins/slash.c @@ -109,6 +109,95 @@ for (u64 v=X[i]; v; v&=v-1) R[j++] = i*64 + CTZ(v); \ } while (0) +static B where(B x, usz xia, u64 s) { + B r; + u64* xp = bitarr_ptr(x); + usz q=xia%64; if (q) xp[xia/64] &= ((u64)1<= xia/16) { + i16* rp = m_tyarrvO(&r, 2, s, t_i16arr, 16); + bmipopc_1slash16(xp, rp, xia); + } + #else + if (s >= xia/2) { + i16* rp = m_tyarrvO(&r, 2, s, t_i16arr, 2); + for (usz i=0; i<(xia+7)/8; i++) { + u8 v = ((u8*)xp)[i]; + for (usz k=0; k<8; k++) { *rp=8*i+k; rp+=v&1; v>>=1; } + } + } + #endif + else { + i16* rp; r=m_i16arrv(&rp,s); WHERE_SPARSE(xp,rp,s); + } + } else { + assert(xia <= (usz)I32_MAX+1); + #if SINGELI && defined(__BMI2__) + i32* rp; r = m_i32arrv(&rp, s); + #else + i32* rp = m_tyarrvO(&r, 4, s, t_i32arr, 4); + #endif + usz b = 1<<11; // Maximum allowed for branchless sparse method + TALLOC(i16, buf, b); + i32* rq=rp; usz i=0; + for (; ixia-i) { + b = xia-i; + bs = s-(rq-rp); + } else { + bs = bit_sum(xp,b); + } + #if SINGELI && defined(__BMI2__) + if (bs >= b/8+b/16) { + bmipopc_1slash16(xp, buf, b); + for (usz j=0; j= b/2) { + for (usz ii=0; ii<(b+7)/8; ii++) { + u8 v = ((u8*)xp)[ii]; + i32* rs=rq; + for (usz k=0; k<8; k++) { *rs=i+8*ii+k; rs+=v&1; v>>=1; } + } + } + #endif + else if (bs >= b/256) { // Branchless sparse + for (usz j=0; j>=24; + p=(u32)u&(top-1); rq[j]+=(3*top)|p; j+=POPC(p); u>>=24; + p=(u32)u ; rq[j]+=(3*top)|p; j+=POPC(p); + } + u64 t=((u64)i<<21)-2*top; + for (usz j=0; j>24) + CTZ((u32)t); + t &= t-1; + } + } else { // Branched very sparse + for (usz ii=i/64, j=0; j 32768) { - usz q=xia%64; if (q) xp[xia/64] &= ((u64)1<xia-i) { - b = xia-i; - bs = s-(rq-rp); - } else { - bs = bit_sum(xp,b); - } - #if SINGELI && defined(__BMI2__) - if (bs >= b/8+b/16) { - bmipopc_1slash16(xp, buf, b); - for (usz j=0; j= b/2) { - for (usz ii=0; ii<(b+7)/8; ii++) { - u8 v = ((u8*)xp)[ii]; - i32* rs=rq; - for (usz k=0; k<8; k++) { *rs=i+8*ii+k; rs+=v&1; v>>=1; } - } - #endif - } else if (bs >= b/256) { // Branchless sparse - for (usz j=0; j>=24; - p=(u32)u&(top-1); rq[j]+=(3*top)|p; j+=POPC(p); u>>=24; - p=(u32)u ; rq[j]+=(3*top)|p; j+=POPC(p); - } - u64 t=((u64)i<<21)-2*top; - for (usz j=0; j>24) + CTZ((u32)t); - t &= t-1; - } - } else { // Branched very sparse - for (usz ii=i/64, j=0; j128 && s < xia/16) { - #else - if (xia<=128 || s < xia/2) { - #endif - usz q=xia%64; if (q) xp[xia/64] &= ((u64)1<>=1; } - } - } - #endif + r = where(x, xia, s); } else { i32* rp; r = m_i32arrv(&rp, s); if (xe==el_i8) { From 4bc5350747709deb0f185c40d25b110d9f836ae7 Mon Sep 17 00:00:00 2001 From: Marshall Lochbaum Date: Sun, 11 Sep 2022 11:00:17 -0400 Subject: [PATCH 11/19] Branchless sparse i16 Where --- src/builtins/slash.c | 96 ++++++++++++++++++++++++++++++-------------- 1 file changed, 66 insertions(+), 30 deletions(-) diff --git a/src/builtins/slash.c b/src/builtins/slash.c index 6bb45fe4..e5b527bd 100644 --- a/src/builtins/slash.c +++ b/src/builtins/slash.c @@ -104,11 +104,58 @@ #endif #endif -#define WHERE_SPARSE(X,R,S) do { \ - for (usz i=0, j=0; j>=24; + p=u&bsp_mask; buf[j]+=(3*bsp_top)|p; j+=POPC(p); u>>=24; + p=u ; buf[j]+=(3*bsp_top)|p; j+=POPC(p); + } + return j; +} +static void bsp_block_u32(u64* src, u32* dst, usz len, usz sum, usz off) { + for (usz j=0; j>24) + CTZ((u32)t); + t &= t-1; + } +} +static void bsp_u16(u64* src, u16* dst, usz len, usz sum) { + usz b = bsp_max; + usz bufsize = b len-i) b = len-i; + usz bs = bsp_fill(src+i/64, buf, b); + u64 t=((u64)i<<21)-2*bsp_top; + for (usz j=0; j>24) + CTZ((u32)t); + t &= t-1; + } + buf[bs] = 0; + dst+= bs; + } + TFREE(buf); +} + static B where(B x, usz xia, u64 s) { B r; u64* xp = bitarr_ptr(x); @@ -118,16 +165,16 @@ static B where(B x, usz xia, u64 s) { i8* rp = m_tyarrvO(&r, 1, s, t_i8arr, 8); bmipopc_1slash8(xp, rp, xia); #else - i8* rp; r=m_i8arrv(&rp,s); WHERE_SPARSE(xp,rp,s); + i8* rp; r=m_i8arrv(&rp,s); WHERE_SPARSE(xp,rp,s,0,); #endif } else if (xia <= 32768) { #if SINGELI && defined(__BMI2__) - if (s >= xia/16) { + if (s >= xia/8) { i16* rp = m_tyarrvO(&r, 2, s, t_i16arr, 16); bmipopc_1slash16(xp, rp, xia); } #else - if (s >= xia/2) { + if (s >= xia/4+xia/8) { i16* rp = m_tyarrvO(&r, 2, s, t_i16arr, 2); for (usz i=0; i<(xia+7)/8; i++) { u8 v = ((u8*)xp)[i]; @@ -136,7 +183,12 @@ static B where(B x, usz xia, u64 s) { } #endif else { - i16* rp; r=m_i16arrv(&rp,s); WHERE_SPARSE(xp,rp,s); + i16* rp; r=m_i16arrv(&rp,s); + if (s >= xia/128) { + bsp_u16(xp, (u16*)rp, xia, s); + } else { + WHERE_SPARSE(xp, rp, s, 0, RARE); + } } } else { assert(xia <= (usz)I32_MAX+1); @@ -145,10 +197,9 @@ static B where(B x, usz xia, u64 s) { #else i32* rp = m_tyarrvO(&r, 4, s, t_i32arr, 4); #endif - usz b = 1<<11; // Maximum allowed for branchless sparse method - TALLOC(i16, buf, b); - i32* rq=rp; usz i=0; - for (; ixia-i) { b = xia-i; @@ -170,25 +221,10 @@ static B where(B x, usz xia, u64 s) { } } #endif - else if (bs >= b/256) { // Branchless sparse - for (usz j=0; j>=24; - p=(u32)u&(top-1); rq[j]+=(3*top)|p; j+=POPC(p); u>>=24; - p=(u32)u ; rq[j]+=(3*top)|p; j+=POPC(p); - } - u64 t=((u64)i<<21)-2*top; - for (usz j=0; j>24) + CTZ((u32)t); - t &= t-1; - } - } else { // Branched very sparse - for (usz ii=i/64, j=0; j= b/256) { + bsp_block_u32(xp, (u32*)rq, b, bs, i); + } else { + WHERE_SPARSE(xp-i/64, rq, bs, i/64, RARE); } rq+= bs; xp+= b/64; From 0d6894ddbc210dab8cb8bf696351b85d14fba7d0 Mon Sep 17 00:00:00 2001 From: Marshall Lochbaum Date: Sun, 11 Sep 2022 15:15:56 -0400 Subject: [PATCH 12/19] Fix broken non-Singeli dense i32 Where --- src/builtins/slash.c | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/src/builtins/slash.c b/src/builtins/slash.c index e5b527bd..8f71d124 100644 --- a/src/builtins/slash.c +++ b/src/builtins/slash.c @@ -214,9 +214,9 @@ static B where(B x, usz xia, u64 s) { } #else if (bs >= b/2) { + i32* rs=rq; for (usz ii=0; ii<(b+7)/8; ii++) { u8 v = ((u8*)xp)[ii]; - i32* rs=rq; for (usz k=0; k<8; k++) { *rs=i+8*ii+k; rs+=v&1; v>>=1; } } } From fda3efe2176a727c5cfa0623d24bda7dbabbf233 Mon Sep 17 00:00:00 2001 From: Marshall Lochbaum Date: Sun, 11 Sep 2022 15:31:23 -0400 Subject: [PATCH 13/19] More macros for Where --- src/builtins/slash.c | 44 +++++++++++++++++++++----------------------- 1 file changed, 21 insertions(+), 23 deletions(-) diff --git a/src/builtins/slash.c b/src/builtins/slash.c index 8f71d124..09b1d246 100644 --- a/src/builtins/slash.c +++ b/src/builtins/slash.c @@ -104,9 +104,18 @@ #endif #endif +// Dense Where, still significantly worse than SIMD +// Assumes modifiable DST +#define WHERE_DENSE(SRC, DST, LEN, OFF) do { \ + for (usz ii=0; ii<(LEN+7)/8; ii++) { \ + u8 v = ((u8*)SRC)[ii]; \ + for (usz k=0; k<8; k++) { *DST=OFF+8*ii+k; DST+=v&1; v>>=1; } \ + } \ + } while (0) + // Sparse Where with branching #define WHERE_SPARSE(X,R,S,I0,COND) do { \ - for (usz ii=I0, j=0; j>24) + CTZ((u32)t); \ + t &= t-1; \ + } static void bsp_block_u32(u64* src, u32* dst, usz len, usz sum, usz off) { for (usz j=0; j>24) + CTZ((u32)t); - t &= t-1; - } + BSP_WRITE(dst, dst, sum, off,); } static void bsp_u16(u64* src, u16* dst, usz len, usz sum) { usz b = bsp_max; @@ -144,13 +155,7 @@ static void bsp_u16(u64* src, u16* dst, usz len, usz sum) { for (usz i=0; i len-i) b = len-i; usz bs = bsp_fill(src+i/64, buf, b); - u64 t=((u64)i<<21)-2*bsp_top; - for (usz j=0; j>24) + CTZ((u32)t); - t &= t-1; - } - buf[bs] = 0; + BSP_WRITE(buf, dst, bs, i, buf[j]=0;); buf[bs]=0; dst+= bs; } TFREE(buf); @@ -176,10 +181,7 @@ static B where(B x, usz xia, u64 s) { #else if (s >= xia/4+xia/8) { i16* rp = m_tyarrvO(&r, 2, s, t_i16arr, 2); - for (usz i=0; i<(xia+7)/8; i++) { - u8 v = ((u8*)xp)[i]; - for (usz k=0; k<8; k++) { *rp=8*i+k; rp+=v&1; v>>=1; } - } + WHERE_DENSE(xp, rp, xia, 0); } #endif else { @@ -214,11 +216,7 @@ static B where(B x, usz xia, u64 s) { } #else if (bs >= b/2) { - i32* rs=rq; - for (usz ii=0; ii<(b+7)/8; ii++) { - u8 v = ((u8*)xp)[ii]; - for (usz k=0; k<8; k++) { *rs=i+8*ii+k; rs+=v&1; v>>=1; } - } + i32* rs=rq; WHERE_DENSE(xp, rs, b, i); } #endif else if (bs >= b/256) { From 819899c74594052da45f20ce287eb10c3d51f03e Mon Sep 17 00:00:00 2001 From: Marshall Lochbaum Date: Sun, 11 Sep 2022 15:42:28 -0400 Subject: [PATCH 14/19] Implement f64 Where with u16 blocks --- src/builtins/slash.c | 51 +++++++++++++++++++++++++++++++++----------- 1 file changed, 38 insertions(+), 13 deletions(-) diff --git a/src/builtins/slash.c b/src/builtins/slash.c index 09b1d246..c8c658c1 100644 --- a/src/builtins/slash.c +++ b/src/builtins/slash.c @@ -161,6 +161,23 @@ static void bsp_u16(u64* src, u16* dst, usz len, usz sum) { TFREE(buf); } +static void where_block_u16(u64* src, u16* dst, usz len, usz sum) { + assert(len <= bsp_max); + #if SINGELI && defined(__BMI2__) + if (sum >= len/8) bmipopc_1slash16(src, (i16*)dst, len); + #else + if (sum >= len/4+len/8) WHERE_DENSE(src, dst, len, 0); + #endif + else if (sum >= len/128) { + u32* buf = (u32*)dst; assert(sum*2 <= len); + for (usz j=0; jxia-i) { b=xia-i; bs=s-(rp-rp0); } else { bs=bit_sum(xp,b); } + where_block_u16(xp, buf, b, bs); + for (usz j=0; j=USZ_MAX) thrOOM(); if (s==0) { decG(x); return emptyIVec(); } usz xia = IA(x); - if (RARE(xia>=I32_MAX)) { - usz xia = IA(x); - SGetU(x) - f64* rp; B r = m_f64arrv(&rp, s); usz ri = 0; - for (usz i = 0; i < xia; i++) { - usz c = o2s(GetU(x, i)); - for (usz j = 0; j < c; j++) rp[ri++] = i; - } - decG(x); - return r; - } B r; u8 xe = TI(x,elType); if (xe==el_bit) { r = where(x, xia, s); + } else if (RARE(xia>=I32_MAX)) { + SGetU(x) + f64* rp; r = m_f64arrv(&rp, s); usz ri = 0; + for (usz i = 0; i < xia; i++) { + usz c = o2s(GetU(x, i)); + for (usz j = 0; j < c; j++) rp[ri++] = i; + } } else { i32* rp; r = m_i32arrv(&rp, s); if (xe==el_i8) { From ab25b08f84ed5c1d9f2a326f58f06f6b1c64ed22 Mon Sep 17 00:00:00 2001 From: Marshall Lochbaum Date: Sun, 11 Sep 2022 15:53:26 -0400 Subject: [PATCH 15/19] Move compress (bool/list) to its own function --- src/builtins/slash.c | 177 ++++++++++++++++++++++--------------------- 1 file changed, 90 insertions(+), 87 deletions(-) diff --git a/src/builtins/slash.c b/src/builtins/slash.c index c8c658c1..c0137fd9 100644 --- a/src/builtins/slash.c +++ b/src/builtins/slash.c @@ -261,6 +261,95 @@ static B where(B x, usz xia, u64 s) { return r; } +static B compress(B w, B x, usz wia, B xf) { + usz xia = wia; + usz ri = 0; + B r; + u64* wp = bitarr_ptr(w); + u8 xe = TI(x,elType); + #ifdef __BMI2__ + usz wsum = bit_sum(wp, wia); + switch(xe) { + case el_bit: { + u64* xp = bitarr_ptr(x); u64* rp; r = m_bitarrv(&rp,wsum+128); a(r)->ia = wsum; + u64 cw = 0; // current word + u64 ro = 0; // offset in word where next bit should be written; never 64 + for (usz i=0; i=64) { + *(rp++) = cw; + cw = ro? v>>(64-ro) : 0; + } + ro = ro2&63; + } + if (ro) *rp = cw; + return r; + } + #if SINGELI + case el_i8: case el_c8: { i8* xp=tyany_ptr(x); i8* rp=m_tyarrvO(&r,1,wsum,el2t(xe), 8); bmipopc_2slash8 (wp, xp, rp, wia); return r; } + case el_i16:case el_c16: { i16* xp=tyany_ptr(x); i16* rp=m_tyarrvO(&r,2,wsum,el2t(xe), 16); bmipopc_2slash16(wp, xp, rp, wia); return r; } + case el_i32:case el_c32: { + i32* xp=tyany_ptr(x); i32* rp=m_tyarrv(&r,4,wsum,el2t(xe)); + usz b = 1<<7; + TALLOC(i8, buf, b); + i32* rq=rp; i32* end=xp+xia-b; + while (xp < end) { + bmipopc_1slash8(wp, buf, b); + usz bs = bit_sum(wp, b); + for (usz j=0; j0 && !bitp_get(wp,wia-1)) wia--; + #ifndef __BMI2__ + usz wsum = bit_sum(wp, wia); + #endif + if (wsum==0) { return q_N(xf)? emptyHVec() : isF64(xf)? emptyIVec() : isC32(xf)? emptyCVec() : m_emptyFVec(xf); } + switch(xe) { default: UD; + + #if !SINGELI + case el_i8: case el_c8: { i8* xp=tyany_ptr(x); i8* rp=m_tyarrv(&r,1,wsum,el2t(xe)); for (usz i=0; iia = wsum; - u64 cw = 0; // current word - u64 ro = 0; // offset in word where next bit should be written; never 64 - for (usz i=0; i=64) { - *(rp++) = cw; - cw = ro? v>>(64-ro) : 0; - } - ro = ro2&63; - } - if (ro) *rp = cw; - goto bit_ret; - } - #if SINGELI - case el_i8: case el_c8: { i8* xp=tyany_ptr(x); i8* rp=m_tyarrvO(&r,1,wsum,el2t(xe), 8); bmipopc_2slash8 (wp, xp, rp, wia); goto bit_ret; } - case el_i16:case el_c16: { i16* xp=tyany_ptr(x); i16* rp=m_tyarrvO(&r,2,wsum,el2t(xe), 16); bmipopc_2slash16(wp, xp, rp, wia); goto bit_ret; } - case el_i32:case el_c32: { - i32* xp=tyany_ptr(x); i32* rp=m_tyarrv(&r,4,wsum,el2t(xe)); - usz b = 1<<7; - TALLOC(i8, buf, b); - i32* rq=rp; i32* end=xp+xia-b; - while (xp < end) { - bmipopc_1slash8(wp, buf, b); - usz bs = bit_sum(wp, b); - for (usz j=0; j0 && !bitp_get(wp,wia-1)) wia--; - #ifndef __BMI2__ - usz wsum = bit_sum(wp, wia); - #endif - if (wsum==0) { decG(w); decG(x); return q_N(xf)? emptyHVec() : isF64(xf)? emptyIVec() : isC32(xf)? emptyCVec() : m_emptyFVec(xf); } - switch(xe) { default: UD; - - #if !SINGELI - case el_i8: case el_c8: { i8* xp=tyany_ptr(x); i8* rp=m_tyarrv(&r,1,wsum,el2t(xe)); for (usz i=0; i Date: Sun, 11 Sep 2022 16:28:06 -0400 Subject: [PATCH 16/19] Compress for 4-byte and 8-byte data using index blocks --- src/builtins/slash.c | 49 +++++++++++++++++++------------------------- 1 file changed, 21 insertions(+), 28 deletions(-) diff --git a/src/builtins/slash.c b/src/builtins/slash.c index c0137fd9..bba78dd2 100644 --- a/src/builtins/slash.c +++ b/src/builtins/slash.c @@ -262,8 +262,6 @@ static B where(B x, usz xia, u64 s) { } static B compress(B w, B x, usz wia, B xf) { - usz xia = wia; - usz ri = 0; B r; u64* wp = bitarr_ptr(w); u8 xe = TI(x,elType); @@ -292,24 +290,6 @@ static B compress(B w, B x, usz wia, B xf) { #if SINGELI case el_i8: case el_c8: { i8* xp=tyany_ptr(x); i8* rp=m_tyarrvO(&r,1,wsum,el2t(xe), 8); bmipopc_2slash8 (wp, xp, rp, wia); return r; } case el_i16:case el_c16: { i16* xp=tyany_ptr(x); i16* rp=m_tyarrvO(&r,2,wsum,el2t(xe), 16); bmipopc_2slash16(wp, xp, rp, wia); return r; } - case el_i32:case el_c32: { - i32* xp=tyany_ptr(x); i32* rp=m_tyarrv(&r,4,wsum,el2t(xe)); - usz b = 1<<7; - TALLOC(i8, buf, b); - i32* rq=rp; i32* end=xp+xia-b; - while (xp < end) { - bmipopc_1slash8(wp, buf, b); - usz bs = bit_sum(wp, b); - for (usz j=0; jwia-i) { b=wia-i; bs=wsum-(rp-rp0); } \ + else { bs=bit_sum(wp,b); } \ + where_block_u16(wp, (u16*)buf, b, bs); \ + for (usz j=0; j Date: Sun, 11 Sep 2022 16:28:06 -0400 Subject: [PATCH 17/19] Faster trailing zero trim, and apply it to all cases in compress() --- src/builtins/slash.c | 15 +++++++++------ 1 file changed, 9 insertions(+), 6 deletions(-) diff --git a/src/builtins/slash.c b/src/builtins/slash.c index bba78dd2..a0dc8ef2 100644 --- a/src/builtins/slash.c +++ b/src/builtins/slash.c @@ -264,9 +264,17 @@ static B where(B x, usz xia, u64 s) { static B compress(B w, B x, usz wia, B xf) { B r; u64* wp = bitarr_ptr(w); + u64 we = 0; + usz ie = wia/64; + usz q=wia%64; if (q) we = wp[ie] &= ((u64)1<ia = wsum; @@ -294,11 +302,6 @@ static B compress(B w, B x, usz wia, B xf) { } #endif - while (wia>0 && !bitp_get(wp,wia-1)) wia--; - #ifndef __BMI2__ - usz wsum = bit_sum(wp, wia); - #endif - if (wsum==0) { return q_N(xf)? emptyHVec() : isF64(xf)? emptyIVec() : isC32(xf)? emptyCVec() : m_emptyFVec(xf); } switch(xe) { default: UD; #if !SINGELI From 9b17701c6684929313e74a4a53dc6aa771881c83 Mon Sep 17 00:00:00 2001 From: Marshall Lochbaum Date: Sun, 11 Sep 2022 20:08:20 -0400 Subject: [PATCH 18/19] Sparse 1- and 2-byte compress, with blocking, if total sum is low enough --- src/builtins/slash.c | 46 +++++++++++++++++++++++--------------------- 1 file changed, 24 insertions(+), 22 deletions(-) diff --git a/src/builtins/slash.c b/src/builtins/slash.c index a0dc8ef2..1dd24b4f 100644 --- a/src/builtins/slash.c +++ b/src/builtins/slash.c @@ -262,7 +262,6 @@ static B where(B x, usz xia, u64 s) { } static B compress(B w, B x, usz wia, B xf) { - B r; u64* wp = bitarr_ptr(w); u64 we = 0; usz ie = wia/64; @@ -273,11 +272,14 @@ static B compress(B w, B x, usz wia, B xf) { } wia = 64*(ie+1) - CLZ(we); usz wsum = bit_sum(wp, wia); + + B r; u8 xe = TI(x,elType); - #ifdef __BMI2__ - switch(xe) { + switch(xe) { default: UD; case el_bit: { - u64* xp = bitarr_ptr(x); u64* rp; r = m_bitarrv(&rp,wsum+128); a(r)->ia = wsum; + u64* xp = bitarr_ptr(x); u64* rp; + #if SINGELI && defined(__BMI2__) + r = m_bitarrv(&rp,wsum+128); a(r)->ia = wsum; u64 cw = 0; // current word u64 ro = 0; // offset in word where next bit should be written; never 64 for (usz i=0; i=wia/CUTOFF) { DENSE; } \ + else { rp=m_tyarrv(&r,W/8,wsum,el2t(xe)); COMPRESS_BLOCK(i##W); } \ + break; } + #if SINGELI + case el_i8: case el_c8: WITH_SPARSE( 8, 32, rp=m_tyarrvO(&r,1,wsum,el2t(xe), 8); bmipopc_2slash8 (wp, xp, rp, wia)) + case el_i16:case el_c16: WITH_SPARSE(16, 16, rp=m_tyarrvO(&r,2,wsum,el2t(xe), 16); bmipopc_2slash16(wp, xp, rp, wia)) + #else + case el_i8: case el_c8: WITH_SPARSE( 8, 2, rp=m_tyarrv(&r,1,wsum,el2t(xe)); for (usz i=0; i Date: Sun, 11 Sep 2022 21:48:31 -0400 Subject: [PATCH 19/19] =?UTF-8?q?Just=20return=20=F0=9D=95=A9=20for=20an?= =?UTF-8?q?=20all-ones=20compress?= MIME-Version: 1.0 Content-Type: text/plain; charset=UTF-8 Content-Transfer-Encoding: 8bit --- src/builtins/slash.c | 2 ++ 1 file changed, 2 insertions(+) diff --git a/src/builtins/slash.c b/src/builtins/slash.c index 1dd24b4f..410740ec 100644 --- a/src/builtins/slash.c +++ b/src/builtins/slash.c @@ -270,8 +270,10 @@ static B compress(B w, B x, usz wia, B xf) { if (RARE(ie==0)) { return q_N(xf)? emptyHVec() : isF64(xf)? emptyIVec() : isC32(xf)? emptyCVec() : m_emptyFVec(xf); } we = wp[--ie]; } + usz wia0 = wia; wia = 64*(ie+1) - CLZ(we); usz wsum = bit_sum(wp, wia); + if (wsum == wia0) { dec(xf); return inc(x); } B r; u8 xe = TI(x,elType);