Your code is fast – if you're lucky

But it's not exactly a cosmetic change. x++ is semantically different from x; x++; I wonder if clang would make it branchless if you instead write

  if (BLQS_CMP(x, piv)) { *lwr = x; ++lwr; }
    else { *rwr = x; --rwr; }

The difference is post-increment has strange semantics. While the compiler should be able to understand that the value wasn't used and post increment and pre increment are the same I wouldn't be surprised if it tracks that it was post increment and misses some optimizations because it's trying to garuntee post increment semantics.

Although it's true compilers can be very sensitive to exact phrasing triggering specific optimization passes. So it still might not give the branchless version by changing it to pre increment (which is the same as a normal +=1).

The only way to really know is to dig into what optimization passes clang took in both cases and analyze the difference.

Coincidentally a few days ago, also for a sorting algorithm, I stumbled over a situation where replacing the branchless cmov with branching instructions actually made the code 30% faster: https://github.com/graphhopper/graphhopper/pull/3380 Or at least the AI explained it this way to me, but I'm unsure if this is correct.

Does anyone know exactly what is going on here to cause this difference? I am extremely puzzled that the "beginner friendly" code is not at some point in the compilation pipeline in EXACTLY the same representation as the non-"beginner friendly" code. I would imagine they'd be in the same form very early on, perhaps even at the point of generating an initial syntax tree. And once they take on the same form in the compilation pipeline, the resulting compiler output should be identical. So what is really going on here?

Is it only me..?

Quicksort is supposed to be an algorithm that has O(n) to O(n²) performance and O(n log n) being only an average performance case. Test was made on random data coming from different archs (so I doubt it's characteristic would be remotely identical).

Given input size of 50M it means that performance could be between 50M (5e7) up to 2.5e15. That's like performance instability of 8 orders of magnitude.

I'm not sure here if we can't write instead that "Your code is fast if you picked fast case for it" especially since fix of 6 OOM is smaller than algorithm's performance range.

    else *rwr-- = x;

No. Make that obvious and the PR can pass. Argue, and you're off the project.

What if you wrote this in a branchless way?

bool v = BLQS_CMP(x, piv);

int* ptr = v ? lwr : rwr;

*ptr = x;

ptr += int(v) * 2 - 1;

I really envy programmers who are so skilled at this kind of low-level optimization.

The same meaning, but different performance based on notation—it's ultimately about entering LLVM's optimization pass, which likely comes down to differences in the internal IR pattern. It almost feels like a difference in innate talent...

I feel like I can build CRUD applications well enough, but I still seem to be weak at low-level processing.

Where can I learn these kinds of techniques? I'd appreciate any book recommendations.

Auto-vectorization is not a programming model!

Don't forget to disable all "spectre and friends" mitigations in your linux kernel, and some workloads will become much faster.

Can you do the same on the windows kernel or apple kernels?

My code is not fast. Writing efficient code takes a lot of brain power. My brain is of the lazy type - it wants the computer (but not AI) to solve things. I only write code so I can be lazier lateron.

I think with this approach, we will only win if a language allows for:

1) ease of writing, and 2) fastness

Right now languages don't really combine both. We have ease of writing e. g. ruby or python, but they are slower than C, our godfather language. So far all languages that try to solve both problems, become mega-verbose and tend to gravitate more towards one than the other - usually e. g. "let's write a replacement for C". I wonder if combining both 1) and 2) is possible, kind of like select on your own what to combine, so if my time is precious, I write a quick prototype. If this must become faster, I write it with more details. That's still not really a language that combines 1) and 2) genuinely but perhaps it is an acceptable trade-off. Right now we kind of mix two languages here, say, ruby+java or python+C or any other similar combination.

Wow that is quite surprising. Almost seems like it could be a compiler bug tbh. Very fragile optimisation if not!

I think with this approach, we will only win if a language allows for:

1) ease of writing, and 2) fastness

Auto-vectorization is not a programming model!

What if you wrote this in a branchless way?

bool v = BLQS_CMP(x, piv);

int* ptr = v ? lwr : rwr;

*ptr = x;

ptr += int(v) * 2 - 1;

The ternary isn't guaranteed to be branchless. In your case it should almost always be on a modern compiler, but it really shouldn't be present in what one would call branchless code.

I really envy programmers who are so skilled at this kind of low-level optimization.

I feel like I can build CRUD applications well enough, but I still seem to be weak at low-level processing.

Where can I learn these kinds of techniques? I'd appreciate any book recommendations.

    else *rwr-- = x;

No. Make that obvious and the PR can pass. Argue, and you're off the project.

*foo++ (and --) is an extremely common C idiom. I'd argue it's clearer than the separated version.

Fine. I'll use inline asm then.

Don't forget to disable all "spectre and friends" mitigations in your linux kernel, and some workloads will become much faster.

Can you do the same on the windows kernel or apple kernels?

Yes you can do the same on Windows.

https://www.grc.com/inspectre.htm

Only really on CPUs with the vulnerability

Interesting, did not know you could do that but it does make sense. Make sure it’s air gapped.

Representing code in a compiler is not precisely trivial, and the two statements are actually quite different from a compiler or AST perspective. Just looking at the first branch:

  *lwr = x; lwr++;

This could be be represented with something like this (and this is a very vague approximation of an AST):

  block
    statement (assignment)
      expression
        operator (dereference)
          variable
      expression
        variable
    statement
      expression
        operator (post-increment)
          variable

The second form, that looks this in the source:

  *lwr++ = x

might look like this in the AST:

  block
    statement (assignment)
      expression
        operator (post-increment)
          operator (dereference)
            variable
      expression
        variable

If these two ASTs are to take the same form in the compilation pipeline, there needs to be some kind of pass that transforms one into the other.

As for why the second form is faster (or rather, why the compiler can generate faster code): There is likely an optimisation pass somewhere in llvm that recognises a pattern that this fits into, which allows it to generate branchless instructions. For instance, in the second form, there is a pattern of:

  operator (post-increment)
    operator (dereference)

that might be recognised by a pass. In the former form, the two operators are far apart in the tree, so a pass would have to "look further" to match them up. A single pass likely won't do this, either for (compiler) performance reasons, or for correctness reasons.

Finding which pass that is can be non-trivial, as it's more than a matter of enabling individual passes until one works. It might be that an earlier pass does some code reshaping that allows the relevant pass to work. My suggestion would be to dump the llvm ir at the end, and find the rough pattern that you're looking for, then re-run the compilation with `-mllvm -print-after-all` to see what the IR looks like after each pass, and then manually "look back" until you can't see the pattern any more.

Have you considered what happens in the presence of multiple threads? I know aarch is weakly memory ordered but the assembly output by those two versions is quite different when multiple threads are involved. There must be a reason spreading the mutations across multiple lines causes the compiler to pessimize to the branch version.

Edit: not saying they are different just that it is harder for the compiler to see the safety of the transform when they are quite different.

Yeah, that's an interesting question. I would expect the two to look the same once converted to SSA.

You can ask the compiler to dump its intermediate representations, so it might not be too hard to answer the question.

Complete speculation, but it could be simply that clang only considers applying the branchless optimisation when the code inside the "if" block "looks like" a single instruction/statement -- and this "looks like" check might be done quite early, before "*x++ = y;" and "*x = y; x++;" are converted to the same thing.

Why restrict the optimisation unnecessarily like this? It's never the intention to artificially restrict optimisations, but they are difficult to test and debug as they can interact with other optimisations, and applying an optimisation too broadly leads to horrible correctness bugs. So if you cannot be certain that you understand all possible interactions now and in the future, it makes some sense to be conservative in applying them.

