7 ms·

Improving performance with SIMD intrinsics in three use cases

58 points by jbredeche 6y ago

grecy 6y ago
I always find code samples of SIMD instructions to be so esoteric it makes it harder to follow that necessary. It's 2020, who uses variable names like 'acc', 'a' and 'b'? There's also no real explanation of what 'a' and 'b' (or 'p1' and 'p2') actually are. It's just magic. I get that it's just an example, but the example feels so useless it doesn't demonstrate why these SIMD calls are helpful at all.
- chrisseaton 6y ago
  > It's 2020, who uses variable names like 'acc', 'a' and 'b'? Mathematicians, scientists, engineers, artists. Have you read any numerical work? They often use one-letter variables, or otherwise very abbreviated variables. Sometimes they differentiate variables with the same letter using by using different fonts. Seems strange to us, but clearly it works for them! And then they keep doing the same thing when they translate to code. When I ask them they say they think our long variable names are unreadable. Maybe they're right and we're wrong? Who knows. They generally don't use tests, linting, CI, anything like that either. But again, clearly it works for them!
  - superjan 6y ago
    If you write code that does basic arithmatic, a and b will do just fine in my experience. If I’m implementing an existing algorithm, I like to stick as close to the algoritm’s symbol use as close as possible. I miss that I can’t use greek symbols like ϕ or σ as variables, but that’s just me probably.
    - chrisseaton 6y ago
      When I try to implement an algorithm from a paper, I try to rename the variables from Greek letters in German fonts to have simple English names, but the scientist will say 'why are you making it more complicated!' I feel like we'll never understand each other!
    - minerjoe 6y ago
      Depending on your programming language, you can definitely use those greek symbols as variables.
      - sgillen 6y ago
        Yup was just about to say this. That being said I usually avoid doing this even though I do a lot of numerical work. They are usually more of a pain to type out, and I don’t trust that they will show up properly on every text editor that it might wind up in.
    - cinntaile 6y ago
      I believe you can do this in Julia and possibly other languages as well.
    - zoomablemind 6y ago
      When transcoding the formulas from CAS (comptuer algebra), the Greeks are often already spelled out in English transliteration (phi, sigma). Even more, sometimes for clarity or performance reasons it may make sense to introduce interim variables, I often name those just as what they'd look like in formula, say, pi/2, would be PI_2; A*B is kept in AB, or AxB depending on context. So many times the encoded formulas end up hiding accidental errors resulting from a desire to use programmer's names. With scientific code, it's safer to stick to formulas, as it's just easier to eyeball them or outright copy and paste to/from CAS for cross-check. Of course, it's paramount to have proper test-cases to make sure the encoded formulas match to the expected results.
  - OnlyOneCannolo 6y ago
    When people do this, they're usually following a convention where those short symbols are understood to mean certain things in that domain. If there's documentation describing what the code is doing and why it's doing it, that can be much more useful than longer variable names.
- The_rationalist 6y ago
  Acc means accumulator and is a standard name when using e.g a fold/reduce functor. That being said, single letters or non standard abbreviative variable names are a cancer to code clarity as the semantic of what they denote become unaccessible.
- creato 6y ago
  There's no context here, why would it matter what the variables are? Do you really need 'a' and 'b' to be named 'add_operand_1' and 'add_operand_2' to understand what they are? Same with 'acc', is it not obvious what it is from the context? Unnecessarily long names can be less readable too. If a variable's lifetime is only a few lines long, and the meaning is obvious from the context, a single letter name might be best. As for the example: if you think a dot product is useless, then maybe SIMD really is not useful for you, and you can forget about it if you want.
- Const-me 6y ago
  Sorry about the confusion. I’m not a professional writer, quite the opposite, I’m a professional programmer. acc stands for “accumulator”. p1 and p2 are two pointers, of type const float*. You’ll find more context on github, here’s the repository with these samples: https://github.com/Const-me/SimdIntroArticle https://github.com/Const-me/SimdIntroArticle
