Compiler Explorer

Source code

#include <stdio.h>
#include <stdint.h>
#define NDEBUG
#include <assert.h>

#ifdef __SSE4_1__
// TODO: MSVC detection of SSE4
#include <immintrin.h>
uint32_t fastpack2_sse4(uint32_t x)
{
        __m128i v = _mm_cvtsi32_si128(x);
        // ppoonnmm...ddccbbaa
        v = _mm_unpacklo_epi8(v, _mm_slli_epi16(v,4));      // movaps copy + psllw + punpcklbw
        // ppoonnmm | nnmmllkk | ... | ddccbbaa | bbaa0000    // in low 64 bits of vector
        v = _mm_cvtepi8_epi16(v);    // pmovsxbw
        // pppppppp | ppoonnmm | ... | dddddddd | ddccbbaa | bbbbbbbb | bbaa0000   // filling 128-bit reg
        v = _mm_slli_epi16(v, 2);    // 0xbba0
        // pppppppp | oonnmm00 | ... | dddddddd | ccbbaa00 | bbbbbbbb | aa000000
        return _mm_movemask_epi8(v); // grab the top bit of each byte.
}

// scalar is better
uint32_t fastpack4_sse4(uint32_t x)
{
        // even if x is being loaded from memory, pmovsxbw would load 64 bits from a memory source
        // which would give us high garbage even if it was known safe (non-faulting and not causing an expensive page-split)
        // uint32_t x = *px;
        __m128i v = _mm_cvtsi32_si128(x);  // movd load or from a GPR
        v = _mm_cvtepi8_epi16(v);    // pmovsxbw : widen each 0xba element to 0xbbba
        v = _mm_slli_epi16(v, 4);    // 0xbba0
        return _mm_movemask_epi8(v); // grab the top bit of each byte.
}

uint32_t fastpack8_sse2(uint32_t x)
{
        __m128i v = _mm_cvtsi32_si128(x);  // movd, or a memory source operand for pmovmskb
        return _mm_movemask_epi8(v);
}
#endif

uint32_t pack_medium_n(uint32_t x, unsigned n)
{
    // at most 4 groups, 2 pairs
    assert(n>=4);   // works(?) but is less efficient for n>10
    // n==4..10

// low pair starting at bit #n-1 (counting from 0)
    uint32_t lowpair = x>>(n-1);
    lowpair &= 3;
    //  high pair starting at bit #3*n - 1
    uint32_t highpair = x>>(3*n - 1 + 2);
    uint32_t highmask = n <= 10 ? 3<<2 : 1<<2;
    highpair &= highmask;      // 1<<2 for n=9 or 10.  FIXME, needs to be 0 for n==11 and higher if there's high garbage.
    return lowpair | highpair;
}

uint32_t pack_large_n(uint32_t x, unsigned n)
{
    assert(n>=11 && n<=32);
    // higher n would be shift-count UB, lower would have a third group not part of the low pair.

// n=11 to 16: one pair : return x>>(n-1) & 3
    // n=17 and up: one bit : return x&1 or x>>(n-1) & 1;
    // if high padding (above the highest group) is known 0, then x>>(n-1) & 3 can be used unconditionally
#if 1
    x >>= n-1;
    uint32_t mask = n>16 ? 1 : 3;  // get a mask ready in parallel with shifting x.  GCC uses SBB/AND/ADD with a possible false-dependency on EAX on uarches other than AMD

// nope, this is backwards and off-by-1, would give n>=16 : 3 : 1;
//    uint32_t mask = (n>>3) | 1;  // n is 1..32.  For n=16..31, the 1<<4 = 16 bit is set, and isn't for smaller n.  And all higher bits are clear.
//                                 // for n=32, n>>4 is 2 not 1, but  x>>31  has no set bits above the LSB
    return x & mask;               // both dependency chains needed only as inputs to AND
#else
// if specializing for a known or invariant n, avoids a shift for n>16
    uint32_t onegroup = x&1;
    uint32_t twogroups = x>>(n-1) & 3;
    return n > 16 ? onegroup : twogroups;
#endif
}