But it's not exactly a cosmetic change. x++ is semantically different from x; x++; I wonder if clang would make it branchless if you instead write

  if (BLQS_CMP(x, piv)) { *lwr = x; ++lwr; }
    else { *rwr = x; --rwr; }

The only way to really know is to dig into what optimization passes clang took in both cases and analyze the difference.

Looking at the optimization pipeline in Compiler Explorer (https://godbolt.org/z/rd3qber3b) after the code is converted to SSA (SROAPass), it seems that the fast version computes the pointer increment (`lwr++`) before storing `x` (`*lwr = x`), whereas the slow version does it the other way around. This happens because Clang chooses a different instruction order when generating a single statement compared to two statements.

The computation is the same, but apparently, this tiny change prevents SimplifyCFGPass from turning the code into the branchless version later on. I'm not sure why this happens, perhaps because it messes something up in the pattern recognition of the pass?

Note: in some cases (not completely sure if it applies in this situation as I haven't looked into the details of this code enough), using pre-increment AND when the value is used is a dependency on that operation, so "+= 1" is actually faster for out of order execution.

So yeah, as you say, compilers can be very sensitive to this, and in the past (although more than 10 years ago now, so might not be relevant now), ICC often used to generate better (faster executing) code when using "val += 1" vs "++val".

Is it only me..?

Given input size of 50M it means that performance could be between 50M (5e7) up to 2.5e15. That's like performance instability of 8 orders of magnitude.

I'm not sure here if we can't write instead that "Your code is fast if you picked fast case for it" especially since fix of 6 OOM is smaller than algorithm's performance range.

You're talking about the complexity of the Quicksort algorithm, whereas the article is about code generation.

Both versions sort the same data using the same algorithm. Just a tiny change in the source code caused Clang to generate different machine code.

Using different seed values - (srand(1), srand(2), srand(time(NULL))) essentially leads to the same result. With a good choice of pivot, Quicksort is very close to O(n log n) in practice, so that’s not the key factor here.

The interesting thing is that the generated machine code changes significantly.

I’m assuming he measured time by averaging on 100s of instances, or he maintained the exact same input for both versions of code. Would be a big oversight if not!

Looking at big-O isn't very informative. We have plenty of statistical tools for telling whether there is an effect even with noisy data.

I hope the test data is the same when comparing the different runs. So the big o notation should be the same across different runs.

if you have decent (randomized) pivoting, you never hit the worst case or anything like it

Wow that is quite surprising. Almost seems like it could be a compiler bug tbh. Very fragile optimisation if not!

A large part of optimization is understanding the hardware architecture in detail and then making your software architecture mirror that hardware architecture as closely as possible. A compiler can't do this for you. Much of "performance engineering" is applying your understanding of how the hardware components work and are connected to the software design.

Most software introduces a large number of unnecessary stalls where one part of the hardware is waiting on or bottlenecked by a different piece of hardware. Optimization is often about removing or minimizing these stalls to the extent possible. But to do this you need to understand why specific code choices cause the hardware to stall in the first place. The most basic form of this is CPU cache locality but there are many levels.

Compilers are great at localized micro-optimizations, except for SIMD. Most idiomatic code can be detected and transformed by the compiler into something that takes advantage of the hardware design. You don't need to do tricky bit-twiddling stuff, the compiler is better at it than you are in most cases (this wasn't always true). I always recommend people interested in optimization experiment running snippets of code through the compiler using godbolt with optimization turned on. You can learn a lot about how the compiler sees and understands your code by studying this. It is also a good way to became familiar with assembly.

Leveraging SIMD and vector ISAs is a mess. Different hardware architectures rely on different idioms and compilers cannot auto-vectorize most code that can be vectorized. You need to learn the idioms of each vector ISA and how to write them using intrinsics. A great way to learn this is reading other peoples' SIMD code. Modern SIMD code is wide in that you need to memorize a lot of things but conceptually pretty shallow. It mostly comes down to learning the idiomatic tricks and gadgets for an ISA -- code is composed from these conceptual primitives.

You could read compiler books, but I would actually recommend reading about CPUs and computer architecture directly. If you understand how the hardware works, then the optimizations are all very natural and fit into the picture perfectly, instead of being some arcane compiler magic that you have to take as a disconnected fact.

Personally I actually haven't read too many books on optimizations, I just absorbed information over years one thing at a time, but something like Computer Organization and Design is a pretty good intro to the low-level picture. If you want to drown in extremely dense technical topics that will give you a lot of jumping off points to search, read Agner Fog's microarchicture optimization guide (https://www.agner.org/optimize/). It won't tell you what LLVM is doing, but it'll tell you why it's doing it. Fair warning, it's dense and pretty dry.

Then it depends how interested you are in doing low-level nonsense. If you spend a lot of time writing performance oriented systems code, you'll come to use profiling tools that show you the assembly. If you stare at it long enough, you sometimes start to question why the compiler wrote it this way. And you're naturally led as you try to optimize your code to wonder how LLVM is coming up with this ASM that it spits out and why it sometimes gets it wrong.

There's nothing magical or that requires innate talent. You can learn all of this very naturally if you work close to the metal and take the time to question how the abstraction layer below you actually works. If you keep doing this, you eventually find out it's not that deep, it's just a lot of stuff accumulated over time, but none of it particularly difficult or inaccessible.

I'm not sure it's a "technique" but the general insight worth taking away from this is that compiler authors often write optimizers to recognize specific patterns so writing your code in a more idiomatic form increases the odds an optimizer will be able to optimize it.

In this specific instance, at the hardware level it helps to understand how the branch predictor works and why quicksort in particular is essentially the worst case for the branch predictor, and then you'll understand why the cmov/csel optimization is such a big win.

Casey Muratori's https://computerenhance.com is a good resource

I don't have a book, but my best recommendation would be: Learn how to measure. Whatever language you are using should have a profiler. Learn how to use it and look at how your code holds up. Look for surprises; those are optimization opportunities and learning opportunities. Make a change that you think has a reasonable shot at optimizing performance, look at how the profile and benchmarks change. If you're using a compiled language, look at what the generated code looks like. Is it what you'd expect, or did the compiler do something that blindsided you? Can you get rid of that overhead? Can you use a faster algorithm? Can you make things faster without making it an unmaintainable mess?

In general, you must expect to try ten things and then one of them will help you; fewer if your ideas are bad. :-) Occasionally, you learn new things (either about your machine or your language or your code base) and then you can try that elsewhere in the code (but don't go overboard, every technique has its limits).

DO NOT FALL PREY TO SUPERSTITION. Always measure in one way or the other. Don't do stuff blindly just because someone on the Internet told you (there's a _lot_ of bad performance advice out there).

First important step is the need to make your application faster.

Second step is to actually make it faster.

And I would say 99% of the time you can make applications faster without going down so deep on this low-level (which is still not easy) but sometimes you can only get faster via low-level optimizations.

In order of priority, I’d say:

1. If you write CRUD apps, make sure that the database does the heavy lifting.

2. Take note of algorithm complexities, use hashtables as appropriate, and write good hash functions when you start using hashtables/dictionaries.

3. Avoid pointer-heavy datastructures. In high level languages like Java, an object reference is a pointer dereference that can stall the CPU waiting for memory. This is sometimes optimized, but you can’t depend on it. The true zealots call this “data oriented design”.

4. If you write C/C++,rust, or the like, you might want to learn to read assembly. Godbolt.com is a fun way to learn. Note that not all instructions are equally fast: Long division and trigonomic functions are slower than integer adds, even when they are both a single instruction.

5. The next level is probably going for vectorized instructions: SIMD (ARM Neon, AVX). The most original applications can be found at lemire.me: a professor exploring optimizing things like JSON parsing using the latest processor features.

I personally learned a lot from cpp conference videos. I would highly recommmend it.

Got to my first job out of college and they gave me core dumps and had me debug the kernel for a year. Not my kind of fun, but definitely got me skilled in the art of low level.

I did Nand2Tetris and then I understod why I need to vectorize. You get a view from nand-gate to software and get to see all the interfaces. Its a wonderful course.

it's arguably more of a cpu bug than a computer bug. The problem is that predictable data determined whether a cmov or a branch is faster. cmov is only faster than if when the branch is unpredictable. Summer the compiler doesn't know what values your program will be called with, it can only pick and hope. To fix this, cpus could have an instruction like cmov but that learns whether speculation would be profitable and converts to treat it like a branch if better.

I did Nand2Tetris and then I understod why I need to vectorize. You get a view from nand-gate to software and get to see all the interfaces. Its a wonderful course.

In order of priority, I’d say:

1. If you write CRUD apps, make sure that the database does the heavy lifting.

2. Take note of algorithm complexities, use hashtables as appropriate, and write good hash functions when you start using hashtables/dictionaries.

Got to my first job out of college and they gave me core dumps and had me debug the kernel for a year. Not my kind of fun, but definitely got me skilled in the art of low level.

You can, in fact, give compiler a hint wether the branch will be predictable or not. E.g. clang has __builtin_unpredictable().

tbh, I'm not sure if we need "more clever" CPUs since they cause more and more unpredictable performance outcome -- and more pressing issues like Spectre. Besides, looking at the evolution of technology tells us that performance increases mainly comes from parallelism, like multicore and GPUs. Which doesn't apply to "annoyingly not parallel" problems of course. But besides "gambling" approaches like branch prediction, there's not much you can do about them in the general case.

The ternary isn't guaranteed to be branchless. In your case it should almost always be on a modern compiler, but it really shouldn't be present in what one would call branchless code.

You could rewrite it to

    int* ptr = ((v != 0) * lwr) + ((v == 0) * rwr);

I doubt it would be faster than the ternary though.

Does any standard library have a glsl-like select(condition, true_val, false_val) that is more or less guaranteed to correspond to a csel?

I also agree that computer architecture is more important - it grounds your understanding of how to write efficient code regardless of platform since most machines today share very similar ideas (OOO execution, caches, NUMA etc).

How ever, I will disagree slightly that all the optimizations compilers do are about optimizing for a given architecture; some transformations are just weird algorithmic black magic about optimizing the underlying code itself. Knowing how to make sure the compiler sees through a given construct to give you the low level expression you want is too much art and randomness; we need better ways to express optimization expectations so that if the compiler fails to match expectations it becomes a loud compiler error.

Thee optimization in question here was not obvious at all. It's a bs clang codegen quirk

Casey Muratori's https://computerenhance.com is a good resource

DO NOT FALL PREY TO SUPERSTITION. Always measure in one way or the other. Don't do stuff blindly just because someone on the Internet told you (there's a _lot_ of bad performance advice out there).

