Efficient way to multiply a vector of size (Nx,1) with a matrix of size (Nx+1,Nx+1)

Question

Matthew Kehoe 2021-7-28

0
链接

此问题的直接链接

https://ww2.mathworks.cn/matlabcentral/answers/887977-efficient-way-to-multiply-a-vector-of-size-nx-1-with-a-matrix-of-size-nx-1-nx-1

评论： the cyclist 2021-7-29

My code repeatedly computes a vector of size (Nx,1) by a matrix of size (Nx+1,Nx+1):

clc; clear all;
Nx = 32;
a = 1.5;
b = 2;
gamma = rand(Nx,1); 
D = rand(Nx+1,Nx+1);
D2 = D*D;
identy = eye(Nx+1,Nx+1);
% Computation performed by my code millions of times
for j=1:Nx
  A = a*D2 + b*D + gamma(j)*identy;
  % Can gamma(j)*identy be avoided through vectorization or a faster
  % technique?
end

Is there an alternative way to compute A without using the for loop (would vectorization or repmat also work)?

2 个评论
显示无隐藏无

Matt J 2021-7-28

在 MATLAB Online 中打开

In your example code, A does not evolve throughout the loop. Rather it is repeatedly over-written resulting in

A = a*D2 + b*D + gamma(Nx);

This is probably not what you want, but the reader cannot tell what it should really be doing.

Matthew Kehoe 2021-7-28

编辑：Matthew Kehoe 2021-7-28

I agree @Matt J. I didn't want to ask about all of my code (as it might have appeared to be a difficult question) so I only asked about a small portion of my code. This might have been a mistake. The answer provided by the cyclist has helped me fix/improve my own code significantly.

请先登录，再进行评论。

请先登录，再回答此问题。

Answer 1

the cyclist 2021-7-28

1
链接

此回答的直接链接

https://ww2.mathworks.cn/matlabcentral/answers/887977-efficient-way-to-multiply-a-vector-of-size-nx-1-with-a-matrix-of-size-nx-1-nx-1#answer_755962

在 MATLAB Online 中打开

If you permute gamma to be a vector along the dimension 3, then you can multiply it by identy, and the result will be a 33x33x32 array, where each "slice" along dimension 3 is the multiple with the corresponding value of gamma.

So, this calculation will do all of the calculations of your for loop (and store them all, rather than overwriting as your code does).

A = a*D2 + b*D + permute(gamma,[3,2,1]) .* identy;

You may now have a memory problem, though, if your arrays are large.

5 个评论
显示 3更早的评论隐藏 3更早的评论

Matthew Kehoe 2021-7-28

在 MATLAB Online 中打开

I see, I think that it depends on how you manipulate A inside the for loop. If you have other actions inside the for loop then computing A inside the for loop appears to be faster.

clc; clear all;
Nx = 32;
a = 1.5;
b = 2;
gamma = rand(Nx,1); 
D = rand(Nx+1,Nx+1);
D2 = D*D;
D_start = D(1,:);
D_end = D(end,:);
identy = eye(Nx+1,Nx+1);
n_min = 1.2;
d_min = 1.6;
n_max = 7;
d_max = 1.0;
tic
% Computation performed by my code millions of times
% Run 1000 computations
for ntests=1:1000
  for j=1:Nx
    A = a*D2 + b*D + gamma(j)*identy;
    A(end,:) = n_min*D_end;
    A(end) = A(end) + d_min;
    A(1,:) = n_max*D_start;
    A(1,1) = A(1,1) + d_max;
  end
end
toc
% Elapsed time is 0.042577 seconds.
tic
% Move A outside of the j loop
A2 = a*D2 + b*D + permute(gamma,[3,2,1]) .* identy;
% Computation performed by my code millions of times
% Run 1000 computations
for ntests=1:1000
  for j=1:Nx
    A2(end,:,j) = n_min*D_end;
    A2(end,end,j) = A2(end,end,j) + d_min;
    A2(1,:,j) = n_max*D_start;
    A2(1,1,j) = A2(1,1,j) + d_max;
    A_orig = A2(:,:,j);
  end
end
toc
% Elapsed time is 0.124658 seconds.

Therefore, if you need to change indices then it might be necessary to do everything inside the for loop.

the cyclist 2021-7-28

在 MATLAB Online 中打开

But those other operations don't need to be done in the loop, either. When I ran this on my machine, the fully vectorized code is about 3 times faster.

NTEST = 10000;
Nx = 32;
a = 1.5;
b = 2;
gamma = rand(Nx,1); 
D = rand(Nx+1,Nx+1);
D2 = D*D;
D_start = D(1,:);
D_end = D(end,:);
identy = eye(Nx+1,Nx+1);
n_min = 1.2;
d_min = 1.6;
n_max = 7;
d_max = 1.0;
tic
% Computation performed by my code millions of times
% Run 1000 computations
for ntests=1:NTEST
  for j=1:Nx
    A = a*D2 + b*D + gamma(j)*identy;
    A(end,:) = n_min*D_end;
    A(end) = A(end) + d_min;
    A(1,:) = n_max*D_start;
    A(1,1) = A(1,1) + d_max;
  end
end
toc
Elapsed time is 0.311290 seconds.
tic
% Move A outside of the j loop
A2 = a*D2 + b*D + permute(gamma,[3,2,1]) .* identy;
% Computation performed by my code millions of times
% Run 1000 computations
for ntests=1:NTEST
  for j=1:Nx
    A2(end,:,j) = n_min*D_end;
    A2(end,end,j) = A2(end,end,j) + d_min;
    A2(1,:,j) = n_max*D_start;
    A2(1,1,j) = A2(1,1,j) + d_max;
    A_orig = A2(:,:,j);
  end
end
toc
Elapsed time is 1.313793 seconds.
tic
% Move A outside of the j loop
A3 = a*D2 + b*D + permute(gamma,[3,2,1]) .* identy;
% Computation performed by my code millions of times
% Run 1000 computations
for ntests=1:NTEST
    A3(end,:,:) = repmat(n_min*D_end,[1,1,Nx]);
    A3(end,end,:) = A3(end,end,:) + d_min;
    A3(1,:,:) = repmat(n_max*D_start,[1,1,Nx]);
    A3(1,1,:) = A3(1,1,:) + d_max;
end
toc
Elapsed time is 0.108381 seconds.

Matthew Kehoe 2021-7-28

Yes, I forgot to put A2 inside the ntests loop. Thanks for all of your help @the cyclist. Would it be alright if I asked you one more question about my code (I almost removed the j loop completely but got stuck on one part)? My code has one more place that I can optimize. Would it be alright if I asked you this in a comment (as opposed to making a new question)?

the cyclist 2021-7-29

That's fine, but if the question is self-contained, it might be better to open a new one, because it might get wider exposure.

You can always tag me as you did in your comment, and I will get a notification.

请先登录，再进行评论。

Efficient way to multiply a vector of size (Nx,1) with a matrix of size (Nx+1,Nx+1)

2 个评论
显示无隐藏无

采纳的回答

5 个评论
显示 3更早的评论隐藏 3更早的评论

更多回答（0 个）

另请参阅

类别

标签

产品

版本

Community Treasure Hunt

Efficient way to multiply a vector of size (Nx,1) with a matrix of size (Nx+1,Nx+1)

2 个评论 显示 无隐藏 无

采纳的回答

5 个评论 显示 3更早的评论隐藏 3更早的评论

更多回答（0 个）

另请参阅

类别

标签

产品

版本

Community Treasure Hunt

2 个评论
显示无隐藏无

5 个评论
显示 3更早的评论隐藏 3更早的评论