function s = matstats(A)
% MATSTATS : Matrix symmetry statistics for CS 240A problem 3.1
%
% s = matstats(A) returns a struct with the following statistics:
%
% # rows
% # columns
% # nonzeros
% structural rank
% structural symmetry (if square)
% numeric symmetry    (if square)
% # blocks in dm decomposition (if square & full struct rank = sfsr)
% fraction of lower triangle area within diagonal blocks   (if sfsr)
% fraction of matrix nonzeros within diagonal blocks       (if sfsr)
% column diagonal dominance                                (if sfsr)
% row diagonal dominance                                   (if sfsr)
%
% John R. Gilbert, 25 April 2013

[s.nrows,s.ncols] = size(A);
s.nonzeros = nnz(A);
s.srank = sprank(A);
if s.nrows ~= s.ncols
    warning('matrix is not square');
    return;
end;
D = diag(diag(A));
AD = A-D;
if nnz(AD) > 0
    s.ssymm = nnz(AD & AD')/nnz(AD);
else
    s.ssymm = 1;
end;
Asymm = (A+A')/2;
Askew = (A-A')/2;
s.nsymm = (norm(Asymm,'fro')-norm(Askew,'fro'))/norm(A,'fro');
if s.srank ~= min(s.nrows,s.ncols)
    warning('matrix does not have full structural rank');
    return;
end;
[p q r rr] = dmperm(A);
assert(all(r==rr));
B = A(p,q);
s.dmblocks = length(r)-1;
nzb = 0;
for i = 1:s.dmblocks
    I = r(i):r(i+1)-1;
    nzb = nzb + nnz(B(I,I));
end;
s.dblockfractbynnz  = nzb/nnz(A);
s.dblockfractbyarea = sum(diff(r).^2)/prod(size(A));
d = abs(diag(A))';
s.colddom = min( d ./ sum(abs(AD)));
s.rowddom = min( d ./ sum(abs(AD')));