These are concepts I've come across in books I've read, and your explanation is easy to understand. But I can tell my understanding is lacking in some parts because I haven't actually tried it in practice. As you said, I think using Godbolt to disassemble things myself is the right approach.

Thanks for the advice.

About that kind of 'technique', I guess I should make it a habit to dig into the compiler, which is still a black box to me. I should study a few techniques myself. Have a good day

First important step is the need to make your application faster.

Second step is to actually make it faster.

I personally learned a lot from cpp conference videos. I would highly recommmend it.

I should study that. Thank you

Fine. I'll use inline asm then.

I sometimes wonder whether that's a better idea for these micro-optimizations, rather than looking at the assembly code and trying to coax the compiler into generating what you want.

That said, I'm not keen on "argue and you're off the project" work environment.

That would unironically be easier to understand than doing *rwr-- = x.

*foo++ (and --) is an extremely common C idiom. I'd argue it's clearer than the separated version.

As I experienced while trying to write out an AST for this pattern, the operator precedence makes it harder to read. I would at least prefer that it's written as *(foo++).

Only really on CPUs with the vulnerability

You can, in fact, give compiler a hint wether the branch will be predictable or not. E.g. clang has __builtin_unpredictable().

You could rewrite it to

    int* ptr = ((v != 0) * lwr) + ((v == 0) * rwr);

I doubt it would be faster than the ternary though.

Yeah, that's an interesting question. I would expect the two to look the same once converted to SSA.

You can ask the compiler to dump its intermediate representations, so it might not be too hard to answer the question.

Thee optimization in question here was not obvious at all. It's a bs clang codegen quirk

Does any standard library have a glsl-like select(condition, true_val, false_val) that is more or less guaranteed to correspond to a csel?

I hope the test data is the same when comparing the different runs. So the big o notation should be the same across different runs.

Yes you can do the same on Windows.

https://www.grc.com/inspectre.htm

Thanks for the advice.

As I experienced while trying to write out an AST for this pattern, the operator precedence makes it harder to read. I would at least prefer that it's written as *(foo++).

I sometimes wonder whether that's a better idea for these micro-optimizations, rather than looking at the assembly code and trying to coax the compiler into generating what you want.

That said, I'm not keen on "argue and you're off the project" work environment.

I should study that. Thank you

Interesting, did not know you could do that but it does make sense. Make sure it’s air gapped.

I’m assuming he measured time by averaging on 100s of instances, or he maintained the exact same input for both versions of code. Would be a big oversight if not!

Both versions use the same input data. I also tried different random initial values and got essentially the same result. I didn't test hundreds of inputs, since that would have been mostly a waste of time in this case. The algorithm and the data distribution remain practically the same. What I'm measuring is the machine code that Clang generates for the hot loop.

You're talking about the complexity of the Quicksort algorithm, whereas the article is about code generation.

Both versions sort the same data using the same algorithm. Just a tiny change in the source code caused Clang to generate different machine code.

The interesting thing is that the generated machine code changes significantly.

Look at the source, Luke!

It's just rand() in test.c

Right. You won't learn from a computer architecture book or a uarch guide about SSA form, or LICM, or other famous compiler principle like the central role of inlining decisions ("the mother of all optimizations"). I don't have a good resource to recommend here, this is where my lack of formal training bites. "Go read a hundred blog articles by compiler experts" doesn't feel like very useful advice.

>Knowing how to make sure the compiler sees through a given construct to give you the low level expression you want is too much art and randomness; we need better ways to express optimization expectations so that if the compiler fails to match expectations it becomes a loud compiler error.

There's a parallel with hardware there. Verilog is a kind of hardware language designed for an abstract simulator, in the same way than C is designed for a standard abstract machine for the sake of portability. You end up with an idea of the assembly/RTL you want the compiler/synthetizer to generate in your head, and then it's a game of writing the right pattern that will be recognized and generate the output you want.

I think this is partially unavoidable, because we're inherently asking the compiler to generate a non-portable target-specific output in what is supposed to be a portable high-level language. If you start injecting compiler hints or requirements in your "portable" code, it all becomes a bit of a mess. Part of the problem is also that the high-level languages we're using were designed at a time were many questions were still unsettled. Things like signed integers being two's complements is a recent change in C and C++. But I think some of it is intrinsic impedance mismatch between high-level code and machine code.