// unsigned long n avoids the need to zero-extend to pointer width, if it didn't inline.
uint32_t pack_LUT(uint32_t x, unsigned long n)
{
    struct LUT {
        // group both arrays into one static variable so compilers don't waste instructions getting the base address separately for each.
        uint8_t shift[32-4];  // counts first so the offset to next array is small, fitting in a disp8
        struct {uint32_t mask, mult;} maskmult[]; // [32-4]
    };
    static const struct LUT lut = {
      .maskmult = {
        {0b00011000000110000001100000011000, ((1uL<<18) | (1u<<12) | (1u<<6) | 1)<<3 }, // n=4 ; 8 groups = 4 pairs
        {0b00000011000000001100000000110000, ((1uL<<16) | (1u<<8)  | 1)<<6 },  // n=5 ; 6 groups = 3 pairs
        {0b00100000000001100000000001100000, ((1uL<<20) | (1u<<10) | 1)<<2 },  // n=6 ; 5 groups = 3ish pairs
        {0b00000000001100000000000011000000, ((1uL<<12) | 1)<<0 }, // n=7 ; 4 groups, 2 pairs // no high garbage even without shifting to the top, shift=18
        {0b00000001100000000000000110000000, ((1uL<<14) | 1)<<0 }, // n=8 ; 4 groups, 2 pairs
        {0b00000100000000000000001100000000, ((1uL<<16) | 1)<<0 }, // n=9 ; 3 groups, 2ish pairs
        {0b00100000000000000000011000000000, ((1uL<<18) | 1)<<0 }, // n=10 ; 3 groups, 2ish pairs
        {0b00000000000000000000110000000000, 1}, // n=11   one pair, don't need to multiply anymore if we want to dispatch to a third version.
        {0b00000000000000000001100000000000, 1}, // n=12
        {3<<12, 1}, // n=13
        {3<<13, 1}, // n=14
        // TODO: finish this, or use a simpler strategy for larger n to save table size at the cost of code size and more instructions for runtime-variable n
      },
      .shift = {
        32-8, 32-6, 32-5,  // n = 4..6 shift from the top to clear high garbage
        18, 21, 24, 27,    // n=7..10 have a second bit-pair (or single bit) at the top
        // n=11 and up just have 2 or fewer groups, trivial shift and keep only 2 or 1 bit
        10, 11, 12, 13, 14, 16, 17, 18, 19, 20, 21, 22, 23, 24, 25, 26, 27, 28, 29, 30, 31
        // could choose multipliers so these are all 30
      }
    };
    x &= lut.maskmult[n-4].mask;
    x *= lut.maskmult[n-4].mult;
    x >>= lut.shift[n-4];
    return x;
}

// untested
uint32_t fastpack2(uint32_t x)
{
        // 0po00nm00lk00ji00hg00fe00dc00ba0   AND mask to keep bit-pairs
    x &= 0b01100110011001100110011001100110;
        // 0ponm0000lkji0000hgfe0000dcba000    x += x<<2 /*LEA*/;  x &= mask;
    x += x<<2;
    x &= 0b01111000011110000111100001111000;
        // 0ponmlkji00000000hgfedcba0000000    x += x<<4;  x &= mask2;
    x += x<<4;
    x &= 0b01111111100000000111111110000000;
        // ponmlkjihgfedcba0000000000000000    x *= (1<<1) | (1<<9)
    x *= (1u<<1)|(1u<<9);
    x >>= 16;
    return x;
}

#define fastpack3 fastpack3_new

