Fast computation of diagonal elements

I have two matrices, A of size m-by-n and B of size m-by-m. I need to quickly compute the digonal elements of (A'*B*A). Here m is of order 100 while n is of order 10000.

How would I do it ? Doing diag(A'*B*A) will be slow.

Neither B nor A is diagonal. they are full matrices.

Answer 1

[Matt J]

 diagonal = sum(A.*(B*A),1);

[Azzi Abdelmalek]

or

sum((A'*B).*A',2)

[Matt J]

It works. But you'll have overhead from the transpositions.

Answer 2

[Sean de Wolski]

m = 100;
n = 10000;
A = rand(m,n);
B = rand(m);

timeit(@()diag(A'*B*A))    % R2013b

ans = 0.8128

So it takes slightly less than a second on my Win7x64 laptop. How many times do you need to do this computation?

[Saurabh]

I need to compute this for possibly 200 iterations.

[Sean de Wolski]

I.e. less than three minutes total? I would just do the above 200x and go get a cup of coffee while MATLAB is crunchin'

Answer 3

[Azzi Abdelmalek]

This will reduce the time

S=A'*B;
d1=arrayfun(@(x) S(x,:)*A(:,x),1:n);

[Saurabh]

This solution is slightly slower than Matt's.