OpenBLAS自动并行DGEMM效率问题：大公共维度性能下降咨询

核心原因：并行开销大于收益 #

你遇到的不是OpenBLAS的局限，而是并行计算的通用问题——当计算任务的粒度太小时，线程调度、数据同步的开销会超过并行带来的性能提升。

具体到你的测试案例：

• 6×100000 与 100000×6 的矩阵乘法，最终结果是6×6的矩阵。每个线程分到的计算量极小，线程创建、切换，以及数据在缓存间的同步开销，完全抵消甚至超过了多线程带来的计算加速。
• 而1000×1000的矩阵乘法计算量足够大，并行的收益能覆盖开销，所以表现出性能提升。

OpenBLAS的并行策略 #

OpenBLAS的DGEMM会根据矩阵尺寸自动决定是否启用并行以及分配线程数。对于这种窄长型矩阵的乘法，它的自动并行判断可能没有最优适配，但这不是Bug，是并行计算的固有特性。

解决方案 #

1. 手动控制线程数 #

针对小尺寸矩阵的计算，直接强制使用单线程，避免并行开销：
在代码中根据矩阵尺寸判断，比如当结果矩阵的行数+列数小于某个阈值时，调用 openblas_set_num_threads(1)，大尺寸计算再启用多线程。

也可以通过环境变量全局限制小任务的线程数：

export OPENBLAS_THREAD_LIMIT=1

2. 转置矩阵优化计算模式 #

对于6×100000 × 100000×6的乘法，先转置B矩阵得到B^{T（6×100000），此时C的每个元素是A的行与B}T的行的点积，这种计算模式能更好利用CPU缓存，同时减少并行调度的无效开销。

3. 是否需要转用MKL？ #

没必要。MKL的自动并行策略可能更精细，但同样会遇到小任务并行开销的问题。除非你的场景中绝大多数是大矩阵计算，否则转用MKL解决不了这个特定问题，手动控制线程数或优化计算模式的性价比更高。

4. 手动实现并行？ #

完全没必要。OpenBLAS和MKL的DGEMM已经是经过高度优化的工业级实现，手动并行不仅很难达到它们的性能水平，还容易引入线程安全、缓存命中率低等问题。

原测试内容 #

Fortran测试代码 #

Program test

    Use, Intrinsic :: iso_fortran_env, Only :  wp => real64, li => int64
    integer, parameter :: dp = selected_real_kind(15, 307)
    
    Real( dp ), Dimension( :, : ), Allocatable :: a
    Real( dp ), Dimension( :, : ), Allocatable :: b
    Real( dp ), Dimension( :, : ), Allocatable :: c

    Integer :: na, nb, nc, m_iter
    Integer :: numthreads
    Integer( li ) :: start, finish, rate
    real(dp) :: sum_time1, alpha, beta

  
    Write( *, * ) 'number of threads'
    Read( *, * ) numthreads

    Write( *, * ) 'na, nb, nc?'
    Read( *, * ) na, nb, nc
    Allocate( a ( 1:na, 1:nb ) ) 
    Allocate( b ( 1:nb, 1:nc ) )  
    Allocate( c( 1:na, 1:nc ) )

  
    Call Random_number( a )
    Call Random_number( b )
    c  = 0.0_dp
 
    alpha = 1.0_dp
    beta = 0.0_dp
    m_iter = 10
  
    write (*,*) 'm_iter average', m_iter 
    !call omp_set_num_threads(1)
    call openblas_set_num_threads(1)
    sum_time1 = 0.0  
    do m = 1, m_iter     
      Call System_clock( start, rate )
      Call dgemm( 'N', 'N', na, nc, nb, alpha, a , na,  b , nb, beta, c, nc )
      Call System_clock( finish, rate )
      sum_time1 = sum_time1 + Real( finish - start, dp ) / rate  
    end do  
    Write( *, * ) 'Time for dgemm-1', sum_time1 / m_iter
      

   ! call omp_set_num_threads(numthreads)  
  !  call mkl_set_num_threads(numthreads)
    call openblas_set_num_threads(numthreads)
    sum_time1 = 0.0  
    do m = 1, m_iter     
      Call System_clock( start, rate )
      Call dgemm( 'N', 'N', na, nc, nb, alpha, a , na,  b , nb, beta, c, nc )
      Call System_clock( finish, rate )
      sum_time1 = sum_time1 + Real( finish - start, dp ) / rate  
    end do  
    Write( *, * ) 'Time for dgemm-auto parallel', sum_time1 / m_iter
 
     deallocate(a, b, c) 

  End

编译命令 #

gfortran-11 -O3 -fopenmp test.f90 -lopenblas

运行结果 #

number of threads

2

           2

 na, nb, nc?

6 100000 6

 m_iter average          10
Time for dgemm-1  0.15203639000000002
 Time for dgemm-auto parallel  0.20141369999999997

内容的提问来源于stack exchange，提问作者user235862