3D matrix multiplication

Hi,

I have two 3D matrices that I need to multiply in a specific way. I need to do for k=1:Z M(:,:,k)=A(:,:,k)*B(:,:,k) end

but I do not want to use for loop because it makes my code run slower and I need to this multiplication for a unsteady flow profile calculation (close to million time steps). Is there a more efficient way to do it?

[Steven Lord]

"but I do not want to use for loop because it makes my code run slower "

I know this question was originally asked several years ago, but I know others have this same belief. Before avoiding the for loop, consider if you think that the loop makes your code run slower or whether you have measured and confirmed that it makes your code run slower. If you haven't checked, check before seeking out an alternate solution.

We have improved the performance of MATLAB in many ways over the past several years and we continue to do so now, so that "slow for loop" may not be so slow any more.

[James Tursa]

In this particular case it really isn't the for-loop per-se that is slowing things down. It is the data copying of extracting the individual 2D slices that is the culprit (which can be avoided in a mex routine).

I'm still hoping that someday TMW will choose to do "implicit 2D slice expansion" for the * matrix multiply operator and \ and / slash operators in the same sense that they now do implicit scalar expansion for the .* and ./ etc operators.

[Steven Lord]

We don't do implicit slice expansion, but as of release R2020b we provide the pagemtimes function for this purpose.

[James Tursa]

Finally! Guess I can now stop worrying about updating MTIMESX.

Answer 1

[James Tursa]

M = mtimesx(A,B);

You can find MTIMESX here:

http://www.mathworks.com/matlabcentral/fileexchange/25977-mtimesx-fast-matrix-multiply-with-multi-dimensional-support

[James Tursa]

P.S. If your 2D slices are small (4x4 or less), be sure to use 'SPEED' or 'SPEEDOMP' mode. This will cause inline code to be used that will run faster than the default BLAS library routine calls.

Answer 2

[Daniel Frisch]

Matrix multiplication can also be expressed using native Matlab code (times and sum):

% A : (a x c x Z)
% B : (c x b x Z)
Ap = permute(A,[2,1,4,3]); % (c x a x 1 x Z)
Bp = permute(B,[1,4,2,3]); % (c x 1 x b x Z)
M = Ap .* Bp;              % (c x a x b x Z)
M = sum(M,1);              % (1 x a x b x Z)
M = permute(M,[2,3,4,1]);  % (a x b x Z)

[James Tursa]

Note that the permute( ) operations will involve a deep data copy, so this method will run slower than the mex routines.

Answer 3

[Walter Roberson]

There is no built-in MATLAB support for 3D multiplications. The program James refers to is probably a good choice.

By the way, especially in the releases of the last few years, "for" loops are sometimes faster than vectorization, especially for large matricies.

Answer 4

[Kamuran]

to use this file I need to install a C compiler ? If so can you suggest me one, Because I never used C :) .

Thank you

[Walter Roberson]

Appropriate C compilers depend upon your OS and MATLAB release. You can start at the following, but be sure to look for information about your particular release, as the list of supported compilers is different between releases. Also note that 64 bit Windows has different compilers than 32 bit Windows.

http://www.mathworks.com/support/compilers/R2011a/win32.html

[James Tursa]

Type the following at the MATLAB prompt:

mex -setup

Then press Enter and see what shows up. If you have a 32-bit system then hopefully a compiler will already be there. If you have a 64-bit system then you will have to install one yourself unfortunately.

[Walter Roberson]

I checked back... we do not have any information about the OS or MATLAB version the poster is using.

Answer 5

[Mr M.]

I got the following error message:

Error using mtimesx_build (line 120) Unable to compile mtimesx.c

Error in mtimesx (line 271) mtimesx_build;

[Mr M.]

However XCode is installed and licence is accepted (mex -setup)

[James Tursa]

MTIMESX is way overdue for an update (TMW changed mex procedures some time ago). In the meantime, try this or some variation of this modified for your system:

lib_blas = fullfile(matlabroot,'extern','lib',computer('arch'),'microsoft','libmwblas.lib');
mex('mtimesx.c',lib_blas,'-largeArrayDims')