gnufx 6y ago
Note that the code GCC generates for reductions (dot product, anyhow) will not be scalar with an option like -Ofast. That enables loop vectorization but also use of non-standard-conforming SIMD for the non-associative operations. It likely does as well as simple use of intrinsics. Level 1 operations are only worth worrying about much with small arrays, as they're ultimately memory-bound.
- jeffbee 6y ago
  It's not even scalar at -O1. https://godbolt.org/z/TRS4W6 https://godbolt.org/z/TRS4W6 It generates vector instructions for all k8-like architectures, with wider instructions starting from -march=sandybridge.
  - Const-me 6y ago
    > It generates vector instructions vmovss, vmulss and vaddss instructions are not vector. They only load/store, multiply or add the lowest lane of their operands. The equivalent vector instructions are vmovaps, vmulps and vaddps respectively, these would load, multiply or add complete vectors, not just the lowest lanes of them.
    - jeffbee 6y ago
      Right you are. I totally stopped reading at `vxorps` and didn't scan all the other `_ss` properly.
      - gnufx 6y ago
        There would be something amiss if it was vectorizing at O1, and if it was disregarding arithmetic rules (as the Intel compiler would by default).
  - ethelward 6y ago
    Although AVX instructions, those are scalar. The SS pattern at the end is for Scalar Single-precision, so only one lane of the whole SIMD register is used.
qppo 6y ago
It would be useful to add a recursively defined algorithm like trapezoidal integration, which is something that basically every autovectorization tool will fall over y[n] = (x[n] + x[n-1]) / 2 + y[n-1] Which will quickly show you that you can't just replace the numeric operations with their intrinsic counterparts and get speed up (and if you did, you may wind up with incorrect code!). It also shows you that sometimes more is less. To vectorize that code you wind up doing more multiplication than you would normally require for the scalar, and depending on the conditions where you execute it, it might not be able to improve over the scalar implementation. There are more formal methods for converting these kinds of operations (linear, time-invariant recursive systems/filters) to vector equivalents using state space matrix forms, but then you start having to benchmark whether or not you get any benefits.
- dragontamer 6y ago
  I got nerd-sniped. Lets see... y[n] = (x[n] + x[n-1]) / 2 y = prefixSum(y); I'm pretty sure prefixsum over Y would solve the problem (see: https://en.wikipedia.org/wiki/Prefix_sum#Algorithm_1:_Shorter_span,_more_parallel https://en.wikipedia.org/wiki/Prefix_sum#Algorithm_1:_Shorte... for the Vectorized solution to prefix sum). Without the addition of y[n-1], the first line is now trivially parallelizable by most auto-vectorizers. The prefix-sum probably needs to be written with intrinsics by hand, but once written it generalizes to many operations. (Prefix-max, prefix-min, prefix-XOR). "Prefix-sum" is sometimes called "scan" in the literature (https://developer.nvidia.com/gpugems/gpugems3/part-vi-gpu-computing/chapter-39-parallel-prefix-sum-scan-cuda https://developer.nvidia.com/gpugems/gpugems3/part-vi-gpu-co...) ----------- EDIT: A surprising number of "recursive" definitions can become a prefix-operation (or scan). Associativity is the main requirement. (Note that floating-points are not associative: error builds up non-associatively). So you'd need to redo your error analysis if these were floats.
  - Const-me 6y ago
    > the first line is now trivially parallelizable by most auto-vectorizers It’s not that trivial to vectorize, would be something like this for SSE: const __m128 x = _mm_loadu_ps( px ); // Rotate the lanes __m128 xPrev = _mm_shuffle_ps( x, x, _MM_SHUFFLE( 2, 1, 0, 3 ) ); // Insert the highest lane of previousStepX into lowest lane of xPrev // The IMM number means extract lane 3, insert into lane 0, don't zero out any lanes xPrev = _mm_insert_ps( xPrev, previousStepX, 0b11000000 ); const __m128 halfSum = _mm_mul_ps( _mm_add_ps( x, xPrev ), _mm_set1_ps( 0.5f ) ); previousStepX = x; Which compiler do you think would do that automatically?
dialecticDolt 6y ago
Does anyone know any more sources that would be useful to learn more about SIMD intrinsics and how to use them? I've tried a few times (with moderate success) to parse through source-code using them but I've never found a good introduction.
- QShift 6y ago
  Handmade Hero is one of the best learning resources on SIMD (and, in my opinion, programming in general). Here's the episode where SIMD is properly introduced: https://guide.handmadehero.org/code/day115 https://guide.handmadehero.org/code/day115
- Const-me 6y ago
  If you’re interested in AMD64 SIMD, deep inside that article on SO there’s a link to my earlier article on the subject: http://const.me/articles/simd/simd.pdf http://const.me/articles/simd/simd.pdf I also made offline docs, there: https://github.com/Const-me/IntelIntrinsics https://github.com/Const-me/IntelIntrinsics
