ove*_*yan 5 fortran simd gfortran avx2 intel-fortran
当 gfortran 对方程（例如 x = y*z ）进行向量化时，编译器将仅使用向量寄存器（例如 YMM）和向量操作码，当且仅当编译器知道所有数据 (x,y,z) 都可以是加载到所述向量寄存器上。然而，如果编译器不知道是否所有数据都可以加载到寄存器中，它将生成冗余代码，以便可以操作剩余的数据。
此冗余代码通常不会使用最高可用的 SIMD 向量寄存器，并且会使可执行文件膨胀。我还没有进行广泛的测试，但常见的推理使我相信它可能会导致 ICACHE 缺失或阻止函数/子例程内联优化。
我想删除这个冗余代码，因为我知道我的所有数据 (x,y,z) 将完美地适合 YMM 向量寄存器。
我之前曾在这里更详细地描述过我的观察结果：https: //www.cfd-online.com/Forums/main/231759-fortran- assembly-how-remove-redundant-non-vectorized-code.html
不过，我想在这里举一些简单的例子来展示我的观点。
对于以下操作：x(:) = y*z其中x，y，z都是可分配数组，编译器不知道该操作是否可以完全矢量化，因为它不知道数组的长度是多少。如果 x,y,z 是整数数组并且长度是 8 的倍数，我们可以推断整个操作可以简单地使用 AVX2 向量化操作来完成。然而，由于数组长度在编译时未知，编译器必须生成可用于任意长度数组的冗余代码。
上帝螺栓链接：https://gcc.godbolt.org/z/1fdzK8
这是为该x(:) = y*z部分生成的程序集：
.L7:
vmovdqu ymm1, YMMWORD PTR [r13+0+rax]
vpmulld ymm0, ymm1, YMMWORD PTR [r12+rax]
vmovdqu YMMWORD PTR [r14+rax], ymm0
add rax, 32
cmp rdx, rax
jne .L7
mov rdx, r15
and rdx, -8
lea rax, [rdx+1]
cmp rdx, r15
je .L22
vzeroupper
.L6:
lea rdx, [rax-1]
mov esi, DWORD PTR [r13+0+rdx*4]
imul esi, DWORD PTR [r12+rdx*4]
mov DWORD PTR [r14+rdx*4], esi
lea rdx, [rax+1]
cmp r15, rdx
jl .L10
mov esi, DWORD PTR [r13+0+rax*4]
imul esi, DWORD PTR [r12+rax*4]
mov DWORD PTR [r14+rax*4], esi
lea rsi, [rax+2]
cmp r15, rsi
jl .L10
mov edi, DWORD PTR [r13+0+rdx*4]
imul edi, DWORD PTR [r12+rdx*4]
mov DWORD PTR [r14+rdx*4], edi
lea rdx, [rax+3]
cmp r15, rdx
jl .L10
mov edi, DWORD PTR [r13+0+rsi*4]
imul edi, DWORD PTR [r12+rsi*4]
mov DWORD PTR [r14+rsi*4], edi
lea rsi, [rax+4]
cmp r15, rsi
jl .L10
mov edi, DWORD PTR [r13+0+rdx*4]
imul edi, DWORD PTR [r12+rdx*4]
mov DWORD PTR [r14+rdx*4], edi
lea rdx, [rax+5]
cmp r15, rdx
jl .L10
mov edi, DWORD PTR [r13+0+rsi*4]
add rax, 6
imul edi, DWORD PTR [r12+rsi*4]
mov DWORD PTR [r14+rsi*4], edi
cmp r15, rax
jl .L10
mov eax, DWORD PTR [r12+rdx*4]
imul eax, DWORD PTR [r13+0+rdx*4]
mov DWORD PTR [r14+rdx*4], eax
我们可以清楚地看到上面代码的.L7部分包含矢量化代码，.L6部分包含冗余代码，以确保任意长度的数组都可以使用相同的代码。
虽然这是一个很好的功能，但我想可靠地关闭它并指示编译器仅在 .L7 部分生成矢量化代码。我还没有进行广泛的性能测试，但是我仍然想删除 .L6 部分。
我能够使英特尔的 Fortran 编译器仅使用编译器指令生成矢量化代码部分，例如：
!DIR$ ASSUME_ALIGNED X: 32
!DIR$ ASSUME (MOD(N,8) .EQ. 0)
do i=1,n
x(i) = y(i)*z(i)
enddo
上帝螺栓链接：https://gcc.godbolt.org/z/PPf1zE
生成的程序集是：
..B1.15: # Preds ..B1.15 ..B1.14
vmovdqu ymm0, YMMWORD PTR [r13+r12*4] #22.12
vpmulld ymm1, ymm0, YMMWORD PTR [r14+r12*4] #22.5
vmovdqu YMMWORD PTR [rax+r12*4], ymm1 #22.5
add r12, 8 #21.1
cmp r12, rcx #21.1
jb ..B1.15 # Prob 82% #21.1
我还能够通过使用一些编译器巫术使 gfortran 编译器仅生成矢量化代码部分。我们只是用来do i=1,AND(N,-8)指示编译器循环以 8 的倍数结束。从数学上讲，它发生的原因是AND(N,-8) == N - MOD(N,8)
do i=1,AND(n,-8)
x(i) = y(i)*z(i)
end do
上帝螺栓链接： https: //gcc.godbolt.org/z/MW7cPG
生成的程序集是：
.L19:
mov rcx, QWORD PTR [rsp+16]
vmovdqu ymm4, YMMWORD PTR [rcx+rax]
mov rcx, QWORD PTR [rsp+8]
vpmulld ymm0, ymm4, YMMWORD PTR [r15+rax]
vmovdqu YMMWORD PTR [rcx+rax], ymm0
add rax, 32
cmp rax, rdx
jne .L19
然而，这对于 gfortran 来说并不是一个可靠的方法，并且无法消除我的一些测试用例中的冗余膨胀。AFAIK gfortran 没有任何像 这样的编译器指令!DIR$ ASSUME (MOD(N,8) .EQ. 0)，所以这个方法对我来说似乎非常不可靠。
我希望以可靠且一致的方式进行此优化。
额外数据：
对看似复杂的代码执行了相同的操作：
vols(:) = 0.5*sqrt((x1*(y2-y3))**2.0 + (x2*(y3-y1))**2.0 + (x3*(y1-y2))**2.0)
上帝螺栓链接： https: //gcc.godbolt.org/z/oxcWzd
在上面的示例中，生成的膨胀相当高。
上述优化适用于 Intel 的 Fortran 编译器。
上帝螺栓链接：https://gcc.godbolt.org/z/Y5orh1
上述优化适用于 gfortran 编译器。
从哲学上来说，Fortran 不希望用户修改低级优化，因为该语言应该能够实现这一点。在您的情况下，可以通过将体积定义为函数来实现完全矢量化elemental：
elemental real function triVolume(x1,x2,x3,y1,y2,y3) result(vol)
real, intent(in) :: x1,x2,x3,y1,y2,y3
vol = 0.5*sqrt((x1*(y2-y3))**2 + (x2*(y3-y1))**2 + (x3*(y1-y2))**2)
end function triVolume
forall然后使用or循环它do concurrent：
! Calculate volume of a triangle using cross product
do concurrent (i=1:n)
vols(i) = triVolume(x1(i),x2(i),x3(i),y1(i),y2(i),y3(i))
end do
这会产生干净的矢量化代码。
请注意，do concurrent此处启用了 CPU 矢量化，但在未来的某些版本中，可能会隐藏有趣的内容，例如卸载到 GPU。