From a0a4ecadc745f0c7b534919a5cc579bb038d59f4 Mon Sep 17 00:00:00 2001 From: Marshall Lochbaum Date: Fri, 26 Aug 2022 16:07:00 -0400 Subject: [PATCH 1/8] 1-byte sorting: insertion, bucket, counting --- src/builtins/sort.c | 79 +++++++++++++++++++++++++++++++++++++++++++++ 1 file changed, 79 insertions(+) diff --git a/src/builtins/sort.c b/src/builtins/sort.c index b5b17ca5..7c274800 100644 --- a/src/builtins/sort.c +++ b/src/builtins/sort.c @@ -27,12 +27,69 @@ typedef struct I32I32p { i32 k; i32 v; } I32I32p; #define SORT_TYPE i32 #include "sortTemplate.h" +#define INSERTION_SORT(T, CMP) \ + rp[0] = xp[0]; \ + for (usz i=0; ie; j--) { i+=c0[j]; rp[i]--; } \ + for (usz i=1; i0; ) for (usz c=c0[j]; c--; ) *rp++ = j-C/2; \ + } + B and_c1(B t, B x) { if (isAtm(x) || RNK(x)==0) thrM("∧: Argument cannot have rank 0"); if (RNK(x)!=1) return bqn_merge(and_c1(t, toCells(x))); usz xia = IA(x); + if (xia <= 1) return x; u8 xe = TI(x,elType); if (xe<=el_i32) { + if (xe==el_i8) { + i8* xp = i8any_ptr(x); usz n=xia; + i8* rp; B r = m_i8arrv(&rp, n); + if (n<16) { + INSERTION_SORT(i8,<); + } else if (n<=256) { + // Radix/bucket sort + TALLOC(u8, c0, 256); u8 *c0o=c0+128; // Offset for signedness + for (usz j=0; j<256; j++) c0[j]=0; + for (usz i=0; i>=56; + } + for (usz i=0; i); + } else if (n<=256) { + // Radix/bucket sort + TALLOC(u8, c0, 256); u8 *c0o=c0+128; // Offset for signedness + for (usz j=0; j<256; j++) c0[j]=0; + for (usz i=0; i>8; s+=s>>16; s+=s>>32; + cw[j] = p|(s>>8); s<<=56; + } + for (usz i=0; i Date: Fri, 26 Aug 2022 17:25:40 -0400 Subject: [PATCH 2/8] Move Sort to template file sort.h (for up/down) like Grade --- src/builtins/sort.c | 151 ++------------------------------------------ src/builtins/sort.h | 103 ++++++++++++++++++++++++++++++ 2 files changed, 107 insertions(+), 147 deletions(-) create mode 100644 src/builtins/sort.h diff --git a/src/builtins/sort.c b/src/builtins/sort.c index 7c274800..d270c073 100644 --- a/src/builtins/sort.c +++ b/src/builtins/sort.c @@ -1,8 +1,6 @@ #include "../core.h" #include "../utils/talloc.h" - - #define CAT0(A,B) A##_##B #define CAT(A,B) CAT0(A,B) typedef struct BI32p { B k; i32 v; } BI32p; @@ -17,148 +15,7 @@ typedef struct I32I32p { i32 k; i32 v; } I32I32p; #define GRADE_CHR "⍒" #include "grade.h" -#define SORT_CMP(W, X) compare(W, X) -#define SORT_NAME bA -#define SORT_TYPE B -#include "sortTemplate.h" - -#define SORT_CMP(W, X) ((W) - (i64)(X)) -#define SORT_NAME iA -#define SORT_TYPE i32 -#include "sortTemplate.h" - -#define INSERTION_SORT(T, CMP) \ - rp[0] = xp[0]; \ - for (usz i=0; ie; j--) { i+=c0[j]; rp[i]--; } \ - for (usz i=1; i0; ) for (usz c=c0[j]; c--; ) *rp++ = j-C/2; \ - } - -B and_c1(B t, B x) { - if (isAtm(x) || RNK(x)==0) thrM("∧: Argument cannot have rank 0"); - if (RNK(x)!=1) return bqn_merge(and_c1(t, toCells(x))); - usz xia = IA(x); - if (xia <= 1) return x; - u8 xe = TI(x,elType); - if (xe<=el_i32) { - if (xe==el_i8) { - i8* xp = i8any_ptr(x); usz n=xia; - i8* rp; B r = m_i8arrv(&rp, n); - if (n<16) { - INSERTION_SORT(i8,<); - } else if (n<=256) { - // Radix/bucket sort - TALLOC(u8, c0, 256); u8 *c0o=c0+128; // Offset for signedness - for (usz j=0; j<256; j++) c0[j]=0; - for (usz i=0; i>=56; - } - for (usz i=0; i); - } else if (n<=256) { - // Radix/bucket sort - TALLOC(u8, c0, 256); u8 *c0o=c0+128; // Offset for signedness - for (usz j=0; j<256; j++) c0[j]=0; - for (usz i=0; i>8; s+=s>>16; s+=s>>32; - cw[j] = p|(s>>8); s<<=56; - } - for (usz i=0; i) + +#define TIM_B SORT_UD(bA,bD) +#define TIM_I SORT_UD(iA,iD) + +#define SORT_CMP(W, X) SORT_UD(compare(W, X), compare(X, W)) +#define SORT_NAME TIM_B +#define SORT_TYPE B +#include "sortTemplate.h" + +#define SORT_CMP(W, X) (SORT_UD((W) - (i64)(X), (X) - (i64)(W))) +#define SORT_NAME TIM_I +#define SORT_TYPE i32 +#include "sortTemplate.h" + +#define INSERTION_SORT(T) \ + rp[0] = xp[0]; \ + for (usz i=0; ie; j--) { i+=c0[j]; rp[i]--; } \ + ) \ + for (usz i=1; i0; )) \ + for (usz c=c0[j]; c--; ) *rp++ = j-C/2; \ + } \ + TFREE(c0) + +B SORT_C1(B t, B x) { + if (isAtm(x) || RNK(x)==0) thrM(SORT_UD("∧","∨")": Argument cannot have rank 0"); + if (RNK(x)!=1) return bqn_merge(SORT_C1(t, toCells(x))); + usz xia = IA(x); + if (xia <= 1) return x; + u8 xe = TI(x,elType); + B r; + if (xe==el_i8) { + i8* xp = i8any_ptr(x); usz n=xia; + i8* rp; r = m_i8arrv(&rp, n); + if (n<16) { + INSERTION_SORT(i8); + } else if (n<=256) { + // Radix/bucket sort + TALLOC(u8, c0, 256); u8 *c0o=c0+128; // Offset for signedness + for (usz j=0; j<256; j++) c0[j]=0; + for (usz i=0; i>=56; + } + , + for (usz j=256/8; j--; ) { + u64 p=s; s+=cw[j]; s+=s>>8; s+=s>>16; s+=s>>32; + cw[j] = p|(s>>8); s<<=56; + } + ) + for (usz i=0; i Date: Fri, 26 Aug 2022 20:22:08 -0400 Subject: [PATCH 3/8] 2-byte sorting: insertion, radix, counting --- src/builtins/sort.h | 64 +++++++++++++++++++++++++++++++++------------ 1 file changed, 47 insertions(+), 17 deletions(-) diff --git a/src/builtins/sort.h b/src/builtins/sort.h index fa08df76..130ae88c 100644 --- a/src/builtins/sort.h +++ b/src/builtins/sort.h @@ -14,6 +14,14 @@ #define SORT_TYPE i32 #include "sortTemplate.h" +#define FOR(I,MAX) SORT_UD(for (usz I=0; I>=), PRE_UD(K,>>,<<=)) + #define INSERTION_SORT(T) \ rp[0] = xp[0]; \ for (usz i=0; i0; )) \ - for (usz c=c0[j]; c--; ) *rp++ = j-C/2; \ + FOR(j,C) for (usz c=c0[j]; c--; ) *rp++ = j-C/2; \ } \ TFREE(c0) @@ -55,30 +62,50 @@ B SORT_C1(B t, B x) { i8* rp; r = m_i8arrv(&rp, n); if (n<16) { INSERTION_SORT(i8); - } else if (n<=256) { - // Radix/bucket sort + } else if (n<=256) { // Radix/bucket sort TALLOC(u8, c0, 256); u8 *c0o=c0+128; // Offset for signedness for (usz j=0; j<256; j++) c0[j]=0; for (usz i=0; i>=56; - } - , - for (usz j=256/8; j--; ) { - u64 p=s; s+=cw[j]; s+=s>>8; s+=s>>16; s+=s>>32; - cw[j] = p|(s>>8); s<<=56; - } - ) + FOR(j, 256/8) { PRE64(); } // Prefix sum for (usz i=0; i>8)]++; } \ + PRE_SUM; \ + i16 *r0 = (i16*)(c0+2*256); \ + for (usz i=0; i>8)]++; rp[c]=v; } \ + TFREE(alloc) + RADIX2(u8, + u64 s0=0; u64 s1=0; u64 *cw0=(u64*)c0; u64 *cw1=(u64*)c1; + FOR(j, 256/8) { PRE64(0); PRE64(1); } + ); + } else if (n < 1<<15) { // Radix sort + RADIX2(u32, + u32 s0=0; u32 s1=0; + FOR(j, 256) { + u32 p0=s0; s0+=c0[j]; c0[j]=p0; + u32 p1=s1; s1+=c1[j]; c1[j]=p1; + } + ); + #undef RADIX2 + } else { + COUNTING_SORT(i16); + } + } else if (xe==el_i32) { i32* xp = i32any_ptr(x); i32* rp; r = m_i32arrv(&rp, xia); memcpy(rp, xp, xia*4); @@ -98,6 +125,9 @@ B SORT_C1(B t, B x) { #undef LT #undef TIM_B #undef TIM_I +#undef FOR +#undef PRE_UD +#undef PRE64 #undef INSERTION_SORT #undef COUNTING_SORT #undef SORT_UD From 9d59dff522bfa66a2f8b80942ad11e165f6bcad5 Mon Sep 17 00:00:00 2001 From: Marshall Lochbaum Date: Fri, 26 Aug 2022 20:54:12 -0400 Subject: [PATCH 4/8] 4-byte sorting: insertion and radix --- src/builtins/sort.h | 67 ++++++++++++++++++++++++++++++--------------- 1 file changed, 45 insertions(+), 22 deletions(-) diff --git a/src/builtins/sort.h b/src/builtins/sort.h index 130ae88c..a2195dba 100644 --- a/src/builtins/sort.h +++ b/src/builtins/sort.h @@ -2,24 +2,18 @@ #define LT SORT_UD(<,>) #define TIM_B SORT_UD(bA,bD) -#define TIM_I SORT_UD(iA,iD) #define SORT_CMP(W, X) SORT_UD(compare(W, X), compare(X, W)) #define SORT_NAME TIM_B #define SORT_TYPE B #include "sortTemplate.h" -#define SORT_CMP(W, X) (SORT_UD((W) - (i64)(X), (X) - (i64)(W))) -#define SORT_NAME TIM_I -#define SORT_TYPE i32 -#include "sortTemplate.h" - #define FOR(I,MAX) SORT_UD(for (usz I=0; I>=), PRE_UD(K,>>,<<=)) #define INSERTION_SORT(T) \ @@ -53,12 +47,12 @@ B SORT_C1(B t, B x) { if (isAtm(x) || RNK(x)==0) thrM(SORT_UD("∧","∨")": Argument cannot have rank 0"); if (RNK(x)!=1) return bqn_merge(SORT_C1(t, toCells(x))); - usz xia = IA(x); - if (xia <= 1) return x; + usz n = IA(x); + if (n <= 1) return x; u8 xe = TI(x,elType); B r; if (xe==el_i8) { - i8* xp = i8any_ptr(x); usz n=xia; + i8* xp = i8any_ptr(x); i8* rp; r = m_i8arrv(&rp, n); if (n<16) { INSERTION_SORT(i8); @@ -66,15 +60,14 @@ B SORT_C1(B t, B x) { TALLOC(u8, c0, 256); u8 *c0o=c0+128; // Offset for signedness for (usz j=0; j<256; j++) c0[j]=0; for (usz i=0; i>8)]++; rp[c]=v; } \ TFREE(alloc) RADIX2(u8, - u64 s0=0; u64 s1=0; u64 *cw0=(u64*)c0; u64 *cw1=(u64*)c1; + u64 s0=0; u64 s1=0; FOR(j, 256/8) { PRE64(0); PRE64(1); } ); } else if (n < 1<<15) { // Radix sort @@ -107,15 +100,46 @@ B SORT_C1(B t, B x) { } } else if (xe==el_i32) { i32* xp = i32any_ptr(x); - i32* rp; r = m_i32arrv(&rp, xia); - memcpy(rp, xp, xia*4); - CAT(TIM_I,tim_sort)(rp, xia); + i32* rp; r = m_i32arrv(&rp, n); + if (n < 40) { + INSERTION_SORT(i32); + } else if (n <= 256) { + #define RADIX4(T, PRE_SUM) \ + TALLOC(u8, alloc, 4*256*sizeof(T) + n*4); \ + T *c0=(T*)alloc, *c1=c0+256, *c2=c1+256, *c3=c2+256, *c3o=c3+128; \ + for (usz j=0; j<4*256; j++) c0[j]=0; \ + for (usz i=0; i> 8)]++; \ + c2 [(u8)(v>>16)]++; c3o[(i8)(v>>24)]++; } \ + PRE_SUM; \ + i32 *r0 = (i32*)(c0+4*256); \ + for (usz i=0; i> 8)]++; rp[c]=v; } \ + for (usz i=0; i>16)]++; r0[c]=v; } \ + for (usz i=0; i>24)]++; rp[c]=v; } \ + TFREE(alloc) + RADIX4(u8, + u64 s0=0; u64 s1=0; u64 s2=0; u64 s3=0; + FOR(j, 256/8) { PRE64(0); PRE64(1); PRE64(2); PRE64(3); } + ); + } else { + RADIX4(usz, + usz s0=0; usz s1=0; usz s2=0; usz s3=0; + FOR(j, 256) { + u32 p0=s0; s0+=c0[j]; c0[j]=p0; + u32 p1=s1; s1+=c1[j]; c1[j]=p1; + u32 p2=s2; s2+=c2[j]; c2[j]=p2; + u32 p3=s3; s3+=c3[j]; c3[j]=p3; + } + ); + #undef RADIX4 + } } else { B xf = getFillQ(x); - HArr_p ra = m_harrUv(xia); + HArr_p ra = m_harrUv(n); SGet(x) - for (usz i = 0; i < xia; i++) ra.a[i] = Get(x,i); - CAT(TIM_B,tim_sort)(ra.a, xia); + for (usz i = 0; i < n; i++) ra.a[i] = Get(x,i); + CAT(TIM_B,tim_sort)(ra.a, n); r = withFill(ra.b,xf); } decG(x); @@ -124,7 +148,6 @@ B SORT_C1(B t, B x) { #undef SORT_C1 #undef LT #undef TIM_B -#undef TIM_I #undef FOR #undef PRE_UD #undef PRE64 From 4ec83e51a69625ae6f36975171cdc2534168f761 Mon Sep 17 00:00:00 2001 From: Marshall Lochbaum Date: Sat, 27 Aug 2022 10:48:21 -0400 Subject: [PATCH 5/8] Merge sort.h into grade.h --- src/builtins/grade.h | 171 +++++++++++++++++++++++++++++++++++++++++-- src/builtins/sort.c | 11 +-- src/builtins/sort.h | 156 --------------------------------------- 3 files changed, 168 insertions(+), 170 deletions(-) delete mode 100644 src/builtins/sort.h diff --git a/src/builtins/grade.h b/src/builtins/grade.h index d0e6a3d4..7d78ea1e 100644 --- a/src/builtins/grade.h +++ b/src/builtins/grade.h @@ -1,4 +1,11 @@ #define GRADE_CAT(N) CAT(GRADE_UD(gradeUp,gradeDown),N) +#define GRADE_NEG GRADE_UD(,-) + +// Timsort +#define SORT_CMP(W, X) GRADE_NEG compare(W, X) +#define SORT_NAME GRADE_UD(bA,bD) +#define SORT_TYPE B +#include "sortTemplate.h" #define SORT_CMP(W, X) GRADE_NEG compare((W).k, (X).k) #define SORT_NAME GRADE_CAT(BP) @@ -10,6 +17,153 @@ #define SORT_TYPE I32I32p #include "sortTemplate.h" + +#define LT GRADE_UD(<,>) +#define FOR(I,MAX) GRADE_UD(for (usz I=0; I>=), PRE_UD(K,>>,<<=)) + +#define INSERTION_SORT(T) \ + rp[0] = xp[0]; \ + for (usz i=0; ie; j--) { i+=c0[j]; rp[i]--; } \ + ) \ + for (usz i=1; i>8)]++; } \ + PRE_SUM; \ + i16 *r0 = (i16*)(c0+2*256); \ + for (usz i=0; i>8)]++; rp[c]=v; } \ + TFREE(alloc) + RADIX2(u8, + u64 s0=0; u64 s1=0; + FOR(j, 256/8) { PRE64(0); PRE64(1); } + ); + } else if (n < 1<<15) { // Radix sort + RADIX2(u32, + u32 s0=0; u32 s1=0; + FOR(j, 256) { + u32 p0=s0; s0+=c0[j]; c0[j]=p0; + u32 p1=s1; s1+=c1[j]; c1[j]=p1; + } + ); + #undef RADIX2 + } else { + COUNTING_SORT(i16); + } + } else if (xe==el_i32) { + i32* xp = i32any_ptr(x); + i32* rp; r = m_i32arrv(&rp, n); + if (n < 40) { + INSERTION_SORT(i32); + } else if (n <= 256) { + #define RADIX4(T, PRE_SUM) \ + TALLOC(u8, alloc, 4*256*sizeof(T) + n*4); \ + T *c0=(T*)alloc, *c1=c0+256, *c2=c1+256, *c3=c2+256, *c3o=c3+128; \ + for (usz j=0; j<4*256; j++) c0[j]=0; \ + for (usz i=0; i> 8)]++; \ + c2 [(u8)(v>>16)]++; c3o[(i8)(v>>24)]++; } \ + PRE_SUM; \ + i32 *r0 = (i32*)(c0+4*256); \ + for (usz i=0; i> 8)]++; rp[c]=v; } \ + for (usz i=0; i>16)]++; r0[c]=v; } \ + for (usz i=0; i>24)]++; rp[c]=v; } \ + TFREE(alloc) + RADIX4(u8, + u64 s0=0; u64 s1=0; u64 s2=0; u64 s3=0; + FOR(j, 256/8) { PRE64(0); PRE64(1); PRE64(2); PRE64(3); } + ); + } else { + RADIX4(usz, + usz s0=0; usz s1=0; usz s2=0; usz s3=0; + FOR(j, 256) { + u32 p0=s0; s0+=c0[j]; c0[j]=p0; + u32 p1=s1; s1+=c1[j]; c1[j]=p1; + u32 p2=s2; s2+=c2[j]; c2[j]=p2; + u32 p3=s3; s3+=c3[j]; c3[j]=p3; + } + ); + #undef RADIX4 + } + } else { + B xf = getFillQ(x); + HArr_p ra = m_harrUv(n); + SGet(x) + for (usz i = 0; i < n; i++) ra.a[i] = Get(x,i); + CAT(GRADE_UD(bA,bD),tim_sort)(ra.a, n); + r = withFill(ra.b,xf); + } + decG(x); + return FL_SET(r, CAT(fl,GRADE_UD(asc,dsc))); +} +#undef SORT_C1 +#undef INSERTION_SORT +#undef COUNTING_SORT + + +#define GRADE_CHR GRADE_UD("⍋","⍒") B GRADE_CAT(c1)(B t, B x) { if (isAtm(x) || RNK(x)==0) thrM(GRADE_CHR": Argument cannot be a unit"); if (RNK(x)>1) x = toCells(x); @@ -82,7 +236,7 @@ B GRADE_CAT(c1)(B t, B x) { return r; } - SLOW1(GRADE_UD("⍋","⍒")"𝕩", x); + SLOW1(GRADE_CHR"𝕩", x); TALLOC(BI32p, tmp, ia); SGetU(x) for (usz i = 0; i < ia; i++) { @@ -94,6 +248,8 @@ B GRADE_CAT(c1)(B t, B x) { TFREE(tmp); decG(x); return r; } + + B GRADE_CAT(c2)(B t, B w, B x) { if (isAtm(w) || RNK(w)==0) thrM(GRADE_CHR": 𝕨 must have rank≥1"); if (isAtm(x)) x = m_atomUnit(x); @@ -128,14 +284,14 @@ B GRADE_CAT(c2)(B t, B w, B x) { usz s = 0, e = wia+1; while (e-s > 1) { usz m = (s+(i64)e)/2; - if (c GRADE_UD(<,>) wi[m-1]) e = m; + if (c LT wi[m-1]) e = m; else s = m; } rp[i] = s; } } else { SGetU(x) - SLOW2("𝕨"GRADE_UD("⍋","⍒")"𝕩", w, x); + SLOW2("𝕨"GRADE_CHR"𝕩", w, x); B* wp = arr_bptr(w); if (wp==NULL) { HArr* a = toHArr(w); @@ -152,7 +308,7 @@ B GRADE_CAT(c2)(B t, B w, B x) { usz s = 0, e = wia+1; while (e-s > 1) { usz m = (s+e) / 2; - if (compare(c, wp[m-1]) GRADE_UD(<,>) 0) e = m; + if (compare(c, wp[m-1]) LT 0) e = m; else s = m; } rp[i] = s; @@ -161,7 +317,12 @@ B GRADE_CAT(c2)(B t, B w, B x) { decG(w);decG(x); return r; } -#undef GRADE_CAT #undef GRADE_CHR + +#undef LT +#undef FOR +#undef PRE_UD +#undef PRE64 +#undef GRADE_CAT #undef GRADE_NEG #undef GRADE_UD diff --git a/src/builtins/sort.c b/src/builtins/sort.c index d270c073..2371dad7 100644 --- a/src/builtins/sort.c +++ b/src/builtins/sort.c @@ -1,21 +1,14 @@ #include "../core.h" #include "../utils/talloc.h" +// Defines Sort, Grade, and Bins + #define CAT0(A,B) A##_##B #define CAT(A,B) CAT0(A,B) typedef struct BI32p { B k; i32 v; } BI32p; typedef struct I32I32p { i32 k; i32 v; } I32I32p; #define GRADE_UD(U,D) U -#define GRADE_NEG -#define GRADE_CHR "⍋" #include "grade.h" #define GRADE_UD(U,D) D -#define GRADE_NEG - -#define GRADE_CHR "⍒" #include "grade.h" - -#define SORT_UD(U,D) U -#include "sort.h" -#define SORT_UD(U,D) D -#include "sort.h" diff --git a/src/builtins/sort.h b/src/builtins/sort.h deleted file mode 100644 index a2195dba..00000000 --- a/src/builtins/sort.h +++ /dev/null @@ -1,156 +0,0 @@ -#define SORT_C1 CAT(SORT_UD(and,or),c1) -#define LT SORT_UD(<,>) - -#define TIM_B SORT_UD(bA,bD) - -#define SORT_CMP(W, X) SORT_UD(compare(W, X), compare(X, W)) -#define SORT_NAME TIM_B -#define SORT_TYPE B -#include "sortTemplate.h" - -#define FOR(I,MAX) SORT_UD(for (usz I=0; I>=), PRE_UD(K,>>,<<=)) - -#define INSERTION_SORT(T) \ - rp[0] = xp[0]; \ - for (usz i=0; ie; j--) { i+=c0[j]; rp[i]--; } \ - ) \ - for (usz i=1; i>8)]++; } \ - PRE_SUM; \ - i16 *r0 = (i16*)(c0+2*256); \ - for (usz i=0; i>8)]++; rp[c]=v; } \ - TFREE(alloc) - RADIX2(u8, - u64 s0=0; u64 s1=0; - FOR(j, 256/8) { PRE64(0); PRE64(1); } - ); - } else if (n < 1<<15) { // Radix sort - RADIX2(u32, - u32 s0=0; u32 s1=0; - FOR(j, 256) { - u32 p0=s0; s0+=c0[j]; c0[j]=p0; - u32 p1=s1; s1+=c1[j]; c1[j]=p1; - } - ); - #undef RADIX2 - } else { - COUNTING_SORT(i16); - } - } else if (xe==el_i32) { - i32* xp = i32any_ptr(x); - i32* rp; r = m_i32arrv(&rp, n); - if (n < 40) { - INSERTION_SORT(i32); - } else if (n <= 256) { - #define RADIX4(T, PRE_SUM) \ - TALLOC(u8, alloc, 4*256*sizeof(T) + n*4); \ - T *c0=(T*)alloc, *c1=c0+256, *c2=c1+256, *c3=c2+256, *c3o=c3+128; \ - for (usz j=0; j<4*256; j++) c0[j]=0; \ - for (usz i=0; i> 8)]++; \ - c2 [(u8)(v>>16)]++; c3o[(i8)(v>>24)]++; } \ - PRE_SUM; \ - i32 *r0 = (i32*)(c0+4*256); \ - for (usz i=0; i> 8)]++; rp[c]=v; } \ - for (usz i=0; i>16)]++; r0[c]=v; } \ - for (usz i=0; i>24)]++; rp[c]=v; } \ - TFREE(alloc) - RADIX4(u8, - u64 s0=0; u64 s1=0; u64 s2=0; u64 s3=0; - FOR(j, 256/8) { PRE64(0); PRE64(1); PRE64(2); PRE64(3); } - ); - } else { - RADIX4(usz, - usz s0=0; usz s1=0; usz s2=0; usz s3=0; - FOR(j, 256) { - u32 p0=s0; s0+=c0[j]; c0[j]=p0; - u32 p1=s1; s1+=c1[j]; c1[j]=p1; - u32 p2=s2; s2+=c2[j]; c2[j]=p2; - u32 p3=s3; s3+=c3[j]; c3[j]=p3; - } - ); - #undef RADIX4 - } - } else { - B xf = getFillQ(x); - HArr_p ra = m_harrUv(n); - SGet(x) - for (usz i = 0; i < n; i++) ra.a[i] = Get(x,i); - CAT(TIM_B,tim_sort)(ra.a, n); - r = withFill(ra.b,xf); - } - decG(x); - return FL_SET(r, CAT(fl,SORT_UD(asc,dsc))); -} -#undef SORT_C1 -#undef LT -#undef TIM_B -#undef FOR -#undef PRE_UD -#undef PRE64 -#undef INSERTION_SORT -#undef COUNTING_SORT -#undef SORT_UD From 36e54d1a392235dcd3737c019b43e10c30ecbf2d Mon Sep 17 00:00:00 2001 From: Marshall Lochbaum Date: Sat, 27 Aug 2022 11:10:17 -0400 Subject: [PATCH 6/8] 1-byte radix with length ==256 doesn't work because of SWAR overflow when descending --- src/builtins/grade.h | 6 +++--- 1 file changed, 3 insertions(+), 3 deletions(-) diff --git a/src/builtins/grade.h b/src/builtins/grade.h index 7d78ea1e..21901e27 100644 --- a/src/builtins/grade.h +++ b/src/builtins/grade.h @@ -69,7 +69,7 @@ B SORT_C1(B t, B x) { i8* rp; r = m_i8arrv(&rp, n); if (n<16) { INSERTION_SORT(i8); - } else if (n<=256) { // Radix/bucket sort + } else if (n<256) { // Radix/bucket sort TALLOC(u8, c0, 256); u8 *c0o=c0+128; // Offset for signedness for (usz j=0; j<256; j++) c0[j]=0; for (usz i=0; i Date: Sat, 27 Aug 2022 15:16:20 -0400 Subject: [PATCH 7/8] Radix grades --- src/builtins/grade.h | 189 ++++++++++++++++++++++++------------------- 1 file changed, 108 insertions(+), 81 deletions(-) diff --git a/src/builtins/grade.h b/src/builtins/grade.h index 21901e27..bcedf0ee 100644 --- a/src/builtins/grade.h +++ b/src/builtins/grade.h @@ -22,12 +22,6 @@ #define FOR(I,MAX) GRADE_UD(for (usz I=0; I>=), PRE_UD(K,>>,<<=)) - #define INSERTION_SORT(T) \ rp[0] = xp[0]; \ for (usz i=0; i>=), PRE_UD(K,>>,<<=)) + +#define CHOOSE_SG_SORT(S,G) S +#define CHOOSE_SG_GRADE(S,G) G + +#define RADIX_SORT_i8(T, TYP) \ + TALLOC(T, c0, 256); T *c0o=c0+128; \ + for (usz j=0; j<256; j++) c0[j]=0; \ + for (usz i=0; i>8)]++; } \ + RADIX_SUM_2_##T; \ + i16 *r0 = (i16*)(c0+2*256); \ + CHOOSE_SG_##TYP( \ + for (usz i=0; i>8)]++]=v; } \ + , \ + I *g0 = (i32*)(r0+n); \ + for (usz i=0; i>8)]++]=g0[i]; } \ + ) \ + TFREE(alloc) +#define RADIX_SUM_2_u8 u64 s0=0, s1=0; FOR(j,256/8) { PRE64(0); PRE64(1); } +#define RADIX_SUM_2(T) T s0=0, s1=0; FOR(j,256) { PRE(T,0); PRE(T,1); } +#define RADIX_SUM_2_usz RADIX_SUM_2(usz) +#define RADIX_SUM_2_u32 RADIX_SUM_2(u32) + +#define RADIX_SORT_i32(T, TYP, I) \ + TALLOC(u8, alloc, 4*256*sizeof(T) + n*(4 + CHOOSE_SG_##TYP(0,4+sizeof(I)))); \ + T *c0=(T*)alloc, *c1=c0+256, *c2=c1+256, *c3=c2+256, *c3o=c3+128; \ + for (usz j=0; j<4*256; j++) c0[j]=0; \ + for (usz i=0; i> 8)]++; \ + c2 [(u8)(v>>16)]++; c3o[(i8)(v>>24)]++; } \ + RADIX_SUM_4_##T; \ + i32 *r0 = (i32*)(c0+4*256); \ + CHOOSE_SG_##TYP( \ + for (usz i=0; i> 8)]++; rp[c]=v; } \ + for (usz i=0; i>16)]++; r0[c]=v; } \ + for (usz i=0; i>24)]++; rp[c]=v; } \ + , \ + i32 *r1 = r0+n; I *g0 = (i32*)(r1+n); \ + for (usz i=0; i> 8)]++; r1[c]=v; rp[c]=g0[i]; } \ + for (usz i=0; i>16)]++; r0[c]=v; g0[c]=rp[i]; } \ + for (usz i=0; i>24)]++; rp[c]=g0[i]; } \ + ) \ + TFREE(alloc) +#define RADIX_SUM_4_u8 u64 s0=0, s1=0, s2=0, s3=0; FOR(j, 256/8) { PRE64(0); PRE64(1); PRE64(2); PRE64(3); } +#define RADIX_SUM_4(T) T s0=0, s1=0, s2=0, s3=0; FOR(j, 256) { PRE(u32,0); PRE(u32,1); PRE(u32,2); PRE(u32,3); } +#define RADIX_SUM_4_usz RADIX_SUM_4(usz) +#define RADIX_SUM_4_u32 RADIX_SUM_4(u32) + #define SORT_C1 CAT(GRADE_UD(and,or),c1) B SORT_C1(B t, B x) { if (isAtm(x) || RNK(x)==0) thrM(GRADE_UD("∧","∨")": Argument cannot have rank 0"); @@ -69,13 +134,8 @@ B SORT_C1(B t, B x) { i8* rp; r = m_i8arrv(&rp, n); if (n<16) { INSERTION_SORT(i8); - } else if (n<256) { // Radix/bucket sort - TALLOC(u8, c0, 256); u8 *c0o=c0+128; // Offset for signedness - for (usz j=0; j<256; j++) c0[j]=0; - for (usz i=0; i>8)]++; } \ - PRE_SUM; \ - i16 *r0 = (i16*)(c0+2*256); \ - for (usz i=0; i>8)]++; rp[c]=v; } \ - TFREE(alloc) - RADIX2(u8, - u64 s0=0; u64 s1=0; - FOR(j, 256/8) { PRE64(0); PRE64(1); } - ); - } else if (n < 1<<15) { // Radix sort - RADIX2(u32, - u32 s0=0; u32 s1=0; - FOR(j, 256) { - u32 p0=s0; s0+=c0[j]; c0[j]=p0; - u32 p1=s1; s1+=c1[j]; c1[j]=p1; - } - ); - #undef RADIX2 + } else if (n < 256) { + RADIX_SORT_i16(u8, SORT,); + } else if (n < 1<<15) { + RADIX_SORT_i16(u32, SORT,); } else { COUNTING_SORT(i16); } @@ -117,35 +157,9 @@ B SORT_C1(B t, B x) { if (n < 40) { INSERTION_SORT(i32); } else if (n < 256) { - #define RADIX4(T, PRE_SUM) \ - TALLOC(u8, alloc, 4*256*sizeof(T) + n*4); \ - T *c0=(T*)alloc, *c1=c0+256, *c2=c1+256, *c3=c2+256, *c3o=c3+128; \ - for (usz j=0; j<4*256; j++) c0[j]=0; \ - for (usz i=0; i> 8)]++; \ - c2 [(u8)(v>>16)]++; c3o[(i8)(v>>24)]++; } \ - PRE_SUM; \ - i32 *r0 = (i32*)(c0+4*256); \ - for (usz i=0; i> 8)]++; rp[c]=v; } \ - for (usz i=0; i>16)]++; r0[c]=v; } \ - for (usz i=0; i>24)]++; rp[c]=v; } \ - TFREE(alloc) - RADIX4(u8, - u64 s0=0; u64 s1=0; u64 s2=0; u64 s3=0; - FOR(j, 256/8) { PRE64(0); PRE64(1); PRE64(2); PRE64(3); } - ); + RADIX_SORT_i32(u8, SORT,); } else { - RADIX4(usz, - usz s0=0; usz s1=0; usz s2=0; usz s3=0; - FOR(j, 256) { - u32 p0=s0; s0+=c0[j]; c0[j]=p0; - u32 p1=s1; s1+=c1[j]; c1[j]=p1; - u32 p2=s2; s2+=c2[j]; c2[j]=p2; - u32 p3=s3; s3+=c3[j]; c3[j]=p3; - } - ); - #undef RADIX4 + RADIX_SORT_i32(u32, SORT,); } } else { B xf = getFillQ(x); @@ -183,22 +197,14 @@ B GRADE_CAT(c1)(B t, B x) { else rp[r0++] = i; } decG(x); return r; - } else if (xe==el_i8) { - i8* xp = i8any_ptr(x); - i32 min=-128, range=256; - TALLOC(usz, tmp, range+1); - for (i64 i = 0; i < range+1; i++) tmp[i] = 0; - GRADE_UD( // i8 range-based - for (usz i = 0; i < ia; i++) (tmp-min+1)[xp[i]]++; - for (i64 i = 1; i < range; i++) tmp[i]+= tmp[i-1]; - for (usz i = 0; i < ia; i++) rp[(tmp-min)[xp[i]]++] = i; - , - for (usz i = 0; i < ia; i++) (tmp-min)[xp[i]]++; - for (i64 i = range-2; i >= 0; i--) tmp[i]+= tmp[i+1]; - for (usz i = 0; i < ia; i++) rp[(tmp-min+1)[xp[i]]++] = i; - ) - TFREE(tmp); decG(x); - return r; + } else if (xe==el_i8 && ia>8) { + i8* xp = i8any_ptr(x); usz n=ia; + RADIX_SORT_i8(usz, GRADE); + decG(x); return r; + } else if (xe==el_i16 && ia>16) { + i16* xp = i16any_ptr(x); usz n = ia; + RADIX_SORT_i16(usz, GRADE, i32); + decG(x); return r; } if (xe==el_i32 || xe==el_c32) { // safe to use the same comparison for i32 & c32 as c32 is 0≤x≤1114111 i32* xp = tyany_ptr(x); @@ -224,6 +230,11 @@ B GRADE_CAT(c1)(B t, B x) { TFREE(tmp); decG(x); return r; } + if (ia > 40) { + usz n=ia; + RADIX_SORT_i32(usz, GRADE, i32); + decG(x); return r; + } TALLOC(I32I32p, tmp, ia); for (usz i = 0; i < ia; i++) { @@ -321,8 +332,24 @@ B GRADE_CAT(c2)(B t, B w, B x) { #undef LT #undef FOR +#undef PRE #undef PRE_UD #undef PRE64 +#undef CHOOSE_SG_SORT +#undef CHOOSE_SG_GRADE +#undef RADIX_SORT_i8 +#undef RADIX_SUM_1_u8 +#undef RADIX_SUM_1_usz +#undef RADIX_SORT_i16 +#undef RADIX_SUM_2_u8 +#undef RADIX_SUM_2 +#undef RADIX_SUM_2_usz +#undef RADIX_SUM_2_u32 +#undef RADIX_SORT_i32 +#undef RADIX_SUM_4_u8 +#undef RADIX_SUM_4 +#undef RADIX_SUM_4_usz +#undef RADIX_SUM_4_u32 #undef GRADE_CAT #undef GRADE_NEG #undef GRADE_UD From 9376e02ecb25834f9590d7b405837319f22f744a Mon Sep 17 00:00:00 2001 From: Marshall Lochbaum Date: Sat, 27 Aug 2022 15:44:25 -0400 Subject: [PATCH 8/8] Boolean sorting --- src/builtins/grade.h | 13 ++++++++++++- 1 file changed, 12 insertions(+), 1 deletion(-) diff --git a/src/builtins/grade.h b/src/builtins/grade.h index bcedf0ee..698e1fb6 100644 --- a/src/builtins/grade.h +++ b/src/builtins/grade.h @@ -129,7 +129,18 @@ B SORT_C1(B t, B x) { if (n <= 1) return x; u8 xe = TI(x,elType); B r; - if (xe==el_i8) { + if (xe==el_bit) { + u64* xp = bitarr_ptr(x); + u64* rp; r = m_bitarrv(&rp, n); + usz sum = bit_sum(xp, n); + u64 n0 = GRADE_UD(n-sum, sum); + u64 ones = -1ull; + u64 v0 = GRADE_UD(0, ones); + usz i=0, e=(n+63)/64; + for (; i