- dragontamer 6y ago
  Honestly, most SIMD intrinsics are pretty simple and straightforward to use. You can probably get 90% of the way there by simply browsing Intel's intrinsic reference manual. https://software.intel.com/sites/landingpage/IntrinsicsGuide/ https://software.intel.com/sites/landingpage/IntrinsicsGuide... Intel's intrinsics guide has an approximate clock-count on all instructions. There are a few obscure details to memorize, but if you already understand super-scalar, pipelined, and OoO operations of modern CPUs, its all straight-forward. --------- The last 10% of using SIMD successfully is obscure and eclectic. Its not "difficult", its just a whole bunch of tricks that you'll probably never figure out on your own. I'm talking about pshufb magic, prefix sum / scan patterns, gather/scatter patterns, and the like. As complicated as pshufb is, the Intel intrinsic documentation does a great job precisely describing the operation: https://software.intel.com/sites/landingpage/IntrinsicsGuide/#text=pshufb&expand=5194,5153,5153 https://software.intel.com/sites/landingpage/IntrinsicsGuide... _mm_shuffle_epi8 gives you a register-to-register (no RAM touched, faster than even L1 cache) 16-byte lookup table. An alternative viewpoint: it gives you a register-to-register arbitrary gather command. I've seen so many tricks from this one instruction alone. --------- Another straightforward guide is Intel's optimization manual. Its straightforward, but a lot of material to read: https://www.intel.com/content/dam/www/public/us/en/documents/manuals/64-ia-32-architectures-optimization-manual.pdf https://www.intel.com/content/dam/www/public/us/en/documents... Chapters 4, 5, and 6 deal with SIMD coding.
mayoff 6y ago
If you have a big need for SIMD on the CPU, you might want to look into ispc, the Intel Implicit SPMD Program Compiler. (It supports ARM NEON too.) https://ispc.github.io/ https://ispc.github.io/
alleycat5000 6y ago
Easy SIMD in the context of image processing is one of the cool things about Halide: https://halide-lang.org https://halide-lang.org It was personally eye opening to see how vectorization can be the gift they keeps on giving when you can sprinkle it throughout a computational pipeline vs spreading work across cores.
SloopJon 6y ago
It's great to see tables comparing the different approaches, but it's also a bit discouraging to think that you might have to solve the problem a dozen different ways. Fortunately, it does seem that naive vectorization, for lack of a better term, is pretty effective. For grayscale, "vector floats SSE2" and "vector floats AVX2" aren't far behind the winner. Similarly, for dot product, I might call it a day after AvxDpPs dropped my time from 193 to 34.

58 points by jbredeche 6y ago

I always find code samples of SIMD instructions to be so esoteric it makes it harder to follow that necessary. It's 2020, who uses variable names like 'acc', 'a' and 'b'? There's also no real explanation of what 'a' and 'b' (or 'p1' and 'p2') actually are. It's just magic. I get that it's just an example, but the example feels so useless it doesn't demonstrate why these SIMD calls are helpful at all.

chrisseaton 6y ago
> It's 2020, who uses variable names like 'acc', 'a' and 'b'? Mathematicians, scientists, engineers, artists. Have you read any numerical work? They often use one-letter variables, or otherwise very abbreviated variables. Sometimes they differentiate variables with the same letter using by using different fonts. Seems strange to us, but clearly it works for them! And then they keep doing the same thing when they translate to code. When I ask them they say they think our long variable names are unreadable. Maybe they're right and we're wrong? Who knows. They generally don't use tests, linting, CI, anything like that either. But again, clearly it works for them!
- superjan 6y ago
  If you write code that does basic arithmatic, a and b will do just fine in my experience. If I’m implementing an existing algorithm, I like to stick as close to the algoritm’s symbol use as close as possible. I miss that I can’t use greek symbols like ϕ or σ as variables, but that’s just me probably.
  - chrisseaton 6y ago
    When I try to implement an algorithm from a paper, I try to rename the variables from Greek letters in German fonts to have simple English names, but the scientist will say 'why are you making it more complicated!' I feel like we'll never understand each other!
  - minerjoe 6y ago
    Depending on your programming language, you can definitely use those greek symbols as variables.
    - sgillen 6y ago
      Yup was just about to say this. That being said I usually avoid doing this even though I do a lot of numerical work. They are usually more of a pain to type out, and I don’t trust that they will show up properly on every text editor that it might wind up in.
  - cinntaile 6y ago
    I believe you can do this in Julia and possibly other languages as well.
  - zoomablemind 6y ago
    When transcoding the formulas from CAS (comptuer algebra), the Greeks are often already spelled out in English transliteration (phi, sigma). Even more, sometimes for clarity or performance reasons it may make sense to introduce interim variables, I often name those just as what they'd look like in formula, say, pi/2, would be PI_2; A*B is kept in AB, or AxB depending on context. So many times the encoded formulas end up hiding accidental errors resulting from a desire to use programmer's names. With scientific code, it's safer to stick to formulas, as it's just easier to eyeball them or outright copy and paste to/from CAS for cross-check. Of course, it's paramount to have proper test-cases to make sure the encoded formulas match to the expected results.
