Base64 encoding with SIMD instructions

Base64 decoding could also be vectorized, although the speedup is not very impressive, merely 35%. Read more ...

Base64 encoding with SIMD instructions

An SSE code is more than 2 times faster on Core i7, and around 70% faster on Core i5. Read more...

Fast conversion of floating-point values to string

The conversion to string could be 15 times faster than sprintf. Read more...

Base64 encoding — implementation study

Although base64 encoding is a very basic algorithm, it could be sped up a little (25% sounds good?) Read more...

Boolean function for the rescue

The problem is defined as follows: a set of features is saved using bit-sets (usually large), and there is a list/map/whatever of sets containing features of different objects. We have to find which features are unique. Read more...

Fast exact summation using small and large superaccumulators

Interesting article by Radford M. Neal:

I present two new methods for exactly summing a set of floating-point numbers, and then correctly rounding to the nearest floating-point number. Higher accuracy than simple summation (rounding after each addition) is important in many applications, such as finding the sample mean of data.

Optimizing Dijkstra for real-world performance

Another interesting paper:

Our experimental results currently put our prototype implementation at about twice as fast as the Boost implementation of the algorithm on both real-world and generated large graphs. Furthermore, this preliminary implementation was written in only a few weeks, by a single programmer. The fact that such an early prototype compares favorably against Boost, a well-known open source library developed by expert programmers, gives us reason to believe our design for the queue is indeed better suited to the problem at hand, and the favorable time measurements are not a product of any specific implementation technique we employed.

Conversion numbers to binary ASCII representation - new method

Recently I've checked different methods to convert numbers to binary representation, including use of new PDEP instruction from BMI2 extension.

Today I've updated the article with new SWAR version 2, a tricky use of multiplication. The method is not faster, but I like the approach---in certain conditions multiplication can be seen as multi-shift/bit-or instruction. I've already use multiplication in this way to emulate instruction pmovmskb.

Speeding up bit-parallel population count

Nearly 50% faster than naive version for large data sets. Discovered by accident. :)

SIMD-ized searching in unique constant dictionary

The problem: there is a ordered dictionary containing only unique keys. Dictionary is read only, and keys are 32-bit (SSE) or 64-bit (AVX2). Read more

SSE: Generating mask where n leading (trailing) bytes are set

Informal specification:

__m128i mask_lower(const unsigned n) {
    
    assert(n < 16);
    switch (n) {
        case 0: return {0x00, 0x00, 0x00, 0x00, 0x00, 0x00, 0x00, 0x00, 0x00, 0x00, 0x00, 0x00, 0x00, 0x00, 0x00, 0x00};
        case 1: return {0xff, 0x00, 0x00, 0x00, 0x00, 0x00, 0x00, 0x00, 0x00, 0x00, 0x00, 0x00, 0x00, 0x00, 0x00, 0x00};
        case 2: return {0xff, 0xff, 0x00, 0x00, 0x00, 0x00, 0x00, 0x00, 0x00, 0x00, 0x00, 0x00, 0x00, 0x00, 0x00, 0x00};
        // ...
        case 14: return {0xff, 0xff, 0xff, 0xff, 0xff, 0xff, 0xff, 0xff, 0xff, 0xff, 0xff, 0xff, 0xff, 0xff, 0xff, 0x00};
        case 15: return {0xff, 0xff, 0xff, 0xff, 0xff, 0xff, 0xff, 0xff, 0xff, 0xff, 0xff, 0xff, 0xff, 0xff, 0xff, 0xff};
    }
}

__m128i mask_higher(const unsigned n) {
    
    assert(n < 16);
    return ~mask_lower(15 - n);
}

Read more ...

Conversion int and double to string - comparison

Milo Yip has compared different itoa and dtoa implementations on Core i7, including my itoa algorithm 2, that use SSE2 instructions.

Results for itoa are interesting: SSE2 version is not as good as it seemed to be. Tricky branchlut algorithm is only 10% slower, moreover is perfectly portable. One obvious drawback of this method is using lookup-table - in real environment where is a big pressure on cache, memory access could be a bottleneck.

Intepolation search revisited

Interpolation search revisited -- theoretical complexity of interpolation search is very promising, but...

Conversion number to hexadecimal representation

