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

1 view (last 30 days)
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 Comments
Matt J
Matt J on 28 Jul 2021
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
Matthew Kehoe on 28 Jul 2021
Edited: Matthew Kehoe on 28 Jul 2021
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.

Sign in to comment.

Accepted Answer

the cyclist
the cyclist on 28 Jul 2021
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 Comments
Matthew Kehoe
Matthew Kehoe on 28 Jul 2021
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
the cyclist on 29 Jul 2021
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.

Sign in to comment.

More Answers (0)

Products


Release

R2020a

Community Treasure Hunt

Find the treasures in MATLAB Central and discover how the community can help you!

Start Hunting!