BLAS 子例程 dgemm、dgemv 和 ddot 不适用于标量？

Question

我有一个 Fortran 子例程，它使用 BLAS 的子例程 dgemm、dgemv 和 ddot，它们计算矩阵 * 矩阵、矩阵 * 向量和向量 * 向量。我有 m * m 个矩阵和 m * 1 个向量。在某些情况下，m=1。看来这些子例程在这些情况下不能很好地工作。他们没有给出错误，但结果似乎存在一些数值不稳定。所以我必须写一些类似的东西：

if(m>1) then 
  vtuni(i,t) = yt(i,t) - ct(i,t) - ddot(m, zt(i,1:m,(t-1)*tvar(3)+1), 1, arec, 1)
else 
   vtuni(i,t) = yt(i,t) - ct(i,t) - zt(i,1,(t-1)*tvar(3)+1)*arec(1)

所以我真正的问题是，当 m=1 时，那些 BLAS 子例程无法正常工作，我是否正确，或者我的代码中是否有问题？编译器会影响这个吗？我正在使用 gfortran。

Answer 1

BLAS 例程对于大小为 1 的对象应该能正确运行。我不认为它可以依赖于编译器，但它可能依赖于您所依赖的 BLAS 的实现（尽管我认为它是实现的错误）。参考（阅读：没有目标优化）BLAS 的实现，可以在 网络库, ，很好地处理这种情况。

我已经对大小为 1 的数组和较大数组的大小为 1 的切片（如您自己的代码中所示）进行了一些测试，它们都工作正常：

 $ cat a.f90 
 implicit none
 double precision :: u(1), v(1)
 double precision, external :: ddot
 u(:) = 2
 v(:) = 3
 print *, ddot(1, u, 1, v, 1)
 end
 $ gfortran a.f90 -lblas && ./a.out
  6.0000000000000000     

 $ cat b.f90                       
 implicit none
 double precision, allocatable :: u(:,:,:), v(:)
 double precision, external :: ddot
 integer :: i, j
 allocate(u(3,1,3),v(1))
 u(:,:,:) = 2
 v(:) = 3
 i = 2
 j = 2
 print *, ddot(1, u(i,1:1,j), 1, v, 1)
 end
 $ gfortran b.f90 -lblas && ./a.out
  6.0000000000000000     

我会考虑进一步调试这个问题：

• 检查您的 ddot 定义是正确的
• 将参考 BLAS 替换为优化后的 BLAS，以检查它是否更改了任何内容（您可以在 ddot.f 我之前在答案中链接到的文件）