%-----------------------------------------------------------------------
function C = MatMatSax(A,B)
% C = MatMatSax(A,B)
% computes matrix*matrix product C = A*B (via saxpys) where
% A is an mXp matrix, B a pXn matrix.
[m,p] = size(A);
[p,n] = size(B);
C = zeros(m,n);
for j = 1:n
   % Compute j-th column of C
   for k = 1:p
      C(:,j) = C(:,j) + A(:,k)*B(k,j);
   end
end
%        6                  8             6
%  /--$---\        / ----$--- \     /--.---\
%  |--$---|   _8_  | ----$--- |     |--.---|
%  |--$---|   \    | ----$--- |     |--.---|
%  |--$---| =  >   | ----$--- |     |--.---|
%  |--$---|   /__  | ----$--- |  * k|--$---|
%  |--$---|   k=1  | ----$--- |     |--.---|
% 7\--$---/       7\ ----$--- /     |--.---|
%     j=3                k=4       8\--.---/  
%                                      j=3
%     C              C              A             B
% here, a column of C is built by adding the contributions from
% all the columns of A, each multiplied by the corresponding element in B
