{"id":3937,"date":"2021-01-14T13:27:49","date_gmt":"2021-01-14T18:27:49","guid":{"rendered":"https:\/\/blogs.mathworks.com\/loren\/?p=3937"},"modified":"2021-03-16T10:21:55","modified_gmt":"2021-03-16T14:21:55","slug":"paged-matrix-functions","status":"publish","type":"post","link":"https:\/\/blogs.mathworks.com\/loren\/2021\/01\/14\/paged-matrix-functions\/","title":{"rendered":"Paged Matrix Functions"},"content":{"rendered":"
\r\n
Today's guest blogger is Mary Fenelon, who is the product marketing manager for Optimization and Math here at MathWorks. In today's post she describes how using the new paged matrix functions can simplify your code and improve performance.<\/div>\r\n
\r\n
Table of Contents<\/span><\/div>\r\n
Paged Matrix Functions\r\n<\/a> Easier\r\n<\/a> Faster\r\n<\/a> More\r\n<\/a>Helper Functions<\/a><\/div>\r\n<\/div>\r\n
A batched, page-wise matrix multiply function pagemtimes<\/a>, along with page-wise transpose pagetranspose<\/a>, and complex conjugate transpose pagectranspose<\/a> functions, were added to MATLAB\u00ae in R2020b. The new functions make operations where an N-D array is a container for dense matrices easier to write and faster to run.<\/div>\r\n
For example, if A and B are 3-D arrays and you want to compute the array C to be the collection of matrices formed by multiplying the matrices contained in the first two dimensions of A and B, you could do this with the statement C = pagemtimes(A,B).<\/span><\/div>\r\n
<\/div>\r\n
\"\"<\/div>\r\n
The alternatives besides for-loops to pagemtimes<\/span> are dlmtimes<\/span> for dlarray<\/span> in the Deep Learning Toolbox and pagefun<\/span> for gpuArray<\/span> and distributed arrays in the Parallel Computing Toolbox. These are still available, but now you can use pagemtimes for regular arrays as well as dlarray<\/span> and gpuArray<\/span>.<\/div>\r\n
Sometimes you only need to transpose a collection of matrices. Now you can use pagetranspose<\/span> instead of permute<\/span>.<\/div>\r\n

Easier<\/h3>\r\n
In this example, the matrices to be multiplied are stored in the first two dimensions of a 4-D array<\/div>\r\n
\r\n
\r\n
A = rand(4,5,6,7);<\/span><\/div>\r\n<\/div>\r\n
\r\n
B = rand(2,5,6,7);<\/span><\/div>\r\n<\/div>\r\n<\/div>\r\n
This doubly nested for-loop computes the matrix-matrix products.<\/div>\r\n
\r\n
\r\n
C1 = zeros(size(A,1),size(B,1),6,7);<\/span><\/div>\r\n<\/div>\r\n
\r\n
for <\/span>p = 1:size(A,3)<\/span><\/div>\r\n<\/div>\r\n
\r\n
for <\/span>q = 1:size(A,4)<\/span><\/div>\r\n<\/div>\r\n
\r\n
C1(:,:,p,q) = A(:,:,p,q)*B(:,:,p,q).';<\/span><\/div>\r\n<\/div>\r\n
\r\n
end<\/span><\/span><\/div>\r\n<\/div>\r\n
\r\n
end<\/span><\/span><\/div>\r\n<\/div>\r\n<\/div>\r\n
Using pagemtimes<\/span>, that code can be simplified to a single function call:<\/div>\r\n
\r\n
\r\n
C2 = pagemtimes(A,'none'<\/span>,B,'transpose'<\/span>);<\/span><\/div>\r\n<\/div>\r\n<\/div>\r\n
That is certainly easier. The string options specify that the page-wise transpose of the second array should be used in the multiplication. pagemtimes<\/span> implicitly expands the dimensions beyond three to multiply all combinations of the paged matrices.<\/div>\r\n

Faster<\/span><\/h3>\r\n
Now that I've made the case for easier, how about for faster? Does the speed-up vary with the size of the N-D arrays? The loop below runs through a set of matrices with varying row sizes of the first matrix but fixed dimensions otherwise. I plot the times and the speed-ups.<\/div>\r\n
\r\n
\r\n
rng('default'<\/span>); % for reproducibility<\/span><\/span><\/div>\r\n<\/div>\r\n
\r\n
n = 100; <\/span><\/div>\r\n<\/div>\r\n
\r\n
k = 2;<\/span><\/div>\r\n<\/div>\r\n
\r\n
p = 10000; <\/span><\/div>\r\n<\/div>\r\n
\r\n
mset = 10:10:100;<\/span><\/div>\r\n<\/div>\r\n
\r\n
mSize = numel(mset);<\/span><\/div>\r\n<\/div>\r\n
\r\n
tLoop = zeros(mSize,1);<\/span><\/div>\r\n<\/div>\r\n
\r\n
tPage = zeros(mSize,1);<\/span><\/div>\r\n<\/div>\r\n
\r\n
<\/div>\r\n<\/div>\r\n
\r\n
for <\/span>i = 1:numel(mset) <\/span><\/div>\r\n<\/div>\r\n
\r\n
m = mset(i);<\/span><\/div>\r\n<\/div>\r\n
\r\n
A = rand(m,n,p);<\/span><\/div>\r\n<\/div>\r\n
\r\n
B = rand(n,k);<\/span><\/div>\r\n<\/div>\r\n
\r\n
<\/div>\r\n<\/div>\r\n
\r\n
CPage = pagemtimes(A,B);<\/span><\/div>\r\n<\/div>\r\n
\r\n
f = @()pagemtimes(A,B);<\/span><\/div>\r\n<\/div>\r\n
\r\n
tPage(i) = timeit(f);<\/span><\/div>\r\n<\/div>\r\n
\r\n
<\/div>\r\n<\/div>\r\n
\r\n
CLoop = loopmtimes(A,B,m,k,p);<\/span><\/div>\r\n<\/div>\r\n
\r\n
f = @()loopmtimes(A,B,m,k,p);<\/span><\/div>\r\n<\/div>\r\n
\r\n
tLoop(i) = timeit(f); <\/span><\/div>\r\n<\/div>\r\n
\r\n
end<\/span><\/span><\/div>\r\n<\/div>\r\n
\r\n
plotTimes(mset,n,p,k,tPage,tLoop);<\/span><\/div>\r\n
\r\n
\r\n
<\/div>\r\n<\/div>\r\n<\/div>\r\n<\/div>\r\n<\/div>\r\n
The Times plot shows that the for-loop takes more time than pagemtimes<\/span>. The Speed-ups plot shows that the benefit of using pagemtimes<\/span> is largest when there are many small matrices but even for the larger values of m<\/span>, the speed-up is still substantial. The improvement for smaller matrices is more significant because the overhead that is being eliminated makes up a larger percentage of the overall time.<\/div>\r\n
I like to use timeit<\/span><\/a> instead of tic<\/span>\/toc<\/span> for timings (see a previous post<\/a> on timeit<\/span>) even though it requires putting the code to be timed into a function. timeit<\/span> runs the code a number of times and takes the median measurement. It's the approach recommended in Measure the Performance of Your Code<\/a>. Even so, with code that takes a very small amount of time to run, such as when the matrix sizes are small, timeit<\/span> will have trouble returning accurate timings. The recommendation in that case is to run the code multiple times in a loop.<\/div>\r\n
If you try different values for the dimensions of A and B, you'll see different speed-ups. For example, when the arrays are large and there are only a few pages, there won't be much since most of the work is in the actual matrix multiplication. You may also see some irregularity due to the difficulty of measuring times for small matrices and the effects of cache and threading for large matrices.<\/div>\r\n

More<\/h3>\r\n
Paged functions can be written for other matrix functions. Which ones should we do first and in what applications will they be used? Let us know here<\/a>.<\/div>\r\n

Helper Functions<\/h2>\r\n
\r\n
\r\n
function <\/span>Cloop = loopmtimes(A,B,m,k,p)<\/span><\/div>\r\n<\/div>\r\n
\r\n
Cloop = zeros(m,k,p);<\/span><\/div>\r\n<\/div>\r\n
\r\n
for <\/span>i = 1:p<\/span><\/div>\r\n<\/div>\r\n
\r\n
Cloop(:,:,i) = A(:,:,i)*B;<\/span><\/div>\r\n<\/div>\r\n
\r\n
end<\/span><\/span><\/div>\r\n<\/div>\r\n
\r\n
end<\/span><\/span><\/div>\r\n<\/div>\r\n<\/div>\r\n
<\/div>\r\n
\r\n
\r\n
function <\/span>plotTimes(mset,n,p,k,tPage,tLoop)<\/span><\/div>\r\n<\/div>\r\n
\r\n
figure('Position'<\/span>,[100,100,1200,450])<\/span><\/div>\r\n<\/div>\r\n
\r\n
t = tiledlayout(1,2);<\/span><\/div>\r\n<\/div>\r\n
\r\n
<\/div>\r\n<\/div>\r\n
\r\n
nexttile<\/span><\/div>\r\n<\/div>\r\n
\r\n
plot(mset,tPage,'*'<\/span>);<\/span><\/div>\r\n<\/div>\r\n
\r\n
hold on<\/span><\/span><\/div>\r\n<\/div>\r\n
\r\n
plot(mset,tLoop,'o'<\/span>);<\/span><\/div>\r\n<\/div>\r\n
\r\n
legend({'pagemtimes'<\/span>,'for-loop'<\/span>},'Location'<\/span>,'northwest'<\/span>);<\/span><\/div>\r\n<\/div>\r\n
\r\n
title('Times'<\/span>)<\/span><\/div>\r\n<\/div>\r\n
\r\n
xlabel('m'<\/span>)<\/span><\/div>\r\n<\/div>\r\n
\r\n
ylabel('Seconds'<\/span>)<\/span><\/div>\r\n<\/div>\r\n
\r\n
hold off<\/span><\/span><\/div>\r\n<\/div>\r\n
\r\n
<\/div>\r\n<\/div>\r\n
\r\n
nexttile<\/span><\/div>\r\n<\/div>\r\n
\r\n
plot(mset,tLoop.\/tPage)<\/span><\/div>\r\n<\/div>\r\n
\r\n
title('Speed-ups of pagemtimes over for-loop'<\/span>)<\/span><\/div>\r\n<\/div>\r\n
\r\n
xlabel('m'<\/span>)<\/span><\/div>\r\n<\/div>\r\n
\r\n
ylim([0,inf])<\/span><\/div>\r\n<\/div>\r\n
\r\n
<\/div>\r\n<\/div>\r\n
\r\n
titleStr = sprintf('Multiply A (m-by-%d-by-%d) by B (%d-by-%d)'<\/span>,n,p,n,k);<\/span><\/div>\r\n<\/div>\r\n
\r\n
title(t,titleStr)<\/span><\/div>\r\n<\/div>\r\n
\r\n
end<\/span><\/span><\/div>\r\n<\/div>\r\n<\/div>\r\n
Copyright 2021 The MathWorks, Inc.<\/span><\/div>\r\n<\/div>\r\n