Creating a matrix from a dyadic product

Question

• Ax, Ay, Az: [N-by-N]
• B=AA (dyadic product)

means

B(i,j)= [Ax(i,j);Ay(i,j);Az(i,j)]*[Ax(i,j) Ay(i,j) Az(i,j)] 

B(I,j): a 3x3 matrix.

One way to construct B is:

N=2;
Ax=rand(N); Ay=rand(N); Az=rand(N);
t=1;
F=zeros(3,3,N^2);
for i=1:N
 for j=1:N
 F(:,:,t)= [Ax(i,j);Ay(i,j);Az(i,j)]*[Ax(i,j) Ay(i,j) Az(i,j)];
 t=t+1; %# t is just a counter
 end
end
%# then we can write
B = mat2cell(F,3,3,ones(N^2,1));
B = reshape(B,N,N)'; 
B = cell2mat(B);

Is there a faster way than this, especially when N is large?

Answer 1

There is a problem with vector multiplication in the second loop. You should transpose the second vector before doing multiplication.

F(:,:,t)= [Ax(i,j);Ay(i,j);Az(i,j)]*[Ax(i,j) Ay(i,j) Az(i,j)]';

One of the ways to speed things up is to apply vectorization instead of two for loops, such as with Dyadics.