- OnlyOneCannolo 6y ago
  When people do this, they're usually following a convention where those short symbols are understood to mean certain things in that domain. If there's documentation describing what the code is doing and why it's doing it, that can be much more useful than longer variable names.
The_rationalist 6y ago
Acc means accumulator and is a standard name when using e.g a fold/reduce functor. That being said, single letters or non standard abbreviative variable names are a cancer to code clarity as the semantic of what they denote become unaccessible.
creato 6y ago
There's no context here, why would it matter what the variables are? Do you really need 'a' and 'b' to be named 'add_operand_1' and 'add_operand_2' to understand what they are? Same with 'acc', is it not obvious what it is from the context? Unnecessarily long names can be less readable too. If a variable's lifetime is only a few lines long, and the meaning is obvious from the context, a single letter name might be best. As for the example: if you think a dot product is useless, then maybe SIMD really is not useful for you, and you can forget about it if you want.
Const-me 6y ago
Sorry about the confusion. I’m not a professional writer, quite the opposite, I’m a professional programmer. acc stands for “accumulator”. p1 and p2 are two pointers, of type const float*. You’ll find more context on github, here’s the repository with these samples: https://github.com/Const-me/SimdIntroArticle https://github.com/Const-me/SimdIntroArticle

gnufx 6y ago

Note that the code GCC generates for reductions (dot product, anyhow) will not be scalar with an option like -Ofast. That enables loop vectorization but also use of non-standard-conforming SIMD for the non-associative operations. It likely does as well as simple use of intrinsics. Level 1 operations are only worth worrying about much with small arrays, as they're ultimately memory-bound.

jeffbee 6y ago
It's not even scalar at -O1. https://godbolt.org/z/TRS4W6 https://godbolt.org/z/TRS4W6 It generates vector instructions for all k8-like architectures, with wider instructions starting from -march=sandybridge.
- Const-me 6y ago
  > It generates vector instructions vmovss, vmulss and vaddss instructions are not vector. They only load/store, multiply or add the lowest lane of their operands. The equivalent vector instructions are vmovaps, vmulps and vaddps respectively, these would load, multiply or add complete vectors, not just the lowest lanes of them.
  - jeffbee 6y ago
    Right you are. I totally stopped reading at `vxorps` and didn't scan all the other `_ss` properly.
    - gnufx 6y ago
      There would be something amiss if it was vectorizing at O1, and if it was disregarding arithmetic rules (as the Intel compiler would by default).
- ethelward 6y ago
  Although AVX instructions, those are scalar. The SS pattern at the end is for Scalar Single-precision, so only one lane of the whole SIMD register is used.

qppo 6y ago

It would be useful to add a recursively defined algorithm like trapezoidal integration, which is something that basically every autovectorization tool will fall over y[n] = (x[n] + x[n-1]) / 2 + y[n-1] Which will quickly show you that you can't just replace the numeric operations with their intrinsic counterparts and get speed up (and if you did, you may wind up with incorrect code!). It also shows you that sometimes more is less. To vectorize that code you wind up doing more multiplication than you would normally require for the scalar, and depending on the conditions where you execute it, it might not be able to improve over the scalar implementation. There are more formal methods for converting these kinds of operations (linear, time-invariant recursive systems/filters) to vector equivalents using state space matrix forms, but then you start having to benchmark whether or not you get any benefits.

