Look into enabling LAPACK for Fortran matmul intrinsic

[dhermes]

See some discussion on StackOverflow.

You don't want all MATMULs to be dgemm, it is not profitable for very small matrices. gfortran does what you want

From the gfortran docs for -fexternal-blas:

generate calls to BLAS functions for some matrix operations like MATMUL ... if the size of the matrices involved is larger than a given limit (see -fblas-matmul-limit) ... The BLAS library will have to be specified at link time.

Options, e.g.:

-fblas-matmul-limit
-fblas-matmul-limit  =30
-finline-matmul-limit=30

From -fblas-matmul-limit docs

If the matrices involved are not square, the size comparison is performed using the geometric mean of the dimensions of the argument and result matrices.

Possible BLAS implementations to install (on Debian):

$ apt install libblas3           # reference BLAS
$ apt install libopenblas-base   # openBLAS
$ apt install libatlas3-base     # ATLAS

[dhermes]

For the "how" of this, I used:

! main.f90
program main

  implicit none

  integer, parameter :: dp = kind(0.d0)

  real(dp) :: left(4, 50)
  real(dp) :: right(50, 2)
  real(dp) :: product_(4, 2)
  integer:: i, j

  do i = 1, 4
     do j = 1, 50
        left(i, j) = i + j
        if (i < 3) then
           right(j, i) = i * i - j
        end if
     end do
  end do

  product_ = matmul(left, right)

  do i = 1, 4
     print *, product_(i, :)
  end do

end program main

compiled with 4 (or 6, depending how you count) sets of options:

$ gfortran -o main1 main.f90
$ gfortran -o main2 main.f90 -llapack
$ gfortran -o main3 main.f90 -lblas
$ gfortran -o main4 main.f90 -fexternal-blas -lblas
$ # Won't compile: gfortran -o main5 main.f90 -fexternal-blas
$ # Won't compile: gfortran -o main6 main.f90 -fexternal-blas -llapack

and compared binaries produced:

$ diff -s main1 main2
Files main1 and main2 are identical
$ diff -s main1 main3
Files main1 and main3 are identical
$ diff -s main1 main4
Binary files main1 and main4 differ