// same output as orig for these test cases
uint32_t fastpack3_new(uint32_t x) {
        // 0000j00i00h00g00f00e00d00c00b00a // incoming x
    x &= 0b00001100001100001100001100001100;
        //     ji    hg    fe    dc    ba
    x += x<<4;      // x *= 17
        // ji  jihg  hgfe  fedc  dcba  ba
    x &= 0b00001111000000001111000000001111;
        //     jihg        fedc        ba
        //         fedc        ba   // partial product from 1<<8
        // fedc        ba           // partial product from 1<<16
        // ((1<<0) | (1<<8) | (1<<16))<<4
    x *= (1<<4) | (1<<12) | (1<<20);
    // with the multiplier left-shifted by 4, result goes to the top, no high garbage after shifting
    x >>= 22;
    // x &= 0x3FF;  // 10-bit.  
    return x;
}

// without shifting the multiplier left by 4, leaving high garbage
uint32_t fastpack3_new_v1(uint32_t x) {
        // 0000j00i00h00g00f00e00d00c00b00a // incoming x
    x &= 0b00001100001100001100001100001100;
        //     ji    hg    fe    dc    ba
    x += x<<4;      // x *= 17
        // ji  jihg  hgfe  fedc  dcba  ba
    x &= 0b00001111000000001111000000001111;
        //     jihg        fedc        ba
        //         fedc             // partial product from 1<<8
        // fedc        ba           // partial product from 1<<16
    x *= (1<<0) | (1<<8) | (1<<16);
    x >>= 18;
    x &= 0x3FF;  // 10-bit
    return x;
}

uint32_t fastpack4(uint32_t x) {
        // hhhhggggffffeeeeddddccccbbbbaaaa
    x &= 0b00011000000110000001100000011000;
        // 000hg000000fe000000dc000000ba000  // result, and 1<<0 partial product
        // 00000fe000000dc000000ba000000000  // 1<<6 partial product
        // 0000000dc000000ba000000000000000  // 1<<12
        // 0dc000000ba000000000000000000000  // 1<<18
        // with an extra shift by 3 to put the result at the top, avoiding high garbage
    x *= (1uL<<21) | (1uL<<15) | (1u<<9) | (1u<<3);
    x >>= (32-8);
    return x;
}

uint32_t fastpack5(uint32_t x) {
        // xxfffffeeeeedddddcccccbbbbbaaaaa
    x &= 0b00000011000000001100000000110000;
        // 000000fe00000000dc00000000ba00000  // result, and 1<<0 partial product
        // 00000000dc00000000ba00000          // 1<<8 partial product
        // dc00000000ba00000                  // 1<<16
    // with an extra <<6 to clear high garbage:
    x *= (1uL<<22) | (1u<<14) | (1u<<6);
    x >>= (32-6);
    return x;
}

uint32_t fastpack8_bitpairs(uint32_t x) {
    
    x &= 0b00000001100000000000000110000000;  // this can be an AARch64 immediate
        // 0000000dc00000000000000ba0000000
        // 000000000ba0000000..............  // x<<14
    x += x<<14;
    x >>= 21;  // (32-4) - 7
    return x;
/*
    x = (x >> 7) | (x>>16);  // more ILP but more x86 instructions, worse throughput if latency isn't a bottleneck
    x &= 0xFF;
    return x;
*/
}   
    
uint32_t fastpack8(uint32_t x) {
    
    // 1. Mul
    // 0000000d0000000c0000000b0000000a
    x &= 0b1000000010000000100000001;    // using an evenly-spaced mask so AArch64 can do it as an immediate, not bit-pairs
    // d0000000c0000000b0000000a         // <<7
    // 0c0000000b0000000a                // <<14
    // 00b0000000a                       // <<21
    // -------++++---------------------
    // dcb0000dcba0000cba00000ba000000a
    
    x *= 0b1000000100000010000001 << 7;  // extra shift of the multiplier to put the result at the top
    x >>= (32-4);
    return x;
}

uint32_t pack_dispatcher(uint32_t x, unsigned n)
{
#ifdef HAVE_PEXT
        return pack_pext(x, n);
#endif
    switch (n){
        default: return pack_large_n(x, n);  // alternative: make the LUTs larger so switch is simpler.

case 10:
        case 9:
        case 8:
        case 7:
        case 6:
        case 5:
        case 4:
            // return pack_medium_n(x, n);  // no LUT just shifts
            return pack_LUT(x, n);

case 3:   return fastpack3(x);
        case 2:
#ifdef __SSE4_1__
            return fastpack2_sse4(x);
#else
            return fastpack2(x);
#endif
        case 1:   return x;  // trivial

}
}

