% Yi-Kai Liu
% 3/23/01; revised 4/28/01; revised 5/9/01

% N * N real symmetric matrices, elements are 
% i.i.d. Gaussian with mean 0 and variance 1

% OK: Mehta normalizes by sqrt(2*N), but I use 2*sqrt(N) ??? 
% Apparently this is due to some notation trouble. In Porter's 
% book, p.197, the semicircle is normalized by 2*sqrt(N). 

% DONE: Exporting data
% FIXME: Examine debugging output

distribution_name = 'Gaussian';

N = input('N * N matrices, choose N: ');
num_matrices = input('Number of matrices: ');

save_to_disk = input('Save output to disk? (1=yes, 0=no) ');
if save_to_disk == 1
   base_filename = input('Filename (without suffix): ', 's');
end;

A = zeros(N);                   % Random matrix
d = zeros(N, num_matrices);     % Eigenvalues (all matrices)
sp = zeros(N-1, num_matrices);  % Spacings (all matrices)


for m = 1 : num_matrices
  
  % Generate random matrix
  % for i = 1 : N
  %   for j = 1 : i-1
  %     A(i,j) = A(j,i);
  %     % Always have j < i 
  %     % (row being read < row being written)
  %   end
  %   for j = i : N
  %     A(i,j) = randn;  % Gaussian distribution
  %   end
  % end
  
  % Generate random matrix--vectorized
  A = triu(randn(N), 0);
  A = A + (triu(A,1))';
  
  % Normalized eigenvalues in interval [-1,1]
  d(:,m) = eig(A)./(2.0*sqrt(N));
  
  sd = sort(d(:,m));   % Sorted eigenvalues

  % Compute spacings
  % for k = 1 : N-1
  %   sp(k,m) = sd(k+1) - sd(k);
  % end
  
  % Compute spacings--vectorized
  sp(:,m) = sd(2:N) - sd(1:N-1);
  
  % Normalized spacings with mean=1
  sp(:,m) = sp(:,m) ./ mean(sp(:,m));
  
end


% Save variables to disk
if save_to_disk == 1
   save([base_filename, '-data'],... 
        'distribution_name', 'save_to_disk', 'base_filename',... 
        'N', 'num_matrices', 'd', 'sp');
end;

% Plot distributions
plotdistrib;