I've done this exact task before. But it was mainly to measure power consumption and CPU temperatures. The following code (which is fairly long) achieves close to optimal on my Core i7 2600K.
The key thing to note here is the massive amount of manual loop-unrolling as well as interleaving of multiplies and adds...
The full project can be found on my GitHub: https://github.com/Mysticial/Flops
Warning:
If you decide to compile and run this, pay attention to your CPU temperatures!!!
Make sure you don't overheat it. And make sure CPU-throttling doesn't affect your results!
Furthermore, I take no responsibility for whatever damage that may result from running this code.
Notes:
- This code is optimized for x64. x86 doesn't have enough registers for this to compile well.
- This code has been tested to work well on Visual Studio 2010/2012 and GCC 4.6.
ICC 11 (Intel Compiler 11) surprisingly has trouble compiling it well.
- These are for pre-FMA processors. In order to achieve peak FLOPS on Intel Haswell and AMD Bulldozer processors (and later), FMA (Fused Multiply Add) instructions will be needed. These are beyond the scope of this benchmark.
#include <emmintrin.h>
#include <omp.h>
#include <iostream>
using namespace std;
typedef unsigned long long uint64;
double test_dp_mac_SSE(double x,double y,uint64 iterations){
register __m128d r0,r1,r2,r3,r4,r5,r6,r7,r8,r9,rA,rB,rC,rD,rE,rF;
// Generate starting data.
r0 = _mm_set1_pd(x);
r1 = _mm_set1_pd(y);
r8 = _mm_set1_pd(-0.0);
r2 = _mm_xor_pd(r0,r8);
r3 = _mm_or_pd(r0,r8);
r4 = _mm_andnot_pd(r8,r0);
r5 = _mm_mul_pd(r1,_mm_set1_pd(0.37796447300922722721));
r6 = _mm_mul_pd(r1,_mm_set1_pd(0.24253562503633297352));
r7 = _mm_mul_pd(r1,_mm_set1_pd(4.1231056256176605498));
r8 = _mm_add_pd(r0,_mm_set1_pd(0.37796447300922722721));
r9 = _mm_add_pd(r1,_mm_set1_pd(0.24253562503633297352));
rA = _mm_sub_pd(r0,_mm_set1_pd(4.1231056256176605498));
rB = _mm_sub_pd(r1,_mm_set1_pd(4.1231056256176605498));
rC = _mm_set1_pd(1.4142135623730950488);
rD = _mm_set1_pd(1.7320508075688772935);
rE = _mm_set1_pd(0.57735026918962576451);
rF = _mm_set1_pd(0.70710678118654752440);
uint64 iMASK = 0x800fffffffffffffull;
__m128d MASK = _mm_set1_pd(*(double*)&iMASK);
__m128d vONE = _mm_set1_pd(1.0);
uint64 c = 0;
while (c < iterations){
size_t i = 0;
while (i < 1000){
// Here's the meat - the part that really matters.
r0 = _mm_mul_pd(r0,rC);
r1 = _mm_add_pd(r1,rD);
r2 = _mm_mul_pd(r2,rE);
r3 = _mm_sub_pd(r3,rF);
r4 = _mm_mul_pd(r4,rC);
r5 = _mm_add_pd(r5,rD);
r6 = _mm_mul_pd(r6,rE);
r7 = _mm_sub_pd(r7,rF);
r8 = _mm_mul_pd(r8,rC);
r9 = _mm_add_pd(r9,rD);
rA = _mm_mul_pd(rA,rE);
rB = _mm_sub_pd(rB,rF);
r0 = _mm_add_pd(r0,rF);
r1 = _mm_mul_pd(r1,rE);
r2 = _mm_sub_pd(r2,rD);
r3 = _mm_mul_pd(r3,rC);
r4 = _mm_add_pd(r4,rF);
r5 = _mm_mul_pd(r5,rE);
r6 = _mm_sub_pd(r6,rD);
r7 = _mm_mul_pd(r7,rC);
r8 = _mm_add_pd(r8,rF);
r9 = _mm_mul_pd(r9,rE);
rA = _mm_sub_pd(rA,rD);
rB = _mm_mul_pd(rB,rC);
r0 = _mm_mul_pd(r0,rC);
r1 = _mm_add_pd(r1,rD);
r2 = _mm_mul_pd(r2,rE);
r3 = _mm_sub_pd(r3,rF);
r4 = _mm_mul_pd(r4,rC);
r5 = _mm_add_pd(r5,rD);
r6 = _mm_mul_pd(r6,rE);
r7 = _mm_sub_pd(r7,rF);
r8 = _mm_mul_pd(r8,rC);
r9 = _mm_add_pd(r9,rD);
rA = _mm_mul_pd(rA,rE);
rB = _mm_sub_pd(rB,rF);
r0 = _mm_add_pd(r0,rF);
r1 = _mm_mul_pd(r1,rE);
r2 = _mm_sub_pd(r2,rD);
r3 = _mm_mul_pd(r3,rC);
r4 = _mm_add_pd(r4,rF);
r5 = _mm_mul_pd(r5,rE);
r6 = _mm_sub_pd(r6,rD);
r7 = _mm_mul_pd(r7,rC);
r8 = _mm_add_pd(r8,rF);
r9 = _mm_mul_pd(r9,rE);
rA = _mm_sub_pd(rA,rD);
rB = _mm_mul_pd(rB,rC);
i++;
}
// Need to renormalize to prevent denormal/overflow.
r0 = _mm_and_pd(r0,MASK);
r1 = _mm_and_pd(r1,MASK);
r2 = _mm_and_pd(r2,MASK);
r3 = _mm_and_pd(r3,MASK);
r4 = _mm_and_pd(r4,MASK);
r5 = _mm_and_pd(r5,MASK);
r6 = _mm_and_pd(r6,MASK);
r7 = _mm_and_pd(r7,MASK);
r8 = _mm_and_pd(r8,MASK);
r9 = _mm_and_pd(r9,MASK);
rA = _mm_and_pd(rA,MASK);
rB = _mm_and_pd(rB,MASK);
r0 = _mm_or_pd(r0,vONE);
r1 = _mm_or_pd(r1,vONE);
r2 = _mm_or_pd(r2,vONE);
r3 = _mm_or_pd(r3,vONE);
r4 = _mm_or_pd(r4,vONE);
r5 = _mm_or_pd(r5,vONE);
r6 = _mm_or_pd(r6,vONE);
r7 = _mm_or_pd(r7,vONE);
r8 = _mm_or_pd(r8,vONE);
r9 = _mm_or_pd(r9,vONE);
rA = _mm_or_pd(rA,vONE);
rB = _mm_or_pd(rB,vONE);
c++;
}
r0 = _mm_add_pd(r0,r1);
r2 = _mm_add_pd(r2,r3);
r4 = _mm_add_pd(r4,r5);
r6 = _mm_add_pd(r6,r7);
r8 = _mm_add_pd(r8,r9);
rA = _mm_add_pd(rA,rB);
r0 = _mm_add_pd(r0,r2);
r4 = _mm_add_pd(r4,r6);
r8 = _mm_add_pd(r8,rA);
r0 = _mm_add_pd(r0,r4);
r0 = _mm_add_pd(r0,r8);
// Prevent Dead Code Elimination
double out = 0;
__m128d temp = r0;
out += ((double*)&temp)[0];
out += ((double*)&temp)[1];
return out;
}
void test_dp_mac_SSE(int tds,uint64 iterations){
double *sum = (double*)malloc(tds * sizeof(double));
double start = omp_get_wtime();
#pragma omp parallel num_threads(tds)
{
double ret = test_dp_mac_SSE(1.1,2.1,iterations);
sum[omp_get_thread_num()] = ret;
}
double secs = omp_get_wtime() - start;
uint64 ops = 48 * 1000 * iterations * tds * 2;
cout << "Seconds = " << secs << endl;
cout << "FP Ops = " << ops << endl;
cout << "FLOPs = " << ops / secs << endl;
double out = 0;
int c = 0;
while (c < tds){
out += sum[c++];
}
cout << "sum = " << out << endl;
cout << endl;
free(sum);
}
int main(){
// (threads, iterations)
test_dp_mac_SSE(8,10000000);
system("pause");
}
Output (1 thread, 10000000 iterations) - Compiled with Visual Studio 2010 SP1 - x64 Release:
Seconds = 55.5104
FP Ops = 960000000000
FLOPs = 1.7294e+010
sum = 2.22652
The machine is a Core i7 2600K @ 4.4 GHz. Theoretical SSE peak is 4 flops * 4.4 GHz = 17.6 GFlops. This code achieves 17.3 GFlops - not bad.
Output (8 threads, 10000000 iterations) - Compiled with Visual Studio 2010 SP1 - x64 Release:
Seconds = 117.202
FP Ops = 7680000000000
FLOPs = 6.55279e+010
sum = 17.8122
Theoretical SSE peak is 4 flops * 4 cores * 4.4 GHz = 70.4 GFlops. Actual is 65.5 GFlops.
Let's take this one step further. AVX...
#include <immintrin.h>
#include <omp.h>
#include <iostream>
using namespace std;
typedef unsigned long long uint64;
double test_dp_mac_AVX(double x,double y,uint64 iterations){
register __m256d r0,r1,r2,r3,r4,r5,r6,r7,r8,r9,rA,rB,rC,rD,rE,rF;
// Generate starting data.
r0 = _mm256_set1_pd(x);
r1 = _mm256_set1_pd(y);
r8 = _mm256_set1_pd(-0.0);
r2 = _mm256_xor_pd(r0,r8);
r3 = _mm256_or_pd(r0,r8);
r4 = _mm256_andnot_pd(r8,r0);
r5 = _mm256_mul_pd(r1,_mm256_set1_pd(0.37796447300922722721));
r6 = _mm256_mul_pd(r1,_mm256_set1_pd(0.24253562503633297352));
r7 = _mm256_mul_pd(r1,_mm256_set1_pd(4.1231056256176605498));
r8 = _mm256_add_pd(r0,_mm256_set1_pd(0.37796447300922722721));
r9 = _mm256_add_pd(r1,_mm256_set1_pd(0.24253562503633297352));
rA = _mm256_sub_pd(r0,_mm256_set1_pd(4.1231056256176605498));
rB = _mm256_sub_pd(r1,_mm256_set1_pd(4.1231056256176605498));
rC = _mm256_set1_pd(1.4142135623730950488);
rD = _mm256_set1_pd(1.7320508075688772935);
rE = _mm256_set1_pd(0.57735026918962576451);
rF = _mm256_set1_pd(0.70710678118654752440);
uint64 iMASK = 0x800fffffffffffffull;
__m256d MASK = _mm256_set1_pd(*(double*)&iMASK);
__m256d vONE = _mm256_set1_pd(1.0);
uint64 c = 0;
while (c < iterations){
size_t i = 0;
while (i < 1000){
// Here's the meat - the part that really matters.
r0 = _mm256_mul_pd(r0,rC);
r1 = _mm256_add_pd(r1,rD);
r2 = _mm256_mul_pd(r2,rE);
r3 = _mm256_sub_pd(r3,rF);
r4 = _mm256_mul_pd(r4,rC);
r5 = _mm256_add_pd(r5,rD);
r6 = _mm256_mul_pd(r6,rE);
r7 = _mm256_sub_pd(r7,rF);
r8 = _mm256_mul_pd(r8,rC);
r9 = _mm256_add_pd(r9,rD);
rA = _mm256_mul_pd(rA,rE);
rB = _mm256_sub_pd(rB,rF);
r0 = _mm256_add_pd(r0,rF);
r1 = _mm256_mul_pd(r1,rE);
r2 = _mm256_sub_pd(r2,rD);
r3 = _mm256_mul_pd(r3,rC);
r4 = _mm256_add_pd(r4,rF);
r5 = _mm256_mul_pd(r5,rE);
r6 = _mm256_sub_pd(r6,rD);
r7 = _mm256_mul_pd(r7,rC);
r8 = _mm256_add_pd(r8,rF);
r9 = _mm256_mul_pd(r9,rE);
rA = _mm256_sub_pd(rA,rD);
rB = _mm256_mul_pd(rB,rC);
r0 = _mm256_mul_pd(r0,rC);
r1 = _mm256_add_pd(r1,rD);
r2 = _mm256_mul_pd(r2,rE);
r3 = _mm256_sub_pd(r3,rF);
r4 = _mm256_mul_pd(r4,rC);
r5 = _mm256_add_pd(r5,rD);
r6 = _mm256_mul_pd(r6,rE);
r7 = _mm256_sub_pd(r7,rF);
r8 = _mm256_mul_pd(r8,rC);
r9 = _mm256_add_pd(r9,rD);
rA = _mm256_mul_pd(rA,rE);
rB = _mm256_sub_pd(rB,rF);
r0 = _mm256_add_pd(r0,rF);
r1 = _mm256_mul_pd(r1,rE);
r2 = _mm256_sub_pd(r2,rD);
r3 = _mm256_mul_pd(r3,rC);
r4 = _mm256_add_pd(r4,rF);
r5 = _mm256_mul_pd(r5,rE);
r6 = _mm256_sub_pd(r6,rD);
r7 = _mm256_mul_pd(r7,rC);
r8 = _mm256_add_pd(r8,rF);
r9 = _mm256_mul_pd(r9,rE);
rA = _mm256_sub_pd(rA,rD);
rB = _mm256_mul_pd(rB,rC);
i++;
}
// Need to renormalize to prevent denormal/overflow.
r0 = _mm256_and_pd(r0,MASK);
r1 = _mm256_and_pd(r1,MASK);
r2 = _mm256_and_pd(r2,MASK);
r3 = _mm256_and_pd(r3,MASK);
r4 = _mm256_and_pd(r4,MASK);
r5 = _mm256_and_pd(r5,MASK);
r6 = _mm256_and_pd(r6,MASK);
r7 = _mm256_and_pd(r7,MASK);
r8 = _mm256_and_pd(r8,MASK);
r9 = _mm256_and_pd(r9,MASK);
rA = _mm256_and_pd(rA,MASK);
rB = _mm256_and_pd(rB,MASK);
r0 = _mm256_or_pd(r0,vONE);
r1 = _mm256_or_pd(r1,vONE);
r2 = _mm256_or_pd(r2,vONE);
r3 = _mm256_or_pd(r3,vONE);
r4 = _mm256_or_pd(r4,vONE);
r5 = _mm256_or_pd(r5,vONE);
r6 = _mm256_or_pd(r6,vONE);
r7 = _mm256_or_pd(r7,vONE);
r8 = _mm256_or_pd(r8,vONE);
r9 = _mm256_or_pd(r9,vONE);
rA = _mm256_or_pd(rA,vONE);
rB = _mm256_or_pd(rB,vONE);
c++;
}
r0 = _mm256_add_pd(r0,r1);
r2 = _mm256_add_pd(r2,r3);
r4 = _mm256_add_pd(r4,r5);
r6 = _mm256_add_pd(r6,r7);
r8 = _mm256_add_pd(r8,r9);
rA = _mm256_add_pd(rA,rB);
r0 = _mm256_add_pd(r0,r2);
r4 = _mm256_add_pd(r4,r6);
r8 = _mm256_add_pd(r8,rA);
r0 = _mm256_add_pd(r0,r4);
r0 = _mm256_add_pd(r0,r8);
// Prevent Dead Code Elimination
double out = 0;
__m256d temp = r0;
out += ((double*)&temp)[0];
out += ((double*)&temp)[1];
out += ((double*)&temp)[2];
out += ((double*)&temp)[3];
return out;
}
void test_dp_mac_AVX(int tds,uint64 iterations){
double *sum = (double*)malloc(tds * sizeof(double));
double start = omp_get_wtime();
#pragma omp parallel num_threads(tds)
{
double ret = test_dp_mac_AVX(1.1,2.1,iterations);
sum[omp_get_thread_num()] = ret;
}
double secs = omp_get_wtime() - start;
uint64 ops = 48 * 1000 * iterations * tds * 4;
cout << "Seconds = " << secs << endl;
cout << "FP Ops = " << ops << endl;
cout << "FLOPs = " << ops / secs << endl;
double out = 0;
int c = 0;
while (c < tds){
out += sum[c++];
}
cout << "sum = " << out << endl;
cout << endl;
free(sum);
}
int main(){
// (threads, iterations)
test_dp_mac_AVX(8,10000000);
system("pause");
}
Output (1 thread, 10000000 iterations) - Compiled with Visual Studio 2010 SP1 - x64 Release:
Seconds = 57.4679
FP Ops = 1920000000000
FLOPs = 3.34099e+010
sum = 4.45305
Theoretical AVX peak is 8 flops * 4.4 GHz = 35.2 GFlops. Actual is 33.4 GFlops.
Output (8 threads, 10000000 iterations) - Compiled with Visual Studio 2010 SP1 - x64 Release:
Seconds = 111.119
FP Ops = 15360000000000
FLOPs = 1.3823e+011
sum = 35.6244
Theoretical AVX peak is 8 flops * 4 cores * 4.4 GHz = 140.8 GFlops. Actual is 138.2 GFlops.
Now for some explanations:
The performance critical part is obviously the 48 instructions inside the inner loop. You'll notice that it's broken into 4 blocks of 12 instructions each. Each of these 12 instructions blocks are completely independent from each other - and take on average 6 cycles to execute.
So there's 12 instructions and 6 cycles between issue-to-use. The latency of multiplication is 5 cycles, so it's just enough to avoid latency stalls.
The normalization step is needed to keep the data from over/underflowing. This is needed since the do-nothing code will slowly increase/decrease the magnitude of the data.
So it's actually possible to do better than this if you just use all zeros and get rid of the normalization step. However, since I wrote the benchmark to measure power consumption and temperature, I had to make sure the flops were on "real" data, rather than zeros - as the execution units may very well have special case-handling for zeros that use less power and produce less heat.
More Results:
- Intel Core i7 920 @ 3.5 GHz
- Windows 7 Ultimate x64
- Visual Studio 2010 SP1 - x64 Release
Threads: 1
Seconds = 72.1116
FP Ops = 960000000000
FLOPs = 1.33127e+010
sum = 2.22652
Theoretical SSE Peak: 4 flops * 3.5 GHz = 14.0 GFlops. Actual is 13.3 GFlops.
Threads: 8
Seconds = 149.576
FP Ops = 7680000000000
FLOPs = 5.13452e+010
sum = 17.8122
Theoretical SSE Peak: 4 flops * 4 cores * 3.5 GHz = 56.0 GFlops. Actual is 51.3 GFlops.
My processor temps hit 76C on the multi-threaded run! If you runs these, be sure the results aren't affected by CPU throttling.
- 2 x Intel Xeon X5482 Harpertown @ 3.2 GHz
- Ubuntu Linux 10 x64
- GCC 4.5.2 x64 - (-O2 -msse3 -fopenmp)
Threads: 1
Seconds = 78.3357
FP Ops = 960000000000
FLOPs = 1.22549e+10
sum = 2.22652
Theoretical SSE Peak: 4 flops * 3.2 GHz = 12.8 GFlops. Actual is 12.3 GFlops.
Threads: 8
Seconds = 78.4733
FP Ops = 7680000000000
FLOPs = 9.78676e+10
sum = 17.8122
Theoretical SSE Peak: 4 flops * 8 cores * 3.2 GHz = 102.4 GFlops. Actual is 97.9 GFlops.
Here are theoretical max FLOPs counts (per core) for a number of recent processor microarchitectures and explanation how to achieve them.
In general, to calculate this look up the throughput of the FMA instruction(s) e.g. on https://agner.org/optimize/ or any other microbenchmark result, and multiply
(FMAs per clock) * (vector elements / instruction) * 2 (FLOPs / FMA)
.
Note that achieving this in real code requires very careful tuning (like loop unrolling), and near-zero cache misses, and no bottlenecks on anything else. Modern CPUs have such high FMA throughput that there isn't much room for other instructions to store the results, or to feed them with input. e.g. 2 SIMD loads per clock is also the limit for most x86 CPUs, so a dot product will bottleneck on 2 loads per 1 FMA. A carefully-tuned dense matrix multiply can come close to achieving these numbers, though.
If your workload includes any ADD/SUB or MUL that can't be contracted into FMAs, the theoretical max numbers aren't an appropriate goal for your workload. Haswell/Broadwell have 2-per-clock SIMD FP multiply (on the FMA units), but only 1 per clock SIMD FP add (on a separate vector FP add unit with lower latency). Skylake dropped the separate SIMD FP adder, running add/mul/fma the same at 4c latency, 2-per-clock throughput, for any vector width.
Intel
Note that Celeron/Pentium versions of recent microarchitectures don't support AVX or FMA instructions, only SSE4.2.
Intel Core 2 and Nehalem (SSE/SSE2):
- 4 DP FLOPs/cycle: 2-wide SSE2 addition + 2-wide SSE2 multiplication
- 8 SP FLOPs/cycle: 4-wide SSE addition + 4-wide SSE multiplication
Intel Sandy Bridge/Ivy Bridge (AVX1):
- 8 DP FLOPs/cycle: 4-wide AVX addition + 4-wide AVX multiplication
- 16 SP FLOPs/cycle: 8-wide AVX addition + 8-wide AVX multiplication
Intel Haswell/Broadwell/Skylake/Kaby Lake/Coffee/... (AVX+FMA3):
- 16 DP FLOPs/cycle: two 4-wide FMA (fused multiply-add) instructions
- 32 SP FLOPs/cycle: two 8-wide FMA (fused multiply-add) instructions
- (Using 256-bit vector instructions can reduce max turbo clock speed on some CPUs.)
Intel Skylake-X/Skylake-EP/Cascade Lake/etc (AVX512F) with 1 FMA units: some Xeon Bronze/Silver
- 16 DP FLOPs/cycle: one 8-wide FMA (fused multiply-add) instruction
- 32 SP FLOPs/cycle: one 16-wide FMA (fused multiply-add) instruction
- Same computation throughput as with narrower 256-bit instructions, but speedups can still be possible with AVX512 for wider loads/stores, a few vector operations that don't run on the FMA units like bitwise operations, and wider shuffles.
- (Having 512-bit vector instructions in flight shuts down the vector ALU on port 1. Also reduces the max turbo clock speed, so "cycles" isn't a constant in your performance calculations.)
Intel Skylake-X/Skylake-EP/Cascade Lake/etc (AVX512F) with 2 FMA units: Xeon Gold/Platinum, and i7/i9 high-end desktop (HEDT) chips.
- 32 DP FLOPs/cycle: two 8-wide FMA (fused multiply-add) instructions
- 64 SP FLOPs/cycle: two 16-wide FMA (fused multiply-add) instructions
- (Having 512-bit vector instructions in flight shuts down the vector ALU on port 1. Also reduces the max turbo clock speed.)
Future: Intel Cooper Lake (successor to Cascade Lake) is expected to introduce Brain Float, a float16 format for neural-network workloads, with support for actual SIMD computation on it, unlike the current F16C extension that only has support for load/store with conversion to float32. This should double the FLOP/cycle throughput vs. single-precision on the same hardware.
Current Intel chips only have actual computation directly on standard float16 in the iGPU.
AMD
AMD K10:
- 4 DP FLOPs/cycle: 2-wide SSE2 addition + 2-wide SSE2 multiplication
- 8 SP FLOPs/cycle: 4-wide SSE addition + 4-wide SSE multiplication
AMD Bulldozer/Piledriver/Steamroller/Excavator, per module (two cores):
- 8 DP FLOPs/cycle: 4-wide FMA
- 16 SP FLOPs/cycle: 8-wide FMA
AMD Ryzen
- 8 DP FLOPs/cycle: 4-wide FMA
- 16 SP FLOPs/cycle: 8-wide FMA
x86 low power
Intel Atom (Bonnell/45nm, Saltwell/32nm, Silvermont/22nm):
- 1.5 DP FLOPs/cycle: scalar SSE2 addition + scalar SSE2 multiplication every other cycle
- 6 SP FLOPs/cycle: 4-wide SSE addition + 4-wide SSE multiplication every other cycle
AMD Bobcat:
- 1.5 DP FLOPs/cycle: scalar SSE2 addition + scalar SSE2 multiplication every other cycle
- 4 SP FLOPs/cycle: 4-wide SSE addition every other cycle + 4-wide SSE multiplication every other cycle
AMD Jaguar:
- 3 DP FLOPs/cycle: 4-wide AVX addition every other cycle + 4-wide AVX multiplication in four cycles
- 8 SP FLOPs/cycle: 8-wide AVX addition every other cycle + 8-wide AVX multiplication every other cycle
ARM
ARM Cortex-A9:
- 1.5 DP FLOPs/cycle: scalar addition + scalar multiplication every other cycle
- 4 SP FLOPs/cycle: 4-wide NEON addition every other cycle + 4-wide NEON multiplication every other cycle
ARM Cortex-A15:
- 2 DP FLOPs/cycle: scalar FMA or scalar multiply-add
- 8 SP FLOPs/cycle: 4-wide NEONv2 FMA or 4-wide NEON multiply-add
Qualcomm Krait:
- 2 DP FLOPs/cycle: scalar FMA or scalar multiply-add
- 8 SP FLOPs/cycle: 4-wide NEONv2 FMA or 4-wide NEON multiply-add
IBM POWER
IBM PowerPC A2 (Blue Gene/Q), per core:
- 8 DP FLOPs/cycle: 4-wide QPX FMA every cycle
- SP elements are extended to DP and processed on the same units
IBM PowerPC A2 (Blue Gene/Q), per thread:
- 4 DP FLOPs/cycle: 4-wide QPX FMA every other cycle
- SP elements are extended to DP and processed on the same units
Intel MIC / Xeon Phi
Intel Xeon Phi (Knights Corner), per core:
- 16 DP FLOPs/cycle: 8-wide FMA every cycle
- 32 SP FLOPs/cycle: 16-wide FMA every cycle
Intel Xeon Phi (Knights Corner), per thread:
- 8 DP FLOPs/cycle: 8-wide FMA every other cycle
- 16 SP FLOPs/cycle: 16-wide FMA every other cycle
Intel Xeon Phi (Knights Landing), per core:
- 32 DP FLOPs/cycle: two 8-wide FMA every cycle
- 64 SP FLOPs/cycle: two 16-wide FMA every cycle
The reason why there are per-thread and per-core datum for IBM Blue Gene/Q and Intel Xeon Phi (Knights Corner) is that these cores have a higher instruction issue rate when running more than one thread per core.
Best Answer
Nehalem is capable of executing 4 DP or 8 SP FLOP/cycle. This is accomplished using SSE, which operates on packed floating point values, 2/register in DP and 4/register in SP. In order to achieve 4 DP FLOP/cycle or 8 SP FLOP/cycle the core has to execute 2 SSE instructions per cycle. This is accomplished by executing a MULDP and an ADDDP (or a MULSP and an ADDSP) per cycle. The reason this is possible is because Nehalem has separate execution units for SSE multiply and SSE add, and these units are pipelined so that the
throughput
is one multiply and one add per cycle. Multiplies are in the multiplier pipeline for 4 cycles in SP and 5 cycles in DP. Adds are in the pipeline for 3 cycles independent of SP/DP. The number of cycles in the pipeline is known as thelatency
. To compute peak FLOP/cycle all you need to know is the throughput. So with a throughput of 1 SSE vector instruction/cycle for both the multiplier and the adder (2 execution units) you have 2 x 2 = 4 FLOP/cycle in DP and 2 x 4 = 8 FLOP/cycle in SP. To actually sustain this peak throughput you need to consider latency (so you have at least as many independent operations in the pipeline as the depth of the pipeline) and you need to consider being able to feed the data fast enough. Nehalem has an integrated memory controller capable of very high bandwidth from memory which it can achieve if the data prefetcher correctly anticipates the access pattern of the data (sequentially loading from memory is a trivial pattern that it can anticipate). Typically there isn't enough memory bandwidth to sustain feeding all cores with data at peak FLOP/cycle, so some amount of reuse of the data from the cache is necessary in order to sustain peak FLOP/cycle.Details on where you can find information on the number of independent execution units and their throughput and latency in cycles follows.
See page 105 8.9 Execution units of this document
http://www.agner.org/optimize/microarchitecture.pdf
It says that for Nehalem
In order to get 8 SP FLOP/cycle you need 4 SP ADD/cycle and 4 SP MUL/cycle. The adder and the multiplier are on separate execution units, and dispatch out of separate ports, each can execute on 4 SP packed operands simultaneously using SSE packed (vector) instructions (4x32bit = 128bits). Both have throughput of 1 operation per clock cycle. In order to get that throughput, you need to consider the latency... how many cycles after the instruction issues before you can use the result.. so you have to issue several independent instructions to cover the latency. The multiplier in single precision has a latency of 4 and the adder of 3.
You can find these same throughput and latency numbers for Nehalem in the Intel Optimization guide, table C-15a
http://www.intel.com/content/www/us/en/architecture-and-technology/64-ia-32-architectures-optimization-manual.html