// Konstantin's version
uint32_t fastpack3_orig(uint32_t x) {

// 0000j00i00h00g00f00e00d00c00b00a

// 1. Mul
    // 0000j00000000000f00000000000b000
    // 00000000f00000000000b000
    // f00000000000b000
    // ----++++++++++------------------
    // f000j000f000b000f000b0000000b000

uint32_t x_0 = x & 0b1000000000001000000000001000;
    x_0 *= 0b10000000100000001;

// 2. Mul
    // 00000i00000000000e00000000000a--
    // 000000000e00000000000a
    // 0e00000000000a
    // ----++++++++++------------------
    // 0e000i000e000a000e000a0000000a00

uint32_t x_1 = x & 0b1000000000001000000000001;
    x_1 *= 0b1000000010000000100;

// 3. Mul
    // 000000h00000000000d000000000----
    // 0000000000d000000000
    // ----++++++++++------------------
    // 000000h000d0000000d0000000000000

uint32_t x_2 = x & 0b1000000000001000000000;
    x_2 *= 0b1000000010000;

// 4. Mul
    // 0000000g00000000000c000000------
    // 00000000000c000000
    // ----++++++++++------------------
    // 0000000g000c0000000c000000000000

uint32_t x_3 = x & 0b1000000000001000000;
    x_3 *= 0b100000001000000;

// 5. Or
    // f000j000f000b000f000b0000000b000
    // 0e000i000e000a000e000a0000000a00
    // 000000h000d0000000d0000000000000
    // 0000000g000c0000000c000000000000
    // ----++++++++++------------------
    // fe00jihgfedcba00fedcba000000ba00

x = x_0 | x_1 | x_2 | x_3;
    x >>= 18;
    x &= 0x3FF;
    return x;
}

uint32_t fastpack4_orig(uint32_t x) {

// 000h000g000f000e000d000c000b000a

// 1. Mul
    // 000h0000000f0000000d0000000b0000
    // 00000f0000000d0000000b0000
    // 0000000d0000000b0000
    // 0d0000000b0000
    // ---++++++++---------------------
    // 0d0h0f0d0b0f0d0b000d0b00000b0000

uint32_t x_odd = x & 0b10000000100000001000000010000;
    x_odd *= 0b1000001000001000001;

// 2. Mul
    // 0000g0000000e0000000c0000000a---
    // 000000e0000000c0000000a
    // e0000000c0000000a
    // 00c0000000a
    // ---++++++++---------------------
    // e0c0g0e0c0a0e0c0a000c0a00000a000

uint32_t x_even = x & 0b1000000010000000100000001;
    x_even *= 0b1000001000001000001000;

// 3. Or
    // 0d0h0f0d0b0f0d0b000d0b00000b0000
    // e0c0g0e0c0a0e0c0a000c0a00000a000
    // ---++++++++---------------------
    // edchgfedcbafedcba00dcba0000ba000

x = x_even | x_odd;
    x >>= 21;
    x &= 0xFF;  // (movzx eax, al)
    return x;
}

uint32_t fastpack5_orig(uint32_t x) {

// 110000f0000e0000d0000c0000b0000a

// 1. Mul
    // 00000000000e000000000c000000000a
    // 000e000000000c000000000a
    // 00000c000000000a
    // 0000000a
    // --++++++------------------------
    // 000e0c0a000e0c0a00000c0a0000000a

uint32_t x_even = x & 0b100000000010000000001;
    x_even *= 0b1000000010000000100000001;

// 2. Mul
    // 00f000000000d000000000b00000----
    // 0000d000000000b00000
    // 000000b00000
    // --++++++------------------------
    // 00f0d0b00000d0b0000000b000000000

uint32_t x_odd = x & 0b10000000001000000000100000;
    x_odd *= 0b100000001000000010000;

// 3. Or
    // 000e0c0a000e0c0a00000c0a0000000a
    // 00f0d0b00000d0b0000000b000000000
    // --++++++------------------------
    // 00fedcba000edcba00000cba0000000a

x = x_even | x_odd;
    x >>= 24;
    return x;
}