I'm not sure I would like a proliferation of annotations that direct exactly how the compiler should optimize (like "must use cmov/csel here"), because if internal optimizer choices become public API, people will rely on internals in their large legacy codebases. I expect this would be a force that ossifies the compiler and prevent optimizations from improving. The "register" and "inline" keywords in C used to mean something to the compiler. But they were misused, having them be a requirement would have held back performance more than anything.

Then again I accepted the same justification against Postgres planner hints, and now that the idea has been recast as a plan stability feature I'm actually very happy with that idea. I'm uncomfortable with letting old calcified codebases hold back compiler internal, but at the same time once you find a way to have the compiler generate what you want, there's a real need to not have it break silently when you upgrade.

Edit: not saying they are different just that it is harder for the compiler to see the safety of the transform when they are quite different.

The compiler doesn’t have to consider how other threads will see this. If you don’t use atomics, the compiler is free to reorder things however it wishes as long as the single threaded behavior is the same in the end.

Representing code in a compiler is not precisely trivial, and the two statements are actually quite different from a compiler or AST perspective. Just looking at the first branch:

  *lwr = x; lwr++;

This could be be represented with something like this (and this is a very vague approximation of an AST):

  block
    statement (assignment)
      expression
        operator (dereference)
          variable
      expression
        variable
    statement
      expression
        operator (post-increment)
          variable

The second form, that looks this in the source:

  *lwr++ = x

might look like this in the AST:

  block
    statement (assignment)
      expression
        operator (post-increment)
          operator (dereference)
            variable
      expression
        variable

If these two ASTs are to take the same form in the compilation pipeline, there needs to be some kind of pass that transforms one into the other.

  operator (post-increment)
    operator (dereference)

  block
    statement (assignment)
      expression
        operator (post-increment)
          operator (dereference)
            variable
      expression
        variable

I don't think this could possibly be a valid AST for '*lwr++ = x' because the increment is not an operation on the dereferenced value, its an operation on the pointer. So in this case I don't see how it could help but be transformed into a form similar to the "beginner friendly" case.

Or perhaps I am wrong and it would generate an AST like you describe and rely on later passes to actually create a proper dependency graph. My mental model of how these kinds of postfix operators work always assumed it must very early on turn it into two separate statements. Thank you for the suggestions.

About that kind of 'technique', I guess I should make it a habit to dig into the compiler, which is still a black box to me. I should study a few techniques myself. Have a good day

Looking at big-O isn't very informative. We have plenty of statistical tools for telling whether there is an effect even with noisy data.

While the article is interesting, it's the case of a missed optimization by the compiler and (possibly) something that a future version of the compiler will catch.

This sort of optimization is one that'd I'd not spend too much time trying to fix or catch, unless you are doing it for fun or you have a specific piece of code in a hot path that needs to go fast.

Where I'd spend time if I were trying to write very fast code which a compiler is unlikely to get right is SIMD optimizations. Specifically with floating point values.

One thing compilers can't and won't do is reorganize floating point optimizations (well, unless you explicitly give them permission to do that). That means the way you write your floating point code can really nerf performance and exclude you from much faster assembly.

The sure fire way to actually make such code faster is learning and using SIMD expressions. The compilers can sometimes get this right, but it's quite fragile.

I ran 10k test locally on 2e5 and I'm seeing 4 orders of magnitude instability, but very high local stability (i.e. runs within specific second are very stable, showing almost no deviation, runs couple second later are the same, but results are within 1 OoM of the prior results (smaller batches, 500 tests).

I'm not saying that optimization isn't valid, what I'm saying is that Quicksort shouldn't be optimized over randomized per run data set.

It's a good rule of thumb that can be quite useful without additional analysis. It's not always the right way to do performance tuning, but I can't count the number of times I've changed an O(n^3) to an O(n) and seen massive performance gains as a result.

if you have decent (randomized) pivoting, you never hit the worst case or anything like it

You don't need randomized pivoting for this, there are deterministic ones like median of median that will also result in a O(nlogn) worst case.

Also note that with a randomized pivoting you _might_ hit a O(n^2) worst case, it's just that it's incredibly rare and cannot be forced by an attacker controlling your input, so for most practical purposes can be ignored.

That would unironically be easier to understand than doing *rwr-- = x.

Less about calcifying the specific optimization and more like “this loop is expected to be vectorized”. That way if someone changes something subtle that prevents it from being vectorized it’s a compiler error. Striking that balance and how to express those constraints in a flexible way is the hard bit.

As you say register and inline were wrong, but we have force inline and force inline so clearly the pendulum swung back a little bit because the compiler completely ignoring is also not good. We have ways to force the compiler to do an unconditional move because source level heuristics are completely incorrect for making such a decision. The die is already cast, we just keep living with a shitty status quo instead of something a bit more robust.

Look at the source, Luke!

It's just rand() in test.c

Therefore it is srand(1).

  block
    statement (assignment)
      expression
        operator (post-increment)
          operator (dereference)
            variable
      expression
        variable

Oops, you're right - I got the operators the wrong way around! I forget the precise precedence of * and ++ in C sometimes. Assuming that it would be bracketed as *(lwr++), it should actually be:

  block
    statement (assignment)
      expression
        operator (dereference)
          operator (post-increment)
            variable
      expression
        variable

I'm not saying that optimization isn't valid, what I'm saying is that Quicksort shouldn't be optimized over randomized per run data set.

While the article is interesting, it's the case of a missed optimization by the compiler and (possibly) something that a future version of the compiler will catch.

This sort of optimization is one that'd I'd not spend too much time trying to fix or catch, unless you are doing it for fun or you have a specific piece of code in a hot path that needs to go fast.

Where I'd spend time if I were trying to write very fast code which a compiler is unlikely to get right is SIMD optimizations. Specifically with floating point values.

The sure fire way to actually make such code faster is learning and using SIMD expressions. The compilers can sometimes get this right, but it's quite fragile.

You don't happen to be running `srand(time(0))`?

> I'm not saying that optimization isn't valid, what I'm saying is that Quicksort shouldn't be optimized over randomized per run data set.

But yeah, this is correct. When optimizing, we want to pin every variable other than the change, as much as possible.

Yeah it helps make awful code decent, and some algorithms are better than others, but in terms of high performance code, locality, vectorization, and branching often matter much more big O.

Thank you for giving me the keywords.

You don't need randomized pivoting for this, there are deterministic ones like median of median that will also result in a O(nlogn) worst case.

median of median does give you guarenteed n*log(n) but it doubles your memory reads per pass making it pretty poor. single random is almost guarenteed to take fewer passes (and median of 3-7 random values can make the number of extra passes over the minimum to be very low)

Therefore it is srand(1).

Oops, you're right - I got the operators the wrong way around! I forget the precise precedence of * and ++ in C sometimes. Assuming that it would be bracketed as *(lwr++), it should actually be:

  block
    statement (assignment)
      expression
        operator (dereference)
          operator (post-increment)
            variable
      expression
        variable

Thank you for giving me the keywords.

You don't happen to be running `srand(time(0))`?

> I'm not saying that optimization isn't valid, what I'm saying is that Quicksort shouldn't be optimized over randomized per run data set.

But yeah, this is correct. When optimizing, we want to pin every variable other than the change, as much as possible.

Then you run the risk of overfitting to whatever seed you picked. I think fixing a seed is probably a good idea most of the time but it’s worth considering

Yeah it helps make awful code decent, and some algorithms are better than others, but in terms of high performance code, locality, vectorization, and branching often matter much more big O.

