I have a matrix
A which is of size
(n2, n1) and I am multiplying it by a matrix,
B, of size
(n1, n0). I have identified this single matrix multiplication as the bottleneck in my
Fortran code. Out of ~2000 lines of code, this single line takes about 77% of the runtime.
A is a double precision matrix with floating point values.
B is, currently, a double precision matrix containing only values
0.0. I can easily make this integer, or even binary, but I was using it as
real so that I could preserver precision in
What is a better way to perform this matrix multiplication to cut down on runtime?
Before anyone suggests it, I am using
DGEMM and compiling with
The largest data I have implemented this program on so far,
N = 5000, results in
n2 = 1668,
n1 = 1701, and
n0 = 1631. This algorithm was implemented in
Matlab and has shorter runtime. Matlab version is about 2.5 seconds, while this fortran program is about 7 seconds. Since this single matrix multiplication is so large, I'm thinking that Matlab is doing something interesting with the variable types.
I have compiled this with
MKL and am current linking against
-lblas and using
-fexternal-blas, relying on
matmul to perform the underlying
BLAS routines. The result of
ldd on my binary executable is:
linux-vdso.so.1 => (0x00002aaaaaacb000) liblapack.so.3 => /usr/lib64/atlas/liblapack.so.3 (0x00002aaaaaccd000) libblas.so.3 => /usr/lib64/libblas.so.3 (0x00002aaaab4f0000) libgfortran.so.3 => /usr/lib64/libgfortran.so.3 (0x00002aaaab747000) libm.so.6 => /lib64/libm.so.6 (0x00002aaaaba39000) libgcc_s.so.1 => /lib64/libgcc_s.so.1 (0x0000003f78c00000) libc.so.6 => /lib64/libc.so.6 (0x00002aaaabcbe000) libf77blas.so.3 => /usr/lib64/atlas/libf77blas.so.3 (0x00002aaaac052000) libcblas.so.3 => /usr/lib64/atlas/libcblas.so.3 (0x00002aaaac272000) /lib64/ld-linux-x86-64.so.2 (0x00002aaaaaaab000) libatlas.so.3 => /usr/lib64/atlas/libatlas.so.3 (0x00002aaaac492000) libpthread.so.0 => /lib64/libpthread.so.0 (0x00002aaaacaee000)
B is structured in the way that it has zeros and ones. The lower left portion (not truly lower triangular) has ones and the upper right portion (not triangular) is zeros.
It appears that the Matlab code is treating the
B matrix as