uint32_t fastpack8_orig(uint32_t x) {

// 1. Mul
    // 0000000d0000000c0000000b0000000a
    // d0000000c0000000b0000000a
    // 0c0000000b0000000a
    // 00b0000000a
    // -------++++---------------------
    // dcb0000dcba0000cba00000ba000000a

x &= 0b1000000010000000100000001;
    x *= 0b1000000100000010000001;
    x >>= 21;
    x &= 0xFF;  // (movzx eax, al)
    return x;
}

//#define fastpack4 fastpack4_sse4
//#define fastpack8 fastpack8_sse2

int main(void) {
    uint32_t x;

printf("pack3 in 6 or 13 operations, or 5 on AArch64 with ORR reg, reg, reg, lsl #4:\n");
    x = 0xC0000007; printf("%08X -> %08X\n", x, fastpack3(x));
    x = 0xC0000038; printf("%08X -> %08X\n", x, fastpack3(x));
    x = 0xC000003F; printf("%08X -> %08X\n", x, fastpack3(x));
    x = 0xC00001C0; printf("%08X -> %08X\n", x, fastpack3(x));
    x = 0xC71C71C7; printf("%08X -> %08X\n", x, fastpack3(x));
    x = 0xF8E38E38; printf("%08X -> %08X orig\n", x, fastpack3_orig(x));
    x = 0xF8E38E38; printf("%08X -> %08X new\n", x, fastpack3(x));
     x = 0xFFFFFFFF; printf("%08X -> %08X\n", x, fastpack3(x));

printf("\npack4 in 7 operations:\n");
    x = 0x0000000F; printf("%08X -> %08X\n", x, fastpack4(x));
    x = 0x000000F0; printf("%08X -> %08X\n", x, fastpack4(x));
    x = 0x000000FF; printf("%08X -> %08X\n", x, fastpack4(x));
    x = 0x00000F00; printf("%08X -> %08X\n", x, fastpack4(x));
    x = 0x0F0F0F0F; printf("%08X -> %08X\n", x, fastpack4(x));
    x = 0xF0F0F0F0; printf("%08X -> %08X\n", x, fastpack4(x));
    x = 0xFFFFFFFF; printf("%08X -> %08X\n", x, fastpack4(x));

printf("\npack5 in 6 operations:\n");
    x = 0xC000001F; printf("%08X -> %08X\n", x, fastpack5(x));
    x = 0xC00003E0; printf("%08X -> %08X\n", x, fastpack5(x));
    x = 0xC00003FF; printf("%08X -> %08X\n", x, fastpack5(x));
    x = 0xC0007C00; printf("%08X -> %08X\n", x, fastpack5(x));
    x = 0xC1F07C1F; printf("%08X -> %08X\n", x, fastpack5(x));
    x = 0xFE0F83E0; printf("%08X -> %08X\n", x, fastpack5(x));
    x = 0xFFFFFFFF; printf("%08X -> %08X\n", x, fastpack5(x));

printf("\npack8 in 4 operations:\n");
    x = 0x000000FF; printf("%08X -> %08X\n", x, fastpack8(x));
    x = 0x0000FF00; printf("%08X -> %08X\n", x, fastpack8(x));
    x = 0x0000FFFF; printf("%08X -> %08X\n", x, fastpack8(x));
    x = 0x00FF0000; printf("%08X -> %08X\n", x, fastpack8(x));
    x = 0x00FF00FF; printf("%08X -> %08X\n", x, fastpack8(x));
    x = 0xFF00FF00; printf("%08X -> %08X\n", x, fastpack8(x));
    x = 0xFFFFFFFF; printf("%08X -> %08X\n", x, fastpack8(x));
}