% Tests for random number generators
% =====
% rand(1,N) returns a row, rand (N,1) returns a column
% note: type "rng" for MATLAB help on details and various mathods.
% Last update January 2022

clear  % Clears command window
clc    % Clears command history
tic    % set timer for the program
%clf    % Removes anything in the figure window before simulation. 

N = 2000;   % size of arrays
k = 5;      % to evaluate correlation between (x(i) and x(i+k))

rng('default')  % initilize RNG (default = "Mersenne Twister")
rng('shuffle')  % seed using current time
%s = rng % (provides information about current rng

% Parking lot test - RNG
x = rand(1,N);
y = rand(1,N);

figure('Name','Parking lot test')
plot(x,y,'or')      % parking lot test
xlabel('xi')
ylabel('yi')
str = sprintf('RNG - parking lot test N = %5i ',N);
title(str);

% Parkin lot test - correlation
cx = zeros(N-k,1);
cy = zeros(N-k,1);
for i=1:N-k
    cx(i) = x(i);
    cy(i) = x(i+k);
end

figure('Name','Correlation test')
plot(cx,cy,'ob')      % parking lot test
xlabel('xi')
ylabel('yi')
str = sprintf('RNG - correlation k = %1i ',k);
title(str);


% set some zeros for sums
Sx = 0.0;
Sy = 0.0;

% test for k-th momentum =========
for i=1:N
    Sx = Sx + x(i)^k;
end
Sx = Sx/N;

fprintf(' %1i',k)
fprintf('-th -momentum = %8.6f \n',Sx)
fprintf(' true momentum  = %8.6f \n',1/(k+1))
fprintf(' \n')
%===================

% test near-neighbor correlation
Sk = 0.0;

for i=1:N-k
    Sk = Sk + x(i)*x(i+k);
end
Sk = Sk/(N-k);
fprintf(' correlation       %8.6f \n',Sk)
fprintf(' true correlation  %8.6f \n',1/4)