That depends on what end ends. For a small N, very bad algorithms can still be plenty fast. Sometimes it can even be faster. Some of your lower O algorithms can have very terrible constant factors, which means they're terrible when N is small, but as N gets large.

Big and small N different for different algorithms and hardware both.

Totally agree, but I'm in Java land and cursed to have about the worst case scenario for those things :D

Then you run the risk of overfitting to whatever seed you picked. I think fixing a seed is probably a good idea most of the time but it’s worth considering

Big and small N different for different algorithms and hardware both.

Totally agree, but I'm in Java land and cursed to have about the worst case scenario for those things :D

Oh yeah I'm also mostly in Java land. These concerns are also concerns in Java, plus whatever C2 gets up to.

Valhalla can't come soon enough.

I've definitely had to change things into SOA in order to eke out performance. My coworkers aren't thrilled seeing `double[] x; double[] y;` but that really is about the only way to get the JVM to play nice.

Valhalla can't come soon enough.

Your code is fast - if you’re lucky

Recently, while I was working on an optimized Quicksort implementation, I came across a rather interesting quirk. Modern compilers (especially Clang) optimize loops using fast, branch-free instructions - provided you use the right programming style.

sort.h - a quicksort with sorting networks

// SPDX-License-Identifier: MIT
// sort.h - Branchless Quicksort
// (c) christof.kaser@gmail.com

#ifndef SORT_H
#define SORT_H

#ifndef BLQS_CMP
#define BLQS_CMP(a, b) ((a) < (b))
#endif

#include <stddef.h>
#include <string.h>
#define min(a, b) (((a) < (b)) ? (a) : (b))

#define SMALLPART 1024
#define SWSZ 512
#define UNROLL 16

#define sort2(a, b) do {  \
    unsigned m = BLQS_CMP(a, b); \
    BLQS_TYPE x = a; \
    a = m ? a : b;  \
    b = m ? b : x;  \
} while(0)

#define sort3(a, b, c) do { \
    sort2(a, b); sort2(b, c); sort2(a, b); \
} while(0)

#define sort4(a, b, c, d) do { \
    sort2(a, b); sort2(c, d); sort2(a, c); \
    sort2(b, d); sort2(b, c); \
} while(0)

#define sort5(a, b, c, d, e) do { \
    sort2(b, c); sort2(d, e); sort2(b, d); \
    sort2(a, c); sort2(a, d); sort2(c, e); \
    sort2(a, b); sort2(c, d); sort2(b, c); \
} while(0)

#define sort6(a, b, c, d, e, f) do { \
    sort2(a, b); sort2(c, d); sort2(e, f); \
    sort2(a, c); sort2(b, d); sort2(e, f); \
    sort2(a, e); sort2(b, f); sort2(c, e); \
    sort2(d, f); sort2(b, c); sort2(d, e); \
    sort2(c, d); \
} while(0)

#define sort7(a, b, c, d, e, f, g) do { \
    sort2(a, b); sort2(c, d); sort2(a, c); \
    sort2(b, d); sort2(b, c); sort2(e, f); \
    sort2(e, g); sort2(f, g); sort2(a, e); \
    sort2(b, f); sort2(c, g); sort2(b, e); \
    sort2(d, g); sort2(c, e); sort2(b, c); \
    sort2(d, f); sort2(d, e); \
} while(0)

#define sort8(a,b,c,d,e,f,g,h) do { \
    sort2(a,b); sort2(c,d); sort2(e,f); sort2(g,h); \
    sort2(a,c); sort2(b,d); sort2(e,g); sort2(f,h); \
    sort2(b,c); sort2(f,g); \
    sort2(a,e); sort2(b,f); sort2(c,g); sort2(d,h); \
    sort2(c,e); sort2(d,f); \
    sort2(b,c); sort2(d,e); sort2(f,g); \
} while (0)

#define sort9(a,b,c,d,e,f,g,h,i) do { \
    sort2(a,d); sort2(b,h); sort2(c,f); sort2(e,i); \
    sort2(a,h); sort2(c,e); sort2(d,i); sort2(f,g); \
    sort2(a,c); sort2(b,d); sort2(e,f); sort2(h,i); \
    sort2(b,e); sort2(d,g); sort2(f,h); \
    sort2(a,b); sort2(c,e); sort2(d,f); sort2(g,i); \
    sort2(c,d); sort2(e,f); sort2(g,h); \
    sort2(b,c); sort2(d,e); sort2(f,g); \
} while (0)

#define sort10(a,b,c,d,e,f,g,h,i,j) do { \
    sort2(a,i); sort2(b,j); sort2(c,h); sort2(d,f); sort2(e,g); \
    sort2(a,c); sort2(b,e); sort2(f,i); sort2(h,j); \
    sort2(a,d); sort2(c,e); sort2(f,h); sort2(g,j); \
    sort2(a,b); sort2(d,g); sort2(i,j); \
    sort2(b,f); sort2(c,d); sort2(e,i); sort2(g,h); \
    sort2(b,c); sort2(d,f); sort2(e,g); sort2(h,i); \
    sort2(c,d); sort2(e,f); sort2(g,h); \
    sort2(d,e); sort2(f,g); \
} while (0)

#define sort11(a,b,c,d,e,f,g,h,i,j,k) do { \
    sort2(a,j); sort2(b,g); sort2(c,e); sort2(d,h); sort2(f,i); \
    sort2(a,b); sort2(d,f); sort2(e,k); sort2(g,j); sort2(h,i); \
    sort2(b,d); sort2(c,f); sort2(e,h); sort2(i,k); \
    sort2(a,e); sort2(b,c); sort2(d,h); sort2(f,j); sort2(g,i); \
    sort2(a,b); sort2(c,g); sort2(e,f); sort2(h,i); sort2(j,k); \
    sort2(c,e); sort2(d,g); sort2(f,h); sort2(i,j); \
    sort2(b,c); sort2(d,e); sort2(f,g); sort2(h,i); \
    sort2(c,d); sort2(e,f); sort2(g,h); \
} while (0)

#define sort12(a,b,c,d,e,f,g,h,i,j,k,l) do { \
    sort2(a,i); sort2(b,h); sort2(c,g); sort2(d,l); sort2(e,k); sort2(f,j); \
    sort2(a,c); sort2(b,e); sort2(d,f); sort2(g,i); sort2(h,k); sort2(j,l); \
    sort2(a,b); sort2(c,j); sort2(e,h); sort2(f,g); sort2(k,l); \
    sort2(b,d); sort2(c,h); sort2(e,j); sort2(i,k); \
    sort2(a,b); sort2(c,d); sort2(e,f); sort2(g,h); sort2(i,j); sort2(k,l); \
    sort2(b,c); sort2(d,f); sort2(g,i); sort2(j,k); \
    sort2(c,e); sort2(d,g); sort2(f,i); sort2(h,j); \
    sort2(b,c); sort2(d,e); sort2(f,g); sort2(h,i); sort2(j,k); \
} while (0)