Conversion numbers to hexadecimal representation - SWAR, plain SSE, and draft of BMI2 implementation.

Article SSSE3: printing hex values describes the same topic but is limited to exploit PSHUFB.

Conversion numbers to binary representation

New article Conversion numbers to binary representation — SIMD & SWAR versions.

Few years ago I've described an MMX variant of SIMD algorithm. But the text was in Polish, so audience was limited.

Scalar version of SSE move mask instruction

SSE instruction PMOVMSKB gathers all most significant bits from bytes and stores them as a single 16-bit value; similar action is performed by MOVMSKPD and MOVMSKPS.

Such operation could be easily done using scalar multiplication. Read more ...

Asymmetric numeral systems

Very interesting idea: Asymmetric numeral systems: entropy coding combining speed of Huffman coding with compression rate of arithmetic coding.

SIMD-friendly Rabin-Karp modification

Rabin-Karp algorithm uses a weak hash function to locate possible substring positions. This modification uses merely equality of first and last char of searched substring, however equality of chars can be done very fast in parallel, even without SIMD instruction. Read more ...

Integer log 10 of unsigned integer - SIMD version

Fast calculate ceil(log10(x)) of some unsigned number is described on Bit Twiddling Hacks, this text show the SIMD solution for 32-bit numbers.

Algorithm:

1. populate value in XMM registers. Since maximum value of this function is 10 we need three registers.

movd   %eax, %xmm0          // xmm0 = packed_dword(0, 0, 0, x)
pshufd $0, %xmm0, %xmm0 \n" // xmm0 = packed_dword(x, x, x, x)
movapd %xmm0, %xmm1
movapd %xmm0, %xmm2

2. compare these numbers with sequence of powers of 10.

// powers_a = packed_dword(10^1 - 1, 10^2 - 1, 10^3 - 1, 10^4 - 1)
// powers_c = packed_dword(10^5 - 1, 10^6 - 1, 10^7 - 1, 10^8 - 1)
// powers_c = packed_dword(10^9 - 1, 0, 0, 0)
pcmpgtd powers_a, %xmm0
pcmpgtd powers_b, %xmm1
pcmpgtd powers_c, %xmm2

result of comparisons are: 0 (false) or -1 (true), for example:

xmm0 = packed_dword(-1, -1, -1, -1)
xmm1 = packed_dword( 0, 0, -1, -1)
xmm2 = packed_dword( 0, 0, 0, 0)

3. calculate sum of all dwords

psrld $31, %xmm0       // xmm0 = packed_dword( 1, 1, 1, 1) - convert -1 to 1
psubd %xmm1, %xmm0     // xmm0 = packed_dword( 1, 1, 2, 2)
psubd %xmm2, %xmm0     // xmm0 = packed_dword( 1, 1, 2, 2)

// convert packed_dword to packed_word
pxor %xmm1, %xmm1
packssdw %xmm1, %xmm0 // xmm0 = packed_word(0, 0, 0, 0, 1, 1, 2, 2)

// max value of word in xmm0 is 3, so higher
// bytes are always zero
psadbw %xmm1, %xmm0   // xmm0 = packded_qword(0, 6)

4. save result, i.e. the lowest dword

movd %xmm0, %eax      // eax = 6

Sample program is available.

Mask for zero/non-zero bytes

The description of Determine if a word has a zero byte from "Bit Twiddling Hacks" says about haszero(v): "the result is the high bits set where the bytes in v were zero".

Unfortunately this is not true. High bits are also set for ones followed zeros, i.e. haszero(0xff010100) = 0x00808080. Of course the result is still valid (non-zero if there were any zero byte), but if we want to iterate over all zeros or find last zero index, this could be a problem.

It's possible to create exact mask:

uint32_t nonzeromask(const uint32_t v) {
    // MSB are set if any of 7 lowest bits are set
    const uint32_t nonzero_7bit = (v & 0x7f7f7f7f) + 0x7f7f7f7f;

    return (v | nonzero_7bit) & 0x80808080;
}

uint32_t zeromask(const uint32_t v) {
     // negate MSBs
     return (nonzeromask(v)) ^ 0x80808080;
}

Function nonzeromask requires 4 simple instructions, and zeromask one additional xor.