L1 memory bandwidth: 50% drop in efficiency using addresses which differ by 4096+64 bytes
I think the gap between a and b does not really matter. After leaving only one gap between b and c I've got the following results on Haswell:
k %
-----
1 48
2 48
3 48
4 48
5 46
6 53
7 59
8 67
9 73
10 81
11 85
12 87
13 87
...
0 86
Since Haswell is known to be free of bank conflicts, the only remaining explanation is false dependence between memory addresses (and you've found proper place in Agner Fog's microarchitecture manual explaining exactly this problem). The difference between bank conflict and false sharing is that bank conflict prevents accessing the same bank twice during the same clock cycle while false sharing prevents reading from some offset in 4K piece of memory just after you've written something to same offset (and not only during the same clock cycle but also for several clock cycles after the write).
Since your code (for k=0) writes to any offset just after doing two reads from the same offset and would not read from it for a very long time, this case should be considered as "best", so I placed k=0 at the end of the table. For k=1 you always read from offset that is very recently overwritten, which means false sharing and therefore performance degradation. With larger k time between write and read increases and CPU core has more chances to pass written data through all memory hierarchy (which means two address translations for read and write, updating cache data and tags and getting data from cache, data synchronization between cores, and probably many more stuff). k=12 or 24 clocks (on my CPU) is enough for every written piece of data to be ready for subsequent read operations, so starting with this value performance gets back to usual. Looks not very different from 20+ clocks on AMD (as said by @Mysticial).
TL;DR: For certain values of k, too many 4K aliasing conditions occur, which is the main cause for the bandwidth degradation. In 4K aliasing, a load is stalled unnecessarily, thereby increasing the effective load latency and stalling all later dependent instructions. This in turn results in reduced L1 bandwidth utilization. For these values of k, most 4K aliasing conditions can be eliminated by splitting the loop as follows:
for(int i=0; i<n/64; i++) {
_mm256_store_ps(z1+64*i+ 0,_mm256_add_ps(_mm256_load_ps(x1+64*i+ 0), _mm256_load_ps(y1+64*i+ 0)));
_mm256_store_ps(z1+64*i+ 8,_mm256_add_ps(_mm256_load_ps(x1+64*i+ 8), _mm256_load_ps(y1+64*i+ 8)));
}
for(int i=0; i<n/64; i++) {
_mm256_store_ps(z1+64*i+ 16,_mm256_add_ps(_mm256_load_ps(x1+64*i+16), _mm256_load_ps(y1+64*i+ 16)));
_mm256_store_ps(z1+64*i+ 24,_mm256_add_ps(_mm256_load_ps(x1+64*i+24), _mm256_load_ps(y1+64*i+ 24)));
}
for(int i=0; i<n/64; i++) {
_mm256_store_ps(z1+64*i+ 32,_mm256_add_ps(_mm256_load_ps(x1+64*i+32), _mm256_load_ps(y1+64*i+ 32)));
_mm256_store_ps(z1+64*i+ 40,_mm256_add_ps(_mm256_load_ps(x1+64*i+40), _mm256_load_ps(y1+64*i+ 40)));
}
for(int i=0; i<n/64; i++) {
_mm256_store_ps(z1+64*i+ 48,_mm256_add_ps(_mm256_load_ps(x1+64*i+48), _mm256_load_ps(y1+64*i+ 48)));
_mm256_store_ps(z1+64*i+ 56,_mm256_add_ps(_mm256_load_ps(x1+64*i+56), _mm256_load_ps(y1+64*i+ 56)));
}
This split eliminates most 4K aliasing for the cases when k is an odd positive integer (such as 1). The achieved L1 bandwidth is improved by about 50% on Haswell. There is still room for improvement, for example, by unrolling the loop and figuring out a way to not use the indexed-addressing mode for loads and stores.
However, this split doesn't eliminate 4K aliasing for even values of k. So a different split needs to be used for even values of k. However, when k is 0, optimal performance can be achieved without splitting the loop. In this case, performance is backend-bound on ports 1, 2, 3, 4, and 7 simultaneously.
There could be a penalty of a few cycles in certain cases when performing a load and store at the same time, but in this particular case, this penalty basically does not exist because there are basically no such conflicts (i.e., the addresses of concurrent loads and stores are sufficiently far apart). In addition, the total working set size fits in the L1 so there is no L1-L2 traffic beyond the first execution of the loop.
The rest of this answer includes a detailed explanation of this summary.
First, observe that the three arrays have a total size of 24KB. In addition, since you're initializing the arrays before executing the main loop, most accesses in the main loop will hit into the L1D, which is 32KB in size and 8-way associative on modern Intel processors. So we don't have to worry about misses or hardware prefetching. The most important performance event in this case is LD_BLOCKS_PARTIAL.ADDRESS_ALIAS, which occurs when a partial address comparison involving a later load results in a match with an earlier store and all of the conditions of store forwarding are satisfied, but the target locations are actually different. Intel refers to this situation as 4K aliasing or false store forwarding. The observable performance penalty of 4K aliasing depends on the surrounding code.
By measuring cycles, LD_BLOCKS_PARTIAL.ADDRESS_ALIAS and MEM_UOPS_RETIRED.ALL_LOADS, we can see that for all values of k where the achieved bandwidth is much smaller than the peak bandwidth, LD_BLOCKS_PARTIAL.ADDRESS_ALIAS and MEM_UOPS_RETIRED.ALL_LOADS are almost equal. Also for all values of k where the achieved bandwidth is close to the peak bandwidth, LD_BLOCKS_PARTIAL.ADDRESS_ALIAS is very small compared to MEM_UOPS_RETIRED.ALL_LOADS. This confirms that bandwidth degradation is occurring due to most loads suffering from 4K aliasing.
The Intel optimization manual Section 12.8 says the following:
4-KByte memory aliasing occurs when the code stores to one memory location and shortly after that it loads from a different memory location with a 4-KByte offset between them. For example, a load to linear address 0x400020 follows a store to linear address 0x401020.
The load and store have the same value for bits 5 - 11 of their addresses and the accessed byte offsets should have partial or complete overlap.
That is, there are two necessary conditions for a later load to alias with an earlier store:
- Bits 5-11 of the two linear addresses must be equal.
- The accessed locations must overlap (so that there can be some data to forward).
On processors that support AVX-512, it seems to me that a single load uop can load up to 64 bytes. So I think the range for the first condition should be 6-11 instead of 5-11.
The following listing shows the AVX-based (32-byte) sequence of memory accesses and the least significant 12 bits of their addresses for two different values of k.
======
k=0
======
load x+(0*64+0)*4 = x+0 where x is 4k aligned 0000 000|0 0000
load y+(0*64+0)*4 = y+0 where y is 4k aligned 0000 000|0 0000
store z+(0*64+0)*4 = z+0 where z is 4k aligned 0000 000|0 0000
load x+(0*64+8)*4 = x+32 where x is 4k aligned 0000 001|0 0000
load y+(0*64+8)*4 = y+32 where y is 4k aligned 0000 001|0 0000
store z+(0*64+8)*4 = z+32 where z is 4k aligned 0000 001|0 0000
load x+(0*64+16)*4 = x+64 where x is 4k aligned 0000 010|0 0000
load y+(0*64+16)*4 = y+64 where y is 4k aligned 0000 010|0 0000
store z+(0*64+16)*4= z+64 where z is 4k aligned 0000 010|0 0000
load x+(0*64+24)*4 = x+96 where x is 4k aligned 0000 011|0 0000
load y+(0*64+24)*4 = y+96 where y is 4k aligned 0000 011|0 0000
store z+(0*64+24)*4 = z+96 where z is 4k aligned 0000 011|0 0000
load x+(0*64+32)*4 = x+128 where x is 4k aligned 0000 100|0 0000
load y+(0*64+32)*4 = y+128 where y is 4k aligned 0000 100|0 0000
store z+(0*64+32)*4= z+128 where z is 4k aligned 0000 100|0 0000
.
.
.
======
k=1
======
load x+(0*64+0)*4 = x+0 where x is 4k aligned 0000 000|0 0000
load y+(0*64+0)*4 = y+0 where y is 4k+64 aligned 0000 010|0 0000
store z+(0*64+0)*4 = z+0 where z is 4k+128 aligned 0000 100|0 0000
load x+(0*64+8)*4 = x+32 where x is 4k aligned 0000 001|0 0000
load y+(0*64+8)*4 = y+32 where y is 4k+64 aligned 0000 011|0 0000
store z+(0*64+8)*4 = z+32 where z is 4k+128 aligned 0000 101|0 0000
load x+(0*64+16)*4 = x+64 where x is 4k aligned 0000 010|0 0000
load y+(0*64+16)*4 = y+64 where y is 4k+64 aligned 0000 100|0 0000
store z+(0*64+16)*4= z+64 where z is 4k+128 aligned 0000 110|0 0000
load x+(0*64+24)*4 = x+96 where x is 4k aligned 0000 011|0 0000
load y+(0*64+24)*4 = y+96 where y is 4k+64 aligned 0000 101|0 0000
store z+(0*64+24)*4 = z+96 where z is 4k+128 aligned 0000 111|0 0000
load x+(0*64+32)*4 = x+128 where x is 4k aligned 0000 100|0 0000
load y+(0*64+32)*4 = y+128 where y is 4k+64 aligned 0000 110|0 0000
store z+(0*64+32)*4= z+128 where z is 4k+128 aligned 0001 000|0 0000
.
.
.
Note that when k=0, no load seem satisfy the two conditions of 4K aliasing. On the other hand, when k=1, all loads seem satisfy the conditions. However, it's tedious to do this manually for all iterations and all values of k. So I wrote a program that basically generates the addresses of the memory accesses and calculates the total number of loads that suffered 4K aliasing for different values of k. One issue I faced was that we don't know, for any given load, the number of stores that are still in the store buffer (have not been committed yet). Therefore, I've designed the simulator so that it can use different store throughputs for different values of k, which seems to better reflect what's actually happening on a real processor. The code can be found here.
The following figure show the number of 4K aliasing cases produced by the simulator compared to the measured number using LD_BLOCKS_PARTIAL.ADDRESS_ALIAS on Haswell. I've tuned the store throughput used in the simulator for each value of k to make the two curves as similar as possible. The second figure shows the inverse store throughput (total cycles divided by total number of stores) used in the simulator and measured on Haswell. Note that the store throughput when k=0 doesn't matter because there is no 4K aliasing anyway. Since there are two loads for each store, the inverse load throughput is half of the inverse store throughput.
Obviously the amount of time each store remains in the store buffer is different on Haswell and the simulator, so I needed to use different throughputs to make the two curves similar. The simulator can be used to show how the store throughput can impact the number of 4K aliases. If the store throughput is very close to 1c/store, then the number of 4K aliasing cases would have been much smaller. 4K aliasing conditions don't result in pipeline flushes, but they may result in uop replays from the RS. In this particular case, I didn't observe any replays though.
I think I can explain these numbers if I assume that for k=1 that writes and reads cannot happen in the same clock cycle.
There is actually a penalty of a few cycles when executing a load and store at the same time, but they can only happen when the addresses of the load and store are within 64 bytes (but not equal) on Haswell or 32 bytes on Ivy Bridge and Sandy Bridge. Weird performance effects from nearby dependent stores in a pointer-chasing loop on IvyBridge. Adding an extra load speeds it up?. In this case, the addresses of all accesses are 32-byte aligned, but, on IvB, the L1 ports are all 16-byte in size, so the penalty can be incurred on Haswell and IvB. In fact, since loads and stores may take more time to retire and since there are more load buffers than store buffers, it's more likely that a later load will false-alias an earlier store. This raises the question, though, how the 4K alias penalty and the L1 access penalty interact with each other and contribute to the overall performance. Using the CYCLE_ACTIVITY.STALLS_LDM_PENDING event and the load latency performance monitoring facility MEM_TRANS_RETIRED.LOAD_LATENCY_GT_*, it seems to me that there is no observable L1 access penalty. This implies that most of the time the addresses of concurrent loads and stores don't induce the penalty. Hence, the 4K aliasing penalty is the main cause for bandwidth degradation.
I've used the following code to make measurements on Haswell. This is essentially the same code emitted by g++ -O3 -mavx.
%define SIZE 64*64*2
%define K_ 10
BITS 64
DEFAULT REL
GLOBAL main
EXTERN printf
EXTERN exit
section .data
align 4096
bufsrc1: times (SIZE+(64*K_)) db 1
bufsrc2: times (SIZE+(64*K_)) db 1
bufdest: times SIZE db 1
section .text
global _start
_start:
mov rax, 1000000
.outer:
mov rbp, SIZE/256
lea rsi, [bufsrc1]
lea rdi, [bufsrc2]
lea r13, [bufdest]
.loop:
vmovaps ymm1, [rsi]
vaddps ymm0, ymm1, [rdi]
add rsi, 256
add rdi, 256
add r13, 256
vmovaps[r13-256], ymm0
vmovaps ymm2, [rsi-224]
vaddps ymm0, ymm2, [rdi-224]
vmovaps [r13-224], ymm0
vmovaps ymm3, [rsi-192]
vaddps ymm0, ymm3, [rdi-192]
vmovaps [r13-192], ymm0
vmovaps ymm4, [rsi-160]
vaddps ymm0, ymm4, [rdi-160]
vmovaps [r13-160], ymm0
vmovaps ymm5, [rsi-128]
vaddps ymm0, ymm5, [rdi-128]
vmovaps [r13-128], ymm0
vmovaps ymm6, [rsi-96]
vaddps ymm0, ymm6, [rdi-96]
vmovaps [r13-96], ymm0
vmovaps ymm7, [rsi-64]
vaddps ymm0, ymm7, [rdi-64]
vmovaps [r13-64], ymm0
vmovaps ymm1, [rsi-32]
vaddps ymm0, ymm1, [rdi-32]
vmovaps [r13-32], ymm0
dec rbp
jg .loop
dec rax
jg .outer
xor edi,edi
mov eax,231
syscall