static void sorting_network(BLQS_TYPE* l, int partszm1_min_1) {
    switch (partszm1_min_1) {
    case 0: break;
    case 1: sort2(l[0],l[1]); break;
    case 2: sort3(l[0],l[1],l[2]); break;
    case 3: sort4(l[0],l[1],l[2],l[3]); break;
    case 4: sort5(l[0],l[1],l[2],l[3],l[4]); break;
    case 5: sort6(l[0],l[1],l[2],l[3],l[4],l[5]); break;
    case 6: sort7(l[0],l[1],l[2],l[3],l[4],l[5],l[6]); break;
    case 7: sort8(l[0],l[1],l[2],l[3],l[4],l[5],l[6],l[7]); break;
    case 8: sort9(l[0],l[1],l[2],l[3],l[4],l[5],l[6],l[7],l[8]); break;
    case 9: sort10(l[0],l[1],l[2],l[3],l[4],l[5],l[6],l[7],l[8],l[9]); break;
    case 10: sort11(l[0],l[1],l[2],l[3],l[4],l[5],l[6],l[7],l[8],l[9],l[10]); break;
    case 11: sort12(l[0],l[1],l[2],l[3],l[4],l[5],l[6],l[7],l[8],l[9],l[10],l[11]); break;
    }
}
#define med5(a,b,c,d,e) do { \
    sort2(a,b); sort2(c,d); sort2(a,c); \
    sort2(b,d); sort2(b,c); sort2(c,e); \
    sort2(b,c); \
} while(0)


static BLQS_TYPE* partition_small(BLQS_TYPE* left, BLQS_TYPE* right) {

    BLQS_TYPE* outerleft = left;
    BLQS_TYPE* pivp = left + 6;

    BLQS_TYPE piv = *pivp;

    BLQS_TYPE l1 = left[1],l2 = left[2];
    BLQS_TYPE r1 = right[-1], r0 = *right;
    med5(l1, l2, piv, r1, r0);
    left[1] = l1; left[2] = l2;
    right[-1] = r1; *right = r0;
    left += 3;
    right -= 2;

    *pivp = *outerleft;

    BLQS_TYPE swbuf[SMALLPART];
    BLQS_TYPE* sw = swbuf;
    BLQS_TYPE* lwr = left;

    while (right - left >= UNROLL) for (int i = UNROLL; i--;) {
        BLQS_TYPE x = *left++;
        if (BLQS_CMP(x, piv)) { *lwr = x; lwr++; }
        else { *sw = x; sw++; }
    }
    while (left <= right) {
        BLQS_TYPE x = *left++;
        if (BLQS_CMP(x, piv)) { *lwr = x; lwr++; }
        else { *sw = x; sw++; }
    }
    memcpy(lwr, swbuf, (sw - swbuf) * sizeof(BLQS_TYPE));
    lwr -= 1;
    *outerleft = *lwr;
    *lwr = piv;
    return lwr;
}

static BLQS_TYPE* partition(BLQS_TYPE* left, BLQS_TYPE* right) {
    BLQS_TYPE* outerleft = left;
    BLQS_TYPE* pivp = left + (right - left) / 2;

    BLQS_TYPE piv = *pivp;

    med5(left[1],left[2],left[3],left[4],left[5]);
    med5(left[11],left[12],left[13],left[14],left[15]);
    med5(pivp[-2], pivp[-1], piv, pivp[1], pivp[2]);
    med5(right[-14], right[-13], right[-12], right[-11], right[-10]);
    med5(right[-4], right[-3], right[-2], right[-1], right[0]);
    med5(left[3], left[13], piv, right[-12], right[-2]);

    left += 1;
    *pivp = *outerleft;

    BLQS_TYPE swbuf[SWSZ];
    BLQS_TYPE *rwr = right, *sw = swbuf;
    BLQS_TYPE *lwr = left;

    while (UNROLL < SWSZ - (sw - swbuf) && left < right - UNROLL) {
        ptrdiff_t avail = min(right - left, SWSZ - (sw - swbuf));
        BLQS_TYPE* endp = right - avail;
        while (right > endp + UNROLL) {
            for (int i = UNROLL; i--;) {
                BLQS_TYPE x = *right--;
                if (BLQS_CMP(x, piv)) { *sw = x; sw++; }
                else { *rwr = x; rwr--; }
            }
        }
    }

    while (right - left >= UNROLL &&
            (rwr - right > UNROLL || left - lwr > UNROLL)) {

        while (rwr - right > UNROLL && right - left >= UNROLL) {
            for (int i = UNROLL; i--;) {
                BLQS_TYPE x = *left++;
                if (BLQS_CMP(x, piv)) { *lwr = x; lwr++; }
                else { *rwr = x; rwr--; }
            }
        }
        while (left - lwr > UNROLL && right - left >= UNROLL) {
            for (int i = UNROLL; i--;) {
                BLQS_TYPE x = *right--;
                if (BLQS_CMP(x, piv)) { *lwr = x; lwr++; }
                else { *rwr = x; rwr--; }
            }
        }
    }
    do {
        while (rwr > right && left <= right) {
            BLQS_TYPE x = *left++;
            if (BLQS_CMP(x, piv)) { *lwr = x; lwr++; }
            else { *rwr = x; rwr--; }
        }
        while (lwr < left && left <= right) {
            BLQS_TYPE x = *right--;
            if (BLQS_CMP(x, piv)) { *lwr = x; lwr++; }
            else { *rwr = x; rwr--; }
        }
    } while ((lwr < left||rwr > right) && left <= right);

    while (left <= right && !BLQS_CMP(*right, piv)) {
        right--;
        rwr--;
    }
    memcpy(lwr, swbuf, (sw - swbuf) * sizeof(BLQS_TYPE));
    *outerleft = *rwr;
    *rwr = piv;
    return rwr;
}

static void smallsort(BLQS_TYPE* left, BLQS_TYPE* right) {
    while (right - left > 11) {
        BLQS_TYPE* mid = partition_small(left, right);
        smallsort(left, mid - 1);
        left = mid + 1;
    }
    sorting_network(left, right - left);
}

static void sortr(BLQS_TYPE* left, BLQS_TYPE* right) {
    while (1) {
        ptrdiff_t partszm1 = right - left;
        if (partszm1 <= SMALLPART) break;
        BLQS_TYPE* mid = partition(left, right);

        if (mid - left < partszm1 / 16) {
            if (mid > left) sortr(left, mid - 1);
            BLQS_TYPE piv = *mid;
            mid += 1;
            // collect duplicates
            for (BLQS_TYPE* p = mid; p <= right; p++) {
                if (!BLQS_CMP(piv, *p)) {
                    BLQS_TYPE h = *mid;
                    *mid = *p;
                    *p = h;
                    mid++;
                }
            }
            left = mid;
            if (right - left < SMALLPART) break;
            mid = partition(left, right);
        }
        if (mid - left < right - mid) {
            sortr(left, mid - 1);
            left = mid + 1;
        }
        else {
            sortr(mid + 1, right);
            right = mid - 1;
        }
    }
    smallsort(left, right);
}


static void sort(BLQS_TYPE* data, int len) {
    if (len < 2) return;
    sortr(data, data + len - 1);
}

#endif

test.c - sorting 50 million doubles

// SPDX-License-Identifier: MIT
#include <stdio.h>
#include <stdlib.h>
#include <time.h>
#include <sys/time.h>

#define BLQS_CMP(a, b) ((a) < (b))
#define BLQS_TYPE double
#include "sort.h"

#define SIZE 50000000
double data[SIZE];

double ts(void) {
    struct timeval tv;
    gettimeofday(&tv, NULL);
    return tv.tv_sec + tv.tv_usec / 1000000.0;
}

int main() {
    double t0;
    for (int i = 0; i < SIZE; i++) data[i] = rand() / 1024.0;
    t0 = ts();
    sort(data, SIZE);
    printf("Time: %.2fs\n", ts() - t0);
}

On macOS/M1 (Clang, -O3):