dragontamer 6y ago
I got nerd-sniped. Lets see... y[n] = (x[n] + x[n-1]) / 2 y = prefixSum(y); I'm pretty sure prefixsum over Y would solve the problem (see: https://en.wikipedia.org/wiki/Prefix_sum#Algorithm_1:_Shorter_span,_more_parallel https://en.wikipedia.org/wiki/Prefix_sum#Algorithm_1:_Shorte... for the Vectorized solution to prefix sum). Without the addition of y[n-1], the first line is now trivially parallelizable by most auto-vectorizers. The prefix-sum probably needs to be written with intrinsics by hand, but once written it generalizes to many operations. (Prefix-max, prefix-min, prefix-XOR). "Prefix-sum" is sometimes called "scan" in the literature (https://developer.nvidia.com/gpugems/gpugems3/part-vi-gpu-computing/chapter-39-parallel-prefix-sum-scan-cuda https://developer.nvidia.com/gpugems/gpugems3/part-vi-gpu-co...) ----------- EDIT: A surprising number of "recursive" definitions can become a prefix-operation (or scan). Associativity is the main requirement. (Note that floating-points are not associative: error builds up non-associatively). So you'd need to redo your error analysis if these were floats.
- Const-me 6y ago
  > the first line is now trivially parallelizable by most auto-vectorizers It’s not that trivial to vectorize, would be something like this for SSE: const __m128 x = _mm_loadu_ps( px ); // Rotate the lanes __m128 xPrev = _mm_shuffle_ps( x, x, _MM_SHUFFLE( 2, 1, 0, 3 ) ); // Insert the highest lane of previousStepX into lowest lane of xPrev // The IMM number means extract lane 3, insert into lane 0, don't zero out any lanes xPrev = _mm_insert_ps( xPrev, previousStepX, 0b11000000 ); const __m128 halfSum = _mm_mul_ps( _mm_add_ps( x, xPrev ), _mm_set1_ps( 0.5f ) ); previousStepX = x; Which compiler do you think would do that automatically?

dialecticDolt 6y ago

Does anyone know any more sources that would be useful to learn more about SIMD intrinsics and how to use them? I've tried a few times (with moderate success) to parse through source-code using them but I've never found a good introduction.

QShift 6y ago
Handmade Hero is one of the best learning resources on SIMD (and, in my opinion, programming in general). Here's the episode where SIMD is properly introduced: https://guide.handmadehero.org/code/day115 https://guide.handmadehero.org/code/day115
Const-me 6y ago
If you’re interested in AMD64 SIMD, deep inside that article on SO there’s a link to my earlier article on the subject: http://const.me/articles/simd/simd.pdf http://const.me/articles/simd/simd.pdf I also made offline docs, there: https://github.com/Const-me/IntelIntrinsics https://github.com/Const-me/IntelIntrinsics
dragontamer 6y ago
Honestly, most SIMD intrinsics are pretty simple and straightforward to use. You can probably get 90% of the way there by simply browsing Intel's intrinsic reference manual. https://software.intel.com/sites/landingpage/IntrinsicsGuide/ https://software.intel.com/sites/landingpage/IntrinsicsGuide... Intel's intrinsics guide has an approximate clock-count on all instructions. There are a few obscure details to memorize, but if you already understand super-scalar, pipelined, and OoO operations of modern CPUs, its all straight-forward. --------- The last 10% of using SIMD successfully is obscure and eclectic. Its not "difficult", its just a whole bunch of tricks that you'll probably never figure out on your own. I'm talking about pshufb magic, prefix sum / scan patterns, gather/scatter patterns, and the like. As complicated as pshufb is, the Intel intrinsic documentation does a great job precisely describing the operation: https://software.intel.com/sites/landingpage/IntrinsicsGuide/#text=pshufb&expand=5194,5153,5153 https://software.intel.com/sites/landingpage/IntrinsicsGuide... _mm_shuffle_epi8 gives you a register-to-register (no RAM touched, faster than even L1 cache) 16-byte lookup table. An alternative viewpoint: it gives you a register-to-register arbitrary gather command. I've seen so many tricks from this one instruction alone. --------- Another straightforward guide is Intel's optimization manual. Its straightforward, but a lot of material to read: https://www.intel.com/content/dam/www/public/us/en/documents/manuals/64-ia-32-architectures-optimization-manual.pdf https://www.intel.com/content/dam/www/public/us/en/documents... Chapters 4, 5, and 6 deal with SIMD coding.

mayoff 6y ago

If you have a big need for SIMD on the CPU, you might want to look into ispc, the Intel Implicit SPMD Program Compiler. (It supports ARM NEON too.) https://ispc.github.io/ https://ispc.github.io/

alleycat5000 6y ago

Easy SIMD in the context of image processing is one of the cool things about Halide: https://halide-lang.org https://halide-lang.org It was personally eye opening to see how vectorization can be the gift they keeps on giving when you can sprinkle it throughout a computational pipeline vs spreading work across cores.

SloopJon 6y ago

It's great to see tables comparing the different approaches, but it's also a bit discouraging to think that you might have to solve the problem a dozen different ways. Fortunately, it does seem that naive vectorization, for lack of a better term, is pretty effective. For grayscale, "vector floats SSE2" and "vector floats AVX2" aren't far behind the winner. Similarly, for dot product, I might call it a day after AvxDpPs dropped my time from 193 to 34.