From 57d33990e36ac0df7f9ef61e1ceb287ca6fb0b38 Mon Sep 17 00:00:00 2001 From: Marshall Lochbaum Date: Wed, 14 Sep 2022 14:48:33 -0400 Subject: [PATCH 01/20] Make self-search special code only depend on cell size --- src/builtins/selfsearch.c | 52 +++++++++++++++++++-------------------- src/core/arrFns.h | 15 ++++++++++- src/core/tyarr.c | 10 +++++++- 3 files changed, 49 insertions(+), 28 deletions(-) diff --git a/src/builtins/selfsearch.c b/src/builtins/selfsearch.c index 5eb430ca..7dcb30cc 100644 --- a/src/builtins/selfsearch.c +++ b/src/builtins/selfsearch.c @@ -7,12 +7,11 @@ B memberOf_c1(B t, B x) { if (isAtm(x) || RNK(x)==0) thrM("∊: Argument cannot have rank 0"); usz n = *SH(x); if (n==0) { decG(x); return emptyIVec(); } - if (RNK(x)>1) x = toCells(x); - u8 xe = TI(x,elType); - if (elChr(xe)) xe -= el_c8-el_i8; + u8 lw = cellWidthLog(x); + void* xv = tyany_ptr(x); #define BRUTE(T) \ - i##T* xp = tyany_ptr(x); \ + i##T* xp = xv; \ u64* rp; B r = m_bitarrv(&rp, n); bitp_set(rp, 0, 1); \ for (usz i=1; i8 && n3*n?tn:3*n)+(2*rx+1)*sizeof(usz))); @@ -81,6 +80,7 @@ B memberOf_c1(B t, B x) { } #undef BRUTE + if (RNK(x)>1) x = toCells(x); u64* rp; B r = m_bitarrv(&rp, n); H_Sb* set = m_Sb(64); SGetU(x) @@ -94,12 +94,11 @@ B count_c1(B t, B x) { usz n = *SH(x); if (n==0) { decG(x); return emptyIVec(); } if (n>(usz)I32_MAX+1) thrM("⊒: Argument length >2⋆31 not supported"); - if (RNK(x)>1) x = toCells(x); - u8 xe = TI(x,elType); - if (elChr(xe)) xe -= el_c8-el_i8; + u8 lw = cellWidthLog(x); + void* xv = tyany_ptr(x); #define BRUTE(T) \ - i##T* xp = tyany_ptr(x); \ + i##T* xp = xv; \ i8* rp; B r = m_i8arrv(&rp, n); rp[0]=0; \ for (usz i=1; i8 && nn?tn:n)+(2*rx+1)*sizeof(usz))); @@ -169,6 +168,7 @@ B count_c1(B t, B x) { } #undef BRUTE + if (RNK(x)>1) x = toCells(x); i32* rp; B r = m_i32arrv(&rp, n); H_b2i* map = m_b2i(64); SGetU(x) @@ -186,12 +186,11 @@ B indexOf_c1(B t, B x) { usz n = *SH(x); if (n==0) { decG(x); return emptyIVec(); } if (n>(usz)I32_MAX+1) thrM("⊐: Argument length >2⋆31 not supported"); - if (RNK(x)>1) x = toCells(x); - u8 xe = TI(x,elType); - if (elChr(xe)) xe -= el_c8-el_i8; + u8 lw = cellWidthLog(x); + void* xv = tyany_ptr(x); #define BRUTE(T) \ - i##T* xp = tyany_ptr(x); \ + i##T* xp = xv; \ i8* rp; B r = m_i8arrv(&rp, n); rp[0]=0; \ TALLOC(i##T, uniq, n); uniq[0]=xp[0]; \ for (usz i=1, u=1; i8 && n1) x = toCells(x); i32* rp; B r = m_i32arrv(&rp, n); H_b2i* map = m_b2i(64); SGetU(x) diff --git a/src/core/arrFns.h b/src/core/arrFns.h index 10c6e109..33f219a6 100644 --- a/src/core/arrFns.h +++ b/src/core/arrFns.h @@ -54,4 +54,17 @@ extern u8 elType2type[]; extern u8 elTypeWidth[]; #define elWidth(X) elTypeWidth[X] extern u8 arrTypeWidthLog[]; -#define arrTypeWidthLog(X) arrTypeWidthLog[X] \ No newline at end of file +#define arrTypeWidthLog(X) arrTypeWidthLog[X] +extern u8 arrTypeBitsLog[]; +#define arrTypeBitsLog(X) arrTypeBitsLog[X] + +// Log of width in bits: max of 7, and also return 7 if not power of 2 +static u8 cellWidthLog(B x) { + assert(isArr(x) && RNK(x)>=1); + u8 lw = arrTypeBitsLog(TY(x)); + if (LIKELY(RNK(x)==1)) return lw; + usz csz = arr_csz(x); + if (csz & (csz-1)) return 7; // Not power of 2 + lw += CTZ(csz); + return lw<7? lw : 7; +} diff --git a/src/core/tyarr.c b/src/core/tyarr.c index ec09a771..dfb11f2c 100644 --- a/src/core/tyarr.c +++ b/src/core/tyarr.c @@ -20,6 +20,14 @@ u8 arrTypeWidthLog[] = { [t_i32arr]=2, [t_i32slice]=2, [t_c32arr]=2, [t_c32slice]=2, [t_f64arr]=3, [t_f64slice]=3 }; +u8 arrTypeBitsLog[] = { + [t_bitarr]=0, + [t_i8arr ]=3, [t_i8slice ]=3, [t_c8arr ]=3, [t_c8slice ]=3, + [t_i16arr]=4, [t_i16slice]=4, [t_c16arr]=4, [t_c16slice]=4, + [t_i32arr]=5, [t_i32slice]=5, [t_c32arr]=5, [t_c32slice]=5, + [t_f64arr]=6, [t_f64slice]=6, + [t_harr ]=6, [t_hslice ]=6, [t_fillarr]=6,[t_fillslice]=6 +}; B m_i8(i8 x) { return m_i32(x); } B m_i16(i16 x) { return m_i32(x); } B m_c8(u8 x) { return m_c32(x); } B m_c16(u16 x) { return m_c32(x); } @@ -87,4 +95,4 @@ void tyarr_init() { Arr* emptySVec = arr_shVec(m_fillarrp(0)); fillarr_setFill(emptySVec, emptyCVec()); bi_emptySVec = taga(emptySVec); gc_add(bi_emptySVec); -} \ No newline at end of file +} From b84a2d4e1f2aa5ecd153e60714e2f0896efbd7b5 Mon Sep 17 00:00:00 2001 From: Marshall Lochbaum Date: Wed, 14 Sep 2022 15:42:21 -0400 Subject: [PATCH 02/20] Special-casing for length 1 in self-search functions --- src/builtins/selfsearch.c | 6 +++--- 1 file changed, 3 insertions(+), 3 deletions(-) diff --git a/src/builtins/selfsearch.c b/src/builtins/selfsearch.c index 7dcb30cc..37e647fd 100644 --- a/src/builtins/selfsearch.c +++ b/src/builtins/selfsearch.c @@ -6,7 +6,7 @@ B memberOf_c1(B t, B x) { if (isAtm(x) || RNK(x)==0) thrM("∊: Argument cannot have rank 0"); usz n = *SH(x); - if (n==0) { decG(x); return emptyIVec(); } + if (n<=1) { decG(x); return n ? taga(arr_shVec(allOnes(1))) : emptyIVec(); } u8 lw = cellWidthLog(x); void* xv = tyany_ptr(x); @@ -92,7 +92,7 @@ B memberOf_c1(B t, B x) { B count_c1(B t, B x) { if (isAtm(x) || RNK(x)==0) thrM("⊒: Argument cannot have rank 0"); usz n = *SH(x); - if (n==0) { decG(x); return emptyIVec(); } + if (n<=1) { decG(x); return n ? taga(arr_shVec(allZeroes(1))) : emptyIVec(); } if (n>(usz)I32_MAX+1) thrM("⊒: Argument length >2⋆31 not supported"); u8 lw = cellWidthLog(x); @@ -184,7 +184,7 @@ extern B rt_indexOf; B indexOf_c1(B t, B x) { if (isAtm(x) || RNK(x)==0) thrM("⊐: 𝕩 cannot have rank 0"); usz n = *SH(x); - if (n==0) { decG(x); return emptyIVec(); } + if (n<=1) { decG(x); return n ? taga(arr_shVec(allZeroes(1))) : emptyIVec(); } if (n>(usz)I32_MAX+1) thrM("⊐: Argument length >2⋆31 not supported"); u8 lw = cellWidthLog(x); From 0b1ded14e0c0e48570dbcf4078fad3c3d5f7f8e4 Mon Sep 17 00:00:00 2001 From: Marshall Lochbaum Date: Wed, 14 Sep 2022 17:45:48 -0400 Subject: [PATCH 03/20] Type-independent Reverse --- src/builtins/sfns.c | 32 +++++++++++++++++--------------- src/core/arrFns.h | 1 + 2 files changed, 18 insertions(+), 15 deletions(-) diff --git a/src/builtins/sfns.c b/src/builtins/sfns.c index ffb26c30..e4de367b 100644 --- a/src/builtins/sfns.c +++ b/src/builtins/sfns.c @@ -1030,22 +1030,24 @@ B group_c1(B t, B x) { extern B rt_reverse; B reverse_c1(B t, B x) { if (isAtm(x) || RNK(x)==0) thrM("⌽: Argument cannot be a unit"); - usz xia = IA(x); - if (xia==0) return x; - u8 xe = TI(x,elType); - if (RNK(x)==1) { + usz n = *SH(x); + if (n==0) return x; + u8 xl = cellWidthLog(x); + u8 xt = arrNewType(TY(x)); + if (xl <= 6 && (xl >=3 || xl == 0)) { + void* xv = tyany_ptr(x); B r; - switch(xe) { default: UD; - case el_bit: { u64* xp = bitarr_ptr(x); u64* rp; r = m_bitarrv(&rp, xia); for (usz i = 0; i < xia; i++) bitp_set(rp, i, bitp_get(xp, xia-i-1)); break; } - case el_i8: case el_c8: { u8* xp = tyany_ptr(x); u8* rp = m_tyarrv(&r, 1, xia, el2t(xe)); for (usz i = 0; i < xia; i++) rp[i] = xp[xia-i-1]; break; } - case el_i16:case el_c16: { u16* xp = tyany_ptr(x); u16* rp = m_tyarrv(&r, 2, xia, el2t(xe)); for (usz i = 0; i < xia; i++) rp[i] = xp[xia-i-1]; break; } - case el_i32:case el_c32: { u32* xp = tyany_ptr(x); u32* rp = m_tyarrv(&r, 4, xia, el2t(xe)); for (usz i = 0; i < xia; i++) rp[i] = xp[xia-i-1]; break; } - case el_f64: { f64* xp = f64any_ptr(x); f64* rp; r = m_f64arrv(&rp, xia); for (usz i = 0; i < xia; i++) rp[i] = xp[xia-i-1]; break; } - case el_B: { + switch(xl) { default: UD; break; + case 0: { u64* xp = xv; u64* rp; r = m_bitarrc(&rp, x); for (usz i = 0; i < n; i++) bitp_set(rp, i, bitp_get(xp, n-i-1)); break; } + case 3: { u8* xp = xv; u8* rp = m_tyarrc(&r, 1, x, xt); for (usz i = 0; i < n; i++) rp[i] = xp[n-i-1]; break; } + case 4: { u16* xp = xv; u16* rp = m_tyarrc(&r, 2, x, xt); for (usz i = 0; i < n; i++) rp[i] = xp[n-i-1]; break; } + case 5: { u32* xp = xv; u32* rp = m_tyarrc(&r, 4, x, xt); for (usz i = 0; i < n; i++) rp[i] = xp[n-i-1]; break; } + case 6: if (TI(x,elType)!=el_B) { u64* xp = xv; u64* rp = m_tyarrc(&r, 8, x, xt); for (usz i = 0; i < n; i++) rp[i] = xp[n-i-1]; break; } + else { HArr_p rp = m_harrUc(x); B* xp = arr_bptr(x); - if (xp!=NULL) for (usz i = 0; i < xia; i++) rp.a[i] = inc(xp[xia-i-1]); - else { SGet(x) for (usz i = 0; i < xia; i++) rp.a[i] = Get(x, xia-i-1); } + if (xp!=NULL) for (usz i = 0; i < n; i++) rp.a[i] = inc(xp[n-i-1]); + else { SGet(x) for (usz i = 0; i < n; i++) rp.a[i] = Get(x, n-i-1); } r = rp.b; B xf = getFillQ(x); decG(x); @@ -1060,8 +1062,8 @@ B reverse_c1(B t, B x) { usz csz = arr_csz(x); usz cam = SH(x)[0]; usz rp = 0; - usz ip = xia; - MAKE_MUT(r, xia); mut_init(r, xe); + usz ip = IA(x); + MAKE_MUT(r, ip); mut_init(r, TI(x,elType)); MUTG_INIT(r); for (usz i = 0; i < cam; i++) { ip-= csz; diff --git a/src/core/arrFns.h b/src/core/arrFns.h index 33f219a6..9ee17682 100644 --- a/src/core/arrFns.h +++ b/src/core/arrFns.h @@ -57,6 +57,7 @@ extern u8 arrTypeWidthLog[]; #define arrTypeWidthLog(X) arrTypeWidthLog[X] extern u8 arrTypeBitsLog[]; #define arrTypeBitsLog(X) arrTypeBitsLog[X] +#define arrNewType(X) el2t(ti_elType[X]) // Log of width in bits: max of 7, and also return 7 if not power of 2 static u8 cellWidthLog(B x) { From 21033fa3551b7b81638d747518ef8628d95094ab Mon Sep 17 00:00:00 2001 From: Marshall Lochbaum Date: Fri, 16 Sep 2022 10:59:03 -0400 Subject: [PATCH 04/20] Move Group to its own file --- makefile | 2 +- src/builtins/group.c | 225 ++++++++++++++++++++++++++++++++++++ src/builtins/sfns.c | 265 ------------------------------------------- src/builtins/slash.c | 51 ++++++++- src/opt/single.c | 1 + 5 files changed, 277 insertions(+), 267 deletions(-) create mode 100644 src/builtins/group.c diff --git a/makefile b/makefile index 6da1bdf2..38b18583 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 slash.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 group.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/group.c b/src/builtins/group.c new file mode 100644 index 00000000..7cca8b1a --- /dev/null +++ b/src/builtins/group.c @@ -0,0 +1,225 @@ +#include "../core.h" +#include "../utils/talloc.h" +#include "../utils/mut.h" +#include "../builtins.h" + +extern B ne_c2(B, B, B); +extern B slash_c1(B, B); +extern B slash_c2(B, B, B); +extern B select_c2(B, B, B); +extern B take_c2(B, B, B); +extern B drop_c2(B, B, B); +extern B join_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) { + usz wia = IA(w); + usz xia = IA(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_bit) w = taga(cpyI8Arr(w)); + i64 ria = 0; + bool bad = false, sort = true; + usz neg = 0, change = 0; + void *wp0 = tyany_ptr(w); + #define CASE(T) case el_##T: { \ + 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; \ + change += prev != n; \ + prev = n; \ + } \ + 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)) 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)); + 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); + 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 + + 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: { + void* xp = tyany_ptr(x); + 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; + case 4: for (usz i = 0; i < xia; i++) { i32 n = wp[i]; if (n>=0) ((u32*)tyarr_ptr(rp[n]))[pos[n]++] = ((u32*)xp)[i]; } break; + case 8: for (usz i = 0; i < xia; i++) { i32 n = wp[i]; if (n>=0) ((f64*)tyarr_ptr(rp[n]))[pos[n]++] = ((f64*)xp)[i]; } break; + } + break; + } + case el_bit: case el_B: { + for (usz i = 0; i < ria; i++) { + Arr* c = m_fillarrp(len[i]); + c->ia = 0; + fillarr_setFill(c, inc(xf)); + arr_shVec(c); + rp[i] = taga(c); + } + SLOW2("𝕨⊔𝕩", w, x); + SGet(x) + for (usz i = 0; i < xia; i++) { + i32 n = wp[i]; + if (n>=0) fillarr_ptr(a(rp[n]))[pos[n]++] = Get(x, i); + } + for (usz i = 0; i < ria; i++) a(rp[i])->ia = len[i]; + break; + } + } + intvec_ret: + fillarr_setFill(rf, xf); + decG(w); decG(x); TFREE(pos); + return taga(r); + } else { + SLOW2("𝕨⊔𝕩", w, x); + SGetU(w) + i64 ria = wia==xia? 0 : o2i64(GetU(w, xia)); + if (ria<-1) thrM("⊔: 𝕨 can't contain elements less than ¯1"); + ria--; + for (usz i = 0; i < xia; i++) { + B cw = GetU(w, i); + if (!q_i64(cw)) goto base; + i64 c = o2i64G(cw); + if (c>ria) ria = c; + if (c<-1) thrM("⊔: 𝕨 can't contain elements less than ¯1"); + } + if (ria > (i64)(USZ_MAX-1)) thrOOM(); + 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++) len[o2i64G(GetU(w, i))]++; + + 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); + + for (usz i = 0; i < ria; i++) { + Arr* c = m_fillarrp(len[i]); + c->ia = 0; + fillarr_setFill(c, inc(xf)); + arr_shVec(c); + rp[i] = taga(c); + } + Arr* rf = m_fillarrp(0); arr_shVec(rf); + fillarr_setFill(rf, xf); + fillarr_setFill(r, taga(rf)); + SGet(x) + for (usz i = 0; i < xia; i++) { + i64 n = o2i64G(GetU(w, i)); + if (n>=0) fillarr_ptr(a(rp[n]))[pos[n]++] = Get(x, i); + } + for (usz i = 0; i < ria; i++) a(rp[i])->ia = len[i]; + decG(w); decG(x); TFREE(lenO); TFREE(pos); + return taga(r); + } + } + base: + return c2(rt_group, w, x); +} +B ud_c1(B, B); +B group_c1(B t, B x) { + if (isArr(x) && RNK(x)==1 && TI(x,arrD1)) { + usz ia = IA(x); + B range = ud_c1(t, m_f64(ia)); + return group_c2(m_f64(0), x, range); + } + return c1(rt_group, x); +} + diff --git a/src/builtins/sfns.c b/src/builtins/sfns.c index e4de367b..d635da65 100644 --- a/src/builtins/sfns.c +++ b/src/builtins/sfns.c @@ -811,222 +811,6 @@ 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 slash_c1(B, B); -extern B slash_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) { - usz wia = IA(w); - usz xia = IA(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_bit) w = taga(cpyI8Arr(w)); - i64 ria = 0; - bool bad = false, sort = true; - usz neg = 0, change = 0; - void *wp0 = tyany_ptr(w); - #define CASE(T) case el_##T: { \ - 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; \ - change += prev != n; \ - prev = n; \ - } \ - 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)) 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)); - 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); - 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 - - 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: { - void* xp = tyany_ptr(x); - 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; - case 4: for (usz i = 0; i < xia; i++) { i32 n = wp[i]; if (n>=0) ((u32*)tyarr_ptr(rp[n]))[pos[n]++] = ((u32*)xp)[i]; } break; - case 8: for (usz i = 0; i < xia; i++) { i32 n = wp[i]; if (n>=0) ((f64*)tyarr_ptr(rp[n]))[pos[n]++] = ((f64*)xp)[i]; } break; - } - break; - } - case el_bit: case el_B: { - for (usz i = 0; i < ria; i++) { - Arr* c = m_fillarrp(len[i]); - c->ia = 0; - fillarr_setFill(c, inc(xf)); - arr_shVec(c); - rp[i] = taga(c); - } - SLOW2("𝕨⊔𝕩", w, x); - SGet(x) - for (usz i = 0; i < xia; i++) { - i32 n = wp[i]; - if (n>=0) fillarr_ptr(a(rp[n]))[pos[n]++] = Get(x, i); - } - for (usz i = 0; i < ria; i++) a(rp[i])->ia = len[i]; - break; - } - } - intvec_ret: - fillarr_setFill(rf, xf); - decG(w); decG(x); TFREE(pos); - return taga(r); - } else { - SLOW2("𝕨⊔𝕩", w, x); - SGetU(w) - i64 ria = wia==xia? 0 : o2i64(GetU(w, xia)); - if (ria<-1) thrM("⊔: 𝕨 can't contain elements less than ¯1"); - ria--; - for (usz i = 0; i < xia; i++) { - B cw = GetU(w, i); - if (!q_i64(cw)) goto base; - i64 c = o2i64G(cw); - if (c>ria) ria = c; - if (c<-1) thrM("⊔: 𝕨 can't contain elements less than ¯1"); - } - if (ria > (i64)(USZ_MAX-1)) thrOOM(); - 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++) len[o2i64G(GetU(w, i))]++; - - 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); - - for (usz i = 0; i < ria; i++) { - Arr* c = m_fillarrp(len[i]); - c->ia = 0; - fillarr_setFill(c, inc(xf)); - arr_shVec(c); - rp[i] = taga(c); - } - Arr* rf = m_fillarrp(0); arr_shVec(rf); - fillarr_setFill(rf, xf); - fillarr_setFill(r, taga(rf)); - SGet(x) - for (usz i = 0; i < xia; i++) { - i64 n = o2i64G(GetU(w, i)); - if (n>=0) fillarr_ptr(a(rp[n]))[pos[n]++] = Get(x, i); - } - for (usz i = 0; i < ria; i++) a(rp[i])->ia = len[i]; - decG(w); decG(x); TFREE(lenO); TFREE(pos); - return taga(r); - } - } - base: - return c2(rt_group, w, x); -} -B ud_c1(B, B); -B group_c1(B t, B x) { - if (isArr(x) && RNK(x)==1 && TI(x,arrD1)) { - usz ia = IA(x); - B range = ud_c1(t, m_f64(ia)); - return group_c2(m_f64(0), x, range); - } - return c1(rt_group, x); -} - extern B rt_reverse; B reverse_c1(B t, B x) { if (isAtm(x) || RNK(x)==0) thrM("⌽: Argument cannot be a unit"); @@ -1194,54 +978,6 @@ B pick_ucw(B t, B o, B w, B x) { return qWithFill(mut_fcd(r, x), xf); } -B slash_ucw(B t, B o, B w, B x) { - if (isAtm(w) || isAtm(x) || RNK(w)!=1 || RNK(x)!=1 || IA(w)!=IA(x)) return def_fn_ucw(t, o, w, x); - usz ia = IA(x); - SGetU(w) - if (!elInt(TI(w,elType))) for (usz i = 0; i < ia; i++) if (!q_i32(GetU(w,i))) return def_fn_ucw(t, o, w, x); - B arg = slash_c2(t, inc(w), inc(x)); - usz argIA = IA(arg); - B rep = c1(o, arg); - if (isAtm(rep) || RNK(rep)!=1 || IA(rep) != argIA) thrF("𝔽⌾(a⊸/)𝕩: Result of 𝔽 must have the same shape as a/𝕩 (expected ⟨%s⟩, got %H)", argIA, rep); - MAKE_MUT(r, ia); mut_init(r, el_or(TI(x,elType), TI(rep,elType))); - SGet(x) - SGet(rep) - usz repI = 0; - if (TY(w) == t_bitarr) { - u64* d = bitarr_ptr(w); - if (elInt(TI(x,elType)) && elInt(TI(rep,elType))) { - if (r->fns->elType!=el_i32) mut_to(r, el_i32); - i32* rp = r->ai32; - x = toI32Any(x); i32* xp = i32any_ptr(x); - rep = toI32Any(rep); i32* np = i32any_ptr(rep); - for (usz i = 0; i < ia; i++) { - bool v = bitp_get(d, i); - i32 nc = np[repI]; - i32 xc = xp[i]; - rp[i] = v? nc : xc; - repI+= v; - } - } else { - MUTG_INIT(r); - for (usz i = 0; i < ia; i++) mut_setG(r, i, bitp_get(d, i)? Get(rep,repI++) : Get(x,i)); - } - } else { - SGetU(rep) - MUTG_INIT(r); - for (usz i = 0; i < ia; i++) { - i32 cw = o2iG(GetU(w, i)); - if (cw) { - B cr = Get(rep,repI); - if (CHECK_VALID) for (i32 j = 1; j < cw; j++) if (!equal(GetU(rep,repI+j), cr)) { mut_pfree(r,i); thrM("𝔽⌾(a⊸/): Incompatible result elements"); } - mut_setG(r, i, cr); - repI+= cw; - } else mut_setG(r, i, Get(x,i)); - } - } - decG(w); decG(rep); - return mut_fcd(r, x); -} - static B takedrop_ucw(i64 wi, B o, u64 am, B x, size_t xr) { usz xia = IA(x); usz csz = arr_csz(x); @@ -1327,7 +1063,6 @@ void sfns_init() { c(BFn,bi_pick)->uc1 = pick_uc1; c(BFn,bi_reverse)->uc1 = reverse_uc1; c(BFn,bi_pick)->ucw = pick_ucw; - c(BFn,bi_slash)->ucw = slash_ucw; c(BFn,bi_select)->ucw = select_ucw; // TODO move to new init fn c(BFn,bi_shape)->uc1 = shape_uc1; c(BFn,bi_transp)->uc1 = transp_uc1; diff --git a/src/builtins/slash.c b/src/builtins/slash.c index 437f9480..208aa15d 100644 --- a/src/builtins/slash.c +++ b/src/builtins/slash.c @@ -617,6 +617,55 @@ B slash_im(B t, B x) { } } +B slash_ucw(B t, B o, B w, B x) { + if (isAtm(w) || isAtm(x) || RNK(w)!=1 || RNK(x)!=1 || IA(w)!=IA(x)) return def_fn_ucw(t, o, w, x); + usz ia = IA(x); + SGetU(w) + if (!elInt(TI(w,elType))) for (usz i = 0; i < ia; i++) if (!q_i32(GetU(w,i))) return def_fn_ucw(t, o, w, x); + B arg = slash_c2(t, inc(w), inc(x)); + usz argIA = IA(arg); + B rep = c1(o, arg); + if (isAtm(rep) || RNK(rep)!=1 || IA(rep) != argIA) thrF("𝔽⌾(a⊸/)𝕩: Result of 𝔽 must have the same shape as a/𝕩 (expected ⟨%s⟩, got %H)", argIA, rep); + MAKE_MUT(r, ia); mut_init(r, el_or(TI(x,elType), TI(rep,elType))); + SGet(x) + SGet(rep) + usz repI = 0; + if (TY(w) == t_bitarr) { + u64* d = bitarr_ptr(w); + if (elInt(TI(x,elType)) && elInt(TI(rep,elType))) { + if (r->fns->elType!=el_i32) mut_to(r, el_i32); + i32* rp = r->ai32; + x = toI32Any(x); i32* xp = i32any_ptr(x); + rep = toI32Any(rep); i32* np = i32any_ptr(rep); + for (usz i = 0; i < ia; i++) { + bool v = bitp_get(d, i); + i32 nc = np[repI]; + i32 xc = xp[i]; + rp[i] = v? nc : xc; + repI+= v; + } + } else { + MUTG_INIT(r); + for (usz i = 0; i < ia; i++) mut_setG(r, i, bitp_get(d, i)? Get(rep,repI++) : Get(x,i)); + } + } else { + SGetU(rep) + MUTG_INIT(r); + for (usz i = 0; i < ia; i++) { + i32 cw = o2iG(GetU(w, i)); + if (cw) { + B cr = Get(rep,repI); + if (CHECK_VALID) for (i32 j = 1; j < cw; j++) if (!equal(GetU(rep,repI+j), cr)) { mut_pfree(r,i); thrM("𝔽⌾(a⊸/): Incompatible result elements"); } + mut_setG(r, i, cr); + repI+= cw; + } else mut_setG(r, i, Get(x,i)); + } + } + decG(w); decG(rep); + return mut_fcd(r, x); +} + void slash_init() { c(BFn,bi_slash)->im = slash_im; -} \ No newline at end of file + c(BFn,bi_slash)->ucw = slash_ucw; +} diff --git a/src/opt/single.c b/src/opt/single.c index c225896d..f17124e3 100644 --- a/src/opt/single.c +++ b/src/opt/single.c @@ -17,6 +17,7 @@ #include "../builtins/sfns.c" #include "../builtins/select.c" #include "../builtins/slash.c" +#include "../builtins/group.c" #include "../builtins/sysfn.c" #include "../builtins/sort.c" #include "../builtins/selfsearch.c" From 3505e1515e8005e3c2794d00c396ea2c063d3967 Mon Sep 17 00:00:00 2001 From: Marshall Lochbaum Date: Fri, 16 Sep 2022 11:42:55 -0400 Subject: [PATCH 05/20] =?UTF-8?q?Extend=20Group=20where=20=F0=9D=95=A8=20i?= =?UTF-8?q?s=20an=20integer=20list=20to=20handle=20higher-rank=20?= =?UTF-8?q?=F0=9D=95=A9?= MIME-Version: 1.0 Content-Type: text/plain; charset=UTF-8 Content-Transfer-Encoding: 8bit --- src/builtins/group.c | 162 +++++++++++++++++++++++++------------------ 1 file changed, 93 insertions(+), 69 deletions(-) diff --git a/src/builtins/group.c b/src/builtins/group.c index 7cca8b1a..aa84021c 100644 --- a/src/builtins/group.c +++ b/src/builtins/group.c @@ -3,6 +3,7 @@ #include "../utils/mut.h" #include "../builtins.h" +extern B ud_c1(B, B); extern B ne_c2(B, B, B); extern B slash_c1(B, B); extern B slash_c2(B, B, B); @@ -11,11 +12,29 @@ extern B take_c2(B, B, B); extern B drop_c2(B, B, B); extern B join_c2(B, B, B); +static Arr* arr_shChangeLen(Arr* a, ur r, usz* xsh, usz len) { + assert(r > 1); + usz* sh = a->sh = m_shArr(r)->a; + SPRNK(a,r); + sh[0] = len; + shcpy(sh+1, xsh+1, r-1); + return a; +} +static B m_shChangeLen(u8 xt, ur xr, usz* xsh, usz l, usz cw, usz csz) { + return taga(arr_shChangeLen(m_arr(offsetof(TyArr, a)+l*cw, xt, l*csz), xr, xsh, l)); +} +static void allocGroups(B* rp, usz ria, B z, u8 xt, ur xr, usz* xsh, i32* len, usz width, usz csz) { + if (xr==1) for (usz j = 0; j < ria; j++) { usz l=len[j]; if (!l) rp[j] = inc(z); else m_tyarrv(rp+j, width, l, xt); } + else for (usz j = 0; j < ria; j++) { usz l=len[j]; rp[j] = !l ? inc(z) : m_shChangeLen(xt, xr, xsh, l, width, csz); } +} + 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) { + ur xr = RNK(x); + if (isArr(w)&isArr(x) && RNK(w)==1 && xr>=1 && depth(w)==1) { usz wia = IA(w); - usz xia = IA(x); + usz* xsh = SH(x); + usz xia = *xsh; if (wia-xia > 1) thrF("⊔: ≠𝕨 must be either ≠𝕩 or one bigger (%s≡≠𝕨, %s≡≠𝕩)", wia, xia); u8 we = TI(w,elType); if (elInt(we)) { @@ -48,7 +67,8 @@ B group_c2(B t, B w, B x) { 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)); + Arr* rf = m_fillarrp(0); if (xr==1) arr_shVec(rf); else arr_shChangeLen(rf, xr, xsh, 0); + fillarr_setFill(rf, m_f64(0)); B z = taga(rf); fillarr_setFill(r, z); @@ -60,10 +80,19 @@ B group_c2(B t, B w, B x) { goto intvec_ret; } - u8 xe = TI(x,elType); - u8 width = elWidth(xe); - u64 xw; - if (xia>64 && (xw=(u64)xia*width)<=I32_MAX && change1)) { + width *= csz = arr_csz(x); + xl += CTZ(csz); + if (bits && xl>=3) { bits=1; width>>=3; } + if ((csz & (csz-1)) || xl>7) xl = 7; + } + if (xia>64 && notB && !bits && change<(xia*width)/32) { #define C1(F,X ) F##_c1(m_f64(0),X ) #define C2(F,X,W) F##_c2(m_f64(0),X,W) if (wia>xia) w = C2(take, m_f64(xia), w); @@ -86,86 +115,83 @@ B group_c2(B t, B w, B x) { for (usz i = 0; i < wia; i++) len[wp[i]]+=ip[i]; void* xp = tyany_ptr(x); - u8 xt = el2t(xe); - for (usz j = 0; j < ria; j++) { - usz l=len[j]; - if (!l) rp[j]=inc(z); else m_tyarrv(rp+j, width, l, xt); - } - for (usz i=0, k=i0*width; 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)); - 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); - 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 - - 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: { - void* xp = tyany_ptr(x); - 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; - case 4: for (usz i = 0; i < xia; i++) { i32 n = wp[i]; if (n>=0) ((u32*)tyarr_ptr(rp[n]))[pos[n]++] = ((u32*)xp)[i]; } break; - case 8: for (usz i = 0; i < xia; i++) { i32 n = wp[i]; if (n>=0) ((f64*)tyarr_ptr(rp[n]))[pos[n]++] = ((f64*)xp)[i]; } break; - } - break; + decG(ind); + } else { + 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); + x = slash_c2(m_f64(0), m, x); xia = IA(x); + neg = 0; } - case el_bit: case el_B: { + if (TI(w,elType)!=el_i32) w = taga(cpyI32Arr(w)); + i32* wp = i32any_ptr(w); + 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 + + u8 xk = xl - 3; + if (notB && !bits && sort) { + void* xp = tyany_ptr(x); + u64 i=neg*width; for (usz j=0; j=0) ((u8* )tyarr_ptr(rp[n]))[pos[n]++] = ((u8* )xp)[i]; } break; + case 1: 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; + case 2: for (usz i = 0; i < xia; i++) { i32 n = wp[i]; if (n>=0) ((u32*)tyarr_ptr(rp[n]))[pos[n]++] = ((u32*)xp)[i]; } break; + case 3: for (usz i = 0; i < xia; i++) { i32 n = wp[i]; if (n>=0) ((f64*)tyarr_ptr(rp[n]))[pos[n]++] = ((f64*)xp)[i]; } break; + } + } else { for (usz i = 0; i < ria; i++) { - Arr* c = m_fillarrp(len[i]); + usz l = len[i]; + Arr* c = m_fillarrp(l*csz); c->ia = 0; fillarr_setFill(c, inc(xf)); - arr_shVec(c); + if (xr==1) arr_shVec(c); else arr_shChangeLen(c, xr, xsh, l); rp[i] = taga(c); } SLOW2("𝕨⊔𝕩", w, x); SGet(x) - for (usz i = 0; i < xia; i++) { - i32 n = wp[i]; - if (n>=0) fillarr_ptr(a(rp[n]))[pos[n]++] = Get(x, i); + if (csz == 1) { + for (usz i = 0; i < xia; i++) { + i32 n = wp[i]; + if (n>=0) fillarr_ptr(a(rp[n]))[pos[n]++] = Get(x, i); + } + } else { + for (usz i = 0; i < xia; i++) { + i32 n = wp[i]; + if (n<0) continue; + usz p = (pos[n]++)*csz; + B* rnp = fillarr_ptr(a(rp[n])) + p; + for (usz j = 0; j < csz; j++) rnp[j] = Get(x, i*csz + j); + } } - for (usz i = 0; i < ria; i++) a(rp[i])->ia = len[i]; - break; + for (usz i = 0; i < ria; i++) a(rp[i])->ia = len[i]*csz; } } intvec_ret: fillarr_setFill(rf, xf); decG(w); decG(x); TFREE(pos); return taga(r); - } else { + } else if (xr==1) { SLOW2("𝕨⊔𝕩", w, x); SGetU(w) i64 ria = wia==xia? 0 : o2i64(GetU(w, xia)); @@ -213,7 +239,6 @@ B group_c2(B t, B w, B x) { base: return c2(rt_group, w, x); } -B ud_c1(B, B); B group_c1(B t, B x) { if (isArr(x) && RNK(x)==1 && TI(x,arrD1)) { usz ia = IA(x); @@ -222,4 +247,3 @@ B group_c1(B t, B x) { } return c1(rt_group, x); } - From 26c3dbfeea7da0c4748bd42f9e03668bd894fabf Mon Sep 17 00:00:00 2001 From: Marshall Lochbaum Date: Fri, 16 Sep 2022 19:26:55 -0400 Subject: [PATCH 06/20] Start using type-independent utilities in compress() --- src/builtins/slash.c | 34 ++++++++++++++++++---------------- 1 file changed, 18 insertions(+), 16 deletions(-) diff --git a/src/builtins/slash.c b/src/builtins/slash.c index 208aa15d..83df1c70 100644 --- a/src/builtins/slash.c +++ b/src/builtins/slash.c @@ -261,24 +261,25 @@ static B where(B x, usz xia, u64 s) { return r; } -static B compress(B w, B x, usz wia, B xf) { +static B compress(B w, B x, usz wia) { 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; @@ -318,19 +319,20 @@ static B compress(B w, B x, usz wia, B xf) { #define WITH_SPARSE(W, CUTOFF, DENSE) { \ i##W *xp=tyany_ptr(x), *rp; \ if (wsum>=wia/CUTOFF) { DENSE; } \ - else { rp=m_tyarrv(&r,W/8,wsum,el2t(xe)); COMPRESS_BLOCK(i##W); } \ + else { rp=m_tyarrv(&r,W/8,wsum,xt); 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)) + case 3: WITH_SPARSE( 8, 32, rp=m_tyarrvO(&r,1,wsum,xt, 8); bmipopc_2slash8 (wp, xp, rp, wia)) + case 4: WITH_SPARSE(16, 16, rp=m_tyarrvO(&r,2,wsum,xt, 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: Fri, 16 Sep 2022 20:02:21 -0400 Subject: [PATCH 07/20] Sparse Indices prototype --- src/builtins/grade.h | 17 +++++++---------- src/builtins/slash.c | 20 +++++++++++++------- 2 files changed, 20 insertions(+), 17 deletions(-) diff --git a/src/builtins/grade.h b/src/builtins/grade.h index 27c7e078..5d7d634c 100644 --- a/src/builtins/grade.h +++ b/src/builtins/grade.h @@ -35,16 +35,13 @@ TALLOC(usz, c0, C); usz *c0o=c0+C/2; \ for (usz j=0; je; j--) { i+=c0[j]; rp[i]--; } \ - ) \ - for (usz i=1; i0 && !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 { + if (s/16 <= xia) { + usz b = 1<<10; + for (usz k=0, j=0, js=0, ij=xp[0]; ; ) { + usz e = b0 && !xp[xia-1]) xia--; + for (u64 i = 0; i < xia; i++) { + i32 c = xp[i]; for (i32 j = 0; j < c; j++) *rp++ = i; } } From 5f196c49d9d87f55c6fe86fe9cf10d9e5cb441f7 Mon Sep 17 00:00:00 2001 From: Marshall Lochbaum Date: Fri, 16 Sep 2022 20:08:23 -0400 Subject: [PATCH 08/20] Extend sparse and dense where to all cases fitting in i32s --- src/builtins/slash.c | 71 ++++++++++++++++++++++---------------------- 1 file changed, 36 insertions(+), 35 deletions(-) diff --git a/src/builtins/slash.c b/src/builtins/slash.c index b7bf7794..3840e786 100644 --- a/src/builtins/slash.c +++ b/src/builtins/slash.c @@ -361,54 +361,55 @@ B slash_c1(B t, B x) { usz xia = IA(x); B r; u8 xe = TI(x,elType); + if (xe!=el_bit && s<=xia) { x = num_squeeze(x); xe = TI(x,elType); } if (xe==el_bit) { r = where(x, xia, s); - } else if (RARE(xia>=I32_MAX)) { + } else if (RARE(xia > (usz)I32_MAX+1)) { 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 { + } else if (RARE(xe > el_i32)) { i32* rp; r = m_i32arrv(&rp, s); - if (xe==el_i8) { - i8* xp = i8any_ptr(x); - if (s/16 <= xia) { - usz b = 1<<10; - for (usz k=0, j=0, js=0, ij=xp[0]; ; ) { - usz e = b0 && !xp[xia-1]) xia--; - for (u64 i = 0; i < xia; i++) { - i32 c = xp[i]; - for (i32 j = 0; j < c; j++) *rp++ = i; + i32* rp; r = m_i32arrv(&rp, s); + if (xe == el_i8 ) { SPARSE_IND(i8 ); } + else if (xe == el_i16) { SPARSE_IND(i16); } + else { SPARSE_IND(i32); } + #undef SPARSE_IND + } else { // Dense case: only result type matters + #define DENSE_IND(T) \ + T* rp; r = m_##T##arrv(&rp, s); \ + for (u64 i = 0; i < xia; i++) { \ + i32 c = xp[i]; \ + for (i32 j = 0; j < c; j++) *rp++ = i; \ } - } - } else if (xe==el_i32) { + if (xe < el_i32) x = taga(cpyI32Arr(x)); 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; - } + if (xia <= 128) { DENSE_IND(i8 ); } + else if (xia <= 32768) { DENSE_IND(i16); } + else { DENSE_IND(i32); } + #undef DENSE_IND } } decG(x); From c72b9b5abd03644d469b36e4c9506498efc191aa Mon Sep 17 00:00:00 2001 From: Marshall Lochbaum Date: Sat, 17 Sep 2022 10:22:34 -0400 Subject: [PATCH 09/20] Redo replicate with sparse and dense cases --- src/builtins/slash.c | 166 +++++++++++++++++++++++++++---------------- 1 file changed, 104 insertions(+), 62 deletions(-) diff --git a/src/builtins/slash.c b/src/builtins/slash.c index 3840e786..c7caa0f4 100644 --- a/src/builtins/slash.c +++ b/src/builtins/slash.c @@ -424,74 +424,116 @@ B slash_c2(B t, B w, B x) { if (wia==0) { decG(w); return x; } thrF("/: Lengths of components of 𝕨 must match 𝕩 (%s ≠ %s)", wia, xia); } - - if (TI(w,elType)==el_bit) { + + u8 we = TI(w,elType); + if (we > el_i32) { w = any_squeeze(w); we = TI(w,elType); } + if (we==el_bit) { + wbool: B r = compress(w, x, wia); decG(w); decG(x); return r; } - B xf = getFillQ(x); - #define CASE(WT,XT) if (TI(x,elType)==el_##XT) { \ - XT* xp = XT##any_ptr(x); \ - XT* rp; B r = m_##XT##arrv(&rp, wsum); \ - if (or<2) for (usz i = 0; i < wia; i++) { \ - *rp = xp[i]; \ - rp+= wp[i]; \ - } else for (usz i = 0; i < wia; i++) { \ - WT cw = wp[i]; XT cx = xp[i]; \ - for (i64 j = 0; j < cw; j++) *rp++ = cx; \ - } \ - decG(w); decG(x); return r; \ + u64 s = usum(w); + if (we!=el_bit && s<=wia) { + w = num_squeeze(w); we = TI(w,elType); + if (we==el_bit) goto wbool; } - #define TYPED(WT,SIGN) { \ - WT* wp = WT##any_ptr(w); \ - while (wia>0 && !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); + B r; + u8 xe = TI(x,elType); + if (RARE(we>el_i32 || xe==el_B)) { // Slow case + SLOW2("𝕨/𝕩", w, x); + B xf = getFillQ(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); + } else if (xe == el_bit) { + u64* xp = bitarr_ptr(x); + u64* rp; r = m_bitarrv(&rp, s); + if (s/256 <= wia) { + #define SPARSE_REP(T) \ + T* wp = T##any_ptr(w); \ + usz b = 1<<12; \ + u64 xx=xp[0], xs=xx>>63, js=-(xx&1); xx^=xx<<1; \ + for (usz k=0, j=0, ij=wp[0]; ; ) { \ + usz e = b>=1; j++; if (j%64==0) { u64 v=xp[j/64]; xx=v^(v<<1)^xs; xs=v>>63; } \ + rp[ij/64]^=(-(xx&1))<<(ij%64); ij+=wp[j]; \ + } \ + for (usz i=k/64; i>63); \ + if (e==s) break; k=e; \ + } + if (we == el_i8 ) { SPARSE_REP(i8 ); } + else if (we == el_i16) { SPARSE_REP(i16); } + else { SPARSE_REP(i32); } + #undef SPARSE_REP + } else { + if (we < el_i32) w = taga(cpyI32Arr(w)); + i32* wp = i32any_ptr(w); + u64 ri=0, rc=0, xc=0; usz j=0; + for (usz i = 0; i < wia; i++) { + u64 v = -(u64)bitp_get(xp,i); + rc ^= (v^xc) << (ri%64); + xc = v; + ri += wp[i]; usz e = ri/64; + if (j < e) { + rp[j++] = rc; + while (j < e) rp[j++] = v; + rc = v; + } + } + if (ri%64) rp[j] = rc; + } + } else { + u8 xt = TY(x); + u8 xl = arrTypeBitsLog(xt)-3; + void* rv = m_tyarrv(&r, 1<0 && !wp[wia-1]) wia--; + switch (xl) { default: UD; CASE(0,u8) CASE(1,u16) CASE(2,u32) CASE(3,u64) } + #undef CASE } } - decG(w); decG(x); - return withFill(HARR_FV(r), xf); + decG(w); decG(x); return r; } if (isArr(x) && RNK(x)==1 && q_i32(w)) { usz xia = IA(x); From ea1367e639a7e39b1b966c4ed4f7cd6917a856e8 Mon Sep 17 00:00:00 2001 From: Marshall Lochbaum Date: Sat, 17 Sep 2022 20:54:27 -0400 Subject: [PATCH 10/20] Handle cell size 0 correctly in cellWidthLog() --- src/core/arrFns.h | 5 ++--- 1 file changed, 2 insertions(+), 3 deletions(-) diff --git a/src/core/arrFns.h b/src/core/arrFns.h index 9ee17682..20652800 100644 --- a/src/core/arrFns.h +++ b/src/core/arrFns.h @@ -65,7 +65,6 @@ static u8 cellWidthLog(B x) { u8 lw = arrTypeBitsLog(TY(x)); if (LIKELY(RNK(x)==1)) return lw; usz csz = arr_csz(x); - if (csz & (csz-1)) return 7; // Not power of 2 - lw += CTZ(csz); - return lw<7? lw : 7; + if (csz & (csz-1)) return 7; // Not power of 2 + return lw + CTZ(csz | 128>>lw); // Max of 7; also handle csz==0 } From cae65947cde1be3e318e9a0a8355484890f4ad8a Mon Sep 17 00:00:00 2001 From: Marshall Lochbaum Date: Sat, 17 Sep 2022 20:36:43 -0400 Subject: [PATCH 11/20] Apply fast Replicate code whenever the cell size is right --- src/builtins/slash.c | 71 ++++++++++++++++++++++++++++++-------------- 1 file changed, 49 insertions(+), 22 deletions(-) diff --git a/src/builtins/slash.c b/src/builtins/slash.c index c7caa0f4..ae9bafe3 100644 --- a/src/builtins/slash.c +++ b/src/builtins/slash.c @@ -261,13 +261,19 @@ static B where(B x, usz xia, u64 s) { return r; } -static B compress(B w, B x, usz wia) { +extern B take_c2(B, B, B); +static B compress(B w, B x, usz wia, u8 xl, u8 xt) { u64* wp = bitarr_ptr(w); u64 we = 0; usz ie = wia/64; usz q=wia%64; if (q) we = wp[ie] &= ((u64)1<1) return take_c2(m_f64(0), m_f64(0), inc(x)); + u8 xe = TI(x,elType); + if (xe != el_B) return elNum(xe)? emptyIVec() : emptyCVec(); + B xf = getFillQ(x); return q_N(xf)? emptyHVec() : m_emptyFVec(xf); + } we = wp[--ie]; } usz wia0 = wia; @@ -276,8 +282,6 @@ static B compress(B w, B x, usz wia) { if (wsum == wia0) return inc(x); B r; - u8 xl = arrTypeBitsLog(TY(x)); - u8 xt = arrNewType(TY(x)); switch(xl) { default: UD; case 0: { u64* xp = bitarr_ptr(x); u64* rp; @@ -330,7 +334,7 @@ static B compress(B w, B x, usz wia) { #endif #undef WITH_SPARSE case 5: { i32* xp= tyany_ptr(x); i32* rp=m_tyarrv(&r,4,wsum,xt); COMPRESS_BLOCK(i32); break; } - case 6: if (TI(x,elType)!=el_B) { f64* xp=f64any_ptr(x); f64* rp; r = m_f64arrv(&rp,wsum); COMPRESS_BLOCK(f64); break; } + case 6: if (TI(x,elType)!=el_B) { u64* xp=tyany_ptr(x); u64* rp=m_tyarrv(&r,8,wsum,xt); COMPRESS_BLOCK(u64); break; } else { B xf = getFillQ(x); B* xp = arr_bptr(x); @@ -349,6 +353,13 @@ static B compress(B w, B x, usz wia) { } #undef COMPRESS_BLOCK } + ur xr = RNK(x); + if (xr > 1) { + Arr* ra=a(r); SPRNK(ra,xr); + usz* sh = PSH(ra) = m_shArr(xr)->a; + sh[0] = PIA(ra); PIA(ra) *= arr_csz(x); + shcpy(sh+1, SH(x)+1, xr-1); + } return r; } @@ -417,19 +428,23 @@ B slash_c1(B t, B x) { } B slash_c2(B t, B w, B x) { - if (isArr(x) && RNK(x)==1 && isArr(w) && RNK(w)==1 && depth(w)==1) { + if (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); - } + if (wia==0) { decG(w); return isArr(x)? x : m_atomUnit(x); } + if (isAtm(x) || RNK(x)==0) thrM("/: 𝕩 must have rank at least 1 for simple 𝕨"); + ur xr = RNK(x); + usz xlen = *SH(x); + if (RARE(wia!=xlen)) thrF("/: Lengths of components of 𝕨 must match 𝕩 (%s ≠ %s)", wia, xlen); + + u8 xl = cellWidthLog(x); + u8 xt = arrNewType(TY(x)); + if (xl > 6 || (xl < 3 && xl != 0)) goto base; u8 we = TI(w,elType); if (we > el_i32) { w = any_squeeze(w); we = TI(w,elType); } if (we==el_bit) { wbool: - B r = compress(w, x, wia); + B r = compress(w, x, wia, xl, xt); decG(w); decG(x); return r; } u64 s = usum(w); @@ -437,9 +452,9 @@ B slash_c2(B t, B w, B x) { w = num_squeeze(w); we = TI(w,elType); if (we==el_bit) goto wbool; } - B r; - u8 xe = TI(x,elType); - if (RARE(we>el_i32 || xe==el_B)) { // Slow case + + if (RARE(we>el_i32 || TI(x,elType)==el_B)) { // Slow case + if (xr > 1) goto base; SLOW2("𝕨/𝕩", w, x); B xf = getFillQ(x); u64 ria = usum(w); @@ -454,9 +469,20 @@ B slash_c2(B t, B w, B x) { } decG(w); decG(x); return withFill(HARR_FV(r), xf); - } else if (xe == el_bit) { + } + + B r; + // Make shape if needed; all cases below use it + usz* rsh = NULL; + if (xr > 1) { + usz* sh = rsh = m_shArr(xr)->a; + sh[0] = s; + shcpy(sh+1, SH(x)+1, xr-1); + } + + if (xl == 0) { u64* xp = bitarr_ptr(x); - u64* rp; r = m_bitarrv(&rp, s); + u64* rp; r = m_bitarrv(&rp, s); if (rsh) { SPRNK(a(r),xr); SH(r) = rsh; } if (s/256 <= wia) { #define SPARSE_REP(T) \ T* wp = T##any_ptr(w); \ @@ -495,9 +521,9 @@ B slash_c2(B t, B w, B x) { if (ri%64) rp[j] = rc; } } else { - u8 xt = TY(x); - u8 xl = arrTypeBitsLog(xt)-3; - void* rv = m_tyarrv(&r, 1<0 && !wp[wia-1]) wia--; - switch (xl) { default: UD; CASE(0,u8) CASE(1,u16) CASE(2,u32) CASE(3,u64) } + switch (xk) { default: UD; CASE(0,u8) CASE(1,u16) CASE(2,u32) CASE(3,u64) } #undef CASE } } @@ -563,6 +589,7 @@ B slash_c2(B t, B w, B x) { return withFill(r.b, xf); } } + base: return c2(rt_slash, w, x); } From 12c4fd0b078fda74afb14e448f1c0011731a9750 Mon Sep 17 00:00:00 2001 From: Marshall Lochbaum Date: Sun, 18 Sep 2022 16:15:20 -0400 Subject: [PATCH 12/20] Vectorizable usum code for i8, i16, i32 --- src/core/numarr.h | 6 +++--- 1 file changed, 3 insertions(+), 3 deletions(-) diff --git a/src/core/numarr.h b/src/core/numarr.h index 1a0f050b..b029ed7b 100644 --- a/src/core/numarr.h +++ b/src/core/numarr.h @@ -101,9 +101,9 @@ static u64 usum(B x) { // doesn't consume; may error usz xia = IA(x); u8 xe = TI(x,elType); if (xe==el_bit) return bit_sum(bitarr_ptr(x), xia); - else if (xe==el_i8 ) { i8* p = i8any_ptr (x); for (usz i = 0; i < xia; i++) { if (RARE(p[i]<0)) goto neg; r+= p[i]; } } - else if (xe==el_i16) { i16* p = i16any_ptr(x); for (usz i = 0; i < xia; i++) { if (RARE(p[i]<0)) goto neg; if (addOn(r,p[i])) goto overflow; } } - else if (xe==el_i32) { i32* p = i32any_ptr(x); for (usz i = 0; i < xia; i++) { if (RARE(p[i]<0)) goto neg; if (addOn(r,p[i])) goto overflow; } } + else if (xe==el_i8 ) { i8* p = i8any_ptr (x); i8 m=0; for (usz i = 0; i < xia; ) { usz b=1<< 8; i16 s=0; for (usz e = xia-i Date: Sun, 18 Sep 2022 16:43:19 -0400 Subject: [PATCH 13/20] Split up partitioning loop in Group for rank 1 and rank>1 cases --- src/builtins/group.c | 22 +++++++++++++--------- 1 file changed, 13 insertions(+), 9 deletions(-) diff --git a/src/builtins/group.c b/src/builtins/group.c index aa84021c..097c5cdd 100644 --- a/src/builtins/group.c +++ b/src/builtins/group.c @@ -141,15 +141,19 @@ B group_c2(B t, B w, B x) { u8 xk = xl - 3; if (notB && !bits && sort) { void* xp = tyany_ptr(x); - u64 i=neg*width; for (usz j=0; j Date: Sun, 18 Sep 2022 21:42:24 -0400 Subject: [PATCH 14/20] Grouped compress for odd cell widths --- src/builtins/slash.c | 57 ++++++++++++++++++++++++++++++++++++++++++-- 1 file changed, 55 insertions(+), 2 deletions(-) diff --git a/src/builtins/slash.c b/src/builtins/slash.c index ae9bafe3..e63c65a7 100644 --- a/src/builtins/slash.c +++ b/src/builtins/slash.c @@ -178,6 +178,57 @@ static void where_block_u16(u64* src, u16* dst, usz len, usz sum) { } } +static B compress_grouped(u64* wp, B x, usz wia, usz wsum, u8 xt) { + B r; + usz csz = arr_csz(x); + u8 xl = arrTypeBitsLog(TY(x)); + #define COMPRESS_GROUP(CPY) \ + u64 ri = 0; \ + u64 wv = wp[0]; usz i = 0, wn = (wia-1)/64+1; \ + for (u64 e=wsum*width; ri < e; ) { \ + while (wv== 0) wv=wp[++i]; usz i0=64*i+CTZ( wv); wv|=wv-1; \ + while (wv==-(u64)1) wv=++i0 || csz%8==0) { // Full bytes + u64 width = xl==0 ? csz/8 : csz << (xl-3); + u8* xp; u8* rp; + bool is_B = TI(x,elType) == el_B; + if (!is_B) { + xp = tyany_ptr(x); + rp = m_tyarrv(&r,width,wsum,xt); + } else { + B xf = getFillQ(x); + xp = (u8*)arr_bptr(x); + usz ria = wsum*csz; + if (xp != NULL) { + HArr_p rh = m_harrUv(ria); + r = withFill(rh.b, xf); + IA(r) = wsum; // Shape-setting code at end of compress expects this + rp = (u8*)rh.a; + } else { + SLOW2("𝕨/𝕩", w, x); + M_HARR(rp, ria) SGet(x) + for (usz i = 0; i < wia; i++) if (bitp_get(wp,i)) { + for (usz j = 0; j < csz; j++) HARR_ADDA(rp, Get(x,i*csz+j)); + } + return withFill(HARR_FV(rp), xf); + } + } + #define MEM_CPY(R,RI,X,XI,L) memcpy(R+RI, X+XI, L) + COMPRESS_GROUP(MEM_CPY) + #undef MEM_CPY + if (is_B) for (usz i = 0; i < wsum*csz; i++) inc(((B*)rp)[i]); + } else { // Bits + usz width = csz; + u64* xp = tyany_ptr(x); + u64* rp; r = m_bitarrv(&rp,wsum*width); IA(r) = wsum; + COMPRESS_GROUP(bit_cpy) + } + return r; +} + static B where(B x, usz xia, u64 s) { B r; u64* xp = bitarr_ptr(x); @@ -282,7 +333,8 @@ static B compress(B w, B x, usz wia, u8 xl, u8 xt) { if (wsum == wia0) return inc(x); B r; - switch(xl) { default: UD; + switch(xl) { + default: r = compress_grouped(wp, x, wia, wsum, xt); break; case 0: { u64* xp = bitarr_ptr(x); u64* rp; #if SINGELI && defined(__BMI2__) @@ -438,15 +490,16 @@ B slash_c2(B t, B w, B x) { u8 xl = cellWidthLog(x); u8 xt = arrNewType(TY(x)); - if (xl > 6 || (xl < 3 && xl != 0)) goto base; u8 we = TI(w,elType); if (we > el_i32) { w = any_squeeze(w); we = TI(w,elType); } if (we==el_bit) { wbool: + if (xl > 6 && TI(x,elType) == el_B) goto base; // TODO: fix bugs, enable B r = compress(w, x, wia, xl, xt); decG(w); decG(x); return r; } + if (xl > 6 || (xl < 3 && xl != 0)) goto base; u64 s = usum(w); if (we!=el_bit && s<=wia) { w = num_squeeze(w); we = TI(w,elType); From 11e3db178705dff4b6eb677ac772cd83dbe61e60 Mon Sep 17 00:00:00 2001 From: Marshall Lochbaum Date: Mon, 19 Sep 2022 11:15:47 -0400 Subject: [PATCH 15/20] Can't call harr withFill before adding data to array --- src/builtins/slash.c | 16 ++++++++-------- 1 file changed, 8 insertions(+), 8 deletions(-) diff --git a/src/builtins/slash.c b/src/builtins/slash.c index e63c65a7..a6b454c7 100644 --- a/src/builtins/slash.c +++ b/src/builtins/slash.c @@ -194,18 +194,15 @@ static B compress_grouped(u64* wp, B x, usz wia, usz wsum, u8 xt) { if (xl>0 || csz%8==0) { // Full bytes u64 width = xl==0 ? csz/8 : csz << (xl-3); u8* xp; u8* rp; - bool is_B = TI(x,elType) == el_B; + bool is_B = TI(x,elType) == el_B; HArr_p rh; if (!is_B) { xp = tyany_ptr(x); rp = m_tyarrv(&r,width,wsum,xt); } else { - B xf = getFillQ(x); xp = (u8*)arr_bptr(x); usz ria = wsum*csz; if (xp != NULL) { - HArr_p rh = m_harrUv(ria); - r = withFill(rh.b, xf); - IA(r) = wsum; // Shape-setting code at end of compress expects this + rh = m_harrUv(ria); rp = (u8*)rh.a; } else { SLOW2("𝕨/𝕩", w, x); @@ -213,13 +210,17 @@ static B compress_grouped(u64* wp, B x, usz wia, usz wsum, u8 xt) { for (usz i = 0; i < wia; i++) if (bitp_get(wp,i)) { for (usz j = 0; j < csz; j++) HARR_ADDA(rp, Get(x,i*csz+j)); } - return withFill(HARR_FV(rp), xf); + return withFill(HARR_FV(rp), getFillQ(x)); } } #define MEM_CPY(R,RI,X,XI,L) memcpy(R+RI, X+XI, L) COMPRESS_GROUP(MEM_CPY) #undef MEM_CPY - if (is_B) for (usz i = 0; i < wsum*csz; i++) inc(((B*)rp)[i]); + if (is_B) { + for (usz i = 0; i < wsum*csz; i++) inc(((B*)rp)[i]); + r = withFill(rh.b, getFillQ(x)); + IA(r) = wsum; // Shape-setting code at end of compress expects this + } } else { // Bits usz width = csz; u64* xp = tyany_ptr(x); @@ -495,7 +496,6 @@ B slash_c2(B t, B w, B x) { if (we > el_i32) { w = any_squeeze(w); we = TI(w,elType); } if (we==el_bit) { wbool: - if (xl > 6 && TI(x,elType) == el_B) goto base; // TODO: fix bugs, enable B r = compress(w, x, wia, xl, xt); decG(w); decG(x); return r; } From ac4729bb29709e6fce4adf89fbeb130df05ce611 Mon Sep 17 00:00:00 2001 From: Marshall Lochbaum Date: Mon, 19 Sep 2022 19:31:44 -0400 Subject: [PATCH 16/20] Use grouped compress for all byte sizes if the number of groups is small enough --- src/builtins/slash.c | 29 +++++++++++++++++++++++++---- 1 file changed, 25 insertions(+), 4 deletions(-) diff --git a/src/builtins/slash.c b/src/builtins/slash.c index a6b454c7..1eda65be 100644 --- a/src/builtins/slash.c +++ b/src/builtins/slash.c @@ -313,6 +313,21 @@ static B where(B x, usz xia, u64 s) { return r; } +// Is the number of values switches in w at most max? +static bool groups_lt(u64* wp, usz len, usz max) { + usz r = 0; + u64 prev = 0; + usz l = (len-1)/64+1; // assume trailing bits are zeroed out + usz b = 1<<8; + for (usz i = 0; i < l; ) { + for (usz e = l-i>63; + } + if (r > max) return 0; + } + return 1; +} + extern B take_c2(B, B, B); static B compress(B w, B x, usz wia, u8 xl, u8 xt) { u64* wp = bitarr_ptr(w); @@ -375,8 +390,9 @@ static B compress(B w, B x, usz wia, u8 xl, u8 xt) { TFREE(buf) #define WITH_SPARSE(W, CUTOFF, DENSE) { \ i##W *xp=tyany_ptr(x), *rp; \ - if (wsum>=wia/CUTOFF) { DENSE; } \ - else { rp=m_tyarrv(&r,W/8,wsum,xt); COMPRESS_BLOCK(i##W); } \ + if (wsum=wia/8 && groups_lt(wp,wia, wia/16)) r = compress_grouped(wp, x, wia, wsum, xt); \ + else { T* xp=tyany_ptr(x); T* rp=m_tyarrv(&r,sizeof(T),wsum,xt); COMPRESS_BLOCK(T); break; } \ + break; + case 5: BLOCK_OR_GROUPED(i32) + case 6: if (TI(x,elType)!=el_B) { BLOCK_OR_GROUPED(u64) } + #undef BLOCK_OR_GROUPED else { B xf = getFillQ(x); B* xp = arr_bptr(x); From 01f73dc83b192981f55b570468b1ab9f30573ba4 Mon Sep 17 00:00:00 2001 From: Marshall Lochbaum Date: Mon, 19 Sep 2022 19:38:26 -0400 Subject: [PATCH 17/20] Move usum and bit_sum from numarr.h to stuff.c --- src/core/numarr.h | 39 ++------------------------------------- src/core/stuff.c | 38 ++++++++++++++++++++++++++++++++++++++ 2 files changed, 40 insertions(+), 37 deletions(-) diff --git a/src/core/numarr.h b/src/core/numarr.h index b029ed7b..ce47a59e 100644 --- a/src/core/numarr.h +++ b/src/core/numarr.h @@ -88,40 +88,5 @@ static B toF64Any(B x) { u8 t=TY(x); return t==t_f64arr || t==t_f64slice? x : ta B m_cai32(usz ia, i32* a); B m_caf64(usz sz, f64* a); -static i64 bit_sum(u64* x, u64 am) { - i64 r = 0; - for (u64 i = 0; i < (am>>6); i++) r+= POPC(x[i]); - if (am&63) r+= POPC(x[am>>6]<<(64-am & 63)); - return r; -} - -static u64 usum(B x) { // doesn't consume; may error - assert(isArr(x)); - u64 r = 0; - usz xia = IA(x); - u8 xe = TI(x,elType); - if (xe==el_bit) return bit_sum(bitarr_ptr(x), xia); - else if (xe==el_i8 ) { i8* p = i8any_ptr (x); i8 m=0; for (usz i = 0; i < xia; ) { usz b=1<< 8; i16 s=0; for (usz e = xia-i>6); i++) r+= POPC(x[i]); + if (am&63) r+= POPC(x[am>>6]<<(64-am & 63)); + return r; +} + +u64 usum(B x) { // doesn't consume; may error + assert(isArr(x)); + u64 r = 0; + usz xia = IA(x); + u8 xe = TI(x,elType); + if (xe==el_bit) return bit_sum(bitarr_ptr(x), xia); + else if (xe==el_i8 ) { i8* p = i8any_ptr (x); i8 m=0; for (usz i = 0; i < xia; ) { usz b=1<< 8; i16 s=0; for (usz e = xia-i Date: Tue, 20 Sep 2022 17:31:39 +0300 Subject: [PATCH 18/20] slash.c cleanup --- src/builtins/slash.c | 81 +++++++++++++++++++++++--------------------- src/core/arrFns.h | 2 +- 2 files changed, 44 insertions(+), 39 deletions(-) diff --git a/src/builtins/slash.c b/src/builtins/slash.c index 1eda65be..8b415676 100644 --- a/src/builtins/slash.c +++ b/src/builtins/slash.c @@ -183,13 +183,13 @@ static B compress_grouped(u64* wp, B x, usz wia, usz wsum, u8 xt) { usz csz = arr_csz(x); u8 xl = arrTypeBitsLog(TY(x)); #define COMPRESS_GROUP(CPY) \ - u64 ri = 0; \ - u64 wv = wp[0]; usz i = 0, wn = (wia-1)/64+1; \ - for (u64 e=wsum*width; ri < e; ) { \ - while (wv== 0) wv=wp[++i]; usz i0=64*i+CTZ( wv); wv|=wv-1; \ - while (wv==-(u64)1) wv=++i0 || csz%8==0) { // Full bytes u64 width = xl==0 ? csz/8 : csz << (xl-3); @@ -474,7 +474,7 @@ B slash_c1(B t, B x) { for (usz i=k; i el_i32) { w = any_squeeze(w); we = TI(w,elType); } + if (!elInt(we)) { + w=any_squeeze(w); we=TI(w,elType); + if (!elInt(we)) goto slow; + } if (we==el_bit) { wbool: - B r = compress(w, x, wia, xl, xt); - decG(w); decG(x); return r; + r = compress(w, x, wia, xl, xt); + goto decWX_ret; } - if (xl > 6 || (xl < 3 && xl != 0)) goto base; + if (xl>6 || (xl<3 && xl!=0)) goto base; u64 s = usum(w); - if (we!=el_bit && s<=wia) { - w = num_squeeze(w); we = TI(w,elType); + if (s<=wia) { + w=num_squeezeChk(w); we=TI(w,elType); if (we==el_bit) goto wbool; } - - if (RARE(we>el_i32 || TI(x,elType)==el_B)) { // Slow case + + if (RARE(TI(x,elType)==el_B)) { // Slow case + slow: if (xr > 1) goto base; SLOW2("𝕨/𝕩", w, x); B xf = getFillQ(x); - u64 ria = usum(w); - if (ria>=USZ_MAX) thrOOM(); - M_HARR(r, ria) SGetU(w) SGetU(x) + MAKE_MUT(r0, s) mut_init(r0, el_B); MUTG_INIT(r0); + SGetU(w) SGetU(x) + usz ri = 0; 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); + mut_fillG(r0, ri, GetU(x, i), c); + ri+= c; } } - decG(w); decG(x); - return withFill(HARR_FV(r), xf); + r = withFill(mut_fv(r0), xf); + decWX_ret: decG(w); + decX_ret: decG(x); + return r; } - - B r; + // Make shape if needed; all cases below use it usz* rsh = NULL; if (xr > 1) { @@ -571,11 +577,11 @@ B slash_c2(B t, B w, B x) { rp[ij/64]^=(-(xx&1))<<(ij%64); ij+=wp[j]; \ } \ for (usz i=k/64; i>63); \ - if (e==s) break; k=e; \ + if (e==s) {break;} k=e; \ } - if (we == el_i8 ) { SPARSE_REP(i8 ); } - else if (we == el_i16) { SPARSE_REP(i16); } - else { SPARSE_REP(i32); } + if (we==el_i8 ) { SPARSE_REP(i8 ); } + else if (we==el_i16) { SPARSE_REP(i16); } + else { SPARSE_REP(i32); } #undef SPARSE_REP } else { if (we < el_i32) w = taga(cpyI32Arr(w)); @@ -609,7 +615,7 @@ B slash_c2(B t, B w, B x) { for (usz i=k; i Date: Tue, 20 Sep 2022 17:36:55 +0300 Subject: [PATCH 19/20] =?UTF-8?q?fix=20=F0=9D=95=A8=E2=8A=940?= MIME-Version: 1.0 Content-Type: text/plain; charset=UTF-8 Content-Transfer-Encoding: 8bit --- src/builtins/group.c | 3 ++- 1 file changed, 2 insertions(+), 1 deletion(-) diff --git a/src/builtins/group.c b/src/builtins/group.c index 097c5cdd..8468a20b 100644 --- a/src/builtins/group.c +++ b/src/builtins/group.c @@ -30,8 +30,9 @@ static void allocGroups(B* rp, usz ria, B z, u8 xt, ur xr, usz* xsh, i32* len, u extern B rt_group; B group_c2(B t, B w, B x) { + if (!isArr(x)) thrM("⊔: 𝕩 must be an array"); ur xr = RNK(x); - if (isArr(w)&isArr(x) && RNK(w)==1 && xr>=1 && depth(w)==1) { + if (isArr(w) && RNK(w)==1 && xr>=1 && depth(w)==1) { usz wia = IA(w); usz* xsh = SH(x); usz xia = *xsh; From 42e36851b9cb60a7c03b2c02c2aa407b0dc1dcd8 Mon Sep 17 00:00:00 2001 From: dzaima Date: Tue, 20 Sep 2022 18:57:19 +0300 Subject: [PATCH 20/20] more cleanup --- src/builtins/group.c | 31 ++++++++++++++--------------- src/builtins/slash.c | 46 ++++++++++++++++++++++---------------------- 2 files changed, 39 insertions(+), 38 deletions(-) diff --git a/src/builtins/group.c b/src/builtins/group.c index 8468a20b..250f791a 100644 --- a/src/builtins/group.c +++ b/src/builtins/group.c @@ -62,25 +62,25 @@ 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 = m_fillarrp(0); if (xr==1) arr_shVec(rf); else arr_shChangeLen(rf, xr, xsh, 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 (ria==0) goto setfill_dec_ret; if (neg==xia) { for (usz i = 0; i < ria; i++) rp[i] = inc(z); - goto intvec_ret; + goto setfill_dec_ret; } - + TALLOC(i32, pos, 2*ria+1); i32* len = pos+ria+1; + bool notB = TI(x,elType) != el_B; u8 xt = arrNewType(TY(x)); u8 xl = arrTypeBitsLog(TY(x)); @@ -90,7 +90,7 @@ B group_c2(B t, B w, B x) { if (RARE(xr>1)) { width *= csz = arr_csz(x); xl += CTZ(csz); - if (bits && xl>=3) { bits=1; width>>=3; } + if (bits && xl>=3) { bits=0; width>>=3; } if ((csz & (csz-1)) || xl>7) xl = 7; } if (xia>64 && notB && !bits && change<(xia*width)/32) { @@ -106,7 +106,7 @@ B group_c2(B t, B w, B x) { 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]; @@ -114,9 +114,9 @@ B group_c2(B t, B w, B x) { ip[wia-1] = xia-ip[wia-1]; for (usz i = 0; i < ria; i++) len[i] = pos[i] = 0; for (usz i = 0; i < wia; i++) len[wp[i]]+=ip[i]; - + void* xp = tyany_ptr(x); - + allocGroups(rp, ria, z, xt, xr, xsh, len, width, csz); for (u64 i=0, k=i0*width; i=0) ((u8* )tyarr_ptr(rp[n]))[pos[n]++] = ((u8* )xp)[i]; } break; case 1: 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; case 2: for (usz i = 0; i < xia; i++) { i32 n = wp[i]; if (n>=0) ((u32*)tyarr_ptr(rp[n]))[pos[n]++] = ((u32*)xp)[i]; } break; - case 3: for (usz i = 0; i < xia; i++) { i32 n = wp[i]; if (n>=0) ((f64*)tyarr_ptr(rp[n]))[pos[n]++] = ((f64*)xp)[i]; } break; + case 3: for (usz i = 0; i < xia; i++) { i32 n = wp[i]; if (n>=0) ((u64*)tyarr_ptr(rp[n]))[pos[n]++] = ((u64*)xp)[i]; } break; } } else { for (usz i = 0; i < ria; i++) { @@ -192,9 +192,10 @@ B group_c2(B t, B w, B x) { for (usz i = 0; i < ria; i++) a(rp[i])->ia = len[i]*csz; } } - intvec_ret: + TFREE(pos); + setfill_dec_ret: fillarr_setFill(rf, xf); - decG(w); decG(x); TFREE(pos); + decG(w); decG(x); return taga(r); } else if (xr==1) { SLOW2("𝕨⊔𝕩", w, x); diff --git a/src/builtins/slash.c b/src/builtins/slash.c index 8b415676..e40dfd09 100644 --- a/src/builtins/slash.c +++ b/src/builtins/slash.c @@ -347,7 +347,7 @@ static B compress(B w, B x, usz wia, u8 xl, u8 xt) { wia = 64*(ie+1) - CLZ(we); usz wsum = bit_sum(wp, wia); if (wsum == wia0) return inc(x); - + B r; switch(xl) { default: r = compress_grouped(wp, x, wia, wsum, xt); break; @@ -404,27 +404,27 @@ static B compress(B w, B x, usz wia, u8 xl, u8 xt) { #undef WITH_SPARSE #define BLOCK_OR_GROUPED(T) \ if (wsum>=wia/8 && groups_lt(wp,wia, wia/16)) r = compress_grouped(wp, x, wia, wsum, xt); \ - else { T* xp=tyany_ptr(x); T* rp=m_tyarrv(&r,sizeof(T),wsum,xt); COMPRESS_BLOCK(T); break; } \ - break; - case 5: BLOCK_OR_GROUPED(i32) - case 6: if (TI(x,elType)!=el_B) { BLOCK_OR_GROUPED(u64) } - #undef BLOCK_OR_GROUPED - else { - B xf = getFillQ(x); - B* xp = arr_bptr(x); - if (xp!=NULL) { - HArr_p rh = m_harrUv(wsum); - B *rp = rh.a; COMPRESS_BLOCK(B); - for (usz i=0; i