Time: 4.39

C++ std::sort needs 1.33 seconds for this.

A few cosmetic changes

It is already micro-optimized using sorting networks and loop unrolling. Only a few cosmetic changes remain.

We rewrite this beginner‑friendly style, which explicitly shows how the pointers are moved:

if (BLQS_CMP(x, piv)) { *lwr = x; lwr++; }
else { *rwr = x; rwr--; }

into a more idiomatic and compact C form:

if (BLQS_CMP(x, piv)) *lwr++ = x;
else *rwr-- = x;

sort.h - rewritten

// SPDX-License-Identifier: MIT
// blqsort.h - Branchless Quicksort
// (c) christof.kaser@gmail.com

#ifndef SORT_H
#define SORT_H

#ifndef BLQS_CMP
#define BLQS_CMP(a, b) ((a) < (b))
#endif

#include <stddef.h>
#include <string.h>
#define min(a, b) (((a) < (b)) ? (a) : (b))

#define SMALLPART 2048
#define SWSZ 512
#define UNROLL 16

#define sort2(a, b) do {  \
    unsigned m = BLQS_CMP(a, b); \
    BLQS_TYPE x = a; \
    a = m ? a : b;  \
    b = m ? b : x;  \
} while(0)

#define sort3(a, b, c) do { \
    sort2(a, b); sort2(b, c); sort2(a, b); \
} while(0)

#define sort4(a, b, c, d) do { \
    sort2(a, b); sort2(c, d); sort2(a, c); \
    sort2(b, d); sort2(b, c); \
} while(0)

#define sort5(a, b, c, d, e) do { \
    sort2(b, c); sort2(d, e); sort2(b, d); \
    sort2(a, c); sort2(a, d); sort2(c, e); \
    sort2(a, b); sort2(c, d); sort2(b, c); \
} while(0)

#define sort6(a, b, c, d, e, f) do { \
    sort2(a, b); sort2(c, d); sort2(e, f); \
    sort2(a, c); sort2(b, d); sort2(e, f); \
    sort2(a, e); sort2(b, f); sort2(c, e); \
    sort2(d, f); sort2(b, c); sort2(d, e); \
    sort2(c, d); \
} while(0)

#define sort7(a, b, c, d, e, f, g) do { \
    sort2(a, b); sort2(c, d); sort2(a, c); \
    sort2(b, d); sort2(b, c); sort2(e, f); \
    sort2(e, g); sort2(f, g); sort2(a, e); \
    sort2(b, f); sort2(c, g); sort2(b, e); \
    sort2(d, g); sort2(c, e); sort2(b, c); \
    sort2(d, f); sort2(d, e); \
} while(0)

#define sort8(a,b,c,d,e,f,g,h) do { \
    sort2(a,b); sort2(c,d); sort2(e,f); sort2(g,h); \
    sort2(a,c); sort2(b,d); sort2(e,g); sort2(f,h); \
    sort2(b,c); sort2(f,g); \
    sort2(a,e); sort2(b,f); sort2(c,g); sort2(d,h); \
    sort2(c,e); sort2(d,f); \
    sort2(b,c); sort2(d,e); sort2(f,g); \
} while (0)

#define sort9(a,b,c,d,e,f,g,h,i) do { \
    sort2(a,d); sort2(b,h); sort2(c,f); sort2(e,i); \
    sort2(a,h); sort2(c,e); sort2(d,i); sort2(f,g); \
    sort2(a,c); sort2(b,d); sort2(e,f); sort2(h,i); \
    sort2(b,e); sort2(d,g); sort2(f,h); \
    sort2(a,b); sort2(c,e); sort2(d,f); sort2(g,i); \
    sort2(c,d); sort2(e,f); sort2(g,h); \
    sort2(b,c); sort2(d,e); sort2(f,g); \
} while (0)

#define sort10(a,b,c,d,e,f,g,h,i,j) do { \
    sort2(a,i); sort2(b,j); sort2(c,h); sort2(d,f); sort2(e,g); \
    sort2(a,c); sort2(b,e); sort2(f,i); sort2(h,j); \
    sort2(a,d); sort2(c,e); sort2(f,h); sort2(g,j); \
    sort2(a,b); sort2(d,g); sort2(i,j); \
    sort2(b,f); sort2(c,d); sort2(e,i); sort2(g,h); \
    sort2(b,c); sort2(d,f); sort2(e,g); sort2(h,i); \
    sort2(c,d); sort2(e,f); sort2(g,h); \
    sort2(d,e); sort2(f,g); \
} while (0)

#define sort11(a,b,c,d,e,f,g,h,i,j,k) do { \
    sort2(a,j); sort2(b,g); sort2(c,e); sort2(d,h); sort2(f,i); \
    sort2(a,b); sort2(d,f); sort2(e,k); sort2(g,j); sort2(h,i); \
    sort2(b,d); sort2(c,f); sort2(e,h); sort2(i,k); \
    sort2(a,e); sort2(b,c); sort2(d,h); sort2(f,j); sort2(g,i); \
    sort2(a,b); sort2(c,g); sort2(e,f); sort2(h,i); sort2(j,k); \
    sort2(c,e); sort2(d,g); sort2(f,h); sort2(i,j); \
    sort2(b,c); sort2(d,e); sort2(f,g); sort2(h,i); \
    sort2(c,d); sort2(e,f); sort2(g,h); \
} while (0)

#define sort12(a,b,c,d,e,f,g,h,i,j,k,l) do { \
    sort2(a,i); sort2(b,h); sort2(c,g); sort2(d,l); sort2(e,k); sort2(f,j); \
    sort2(a,c); sort2(b,e); sort2(d,f); sort2(g,i); sort2(h,k); sort2(j,l); \
    sort2(a,b); sort2(c,j); sort2(e,h); sort2(f,g); sort2(k,l); \
    sort2(b,d); sort2(c,h); sort2(e,j); sort2(i,k); \
    sort2(a,b); sort2(c,d); sort2(e,f); sort2(g,h); sort2(i,j); sort2(k,l); \
    sort2(b,c); sort2(d,f); sort2(g,i); sort2(j,k); \
    sort2(c,e); sort2(d,g); sort2(f,i); sort2(h,j); \
    sort2(b,c); sort2(d,e); sort2(f,g); sort2(h,i); sort2(j,k); \
} while (0)

static void sorting_network(BLQS_TYPE* l, int partszm1_min_1) {
    switch (partszm1_min_1) {
    case 0: break;
    case 1: sort2(l[0],l[1]); break;
    case 2: sort3(l[0],l[1],l[2]); break;
    case 3: sort4(l[0],l[1],l[2],l[3]); break;
    case 4: sort5(l[0],l[1],l[2],l[3],l[4]); break;
    case 5: sort6(l[0],l[1],l[2],l[3],l[4],l[5]); break;
    case 6: sort7(l[0],l[1],l[2],l[3],l[4],l[5],l[6]); break;
    case 7: sort8(l[0],l[1],l[2],l[3],l[4],l[5],l[6],l[7]); break;
    case 8: sort9(l[0],l[1],l[2],l[3],l[4],l[5],l[6],l[7],l[8]); break;
    case 9: sort10(l[0],l[1],l[2],l[3],l[4],l[5],l[6],l[7],l[8],l[9]); break;
    case 10: sort11(l[0],l[1],l[2],l[3],l[4],l[5],l[6],l[7],l[8],l[9],l[10]); break;
    case 11: sort12(l[0],l[1],l[2],l[3],l[4],l[5],l[6],l[7],l[8],l[9],l[10],l[11]); break;
    }
}
#define med5(a,b,c,d,e) do { \
    sort2(a,b); sort2(c,d); sort2(a,c); \
    sort2(b,d); sort2(b,c); sort2(c,e); \
    sort2(b,c); \
} while(0)

