I have a matrix multiplication in Matlab that goes as follows

$$\hat{W} = N W N^{T},$$ where $^T$ means a transposition. $N$ is an incidence matrix with the dimensions m x n and W = diag(G), where G is a 1 x n row vector, which makes diag(G) n x n. So basically, when I have a large vector G, e.g. more than 1000 elements I run out of memory because I have a diagonal matrix with n x n elements and most of them are zero. Is there a way to write such multiplication in a more efficient way in Matlab?

N * (G' .* N')

A few releases ago Matlab introduced singleton expansion: in expressions like the one in parentheses the $n\times 1$ matrix G' is "upgraded" to an $n\times m$ matrix (with all columns equal) before the elementwise product.

This new feature gives a cleaner way to implement some tricks that previously required bsxfun.

  • $\begingroup$ This is a nice option as well. Is there a way to represent my other comment. I thank you for your reply. $\endgroup$ – Bakr Jun 5 '20 at 19:27
  • $\begingroup$ @Bakr The one with A + diag(G)? There would be a couple of minor improvements, but in that case the addition is the least of your problems; it is the linear system solution that will take the most CPU time. $\endgroup$ – Federico Poloni Jun 5 '20 at 20:41
  • $\begingroup$ Alright then. I thank you $\endgroup$ – Bakr Jun 5 '20 at 20:49

Solution slightly modified from here

W2 = N * bsxfun(@times, N, G).';

This works for N of size m x n, and G of size 1 x n.

  • $\begingroup$ This is perfect thank you. I searched so much for an answer but I didn't find anything $\endgroup$ – Bakr Jun 5 '20 at 19:00
  • $\begingroup$ Actually I have one more question. I have another formula with $[A - diag(G)]^{-1}[A + diag(G)]$. Is there an alternative for this formula as well? $\endgroup$ – Bakr Jun 5 '20 at 19:23
  • $\begingroup$ Forgot to say that A is a nxn matrix $\endgroup$ – Bakr Jun 5 '20 at 19:42
  • $\begingroup$ @Bakr No, I don't think there's much you can do with that. $\endgroup$ – LedHead Jun 5 '20 at 20:08

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service, privacy policy and cookie policy

Not the answer you're looking for? Browse other questions tagged or ask your own question.