# Efficient method to multiply floating point matrix with binary matrix and get double precision results

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 1.0 and 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 matmul(A,B).

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 -O3 and -mavx for gfortran, and -O3 with -xhost on ifort.

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 ifort using 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)


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 logical.

• The right way to link to the MKL these days with ifort is to use -mkl. -lblas may find a system library that is not optimized. – Bill Barth Jul 28 '15 at 2:38
• I only use -lblas when compiling with gfortran. I use specific Intel compilation flags when linking with ifort – drjrm3 Jul 28 '15 at 2:40
• a) you can use those same options with gfortran to get the MKL, and b) you might try -mkl with ifort just to be sure. – Bill Barth Jul 28 '15 at 2:43
• I don't see how using real for the integer matrix would give you better precision. You can check for your BLAS using ldd on the binary you compile. – AlexE Jul 28 '15 at 11:47
• Your latest edit shows you are actually not linking against the MKL, but rather against ATLAS (which I'd personally drop any time in favour of OpenBLAS). You may try to leave out -lblas in linking. – AlexE Jul 28 '15 at 13:26

What is the percentage of nonzero entries in $B$? If a high percentage of the entries are 0's, then you might well be better off treating $B$ as a sparse matrix in the multiplication.