SSE Assembly vs GCC Compiler - Dot Product

Question

I am currently taking an introductory course in computer architecture. Our goal was to write a dot-product function in x86 Assembly which would use SSE and SIMD (without AVX).

I am not to that confident with my solution:

#include <time.h>
#define NUM 25600
//Example: scalarProduct using SIMD
extern float scalarProduct(float *, float *, int size);

float vsC(float * a, float * b, int size){
    float sum = 0;
    for(int i = 0; i < size; i++){
        sum += b[i]*a[i];
    }
    return sum;
}
int main(int argc, char ** argv){
    float * a = malloc(NUM * sizeof(double));
    float * b = malloc(NUM * sizeof(double));
    for(int i = 0; i < NUM; i++){
        a[i] = 1;
        b[i] = 1.0/(i+1);
    }
    clock_t start, end;
        double cpu_time_used;

         start = clock();
     float sum = scalarProduct(a,b,NUM);
         end = clock();
     cpu_time_used = ((double) (end - start))/CLOCKS_PER_SEC;
     printf("%.15f\n", cpu_time_used);
     printf("Solution %.15f\n", (double)(sum));
     start = clock();
     sum = vsC(a,b,NUM);
     end = clock();
     cpu_time_used = ((double) (end - start))/CLOCKS_PER_SEC;
     printf("%.15f\n", cpu_time_used);
     printf("Solution %.15f\n", (double)(sum));
}

Assembly File

.intel_syntax noprefix
.text
.global scalarProduct

scalarProduct:
        mov rax, rdx
        xorps xmm0, xmm0
        mov rcx, 0
start:
        cmp rax, 4
        jl end
        movdqu xmm3, [rsi + rcx]
        movdqu xmm2, [rdi + rcx]
        vmulps xmm1, xmm2, xmm3
        haddps  xmm7, xmm1
        haddps xmm7, xmm7
        psrldq xmm7, 4 //Shift to pos0
        addss xmm0, xmm7
        xorps xmm7, xmm7
        sub rax, 4
        add rcx, 16
        jmp start
end:
        cmp rax, 0
        je ret
        dec rax
        movss xmm1, [rsi + rcx]
        movss xmm2, [rdi + rcx]
        mulss xmm1, xmm2
        addss xmm0, xmm1
        add rcx, 4
        jmp end
ret:
        ret

Obviously, this Assembly - Code is far from perfect. How could I do better using basic SIMD and SSE?

The second thing which made me wonder is, I indeed outperformed GCC on a Xeon processor, which is irritating.

Compiling the code with:

gcc -o main -O7 main.c scalarProduct.S 

Shows the following result:

./main
0.000015000000000
Solution 10.727574348449707
0.000026000000000
Solution 10.727569580078125

How would I have to improve my C Code so that GCC can step up?

DISCLAIMER:

My homework does not affect my grades and its editing is optional.

Answer 1

While harold reviewed your assembly code, I'll just comment on how you are compiling your C code:

Increase the duration of the benchmark

Your code runs for a very short amount of time, only tens of microseconds. This is too short to get accurate measurements:

• clock() only has a resolution of a microsecond, so this is not insignificant compared to the duration.
• the processor's caches, branch predictors and prefetchers might not have had time to fully warm up.
• interrupts and other processes running on your system add noise.

Aim to run your benchmarks for at least a second. That might be easier said than done: if you increase the size of the arrays you might become memory bandwidth bound instead of CPU bound. And if you naively just repeat the computations with the small array, the compiler might notice and just calculate the sum once and the multiply it by the number of repeats, depending on the level of optimization used.

Use a proper optimization level

There is no such thing as -O7. The highest supported optimization level for GCC is -O3. However, even that does not enable all possible optimzations. In particular, GCC is a bit careful when it comes to floating point math, and tries to ensure the code is correct even if there are infinities, NaNs and denormals. It also knows that floating point math is not strictly associative and commutative, and so will try to keep operations in the same order as you specified them, which prevents it from using some vectorization tricks. If you don't care about that, you can enable -ffast-math, or use -Ofast.

Consider using -mtune=... and/or -march=...

If you don't specify any specific CPU, then on an x86-64 platform, GCC will output code that can run on any 64-bit Intel or AMD CPU, and might not be able to use certain SSE instructions that are not available in the x86-64 baseline. Also, the compiler will assume a certain CPU for instruction timing, delay slots, and other micro-architecural optimizations, which might not be ideal for the Xeon CPU you are running the code on. Typically you would use -march=native to ensure the compiler will provide code using all features of the CPU you are compiling on, but that might cause it to use AVX instructions if your CPU supports those.

Consider using SSE intrinsics

Instead of having a pure assembly version and a pure C version, you can have something inbetween by using SSE intrinsics. These are functions that are compiled into specific CPU instructions. However, the function arguments and return values are just variables (either regular ones or special vector type variables), not registers. The compiler will pick registers as it sees fit, and will also be able to reorder the intrinsics calls if possible, using its knowledge about the CPU's micro-architecture, and if you use the intrinsics in a loop then the compiler can unroll the loop for you.

While GCC and Clang are able to vectorize certain loops, they are quite bad at using horizontal operations such as haddps, so if you help it by using intrinsics you might get code comparable to the best hand-optimized assembly.

Some results

If I increase the size of the array to 25600000, and compile with -O7 (which will effectively be -O3), I get the following results on an AMD Ryzen 9 3900X:

0.008444000000000
Solution 16.000000000000000
0.018092000000000
Solution 15.403682708740234

With -Ofast I get:

0.008399000000000
Solution 16.000000000000000
0.006617000000000
Solution 16.419670104980469

So clearly there is a speed-up when going to -Ofast, but the resulting solution is also different.

Answer 2

How could I do better using basic SIMD and SSE?

The most important things are:

Delay horizontal addition as long as possible

haddps costs two shuffles and a normal addition. It is used twice, plus a shift and scalar-add. That's a lot of cost, and none of it is necessary: the main loop can just use addps. When the main loop is done, then you still need horizontal addition, but that cost is only paid once so it's not bad.

By the way, you can horizontally sum xmm1 like this:

haddps  xmm1, xmm1
haddps  xmm1, xmm1

No pre-zeroed register required, and no shift. It wastes a lot of work that haddps does, but it's short and simple.

Use multiple accumulators

When accumulating via addps, its latency becomes a bottleneck. mulps can be executed once or twice per cycle depending on the architecture, while addps has a latency of 3 to 4. Two mulps per cycles is not a reachable goal for a dot product (too much data needs to be loaded), but one per cycle is. Using a single accumulator means the loop is (eventually) limited to 1 iteration every 3 (or 4) cycles, it can get started quicker but a backlog of dependent addps builds up until it start to block progress.

Using multiple accumulators fixes that issue by spreading the work across some independent addps, so progress can be made faster.

As a bonus, the haddps after the loop has more useful work to do.

Use a 1-jump loop

start:
    cmp rax, 4
    jl end
    ...
    jmp start

Is a 2-jump loop, but you can make it a 1-jump loop like this:

    jmp loopentry
start:
    ...
loopentry:
    cmp rax, 4
    jnl start

There are still two jumps, but one of them is not in the loop anymore.