% Experiment with an iterative map.
% The continuous variable x is a number between -1 and +1 that is changed with each
% iteration. x is analogous to the position of a card in a deck that is shuffled
% once at each iteration.
%
% The behavior of the iterative map is controlled by a dimensionless parameter delta,
% which controls the overlap of the two halves of the deck that will be shuffled
% together.
%
% When delta=0, this map becomes the famous Dyadic Map.
%
% The iterative map is illustrated in the document shuffle.pdf.

Npoints = 500;  % Number of iterations in each trial.
Ndelt = 1000;   % Do this many trials, each with a different value of delta.
x = zeros(Ndelt, Npoints);  % Storage for x values.
delta = 2*(0:Ndelt-1)/(Ndelt-1); % Array of delta values, 0 <= delta <= 2.
t0 = cputime;
for d = 1:Ndelt  % Loop through Ndelt trials.
   x(d,1) = 2*rand - 1;  % Random initial value of x.
   for n=2:Npoints
      x(d,n) = 2*x(d,n-1)*(1 - delta(d)/2) - sign(x(d,n-1)); % Iterative map.
   end
end
seconds_calculating = cputime - t0

% Plot of x values occupied during Npoints iterations for each delta value.
hold off
plot(delta'*ones(1,Npoints), x, '.', 'markersize', 1)
title('Card Shuffling: Bifurcation Diagram')
xlabel('\delta')
ylabel('x_n')
print('bifurcation_diagram.pdf','-dpdfwrite')
seconds_plotting = cputime - seconds_calculating
