Fortran 中具有对称性的高阶张量的 BLAS

Question

如果我有张量收缩 A[a,b] * B[b,c,d] = C[a,c,d] 它具有 B[b,c,d] = B[b,d,c] 和 C[a,c,d] = C[a,d,c] 属性，如何设置 BLAS 以利用这种对称性？

这里假设爱因斯坦求和符号，即重复索引表示求和。

我可以尝试将 flat/reshape 张量 B 转换为低维数组，但似乎 flat/reshape 张量也需要时间，至少在 Fortran 中是这样。

Answer 1

矩阵运算 C_{acd} = A_{ab} 。 B_{bcd} 可以编程方式编写为 matrix * vector 操作的双循环（为清楚起见，使用 matmul；根据需要替换为 BLAS）：

n = size(B,3) ! = size(B,2)
do d=1,n
  do c=1,n
    C(:,c,d) = matmul(A(:,:), B(:,c,d))
  enddo
enddo

因为“C[a,d,c]=C[a,c,d]”，所以matmul的方环可以用matmul的三角环和刚复制的三角环代替，如：

n = size(B,3) ! = size(B,2)
do d=1,n
  do c=1,d
    C(:,c,d) = matmul(A(:,:), B(:,c,d))
  enddo
  
  do c=d+1,n
    C(:,c,d) = C(:,d,c)
  enddo
enddo

这利用对称性来减少 BLAS 操作的数量，提高性能，但必须进行大量 matrix * vector 乘法而不是一次大的 matrix * matrix 乘法会降低性能。这种方法会整体提高还是降低性能？找出答案的最好方法可能是尝试看看。

BLAS with symmetry in higher order tensor in Fortran