static BLQS_TYPE* partition_small(BLQS_TYPE* left, BLQS_TYPE* right) {

    BLQS_TYPE* outerleft = left;
    BLQS_TYPE* pivp = left + 6;

    BLQS_TYPE piv = *pivp;

    BLQS_TYPE l1 = left[1],l2 = left[2];
    BLQS_TYPE r1 = right[-1], r0 = *right;
    med5(l1, l2, piv, r1, r0);
    left[1] = l1; left[2] = l2;
    right[-1] = r1; *right = r0;
    left += 3;
    right -= 2;

    *pivp = *outerleft;

    BLQS_TYPE swbuf[SMALLPART];
    BLQS_TYPE* sw = swbuf;
    BLQS_TYPE* lwr = left;

    while (right - left >= UNROLL) for (int i = UNROLL; i--;) {
        BLQS_TYPE x = *left++;
        if (BLQS_CMP(x, piv)) *lwr++ = x;
        else *sw++ = x;
    }
    while (left <= right) {
        BLQS_TYPE x = *left++;
        if (BLQS_CMP(x, piv)) *lwr++ = x;
        else *sw++ = x;
    }
    memcpy(lwr, swbuf, (sw - swbuf) * sizeof(BLQS_TYPE));
    lwr -= 1;
    *outerleft = *lwr;
    *lwr = piv;
    return lwr;
}

static BLQS_TYPE* partition(BLQS_TYPE* left, BLQS_TYPE* right) {
    BLQS_TYPE* outerleft = left;
    BLQS_TYPE* pivp = left + (right - left) / 2;

    BLQS_TYPE piv = *pivp;

    med5(left[1],left[2],left[3],left[4],left[5]);
    med5(left[11],left[12],left[13],left[14],left[15]);
    med5(pivp[-2], pivp[-1], piv, pivp[1], pivp[2]);
    med5(right[-14], right[-13], right[-12], right[-11], right[-10]);
    med5(right[-4], right[-3], right[-2], right[-1], right[0]);
    med5(left[3], left[13], piv, right[-12], right[-2]);

    left += 1;
    *pivp = *outerleft;

    BLQS_TYPE swbuf[SWSZ];
    BLQS_TYPE *rwr = right, *sw = swbuf;
    BLQS_TYPE *lwr = left;

    while (UNROLL < SWSZ - (sw - swbuf) && left < right - UNROLL) {
        ptrdiff_t avail = min(right - left, SWSZ - (sw - swbuf));
        BLQS_TYPE* endp = right - avail;
        while (right > endp + UNROLL) {
            for (int i = UNROLL; i--;) {
                BLQS_TYPE x = *right--;
                if (BLQS_CMP(x, piv)) *sw++ = x;
                else *rwr-- = x;

            }
        }
    }

    while (right - left >= UNROLL &&
            (rwr - right > UNROLL || left - lwr > UNROLL)) {

        while (rwr - right > UNROLL && right - left >= UNROLL) {
            for (int i = UNROLL; i--;) {
                BLQS_TYPE x = *left++;
                if (BLQS_CMP(x, piv)) *lwr++ = x;
                else *rwr-- = x;
            }
        }
        while (left - lwr > UNROLL && right - left >= UNROLL) {
            for (int i = UNROLL; i--;) {
                BLQS_TYPE x = *right--;
                if (BLQS_CMP(x, piv)) *lwr++ = x;
                else *rwr-- = x;
            }
        }
    }

    do {
        while (rwr > right && left <= right) {
            BLQS_TYPE x = *left++;
            if (BLQS_CMP(x, piv)) *lwr++ = x;
            else *rwr-- = x;
        }
        while (lwr < left && left <= right) {
            BLQS_TYPE x = *right--;
            if (BLQS_CMP(x, piv)) *lwr++ = x;
            else *rwr-- = x;
        }
    } while ((lwr < left||rwr > right) && left <= right);

    while (left <= right && !BLQS_CMP(*right, piv)) {
        right--;
        rwr--;
    }
    memcpy(lwr, swbuf, (sw - swbuf) * sizeof(BLQS_TYPE));
    *outerleft = *rwr;
    *rwr = piv;
    return rwr;
}

static void smallsort(BLQS_TYPE* left, BLQS_TYPE* right) {
    while (right - left > 11) {
        BLQS_TYPE* mid = partition_small(left, right);
        smallsort(left, mid - 1);
        left = mid + 1;
    }
    sorting_network(left, right - left);
}

static void sortr(BLQS_TYPE* left, BLQS_TYPE* right) {
    while (1) {
        ptrdiff_t partszm1 = right - left;
        if (partszm1 <= SMALLPART) break;
        BLQS_TYPE* mid = partition(left, right);

        if (mid - left < partszm1 / 16) {
            if (mid > left) sortr(left, mid - 1);
            BLQS_TYPE piv = *mid;
            mid += 1;
            // collect duplicates
            for (BLQS_TYPE* p = mid; p <= right; p++) {
                if (!BLQS_CMP(piv, *p)) {
                    BLQS_TYPE h = *mid;
                    *mid = *p;
                    *p = h;
                    mid++;
                }
            }
            left = mid;
            if (right - left < SMALLPART) break;
            mid = partition(left, right);
        }
        if (mid - left < right - mid) {
            sortr(left, mid - 1);
            left = mid + 1;
        }
        else {
            sortr(mid + 1, right);
            right = mid - 1;
        }
    }
    smallsort(left, right);
}

static void sort(BLQS_TYPE* data, int len) {
    if (len < 2) return;
    sortr(data, data + len - 1);
}

#endif

A quick test to make sure everything still works.

Time: 0.70

More than 6 times faster than before, and nearly twice as fast as std::sort. That’s quite something.

But what actually happened?

This “small cosmetic” change causes Clang to replace branches with csel.

With branches

; x20 = left, x9 = right, d8 = pivot

loop:
    ldr  d0, [x12], #8
    fcmp d0, d8         ; compare value against pivot
    b.pl ge_case
    str  d0, [x20], #8  ; left++
    b    next
ge_case:
    str  d0, [x9], #-8  ; right--
next:
    cmp  x12, x_end
    b.lt loop

Fast with csel (branchless)

; x9 = current element pointer ; x0 = left destination (grows forwards)
; x20 = right destination (grows backwards); d8  = pivot
; x8  = 8  (positive increment) ; x10 = -8 (negative increment)

loop:
    ldr d0, [x9], #8
    fcmp  d0, d8            ; compare value against pivot
    csel x11, x0, x20, mi   ; target is left (x0) or right (x20)
    str d0, [x11]           ; store value to the selected destination
    csel  x12, x8, xzr, mi  ; left_step is 8 or 0
    csel  x13, xzr, x10, mi ; right_step is 0 or -8
    add   x0, x0, x12       ; left += left_step
    add   x20, x20, x13     ; right += right_step
    b.ls  loop

On x86, Clang behaves similarly: with the compact if, it generates branchless code using cmov (conditional move).

GCC does not exhibit this “quirk” (different code generation for logically equivalent source). It consistently generates the slower branch-based version.

Links

blqsort - Fast Quicksort with C and C++ Interface

When ‘if’ slows you down, avoid it

Interactive sorting demo

christof.kaser@gmail.com

Hacker Times